Re: Triangulation Problem

*Subject*: [mg3083] Re: [mg2968] Triangulation Problem*From*: danl (Daniel Lichtblau)*Date*: 31 Jan 1996 04:39:42 -0600*Approved*: usenet@wri.com*Distribution*: local*Newsgroups*: wri.mathgroup*Organization*: Wolfram Research, Inc.*Sender*: daemon at wri.com

In article <4ekebp$5o1 at dragonfly.wri.com> Preston Nichols <nichols at godel.math.cmu.edu> writes: > David Hoare wrote: > <<<< > I am looking for the trigonometric / algebreic formulae one would > use to triangulate an unknown position, knowing 3 fixed points and > the respective distances to the unknown point. (make sense?) Oh - in > 3 dimensional space. > > [....] > > It's for a program I'm writing, and it's been a while since I've done > this sort of thing - - my math has left me! > Any Help would be GREATLY appreciated!! > >>>> > > I have seen essentially two responses to this, and I have a few > comments on each. > > 1. Anthony D. Gleckler's approach [mg3007], slightly paraphrased > (and the approach I used): > > If a and b are points (given as lists of coordinates), then the > distance between them is, in Mathematica, > > Sqrt[(a-b).(a-b)] > > Let the three given points be > a = {a1,a2,a3}; > b = {b1,b2,b3}; > c = {c1,c2,c3}; > > let the associated distances be A,B,C, and let the unknown point be > > X = {x,y,z}; > > Then the desired point is the simultaneous solution of these three > equations: > > spherea := ( (X-a).(X-a) == A^2 ); > sphereb := ( (X-b).(X-b) == B^2 ); > spherec := ( (X-c).(X-c) == C^2 ); > > each of which just states the desired equality of distances, but > with both sides squared. This is three quadratic equations in three > unknowns; > > sols = Solve[{eqa,eqb,eqc}, X] > > will return the solution as replacement rules, and > > soughtpoints = X/.sols > > will give the unknown point as a list of coordinates -- but not > really! First of all, the Mathematica commands given above work > fine for numerical values for a, b, c, A, B, and C. If you try to > get a general symbolic solution, however, the expressions get > unmanageably big. > > But there's another issue where writing the necessary code is not > hard, but deciding just what you want the code to do requires some > reflection. Unless the measurements of the given data (a, b, c, A, > B, and C) are absolutely perfect, the system of equations will have > TWO DISTINCT SOLUTIONS, possibly complex (having equal and opposite > imaginary parts), and the soughtpoints as computed above will > actually be a list of TWO DISTINCT POINTS, possibly with complex > coordinates. Each of the equations shperea, sphereb, and spherec > define a sphere; two of the spheres generically intersect in a > circle, which will generically intersect the third sphere in two > points. These two points (if real) determine a line perpendicular > to the plane containing the given points a,b,c, and the two points > are equidistant from that plane. > > If the measurements of a, b, c, A, B, and C are only slightly off, > these two points will be very close together, and simply averaging > them will give a reasonable, usable result: > > thepoint = Apply[Plus, soughtpoints]/2 > > (Conveniently, this will give a real result even if complex > solutions are involved!) > > You might, however, wish to use the distance between the two > solutions as an indicator of the quality of your data (i.e. the > given points and distances). At the simplest level, a relatively > "large" distance between the two solution points indicates that your > values for A, B, and C *can't* be accepted as distances of a, b, > and c from a common fourth point; you could use that "largeness" > simply as a tipoff that something's wrong with the data. With some > more effort, you could also convert the size of that distance into a > numerical measurement of how well all your pieces are fitting > together. > > --------------------------------------------------- > 2. Dave Rusin's approach [mg3015] > > I really enjoyed Dave's method, because it so closely follows what > we might do if we could directly manipulate abstract geometrical > objects in 3-space with our hands. He suggests locating and > explicitly parametrizing the circle where the first two spheres > intersect. (I also think parametrizing curves is underappreciated > in general.) This can be done in Mathematica by using some tricks > with NullSpace and LinearAlgebra`Orthogonalization`GramSchmidt. > > Dave's method will lead to the same two solutions for less than > perfect data, but this time they appear in the (ususally) multiple > solutions of the trigonometric equation > > > D3^2 = (r cos(t)-a)^2 + (r sin(t)-b)^2 + (x-c)^2 > > the solving of which involves taking multi-valued arcsines. (In 1, > we might have been alerted to multiple solutions by the fact that > we were solving quadratic equations.) > > --------------------------------------------------- > > Preston Nichols > Visiting Assistant Professor > Department of Mathematics > Carnegie Mellon University > (seeking employment for Fall 1996) > > If the system has any real solution, then it will have two such. Averaging them will not yield a solution; it will give a fourth point in the same plane as the first three. But there is no reason for the solution to lie in that plane, so far as the original problem stated. I am guessing that the original user had the following type of scenario in mind: fix three points on the surface of the earth, sufficiently close that curvature of the surface does not matter (or else, account for it at the end of the problem), and then use the measured distances to locate, say, the tip of Mt Everest, with respect to the coordinate system. Of course the measured distances give two possible solutions, but one of them will be several miles BELOW the surface, hence is unimportant. While I also found the cited methods enlightening, I'd say this problem is best handled by FindRoot (which may well also have been proposed; I do not recall). The input is inexact and with a bit of thought one can guess at an approximation of the solution. Exactly what one wants for FindRoot. This will avoid the extraneous solution and will also be faster for this sort of problem. Daniel Lichtblau Wolfram Research, Inc. danl at wri.com