OpenCV模板匹配算法详解
1 理论介绍
模板匹配是在⼀幅图像中寻⼀个特定⽬标的⽅法之⼀,这种⽅法的原理⾮常简单,遍历图像中的每⼀个可能的位置,⽐较各处与模板是否“相似”,当相似度⾜够⾼时,就认为到了我们的⽬标。OpenCV提供了6种模板匹配算法:
1. 平⽅差匹配法CV_TM_SQDIFF
2. 归⼀化平⽅差匹配法CV_TM_SQDIFF_NORMED
3. 相关匹配法CV_TM_CCORR
4. 归⼀化相关匹配法CV_TM_CCORR_NORMED
5. 相关系数匹配法CV_TM_CCOEFF
6. 归⼀化相关系数匹配法CV_TM_CCOEFF_NORMED
⽤T表⽰模板图像,I表⽰待匹配图像,切模板图像的宽为w⾼为h,⽤R表⽰匹配结果,匹配过程如下图
所⽰:
上述6中匹配⽅法可⽤以下公式进⾏描述:
2 ⽰例代码
下⾯给出⽅法6的python代码
1import numpy as np
2import cv2
3
4def EM(pModel, width, height):
5    sum = np.double(0.0)
6for i in range(0,height):
7for j in range(0,width):
8            sum += pModel[i][j]
9return sum
10
11def EM2(pModel, width, height):
12    sum = np.double(0.0)
13for i in range(0,height):
14for j in range(0,width):
15            sum += pModel[i][j]*1.0*pModel[i][j]
16return sum
17
18def EI(pToSearch, l, h, u, v, pModel, width, height):
19    sum = np.double(0.0)
20    roi = pToSearch[v:v+height, u:u+width]
21for i in range(0,height):
22for j in range(0,width):
23            sum += roi[i][j]
24return sum
25
26def EI2(pToSearch, l, h, u, v, pModel, width, height):
27    sum = np.double(0.0)
28    roi = pToSearch[v:v+height, u:u+width]
29for i in range(0,height):
30for j in range(0,width):
31            sum += roi[i][j]*1.0*roi[i][j]
32return sum
33
34def EIM(pToSearch, l, h, u, v, pModel, width, height):
35    sum = np.double(0.0)
36    roi = pToSearch[v:v+height, u:u+width]
37for i in range(0,height):
38for j in range(0,width):
39            sum += pModel[i][j]*1.0*roi[i][j]
40return sum
41
42def Match(pToSearch, l, h, pModel, width, height):
43    uMax = l-width
44    vMax = h-height
45    N = width*height
46    len = (uMax+1)*(vMax+1)
47    MatchRec = [0.0 for x in range(0, len)]
48    k = 0
49
50    M = EM(pModel,width,height)
51    M2 = EM2(pModel,width,height)
52for p in range(0, uMax+1):
53for q in range(0, vMax+1):
54            I = EI(pToSearch,l,h,p,q,pModel,width,height)
55            I2 = EI2(pToSearch,l,h,p,q,pModel,width,height)
56            IM = EIM(pToSearch,l,h,p,q,pModel,width,height)
57
58            numerator=(N*IM-I*M)*(N*IM-I*M)
59            denominator=(N*I2-I*I)*(N*M2-M*M)
60
61            ret = numerator/denominator
rectangle函数opencv62            MatchRec[k]=ret
63            k+=1
64
65    val = 0
66    k = 0
67    x = y = 0
68for p in range(0, uMax+1):
69for q in range(0, vMax+1):
70if MatchRec[k] > val:
71                val = MatchRec[k]
72                x = p
73                y = q
74            k+=1
75print"val: %f"%val
76return (x, y)
77
78def main():
79    img = cv2.imread('niu.jpg', cv2.IMREAD_GRAYSCALE)
80    temp = cv2.imread('temp.png', cv2.IMREAD_GRAYSCALE)
81
82print temp.shape
83    imgHt, imgWd = img.shape
84    tempHt, tempWd = temp.shape
85#print EM(temp, tempWd, tempHt)
86    (x, y) = Match(img, imgWd, imgHt, temp, tempWd, tempHt)
87    angle(img, (x, y), (x+tempWd, y+tempHt), (0,0,0), 2)
88    cv2.imshow("temp", temp)
89    cv2.imshow("result", img)
90    cv2.waitKey(0)
91    cv2.destroyAllWindows()
92
93if__name__ == '__main__':
94    main()
归⼀化相关系数匹配法
代码58⾏中的N就是公式(6)中的w*h,由于python代码运⾏速度⽐较慢,代码的58、59⾏相当于对公式(6)的分⼦分母都进⾏了平⽅操作,并且分⼦
分母都乘以了N⽅,以减⼩计算量,所以代码61⾏的ret相当于公式(6)中的R(x,y)的平⽅,
为了更快的进⾏算法验证,⽤上述代码进⾏验证时请尽量选⽤较⼩的匹配图像及模板图像,下图显⽰了我的匹配结果(待匹配图像295x184模板
69x46⽤了⼗⼏分钟):
3 OpenCV源码
较新版本的OpenCV库中的模板匹配已经进⾏了较多的算法改进,直接看新版本中的算法需要了解很多相关理论知识,所以我们结合OpenCV0.9.5的源码进⾏讲解,该版本的源码基本上是C风格代码更容易进⾏理解(如果要对
OpenCV源码进⾏研究,建议⽤该版本进⾏⼊门),仍以归⼀化相关系数匹配法为例进⾏分析。
1/*
2* pImage: 待匹配图像
3* image: 待匹配图像宽(width*depth并已4字节对齐)
4* roiSize: 待匹配图像尺⼨
5* pTemplate: 模板图像
6* templStep: 模板图像宽
7* templSize: 模板图像尺⼨
8* pResult: 匹配结果
9* resultStep: 匹配结果宽
10* pBuffer: 中间结果数据缓存
11*/
12 IPCVAPI_IMPL( CvStatus, icvMatchTemplate_CoeffNormed_32f_C1R,
13              (const float *pImage, int imageStep, CvSize roiSize,
14const float *pTemplate, int templStep, CvSize templSize,
15float *pResult, int resultStep, void *pBuffer) )
16 {
17float *imgBuf = 0;              // 待匹配图像相关数据
18float *templBuf = 0;            // 模板图像数据
19double *sumBuf = 0;            // 待匹配图像遍历块单⾏和
20double *sqsumBuf = 0;          // 待匹配图像遍历块单⾏平⽅和
21double *resNum = 0;            // 模板图像和待匹配图像遍历块内积
22double *resDenom = 0;          // 待匹配图像遍历块累加和及待匹配图像遍历块平⽅累加和
23double templCoeff = 0;          // 模板图像均分差倒数
24double templSum = 0;            // 模板图像累加和
25
26int winLen = templSize.width * templSize.height;
27double winCoeff = 1. / (winLen + DBL_EPSILON);          // + DBL_EPSILON 加⼀个⼩整数防⽌分母为零
28
29    CvSize resultSize = cvSize( roiSize.width - templSize.width + 1,
30                                roiSize.height - templSize.height + 1 );
31int x, y;
32
33// 计算并为imgBuf、templBuf、sumBuf、sqsumBuf、resNum、resDenom分配存储空间
34    CvStatus result = icvMatchTemplateEntry( pImage, imageStep, roiSize,
35                                              pTemplate, templStep, templSize,
36                                              pResult, resultStep, pBuffer,
37                                              cv32f, 1, 1,
38                                              (void **) &imgBuf, (void **) &templBuf,
39                                              (void **) &sumBuf, (void **) &sqsumBuf,
40                                              (void **) &resNum, (void **) &resDenom );
41
42if( result != CV_OK )
43return result;
44
45    imageStep /= sizeof_float;
46    templStep /= sizeof_float;
47    resultStep /= sizeof_float;
48
49/* calc common statistics for template and image */
50    {
51const float *rowPtr = (const float *) imgBuf;
52double templSqsum = icvCrossCorr_32f_C1( templBuf, templBuf, winLen );          // 模板图像平⽅累加和
53
54        templSum = icvSumPixels_32f_C1( templBuf, winLen );                            // 模板图像累加和
55        templCoeff = (double) templSqsum - ((double) templSum) * templSum * winCoeff;  // 模板图像均⽅差的平⽅
56        templCoeff = icvInvSqrt64d( fabs( templCoeff ) + FLT_EPSILON );                // 模板图像均⽅差倒数
57
58for( y = 0; y < roiSize.height; y++, rowPtr += templSize.width )
59        {
60            sumBuf[y] = icvSumPixels_32f_C1( rowPtr, templSize.width );                // 待匹配图像按模板图像宽度求每⾏之和(遍历位置第⼀列)
61            sqsumBuf[y] = icvCrossCorr_32f_C1( rowPtr, rowPtr, templSize.width );      // 待匹配图像按模板图像宽度求每⾏平⽅之和(遍历位置第⼀列)
62        }
63    }
64
65/* main loop - through x coordinate of the result */
66for( x = 0; x < resultSize.width; x++ )
67    {
68double sum = 0;
69double sqsum = 0;
70float *imgPtr = imgBuf + x;                                                      // 待匹配图像起始位置
71
72/* update sums and image band buffer */// 如果不是第1列需重新更新sumBuf,更新后sumBuf为遍历位置第x列每⾏之和(⾏宽为模板图像宽) 73if( x > 0 )
74        {
75const float *src = pImage + x + templSize.width - 1;
76float *dst = imgPtr - 1;
77float out_val = dst[0];
78
79            dst += templSize.width;
80
81for( y = 0; y < roiSize.height; y++, src += imageStep, dst += templSize.width )
82            {
83float in_val = src[0];
84
85                sumBuf[y] += in_val - out_val;
86                sqsumBuf[y] += (in_val - out_val) * (in_val + out_val);
87                out_val = dst[0];
88                dst[0] = (float) in_val;
89            }
90        }
91
92for( y = 0; y < templSize.height; y++ )                                          // 求遍历位置第x列,第1⾏处遍历块累加和sum及平⽅累加和sqsum
93        {
94            sum += sumBuf[y];
95            sqsum += sqsumBuf[y];
96        }
97
98for( y = 0; y < resultSize.height; y++, imgPtr += templSize.width )
99        {
100double res = icvCrossCorr_32f_C1( imgPtr, templBuf, winLen );              // 求模板图像和待匹配图像y⾏x列处遍历块的内积
101
102if( y > 0 )                                                                // 如果不是第1⾏需更新遍历块累加和sum及平⽅累加和sqsum
103            {
104                sum -= sumBuf[y - 1];
105                sum += sumBuf[y + templSize.height - 1];
106                sqsum -= sqsumBuf[y - 1];
107                sqsum += sqsumBuf[y + templSize.height - 1];
108            }
109            resNum[y] = res;
110            resDenom[y] = sum;
111            resDenom[y + resultSize.height] = sqsum;
112        }
113
114for( y = 0; y < resultSize.height; y++ )
115        {
116double sum = ((double) resDenom[y]);
117double wsum = winCoeff * sum;
118double res = ((double) resNum[y]) - wsum * templSum;
119double nrm_s = ((double) resDenom[y + resultSize.height]) - wsum * sum;
120
121            res *= templCoeff * icvInvSqrt64d( fabs( nrm_s ) + FLT_EPSILON );
122            pResult[x + y * resultStep] = (float) res;
123        }
124    }
125
126return CV_OK;
127 }
以上代码是归⼀化相关系数法核⼼函数icvMatchTemplate_CoeffNormed_32f_C1R的源码,我已经在源码中进⾏了详细的注释,读者需⾃⼰再进⾏理解,需要进⼀步说明的是:
代码118⾏res就是计算公式(6)的分⼦部分,代码56⾏templCoeff就是计算公式(6)分母的左半部分,代码121⾏icvInvSqrt64d函数就是在计算公式(6)分母的右半部分,该⾏res的最终结果正是公式(6)中的R(x,y)。
4 结束语
OpenCV0.9.5源码下载:
参考⽂章:

版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198,我们将在24小时内删除。