From e8d637f2bb511b575d91ee8bf6328683f278ff0c Mon Sep 17 00:00:00 2001 From: tyr Date: Thu, 12 Sep 2013 13:55:07 +0200 Subject: [PATCH] Add iD.geo.sphericalDistance iD.geo.euclideanDistance should only be used for calculations of projected coordinates or display (pixel) coordinates. iD.geo.sphericalDistance calculates approximate geographical distances, accounting for distortions at higher latitudes. This can be used for determining the nearest node (operations.Delete, actions.Circularize) or relative length comparisons (actions.Split). --- js/id/actions/circularize.js | 4 ++-- js/id/actions/orthogonalize.js | 2 +- js/id/actions/split.js | 2 +- js/id/behavior/draw.js | 4 ++-- js/id/geo.js | 10 +++++++-- js/id/operations/delete.js | 4 ++-- js/id/svg.js | 2 +- js/id/svg/midpoints.js | 2 +- test/spec/actions/circularize.js | 6 ++--- test/spec/actions/orthogonalize.js | 4 ++-- test/spec/geo.js | 36 ++++++++++++++++++++++++++---- 11 files changed, 55 insertions(+), 21 deletions(-) diff --git a/js/id/actions/circularize.js b/js/id/actions/circularize.js index d1c973723..10a8f0f7c 100644 --- a/js/id/actions/circularize.js +++ b/js/id/actions/circularize.js @@ -8,7 +8,7 @@ iD.actions.Circularize = function(wayId, projection, maxAngle) { points = nodes.map(function(n) { return projection(n.loc); }), keyPoints = keyNodes.map(function(n) { return projection(n.loc); }), centroid = d3.geom.polygon(points).centroid(), - radius = d3.median(points, function(p) { return iD.geo.dist(centroid, p); }), + radius = d3.median(points, function(p) { return iD.geo.euclideanDistance(centroid, p); }), sign = d3.geom.polygon(points).area() > 0 ? 1 : -1, ids; @@ -44,7 +44,7 @@ iD.actions.Circularize = function(wayId, projection, maxAngle) { } // position this key node - distance = iD.geo.dist(centroid, keyPoints[i]); + distance = iD.geo.euclideanDistance(centroid, keyPoints[i]); keyPoints[i] = [ centroid[0] + (keyPoints[i][0] - centroid[0]) / distance * radius, centroid[1] + (keyPoints[i][1] - centroid[1]) / distance * radius]; diff --git a/js/id/actions/orthogonalize.js b/js/id/actions/orthogonalize.js index 6789768c8..d21cd5372 100644 --- a/js/id/actions/orthogonalize.js +++ b/js/id/actions/orthogonalize.js @@ -84,7 +84,7 @@ iD.actions.Orthogonalize = function(wayId, projection) { p = subtractPoints(a, b), q = subtractPoints(c, b); - var scale = 2*Math.min(iD.geo.dist(p, [0, 0]), iD.geo.dist(q, [0, 0])); + var scale = 2*Math.min(iD.geo.euclideanDistance(p, [0, 0]), iD.geo.euclideanDistance(q, [0, 0])); p = normalizePoint(p, 1.0); q = normalizePoint(q, 1.0); diff --git a/js/id/actions/split.js b/js/id/actions/split.js index 8875c1e1d..be1a9153f 100644 --- a/js/id/actions/split.js +++ b/js/id/actions/split.js @@ -44,7 +44,7 @@ iD.actions.Split = function(nodeId, newWayIds) { return iD.util.wrap(index,nodes.length); } function _dist(nA, nB) { - return iD.geo.dist(graph.entity(nA).loc, graph.entity(nB).loc); + return iD.geo.sphericalDistance(graph.entity(nA).loc, graph.entity(nB).loc); } // calculate lengths diff --git a/js/id/behavior/draw.js b/js/id/behavior/draw.js index e202030a2..d305d944f 100644 --- a/js/id/behavior/draw.js +++ b/js/id/behavior/draw.js @@ -34,8 +34,8 @@ iD.behavior.Draw = function(context) { d3.select(window).on('mouseup.draw', function() { element.on('mousemove.draw', mousemove); - if (iD.geo.dist(pos, point()) < closeTolerance || - (iD.geo.dist(pos, point()) < tolerance && + if (iD.geo.euclideanDistance(pos, point()) < closeTolerance || + (iD.geo.euclideanDistance(pos, point()) < tolerance && (+new Date() - time) < 500)) { // Prevent a quick second click diff --git a/js/id/geo.js b/js/id/geo.js index ff6b060c8..3b901d263 100644 --- a/js/id/geo.js +++ b/js/id/geo.js @@ -10,17 +10,23 @@ iD.geo.interp = function(p1, p2, t) { }; // http://jsperf.com/id-dist-optimization -iD.geo.dist = function(a, b) { +iD.geo.euclideanDistance = function(a, b) { var x = a[0] - b[0], y = a[1] - b[1]; return Math.sqrt((x * x) + (y * y)); }; +// Equirectangular approximation of spherical distances on Earth +iD.geo.sphericalDistance = function(a, b) { + var x = Math.cos(a[1]*Math.PI/180) * (a[0] - b[0]), + y = a[1] - b[1]; + return 6.3710E6 * Math.sqrt((x * x) + (y * y)) * Math.PI/180; +}; // Choose the edge with the minimal distance from `point` to its orthogonal // projection onto that edge, if such a projection exists, or the distance to // the closest vertex on that edge. Returns an object with the `index` of the // chosen edge, the chosen `loc` on that edge, and the `distance` to to it. iD.geo.chooseEdge = function(nodes, point, projection) { - var dist = iD.geo.dist, + var dist = iD.geo.euclideanDistance, points = nodes.map(function(n) { return projection(n.loc); }), min = Infinity, idx, loc; diff --git a/js/id/operations/delete.js b/js/id/operations/delete.js index c4863cdbd..d81d84a0d 100644 --- a/js/id/operations/delete.js +++ b/js/id/operations/delete.js @@ -27,8 +27,8 @@ iD.operations.Delete = function(selectedIDs, context) { } else if (i === nodes.length - 1) { i--; } else { - var a = iD.geo.dist(entity.loc, context.entity(nodes[i - 1]).loc), - b = iD.geo.dist(entity.loc, context.entity(nodes[i + 1]).loc); + var a = iD.geo.sphericalDistance(entity.loc, context.entity(nodes[i - 1]).loc), + b = iD.geo.sphericalDistance(entity.loc, context.entity(nodes[i + 1]).loc); i = a < b ? i - 1 : i + 1; } diff --git a/js/id/svg.js b/js/id/svg.js index a74b56986..00e6a8c2b 100644 --- a/js/id/svg.js +++ b/js/id/svg.js @@ -63,7 +63,7 @@ iD.svg = { b = [x, y]; if (a) { - var span = iD.geo.dist(a, b) - offset; + var span = iD.geo.euclideanDistance(a, b) - offset; if (span >= 0) { var angle = Math.atan2(b[1] - a[1], b[0] - a[0]), diff --git a/js/id/svg/midpoints.js b/js/id/svg/midpoints.js index 892efd141..55a069172 100644 --- a/js/id/svg/midpoints.js +++ b/js/id/svg/midpoints.js @@ -20,7 +20,7 @@ iD.svg.Midpoints = function(projection, context) { // If neither of the nodes changed, no need to redraw midpoint if (!midpoints[id] && (filter(a) || filter(b))) { var loc = iD.geo.interp(a.loc, b.loc, 0.5); - if (extent.intersects(loc) && iD.geo.dist(projection(a.loc), projection(b.loc)) > 40) { + if (extent.intersects(loc) && iD.geo.euclideanDistance(projection(a.loc), projection(b.loc)) > 40) { midpoints[id] = { type: 'midpoint', id: id, diff --git a/test/spec/actions/circularize.js b/test/spec/actions/circularize.js index 9456f109e..282310796 100644 --- a/test/spec/actions/circularize.js +++ b/test/spec/actions/circularize.js @@ -48,7 +48,7 @@ describe("iD.actions.Circularize", function () { graph = iD.actions.Circularize('-', projection)(graph); - expect(iD.geo.dist(graph.entity('d').loc, [2, -2])).to.be.lt(0.5); + expect(iD.geo.euclideanDistance(graph.entity('d').loc, [2, -2])).to.be.lt(0.5); }); function angle(point1, point2, center) { @@ -56,10 +56,10 @@ describe("iD.actions.Circularize", function () { vector2 = [point2[0] - center[0], point2[1] - center[1]], distance; - distance = iD.geo.dist(vector1, [0, 0]); + distance = iD.geo.euclideanDistance(vector1, [0, 0]); vector1 = [vector1[0] / distance, vector1[1] / distance]; - distance = iD.geo.dist(vector2, [0, 0]); + distance = iD.geo.euclideanDistance(vector2, [0, 0]); vector2 = [vector2[0] / distance, vector2[1] / distance]; return 180 / Math.PI * Math.acos(vector1[0] * vector2[0] + vector1[1] * vector2[1]); diff --git a/test/spec/actions/orthogonalize.js b/test/spec/actions/orthogonalize.js index 396403749..ed889bda0 100644 --- a/test/spec/actions/orthogonalize.js +++ b/test/spec/actions/orthogonalize.js @@ -97,12 +97,12 @@ describe("iD.actions.Orthogonalize", function () { 'd': iD.Node({id: 'd', loc: tests[i][3]}), '-': iD.Way({id: '-', nodes: ['a', 'b', 'c', 'd', 'a']}) }), - initialWidth = iD.geo.dist(graph.entity('a').loc, graph.entity('b').loc), + initialWidth = iD.geo.sphericalDistance(graph.entity('a').loc, graph.entity('b').loc), finalWidth; graph = iD.actions.Orthogonalize('-', projection)(graph); - finalWidth = iD.geo.dist(graph.entity('a').loc, graph.entity('b').loc); + finalWidth = iD.geo.sphericalDistance(graph.entity('a').loc, graph.entity('b').loc); expect(finalWidth / initialWidth).within(0.90, 1.10); } }); diff --git a/test/spec/geo.js b/test/spec/geo.js index ba2f401a9..7592782b7 100644 --- a/test/spec/geo.js +++ b/test/spec/geo.js @@ -18,21 +18,49 @@ describe('iD.geo', function() { }); }); - describe('.dist', function() { + describe('.euclideanDistance', function() { it('distance between two same points is zero', function() { var a = [0, 0], b = [0, 0]; - expect(iD.geo.dist(a, b)).to.eql(0); + expect(iD.geo.euclideanDistance(a, b)).to.eql(0); }); it('a straight 10 unit line is 10', function() { var a = [0, 0], b = [10, 0]; - expect(iD.geo.dist(a, b)).to.eql(10); + expect(iD.geo.euclideanDistance(a, b)).to.eql(10); }); it('a pythagorean triangle is right', function() { var a = [0, 0], b = [4, 3]; - expect(iD.geo.dist(a, b)).to.eql(5); + expect(iD.geo.euclideanDistance(a, b)).to.eql(5); + }); + }); + + describe('.sphericalDistance', function() { + it('distance between two same points is zero', function() { + var a = [0, 0], + b = [0, 0]; + expect(iD.geo.sphericalDistance(a, b)).to.eql(0); + }); + it('a straight 1 degree line at the equator is aproximately 111 km', function() { + var a = [0, 0], + b = [1, 0]; + expect(iD.geo.sphericalDistance(a, b)).to.be.within(100E3,120E3); + }); + it('a pythagorean triangle is right', function() { + var a = [0, 0], + b = [4, 3]; + expect(iD.geo.sphericalDistance(a, b)).to.be.within(500E3,600E3); + }); + it('east-west distances at high latitude are shorter', function() { + var a = [0, 60], + b = [1, 60]; + expect(iD.geo.sphericalDistance(a, b)).to.be.within(50E3,60E3); + }); + it('north-south distances at high latitude are not shorter', function() { + var a = [0, 60], + b = [0, 61]; + expect(iD.geo.sphericalDistance(a, b)).to.be.within(100E3,120E3); }); });