本文側重于如何使用Python語言實作SIFT演算法

所有程式已打包：基于OpenCV-Python的SIFT演算法的實作

一、什么是SIFT演算法

??SIFT，即尺度不變特征變換（Scale-invariant feature transform，SIFT），是用于影像處理領域的一種描述，這種描述具有尺度不變性，可在影像中檢測出關鍵點，是一種區域特征描述子，

二、準備作業

2.1 實驗設備

??本文在Windows10系統上，使用pycharm軟體完成所有實驗，

2.2 OpenCV安裝

??我們可以使用OpenCV庫中的cv2.xfeatures2d.SIFT_create()函式實作SIFT，但由于專利保護，很多版本的OpenCV庫已無法提供該函式，目前僅3.4.2.16版本的OpenCV庫可使用此函式，

安裝教程：
??（1）查看當前版本opencv：進入cmd（組合鍵win+R，輸入cmd），輸入conda list，查看當前pycharm所有庫并找到opencv-python，若找不到庫，說明沒有安裝，
??（2）卸載原版本（在cmd中輸入：pip uninstall opencv）
??（3）安裝新版本（在cmd中輸入：pip install opencv-python==3.4.2.16 -i "https://pypi.doubanio.com/simple/"
??（4）安裝附屬庫（在cmd中輸入：pip install opencv-contrib-python==3.4.2.16 -i "https://pypi.doubanio.com/simple/"）

三、實驗作業

3.1 影像選擇

??這里選擇經典的lena影像作為實驗物件，為了選擇一個待匹配影像，本文使用如下代碼對lena影像進行逆時針45°旋轉，

from PIL import Image

img = Image.open('lena.png')
img2 = img.rotate(45)       # 逆時針旋轉45°
img2.save("lena_rot45.png")
img2.show()

參考影像與待匹配影像（即旋轉影像）如下圖所示：

3.2 程式實作

"""
影像匹配——SIFT點特征匹配實作步驟：
    （1）讀取影像；
    （2）定義sift算子；
    （3）通過sift算子對需要匹配的影像進行特征點獲取；
        a.可獲取各匹配影像經過sift算子的特征點數目
    （4）可視化特征點（在原圖中標記為圓圈）；
        a.為方便觀察，可將匹配影像橫向拼接
    （5）影像匹配（特征點匹配）；
        a.通過調整ratio獲取需要進行影像匹配的特征點數量（ratio值越大，匹配的線條越密集，但錯誤匹配點也會增多）
        b.通過索引ratio選擇固定的特征點進行影像匹配
    （6）將待匹配影像通過旋轉、變換等方式將其與目標影像對齊
"""

import cv2              # opencv版本需為3.4.2.16
import numpy as np      # 矩陣運算庫
import time             # 時間庫

original_lena = cv2.imread('lena.png')          # 讀取lena原圖
lena_rot45 = cv2.imread('lena_rot45.png')       # 讀取lena旋轉45°圖

sift = cv2.xfeatures2d.SIFT_create()

# 獲取各個影像的特征點及sift特征向量
# 回傳值kp包含sift特征的方向、位置、大小等資訊；des的shape為（sift_num， 128）， sift_num表示影像檢測到的sift特征數量
(kp1, des1) = sift.detectAndCompute(original_lena, None)
(kp2, des2) = sift.detectAndCompute(lena_rot45, None)

# 特征點數目顯示
print("=========================================")
print("=========================================")
print('lena 原圖  特征點數目：', des1.shape[0])
print('lena 旋轉圖 特征點數目：', des2.shape[0])
print("=========================================")
print("=========================================")

# 舉例說明kp中的引數資訊
for i in range(2):
    print("關鍵點", i)
    print("資料型別:", type(kp1[i]))
    print("關鍵點坐標:", kp1[i].pt)
    print("鄰域直徑:", kp1[i].size)
    print("方向:", kp1[i].angle)
    print("所在的影像金字塔的組:", kp1[i].octave)

print("=========================================")
print("=========================================")
"""
首先對原圖和旋轉圖進行特征匹配，即圖original_lena和圖lena_rot45
"""
# 繪制特征點，并顯示為紅色圓圈
sift_original_lena = cv2.drawKeypoints(original_lena, kp1, original_lena, color=(255, 0, 255))
sift_lena_rot45 = cv2.drawKeypoints(lena_rot45, kp2, lena_rot45, color=(255, 0, 255))

sift_cat1 = np.hstack((sift_original_lena, sift_lena_rot45))        # 對提取特征點后的影像進行橫向拼接
cv2.imwrite("sift_cat1.png", sift_cat1)
print('原圖與旋轉圖 特征點繪制影像已保存')
cv2.imshow("sift_point1", sift_cat1)
cv2.waitKey()

# 特征點匹配
# K近鄰演算法求取在空間中距離最近的K個資料點，并將這些資料點歸為一類
start = time.time()     # 計算匹配點匹配時間
bf = cv2.BFMatcher()
matches1 = bf.knnMatch(des1, des2, k=2)
print('用于 原圖和旋轉圖 影像匹配的所有特征點數目：', len(matches1))

# 調整ratio
# ratio=0.4：對于準確度要求高的匹配；
# ratio=0.6：對于匹配點數目要求比較多的匹配；
# ratio=0.5：一般情況下，
ratio1 = 0.5
good1 = []

for m1, n1 in matches1:
    # 如果最接近和次接近的比值大于一個既定的值，那么我們保留這個最接近的值，認為它和其匹配的點為good_match
    if m1.distance < ratio1 * n1.distance:
        good1.append([m1])

end = time.time()
print("匹配點匹配運行時間:%.4f秒" % (end-start))

# 通過對good值進行索引，可以指定固定數目的特征點進行匹配，如good[:20]表示對前20個特征點進行匹配
match_result1 = cv2.drawMatchesKnn(original_lena, kp1, lena_rot45, kp2, good1, None, flags=2)
cv2.imwrite("match_result1.png", match_result1)

print('原圖與旋轉圖 特征點匹配影像已保存')
print("=========================================")
print("=========================================")
print("原圖與旋轉圖匹配對的數目:", len(good1))

for i in range(2):
    print("匹配", i)
    print("資料型別:", type(good1[i][0]))
    print("描述符之間的距離:", good1[i][0].distance)
    print("查詢影像中描述符的索引:", good1[i][0].queryIdx)
    print("目標影像中描述符的索引:", good1[i][0].trainIdx)

print("=========================================")
print("=========================================")
cv2.imshow("original_lena and lena_rot45 feature matching result", match_result1)
cv2.waitKey()

# 將待匹配影像通過旋轉、變換等方式將其與目標影像對齊，這里使用單應性矩陣，
# 單應性矩陣有八個引數，如果要解這八個引數的話，需要八個方程，由于每一個對應的像素點可以產生2個方程(x一個，y一個)，那么總共只需要四個像素點就能解出這個單應性矩陣，
if len(good1) > 4:
    ptsA = np.float32([kp1[m[0].queryIdx].pt for m in good1]).reshape(-1, 1, 2)
    ptsB = np.float32([kp2[m[0].trainIdx].pt for m in good1]).reshape(-1, 1, 2)
    ransacReprojThreshold = 4
    # RANSAC演算法選擇其中最優的四個點
    H, status =cv2.findHomography(ptsA, ptsB, cv2.RANSAC, ransacReprojThreshold)
    imgout = cv2.warpPerspective(lena_rot45, H, (original_lena.shape[1], original_lena.shape[0]),
                                 flags=cv2.INTER_LINEAR + cv2.WARP_INVERSE_MAP)

    cv2.imwrite("imgout.png", imgout)
    cv2.imshow("lena_rot45's result after transformation", imgout)
    cv2.waitKey()

3.3 程式結果