Changeset 1085
- Timestamp:
- Jun 24, 2004, 2:31:12 PM (22 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/psLib/src/dataManip/psMatrixVectorArithmetic.c
r1082 r1085 29 29 * @author Ross Harman, MHPCC 30 30 * 31 * @version $Revision: 1. 5$ $Name: not supported by cvs2svn $32 * @date $Date: 2004-06-2 4 18:53:42 $31 * @version $Revision: 1.6 $ $Name: not supported by cvs2svn $ 32 * @date $Date: 2004-06-25 00:31:12 $ 33 33 * 34 34 * Copyright 2004 Maui High Performance Computing Center, University of Hawaii … … 76 76 /*****************************************************************************/ 77 77 78 // Local version of min function (cast to double for complex numbers) 78 79 #define MIN(A,B)(((double)(A)<(double)(B))?(A):(B)) 79 80 81 // Local version of max function (cast to double for complex numbers) 80 82 #define MAX(A,B)(((double)(A)>(double)(B))?(A):(B)) 81 83 82 83 // VECTOR_XXXX operations 84 // Conversion for degrees to radians 85 #define D2R M_PI/360.0 86 87 // Conversion for radians to degrees 88 #define R2D 360.0/M_PI 89 90 // Binary SCALAR_XXXX operations 91 #define SCALAR_SCALAR(OUT,IN1,OP,IN2,TYPE) \ 92 { \ 93 ps##TYPE *o = NULL; \ 94 ps##TYPE *i1 = NULL; \ 95 ps##TYPE *i2 = NULL; \ 96 o = &((psScalar*)OUT)->data.TYPE; \ 97 i1 = &((psScalar*)IN1)->data.TYPE; \ 98 i2 = &((psScalar*)IN2)->data.TYPE; \ 99 *o = OP; \ 100 } 101 102 #define SCALAR_VECTOR(OUT,IN1,OP,IN2,TYPE) \ 103 { \ 104 int i = 0; \ 105 int npt = 0; \ 106 ps##TYPE *o = NULL; \ 107 ps##TYPE *i1 = NULL; \ 108 ps##TYPE *i2 = NULL; \ 109 npt = ((psVector*)IN1)->n; \ 110 o = ((psVector*)OUT)->data.TYPE; \ 111 i1 = &((psScalar*)IN1)->data.TYPE; \ 112 i2 = ((psVector*)IN2)->data.TYPE; \ 113 for (i=0; i < npt; i++, o++, i2++) { \ 114 *o = OP; \ 115 } \ 116 } 117 118 #define SCALAR_IMAGE(OUT,IN1,OP,IN2,TYPE) \ 119 { \ 120 int i = 0; \ 121 int j = 0; \ 122 int numRows = 0; \ 123 int numCols = 0; \ 124 ps##TYPE *o = NULL; \ 125 ps##TYPE *i1 = NULL; \ 126 ps##TYPE *i2 = NULL; \ 127 numRows = ((psImage*)IN1)->numRows; \ 128 numCols = ((psImage*)IN1)->numCols; \ 129 for(j = 0; j < numCols; j++) { \ 130 o = ((psImage*)OUT)->data.TYPE[j]; \ 131 i1 = &((psScalar*)IN1)->data.TYPE; \ 132 i2 = ((psImage*)IN2)->data.TYPE[j]; \ 133 for(i = 0; i < numRows; i++, o++, i2++) { \ 134 *o = OP; \ 135 } \ 136 } \ 137 } 138 139 // Binary VECTOR_XXXX operations 84 140 #define VECTOR_SCALAR(OUT,IN1,OP,IN2,TYPE) \ 85 141 { \ … … 109 165 n2 = ((psVector*)IN2)->n; \ 110 166 if(n1 != n2) { \ 111 psError(__func__, ": Inconsistent element count: %d vs %d", n1, n2); \112 return OUT; \167 psError(__func__, ": Inconsistent element count: %d vs %d", n1, n2); \ 168 return OUT; \ 113 169 } \ 114 170 o = ((psVector*)OUT)->data.TYPE; \ … … 138 194 if(dim1 == PS_DIMEN_VECTOR) { /* Regular vectors */ \ 139 195 if(n1!=numRows2) { \ 140 psError(__func__, ": Inconsistent element count: %d vs %d", n1, numRows2); \141 return OUT; \196 psError(__func__, ": Inconsistent element count: %d vs %d", n1, numRows2); \ 197 return OUT; \ 142 198 } \ 143 199 \ … … 156 212 } \ 157 213 \ 158 for(j = 0; j < numRows2; j++) { \159 o = ((psImage*)OUT)->data.TYPE[j]; \160 i1 = ((psVector*)IN1)->data.TYPE; \161 i2 = ((psImage*)IN2)->data.TYPE[j]; \162 for(i = 0; i < numCols2; i++, o++, i1++, i2++) { \163 *o = OP; \164 } \165 } \166 } \ 167 } 168 169 // IMAGE_XXXX operations214 for(j = 0; j < numRows2; j++) { \ 215 o = ((psImage*)OUT)->data.TYPE[j]; \ 216 i1 = ((psVector*)IN1)->data.TYPE; \ 217 i2 = ((psImage*)IN2)->data.TYPE[j]; \ 218 for(i = 0; i < numCols2; i++, o++, i1++, i2++) { \ 219 *o = OP; \ 220 } \ 221 } \ 222 } \ 223 } 224 225 // Binary IMAGE_XXXX operations 170 226 #define IMAGE_SCALAR(OUT,IN1,OP,IN2,TYPE) \ 171 227 { \ … … 207 263 if(dim2 == PS_DIMEN_VECTOR) { /* Regular vectors */ \ 208 264 if(n2!=numRows1) { \ 209 psError(__func__, ": Inconsistent element count: %d vs %d", n2, numRows1); \210 return OUT; \265 psError(__func__, ": Inconsistent element count: %d vs %d", n2, numRows1); \ 266 return OUT; \ 211 267 } \ 212 268 \ … … 225 281 } \ 226 282 \ 227 for(j = 0; j < numRows1; j++) { \228 o = ((psImage*)OUT)->data.TYPE[j]; \229 i1 = ((psVector*)IN2)->data.TYPE; \230 i2 = ((psImage*)IN1)->data.TYPE[j]; \231 for(i = 0; i < numCols1; i++, o++, i2++, i1++) { \232 *o = OP; \233 } \234 } \283 for(j = 0; j < numRows1; j++) { \ 284 o = ((psImage*)OUT)->data.TYPE[j]; \ 285 i1 = ((psVector*)IN2)->data.TYPE; \ 286 i2 = ((psImage*)IN1)->data.TYPE[j]; \ 287 for(i = 0; i < numCols1; i++, o++, i2++, i1++) { \ 288 *o = OP; \ 289 } \ 290 } \ 235 291 } \ 236 292 } … … 252 308 numCols2 = ((psImage*)IN2)->numCols; \ 253 309 if(numRows1!=numRows2 || numCols1!=numCols2) { \ 254 psError(__func__, ": Inconsistent element count: numRows: %d vs %d numCols: %d vs %d", numRows1, \255 numRows2, numCols1, numCols2); 256 return OUT; \310 psError(__func__, ": Inconsistent element count: numRows: %d vs %d numCols: %d vs %d", numRows1, \ 311 numRows2, numCols1, numCols2); \ 312 return OUT; \ 257 313 } \ 258 314 for(j = 0; j < numRows1; j++) { \ … … 266 322 } 267 323 268 // SCALAR_XXXX operations 269 #define SCALAR_SCALAR(OUT,IN1,OP,IN2,TYPE) \ 270 { \ 271 ps##TYPE *o = NULL; \ 272 ps##TYPE *i1 = NULL; \ 273 ps##TYPE *i2 = NULL; \ 274 o = &((psScalar*)OUT)->data.TYPE; \ 275 i1 = &((psScalar*)IN1)->data.TYPE; \ 276 i2 = &((psScalar*)IN2)->data.TYPE; \ 277 *o = OP; \ 278 } 279 280 #define SCALAR_VECTOR(OUT,IN1,OP,IN2,TYPE) \ 281 { \ 282 int i = 0; \ 283 int npt = 0; \ 284 ps##TYPE *o = NULL; \ 285 ps##TYPE *i1 = NULL; \ 286 ps##TYPE *i2 = NULL; \ 287 npt = ((psVector*)IN1)->n; \ 288 o = ((psVector*)OUT)->data.TYPE; \ 289 i1 = &((psScalar*)IN1)->data.TYPE; \ 290 i2 = ((psVector*)IN2)->data.TYPE; \ 291 for (i=0; i < npt; i++, o++, i2++) { \ 292 *o = OP; \ 293 } \ 294 } 295 296 #define SCALAR_IMAGE(OUT,IN1,OP,IN2,TYPE) \ 297 { \ 298 int i = 0; \ 299 int j = 0; \ 300 int numRows = 0; \ 301 int numCols = 0; \ 302 ps##TYPE *o = NULL; \ 303 ps##TYPE *i1 = NULL; \ 304 ps##TYPE *i2 = NULL; \ 305 numRows = ((psImage*)IN1)->numRows; \ 306 numCols = ((psImage*)IN1)->numCols; \ 307 for(j = 0; j < numCols; j++) { \ 308 o = ((psImage*)OUT)->data.TYPE[j]; \ 309 i1 = &((psScalar*)IN1)->data.TYPE; \ 310 i2 = ((psImage*)IN2)->data.TYPE[j]; \ 311 for(i = 0; i < numRows; i++, o++, i2++) { \ 312 *o = OP; \ 313 } \ 314 } \ 315 } 316 324 // Preprocessor macro function to create arithmetic function based on input type 317 325 #define LEVEL1(DIM1,DIM2,OUT,IN1,OP,IN2) \ 318 326 switch (IN1->type) { \ … … 357 365 } 358 366 367 // Preprocessor macro function to create arithmetic function operation name 359 368 #define LEVEL2(DIM1,DIM2,OUT,IN1,OP,IN2) \ 360 369 if(!strncmp(OP, "+", 1)) { \ … … 373 382 } \ 374 383 } else if(!strncmp(OP, "min", 3)) { \ 375 LEVEL1(DIM1,DIM2,OUT,IN1,MIN(*i1,*i2),IN2); \384 LEVEL1(DIM1,DIM2,OUT,IN1,MIN(*i1,*i2),IN2); \ 376 385 } else if(!strncmp(OP, "max", 3)) { \ 377 LEVEL1(DIM1,DIM2,OUT,IN1,MAX(*i1,*i2),IN2); \386 LEVEL1(DIM1,DIM2,OUT,IN1,MAX(*i1,*i2),IN2); \ 378 387 } else { \ 379 388 psError(__func__, ": Line %d - Invalid operator: %s", __LINE__, #OP); \ … … 407 416 if(psType2 == NULL) { 408 417 psError(__func__, ": Line %d - Null in2 argument", __LINE__); 418 return out; 419 } 420 421 if(op == NULL) { 422 psError(__func__, ": Line %d - Null op argument", __LINE__); 409 423 return out; 410 424 } … … 465 479 return out; 466 480 } 481 482 // Unary SCALAR operations 483 #define SCALAR(OUT,IN,OP,TYPE) \ 484 { \ 485 ps##TYPE *o = NULL; \ 486 ps##TYPE *i1 = NULL; \ 487 o = &((psScalar*)OUT)->data.TYPE; \ 488 i1 = &((psScalar*)IN)->data.TYPE; \ 489 *o = OP; \ 490 } 491 492 // Unary IMAGE operations 493 #define VECTOR(OUT,IN,OP,TYPE) \ 494 { \ 495 int i = 0; \ 496 int nIn = 0; \ 497 int nOut = 0; \ 498 ps##TYPE *o = NULL; \ 499 ps##TYPE *i1 = NULL; \ 500 nIn = ((psVector*)IN)->n; \ 501 nOut = ((psVector*)OUT)->n; \ 502 if(nIn != nOut) { \ 503 psError(__func__, ": Inconsistent element count: %d vs %d", nIn, nOut); \ 504 return OUT; \ 505 } \ 506 o = ((psVector*)OUT)->data.TYPE; \ 507 i1 = ((psVector*)IN)->data.TYPE; \ 508 for(i = 0; i < nIn; i++, o++, i1++) { \ 509 *o = OP; \ 510 } \ 511 } 512 513 // Unary IMAGE operations 514 #define IMAGE(OUT,IN,OP,TYPE) \ 515 { \ 516 int i = 0; \ 517 int j = 0; \ 518 int numRowsIn = 0; \ 519 int numColsIn = 0; \ 520 int numRowsOut = 0; \ 521 int numColsOut = 0; \ 522 ps##TYPE *o = NULL; \ 523 ps##TYPE *i1 = NULL; \ 524 numRowsIn = ((psImage*)IN)->numRows; \ 525 numColsIn = ((psImage*)IN)->numCols; \ 526 numRowsOut = ((psImage*)OUT)->numRows; \ 527 numColsOut = ((psImage*)OUT)->numCols; \ 528 if(numRowsIn!=numRowsOut || numColsIn!=numColsOut) { \ 529 psError(__func__, ": Inconsistent element count: numRows: %d vs %d numCols: %d vs %d", numRowsIn, \ 530 numRowsOut, numColsIn, numColsOut); \ 531 return OUT; \ 532 } \ 533 for(j = 0; j < numRowsIn; j++) { \ 534 o = ((psImage*)OUT)->data.TYPE[j]; \ 535 i1 = ((psImage*)IN)->data.TYPE[j]; \ 536 for(i = 0; i < numColsIn; i++, o++, i1++) { \ 537 *o = OP; \ 538 } \ 539 } \ 540 } 541 542 543 // Preprocessor macro function to create arithmetic function based on input type 544 #define UNARY_TYPE(DIM,OUT,IN,OP) \ 545 switch (IN->type) { \ 546 case PS_TYPE_S8: \ 547 DIM(OUT,IN,OP,S8); \ 548 break; \ 549 case PS_TYPE_U8: \ 550 DIM(OUT,IN,OP,U8); \ 551 break; \ 552 case PS_TYPE_S16: \ 553 DIM(OUT,IN,OP,S16); \ 554 break; \ 555 case PS_TYPE_U16: \ 556 DIM(OUT,IN,OP,U16); \ 557 break; \ 558 case PS_TYPE_S32: \ 559 DIM(OUT,IN,OP,S32); \ 560 break; \ 561 case PS_TYPE_U32: \ 562 DIM(OUT,IN,OP,U32); \ 563 break; \ 564 case PS_TYPE_S64: \ 565 DIM(OUT,IN,OP,S64); \ 566 break; \ 567 case PS_TYPE_U64: \ 568 DIM(OUT,IN,OP,U64); \ 569 break; \ 570 case PS_TYPE_F32: \ 571 DIM(OUT,IN,OP,F32); \ 572 break; \ 573 case PS_TYPE_F64: \ 574 DIM(OUT,IN,OP,F64); \ 575 break; \ 576 case PS_TYPE_C32: \ 577 DIM(OUT,IN,OP,C32); \ 578 break; \ 579 case PS_TYPE_C64: \ 580 DIM(OUT,IN,OP,C64); \ 581 break; \ 582 default: \ 583 psError(__func__, ": Line %d - Invalid PS_TYPE: %d", __LINE__, IN->type); \ 584 } 585 586 587 // Preprocessor macro function to create arithmetic function operation name 588 #define UNARY_OP(DIM,OUT,IN,OP) \ 589 if(!strncmp(OP, "abs", 3)) { \ 590 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 591 UNARY_TYPE(DIM,OUT,IN,cabs(*i1)); \ 592 } else { \ 593 UNARY_TYPE(DIM,OUT,IN,fabs(*i1)); \ 594 } \ 595 } else if(!strncmp(OP, "exp", 3)) { \ 596 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 597 UNARY_TYPE(DIM,OUT,IN,exp(*i1)); \ 598 } else { \ 599 UNARY_TYPE(DIM,OUT,IN,cexp(*i1)); \ 600 } \ 601 } else if(!strncmp(OP, "ln", 2)) { \ 602 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 603 UNARY_TYPE(DIM,OUT,IN,log(*i1)); \ 604 } else { \ 605 UNARY_TYPE(DIM,OUT,IN,clog(*i1)); \ 606 } \ 607 } else if(!strncmp(OP, "ten", 3)) { \ 608 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 609 UNARY_TYPE(DIM,OUT,IN,exp10(*i1)); \ 610 } else { \ 611 UNARY_TYPE(DIM,OUT,IN,cpow(*i1,10.0)); \ 612 } \ 613 } else if(!strncmp(OP, "log", 3)) { \ 614 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 615 UNARY_TYPE(DIM,OUT,IN,log10(*i1)); \ 616 } else { \ 617 UNARY_TYPE(DIM,OUT,IN,clog10(*i1)); \ 618 } \ 619 } else if(!strncmp(OP, "sin", 3)) { \ 620 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 621 UNARY_TYPE(DIM,OUT,IN,sin(*i1)); \ 622 } else { \ 623 UNARY_TYPE(DIM,OUT,IN,csin(*i1)); \ 624 } \ 625 } else if(!strncmp(OP, "dsin", 4)) { \ 626 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 627 UNARY_TYPE(DIM,OUT,IN,sin(*i1*D2R)); \ 628 } else { \ 629 UNARY_TYPE(DIM,OUT,IN,csin(*i1*D2R)); \ 630 } \ 631 } else if(!strncmp(OP, "cos", 3)) { \ 632 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 633 UNARY_TYPE(DIM,OUT,IN,cos(*i1)); \ 634 } else { \ 635 UNARY_TYPE(DIM,OUT,IN,ccos(*i1)); \ 636 } \ 637 } else if(!strncmp(OP, "dcos", 4)) { \ 638 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 639 UNARY_TYPE(DIM,OUT,IN,cos(*i1*D2R)); \ 640 } else { \ 641 UNARY_TYPE(DIM,OUT,IN,ccos(*i1*D2R)); \ 642 } \ 643 } else if(!strncmp(OP, "tan", 3)) { \ 644 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 645 UNARY_TYPE(DIM,OUT,IN,tan(*i1)); \ 646 } else { \ 647 UNARY_TYPE(DIM,OUT,IN,ctan(*i1)); \ 648 } \ 649 } else if(!strncmp(OP, "dtan", 4)) { \ 650 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 651 UNARY_TYPE(DIM,OUT,IN,tan(*i1*D2R)); \ 652 } else { \ 653 UNARY_TYPE(DIM,OUT,IN,ctan(*i1*D2R)); \ 654 } \ 655 } else if(!strncmp(OP, "asin", 4)) { \ 656 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 657 UNARY_TYPE(DIM,OUT,IN,asin(*i1)); \ 658 } else { \ 659 UNARY_TYPE(DIM,OUT,IN,casin(*i1)); \ 660 } \ 661 } else if(!strncmp(OP, "dasin", 5)) { \ 662 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 663 UNARY_TYPE(DIM,OUT,IN,R2D*asin(*i1)); \ 664 } else { \ 665 UNARY_TYPE(DIM,OUT,IN,R2D*casin(*i1)); \ 666 } \ 667 } else if(!strncmp(OP, "acos", 4)) { \ 668 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 669 UNARY_TYPE(DIM,OUT,IN,acos(*i1)); \ 670 } else { \ 671 UNARY_TYPE(DIM,OUT,IN,cacos(*i1)); \ 672 } \ 673 } else if(!strncmp(OP, "dacos", 5)) { \ 674 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 675 UNARY_TYPE(DIM,OUT,IN,R2D*acos(*i1)); \ 676 } else { \ 677 UNARY_TYPE(DIM,OUT,IN,R2D*cacos(*i1)); \ 678 } \ 679 } else if(!strncmp(OP, "atan", 4)) { \ 680 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 681 UNARY_TYPE(DIM,OUT,IN,atan(*i1)); \ 682 } else { \ 683 UNARY_TYPE(DIM,OUT,IN,catan(*i1)); \ 684 } \ 685 } else if(!strncmp(OP, "datan", 5)) { \ 686 if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ 687 UNARY_TYPE(DIM,OUT,IN,R2D*atan(*i1)); \ 688 } else { \ 689 UNARY_TYPE(DIM,OUT,IN,R2D*catan(*i1)); \ 690 } \ 691 } else { \ 692 psError(__func__, ": Line %d - Invalid operator: %s", __LINE__, #OP); \ 693 } 694 695 void *psUnaryOp(void *out, void *in, char *op) 696 { 697 psDimen dimIn = 0; 698 psDimen dimOut = 0; 699 psElemType elTypeIn = 0; 700 psElemType elTypeOut = 0; 701 psType *psTypeIn = NULL; 702 psType *psTypeOut = NULL; 703 704 psTypeOut = (psType*)out; 705 if(psTypeOut == NULL) { 706 psError(__func__, ": Line %d - Null out argument", __LINE__); 707 return out; 708 } 709 710 psTypeIn = (psType*)in; 711 if(psTypeIn == NULL) { 712 psError(__func__, ": Line %d - Null in argument", __LINE__); 713 return out; 714 } 715 716 if(op == NULL) { 717 psError(__func__, ": Line %d - Null op argument", __LINE__); 718 return out; 719 } 720 721 dimIn = psTypeIn->dimen; 722 dimOut = psTypeOut->dimen; 723 elTypeIn = psTypeIn->type; 724 elTypeOut = psTypeOut->type; 725 726 if(elTypeIn!=elTypeOut) { 727 psError(__func__, ": Line %d - Element types for arguments inconsistent: (%d, %d)", __LINE__, 728 elTypeIn, elTypeOut); 729 return out; 730 } 731 732 if(dimIn==PS_DIMEN_OTHER || dimOut==PS_DIMEN_OTHER) { 733 psError(__func__, ": Line %d - PS_DIMEN_OTHER not allowed for arguments: (%d, %d)", __LINE__, 734 dimIn, dimOut); 735 return out; 736 } 737 738 if(dimIn == PS_DIMEN_SCALAR) { 739 UNARY_OP(SCALAR,out,psTypeIn,op); // scalar 740 } else if(dimIn==PS_DIMEN_VECTOR || dimIn==PS_DIMEN_TRANSV) { 741 UNARY_OP(VECTOR,out,psTypeIn,op); // vector 742 } else if(dimIn == PS_DIMEN_IMAGE) { 743 UNARY_OP(IMAGE,out,psTypeIn,op); // image 744 } else { 745 psError(__func__, ": Line %d - Invalid dimensionality for in arg: %d", __LINE__, dimIn); 746 } 747 748 return out; 749 }
Note:
See TracChangeset
for help on using the changeset viewer.
