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 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 \* ************************************************************************** */ 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 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 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 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 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 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 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 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
