Vincenty Sphere Distance

[GilChrist19]

Hello,

I was testing the distance functions in src/dist.c. The Vincenty sphere distance currently reads:

double distVinSph(double lon1, double lat1, double lon2, double lat2, double r) {

    double x, x1, x2, y;

    lon1 = toRad(lon1);
    lon2 = toRad(lon2);
    lat1 = toRad(lat1);
    lat2 = toRad(lat2);

    x1 = sqrt(cos(lat2) * sin(lon1-lon2));
    x2 = cos(lat1) * sin(lat2) - sin(lat1) * cos(lat2) * cos(lon1-lon2);
    x = x1*x1 + x2*x2;
    y = sin(lat1) * sin(lat2) + cos(lat1) * cos(lat2) * cos(lon1-lon2);

    return r * atan2(x, y);
}

However, this doesn't match the corresponding R code. I propose an updated version read:

double distVinSph(double lon1, double lat1, double lon2, double lat2, double r) {

    double x, x1, x2, y;

    lon1 = toRad(lon1);
    lon2 = toRad(lon2);
    lat1 = toRad(lat1);
    lat2 = toRad(lat2);

        // changed region
        x1 = cos(lat2) * sin(lon1-lon2);
        x2 = cos(lat1) * sin(lat2) - sin(lat1) * cos(lat2) * cos(lon1-lon2);
        x = sqrt(x1*x1 + x2*x2);
        // end change
    y = sin(lat1) * sin(lat2) + cos(lat1) * cos(lat2) * cos(lon1-lon2);

    return r * atan2(x, y);
}

This brings it in line with the R code.

[edzer]

Thanks! This is only a fork of the package, I'm not even one of its authors. Sharing this with @rhijmans, here. Robert, will you take this from here?

[GilChrist19]

Ah, I'm sorry! I thought this was the development version of the CRAN branch.

[rhijmans]

Thank you.