[PROJ] ISEA projection: proposed parameters deprecation and input latitude interpretation

Tue Aug 27 02:47:49 PDT 2024

Hi Luís,

Thanks for the input.

> adding azi= instead.
To be clear, the +azi= parameter is already there for the forward 
projection, just not in the new inverse yet.

> I wonder whether the angular coordinates interpretation shouldn't also 
> be left to the user to decide. Eventually with the authalic 
> transformation as default behaviour.
Do you mean by introducing an additional parameter? I don't think this 
should be necessary, and I think it would go against how ellipsoids are 
generally handled for all projections in PROJ.

If using a perfectly spherical ellipsoid, all types of latitudes are the 
same.

The bold title I used for that section of my last e-mail was confusing, 
sorry.

What I meant is that if an oblate spheroid ellipsoid is used (which is 
the default if you don't specify any ellipsoid parameter e.g., simply 
saying +proj=isea), then the latitude should be interpreted as geodetic 
latitude, but it should be transformed to an authalic latitude prior to 
applying the ISEA projection which is only defined on a sphere, just 
like it is already being done for Lambert Azimuthal Equal-Area. This is 
the only correct thing to do I believe.

Some systems do, arguably incorrectly, take the geodetic latitude, and 
directly apply the ISEA spherical projection pretending that latitude is 
already the authatlic latitude on the sphere.
This could still be done with PROJ simply by specifying 
+R=6371007.18091875, and still giving the same geodetic latitude as input.

Other systems (still arguably incorrectly, like ours currently) convert 
the geodetic latitude to geocentric latitude before applying the ISEA 
projection.
This could still be achieved by specifying +proj=pipeline +step 
+proj=geoc +step +proj=isea +R=6371007.18091875

Therefore I think even these incorrect interpretations could still be 
explicitly specified by the user, using the non-projection-specific PROJ 
mechanisms, without the introduction of a projection-specific parameter.

But the default when using a non-spherical ellipsoid would be the 
correct way of first transforming geoetic latitude to authalic latitude.

I hope that makes sense.

Kind regards,

-Jerome

On 8/27/24 5:20 AM, Luí s Moreira de Sousa wrote:
> Hi Jérôme,
>
> this is all very logical, dropping those DGGS paramaters and adding 
> azi= instead.
>
> I wonder whether the angular coordinates interpretation shouldn't also 
> be left to the user to decide. Eventually with the authalic 
> transformation as default behaviour.
>
> Regards.
>
> -- 
> Luís Moreira de Sousa
> Mastodon: https://mastodon.social/@luis_de_sousa
> URL: https://ldesousa.codeberg.page
>
> Sent with Proton Mail <https://proton.me/> secure email.
>
> On Saturday, 24 August 2024 at 15:36, Jérôme St-Louis via PROJ 
> <proj at lists.osgeo.org> wrote:
>>
>> Dear PROJ community,
>>
>> I recently had the chance to integrate an implementation of the 
>> inverse transformation for the Icosahedral Snyder Equal-Area (ISEA) 
>> projection <https://proj.org/en/9.4/operations/projections/isea.html> 
>> (PR #4211 <https://github.com/OSGeo/PROJ/pull/4211> which was 
>> recently merged).
>>
>> Doing so resulted in some initial re-factoring to share code between 
>> the forward and inverse projection.
>>
>> *Deprecating grid parameters*
>>
>> During this process, I realized that several parameters for the 
>> projection (+proj=isea) really have nothing to do with the projection 
>> itself, but are parameters for constructing discrete global grids on 
>> top of an ISEA projection.
>>
>> In fact, the code applying these parameters 
>> <https://github.com/OSGeo/PROJ/blob/master/src/projections/isea.cpp#L931> 
>> happens as a later step after the default planar projection (with the 
>> exception of the trivial isea_tri_plane() 
>> <https://github.com/OSGeo/PROJ/blob/master/src/projections/isea.cpp#L656>).
>>
>> Therefore, I suggest this code really does not belong in the PROJ 
>> library, but outside, in the implementation of a discrete global grid 
>> system, such as in libraries like DGGRID 
>> <https://github.com/sahrk/DGGRID> and geogrid 
>> <https://github.com/mocnik-science/geogrid/>, and these parameters 
>> should be deprecated.
>>
>> In addition:
>>
>> - there is no documentation at all for the parameters and they are 
>> very poorly understood,
>> - there are no existing tests for these parameters,
>> - the new inverse transformation does not support these parameters.
>>
>> The three parameters proposed for deprecation are:
>>
>> *+mode* which can be set to /plane/ (the default), /*di*/, */dd/*, or 
>> /*hex*/
>> The default mode is /plane/ which would become the only thing supported.
>>
>> /hex/ referred to returning coordinates related to a hexagonal grid, 
>> but it is really not clear.
>> /di/ and /dd/ seemed to refer to a rhombus grid, in a way very 
>> similar to the approach we use for an ISEA9R DGGRS 
>> <https://docs.ogc.org/DRAFTS/21-038.html#isea9r-dggrs> (or ISEA4R for 
>> a quad-tree).
>> The ISEA planar projection can be very easily transformed to a space 
>> where icosahedron triangles merge into 10 axis-aligned rhombuses by 
>> applying a 60 degrees rotation and 30 degrees horizontal shearing, 
>> which can be done with +proj=affine.
>> This seems to be what isea_ptdd() 
>> <https://github.com/OSGeo/PROJ/blob/master/src/projections/isea.cpp#L669> 
>> is doing for /dd/.
>> /di/ seems to be referring to an extra step done in isea_ptdi() 
>> <https://github.com/OSGeo/PROJ/blob/master/src/projections/isea.cpp#L849>. 
>> I have no idea what this is doing.
>>
>> *+aperture* which is the aperture (refinement ratio) of the discrete 
>> global grid for those di, dd and hex modes
>>
>> *+resolution* which I am not sure how it affects the di, dd, and hex 
>> modes, but is only used with those modes proposed for deprecation
>>
>> Systems using these parameters could still extract the code for these 
>> extra steps in its current state 
>> <https://github.com/OSGeo/PROJ/blob/2ce6c7a9b1ca44b8a7cfb1acca7ff10c8c616774/src/projections/isea.cpp>, 
>> and invoke as an extra step after the planar projection (after 
>> undoing the trivial 180 deg rotation and translation that 
>> isea_tri_plane() 
>> <https://github.com/OSGeo/PROJ/blob/master/src/projections/isea.cpp#L656> 
>> is doing).
>>
>> If you make use of these parameters and object to their deprecation, 
>> please reply with feedback.
>> I plan to submit a pull request for a staged deprecation, which first 
>> mentions in the documentation that these parameters will be 
>> deprecated in a future version.
>>
>> On the flip side, I hope to eventually contribute support for the 
>> remaining *+azi*= parameter not yet supported by the inverse 
>> projection (which allows support for orientations such as the R. 
>> Buckminster Fuller's Dymaxion Orientation).
>> Possibly also eventually allowing to specify an arbitrary latitude 
>> and longitude for the first icosahedron vertex, as opposed to the two 
>> available options of +orient=pole and +orient=isea with fixed values.
>>
>> *
>> Proposing to interpret input latitude as authalic latitude*
>>
>> As discussed in the PR 
>> <https://github.com/OSGeo/PROJ/pull/4211#issuecomment-2275206640>, 
>> the ISEA projection is only defined on a sphere.
>>
>> Although I clarified this in the documentation 
>> <https://osgeo-proj--4211.org.readthedocs.build/en/4211/operations/projections/isea.html>, 
>> the most common use case is having geodetic WGS84 latitude and 
>> longitude coordinates which you want to project to the ISEA planar 
>> projection.
>>
>> This includes the use case of integrating data into a discrete global 
>> grid system using a DGGRS based on the ISEA projection, such as 
>> ISEA3H or ISEA9R.
>>
>> Before applying the spherical projection with a sphere whose surface 
>> area is the same as the GRS80 or WGS84 ellipsoid (authalic sphere), 
>> the geodetic latitude should be converted to an authalic latitude 
>> <https://en.wikipedia.org/wiki/Latitude#Authalic_latitude>.
>>
>> This is already done by PROJ in other projections, including the 
>> Lamber Azimuthal Equal-Area projection of which ISEA is modified 
>> form, if I understand correctly (see 
>> https://github.com/OSGeo/PROJ/blob/master/src/projections/laea.cpp#L38 
>> and https://github.com/OSGeo/PROJ/pull/4211#issuecomment-2284383509 ).
>>
>> I propose that the geodetic latitude be automatically converted to an 
>> authalic latitude when using the projection with a non-spherical 
>> ellipsoid (which is the default, and therefore currently quite prone 
>> to misuse).
>>
>> Please reply if you see any issue with this change, for which I also 
>> plan to submit a pull request.
>>
>> Thank you very much!
>>
>> Kind regards,
>>
>> -Jerome
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/proj/attachments/20240827/db9ba1a2/attachment.htm>