Index: trunk/psLib/src/dataManip/psMatrixVectorArithmetic.c =================================================================== --- trunk/psLib/src/dataManip/psMatrixVectorArithmetic.c (revision 1648) +++ trunk/psLib/src/dataManip/psMatrixVectorArithmetic.c (revision 2021) @@ -1,3 +1,2 @@ - /** @file psMatrixVectorArithmetic.c * @@ -30,15 +29,13 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.27 $ $Name: not supported by cvs2svn $ - * @date $Date: 2004-08-27 23:33:31 $ + * @version $Revision: 1.28 $ $Name: not supported by cvs2svn $ + * @date $Date: 2004-10-08 02:48:01 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii */ -/******************************************************************************/ - -/* INCLUDE FILES */ - -/******************************************************************************/ +/****************************************************************************** + * INCLUDE FILES * + ******************************************************************************/ #include #include @@ -50,42 +47,9 @@ #include "psVector.h" #include "psScalar.h" - -/******************************************************************************/ - -/* DEFINE STATEMENTS */ - -/******************************************************************************/ - -// None - -/******************************************************************************/ - -/* TYPE DEFINITIONS */ - -/******************************************************************************/ - -// None - -/*****************************************************************************/ - -/* GLOBAL VARIABLES */ - -/*****************************************************************************/ - -// None - -/*****************************************************************************/ - -/* FILE STATIC VARIABLES */ - -/*****************************************************************************/ - -// None - -/*****************************************************************************/ - -/* FUNCTION IMPLEMENTATION - LOCAL */ - -/*****************************************************************************/ +#include "psLogMsg.h" + +/***************************************************************************** + * FUNCTION IMPLEMENTATION - LOCAL * + *****************************************************************************/ // Conversion for degrees to radians @@ -95,16 +59,18 @@ #define R2D 57.29577950924861 /* 180.0/PI */ +/* IEEE Standards on floating-point makes this function pointless, really. // Division with NAN checking static complex double psNanDiv(complex double a, complex double b) { complex double out = 0 + 0i; - + out = a / b; if (isnan(creal(out)) || isnan(cimag(out))) { psError(__func__, ": Divide by zero"); } - + return out; } +*/ // Binary SCALAR_XXXX operations @@ -186,5 +152,8 @@ if(n1 != n2) { \ psError(__func__, ": Inconsistent element count: %d vs %d", n1, n2); \ - return OUT; \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ } \ o = ((psVector* )OUT)->data.TYPE; \ @@ -215,5 +184,8 @@ if(n1!=numRows2) { \ psError(__func__, ": Inconsistent element count: %d vs %d", n1, numRows2); \ - return OUT; \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ } \ \ @@ -229,5 +201,8 @@ if(n1!=numCols2) { \ psError(__func__, ": Inconsistent element count: %d vs %d", n1, numCols2); \ - return OUT; \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ } \ \ @@ -284,5 +259,8 @@ if(n2!=numRows1) { \ psError(__func__, ": Inconsistent element count: %d vs %d", n2, numRows1); \ - return OUT; \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ } \ \ @@ -298,5 +276,8 @@ if(n2!=numCols1) { \ psError(__func__, ": Inconsistent element count: %d vs %d", n2, numCols1); \ - return OUT; \ + if (OUT != IN1) { \ + psFree(OUT); \ + } \ + return NULL; \ } \ \ @@ -330,5 +311,8 @@ psError(__func__, ": Inconsistent element count: numRows: %d vs %d numCols: %d vs %d", numRows1, \ numRows2, numCols1, numCols2); \ - return OUT; \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ } \ for(j = 0; j < numRows1; j++) { \ @@ -345,7 +329,28 @@ #define BINARY_TYPE(DIM1,DIM2,OUT,IN1,OP,IN2) \ switch (IN1->type) { \ +case PS_TYPE_U8: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,U8); \ + break; \ +case PS_TYPE_U16: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,U16); \ + break; \ +case PS_TYPE_U32: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,U32); \ + break; \ +case PS_TYPE_U64: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,U64); \ + break; \ +case PS_TYPE_S8: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,S8); \ + break; \ +case PS_TYPE_S16: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,S16); \ + break; \ case PS_TYPE_S32: \ DIM1##_##DIM2(OUT,IN1,OP,IN2,S32); \ break; \ +case PS_TYPE_S64: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,S64); \ + break; \ case PS_TYPE_F32: \ DIM1##_##DIM2(OUT,IN1,OP,IN2,F32); \ @@ -357,14 +362,13 @@ DIM1##_##DIM2(OUT,IN1,OP,IN2,C32); \ break; \ +case PS_TYPE_C64: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,C64); \ + break; \ default: \ -case PS_TYPE_S8: \ -case PS_TYPE_U8: \ -case PS_TYPE_S16: \ -case PS_TYPE_U16: \ -case PS_TYPE_U32: \ -case PS_TYPE_S64: \ -case PS_TYPE_U64: \ -case PS_TYPE_C64: \ psError(__func__, ": Invalid PS_TYPE: %d", IN1->type); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ } @@ -378,5 +382,5 @@ BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 * *i2,IN2); \ } else if(!strncmp(OP, "/", 1)) { \ - BINARY_TYPE(DIM1,DIM2,OUT,IN1,psNanDiv(*i1,*i2),IN2); \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 / *i2,IN2); \ } else if(!strncmp(OP, "^", 1)) { \ if(PS_IS_PSELEMTYPE_COMPLEX(IN1->type)) { \ @@ -388,4 +392,8 @@ if(PS_IS_PSELEMTYPE_COMPLEX(IN1->type)) { \ psError(__func__, ": Minimum operation not supported for complex numbers"); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ } else { \ BINARY_TYPE(DIM1,DIM2,OUT,IN1,fmin(*i1,*i2),IN2); \ @@ -394,4 +402,8 @@ if(PS_IS_PSELEMTYPE_COMPLEX(IN1->type)) { \ psError(__func__, ": Maximum operation not supported for complex numbers"); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ } else { \ BINARY_TYPE(DIM1,DIM2,OUT,IN1,fmax(*i1,*i2),IN2); \ @@ -399,4 +411,8 @@ } else { \ psError(__func__, ": Invalid operation: %s", OP); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ } @@ -405,22 +421,16 @@ psDimen dim1 = 0; psDimen dim2 = 0; - psDimen dimOut = 0; psElemType elType1 = 0; psElemType elType2 = 0; - psElemType elTypeOut = 0; psType* psType1 = NULL; psType* psType2 = NULL; - psType* psTypeOut = NULL; - - psTypeOut = (psType* ) out; - if (psTypeOut == NULL) { - psError(__func__, ": Line %d - Null out argument", __LINE__); - return out; - } psType1 = (psType* ) in1; if (psType1 == NULL) { psError(__func__, ": Line %d - Null in1 argument", __LINE__); - return out; + if (out != in1 && out != in2) { + psFree(out); + } + return NULL; } @@ -428,57 +438,89 @@ if (psType2 == NULL) { psError(__func__, ": Line %d - Null in2 argument", __LINE__); - return out; + if (out != in1 && out != in2) { + psFree(out); + } + return NULL; } if (op == NULL) { psError(__func__, ": Line %d - Null op argument", __LINE__); - return out; + if (out != in1 && out != in2) { + psFree(out); + } + return NULL; } dim1 = psType1->dimen; dim2 = psType2->dimen; - dimOut = psTypeOut->dimen; elType1 = psType1->type; elType2 = psType2->type; - elTypeOut = psTypeOut->type; - - if (elType1 != elType2 && elType1 != elTypeOut) { - psError(__func__, ": Line %d - Element types for arguments inconsistent: (%d, %d, %d)", __LINE__, - elType1, elType2, elTypeOut); - return out; - } - - if (dim1 == PS_DIMEN_OTHER || dim2 == PS_DIMEN_OTHER || dimOut == PS_DIMEN_OTHER) { - psError(__func__, ": Line %d - PS_DIMEN_OTHER not allowed for arguments: (%d, %d, %d)", __LINE__, - dim1, dim2, dimOut); - return out; - } - - if (dim1 == PS_DIMEN_VECTOR) { + + if (elType1 != elType2) { + psError(__func__, ": Line %d - Element types for arguments inconsistent: (%d, %d)", __LINE__, + elType1, elType2); + if (out != in1 && out != in2) { + psFree(out); + } + return NULL; + } + + if (dim1 == PS_DIMEN_OTHER || dim2 == PS_DIMEN_OTHER) { + psError(__func__, ": Line %d - PS_DIMEN_OTHER not allowed for arguments: (%d, %d)", __LINE__, + dim1, dim2); + if (out != in1 && out != in2) { + psFree(out); + } + return NULL; + } + + if (dim1 == PS_DIMEN_VECTOR || dim1 == PS_DIMEN_TRANSV) { if (((psVector* ) in1)->n == 0) { - psError(__func__, ": Line %d - Vector contains zero elements"); + psLogMsg(__func__, PS_LOG_WARN, ": Line %d - Vector contains zero elements"); } } else if (dim1 == PS_DIMEN_IMAGE) { if (((psImage* ) in1)->numCols == 0 || ((psImage* ) in1)->numRows == 0) { - psError(__func__, ": Line %d - Image contains zero length row or cols"); - } - } - - if (dim2 == PS_DIMEN_VECTOR) { + psLogMsg(__func__, PS_LOG_WARN, ": Line %d - Image contains zero length row or cols"); + } + } + + if (dim2 == PS_DIMEN_VECTOR || dim2 == PS_DIMEN_TRANSV) { if (((psVector* ) in2)->n == 0) { - psError(__func__, ": Line %d - Vector contains zero elements"); + psLogMsg(__func__, PS_LOG_WARN, ": Line %d - Vector contains zero elements"); } } else if (dim2 == PS_DIMEN_IMAGE) { if (((psImage* ) in2)->numCols == 0 || ((psImage* ) in2)->numRows == 0) { - psError(__func__, ": Line %d - Image contains zero length row or cols"); + psLogMsg(__func__, PS_LOG_WARN, ": Line %d - Image contains zero length row or cols"); } } if (dim1 == PS_DIMEN_SCALAR) { + if ( out != NULL && ((psType*)out)->dimen != dim2) { + if (out != in1 && out != in2) { + psFree(out); + } + out = NULL; + } if (dim2 == PS_DIMEN_SCALAR) { + if (out == NULL || ((psScalar*)out)->type.type != elType1) { + if (out != in1 && out != in2) { + psFree(out); + } + out = psScalarAlloc(0.0,elType1); + } BINARY_OP(SCALAR, SCALAR, out, psType1, op, psType2); // scalar op scalar } else if (dim2 == PS_DIMEN_VECTOR || dim2 == PS_DIMEN_TRANSV) { + out = psVectorRecycle(out,((psVector*)in2)->n,elType1); + if (out == NULL) { + psError(__func__, "Couldn't create a proper output psVector."); + return NULL; + } BINARY_OP(SCALAR, VECTOR, out, psType1, op, psType2); // scalar op vector } else if (dim2 == PS_DIMEN_IMAGE) { + out = psImageRecycle(out, ((psImage* ) in2)->numCols, ((psImage* ) in2)->numRows,elType1); + if (out == NULL) { + psError(__func__, "Couldn't create a proper output psImage."); + return NULL; + } BINARY_OP(SCALAR, IMAGE, out, psType1, op, psType2); // scalar op image } else { @@ -487,8 +529,23 @@ } else if (dim1 == PS_DIMEN_VECTOR || dim1 == PS_DIMEN_TRANSV) { if (dim2 == PS_DIMEN_SCALAR) { + out = psVectorRecycle(out,((psVector*)in1)->n,elType1); + if (out == NULL) { + psError(__func__, "Couldn't create a proper output psVector."); + return NULL; + } BINARY_OP(VECTOR, SCALAR, out, psType1, op, psType2); // vector op scalar } else if (dim2 == PS_DIMEN_VECTOR || dim2 == PS_DIMEN_TRANSV) { + out = psVectorRecycle(out,((psVector*)in2)->n,elType2); + if (out == NULL) { + psError(__func__, "Couldn't create a proper output psVector."); + return NULL; + } BINARY_OP(VECTOR, VECTOR, out, psType1, op, psType2); // vector op vector } else if (dim2 == PS_DIMEN_IMAGE) { + out = psImageRecycle(out, ((psImage* ) in2)->numCols, ((psImage* ) in2)->numRows, elType2); + if (out == NULL) { + psError(__func__, "Couldn't create a proper output psImage."); + return NULL; + } BINARY_OP(VECTOR, IMAGE, out, psType1, op, psType2); // vector op image } else { @@ -496,4 +553,9 @@ } } else if (dim1 == PS_DIMEN_IMAGE) { + out = psImageRecycle(out, ((psImage*)in1)->numCols, ((psImage*)in1)->numRows, elType1); + if (out == NULL) { + psError(__func__, "Couldn't create a proper output psImage."); + return NULL; + } if (dim2 == PS_DIMEN_SCALAR) { BINARY_OP(IMAGE, SCALAR, out, psType1, op, psType2); // image op scalar @@ -503,8 +565,16 @@ BINARY_OP(IMAGE, IMAGE, out, psType1, op, psType2); // image op image } else { + if (out != in1 && out != in2) { + psFree(out); + } psError(__func__, ": Line %d - Invalid dimensionality for in2 arg: %d", __LINE__, dim2); + return NULL; } } else { + if (out != in1 && out != in2) { + psFree(out); + } psError(__func__, ": Line %d - Invalid dimensionality for in1 arg: %d", __LINE__, dim1); + return NULL; } @@ -544,5 +614,8 @@ if(nIn != nOut) { \ psError(__func__, ": Inconsistent element count: %d vs %d", nIn, nOut); \ - return OUT; \ + if (OUT != IN) { \ + psFree(OUT); \ + } \ + return NULL; \ } \ o = ((psVector* )OUT)->data.TYPE; \ @@ -571,5 +644,8 @@ psError(__func__, ": Inconsistent element count: numRows: %d vs %d numCols: %d vs %d", numRowsIn, \ numRowsOut, numColsIn, numColsOut); \ - return OUT; \ + if (OUT != IN) { \ + psFree(OUT); \ + } \ + return NULL; \ } \ for(j = 0; j < numRowsIn; j++) { \ @@ -597,14 +673,34 @@ DIM(OUT,IN,OP,C32); \ break; \ +case PS_TYPE_S8: \ + DIM(OUT,IN,OP,C32); \ + break; \ +case PS_TYPE_U8: \ + DIM(OUT,IN,OP,C32); \ + break; \ +case PS_TYPE_S16: \ + DIM(OUT,IN,OP,C32); \ + break; \ +case PS_TYPE_U16: \ + DIM(OUT,IN,OP,C32); \ + break; \ +case PS_TYPE_U32: \ + DIM(OUT,IN,OP,C32); \ + break; \ +case PS_TYPE_S64: \ + DIM(OUT,IN,OP,C32); \ + break; \ +case PS_TYPE_U64: \ + DIM(OUT,IN,OP,C32); \ + break; \ +case PS_TYPE_C64: \ + DIM(OUT,IN,OP,C32); \ + break; \ default: \ -case PS_TYPE_S8: \ -case PS_TYPE_U8: \ -case PS_TYPE_S16: \ -case PS_TYPE_U16: \ -case PS_TYPE_U32: \ -case PS_TYPE_S64: \ -case PS_TYPE_U64: \ -case PS_TYPE_C64: \ psError(__func__, ": Invalid PS_TYPE: %d", IN->type); \ + if (OUT != IN) { \ + psFree(OUT); \ + } \ + return NULL; \ } @@ -725,78 +821,83 @@ { psDimen dimIn = 0; - psDimen dimOut = 0; psElemType elTypeIn = 0; - psElemType elTypeOut = 0; psType* psTypeIn = NULL; - psType* psTypeOut = NULL; - - psTypeOut = (psType* ) out; - if (psTypeOut == NULL) { - psError(__func__, ": Line %d - Null out argument", __LINE__); - return out; - } psTypeIn = (psType* ) in; if (psTypeIn == NULL) { psError(__func__, ": Line %d - Null in argument", __LINE__); - return out; + if (out != in) { + psFree(out); + } + return NULL; } if (op == NULL) { psError(__func__, ": Line %d - Null op argument", __LINE__); - return out; + if (out != in) { + psFree(out); + } + return NULL; } dimIn = psTypeIn->dimen; - dimOut = psTypeOut->dimen; elTypeIn = psTypeIn->type; - elTypeOut = psTypeOut->type; - - if (elTypeIn != elTypeOut) { - psError(__func__, ": Line %d - Element types for arguments inconsistent: (%d, %d)", __LINE__, - elTypeIn, elTypeOut); - return out; - } - - if (dimIn != dimOut) { - psError(__func__, ": Line %d - Dimensionality for arguments inconsistent: (%d, %d)", __LINE__, - dimIn, dimOut); - return out; - } - - if (dimIn == PS_DIMEN_OTHER || dimOut == PS_DIMEN_OTHER) { - psError(__func__, ": Line %d - PS_DIMEN_OTHER not allowed for arguments: (%d, %d)", __LINE__, - dimIn, dimOut); - return out; - } - - if (dimIn == PS_DIMEN_VECTOR) { - if (((psVector* ) in)->n == 0) { + + switch (dimIn) { + case PS_DIMEN_SCALAR: + if (out == NULL || + ((psType*)out)->dimen != PS_DIMEN_SCALAR || + ((psScalar*)out)->type.type != elTypeIn) { + psFree(out); + out = psScalarAlloc(0.0,elTypeIn); + } + UNARY_OP(SCALAR, out, psTypeIn, op); // scalar + break; + case PS_DIMEN_VECTOR: + case PS_DIMEN_TRANSV: + if (((psVector*)in)->n == 0) { + if (out != in) { + psFree(out); + } psError(__func__, ": Line %d - Vector contains zero elements"); - } - } else if (dimIn == PS_DIMEN_IMAGE) { + return NULL; + } + + out = psVectorRecycle(out, + ((psVector*)in)->n, + elTypeIn); + if (out == NULL) { + psError(__func__, "Couldn't create a proper output psVector."); + return NULL; + } + + UNARY_OP(VECTOR, out, psTypeIn, op); // vector + break; + case PS_DIMEN_IMAGE: if (((psImage* ) in)->numCols == 0 || ((psImage* ) in)->numRows == 0) { + if (out != in) { + psFree(out); + } psError(__func__, ": Line %d - Image contains zero length row or cols"); - } - } - - if (dimOut == PS_DIMEN_VECTOR) { - if (((psVector* ) out)->n == 0) { - psError(__func__, ": Line %d - Vector contains zero elements"); - } - } else if (dimOut == PS_DIMEN_IMAGE) { - if (((psImage* ) out)->numCols == 0 || ((psImage* ) out)->numRows == 0) { - psError(__func__, ": Line %d - Image contains zero length row or cols"); - } - } - - if (dimIn == PS_DIMEN_SCALAR) { - UNARY_OP(SCALAR, out, psTypeIn, op); // scalar - } else if (dimIn == PS_DIMEN_VECTOR || dimIn == PS_DIMEN_TRANSV) { - UNARY_OP(VECTOR, out, psTypeIn, op); // vector - } else if (dimIn == PS_DIMEN_IMAGE) { + return NULL; + } + + out = psImageRecycle(out, + ((psImage*)in)->numCols, + ((psImage*)in)->numRows, + elTypeIn); + if (out == NULL) { + psError(__func__, "Couldn't create a proper output psImage."); + return NULL; + } + UNARY_OP(IMAGE, out, psTypeIn, op); // image - } else { + break; + default: + if (out != in) { + psFree(out); + } psError(__func__, ": Line %d - Invalid dimensionality for in arg: %d", __LINE__, dimIn); + return NULL; }