Internal API
frontend
Frontend for xESMF, exposed to users.
- xesmf.frontend.subset_regridder(ds_out, ds_in, method, in_dims, out_dims, locstream_in, locstream_out, periodic, **kwargs)
Compute subset of weights
- xesmf.frontend.ds_to_ESMFgrid(ds, need_bounds=False, periodic=None, append=None)
Convert xarray DataSet or dictionary to ESMF.Grid object.
- Parameters:
- dsxarray DataSet or dictionary
Contains variables
lon
,lat
, and optionallylon_b
,lat_b
if need_bounds=True.Shape should be
(n_lat, n_lon)
or(n_y, n_x)
, as normal C or Python ordering. Will be then tranposed to F-ordered.- need_boundsbool, optional
Need cell boundary values?
- periodicbool, optional
Periodic in longitude?
- Returns:
- grid
ESMF.Grid object
- shape
Shape of the grid
- dim_names
Dimension names of the grid
- xesmf.frontend.ds_to_ESMFlocstream(ds)
Convert xarray DataSet or dictionary to ESMF.LocStream object.
- Parameters:
- dsxarray DataSet or dictionary
Contains variables
lon
,lat
.
- Returns:
- locstreamESMF.LocStream object
- xesmf.frontend.polys_to_ESMFmesh(polys)
Convert a sequence of shapely Polygons to a ESMF.Mesh object.
MultiPolygons are split in their polygon parts and holes are ignored.
- Parameters:
- polyssequence of shapely Polygon or MultiPolygon
- Returns:
- exteriorESMF.Mesh
A mesh where elements are the exterior rings of the polygons
- tuple
The shape of the mesh : (1, N_elements)
backend
Backend for xESMF. This module wraps ESMPy’s complicated API and can create ESMF Grid and Regrid objects only using basic numpy arrays.
General idea:
1) Only use pure numpy array in this low-level backend. xarray should only be used in higher-level APIs which interface with this low-level backend.
2) Use simple, procedural programming here. Because ESMPy Classes are complicated enough, building new Classes will make debugging very difficult.
3) Add some basic error checking in this wrapper level. ESMPy is hard to debug because the program often dies in the Fortran level. So it would be helpful to catch some common mistakes in Python level.
- xesmf.backend.warn_f_contiguous(a)
Give a warning if input array if not Fortran-ordered.
ESMPy expects Fortran-ordered array. Passing C-ordered array will slow down performance due to memory rearrangement.
- Parameters:
- anumpy array
- xesmf.backend.warn_lat_range(lat)
Give a warning if latitude is outside of [-90, 90]
Longitute, on the other hand, can be in any range, since the it the transform is done in (x, y, z) space.
- Parameters:
- latnumpy array
- class xesmf.backend.Grid(*args: Any, **kwargs: Any)
- classmethod from_xarray(lon, lat, periodic=False, mask=None)
Create an ESMF.Grid object, for constructing ESMF.Field and ESMF.Regrid.
- Parameters:
- lon, lat2D numpy array
Longitute/Latitude of cell centers.
Recommend Fortran-ordering to match ESMPy internal.
Shape should be
(Nlon, Nlat)
for rectilinear grid, or(Nx, Ny)
for general quadrilateral grid.- periodicbool, optional
Periodic in longitude? Default to False. Only useful for source grid.
- mask2D numpy array, optional
Grid mask. According to the ESMF convention, masked cells are set to 0 and unmasked cells to 1.
Shape should be
(Nlon, Nlat)
for rectilinear grid, or(Nx, Ny)
for general quadrilateral grid.
- Returns:
- gridESMF.Grid object
- get_shape(loc=esmpy.StaggerLoc.CENTER)
Return shape of grid for specified StaggerLoc
- class xesmf.backend.LocStream(*args: Any, **kwargs: Any)
- classmethod from_xarray(lon, lat)
Create an ESMF.LocStream object, for contrusting ESMF.Field and ESMF.Regrid
- Parameters:
- lon, lat1D numpy array
Longitute/Latitude of cell centers.
- Returns:
- locstreamESMF.LocStream object
- get_shape()
Return LocStream shape.
- xesmf.backend.add_corner(grid, lon_b, lat_b)
Add corner information to ESMF.Grid for conservative regridding.
Not needed for other methods like bilinear or nearest neighbour.
- Parameters:
- gridESMF.Grid object
Generated by
Grid.from_xarray()
. Will be modified in-place.- lon_b, lat_b2D numpy array
Longitute/Latitude of cell corner Recommend Fortran-ordering to match ESMPy internal. Shape should be
(Nlon+1, Nlat+1)
, or(Nx+1, Ny+1)
- class xesmf.backend.Mesh(*args: Any, **kwargs: Any)
- classmethod from_polygons(polys, element_coords='centroid')
Create an ESMF.Mesh object from a list of polygons.
All exterior ring points are added to the mesh as nodes and each polygon is added as an element, with the polygon centroid as the element’s coordinates.
- Parameters:
- polyssequence of shapely Polygon
Holes are not represented by the Mesh.
- element_coordsarray or “centroid”, optional
If “centroid”, the polygon centroids will be used (default) If an array of shape (len(polys), 2) : the element coordinates of the mesh. If None, the Mesh’s elements will not have coordinates.
- Returns:
- meshESMF.Mesh
A mesh where each polygon is represented as an Element.
- get_shape(loc=esmpy.MeshLoc.ELEMENT)
Return the shape of the Mesh at specified MeshLoc location.
- xesmf.backend.esmf_regrid_build(sourcegrid, destgrid, method, filename=None, extra_dims=None, extrap_method=None, extrap_dist_exponent=None, extrap_num_src_pnts=None, ignore_degenerate=None)
Create an ESMF.Regrid object, containing regridding weights.
- Parameters:
- sourcegrid, destgridESMF.Grid or ESMF.Mesh object
Source and destination grids.
Should create them by
Grid.from_xarray()
(with optionallyadd_corner()
), instead of ESMPy’s original API.- methodstr
Regridding method. Options are
‘bilinear’
‘conservative’, need grid corner information
‘conservative_normed’, need grid corner information
‘patch’
‘nearest_s2d’
‘nearest_d2s’
- filenamestr, optional
Offline weight file. Require ESMPy 7.1.0.dev38 or newer. With the weights available, we can use Scipy’s sparse matrix multiplication to apply weights, which is faster and more Pythonic than ESMPy’s online regridding. If None, weights are stored in memory only.
- extra_dimsa list of integers, optional
Extra dimensions (e.g. time or levels) in the data field
This does NOT affect offline weight file, only affects online regrid.
Extra dimensions will be stacked to the fastest-changing dimensions, i.e. following Fortran-like instead of C-like conventions. For example, if extra_dims=[Nlev, Ntime], then the data field dimension will be [Nlon, Nlat, Nlev, Ntime]
- extrap_methodstr, optional
Extrapolation method. Options are
‘inverse_dist’
‘nearest_s2d’
- extrap_dist_exponentfloat, optional
The exponent to raise the distance to when calculating weights for the extrapolation method. If none are specified, defaults to 2.0
- extrap_num_src_pntsint, optional
The number of source points to use for the extrapolation methods that use more than one source point. If none are specified, defaults to 8
- ignore_degeneratebool, optional
If False (default), raise error if grids contain degenerated cells (i.e. triangles or lines, instead of quadrilaterals)
- Returns:
- gridESMF.Grid object
- xesmf.backend.esmf_regrid_apply(regrid, indata)
Apply existing regridding weights to the data field, using ESMPy’s built-in functionality.
xESMF use Scipy to apply weights instead of this. This is only for benchmarking Scipy’s result and performance.
- Parameters:
- regridESMF.Regrid object
Contains the mapping from the source grid to the destination grid.
Users should create them by esmf_regrid_build(), instead of ESMPy’s original API.
- indatanumpy array of shape
(Nlon, Nlat, N1, N2, ...)
Extra dimensions
(N1, N2, ...)
are specified inesmf_regrid_build()
.Recommend Fortran-ordering to match ESMPy internal.
- Returns:
- outdatanumpy array of shape
(Nlon_out, Nlat_out, N1, N2, ...)
- outdatanumpy array of shape
- xesmf.backend.esmf_regrid_finalize(regrid)
Free the underlying Fortran array to avoid memory leak.
After calling
destroy()
on regrid or its fields, we cannot use the regrid method anymore, but the input and output data still exist.- Parameters:
- regridESMF.Regrid object
smm
Sparse matrix multiplication (SMM) using scipy.sparse library.
- xesmf.smm.read_weights(weights, n_in, n_out)
Read regridding weights into a DataArray (sparse COO matrix).
- Parameters:
- weightsstr, Path, xr.Dataset, xr.DataArray, sparse.COO
Weights generated by ESMF. Can be a path to a netCDF file generated by ESMF, an xarray.Dataset, a dictionary created by ESMPy.api.Regrid.get_weights_dict or directly the sparse array as returned by this function.
- n_in, n_outintegers
(N_out, N_in)
will be the shape of the returning sparse matrix. They are the total number of grid boxes in input and output grids:N_in = Nx_in * Ny_in N_out = Nx_out * Ny_out
We need them because the shape cannot always be inferred from the largest column and row indices, due to unmapped grid boxes.
- Returns:
- xr.DataArray
A DataArray backed by a sparse.COO array, with dims (‘out_dim’, ‘in_dim’) and size (n_out, n_in).
- xesmf.smm.check_shapes(indata, weights, shape_in, shape_out)
Compare the shapes of the input array, the weights and the regridder and raises potential errors.
- Parameters:
- indataarray
Input array with the two spatial dimensions at the end, which should fit shape_in.
- weightsarray
Weights 2D array of shape (out_dim, in_dim). First element should be the product of shape_out. Second element should be the product of shape_in.
- shape_in2-tuple of int
Shape of the input of the Regridder.
- shape_out2-tuple of int
Shape of the output of the Regridder.
- Raises:
- ValueError
If any of the conditions is not respected.
- xesmf.smm.apply_weights(weights, indata, shape_in, shape_out)
Apply regridding weights to data.
- Parameters:
- weightssparse COO matrix
Regridding weights.
- indatanumpy array of shape
(..., n_lat, n_lon)
or(..., n_y, n_x)
. Should be C-ordered. Will be then tranposed to F-ordered.
- shape_in, shape_outtuple of two integers
Input/output data shape. For rectilinear grid, it is just
(n_lat, n_lon)
.
- Returns:
- outdatanumpy array of shape
(..., shape_out[0], shape_out[1])
. Extra dimensions are the same as indata. If input data is C-ordered, output will also be C-ordered.
- outdatanumpy array of shape
- xesmf.smm.add_nans_to_weights(weights)
Add NaN in empty rows of the regridding weights sparse matrix.
By default, empty rows in the weights sparse matrix are interpreted as zeroes. This can become problematic when the field being interpreted has legitimate null values. This function inserts NaN values in each row to make sure empty weights are propagated as NaNs instead of zeros.
- Parameters:
- weightsDataArray backed by a sparse.COO array
Sparse weights matrix.
- Returns:
- DataArray backed by a sparse.COO array
Sparse weights matrix.