使用 numpy 获得二维多项式变换的逆（用于图像或光栅图像变形/采样）

Question

如果我有一个一阶/仿射、二阶或三阶的二维（x 和 y 坐标）多项式变换函数（即我有系数/变换矩阵 A），那么数学或编程方法是什么得到这个函数的精确倒数？理想情况下，我将如何在 Numpy 中实现这一点？这是在图像变形或地图地理配准的背景下，即在新的变形坐标系中将坐标从输入图像转换或变形到输出图像。

尝试的解决方案

为了解决这个问题，我尝试了一种矩阵代数方法来求解方程组。在数学上，转换过程表示为Au = v。前向变换很容易，您可以根据输入坐标将u 计算为包含多项式方程项的列矩阵，然后将u 与变换矩阵A 相乘，以获得转换后的输出列矩阵v 包含输出坐标。另一方面，反向变换意味着我们知道输出坐标v，并希望找到输入坐标u，因此我们需要将方程重新洗牌为u = Av。根据矩阵代数的规则，A 矩阵在移动时必须反转。在 Numpy 中为二阶多项式变换实现这一点，它似乎确实有效：

import numpy as np

# input coords
x = np.array([13])
y = np.array([13])
    
# terms of the 2nd order polynomial equation
x = x
y = y
xx = x*x
xy = x*y
yy = y*y
ones = np.ones(x.shape)

# u consists of each term in 2nd order polynomial equation
# with each term being array if want to transform multiple
u = np.array([xx,xy,yy,x,y,ones])
print('original input u', u)

## output:
## ('original input u', array([[169.],
##       [169.],
##       [169.],
##       [ 13.],
##       [ 13.],
##       [  1.]]))

# forward transform matrix
A = np.array([[1,2,3,1,6,8],
              [5,2,9,2,0,1],
              [8,1,5,8,4,3],
              [1,4,8,2,3,9],
              [9,3,2,1,9,5],
              [4,2,5,6,2,1]])

# get forward coords
v = A.dot(u)
print('output v', v)

## output:
## ('output v', array([[1113.],
##       [2731.],
##       [2525.],
##       [2271.],
##       [2501.],
##       [1964.]]))

# get backward coords (should exactly reproduce the input coords)
Ainv = np.linalg.inv(A)
u_pred = Ainv.dot(v)
print('backwards predicted input u', u_pred)

## output:
## ('backwards predicted input u', array([[169.],
##       [169.],
##       [169.],
##       [ 13.],
##       [ 13.],
##       [  1.]]))

在上面的示例中，输出 v 实际上是一个 1x6 矩阵，其中只有前两行/值表示转换后的 x 和 y 坐标。问题变成了我们需要所有 v 中的附加值，以便精确地反转坐标。但在实际场景中，我们只知道转换后的 x 和 y 值（即v 的前两行/值），我们不知道完整的 1x6 v 矩阵。

也许我在想这个错误，或者这个矩阵代数方法不是正确的方法，因为二阶多项式和更高阶的多项式不再是线性的？任何用于反转多项式变换的替代编程/numpy 方法？

一些上下文

我查找了许多类似的问题和网站以及诸如numpy.polynomial.Polynomial.fit 之类的 numpy 函数，但它们中的大多数仅与反一维多项式变换有关。我发现的几个关于二维变换的链接说没有确切的方法来反转它，这是没有意义的，因为这是图像变形/重采样和地图地理配准中非常常见的操作。例如，扭曲图像的步骤通常分解为：

1. 使用变换函数/矩阵A 将所有原始像素（列-行）坐标u 前向投影，以便找到变换后的坐标空间v 的边界。
2. 然后，对于在变换后的坐标空间边界（在步骤 1 中找到）中以规则间隔采样的每个坐标，在变换后的坐标系中反向采样这些 v 坐标以找到它们的原始坐标 u。这决定了为变换图像中的每个位置采样的原始像素。

然后我的问题是我有第 1 步所需的前向变换，但我需要找到第 2 步中反向采样所需的变换的精确逆。数学答案或 numpy 解决方案都可以。

Answer 1

2D 仿射函数的反演非常简单。它需要 2x2 线性方程组的分辨率。

二次多项式和三次多项式的情况更成问题。如果我是对的，两个未知数的系统等效于单个四次或非二次（9 次）多项式方程。四次情况下存在显式（虽然很复杂）公式，但非二次情况下没有公式，您将不得不求助于数值方法（牛顿迭代法）。

此外，这些非线性方程的解不是唯一的（您可以有 4 或 9 个解），您需要保留正确的解。

---

如果您的变换仍然接近仿射（例如在校正图像失真时），我建议选择一个近似完整方程的仿射变换，使用反向变换找到初始近似值，然后使用牛顿进行细化。