2016-09-07 41 views
-1

私は[0、a] x [0、b]平面のあらゆる点(x、y)に(x、y)\に不規則に分布しているのに対し、約80.000点x、y)の場合、物理量zはある値をとる。データをさらに評価するために、グリッド上にそれを補間したい。Pythonの非線形グリッド上で不規則に分布したデータを補間する方法は?

以前は、scipy.interpolate.griddataを使用して、正則な2次の2次元グリッド上に補間しました。しかし、この規則的なグリッドは、zが劇的に変化する領域を適切にモデル化することはできず、zのわずかな変化しかない領域には多くのデータ点が存在するという欠点がある。

私は、zの劇的な変化の領域に多くのグリッド点をもち、zのわずかな変化の領域にはデータ点が少ない、非線形(好ましくはまだ2次的ですが、可変メッシュサイズ)のグリッドを望みます。

答えて

0

あなたはそれを後ろに持っていると思います:あなたのグリッドはできるだけ定期的にすることができますが、各サンプルポイントの同じ数を使って各グリッドポイントを評価する必要があります。データスパース領域での滑らかさ。

私はそのために逆距離加重木を使用します。 私はpythonで漂っている実装:

class invdisttree(object): 
    """ 
    Compute the score of query points based on the scores of their k-nearest neighbours, 
    weighted by the inverse of their distances. 

    @reference: 
    https://en.wikipedia.org/wiki/Inverse_distance_weighting 

    Example: 
    -------- 

    import numpy as np 
    import matplotlib.pyplot as plt 
    from invdisttree import invdisttree 

    import matplotlib.pyplot as plt 

    # create sample points with structured scores 
    X1 = 10 * np.random.rand(1000, 2) -5 

    def func(x, y): 
     return np.sin(x**2 + y**2)/(x**2 + y**2) 

    z1 = func(X1[:,0], X1[:,1]) 

    # 'train' 
    tree = invdisttree(X1, z1) 

    # 'test' 
    spacing = np.linspace(-5., 5., 100) 
    X2 = np.meshgrid(spacing, spacing) 
    grid_shape = X2[0].shape 
    X2 = np.reshape(X2, (2, -1)).T 
    z2 = tree(X2) 

    fig, (ax1, ax2, ax3) = plt.subplots(1,3, sharex=True, sharey=True, figsize=(10,3)) 
    ax1.contourf(spacing, spacing, func(*np.meshgrid(spacing, spacing))) 
    ax1.set_title('Ground truth') 
    ax2.scatter(X1[:,0], X1[:,1], c=z1, linewidths=0) 
    ax2.set_title('Samples') 
    ax3.contourf(spacing, spacing, z2.reshape(grid_shape)) 
    ax3.set_title('Reconstruction') 
    plt.show() 

    """ 
    def __init__(self, X=None, z=None, leafsize=10): 
     if not X is None: 
      self.tree = cKDTree(X, leafsize=leafsize) 
     if not z is None: 
      self.z = z 

    def fit(self, X=None, z=None, leafsize=10): 
     """ 
     Arguments: 
     ---------- 
      X: (N, d) ndarray 
       Coordinates of N sample points in a d-dimensional space. 
      z: (N,) ndarray 
       Corresponding scores. 
      leafsize: int (default 10) 
       Leafsize of KD-tree data structure; 
       should be less than 20. 

     Returns: 
     -------- 
      invdisttree instance: object 
     """ 
     return self.__init__(X, z, leafsize) 

    def __call__(self, X, k=6, eps=1e-6, p=2, regularize_by=1e-9): 
     self.distances, self.idx = self.tree.query(X, k, eps=eps, p=p) 
     self.distances += regularize_by 
     weights = self.z[self.idx.ravel()].reshape(self.idx.shape) 
     mw = np.sum(weights/self.distances, axis=1)/np.sum(1./self.distances, axis=1) 
     return mw 

    def transform(self, X, k=6, p=2, eps=1e-6, regularize_by=1e-9): 
     """ 
     Arguments: 
     ---------- 
      X: (N, d) ndarray 
       Coordinates of N query points in a d-dimensional space. 

      k: int (default 6) 
       Number of nearest neighbours to use. 

      p: int or inf 
       Which Minkowski p-norm to use. 
       1 is the sum-of-absolute-values "Manhattan" distance 
       2 is the usual Euclidean distance 
       infinity is the maximum-coordinate-difference distance 

      eps: float (default 1e-6) 
       Return approximate nearest neighbors; the k-th returned value 
       is guaranteed to be no further than (1+eps) times the 
       distance to the real k-th nearest neighbor. 

      regularise_by: float (default 1e-9) 
       Regularise distances to prevent division by zero 
       for sample points with the same location as query points. 

     Returns: 
     -------- 
      z: (N,) ndarray 
       Corresponding scores. 
     """ 
     return self.__call__(X, k, eps, p, regularize_by) 
+0

おかげであなたの答えのためのscipy.spatial輸入cKDTreeからNP として 輸入numpyの。私は持っていない、なぜ同じポイントの点ですべてのグリッドポイントを評価すると、この場合私を助けるだろう。最後に、補間されたデータを有限差分方程式のシステムに入れたいと思います。グリッドはまだ規則正しいと言われているので、これがどのようにして有限差分方程式(それ自体近似)の精度を向上させるのか分かりません。つまり、データポイントに関しては、地域が人口密度の高い地域や人口の少ない地域ではないということです。一部の地域で物理量zが急激に変化したことについて、補間したいことです。 – Ben

+0

あなたは有限差分方程式で結果を使用したいとは言及していませんでした。 scipy.interpolateはスプラインに依存しています。スプラインは一般に「グローバル」な滑らかさのパラメータ(AFAIK;私はscipy.interpolateの実装をチェックしていません)を持っています。 – Paul

+0

また、滑らかに変化するメッシュサイズのグリッドを四角形で作ることはできません。三角形はあなたが見ているものです。 – Paul

関連する問題