In-place Ramer–Douglas–Peucker for Go & TypeScript

rdp.go

package rdp

import "math"

type point [2]float64 // [x, y]

// RDPSimplify is an in-place implementation of Ramer–Douglas–Peucker.
func RDPSimplify(points []point, epsilon float64) []point {
	return points[:rdpCompress(points, epsilon)]
}

func rdpCompress(points []point, epsilon float64) int {
	end := len(points)

	if end < 3 {
		// return points
		return end
	}

	// Find the point with the maximum distance
	var (
		first = points[0]
		last  = points[end-1]

		flDist  = distance(first, last)
		flCross = first[0]*last[1] - first[1]*last[0]

		dmax  float64
		index int
	)

	for i := 2; i < end-1; i++ {
		d := perpendicularDistance(points[i], first, last, flDist, flCross)
		if d > dmax {
			dmax, index = d, i
		}
	}

	// If max distance is lte to epsilon, return segment containing
	// the first and last points.
	if dmax <= epsilon {
		// return []point{first, last}
		points[1] = last
		return 2
	}

	// Recursive call
	r1 := rdpCompress(points[:index+1], epsilon)
	r2 := rdpCompress(points[index:], epsilon)

	// Build the result list
	// return append(r1[:len(r1)-1], r2...)
	x := r1 - 1
	n := copy(points[x:], points[index:index+r2])

	return x + n
}

func distance(a, b point) float64 {
	x := a[0] - b[0]
	y := a[1] - b[1]
	return math.Sqrt(x*x + y*y)
}

func perpendicularDistance(p, fp, lp point, dist, cross float64) float64 {
	return math.Abs(cross+
		lp[0]*p[1]+p[0]*fp[1]-
		p[0]*lp[1]-fp[0]*p[1]) / dist
}

rdp.ts

module RDP {
	export interface Point {
		x: number;
		y: number;
	}

	// Simplify is an in-place* implementation of Ramer–Douglas–Peucker.
	export function Simplify(points: Point[], epsilon: number): Point[] {
		epsilon = Math.abs(epsilon); // errors are absolute

		// *makes copy
		return points.slice(0, compress(points, 0, points.length, epsilon));
	}

	function compress(points: Point[], start: number, stop: number, epsilon: number): number {
		const end: number = stop - start;

		if(end < 3) {
			// return points
			return end;
		}

		// Find the point with the maximum distance
		const first: Point = points[start];
		const last: Point = points[stop - 1];

		const flDist: number = distance(first, last);
		const flCross: number = first.x * last.y - first.y * last.x;

		let dmax: number = 0.0;
		let index: number = 0;

		for(let i = 2; i < end - 1; i++) {
			const d: number = perpendicularDistance(points[start + i], first, last, flDist, flCross);
			if(dmax < d) {
				dmax = d;
				index = i;
			}
		}

		// If max distance is lte to epsilon, return segment containing
		// the first and last points.
		if(dmax <= epsilon) {
			points[start + 1] = last;
			return 2;
		}

		index += start;

		// Recursive call
		const r1: number = compress(points, start, index + 1, epsilon);
		const r2: number = compress(points, index, stop, epsilon);

		// Build the result list
		const x: number = r1 - 1;
		for(let i = 0; i < r2; i++) {
			points[start + x + i] = points[index + i];
		}

		return x + r2;
	}

	function distance(a: Point, b: Point): number {
		const x: number = a.x - b.x;
		const y: number = a.y - b.y;
		return Math.sqrt(x * x + y * y);
	}

	function perpendicularDistance(p: Point, fp: Point, lp: Point, dist: number, cross: number): number {
		return Math.abs(cross +
			lp.x * p.y + p.x * fp.y -
			p.x * lp.y - fp.x * p.y) / dist;
	}
}

Sample run on data generated by a sine function.

func main() {
    ps := make([][2]float64, 1000)
    step := 0.04

    var x float64
    for i := range ps {
        ps[i] = point{x * 10, math.Sin(x) * float64(i)}
        x += step
    }

    xysBefore := clonePoints(ps) // simple copy fn
    xysAfter := RDPSimplify(ps, 0.5)

    plot(xysBefore, xysAfter) // uses github.com/gonum/plot
}