2012-01-09 41 views
5

我有一個大型的網格數據點,我已經從模擬產生了,並且與xy平面中的每個點相關聯的是z值(模擬結果)。由數據點定義的「平面」下的體積-python

我將x,y,z值轉儲爲純文本文件,我想要做的是測量xy平面(即z = 0)和由數據點。數據點目前並不是均勻間隔的,儘管它們應該在模擬運行完成之後。

我一直在尋找scipy文檔,我不確定scipy.integrate是否提供了我需要的功能 - 似乎只有在2d中才能做到這一點,而不是我需要的功能。首先,除非有必要,我可以在沒有內插的情況下完全基於「梯形法則」或類似的近似來進行積分,這是一個很好的基礎。

任何幫助表示讚賞。

感謝

編輯:雙雙跌破工作描述以及解決方案。就我而言,事實證明使用樣條可能會在飛機的尖銳最大值周圍產生「漣漪」,因此Delaunay方法效果更好,但我建議人們檢查兩者。

回答

3

如果要嚴格堅持梯形規則,你可以做一些與此類似:

import numpy as np 
import scipy.spatial 

def main(): 
    xyz = np.random.random((100, 3)) 
    area_underneath = trapezoidal_area(xyz) 
    print area_underneath 

def trapezoidal_area(xyz): 
    """Calculate volume under a surface defined by irregularly spaced points 
    using delaunay triangulation. "x,y,z" is a <numpoints x 3> shaped ndarray.""" 
    d = scipy.spatial.Delaunay(xyz[:,:2]) 
    tri = xyz[d.vertices] 

    a = tri[:,0,:2] - tri[:,1,:2] 
    b = tri[:,0,:2] - tri[:,2,:2] 
    proj_area = np.cross(a, b).sum(axis=-1) 
    zavg = tri[:,:,2].sum(axis=1) 
    vol = zavg * np.abs(proj_area)/6.0 
    return vol.sum() 

main() 

無論花鍵或線性(trapezodial)插值是一個更適合將在很大程度上取決於你的問題。

+0

謝謝,我也會看看這個選擇。 – samb8s 2012-01-10 17:05:47