Reverse use of Haversine formula

10.1k Views Asked by At

Alright the title is not the best. What I want to do is to change the given parameters in Haversine's formula.

If we know the lat,lng of two points we can calculate their distance. I found the following formulas.

dlon = lon2 - lon1
dlat = lat2 - lat1
a = (sin(dlat/2))^2 + cos(lat1) * cos(lat2) * (sin(dlon/2))^2
c = 2 * atan2(sqrt(a), sqrt(1-a)) 
d = R * c

Let's assume R = 6371 #km

The change I want to make is for only one given point and a given distance to find a point. Now for only those two parameters we create a circle with ofc infinite points around the given point.

Additionally I will create my program to start from the point -90,-180 up to 90,180. So we need to find every time either:

  1. For the same latitude from our circle we have two points in the line. We want to find the one with the bigger longitude. So for the starting point of -70,-170 and a given distance of 100km we can have two points that have the same latitude (approximately): -70,-167.37029 and -70,-172.6297. So I will pick the -70,-167.37029

  2. For the same longitude, I do the same.

Any possible ways to reverse the formula? It has been some time since my last trigonometric problem... Thanks

2

There are 2 best solutions below

0
On

Firstly let us use these variables:

  • $\phi_1, \phi_2$: The latitude of points 1 and 2
  • $\lambda_1, \lambda_2$: The longitude of points 1 and 2
  • $r$: the radius of the sphere (in this case $6371$km)
  • $d$: the distance that we want to find

Now the haversine formula is defined as: $$\text{haversin}\left(\frac dr\right)=\text{haversin}(\phi_2-\phi_1)+\cos(\phi_1)\cos(\phi_2)\text{haversin}(\lambda_2-\lambda_1)$$ and the $\text{haversin}$ function is defined as: $$\text{haversin}(\theta)=\sin^2\left(\frac\theta2\right)=\frac{1-\cos(\theta)}2$$ where $\theta$ is some angle in radians.


Note: your angles are presumably in degrees. This may cause problems in the forumla. To fix this just multiply the degrees angle by the conversion factor of $\frac\pi{180}$. If you want to convert back just multiply by $\frac{180}\pi$.


As you may have noticed $d$ is nested within the function so to extract it we need to use the inverse haversin function given by: $$\text{haversin}^{-1}(\theta)=2\sin^{-1}(\sqrt \theta)$$ Now applying this knowledge we can give a full formula for $d$: $$d=2r\sin^{-1}\left(\sqrt{\text{haversin}(\phi_2-\phi_1)+\cos(\phi_1)\cos(\phi_2)\text{haversin}(\lambda_2-\lambda_1)}\right)$$


On the topic of your "circle" description. You want to have to a point on a sphere and connect it to another point finding $d$. Then you want find all points at distance $d$ from your point. To do this you would need to know $d$ by processing your two points using the haversin formula above. Once $d$ has been found you will be able to rearrange the equation such you can find a point. I am warning you though. This is a pretty hard rearrangement.

0
On

Here is JavaScript that implements the rearrangement of the formula mentioned in the answer.

The "reverse use of the Haversine formula" is in function calculateOuterBounds, which calculates the adjustment needed to latitude or longitude to find the points due north, south, east or west on the circle on the surface of the sphere that has radius centered at given location.

Also included is an implementation of the Haversine formula, distance. Plus some sample code to show how to use the reverse form.

function circle (location, radius) {
    const R = 6371; // Radius of the earth in km

    function distanceTo(point) {
        return distance(location, point)
    }

    function distance(a, b) {
        const lat1 = deg2rad(a.latitude)
        const lat2 = deg2rad(b.latitude)
        const deltaLatitude = deg2rad(a.latitude - b.latitude);
        const deltaLongitude = deg2rad(a.longitude - b.longitude); 
        const x = 
            haversine(deltaLatitude) +
            Math.cos(lat1) * Math.cos(lat2) * haversine(deltaLongitude)
        const c = 2 * Math.atan2(Math.sqrt(x), Math.sqrt(1-x))
        return R * c; // Distance in km
    }

    function calculateOuterBounds() {
        const deltaLatitude = radius / R
        const deltaLongitude = 2 * Math.asin(Math.sqrt(haversine(radius / R) / (Math.cos(deg2rad(location.latitude)) ** 2)))

        return {
            deltaLatitude: rad2deg(deltaLatitude), 
            deltaLongitude: rad2deg(deltaLongitude)
        };
    }

    function haversine(theta) {
        return Math.sin(theta/2) ** 2
    }

    function deg2rad(deg) {
        return deg * (Math.PI/180)
    }

    function rad2deg(rad) {
        return rad * (180/Math.PI)
    }

    return {distanceTo, calculateOuterBounds}
}


/////////////////////////////////////////
// sample code to show how to use

function adjust(point, by) {
    var adjusted = {...point}
    adjusted.latitude += by[0]
    adjusted.longitude += by[1]
    return `${adjusted.latitude}, ${adjusted.longitude}`
}

const theSpire = {latitude: 53.34982169772085, longitude: -6.260252860015207}
const radius = 2

const outerBounds = circle(theSpire, radius).calculateOuterBounds() 

const bounds = {
    north: adjust(theSpire, [outerBounds.deltaLatitude, 0]),
    south: adjust(theSpire, [-outerBounds.deltaLatitude, 0]),
    east: adjust(theSpire, [0, outerBounds.deltaLongitude]),
    west: adjust(theSpire, [0, -outerBounds.deltaLongitude])
}

console.log(`the point ${radius}km north of the spire is ${bounds.north}, and the point the same distance to the west is ${bounds.west}`)