Get Lat/Long Given Current Point, Distance and Bearing

Get lat/long given current point, distance and bearing

Needed to convert answers from radians back to degrees. Working code below:

import math

R = 6378.1 #Radius of the Earth
brng = 1.57 #Bearing is 90 degrees converted to radians.
d = 15 #Distance in km

#lat2  52.20444 - the lat result I'm hoping for
#lon2  0.36056 - the long result I'm hoping for.

lat1 = math.radians(52.20472) #Current lat point converted to radians
lon1 = math.radians(0.14056) #Current long point converted to radians

lat2 = math.asin( math.sin(lat1)*math.cos(d/R) +
     math.cos(lat1)*math.sin(d/R)*math.cos(brng))

lon2 = lon1 + math.atan2(math.sin(brng)*math.sin(d/R)*math.cos(lat1),
             math.cos(d/R)-math.sin(lat1)*math.sin(lat2))

lat2 = math.degrees(lat2)
lon2 = math.degrees(lon2)

print(lat2)
print(lon2)

Find destination coordinates given starting coordinates, bearing, and distance

According to this page

math.sin(x) ... Return the sine of x radians.

So, you need to convert lat1 and lon1 to radians before your equation and then later you can convert lat2 and lon2 back into degrees.

Get accurate lat/long given current point, distance and bearing

You can use the spherical geometry utilities of the Google Maps SDK for Android utility library.

The specific method that you can use will be the computeOffset method. You can see it when you search the SphericalUtil in this link.

ComputeOffset returns the LatLng resulting from moving a distance from an origin in the specified heading (expressed in degrees clockwise from north). Here are the parameters that you need to specify:

from - The LatLng from which to start.
distance - The distance to travel.
heading - The heading in degrees clockwise from north.

To use this, add the utility library to your dependencies:

implementation 'com.google.maps.android:android-maps-utils:0.5'

This is the sample of a line of code that I use:

LatLng newCoordinates = SphericalUtil.computeOffset(new LatLng(-25.363, 131.044),30000, 90);

Hope this helps!

calculating Lat and Long from Bearing and Distance

Your problem is on your first line.

Try

double dist = 150.0 / 6371.0;

The reason is that 150/6371 gets calculated as 0, because it performs integer division (rather than floating point division). This is true even though the result is being stored in a double. You can force floating point division by making one of the two numbers a floating point literal.

Attempting to calculate lat and long given distance and bearing from start point

Here is my code that solved the problem:

public static (double Lat, double Lon) Destination((double Lat, double Lon) startPoint, double distance, double bearing)
    {
        double φ1 = startPoint.Lat * (Math.PI / 180);
        double λ1 = startPoint.Lon * (Math.PI / 180);
        double brng = bearing * (Math.PI / 180);
        double φ2 = Math.Asin(Math.Sin(φ1) * Math.Cos(distance / radius) + Math.Cos(φ1) * Math.Sin(distance / radius) * Math.Cos(brng));
        double λ2 = λ1 + Math.Atan2(Math.Sin(brng) * Math.Sin(distance / radius) * Math.Cos(φ1), Math.Cos(distance / radius) - Math.Sin(φ1) * Math.Sin(φ2));
        return (φ2 * (180 / Math.PI), λ2 * (180 / Math.PI));
    }

radius is the radius of the Earth in meters.

Calculating coordinates given a bearing and a distance

It seems like these are the issues in your code:

1. You need to convert lat1 and lon1 to radians before calling your function.
2. You may be scaling radialDistance incorrectly.
3. Testing a floating-point number for equality is dangerous. Two numbers that are equal after exact arithmetic might not be exactly equal after floating-point arithmetic. Thus abs(x-y) < threshold is safer than x == y for testing two floating-point numbers x and y for equality.
4. I think you want to convert lat and lon from radians to degrees.

Here is my implementation of your code in Python:

#!/usr/bin/env python

from math import asin,cos,pi,sin

rEarth = 6371.01 # Earth's average radius in km
epsilon = 0.000001 # threshold for floating-point equality

def deg2rad(angle):
    return angle*pi/180

def rad2deg(angle):
    return angle*180/pi

def pointRadialDistance(lat1, lon1, bearing, distance):
    """
    Return final coordinates (lat2,lon2) [in degrees] given initial coordinates
    (lat1,lon1) [in degrees] and a bearing [in degrees] and distance [in km]
    """
    rlat1 = deg2rad(lat1)
    rlon1 = deg2rad(lon1)
    rbearing = deg2rad(bearing)
    rdistance = distance / rEarth # normalize linear distance to radian angle

    rlat = asin( sin(rlat1) * cos(rdistance) + cos(rlat1) * sin(rdistance) * cos(rbearing) )

    if cos(rlat) == 0 or abs(cos(rlat)) < epsilon: # Endpoint a pole
        rlon=rlon1
    else:
        rlon = ( (rlon1 - asin( sin(rbearing)* sin(rdistance) / cos(rlat) ) + pi ) % (2*pi) ) - pi

    lat = rad2deg(rlat)
    lon = rad2deg(rlon)
    return (lat, lon)

def main():
    print "lat1 \t lon1 \t\t bear \t dist \t\t lat2 \t\t lon2"
    testcases = []
    testcases.append((0,0,0,1))
    testcases.append((0,0,90,1))
    testcases.append((0,0,0,100))
    testcases.append((0,0,90,100))
    testcases.append((49.25705,-123.140259,225,1))
    testcases.append((49.25705,-123.140259,225,100))
    testcases.append((49.25705,-123.140259,225,1000))
    for lat1, lon1, bear, dist in testcases:
        (lat,lon) = pointRadialDistance(lat1,lon1,bear,dist)
        print "%6.2f \t %6.2f \t %4.1f \t %6.1f \t %6.2f \t %6.2f" % (lat1,lon1,bear,dist,lat,lon)

if __name__ == "__main__":
    main()

Here is the output:

lat1     lon1        bear    dist        lat2        lon2
  0.00     0.00       0.0       1.0        0.01        0.00
  0.00     0.00      90.0       1.0        0.00       -0.01
  0.00     0.00       0.0     100.0        0.90        0.00
  0.00     0.00      90.0     100.0        0.00       -0.90
 49.26   -123.14     225.0      1.0       49.25      -123.13
 49.26   -123.14     225.0    100.0       48.62      -122.18
 49.26   -123.14     225.0   1000.0       42.55      -114.51

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

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
 */

Calculate second point knowing the starting point and distance

It seems you are measuring distance (R) in meters, and bearing (theta) counterclockwise from due east. And for your purposes (hundereds of meters), plane geometry should be accurate enough. In that case,

dx = R*cos(theta) ; theta measured counterclockwise from due east
dy = R*sin(theta) ; dx, dy same units as R

If theta is measured clockwise from due north (for example, compass bearings),
the calculation for dx and dy is slightly different:

dx = R*sin(theta)  ; theta measured clockwise from due north
dy = R*cos(theta)  ; dx, dy same units as R

In either case, the change in degrees longitude and latitude is:

delta_longitude = dx/(111320*cos(latitude))  ; dx, dy in meters
delta_latitude = dy/110540                   ; result in degrees long/lat

The difference between the constants 110540 and 111320 is due to the earth's oblateness
(polar and equatorial circumferences are different).

Here's a worked example, using the parameters from a later question of yours:

Given a start location at longitude -87.62788 degrees, latitude 41.88592 degrees,
find the coordinates of the point 500 meters northwest from the start location.

If we're measuring angles counterclockwise from due east, "northwest" corresponds
to theta=135 degrees. R is 500 meters.

dx = R*cos(theta) 
   = 500 * cos(135 deg) 
   = -353.55 meters

dy = R*sin(theta) 
   = 500 * sin(135 deg) 
   = +353.55 meters

delta_longitude = dx/(111320*cos(latitude)) 
                = -353.55/(111320*cos(41.88592 deg))
                = -.004266 deg (approx -15.36 arcsec)

delta_latitude = dy/110540
               = 353.55/110540
               =  .003198 deg (approx 11.51 arcsec)

Final longitude = start_longitude + delta_longitude
                = -87.62788 - .004266
                = -87.632146

Final latitude = start_latitude + delta_latitude
               = 41.88592 + .003198
               = 41.889118

---

Related Topics