/**
 * Function to calculate distance between two coordinates using the haversine algorithm (less precise).
 * @param coord1    first set of coordinates.
 * @param coord2    second set of coordinates.
 * @returns {number}    distance between the two points in meters.
 */
function distanceInMetersBetweenEarthCoordinates(coord1, coord2) {
    let long1 = coord1[0]
    let lat1 = coord1[1]
    let long2 = coord2[0]
    let lat2 = coord2[1]

    let earthRadiusM = 6371000

    let phi1 = lat1 * Math.PI / 180
    let phi2 = lat2 * Math.PI / 180

    let dlat = (lat2-lat1) * Math.PI / 180
    let dlong = (long2 - long1) * Math.PI / 180

    let a = Math.sin(dlat/2) * Math.sin(dlat/2) +
        Math.cos(phi1) * Math.cos(phi2) *
        Math.sin(dlong/2) * Math.sin(dlong/2)
    let c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1-a))
    return earthRadiusM * c
}


Number.prototype.toRad = function () { return this * Math.PI / 180; }

/**
 * Function to calculate distance between two coordinates using the vincenty algorithm (more precise).
 * @param coord1    first set of coordinates.
 * @param coord2    second set of coordinates.
 * @returns {string|number} distance between the two point is meters.
 */
function distVincenty(coord1, coord2) {
    const lon1 = coord1[0]
    const lat1 = coord1[1]
    const lon2 = coord2[0]
    const lat2 = coord2[1]
    var a = 6378137, b = 6356752.314245,  f = 1/298.257223563;  // WGS-84 ellipsoid params
    var L = (lon2-lon1).toRad()
    var U1 = Math.atan((1-f) * Math.tan(lat1.toRad()));
    var U2 = Math.atan((1-f) * Math.tan(lat2.toRad()));
    var sinU1 = Math.sin(U1), cosU1 = Math.cos(U1);
    var sinU2 = Math.sin(U2), cosU2 = Math.cos(U2);

    var lambda = L, lambdaP, iterLimit = 100;
    do {
        var sinLambda = Math.sin(lambda), cosLambda = Math.cos(lambda);
        var sinSigma = Math.sqrt((cosU2*sinLambda) * (cosU2*sinLambda) +
            (cosU1*sinU2-sinU1*cosU2*cosLambda) * (cosU1*sinU2-sinU1*cosU2*cosLambda));
        if (sinSigma===0) return 0;  // co-incident points
        var cosSigma = sinU1*sinU2 + cosU1*cosU2*cosLambda;
        var sigma = Math.atan2(sinSigma, cosSigma);
        var sinAlpha = cosU1 * cosU2 * sinLambda / sinSigma;
        var cosSqAlpha = 1 - sinAlpha*sinAlpha;
        var cos2SigmaM = cosSigma - 2*sinU1*sinU2/cosSqAlpha;
        if (isNaN(cos2SigmaM)) cos2SigmaM = 0;  // equatorial line: cosSqAlpha=0 (§6)
        var C = f/16*cosSqAlpha*(4+f*(4-3*cosSqAlpha));
        lambdaP = lambda;
        lambda = L + (1-C) * f * sinAlpha *
            (sigma + C*sinSigma*(cos2SigmaM+C*cosSigma*(-1+2*cos2SigmaM*cos2SigmaM)));
    } while (Math.abs(lambda-lambdaP) > 1e-12 && --iterLimit>0);

    if (iterLimit===0) return NaN  // formula failed to converge

    var uSq = cosSqAlpha * (a*a - b*b) / (b*b);
    var A = 1 + uSq/16384*(4096+uSq*(-768+uSq*(320-175*uSq)));
    var B = uSq/1024 * (256+uSq*(-128+uSq*(74-47*uSq)));
    var deltaSigma = B*sinSigma*(cos2SigmaM+B/4*(cosSigma*(-1+2*cos2SigmaM*cos2SigmaM)-
        B/6*cos2SigmaM*(-3+4*sinSigma*sinSigma)*(-3+4*cos2SigmaM*cos2SigmaM)));
    var s = b*A*(sigma-deltaSigma);

    s = s.toFixed(3); // round to 1mm precision
    return s;
}