geocat.comp.triple_to_grid
geocat.comp.triple_to_grid#
- geocat.comp.triple_to_grid(data, x_in, y_in, x_out, y_out, method=1, domain=1.0, distmx=None, missing_value=None, meta=False)#
Places unstructured (randomly-spaced) data onto the nearest locations of a rectilinear grid.
- Parameters
data (
xarray.DataArray
:,numpy.ndarray
:) – A multi-dimensional array, whose rightmost dimension is the same length as x_in and y_in, containing the values associated with the “x” and “y” coordinates. Missing values may be present but will be ignored.x_in (
xarray.DataArray
:,numpy.ndarray
:) – A one-dimensional array that specifies the x-coordinate associated with the input (data).y_in (
xarray.DataArray
:,numpy.ndarray
:) – A one-dimensional array that specifies the y-coordinate associated with the input (data).x_out (
xarray.DataArray
:,numpy.ndarray
:) – A one-dimensional array of length M containing the x-coordinates associated with the returned two-dimensional grid. The coordinate values must be monotonically increasing.y_out (
xarray.DataArray
:or :class:`numpy.ndarray
:`) – A one-dimensional array of length N containing the y-coordinates associated with the returned two-dimensional grid. The coordinate values must be monotonically increasing.Optional Parameters
——————-
method :obj:`int` – An integer value that can be 0 or 1. The default value is 1. A value of 1 means to use the great circle distance formula for distance calculations. Warning: method = 0, together with domain = 1.0, could result in many of the target grid points to be set to the missing value if the number of grid points is large (ie: a high resolution grid) and the number of observations relatively small.
domain :obj:`float` – A float value that should be set to a value >= 0. The default value is 1.0. If present, the larger this factor, the wider the spatial domain allowed to influence grid boundary points. Typically, domain is 1.0 or 2.0. If domain <= 0.0, then values located outside the grid domain specified by x_out and y_out arguments will not be used.
distmx :obj:`float` – Setting distmx allows the user to specify a search radius (km) beyond which observations are not considered for nearest neighbor. Only applicable when method = 1. The default distmx`=1e20 (km) means that every grid point will have a nearest neighbor. It is suggested that users specify a reasonable value for `distmx.
missing_value (
numpy.number
:) – A numpy scalar value that represent a missing value in data. The default value is np.nan. If specified explicitly, this argument allows the user to use a missing value scheme other than NaN or masked arrays.meta (
bool
:) – If set to True and the input array is an Xarray, the metadata from the input array will be copied to the output array; default is False. Warning: This option is not yet supported for this function.
- Returns
grid (
xarray.DataArray
,numpy.ndarray
:) – The returned array will be K x N x M, where K represents the leftmost dimensions of data, N represent the size of y_out, and M represent the size of x_out coordinate vectors.Description
-----------
– This function puts unstructured data (randomly-spaced) onto the nearest locations of a rectilinear grid. A default value of domain option is now set to 1.0 instead of 0.0.This function does not perform interpolation; rather, each individual data point is assigned to the nearest grid point. It is possible that upon return, grid will contain grid points set to missing value if no x_in(n), y_in(n) are nearby.
Examples
Example 1: Using triple_to_grid with
xarray.DataArray
inputimport numpy as np import xarray as xr import geocat.comp # Open a netCDF data file using xarray default engine and load the data stream ds = xr.open_dataset("./ruc.nc") # [INPUT] Grid & data info on the source curvilinear data = ds.DIST_236_CBL[:] x_in = ds.gridlat_236[:] y_in = ds.gridlon_236[:] x_out = ds.gridlat_236[:] y_out = ds.gridlon_236[:] # [OUTPUT] Grid on destination points grid (or read the 1D lat and lon from # an other .nc file. newlat1D_points=np.linspace(lat2D_curv.min(), lat2D_curv.max(), 100) newlon1D_points=np.linspace(lon2D_curv.min(), lon2D_curv.max(), 100) output = geocat.comp.triple_to_grid(data, x_out, y_out, x_in, y_in)