Disappointment

[quadibloc]

In doing a Google search for items related to Athelstan Spilhaus, one image result I obtained was of an interesting-looking map projection.
Following the link, I found a paper about equal-area map projections that were kinder to the polar regions:
https://www.jstor.org/stable/986728
I did not immediately log in to access the paper, instead I tried thinking about the matter, to figure out how one might devise a projection that resembled what I saw and was equal-area. But I eventually decided that there was no obvious answer.
Then I did read the paper... and I found that the image I saw was of his "Equidistant" projection, and none of the equal-area projections he presented looked quite so nice as that one!

[PeteD]

Is this his equidistant projection?

spilhaus.jpg (12.16 KiB) Viewed 1722 times

Apart from Australia, I like it. Any chance you could post its formula or a quick summary of its derivation for those of us who don't have a JSTOR account, please?

[quadibloc]

I'll note simply that a JSTOR account to view papers like that is available to anyone free of charge.
I liked it too, apart from Australia, but since it wasn't equal area, I hadn't seen a point. I'll have to look at the paper again carefully; it appeared, on my initial superficial examination, that the details were expressed in a manner that was difficult for me to understand.
The corresponding equal-area projection comes to a point at both ends of the Equator. I now see the reason why.
The basis of the projection was derived by this principle:
Start with the Werner Cardioid projection. Modify it so that an area around the pole, to a target latitude, is on the polar case of the Azimuthal Equidistant projection. Below that latitude, the lengths of the parallels shorten, so that the space between them doesn't have to shrink as rapidly.
If you take the hemisphere from North Pole to South Pole of this projection with the central meridian in the middle, the boundary is a straight line in the azimuthal portion, so two hemispheres in this projection can be abutted.
Transverse the result, and you get his equal-area projection.
Making it equidistant instead of equal-area gets rid of the pointy ends.
But the formulae used are not presented in a way I find clear.

[PeteD]

I'll note simply that a JSTOR account to view papers like that is available to anyone free of charge.

Thanks for pointing that out. I saw "Log in through your school or library" and assumed that was the only way to log in.

I'll have to look at the paper again carefully; it appeared, on my initial superficial examination, that the details were expressed in a manner that was difficult for me to understand.

Six weeks later, I've finally got round to having a proper look at the paper, and I also found some parts of it rather unclear. The part about the equal-area projection isn't too hard to follow, and I've managed to reproduce it as it appears on plate VI, i.e. without the extra interruptions along coastlines:

spilhaus_equal_area.png (236.52 KiB) Viewed 1665 times

In case my code is easier to follow than the equations in the paper, I've pasted it below. The equations in the paper describe the polar aspect of the projection with a single cusp and a single point as shown in figures 2 and 3, so I implemented that as a first step:

function spilhausEqualAreaPolarRaw(lambda, phi) {
var p = halfPi - phi,
p0 = 115 * radians,
CON = 1.145;

if (p <= p0) {
var rho = sqrt(2 * (1-cos(p))),
theta = lambda;
} else {
var rho0 = sqrt(2 * (1-cos(p0))),
drho0 = sin(p0) / rho0,
d2rho0 = 1/rho0 * (cos(p0) - pow(sin(p0)/rho0, 2)),
dp = p - p0,
rho = rho0 + drho0 * dp + d2rho0/2 * pow(dp, 2) + CON/6 * pow(dp, 3),
drho = drho0 + d2rho0 * dp + CON/2 * pow(dp, 2),
theta = sin(p)/(rho*drho) * lambda;
}

var x = rho * sin(theta),
y = -rho * cos(theta);

return[x, y];
}

then as a second step, I converted to the version in equatorial aspect with two cusps and two points as shown on plate VI:

function spilhausEqualAreaEquatorialRaw(lambda, phi) {
var p = acos(cos(phi)*cos(lambda)),
p0 = 115 * radians,
CON = 1.145;

if (p <= p0) {
var rho = sqrt(2 * (1-cos(p))),
theta = atan2(sin(lambda), -tan(phi));
} else {
var rho0 = sqrt(2 * (1-cos(p0))),
drho0 = sin(p0) / rho0,
d2rho0 = 1/rho0 * (cos(p0) - pow(sin(p0)/rho0, 2)),
dp = p - p0,
rho = rho0 + drho0 * dp + d2rho0/2 * pow(dp, 2) + CON/6 * pow(dp, 3),
drho = drho0 + d2rho0 * dp + CON/2 * pow(dp, 2),
theta = sin(p)/(rho*drho) * (atan2(sin(lambda), -tan(phi)) - sign(lambda)*halfPi) + sign(lambda)*halfPi;
}

var x = rho * sin(theta),
y = -rho * cos(theta);

return[x, y];
}

The corresponding equal-area projection comes to a point at both ends of the Equator.

You can make it a lot less pointy by reducing the value of the parameter CON. This also straightens out the weird squiggles in the parallels. Here it is with CON = 0:

spilhaus_equal_area_CON=0.png (315.34 KiB) Viewed 1665 times

Then I did read the paper... and I found that the image I saw was of his "Equidistant" projection, and none of the equal-area projections he presented looked quite so nice as that one!

Personally, I think equidistant projections generally look nicer than otherwise similar equal-area ones, but if it was the points and/or the squiggles that were putting you off his equal-area projection, now you know how to get rid of them.

Speaking of his equidistant projection, this is where I'm obviously misunderstanding something in the paper because when I try to plot it, I get this:

spilhaus_equidistant.png (232.31 KiB) Viewed 1665 times

I've posted my code below in the hope that someone can see what I'm doing wrong. I've only posted the code for the version in polar aspect with a single cusp and a single point because it corresponds to the equations in the paper and because the process of converting it to the version in equatorial aspect with two cusps and two points should be the same as for the equal-area projection, so I don't think that's where the mistake is.

function depressedCubic(p, q) { // finds a real root of a cubic equation of the form x^3 + p*x + q = 0
var Delta = -(4 * pow(p, 3) + 27 * pow(q, 2));

if (Delta < 0) { // 1 real root, 2 complex roots
var C = pow(-q/2 + sqrt(pow(q/2, 2) + pow(p/3, 3)), 1/3),
x = C - p/(3*C); // real root
} else if (Delta > 0) { // 3 real roots
var k = 1, // choose which root (0, 1 or 2)
x = 2 * sqrt(-p/3) * cos(1/3 * acos(3*q/(2*p) * sqrt(-3/p)) - 2*pi/3 * k);
}

return x;
};

function spilhausEquidistantPolarRaw(lambda, phi) {
var P_D = halfPi - phi,
p0 = 115 * radians,
rho0 = sqrt(2 * (1-cos(p0))),
K = rho0 / p0,
rho = K * P_D;

if (P_D <= p0) {
var theta = lambda;
} else {
var drho0 = sin(p0) / rho0,
d2rho0 = 1/rho0 * (cos(p0) - pow(sin(p0)/rho0, 2)),
CON = 6/pow(pi-p0, 3) * ((rho0*pi)/p0 - rho0 - drho0 * (pi-p0) - d2rho0/2 * pow(pi-p0, 2));

// terms of cubic equation a*dp^3 + b*dp^2 + c*dp + d = 0

var a = CON/6,
b = d2rho0/2,
c = drho0,
d = rho0-rho;

// convert to depressed cubic t^3 + r*t + s = 0 and solve for t

var r = (3*a*c - pow(b, 2)) / (3*pow(a, 2)),
s = (2*pow(b, 3) - 9*a*b*c + 27*pow(a, 2)*d) / (27*pow(a, 3)),
t = depressedCubic(r, s);

// convert back to dp in order to calculate p and hence theta

var dp = t - b/(3*a),
p = p0 + dp,
theta = sin(p)/(rho*drho0) * lambda;

// check dp is really a root of cubic equation

console.log(a*pow(dp, 3) + b*pow(dp, 2) + c*dp + d)
}

var x = rho * sin(theta),
y = -rho * cos(theta);

return[x, y];
}

I liked it too, apart from Australia, but since it wasn't equal area, I hadn't seen a point.

Why is there no point in equidistant projections?

[PeteD]

I still haven't managed to work what's wrong with my implementation of Spilhaus's equidistant projection, but the more I think about it, the more I think that Spilhaus made it unnecessarily complicated.

In polar aspect, for values of the colatitude p less than a threshold value p0, it's just the equidistant azimuthal projection, which in polar coordinates is normally given by (rho, theta) = (p, lambda). Spilhaus scales rho by a constant factor K, but this doesn't really matter as it only affects the nominal scale of the projection, not its appearance.

For p > p0, rho must still be equal to p (or K*p) for the projection to be equidistant, but theta starts to decrease along each meridian and becomes zero at the boundary so that the projection tapers to a point. It does this in such a way that dtheta/dp = 0 at p = p0 in order to avoid a kink in the meridians at the boundary between the two regions.

This last point is where my attempt at implementing Spilhaus's formula clearly fails. However, this criterion is easy to fulfil, and as far as I can tell, the projection isn't supposed to fulfil any further criteria. The most obvious way is simply to specify that theta be a parabolic, elliptic or sinusoidal function of p, then converting to the version in equatorial aspect with two cusps and two points gives the following three projections:

spilhaus_denner_parabolic.png (228.91 KiB) Viewed 1651 times

spilhaus_denner_elliptic.png (238.16 KiB) Viewed 1651 times

spilhaus_denner_sinusoidal.png (227.31 KiB) Viewed 1651 times

Although perhaps the increased complexity of Spilhaus's formula is justified on the grounds that his projection is slightly better than these three -- it avoids both the pointy sides of my parabolic and sinusoidal versions and the large inflation towards the edge of my elliptic version.

[Milo]

I have to wonder about the point of these projections. If I just run the azimuthal equidistant projection with the same center point, I get a projection that looks pretty much the same. The only place that really looks noticeably different is New Zealand, and it's not like any of these projections actually do a good job there. Continents like North America and Australia are visibly distorted from their proper shapes, but with pretty much the same distortions either way.

Do these projections have any actual advantage that justifies their more complicated math?

I liked it too, apart from Australia, but since it wasn't equal area, I hadn't seen a point.

Why is there no point in equidistant projections?

I do consider the azimuthal equidistant and cylindrical equidistant projections to be useful, because they work well as "compromise projections" between equal-area and conformal properties, while also being very simple and easy to work with.

However, the equidistant property by itself is not too useful, especially if you're not making it clear which distances are being preserved (since no projection can preserve all distances). If you do care about measuring some specific distances, then the azimuthal equidistant or two-point equidistant projection is probably the one you want.

[quadibloc]

Personally, I think equidistant projections generally look nicer than otherwise similar equal-area ones,

Well, of course they do, and conformal projections look nicer yet, in general, but sometimes the equal-area property is needed.
However, that is only a general rule, and the Orthographic projection, which is in the direction from equal-area projections away from equidistant and conformal ones, looks nicer yet.
So it is nice to know how to get rid of the points and squiggles.

Why is there no point in equidistant projections?

I hadn't meant that in so general a sense. The azimuthal equidistant projection is certainly useful and valuable. So, for that matter, is the simple conic, for example.