Invalid distance calculation using Haversine formula

Question

DISCLAIMER This is a hometask. Please, don't take it serious or think about this problem for more than 2 seconds. It's just not worth it. I'm sorry.

I am trying to calculate the approximate distance of the route. In the getPoints method I am requesting all the route points from the database and writing them into the array. This part works well and as expected. Here you can see the code snippet:

public ArrayList<Double> getPoints() {
    ArrayList<Double> location = new ArrayList<>();
    SQLiteDatabase db = this.getReadableDatabase();
    Cursor cursor = db.rawQuery("select latitude,longitude from " + Table_Name_Location, null);
    if (cursor.getCount() > 0) {
        while (cursor.moveToNext()) {
            Double latitude = cursor.getDouble(cursor.getColumnIndex("Lat"));
            Double longitude = cursor.getDouble(cursor.getColumnIndex("Longi"));
            location.add(latitude);
            location.add(longitude);
        }
    }
    cursor.close();
    return location;
}

private double distance(double lat1, double lon1, double lat2, double lon2) {
    double theta = lon1 - lon2;
    double dist = Math.sin(deg2rad(lat1))
            * Math.sin(deg2rad(lat2))
            + Math.cos(deg2rad(lat1))
            * Math.cos(deg2rad(lat2))
            * Math.cos(deg2rad(theta));
    dist = Math.acos(dist);
    dist = rad2deg(dist);
    dist = dist * 60 * 1.1515;
    return (dist);
}

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

private double rad2deg(double rad) {
    return (rad * 180.0 / Math.PI);
}

However, when the time comes to the calculating the distance itself the problem occurs. For some reason the Haversine formula calculates the distances incorrectly. I recon, that the problem somewhere in the distance method, but cannot identify it myself. I've compared my implementation of Haversine formula, but haven’t found any differences with the original one.

Answer 1

The Android SDK includes a Location class, which computes the Haversine formula.

Here is the source code.

• Google – java/android/location/Location.java
• Android Code Search – Location.java

The method is titled, "distanceBetween", and utilizes a private "computeDistanceAndBearing" method.
And, the computed values are, distance in meters, initial bearing, and final bearing.

Here is a ported version.
I've removed the BearingDistanceCache argument, and adjusted it to return a double array.

double[] computeDistanceAndBearing(double lat1, double lon1, double lat2, double lon2) {
    // Based on http://www.ngs.noaa.gov/PUBS_LIB/inverse.pdf
    // using the "Inverse Formula" (section 4)

    // Convert lat/long to radians
    lat1 *= Math.PI / 180.0;
    lat2 *= Math.PI / 180.0;
    lon1 *= Math.PI / 180.0;
    lon2 *= Math.PI / 180.0;

    double a = 6378137.0; // WGS84 major axis
    double b = 6356752.3142; // WGS84 semi-major axis
    double f = (a - b) / a;
    double aSqMinusBSqOverBSq = (a * a - b * b) / (b * b);

    double l = lon2 - lon1;
    double aA = 0.0;
    double u1 = Math.atan((1.0 - f) * Math.tan(lat1));
    double u2 = Math.atan((1.0 - f) * Math.tan(lat2));

    double cosU1 = Math.cos(u1);
    double cosU2 = Math.cos(u2);
    double sinU1 = Math.sin(u1);
    double sinU2 = Math.sin(u2);
    double cosU1cosU2 = cosU1 * cosU2;
    double sinU1sinU2 = sinU1 * sinU2;

    double sigma = 0.0;
    double deltaSigma = 0.0;
    double cosSqAlpha;
    double cos2SM;
    double cosSigma;
    double sinSigma;
    double cosLambda = 0.0;
    double sinLambda = 0.0;

    double lambda = l; // initial guess
    for (int iter = 0; iter < 20; iter++) {
        double lambdaOrig = lambda;
        cosLambda = Math.cos(lambda);
        sinLambda = Math.sin(lambda);
        double t1 = cosU2 * sinLambda;
        double t2 = cosU1 * sinU2 - sinU1 * cosU2 * cosLambda;
        double sinSqSigma = t1 * t1 + t2 * t2;
        sinSigma = Math.sqrt(sinSqSigma);
        cosSigma = sinU1sinU2 + cosU1cosU2 * cosLambda;
        sigma = Math.atan2(sinSigma, cosSigma);
        double sinAlpha = (sinSigma == 0) ? 0.0 :
            cosU1cosU2 * sinLambda / sinSigma;
        cosSqAlpha = 1.0 - sinAlpha * sinAlpha;
        cos2SM = (cosSqAlpha == 0) ? 0.0 : cosSigma - 2.0 * sinU1sinU2 / cosSqAlpha;

        double uSquared = cosSqAlpha * aSqMinusBSqOverBSq;
        aA = 1 + (uSquared / 16384.0) * (4096.0 + uSquared * (-768 + uSquared * (320.0
            - 175.0 * uSquared)));
        double bB = (uSquared / 1024.0) * (256.0 + uSquared * (-128.0 + uSquared * (74.0
            - 47.0 * uSquared)));
        double cC = (f / 16.0) * cosSqAlpha * (4.0 + f * (4.0 - 3.0 * cosSqAlpha));
        double cos2SMSq = cos2SM * cos2SM;
        deltaSigma = bB * sinSigma * (cos2SM + (bB / 4.0) * (cosSigma * (-1.0 + 2.0 * cos2SMSq)
            - (bB / 6.0) * cos2SM * (-3.0 + 4.0 * sinSigma * sinSigma) * (-3.0
            + 4.0 * cos2SMSq)));

        lambda = l + (1.0 - cC) * f * sinAlpha * (sigma + cC * sinSigma * (cos2SM
            + cC * cosSigma * (-1.0 + 2.0 * cos2SM * cos2SM)));

        double delta = (lambda - lambdaOrig) / lambda;
        if (Math.abs(delta) < 1.0e-12) {
            break;
        }
    }

    float mDistance, mInitialBearing, mFinalBearing;
    mDistance = (float) (b * aA * (sigma - deltaSigma));
    float initialBearing = (float) Math.atan2(cosU2 * sinLambda,
        cosU1 * sinU2 - sinU1 * cosU2 * cosLambda);
    initialBearing = (float) (initialBearing * (180.0 / Math.PI));
    mInitialBearing = initialBearing;
    float finalBearing = (float) Math.atan2(cosU1 * sinLambda,
        -sinU1 * cosU2 + cosU1 * sinU2 * cosLambda);
    finalBearing = (float) (finalBearing * (180.0 / Math.PI));
    mFinalBearing = finalBearing;

    return new double[] { mDistance, mInitialBearing, mFinalBearing };
}

Here is an example usage.
According to a few sources, the distance from New York City to Los Angeles is approximately, 3,934 km.

double[] nyc = { 40.712778, -74.006111 };
double[] la = { 34.05, -118.25 };
double[] x = computeDistanceAndBearing(nyc[0], nyc[1], la[0], la[1]);

Output

[3945043.75, -86.26732635498047, -114.04286193847656]