[Qgis-user] Does QGIS require units "m" for the projection coordinate variables in a netCDF/CF file?

[Thomas Lavergne]

Thank you Even and Richard,

Since Even opened a pull request on gdal's github (
https://github.com/OSGeo/gdal/pull/8407), I will continue the discussion
over there.

Thomas

lør. 16. sep. 2023 kl. 17:10 skrev Even Rouault <even.rouault at spatialys.com
>:

> Thomas,
> >
> > I am not from the GIS world and not even a QGIS user, so bear with me
> > if my answers do not use the expected vocabulary.
> It takes years to get fully familar with GIS oddities, and netCDF is an
> area where even a life long of experience will not be enough to fully
> overcome the "creativity" of the netCDF community...
> >
> > From my point of view, my netCDF files describe their geolocation in a
> > CF-compliant manner, with a grid_mapping CRS variable, and the x and y
> > projection coordinate variables. But QGIS does not recognize the CRS,
> > seemingly because I have my x and y as km (and I do specify
> > :units="km", I am not hiding this).
> >
> > Since I am also a programmer, I tried to find where in the QGIS code
> > the netCDF/CF geolocation is read and decoded, to see if there was a
> > strict test on :units="m". But I failed to find this in the code.
>
> This is actually not really done in QGIS, but in PROJ (more exactly
> QGIS's QgsProjUtils::identifyCrs() delegates to PROJ's proj_identify()
> in
>
> https://github.com/qgis/QGIS/blob/8723b82ec19a3bec7f39b46128c798c7a2ee4230/src/core/proj/qgsprojutils.cpp#L271
> to try to correlate the CRS definition from the netCDF file to one in
> the database).
>
> As you likely open it as a raster (and not a mesh as Richard
> mentionned), the GDAL netCDF driver is also involved since it is it that
> reconstructs a CRS object from the attributes of the netCDF CF conventions.
>
>  From what I can tell, both PROJ and GDAL do the job "as expected". It
> is just that the way this CRS is encoded is not going to find any match
> with a known CRS.  To get what you want, the GDAL netCDF driver should
> both modify the CRS definition to change the km unit to m, and alter the
> geotransform matrix (which gives the coordinate of the top left corner
> and pixel size) to multiply their value by 1000. This could be done, but
> should it be done...? I'm not so sure. The issue is more than the
> "philosophy" of netCDF data producers is sometimes at odds with the
> usual CRS definitions.
>
> That said, for that very particular case, there's a proj4_string
> attribute in Lambert_Azimuthal_Equal_Area variable, with value
> "+proj=laea +lat_0=90 +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m
> +no_defs +type=crs", which indicates the intent to have a metre unit.
> There's also a global attribute geospatial_bounds_crs = "EPSG:6931"
> which further correlate this. Hence this enhancement to the GDAL netCDF
> driver to renormalize to metre: https://github.com/OSGeo/gdal/pull/8407
>
> Even
>
> --
> 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/qgis-user/attachments/20230918/945f1ddc/attachment.htm>