I'm trying to implement BLUE estimator in MATLAB for source localization and after my research I've come up with a theoretical example in Steven Kay's "Fundamentals of Statistical Signal Processing: Estimation Theory" book (Example 6.3).

In this solution, while estimating $\theta$, the Gauss-Markov theorem is used and a point called as nominal source position($R_{n_i}$) is defined rather than real source position($R_i$). There is equation which shows the relation between $R_{n_i}$ and $R_i$ found by first order Taylor series expansion.

I'm confused about finding the value of $R_{n_i}$ variable, can anyone explain me this solution?

This could be a very simple question, but I want to understand the theory and my mind is not clear now. When you look at the solution in the book, I think you will understand better than I wrote.

Figure 6.3 from Kay.


N = 8; %total number of receivers
var = 1;
c = 3*10^5; % speed of electromagnetic waves
x = [0 8 7 1 5 9 3 1]; %xi receiver location
y = [0 1 9 9 2 5 8 4]; %yi receiver location
xs = 3; ys = 6; %real source location
xn = 6; yn = 10; %nominal source location
Sxs = xs-xn; Sys = ys-yn;
for i=1:N
    R(i) = sqrt((xs-x(i))^2*(ys-y(i))^2);
    Rn(i) = sqrt((xn-x(i))^2*(yn-y(i))^2);
    co(i) = (xn-x(i))/Rn(i); %cosai
    si(i) = (yn-y(i))/Rn(i); %sinai
    noise(i) = var*randn(1,1); %e1,e2,..
for i=1:N-1
    e(i) = 1/c*(co(i+1)-co(i))*Sxs+1/c*(si(i+1)-si(i))*Sys+noise(i+1)-noise(i); %linear model (6.26)
H = 1/c*[co(2)-co(1) si(2)-si(1);
         co(3)-co(2) si(3)-si(2);
         co(4)-co(3) si(4)-si(3);
         co(5)-co(4) si(5)-si(4);
         co(6)-co(5) si(6)-si(5);
         co(7)-co(6) si(7)-si(6);
         co(8)-co(7) si(8)-si(7)];

A = [-1 1 0 0 0 0 0 0;
     0 -1 1 0 0 0 0 0;
     0 0 -1 1 0 0 0 0;
     0 0 0 -1 1 0 0 0;
     0 0 0 0 -1 1 0 0;
     0 0 0 0 0 -1 1 0;
     0 0 0 0 0 0 -1 1];
inaaT = inv(A*transpose(A));
estTheta = inv(transpose(H)*inaaT*H)*transpose(H)*inaaT*transpose(e); % eq. 6.27

This is not a modular code, I just want to be sure about not missing any point. However, when I run this code I get useless estTheta values. I think again I leave out something in the theory.


1 Answer 1


One problem with your code is that the range calculations should be

R(i) = sqrt((xs-x(i))^2$\color{red}{\bf +}$(ys-y(i))^2);

Rn(i) = sqrt((xn-x(i))^2$\color{red}{\bf +}$(yn-y(i))^2);

But that doesn't appear to help.

The thing that seems to be the issue is that your noise variance is way too high. If I set it to var = 0.00000001 then I get good estimates.

The reason is because you're using this variance (or a standard deviation of 0.0001) to add noise to your H matrix. Looking at the values there:

H matrix values

you can see you're adding something with a variance of 1 (standard deviation of 1), to something of the order of $10^{-5}$. That means your SNR (signal-to-noise ratio) is about -100dB. I know of no algorithm that will help you at SNRs that low.

You left out that Kay defines $$ R_i(x_s,y_s) = \sqrt{(x_s - x_i)^2 + (y_s - y_i)^2} $$ and $$ R_{n_i}(x_n,y_n) = \sqrt{(x_n - x_i)^2 + (y_n - y_i)^2} $$ where $x_s = x_n + \delta x_s$ and $y_s = y_n + \delta y_s$.

Then the aim is to write $R_i$ in terms of $R_{n_i}$.

Wikipedia has an example of Taylor series in two variables, and the first order approximation is just:

$$ f(x,y) \approx f(a,b) + (x-a) \frac{df}{dx} (a,b) + (y-b) \frac{df}{dy} (a,b) + \ldots $$

Then all you need to do is apply Wikipedia's expression to Kay's to get

$$ R_i \approx R_{n_i} + \frac{x_n-x_i}{R_{n_i}} \delta x_s + \frac{y_n - y_i}{R_{n_i}} \delta y_s. $$

  • $\begingroup$ Thanks for the explanation, but how do I obtain xn and yn? I should use the xs=xn+δxs equation, then finding δxs is another blur thing for me. $\endgroup$ Commented Apr 30, 2020 at 20:28
  • $\begingroup$ Kay says: This nominal position, which is close to the true source position, could have been obtained from previous measurements.. The aim is to estimate $\delta x_s$ and $\delta y_s$. $\endgroup$
    – Peter K.
    Commented Apr 30, 2020 at 21:40
  • $\begingroup$ I see, however, in estimation equation (6.27), there is a H matrix which consists of cos and sin values depends on xn, xi and Rni. Here, isn't xn independent of x(n-2), x(n-3) etc. or is xn a constant point? So, what refers to previous measurements exactly? Maybe, if it is possible you can give me the steps of a few iterations for calculating cos and sin values. Some further elaboration could be awesome for my side. $\endgroup$ Commented May 1, 2020 at 1:45
  • $\begingroup$ @sfwmonster : I'm not quite sure what you're asking. What is x(n-2) and x(n-2) ? $\endgroup$
    – Peter K.
    Commented May 1, 2020 at 1:57
  • $\begingroup$ Sorry, let me clear that if the nominal position point is not a constant, I mean if we are sliding it to the real source position, then there will be more than one nominal position points which are at x_(n-1) point, x_(n+1) point etc. It's a better idea to ask only my question that I just want to construct the H matrix depends on cos and sin values. For example if I calculate/estimate cosα_1 = (x_n-x_i)/R_ni, what values should I put instead of x_n and R_ni ? $\endgroup$ Commented May 1, 2020 at 6:44

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.