busstoptaktik / geodesy

Rust geodesy
Apache License 2.0
66 stars 6 forks source link

Support translation from PROJ strings to Geodesy format #56

Closed busstoptaktik closed 1 year ago

busstoptaktik commented 1 year ago

Over at geodesy-wasm, @Rennzie identifies the need for a

Proj string to geodesy_rs parser so that proj strings can be used with geodesy_wasm

Given that the Geodesy crate contains all necessary material for doing this, we really should support it directly in the library. Whether it should be directly supported as input format is less evident.

busstoptaktik commented 1 year ago

@Rennzie: I had forgotten how cursed the PROJ pipeline case really is (which is odd, since I introduced and implemented it in PROJ, some 5 years ago).

I think I have most cases correct now, but still missing the pathological case of inverted pipelines:

# Inversion of    a pipeline converting UTM zone 32 to UTM zone 33,
# hence obtaining a pipeline converting UTM zone 33 to UTM zone 32,
pipeline inv step proj=utm inv zone=32 step proj=utm zone=33

where

  1. The global inv should be inserted into the second step, but removed from the first, due to inv inv
  2. The steps should be collected in reverse order for output
Rennzie commented 1 year ago

@busstoptaktik thank you for implementing this.

It will be very helpful in geodesy-wasm. The pathological case you've now implemented is the one example I could think of after reading your e-mail (I'll be able to reply to the rest tomorrow).

I had a couple of other questions:

  1. Should the parser handle differences in operator keys? For example the gridshift operator has the key grid but PROJ has the keys hgridshift & vgridshift as documented in Rumination-002? A brief scan suggests omerc is the other operator that has differences.

  2. In the same vein, should the parser handle adapt. For example how would it handle a PROJ pipeline such as the below?

+proj=pipeline
  +step +proj=axisswap +order=2,1
  +step +proj=unitconvert +xy_in=deg +xy_out=rad
  +step +proj=webmerc +lat_0=0 +lon_0=0 +x_0=0 +y_0=0 +ellps=WGS84

I've not read the adapt code in depth but as I understand it so far it would handle both the unitconvert and axisswap steps of the above pipeline?

Should these be handled by the parser or should the callers be responsible for making these changes?

busstoptaktik commented 1 year ago

@Rennzie

  1. PROJ also has the gridshift [h+v] operator. In reverse, I intend to implement specializations to hgridshift/vgridshift as well. grid could change to grids once multigrid support is available.
    • omerc implements the options described in EPSG Geodetic guidance note 7-2, so it covers everything needed to handle transformations and systems represented in the EPSG registry. I would happily accept a pull request making it more PROJ compatible, but I do not intend to do any further work on omerc myself.
  2. adapt should be supplemented with the axisswap and unitconvert operators. adapt is an experiment showing how these operations can be described in declarative, rather than imperative terms. It was not intended as a replacement for the two PROJ operators, but as a supplement - I just happen to not have had any need for those two myself yet, so I have not gotten around to implement them.

So in brief, the pipeline structure is a means for stringing together existing operators, not a way of mediating between two different parametrisations of partially overlapping functionality (that would be a recipe for disaster 😊)

busstoptaktik commented 1 year ago

The pathological case you've now implemented is the one example I could think of after reading your e-mail

Don't worry - I have plenty of others (cf the test code) 😊