transforms.aggregate.spatial ============================ .. py:module:: transforms.aggregate.spatial Attributes ---------- .. autoapisummary:: transforms.aggregate.spatial.logger Functions --------- .. autoapisummary:: transforms.aggregate.spatial.rasterize transforms.aggregate.spatial.mask_contains_points transforms.aggregate.spatial.shapes_to_masks transforms.aggregate.spatial.shapes_to_mask transforms.aggregate.spatial.get_mask_dim_index transforms.aggregate.spatial.masks transforms.aggregate.spatial.mask transforms.aggregate.spatial.reduce Module Contents --------------- .. py:data:: logger .. py:function:: rasterize(shape_list, coords, lat_key = 'latitude', lon_key = 'longitude', **kwargs) Rasterize a list of geometries onto the given xarray coordinates. This only works for regular and contiguous latitude and longitude grids. :param shape_list: List of geometries :type shape_list: :class:`affine.Affine` :param coords: Coordinates of dataarray to be masked :type coords: :class:`xarray.coords` :param lat_key: name of the latitude variables in the coordinates object :param lon_key: name of the longitude variables in the coordinates object :param dtype: datatype of the returned mask, default is `int` :param kwargs: Any other kwargs accepted by rasterio.features.rasterize :returns: A mask where points not inside the shape_list are set to `fill` value :rtype: :class:`xr.DataArray` .. py:function:: mask_contains_points(shape_list, coords, lat_key = 'latitude', lon_key = 'longitude', **_kwargs) Return a mask array for the spatial points of data that lie within shapes in shape_list. Function uses matplotlib.Path so can accept a list of points, this is much faster than shapely. It was initially included for use with irregular data but has been constructed to also accept regular data and return in the same format as the rasterize function. :param shape_list: List of geometries :param coords: Coordinates of dataarray to be masked :type coords: :class:`xarray.coords` :param lat_key: name of the latitude variables in the coordinates object :param lon_key: name of the longitude variables in the coordinates object :param dtype: datatype of the returned mask, default is `int` :returns: A mask where points not inside the shape_list are set to `fill` value :rtype: :class:`xr.DataArray` .. py:function:: shapes_to_masks(shapes, target, regular=True, **kwargs) Create a list of masked dataarrays, if possible use the shape_mask_iterator. :param shapes: A geodataframe or list of geodataframes containing the polygons for masks :param target: A dataarray to to create a mask for, only the geospatial coordinates are used :param regular: If True, data is on a regular grid so use rasterize method, if False use mask_contains_points :param all_touched: If True, all pixels touched by geometries will be considered in, if False, only pixels whose center is within. Default is False. Only valid for regular data. :param kwargs: kwargs accepted by the masking methods, rasterize or mask_contains_points :returns: A list of masks where points inside each geometry are 1, and those outside are np.nan :rtype: :class:`list[xr.DataArray]` .. py:function:: shapes_to_mask(shapes, target, regular=True, **kwargs) Create a single masked dataarray based on all features in shapes. If possible use the shape_mask_iterator. :param shapes: A geodataframe or list of geodataframes containing the polygons for masks :param target: A dataarray to to create a mask for, only the geospatial coordinates are used :param regular: If True, data is on a regular grid so use rasterize method, if False use mask_contains_points :param all_touched: If True, all pixels touched by geometries will be considered in, if False, only pixels whose center is within. Default is False. Only valid for regular data. :param kwargs: kwargs accepted by the masking methods, rasterize or mask_contains_points :returns: A mask where points inside any geometry are 1, and those outside are np.nan :rtype: :class:`xr.DataArray` .. py:function:: get_mask_dim_index(mask_dim, geodataframe, default_index_name = 'index') .. py:function:: masks(dataarray, geodataframe, *_args, **_kwargs) .. py:function:: mask(dataarray, geodataframe, mask_dim = None, lat_key = None, lon_key = None, chunk = True, union_geometries = False, **mask_kwargs) Apply multiple shape masks to some gridded data. Each feature in shape is treated as an individual mask to apply to data. The data provided is returned with an additional dimension equal in length to the number of features in the shape object, this can result in very large files which will slow down your script. It may be better to loop over individual features, or directly apply the mask with the shapes.reduce. :param dataarray: Xarray data object (must have geospatial coordinates). :param geodataframe: Geopandas Dataframe containing the polygons for aggregations :param mask_dim: dimension that will be created to accomodate the masked arrays, default is the index of the geodataframe :param all_touched: If True, all pixels touched by geometries will be considered in, if False, only pixels whose center is within. Default is False. Only valid for regular data. :param lat_key: key for latitude variable, default behaviour is to detect variable keys. :param lon_key: key for longitude variable, default behaviour is to detect variable keys. :param chunk: Boolean to indicate whether to use chunking, default = `True`. This is advised as spatial.masks can create large results. If you are working with small arrays, or you have implemented you own chunking rules you may wish to disable it. :type chunk: :class:`bool` :param union_geometries: Boolean to indicate whether to union all geometries before masking. Default is `False`, which will apply each geometry in the geodataframe as a separate mask. :type union_geometries: :class:`bool` :param mask_kwargs: Any kwargs to pass into the mask method :returns: A masked data array with dimensions [feautre_id] + [data.dims]. Each slice of layer corresponds to a feature in layer. :rtype: :class:`xr.Dataset | xr.DataArray` .. py:function:: reduce(dataarray, geodataframe = None, mask_arrays = None, **kwargs) Apply a shape object to an xarray.DataArray object using the specified 'how' method. Geospatial coordinates are reduced to a dimension representing the list of features in the shape object. :param dataarray: Xarray data object (must have geospatial coordinates). :param geodataframe: Geopandas Dataframe containing the polygons for aggregations :param mask_arrays: precomputed mask array[s], if provided this will be used instead of creating a new mask. They must be on the same spatial grid as the dataarray. :param how: method used to apply mask. Default='mean', which calls np.nanmean :param weights: Provide weights for aggregation, also accepts recognised keys for weights, e.g. 'latitude' :param lat_key/lon_key: key for latitude/longitude variable, default behaviour is to detect variable keys. :param extra_reduce_dims: any additional dimensions to aggregate over when reducing over spatial dimensions :param mask_dim: dimension that will be created after the reduction of the spatial dimensions, default is the index of the dataframe :param all_touched: If True, all pixels touched by geometries will be considered in, if False, only pixels whose center is within. Default is False. Only valid for regular data. :param mask_kwargs: Any kwargs to pass into the mask method :param mask_arrays: precomputed mask array[s], if provided this will be used instead of creating a new mask. They must be on the same spatial grid as the dataarray. :param return_as: what format to return the data object, `pandas` or `xarray`. Work In Progress :param compact: If True, return a compact pandas.DataFrame with the reduced data as a new column. If False, return a fully expanded pandas.DataFrame. Only valid if return_as is `pandas` :param how_label: label to append to variable name in returned object, default is not to append :param kwargs: kwargs recognised by the how function :returns: A data array with dimensions `features` + `data.dims not in 'lat','lon'`. Each slice of layer corresponds to a feature in layer. :rtype: :class:`xr.Dataset | xr.DataArray`