Data & configuration¶

Reading, preparing and writing rasters, vectors and tables, saving and plotting results, colour tables, collections of files and the thread settings shared by every tool. The page Loading and saving data shows them together.

Every example runs

The examples on these pages use the maps that ship with Landschaft (ls.example_data) and are run as tests, outputs included.

Rasters¶

load_raster ¶

load_raster(
    source,
    *,
    band: int | Iterable[int] | None = None,
    masked: bool | str = "auto",
    chunks: Any = None,
    validate: str | None = None,
    like=None,
    resampling: str = "auto",
    clip=None,
    stack: bool = False,
    errors: str = "raise",
    recursive: bool = False,
)

Read one raster, or collect many.

One function for every raster landschaft reads: habitat and land-cover maps, continuous surfaces, species stacks, costs. The result is a georeferenced xarray.DataArray (CRS, transform and NoData via .rio) that write_raster() writes back as it was.

Parameters:

Name	Type	Description	Default
`source`		One raster (path or GDAL URL), a directory (its `.tif`/`.tiff` files), a glob pattern (`"tiles/*.tif"`) or a list of any of these.	required
`band`	`int \| Iterable[int] \| None`	Bands to read (1-based): `None` (default) reads them all, an integer one band, a list several. One band gives a `(y, x)` array, several a `(band, y, x)` array whose `band` coordinate holds the band descriptions when every band has a distinct one (else the band numbers). Several rasters: each file must have one band, or give `band=N`.	`None`
`masked`	`bool \| str`	NoData handling. `"auto"` (default): in floating-point rasters NoData becomes NaN; integer rasters (class maps) keep their codes, with the NoData value in `.rio.nodata`, so they are not turned into floats. `True`: NoData becomes NaN in every raster (integers become `float32`, or `float64` for 32- and 64-bit integers). `False`: the values as stored.	`'auto'`
`chunks`	`Any`	Single raster only: `None` (default) loads the data into memory; anything else (e.g. `"auto"`) returns a lazy dask-backed array.	`None`
`validate`	`str \| None`	`"mspa"` checks that the map is an MSPA input: `uint8` coded 0 = missing, 1 = background, 2 = foreground, with some foreground and some background. `None` (default) reads any raster; the habitat-pattern tools check their input themselves.	`None`
`like`		A raster (DataArray or path) whose grid the result is put on, reprojected and resampled when the CRS, resolution or extent differ, with the coordinates of `like` (so that arithmetic between the two aligns cell by cell). Cells outside the source are NoData. Several rasters (a list, a folder or a pattern) with `like` give one `(band, y, x)` stack named by the files, whatever their grids.	`None`
`resampling`	`str`	With `like`: `"auto"` (default) takes the nearest cell for integer rasters (classes) and interpolates floats bilinearly; or any GDAL method: `"nearest"`, `"bilinear"`, `"cubic"`, `"average"`, `"mode"`, `"sum"`, `"min"`, `"max"`...	`'auto'`
`clip`		A region (polygons as a GeoDataFrame, GeoSeries or shapely geometries, or bounds `(xmin, ymin, xmax, ymax)` in the raster's CRS): only the window of the region is read, the cells outside the polygons are NoData (see `clip()`). For a small area of a large raster, reading is much faster.	`None`
`stack`	`bool`	Several rasters only: return a lazy `(file, y, x)` DataArray (one dask chunk per file) instead of a `RasterCollection`. All rasters must share the same grid.	`False`
`errors`	`str`	Several rasters only: `"raise"` (default) if any file is unusable, `"skip"` to leave such files out (listed in `.skipped`).	`'raise'`
`recursive`	`bool`	Several rasters only: also search sub-directories.	`False`

Returns:

Type	Description
`DataArray or RasterCollection`	For one raster, a georeferenced `(y, x)` or `(band, y, x)` array. For several, a `RasterCollection` that reads only the file headers now; pass it to `mspa()` to process every raster in parallel. With `stack=True`, a lazy DataArray.

Raises:

Type	Description
`InvalidInputError`	If a file does not satisfy the requested checks.

Examples:

A habitat map (0 missing, 1 background, 2 habitat) and a land-cover map:

>>> import landschaft as ls
>>> forest = ls.load_raster(ls.example_data("input"))
>>> forest.shape, str(forest.dtype)
((1000, 1000), 'uint8')
>>> landcover = ls.load_raster(ls.example_data("clc"))
>>> landcover.rio.crs.to_epsg()
3035

On the grid of another raster (here a coarser one in another CRS):

>>> coarse = ls.write_raster(landcover.rio.reproject("EPSG:4326")[::4, ::4], "coarse.tif")
>>> ls.load_raster(coarse, like=landcover).rio.crs == landcover.rio.crs
True

Only a study area, read from a large raster:

>>> x0, y0, x1, y1 = landcover.rio.bounds()
>>> ls.load_raster(ls.example_data("clc"), clip=(x0, y0, x0 + 10000, y0 + 5000)).shape
(200, 400)

Checked as an MSPA input on reading:

>>> forest = ls.load_raster(ls.example_data("input"), validate="mspa")

A multi-band stack keeps its band names; NoData of floats is NaN:

>>> import numpy as np
>>> import xarray as xr
>>> habitats = xr.concat([landcover.isin([23, 24, 25]), landcover == 41], dim="band").astype("float32")
>>> habitats = habitats.assign_coords(band=["forest", "water"]).where(landcover != 1)  # no data in towns
>>> _ = ls.write_raster(habitats, "stacks/habitats.tif")   # creates the folder
>>> stack = ls.load_raster("stacks/habitats.tif")
>>> stack.dims, stack.band.values.tolist()
(('band', 'y', 'x'), ['forest', 'water'])
>>> bool(stack.sel(band="water").isnull().any()), ls.load_raster("stacks/habitats.tif", band=2).dims
(True, ('y', 'x'))

Many rasters at once (a folder, a pattern or a list) give a collection:

>>> _ = ls.write_raster(forest[:500], "north.tif")
>>> _ = ls.write_raster(forest[500:], "south.tif")
>>> maps = ls.load_raster(["north.tif", "south.tif"])
>>> len(maps)
2

write_raster ¶

write_raster(
    data,
    path: PathLike,
    *,
    like=None,
    crs=None,
    transform=None,
    nodata=None,
    dtype=None,
    band_names=None,
    colormap: str | dict | None = "auto",
    compress: str = "lzw",
    description: str | None = None,
    driver: str | None = None,
) -> Path

Write a raster: one band or several, any format GDAL writes.

One function for every raster landschaft writes. It keeps what the data says about itself (georeferencing, NoData, band names, the colours and parameters of landschaft's products) and stores it so that GIS software and load_raster() read it back as it was:

floating-point data with NaN get NaN as NoData (or nodata, which then replaces the NaN); integers keep the .rio.nodata of the data;
bool is written as uint8 and 64-bit integers as the smallest integer type that holds them;
GeoTIFFs are compressed with the predictor that suits the type, tiled when larger than 512 x 512 and BigTIFF when needed; the format follows the file extension (.tif, .asc, .img, .nc...) or driver ("COG": Cloud-Optimised GeoTIFF);
dask arrays are written block by block, never whole in memory;
missing folders of path are created.

Parameters:

Name	Type	Description	Default
`data`		2-D `(y, x)` or 3-D `(band, y, x)` `xarray.DataArray` (georeferencing taken from `.rio`), numpy or dask array.	required
`path`	`PathLike`	Output file.	required
`like`		A raster (path or DataArray) on the same grid whose CRS and transform georeference `data`, e.g. a numpy result computed from an input.	`None`
`crs`		Coordinate reference system and affine transform; they override those of `data` and `like`.	`None`
`transform`		Coordinate reference system and affine transform; they override those of `data` and `like`.	`None`
`nodata`		NoData value; overrides the one of `data`. NaN in the data are written as this value.	`None`
`dtype`		Data type to write (default: the data's, see above). Floats cast to integers get `nodata` (or the type's maximum, minimum if signed) where they are NaN.	`None`
`band_names`		Band descriptions (default: a `band` coordinate of strings, e.g. from `load_raster()`).	`None`
`colormap`	`str \| dict \| None`	`"auto"` (default) embeds the colour table of the product when `data` is a class map made by landschaft (MSPA, fragmentation, change...), `"mspa"` forces the MSPA colours, a `{value: colour}` mapping (colours as `"#rrggbb"`, names or `(r, g, b[, a])`) embeds your own and `None` none. One-band `uint8`/`uint16` only.	`'auto'`
`compress`	`str`	GeoTIFF compression: `"lzw"` (default), `"deflate"`, `"zstd"`, `None`...	`'lzw'`
`description`	`str \| None`	TIFF image description (MSPA results record their parameters there).	`None`
`driver`	`str \| None`	GDAL format (default: from the extension of `path`).	`None`

Returns:

Type	Description
`Path`	The written file.

Examples:

A landschaft result keeps its colours and parameters:

>>> import landschaft as ls
>>> forest = ls.load_raster(ls.example_data("input"))
>>> ls.write_raster(ls.mspa(forest), "forest_mspa.tif").name
'forest_mspa.tif'

A numpy array georeferenced like the input, with your own colours:

>>> import numpy as np
>>> mask = (forest.values == 2).astype("uint8")
>>> path = ls.write_raster(mask, "mask.tif", like=forest,
...                        colormap={0: (255, 255, 255), 1: (21, 128, 61)})

Floats with NaN: NaN is the NoData, or the value you choose:

>>> share = forest.where(forest > 0) / 2.0
>>> float(ls.load_raster(ls.write_raster(share, "share.tif")).rio.nodata)
nan
>>> ls.write_raster(share, "share.asc", nodata=-9999).name
'share.asc'

Several bands, named:

>>> two = ls.write_raster(np.stack([mask, 1 - mask]), "stacks/two.tif", like=forest, band_names=["forest", "open"])
>>> ls.load_raster(two).band.values.tolist()
['forest', 'open']

match_grid ¶

match_grid(data: DataArray, like, *, resampling: str = 'auto') -> DataArray

Put a raster on the grid of another: same CRS, cells and coordinates.

Rasters from different sources (a species model, a cost layer, a land-cover map) rarely share their grid; most tools need them cell by cell. match_grid reprojects and resamples data onto the grid of like with GDAL's resampling (in parallel), and gives it the coordinates of like. load_raster() does the same on reading with like=, and reproject() with like= (which also reprojects to a CRS, and vectors).

Parameters:

Name	Type	Description	Default
`data`	`DataArray`	Georeferenced `(y, x)` or `(band, y, x)` DataArray.	required
`like`		The target grid: a DataArray or a raster path.	required
`resampling`	`str`	`"auto"` (default): nearest cell for integer rasters, bilinear for floats; or a GDAL method (`"nearest"`, `"bilinear"`, `"cubic"`, `"average"`, `"mode"`, `"sum"`...).	`'auto'`

Returns:

Type	Description
`DataArray`	`data` on the grid of `like`; cells outside `data` are NoData (NaN for floats). A lazy (dask) raster stays lazy.

Examples:

>>> import landschaft as ls
>>> landcover = ls.load_raster(ls.example_data("clc"))
>>> coarse = landcover.coarsen(x=4, y=4).mean()          # 100 m instead of 25 m
>>> back = ls.match_grid(coarse, landcover)
>>> back.shape == landcover.shape, bool((back.x == landcover.x).all())
(True, True)

reproject ¶

reproject(
    data,
    crs=None,
    *,
    like=None,
    res=None,
    resampling: str = "auto",
    tolerance: float = 0.125,
    densify: float | None = None,
    n_threads: int | None = None,
)

Reproject a raster or vector to another CRS, or onto the grid of another raster.

Like terra's project(): a raster goes to a CRS (on the grid GDAL suggests, or with the cell size res), or onto the grid of a reference raster (like=: its CRS, cells and extent, so that the two combine cell by cell); a vector goes to a CRS or to the CRS of a reference. With res and no CRS, a raster is resampled in its own CRS.

Rasters are warped by GDAL's resampling kernels (ported to C++, run in parallel), on the source coordinates GDAL's own transformation gives for a fixed lattice of points of the result, interpolated in between: the values do not depend on the number of threads, nor on whether the raster is in memory, lazy (dask, load_raster(..., chunks=...)) or a file, so a raster larger than memory gets the values it would get in memory. With tolerance=0 they are GDAL's (gdalwarp -et 0), bit for bit. A lazy raster stays lazy; a path is read window by window, never whole.

Parameters:

Name	Type	Description	Default
`data`		A georeferenced `(y, x)`, `(band, y, x)`... DataArray (any dimensions before `y` and `x`), a GeoDataFrame or GeoSeries, or the path of a raster or vector file.	required
`crs`		The target CRS: `"EPSG:4326"`, `4326`, a WKT or PROJ string, a pyproj or rasterio CRS, or a raster or vector (object or path) whose CRS is taken (and only its CRS).	`None`
`like`		A reference: a raster (DataArray or path) whose grid a raster is put on (CRS, cells, extent and coordinates), or whose CRS a vector takes; or a vector (object or path) whose CRS is taken. Not with `crs`.	`None`
`res`		Rasters: the cell size of the result, one number or `(xres, yres)`, in the units of the target CRS (default: about the source's, as GDAL suggests). Not with a raster `like`.	`None`
`resampling`	`str`	Rasters: `"auto"` (default) takes the nearest cell for integer rasters (classes) and interpolates floats bilinearly; or any GDAL method: `"nearest"`, `"bilinear"`, `"cubic"`, `"cubic_spline"`, `"lanczos"`, `"average"`, `"mode"`, `"min"`, `"max"`, `"med"`, `"q1"`, `"q3"`, `"sum"`, `"rms"`.	`'auto'`
`tolerance`	`float`	Rasters: the largest error, in source cells, of the coordinates interpolated between the points transformed exactly (default 0.125, as `gdalwarp`; checked midway between them, finer where the projection bends, as GDAL's approximate transformation does). `0` transforms every cell exactly, four to seven times slower: GDAL's values bit for bit.	`0.125`
`densify`	`float \| None`	Vectors: add vertices so that no segment is longer than this (in the units of the source CRS) before reprojecting, so that long straight edges follow the curvature of the target CRS. Default: the vertices as they are.	`None`
`n_threads`	`int \| None`	Rasters in memory or from a file: threads (default: all cores but one); the result does not depend on it. A lazy raster stays lazy (its tasks run when computed, with dask's threads).	`None`

Returns:

Type	Description
`(DataArray, GeoDataFrame or GeoSeries)`	Rasters: the same type and dimensions; cells outside the source are NoData (NaN for floats, the NoData value of integer rasters, else the largest value of unsigned types and the smallest of signed ones). Boolean rasters become `uint8`. Vectors: the same type, attributes kept.

Examples:

>>> import landschaft as ls
>>> landcover = ls.load_raster(ls.example_data("clc"))            # 25 m, EPSG:3035
>>> lonlat = ls.reproject(landcover, "EPSG:4326")
>>> lonlat.rio.crs.to_epsg(), str(lonlat.dtype), int(lonlat.rio.nodata)   # classes kept
(4326, 'uint8', 255)

With a cell size, in metres of UTM 32N:

>>> utm = ls.reproject(landcover, 32632, res=100)
>>> utm.rio.resolution()
(100.0, -100.0)

Onto the grid of a reference raster, cell by cell:

>>> back = ls.reproject(lonlat, like=landcover)
>>> back.shape == landcover.shape, bool((back.x == landcover.x).all())
(True, True)

Vectors go to the CRS of a reference (raster or vector):

>>> import geopandas as gpd
>>> from shapely.geometry import box
>>> x0, y0, x1, y1 = landcover.rio.bounds()
>>> area = gpd.GeoDataFrame({"name": ["study"]}, geometry=[box(x0, y0, x1, y1)], crs=landcover.rio.crs)
>>> ls.reproject(area, like=lonlat).crs.to_epsg()
4326

clip ¶

clip(
    raster,
    region,
    *,
    crop: bool = True,
    all_touched: bool = False,
    invert: bool = False,
) -> DataArray

Cut a raster to a region: the study area, a biome, a municipality.

The raster is cropped to the region's bounding box and the cells outside its polygons become NoData (NaN for floats, the NoData value for integer maps, which keep their type). The polygons are reprojected to the raster's CRS; the cells kept are those rasterize() burns. On a lazy raster (load_raster(..., chunks=...)) the crop is lazy too; to read only the region from a file, use load_raster(path, clip=region).

Parameters:

Name	Type	Description	Default
`raster`		Georeferenced `(y, x)` or `(band, y, x)` DataArray.	required
`region`		Polygons (GeoDataFrame, GeoSeries, shapely geometries), or bounds `(xmin, ymin, xmax, ymax)` in the raster's CRS (a crop only).	required
`crop`	`bool`	Crop to the region's bounding box (default); `False` keeps the extent and only masks.	`True`
`all_touched`	`bool`	Keep every cell the polygons touch (default: the cells whose centre is inside).	`False`
`invert`	`bool`	Mask the inside of the region instead (e.g. to leave urban areas out); the extent is kept.	`False`

Returns:

Type	Description
`DataArray`

Examples:

>>> import geopandas as gpd
>>> from shapely.geometry import Point
>>> import landschaft as ls
>>> landcover = ls.load_raster(ls.example_data("clc"))
>>> x0, y0, x1, y1 = landcover.rio.bounds()
>>> area = gpd.GeoDataFrame(geometry=[Point(x0 + 10000, y0 + 10000).buffer(5000)], crs=landcover.rio.crs)
>>> study = ls.clip(landcover, area)
>>> study.shape, int(study.rio.nodata), str(study.dtype)              # 10 km box, classes kept
((400, 400), 255, 'uint8')
>>> round(int((study != 255).sum()) * 25 * 25 / 1e6, 1)              # km2 inside the circle
78.4

reclassify ¶

reclassify(
    data,
    table=None,
    *,
    foreground: Iterable[int] | int | Callable[[ndarray], ndarray] | None = None,
    missing: Iterable[int] | int | None = None,
    background: Iterable[int] | int | None = None,
    others="keep",
    dtype=None,
    n_threads: int | None = None,
)

Give the values of a raster new values: from a table, or into a habitat map.

With a table, every value (or range of values) becomes the value the table gives it: land-cover classes into resistances or costs, classes grouped, a continuous map cut into classes. The table is

a dict {old: new}; a key may be a range (low, high);
a list of pairs (old, new) or of ranges (low, high, new), which hold the values low <= v < high (-np.inf / np.inf for open ends); where ranges overlap, the first one wins;
a DataFrame (e.g. from load_table()) of two columns (old, new) or three (low, high, new), or a 2-D array of the same.

With foreground (or a value, a list of values or a function in place of the table), the result is a habitat map (0 = missing, 1 = background, 2 = habitat), the input of mspa(), fad() and the other habitat-pattern tools.

Parameters:

Name	Type	Description	Default
`data`		Raster: numpy array, DataArray (dask-backed ones stay lazy), any type.	required
`table`		The new values (see above).	`None`
`foreground`	`Iterable[int] \| int \| Callable[[ndarray], ndarray] \| None`	Habitat map: value(s) that become habitat (2), or a function of the values returning a boolean mask (e.g. `lambda a: a >= 30` on a tree-cover percentage).	`None`
`missing`	`Iterable[int] \| int \| None`	Habitat map: value(s) that become missing (0). NaN and the raster's NoData value are always missing.	`None`
`background`	`Iterable[int] \| int \| None`	Habitat map: value(s) that become background (1). By default every value that is neither habitat nor missing; when given, the remaining values become missing.	`None`
`others`		Table: what values the table does not hold become: `"keep"` (default, their own value), `"nodata"`, or a number. NoData stays NoData.	`'keep'`
`dtype`		Table: output type (default: the smallest that holds the new values, and the kept ones).	`None`
`n_threads`	`int \| None`	Threads (default: all cores but one).	`None`

Returns:

Type	Description
`DataArray or ndarray`	Same type and georeferencing as `data`.

Examples:

Land-cover classes into resistance to movement (the other classes: 50):

>>> import landschaft as ls
>>> landcover = ls.load_raster(ls.example_data("clc"))
>>> resistance = ls.reclassify(landcover, {23: 1, 24: 1, 25: 1, 18: 10, 1: 1000}, others=50)
>>> str(resistance.dtype), sorted(int(v) for v in set(resistance.values.ravel()))
('uint16', [1, 10, 50, 1000])

Ranges of a continuous map into classes, [low, high):

>>> import numpy as np
>>> footprint = np.array([[0.0, 3.5], [12.0, np.nan]])
>>> ls.reclassify(footprint, [(0, 4, 1), (4, 10, 2), (10, np.inf, 3)])
array([[ 1.,  1.],
       [ 3., nan]])

Forests of a CORINE land-cover map (classes 23, 24 and 25) as a habitat map:

>>> forest = ls.reclassify(landcover, foreground=[23, 24, 25])
>>> sorted(int(v) for v in set(forest.values.ravel()))
[1, 2]
>>> morph = ls.mspa(forest)

A rule on a continuous map, water (class 41) left out as missing:

>>> tall = ls.reclassify(landcover, foreground=lambda a: (a >= 23) & (a <= 25), missing=[41])

Vectors¶

load_vector ¶

load_vector(
    source,
    *,
    layer=None,
    like=None,
    x: str | None = None,
    y: str | None = None,
    crs=None,
    columns=None,
    bbox=None,
)

Read a vector layer (polygons, lines, points) as a GeoDataFrame.

Any format GDAL reads (GeoPackage, Shapefile, GeoJSON, FlatGeobuf, KML...), GeoParquet, or a table of point coordinates (CSV, Excel, Parquet): the coordinate columns are found by name (x/y, lon/lat, longitude/latitude, Darwin Core's decimalLongitude...).

Parameters:

Name	Type	Description	Default
`source`		File path (or a GDAL URL).	required
`layer`		Layer of a multi-layer file (GeoPackage...), by name or number.	`None`
`like`		A raster or vector (object or path) whose CRS the result is reprojected to, so that it overlays the rasters it will be used with.	`None`
`x`	`str \| None`	Coordinate columns of a table (default: found by name).	`None`
`y`	`str \| None`	Coordinate columns of a table (default: found by name).	`None`
`crs`		CRS of the coordinates of a table (default: longitude/latitude, EPSG:4326, for columns named like longitude and latitude, else none).	`None`
`columns`		Attribute columns to read (default: all).	`None`
`bbox`		`(xmin, ymin, xmax, ymax)`: only the features intersecting it, in the CRS of the file.	`None`

Returns:

Type	Description
`GeoDataFrame`

Examples:

>>> import pandas as pd
>>> import landschaft as ls
>>> sites = pd.DataFrame({"name": ["a", "b"], "longitude": [10.1, 10.4], "latitude": [52.3, 52.6]})
>>> _ = ls.write_table(sites, "sites.csv")
>>> points = ls.load_vector("sites.csv")                  # columns found by name, EPSG:4326
>>> points.crs.to_epsg(), len(points)
(4326, 2)
>>> landcover = ls.load_raster(ls.example_data("clc"))
>>> ls.load_vector("sites.csv", like=landcover).crs.to_epsg()    # on the raster's CRS
3035

write_vector ¶

write_vector(data, path, **kwargs) -> Path

Write a vector layer, the format given by the file extension.

.gpkg (GeoPackage), .shp, .geojson, .fgb (FlatGeobuf), .kml, .parquet (GeoParquet), .csv (geometry as WKT). Missing folders are created.

Parameters:

Name	Description	Default
`data`	A GeoDataFrame or GeoSeries (e.g. from `load_vector()` or `zonal_metrics()`).	required
`path`	Output file.	required
`**kwargs`	Passed to geopandas (`layer=` for a GeoPackage...).	`{}`

Returns:

Type	Description
`Path`

Examples:

>>> import geopandas as gpd
>>> from shapely.geometry import box
>>> import landschaft as ls
>>> reserves = gpd.GeoDataFrame({"name": ["north"]}, geometry=[box(0, 0, 10, 10)], crs="EPSG:3035")
>>> ls.write_vector(reserves, "out/reserves.gpkg").name
'reserves.gpkg'
>>> ls.load_vector("out/reserves.gpkg")["name"].tolist()
['north']

rasterize ¶

rasterize(
    vector,
    like,
    *,
    field: str | None = None,
    value=1,
    fill=0,
    all_touched: bool = False,
    merge: str = "replace",
    dtype=None,
    name: str | None = None,
) -> DataArray

Burn vector features into a raster on the grid of another raster.

Protected areas, municipalities, roads or sample sites become a raster that lines up cell by cell with like: a mask (1 inside, 0 outside), the values of an attribute, or the number of overlapping features. The features are reprojected to the CRS of like first.

Parameters:

Name	Type	Description	Default
`vector`		GeoDataFrame, GeoSeries, a list of shapely geometries or a vector path.	required
`like`		The grid: a DataArray or a raster path.	required
`field`	`str \| None`	Attribute whose values are burnt (default: `value` for every feature).	`None`
`value`		The value burnt when `field` is not given (default 1).	`1`
`fill`		The value of the cells no feature covers (default 0; `np.nan` for a float raster with NoData outside).	`0`
`all_touched`	`bool`	Burn every cell a feature touches (default: the cells whose centre is inside, which does not overestimate small polygons).	`False`
`merge`	`str`	Where features overlap: `"replace"` (the last one wins) or `"add"` (the values are summed, e.g. to count overlapping features).	`'replace'`
`dtype`		Output type (default: the smallest that holds the values and `fill`).	`None`
`name`	`str \| None`	Name of the result (default: `field`, else `"rasterized"`).	`None`

Returns:

Type	Description
`DataArray`	On the grid of `like`, with its coordinates and CRS.

Examples:

>>> import geopandas as gpd
>>> from shapely.geometry import box
>>> import landschaft as ls
>>> landcover = ls.load_raster(ls.example_data("clc"))
>>> x0, y0 = float(landcover.x[0]), float(landcover.y[-1])
>>> reserves = gpd.GeoDataFrame({"id": [1, 2]}, crs=landcover.rio.crs,
...                             geometry=[box(x0, y0, x0 + 5000, y0 + 5000), box(x0 + 20000, y0, x0 + 30000, y0 + 8000)])
>>> mask = ls.rasterize(reserves, like=landcover)          # 1 inside the reserves, 0 outside
>>> str(mask.dtype), int(mask.sum()) * 25 * 25 / 1e6       # km2
('uint8', 105.375)
>>> ls.rasterize(reserves, like=landcover, field="id").max().item()
2

Values by zone and at points¶

zonal_stats ¶

zonal_stats(
    raster,
    zones,
    stats=("count", "mean", "min", "max", "std"),
    *,
    all_touched: bool = False,
    n_threads: int | None = None,
)

Statistics of a raster in every zone: polygons, or the labels of a raster.

Mean NDVI per municipality, total population per watershed, the share of each land-cover class per property, the median priority rank of every protected area: computed in C++ in parallel, on every band of a stack at once. Polygons may overlap (e.g. buffers around points): each zone gets its own cells, the cells GDAL burns for it (rasterize()), reprojected to the raster's CRS first.

Parameters:

Name	Type	Description	Default
`raster`		Georeferenced `(y, x)` or `(band, y, x)` DataArray (e.g. from `load_raster()`); NaN and its NoData are no data.	required
`zones`		Polygons (GeoDataFrame, GeoSeries, shapely geometries), or a raster of zone labels on the same grid (integers; its NoData and NaN are no zone), e.g. municipalities or land-cover classes.	required
`stats`		Any of `"count"` (cells with data), `"nodata"` (cells without), `"sum"`, `"mean"`, `"min"`, `"max"`, `"range"`, `"std"`, `"var"` (population: over N), `"median"`, `"majority"` and `"minority"` (the most and least frequent value; ties: the smallest), `"unique"` (distinct values) and percentiles `"p10"`, `"p90"`... (linear interpolation, as NumPy's).	`('count', 'mean', 'min', 'max', 'std')`
`all_touched`	`bool`	Polygons: every cell they touch (default: the cells whose centre is inside, as `rasterize()`).	`False`
`n_threads`	`int \| None`	Threads (default: all cores but one); the result does not depend on it.	`None`

Returns:

Type	Description
`GeoDataFrame or DataFrame`	For polygons, the zones with one column per statistic (per band and statistic, `<band>_<stat>`, for a stack); for a raster of zones, one row per label (column `zone`). Zones without data: count 0, NaN.

Examples:

>>> import geopandas as gpd
>>> from shapely.geometry import box
>>> import landschaft as ls
>>> dem = ls.load_raster(ls.example_data("orelevation"))          # elevation (m)
>>> x0, y0, x1, y1 = dem.rio.bounds()
>>> halves = gpd.GeoDataFrame({"part": ["west", "east"]}, crs=dem.rio.crs,
...                           geometry=[box(x0, y0, (x0 + x1) / 2, y1), box((x0 + x1) / 2, y0, x1, y1)])
>>> ls.zonal_stats(dem, halves, ["mean", "max", "p90"]).round(1)[["part", "mean", "max", "p90"]]
   part   mean     max     p90
0  west  894.5  1252.0  1112.0
1  east  903.5  1390.0  1228.0

Zones from a raster: the cells of three elevation classes.

>>> classes = ls.reclassify(dem, [(0, 800, 1), (800, 1200, 2), (1200, 3000, 3)])
>>> ls.zonal_stats(dem, classes, ["count", "min", "max"])
   zone  count     min     max
0     1  46678   433.0   799.0
1     2  80653   800.0  1199.0
2     3  10310  1200.0  1390.0

extract ¶

extract(
    raster,
    points,
    *,
    method: str = "nearest",
    buffer: float | None = None,
    stats="mean",
    all_touched: bool = False,
    n_threads: int | None = None,
)

The values of a raster at points.

Habitat suitability at the occurrence records, the resistance at the sampling sites, every band of a stack at a million points: computed in C++ in parallel. The points are reprojected to the raster's CRS. A dask-backed raster (load_raster(..., chunks="auto")) is read only in the blocks of rows that hold points.

Parameters:

Name	Type	Description	Default
`raster`		Georeferenced `(y, x)` or `(band, y, x)` DataArray; NaN and its NoData are no data.	required
`points`		GeoDataFrame or GeoSeries of points, shapely points, an `(n, 2)` array or a pair of `x`, `y` arrays (in the raster's CRS).	required
`method`	`str`	`"nearest"` (default): the value of the cell holding the point; `"bilinear"`: interpolated from the four nearest cell centres (weights renormalised over those with data).	`'nearest'`
`buffer`	`float \| None`	A radius in the units of the raster's CRS: the statistics `stats` of the cells around every point (see `zonal_stats()`) instead of one value.	`None`
`stats`		With `buffer`: the statistics, as `zonal_stats()` names them.	`'mean'`
`all_touched`	`bool`	With `buffer`: every cell the circle touches.	`False`
`n_threads`	`int \| None`	Threads (default: all cores but one).	`None`

Returns:

Type	Description
`GeoDataFrame or DataFrame`	The points with one column per band (named by the band, or `value`), or per statistic with `buffer`; NaN outside the raster or on no data. Coordinates give a DataFrame with columns `x` and `y`.

Examples:

>>> import pandas as pd
>>> import landschaft as ls
>>> dem = ls.load_raster(ls.example_data("maungawhau"))       # elevation (m), 10 m cells
>>> table = pd.DataFrame({"site": ["summit", "slope"], "x": [2667705.0, 2667812.0], "y": [6478895.0, 6479004.0]})
>>> sites = ls.load_vector(ls.write_table(table, "sites.csv"), crs=dem.rio.crs)
>>> ls.extract(dem, sites)[["site", "maungawhau"]]
     site  maungawhau
0  summit       195.0
1   slope       168.0
>>> ls.extract(dem, sites, method="bilinear")["maungawhau"].round(1).tolist()
[195.0, 169.2]

The cells within 50 m of every site:

>>> ls.extract(dem, sites, buffer=50, stats=["mean", "max"])[["site", "mean", "max"]].round(1)
     site   mean    max
0  summit  186.5  195.0
1   slope  168.2  180.0

Tables¶

load_table ¶

load_table(
    source,
    *,
    sheet=None,
    sep: str | None = None,
    decimal: str | None = None,
    columns=None,
    encoding: str | None = None,
    **kwargs,
) -> DataFrame

Read a table, the format given by the file extension.

.csv / .txt / .tsv (also .gz; the separator , ; tab or |, a decimal comma and a Windows / Latin-1 encoding are detected), .xlsx / .xls (needs openpyxl), .parquet, .feather, .json, and the attribute table of a vector file (.gpkg, .shp, .dbf...).

Parameters:

Name	Type	Description	Default
`source`		File path.	required
`sheet`		Excel: the sheet, by name or number (default: the first).	`None`
`sep`	`str \| None`	Delimited text: the column separator and the decimal mark (default: detected).	`None`
`decimal`	`str \| None`	Delimited text: the column separator and the decimal mark (default: detected).	`None`
`columns`		Columns to read (default: all).	`None`
`encoding`	`str \| None`	Delimited text: the character encoding (default: UTF-8, or Windows-1252 / Latin-1 when the file is not valid UTF-8, as tables saved by Excel on Windows often are).	`None`
`**kwargs`		Passed to the pandas reader.	`{}`

Returns:

Type	Description
`DataFrame`

Examples:

>>> from pathlib import Path
>>> import landschaft as ls
>>> _ = Path("species.csv").write_text("especie;peso\nonca;2,5\nanta;1\n", encoding="utf-8")
>>> ls.load_table("species.csv")                       # ';' and the decimal comma detected
  especie  peso
0    onca   2.5
1    anta   1.0

write_table ¶

write_table(table, path, **kwargs) -> Path

Write a table, the format given by the file extension.

.csv (.csv.gz compressed), .tsv / .txt (tab-separated), .xlsx (needs openpyxl), .parquet (needs pyarrow; keeps the categorical and nullable integer columns, compact and fast for tables of millions of rows), .feather, .json (one record per row), or, for a GeoDataFrame, a vector format (.gpkg, .geojson, .shp, .fgb; see write_vector()). A named index (e.g. the metric names of a Series) is written as columns; missing folders are created.

Parameters:

Name	Description	Default
`table`	A DataFrame, Series or GeoDataFrame.	required
`path`	Output file; its extension selects the format.	required
`**kwargs`	Passed to the pandas / geopandas writer.	`{}`

Returns:

Type	Description
`Path`	The file written.

Examples:

>>> import landschaft as ls
>>> table = ls.landscape_metrics(ls.load_raster(ls.example_data("clc")), level="class", metric="ca")
>>> ls.write_table(table, "class_area.csv").name
'class_area.csv'
>>> ls.load_table("class_area.csv").columns.tolist()[:4]
['level', 'class', 'id', 'metric']

Saving results¶

save ¶

save(obj, path, *, like=None, table_format: str = 'csv', **kwargs) -> Path

Save data or a result, each part in the format that suits it.

One call for everything landschaft produces:

a raster (DataArray, numpy or dask array of 2 or 3 dimensions) becomes a GeoTIFF, as write_raster() writes it;
a table (DataFrame, Series) a CSV (or the format of the extension, see write_table()), a GeoDataFrame a GeoPackage (write_vector());
a tree a Newick file;
a result (of fos(), restoration_status(), zonation()... or a dict of results) a folder: one GeoTIFF per map, one table per table, and the numbers and settings in summary.json. Maps that are numpy arrays get the georeferencing of a georeferenced map of the same result, or of like. Results that have a format of their own (Zonation's output files, Biodiverse's) are written in it.

Without an extension, path gets the usual one (.tif, .csv, .gpkg) or is a folder for a result. Missing folders are created.

Parameters:

Name	Type	Description	Default
`obj`		What to save.	required
`path`		File, or folder for a result.	required
`like`		A raster (DataArray or path) that georeferences numpy maps.	`None`
`table_format`	`str`	Extension of the tables of a result: `"csv"` (default), `"parquet"`, `"xlsx"`...	`'csv'`
`**kwargs`		For a single raster, table or vector: passed to its writer (`nodata=`, `compress=`...).	`{}`

Returns:

Type	Description
`Path`	The file or folder written.

Examples:

>>> import landschaft as ls
>>> forest = ls.load_raster(ls.example_data("input"))
>>> ls.save(ls.mspa(forest), "results/mspa").name           # a raster: GeoTIFF
'mspa.tif'
>>> status = ls.restoration_status(forest)
>>> folder = ls.save(status, "results/restoration")         # a result: a folder
>>> sorted(p.name for p in folder.iterdir())
['summary.json']
>>> import json
>>> json.loads((folder / "summary.json").read_text())["nr_comp"]
2850

load_result ¶

load_result(path) -> dict

Read back a result folder written by save() (or result.save).

Every map comes back as a DataArray (load_raster()), every table as a DataFrame (load_table()), a vector layer as a GeoDataFrame, a .npy array as numpy, a sub-folder (a nested result) as a dict of its own, and the numbers and settings of summary.json as they were written ("result": the name of the result's class). The result's class itself is not rebuilt: its parts are.

Parameters:

Name	Type	Description	Default
`path`		The folder.	required

Returns:

Type	Description
`dict`	Each part by its name (the file name without its extension).

Examples:

>>> import landschaft as ls
>>> forest = ls.load_raster(ls.example_data("input"))
>>> folder = ls.fos(forest).save("fos_result")
>>> parts = ls.load_result(folder)
>>> parts["result"], type(parts["map"]).__name__, type(parts["stats"]).__name__
('FOSResult', 'DataArray', 'DataFrame')
>>> parts["n_patches"] == ls.fos(forest).n_patches
True

Every result object (of fos, restoration_status, zonation...) also has a save(folder) method that calls save.

Plotting¶

plot ¶

plot(
    obj,
    *,
    column=None,
    cmap=None,
    labels=None,
    categorical=None,
    vmin=None,
    vmax=None,
    log=None,
    opacity=None,
    bands=None,
    title=None,
    basemap="satellite",
    compare=False,
    static=False,
    max_size: int = 2048,
    show=None,
    **kwargs,
) -> Plot

Plot anything: a raster, a vector layer, a table or the result of a tool, as it calls for.

The form follows the data. A raster (any, or a landschaft product) goes on an interactive map over satellite imagery: a class map (integers, a product such as MSPA or FAD) in its classes' colours with their names, continuous values on a colour ramp that suits them (elevation, NDVI, priority ranks, currents on a log scale, values around zero on a diverging ramp), a stack band by band (or as true colour), a vector layer as polygons, lines or points, coloured by a column. A table becomes the chart it calls for; the result of a tool, its maps, charts and numbers together. Several rasters and vectors in a list share one map.

In a notebook the plot shows inline; elsewhere it opens in a window of its own (pywebview, else the web browser). Rasters are reduced to the screen's resolution before anything else (a lazy raster is read reduced), so a map of a whole country draws in about the time of a small one. Vector layers of over 20 000 features (up to 2 million) are drawn by the graphics card (WebGL), sharp at every zoom, the feature under the cursor read as for rasters. While Python runs (a notebook's kernel, a script's window), zooming in past the overview draws the part in view in full detail (LANDSCHAFT_PLOT_LIVE=0: off).

Parameters:

Name	Type	Description	Default
`obj`		A DataArray (2-D, or a `(band, y, x)` stack), GeoDataFrame, GeoSeries, shapely geometry, DataFrame, Series, numpy array, a result of a landschaft tool, a dict of these, or a list of rasters and vectors.	required
`column`		Vectors: the column that colours the features.	`None`
`cmap`		The colours: a palette name (`"viridis"`, `"magma"`, `"terrain"`, `"rdylgn"`, `"rdbu"`, `"blues"`, `"greens"`, `"ylorrd"`, `"rank"`, `"current"`, or any matplotlib colormap; `"_r"` reverses it), a list of colours, or a dict `{value: colour}` for classes.	`None`
`labels`		`{value: name}`: the names of the classes in the legend.	`None`
`categorical`		`True` draws the values as classes, `False` on a ramp (default: detected: products, booleans and integers with at most 64 values are classes).	`None`
`vmin`		The ramp's range (default: the 2nd and 98th percentiles) and a log scale.	`None`
`vmax`		The ramp's range (default: the 2nd and 98th percentiles) and a log scale.	`None`
`log`		The ramp's range (default: the 2nd and 98th percentiles) and a log scale.	`None`
`opacity`		Opacity of the layers.	`None`
`bands`		A stack: the bands to draw, by their labels (names, or the band numbers of a file read by `load_raster()`, from 1), else by position (from 0).	`None`
`title`		Title of the map or chart.	`None`
`basemap`		`"satellite"` (default, Esri World Imagery), `"osm"`, `"topo"`, `"light"`, `"relief"` or `None`; the others stay selectable.	`'satellite'`
`compare`		`obj` a list of two rasters or vector layers `[before, after]`: the first on the left of a line dragged across the map, the second on its right.	`False`
`static`		A static matplotlib figure instead (for publications and offline use; the satellite imagery needs contextily): a stack's bands, a result's rasters or two dates compared side by side as small multiples, with a scale bar, a north arrow and, on a single map, longitudes and latitudes. Its options, given as keywords: `figsize=(w, h)` (inches), `ncols=`, `scalebar=False`, `north_arrow=False`.	`False`
`max_size`	`int`	Pixels on the long side of every raster drawn (default 2048).	`2048`
`**kwargs`		Options of a tool's result: `env=` (a species distribution model: the predictor rasters, to map its suitability), `scale=` (a FAD change: the window size, default the multiscale map), `distance=` (a habitat graph: the links up to that distance); of a static figure (see `static`).	`{}`
`show`		Show the plot now (default: outside notebooks, unless the environment variable `LANDSCHAFT_PLOT_SHOW` is `0`); in a notebook, the returned Plot shows itself.	`None`

Returns:

Type	Description
`Plot`	Shows itself in a notebook; `.save("map.html")` writes it.

Examples:

>>> import landschaft as ls
>>> forest = ls.load_raster(ls.example_data("input"))
>>> p = ls.plot(ls.mspa(forest))                      # GuidosToolbox's colours and class names
>>> p
Plot(map)
>>> dem = ls.load_raster(ls.example_data("orelevation"))
>>> p = ls.plot(dem, cmap="terrain", title="Elevation")
>>> p.save("elevation.html").name
'elevation.html'
>>> after = forest.where(forest.x < forest.x[400], 1)      # the east cut
>>> p = ls.plot([forest, after], compare=True)            # drag the line to compare

Plot ¶

Plot(panels: list[_Panel], title: str | None = None, static_recipe=None)

A plot made by plot(): interactive maps and charts (or a static figure).

It shows itself inline in a notebook (Jupyter, JupyterLab, VS Code, Colab); elsewhere plot() opens it in a window of its own. save() writes a self-contained HTML page (or, for a static plot, PNG, PDF, SVG).

maps `property` ¶

maps: list

The folium maps (to customise them further).

charts `property` ¶

charts: list

The plotly figures.

figure `property` ¶

figure

The matplotlib figure of a static plot.

show ¶

show() -> None

Show the plot: inline in a notebook, else in a window of its own.

save ¶

save(path, *, dpi: int = 200) -> Path

Write the plot: .html (interactive, self-contained), or .png, .pdf, .svg, .jpg (a static rendering). Missing folders are created.

to_html ¶

to_html() -> str

A self-contained HTML page of the plot: it opens without a network connection (but for the satellite imagery and other base maps).

The page Plotting shows the plots of every kind of data, and how to choose their colours.

Example data¶

example_data ¶

example_data(name: str = 'input') -> Path

Path of an example data set shipped with the package.

One function for the examples of every module: open them with load_raster(), load_vector() or load_table().

Parameters:

Name	Type	Description	Default
`name`	`str`	Rasters of habitat and land cover (`EXAMPLE_DATA`): `"input"` (default): habitat map, 1000 x 1000 cells, coded 0 missing, 1 background, 2 forest; `"clc"`: CORINE land cover, 2000 x 2000 cells of 25 m (EPSG:3035), classes 1-44 (forests: 23, 24, 25); `"gscinput"`: percentage map (0-100, 255 missing), 1990 x 1289 cells. Surfaces (see `example_data()`): `"normforest"`, `"orforest"` (NDVI), `"orelevation"` (elevation), 371 x 371 cells. Grid distances (see `example_data()`): `"maungawhau"` (elevation of a volcano), `"europe"` (land and sea), `"gen_dist"` and `"pop_coord"` (tables of genetic distances and coordinates of European populations). Biodiversity (see `example_data()`): `"sites"` (species records, CSV) and `"tree"` (their phylogeny, Nexus). Species distributions (see `example_data()`): `"bradypus"` and `"acaule"` (occurrence records, CSV), `"bioclim"` (nine bioclimatic variables) and `"biome"` (biomes), 0.5 degree over Central and South America. Watershed model (see `load_txtinout()`): `"swat_hubbard_brook"` (the SWAT+ project of Hubbard Brook W2, a zip of its TxtInOut files, 2000-2002; its parameters.bsn, written for an earlier SWAT+, has the peak rate factor SWAT+ 62 expects, 484).	`'input'`

Returns:

Type	Description
`Path`	Path of the file.

Examples:

>>> import landschaft as ls
>>> ls.example_data("clc").name
'clc.tif'
>>> forest = ls.load_raster(ls.example_data("input"))
>>> dem = ls.load_raster(ls.example_data("orelevation"))
>>> sites = ls.load_table(ls.example_data("sites"))
>>> dem.shape, len(sites)
((371, 371), 1795)

Example data¶

ls.EXAMPLE_DATA describes the habitat and land-cover rasters; the examples of the other modules are listed under example_data above:

Name	Content
`"input"`	habitat map, 1000 × 1000 cells, 0 missing, 1 background, 2 forest (EPSG:4326)
`"clc"`	CORINE land cover excerpt, 2000 × 2000 cells of 25 m (EPSG:3035), classes 1–44
`"gscinput"`	percentage map, 1990 × 1289 cells of 2430 m (EPSG:5070), 0–100, 255 missing

Collections of rasters¶

RasterCollection ¶

RasterCollection(
    infos: Iterable[RasterInfo],
    *,
    band: int | None = None,
    validate: str | None = None,
    masked: bool | str = "auto",
    skipped: dict[Path, str] | None = None,
)

Bases: Sequence

Many rasters handled as one object; pixels are read only when processed.

Create it with load_raster() given a list of paths, a folder or a glob pattern. Pass it to mspa() or batch() to process every raster in parallel, or turn it into a lazy (file, y, x) stack with to_stack() when all rasters share a grid.

Attributes:

Name	Type	Description
`paths`	`list of pathlib.Path`	The files, in order.
`total_pixels`	`int`	Pixels of all rasters together.
`same_grid`	`bool`	Whether all rasters share shape, CRS and transform.
`skipped`	`dict`	Files left out by `errors="skip"`, with the reason.

Examples:

>>> import landschaft as ls
>>> forest = ls.load_raster(ls.example_data("input"))
>>> _ = ls.write_raster(forest[:500], "north.tif")
>>> _ = ls.write_raster(forest[500:], "south.tif")
>>> maps = ls.load_raster("*.tif")
>>> len(maps), maps.same_grid          # two different places: not the same grid
(2, False)
>>> first = maps.load(0)                       # one raster, as a DataArray
>>> results = ls.mspa(maps, edge_width=2)      # every raster, in parallel

load ¶

load(index: int) -> DataArray

Read one raster into memory (with the collection's checks).

Examples:

>>> import landschaft as ls
>>> forest = ls.load_raster(ls.example_data("input"))
>>> _ = ls.write_raster(forest[:500], "a.tif")
>>> _ = ls.write_raster(forest[500:], "b.tif")
>>> ls.load_raster(["a.tif", "b.tif"]).load(1).shape
(500, 1000)

to_stack ¶

to_stack() -> DataArray

Lazy (file, y, x) DataArray, one dask chunk per raster.

Requires a shared grid (see same_grid). Pixel values are not checked here; mspa() checks them when it runs.

Examples:

>>> import landschaft as ls
>>> forest = ls.load_raster(ls.example_data("input"))
>>> _ = ls.write_raster(forest, "y2010.tif")
>>> _ = ls.write_raster(forest, "y2020.tif")
>>> stack = ls.load_raster(["y2010.tif", "y2020.tif"]).to_stack()
>>> stack.dims
('file', 'y', 'x')

RasterInfo `dataclass` ¶

RasterInfo(
    path: Path,
    height: int,
    width: int,
    count: int,
    dtype: str,
    crs: Any,
    transform: Any,
    nodata: float | None,
)

Header of one raster of a collection (the pixels are not read).

Attributes:

Name	Type	Description
`path`	`Path`	The file.
`height, width, count`	`int`	Rows, columns and bands.
`dtype`	`str`	Data type of the first band.
`crs, transform`		Coordinate reference system and affine transform (rasterio objects).
`nodata`	`float or None`	NoData value.

Examples:

>>> import landschaft as ls
>>> info = ls.RasterInfo.read(ls.example_data("input"))
>>> info.height, info.width, info.pixels
(1000, 1000, 1000000)

Checking inputs¶

validate_mspa_input ¶

validate_mspa_input(data, *, spatial_only: bool = True) -> None

Check that data is a habitat map MSPA can process.

MSPA expects a single-band uint8 map coded 0 = missing, 1 = background, 2 = habitat, with at least one habitat pixel and at least one pixel that is not habitat.

Parameters:

Name	Type	Description	Default
`data`		numpy array, dask array or `xarray.DataArray`.	required
`spatial_only`	`bool`	Require a 2-D map (default); otherwise leading dimensions (e.g. time) are allowed and all slices are checked together.	`True`

Raises:

Type	Description
`InvalidInputError`	With a message that says what is wrong and how to fix it.

Examples:

>>> import landschaft as ls
>>> ls.validate_mspa_input(ls.load_raster(ls.example_data("input")))   # passes silently

validate_fad_input ¶

validate_fad_input(data) -> None

Check that data is a habitat map the density and object tools accept.

A single-band uint8 map coded 0 = missing, 1 = background, 2 = habitat and 3/4 = special background, with at least one habitat pixel (used by accounting() and fos()).

Parameters:

Name	Type	Description	Default
`data`		numpy array or `xarray.DataArray`, 2-D.	required

Raises:

Type	Description
`InvalidInputError`	With a message that says what is wrong and how to fix it.

Examples:

>>> import landschaft as ls
>>> ls.validate_fad_input(ls.load_raster(ls.example_data("input")))

InvalidInputError ¶

Bases: ValueError

The raster cannot be processed by the requested tool.

Raised by load_raster() and the habitat-pattern tools when a map is not coded as they expect; the message says what is wrong and how to fix it. It is a ValueError.

Examples:

>>> import numpy as np
>>> import landschaft as ls
>>> try:
...     ls.validate_mspa_input(np.array([[5, 1], [2, 2]], dtype="uint8"))
... except ls.InvalidInputError as err:
...     message = str(err)
>>> "larger than 2" in message
True

Colour tables¶

Every class map has a colour table: a {value: (red, green, blue, alpha)} dictionary with 256 entries, embedded by write_raster. Each is also available as a constant (ls.MSPA_CLASSES names the MSPA codes; ls.FAD_COLORMAP, ls.LM_COLORMAP... are the tables themselves).

mspa_colormap ¶

mspa_colormap(transition: bool = True) -> dict[int, tuple[int, int, int, int]]

Colour table of MSPA class maps (mspa()).

The standard colours of the product, as embedded by write_raster().

Parameters:

Name	Type	Description	Default
`transition`	`bool`	Colour bridges and loops that cross an edge or a perforation in their own colours (default) or as the class they cross.	`True`

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.mspa_colormap()
>>> len(colours), len(colours[0])
(256, 4)

fad_colormap ¶

fad_colormap() -> dict[int, tuple[int, int, int, int]]

Colour table of density maps of fad(), fad_multiscale() and fos() (6 classes).

The standard colours of the product, as embedded by write_raster().

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.fad_colormap()
>>> len(colours), len(colours[0])
(256, 4)

fad5_colormap ¶

fad5_colormap() -> dict[int, tuple[int, int, int, int]]

Colour table of 5-class fragmentation maps (fos() with classes=5).

The standard colours of the product, as embedded by write_raster().

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.fad5_colormap()
>>> len(colours), len(colours[0])
(256, 4)

lm_colormap ¶

lm_colormap() -> dict[int, tuple[int, int, int, int]]

Colour table of landscape mosaic maps (landscape_mosaic()).

The standard colours of the product, as embedded by write_raster().

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.lm_colormap()
>>> len(colours), len(colours[0])
(256, 4)

accounting_colormap ¶

accounting_colormap() -> dict[int, tuple[int, int, int, int]]

Colour table of object size-class maps (accounting()).

The standard colours of the product, as embedded by write_raster().

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.accounting_colormap()
>>> len(colours), len(colours[0])
(256, 4)

change_colormap ¶

change_colormap() -> dict[int, tuple[int, int, int, int]]

Colour table of change maps (simple_change(), morphological_change()).

The standard colours of the product, as embedded by write_raster().

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.change_colormap()
>>> len(colours), len(colours[0])
(256, 4)

foschange_colormap ¶

foschange_colormap() -> dict[int, tuple[int, int, int, int]]

Colour table of fragmentation change maps (fos_change()).

The standard colours of the product, as embedded by write_raster().

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.foschange_colormap()
>>> len(colours), len(colours[0])
(256, 4)

dist_colormap ¶

dist_colormap() -> dict[int, tuple[int, int, int, int]]

Colour table of distance maps (euclidean_distance()).

The standard colours of the product, as embedded by write_raster().

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.dist_colormap()
>>> len(colours), len(colours[0])
(256, 4)

entropy_colormap ¶

entropy_colormap() -> dict[int, tuple[int, int, int, int]]

Colour table of spatial entropy maps (spatial_entropy_map()).

The standard colours of the product, as embedded by write_raster().

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.entropy_colormap()
>>> len(colours), len(colours[0])
(256, 4)

contortion_colormap ¶

contortion_colormap() -> dict[int, tuple[int, int, int, int]]

Colour table of contortion maps (contortion()).

The standard colours of the product, as embedded by write_raster().

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.contortion_colormap()
>>> len(colours), len(colours[0])
(256, 4)

resistance_colormap ¶

resistance_colormap() -> dict[int, tuple[int, int, int, int]]

Colour table of resistance and isochrone maps (distance_to_resistance(), isochrone_map()).

The standard colours of the product, as embedded by write_raster().

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.resistance_colormap()
>>> len(colours), len(colours[0])
(256, 4)

gtn_components_colormap ¶

gtn_components_colormap() -> dict[int, tuple[int, int, int, int]]

Colour table of network component maps (gtn_components()).

The standard colours of the product, as embedded by write_raster().

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.gtn_components_colormap()
>>> len(colours), len(colours[0])
(256, 4)

gtn_importance_colormap ¶

gtn_importance_colormap() -> dict[int, tuple[int, int, int, int]]

Colour table of network importance maps (gtn_importance()).

The standard colours of the product, as embedded by write_raster().

Returns:

Type	Description
`dict`	`{value: (red, green, blue, alpha)}` for the 256 values 0-255.

Examples:

>>> import landschaft as ls
>>> colours = ls.gtn_importance_colormap()
>>> len(colours), len(colours[0])
(256, 4)

Threads¶

default_threads ¶

default_threads() -> int

The number of threads a parallel function uses when not given one.

All the cores the process may use but one (at least 1), so the computer stays responsive, or the value set with set_default_threads().

Returns:

Type	Description
`int`

Examples:

>>> import landschaft as ls
>>> n = ls.default_threads()
>>> n >= 1
True

set_default_threads ¶

set_default_threads(n: int | None = None) -> int

Set the default thread count of every parallel function for this process.

Parameters:

Name	Type	Description	Default
`n`	`int \| None`	Threads to use when a function is not given `n_threads`; `None` goes back to all cores but one.	`None`

Returns:

Type	Description
`int`	The new default.

Examples:

>>> import landschaft as ls
>>> ls.set_default_threads(2)
2
>>> result = ls.fad(ls.load_raster(ls.example_data("input")))   # 2 threads
>>> _ = ls.set_default_threads(None)                            # all cores but one again

max_threads ¶

max_threads() -> int

Number of processors the compiled core can use.

1 when landschaft was built without OpenMP. Parallel functions use default_threads() of them unless told otherwise.

Returns:

Type	Description
`int`

Examples:

>>> import landschaft as ls
>>> ls.max_threads() >= ls.default_threads() >= 1
True

Data & configuration¶

Rasters¶

load_raster ¶

write_raster ¶

match_grid ¶

reproject ¶

clip ¶

reclassify ¶

Vectors¶

load_vector ¶

write_vector ¶

rasterize ¶

Values by zone and at points¶

zonal_stats ¶

extract ¶

Tables¶

load_table ¶

write_table ¶

Saving results¶

save ¶

load_result ¶

Plotting¶

plot ¶

Plot ¶

maps property ¶

charts property ¶

figure property ¶

show ¶

save ¶

to_html ¶

Example data¶

example_data ¶

Example data¶

Collections of rasters¶

RasterCollection ¶

load ¶

to_stack ¶

RasterInfo dataclass ¶

Checking inputs¶

validate_mspa_input ¶

validate_fad_input ¶

InvalidInputError ¶

Colour tables¶

mspa_colormap ¶

fad_colormap ¶

fad5_colormap ¶

lm_colormap ¶

accounting_colormap ¶

change_colormap ¶

foschange_colormap ¶

dist_colormap ¶

entropy_colormap ¶

contortion_colormap ¶

resistance_colormap ¶

gtn_components_colormap ¶

gtn_importance_colormap ¶

Threads¶

default_threads ¶

set_default_threads ¶

max_threads ¶

maps `property` ¶

charts `property` ¶

figure `property` ¶

RasterInfo `dataclass` ¶