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 // Intel License Agreement
11 // For Open Source Computer Vision Library
12 //
13 // Copyright (C) 2000, Intel Corporation, all rights reserved.
14 // Third party copyrights are property of their respective owners.
15 //
16 // Redistribution and use in source and binary forms, with or without modification,
17 // are permitted provided that the following conditions are met:
18 //
19 // * Redistribution's of source code must retain the above copyright notice,
20 // this list of conditions and the following disclaimer.
21 //
22 // * Redistribution's in binary form must reproduce the above copyright notice,
23 // this list of conditions and the following disclaimer in the documentation
24 // and/or other materials provided with the distribution.
25 //
26 // * The name of Intel Corporation may not be used to endorse or promote products
27 // derived from this software without specific prior written permission.
28 //
29 // This software is provided by the copyright holders and contributors "as is" and
30 // any express or implied warranties, including, but not limited to, the implied
31 // warranties of merchantability and fitness for a particular purpose are disclaimed.
32 // In no event shall the Intel Corporation or contributors be liable for any direct,
33 // indirect, incidental, special, exemplary, or consequential damages
34 // (including, but not limited to, procurement of substitute goods or services;
35 // loss of use, data, or profits; or business interruption) however caused
36 // and on any theory of liability, whether in contract, strict liability,
37 // or tort (including negligence or otherwise) arising in any way out of
38 // the use of this software, even if advised of the possibility of such damage.
39 //
40 //M*/
41
42 #include "_cv.h"
43
44 #define ICV_DEF_INTEGRAL_OP_C1( flavor, arrtype, sumtype, sqsumtype, worktype, \
45 cast_macro, cast_sqr_macro ) \
46 static CvStatus CV_STDCALL \
47 icvIntegralImage_##flavor##_C1R( const arrtype* src, int srcstep,\
48 sumtype* sum, int sumstep, \
49 sqsumtype* sqsum, int sqsumstep,\
50 sumtype* tilted, int tiltedstep,\
51 CvSize size ) \
52 { \
53 int x, y; \
54 sumtype s; \
55 sqsumtype sq; \
56 sumtype* buf = 0; \
57 \
58 srcstep /= sizeof(src[0]); \
59 \
60 memset( sum, 0, (size.width+1)*sizeof(sum[0])); \
61 sumstep /= sizeof(sum[0]); \
62 sum += sumstep + 1; \
63 \
64 if( sqsum ) \
65 { \
66 memset( sqsum, 0, (size.width+1)*sizeof(sqsum[0])); \
67 sqsumstep /= sizeof(sqsum[0]); \
68 sqsum += sqsumstep + 1; \
69 } \
70 \
71 if( tilted ) \
72 { \
73 memset( tilted, 0, (size.width+1)*sizeof(tilted[0])); \
74 tiltedstep /= sizeof(tilted[0]); \
75 tilted += tiltedstep + 1; \
76 } \
77 \
78 if( sqsum == 0 && tilted == 0 ) \
79 { \
80 for( y = 0; y < size.height; y++, src += srcstep, \
81 sum += sumstep ) \
82 { \
83 sum[-1] = 0; \
84 for( x = 0, s = 0; x < size.width; x++ ) \
85 { \
86 sumtype t = cast_macro(src[x]); \
87 s += t; \
88 sum[x] = sum[x - sumstep] + s; \
89 } \
90 } \
91 } \
92 else if( tilted == 0 ) \
93 { \
94 for( y = 0; y < size.height; y++, src += srcstep, \
95 sum += sumstep, sqsum += sqsumstep ) \
96 { \
97 sum[-1] = 0; \
98 sqsum[-1] = 0; \
99 \
100 for( x = 0, s = 0, sq = 0; x < size.width; x++ ) \
101 { \
102 worktype it = src[x]; \
103 sumtype t = cast_macro(it); \
104 sqsumtype tq = cast_sqr_macro(it); \
105 s += t; \
106 sq += tq; \
107 t = sum[x - sumstep] + s; \
108 tq = sqsum[x - sqsumstep] + sq; \
109 sum[x] = t; \
110 sqsum[x] = tq; \
111 } \
112 } \
113 } \
114 else \
115 { \
116 if( sqsum == 0 ) \
117 { \
118 assert(0); \
119 return CV_NULLPTR_ERR; \
120 } \
121 \
122 buf = (sumtype*)cvStackAlloc((size.width + 1 )* sizeof(buf[0]));\
123 sum[-1] = tilted[-1] = 0; \
124 sqsum[-1] = 0; \
125 \
126 for( x = 0, s = 0, sq = 0; x < size.width; x++ ) \
127 { \
128 worktype it = src[x]; \
129 sumtype t = cast_macro(it); \
130 sqsumtype tq = cast_sqr_macro(it); \
131 buf[x] = tilted[x] = t; \
132 s += t; \
133 sq += tq; \
134 sum[x] = s; \
135 sqsum[x] = sq; \
136 } \
137 \
138 if( size.width == 1 ) \
139 buf[1] = 0; \
140 \
141 for( y = 1; y < size.height; y++ ) \
142 { \
143 worktype it; \
144 sumtype t0; \
145 sqsumtype tq0; \
146 \
147 src += srcstep; \
148 sum += sumstep; \
149 sqsum += sqsumstep; \
150 tilted += tiltedstep; \
151 \
152 it = src[0/*x*/]; \
153 s = t0 = cast_macro(it); \
154 sq = tq0 = cast_sqr_macro(it); \
155 \
156 sum[-1] = 0; \
157 sqsum[-1] = 0; \
158 /*tilted[-1] = buf[0];*/ \
159 tilted[-1] = tilted[-tiltedstep]; \
160 \
161 sum[0] = sum[-sumstep] + t0; \
162 sqsum[0] = sqsum[-sqsumstep] + tq0; \
163 tilted[0] = tilted[-tiltedstep] + t0 + buf[1]; \
164 \
165 for( x = 1; x < size.width - 1; x++ ) \
166 { \
167 sumtype t1 = buf[x]; \
168 buf[x-1] = t1 + t0; \
169 it = src[x]; \
170 t0 = cast_macro(it); \
171 tq0 = cast_sqr_macro(it); \
172 s += t0; \
173 sq += tq0; \
174 sum[x] = sum[x - sumstep] + s; \
175 sqsum[x] = sqsum[x - sqsumstep] + sq; \
176 t1 += buf[x+1] + t0 + tilted[x - tiltedstep - 1];\
177 tilted[x] = t1; \
178 } \
179 \
180 if( size.width > 1 ) \
181 { \
182 sumtype t1 = buf[x]; \
183 buf[x-1] = t1 + t0; \
184 it = src[x]; /*+*/ \
185 t0 = cast_macro(it); \
186 tq0 = cast_sqr_macro(it); \
187 s += t0; \
188 sq += tq0; \
189 sum[x] = sum[x - sumstep] + s; \
190 sqsum[x] = sqsum[x - sqsumstep] + sq; \
191 tilted[x] = t0 + t1 + tilted[x - tiltedstep - 1];\
192 buf[x] = t0; \
193 } \
194 } \
195 } \
196 \
197 return CV_OK; \
198 }
199
200
201 ICV_DEF_INTEGRAL_OP_C1( 8u32s, uchar, int, double, int, CV_NOP, CV_8TO32F_SQR )
202 ICV_DEF_INTEGRAL_OP_C1( 8u64f, uchar, double, double, int, CV_8TO32F, CV_8TO32F_SQR )
203 ICV_DEF_INTEGRAL_OP_C1( 32f64f, float, double, double, double, CV_NOP, CV_SQR )
204 ICV_DEF_INTEGRAL_OP_C1( 64f, double, double, double, double, CV_NOP, CV_SQR )
205
206
207 #define ICV_DEF_INTEGRAL_OP_CN( flavor, arrtype, sumtype, sqsumtype, \
208 worktype, cast_macro, cast_sqr_macro ) \
209 static CvStatus CV_STDCALL \
210 icvIntegralImage_##flavor##_CnR( const arrtype* src, int srcstep,\
211 sumtype* sum, int sumstep, \
212 sqsumtype* sqsum, int sqsumstep,\
213 CvSize size, int cn ) \
214 { \
215 int x, y; \
216 srcstep /= sizeof(src[0]); \
217 \
218 memset( sum, 0, (size.width+1)*cn*sizeof(sum[0])); \
219 sumstep /= sizeof(sum[0]); \
220 sum += sumstep + cn; \
221 \
222 if( sqsum ) \
223 { \
224 memset( sqsum, 0, (size.width+1)*cn*sizeof(sqsum[0])); \
225 sqsumstep /= sizeof(sqsum[0]); \
226 sqsum += sqsumstep + cn; \
227 } \
228 \
229 size.width *= cn; \
230 \
231 if( sqsum == 0 ) \
232 { \
233 for( y = 0; y < size.height; y++, src += srcstep, \
234 sum += sumstep ) \
235 { \
236 for( x = -cn; x < 0; x++ ) \
237 sum[x] = 0; \
238 \
239 for( x = 0; x < size.width; x++ ) \
240 sum[x] = cast_macro(src[x]) + sum[x - cn]; \
241 \
242 for( x = 0; x < size.width; x++ ) \
243 sum[x] = sum[x] + sum[x - sumstep]; \
244 } \
245 } \
246 else \
247 { \
248 for( y = 0; y < size.height; y++, src += srcstep, \
249 sum += sumstep, sqsum += sqsumstep ) \
250 { \
251 for( x = -cn; x < 0; x++ ) \
252 { \
253 sum[x] = 0; \
254 sqsum[x] = 0; \
255 } \
256 \
257 for( x = 0; x < size.width; x++ ) \
258 { \
259 worktype it = src[x]; \
260 sumtype t = cast_macro(it) + sum[x-cn]; \
261 sqsumtype tq = cast_sqr_macro(it) + sqsum[x-cn];\
262 sum[x] = t; \
263 sqsum[x] = tq; \
264 } \
265 \
266 for( x = 0; x < size.width; x++ ) \
267 { \
268 sumtype t = sum[x] + sum[x - sumstep]; \
269 sqsumtype tq = sqsum[x] + sqsum[x - sqsumstep]; \
270 sum[x] = t; \
271 sqsum[x] = tq; \
272 } \
273 } \
274 } \
275 \
276 return CV_OK; \
277 }
278
279
280 ICV_DEF_INTEGRAL_OP_CN( 8u32s, uchar, int, double, int, CV_NOP, CV_8TO32F_SQR )
281 ICV_DEF_INTEGRAL_OP_CN( 8u64f, uchar, double, double, int, CV_8TO32F, CV_8TO32F_SQR )
282 ICV_DEF_INTEGRAL_OP_CN( 32f64f, float, double, double, double, CV_NOP, CV_SQR )
283 ICV_DEF_INTEGRAL_OP_CN( 64f, double, double, double, double, CV_NOP, CV_SQR )
284
285
icvInitIntegralImageTable(CvFuncTable * table_c1,CvFuncTable * table_cn)286 static void icvInitIntegralImageTable( CvFuncTable* table_c1, CvFuncTable* table_cn )
287 {
288 table_c1->fn_2d[CV_8U] = (void*)icvIntegralImage_8u64f_C1R;
289 table_c1->fn_2d[CV_32F] = (void*)icvIntegralImage_32f64f_C1R;
290 table_c1->fn_2d[CV_64F] = (void*)icvIntegralImage_64f_C1R;
291
292 table_cn->fn_2d[CV_8U] = (void*)icvIntegralImage_8u64f_CnR;
293 table_cn->fn_2d[CV_32F] = (void*)icvIntegralImage_32f64f_CnR;
294 table_cn->fn_2d[CV_64F] = (void*)icvIntegralImage_64f_CnR;
295 }
296
297
298 typedef CvStatus (CV_STDCALL * CvIntegralImageFuncC1)(
299 const void* src, int srcstep, void* sum, int sumstep,
300 void* sqsum, int sqsumstep, void* tilted, int tiltedstep,
301 CvSize size );
302
303 typedef CvStatus (CV_STDCALL * CvIntegralImageFuncCn)(
304 const void* src, int srcstep, void* sum, int sumstep,
305 void* sqsum, int sqsumstep, CvSize size, int cn );
306
307 icvIntegral_8u32s_C1R_t icvIntegral_8u32s_C1R_p = 0;
308 icvSqrIntegral_8u32s64f_C1R_t icvSqrIntegral_8u32s64f_C1R_p = 0;
309
310 CV_IMPL void
cvIntegral(const CvArr * image,CvArr * sumImage,CvArr * sumSqImage,CvArr * tiltedSumImage)311 cvIntegral( const CvArr* image, CvArr* sumImage,
312 CvArr* sumSqImage, CvArr* tiltedSumImage )
313 {
314 static CvFuncTable tab_c1, tab_cn;
315 static int inittab = 0;
316
317 CV_FUNCNAME( "cvIntegralImage" );
318
319 __BEGIN__;
320
321 CvMat src_stub, *src = (CvMat*)image;
322 CvMat sum_stub, *sum = (CvMat*)sumImage;
323 CvMat sqsum_stub, *sqsum = (CvMat*)sumSqImage;
324 CvMat tilted_stub, *tilted = (CvMat*)tiltedSumImage;
325 int coi0 = 0, coi1 = 0, coi2 = 0, coi3 = 0;
326 int depth, cn;
327 int src_step, sum_step, sqsum_step, tilted_step;
328 CvIntegralImageFuncC1 func_c1 = 0;
329 CvIntegralImageFuncCn func_cn = 0;
330 CvSize size;
331
332 if( !inittab )
333 {
334 icvInitIntegralImageTable( &tab_c1, &tab_cn );
335 inittab = 1;
336 }
337
338 CV_CALL( src = cvGetMat( src, &src_stub, &coi0 ));
339 CV_CALL( sum = cvGetMat( sum, &sum_stub, &coi1 ));
340
341 if( sum->width != src->width + 1 ||
342 sum->height != src->height + 1 )
343 CV_ERROR( CV_StsUnmatchedSizes, "" );
344
345 if( (CV_MAT_DEPTH( sum->type ) != CV_64F &&
346 (CV_MAT_DEPTH( src->type ) != CV_8U ||
347 CV_MAT_DEPTH( sum->type ) != CV_32S )) ||
348 !CV_ARE_CNS_EQ( src, sum ))
349 CV_ERROR( CV_StsUnsupportedFormat,
350 "Sum array must have 64f type (or 32s type in case of 8u source array) "
351 "and the same number of channels as the source array" );
352
353 if( sqsum )
354 {
355 CV_CALL( sqsum = cvGetMat( sqsum, &sqsum_stub, &coi2 ));
356 if( !CV_ARE_SIZES_EQ( sum, sqsum ) )
357 CV_ERROR( CV_StsUnmatchedSizes, "" );
358 if( CV_MAT_DEPTH( sqsum->type ) != CV_64F || !CV_ARE_CNS_EQ( src, sqsum ))
359 CV_ERROR( CV_StsUnsupportedFormat,
360 "Squares sum array must be 64f "
361 "and the same number of channels as the source array" );
362 }
363
364 if( tilted )
365 {
366 if( !sqsum )
367 CV_ERROR( CV_StsNullPtr,
368 "Squared sum array must be passed if tilted sum array is passed" );
369
370 CV_CALL( tilted = cvGetMat( tilted, &tilted_stub, &coi3 ));
371 if( !CV_ARE_SIZES_EQ( sum, tilted ) )
372 CV_ERROR( CV_StsUnmatchedSizes, "" );
373 if( !CV_ARE_TYPES_EQ( sum, tilted ) )
374 CV_ERROR( CV_StsUnmatchedFormats,
375 "Sum and tilted sum must have the same types" );
376 if( CV_MAT_CN(tilted->type) != 1 )
377 CV_ERROR( CV_StsNotImplemented,
378 "Tilted sum can not be computed for multi-channel arrays" );
379 }
380
381 if( coi0 || coi1 || coi2 || coi3 )
382 CV_ERROR( CV_BadCOI, "COI is not supported by the function" );
383
384 depth = CV_MAT_DEPTH(src->type);
385 cn = CV_MAT_CN(src->type);
386
387 if( CV_MAT_DEPTH( sum->type ) == CV_32S )
388 {
389 func_c1 = (CvIntegralImageFuncC1)icvIntegralImage_8u32s_C1R;
390 func_cn = (CvIntegralImageFuncCn)icvIntegralImage_8u32s_CnR;
391 }
392 else
393 {
394 func_c1 = (CvIntegralImageFuncC1)tab_c1.fn_2d[depth];
395 func_cn = (CvIntegralImageFuncCn)tab_cn.fn_2d[depth];
396 if( !func_c1 && !func_cn )
397 CV_ERROR( CV_StsUnsupportedFormat, "This source image format is unsupported" );
398 }
399
400 size = cvGetMatSize(src);
401 src_step = src->step ? src->step : CV_STUB_STEP;
402 sum_step = sum->step ? sum->step : CV_STUB_STEP;
403 sqsum_step = !sqsum ? 0 : sqsum->step ? sqsum->step : CV_STUB_STEP;
404 tilted_step = !tilted ? 0 : tilted->step ? tilted->step : CV_STUB_STEP;
405
406 if( cn == 1 )
407 {
408 if( depth == CV_8U && !tilted && CV_MAT_DEPTH(sum->type) == CV_32S )
409 {
410 if( !sqsum && icvIntegral_8u32s_C1R_p &&
411 icvIntegral_8u32s_C1R_p( src->data.ptr, src_step,
412 sum->data.i, sum_step, size, 0 ) >= 0 )
413 EXIT;
414
415 if( sqsum && icvSqrIntegral_8u32s64f_C1R_p &&
416 icvSqrIntegral_8u32s64f_C1R_p( src->data.ptr, src_step, sum->data.i,
417 sum_step, sqsum->data.db, sqsum_step, size, 0, 0 ) >= 0 )
418 EXIT;
419 }
420
421 IPPI_CALL( func_c1( src->data.ptr, src_step, sum->data.ptr, sum_step,
422 sqsum ? sqsum->data.ptr : 0, sqsum_step,
423 tilted ? tilted->data.ptr : 0, tilted_step, size ));
424 }
425 else
426 {
427 IPPI_CALL( func_cn( src->data.ptr, src_step, sum->data.ptr, sum_step,
428 sqsum ? sqsum->data.ptr : 0, sqsum_step, size, cn ));
429 }
430
431 __END__;
432 }
433
434
435 /* End of file. */
436