Metadata-Version: 2.5
Name: gri-convolve
Version: 0.6.0
Summary: Ellipsoid convolution for combining geolocation estimates with outlier detection
Project-URL: Homepage, https://geosolresearch.com
Project-URL: Repository, https://gitlab.com/geosol-foss/python/gri-convolve
Project-URL: Issues, https://gitlab.com/geosol-foss/python/gri-convolve/-/issues
Project-URL: Changelog, https://gitlab.com/geosol-foss/python/gri-convolve/-/releases
Author-email: GeoSol Research Inc <contact@geosolresearch.com>
License-Expression: MIT
License-File: LICENSE
Keywords: clustering,convolution,ellipsoid,fusion,geolocation
Classifier: Development Status :: 3 - Alpha
Classifier: Intended Audience :: Science/Research
Classifier: License :: OSI Approved :: MIT License
Classifier: Programming Language :: Python :: 3
Classifier: Topic :: Scientific/Engineering
Classifier: Topic :: Scientific/Engineering :: GIS
Requires-Python: >=3.12
Requires-Dist: gri-ell>=0.2.0
Requires-Dist: gri-pos>=0.2.0
Requires-Dist: gri-utils>=0.5.1
Requires-Dist: numpy>=2.3.3
Requires-Dist: scipy>=1.16.3
Description-Content-Type: text/markdown

[![GeoSol Research Logo](https://geosolresearch.com/logos/foss_logo.png "GeoSol Research")](https://geosolresearch.com)

# Convolve (Ellipsoid Fusion)

Ellipsoid convolution functions for combining geolocation estimates with outlier detection and multi-cluster support.

## Overview

gri-convolve provides three functions of increasing sophistication for fusing collections of `Ell` (ellipsoid) objects into combined position estimates:

- **`convolve`** -- combine all input ellipsoids into a single fused result with no outlier rejection
- **`smart_convolve`** -- iteratively remove outliers by Mahalanobis distance before fusing
- **`cluster_convolve`** -- find multiple clusters within a dataset and fuse each independently

Each function operates on `Ell` objects from gri-ell, which pair a 3D position with a statistical covariance (or information matrix). The output is one or more fused `Ell` objects representing the combined position estimate and its uncertainty.

Recursive single-target tracking -- the `IMM` / `SmartSegmentedIMM` filters, the
motion-model bank, EKF/UKF observable updates, and RTS smoothing -- is provided
by the companion **`gri-kalman`** package. Convolution is the batch face and the
filter is the recursive face of the same estimation problem: the convolver seeds
and refines tracks, the filter does the online sequential update. See the
`gri-kalman` README for the tracker API and EKF-vs-UKF guidance.

Requires Python 3.12+.

## Mathematical Background

**Information matrix fusion.** Given N ellipsoids, each with position `x_k` and information matrix `I_k` (the inverse of the covariance matrix, in XYZ coordinates, 1/m^2, 1-sigma), the fused position and information matrix are:

    S = sum(I_k)              (combined information matrix)
    x = S^{-1} sum(I_k x_k)  (fused position)

This is the maximum-likelihood estimator under Gaussian assumptions.

**Inflation methods.** The raw fusion above underestimates uncertainty when inputs are inconsistent. Three modes control how the output covariance is inflated:

- `"none"` -- strict information matrix combination (no inflation)
- `"std"` -- inflate by the sample scatter of input positions in XYZ
- `"bart"` -- inflate along the semi-major axis direction in ENU (recommended default)

**Outlier detection.** `smart_convolve` computes the Mahalanobis distance from each input to the fused point:

    d_M = sqrt((x - mu)^T I (x - mu))

where `mu` is the fused position and `I` is its information matrix. Distances are normalized to 95% confidence scale. Inputs exceeding `max_norm` are iteratively removed, worst first.

Reference: Mahalanobis, P.C. (1936). "On the generalized distance in statistics."

## Installation

```bash
pip install gri-convolve
```

For development:

```bash
git clone https://gitlab.com/geosol-foss/python/gri-convolve.git
cd gri-convolve
uv sync
```

## Quick Start

```python
from gri_convolve import convolve, smart_convolve, cluster_convolve
from gri_ell import Ell
from gri_pos import Pos
import numpy as np

# Create some ellipsoids at nearby positions
e1 = Ell.from_2d(Pos.LLA(40.0, -105.0, 1600), 100, 50, 45)
e2 = Ell.from_2d(Pos.LLA(40.001, -105.001, 1610), 120, 60, 30)
e3 = Ell.from_2d(Pos.LLA(40.0005, -104.999, 1605), 90, 45, 50)

# Simple fusion
fused = convolve([e1, e2, e3])
print(fused.lla)            # Fused position
print(fused.ellipse.sma_95) # Fused semi-major axis (95%, meters)
```

## `convolve()`

Fuses all input ellipsoids into a single result. No outlier detection.

```python
fused = convolve(ells, inflation="bart")
```

**Parameters:**

- `ells` -- sequence or generator of `Ell` objects
- `inflation` -- `"none"`, `"std"`, or `"bart"` (default: `"bart"`)

**Returns:** A single fused `Ell`.

## `smart_convolve()`

Fuses with iterative outlier rejection. Computes the fused point, finds the input with the largest normalized Mahalanobis distance, and removes it if it exceeds `max_norm`. Repeats until all remaining inputs are within tolerance or fewer than `min_pts` remain.

```python
result = smart_convolve(ells, max_norm=2.0, min_pts=3)
if result is not None:
    fused_ell, used_indices, discarded_indices = result
```

**Parameters:**

- `ells` -- sequence or generator of `Ell` objects
- `max_norm` -- maximum allowed normalized distance (default: 2.0)
- `min_pts` -- minimum inputs required for a valid result (default: 3)
- `check_scale` -- warn on a likely unit mismatch or near-total rejection (default: True; see [Unit sanity checks](#unit-sanity-checks))

**Returns:** `(Ell, list[int], list[int])` or `None` if no valid cluster is found.

Pre-cluster your data before calling `smart_convolve`. Without pre-clustering, a large group of scattered noise points can cause valid clusters to be discarded first.

## `cluster_convolve()`

Finds multiple clusters within a dataset by iteratively applying `smart_convolve`. After finding the largest valid cluster, the discarded points are passed back in to find additional clusters.

```python
locations, used_per_location, discarded = cluster_convolve(
    ells,
    max_norm=2.0,
    min_pts=3,
    max_pts=10,
    min_sma_m=50.0,
)

for loc, indices in zip(locations, used_per_location):
    print(f"Cluster at {loc.lla} using {len(indices)} inputs")
```

**Parameters:**

- `ells` -- sequence or generator of `Ell` objects
- `max_norm` -- maximum normalized Mahalanobis distance (default: 2.0)
- `min_pts` -- minimum inputs per cluster (default: 3)
- `max_pts` -- maximum inputs per cluster; splits larger groups (default: None)
- `min_sma_m` -- minimum semi-major axis for output ellipsoids in meters (default: 0)
- `max_ori_spread` -- sort by orientation before splitting for diversity (default: True)
- `check_scale` -- warn on a likely unit mismatch (default: True; see [Unit sanity checks](#unit-sanity-checks))
- `alt_post_process` -- callback for altitude correction (e.g., snap to terrain)

**Returns:** `(list[Ell], list[list[int]], list[int])`

## Unit sanity checks

The convolution math is scale-invariant -- a Mahalanobis distance is dimensionless, so entering every quantity in kilometers instead of meters changes nothing. Trouble appears only when the *position* scale and the *covariance* scale disagree. The common mistake is positions in meters (positions always resolve to ECEF meters) paired with covariance sigmas typed in kilometers: the information matrices come out about a million times too large, every residual distance is about a thousand times too large, so `smart_convolve` and `cluster_convolve` reject nearly every point one at a time -- slow, and with an almost-empty result.

To catch this, `smart_convolve` and `cluster_convolve` run two cheap checks on the full input (controlled by `check_scale`, default `True`) and emit a `ConvolveScaleWarning` when either trips:

- **Scale pre-flight** (before the rejection loop): if the tightest-packed inputs still sit many 1-sigma widths from their nearest neighbor -- meaning no cluster can form -- it warns that the sigmas look too small for the point spacing, and flags a ratio near 1000 as a likely kilometer/meter mix-up. Because it measures nearest-neighbor spacing (a local quantity), legitimate multi-site data with many separate clusters does not trip it.
- **Utilization backstop** (after): if nearly all inputs land in the discard list, it warns that most of the data was rejected. (`cluster_convolve` suppresses this when a cluster size floor is set, since trimming small clusters discards points by design.)

The checks never change the result -- they only warn. Filter or catch them with the exported `ConvolveScaleWarning` category, or pass `check_scale=False` to silence them (e.g. for data that is intentionally spread far in sigma units):

```python
import warnings
from gri_convolve import ConvolveScaleWarning, smart_convolve

with warnings.catch_warnings():
    warnings.simplefilter("error", ConvolveScaleWarning)  # promote to an exception
    result = smart_convolve(ells)
```

## Degenerate inputs are rejected

Separate from the warnings above, and not optional. Every semi-axis of an input ellipsoid must be finite and greater than zero. A zero-width axis -- most often `alt_95_m=0`, written when only a 2D ellipse is available -- makes the covariance rank deficient, and there is no fused answer to give.

This used to fail in whichever of four ways the arithmetic happened to land in, decided by float noise in the ENU-to-ECEF rotation at that latitude and longitude: a bare `LinAlgError: Singular matrix` at build time, the same error much later from the inflation step, a silent `None` with every point discarded, or a silently wrong answer built from a subset. It now raises `ConvolveInputError` up front, naming the offending rows:

```python
from gri_convolve import ConvolveInputError, smart_convolve

try:
    result = smart_convolve(ells)
except ConvolveInputError as exc:
    print(exc)
    # 1 of 5000 input covariances have a variance of zero or less on the
    # diagonal: rows [4999]. A zero-width axis (commonly alt_95_m=0) has no
    # invertible covariance; the convolve solves in 3D, so even a pinned
    # altitude needs a real uncertainty.
```

`ConvolveInputError` subclasses `numpy.linalg.LinAlgError` (itself a `ValueError`), so code already guarding a convolve call against a singular matrix keeps working and simply gets a message that says which input is at fault.

No default is substituted for a zero altitude uncertainty. The solver works in 3D, and whether the right value is a meter, ten meters, or a fraction of the semi-major axis is a modeling decision that belongs to the caller. A 2D-only `Ell` built without `alt_95_m` is rejected with its own message for the same reason.

Eccentricity is not degeneracy: a semi-major/semi-minor ratio of 1e7 is accepted. `params_to_arrays` builds the information matrix directly from the ellipse parameters instead of inverting a covariance, so it keeps full float64 precision at any ratio, where inverting loses accuracy as the square of the ratio.

## Units and Conventions

- Positions are in ECEF XYZ (meters) internally
- Information matrices are in XYZ, 1/m^2, 1-sigma
- Covariance matrices are in ENU, m^2, 1-sigma
- Output ellipse parameters (SMA, SMI, orientation) are at 95% confidence
- Mahalanobis distances are normalized to 95% scale for `max_norm` comparisons

## Dependencies

- **gri-ell**: Ellipsoid objects with position and covariance
- **gri-pos**: Position objects (XYZ, LLA coordinates)
- **gri-utils**: Coordinate conversions and constants
- **numpy**: Array operations
- **scipy**: Nearest-neighbor search (`KDTree`) for the unit sanity check


## Other Projects

Current list of other [GRI FOSS Projects](https://gitlab.com/geosol-foss/python/gri-convolve/-/blob/main/.docs_other_projects.md) we are building and maintaining.

## License

MIT License. See [LICENSE](https://gitlab.com/geosol-foss/python/gri-convolve/-/blob/main/LICENSE) for details.
