combineSimilarRasters

Modules:

• GEOM –
• RASTER –
• SRS –
• UTIL –

  The Util sub-module contains a number of generally helpful utility functions, classes, and constants.

• VECTOR –

Classes:

• Extent –

  Geographic extent.

• RegionMask –

  The RegionMask object represents a given region and exposes methods allowing

Functions:

• checkSimilarRasters –

  Parameters

• combineSimilarRasters –

  Combines several similar raster files into one single raster file.

Extent

Extent(*args, srs: srs_input = 'latlon')

Geographic extent.

The Extent object represents geographic extents of an area and exposes useful methods which depend on those extents. This includes: - Easily representing the boundaries as (xMin, xMax, yMin, yMax) or (xMin, yMin, xMax, yMax) - Casting to another projection system - Padding and shifting the boundaries - "Fitting" the boundaries onto a given resolution - Clipping a given raster file

Initialization:

• Extent(xMin, yMin, xMax, yMax [, srs])
• Extent.from_xyXY( (xMin, yMin, xMax, yMax) [, srs])
• Extent.from_xXyY( (xMin, xMax, yMin, yMax) [, srs])
• Extent.fromGeom( geom [, srs] )
• Extent.fromVector( vector-file-path )
• Extent.fromRaster( raster-file-path )
• Extent.load( args )

Create extent from explicitly defined boundaries.

Usage:

Extent(xMin, yMin, xMax, yMax [, srs=]) Extent( (xMin, yMin, xMax, yMax) [, srs=])

Where: xMin - The minimal x value in the respective SRS yMin - The minimal y value in the respective SRS xMax - The maximal x value in the respective SRS yMax - The maximal y value in the respective SRS srs - The Spatial Reference system to use

Methods:

• castTo –

  Creates a new Extent by transforming an extent from the original Extent's

• center –

  Get the Extent's center.

• clipRaster –

  Clip a given raster source to the calling Extent.

• computePixelSize –

  Finds the pixel resolution which fits to the Extent for a given pixel count.

• contains –

  Tests if the extent contains another given extent.

• containsLoc –

  Test if the extent contains a location or an iterable of locations.

• contoursFromRaster –

  Convenience wrapper for geokit.raster.contours which automatically

• corners –

  Returns the four corners of the extent as ogr.gGometry points or as (x,y)

• createRaster –

  Convenience function for geokit.raster.createRaster which sets 'bounds'

• drawSmopyMap –

  Draws a basemap using the "smopy" python package.

• exportWKT –

  Export the extent to a Well-Known_Text string.

• extractFeatures –

  Convenience wrapper for geokit.vector.extractFeatures() by setting the

• extractMatrix –

  Convenience wrapper around geokit.raster.extractMatrix(). Extracts the

• filterSources –

  Filter a list of sources by those whose's envelope overlaps the Extent.

• findWithin –

  Finds the indexes of the given extent within the main extent according

• fit –

  Fit the extent to a given pixel resolution.

• fitsResolution –

  Test if calling Extent first around the given unit(s) (at least within

• fromGeom –

  Create extent around a given geometry.

• fromLocationSet –

  Create extent around the contents of a LocationSet object.

• fromRaster –

  Create extent around the contents of a raster source.

• fromTile –

  Generates an Extent corresponding to tiles used for "slippery maps".

• fromTileAt –

  Generates an Extent corresponding to tiles used for "slippery maps"

• fromVector –

  Create extent around the contemts of a vector source.

• fromWKT –

  Create extent from a Well-Known_Text string.

• from_xXyY –

  Create an Extent from explicitly defined boundaries.

• inSourceExtent –

  Tests if the extent box is at least partially contained in the extent-box

• load –

  Attempts to load an Extent from a variety of inputs in the most

• mutateRaster –

  Convenience function for geokit.raster.mutateRaster which automatically

• mutateVector –

  Convenience function for geokit.vector.mutateVector which automatically

• overlaps –

  Tests if the extent overlaps with another given extent.

• pad –

  Pad the extent in all directions.

• rasterMosaic –

  Create a raster source surrounding the Extent from a collection of other rasters.

• rasterize –

  Convenience function for geokit.vector.rasterize() which automatically

• shift –

  Shift the extent in the X and/or Y dimensions.

• subTiles –

  Generates tile Extents at a given zoom level which encompass the invoking Extent.

• tileBox –

  Determine the tile Extent at a given zoom level which surround the invoked Extent.

• tileIndexBox –

  Determine the tile indexes at a given zoom level which surround the invoked Extent.

• tileMosaic –

  Create a raster source surrounding the Extent from a collection of tiles.

• tileSources –

  Get the tiles sources which contribute to the invoking Extent.

• warp –

  Convenience function for geokit.raster.warp() which automatically sets the

Attributes:

• YxyX (tuple[float, float, float, float]) –

  Returns a tuple of the extent boundaries in order:

• box –

  Returns a rectangular ogr.Geometry object representing the extent.

• xXyY (tuple[float, float, float, float]) –

  Returns a tuple of the extent boundaries in order:

• xYXy (tuple[float, float, float, float]) –

  Returns a tuple of the extent boundaries in order:

• xlim (tuple[float, float]) –

  Returns a tuple of the x-axis extent boundaries in order:

• xyXY (tuple[float, float, float, float]) –

  Returns a tuple of the extent boundaries in order:

• ylim (tuple[float, float]) –

  Returns a tuple of the y-axis extent boundaries in order:

• yxYX (tuple[float, float, float, float]) –

  Returns a tuple of the extent boundaries in order:

Source code in geokit/core/extent.py

def __init__(self, *args, srs: srs_input = "latlon"):
    """Create extent from explicitly defined boundaries.

    Usage:
    ------
    Extent(xMin, yMin, xMax, yMax [, srs=<srs>])
    Extent( (xMin, yMin, xMax, yMax) [, srs=<srs>])

    Where:
        xMin - The minimal x value in the respective SRS
        yMin - The minimal y value in the respective SRS
        xMax - The maximal x value in the respective SRS
        yMax - The maximal y value in the respective SRS
        srs - The Spatial Reference system to use
    """
    # Unpack args
    if len(args) == 1:
        xMin, yMin, xMax, yMax = args[0]
    elif len(args) == 4:
        xMin, yMin, xMax, yMax = args
    else:
        raise GeoKitExtentError(
            "Incorrect number of positional arguments given in init (accepts 1 or 4). Is an srs given as 'srs=...'?"
        )

    # Ensure good inputs
    xMin, xMax = min(xMin, xMax), max(xMin, xMax)
    yMin, yMax = min(yMin, yMax), max(yMin, yMax)

    self.xMin = xMin
    self.xMax = xMax
    self.yMin = yMin
    self.yMax = yMax
    self.srs = SRS.loadSRS(srs)

    self._box = GEOM.box(self.xMin, self.yMin, self.xMax, self.yMax, srs=self.srs)

YxyX property

YxyX: tuple[float, float, float, float]

Returns a tuple of the extent boundaries in order: yMax, xMin, yMin, xMax.

box property

box

Returns a rectangular ogr.Geometry object representing the extent.

xXyY property

xXyY: tuple[float, float, float, float]

Returns a tuple of the extent boundaries in order: xMin, xMax, yMin, yMax.

xYXy property

xYXy: tuple[float, float, float, float]

Returns a tuple of the extent boundaries in order: xMin, yMax, xMax, yMin.

xlim property

xlim: tuple[float, float]

Returns a tuple of the x-axis extent boundaries in order: xMin, xMax.

xyXY property

xyXY: tuple[float, float, float, float]

Returns a tuple of the extent boundaries in order: xMin, yMin, xMax, yMax.

ylim property

ylim: tuple[float, float]

Returns a tuple of the y-axis extent boundaries in order: yMin, yMax.

yxYX property

yxYX: tuple[float, float, float, float]

Returns a tuple of the extent boundaries in order: yMin, xMin, yMax, xMax.

castTo

castTo(srs: srs_input, segments: numeric = 100)

Creates a new Extent by transforming an extent from the original Extent's srs to a target SRS.

Note:

The resulting region spanned by the extent will be equal-to or (almost certainly) larger than the original

Parameters:

• srs

  (Anything acceptable to geokit.srs.loadSRS()) –

  The srs to cast the Extent object to

Returns:

• Extent –

Source code in geokit/core/extent.py

def castTo(self, srs: srs_input, segments: numeric = 100):
    """
    Creates a new Extent by transforming an extent from the original Extent's
    srs to a target SRS.

    Note:
    -----
    The resulting region spanned by the extent will be equal-to or (almost
    certainly) larger than the original

    Parameters
    ----------
    srs : Anything acceptable to geokit.srs.loadSRS()
        The srs to cast the Extent object to

    Returns
    -------
    Extent
    """
    srs = SRS.loadSRS(srs)

    if srs.IsSame(self.srs):
        return self

    segment_size = min(self.xMax - self.xMin, self.yMax - self.yMin) / segments
    box = self.box
    box.Segmentize(segment_size)
    box = GEOM.transform(box, toSRS=srs)
    xMin, xMax, yMin, yMax = box.GetEnvelope()
    return Extent(xMin, yMin, xMax, yMax, srs=srs)

center

center(srs=None)

Get the Extent's center.

Source code in geokit/core/extent.py

def center(self, srs=None):
    """Get the Extent's center."""
    x, y = (self.xMax + self.xMin) / 2, (self.yMax + self.yMin) / 2
    if not srs is None:
        srs = SRS.loadSRS(srs)
        if not srs.IsSame(self.srs):
            xy = SRS.xyTransform(x, y, fromSRS=self.srs, toSRS=srs, outputFormat="xy")

            x = xy.x
            y = xy.y
    return x, y

clipRaster

clipRaster(
    source, output: None | str = None, **kwargs
) -> None | Dataset

Clip a given raster source to the calling Extent.

Parameters:

• source

  (Anything acceptable to geokit.raster.loadRaster()) –

  The source to clip

• **kwargs

  –

  All other keyword arguments are passed to gdal.Translate

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/extent.py

def clipRaster(self, source, output: None | str = None, **kwargs) -> None | gdal.Dataset:
    """Clip a given raster source to the calling Extent.

    Parameters
    ----------
    source : Anything acceptable to geokit.raster.loadRaster()
        The source to clip

    **kwargs:
        All other keyword arguments are passed to gdal.Translate

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    from time import time_ns

    opts = gdal.TranslateOptions(
        projWin=[self.xMin, self.yMax, self.xMax, self.yMin],
        projWinSRS=self.srs,
        **kwargs,
    )

    if output is None:
        fname = "/vsimem/clip_{}.tif".format(time_ns())
    else:
        fname = output
    ds = gdal.Translate(fname, source, options=opts)

    if output is None:
        return ds
    else:
        return output

computePixelSize

computePixelSize(*args)

Finds the pixel resolution which fits to the Extent for a given pixel count.

Note:

• If only one integer argument is given, it is assumed to fit to both the X and Y dimensions
• If two integer arguments are given, it is assumed to be in the order X then Y

Returns:

• tuple -> (pixelWidth, pixelHeight) –

Source code in geokit/core/extent.py

def computePixelSize(self, *args):
    """Finds the pixel resolution which fits to the Extent for a given pixel count.

    Note:
    -----
    * If only one integer argument is given, it is assumed to fit to both the X and Y dimensions
    * If two integer arguments are given, it is assumed to be in the order X then Y

    Returns
    -------
        tuple -> (pixelWidth, pixelHeight)
    """
    if len(args) == 1:
        pixels_x = args[0]
        pixels_y = args[0]
    else:
        pixels_x, pixels_y = args

    pixelWidth = (self.xMax - self.xMin) / pixels_x
    pixelHeight = (self.yMax - self.yMin) / pixels_y

    return pixelWidth, pixelHeight

contains

contains(extent, res=None)

Tests if the extent contains another given extent.

Note:

If an optional resolution ('res') is given, the containment value is also dependent on whether or not the given extent fits within the larger extent AND is situated along the given resolution

Parameters:

• extent

  (Extent) –

  The Extent object to test for containment

• res

  (numeric or tuple, default: None ) –

  The X & Y resolution to enforce

Returns:

• bool –

Source code in geokit/core/extent.py

def contains(self, extent, res=None):
    """Tests if the extent contains another given extent.

    Note:
    -----
    If an optional resolution ('res') is given, the containment value is also
    dependent on whether or not the given extent fits within the larger extent
    AND is situated along the given resolution

    Parameters
    ----------
    extent : Extent
        The Extent object to test for containment

    res : numeric or tuple
        The X & Y resolution to enforce

    Returns
    -------
    bool
    """
    # test raw bounds
    if (
        not extent.srs.IsSame(self.srs)
        or extent.xMin < self.xMin
        or extent.yMin < self.yMin
        or extent.xMax > self.xMax
        or extent.yMax > self.yMax
    ):
        return False

    if res:
        # unpack resolution
        try:
            dx, dy = res
        except:
            dx, dy = res, res

        # Test for factor of resolutions
        thresh = dx / 1000
        if (
            (extent.xMin - self.xMin) % dx > thresh
            or (extent.yMin - self.yMin) % dy > thresh
            or (self.xMax - extent.xMax) % dx > thresh
            or (self.yMax - extent.yMax) % dy > thresh
        ):
            return False
    return True

containsLoc

containsLoc(locs)

Test if the extent contains a location or an iterable of locations.

Parameters:

• locs

  (Anything acceptable to LocationSet()) –

  The locations to be checked

Returns:

• * If a single location is checker: bool –
• * If multiple locations are checked: numpy.ndarray –

Source code in geokit/core/extent.py

def containsLoc(
    self,
    locs,
):
    """Test if the extent contains a location or an iterable of locations.

    Parameters
    ----------
    locs : Anything acceptable to LocationSet()
        The locations to be checked

    Returns
    -------
    * If a single location is checker: bool
    * If multiple locations are checked: numpy.ndarray
    """
    self.box  # initialize the box

    # Normalize locations
    locs = LocationSet(locs).asXY(self.srs)

    # Do tests
    sel = np.ones(locs.shape[0], dtype=bool)
    sel *= locs[:, 0] >= self.xMin
    sel *= locs[:, 0] <= self.xMax
    sel *= locs[:, 1] >= self.yMin
    sel *= locs[:, 1] <= self.yMax

    # Done!
    if sel.size == 1:
        return sel[0]
    else:
        return sel

contoursFromRaster

contoursFromRaster(
    raster,
    contourEdges: list[float],
    transformGeoms: bool = True,
    **kwargs,
) -> DataFrame

Convenience wrapper for geokit.raster.contours which automatically clips a raster to the invoked Extent.

Parameters:

• raster

  (The raster datasource to warp from) –

• contourEdges

  (list[float]) –

  The edges to search for within the raster dataset * This parameter can be set as "None", in which case an additional argument should be given to specify how the edges should be determined - See the documentation of "GDALContourGenerateEx" - Ex. "LEVEL_INTERVAL=10", contourEdges=None

• transformGeoms

  (bool, default: True ) –

  If True, geometries are transformed to the Extent's SRS, otherwise they are left in their native SRS

• kwargs

  –

  Keyword arguments to pass on to the contours function * See geokit.raster.contours

Returns:

• DataFrame –
• With columns: –

  'geom' -> The contiguous-valued geometries 'ID' -> The associated contour edge for each object

Source code in geokit/core/extent.py

def contoursFromRaster(
    self, raster, contourEdges: list[float], transformGeoms: bool = True, **kwargs
) -> pd.DataFrame:
    """Convenience wrapper for geokit.raster.contours which automatically
    clips a raster to the invoked Extent.

    Parameters
    ----------
    raster : The raster datasource to warp from

    contourEdges : list[float]
        The edges to search for within the raster dataset
        * This parameter can be set as "None", in which case an additional
            argument should be given to specify how the edges should be determined
            - See the documentation of "GDALContourGenerateEx"
            - Ex. "LEVEL_INTERVAL=10", contourEdges=None

    transformGeoms : bool
        If True, geometries are transformed to the Extent's SRS, otherwise they
        are left in their native SRS

    kwargs
        Keyword arguments to pass on to the contours function
        * See geokit.raster.contours

    Returns
    -------
    pandas.DataFrame

    With columns:
        'geom' -> The contiguous-valued geometries
        'ID' -> The associated contour edge for each object
    """
    pass
    raster = self.clipRaster(raster)
    geoms = RASTER.contours(raster, contourEdges, **kwargs)

    if transformGeoms:
        geoms.geom = GEOM.transform(geoms.geom, toSRS=self.srs)

    return geoms

corners

corners(asPoints=False)

Returns the four corners of the extent as ogr.gGometry points or as (x,y) coordinates in the extent's srs.

Source code in geokit/core/extent.py

def corners(self, asPoints=False):
    """Returns the four corners of the extent as ogr.gGometry points or as (x,y)
    coordinates in the extent's srs.
    """
    if asPoints:
        # Make corner points
        bl = GEOM.point(self.xMin, self.yMin, srs=self.srs)
        br = GEOM.point(self.xMax, self.yMin, srs=self.srs)
        tl = GEOM.point(self.xMin, self.yMax, srs=self.srs)
        tr = GEOM.point(self.xMax, self.yMax, srs=self.srs)

    else:
        # Make corner points
        bl = (self.xMin, self.yMin)
        br = (self.xMax, self.yMin)
        tl = (self.xMin, self.yMax)
        tr = (self.xMax, self.yMax)

    return (bl, br, tl, tr)

createRaster

createRaster(pixelWidth, pixelHeight, **kwargs)

Convenience function for geokit.raster.createRaster which sets 'bounds' and 'srs' inputs.

• The input resolution MUST fit within the extent

Parameters:

• pixelWidth

  (numeric) –

  The pixel width of the raster in units of the input srs * The keyword 'dx' can be used as well and will override anything given assigned to 'pixelWidth'

• pixelHeight

  (numeric) –

  The pixel height of the raster in units of the input srs * The keyword 'dy' can be used as well and will override anything given assigned to 'pixelHeight'

• **kwargs

  –

  All other keyword arguments are passed on to geokit.raster.createRaster()

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/extent.py

def createRaster(self, pixelWidth, pixelHeight, **kwargs):
    """Convenience function for geokit.raster.createRaster which sets 'bounds'
    and 'srs' inputs.

    * The input resolution MUST fit within the extent

    Parameters
    ----------
    pixelWidth : numeric
        The pixel width of the raster in units of the input srs
        * The keyword 'dx' can be used as well and will override anything given
        assigned to 'pixelWidth'

    pixelHeight : numeric
        The pixel height of the raster in units of the input srs
        * The keyword 'dy' can be used as well and will override anything given
          assigned to 'pixelHeight'

    **kwargs:
        All other keyword arguments are passed on to geokit.raster.createRaster()

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    if not self.fitsResolution((pixelWidth, pixelHeight)):
        raise GeoKitExtentError("The given resolution does not fit to the Extent boundaries")
    return RASTER.createRaster(
        bounds=self.xyXY,
        pixelWidth=pixelWidth,
        pixelHeight=pixelHeight,
        srs=self.srs,
        **kwargs,
    )

drawSmopyMap

drawSmopyMap(
    zoom,
    tileserver="https://a.tile.openstreetmap.org/{z}/{x}/{y}.png",
    tilesize=256,
    maxtiles=100,
    ax=None,
    **kwargs,
)

Draws a basemap using the "smopy" python package.

• See more details about smopy here: https://github.com/rossant/smopy

Returns:

• namedtuple –
  • .ax -> The axes draw on
  • .srs -> The SRS used when drawing (will always be EPSG 3857)
  • .bounds -> The boundaries of the drawn map

Source code in geokit/core/extent.py

def drawSmopyMap(
    self,
    zoom,
    tileserver="https://a.tile.openstreetmap.org/{z}/{x}/{y}.png",
    tilesize=256,
    maxtiles=100,
    ax=None,
    **kwargs,
):
    """
    Draws a basemap using the "smopy" python package.

    * See more details about smopy here: https://github.com/rossant/smopy

    Parameters
    ----------
        zoom : int
            The zoom level to draw (between 1-20)
            * I suggest starting low (e.g. 4), and zooming in until you find a level that suits your needs

        tileserver : string
            The tile server to use

        tilesize : int
            The pixel size of the tiles from 'tileserver'

        maxtiles : int
            The maximum tiles to use when drawing an image
            * Be careful to adhere to the usage conditions stated by your selected tileserver!

        ax : matplotlib.axes
            The matplotlib axes to draw on
            * If 'None', then one will be generated automatically

        kwargs
            All extra keyword arguments are passed on to matplotlib.ax.imshow

    Returns
    -------
        namedtuple
            * .ax     -> The axes draw on
            * .srs    -> The SRS used when drawing (will always be EPSG 3857)
            * .bounds -> The boundaries of the drawn map
    """
    return RASTER.drawSmopyMap(
        bounds=self.castTo(SRS.EPSG4326).xyXY,
        zoom=zoom,
        tileserver=tileserver,
        tilesize=tilesize,
        maxtiles=maxtiles,
        ax=ax,
        **kwargs,
    )

exportWKT

exportWKT(delimiter='|') -> str

Export the extent to a Well-Known_Text string.

• Actually the will be two WKT strings separated by a "|" character
• These correspond to "|"

Parameters:

• delimiter

  (The delimiter which separates the two WKT sections, default: '|' ) –

Returns:

• string –

Source code in geokit/core/extent.py

def exportWKT(self, delimiter="|") -> str:
    """Export the extent to a Well-Known_Text string.

    * Actually the will be two WKT strings separated by a "|" character
    * These correspond to "<A Geometry WKT>|<an SRS WKT>"

    Parameters
    ----------
    delimiter : The delimiter which separates the two WKT sections

    Returns
    -------
    string
    """
    return "{}{}{}".format(self.box.ExportToWkt(), delimiter, self.srs.ExportToWkt())

extractFeatures

extractFeatures(source, **kwargs)

Convenience wrapper for geokit.vector.extractFeatures() by setting the 'geom' input to the extent's box.

Parameters:

• source

  (str) –

  The path to the vector file to load

• **kwargs

  –

  All other keyword arguments are passed on to vector.extractFeatures()

Returns:

• * If asPandas is True: pandas.DataFrame or pandas.Series –
• * If asPandas is False: generator –

Source code in geokit/core/extent.py

def extractFeatures(self, source, **kwargs):
    """Convenience wrapper for geokit.vector.extractFeatures() by setting the
    'geom' input to the extent's box.

    Parameters
    ----------
    source : str
        The path to the vector file to load

    **kwargs:
        All other keyword arguments are passed on to vector.extractFeatures()

    Returns
    -------
    * If asPandas is True: pandas.DataFrame or pandas.Series
    * If asPandas is False: generator
    """
    return VECTOR.extractFeatures(source=source, geom=self._box, **kwargs)

extractMatrix

extractMatrix(source, strict=True, **kwargs)

Convenience wrapper around geokit.raster.extractMatrix(). Extracts the extent directly from the given raster source as a matrix around the Extent.

Note:

The called extent must fit somewhere within the raster's grid

Parameters:

• source

  –

  The raster source to be read

• strict

  –

  Whether or not to allow a returned value which does not fit to the given extent !! If this is set to False, it is STRONGLY recommended to also set the argument 'returnBounds' as True so that the new computed boundary can be known

• **kwargs

  –

  All keyword arguments are passed to geokit.raster.extractMatrix

Returns:

• ndarray or tuple –
  • See geokit.raster.extractMatrix

Source code in geokit/core/extent.py

def extractMatrix(self, source, strict=True, **kwargs):
    """Convenience wrapper around geokit.raster.extractMatrix(). Extracts the
    extent directly from the given raster source as a matrix around the Extent.

    Note:
    -----
    The called extent must fit somewhere within the raster's grid

    Parameters
    ----------
    source: gdal.Dataset or str
        The raster source to be read

    strict: bool; optional
        Whether or not to allow a returned value which does not fit to the
        given extent
        !! If this is set to False, it is STRONGLY recommended to also set the
           argument 'returnBounds' as True so that the new computed boundary
           can be known

    **kwargs
        All keyword arguments are passed to geokit.raster.extractMatrix

    Returns
    -------
    numpy.ndarray or tuple
        * See geokit.raster.extractMatrix
    """
    if strict:
        ri = RASTER.rasterInfo(source)
        if not self.srs.IsSame(ri.srs):
            raise GeoKitExtentError("Extent and source do not share an srs")
        if not Extent._fromInfo(ri).contains(self, (ri.dx, ri.dy)):
            raise GeoKitExtentError("Extent does not fit the raster's resolution")

    return RASTER.extractMatrix(source, bounds=self.xyXY, boundsSRS=self.srs, **kwargs)

filterSources

filterSources(sources, error_on_missing=True)

Filter a list of sources by those whose's envelope overlaps the Extent.

Note:

Creates a filter object which can be immediately iterated over, or else can be cast as a list

Parameters:

• sources

  (list or str) –

  The sources to filter * An iterable of vector/raster sources * An iterable of paths pointing to vector/raster sources * A glob string which will generate a list of source paths - see glob.glob for more info

• error_on_missing

  (bool, default: True ) –

  If True, then if a file path is given which does not exist, a RunTime error is raised. Otherwise a warning is given Only performs check when input is a string

Returns:

• filter –

Source code in geokit/core/extent.py

def filterSources(self, sources, error_on_missing=True):
    """Filter a list of sources by those whose's envelope overlaps the Extent.

    Note:
    -----
    Creates a filter object which can be immediately iterated over, or else
    can be cast as a list

    Parameters
    ----------
    sources : list or str
        The sources to filter
        * An iterable of vector/raster sources
        * An iterable of paths pointing to vector/raster sources
        * A glob string which will generate a list of source paths
            - see glob.glob for more info

    error_on_missing : bool, optional
        If True, then if a file path is given which does not exist, a RunTime
            error is raised. Otherwise a warning is given
        Only performs check when input is a string

    Returns
    -------
    filter
    """
    # create list of searchable files
    if isinstance(sources, str):
        _directoryList = glob(sources)
    else:
        _directoryList = sources

    # Ensure all files exist (only check if input is a string)
    directoryList = []
    for source in _directoryList:
        if isinstance(source, str) and not isfile(source):
            if error_on_missing:
                raise RuntimeError("Cannot find file: " + source)
            else:
                warnings.warn("Skipping missing file: " + source)
        else:
            directoryList.append(source)

    return filter(self.inSourceExtent, directoryList)

findWithin

findWithin(extent, res=100, yAtTop=True)

Finds the indexes of the given extent within the main extent according to the given resolution.

Note:

• Use this to compute the index offsets and window sizes of a window within a raster dataset
• The two extents MUST share the same SRS

Parameters:

• extent

  (Extent) –

  The extent to find within the calling extent

• res

  (numeric or tuple, default: 100 ) –

  A resolution to check containment on

• yAtTop

  (bool; optional, default: True ) –

  Instructs the offsetting to begin from yMax instead of from yMin

Returns:

• tuple -> (xOffset, yOffset, xWindowSize, yWindowSize) –

Source code in geokit/core/extent.py

def findWithin(self, extent, res=100, yAtTop=True):
    """Finds the indexes of the given extent within the main extent according
    to the given resolution.

    Note:
    -----
    * Use this to compute the index offsets and window sizes of a window
      within a raster dataset
    * The two extents MUST share the same SRS

    Parameters
    ----------
    extent : Extent
        The extent to find within the calling extent

    res : numeric or tuple
        A resolution to check containment on

    yAtTop : bool; optional
        Instructs the offsetting to begin from yMax instead of from yMin

    Returns
    -------
        tuple -> (xOffset, yOffset, xWindowSize, yWindowSize)
    """
    # test srs
    if not self.srs.IsSame(extent.srs):
        raise GeoKitExtentError("extents are not of the same srs")

    # try to unpack the resolution
    try:
        dx, dy = res
    except:
        dx, dy = res, res

    # Get offsets
    tmpX = (extent.xMin - self.xMin) / dx
    xOff = int(np.round(tmpX))

    if yAtTop:
        tmpY = (self.yMax - extent.yMax) / dy
    else:
        tmpY = (extent.yMin - self.yMin) / dy
    yOff = int(np.round(tmpY))

    if not (np.isclose(xOff, tmpX) and np.isclose(yOff, tmpY)):
        raise GeoKitExtentError("The extents are not relatable on the given resolution")

    # Get window sizes
    tmpX = (extent.xMax - extent.xMin) / dx
    xWin = int(np.round(tmpX))

    tmpY = (extent.yMax - extent.yMin) / dy
    yWin = int(np.round(tmpY))

    if not (np.isclose(xWin, tmpX) and np.isclose(yWin, tmpY)):
        raise GeoKitExtentError("The extents are not relatable on the given resolution")

    # Done!
    return IndexSet(xOff, yOff, xWin, yWin, xOff + xWin, yOff + yWin)

fit

fit(unit, dtype=None, start_raster=None) -> Extent

Fit the extent to a given pixel resolution.

Note:

The extent is always expanded to fit onto the given unit

Parameters:

• unit

  (numeric or tuple) –

  The unit value(s) to check * If numeric, a single value is assumed for both X and Y dim * If tuple, resolutions for both dimensions (x, y)

• start_raster

  (str(left or right), default: None ) –

  If None passed, the extent will be centered on the shape with overlaps left and right depending on individual shape width and raster cell width. If 'left'/'right' passed, extent edges will be matching min/max longitude of shape.

• dtype

  (Type or dtype, default: None ) –

  The final data type of the boundary values

Returns:

• Extent –

Source code in geokit/core/extent.py

def fit(self, unit, dtype=None, start_raster=None) -> "Extent":
    """Fit the extent to a given pixel resolution.

    Note:
    -----
    The extent is always expanded to fit onto the given unit

    Parameters
    ----------
    unit : numeric or tuple
        The unit value(s) to check
        * If numeric, a single value is assumed for both X and Y dim
        * If tuple, resolutions for both dimensions (x, y)

    start_raster : str ('left' or 'right')
        If None passed, the extent will be centered on the shape with overlaps left and right
        depending on individual shape width and raster cell width. If 'left'/'right' passed,
        extent edges will be matching min/max longitude of shape.

    dtype : Type or np.dtype
        The final data type of the boundary values

    Returns
    -------
    Extent
    """
    try:
        unitX, unitY = unit
    except:
        unitX, unitY = unit, unit

    # Look for bad sizes
    if unitX > self.xMax - self.xMin:
        raise GeoKitExtentError("Unit size is larger than extent width")
    if unitY > self.yMax - self.yMin:
        raise GeoKitExtentError("Unit size is larger than extent width")

    # Calculate new extent
    if start_raster == "left":
        newXMin = self.xMin
        newYMin = np.floor(self.yMin / unitY) * unitY
        newXMax = self.xMin + np.ceil((self.xMax - self.xMin) / unitX) * unitX
        newYMax = np.ceil(self.yMax / unitY) * unitY
    elif start_raster == "right":
        newXMin = self.xMax - np.ceil((self.xMax - self.xMin) / unitX) * unitX
        newYMin = np.floor(self.yMin / unitY) * unitY
        newXMax = self.xMax
        newYMax = np.ceil(self.yMax / unitY) * unitY
    elif start_raster == None:
        newXMin = np.floor(self.xMin / unitX) * unitX
        newYMin = np.floor(self.yMin / unitY) * unitY
        newXMax = np.ceil(self.xMax / unitX) * unitX
        newYMax = np.ceil(self.yMax / unitY) * unitY
    else:
        raise GeoKitExtentError(
            "start_raster parameter must be either 'left' or 'right' (or None). Process terminated."
        )

    # Done!
    if dtype is None or isinstance(unitX, dtype):
        return Extent(newXMin, newYMin, newXMax, newYMax, srs=self.srs)
    else:
        return Extent(
            dtype(newXMin),
            dtype(newYMin),
            dtype(newXMax),
            dtype(newYMax),
            srs=self.srs,
        )

fitsResolution

fitsResolution(unit, tolerance=1e-06)

Test if calling Extent first around the given unit(s) (at least within an error defined by 'tolerance').

Parameters:

• unit

  (numeric or tuple) –

  The unit value(s) to check * If float, a single resolution value is assumed for both X and Y dim * If tuple, resolutions for both dimensions (x, y)

• tolerance

  (float, default: 1e-06 ) –

  The tolerance to allow when comparing float values

Returns:

• Extent –

Examples:

>>> ex = Extent( 100, 100, 300, 500)
>>> ex.fitsResolution(25) # True!
>>> ex.fitsResolution( (25, 10) ) # True!
>>> ex.fitsResolution(33) # False!
>>> ex.fitsResolution( (25, 33) ) # False!

Source code in geokit/core/extent.py

def fitsResolution(self, unit, tolerance=1e-6):
    """Test if calling Extent first around the given unit(s) (at least within
    an error defined by 'tolerance').

    Parameters
    ----------
    unit : numeric or tuple
        The unit value(s) to check
        * If float, a single resolution value is assumed for both X and Y dim
        * If tuple, resolutions for both dimensions (x, y)

    tolerance : float
        The tolerance to allow when comparing float values

    Returns
    -------
    Extent

    Examples
    --------
    >>> ex = Extent( 100, 100, 300, 500)
    >>> ex.fitsResolution(25) # True!
    >>> ex.fitsResolution( (25, 10) ) # True!
    >>> ex.fitsResolution(33) # False!
    >>> ex.fitsResolution( (25, 33) ) # False!
    """
    try:
        unitX, unitY = unit
    except:
        unitX, unitY = unit, unit

    xSteps = (self.xMax - self.xMin) / unitX
    xResidual = abs(xSteps - np.round(xSteps))
    if xResidual > tolerance:
        return False

    ySteps = (self.yMax - self.yMin) / unitY
    yResidual = abs(ySteps - np.round(ySteps))
    if yResidual > tolerance:
        return False

    return True

fromGeom staticmethod

fromGeom(geom) -> Extent

Create extent around a given geometry.

Parameters:

• geom

  (Geometry) –

  The geometry from which to extract the extent

Returns:

• Extent –

Source code in geokit/core/extent.py

@staticmethod
def fromGeom(geom) -> "Extent":
    """Create extent around a given geometry.

    Parameters
    ----------
    geom : ogr.Geometry
        The geometry from which to extract the extent

    Returns
    -------
    Extent
    """
    # Read Envelope
    xMin, xMax, yMin, yMax = geom.GetEnvelope()

    # Done!
    return Extent(xMin, yMin, xMax, yMax, srs=geom.GetSpatialReference())

fromLocationSet staticmethod

fromLocationSet(locs) -> Extent

Create extent around the contents of a LocationSet object.

Parameters:

• locs

  (LocationSet) –

Returns:

• Extent –

Source code in geokit/core/extent.py

@staticmethod
def fromLocationSet(locs) -> "Extent":
    """Create extent around the contents of a LocationSet object.

    Parameters
    ----------
    locs : LocationSet

    Returns
    -------
    Extent
    """
    lonMin, latMin, lonMax, latMax = locs.getBounds()
    return Extent(lonMin, latMin, lonMax, latMax, srs=SRS.EPSG4326)

fromRaster staticmethod

fromRaster(source)

Create extent around the contents of a raster source.

Parameters:

• source

  (Anything acceptable by loadRaster()) –

  The vector datasource to read from

Returns:

• Extent –

Source code in geokit/core/extent.py

@staticmethod
def fromRaster(source):
    """Create extent around the contents of a raster source.

    Parameters
    ----------
    source : Anything acceptable by loadRaster()
        The vector datasource to read from

    Returns
    -------
    Extent
    """
    dsInfo = RASTER.rasterInfo(source)

    xMin, yMin, xMax, yMax = dsInfo.bounds

    return Extent(xMin, yMin, xMax, yMax, srs=dsInfo.srs)

fromTile staticmethod

fromTile(xi: int, yi: int, zoom: int) -> Extent

Generates an Extent corresponding to tiles used for "slippery maps".

Parameters:

• xi

  (int) –

  The tile's X-index - Range depends on zoom value

• yi

  (int) –

  The tile's Y-index - Range depends on zoom value

• zoom

  (int) –

  The tile's zoom index - Range is between 0 and 18

Returns:

• Extent –

Source code in geokit/core/extent.py

@staticmethod
def fromTile(xi: int, yi: int, zoom: int) -> "Extent":
    """Generates an Extent corresponding to tiles used for "slippery maps".

    Parameters
    ----------
    xi : int
        The tile's X-index
        - Range depends on zoom value

    yi : int
        The tile's Y-index
        - Range depends on zoom value

    zoom : int
        The tile's zoom index
        - Range is between 0 and 18

    Returns
    -------
    geokit.Extent
    """
    tl = smopy.num2deg(xi - 0.0, yi + 1.0, zoom)[::-1]
    br = smopy.num2deg(xi + 1.0, yi - 0.0, zoom)[::-1]

    o = SRS.xyTransform([tl, br], fromSRS=SRS.EPSG4326, toSRS=SRS.EPSG3857, outputFormat="xy")
    assert isinstance(o, TransformedPointsXY)
    return Extent(o.x.min(), o.y.min(), o.x.max(), o.y.max(), srs=SRS.EPSG3857)

fromTileAt staticmethod

fromTileAt(x, y, zoom, srs) -> Extent

Generates an Extent corresponding to tiles used for "slippery maps" at the coordinates ('x','y') in the 'srs' reference system.

Parameters:

• x

  (float) –

  The X coordinate to search for a tile around

• y

  (float) –

  The Y coordinate to search for a tile around

• zoom

  (int) –

  The tile's zoom index - Range is between 0 and 18

• srs

  (anything acceptable to SRS.loadSRS) –

  The SRS of the given 'x' & 'y' coordinates

Returns:

• Extent –

Source code in geokit/core/extent.py

@staticmethod
def fromTileAt(x, y, zoom, srs) -> "Extent":
    """Generates an Extent corresponding to tiles used for "slippery maps"
    at the coordinates ('x','y') in the 'srs' reference system.

    Parameters
    ----------
    x : float
        The X coordinate to search for a tile around

    y : float
        The Y coordinate to search for a tile around

    zoom : int
        The tile's zoom index
        - Range is between 0 and 18

    srs : anything acceptable to SRS.loadSRS
        The SRS of the given 'x' & 'y' coordinates

    Returns
    -------
    geokit.Extent
    """
    t = SRS.tileIndexAt(x=x, y=y, zoom=zoom, srs=srs)

    return Extent.fromTile(t.xi, t.yi, t.zoom)

fromVector staticmethod

fromVector(
    source,
    where: str | None = None,
    geom: Geometry | None = None,
) -> Extent

Create extent around the contemts of a vector source.

Parameters:

• source

  (Anything acceptable by loadVector()) –

  The vector datasource to read from

• where

  (str; optional, default: None ) –

  An SQL-style filtering string * Can be used to filter the input source according to their attributes * For tips, see "http://www.gdal.org/ogr_sql.html" Ex: where="eye_color='Green' AND IQ>90"

• geom

  (ogr.Geometry; optional, default: None ) –

  The geometry to search within * All features are extracted which touch this Geometry

Returns:

• Extent –

Source code in geokit/core/extent.py

@staticmethod
def fromVector(source, where: str | None = None, geom: ogr.Geometry | None = None) -> "Extent":
    """Create extent around the contemts of a vector source.

    Parameters
    ----------
    source : Anything acceptable by loadVector()
        The vector datasource to read from

    where : str; optional
        An SQL-style filtering string
        * Can be used to filter the input source according to their attributes
        * For tips, see "http://www.gdal.org/ogr_sql.html"
        Ex:
          where="eye_color='Green' AND IQ>90"

    geom : ogr.Geometry; optional
        The geometry to search within
        * All features are extracted which touch this Geometry

    Returns
    -------
    Extent
    """
    if where is None and geom is None:
        shapeDS = VECTOR.loadVector(source)

        shapeLayer = shapeDS.GetLayer()
        shapeSRS = shapeLayer.GetSpatialRef()

        xMin, xMax, yMin, yMax = shapeLayer.GetExtent()

        return Extent(xMin, yMin, xMax, yMax, srs=shapeSRS)
    else:
        geom = VECTOR.extractFeature(source=source, where=where, geom=geom, onlyGeom=True)
        return Extent.fromGeom(geom)

fromWKT staticmethod

fromWKT(wkt, delimiter='|') -> Extent

Create extent from a Well-Known_Text string.

• Actually the input should be two WKT strings separated by a "|" character
• These correspond to "|"

Parameters:

• wkt

  (The string to be processed) –

• delimiter

  (The delimiter which separates the two WKT sections, default: '|' ) –

Returns:

• Extent –

Source code in geokit/core/extent.py

@staticmethod
def fromWKT(wkt, delimiter="|") -> "Extent":
    """Create extent from a Well-Known_Text string.

    * Actually the input should be two WKT strings separated by a "|" character
    * These correspond to "<A Geometry WKT>|<an SRS WKT>"

    Parameters
    ----------
    wkt : The string to be processed

    delimiter : The delimiter which separates the two WKT sections

    Returns
    -------
    Extent
    """
    geomWKT, srsWKT = wkt.split(delimiter)
    srs = SRS.loadSRS(srsWKT)
    geom = GEOM.convertWKT(geomWKT, srs=srs)

    return Extent.fromGeom(geom)

from_xXyY staticmethod

from_xXyY(bounds, srs='latlon') -> Extent

Create an Extent from explicitly defined boundaries.

Parameters:

• bounds

  (tuple) –

  The (xMin, xMax, yMin, yMax) values for the extent

• srs

  (Anything acceptable to geokit.srs.loadSRS(); optional, default: 'latlon' ) –

  The srs of the input coordinates * If not given, lat/lon coordinates are assumed

Returns:

• Extent –

Source code in geokit/core/extent.py

@staticmethod
def from_xXyY(bounds, srs="latlon") -> "Extent":
    """Create an Extent from explicitly defined boundaries.

    Parameters
    ----------
    bounds : tuple
        The (xMin, xMax, yMin, yMax) values for the extent

    srs : Anything acceptable to geokit.srs.loadSRS(); optional
        The srs of the input coordinates
          * If not given, lat/lon coordinates are assumed

    Returns
    -------
    Extent
    """
    return Extent(bounds[0], bounds[2], bounds[1], bounds[3], srs=srs)

inSourceExtent

inSourceExtent(source)

Tests if the extent box is at least partially contained in the extent-box of the given vector or raster source.

Parameters:

• sources

  (str) –

  The sources to test

Source code in geokit/core/extent.py

def inSourceExtent(self, source):
    """Tests if the extent box is at least partially contained in the extent-box
    of the given vector or raster source.

    Parameters
    ----------
    sources : str
        The sources to test
    """
    try:
        # cover the case that the source file is an empty or single-point vector = no extent possible
        _tmp = VECTOR.extractFeatures(source)
        if len(_tmp) == 1 and _tmp.geom.iloc[0].GetGeometryName() == "POINT":
            # we have a single point, add minimal tolerance in each direction to enable an extent
            tol = 0.00001
            sourceExtent = Extent(
                _tmp.geom.iloc[0].GetX() * (1 - tol),
                _tmp.geom.iloc[0].GetY() * (1 - tol),
                _tmp.geom.iloc[0].GetX() * (1 + tol),
                _tmp.geom.iloc[0].GetY() * (1 + tol),
            )
            return self.overlaps(sourceExtent)
        elif len(_tmp) == 0:
            # source is an empty vector file, overlap is not possible
            return False
        else:
            # fall back to normal path below in except
            raise Exception()
    except:
        sourceExtent = Extent.load(source)
        return self.overlaps(sourceExtent)

load staticmethod

load(source, **kwargs) -> Extent

Attempts to load an Extent from a variety of inputs in the most appropriate manner.

One Extent initializer (.fromXXX) is called depending on the inputs

source is LocationSet -> Extent.fromLocationSet( source )
source is ogr.Geometry -> Extent.fromGeom( source )
source is not a string: -> Extent(*source, **kwargs)
source is a string:
    First try: Extent.fromVector(source)
    Then try: Extent.fromRaster(source)

If none of the above works, an error is raised

Returns:

• Extent –

Source code in geokit/core/extent.py

@staticmethod
def load(source, **kwargs) -> "Extent":
    """Attempts to load an Extent from a variety of inputs in the most
    appropriate manner.

    One Extent initializer (.fromXXX) is called depending on the inputs

        source is LocationSet -> Extent.fromLocationSet( source )
        source is ogr.Geometry -> Extent.fromGeom( source )
        source is not a string: -> Extent(*source, **kwargs)
        source is a string:
            First try: Extent.fromVector(source)
            Then try: Extent.fromRaster(source)

    If none of the above works, an error is raised

    Returns
    -------
    Extent
    """
    if isinstance(source, Extent):
        return source
    elif isinstance(source, LocationSet):
        return Extent.fromLocationSet(source)
    elif isinstance(source, ogr.Geometry):
        return Extent.fromGeom(source)
    elif UTIL.isVector(source):
        return Extent.fromVector(source)
    elif UTIL.isRaster(source):
        return Extent.fromRaster(source)
    elif isinstance(source, str):
        return Extent.fromWKT(source)

    try:  # Maybe the source is an iterable giving xyXY
        vals = list(source)
        if len(vals) == 4:
            return Extent(vals, **kwargs)
    except:
        pass

    raise GeoKitExtentError("Could not load the source")

mutateRaster

mutateRaster(
    source: load_raster_input,
    pixelWidth: numeric | None = None,
    pixelHeight: numeric | None = None,
    matchContext: bool = False,
    warpArgs: dict | None = None,
    processor: Callable | None = None,
    resampleAlg: Literal[
        "near", "bilinear", "cubic", "average"
    ] = "bilinear",
    **mutateArgs,
)

Convenience function for geokit.raster.mutateRaster which automatically warps the raster to the extent's area and srs before mutating.

Note:

If this is called without any arguments except for a source, it serves to clip the raster source around the Extent, therefore performing the same function as Extent.warp(...) on an Extent which has been cast to the source's srs

Parameters:

• source

  (Anything acceptable to geokit.raster.loadRaster()) –

  The source to mutate

• pixelHeight

  (numeric, default: None ) –

  The pixel height (y-resolution) of the output raster

• pixelWidth

  (numeric, default: None ) –

  The pixel width (x-resolution) of the output raster

• matchContext

  (bool; optional, default: False ) –
  • If True, Warp to the Extent's boundaries and srs before mutating
    • pixelHeight and pixelWidth MUST be provided in this case
  • If False, only warp to the Extent's boundaries, but keep its srs and resolution intact

• warpArgs

  (dict; optional, default: None ) –

  Arguments to apply to the warping step * See geokit.raster.warp()

• processor

  (Callable | None, default: None ) –

  The function performing the mutation of the raster's data * The function will take single argument (a 2D numpy.ndarray) * The function must return a numpy.ndarray of the same size as the input * The return type must also be containable within a Float32 (int and boolean is okay) * See example in geokit.raster.mutateRaster for more info

• resampleAlg

  (str; optional, default: 'bilinear' ) –

  The resampling algorithm to use while warping * Knowing which option to use can have significant impacts! * Options are: 'near', 'bilinear', 'cubic', 'average'

• **kwargs

  –

  All other keyword arguments are passed to geokit.vector.mutateVector

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/extent.py

def mutateRaster(
    self,
    source: load_raster_input,
    pixelWidth: numeric | None = None,
    pixelHeight: numeric | None = None,
    matchContext: bool = False,
    warpArgs: dict | None = None,
    processor: Callable | None = None,
    resampleAlg: Literal["near", "bilinear", "cubic", "average"] = "bilinear",
    **mutateArgs,
):
    """Convenience function for geokit.raster.mutateRaster which automatically
    warps the raster to the extent's area and srs before mutating.

    Note:
    -----
    If this is called without any arguments except for a source, it serves
    to clip the raster source around the Extent, therefore performing
    the same function as Extent.warp(...) on an Extent which has been cast
    to the source's srs

    Parameters
    ----------
    source : Anything acceptable to geokit.raster.loadRaster()
        The source to mutate

    pixelHeight : numeric
        The pixel height (y-resolution) of the output raster

    pixelWidth : numeric
        The pixel width (x-resolution) of the output raster

    matchContext : bool; optional
        * If True, Warp to the Extent's boundaries and srs before mutating
            - pixelHeight and pixelWidth MUST be provided in this case
        * If False, only warp to the Extent's boundaries, but keep its
          srs and resolution intact

    warpArgs : dict; optional
        Arguments to apply to the warping step
        * See geokit.raster.warp()

    processor - function; optional
        The function performing the mutation of the raster's data
        * The function will take single argument (a 2D numpy.ndarray)
        * The function must return a numpy.ndarray of the same size as the input
        * The return type must also be containable within a Float32 (int and
          boolean is okay)
        * See example in geokit.raster.mutateRaster for more info

    resampleAlg : str; optional
        The resampling algorithm to use while warping
        * Knowing which option to use can have significant impacts!
        * Options are: 'near', 'bilinear', 'cubic', 'average'

    **kwargs:
        All other keyword arguments are passed to geokit.vector.mutateVector

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    if warpArgs is None:
        warpArgs = {}

    if processor is None:  # We won't do a mutation without a processor, since everything else
        # can be handled by Warp. Therefore we pass on any 'output' that is
        # given to the warping stage, unless one was already given
        warpArgs["output"] = warpArgs.get("output", mutateArgs.get("output", None))

    # Warp the source
    # TODO: Should the warping be updated to use Extent.clipRaster???
    if matchContext:
        if pixelWidth is None or pixelHeight is None:
            raise GeoKitExtentError("pixelWidth and pixelHeight must be provided when matchContext is True")

        source = self.warp(
            source=source,
            resampleAlg=resampleAlg,
            pixelWidth=pixelWidth,
            pixelHeight=pixelWidth,
            strict=True,
            **warpArgs,
        )
    else:
        if not "srs" in mutateArgs:
            source = RASTER.loadRaster(source)
            srs = source.GetProjectionRef()

        ext = self.castTo(srs)
        source = ext.warp(
            source=source,
            resampleAlg=resampleAlg,
            pixelWidth=pixelWidth,
            pixelHeight=pixelWidth,
            strict=False,
            **warpArgs,
        )

    # mutate the source
    if not processor is None:
        return RASTER.mutateRaster(source, processor=processor, **mutateArgs)
    else:
        return source

mutateVector

mutateVector(source, matchContext=False, **kwargs)

Convenience function for geokit.vector.mutateVector which automatically sets 'srs' and 'geom' input to the Extent's srs and geometry.

Note:

If this is called without any arguments except for a source, it serves to clip the vector source around the extent

Parameters:

• source

  (Anything acceptable to geokit.vector.loadVector()) –

  The source to clip

• matchContext

  (bool; optional, default: False ) –
  • If True, transforms all geometries to the Extent's srs before mutating
  • If False, the Extent is cast to the source's srs, and all filtering and mutating happens in that context

• **kwargs

  –

  All other keyword arguments are passed to geokit.vector.mutateVector

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/extent.py

def mutateVector(self, source, matchContext=False, **kwargs):
    """Convenience function for geokit.vector.mutateVector which automatically
    sets 'srs' and 'geom' input to the Extent's srs and geometry.

    Note:
    -----
    If this is called without any arguments except for a source, it serves
    to clip the vector source around the extent

    Parameters
    ----------
    source : Anything acceptable to geokit.vector.loadVector()
        The source to clip

    matchContext : bool; optional
        * If True, transforms all geometries to the Extent's srs before
          mutating
        * If False, the Extent is cast to the source's srs, and all filtering
          and mutating happens in that context

    **kwargs:
        All other keyword arguments are passed to geokit.vector.mutateVector

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    # Get the working srs
    if not matchContext:
        vinfo = VECTOR.vectorInfo(source)
        ext = self.castTo(vinfo.srs)
    else:
        ext = self

    # mutate the source
    return VECTOR.mutateVector(source, srs=ext.srs, geom=ext._box, **kwargs)

overlaps

overlaps(extent, referenceSRS=EPSG4326)

Tests if the extent overlaps with another given extent.

Note:

If an optional resolution ('res') is given, the containment value is also dependent on whether or not the given extent fits within the larger extent AND is situated along the given resolution

Parameters:

• extent

  (Extent) –

  The Extent object to test for containment

• referenceSRS

  –

  The spatial reference frame to do the comparison in * Can be 'self'

Returns:

• bool –

Source code in geokit/core/extent.py

def overlaps(self, extent, referenceSRS=SRS.EPSG4326):
    """Tests if the extent overlaps with another given extent.

    Note:
    -----
    If an optional resolution ('res') is given, the containment value is also
    dependent on whether or not the given extent fits within the larger extent
    AND is situated along the given resolution

    Parameters
    ----------
    extent : Extent
        The Extent object to test for containment

    referenceSRS
        The spatial reference frame to do the comparison in
        * Can be 'self'

    Returns
    -------
    bool
    """
    if referenceSRS != "self":
        self = self.castTo(referenceSRS)
    else:
        referenceSRS = self.srs
    extent = extent.castTo(referenceSRS)

    if self.box.Intersects(extent.box):
        return True
    if extent.box.Intersects(self.box):
        return True

    return False

pad

pad(pad, percent=False) -> Extent

Pad the extent in all directions.

Parameters:

• pad

  (float) –

  The amount to pad in all directions * In units of the extent's srs * Can also accept a negative padding

• percent

  (bool, default: False ) –

  If True, the padding values are understood to be a percentage of the unpadded extent

Returns:

• Extent –

Source code in geokit/core/extent.py

def pad(self, pad, percent=False) -> "Extent":
    """Pad the extent in all directions.

    Parameters
    ----------
    pad : float
        The amount to pad in all directions
        * In units of the extent's srs
        * Can also accept a negative padding

    percent : bool, optional
        If True, the padding values are understood to be a percentage of the
        unpadded extent

    Returns
    -------
    Extent
    """
    # Check for no input pads
    if pad is None:
        return self

    # try breaking apart by x and y component
    try:
        xpad, ypad = pad
    except:
        xpad = pad
        ypad = pad

    if percent:
        xpad = xpad / 100 * (self.xMax - self.xMin) / 2
        ypad = ypad / 100 * (self.yMax - self.yMin) / 2

    # Pad the extent
    return Extent(
        self.xMin - xpad,
        self.yMin - ypad,
        self.xMax + xpad,
        self.yMax + ypad,
        srs=self.srs,
    )

rasterMosaic

rasterMosaic(
    sources: list[load_raster_input],
    resampleAlg: Literal[
        "near", "bilinear", "cubic", "average"
    ] = "near",
    _warpKwargs={},
    _skipFiltering: bool = False,
)

Create a raster source surrounding the Extent from a collection of other rasters.

Parameters:

• sources

  (list, or something acceptable to gk.Extent.filterSources) –

  The sources to add together over the invoking Extent

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/extent.py

def rasterMosaic(
    self,
    sources: list[load_raster_input],
    resampleAlg: Literal["near", "bilinear", "cubic", "average"] = "near",
    _warpKwargs={},
    _skipFiltering: bool = False,
):
    """Create a raster source surrounding the Extent from a collection of other rasters.

    Parameters
    ----------
    sources : list, or something acceptable to gk.Extent.filterSources
        The sources to add together over the invoking Extent

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    if _skipFiltering:
        sources = sorted(list(sources))
    else:
        sources = sorted(list(self.filterSources(sources)))

    if len(sources) == 0:
        warnings.warn("No suitable sources found")
        return None

    raster_info = RASTER.rasterInfo(sources[0])

    ext = self.castTo(raster_info.srs).fit((raster_info.dx, raster_info.dy))

    master_raster = ext._quickRaster(
        dx=raster_info.pixelWidth,
        dy=raster_info.pixelHeight,
        noData=raster_info.noData,
        scale=raster_info.scale,
        offset=raster_info.offset,
        dtype=raster_info.data_type_name_str,
    )
    gdal.Warp(
        master_raster,
        sources,
        resampleAlg=resampleAlg,
        **_warpKwargs,
    )

    return master_raster

rasterize

rasterize(
    source, pixelWidth, pixelHeight, strict=True, **kwargs
)

Convenience function for geokit.vector.rasterize() which automatically sets the 'srs' and 'bounds' input.

Note:

When creating an 'in memory' raster vs one which is saved to disk, a slightly different algorithm is used which can sometimes add an extra row of pixels. Be aware of this if you intend to compare value-matricies directly from rasters generated with this function.

Parameters:

• source

  (str) –

  The path to the vector file to load

• pixelHeight

  (numeric; optional) –

  The pixel height (y-resolution) of the output raster * Only required if this value should be changed

• pixelWidth

  (numeric; optional) –

  The pixel width (x-resolution) of the output raster * Only required if this value should be changed

• strict

  (bool, default: True ) –

  If True, raise an error if trying to rasterize to a pixelWidth and pixelHeight which does not fit into the Extent

• **kwargs

  –

  All other keyword arguments are passed on to geokit.raster.warp()

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/extent.py

def rasterize(self, source, pixelWidth, pixelHeight, strict=True, **kwargs):
    """Convenience function for geokit.vector.rasterize() which automatically
    sets the 'srs' and 'bounds' input.

    Note:
    -----
    When creating an 'in memory' raster vs one which is saved to disk, a slightly
    different algorithm is used which can sometimes add an extra row of pixels. Be
    aware of this if you intend to compare value-matricies directly from rasters
    generated with this function.

    Parameters
    ----------
    source : str
        The path to the vector file to load

    pixelHeight : numeric; optional
        The pixel height (y-resolution) of the output raster
        * Only required if this value should be changed

    pixelWidth : numeric; optional
        The pixel width (x-resolution) of the output raster
        * Only required if this value should be changed

    strict : bool
        If True, raise an error if trying to rasterize to a pixelWidth and
        pixelHeight which does not fit into the Extent

    **kwargs:
        All other keyword arguments are passed on to geokit.raster.warp()

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    if strict and not self.fitsResolution((pixelWidth, pixelHeight)):
        raise GeoKitExtentError("The given resolution does not fit to the Extent boundaries")
    return VECTOR.rasterize(
        source=source,
        pixelWidth=pixelWidth,
        pixelHeight=pixelHeight,
        srs=self.srs,
        bounds=self.xyXY,
        **kwargs,
    )

shift

shift(dx=0, dy=0) -> Extent

Shift the extent in the X and/or Y dimensions.

Parameters:

• dx

  (float, default: 0 ) –

  The amount to shift in the x dimension * In units of the extent's srs

• dy

  (float, default: 0 ) –

  The amount to shift in the y dimension * In units of the extent's srs

Returns:

• Extent –

Source code in geokit/core/extent.py

def shift(self, dx=0, dy=0) -> "Extent":
    """Shift the extent in the X and/or Y dimensions.

    Parameters
    ----------
    dx : float
        The amount to shift in the x dimension
        * In units of the extent's srs

    dy : float
        The amount to shift in the y dimension
        * In units of the extent's srs

    Returns
    -------
    Extent
    """
    return Extent(self.xMin + dx, self.yMin + dy, self.xMax + dx, self.yMax + dy, srs=self.srs)

subTiles

subTiles(zoom, asGeom=False)

Generates tile Extents at a given zoom level which encompass the invoking Extent.

Parameters:

• zoom

  (int) –

  The zoom level of the expected tile source

• asGeom

  (bool, default: False ) –

  If True, returns tuple of ogr.Geometries in stead of (xi,yi,zoom) tuples

Returns:

• Generator of Geometries or (xi,yi,zoom) tuples –

Source code in geokit/core/extent.py

def subTiles(self, zoom, asGeom=False):
    """Generates tile Extents at a given zoom level which encompass the invoking Extent.

    Parameters
    ----------
    zoom : int
        The zoom level of the expected tile source

    asGeom : bool
        If True, returns tuple of ogr.Geometries in stead of (xi,yi,zoom) tuples

    Returns
    -------
    Generator of Geometries or (xi,yi,zoom) tuples
    """
    yield from GEOM.subTiles(self.box, zoom, checkIntersect=False, asGeom=asGeom)

tileBox

tileBox(zoom: int, return_index_box: bool = False)

Determine the tile Extent at a given zoom level which surround the invoked Extent.

Parameters:

• zoom

  (int) –

  The zoom level of the expected tile source

• return_index_box

  (bool, default: False ) –

  If true, also return the index box at the specified zoom level (from self.tileIndexBox)

Returns:

• if return_index_box is False: geokit.Extent –
• if return_index_box is True: Tuple –
  • Item 0: geokit.Extent
  • Item 1: namedtuple(xi_start, xi_stop, yi_start, yi_stop)

Source code in geokit/core/extent.py

def tileBox(self, zoom: int, return_index_box: bool = False):
    """Determine the tile Extent at a given zoom level which surround the invoked Extent.

    Parameters
    ----------
    zoom : int
        The zoom level of the expected tile source

    return_index_box : bool
        If true, also return the index box at the specified zoom level (from self.tileIndexBox)

    Returns
    -------
    if return_index_box is False: geokit.Extent

    if return_index_box is True: Tuple
        - Item 0: geokit.Extent
        - Item 1: namedtuple(xi_start, xi_stop, yi_start, yi_stop)
    """
    # Get Bounds of new raster in EPSG3857
    tb = self.tileIndexBox(zoom)

    tl_tile_xi = tb.xi_start
    tl_tile_yi = tb.yi_start
    br_tile_xi = tb.xi_stop
    br_tile_yi = tb.yi_stop

    tl_lat, tl_lon = smopy.num2deg(tl_tile_xi, tl_tile_yi, zoom)
    br_lat, br_lon = smopy.num2deg(br_tile_xi + 1, br_tile_yi + 1, zoom)

    coords3857 = SRS.xyTransform(
        [(tl_lon, tl_lat), (br_lon, br_lat)],
        fromSRS=SRS.EPSG4326,
        toSRS=SRS.EPSG3857,
        outputFormat="xy",
    )

    ext = Extent(
        coords3857.x.min(),
        coords3857.y.min(),
        coords3857.x.max(),
        coords3857.y.max(),
        srs=SRS.EPSG3857,
    )

    if return_index_box:
        return ext, tb
    else:
        return ext

tileIndexBox

tileIndexBox(zoom: int)

Determine the tile indexes at a given zoom level which surround the invoked Extent.

Parameters:

• zoom

  (int) –

  The zoom level of the expected tile source

Returns:

• namedtuple –
  • xi_start: int - The starting x index
  • xi_stop: int - The ending x index
  • yi_start: int - The starting y index
  • yi_stop: int - The ending y index

Source code in geokit/core/extent.py

def tileIndexBox(self, zoom: int):
    """Determine the tile indexes at a given zoom level which surround the invoked Extent.

    Parameters
    ----------
    zoom : int
        The zoom level of the expected tile source

    Returns
    -------
    namedtuple:
        - xi_start: int - The starting x index
        - xi_stop:  int - The ending x index
        - yi_start: int - The starting y index
        - yi_stop:  int - The ending y index
    """
    ext4326 = self.castTo(SRS.EPSG4326)

    tl_tile_xi, tl_tile_yi = smopy.deg2num(ext4326.yMax, ext4326.xMin, zoom)
    br_tile_xi, br_tile_yi = smopy.deg2num(ext4326.yMin, ext4326.xMax, zoom)

    return TileIndexBox(
        xi_start=tl_tile_xi,
        xi_stop=br_tile_xi,
        yi_start=tl_tile_yi,
        yi_stop=br_tile_yi,
        zoom=zoom,
    )

tileMosaic

tileMosaic(source: str, zoom: int, **kwargs)

Create a raster source surrounding the Extent from a collection of tiles.

Parameters:

• source

  (str) –

  The source to fetch tiles from * Must include indicators for: {z} -> The tile's zoom level {x} -> The tile's x-index {y} -> The tile's y-index * Ex: File on disk : "/path/to/tile/directory/{z}/{x}/{y}/filename.tif" Remote HTTP file : "/vsicurl_streaming/http://path/to/resource/{z}/{x}/{y}/filename.tif" * Find more info at https://gdal.org/user/virtual_file_systems.html

• zoom

  (int) –

  The zoom level of the expected tile source

• pixelsPerTile

  ((int, (int, int))) –

  The number of pixels found in each tile

• workingType

  (dtype) –

  The datatype of the working matrix (should match the raster source)

• noData

  (numeric) –

  The value to treat as 'no data'

• output

  (str) –

  An optional path for an output raster (.tif) file

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/extent.py

def tileMosaic(self, source: str, zoom: int, **kwargs):
    """Create a raster source surrounding the Extent from a collection of tiles.

    Parameters
    ----------
    source : str
        The source to fetch tiles from
        * Must include indicators for:
          {z} -> The tile's zoom level
          {x} -> The tile's x-index
          {y} -> The tile's y-index
        * Ex:
          File on disk     : "/path/to/tile/directory/{z}/{x}/{y}/filename.tif"
          Remote HTTP file : "/vsicurl_streaming/http://path/to/resource/{z}/{x}/{y}/filename.tif"
        * Find more info at https://gdal.org/user/virtual_file_systems.html

    zoom : int
        The zoom level of the expected tile source

    pixelsPerTile : int, (int,int)
        The number of pixels found in each tile

    workingType : np.dtype
        The datatype of the working matrix (should match the raster source)

    noData : numeric
        The value to treat as 'no data'

    output : str
        An optional path for an output raster (.tif) file

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    sources = list(self.tileSources(zoom=zoom, source=source))
    return self.rasterMosaic(sources, _skipFiltering=True, **kwargs)

tileSources

tileSources(zoom: int, source: str | None = None)

Get the tiles sources which contribute to the invoking Extent.

Parameters:

• zoom

  (int) –

  The zoom level of the expected tile source

• source

  (str, default: None ) –

  The source to fetch tiles from * Must include indicators for: {z} -> The tile's zoom level {x} -> The tile's x-index {y} -> The tile's y-index * Ex: File on disk : "/path/to/tile/directory/{z}/{x}/{y}/filename.tif" Remote HTTP file : "/vsicurl_streaming/http://path/to/resource/{z}/{x}/{y}/filename.tif" * Find more info at https://gdal.org/user/virtual_file_systems.html

Yields:

• if source is given: str –
• if source is not given: (xi,yi,zoom) –

Source code in geokit/core/extent.py

def tileSources(self, zoom: int, source: str | None = None):
    """Get the tiles sources which contribute to the invoking Extent.

    Parameters
    ----------
    zoom : int
        The zoom level of the expected tile source

    source : str
        The source to fetch tiles from
        * Must include indicators for:
          {z} -> The tile's zoom level
          {x} -> The tile's x-index
          {y} -> The tile's y-index
        * Ex:
          File on disk     : "/path/to/tile/directory/{z}/{x}/{y}/filename.tif"
          Remote HTTP file : "/vsicurl_streaming/http://path/to/resource/{z}/{x}/{y}/filename.tif"
        * Find more info at https://gdal.org/user/virtual_file_systems.html

    Yields
    ------
    if source is given:     str
    if source is not given: (xi,yi,zoom)
    """
    tb = self.tileIndexBox(zoom=zoom)
    for xi in range(tb.xi_start, tb.xi_stop + 1):
        for yi in range(tb.yi_start, tb.yi_stop + 1):
            if source is None:
                yield (xi, yi, zoom)
            else:
                yield (source.replace("{z}", str(zoom)).replace("{x}", str(xi)).replace("{y}", str(yi)))

warp

warp(
    source: load_raster_input,
    pixelWidth,
    pixelHeight,
    strict=True,
    **kwargs,
)

Convenience function for geokit.raster.warp() which automatically sets the 'srs' and 'bounds' input.

Note:

When creating an 'in memory' raster vs one which is saved to disk, a slightly different algorithm is used which can sometimes add an extra row of pixels. Be aware of this if you intend to compare value-matricies directly from rasters generated with this function.

Parameters:

• source

  (str) –

  The path to the vector file to load

• pixelHeight

  (numeric; optional) –

  The pixel height (y-resolution) of the output raster * Only required if this value should be changed

• pixelWidth

  (numeric; optional) –

  The pixel width (x-resolution) of the output raster * Only required if this value should be changed

• strict

  (bool, default: True ) –

  If True, raise an error if trying to warp to a pixelWidth and pixelHeight which does not fit into the Extent

• **kwargs

  –

  All other keyword arguments are passed on to geokit.raster.warp()

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/extent.py

def warp(self, source: load_raster_input, pixelWidth, pixelHeight, strict=True, **kwargs):
    """Convenience function for geokit.raster.warp() which automatically sets the
    'srs' and 'bounds' input.

    Note:
    -----
    When creating an 'in memory' raster vs one which is saved to disk, a slightly
    different algorithm is used which can sometimes add an extra row of pixels. Be
    aware of this if you intend to compare value-matricies directly from rasters
    generated with this function.

    Parameters
    ----------
    source : str
        The path to the vector file to load

    pixelHeight : numeric; optional
        The pixel height (y-resolution) of the output raster
        * Only required if this value should be changed

    pixelWidth : numeric; optional
        The pixel width (x-resolution) of the output raster
        * Only required if this value should be changed

    strict : bool
        If True, raise an error if trying to warp to a pixelWidth and
        pixelHeight which does not fit into the Extent

    **kwargs:
        All other keyword arguments are passed on to geokit.raster.warp()

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    if strict and not self.fitsResolution((pixelWidth, pixelHeight)):
        raise GeoKitExtentError("The given resolution does not fit to the Extent boundaries")
    return RASTER.warp(
        source=source,
        pixelWidth=pixelWidth,
        pixelHeight=pixelHeight,
        srs=self.srs,
        bounds=self.xyXY,
        **kwargs,
    )

RegionMask

RegionMask(
    extent: Extent,
    pixelRes,
    mask=None,
    geom=None,
    attributes=None,
    **kwargs,
)

The RegionMask object represents a given region and exposes methods allowing for easy manipulation of geospatial data around that region.

RegionMask objects are defined by providing a polygon (either via a vector file, an ogr.Geometry, or a Well-Known-Text (WKT) string), a projection system to work in, and an extent and pixel resolution to create a matrix mask (i.e. boolean values) of.

• The extent of the generated mask matrix is the tightest fit around the region in units of the pixel resolution. However, the extent can be defined explicitly if desired
• The region can be manipulated as a vector polygon via the ".geometry" attribute, which exposes the geometry as an ogr.Geometry. To incorporate this into other vector-handeling libraries it is suggested to use the ".ExportToWkt()" method available via OGR.
• The region can be manipulated as a raster matrix via the ".mask" attribute which exposes the mask as a boolean numpy.ndarray
• Any raster source can be easily warped onto the region-mask's extent, projection, and resolution via the ".warp" method
• Any vector source can be rasterized onto the region-mask's extent, projection, and resolution via the ".rasterize" method
• The default mask set-up is defined by the constant members: DEFAULT_SRS, DEFAULT_RES, and DEFAULT_PAD

Initializers:

• RegionMask(...)

  • This is not the preferred way

• RegionMask.fromVector( ... )

• RegionMask.fromVectorFeature( ... )

• RegionMask.fromGeom( ... )

• RegionMask.fromMask( ... )

• RegionMask.load( ... )

  • This function tries to determine which of the other initializers should be used based off the input

The default constructor for RegionMask objects. Creates a RegionMask directly from a matrix mask and a given extent (and optionally a geometry). Pixel resolution is calculated in accordance with the shape of the mask mask and the provided extent.

• Generally one should use the '.load' or else one of the '.fromXXX' methods to create RegionMasks

Parameters:

• extent

  (Extent object) –

  The geospatial context of the region mask * The extent must fit the given pixel sizes * All computations using the RegionMask will be evaluated within this spatial context

• pixelRes

  (float or tuple) –

  The RegionMask's native pixel size(s) * If float : A pixel size to apply to both the X and Y dimension * If (float float) : An X-dimension and Y-dimension pixel size * All computations using the RegionMask will generate results in reference to these pixel sizes (i.e. either at this resolution or at some scaling of this resolution)

• mask

  (numpy - ndarray, default: None ) –

  A mask over the context area defining which pixel as inside the region and which are outside * Must be a 2-Dimensional bool-matrix describing the region, where: - 0/False -> "not in the region" - 1/True -> "Inside the region" * Either a mask or a geometry must be given, but not both

• geom

  (ogr - Geomertry, default: None ) –

  A geometric representation of the RegionMask's region * Either a mask or a geometry must be given, but not both

• attributes

  (dict, default: None ) –

  Keyword attributes and values to carry along with the RegionMask

Methods:

• applyMask –

  Shortcut to apply the RegionMask's mask to an array. Mainly intended

• buildGeometry –

  Explicitly build the RegionMask's geometry.

• buildMask –

  Explicitly build the RegionMask's mask matrix.

• contoursFromMask –

  Convenience wrapper for geokit.raster.contours which automatically

• contoursFromMatrix –

  Convenience wrapper for geokit.raster.contours which automatically

• contoursFromRaster –

  Convenience wrapper for geokit.raster.contours which automatically

• createRaster –

  Convenience wrapper for geokit.raster.createRaster which sets 'srs',

• drawGeoms –

  Convenience wrapper around geokit.geom.drawGeoms which plots geometries

• drawImage –

  Convenience wrapper around geokit.util.drawImage which plots matrix data

• drawMask –

  Convenience wrapper around geokit.util.drawImage which plots the

• drawRaster –

  Convenience wrapper around geokit.raster.drawRaster which plots a raster

• drawSelf –

  Convenience wrapper around geokit.geom.drawGeoms which plots the

• extractFeatures –

  Convenience wrapper for geokit.vector.extractFeatures() by setting the

• fromGeom –

  Make a RasterMask from a given geometry.

• fromMask –

  Make a RegionMask directly from a mask matrix and extent.

• fromVector –

  Make a RasterMask from a given vector source.

• indicateFeatures –

  Indicates the RegionMask pixels which are found within the features (or

• indicateGeoms –

  Convenience wrapper to indicate values found within a geometry (or a

• indicateValueToGeoms –

  TODO: UPDATE ME.

• indicateValues –

  Indicates those pixels in the RegionMask which correspond to a particular

• load –

  Tries to initialize and return a RegionMask in the most appropriate way.

• mutateRaster –

  Convenience wrapper for geokit.vector.mutateRaster which automatically

• mutateVector –

  Convenience wrapper for geokit.vector.mutateVector which automatically

• polygonizeMask –

  Convenience wrapper for geokit.geom.polygonizeMask which automatically

• polygonizeMatrix –

  Convenience wrapper for geokit.geom.polygonizeMatrix which automatically

• rasterize –

  Convenience wrapper for geokit.vector.rasterize() which automatically

• subRegions –

  Generate a number of sub regions on a grid which combine into the total

• subTiles –

  Generates tile Extents at a given zoom level which encompass the invoking Regionmask.

• warp –

  Convenience wrapper for geokit.raster.warp() which automatically sets

Attributes:

• geometry –

  Fetches a clone of the RegionMask's geometry as an OGR Geometry object.

• mask –

  The RegionMask's mask array as an 2-dimensional boolean numpy array.

• pixelRes –

  The RegionMask's pixel size.

• vector –

  Returns a vector saved in memory which is built only once.

• vectorPath –

  Returns a path to a vector path on disc which is built only once.

Source code in geokit/core/regionmask.py

def __init__(self, extent: Extent, pixelRes, mask=None, geom=None, attributes=None, **kwargs):
    """The default constructor for RegionMask objects. Creates a RegionMask
    directly from a matrix mask and a given extent (and optionally a geometry).
    Pixel resolution is calculated in accordance with the shape of the mask
    mask and the provided extent.

    * Generally one should use the '.load' or else one of the '.fromXXX'
      methods to create RegionMasks

    Parameters
    ----------
    extent : Extent object
        The geospatial context of the region mask
        * The extent must fit the given pixel sizes
        * All computations using the RegionMask will be evaluated within this
          spatial context

    pixelRes : float or tuple
        The RegionMask's native pixel size(s)
        * If float : A pixel size to apply to both the X and Y dimension
        * If (float float) : An X-dimension and Y-dimension pixel size
        * All computations using the RegionMask will generate results in
          reference to these pixel sizes (i.e. either at this resolution or
          at some scaling of this resolution)

    mask : numpy-ndarray
        A mask over the context area defining which pixel as inside the region
        and which are outside
        * Must be a 2-Dimensional bool-matrix describing the region, where:
            - 0/False -> "not in the region"
            - 1/True  -> "Inside the region"
        * Either a mask or a geometry must be given, but not both

    geom : ogr-Geomertry
        A geometric representation of the RegionMask's region
        * Either a mask or a geometry must be given, but not both

    attributes : dict
        Keyword attributes and values to carry along with the RegionMask
    """
    # Check for bad input
    if mask is None and geom is None:
        raise GeoKitRegionMaskError("Either mask or geom should be defined")
    if not kwargs.get("mask_plus_geom_is_okay", False):
        if not mask is None and not geom is None:
            raise GeoKitRegionMaskError("mask and geom cannot be defined simultaneously")

    # Set basic values
    self.extent: Extent = extent
    self.srs: osr.SpatialReference = extent.srs

    if self.srs is None:
        raise GeoKitRegionMaskError("Extent SRS cannot be None")

    # Set Pixel Size)
    if not extent.fitsResolution(pixelRes):
        raise GeoKitRegionMaskError("The given extent does not fit the given pixelRes")

    try:
        pixelWidth, pixelHeight = pixelRes
    except:
        pixelWidth, pixelHeight = pixelRes, pixelRes

    self.pixelWidth = abs(pixelWidth)
    self.pixelHeight = abs(pixelHeight)

    if self.pixelHeight == self.pixelWidth:
        self._pixelRes = self.pixelHeight
    else:
        self._pixelRes = None

    # set height and width
    # It turns out that I can't set these values here, since sometimes gdal
    # functions can add an extra row (due to float comparison issues?) when
    ## warping and rasterizing
    # int(np.round((self.extent.xMax-s.extent.xMin)/s.pixelWidth))
    self.width = None
    # int(np.round((self.extent.yMax-s.extent.yMin)/s.pixelHeight))
    self.height = None

    # Set mask
    self._mask = mask
    if not mask is None:  # test the mask
        # test type
        if mask.dtype != "bool" and mask.dtype != "uint8":
            raise GeoKitRegionMaskError("Mask must be bool type")
        if mask.dtype == "uint8":
            mask = mask.astype("bool")

        if not np.isclose(extent.xMin + pixelWidth * mask.shape[1], extent.xMax) or not np.isclose(
            extent.yMin + pixelHeight * mask.shape[0], extent.yMax
        ):
            raise GeoKitRegionMaskError("Extent and pixels sizes do not correspond to mask shape")

    # Set geometry
    if not geom is None:  # test the geometry
        if not isinstance(geom, ogr.Geometry):
            raise GeoKitRegionMaskError("geom is not an ogr.Geometry object")

        self._geometry = geom.Clone()
        gSRS = geom.GetSpatialReference()
        if gSRS is None:
            raise GeoKitRegionMaskError("geom does not have an srs")

        if not gSRS.IsSame(self.srs):
            self._geometry = GEOM.transform(self._geometry, toSRS=self.srs, fromSRS=gSRS)
    else:
        self._geometry = None

    # Set other containers
    self._vector = None
    self._vectorPath = None

    # set attributes
    self.attributes = {} if attributes is None else attributes

geometry property

geometry

Fetches a clone of the RegionMask's geometry as an OGR Geometry object.

• If a geometry was not provided when the RegionMask was initialized, then one will be generated from the RegionMask's mask matrix in the RegionMask's extent
• The geometry can always be deleted and rebuild using the RegionMask.rebuildGeometry() function

mask property

mask

The RegionMask's mask array as an 2-dimensional boolean numpy array.

• If no mask was given at the time of the RegionMask's creation, then a mask will be generated on first access to the 'mask' property
• The mask can be rebuilt in a customized way using the RegionMask.buildMask() function

pixelRes property

pixelRes

The RegionMask's pixel size.

!!Only available when pixelWidth equals pixelHeight!!

vector property

vector

Returns a vector saved in memory which is built only once.

vectorPath property

vectorPath

Returns a path to a vector path on disc which is built only once.

applyMask

applyMask(mat, noData=0)

Shortcut to apply the RegionMask's mask to an array. Mainly intended for internal use.

• When the passed matrix does not have the same extents of the given matrix, it is assumed that the RegionMask's mask needs to be scaled so that the matrix dimensions match

• The RM's mask can only be scaled UP, and the given matrix's dimensions must be multiples of the mask's dimensions

Parameters:

• mat

  (ndarray) –

  The matrix to apply the mask to * Must have dimensions equal, or are multiples of, the mask's

• noData

  (float, default: 0 ) –

  The no-data value to set into matrix's values which are not within the region

Returns:

• ndarray –

Source code in geokit/core/regionmask.py

def applyMask(self, mat, noData=0):
    """Shortcut to apply the RegionMask's mask to an array. Mainly intended
    for internal use.

    * When the passed matrix does not have the same extents of the given matrix,
      it is assumed that the RegionMask's mask needs to be scaled so that the
      matrix dimensions match

    * The RM's mask can only be scaled UP, and the given matrix's dimensions
      must be multiples of the mask's dimensions

    Parameters
    ----------
    mat : np.ndarray
        The matrix to apply the mask to
        * Must have dimensions equal, or are multiples of, the mask's

    noData : float
        The no-data value to set into matrix's values which are not within
        the region

    Returns
    -------
    numpy.ndarray
    """
    if noData is None:
        noData = 0

    # set matrix datatype to float if float noData value (like e.g. nan) is passed
    if isinstance(noData, float):
        mat = mat.astype(np.float64)

    # Get size
    Y, X = mat.shape

    # make output array
    out = np.array(mat)

    # Apply mask
    if self.mask.shape == mat.shape:  # matrix dimensions coincide with mask's data
        out[~self.mask] = noData

    elif Y > self.height and X > self.width:
        if not Y % self.height == 0 or not X % self.width == 0:
            raise GeoKitRegionMaskError("Matrix dimensions must be multiples of mask dimensions")

        yScale = Y // self.height
        xScale = X // self.width

        scaledMask = UTIL.scaleMatrix(self.mask, (yScale, xScale))
        sel = np.where(~scaledMask)
        out[sel] = noData

    else:
        raise GeoKitRegionMaskError("Could not map mask onto matrix")

    return out

buildGeometry

buildGeometry()

Explicitly build the RegionMask's geometry.

Source code in geokit/core/regionmask.py

def buildGeometry(self):
    """Explicitly build the RegionMask's geometry."""
    if self._mask is None:
        raise GeoKitRegionMaskError("Cannot build geometry when mask is None")

    self._geometry = None
    self._geometry = GEOM.polygonizeMask(self.mask, bounds=self.extent.xyXY, srs=self.extent.srs, flat=True)

buildMask

buildMask(**kwargs)

Explicitly build the RegionMask's mask matrix.

• The 'width' and 'height' attributes for the RegionMask are also set when this function is called
• All kwargs are passed on to a call to geokit.vector.rasterize()

Source code in geokit/core/regionmask.py

def buildMask(self, **kwargs):
    """Explicitly build the RegionMask's mask matrix.

    * The 'width' and 'height' attributes for the RegionMask are also set
      when this function is called
    * All kwargs are passed on to a call to geokit.vector.rasterize()
    """
    if self._geometry is None:
        raise GeoKitRegionMaskError("Cannot build mask when geometry is None")

    self._mask = None
    self._mask = self.rasterize(self.vectorPath, applyMask=False, **kwargs).astype(bool)
    self.height, self.width = self._mask.shape

contoursFromMask

contoursFromMask(
    mask,
    truthThreshold=0.5,
    trueAboveThreshold=True,
    contoursKwargs={},
    createRasterKwargs={},
)

Convenience wrapper for geokit.raster.contours which automatically creates a raster for the given mask (which is assumed to match the domain of the RegionMask), and extracts the geometries which are indicated in the mask as "True".

Parameters:

• mask

  (matrix_like) –

  The mask which will be turned into a geometry set * Must be 2 dimensional

• truthThreshold

  ([float], default: 0.5 ) –

  The value which separates "True" from "False" values * Values are True when they are above the threshold unless trueAboveThreshold is set as False

• trueAboveThreshold

  –

  If true, then pixels with values above the threshold are identified as "True"

• contoursKwargs

  (dict, default: {} ) –

  Keyword arguments to pass on to the contours function * See geokit.raster.contours

• createRasterKwargs

  (dict, default: {} ) –

  Keyword arguments to pass on to the raster creation function * See geokit.RegionMask.createRaster

Returns:

• DataFrame –
• With columns: –

  'geom' -> The contiguous-valued geometries 'ID' -> The associated contour edge for each object

Source code in geokit/core/regionmask.py

def contoursFromMask(
    self,
    mask,
    truthThreshold=0.5,
    trueAboveThreshold=True,
    contoursKwargs={},
    createRasterKwargs={},
):
    """Convenience wrapper for geokit.raster.contours which automatically
    creates a raster for the given mask (which is assumed to match the
    domain of the RegionMask), and extracts the geometries which are indicated
    in the mask as "True".

    Parameters
    ----------
    mask : matrix_like
        The mask which will be turned into a geometry set
          * Must be 2 dimensional

    truthThreshold : [float,]
        The value which separates "True" from "False" values
        * Values are True when they are above the threshold unless
          trueAboveThreshold is set as False

    trueAboveThreshold: bool
        If true, then pixels with values above the threshold are identified
        as "True"

    contoursKwargs : dict
        Keyword arguments to pass on to the contours function
        * See geokit.raster.contours

    createRasterKwargs : dict
        Keyword arguments to pass on to the raster creation function
        * See geokit.RegionMask.createRaster

    Returns
    -------
    pandas.DataFrame

    With columns:
        'geom' -> The contiguous-valued geometries
        'ID' -> The associated contour edge for each object
    """
    geomDF = self.contoursFromMatrix(
        matrix=mask,
        contourEdges=[
            truthThreshold,
        ],
        createRasterKwargs=createRasterKwargs,
        contoursKwargs=contoursKwargs,
    )

    geoms = geomDF[geomDF.ID == (1 if trueAboveThreshold else 0)].geom

    return geoms

contoursFromMatrix

contoursFromMatrix(
    matrix,
    contourEdges,
    contoursKwargs={},
    createRasterKwargs={},
)

Convenience wrapper for geokit.raster.contours which automatically creates a raster for the given matrix (which is assumed to match the domain of the RegionMask).

Parameters:

• matrix

  (matrix_like) –

  The matrix which will be turned into a geometry set * Must be 2 dimensional

• contourEdges

  ([float]) –

  The edges to search for within the raster dataset * This parameter can be set as "None", in which case an additional argument should be given to specify how the edges should be determined - See the documentation of "GDALContourGenerateEx" - Ex. "LEVEL_INTERVAL=10", contourEdges=None

• contoursKwargs

  (dict, default: {} ) –

  Keyword arguments to pass on to the contours function * See geokit.raster.contours

• createRasterKwargs

  (dict, default: {} ) –

  Keyword arguments to pass on to the raster creation function * See geokit.RegionMask.createRaster

Returns:

• DataFrame –
• With columns: –

  'geom' -> The contiguous-valued geometries 'ID' -> The associated contour edge for each object

Source code in geokit/core/regionmask.py

def contoursFromMatrix(self, matrix, contourEdges, contoursKwargs={}, createRasterKwargs={}):
    """Convenience wrapper for geokit.raster.contours which automatically
    creates a raster for the given matrix (which is assumed to match the
    domain of the RegionMask).

    Parameters
    ----------
    matrix : matrix_like
        The matrix which will be turned into a geometry set
          * Must be 2 dimensional

    contourEdges : [float,]
        The edges to search for within the raster dataset
        * This parameter can be set as "None", in which case an additional
            argument should be given to specify how the edges should be determined
            - See the documentation of "GDALContourGenerateEx"
            - Ex. "LEVEL_INTERVAL=10", contourEdges=None

    contoursKwargs : dict
        Keyword arguments to pass on to the contours function
        * See geokit.raster.contours

    createRasterKwargs : dict
        Keyword arguments to pass on to the raster creation function
        * See geokit.RegionMask.createRaster

    Returns
    -------
    pandas.DataFrame

    With columns:
        'geom' -> The contiguous-valued geometries
        'ID' -> The associated contour edge for each object
    """
    raster = self.createRaster(data=matrix, **createRasterKwargs)
    geoms = RASTER.contours(raster, contourEdges, **contoursKwargs)

    return geoms

contoursFromRaster

contoursFromRaster(
    raster,
    contourEdges,
    applyMask=True,
    contoursKwargs={},
    warpKwargs={},
)

Convenience wrapper for geokit.raster.contours which automatically warps a raster to the invoking RegioNmask.

NOTE:

• The raster is first warped to the RegionMask before the contours are determined. If this behavior is not desired, consider using the function Extent.contoursFromRaster

Parameters:

• raster

  (The raster datasource to warp from) –

• contourEdges

  ([float]) –

  The edges to search for within the raster dataset * This parameter can be set as "None", in which case an additional argument should be given to specify how the edges should be determined - See the documentation of "GDALContourGenerateEx" - Ex. "LEVEL_INTERVAL=10", contourEdges=None

• contoursKwargs

  (dict, default: {} ) –

  Keyword arguments to pass on to the contours function * See geokit.raster.contours

• warpKwargs

  (dict, default: {} ) –

  Keyword arguments to pass on to the raster warp function * See geokit.RegionMask.warp

Returns:

• DataFrame –
• With columns: –

  'geom' -> The contiguous-valued geometries 'ID' -> The associated contour edge for each object

Source code in geokit/core/regionmask.py

def contoursFromRaster(self, raster, contourEdges, applyMask=True, contoursKwargs={}, warpKwargs={}):
    """Convenience wrapper for geokit.raster.contours which automatically
    warps a raster to the invoking RegioNmask.

    NOTE:
    -----
    * The raster is first warped to the RegionMask before the contours are
      determined. If this behavior is not desired, consider using the function
      Extent.contoursFromRaster

    Parameters
    ----------
    raster : The raster datasource to warp from

    contourEdges : [float,]
        The edges to search for within the raster dataset
        * This parameter can be set as "None", in which case an additional
            argument should be given to specify how the edges should be determined
            - See the documentation of "GDALContourGenerateEx"
            - Ex. "LEVEL_INTERVAL=10", contourEdges=None

    contoursKwargs : dict
        Keyword arguments to pass on to the contours function
        * See geokit.raster.contours

    warpKwargs : dict
        Keyword arguments to pass on to the raster warp function
        * See geokit.RegionMask.warp

    Returns
    -------
    pandas.DataFrame

    With columns:
        'geom' -> The contiguous-valued geometries
        'ID' -> The associated contour edge for each object
    """
    raster = self.warp(raster, applyMask=applyMask, returnMatrix=False, **warpKwargs)
    geoms = RASTER.contours(raster, contourEdges, **contoursKwargs)

    return geoms

createRaster

createRaster(output=None, resolutionDiv=1, **kwargs)

Convenience wrapper for geokit.raster.createRaster which sets 'srs', 'bounds', 'pixelWidth', and 'pixelHeight' inputs.

Parameters:

• output

  (str; optional, default: None ) –

  A path to an output file to write to

• resolutionDiv

  (int, default: 1 ) –

  The factor by which to divide the RegionMask's native resolution * This is useful if you need to represent very fine details

• **kwargs

  –

  All other keywargs are passed on to geokit.raster.createRaster() * See below for argument descriptions

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/regionmask.py

def createRaster(self, output=None, resolutionDiv=1, **kwargs):
    """Convenience wrapper for geokit.raster.createRaster which sets 'srs',
    'bounds', 'pixelWidth', and 'pixelHeight' inputs.

    Parameters
    ----------
    output : str; optional
        A path to an output file to write to

    resolutionDiv : int
        The factor by which to divide the RegionMask's native resolution
        * This is useful if you need to represent very fine details

    **kwargs:
        All other keywargs are passed on to geokit.raster.createRaster()
        * See below for argument descriptions

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    pW, pH = self._resolve(resolutionDiv)
    return self.extent.createRaster(pixelWidth=pW, pixelHeight=pH, output=output, **kwargs)

drawGeoms

drawGeoms(
    geoms, ax=None, drawSelf=True, **kwargs
) -> AxHands

Convenience wrapper around geokit.geom.drawGeoms which plots geometries which are then plotted within the context of the RegionMask.

• See geokit.geom.drawGeoms for more info on argument options
• Geometries are always plotted in the RegionMask's SRS
• Unless specified, x and y limits are set to the RegionMask's extent
  • This only plays a role when generating a new axis

Source code in geokit/core/regionmask.py

def drawGeoms(self, geoms, ax=None, drawSelf=True, **kwargs) -> AxHands:
    """Convenience wrapper around geokit.geom.drawGeoms which plots geometries
    which are then plotted within the context of the RegionMask.

    * See geokit.geom.drawGeoms for more info on argument options
    * Geometries are always plotted in the RegionMask's SRS
    * Unless specified, x and y limits are set to the RegionMask's extent
        - This only plays a role when generating a new axis
    """
    xlim = kwargs.pop("xlim", (self.extent.xMin, self.extent.xMax))
    ylim = kwargs.pop("ylim", (self.extent.yMin, self.extent.yMax))
    ax = GEOM.drawGeoms(geoms, ax=ax, srs=self.srs, xlim=xlim, ylim=ylim, **kwargs)

    if drawSelf:
        self.drawSelf(ax=ax, fc="None", ec="k", linewidth=2)
    return ax

drawImage

drawImage(
    matrix, ax_hands=None, drawSelf=True, **kwargs
) -> AxHands

Convenience wrapper around geokit.util.drawImage which plots matrix data which is assumed to match the boundaries of the RegionMask.

• See geokit.util.drawImage for more info on argument options
• Unless specified, the plotting extent is set to the RegionMask's extent
  • This only plays a role when generating a new axis

Source code in geokit/core/regionmask.py

def drawImage(self, matrix, ax_hands=None, drawSelf=True, **kwargs) -> AxHands:
    """Convenience wrapper around geokit.util.drawImage which plots matrix data
    which is assumed to match the boundaries of the RegionMask.

    * See geokit.util.drawImage for more info on argument options
    * Unless specified, the plotting extent is set to the RegionMask's extent
        - This only plays a role when generating a new axis
    """
    xlim = kwargs.pop("xlim", (self.extent.xMin, self.extent.xMax))
    ylim = kwargs.pop("ylim", (self.extent.yMin, self.extent.yMax))

    ax_hands = UTIL.drawImage(matrix, ax=ax_hands, xlim=xlim, ylim=ylim, **kwargs)

    if drawSelf:
        self.drawSelf(ax=ax_hands, fc="None", ec="k", linewidth=2)
    return ax_hands

drawMask

drawMask(ax=None, **kwargs)

Convenience wrapper around geokit.util.drawImage which plots the RegionMask's mask over the RegionMask's context.

• See geokit.util.drawImage for more info on argument options
• Unless specified, the plotting extent is set to the RegionMask's extent
  • This only plays a role when generating a new axis

Source code in geokit/core/regionmask.py

def drawMask(self, ax=None, **kwargs):
    """Convenience wrapper around geokit.util.drawImage which plots the
    RegionMask's mask over the RegionMask's context.

    * See geokit.util.drawImage for more info on argument options
    * Unless specified, the plotting extent is set to the RegionMask's extent
        - This only plays a role when generating a new axis
    """
    xlim = kwargs.pop("xlim", (self.extent.xMin, self.extent.xMax))
    ylim = kwargs.pop("ylim", (self.extent.yMin, self.extent.yMax))
    return UTIL.drawImage(self.mask, ax=ax, xlim=xlim, ylim=ylim, **kwargs)

drawRaster

drawRaster(ax=None, drawSelf=True, **kwargs)

Convenience wrapper around geokit.raster.drawRaster which plots a raster dataset within the context of the RegionMask.

• See geokit.raster.drawRaster for more info on argument options
• The raster is always warped to the RegionMask's SRS
• Unless specified, x and y limits are set to the RegionMask's extent
  • This only plays a role when generating a new axis

Source code in geokit/core/regionmask.py

def drawRaster(self, ax=None, drawSelf=True, **kwargs):
    """Convenience wrapper around geokit.raster.drawRaster which plots a raster
    dataset within the context of the RegionMask.

    * See geokit.raster.drawRaster for more info on argument options
    * The raster is always warped to the RegionMask's SRS
    * Unless specified, x and y limits are set to the RegionMask's extent
        - This only plays a role when generating a new axis
    """
    xlim = kwargs.pop("xlim", (self.extent.xMin, self.extent.xMax))
    ylim = kwargs.pop("ylim", (self.extent.yMin, self.extent.yMax))
    ax = GEOM.drawGeoms(self.geometry, ax=ax, srs=self.srs, xlim=xlim, ylim=ylim, **kwargs)

    if drawSelf:
        self.drawSelf(ax=ax, fc="None", ec="k", linewidth=2)
    return ax

drawSelf

drawSelf(ax=None, **kwargs)

Convenience wrapper around geokit.geom.drawGeoms which plots the RegionMask's geometry.

• See geokit.geom.drawGeoms for more info on argument options
• Geometry are always plotted in the RegionMask's SRS
• Unless specified, x and y limits are set to the RegionMask's extent
  • This only plays a role when generating a new axis

Source code in geokit/core/regionmask.py

def drawSelf(self, ax=None, **kwargs):
    """Convenience wrapper around geokit.geom.drawGeoms which plots the
    RegionMask's geometry.

    * See geokit.geom.drawGeoms for more info on argument options
    * Geometry are always plotted in the RegionMask's SRS
    * Unless specified, x and y limits are set to the RegionMask's extent
        - This only plays a role when generating a new axis
    """
    xlim = kwargs.pop("xlim", (self.extent.xMin, self.extent.xMax))
    ylim = kwargs.pop("ylim", (self.extent.yMin, self.extent.yMax))
    return GEOM.drawGeoms(self.geometry, ax=ax, srs=self.srs, xlim=xlim, ylim=ylim, **kwargs)

extractFeatures

extractFeatures(source, **kwargs)

Convenience wrapper for geokit.vector.extractFeatures() by setting the 'geom' input to the RegionMask's geometry.

Parameters:

• source

  (str) –

  The path to the vector file to load

• **kwargs

  –

  All other keyword arguments are passed on to vector.extractFeatures()

Returns:

• * If asPandas is True: pandas.DataFrame or pandas.Series –
• * If asPandas is False: generator –

Source code in geokit/core/regionmask.py

def extractFeatures(self, source, **kwargs):
    """Convenience wrapper for geokit.vector.extractFeatures() by setting the
    'geom' input to the RegionMask's geometry.

    Parameters
    ----------
    source : str
        The path to the vector file to load

    **kwargs:
        All other keyword arguments are passed on to vector.extractFeatures()

    Returns
    -------
    * If asPandas is True: pandas.DataFrame or pandas.Series
    * If asPandas is False: generator
    """
    return VECTOR.extractFeatures(source=source, geom=self.geometry, **kwargs)

fromGeom staticmethod

fromGeom(
    geom,
    srs,
    pixelRes=DEFAULT_RES,
    start_raster=None,
    extent=None,
    padExtent=DEFAULT_PAD,
    attributes=None,
    **region_mask_key_worded_arugments,
) -> RegionMask

Make a RasterMask from a given geometry.

Parameters:

• geom

  (ogr - Geomertry or str) –

  A geometric representation of the RegionMask's region * If a string is given, geokit.geom.convertWKT(geom, srs) is called to convert it to an ogr.Geometry

• pixelRes

  (float or tuple, default: DEFAULT_RES ) –

  The RegionMask's native pixel resolution(s) * If float : A pixel size to apply to both the X and Y dimension * If (float float) : An X-dimension and Y-dimension pixel size

• srs

  (Anything acceptable to geokit.srs.loadSRS()) –

  The srs context of the generated RegionMask object * This srs is superseded by the srs in an explicitly defined extent

• extent

  (Extent object, default: None ) –

  The geospatial context of the generated region mask * The extent must fit the given pixel sizes

• padExtent

  (float; optional, default: DEFAULT_PAD ) –

  An amount by which to pad the extent before generating the RegionMask

• attributes

  (dict, default: None ) –

  Keyword attributes and values to carry along with the RegionMask

Returns:

• RegionMask –

Source code in geokit/core/regionmask.py

@staticmethod
def fromGeom(
    geom,
    srs,
    pixelRes=DEFAULT_RES,
    start_raster=None,
    extent=None,
    padExtent=DEFAULT_PAD,
    attributes=None,
    **region_mask_key_worded_arugments,
) -> "RegionMask":
    """Make a RasterMask from a given geometry.

    Parameters
    ----------
    geom : ogr-Geomertry or str
        A geometric representation of the RegionMask's region
        * If a string is given, geokit.geom.convertWKT(geom, srs) is called
          to convert it to an ogr.Geometry

    pixelRes : float or tuple
        The RegionMask's native pixel resolution(s)
        * If float : A pixel size to apply to both the X and Y dimension
        * If (float float) : An X-dimension and Y-dimension pixel size

    srs : Anything acceptable to geokit.srs.loadSRS()
        The srs context of the generated RegionMask object
        * This srs is superseded by the srs in an explicitly defined extent

    extent : Extent object
        The geospatial context of the generated region mask
        * The extent must fit the given pixel sizes

    padExtent : float; optional
        An amount by which to pad the extent before generating the RegionMask

    attributes : dict
        Keyword attributes and values to carry along with the RegionMask

    Returns
    -------
    RegionMask
    """
    # make sure we have a geometry with an srs
    if isinstance(geom, str):
        assert not (isinstance(srs, str) and srs.upper() == "LAEA"), "srs cannot be LAEA when geom is WKT"
        srs = SRS.loadSRS(srs)
        geom = GEOM.convertWKT(geom, srs)
    elif isinstance(srs, str) and srs.upper() == "LAEA":
        # First convert to epsg 4326 to get lat lon centroid then generate centered LAEA
        geom = GEOM.transform(geoms=geom, toSRS=4326)
        srs = SRS.loadSRS(source="laea", geom=geom)
    else:
        srs = SRS.loadSRS(srs)

    # clone to make sure we're free of outside dependencies
    # convert to regionmask srs to ensure that rm.geometry.GetSpatialReference() is equal to rm.srs
    geom = GEOM.transform(geom.Clone(), toSRS=srs)

    # set extent (if not given)
    if extent is None:
        extent = Extent.fromGeom(geom).castTo(srs).pad(padExtent).fit(unit=pixelRes, start_raster=start_raster)
    else:
        if not extent.srs.IsSame(srs):
            raise GeoKitRegionMaskError("The given srs does not match the extent's srs")
        # extent = extent.pad(padExtent)

    # make a RegionMask object
    return RegionMask(
        extent=extent, pixelRes=pixelRes, geom=geom, attributes=attributes, **region_mask_key_worded_arugments
    )

fromMask staticmethod

fromMask(extent, mask, attributes=None) -> RegionMask

Make a RegionMask directly from a mask matrix and extent.

Note:

Pixel sizes are calculated from the extent boundaries and mask dimensional sizes

Parameters:

• extent

  (Extent object) –

  The geospatial context of the region mask * The extent must fit the given pixel sizes * All computations using the RegionMask will be evaluated within this spatial context

• mask

  (numpy - ndarray) –

  A mask over the context area defining which pixel as inside the region and which are outside * Must be a 2-Dimensional bool-matrix describing the region, where: - 0/False -> "not in the region" - 1/True -> "Inside the region"

• attributes

  (dict, default: None ) –

  Keyword attributes and values to carry along with the RegionMask

Returns:

• RegionMask –

Source code in geokit/core/regionmask.py

@staticmethod
def fromMask(extent, mask, attributes=None) -> "RegionMask":
    """Make a RegionMask directly from a mask matrix and extent.

    Note:
    -----
    Pixel sizes are calculated from the extent boundaries and mask dimensional
    sizes

    Parameters
    ----------
    extent : Extent object
        The geospatial context of the region mask
        * The extent must fit the given pixel sizes
        * All computations using the RegionMask will be evaluated within this
          spatial context

    mask : numpy-ndarray
        A mask over the context area defining which pixel as inside the region
        and which are outside
        * Must be a 2-Dimensional bool-matrix describing the region, where:
            - 0/False -> "not in the region"
            - 1/True  -> "Inside the region"

    attributes : dict
        Keyword attributes and values to carry along with the RegionMask

    Returns
    -------
    RegionMask
    """
    # get pixelWidth and pixelHeight
    pixelWidth = (extent.xMax - extent.xMin) / (mask.shape[1])
    pixelHeight = (extent.yMax - extent.yMin) / (mask.shape[0])

    return RegionMask(
        extent=extent,
        pixelRes=(pixelWidth, pixelHeight),
        mask=mask,
        attributes=attributes,
    )

fromVector staticmethod

fromVector(
    source,
    srs,
    where=None,
    geom=None,
    start_raster=None,
    pixelRes=DEFAULT_RES,
    extent=None,
    padExtent=DEFAULT_PAD,
    limitOne=True,
    **kwargs,
) -> RegionMask

Make a RasterMask from a given vector source.

Note:

Be careful when creating a RegionMask over a large area (such as a country)! Using the default pixel size for a large area (such as a country) can easily consume your system's memory

Parameters:

• source

  (Anything acceptable by loadVector()) –

  The vector data source to read from

• where

  (str, int; optional, default: None ) –

  If string -> An SQL-like where statement to apply to the source If int -> The feature's ID within the dataset * Feature attribute name do not need quotes * String values should be wrapped in 'single quotes' Example: If the source vector has a string attribute called "ISO" and a integer attribute called "POP", you could use....

  where = "ISO='DEU' AND POP>1000"
  

• geom

  (ogr.Geometry; optional, default: None ) –

  The geometry to search with * All features are extracted which touch this geometry

• pixelRes

  (float or tuple, default: DEFAULT_RES ) –

  The RegionMask's native pixel resolution(s) * If float : A pixel size to apply to both the X and Y dimension * If (float float) : An X-dimension and Y-dimension pixel size

• srs

  (Anything acceptable to geokit.srs.loadSRS()) –

  The srs context of the generated RegionMask object * This srs is superseded by the srs in an explicitly defined extent

• extent

  (Extent object, default: None ) –

  The geospatial context of the generated region mask * The extent must fit the given pixel sizes * If not specified, the entire extent of the vector file is assumed

• padExtent

  (float; optional, default: DEFAULT_PAD ) –

  An amount by which to pad the extent before generating the RegionMask

• limitOne

  (bool; optional, default: True ) –

  Whether or not to allow more than one feature to be extracted

Returns:

• RegionMask –

Source code in geokit/core/regionmask.py

@staticmethod
def fromVector(
    source,
    srs,
    where=None,
    geom=None,
    start_raster=None,
    pixelRes=DEFAULT_RES,
    extent=None,
    padExtent=DEFAULT_PAD,
    limitOne=True,
    **kwargs,
) -> "RegionMask":
    """Make a RasterMask from a given vector source.

    Note:
    -----
    Be careful when creating a RegionMask over a large area (such as a country)!
    Using the default pixel size for a large area (such as a country) can
    easily consume your system's memory

    Parameters
    ----------
    source : Anything acceptable by loadVector()
        The vector data source to read from

    where : str, int; optional
        If string -> An SQL-like where statement to apply to the source
        If int -> The feature's ID within the dataset
        * Feature attribute name do not need quotes
        * String values should be wrapped in 'single quotes'
        Example: If the source vector has a string attribute called "ISO" and
                 a integer attribute called "POP", you could use....

            where = "ISO='DEU' AND POP>1000"

    geom : ogr.Geometry; optional
        The geometry to search with
        * All features are extracted which touch this geometry

    pixelRes : float or tuple
        The RegionMask's native pixel resolution(s)
        * If float : A pixel size to apply to both the X and Y dimension
        * If (float float) : An X-dimension and Y-dimension pixel size

    srs : Anything acceptable to geokit.srs.loadSRS()
        The srs context of the generated RegionMask object
        * This srs is superseded by the srs in an explicitly defined extent

    extent : Extent object
        The geospatial context of the generated region mask
        * The extent must fit the given pixel sizes
        * If not specified, the entire extent of the vector file is assumed

    padExtent : float; optional
        An amount by which to pad the extent before generating the RegionMask

    limitOne : bool; optional
        Whether or not to allow more than one feature to be extracted

    Returns
    -------
    RegionMask
    """
    # Get all geoms which fit the search criteria
    if isinstance(where, int):
        geom, attr = VECTOR.extractFeature(source=source, where=where, srs=4326)
    else:
        ftrs = list(VECTOR.extractFeatures(source=source, where=where, srs=4326, asPandas=False))

        if len(ftrs) == 0:
            raise GeoKitRegionMaskError("Zero features found")
        elif len(ftrs) == 1:
            geom = ftrs[0].geom
            attr = ftrs[0].attr
        else:
            if limitOne:
                raise GeoKitRegionMaskError(
                    "Multiple features found. If you are okay with this, set 'limitOne' to False"
                )
            geom = GEOM.flatten([f.geom for f in ftrs])
            attr = None
    srs = SRS.loadSRS(srs, geom=geom)
    geom = GEOM.transform(geoms=geom, toSRS=srs)

    # Done!
    return RegionMask.fromGeom(
        geom,
        extent=extent,
        start_raster=start_raster,
        pixelRes=pixelRes,
        attributes=attr,
        padExtent=padExtent,
        srs=srs,
        **kwargs,
    )

indicateFeatures

indicateFeatures(
    source: load_vector_input,
    where: str | None = None,
    buffer: numeric | None = None,
    bufferMethod: Literal[
        "geom", "area", "contour"
    ] = "geom",
    resolutionDiv: int = 1,
    forceMaskShape: bool = False,
    applyMask: bool = True,
    noData: numeric = 0,
    preBufferSimplification: numeric | None = None,
    multiProcess: bool = True,
    **kwargs,
)

Indicates the RegionMask pixels which are found within the features (or a subset of the features) contained in a given vector datasource.

• A Rasterization is performed from the input data set to the RegionMask's mask. -See geokit.vector.rasterize or, more specifically gdal.RasterizeOptions kwargs for more info on how to control the rasterization step

Parameters:

• source

  (str or Dataset) –

  The vector datasource to indicate from

• where

  (str; optional, default: None ) –

  An SQL-style filtering string * Can be used to filter the input source according to their attributes * For tips, see "http://www.gdal.org/ogr_sql.html" Ex: where="eye_color='Green' AND IQ>90"

• buffer

  (float; optional, default: None ) –

  A buffer region to add around the indicated pixels * Units are in the RegionMask's srs

• bufferMethod

  (str; optional, default: 'geom' ) –

  An indicator determining the method to use when buffering * Options are: 'geom', 'area', and 'contour' * If 'geom', the function will attempt to grow each of the geometries directly using the ogr library - This can fail sometimes when the geometries are particularly complex or if some of the geometries are not valid (as in, they have self-intersections) * If 'area', the function will first rasterize the raw geometries and will then apply the buffer to the indicated pixels - This is the safer option although is not as accurate as the 'geom' option since it does not capture the exact edges of the geometries - This method can be made more accurate by increasing the 'resolutionDiv' input * If 'contour', the function will still rasterize the raw geometries, but will then create geometries via mask contours (not the explicit pixel edges) - This option will recreate geometries which are more similar to the original geometries compared to the 'area' method - This method can be made more accurate by increasing the 'resolutionDiv' input

• resolutionDiv

  (int; optional, default: 1 ) –

  The factor by which to divide the RegionMask's native resolution * This is useful if you need to represent very fine details

• forceMaskShape

  (bool; optional, default: False ) –

  If True, forces the returned matrix to have the same dimension as the RegionMask's mask regardless of the 'resolutionDiv' argument

• applyMask

  (bool; optional, default: True ) –

  When True, the RegionMask's mask will be applied to the outputData as described by RegionMask.applyMask

• noData

  (numeric, default: 0 ) –

  The noData value to use when applying the mask

• preBufferSimplification

  (numeric | None, default: None ) –

  If given, then geometries will be simplified (using ogr.Geometry.Simplify) using the specified value before being buffered - Using this can drastically decrease the time it takes to perform the buffering procedure, but can decrease accuracy if it is too high

• multiProcess

  (bool, default: True ) –

  If True, multiple parallel processes will be spawned within the function to improve RAM efficiency, else it will fall back on linear execution. By default True.

• kwargs

  –
  • Most notably: 'allTouched'

Returns:

• ndarray –

Source code in geokit/core/regionmask.py

def indicateFeatures(
    self,
    source: load_vector_input,
    where: str | None = None,
    buffer: numeric | None = None,
    bufferMethod: Literal["geom", "area", "contour"] = "geom",
    resolutionDiv: int = 1,
    forceMaskShape: bool = False,
    applyMask: bool = True,
    noData: numeric = 0,
    preBufferSimplification: numeric | None = None,
    multiProcess: bool = True,
    **kwargs,
):
    """
    Indicates the RegionMask pixels which are found within the features (or
    a subset of the features) contained in a given vector datasource.

    * A Rasterization is performed from the input data set to the
    RegionMask's mask.
    -See geokit.vector.rasterize or, more specifically gdal.RasterizeOptions
    kwargs for more info on how to control the rasterization step

    Parameters
    ----------
    source : str or gdal.Dataset
        The vector datasource to indicate from

    where : str; optional
        An SQL-style filtering string
        * Can be used to filter the input source according to their attributes
        * For tips, see "http://www.gdal.org/ogr_sql.html"
        Ex:
        where="eye_color='Green' AND IQ>90"

    buffer : float; optional
        A buffer region to add around the indicated pixels
        * Units are in the RegionMask's srs

    bufferMethod : str; optional
        An indicator determining the method to use when buffering
        * Options are: 'geom', 'area', and 'contour'
        * If 'geom', the function will attempt to grow each of the geometries
        directly using the ogr library
        - This can fail sometimes when the geometries are particularly
            complex or if some of the geometries are not valid (as in, they
            have self-intersections)
        * If 'area', the function will first rasterize the raw geometries and
        will then apply the buffer to the indicated pixels
        - This is the safer option although is not as accurate as the 'geom'
            option since it does not capture the exact edges of the geometries
        - This method can be made more accurate by increasing the
            'resolutionDiv' input
        * If 'contour', the function will still rasterize the raw geometries,
        but will then create geometries via mask contours (not the explicit
        pixel edges)
        - This option will recreate geometries which are more similar to the
            original geometries compared to the 'area' method
        - This method can be made more accurate by increasing the
            'resolutionDiv' input

    resolutionDiv : int; optional
        The factor by which to divide the RegionMask's native resolution
        * This is useful if you need to represent very fine details

    forceMaskShape : bool; optional
        If True, forces the returned matrix to have the same dimension as
        the RegionMask's mask regardless of the 'resolutionDiv' argument

    applyMask : bool; optional
        When True, the RegionMask's mask will be applied to the outputData
        as described by RegionMask.applyMask

    noData : numeric
        The noData value to use when applying the mask

    preBufferSimplification: numeric
        If given, then geometries will be simplified (using ogr.Geometry.Simplify)
        using the specified value before being buffered
        - Using this can drastically decrease the time it takes to perform the
        buffering procedure, but can decrease accuracy if it is too high

    multiProcess: boolean, optional
        If True, multiple parallel processes will be spawned within the function to
        improve RAM efficiency, else it will fall back on linear execution. By default True.

    kwargs -- Passed on to RegionMask.rasterize()
        * Most notably: 'allTouched'

    Returns
    -------
    numpy.ndarray
    """

    def _indicateFeatures(
        source: load_vector_input,
        where: str | None = None,
        buffer: numeric | None = None,
        bufferMethod: Literal["geom", "area", "contour"] = "geom",
        resolutionDiv: numeric = 1,
        forceMaskShape: bool = False,
        applyMask: bool = True,
        noData: numeric = 0,
        preBufferSimplification: numeric | None = None,
        resultsCollector: dict | multiprocessing.managers.DictProxy | None = None,
        **kwargs,
    ):
        """The core method for _indicateFeatures, can be called either directly or via multiprocessing."""
        assert bufferMethod in ["geom", "area", "contour"]
        # Ensure path to dataSet exists
        source = VECTOR.loadVector(source)  # small ram increase

        # Do we need to buffer?
        if buffer == 0:
            buffer = None
        if buffer is not None and bufferMethod == "geom":

            def doBuffer(ftr):
                if preBufferSimplification is not None:
                    geom = ftr.geom.Simplify(preBufferSimplification)
                else:
                    geom = ftr.geom
                return {"geom": geom.Buffer(float(buffer))}

            source = self.mutateVector(
                source,
                where=where,
                processor=doBuffer,
                matchContext=True,
                keepAttributes=False,
                _slim=True,
                **kwargs,
            )  # big ram increase, much bigger than actual file size. won't bereleased either!!!
            print(f"Memory usage during calc:", str(usage()), "MB")

            where = None  # Set where to None since the filtering has already been done

            if source is None:  # this happens when the returned dataset is empty
                resultsCollector["indications"] = self._returnBlank(
                    resolutionDiv=resolutionDiv,
                    forceMaskShape=forceMaskShape,
                    applyMask=applyMask,
                    noData=noData,
                )
                return

        # Do rasterize
        final = self.rasterize(
            source,
            dtype="float32",
            value=1,
            where=where,
            resolutionDiv=resolutionDiv,
            applyMask=False,
            noData=noData,
        )  # no ram incrreas. final is a "small" np.array
        # Check for results
        if not (final > 0).any():
            # no results were found
            resultsCollector["indications"] = self._returnBlank(
                resolutionDiv=resolutionDiv,
                forceMaskShape=forceMaskShape,
                applyMask=applyMask,
                noData=noData,
            )
            return

        # maybe we want to do the other buffer method
        if buffer is not None and (bufferMethod == "area" or bufferMethod == "contour"):
            if bufferMethod == "area":
                geoms = self.polygonizeMask(final > 0.5, flat=False)
            elif bufferMethod == "contour":
                geoms = self.contoursFromMask(final)

            if preBufferSimplification is not None:
                geoms = [g.Simplify(preBufferSimplification) for g in geoms]

            if len(geoms) > 0:
                geoms = [g.Buffer(float(buffer)) for g in geoms]
                dataSet = VECTOR.createVector(geoms)
                final = self.rasterize(
                    dataSet,
                    dtype="float32",
                    bands=[1],
                    burnValues=[1],
                    resolutionDiv=resolutionDiv,
                    applyMask=False,
                    noData=noData,
                )
            else:
                # no results were found
                resultsCollector["indications"] = self._returnBlank(
                    resolutionDiv=resolutionDiv,
                    forceMaskShape=forceMaskShape,
                    applyMask=applyMask,
                    noData=noData,
                )
                return

        # Make sure we have the mask's shape
        if forceMaskShape:
            if resolutionDiv > 1:
                final = UTIL.scaleMatrix(final, -1 * resolutionDiv)

        # Apply mask?
        if applyMask:
            final = self.applyMask(final, noData)

        # Write results into the resultsCollector dict and exit
        resultsCollector["indications"] = final

    ####################
    # EXECUTE FUNCTION #
    ####################

    if multiProcess:
        # multiprocessing is requested, try to execute
        with multiprocessing.Manager() as manager:
            # initialize a dict to combine the results from parallel processes
            resultsCollector = manager.dict()
            p = multiprocessing.Process(
                target=_indicateFeatures,
                args=(source,),
                kwargs={
                    **{
                        "where": where,
                        "buffer": buffer,
                        "bufferMethod": bufferMethod,
                        "resolutionDiv": resolutionDiv,
                        "forceMaskShape": forceMaskShape,
                        "applyMask": applyMask,
                        "noData": noData,
                        "preBufferSimplification": preBufferSimplification,
                        "resultsCollector": resultsCollector,
                    },
                    **kwargs,
                },
            )
            p.start()
            p.join()

            if "indications" not in resultsCollector.keys():
                error_msg = f"'indications' key not in resultsCollector dict after feature indication."
                # print error statement before raising Error to show it even in try/except loop
                print(error_msg, flush=True)
                raise KeyError(error_msg)

            # extract results from multiprocessing manager
            result = resultsCollector["indications"]

    else:
        # set up a clean resultsCollector in case of multiprocess False or partial results from a multiprocessing failed half-way
        resultsCollector = dict()  # manager.dict()
        # run _indicateFeatures in a single process
        _indicateFeatures(
            source=source,
            where=where,
            buffer=buffer,
            bufferMethod=bufferMethod,
            resolutionDiv=resolutionDiv,
            forceMaskShape=forceMaskShape,
            applyMask=applyMask,
            noData=noData,
            preBufferSimplification=preBufferSimplification,
            resultsCollector=resultsCollector,
            **kwargs,
        )

        # write into return variable
        result = resultsCollector["indications"]

    return result

indicateGeoms

indicateGeoms(geom, **kwargs)

Convenience wrapper to indicate values found within a geometry (or a list of geometries).

• Simply creates a new vector source from the given geometry and then calls RegionMask.indicateFeatures
• All keywords are passed on to RegionMask.indicateFeatures

Source code in geokit/core/regionmask.py

def indicateGeoms(self, geom, **kwargs):
    """
    Convenience wrapper to indicate values found within a geometry (or a
    list of geometries).

    * Simply creates a new vector source from the given geometry and then
      calls RegionMask.indicateFeatures
    * All keywords are passed on to RegionMask.indicateFeatures
    """
    # Make a vector dataset
    ds = UTIL.quickVector(geom)

    # Indicate features
    return self.indicateFeatures(ds, **kwargs)

indicateValueToGeoms

indicateValueToGeoms(
    source, value, contours=False, transformGeoms=True
)

TODO: UPDATE ME.

Source code in geokit/core/regionmask.py

def indicateValueToGeoms(self, source, value, contours=False, transformGeoms=True):
    """TODO: UPDATE ME."""
    # Unpack value
    if isinstance(value, tuple):
        valueMin, valueMax = value
        valueEquals = None
    else:
        valueMin, valueMax = None, None
        valueEquals = value

    # make processor
    def processor(data):
        # Do processing
        if not valueEquals is None:
            output = data == valueEquals
        else:
            output = np.ones(data.shape, dtype="bool")

            if not valueMin is None:
                np.logical_and(data >= valueMin, output, output)
            if not valueMax is None:
                np.logical_and(data <= valueMax, output, output)

        # Done!
        return output

    # Do processing
    newDS = self.extent.mutateRaster(source, processor=processor, dtype="bool", matchContext=False)

    # Make contours
    if contours:
        geomDF = self.extent.contoursFromRaster(newDS, [1], transformGeoms=False)
        geoms = geomDF[geomDF.ID == 1].geom
    else:
        geomDF = RASTER.polygonizeRaster(newDS)
        geoms = geomDF[geomDF.value == 1].geom

    if len(geoms) == 0:
        return []
    else:
        geoms = list(geoms)

    if transformGeoms:
        geoms = GEOM.transform(geoms, toSRS=self.srs)

    return geoms

indicateValues

indicateValues(
    source: load_raster_input,
    value: numeric | tuple | list | str,
    buffer: numeric | None = None,
    resolutionDiv: int = 1,
    forceMaskShape: bool = False,
    applyMask: bool = True,
    noData: numeric | None = None,
    resampleAlg: Literal[
        "near",
        "bilinear",
        "cubic",
        "average",
        "mode",
        "max",
        "min",
    ] = "bilinear",
    bufferMethod: Literal["area", "contour"] = "area",
    preBufferSimplification=None,
    warpDType: geokit_c_data_types_literal | None = None,
    prunePatchSize: numeric = 0,
    threshold: numeric = 0.5,
    multiProcess: bool = True,
    **kwargs,
)

Indicates those pixels in the RegionMask which correspond to a particular value, or range of values, from a given raster datasource.

Returns a matrix matching the RegionMask's mask dimensions wherein 0 means the pixels is not included in the indicated set, and 1 meaning the pixel is included in the indicated set. Intermediate values are also possible. This results from a scenario when the datasource's resolution does not line up perfectly with the RegionMask's resolution and, as a result, a RegionMask pixel overlaps multiple datasource pixels which are not all indicated (or not-indicated).

• Value processing is performed BEFORE a warp takes place
• Output from the warp is clipped to values between 0 and 1
• If a boolean matrix is desired of the result, use "result > 0.5"

Parameters:

• source

  (str or Dataset) –

  The raster datasource to indicate from

• value

  ((numeric, tuple, iterable or str)) –

  The value, range, or set of values to indicate on * If float : The exact value to accept * If tuple : The inclusive range to accept. Given as (low,high) - Assumes exactly 2 values are present - If either value is "None", then the range is assumed to be unbounded on that side * If any other iterable : The list of exact values to accept * If str : The formatted set of elements to accept - Each element in the set is separated by a "," - Each element must be either a singular numeric value, or a range - A range element begins with either "[" or "(", and ends with either "]" or ")" and should have an '-' in between - "[" and "]" imply inclusivity - "(" and ")" imply exclusivity - Numbers on either side can be omitted, impling no limit on that side - Examples: - "[1-5]" -> Indicate values from 1 up to 5, inclusively - "[1-5)" -> Indicate values from 1 up to 5, but not including 5 - "(1-]" -> Indicate values above 1 (but not including 1) up to infinity - "[-5]" -> Indicate values from negative infinity up to and including 5 - "[-]" -> Indicate values from negative infinity to positive infinity (dont do this..) - All whitespaces will be ignored (so feel free to use them as you wish) - Example: - "[-2),[5-7),12,(22-26],29,33,[40-]" will indicate all of the following: - Everything below 2, but not including 2 - Values between 5 up to 7, but not including 7 - 12 - Values above 22 up to and including 26 - 29 - 33 - Everything above 40, including 40

• buffer

  (float; optional, default: None ) –

  A buffer region to add around the indicated pixels * Units are in the RegionMask's srs * The buffering occurs AFTER the indication and warping step and so it may not represent the original dataset exactly - Buffering can be made more accurate by increasing the 'resolutionDiv' input

• resolutionDiv

  (int, default: 1 ) –

  The factor by which to divide the RegionMask's native resolution * This is useful if you need to represent very fine details

• resampleAlg

  (str; optional, default: 'bilinear' ) –

  The resampling algorithm to use when warping values * Options are: 'near', 'bilinear', 'cubic', 'average', 'mode', 'max', 'min' * Knowing which option to use can have significant impacts! When indicating from a low resolution raster (relative to the region mask), then it is best to use one of 'near', 'bilinear', or 'cubic'. However, when indicating from a high resolution raster file (again, relative to the region mask) then one of 'average', 'mode', 'max', or 'min' is likely better.

• warpDType

  (str or None; optional, default: None ) –

  If given, this controls the raster datatype of the warped indication matrix. If not given, then a default datatype is assumed based off resampleAlg:

  reampleAlg : assumed dtype

  ---

     'near' : 'uint8'
  

  'bilinear' : 'float32' 'cubic' : 'float32' 'average' : 'float32' 'mode' : 'uint8' 'max' : 'uint8' 'min' : 'uint8'

• forceMaskShape

  (bool, default: False ) –

  If True, forces the returned matrix to have the same dimension as the RegionMask's mask regardless of the 'resolutionDiv' argument

• applyMask

  (bool, default: True ) –

  When True, the RegionMask's mask will be applied to the outputData as described by RegionMask.applyMask

• noData

  (numeric, default: None ) –

  The noData value to use when applying the mask

• geomsFromContours

  –

  If True, then geometries will be constructed from the function geokit.RegionMask.contoursFromMatrix, as opposed to using geokit.RegionMask.polygonizeMask. - This will result in simpler geometries which are easier to grow, but which do not strictly follow the edges of the indicated pixels

• bufferMethod

  (str; optional, default: 'area' ) –

  An indicator determining the method to use when buffereing * Options are: 'area' and 'contour' * If 'area', the function will first rasterize the raw geometries and will then apply the buffer to the indicated pixels - Uses geokit.RegionMask.polygonizeMask - This is the safer option although is not as accurate as the 'geom' option since it does not capture the exact edges of the geometries - This method can be made more accurate by increasing the 'resolutionDiv' input * If 'contour', the function will still rasterize the raw geometries, but will then create geometries via mask contours (not the explicit pixel edges) - Uses geokit.RegionMask.contoursFromMatrix - This option will recreate geometries which are more similar to the original geometries compared to the 'area' method - This method can be made more accurate by increasing the 'resolutionDiv' input

• preBufferSimplification

  –

  If given, then geometries will be simplified (using ogr.Geometry.Simplify) using the specified value before being buffered - Using this can drastically decrease the time it takes to perform the buffering procedure, but can decrease accuracy if it is too high

• prunePatchSize

  (numeric, default: 0 ) –

  If given, then isolated non-indicated patches below the given size will be removed. The given value corresponds to the minimum area in the unit of the regionmask SRS that will not be removed. Defaults to 0, i.e. no patches will be removed. Note: This is applied to the geoms after buffer application and can deviate from the patch size after final rasterization.

• threshold

  (numeric, default: 0.5 ) –

  The cell value ABOVE which cells count as positively indicated, relevant for partial overlaps with buffer method 'area'. Defaults to 0.5.

• multiProcess

  (bool, default: True ) –

  If True, multiple parallel processes will be spawned within the function to improve RAM efficiency, else it will fall back on linear execution. By default True. Only works on Linux and will be deactivated on Windows and Mac.

• kwargs

  –
  • Most notably: 'resampleAlg'

Returns:

• ndarray –

Source code in geokit/core/regionmask.py

def indicateValues(
    self,
    source: load_raster_input,
    value: numeric | tuple | list | str,
    buffer: numeric | None = None,
    resolutionDiv: int = 1,
    forceMaskShape: bool = False,
    applyMask: bool = True,
    noData: numeric | None = None,
    resampleAlg: Literal["near", "bilinear", "cubic", "average", "mode", "max", "min"] = "bilinear",
    bufferMethod: Literal["area", "contour"] = "area",
    preBufferSimplification=None,
    warpDType: geokit_c_data_types_literal | None = None,
    prunePatchSize: numeric = 0,
    threshold: numeric = 0.5,
    multiProcess: bool = True,
    **kwargs,
):
    """
    Indicates those pixels in the RegionMask which correspond to a particular
    value, or range of values, from a given raster datasource.

    Returns a matrix matching the RegionMask's mask dimensions wherein 0 means
    the pixels is not included in the indicated set, and 1 meaning the pixel
    is included in the indicated set. Intermediate values are also possible.
    This results from a scenario when the datasource's resolution does not
    line up perfectly with the RegionMask's resolution and, as a result, a
    RegionMask pixel overlaps multiple datasource pixels which are not all
    indicated (or not-indicated).

    * Value processing is performed BEFORE a warp takes place
    * Output from the warp is clipped to values between 0 and 1
    * If a boolean matrix is desired of the result, use "result > 0.5"

    Parameters
    ----------
    source : str or gdal.Dataset
        The raster datasource to indicate from

    value : numeric, tuple, iterable or str
        The value, range, or set of values to indicate on
        * If float : The exact value to accept
        * If tuple : The inclusive range to accept. Given as (low,high)
          - Assumes exactly 2 values are present
          - If either value is "None", then the range is assumed to be unbounded on that side
        * If any other iterable : The list of exact values to accept
        * If str : The formatted set of elements to accept
          - Each element in the set is separated by a ","
          - Each element must be either a singular numeric value, or a range
          - A range element begins with either "[" or "(", and ends with either "]" or ")"
            and should have an '-' in between
            - "[" and "]" imply inclusivity
            - "(" and ")" imply exclusivity
            - Numbers on either side can be omitted, impling no limit on that side
            - Examples:
              - "[1-5]" -> Indicate values from 1 up to 5, inclusively
              - "[1-5)" -> Indicate values from 1 up to 5, but not including 5
              - "(1-]"  -> Indicate values above 1 (but not including 1) up to infinity
              - "[-5]"  -> Indicate values from negative infinity up to and including 5
              - "[-]"   -> Indicate values from negative infinity to positive infinity (dont do this..)
          - All whitespaces will be ignored (so feel free to use them as you wish)
          - Example:
            - "[-2),[5-7),12,(22-26],29,33,[40-]" will indicate all of the following:
              - Everything below 2, but not including 2
              - Values between 5 up to 7, but not including 7
              - 12
              - Values above 22 up to and including 26
              - 29
              - 33
              - Everything above 40, including 40

    buffer : float; optional
        A buffer region to add around the indicated pixels
        * Units are in the RegionMask's srs
        * The buffering occurs AFTER the indication and warping step and
          so it may not represent the original dataset exactly
          - Buffering can be made more accurate by increasing the
            'resolutionDiv' input

    resolutionDiv : int
        The factor by which to divide the RegionMask's native resolution
        * This is useful if you need to represent very fine details

    resampleAlg : str; optional
        The resampling algorithm to use when warping values
        * Options are: 'near', 'bilinear', 'cubic', 'average', 'mode', 'max', 'min'
        * Knowing which option to use can have significant impacts!
            When indicating from a low resolution raster (relative to the region mask),
            then it is best to use one of 'near', 'bilinear', or 'cubic'. However,
            when indicating from a high resolution raster file (again, relative to the region
            mask) then one of 'average', 'mode', 'max', or 'min' is likely better.

    warpDType : str or None; optional
        If given, this controls the raster datatype of the warped indication matrix.
        If not given, then a default datatype is assumed based off `resampleAlg`:

           reampleAlg : assumed dtype
           --------------------------
               'near' : 'uint8'
           'bilinear' : 'float32'
              'cubic' : 'float32'
            'average' : 'float32'
               'mode' : 'uint8'
                'max' : 'uint8'
                'min' : 'uint8'

    forceMaskShape : bool
        If True, forces the returned matrix to have the same dimension as
        the RegionMask's mask regardless of the 'resolutionDiv' argument

    applyMask : bool
        When True, the RegionMask's mask will be applied to the outputData
        as described by RegionMask.applyMask

    noData : numeric
        The noData value to use when applying the mask

    geomsFromContours: bool
        If True, then geometries will be constructed from the function
        geokit.RegionMask.contoursFromMatrix, as opposed to using
        geokit.RegionMask.polygonizeMask.
        - This will result in simpler geometries which are easier to grow,
          but which do not strictly follow the edges of the indicated pixels

    bufferMethod : str; optional
        An indicator determining the method to use when buffereing
        * Options are: 'area' and 'contour'
        * If 'area', the function will first rasterize the raw geometries and
          will then apply the buffer to the indicated pixels
          - Uses geokit.RegionMask.polygonizeMask
          - This is the safer option although is not as accurate as the 'geom'
            option since it does not capture the exact edges of the geometries
          - This method can be made more accurate by increasing the
            'resolutionDiv' input
        * If 'contour', the function will still rasterize the raw geometries,
          but will then create geometries via mask contours (not the explicit
          pixel edges)
          - Uses geokit.RegionMask.contoursFromMatrix
          - This option will recreate geometries which are more similar to the
            original geometries compared to the 'area' method
          - This method can be made more accurate by increasing the
            'resolutionDiv' input

    preBufferSimplification: numeric; optional
        If given, then geometries will be simplified (using ogr.Geometry.Simplify)
        using the specified value before being buffered
        - Using this can drastically decrease the time it takes to perform the
          buffering procedure, but can decrease accuracy if it is too high

    prunePatchSize: numeric; optional
        If given, then isolated non-indicated patches below the given size
        will be removed. The given value corresponds to the minimum area in
        the unit of the regionmask SRS that will not be removed. Defaults
        to 0, i.e. no patches will be removed.
        Note: This is applied to the geoms after buffer application and can
        deviate from the patch size after final rasterization.

    threshold: float; optional
        The cell value ABOVE which cells count as positively indicated,
        relevant for partial overlaps with buffer method 'area'. Defaults to 0.5.

    multiProcess: boolean, optional
        If True, multiple parallel processes will be spawned within the function to
        improve RAM efficiency, else it will fall back on linear execution. By default True.
        Only works on Linux and will be deactivated on Windows and Mac.

    kwargs -- Passed on to RegionMask.warp()
        * Most notably: 'resampleAlg'

    Returns
    -------
    numpy.ndarray
    """
    multi_processing_warning_message = (
        "Multiprocessing has been set to 'False' because it is not available for Windows or Mac."
        " To deactivate this warning, please set the multiProcess variable to False. On Windows and "
        "Mac, new processes must be spawned, which requires the serialisation of the method to be "
        "executed via multiprocessing. However, Geokit contains objects that cannot be serialised by Pickle. "
        "On Linux, however, new processes are inherited and no serialisation is required."
    )
    if platform == "linux" or platform == "linux2":
        pass
    elif platform == "darwin" and multiProcess is True:
        multiProcess = False
        warn(multi_processing_warning_message)
    elif platform == "win32" and multiProcess is True:
        multiProcess = False
        warn(multi_processing_warning_message)
    elif platform == "darwin" and multiProcess is False:
        pass
    elif platform == "win32" and multiProcess is False:
        pass
    else:
        multiProcess = False
        warn(multi_processing_warning_message)

    def _indicateValues(
        source: load_raster_input,
        value,
        buffer=None,
        resolutionDiv=1,
        forceMaskShape=False,
        applyMask=True,
        noData=None,
        resampleAlg: Literal["near", "bilinear", "cubic", "average", "mode", "max", "min"] = "bilinear",
        bufferMethod="area",
        preBufferSimplification=None,
        warpDType=None,
        prunePatchSize=0,
        threshold=0.5,
        resultsCollector=None,
        **kwargs,
    ):
        """The core method for _indicateValues, can be called either directly or via multiprocessing."""
        assert bufferMethod in ["area", "contour"]
        assert isinstance(prunePatchSize, int)

        # format value input
        if isinstance(value, str):
            pass
        elif isinstance(value, tuple):  # Assume a range in implied
            value = "[{}-{}]".format(
                "" if value[0] is None else value[0],
                "" if value[1] is None else value[1],
            )
        else:
            try:  # Try treating value as an iterable
                _value = ""
                for v in value:
                    _value += "{},".format(v)
                value = _value[:-1]
            except TypeError:  # Value should be just a number
                value = str(value)

        # make processor
        def processor(data):
            # Find nan values, maybe
            if not noData is None:
                nodat = np.isnan(data)

            # Indicate value elements
            output = np.zeros(data.shape, dtype="bool")
            value_re = re.compile(
                r"(?P<range>(?P<open>[\[\(])(?P<low>[-+]?(\d*\.\d+|\d+\.?))?-(?P<high>[-+]?(\d*\.\d+|\d+\.?))?(?P<close>[\]\)]))|(?P<value>[-+]?(\d*\.\d+|\d+\.?))"
            )
            for element in value.split(","):
                element = element.replace(" ", "")
                if element == "":
                    continue

                m = value_re.match(element)
                if m is None or (m["value"] is None and m["range"] is None):
                    raise RuntimeError('The element "{}" does not match an expected format'.format(element))

                if m["value"] is not None:
                    update_sel = data == float(m["value"])

                else:  # We are dealing with a range
                    update_sel = np.ones(data.shape, dtype="bool")

                    if m["low"] is not None and m["open"] == "[":
                        np.logical_and(data >= float(m["low"]), update_sel, update_sel)

                    elif m["low"] is not None and m["open"] == "(":
                        np.logical_and(data > float(m["low"]), update_sel, update_sel)

                    if m["high"] is not None and m["close"] == "]":
                        np.logical_and(data <= float(m["high"]), update_sel, update_sel)

                    elif m["high"] is not None and m["close"] == ")":
                        np.logical_and(data < float(m["high"]), update_sel, update_sel)

                np.logical_or(update_sel, output, output)

            # Fill nan values, maybe
            if noData is not None:
                output[nodat] = noData

            # Done!
            return output

        # Do processing
        newDS = self.extent.mutateRaster(
            source,
            processor=processor,
            dtype="uint8",
            noData=noData,
            matchContext=False,
        )
        print(f"Memory usage during calc:", str(usage()), "MB")

        # Warp onto region
        if warpDType is None:
            if resampleAlg in ["bilinear", "cubic", "average"]:
                warpDType = "float32"
            elif resampleAlg in ["near", "mode", "max", "min"]:
                warpDType = "uint8"
            else:
                warpDType = "float32"

        final = self.warp(
            newDS,
            dtype=warpDType,
            resolutionDiv=resolutionDiv,
            resampleAlg=resampleAlg,
            applyMask=False,
            noData=noData,
            returnMatrix=True,
            **kwargs,
        )

        # Check for results
        if not (final > 0).any():
            # no results were found
            resultsCollector["indications"] = self._returnBlank(
                resolutionDiv=resolutionDiv,
                forceMaskShape=forceMaskShape,
                applyMask=applyMask,
                noData=noData,
            )
            return

        # Apply a buffer if requested
        if not buffer is None:
            if bufferMethod == "contour":
                geoms = self.contoursFromMask(final)
            elif bufferMethod == "area":
                # Check for results
                if not (final > threshold).any():
                    # no cells > threshold exist that would be buffered, so return empty
                    resultsCollector["indications"] = self._returnBlank(
                        resolutionDiv=resolutionDiv,
                        forceMaskShape=forceMaskShape,
                        applyMask=applyMask,
                        noData=noData,
                    )
                    return
                geoms = self.polygonizeMask(final > threshold, flat=False)

            if preBufferSimplification is not None:
                geoms = [g.Simplify(preBufferSimplification) for g in geoms]

            if len(geoms) > 0:
                geoms = [g.Buffer(float(buffer)) for g in geoms]
                if not prunePatchSize == 0:
                    # create Union of all geoms, this will merge overlapping
                    # TODO speed up via cascaded union or better sieve raster or even inverted raster (approximate exclusions based on geoms directly)!
                    union = geoms[0]
                    if len(geoms) > 1:
                        for g in geoms:
                            union = union.Union(g)

                    # create a new empty list of geoms that will be filled with manipulated shapes partly without holes
                    geoms = []
                    # if regionmask extent contains only a single polygon, make it a single-entry list to make it iterable
                    if union.GetGeometryName() == "POLYGON":
                        union = [union]
                    # iterate over union first to separate independent polygons
                    for polygon in union:
                        # check if the independent iteration polygon has more than one geometry, if so has holes/inner rings
                        NoIndivPol = polygon.GetGeometryCount()
                        if NoIndivPol > 1:
                            print("inner polygons found")
                            # create a new base polygon containing the outer ring
                            newpolygon = ogr.Geometry(ogr.wkbPolygon)
                            newpolygon.AddGeometry(polygon.GetGeometryRef(0))
                            # then iterate through the other geometries = inner rings
                            for i in range(1, NoIndivPol):
                                # only if the area of the inner ring exceeds the prunePatchSize criterion, add back as inner ring
                                # else do nothing i.e. drop the inner ring
                                if polygon.GetGeometryRef(i).Area() >= prunePatchSize:
                                    # TODO this part might be faster with polygon.RemoveGeometry() without creating a 'newpolygon' above
                                    newpolygon.AddGeometry(polygon.GetGeometryRef(i))
                                else:
                                    print("a innerring was dropped!")
                            # add the new polygon with possibly less holes back to geoms
                            geoms.append(newpolygon)

                        # if the polygon has only one geometry, it has no inner ring and can be added to geom list as such
                        else:
                            geoms.append(polygon)

                areaDS = VECTOR.createVector(geoms)
                final = self.rasterize(
                    areaDS,
                    dtype="float32",
                    bands=[1],
                    burnValues=[1],
                    resolutionDiv=resolutionDiv,
                    applyMask=False,
                    noData=noData,
                )
            else:
                # no results were found
                resultsCollector["indications"] = self._returnBlank(
                    resolutionDiv=resolutionDiv,
                    forceMaskShape=forceMaskShape,
                    applyMask=applyMask,
                    noData=noData,
                )
                return

        # apply a threshold in case of funky warping issues
        final[final > 1.0] = 1
        final[final < 0.0] = 0

        # Make sure we have the mask's shape

        if forceMaskShape:
            if resolutionDiv > 1:
                final = UTIL.scaleMatrix(final, -1 * resolutionDiv)

        # Apply mask?
        if applyMask:
            final = self.applyMask(final, noData)

        # Write results into the resultsCollector dict and exit
        resultsCollector["indications"] = final
        return

    ####################
    # EXECUTE FUNCTION #
    ####################

    # try multiprocessing or fall back to linear processing

    if multiProcess:
        # multiprocessing is requested, try to execute
        with multiprocessing.Manager() as manager:
            # initialize a dict to combine the results from parallel processes
            resultsCollector = manager.dict()
            p = multiprocessing.Process(
                target=_indicateValues,
                args=(source,),
                kwargs={
                    **{
                        "value": value,
                        "buffer": buffer,
                        "resolutionDiv": resolutionDiv,
                        "forceMaskShape": forceMaskShape,
                        "applyMask": applyMask,
                        "noData": noData,
                        "resampleAlg": resampleAlg,
                        "bufferMethod": bufferMethod,
                        "preBufferSimplification": preBufferSimplification,
                        "warpDType": warpDType,
                        "prunePatchSize": prunePatchSize,
                        "threshold": threshold,
                        "resultsCollector": resultsCollector,
                    },
                    **kwargs,
                },
            )
            p.start()
            p.join()

            if not "indications" in resultsCollector.keys():
                error_msg = f"'indications' key not in resultsCollector dict after value indication."
                # print error statement before raising Error to show it even in try/except loop
                print(error_msg, flush=True)
                raise KeyError(error_msg)

            # extract results from multiprocessing manager
            result = resultsCollector["indications"]

    else:
        # set up a clean resultsCollector in case of multiprocess False or partial results from a multiprocessing failed half-way
        resultsCollector = dict()  # manager.dict()
        # run _indicateValues in a single process
        _indicateValues(
            source=source,
            value=value,
            buffer=buffer,
            resolutionDiv=resolutionDiv,
            forceMaskShape=forceMaskShape,
            applyMask=applyMask,
            noData=noData,
            resampleAlg=resampleAlg,
            bufferMethod=bufferMethod,
            preBufferSimplification=preBufferSimplification,
            warpDType=warpDType,
            prunePatchSize=prunePatchSize,
            threshold=threshold,
            resultsCollector=resultsCollector,
            **kwargs,
        )

        # write into return variable
        result = resultsCollector["indications"]

    return result

load staticmethod

load(region, start_raster=None, **kwargs) -> RegionMask

Tries to initialize and return a RegionMask in the most appropriate way.

Note:

If 'region' input is... * Already a RegionMask, simply return it * A file path, use RegionMask.fromVector * An OGR Geometry object, assume is it to be loaded by RegionMask.fromGeom * A NumPy array, assume is it to be loaded by RegionMask.fromMask - An 'extent' input must also be given

Parameters:

• region

  (Can be RegionMask, str, ogr.Geometry, numpy.ndarray) –

  The shape defining the region over which to build the RegionMask * See the note above

• where

  (str, int; optional) –

  If string -> An SQL-like where statement to apply to the source If int -> The feature's ID within the dataset * Feature attribute name do not need quotes * String values should be wrapped in 'single quotes' Example: If the source vector has a string attribute called "ISO" and a integer attribute called "POP", you could use....

  where = "ISO='DEU' AND POP>1000"
  

• geom

  (ogr.Geometry; optional) –

  The geometry to search with * All features are extracted which touch this geometry

• pixelRes

  (float or tuple) –

  The RegionMask's native pixel resolution(s) * If float : A pixel size to apply to both the X and Y dimension * If (float float) : An X-dimension and Y-dimension pixel size

• srs

  (Anything acceptable to geokit.srs.loadSRS()) –

  The srs context of the generated RegionMask object * This srs is superseded by the srs in an explicitly defined extent

• extent

  (Extent object) –

  The geospatial context of the generated region mask * The extent must fit the given pixel sizes * If not specified, the entire extent of the vector file is assumed

• padExtent

  (float; optional) –

  An amount by which to pad the extent before generating the RegionMask

Source code in geokit/core/regionmask.py

@staticmethod
def load(region, start_raster=None, **kwargs) -> "RegionMask":
    """Tries to initialize and return a RegionMask in the most appropriate way.

    Note:
    -----
    If 'region' input is...
        * Already a RegionMask, simply return it
        * A file path, use RegionMask.fromVector
        * An OGR Geometry object, assume is it to be loaded by RegionMask.fromGeom
        * A NumPy array, assume is it to be loaded by RegionMask.fromMask
            - An 'extent' input must also be given

    Parameters
    ----------
    region : Can be RegionMask, str, ogr.Geometry, numpy.ndarray
        The shape  defining the region over which to build the RegionMask
        * See the note above

    where : str, int; optional
        If string -> An SQL-like where statement to apply to the source
        If int -> The feature's ID within the dataset
        * Feature attribute name do not need quotes
        * String values should be wrapped in 'single quotes'
        Example: If the source vector has a string attribute called "ISO" and
                 a integer attribute called "POP", you could use....

            where = "ISO='DEU' AND POP>1000"

    geom : ogr.Geometry; optional
        The geometry to search with
        * All features are extracted which touch this geometry

    pixelRes : float or tuple
        The RegionMask's native pixel resolution(s)
        * If float : A pixel size to apply to both the X and Y dimension
        * If (float float) : An X-dimension and Y-dimension pixel size

    srs : Anything acceptable to geokit.srs.loadSRS()
        The srs context of the generated RegionMask object
        * This srs is superseded by the srs in an explicitly defined extent

    extent : Extent object
        The geospatial context of the generated region mask
        * The extent must fit the given pixel sizes
        * If not specified, the entire extent of the vector file is assumed

    padExtent : float; optional
        An amount by which to pad the extent before generating the RegionMask
    """
    if isinstance(region, RegionMask):
        return region
    elif isinstance(region, str):
        return RegionMask.fromVector(region, start_raster=start_raster, **kwargs)
    elif isinstance(region, ogr.Geometry):
        return RegionMask.fromGeom(region, start_raster=start_raster, **kwargs)
    elif isinstance(region, np.ndarray):
        return RegionMask.fromMask(region, **kwargs)
    else:
        raise GeoKitRegionMaskError("Could not understand region input")

mutateRaster

mutateRaster(
    source: load_raster_input,
    matchContext: bool = True,
    warpArgs: dict | None = None,
    applyMask: bool = True,
    processor: Callable | None = None,
    resampleAlg: Literal[
        "near", "bilinear", "cubic", "average"
    ] = "bilinear",
    **mutateArgs,
)

Convenience wrapper for geokit.vector.mutateRaster which automatically sets 'bounds'. It also warps the raster to the RegionMask's area and srs before mutating.

Note:

If this is called without any arguments except for a source, it serves to clip the raster source around the RegionMask, therefore performing the same function as RegionMask.warp(..., returnMatrix=False)

Parameters:

• source

  (Anything acceptable to geokit.raster.loadRaster()) –

  The source to mutate

• matchContext

  (bool; optional, default: True ) –
  • If True, Warp to the RegionMask's boundaries, srs and pixel size before mutating
  • If False, only warp to the RegionMask's boundaries, but keep its srs and resolution intact

• resampleAlg

  (str; optional, default: 'bilinear' ) –

  The resampling algorithm to use when warping values * Knowing which option to use can have significant impacts! * Options are: 'nearesampleAlg=resampleAlg, r', 'bilinear', 'cubic', 'average'

• warpArgs

  (dict; optional, default: None ) –

  Arguments to apply to the warping step * See geokit.raster.warp()

• processor

  (Callable | None, default: None ) –

  The function performing the mutation of the raster's data * The function will take single argument (a 2D numpy.ndarray) * The function must return a numpy.ndarray of the same size as the input * The return type must also be containable within a Float32 (int and boolean is okay) * See example in geokit.raster.mutateRaster for more info

• applyMask

  (bool; optional, default: True ) –

  When True, the RegionMask's mask will be applied to the outputData as described by RegionMask.applyMask

• **mutateArgs

  –

  All other keyword arguments are passed to geokit.vector.mutateVector

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/regionmask.py

def mutateRaster(
    self,
    source: load_raster_input,
    matchContext: bool = True,
    warpArgs: dict | None = None,
    applyMask: bool = True,
    processor: Callable | None = None,
    resampleAlg: Literal["near", "bilinear", "cubic", "average"] = "bilinear",
    **mutateArgs,
):
    """Convenience wrapper for geokit.vector.mutateRaster which automatically
    sets 'bounds'. It also warps the raster to the RegionMask's area
    and srs before mutating.

    Note:
    -----
    If this is called without any arguments except for a source, it serves
    to clip the raster source around the RegionMask, therefore performing
    the same function as RegionMask.warp(..., returnMatrix=False)

    Parameters
    ----------
    source : Anything acceptable to geokit.raster.loadRaster()
        The source to mutate

    matchContext : bool; optional
        * If True, Warp to the RegionMask's boundaries, srs and pixel size
          before mutating
        * If False, only warp to the RegionMask's boundaries, but keep its
          srs and resolution intact

    resampleAlg : str; optional
        The resampling algorithm to use when warping values
        * Knowing which option to use can have significant impacts!
        * Options are: 'nearesampleAlg=resampleAlg, r', 'bilinear', 'cubic',
          'average'

    warpArgs : dict; optional
        Arguments to apply to the warping step
        * See geokit.raster.warp()

    processor - function; optional
        The function performing the mutation of the raster's data
        * The function will take single argument (a 2D numpy.ndarray)
        * The function must return a numpy.ndarray of the same size as the input
        * The return type must also be containable within a Float32 (int and
          boolean is okay)
        * See example in geokit.raster.mutateRaster for more info

    applyMask : bool; optional
        When True, the RegionMask's mask will be applied to the outputData
        as described by RegionMask.applyMask

    **mutateArgs:
        All other keyword arguments are passed to geokit.vector.mutateVector

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    output = mutateArgs.pop("output", None)
    if warpArgs is None:
        warpArgs = {}

    # Do the warp and mutation
    if matchContext:
        source = self.warp(
            source,
            returnMatrix=False,
            applyMask=applyMask,
            resampleAlg=resampleAlg,
            **warpArgs,
        )

        if processor is None:
            return source
        else:
            return RASTER.mutateRaster(source, output=output, **mutateArgs)

    else:
        if applyMask:
            if "cutline" in warpArgs:
                raise GeoKitRegionMaskError("Cannot apply both a cutline and the mask during prewarping")
            warpArgs["cutline"] = self.vector

        return self.extent.mutateRaster(
            source,
            matchContext=False,
            warpArgs=warpArgs,
            resampleAlg=resampleAlg,
            **mutateArgs,
        )

mutateVector

mutateVector(
    source, matchContext=False, regionPad=0, **kwargs
)

Convenience wrapper for geokit.vector.mutateVector which automatically sets 'srs' and 'geom' inputs to the RegionMask's srs and geometry.

• The RegionMask's geometry is always used to select features within the source. If you need a broader scope, try using the RegionMask's extent's version of this function

Note:

If this is called without any arguments except for a source, it serves to clip the vector source around the RegionMask

Parameters:

• source

  (Anything acceptable to geokit.vector.loadVector()) –

  The source to clip

• matchContext

  (bool; optional, default: False ) –
  • If True, transforms all geometries to the RegionMask's srs before mutating
  • If False, only selects the geometries which touch the RegionMask

• regionPad

  –

  Will buffer the regionmask geometry by a value (in the unit of the regionmask srs) before mutating the source vector onto it. Defaults to 0.

• **kwargs

  –

  All other keyword arguments are passed to geokit.vector.mutateVector

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/regionmask.py

def mutateVector(self, source, matchContext=False, regionPad=0, **kwargs):
    """Convenience wrapper for geokit.vector.mutateVector which automatically
    sets 'srs' and 'geom' inputs to the RegionMask's srs and geometry.

    * The RegionMask's geometry is always used to select features within the
    source. If you need a broader scope, try using the RegionMask's extent's
    version of this function

    Note:
    -----
    If this is called without any arguments except for a source, it serves
    to clip the vector source around the RegionMask

    Parameters
    ----------
    source : Anything acceptable to geokit.vector.loadVector()
        The source to clip

    matchContext : bool; optional
        * If True, transforms all geometries to the RegionMask's srs before
          mutating
        * If False, only selects the geometries which touch the RegionMask

    regionPad: int, optional
        Will buffer the regionmask geometry by a value (in the unit of the
        regionmask srs) before mutating the source vector onto it. Defaults
        to 0.

    **kwargs:
        All other keyword arguments are passed to geokit.vector.mutateVector

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    # Get the working srs
    if not matchContext:
        vinfo = VECTOR.vectorInfo(source)
        ext = self.extent.castTo(vinfo.srs)
    else:
        ext = self.extent

    # mutate the source
    if regionPad == None:  # cannot calculate with none, so use 0
        regionPad = 0
    return VECTOR.mutateVector(source, srs=ext.srs, geom=self.geometry.Buffer(float(regionPad)), **kwargs)

polygonizeMask

polygonizeMask(mask, flat=True, shrink=True)

Convenience wrapper for geokit.geom.polygonizeMask which automatically sets the 'bounds' and 'srs' inputs. The mask data is assumed to span the RegionMask exactly.

Each True-valued group of pixels will be converted to a geometry

Parameters:

• mask

  (matrix_like) –

  The mask which will be turned into a geometry set * Must be 2 dimensional * Must be boolean type * True values are interpreted as 'in the geometry'

• flat

  (bool, default: True ) –

  If True, flattens the resulting geometries into a single geometry

• shrink

  (bool, default: True ) –

  If True, shrink all geoms by a tiny amount in order to avoid geometry overlapping issues * The total amount shrunk should be very very small * Generally this should be left as True unless it is ABSOLUTELY necessary to maintain the same area

Returns:

• If 'flat' is True: ogr.Geometry –
• else ( [Geometry] ) –

Source code in geokit/core/regionmask.py

def polygonizeMask(self, mask, flat=True, shrink=True):
    """Convenience wrapper for geokit.geom.polygonizeMask which automatically
    sets the 'bounds' and 'srs' inputs. The mask data is assumed to span the
    RegionMask exactly.

    Each True-valued group of pixels will be converted to a geometry

    Parameters
    ----------
    mask : matrix_like
        The mask which will be turned into a geometry set
          * Must be 2 dimensional
          * Must be boolean type
          * True values are interpreted as 'in the geometry'

    flat : bool
        If True, flattens the resulting geometries into a single geometry

    shrink : bool
        If True, shrink all geoms by a tiny amount in order to avoid geometry
        overlapping issues
          * The total amount shrunk should be very very small
          * Generally this should be left as True unless it is ABSOLUTELY
            necessary to maintain the same area

    Returns
    -------
    If 'flat' is True: ogr.Geometry
    else: [ogr.Geometry,  ]
    """
    return GEOM.polygonizeMask(mask, bounds=self.extent.xyXY, srs=self.srs, flat=flat, shrink=shrink)

polygonizeMatrix

polygonizeMatrix(
    matrix, flat=False, shrink=True, _raw=False
)

Convenience wrapper for geokit.geom.polygonizeMatrix which automatically sets the 'bounds' and 'srs' inputs. The matrix data is assumed to span the RegionMask exactly.

Each unique-valued group of pixels will be converted to a geometry

Parameters:

• matrix

  (matrix_like) –

  The matrix which will be turned into a geometry set * Must be 2 dimensional * Must be integer or boolean type

• flat

  (bool, default: False ) –

  If True, flattens the resulting geometries which share a contiguous matrix value into a single geometry object

• shrink

  (bool, default: True ) –

  If True, shrink all geoms by a tiny amount in order to avoid geometry overlapping issues * The total amount shrunk should be very very small * Generally this should be left as True unless it is ABSOLUTELY necessary to maintain the same area

Returns:

• pandas.DataFrame -> With columns: –

  'geom' -> The contiguous-valued geometries 'value' -> The value for each geometry

Source code in geokit/core/regionmask.py

def polygonizeMatrix(self, matrix, flat=False, shrink=True, _raw=False):
    """Convenience wrapper for geokit.geom.polygonizeMatrix which automatically
    sets the 'bounds' and 'srs' inputs. The matrix data is assumed to span the
    RegionMask exactly.

    Each unique-valued group of pixels will be converted to a geometry

    Parameters
    ----------
    matrix : matrix_like
        The matrix which will be turned into a geometry set
          * Must be 2 dimensional
          * Must be integer or boolean type

    flat : bool
        If True, flattens the resulting geometries which share a contiguous matrix
        value into a single geometry object

    shrink : bool
        If True, shrink all geoms by a tiny amount in order to avoid geometry
        overlapping issues
          * The total amount shrunk should be very very small
          * Generally this should be left as True unless it is ABSOLUTELY
            necessary to maintain the same area

    Returns
    -------
    pandas.DataFrame -> With columns:
                            'geom' -> The contiguous-valued geometries
                            'value' -> The value for each geometry
    """
    return GEOM.polygonizeMatrix(
        matrix,
        bounds=self.extent.xyXY,
        srs=self.srs,
        flat=flat,
        shrink=shrink,
        _raw=_raw,
    )

rasterize

rasterize(
    source,
    output=None,
    resolutionDiv=1,
    returnMatrix=True,
    applyMask=True,
    noData=None,
    **kwargs,
)

Convenience wrapper for geokit.vector.rasterize() which automatically sets the 'srs', 'bounds', 'pixelWidth', and 'pixelHeight' inputs.

Note:

When creating an 'in memory' raster vs one which is saved to disk, a slightly different algorithm is used which can sometimes add an extra row of pixels. Be aware of this if you intend to compare value-matricies directly from rasters generated with this function.

Parameters:

• source

  (str) –

  The path to the vector file to load

• output

  (str; optional, default: None ) –

  A path to an output file to write to

• resolutionDiv

  (int; optional, default: 1 ) –

  The factor by which to divide the RegionMask's native resolution * This is useful if you need to represent very fine details

• returnAsMatrix

  (bool; optional) –

  When True, the resulting raster's matrix is return * Should have the same dimensions as the RegionMask's mask matrix

• applyMask

  (bool; optional, default: True ) –

  When True, the RegionMask's mask will be applied to the outputData as described by RegionMask.applyMask

• noData

  (numeric; optional, default: None ) –

  The noData value to use when applying the mask

• **kwargs

  –

  All other keywargs are passed on to geokit.vector.rasterize()

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/regionmask.py

def rasterize(
    self,
    source,
    output=None,
    resolutionDiv=1,
    returnMatrix=True,
    applyMask=True,
    noData=None,
    **kwargs,
):
    """Convenience wrapper for geokit.vector.rasterize() which automatically
    sets the 'srs', 'bounds', 'pixelWidth', and 'pixelHeight' inputs.

    Note:
    -----
    When creating an 'in memory' raster vs one which is saved to disk, a slightly
    different algorithm is used which can sometimes add an extra row of pixels. Be
    aware of this if you intend to compare value-matricies directly from rasters
    generated with this function.

    Parameters
    ----------
    source : str
        The path to the vector file to load

    output : str; optional
        A path to an output file to write to

    resolutionDiv : int; optional
        The factor by which to divide the RegionMask's native resolution
        * This is useful if you need to represent very fine details

    returnAsMatrix : bool; optional
        When True, the resulting raster's matrix is return
        * Should have the same dimensions as the RegionMask's mask matrix

    applyMask : bool; optional
        When True, the RegionMask's mask will be applied to the outputData
        as described by RegionMask.applyMask

    noData : numeric; optional
        The noData value to use when applying the mask

    **kwargs:
        All other keywargs are passed on to geokit.vector.rasterize()

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    pW, pH = self._resolve(resolutionDiv)

    # do rasterization
    if returnMatrix:
        newDS = self.extent.rasterize(
            source=source,
            pixelWidth=pW,
            pixelHeight=pH,
            output=output,
            noData=noData,
            **kwargs,
        )
    else:
        if applyMask:
            if "cutline" in kwargs:
                raise GeoKitRegionMaskError(
                    "Cannot apply both a cutline and the mask when returning the rasterized dataset"
                )
            newDS = self.extent.rasterize(
                source=source,
                pixelWidth=pW,
                pixelHeight=pH,
                cutline=self.vectorPath,
                output=output,
                noData=noData,
                **kwargs,
            )
        else:
            newDS = self.extent.rasterize(
                source=source,
                pixelWidth=pW,
                pixelHeight=pH,
                output=output,
                noData=noData,
                **kwargs,
            )

    if not returnMatrix:
        return newDS

    # Read array
    if newDS is None:
        newDS = output
    final = RASTER.extractMatrix(newDS)

    # clean up
    del newDS

    # Apply mask, maybe
    if applyMask:
        final = self.applyMask(final, noData)

    # Return
    return final

subRegions

subRegions(gridSize, asMaskAndExtent=False)

Generate a number of sub regions on a grid which combine into the total RegionMask area.

Source code in geokit/core/regionmask.py

def subRegions(self, gridSize, asMaskAndExtent=False):
    """Generate a number of sub regions on a grid which combine into the total
    RegionMask area.
    """
    # get useful matrix info
    yN, xN = self.mask.shape
    pixelGridSize = int(gridSize / min(self.pixelWidth, self.pixelHeight))

    # Make grid areas
    count = 0
    for ys in range(0, yN, pixelGridSize):
        yn = min(yN, ys + pixelGridSize)
        yMax = self.extent.yMax - ys * self.pixelHeight
        yMin = self.extent.yMax - yn * self.pixelHeight

        for xs in range(0, xN, pixelGridSize):
            xn = min(xN, xs + pixelGridSize)
            xMin = self.extent.xMin + xs * self.pixelWidth
            xMax = self.extent.xMin + xn * self.pixelWidth

            sectionMask = self.mask[ys:yn, xs:xn]
            if not sectionMask.any():
                continue

            sectionExtent = Extent(xMin, yMin, xMax, yMax, srs=self.srs).fit((self.pixelWidth, self.pixelHeight))

            if asMaskAndExtent:
                yield MaskAndExtent(sectionMask, sectionExtent, count)
            else:
                yield RegionMask.fromMask(sectionExtent, sectionMask, dict(id=count))

            count += 1

subTiles

subTiles(zoom, checkIntersect=True, asGeom=False)

Generates tile Extents at a given zoom level which encompass the invoking Regionmask.

Parameters:

• zoom

  (int) –

  The zoom level of the expected tile source

• checkIntersect

  (bool, default: True ) –

  If True, exclude tiles which do not intersect with the RegionMask's geometry

• asGeom

  (bool, default: False ) –

  If True, returns tuple of ogr.Geometries in stead of (xi,yi,zoom) tuples

Returns:

• Generator of Geometries or (xi,yi,zoom) tuples –

Source code in geokit/core/regionmask.py

def subTiles(self, zoom, checkIntersect=True, asGeom=False):
    """Generates tile Extents at a given zoom level which encompass the invoking Regionmask.

    Parameters
    ----------
    zoom : int
        The zoom level of the expected tile source

    checkIntersect : bool
        If True, exclude tiles which do not intersect with the RegionMask's geometry

    asGeom : bool
        If True, returns tuple of ogr.Geometries in stead of (xi,yi,zoom) tuples

    Returns
    -------
    Generator of Geometries or (xi,yi,zoom) tuples
    """
    yield from GEOM.subTiles(self.geometry, zoom, checkIntersect=checkIntersect, asGeom=asGeom)

warp

warp(
    source,
    output: str | None = None,
    resolutionDiv=1,
    returnMatrix=True,
    applyMask=True,
    noData=None,
    resampleAlg="bilinear",
    **kwargs,
) -> Dataset | ndarray | None

Convenience wrapper for geokit.raster.warp() which automatically sets 'srs', 'bounds', 'pixelWidth', and 'pixelHeight' inputs.

Note:

When creating an 'in memory' raster vs one which is saved to disk, a slightly different algorithm is used which can sometimes add an extra row of pixels. Be aware of this if you intend to compare value-matricies directly from rasters generated with this function.

Parameters:

• source

  (str) –

  The path to the raster file to warp

• output

  (str; optional, default: None ) –

  A path to an output file to write to

• resampleAlg

  (str; optional, default: 'bilinear' ) –

  The resampling algorithm to use when warping values * Knowing which option to use can have significant impacts! * Options are: 'nearesampleAlg=resampleAlg, r', 'bilinear', 'cubic', 'average'

• resolutionDiv

  (int, default: 1 ) –

  The factor by which to divide the RegionMask's native resolution * This is useful if you need to represent very fine details

• returnAsMatrix

  (bool) –

  When True, the resulting raster's matrix is return * Should have the same dimensions as the RegionMask's mask matrix

• applyMask

  (bool, default: True ) –

  When True, the RegionMask's mask will be applied to the outputData as described by RegionMask.applyMask

• noData

  (numeric, default: None ) –

  The noData value to use when applying the mask

• **kwargs

  –

  All other keywargs are passed on to geokit.raster.warp()

Returns:

• * If 'output' is None: gdal.Dataset –
• * If 'output' is a string: None –

Source code in geokit/core/regionmask.py

def warp(
    self,
    source,
    output: str | None = None,
    resolutionDiv=1,
    returnMatrix=True,
    applyMask=True,
    noData=None,
    resampleAlg="bilinear",
    **kwargs,
) -> gdal.Dataset | np.ndarray | None:
    """Convenience wrapper for geokit.raster.warp() which automatically sets
    'srs', 'bounds', 'pixelWidth', and 'pixelHeight' inputs.

    Note:
    -----
    When creating an 'in memory' raster vs one which is saved to disk, a slightly
    different algorithm is used which can sometimes add an extra row of pixels. Be
    aware of this if you intend to compare value-matricies directly from rasters
    generated with this function.

    Parameters
    ----------
    source : str
        The path to the raster file to warp

    output : str; optional
        A path to an output file to write to

    resampleAlg : str; optional
        The resampling algorithm to use when warping values
        * Knowing which option to use can have significant impacts!
        * Options are: 'nearesampleAlg=resampleAlg, r', 'bilinear', 'cubic',
          'average'

    resolutionDiv : int
        The factor by which to divide the RegionMask's native resolution
        * This is useful if you need to represent very fine details

    returnAsMatrix : bool
        When True, the resulting raster's matrix is return
        * Should have the same dimensions as the RegionMask's mask matrix

    applyMask : bool
        When True, the RegionMask's mask will be applied to the outputData
        as described by RegionMask.applyMask

    noData : numeric
        The noData value to use when applying the mask

    **kwargs:
        All other keywargs are passed on to geokit.raster.warp()

    Returns
    -------
    * If 'output' is None: gdal.Dataset
    * If 'output' is a string: None
    """
    pW, pH = self._resolve(resolutionDiv)

    # do warp
    if returnMatrix:
        newDS = self.extent.warp(
            source=source,
            pixelWidth=pW,
            pixelHeight=pH,
            resampleAlg=resampleAlg,
            output=output,
            noData=noData,
            **kwargs,
        )
    else:
        if applyMask:
            if "cutline" in kwargs:
                raise GeoKitRegionMaskError(
                    "Cannot apply both a cutline and the mask when returning the warped dataset"
                )
            newDS = self.extent.warp(
                source=source,
                pixelWidth=pW,
                pixelHeight=pH,
                resampleAlg=resampleAlg,
                cutline=self.vectorPath,
                output=output,
                noData=noData,
                **kwargs,
            )
        else:
            newDS = self.extent.warp(
                source=source,
                pixelWidth=pW,
                pixelHeight=pH,
                resampleAlg=resampleAlg,
                output=output,
                noData=noData,
                **kwargs,
            )

    if not returnMatrix:
        return newDS

    # Read array
    if newDS is None:
        newDS = output
    final = RASTER.extractMatrix(newDS)

    # clean up
    del newDS

    # Apply mask, maybe
    if applyMask:
        final = self.applyMask(final, noData)

    # Return
    return final

checkSimilarRasters

checkSimilarRasters(datasets, rtol=0)

Parameters:

• datasets

  (string or list) –

  glob string path describing datasets to combine, alternatively list of gdal.Datasets or iterable object with paths.

• rtol

  ((int, float), default: 0 ) –

  The relative tolerance that is allowed for numeric deviations. By default 0, i.e. an exact match (within data type accuracy) is required.

Returns:

• output datasets: list –

  List of osgeo.gdal.Datasets with similar contexts.

Source code in geokit/_algorithms/combineSimilarRasters.py

def checkSimilarRasters(
    datasets,
    rtol=0,
):
    """
    Parameters
    ----------
    datasets : string or list
        glob string path describing datasets to combine, alternatively list of
        gdal.Datasets or iterable object with paths.
    rtol : int, float, optional
        The relative tolerance that is allowed for numeric deviations. By
        default 0, i.e. an exact match (within data type accuracy) is required.

    Returns
    -------
    output datasets: list
        List of osgeo.gdal.Datasets with similar contexts.
    """
    assert isinstance(rtol, (int, float)) and rtol >= 0, f"rtol must be a float or int >= 0"
    # ensure we have a list of raster datasets
    if isinstance(datasets, str):
        datasets = glob(datasets)
        if len(datasets) == 0:
            raise FileNotFoundError(f"datasets given as a string but does not lead to any existing files: '{datasets}'")
        datasets.sort()
    if not isinstance(datasets, list):
        raise TypeError(f"datasets must be a list")

    # check and load all datasets
    _datasets = []
    for dataset in datasets:
        if isinstance(dataset, str):
            if not os.path.isfile(dataset):
                raise FileNotFoundError(f"datasets string entry is not an existing file: '{dataset}'")
            _datasets.append(loadRaster(dataset))
        elif isinstance(dataset, gdal.Dataset):
            _datasets.append(dataset)
        else:
            raise TypeError(f"datasets must contain only string or osgeo.gdal.Dataset entries.")
    datasets = _datasets

    # get all raster infos
    infoDataset = [rasterInfo(sourceDS=d, compute_statistics=True) for d in datasets]
    # check all relevant variables
    for rInfo in infoDataset[1:]:
        # same srs is required
        if not infoDataset[0].srs.IsSame(rInfo.srs):
            raise GeoKitError(f"SRS mismatch between datasets.")
        # pixel width and height must be same/similar
        if not np.isclose(infoDataset[0].dx, rInfo.dx, rtol=rtol, atol=0):
            raise GeoKitError(f"dx mismatch between datasets.")
        if not np.isclose(infoDataset[0].dy, rInfo.dy, rtol=rtol, atol=0):
            raise GeoKitError(f"dy mismatch between datasets.")
        # bounds shift must be an exact/close match to a multiple of dx/dy
        diffx = infoDataset[0].bounds[0] - rInfo.bounds[0]
        if not (
            round(diffx / rInfo.dx, 0) == 0
            or np.isclose(diffx / round(diffx / rInfo.dx, 0), rInfo.dx, rtol=rtol, atol=0)
        ):
            raise GeoKitError(f"horizontal bounds shift between datasets is not a multiple of dx.")
        diffy = infoDataset[0].bounds[1] - rInfo.bounds[1]
        if not (
            round(diffy / rInfo.dy, 0) == 0
            or np.isclose(diffy / round(diffy / rInfo.dy, 0), rInfo.dy, rtol=rtol, atol=0)
        ):
            raise GeoKitError(f"vertical bounds shift between datasets is not a multiple of dy.")
        # # noData should be the same
        # if not infoDataset[0].noData == rInfo.noData:
        #     raise GeoKitError(f"noData mismatch between datasets.")

        # noData equality
        if not nodata_equal(infoDataset[0].noData, rInfo.noData):
            raise GeoKitError("noData mismatch between datasets.")

    # make sure the datatypes are the same or can be combined
    list_of_datatypes_in_raster = [rInfo.data_type_name_str for rInfo in infoDataset]
    list_of_minimum_values_in_raster = [rInfo.minimum_value for rInfo in infoDataset]
    list_of_maximum_values_in_raster = [rInfo.maximum_value for rInfo in infoDataset]
    if rtol == 0 and not len(set(list_of_datatypes_in_raster)):
        # no tolerance allowed - assume dtypes must also match exactly
        raise TypeError(f"dtypes or rasters differ but rtol is zero.")
    elif rtol > 0:
        # accept different dtypes as long as they can be combined into one
        list_of_numbers_to_consider = [*list_of_minimum_values_in_raster, *list_of_maximum_values_in_raster]
        MinimumCDataTypeHandler.get_valid_gdal_data_type_as_string(
            list_of_numbers=list_of_numbers_to_consider, minimum_gdal_type_list=list_of_datatypes_in_raster
        )  # fail if no common dtype
    # return list of preloaded, similar datasets
    return datasets

combineSimilarRasters

combineSimilarRasters(
    datasets,
    output=None,
    combiningFunc=None,
    verbose=True,
    updateMeta=False,
    allowNumericMismatch=False,
    **kwargs,
)

Combines several similar raster files into one single raster file.

Parameters:

• datasets

  (string or list) –

  glob string path describing datasets to combine, alternatively list of gdal.Datasets or iterable object with paths.

• output

  (string, default: None ) –

  Filepath to output raster file. If it is an existing file, datasets will be added to output. Recommended to create a new file every time though. If None, no output dataset will be loaded or created on disk and output dataset kept in memory only, by default None

• combiningFunc

  ([type], default: None ) –

  Allows customized functions to combine matrices, by default None

• verbose

  (bool, default: True ) –

  If True, additional status print stamenets will be issued, by default True

• updateMeta

  (bool, default: False ) –

  If True, metadata of output dataset will be a combination of all input rasters, by default False. NOTE: In the case of multiple values for the metadata keys, the last dataset metadata will take precedence.

• allowNumericMismatch

  (bool, default: False ) –

  If True, minor deviations in raster context will be ignored/corrected. By default False, i.e. only exactly similar rasters will be combined.

• **kwargs

  –

  Will be passed on to geokit.raster.createRaster().

Returns:

• output dataset: osgeo.gdal.Dataset –

  Raster file containing the combined matrices of all input datasets.

Source code in geokit/_algorithms/combineSimilarRasters.py

def combineSimilarRasters(
    datasets,
    output=None,
    combiningFunc=None,
    verbose=True,
    updateMeta=False,
    allowNumericMismatch=False,
    **kwargs,
):
    """
    Combines several similar raster files into one single raster file.

    Parameters
    ----------
    datasets : string or list
        glob string path describing datasets to combine, alternatively list of
        gdal.Datasets or iterable object with paths.
    output : string, optional
        Filepath to output raster file. If it is an existing file, datasets will
        be added to output. Recommended to create a new file every time though.
        If None, no output dataset will be loaded or created on disk and output
        dataset kept in memory only, by default None
    combiningFunc : [type], optional
        Allows customized functions to combine matrices, by default None
    verbose : bool, optional
        If True, additional status print stamenets will be issued, by default
        True
    updateMeta : bool, optional
        If True, metadata of output dataset will be a combination of all input
        rasters, by default False.
        NOTE: In the case of multiple values for the metadata keys, the last
        dataset metadata will take precedence.
    allowNumericMismatch : bool, optional
        If True, minor deviations in raster context will be ignored/corrected.
        By default False, i.e. only exactly similar rasters will be combined.
    **kwargs
        Will be passed on to geokit.raster.createRaster().

    Returns
    -------
    output dataset: osgeo.gdal.Dataset
        Raster file containing the combined matrices of all input datasets.
    """
    if not isinstance(allowNumericMismatch, bool):
        raise TypeError(f"allowNumericMismatch must be boolean.")

    # CHECK AND PREPROCESS INPUT DATASETS

    # Ensure we have a list of raster datasets
    if isinstance(datasets, str):
        datasets = glob(datasets)
        datasets.sort()
    elif isinstance(datasets, gdal.Dataset):
        datasets = [
            datasets,
        ]
    else:  # assume datasets is iterable
        datasets = list(datasets)

    # make sure we do actually have rasters and they are all indeed "similar"
    if len(datasets) == 0:
        raise GeoKitError("No datasets given")
    rtol = 0
    if allowNumericMismatch:
        rtol = 0.0001  # allow a numeric deviation of 0.01%
    datasets = checkSimilarRasters(datasets=datasets, rtol=rtol)

    # GET REFERENCE CONTEXT FOR THE OUTPUT RASTER

    # determine info for all datasets
    raster_info_list = [rasterInfo(sourceDS=d, compute_statistics=True) for d in datasets]

    # get reference srs - are all the same thanks to checkSimilarRasters
    srs_ref = raster_info_list[0].srs

    # get the unique actual dtypes in input rasters
    data_type_list = []
    minimum_and_maximum_values = []
    for current_raster_info in raster_info_list:
        if current_raster_info.dtype is None:
            raise GeoKitError("Input raster has no dtype.")
        data_type_list.append(current_raster_info.data_type_name_str)
        minimum_and_maximum_values.append(current_raster_info.minimum_value)
        minimum_and_maximum_values.append(current_raster_info.maximum_value)

    # now get the most lightweight commonly usable dtype
    gdal_data_type_as_string = MinimumCDataTypeHandler.get_valid_gdal_data_type_as_string(
        list_of_numbers=minimum_and_maximum_values, minimum_gdal_type_list=data_type_list
    )
    # get_common_dtype(dtypes=dtypes, fallback=None)

    # get the reference resolution in x and y dir as the most commonly used value
    dx_ref = statistics.mode([_i.pixelWidth for _i in raster_info_list])
    dy_ref = statistics.mode([_i.pixelHeight for _i in raster_info_list])

    # try to align all bounds to the first matching raster which has correct x_ref and y_ref resolution
    i_match = next(
        (
            i
            for i, (x, y) in enumerate(
                zip(
                    [_i.pixelWidth for _i in raster_info_list],
                    [_i.pixelHeight for _i in raster_info_list],
                )
            )
            if x == dx_ref and y == dy_ref
        ),
        None,
    )
    if i_match is not None:
        # we have a "perfect" raster, use the min. bounds of that one as reference for all other rasters
        boundsXmin_ref = raster_info_list[i_match].bounds[0]
        boundsYmin_ref = raster_info_list[i_match].bounds[2]
    else:
        # we do not have any raster which matches both r_ref any y_ref in its resolution. Simply use the first raster.
        boundsXmin_ref = raster_info_list[0].bounds[0]
        boundsYmin_ref = raster_info_list[0].bounds[2]
        if verbose:
            print(
                datetime.datetime.now(),
                f"NOTE: None of the rasters matches both reference resolutions in x and y direction. Use first raster bounds as reference.",
                flush=True,
            )

    # calculate the possibly adapted bounds for all datasets
    boundsSet = []
    for _info in raster_info_list:
        # calculate the new bounds by aligning bottom left corner with boundsXmin_ref/boundsYmin_ref + multiple of cell size
        _bounds_Xmin = boundsXmin_ref + round((_info.bounds[0] - boundsXmin_ref) / dx_ref) * dx_ref
        _bounds_Ymin = boundsYmin_ref + round((_info.bounds[1] - boundsYmin_ref) / dy_ref) * dy_ref
        _bounds = (
            _bounds_Xmin,
            _bounds_Ymin,
            _bounds_Xmin + _info.xWinSize * dx_ref,
            _bounds_Ymin + _info.yWinSize * dy_ref,
        )
        # make sure the number of cells would remain the same in the new bounds
        assert round(_info.xWinSize * dx_ref / _info.dx, 0) == _info.xWinSize, (
            f"The change in bounds width would lead to a different amount of cell columns in the original raster."
        )
        assert round(_info.yWinSize * dy_ref / _info.dy, 0) == _info.yWinSize, (
            f"The change in bounds height would lead to a different amount of cell rows in the original raster."
        )
        boundsSet.append(_bounds)

    # get summary info about the whole dataset group
    dataXMin = min([i[0] for i in boundsSet])
    dataXMax = max([i[2] for i in boundsSet])
    dataYMin = min([i[1] for i in boundsSet])
    dataYMax = max([i[3] for i in boundsSet])

    # get noData value from kwargs, else take from dataset infos
    noData_ref = kwargs.pop("noData", None)
    if noData_ref is None:
        noDataSet = set([i.noData for i in raster_info_list])
        assert len(noDataSet) == 1  # make sure, is enforced by checkSimilarRasters
        noData_ref = noDataSet.pop()

    # Maybe create a new output dataset
    if isinstance(output, str):
        if not os.path.isfile(output):
            # we will need to create a output source
            createRaster(
                bounds=(dataXMin, dataYMin, dataXMax, dataYMax),
                output=output,
                dtype=gdal_data_type_as_string,
                pixelWidth=dx_ref,
                pixelHeight=dy_ref,
                noData=noData_ref,
                srs=srs_ref,
                fill=noData_ref,
                **kwargs,
            )
        else:
            warn(
                "WARNING: Overwriting existing output file. Sometimes writing to an non empty output fails. Recommended to write to a non existing location instead and include maser into datasets."
            )
    elif output is None:
        # create raster in memory
        outputDS = createRaster(
            bounds=(dataXMin, dataYMin, dataXMax, dataYMax),
            dtype=gdal_data_type_as_string,
            pixelWidth=dx_ref,
            pixelHeight=dy_ref,
            noData=noData_ref,
            srs=srs_ref,
            fill=noData_ref,
            **kwargs,
        )
    else:
        raise TypeError(
            "output must be None or a str formatted file path to an existing output file or a file to be created."
        )

    # open output dataset if required and check parameters
    if not output is None:
        outputDS = gdal.Open(output, gdal.GA_Update)
    mInfo = rasterInfo(outputDS)
    mExtent = Extent(mInfo.bounds, srs=mInfo.srs)

    outputBand = outputDS.GetRasterBand(1)

    # Make a meta container
    if updateMeta:
        meta = outputDS.GetMetadata_Dict()

    # Add each dataset to output
    for i in range(len(datasets)):
        if verbose:
            print(
                datetime.datetime.now(),
                f"Now adding raster No. {i + 1}/{len(datasets)}",
            )
        # create dataset extent
        dExtent = Extent(boundsSet[i], srs=srs_ref)

        # extract the dataset's matrix
        dMatrix = extractMatrix(datasets[i])
        if not raster_info_list[i].yAtTop:
            dMatrix = dMatrix[::-1, :]

        # Calculate starting indices
        idx = mExtent.findWithin(dExtent, (mInfo.dx, mInfo.dy), yAtTop=mInfo.yAtTop)

        # Get output data
        mMatrix = outputBand.ReadAsArray(xoff=idx.xStart, yoff=idx.yStart, win_xsize=idx.xWin, win_ysize=idx.yWin)
        if mMatrix is None:
            raise GeoKitError("mMatrix is None")

        # create selector
        if not combiningFunc is None:
            # update rasterInfo since we might have slightly changed cell/bounds info
            rInfo_dict = raster_info_list[i]._asdict()
            rInfo_dict["bounds"] = boundsSet[i]
            rInfo_dict["pixelWidth"] = dx_ref
            rInfo_dict["pixelHeight"] = dy_ref
            rInfo_dict["dx"] = dx_ref
            rInfo_dict["dy"] = dy_ref
            rInfo_dict["xMin"] = boundsSet[i][0]
            rInfo_dict["yMin"] = boundsSet[i][1]
            rInfo_dict["xMax"] = boundsSet[i][2]
            rInfo_dict["yMax"] = boundsSet[i][3]
            rInfo_upd = RasterInfo(**rInfo_dict)
            writeMatrix = combiningFunc(mMatrix=mMatrix, mInfo=mInfo, dMatrix=dMatrix, dInfo=rInfo_upd)
        elif not raster_info_list[i].noData is None:
            sel = dMatrix != raster_info_list[i].noData
            mMatrix[sel] = dMatrix[sel]
            writeMatrix = mMatrix
        else:
            writeMatrix = dMatrix

        # Add to output
        outputBand.WriteArray(writeMatrix, idx.xStart, idx.yStart)
        outputBand.FlushCache()

        # update metaData, maybe
        if updateMeta:
            meta.update(raster_info_list[i].meta)

    if updateMeta:
        outputDS.SetMetadata(meta)

    # Write final raster
    outputDS.FlushCache()
    outputBand.ComputeRasterMinMax(0)
    outputBand.ComputeBandStats(0)

    return outputDS