算法：如何计算双线性插值的逆？ [英] Algorithm: how calculate INVERSE of bilinear interpolation?

问题描述

双线性插值法是微不足道的计算。 但我需要一个算法，做相反的操作。 （算法将有利于我的假code，或任何广泛使用的计算机语言）

例如，这里是一个Visual Basic实现双线性插值的。

 xyWgt范围（0,1）在x和y。 （0,0）将返回X0Y0，
（0,1）将返回X0Y1等
'例如，如果xyWgt是图像内的相对位置，
'和XnYn值的GPS的coords在4个角的图像，
其结果是对应于xyWgt GPS坐标。
' 例如。定（0.5，0.5），则结果将是在图像的中心的GPS坐标。
公共功能Lerp2D（xyWgt作为的Point2D，X0Y0作为的Point2D，X1Y0作为的Point2D，X0Y1作为的Point2D，X1Y1作为的Point2D）作为的Point2D
    昏暗xY0作为的Point2D =线性插值（X0Y0，X1Y0，xyWgt.X）
    昏暗XY1作为的Point2D =线性插值（X0Y1，X1Y1，xyWgt.X）

    昏暗的XY作为的Point2D =线性插值（xY0，XY1，xyWgt.Y）
    返回XY
端功能
 

其中

 加权两点平均。
公共职能线性插值（BYVAL一个作为的Point2D，BYVAL b以的Point2D，BYVAL wgtB作为双人间）作为的Point2D
    返回新的Point2D（线性插值（AX，BX，wgtB），线性插值（AY，BY，wgtB））
端功能
 

和

 加权两个数字的平均值。
当wgtB == 0，返回，当wgtB == 1，返回湾
含蓄，wgtA = 1  -  wgtB。也就是说，权重进行归一化。
公共职能线性插值（BYVAL一个为DOUBLE，BYVAL b以双人间，BYVAL wgtB作为双人间）作为双
    返回一个+（wgtB *（B  -  A））
端功能
 

---

在1-D，I已经确定线性插值的反函数：

 计算wgtB，将返回的结果，如果没有线性插值（A，B，wgtB）。
'即，在结果是，w.r.t. a和b。
'＆LT; 0是一个，＆GT之前; 1为b之后。
公共职能WgtFromResult（BYVAL一个为DOUBLE，BYVAL b以双人间，BYVAL结果作为双人间）作为双

    昏暗的分母为Double = B  - 一个

    如果Math.Abs​​（分母）＆LT; 0.00000001然后
        '避免除以零（一个和b几乎相等）。
        如果Math.Abs​​（结果 - 一）＆LT; 0.00000001然后
            结果是接近（但也到b）：给最简单的答案是：平均它们。
            返回0.5
        结束如果
        '无法计算。
        返回Double.NaN
    结束如果

    结果= A +（WGT *（B  -  A））=＆GT;
    WGT *（B  -  A）=（结果 - 一）=＆GT;
    昏暗WGT作为双=（结果 - 一）/分母

    昏暗的验证作为双=线性插值（A，B，WGT）
    如果不NearlyEqual（结果验证），然后
    昏暗的测试= 0'测试
    结束如果

    返回WGT
端功能
 

现在我需要在2-D做相同的：

 '返回xyWgt，这要是给Lerp2D，将返回这个XY。
所以，如果XY = X0Y0，将返回（0，0）。如果XY = X1Y0，将返回（1,0）等。
'例如，如果4个角是GPS坐标中的图像的角部，
，并通过在GPS坐标，
'图像内返回相对位置。
公共功能InverseLerp2D（XY作为的Point2D，X0Y0作为的Point2D，X1Y0作为的Point2D，X0Y1作为的Point2D，X1Y1作为的Point2D）作为的Point2D
    ' 去做 ？？？？
端功能
 

解决方案

要简化，让我们首先考虑只是一个单一的intepolated值的以Z 的。照片 假设四个值的以Z 的<分> 00 ，以Z 的<分> 01 ，以Z 的<分> 10 < /分>，以Z 的<分> 10 ，和两个权重的是W 的<子> 0 和是W 的<子> 1 施加到第一和第二索引，给

以Z 的<子> 0 = 以Z 的<分> 00 + 是W 的<子> 0 < /分>＆倍; （以Z 的<分> 10 - 以Z 的<分> 00 ）
 以Z 的<分> 1 = 以Z 的<分> 01 + 是W 的<子> 0 ＆倍; （以Z 的<分> 11 - 以Z 的<分> 01 ）

最后

以Z = 以Z 的<子> 0 + 是W 的<分> 1 ＆倍; （以Z 的<分> 1 - 以Z 的<子> 0 ）
 ＆NBSP;＆NBSP; = 以Z 的<分> 00 + 是W 的<子> 0 ＆倍; （以Z 的<分> 10 - 以Z 的<分> 00 ）+ 是W 的<分> 1 ＆倍; （以Z 的<分> 01 - 以Z 的<分> 00 ）+ 是W 的<分> 1 ＆倍; 是W 的<子> 0 ＆倍; （以Z 的<分> 11 - 以Z 的<分> 10 - 以Z 的<分> 01 + 以Z 的<分> 00 ）

因此​​，对于您的问题，您将有反转二维二次方程

X = X 的<分> 00 + 是W 的<子> 0 ＆倍; （ X 的<分> 10 - X 的<分> 00 ）+ 是W 的<分> 1 ＆倍; （ X 的<分> 01 - X 的<分> 00 ）+ 是W 的<分> 1 ＆倍; 是W 的<子> 0 ＆倍; （ X 的<分> 11 - X 的<分> 10 - X 的<分> 01 + X 的<分> 00 ）
 是 = 是的<分> 00 + 是W 的<子> 0 ＆倍; （是的<分> 10 - 是的<分> 00 ）+ 是W 的<分> 1 ＆倍; （是的<分> 01 - 是的<分> 00 ）+ 是W 的<分> 1 ＆倍; 是W 的<子> 0 ＆倍; （是的<分> 11 - 是的<分> 10 - 是的<分> 01 + 是的<分> 00 ）

不幸的是，没有一个简单的公式来恢复的是W 的<子> 0 和是W 的<分> 1 的 X 和是的。你可以，但是，把它当作一个非线性最小二乘问题，并尽量减少

（ X _是W 的（是W 的<子> 0 ，是W 的<分> 1 ） - X ） ² +（是_是W 的（是W 的<子> 0 ，是W 的<分> 1 ） - 是的） ²

，你可以有效地与 Levenberg-Marquardt算法做。

修改：再思考

它发生，我认为你可能会得到满意的插值（ X 的是的）到（是W 的<子> 0 ，瓦特的<子> 1 ），而不是实际的逆。这将是在这个意义上不太准确该转（FWD（的瓦特的<子> ，瓦特的<子> ）0 1）将可能进一步从（是W 的<子> 0 ，是W 的<分> 1 ）比实际的倒数。

这是你插了一个不规则的网状而不是一个常规电网的事实将会使这是一个棘手的命题。理想情况下，你应该加入你的（ X 的是的）点不重叠的三角形和使用的重心坐标线性插值。
对于数值稳定性应避免浅，尖尖的三角形。幸运的是， Delaunay三角 satifies这一要求，而不是的太的困难构建在两个维度

如果您希望您的反向插补采取类似的形式，以你的前进插值可以使用基函数

1
X 的
是的
X 的＆倍; 是的

和计算系数的在的<子>的我的， B 的<子>的我的< /分>， C 的<子>的我的和ð的<子>的我的（我的等于0或1），使得

是W 的<子> 0 = 在的<子> 0 + B 的<子> 0 < /分>＆倍; X + C 的<子> 0 ＆倍; 是 + ð的<子> 0 ＆倍; X 的＆倍; 是的
是W 的<分> 1 = 在的<分> 1 + B 的<分> 1 ＆倍; X + C 的<分> 1 ＆倍; 是 + ð的<分> 1 ＆倍; X 的＆倍; 是的

通过替代的的 X 的是的是W 的<子> 0 和 W精品相关已知值 <子> 1 你会得到四个同步的线性的方程为每个的是W 的，你可以解决得到其系数。
理想情况下，你应该用一个数值稳定的矩阵求逆算法，可以应付近奇异矩阵（如 SVD ） ，但你的可以的能逃脱高斯消元如果你再小心。

抱歉，我不能给你任何简单的选择，但我恐怕，这确实是一个比较棘手的问题！

Bilinear interpolation is trivial to compute. But I need an algorithm that does the INVERSE operation. (algorithm will be useful to me in pseudo-code, or any widely-used computer language)

For example, here is a Visual Basic implementation of bilinear interpolation.

' xyWgt ranges (0..1) in x and y. (0,0) will return X0Y0,
(0,1) will return X0Y1, etc.
' For example, if xyWgt is relative location within an image,
' and the XnYn values are GPS coords at the 4 corners of the image,
' The result is GPS coord corresponding to xyWgt.
' E.g. given (0.5, 0.5), the result will be the GPS coord at center of image.
Public Function Lerp2D(xyWgt As Point2D, X0Y0 As Point2D, X1Y0 As Point2D, X0Y1 As Point2D, X1Y1 As Point2D) As Point2D
    Dim xY0 As Point2D = Lerp(X0Y0, X1Y0, xyWgt.X)
    Dim xY1 As Point2D = Lerp(X0Y1, X1Y1, xyWgt.X)

    Dim xy As Point2D = Lerp(xY0, xY1, xyWgt.Y)
    Return xy
End Function

where

' Weighted Average of two points.
Public Function Lerp(ByVal a As Point2D, ByVal b As Point2D, ByVal wgtB As Double) As Point2D
    Return New Point2D(Lerp(a.X, b.X, wgtB), Lerp(a.Y, b.Y, wgtB))
End Function

and

' Weighted Average of two numbers.
' When wgtB==0, returns a, when wgtB==1, returns b.
' Implicitly, wgtA = 1 - wgtB. That is, the weights are normalized.
Public Function Lerp(ByVal a As Double, ByVal b As Double, ByVal wgtB As Double) As Double
    Return a + (wgtB * (b - a))
End Function

---

In 1-D, I have determined the inverse function of Lerp:

' Calculate wgtB that would return result, if did Lerp(a, b, wgtB).
' That is, where result is, w.r.t. a and b.
' < 0 is before a, > 1 is after b.
Public Function WgtFromResult(ByVal a As Double, ByVal b As Double, ByVal result As Double) As Double

    Dim denominator As Double = b - a

    If Math.Abs(denominator) < 0.00000001 Then
        ' Avoid divide-by-zero (a & b are nearly equal).
        If Math.Abs(result - a) < 0.00000001 Then
            ' Result is close to a (but also to b): Give simplest answer: average them.
            Return 0.5
        End If
        ' Cannot compute.
        Return Double.NaN
    End If

    ' result = a + (wgt * (b - a))   =>
    ' wgt * (b - a) = (result - a)   =>
    Dim wgt As Double = (result - a) / denominator

    'Dim verify As Double = Lerp(a, b, wgt)
    'If Not NearlyEqual(result, verify) Then
    '    Dim test = 0    ' test
    'End If

    Return wgt
End Function

Now I need to do the same in 2-D:

' Returns xyWgt, which if given to Lerp2D, would return this "xy".
' So if xy = X0Y0, will return (0, 0). if xy = X1Y0, will return (1, 0), etc.
' For example, if 4 corners are GPS coordinates in corners of an image,
' and pass in a GPS coordinate,
' returns relative location within the image.
Public Function InverseLerp2D(xy As Point2D, X0Y0 As Point2D, X1Y0 As Point2D, X0Y1 As Point2D, X1Y1 As Point2D) As Point2D
    ' TODO ????
End Function

解决方案

To simplify, let's begin by just considering a single intepolated value z.
Assume four values z₀₀, z₀₁, z₁₀, z₁₀, and two weights w₀ and w₁ applied to the first and second index, giving

z₀ = z₀₀ + w₀ × (z₁₀ - z₀₀)
z₁ = z₀₁ + w₀ × (z₁₁ - z₀₁)

and finally

z = z₀ + w₁ × (z₁ - z₀)
   = z₀₀ + w₀ × (z₁₀ - z₀₀) + w₁ × (z₀₁ - z₀₀) + w₁ × w₀ × (z₁₁ - z₁₀ - z₀₁ + z₀₀)

So, for your problem you will have to invert a two dimensional quadratic equation

x = x₀₀ + w₀ × (x₁₀ - x₀₀) + w₁ × (x₀₁ - x₀₀) + w₁ × w₀ × (x₁₁ - x₁₀ - x₀₁ + x₀₀)
y = y₀₀ + w₀ × (y₁₀ - y₀₀) + w₁ × (y₀₁ - y₀₀) + w₁ × w₀ × (y₁₁ - y₁₀ - y₀₁ + y₀₀)

Unfortunately, there isn't a simple formula to recover w₀ and w₁ from x and y. You can, however, treat it as a non-linear least squares problem and minimise

(x_w(w₀,w₁) - x)² + (y_w(w₀,w₁) - y)²

which you can do efficiently with the Levenberg–Marquardt algorithm.

Edit: Further Thoughts

It has occurred to me that you might be satisfied with an interpolation from (x, y) to (w₀, w₁) rather than the actual inverse. This will be less accurate in the sense that rev(fwd(w₀, w₁)) will likely be further from (w₀, w₁) than the actual inverse.

The fact that you're interpolating over an irregular mesh rather than a regular grid is going to make this a trickier proposition. Ideally you should join up your (x, y) points with non-overlapping triangles and use barycentric coordinates to linearly interpolate.
For numerical stability you should avoid shallow, pointy triangles. Fortunately, the Delaunay triangulation satifies this requirement and isn't too difficult to construct in two dimensions.

If you would like your reverse interpolation to take a similar form to your forward interpolation you can use the basis functions

1
x
y
x × y

and compute coefficients a_i, b_i, c_i and d_i (i equal to 0 or 1) such that

w₀ = a₀ + b₀ × x + c₀ × y + d₀ × x × y
w₁ = a₁ + b₁ × x + c₁ × y + d₁ × x × y

By substituting the relevant known values of x, y, w₀ and w₁ you'll get four simultaneous linear equations for each w that you can solve to get its coefficients.
Ideally you should use a numerically stable matrix inversion algorithm that can cope with near singular matrices (e.g. SVD), but you may be able to get away with Gaussian elimination if you're careful.

Sorry I can't give you any simpler options, but I'm afraid that this really is a rather tricky problem!

这篇关于算法：如何计算双线性插值的逆？的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！