Mapping (pyart.map)

Py-ART has a robust function for mapping radar data from the collected radar coordinates to Cartesian coordinates.

grid_from_radars(radars, grid_shape, grid_limits) Map one or more radars to a Cartesian grid returning a Grid object.
map_to_grid(radars, grid_shape, grid_limits) Map one or more radars to a Cartesian grid.
map_gates_to_grid(radars, grid_shape, …[, …]) Map gates from one or more radars to a Cartesian grid.
example_roi_func_constant(zg, yg, xg) Example RoI function which returns a constant radius.
example_roi_func_dist(zg, yg, xg) Example RoI function which returns a radius which grows with distance.
example_roi_func_dist_beam(zg, yg, xg) Example RoI function which returns a radius which grows with distance and whose parameters are based on virtual beam size.
pyart.map.example_roi_func_constant(zg, yg, xg)[source]

Example RoI function which returns a constant radius.

Parameters:
zg, yg, xg : float

Distance from the grid center in meters for the x, y and z axes.

Returns:
roi : float

Radius of influence in meters

pyart.map.example_roi_func_dist(zg, yg, xg)[source]

Example RoI function which returns a radius which grows with distance.

Parameters:
zg, yg, xg : float

Distance from the grid center in meters for the x, y and z axes.

Returns:
roi : float
pyart.map.example_roi_func_dist_beam(zg, yg, xg)[source]

Example RoI function which returns a radius which grows with distance and whose parameters are based on virtual beam size.

Parameters:
zg, yg, xg : float

Distance from the grid center in meters for the x, y and z axes.

Returns:
roi : float
pyart.map.get_earth_radius(latitude)[source]

Computes the earth radius for a given latitude

Parameters:
latitude: latitude in degrees (WGS84)
Returns:
earth_radius : the radius of the earth at the given latitude
pyart.map.grid_from_radars(radars, grid_shape, grid_limits, gridding_algo='map_gates_to_grid', **kwargs)[source]

Map one or more radars to a Cartesian grid returning a Grid object.

Additional arguments are passed to map_to_grid() or map_gates_to_grid().

Parameters:
radars : Radar or tuple of Radar objects.

Radar objects which will be mapped to the Cartesian grid.

grid_shape : 3-tuple of floats

Number of points in the grid (z, y, x).

grid_limits : 3-tuple of 2-tuples

Minimum and maximum grid location (inclusive) in meters for the z, y, x coordinates.

gridding_algo : ‘map_to_grid’ or ‘map_gates_to_grid’

Algorithm to use for gridding. ‘map_to_grid’ finds all gates within a radius of influence for each grid point, ‘map_gates_to_grid’ maps each radar gate onto the grid using a radius of influence and is typically significantly faster.

Returns:
grid : Grid

A pyart.io.Grid object containing the gridded radar data.

See also

map_to_grid
Map to grid and return a dictionary of radar fields.
map_gates_to_grid
Map each gate onto a grid returning a dictionary of radar fields.

References

Barnes S., 1964: A Technique for Maximizing Details in Numerical Weather Map Analysis. Journal of Applied Meteorology and Climatology, 3(4), 396-409.

Cressman G., 1959: An operational objective analysis system. Monthly Weather Review, 87(10), 367-374.

Pauley, P. M. and X. Wu, 1990: The theoretical, discrete, and actual response of the Barnes objective analysis scheme for one- and two-dimensional fields. Monthly Weather Review, 118, 1145-1164

pyart.map.map_gates_to_grid(radars, grid_shape, grid_limits, grid_origin=None, grid_origin_alt=None, grid_projection=None, fields=None, gatefilters=False, map_roi=True, weighting_function='Barnes', toa=17000.0, roi_func='dist_beam', constant_roi=None, z_factor=0.05, xy_factor=0.02, min_radius=500.0, h_factor=1.0, nb=1.5, bsp=1.0, **kwargs)[source]

Map gates from one or more radars to a Cartesian grid.

Generate a Cartesian grid of points for the requested fields from the collected points from one or more radars. For each radar gate that is not filtered a radius of influence is calculated. The weighted field values for that gate are added to all grid points within that radius. This routine scaled linearly with the number of radar gates and the effective grid size.

Parameters not defined below are identical to those in map_to_grid().

Parameters:
roi_func : str or RoIFunction

Radius of influence function. A function which takes an z, y, x grid location, in meters, and returns a radius (in meters) within which all collected points will be included in the weighting for that grid points. Examples can be found in the Typically following strings can use to specify a built in radius of influence function:

  • constant: constant radius of influence.
  • dist: radius grows with the distance from each radar.
  • dist_beam: radius grows with the distance from each radar and parameter are based of virtual beam sizes.

A custom RoIFunction can be defined using the RoIFunction class and defining a get_roi method which returns the radius. For efficient mapping this class should be implemented in Cython.

Returns:
grids : dict

Dictionary of mapped fields. The keys of the dictionary are given by parameter fields. Each elements is a grid_size float64 array containing the interpolated grid for that field.

See also

grid_from_radars
Map to a grid and return a Grid object
map_to_grid
Create grid by finding the radius of influence around each grid point.
pyart.map.map_to_grid(radars, grid_shape, grid_limits, grid_origin=None, grid_origin_alt=None, grid_projection=None, fields=None, gatefilters=False, map_roi=True, weighting_function='Barnes', toa=17000.0, copy_field_data=True, algorithm='kd_tree', leafsize=10.0, roi_func='dist_beam', constant_roi=None, z_factor=0.05, xy_factor=0.02, min_radius=500.0, h_factor=1.0, nb=1.5, bsp=1.0, **kwargs)[source]

Map one or more radars to a Cartesian grid.

Generate a Cartesian grid of points for the requested fields from the collected points from one or more radars. The field value for a grid point is found by interpolating from the collected points within a given radius of influence and weighting these nearby points according to their distance from the grid points. Collected points are filtered according to a number of criteria so that undesired points are not included in the interpolation.

Parameters:
radars : Radar or tuple of Radar objects.

Radar objects which will be mapped to the Cartesian grid.

grid_shape : 3-tuple of floats

Number of points in the grid (z, y, x).

grid_limits : 3-tuple of 2-tuples

Minimum and maximum grid location (inclusive) in meters for the z, y, x coordinates.

grid_origin : (float, float) or None

Latitude and longitude of grid origin. None sets the origin to the location of the first radar.

grid_origin_alt: float or None

Altitude of grid origin, in meters. None sets the origin to the location of the first radar.

grid_projection : dic or str

Projection parameters defining the map projection used to transform the locations of the radar gates in geographic coordinate to Cartesian coodinates. None will use the default dictionary which uses a native azimutal equidistance projection. See pyart.core.Grid() for additional details on this parameter. The geographic coordinates of the radar gates are calculated using the projection defined for each radar. No transformation is used if a grid_origin and grid_origin_alt are None and a single radar is specified.

fields : list or None

List of fields within the radar objects which will be mapped to the cartesian grid. None, the default, will map the fields which are present in all the radar objects.

gatefilters : GateFilter, tuple of GateFilter objects, optional

Specify what gates from each radar will be included in the interpolation onto the grid. Only gates specified in each gatefilters will be included in the mapping to the grid. A single GateFilter can be used if a single Radar is being mapped. A value of False for a specific element or the entire parameter will apply no filtering of gates for a specific radar or all radars (the default). Similarily a value of None will create a GateFilter from the radar moments using any additional arguments by passing them to moment_based_gate_filter().

roi_func : str or function

Radius of influence function. A functions which takes an z, y, x grid location, in meters, and returns a radius (in meters) within which all collected points will be included in the weighting for that grid points. Examples can be found in the example_roi_func_constant(), example_roi_func_dist(), and example_roi_func_dist_beam(). Alternatively the following strings can use to specify a built in radius of influence function:

  • constant: constant radius of influence.
  • dist: radius grows with the distance from each radar.
  • dist_beam: radius grows with the distance from each radar and parameter are based of virtual beam sizes.

The parameters which control these functions are listed in the Other Parameters section below.

map_roi : bool

True to include a radius of influence field in the returned dictionary under the ‘ROI’ key. This is the value of roi_func at all grid points.

weighting_function : ‘Barnes’ or ‘Barnes2’ or ‘Cressman’ or ‘Nearest’

Functions used to weight nearby collected points when interpolating a grid point.

toa : float

Top of atmosphere in meters. Collected points above this height are not included in the interpolation.

Returns:
grids : dict

Dictionary of mapped fields. The keys of the dictionary are given by parameter fields. Each elements is a grid_size float64 array containing the interpolated grid for that field.

Other Parameters:
 
constant_roi : float

Radius of influence parameter for the built in ‘constant’ function. This parameter is the constant radius in meter for all grid points. This parameter is used when roi_func is constant or constant_roi is not None. If constant_roi is not None, the constant roi_func is used automatically.

z_factor, xy_factor, min_radius : float

Radius of influence parameters for the built in ‘dist’ function. The parameter correspond to the radius size increase, in meters, per meter increase in the z-dimension from the nearest radar, the same foreach meteter in the xy-distance from the nearest radar, and the minimum radius of influence in meters. These parameters are only used when roi_func is ‘dist’.

h_factor, nb, bsp, min_radius : float

Radius of influence parameters for the built in ‘dist_beam’ function. The parameter correspond to the height scaling, virtual beam width, virtual beam spacing, and minimum radius of influence. These parameters are only used when roi_func is ‘dist_mean’.

copy_field_data : bool

True to copy the data within the radar fields for faster gridding, the dtype for all fields in the grid will be float64. False will not copy the data which preserves the dtype of the fields in the grid, may use less memory but results in significantly slower gridding times. When False gates which are masked in a particular field but are not masked in the refl_field field will still be included in the interpolation. This can be prevented by setting this parameter to True or by gridding each field individually setting the refl_field parameter and the fields parameter to the field in question. It is recommended to set this parameter to True.

algorithm : ‘kd_tree’.

Algorithms to use for finding the nearest neighbors. ‘kd_tree’ is the only valid option.

leafsize : int

Leaf size passed to the neighbor lookup tree. This can affect the speed of the construction and query, as well as the memory required to store the tree. The optimal value depends on the nature of the problem. This value should only effect the speed of the gridding, not the results.

See also

grid_from_radars
Map to grid and return a Grid object.
pyart.map.polar_to_cartesian(radar_sweep, field_name, cart_res=75, max_range=None, mapping=None)[source]

Interpolates a PPI or RHI scan in polar coordinates to a regular cartesian grid of South-North and West-East coordinates (for PPI) or distance at ground and altitude coordinates (for RHI)

Parameters:
radar : Radar

Radar instance as generated py pyart

sweep : int

Sweep number to project to cartesian coordinates.

field_name : str

Name of the radar field to be interpolated

cart_res : int, optional

Resolution (in m.) of the cartesian grid to which polar data is interpolated

max_range : int, optional

Maximal allowed range (in m.) from radar for gates to be interpolated

mapping : dict, optional

Dictionnary of mapping indexes (from polar to cartesian), gets returned by the function (see below). Can be used as input when interpolating sequentially several variables for the same scan, to save significant time

Returns:
coords : tuple of 2 arrays

2D coordinates of the cartesian grid

cart_data : 2D array

Interpolated radar measurements (on the cartesian grid)

mapping,: dict

Dictionnary of mapping indexes (from polar to cartesian),which contains the indexes mapping the polar grid to the cartesian grid as well as some metadata.