4
$\begingroup$

I want to regrid GCM data available on Gaussian grids to regular grids, say at 1 deg X 1 deg. How to do this using CDO or Python or MATLAB.

Edit 1: Data format is netcdf.

$\endgroup$
  • 1
    $\begingroup$ For the usage of CDO: Do you have your data available as netCDF files? How is your output grid exactly defined? The easiest way with cdo: create a grid description file of the target grid (with griddes from an existing target grid netCDF file) and apply an interpolation routine (e.g. remapbil) on your files that should be converted. $\endgroup$ – daniel.neumann May 30 '16 at 12:45
  • 1
    $\begingroup$ With respect to Matlab: The second and third search result (when I just searched it) provide the interp2 and griddedInterpolant functions which seem to be what you need. Please provide more information if you need further help. $\endgroup$ – daniel.neumann May 30 '16 at 15:53
  • $\begingroup$ @daniel.neumann thanks for your valuable suggestion, my data is in netcdf format. I am not sure whether interp2 can be applied on Gussian Grids though I have used this function heavily for rectilinear grids. $\endgroup$ – Mario May 31 '16 at 6:39
9
$\begingroup$

climate data operators (CDO)

define grid

We define a lat-lon target grid with 1°x1° grid cell size 30x30 grid cells starting at 40°N and -10°E (=10°W):

gridtype = lonlat 
xsize    =  30 
ysize    =  30 
xfirst   = -10 
xinc     =   1 
yfirst   =  40 
yinc     =   1

This text is written into a text file. See section 1.3.2 CDO Manual for details and further examples.

If you already have a netCDF file which has data on your target grid, you can also extract the grid definition from that file:

cdo griddes FILE_WITH_TARGET_GRID.nc > myGridDef

myGridDef is a text file.

interpolate data

I assume that we do a bilinear interpolation. This is done with the remapbil operator via:

cdo remapbil,myGridDef INPUT_FILE.nc OUTPUT_FILE.nc

The grid of your input data needs to be properly defined in your input file (that is often a problem).

'properly defined' grid

Alternative 1

The spatial dimensions of the data variable (SST in the example below) need to be named lon and lat. Additionally, coordinate variables lon and lat have to exist (time independent!). These need to have correct coordinate variable attributes (see example below).

Example:

dimensions:
    lon = 30 ;
    lat = 30 ;
    time = UNLIMITED ; // (24 currently)
variables:
    double lon(lon) ;
            lon:standard_name = "longitude" ;
            lon:long_name = "longitude" ;
            lon:units = "degrees_east" ;
            lon:axis = "X" ;
    double lat(lat) ;
            lat:standard_name = "latitude" ;
            lat:long_name = "latitude" ;
            lat:units = "degrees_north" ;
            lat:axis = "Y" ;
    double time(time) ;
            time:standard_name = "time" ;
            time:long_name = "time" ;
            time:units = "seconds since 1900-01-01 00:00:00" ;
            time:calendar = "standard" ;
            time:axis = "T" ;
    float SST(time, lat, lon) ;
            SST:long_name = "sea_surface_temperature" ;
            SST:units = "degree celsius" ;
            SST:_FillValue = NaNf ;
            SST:missing_value = NaNf ;
            SST:var_desc = "sea surface temperature" ;

Alternative 2

The spatial dimensions and corresponding variables exist in the source file but are not denoted as lon and lat. Then each data variable (SST, in this example) needs to have an attribute coordinates that holds the names of the coordinate variables.

SST:coordinates = "lon lat" ;

Example:

dimensions:
    TSTEP = UNLIMITED ; // (24 currently)
    COL = 112 ;
    ROW = 106 ;
variables:                     " ;
    double lon(ROW, COL) ;
        lon:standard_name = "longitude" ;
        lon:long_name = "longitude coordinate" ;
        lon:units = "degrees_east" ;
        lon:_CoordinateAxisType = "Lon" ;
    double lat(ROW, COL) ;
        lat:standard_name = "latitude" ;
        lat:long_name = "latitude coordinate" ;
        lat:units = "degrees_north" ;
        lat:_CoordinateAxisType = "Lat" ;
    double time(TSTEP) ;
        time:standard_name = "time" ;
        time:long_name = "time" ;
        time:units = "seconds since 1900-01-01 00:00:00" ;
        time:calendar = "standard" ;
    float SST(TSTEP, ROW, COL) ;
        SST:long_name = "sea_surface_temperature" ;
        SST:units = "degree celsius" ;
        SST:coordinates = "lon lat" ;
        SST:var_desc = "sea surface temperature" ;

Note:

If the coordinate variables are missing, then they can be created via the setgrid operator. If we want to add a grid definition to the source file then we need the source file's grid definition in the text file (see 'define grid').

cdo setgrid,mySourceGridDef INFILE_NO_COORDS.nc INFILE_WI_COORDS.nc
$\endgroup$
  • $\begingroup$ Is there an equally simple way to proceed when setting coordinates for a projected coordinate system? I didn't understand much from the documentation. $\endgroup$ – Janina Nov 2 '17 at 9:18
  • $\begingroup$ @Janina You mean when you have some projection, which does not project your model grid on a regular lat-lon grid? Do you mean interpolation or just setting a grid definition? Better you ask a new question (refer to this one) including an example of your data set. $\endgroup$ – daniel.neumann Nov 2 '17 at 9:48
  • $\begingroup$ yes, exactly, something like Lambert Conformal Conic etc. Essentially, I wanted to set the grid definition, so that for example in Panoply it would plot the data in the correct geographic location. And yes, ultimately I'm interested in regridding. Ok, I guess I'll turn this into a question... :-) $\endgroup$ – Janina Nov 2 '17 at 10:19
  • $\begingroup$ @Janina I just re-read my answer. It actually solves your issue if you already have your coordinates available in one netCDF file. $\endgroup$ – daniel.neumann Nov 2 '17 at 10:31
  • 1
    $\begingroup$ You are right, actually. After some trial-and-error progress I was able to do the regridding successfully - thank you for your detailed description! But the validity of my workflow is questionable. I created and saved the data in R with original coordinates in EPSG:3333. I solved the georeferencing problem by adding projection metadata. I did the transformation to WGS84 to get the cell centre coordinates in lat-long and saved them to variables lon and lat (just like your alternative 2). The output of griddes for this file was a curvilinear grid. $\endgroup$ – Janina Nov 3 '17 at 10:01
0
$\begingroup$

Using Python, the best solution I've been able to find is pyresample. See this blog post.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.