Python, GIS and Fortran: Trying to create multiple polygons from xy point data

I've been working on this problem for awhile and have found no joy on the ESRI forum page or with some FORTRAN triangulation script I wrote.

I have two .csv files with h开发者_StackOverflow中文版undreds of xy point data in. These points represent the high and low end of an intertidal range. The high and low points run parallel to each other and I want to create polygon slivers that connect four of those points every foot into separate polygons. The height of the polygon would be x depending upon the distance between the high and low points. The below link shows two images that illustrate what I mean:

http://forums.arcgis.com/threads/39757-Feature-to-Line..?p=135880&posted=1#post135880

The main problem has been scripting the polygons to form correctly in the corners. I understand that you cannot have a polygon that is 1ft in diameter at the bottom and 1ft in diameter at the top when moving round a bend. But this is just one of many problems I've encountered in trying to solve this...

Any help would be greatly appreciated, Thanks

This should work for the interpolation (say the tide line is x, y distance from some control point on shore):

import math

example_triangle = [[0,0], [5,5], [0,5], [0,0]]
high_tide_line = [[0, 0], [5., 1.5], [10., 3.2], [20., 1.], [30., 4.5], [80.,2.], [80,0]]
low_tide_line = [[0, 10.], [5., 11.5], [10., 13.2], [20., 11.], [30., 14.5], [80., 12.], [80, 10]]

def points_from_geom(geom):
    idx = 0
    line_lengths = []
    unit_vectors = []
    interpolated_points = []
    while idx < (len(geom) - 1):
        dy, dx = ((geom[idx+1][1] - geom[idx][1]), (geom[idx+1][0] - geom[idx][0]))
        line_lengths.append(math.sqrt(dy**2 + dx**2))
        try:
            angle = math.atan(dy/dx)
            unit_vectors.append([math.cos(angle)*cmp(dx, 0),
                math.sin(angle)*cmp(dy, 0)])
        except ZeroDivisionError:
            if geom[idx+1][1] < geom[idx][1]:
                direction = [0, -1]
            else:
                direction = [0, 1]
            unit_vectors.append(direction)
        idx += 1

    for i, length in enumerate(line_lengths):
        inter = 0
        while inter <= length:
            interpolated_points.append([geom[i][0] + unit_vectors[i][0]*inter,\
                geom[i][1] + unit_vectors[i][1]*inter])
            inter += .3048 # a ft in proper units ;)

    return interpolated_points

ln1 = points_from_geom(example_triangle)
ln2 = points_from_geom(high_tide_line)
ln3 = points_from_geom(low_tide_line)

print ln1, ln2, ln3

There are a lot of different ways to do the polygons. What I might try is determine a reference shoreline, then take perpendicular lines to it at fixed intervals, or in the middle of each line segment, then make a polygon or bounding box with the intersection on the adjacent shoreline. BTW in a real application you would need to make sure cmp is not operating on (0,0).