Improve support for dynamic crs -> dynamic CRS transform at different epochs (since PROJ 9.4.0)

[agiudiceandrea]

Feature description

Since PROJ 9.4.0, the PROJ library can perform transformations involving coordinate epoch changes https://github.com/OSGeo/PROJ/pull/3884, https://github.com/OSGeo/gdal/pull/8340.

Refs.: https://github.com/qgis/QGIS/pull/43143, https://github.com/qgis/QGIS/pull/43188.

Additional context

No response

[nyalldawson]

There's no reason this shouldn't already be working in QGIS now -- the limitation was entirely on the PROJ side

[agiudiceandrea]

Even using QGIS 3.36.3 with PRJ 9.4.0 (OSGeo4W), the following error message is sometimes displayed, e.g.:

CRITICAL    Cannot transform between dynamic CRS at difference coordinate epochs

Unsupported Transformation
Transformation between EPSG:9989 - ITRF2020 @ 2024.25 and EPSG:4326 - WGS 84 is not currently supported.
The results will be unpredictable and should not be used for high accuracy work.

Using the GUI, it is not possible to set or change the Epoch for a layer or for the project.

It is not possible to set the Epoch in the "Reproject layer" or "Assign projection" processing algorithm's dialog windows or in the "Save Vector Layer as..." dialog window.

It is not possible to to set the Epoch in the translate function https://docs.qgis.org/testing/en/docs/user_manual/expressions/functions_list.html#transform.

Setting the Epoch for layer and project using PyQGIS, the coordinate transformation between different CRS @ Epoch doesn't give the same result as using ogr2ogr or cs2cs in some cases.

For example, Pointz_9989-2024_249.zip is a PointZ layer "EPSG:9989 - ITRF2020 @ 2024.25" (actually 2024.249: Apr 1st, 2024) with 1 PointZ having coordinates: -80.09282164175 43.918799362083334 387.746

Adding such layer to a new project in QGIS, the error message "CRITICAL Cannot transform between dynamic CRS at difference coordinate epochs" is displayed.

Reprojecting such layer using the "Convert format" processing algorithm with the additional options -s_srs EPSG:9989 -s_coord_epoch 2024.249 -t_srs EPSG:22717 -t_coord_epoch 2010 (or using the corresponding ogr2ogr command in the OSGeo4W Shell) produces a PointZ layer "EPSG:22717 - NAD83(CSRS)v7 / UTM zone 17N @ 2010" (Jan 1st, 2010) with 1 PointZ having coordinates: 572832.228 4863253.064 388.857

The same as using the cs2cs tool or the Natural Resources Canada TRX tool:

cs2cs -d 3 --s_epoch 2024.249 --t_epoch 2010 EPSG:9989 EPSG:22717
43.918799362083334 -80.09282164175 387.746
572832.228      4863253.064 388.857

Using QGIS, either with the on-the-fly reprojection or using the "Reproject layer" processing algorithm (whichever transformation is selected) or the "Save Vector Layer as..." dialog window (whichever transformation is set in the project between such CRS's), the reprojected coordinates are different: 572832.250 4863253.041 387.746

The same occurs using QgsCoordinateTransform:

sourceCrs = QgsCoordinateReferenceSystem('EPSG:9989')
sourceCrs.setCoordinateEpoch(2024.249)
destCrs = QgsCoordinateReferenceSystem('EPSG:22717')
destCrs.setCoordinateEpoch(2010)
tr = QgsCoordinateTransform(sourceCrs, destCrs, QgsProject.instance())
geom = QgsGeometry(QgsPoint(-80.09282164175, 43.918799362083334, 387.746))
geom.transform(tr)
print(geom)
<QgsGeometry: PointZ (572832.25049596535973251 4863253.04105716664344072 387.7459999999999809)>

In the QgsCoordinateTransoform API docs:

Since QGIS 3.20 The QgsCoordinateTransform class can perform time-dependent transformations between a static and dynamic CRS based on either the source OR destination CRS coordinate epoch, however dynamic CRS to dynamic CRS transformations are not currently supported.