Sunday, 2 December 2018

coordinate system - How accurate is approximating the Earth as a sphere?


What level of error do I encounter when approximating the earth as a sphere? Specifically, when dealing with the location of points and, for example, the great circle distances between them.


Are there any studies on the average and worst case error compared to an ellipsoid? I'm wondering how much accuracy I'd be sacrificing if I go with a sphere for the sake of easier calculations.


My particular scenario involves directly mapping WGS84 coordinates as if they were coordinates on a perfect sphere (with the mean radius defined by the IUGG) without any transformation.



Answer




In short, the distance can be in error up to roughly 22km or 0.3%, depending on the points in question. That is:




  • The error can be expressed in several natural, useful ways, such as (i) (residual) error, equal to the difference between the two calculated distances (in kilometers), and (ii) relative error, equal to the difference divided by the "correct" (ellipsoidal) value. To produce numbers convenient to work with, I multiply these ratios by 1000 to express the relative error in parts per thousand.




  • The errors depend on the endpoints. Due to the rotational symmetry of the ellipsoid and sphere and their bilateral (north-south and east-west) symmetries, we may place one of the endpoints somewhere along the prime meridian (longitude 0) in the northern hemisphere (latitude between 0 and 90) and the other endpoint in the eastern hemisphere (longitude between 0 and 180).




To explore these dependencies, I have plotted the errors between endpoints at (lat,lon) = (mu,0) and (x,lambda) as a function of latitude x between -90 and 90 degrees. (All points are nominally at an ellipsoid height of zero.) In the figures, rows correspond to values of mu at {0, 22.5, 45, 67.5} degrees and columns to values of lambda at {0, 45, 90, 180} degrees. This gives us a good view of the spectrum of possibilities. As expected, their maximum sizes are approximately the flattening (around 1/300) times the major axis (around 6700 km), or about 22 km.



Errors


Residual errors


Relative errors


Relative errors


Contour plot


Another way to visualize the errors is to fix one endpoint and let the other vary, contouring the errors that arise. Here, for example, is a contour plot where the first endpoint is at 45 degrees north latitude, 0 degrees longitude. As before, error values are in kilometers and positive errors mean the spherical calculation is too large:


Contour plot


It might be easier to read when wrapped around the globe:


Globe plot


The red dot in the south of France shows the location of the first endpoint.



For the record, here is the Mathematica 8 code used for the calculations:


WGS84[x_, y_] := GeoDistance @@ (GeoPosition[Append[#, 0], "WGS84"] & /@ {x, y});
sphere[x_, y_] := GeoDistance @@
(GeoPosition[{GeodesyData["WGS84", {"ReducedLatitude", #[[1]]}], #[[2]], 0}, "WGS84"] & /@ {x, y});

And one of the plotting commands:


With[{mu = 45}, ContourPlot[(sphere[{mu, 0}, {x, y}] - WGS84[{mu, 0}, {x, y}]) / 1000, 
{y, 0, 180}, {x, -90, 90}, ContourLabels -> True]]

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...