Friday 20 March 2015

coordinate system - Long range Trilateration from 3 lat, long, range


I've been looking at the following solution to a similar problem Trilateration using 3 latitude and longitude points, and 3 distances However with my problem the 3 points are spread over a large distance and the distance is a distance-over-ground (ignoring altitude). So the curvature of the earth comes into play.



Example of one of the data points is...[lat,long,km]


51.505348978413,-0.11270052661132,75.639168 51.4845086249581,-3.16947466601561,173.809152 55.9568033803799,-3.20266968478392,465.100416


This can be visualised at http://www.freemaptools.com/radius-around-point.htm by pasting it into the section "CSV Upload - [latitude,longitude,radius(km)] per line"


Using the solution in Trilateration using 3 latitude and longitude points, and 3 distances I get the result..


51.8247121391 -0.802763718095


which isn't bad but it's about 14 miles out. I'm really wanting this down to ~5 miles of error.


I changed the ECEF calculation to one from here http://www.mathworks.co.uk/matlabcentral/fileexchange/7942-covert-lat-lon-alt-to-ecef-cartesian/content/lla2ecef.m and used an altitude of 0.


However, this is way way off the mark with...


56.2518456429 -0.813853131422


Intermediate results for xA,yA,zA...etc. are



lla2ecef


xA = [3978.79751502] yA = [-7.82628596545] zA = [5342.89431232]


xB = [3974.53105834] yB = [-220.086727248] zB = [5342.88794881]


xC = [3573.7981748] yC = [-199.973379186] zC = [5344.22411569]


original method


xA = [3965.56758195] yA = [-7.8002627162] zA = [4986.36680448]


xB = [3961.32002974] yB = [-219.355176278] zB = [4984.92406448]


xC = [3561.02859617] yC = [-199.258852046] zC = [5279.1109334]


So... the Z values for lla2ecef method are quite different. I'm guessing this is what's causing the issue, something to do with fixing the altitude to 0 and the fact my distances are distance-over-ground.


Am I going to need to do some projections instead?



My brain is about melted for today, can anyone shed light on ways to tackle this problem?



Answer



Using the Non Linear Fitting method as provided by whuber in (Trilateration algorithm for n amount of points) I used Matlab to solve the problem.


P = [51.505348978413,-0.11270052661132;
51.4845086249581,-3.16947466601561;
55.9568033803799,-3.20266968478392];

R = [75.639168, 173.809152, 465.100416];
R = R*1000; %convert km to m


opts = statset('MaxIter',1000);
beta0 = [53.374, -1.394];
BETA=nlinfit(P,R,'myfun',beta0,opts)


  • P are the 3 trilateration points (lat, lon)

  • R are the ranges from my sample point to points in P distance-over-ground

  • beta0 are the initial guesses for a position that meets the condition defined by P and R


myfun is in a separate file



function result = myfun(beta, X)

ellipsoid = [6378137, 8.1819190842622e-2];
result(1) = distance(51.505348978413, -0.11270052661132, beta(1), beta(2), ellipsoid);
result(2) = distance(51.4845086249581, -3.16947466601561, beta(1), beta(2), ellipsoid);
result(3) = distance(55.9568033803799, -3.20266968478392, beta(1), beta(2), ellipsoid);

I'm not sure what X does in this case, but myfun needs to have two inputs (beta & X). The initial guess returns three distances to my points P. result is compared to R until a maximum of 1000 iterations (or some error is reduced below a limit - not sure what this limit is or where to set it yet).


Once the values for beta have been found that results in the same values of R the final values are stored in BETA


No comments:

Post a Comment

arcpy - Changing output name when exporting data driven pages to JPG?

Is there a way to save the output JPG, changing the output file name to the page name, instead of page number? I mean changing the script fo...