Calculate the latitude / longitude distance between two points

Background

I used to work with maps when I was in vocational school, and I thought I would take this opportunity to keep the functions I made at that time.

Commentary

I used the Hubeni formula explained in Mountains. ceremony

\sqrt{(d_yM)^2+(d_zN \cos µ_y)^2}

is.

If you translate it into words in an easy-to-understand manner

\sqrt{(Latitude difference x meridian curve radius)^2 + (Longitude difference x Prime vertical radius of curvature x cos Latitude mean)^2}

I wondered what I was doing. If you look closely

\sqrt{(c−a)^2+(d−b)^2}

It is a formula to find the distance between two points. It seems that the formula corrects the error when assuming that the earth is an ellipse.

code

HubenyDistance.java

public class HubenyDistance {

    //World observation system
    public static final double GRS80_A = 6378137.000;//Semimajor axis a(m)
    public static final double GRS80_E2 = 0.00669438002301188;//First centrifuge e squared

    public static double deg2rad(double deg){
        return deg * Math.PI / 180.0;
    }

    public static double calcDistance(double lat1, double lng1, double lat2, double lng2){
        double my = deg2rad((lat1 + lat2) / 2.0); //Mean of latitude
        double dy = deg2rad(lat1 - lat2); //Latitude difference
        double dx = deg2rad(lng1 - lng2); //Longitude difference

        //Find the radius of curvature of the Prime Vertical line(Radius of the line connecting east and west)
        double sinMy = Math.sin(my);
        double w = Math.sqrt(1.0 - GRS80_E2 * sinMy * sinMy);
        double n = GRS80_A / w;

        //Find the meridian curve radius(Radius of the line connecting north and south)
        double mnum = GRS80_A * (1 - GRS80_E2);
        double m = mnum / (w * w * w);

        //Hubeni's official
        double dym = dy * m;
        double dxncos = dx * n * Math.cos(my);
        return Math.sqrt(dym * dym + dxncos * dxncos);
    }
}