Tobler’s hyperelliptical projection

[rschmunk]

My original goal in writing that program was to explore the Tobler hyperelliptical projection better than available tools let me do. NASA's G.Projector, which I usually use, does include the Tobler projection, but does not allow variation of the exponent (presumably, it's hardcoded to use the 2.5 value that Wikipedia lists as Tobler's favorite), which I wanted to experiment with.

I have worked on the general case of Tobler's Hyperelliptical, allowing the user to specify all three parameters, and recently. But the results have been unsatisfactory. Even when the output map comes out looking good, rendering time can be sluggish. Even rendering Tobler's preferred case is a bit slow.

[Milo]

I have worked on the general case of Tobler's Hyperelliptical, allowing the user to specify all three parameters, and recently. But the results have been unsatisfactory. Even when the output map comes out looking good, rendering time can be sluggish. Even rendering Tobler's preferred case is a bit slow.

How sluggish are you talking about? A quick test shows mine can render a fullscreen-ish map in about three-quarters of a second, which is short enough to not cause impatience when I'm rendering one map at a time from the command line, though it could be a problem if you're trying to render maps in real-time in response to a user shifting a slider in a GUI. Then again, I suspect that the majority of that time is spent on the PNG compression algorithm when saving the results (I have it set to maximum compression quality), rather than the actual map-rendering code itself.

[daan]

How sluggish are you talking about? A quick test shows mine can render a fullscreen-ish map in about three-quarters of a second…

That’s not bad, and possibly exceptional, depending on your hardware. My original implementation in the early 90s was unusably slow due to the quadrature. I made it usable by caching fixed values and interpolating. However, with improved skills and tools, I find now that the integral can be expressed in terms of the incomplete Beta function, which I have an implementation for, or the Gaussian hypergeometric series, which I also have an implementation for. Those would yield good speeds without caching and interpolation, and, since I rely on the projection inverse for proper “antialiasing”, I should also replace the generic projection root-finder that I use for inverse in most cases with one properly tailored to the projection.

The projection also has the problem that it’s not monotonically increasing in y with respect to latitude for many choices of parameters. That obliges my implementation to check for that, since bad things happen without a strict-ish inverse.

Cheers,
— daan

[Milo]

That’s not bad, and possibly exceptional, depending on your hardware. My original implementation in the early 90s was unusably slow due to the quadrature. I made it usable by caching fixed values and interpolating. However, with improved skills and tools, I find now that the integral can be expressed in terms of the incomplete Beta function, which I have an implementation for, or the Gaussian hypergeometric series, which I also have an implementation for. Those would yield good speeds without caching and interpolation, and, since I rely on the projection inverse for proper “antialiasing”, I should also replace the generic projection root-finder that I use for inverse in most cases with one properly tailored to the projection.

As I noted in my earlier post on the subject, my implementation does not actually perform any integration at all. Rather, it computes which pixels (or anti-aliased subpixels) do or don't belong to the hyperellipse using the |2x/w|^a + |2y/h|^b ≤ 1 formula, and then counts the pixels.

Technically, computing the hyperellipse is done by a different program as a sort of "preprocessing" step that wasn't included in the measurement above. However, I tested it just now, and it takes a small enough fraction of a second that both steps put together would still take less than a second.

My current implementation only handles the "inverse" projection (which is what you need to render raster data), but in theory the forward projection (which is what you need to render vector data, including graticules) could be handled using the same principle.

The projection also has the problem that it’s not monotonically increasing in y with respect to latitude for many choices of parameters. That obliges my implementation to check for that, since bad things happen without a strict-ish inverse.

I don't understand what you're talking about here?

The relationship between y and latitude is always monotonic for any pseudocylindrical projection, obviously including Tobler (even in silly cases, like when the exponent lies between 0 and 1). What are you claiming isn't monotonic, exactly?

[daan]

I don't understand what you're talking about here?

The relationship between y and latitude is always monotonic for any pseudocylindrical projection, obviously including Tobler (even in silly cases, like when the exponent lies between 0 and 1). What are you claiming isn't monotonic, exactly?

Okay. Thanks for the conversation. That sent me down a rabbit-hole from which I have emerged with a vastly improved (actually rigorous, both forward and inverse) implementation. My 25-year old code and memory were both decrepit. Interpolation, gone. Quadrature, gone. Failed combinations of parameters, gone.

Cheers,
— daan

[Milo]

Glad to help, I guess.

Did you imitate my method, or did I somehow inspire you to come up with a third approach?

I'm also still curious about what method rschmunk (R. Schmunk?) was using...

[daan]

Did you imitate my method, or did I somehow inspire you to come up with a third approach?

The help came in the form of prompting me to actually look at what I had and figure out why I would claim, from ancient memory, that there were incompatible mixes of parameters. The argument, “The relationship between y and latitude is always monotonic for any pseudocylindrical projection,” while true by definition, isn’t always true by formulation. Hufnagel’s projection, for example, has large domains in which the result isn’t a pseudocylindrical projection due to the failure of monotonicity in those domains, but when monoticity does not fail, the result is a pseudocylindrical projection. My ancient memories had me thinking I was having the same problem with Tobler’s hyperelliptic. What I found is that my janky numerical techniques were at fault. A generalized midpoint quadrature just isn’t a good idea here, among other messes in the pile of code.

I go off of the original formulation for both forward and inverse. To solve y, I use a Newton’s iteration. I use Gauss’s hypergeometric series to compute the integral. Seed value for Newton’s iteration can get tricky, though.

— daan

[rschmunk]

I have worked on the general case of Tobler's Hyperelliptical, allowing the user to specify all three parameters, and recently. But the results have been unsatisfactory. Even when the output map comes out looking good, rendering time can be sluggish. Even rendering Tobler's preferred case is a bit slow.

How sluggish are you talking about? A quick test shows mine can render a fullscreen-ish map in about three-quarters of a second, which is short enough to not cause impatience when I'm rendering one map at a time from the command line, though it could be a problem if you're trying to render maps in real-time in response to a user shifting a slider in a GUI. Then again, I suspect that the majority of that time is spent on the PNG compression algorithm when saving the results (I have it set to maximum compression quality), rather than the actual map-rendering code itself.

It varies. Generally ~3 seconds to re-render on screen, but for some parameter choices, it goes all to hell and could be 20 seconds.

But the problem is in the forward equations, used to render coastlines and such. Turn that off, and maps are rendered as soon as you hit return.

[Atarimaster]

Failed combinations of parameters, gone.

Does that mean we won’t get the “Projection parameters are not compatible with each other” message anymore on the hyperelliptical projection?

[mapnerd2022]

Speaking of which, which is your favourite configuration of the Hyperelliptical projection?