Merge branch 'master' into python2and3
[profile/ivi/opencv.git] / modules / imgproc / src / corner.cpp
1 /*M///////////////////////////////////////////////////////////////////////////////////////
2 //
3 //  IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
4 //
5 //  By downloading, copying, installing or using the software you agree to this license.
6 //  If you do not agree to this license, do not download, install,
7 //  copy or use the software.
8 //
9 //
10 //                           License Agreement
11 //                For Open Source Computer Vision Library
12 //
13 // Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
14 // Copyright (C) 2009, Willow Garage Inc., all rights reserved.
15 // Third party copyrights are property of their respective owners.
16 //
17 // Redistribution and use in source and binary forms, with or without modification,
18 // are permitted provided that the following conditions are met:
19 //
20 //   * Redistribution's of source code must retain the above copyright notice,
21 //     this list of conditions and the following disclaimer.
22 //
23 //   * Redistribution's in binary form must reproduce the above copyright notice,
24 //     this list of conditions and the following disclaimer in the documentation
25 //     and/or other materials provided with the distribution.
26 //
27 //   * The name of the copyright holders may not be used to endorse or promote products
28 //     derived from this software without specific prior written permission.
29 //
30 // This software is provided by the copyright holders and contributors "as is" and
31 // any express or implied warranties, including, but not limited to, the implied
32 // warranties of merchantability and fitness for a particular purpose are disclaimed.
33 // In no event shall the Intel Corporation or contributors be liable for any direct,
34 // indirect, incidental, special, exemplary, or consequential damages
35 // (including, but not limited to, procurement of substitute goods or services;
36 // loss of use, data, or profits; or business interruption) however caused
37 // and on any theory of liability, whether in contract, strict liability,
38 // or tort (including negligence or otherwise) arising in any way out of
39 // the use of this software, even if advised of the possibility of such damage.
40 //
41 //M*/
42
43 #include "precomp.hpp"
44 #include "opencl_kernels.hpp"
45
46 namespace cv
47 {
48
49 static void calcMinEigenVal( const Mat& _cov, Mat& _dst )
50 {
51     int i, j;
52     Size size = _cov.size();
53 #if CV_SSE
54     volatile bool simd = checkHardwareSupport(CV_CPU_SSE);
55 #endif
56
57     if( _cov.isContinuous() && _dst.isContinuous() )
58     {
59         size.width *= size.height;
60         size.height = 1;
61     }
62
63     for( i = 0; i < size.height; i++ )
64     {
65         const float* cov = (const float*)(_cov.data + _cov.step*i);
66         float* dst = (float*)(_dst.data + _dst.step*i);
67         j = 0;
68     #if CV_SSE
69         if( simd )
70         {
71             __m128 half = _mm_set1_ps(0.5f);
72             for( ; j <= size.width - 5; j += 4 )
73             {
74                 __m128 t0 = _mm_loadu_ps(cov + j*3); // a0 b0 c0 x
75                 __m128 t1 = _mm_loadu_ps(cov + j*3 + 3); // a1 b1 c1 x
76                 __m128 t2 = _mm_loadu_ps(cov + j*3 + 6); // a2 b2 c2 x
77                 __m128 t3 = _mm_loadu_ps(cov + j*3 + 9); // a3 b3 c3 x
78                 __m128 a, b, c, t;
79                 t = _mm_unpacklo_ps(t0, t1); // a0 a1 b0 b1
80                 c = _mm_unpackhi_ps(t0, t1); // c0 c1 x x
81                 b = _mm_unpacklo_ps(t2, t3); // a2 a3 b2 b3
82                 c = _mm_movelh_ps(c, _mm_unpackhi_ps(t2, t3)); // c0 c1 c2 c3
83                 a = _mm_movelh_ps(t, b);
84                 b = _mm_movehl_ps(b, t);
85                 a = _mm_mul_ps(a, half);
86                 c = _mm_mul_ps(c, half);
87                 t = _mm_sub_ps(a, c);
88                 t = _mm_add_ps(_mm_mul_ps(t, t), _mm_mul_ps(b,b));
89                 a = _mm_sub_ps(_mm_add_ps(a, c), _mm_sqrt_ps(t));
90                 _mm_storeu_ps(dst + j, a);
91             }
92         }
93     #endif
94         for( ; j < size.width; j++ )
95         {
96             float a = cov[j*3]*0.5f;
97             float b = cov[j*3+1];
98             float c = cov[j*3+2]*0.5f;
99             dst[j] = (float)((a + c) - std::sqrt((a - c)*(a - c) + b*b));
100         }
101     }
102 }
103
104
105 static void calcHarris( const Mat& _cov, Mat& _dst, double k )
106 {
107     int i, j;
108     Size size = _cov.size();
109 #if CV_SSE
110     volatile bool simd = checkHardwareSupport(CV_CPU_SSE);
111 #endif
112
113     if( _cov.isContinuous() && _dst.isContinuous() )
114     {
115         size.width *= size.height;
116         size.height = 1;
117     }
118
119     for( i = 0; i < size.height; i++ )
120     {
121         const float* cov = (const float*)(_cov.data + _cov.step*i);
122         float* dst = (float*)(_dst.data + _dst.step*i);
123         j = 0;
124
125     #if CV_SSE
126         if( simd )
127         {
128             __m128 k4 = _mm_set1_ps((float)k);
129             for( ; j <= size.width - 5; j += 4 )
130             {
131                 __m128 t0 = _mm_loadu_ps(cov + j*3); // a0 b0 c0 x
132                 __m128 t1 = _mm_loadu_ps(cov + j*3 + 3); // a1 b1 c1 x
133                 __m128 t2 = _mm_loadu_ps(cov + j*3 + 6); // a2 b2 c2 x
134                 __m128 t3 = _mm_loadu_ps(cov + j*3 + 9); // a3 b3 c3 x
135                 __m128 a, b, c, t;
136                 t = _mm_unpacklo_ps(t0, t1); // a0 a1 b0 b1
137                 c = _mm_unpackhi_ps(t0, t1); // c0 c1 x x
138                 b = _mm_unpacklo_ps(t2, t3); // a2 a3 b2 b3
139                 c = _mm_movelh_ps(c, _mm_unpackhi_ps(t2, t3)); // c0 c1 c2 c3
140                 a = _mm_movelh_ps(t, b);
141                 b = _mm_movehl_ps(b, t);
142                 t = _mm_add_ps(a, c);
143                 a = _mm_sub_ps(_mm_mul_ps(a, c), _mm_mul_ps(b, b));
144                 t = _mm_mul_ps(_mm_mul_ps(k4, t), t);
145                 a = _mm_sub_ps(a, t);
146                 _mm_storeu_ps(dst + j, a);
147             }
148         }
149     #endif
150
151         for( ; j < size.width; j++ )
152         {
153             float a = cov[j*3];
154             float b = cov[j*3+1];
155             float c = cov[j*3+2];
156             dst[j] = (float)(a*c - b*b - k*(a + c)*(a + c));
157         }
158     }
159 }
160
161
162 static void eigen2x2( const float* cov, float* dst, int n )
163 {
164     for( int j = 0; j < n; j++ )
165     {
166         double a = cov[j*3];
167         double b = cov[j*3+1];
168         double c = cov[j*3+2];
169
170         double u = (a + c)*0.5;
171         double v = std::sqrt((a - c)*(a - c)*0.25 + b*b);
172         double l1 = u + v;
173         double l2 = u - v;
174
175         double x = b;
176         double y = l1 - a;
177         double e = fabs(x);
178
179         if( e + fabs(y) < 1e-4 )
180         {
181             y = b;
182             x = l1 - c;
183             e = fabs(x);
184             if( e + fabs(y) < 1e-4 )
185             {
186                 e = 1./(e + fabs(y) + FLT_EPSILON);
187                 x *= e, y *= e;
188             }
189         }
190
191         double d = 1./std::sqrt(x*x + y*y + DBL_EPSILON);
192         dst[6*j] = (float)l1;
193         dst[6*j + 2] = (float)(x*d);
194         dst[6*j + 3] = (float)(y*d);
195
196         x = b;
197         y = l2 - a;
198         e = fabs(x);
199
200         if( e + fabs(y) < 1e-4 )
201         {
202             y = b;
203             x = l2 - c;
204             e = fabs(x);
205             if( e + fabs(y) < 1e-4 )
206             {
207                 e = 1./(e + fabs(y) + FLT_EPSILON);
208                 x *= e, y *= e;
209             }
210         }
211
212         d = 1./std::sqrt(x*x + y*y + DBL_EPSILON);
213         dst[6*j + 1] = (float)l2;
214         dst[6*j + 4] = (float)(x*d);
215         dst[6*j + 5] = (float)(y*d);
216     }
217 }
218
219 static void calcEigenValsVecs( const Mat& _cov, Mat& _dst )
220 {
221     Size size = _cov.size();
222     if( _cov.isContinuous() && _dst.isContinuous() )
223     {
224         size.width *= size.height;
225         size.height = 1;
226     }
227
228     for( int i = 0; i < size.height; i++ )
229     {
230         const float* cov = (const float*)(_cov.data + _cov.step*i);
231         float* dst = (float*)(_dst.data + _dst.step*i);
232
233         eigen2x2(cov, dst, size.width);
234     }
235 }
236
237
238 enum { MINEIGENVAL=0, HARRIS=1, EIGENVALSVECS=2 };
239
240
241 static void
242 cornerEigenValsVecs( const Mat& src, Mat& eigenv, int block_size,
243                      int aperture_size, int op_type, double k=0.,
244                      int borderType=BORDER_DEFAULT )
245 {
246 #ifdef HAVE_TEGRA_OPTIMIZATION
247     if (tegra::cornerEigenValsVecs(src, eigenv, block_size, aperture_size, op_type, k, borderType))
248         return;
249 #endif
250
251     int depth = src.depth();
252     double scale = (double)(1 << ((aperture_size > 0 ? aperture_size : 3) - 1)) * block_size;
253     if( aperture_size < 0 )
254         scale *= 2.0;
255     if( depth == CV_8U )
256         scale *= 255.0;
257     scale = 1.0/scale;
258
259     CV_Assert( src.type() == CV_8UC1 || src.type() == CV_32FC1 );
260
261     Mat Dx, Dy;
262     if( aperture_size > 0 )
263     {
264         Sobel( src, Dx, CV_32F, 1, 0, aperture_size, scale, 0, borderType );
265         Sobel( src, Dy, CV_32F, 0, 1, aperture_size, scale, 0, borderType );
266     }
267     else
268     {
269         Scharr( src, Dx, CV_32F, 1, 0, scale, 0, borderType );
270         Scharr( src, Dy, CV_32F, 0, 1, scale, 0, borderType );
271     }
272
273     Size size = src.size();
274     Mat cov( size, CV_32FC3 );
275     int i, j;
276
277     for( i = 0; i < size.height; i++ )
278     {
279         float* cov_data = (float*)(cov.data + i*cov.step);
280         const float* dxdata = (const float*)(Dx.data + i*Dx.step);
281         const float* dydata = (const float*)(Dy.data + i*Dy.step);
282
283         for( j = 0; j < size.width; j++ )
284         {
285             float dx = dxdata[j];
286             float dy = dydata[j];
287
288             cov_data[j*3] = dx*dx;
289             cov_data[j*3+1] = dx*dy;
290             cov_data[j*3+2] = dy*dy;
291         }
292     }
293
294     boxFilter(cov, cov, cov.depth(), Size(block_size, block_size),
295         Point(-1,-1), false, borderType );
296
297     if( op_type == MINEIGENVAL )
298         calcMinEigenVal( cov, eigenv );
299     else if( op_type == HARRIS )
300         calcHarris( cov, eigenv, k );
301     else if( op_type == EIGENVALSVECS )
302         calcEigenValsVecs( cov, eigenv );
303 }
304
305 #ifdef HAVE_OPENCL
306
307 static bool extractCovData(InputArray _src, UMat & Dx, UMat & Dy, int depth,
308                            float scale, int aperture_size, int borderType)
309 {
310     UMat src = _src.getUMat();
311
312     Size wholeSize;
313     Point ofs;
314     src.locateROI(wholeSize, ofs);
315
316     const int sobel_lsz = 16;
317     if ((aperture_size == 3 || aperture_size == 5 || aperture_size == 7 || aperture_size == -1) &&
318         wholeSize.height > sobel_lsz + (aperture_size >> 1) &&
319         wholeSize.width > sobel_lsz + (aperture_size >> 1))
320     {
321         CV_Assert(depth == CV_8U || depth == CV_32F);
322
323         Dx.create(src.size(), CV_32FC1);
324         Dy.create(src.size(), CV_32FC1);
325
326         size_t localsize[2] = { sobel_lsz, sobel_lsz };
327         size_t globalsize[2] = { localsize[0] * (1 + (src.cols - 1) / localsize[0]),
328                                  localsize[1] * (1 + (src.rows - 1) / localsize[1]) };
329
330         int src_offset_x = (int)((src.offset % src.step) / src.elemSize());
331         int src_offset_y = (int)(src.offset / src.step);
332
333         const char * const borderTypes[] = { "BORDER_CONSTANT", "BORDER_REPLICATE", "BORDER_REFLECT",
334                                              "BORDER_WRAP", "BORDER_REFLECT101" };
335
336         ocl::Kernel k(format("sobel%d", aperture_size).c_str(), ocl::imgproc::covardata_oclsrc,
337                       cv::format("-D BLK_X=%d -D BLK_Y=%d -D %s -D SRCTYPE=%s%s",
338                                  (int)localsize[0], (int)localsize[1], borderTypes[borderType], ocl::typeToStr(depth),
339                                  aperture_size < 0 ? " -D SCHARR" : ""));
340         if (k.empty())
341             return false;
342
343         k.args(ocl::KernelArg::PtrReadOnly(src), (int)src.step, src_offset_x, src_offset_y,
344                ocl::KernelArg::WriteOnlyNoSize(Dx), ocl::KernelArg::WriteOnly(Dy),
345                wholeSize.height, wholeSize.width, scale);
346
347         return k.run(2, globalsize, localsize, false);
348     }
349     else
350     {
351         if (aperture_size > 0)
352         {
353             Sobel(_src, Dx, CV_32F, 1, 0, aperture_size, scale, 0, borderType);
354             Sobel(_src, Dy, CV_32F, 0, 1, aperture_size, scale, 0, borderType);
355         }
356         else
357         {
358             Scharr(_src, Dx, CV_32F, 1, 0, scale, 0, borderType);
359             Scharr(_src, Dy, CV_32F, 0, 1, scale, 0, borderType);
360         }
361     }
362
363     return true;
364 }
365
366 static bool ocl_cornerMinEigenValVecs(InputArray _src, OutputArray _dst, int block_size,
367                                       int aperture_size, double k, int borderType, int op_type)
368 {
369     CV_Assert(op_type == HARRIS || op_type == MINEIGENVAL);
370
371     if ( !(borderType == BORDER_CONSTANT || borderType == BORDER_REPLICATE ||
372            borderType == BORDER_REFLECT || borderType == BORDER_REFLECT_101) )
373         return false;
374
375     int type = _src.type(), depth = CV_MAT_DEPTH(type);
376     if ( !(type == CV_8UC1 || type == CV_32FC1) )
377         return false;
378
379     const char * const borderTypes[] = { "BORDER_CONSTANT", "BORDER_REPLICATE", "BORDER_REFLECT",
380                                          "BORDER_WRAP", "BORDER_REFLECT101" };
381     const char * const cornerType[] = { "CORNER_MINEIGENVAL", "CORNER_HARRIS", 0 };
382
383
384     double scale = (double)(1 << ((aperture_size > 0 ? aperture_size : 3) - 1)) * block_size;
385     if (aperture_size < 0)
386         scale *= 2.0;
387     if (depth == CV_8U)
388         scale *= 255.0;
389     scale = 1.0 / scale;
390
391     UMat Dx, Dy;
392     if (!extractCovData(_src, Dx, Dy, depth, (float)scale, aperture_size, borderType))
393         return false;
394
395     ocl::Kernel cornelKernel("corner", ocl::imgproc::corner_oclsrc,
396                              format("-D anX=%d -D anY=%d -D ksX=%d -D ksY=%d -D %s -D %s",
397                                     block_size / 2, block_size / 2, block_size, block_size,
398                                     borderTypes[borderType], cornerType[op_type]));
399     if (cornelKernel.empty())
400         return false;
401
402     _dst.createSameSize(_src, CV_32FC1);
403     UMat dst = _dst.getUMat();
404
405     cornelKernel.args(ocl::KernelArg::ReadOnly(Dx), ocl::KernelArg::ReadOnly(Dy),
406                       ocl::KernelArg::WriteOnly(dst), (float)k);
407
408     size_t blockSizeX = 256, blockSizeY = 1;
409     size_t gSize = blockSizeX - block_size / 2 * 2;
410     size_t globalSizeX = (Dx.cols) % gSize == 0 ? Dx.cols / gSize * blockSizeX : (Dx.cols / gSize + 1) * blockSizeX;
411     size_t rows_per_thread = 2;
412     size_t globalSizeY = ((Dx.rows + rows_per_thread - 1) / rows_per_thread) % blockSizeY == 0 ?
413                          ((Dx.rows + rows_per_thread - 1) / rows_per_thread) :
414                          (((Dx.rows + rows_per_thread - 1) / rows_per_thread) / blockSizeY + 1) * blockSizeY;
415
416     size_t globalsize[2] = { globalSizeX, globalSizeY }, localsize[2] = { blockSizeX, blockSizeY };
417     return cornelKernel.run(2, globalsize, localsize, false);
418 }
419
420 static bool ocl_preCornerDetect( InputArray _src, OutputArray _dst, int ksize, int borderType, int depth )
421 {
422     UMat Dx, Dy, D2x, D2y, Dxy;
423
424     if (!extractCovData(_src, Dx, Dy, depth, 1, ksize, borderType))
425         return false;
426
427     Sobel( _src, D2x, CV_32F, 2, 0, ksize, 1, 0, borderType );
428     Sobel( _src, D2y, CV_32F, 0, 2, ksize, 1, 0, borderType );
429     Sobel( _src, Dxy, CV_32F, 1, 1, ksize, 1, 0, borderType );
430
431     _dst.create( _src.size(), CV_32FC1 );
432     UMat dst = _dst.getUMat();
433
434     double factor = 1 << (ksize - 1);
435     if( depth == CV_8U )
436         factor *= 255;
437     factor = 1./(factor * factor * factor);
438
439     ocl::Kernel k("preCornerDetect", ocl::imgproc::precornerdetect_oclsrc);
440     if (k.empty())
441         return false;
442
443     k.args(ocl::KernelArg::ReadOnlyNoSize(Dx), ocl::KernelArg::ReadOnlyNoSize(Dy),
444            ocl::KernelArg::ReadOnlyNoSize(D2x), ocl::KernelArg::ReadOnlyNoSize(D2y),
445            ocl::KernelArg::ReadOnlyNoSize(Dxy), ocl::KernelArg::WriteOnly(dst), (float)factor);
446
447     size_t globalsize[2] = { dst.cols, dst.rows };
448     return k.run(2, globalsize, NULL, false);
449 }
450
451 #endif
452
453 }
454
455 void cv::cornerMinEigenVal( InputArray _src, OutputArray _dst, int blockSize, int ksize, int borderType )
456 {
457     CV_OCL_RUN(_src.dims() <= 2 && _dst.isUMat(),
458                ocl_cornerMinEigenValVecs(_src, _dst, blockSize, ksize, 0.0, borderType, MINEIGENVAL))
459
460     Mat src = _src.getMat();
461     _dst.create( src.size(), CV_32FC1 );
462     Mat dst = _dst.getMat();
463 #if defined(HAVE_IPP) && (IPP_VERSION_MAJOR >= 8)
464     typedef IppStatus (CV_STDCALL * ippiMinEigenValGetBufferSize)(IppiSize, int, int, int*);
465     typedef IppStatus (CV_STDCALL * ippiMinEigenVal)(const void*, int, Ipp32f*, int, IppiSize, IppiKernelType, int, int, Ipp8u*);
466     IppiKernelType kerType;
467     int kerSize = ksize;
468     if (ksize < 0)
469     {
470         kerType = ippKernelScharr;
471         kerSize = 3;
472     } else
473     {
474         kerType = ippKernelSobel;
475     }
476     bool isolated = (borderType & BORDER_ISOLATED) != 0;
477     int borderTypeNI = borderType & ~BORDER_ISOLATED;
478     if ((borderTypeNI == BORDER_REPLICATE && (!src.isSubmatrix() || isolated)) &&
479         (kerSize == 3 || kerSize == 5) && (blockSize == 3 || blockSize == 5))
480     {
481         ippiMinEigenValGetBufferSize getBufferSizeFunc = 0;
482         ippiMinEigenVal minEigenValFunc = 0;
483         float norm_coef = 0.f;
484
485         if (src.type() == CV_8UC1)
486         {
487             getBufferSizeFunc = (ippiMinEigenValGetBufferSize) ippiMinEigenValGetBufferSize_8u32f_C1R;
488             minEigenValFunc = (ippiMinEigenVal) ippiMinEigenVal_8u32f_C1R;
489             norm_coef = 1.f / 255.f;
490         } else if (src.type() == CV_32FC1)
491         {
492             getBufferSizeFunc = (ippiMinEigenValGetBufferSize) ippiMinEigenValGetBufferSize_32f_C1R;
493             minEigenValFunc = (ippiMinEigenVal) ippiMinEigenVal_32f_C1R;
494             norm_coef = 255.f;
495         }
496         norm_coef = kerType == ippKernelSobel ? norm_coef : norm_coef / 2.45f;
497
498         if (getBufferSizeFunc && minEigenValFunc)
499         {
500             int bufferSize;
501             IppiSize srcRoi = { src.cols, src.rows };
502             IppStatus ok = getBufferSizeFunc(srcRoi, kerSize, blockSize, &bufferSize);
503             if (ok >= 0)
504             {
505                 AutoBuffer<uchar> buffer(bufferSize);
506                 ok = minEigenValFunc(src.data, (int) src.step, (Ipp32f*) dst.data, (int) dst.step, srcRoi, kerType, kerSize, blockSize, buffer);
507                 CV_SUPPRESS_DEPRECATED_START
508                 if (ok >= 0) ok = ippiMulC_32f_C1IR(norm_coef, (Ipp32f*) dst.data, (int) dst.step, srcRoi);
509                 CV_SUPPRESS_DEPRECATED_END
510                 if (ok >= 0)
511                     return;
512             }
513             setIppErrorStatus();
514         }
515     }
516 #endif
517     cornerEigenValsVecs( src, dst, blockSize, ksize, MINEIGENVAL, 0, borderType );
518 }
519
520 void cv::cornerHarris( InputArray _src, OutputArray _dst, int blockSize, int ksize, double k, int borderType )
521 {
522     CV_OCL_RUN(_src.dims() <= 2 && _dst.isUMat(),
523                ocl_cornerMinEigenValVecs(_src, _dst, blockSize, ksize, k, borderType, HARRIS))
524
525     Mat src = _src.getMat();
526     _dst.create( src.size(), CV_32FC1 );
527     Mat dst = _dst.getMat();
528
529 #if IPP_VERSION_X100 >= 801 && 0
530     int type = src.type(), depth = CV_MAT_DEPTH(type), cn = CV_MAT_CN(type);
531     int borderTypeNI = borderType & ~BORDER_ISOLATED;
532     bool isolated = (borderType & BORDER_ISOLATED) != 0;
533
534     if ( (ksize == 3 || ksize == 5) && (type == CV_8UC1 || type == CV_32FC1) &&
535         (borderTypeNI == BORDER_CONSTANT || borderTypeNI == BORDER_REPLICATE) && cn == 1 && (!src.isSubmatrix() || isolated) )
536     {
537         IppiSize roisize = { src.cols, src.rows };
538         IppiMaskSize masksize = ksize == 5 ? ippMskSize5x5 : ippMskSize3x3;
539         IppDataType datatype = type == CV_8UC1 ? ipp8u : ipp32f;
540         Ipp32s bufsize = 0;
541
542         double scale = (double)(1 << ((ksize > 0 ? ksize : 3) - 1)) * blockSize;
543         if (ksize < 0)
544             scale *= 2.0;
545         if (depth == CV_8U)
546             scale *= 255.0;
547         scale = std::pow(scale, -4.0);
548
549         if (ippiHarrisCornerGetBufferSize(roisize, masksize, blockSize, datatype, cn, &bufsize) >= 0)
550         {
551             Ipp8u * buffer = ippsMalloc_8u(bufsize);
552             IppiDifferentialKernel filterType = ksize > 0 ? ippFilterSobel : ippFilterScharr;
553             IppiBorderType borderTypeIpp = borderTypeNI == BORDER_CONSTANT ? ippBorderConst : ippBorderRepl;
554             IppStatus status = (IppStatus)-1;
555
556             if (depth == CV_8U)
557                 status = ippiHarrisCorner_8u32f_C1R((const Ipp8u *)src.data, (int)src.step, (Ipp32f *)dst.data, (int)dst.step, roisize,
558                                                     filterType, masksize, blockSize, (Ipp32f)k, (Ipp32f)scale, borderTypeIpp, 0, buffer);
559             else if (depth == CV_32F)
560                 status = ippiHarrisCorner_32f_C1R((const Ipp32f *)src.data, (int)src.step, (Ipp32f *)dst.data, (int)dst.step, roisize,
561                                                   filterType, masksize, blockSize, (Ipp32f)k, (Ipp32f)scale, borderTypeIpp, 0, buffer);
562             ippsFree(buffer);
563
564             if (status >= 0)
565                 return;
566         }
567         setIppErrorStatus();
568     }
569 #endif
570
571     cornerEigenValsVecs( src, dst, blockSize, ksize, HARRIS, k, borderType );
572 }
573
574
575 void cv::cornerEigenValsAndVecs( InputArray _src, OutputArray _dst, int blockSize, int ksize, int borderType )
576 {
577     Mat src = _src.getMat();
578     Size dsz = _dst.size();
579     int dtype = _dst.type();
580
581     if( dsz.height != src.rows || dsz.width*CV_MAT_CN(dtype) != src.cols*6 || CV_MAT_DEPTH(dtype) != CV_32F )
582         _dst.create( src.size(), CV_32FC(6) );
583     Mat dst = _dst.getMat();
584     cornerEigenValsVecs( src, dst, blockSize, ksize, EIGENVALSVECS, 0, borderType );
585 }
586
587
588 void cv::preCornerDetect( InputArray _src, OutputArray _dst, int ksize, int borderType )
589 {
590     int type = _src.type();
591     CV_Assert( type == CV_8UC1 || type == CV_32FC1 );
592
593     CV_OCL_RUN( _src.dims() <= 2 && _dst.isUMat(),
594                 ocl_preCornerDetect(_src, _dst, ksize, borderType, CV_MAT_DEPTH(type)))
595
596     Mat Dx, Dy, D2x, D2y, Dxy, src = _src.getMat();
597     _dst.create( src.size(), CV_32FC1 );
598     Mat dst = _dst.getMat();
599
600     Sobel( src, Dx, CV_32F, 1, 0, ksize, 1, 0, borderType );
601     Sobel( src, Dy, CV_32F, 0, 1, ksize, 1, 0, borderType );
602     Sobel( src, D2x, CV_32F, 2, 0, ksize, 1, 0, borderType );
603     Sobel( src, D2y, CV_32F, 0, 2, ksize, 1, 0, borderType );
604     Sobel( src, Dxy, CV_32F, 1, 1, ksize, 1, 0, borderType );
605
606     double factor = 1 << (ksize - 1);
607     if( src.depth() == CV_8U )
608         factor *= 255;
609     factor = 1./(factor * factor * factor);
610
611 #if CV_SSE2
612     volatile bool haveSSE2 = cv::checkHardwareSupport(CV_CPU_SSE2);
613     __m128 v_factor = _mm_set1_ps((float)factor), v_m2 = _mm_set1_ps(-2.0f);
614 #endif
615
616     Size size = src.size();
617     int i, j;
618     for( i = 0; i < size.height; i++ )
619     {
620         float* dstdata = (float*)(dst.data + i*dst.step);
621         const float* dxdata = (const float*)(Dx.data + i*Dx.step);
622         const float* dydata = (const float*)(Dy.data + i*Dy.step);
623         const float* d2xdata = (const float*)(D2x.data + i*D2x.step);
624         const float* d2ydata = (const float*)(D2y.data + i*D2y.step);
625         const float* dxydata = (const float*)(Dxy.data + i*Dxy.step);
626
627         j = 0;
628
629 #if CV_SSE2
630         if (haveSSE2)
631         {
632             for( ; j <= size.width - 4; j += 4 )
633             {
634                 __m128 v_dx = _mm_loadu_ps((const float *)(dxdata + j));
635                 __m128 v_dy = _mm_loadu_ps((const float *)(dydata + j));
636
637                 __m128 v_s1 = _mm_mul_ps(_mm_mul_ps(v_dx, v_dx), _mm_loadu_ps((const float *)(d2ydata + j)));
638                 __m128 v_s2 = _mm_mul_ps(_mm_mul_ps(v_dy, v_dy), _mm_loadu_ps((const float *)(d2xdata + j)));
639                 __m128 v_s3 = _mm_mul_ps(_mm_mul_ps(v_dx, v_dy), _mm_loadu_ps((const float *)(dxydata + j)));
640                 v_s1 = _mm_mul_ps(v_factor, _mm_add_ps(v_s1, _mm_add_ps(v_s2, _mm_mul_ps(v_s3, v_m2))));
641                 _mm_storeu_ps(dstdata + j, v_s1);
642             }
643         }
644 #endif
645
646         for( ; j < size.width; j++ )
647         {
648             float dx = dxdata[j];
649             float dy = dydata[j];
650             dstdata[j] = (float)(factor*(dx*dx*d2ydata[j] + dy*dy*d2xdata[j] - 2*dx*dy*dxydata[j]));
651         }
652     }
653 }
654
655 CV_IMPL void
656 cvCornerMinEigenVal( const CvArr* srcarr, CvArr* dstarr,
657                      int block_size, int aperture_size )
658 {
659     cv::Mat src = cv::cvarrToMat(srcarr), dst = cv::cvarrToMat(dstarr);
660
661     CV_Assert( src.size() == dst.size() && dst.type() == CV_32FC1 );
662     cv::cornerMinEigenVal( src, dst, block_size, aperture_size, cv::BORDER_REPLICATE );
663 }
664
665 CV_IMPL void
666 cvCornerHarris( const CvArr* srcarr, CvArr* dstarr,
667                 int block_size, int aperture_size, double k )
668 {
669     cv::Mat src = cv::cvarrToMat(srcarr), dst = cv::cvarrToMat(dstarr);
670
671     CV_Assert( src.size() == dst.size() && dst.type() == CV_32FC1 );
672     cv::cornerHarris( src, dst, block_size, aperture_size, k, cv::BORDER_REPLICATE );
673 }
674
675
676 CV_IMPL void
677 cvCornerEigenValsAndVecs( const void* srcarr, void* dstarr,
678                           int block_size, int aperture_size )
679 {
680     cv::Mat src = cv::cvarrToMat(srcarr), dst = cv::cvarrToMat(dstarr);
681
682     CV_Assert( src.rows == dst.rows && src.cols*6 == dst.cols*dst.channels() && dst.depth() == CV_32F );
683     cv::cornerEigenValsAndVecs( src, dst, block_size, aperture_size, cv::BORDER_REPLICATE );
684 }
685
686
687 CV_IMPL void
688 cvPreCornerDetect( const void* srcarr, void* dstarr, int aperture_size )
689 {
690     cv::Mat src = cv::cvarrToMat(srcarr), dst = cv::cvarrToMat(dstarr);
691
692     CV_Assert( src.size() == dst.size() && dst.type() == CV_32FC1 );
693     cv::preCornerDetect( src, dst, aperture_size, cv::BORDER_REPLICATE );
694 }
695
696 /* End of file */