ramer_douglas_peucker_funneling.py

"""
This is based on https://en.wikipedia.org/wiki/Ramer%E2%80%93Douglas%E2%80%93Peucker_algorithm,
with the requirement that the maximum norm is used for distance calculations.

The algorithm is as follows, described for R^2, but it is easily extensible to
arbitrary dimensions >2:
- given points p_1, ..., p_n, we start by constructing our funnel
  by calculating
    dp_1 = (p_2 - p_1)
    x_1_top    = (dp_1_y + epsilon) / dp_1_x
    x_1_bottom = (dp_1_y - epsilon) / dp_1_x
- we then iterate over all points p_i with i > 1 and calculate the funnel
  values x_i_top_current and x_i_bottom_current of p_1 and p_i. then:
  x_i_top = min(x_(i-1)_top, x_i_top_current)
  x_i_bottom = max(x_(i-1)_bottom, x_i_bottom_current)
- as soon as we have x_i_top < x_i_bottom, we know we need to retain
  p_(i-1), as p_i is the first outlier, escaping the funnel.
- repeat the whole stuff from p_(i-1) onwards
- add the last point p_n

tadaa
"""

from typing import Tuple, Sequence, List

Sample = Tuple[float, float]


def iterative_split(
    samples: Sequence[Sample], p0: int, epsilon: float
) -> Optional[int]:
    if p0 >= len(samples) - 1:
        return None

    x0, y0 = samples[p0]

    def calc_xi(sample: Sample) -> Tuple[float, float]:
        x, y = sample
        dx = x - x0
        dy = y - y0
        return (
            (dy + epsilon) / dx,
            (dy - epsilon) / dx,
        )

    xi_top, xi_bottom = calc_xi(samples[p0 + 1])

    for i, sample in enumerate(samples[p0 + 2 :], p0 + 2):
        current_xi_top, current_xi_bottom = calc_xi(sample)
        xi_top = min(xi_top, current_xi_top)
        xi_bottom = max(xi_bottom, current_xi_bottom)
        if xi_top < xi_bottom:
            return i - 1
    return len(samples) - 1


def simplify_iterative(all_samples: Sequence[Sample], epsilon: float) -> List[int]:
    keep = [0]
    while True:
        split = iterative_split(all_samples, keep[-1], error_threshold)
        if split is None:
            return keep
        keep.append(split)