CHAPTER 4 Regridding Data

Overview

CDMS has functions to interpolate gridded data:

from one horizontal (lat/lon) grid to another
from one set of pressure levels to another
from one vertical (lat/level) cross-section to another vertical cross-section.

Horizontal regridder

The simplest method to regrid a variable from one horizontal, lat/lon grid to another is to use the regrid function defined for variables. This function takes the target grid as an argument, and returns the variable regridded to the target grid:

>>> import cdms

>>> f = cdms.open('/pcmdi/cdms/exp/cmip2/ccc/perturb.xml')

>>> rlsf = f('rls') # Read the data

>>> rlsf.shape

(4, 48, 96)

>>> g = cdms.open('/pcmdi/cdms/exp/cmip2/mri/perturb.xml')

>>> rlsg = g['rls'] # Get the file variable (no data read)

>>> outgrid = rlsg.getGrid() # Get the target grid

>>> rlsnew = rlsf.regrid(outgrid)

>>> rlsnew.shape

(4, 46, 72)

>>> outgrid.shape

(46, 72)

A somewhat more efficient method is to create a regridder function. This has the advantage that the mapping is created only once and can be used for multiple arrays. Also, this method can be used with data in the form of an MA.MaskedArray or Numeric array. The steps in this process are:

Given an input grid and output grid, generate a regridder function.
Call the regridder function on a Numeric array, resulting in an array defined on the output grid. The regridder function can be called with any array or variable defined on the input grid.

The following example illustrates this process. The regridder function is generated at line 9, and the regridding is performed at line 10:

1 #!/usr/bin/env python

2 import cdms

3 from regrid import Regridder

4 f = cdms.open('/pcmdi/cdms/exp/cmip2/ccc/perturb.xml')

5 rlsf = f['rls']

6 ingrid = rlsf.getGrid()

7 g = cdms.open('/pcmdi/cdms/exp/cmip2/mri/perturb.xml')

8 outgrid = g['rls'].getGrid()

9 regridfunc = Regridder(ingrid, outgrid)

10 rlsnew = regridfunc(rlsf)

11 f.close()

12 g.close()

Line	Notes
2	Makes the CDMS module available.
3	Makes the Regridder class available from the regrid module.
4	Opens the input dataset.
5	Gets the variable object named `rls' . No data is read.
6	Gets the input grid.
7	Opens a dataset to retrieve the output grid.
8	The output grid is the grid associated with the variable named `rls' in dataset g . Just the grid is retrieved, not the data.
9	Generates a regridder function regridfunc .
10	Reads all data for variable rlsf , and calls the regridder function on that data, resulting in a transient variable rlsnew .

Pressure-level regridder

To regrid a variable which is a function of latitude, longitude, pressure level, and (optionally) time to a new set of pressure levels, use the pressureRegrid function defined for variables. This function takes an axis representing the target set of pressure levels, and returns a new variable d regridded to that dimension.

>>> var.shape

(3, 16, 32)

>>> var.getAxisIds()

['level', 'latitude', 'longitude']

>>> len(levout)

>>> result = var.pressureRegrid(levout)

>>> result.shape

(2, 16, 32)

Cross-section regridder

To regrid a variable which is a function of latitude, height, and (optionally) time to a new latitude/height cross-section, use the crossSectionRegridder defined for variables. This function takes as arguments the new latitudes and heights, and returns the variable regridded to those axes.

>>> varin.shape

(11, 46)

>>> varin.getAxisIds()

['level', 'latitude']

>>> levOut[:]

[ 10., 30., 50., 70., 100., 200., 300., 400., 500., 700., 850.,

1000.,]

>>> varout = varin.crossSectionRegrid(levOut, latOut)

>>> varout.shape

(12, 64)

Regridder Constructor

regridFunction = Regridder(inputGrid, outputGrid) Create a regridder function which interpolates a data array from input to output grid. Table 4.2 describes the calling sequence of this function. inputGrid and outputGrid are CDMS grid objects. Note: To set the mask associated with inputGrid or outputGrid, use the grid setMask function.

regridder functions

A regridder function is an instance of the Regridder class. The function is associated with an input and output grid. Typically its use is straightforward: the function is passed an input array and returns the regridded array. However, when the array has missing data, or the input and/or output grids are masked, the logic becomes more complicated.

Step 1: The regridder function first forms an input mask. This mask is either two-dimensional or `n-dimensional', depending on the rank of the user-supplied mask. If no mask or missing value is specified, the mask is obtained from the data array mask if present.

Two-dimensional case:

Let mask_1 be the two-dimensional user mask supplied via the mask argument, or the mask of the input grid if no user mask is specified.
If a missing-data value is specified via the missing argument, let the implicit_mask be the two-dimensional mask defined as 0 where the first horizontal slice of the input array is missing, 1 elsewhere.
The input mask is the logical AND(mask_1, implicit_mask)

N-dimensional case: If the user mask is 3 or 4-dimensional with the same shape as the input array, it is used as the input mask.

Step 2: The data is then regridded. In the two-dimensional case, the input mask is `broadcast' across the other dimensions of the array. In other words, it assumes that all horizontal slices of the array have the same mask. The result is a new array, defined on the output grid. Optionally, the regridder function can also return an array having the same shape as the output array, defining the fractional area of the output array which overlaps a non-missing input grid cell. This is useful for calculating area-weighted means of masked data.

Step 3: Finally, if the output grid has a mask, it is applied to the result array. Where the output mask is 0, data values are set to the missing data value, or 1.0e20 if undefined. The result array or transient variable will have a mask value of 1 (invalid value) for those output grid cells which completely overlap input grid cells with missing values.

Regridder function
Type
Array or TransientVariable	regridFunction(array, missing=None, order=None, mask=None) Interpolate a gridded data array to a new grid. The interpolation preserves the area-weighted mean on each horizontal slice. If array is a Variable, a TransientVariable of the same rank as the input array is returned, otherwise a masked array is returned. array is a Variable, masked array, or Numeric array of rank 2, 3, or 4. missing is a Float specifying the missing data value. The default is 1.0e20. order is a string indicating the order of dimensions of the array. It has the form returned from variable.getOrder(). For example, the string "tzyx" indicates that the dimension order of array is (time, level, latitude, longitude). If unspecified, the function assumes that the last two dimensions of array match the input grid. mask is a Numeric array, of datatype Integer or Float, consisting of ones and zeros. A value of 0 or 0.0 indicates that the corresponding data value is to be ignored for purposes of regridding. If mask is two-dimensional of the same shape as the input grid, it overrides the mask of the input grid. If the mask has more than two dimensions, it must have the same shape as array. In this case, the missing data value is also ignored. Such an n-dimensional mask is useful if the pattern of missing data varies with level (e.g., ocean data) or time. Note: If neither missing or mask is set, the default mask is obtained from the mask of the array if any.
Array, Array	regridFunction(ar, missing=None, order=None, mask=None, returnTuple=1) If called with the optional returnTuple argument equal to 1, the function returns a tuple (dataArray, maskArray). dataArray is the result data array. maskArray is a Float32 array of the same shape as dataArray, such that maskArray[i,j] is fraction of the output grid cell [i,j] overlapping a non-missing cell of the input grid.

Examples

Example: Regrid data to a uniform output grid.

1 #!/usr/local/bin/python

2 import cdms

3 from regrid import Regridder

4 f = cdms.open('rls_ccc_per.nc')

5 rlsf = f.variables['rls']

6 ingrid = rlsf.getGrid()

7 outgrid = cdms.createUniformGrid(90.0, 46, -4.0, 0.0, 72, 5.0)

8 regridFunc = Regridder(ingrid, outgrid)

9 newrls = regridFunc(rlsf)

10 f.close()

Line	Notes
4	Open a netCDF file for input.
7	Create a 4 x 5 degree output grid. Note that this grid is not associated with a file or dataset
8	Create the regridder function
9	Read all data and regrid. The missing data value is obtained from variable rlsf.

Example: Get a mask from a separate file, and set as the input grid mask.

1 import cdms

2 from regrid import Regridder

3 #

4 f = cdms.open('so_ccc_per.nc')

5 sof = f.variables['so']

6 ingrid = sof.getGrid()

7 g = cdms.open('rls_mri_per.nc')

8 rlsg = g.variables['rls']

9 outgrid = rlsg.getGrid()

10 regridFunc = Regridder(ingrid,outgrid)

11 h = cdms.open('sft_ccc.nc')

12 sfmaskvar = h.variables['sfmask']

13 sfmask = sfmaskvar[:]

14 outArray = regridFunc(sof.subSlice(time=0),mask=sfmask)

15 f.close()

16 g.close()

17 h.close()

Line	Notes
6	Get the input grid.
9	Get the output grid
10	Create the regridder function.
13	Get the mask.
14	Regrid with a user mask. The subslice call returns a transient variable corresponding to variable sof at time 0. Note: Although it cannot be determined from the code, both mask and the input array sof are four-dimensional. This is the `n-dimensional' case.

Example: Generate an array of zonal mean values.

1 f = cdms.open('rls_ccc_per.nc')

2 rlsf = f.variables['rls']

3 ingrid = rlsf.getGrid()

4 outgrid = cdms.createZonalGrid(ingrid)

5 regridFunc = Regridder(ingrid,outgrid)

6 mean = regridFunc(rlsf)

7 f.close()

Line	Notes
3	Get the input grid.
4	Create a zonal grid. outgrid has the same latitudes as ingrid , and a singleton longitude dimension. createGlobalMeanGrid could be used here to generate a global mean array.
5	Generate the regridder function.
6	Generate the zonal mean array.

Example: Regrid an array with missing data, and calculate the area-weighted mean of the result.

1 from cdms.MV import *

...

2 outgrid = cdms.createUniformGrid(90.0, 46, -4.0, 0.0, 72, 5.0)

3 outlatw, outlonw = outgrid.getWeights()

4 outweights = outerproduct(outlatw, outlonw)

5 grid = var.getGrid()

6 sample = var[0,0]

7 latw, lonw = grid.getWeights()

8 weights = outerproduct(latw, lonw)

9 inmask = where(greater(absolute(sample),1.e15),0,1)

10 mean = add.reduce(ravel(inmask*weights*sample))/add.reduce(ravel(inmask*weights))

11 regridFunc = Regridder(grid, outgrid)

12 outsample, outmask = regridFunc(sample, mask=inmask, returnTuple=1)

13 outmean = add.reduce(ravel(outmask*outweights*outsample))/add.reduce(ravel(outmask*outweights))

Line	Notes
2	Create a uniform target grid.
3	Get the latitude and longitude weights.
4	Generate a 2-D weights array.
5	Get the input grid. var is a 4-D variable.
6	Get the first horizontal slice from var .
7-8	Get the input weights, and generate a 2-D weights array.
9	Set the 2-D input mask.
10	Calculate the input array area-weighted mean.
11	Create the regridder function.
12	Regrid. Because returnTuple is set to 1, the result is a tuple (dataArray, maskArray).
13	Calculate the area-weighted mean of the regridded data. mean and outmean should be approximately equal.