Increase MAX_ITER so Mollweide forward projection works near the poles.

[erykoff]

• Closes #xxxx
• Tests added
• Added clear title that can be used to generate release notes
• Fully documented, including updating docs/source/*.rst for new API

As has been reported on cartopy here, and here, the PROJ implementation of the Mollweide projection fails near the poles. In particular, with all versions up to and including 9.0.0, you can see that all declinations north of ~89.15 map to the same x/y point:

$ echo 0.0 89.1 | proj +proj=moll +ellps=GRS80
0.00	9000253.15
$ echo 0.0 89.2 | proj +proj=moll +ellps=GRS80
0.00	9020047.85
$ echo 0.0 89.3 | proj +proj=moll +ellps=GRS80
0.00	9020047.85
$ echo 0.0 89.4 | proj +proj=moll +ellps=GRS80
0.00	9020047.85

It turns out the problem is simply that the MAX_ITER value is too low. Setting MAX_ITER = 30 as in this PR allows the projection to work properly up to (at least) 89.99999. There should be no penalty for performance at other declinations since the convergence tolerance has not been changed.

With this PR:

$ echo 0.0 89.1 | proj +proj=moll +ellps=GRS80
0.00	9000253.15
$ echo 0.0 89.2 | proj +proj=moll +ellps=GRS80
0.00	9003130.48
$ echo 0.0 89.3 | proj +proj=moll +ellps=GRS80
0.00	9005889.97
$ echo 0.0 89.4 | proj +proj=moll +ellps=GRS80
0.00	9008520.64
$ echo 0.0 89.999 | proj +proj=moll +ellps=GRS80
0.00	9020045.57
$ echo 0.0 89.9999 | proj +proj=moll +ellps=GRS80 
0.00	9020047.74
$ echo 0.0 89.99999 | proj +proj=moll +ellps=GRS80
0.00	9020047.84

For further information, setting MAX_ITER = 20 works up to 89.999.

[kbevers]

Seems reasonable to me. Please add your test points to builtins.gie:

PROJ/test/gie/builtins.gie

Lines 3581 to 3612 in 6243d11

	===============================================================================
	# Mollweide
	# PCyl., Sph.
	===============================================================================
	-------------------------------------------------------------------------------
	operation +proj=moll +a=6400000
	-------------------------------------------------------------------------------
	tolerance 0.1 mm
	accept 2 1
	expect 201113.698641813 124066.283433860
	roundtrip 1
	accept 2 -1
	expect 201113.698641813 -124066.283433860
	roundtrip 1
	accept -2 1
	expect -201113.698641813 124066.283433860
	roundtrip 1
	accept -2 -1
	expect -201113.698641813 -124066.283433860
	roundtrip 1
	direction inverse
	accept 200 100
	expect 0.001988738 0.000806005
	accept 200 -100
	expect 0.001988738 -0.000806005
	accept -200 100
	expect -0.001988738 0.000806005
	accept -200 -100
	expect -0.001988738 -0.000806005

Also, is Mollweide undefined at the poles or can latitudes very close to 90 be special cased to return the coordinates at the poles?

[erykoff]

The projection is defined up to +/-90.0, and I added the test to confirm that works. The problem before was that everything within +/-0.85 degrees snapped to the value at the pole. The additional test ensures that +/-89.99 is indeed different from the pole. If you wish I could change this to more 9s, but I thought this would be sufficient.

[erykoff]

I can't figure out from the logs what the test failures are, and building locally seems to work for me.

[kbevers]

I would expect the logs to give more information than basically "test failed". @rouault @mwtoews is this normal behaviour?

[mwtoews]

I'll make a quick PR to use ctest --output-on-failure which will show these errors. For this one, I get it locally:

$ ctest --output-on-failure
Test project /tmp/build
      Start  1: geodesic-test
 1/57 Test  #1: geodesic-test ....................   Passed    0.01 sec
      Start  2: Builtins
 2/57 Test  #2: Builtins .........................***Failed    0.04 sec
proj_create: Error 1027 (Invalid value for an argument): aea: Invalid value for lat_1: |lat_1| should be <= 90°
proj_create: Error 1027 (Invalid value for an argument): aea: Invalid value for lat_2: |lat_2| should be <= 90°
proj_create: Error 1027 (Invalid value for an argument): aea: Invalid value for lat_1 and lat_2: |lat_1 + lat_2| should be > 0
proj_create: Error 1027 (Invalid value for an argument): aea: Invalid value for eccentricity
eqdc: Invalid latitude
proj_create: Error 1027 (Invalid value for an argument): eqdc: Invalid value for lat_1 and lat_2: |lat_1 + lat_2| should be > 0
proj_create: Error 1027 (Invalid value for an argument): eqdc: Invalid value for lat_1: |lat_1| should be <= 90°
proj_create: Error 1027 (Invalid value for an argument): eqdc: Invalid value for lat_2: |lat_2| should be <= 90°
proj_create: Error 1027 (Invalid value for an argument): eqdc: Invalid value for lat_1 and lat_2: |lat_1 + lat_2| should be > 0
proj_create: Error 1027 (Invalid value for an argument): geos: Invalid value for h.
proj_create: Error 1027 (Invalid value for an argument): geos: Invalid value for h.
proj_create: Error 1027 (Invalid value for an argument): krovak: Invalid value for lat_0: lat_0 + PI/4 should be different from 0
proj_create: Error 1027 (Invalid value for an argument): labrd: Invalid value for lat_0: lat_0 should be different from 0
proj_create: Error 1027 (Invalid value for an argument): laea: Invalid value for lat_0: |lat_0| should be <= 90°
proj_create: Error 1027 (Invalid value for an argument): lagrng: Invalid value for W: it should be > 0
proj_create: Error 1027 (Invalid value for an argument): lagrng: Invalid value for lat_1: |lat_1| should be < 90°
proj_create: Error 1027 (Invalid value for an argument): lcc: Invalid value for eccentricity
proj_create: Error 1027 (Invalid value for an argument): lcc: Invalid value for eccentricity
proj_create: Error 1027 (Invalid value for an argument): lcc: Invalid value for lat_2: |lat_2| should be < 90°
proj_create: Error 1027 (Invalid value for an argument): lcc: Invalid value for lat_1: |lat_1| should be < 90°
proj_create: Error 1027 (Invalid value for an argument): lcc: Invalid value for lat_1: |lat_1| should be < 90°
proj_create: Error 1027 (Invalid value for an argument): lcc: Invalid value for lat_2: |lat_2| should be < 90°
proj_create: Error 1027 (Invalid value for an argument): lcc: Invalid value for lat_1: |lat_1| should be < 90°
proj_create: Error 1027 (Invalid value for an argument): lcc: Invalid value for lat_1: |lat_1| should be < 90°
proj_create: Error 1027 (Invalid value for an argument): lcc: Invalid value for lat_2: |lat_2| should be < 90°
proj_create: Error 1027 (Invalid value for an argument): nsper: Invalid value for h
proj_create: Error 1027 (Invalid value for an argument): nsper: Invalid value for h
proj_create: Error 1026 (Missing argument): ob_tran: Failed to find projection to be rotated
proj_create: Error 1027 (Invalid value for an argument): omerc: Invalid value for lat_01: |lat_0| should be < 90°
proj_create: Error 1027 (Invalid value for an argument): omerc: Invalid value for eccentricity
proj_create: Error 1027 (Invalid value for an argument): omerc: Invalid value for lat_1: |lat_1| should be < 90°
proj_create: Error 1027 (Invalid value for an argument): omerc: Invalid value for lat_2: |lat_2| should be < 90°
proj_create: Error 1027 (Invalid value for an argument): s2: Invalid value for s2 parameter: should be linear, quadratic, tangent, or none.
proj_create: Error 4096 (Unknown error (code 4096))
proj_create: Error 1027 (Invalid value for an argument): tpeqd: Invalid value for lat_1 and lat_2: their absolute value should be < 90°.
proj_create: Error 1027 (Invalid value for an argument): ups: Invalid value for es: only ellipsoidal formulation supported
proj_create: Error 1027 (Invalid value for an argument): urm5: Invalid value for n / alpha: n * sin(|alpha|) should be < 1.
proj_create: Error 1027 (Invalid value for an argument): utm: Invalid value for eccentricity: it should not be zero
proj_create: Error 1027 (Invalid value for an argument): utm: Invalid value for zone
proj_create: Error 1026 (Missing argument): topocentric: missing X_0 or lon_0
proj_create: Error 1026 (Missing argument): topocentric: missing Y_0 and/or Z_0
proj_create: Error 1026 (Missing argument): topocentric: missing lat_0
proj_create: Error 1028 (Mutually exclusive arguments): topocentric: (X_0,Y_0,Z_0) and (lon_0,lat_0,h_0) are mutually exclusive
-------------------------------------------------------------------------------
Reading file '/io/test/gie/builtins.gie'
     -----
     FAILURE in builtins.gie(3606):
     expected: 0.00 9050966.799011318013
     got:      0.000000000000   9050966.799187807366
     deviation:  0.176489 mm,  expected:  0.100000 mm
     -----
     FAILURE in builtins.gie(3612):
     expected: 0.00 -9050966.799011318013
     got:      0.000000000000   -9050966.799187807366
     deviation:  0.176489 mm,  expected:  0.100000 mm
-------------------------------------------------------------------------------
total: 2140 tests succeeded,  0 tests skipped,  2 tests FAILED!
-------------------------------------------------------------------------------

[kbevers]

I'll make a quick PR to use ctest --output-on-failure which will show these errors

Awesome. Thanks!

I guess the output you get here also shows all the noise made from tests that expect an error code. I recall having been annoyed by that previously without finding a good solution.

[erykoff]

Oh wonderful, it looks like the test is machine dependent at the pole (and I don't believe that part of the code was touched by my PR). I've changed the test to look very near the pole which is the primary problem I'm trying to solve here.