сплайн-интерполяция поверхности - PullRequest
Новая
       11

сплайн-интерполяция поверхности

8 голосов
/ 11 марта 2012

Допустим, у меня есть 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). Я ищу алгоритм для линейного и особенно сплайн-приближения. Пример алгоритмов или хотя бы несколько указателей было бы здорово.

Ответы [ 3 ]

4 голосов
/ 06 мая 2014

Это расплывчатое описание подхода к линейному приближению.

  1. Определите диаграмму Вороного ваших точек (для каждой точки на плоскости найдите ближайшую(x_i,y_i))
  2. Возьмите двойное из этого, чтобы получить триангуляцию Делоне : соедините (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)
2 голосов
/ 06 мая 2014

Интерполяция нерегулярных 2D-данных не так проста.Я не знаю ни одного истинного обобщения сплайна для нерегулярного 2D.

Помимо подходов, основанных на триангуляции, вы можете взглянуть на Барнса (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 голос
/ 11 марта 2012

Вы можете использовать свои точки в качестве контрольных точек поверхности Безье (или Bspline), особенно если (xi, yi) выбирает прямоугольник в плоскости XY.В этом отношении не требуется подгонка.

Поверхность, которую вы получите, будет находиться в выпуклой оболочке ваших точек и будет пересекать (интерполировать) точки на границе {xi, yi}.

Если вы хотите поэкспериментировать, На этом форуме , похоже, содержится простой код в Matlab, и вы можете использовать GUIRIT , чтобы сделать то же самое, если у вас нетMatlab (хотя это требует определения формата файла программы).

...