void genPyr(const Mat& imgSrc, vector<Mat>& outPutArray, int TYPE, int level)
{outPutArray.assign(level + 1, Mat());outPutArray[0] = imgSrc.clone(); // the 0 level is the image. for (int i = 0; i != level; i++){pyrDown(outPutArray[i], outPutArray[i + 1]);}if (PYR_GUASS == TYPE){return;}for (int i = 0; i != level; i++){Mat UpSampleImg;pyrUp(outPutArray[i + 1], UpSampleImg, outPutArray[i].size());outPutArray[i] -= UpSampleImg;}
}

3 為什么用高斯核

圖像的金字塔化能高效地（計算效率也較高）對圖像進行多尺度的表達，但它缺乏堅實的理論基礎，不能分析圖像中物體的各種尺度（雖然我們有小貓的金字塔圖像，我們還是不知道原圖像內小貓的大小）。

信號的尺度空間剛提出是就是通過一系列單參數、寬度遞增的高斯濾波器將原始信號濾波得到到組低頻信號。那么一個很明顯的疑問是：除了高斯濾波之外，其他帶有參數t的低通濾波器是否也可以用來生成一個尺度空間。

后來Koenerink、Lindeberg[Scale-space theory in computer vision]、Florack等人用精確的數學形式通過不同的途徑都證明了高斯核是實現尺度變換的唯一變換核。

雖然很多研究者從可分性、旋轉不變性、因果性等特性推出高斯濾波器是建立線性尺度空間的最優濾波器。然后在數字圖像處理中，需要對核函數進行采樣，離散的高斯函數并不滿足連續高斯函數的的一些優良的性質。所以后來出現了一些非線性的濾波器組來建立尺度空間，如B樣條核函數。

使用高斯濾波器對圖像進行尺度空間金塔塔圖的構建，讓這個尺度空間具有下面的性質：

1）加權平均和有限孔徑效應

信號在尺度t上的表達可以看成是原信號在空間上的一系列加權平均，權重就是具有不同尺度參數的高斯核。

信號在尺度t上的表達也對應于用一個無方向性的孔徑函數（特征長度為σ=t√）來觀測信號的結果。這時候信號中特征長度小于σ的精細結構會被抑制[理解為一維信號上小于σ的波動會被平滑掉。]。

2）層疊平滑

也叫高斯核族的半群（Semi-Group）性質：兩個高斯核的卷積等同于另外一個不同核參數的高斯核卷積。

這個性質的意思就是說不同的高斯核對圖像的平滑是連續的。

3）局部極值遞性

這個特征可以從人眼的視覺原理去理解，人在看一件物體時，離得越遠，物體的細節看到的越少，細節特征是在減少的。

高斯核對圖像進行濾波具有壓制局部細節的性質。

4）尺度伸縮不變性。

這里只是一個公式推導的問題，對原來的信號加一個變換函數，對變換后的信號再進行高斯核的尺度空間生成，新的信號的極值點等特征是不變的。

Young對經生理學的研究中發現，哺乳動物的視網膜和視覺皮層的感受區域可以很好地用4階以內的高斯微分來建模。

4 尺度的選擇

一般我們采集到的圖像中，我們并不知道我們感興趣的目標在圖像中的尺度，在這樣的情況下，我們對圖像進行分析時就無法選擇合適的參數，比如邊緣檢測，可能由于參數不當，而造成過多的局部細節。

如下圖所示：紅色圓圈內的斑點的大小（直徑）比例對應著兩幅圖像之間尺度比例（scale ratio）。如果對兩幅圖像采用相同的固定尺度的LoG檢測器檢測，很難將這兩個斑點檢測出來。LoG檢測器相當于一個匹配濾波器，只有當LoG的尺度與圖片中斑點結構尺度相當時才會有較強的響應。如果用與左圖中斑點結構相當大小尺度LoG算子，在中的大斑點的對應的LoG響應很小不能被檢測出來，反之亦然。因此固定尺度的LoG斑點檢測器不具有尺度不變性。使用尺度空間進行多尺度檢測可以將兩幅圖像中不同尺度的斑點檢測出來。但是由于斑點結構是在一定尺度范圍之內存在的，比如用5~8尺度的LoG可能都能檢測出來右邊圖像中的斑點結構，所以在尺度空間中進行斑點檢測會有重復檢測的缺點。

?
在實際操作中，我們需要定義一個特征響應函數，在不同的尺度空間上尋找一個極值點。比如小貓的金字塔圖像分析時，我們定義了一個大小為[w,h]的小貓的模板，用這個模板去與金字塔系列圖像匹配，一定有匹配度最佳（即特征響應最強）。

需要注意的是，圖像結構往往是在粗糙的尺度上被檢測到，此時位置信息未必是最準確的，因此通常圖像的尺度分析包含兩個階段：首先在粗尺度上進行特征（結構）檢測，然后再在細尺度上進行精確定位。

5.高斯金字塔

Gaussian Pyramid

• Imagine the pyramid as a set of layers in which the higher the layer, the smaller the size.

• Every layer is numbered from bottom to top, so layer??(denoted as??is smaller than layer??().

• To produce layer??in the Gaussian pyramid, we do the following:

  • Convolve??with a Gaussian kernel:

  • Remove every even-numbered row and column.

• You can easily notice that the resulting image will be exactly one-quarter the area of its predecessor. Iterating this process on the input image??(original image) produces the entire pyramid.

• The procedure above was useful to downsample an image. What if we want to make it bigger?:

  • First, upsize the image to twice the original in each dimension, wit the new even rows and columns filled with zeros ()
  • Perform a convolution with the same kernel shown above (multiplied by 4) to approximate the values of the “missing pixels”

• These two procedures (downsampling and upsampling as explained above) are implemented by the OpenCV functions?pyrUp?and?pyrDown, as we will see in an example with the code below:

Note

When we reduce the size of an image, we are actually?losing?information of the image.

Code

This tutorial code’s is shown lines below. You can also download it from?here

#include "opencv2/imgproc/imgproc.hpp"
#include "opencv2/highgui/highgui.hpp"
#include <math.h>
#include <stdlib.h>
#include <stdio.h>using namespace cv;/// Global variables
Mat src, dst, tmp;
char* window_name = "Pyramids Demo";/**
 * @function main
 */
int main( int argc, char** argv )
{/// General instructionsprintf( "\n Zoom In-Out demo  \n " );printf( "------------------ \n" );printf( " * [u] -> Zoom in  \n" );printf( " * [d] -> Zoom out \n" );printf( " * [ESC] -> Close program \n \n" );/// Test image - Make sure it s divisible by 2^{n}src = imread( "../images/chicky_512.jpg" );if( !src.data ){ printf(" No data! -- Exiting the program \n");return -1; }tmp = src;dst = tmp;/// Create windownamedWindow( window_name, CV_WINDOW_AUTOSIZE );imshow( window_name, dst );/// Loopwhile( true ){int c;c = waitKey(10);if( (char)c == 27 ){ break; }if( (char)c == 'u' ){ pyrUp( tmp, dst, Size( tmp.cols*2, tmp.rows*2 ) );printf( "** Zoom In: Image x 2 \n" );}else if( (char)c == 'd' ){ pyrDown( tmp, dst, Size( tmp.cols/2, tmp.rows/2 ) );printf( "** Zoom Out: Image / 2 \n" );}imshow( window_name, dst );tmp = dst;}return 0;
}

Explanation

1. Let’s check the general structure of the program:

   • Load an image (in this case it is defined in the program, the user does not have to enter it as an argument)

     /// Test image - Make sure it s divisible by 2^{n}
     src = imread( "../images/chicky_512.jpg" );
     if( !src.data ){ printf(" No data! -- Exiting the program \n");return -1; }
     

   • Create a Mat object to store the result of the operations (dst) and one to save temporal results (tmp).

     Mat src, dst, tmp;
     /* ... */
     tmp = src;
     dst = tmp;
     

   • Create a window to display the result

     namedWindow( window_name, CV_WINDOW_AUTOSIZE );
     imshow( window_name, dst );
     

   • Perform an infinite loop waiting for user input.

     while( true )
     {int c;c = waitKey(10);if( (char)c == 27 ){ break; }if( (char)c == 'u' ){ pyrUp( tmp, dst, Size( tmp.cols*2, tmp.rows*2 ) );printf( "** Zoom In: Image x 2 \n" );}else if( (char)c == 'd' ){ pyrDown( tmp, dst, Size( tmp.cols/2, tmp.rows/2 ) );printf( "** Zoom Out: Image / 2 \n" );}imshow( window_name, dst );tmp = dst;
     }
     

     Our program exits if the user presses?ESC. Besides, it has two options:

     • Perform upsampling (after pressing ‘u’)

       pyrUp( tmp, dst, Size( tmp.cols*2, tmp.rows*2 )
       

       We use the function?pyrUp?with 03 arguments:

       • tmp: The current image, it is initialized with the?src?original image.
       • dst: The destination image (to be shown on screen, supposedly the double of the input image)
       • Size( tmp.cols*2, tmp.rows*2 )?: The destination size. Since we are upsampling,?pyrUp?expects a size double than the input image (in this case?tmp).

     • Perform downsampling (after pressing ‘d’)

       pyrDown( tmp, dst, Size( tmp.cols/2, tmp.rows/2 )
       

       Similarly as with?pyrUp, we use the function?pyrDown?with 03 arguments:

       • tmp: The current image, it is initialized with the?src?original image.
       • dst: The destination image (to be shown on screen, supposedly half the input image)
       • Size( tmp.cols/2, tmp.rows/2 )?: The destination size. Since we are upsampling,?pyrDown?expects half the size the input image (in this case?tmp).

     • Notice that it is important that the input image can be divided by a factor of two (in both dimensions). Otherwise, an error will be shown.

     • Finally, we update the input image?tmp?with the current image displayed, so the subsequent operations are performed on it.

       tmp = dst;
       

Results

• After compiling the code above we can test it. The program calls an image?chicky_512.jpg?that comes in the?tutorial_code/image?folder. Notice that this image is?, hence a downsample won’t generate any error (). The original image is shown below:

• First we apply two successive?pyrDown?operations by pressing ‘d’. Our output is:

• Note that we should have lost some resolution due to the fact that we are diminishing the size of the image. This is evident after we apply?pyrUp?twice (by pressing ‘u’). Our output is now:

六：補充說明【OpenCV】SIFT原理與源碼分析：DoG尺度空間構造

《SIFT原理與源碼分析》系列文章索引：http://blog.csdn.net/xiaowei_cqu/article/details/8069548

尺度空間理論

為什么要討論尺度空間？

圖像的尺度空間表達就是圖像在所有尺度下的描述。

尺度空間表達與金字塔多分辨率表達

高斯模糊

高斯核是唯一可以產生多尺度空間的核（《Scale-space theory: A basic tool for analysing structures at different scales》）。一個圖像的尺度空間L(x,y,σ) ,定義為原始圖像I(x,y)與一個可變尺度的2維高斯函數G(x,y,σ)卷積運算。?

二維空間高斯函數：

尺度空間：

尺度是自然客觀存在的，不是主觀創造的。高斯卷積只是表現尺度空間的一種形式。

二維空間高斯函數是等高線從中心成正太分布的同心圓：

分布不為零的點組成卷積陣與原始圖像做變換，即每個像素值是周圍相鄰像素值的高斯平均。一個5*5的高斯模版如下所示：

高斯模糊另一個很厲害的性質就是線性可分：使用二維矩陣變換的高斯模糊可以通過在水平和豎直方向各進行一維高斯矩陣變換相加得到。

O(N^2*m*n)次乘法就縮減成了O(N*m*n)+O(N*m*n)次乘法。（N為高斯核大小，m,n為二維圖像高和寬）

其實高斯這一部分只需要簡單了解就可以了，在OpenCV也只需要一句代碼：

1. GaussianBlur(dbl,?dbl,?Size(),?sig_diff,?sig_diff);??

金字塔多分辨率

金字塔是早期圖像多尺度的表示形式。圖像金字塔化一般包括兩個步驟：使用低通濾波器平滑圖像；對平滑圖像進行降采樣（通常是水平，豎直方向1/2），從而得到一系列尺寸縮小的圖像。

上圖中（a）是對原始信號進行低通濾波，（b）是降采樣得到的信號。

而對于二維圖像，一個傳統的金字塔中，每一層圖像由上一層分辨率的長、寬各一半，也就是四分之一的像素組成：

多尺度和多分辨率

尺度空間表達和金字塔多分辨率表達之間最大的不同是：

• 尺度空間表達是由不同高斯核平滑卷積得到，在所有尺度上有相同的分辨率；
• 而金字塔多分辨率表達每層分辨率減少固定比率。

DoG（Difference of Gaussian）

高斯拉普拉斯LoG金字塔

高斯差分DoG金字塔

金字塔構建

構建高斯金字塔

1. //?構建nOctaves組（每組nOctaves+3層）高斯金字塔??
2. void?SIFT::buildGaussianPyramid(?const?Mat&?base,?vector<Mat>&?pyr,?int?nOctaves?)?const??
3. {??
4. ????vector<double>?sig(nOctaveLayers?+?3);??
5. ????pyr.resize(nOctaves*(nOctaveLayers?+?3));??
6. ??
7. ????//?precompute?Gaussian?sigmas?using?the?following?formula:??
8. ????//??\sigma_{total}^2?=?\sigma_{i}^2?+?\sigma_{i-1}^2、??
9. ????//?計算對圖像做不同尺度高斯模糊的尺度因子??
10. ????sig[0]?=?sigma;??
11. ????double?k?=?pow(?2.,?1.?/?nOctaveLayers?);??
12. ????for(?int?i?=?1;?i?<?nOctaveLayers?+?3;?i++?)??
13. ????{??
14. ????????double?sig_prev?=?pow(k,?(double)(i-1))*sigma;??
15. ????????double?sig_total?=?sig_prev*k;??
16. ????????sig[i]?=?std::sqrt(sig_total*sig_total?-?sig_prev*sig_prev);??
17. ????}??
18. ??
19. ????for(?int?o?=?0;?o?<?nOctaves;?o++?)??
20. ????{??
21. ????????//?DoG金子塔需要nOctaveLayers+2層圖像來檢測nOctaves層尺度??
22. ????????//?所以高斯金字塔需要nOctaveLayers+3層圖像得到nOctaveLayers+2層DoG金字塔??
23. ????????for(?int?i?=?0;?i?<?nOctaveLayers?+?3;?i++?)??
24. ????????{??
25. ????????????//?dst為第o組（Octave）金字塔??
26. ????????????Mat&?dst?=?pyr[o*(nOctaveLayers?+?3)?+?i];??
27. ????????????//?第0組第0層為原始圖像??
28. ????????????if(?o?==?0??&&??i?==?0?)??
29. ????????????????dst?=?base;??
30. ??????????????
31. ????????????//?base?of?new?octave?is?halved?image?from?end?of?previous?octave??
32. ????????????//?每一組第0副圖像時上一組倒數第三幅圖像隔點采樣得到??
33. ????????????else?if(?i?==?0?)??
34. ????????????{??
35. ????????????????const?Mat&?src?=?pyr[(o-1)*(nOctaveLayers?+?3)?+?nOctaveLayers];??
36. ????????????????resize(src,?dst,?Size(src.cols/2,?src.rows/2),??
37. ???????????????????????0,?0,?INTER_NEAREST);??
38. ????????????}??
39. ????????????//?每一組第i副圖像是由第i-1副圖像進行sig[i]的高斯模糊得到??
40. ????????????//?也就是本組圖像在sig[i]的尺度空間下的圖像??
41. ????????????else??
42. ????????????{??
43. ????????????????const?Mat&?src?=?pyr[o*(nOctaveLayers?+?3)?+?i-1];??
44. ????????????????GaussianBlur(src,?dst,?Size(),?sig[i],?sig[i]);??
45. ????????????}??
46. ????????}??
47. ????}??
48. }??

構建DoG金字塔

1. //?構建nOctaves組（每組nOctaves+2層）高斯差分金字塔??
2. void?SIFT::buildDoGPyramid(?const?vector<Mat>&?gpyr,?vector<Mat>&?dogpyr?)?const??
3. {??
4. ????int?nOctaves?=?(int)gpyr.size()/(nOctaveLayers?+?3);??
5. ????dogpyr.resize(?nOctaves*(nOctaveLayers?+?2)?);??
6. ??
7. ????for(?int?o?=?0;?o?<?nOctaves;?o++?)??
8. ????{??
9. ????????for(?int?i?=?0;?i?<?nOctaveLayers?+?2;?i++?)??
10. ????????{??
11. ????????????//?第o組第i副圖像為高斯金字塔中第o組第i+1和i組圖像相減得到??
12. ????????????const?Mat&?src1?=?gpyr[o*(nOctaveLayers?+?3)?+?i];??
13. ????????????const?Mat&?src2?=?gpyr[o*(nOctaveLayers?+?3)?+?i?+?1];??
14. ????????????Mat&?dst?=?dogpyr[o*(nOctaveLayers?+?2)?+?i];??
15. ????????????subtract(src2,?src1,?dst,?noArray(),?CV_16S);??
16. ????????}??
17. ????}??
18. } ?

SIFT簡介

Scale Invariant Feature Transform，尺度不變特征變換匹配算法，是由David G. Lowe在1999年（《Object Recognition from Local Scale-Invariant Features》）提出的高效區域檢測算法，在2004年（《Distinctive Image Features from Scale-Invariant Keypoints》）得以完善。

SIFT特征對旋轉、尺度縮放、亮度變化等保持不變性，是非常穩定的局部特征，現在應用很廣泛。而SIFT算法是將Blob檢測，特征矢量生成，特征匹配搜索等步驟結合在一起優化。我會更新一系列文章，分析SIFT算法原理及OpenCV 2.4.2實現的SIFT源碼：

1. DoG尺度空間構造（Scale-space extrema detection）http://blog.csdn.net/xiaowei_cqu/article/details/8067881
2. 關鍵點搜索與定位（Keypoint localization）http://blog.csdn.net/xiaowei_cqu/article/details/8087239
3. 方向賦值（Orientation assignment）http://blog.csdn.net/xiaowei_cqu/article/details/8096072
4. 關鍵點描述（Keypoint descriptor）http://blog.csdn.net/xiaowei_cqu/article/details/8113565
5. OpenCV實現：特征檢測器FeatureDetectorhttp://blog.csdn.net/xiaowei_cqu/article/details/8652096
   
6. SIFT中LoG和DoG的比較http://blog.csdn.net/xiaowei_cqu/article/details/27692123

SIFT in OpenCV

構造函數：

1. SIFT::SIFT(int?nfeatures=0,?int?nOctaveLayers=3,?double?contrastThreshold=0.04,?double?edgeThreshold=??
2. 10,?double?sigma=1.6)??

重載操作符：

1. void?SIFT::operator()(InputArray?img,?InputArray?mask,?vector<KeyPoint>&?keypoints,?OutputArray??
2. descriptors,?bool?useProvidedKeypoints=false)??

函數源碼

1. SIFT::SIFT(?int?_nfeatures,?int?_nOctaveLayers,??
2. ???????????double?_contrastThreshold,?double?_edgeThreshold,?double?_sigma?)??
3. ????:?nfeatures(_nfeatures),?nOctaveLayers(_nOctaveLayers),??
4. ????contrastThreshold(_contrastThreshold),?edgeThreshold(_edgeThreshold),?sigma(_sigma)??
5. ????//?sigma：對第0層進行高斯模糊的尺度空間因子。??
6. ????//?默認為1.6（如果是軟鏡攝像頭捕獲的圖像，可以適當減小此值）??
7. {??
8. }??

1. void?SIFT::operator()(InputArray?_image,?InputArray?_mask,??
2. ??????????????????????vector<KeyPoint>&?keypoints,??
3. ??????????????????????OutputArray?_descriptors,??
4. ??????????????????????bool?useProvidedKeypoints)?const??
5. //?mask?：Optional?input?mask?that?marks?the?regions?where?we?should?detect?features.??
6. //?Boolean?flag.?If?it?is?true,?the?keypoint?detector?is?not?run.?Instead,??
7. //?the?provided?vector?of?keypoints?is?used?and?the?algorithm?just?computes?their?descriptors.??
8. //?descriptors?–?The?output?matrix?of?descriptors.??
9. //?Pass?cv::noArray()?if?you?do?not?need?them.??????????????
10. {??
11. ????Mat?image?=?_image.getMat(),?mask?=?_mask.getMat();??
12. ??
13. ????if(?image.empty()?||?image.depth()?!=?CV_8U?)??
14. ????????CV_Error(?CV_StsBadArg,?"image?is?empty?or?has?incorrect?depth?(!=CV_8U)"?);??
15. ??
16. ????if(?!mask.empty()?&&?mask.type()?!=?CV_8UC1?)??
17. ????????CV_Error(?CV_StsBadArg,?"mask?has?incorrect?type?(!=CV_8UC1)"?);??
18. ??
19. ??????????
20. ????//?得到第1組（Octave）圖像??
21. ????Mat?base?=?createInitialImage(image,?false,?(float)sigma);??
22. ????vector<Mat>?gpyr,?dogpyr;??
23. ????//?每層金字塔圖像的組數（Octave）??
24. ????int?nOctaves?=?cvRound(log(?(double)std::min(?base.cols,?base.rows?)?)?/?log(2.)?-?2);??
25. ??
26. ????//?double?t,?tf?=?getTickFrequency();??
27. ????//?t?=?(double)getTickCount();??
28. ??????
29. ????//?構建金字塔（金字塔層數和組數相等）??
30. ????buildGaussianPyramid(base,?gpyr,?nOctaves);??
31. ????//?構建高斯差分金字塔??
32. ????buildDoGPyramid(gpyr,?dogpyr);??
33. ??
34. ????//t?=?(double)getTickCount()?-?t;??
35. ????//printf("pyramid?construction?time:?%g\n",?t*1000./tf);??
36. ??????
37. ????//?useProvidedKeypoints默認為false??
38. ????//?使用keypoints并計算特征點的描述符??
39. ????if(?!useProvidedKeypoints?)??
40. ????{??
41. ????????//t?=?(double)getTickCount();??
42. ????????findScaleSpaceExtrema(gpyr,?dogpyr,?keypoints);??
43. ????????//除去重復特征點??
44. ????????KeyPointsFilter::removeDuplicated(?keypoints?);???
45. ??
46. ????????//?mask標記檢測區域（可選）??
47. ????????if(?!mask.empty()?)??
48. ????????????KeyPointsFilter::runByPixelsMask(?keypoints,?mask?);??
49. ??
50. ????????//?retainBest:根據相應保留指定數目的特征點（features2d.hpp）??
51. ????????if(?nfeatures?>?0?)??
52. ????????????KeyPointsFilter::retainBest(keypoints,?nfeatures);??
53. ????????//t?=?(double)getTickCount()?-?t;??
54. ????????//printf("keypoint?detection?time:?%g\n",?t*1000./tf);??
55. ????}??
56. ????else??
57. ????{??
58. ????????//?filter?keypoints?by?mask??
59. ????????//?KeyPointsFilter::runByPixelsMask(?keypoints,?mask?);??
60. ????}??
61. ??
62. ????//?特征點輸出數組??
63. ????if(?_descriptors.needed()?)??
64. ????{??
65. ????????//t?=?(double)getTickCount();??
66. ????????int?dsize?=?descriptorSize();??
67. ????????_descriptors.create((int)keypoints.size(),?dsize,?CV_32F);??
68. ????????Mat?descriptors?=?_descriptors.getMat();??
69. ??
70. ????????calcDescriptors(gpyr,?keypoints,?descriptors,?nOctaveLayers);??
71. ????????//t?=?(double)getTickCount()?-?t;??
72. ????????//printf("descriptor?extraction?time:?%g\n",?t*1000./tf);??
73. ????}??
74. }??

[8]

總結

以上是生活随笔為你收集整理的尺度空间理论与图像金字塔的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。