Skip to content

Commit e1beaab

Browse files
jasondaviesjasondavies
committed
d3.geo.bounds: antimeridian and inflection points.
The behaviour of d3.geo.bounds has been modified slightly, such that the minimum and maximum longitudes denote minimum and maximum meridians when going from left to right in a standard equirectangular projection. If the maximum longitude is in fact numerically smaller than the meridian, this simply means that the bounds cross the antimeridian. This fixes several issues: * The bounding box is now a true minimum-width bounding box. * Line segments that cross the antimeridian are now correctly handled. * Similarly, polygons that cross the antimeridian are now correctly handled. Secondly, inflection points are now detected, so that the latitudinal range is extended in the case where the great arc contains the great circle’s inflection point. Lastly, support for polygons has been improved, e.g. those that wind around the poles, and those that have a counter-clockwise winding order. Note that as part of the winding order detection, the performance of d3.geo.area had to be reduced somewhat (by as much as 25%). This is because the trick to avoid atan2 calls fails for winding order detection in the case where a counter-clockwise area strays into a positive quadrant. We can reinstate the trick for area-only calculations in future; I think correct winding order detection is more important at the moment. See d3#1154.
1 parent 0c7d186 commit e1beaab

7 files changed

Lines changed: 486 additions & 111 deletions

File tree

‎d3.js‎

Lines changed: 168 additions & 57 deletions
Original file line numberDiff line numberDiff line change
@@ -2076,7 +2076,7 @@ d3 = function() {
20762076
d3.geo.stream(object, d3_geo_area);
20772077
return d3_geo_areaSum;
20782078
};
2079-
var d3_geo_areaSum, d3_geo_areaRingU, d3_geo_areaRingV;
2079+
var d3_geo_areaSum, d3_geo_areaRingSum;
20802080
var d3_geo_area = {
20812081
sphere: function() {
20822082
d3_geo_areaSum += 4 * π;
@@ -2085,11 +2085,11 @@ d3 = function() {
20852085
lineStart: d3_noop,
20862086
lineEnd: d3_noop,
20872087
polygonStart: function() {
2088-
d3_geo_areaRingU = 1, d3_geo_areaRingV = 0;
2088+
d3_geo_areaRingSum = 0;
20892089
d3_geo_area.lineStart = d3_geo_areaRingStart;
20902090
},
20912091
polygonEnd: function() {
2092-
var area = 2 * Math.atan2(d3_geo_areaRingV, d3_geo_areaRingU);
2092+
var area = 2 * d3_geo_areaRingSum;
20932093
d3_geo_areaSum += area < 0 ? 4 * π + area : area;
20942094
d3_geo_area.lineStart = d3_geo_area.lineEnd = d3_geo_area.point = d3_noop;
20952095
}
@@ -2104,44 +2104,169 @@ d3 = function() {
21042104
function nextPoint(λ, φ) {
21052105
λ *= d3_radians;
21062106
φ = φ * d3_radians / 2 + π / 4;
2107-
var dλ = λ - λ0, cosφ = Math.cos(φ), sinφ = Math.sin(φ), k = sinφ0 * sinφ, u0 = d3_geo_areaRingU, v0 = d3_geo_areaRingV, u = cosφ0 * cosφ + k * Math.cos(dλ), v = k * Math.sin(dλ);
2108-
d3_geo_areaRingU = u0 * u - v0 * v;
2109-
d3_geo_areaRingV = v0 * u + u0 * v;
2107+
var dλ = λ - λ0, cosφ = Math.cos(φ), sinφ = Math.sin(φ), k = sinφ0 * sinφ, u = cosφ0 * cosφ + k * Math.cos(dλ), v = k * Math.sin(dλ);
2108+
d3_geo_areaRingSum += Math.atan2(v, u);
21102109
λ0 = λ, cosφ0 = cosφ, sinφ0 = sinφ;
21112110
}
21122111
d3_geo_area.lineEnd = function() {
21132112
nextPoint(λ00, φ00);
21142113
};
21152114
}
2116-
d3.geo.bounds = d3_geo_bounds(d3_identity);
2117-
function d3_geo_bounds(projectStream) {
2118-
var x0, y0, x1, y1;
2115+
function d3_geo_cartesian(spherical) {
2116+
var λ = spherical[0], φ = spherical[1], cosφ = Math.cos(φ);
2117+
return [ cosφ * Math.cos(λ), cosφ * Math.sin(λ), Math.sin(φ) ];
2118+
}
2119+
function d3_geo_cartesianDot(a, b) {
2120+
return a[0] * b[0] + a[1] * b[1] + a[2] * b[2];
2121+
}
2122+
function d3_geo_cartesianCross(a, b) {
2123+
return [ a[1] * b[2] - a[2] * b[1], a[2] * b[0] - a[0] * b[2], a[0] * b[1] - a[1] * b[0] ];
2124+
}
2125+
function d3_geo_cartesianAdd(a, b) {
2126+
a[0] += b[0];
2127+
a[1] += b[1];
2128+
a[2] += b[2];
2129+
}
2130+
function d3_geo_cartesianScale(vector, k) {
2131+
return [ vector[0] * k, vector[1] * k, vector[2] * k ];
2132+
}
2133+
function d3_geo_cartesianNormalize(d) {
2134+
var l = Math.sqrt(d[0] * d[0] + d[1] * d[1] + d[2] * d[2]);
2135+
d[0] /= l;
2136+
d[1] /= l;
2137+
d[2] /= l;
2138+
}
2139+
function d3_geo_spherical(cartesian) {
2140+
return [ Math.atan2(cartesian[1], cartesian[0]), Math.asin(Math.max(-1, Math.min(1, cartesian[2]))) ];
2141+
}
2142+
function d3_geo_sphericalEqual(a, b) {
2143+
return Math.abs(a[0] - b[0]) < ε && Math.abs(a[1] - b[1]) < ε;
2144+
}
2145+
d3.geo.bounds = function() {
2146+
var λ0, φ0, λ1, φ1, λ_, λ__, φ__, p0, dλSum, ranges, range;
21192147
var bound = {
2120-
point: boundPoint,
2121-
lineStart: d3_noop,
2122-
lineEnd: d3_noop,
2148+
point: point,
2149+
lineStart: lineStart,
2150+
lineEnd: lineEnd,
21232151
polygonStart: function() {
2124-
bound.lineEnd = boundPolygonLineEnd;
2152+
bound.point = ringPoint;
2153+
bound.lineStart = ringStart;
2154+
bound.lineEnd = ringEnd;
2155+
dλSum = 0;
2156+
d3_geo_area.polygonStart();
21252157
},
21262158
polygonEnd: function() {
2127-
bound.point = boundPoint;
2159+
d3_geo_area.polygonEnd();
2160+
bound.point = point;
2161+
bound.lineStart = lineStart;
2162+
bound.lineEnd = lineEnd;
2163+
if (d3_geo_areaRingSum < 0) λ0 = -(λ1 = 180), φ0 = -(φ1 = 90); else if (dλSum > ε) φ1 = 90; else if (dλSum < -ε) φ0 = -90;
2164+
range[0] = λ0, range[1] = λ1;
2165+
}
2166+
};
2167+
function point(λ, φ) {
2168+
ranges.push(range = [ λ0 = λ, λ1 = λ ]);
2169+
if (φ < φ0) φ0 = φ;
2170+
if (φ > φ1) φ1 = φ;
2171+
}
2172+
function linePoint(λ, φ) {
2173+
var p = d3_geo_cartesian([ λ * d3_radians, φ * d3_radians ]);
2174+
if (p0) {
2175+
var normal = d3_geo_cartesianCross(p0, p), equatorial = [ normal[1], -normal[0], 0 ], inflection = d3_geo_cartesianCross(equatorial, normal);
2176+
d3_geo_cartesianNormalize(inflection);
2177+
inflection = d3_geo_spherical(inflection);
2178+
var dλ = λ - λ_, s = dλ > 0 ? 1 : -1, λi = inflection[0] * d3_degrees * s, antimeridian = Math.abs(dλ) > 180;
2179+
if (antimeridian ^ (s * λ_ < λi && λi < s * λ)) {
2180+
var φi = inflection[1] * d3_degrees;
2181+
if (φi > φ1) φ1 = φi;
2182+
} else if (λi = (λi + 360) % 360 - 180, antimeridian ^ (s * λ_ < λi && λi < s * λ)) {
2183+
var φi = -inflection[1] * d3_degrees;
2184+
if (φi < φ0) φ0 = φi;
2185+
} else {
2186+
if (φ < φ0) φ0 = φ;
2187+
if (φ > φ1) φ1 = φ;
2188+
}
2189+
if (antimeridian) {
2190+
if (λ < λ_) {
2191+
if (angle(λ0, λ) > angle(λ0, λ1)) λ1 = λ;
2192+
} else {
2193+
if (angle(λ, λ1) > angle(λ0, λ1)) λ0 = λ;
2194+
}
2195+
} else {
2196+
if (λ1 >= λ0) {
2197+
if (λ < λ0) λ0 = λ;
2198+
if (λ > λ1) λ1 = λ;
2199+
} else {
2200+
if (λ > λ_) {
2201+
if (angle(λ0, λ) > angle(λ0, λ1)) λ1 = λ;
2202+
} else {
2203+
if (angle(λ, λ1) > angle(λ0, λ1)) λ0 = λ;
2204+
}
2205+
}
2206+
}
2207+
} else {
2208+
point(λ, φ);
21282209
}
2129-
};
2130-
function boundPoint(x, y) {
2131-
if (x < x0) x0 = x;
2132-
if (x > x1) x1 = x;
2133-
if (y < y0) y0 = y;
2134-
if (y > y1) y1 = y;
2210+
p0 = p, λ_ = λ;
2211+
}
2212+
function lineStart() {
2213+
bound.point = linePoint;
2214+
}
2215+
function lineEnd() {
2216+
range[0] = λ0, range[1] = λ1;
2217+
bound.point = point;
2218+
p0 = null;
2219+
}
2220+
function ringPoint(λ, φ) {
2221+
if (p0) {
2222+
var dλ = λ - λ_;
2223+
dλSum += Math.abs(dλ) > 180 ? dλ + (dλ > 0 ? 360 : -360) : dλ;
2224+
} else λ__ = λ, φ__ = φ;
2225+
d3_geo_area.point(λ, φ);
2226+
linePoint(λ, φ);
2227+
}
2228+
function ringStart() {
2229+
d3_geo_area.lineStart();
2230+
}
2231+
function ringEnd() {
2232+
ringPoint(λ__, φ__);
2233+
d3_geo_area.lineEnd();
2234+
if (Math.abs(dλSum) > ε) λ0 = -(λ1 = 180);
2235+
range[0] = λ0, range[1] = λ1;
2236+
p0 = null;
2237+
}
2238+
function angle(λ0, λ1) {
2239+
return (λ1 -= λ0) < 0 ? λ1 + 360 : λ1;
2240+
}
2241+
function compareRanges(a, b) {
2242+
return a[0] - b[0];
21352243
}
2136-
function boundPolygonLineEnd() {
2137-
bound.point = bound.lineEnd = d3_noop;
2244+
function withinRange(x, range) {
2245+
return range[0] <= range[1] ? range[0] <= x && x <= range[1] : x < range[0] || range[1] < x;
21382246
}
21392247
return function(feature) {
2140-
y1 = x1 = -(x0 = y0 = Infinity);
2141-
d3.geo.stream(feature, projectStream(bound));
2142-
return [ [ x0, y0 ], [ x1, y1 ] ];
2248+
φ1 = λ1 = -(λ0 = φ0 = Infinity);
2249+
ranges = [];
2250+
d3.geo.stream(feature, bound);
2251+
ranges.sort(compareRanges);
2252+
for (var i = 1, n = ranges.length, a = ranges[0], b, merged = [ a ]; i < n; ++i) {
2253+
b = ranges[i];
2254+
if (withinRange(b[0], a) || withinRange(b[1], a)) {
2255+
if (angle(a[0], b[1]) > angle(a[0], a[1])) a[1] = b[1];
2256+
if (angle(b[0], a[1]) > angle(a[0], a[1])) a[0] = b[0];
2257+
} else {
2258+
merged.push(a = b);
2259+
}
2260+
}
2261+
var best = -Infinity, dλ;
2262+
for (var n = merged.length - 1, i = 0, a = merged[n], b; i <= n; a = b, ++i) {
2263+
b = merged[i];
2264+
if ((dλ = angle(a[1], b[0])) > best) best = dλ, λ0 = b[0], λ1 = a[1];
2265+
}
2266+
ranges = range = null;
2267+
return [ [ λ0, φ0 ], [ λ1, φ1 ] ];
21432268
};
2144-
}
2269+
}();
21452270
d3.geo.centroid = function(object) {
21462271
d3_geo_centroidDimension = d3_geo_centroidW = d3_geo_centroidX = d3_geo_centroidY = d3_geo_centroidZ = 0;
21472272
d3.geo.stream(object, d3_geo_centroid);
@@ -2223,39 +2348,9 @@ d3 = function() {
22232348
function d3_geo_centroidLineEnd() {
22242349
d3_geo_centroid.point = d3_geo_centroidPoint;
22252350
}
2226-
function d3_geo_cartesian(spherical) {
2227-
var λ = spherical[0], φ = spherical[1], cosφ = Math.cos(φ);
2228-
return [ cosφ * Math.cos(λ), cosφ * Math.sin(λ), Math.sin(φ) ];
2229-
}
2230-
function d3_geo_cartesianDot(a, b) {
2231-
return a[0] * b[0] + a[1] * b[1] + a[2] * b[2];
2232-
}
2233-
function d3_geo_cartesianCross(a, b) {
2234-
return [ a[1] * b[2] - a[2] * b[1], a[2] * b[0] - a[0] * b[2], a[0] * b[1] - a[1] * b[0] ];
2235-
}
2236-
function d3_geo_cartesianAdd(a, b) {
2237-
a[0] += b[0];
2238-
a[1] += b[1];
2239-
a[2] += b[2];
2240-
}
2241-
function d3_geo_cartesianScale(vector, k) {
2242-
return [ vector[0] * k, vector[1] * k, vector[2] * k ];
2243-
}
2244-
function d3_geo_cartesianNormalize(d) {
2245-
var l = Math.sqrt(d[0] * d[0] + d[1] * d[1] + d[2] * d[2]);
2246-
d[0] /= l;
2247-
d[1] /= l;
2248-
d[2] /= l;
2249-
}
22502351
function d3_true() {
22512352
return true;
22522353
}
2253-
function d3_geo_spherical(cartesian) {
2254-
return [ Math.atan2(cartesian[1], cartesian[0]), Math.asin(Math.max(-1, Math.min(1, cartesian[2]))) ];
2255-
}
2256-
function d3_geo_sphericalEqual(a, b) {
2257-
return Math.abs(a[0] - b[0]) < ε && Math.abs(a[1] - b[1]) < ε;
2258-
}
22592354
function d3_geo_clipPolygon(segments, compare, inside, interpolate, listener) {
22602355
var subject = [], clip = [];
22612356
segments.forEach(function(segment) {
@@ -3341,6 +3436,20 @@ d3 = function() {
33413436
nextPoint(x00, y00);
33423437
};
33433438
}
3439+
var d3_geo_pathBoundsX0, d3_geo_pathBoundsY0, d3_geo_pathBoundsX1, d3_geo_pathBoundsY1;
3440+
var d3_geo_pathBounds = {
3441+
point: d3_geo_pathBoundsPoint,
3442+
lineStart: d3_noop,
3443+
lineEnd: d3_noop,
3444+
polygonStart: d3_noop,
3445+
polygonEnd: d3_noop
3446+
};
3447+
function d3_geo_pathBoundsPoint(x, y) {
3448+
if (x < d3_geo_pathBoundsX0) d3_geo_pathBoundsX0 = x;
3449+
if (x > d3_geo_pathBoundsX1) d3_geo_pathBoundsX1 = x;
3450+
if (y < d3_geo_pathBoundsY0) d3_geo_pathBoundsY0 = y;
3451+
if (y > d3_geo_pathBoundsY1) d3_geo_pathBoundsY1 = y;
3452+
}
33443453
function d3_geo_pathBuffer() {
33453454
var pointCircle = d3_geo_pathCircle(4.5), buffer = [];
33463455
var stream = {
@@ -3506,7 +3615,9 @@ d3 = function() {
35063615
return d3_geo_centroidZ ? [ d3_geo_centroidX / d3_geo_centroidZ, d3_geo_centroidY / d3_geo_centroidZ ] : undefined;
35073616
};
35083617
path.bounds = function(object) {
3509-
return d3_geo_bounds(projectStream)(object);
3618+
d3_geo_pathBoundsX1 = d3_geo_pathBoundsY1 = -(d3_geo_pathBoundsX0 = d3_geo_pathBoundsY0 = Infinity);
3619+
d3.geo.stream(object, projectStream(d3_geo_pathBounds));
3620+
return [ [ d3_geo_pathBoundsX0, d3_geo_pathBoundsY0 ], [ d3_geo_pathBoundsX1, d3_geo_pathBoundsY1 ] ];
35103621
};
35113622
path.projection = function(_) {
35123623
if (!arguments.length) return projection;

‎d3.min.js‎

Lines changed: 5 additions & 5 deletions
Some generated files are not rendered by default. Learn more about customizing how changed files appear on GitHub.

‎src/geo/area.js‎

Lines changed: 4 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -10,8 +10,7 @@ d3.geo.area = function(object) {
1010
};
1111

1212
var d3_geo_areaSum,
13-
d3_geo_areaRingU,
14-
d3_geo_areaRingV;
13+
d3_geo_areaRingSum;
1514

1615
var d3_geo_area = {
1716
sphere: function() { d3_geo_areaSum += 4 * π; },
@@ -21,11 +20,11 @@ var d3_geo_area = {
2120

2221
// Only count area for polygon rings.
2322
polygonStart: function() {
24-
d3_geo_areaRingU = 1, d3_geo_areaRingV = 0;
23+
d3_geo_areaRingSum = 0;
2524
d3_geo_area.lineStart = d3_geo_areaRingStart;
2625
},
2726
polygonEnd: function() {
28-
var area = 2 * Math.atan2(d3_geo_areaRingV, d3_geo_areaRingU);
27+
var area = 2 * d3_geo_areaRingSum;
2928
d3_geo_areaSum += area < 0 ? 4 * π + area : area;
3029
d3_geo_area.lineStart = d3_geo_area.lineEnd = d3_geo_area.point = d3_noop;
3130
}
@@ -52,13 +51,9 @@ function d3_geo_areaRingStart() {
5251
cosφ = Math.cos(φ),
5352
sinφ = Math.sin(φ),
5453
k = sinφ0 * sinφ,
55-
u0 = d3_geo_areaRingU,
56-
v0 = d3_geo_areaRingV,
5754
u = cosφ0 * cosφ + k * Math.cos(dλ),
5855
v = k * Math.sin(dλ);
59-
// ∑ arg(z) = arg(∏ z), where z = u + iv.
60-
d3_geo_areaRingU = u0 * u - v0 * v;
61-
d3_geo_areaRingV = v0 * u + u0 * v;
56+
d3_geo_areaRingSum += Math.atan2(v, u);
6257

6358
// Advance the previous points.
6459
λ0 = λ, cosφ0 = cosφ, sinφ0 = sinφ;

0 commit comments

Comments
 (0)