【问题标题】:spline surface interpolation样条曲面插值
【发布时间】:2012-03-28 05:23:21
【问题描述】:

假设我有 n 个点在 z 轴上定义一个曲面

f(x1,y1) = 10
f(x2,y2) = 12
f(x3,y3) = 5
f(x4,y4) = 2
...
f(xn,yn) = 21

现在我希望能够近似 f(x,y)。我正在寻找一种线性算法,尤其是样条近似算法。一个示例算法或至少一些指针会很棒。

【问题讨论】:

  • [wikipedia][1] 文章有点令人生畏,但至少可以尝试查看示例部分。 [1]:en.wikipedia.org/wiki/Spline_interpolation
  • 你是规则网格上的 x,y 控制点吗?
  • 对于 f(x,y) 形式的函数,更常见的是先假设基础数据的形式(K 次多项式、N 高斯之和等),然后通过最小二乘法确定系数。这在这里行得通吗?你知道数据代表什么吗?如果您真的想要样条曲线,那么 NURBS en.wikipedia.org/wiki/NURBS 值得一看。它们具有很好的渲染属性。构造 (x,y) 点的 Delaunay 三角剖分以获得基础,除非它们位于规则网格上。对于平面拟合,您需要标准最小二乘拟合。
  • 我不一定需要成为样条线。线性就足够了。关于平面拟合我不能说太多,但数据点在规则网格上。只是缺少一些数据点。

标签: algorithm math 3d spline


【解决方案1】:

您可以将您的点用作 Bezier(或 Bspline)曲面的控制点,尤其是当 (xi, yi)XY 平面中采样一个矩形时。在这方面,不涉及任何拟合。

您将获得的表面将位于您的点的凸包中,并将与{xi, yi} 边界处的点相交(插值)。

如果你想尝试一下,This forums posting 似乎在Matlab 中包含简单的代码,如果你没有Matlab,你可以使用GuIRIT 来做同样的事情(尽管它需要弄清楚程序的文件格式)。

【讨论】:

  • 最终的实现必须在 ruby​​ 中——所以 Matlab 并不是一个真正的选择。但问题确实是关于 XY 平面上的一个矩形。
  • 我从来没有用过 Ruby,但我确信有一个 Bezier/Bspline 包。
【解决方案2】:

这是对进行线性近似的方法的模糊描述。

  1. 确定您的点的Voronoi diagram(对于平面中的每个点,找到最近的(x_i,y_i)
  2. 取这个的对偶得到Delaunay triangulation:如果有一条线段点连接(x_i,y_i)(x_j,y_j),那么(x_i,y_i)(x_j,y_j)是等距的(并且比任何其他对更近) .
  3. 在每个三角形上,找到通过三个顶点的平面。这是您需要的线性插值。

以下实现了 Python 中的前两个步骤。你的规律性 网格可以让你加快速度(它也可能会弄乱三角测量)。

import itertools

""" Based on http://stackoverflow.com/a/1165943/2336725 """
def is_ccw(tri):
    return ( ( tri[1][0]-tri[0][0] )*( tri[1][1]+tri[0][1] )
            + ( tri[2][0]-tri[1][0] )*( tri[2][1]+tri[1][1] )
            + ( tri[0][0]-tri[2][0] )*( tri[0][1]+tri[2][1] ) ) < 0

def circumcircle_contains_point(triangle,point):
    import numpy as np
    matrix = np.array( [ [p[0],p[1],p[0]**2+p[1]**2,1] for p in triangle+point ] )
    return is_ccw(triangle) == ( np.linalg.det(matrix) > 0 )

triangulation = set()
"""
A triangle is in the Delaunay triangulation if and only if its circumscribing
circle contains no other point.  This implementation is O(n^4).  Faster methods
are at http://en.wikipedia.org/wiki/Delaunay_triangulation
"""
for triangle in itertools.combinations(points,3):
    triangle_empty = True
    for point in points:
        if point in triangle:
            next
        if circumcircle_contains_point(triangle,point):
            triangle_empty = False
            break
    if triangle_empty:
        triangulation.add(triangle)

【讨论】:

    【解决方案3】:

    对不规则的 2D 数据进行插值并不容易。我知道没有真正的样条泛化到不规则 2D。

    除了基于三角测量的方法之外,您还可以查看 Barnes (http://en.wikipedia.org/wiki/Barnes_interpolation) 和反向距离加权 (http://en.wikipedia.org/wiki/Inverse_distance_weighting),或者更普遍的 RBF (http://en.wikipedia.org/wiki/Radial_basis_functions)。

    如果您的点非常不均匀地分布(密集簇),则可能需要使函数的大小具有自适应性,或者采用近似而不是插值。

    【讨论】:

    • 实际上它们的传播范围并不大。在最坏的情况下,我什至可以使用线性插值,但更喜欢更平滑的表面。
    • 近似听起来也是一个有趣的角度。
    • +1 表示径向基函数。几年前,我写了一篇论文,使用这些对象泛化为流形上的函数。只要您没有密集的簇并且点数 n 不会变得太大,它们就可以很好地工作。 (如果 n 确实变大了,你希望你的 RBF 有紧凑的支持,这样所涉及的矩阵是稀疏的。但是你需要使用稀疏线性代数来保持它可以接受的缩放。)很好,平滑的插值。
    猜你喜欢
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2014-08-30
    • 1970-01-01
    • 1970-01-01
    相关资源
    最近更新 更多