How do I solve for the object distance to each receiver for three radar receivers on the ground, each the same distance from the other, and each receiving echoes, reflected from an object overhead, of a signal pulse from a single transmitter located on the ground at the exact center of the receivers?

Transmitter is located at the exact center of an equlaternal triangle bounded by the the three receivers. The transmitter is t. The three receivers are a, b, c. The overhead object is o.

The physical distance between each pair of receivers is known. (a—b, b—c, c—a ) The time it takes the signal from the transmitter to the object and echoed back to each receiver are known (t—o—a, t—o—b, t—o—c). These three times will be equivalent if the object is directly overhead the transmitter. The three times will be different if the object is at a different distance from each receiver.

The signal from transmitter, reflected from object, and received by receivers, always travels at a constant speed 's'.

Ground is a plane. Any environmental factor (such as air density) that affects speed of transmission pulse is a constant.

=================================================== UPDATE:

I implemented the solution in the answer below into 'C' to test it. The x,y,z coordinate solution for case 1 & 2 looks correct. But the x coordinate is incorrectly 0 for case 3 & 4. It should be < 0 for case 3, and > 0 for case 4.

int main()


{ printf( "\n\n CASE 1: object directly obove TX:"); generate_3D_vector( 10, 10, 10 );

printf( "\n\n CASE 2: object closer to RX-a:");
generate_3D_vector( 9, 10, 10 );

printf( "\n\n CASE 3: object closer to RX-b:");
generate_3D_vector( 10, 9, 10 );

printf( "\n\n CASE 4: object closer to RX-c:");
generate_3D_vector( 10, 10, 9 );


}

void generate_3D_vector( float ra, float rb, float rc )

{ float x, y, z;

x = sqrt( 3.0 ) * ( rb - rc )*( 15 / 24 - 1 / 3 * ra*( ra - rb - rc ) + 2 / 3 * rb*rc ) / ( ra + rb + rc );
y = ( 1.0/2 *ra*( 1 / 8 - ra*( rb + rc ) + rb*rb + rc*rc ) - 1 / 2 * ( ra - rb - rc ) ) / ( ra + rb + rc );
z = sqrt( ra*ra*( ra*ra - 2 - 4 * ( x*x + y*y - y ) ) + ( pow( 2 * y - 1, 2 ))) / ( 2 * ra );

printf( "\n ra rb rc = %.f %.f %.f  ==>  x y z1 = %f %f %f ", ra, rb, rc,    x, y, z );


}

=================================================== UPDATE 4/18/2016:

I implemented the latest x,y,z solution from the formulae, but it didn't work. Next, I copy and pasted, as is, the three AWK equations for x,y,z and got the results below, with problems: x looks ok. Y looks ok except for the 9,10,10 case in which y at -8 seems incorrect. z in each case crashes.

Here is my 'C' source code:

int main()


{

printf( "\n\n CASE 1: object directly obove TX:" );
generate_3D_vector( 10, 10, 10 );

printf( "\n\n CASE 2: object closer to RX-a:" );
generate_3D_vector( 9, 10, 10 );

printf( "\n\n CASE 3: object closer to RX-b:" );
generate_3D_vector( 10, 9, 10 );

printf( "\n\n CASE 4: object closer to RX-c:" );
generate_3D_vector( 10, 10, 9 );


}

void generate_3D_vector( float ra, float rb, float rc ) { float x, y, z; float ra2 = ra*ra; float ra4 = pow( ra, 4 );

float C1 = -sqrt( 3.0 ) / 6.0;
float C2 = -0.5;
float r_a = ra;
float r_b = rb;
float r_c = rc;

x = C1 * ( r_b - r_c ) * ( r_a*( r_a - r_b - r_c ) + 2.0 * r_b*r_c + 1 ) / ( r_a - r_b - r_c );
y = C2 * ( r_a*( r_a*( r_b + r_c ) - r_b*r_b - r_c*r_c + 2.0 ) - r_b - r_c ) / ( r_a - r_b - r_c );
z = sqrt( r_a*r_a * ( r_a*r_a - 4.0 * ( x*x + y*y + y ) - 2 ) + 4 * ( y*y + y ) + 1 ) / ( 2 * r_a );

printf( "\n ra rb rc = %.f %.f %.f  ==>  x y z1 = %f %f %f ", ra, rb, rc, x, y, z );


}

AND HERE IS MY 'C' OUTPUT........................

CASE 1: object directly obove TX: ra rb rc = 10 10 10 ==> x y z1 = 0.000000 0.000000 4.950000

CASE 2: object closer to RX-a: ra rb rc = 9 10 10 ==> x y z1 = 0.000000 -8.272727 -nan(ind) <<<<<<<<<<<<<<<<<<<<< I would expect y to be greater (closer to 0).

CASE 3: object closer to RX-b: ra rb rc = 10 9 10 ==> x y z1 = -2.918826 5.055555 -nan(ind)

CASE 4: object closer to RX-c: ra rb rc = 10 10 9 ==> x y z1 = 2.918826 5.055555 -nan(ind)

• It seems there was a typo in the formula in my answer I didn't notice. I apologize. I'll fix my answer, and include test bench scripts. – Nominal Animal Apr 14 '16 at 17:46

Rewritten on 2016-04-19. The situation is as follows:

(source: nominal-animal.net)

Let us use a coordinate system where the transmitter $$t$$ is at origin, and each receiver is at unit distance from the transmitter, $$d_t = 1$$: \begin{align} d_t &= 1 \\ d_e &= \sqrt{3} \\ \vec{a} = (x_a, y_a, z_a) &= (0, 1, 0) \\ \vec{b} = (x_b, y_b, z_b) &= (-\sqrt{3/4}, -1/2, 0) \approx (-0.8660254, -0.5, 0 ) \\ \vec{c} = (x_c, y_c, z_c) &= (\sqrt{3/4}, -1/2, 0) \approx (0.8660254, -0.5, 0 ) \end{align} If $$\tau$$ is the duration a signal takes from the transmitter to any of the receivers (when there are no reflections; i.e. the signal travels the unit distance), we can derive the signal path lengths $$r$$ from the time it takes for a signal to reflect from the object and reach a receiver: \begin{align} r_a &= \frac{d_a}{\tau} \ge 1 \\ r_b &= \frac{d_b}{\tau} \ge 1 \\ r_c &= \frac{d_c}{\tau} \ge 1 \end{align} They must all be at least one, since there is no shorter path from the transmitter to each receiver. We shall use these signal path lengths to determine the location of the object reflecting the radar signal.

The signal path lengths can be computed using Pythagorean theorem. The object is at $$(x, y, z)$$, and the radar signal must first travel from the transmitter to the object. From the object, the signal must travel to each receiver. Thus, we have a set of three equations: \left \lbrace \begin{align} r_a &= \sqrt{x^2 + y^2 + z^2} + \sqrt{(x - x_a)^2 + (y - y_a)^2 + (z - z_a)^2 } \\ r_b &= \sqrt{x^2 + y^2 + z^2} + \sqrt{(x - x_b)^2 + (y - y_c)^2 + (z - z_b)^2 } \\ r_c &= \sqrt{x^2 + y^2 + z^2} + \sqrt{(x - x_c)^2 + (y - y_b)^2 + (z - z_c)^2 } \end{align} \right. Using Maple to solve the above for $$(x, y, z)$$ we get two solutions, that only differ by the sign of $$z$$. Since our transmitter and receivers are on the same plane, we can only detect objects in one half-space (either nonnegative $$z$$ or nonpositive $$z$$). For simplicity, we select $$z \ge 0$$: \left \lbrace \begin{align} x &= \frac{ -\sqrt{1/12} (r_b - r_c) (r_a (r_a - r_bc) - 2 r_b r_c - 3)}{r_a + r_b + r_c} \\ y &= \frac{ -1/2 (r_a (r_a r_b r_c - r_b^2 - r_c^2 + 2) - r_b r_c)}{r_a + r_b + r_c} \\ z &= \sqrt{ 1/4 r_a^2 - x^2 - y (y - 1) + (y (y - 1) + 1/4) / r_a^2 - 1/2 } \end{align} \right. Note that when the object is very close to the ground, the argument for the square root above may be negative, due to numerical errors. If that occurs, just use $$z = 0$$ instead. (If it is negative but not close to zero, there must be something wrong with the assumed signal path lengths, perhaps an incorrect signal.)

To be exact, the solution Maple finds for $$z$$ is $$z = \frac{\sqrt{terms \; / \; 12}}{r_a + r_b + r_c}$$ where $$terms = - 4 r_a^4 r_b^2 - 4 r_a^4 r_b r_c - 4 r_a^4 r_c^2 + 8 r_a^3 r_b^3 + 4 r_a^3 r_b^2 r_c + 4 r_a^3 r_b r_c^2 + 8 r_a^3 r_c^3 - 4 r_a^2 r_b^4 + 4 r_a^2 r_b^3 r_c - 12 r_a^2 r_b^2 r_c^2 + 4 r_a^2 r_b r_c^3 - 4 r_a^2 r_c^4 - 4 r_a r_b^4 r_c + 4 r_a r_b^3 r_c^2 + 4 r_a r_b^2 r_c^3 - 4 r_a r_b r_c^4 - 4 r_b^4 r_c^2 + 8 r_b^3 r_c^3 - 4 r_b^2 r_c^4 + 3 r_a^4 - 12 r_a^3 r_b - 12 r_a^3 r_c + 30 r_a^2 r_b^2 + 30 r_a^2 r_c^2 - 12 r_a r_b^3 - 12 r_a r_c^3 + 3 r_b^4 - 12 r_b^3 r_c + 30 r_b^2 r_c^2 - 12 r_b r_c^3 + 3 r_c^4 - 30 r_a^2 + 12 r_a r_b + 12 r_a r_c - 30 r_b^2 + 12 r_b r_c - 30 r_c^2 + 27$$ Substituting earlier $$x$$ and $$y$$ into the earlier formula for $$z$$ yields this same expression, so the earlier formula for $$z$$ is mathematically correct, too. (Note that I moved all terms inside the square root in the earlier formula, including the divisor. This makes the case where the argument of the square root is negative much easier to manage, in cases where you cannot fully trust the radar signals. In other words, if the argument is negative and not close to zero, there was something wrong with the radar signal measurements.)

It might be possible to simplify and reorder this latter expression, minimizing numerical errors like cancellation, but I did not find a way to do so in Maple. Testing indicates it is not that important, either: numerical errors in $$z$$ are larger than in $$x$$ or $$y$$, but quite acceptable.

Here are three awk scripts that can be used for testing.

To create random points, random.awk:

#!/usr/bin/awk -f
BEGIN {
if (ARGV[1] == "-h" || ARGV[1] == "--help" || ARGV[1] == "") {
printf "\n" > "/dev/stderr"
printf "Usage: %s [ -h | --help ]\n", ARGV[0] > "/dev/stderr"
printf "       %s POINTS [ RANGE [ SEED ]]\n", ARGV[0] > "/dev/stderr"
printf "\n" > "/dev/stderr"
printf "Default range is 10, and seed is automatically randomized.\n" > "/dev/stderr"
printf "\n" > "/dev/stderr"
exit(0)
}

N = int(ARGV[1])
if (N < 1) {
printf "%s: Invalid number of points.\n", ARGV[1] > "/dev/stderr"
exit(1)
}

if (length("" ARGV[2]) > 0) {
range = 1.0 * ARGV[2]
if (range <= 0) {
printf "%s: Invalid range.\n", ARGV[2] > "/dev/stderr"
exit(1)
}
} else
range = 10.0

if (length("" ARGV[3]) > 0)
srand(ARGV[3]);
else
srand();

for (i = 0; i < N; i++)
printf "%12.6f %12.6f %12.6f\n", range*(2*rand()-1), range*(2*rand()-1), range*rand()

exit(0)
}


#!/usr/bin/awk -f
BEGIN {
x_a =  0.0
y_a =  1.0
z_a =  0.0

x_b = -sqrt(0.75)
y_b = -0.5
z_b =  0.0

x_c =  sqrt(0.75)
y_c = -0.5
z_c =  0.0
}

NF >= 3 {

x = $$1 * 1.0 y =$$2 * 1.0
z = $3 * 1.0 r_a = sqrt(x*x + y*y + z*z) + sqrt( (x - x_a)*(x - x_a) + (y - y_a)*(y - y_a) + (z - z_a)*(z - z_a) ) r_b = sqrt(x*x + y*y + z*z) + sqrt( (x - x_b)*(x - x_b) + (y - y_b)*(y - y_b) + (z - z_b)*(z - z_b) ) r_c = sqrt(x*x + y*y + z*z) + sqrt( (x - x_c)*(x - x_c) + (y - y_c)*(y - y_c) + (z - z_c)*(z - z_c) ) printf "%12.6f %12.6f %12.6f %12.6f %12.6f %12.6f\n", x, y, z, r_a, r_b, r_c }  To estimate the object position based on the radar signal paths, and to calculate the error ranges, radars.awk: #!/usr/bin/awk -f BEGIN { C1 = -sqrt(1.0/12) x_a = 0.0 y_a = 1.0 z_a = 0.0 x_b = -sqrt(0.75) y_b = -0.5 z_b = 0.0 x_c = sqrt(0.75) y_c = -0.5 z_c = 0.0 } NF >= 6 { r_a = $$4 * 1.0 r_b =$$5 * 1.0 r_c =$6 * 1.0

r_abc = r_a + r_b + r_c
r_bc = r_b + r_c

x = C1 * (r_b - r_c) * (r_a * (r_a - r_bc) - 2.0 * r_b * r_c - 3.0) / r_abc;
y = -0.5 * (r_a * (r_a * r_bc - r_b*r_b - r_c*r_c + 2.0) - r_bc) / r_abc;

t1 = r_a * r_a
t2 = y * (y - 1.0)
s = 0.25 * t1 - x*x - t2 + (t2 + 0.25) / t1 - 0.5
if (s > 0.0)
z = sqrt(s)
else
z = 0

err_x = $$1 - x err_y =$$2 - y
err_z = $3 - z if (minerr_x > err_x) minerr_x = err_x if (minerr_y > err_y) minerr_y = err_y if (minerr_z > err_z) minerr_z = err_z if (maxerr_x < err_x) maxerr_x = err_x if (maxerr_y < err_y) maxerr_y = err_y if (maxerr_z < err_z) maxerr_z = err_z err_a = sqrt(x*x + y*y + z*z) + sqrt( (x - x_a)*(x - x_a) + (y - y_a)*(y - y_a) + (z - z_a)*(z - z_a) ) - r_a err_b = sqrt(x*x + y*y + z*z) + sqrt( (x - x_b)*(x - x_b) + (y - y_b)*(y - y_b) + (z - z_b)*(z - z_b) ) - r_b err_c = sqrt(x*x + y*y + z*z) + sqrt( (x - x_c)*(x - x_c) + (y - y_c)*(y - y_c) + (z - z_c)*(z - z_c) ) - r_c if (minerr_a > err_a) minerr_a = err_a if (minerr_b > err_b) minerr_b = err_b if (minerr_c > err_c) minerr_c = err_c if (maxerr_a < err_a) maxerr_a = err_a if (maxerr_b < err_b) maxerr_b = err_b if (maxerr_c < err_c) maxerr_c = err_c printf "%12.6f %12.6f %12.6f %12.6f %12.6f %12.6f %12.6f %12.6f %12.6f %+.6f %+.6f %+.6f %+.6f %+.6f %+.6f\n", $$1,$$2,$3, r_a, r_b, r_c, x, y, z, err_x, err_y, err_z, err_a, err_b, err_c
}

END {
printf "Errors in object coordinates:\n" > "/dev/stderr"
printf "       x: %12.6f .. %12.6f\n", minerr_x, maxerr_x > "/dev/stderr"
printf "       y: %12.6f .. %12.6f\n", minerr_y, maxerr_y > "/dev/stderr"
printf "       z: %12.6f .. %12.6f\n", minerr_z, maxerr_z > "/dev/stderr"
printf "Resulting errors in radar signal path lengths:\n" > "/dev/stderr"
printf " Radar 1: %12.6f .. %12.6f\n", minerr_a, maxerr_a > "/dev/stderr"
printf " Radar 2: %12.6f .. %12.6f\n", minerr_b, maxerr_b > "/dev/stderr"
printf " Radar 3: %12.6f .. %12.6f\n", minerr_c, maxerr_c > "/dev/stderr"
}


These can be chained to run a set of tests. For example, to generate 20 test points within -500..500, -500..500, 0..500, run

./random.awk 20 500 | ./blip.awk | ./radars.awk


which outputs something like

  -42.470214  -460.120727   181.664547     993.936046   992.472463   992.620758     -42.470204  -460.120604   181.664860   -0.000010 -0.000123 -0.000313   +0.000000 +0.000000 +0.000000
-217.450030    -2.192033    96.677738     475.977300   475.170850   476.753487    -217.450114    -2.192003    96.677549   +0.000084 -0.000030 +0.000189   +0.000000 +0.000000 +0.000000
-468.288971  -111.369017   210.325364    1050.802046  1049.711293  1051.255690    -468.288923  -111.369260   210.325343   -0.000048 +0.000243 +0.000021   -0.000000 -0.000000 -0.000000
276.991607   482.728213   344.960062    1308.840645  1310.312838  1309.580551     276.991665   482.728308   344.959883   -0.000058 -0.000095 +0.000179   +0.000000 +0.000000 +0.000000
-202.699862  -150.804663   180.859238     621.902211   620.608527   621.739364    -202.699938  -150.804740   180.859089   +0.000076 +0.000077 +0.000149   +0.000000 +0.000000 +0.000000
-167.340633  -114.491254    28.985180     410.200868   408.653227   410.070261    -167.340621  -114.491243    28.985295   -0.000012 -0.000011 -0.000115   +0.000000 +0.000000 +0.000000
-321.876870   435.167162   392.943801    1337.079259  1337.638721  1338.471828    -321.877012   435.167208   392.943634   +0.000142 -0.000046 +0.000167   +0.000000 +0.000000 -0.000000
29.890078    79.345123   265.040198     556.260754   556.781522   556.595574      29.889988    79.345131   265.040206   +0.000090 -0.000008 -0.000008   -0.000000 -0.000000 -0.000000
-250.865895   230.421345   453.245799    1133.542210  1133.768749  1134.534845    -250.866017   230.421209   453.245800   +0.000122 +0.000136 -0.000001   +0.000000 -0.000000 +0.000000
190.106943  -261.430075   302.619440     886.175382   885.661968   884.917840     190.107006  -261.429917   302.619537   -0.000063 -0.000158 -0.000097   +0.000000 +0.000000 +0.000000
10.765294  -241.908998   255.620615     704.897282   703.893734   703.840727      10.765200  -241.909083   255.620539   +0.000094 +0.000085 +0.000076   -0.000000 +0.000000 +0.000000
152.144771  -485.665604   403.917929    1300.238177  1299.320168  1298.914357     152.144690  -485.665479   403.918110   +0.000081 -0.000125 -0.000181   +0.000000 +0.000000 +0.000000
54.557318  -490.596258   340.339875    1199.972247  1198.824162  1198.666450      54.557322  -490.595937   340.340337   -0.000004 -0.000321 -0.000462   +0.000000 +0.000000 +0.000000
-165.436771   283.870165   372.324676     992.559940   993.129242   993.705963    -165.436869   283.870154   372.324641   +0.000098 +0.000011 +0.000035   +0.000000 +0.000000 +0.000000
402.780470   490.265358   162.968235    1309.444109  1311.098958  1310.034629     402.780615   490.265466   162.967552   -0.000145 -0.000108 +0.000683   +0.000000 +0.000000 -0.000000
246.386135   475.243112   130.961682    1101.341321  1103.022096  1102.248338     246.386223   475.243160   130.961342   -0.000088 -0.000048 +0.000340   +0.000000 -0.000000 -0.000000
467.384508  -385.128184   299.503320    1351.828069  1351.572289  1350.373591     467.384503  -385.128265   299.503224   +0.000005 +0.000081 +0.000096   +0.000000 -0.000000 +0.000000
237.992954  -186.075437   391.296435     989.068560   988.920971   988.086793     237.992979  -186.075157   391.296553   -0.000025 -0.000280 -0.000118   +0.000000 +0.000000 +0.000000
337.241040  -420.753256   295.082717    1230.056531  1229.505282  1228.554482     337.240889  -420.753160   295.083027   +0.000151 -0.000096 -0.000310   +0.000000 +0.000000 +0.000000
117.893074   138.442193    85.819730     401.457263   402.997179   401.983383     117.893083   138.442251    85.819625   -0.000009 -0.000058 +0.000105   +0.000000 -0.000000 +0.000000
Errors in object coordinates:
x:    -0.000145 ..     0.000151
y:    -0.000321 ..     0.000243
z:    -0.000462 ..     0.000683
Resulting errors in radar signal path lengths:


To only look at the errors, testing say a million points, you can run

./random.awk 1e6 500 | ./blip.awk | ./radars.awk >/dev/null


which outputs something like

Errors in object coordinates:
x:    -0.000470 ..     0.000476
y:    -0.000497 ..     0.000509
z:    -0.549315 ..     0.544149
Resulting errors in radar signal path lengths:


As you can see, the $$z$$ coordinate has larger errors compared to $$x$$ and $$y$$. The small errors in radar signal path lengths shows that the calculated coordinates do yield the correct signal path lengths; i.e., that the results should be correct.

In C, the function to compute the object position given radar signal path lengths is for example

typedef struct {
double x;
double y;
double z;
} vec3d;

#define SQRT1OF12 0.2886751345948128822545743902509787278239

vec3d detect(const double r_a, const double r_b, const double r_c)
{
const double r_abc = r_a + r_b + r_c;
const double t1 = r_b + r_c;
vec3d p;

if (r_a < 1.0 || r_b < 1.0 || r_c < 1.0) {
/* Impossible signal path lengths! */
p.x =  0.0;
p.y =  0.0;
p.z = -1.0;
return p;
}

p.x = -SQRT1OF12 * (r_b - r_c) * (r_a * (r_a - t1) - 2.0 * r_b * r_c - 3.0) / r_abc;
p.y = -0.5 * (r_a * (r_a * t1 - r_b*r_b - r_c*r_c + 2.0) - t1) / r_abc;

{
const double t2 = r_a * r_a;
const double t3 = p.y * (p.y - 1.0);
const double s = 0.25 * t2 - p.x * p.x - t3 + (t3 + 0.25) / t2 - 0.5;
if (s > 0.0)
p.z = sqrt(s);
else
p.z = 0.0;
}

return p;
}


I am quite puzzled about the differences to earlier versions to this answer with regards to the actual solution formula. I'm quite sure I used the same Maple recipe, but on other days it yields a quite different solution, one with $$r_a - r_b - r_c$$ in the divisor for $$x$$ and $$y$$, instead of $$r_a + r_b + r_c$$ as here.

These results have been tested for several million random points $$(-1000\dots1000, -1000\dots1000, 0\dots1000)$$, and the errors in both estimated object coordinates $$(-0.001\dots0.001, -0.001\dots0.001, -1\dots1)$$ and in the signal path lengths using the estimated object coordinates ($$0\dots0.001$$) seem acceptable.

• I am pretty sure the question is assuming a 3D case. It mentions the object is overhead, while the transmitter and receivers are on the ground. – spektr Mar 22 '16 at 17:02
• @choward: Hmm, good point. (English is not my native language, and I assumed the flatness referred to the entire problem being a 2D case.) I shall amend my answer accordingly. – Nominal Animal Mar 23 '16 at 8:18
• 3D case, indeed: "overhead" object and "ground based" receivers. – Doug Null Mar 31 '16 at 19:59
• @DougNull: Right. So, the set of $(x,y,z)$ with $\sqrt{terms}$ is then the one you need to apply. If there is something you need clarification on, let me know. – Nominal Animal Mar 31 '16 at 21:55
• Nom.An.: Awesome! How did you derive 'terms'? Where did all those constant numbers come from? (I'm an engineer) – Doug Null Apr 7 '16 at 14:55

For each receiver, the possible locations of the target are on the surface of an ellipsoid of revolution that has the transmitter as one focus and the receiver as the other. The target will be at the intersection of these three surfaces.