devbox
你好赵伟
这个屌丝很懒，什么也没留下！
关注作者
热门标签
热门文章
当前位置:   article > 正文

适合新手--特征点匹配及消除误匹配点-python语言-适用于深度学习提取的特征点和描述符（超简单）_如何剔除误匹配

作者：你好赵伟 | 2024-08-06 07:51:56

赞

踩

如何剔除误匹配

这里用RANSAC来过滤误匹配

本代码主要是实现图像误匹配，理论部分及参考文章移步：

计算机视觉：RANSAC剔除基础矩阵F错误匹配(Python实现)_Sunrise永不言弃的博客-CSDN博客

以下为代码部分，只需要修改350-366行就可以了，自己根据实际情况修改

  1.  
  2. import os
  3. import cv2
  4. import numpy as np
  5. import random
  6. os.environ["KMP_DUPLICATE_LIB_OK"] = "TRUE"
  7. if int(cv2.__version__[0]) < 3:  # pragma: no cover
  8.     print('Warning: OpenCV 3 is not installed')
  9.  
  10. #ransac
  11. def compute_fundamental(x1, x2):
  12.     n = x1.shape[1]
  13.     if x2.shape[1] != n:
  14.         raise ValueError("Number of points don't match.")
  15.  
  16.     # build matrix for equations
  17.     A = np.zeros((n, 9))
  18.     for i in range(n):
  19.         A[i] = [x1[0, i] * x2[0, i], x1[0, i] * x2[1, i], x1[0, i] * x2[2, i],
  20.                 x1[1, i] * x2[0, i], x1[1, i] * x2[1, i], x1[1, i] * x2[2, i],
  21.                 x1[2, i] * x2[0, i], x1[2, i] * x2[1, i], x1[2, i] * x2[2, i]]
  22.  
  23.     # compute linear least square solution
  24.     U, S, V = np.linalg.svd(A)
  25.     F = V[-1].reshape(3, 3)
  26.  
  27.     # constrain F
  28.     # make rank 2 by zeroing out last singular value
  29.     U, S, V = np.linalg.svd(F)
  30.     S[2] = 0
  31.     F = np.dot(U, np.dot(np.diag(S), V))
  32.  
  33.     return F / F[2, 2]
  34.  
  35.  
  36. def compute_fundamental_normalized(x1, x2):
  37.     """    Computes the fundamental matrix from corresponding points
  38.         (x1,x2 3*n arrays) using the normalized 8 point algorithm.
  39.         从相应点计算基本矩阵
  40.          （x1，x2 3*n数组）使用标准化的8点算法
  42.         """
  43.  
  44.     n = x1.shape[1]
  45.     if x2.shape[1] != n:
  46.         raise ValueError("Number of points don't match.")
  47.  
  48.     # normalize image coordinates
  49.     #标准化图像坐标
  50.     x1 = x1 / x1[2]
  51.     mean_1 = np.mean(x1[:2], axis=1)
  52.     S1 = np.sqrt(2) / np.std(x1[:2])
  53.     T1 = np.array([[S1, 0, -S1 * mean_1[0]], [0, S1, -S1 * mean_1[1]], [0, 0, 1]])
  54.     x1 = np.dot(T1, x1)
  55.  
  56.     x2 = x2 / x2[2]
  57.     mean_2 = np.mean(x2[:2], axis=1)
  58.     S2 = np.sqrt(2) / np.std(x2[:2])
  59.     T2 = np.array([[S2, 0, -S2 * mean_2[0]], [0, S2, -S2 * mean_2[1]], [0, 0, 1]])
  60.     x2 = np.dot(T2, x2)
  61.  
  62.     # compute F with the normalized coordinates
  63.     #用归一化坐标计算F
  64.     F = compute_fundamental(x1, x2)
  65.     # print (F)
  66.     # reverse normalization
  67.     F = np.dot(T1.T, np.dot(F, T2))
  68.  
  69.     return F / F[2, 2]
  70.  
  71. def randSeed(good, num = 8):
  72.     '''
  73.     :param good: 初始的匹配点对
  74.     :param num: 选择随机选取的点对数量
  75.     :return: 8个点对list
  76.     '''
  77.     eight_point = random.sample(good, num)
  78.     return eight_point
  79.  
  80. def PointCoordinates(eight_points, keypoints1, keypoints2):  #随机八点求  eight_points是匹配好的点
  81.     '''
  82.     :param eight_points: 随机八点
  83.     :param keypoints1: 点坐标
  84.     :param keypoints2: 点坐标
  85.     :return:8个点
  86.     '''
  87.     x1 = []
  88.     x2 = []
  89.     tuple_dim = [1.,]             #这里改成了np类型类型
  90.     #print(type(eight_points))   类型是list
  91.     #print(eight_points)
  92.  
  93.     fin_list_1=[]           #初始化保存8点的坐标
  94.     fin_list_2 = []
  95.     for i in eight_points:
  96.         #        pts[0, :] = ys
  97.         #        pts[1, :] = xs
  98.         #print(type(keypoints1))   #元组
  99.  
  100.         #tuple_x1 = keypoints1[i][:2] + tuple_dim
  101.         #tuple_x2 = keypoints2[i][:2] + tuple_dim
  102.         idex_0=i[0].queryIdx                  #查询的描述符索引
  103.         idex_1 = i[0].trainIdx                #目标的描述符索引
  104.         #print("xxxxxxx",idex_0)
  105.         tuple_x1 = keypoints1[idex_0][:2]           #得到关键点1的坐标  关键点转置后是N*3的矩阵，取前2维度
  106.         tuple_x2 = keypoints2[idex_1][:2]           #得到关键点2的坐标
  107.         #print("坐标是",tuple_x1,tuple_x2)           目前输出[257. 418.] [220. 464.]
  108.  
  109.         # 转化为其次，最后+1个维度的那种
  110.         new_tuple_x1=np.append(tuple_x1,tuple_dim)
  111.         new_tuple_x2 = np.append(tuple_x2, tuple_dim)
  112.         #print("xxxxxxx",new_tuple_x1,new_tuple_x2)    #  [663. 445.   1.] [269. 246.   1.]
  113.  
  114.         fin_list_1=np.append(fin_list_1,new_tuple_x1)
  115.         fin_list_2 = np.append(fin_list_2, new_tuple_x2)
  116.         #print(type(fin_list_1))   #fin_list_1 np类型的，24个数的集合
  117.         #arr = np.arange(24).reshape(3, 8)
  118.  
  119.         #原来的
  120.         #tuple_x1 = keypoints1[i[0].queryIdx].pt + tuple_dim
  121.         #tuple_x2 = keypoints2[i[0].trainIdx].pt + tuple_dim
  122.         #x1.append(tuple_x1)
  123.         #x2.append(tuple_x2)
  124.     #print("xxxxxxx",fin_list_1.shape)
  125.     x1 = fin_list_1.reshape(-1, 3)     #列是3 ，行数计算
  126.     x2 = fin_list_2.reshape(-1, 3)     #列是3 ，行数计算
  127.     #print(x1)
  128.     #print("xxxxxxxx",x1.shape)
  129.     return np.array(x1, dtype=float), np.array(x2, dtype=float)
  130.  
  131.  
  132. def ransac(good, keypoints1, keypoints2, confidence,iter_num):
  133.     Max_num = 0
  134.     good_F = np.zeros([3,3])
  135.     inlier_points = []
  136.     for i in range(iter_num):
  137.         eight_points = randSeed(good)
  138.         x1,x2 = PointCoordinates(eight_points, keypoints1, keypoints2)
  139.         F = compute_fundamental_normalized(x1.T, x2.T)
  140.         num, ransac_good = inlier(F, good, keypoints1, keypoints2, confidence)
  141.         if num > Max_num:
  142.             Max_num = num
  143.             good_F = F
  144.             inlier_points = ransac_good
  145.     print(Max_num, good_F)
  146.     return Max_num, good_F, inlier_points
  147.  
  148.  
  149. def computeReprojError(x1, x2, F):
  150.     """
  151.     计算投影误差
  152.     """
  153.     ww = 1.0/(F[2,0]*x1[0]+F[2,1]*x1[1]+F[2,2])
  154.     dx = (F[0,0]*x1[0]+F[0,1]*x1[1]+F[0,2])*ww - x2[0]
  155.     dy = (F[1,0]*x1[0]+F[1,1]*x1[1]+F[1,2])*ww - x2[1]
  156.     return dx*dx + dy*dy
  157.  
  158. def inlier(F,good, keypoints1,keypoints2,confidence):   #找到内点
  159.     num = 0
  160.     ransac_good = []
  161.     x1, x2 = PointCoordinates(good, keypoints1, keypoints2)
  162.     for i in range(len(x2)):
  163.         line = F.dot(x1[i].T)
  164.         #在对极几何中极线表达式为[A B C],Ax+By+C=0,  方向向量可以表示为[-B,A]
  165.         line_v = np.array([-line[1], line[0]])
  166.         err = h = np.linalg.norm(np.cross(x2[i,:2], line_v)/np.linalg.norm(line_v))
  167.         # err = computeReprojError(x1[i], x2[i], F)
  168.         if abs(err) < confidence:
  169.             ransac_good.append(good[i])
  170.             num += 1
  171.     return num, ransac_good
  172.  
  173.  
  174.  
  175.  
  176.  
  177.  
  178.  
  179. def nn_match_two_way(desc1, desc2, nn_thresh):
  180.     assert desc1.shape[0] == desc2.shape[0]
  181.     if desc1.shape[1] == 0 or desc2.shape[1] == 0:
  182.         return np.zeros((3, 0))
  183.     if nn_thresh < 0.0:
  184.         raise ValueError('\'nn_thresh\' should be non-negative')
  185.     # Compute L2 distance. Easy since vectors are unit normalized.
  186.     #计算L2距离容易，因为向量已归一化
  187.     dmat = np.dot(desc1.T, desc2)
  188.     dmat = np.sqrt(2 - 2 * np.clip(dmat, -1, 1))
  189.  
  190.     # Get NN indices and scores.
  191.     #获取NN索引和分数。
  192.     idx = np.argmin(dmat, axis=1)                 #最近邻索引   argmin返回这个数组中最小的索引
  193.     scores = dmat[np.arange(dmat.shape[0]), idx]  #最近邻分数
  194.     # Threshold the NN matches.
  195.     #NN匹配的阈值
  196.     keep = scores < nn_thresh
  197.     # Check if nearest neighbor goes both directions and keep those.
  198.     #检查最近的邻居是否两个方向并保留这些方向
  199.     idx2 = np.argmin(dmat, axis=0)
  200.     keep_bi = np.arange(len(idx)) == idx2[idx]
  201.     keep = np.logical_and(keep, keep_bi)
  202.     idx = idx[keep]
  203.     scores = scores[keep]
  204.     # Get the surviving point indices.
  205.     #获取幸存点索引
  206.     m_idx1 = np.arange(desc1.shape[1])[keep]  #keep是布尔类型的 ture或者false  ,输出这样的[[0 1 2 3 4]]
  207.     m_idx2 = idx
  208.     # Populate the final 3xN match data structure.
  209.     #填充最终的3XN匹配数据结构
  210.     matches = np.zeros((3, int(keep.sum())))
  211.     matches[0, :] = m_idx1
  212.     #print(matches[0])  # 5.  22.  40.  68.  85.  94.  96. 109. 110. 126. 132. 135
  213.     matches[1, :] = m_idx2
  214.     #print(matches[1])  #56.  11. 369. 286.  94
  215.     matches[2, :] = scores   #最近邻分数
  216.     #print(matches[2])  #0.64095849 0.58723611 0.68436664 0.55057442
  217.     return matches    #matches保存的就是点对的对应关系
  218.  
  219.  
  220. def match_descriptors(kp1, desc1, kp2, desc2):
  221.     # Match the keypoints with the warped_keypoints with nearest neighbor search
  222.     bf = cv2.BFMatcher(cv2.NORM_L2, crossCheck=True)
  223.     ##是一个暴力匹配的对象，取desc1中一个描述子，再与desc2中的所有描述子计算欧式距离
  224.     matches = bf.match(desc1, desc2)
  225.  
  226.     ##上一行的返回值类似元组的集合（i，j）代表第一个集合的第i个点的最佳匹配是第二个集合的第j个点
  227.     matches_idx = np.array([m.queryIdx for m in matches])
  228.     m_kp1 = [kp1[idx] for idx in matches_idx]
  229.     matches_idx = np.array([m.trainIdx for m in matches])
  230.     m_kp2 = [kp2[idx] for idx in matches_idx]
  231.  
  232.     ####m_kp1是第一张图片的特征点，m_kp2是第二张图片的特征点，此时它们已经一一对应了（至少是最对应的，距离最小的
  233.     return m_kp1, m_kp2, matches
  234.  
  235.  
  236. def showpoint(img, ptx):
  237.     for i in range(ptx.shape[1]):
  238.         x = int(round(ptx[0, i]))
  239.         y = int(round(ptx[1, i]))
  240.         # if x>20 and y>20 and x<640 and y <450:
  241.         #   None
  242.         cv2.circle(img, (x, y), 3, color=(255, 0, 0))
  243.     return img
  244. ## test def drawMatches(img1, kp1, img2, kp2, matches)
  245. def drawMatches(img1, kp1, img2, kp2, matches):
  246.     """
  247.     My own implementation of cv2.drawMatches as OpenCV 2.4.9
  248.     does not have this function available but it's supported in
  249.     OpenCV 3.0.0
  250.     This function takes in two images with their associated
  251.     keypoints, as well as a list of DMatch data structure (matches)
  252.     that contains which keypoints matched in which images.
  253.     An image will be produced where a montage is shown with
  254.     the first image followed by the second image beside it.
  255.     Keypoints are delineated with circles, while lines are connected
  256.     between matching keypoints.
  257.     img1,img2 - Grayscale images #灰度图像
  258.     kp1,kp2 - Detected list of keypoints through any of the OpenCV keypoint
  259.               detection algorithms
  260.     #通过任何OPENCV关键点检测到的关键点列表
  261.                检测算法
  262.     matches - A list of matches of corresponding keypoints through any
  263.               OpenCV keypoint matching algorithm
  264.     #通过任何OPENCV关键点匹配算法
  265.     """
  266.  
  267.     # Create a new output image that concatenates the two images together
  268.     #创建一个新的输出图像，将两个图像串联在一起
  269.     # (a.k.a) a montage #剪辑
  270.     rows1 = img1.shape[0]
  271.     cols1 = img1.shape[1]
  272.     rows2 = img2.shape[0]
  273.     cols2 = img2.shape[1]
  274.  
  275.  
  276.     #将第一个图像放在左侧
  277.     #i1 = np.dstack([img1, img1, img1])
  278.     #i2 = np.dstack([img2, img2, img2])
  279.     #i1 = cv2.imread("D:\essay\SuperPoint\SuperPointPretrainedNetwork-master/assets/test/1.png")
  280.     #i2 = cv2.imread("D:\essay\SuperPoint\SuperPointPretrainedNetwork-master/assets/test/2.png")
  281.  
  282.     #cv2.imshow("image1", i1)
  283.     #cv2.imshow("image2", i1)
  284.     #np.hstack将参数元组的元素数组按水平方向进行叠加
  285.     #就是图像并排放一起
  286.     out = np.hstack([img1, img2])
  287.     #cv2.imshow("image2", out)
  288.     print("输出尺寸", out.shape)
  289.     # Place the next image to the right of it
  290.     # out[0:480,640:1280] = np.dstack([img2, img2, img2])
  291.  
  292.     # For each pair of points we have between both images
  293.     # draw circles, then connect a line between them
  294.     #对于每对点，我们在这两个图像之间都有
  295.     #绘制圆圈，然后在它们之间连接一条线
  296.     for i in range(matches.shape[1]):   #[1]点的个数
  297.         # Get the matching keypoints for each of the images
  298.         #获取每个图像的匹配关键点
  299.  
  300.         img1_idx = matches[0, i]
  301.         img2_idx = matches[1, i]
  302.         x11 = int(img1_idx)
  303.         y11 = int(img1_idx)
  304.         x22 = int(img2_idx)
  305.         y22 = int(img2_idx)
  306.  
  307.         # x - columns    x列
  308.         # y - rows       y行
  309.         x1 = kp1[0, x11]
  310.         y1 = kp1[1, y11]
  311.         x2 = kp2[0, x22]
  312.         y2 = kp2[1, y22]
  313.  
  314.         # Draw a small circle at both co-ordinates
  315.         # radius 4
  316.         # colour blue
  317.         # thickness = 1
  318.         # a = np.random.randint(0, 256)  #返回一个随机整数值
  319.         # b = np.random.randint(0, 256)
  320.         # c = np.random.randint(0, 256)
  321.  
  322.         a = np.random.randint(5, 7)     #B 蓝色返回一个随机整数值
  323.         b = np.random.randint(5, 7)     #G 绿色
  324.         c = np.random.randint(244, 246) #R 红色
  325.  
  326.         cv2.circle(out, (int(np.round(x1)), int(np.round(y1))), 2, (a, b, c), 1)  # 画圆，cv2.circle()参考官方文档
  327.         cv2.circle(out, (int(np.round(x2) + cols1), int(np.round(y2))), 2, (a, b, c), 1)
  328.  
  329.         # Draw a line in between the two points
  330.         # thickness = 1
  331.         # colour blue
  332.         cv2.line(out, (int(np.round(x1)), int(np.round(y1))), (int(np.round(x2) + cols1), int(np.round(y2))),
  333.                  (a, b, c),2, shift=0)  # 画线，cv2.line()参考官方文档
  334.  
  335.     # Also return the image if you'd like a copy
  336.     return out
  337.  
  338.  
  339.  
  340. if __name__ == '__main__':
  341.     font = cv2.FONT_HERSHEY_DUPLEX  # 设置可视化字体
  342.     font_clr = (255, 255, 255)
  343.     font_pt = (4, 12)
  344.     font_sc = 0.4
  345.  
  346.  
  347.  
  348.     #读图
  349.     #第一，二图的特征点和描述符
  350.     #.d2-net是NPZ文件，你换成自己的就就好了
  351.     data_1 = np.load("D:\essay\SuperPoint\superpoint\pictuer\style_image\match/11.jpg.d2-net")
  352.     data_2 = np.load("D:\essay\SuperPoint\superpoint\pictuer\style_image\match/12.jpg.d2-net")
  353.  
  354.     #读入2张图像
  355.     imgx_1 = cv2.imread('D:\essay\SuperPoint\superpoint\pictuer\style_image\match/11.jpg', 1)
  356.     imgx_2 = cv2.imread('D:\essay\SuperPoint\superpoint\pictuer\style_image\match/12.jpg', 1)
  357.  
  358.     #可以用这个看下你的npz文件的属性
  359.     # print(data.files)
  360.  
  361.     #我这里是(3473, 3)维，N是3473，不同的网络提取的特征点不一样。如果你的不是这样，你可以转置下  例如pts=pts.T
  362.     #如果你的是(3473, 2)也是没问题的
  363.     pts = data_1['keypoints']  # np类型  (3473, 3)  3473是N    pts[0]是[231.68674 271.08926   1.     ]  齐次坐标都有了
  364.     # print(pts[0])
  365.     desc = data_1['descriptors']  # np类型  (3473, 512)
  366.     print("图1提取特征点数目", pts.shape[0])
  367.  
  368.  
  369.     #imgx = img1.copy()  #不要copy就不是灰色的
  370.     pts=pts.T                               # 因为其他函数写的时候都是对应3*N的来写的，转置以后成为3*N
  371.     desc=desc.T                     #转职后成为512*N，本身是(3473, 512)
  372.     img11 = showpoint(imgx_1, pts)
  373.     cv2.namedWindow('imgone', 0)
  374.     cv2.imshow("imgone", img11)  #第一张图显示特征点
  375.     #cv2.waitKey(0)   有这个才能显示窗口
  376.     ##第2图的特征点和描述符
  377.  
  378.     pts1 = data_2['keypoints']  # np类型  (3473, 3)  3473是N    pts[0]是[231.68674 271.08926   1.     ]  齐次坐标都有了
  379.     desc1 = data_2['descriptors']  # np类型  (3473, 512)
  380.     desc1=desc1.T
  381.     print("图1提取特征点数目", pts1.shape[0])
  382.     pts1=pts1.T
  383.     img22 = showpoint(imgx_2, pts1)
  384.     cv2.namedWindow('imgtwo', 0)
  385.     cv2.imshow("imgtwo", img22)  #第一张图显示特征点
  386.     #cv2.waitKey(0)   #有这个才能显示窗口
  387.  
  388.  
  389.  
  390.  
  391.     #下面都是共同使用的东西
  392.     # 尝试暴力匹配+KNN算法实现
  393.     desc_test_0=desc.T    #矩阵转置后 #可以顺利输出 这里必须(1099, 128)N在前
  394.     desc_test_1 = desc1.T
  395.     pts_test_0 =pts.T       #特征点也转置
  396.     pts_test_1 =pts1.T
  397.     #特征点转为元组
  398.     pts_test_0 = tuple(pts_test_0)
  399.     pts_test_1 = tuple(pts_test_1)
  400.  
  401.     match_test = cv2.BFMatcher()
  402.     matches = match_test.knnMatch(desc_test_0, desc_test_1, k=2)
  403.     #print("这里的长度是",len(matches))  #tuple  长度893
  404.  
  405.     good = []
  406.     for m, n in matches:
  407.         if m.distance < 0.8 * n.distance:   #可以自己改
  408.             good.append([m])
  409.     print("number of feature points:",pts.shape[1], pts1.shape[1])  #特征点个数，调试正确
  410.     print("good match num:{} good match points:".format(len(good))) #正确  knn匹配到25个
  411.  
  412.  
  413.  
  414.     #
  415.     #这里只能将pts, pts1抓换后才能用
  416.     Max_num, good_F, inlier_points = ransac(good, pts_test_0, pts_test_1, confidence=30, iter_num=500)
  417.     print(pts.shape[1], "/", pts1.shape[1], "/", len(good), "/", Max_num, "/", '%.2f' % (Max_num / len(good)))
  418.     #print(inlier_points)   #这个就是得出的ransac后的匹配数
  419.  
  420.     #inlier_points 转化为match的坐标  #match是一个3*67的数组，这里转换后还要转置下T
  421.     #print("xxxxxxxxx",len(inlier_points))
  422.     #print(inlier_points)
  423.     tuple_x = []
  424.     tuple_y = []
  425.     #matche = np.zeros((2,len(inlier_points)))   #初始化matche， 0填充 M行 2列，M为匹配到的点数量  最后转置下T
  426.     #print(type(matche))  #np类型的
  427.     for i in inlier_points:
  428.         #print(i)
  429.         idex_x=i[0].queryIdx                  #查询的描述符索引 845个
  430.         idex_y = i[0].trainIdx                #目标的描述符索引 对应555个
  431.         #print("xxxxxxx",idex_0)
  432.         tuple_x=np.append(tuple_x,idex_x)       #  查询的描述符的索引弄到一起
  433.         #print(tuple_x)
  434.         tuple_y = np.append(tuple_y, idex_y )   #  目标的描述符的索引弄到一起
  435.  
  436.  
  437.     match=np.vstack((tuple_x,tuple_y))
  438.  
  439.  
  440.  
  441.  
  442.  
  443.  
  444.  
  445.     #画图  test
  446.     out =drawMatches(imgx_1, pts, imgx_2,pts1, match)   #传入了彩色图
  447.     cv2.namedWindow("matcher", 0)
  448.     cv2.imshow("matcher", out)
  449.     cv2.waitKey(0) # 函数的功能是不断刷新图像 , 频率时间为delay , 单位为ms ,处理视频的时候用的
  450.     print('==> Finshed ')
  451.  
  452.     #
  453.  
  454.  
  455.

最终显示图片    D2-net实现

 这个是Superpoint实现的

这里显示的是 ：特征点数量 A，B/匹配对数/正确对数/准确率

声明：本文内容由网友自发贡献，不代表【wpsshop博客】立场，版权归原作者所有，本站不承担相应法律责任。如您发现有侵权的内容，请联系我们。转载请注明出处：https://www.wpsshop.cn/w/你好赵伟/article/detail/936508
推荐阅读
相关标签

Copyright © 2003-2013 www.wpsshop.cn 版权所有，并保留所有权利。

  

闽ICP备14008679号