GDAL NetCDF corruption/compatability

HDF5 Library h5py
1

Hi there, i’m attempting to create a netcdf for a specifc data type (a digital elevation model) using a converted GeoTIFF. I’m converting this geotiff using GDAL. When I try to open the netcdf using h5py (which is the library i’m reliant on for what I need), it’s saying it can’t read it. Yet I can open the file absolutely fine using ncdump or xarray. Does anyone have any ideas as to how I can resolve this to get the file open in h5py? I can’t attach the file itself but here’s a copy of the ncdump -h output:

*ncdump -h Copernicus_DEM_test2.nc *
netcdf Copernicus_DEM_test2 {
dimensions:

  • lon = 15900 ;*
  • lat = 12609 ;*
    variables:
  • char crs ;*
  • crs:grid_mapping_name = "latitude_longitude" ;*

  • crs:long_name = "CRS definition" ;*

  • crs:longitude_of_prime_meridian = 0. ;*

  • crs:semi_major_axis = 6378137. ;*

  • crs:inverse_flattening = 298.257223563 ;*

  • crs:spatial_ref = "GEOGCS[\"WGS 84\",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[\"Latitude\",NORTH],AXIS[\"Longitude\",EAST],AUTHORITY[\"EPSG\",\"4326\"]]" ;*

  • crs:crs_wkt = "GEOGCS[\"WGS 84\",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[\"Latitude\",NORTH],AXIS[\"Longitude\",EAST],AUTHORITY[\"EPSG\",\"4326\"]]" ;*

  • crs:GeoTransform = "24.17635924150868 0.0002694945852358565 0 -28.51387459087979 0 -0.0002694945852358565 " ;*
  • double lat(lat) ;*
  • lat:standard_name = "latitude" ;*

  • lat:long_name = "latitude" ;*

  • lat:units = "degrees_north" ;*
  • double lon(lon) ;*
  • lon:standard_name = "longitude" ;*

  • lon:long_name = "longitude" ;*

  • lon:units = "degrees_east" ;*
  • float Band1(lat, lon) ;*
  • Band1:long_name = "GDAL Band Number 1" ;*

  • Band1:_FillValue = 9.96921e+36f ;*

  • Band1:grid_mapping = "crs" ;*

// global attributes:

  • :GDAL_AREA_OR_POINT = "Area" ;*

  • :Conventions = "CF-1.5" ;*

  • :GDAL = "GDAL 3.6.3, released 2023/03/07" ;*

  • :history = "Wed Nov 15 13:48:31 2023: GDAL CreateCopy( /home/bwright/Documents/active/insar/data/jagers_CGG_test/Copernicus_DEM_test2.nc, ... )" ;*

}

Thanks in advance!

2

Is there a chance that the converted file is netCDF-3, not netCDF-4?
Only netCDF-4 is based on HDF5.
h5py can open only netCDF-4, not netCDF-3.
You can use h5dump on your .nc file to make sure that the converted file is netCDF-4.

3

I don’t believe so. I seemed to resolve my issue by extracting ‘Band1’ from the gdal output using xarray using xds[‘Band1’].to_netcdf(‘path/to/output’). The other band appears to be a data array labeled ‘crs’ which contains some of the reference system information.

4

Diagnostics please. Show the output from these two shell commands.

ncdump -k Copernicus_DEM_test2.nc
od -N16 -c Copernicus_DEM_test2.nc

1 Like
5

Both tools report the file as “classic”. Apologies, I think you were right the first time around! I’ve done some further reading and it appears that GDAL does not default to NetCDF4 and the additional command of “-co FORMAT=NC4” is required when using gdal_translate.

2 Likes