Establecer la proyección de NetCDF4 en Python

Question

Estoy utilizando la librería NetCDF4 para Python, pero no encuentro cómo establecer la proyección a WGS84.

import netCDF4 as nc

dsout = nc.Dataset(outfile, 'r')

time = dsout.createDimension('time', 0)
lat = dsout.createDimension('lat', rows)
lat = dsout.createVariable('lat', 'f4', ('lat',))
lat.standard_name = 'latitude'
lat.units = 'degrees_north'
lat.axis = "Y"
lat[:] = lats
lon = dsout.createDimension('lon', cols)
lon = dsout.createVariable('lon', 'f4', ('lon',))
lon.standard_name = 'longitude'
lon.units = 'degrees_east'
lon.axis = "X"
lon[:] = lons
times = dsout.createVariable('time', 'f4', ('time',))
times.standard_name = 'time'
times.long_name = 'time'
times.units = 'hours since 1970-01-01 00:00:00'
times.calendar = 'gregorian'
acc_precip = dsout.createVariable(
    'acc_precipitation_amount',
    'f4',
    ('time', 'lat', 'lon'),
    zlib=True,
    complevel=9,
    least_significant_digit=1,
    fill_value=-9999
)
acc_precip.standard_name = 'acc_precipitation_amount'
acc_precip.units = 'mm'

Answer 1

No sé mucho sobre la escritura de archivos NetCDF, pero lo siguiente me permitió escribir un archivo NetCDF que GDAL reconoció el CRS.

import numpy as np
import netCDF4 as nc

outfile = r'test.nc'
lats = tuple(range(-90, 90))
lons = tuple(range(0, 180))
cols = len(lons)
rows = len(lats)

dsout = nc.Dataset(outfile, 'w', clobber=True)

time = dsout.createDimension('time', 0)

lat = dsout.createDimension('lat', rows)
lat = dsout.createVariable('lat', 'f4', ('lat',))
lat.standard_name = 'latitude'
lat.units = 'degrees_north'
lat.axis = "Y"
lat[:] = lats

lon = dsout.createDimension('lon', cols)
lon = dsout.createVariable('lon', 'f4', ('lon',))
lon.standard_name = 'longitude'
lon.units = 'degrees_east'
lon.axis = "X"
lon[:] = lons

times = dsout.createVariable('time', 'f4', ('time',))
times.standard_name = 'time'
times.long_name = 'time'
times.units = 'hours since 1970-01-01 00:00:00'
times.calendar = 'gregorian'

acc_precip = dsout.createVariable(
    'acc_precipitation_amount',
    'f4',
    ('time', 'lat', 'lon'),
    zlib=True,
    complevel=9,
    least_significant_digit=1,
    fill_value=-9999
)
acc_precip[:] = np.ones((1,180,180))
acc_precip.standard_name = 'acc_precipitation_amount'
acc_precip.units = 'mm'
acc_precip.setncattr('grid_mapping', 'spatial_ref')

crs = dsout.createVariable('spatial_ref', 'i4')
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"]],AUTHORITY["EPSG","4326"]]'

gdalinfo de salida:

$ gdalinfo netcdf:test.nc
Warning 1: NetCDF driver detected file type=5, but libnetcdf detected type=3
Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute
Driver: netCDF/Network Common Data Format
Files: none associated
Size is 180, 180
Coordinate System is:
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"]],
    AUTHORITY["EPSG","4326"]]
Origin = (-0.500000000000000,89.500000000000000)
Pixel Size = (1.000000000000000,-1.000000000000000)
Metadata:
  acc_precipitation_amount#grid_mapping=spatial_ref
  acc_precipitation_amount#least_significant_digit=1
  acc_precipitation_amount#standard_name=acc_precipitation_amount
  acc_precipitation_amount#units=mm
  acc_precipitation_amount#_FillValue=-9999
  lat#axis=Y
  lat#standard_name=latitude
  lat#units=degrees_north
  lon#axis=X
  lon#standard_name=longitude
  lon#units=degrees_east
  NETCDF_DIM_EXTRA={time}
  NETCDF_DIM_time_DEF={1,5}
  NETCDF_DIM_time_VALUES=9.96921e+36
  spatial_ref#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"]],AUTHORITY["EPSG","4326"]]

  time#calendar=gregorian
  time#long_name=time
  time#standard_name=time
  time#units=hours since 1970-01-01 00:00:00
Corner Coordinates:
Upper Left  (  -0.5000000,  89.5000000) (  0d30' 0.00"W, 89d30' 0.00"N)
Lower Left  (  -0.5000000, -90.5000000) (  0d30' 0.00"W, 90d30' 0.00"S)
Upper Right ( 179.5000000,  89.5000000) (179d30' 0.00"E, 89d30' 0.00"N)
Lower Right ( 179.5000000, -90.5000000) (179d30' 0.00"E, 90d30' 0.00"S)
Center      (  89.5000000,  -0.5000000) ( 89d30' 0.00"E,  0d30' 0.00"S)
Band 1 Block=180x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
  Unit Type: mm
  Metadata:
    grid_mapping=spatial_ref
    least_significant_digit=1
    NETCDF_DIM_time=9.96921e+36
    NETCDF_VARNAME=acc_precipitation_amount
    standard_name=acc_precipitation_amount
    units=mm
    _FillValue=-9999