[gdal-dev] Where is the projection info stored in a GDAL dataset?

[Joaquim Manuel Freire Luís]

Thanks Even,

That led me to OGR_L_GetSpatialRef()

julia> Gdal.OGR_L_GetSpatialRef(getlayer(ds2,0).ptr)
Ptr{Nothing} @0x0000000067ce4610

and

julia> toWKT(Gdal.SpatialRef(p))
"GEOGCS[\"unknown\",DATUM[\"WGS_1984\",SPHEROID[\"WGS 84\",6378137,298.257223563,AUTHORITY[\"EPSG\",\"7030\"]],AUTHORITY[\"EPSG\",\"6326\"]],PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]],UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AXIS[\"Longitude\",EAST],AXIS[\"Latitude\",NORTH]]"

Just have to make this a bit more smooth.


From: Even Rouault <even.rouault at spatialys.com>
Sent: Sunday, March 28, 2021 7:14 PM
To: Joaquim Manuel Freire Luís <jluis at ualg.pt>; gdal-dev at lists.osgeo.org
Subject: Re: [gdal-dev] Where is the projection info stored in a GDAL dataset?


Joaquim,

ogr2ogr -t_srs passes the CRS object as an argument of the GDALDataset::CreateLayer() API

GDALDataset::SetProjection() / GDALSetProjection() is (at least in all implementations of in-tree drivers) reserved for the raster API, and raster drivers that support the Create() API.

Similarly GDALDataset::GetSpatialRef() / GDALGetSpatialRef() is for the raster API

For the vector API, use OGRLayer::GetSpatialRef() (different layers in a same dataset could have different CRS)

Best regards,

Even
Le 28/03/2021 à 19:57, Joaquim Manuel Freire Luís a écrit :
Hi,

I am trying to use the GDAL shared lib directly from Julia and there are things that keep breaking my head. For example, the lines below use a mix of GMT and GDAL wrapper functions calls that WORK. First command creates a GMT dataset with two UTM points type with an assigned coordinate system. Next, in 2), I create a GDAL dataset of a LineString, convert it to geog in 3) and save it to file using the OGR_GMT driver in 4)

1) D2 = mat2ds([588977.324434907 4095339.69117669; 579551.921920776 4150721.71333039], proj="+proj=utm +zone=29");
2) ds = gmt2gd(D2);
3) ds2=ogr2ogr(ds, ["-t_srs", "+proj=longlat", "-overwrite"])
4) ogr2ogr(ds2, dest="lixo2.gmt")

Now the troubles. When I tried to use GDALSetProjection() I got an error like this

ERROR 6: /vsimem/##269: Dataset does not support the SetSpatialRef() method.

So I found that I could set the CRS with GDALDatasetCreateLayer() while creating a layer on the dataset (using the Memory driver)
(https://github.com/GenericMappingTools/GMT.jl/blob/master/src/gdal_utils.jl#L254)

But oddly trying to inquiry the *ds* with GDALGetSpatialRef() returns an empty pointer. I.e.:

julia> Gdal.GDALGetSpatialRef(ds.ptr)
Ptr{Nothing} @0x0000000000000000

However the CSR was well sored in some place because the ="lixo2.gmt" is correct (it has this)

# @VGMT1.0 @GLINESTRING
# @R-8.1/-8/37/37.5
# @Jp"+proj=longlat +datum=WGS84 +no_defs"
# @Jw"GEOGCS[\"unknown\",DATUM[\"WGS_1984\",SPHEROID[\"WGS 84\",6378137,298.257223563,AUTHORITY[\"EPSG\",\"7030\"]],AUTHORITY[\"EPSG\",\"6326\"]],PRIMEM[\"Greenwich\",0,AUTHORITY[\"EPSG\",\"8901\"]],UNIT[\"degree\",0.0174532925199433,AUTHORITY[\"EPSG\",\"9122\"]],AXIS[\"Longitude\",EAST],AXIS[\"Latitude\",NORTH]]"
# FEATURE_DATA
>
-8 37.0
-8.1 37.5

But my trouble is that I want to fish the CRS from the ogr2ogr dataset (the ds2 above) so that I can convert back from GDAL to GMT (which I already can) without having to save to a file in disk.

So basically my question is where is that CRS stored and how can I access it?

Thanks

Joaquim



_______________________________________________

gdal-dev mailing list

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

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

--

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/gdal-dev/attachments/20210328/a1ee5562/attachment-0001.html>