Python映像特徵檢測演算法（1）：Python實現SIFT和Harris

Python映像特徵檢測演算法（1）：Python實現SIFT和Harris

本文將介紹用於映像匹配的兩種局部描述子演算法，SIFT[論文連結]和Harris，它們在很多應用中都有比較重要的作用，比如目標匹配、目標跟蹤、建立全景、增強現實技術以及計算映像的三維重建等，而常用的特徵有顏色、角點、特徵點、輪廓、紋理等，很多內容中都會用到這些特徵。Harris角點檢測演算法（也稱Harris&Stephens角點檢測器）是特徵點檢測的基礎，是一個極為簡單的角點檢測演算法。Harris提出了應用鄰近像素點灰階差值的概念，從而進行判斷是否為角點（興趣點）、邊緣、平滑地區。Harris角點檢測原理是利用移動的視窗在映像中計算灰階變化值，其中關鍵流程包括轉化為灰階映像、計算差分映像、高斯平滑、計算局部極值、確認角點（關鍵點）。SIFT[首頁]，即尺度不變特徵變換（Scale-invariantfeature transform，SIFT），是用於影像處理領域的一種描述。這種描述具有尺度不變性，可在映像中檢測出關鍵點，是一種局部特徵描述子。該方法於1999年由David Lowe[首頁]首先發表於電腦視覺國際會議（InternationalConference on Computer Vision，ICCV），2004年再次經David Lowe整理完善後發表於International journal of computer vision(IJCV)。截止2014年8月，該論文單篇被引次數達25000餘次。

一、前期準備

1.安裝配置Python及其相關軟體包（推薦使用Anaconda，其內建完整的庫，無需費時費力安裝）

2.安裝配置VLFeat開源工具包（windows10/64位系統為例）

我們使開源工具包VLFeat提供的二進位檔案來計算映像的SIFT特徵，VLFeat工具包可以從http://www.vlfeat.org/下載，二進位檔案可以在所有的主要平台上運行。VLFeat的首頁如下圖所示，可以看出VLFeat工具包中含有大量的完整的代碼實現。VLFeat庫是用C語言來寫的，但是我們可以使用該庫提供的命令列介面。

在此介面上下載檔案vlfeat-0.9.20-bin，解壓縮後，將vlfeat-0.9.20/bin/win32（win64下的測試不相容）檔案夾下的sift.exe和vl.dll拷貝到當前工程目錄下。這樣，前期工作就已完成。

二、編寫代碼

Python實現SIFT和Harris的具體代碼如下：

# -*- coding: utf-8 -*-# Yan Zhenguofrom PIL import Imagefrom pylab import *from numpy import *import osfrom scipy.ndimage import filtersdef process_image(imagename, resultname, params="--edge-thresh 10 --peak-thresh 5"):    """ Process an image and save the results in a file. """    if imagename[-3:] != 'pgm':        # create a pgm file        im = Image.open(imagename).convert('L')        im.save('tmp.pgm')        imagename = 'tmp.pgm'    cmmd = str("sift " + imagename + " --output=" + resultname + " " + params)    os.system(cmmd)    print('processed', imagename, 'to', resultname)def read_features_from_file(filename):    """ Read feature properties and return in matrix form. """    f = loadtxt(filename)    return f[:, :4], f[:, 4:]  # feature locations, descriptorsdef write_features_to_file(filename, locs, desc):    """ Save feature location and descriptor to file. """    savetxt(filename, hstack((locs, desc)))def plot_features(im, locs, circle=False):    """ Show image with features. input: im (image as array),        locs (row, col, scale, orientation of each feature). """    def draw_circle(c, r):        t = arange(0, 1.01, .01) * 2 * pi        x = r * cos(t) + c[0]        y = r * sin(t) + c[1]        plot(x, y, 'b', linewidth=2)    imshow(im)    if circle:        for p in locs:            draw_circle(p[:2], p[2])    else:        plot(locs[:, 0], locs[:, 1], 'ob')    axis('off')def compute_harris_response(im, sigma=3):    """ Compute the Harris corner detector response function        for each pixel in a graylevel image. """    # derivatives    imx = zeros(im.shape)    filters.gaussian_filter(im, (sigma, sigma), (0, 1), imx)    imy = zeros(im.shape)    filters.gaussian_filter(im, (sigma, sigma), (1, 0), imy)    # compute components of the Harris matrix    Wxx = filters.gaussian_filter(imx * imx, sigma)    Wxy = filters.gaussian_filter(imx * imy, sigma)    Wyy = filters.gaussian_filter(imy * imy, sigma)    # determinant and trace    Wdet = Wxx * Wyy - Wxy ** 2    Wtr = Wxx + Wyy    return Wdet / Wtrdef get_harris_points(harrisim, min_dist=10, threshold=0.1):    """ Return corners from a Harris response image        min_dist is the minimum number of pixels separating        corners and image boundary. """    # find top corner candidates above a threshold    corner_threshold = harrisim.max() * threshold    harrisim_t = (harrisim > corner_threshold) * 1    # get coordinates of candidates    coords = array(harrisim_t.nonzero()).T    # ...and their values    candidate_values = [harrisim[c[0], c[1]] for c in coords]    # sort candidates (reverse to get descending order)    index = argsort(candidate_values)[::-1]    # store allowed point locations in array    allowed_locations = zeros(harrisim.shape)    allowed_locations[min_dist:-min_dist, min_dist:-min_dist] = 1    # select the best points taking min_distance into account    filtered_coords = []    for i in index:        if allowed_locations[coords[i, 0], coords[i, 1]] == 1:            filtered_coords.append(coords[i])            allowed_locations[(coords[i, 0] - min_dist):(coords[i, 0] + min_dist),            (coords[i, 1] - min_dist):(coords[i, 1] + min_dist)] = 0    return filtered_coordsdef plot_harris_points(image, filtered_coords):    """ Plots corners found in image. """    figure()    gray()    imshow(image)    plot([p[1] for p in filtered_coords],[p[0] for p in filtered_coords], '*')    axis('off')    title('Harris-Features')    show()imname = 'building.jpg'im = array(Image.open(imname).convert('L'))process_image(imname, 'building.sift')l1, d1 = read_features_from_file('building.sift')figure()gray()"""Figure1:SIFT特徵"""plot_features(im, l1, circle=False)title('SIFT-Features')"""Figure2:使用圓圈表示特徵尺度的SIFT特徵"""figure()gray()plot_features(im, l1, circle=True)title('Detect-SIFT-Features')"""Figure3:Harris角點檢測的結果"""harrisim = compute_harris_response(im)filtered_coords = get_harris_points(harrisim, 6, 0.05)plot_harris_points(im, filtered_coords)

三、運行結果

1、SIFT特徵

2、使用圓圈表示特徵尺度的SIFT特徵

3、Harris角點

在後續工作中，我將繼續為大家展現映像特徵檢測相關得到工作和深度學習網路帶來的無盡樂趣，我將和大家一起探討映像世界和深度學習的奧秘。當然，如果你感興趣，我的Weibo將與你一起分享最前沿的人工智慧、機器學習、深度學習與電腦視覺方面的技術。

四、參考文獻

1.Python電腦視覺編程

2.百度百科SIFT：https://baike.baidu.com/item/SIFT/1396275?fr=aladdin