Rework DEM.vcrs to support any 3D transformation (and fix with PROJ update)

[rhugonnet]

This PR reworks the vertical CRS attribute of DEM from .vref into .vcrs and uses more rigorously pyproj to define .vcrs as a vertical CRS (1D vertical axis) and ccrs as a 3D CRS (compound 2D + 1D, or promoted 2D to 3D CRS) required for vertical transformation.

The previous implementation stopped working because, since PROJ6, having 3D CRS as both input and output of a transformation is necessary to transform z (otherwise runs but without transforming). See: https://pyproj4.github.io/pyproj/stable/advanced_examples.html#promote-crs-to-3d.

Now, two methods are implemented in DEM: to_vcrs() and set_vcrs (very much like to_crs() and set_crs() in GeoPandas). The functionalities resides in a separate module vcrs.py.
I ran into some difficulties implementing the functionality because the function pyproj.CRS.to_2d() did not exist while pyproj.CRS.to_3d() does. This is now solved as it existed in PROJ and just needed to be wrapped, see related discussion: pyproj4/pyproj#1265 and related PR: pyproj4/pyproj#1267.

Also removes one time test that still failed randomly (forgotten in #355).

Details

This PR:

• Moves all vertical CRS operations in a vcrs.py module with his own tests test_vcrs.py,
• Adds a _vcrs_from_user_input() function to parse any vertical CRS from either common names ("EGM96", "EGM08", "Ellipsoid"), or a PROJ grid path, or an EPSG code or a CRS object,
• Adds a _build_vcrs_from_grid() function to consistently build a vertical CRS from a grid. If the grid does not exist locally, it attempts to download it from the PROJ fixed repo: https://cdn.proj.org/.
• Adds a _build_ccrs_from_crs_and_vcrs() function to consistently build a 3D CRS from a 2D/3D CRS stored in .crs by Raster (dropping the 3D if needed), and a vertical CRS in .vcrs. Additionally, if the target is the ellipsoid, the function casts to 3D.
• Adds a _transform_zz function that creates a TransformerGroup to list all possible transformers and, if the best_available is not there, downloads the grids automatically,
• Re-works DEM.set_vcrs and DEM.to_vcrs using the previous functions of vcrs.py,
• Adds an overridden DEM.from_array() to set vcrs during DEM creation from array,
• Adds a detection of 3D CRS in DEM.__init__ to set a vcrs if the .crs is 3D,
• Removes proj-data as a dependency,
• Forces pyproj>=3.4 to ensure that Transformer.transform() returns a np.ma.MaskedArray,
• Writes the documentation page on "Vertical referencing".

Issues

Resolves #344
Resolves #345
Resolves #331
Resolves #262
Resolves #311
Resolves #223

[erikmannerfelt]

Wow, what a big change with so few things for me to say! This looks like a great modification of the current behaviour and I only have nitpicks.

[rhugonnet]

Thanks a lot for the review @erikmannerfelt! 😄

I realized there was still a bit more to do to have everything practical and consistent!
What's new:

• Added the automatic downloading of PROJ grids and removed the proj-data dependency (fixes PROJ grids not found on Windows #311 and Remove proj-data dependency #345 in one go),
• Now, before transformation, pyproj.transformer.TransformerGroup is called to get all existing transformations for the area-of-interest sorted by best accuracy. If the best one is not available, it is downloaded automatically. Then it is used,
• I added an overridden from_array that accepts vcrs,
• I added a loop during DEM.__init__ that tries to estimate the vertical CRS from DEM.crs (in case it is 3D, which won't be the case for old archives, but is now getting more and more common! 😄).

And, most importantly, wrote the documentation page! Here is it: https://xdem-rhugonnet.readthedocs.io/en/fix_vref/vertical_ref.html

[rhugonnet]

Finally fixed the last little bit! All yours to review @erikmannerfelt!

One big question remains for me: Should we integrate the vcrs into the crs (overridding Raster.crs in DEM.crs)? This way the DEM would be written to file/read with the correct 3D CRS! 😄
(only small downside I see, is that people unaware of 3D CRS might not understand why their DEM CRS has changed, but that should become obvious from the doc)

[erikmannerfelt]

I only have documentation comments now!

Regarding your suggestion to override crs with vcrs, do you mean to override crs with ccrs, or am I missing something? I'm still a bit new to vcrs so it might be a dumb question!

[erikmannerfelt]

"[...] in most installations of proj [...]". On some installations like in nix, the proj folders are readonly.

[rhugonnet]

I'm not sure where you would place the addition?

[rhugonnet]

Thanks for the comments!

Regarding your suggestion to override crs with vcrs, do you mean to override crs with ccrs, or am I missing something? I'm still a bit new to vcrs so it might be a dumb question!

Yes, exactly, overriding it with ccrs so that it's always 3D! (then I guess we could remove ccrs and have it just be crs).

@adehecq: Curious about your feedback on this PR as well!
FYI: It's a problem also of interest to the ICESat-2 SlideRule team that wants to integrate vertical transforms in their cloud processing, so we're discussing it a bit there too in parallel.

[rhugonnet]

Merging as @erikmannerfelt is on fieldwork so I assume the review is finished! 😄
Opening an issue for overridding crs.