NOTICE: This version of the NSF Unidata web site (archive.unidata.ucar.edu) is no longer being updated.
Current content can be found at unidata.ucar.edu.
To learn about what's going on, see About the Archive Site.
NOTICE: This version of the NSF Unidata web site (archive.unidata.ucar.edu) is no longer being updated.
Current content can be found at unidata.ucar.edu.
To learn about what's going on, see About the Archive Site.
Since your question (and possible answers) involve code, I suggest
submitting this question to stackoverflow.com or gis.stackexchange.com
and tagging with "netcdf" and "python". It typically works a lot
better than e-mail for these type of questions.
On Thu, Feb 27, 2014 at 1:18 PM, Michael Nolde
<nolde@xxxxxxxxxxxxxxxxxxxxxx> wrote:
> Dear list members,
>
> I have been trying to write a netcdf with an irregular grid for weeks now,
> but failed. I used 'Python' with 'python-netcdf4' package, as well as 'R'
> with
> 'ncdf4'-package. I also installed 'cdo' with netcdf-support and tried to
> reproject
> the data in my output file to my input file (with correct irregular grid).
> My current approch is to create a regular grid within Python and assign the
> values from the irregular grid to it, which costs resolution and takes
> considerably
> more disk space (source code attached below). Of course, I searched the web
> up
> and down, but could not find anything related to that problem. Could anyone
> please help me?
> Any way to create a netcdf file with an irregular grid would be highly
> appreciated.
>
>
>
> Here's the details:
>
> I've got a netCDF file which I want to read, modify, and save, which works
> well with R ('ncdf4') and Python ('python-netcdf').
> The problem is, that the grid representing the data (temperature values) is
> irregular, and I am not able to output an irregular grid. In the original
> file I have two variables ('lat' and 'lon'), which are both two-dimensional
> (356 x 236), with each entry representing the latitude or longitude,
> respectively, of each cell in the original grid (356 x 236 x 365 [days]).
>
> In the input file, the data grid (TMAX_2M, please see blow) references the
> 'lat' and 'lon' variables by the argument 'coordinates = "lon lat"', as far
> as I understood, but I have been unable to produce the desired behaviour
> even if I hand the datagrid with the 'coordinates'-flag to the output file.
>
> I would love to know how to make the netCDF file aware that the
> latitude/longitude coordinates of the data grid cells are stored in the
> "external" 'lat' and 'lon' variables, and not in the grid object itself.
> Below I'm posting the structure of my input file, which I hope to reproduce,
> as well as my source code (in 'R') so far.
> Any help is highly appreciated, thank you very much in advance.
>
>
>
> Best regards, Michael
>
>
> this is the structure of the source netCDF file:
> -----------------------------------------------------
>
> dimensions:
> x = 356;
> y = 236;
> height_2m = 1;
> time = UNLIMITED; // (365 currently
> nb2 = 2;
>
> variables:
> float lon(y=236, x=356);
> :standard_name = "longitude";
> :long_name = "longitude";
> :units = "degrees_east";
> :_CoordinateAxisType = "Lon";
>
> float lat(y=236, x=356);
> :standard_name = "latitude";
> :long_name = "latitude";
> :units = "degrees_north";
> :_CoordinateAxisType = "Lat";
>
> float height_2m(height_2m=1);
> :standard_name = "height";
> :long_name = "height above the surface";
> :units = "m";
> :positive = "up";
> :axis = "Z";
>
> double time(time=365);
> :standard_name = "time";
> :long_name = "time";
> :bounds = "time_bnds";
> :units = "seconds since 1979-01-01 00:00:00";
> :calendar = "proleptic_gregorian";
>
> double time_bnds(time=365, nb2=2);
> :units = "seconds since 1979-01-01 00:00:00";
> :calendar = "proleptic_gregorian";
>
> float TMAX_2M(time=365, height_2m=1, y=236, x=356);
> :standard_name = "air_temperature";
> :long_name = "2m maximum temperature";
> :units = "K";
> :coordinates = "lon lat";
> :cell_methods = "time: maximum";
> :_FillValue = -9.0E33f; // float
>
>
>
> this is my 'Python' - code so far,
> based on an example by Jeff Whitaker (found on the Unidata-website)
> ----------------------------------------------------------------------------------------------
> # the Scientific Python netCDF 3 interface
> # http://dirac.cnrs-orleans.fr/ScientificPython/
> from Scientific.IO.NetCDF import NetCDFFile as Dataset
>
> # the 'classic' version of the netCDF4 python interface
> # http://code.google.com/p/netcdf4-python/
> #from netCDF4_classic import Dataset
> from numpy import arange, dtype # array module from http://numpy.scipy.org
>
>
> import netCDF4
> import numpy as np
> import sys
>
> f = netCDF4.Dataset('TMAX_2M_daily1981_box.nc', 'r')
> for v in f.variables:
> print(v),
>
> type(f.variables) # f.variables is a dictionary data structure
>
> for d in f.dimensions:
> print(d),
> print
>
>
> x = f.dimensions['x']
> y = f.dimensions['y']
> #time = f.dimensions['time']
>
> temp = f.variables['TMAX_2M']
> lon = f.variables['lon']
> lat = f.variables['lat']
> time = f.variables['time']
> height = f.variables['height_2m']
>
>
> npLon = np.array(lon)
> npLat = np.array(lat)
>
> npLonMin = npLon.min()
> npLonMax = npLon.max()
> npLatMin = npLat.min()
> npLatMax = npLat.max()
>
>
> lonOut = np.arange(-18.75, 59.00, 0.25)
> latOut = np.arange(24.75, 60.25, 0.25)
>
>
> tempOut = np.zeros(shape=(2,len(latOut),len(lonOut)),dtype='float32')
> tempOut.fill(np.nan)
>
>
>
> iTime = 0
> iHeight = 0
> tempOutOld = 0
>
> for iY in range(0,236):
> for iX in range(0,356):
> cellX = int(round((lon[iY,iX] - npLonMin) * 4,0))
> cellY = int(round((lat[iY,iX] - npLatMin) * 4,0))
> tempOut[iTime,cellY,cellX] = temp[iTime,iHeight,iY,iX]
>
>
>
> """
> This is an example program which writes some 4D pressure and
> temperatures.
> The companion program pres_temp_4D_rd.py shows how
> to read the netCDF data file created by this program.
>
> This example demonstrates the netCDF Python API.
> It will work either with the Scientific Python NetCDF version 3 interface
> (http://dirac.cnrs-orleans.fr/ScientificPython/)
> of the 'classic' version of the netCDF4 interface.
> (http://netcdf4-python.googlecode.com/svn/trunk/docs/netCDF4_classic-module.html)
> To switch from one to another, just comment/uncomment the appropriate
> import statements at the beginning of this file.
>
> Jeff Whitaker <jeffrey.s.whitaker@xxxxxxxx> 20070202
> """
>
> # the netCDF variable will be nrecs x nlevs x nlats x nlons.
> nrecs = 2; nlevs = 2; nlats = 142; nlons = 311
> #1*284*622
>
> #lonOut = np.arange(-18.75, 59.00, 0.125)
> #latOut = np.arange(24.75, 60.25, 0.125)
>
> # open a new netCDF file for writing.
> ncfile = Dataset('pres_temp_4D.nc','w')
>
> # latitudes and longitudes of grid
> #lats_out = 0.0 + 5.0*arange(nlats,dtype='float32')
> #lons_out = 0.0 + 5.0*arange(nlons,dtype='float32')
> lats_out = arange(start=24.75, stop=60.25, step=0.25,dtype='float32')
> lons_out = arange(start=-18.75, stop=59.0, step=0.25,dtype='float32')
>
> print len(lats_out)
> print len(lons_out)
>
> # output data.
> press_out = 900. + arange(nlevs*nlats*nlons,dtype='float32') # 1d array
> press_out.shape = (nlevs,nlats,nlons) # reshape to 2d array
> #temp_out = 9. + arange(nlevs*nlats*nlons,dtype='float32') # 1d array
> temp_out = tempOut
> temp_out.shape = (nlevs,nlats,nlons) # reshape to 2d array
>
> # create the lat and lon dimensions.
> ncfile.createDimension('latitude',nlats)
> ncfile.createDimension('longitude',nlons)
>
> # create level dimension.
> ncfile.createDimension('level',nlevs)
>
> # create time dimension (record, or unlimited dimension)
> ncfile.createDimension('time',None)
>
> # Define the coordinate variables. They will hold the coordinate
> # information, that is, the latitudes and longitudes.
> # Coordinate variables only given for lat and lon.
> lats = ncfile.createVariable('latitude',dtype('float32').char,('latitude',))
> lons =
> ncfile.createVariable('longitude',dtype('float32').char,('longitude',))
>
> # Assign units attributes to coordinate var data. This attaches a
> # text attribute to each of the coordinate variables, containing the
> # units.
> lats.units = 'degrees_north'
> lons.units = 'degrees_east'
>
> # write data to coordinate vars.
> lats[:] = lats_out
> lons[:] = lons_out
>
> # create the pressure and temperature variables
> press =
> ncfile.createVariable('pressure',dtype('float32').char,('time','level','latitude','longitude'))
> temp =
> ncfile.createVariable('temperature',dtype('float32').char,('time','level','latitude','longitude'))
>
> # set the units attribute.
> press.units = 'hPa'
> temp.units = 'celsius'
>
> # write data to variables along record (unlimited) dimension.
> # same data is written for each record.
> for nrec in range(nrecs):
> press[nrec,:,::] = press_out
> temp[nrec,:,::] = temp_out
>
> # close the file.
> ncfile.close()
> print '*** SUCCESS writing example file pres_temp_4D.nc'
>
> _______________________________________________
> netcdfgroup mailing list
> netcdfgroup@xxxxxxxxxxxxxxxx
> For list information or to unsubscribe, visit:
> http://www.unidata.ucar.edu/mailing_lists/
--
Dr. Richard P. Signell (508) 457-2229
USGS, 384 Woods Hole Rd.
Woods Hole, MA 02543-1598
netcdfgroup archives: