图片拼接是一个很常见的应用，本文总结一下通过OpenCV实现多张图的全景拼接

1.获取特征点

  • 方法一： SIFT算法
1
2

sift = cv2.xfeatures2d.SIFT_create(nfeatures=999999, nOctaveLayers=9 ,contrastThreshold=0.0001, edgeThreshold=2)
kp1, des1 = sift.detectAndCompute(img1, None)
  • 方法二： SURF算法
1
2

surf = cv2.xfeatures2d.SURF_create(hessianThreshold=100, nOctaves=1, nOctaveLayers=1, extended=False, upright=True)
kp1, des1 = surf.detectAndCompute(img1, None)

参数可调。 实践证明SIFT的准确度明显高于SURF。

2.计算单应性矩阵

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25

matcher = cv2.FlannBasedMatcher()
# matcher = cv2.BFMatcher()

raw_matches = matcher.knnMatch(np.asarray(des1, np.float32), np.asarray(des2, np.float32), k=2)

good_points = []
good_matches = []
for m1, m2 in raw_matches:
    if m1.distance < 0.8 * m2.distance:
        good_points.append((m1.trainIdx, m1.queryIdx))
        good_matches.append([m1])

if len(good_points) > 12:
    image1_kp = np.float32([kp1[i].pt for (_, i) in good_points])
    image2_kp = np.float32([kp2[i].pt for (i, _) in good_points])

    H21, status21 = cv2.findHomography(image2_kp, image1_kp, cv2.RANSAC, 5.0)
    lastmatch21 = len([i[0] for i in status21 if i[0] == 1])

    H12, status12 = cv2.findHomography(image1_kp, image2_kp, cv2.RANSAC, 5.0)
    lastmatch12 = len([i[0] for i in status12 if i[0] == 1])

    if lastmatch21 < 10 or lastmatch12 < 10:
        raise Exception("lastmatch less than 10.")

H21表示将img2拼接到img1时img2的变换矩阵 lastmatch是最终匹配特征点点的个数，如果特征点不够将无法拼接。

上述示例用到了findHomography进行透视变换，也可以用estimateAffine2D进行仿射变换，参数一样，不同的是estimateAffine2D得到的是23的矩阵， 可填充为33（H[2, 0] = 0.0, H[2, 1] = 0, H[2, 2] = 1.0)

3.拼接

两张图的场景

完整代码参考: 两张图的拼接

1
2
3
4
5
6
7
8

h, w, _ = img1.shape
panorama = np.zeros((h, w * 2, 3))
panorama[0:h, 0:w, :] = img1
wrap = cv2.warpPerspective(img2, H21, (w * 2, h))
label = (cv2.cvtColor(wrap[:h, :w * 2, :], cv2.COLOR_BGR2GRAY) > 0)
for i in range(panorama.shape[2]):
    panorama[:h, :w * 2, i] = wrap[:h, :w * 2, i] * (label > 0) + panorama[:h, :w * 2, i] * (label < 1)
cv2.imwrite('stitch_im.jpg', panorama)

三张以上图的场景

1. 首先选定基准图，然后把剩余的图的homography矩阵变换

有图集[img1, img2, img3, img4, img5]

单应性矩阵： H1: [H12, H23, H34, H45] H2: [H21, H32, H43, H54]

选定一张基准图img2, 需要计算变换矩阵列表： [H13, H23, H, H43, H53]

其中: H是单位矩阵 H13 = np.dot(H12, H23)

其余类似，如有更多图，可多级连乘（e.g.: H85 = H87 · H76 · H65）

矩阵累积可以用递归的函数实现，示例：

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16

def multiply_homos12(homos):
    if len(homos) == 1:
        return homos[0]
    elif len(homos) == 2:
        return np.dot(homos[0], homos[1])
    else:
        return np.dot(homos[0], multiplyhomos1_2(homos[1:]))


def multiply_homos21(homos):
    if len(homos) == 1:
        return homos[0]
    elif len(homos) == 2:
        return np.dot(homos[1], homos[0])
    else:
        return np.dot(multiplyhomos2_1(homos[1:]), homos[0])

然后每个矩阵都要点乘Hpivot

1

homos = [np.dot(Hpivot, X) for X in homos_]

Hpivot是基准图对应的变换矩阵，只需要计算基准图距离最左边的距离，方法如下

1
2
3
4
5
6
7

H_ = multiply_homos21(homos12[:mid])
a = np.dot(H_, np.array([0.0, 0.0, 1.0]).reshape(3, 1)) # 基准图左上点的右移距离
b = np.dot(H_,np.array([0.0, float(h), 1.0]).reshape(3, 1)) # 基准图左下点的右移距离
distance = abs(min(a[0, 0] / a[2, 0], b[0, 0] / b[2, 0]))
Hpivot = np.array([[1.0, 0.0, distance], 
                   [0.0, 1.0, 0.0],
                   [0.0, 0.0, 1.0]])

2. 计算变换后的图 （参考两张图的场景）

3. 将所有变换后的图拼接在一起 （参考两张图的场景） 边界处可进行融合算法

融合渐变函数示例

1
2
3
4
5
6
7

# h是图片高度，w_是融合区域左右宽度
temp1 = np.ones((h, w_), dtype=float)
temp2 = np.ones((h, w_), dtype=float)
for i in range(w_):
    ex = 1 - 2.0 * float(i) / float(w_ - 1)
    temp1[:, i] = 1.0 / (1.0 + math.pow(math.e, -40 * ex))
    temp2[:, i] = 1.0 - temp1[:, i]