[PROJ] Proj transformation for arbitrary (quaternion/Euler) spatial rotation

[Even Rouault]

Daven,

finding the right parameters is always a bit tricky. For example, for 
meteorological data, there are at least 2 different conventions to 
express angles.

See PROJ_WKT2_NAME_METHOD_POLE_ROTATION_GRIB_CONVENTION and 
PROJ_WKT2_NAME_METHOD_POLE_ROTATION_NETCDF_CF_CONVENTION at 
https://github.com/OSGeo/PROJ/blob/ae3180303d974653b126bc65dcfdad3cecb5a313/src/iso19111/operation/conversion.cpp#L3783 
where there's a mapping from south pole / north pole conventions to 
+proj=ob_tran. Setting lon_0 is often necessary.

Even

Le 09/11/2022 à 13:00, DAVEN P QUINN via PROJ a écrit :
> Hello all,
>
> I’m a geologist who is working with paleogeographic reconstructions 
> (see example here 
> <https://davenquinn.com/viz/corelle-demo-pbdb/?time=295>). These 
> require composite reproductions where different parts of a feature 
> dataset (continents) are rotated with different axis-angle 
> transformations. I’ve had good luck doing these rotations with 
> quaternion math in the browser/Python environments, but I am now 
> trying to use Proj transformations in order to apply rotations 
> directly within PostGIS queries.
>
> This linked image 
> <https://pbs.twimg.com/media/Fgpx4UGXEAQmPOt?format=jpg&name=4096x4096> 
> shows the desired result, a plate reconstruction to 250 Ma with a 
> different rotation applied to each plate. This was produced by 
> applying the desired quaternions through pl/pgsql math. Moving this 
> math to Proj internals would result in a ~50-100x speedup.
>
> The tool that seems most fitting is the `ob_trans` family of 
> projections, but I have been unable to define a rotation that can 
> handle my preferred representation (a pole defined in Lon-lat 
> coordinates and an associated angle of rotation around it). From my 
> reading of the docs, I believe the `o_lon_c`, `o_lat_c` and `o_alpha` 
> parameters should do this, but I cannot get them to work reliably. In 
> fact, I can’t even reliably define a ’no-op’ transformation that 
> leaves coordinates unchanged. Perhaps I have the angular coordinate 
> system wrong, or the rotation is being done in Cartesian space even 
> though `proj=lonlat` is used.
>
> Is it possible to define an arbitrary spherical spatial rotation in 
> Proj transformations? Maybe I need to use pipelines instead? I’d 
> appreciate any guidance! More details below the fold...
>
> Regards,
>
> *Daven P. Quinn*
> Research scientist II · /U of Wisconsin Madison/
> PhD · structural geology · /Caltech/ ‘18
> https://davenquinn.com
> +1 704 920 8487
>
> -------------
>
> Here is an example of the pl/pgsql I am currently using to assemble a 
> projection (I have gone through many iterations testing different 
> offsets and the different ways to specify the transformation):
> ```
> RETURN '+proj=ob_tran +o_proj=longlat +o_alpha=' || pi()/2+angle || 'r 
> +o_lon_c=' || pi()/2+lon || 'r +o_lat_c=' || lat || 'r' proj
> ```
> where (lat, lon, angle) defines a rotation pole.
>
> This results in some geometries (with a fortuitous set of poles, I 
> guess) attaining ballpark-correct transformations while other features 
> are shifted far outside of their origin tiles.
>
> _______________________________________________
> PROJ mailing list
> PROJ at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/proj

-- 
http://www.spatialys.com
My software is free, but my time generally not.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/proj/attachments/20221109/e029b128/attachment.htm>