Problems in regridding from high resolution (0.5x0.625) to low resolution (4x5) for rectilinear grids

[uzzzaldash]

Hi, I'm new to xESMF, trying to regridding from high resolution to low resolution for rectilinear grids following https://xesmf.readthedocs.io/en/latest/notebooks/Rectilinear_grid.html. First, I tried to regrid from GEOS-Chem model data (Res:4x5) to CrIS Satellite data (Res:0.5x0.625). And it seems worked perfectly as the spatial maps for original data and regridded data was same. Then I tried to do opposite which is regidded from CrIS Satellite data (Res:0.5x0.625) to GEOS-Chem model data (Res:4x5). In this case, I just changed the input and output under the xe.Regridder function and plotted the spatial maps. However the maps generated from original and regridded data are not look same in this case. I am giving my code and figures below:

# Import required libraries
%matplotlib inline
import numpy as np
import xarray as xr
import matplotlib.pyplot as plt
import cartopy.crs as ccrs  #Coordinate reference system
import pandas as pd
import xesmf as xe

Import in a Sample GEOS-Chem 4x5 Model Run

GC_path = '.../OutputDir/GEOSChem.SatDiagn.20190701_0000z.nc4'
GC4x5 = xr.open_dataset(GC_path)

# Let's just take the surface layer of data. We're going to regrid surface level isoprene mixing ratio
GC4x5_ISOP_surflev = GC4x5['SatDiagnConc_ISOP'].isel(lev=0)  #time: 7 lat: 46 lon: 72

Import in Daily CrIS Data

CrIS_ISOP_path = '.../CrIS_ISOP_data/201901_daily_CrIS_Isoprene_screened.nc'
CrIS_ISOP = xr.open_dataset(CrIS_ISOP_path, decode_times=True)

# Need to get time as an xarray coordinate
# Bcz GC dim [time, lat, lon] but CrIS dim [day, lat, lon] 
temp_staging_area = pd.to_datetime(CrIS_ISOP['day']+1,format='%d') + pd.DateOffset(hour=12, month=1, years=2019-1900)
CrIS_ISOP['day'] = temp_staging_area

# Rename CrIS coordinate of day to time to match with GEOS-Chem naming convention
CrIS_ISOP = CrIS_ISOP.rename({'day':'time'})
CrIS_ISOP['time'].attrs['long_name'] = 'Time'
CrIS_ISOP['time'].attrs['axis'] = 'T'
CrIS_ISOP    #time: 31 lat: 361 lon: 576

Setting up the CrIS grid (0.5x0.625)

#CrIS_ISOP['lat'] #[-89.875, -89.5, -89., ...,89., 89.5, 89.875] #361
#CrIS_ISOP['lon'] # [-180.   , -179.375, -178.75 , ...,  178.125,  178.75 ,  179.375] #576

# Midpoints of latitude coordinate in CrIS data
a = np.arange(-89.5, 90, 0.5) # [-89.5, -89. , -88.5,..89.5]
b = np.append(a, 89.875)
c = np.insert(b, 0, -89.875) #[-89.875, -89.5  , -89.   , -88.5,...88.5  ,  89.   ,  89.5  ,  89.875]

# Boundaries of latitude coordinate in CrIS data
d = np.linspace(-89.5-0.5/2, 89.5+0.5/2, 360) #[-89.75, -89.25, -88.75,...,88.75,  89.25,  89.75]]
e = np.append(d, 90)
f = np.insert(e, 0, -90) # [-90.  , -89.75, -89.25,...89.25, 89.75,  90.]

CrIS_grid_with_bounds = {'lon': np.linspace(-180,179.375,576),
                         'lat': c,
                         'lon_b': np.linspace(-180-0.625/2,179.375+0.625/2,577),
                         'lat_b':f
                        }

# lon[-180., -179.375, -178.75, -178.125,...,177.5, 178.125, 178.75 , 179.375] #576
# lon_b[-180.3125, -179.6875, -179.0625, -178.4375,...,177.8125, 178.4375, 179.0625, 179.6875] #577
# lat[-89.875, -89.5, -89., -88.5,...,88.5, 89., 89.5, 89.875] #361
# lat_b[-90., -89.75, -89.25, -88.75,...,88.75, 89.25,89.75, 90.] #362

Setting up Grid for GEOS-Chem 4x5

#GC4x5['lat'] # [-89., -86., -82., -78.,..., 78., 82., 86., 89.] #46
#GC4x5['lon'] # [-180., -175., -170.,..., 165., 170., 175.] #72

# Midpoints of latitude coordinate in GC 4x5
a = np.arange(-86, 90, 4)
b = np.append(a, 89)
c = np.insert(b, 0, -89) #Use this variable

# Boundaries of latitude coordinate in GC 4x5
d = np.arange(-88, 90, 4)
e = np.append(d,90)
f = np.insert(e,0,-90) #Use this variable

# Only valid for HEMCO standalone grid
GC_standalone_grid_with_bounds = {'lon': np.arange(-180,180,5),
                                  'lat': c,
                                  'lon_b': np.arange(-182.5, 180,5),
                                  'lat_b':f
                                 }

#lon[-180, -175, -170, -165,...,  160,  165,  170,  175] 
#lat[-89, -86, -82, -78, -74,...,  78,  82,  86,  89]
#lat_b[-90, -88, -84, -80,...,  76,  80,  84,  88,  90]
#lon_b[-182.5, -177.5, -172.5, -167.5,...,  162.5,  167.5,  172.5,177.5]

Make a Regridding Object

regridder = xe.Regridder(CrIS_grid_with_bounds, GC_standalone_grid_with_bounds, 'conservative_normed', periodic=True)
regridder

Output:
xESMF Regridder 
Regridding algorithm:       conservative_normed 
Weight filename:            conservative_normed_361x576_46x72.nc 
Reuse pre-computed weights? False 
Input grid shape:           (361, 576) 
Output grid shape:          (46, 72) 
Periodic in longitude?      False

Perform the Regridding from 0.5x0.625 to 4x5

CrIS_ISOP_regridded = regridder(CrIS_ISOP, keep_attrs=True)

Plot rigridded data

CrIS_ISOP_regridded.isel(time=0).plot()

Plot Original data

CrIS_ISOP['Isop'].isel(time=0).plot()

Here, the first figure is for regridded data and the second one is for the original data. But they are not same. But when I regidded from low to high resolution (4x5 to 0.5x0.625) in the same way, I got the exactly same figure for original and regridded data. What I am missing here then?

[raphaeldussin]

hi @uzzzaldash

this is an interesting application! My understanding is that your input data is not defined everywhere and the resulting field looks like a dilution of the original. Have you tried to define a mask for your input data? you can add a variable "mask" in your dataset that you could define like ds["mask"] = xr.where(np.isnan(data), 0, 1)

[uzzzaldash]

Hi @raphaeldussin , Thanks for the response. i'm not fully sure what do you mean by 'defining my input data'. So, I am telling you what did after getting your response.

CrIS_ISOP is my input dataset. It has actually two variables. 1. Isop 2.Ndata
So, modified my code like:

CrIS_ISOP_path = '.../CrIS_ISOP_data/201907_daily_CrIS_Isoprene_screened.nc'
CrIS_ISOP_data = xr.open_dataset(CrIS_ISOP_path, decode_times=True) 
CrIS_ISOP = CrIS_ISOP_data['Isop']

This what I understood by defining input data.

Then I added your suggested line at the end of section named Import in Daily CrIS Data which is 3rd paragraph in my code. So, the modified code is like:

Import in Daily CrIS Data

CrIS_ISOP_path = '.../CrIS_ISOP_data/201901_daily_CrIS_Isoprene_screened.nc'
CrIS_ISOP = xr.open_dataset(CrIS_ISOP_path, decode_times=True)

# Need to get time as an xarray coordinate
# Bcz GC dim [time, lat, lon] but CrIS dim [day, lat, lon] 
temp_staging_area = pd.to_datetime(CrIS_ISOP['day']+1,format='%d') + pd.DateOffset(hour=12, month=1, years=2019-1900)
CrIS_ISOP['day'] = temp_staging_area

# Rename CrIS coordinate of day to time to match with GEOS-Chem naming convention
CrIS_ISOP = CrIS_ISOP.rename({'day':'time'})
CrIS_ISOP['time'].attrs['long_name'] = 'Time'
CrIS_ISOP['time'].attrs['axis'] = 'T'
CrIS_ISOP    #time: 31 lat: 361 lon: 576

CrIS_ISOP["mask"] = xr.where(np.isnan(CrIS_ISOP),0,1)

After doing these, The figures I am getting still the same. Could you please tell me little detail about it?