Interpolation between two xarray datasets with Basemap

Question

I have two different xarray datasets that have different latitude/longitude grid resolutions. I want to regrid the one xarray with lower resolution to the same resolution as the one xarray with higher resolution. I found some examples (e.g., http://earthpy.org/interpolation_between_grids_with_basemap.html), but it does not work for me. Here is one example that I made for testing:

import numpy as np
import xarray as xray
import mpl_toolkits.basemap

var1=xray.DataArray(np.random.randn(len(np.linspace(40.5,49.5,10)),len(np.linspace(-39.5,-20.5,20))),coords=[np.linspace(40.5,49.5,10), np.linspace(-39.5,-20.5,20)],dims=['lat','lon'])

(xlon, xlat)=np.meshgrid(np.linspace(-39.875,-20.125,80),np.linspace(40.125,49.875,40))
var2=xray.DataArray(-xlon**2+xlat**2,coords=[np.linspace(40.125,49.875,40),np.linspace(-39.875,-20.125,80)],dims=['lat','lon'])

mpl_toolkits.basemap.interp(var1,var1.lon,var1.lat,var2.lon,var2.lat,checkbounds=False,masked=False,order=0)

I get following error:

ValueError: xout and yout must have same shape!

Screenshot: see screenshot

Does basemap.interp() require xout and yout to be the same shape? So var2 needs to be a square? This is almost never the case with any of my datasets! How can I regrid var1 to be the same resolution as var2?

Note: After regridding, I want to subsample var1 given some condition related to var2. For example:

var1_subset = var1.where(var2>1000)

So I want to minimize any loss of grid points during the interpolation.

suvy · Accepted Answer · 2017-08-10 21:18:56Z

basemap.interp will work only when xout and yout are same in number or number of output nlons and nlats are same,

why not generate same length output nlats and nlons and subset it later.

For example:

import numpy as np
import xarray as xray 
import mpl_toolkits.basemap
var1=xray.DataArray(np.random.randn(len(np.linspace(40.5,49.5,10)),len(np.linspace(-39.5,-20.5,20))),coords=[np.linspace(40.5,49.5,10), np.linspace(-39.5,-20.5,20)],dims=['lat','lon'])

(xlon,xlat)=np.meshgrid(np.linspace(-39.875,20.125,80),np.linspace(40.125,49.875,80))
var2=xray.DataArray(-xlon**2+xlat**2,coords[np.linspace(40.125,49.875,80),np.linspace(-39.875,-20.125,80)],dims=['lat','lon'])
mpl_toolkits.basemap.interp(var1,var1.lon,var1.lat,var2.lon,var2.lat,checkbounds=False,masked=False,order=0)

Here is another cool trick with xarray.

lonreg=var1.groupby_bins('lon',np.linspace(-39.875,20.125,80)).mean(dim='lon')

regridded=lonreg.groupby_bins('lat',np.linspace(-39.5,20.5,20)).mean(dim='lat')

if you want weighted averaged regridding, it is easy to extend this for area averaged regridding by using weights and sum function on groupby object.

Collectives™ on Stack Overflow

Interpolation between two xarray datasets with Basemap

1 Answer 1

Comments

Your Answer

Hot Network Questions

Collectives™ on Stack Overflow

1 Answer 1

Comments

Your Answer

Sign up or log in

Post as a guest

Related