Skip to main content
Announcements
Qlik Connect 2024! Seize endless possibilities! LEARN MORE
cancel
Showing results for 
Search instead for 
Did you mean: 
Not applicable

WGS84 Compliant Great Circle Distance Calculation

Hi,

I am trying to create a system that shows the distance between two coordinates. So far I've managed to get part way there using the Haversine formula to calculate great circle distance based on a spherical earth model. However, the requirements have just been tightened up and I now need to try and find a way of using a WGS84 compliant method (oblate spheroid model) which is much more complex.

I've found a suitable formula (the Vincenty algorithm http://www.movable-type.co.uk/scripts/latlong-vincenty.html ) but have yet to figure out how to get this into QlikView! The database the information is coming out of is clunky, in-house written (and those that wrote it no longer work for the company) and not Oracle based so ideally I need to figure this out without involving the source database. I have a a list of journeys in a qvd file with to and from coordinates in them and I need to calculate the distance ... any ideas? 🙂

Ever hopeful!

Emma

1 Solution

Accepted Solutions
Not applicable
Author

hi,

Take a look at the source of the linked HTML page == there is Javascript with the formula 🙂 !

Makea QlikView file with a module written in Jscript. Paste the code there.


/**
* Calculates geodetic distance between two points specified by latitude/longitude using
* Vincenty inverse formula for ellipsoids
*
* @param {Number} lat1, lon1: first point in decimal degrees
* @param {Number} lat2, lon2: second point in decimal degrees
* @returns (Number} distance in metres between points
*/
function distVincenty(lat1, lon1, lat2, lon2) {
var a = 6378137, b = 6356752.3142, 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;

// note: to return initial/final bearings in addition to distance, use something like:
var fwdAz = Math.atan2(cosU2*sinLambda, cosU1*sinU2-sinU1*cosU2*cosLambda);
var revAz = Math.atan2(cosU1*sinLambda, -sinU1*cosU2+cosU1*sinU2*cosLambda);
return { distance: s, initialBearing: fwdAz.toDeg(), finalBearing: revAz.toDeg() };
}

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

/** extend Number object with methods for converting degrees/radians */
toRad = function() { // convert degrees to radians
return this * Math.PI / 180;
}
toDeg = function() { // convert radians to degrees
return this * 180 / Math.PI;
}

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
/* Latitude/longitude spherical geodesy formulae & scripts (c) Chris Veness 2002-2010 */
/* - www.movable-type.co.uk/scripts/latlong.html */
/* */
/* Sample usage: */
/* var p1 = new LatLon(51.5136, -0.0983); */
/* var p2 = new LatLon(51.4778, -0.0015); */
/* var dist = p1.distanceTo(p2); // in km */
/* var brng = p1.bearingTo(p2); // in degrees clockwise from north */
/* ... etc */
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
/* Note that minimal error checking is performed in this example code! */
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */


/**
* Creates a point on the earth's surface at the supplied latitude / longitude
*
* @constructor
* @param {Number} lat: latitude in numeric degrees
* @param {Number} lon: longitude in numeric degrees
* @param {Number} [rad=6371]: radius of earth if different value is required from standard 6,371km
*/
function LatLon(lat, lon, rad) {
if (typeof(rad) == 'undefined') rad = 6371; // earth's mean radius in km
// only accept numbers or valid numeric strings
this._lat = typeof(lat)=='number' ? lat : typeof(lat)=='string' && lat.trim()!='' ? +lat : NaN;
this._lon = typeof(lat)=='number' ? lon : typeof(lon)=='string' && lon.trim()!='' ? +lon : NaN;
this._radius = typeof(rad)=='number' ? rad : typeof(rad)=='string' && trim(lon)!='' ? +rad : NaN;
}


/**
* Returns the distance from this point to the supplied point, in km
* (using Haversine formula)
*
* from: Haversine formula - R. W. Sinnott, "Virtues of the Haversine",
* Sky and Telescope, vol 68, no 2, 1984
*
* @param {LatLon} point: Latitude/longitude of destination point
* @param {Number} [precision=4]: no of significant digits to use for returned value
* @returns {Number} Distance in km between this point and destination point
*/
distanceTo = function(point, precision) {
// default 4 sig figs reflects typical 0.3% accuracy of spherical model
if (typeof precision == 'undefined') precision = 4;

var R = this._radius;
var lat1 = this._lat.toRad(), lon1 = this._lon.toRad();
var lat2 = point._lat.toRad(), lon2 = point._lon.toRad();
var dLat = lat2 - lat1;
var dLon = lon2 - lon1;

var a = Math.sin(dLat/2) * Math.sin(dLat/2) +
Math.cos(lat1) * Math.cos(lat2) *
Math.sin(dLon/2) * Math.sin(dLon/2);
var c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1-a));
var d = R * c;
return d.toPrecisionFixed(precision);
}


/**
* Returns the (initial) bearing from this point to the supplied point, in degrees
* see http://williams.best.vwh.net/avform.htm#Crs
*
* @param {LatLon} point: Latitude/longitude of destination point
* @returns {Number} Initial bearing in degrees from North
*/
bearingTo = function(point) {
var lat1 = this._lat.toRad(), lat2 = point._lat.toRad();
var dLon = (point._lon-this._lon).toRad();

var y = Math.sin(dLon) * Math.cos(lat2);
var x = Math.cos(lat1)*Math.sin(lat2) -
Math.sin(lat1)*Math.cos(lat2)*Math.cos(dLon);
var brng = Math.atan2(y, x);

return (brng.toDeg()+360) % 360;
}


/**
* Returns final bearing arriving at supplied destination point from this point; the final bearing
* will differ from the initial bearing by varying degrees according to distance and latitude
*
* @param {LatLon} point: Latitude/longitude of destination point
* @returns {Number} Final bearing in degrees from North
*/
finalBearingTo = function(point) {
// get initial bearing from supplied point back to this point...
var lat1 = point._lat.toRad(), lat2 = this._lat.toRad();
var dLon = (this._lon-point._lon).toRad();

var y = Math.sin(dLon) * Math.cos(lat2);
var x = Math.cos(lat1)*Math.sin(lat2) -
Math.sin(lat1)*Math.cos(lat2)*Math.cos(dLon);
var brng = Math.atan2(y, x);

// ... & reverse it by adding 180°
return (brng.toDeg()+180) % 360;
}


/**
* Returns the midpoint between this point and the supplied point.
* see http://mathforum.org/library/drmath/view/51822.html for derivation
*
* @param {LatLon} point: Latitude/longitude of destination point
* @returns {LatLon} Midpoint between this point and the supplied point
*/
midpointTo = function(point) {
lat1 = this._lat.toRad(), lon1 = this._lon.toRad();
lat2 = point._lat.toRad();
var dLon = (point._lon-this._lon).toRad();

var Bx = Math.cos(lat2) * Math.cos(dLon);
var By = Math.cos(lat2) * Math.sin(dLon);

lat3 = Math.atan2(Math.sin(lat1)+Math.sin(lat2),
Math.sqrt( (Math.cos(lat1)+Bx)*(Math.cos(lat1)+Bx) + By*By) );
lon3 = lon1 + Math.atan2(By, Math.cos(lat1) + Bx);

return new LatLon(lat3.toDeg(), lon3.toDeg());
}


/**
* Returns the destination point from this point having travelled the given distance (in km) on the
* given initial bearing (bearing may vary before destination is reached)
*
* see http://williams.best.vwh.net/avform.htm#LL
*
* @param {Number} brng: Initial bearing in degrees
* @param {Number} dist: Distance in km
* @returns {LatLon} Destination point
*/
destinationPoint = function(brng, dist) {
dist = typeof(dist)=='number' ? dist : typeof(dist)=='string' && dist.trim()!='' ? +dist : NaN;
dist = dist/this._radius; // convert dist to angular distance in radians
brng = brng.toRad(); //
var lat1 = this._lat.toRad(), lon1 = this._lon.toRad();

var lat2 = Math.asin( Math.sin(lat1)*Math.cos(dist) +
Math.cos(lat1)*Math.sin(dist)*Math.cos(brng) );
var lon2 = lon1 + Math.atan2(Math.sin(brng)*Math.sin(dist)*Math.cos(lat1),
Math.cos(dist)-Math.sin(lat1)*Math.sin(lat2));
lon2 = (lon2+3*Math.PI)%(2*Math.PI) - Math.PI; // normalise to -180...+180

return new LatLon(lat2.toDeg(), lon2.toDeg());
}


/**
* Returns the point of intersection of two paths defined by point and bearing
*
* see http://williams.best.vwh.net/avform.htm#Intersection
*
* @param {LatLon} p1: First point
* @param {Number} brng1: Initial bearing from first point
* @param {LatLon} p2: Second point
* @param {Number} brng2: Initial bearing from second point
* @returns {LatLon} Destination point (null if no unique intersection defined)
*/
intersection = function(p1, brng1, p2, brng2) {
brng1 = typeof brng1 == 'number' ? brng1 : typeof brng1 == 'string' && trim(brng1)!='' ? +brng1 : NaN;
brng2 = typeof brng2 == 'number' ? brng2 : typeof brng2 == 'string' && trim(brng2)!='' ? +brng2 : NaN;
lat1 = p1._lat.toRad(), lon1 = p1._lon.toRad();
lat2 = p2._lat.toRad(), lon2 = p2._lon.toRad();
brng13 = brng1.toRad(), brng23 = brng2.toRad();
dLat = lat2-lat1, dLon = lon2-lon1;

dist12 = 2*Math.asin( Math.sqrt( Math.sin(dLat/2)*Math.sin(dLat/2) +
Math.cos(lat1)*Math.cos(lat2)*Math.sin(dLon/2)*Math.sin(dLon/2) ) );
if (dist12 == 0) return null;

// initial/final bearings between points
brngA = Math.acos( ( Math.sin(lat2) - Math.sin(lat1)*Math.cos(dist12) ) /
( Math.sin(dist12)*Math.cos(lat1) ) );
if (isNaN(brngA)) brngA = 0; // protect against rounding
brngB = Math.acos( ( Math.sin(lat1) - Math.sin(lat2)*Math.cos(dist12) ) /
( Math.sin(dist12)*Math.cos(lat2) ) );

if (Math.sin(lon2-lon1) > 0) {
brng12 = brngA;
brng21 = 2*Math.PI - brngB;
} else {
brng12 = 2*Math.PI - brngA;
brng21 = brngB;
}

alpha1 = (brng13 - brng12 + Math.PI) % (2*Math.PI) - Math.PI; // angle 2-1-3
alpha2 = (brng21 - brng23 + Math.PI) % (2*Math.PI) - Math.PI; // angle 1-2-3

if (Math.sin(alpha1)==0 && Math.sin(alpha2)==0) return null; // infinite intersections
if (Math.sin(alpha1)*Math.sin(alpha2) < 0) return null; // ambiguous intersection

//alpha1 = Math.abs(alpha1);
//alpha2 = Math.abs(alpha2);
// ... Ed Williams takes abs of alpha1/alpha2, but seems to break calculation?

alpha3 = Math.acos( -Math.cos(alpha1)*Math.cos(alpha2) +
Math.sin(alpha1)*Math.sin(alpha2)*Math.cos(dist12) );
dist13 = Math.atan2( Math.sin(dist12)*Math.sin(alpha1)*Math.sin(alpha2),
Math.cos(alpha2)+Math.cos(alpha1)*Math.cos(alpha3) )
lat3 = Math.asin( Math.sin(lat1)*Math.cos(dist13) +
Math.cos(lat1)*Math.sin(dist13)*Math.cos(brng13) );
dLon13 = Math.atan2( Math.sin(brng13)*Math.sin(dist13)*Math.cos(lat1),
Math.cos(dist13)-Math.sin(lat1)*Math.sin(lat3) );
lon3 = lon1+dLon13;
lon3 = (lon3+Math.PI) % (2*Math.PI) - Math.PI; // normalise to -180..180º

return new LatLon(lat3.toDeg(), lon3.toDeg());
}


/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

/**
* Returns the distance from this point to the supplied point, in km, travelling along a rhumb line
*
* see http://williams.best.vwh.net/avform.htm#Rhumb
*
* @param {LatLon} point: Latitude/longitude of destination point
* @returns {Number} Distance in km between this point and destination point
*/
rhumbDistanceTo = function(point) {
var R = this._radius;
var lat1 = this._lat.toRad(), lat2 = point._lat.toRad();
var dLat = (point._lat-this._lat).toRad();
var dLon = Math.abs(point._lon-this._lon).toRad();

var dPhi = Math.log(Math.tan(lat2/2+Math.PI/4)/Math.tan(lat1/2+Math.PI/4));
var q = (!isNaN(dLat/dPhi)) ? dLat/dPhi : Math.cos(lat1); // E-W line gives dPhi=0
// if dLon over 180° take shorter rhumb across 180° meridian:
if (dLon > Math.PI) dLon = 2*Math.PI - dLon;
var dist = Math.sqrt(dLat*dLat + q*q*dLon*dLon) * R;

return dist.toPrecisionFixed(4); // 4 sig figs reflects typical 0.3% accuracy of spherical model
}

/**
* Returns the bearing from this point to the supplied point along a rhumb line, in degrees
*
* @param {LatLon} point: Latitude/longitude of destination point
* @returns {Number} Bearing in degrees from North
*/
LatLon.prototype.rhumbBearingTo = function(point) {
var lat1 = this._lat.toRad(), lat2 = point._lat.toRad();
var dLon = (point._lon-this._lon).toRad();

var dPhi = Math.log(Math.tan(lat2/2+Math.PI/4)/Math.tan(lat1/2+Math.PI/4));
if (Math.abs(dLon) > Math.PI) dLon = dLon>0 ? -(2*Math.PI-dLon) : (2*Math.PI+dLon);
var brng = Math.atan2(dLon, dPhi);

return (brng.toDeg()+360) % 360;
}

/**
* Returns the destination point from this point having travelled the given distance (in km) on the
* given bearing along a rhumb line
*
* @param {Number} brng: Bearing in degrees from North
* @param {Number} dist: Distance in km
* @returns {LatLon} Destination point
*/
rhumbDestinationPoint = function(brng, dist) {
var R = this._radius;
var d = parseFloat(dist)/R; // d = angular distance covered on earth's surface
var lat1 = this._lat.toRad(), lon1 = this._lon.toRad();
brng = brng.toRad();

var lat2 = lat1 + d*Math.cos(brng);
var dLat = lat2-lat1;
var dPhi = Math.log(Math.tan(lat2/2+Math.PI/4)/Math.tan(lat1/2+Math.PI/4));
var q = (!isNaN(dLat/dPhi)) ? dLat/dPhi : Math.cos(lat1); // E-W line gives dPhi=0
var dLon = d*Math.sin(brng)/q;
// check for some daft bugger going past the pole
if (Math.abs(lat2) > Math.PI/2) lat2 = lat2>0 ? Math.PI-lat2 : -(Math.PI-lat2);
lon2 = (lon1+dLon+3*Math.PI)%(2*Math.PI) - Math.PI;

return new LatLon(lat2.toDeg(), lon2.toDeg());
}

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */


/**
* Returns the latitude of this point; signed numeric degrees if no format, otherwise format & dp
* as per Geo.toLat()
*
* @param {String} [format]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to display
* @returns {Number|String} Numeric degrees if no format specified, otherwise deg/min/sec
*
* @requires Geo
*/
lat = function(format, dp) {
if (typeof format == 'undefined') return this._lat;

return Geo.toLat(this._lat, format, dp);
}

/**
* Returns the longitude of this point; signed numeric degrees if no format, otherwise format & dp
* as per Geo.toLon()
*
* @param {String} [format]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to display
* @returns {Number|String} Numeric degrees if no format specified, otherwise deg/min/sec
*
* @requires Geo
*/
lon = function(format, dp) {
if (typeof format == 'undefined') return this._lon;

return Geo.toLon(this._lon, format, dp);
}

/**
* Returns a string representation of this point; format and dp as per lat()/lon()
*
* @param {String} [format]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to display
* @returns {String} Comma-separated latitude/longitude
*
* @requires Geo
*/
toString = function(format, dp) {
if (typeof format == 'undefined') format = 'dms';

if (isNaN(this._lat) || isNaN(this._lon)) return '-,-';

return Geo.toLat(this._lat, format, dp) + ', ' + Geo.toLon(this._lon, format, dp);
}

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

// ---- extend Number object with methods for converting degrees/radians

/** Converts numeric degrees to radians */
if (typeof(Number.prototype.toRad) === "undefined") {
Number.prototype.toRad = function() {
return this * Math.PI / 180;
}
}

/** Converts radians to numeric (signed) degrees */
if (typeof(Number.prototype.toDeg) === "undefined") {
Number.prototype.toDeg = function() {
return this * 180 / Math.PI;
}
}

/**
* Formats the significant digits of a number, using only fixed-point notation (no exponential)
*
* @param {Number} precision: Number of significant digits to appear in the returned string
* @returns {String} A string representation of number which contains precision significant digits
*/
if (typeof(Number.prototype.toPrecisionFixed) === "undefined") {
Number.prototype.toPrecisionFixed = function(precision) {
if (isNaN(this)) return 'NaN';
var numb = this < 0 ? -this : this; // can't take log of -ve number...
var sign = this < 0 ? '-' : '';

if (numb == 0) { n = '0.'; while (precision--) n += '0'; return n }; // can't take log of zero

var scale = Math.ceil(Math.log(numb)*Math.LOG10E); // no of digits before decimal
var n = String(Math.round(numb * Math.pow(10, precision-scale)));
if (scale > 0) { // add trailing zeros & insert decimal as required
l = scale - n.length;
while (l-- > 0) n = n + '0';
if (scale < n.length) n = n.slice(0,scale) + '.' + n.slice(scale);
} else { // prefix decimal and leading zeros if required
while (scale++ < 0) n = '0' + n;
n = '0.' + n;
}
return sign + n;
}
}

/** Trims whitespace from string (q.v. blog.stevenlevithan.com/archives/faster-trim-javascript) */
if (typeof(String.prototype.trim) === "undefined") {
String.prototype.trim = function() {
return String(this).replace(/^\s\s*/, '').replace(/\s\s*$/, '');
}
}


/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
/* Geodesy representation conversion functions (c) Chris Veness 2002-2010 */
/* - www.movable-type.co.uk/scripts/latlong.html */
/* */
/* Sample usage: */
/* var lat = Geo.parseDMS('51° 28′ 40.12″ N'); */
/* var lon = Geo.parseDMS('000° 00′ 05.31″ W'); */
/* var p1 = new LatLon(lat, lon); */
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

var Geo = {}; // Geo namespace, representing static class

/**
* Parses string representing degrees/minutes/seconds into numeric degrees
*
* This is very flexible on formats, allowing signed decimal degrees, or deg-min-sec optionally
* suffixed by compass direction (NSEW). A variety of separators are accepted (eg 3º 37' 09"W)
* or fixed-width format without separators (eg 0033709W). Seconds and minutes may be omitted.
* (Note minimal validation is done).
*
* @param {String|Number} dmsStr: Degrees or deg/min/sec in variety of formats
* @returns {Number} Degrees as decimal number
* @throws {TypeError} dmsStr is an object, perhaps DOM object without .value?
*/
parseDMS = function(dmsStr) {
if (typeof deg == 'object') throw new TypeError('Geo.parseDMS - dmsStr is [DOM?] object');

// check for signed decimal degrees without NSEW, if so return it directly
if (typeof dmsStr === 'number' && isFinite(dmsStr)) return Number(dmsStr);

// strip off any sign or compass dir'n & split out separate d/m/s
var dms = String(dmsStr).trim().replace(/^-/,'').replace(/[NSEW]$/i,'').split(/[^0-9.,]+/);
if (dms[dms.length-1]=='') dms.splice(dms.length-1); // from trailing symbol

if (dms == '') return NaN;

// and convert to decimal degrees...
switch (dms.length) {
case 3: // interpret 3-part result as d/m/s
var deg = dms[0]/1 + dms[1]/60 + dms[2]/3600;
break;
case 2: // interpret 2-part result as d/m
var deg = dms[0]/1 + dms[1]/60;
break;
case 1: // just d (possibly decimal) or non-separated dddmmss
var deg = dms[0];
// check for fixed-width unseparated format eg 0033709W
if (/[NS]/i.test(dmsStr)) deg = '0' + deg; // - normalise N/S to 3-digit degrees
if (/[0-9]{7}/.test(deg)) deg = deg.slice(0,3)/1 + deg.slice(3,5)/60 + deg.slice(5)/3600;
break;
default:
return NaN;
}
if (/^-|[WS]$/i.test(dmsStr.trim())) deg = -deg; // take '-', west and south as -ve
return Number(deg);
}

/**
* Convert decimal degrees to deg/min/sec format
* - degree, prime, double-prime symbols are added, but sign is discarded, though no compass
* direction is added
*
* @private
* @param {Number} deg: Degrees
* @param {String} [format=dms]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to use - default 0 for dms, 2 for dm, 4 for d
* @returns {String} deg formatted as deg/min/secs according to specified format
* @throws {TypeError} deg is an object, perhaps DOM object without .value?
*/
toDMS = function(deg, format, dp) {
if (typeof deg == 'object') throw new TypeError('Geo.toDMS - deg is [DOM?] object');
if (isNaN(deg)) return 'NaN'; // give up here if we can't make a number from deg

// default values
if (typeof format == 'undefined') format = 'dms';
if (typeof dp == 'undefined') {
switch (format) {
case 'd': dp = 4; break;
case 'dm': dp = 2; break;
case 'dms': dp = 0; break;
default: format = 'dms'; dp = 0; // be forgiving on invalid format
}
}

deg = Math.abs(deg); // (unsigned result ready for appending compass dir'n)

switch (format) {
case 'd':
d = deg.toFixed(dp); // round degrees
if (d<100) d = '0' + d; // pad with leading zeros
if (d<10) d = '0' + d;
dms = d + '\u00B0'; // add º symbol
break;
case 'dm':
var min = (deg*60).toFixed(dp); // convert degrees to minutes & round
var d = Math.floor(min / 60); // get component deg/min
var m = (min % 60).toFixed(dp); // pad with trailing zeros
if (d<100) d = '0' + d; // pad with leading zeros
if (d<10) d = '0' + d;
if (m<10) m = '0' + m;
dms = d + '\u00B0' + m + '\u2032'; // add º, ' symbols
break;
case 'dms':
var sec = (deg*3600).toFixed(dp); // convert degrees to seconds & round
var d = Math.floor(sec / 3600); // get component deg/min/sec
var m = Math.floor(sec/60) % 60;
var s = (sec % 60).toFixed(dp); // pad with trailing zeros
if (d<100) d = '0' + d; // pad with leading zeros
if (d<10) d = '0' + d;
if (m<10) m = '0' + m;
if (s<10) s = '0' + s;
dms = d + '\u00B0' + m + '\u2032' + s + '\u2033'; // add º, ', " symbols
break;
}

return dms;
}

/**
* Convert numeric degrees to deg/min/sec latitude (suffixed with N/S)
*
* @param {Number} deg: Degrees
* @param {String} [format=dms]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to use - default 0 for dms, 2 for dm, 4 for d
* @returns {String} Deg/min/seconds
*/
toLat = function(deg, format, dp) {
var lat = Geo.toDMS(deg, format, dp);
return lat=='' ? '' : lat.slice(1) + (deg<0 ? 'S' : 'N'); // knock off initial '0' for lat!
}

/**
* Convert numeric degrees to deg/min/sec longitude (suffixed with E/W)
*
* @param {Number} deg: Degrees
* @param {String} [format=dms]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to use - default 0 for dms, 2 for dm, 4 for d
* @returns {String} Deg/min/seconds
*/
toLon = function(deg, format, dp) {
var lon = Geo.toDMS(deg, format, dp);
return lon=='' ? '' : lon + (deg<0 ? 'W' : 'E');
}

/**
* Convert numeric degrees to deg/min/sec as a bearing (0º..360º)
*
* @param {Number} deg: Degrees
* @param {String} [format=dms]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to use - default 0 for dms, 2 for dm, 4 for d
* @returns {String} Deg/min/seconds
*/
toBrng = function(deg, format, dp) {
deg = (Number(deg)+360) % 360; // normalise -ve values to 180º..360º
var brng = Geo.toDMS(deg, format, dp);
return brng.replace('360', '0'); // just in case rounding took us up to 360º!
}


Load some example data:


GEO:
LOAD * INLINE [
ROUTE, LAT_1, LONG_1, LAT_2, LONG_2
R1, 53 09 02N, 001 50 40W, 52 12 19N, 000 08 33W
];


DIST:
LOAD
*,
distVincenty(parseDMS(LAT_1), parseDMS(LONG_1), parseDMS(LAT_2), parseDMS(LONG_2)) as DIST
resident GEO;

drop table GEO;


-Alex

View solution in original post

5 Replies
Not applicable
Author

hi,

Take a look at the source of the linked HTML page == there is Javascript with the formula 🙂 !

Makea QlikView file with a module written in Jscript. Paste the code there.


/**
* Calculates geodetic distance between two points specified by latitude/longitude using
* Vincenty inverse formula for ellipsoids
*
* @param {Number} lat1, lon1: first point in decimal degrees
* @param {Number} lat2, lon2: second point in decimal degrees
* @returns (Number} distance in metres between points
*/
function distVincenty(lat1, lon1, lat2, lon2) {
var a = 6378137, b = 6356752.3142, 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;

// note: to return initial/final bearings in addition to distance, use something like:
var fwdAz = Math.atan2(cosU2*sinLambda, cosU1*sinU2-sinU1*cosU2*cosLambda);
var revAz = Math.atan2(cosU1*sinLambda, -sinU1*cosU2+cosU1*sinU2*cosLambda);
return { distance: s, initialBearing: fwdAz.toDeg(), finalBearing: revAz.toDeg() };
}

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

/** extend Number object with methods for converting degrees/radians */
toRad = function() { // convert degrees to radians
return this * Math.PI / 180;
}
toDeg = function() { // convert radians to degrees
return this * 180 / Math.PI;
}

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
/* Latitude/longitude spherical geodesy formulae & scripts (c) Chris Veness 2002-2010 */
/* - www.movable-type.co.uk/scripts/latlong.html */
/* */
/* Sample usage: */
/* var p1 = new LatLon(51.5136, -0.0983); */
/* var p2 = new LatLon(51.4778, -0.0015); */
/* var dist = p1.distanceTo(p2); // in km */
/* var brng = p1.bearingTo(p2); // in degrees clockwise from north */
/* ... etc */
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
/* Note that minimal error checking is performed in this example code! */
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */


/**
* Creates a point on the earth's surface at the supplied latitude / longitude
*
* @constructor
* @param {Number} lat: latitude in numeric degrees
* @param {Number} lon: longitude in numeric degrees
* @param {Number} [rad=6371]: radius of earth if different value is required from standard 6,371km
*/
function LatLon(lat, lon, rad) {
if (typeof(rad) == 'undefined') rad = 6371; // earth's mean radius in km
// only accept numbers or valid numeric strings
this._lat = typeof(lat)=='number' ? lat : typeof(lat)=='string' && lat.trim()!='' ? +lat : NaN;
this._lon = typeof(lat)=='number' ? lon : typeof(lon)=='string' && lon.trim()!='' ? +lon : NaN;
this._radius = typeof(rad)=='number' ? rad : typeof(rad)=='string' && trim(lon)!='' ? +rad : NaN;
}


/**
* Returns the distance from this point to the supplied point, in km
* (using Haversine formula)
*
* from: Haversine formula - R. W. Sinnott, "Virtues of the Haversine",
* Sky and Telescope, vol 68, no 2, 1984
*
* @param {LatLon} point: Latitude/longitude of destination point
* @param {Number} [precision=4]: no of significant digits to use for returned value
* @returns {Number} Distance in km between this point and destination point
*/
distanceTo = function(point, precision) {
// default 4 sig figs reflects typical 0.3% accuracy of spherical model
if (typeof precision == 'undefined') precision = 4;

var R = this._radius;
var lat1 = this._lat.toRad(), lon1 = this._lon.toRad();
var lat2 = point._lat.toRad(), lon2 = point._lon.toRad();
var dLat = lat2 - lat1;
var dLon = lon2 - lon1;

var a = Math.sin(dLat/2) * Math.sin(dLat/2) +
Math.cos(lat1) * Math.cos(lat2) *
Math.sin(dLon/2) * Math.sin(dLon/2);
var c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1-a));
var d = R * c;
return d.toPrecisionFixed(precision);
}


/**
* Returns the (initial) bearing from this point to the supplied point, in degrees
* see http://williams.best.vwh.net/avform.htm#Crs
*
* @param {LatLon} point: Latitude/longitude of destination point
* @returns {Number} Initial bearing in degrees from North
*/
bearingTo = function(point) {
var lat1 = this._lat.toRad(), lat2 = point._lat.toRad();
var dLon = (point._lon-this._lon).toRad();

var y = Math.sin(dLon) * Math.cos(lat2);
var x = Math.cos(lat1)*Math.sin(lat2) -
Math.sin(lat1)*Math.cos(lat2)*Math.cos(dLon);
var brng = Math.atan2(y, x);

return (brng.toDeg()+360) % 360;
}


/**
* Returns final bearing arriving at supplied destination point from this point; the final bearing
* will differ from the initial bearing by varying degrees according to distance and latitude
*
* @param {LatLon} point: Latitude/longitude of destination point
* @returns {Number} Final bearing in degrees from North
*/
finalBearingTo = function(point) {
// get initial bearing from supplied point back to this point...
var lat1 = point._lat.toRad(), lat2 = this._lat.toRad();
var dLon = (this._lon-point._lon).toRad();

var y = Math.sin(dLon) * Math.cos(lat2);
var x = Math.cos(lat1)*Math.sin(lat2) -
Math.sin(lat1)*Math.cos(lat2)*Math.cos(dLon);
var brng = Math.atan2(y, x);

// ... & reverse it by adding 180°
return (brng.toDeg()+180) % 360;
}


/**
* Returns the midpoint between this point and the supplied point.
* see http://mathforum.org/library/drmath/view/51822.html for derivation
*
* @param {LatLon} point: Latitude/longitude of destination point
* @returns {LatLon} Midpoint between this point and the supplied point
*/
midpointTo = function(point) {
lat1 = this._lat.toRad(), lon1 = this._lon.toRad();
lat2 = point._lat.toRad();
var dLon = (point._lon-this._lon).toRad();

var Bx = Math.cos(lat2) * Math.cos(dLon);
var By = Math.cos(lat2) * Math.sin(dLon);

lat3 = Math.atan2(Math.sin(lat1)+Math.sin(lat2),
Math.sqrt( (Math.cos(lat1)+Bx)*(Math.cos(lat1)+Bx) + By*By) );
lon3 = lon1 + Math.atan2(By, Math.cos(lat1) + Bx);

return new LatLon(lat3.toDeg(), lon3.toDeg());
}


/**
* Returns the destination point from this point having travelled the given distance (in km) on the
* given initial bearing (bearing may vary before destination is reached)
*
* see http://williams.best.vwh.net/avform.htm#LL
*
* @param {Number} brng: Initial bearing in degrees
* @param {Number} dist: Distance in km
* @returns {LatLon} Destination point
*/
destinationPoint = function(brng, dist) {
dist = typeof(dist)=='number' ? dist : typeof(dist)=='string' && dist.trim()!='' ? +dist : NaN;
dist = dist/this._radius; // convert dist to angular distance in radians
brng = brng.toRad(); //
var lat1 = this._lat.toRad(), lon1 = this._lon.toRad();

var lat2 = Math.asin( Math.sin(lat1)*Math.cos(dist) +
Math.cos(lat1)*Math.sin(dist)*Math.cos(brng) );
var lon2 = lon1 + Math.atan2(Math.sin(brng)*Math.sin(dist)*Math.cos(lat1),
Math.cos(dist)-Math.sin(lat1)*Math.sin(lat2));
lon2 = (lon2+3*Math.PI)%(2*Math.PI) - Math.PI; // normalise to -180...+180

return new LatLon(lat2.toDeg(), lon2.toDeg());
}


/**
* Returns the point of intersection of two paths defined by point and bearing
*
* see http://williams.best.vwh.net/avform.htm#Intersection
*
* @param {LatLon} p1: First point
* @param {Number} brng1: Initial bearing from first point
* @param {LatLon} p2: Second point
* @param {Number} brng2: Initial bearing from second point
* @returns {LatLon} Destination point (null if no unique intersection defined)
*/
intersection = function(p1, brng1, p2, brng2) {
brng1 = typeof brng1 == 'number' ? brng1 : typeof brng1 == 'string' && trim(brng1)!='' ? +brng1 : NaN;
brng2 = typeof brng2 == 'number' ? brng2 : typeof brng2 == 'string' && trim(brng2)!='' ? +brng2 : NaN;
lat1 = p1._lat.toRad(), lon1 = p1._lon.toRad();
lat2 = p2._lat.toRad(), lon2 = p2._lon.toRad();
brng13 = brng1.toRad(), brng23 = brng2.toRad();
dLat = lat2-lat1, dLon = lon2-lon1;

dist12 = 2*Math.asin( Math.sqrt( Math.sin(dLat/2)*Math.sin(dLat/2) +
Math.cos(lat1)*Math.cos(lat2)*Math.sin(dLon/2)*Math.sin(dLon/2) ) );
if (dist12 == 0) return null;

// initial/final bearings between points
brngA = Math.acos( ( Math.sin(lat2) - Math.sin(lat1)*Math.cos(dist12) ) /
( Math.sin(dist12)*Math.cos(lat1) ) );
if (isNaN(brngA)) brngA = 0; // protect against rounding
brngB = Math.acos( ( Math.sin(lat1) - Math.sin(lat2)*Math.cos(dist12) ) /
( Math.sin(dist12)*Math.cos(lat2) ) );

if (Math.sin(lon2-lon1) > 0) {
brng12 = brngA;
brng21 = 2*Math.PI - brngB;
} else {
brng12 = 2*Math.PI - brngA;
brng21 = brngB;
}

alpha1 = (brng13 - brng12 + Math.PI) % (2*Math.PI) - Math.PI; // angle 2-1-3
alpha2 = (brng21 - brng23 + Math.PI) % (2*Math.PI) - Math.PI; // angle 1-2-3

if (Math.sin(alpha1)==0 && Math.sin(alpha2)==0) return null; // infinite intersections
if (Math.sin(alpha1)*Math.sin(alpha2) < 0) return null; // ambiguous intersection

//alpha1 = Math.abs(alpha1);
//alpha2 = Math.abs(alpha2);
// ... Ed Williams takes abs of alpha1/alpha2, but seems to break calculation?

alpha3 = Math.acos( -Math.cos(alpha1)*Math.cos(alpha2) +
Math.sin(alpha1)*Math.sin(alpha2)*Math.cos(dist12) );
dist13 = Math.atan2( Math.sin(dist12)*Math.sin(alpha1)*Math.sin(alpha2),
Math.cos(alpha2)+Math.cos(alpha1)*Math.cos(alpha3) )
lat3 = Math.asin( Math.sin(lat1)*Math.cos(dist13) +
Math.cos(lat1)*Math.sin(dist13)*Math.cos(brng13) );
dLon13 = Math.atan2( Math.sin(brng13)*Math.sin(dist13)*Math.cos(lat1),
Math.cos(dist13)-Math.sin(lat1)*Math.sin(lat3) );
lon3 = lon1+dLon13;
lon3 = (lon3+Math.PI) % (2*Math.PI) - Math.PI; // normalise to -180..180º

return new LatLon(lat3.toDeg(), lon3.toDeg());
}


/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

/**
* Returns the distance from this point to the supplied point, in km, travelling along a rhumb line
*
* see http://williams.best.vwh.net/avform.htm#Rhumb
*
* @param {LatLon} point: Latitude/longitude of destination point
* @returns {Number} Distance in km between this point and destination point
*/
rhumbDistanceTo = function(point) {
var R = this._radius;
var lat1 = this._lat.toRad(), lat2 = point._lat.toRad();
var dLat = (point._lat-this._lat).toRad();
var dLon = Math.abs(point._lon-this._lon).toRad();

var dPhi = Math.log(Math.tan(lat2/2+Math.PI/4)/Math.tan(lat1/2+Math.PI/4));
var q = (!isNaN(dLat/dPhi)) ? dLat/dPhi : Math.cos(lat1); // E-W line gives dPhi=0
// if dLon over 180° take shorter rhumb across 180° meridian:
if (dLon > Math.PI) dLon = 2*Math.PI - dLon;
var dist = Math.sqrt(dLat*dLat + q*q*dLon*dLon) * R;

return dist.toPrecisionFixed(4); // 4 sig figs reflects typical 0.3% accuracy of spherical model
}

/**
* Returns the bearing from this point to the supplied point along a rhumb line, in degrees
*
* @param {LatLon} point: Latitude/longitude of destination point
* @returns {Number} Bearing in degrees from North
*/
LatLon.prototype.rhumbBearingTo = function(point) {
var lat1 = this._lat.toRad(), lat2 = point._lat.toRad();
var dLon = (point._lon-this._lon).toRad();

var dPhi = Math.log(Math.tan(lat2/2+Math.PI/4)/Math.tan(lat1/2+Math.PI/4));
if (Math.abs(dLon) > Math.PI) dLon = dLon>0 ? -(2*Math.PI-dLon) : (2*Math.PI+dLon);
var brng = Math.atan2(dLon, dPhi);

return (brng.toDeg()+360) % 360;
}

/**
* Returns the destination point from this point having travelled the given distance (in km) on the
* given bearing along a rhumb line
*
* @param {Number} brng: Bearing in degrees from North
* @param {Number} dist: Distance in km
* @returns {LatLon} Destination point
*/
rhumbDestinationPoint = function(brng, dist) {
var R = this._radius;
var d = parseFloat(dist)/R; // d = angular distance covered on earth's surface
var lat1 = this._lat.toRad(), lon1 = this._lon.toRad();
brng = brng.toRad();

var lat2 = lat1 + d*Math.cos(brng);
var dLat = lat2-lat1;
var dPhi = Math.log(Math.tan(lat2/2+Math.PI/4)/Math.tan(lat1/2+Math.PI/4));
var q = (!isNaN(dLat/dPhi)) ? dLat/dPhi : Math.cos(lat1); // E-W line gives dPhi=0
var dLon = d*Math.sin(brng)/q;
// check for some daft bugger going past the pole
if (Math.abs(lat2) > Math.PI/2) lat2 = lat2>0 ? Math.PI-lat2 : -(Math.PI-lat2);
lon2 = (lon1+dLon+3*Math.PI)%(2*Math.PI) - Math.PI;

return new LatLon(lat2.toDeg(), lon2.toDeg());
}

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */


/**
* Returns the latitude of this point; signed numeric degrees if no format, otherwise format & dp
* as per Geo.toLat()
*
* @param {String} [format]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to display
* @returns {Number|String} Numeric degrees if no format specified, otherwise deg/min/sec
*
* @requires Geo
*/
lat = function(format, dp) {
if (typeof format == 'undefined') return this._lat;

return Geo.toLat(this._lat, format, dp);
}

/**
* Returns the longitude of this point; signed numeric degrees if no format, otherwise format & dp
* as per Geo.toLon()
*
* @param {String} [format]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to display
* @returns {Number|String} Numeric degrees if no format specified, otherwise deg/min/sec
*
* @requires Geo
*/
lon = function(format, dp) {
if (typeof format == 'undefined') return this._lon;

return Geo.toLon(this._lon, format, dp);
}

/**
* Returns a string representation of this point; format and dp as per lat()/lon()
*
* @param {String} [format]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to display
* @returns {String} Comma-separated latitude/longitude
*
* @requires Geo
*/
toString = function(format, dp) {
if (typeof format == 'undefined') format = 'dms';

if (isNaN(this._lat) || isNaN(this._lon)) return '-,-';

return Geo.toLat(this._lat, format, dp) + ', ' + Geo.toLon(this._lon, format, dp);
}

/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

// ---- extend Number object with methods for converting degrees/radians

/** Converts numeric degrees to radians */
if (typeof(Number.prototype.toRad) === "undefined") {
Number.prototype.toRad = function() {
return this * Math.PI / 180;
}
}

/** Converts radians to numeric (signed) degrees */
if (typeof(Number.prototype.toDeg) === "undefined") {
Number.prototype.toDeg = function() {
return this * 180 / Math.PI;
}
}

/**
* Formats the significant digits of a number, using only fixed-point notation (no exponential)
*
* @param {Number} precision: Number of significant digits to appear in the returned string
* @returns {String} A string representation of number which contains precision significant digits
*/
if (typeof(Number.prototype.toPrecisionFixed) === "undefined") {
Number.prototype.toPrecisionFixed = function(precision) {
if (isNaN(this)) return 'NaN';
var numb = this < 0 ? -this : this; // can't take log of -ve number...
var sign = this < 0 ? '-' : '';

if (numb == 0) { n = '0.'; while (precision--) n += '0'; return n }; // can't take log of zero

var scale = Math.ceil(Math.log(numb)*Math.LOG10E); // no of digits before decimal
var n = String(Math.round(numb * Math.pow(10, precision-scale)));
if (scale > 0) { // add trailing zeros & insert decimal as required
l = scale - n.length;
while (l-- > 0) n = n + '0';
if (scale < n.length) n = n.slice(0,scale) + '.' + n.slice(scale);
} else { // prefix decimal and leading zeros if required
while (scale++ < 0) n = '0' + n;
n = '0.' + n;
}
return sign + n;
}
}

/** Trims whitespace from string (q.v. blog.stevenlevithan.com/archives/faster-trim-javascript) */
if (typeof(String.prototype.trim) === "undefined") {
String.prototype.trim = function() {
return String(this).replace(/^\s\s*/, '').replace(/\s\s*$/, '');
}
}


/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
/* Geodesy representation conversion functions (c) Chris Veness 2002-2010 */
/* - www.movable-type.co.uk/scripts/latlong.html */
/* */
/* Sample usage: */
/* var lat = Geo.parseDMS('51° 28′ 40.12″ N'); */
/* var lon = Geo.parseDMS('000° 00′ 05.31″ W'); */
/* var p1 = new LatLon(lat, lon); */
/* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

var Geo = {}; // Geo namespace, representing static class

/**
* Parses string representing degrees/minutes/seconds into numeric degrees
*
* This is very flexible on formats, allowing signed decimal degrees, or deg-min-sec optionally
* suffixed by compass direction (NSEW). A variety of separators are accepted (eg 3º 37' 09"W)
* or fixed-width format without separators (eg 0033709W). Seconds and minutes may be omitted.
* (Note minimal validation is done).
*
* @param {String|Number} dmsStr: Degrees or deg/min/sec in variety of formats
* @returns {Number} Degrees as decimal number
* @throws {TypeError} dmsStr is an object, perhaps DOM object without .value?
*/
parseDMS = function(dmsStr) {
if (typeof deg == 'object') throw new TypeError('Geo.parseDMS - dmsStr is [DOM?] object');

// check for signed decimal degrees without NSEW, if so return it directly
if (typeof dmsStr === 'number' && isFinite(dmsStr)) return Number(dmsStr);

// strip off any sign or compass dir'n & split out separate d/m/s
var dms = String(dmsStr).trim().replace(/^-/,'').replace(/[NSEW]$/i,'').split(/[^0-9.,]+/);
if (dms[dms.length-1]=='') dms.splice(dms.length-1); // from trailing symbol

if (dms == '') return NaN;

// and convert to decimal degrees...
switch (dms.length) {
case 3: // interpret 3-part result as d/m/s
var deg = dms[0]/1 + dms[1]/60 + dms[2]/3600;
break;
case 2: // interpret 2-part result as d/m
var deg = dms[0]/1 + dms[1]/60;
break;
case 1: // just d (possibly decimal) or non-separated dddmmss
var deg = dms[0];
// check for fixed-width unseparated format eg 0033709W
if (/[NS]/i.test(dmsStr)) deg = '0' + deg; // - normalise N/S to 3-digit degrees
if (/[0-9]{7}/.test(deg)) deg = deg.slice(0,3)/1 + deg.slice(3,5)/60 + deg.slice(5)/3600;
break;
default:
return NaN;
}
if (/^-|[WS]$/i.test(dmsStr.trim())) deg = -deg; // take '-', west and south as -ve
return Number(deg);
}

/**
* Convert decimal degrees to deg/min/sec format
* - degree, prime, double-prime symbols are added, but sign is discarded, though no compass
* direction is added
*
* @private
* @param {Number} deg: Degrees
* @param {String} [format=dms]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to use - default 0 for dms, 2 for dm, 4 for d
* @returns {String} deg formatted as deg/min/secs according to specified format
* @throws {TypeError} deg is an object, perhaps DOM object without .value?
*/
toDMS = function(deg, format, dp) {
if (typeof deg == 'object') throw new TypeError('Geo.toDMS - deg is [DOM?] object');
if (isNaN(deg)) return 'NaN'; // give up here if we can't make a number from deg

// default values
if (typeof format == 'undefined') format = 'dms';
if (typeof dp == 'undefined') {
switch (format) {
case 'd': dp = 4; break;
case 'dm': dp = 2; break;
case 'dms': dp = 0; break;
default: format = 'dms'; dp = 0; // be forgiving on invalid format
}
}

deg = Math.abs(deg); // (unsigned result ready for appending compass dir'n)

switch (format) {
case 'd':
d = deg.toFixed(dp); // round degrees
if (d<100) d = '0' + d; // pad with leading zeros
if (d<10) d = '0' + d;
dms = d + '\u00B0'; // add º symbol
break;
case 'dm':
var min = (deg*60).toFixed(dp); // convert degrees to minutes & round
var d = Math.floor(min / 60); // get component deg/min
var m = (min % 60).toFixed(dp); // pad with trailing zeros
if (d<100) d = '0' + d; // pad with leading zeros
if (d<10) d = '0' + d;
if (m<10) m = '0' + m;
dms = d + '\u00B0' + m + '\u2032'; // add º, ' symbols
break;
case 'dms':
var sec = (deg*3600).toFixed(dp); // convert degrees to seconds & round
var d = Math.floor(sec / 3600); // get component deg/min/sec
var m = Math.floor(sec/60) % 60;
var s = (sec % 60).toFixed(dp); // pad with trailing zeros
if (d<100) d = '0' + d; // pad with leading zeros
if (d<10) d = '0' + d;
if (m<10) m = '0' + m;
if (s<10) s = '0' + s;
dms = d + '\u00B0' + m + '\u2032' + s + '\u2033'; // add º, ', " symbols
break;
}

return dms;
}

/**
* Convert numeric degrees to deg/min/sec latitude (suffixed with N/S)
*
* @param {Number} deg: Degrees
* @param {String} [format=dms]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to use - default 0 for dms, 2 for dm, 4 for d
* @returns {String} Deg/min/seconds
*/
toLat = function(deg, format, dp) {
var lat = Geo.toDMS(deg, format, dp);
return lat=='' ? '' : lat.slice(1) + (deg<0 ? 'S' : 'N'); // knock off initial '0' for lat!
}

/**
* Convert numeric degrees to deg/min/sec longitude (suffixed with E/W)
*
* @param {Number} deg: Degrees
* @param {String} [format=dms]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to use - default 0 for dms, 2 for dm, 4 for d
* @returns {String} Deg/min/seconds
*/
toLon = function(deg, format, dp) {
var lon = Geo.toDMS(deg, format, dp);
return lon=='' ? '' : lon + (deg<0 ? 'W' : 'E');
}

/**
* Convert numeric degrees to deg/min/sec as a bearing (0º..360º)
*
* @param {Number} deg: Degrees
* @param {String} [format=dms]: Return value as 'd', 'dm', 'dms'
* @param {Number} [dp=0|2|4]: No of decimal places to use - default 0 for dms, 2 for dm, 4 for d
* @returns {String} Deg/min/seconds
*/
toBrng = function(deg, format, dp) {
deg = (Number(deg)+360) % 360; // normalise -ve values to 180º..360º
var brng = Geo.toDMS(deg, format, dp);
return brng.replace('360', '0'); // just in case rounding took us up to 360º!
}


Load some example data:


GEO:
LOAD * INLINE [
ROUTE, LAT_1, LONG_1, LAT_2, LONG_2
R1, 53 09 02N, 001 50 40W, 52 12 19N, 000 08 33W
];


DIST:
LOAD
*,
distVincenty(parseDMS(LAT_1), parseDMS(LONG_1), parseDMS(LAT_2), parseDMS(LONG_2)) as DIST
resident GEO;

drop table GEO;


-Alex

Not applicable
Author

See file from my profile: http://community.qlik.com/members/alxtoth/files/Attached+Files/geo.qvw.aspx

-Alex

Not applicable
Author

Excellent ... many thanks ... now all I need to do is figure out how to use Java script stuff! 😮

Keep watching this space ...

Cheers,

Emma

Not applicable
Author

Awsome! Many thanks - I've got it working 🙂

Now, all I need to do is get approval from the verifier and you'll have made us WGS84 compliant - as everyone should be of course! Thank you so much!!

🙂

Emma

Not applicable
Author

Hi Alex,

Many thanks again for the script and instructions, I had a small (well ... quite large) panic yesterday when I compared the results I got using this to the results I had using the Haversine (spherical model) formula. Luckily, after a bit of investigation I found that all it was, was that the Lat and Long format in my file has the NSEW indicator at the front not the end so although the formula will still run (and give valid results on N and E coordinates which was confusing) it didn't seem to spot the indicator and do anything with it for the S and W directions.

I'm just really putting this on here for completeness in case anyone else uses the script ...

I have found someone here with more knowledge of Javascipt than me (not difficult!!) and there is one line of code we've amended to give a little more flexibility and allow for either a leading or trailing compass letter:

if (/^-|[WS]$/i.test(dmsStr.trim())) deg = -deg; // take '-', west and south as -ve

we have changed to

if (/^-|[WS]|[WS]$/i.test(dmsStr.trim())) deg = -deg; // take '-', west and south as -ve

A small change but means that I can plug it straight into our data without needing to convert to degrees first or try and move the compass letter from the beginning to the end of the coordinates!

Many thanks once again.

Cheers,

Emma