Raster Features API

obia.features.extract.extract_features(segments, image, *, bands=None, spectral=SPECTRAL_FEATURES, percentiles=(10, 90), indices=None, shape=(), shape_units='map', texture=None, mask=None)

Add selected features to a copy of an existing segment table.

Parameters:

Name	Type	Description	Default
segments	GeoDataFrame	Polygon GeoDataFrame in the raster CRS. Its IDs, index, geometry, and existing attributes are preserved. Requested feature columns replace columns with the same names when recalculating.	required
image	Image	OBIA Image; calculations use its loaded (height, width, bands) values, including their selected order, rather than rereading files.	required
bands	tuple[int, ...] | None	Zero-based band indices for spectral statistics and, by default, texture. None selects all loaded bands; () selects none.	None
spectral	tuple[str, ...]	Names in SPECTRAL_FEATURES. Use () to disable these summaries.	SPECTRAL_FEATURES
percentiles	tuple[float, ...]	Additional percentiles in [0, 100], using linear interpolation. Use () to disable them. Applied to bands and calculated indices.	(10, 90)
indices	dict[str, tuple[int, ...]] | None	Index names mapped to ordered numeric band indices, e.g. {'ndvi': (3, 0)}. Indices are computed per pixel before summarization.	None
shape	tuple[str, ...] | bool	Names in SHAPE_FEATURES, True for all, or False/() for none.	()
shape_units	str	'map' requires a projected CRS and uses its coordinate units; 'pixel' uses the inverse raster transform, including rotation/shear.	'map'
texture	TextureConfig | None	Optional texture configuration. Ranges are fitted once per band over finite included image pixels, or supplied explicitly.	None
mask	ndarray | None	Optional boolean (height, width) mask. False pixels contribute to neither distributions nor texture. Geometry itself is unchanged.	None

Returns:

Type	Description
GeoDataFrame	gpd.GeoDataFrame: Feature columns prefixed bN_, index_METHOD_bN_..., or shape_. Undefined measurements are NaN, valid counts are zero for objects without samples. Configuration and fitted texture ranges are recorded in attrs['obia_features'] (in-memory metadata).

Source code in obia/features/extract.py

def extract_features(
    segments: gpd.GeoDataFrame, image: "Image", *, bands: tuple[int, ...] | None = None,
    spectral: tuple[str, ...] = SPECTRAL_FEATURES, percentiles: tuple[float, ...] = (10, 90),
    indices: dict[str, tuple[int, ...]] | None = None,
    shape: tuple[str, ...] | bool = (), shape_units: str = "map",
    texture: TextureConfig | None = None, mask: np.ndarray | None = None,
) -> gpd.GeoDataFrame:
    """Add selected features to a copy of an existing segment table.

    Args:
        segments: Polygon GeoDataFrame in the raster CRS. Its IDs, index,
            geometry, and existing attributes are preserved. Requested feature
            columns replace columns with the same names when recalculating.
        image: OBIA Image; calculations use its loaded (height, width, bands)
            values, including their selected order, rather than rereading files.
        bands: Zero-based band indices for spectral statistics and, by default,
            texture. None selects all loaded bands; () selects none.
        spectral: Names in ``SPECTRAL_FEATURES``. Use () to disable these summaries.
        percentiles: Additional percentiles in [0, 100], using linear interpolation.
            Use () to disable them. Applied to bands and calculated indices.
        indices: Index names mapped to ordered numeric band indices, e.g.
            {'ndvi': (3, 0)}. Indices are computed per pixel before summarization.
        shape: Names in ``SHAPE_FEATURES``, True for all, or False/() for none.
        shape_units: 'map' requires a projected CRS and uses its coordinate units;
            'pixel' uses the inverse raster transform, including rotation/shear.
        texture: Optional texture configuration. Ranges are fitted once per band
            over finite included image pixels, or supplied explicitly.
        mask: Optional boolean (height, width) mask. False pixels contribute to
            neither distributions nor texture. Geometry itself is unchanged.

    Returns:
        gpd.GeoDataFrame: Feature columns prefixed bN_, index_METHOD_bN_..., or
            shape_. Undefined measurements are NaN, valid counts are zero for
            objects without samples. Configuration and fitted texture ranges are
            recorded in attrs['obia_features'] (in-memory metadata).
    """
    data = np.ma.asarray(image.img_data)
    if data.ndim != 3 or any(size == 0 for size in data.shape):
        raise ValueError("image.img_data must have nonempty shape (height, width, bands).")
    if not isinstance(segments, gpd.GeoDataFrame):
        raise TypeError("segments must be a GeoDataFrame.")
    if (segments.crs is None) != (image.crs is None) or (
        segments.crs is not None and CRS.from_user_input(segments.crs) != CRS.from_user_input(image.crs)
    ):
        raise ValueError("segments and image must have the same CRS; reproject segments first.")
    for geometry in segments.geometry:
        if geometry is None or geometry.geom_type not in {"Polygon", "MultiPolygon"} or not geometry.is_valid:
            raise ValueError("segments must contain valid Polygon or MultiPolygon geometries.")
    selected = band_indices(bands, data.shape[2])
    spectral, percentiles = validate_statistics(spectral, percentiles)
    if shape is True:
        shape = SHAPE_FEATURES
    elif shape is False:
        shape = ()
    shape = choices(shape, SHAPE_FEATURES, "shape features")
    if shape_units not in {"map", "pixel"}:
        raise ValueError("shape_units must be 'map' or 'pixel'.")
    if shape and shape_units == "map" and (image.crs is None or not CRS.from_user_input(image.crs).is_projected):
        raise ValueError("Map-unit shape measurements require a projected CRS; use shape_units='pixel' instead.")
    included = np.ones(data.shape[:2], dtype=bool) if mask is None else np.asarray(mask, dtype=bool)
    if included.shape != data.shape[:2]:
        raise ValueError("mask must have shape (height, width).")
    image_mask = getattr(image, "valid_mask", None)
    if image_mask is not None and np.shape(image_mask) != data.shape:
        raise ValueError("image.valid_mask must match image.img_data.")
    if texture is not None and not isinstance(texture, TextureConfig):
        raise TypeError("texture must be a TextureConfig or None.")
    texture_bands = () if texture is None else band_indices(
        selected if texture.bands is None else texture.bands, data.shape[2])
    indices = {} if indices is None else dict(indices)
    indices = {name: validate_index(name, values, data.shape[2]) for name, values in indices.items()}

    cached_bands = {}

    def load_band(index):
        if index not in cached_bands:
            values = np.ma.asarray(data[:, :, index], dtype=float).filled(np.nan)
            valid = included & np.isfinite(values)
            if image_mask is not None:
                valid &= np.asarray(image_mask[:, :, index], dtype=bool)
            cached_bands[index] = np.where(valid, values, np.nan)
        return cached_bands[index]

    layers = {}
    if spectral or percentiles:
        layers.update({f"b{b}": load_band(b) for b in selected})
        for name, inputs in indices.items():
            values = np.stack([load_band(b) for b in inputs], axis=-1)
            prefix = f"index_{name}_" + "_".join(f"b{b}" for b in inputs)
            layers[prefix] = compute_index(values, name, tuple(range(len(inputs))))
    ranges = {}
    for band in texture_bands:
        values = load_band(band)
        valid = values[np.isfinite(values)]
        ranges[band] = texture.value_range or (
            (float(valid.min()), float(valid.max())) if valid.size else (0.0, 0.0))

    stats_columns = tuple(summarize([], 0, spectral, percentiles))
    columns = [f"{prefix}_{stat}" for prefix in layers for stat in stats_columns]
    columns += [f"b{b}_{name}" for b in texture_bands for name in texture_columns(texture)]
    columns += [f"shape_{name}" for name in shape]
    grid = SimpleNamespace(transform=image.transform, width=data.shape[1], height=data.shape[0])
    rows = []
    for geometry in segments.geometry:
        row = {}
        window = _object_window(geometry, grid) if layers or texture_bands else None
        if window is None:
            for prefix in layers:
                row.update({f"{prefix}_{name}": value for name, value in
                            summarize([], 0, spectral, percentiles).items()})
            for band in texture_bands:
                row.update({f"b{band}_{name}": np.nan for name in texture_columns(texture)})
        else:
            r, c, inside = window
            total = int(inside.sum())
            for prefix, values in layers.items():
                summary = summarize(values[r, c][inside], total, spectral, percentiles)
                row.update({f"{prefix}_{name}": value for name, value in summary.items()})
            for band in texture_bands:
                measured = measure_texture(load_band(band)[r, c], inside, texture, ranges[band])
                row.update({f"b{band}_{name}": value for name, value in measured.items()})
        if shape:
            measured = shape_features(geometry, shape, image.transform if shape_units == "pixel" else None)
            row.update({f"shape_{name}": value for name, value in measured.items()})
        rows.append(row)
    result = segments.copy()
    features = pd.DataFrame(rows, index=segments.index, columns=columns)
    # Assign together to avoid fragmenting tables with many bands and scales.
    result = gpd.GeoDataFrame(pd.concat([result.drop(columns=columns, errors="ignore"), features], axis=1),
                             geometry=segments.geometry.name, crs=segments.crs)
    result.attrs = dict(segments.attrs)
    result.attrs["obia_features"] = {
        "bands": selected, "spectral": spectral, "percentiles": percentiles,
        "indices": indices, "shape": shape, "shape_units": shape_units,
        "texture": asdict(texture) if texture is not None else None,
        "texture_ranges": ranges,
    }
    return result

obia.features.spectral.compute_index(data, method, bands)

Calculate a spectral index from explicit zero-based band indices.

Parameters:

Name	Type	Description	Default
data	ndarray	Raster array shaped (height, width, bands), optionally masked.	required
method	str	Name in SPECTRAL_INDICES. Normalized-difference methods use (a-b)/(a+b), ratio uses a/b, and difference uses a-b. exg uses (2*b-a-c)/(a+b+c), vari uses (b-a)/(b+a-c), and brightness is the arithmetic mean. Here a, b, c are the bands in supplied order.	required
bands	tuple[int, ...]	Ordered numeric band indices. No wavelengths are inferred. For example, ndvi with (3, 0) calculates (b3-b0)/(b3+b0).	required

Returns:

Type	Description
ndarray	np.ndarray: A float64 index image. Masked/nonfinite operands, zero denominators, and nonfinite results become NaN. Input is unchanged.

Source code in obia/features/spectral.py

def compute_index(data: np.ndarray, method: str, bands: tuple[int, ...]) -> np.ndarray:
    """Calculate a spectral index from explicit zero-based band indices.

    Args:
        data: Raster array shaped (height, width, bands), optionally masked.
        method: Name in ``SPECTRAL_INDICES``. Normalized-difference methods use
            (a-b)/(a+b), ratio uses a/b, and difference uses a-b. ``exg`` uses
            (2*b-a-c)/(a+b+c), ``vari`` uses (b-a)/(b+a-c), and brightness is
            the arithmetic mean. Here a, b, c are the bands in supplied order.
        bands: Ordered numeric band indices. No wavelengths are inferred.
            For example, ndvi with (3, 0) calculates (b3-b0)/(b3+b0).

    Returns:
        np.ndarray: A float64 index image. Masked/nonfinite operands, zero
            denominators, and nonfinite results become NaN. Input is unchanged.
    """
    data = np.ma.asarray(data, dtype=float).filled(np.nan)
    if data.ndim != 3:
        raise ValueError("data must have shape (height, width, bands).")
    bands = validate_index(method, bands, data.shape[2])
    selected = data[:, :, bands]
    valid = np.isfinite(selected).all(axis=2)
    a = selected[:, :, 0]
    with np.errstate(divide="ignore", invalid="ignore", over="ignore"):
        if method == "brightness":
            result = selected.mean(axis=2)
        else:
            b = selected[:, :, 1]
            if method in _NORMALIZED_DIFFERENCE:
                result = (a - b) / (a + b)
            elif method == "ratio":
                result = a / b
            elif method == "difference":
                result = a - b
            elif method == "exg":
                c = selected[:, :, 2]
                result = (2 * b - a - c) / (a + b + c)
            else:  # VARI
                c = selected[:, :, 2]
                result = (b - a) / (b + a - c)
    return np.where(valid & np.isfinite(result), result, np.nan)

obia.features.texture.TextureConfig dataclass

Choose texture bands, spatial scales, directions, and quantization.

Attributes:

Name	Type	Description
bands	tuple[int, ...] | None	Zero-based band indices; None uses the extraction's selected bands.
features	tuple[str, ...]	GLCM feature names; an empty tuple disables GLCM calculations.
distances	tuple[int, ...]	Positive integer pixel distances, reported separately.
angles	tuple[float, ...]	Directions in radians, measured from increasing columns toward increasing rows. With angle_aggregation='none', output a0, a1, ... refer to positions in this tuple.
levels	int	Number of uniformly spaced gray levels, between 2 and 256.
value_range	tuple[float, float] | None	Fixed (minimum, maximum), or None for each band's finite range over the entire included image. Never fitted per object.
symmetric	bool	Include both pair orders in each co-occurrence matrix.
angle_aggregation	str	'mean' averages available directions; 'none' retains each direction. A direction with no valid pixel pairs is NaN.
entropy_radii	tuple[int, ...]	Positive integer radii for disk-shaped local entropy windows. Their mean and population standard deviation are reported.

Source code in obia/features/texture.py

@dataclass(frozen=True)
class TextureConfig:
    """Choose texture bands, spatial scales, directions, and quantization.

    Attributes:
        bands: Zero-based band indices; None uses the extraction's selected bands.
        features: GLCM feature names; an empty tuple disables GLCM calculations.
        distances: Positive integer pixel distances, reported separately.
        angles: Directions in radians, measured from increasing columns toward
            increasing rows. With angle_aggregation='none', output a0, a1, ...
            refer to positions in this tuple.
        levels: Number of uniformly spaced gray levels, between 2 and 256.
        value_range: Fixed (minimum, maximum), or None for each band's finite
            range over the entire included image. Never fitted per object.
        symmetric: Include both pair orders in each co-occurrence matrix.
        angle_aggregation: 'mean' averages available directions; 'none' retains
            each direction. A direction with no valid pixel pairs is NaN.
        entropy_radii: Positive integer radii for disk-shaped local entropy
            windows. Their mean and population standard deviation are reported.
    """

    bands: tuple[int, ...] | None = None
    features: tuple[str, ...] = GLCM_FEATURES[:7]
    distances: tuple[int, ...] = (1,)
    angles: tuple[float, ...] = (0, np.pi / 4, np.pi / 2, 3 * np.pi / 4)
    levels: int = 32
    value_range: tuple[float, float] | None = None
    symmetric: bool = True
    angle_aggregation: str = "mean"
    entropy_radii: tuple[int, ...] = ()

    def __post_init__(self):
        if self.bands is not None:
            object.__setattr__(self, "bands", tuple(self.bands))
        object.__setattr__(self, "features", choices(self.features, GLCM_FEATURES, "GLCM features"))
        for name in ("distances", "entropy_radii"):
            values = tuple(getattr(self, name))
            if any(isinstance(v, (bool, np.bool_)) or not isinstance(v, (int, np.integer)) or v <= 0
                   for v in values):
                raise ValueError(f"{name} must contain positive integers in pixels.")
            if len(set(values)) != len(values):
                raise ValueError(f"{name} must not contain duplicates.")
            object.__setattr__(self, name, values)
        angles = tuple(float(a) for a in self.angles)
        if not angles or any(not np.isfinite(a) for a in angles) or len(set(angles)) != len(angles):
            raise ValueError("angles must contain distinct finite radians.")
        object.__setattr__(self, "angles", angles)
        if self.features and not self.distances:
            raise ValueError("GLCM features require at least one distance.")
        if isinstance(self.levels, (bool, np.bool_)) or not isinstance(self.levels, (int, np.integer)) or not 2 <= self.levels <= 256:
            raise ValueError("levels must be an integer between 2 and 256.")
        if self.angle_aggregation not in {"mean", "none"}:
            raise ValueError("angle_aggregation must be 'mean' or 'none'.")
        if not isinstance(self.symmetric, (bool, np.bool_)):
            raise TypeError("symmetric must be a boolean.")
        if self.value_range is not None:
            bounds = tuple(float(v) for v in self.value_range)
            if len(bounds) != 2 or not np.isfinite(bounds).all() or bounds[0] >= bounds[1]:
                raise ValueError("value_range must contain finite minimum < maximum.")
            object.__setattr__(self, "value_range", bounds)

Available feature names can be inspected with SPECTRAL_FEATURES, SPECTRAL_INDICES, SHAPE_FEATURES, and GLCM_FEATURES from obia.features. See the raster feature guide for formulas, units, mask handling, and configuration examples.