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 "_cxcore.h"
43
44 #ifdef HAVE_CONFIG_H
45 #include <cvconfig.h>
46 #endif
47
48 #define ICV_MATH_BLOCK_SIZE 256
49
50 #define _CV_SQRT_MAGIC 0xbe6f0000
51
52 #define _CV_SQRT_MAGIC_DBL CV_BIG_UINT(0xbfcd460000000000)
53
54 #define _CV_ATAN_CF0 (-15.8131890796f)
55 #define _CV_ATAN_CF1 (61.0941945596f)
56 #define _CV_ATAN_CF2 0.f /*(-0.140500406322f)*/
57
58 static const float icvAtanTab[8] = { 0.f + _CV_ATAN_CF2, 90.f - _CV_ATAN_CF2,
59 180.f - _CV_ATAN_CF2, 90.f + _CV_ATAN_CF2,
60 360.f - _CV_ATAN_CF2, 270.f + _CV_ATAN_CF2,
61 180.f + _CV_ATAN_CF2, 270.f - _CV_ATAN_CF2
62 };
63
64 static const int icvAtanSign[8] =
65 { 0, 0x80000000, 0x80000000, 0, 0x80000000, 0, 0, 0x80000000 };
66
67 CV_IMPL float
cvFastArctan(float y,float x)68 cvFastArctan( float y, float x )
69 {
70 Cv32suf _x, _y;
71 int ix, iy, ygx, idx;
72 double z;
73
74 _x.f = x; _y.f = y;
75 ix = _x.i; iy = _y.i;
76 idx = (ix < 0) * 2 + (iy < 0) * 4;
77
78 ix &= 0x7fffffff;
79 iy &= 0x7fffffff;
80
81 ygx = (iy <= ix) - 1;
82 idx -= ygx;
83
84 idx &= ((ix == 0) - 1) | ((iy == 0) - 1);
85
86 /* swap ix and iy if ix < iy */
87 ix ^= iy & ygx;
88 iy ^= ix & ygx;
89 ix ^= iy & ygx;
90
91 _y.i = iy ^ icvAtanSign[idx];
92
93 /* ix = ix != 0 ? ix : 1.f */
94 _x.i = ((ix ^ CV_1F) & ((ix == 0) - 1)) ^ CV_1F;
95
96 z = _y.f / _x.f;
97 return (float)((_CV_ATAN_CF0*fabs(z) + _CV_ATAN_CF1)*z + icvAtanTab[idx]);
98 }
99
100
101 IPCVAPI_IMPL( CvStatus, icvFastArctan_32f,
102 (const float *__y, const float *__x, float *angle, int len ), (__y, __x, angle, len) )
103 {
104 int i = 0;
105 const int *y = (const int*)__y, *x = (const int*)__x;
106
107 if( !(y && x && angle && len >= 0) )
108 return CV_BADFACTOR_ERR;
109
110 /* unrolled by 4 loop */
111 for( ; i <= len - 4; i += 4 )
112 {
113 int j, idx[4];
114 float xf[4], yf[4];
115 double d = 1.;
116
117 /* calc numerators and denominators */
118 for( j = 0; j < 4; j++ )
119 {
120 int ix = x[i + j], iy = y[i + j];
121 int ygx, k = (ix < 0) * 2 + (iy < 0) * 4;
122 Cv32suf _x, _y;
123
124 ix &= 0x7fffffff;
125 iy &= 0x7fffffff;
126
127 ygx = (iy <= ix) - 1;
128 k -= ygx;
129
130 k &= ((ix == 0) - 1) | ((iy == 0) - 1);
131
132 /* swap ix and iy if ix < iy */
133 ix ^= iy & ygx;
134 iy ^= ix & ygx;
135 ix ^= iy & ygx;
136
137 _y.i = iy ^ icvAtanSign[k];
138
139 /* ix = ix != 0 ? ix : 1.f */
140 _x.i = ((ix ^ CV_1F) & ((ix == 0) - 1)) ^ CV_1F;
141 idx[j] = k;
142 yf[j] = _y.f;
143 d *= (xf[j] = _x.f);
144 }
145
146 d = 1. / d;
147
148 {
149 double b = xf[2] * xf[3], a = xf[0] * xf[1];
150
151 float z0 = (float) (yf[0] * xf[1] * b * d);
152 float z1 = (float) (yf[1] * xf[0] * b * d);
153 float z2 = (float) (yf[2] * xf[3] * a * d);
154 float z3 = (float) (yf[3] * xf[2] * a * d);
155
156 z0 = (float)((_CV_ATAN_CF0*fabs(z0) + _CV_ATAN_CF1)*z0 + icvAtanTab[idx[0]]);
157 z1 = (float)((_CV_ATAN_CF0*fabs(z1) + _CV_ATAN_CF1)*z1 + icvAtanTab[idx[1]]);
158 z2 = (float)((_CV_ATAN_CF0*fabs(z2) + _CV_ATAN_CF1)*z2 + icvAtanTab[idx[2]]);
159 z3 = (float)((_CV_ATAN_CF0*fabs(z3) + _CV_ATAN_CF1)*z3 + icvAtanTab[idx[3]]);
160
161 angle[i] = z0;
162 angle[i+1] = z1;
163 angle[i+2] = z2;
164 angle[i+3] = z3;
165 }
166 }
167
168 /* process the rest */
169 for( ; i < len; i++ )
170 angle[i] = cvFastArctan( __y[i], __x[i] );
171
172 return CV_OK;
173 }
174
175
176 /* ************************************************************************** *\
177 Fast cube root by Ken Turkowski
178 (http://www.worldserver.com/turk/computergraphics/papers.html)
179 \* ************************************************************************** */
cvCbrt(float value)180 CV_IMPL float cvCbrt( float value )
181 {
182 float fr;
183 Cv32suf v, m;
184 int ix, s;
185 int ex, shx;
186
187 v.f = value;
188 ix = v.i & 0x7fffffff;
189 s = v.i & 0x80000000;
190 ex = (ix >> 23) - 127;
191 shx = ex % 3;
192 shx -= shx >= 0 ? 3 : 0;
193 ex = (ex - shx) / 3; /* exponent of cube root */
194 v.i = (ix & ((1<<23)-1)) | ((shx + 127)<<23);
195 fr = v.f;
196
197 /* 0.125 <= fr < 1.0 */
198 /* Use quartic rational polynomial with error < 2^(-24) */
199 fr = (float)(((((45.2548339756803022511987494 * fr +
200 192.2798368355061050458134625) * fr +
201 119.1654824285581628956914143) * fr +
202 13.43250139086239872172837314) * fr +
203 0.1636161226585754240958355063)/
204 ((((14.80884093219134573786480845 * fr +
205 151.9714051044435648658557668) * fr +
206 168.5254414101568283957668343) * fr +
207 33.9905941350215598754191872) * fr +
208 1.0));
209
210 /* fr *= 2^ex * sign */
211 m.f = value;
212 v.f = fr;
213 v.i = (v.i + (ex << 23) + s) & (m.i*2 != 0 ? -1 : 0);
214 return v.f;
215 }
216
217 //static const double _0_5 = 0.5, _1_5 = 1.5;
218
219 IPCVAPI_IMPL( CvStatus, icvInvSqrt_32f, (const float *src, float *dst, int len), (src, dst, len) )
220 {
221 int i = 0;
222
223 if( !(src && dst && len >= 0) )
224 return CV_BADFACTOR_ERR;
225
226 for( ; i < len; i++ )
227 dst[i] = (float)(1.f/sqrt(src[i]));
228
229 return CV_OK;
230 }
231
232
233 IPCVAPI_IMPL( CvStatus, icvSqrt_32f, (const float *src, float *dst, int len), (src, dst, len) )
234 {
235 int i = 0;
236
237 if( !(src && dst && len >= 0) )
238 return CV_BADFACTOR_ERR;
239
240 for( ; i < len; i++ )
241 dst[i] = (float)sqrt(src[i]);
242
243 return CV_OK;
244 }
245
246
247 IPCVAPI_IMPL( CvStatus, icvSqrt_64f, (const double *src, double *dst, int len), (src, dst, len) )
248 {
249 int i = 0;
250
251 if( !(src && dst && len >= 0) )
252 return CV_BADFACTOR_ERR;
253
254 for( ; i < len; i++ )
255 dst[i] = sqrt(src[i]);
256
257 return CV_OK;
258 }
259
260
261 IPCVAPI_IMPL( CvStatus, icvInvSqrt_64f, (const double *src, double *dst, int len), (src, dst, len) )
262 {
263 int i = 0;
264
265 if( !(src && dst && len >= 0) )
266 return CV_BADFACTOR_ERR;
267
268 for( ; i < len; i++ )
269 dst[i] = 1./sqrt(src[i]);
270
271 return CV_OK;
272 }
273
274
275 #define ICV_DEF_SQR_MAGNITUDE_FUNC(flavor, arrtype, magtype)\
276 static CvStatus CV_STDCALL \
277 icvSqrMagnitude_##flavor(const arrtype* x, const arrtype* y,\
278 magtype* mag, int len) \
279 { \
280 int i; \
281 \
282 for( i = 0; i <= len - 4; i += 4 ) \
283 { \
284 magtype x0 = (magtype)x[i], y0 = (magtype)y[i]; \
285 magtype x1 = (magtype)x[i+1], y1 = (magtype)y[i+1]; \
286 \
287 x0 = x0*x0 + y0*y0; \
288 x1 = x1*x1 + y1*y1; \
289 mag[i] = x0; \
290 mag[i+1] = x1; \
291 x0 = (magtype)x[i+2], y0 = (magtype)y[i+2]; \
292 x1 = (magtype)x[i+3], y1 = (magtype)y[i+3]; \
293 x0 = x0*x0 + y0*y0; \
294 x1 = x1*x1 + y1*y1; \
295 mag[i+2] = x0; \
296 mag[i+3] = x1; \
297 } \
298 \
299 for( ; i < len; i++ ) \
300 { \
301 magtype x0 = (magtype)x[i], y0 = (magtype)y[i]; \
302 mag[i] = x0*x0 + y0*y0; \
303 } \
304 \
305 return CV_OK; \
306 }
307
308
309 ICV_DEF_SQR_MAGNITUDE_FUNC( 32f, float, float )
310 ICV_DEF_SQR_MAGNITUDE_FUNC( 64f, double, double )
311
312 /****************************************************************************************\
313 * Cartezian -> Polar *
314 \****************************************************************************************/
315
316 CV_IMPL void
cvCartToPolar(const CvArr * xarr,const CvArr * yarr,CvArr * magarr,CvArr * anglearr,int angle_in_degrees)317 cvCartToPolar( const CvArr* xarr, const CvArr* yarr,
318 CvArr* magarr, CvArr* anglearr,
319 int angle_in_degrees )
320 {
321 CV_FUNCNAME( "cvCartToPolar" );
322
323 __BEGIN__;
324
325 float* mag_buffer = 0;
326 float* x_buffer = 0;
327 float* y_buffer = 0;
328 int block_size = 0;
329 CvMat xstub, *xmat = (CvMat*)xarr;
330 CvMat ystub, *ymat = (CvMat*)yarr;
331 CvMat magstub, *mag = (CvMat*)magarr;
332 CvMat anglestub, *angle = (CvMat*)anglearr;
333 int coi1 = 0, coi2 = 0, coi3 = 0, coi4 = 0;
334 int depth;
335 CvSize size;
336 int x, y;
337 int cont_flag = CV_MAT_CONT_FLAG;
338
339 if( !CV_IS_MAT(xmat))
340 CV_CALL( xmat = cvGetMat( xmat, &xstub, &coi1 ));
341
342 if( !CV_IS_MAT(ymat))
343 CV_CALL( ymat = cvGetMat( ymat, &ystub, &coi2 ));
344
345 if( !CV_ARE_TYPES_EQ( xmat, ymat ) )
346 CV_ERROR_FROM_CODE( CV_StsUnmatchedFormats );
347
348 if( !CV_ARE_SIZES_EQ( xmat, ymat ) )
349 CV_ERROR_FROM_CODE( CV_StsUnmatchedSizes );
350
351 depth = CV_MAT_DEPTH( xmat->type );
352 if( depth < CV_32F )
353 CV_ERROR( CV_StsUnsupportedFormat, "" );
354
355 if( mag )
356 {
357 CV_CALL( mag = cvGetMat( mag, &magstub, &coi3 ));
358
359 if( !CV_ARE_TYPES_EQ( mag, xmat ) )
360 CV_ERROR_FROM_CODE( CV_StsUnmatchedFormats );
361
362 if( !CV_ARE_SIZES_EQ( mag, xmat ) )
363 CV_ERROR_FROM_CODE( CV_StsUnmatchedSizes );
364 cont_flag = mag->type;
365 }
366
367 if( angle )
368 {
369 CV_CALL( angle = cvGetMat( angle, &anglestub, &coi4 ));
370
371 if( !CV_ARE_TYPES_EQ( angle, xmat ) )
372 CV_ERROR_FROM_CODE( CV_StsUnmatchedFormats );
373
374 if( !CV_ARE_SIZES_EQ( angle, xmat ) )
375 CV_ERROR_FROM_CODE( CV_StsUnmatchedSizes );
376 cont_flag &= angle->type;
377 }
378
379 if( coi1 != 0 || coi2 != 0 || coi3 != 0 || coi4 != 0 )
380 CV_ERROR( CV_BadCOI, "" );
381
382 size = cvGetMatSize(xmat);
383 size.width *= CV_MAT_CN(xmat->type);
384
385 if( CV_IS_MAT_CONT( xmat->type & ymat->type & cont_flag ))
386 {
387 size.width *= size.height;
388 size.height = 1;
389 }
390
391 block_size = MIN( size.width, ICV_MATH_BLOCK_SIZE );
392 if( depth == CV_64F && angle )
393 {
394 x_buffer = (float*)cvStackAlloc( block_size*sizeof(float));
395 y_buffer = (float*)cvStackAlloc( block_size*sizeof(float));
396 }
397 else if( depth == CV_32F && mag )
398 {
399 mag_buffer = (float*)cvStackAlloc( block_size*sizeof(float));
400 }
401
402 if( depth == CV_32F )
403 {
404 for( y = 0; y < size.height; y++ )
405 {
406 float* x_data = (float*)(xmat->data.ptr + xmat->step*y);
407 float* y_data = (float*)(ymat->data.ptr + ymat->step*y);
408 float* mag_data = mag ? (float*)(mag->data.ptr + mag->step*y) : 0;
409 float* angle_data = angle ? (float*)(angle->data.ptr + angle->step*y) : 0;
410
411 for( x = 0; x < size.width; x += block_size )
412 {
413 int len = MIN( size.width - x, block_size );
414
415 if( mag )
416 icvSqrMagnitude_32f( x_data + x, y_data + x, mag_buffer, len );
417
418 if( angle )
419 {
420 icvFastArctan_32f( y_data + x, x_data + x, angle_data + x, len );
421 if( !angle_in_degrees )
422 icvScale_32f( angle_data + x, angle_data + x,
423 len, (float)(CV_PI/180.), 0 );
424 }
425
426 if( mag )
427 icvSqrt_32f( mag_buffer, mag_data + x, len );
428 }
429 }
430 }
431 else
432 {
433 for( y = 0; y < size.height; y++ )
434 {
435 double* x_data = (double*)(xmat->data.ptr + xmat->step*y);
436 double* y_data = (double*)(ymat->data.ptr + ymat->step*y);
437 double* mag_data = mag ? (double*)(mag->data.ptr + mag->step*y) : 0;
438 double* angle_data = angle ? (double*)(angle->data.ptr + angle->step*y) : 0;
439
440 for( x = 0; x < size.width; x += block_size )
441 {
442 int len = MIN( size.width - x, block_size );
443
444 if( angle )
445 {
446 icvCvt_64f32f( x_data + x, x_buffer, len );
447 icvCvt_64f32f( y_data + x, y_buffer, len );
448 }
449
450 if( mag )
451 {
452 icvSqrMagnitude_64f( x_data + x, y_data + x, mag_data + x, len );
453 icvSqrt_64f( mag_data + x, mag_data + x, len );
454 }
455
456 if( angle )
457 {
458 icvFastArctan_32f( y_buffer, x_buffer, x_buffer, len );
459 if( !angle_in_degrees )
460 icvScale_32f( x_buffer, x_buffer, len, (float)(CV_PI/180.), 0 );
461 icvCvt_32f64f( x_buffer, angle_data + x, len );
462 }
463 }
464 }
465 }
466
467 __END__;
468 }
469
470
471 /****************************************************************************************\
472 * Polar -> Cartezian *
473 \****************************************************************************************/
474
475 static CvStatus CV_STDCALL
icvSinCos_32f(const float * angle,float * sinval,float * cosval,int len,int angle_in_degrees)476 icvSinCos_32f( const float *angle,float *sinval, float* cosval,
477 int len, int angle_in_degrees )
478 {
479 const int N = 64;
480
481 static const double sin_table[] =
482 {
483 0.00000000000000000000, 0.09801714032956060400,
484 0.19509032201612825000, 0.29028467725446233000,
485 0.38268343236508978000, 0.47139673682599764000,
486 0.55557023301960218000, 0.63439328416364549000,
487 0.70710678118654746000, 0.77301045336273699000,
488 0.83146961230254524000, 0.88192126434835494000,
489 0.92387953251128674000, 0.95694033573220894000,
490 0.98078528040323043000, 0.99518472667219682000,
491 1.00000000000000000000, 0.99518472667219693000,
492 0.98078528040323043000, 0.95694033573220894000,
493 0.92387953251128674000, 0.88192126434835505000,
494 0.83146961230254546000, 0.77301045336273710000,
495 0.70710678118654757000, 0.63439328416364549000,
496 0.55557023301960218000, 0.47139673682599786000,
497 0.38268343236508989000, 0.29028467725446239000,
498 0.19509032201612861000, 0.09801714032956082600,
499 0.00000000000000012246, -0.09801714032956059000,
500 -0.19509032201612836000, -0.29028467725446211000,
501 -0.38268343236508967000, -0.47139673682599764000,
502 -0.55557023301960196000, -0.63439328416364527000,
503 -0.70710678118654746000, -0.77301045336273666000,
504 -0.83146961230254524000, -0.88192126434835494000,
505 -0.92387953251128652000, -0.95694033573220882000,
506 -0.98078528040323032000, -0.99518472667219693000,
507 -1.00000000000000000000, -0.99518472667219693000,
508 -0.98078528040323043000, -0.95694033573220894000,
509 -0.92387953251128663000, -0.88192126434835505000,
510 -0.83146961230254546000, -0.77301045336273688000,
511 -0.70710678118654768000, -0.63439328416364593000,
512 -0.55557023301960218000, -0.47139673682599792000,
513 -0.38268343236509039000, -0.29028467725446250000,
514 -0.19509032201612872000, -0.09801714032956050600,
515 };
516
517 static const double k2 = (2*CV_PI)/N;
518
519 static const double sin_a0 = -0.166630293345647*k2*k2*k2;
520 static const double sin_a2 = k2;
521
522 static const double cos_a0 = -0.499818138450326*k2*k2;
523 /*static const double cos_a2 = 1;*/
524
525 double k1;
526 int i;
527
528 if( !angle_in_degrees )
529 k1 = N/(2*CV_PI);
530 else
531 k1 = N/360.;
532
533 for( i = 0; i < len; i++ )
534 {
535 double t = angle[i]*k1;
536 int it = cvRound(t);
537 t -= it;
538 int sin_idx = it & (N - 1);
539 int cos_idx = (N/4 - sin_idx) & (N - 1);
540
541 double sin_b = (sin_a0*t*t + sin_a2)*t;
542 double cos_b = cos_a0*t*t + 1;
543
544 double sin_a = sin_table[sin_idx];
545 double cos_a = sin_table[cos_idx];
546
547 double sin_val = sin_a*cos_b + cos_a*sin_b;
548 double cos_val = cos_a*cos_b - sin_a*sin_b;
549
550 sinval[i] = (float)sin_val;
551 cosval[i] = (float)cos_val;
552 }
553
554 return CV_OK;
555 }
556
557
558 CV_IMPL void
cvPolarToCart(const CvArr * magarr,const CvArr * anglearr,CvArr * xarr,CvArr * yarr,int angle_in_degrees)559 cvPolarToCart( const CvArr* magarr, const CvArr* anglearr,
560 CvArr* xarr, CvArr* yarr, int angle_in_degrees )
561 {
562 CV_FUNCNAME( "cvPolarToCart" );
563
564 __BEGIN__;
565
566 float* x_buffer = 0;
567 float* y_buffer = 0;
568 int block_size = 0;
569 CvMat xstub, *xmat = (CvMat*)xarr;
570 CvMat ystub, *ymat = (CvMat*)yarr;
571 CvMat magstub, *mag = (CvMat*)magarr;
572 CvMat anglestub, *angle = (CvMat*)anglearr;
573 int coi1 = 0, coi2 = 0, coi3 = 0, coi4 = 0;
574 int depth;
575 CvSize size;
576 int x, y;
577 int cont_flag;
578
579 if( !CV_IS_MAT(angle))
580 CV_CALL( angle = cvGetMat( angle, &anglestub, &coi4 ));
581
582 depth = CV_MAT_DEPTH( angle->type );
583 if( depth < CV_32F )
584 CV_ERROR( CV_StsUnsupportedFormat, "" );
585 cont_flag = angle->type;
586
587 if( mag )
588 {
589 if( !CV_IS_MAT(mag))
590 CV_CALL( mag = cvGetMat( mag, &magstub, &coi3 ));
591
592 if( !CV_ARE_TYPES_EQ( angle, mag ) )
593 CV_ERROR_FROM_CODE( CV_StsUnmatchedFormats );
594
595 if( !CV_ARE_SIZES_EQ( angle, mag ) )
596 CV_ERROR_FROM_CODE( CV_StsUnmatchedSizes );
597
598 cont_flag &= mag->type;
599 }
600
601 if( xmat )
602 {
603 if( !CV_IS_MAT(xmat))
604 CV_CALL( xmat = cvGetMat( xmat, &xstub, &coi1 ));
605
606 if( !CV_ARE_TYPES_EQ( angle, xmat ) )
607 CV_ERROR_FROM_CODE( CV_StsUnmatchedFormats );
608
609 if( !CV_ARE_SIZES_EQ( angle, xmat ) )
610 CV_ERROR_FROM_CODE( CV_StsUnmatchedSizes );
611
612 cont_flag &= xmat->type;
613 }
614
615 if( ymat )
616 {
617 if( !CV_IS_MAT(ymat))
618 CV_CALL( ymat = cvGetMat( ymat, &ystub, &coi2 ));
619
620 if( !CV_ARE_TYPES_EQ( angle, ymat ) )
621 CV_ERROR_FROM_CODE( CV_StsUnmatchedFormats );
622
623 if( !CV_ARE_SIZES_EQ( angle, ymat ) )
624 CV_ERROR_FROM_CODE( CV_StsUnmatchedSizes );
625
626 cont_flag &= ymat->type;
627 }
628
629 if( coi1 != 0 || coi2 != 0 || coi3 != 0 || coi4 != 0 )
630 CV_ERROR( CV_BadCOI, "" );
631
632 size = cvGetMatSize(angle);
633 size.width *= CV_MAT_CN(angle->type);
634
635 if( CV_IS_MAT_CONT( cont_flag ))
636 {
637 size.width *= size.height;
638 size.height = 1;
639 }
640
641 block_size = MIN( size.width, ICV_MATH_BLOCK_SIZE );
642 x_buffer = (float*)cvStackAlloc( block_size*sizeof(float));
643 y_buffer = (float*)cvStackAlloc( block_size*sizeof(float));
644
645 if( depth == CV_32F )
646 {
647 for( y = 0; y < size.height; y++ )
648 {
649 float* x_data = (float*)(xmat ? xmat->data.ptr + xmat->step*y : 0);
650 float* y_data = (float*)(ymat ? ymat->data.ptr + ymat->step*y : 0);
651 float* mag_data = (float*)(mag ? mag->data.ptr + mag->step*y : 0);
652 float* angle_data = (float*)(angle->data.ptr + angle->step*y);
653
654 for( x = 0; x < size.width; x += block_size )
655 {
656 int i, len = MIN( size.width - x, block_size );
657
658 icvSinCos_32f( angle_data+x, y_buffer, x_buffer, len, angle_in_degrees );
659
660 for( i = 0; i < len; i++ )
661 {
662 float tx = x_buffer[i];
663 float ty = y_buffer[i];
664
665 if( mag_data )
666 {
667 float magval = mag_data[x+i];
668 tx *= magval;
669 ty *= magval;
670 }
671
672 if( xmat )
673 x_data[x+i] = tx;
674 if( ymat )
675 y_data[x+i] = ty;
676 }
677 }
678 }
679 }
680 else
681 {
682 for( y = 0; y < size.height; y++ )
683 {
684 double* x_data = (double*)(xmat ? xmat->data.ptr + xmat->step*y : 0);
685 double* y_data = (double*)(ymat ? ymat->data.ptr + ymat->step*y : 0);
686 double* mag_data = (double*)(mag ? mag->data.ptr + mag->step*y : 0);
687 double* angle_data = (double*)(angle->data.ptr + angle->step*y);
688 double C = angle_in_degrees ? CV_PI/180. : 1;
689
690 for( x = 0; x < size.width; x++ )
691 {
692 double phi = angle_data[x]*C;
693 double magval = mag_data ? mag_data[x] : 1.;
694 if( xmat )
695 x_data[x] = cos(phi)*magval;
696 if( ymat )
697 y_data[x] = sin(phi)*magval;
698 }
699 }
700 }
701
702 __END__;
703 }
704
705 /****************************************************************************************\
706 * E X P *
707 \****************************************************************************************/
708
709 typedef union
710 {
711 struct {
712 #if ( defined( WORDS_BIGENDIAN ) && !defined( OPENCV_UNIVERSAL_BUILD ) ) || defined( __BIG_ENDIAN__ )
713 int hi;
714 int lo;
715 #else
716 int lo;
717 int hi;
718 #endif
719 } i;
720 double d;
721 }
722 DBLINT;
723
724 #define EXPTAB_SCALE 6
725 #define EXPTAB_MASK ((1 << EXPTAB_SCALE) - 1)
726
727 #define EXPPOLY_32F_A0 .9670371139572337719125840413672004409288e-2
728
729 static const double icvExpTab[] = {
730 1.0 * EXPPOLY_32F_A0,
731 1.0108892860517004600204097905619 * EXPPOLY_32F_A0,
732 1.0218971486541166782344801347833 * EXPPOLY_32F_A0,
733 1.0330248790212284225001082839705 * EXPPOLY_32F_A0,
734 1.0442737824274138403219664787399 * EXPPOLY_32F_A0,
735 1.0556451783605571588083413251529 * EXPPOLY_32F_A0,
736 1.0671404006768236181695211209928 * EXPPOLY_32F_A0,
737 1.0787607977571197937406800374385 * EXPPOLY_32F_A0,
738 1.0905077326652576592070106557607 * EXPPOLY_32F_A0,
739 1.1023825833078409435564142094256 * EXPPOLY_32F_A0,
740 1.1143867425958925363088129569196 * EXPPOLY_32F_A0,
741 1.126521618608241899794798643787 * EXPPOLY_32F_A0,
742 1.1387886347566916537038302838415 * EXPPOLY_32F_A0,
743 1.151189229952982705817759635202 * EXPPOLY_32F_A0,
744 1.1637248587775775138135735990922 * EXPPOLY_32F_A0,
745 1.1763969916502812762846457284838 * EXPPOLY_32F_A0,
746 1.1892071150027210667174999705605 * EXPPOLY_32F_A0,
747 1.2021567314527031420963969574978 * EXPPOLY_32F_A0,
748 1.2152473599804688781165202513388 * EXPPOLY_32F_A0,
749 1.2284805361068700056940089577928 * EXPPOLY_32F_A0,
750 1.2418578120734840485936774687266 * EXPPOLY_32F_A0,
751 1.2553807570246910895793906574423 * EXPPOLY_32F_A0,
752 1.2690509571917332225544190810323 * EXPPOLY_32F_A0,
753 1.2828700160787782807266697810215 * EXPPOLY_32F_A0,
754 1.2968395546510096659337541177925 * EXPPOLY_32F_A0,
755 1.3109612115247643419229917863308 * EXPPOLY_32F_A0,
756 1.3252366431597412946295370954987 * EXPPOLY_32F_A0,
757 1.3396675240533030053600306697244 * EXPPOLY_32F_A0,
758 1.3542555469368927282980147401407 * EXPPOLY_32F_A0,
759 1.3690024229745906119296011329822 * EXPPOLY_32F_A0,
760 1.3839098819638319548726595272652 * EXPPOLY_32F_A0,
761 1.3989796725383111402095281367152 * EXPPOLY_32F_A0,
762 1.4142135623730950488016887242097 * EXPPOLY_32F_A0,
763 1.4296133383919700112350657782751 * EXPPOLY_32F_A0,
764 1.4451808069770466200370062414717 * EXPPOLY_32F_A0,
765 1.4609177941806469886513028903106 * EXPPOLY_32F_A0,
766 1.476826145939499311386907480374 * EXPPOLY_32F_A0,
767 1.4929077282912648492006435314867 * EXPPOLY_32F_A0,
768 1.5091644275934227397660195510332 * EXPPOLY_32F_A0,
769 1.5255981507445383068512536895169 * EXPPOLY_32F_A0,
770 1.5422108254079408236122918620907 * EXPPOLY_32F_A0,
771 1.5590044002378369670337280894749 * EXPPOLY_32F_A0,
772 1.5759808451078864864552701601819 * EXPPOLY_32F_A0,
773 1.5931421513422668979372486431191 * EXPPOLY_32F_A0,
774 1.6104903319492543081795206673574 * EXPPOLY_32F_A0,
775 1.628027421857347766848218522014 * EXPPOLY_32F_A0,
776 1.6457554781539648445187567247258 * EXPPOLY_32F_A0,
777 1.6636765803267364350463364569764 * EXPPOLY_32F_A0,
778 1.6817928305074290860622509524664 * EXPPOLY_32F_A0,
779 1.7001063537185234695013625734975 * EXPPOLY_32F_A0,
780 1.7186192981224779156293443764563 * EXPPOLY_32F_A0,
781 1.7373338352737062489942020818722 * EXPPOLY_32F_A0,
782 1.7562521603732994831121606193753 * EXPPOLY_32F_A0,
783 1.7753764925265212525505592001993 * EXPPOLY_32F_A0,
784 1.7947090750031071864277032421278 * EXPPOLY_32F_A0,
785 1.8142521755003987562498346003623 * EXPPOLY_32F_A0,
786 1.8340080864093424634870831895883 * EXPPOLY_32F_A0,
787 1.8539791250833855683924530703377 * EXPPOLY_32F_A0,
788 1.8741676341102999013299989499544 * EXPPOLY_32F_A0,
789 1.8945759815869656413402186534269 * EXPPOLY_32F_A0,
790 1.9152065613971472938726112702958 * EXPPOLY_32F_A0,
791 1.9360617934922944505980559045667 * EXPPOLY_32F_A0,
792 1.9571441241754002690183222516269 * EXPPOLY_32F_A0,
793 1.9784560263879509682582499181312 * EXPPOLY_32F_A0,
794 };
795
796 static const double exp_prescale = 1.4426950408889634073599246810019 * (1 << EXPTAB_SCALE);
797 static const double exp_postscale = 1./(1 << EXPTAB_SCALE);
798 static const double exp_max_val = 3000.*(1 << EXPTAB_SCALE); // log10(DBL_MAX) < 3000
799
800 IPCVAPI_IMPL( CvStatus, icvExp_32f, ( const float *_x, float *y, int n ), (_x, y, n) )
801 {
802 static const double
803 EXPPOLY_32F_A4 = 1.000000000000002438532970795181890933776 / EXPPOLY_32F_A0,
804 EXPPOLY_32F_A3 = .6931471805521448196800669615864773144641 / EXPPOLY_32F_A0,
805 EXPPOLY_32F_A2 = .2402265109513301490103372422686535526573 / EXPPOLY_32F_A0,
806 EXPPOLY_32F_A1 = .5550339366753125211915322047004666939128e-1 / EXPPOLY_32F_A0;
807
808 #undef EXPPOLY
809 #define EXPPOLY(x) \
810 (((((x) + EXPPOLY_32F_A1)*(x) + EXPPOLY_32F_A2)*(x) + EXPPOLY_32F_A3)*(x) + EXPPOLY_32F_A4)
811
812 int i = 0;
813 DBLINT buf[4];
814 const Cv32suf* x = (const Cv32suf*)_x;
815
816 if( !x || !y )
817 return CV_NULLPTR_ERR;
818 if( n <= 0 )
819 return CV_BADSIZE_ERR;
820
821 buf[0].i.lo = buf[1].i.lo = buf[2].i.lo = buf[3].i.lo = 0;
822
823 for( ; i <= n - 4; i += 4 )
824 {
825 double x0 = x[i].f * exp_prescale;
826 double x1 = x[i + 1].f * exp_prescale;
827 double x2 = x[i + 2].f * exp_prescale;
828 double x3 = x[i + 3].f * exp_prescale;
829 int val0, val1, val2, val3, t;
830
831 if( ((x[i].i >> 23) & 255) > 127 + 10 )
832 x0 = x[i].i < 0 ? -exp_max_val : exp_max_val;
833
834 if( ((x[i+1].i >> 23) & 255) > 127 + 10 )
835 x1 = x[i+1].i < 0 ? -exp_max_val : exp_max_val;
836
837 if( ((x[i+2].i >> 23) & 255) > 127 + 10 )
838 x2 = x[i+2].i < 0 ? -exp_max_val : exp_max_val;
839
840 if( ((x[i+3].i >> 23) & 255) > 127 + 10 )
841 x3 = x[i+3].i < 0 ? -exp_max_val : exp_max_val;
842
843 val0 = cvRound(x0);
844 val1 = cvRound(x1);
845 val2 = cvRound(x2);
846 val3 = cvRound(x3);
847
848 x0 = (x0 - val0)*exp_postscale;
849 x1 = (x1 - val1)*exp_postscale;
850 x2 = (x2 - val2)*exp_postscale;
851 x3 = (x3 - val3)*exp_postscale;
852
853 t = (val0 >> EXPTAB_SCALE) + 1023;
854 t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
855 buf[0].i.hi = t << 20;
856
857 t = (val1 >> EXPTAB_SCALE) + 1023;
858 t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
859 buf[1].i.hi = t << 20;
860
861 t = (val2 >> EXPTAB_SCALE) + 1023;
862 t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
863 buf[2].i.hi = t << 20;
864
865 t = (val3 >> EXPTAB_SCALE) + 1023;
866 t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
867 buf[3].i.hi = t << 20;
868
869 x0 = buf[0].d * icvExpTab[val0 & EXPTAB_MASK] * EXPPOLY( x0 );
870 x1 = buf[1].d * icvExpTab[val1 & EXPTAB_MASK] * EXPPOLY( x1 );
871
872 y[i] = (float)x0;
873 y[i + 1] = (float)x1;
874
875 x2 = buf[2].d * icvExpTab[val2 & EXPTAB_MASK] * EXPPOLY( x2 );
876 x3 = buf[3].d * icvExpTab[val3 & EXPTAB_MASK] * EXPPOLY( x3 );
877
878 y[i + 2] = (float)x2;
879 y[i + 3] = (float)x3;
880 }
881
882 for( ; i < n; i++ )
883 {
884 double x0 = x[i].f * exp_prescale;
885 int val0, t;
886
887 if( ((x[i].i >> 23) & 255) > 127 + 10 )
888 x0 = x[i].i < 0 ? -exp_max_val : exp_max_val;
889
890 val0 = cvRound(x0);
891 t = (val0 >> EXPTAB_SCALE) + 1023;
892 t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
893
894 buf[0].i.hi = t << 20;
895 x0 = (x0 - val0)*exp_postscale;
896
897 y[i] = (float)(buf[0].d * icvExpTab[val0 & EXPTAB_MASK] * EXPPOLY(x0));
898 }
899
900 return CV_OK;
901 }
902
903
904 IPCVAPI_IMPL( CvStatus, icvExp_64f, ( const double *_x, double *y, int n ), (_x, y, n) )
905 {
906 static const double
907 A5 = .99999999999999999998285227504999 / EXPPOLY_32F_A0,
908 A4 = .69314718055994546743029643825322 / EXPPOLY_32F_A0,
909 A3 = .24022650695886477918181338054308 / EXPPOLY_32F_A0,
910 A2 = .55504108793649567998466049042729e-1 / EXPPOLY_32F_A0,
911 A1 = .96180973140732918010002372686186e-2 / EXPPOLY_32F_A0,
912 A0 = .13369713757180123244806654839424e-2 / EXPPOLY_32F_A0;
913
914 #undef EXPPOLY
915 #define EXPPOLY(x) (((((A0*(x) + A1)*(x) + A2)*(x) + A3)*(x) + A4)*(x) + A5)
916
917 int i = 0;
918 DBLINT buf[4];
919 const Cv64suf* x = (const Cv64suf*)_x;
920
921 if( !x || !y )
922 return CV_NULLPTR_ERR;
923 if( n <= 0 )
924 return CV_BADSIZE_ERR;
925
926 buf[0].i.lo = buf[1].i.lo = buf[2].i.lo = buf[3].i.lo = 0;
927
928 for( ; i <= n - 4; i += 4 )
929 {
930 double x0 = x[i].f * exp_prescale;
931 double x1 = x[i + 1].f * exp_prescale;
932 double x2 = x[i + 2].f * exp_prescale;
933 double x3 = x[i + 3].f * exp_prescale;
934
935 double y0, y1, y2, y3;
936 int val0, val1, val2, val3, t;
937
938 t = (int)(x[i].i >> 52);
939 if( (t & 2047) > 1023 + 10 )
940 x0 = t < 0 ? -exp_max_val : exp_max_val;
941
942 t = (int)(x[i+1].i >> 52);
943 if( (t & 2047) > 1023 + 10 )
944 x1 = t < 0 ? -exp_max_val : exp_max_val;
945
946 t = (int)(x[i+2].i >> 52);
947 if( (t & 2047) > 1023 + 10 )
948 x2 = t < 0 ? -exp_max_val : exp_max_val;
949
950 t = (int)(x[i+3].i >> 52);
951 if( (t & 2047) > 1023 + 10 )
952 x3 = t < 0 ? -exp_max_val : exp_max_val;
953
954 val0 = cvRound(x0);
955 val1 = cvRound(x1);
956 val2 = cvRound(x2);
957 val3 = cvRound(x3);
958
959 x0 = (x0 - val0)*exp_postscale;
960 x1 = (x1 - val1)*exp_postscale;
961 x2 = (x2 - val2)*exp_postscale;
962 x3 = (x3 - val3)*exp_postscale;
963
964 t = (val0 >> EXPTAB_SCALE) + 1023;
965 t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
966 buf[0].i.hi = t << 20;
967
968 t = (val1 >> EXPTAB_SCALE) + 1023;
969 t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
970 buf[1].i.hi = t << 20;
971
972 t = (val2 >> EXPTAB_SCALE) + 1023;
973 t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
974 buf[2].i.hi = t << 20;
975
976 t = (val3 >> EXPTAB_SCALE) + 1023;
977 t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
978 buf[3].i.hi = t << 20;
979
980 y0 = buf[0].d * icvExpTab[val0 & EXPTAB_MASK] * EXPPOLY( x0 );
981 y1 = buf[1].d * icvExpTab[val1 & EXPTAB_MASK] * EXPPOLY( x1 );
982
983 y[i] = y0;
984 y[i + 1] = y1;
985
986 y2 = buf[2].d * icvExpTab[val2 & EXPTAB_MASK] * EXPPOLY( x2 );
987 y3 = buf[3].d * icvExpTab[val3 & EXPTAB_MASK] * EXPPOLY( x3 );
988
989 y[i + 2] = y2;
990 y[i + 3] = y3;
991 }
992
993 for( ; i < n; i++ )
994 {
995 double x0 = x[i].f * exp_prescale;
996 int val0, t;
997
998 t = (int)(x[i].i >> 52);
999 if( (t & 2047) > 1023 + 10 )
1000 x0 = t < 0 ? -exp_max_val : exp_max_val;
1001
1002 val0 = cvRound(x0);
1003 t = (val0 >> EXPTAB_SCALE) + 1023;
1004 t = (t | ((t < 2047) - 1)) & (((t < 0) - 1) & 2047);
1005
1006 buf[0].i.hi = t << 20;
1007 x0 = (x0 - val0)*exp_postscale;
1008
1009 y[i] = buf[0].d * icvExpTab[val0 & EXPTAB_MASK] * EXPPOLY( x0 );
1010 }
1011
1012 return CV_OK;
1013 }
1014
1015 #undef EXPTAB_SCALE
1016 #undef EXPTAB_MASK
1017 #undef EXPPOLY_32F_A0
1018
cvExp(const CvArr * srcarr,CvArr * dstarr)1019 CV_IMPL void cvExp( const CvArr* srcarr, CvArr* dstarr )
1020 {
1021 CV_FUNCNAME( "cvExp" );
1022
1023 __BEGIN__;
1024
1025 CvMat srcstub, *src = (CvMat*)srcarr;
1026 CvMat dststub, *dst = (CvMat*)dstarr;
1027 int coi1 = 0, coi2 = 0, src_depth, dst_depth;
1028 double* buffer = 0;
1029 CvSize size;
1030 int x, y, dx = 0;
1031
1032 if( !CV_IS_MAT(src))
1033 CV_CALL( src = cvGetMat( src, &srcstub, &coi1 ));
1034
1035 if( !CV_IS_MAT(dst))
1036 CV_CALL( dst = cvGetMat( dst, &dststub, &coi2 ));
1037
1038 if( coi1 != 0 || coi2 != 0 )
1039 CV_ERROR( CV_BadCOI, "" );
1040
1041 src_depth = CV_MAT_DEPTH(src->type);
1042 dst_depth = CV_MAT_DEPTH(dst->type);
1043
1044 if( !CV_ARE_CNS_EQ( src, dst ) || src_depth < CV_32F || dst_depth < src_depth )
1045 CV_ERROR_FROM_CODE( CV_StsUnmatchedFormats );
1046
1047 if( !CV_ARE_SIZES_EQ( src, dst ) )
1048 CV_ERROR_FROM_CODE( CV_StsUnmatchedSizes );
1049
1050 size = cvGetMatSize(src);
1051 size.width *= CV_MAT_CN(src->type);
1052
1053 if( CV_IS_MAT_CONT( src->type & dst->type ))
1054 {
1055 size.width *= size.height;
1056 size.height = 1;
1057 }
1058
1059 if( !CV_ARE_DEPTHS_EQ( src, dst ))
1060 {
1061 dx = MIN( 1024, size.width );
1062 buffer = (double*)cvStackAlloc( dx*sizeof(buffer[0]) );
1063 }
1064
1065 for( y = 0; y < size.height; y++ )
1066 {
1067 uchar* src_data = src->data.ptr + src->step*y;
1068 uchar* dst_data = dst->data.ptr + dst->step*y;
1069
1070 if( src_depth == CV_64F )
1071 {
1072 icvExp_64f( (double*)src_data, (double*)dst_data, size.width );
1073 }
1074 else if( src_depth == dst_depth )
1075 {
1076 icvExp_32f( (float*)src_data, (float*)dst_data, size.width );
1077 }
1078 else
1079 {
1080 for( x = 0; x < size.width; x += dx )
1081 {
1082 int len = dx;
1083 if( x + len > size.width )
1084 len = size.width - x;
1085 icvCvt_32f64f( (float*)src_data + x, buffer, len );
1086 icvExp_64f( buffer, (double*)dst_data + x, len );
1087 }
1088 }
1089 }
1090
1091 __END__;
1092 }
1093
1094
1095 /****************************************************************************************\
1096 * L O G *
1097 \****************************************************************************************/
1098
1099 #define LOGTAB_SCALE 8
1100 #define LOGTAB_MASK ((1 << LOGTAB_SCALE) - 1)
1101 #define LOGTAB_MASK2 ((1 << (20 - LOGTAB_SCALE)) - 1)
1102 #define LOGTAB_MASK2_32F ((1 << (23 - LOGTAB_SCALE)) - 1)
1103
1104 static const double icvLogTab[] = {
1105 0.0000000000000000000000000000000000000000, 1.000000000000000000000000000000000000000,
1106 .00389864041565732288852075271279318258166, .9961089494163424124513618677042801556420,
1107 .00778214044205494809292034119607706088573, .9922480620155038759689922480620155038760,
1108 .01165061721997527263705585198749759001657, .9884169884169884169884169884169884169884,
1109 .01550418653596525274396267235488267033361, .9846153846153846153846153846153846153846,
1110 .01934296284313093139406447562578250654042, .9808429118773946360153256704980842911877,
1111 .02316705928153437593630670221500622574241, .9770992366412213740458015267175572519084,
1112 .02697658769820207233514075539915211265906, .9733840304182509505703422053231939163498,
1113 .03077165866675368732785500469617545604706, .9696969696969696969696969696969696969697,
1114 .03455238150665972812758397481047722976656, .9660377358490566037735849056603773584906,
1115 .03831886430213659461285757856785494368522, .9624060150375939849624060150375939849624,
1116 .04207121392068705056921373852674150839447, .9588014981273408239700374531835205992509,
1117 .04580953603129420126371940114040626212953, .9552238805970149253731343283582089552239,
1118 .04953393512227662748292900118940451648088, .9516728624535315985130111524163568773234,
1119 .05324451451881227759255210685296333394944, .9481481481481481481481481481481481481481,
1120 .05694137640013842427411105973078520037234, .9446494464944649446494464944649446494465,
1121 .06062462181643483993820353816772694699466, .9411764705882352941176470588235294117647,
1122 .06429435070539725460836422143984236754475, .9377289377289377289377289377289377289377,
1123 .06795066190850773679699159401934593915938, .9343065693430656934306569343065693430657,
1124 .07159365318700880442825962290953611955044, .9309090909090909090909090909090909090909,
1125 .07522342123758751775142172846244648098944, .9275362318840579710144927536231884057971,
1126 .07884006170777602129362549021607264876369, .9241877256317689530685920577617328519856,
1127 .08244366921107458556772229485432035289706, .9208633093525179856115107913669064748201,
1128 .08603433734180314373940490213499288074675, .9175627240143369175627240143369175627240,
1129 .08961215868968712416897659522874164395031, .9142857142857142857142857142857142857143,
1130 .09317722485418328259854092721070628613231, .9110320284697508896797153024911032028470,
1131 .09672962645855109897752299730200320482256, .9078014184397163120567375886524822695035,
1132 .10026945316367513738597949668474029749630, .9045936395759717314487632508833922261484,
1133 .10379679368164355934833764649738441221420, .9014084507042253521126760563380281690141,
1134 .10731173578908805021914218968959175981580, .8982456140350877192982456140350877192982,
1135 .11081436634029011301105782649756292812530, .8951048951048951048951048951048951048951,
1136 .11430477128005862852422325204315711744130, .8919860627177700348432055749128919860627,
1137 .11778303565638344185817487641543266363440, .8888888888888888888888888888888888888889,
1138 .12124924363286967987640707633545389398930, .8858131487889273356401384083044982698962,
1139 .12470347850095722663787967121606925502420, .8827586206896551724137931034482758620690,
1140 .12814582269193003360996385708858724683530, .8797250859106529209621993127147766323024,
1141 .13157635778871926146571524895989568904040, .8767123287671232876712328767123287671233,
1142 .13499516453750481925766280255629681050780, .8737201365187713310580204778156996587031,
1143 .13840232285911913123754857224412262439730, .8707482993197278911564625850340136054422,
1144 .14179791186025733629172407290752744302150, .8677966101694915254237288135593220338983,
1145 .14518200984449788903951628071808954700830, .8648648648648648648648648648648648648649,
1146 .14855469432313711530824207329715136438610, .8619528619528619528619528619528619528620,
1147 .15191604202584196858794030049466527998450, .8590604026845637583892617449664429530201,
1148 .15526612891112392955683674244937719777230, .8561872909698996655518394648829431438127,
1149 .15860503017663857283636730244325008243330, .8533333333333333333333333333333333333333,
1150 .16193282026931324346641360989451641216880, .8504983388704318936877076411960132890365,
1151 .16524957289530714521497145597095368430010, .8476821192052980132450331125827814569536,
1152 .16855536102980664403538924034364754334090, .8448844884488448844884488448844884488449,
1153 .17185025692665920060697715143760433420540, .8421052631578947368421052631578947368421,
1154 .17513433212784912385018287750426679849630, .8393442622950819672131147540983606557377,
1155 .17840765747281828179637841458315961062910, .8366013071895424836601307189542483660131,
1156 .18167030310763465639212199675966985523700, .8338762214983713355048859934853420195440,
1157 .18492233849401198964024217730184318497780, .8311688311688311688311688311688311688312,
1158 .18816383241818296356839823602058459073300, .8284789644012944983818770226537216828479,
1159 .19139485299962943898322009772527962923050, .8258064516129032258064516129032258064516,
1160 .19461546769967164038916962454095482826240, .8231511254019292604501607717041800643087,
1161 .19782574332991986754137769821682013571260, .8205128205128205128205128205128205128205,
1162 .20102574606059073203390141770796617493040, .8178913738019169329073482428115015974441,
1163 .20421554142869088876999228432396193966280, .8152866242038216560509554140127388535032,
1164 .20739519434607056602715147164417430758480, .8126984126984126984126984126984126984127,
1165 .21056476910734961416338251183333341032260, .8101265822784810126582278481012658227848,
1166 .21372432939771812687723695489694364368910, .8075709779179810725552050473186119873817,
1167 .21687393830061435506806333251006435602900, .8050314465408805031446540880503144654088,
1168 .22001365830528207823135744547471404075630, .8025078369905956112852664576802507836991,
1169 .22314355131420973710199007200571941211830, .8000000000000000000000000000000000000000,
1170 .22626367865045338145790765338460914790630, .7975077881619937694704049844236760124611,
1171 .22937410106484582006380890106811420992010, .7950310559006211180124223602484472049689,
1172 .23247487874309405442296849741978803649550, .7925696594427244582043343653250773993808,
1173 .23556607131276688371634975283086532726890, .7901234567901234567901234567901234567901,
1174 .23864773785017498464178231643018079921600, .7876923076923076923076923076923076923077,
1175 .24171993688714515924331749374687206000090, .7852760736196319018404907975460122699387,
1176 .24478272641769091566565919038112042471760, .7828746177370030581039755351681957186544,
1177 .24783616390458124145723672882013488560910, .7804878048780487804878048780487804878049,
1178 .25088030628580937353433455427875742316250, .7781155015197568389057750759878419452888,
1179 .25391520998096339667426946107298135757450, .7757575757575757575757575757575757575758,
1180 .25694093089750041913887912414793390780680, .7734138972809667673716012084592145015106,
1181 .25995752443692604627401010475296061486000, .7710843373493975903614457831325301204819,
1182 .26296504550088134477547896494797896593800, .7687687687687687687687687687687687687688,
1183 .26596354849713793599974565040611196309330, .7664670658682634730538922155688622754491,
1184 .26895308734550393836570947314612567424780, .7641791044776119402985074626865671641791,
1185 .27193371548364175804834985683555714786050, .7619047619047619047619047619047619047619,
1186 .27490548587279922676529508862586226314300, .7596439169139465875370919881305637982196,
1187 .27786845100345625159121709657483734190480, .7573964497041420118343195266272189349112,
1188 .28082266290088775395616949026589281857030, .7551622418879056047197640117994100294985,
1189 .28376817313064456316240580235898960381750, .7529411764705882352941176470588235294118,
1190 .28670503280395426282112225635501090437180, .7507331378299120234604105571847507331378,
1191 .28963329258304265634293983566749375313530, .7485380116959064327485380116959064327485,
1192 .29255300268637740579436012922087684273730, .7463556851311953352769679300291545189504,
1193 .29546421289383584252163927885703742504130, .7441860465116279069767441860465116279070,
1194 .29836697255179722709783618483925238251680, .7420289855072463768115942028985507246377,
1195 .30126133057816173455023545102449133992200, .7398843930635838150289017341040462427746,
1196 .30414733546729666446850615102448500692850, .7377521613832853025936599423631123919308,
1197 .30702503529491181888388950937951449304830, .7356321839080459770114942528735632183908,
1198 .30989447772286465854207904158101882785550, .7335243553008595988538681948424068767908,
1199 .31275571000389684739317885942000430077330, .7314285714285714285714285714285714285714,
1200 .31560877898630329552176476681779604405180, .7293447293447293447293447293447293447293,
1201 .31845373111853458869546784626436419785030, .7272727272727272727272727272727272727273,
1202 .32129061245373424782201254856772720813750, .7252124645892351274787535410764872521246,
1203 .32411946865421192853773391107097268104550, .7231638418079096045197740112994350282486,
1204 .32694034499585328257253991068864706903700, .7211267605633802816901408450704225352113,
1205 .32975328637246797969240219572384376078850, .7191011235955056179775280898876404494382,
1206 .33255833730007655635318997155991382896900, .7170868347338935574229691876750700280112,
1207 .33535554192113781191153520921943709254280, .7150837988826815642458100558659217877095,
1208 .33814494400871636381467055798566434532400, .7130919220055710306406685236768802228412,
1209 .34092658697059319283795275623560883104800, .7111111111111111111111111111111111111111,
1210 .34370051385331840121395430287520866841080, .7091412742382271468144044321329639889197,
1211 .34646676734620857063262633346312213689100, .7071823204419889502762430939226519337017,
1212 .34922538978528827602332285096053965389730, .7052341597796143250688705234159779614325,
1213 .35197642315717814209818925519357435405250, .7032967032967032967032967032967032967033,
1214 .35471990910292899856770532096561510115850, .7013698630136986301369863013698630136986,
1215 .35745588892180374385176833129662554711100, .6994535519125683060109289617486338797814,
1216 .36018440357500774995358483465679455548530, .6975476839237057220708446866485013623978,
1217 .36290549368936841911903457003063522279280, .6956521739130434782608695652173913043478,
1218 .36561919956096466943762379742111079394830, .6937669376693766937669376693766937669377,
1219 .36832556115870762614150635272380895912650, .6918918918918918918918918918918918918919,
1220 .37102461812787262962487488948681857436900, .6900269541778975741239892183288409703504,
1221 .37371640979358405898480555151763837784530, .6881720430107526881720430107526881720430,
1222 .37640097516425302659470730759494472295050, .6863270777479892761394101876675603217158,
1223 .37907835293496944251145919224654790014030, .6844919786096256684491978609625668449198,
1224 .38174858149084833769393299007788300514230, .6826666666666666666666666666666666666667,
1225 .38441169891033200034513583887019194662580, .6808510638297872340425531914893617021277,
1226 .38706774296844825844488013899535872042180, .6790450928381962864721485411140583554377,
1227 .38971675114002518602873692543653305619950, .6772486772486772486772486772486772486772,
1228 .39235876060286384303665840889152605086580, .6754617414248021108179419525065963060686,
1229 .39499380824086893770896722344332374632350, .6736842105263157894736842105263157894737,
1230 .39762193064713846624158577469643205404280, .6719160104986876640419947506561679790026,
1231 .40024316412701266276741307592601515352730, .6701570680628272251308900523560209424084,
1232 .40285754470108348090917615991202183067800, .6684073107049608355091383812010443864230,
1233 .40546510810816432934799991016916465014230, .6666666666666666666666666666666666666667,
1234 .40806588980822172674223224930756259709600, .6649350649350649350649350649350649350649,
1235 .41065992498526837639616360320360399782650, .6632124352331606217616580310880829015544,
1236 .41324724855021932601317757871584035456180, .6614987080103359173126614987080103359173,
1237 .41582789514371093497757669865677598863850, .6597938144329896907216494845360824742268,
1238 .41840189913888381489925905043492093682300, .6580976863753213367609254498714652956298,
1239 .42096929464412963239894338585145305842150, .6564102564102564102564102564102564102564,
1240 .42353011550580327293502591601281892508280, .6547314578005115089514066496163682864450,
1241 .42608439531090003260516141381231136620050, .6530612244897959183673469387755102040816,
1242 .42863216738969872610098832410585600882780, .6513994910941475826972010178117048346056,
1243 .43117346481837132143866142541810404509300, .6497461928934010152284263959390862944162,
1244 .43370832042155937902094819946796633303180, .6481012658227848101265822784810126582278,
1245 .43623676677491801667585491486534010618930, .6464646464646464646464646464646464646465,
1246 .43875883620762790027214350629947148263450, .6448362720403022670025188916876574307305,
1247 .44127456080487520440058801796112675219780, .6432160804020100502512562814070351758794,
1248 .44378397241030093089975139264424797147500, .6416040100250626566416040100250626566416,
1249 .44628710262841947420398014401143882423650, .6400000000000000000000000000000000000000,
1250 .44878398282700665555822183705458883196130, .6384039900249376558603491271820448877805,
1251 .45127464413945855836729492693848442286250, .6368159203980099502487562189054726368159,
1252 .45375911746712049854579618113348260521900, .6352357320099255583126550868486352357320,
1253 .45623743348158757315857769754074979573500, .6336633663366336633663366336633663366337,
1254 .45870962262697662081833982483658473938700, .6320987654320987654320987654320987654321,
1255 .46117571512217014895185229761409573256980, .6305418719211822660098522167487684729064,
1256 .46363574096303250549055974261136725544930, .6289926289926289926289926289926289926290,
1257 .46608972992459918316399125615134835243230, .6274509803921568627450980392156862745098,
1258 .46853771156323925639597405279346276074650, .6259168704156479217603911980440097799511,
1259 .47097971521879100631480241645476780831830, .6243902439024390243902439024390243902439,
1260 .47341577001667212165614273544633761048330, .6228710462287104622871046228710462287105,
1261 .47584590486996386493601107758877333253630, .6213592233009708737864077669902912621359,
1262 .47827014848147025860569669930555392056700, .6198547215496368038740920096852300242131,
1263 .48068852934575190261057286988943815231330, .6183574879227053140096618357487922705314,
1264 .48310107575113581113157579238759353756900, .6168674698795180722891566265060240963855,
1265 .48550781578170076890899053978500887751580, .6153846153846153846153846153846153846154,
1266 .48790877731923892879351001283794175833480, .6139088729016786570743405275779376498801,
1267 .49030398804519381705802061333088204264650, .6124401913875598086124401913875598086124,
1268 .49269347544257524607047571407747454941280, .6109785202863961813842482100238663484487,
1269 .49507726679785146739476431321236304938800, .6095238095238095238095238095238095238095,
1270 .49745538920281889838648226032091770321130, .6080760095011876484560570071258907363420,
1271 .49982786955644931126130359189119189977650, .6066350710900473933649289099526066350711,
1272 .50219473456671548383667413872899487614650, .6052009456264775413711583924349881796690,
1273 .50455601075239520092452494282042607665050, .6037735849056603773584905660377358490566,
1274 .50691172444485432801997148999362252652650, .6023529411764705882352941176470588235294,
1275 .50926190178980790257412536448100581765150, .6009389671361502347417840375586854460094,
1276 .51160656874906207391973111953120678663250, .5995316159250585480093676814988290398126,
1277 .51394575110223428282552049495279788970950, .5981308411214953271028037383177570093458,
1278 .51627947444845445623684554448118433356300, .5967365967365967365967365967365967365967,
1279 .51860776420804555186805373523384332656850, .5953488372093023255813953488372093023256,
1280 .52093064562418522900344441950437612831600, .5939675174013921113689095127610208816705,
1281 .52324814376454775732838697877014055848100, .5925925925925925925925925925925925925926,
1282 .52556028352292727401362526507000438869000, .5912240184757505773672055427251732101617,
1283 .52786708962084227803046587723656557500350, .5898617511520737327188940092165898617512,
1284 .53016858660912158374145519701414741575700, .5885057471264367816091954022988505747126,
1285 .53246479886947173376654518506256863474850, .5871559633027522935779816513761467889908,
1286 .53475575061602764748158733709715306758900, .5858123569794050343249427917620137299771,
1287 .53704146589688361856929077475797384977350, .5844748858447488584474885844748858447489,
1288 .53932196859560876944783558428753167390800, .5831435079726651480637813211845102505695,
1289 .54159728243274429804188230264117009937750, .5818181818181818181818181818181818181818,
1290 .54386743096728351609669971367111429572100, .5804988662131519274376417233560090702948,
1291 .54613243759813556721383065450936555862450, .5791855203619909502262443438914027149321,
1292 .54839232556557315767520321969641372561450, .5778781038374717832957110609480812641084,
1293 .55064711795266219063194057525834068655950, .5765765765765765765765765765765765765766,
1294 .55289683768667763352766542084282264113450, .5752808988764044943820224719101123595506,
1295 .55514150754050151093110798683483153581600, .5739910313901345291479820627802690582960,
1296 .55738115013400635344709144192165695130850, .5727069351230425055928411633109619686801,
1297 .55961578793542265941596269840374588966350, .5714285714285714285714285714285714285714,
1298 .56184544326269181269140062795486301183700, .5701559020044543429844097995545657015590,
1299 .56407013828480290218436721261241473257550, .5688888888888888888888888888888888888889,
1300 .56628989502311577464155334382667206227800, .5676274944567627494456762749445676274945,
1301 .56850473535266865532378233183408156037350, .5663716814159292035398230088495575221239,
1302 .57071468100347144680739575051120482385150, .5651214128035320088300220750551876379691,
1303 .57291975356178548306473885531886480748650, .5638766519823788546255506607929515418502,
1304 .57511997447138785144460371157038025558000, .5626373626373626373626373626373626373626,
1305 .57731536503482350219940144597785547375700, .5614035087719298245614035087719298245614,
1306 .57950594641464214795689713355386629700650, .5601750547045951859956236323851203501094,
1307 .58169173963462239562716149521293118596100, .5589519650655021834061135371179039301310,
1308 .58387276558098266665552955601015128195300, .5577342047930283224400871459694989106754,
1309 .58604904500357812846544902640744112432000, .5565217391304347826086956521739130434783,
1310 .58822059851708596855957011939608491957200, .5553145336225596529284164859002169197397,
1311 .59038744660217634674381770309992134571100, .5541125541125541125541125541125541125541,
1312 .59254960960667157898740242671919986605650, .5529157667386609071274298056155507559395,
1313 .59470710774669277576265358220553025603300, .5517241379310344827586206896551724137931,
1314 .59685996110779382384237123915227130055450, .5505376344086021505376344086021505376344,
1315 .59900818964608337768851242799428291618800, .5493562231759656652360515021459227467811,
1316 .60115181318933474940990890900138765573500, .5481798715203426124197002141327623126338,
1317 .60329085143808425240052883964381180703650, .5470085470085470085470085470085470085470,
1318 .60542532396671688843525771517306566238400, .5458422174840085287846481876332622601279,
1319 .60755525022454170969155029524699784815300, .5446808510638297872340425531914893617021,
1320 .60968064953685519036241657886421307921400, .5435244161358811040339702760084925690021,
1321 .61180154110599282990534675263916142284850, .5423728813559322033898305084745762711864,
1322 .61391794401237043121710712512140162289150, .5412262156448202959830866807610993657505,
1323 .61602987721551394351138242200249806046500, .5400843881856540084388185654008438818565,
1324 .61813735955507864705538167982012964785100, .5389473684210526315789473684210526315789,
1325 .62024040975185745772080281312810257077200, .5378151260504201680672268907563025210084,
1326 .62233904640877868441606324267922900617100, .5366876310272536687631027253668763102725,
1327 .62443328801189346144440150965237990021700, .5355648535564853556485355648535564853556,
1328 .62652315293135274476554741340805776417250, .5344467640918580375782881002087682672234,
1329 .62860865942237409420556559780379757285100, .5333333333333333333333333333333333333333,
1330 .63068982562619868570408243613201193511500, .5322245322245322245322245322245322245322,
1331 .63276666957103777644277897707070223987100, .5311203319502074688796680497925311203320,
1332 .63483920917301017716738442686619237065300, .5300207039337474120082815734989648033126,
1333 .63690746223706917739093569252872839570050, .5289256198347107438016528925619834710744,
1334 .63897144645792069983514238629140891134750, .5278350515463917525773195876288659793814,
1335 .64103117942093124081992527862894348800200, .5267489711934156378600823045267489711934,
1336 .64308667860302726193566513757104985415950, .5256673511293634496919917864476386036961,
1337 .64513796137358470073053240412264131009600, .5245901639344262295081967213114754098361,
1338 .64718504499530948859131740391603671014300, .5235173824130879345603271983640081799591,
1339 .64922794662510974195157587018911726772800, .5224489795918367346938775510204081632653,
1340 .65126668331495807251485530287027359008800, .5213849287169042769857433808553971486762,
1341 .65330127201274557080523663898929953575150, .5203252032520325203252032520325203252033,
1342 .65533172956312757406749369692988693714150, .5192697768762677484787018255578093306288,
1343 .65735807270835999727154330685152672231200, .5182186234817813765182186234817813765182,
1344 .65938031808912778153342060249997302889800, .5171717171717171717171717171717171717172,
1345 .66139848224536490484126716182800009846700, .5161290322580645161290322580645161290323,
1346 .66341258161706617713093692145776003599150, .5150905432595573440643863179074446680080,
1347 .66542263254509037562201001492212526500250, .5140562248995983935742971887550200803213,
1348 .66742865127195616370414654738851822912700, .5130260521042084168336673346693386773547,
1349 .66943065394262923906154583164607174694550, .5120000000000000000000000000000000000000,
1350 .67142865660530226534774556057527661323550, .5109780439121756487025948103792415169661,
1351 .67342267521216669923234121597488410770900, .5099601593625498007968127490039840637450,
1352 .67541272562017662384192817626171745359900, .5089463220675944333996023856858846918489,
1353 .67739882359180603188519853574689477682100, .5079365079365079365079365079365079365079,
1354 .67938098479579733801614338517538271844400, .5069306930693069306930693069306930693069,
1355 .68135922480790300781450241629499942064300, .5059288537549407114624505928853754940711,
1356 .68333355911162063645036823800182901322850, .5049309664694280078895463510848126232742,
1357 .68530400309891936760919861626462079584600, .5039370078740157480314960629921259842520,
1358 .68727057207096020619019327568821609020250, .5029469548133595284872298624754420432220,
1359 .68923328123880889251040571252815425395950, .5019607843137254901960784313725490196078,
1360 .69314718055994530941723212145818, 5.0e-01,
1361 };
1362
1363
1364
1365 #define LOGTAB_TRANSLATE(x,h) (((x) - 1.)*icvLogTab[(h)+1])
1366 static const double ln_2 = 0.69314718055994530941723212145818;
1367
1368 IPCVAPI_IMPL( CvStatus, icvLog_32f, ( const float *_x, float *y, int n ), (_x, y, n) )
1369 {
1370 static const double shift[] = { 0, -1./512 };
1371 static const double
1372 A0 = 0.3333333333333333333333333,
1373 A1 = -0.5,
1374 A2 = 1;
1375
1376 #undef LOGPOLY
1377 #define LOGPOLY(x,k) ((x)+=shift[k],((A0*(x) + A1)*(x) + A2)*(x))
1378
1379 int i = 0;
1380 union
1381 {
1382 int i;
1383 float f;
1384 }
1385 buf[4];
1386
1387 const int* x = (const int*)_x;
1388
1389 if( !x || !y )
1390 return CV_NULLPTR_ERR;
1391 if( n <= 0 )
1392 return CV_BADSIZE_ERR;
1393
1394 for( i = 0; i <= n - 4; i += 4 )
1395 {
1396 double x0, x1, x2, x3;
1397 double y0, y1, y2, y3;
1398 int h0, h1, h2, h3;
1399
1400 h0 = x[i];
1401 h1 = x[i+1];
1402 buf[0].i = (h0 & LOGTAB_MASK2_32F) | (127 << 23);
1403 buf[1].i = (h1 & LOGTAB_MASK2_32F) | (127 << 23);
1404
1405 y0 = (((h0 >> 23) & 0xff) - 127) * ln_2;
1406 y1 = (((h1 >> 23) & 0xff) - 127) * ln_2;
1407
1408 h0 = (h0 >> (23 - LOGTAB_SCALE - 1)) & LOGTAB_MASK * 2;
1409 h1 = (h1 >> (23 - LOGTAB_SCALE - 1)) & LOGTAB_MASK * 2;
1410
1411 y0 += icvLogTab[h0];
1412 y1 += icvLogTab[h1];
1413
1414 h2 = x[i+2];
1415 h3 = x[i+3];
1416
1417 x0 = LOGTAB_TRANSLATE( buf[0].f, h0 );
1418 x1 = LOGTAB_TRANSLATE( buf[1].f, h1 );
1419
1420 buf[2].i = (h2 & LOGTAB_MASK2_32F) | (127 << 23);
1421 buf[3].i = (h3 & LOGTAB_MASK2_32F) | (127 << 23);
1422
1423 y2 = (((h2 >> 23) & 0xff) - 127) * ln_2;
1424 y3 = (((h3 >> 23) & 0xff) - 127) * ln_2;
1425
1426 h2 = (h2 >> (23 - LOGTAB_SCALE - 1)) & LOGTAB_MASK * 2;
1427 h3 = (h3 >> (23 - LOGTAB_SCALE - 1)) & LOGTAB_MASK * 2;
1428
1429 y2 += icvLogTab[h2];
1430 y3 += icvLogTab[h3];
1431
1432 x2 = LOGTAB_TRANSLATE( buf[2].f, h2 );
1433 x3 = LOGTAB_TRANSLATE( buf[3].f, h3 );
1434
1435 y0 += LOGPOLY( x0, h0 == 510 );
1436 y1 += LOGPOLY( x1, h1 == 510 );
1437
1438 y[i] = (float) y0;
1439 y[i + 1] = (float) y1;
1440
1441 y2 += LOGPOLY( x2, h2 == 510 );
1442 y3 += LOGPOLY( x3, h3 == 510 );
1443
1444 y[i + 2] = (float) y2;
1445 y[i + 3] = (float) y3;
1446 }
1447
1448 for( ; i < n; i++ )
1449 {
1450 int h0 = x[i];
1451 double x0, y0;
1452
1453 y0 = (((h0 >> 23) & 0xff) - 127) * ln_2;
1454
1455 buf[0].i = (h0 & LOGTAB_MASK2_32F) | (127 << 23);
1456 h0 = (h0 >> (23 - LOGTAB_SCALE - 1)) & LOGTAB_MASK * 2;
1457
1458 y0 += icvLogTab[h0];
1459 x0 = LOGTAB_TRANSLATE( buf[0].f, h0 );
1460 y0 += LOGPOLY( x0, h0 == 510 );
1461
1462 y[i] = (float)y0;
1463 }
1464
1465 return CV_OK;
1466 }
1467
1468
1469 IPCVAPI_IMPL( CvStatus, icvLog_64f, ( const double *x, double *y, int n ), (x, y, n) )
1470 {
1471 static const double shift[] = { 0, -1./512 };
1472 static const double
1473 A0 = -.1666666666666666666666666666666666666666,
1474 A1 = +0.2,
1475 A2 = -0.25,
1476 A3 = +0.3333333333333333333333333333333333333333,
1477 A4 = -0.5,
1478 A5 = +1.0;
1479
1480 #undef LOGPOLY
1481 #define LOGPOLY(x,k) ((x)+=shift[k], (xq) = (x)*(x),\
1482 ((A0*(xq) + A2)*(xq) + A4)*(xq) + ((A1*(xq) + A3)*(xq) + A5)*(x))
1483
1484 int i = 0;
1485 DBLINT buf[4];
1486 DBLINT *X = (DBLINT *) x;
1487
1488 if( !x || !y )
1489 return CV_NULLPTR_ERR;
1490 if( n <= 0 )
1491 return CV_BADSIZE_ERR;
1492
1493 for( ; i <= n - 4; i += 4 )
1494 {
1495 double xq;
1496 double x0, x1, x2, x3;
1497 double y0, y1, y2, y3;
1498 int h0, h1, h2, h3;
1499
1500 h0 = X[i].i.lo;
1501 h1 = X[i + 1].i.lo;
1502 buf[0].i.lo = h0;
1503 buf[1].i.lo = h1;
1504
1505 h0 = X[i].i.hi;
1506 h1 = X[i + 1].i.hi;
1507 buf[0].i.hi = (h0 & LOGTAB_MASK2) | (1023 << 20);
1508 buf[1].i.hi = (h1 & LOGTAB_MASK2) | (1023 << 20);
1509
1510 y0 = (((h0 >> 20) & 0x7ff) - 1023) * ln_2;
1511 y1 = (((h1 >> 20) & 0x7ff) - 1023) * ln_2;
1512
1513 h2 = X[i + 2].i.lo;
1514 h3 = X[i + 3].i.lo;
1515 buf[2].i.lo = h2;
1516 buf[3].i.lo = h3;
1517
1518 h0 = (h0 >> (20 - LOGTAB_SCALE - 1)) & LOGTAB_MASK * 2;
1519 h1 = (h1 >> (20 - LOGTAB_SCALE - 1)) & LOGTAB_MASK * 2;
1520
1521 y0 += icvLogTab[h0];
1522 y1 += icvLogTab[h1];
1523
1524 h2 = X[i + 2].i.hi;
1525 h3 = X[i + 3].i.hi;
1526
1527 x0 = LOGTAB_TRANSLATE( buf[0].d, h0 );
1528 x1 = LOGTAB_TRANSLATE( buf[1].d, h1 );
1529
1530 buf[2].i.hi = (h2 & LOGTAB_MASK2) | (1023 << 20);
1531 buf[3].i.hi = (h3 & LOGTAB_MASK2) | (1023 << 20);
1532
1533 y2 = (((h2 >> 20) & 0x7ff) - 1023) * ln_2;
1534 y3 = (((h3 >> 20) & 0x7ff) - 1023) * ln_2;
1535
1536 h2 = (h2 >> (20 - LOGTAB_SCALE - 1)) & LOGTAB_MASK * 2;
1537 h3 = (h3 >> (20 - LOGTAB_SCALE - 1)) & LOGTAB_MASK * 2;
1538
1539 y2 += icvLogTab[h2];
1540 y3 += icvLogTab[h3];
1541
1542 x2 = LOGTAB_TRANSLATE( buf[2].d, h2 );
1543 x3 = LOGTAB_TRANSLATE( buf[3].d, h3 );
1544
1545 y0 += LOGPOLY( x0, h0 == 510 );
1546 y1 += LOGPOLY( x1, h1 == 510 );
1547
1548 y[i] = y0;
1549 y[i + 1] = y1;
1550
1551 y2 += LOGPOLY( x2, h2 == 510 );
1552 y3 += LOGPOLY( x3, h3 == 510 );
1553
1554 y[i + 2] = y2;
1555 y[i + 3] = y3;
1556 }
1557
1558 for( ; i < n; i++ )
1559 {
1560 int h0 = X[i].i.hi;
1561 double xq;
1562 double x0, y0 = (((h0 >> 20) & 0x7ff) - 1023) * ln_2;
1563
1564 buf[0].i.hi = (h0 & LOGTAB_MASK2) | (1023 << 20);
1565 buf[0].i.lo = X[i].i.lo;
1566 h0 = (h0 >> (20 - LOGTAB_SCALE - 1)) & LOGTAB_MASK * 2;
1567
1568 y0 += icvLogTab[h0];
1569 x0 = LOGTAB_TRANSLATE( buf[0].d, h0 );
1570 y0 += LOGPOLY( x0, h0 == 510 );
1571 y[i] = y0;
1572 }
1573
1574 return CV_OK;
1575 }
1576
1577
cvLog(const CvArr * srcarr,CvArr * dstarr)1578 CV_IMPL void cvLog( const CvArr* srcarr, CvArr* dstarr )
1579 {
1580 CV_FUNCNAME( "cvLog" );
1581
1582 __BEGIN__;
1583
1584 CvMat srcstub, *src = (CvMat*)srcarr;
1585 CvMat dststub, *dst = (CvMat*)dstarr;
1586 int coi1 = 0, coi2 = 0, src_depth, dst_depth;
1587 double* buffer = 0;
1588 CvSize size;
1589 int x, y, dx = 0;
1590
1591 if( !CV_IS_MAT(src))
1592 CV_CALL( src = cvGetMat( src, &srcstub, &coi1 ));
1593
1594 if( !CV_IS_MAT(dst))
1595 CV_CALL( dst = cvGetMat( dst, &dststub, &coi2 ));
1596
1597 if( coi1 != 0 || coi2 != 0 )
1598 CV_ERROR( CV_BadCOI, "" );
1599
1600 src_depth = CV_MAT_DEPTH(src->type);
1601 dst_depth = CV_MAT_DEPTH(dst->type);
1602
1603 if( !CV_ARE_CNS_EQ( src, dst ) || dst_depth < CV_32F || src_depth < dst_depth )
1604 CV_ERROR_FROM_CODE( CV_StsUnmatchedFormats );
1605
1606 if( !CV_ARE_SIZES_EQ( src, dst ) )
1607 CV_ERROR_FROM_CODE( CV_StsUnmatchedSizes );
1608
1609 size = cvGetMatSize(src);
1610 size.width *= CV_MAT_CN(src->type);
1611
1612 if( CV_IS_MAT_CONT( src->type & dst->type ))
1613 {
1614 size.width *= size.height;
1615 size.height = 1;
1616 }
1617
1618 if( !CV_ARE_DEPTHS_EQ( src, dst ))
1619 {
1620 dx = MIN( 1024, size.width );
1621 buffer = (double*)cvStackAlloc( dx*sizeof(buffer[0]) );
1622 }
1623
1624 for( y = 0; y < size.height; y++ )
1625 {
1626 uchar* src_data = src->data.ptr + src->step*y;
1627 uchar* dst_data = dst->data.ptr + dst->step*y;
1628
1629 if( dst_depth == CV_64F )
1630 {
1631 icvLog_64f( (double*)src_data, (double*)dst_data, size.width );
1632 }
1633 else if( src_depth == dst_depth )
1634 {
1635 icvLog_32f( (float*)src_data, (float*)dst_data, size.width );
1636 }
1637 else
1638 {
1639 for( x = 0; x < size.width; x += dx )
1640 {
1641 int len = dx;
1642 if( x + len > size.width )
1643 len = size.width - x;
1644 icvLog_64f( (double*)src_data + x, buffer, len );
1645 icvCvt_64f32f( buffer, (float*)dst_data + x, len );
1646 }
1647 }
1648 }
1649
1650 __END__;
1651 }
1652
1653
1654 /****************************************************************************************\
1655 * P O W E R *
1656 \****************************************************************************************/
1657
1658 #define ICV_DEF_IPOW_OP( flavor, arrtype, worktype, cast_macro ) \
1659 static CvStatus CV_STDCALL \
1660 icvIPow_##flavor( const arrtype* src, arrtype* dst, int len, int power ) \
1661 { \
1662 int i; \
1663 \
1664 for( i = 0; i < len; i++ ) \
1665 { \
1666 worktype a = 1, b = src[i]; \
1667 int p = power; \
1668 while( p > 1 ) \
1669 { \
1670 if( p & 1 ) \
1671 a *= b; \
1672 b *= b; \
1673 p >>= 1; \
1674 } \
1675 \
1676 a *= b; \
1677 dst[i] = cast_macro(a); \
1678 } \
1679 \
1680 return CV_OK; \
1681 }
1682
1683
1684 ICV_DEF_IPOW_OP( 8u, uchar, int, CV_CAST_8U )
1685 ICV_DEF_IPOW_OP( 16u, ushort, int, CV_CAST_16U )
1686 ICV_DEF_IPOW_OP( 16s, short, int, CV_CAST_16S )
1687 ICV_DEF_IPOW_OP( 32s, int, int, CV_CAST_32S )
1688 ICV_DEF_IPOW_OP( 32f, float, double, CV_CAST_32F )
1689 ICV_DEF_IPOW_OP( 64f, double, double, CV_CAST_64F )
1690
1691 #define icvIPow_8s 0
1692
1693 CV_DEF_INIT_FUNC_TAB_1D( IPow )
1694
1695 typedef CvStatus (CV_STDCALL * CvIPowFunc)( const void* src, void* dst, int len, int power );
1696 typedef CvStatus (CV_STDCALL * CvSqrtFunc)( const void* src, void* dst, int len );
1697
cvPow(const CvArr * srcarr,CvArr * dstarr,double power)1698 CV_IMPL void cvPow( const CvArr* srcarr, CvArr* dstarr, double power )
1699 {
1700 static CvFuncTable ipow_tab;
1701 static int inittab = 0;
1702
1703 CV_FUNCNAME( "cvPow" );
1704
1705 __BEGIN__;
1706
1707 void* temp_buffer = 0;
1708 int block_size = 0;
1709 CvMat srcstub, *src = (CvMat*)srcarr;
1710 CvMat dststub, *dst = (CvMat*)dstarr;
1711 int coi1 = 0, coi2 = 0;
1712 int depth;
1713 CvSize size;
1714 int x, y;
1715 int ipower = cvRound( power );
1716 int is_ipower = 0;
1717
1718 if( !CV_IS_MAT(src))
1719 CV_CALL( src = cvGetMat( src, &srcstub, &coi1 ));
1720
1721 if( !CV_IS_MAT(dst))
1722 CV_CALL( dst = cvGetMat( dst, &dststub, &coi2 ));
1723
1724 if( coi1 != 0 || coi2 != 0 )
1725 CV_ERROR( CV_BadCOI, "" );
1726
1727 if( !CV_ARE_TYPES_EQ( src, dst ))
1728 CV_ERROR_FROM_CODE( CV_StsUnmatchedFormats );
1729
1730 if( !CV_ARE_SIZES_EQ( src, dst ) )
1731 CV_ERROR_FROM_CODE( CV_StsUnmatchedSizes );
1732
1733 depth = CV_MAT_DEPTH( src->type );
1734
1735 if( fabs(ipower - power) < DBL_EPSILON )
1736 {
1737 if( !inittab )
1738 {
1739 icvInitIPowTable( &ipow_tab );
1740 inittab = 1;
1741 }
1742
1743 if( ipower < 0 )
1744 {
1745 CV_CALL( cvDiv( 0, src, dst ));
1746
1747 if( ipower == -1 )
1748 EXIT;
1749 ipower = -ipower;
1750 src = dst;
1751 }
1752
1753 switch( ipower )
1754 {
1755 case 0:
1756 cvSet( dst, cvScalarAll(1));
1757 EXIT;
1758 case 1:
1759 cvCopy( src, dst );
1760 EXIT;
1761 case 2:
1762 cvMul( src, src, dst );
1763 EXIT;
1764 default:
1765 is_ipower = 1;
1766 }
1767 }
1768 else if( depth < CV_32F )
1769 CV_ERROR( CV_StsUnsupportedFormat,
1770 "Fractional or negative integer power factor can be used "
1771 "with floating-point types only");
1772
1773 size = cvGetMatSize(src);
1774 size.width *= CV_MAT_CN(src->type);
1775
1776 if( CV_IS_MAT_CONT( src->type & dst->type ))
1777 {
1778 size.width *= size.height;
1779 size.height = 1;
1780 }
1781
1782 if( is_ipower )
1783 {
1784 CvIPowFunc pow_func = (CvIPowFunc)ipow_tab.fn_2d[depth];
1785 if( !pow_func )
1786 CV_ERROR( CV_StsUnsupportedFormat, "The data type is not supported" );
1787
1788 for( y = 0; y < size.height; y++ )
1789 {
1790 uchar* src_data = src->data.ptr + src->step*y;
1791 uchar* dst_data = dst->data.ptr + dst->step*y;
1792
1793 pow_func( src_data, dst_data, size.width, ipower );
1794 }
1795 }
1796 else if( fabs(fabs(power) - 0.5) < DBL_EPSILON )
1797 {
1798 CvSqrtFunc sqrt_func = power < 0 ?
1799 (depth == CV_32F ? (CvSqrtFunc)icvInvSqrt_32f : (CvSqrtFunc)icvInvSqrt_64f) :
1800 (depth == CV_32F ? (CvSqrtFunc)icvSqrt_32f : (CvSqrtFunc)icvSqrt_64f);
1801
1802 for( y = 0; y < size.height; y++ )
1803 {
1804 uchar* src_data = src->data.ptr + src->step*y;
1805 uchar* dst_data = dst->data.ptr + dst->step*y;
1806
1807 sqrt_func( src_data, dst_data, size.width );
1808 }
1809 }
1810 else
1811 {
1812 block_size = MIN( size.width, ICV_MATH_BLOCK_SIZE );
1813 temp_buffer = cvStackAlloc( block_size*CV_ELEM_SIZE(depth) );
1814
1815 for( y = 0; y < size.height; y++ )
1816 {
1817 uchar* src_data = src->data.ptr + src->step*y;
1818 uchar* dst_data = dst->data.ptr + dst->step*y;
1819
1820 for( x = 0; x < size.width; x += block_size )
1821 {
1822 int len = MIN( size.width - x, block_size );
1823 if( depth == CV_32F )
1824 {
1825 icvLog_32f( (float*)src_data + x, (float*)temp_buffer, len );
1826 icvScale_32f( (float*)temp_buffer, (float*)temp_buffer, len, (float)power, 0 );
1827 icvExp_32f( (float*)temp_buffer, (float*)dst_data + x, len );
1828 }
1829 else
1830 {
1831 icvLog_64f( (double*)src_data + x, (double*)temp_buffer, len );
1832 icvScale_64f( (double*)temp_buffer, (double*)temp_buffer, len, power, 0 );
1833 icvExp_64f( (double*)temp_buffer, (double*)dst_data + x, len );
1834 }
1835 }
1836 }
1837 }
1838
1839 __END__;
1840 }
1841
1842
1843 /************************** CheckArray for NaN's, Inf's *********************************/
1844
1845 IPCVAPI_IMPL( CvStatus, icvCheckArray_32f_C1R,
1846 ( const float* src, int srcstep, CvSize size, int flags, double min_val, double max_val ),
1847 (src, srcstep, size, flags, min_val, max_val) )
1848 {
1849 Cv32suf a, b;
1850 int ia, ib;
1851 const int* isrc = (const int*)src;
1852
1853 if( !src )
1854 return CV_NULLPTR_ERR;
1855
1856 if( size.width <= 0 || size.height <= 0 )
1857 return CV_BADSIZE_ERR;
1858
1859 if( flags & CV_CHECK_RANGE )
1860 {
1861 a.f = (float)min_val;
1862 b.f = (float)max_val;
1863 }
1864 else
1865 {
1866 a.f = -FLT_MAX;
1867 b.f = FLT_MAX;
1868 }
1869
1870 ia = CV_TOGGLE_FLT(a.i);
1871 ib = CV_TOGGLE_FLT(b.i);
1872
1873 srcstep /= sizeof(isrc[0]);
1874 for( ; size.height--; isrc += srcstep )
1875 {
1876 int i;
1877 for( i = 0; i < size.width; i++ )
1878 {
1879 int val = isrc[i];
1880 val = CV_TOGGLE_FLT(val);
1881
1882 if( val < ia || val >= ib )
1883 return CV_BADRANGE_ERR;
1884 }
1885 }
1886
1887 return CV_OK;
1888 }
1889
1890
1891 IPCVAPI_IMPL( CvStatus, icvCheckArray_64f_C1R,
1892 ( const double* src, int srcstep, CvSize size, int flags, double min_val, double max_val ),
1893 (src, srcstep, size, flags, min_val, max_val) )
1894 {
1895 Cv64suf a, b;
1896 int64 ia, ib;
1897 const int64* isrc = (const int64*)src;
1898
1899 if( !src )
1900 return CV_NULLPTR_ERR;
1901
1902 if( size.width <= 0 || size.height <= 0 )
1903 return CV_BADSIZE_ERR;
1904
1905 if( flags & CV_CHECK_RANGE )
1906 {
1907 a.f = min_val;
1908 b.f = max_val;
1909 }
1910 else
1911 {
1912 a.f = -DBL_MAX;
1913 b.f = DBL_MAX;
1914 }
1915
1916 ia = CV_TOGGLE_DBL(a.i);
1917 ib = CV_TOGGLE_DBL(b.i);
1918
1919 srcstep /= sizeof(isrc[0]);
1920 for( ; size.height--; isrc += srcstep )
1921 {
1922 int i;
1923 for( i = 0; i < size.width; i++ )
1924 {
1925 int64 val = isrc[i];
1926 val = CV_TOGGLE_DBL(val);
1927
1928 if( val < ia || val >= ib )
1929 return CV_BADRANGE_ERR;
1930 }
1931 }
1932
1933 return CV_OK;
1934 }
1935
1936
cvCheckArr(const CvArr * arr,int flags,double minVal,double maxVal)1937 CV_IMPL int cvCheckArr( const CvArr* arr, int flags,
1938 double minVal, double maxVal )
1939 {
1940 int result = 0;
1941
1942 CV_FUNCNAME( "cvCheckArr" );
1943
1944 __BEGIN__;
1945
1946 if( arr )
1947 {
1948 CvStatus status = CV_OK;
1949 CvMat stub, *mat = (CvMat*)arr;
1950 int type;
1951 CvSize size;
1952
1953 if( !CV_IS_MAT( mat ))
1954 CV_CALL( mat = cvGetMat( mat, &stub, 0, 1 ));
1955
1956 type = CV_MAT_TYPE( mat->type );
1957 size = cvGetMatSize( mat );
1958
1959 size.width *= CV_MAT_CN( type );
1960
1961 if( CV_IS_MAT_CONT( mat->type ))
1962 {
1963 size.width *= size.height;
1964 size.height = 1;
1965 }
1966
1967 if( CV_MAT_DEPTH(type) == CV_32F )
1968 {
1969 status = icvCheckArray_32f_C1R( mat->data.fl, mat->step, size,
1970 flags, minVal, maxVal );
1971 }
1972 else if( CV_MAT_DEPTH(type) == CV_64F )
1973 {
1974 status = icvCheckArray_64f_C1R( mat->data.db, mat->step, size,
1975 flags, minVal, maxVal );
1976 }
1977 else
1978 {
1979 CV_ERROR( CV_StsUnsupportedFormat, "" );
1980 }
1981
1982 if( status < 0 )
1983 {
1984 if( status != CV_BADRANGE_ERR || !(flags & CV_CHECK_QUIET))
1985 CV_ERROR( CV_StsOutOfRange, "CheckArray failed" );
1986 EXIT;
1987 }
1988 }
1989
1990 result = 1;
1991
1992 __END__;
1993
1994 return result;
1995 }
1996
1997
1998 /* End of file. */
1999