反转实数索引网格

Question

OpenCV的remap()使用实值索引网格，使用双线性插值对图像中的值网格进行采样，并将样本网格作为新图像返回。

准确地说，让：

A = an image 
X = a grid of real-valued X coords into the image. 
Y = a grid of real-valued Y coords into the image. 
B = remap(A, X, Y) 

然后，对于所有像素坐标I，J，

B[i, j] = A(X[i, j], Y[i, j]) 

当圆形括号符号A(x, y)表示使用双线性插值来求解的像素值图像A使用浮动值的坐标x和y。

我的问题是：给出一个索引网格X，Y，我怎么能产生 “逆格” X^-1，Y^-1使得：

X(X^-1[i, j], Y^-1[i, j]) = i 
Y(X^-1[i, j], Y^-1[i, j]) = j 

而且

X^-1(X[i, j], Y[i, j]) = i 
Y^-1(X[i, j], Y[i, j]) = j 

对于所有的整数像素坐标i, j？

FWIW，图像和索引图X和Y是相同的形状。但是，索引图X和Y没有先验结构。例如，它们不一定是仿射或刚性变换。它们甚至可能是不可逆的，例如如果X, Y将A中的多个像素映射到B中相同的精确像素坐标。我正在寻找一种方法的想法，如果存在的话，它会找到合理的逆映射。

该解决方案不需要基于OpenCV，因为我没有使用OpenCV，而是另一个具有remap()实现的库。尽管欢迎任何建议，但我特别热衷于“数学上正确”的事情，即如果我的映射M是完全可逆的，则该方法应该在机器精度的某个小幅度内找到完美的逆。

Answer 1

那么我只需要自己解决这个重新映射问题，我将概述我的解决方案。

鉴于X，Y为remap()函数执行以下操作：

B[i, j] = A(X[i, j], Y[i, j]) 

我计算Xinv，Yinv可由remap()函数反相的方法中使用：

A[x, y] = B(Xinv[x,y],Yinv[x,y]) 

首先我建立一个KD-Tree为2D点集{(X[i,j],Y[i,j]}所以我可以高效地找到N最近的邻居到一个给定的点(x,y).我使用欧几里德距离作为我的距离度量。我在GitHub上发现了很棒的C++ header lib for KD-Trees。

然后我通过A的网格中的所有(x,y)值循环并找到N = 5最近的邻居{(X[i_k,j_k],Y[i_k,j_k]) | k = 0 .. N-1}在我的点集中。

• 如果距离d_k == 0一些k然后Xinv[x,y] = i_k和Yinv[x,y] = j_k，否则......

• 使用Inverse Distance Weighting (IDW)计算的插补值：

  • 让重量w_k = 1/pow(d_k, p)（我用的p = 2）
  • Xinv[x,y] = (sum_k w_k * i_k)/(sum_k w_k)
  • Yinv[x,y] = (sum_k w_k * j_k)/(sum_k w_k)

注意，如果B是W x H图像然后X和Y是浮筒的W x H阵列。如果A是w x h图像，则Xinv和Yinv是w x h浮点数组。与图像和地图大小一致很重要。

工程就像一个魅力！我的第一个版本，我试图蛮力强制搜索，我甚至从来没有等待它完成。我切换到KD树，然后我开始获得合理的运行时间。如果我想花时间将它添加到OpenCV中，

以下第二张图片使用remap()来消除第一张图像中的镜头失真。第三张图像是反转过程的结果。

Answer 2

没有任何标准的方法来做到这一点OpenCV。

如果您正在寻找一个完整的现成解决方案，我不确定我能否提供帮助，但我至少可以描述一种几年前用于完成此任务的方法。

首先，您应该创建与源图像尺寸相同的重映射贴图。我为更简单的插值创建了更大尺寸的地图，并在最后一步将其裁剪为适当的尺寸。然后，你应该用先前重映射图中存在的值填充它们（并不那么难：只需迭代它们，如果地图坐标x和y放在图像的限制范围内，则将它们的行和列作为新的y和x，并放入旧的x和y列和新地图的行）。这是相当简单的解决方案，但它给出了相当好的结果。对于完美的一个，您应该使用插值方法和相邻像素将旧的x和y插值为整数值。

之后，您应该手动重新映射像素颜色，或者完全使用像素坐标填充重新映射图并使用OpenCV中的版本。

您将遇到一个颇具挑战性的任务：您应该在空白区域内插像素。换句话说，您应该根据这些距离将距离与最接近的非零像素坐标和混合颜色（如果重新映射颜色）或坐标（如果继续进行完整的地图计算）分数相乘。实际上，对于线性插值也并不困难，并且您甚至可以在OpenCV github page中查看remap()实现。对于NN插值，我会简单得多 - 只需要最近邻居的颜色/坐标。

最后的任务是重新映射像素区域的边界外的区域。 OpenCV中的算法也可以用作参考。

Answer 3

从我所了解的你有一个原始图像和一个变换后的图像，你希望恢复已经应用的变换的本质而不知道它，但假设它是明智的，比如旋转或鱼 - 眼睛扭曲。

我会尝试的是阈值图像将其转换为二进制，在索引图像和普通图像。然后尝试识别对象。大多数映射至少会保留连通性和欧拉数，大多数指数中最大的对象仍然是平原上最大的对象。

然后花点时间为您匹配的图片/索引对，看看您是否可以移除平移，旋转和缩放。这给你几张反向贴图，然后你可以尝试缝合在一起。 （如果变换不是简单的，那么很难，但重构任何变换的一般问题无法解决）。

Answer 4

如果您的映射来自单应H，您可以反转H并直接使用cv :: initUndistortRectifyMap（）创建反向映射。

例如在python：

约initUndistortRectifyMap（）

import numpy as np. 
map_size =() # fill in your map size 
H_inv = np.linalg.inv(H) 
map1, map2 = cv2.initUndistortRectifyMap(cameraMatrix=np.eye(3), distCoeffs=np.zeros(5), R=H_inv, newCameraMatrix=np.eye(3), size=map_size, m1type=cv2.CV_32FC1) 

在OpenCV文档状态： “功能实际上建立用于所使用的重映射（逆映射算法）地图也就是说，对于每个像素（U，V。 ）在目标图像中，函数计算源图像中的相应坐标。“

如果您刚给出地图，您必须自己做。 Hoewever，新地图坐标的插值并不是微不足道的，因为一个像素的支持区域可能非常大。

这是一个简单的python解决方案，通过点对点映射来反转地图。这可能会留下一些未指定的坐标，而其他坐标则会多次更新。所以地图上可能会出现漏洞。

这里是一个小的Python程序演示了这两种方法：这里

import cv2 
import numpy as np 


def invert_maps(map_x, map_y): 
    assert(map_x.shape == map_y.shape) 
    rows = map_x.shape[0] 
    cols = map_x.shape[1] 
    m_x = np.ones(map_x.shape, dtype=map_x.dtype) * -1 
    m_y = np.ones(map_y.shape, dtype=map_y.dtype) * -1 
    for i in range(rows): 
     for j in range(cols): 
      i_ = round(map_y[i, j]) 
      j_ = round(map_x[i, j]) 
      if 0 <= i_ < rows and 0 <= j_ < cols: 
       m_x[i_, j_] = j 
       m_y[i_, j_] = i 
    return m_x, m_y 


def main(): 
    img = cv2.imread("pigeon.png", cv2.IMREAD_GRAYSCALE) 

    # a simply rotation by 45 degrees 
    H = np.array([np.sin(np.pi/4), -np.cos(np.pi/4), 0, np.cos(np.pi/4), np.sin(np.pi/4), 0, 0, 0, 1]).reshape((3,3)) 
    H_inv = np.linalg.inv(H) 
    map_size = (img.shape[1], img.shape[0]) 

    map1, map2 = cv2.initUndistortRectifyMap(cameraMatrix=np.eye(3), distCoeffs=np.zeros(5), R=H, newCameraMatrix=np.eye(3), size=map_size, m1type=cv2.CV_32FC1) 
    map1_inv, map2_inv = cv2.initUndistortRectifyMap(cameraMatrix=np.eye(3), distCoeffs=np.zeros(5), R=H_inv, newCameraMatrix=np.eye(3), size=map_size, m1type=cv2.CV_32FC1) 
    map1_simple_inv, map2_simple_inv = invert_maps(map1, map2) 

    img1 = cv2.remap(src=img, map1=map1, map2=map2, interpolation=cv2.INTER_LINEAR) 
    img2 = cv2.remap(src=img1, map1=map1_inv, map2=map2_inv, interpolation=cv2.INTER_LINEAR) 
    img3 = cv2.remap(src=img1, map1=map1_simple_inv, map2=map2_simple_inv, 
           interpolation=cv2.INTER_LINEAR) 

    cv2.imshow("Original image", img) 
    cv2.imshow("Mapped image", img1) 
    cv2.imshow("Mapping forth and back with H_inv", img2) 
    cv2.imshow("Mapping forth and back with invert_maps()", img3) 
    cv2.waitKey(0) 


if __name__ == '__main__': 
    main() 

Answer 5

OP。我想我找到了答案。我还没有实现它，如果有人提出了一个不太方便的解决方案（或者发现这个问题有什么问题），我会选择他们的答案。

问题陈述

设A为源图像，B是目标图像，和M是从A的坐标中映射到B的COORDS，即：

B[k, l, :] == A(M[k, l, 0], M[k, l, 1], :) 
for all k, l in B's coords. 

...其中方括号表示使用整数索引进行数组查找，圆括号表示使用浮点索引进行双线性插值查找。我们重述以上使用更经济的表示法：

B = A(M) 

我们希望找到映射乙回到A为最好的，因为有可能一个逆映射N：

Find N s.t. A \approx B(N) 

问题可以不用说参考A或B：

Find N = argmin_N || M(N) - I_n || 

...其中||*||表示Frobenius范数，和I_n是恒等映射具有相同的尺寸N，即图中：

I_n[i, j, :] == [i, j] 
for all i, j 

天真的

如果M的值都是整数，M是同构的，那么你就可以构建ñ直接为：

N[M[k, l, 0], M[k, l, 1], :] = [k, l] 
for all k, l 

或在我们的简化表示法：

N[M] = I_m 

...其中I_M是具有相同尺寸的恒等映射为M.

有两个问题：

1. M不是一个同构，所以上述将在n在N [I，J，：]离开的 “洞” 为任何[I，J]不间M中的值。
2. M的值是浮点坐标[i，j]，而不是整数坐标。对于浮点值i，j，我们不能简单地为双线性内插数量N（i，j，:)分配一个值。为了达到等效效果，我们必须设置[i，j]的四个边角N [floor（i），floor（j），]，N [floor（i），ceil（j）， ：]，N [ceil（i），floor（j），]，N [ceil（i），ceil（j），]：使得内插值N（i，j，:)等于期望值[ K，L]，对于所有像素映射[I，J] - 在M.
> [K，1]

解

构建空N作为漂浮的3D张量：

N = zeros(size=(A.shape[0], A.shape[1], 2)) 

对于A坐标空间中的每个坐标[i，j]，执行：

1. 查找M中[i，j]所在的A坐标的2×2网格。 计算单应矩阵H，将这些A坐标映射到它们对应的B坐标（由2x2网格的像素索引给出）。
2. 集N [I，J，：] = MATMUL（H，[I，J]）

这里的潜在昂贵的步骤将是在A-坐标的2×2网格步骤1中的搜索包围[i，j]的M。蛮力搜索将使得整个算法O（n * m）其中n是A中像素的数量，m是B中像素的数量。

为了将其减少到O（n）而是在每个A坐标四边形内运行扫描线算法以识别其包含的所有整数值坐标[i，j]。这可以预先计算为将整数值的A坐标[i，j]映射到其环绕的四边形B坐标[k，l]的左上角的散列图。