Planar Approximation: Latitude and Longitude
Date: 04/18/2003 at 11:01:48 From: Jessica Subject: When can two points (lat and long) be treated as planar How can I calculate the distance between two points below which they can be treated as if they were in a plane rather than on a sphere? I am trying to calculate the midpoint between cases of legionella and their nearest neighbor. The cases are usually within 15 miles, but how can I calculate the maximum distance at which I need to treat the two cases as if they were on a sphere (i.e. account for the curvature of the earth)? I've seen many people mention that below 20 km or below 10 miles the earth's curvature doesn't matter - how do they calculate this?
Date: 04/18/2003 at 11:38:18 From: Doctor Rick Subject: Re: When can two points (lat and long) be treated as planar Hi, Jessica. It's a complicated matter. How much error the planar approximation introduces in the distance between points depends not only on the distance but on the latitude and on the direction between the points. For instance, if one point is on the equator, distances north and south, and distances east and west, are both exact regardless of the distance (ignoring effects of nonsphericity). Errors do enter when the direction is at an angle. At the poles, the error is extreme. In fact, if one point is at a pole, the correct formula for distance is (90-lat)*pi*R/180 where lat is the latitude of the second point and R is the radius of the earth. This is very different from the Pythagorean formula! If some distances will be great enough to require the spherical formula, why not just use it all the time? I know it takes more trig calculations, and if one point is fixed, the one trig function in the planar calculation can be pre-calculated, so there may be very good reasons for using the planar calculation whenever possible. But is that the case here? I would use the Haversine formula for distance when distances may be quite small: Deriving the Haversine Formula http://mathforum.org/library/drmath/view/51879.html If you'd like, I could generate a table showing the errors for various distances and directions from a given point. Let me know if this will help, and tell me what latitude to use for the reference point. - Doctor Rick, The Math Forum http://mathforum.org/dr.math/
Date: 04/18/2003 at 12:41:12 From: Jessica Subject: When can two points (lat & long) be treated as planar Thanks for your quick response Dr. Rick! You've convinced me to use the Haversine formula. I am programming this in SAS and need to calculate the distance 1.3 billion times in the program (this simulation is run every day too), so efficiency in calculations will factor in. All of the cases are in NYC (Lat range is 40.5 to ~41 dd and Long is -73.7 to -74.3 dd). I think an error table will be very useful; thank you for the offer. I have been using the law of cosines for spherical trig formula to calculate the distance between two points, but will replace it with the Haversine formula for all distance calculations. After I've found each case's two closest neighbors using the Haversine formula, I need to calculate the midpoint. If the distance is small (i.e. less than 20 km) can I use a simple average of the two lats or longs (see below)? midptlong = ((x1 + x2)/2) midptlat = ((y1 + y2)/2) How would I calculate the midpoint if the distance were greater than 20 km? Many thanks again! Jessica
Date: 04/18/2003 at 20:36:17 From: Doctor Rick Subject: Re: When can two points (lat and long) be treated as planar Hi, Jessica. See: Latitude and Longitude of a Point Halfway between Two Points http://mathforum.org/library/drmath/view/51822.html You just need the formulas for dlon3 and lat3; to find the longitude of the midpoint, use lon3 = lon1 + dlon3 Since you need an efficient algorithm and you are (in the first phase) only interested in *comparing* distances, I suggest that you only calculate 'a' in the Haversine algorithm. It is the square of half the chord length between the two points, so of two choices for point B, the one that has a smaller 'a' value is closer to point A. This saves you two square roots and an atan2 for each point in your search. Of course, you can pre-calculate cos(lat1) for the given point A, leaving only three trig calculations per point in the search. For the error table: I set up a spreadsheet to compute a table of relative errors due to the planar approximation for points at various distances and bearings. The table has a column for each distance (in kilometers) and a row for each bearing (in degrees from north), from a point at a given latitude (41 degrees). For each distance and bearing, I computed the latitude and longitude of the second point using the exact spherical formula: LAT, LON GIVEN DISTANCE (d*R) AND AZIMUTH (tc) lat2=asin(sin(lat1)*cos(d)+cos(lat1)*sin(d)*cos(tc)) IF (cos(lat)=0) lon2=lon1 // endpoint a pole ELSE lon2=mod(lon1-asin(sin(tc)*sin(d)/cos(lat1))+pi,2*pi)-pi ENDIF Then I computed the distance using the approximate formula: PLANAR APPROXIMATION x = (lon2-lon1)*cos(lat1)*pi/180 y = (lat2-lat1)*pi/180 d = R * sqrt(x^2 + y^2) The distances I found were: CALCULATED DISTANCE BY PLANAR APPROXIMATION distance (km) 1 2 5 10 50 0 1 2 5 10 50 15 1.000004416 2.000017668 5.000110532 10.00044283 50.01121206 30 1.000014774 2.000059113 5.000369752 10.00148097 50.03741972 45 1.000024125 2.000096516 5.000603557 10.00241644 50.06085552 60 1.000025585 2.000102347 5.000639801 10.00256008 50.06417495 75 1.000016472 2.000065881 5.000411598 10.00164534 50.0409209 90 0.999999996 1.999999969 4.999999515 9.999996121 49.99951521 105 0.999983524 1.999934086 4.999587879 9.998350475 49.95855614 120 0.999974418 1.999897681 4.999360636 9.997443424 49.9362629 135 0.999975884 1.999903555 4.999397548 9.997592399 49.94024953 150 0.999985233 1.999940949 4.999631225 9.998526849 49.96355763 165 0.999995587 1.999982354 4.999889816 9.99955996 49.98913625 180 1 2 5 10 50 Finally I computed the relative error, that is, (approx. distance - actual distance)/(actual distance) RELATIVE ERROR distance (km) 1 2 5 10 50 0 1.40554E-13 5.3535E-13 1.4051E-13 -1.74083E-14 1.42109E-14 15 4.41569E-06 8.83417E-06 2.21064E-05 4.42826E-05 0.000224241 30 1.47744E-05 2.95567E-05 7.39505E-05 0.000148097 0.000748394 45 2.41246E-05 4.8258E-05 0.000120711 0.000241644 0.00121711 60 2.55851E-05 5.11736E-05 0.00012796 0.000256008 0.001283499 75 1.64723E-05 3.29404E-05 8.23196E-05 0.000164534 0.000818418 90 -3.87847E-09 -1.55143E-08 -9.6964E-08 -3.87854E-07 -9.69588E-06 105 -1.64765E-05 -3.29571E-05 -8.24242E-05 -0.000164952 -0.000828877 120 -2.55816E-05 -5.11596E-05 -0.000127873 -0.000255658 -0.001274742 135 -2.41158E-05 -4.82227E-05 -0.00012049 -0.00024076 -0.001195009 150 -1.47666E-05 -2.95254E-05 -7.3755E-05 -0.000147315 -0.000728847 165 -4.41291E-06 -8.82303E-06 -2.20367E-05 -4.4004E-05 -0.000217275 180 1.40554E-13 1.40554E-13 1.4051E-13 -1.74083E-14 1.42109E-14 You can see that, at latitude 41 degrees, the greatest error at 10 km is 0.0256%, or 2.56 m. Not bad, I'd say. At 50 km the greatest error is 0.128%, or 64 m. That's better than I expected, frankly. At 100 km the error is still under 258 meters (0.26%). - Doctor Rick, The Math Forum http://mathforum.org/dr.math/
Date: 04/23/2003 at 10:13:18 From: Jessica Subject: Thank you (When can two points (lat and long) be treated as planar) Many thanks to you! That was brilliant and very useful! One final question - why is the cos(lat1) term added to the x term? PLANAR APPROXIMATION x = (lon2-lon1)*cos(lat1)*pi/180 y = (lat2-lat1)*pi/180 d = R * sqrt(x^2 + y^2) Thanks again, Jessica
Date: 04/23/2003 at 10:43:08 From: Doctor Rick Subject: Re: When can two points (lat and long) be treated as planar Hi, Jessica. The lines of longitude get closer together as they head toward the poles. Therefore a degree of longitude corresponds to a shorter distance as the latitude increases. To make this quantitative, consider a line of latitude, which is a circle. Looking at a side view of the earth (a cross-section along a line of longitude), we see that the radius of the latitude circle, r, and the radius of the earth, R, form a right triangle with an angle that is the latitude: ********* ****** | ****** *** +----------*** ** | r /** ** | / | ** ** | / | ** * | / | * * | / | * * | / | * * | /R | * * | / | * * | / | * * | / | * * | / | * * | / | * * |/ lat | * *-------------------------*-------------+-----------* r Thus we see that cos(lat) = r/R, so r = R*cos(lat). Now, the difference in longitude between two points at the same latitude is a central angle of the latitude circle. The distance between the points is the radius of the circle times the angle in radians. To convert the longitude difference to radians, we multiply it by pi/180. Then we multiply by r = R*cos(lat) to get the distance x: x = (lon2-lon1)*(pi/180)*R*cos(lat) - Doctor Rick, The Math Forum http://mathforum.org/dr.math/
Search the Dr. Math Library:
Ask Dr. MathTM
© 1994- The Math Forum at NCTM. All rights reserved.