git.maemo.org Git - opencv/blob - src/cv/cvtemplmatch.cpp

   1 /*M///////////////////////////////////////////////////////////////////////////////////////\r
   2 //\r
   3 //  IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.\r
   4 //\r
   5 //  By downloading, copying, installing or using the software you agree to this license.\r
   6 //  If you do not agree to this license, do not download, install,\r
   7 //  copy or use the software.\r
   8 //\r
   9 //\r
  10 //                        Intel License Agreement\r
  11 //                For Open Source Computer Vision Library\r
  12 //\r
  13 // Copyright (C) 2000, Intel Corporation, all rights reserved.\r
  14 // Third party copyrights are property of their respective owners.\r
  15 //\r
  16 // Redistribution and use in source and binary forms, with or without modification,\r
  17 // are permitted provided that the following conditions are met:\r
  18 //\r
  19 //   * Redistribution's of source code must retain the above copyright notice,\r
  20 //     this list of conditions and the following disclaimer.\r
  21 //\r
  22 //   * Redistribution's in binary form must reproduce the above copyright notice,\r
  23 //     this list of conditions and the following disclaimer in the documentation\r
  24 //     and/or other materials provided with the distribution.\r
  25 //\r
  26 //   * The name of Intel Corporation may not be used to endorse or promote products\r
  27 //     derived from this software without specific prior written permission.\r
  28 //\r
  29 // This software is provided by the copyright holders and contributors "as is" and\r
  30 // any express or implied warranties, including, but not limited to, the implied\r
  31 // warranties of merchantability and fitness for a particular purpose are disclaimed.\r
  32 // In no event shall the Intel Corporation or contributors be liable for any direct,\r
  33 // indirect, incidental, special, exemplary, or consequential damages\r
  34 // (including, but not limited to, procurement of substitute goods or services;\r
  35 // loss of use, data, or profits; or business interruption) however caused\r
  36 // and on any theory of liability, whether in contract, strict liability,\r
  37 // or tort (including negligence or otherwise) arising in any way out of\r
  38 // the use of this software, even if advised of the possibility of such damage.\r
  39 //\r
  40 //M*/\r
  41 \r
  42 #include "_cv.h"\r
  43 \r
  44 void\r
  45 icvCrossCorr( const CvArr* _img, const CvArr* _templ, CvArr* _corr,\r
  46               CvPoint anchor, double delta, int borderType )\r
  47 {\r
  48     const double block_scale = 4.5;\r
  49     const int min_block_size = 256;\r
  50     CvMat* dft_img[CV_MAX_THREADS] = {0};\r
  51     CvMat* dft_templ = 0;\r
  52     void* buf[CV_MAX_THREADS] = {0};\r
  53     int k, num_threads = 0;\r
  54     \r
  55     CV_FUNCNAME( "icvCrossCorr" );\r
  56     \r
  57     __BEGIN__;\r
  58 \r
  59     CvMat istub, *img = (CvMat*)_img;\r
  60     CvMat tstub, *templ = (CvMat*)_templ;\r
  61     CvMat cstub, *corr = (CvMat*)_corr;\r
  62     CvSize dftsize, blocksize;\r
  63     int depth, templ_depth, corr_depth, max_depth = CV_32F,\r
  64         cn, templ_cn, corr_cn, buf_size = 0,\r
  65         tile_count_x, tile_count_y, tile_count;\r
  66 \r
  67     CV_CALL( img = cvGetMat( img, &istub ));\r
  68     CV_CALL( templ = cvGetMat( templ, &tstub ));\r
  69     CV_CALL( corr = cvGetMat( corr, &cstub ));\r
  70 \r
  71     if( CV_MAT_DEPTH( img->type ) != CV_8U &&\r
  72         CV_MAT_DEPTH( img->type ) != CV_16U &&\r
  73         CV_MAT_DEPTH( img->type ) != CV_32F )\r
  74         CV_ERROR( CV_StsUnsupportedFormat,\r
  75         "The function supports only 8u, 16u and 32f data types" );\r
  76 \r
  77     if( !CV_ARE_DEPTHS_EQ( img, templ ) && CV_MAT_DEPTH( templ->type ) != CV_32F )\r
  78         CV_ERROR( CV_StsUnsupportedFormat,\r
  79         "Template (kernel) must be of the same depth as the input image, or be 32f" );\r
  80     \r
  81     if( !CV_ARE_DEPTHS_EQ( img, corr ) && CV_MAT_DEPTH( corr->type ) != CV_32F &&\r
  82         CV_MAT_DEPTH( corr->type ) != CV_64F )\r
  83         CV_ERROR( CV_StsUnsupportedFormat,\r
  84         "The output image must have the same depth as the input image, or be 32f/64f" );\r
  85 \r
  86     if( (!CV_ARE_CNS_EQ( img, corr ) || CV_MAT_CN(templ->type) > 1) &&\r
  87         (CV_MAT_CN( corr->type ) > 1 || !CV_ARE_CNS_EQ( img, templ)) )\r
  88         CV_ERROR( CV_StsUnsupportedFormat,\r
  89         "The output must have the same number of channels as the input (when the template has 1 channel), "\r
  90         "or the output must have 1 channel when the input and the template have the same number of channels" );\r
  91 \r
  92     depth = CV_MAT_DEPTH(img->type);\r
  93     cn = CV_MAT_CN(img->type);\r
  94     templ_depth = CV_MAT_DEPTH(templ->type);\r
  95     templ_cn = CV_MAT_CN(templ->type);\r
  96     corr_depth = CV_MAT_DEPTH(corr->type);\r
  97     corr_cn = CV_MAT_CN(corr->type);\r
  98 \r
  99     CV_Assert( corr_cn == 1 || delta == 0 );\r
 100 \r
 101     max_depth = MAX( max_depth, templ_depth );\r
 102     max_depth = MAX( max_depth, depth );\r
 103     max_depth = MAX( max_depth, corr_depth );\r
 104     if( depth > CV_8U )\r
 105         max_depth = CV_64F;\r
 106 \r
 107     if( img->cols < templ->cols || img->rows < templ->rows )\r
 108         CV_ERROR( CV_StsUnmatchedSizes,\r
 109         "Such a combination of image and template/filter size is not supported" );\r
 110 \r
 111     if( corr->rows > img->rows + templ->rows - 1 ||\r
 112         corr->cols > img->cols + templ->cols - 1 )\r
 113         CV_ERROR( CV_StsUnmatchedSizes,\r
 114         "output image should not be greater than (W + w - 1)x(H + h - 1)" );\r
 115 \r
 116     blocksize.width = cvRound(templ->cols*block_scale);\r
 117     blocksize.width = MAX( blocksize.width, min_block_size - templ->cols + 1 );\r
 118     blocksize.width = MIN( blocksize.width, corr->cols );\r
 119     blocksize.height = cvRound(templ->rows*block_scale);\r
 120     blocksize.height = MAX( blocksize.height, min_block_size - templ->rows + 1 );\r
 121     blocksize.height = MIN( blocksize.height, corr->rows );\r
 122 \r
 123     dftsize.width = cvGetOptimalDFTSize(blocksize.width + templ->cols - 1);\r
 124     if( dftsize.width == 1 )\r
 125         dftsize.width = 2;\r
 126     dftsize.height = cvGetOptimalDFTSize(blocksize.height + templ->rows - 1);\r
 127     if( dftsize.width <= 0 || dftsize.height <= 0 )\r
 128         CV_ERROR( CV_StsOutOfRange, "the input arrays are too big" );\r
 129 \r
 130     // recompute block size\r
 131     blocksize.width = dftsize.width - templ->cols + 1;\r
 132     blocksize.width = MIN( blocksize.width, corr->cols );\r
 133     blocksize.height = dftsize.height - templ->rows + 1;\r
 134     blocksize.height = MIN( blocksize.height, corr->rows );\r
 135 \r
 136     CV_CALL( dft_templ = cvCreateMat( dftsize.height*templ_cn, dftsize.width, max_depth ));\r
 137 \r
 138     num_threads = cvGetNumThreads();\r
 139 \r
 140     for( k = 0; k < num_threads; k++ )\r
 141         CV_CALL( dft_img[k] = cvCreateMat( dftsize.height, dftsize.width, max_depth ));\r
 142 \r
 143     if( templ_cn > 1 && templ_depth != max_depth )\r
 144         buf_size = templ->cols*templ->rows*CV_ELEM_SIZE(templ_depth);\r
 145 \r
 146     if( cn > 1 && depth != max_depth )\r
 147         buf_size = MAX( buf_size, (blocksize.width + templ->cols - 1)*\r
 148             (blocksize.height + templ->rows - 1)*CV_ELEM_SIZE(depth));\r
 149 \r
 150     if( (corr_cn > 1 || cn > 1) && corr_depth != max_depth )\r
 151         buf_size = MAX( buf_size, blocksize.width*blocksize.height*CV_ELEM_SIZE(corr_depth));\r
 152 \r
 153     if( buf_size > 0 )\r
 154     {\r
 155         for( k = 0; k < num_threads; k++ )\r
 156             CV_CALL( buf[k] = cvAlloc(buf_size) );\r
 157     }\r
 158 \r
 159     // compute DFT of each template plane\r
 160     for( k = 0; k < templ_cn; k++ )\r
 161     {\r
 162         CvMat dstub, *src, *dst, temp;\r
 163         CvMat* planes[] = { 0, 0, 0, 0 };\r
 164         int yofs = k*dftsize.height;\r
 165 \r
 166         src = templ;\r
 167         dst = cvGetSubRect( dft_templ, &dstub, cvRect(0,yofs,templ->cols,templ->rows));\r
 168     \r
 169         if( templ_cn > 1 )\r
 170         {\r
 171             planes[k] = templ_depth == max_depth ? dst :\r
 172                 cvInitMatHeader( &temp, templ->rows, templ->cols, templ_depth, buf[0] );\r
 173             cvSplit( templ, planes[0], planes[1], planes[2], planes[3] );\r
 174             src = planes[k];\r
 175             planes[k] = 0;\r
 176         }\r
 177 \r
 178         if( dst != src )\r
 179             cvConvert( src, dst );\r
 180 \r
 181         if( dft_templ->cols > templ->cols )\r
 182         {\r
 183             cvGetSubRect( dft_templ, dst, cvRect(templ->cols, yofs,\r
 184                           dft_templ->cols - templ->cols, templ->rows) );\r
 185             cvZero( dst );\r
 186         }\r
 187         cvGetSubRect( dft_templ, dst, cvRect(0,yofs,dftsize.width,dftsize.height) );\r
 188         cvDFT( dst, dst, CV_DXT_FORWARD + CV_DXT_SCALE, templ->rows );\r
 189     }\r
 190 \r
 191     tile_count_x = (corr->cols + blocksize.width - 1)/blocksize.width;\r
 192     tile_count_y = (corr->rows + blocksize.height - 1)/blocksize.height;\r
 193     tile_count = tile_count_x*tile_count_y;\r
 194 \r
 195     {\r
 196 #ifdef _OPENMP\r
 197     #pragma omp parallel for num_threads(num_threads) schedule(dynamic)\r
 198 #endif\r
 199     // calculate correlation by blocks\r
 200     for( k = 0; k < tile_count; k++ )\r
 201     {\r
 202         int thread_idx = cvGetThreadNum();\r
 203         int x = (k%tile_count_x)*blocksize.width;\r
 204         int y = (k/tile_count_x)*blocksize.height;\r
 205         int i, yofs;\r
 206         CvMat sstub, dstub, *src, *dst, temp;\r
 207         CvMat* planes[] = { 0, 0, 0, 0 };\r
 208         CvMat* _dft_img = dft_img[thread_idx];\r
 209         void* _buf = buf[thread_idx];\r
 210         CvSize csz = { blocksize.width, blocksize.height }, isz;\r
 211         int x0 = x - anchor.x, y0 = y - anchor.y;\r
 212         int x1 = MAX( 0, x0 ), y1 = MAX( 0, y0 ), x2, y2;\r
 213         csz.width = MIN( csz.width, corr->cols - x );\r
 214         csz.height = MIN( csz.height, corr->rows - y );\r
 215         isz.width = csz.width + templ->cols - 1;\r
 216         isz.height = csz.height + templ->rows - 1;\r
 217         x2 = MIN( img->cols, x0 + isz.width );\r
 218         y2 = MIN( img->rows, y0 + isz.height );\r
 219         \r
 220         for( i = 0; i < cn; i++ )\r
 221         {\r
 222             CvMat dstub1, *dst1;\r
 223             yofs = i*dftsize.height;\r
 224 \r
 225             src = cvGetSubRect( img, &sstub, cvRect(x1,y1,x2-x1,y2-y1) );\r
 226             dst = cvGetSubRect( _dft_img, &dstub,\r
 227                 cvRect(0,0,isz.width,isz.height) );\r
 228             dst1 = dst;\r
 229             \r
 230             if( x2 - x1 < isz.width || y2 - y1 < isz.height )\r
 231                 dst1 = cvGetSubRect( _dft_img, &dstub1,\r
 232                     cvRect( x1 - x0, y1 - y0, x2 - x1, y2 - y1 ));\r
 233 \r
 234             if( cn > 1 )\r
 235             {\r
 236                 planes[i] = dst1;\r
 237                 if( depth != max_depth )\r
 238                     planes[i] = cvInitMatHeader( &temp, y2 - y1, x2 - x1, depth, _buf );\r
 239                 cvSplit( src, planes[0], planes[1], planes[2], planes[3] );\r
 240                 src = planes[i];\r
 241                 planes[i] = 0;\r
 242             }\r
 243 \r
 244             if( dst1 != src )\r
 245                 cvConvert( src, dst1 );\r
 246 \r
 247             if( dst != dst1 )\r
 248                 cvCopyMakeBorder( dst1, dst, cvPoint(x1 - x0, y1 - y0), borderType );\r
 249 \r
 250             if( dftsize.width > isz.width )\r
 251             {\r
 252                 cvGetSubRect( _dft_img, dst, cvRect(isz.width, 0,\r
 253                       dftsize.width - isz.width,dftsize.height) );\r
 254                 cvZero( dst );\r
 255             }\r
 256 \r
 257             cvDFT( _dft_img, _dft_img, CV_DXT_FORWARD, isz.height );\r
 258             cvGetSubRect( dft_templ, dst,\r
 259                 cvRect(0,(templ_cn>1?yofs:0),dftsize.width,dftsize.height) );\r
 260 \r
 261             cvMulSpectrums( _dft_img, dst, _dft_img, CV_DXT_MUL_CONJ );\r
 262             cvDFT( _dft_img, _dft_img, CV_DXT_INVERSE, csz.height );\r
 263 \r
 264             src = cvGetSubRect( _dft_img, &sstub, cvRect(0,0,csz.width,csz.height) );\r
 265             dst = cvGetSubRect( corr, &dstub, cvRect(x,y,csz.width,csz.height) );\r
 266 \r
 267             if( corr_cn > 1 )\r
 268             {\r
 269                 planes[i] = src;\r
 270                 if( corr_depth != max_depth )\r
 271                 {\r
 272                     planes[i] = cvInitMatHeader( &temp, csz.height, csz.width,\r
 273                                                  corr_depth, _buf );\r
 274                     cvConvertScale( src, planes[i], 1, delta );\r
 275                 }\r
 276                 cvMerge( planes[0], planes[1], planes[2], planes[3], dst );\r
 277                 planes[i] = 0;                    \r
 278             }\r
 279             else\r
 280             {\r
 281                 if( i == 0 )\r
 282                     cvConvertScale( src, dst, 1, delta );\r
 283                 else\r
 284                 {\r
 285                     if( max_depth > corr_depth )\r
 286                     {\r
 287                         cvInitMatHeader( &temp, csz.height, csz.width,\r
 288                                          corr_depth, _buf );\r
 289                         cvConvert( src, &temp );\r
 290                         src = &temp;\r
 291                     }\r
 292                     cvAcc( src, dst );\r
 293                 }\r
 294             }\r
 295         }\r
 296     }\r
 297     }\r
 298 \r
 299     __END__;\r
 300 \r
 301     cvReleaseMat( &dft_templ );\r
 302 \r
 303     for( k = 0; k < num_threads; k++ )\r
 304     {\r
 305         cvReleaseMat( &dft_img[k] );\r
 306         cvFree( &buf[k] );\r
 307     }\r
 308 }\r
 309 \r
 310 void\r
 311 cv::crossCorr( const Mat& img, const Mat& templ, Mat& corr,\r
 312            Point anchor, double delta, int borderType )\r
 313 {\r
 314     CvMat _img = img, _templ = templ, _corr = corr;\r
 315     icvCrossCorr( &_img, &_templ, &_corr, anchor, delta, borderType );\r
 316 }\r
 317 \r
 318 \r
 319 /*****************************************************************************************/\r
 320 \r
 321 CV_IMPL void\r
 322 cvMatchTemplate( const CvArr* _img, const CvArr* _templ, CvArr* _result, int method )\r
 323 {\r
 324     CvMat* sum = 0;\r
 325     CvMat* sqsum = 0;\r
 326     \r
 327     CV_FUNCNAME( "cvMatchTemplate" );\r
 328 \r
 329     __BEGIN__;\r
 330 \r
 331     int coi1 = 0, coi2 = 0;\r
 332     int depth, cn;\r
 333     int i, j, k;\r
 334     CvMat stub, *img = (CvMat*)_img;\r
 335     CvMat tstub, *templ = (CvMat*)_templ;\r
 336     CvMat rstub, *result = (CvMat*)_result;\r
 337     CvScalar templ_mean = cvScalarAll(0);\r
 338     double templ_norm = 0, templ_sum2 = 0;\r
 339     \r
 340     int idx = 0, idx2 = 0;\r
 341     double *p0, *p1, *p2, *p3;\r
 342     double *q0, *q1, *q2, *q3;\r
 343     double inv_area;\r
 344     int sum_step, sqsum_step;\r
 345     int num_type = method == CV_TM_CCORR || method == CV_TM_CCORR_NORMED ? 0 :\r
 346                    method == CV_TM_CCOEFF || method == CV_TM_CCOEFF_NORMED ? 1 : 2;\r
 347     int is_normed = method == CV_TM_CCORR_NORMED ||\r
 348                     method == CV_TM_SQDIFF_NORMED ||\r
 349                     method == CV_TM_CCOEFF_NORMED;\r
 350 \r
 351     CV_CALL( img = cvGetMat( img, &stub, &coi1 ));\r
 352     CV_CALL( templ = cvGetMat( templ, &tstub, &coi2 ));\r
 353     CV_CALL( result = cvGetMat( result, &rstub ));\r
 354 \r
 355     if( CV_MAT_DEPTH( img->type ) != CV_8U &&\r
 356         CV_MAT_DEPTH( img->type ) != CV_32F )\r
 357         CV_ERROR( CV_StsUnsupportedFormat,\r
 358         "The function supports only 8u and 32f data types" );\r
 359 \r
 360     if( !CV_ARE_TYPES_EQ( img, templ ))\r
 361         CV_ERROR( CV_StsUnmatchedSizes, "image and template should have the same type" );\r
 362 \r
 363     if( CV_MAT_TYPE( result->type ) != CV_32FC1 )\r
 364         CV_ERROR( CV_StsUnsupportedFormat, "output image should have 32f type" );\r
 365 \r
 366     if( img->rows < templ->rows || img->cols < templ->cols )\r
 367     {\r
 368         CvMat* t;\r
 369         CV_SWAP( img, templ, t );\r
 370     }\r
 371 \r
 372     if( result->rows != img->rows - templ->rows + 1 ||\r
 373         result->cols != img->cols - templ->cols + 1 )\r
 374         CV_ERROR( CV_StsUnmatchedSizes, "output image should be (W - w + 1)x(H - h + 1)" );\r
 375 \r
 376     if( method < CV_TM_SQDIFF || method > CV_TM_CCOEFF_NORMED )\r
 377         CV_ERROR( CV_StsBadArg, "unknown comparison method" );\r
 378 \r
 379     depth = CV_MAT_DEPTH(img->type);\r
 380     cn = CV_MAT_CN(img->type);\r
 381 \r
 382     /*if( is_normed && cn == 1 && templ->rows > 8 && templ->cols > 8 &&\r
 383         img->rows > templ->cols && img->cols > templ->cols )\r
 384     {\r
 385         CvTemplMatchIPPFunc ipp_func =\r
 386             depth == CV_8U ?\r
 387             (method == CV_TM_SQDIFF_NORMED ? (CvTemplMatchIPPFunc)icvSqrDistanceValid_Norm_8u32f_C1R_p :\r
 388             method == CV_TM_CCORR_NORMED ? (CvTemplMatchIPPFunc)icvCrossCorrValid_Norm_8u32f_C1R_p :\r
 389             (CvTemplMatchIPPFunc)icvCrossCorrValid_NormLevel_8u32f_C1R_p) :\r
 390             (method == CV_TM_SQDIFF_NORMED ? (CvTemplMatchIPPFunc)icvSqrDistanceValid_Norm_32f_C1R_p :\r
 391             method == CV_TM_CCORR_NORMED ? (CvTemplMatchIPPFunc)icvCrossCorrValid_Norm_32f_C1R_p :\r
 392             (CvTemplMatchIPPFunc)icvCrossCorrValid_NormLevel_32f_C1R_p);\r
 393 \r
 394         if( ipp_func )\r
 395         {\r
 396             CvSize img_size = cvGetMatSize(img), templ_size = cvGetMatSize(templ);\r
 397 \r
 398             IPPI_CALL( ipp_func( img->data.ptr, img->step ? img->step : CV_STUB_STEP,\r
 399                                  img_size, templ->data.ptr,\r
 400                                  templ->step ? templ->step : CV_STUB_STEP,\r
 401                                  templ_size, result->data.ptr,\r
 402                                  result->step ? result->step : CV_STUB_STEP ));\r
 403             for( i = 0; i < result->rows; i++ )\r
 404             {\r
 405                 float* rrow = (float*)(result->data.ptr + i*result->step);\r
 406                 for( j = 0; j < result->cols; j++ )\r
 407                 {\r
 408                     if( fabs(rrow[j]) > 1. )\r
 409                         rrow[j] = rrow[j] < 0 ? -1.f : 1.f;\r
 410                 }\r
 411             }\r
 412             EXIT;\r
 413         }\r
 414     }*/\r
 415 \r
 416     CV_CALL( icvCrossCorr( img, templ, result ));\r
 417 \r
 418     if( method == CV_TM_CCORR )\r
 419         EXIT;\r
 420 \r
 421     inv_area = 1./((double)templ->rows * templ->cols);\r
 422 \r
 423     CV_CALL( sum = cvCreateMat( img->rows + 1, img->cols + 1,\r
 424                                 CV_MAKETYPE( CV_64F, cn )));\r
 425     if( method == CV_TM_CCOEFF )\r
 426     {\r
 427         CV_CALL( cvIntegral( img, sum, 0, 0 ));\r
 428         CV_CALL( templ_mean = cvAvg( templ ));\r
 429         q0 = q1 = q2 = q3 = 0;\r
 430     }\r
 431     else\r
 432     {\r
 433         CvScalar _templ_sdv = cvScalarAll(0);\r
 434         CV_CALL( sqsum = cvCreateMat( img->rows + 1, img->cols + 1,\r
 435                                       CV_MAKETYPE( CV_64F, cn )));\r
 436         CV_CALL( cvIntegral( img, sum, sqsum, 0 ));\r
 437         CV_CALL( cvAvgSdv( templ, &templ_mean, &_templ_sdv ));\r
 438 \r
 439         templ_norm = CV_SQR(_templ_sdv.val[0]) + CV_SQR(_templ_sdv.val[1]) +\r
 440                     CV_SQR(_templ_sdv.val[2]) + CV_SQR(_templ_sdv.val[3]);\r
 441 \r
 442         if( templ_norm < DBL_EPSILON && method == CV_TM_CCOEFF_NORMED )\r
 443         {\r
 444             cvSet( result, cvScalarAll(1.) );\r
 445             EXIT;\r
 446         }\r
 447         \r
 448         templ_sum2 = templ_norm +\r
 449                      CV_SQR(templ_mean.val[0]) + CV_SQR(templ_mean.val[1]) +\r
 450                      CV_SQR(templ_mean.val[2]) + CV_SQR(templ_mean.val[3]);\r
 451 \r
 452         if( num_type != 1 )\r
 453         {\r
 454             templ_mean = cvScalarAll(0);\r
 455             templ_norm = templ_sum2;\r
 456         }\r
 457         \r
 458         templ_sum2 /= inv_area;\r
 459         templ_norm = sqrt(templ_norm);\r
 460         templ_norm /= sqrt(inv_area); // care of accuracy here\r
 461 \r
 462         q0 = (double*)sqsum->data.ptr;\r
 463         q1 = q0 + templ->cols*cn;\r
 464         q2 = (double*)(sqsum->data.ptr + templ->rows*sqsum->step);\r
 465         q3 = q2 + templ->cols*cn;\r
 466     }\r
 467 \r
 468     p0 = (double*)sum->data.ptr;\r
 469     p1 = p0 + templ->cols*cn;\r
 470     p2 = (double*)(sum->data.ptr + templ->rows*sum->step);\r
 471     p3 = p2 + templ->cols*cn;\r
 472 \r
 473     sum_step = sum ? sum->step / sizeof(double) : 0;\r
 474     sqsum_step = sqsum ? sqsum->step / sizeof(double) : 0;\r
 475 \r
 476     for( i = 0; i < result->rows; i++ )\r
 477     {\r
 478         float* rrow = (float*)(result->data.ptr + i*result->step);\r
 479         idx = i * sum_step;\r
 480         idx2 = i * sqsum_step;\r
 481 \r
 482         for( j = 0; j < result->cols; j++, idx += cn, idx2 += cn )\r
 483         {\r
 484             double num = rrow[j], t;\r
 485             double wnd_mean2 = 0, wnd_sum2 = 0;\r
 486             \r
 487             if( num_type == 1 )\r
 488             {\r
 489                 for( k = 0; k < cn; k++ )\r
 490                 {\r
 491                     t = p0[idx+k] - p1[idx+k] - p2[idx+k] + p3[idx+k];\r
 492                     wnd_mean2 += CV_SQR(t);\r
 493                     num -= t*templ_mean.val[k];\r
 494                 }\r
 495 \r
 496                 wnd_mean2 *= inv_area;\r
 497             }\r
 498 \r
 499             if( is_normed || num_type == 2 )\r
 500             {\r
 501                 for( k = 0; k < cn; k++ )\r
 502                 {\r
 503                     t = q0[idx2+k] - q1[idx2+k] - q2[idx2+k] + q3[idx2+k];\r
 504                     wnd_sum2 += t;\r
 505                 }\r
 506 \r
 507                 if( num_type == 2 )\r
 508                     num = wnd_sum2 - 2*num + templ_sum2;\r
 509             }\r
 510 \r
 511             if( is_normed )\r
 512             {\r
 513                 t = sqrt(MAX(wnd_sum2 - wnd_mean2,0))*templ_norm;\r
 514                 if( t > DBL_EPSILON )\r
 515                 {\r
 516                     num /= t;\r
 517                     if( fabs(num) > 1. )\r
 518                         num = num > 0 ? 1 : -1;\r
 519                 }\r
 520                 else\r
 521                     num = method != CV_TM_SQDIFF_NORMED || num < DBL_EPSILON ? 0 : 1;\r
 522             }\r
 523 \r
 524             rrow[j] = (float)num;\r
 525         }\r
 526     }\r
 527         \r
 528     __END__;\r
 529 \r
 530     cvReleaseMat( &sum );\r
 531     cvReleaseMat( &sqsum );\r
 532 }\r
 533 \r
 534 void cv::matchTemplate( const Mat& image, const Mat& templ, Mat& result, int method )\r
 535 {\r
 536     result.create( std::abs(image.rows - templ.rows) + 1,\r
 537                    std::abs(image.cols - templ.cols) + 1, CV_32F );\r
 538     CvMat _image = image, _templ = templ, _result = result;\r
 539     cvMatchTemplate( &_image, &_templ, &_result, method );    \r
 540 }\r
 541 \r
 542 /* End of file. */\r