Skip to content

Instantly share code, notes, and snippets.

@Fil
Last active Aug 19, 2017
Embed
What would you like to do?
L’Équateur de Marseille
license: gpl-3.0
height: 600
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.

L'équateur de Marseille s'appuie sur la Canebière et la prolonge dans les deux directions, jusqu'à faire le tour du monde.

Cette approche a été conçue et développée par Jean-Luc Arnaud dans « Marseille et la Canebière - Echelles du monde », 2017. <halshs-01484193>

Visualisation interactive par Philippe Rivière -- visionscarto.net


Forked from mbostock's block: Versor Dragging II

// https://d3js.org/d3-geo-projection/ Version 2.1.0. Copyright 2017 Mike Bostock.
(function (global, factory) {
typeof exports === 'object' && typeof module !== 'undefined' ? factory(exports, require('d3-geo'), require('d3-array')) :
typeof define === 'function' && define.amd ? define(['exports', 'd3-geo', 'd3-array'], factory) :
(factory((global.d3 = global.d3 || {}),global.d3,global.d3));
}(this, (function (exports,d3Geo,d3Array) { 'use strict';
var abs = Math.abs;
var atan = Math.atan;
var atan2 = Math.atan2;
var cos = Math.cos;
var exp = Math.exp;
var floor = Math.floor;
var log = Math.log;
var max = Math.max;
var min = Math.min;
var pow = Math.pow;
var round = Math.round;
var sign = Math.sign || function(x) { return x > 0 ? 1 : x < 0 ? -1 : 0; };
var sin = Math.sin;
var tan = Math.tan;
var epsilon = 1e-6;
var epsilon2 = 1e-12;
var pi = Math.PI;
var halfPi = pi / 2;
var quarterPi = pi / 4;
var sqrt1_2 = Math.SQRT1_2;
var sqrt2 = sqrt(2);
var sqrtPi = sqrt(pi);
var tau = pi * 2;
var degrees = 180 / pi;
var radians = pi / 180;
function sinci(x) {
return x ? x / Math.sin(x) : 1;
}
function asin(x) {
return x > 1 ? halfPi : x < -1 ? -halfPi : Math.asin(x);
}
function acos(x) {
return x > 1 ? 0 : x < -1 ? pi : Math.acos(x);
}
function sqrt(x) {
return x > 0 ? Math.sqrt(x) : 0;
}
function tanh(x) {
x = exp(2 * x);
return (x - 1) / (x + 1);
}
function sinh(x) {
return (exp(x) - exp(-x)) / 2;
}
function cosh(x) {
return (exp(x) + exp(-x)) / 2;
}
function arsinh(x) {
return log(x + sqrt(x * x + 1));
}
function arcosh(x) {
return log(x + sqrt(x * x - 1));
}
function airyRaw(beta) {
var tanBeta_2 = tan(beta / 2),
b = 2 * log(cos(beta / 2)) / (tanBeta_2 * tanBeta_2);
function forward(x, y) {
var cosx = cos(x),
cosy = cos(y),
siny = sin(y),
cosz = cosy * cosx,
k = -((1 - cosz ? log((1 + cosz) / 2) / (1 - cosz) : -0.5) + b / (1 + cosz));
return [k * cosy * sin(x), k * siny];
}
forward.invert = function(x, y) {
var r = sqrt(x * x + y * y),
z = -beta / 2,
i = 50, delta;
if (!r) return [0, 0];
do {
var z_2 = z / 2,
cosz_2 = cos(z_2),
sinz_2 = sin(z_2),
tanz_2 = tan(z_2),
lnsecz_2 = log(1 / cosz_2);
z -= delta = (2 / tanz_2 * lnsecz_2 - b * tanz_2 - r) / (-lnsecz_2 / (sinz_2 * sinz_2) + 1 - b / (2 * cosz_2 * cosz_2));
} while (abs(delta) > epsilon && --i > 0);
var sinz = sin(z);
return [atan2(x * sinz, r * cos(z)), asin(y * sinz / r)];
};
return forward;
}
var airy = function() {
var beta = halfPi,
m = d3Geo.geoProjectionMutator(airyRaw),
p = m(beta);
p.radius = function(_) {
return arguments.length ? m(beta = _ * radians) : beta * degrees;
};
return p
.scale(179.976)
.clipAngle(147);
};
function aitoffRaw(x, y) {
var cosy = cos(y), sincia = sinci(acos(cosy * cos(x /= 2)));
return [2 * cosy * sin(x) * sincia, sin(y) * sincia];
}
// Abort if [x, y] is not within an ellipse centered at [0, 0] with
// semi-major axis pi and semi-minor axis pi/2.
aitoffRaw.invert = function(x, y) {
if (x * x + 4 * y * y > pi * pi + epsilon) return;
var x1 = x, y1 = y, i = 25;
do {
var sinx = sin(x1),
sinx_2 = sin(x1 / 2),
cosx_2 = cos(x1 / 2),
siny = sin(y1),
cosy = cos(y1),
sin_2y = sin(2 * y1),
sin2y = siny * siny,
cos2y = cosy * cosy,
sin2x_2 = sinx_2 * sinx_2,
c = 1 - cos2y * cosx_2 * cosx_2,
e = c ? acos(cosy * cosx_2) * sqrt(f = 1 / c) : f = 0,
f,
fx = 2 * e * cosy * sinx_2 - x,
fy = e * siny - y,
dxdx = f * (cos2y * sin2x_2 + e * cosy * cosx_2 * sin2y),
dxdy = f * (0.5 * sinx * sin_2y - e * 2 * siny * sinx_2),
dydx = f * 0.25 * (sin_2y * sinx_2 - e * siny * cos2y * sinx),
dydy = f * (sin2y * cosx_2 + e * sin2x_2 * cosy),
z = dxdy * dydx - dydy * dxdx;
if (!z) break;
var dx = (fy * dxdy - fx * dydy) / z,
dy = (fx * dydx - fy * dxdx) / z;
x1 -= dx, y1 -= dy;
} while ((abs(dx) > epsilon || abs(dy) > epsilon) && --i > 0);
return [x1, y1];
};
var aitoff = function() {
return d3Geo.geoProjection(aitoffRaw)
.scale(152.63);
};
function armadilloRaw(phi0) {
var sinPhi0 = sin(phi0),
cosPhi0 = cos(phi0),
sPhi0 = phi0 >= 0 ? 1 : -1,
tanPhi0 = tan(sPhi0 * phi0),
k = (1 + sinPhi0 - cosPhi0) / 2;
function forward(lambda, phi) {
var cosPhi = cos(phi),
cosLambda = cos(lambda /= 2);
return [
(1 + cosPhi) * sin(lambda),
(sPhi0 * phi > -atan2(cosLambda, tanPhi0) - 1e-3 ? 0 : -sPhi0 * 10) + k + sin(phi) * cosPhi0 - (1 + cosPhi) * sinPhi0 * cosLambda // TODO D3 core should allow null or [NaN, NaN] to be returned.
];
}
forward.invert = function(x, y) {
var lambda = 0,
phi = 0,
i = 50;
do {
var cosLambda = cos(lambda),
sinLambda = sin(lambda),
cosPhi = cos(phi),
sinPhi = sin(phi),
A = 1 + cosPhi,
fx = A * sinLambda - x,
fy = k + sinPhi * cosPhi0 - A * sinPhi0 * cosLambda - y,
dxdLambda = A * cosLambda / 2,
dxdPhi = -sinLambda * sinPhi,
dydLambda = sinPhi0 * A * sinLambda / 2,
dydPhi = cosPhi0 * cosPhi + sinPhi0 * cosLambda * sinPhi,
denominator = dxdPhi * dydLambda - dydPhi * dxdLambda,
dLambda = (fy * dxdPhi - fx * dydPhi) / denominator / 2,
dPhi = (fx * dydLambda - fy * dxdLambda) / denominator;
lambda -= dLambda, phi -= dPhi;
} while ((abs(dLambda) > epsilon || abs(dPhi) > epsilon) && --i > 0);
return sPhi0 * phi > -atan2(cos(lambda), tanPhi0) - 1e-3 ? [lambda * 2, phi] : null;
};
return forward;
}
var armadillo = function() {
var phi0 = 20 * radians,
sPhi0 = phi0 >= 0 ? 1 : -1,
tanPhi0 = tan(sPhi0 * phi0),
m = d3Geo.geoProjectionMutator(armadilloRaw),
p = m(phi0),
stream_ = p.stream;
p.parallel = function(_) {
if (!arguments.length) return phi0 * degrees;
tanPhi0 = tan((sPhi0 = (phi0 = _ * radians) >= 0 ? 1 : -1) * phi0);
return m(phi0);
};
p.stream = function(stream) {
var rotate = p.rotate(),
rotateStream = stream_(stream),
sphereStream = (p.rotate([0, 0]), stream_(stream));
p.rotate(rotate);
rotateStream.sphere = function() {
sphereStream.polygonStart(), sphereStream.lineStart();
for (var lambda = sPhi0 * -180; sPhi0 * lambda < 180; lambda += sPhi0 * 90) sphereStream.point(lambda, sPhi0 * 90);
while (sPhi0 * (lambda -= phi0) >= -180) { // TODO precision?
sphereStream.point(lambda, sPhi0 * -atan2(cos(lambda * radians / 2), tanPhi0) * degrees);
}
sphereStream.lineEnd(), sphereStream.polygonEnd();
};
return rotateStream;
};
return p
.scale(218.695)
.center([0, 28.0974]);
};
function augustRaw(lambda, phi) {
var tanPhi = tan(phi / 2),
k = sqrt(1 - tanPhi * tanPhi),
c = 1 + k * cos(lambda /= 2),
x = sin(lambda) * k / c,
y = tanPhi / c,
x2 = x * x,
y2 = y * y;
return [
4 / 3 * x * (3 + x2 - 3 * y2),
4 / 3 * y * (3 + 3 * x2 - y2)
];
}
augustRaw.invert = function(x, y) {
x *= 3 / 8, y *= 3 / 8;
if (!x && abs(y) > 1) return null;
var x2 = x * x,
y2 = y * y,
s = 1 + x2 + y2,
sin3Eta = sqrt((s - sqrt(s * s - 4 * y * y)) / 2),
eta = asin(sin3Eta) / 3,
xi = sin3Eta ? arcosh(abs(y / sin3Eta)) / 3 : arsinh(abs(x)) / 3,
cosEta = cos(eta),
coshXi = cosh(xi),
d = coshXi * coshXi - cosEta * cosEta;
return [
sign(x) * 2 * atan2(sinh(xi) * cosEta, 0.25 - d),
sign(y) * 2 * atan2(coshXi * sin(eta), 0.25 + d)
];
};
var august = function() {
return d3Geo.geoProjection(augustRaw)
.scale(66.1603);
};
var sqrt8 = sqrt(8);
var phi0 = log(1 + sqrt2);
function bakerRaw(lambda, phi) {
var phi0 = abs(phi);
return phi0 < quarterPi
? [lambda, log(tan(quarterPi + phi / 2))]
: [lambda * cos(phi0) * (2 * sqrt2 - 1 / sin(phi0)), sign(phi) * (2 * sqrt2 * (phi0 - quarterPi) - log(tan(phi0 / 2)))];
}
bakerRaw.invert = function(x, y) {
if ((y0 = abs(y)) < phi0) return [x, 2 * atan(exp(y)) - halfPi];
var phi = quarterPi, i = 25, delta, y0;
do {
var cosPhi_2 = cos(phi / 2), tanPhi_2 = tan(phi / 2);
phi -= delta = (sqrt8 * (phi - quarterPi) - log(tanPhi_2) - y0) / (sqrt8 - cosPhi_2 * cosPhi_2 / (2 * tanPhi_2));
} while (abs(delta) > epsilon2 && --i > 0);
return [x / (cos(phi) * (sqrt8 - 1 / sin(phi))), sign(y) * phi];
};
var baker = function() {
return d3Geo.geoProjection(bakerRaw)
.scale(112.314);
};
function berghausRaw(lobes) {
var k = 2 * pi / lobes;
function forward(lambda, phi) {
var p = d3Geo.geoAzimuthalEquidistantRaw(lambda, phi);
if (abs(lambda) > halfPi) { // back hemisphere
var theta = atan2(p[1], p[0]),
r = sqrt(p[0] * p[0] + p[1] * p[1]),
theta0 = k * round((theta - halfPi) / k) + halfPi,
alpha = atan2(sin(theta -= theta0), 2 - cos(theta)); // angle relative to lobe end
theta = theta0 + asin(pi / r * sin(alpha)) - alpha;
p[0] = r * cos(theta);
p[1] = r * sin(theta);
}
return p;
}
forward.invert = function(x, y) {
var r = sqrt(x * x + y * y);
if (r > halfPi) {
var theta = atan2(y, x),
theta0 = k * round((theta - halfPi) / k) + halfPi,
s = theta > theta0 ? -1 : 1,
A = r * cos(theta0 - theta),
cotAlpha = 1 / tan(s * acos((A - pi) / sqrt(pi * (pi - 2 * A) + r * r)));
theta = theta0 + 2 * atan((cotAlpha + s * sqrt(cotAlpha * cotAlpha - 3)) / 3);
x = r * cos(theta), y = r * sin(theta);
}
return d3Geo.geoAzimuthalEquidistantRaw.invert(x, y);
};
return forward;
}
var berghaus = function() {
var lobes = 5,
m = d3Geo.geoProjectionMutator(berghausRaw),
p = m(lobes),
projectionStream = p.stream,
epsilon$$1 = 1e-2,
cr = -cos(epsilon$$1 * radians),
sr = sin(epsilon$$1 * radians);
p.lobes = function(_) {
return arguments.length ? m(lobes = +_) : lobes;
};
p.stream = function(stream) {
var rotate = p.rotate(),
rotateStream = projectionStream(stream),
sphereStream = (p.rotate([0, 0]), projectionStream(stream));
p.rotate(rotate);
rotateStream.sphere = function() {
sphereStream.polygonStart(), sphereStream.lineStart();
for (var i = 0, delta = 360 / lobes, delta0 = 2 * pi / lobes, phi = 90 - 180 / lobes, phi0 = halfPi; i < lobes; ++i, phi -= delta, phi0 -= delta0) {
sphereStream.point(atan2(sr * cos(phi0), cr) * degrees, asin(sr * sin(phi0)) * degrees);
if (phi < -90) {
sphereStream.point(-90, -180 - phi - epsilon$$1);
sphereStream.point(-90, -180 - phi + epsilon$$1);
} else {
sphereStream.point(90, phi + epsilon$$1);
sphereStream.point(90, phi - epsilon$$1);
}
}
sphereStream.lineEnd(), sphereStream.polygonEnd();
};
return rotateStream;
};
return p
.scale(87.8076)
.center([0, 17.1875])
.clipAngle(180 - 1e-3);
};
function mollweideBromleyTheta(cp, phi) {
var cpsinPhi = cp * sin(phi), i = 30, delta;
do phi -= delta = (phi + sin(phi) - cpsinPhi) / (1 + cos(phi));
while (abs(delta) > epsilon && --i > 0);
return phi / 2;
}
function mollweideBromleyRaw(cx, cy, cp) {
function forward(lambda, phi) {
return [cx * lambda * cos(phi = mollweideBromleyTheta(cp, phi)), cy * sin(phi)];
}
forward.invert = function(x, y) {
return y = asin(y / cy), [x / (cx * cos(y)), asin((2 * y + sin(2 * y)) / cp)];
};
return forward;
}
var mollweideRaw = mollweideBromleyRaw(sqrt2 / halfPi, sqrt2, pi);
var mollweide = function() {
return d3Geo.geoProjection(mollweideRaw)
.scale(169.529);
};
var k = 2.00276;
var w = 1.11072;
function boggsRaw(lambda, phi) {
var theta = mollweideBromleyTheta(pi, phi);
return [k * lambda / (1 / cos(phi) + w / cos(theta)), (phi + sqrt2 * sin(theta)) / k];
}
boggsRaw.invert = function(x, y) {
var ky = k * y, theta = y < 0 ? -quarterPi : quarterPi, i = 25, delta, phi;
do {
phi = ky - sqrt2 * sin(theta);
theta -= delta = (sin(2 * theta) + 2 * theta - pi * sin(phi)) / (2 * cos(2 * theta) + 2 + pi * cos(phi) * sqrt2 * cos(theta));
} while (abs(delta) > epsilon && --i > 0);
phi = ky - sqrt2 * sin(theta);
return [x * (1 / cos(phi) + w / cos(theta)) / k, phi];
};
var boggs = function() {
return d3Geo.geoProjection(boggsRaw)
.scale(160.857);
};
var parallel1 = function(projectAt) {
var phi0 = 0,
m = d3Geo.geoProjectionMutator(projectAt),
p = m(phi0);
p.parallel = function(_) {
return arguments.length ? m(phi0 = _ * radians) : phi0 * degrees;
};
return p;
};
function sinusoidalRaw(lambda, phi) {
return [lambda * cos(phi), phi];
}
sinusoidalRaw.invert = function(x, y) {
return [x / cos(y), y];
};
var sinusoidal = function() {
return d3Geo.geoProjection(sinusoidalRaw)
.scale(152.63);
};
function bonneRaw(phi0) {
if (!phi0) return sinusoidalRaw;
var cotPhi0 = 1 / tan(phi0);
function forward(lambda, phi) {
var rho = cotPhi0 + phi0 - phi,
e = rho ? lambda * cos(phi) / rho : rho;
return [rho * sin(e), cotPhi0 - rho * cos(e)];
}
forward.invert = function(x, y) {
var rho = sqrt(x * x + (y = cotPhi0 - y) * y),
phi = cotPhi0 + phi0 - rho;
return [rho / cos(phi) * atan2(x, y), phi];
};
return forward;
}
var bonne = function() {
return parallel1(bonneRaw)
.scale(123.082)
.center([0, 26.1441])
.parallel(45);
};
function bottomleyRaw(sinPsi) {
function forward(lambda, phi) {
var rho = halfPi - phi,
eta = rho ? lambda * sinPsi * sin(rho) / rho : rho;
return [rho * sin(eta) / sinPsi, halfPi - rho * cos(eta)];
}
forward.invert = function(x, y) {
var x1 = x * sinPsi,
y1 = halfPi - y,
rho = sqrt(x1 * x1 + y1 * y1),
eta = atan2(x1, y1);
return [(rho ? rho / sin(rho) : 1) * eta / sinPsi, halfPi - rho];
};
return forward;
}
var bottomley = function() {
var sinPsi = 0.5,
m = d3Geo.geoProjectionMutator(bottomleyRaw),
p = m(sinPsi);
p.fraction = function(_) {
return arguments.length ? m(sinPsi = +_) : sinPsi;
};
return p
.scale(158.837);
};
var bromleyRaw = mollweideBromleyRaw(1, 4 / pi, pi);
var bromley = function() {
return d3Geo.geoProjection(bromleyRaw)
.scale(152.63);
};
// Azimuthal distance.
function distance(dPhi, c1, s1, c2, s2, dLambda) {
var cosdLambda = cos(dLambda), r;
if (abs(dPhi) > 1 || abs(dLambda) > 1) {
r = acos(s1 * s2 + c1 * c2 * cosdLambda);
} else {
var sindPhi = sin(dPhi / 2), sindLambda = sin(dLambda / 2);
r = 2 * asin(sqrt(sindPhi * sindPhi + c1 * c2 * sindLambda * sindLambda));
}
return abs(r) > epsilon ? [r, atan2(c2 * sin(dLambda), c1 * s2 - s1 * c2 * cosdLambda)] : [0, 0];
}
// Angle opposite a, and contained between sides of lengths b and c.
function angle(b, c, a) {
return acos((b * b + c * c - a * a) / (2 * b * c));
}
// Normalize longitude.
function longitude(lambda) {
return lambda - 2 * pi * floor((lambda + pi) / (2 * pi));
}
function chamberlinRaw(p0, p1, p2) {
var points = [
[p0[0], p0[1], sin(p0[1]), cos(p0[1])],
[p1[0], p1[1], sin(p1[1]), cos(p1[1])],
[p2[0], p2[1], sin(p2[1]), cos(p2[1])]
];
for (var a = points[2], b, i = 0; i < 3; ++i, a = b) {
b = points[i];
a.v = distance(b[1] - a[1], a[3], a[2], b[3], b[2], b[0] - a[0]);
a.point = [0, 0];
}
var beta0 = angle(points[0].v[0], points[2].v[0], points[1].v[0]),
beta1 = angle(points[0].v[0], points[1].v[0], points[2].v[0]),
beta2 = pi - beta0;
points[2].point[1] = 0;
points[0].point[0] = -(points[1].point[0] = points[0].v[0] / 2);
var mean = [
points[2].point[0] = points[0].point[0] + points[2].v[0] * cos(beta0),
2 * (points[0].point[1] = points[1].point[1] = points[2].v[0] * sin(beta0))
];
function forward(lambda, phi) {
var sinPhi = sin(phi),
cosPhi = cos(phi),
v = new Array(3), i;
// Compute distance and azimuth from control points.
for (i = 0; i < 3; ++i) {
var p = points[i];
v[i] = distance(phi - p[1], p[3], p[2], cosPhi, sinPhi, lambda - p[0]);
if (!v[i][0]) return p.point;
v[i][1] = longitude(v[i][1] - p.v[1]);
}
// Arithmetic mean of interception points.
var point = mean.slice();
for (i = 0; i < 3; ++i) {
var j = i == 2 ? 0 : i + 1;
var a = angle(points[i].v[0], v[i][0], v[j][0]);
if (v[i][1] < 0) a = -a;
if (!i) {
point[0] += v[i][0] * cos(a);
point[1] -= v[i][0] * sin(a);
} else if (i == 1) {
a = beta1 - a;
point[0] -= v[i][0] * cos(a);
point[1] -= v[i][0] * sin(a);
} else {
a = beta2 - a;
point[0] += v[i][0] * cos(a);
point[1] += v[i][0] * sin(a);
}
}
point[0] /= 3, point[1] /= 3;
return point;
}
return forward;
}
function pointRadians(p) {
return p[0] *= radians, p[1] *= radians, p;
}
function chamberlinAfrica() {
return chamberlin([0, 22], [45, 22], [22.5, -22])
.scale(380)
.center([22.5, 2]);
}
function chamberlin(p0, p1, p2) { // TODO order matters!
var c = d3Geo.geoCentroid({type: "MultiPoint", coordinates: [p0, p1, p2]}),
R = [-c[0], -c[1]],
r = d3Geo.geoRotation(R),
p = d3Geo.geoProjection(chamberlinRaw(pointRadians(r(p0)), pointRadians(r(p1)), pointRadians(r(p2)))).rotate(R),
center = p.center;
delete p.rotate;
p.center = function(_) {
return arguments.length ? center(r(_)) : r.invert(center());
};
return p
.clipAngle(90);
}
function collignonRaw(lambda, phi) {
var alpha = sqrt(1 - sin(phi));
return [(2 / sqrtPi) * lambda * alpha, sqrtPi * (1 - alpha)];
}
collignonRaw.invert = function(x, y) {
var lambda = (lambda = y / sqrtPi - 1) * lambda;
return [lambda > 0 ? x * sqrt(pi / lambda) / 2 : 0, asin(1 - lambda)];
};
var collignon = function() {
return d3Geo.geoProjection(collignonRaw)
.scale(95.6464)
.center([0, 30]);
};
function craigRaw(phi0) {
var tanPhi0 = tan(phi0);
function forward(lambda, phi) {
return [lambda, (lambda ? lambda / sin(lambda) : 1) * (sin(phi) * cos(lambda) - tanPhi0 * cos(phi))];
}
forward.invert = tanPhi0 ? function(x, y) {
if (x) y *= sin(x) / x;
var cosLambda = cos(x);
return [x, 2 * atan2(sqrt(cosLambda * cosLambda + tanPhi0 * tanPhi0 - y * y) - cosLambda, tanPhi0 - y)];
} : function(x, y) {
return [x, asin(x ? y * tan(x) / x : y)];
};
return forward;
}
var craig = function() {
return parallel1(craigRaw)
.scale(249.828)
.clipAngle(90);
};
var sqrt3 = sqrt(3);
function crasterRaw(lambda, phi) {
return [sqrt3 * lambda * (2 * cos(2 * phi / 3) - 1) / sqrtPi, sqrt3 * sqrtPi * sin(phi / 3)];
}
crasterRaw.invert = function(x, y) {
var phi = 3 * asin(y / (sqrt3 * sqrtPi));
return [sqrtPi * x / (sqrt3 * (2 * cos(2 * phi / 3) - 1)), phi];
};
var craster = function() {
return d3Geo.geoProjection(crasterRaw)
.scale(156.19);
};
function cylindricalEqualAreaRaw(phi0) {
var cosPhi0 = cos(phi0);
function forward(lambda, phi) {
return [lambda * cosPhi0, sin(phi) / cosPhi0];
}
forward.invert = function(x, y) {
return [x / cosPhi0, asin(y * cosPhi0)];
};
return forward;
}
var cylindricalEqualArea = function() {
return parallel1(cylindricalEqualAreaRaw)
.parallel(38.58) // acos(sqrt(width / height / pi)) * radians
.scale(195.044); // width / (sqrt(width / height / pi) * 2 * pi)
};
function cylindricalStereographicRaw(phi0) {
var cosPhi0 = cos(phi0);
function forward(lambda, phi) {
return [lambda * cosPhi0, (1 + cosPhi0) * tan(phi / 2)];
}
forward.invert = function(x, y) {
return [x / cosPhi0, atan(y / (1 + cosPhi0)) * 2];
};
return forward;
}
var cylindricalStereographic = function() {
return parallel1(cylindricalStereographicRaw)
.scale(124.75);
};
function eckert1Raw(lambda, phi) {
var alpha = sqrt(8 / (3 * pi));
return [
alpha * lambda * (1 - abs(phi) / pi),
alpha * phi
];
}
eckert1Raw.invert = function(x, y) {
var alpha = sqrt(8 / (3 * pi)),
phi = y / alpha;
return [
x / (alpha * (1 - abs(phi) / pi)),
phi
];
};
var eckert1 = function() {
return d3Geo.geoProjection(eckert1Raw)
.scale(165.664);
};
function eckert2Raw(lambda, phi) {
var alpha = sqrt(4 - 3 * sin(abs(phi)));
return [
2 / sqrt(6 * pi) * lambda * alpha,
sign(phi) * sqrt(2 * pi / 3) * (2 - alpha)
];
}
eckert2Raw.invert = function(x, y) {
var alpha = 2 - abs(y) / sqrt(2 * pi / 3);
return [
x * sqrt(6 * pi) / (2 * alpha),
sign(y) * asin((4 - alpha * alpha) / 3)
];
};
var eckert2 = function() {
return d3Geo.geoProjection(eckert2Raw)
.scale(165.664);
};
function eckert3Raw(lambda, phi) {
var k = sqrt(pi * (4 + pi));
return [
2 / k * lambda * (1 + sqrt(1 - 4 * phi * phi / (pi * pi))),
4 / k * phi
];
}
eckert3Raw.invert = function(x, y) {
var k = sqrt(pi * (4 + pi)) / 2;
return [
x * k / (1 + sqrt(1 - y * y * (4 + pi) / (4 * pi))),
y * k / 2
];
};
var eckert3 = function() {
return d3Geo.geoProjection(eckert3Raw)
.scale(180.739);
};
function eckert4Raw(lambda, phi) {
var k = (2 + halfPi) * sin(phi);
phi /= 2;
for (var i = 0, delta = Infinity; i < 10 && abs(delta) > epsilon; i++) {
var cosPhi = cos(phi);
phi -= delta = (phi + sin(phi) * (cosPhi + 2) - k) / (2 * cosPhi * (1 + cosPhi));
}
return [
2 / sqrt(pi * (4 + pi)) * lambda * (1 + cos(phi)),
2 * sqrt(pi / (4 + pi)) * sin(phi)
];
}
eckert4Raw.invert = function(x, y) {
var A = y * sqrt((4 + pi) / pi) / 2,
k = asin(A),
c = cos(k);
return [
x / (2 / sqrt(pi * (4 + pi)) * (1 + c)),
asin((k + A * (c + 2)) / (2 + halfPi))
];
};
var eckert4 = function() {
return d3Geo.geoProjection(eckert4Raw)
.scale(180.739);
};
function eckert5Raw(lambda, phi) {
return [
lambda * (1 + cos(phi)) / sqrt(2 + pi),
2 * phi / sqrt(2 + pi)
];
}
eckert5Raw.invert = function(x, y) {
var k = sqrt(2 + pi),
phi = y * k / 2;
return [
k * x / (1 + cos(phi)),
phi
];
};
var eckert5 = function() {
return d3Geo.geoProjection(eckert5Raw)
.scale(173.044);
};
function eckert6Raw(lambda, phi) {
var k = (1 + halfPi) * sin(phi);
for (var i = 0, delta = Infinity; i < 10 && abs(delta) > epsilon; i++) {
phi -= delta = (phi + sin(phi) - k) / (1 + cos(phi));
}
k = sqrt(2 + pi);
return [
lambda * (1 + cos(phi)) / k,
2 * phi / k
];
}
eckert6Raw.invert = function(x, y) {
var j = 1 + halfPi,
k = sqrt(j / 2);
return [
x * 2 * k / (1 + cos(y *= k)),
asin((y + sin(y)) / j)
];
};
var eckert6 = function() {
return d3Geo.geoProjection(eckert6Raw)
.scale(173.044);
};
var eisenlohrK = 3 + 2 * sqrt2;
function eisenlohrRaw(lambda, phi) {
var s0 = sin(lambda /= 2),
c0 = cos(lambda),
k = sqrt(cos(phi)),
c1 = cos(phi /= 2),
t = sin(phi) / (c1 + sqrt2 * c0 * k),
c = sqrt(2 / (1 + t * t)),
v = sqrt((sqrt2 * c1 + (c0 + s0) * k) / (sqrt2 * c1 + (c0 - s0) * k));
return [
eisenlohrK * (c * (v - 1 / v) - 2 * log(v)),
eisenlohrK * (c * t * (v + 1 / v) - 2 * atan(t))
];
}
eisenlohrRaw.invert = function(x, y) {
if (!(p = augustRaw.invert(x / 1.2, y * 1.065))) return null;
var lambda = p[0], phi = p[1], i = 20, p;
x /= eisenlohrK, y /= eisenlohrK;
do {
var _0 = lambda / 2,
_1 = phi / 2,
s0 = sin(_0),
c0 = cos(_0),
s1 = sin(_1),
c1 = cos(_1),
cos1 = cos(phi),
k = sqrt(cos1),
t = s1 / (c1 + sqrt2 * c0 * k),
t2 = t * t,
c = sqrt(2 / (1 + t2)),
v0 = (sqrt2 * c1 + (c0 + s0) * k),
v1 = (sqrt2 * c1 + (c0 - s0) * k),
v2 = v0 / v1,
v = sqrt(v2),
vm1v = v - 1 / v,
vp1v = v + 1 / v,
fx = c * vm1v - 2 * log(v) - x,
fy = c * t * vp1v - 2 * atan(t) - y,
deltatDeltaLambda = s1 && sqrt1_2 * k * s0 * t2 / s1,
deltatDeltaPhi = (sqrt2 * c0 * c1 + k) / (2 * (c1 + sqrt2 * c0 * k) * (c1 + sqrt2 * c0 * k) * k),
deltacDeltat = -0.5 * t * c * c * c,
deltacDeltaLambda = deltacDeltat * deltatDeltaLambda,
deltacDeltaPhi = deltacDeltat * deltatDeltaPhi,
A = (A = 2 * c1 + sqrt2 * k * (c0 - s0)) * A * v,
deltavDeltaLambda = (sqrt2 * c0 * c1 * k + cos1) / A,
deltavDeltaPhi = -(sqrt2 * s0 * s1) / (k * A),
deltaxDeltaLambda = vm1v * deltacDeltaLambda - 2 * deltavDeltaLambda / v + c * (deltavDeltaLambda + deltavDeltaLambda / v2),
deltaxDeltaPhi = vm1v * deltacDeltaPhi - 2 * deltavDeltaPhi / v + c * (deltavDeltaPhi + deltavDeltaPhi / v2),
deltayDeltaLambda = t * vp1v * deltacDeltaLambda - 2 * deltatDeltaLambda / (1 + t2) + c * vp1v * deltatDeltaLambda + c * t * (deltavDeltaLambda - deltavDeltaLambda / v2),
deltayDeltaPhi = t * vp1v * deltacDeltaPhi - 2 * deltatDeltaPhi / (1 + t2) + c * vp1v * deltatDeltaPhi + c * t * (deltavDeltaPhi - deltavDeltaPhi / v2),
denominator = deltaxDeltaPhi * deltayDeltaLambda - deltayDeltaPhi * deltaxDeltaLambda;
if (!denominator) break;
var deltaLambda = (fy * deltaxDeltaPhi - fx * deltayDeltaPhi) / denominator,
deltaPhi = (fx * deltayDeltaLambda - fy * deltaxDeltaLambda) / denominator;
lambda -= deltaLambda;
phi = max(-halfPi, min(halfPi, phi - deltaPhi));
} while ((abs(deltaLambda) > epsilon || abs(deltaPhi) > epsilon) && --i > 0);
return abs(abs(phi) - halfPi) < epsilon ? [0, phi] : i && [lambda, phi];
};
var eisenlohr = function() {
return d3Geo.geoProjection(eisenlohrRaw)
.scale(62.5271);
};
var faheyK = cos(35 * radians);
function faheyRaw(lambda, phi) {
var t = tan(phi / 2);
return [lambda * faheyK * sqrt(1 - t * t), (1 + faheyK) * t];
}
faheyRaw.invert = function(x, y) {
var t = y / (1 + faheyK);
return [x && x / (faheyK * sqrt(1 - t * t)), 2 * atan(t)];
};
var fahey = function() {
return d3Geo.geoProjection(faheyRaw)
.scale(137.152);
};
function foucautRaw(lambda, phi) {
var k = phi / 2, cosk = cos(k);
return [ 2 * lambda / sqrtPi * cos(phi) * cosk * cosk, sqrtPi * tan(k)];
}
foucautRaw.invert = function(x, y) {
var k = atan(y / sqrtPi), cosk = cos(k), phi = 2 * k;
return [x * sqrtPi / 2 / (cos(phi) * cosk * cosk), phi];
};
var foucaut = function() {
return d3Geo.geoProjection(foucautRaw)
.scale(135.264);
};
function gilbertForward(point) {
return [point[0] / 2, asin(tan(point[1] / 2 * radians)) * degrees];
}
function gilbertInvert(point) {
return [point[0] * 2, 2 * atan(sin(point[1] * radians)) * degrees];
}
var gilbert = function(projectionType) {
if (projectionType == null) projectionType = d3Geo.geoOrthographic;
var projection = projectionType(),
equirectangular = d3Geo.geoEquirectangular().scale(degrees).precision(0).clipAngle(null).translate([0, 0]); // antimeridian cutting
function gilbert(point) {
return projection(gilbertForward(point));
}
if (projection.invert) gilbert.invert = function(point) {
return gilbertInvert(projection.invert(point));
};
gilbert.stream = function(stream) {
var s1 = projection.stream(stream), s0 = equirectangular.stream({
point: function(lambda, phi) { s1.point(lambda / 2, asin(tan(-phi / 2 * radians)) * degrees); },
lineStart: function() { s1.lineStart(); },
lineEnd: function() { s1.lineEnd(); },
polygonStart: function() { s1.polygonStart(); },
polygonEnd: function() { s1.polygonEnd(); }
});
s0.sphere = s1.sphere;
return s0;
};
function property(name) {
gilbert[name] = function(_) {
return arguments.length ? (projection[name](_), gilbert) : projection[name]();
};
}
gilbert.rotate = function(_) {
return arguments.length ? (equirectangular.rotate(_), gilbert) : equirectangular.rotate();
};
gilbert.center = function(_) {
return arguments.length ? (projection.center(gilbertForward(_)), gilbert) : gilbertInvert(projection.center());
};
property("clipAngle");
property("clipExtent");
property("scale");
property("translate");
property("precision");
return gilbert
.scale(249.5);
};
function gingeryRaw(rho, n) {
var k = 2 * pi / n,
rho2 = rho * rho;
function forward(lambda, phi) {
var p = d3Geo.geoAzimuthalEquidistantRaw(lambda, phi),
x = p[0],
y = p[1],
r2 = x * x + y * y;
if (r2 > rho2) {
var r = sqrt(r2),
theta = atan2(y, x),
theta0 = k * round(theta / k),
alpha = theta - theta0,
rhoCosAlpha = rho * cos(alpha),
k_ = (rho * sin(alpha) - alpha * sin(rhoCosAlpha)) / (halfPi - rhoCosAlpha),
s_ = gingeryLength(alpha, k_),
e = (pi - rho) / gingeryIntegrate(s_, rhoCosAlpha, pi);
x = r;
var i = 50, delta;
do {
x -= delta = (rho + gingeryIntegrate(s_, rhoCosAlpha, x) * e - r) / (s_(x) * e);
} while (abs(delta) > epsilon && --i > 0);
y = alpha * sin(x);
if (x < halfPi) y -= k_ * (x - halfPi);
var s = sin(theta0),
c = cos(theta0);
p[0] = x * c - y * s;
p[1] = x * s + y * c;
}
return p;
}
forward.invert = function(x, y) {
var r2 = x * x + y * y;
if (r2 > rho2) {
var r = sqrt(r2),
theta = atan2(y, x),
theta0 = k * round(theta / k),
dTheta = theta - theta0;
x = r * cos(dTheta);
y = r * sin(dTheta);
var x_halfPi = x - halfPi,
sinx = sin(x),
alpha = y / sinx,
delta = x < halfPi ? Infinity : 0,
i = 10;
while (true) {
var rhosinAlpha = rho * sin(alpha),
rhoCosAlpha = rho * cos(alpha),
sinRhoCosAlpha = sin(rhoCosAlpha),
halfPi_RhoCosAlpha = halfPi - rhoCosAlpha,
k_ = (rhosinAlpha - alpha * sinRhoCosAlpha) / halfPi_RhoCosAlpha,
s_ = gingeryLength(alpha, k_);
if (abs(delta) < epsilon2 || !--i) break;
alpha -= delta = (alpha * sinx - k_ * x_halfPi - y) / (
sinx - x_halfPi * 2 * (
halfPi_RhoCosAlpha * (rhoCosAlpha + alpha * rhosinAlpha * cos(rhoCosAlpha) - sinRhoCosAlpha) -
rhosinAlpha * (rhosinAlpha - alpha * sinRhoCosAlpha)
) / (halfPi_RhoCosAlpha * halfPi_RhoCosAlpha));
}
r = rho + gingeryIntegrate(s_, rhoCosAlpha, x) * (pi - rho) / gingeryIntegrate(s_, rhoCosAlpha, pi);
theta = theta0 + alpha;
x = r * cos(theta);
y = r * sin(theta);
}
return d3Geo.geoAzimuthalEquidistantRaw.invert(x, y);
};
return forward;
}
function gingeryLength(alpha, k) {
return function(x) {
var y_ = alpha * cos(x);
if (x < halfPi) y_ -= k;
return sqrt(1 + y_ * y_);
};
}
// Numerical integration: trapezoidal rule.
function gingeryIntegrate(f, a, b) {
var n = 50,
h = (b - a) / n,
s = f(a) + f(b);
for (var i = 1, x = a; i < n; ++i) s += 2 * f(x += h);
return s * 0.5 * h;
}
var gingery = function() {
var n = 6,
rho = 30 * radians,
cRho = cos(rho),
sRho = sin(rho),
m = d3Geo.geoProjectionMutator(gingeryRaw),
p = m(rho, n),
stream_ = p.stream,
epsilon$$1 = 1e-2,
cr = -cos(epsilon$$1 * radians),
sr = sin(epsilon$$1 * radians);
p.radius = function(_) {
if (!arguments.length) return rho * degrees;
cRho = cos(rho = _ * radians);
sRho = sin(rho);
return m(rho, n);
};
p.lobes = function(_) {
if (!arguments.length) return n;
return m(rho, n = +_);
};
p.stream = function(stream) {
var rotate = p.rotate(),
rotateStream = stream_(stream),
sphereStream = (p.rotate([0, 0]), stream_(stream));
p.rotate(rotate);
rotateStream.sphere = function() {
sphereStream.polygonStart(), sphereStream.lineStart();
for (var i = 0, delta = 2 * pi / n, phi = 0; i < n; ++i, phi -= delta) {
sphereStream.point(atan2(sr * cos(phi), cr) * degrees, asin(sr * sin(phi)) * degrees);
sphereStream.point(atan2(sRho * cos(phi - delta / 2), cRho) * degrees, asin(sRho * sin(phi - delta / 2)) * degrees);
}
sphereStream.lineEnd(), sphereStream.polygonEnd();
};
return rotateStream;
};
return p
.rotate([90, -40])
.scale(91.7095)
.clipAngle(180 - 1e-3);
};
var ginzburgPolyconicRaw = function(a, b, c, d, e, f, g, h) {
if (arguments.length < 8) h = 0;
function forward(lambda, phi) {
if (!phi) return [a * lambda / pi, 0];
var phi2 = phi * phi,
xB = a + phi2 * (b + phi2 * (c + phi2 * d)),
yB = phi * (e - 1 + phi2 * (f - h + phi2 * g)),
m = (xB * xB + yB * yB) / (2 * yB),
alpha = lambda * asin(xB / m) / pi;
return [m * sin(alpha), phi * (1 + phi2 * h) + m * (1 - cos(alpha))];
}
forward.invert = function(x, y) {
var lambda = pi * x / a,
phi = y,
deltaLambda, deltaPhi, i = 50;
do {
var phi2 = phi * phi,
xB = a + phi2 * (b + phi2 * (c + phi2 * d)),
yB = phi * (e - 1 + phi2 * (f - h + phi2 * g)),
p = xB * xB + yB * yB,
q = 2 * yB,
m = p / q,
m2 = m * m,
dAlphadLambda = asin(xB / m) / pi,
alpha = lambda * dAlphadLambda,
xB2 = xB * xB,
dxBdPhi = (2 * b + phi2 * (4 * c + phi2 * 6 * d)) * phi,
dyBdPhi = e + phi2 * (3 * f + phi2 * 5 * g),
dpdPhi = 2 * (xB * dxBdPhi + yB * (dyBdPhi - 1)),
dqdPhi = 2 * (dyBdPhi - 1),
dmdPhi = (dpdPhi * q - p * dqdPhi) / (q * q),
cosAlpha = cos(alpha),
sinAlpha = sin(alpha),
mcosAlpha = m * cosAlpha,
msinAlpha = m * sinAlpha,
dAlphadPhi = ((lambda / pi) * (1 / sqrt(1 - xB2 / m2)) * (dxBdPhi * m - xB * dmdPhi)) / m2,
fx = msinAlpha - x,
fy = phi * (1 + phi2 * h) + m - mcosAlpha - y,
deltaxDeltaPhi = dmdPhi * sinAlpha + mcosAlpha * dAlphadPhi,
deltaxDeltaLambda = mcosAlpha * dAlphadLambda,
deltayDeltaPhi = 1 + dmdPhi - (dmdPhi * cosAlpha - msinAlpha * dAlphadPhi),
deltayDeltaLambda = msinAlpha * dAlphadLambda,
denominator = deltaxDeltaPhi * deltayDeltaLambda - deltayDeltaPhi * deltaxDeltaLambda;
if (!denominator) break;
lambda -= deltaLambda = (fy * deltaxDeltaPhi - fx * deltayDeltaPhi) / denominator;
phi -= deltaPhi = (fx * deltayDeltaLambda - fy * deltaxDeltaLambda) / denominator;
} while ((abs(deltaLambda) > epsilon || abs(deltaPhi) > epsilon) && --i > 0);
return [lambda, phi];
};
return forward;
};
var ginzburg4Raw = ginzburgPolyconicRaw(2.8284, -1.6988, 0.75432, -0.18071, 1.76003, -0.38914, 0.042555);
var ginzburg4 = function() {
return d3Geo.geoProjection(ginzburg4Raw)
.scale(149.995);
};
var ginzburg5Raw = ginzburgPolyconicRaw(2.583819, -0.835827, 0.170354, -0.038094, 1.543313, -0.411435,0.082742);
var ginzburg5 = function() {
return d3Geo.geoProjection(ginzburg5Raw)
.scale(153.93);
};
var ginzburg6Raw = ginzburgPolyconicRaw(5 / 6 * pi, -0.62636, -0.0344, 0, 1.3493, -0.05524, 0, 0.045);
var ginzburg6 = function() {
return d3Geo.geoProjection(ginzburg6Raw)
.scale(130.945);
};
function ginzburg8Raw(lambda, phi) {
var lambda2 = lambda * lambda,
phi2 = phi * phi;
return [
lambda * (1 - 0.162388 * phi2) * (0.87 - 0.000952426 * lambda2 * lambda2),
phi * (1 + phi2 / 12)
];
}
ginzburg8Raw.invert = function(x, y) {
var lambda = x,
phi = y,
i = 50, delta;
do {
var phi2 = phi * phi;
phi -= delta = (phi * (1 + phi2 / 12) - y) / (1 + phi2 / 4);
} while (abs(delta) > epsilon && --i > 0);
i = 50;
x /= 1 -0.162388 * phi2;
do {
var lambda4 = (lambda4 = lambda * lambda) * lambda4;
lambda -= delta = (lambda * (0.87 - 0.000952426 * lambda4) - x) / (0.87 - 0.00476213 * lambda4);
} while (abs(delta) > epsilon && --i > 0);
return [lambda, phi];
};
var ginzburg8 = function() {
return d3Geo.geoProjection(ginzburg8Raw)
.scale(131.747);
};
var ginzburg9Raw = ginzburgPolyconicRaw(2.6516, -0.76534, 0.19123, -0.047094, 1.36289, -0.13965,0.031762);
var ginzburg9 = function() {
return d3Geo.geoProjection(ginzburg9Raw)
.scale(131.087);
};
var squareRaw = function(project) {
var dx = project(halfPi, 0)[0] - project(-halfPi, 0)[0];
function projectSquare(lambda, phi) {
var s = lambda > 0 ? -0.5 : 0.5,
point = project(lambda + s * pi, phi);
point[0] -= s * dx;
return point;
}
if (project.invert) projectSquare.invert = function(x, y) {
var s = x > 0 ? -0.5 : 0.5,
location = project.invert(x + s * dx, y),
lambda = location[0] - s * pi;
if (lambda < -pi) lambda += 2 * pi;
else if (lambda > pi) lambda -= 2 * pi;
location[0] = lambda;
return location;
};
return projectSquare;
};
function gringortenRaw(lambda, phi) {
var sLambda = sign(lambda),
sPhi = sign(phi),
cosPhi = cos(phi),
x = cos(lambda) * cosPhi,
y = sin(lambda) * cosPhi,
z = sin(sPhi * phi);
lambda = abs(atan2(y, z));
phi = asin(x);
if (abs(lambda - halfPi) > epsilon) lambda %= halfPi;
var point = gringortenHexadecant(lambda > pi / 4 ? halfPi - lambda : lambda, phi);
if (lambda > pi / 4) z = point[0], point[0] = -point[1], point[1] = -z;
return (point[0] *= sLambda, point[1] *= -sPhi, point);
}
gringortenRaw.invert = function(x, y) {
if (abs(x) > 1) x = sign(x) * 2 - x;
if (abs(y) > 1) y = sign(y) * 2 - y;
var sx = sign(x),
sy = sign(y),
x0 = -sx * x,
y0 = -sy * y,
t = y0 / x0 < 1,
p = gringortenHexadecantInvert(t ? y0 : x0, t ? x0 : y0),
lambda = p[0],
phi = p[1],
cosPhi = cos(phi);
if (t) lambda = -halfPi - lambda;
return [sx * (atan2(sin(lambda) * cosPhi, -sin(phi)) + pi), sy * asin(cos(lambda) * cosPhi)];
};
function gringortenHexadecant(lambda, phi) {
if (phi === halfPi) return [0, 0];
var sinPhi = sin(phi),
r = sinPhi * sinPhi,
r2 = r * r,
j = 1 + r2,
k = 1 + 3 * r2,
q = 1 - r2,
z = asin(1 / sqrt(j)),
v = q + r * j * z,
p2 = (1 - sinPhi) / v,
p = sqrt(p2),
a2 = p2 * j,
a = sqrt(a2),
h = p * q,
x,
i;
if (lambda === 0) return [0, -(h + r * a)];
var cosPhi = cos(phi),
secPhi = 1 / cosPhi,
drdPhi = 2 * sinPhi * cosPhi,
dvdPhi = (-3 * r + z * k) * drdPhi,
dp2dPhi = (-v * cosPhi - (1 - sinPhi) * dvdPhi) / (v * v),
dpdPhi = (0.5 * dp2dPhi) / p,
dhdPhi = q * dpdPhi - 2 * r * p * drdPhi,
dra2dPhi = r * j * dp2dPhi + p2 * k * drdPhi,
mu = -secPhi * drdPhi,
nu = -secPhi * dra2dPhi,
zeta = -2 * secPhi * dhdPhi,
lambda1 = 4 * lambda / pi,
delta;
// Slower but accurate bisection method.
if (lambda > 0.222 * pi || phi < pi / 4 && lambda > 0.175 * pi) {
x = (h + r * sqrt(a2 * (1 + r2) - h * h)) / (1 + r2);
if (lambda > pi / 4) return [x, x];
var x1 = x, x0 = 0.5 * x;
x = 0.5 * (x0 + x1), i = 50;
do {
var g = sqrt(a2 - x * x),
f = (x * (zeta + mu * g) + nu * asin(x / a)) - lambda1;
if (!f) break;
if (f < 0) x0 = x;
else x1 = x;
x = 0.5 * (x0 + x1);
} while (abs(x1 - x0) > epsilon && --i > 0);
}
// Newton-Raphson.
else {
x = epsilon, i = 25;
do {
var x2 = x * x,
g2 = sqrt(a2 - x2),
zetaMug = zeta + mu * g2,
f2 = x * zetaMug + nu * asin(x / a) - lambda1,
df = zetaMug + (nu - mu * x2) / g2;
x -= delta = g2 ? f2 / df : 0;
} while (abs(delta) > epsilon && --i > 0);
}
return [x, -h - r * sqrt(a2 - x * x)];
}
function gringortenHexadecantInvert(x, y) {
var x0 = 0,
x1 = 1,
r = 0.5,
i = 50;
while (true) {
var r2 = r * r,
sinPhi = sqrt(r),
z = asin(1 / sqrt(1 + r2)),
v = (1 - r2) + r * (1 + r2) * z,
p2 = (1 - sinPhi) / v,
p = sqrt(p2),
a2 = p2 * (1 + r2),
h = p * (1 - r2),
g2 = a2 - x * x,
g = sqrt(g2),
y0 = y + h + r * g;
if (abs(x1 - x0) < epsilon2 || --i === 0 || y0 === 0) break;
if (y0 > 0) x0 = r;
else x1 = r;
r = 0.5 * (x0 + x1);
}
if (!i) return null;
var phi = asin(sinPhi),
cosPhi = cos(phi),
secPhi = 1 / cosPhi,
drdPhi = 2 * sinPhi * cosPhi,
dvdPhi = (-3 * r + z * (1 + 3 * r2)) * drdPhi,
dp2dPhi = (-v * cosPhi - (1 - sinPhi) * dvdPhi) / (v * v),
dpdPhi = 0.5 * dp2dPhi / p,
dhdPhi = (1 - r2) * dpdPhi - 2 * r * p * drdPhi,
zeta = -2 * secPhi * dhdPhi,
mu = -secPhi * drdPhi,
nu = -secPhi * (r * (1 + r2) * dp2dPhi + p2 * (1 + 3 * r2) * drdPhi);
return [pi / 4 * (x * (zeta + mu * g) + nu * asin(x / sqrt(a2))), phi];
}
var gringorten = function() {
return d3Geo.geoProjection(squareRaw(gringortenRaw))
.scale(239.75);
};
// Returns [sn, cn, dn](u + iv|m).
function ellipticJi(u, v, m) {
var a, b, c;
if (!u) {
b = ellipticJ(v, 1 - m);
return [
[0, b[0] / b[1]],
[1 / b[1], 0],
[b[2] / b[1], 0]
];
}
a = ellipticJ(u, m);
if (!v) return [[a[0], 0], [a[1], 0], [a[2], 0]];
b = ellipticJ(v, 1 - m);
c = b[1] * b[1] + m * a[0] * a[0] * b[0] * b[0];
return [
[a[0] * b[2] / c, a[1] * a[2] * b[0] * b[1] / c],
[a[1] * b[1] / c, -a[0] * a[2] * b[0] * b[2] / c],
[a[2] * b[1] * b[2] / c, -m * a[0] * a[1] * b[0] / c]
];
}
// Returns [sn, cn, dn, ph](u|m).
function ellipticJ(u, m) {
var ai, b, phi, t, twon;
if (m < epsilon) {
t = sin(u);
b = cos(u);
ai = m * (u - t * b) / 4;
return [
t - ai * b,
b + ai * t,
1 - m * t * t / 2,
u - ai
];
}
if (m >= 1 - epsilon) {
ai = (1 - m) / 4;
b = cosh(u);
t = tanh(u);
phi = 1 / b;
twon = b * sinh(u);
return [
t + ai * (twon - u) / (b * b),
phi - ai * t * phi * (twon - u),
phi + ai * t * phi * (twon + u),
2 * atan(exp(u)) - halfPi + ai * (twon - u) / b
];
}
var a = [1, 0, 0, 0, 0, 0, 0, 0, 0],
c = [sqrt(m), 0, 0, 0, 0, 0, 0, 0, 0],
i = 0;
b = sqrt(1 - m);
twon = 1;
while (abs(c[i] / a[i]) > epsilon && i < 8) {
ai = a[i++];
c[i] = (ai - b) / 2;
a[i] = (ai + b) / 2;
b = sqrt(ai * b);
twon *= 2;
}
phi = twon * a[i] * u;
do {
t = c[i] * sin(b = phi) / a[i];
phi = (asin(t) + phi) / 2;
} while (--i);
return [sin(phi), t = cos(phi), t / cos(phi - b), phi];
}
// Calculate F(phi+iPsi|m).
// See Abramowitz and Stegun, 17.4.11.
function ellipticFi(phi, psi, m) {
var r = abs(phi),
i = abs(psi),
sinhPsi = sinh(i);
if (r) {
var cscPhi = 1 / sin(r),
cotPhi2 = 1 / (tan(r) * tan(r)),
b = -(cotPhi2 + m * (sinhPsi * sinhPsi * cscPhi * cscPhi) - 1 + m),
c = (m - 1) * cotPhi2,
cotLambda2 = (-b + sqrt(b * b - 4 * c)) / 2;
return [
ellipticF(atan(1 / sqrt(cotLambda2)), m) * sign(phi),
ellipticF(atan(sqrt((cotLambda2 / cotPhi2 - 1) / m)), 1 - m) * sign(psi)
];
}
return [
0,
ellipticF(atan(sinhPsi), 1 - m) * sign(psi)
];
}
// Calculate F(phi|m) where m = k² = sin²α.
// See Abramowitz and Stegun, 17.6.7.
function ellipticF(phi, m) {
if (!m) return phi;
if (m === 1) return log(tan(phi / 2 + quarterPi));
var a = 1,
b = sqrt(1 - m),
c = sqrt(m);
for (var i = 0; abs(c) > epsilon; i++) {
if (phi % pi) {
var dPhi = atan(b * tan(phi) / a);
if (dPhi < 0) dPhi += pi;
phi += dPhi + ~~(phi / pi) * pi;
} else phi += phi;
c = (a + b) / 2;
b = sqrt(a * b);
c = ((a = c) - b) / 2;
}
return phi / (pow(2, i) * a);
}
function guyouRaw(lambda, phi) {
var k_ = (sqrt2 - 1) / (sqrt2 + 1),
k = sqrt(1 - k_ * k_),
K = ellipticF(halfPi, k * k),
f = -1,
psi = log(tan(pi / 4 + abs(phi) / 2)),
r = exp(f * psi) / sqrt(k_),
at = guyouComplexAtan(r * cos(f * lambda), r * sin(f * lambda)),
t = ellipticFi(at[0], at[1], k * k);
return [-t[1], (phi >= 0 ? 1 : -1) * (0.5 * K - t[0])];
}
function guyouComplexAtan(x, y) {
var x2 = x * x,
y_1 = y + 1,
t = 1 - x2 - y * y;
return [
0.5 * ((x >= 0 ? halfPi : -halfPi) - atan2(t, 2 * x)),
-0.25 * log(t * t + 4 * x2) +0.5 * log(y_1 * y_1 + x2)
];
}
function guyouComplexDivide(a, b) {
var denominator = b[0] * b[0] + b[1] * b[1];
return [
(a[0] * b[0] + a[1] * b[1]) / denominator,
(a[1] * b[0] - a[0] * b[1]) / denominator
];
}
guyouRaw.invert = function(x, y) {
var k_ = (sqrt2 - 1) / (sqrt2 + 1),
k = sqrt(1 - k_ * k_),
K = ellipticF(halfPi, k * k),
f = -1,
j = ellipticJi(0.5 * K - y, -x, k * k),
tn = guyouComplexDivide(j[0], j[1]),
lambda = atan2(tn[1], tn[0]) / f;
return [
lambda,
2 * atan(exp(0.5 / f * log(k_ * tn[0] * tn[0] + k_ * tn[1] * tn[1]))) - halfPi
];
};
var guyou = function() {
return d3Geo.geoProjection(squareRaw(guyouRaw))
.scale(151.496);
};
function hammerRaw(A, B) {
if (arguments.length < 2) B = A;
if (B === 1) return d3Geo.geoAzimuthalEqualAreaRaw;
if (B === Infinity) return hammerQuarticAuthalicRaw;
function forward(lambda, phi) {
var coordinates = d3Geo.geoAzimuthalEqualAreaRaw(lambda / B, phi);
coordinates[0] *= A;
return coordinates;
}
forward.invert = function(x, y) {
var coordinates = d3Geo.geoAzimuthalEqualAreaRaw.invert(x / A, y);
coordinates[0] *= B;
return coordinates;
};
return forward;
}
function hammerQuarticAuthalicRaw(lambda, phi) {
return [
lambda * cos(phi) / cos(phi /= 2),
2 * sin(phi)
];
}
hammerQuarticAuthalicRaw.invert = function(x, y) {
var phi = 2 * asin(y / 2);
return [
x * cos(phi / 2) / cos(phi),
phi
];
};
var hammer = function() {
var B = 2,
m = d3Geo.geoProjectionMutator(hammerRaw),
p = m(B);
p.coefficient = function(_) {
if (!arguments.length) return B;
return m(B = +_);
};
return p
.scale(169.529);
};
function hammerRetroazimuthalRaw(phi0) {
var sinPhi0 = sin(phi0),
cosPhi0 = cos(phi0),
rotate = hammerRetroazimuthalRotation(phi0);
rotate.invert = hammerRetroazimuthalRotation(-phi0);
function forward(lambda, phi) {
var p = rotate(lambda, phi);
lambda = p[0], phi = p[1];
var sinPhi = sin(phi),
cosPhi = cos(phi),
cosLambda = cos(lambda),
z = acos(sinPhi0 * sinPhi + cosPhi0 * cosPhi * cosLambda),
sinz = sin(z),
K = abs(sinz) > epsilon ? z / sinz : 1;
return [
K * cosPhi0 * sin(lambda),
(abs(lambda) > halfPi ? K : -K) // rotate for back hemisphere
* (sinPhi0 * cosPhi - cosPhi0 * sinPhi * cosLambda)
];
}
forward.invert = function(x, y) {
var rho = sqrt(x * x + y * y),
sinz = -sin(rho),
cosz = cos(rho),
a = rho * cosz,
b = -y * sinz,
c = rho * sinPhi0,
d = sqrt(a * a + b * b - c * c),
phi = atan2(a * c + b * d, b * c - a * d),
lambda = (rho > halfPi ? -1 : 1) * atan2(x * sinz, rho * cos(phi) * cosz + y * sin(phi) * sinz);
return rotate.invert(lambda, phi);
};
return forward;
}
// Latitudinal rotation by phi0.
// Temporary hack until D3 supports arbitrary small-circle clipping origins.
function hammerRetroazimuthalRotation(phi0) {
var sinPhi0 = sin(phi0),
cosPhi0 = cos(phi0);
return function(lambda, phi) {
var cosPhi = cos(phi),
x = cos(lambda) * cosPhi,
y = sin(lambda) * cosPhi,
z = sin(phi);
return [
atan2(y, x * cosPhi0 - z * sinPhi0),
asin(z * cosPhi0 + x * sinPhi0)
];
};
}
var hammerRetroazimuthal = function() {
var phi0 = 0,
m = d3Geo.geoProjectionMutator(hammerRetroazimuthalRaw),
p = m(phi0),
rotate_ = p.rotate,
stream_ = p.stream,
circle = d3Geo.geoCircle();
p.parallel = function(_) {
if (!arguments.length) return phi0 * degrees;
var r = p.rotate();
return m(phi0 = _ * radians).rotate(r);
};
// Temporary hack; see hammerRetroazimuthalRotation.
p.rotate = function(_) {
if (!arguments.length) return (_ = rotate_.call(p), _[1] += phi0 * degrees, _);
rotate_.call(p, [_[0], _[1] - phi0 * degrees]);
circle.center([-_[0], -_[1]]);
return p;
};
p.stream = function(stream) {
stream = stream_(stream);
stream.sphere = function() {
stream.polygonStart();
var epsilon$$1 = 1e-2,
ring = circle.radius(90 - epsilon$$1)().coordinates[0],
n = ring.length - 1,
i = -1,
p;
stream.lineStart();
while (++i < n) stream.point((p = ring[i])[0], p[1]);
stream.lineEnd();
ring = circle.radius(90 + epsilon$$1)().coordinates[0];
n = ring.length - 1;
stream.lineStart();
while (--i >= 0) stream.point((p = ring[i])[0], p[1]);
stream.lineEnd();
stream.polygonEnd();
};
return stream;
};
return p
.scale(79.4187)
.parallel(45)
.clipAngle(180 - 1e-3);
};
var healpixParallel = 41 + 48 / 36 + 37 / 3600;
var healpixLambert = cylindricalEqualAreaRaw(0);
function healpixRaw(H) {
var phi0 = healpixParallel * radians,
dx = collignonRaw(pi, phi0)[0] - collignonRaw(-pi, phi0)[0],
y0 = healpixLambert(0, phi0)[1],
y1 = collignonRaw(0, phi0)[1],
dy1 = sqrtPi - y1,
k = tau / H,
w = 4 / tau,
h = y0 + (dy1 * dy1 * 4) / tau;
function forward(lambda, phi) {
var point,
phi2 = abs(phi);
if (phi2 > phi0) {
var i = min(H - 1, max(0, floor((lambda + pi) / k)));
lambda += pi * (H - 1) / H - i * k;
point = collignonRaw(lambda, phi2);
point[0] = point[0] * tau / dx - tau * (H - 1) / (2 * H) + i * tau / H;
point[1] = y0 + (point[1] - y1) * 4 * dy1 / tau;
if (phi < 0) point[1] = -point[1];
} else {
point = healpixLambert(lambda, phi);
}
point[0] *= w, point[1] /= h;
return point;
}
forward.invert = function(x, y) {
x /= w, y *= h;
var y2 = abs(y);
if (y2 > y0) {
var i = min(H - 1, max(0, floor((x + pi) / k)));
x = (x + pi * (H - 1) / H - i * k) * dx / tau;
var point = collignonRaw.invert(x, 0.25 * (y2 - y0) * tau / dy1 + y1);
point[0] -= pi * (H - 1) / H - i * k;
if (y < 0) point[1] = -point[1];
return point;
}
return healpixLambert.invert(x, y);
};
return forward;
}
function sphere(step) {
return {
type: "Polygon",
coordinates: [
d3Array.range(-180, 180 + step / 2, step).map(function(x, i) { return [x, i & 1 ? 90 - 1e-6 : healpixParallel]; })
.concat(d3Array.range(180, -180 - step / 2, -step).map(function(x, i) { return [x, i & 1 ? -90 + 1e-6 : -healpixParallel]; }))
]
};
}
var healpix = function() {
var H = 4,
m = d3Geo.geoProjectionMutator(healpixRaw),
p = m(H),
stream_ = p.stream;
p.lobes = function(_) {
return arguments.length ? m(H = +_) : H;
};
p.stream = function(stream) {
var rotate = p.rotate(),
rotateStream = stream_(stream),
sphereStream = (p.rotate([0, 0]), stream_(stream));
p.rotate(rotate);
rotateStream.sphere = function() { d3Geo.geoStream(sphere(180 / H), sphereStream); };
return rotateStream;
};
return p
.scale(239.75);
};
function hillRaw(K) {
var L = 1 + K,
sinBt = sin(1 / L),
Bt = asin(sinBt),
A = 2 * sqrt(pi / (B = pi + 4 * Bt * L)),
B,
rho0 = 0.5 * A * (L + sqrt(K * (2 + K))),
K2 = K * K,
L2 = L * L;
function forward(lambda, phi) {
var t = 1 - sin(phi),
rho,
omega;
if (t && t < 2) {
var theta = halfPi - phi, i = 25, delta;
do {
var sinTheta = sin(theta),
cosTheta = cos(theta),
Bt_Bt1 = Bt + atan2(sinTheta, L - cosTheta),
C = 1 + L2 - 2 * L * cosTheta;
theta -= delta = (theta - K2 * Bt - L * sinTheta + C * Bt_Bt1 -0.5 * t * B) / (2 * L * sinTheta * Bt_Bt1);
} while (abs(delta) > epsilon2 && --i > 0);
rho = A * sqrt(C);
omega = lambda * Bt_Bt1 / pi;
} else {
rho = A * (K + t);
omega = lambda * Bt / pi;
}
return [
rho * sin(omega),
rho0 - rho * cos(omega)
];
}
forward.invert = function(x, y) {
var rho2 = x * x + (y -= rho0) * y,
cosTheta = (1 + L2 - rho2 / (A * A)) / (2 * L),
theta = acos(cosTheta),
sinTheta = sin(theta),
Bt_Bt1 = Bt + atan2(sinTheta, L - cosTheta);
return [
asin(x / sqrt(rho2)) * pi / Bt_Bt1,
asin(1 - 2 * (theta - K2 * Bt - L * sinTheta + (1 + L2 - 2 * L * cosTheta) * Bt_Bt1) / B)
];
};
return forward;
}
var hill = function