Latest Release		Downloads
License		People
Build Status		Coverage

---

📍 Fast, Accurate Python library for Raster Operations

⚡ Extensible with Numba

⏩ Scalable with Dask

🎊 Free of GDAL / GEOS Dependencies

🌍 General-Purpose Spatial Processing, Geared Towards GIS Professionals

---

Xarray-Spatial implements common raster analysis functions using Numba and provides an easy-to-install, easy-to-extend codebase for raster analysis.

Installation

# via pip
pip install xarray-spatial

# via conda
conda install -c conda-forge xarray-spatial

Downloading our starter examples and data

Once you have xarray-spatial installed in your environment, you can use one of the following in your terminal (with the environment active) to download our examples and/or sample data into your local directory.

xrspatial examples : Download the examples notebooks and the data used.

xrspatial copy-examples : Download the examples notebooks but not the data. Note: you won't be able to run many of the examples.

xrspatial fetch-data : Download just the data and not the notebooks.

In all the above, the command will download and store the files into your current directory inside a folder named 'xrspatial-examples'.

xarray-spatial grew out of the Datashader project, which provides fast rasterization of vector data (points, lines, polygons, meshes, and rasters) for use with xarray-spatial.

xarray-spatial does not depend on GDAL / GEOS, which makes it fully extensible in Python but does limit the breadth of operations that can be covered. xarray-spatial is meant to include the core raster-analysis functions needed for GIS developers / analysts, implemented independently of the non-Python geo stack.

Our documentation is still under construction, but docs can be found here.

Raster-huh?

Rasters are regularly gridded datasets like GeoTIFFs, JPGs, and PNGs.

In the GIS world, rasters are used for representing continuous phenomena (e.g. elevation, rainfall, distance), either directly as numerical values, or as RGB images created for humans to view. Rasters typically have two spatial dimensions, but may have any number of other dimensions (time, type of measurement, etc.)

Supported Spatial Functions with Supported Inputs

✅ = native backend    🔄 = accepted (CPU fallback)

---

Classification

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
Box Plot	Classifies values into bins based on box plot quartile boundaries	PySAL mapclassify	✅️	✅	✅	🔄
Equal Interval	Divides the value range into equal-width bins	PySAL mapclassify	✅️	✅	✅	✅
Head/Tail Breaks	Classifies heavy-tailed distributions using recursive mean splitting	PySAL mapclassify	✅️	✅	🔄	🔄
Maximum Breaks	Finds natural groupings by maximizing differences between sorted values	PySAL mapclassify	✅️	✅	🔄	🔄
Natural Breaks	Optimizes class boundaries to minimize within-class variance (Jenks)	Jenks 1967, PySAL	✅️	✅	🔄	🔄
Percentiles	Assigns classes based on user-defined percentile breakpoints	PySAL mapclassify	✅️	✅	✅	🔄
Quantile	Distributes values into classes with equal observation counts	PySAL mapclassify	✅️	✅	✅	🔄
Reclassify	Remaps pixel values to new classes using a user-defined lookup	PySAL mapclassify	✅️	✅	✅	✅
Std Mean	Classifies values by standard deviation intervals from the mean	PySAL mapclassify	✅️	✅	✅	✅

---

Diffusion

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
Diffuse	Runs explicit forward-Euler diffusion on a 2D scalar field	Standard (heat equation)	✅️	✅️	✅️	✅️

---

Focal

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
Apply	Applies a custom function over a sliding neighborhood window	Standard	✅️	✅️	✅️	✅️
Hotspots	Identifies statistically significant spatial clusters using Getis-Ord Gi*	Getis & Ord 1992	✅️	✅️	✅️	✅️
Emerging Hotspots	Classifies time-series hot/cold spot trends using Gi* and Mann-Kendall	Getis & Ord 1992, Mann 1945	✅️	✅️	✅️	✅️
Mean	Computes the mean value within a sliding neighborhood window	Standard	✅️	✅️	✅️	✅️
Focal Statistics	Computes summary statistics over a sliding neighborhood window	Standard	✅️	✅️	✅️	✅️
Bilateral	Feature-preserving smoothing via bilateral filtering	Tomasi & Manduchi 1998	✅️	✅️	✅️	✅️
GLCM Texture	Computes Haralick GLCM texture features over a sliding window	Haralick et al. 1973	✅️	✅️	🔄	🔄

---

Morphological

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
Erode	Morphological erosion (local minimum over structuring element)	Standard (morphology)	✅️	✅️	✅️	✅️
Dilate	Morphological dilation (local maximum over structuring element)	Standard (morphology)	✅️	✅️	✅️	✅️
Opening	Erosion then dilation (removes small bright features)	Standard (morphology)	✅️	✅️	✅️	✅️
Closing	Dilation then erosion (fills small dark gaps)	Standard (morphology)	✅️	✅️	✅️	✅️

---

Fire

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
dNBR	Differenced Normalized Burn Ratio (pre minus post NBR)	USGS	✅️	✅️	✅️	✅️
RdNBR	Relative dNBR normalized by pre-fire vegetation density	USGS	✅️	✅️	✅️	✅️
Burn Severity Class	USGS 7-class burn severity from dNBR thresholds	USGS	✅️	✅️	✅️	✅️
Fireline Intensity	Byram's fireline intensity from fuel load and spread rate (kW/m)	Byram 1959	✅️	✅️	✅️	✅️
Flame Length	Flame length derived from fireline intensity (m)	Byram 1959	✅️	✅️	✅️	✅️
Rate of Spread	Simplified Rothermel spread rate with Anderson 13 fuel models (m/min)	Rothermel 1972, Anderson 1982	✅️	✅️	✅️	✅️
KBDI	Keetch-Byram Drought Index single time-step update (0-800 mm)	Keetch & Byram 1968	✅️	✅️	✅️	✅️

---

Multispectral

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
Atmospherically Resistant Vegetation Index (ARVI)	Vegetation index resistant to atmospheric effects using blue band correction	Kaufman & Tanre 1992	✅️	✅️	✅️	✅️
Enhanced Built-Up and Bareness Index (EBBI)	Highlights built-up areas and barren land from thermal and SWIR bands	As-syakur et al. 2012	✅️	✅️	✅️	✅️
Enhanced Vegetation Index (EVI)	Enhanced vegetation index reducing soil and atmospheric noise	Huete et al. 2002	✅️	✅️	✅️	✅️
Green Chlorophyll Index (GCI)	Estimates leaf chlorophyll content from green and NIR reflectance	Gitelson et al. 2003	✅️	✅️	✅️	✅️
Normalized Burn Ratio (NBR)	Measures burn severity using NIR and SWIR band difference	USGS Landsat	✅️	✅️	✅️	✅️
Normalized Burn Ratio 2 (NBR2)	Refines burn severity mapping using two SWIR bands	USGS Landsat	✅️	✅️	✅️	✅️
Normalized Difference Moisture Index (NDMI)	Detects vegetation moisture stress from NIR and SWIR reflectance	USGS Landsat	✅️	✅️	✅️	✅️
Normalized Difference Water Index (NDWI)	Maps open water bodies using green and NIR band difference	McFeeters 1996	✅️	✅️	✅️	✅️
Modified Normalized Difference Water Index (MNDWI)	Detects water in urban areas using green and SWIR bands	Xu 2006	✅️	✅️	✅️	✅️
Normalized Difference Vegetation Index (NDVI)	Quantifies vegetation density from red and NIR band difference	Rouse et al. 1973	✅️	✅️	✅️	✅️
Soil Adjusted Vegetation Index (SAVI)	Vegetation index with soil brightness correction factor	Huete 1988	✅️	✅️	✅️	✅️
Structure Insensitive Pigment Index (SIPI)	Estimates carotenoid-to-chlorophyll ratio for plant stress detection	Penuelas et al. 1995	✅️	✅️	✅️	✅️
True Color	Composites red, green, and blue bands into a natural color image	Standard	✅️	✅️	✅️	✅️

---

Multivariate

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
Mahalanobis Distance	Measures statistical distance from a multi-band reference distribution, accounting for band correlations	Mahalanobis 1936	✅️	✅️	✅️	✅️

---

Pathfinding

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
A* Pathfinding	Finds the least-cost path between two cells on a cost surface	Hart et al. 1968	✅️	✅	🔄	🔄
Multi-Stop Search	Routes through N waypoints in sequence, with optional TSP reordering	Custom	✅️	✅	🔄	🔄

---

Proximity

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
Allocation	Assigns each cell to the identity of the nearest source feature	Standard (Dijkstra)	✅️	✅	✅️	✅️
Balanced Allocation	Partitions a cost surface into territories of roughly equal cost-weighted area	Custom	✅️	✅	✅️	✅️
Cost Distance	Computes minimum accumulated cost to the nearest source through a friction surface	Standard (Dijkstra)	✅️	✅	✅️	✅️
Least-Cost Corridor	Identifies zones of low cumulative cost between two source locations	Standard (Dijkstra)	✅️	✅	✅️	✅️
Direction	Computes the direction from each cell to the nearest source feature	Standard	✅️	✅	✅️	✅️
Proximity	Computes the distance from each cell to the nearest source feature	Standard	✅️	✅	✅️	✅️
Surface Distance	Computes distance along the 3D terrain surface to the nearest source	Standard (Dijkstra)	✅️	✅	✅️	✅️
Surface Allocation	Assigns each cell to the nearest source by terrain surface distance	Standard (Dijkstra)	✅️	✅	✅️	✅️
Surface Direction	Computes compass direction to the nearest source by terrain surface distance	Standard (Dijkstra)	✅️	✅	✅️	✅️

---

Raster to vector

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
Polygonize	Converts contiguous regions of equal value into vector polygons	Standard (CCL)	✅️	✅️	✅️	🔄
Contours	Extracts elevation contour lines (isolines) from a raster surface	Standard (marching squares)	✅️	✅️	🔄	🔄

---

Surface

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
Aspect	Computes downslope direction of each cell in degrees	Horn 1981	✅️	✅️	✅️	✅️
Curvature	Measures rate of slope change (concavity/convexity) at each cell	Zevenbergen & Thorne 1987	✅️	✅️	✅️	✅️
Hillshade	Simulates terrain illumination from a given sun angle and azimuth	GDAL gdaldem	✅️	✅️	✅️	✅️
Roughness	Computes local relief as max minus min elevation in a 3×3 window	GDAL gdaldem	✅️	✅️	✅️	✅️
Sky-View Factor	Measures the fraction of visible sky hemisphere at each cell	Zakek et al. 2011	✅️	✅️	✅️	✅️
Slope	Computes terrain gradient steepness at each cell in degrees	Horn 1981	✅️	✅️	✅️	✅️
Terrain Generation	Generates synthetic terrain from fBm or ridged fractal noise with optional domain warping, Worley blending, and hydraulic erosion	Custom (fBm)	✅️	✅️	✅️	✅️
TPI	Computes Topographic Position Index (center minus mean of neighbors)	Weiss 2001	✅️	✅️	✅️	✅️
TRI	Computes Terrain Ruggedness Index (local elevation variation)	Riley et al. 1999	✅️	✅️	✅️	✅️
Landforms	Classifies terrain into 10 landform types using the Weiss (2001) TPI scheme	Weiss 2001	✅️	✅️	✅️	✅️
Viewshed	Determines visible cells from a given observer point on terrain	GRASS GIS r.viewshed	✅️	✅️	✅️	✅️
Min Observable Height	Finds the minimum observer height needed to see each cell (experimental)	Custom	✅️
Perlin Noise	Generates smooth continuous random noise for procedural textures	Perlin 1985	✅️	✅️	✅️	✅️
Worley Noise	Generates cellular (Voronoi) noise returning distance to the nearest feature point	Worley 1996	✅️	✅️	✅️	✅️
Hydraulic Erosion	Simulates particle-based water erosion to carve valleys and deposit sediment	Custom	✅️	✅️	✅️	✅️
Bump Mapping	Adds randomized bump features to simulate natural terrain variation	Custom	✅️	✅️	✅️	✅️

---

Hydrology

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
Flow Direction (D8)	Computes D8 flow direction from each cell toward the steepest downhill neighbor	O'Callaghan & Mark 1984	✅️	✅️	✅️	✅️
Flow Direction (Dinf)	Computes D-infinity flow direction as a continuous angle toward the steepest downslope facet	Tarboton 1997	✅️	✅️	✅️	✅️
Flow Direction (MFD)	Partitions flow to all downslope neighbors with an adaptive exponent (Qin et al. 2007)	Qin et al. 2007	✅️	✅️	✅️	✅️
Flow Accumulation (D8)	Counts upstream cells draining through each cell in a D8 flow direction grid	Jenson & Domingue 1988	✅️	✅️	✅️	✅️
Flow Accumulation (MFD)	Accumulates upstream area through all MFD flow paths weighted by directional fractions	Qin et al. 2007	✅️	✅️	✅️	🔄
Watershed	Labels each cell with the pour point it drains to via D8 flow direction	Standard (D8 tracing)	✅️	✅️	✅️	✅️
Basins	Delineates drainage basins by labeling each cell with its outlet ID	Standard (D8 tracing)	✅️	✅️	✅️	✅️
Stream Order	Assigns Strahler or Shreve stream order to cells in a drainage network	Strahler 1957, Shreve 1966	✅️	✅️	✅️	✅️
Stream Link	Assigns unique IDs to each stream segment between junctions	Standard	✅️	✅️	✅️	✅️
Snap Pour Point	Snaps pour points to the highest-accumulation cell within a search radius	Custom	✅️	✅️	✅️	✅️
Flow Path	Traces downstream flow paths from start points through a D8 direction grid	Standard (D8 tracing)	✅️	✅️	🔄	🔄
HAND	Computes Height Above Nearest Drainage by tracing D8 flow to the nearest stream cell	Nobre et al. 2011	✅️	✅️	🔄	🔄
TWI	Topographic Wetness Index: ln(specific catchment area / tan(slope))	Beven & Kirkby 1979	✅️	✅️	✅️	🔄

---

Flood

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
Flood Depth	Computes water depth above terrain from a HAND raster and water level	Standard (HAND-based)	✅️	✅️	✅️	✅️
Inundation	Produces a binary flood/no-flood mask from a HAND raster and water level	Standard (HAND-based)	✅️	✅️	✅️	✅️
Curve Number Runoff	Estimates runoff depth from rainfall using the SCS/NRCS curve number method	SCS/NRCS	✅️	✅️	✅️	✅️
Travel Time	Estimates overland flow travel time via simplified Manning's equation	Manning 1891	✅️	✅️	✅️	✅️

---

Interpolation

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
IDW	Inverse Distance Weighting from scattered points to a raster grid	Standard (IDW)	✅️	✅️	✅️	✅️
Kriging	Ordinary Kriging with automatic variogram fitting (spherical, exponential, gaussian)	Standard (ordinary kriging)	✅️	✅️	✅️	✅️
Spline	Thin Plate Spline interpolation with optional smoothing	Standard (TPS)	✅️	✅️	✅️	✅️

---

Zonal

Name	Description	Source	NumPy xr.DataArray	Dask xr.DataArray	CuPy GPU xr.DataArray	Dask GPU xr.DataArray
Apply	Applies a custom function to each zone in a classified raster	Standard	✅️	✅️	✅️	✅️
Crop	Extracts the bounding rectangle of a specific zone	Standard	✅️	✅️	✅️	✅️
Regions	Identifies connected regions of non-zero cells	Standard (CCL)	✅️	✅️	✅️	✅️
Trim	Removes nodata border rows and columns from a raster	Standard	✅️	✅️	✅️	✅️
Zonal Statistics	Computes summary statistics for a value raster within each zone	Standard	✅️	✅️	✅️	🔄
Zonal Cross Tabulate	Cross-tabulates agreement between two categorical rasters	Standard	✅️	✅️	🔄	🔄

Usage

Quick Start

Importing xrspatial registers an .xrs accessor on DataArrays and Datasets, giving you tab-completable access to every spatial operation:

import numpy as np
import xarray as xr
import xrspatial

# Create or load a raster
elevation = xr.DataArray(np.random.rand(100, 100) * 1000, dims=['y', 'x'])

# Surface analysis — call operations directly on the DataArray
slope = elevation.xrs.slope()
hillshaded = elevation.xrs.hillshade(azimuth=315, angle_altitude=45)
aspect = elevation.xrs.aspect()

# Classification
classes = elevation.xrs.equal_interval(k=5)
breaks = elevation.xrs.natural_breaks(k=10)

# Proximity
distance = elevation.xrs.proximity(target_values=[1])

# Multispectral — call on the NIR band, pass other bands as arguments
nir = xr.DataArray(np.random.rand(100, 100), dims=['y', 'x'])
red = xr.DataArray(np.random.rand(100, 100), dims=['y', 'x'])
blue = xr.DataArray(np.random.rand(100, 100), dims=['y', 'x'])

vegetation = nir.xrs.ndvi(red)
enhanced_vi = nir.xrs.evi(red, blue)

Dataset Support

The .xrs accessor works on Datasets too. Single-input functions apply the operation to each data variable. Multi-input functions (multispectral indices) accept string kwargs that map band aliases to variable names:

ds = xr.Dataset({'band_4': red, 'band_5': nir})

# Single-input: slope computed for each variable
slope_ds = ds.xrs.slope()

# Multi-input: map variable names to band parameters
ndvi_result = ds.xrs.ndvi(nir='band_5', red='band_4')

Function Import Style

All operations are also available as standalone functions if you prefer explicit imports:

from xrspatial import hillshade, slope, ndvi

hillshaded = hillshade(elevation)
slope_result = slope(elevation)
vegetation = ndvi(nir, red)

Check out the user guide here.

---

Dependencies

xarray-spatial currently depends on Datashader, but will soon be updated to depend only on xarray and numba, while still being able to make use of Datashader output when available.

Notes on GDAL

Within the Python ecosystem, many geospatial libraries interface with the GDAL C++ library for raster and vector input, output, and analysis (e.g. rasterio, rasterstats, geopandas). GDAL is robust, performant, and has decades of great work behind it. For years, off-loading expensive computations to the C/C++ level in this way has been a key performance strategy for Python libraries (obviously...Python itself is implemented in C!).

However, wrapping GDAL has a few drawbacks for Python developers and data scientists:

• GDAL can be a pain to build / install.
• GDAL is hard for Python developers/analysts to extend, because it requires understanding multiple languages.
• GDAL's data structures are defined at the C/C++ level, which constrains how they can be accessed from Python.

With the introduction of projects like Numba, Python gained new ways to provide high-performance code directly in Python, without depending on or being constrained by separate C/C++ extensions. xarray-spatial implements algorithms using Numba and Dask, making all of its source code available as pure Python without any "black box" barriers that obscure what is going on and prevent full optimization. Projects can make use of the functionality provided by xarray-spatial where available, while still using GDAL where required for other tasks.

Citation

Cite this code:

xarray-contrib/xarray-spatial, https://github.com/xarray-contrib/xarray-spatial, ©2020-2026.