OSGeo / PROJ-data

Repository for proj datum grids (for use by PROJ 7 or later)
Other
72 stars 33 forks source link

Add Danish grids included in EPSG v. 10.088 #101

Closed kbevers closed 1 year ago

kbevers commented 1 year ago

Adds the grids

In response to https://github.com/OSGeo/PROJ/pull/3731

kbevers commented 1 year ago

Thanks, @rouault. The checks are failing because the grids reference EPSG codes that are still not available in a proper PROJ release. I assume we can just ignore that but it would be nice with a second pair of eyes on the grid metadata to confirm that I put in the right information.

rouault commented 1 year ago

The checks are failing because the grids reference EPSG codes that are still not available in a proper PROJ release

I've regenerated the Docker image against latest PROJ master, and restarted the checks. They still fail on "target_crs_epsg_code=1347 is not a valid EPSG code" which indeed appears not to be a valid CRS code. The Computed Min/Max=-8892.409,516.196 on the latitude (and longitude) as well of dk_sdfi_gs_2022.tif looks very suspicious. -8892 arcseconds is - 2.47 degrees. I assume this "GS intermediate datum" is not a real datum (even ancient geodesists couldn't be wrong by one degree!) but a "fictional" one ?

kbevers commented 1 year ago

These transformations and CRS's were a bit tricky to mangle into the EPSG database. They are all quite old CRS's and the original transformations were using the Horner polynomium operation. Those were an even bigger faf to get into the EPSG db so @busstoptaktik came up with a simpler solution. The idea is to approximate the transformation from ETRS89 to the old CRS's using a projection (e.g. tmerc) and then accounting for the residual errors using a horizontal grid adjustment. Since the EPSG doesn't do pipelines we had to introduce the intermediate datums to allow for a two-step transformation like this. The intermediate datums are then based on the adjustment grids.

EPSG:1347 is a Datum code which refer to the GS intermediate datum (GS-IRF). The grid is an transformation from ETRS89 to GS-IRF so I figured it was best to add the codes for those. I don't think that is wrong but it does collide with the check in check_gtiff_grid.py. Perhaps it's a bit too limited to only allow codes from the CRS category? Alternatively I can put in the code for the actual GS CRS (EPSG:10258) but that wouldn't be entirely correct.

The suspicious min/max values are only present at the limits of the grids and well outside the defined usage area of the transformations. They aren't causing any problems and are purely an artifact of the procedure used to generate the grid.

rouault commented 1 year ago

Alternatively I can put in the code for the actual GS CRS (EPSG:10258) but that wouldn't be entirely correct.

Not correct because of what ? Normally a geodetic CRS and its datum are assumed to be somewhat equivalent, once you've defined the prime meridian and coordinate system (at least in the ISO-19111 model)

Perhaps I'd get a better understanding if you could give the full pipeline with the starting and intermediate CRS/datum.

I think it is better to stick with CRS as the source and target of those grids, at least for consistency. If for some reason using EPSG:10258 wouldn't be correct in that situation, perhaps that could just be captured as a comment in TIFFTAG_IMAGEDESCRIPTION ?

kbevers commented 1 year ago

Not correct because of what ?

Because the grid doesn't provide a transformation from ETRS89 (EPSG:4258) to GS (EPSG:10258). The grid transforms from ETRS89 (EPSG:4258) to GS-IRF (EPSG:1347).

Perhaps I'd get a better understanding if you could give the full pipeline with the starting and intermediate CRS/datum.

Here you go:

$ ./bin/projinfo -s EPSG:4258 -t EPSG:10258 -o PROJ
Candidate operations found: 2
-------------------------------------
Operation No. 1:

unknown id, ETRS89 to GS-IRF (1) + Generalstabens System LCC, 0.5 m, Denmark - onshore Jutland, Funen, Zealand and Lolland.

PROJ string:
+proj=pipeline
  +step +proj=axisswap +order=2,1
  +step +proj=unitconvert +xy_in=deg +xy_out=rad
  +step +proj=hgridshift +grids=dk_sdfi_gs_2022.tif
  +step +proj=lcc +lat_1=56 +k_0=1 +lat_0=55 +lon_0=10.37775 +x_0=0 +y_0=0
        +a=6377019.27 +rf=300

-------------------------------------
Operation No. 2:

unknown id, Ballpark geographic offset from ETRS89 to GS-IRF + Generalstabens System LCC, unknown accuracy, World, has ballpark transformation

PROJ string:
+proj=pipeline
  +step +proj=axisswap +order=2,1
  +step +proj=unitconvert +xy_in=deg +xy_out=rad
  +step +proj=lcc +lat_1=56 +k_0=1 +lat_0=55 +lon_0=10.37775 +x_0=0 +y_0=0
        +a=6377019.27 +rf=300

For this we are only concerned with the first one. The gridshift step is from ETRS89 to GS-IRF and the LCC step is from GS-IRF (EPSG Datum) to GS (EPSG CRS based on GS-IRF).

I think it is better to stick with CRS as the source and target of those grids, at least for consistency. If for some reason using EPSG:10258 wouldn't be correct in that situation, perhaps that could just be captured as a comment in TIFFTAG_IMAGEDESCRIPTION ?

Sure. I don't really have an opinion and as far as I am aware PROJ isn't using the metadata for anything.

rouault commented 1 year ago

The gridshift step is from ETRS89 to GS-IRF and the LCC step is from GS-IRF (EPSG Datum) to GS (EPSG CRS based on GS-IRF).

hum, seeing that EPSG:10258 is a projected CRS, which uses EPSG:10256 "GS-IRF" as its base geodetic CRS. Shouldn't EPSG:10256 be used then as the target CRS ? From your explanation I kind of understand that not really, but I believe the grid_check tool will reject providing a projected CRS as the target. So perhaps using EPSG:10256 + a comment that an extra further LCC step is actually needed to really get into the datum ?

kbevers commented 1 year ago

hum, seeing that EPSG:10258 is a projected CRS, which uses EPSG:10256 "GS-IRF" as its base geodetic CRS. Shouldn't EPSG:10256 be used then as the target CRS ?

Yeah, that sounds correct. I must have got the concepts confused when looking trying to find the correct codes to use. I'll make the change and see how it goes. I believe a similar thing is needed for the rest of grids as well.

rouault commented 1 year ago

"TIFFTAG_IMAGEDESCRIPTION=Transformation grid from ETRS89 to the GS intermediate datum" for dk_sdfi_gsb_2022.tif lacks the B of GSB.

kbevers commented 1 year ago

"TIFFTAG_IMAGEDESCRIPTION=Transformation grid from ETRS89 to the GS intermediate datum" for dk_sdfi_gsb_2022.tif lacks the B of GSB.

Thanks. Should be correct now.