tin_reducer_gml 0.2

import time
from osgeo import ogr
from pywps.Process.Process import WPSProcess


def iter_points(layer):
    for feature in layer:
        exterior_ring = feature.GetGeometryRef().GetGeometryRef(0)
        for i in xrange(exterior_ring.GetPointCount()):
            yield exterior_ring.GetPoint(i)


def clean_points(points):
    #
    # Delete the last point of each triangle (we only need three point for our
    # calculations)
    #
    # TODO: The code below has the effect of the original code but it doesn't
    #       match the description above and it repeats points which will get
    #       filtered out in a later processing step anyway.  I guess the
    #       `drop_every_nth()` function below is actually wanted here instead.
    #
    points = iter(points)
    try:
        while True:
            yield next(points)
            next(points)
            point = next(points)
            yield point
            yield point
            next(points)
    except StopIteration:
        pass  # Ignored intentionally.


def drop_every_nth(nth, iterable):
    return (item for i, item in enumerate(iterable, 1) if i % nth != 0)


def distance(point_a, point_b):
    #
    # TODO: To be the usual euklidean distance formula it lacks a `math.sqrt()`
    #       call!
    #
    return sum((a - b)**2 for a, b in zip(point_a, point_b))


def filter_close_points(minimal_distance, points):
    #
    # TODO: The reversing is only neccessary if the points have to be looked
    #       at in the order they are given to this function.  If going through
    #       them in reverse would be okay we could spare this step.
    #
    points = list(reversed(points))
    while points:
        point = points.pop()
        point_count = len(points)
        points = [p for p in points if distance(point, p) >= minimal_distance]
        if len(points) == point_count:
            yield point


def create_geometry(points):
    geometry = ogr.Geometry(ogr.wkbMultiPoint)
    wkb_point = ogr.Geometry(ogr.wkbPoint)
    for point in points:
        wkb_point.AddPoint(*point)
        geometry.AddGeometry(wkb_point)
    return geometry


class Process(WPSProcess):
    def __init__(self):
        WPSProcess.__init__(
            self,
            identifier=__name__,
            title='Pollmueller Tin-Reducer',
            version = '0.1',
            storeSupported='true',
            statusSupported='true',
            abstract=(
                'reduces the points of a given triangulated irregular network'
                ' by a variable factor'
            ),
            grassLocation=False
        )
        self.data_output = self.addComplexOutput(
            identifier='output',
            title='Output vector data',
            formats=[{'mimeType': 'text/xml'}]
        )
        self.time_output = self.addLiteralOutput(
            identifier='output2',
            title='computation time (to compare different approaches)'
        )
        self.data_filename_input = self.addComplexInput(
            identifier='data',
            title='Input vector data',
            formats=[{'mimeType': 'text/xml'}]
        )
        self.minimal_distance_input = self.addLiteralInput(
            identifier='input2', title='minimal distance between points'
        )

    def execute(self):
        #
        # TODO: Different operating systems have different semantics and
        #       resolutions for `time()` and `clock()`.  Make sure `clock()` is
        #       the correct one here.
        #
        start_timestamp = time.clock()

        shape_file = ogr.Open(self.data_filename_input.getValue())
        if shape_file is None:
            raise Exception('failed to open shape file.')
        #
        # Get the first layer.
        #
        layer = shape_file.GetLayer(0)
        layer.ResetReading()
        points = filter_close_points(
            self.minimal_distance_input.getValue(),
            clean_points(iter_points(layer))
        )
        gml = create_geometry(points).ExportToGML()

        elapsed_time = time.clock() - start_timestamp
        self.data_output.setValue(gml)
        self.time_output.setValue(elapsed_time)