Given point of (latitude,longitude), distance and bearing, How to get the new latitude and longitude

Question

I found a piece of code on web. It calculates the Minimum bounding rectangle by a given lat/lon point and a distance.

private static void GetlatLon(double LAT, double LON, double distance, double angle, out double newLon, out double newLat)
        {     
            double dx = distance * 1000 * Math.Sin(angle * Math.PI / 180.0);
            double dy = distance * 1000 * Math.Cos(angle * Math.PI / 180.0);
            double ec = 6356725 + 21412 * (90.0 - LAT) / 90.0;
            double ed = ec * Math.Cos(LAT * Math.PI / 180);
            newLon = (dx / ed + LON * Math.PI / 180.0) * 180.0 / Math.PI;
            newLat = (dy / ec + LAT * Math.PI / 180.0) * 180.0 / Math.PI;
            
            }
            
    public static void GetRectRange(double centorlatitude, double centorLogitude, double distance,
                                  out double maxLatitude, out double minLatitude, out double maxLongitude, out double minLongitude)
{       
            GetlatLon(centorlatitude, centorLogitude, distance, 0, out temp, out maxLatitude);
            GetlatLon(centorlatitude, centorLogitude, distance, 180, out temp, out minLatitude);
            GetlatLon(centorlatitude, centorLogitude, distance, 90, out minLongitude, out temp);
            GetlatLon(centorlatitude, centorLogitude, distance, 270, out maxLongitude, out temp);
}



double ec = 6356725 + 21412 * (90.0 - LAT) / 90.0; //why?
double ed = ec * Math.Cos(LAT * Math.PI / 180);    // why?
dx / ed                                            //why?
dy / ec                                            //why?

6378137 is the equator radius, 6356725 is polar radius, 21412 =6378137 -6356725. from the link, I know a little of the meanings. But these four lines, I don't know why. Could you please help to give more information? Could you please help to let me know the derivation of the formula?

From the link, in the section "Destination point given distance and bearing from start point", it gives another formula to get the result. What is the derivation of the formula?

From this link , I know the derivation of the Haversine Formula, it's very informative. I don't think the formula in the section of "Destination point given distance and bearing from start point" is just a simple reversion of Haversine.

Thanks a lot!

`double ec = 6356725 + 21412 * (90.0 - LAT) / 90.0` looks like a weighted average to estimate the earth's radius between the equator and poles. So if you were on a pole you would just use radius=6356725, and if you were on the equator radius=6378137 (notice it is bigger at the equator). So when you are in between you need to estimate the radius. — Scott, May 03 '15 at 03:41
`double ed = ec * Math.Cos(LAT * Math.PI / 180)` multiplies the estimated earth's radius by the given latitude converted from degrees to radians. So it is just a correction term. Since `dx` and `dy` are calculated distances on a plane `ec` is meant to correct for the curvature of the earth. — Scott, May 03 '15 at 03:55

David C. Rankin · Accepted Answer · 2015-05-03T17:00:53.777

This is a prime example of why commenting your code makes it more readable and maintainable. Mathematically you are looking at the following:

double ec = 6356725 + 21412 * (90.0 - LAT) / 90.0; //why?

This is a measure of eccentricity to account for the equatorial bulge in some fashion. 21412 is, as you know, the difference in earth radius between the equator and pole. 6356725 is the polar radius. (90.0 - LAT) / 90.0 is 1 at the equator, and 0 at the pole. The formula simply estimates how much bulge is present at any given latitude.

double ed = ec * Math.Cos(LAT * Math.PI / 180); // why?

(LAT * Math.PI / 180) is a conversion of latitude from degrees to radians. cos (0) = 1 and cos(1) = 0, so at the equator, you are applying the full amount of the eccentricity while at the pole you are applying none. Similar to the preceding line.

dx / ed //why?

dy / ec //why?

The above seems to be the fractional additions to distance in both the x and y directions attributable to the bulge at any given lat/lon used in the newLon newLat computation to arrive at the new location.

I haven't done any research into the code snippet you found, but mathematically, this is what is taking place. Hopefully that will steer you in the right direction.

Haversine Example in C

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

double m2ft (double l) {            /* convert meters to feet       */
    return l/(1200.0/3937.0);
}

double ft2smi (double l) {          /* convert feet to statute miles*/
    return l/5280.0;
}

double km2smi (double l) {          /* convert km to statute mi.    */
    return ft2smi(m2ft( l * 1000.0 ));
}

static const double deg2rad = 0.017453292519943295769236907684886;
static const double earth_rad_m = 6372797.560856;

typedef struct pointd {
    double lat;
    double lon;
} pointd;

/* Computes the arc, in radian, between two WGS-84 positions.
   The result is equal to Distance(from,to)/earth_rad_m
      = 2*asin(sqrt(h(d/earth_rad_m )))
   where:
      d is the distance in meters between 'from' and 'to' positions.
      h is the haversine function: h(x)=sin²(x/2)
   The haversine formula gives:
      h(d/R) = h(from.lat-to.lat)+h(from.lon-to.lon)+cos(from.lat)*cos(to.lat)
   http://en.wikipedia.org/wiki/Law_of_haversines
 */
double arcradians (const pointd *from, const pointd *to)
{
    double latitudeArc  = (from-> lat - to-> lat) * deg2rad;
    double longitudeArc = (from-> lon - to-> lon) * deg2rad;

    double latitudeH = sin (latitudeArc * 0.5);
    latitudeH *= latitudeH;

    double lontitudeH = sin (longitudeArc * 0.5);
    lontitudeH *= lontitudeH;

    double tmp = cos (from-> lat * deg2rad) * cos (to-> lat * deg2rad);

    return 2.0 * asin (sqrt (latitudeH + tmp*lontitudeH));
}

/* Computes the distance, in meters, between two WGS-84 positions.
   The result is equal to earth_rad_m*ArcInRadians(from,to)
 */
double dist_m (const pointd *from, const pointd *to) {
    return earth_rad_m * arcradians (from, to);
}

int main (int argc, char **argv) {

    if (argc < 5 ) {
        fprintf (stderr, "Error: insufficient input, usage: %s (lat,lon) (lat,lon)\n", argv[0]);
        return 1;
    }

    pointd points[2];

    points[0].lat = strtod (argv[1], NULL);
    points[0].lon = strtod (argv[2], NULL);

    points[1].lat = strtod (argv[3], NULL);
    points[1].lon = strtod (argv[4], NULL);

    printf ("\nThe distance in meters from 1 to 2 (smi): %lf\n\n", km2smi (dist_m (&points[0], &points[1])/1000.0) );

    return 0;
}

/* Results/Example.
    ./bin/gce 31.77 -94.61 31.44 -94.698
    The distance in miles from Nacogdoches to Lufkin, Texas (smi): 23.387997 miles
 */

Thank you David. From the link(http://www.movable-type.co.uk/scripts/latlong.html), in the section "Destination point given distance and bearing from start point", it gives another formula to get the result. What is the derivation of the formula? Are there any web page like the link(http://mathforum.org/library/drmath/view/51879.html), which demonstrate the Haversine Formula very clearly. — gfan, May 03 '15 at 14:29
Haversine formula isn't bad. I think I have a port to C. Basically in all lat/lon problems you are dealing with **great circle distances** (which is just a fancy term for a spherical system of coordinates over the earth's surface.) It all stems from the fact that a line between point `a` and point `b` on a globe isn't a **line** at all, it is an **arc**. Therefore the distance you are finding is the **arc** distance from `a` to `b`. `archlength = radius * angle`. That's why you see the trig, and `r` is changing due to the bulge, so that is why you see the `ec`, `ed` tweaks. — David C. Rankin, May 03 '15 at 16:25
Also **"Destination point given distance and bearing from start point"** equation is the same formula, but instead of having `point1` & `point2`, you have `point1` and the `angle from north & arc length` as your initial conditions. Here you already have the distance, so instead of finding the distance, you use the distance to find `point2`. (been many years since I derived spherical coordinate equations -- and it will be many more until I do it again -- good luck `:p` ) — David C. Rankin, May 03 '15 at 16:54

score 0 · Answer 2 · edited May 23 '17 at 11:58

0

I assume 6356725 has something to do with the radius of the earth. Check out this answer, and also take a look at the Haversine Formula.

edited May 23 '17 at 11:58

Community

1
1

answered May 03 '15 at 03:03

tmthydvnprt

10,398
8
52
72

Given point of (latitude,longitude), distance and bearing, How to get the new latitude and longitude

2 Answers2

Haversine Example in C

Linked