Le 25/01/2024 à 17:08, Kirk Waters - NOAA Federal a écrit :
Even,Thanks for the explanation and tips. I'm not clear on how setting to 26988 (Michigan North meters) would help.
Hopefully forcing 1.0 will fix the issue.
ah, I didn't realize that the "NAD83 / Michigan North" in your WKT was *not* official EPSG:26988. I'd strongly suggest you modify your WKT to use "NAD83 / Michigan North (ftuS)" instead to avoid any confusion!
Actually your WKT is a bit weird as it lacks the name of the CompoundCRS (quite surprise that PROJ manages to parse it) before the PROJCS[]. Try the following instead
COMPD_CS["NAD83 / Michigan North (ftUS) + NAVD88 height (ftUS)",PROJCS["NAD83 / Michigan North (ftUS)",GEOGCS["NAD83",DATUM["North_American_Datum_1983",SPHEROID["GRS 1980",6378137,298.257222101,AUTHORITY["EPSG","7019"]],AUTHORITY["EPSG","6269"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4269"]],PROJECTION["Lambert_Conformal_Conic_2SP"],PARAMETER["standard_parallel_1",47.0833333333333],PARAMETER["standard_parallel_2",45.4833333333333],PARAMETER["latitude_of_origin",44.7833333333333],PARAMETER["central_meridian",-87],PARAMETER["false_easting",26246666.6666667],PARAMETER["false_northing",0],UNIT["US survey foot",0.304800609601219],AXIS["X",EAST],AXIS["Y",NORTH]],VERT_CS["NAVD88 height (ftUS)",VERT_DATUM["North American Vertical Datum 1988",2005],UNIT["US survey foot",0.304800609601219],AXIS["Gravity-related height",UP],AUTHORITY["EPSG","6360"]]]
By itself, that will not improve the interoperability, but at least things will be better labelled, with a GTCitationGeoKey (Ascii,53): "NAD83 / Michigan North (ftUS) + NAVD88 height (ftUS)"
Even
If I do a gdal_translate with that, it shows as being in New Zealand, just as Arc had done. It would certainly make lots more sense to actually put the data in a standard projection that has an EPSG code (meters or international feet), but that gets back to the embarrassing tale. It may well be that I'll have to force GeoTiff 1.0. I think there are other software packages that haven't caught up. Even when they do catch up, I've come across people running 15 year old CAD software and complaining about the DXF version GDAL produces.Thanks to the whole team for the work you do to maintain GDAL. KirkOn Thu, Jan 25, 2024 at 10:47 AM Even Rouault <[email protected]> wrote:Kirk, perhaps you could try using -a_srs EPSG:26988+6360 instead of a full WKT. That way GDAL will *not* write the projected parameter definitions, but just reference the horizontal and vertical CRS codes: Geotiff_Information: Version: 1 Key_Revision: 1.1 Tagged_Information: ModelTiepointTag (2,3): 0 0 0 440720 3751320 0 ModelPixelScaleTag (1,3): 60 60 1 End_Of_Tags. Keyed_Information: GTModelTypeGeoKey (Short,1): ModelTypeProjected GTRasterTypeGeoKey (Short,1): RasterPixelIsArea GTCitationGeoKey (Ascii,46): "NAD83 / Michigan North + NAVD88 height (ftUS)" ProjectedCSTypeGeoKey (Short,1): PCS_NAD83_Michigan_North VerticalCSTypeGeoKey (Short,1): Code-6360 (NAVD88 height (ftUS)) End_Of_Keys. End_Of_Geotiff. Even Le 25/01/2024 à 16:17, Even Rouault via gdal-dev a écrit :Kirk, this is a complicated story indeed... What changed since https://trac.osgeo.org/gdal/changeset/31405 is that OGC GeoTIFF 1.1 was published, which mostly fixes issues with vertical CRS, which your CRS uses. So by default GDAL 3.1 or later, when seeing a CompoundCRS write it as GeoTIFF 1.1, and this new GeoTIFF version was the opportunity to remove writing the ProjLinearUnitsInterpCorrectGeoKey=3059 hack, since using GeoTIFF 1.1 is a way to indicate that you're implementing correctly the standard. If you force writing GeoTIFF 1.0 by adding -co GEOTIFF_VERSION=1.0, GDAL will write the ProjLinearUnitsInterpCorrectGeoKey when needed, and listgeo will show the "Unknown-3059 (Short,1): Unknown-1" geotiff key. To know if a Geotiff is 1.0 or 1.1 with listgeo , look at the Key_Revision in the 3rd line which will be 1.0 or 1.1 So I assume ArcPro doesn't specifically implement 1.1. Either they implement GeoTIFF themselves, or they use GDAL < 3.1 that hasn't GeoTIFF 1.1 explicit support Even Le 25/01/2024 à 15:24, Kirk Waters - NOAA Federal via gdal-dev a écrit :-- http://www.spatialys.comHi, I've run across an odd issue with GeoTIFFS using custom projections written by GDAL and their interpretation in ArcPro. The specific case is doing a projection for Michigan State Plane North using U.S. survey feet. [insert embarrassing tale of another US agency wanting to always use survey feet regardless of state legislation or that the federal government has been officially metric for a long time]. For a file that should land in the Michigan upper peninsula, ArcPro puts it in New Zealand. It appears the issue is that Arc is interpreting the false easting as meters instead of survey feet. If I do a define projection in ArcPro and make a custom projection, it lands in the right place. I used an old listgeo from the geotiff library to look at the tags both in the original file and a copy that had define projection applied to see what was different. They looked essentially the same, including having the same value for the false easting. However, the Arc version had an extra tag (3059) that wasn't in the GDAL produced one. An internet search turned up this 13-year old GDAL issue <https://trac.osgeo.org/gdal/ticket/3901>. If I understand it correctly, a bug in GDAL was found that always used meters for the false easting and northing. In addition to fixing that, a tag (3059) was added to indicate this was fixed so software could differentiate broken GDAL output from new output. It looks like Arc is looking for this tag and if it's not there, false easting is assumed to be meters. A few questions, assuming I interpreted that issue correctly: 1) Why is gdal not adding that 3059 tag? (GDAL command used is below) 2) If it's only supposed to be relevant to GDAL written files, why is Arc writing it? 3) How do I get GDAL to add the tag? Command used in file creation: gdal_translate -a_srs 'COMPD_CS[PROJCS["NAD83 / Michigan North",GEOGCS["NAD83",DATUM["North_American_Datum_1983",SPHEROID["GRS 1980",6378137,298.257222101,AUTHORITY["EPSG","7019"]],TOWGS84[0,0,0,0,0,0,0],AUTHORITY["EPSG","6269"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4269"]],PROJECTION["Lambert_Conformal_Conic_2SP"],PARAMETER["standard_parallel_1",47.08333333333334],PARAMETER["standard_parallel_2",45.48333333333333],PARAMETER["latitude_of_origin",44.78333333333333],PARAMETER["central_meridian",-87],PARAMETER["false_easting",26246666.6666667],PARAMETER["false_northing",0],UNIT["US survey foot",0.304800609601219],AXIS["X",EAST],AXIS["Y",NORTH]],VERTCRS["NAVD88 height (ftUS)", VDATUM["North American Vertical Datum 1988"], CS[vertical,1], AXIS["gravity-related height (H)",up], LENGTHUNIT["US survey foot",0.304800609601], ID["EPSG",6360]]]' -co TILED=YES -co COMPRESS=DEFLATE -co PREDICTOR=3 -co BIGTIFF=IF_SAFER input.tif output.tif Kirk Waters, PhD NOAA Office for Coastal Management Applied Sciences Program coast.noaa.gov/digitalcoast <http://coast.noaa.gov/digitalcoast> _______________________________________________ gdal-dev mailing list [email protected] https://lists.osgeo.org/mailman/listinfo/gdal-devMy software is free, but my time generally not. _______________________________________________ gdal-dev mailing list [email protected] https://lists.osgeo.org/mailman/listinfo/gdal-dev-- http://www.spatialys.comMy software is free, but my time generally not.
-- http://www.spatialys.com My software is free, but my time generally not.
_______________________________________________ gdal-dev mailing list [email protected] https://lists.osgeo.org/mailman/listinfo/gdal-dev
