### Install on Windows

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/installation.md

Install the package via pip after setting up rasterio dependencies.

```bash
pip install rasterstats
```

--------------------------------

### Install on OS X

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/installation.md

Install GDAL via Homebrew before installing the python package.

```bash
brew install gdal
pip install rasterstats
```

--------------------------------

### Install on Ubuntu 14.04

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/installation.md

Install system dependencies via apt-get followed by the python package.

```bash
sudo apt-get install python-numpy libgdal1h gdal-bin libgdal-dev
pip install rasterstats
```

--------------------------------

### Calculate Zonal Statistics with Specific Stats

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/manual.md

Specify the statistics to calculate using the 'stats' argument as a list of strings. This example calculates min, max, median, majority, and sum.

```python
>>> stats = zonal_stats("tests/data/polygons.shp",
...                     "tests/data/slope.tif",
...                     stats=['min', 'max', 'median', 'majority', 'sum'])
```

--------------------------------

### Define and Use Custom Mean Statistic

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/manual.md

Define a custom aggregate function (e.g., 'mymean') using NumPy and register it with the 'add_stats' argument. This example calculates the count and the custom mean.

```python
>>> from __future__ import division
>>> import numpy as np

>>> def mymean(x):
...     return np.ma.mean(x)
```

```python
>>> zonal_stats("tests/data/polygons.shp",
...             "tests/data/slope.tif",
...             stats="count",
...             add_stats={'mymean':mymean})
```

--------------------------------

### Define and Use Custom Statistic with Geometry Properties

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/manual.md

Define a custom aggregate function that accepts geometry properties. Pass the desired properties using the 'properties' argument. This example calculates count and a custom mean multiplied by the 'id' property.

```python
>>> def mymean_prop(x,prop):
...     return np.ma.mean(x) * prop['id']
```

```python
>>> zonal_stats("tests/data/polygons.shp",
...             "tests/data/slope.tif",
...             stats="count",
...             add_stats={'mymean_prop':mymean_prop},
...             properties=['id'])
```

--------------------------------

### Point Query Python API

Source: https://github.com/perrygeo/python-rasterstats/blob/master/README.rst

Perform point queries using the `point_query` function in Python. Pass a GeoJSON-like point dictionary and the raster file path to get the raster value at that point.

```python
from rasterstats import point_query
point = {'type': 'Point', 'coordinates': (245309.0, 1000064.0)}
point_query(point, "tests/data/slope.tif")
```

--------------------------------

### Specify Raster Data Sources

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/manual.md

Verify raster support and load data from files, specific bands, or numpy arrays.

```bash
$ rio info raster.tif
```

```python
>>> zs = zonal_stats('tests/data/polygons.shp', 'tests/data/slope.tif')
```

```python
>>> zs = zonal_stats('tests/data/polygons.shp', 'tests/data/slope.tif', band=1)
```

```python
>>> import rasterio
>>> with rasterio.open('tests/data/slope.tif') as src:
...     affine = src.transform
...     array = src.read(1)
>>> zs = zonal_stats('tests/data/polygons.shp', array, affine=affine)
```

--------------------------------

### Execute point queries via CLI

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Use the rio pointquery plugin to extract raster values at point locations with configurable interpolation and nodata handling.

```bash
# Basic point query with bilinear interpolation (default)
fio cat sample_points.shp | rio pointquery -r elevation.tif > points_with_elevation.geojson

# Use nearest neighbor interpolation
fio cat points.shp | rio pointquery -r elevation.tif \
    --interpolate nearest > output.geojson

# Specify property name for the queried value
fio cat weather_stations.shp | rio pointquery -r temperature.tif \
    --property-name "temperature" > stations_with_temp.geojson

# Query specific band
fio cat points.shp | rio pointquery -r multiband.tif \
    --band 2 > output.geojson

# Override nodata value
fio cat points.shp | rio pointquery -r elevation.tif \
    --nodata -9999 > output.geojson
```

--------------------------------

### Point Query CLI

Source: https://github.com/perrygeo/python-rasterstats/blob/master/README.rst

Utilize the command line interface for point queries by piping point vector data to the `rio pointquery` subcommand.

```bash
$ fio cat points.shp | rio pointquery -r elevation.tif
```

--------------------------------

### Saving Output to File

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/cli.md

Redirect the command output to a file to save the resulting GeoJSON.

```default
fio cat countries.shp | rio zonalstats -r dem.tif --prefix "elevation_" > countries_with_elevation.geojson
```

--------------------------------

### Specifying Custom Statistics

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/cli.md

Use the --stats option to specify a space-delimited string of desired summary statistics.

```default
$ fio cat countries.shp \
    | rio zonalstats -r dem.tif \
          --prefix "elevation_" \
          --stats "min max median percentile_95"
...
{
  "type": "Feature",
  "geometry": {...} ,
  "properties": {
    "country_name": "Grenada",
    "elevation_min": 0.0,
    "elevation_median": 161.33
    "elevation_max": 840.33,
    "elevation_percentile_95": 533.6
  }
}
```

--------------------------------

### Execute zonal statistics via CLI

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Use the rio zonalstats plugin to process vector data from stdin and output GeoJSON statistics.

```bash
# Basic zonal statistics
fio cat polygons.shp | rio zonalstats -r elevation.tif > output.geojson

# Specify statistics to calculate
fio cat countries.shp | rio zonalstats -r dem.tif \
    --stats "min max mean median percentile_95"

# Add prefix to property names
fio cat regions.shp | rio zonalstats -r temperature.tif \
    --prefix "temp_" > regions_with_temp.geojson

# Include all cells touched by geometry
fio cat small_parcels.shp | rio zonalstats -r coarse_dem.tif \
    --all-touched > parcels_with_elevation.geojson

# Categorical raster
fio cat watersheds.shp | rio zonalstats -r landcover.tif \
    --categorical > watersheds_landcover.geojson

# Specify raster band
fio cat polygons.shp | rio zonalstats -r multiband.tif \
    --band 3 > output.geojson

# Stream processing with sequences (memory efficient for large files)
fio cat large.shp | rio zonalstats -r elevation.tif \
    --sequence | some-other-process
```

--------------------------------

### Use memory-efficient point query generators

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Stream results for large point datasets to avoid loading entire files into memory.

```python
from rasterstats import gen_point_query

# Stream results for large point datasets
for value in gen_point_query("millions_of_points.shp", "elevation.tif"):
    if value is not None:
        # Process each value
        pass
```

--------------------------------

### Enable 'all_touched' Rasterization Strategy

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/manual.md

Enable the 'all_touched' rasterization strategy by setting 'all_touched=True'. This includes all pixels that the geometry touches, as opposed to the default which only includes pixels whose center point is within the geometry.

```python
>>> zs = zonal_stats("tests/data/polygons.shp",
...                  "tests/data/slope.tif",
...                  all_touched=True)
```

--------------------------------

### Basic Zonal Statistics with File Paths

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Calculate default summary statistics (min, max, mean, count) for raster values within vector polygons using file paths. Ensure polygons and raster files are accessible.

```python
from rasterstats import zonal_stats
from pprint import pprint

# Basic usage with file paths
stats = zonal_stats("polygons.shp", "elevation.tif")
pprint(stats)
# [{'count': 75, 'max': 22.27, 'mean': 14.66, 'min': 6.58},
#  {'count': 50, 'max': 82.69, 'mean': 56.61, 'min': 16.94}]
```

--------------------------------

### Customizing Output Prefixes

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/cli.md

Use the --prefix option to customize the property names in the output GeoJSON.

```default
$ fio cat countries.shp | rio zonalstats -r dem.tif --prefix "elevation_"
...
{
  "type": "Feature",
  "geometry": {...} ,
  "properties": {
    "country_name": "Grenada",
    "elevation_min": 0.0,
    "elevation_mean": 210.47,
    "elevation_max": 840.33,
    "elevation_count": 94
  }
}
```

--------------------------------

### Specify Vector Data Sources

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/manual.md

Load vector data from files, layers, or objects supporting the __geo_interface__.

```python
>>> zs = zonal_stats('tests/data/polygons.shp', 'tests/data/slope.tif')
```

```python
>>> zs = zonal_stats('tests/data', 'tests/data/slope.tif', layer="polygons")
```

```python
>>> import fiona
>>> with fiona.open('tests/data/polygons.shp') as src:
...    zs = zonal_stats(src, 'tests/data/slope.tif')
```

```python
>>> from shapely.geometry import Point
>>> pt = Point(245000, 1000000)
>>> pt.__geo_interface__
{'type': 'Point', 'coordinates': (245000.0, 1000000.0)}
>>> point_query([pt], 'tests/data/slope.tif')
[21.32739672330894]
```

```python
>>> pt.wkt
'POINT (245000 1000000)'
>>> point_query([pt], 'tests/data/slope.tif')
[21.32739672330894]

>>> pt.wkb
'\x01\x01\x00\x00\x00\x00\x00\x00\x00@\xe8\rA\x00\x00\x00\x00\x80\x84.A'
>>> point_query([pt], 'tests/data/slope.tif')
[21.32739672330894]
```

--------------------------------

### Point Query with Interpolation

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Choose between bilinear interpolation for smoother values or nearest neighbor for exact cell values.

```python
from rasterstats import point_query
from shapely.geometry import Point

point = Point(245309.5, 1000064.5)

# Bilinear interpolation (default) - weighted average of 4 nearest cells
value_bilinear = point_query(
    [point],
    "elevation.tif",
    interpolate="bilinear"
)
print(f"Bilinear: {value_bilinear[0]:.2f}")
# Bilinear: 156.78

# Nearest neighbor - value of closest cell
value_nearest = point_query(
    [point],
    "elevation.tif",
    interpolate="nearest"
)
print(f"Nearest: {value_nearest[0]:.2f}")
# Nearest: 158.00
```

--------------------------------

### Streaming Large Datasets

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/cli.md

Use --sequence to stream line-delimited features for memory-efficient processing of large datasets.

```default
fio cat large.shp | rio zonalstats -r elevation.tif --sequence | some-other-process
```

--------------------------------

### Migrate GeoRaster Output

Source: https://github.com/perrygeo/python-rasterstats/wiki/Migrating-to-1.0

Manual implementation of GeoRaster output following the removal of the GeoRaster dependency.

```python
# before 
zs = zonal_stats(..., raster_out=True, opt_georaster=True)
```

```python
# after
zs = zonal_stats(...)
newstats = []
for s in zs:
    s['georaster'] = georasters.GeoRaster(
        s['mini_raster_array'],
        s['mini_raster_affine'].to_gdal(),
        nodata_value=s['mini_raster_nodata'],
        projection=spatial_ref)
    newstats.append(s)
```

--------------------------------

### Zonal Statistics CLI

Source: https://github.com/perrygeo/python-rasterstats/blob/master/README.rst

Use the command line interface to calculate zonal statistics by piping vector data to the `rio zonalstats` subcommand.

```bash
$ fio cat polygon.shp | rio zonalstats -r elevation.tif
```

--------------------------------

### Replace global_src_extent with Rasterio

Source: https://github.com/perrygeo/python-rasterstats/wiki/Migrating-to-1.0

Manually read raster data using rasterio instead of using the deprecated global_src_extent argument.

```python
# before
zonal_stats(raster="a.tif", ..., global_src_extent=True)
```

```python
# after
with rasterio.open("a.tif) as src:
    arr = src.read(1)
    aff = src.affine
    nodata = src.nodata
zonal_stats(raster=arr, affine=aff, nodata=nodata)
```

--------------------------------

### Extract mini-rasters for features

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/manual.md

Set raster_out=True to include the clipped and masked numpy array, affine transformation, and nodata value in the output dictionary.

```python
>>> zonal_stats('tests/data/polygons.shp',
...             'tests/data/slope_classes.tif',
...             stats="count",
...             raster_out=True)[0].keys()
['count', 'mini_raster_affine', 'mini_raster_array', 'mini_raster_nodata']
```

--------------------------------

### Control Rasterization with all_touched

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Adjust the all_touched parameter to include cells that touch the geometry boundary, which is useful for small polygons.

```python
from rasterstats import zonal_stats

# Default: only cells with centroid inside polygon
stats_default = zonal_stats(
    "small_parcels.shp",
    "coarse_dem.tif",
    all_touched=False  # default
)

# Include all cells that touch the geometry boundary
stats_touched = zonal_stats(
    "small_parcels.shp",
    "coarse_dem.tif",
    all_touched=True
)

# Useful for small polygons that might miss cells with default rasterization
print(f"Default count: {stats_default[0]['count']}")
# Default count: 2
print(f"All touched count: {stats_touched[0]['count']}")
# All touched count: 8
```

--------------------------------

### Perform Basic Zonal Statistics and Point Queries

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/manual.md

Calculate zonal statistics for polygons and query raster values at specific points.

```python
>>> from rasterstats import zonal_stats, point_query
>>> stats = zonal_stats('tests/data/polygons.shp', 'tests/data/slope.tif')
>>> pts = point_query('tests/data/points.shp', 'tests/data/slope.tif')
```

```python
>>> from pprint import pprint
>>> pprint(stats)
[{'count': 75,
  'max': 22.273418426513672,
  'mean': 14.660084635416666,
  'min': 6.575114727020264},
 {'count': 50,
  'max': 82.69043731689453,
  'mean': 56.60576171875,
  'min': 16.940950393676758}]
```

```python
>>> pts
[14.037668283186257, 33.1370268256543, 36.46848854950241]
```

--------------------------------

### Reference Available Zonal Statistics

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

List of all supported statistics for zonal_stats, including percentile support and wildcard usage.

```python
from rasterstats import zonal_stats

# Default statistics (returned if stats not specified)
# count, min, max, mean

# All available statistics
stats = zonal_stats(
    "polygons.shp",
    "raster.tif",
    stats=[
        "count",       # Number of valid (non-nodata) pixels
        "min",         # Minimum value
        "max",         # Maximum value
        "mean",        # Arithmetic mean
        "sum",         # Sum of all values
        "std",         # Standard deviation
        "median",      # Median value
        "majority",    # Most common value (mode)
        "minority",    # Least common value
        "unique",      # Count of unique values
        "range",       # max - min
        "nodata",      # Count of nodata pixels
        "percentile_10",   # 10th percentile
        "percentile_50",   # 50th percentile (same as median)
        "percentile_90",   # 90th percentile
        "percentile_99.5", # Decimal percentiles supported
    ]
)

# Use "*" or "ALL" to get all standard statistics
stats = zonal_stats("polygons.shp", "raster.tif", stats="*")
```

--------------------------------

### Point Query with GeoJSON Output

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Preserve feature properties and geometry while appending queried raster values as a new property.

```python
from rasterstats import point_query

# Return GeoJSON features with raster values as properties
features = point_query(
    "weather_stations.shp",
    "temperature.tif",
    geojson_out=True,
    property_name="temperature"
)

print(features[0])
```

--------------------------------

### Use memory-efficient zonal statistics generators

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Process large datasets feature-by-feature using generators to minimize memory usage. Includes options for progress tracking.

```python
from rasterstats import gen_zonal_stats

# Process features one at a time without loading all results into memory
for stat in gen_zonal_stats("large_dataset.shp", "raster.tif"):
    if stat["mean"] is not None and stat["mean"] > 100:
        print(f"High mean value: {stat['mean']}")
    # Process or write each result individually

# Combine with progress bar (requires tqdm)
stats = zonal_stats(
    "large_dataset.shp",
    "raster.tif",
    progress=True  # Shows tqdm progress bar
)
```

--------------------------------

### Pipe GeoJSON to Zonal Statistics

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/cli.md

Pipe GeoJSON features from `fio cat` directly into `rio zonalstats` to calculate statistics without an intermediate file. The raster is specified using the -r option.

```default
fio cat countries.shp | rio zonalstats -r dem.tif
```

--------------------------------

### CLI: rio zonalstats

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/cli.md

Generates summary statistics of geospatial raster datasets based on vector features.

```APIDOC
## CLI: rio zonalstats

### Description
Generates summary statistics (min, max, mean, etc.) of geospatial raster datasets based on input vector features. The output is a GeoJSON with additional properties per feature.

### Parameters
#### Path Parameters
- **FEATURES** (GeoJSON) - Required - Valid GeoJSON Features to process.

#### Options
- **-r, --raster** (PATH) - Required - File path to the raster dataset.
- **--all-touched** (boolean) - Optional - Include all pixels touched by the geometry.
- **--band** (INTEGER) - Optional - Raster band to process.
- **--categorical** (boolean) - Optional - Treat raster as categorical data.
- **--stats** (TEXT) - Optional - Specific statistics to calculate.

### Request Example
rio zonalstats states.geojson -r rainfall.tif > mean_rainfall_by_state.geojson
```

--------------------------------

### Query Raster Values at Points

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Extract raster values at specific point locations using file paths, Shapely geometries, or GeoJSON dictionaries.

```python
from rasterstats import point_query
from shapely.geometry import Point

# Basic usage with file path
values = point_query("sample_points.shp", "elevation.tif")
print(values)
# [145.2, 203.8, 89.4, 167.1]

# Using Shapely geometries
points = [
    Point(245309.0, 1000064.0),
    Point(245500.0, 1000100.0),
    Point(245100.0, 999900.0)
]
values = point_query(points, "elevation.tif")
print(values)
# [74.09, 82.34, 65.21]

# GeoJSON-like dictionary
point_geojson = {"type": "Point", "coordinates": (245309.0, 1000064.0)}
values = point_query([point_geojson], "slope.tif")
print(values)
# [74.09817594635244]
```

--------------------------------

### Handle diverse vector data input formats

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Rasterstats supports various input types including file paths, fiona sources, GeoJSON, Shapely geometries, WKT, and WKB.

```python
import fiona
from shapely.geometry import shape, mapping, Point, Polygon
from rasterstats import zonal_stats, point_query

# 1. File path (any format supported by fiona)
stats = zonal_stats("data/polygons.shp", "elevation.tif")
stats = zonal_stats("data/regions.geojson", "elevation.tif")
stats = zonal_stats("data/zones.gpkg", "elevation.tif", layer="boundaries")

# 2. Open fiona source
with fiona.open("polygons.shp") as src:
    stats = zonal_stats(src, "elevation.tif")

# 3. GeoJSON-like dictionary
geojson_feature = {
    "type": "Feature",
    "geometry": {
        "type": "Polygon",
        "coordinates": [[[0, 0], [1, 0], [1, 1], [0, 1], [0, 0]]]
    },
    "properties": {"id": 1}
}
stats = zonal_stats(geojson_feature, "elevation.tif")

# 4. GeoJSON FeatureCollection
fc = {
    "type": "FeatureCollection",
    "features": [geojson_feature]
}
stats = zonal_stats(fc, "elevation.tif")

# 5. Shapely geometry
polygon = Polygon([(0, 0), (1, 0), (1, 1), (0, 1)])
stats = zonal_stats([polygon], "elevation.tif")

# 6. WKT string
wkt_geom = "POLYGON ((0 0, 1 0, 1 1, 0 1, 0 0))"
stats = zonal_stats([wkt_geom], "elevation.tif")

# 7. WKB bytes
wkb_geom = polygon.wkb
stats = zonal_stats([wkb_geom], "elevation.tif")
```

--------------------------------

### Zonal Statistics with ArcPy

Source: https://github.com/perrygeo/python-rasterstats/wiki/Alternative-raster-stats-implementations

Use ArcPy's ZonalStatisticsAsTable for calculating statistics. Outputs to DBF format and does not provide accessible Python variables.

```python
import arcpy
arcpy.sa.ZonalStatisticsAsTable("fcid", "FID", "dem_aea2_feet.tif", "E:\\workspace\\stats1.dbf")
```

--------------------------------

### Update Geotransform to Affine

Source: https://github.com/perrygeo/python-rasterstats/wiki/Migrating-to-1.0

Replace GDAL-style geotransform tuples with Affine objects in zonal_stats calls.

```python
# before 
gt = (244300.6, 25.5, 0.0, 1000868.8, 0.0, -25.5)
zonal_stats(..., transform=gt)
```

```python
# after
zonal_stats(..., affine=Affine.from_gdal(*gt))
```

--------------------------------

### Extract Mini-Raster Output

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Use raster_out=True to retrieve the clipped and masked raster array for each feature.

```python
from rasterstats import zonal_stats
import matplotlib.pyplot as plt

stats = zonal_stats(
    "polygons.shp",
    "elevation.tif",
    stats="count mean",
    raster_out=True
)

# Access the mini-raster for first feature
feature_data = stats[0]
mini_array = feature_data["mini_raster_array"]      # Masked numpy array
mini_affine = feature_data["mini_raster_affine"]    # Affine transform
mini_nodata = feature_data["mini_raster_nodata"]    # NoData value

print(f"Array shape: {mini_array.shape}")
print(f"Valid pixels: {mini_array.count()}")
print(f"Affine: {mini_affine}")

# Visualize the clipped raster
plt.imshow(mini_array)
plt.colorbar(label="Elevation (m)")
plt.title("Elevation within polygon")
plt.savefig("mini_raster.png")
```

--------------------------------

### Output Zonal Statistics as GeoJSON Features

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/manual.md

Use the 'geojson_out=True' argument to return a list of GeoJSON features, where each feature includes the calculated zonal statistics as properties. This preserves the original geometries and properties.

```python
>>> stats = zonal_stats("tests/data/polygons.shp",
...                     "tests/data/slope.tif",
...                     geojson_out=True)
```

--------------------------------

### Calculate Zonal Statistics

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/index.md

Compute summary statistics for elevation data within polygon boundaries.

```python
from rasterstats import zonal_stats
zonal_stats("polygons.shp", "elevation.tif",
            stats="count min mean max median")
```

--------------------------------

### Point Query CLI Usage

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/cli.md

Use `rio pointquery` to query raster values at the points of input GeoJSON Features. Raster values are added to feature properties and output as GeoJSON Feature Collection. Supports bilinear (default) or nearest neighbor interpolation.

```console
$ rio pointquery --help
Usage: rio pointquery [OPTIONS] FEATURES...

Queries the raster values at the points of the input GeoJSON Features. The
raster values are added to the features properties and output as GeoJSON
Feature Collection.

If the Features are Points, the point geometery is used. For other
Feauture types, all of the verticies of the geometry will be queried. For
example, you can provide a linestring and get the profile along the line
if the verticies are spaced properly.

You can use either bilinear (default) or nearest neighbor interpolation.

Options:
--version                   Show the version and exit.
-r, --raster PATH           [required]
--band INTEGER
--nodata INTEGER
--indent INTEGER
--interpolate TEXT
--property-name TEXT
--sequence / --no-sequence  Write a LF-delimited sequence of texts
                            containing individual objects or write a single
                            JSON text containing a feature collection object
                            (the default).
--rs / --no-rs              Use RS (0x1E) as a prefix for individual texts
                            in a sequence as per http://tools.ietf.org/html
                            /draft-ietf-json-text-sequence-13 (default is
                            False).
-h, --help                  Show this message and exit.
```

--------------------------------

### Point Query with Non-Point Geometries

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Extract values at all vertices for complex geometry types like LineStrings or Polygons.

```python
from rasterstats import point_query
from shapely.geometry import LineString

# Create a transect line
transect = LineString([
    (245000, 1000000),
    (245100, 1000050),
    (245200, 1000100),
    (245300, 1000150)
])
```

--------------------------------

### Zonal Statistics with QGIS

Source: https://github.com/perrygeo/python-rasterstats/wiki/Alternative-raster-stats-implementations

QGIS's QgsZonalStatistics can calculate statistics but is slower and offers limited stats. Variables are not directly accessible, and it may alter the source dataset.

```python
from qgis.analysis import QgsZonalStatistics

#specify polygon shapefile vector
polygonLayer = QgsVectorLayer('F:/temp/zonalstat/zonePoly.shp', 'zonepolygons', "ogr") 

# specify raster filename
rasterFilePath = 'F:/temp/zonalstat/raster1.tif'

# usage - QgsZonalStatistics (QgsVectorLayer *polygonLayer, const QString &rasterFile, const QString &attributePrefix="", int rasterBand=1)
zoneStat = QgsZonalStatistics (polygonLayer, rasterFilePath, 'pre-', 1)
zoneStat.calculateStatistics(None)
```

--------------------------------

### Zonal Statistics with Numpy Arrays

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Perform zonal statistics using numpy arrays directly. Requires providing the affine transform mapping array coordinates to geographic coordinates. Ensure rasterio is used to read the data and obtain the transform.

```python
import rasterio
import numpy as np
from rasterstats import zonal_stats

# Read raster data with rasterio
with rasterio.open("elevation.tif") as src:
    affine = src.transform
    array = src.read(1)
    nodata = src.nodata

# Pass numpy array with affine transform
stats = zonal_stats(
    "polygons.shp",
    array,
    affine=affine,
    nodata=nodata
)
print(stats)
# [{'count': 75, 'max': 22.27, 'mean': 14.66, 'min': 6.58}]
```

--------------------------------

### Convert Shapefile to GeoJSON Features

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/cli.md

Use `fio cat` to convert a shapefile into GeoJSON Features printed to standard output. This is useful for piping data into other commands.

```default
fio cat countries.shp
```

--------------------------------

### CLI: rio pointquery

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/cli.md

Queries raster values at the points of input GeoJSON features and adds them to feature properties.

```APIDOC
## CLI: rio pointquery

### Description
Queries the raster values at the points of the input GeoJSON Features. If features are not points, all vertices of the geometry are queried.

### Parameters
#### Path Parameters
- **FEATURES** (GeoJSON) - Required - Valid GeoJSON Features to process.

#### Options
- **-r, --raster** (PATH) - Required - File path to the raster dataset.
- **--band** (INTEGER) - Optional - Raster band to process.
- **--interpolate** (TEXT) - Optional - Interpolation method (bilinear or nearest neighbor).
- **--property-name** (TEXT) - Optional - Name for the property to store the queried value.
```

--------------------------------

### Perform point queries on vector data

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Extract raster values at point locations or along transects. Returns a list of values for each feature.

```python
profile = point_query([transect], "elevation.tif")
print(profile[0])  # List of values at each vertex
# [145.2, 167.8, 189.3, 212.1]

# Works with polygons too - returns values at all boundary vertices
polygon_values = point_query("boundaries.shp", "elevation.tif")
# Each feature returns a list of values at its vertices
```

--------------------------------

### Zonal Statistics CLI Usage

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/cli.md

Use `rio zonalstats` to generate summary statistics of geospatial raster datasets based on vector features. Input features should be valid GeoJSON. The raster is specified with the -r/--raster argument.

```console
Usage: rio zonalstats [OPTIONS] FEATURES...

zonalstats generates summary statistics of geospatial raster datasets
based on vector features.

The input arguments to zonalstats should be valid GeoJSON Features. (see
cligj)

The output GeoJSON will be mostly unchanged but have additional properties
per feature describing the summary statistics (min, max, mean, etc.) of
the underlying raster dataset.

The raster is specified by the required -r/--raster argument.

Example, calculate rainfall stats for each state and output to file:

  rio zonalstats states.geojson -r rainfall.tif > mean_rainfall_by_state.geojson

Options:
--version                       Show the version and exit.
-r, --raster PATH               [required]
--all-touched / --no-all-touched
--band INTEGER
--categorical / --no-categorical
--indent INTEGER
--info / --no-info
--nodata INTEGER
--prefix TEXT
--stats TEXT
--sequence / --no-sequence      Write a LF-delimited sequence of texts
                                containing individual objects or write a
                                single JSON text containing a feature
                                collection object (the default).
--rs / --no-rs                  Use RS (0x1E) as a prefix for individual
                                texts in a sequence as per
                                http://tools.ietf.org/html/draft-ietf-json-
                                text-sequence-13 (default is False).
-h, --help                      Show this message and exit.
```

--------------------------------

### Specify Zonal Statistics

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Calculate specific summary statistics by providing a list or space-delimited string of desired statistics. Supports min, max, median, majority, sum, std, range, and percentiles.

```python
# Specify which statistics to calculate
stats = zonal_stats(
    "polygons.shp",
    "elevation.tif",
    stats=["min", "max", "median", "majority", "sum", "std", "range", "percentile_95"]
)
```

```python
# Alternative: space-delimited string
stats = zonal_stats(
    "polygons.shp",
    "slope.tif",
    stats="min max median sum percentile_90"
)
```

```python
# Use all available statistics
stats = zonal_stats("polygons.shp", "elevation.tif", stats="*")
```

--------------------------------

### Calculate Zonal Statistics with Space-Delimited Stats String

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/manual.md

Alternatively, specify statistics as a space-delimited string. This method is equivalent to using a list of strings for the 'stats' argument.

```python
>>> stats = zonal_stats("tests/data/polygons.shp",
...                     "tests/data/slope.tif",
...                     stats="min max median majority sum")
```

--------------------------------

### Zonal Statistics with Categorical Rasters

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Process categorical rasters by setting `categorical=True`. This counts pixels for each unique category value. Optionally, use `category_map` to map numerical category values to meaningful names.

```python
from rasterstats import zonal_stats

# Count pixels by category
stats = zonal_stats(
    "watersheds.shp",
    "landcover.tif",
    categorical=True
)
print(stats[0])
# {1.0: 542, 2.0: 1203, 3.0: 89, 5.0: 2341}
```

```python
# Map category values to meaningful names
category_map = {
    1.0: "water",
    2.0: "forest",
    3.0: "urban",
    5.0: "agriculture"
}

stats = zonal_stats(
    "watersheds.shp",
    "landcover.tif",
    categorical=True,
    category_map=category_map
)
print(stats[0])
# {'water': 542, 'forest': 1203, 'urban': 89, 'agriculture': 2341}
```

--------------------------------

### Zonal Statistics Output Format

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/index.md

The function returns a list of dictionaries containing the calculated statistics for each feature.

```python
[...,
 {'count': 89,
  'max': 69.52958679199219,
  'mean': 20.08093536034059,
  'median': 19.33736801147461,
  'min': 1.5106816291809082},
]
```

--------------------------------

### Zonal Statistics with GeoJSON Output

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Generate GeoJSON output where zonal statistics are appended as properties to the original features. Use `geojson_out=True` and optionally provide a `prefix` for the new property names.

```python
from rasterstats import zonal_stats

# Return GeoJSON features with stats as properties
stats = zonal_stats(
    "polygons.shp",
    "elevation.tif",
    geojson_out=True,
    prefix="elev_"  # prefix for stat property names
)

# Each result is a GeoJSON-like feature
print(stats[0]["type"])
# 'Feature'

print(stats[0]["properties"])
# {'id': 1, 'name': 'Zone A', 'elev_min': 6.58, 'elev_max': 22.27,
#  'elev_mean': 14.66, 'elev_count': 75}

print(stats[0]["geometry"]["type"])
# 'Polygon'
```

--------------------------------

### Calculate categorical zonal statistics

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/manual.md

Use categorical=True to return a dictionary of pixel value counts. Map these values to labels using the category_map argument.

```python
>>> zonal_stats('tests/data/polygons.shp',
...             'tests/data/slope_classes.tif',
...             categorical=True)[1]
{1.0: 1, 2.0: 9, 5.0: 40}
```

```python
>>> cmap = {1.0: 'low', 2.0: 'med', 5.0: 'high'}
>>> zonal_stats('tests/data/polygons.shp',
...             'tests/data/slope_classes.tif',
...             categorical=True, category_map=cmap)[1]
{'high': 40, 'med': 9, 'low': 1}
```

--------------------------------

### Zonal Statistics with Custom Statistics Functions

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Define and use custom aggregate functions for zonal statistics via the `add_stats` argument. Functions can accept the masked numpy array and optionally feature properties.

```python
import numpy as np
from rasterstats import zonal_stats

# Single-argument custom function
def coefficient_of_variation(x):
    return np.ma.std(x) / np.ma.mean(x) * 100

# Two-argument function with access to feature properties
def weighted_mean(x, props):
    return np.ma.mean(x) * props.get("weight", 1.0)

```

--------------------------------

### Zonal Statistics with R

Source: https://github.com/perrygeo/python-rasterstats/wiki/Alternative-raster-stats-implementations

This R script uses the raster and rgdal libraries to extract raster data from polygons and append results to a DBF file.

```r
library(raster)
library(rgdal)
library(foreign)

# extract raster data from overlying polygon
rast <- raster('raster.tif')
poly <- readOGR(dsn="polygons", layer="polygon")
ext.poly <- extract(rast, poly, fun = mean, na.rm=TRUE, df=TRUE)

# and then append the result to the .dbf file of the original polygon
poly.dbf <- read.dbf("polygons/polygon.dbf", as.is = TRUE)
poly.dbf$result <- ext.poly$extraction
write.dbf(poly.dbf, "polygon_with_extraction.dbf")
```

--------------------------------

### Integrate Rasterstats with GeoPandas

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Use GeoPandas to read vector data and join calculated zonal statistics back to the original GeoDataFrame.

```python
import geopandas as gpd
import pandas as pd
from rasterstats import zonal_stats

# Read vector data with GeoPandas
gdf = gpd.read_file("regions.shp")

# Calculate zonal statistics
stats = zonal_stats(gdf, "elevation.tif", stats="min mean max std")

# Convert stats to DataFrame and join with original GeoDataFrame
stats_df = pd.DataFrame(stats)
result = gdf.join(stats_df)

# Now you have a GeoDataFrame with both geometry and statistics
print(result[["name", "min", "mean", "max", "std"]].head())
#          name    min    mean     max    std
# 0    Region A  123.4  456.7   789.0   45.2
# 1    Region B   98.2  234.5   567.8   67.3

# Save to file
result.to_file("regions_with_stats.gpkg", driver="GPKG")
```

--------------------------------

### GeoJSON Feature Output Structure

Source: https://github.com/perrygeo/python-rasterstats/blob/master/docs/cli.md

The default structure of a GeoJSON feature after zonal statistics processing.

```default
{
  "type": "Feature",
  "geometry": {...} ,
  "properties": {
    "country_name": "Grenada",
    "_min": 0.0,
    "_mean": 210.47,
    "_max": 840.33,
    "_count": 94
  }
}
```

--------------------------------

### Custom Zonal Statistics Function

Source: https://context7.com/perrygeo/python-rasterstats/llms.txt

Define a custom function to calculate statistics based on the rasterized geometry mask.

```python
def edge_mean(x, props, rv_array):
    # Calculate mean only at the edges of the geometry
    from scipy import ndimage
    eroded = ndimage.binary_erosion(rv_array)
    edge_mask = rv_array & ~eroded
    return float(np.ma.mean(x[edge_mask])) if edge_mask.any() else None

stats = zonal_stats(
    "polygons.shp",
    "elevation.tif",
    stats="count mean",
    add_stats={
        "cv": coefficient_of_variation,
        "weighted_mean": weighted_mean
    }
)
print(stats[0])
```

--------------------------------

### Zonal Statistics Python API

Source: https://github.com/perrygeo/python-rasterstats/blob/master/README.rst

Calculate zonal statistics in Python by providing a shapefile and a raster file to the `zonal_stats` function. The function returns a list of dictionaries, each containing statistics for a feature.

```python
from rasterstats import zonal_stats
stats = zonal_stats("tests/data/polygons.shp", "tests/data/slope.tif")
stats[0].keys()
[f['mean'] for f in stats]
```

=== COMPLETE CONTENT === This response contains all available snippets from this library. No additional content exists. Do not make further requests.