写入netCDF文件时设置geotransform

如何解决写入netCDF文件时设置geotransform

我一直在使用netcdf4写入netCDF文件。文件写入本身可以正常工作，但是当我使用gdal加载文件并查看geotransform时，它应为（-180.0、0.5、0.0、90.0、0.0，-时，会产生（0.0、1.0、0.0、0.0、0.0、1.0）- 0.5）。写入文件时，我一定会丢失一些东西，因为当我从另一个机构加载.nc文件时，它包含转换。因此它必须作为变量存在于文件中的某个位置...

以下是创建.nc文件的一些代码（python 3.6.8）：

import numpy as np
import netCDF4 as nc

# make data
lon=np.arange(-180,181,5)
lat=np.arange(-90,91,5)
dummy_data=np.ones((len(lat),len(lon)))

# write netcdf file
grp = nc.Dataset("example_file.nc",'w')
grp.createDimension('lat',len(lat))

# set latitude values
latVar = grp.createVariable('lat','f8',('lat',),zlib = True,complevel =6,shuffle = True)
latVar[:] = lat
latVar.units="degrees_north"
latVar.long_name = 'latitude'
latVar.axis='Y'

# set longitude
grp.createDimension('lon',len(lon))
lonVar = grp.createVariable('lon',('lon',complevel = 6,shuffle = True)
lonVar[:] = lon
lonVar.units = 'degrees_east'
lonVar.long_name = 'longitude'
lonVar.axis='X'

crs_var = grp.createVariable('crs',np.int8,())
crs_var.standard_name = 'crs'
crs_var.grid_mapping_name = 'latitude_longitude'
crs_var.crs_wkt = 'GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],"6326"]],PRIMEM["Greenwich","8901"]],UNIT["degree",0.0174532925199433,"9122"]],"4326"]]'

var = grp.createVariable('dummy_data','f4','lon'),shuffle = True,fill_value=-9999)
var.missing_value=-9999.0 # set null value
var.unit='degrees'
var.grid_mapping = 'crs'
var[:] = dummy_data

grp.close()

使用gdal加载.nc文件时，我们可以看到geotransform不是我期望的...

from osgeo import gdal
datafile = gdal.Open("example_file.nc")
print( "Driver: ",datafile.GetDriver().ShortName,datafile.GetDriver().LongName)
print( "Size is ",datafile.RasterXSize,datafile.RasterYSize)
print( "Bands = ",datafile.RasterCount)
print( "Coordinate System is:",datafile.GetProjectionRef ())
print( "GetGeoTransform() = ",datafile.GetGeoTransform ())
print( "GetMetadata() = ",datafile.GetMetadata ())

这会产生...

Driver:  HDF5Image HDF5 Dataset
Size is  73 37
Bands =  1
Coordinate System is:
GetGeoTransform() =  (0.0,1.0,0.0,1.0)
GetMetadata() =  {'_NCProperties': 'version=2,netcdf=4.7.3,hdf5=1.10.6'}

当我加载其他机构生成的.nc文件时，将设置geotransform ...

Driver:  netCDF Network Common Data Format
Size is  720 360
Bands =  1
Coordinate System is:
GetGeoTransform() =  (-180.0,0.5,90.0,-0.5)

我已经检查了以下内容，但到目前为止还没有运气。

• https://www.esri.com/arcgis-blog/products/arcgis/data-management/creating-netcdf-files-for-analysis-and-visualization-in-arcgis/
• https://gis.stackexchange.com/questions/298211/python-netcdf4-writing-file-leaves-crs-as-undefined-and-extent-not-correct
• https://gis.stackexchange.com/questions/230093/set-projection-for-netcdf4-in-python

对于解决此问题的任何帮助将不胜感激。