gcpy.grid
Module Contents
Classes
Generator for cubed-sphere grid geometries. |
Functions
Returns a mask array for picking out the tropospheric grid boxes. |
|
|
Returns resolution of dataset passed to compare_single_level or compare_zonal_means |
|
Create a mask with NaN values removed from an input array |
|
Get min and max lat and lon from an input GEOS-Chem xarray dataset or grid dict |
|
Determine vertical grid of input dataset |
|
Get indices where edge pressure values are within a given pressure range |
|
Add outer indices to edge pressure index list |
|
Get index of pressure level that contains the requested pressure value. |
|
Convert lev dimension to pressure in a GEOS-Chem dataset |
|
Creates a lat/lon grid description. |
|
Creates a cubed-sphere grid description. |
|
Creates a stretched-grid grid description. |
|
Compute longitude edge vector for a rectilinear grid. |
|
Compute latitude edge vector for a rectilinear grid. |
|
Compute grid cell areas (in m2) for a rectilinear grid. |
|
Compute grid cell longitude widths from an edge vector. |
|
Return cubedsphere coordinates with GMAO face orientation |
|
Convert latitude/longitude coordinates along the unit sphere to cartesian |
|
Convert a cartesian coordinate to latitude/longitude coordinates. |
|
Convert spherical coordinates in the form (theta, phi[, r]) to |
|
Convert cartesian coordinates to spherical in the form |
|
Rotate a spherical coordinate in the form (theta, phi[, r]) |
Attributes
- gcpy.grid.get_troposphere_mask(ds)
Returns a mask array for picking out the tropospheric grid boxes.
- Args:
- ds: xarray Dataset
Dataset containing certain met field variables (i.e. Met_TropLev, Met_BXHEIGHT).
- Returns:
- tropmask: numpy ndarray
Tropospheric mask. False denotes grid boxes that are in the troposphere and True in the stratosphere (as per Python masking logic).
- gcpy.grid.get_input_res(data)
Returns resolution of dataset passed to compare_single_level or compare_zonal_means
- Args:
- data: xarray Dataset
Input GEOS-Chem dataset
- Returns:
- res: str or int
Lat/lon res of the form ‘latresxlonres’ or cubed-sphere resolution
- gridtype: str
‘ll’ for lat/lon or ‘cs’ for cubed-sphere
- gcpy.grid.call_make_grid(res, gridtype, in_extent=[-180, 180, -90, 90], out_extent=[-180, 180, -90, 90], sg_params=[1, 170, -90])
Create a mask with NaN values removed from an input array
- Args:
- res: str or int
Resolution of grid (format ‘latxlon’ or csres)
- gridtype: str
‘ll’ for lat/lon or ‘cs’ for cubed-sphere
- Keyword Args (optional):
- in_extent: list[float, float, float, float]
Describes minimum and maximum latitude and longitude of input data in the format [minlon, maxlon, minlat, maxlat] Default value: [-180, 180, -90, 90]
- out_extent: list[float, float, float, float]
Desired minimum and maximum latitude and longitude of output grid in the format [minlon, maxlon, minlat, maxlat] Default value: [-180, 180, -90, 90]
- sg_params: list[float, float, float] (stretch_factor, target_longitude, target_latitude)
Desired stretched-grid parameters in the format [stretch_factor, target_longitude, target_latitude]. Will trigger stretched-grid creation if not default values. Default value: [1, 170, -90] (no stretching)
- Returns:
- [grid, grid_list]: list(dict, list(dict))
Returns the created grid. grid_list is a list of grids if gridtype is ‘cs’, else it is None
- gcpy.grid.get_grid_extents(data, edges=True)
Get min and max lat and lon from an input GEOS-Chem xarray dataset or grid dict
- Args:
- data: xarray Dataset or dict
A GEOS-Chem dataset or a grid dict
- edges (optional): bool
Whether grid extents should use cell edges instead of centers Default value: True
- Returns:
- minlon: float
Minimum longitude of data grid
- maxlon: float
Maximum longitude of data grid
- minlat: float
Minimum latitude of data grid
- maxlat: float
Maximum latitude of data grid
- gcpy.grid.get_vert_grid(dataset, AP=[], BP=[])
Determine vertical grid of input dataset
- Args:
- dataset: xarray Dataset
A GEOS-Chem output dataset
- Keyword Args (optional):
- AP: list-like type
Hybrid grid parameter A in hPa Default value: []
- BP: list-like type
Hybrid grid parameter B (unitless) Default value: []
- Returns:
- p_edge: numpy array
Edge pressure values for vertical grid
- p_mid: numpy array
Midpoint pressure values for vertical grid
- nlev: int
Number of levels in vertical grid
- gcpy.grid.get_pressure_indices(pedge, pres_range)
Get indices where edge pressure values are within a given pressure range
- Args:
- pedge: numpy array
A GEOS-Chem output dataset
- pres_range: list(float, float)
Contains minimum and maximum pressure
- Returns:
- numpy array
Indices where edge pressure values are within a given pressure range
- gcpy.grid.pad_pressure_edges(pedge_ind, max_ind, pmid_len)
Add outer indices to edge pressure index list
- Args:
- pedge_ind: list
List of edge pressure indices
- max_ind: int
Maximum index
- pmid_len: int
Length of pmid which should not be exceeded by indices
- Returns:
- pedge_ind: list
List of edge pressure indices, possibly with new minimum and maximum indices
- gcpy.grid.get_ind_of_pres(dataset, pres)
Get index of pressure level that contains the requested pressure value.
- Args:
- dataset: xarray Dataset
GEOS-Chem dataset
- pres: int or float
Desired pressure value
- Returns:
- index: int
Index of level in dataset that corresponds to requested pressure
- gcpy.grid.convert_lev_to_pres(dataset, pmid, pedge, lev_type='pmid')
Convert lev dimension to pressure in a GEOS-Chem dataset
- Args:
- dataset: xarray Dataset
GEOS-Chem dataset
- pmid: np.array
Midpoint pressure values
- pedge: np.array
Edge pressure values
- lev_type (optional): str
Denote whether lev is ‘pedge’ or ‘pmid’ if grid is not 72/73 or 47/48 levels Default value: ‘pmid’
- Returns:
- dataset: xarray Dataset
Input dataset with “lev” dimension values replaced with pressure values
- gcpy.grid._GEOS_72L_AP
- gcpy.grid._GEOS_72L_BP
- gcpy.grid.GEOS_72L_grid
- gcpy.grid._GEOS_47L_AP
- gcpy.grid._GEOS_47L_BP
- gcpy.grid._xmat_i
- gcpy.grid._xmat_j
- gcpy.grid._xmat_s
- gcpy.grid._x_lev
- gcpy.grid._skip_size_vec = [2, 4]
- gcpy.grid._number_lumped = [4, 7]
- gcpy.grid._i_lev = 36
- gcpy.grid._i_lev_72 = 36
- gcpy.grid._skip_size
- gcpy.grid._xmat_72to47
- gcpy.grid.GEOS_47L_grid
- gcpy.grid._CAM_26L_AP
- gcpy.grid._CAM_26L_BP
- gcpy.grid.CAM_26L_grid
- gcpy.grid.make_grid_LL(llres, in_extent=[-180, 180, -90, 90], out_extent=[])
Creates a lat/lon grid description.
- Args:
- llres: str
lat/lon resolution in ‘latxlon’ format (e.g. ‘4x5’)
- Keyword Args (optional):
- in_extent: list[float, float, float, float]
Describes minimum and maximum latitude and longitude of initial grid in the format [minlon, maxlon, minlat, maxlat] Default value: [-180, 180, -90, 90]
- out_extent: list[float, float, float, float]
Describes minimum and maximum latitude and longitude of target grid in the format [minlon, maxlon, minlat, maxlat]. Needed when intending to use grid to trim extent of input data Default value: [] (assumes value of in_extent)
- Returns:
- llgrid: dict
- dict grid description of format {‘lat’lat midpoints,
‘lon’ : lon midpoints, ‘lat_b’ : lat edges, ‘lon_b’ : lon edges}
- gcpy.grid.make_grid_CS(csres)
Creates a cubed-sphere grid description.
- Args:
- csres: int
cubed-sphere resolution of target grid
- Returns:
- [csgrid, csgrid_list]: list[dict, list[dict]]
- csgrid is a dict of format {‘lat’lat midpoints,
‘lon’ : lon midpoints, ‘lat_b’ : lat edges, ‘lon_b’ : lon edges}
where each value has an extra face dimension of length 6. csgrid_list is a list of dicts separated by face index
- gcpy.grid.make_grid_SG(csres, stretch_factor, target_lon, target_lat)
Creates a stretched-grid grid description.
- Args:
- csres: int
cubed-sphere resolution of target grid
- stretch_factor: float
stretch factor of target grid
- target_lon: float
target stretching longitude of target grid
- target_lon: float
target stretching latitude of target grid
- Returns:
- [csgrid, csgrid_list]: list[dict, list[dict]]
- csgrid is a dict of format {‘lat’lat midpoints,
‘lon’ : lon midpoints, ‘lat_b’ : lat edges, ‘lon_b’ : lon edges}
where each value has an extra face dimension of length 6. csgrid_list is a list of dicts separated by face index
- gcpy.grid.calc_rectilinear_lon_edge(lon_stride, center_at_180)
Compute longitude edge vector for a rectilinear grid. Parameters ———- lon_stride: float
Stride length in degrees. For example, for a standard GEOS-Chem Classic 4x5 grid, lon_stride would be 5.
- center_at_180: bool
Whether or not the grid should have a cell center at 180 degrees (i.e. on the date line). If true, the first grid cell is centered on the date line; if false, the first grid edge is on the date line.
Returns
Longitudes of cell edges in degrees East. Notes —– All values are forced to be between [-180,180]. For a grid with N cells in each band, N+1 edges will be returned, with the first and last value being duplicates. Examples ——– >>> from gcpy.grid.horiz import calc_rectilinear_lon_edge >>> calc_rectilinear_lon_edge(5.0,true) np.array([177.5,-177.5,-172.5,…,177.5]) See Also ——– [NONE]
- gcpy.grid.calc_rectilinear_lat_edge(lat_stride, half_polar_grid)
Compute latitude edge vector for a rectilinear grid. Parameters ———- lat_stride: float
Stride length in degrees. For example, for a standard GEOS-Chem Classic 4x5 grid, lat_stride would be 4.
- half_polar_grid: bool
Whether or not the grid should be “half-polar” (i.e. bands at poles are half the size). In either case the grid will start and end at -/+ 90, but when half_polar_grid is True, the first and last bands will have a width of 1/2 the normal lat_stride.
Returns
Latitudes of cell edges in degrees North. Notes —– All values are forced to be between [-90,90]. For a grid with N cells in each band, N+1 edges will be returned, with the first and last value being duplicates. Examples ——– >>> from gcpy.grid.horiz import calc_rectilinear_lat_edge >>> calc_rectilinear_lat_edge(4.0,true) np.array([-90,-88,-84,-80,…,84,88,90]) See Also ——– [NONE]
- gcpy.grid.calc_rectilinear_grid_area(lon_edge, lat_edge)
Compute grid cell areas (in m2) for a rectilinear grid. Parameters ———- #TODO Returns ——- #TODO Notes —– #TODO Examples ——– #TODO See Also ——– [NONE]
- gcpy.grid.calc_delta_lon(lon_edge)
Compute grid cell longitude widths from an edge vector. Parameters ———- lon_edge: float
Vector of longitude edges, in degrees East.
Returns
Width of each cell, degrees East Notes —– Accounts for looping over the domain. Examples ——– #TODO
- gcpy.grid.csgrid_GMAO(res, offset=-10)
Return cubedsphere coordinates with GMAO face orientation Parameters ———- res: cubed-sphere Resolution This function was originally written by Jiawei Zhuange and included in package cubedsphere: https://github.com/JiaweiZhuang/cubedsphere
- gcpy.grid._INV_SQRT_3
- gcpy.grid._ASIN_INV_SQRT_3
- class gcpy.grid.CSGrid(c, offset=None)
Bases:
objectGenerator for cubed-sphere grid geometries. CSGrid computes the latitutde and longitudes of cell centers and edges on a cubed-sphere grid, providing a way to retrieve these geometries on-the-fly if your model output data does not include them. Attributes ———- {lon,lat}_center: np.ndarray
lat/lon coordinates for each cell center along the cubed-sphere mesh
- {lon,lat}_edge: np.ndarray
lat/lon coordinates for the midpoint of the edges separating each element on the cubed-sphere mesh.
- xyz_{center,edge}: np.ndarray
As above, except coordinates are projected into a 3D cartesian space with common origin to the original lat/lon coordinate system, assuming a unit sphere.
This class was originally written by Jiawei Zhuange and included in package cubedsphere: https://github.com/JiaweiZhuang/cubedsphere
- _initialize()
- gcpy.grid.latlon_to_cartesian(lon, lat)
Convert latitude/longitude coordinates along the unit sphere to cartesian coordinates defined by a vector pointing from the sphere’s center to its surface. This function was originally written by Jiawei Zhuange and included in package cubedsphere: https://github.com/JiaweiZhuang/cubedsphere
- gcpy.grid.vec_latlon_to_cartesian
- gcpy.grid.cartesian_to_latlon(x, y, z, ret_xyz=False)
Convert a cartesian coordinate to latitude/longitude coordinates. Optionally return the original cartesian coordinate as a tuple. This function was originally written by Jiawei Zhuange and included in package cubedsphere: https://github.com/JiaweiZhuang/cubedsphere
- gcpy.grid.vec_cartesian_to_latlon
- gcpy.grid.spherical_to_cartesian(theta, phi, r=1)
Convert spherical coordinates in the form (theta, phi[, r]) to cartesian, with the origin at the center of the original spherical coordinate system. This function was originally written by Jiawei Zhuange and included in package cubedsphere: https://github.com/JiaweiZhuang/cubedsphere
- gcpy.grid.vec_spherical_to_cartesian
- gcpy.grid.cartesian_to_spherical(x, y, z)
Convert cartesian coordinates to spherical in the form (theta, phi[, r]) with the origin remaining at the center of the original spherical coordinate system. This function was originally written by Jiawei Zhuange and included in package cubedsphere: https://github.com/JiaweiZhuang/cubedsphere
- gcpy.grid.vec_cartesian_to_spherical
- gcpy.grid.rotate_sphere_3D(theta, phi, r, rot_ang, rot_axis='x')
Rotate a spherical coordinate in the form (theta, phi[, r]) about the indicating axis, ‘rot_axis’. This method accomplishes the rotation by projecting to a cartesian coordinate system and performing a solid body rotation around the requested axis. This function was originally written by Jiawei Zhuange and included in package cubedsphere: https://github.com/JiaweiZhuang/cubedsphere