Yeah, I agree with your assessment. I hadn’t realized that you could set the 
option externally when I started the thread. I agree with you that the user 
likely set it themselves given what you have shown thus far, so I’ll sort it 
out with them. Apologies for taking some of your time.

Thanks again,
André

From: Even Rouault <even.roua...@spatialys.com>
Sent: April 22, 2022 13:47
To: Vautour, André (INT) <andre.vaut...@teledyne.com>; gdal-dev@lists.osgeo.org
Subject: Re: [gdal-dev] BAG CRS

External Email

André,

hum ok. it might perhaps be possible that WKT2 could come when extract the 
vertical part but I still can't reproduce that trying:

$ gdal_translate autotest/gcore/data/byte.tif byte.bag -a_srs 
'COMPOUNDCRS["foo",PROJCRS["NAD27 / UTM zone 
11N",BASEGEOGCRS["NAD27",DATUM["North American Datum 1927",ELLIPSOID["Clarke 
1866",6378206.4,294.978698213898,LENGTHUNIT["metre",1]]],PRIMEM["Greenwich",0,ANGLEUNIT["degree",0.0174532925199433]],ID["EPSG",4267]],CONVERSION["UTM
 zone 11N",METHOD["Transverse Mercator",ID["EPSG",9807]],PARAMETER["Latitude of 
natural 
origin",0,ANGLEUNIT["degree",0.0174532925199433],ID["EPSG",8801]],PARAMETER["Longitude
 of natural 
origin",-117,ANGLEUNIT["degree",0.0174532925199433],ID["EPSG",8802]],PARAMETER["Scale
 factor at natural 
origin",0.9996,SCALEUNIT["unity",1],ID["EPSG",8805]],PARAMETER["False 
easting",500000,LENGTHUNIT["metre",1],ID["EPSG",8806]],PARAMETER["False 
northing",0,LENGTHUNIT["metre",1],ID["EPSG",8807]]],CS[Cartesian,2],AXIS["(E)",east,ORDER[1],LENGTHUNIT["metre",1]],AXIS["(N)",north,ORDER[2],LENGTHUNIT["metre",1]],USAGE[SCOPE["Engineering
 survey, topographic mapping."],AREA["North America - between 120°W and 114°W - 
onshore. Canada - Alberta; British Columbia; Northwest Territories; Nunavut. 
Mexico. United States (USA) - California; Idaho; Nevada; Oregon; 
Washington."],BBOX[26.93,-120,78.13,-114]],ID["EPSG",26711]],VERTCRS["ellipse",VDATUM["NAD83(2011)
 / UTM zone 10N + ellipse"],CS[vertical,1],AXIS["ellipsoid height 
(h)",up,LENGTHUNIT["metre",1]]]]'

$ gdalinfo byte.bag -mdd xml:BAG | grep VERT_CS
            <gco:CharacterString>VERT_CS["ellipse",VERT_DATUM["NAD83(2011) / 
UTM zone 10N + ellipse",2005],UNIT["metre",1],AXIS["Ellipsoid 
height",UP]]</gco:CharacterString>

(the ,2005 here isn't really appropriate. it should rather be ,2002,  but such 
formulations of CompoundCRS with a VerticalCRS that expressed an ellipsoidal 
height are not ISO-19111 compliant)

Anyway if you can find a way to reproduce the issue with a CRS provided through 
the GetProjection()/SetProjection() API of GDAL, a fix would be welcome. If it 
is the user who provides WKT2 with a manual VAR_VERT_WKT option, then I'd say 
it is their responsibility to provide the right value.

Even
Le 22/04/2022 à 18:35, Vautour, André (INT) a écrit :
Sorry Even,

I should have done more research on my end up front before sending this. A 
client told us his WKT 2.0 BAG was created with GDAL. Looking at it some more, 
it was the vertical CRS that was WKT 2.0.

Their vertical CRS looked as follows:
<gco:CharacterString>VERTCRS["ellipse",VDATUM["NAD83(2011) / UTM zone 10N + 
ellipse"],CS[vertical,1],AXIS["ellipsoid height 
(h)",up,LENGTHUNIT["metre",1]]]</gco:CharacterString>

Cheers,
André

From: Even Rouault 
<even.roua...@spatialys.com><mailto:even.roua...@spatialys.com>
Sent: April 22, 2022 13:22
To: Vautour, André (INT) 
<andre.vaut...@teledyne.com><mailto:andre.vaut...@teledyne.com>; 
gdal-dev@lists.osgeo.org<mailto:gdal-dev@lists.osgeo.org>
Subject: Re: [gdal-dev] BAG CRS

External Email

André,

I don't confirm this trying:

$ gdal_translate autotest/gcore/data/byte.tif byte.bag

$ gdalinfo byte.bag -mdd xml:BAG | grep PROJCS
            <gco:CharacterString>PROJCS["NAD27 / UTM zone 
11N",GEOGCS["NAD27",DATUM["North_American_Datum_1927",SPHEROID["Clarke 
1866",6378206.4,294.978698213898,AUTHORITY["EPSG","7008"]],AUTHORITY["EPSG","6267"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4267"]],PROJECTION["Transverse_Mercator"],PARAMETER["latitude_of_origin",0],PARAMETER["central_meridian",-117],PARAMETER["scale_factor",0.9996],PARAMETER["false_easting",500000],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["Easting",EAST],AXIS["Northing",NORTH],AUTHORITY["EPSG","26711"]]</gco:CharacterString>

This is WKT1

I'm wondering if the confusion might come from looking at the CRS reported by 
gdalinfo, which will always be WKT2 in recent GDAL versions whatever the 
encoding in the source dataset itself ?

Or perhaps you try to write a Geographic 3D CRS ? in which case WKT2 will be 
used since there's no WKT1 representation for this

$ gdalwarp -overwrite autotest/gcore/data/byte.tif byte.bag -t_srs EPSG:4979 
-ot float32

$ gdalinfo byte.bag -mdd xml:BAG | grep GEOGCRS
            <gco:CharacterString>GEOGCRS["WGS 84",ENSEMBLE["World Geodetic 
System 1984 ensemble",MEMBER["World Geodetic System 1984 
(Transit)"],MEMBER["World Geodetic System 1984 (G730)"],MEMBER["World Geodetic 
System 1984 (G873)"],MEMBER["World Geodetic System 1984 (G1150)"],MEMBER["World 
Geodetic System 1984 (G1674)"],MEMBER["World Geodetic System 1984 
(G1762)"],MEMBER["World Geodetic System 1984 (G2139)"],ELLIPSOID["WGS 
84",6378137,298.257223563,LENGTHUNIT["metre",1]],ENSEMBLEACCURACY[2.0]],PRIMEM["Greenwich",0,ANGLEUNIT["degree",0.0174532925199433]],CS[ellipsoidal,3],AXIS["geodetic
 latitude 
(Lat)",north,ORDER[1],ANGLEUNIT["degree",0.0174532925199433]],AXIS["geodetic 
longitude 
(Lon)",east,ORDER[2],ANGLEUNIT["degree",0.0174532925199433]],AXIS["ellipsoidal 
height (h)",up,ORDER[3],LENGTHUNIT["metre",1]], [... snip ....] 
,ID["EPSG",4979]]</gco:CharacterString>

(I see that gdalinfo also issues a "ERROR 1: PROJ: proj_create_compound_crs: 
components of the compound CRS do not belong to one of the allowed combinations 
of http://docs.opengeospatial.org/as/18-005r4/18-005r4.html#34"; when reading 
this, which I suspect comes from the fact that it tries to create a CompoundCRS 
from this 3D GeogCRS and the "VERT_CS["unknown", VERT_DATUM["unknown", 2000]]", 
so writing 3D GeogCRS isn't really appropriate for BAG)

Even
Le 22/04/2022 à 18:01, Vautour, André (INT) a écrit :
Hi all,

It has just come to my attention that the BAG driver in GDAL is writing the CRS 
with a WKT codeSpace and a WKT 2.0 string as the code. While I am fairly sure 
we started that incorrect practice here at CARIS with a WKT 1.0 string, this is 
the first time I’ve seen a WKT 2.0 string written in there, and I am therefore 
concerned about interoperability.

I think we can all agree that what makes the most sense is to write an EPSG 
codebase with the EPSG code if we have it, and I think falling back to a WKT 
1.0 string if it is not EPSG is a reasonable fallback if it is not EPSG, given 
how it has been used in the past, but I am worried that older software won’t 
recognize the WKT 2.0 string.

Would you be willing to take a pull request that does that? Worst-case, would 
you be fine with changing that WKT 2.0 string to a WKT 1.0 string?

Regards,
André




_______________________________________________

gdal-dev mailing list

gdal-dev@lists.osgeo.org<mailto:gdal-dev@lists.osgeo.org>

https://lists.osgeo.org/mailman/listinfo/gdal-dev

--

http://www.spatialys.com

My 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
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev

Reply via email to