Regridding from a curvilinear to rectilinear grid

756 views
Skip to first unread message

Damien Irving

unread,
May 5, 2016, 8:00:15 PM5/5/16
to Iris
I'm working with CMIP5 ocean data, and many of the models use a curvilinear grid (i.e. where the horizontal spatial coordinates are two dimensional). For example:

print cube

sea_water_potential_temperature                  
(time: 60; depth: 50; cell index along second dimension: 300; cell index along first dimension: 360)
     
Dimension coordinates:
          time                                         x          
-                                      -                                      -
          depth                                        
-          x                                      -                                      -
          cell index along second dimension            
-          -                                      x                                      -
          cell index along first dimension            
-          -                                      -                                      x
     
Auxiliary coordinates:
          latitude                                    
-          -                                      x                                      x
          longitude                                    
-          -                                      x                                      x


I want to regrid/interpolate cubes like this to a regular rectilinear lat/lon grid (i.e. where the horizontal spatial coordinates are one dimensional). The regular regrid and interpolate methods in Iris don't seem to be able to handle 2D coordinates, so it looks like the following experimental function is the way to go:
 
import iris, numpy
from iris.experimental.regrid import regrid_weighted_curvilinear_to_rectilinear

def make_grid(lat_values, lon_values):
   
"""Make a dummy cube with desired grid."""
       
    latitude
= iris.coords.DimCoord(lat_values,
                                    standard_name
='latitude',
                                    units
='degrees_north',
                                    coord_system
=None)
    longitude
= iris.coords.DimCoord(lon_values,                    
                                     standard_name
='longitude',
                                     units
='degrees_east',
                                     coord_system
=None)

    dummy_data
= numpy.zeros((len(lat_values), len(lon_values)))
    new_cube
= iris.cube.Cube(dummy_data, dim_coords_and_dims=[(latitude, 0), (longitude, 1)])
   
   
return new_cube

src_cube
= cube[0, 0]
weights
= numpy.ones(src_cube.shape)

lats
= numpy.arange(-90, 91, 1)
lons
= numpy.arange(0, 360, 1)
target_grid_cube
= make_grid(lats, lons)
target_grid_cube
.coord('longitude').guess_bounds()
target_grid_cube
.coord('latitude').guess_bounds()

regridded_cube
= regrid_weighted_curvilinear_to_rectilinear(src_cube, weights, target_grid_cube)

Since regrid_weighted_curvilinear_to_rectilinear only handles 2D fields, I'd then have to iterate over each ['time', 'depth'] slice to regrid the entire cube.

From plotting the resulting the regridded cube it looks like everything worked just fine, however since that function is experimental I just wanted to check I was implementing it correctly? Or do others take a different approach to regridding curvilinear grids?




Denis Sergeev

unread,
Jun 23, 2017, 10:20:13 AM6/23/17
to Iris

Not actually answering your question, but on a related note:


I'm currently trying to regrid Arctic System Reanalysis from North Polar Stereographic coordinates to a regular lat/lon grid.


The subset of original data looks like this:

Temperature / (K)                   (Time: 128; pressure: 10; y_coordinate: 180; x_coordinate: 200)
     
Dimension coordinates:
         
Time                           x              -                 -                  -
          pressure                      
-              x                 -                  -
          y_coordinate                  
-              -                 x                  -
          x_coordinate                  
-              -                 -                  x
     
Auxiliary coordinates:

          latitude                      
-              -                 x                  x
          longitude                      
-              -                 x                  x


Since the coordinate system is not recognised automatically, I define it as following:

asr_coord_system = iris.coord_systems.Stereographic(90, -175)


Then, I'm using almost the same approach to interpolate the data to rectilinear grid:

def make_grid(lat_values, lon_values):
   
"""Make a dummy cube with desired grid."""

    geogcs
= iris.coord_systems.GeogCS(iris.fileformats.pp.EARTH_RADIUS)


    latitude
= iris.coords.DimCoord(lat_values,

                                    standard_name
='latitude',

                                    units
='degrees_north',


                                    coord_system
=geogcs)


    longitude
= iris.coords.DimCoord(lon_values,
                                     standard_name
='longitude',
                                     units
='degrees_east',

                                     coord_system
=geogcs)

    dummy_data
= np.zeros((len(lat_values), len(lon_values)))

    new_cube
= iris.cube.Cube(dummy_data, dim_coords_and_dims=[(latitude, 0), (longitude, 1)])

   
return new_cube

src_cube
= cube[0, 0]

weights
= np.ones(src_cube.shape)
lats
= np.arange(65, 91, 1)
lons
= np.arange(-20, 25, 1)

target_grid_cube
= make_grid(lats, lons)
target_grid_cube
.coord('longitude').guess_bounds()
target_grid_cube
.coord('latitude').guess_bounds()

regridded_cube
= regrid_weighted_curvilinear_to_rectilinear(src_cube, weights, target_grid_cube)


However, I get strange results with only one row of data, while the rest of the array is masked.


Does anyone know what I am doing wrong? Can it be done in iris at all?


Cheers,

Denis

Holly A

unread,
Feb 19, 2018, 11:53:25 AM2/19/18
to SciTools (iris, cartopy, cf_units, etc.) - https://github.com/scitools

Hi,
 
So not actually an answer to this question, but I am having the exact same problems as I need to do multi-model analysis of sea ice and ocean CMIP5 data.

I was wondering if you have found a way to do this as it seems you have had a similar problem?

Sorry I cannot help you with this!

Thank you

Damien Irving

unread,
Feb 22, 2018, 3:27:41 PM2/22/18
to SciTools (iris, cartopy, cf_units, etc.) - https://github.com/scitools
Hi Holly,

I've been keeping track of what I've found out about regridding vector quantities at the following GitHub issue:

So far I haven't implemented any of the ideas listed (I've just put it in the too hard basket for now), but hopefully they are helpful for you :-)


Cheers,
Damien
Reply all
Reply to author
Forward
0 new messages