2011-09-22 71 views
1

我一直在研究這個問題一段時間,並沒有發現ESRI論壇頁面或我寫的一些FORTRAN三角測試腳本沒有喜悅。Python,GIS和Fortran:試圖從xy點數據創建多個多邊形

我有兩個.csv文件,其中包含數百個xy點數據。這些點表示潮間帶的高端和低端。高點和低點相互平行,我想創建多邊形條,將這些點中的四個連接到每個單獨的多邊形中。多邊形的高度取決於高點和低點之間的距離。下面的鏈接顯示,說明我的意思兩個圖像:

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

的主要問題已腳本的多邊形正確地在角落形成。我知道在彎道上移動時,不能在底部直徑爲1英尺,頂部直徑爲1英尺的多邊形。但是,這僅僅是我在試圖解決此遇到了很多的問題之一......

任何幫助將不勝感激,謝謝

回答

0

這應該插值工作(說潮線爲x ,y距岸上某個控制點的距離):

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 

有很多不同的方法來做多邊形。我可能會嘗試的是確定一個參考海岸線,然後以固定間隔或在每個線段中間取垂直線,然後在相鄰海岸線上的相交處製作一個多邊形或邊界框。順便說一句,在真實應用程序中,您需要確保cmp未在(0,0)上運行。