A dihedral projection

[Milo]

Ah, if we say trilinear coordinates, it will be interpreted as that kind of normalization. Not a good name...

Name isn't too important, just explain the math behind what you mean

Your original post said "based on the Chamberlin trimetric, but...", which I interpreted as meaning that you're doing something which is similar to the Chamberlin trimetric but not identical (a description that would indeed apply to trilinear coordinates as defined by Wikipedia).

My intention was Chamberlin's way. The three lines do not intersect at a single point; therefore, find the center of the triangle formed by them.

Center defined how? There are many different ways of determining the center of a non-equilateral triangle, further complicated by the fact that we're properly talking about a circular triangle (convex type) rather than the normal kind.

The way I used trilinear coordinates above will always produce some point lying inside this triangle, and so can be thought of as a type of center, and in any case won't be too far off.

Presumably, if you can define the center in just the right way, you won't need the second tweaking step.

Inverse functions of added things are difficult (recalling Mollweide's equation).

Tedious to calculate in practice, but not difficult to define in the mathematical sense. Sometimes it's just the function you want.

Figuring out if this method actually produces the right results would be a good start, optimizing the calculation can come later.

• Stereographic = exp(Mercator), Mercator = log(Stereographic) (where Mercator is directed N as plus, Stereographic is south-polar view). Translating the Mercator projection is equivalent to scaling the Stereographic projection.
• In stereographic projection, great-circle through the tangent point is projected to a line. scaling the projected image results the line with same appearance, different mapping. Then inverse stereographic results same appearance different scale great-circle.

Yup, you're completely right! I hadn't thought to try a stereographic projection whose center lies on the edge of the hemisphere.

As may be well known, it is possible to convert between a hemisphere and a octant without elliptic integrals.

I didn't know that for a fact, but I'm not at all surprised that it is. What would surprise me greatly was if such a conversion were conformal.

The conversion dummy_index describes is a composition of two conformal transformations, and so is still conformal.

Also, the Riemann mapping theorem says, essentially, that any shape can be conformally mapped onto any other shape (terms and conditions apply), so the existence of a conformal mapping between the hemisphere and octant should be not at all surprising, although the fact that it doesn't need elliptic integrals might be.

[dummy_index]

Center defined how? There are many different ways of determining the center of a non-equilateral triangle, further complicated by the fact that we're properly talking about a circular triangle (convex type) rather than the normal kind.

Sorry, I don't use arcs in the x-y plane at all.

(The story is much the same whether the distance is from each great circle or from the vertex.)
a = latitude
b = asin(cos(latitude)*sin(longitude))
c = asin(cos(latitude)*cos(longitude))
d = a + b + c - pi/2
∴(a - d/3) + (b - d/3) + (c - d/3) = pi/2
x = sqrt(1/3)*((b - d/3)-(c - d/3)) = sqrt(1/3)*(b-c)
y = (a - d/3)

(Ternary plot is the better image. It is not appropriate to borrow Chamberlin's name because it is much simpler...)

Presumably, if you can define the center in just the right way, you won't need the second tweaking step.

I did this because I thought about using the LatTLat as a starting point to keep orthogonality at the edges, and it was obvious that a method like a/(a+b+c) would break the orthogonality...

[Milo]

Ah, I see. "Normalizing" additively instead of multiplicatively.

Meanwhile, my own attempt to find a formula to "flatten the apothem" is failing miserably. Trying to find a match for these formulae:
d = a + b + c - pi/2
x = sqrt(1/3) * (f(b-d/3) - f(c-d/3))
y = f(a-d/3) = f((2*a - b - c + pi/2) / 3)
doesn't work because, while (a-d/3) + (b-d/3) + (c-d/3) is constant, there's no guarantee that (a-d/3) + (b-d/3) + (c-d/3) will be. To solve this, the tweaking function can be applied to a/b/c before you subtract d/3, like this:
d = f(a) + f(b) + f(c) - pi/2
x = sqrt(1/3) * (f(b) - f(c))
y = f(a) - d/3 = (2*f(a) - f(b) - f(c) + pi/2) / 3
...But when you do this, it's far too hard to figure out any meaningful characterization of f (even in terms of its inverse function), even in the special case where b = c = asin(cos(a)/sqrt(2)).

In fact, even the inverse formulae for the basic projection is hard to figure out, because when you use:
x = sqrt(1/3) * (b - c)
y = (a - d/3)
there's no easy way to determine d/3 from x and y. If you could, it would then be easy enough to compute that:
a = y + d/3
b = (pi/2 + sqrt(3)*x - y) / 2 + d/3
c = (pi/2 - sqrt(3)*x - y) / 2 + d/3
latitude = y + d/3
longitude = atan2(sin(b), sin(c)) (or any of several different approaches)

What I can observe, though, is that:

The center of the octant lies at latitude = a = b = c = acos(1/3)/2 ≈ 35.264°.

The center of the planar triangle, though, lies at y = pi/6 = 30°.

Therefore, for resolution-efficiency, the line segment at longitude = pi/4, latitude ∈ [0,acos(1/3)/2] needs to be mapped linearly onto x = 0, y = latitude * pi/3 / acos(1/3). However, the remaining line segment at longitude = pi/4, latitude ∈ [acos(1/3)/2,pi/2] then needs to be mapped differently in order to "catch up" so y = pi/2 when latitude = pi/2. So you need a piecewise function, no matter what.

In theory, it is possible to do this in an infinitely-differentiable way, if you use a flat function, but it still won't be analytic at the center.

[dummy_index]

d = f(a) + f(b) + f(c) - pi/2
x = sqrt(1/3) * (f(b) - f(c))
y = f(a) - d/3 = (2*f(a) - f(b) - f(c) + pi/2) / 3
...But when you do this, it's far too hard to figure out any meaningful characterization of f (even in terms of its inverse function), even in the special case where b = c = asin(cos(a)/sqrt(2)).

(snip)

Therefore, for resolution-efficiency, the line segment at longitude = pi/4, latitude ∈ [0,acos(1/3)/2] needs to be mapped linearly onto x = 0, y = latitude * pi/3 / acos(1/3). However, the remaining line segment at longitude = pi/4, latitude ∈ [acos(1/3)/2,pi/2] then needs to be mapped differently in order to "catch up" so y = pi/2 when longitude = pi/2. So you need a piecewise function, no matter what.

Indeed. Then, where b = c = asin(cos(a)/sqrt(2)) and a ∈ [0,acos(1/3)/2], b = c ∈ [acos(1/3)/2,pi/4]. No overlap. And if we tweak in [acos(1/3)/2,pi/4], it only affect inside the inverted triangle. As you can see, everywhere near the equator needs to be stretched with the goal of improving resolution-efficiency. Hence, f(a) must be tweaked.

Trilinear_MinAxis.png (64.56 KiB) Viewed 7877 times

If f(b) and f(c) can be known within this interval, then

f(b) = f(c) = f(asin(cos(a)/sqrt(2)))
y = f(a) - d/3 = (2*f(a) - f(b) - f(c) + pi/2) / 3 = (2*f(a) - 2*f(asin(cos(a)/sqrt(2))) + pi/2) / 3

It is sufficient if f(a) cancels out f(asin(cos(a)/sqrt(2))) so that a linear term remains. That is, EDIT
f(a) = constant*a + f(asin(cos(a)/sqrt(2))) - pi/4 when a ∈ [0,acos(1/3)/2].

[Milo]

f(a) = constant*a - f(asin(cos(a)/sqrt(2)))

I believe you mean plus, not minus? Or f(a) - f(asin(cos(a)/sqrt(2))) = constant*a.

Other than that, yeah, makes sense.

It's not quite the same thing as Mollweide's equation, which is more of the form f(y) + g(f(y)) = h(x), solved for y given x (with f/g/h all known functions), rather than f(x) + f(g(x)) = h(x), solved for f(x) given x (with g/h known functions but f not). I know how to solve the former kind of equation numerically, or plot it by calculating the closed-form inverse function (solve for x given y) and mirroring it across the x=y axis, but I don't see how to do either of those for the latter kind of equation.

[dummy_index]

f(a) = constant*a - f(asin(cos(a)/sqrt(2)))

I believe you mean plus, not minus? Or f(a) - f(asin(cos(a)/sqrt(2))) = constant*a.

My apologies. Corrected.

And calm down. It was meant to be a definition, rather than an equation...

f(x): tweaking function,
f(x) := f_1(x) when x ∈ [0,acos(1/3)/2],
f(x) := f_2(x) when x ∈ [acos(1/3)/2,pi/2-acos(1/3)/2],
f(x) := f_3(x) when x ∈ [pi/2-acos(1/3)/2,pi/2],
f_1(0) = 0,
f_3(pi/2) = pi/2,
f_2(pi/4) = pi/4,
f_2(x)+f_2(pi/2-x) = pi/2 (for straight edges as the basic projection),
f_1(x)+f_3(pi/2-x) = pi/2 (for straight edges as the basic projection),
f_1(acos(1/3)/2) = f_2(acos(1/3)/2) (for C^0 condition),
f_2(pi/2-acos(1/3)/2) = f_3(pi/2-acos(1/3)/2) (for C^0 condition),
f_2 is given in some way. (simplest case: f_2(x) := x).
y = f(a) - d/3 = (2*f(a) - f(b) - f(c) + pi/2) / 3.

In the case of b = c = asin(cos(a)/sqrt(2)) and a ∈ [0,acos(1/3)/2]:
y = (2*f_1(a) - f_2(b) - f_2(c) + pi/2) / 3
= (2*f_1(a) - 2*f_2(asin(cos(a)/sqrt(2))) + pi/2) / 3
And we want y = a * (pi/3) / acos(1/3).
(2*f_1(a) - 2*f_2(asin(cos(a)/sqrt(2))) + pi/2) / 3 = a * (pi/3) / acos(1/3)
2*f_1(a) - 2*f_2(asin(cos(a)/sqrt(2))) + pi/2 = a * pi / acos(1/3)
2*f_1(a) = a * pi / acos(1/3) + 2*f_2(asin(cos(a)/sqrt(2))) - pi/2
f_1(a) = a * pi / (2*acos(1/3)) + f_2(asin(cos(a)/sqrt(2))) - pi/4
(where f_2(***) - pi/4 = 0 if a = 0)

Thus, given f_2, f_1 is uniquely determined.

[Milo]

That took some effort to unpack, but I think I see the gist of it. The key observation is that when a ∈ [0,acos(1/3)/2], then asin(cos(a)/sqrt(2)) ∈ [acos(1/3)/2,pi/4] always, so the resolution-efficiency constraint only ever compares f_1 against f_2, rather than needing to compare f_1 against itself evaluated at a different point.

f_2(pi/4) = pi/4 follows logically from f_1(0) = 0 and the resolution-efficiency constraint. I had a bit more trouble with f_2(x)+f_2(pi/2-x) = pi/2 and f_1(x)+f_3(pi/2-x) = pi/2, but then I figured out from the behavior of Euclidean ternary plots that if you assume the untweaked projection satisfies perpendicularity at the edges (which I still haven't worked out a theoretical proof for, but is empirically true), then this behavior ensures the tweaked projection will still do so.

The constraint becomes:
f_1(a) = a * pi/2 / acos(1/3) + f_2(asin(cos(a)/sqrt(2))) - pi/4
We can also derive this to get:
f_1'(a) = pi/2 / acos(1/3) - sin(a) / sqrt(2 - cos(a)^2) * f_2'(asin(cos(a)/sqrt(2)))
Plugging in a = acos(1/3)/2 gives:
f_1(acos(1/3)/2) = f_2(acos(1/3)/2)
f_1'(acos(1/3)/2) = pi/2 / acos(1/3) - 1/2 * f_2'(acos(1/3)/2)
The first is something we obviously wanted anyway. The second, though, can be combined with the requirement that f_1'(acos(1/3)/2) = f_2'(acos(1/3)/2) (for smoothness), to give:
f_1'(acos(1/3)/2) = f_2'(acos(1/3)/2) = pi/3 / acos(1/3)
Alternatively, plugging in a = 0 gives:
f_1'(0) = pi/2 / acos(1/3)
So, f_1 at least can't be linear. We don't have a constraint for f_2'(pi/4), though, because it gets multiplied by zero.

f_2 could conceivably be linear, but then it would have to be:
f_2(a) = (a - pi/4) * pi/3 / acos(1/3) + pi/4
In which case:
f_1(a) = (a/2 - pi/12 + asin(cos(a)/sqrt(2))/3) * pi / acos(1/3)
f_1(acos(1/3)/2) = f_2(acos(1/3)/2) = (5 - pi/acos(1/3)) * pi/12
f_1''(acos(1/3)/2) = -pi/acos(1/3)/4/sqrt(2) ≈ -0.45116 ≠ 0 = f_2''(acos(1/3)/2)

plot [a=0:pi/2] [0:pi/2] a <= acos(1.0/3)/2 ? (a/2 - pi/12 + asin(cos(a)/sqrt(2))/3) * pi / acos(1.0/3) : a >= pi/2 - acos(1.0/3)/2 ? pi/2 - (pi/6 - a/2 + asin(sin(a)/sqrt(2))/3) * pi / acos(1.0/3) : (a - pi/4) * pi/3 / acos(1.0/3) + pi/4

tweak.png (2.18 KiB) Viewed 7863 times

First-differentiability is enough to ensure a smooth appearance, but we could push farther and attempt to derive constraints on the second derivative too. The complete formula is messy, but if I have Maxima work it out and then plug in the relevant boundary conditions, then, assuming I didn't make any mistakes (it's pretty likely I did), I arrive at:
f_1''(0) = -f_2'(pi/4)
f_1''(acos(1/3)/2) = -3/4/sqrt(2) * f_2'(acos(1/3)/2) + 1/4 * f_2''(acos(1/3)/2)
f_1''(acos(1/3)/2) = f_2''(acos(1/3)/2) = -1/sqrt(2) * f_2'(acos(1/3)/2) = -pi/acos(1/3)/3/sqrt(2)

[dummy_index]

I never imagined that differential coefficients could be so reusable. Great!

Class C^1 (First-differentiability): confirmed.
f_2(x) = (x-pi/4) * pi/3/acos(1.0/3) + pi/4

Class C^2: confirmed.
I now assume
f_2(x) = ( x+g*(3*(x-pi/4)-(x-pi/4)**3/(pi/4-acos(1.0/3)/2)**2) - pi/4 )*pi/3/acos(1.0/3)+pi/4
and
f_2''(acos(1/3)/2) = -pi/acos(1/3)/3/sqrt(2)
therefore
g = -(pi/4-acos(1.0/3)/2)/6/2**0.5.
(When the second-order central differences are plotted and compared, my hand-adjusted result looks a little better. Nevertheless, it is a relief to have the evidence.)

[Milo]

I arrive at:
f_2(a) = (a-pi/4)*u + (a-pi/4)^3*v + pi/4
f_2'(a) = u + 3*(a-pi/4)^2*v
f_2''(a) = 6*(a-pi/4)*v
f_2'(acos(1/3)/2) = u + 3*(acos(1/3)/2-pi/4)^2*v = pi/3 / acos(1/3)
f_2''(acos(1/3)/2) = 6*(acos(1/3)/2-pi/4)*v = -pi/acos(1/3)/3/sqrt(2)
u = pi/3 * ((1-pi/8/sqrt(2))/acos(1/3) + 1/4/sqrt(2)) = pi/3/acos(1/3) * (1 + (acos(1/3)-pi/2)/4/sqrt(2))
v = -pi/9*sqrt(2) / acos(1/3) / (2*acos(1/3)-pi)
I don't even want to figure out f_1 from that.

Plotting suggests this is the same f_2 as yours, so I guess we did it right.

In theory, it should be possible to inductively define a C^∞ Taylor series, but you'd have to find some way around the coefficients becoming exponentially more complicated with each iteration.

[dummy_index]

Please wait a while as I would like to derive the C3 class solution on my own. From numerical observation, it seems that the coefficient of the third-order term at the endpoints (actually the fifth-order term due to symmetry) can be left as is or halved to achieve optimal resolution-efficiency.

g3 = -(pi/4-acos(1.0/3)/2)/6/2**0.5
g5 = 0.00276
f_2(x) = ( x+g3*(3*(x-pi/4)-(x-pi/4)**3/(pi/4-acos(1.0/3)/2)**2)+g5*(5*(x-pi/4)-10.0/3*(x-pi/4)**3/(pi/4-acos(1.0/3)/2)**2+(x-pi/4)**5/(pi/4-acos(1.0/3)/2)**4)) - pi/4 )*pi/3/acos(1.0/3)+pi/4

By the way, I would like to change the name of the basic projection. In addition, I will modify the scale of my equation to match Milo's equation (so that the scale at the center point of the equilateral triangle is 1).
Name: Isometric ternary projection (octant)
Description: Three orthogonal great circles are represented by an equilateral triangle. If the distances from a point inside the equilateral triangle to each side are a', b', and c', then a'-b', b'-c', and c'-a' correctly represent the corresponding a-b, b-c, and c-a on the sphere. (An example of properties emerging from additive normalization. Other suggestions welcome)
Distortion: scale on center is correct. along the edges are constant. Line segments orthogonal to the great circles are mapped orthogonally to the corresponding edges.