Douglas Peucker

gistfile1.py

from math import sqrt
from itertools import islice, chain

def subtract_coords(a, b):
    '''Returns a - b'''

    return (a[0] - b[0], a[1] - b[1])

def dot_product(a, b):
    return a[0] * b[0] + a[1] * b[1]

def magnitude(a):
    return sqrt(a[0] * a[0] + a[1] * a[1])

def douglas_peucker(coords, epsilon):
    result = _douglas_peucker(coords, epsilon, 0, len(coords))
    return list(result)

def _douglas_peucker(coords, epsilon, start_at, end_at):
    head_coord = coords[start_at]
    tail_coord = coords[end_at - 1]

    dmax = 0
    index = 0

    u = subtract_coords(tail_coord, head_coord)
    # the distance from the first to the last coord
    total_distance = magnitude(u)

    if total_distance > 0:
        u = (u[0] / total_distance, u[1] / total_distance)

    for i in xrange(start_at + 1, end_at - 1):
        v = subtract_coords(coords[i],  head_coord)
        dist_from_head = magnitude(v) # the distance from the head_coord to coord[i]

        d = dist_from_head
        if total_distance > 0.0:
            t = dot_product(u, v) /  dist_from_head
            d *= sqrt(1 - t * t)

        if d > dmax:
            index = i
            dmax = d

    # if max distance is greater than epsilon, recursively simplify
    if dmax >= epsilon:
        left_half =  _douglas_peucker(coords, epsilon, start_at, index + 1)
        right_half = _douglas_peucker(coords, epsilon, index, end_at)

        i = iter(right_half)
        i.next()

        for coord in left_half:
            yield coord
        for coord in i:
            yield coord

    else:
        yield head_coord
        yield tail_coord

if __name__ == '__main__':

    less_complicated = douglas_peucker(complicated_linestring, 0.0005)