busstoptaktik / geodesy

Rust geodesy
Apache License 2.0
66 stars 6 forks source link

Support oblique mercator #41

Closed busstoptaktik closed 1 year ago

busstoptaktik commented 1 year ago

(As per request from @pka, on GeoRust Discord, for use in tile-grid):

@pka:

Non-mercator transformations are based on Proj. But would be great, if we could cover selected projections with Geodesy!

@busstoptaktik:

That would be awesome, although the selection of projections is very limited compared to PROJ's huge gamut

@pka:

I'm mainly interested in ETRS89 and https://spatialreference.org/ref/epsg/2056/ - would these be possible?

@busstoptaktik:

epsg:4258 to epsg:2056 is a two step process (cf. projinfo -s epsg:4258 -t epsg:2056), involving first a datum shift from ETRS89 to CH1903+, then a projection to "Hotine oblique mercator variant B".

The datum shift is registered in EPSG as a 3 parameter Helmert shift, which is supported by Geodesy, while the second step, the projection to Oblique Mercator is not (yet).

It is, however, documented in "EPSG Guidance Note 7 part 2", page 64ff, and should not be too hard to implement in Rust - I would happily accept a pull request 🙂

pka commented 1 year ago

To be more precise - for tile calculations using https://github.com/pka/tile-grid/tree/main/src/transform I only need WGS84 -> ETRS89 and WGS84 -> EPSG:2056. Transformation back to WGS84 would be nice to have.

busstoptaktik commented 1 year ago

@pka said:

I only need WGS84 -> ETRS89 and WGS84 -> EPSG:2056. Transformation back to WGS84 would be nice to have.

WGS84 (2D/3D: EPSG:4326/4327) and ETRS89 (2D/3D: EPSG:4258/4937) are both datum ensembles. WGS84 with an ensemble-internal consistency of approx 2 m, ETRS89 with an ensemble-internal consistency of approx. 0.1 m.

The WGS84 ensemble comprises all realizations from WGS84(Transit) to WGS84(G2139). The ETRS89 ensemble comprises all the national realizations of ETRS89, related to a munber of ETRFs.

At the 2 m accuracy assumed by using the WGS84 ensemble, WGS84 and ETRS89 are identical, as also indicated in the official remark for the 3 parameter Helmert transformation from ETRS89 to CH1903+:

This transformation is also given as CH1903+ to CHTRS95 (1) (code 1509). CHTRS95 is a national realization of ETRS89. May be taken as approximate transformation CH1903+ to WGS 84 - see code 1676.

So unless you know your coordinates are related to a specific WGS84 realization, you can safely assert WGS84==ETRS89, since your overall accuracy is not better than 2 m anyway. So WGS84->ETRS89->WGS84 is available out of the box by doing nothing.

The datum shift comes into play when sailing into EPSG:2056, which is a projected system on the CH1903+ datum. Since that datum shift is based on a 3 parameter Helmert opertion, it is already supported by Rust Geodesy.

The missing link is the implementation of the Hotine Oblique Mercator projection step. Which EDIT: ~I may take a look at soon - and would happily accept a pull request if anyone else should feel the urge to implement it, before I get around to doing so~ I'm currently working on, and have so far completed the forward case.

pka commented 1 year ago

Reading your explanation, I realized that I forgot the more important second part of the transformation - I meant ETRS89 / UTM32: https://spatialreference.org/ref/epsg/25832/

For tile calculations I have to translate from geographic coordinates (which I think are always WGS84 in general web mapping applications) to the coordinate reference system of the tiles, which is usually EPSG:3857. But I want also support other systems like EPSG:25832 and EPSG:2056.

busstoptaktik commented 1 year ago

EPSG:3857 "Web Mercator", is constructed by the geodetically dubious practice of converting ellipsoidal geographical coordinates (WGS84-ish) to Mercator-ish coordinates by applying a spherical implementation of the Mercator projection - hence throwing away the conformal properties of the projection.

At least, that's how it's implemented in PROJ, (here).

As Geodesy in general does not support spherical versions of projections, this will not work here. It is, however, trivial to implement, and I will give it a shot once the omerc implementation is merged.

busstoptaktik commented 1 year ago

@pka I have merged the support for Web Mercator and Oblique Mercator now. It would, however, be useful with some additional test points. Do you happen to have some samples of (say) homologous EPSG:25832 and EPSG:2056 coordinates?

pka commented 1 year ago

Thanks a lot! I will test this this week. I also looked for reference data, but didn't find anything with both geographic and projected coordinates.