Python 3D 多项式曲面拟合，顺序相关

Question

我目前正在处理天文数据，其中有 cometd 图像。由于拍摄时间(黄昏)，我想删除这些图像中的背景天空渐变。我为此开发的第一个程序从 Matplotlib 的“ginput”(x,y) 中获取用户选择的点，提取每个坐标 (z) 的数据，然后使用 SciPy 的“griddata”将数据网格化到一个新数组中。

由于假定背景仅略有变化，因此我想将 3d 低阶多项式拟合到这组 (x,y,z) 点。但是，“griddata”不允许输入顺序:

griddata(points,values, (dimension_x,dimension_y), method='nearest/linear/cubic')

关于可以使用的其他函数或开发允许我控制顺序的最小二乘法拟合的方法有什么想法吗？

Answer 1

Griddata 使用样条拟合。三阶样条与三阶多项式不同(相反，它在每一点都是不同的三阶多项式)。

如果您只想将 2D、3 阶多项式拟合到您的数据，则执行类似以下操作以使用所有您的数据点来估计 16 个系数。

import itertools
import numpy as np
import matplotlib.pyplot as plt

def main():
    # Generate Data...
    numdata = 100
    x = np.random.random(numdata)
    y = np.random.random(numdata)
    z = x**2 + y**2 + 3*x**3 + y + np.random.random(numdata)

    # Fit a 3rd order, 2d polynomial
    m = polyfit2d(x,y,z)

    # Evaluate it on a grid...
    nx, ny = 20, 20
    xx, yy = np.meshgrid(np.linspace(x.min(), x.max(), nx), 
                         np.linspace(y.min(), y.max(), ny))
    zz = polyval2d(xx, yy, m)

    # Plot
    plt.imshow(zz, extent=(x.min(), y.max(), x.max(), y.min()))
    plt.scatter(x, y, c=z)
    plt.show()

def polyfit2d(x, y, z, order=3):
    ncols = (order + 1)**2
    G = np.zeros((x.size, ncols))
    ij = itertools.product(range(order+1), range(order+1))
    for k, (i,j) in enumerate(ij):
        G[:,k] = x**i * y**j
    m, _, _, _ = np.linalg.lstsq(G, z)
    return m

def polyval2d(x, y, m):
    order = int(np.sqrt(len(m))) - 1
    ij = itertools.product(range(order+1), range(order+1))
    z = np.zeros_like(x)
    for a, (i,j) in zip(m, ij):
        z += a * x**i * y**j
    return z

main()