[PROJ] Correcting metadata in GeoTIFF files

Tue Jan 25 06:58:05 PST 2022

Le 25/01/2022 à 15:43, Lesparre, Jochem a écrit :
>
> Thanks.
>
> I encountered a new problem:
>
> We use a horizontal grid as well as a vertical grid for the 
> transformation from ETRS89 (EPSG:4937) to our national CRS: RD 
> (EPSG:28992) with NAP height (EPSG:5709). With PROJ4-style strings 
> [1],an interpolation CRS (EPSG:4289) is needed for the Geodetic TIFF 
> grid of the geoid, next to the source (EPSG:4937) and target 
> (EPSG:5709) CRS. To specify this, I have to use --type 
> VERTICAL_TO_VERTICALin the Python script for the conversion from the 
> VDatum grid. This results in the description “vertical offset” instead 
> of “geoid undulation”. This is a bit strange, as the grid is in fact a 
> geoid. Problematic is that (according to the documentation [2]) the 
> source and target CRS both must be a vertical CRS, while one is our 
> national levelling height system NAP, but the other is ellipsoidal 
> ETRS89 height.
>
You should use --type GEOGRAPHIC_TO_VERTICAL if it is a geoid file (the 
name is perhaps a bit misleading as the values are to go from vertical 
CRS to geographic CRS but this is the EPSG naming convention for the 
transformation)
>
> The latter doesn’t have a vertical EPSG code. What is the WKT string 
> of ETRS89 ellipsoidal height?
>
You should use EPSG:4937 for ETRS89 geographic 3D
>
> Next to this, I have some new questions on the use of WKT strings 
> instead of EPSG codes for thespecification of a CRS in a Geodetic TIFF 
> grid file:
>
>   * Are the simpler EPSG codes preferred over WKT strings?
>
Yes
>
>   * Is it OK to use a WKT string for the source and a EPSG code for
>     the target,
>
Yes
>
>   * Can I use the WKT string given by EPSG.org for a CRS with an EPSG
>     code and edit the CRS name in this WKT string to use an alias
>     instead of the official EPSG name (see example [3])?
>
You likely need to remove the ID["EPSG",4289] node. And it is a bit 
weird to still use Amesfoort as the datum, since, as far as I 
understand, if 2 geographic 2D CRS that have same axis, units, datum and 
prime meridian, they should be the same... That said that will have 
little consequences since PROJ ignores that anyway.
>
>  *
>
> Regards, Jochem
>
> [1] PROJ4-style transformation from ETRS89 (EPSG:4937) to national 
> horizontal and vertical CRS of the Netherlands:
>
> cs2cs +init=epsg:4937 +to +proj=sterea +lat_0=52.156160556 
> +lon_0=5.387638889 +k=0.9999079 +x_0=155000 +y_0=463000 +ellps=bessel 
> +nadgrids=rdtrans2018.gsb +units=m +geoidgrids=naptrans2018.gtx 
> +vunits=m+no_defs
>
> [2] 
> https://github.com/OSGeo/PROJ-data/blob/master/grid_tools/README.md 
> <https://github.com/OSGeo/PROJ-data/blob/master/grid_tools/README.md>
>
> [3] WKT string of EPSG:4289 where I replaced the official EPSG name 
> “Amersfoort” to the official Dutch name “RD Bessel”:
>
> GEOGCRS["_RD Bessel_",
>
> DATUM["Amersfoort",
>
> ELLIPSOID["Bessel 1841",6377397.155,299.1528128,
>
> LENGTHUNIT["metre",1,
>
> ID["EPSG",9001]
>
> ],
>
> ID["EPSG",7004]
>
> ],
>
> ID["EPSG",6289]],
>
> CS[ellipsoidal,2,
>
> ID["EPSG",6422]
>
> ],
>
> AXIS["latitude (Lat)",north],
>
> AXIS["longitude (Lon)",east],
>
> ANGLEUNIT["degree",0.0174532925199433,
>
> ID["EPSG",9102]
>
> ],
>
> ID["EPSG",4289]
>
> ]
>
> *From:*Even Rouault <even.rouault at spatialys.com>
> *Sent:* maandag 24 januari 2022 21:17
> *To:* Lesparre, Jochem <Jochem.Lesparre at kadaster.nl>; PROJ at lists.osgeo.org
> *Subject:* Re: [PROJ] Correcting metadata in GeoTIFF files
>
> Jochem,
>
> Le 24/01/2022 à 20:09, Lesparre, Jochem via PROJ a écrit :
>
>     Dear list,
>
>     I'm from Netherlands Partnership Geodetic Infrastructure (NSGI).
>     We published NTv2 (.gsb) and VDatum (.gtx) grid files for the
>     national CRS of the Netherlands (RD coordinates with NAP height)
>     in 2019. PROJ converted these to GeoTIFF (.tif) in 2020 [1].
>     Unfortunately, some wrong information was stored in the metadata
>     of the GeoTIFF files in the conversion. We are now creating new
>     GeoTIFF files with corrected metadata, using the Python scrips
>     ntv2_to_gtiff.py [2] and vertoffset_grid_to_gtiff.py [3].
>
>     For most changes in the metadata this is straightforward, but
>     there is one more difficult issue:
>
>     We have two variants of the horizontal transformation:
>
>     1.Conventional 1-step transformation (variant 2): geographic
>     coordinates of national CRS ---[nl_nsgi_rdtrans2018.tif]---> ETRS89
>
>     2.Better 2-step transformation (variant 1): geographic coordinates
>     of national CRS ---[nl_nsgi_rdcorr2018.tif]---> corrected
>     geographic coordinates of national CRS
>     ---[7_parameter_transformation]---> ETRS89
>
>     Since there is no separate EPSG code for the corrected geographic
>     coordinates of the national CRS, we want to use the same EPSG code
>     for both corrected and uncorrected coordinates. Or will it give
>     problems when the source and target CRS of a GeoTIFF file are the
>     same?
>
> Instead of the target_crs_code metadata item, you could include a 
> target_crw_wkt with a WKT2 string. That can be done through the 
> ntv2_to_gtiff.py script as it can accept a WKT2 CRS string as the 
> value for --target-crs.
>
>     Next to this, I have some other questions:
>
>      1. Shouldn't a GeoTIFF grid file for a vertical transformation
>         have an accuracy band like a grid for a horizontal
>         transformation? Is it possible with the Python script to
>         create an accuracy band in the GeoTIFF from a VDatum grid file?
>
> You could possibly add with GDAL standard tools (gdal_translate, etc) 
> a band with the accuracy and a description of 
> "geoid_undulation_accuracy" (you may need to go through a VRT to 
> manually add the Units to the band). The 
> validate_vertical_offset_geographic_to_vertical() method of the 
> check_gtiff_grid.py script would likely have to be updated so that it 
> doesn't emit an information message about the new band not being 
> recognized. And the spec at 
> https://proj.org/specifications/geodetictiffgrids.html 
> <https://secure-web.cisco.com/1IPsQrcS7XiA25_SIbwFzb-QZAKsmLdiZCviwwLfntX_LP5RL4MhT86HqkIDLHoWnxjGPsOj9qatKSQRQ6s-3OgTu2NWmUxO08eXFQDCyBNoiM6PgcO9f94P0uPY5Tod4uKHZ2DYmeBgeH9t6R_uf9AxSLI6uZ7h4ReNQ8fupVoj7z9uP0to6Nlr1Kti7lJycLJvZgYL1hGJywSBwFWky-TKX6tQWP90b26Q0Iub853RFhi4Meej8qTYaHcMfPtczIIxD62IAkCtQyIv3XMdTneUYyFMRszve23cNIudpKQzeuZzt3yS_yl6bQIPnPMO6/https%3A%2F%2Fproj.org%2Fspecifications%2Fgeodetictiffgrids.html>as 
> well
>
>      2. Can I add a recommended_interpolation_method with the Python
>         scripts?
>
> Not currently.  Either enhance them or add the metadata item manually 
> through GDAL tools (gdal_translate -mo 
> recommended_interpolation_method=foo in.tif out.tif -co 
> COMPRESS=DEFLATE -co INTERLEAVE=BAND -co PREDICTOR=3). Note that PROJ 
> itself will ignore that item and always apply bilinear interpolation.
>
>      3. How should I supply the corrected GeoTIFF files to PROJ?
>
> Pull request against https://github.com/OSGeo/PROJ-data 
> <https://github.com/OSGeo/PROJ-data/tree/master/nl_nsgi>
>
>      4. Is it necessary for PROJ to use new files names to distinguish
>         them from the old version with incorrect metadata?
>
> No
>
>     5.
>
>     Regards, Jochem
>
>     [1] https://github.com/OSGeo/PROJ-data/tree/master/nl_nsgi
>     <https://github.com/OSGeo/PROJ-data/tree/master/nl_nsgi>
>
>     [2]
>     https://github.com/OSGeo/PROJ-data/blob/master/grid_tools/ntv2_to_gtiff.py
>     <https://github.com/OSGeo/PROJ-data/blob/master/grid_tools/ntv2_to_gtiff.py>
>
>     [3]
>     https://github.com/OSGeo/PROJ-data/blob/master/grid_tools/vertoffset_grid_to_gtiff.py
>     <https://github.com/OSGeo/PROJ-data/blob/master/grid_tools/vertoffset_grid_to_gtiff.py>
>
>     **
>
>     J. Lesparre
>
>     Netherlands Partnership Geodetic Infrastructure (NSGI.nl)
>
>
>
>     Disclaimer:
>     De inhoud van deze e-mail is vertrouwelijk en uitsluitend bestemd
>     voor de geadresseerde(n).
>     Gebruik, openbaarmaking, vermenigvuldiging, verspreiding en/of
>     verstrekking van deze informatie aan derden is niet toegestaan.
>     Op al onze producten en diensten zijn onze algemene
>     leveringsvoorwaarden van toepassing
>     [https://www.kadaster.nl/algemene-leveringsvoorwaarden
>     <https://www.kadaster.nl/algemene-leveringsvoorwaarden>].
>
>     Disclaimer:
>     This email and any files transmitted with it are confidential and
>     intended solely for the use of the individual or entity to whom
>     they are addressed.
>     If you are not the intended recipient, you are notified that
>     disclosing, copying, distributing or taking any action in reliance
>     on the contents of this information is strictly prohibited.
>     Our general terms and conditions of delivery apply to all our
>     products and services
>     [https://www.kadaster.com/general-terms-and-conditions
>     <https://secure-web.cisco.com/1fU8sg4FrMFO4TyIW9ZofFSwdqxB0Kag8iZvKOseM_86a_Jx5Z2YV1Hrtt_e8mtlHLMQ0H8zu1i3lanHmsFjggZd1QXIC7i_zohRumwvppYst-xwoiyIzVpRIaUe7MImO2MSCS5aoznYYHb6E-LdzL2tarUM3ByI1fBzPf6EX-yQBpfd_vLLLUaYH2he7-jy3HzpLC2iu2ee6zOAfHGtKhEMeppr62m_f186VKNfXWlLwAxhGXGppzPfI_g4RpQ9OMJiT8koLCKdzzpwTMdFYDY2DaE08KWCvO9WNaB0ciQrOfxgdxfjiQ802OXi4mhKe/https%3A%2F%2Fwww.kadaster.com%2Fgeneral-terms-and-conditions>].
>
>
>
>     _______________________________________________
>
>     PROJ mailing list
>
>     PROJ at lists.osgeo.org  <mailto:PROJ at lists.osgeo.org>
>
>     https://lists.osgeo.org/mailman/listinfo/proj  <https://lists.osgeo.org/mailman/listinfo/proj>
>
> -- 
> http://www.spatialys.com  <http://secure-web.cisco.com/1fXZS9zkl1X2JBd5DEK4YAq1gPWoDfSCpDSaFlOf5HpAUYrV7SIEw12ok2xJ1dLlvazjXyElbKEjEJLnTBo5NsMUQgMzUh1-_47mgVQpSA2Wih89heEKMsI5iXG5QbxXLAERdl1W-2up80jCvknOD8vY83_Wh-5g43WZIwWQNwaJ38eTscgBCOajh8Bgq8Xmw06STQekuuToXJWmYwQ6woQiNYlE3tHnjcvuhJOGF6aedm4xTRlC4hjYeQ9b8T3HPmwawXXFsXFnlPiuOSR06jIeZ0FCvBXtbexj0ZloQCmzM8W3sBGc11dNGXri3ia7w/http%3A%2F%2Fwww.spatialys.com>
> My software is free, but my time generally not.
>
>
> Disclaimer:
> De inhoud van deze e-mail is vertrouwelijk en uitsluitend bestemd voor 
> de geadresseerde(n).
> Gebruik, openbaarmaking, vermenigvuldiging, verspreiding en/of 
> verstrekking van deze informatie aan derden is niet toegestaan.
> Op al onze producten en diensten zijn onze algemene 
> leveringsvoorwaarden van toepassing
> [https://www.kadaster.nl/algemene-leveringsvoorwaarden].
>
> Disclaimer:
> This email and any files transmitted with it are confidential and 
> intended solely for the use of the individual or entity to whom they 
> are addressed.
> If you are not the intended recipient, you are notified that 
> disclosing, copying, distributing or taking any action in reliance on 
> the contents of this information is strictly prohibited.
> Our general terms and conditions of delivery apply to all our products 
> and services
> [https://www.kadaster.com/general-terms-and-conditions].

-- 
http://www.spatialys.com
My software is free, but my time generally not.

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/proj/attachments/20220125/610182ed/attachment-0001.html>