Index: trunk/psLib/src/astro/psCoord.c =================================================================== --- trunk/psLib/src/astro/psCoord.c (revision 3171) +++ trunk/psLib/src/astro/psCoord.c (revision 3182) @@ -10,6 +10,6 @@ * @author GLG, MHPCC * -* @version $Revision: 1.49 $ $Name: not supported by cvs2svn $ -* @date $Date: 2005-02-03 00:45:06 $ +* @version $Revision: 1.50 $ $Name: not supported by cvs2svn $ +* @date $Date: 2005-02-10 02:36:41 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -267,9 +267,9 @@ { while (angle < FLT_EPSILON) { - angle+=M_PI; - } - - while (angle >= M_PI) { - angle-=M_PI; + angle+=PS_PI; + } + + while (angle >= PS_PI) { + angle-=PS_PI; } return(angle); @@ -547,5 +547,5 @@ R = PS_SQRT_F32((coord->x * coord->x) + (coord->y * coord->y)); tmp->d = atan2(-coord->y, coord->x); - tmp->r = acos((R * M_PI) / 180.0); + tmp->r = acos((R * PS_PI) / 180.0); } else if (projection->type == PS_PROJ_CAR) { @@ -555,14 +555,14 @@ } else if (projection->type == PS_PROJ_MER) { tmp->d = coord->x; - tmp->r = (2.0 * atan(exp((coord->y * M_PI / 180.0)))) - 180.0; + tmp->r = (2.0 * atan(exp((coord->y * PS_PI / 180.0)))) - 180.0; } else if (projection->type == PS_PROJ_AIT) { - chu1 = (coord->x * M_PI) / 720.0; + chu1 = (coord->x * PS_PI) / 720.0; chu1 *= chu1; - chu2 = (coord->y * M_PI) / 360.0; + chu2 = (coord->y * PS_PI) / 360.0; chu2 *= chu2; chu = PS_SQRT_F32(1.0 - chu1 - chu2); - tmp->d = 2.0 * atan2((2.0 * chu * chu) - 1.0, (coord->x * chu * M_PI) / 360.0); - tmp->r = asin((coord->y * chu * M_PI) / 180.0); + tmp->d = 2.0 * atan2((2.0 * chu * chu) - 1.0, (coord->x * chu * PS_PI) / 360.0); + tmp->r = asin((coord->y * chu * PS_PI) / 180.0); } else if (projection->type == PS_PROJ_PAR) { @@ -640,6 +640,6 @@ // Wrap these to an acceptable range. This assumes that all // angles are in radians. - tmp->r = fmod(tmp->r, 2*M_PI); - tmp->d = fmod(tmp->d, 2*M_PI); + tmp->r = fmod(tmp->r, 2*PS_PI); + tmp->d = fmod(tmp->d, 2*PS_PI); tmp->rErr = 0.0; tmp->dErr = 0.0; @@ -731,7 +731,7 @@ tmp = (psSphere* ) psAlloc(sizeof(psSphere)); tmp->r = position->r + tmpR; - tmp->r = fmod(tmp->r, 2.0*M_PI); + tmp->r = fmod(tmp->r, 2.0*PS_PI); tmp->d = position->d + tmpD; - tmp->d = fmod(tmp->d, 2.0*M_PI); + tmp->d = fmod(tmp->d, 2.0*PS_PI); tmp->rErr = 0.0; tmp->dErr = 0.0; Index: trunk/psLib/src/astro/psTime.c =================================================================== --- trunk/psLib/src/astro/psTime.c (revision 3171) +++ trunk/psLib/src/astro/psTime.c (revision 3182) @@ -10,6 +10,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.50 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-02-03 00:45:06 $ + * @version $Revision: 1.51 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -43,8 +43,8 @@ /** Two times pi with double precision accuracy */ -#define TWOPI (2.0*M_PI) +#define TWOPI (2.0*PS_PI) /** Conversion from radians to degrees */ -#define R2DEG = (180.0/M_PI) +#define R2DEG = (180.0/PS_PI) /** Maximum length of time string */ @@ -548,5 +548,5 @@ t = 2000.0 + (mjd - 51544.03)/365.2422; - dut2ut1 = 0.022*sin(TWOPI*t) - 0.012*cos(TWOPI*t) - 0.006*sin(4.0*M_PI*t) + 0.007*cos(4.0*M_PI*t); + dut2ut1 = 0.022*sin(TWOPI*t) - 0.012*cos(TWOPI*t) - 0.006*sin(4.0*PS_PI*t) + 0.007*cos(4.0*PS_PI*t); result = dut->data.F64[0] + dut->data.F64[1]*(mjd - dut->data.F64[2]) - dut2ut1; @@ -677,8 +677,8 @@ } - // Create output sphere and convert arcsec to radians (i.e. x/60/60*M_PI/180) + // Create output sphere and convert arcsec to radians (i.e. x/60/60*PS_PI/180) output = psAlloc(sizeof(psSphere)); - output->r = x * M_PI / 648000.0; - output->d = y * M_PI / 648000.0; + output->r = x * PS_PI / 648000.0; + output->d = y * PS_PI / 648000.0; return output; Index: trunk/psLib/src/astronomy/astronomy.i =================================================================== --- trunk/psLib/src/astronomy/astronomy.i (revision 3171) +++ trunk/psLib/src/astronomy/astronomy.i (revision 3182) @@ -17,27 +17,60 @@ // add a method to add primitives to metadata %inline %{ -psBool psMetadataAddStr(psMetadata* md, psS32 where, const char* name, const char* comment, const char* value) { +psBool psMetadataAddStr(psMetadata* md, psS32 where, const char* name, +const char* comment, const char* value) { return psMetadataAdd(md,where,name,PS_TYPE_PTR, PS_META_STR,comment,value); } -psBool psMetadataAddF32(psMetadata* md, psS32 where, const char* name, const char* comment, psF32 value) { +psBool psMetadataAddF32(psMetadata* md, + psS32 where, + const char* name, + const char* comment, + psF32 value) +{ return psMetadataAdd(md,where,name,PS_TYPE_F32, PS_META_PRIMITIVE,comment,value); } -psBool psMetadataAddF64(psMetadata* md, psS32 where, const char* name, const char* comment, psF64 value) { + +psBool psMetadataAddF64(psMetadata* md, + psS32 where, + const char* name, + const char* comment, + psF64 value) +{ return psMetadataAdd(md,where,name,PS_TYPE_F64, PS_META_PRIMITIVE,comment,value); } -psBool psMetadataAddS32(psMetadata* md, psS32 where, const char* name, const char* comment, psS32 value) { + +psBool psMetadataAddS32(psMetadata* md, + psS32 where, + const char* name, + const char* comment, + psS32 value) +{ return psMetadataAdd(md,where,name,PS_TYPE_S32, PS_META_PRIMITIVE,comment,value); } -psMetadataItem* psMetadataItemAllocStr(const char* name, const char* comment, const char* value) { +psMetadataItem* psMetadataItemAllocStr(const char* name, + const char* comment, + const char* value) +{ return psMetadataItemAlloc(name, PS_TYPE_PTR, PS_META_STR, comment, value); } -psMetadataItem* psMetadataItemAllocF32(const char* name, const char* comment, psF32 value) { + +psMetadataItem* psMetadataItemAllocF32(const char* name, + const char* comment, + psF32 value) +{ return psMetadataItemAlloc(name, PS_TYPE_F32, PS_META_PRIMITIVE, comment, value); } -psMetadataItem* psMetadataItemAllocF64(const char* name, const char* comment, psF64 value) { + +psMetadataItem* psMetadataItemAllocF64(const char* name, + const char* comment, + psF64 value) +{ return psMetadataItemAlloc(name, PS_TYPE_F64, PS_META_PRIMITIVE, comment, value); } -psMetadataItem* psMetadataItemAllocS32(const char* name, const char* comment, psS32 value) { + +psMetadataItem* psMetadataItemAllocS32(const char* name, + const char* comment, + psS32 value) +{ return psMetadataItemAlloc(name, PS_TYPE_S32, PS_META_PRIMITIVE, comment, value); } Index: trunk/psLib/src/astronomy/psCoord.c =================================================================== --- trunk/psLib/src/astronomy/psCoord.c (revision 3171) +++ trunk/psLib/src/astronomy/psCoord.c (revision 3182) @@ -10,6 +10,6 @@ * @author GLG, MHPCC * -* @version $Revision: 1.49 $ $Name: not supported by cvs2svn $ -* @date $Date: 2005-02-03 00:45:06 $ +* @version $Revision: 1.50 $ $Name: not supported by cvs2svn $ +* @date $Date: 2005-02-10 02:36:41 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -267,9 +267,9 @@ { while (angle < FLT_EPSILON) { - angle+=M_PI; - } - - while (angle >= M_PI) { - angle-=M_PI; + angle+=PS_PI; + } + + while (angle >= PS_PI) { + angle-=PS_PI; } return(angle); @@ -547,5 +547,5 @@ R = PS_SQRT_F32((coord->x * coord->x) + (coord->y * coord->y)); tmp->d = atan2(-coord->y, coord->x); - tmp->r = acos((R * M_PI) / 180.0); + tmp->r = acos((R * PS_PI) / 180.0); } else if (projection->type == PS_PROJ_CAR) { @@ -555,14 +555,14 @@ } else if (projection->type == PS_PROJ_MER) { tmp->d = coord->x; - tmp->r = (2.0 * atan(exp((coord->y * M_PI / 180.0)))) - 180.0; + tmp->r = (2.0 * atan(exp((coord->y * PS_PI / 180.0)))) - 180.0; } else if (projection->type == PS_PROJ_AIT) { - chu1 = (coord->x * M_PI) / 720.0; + chu1 = (coord->x * PS_PI) / 720.0; chu1 *= chu1; - chu2 = (coord->y * M_PI) / 360.0; + chu2 = (coord->y * PS_PI) / 360.0; chu2 *= chu2; chu = PS_SQRT_F32(1.0 - chu1 - chu2); - tmp->d = 2.0 * atan2((2.0 * chu * chu) - 1.0, (coord->x * chu * M_PI) / 360.0); - tmp->r = asin((coord->y * chu * M_PI) / 180.0); + tmp->d = 2.0 * atan2((2.0 * chu * chu) - 1.0, (coord->x * chu * PS_PI) / 360.0); + tmp->r = asin((coord->y * chu * PS_PI) / 180.0); } else if (projection->type == PS_PROJ_PAR) { @@ -640,6 +640,6 @@ // Wrap these to an acceptable range. This assumes that all // angles are in radians. - tmp->r = fmod(tmp->r, 2*M_PI); - tmp->d = fmod(tmp->d, 2*M_PI); + tmp->r = fmod(tmp->r, 2*PS_PI); + tmp->d = fmod(tmp->d, 2*PS_PI); tmp->rErr = 0.0; tmp->dErr = 0.0; @@ -731,7 +731,7 @@ tmp = (psSphere* ) psAlloc(sizeof(psSphere)); tmp->r = position->r + tmpR; - tmp->r = fmod(tmp->r, 2.0*M_PI); + tmp->r = fmod(tmp->r, 2.0*PS_PI); tmp->d = position->d + tmpD; - tmp->d = fmod(tmp->d, 2.0*M_PI); + tmp->d = fmod(tmp->d, 2.0*PS_PI); tmp->rErr = 0.0; tmp->dErr = 0.0; Index: trunk/psLib/src/astronomy/psTime.c =================================================================== --- trunk/psLib/src/astronomy/psTime.c (revision 3171) +++ trunk/psLib/src/astronomy/psTime.c (revision 3182) @@ -10,6 +10,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.50 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-02-03 00:45:06 $ + * @version $Revision: 1.51 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -43,8 +43,8 @@ /** Two times pi with double precision accuracy */ -#define TWOPI (2.0*M_PI) +#define TWOPI (2.0*PS_PI) /** Conversion from radians to degrees */ -#define R2DEG = (180.0/M_PI) +#define R2DEG = (180.0/PS_PI) /** Maximum length of time string */ @@ -548,5 +548,5 @@ t = 2000.0 + (mjd - 51544.03)/365.2422; - dut2ut1 = 0.022*sin(TWOPI*t) - 0.012*cos(TWOPI*t) - 0.006*sin(4.0*M_PI*t) + 0.007*cos(4.0*M_PI*t); + dut2ut1 = 0.022*sin(TWOPI*t) - 0.012*cos(TWOPI*t) - 0.006*sin(4.0*PS_PI*t) + 0.007*cos(4.0*PS_PI*t); result = dut->data.F64[0] + dut->data.F64[1]*(mjd - dut->data.F64[2]) - dut2ut1; @@ -677,8 +677,8 @@ } - // Create output sphere and convert arcsec to radians (i.e. x/60/60*M_PI/180) + // Create output sphere and convert arcsec to radians (i.e. x/60/60*PS_PI/180) output = psAlloc(sizeof(psSphere)); - output->r = x * M_PI / 648000.0; - output->d = y * M_PI / 648000.0; + output->r = x * PS_PI / 648000.0; + output->d = y * PS_PI / 648000.0; return output; Index: trunk/psLib/src/dataManip/Makefile.am =================================================================== --- trunk/psLib/src/dataManip/Makefile.am (revision 3171) +++ trunk/psLib/src/dataManip/Makefile.am (revision 3182) @@ -11,12 +11,6 @@ noinst_LTLIBRARIES = libpslibdataManip.la -libpslibdataManip_la_SOURCES = \ - psMatrixVectorArithmetic.c \ - psStats.c \ - psFunctions.c \ - psMatrix.c \ - psVectorFFT.c \ - psMinimize.c \ - psRandom.c +libpslibdataManip_la_SOURCES = psUnaryOp.c psBinaryOp.c psStats.c \ + psFunctions.c psMatrix.c psVectorFFT.c psMinimize.c psRandom.c BUILT_SOURCES = psDataManipErrors.h @@ -32,5 +26,6 @@ psFunctions.h \ psMatrix.h \ - psMatrixVectorArithmetic.h \ + psBinaryOp.h \ + psUnaryOp.h \ psVectorFFT.h \ psMinimize.h \ Index: trunk/psLib/src/dataManip/psBinaryOp.c =================================================================== --- trunk/psLib/src/dataManip/psBinaryOp.c (revision 3182) +++ trunk/psLib/src/dataManip/psBinaryOp.c (revision 3182) @@ -0,0 +1,582 @@ +/** @file psBinaryOp.c + * + * @brief Provides binary functions for simple matrix and vector element operations. Functions + * include: + * + * Addition (+) + * Subtraction (-) + * Multiplication (*) + * Division (/) + * Power (^) + * Minimum (min) + * Maximum (max) + * Absolute value (abs) + * Exponent (exp) + * Natural Log (ln) + * Power of 10 (ten) + * Log (log) + * Sine (sin or dsin) + * Cosine (cos or dcos) + * Tangent (tan or dtan) + * Arcsine (asin or dasin) + * Arccosine (acos or dacos) + * Arctan (atan or datan) + * + * Currently only vector-vector and image-image binary operations are supported. + * + * @ingroup MatrixArithmetic + * + * @author Ross Harman, MHPCC + * @author Robert DeSonia, MHPCC + * + * @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ + * + * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii + */ + +/****************************************************************************** + * INCLUDE FILES * + ******************************************************************************/ +#include +#include +#include + +#include "psMemory.h" +#include "psError.h" +#include "psImage.h" +#include "psVector.h" +#include "psScalar.h" +#include "psLogMsg.h" +#include "psConstants.h" +#include "psDataManipErrors.h" + +/***************************************************************************** + * FUNCTION IMPLEMENTATION - LOCAL * + *****************************************************************************/ + +// Conversion for degrees to radians +#define D2R 0.01745329252111111 /* PI/180 */ + +// Conversion for radians to degrees +#define R2D 57.29577950924861 /* 180.0/PI */ + +// Binary SCALAR_XXXX operations +#define SCALAR_SCALAR(OUT,IN1,OP,IN2,TYPE) \ +{ \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + o = &((psScalar* )OUT)->data.TYPE; \ + i1 = &((psScalar* )IN1)->data.TYPE; \ + i2 = &((psScalar* )IN2)->data.TYPE; \ + *o = OP; \ +} + +#define SCALAR_VECTOR(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 npt = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + npt = ((psVector* )IN2)->n; \ + o = ((psVector* )OUT)->data.TYPE; \ + i1 = &((psScalar* )IN1)->data.TYPE; \ + i2 = ((psVector* )IN2)->data.TYPE; \ + for (i=0; i < npt; i++, o++, i2++) { \ + *o = OP; \ + } \ +} + +#define SCALAR_IMAGE(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 j = 0; \ + psS32 numRows = 0; \ + psS32 numCols = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + numRows = ((psImage* )IN2)->numRows; \ + numCols = ((psImage* )IN2)->numCols; \ + for(j = 0; j < numCols; j++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = &((psScalar* )IN1)->data.TYPE; \ + i2 = ((psImage* )IN2)->data.TYPE[j]; \ + for(i = 0; i < numRows; i++, o++, i2++) { \ + *o = OP; \ + } \ + } \ +} + +// Binary VECTOR_XXXX operations +#define VECTOR_SCALAR(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 n1 = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + n1 = ((psVector* )IN1)->n; \ + o = ((psVector* )OUT)->data.TYPE; \ + i1 = ((psVector* )IN1)->data.TYPE; \ + i2 = &((psScalar* )IN2)->data.TYPE; \ + for (i=0; i < n1; i++, o++, i1++) { \ + *o = OP; \ + } \ +} + +#define VECTOR_VECTOR(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 n1 = 0; \ + psS32 n2 = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + n1 = ((psVector* )IN1)->n; \ + n2 = ((psVector* )IN2)->n; \ + if(n1 != n2) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ + n1, n2); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + o = ((psVector* )OUT)->data.TYPE; \ + i1 = ((psVector* )IN1)->data.TYPE; \ + i2 = ((psVector* )IN2)->data.TYPE; \ + for (i=0; i < n1; i++, o++, i1++, i2++) { \ + *o = OP; \ + } \ +} + +#define VECTOR_IMAGE(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 j = 0; \ + psS32 n1 = 0; \ + psS32 numRows2 = 0; \ + psS32 numCols2 = 0; \ + psDimen dim1 = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + n1 = ((psVector* )IN1)->n; \ + dim1 = ((psVector* )IN1)->type.dimen; \ + numRows2 = ((psImage* )IN2)->numRows; \ + numCols2 = ((psImage* )IN2)->numCols; \ + \ + if(dim1 == PS_DIMEN_VECTOR) { /* Regular vectors */ \ + if(n1!=numRows2) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ + n1, numRows2); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + \ + i1 = ((psVector* )IN1)->data.TYPE; \ + for(j = 0; j < numRows2; j++, i1++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i2 = ((psImage* )IN2)->data.TYPE[j]; \ + for(i = 0; i < numCols2; i++, o++, i2++) { \ + *o = OP; \ + } \ + } \ + } else { /* Transposed vectors */ \ + if(n1!=numCols2) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ + n1, numCols2); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + \ + for(j = 0; j < numRows2; j++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = ((psVector* )IN1)->data.TYPE; \ + i2 = ((psImage* )IN2)->data.TYPE[j]; \ + for(i = 0; i < numCols2; i++, o++, i1++, i2++) { \ + *o = OP; \ + } \ + } \ + } \ +} + +// Binary IMAGE_XXXX operations +#define IMAGE_SCALAR(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 j = 0; \ + psS32 numRows = 0; \ + psS32 numCols = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + numRows = ((psImage* )IN1)->numRows; \ + numCols = ((psImage* )IN1)->numCols; \ + for(j = 0; j < numRows; j++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = ((psImage* )IN1)->data.TYPE[j]; \ + i2 = &((psScalar* )IN2)->data.TYPE; \ + for(i = 0; i < numCols; i++, o++, i1++) { \ + *o = OP; \ + } \ + } \ +} + +#define IMAGE_VECTOR(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 j = 0; \ + psS32 n2 = 0; \ + psS32 numRows1 = 0; \ + psS32 numCols1 = 0; \ + psDimen dim2 = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + n2 = ((psVector* )IN2)->n; \ + dim2 = ((psVector* )IN2)->type.dimen; \ + numRows1 = ((psImage* )IN1)->numRows; \ + numCols1 = ((psImage* )IN1)->numCols; \ + \ + if(dim2 == PS_DIMEN_VECTOR) { /* Regular vectors */ \ + if(n2!=numRows1) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ + n2, numRows1); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + \ + i2 = ((psVector* )IN2)->data.TYPE; \ + for(j = 0; j < numRows1; j++, i1++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = ((psImage* )IN1)->data.TYPE[j]; \ + for(i = 0; i < numCols1; i++, o++, i1++) { \ + *o = OP; \ + } \ + } \ + } else { /* Transposed vectors */ \ + if(n2!=numCols1) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ + n2, numCols1); \ + if (OUT != IN1) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + \ + for(j = 0; j < numRows1; j++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = ((psVector* )IN2)->data.TYPE; \ + i2 = ((psImage* )IN1)->data.TYPE[j]; \ + for(i = 0; i < numCols1; i++, o++, i2++, i1++) { \ + *o = OP; \ + } \ + } \ + } \ +} + +#define IMAGE_IMAGE(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 j = 0; \ + psS32 numRows1 = 0; \ + psS32 numCols1 = 0; \ + psS32 numRows2 = 0; \ + psS32 numCols2 = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + numRows1 = ((psImage* )IN1)->numRows; \ + numCols1 = ((psImage* )IN1)->numCols; \ + numRows2 = ((psImage* )IN2)->numRows; \ + numCols2 = ((psImage* )IN2)->numCols; \ + if(numRows1!=numRows2 || numCols1!=numCols2) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_IMAGE_SIZE_DIFFERS, \ + numCols1, numRows1, numCols2, numRows2); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + for(j = 0; j < numRows1; j++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = ((psImage* )IN1)->data.TYPE[j]; \ + i2 = ((psImage* )IN2)->data.TYPE[j]; \ + for(i = 0; i < numCols1; i++, o++, i1++, i2++) { \ + *o = OP; \ + } \ + } \ +} + +// Preprocessor macro function to create arithmetic function based on input type +#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); \ + break; \ +case PS_TYPE_F64: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,F64); \ + break; \ +case PS_TYPE_C32: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,C32); \ + break; \ +case PS_TYPE_C64: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,C64); \ + break; \ +default: \ + /* char* strType; \ + PS_TYPE_NAME(strType,IN1->type); \ + psError(PS_ERR_BAD_PARAMETER_TYPE, true, \ + PS_ERRORTEXT_psMatrix_TYPE_MISMATCH, \ + strType); */ \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ +} + +// Preprocessor macro function to create arithmetic function operation name +#define BINARY_OP(DIM1,DIM2,OUT,IN1,OP,IN2) \ +if(!strncmp(OP, "=", 1)) { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i2,IN2); \ +} else if(!strncmp(OP, "+", 1)) { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 + *i2,IN2); \ +} else if(!strncmp(OP, "-", 1)) { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 - *i2,IN2); \ +} else if(!strncmp(OP, "*", 1)) { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 * *i2,IN2); \ +} else if(!strncmp(OP, "/", 1)) { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 / *i2,IN2); \ +} else if(!strncmp(OP, "^", 1)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN1->type)) { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,cpow(*i1,*i2),IN2); \ + } else { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,pow(*i1,*i2),IN2); \ + } \ +} else if(!strncmp(OP, "min", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN1->type)) { \ + psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ + PS_ERRORTEXT_psMatrix_MIN_COMPLEX_SUPPORT); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } else { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,fmin(*i1,*i2),IN2); \ + } \ +} else if(!strncmp(OP, "max", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN1->type)) { \ + psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ + PS_ERRORTEXT_psMatrix_MAX_COMPLEX_SUPPORT); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } else { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,fmax(*i1,*i2),IN2); \ + } \ +} else { \ + psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ + PS_ERRORTEXT_psMatrix_OPERATION_UNSUPPORTED, \ + OP); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ +} + +psPtr psBinaryOp(psPtr out, psPtr in1, char *op, psPtr in2) +{ + + psVector* input1 = (psVector* ) in1; + psVector* input2 = (psVector* ) in2; + + #define psBinaryOp_EXIT { \ + if (out != in1 && out != in2) { \ + psFree(out); \ + } \ + return NULL; \ + } + + PS_PTR_CHECK_NULL_GENERAL(input1, psBinaryOp_EXIT); + PS_PTR_CHECK_NULL_GENERAL(input2, psBinaryOp_EXIT); + PS_PTR_CHECK_NULL_GENERAL(op, psBinaryOp_EXIT); + + PS_PTR_CHECK_TYPE_EQUAL_GENERAL(input1,input2, psBinaryOp_EXIT); + + PS_PTR_CHECK_DIMEN_GENERAL_NOT(input1, PS_DIMEN_OTHER, psBinaryOp_EXIT); + PS_PTR_CHECK_DIMEN_GENERAL_NOT(input2, PS_DIMEN_OTHER, psBinaryOp_EXIT); + + psType* psType1 = (psType*)in1; + psType* psType2 = (psType*)in2; + psDimen dim1 = psType1->dimen; + psDimen dim2 = psType2->dimen; + psElemType elType1 = psType1->type; + psElemType elType2 = psType2->type; + + if (dim1 == PS_DIMEN_VECTOR || dim1 == PS_DIMEN_TRANSV) { + if (((psVector* ) in1)->n == 0) { + psLogMsg(__func__, PS_LOG_WARN, "Vector contains zero elements"); + } + } else if (dim1 == PS_DIMEN_IMAGE) { + if (((psImage* ) in1)->numCols == 0 || ((psImage* ) in1)->numRows == 0) { + psLogMsg(__func__, PS_LOG_WARN, "Image contains zero length row or cols"); + } + } + + if (dim2 == PS_DIMEN_VECTOR || dim2 == PS_DIMEN_TRANSV) { + if (((psVector* ) in2)->n == 0) { + psLogMsg(__func__, PS_LOG_WARN, "Vector contains zero elements"); + } + } else if (dim2 == PS_DIMEN_IMAGE) { + if (((psImage* ) in2)->numCols == 0 || ((psImage* ) in2)->numRows == 0) { + psLogMsg(__func__, PS_LOG_WARN, "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(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_VECTOR_NOT_CREATED); + 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(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_IMAGE_NOT_CREATED); + return NULL; + } + BINARY_OP(SCALAR, IMAGE, out, psType1, op, psType2); // scalar op image + } else { + psError(PS_ERR_BAD_PARAMETER_TYPE, true, + PS_ERRORTEXT_psMatrix_DIMEN_INVALID, + "in2",dim2); + psBinaryOp_EXIT; + } + } 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(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_VECTOR_NOT_CREATED); + 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(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_VECTOR_NOT_CREATED); + 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(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_IMAGE_NOT_CREATED); + return NULL; + } + BINARY_OP(VECTOR, IMAGE, out, psType1, op, psType2); // vector op image + } else { + psError(PS_ERR_BAD_PARAMETER_TYPE, true, + PS_ERRORTEXT_psMatrix_DIMEN_INVALID, + "in2",dim2); + psBinaryOp_EXIT; + } + } else if (dim1 == PS_DIMEN_IMAGE) { + out = psImageRecycle(out, ((psImage*)in1)->numCols, ((psImage*)in1)->numRows, elType1); + if (out == NULL) { + psError(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_IMAGE_NOT_CREATED); + return NULL; + } + if (dim2 == PS_DIMEN_SCALAR) { + BINARY_OP(IMAGE, SCALAR, out, psType1, op, psType2); // image op scalar + } else if (dim2 == PS_DIMEN_VECTOR || dim2 == PS_DIMEN_TRANSV) { + BINARY_OP(IMAGE, VECTOR, out, psType1, op, psType2); // image op vector + } else if (dim2 == PS_DIMEN_IMAGE) { + BINARY_OP(IMAGE, IMAGE, out, psType1, op, psType2); // image op image + } else { + psError(PS_ERR_BAD_PARAMETER_TYPE, true, + PS_ERRORTEXT_psMatrix_DIMEN_INVALID, + "in2",dim2); + psBinaryOp_EXIT; + } + } else { + psError(PS_ERR_BAD_PARAMETER_TYPE, true, + PS_ERRORTEXT_psMatrix_DIMEN_INVALID, + "in1",dim1); + psBinaryOp_EXIT; + } + + // Automtically free psScalar types, since they are usually allocated in the argument list when this + // function is called, provided that the input is not the output. + if(psType1->dimen==PS_DIMEN_SCALAR && in1!=out) { + psFree(in1); + } + + if(psType2->dimen==PS_DIMEN_SCALAR && in2!=out) { + psFree(in2); + } + + return out; +} Index: trunk/psLib/src/dataManip/psBinaryOp.h =================================================================== --- trunk/psLib/src/dataManip/psBinaryOp.h (revision 3182) +++ trunk/psLib/src/dataManip/psBinaryOp.h (revision 3182) @@ -0,0 +1,67 @@ + +/** @file psBinaryOp.h + * + * @brief Provides binary functions for simple matrix and vector element operations. Functions + * include: + * + * Addition (+) + * Subtraction (-) + * Multiplication (*) + * Division (/) + * Power (^) + * Minimum (min) + * Maximum (max) + * Absolute value (abs) + * Exponent (exp) + * Natural Log (ln) + * Power of 10 (ten) + * Log (log) + * Sine (sin or dsin) + * Cosine (cos or dcos) + * Tangent (tan or dtan) + * Arcsine (asin or dasin) + * Arccosine (acos or dacos) + * Arctan (atan or datan) + * + * Currently only vector-vector and image-image binary operations are supported. + * + * @ingroup MatrixArithmetic + * + * @author Ross Harman, MHPCC + * @author Robert DeSonia, MHPCC + * + * @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ + * + * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii + */ + +#ifndef PSBINARY_OP_H +#define PSBINARY_OP_H + +/// @addtogroup MatrixArithmetic +/// @{ + +/** Perform simple binary arithmetic with images or vectors + * + * Performs addition, subtraction, multiplication, division, power, minumum, and maximum arithmetic + * operations with images and vectors. Uses the form: + * + * out = in1 op in2, + * + * Where op is: "=", "+", "-", "*", "/", "^", "min", or "max" + * + * This function only supports vector-vector or image-image operations. + * + * @return psType* : Pointer to either psImage or psVector. + */ +psType* psBinaryOp( + psPtr out, ///< Output type, either psImage or psVector. + psPtr in1, ///< First input, either psImage or psVector. + char *op, ///< Operator. + psPtr in2 ///< Second input, either psImage or psVector. +); + +/// @} + +#endif Index: trunk/psLib/src/dataManip/psConstants.h =================================================================== --- trunk/psLib/src/dataManip/psConstants.h (revision 3171) +++ trunk/psLib/src/dataManip/psConstants.h (revision 3182) @@ -6,6 +6,6 @@ * @author GLG, MHPCC * - * @version $Revision: 1.53 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-01-10 19:47:11 $ + * @version $Revision: 1.54 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -37,10 +37,10 @@ #define PS_COT(X) (1.0 / atan(X)) -#define DEG_TO_RAD(DEGREES) ((DEGREES) * M_PI / 180.0) -#define MIN_TO_RAD(MINUTES) ((MINUTES) * M_PI / (180.0 * 60.0)) -#define SEC_TO_RAD(SECONDS) ((SECONDS) * M_PI / (180.0 * 60.0 * 60.0)) -#define RAD_TO_DEG(RADIANS) ((RADIANS) * 180.0 / M_PI) -#define RAD_TO_MIN(RADIANS) ((RADIANS) * 180.0 * 60.0 / M_PI) -#define RAD_TO_SEC(RADIANS) ((RADIANS) * 180.0 * 60.0 * 60.0 / M_PI) +#define DEG_TO_RAD(DEGREES) ((DEGREES) * PS_PI / 180.0) +#define MIN_TO_RAD(MINUTES) ((MINUTES) * PS_PI / (180.0 * 60.0)) +#define SEC_TO_RAD(SECONDS) ((SECONDS) * PS_PI / (180.0 * 60.0 * 60.0)) +#define RAD_TO_DEG(RADIANS) ((RADIANS) * 180.0 / PS_PI) +#define RAD_TO_MIN(RADIANS) ((RADIANS) * 180.0 * 60.0 / PS_PI) +#define RAD_TO_SEC(RADIANS) ((RADIANS) * 180.0 * 60.0 * 60.0 / PS_PI) /***************************************************************************** Index: trunk/psLib/src/dataManip/psFunctions.c =================================================================== --- trunk/psLib/src/dataManip/psFunctions.c (revision 3171) +++ trunk/psLib/src/dataManip/psFunctions.c (revision 3182) @@ -7,6 +7,6 @@ * polynomials. It also contains a Gaussian functions. * - * @version $Revision: 1.87 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-02-09 20:56:17 $ + * @version $Revision: 1.88 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -975,5 +975,5 @@ if (normal == true) { - tmp = 1.0 / PS_SQRT_F32(2.0 * M_PI * (sigma * sigma)); + tmp = 1.0 / PS_SQRT_F32(2.0 * PS_PI * (sigma * sigma)); } Index: trunk/psLib/src/dataManip/psMatrixVectorArithmetic.c =================================================================== --- trunk/psLib/src/dataManip/psMatrixVectorArithmetic.c (revision 3171) +++ (revision ) @@ -1,905 +1,0 @@ -/** @file psMatrixVectorArithmetic.c - * - * @brief Provides unary and binary functions for simple matrix and vector element operations. Functions - * include: - * - * Addition (+) - * Subtraction (-) - * Multiplication (*) - * Division (/) - * Power (^) - * Minimum (min) - * Maximum (max) - * Absolute value (abs) - * Exponent (exp) - * Natural Log (ln) - * Power of 10 (ten) - * Log (log) - * Sine (sin or dsin) - * Cosine (cos or dcos) - * Tangent (tan or dtan) - * Arcsine (asin or dasin) - * Arccosine (acos or dacos) - * Arctan (atan or datan) - * - * Currently only vector-vector and image-image binary operations are supported. - * - * @ingroup MatrixVectorArithmetic - * - * @author Ross Harman, MHPCC - * - * @version $Revision: 1.34 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-01-12 22:17:01 $ - * - * Copyright 2004 Maui High Performance Computing Center, University of Hawaii - */ - -/****************************************************************************** - * INCLUDE FILES * - ******************************************************************************/ -#include -#include -#include - -#include "psMemory.h" -#include "psError.h" -#include "psImage.h" -#include "psVector.h" -#include "psScalar.h" -#include "psLogMsg.h" -#include "psConstants.h" -#include "psDataManipErrors.h" - -/***************************************************************************** - * FUNCTION IMPLEMENTATION - LOCAL * - *****************************************************************************/ - -// Conversion for degrees to radians -#define D2R 0.01745329252111111 /* PI/180 */ - -// Conversion for radians to degrees -#define R2D 57.29577950924861 /* 180.0/PI */ - -// Binary SCALAR_XXXX operations -#define SCALAR_SCALAR(OUT,IN1,OP,IN2,TYPE) \ -{ \ - ps##TYPE *o = NULL; \ - ps##TYPE *i1 = NULL; \ - ps##TYPE *i2 = NULL; \ - o = &((psScalar* )OUT)->data.TYPE; \ - i1 = &((psScalar* )IN1)->data.TYPE; \ - i2 = &((psScalar* )IN2)->data.TYPE; \ - *o = OP; \ -} - -#define SCALAR_VECTOR(OUT,IN1,OP,IN2,TYPE) \ -{ \ - psS32 i = 0; \ - psS32 npt = 0; \ - ps##TYPE *o = NULL; \ - ps##TYPE *i1 = NULL; \ - ps##TYPE *i2 = NULL; \ - npt = ((psVector* )IN2)->n; \ - o = ((psVector* )OUT)->data.TYPE; \ - i1 = &((psScalar* )IN1)->data.TYPE; \ - i2 = ((psVector* )IN2)->data.TYPE; \ - for (i=0; i < npt; i++, o++, i2++) { \ - *o = OP; \ - } \ -} - -#define SCALAR_IMAGE(OUT,IN1,OP,IN2,TYPE) \ -{ \ - psS32 i = 0; \ - psS32 j = 0; \ - psS32 numRows = 0; \ - psS32 numCols = 0; \ - ps##TYPE *o = NULL; \ - ps##TYPE *i1 = NULL; \ - ps##TYPE *i2 = NULL; \ - numRows = ((psImage* )IN2)->numRows; \ - numCols = ((psImage* )IN2)->numCols; \ - for(j = 0; j < numCols; j++) { \ - o = ((psImage* )OUT)->data.TYPE[j]; \ - i1 = &((psScalar* )IN1)->data.TYPE; \ - i2 = ((psImage* )IN2)->data.TYPE[j]; \ - for(i = 0; i < numRows; i++, o++, i2++) { \ - *o = OP; \ - } \ - } \ -} - -// Binary VECTOR_XXXX operations -#define VECTOR_SCALAR(OUT,IN1,OP,IN2,TYPE) \ -{ \ - psS32 i = 0; \ - psS32 n1 = 0; \ - ps##TYPE *o = NULL; \ - ps##TYPE *i1 = NULL; \ - ps##TYPE *i2 = NULL; \ - n1 = ((psVector* )IN1)->n; \ - o = ((psVector* )OUT)->data.TYPE; \ - i1 = ((psVector* )IN1)->data.TYPE; \ - i2 = &((psScalar* )IN2)->data.TYPE; \ - for (i=0; i < n1; i++, o++, i1++) { \ - *o = OP; \ - } \ -} - -#define VECTOR_VECTOR(OUT,IN1,OP,IN2,TYPE) \ -{ \ - psS32 i = 0; \ - psS32 n1 = 0; \ - psS32 n2 = 0; \ - ps##TYPE *o = NULL; \ - ps##TYPE *i1 = NULL; \ - ps##TYPE *i2 = NULL; \ - n1 = ((psVector* )IN1)->n; \ - n2 = ((psVector* )IN2)->n; \ - if(n1 != n2) { \ - psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ - PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ - n1, n2); \ - if (OUT != IN1 && OUT != IN2) { \ - psFree(OUT); \ - } \ - return NULL; \ - } \ - o = ((psVector* )OUT)->data.TYPE; \ - i1 = ((psVector* )IN1)->data.TYPE; \ - i2 = ((psVector* )IN2)->data.TYPE; \ - for (i=0; i < n1; i++, o++, i1++, i2++) { \ - *o = OP; \ - } \ -} - -#define VECTOR_IMAGE(OUT,IN1,OP,IN2,TYPE) \ -{ \ - psS32 i = 0; \ - psS32 j = 0; \ - psS32 n1 = 0; \ - psS32 numRows2 = 0; \ - psS32 numCols2 = 0; \ - psDimen dim1 = 0; \ - ps##TYPE *o = NULL; \ - ps##TYPE *i1 = NULL; \ - ps##TYPE *i2 = NULL; \ - n1 = ((psVector* )IN1)->n; \ - dim1 = ((psVector* )IN1)->type.dimen; \ - numRows2 = ((psImage* )IN2)->numRows; \ - numCols2 = ((psImage* )IN2)->numCols; \ - \ - if(dim1 == PS_DIMEN_VECTOR) { /* Regular vectors */ \ - if(n1!=numRows2) { \ - psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ - PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ - n1, numRows2); \ - if (OUT != IN1 && OUT != IN2) { \ - psFree(OUT); \ - } \ - return NULL; \ - } \ - \ - i1 = ((psVector* )IN1)->data.TYPE; \ - for(j = 0; j < numRows2; j++, i1++) { \ - o = ((psImage* )OUT)->data.TYPE[j]; \ - i2 = ((psImage* )IN2)->data.TYPE[j]; \ - for(i = 0; i < numCols2; i++, o++, i2++) { \ - *o = OP; \ - } \ - } \ - } else { /* Transposed vectors */ \ - if(n1!=numCols2) { \ - psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ - PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ - n1, numCols2); \ - if (OUT != IN1 && OUT != IN2) { \ - psFree(OUT); \ - } \ - return NULL; \ - } \ - \ - for(j = 0; j < numRows2; j++) { \ - o = ((psImage* )OUT)->data.TYPE[j]; \ - i1 = ((psVector* )IN1)->data.TYPE; \ - i2 = ((psImage* )IN2)->data.TYPE[j]; \ - for(i = 0; i < numCols2; i++, o++, i1++, i2++) { \ - *o = OP; \ - } \ - } \ - } \ -} - -// Binary IMAGE_XXXX operations -#define IMAGE_SCALAR(OUT,IN1,OP,IN2,TYPE) \ -{ \ - psS32 i = 0; \ - psS32 j = 0; \ - psS32 numRows = 0; \ - psS32 numCols = 0; \ - ps##TYPE *o = NULL; \ - ps##TYPE *i1 = NULL; \ - ps##TYPE *i2 = NULL; \ - numRows = ((psImage* )IN1)->numRows; \ - numCols = ((psImage* )IN1)->numCols; \ - for(j = 0; j < numRows; j++) { \ - o = ((psImage* )OUT)->data.TYPE[j]; \ - i1 = ((psImage* )IN1)->data.TYPE[j]; \ - i2 = &((psScalar* )IN2)->data.TYPE; \ - for(i = 0; i < numCols; i++, o++, i1++) { \ - *o = OP; \ - } \ - } \ -} - -#define IMAGE_VECTOR(OUT,IN1,OP,IN2,TYPE) \ -{ \ - psS32 i = 0; \ - psS32 j = 0; \ - psS32 n2 = 0; \ - psS32 numRows1 = 0; \ - psS32 numCols1 = 0; \ - psDimen dim2 = 0; \ - ps##TYPE *o = NULL; \ - ps##TYPE *i1 = NULL; \ - ps##TYPE *i2 = NULL; \ - n2 = ((psVector* )IN2)->n; \ - dim2 = ((psVector* )IN2)->type.dimen; \ - numRows1 = ((psImage* )IN1)->numRows; \ - numCols1 = ((psImage* )IN1)->numCols; \ - \ - if(dim2 == PS_DIMEN_VECTOR) { /* Regular vectors */ \ - if(n2!=numRows1) { \ - psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ - PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ - n2, numRows1); \ - if (OUT != IN1 && OUT != IN2) { \ - psFree(OUT); \ - } \ - return NULL; \ - } \ - \ - i2 = ((psVector* )IN2)->data.TYPE; \ - for(j = 0; j < numRows1; j++, i1++) { \ - o = ((psImage* )OUT)->data.TYPE[j]; \ - i1 = ((psImage* )IN1)->data.TYPE[j]; \ - for(i = 0; i < numCols1; i++, o++, i1++) { \ - *o = OP; \ - } \ - } \ - } else { /* Transposed vectors */ \ - if(n2!=numCols1) { \ - psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ - PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ - n2, numCols1); \ - if (OUT != IN1) { \ - psFree(OUT); \ - } \ - return NULL; \ - } \ - \ - for(j = 0; j < numRows1; j++) { \ - o = ((psImage* )OUT)->data.TYPE[j]; \ - i1 = ((psVector* )IN2)->data.TYPE; \ - i2 = ((psImage* )IN1)->data.TYPE[j]; \ - for(i = 0; i < numCols1; i++, o++, i2++, i1++) { \ - *o = OP; \ - } \ - } \ - } \ -} - -#define IMAGE_IMAGE(OUT,IN1,OP,IN2,TYPE) \ -{ \ - psS32 i = 0; \ - psS32 j = 0; \ - psS32 numRows1 = 0; \ - psS32 numCols1 = 0; \ - psS32 numRows2 = 0; \ - psS32 numCols2 = 0; \ - ps##TYPE *o = NULL; \ - ps##TYPE *i1 = NULL; \ - ps##TYPE *i2 = NULL; \ - numRows1 = ((psImage* )IN1)->numRows; \ - numCols1 = ((psImage* )IN1)->numCols; \ - numRows2 = ((psImage* )IN2)->numRows; \ - numCols2 = ((psImage* )IN2)->numCols; \ - if(numRows1!=numRows2 || numCols1!=numCols2) { \ - psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ - PS_ERRORTEXT_psMatrix_IMAGE_SIZE_DIFFERS, \ - numCols1, numRows1, numCols2, numRows2); \ - if (OUT != IN1 && OUT != IN2) { \ - psFree(OUT); \ - } \ - return NULL; \ - } \ - for(j = 0; j < numRows1; j++) { \ - o = ((psImage* )OUT)->data.TYPE[j]; \ - i1 = ((psImage* )IN1)->data.TYPE[j]; \ - i2 = ((psImage* )IN2)->data.TYPE[j]; \ - for(i = 0; i < numCols1; i++, o++, i1++, i2++) { \ - *o = OP; \ - } \ - } \ -} - -// Preprocessor macro function to create arithmetic function based on input type -#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); \ - break; \ -case PS_TYPE_F64: \ - DIM1##_##DIM2(OUT,IN1,OP,IN2,F64); \ - break; \ -case PS_TYPE_C32: \ - DIM1##_##DIM2(OUT,IN1,OP,IN2,C32); \ - break; \ -case PS_TYPE_C64: \ - DIM1##_##DIM2(OUT,IN1,OP,IN2,C64); \ - break; \ -default: \ - /* char* strType; \ - PS_TYPE_NAME(strType,IN1->type); \ - psError(PS_ERR_BAD_PARAMETER_TYPE, true, \ - PS_ERRORTEXT_psMatrix_TYPE_MISMATCH, \ - strType); */ \ - if (OUT != IN1 && OUT != IN2) { \ - psFree(OUT); \ - } \ - return NULL; \ -} - -// Preprocessor macro function to create arithmetic function operation name -#define BINARY_OP(DIM1,DIM2,OUT,IN1,OP,IN2) \ -if(!strncmp(OP, "=", 1)) { \ - BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i2,IN2); \ -} else if(!strncmp(OP, "+", 1)) { \ - BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 + *i2,IN2); \ -} else if(!strncmp(OP, "-", 1)) { \ - BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 - *i2,IN2); \ -} else if(!strncmp(OP, "*", 1)) { \ - BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 * *i2,IN2); \ -} else if(!strncmp(OP, "/", 1)) { \ - BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 / *i2,IN2); \ -} else if(!strncmp(OP, "^", 1)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN1->type)) { \ - BINARY_TYPE(DIM1,DIM2,OUT,IN1,cpow(*i1,*i2),IN2); \ - } else { \ - BINARY_TYPE(DIM1,DIM2,OUT,IN1,pow(*i1,*i2),IN2); \ - } \ -} else if(!strncmp(OP, "min", 3)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN1->type)) { \ - psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ - PS_ERRORTEXT_psMatrix_MIN_COMPLEX_SUPPORT); \ - if (OUT != IN1 && OUT != IN2) { \ - psFree(OUT); \ - } \ - return NULL; \ - } else { \ - BINARY_TYPE(DIM1,DIM2,OUT,IN1,fmin(*i1,*i2),IN2); \ - } \ -} else if(!strncmp(OP, "max", 3)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN1->type)) { \ - psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ - PS_ERRORTEXT_psMatrix_MAX_COMPLEX_SUPPORT); \ - if (OUT != IN1 && OUT != IN2) { \ - psFree(OUT); \ - } \ - return NULL; \ - } else { \ - BINARY_TYPE(DIM1,DIM2,OUT,IN1,fmax(*i1,*i2),IN2); \ - } \ -} else { \ - psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ - PS_ERRORTEXT_psMatrix_OPERATION_UNSUPPORTED, \ - OP); \ - if (OUT != IN1 && OUT != IN2) { \ - psFree(OUT); \ - } \ - return NULL; \ -} - -psPtr psBinaryOp(psPtr out, psPtr in1, char *op, psPtr in2) -{ - - psVector* input1 = (psVector* ) in1; - psVector* input2 = (psVector* ) in2; - - #define psBinaryOp_EXIT { \ - if (out != in1 && out != in2) { \ - psFree(out); \ - } \ - return NULL; \ - } - - PS_PTR_CHECK_NULL_GENERAL(input1, psBinaryOp_EXIT); - PS_PTR_CHECK_NULL_GENERAL(input2, psBinaryOp_EXIT); - PS_PTR_CHECK_NULL_GENERAL(op, psBinaryOp_EXIT); - - PS_PTR_CHECK_TYPE_EQUAL_GENERAL(input1,input2, psBinaryOp_EXIT); - - PS_PTR_CHECK_DIMEN_GENERAL_NOT(input1, PS_DIMEN_OTHER, psBinaryOp_EXIT); - PS_PTR_CHECK_DIMEN_GENERAL_NOT(input2, PS_DIMEN_OTHER, psBinaryOp_EXIT); - - psType* psType1 = (psType*)in1; - psType* psType2 = (psType*)in2; - psDimen dim1 = psType1->dimen; - psDimen dim2 = psType2->dimen; - psElemType elType1 = psType1->type; - psElemType elType2 = psType2->type; - - if (dim1 == PS_DIMEN_VECTOR || dim1 == PS_DIMEN_TRANSV) { - if (((psVector* ) in1)->n == 0) { - psLogMsg(__func__, PS_LOG_WARN, "Vector contains zero elements"); - } - } else if (dim1 == PS_DIMEN_IMAGE) { - if (((psImage* ) in1)->numCols == 0 || ((psImage* ) in1)->numRows == 0) { - psLogMsg(__func__, PS_LOG_WARN, "Image contains zero length row or cols"); - } - } - - if (dim2 == PS_DIMEN_VECTOR || dim2 == PS_DIMEN_TRANSV) { - if (((psVector* ) in2)->n == 0) { - psLogMsg(__func__, PS_LOG_WARN, "Vector contains zero elements"); - } - } else if (dim2 == PS_DIMEN_IMAGE) { - if (((psImage* ) in2)->numCols == 0 || ((psImage* ) in2)->numRows == 0) { - psLogMsg(__func__, PS_LOG_WARN, "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(PS_ERR_UNKNOWN, false, - PS_ERRORTEXT_psMatrix_OUTPUT_VECTOR_NOT_CREATED); - 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(PS_ERR_UNKNOWN, false, - PS_ERRORTEXT_psMatrix_OUTPUT_IMAGE_NOT_CREATED); - return NULL; - } - BINARY_OP(SCALAR, IMAGE, out, psType1, op, psType2); // scalar op image - } else { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psMatrix_DIMEN_INVALID, - "in2",dim2); - psBinaryOp_EXIT; - } - } 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(PS_ERR_UNKNOWN, false, - PS_ERRORTEXT_psMatrix_OUTPUT_VECTOR_NOT_CREATED); - 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(PS_ERR_UNKNOWN, false, - PS_ERRORTEXT_psMatrix_OUTPUT_VECTOR_NOT_CREATED); - 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(PS_ERR_UNKNOWN, false, - PS_ERRORTEXT_psMatrix_OUTPUT_IMAGE_NOT_CREATED); - return NULL; - } - BINARY_OP(VECTOR, IMAGE, out, psType1, op, psType2); // vector op image - } else { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psMatrix_DIMEN_INVALID, - "in2",dim2); - psBinaryOp_EXIT; - } - } else if (dim1 == PS_DIMEN_IMAGE) { - out = psImageRecycle(out, ((psImage*)in1)->numCols, ((psImage*)in1)->numRows, elType1); - if (out == NULL) { - psError(PS_ERR_UNKNOWN, false, - PS_ERRORTEXT_psMatrix_OUTPUT_IMAGE_NOT_CREATED); - return NULL; - } - if (dim2 == PS_DIMEN_SCALAR) { - BINARY_OP(IMAGE, SCALAR, out, psType1, op, psType2); // image op scalar - } else if (dim2 == PS_DIMEN_VECTOR || dim2 == PS_DIMEN_TRANSV) { - BINARY_OP(IMAGE, VECTOR, out, psType1, op, psType2); // image op vector - } else if (dim2 == PS_DIMEN_IMAGE) { - BINARY_OP(IMAGE, IMAGE, out, psType1, op, psType2); // image op image - } else { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psMatrix_DIMEN_INVALID, - "in2",dim2); - psBinaryOp_EXIT; - } - } else { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psMatrix_DIMEN_INVALID, - "in1",dim1); - psBinaryOp_EXIT; - } - - // Automtically free psScalar types, since they are usually allocated in the argument list when this - // function is called, provided that the input is not the output. - if(psType1->dimen==PS_DIMEN_SCALAR && in1!=out) { - psFree(in1); - } - - if(psType2->dimen==PS_DIMEN_SCALAR && in2!=out) { - psFree(in2); - } - - return out; -} - -// Unary SCALAR operations -#define SCALAR(OUT,IN,OP,TYPE) \ -{ \ - ps##TYPE *o = NULL; \ - ps##TYPE *i1 = NULL; \ - o = &((psScalar* )OUT)->data.TYPE; \ - i1 = &((psScalar* )IN)->data.TYPE; \ - *o = OP; \ -} - -// Unary IMAGE operations -#define VECTOR(OUT,IN,OP,TYPE) \ -{ \ - psS32 i = 0; \ - psS32 nIn = 0; \ - psS32 nOut = 0; \ - ps##TYPE *o = NULL; \ - ps##TYPE *i1 = NULL; \ - nIn = ((psVector* )IN)->n; \ - nOut = ((psVector* )OUT)->n; \ - if(nIn != nOut) { \ - psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ - PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ - nIn, nOut); \ - if (OUT != IN) { \ - psFree(OUT); \ - } \ - return NULL; \ - } \ - o = ((psVector* )OUT)->data.TYPE; \ - i1 = ((psVector* )IN)->data.TYPE; \ - for(i = 0; i < nIn; i++, o++, i1++) { \ - *o = OP; \ - } \ -} - -// Unary IMAGE operations -#define IMAGE(OUT,IN,OP,TYPE) \ -{ \ - psS32 i = 0; \ - psS32 j = 0; \ - psS32 numRowsIn = 0; \ - psS32 numColsIn = 0; \ - psS32 numRowsOut = 0; \ - psS32 numColsOut = 0; \ - ps##TYPE *o = NULL; \ - ps##TYPE *i1 = NULL; \ - numRowsIn = ((psImage* )IN)->numRows; \ - numColsIn = ((psImage* )IN)->numCols; \ - numRowsOut = ((psImage* )OUT)->numRows; \ - numColsOut = ((psImage* )OUT)->numCols; \ - if(numRowsIn!=numRowsOut || numColsIn!=numColsOut) { \ - psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ - PS_ERRORTEXT_psMatrix_IMAGE_SIZE_DIFFERS, \ - numColsIn, numRowsIn, numColsOut, numRowsOut); \ - if (OUT != IN) { \ - psFree(OUT); \ - } \ - return NULL; \ - } \ - for(j = 0; j < numRowsIn; j++) { \ - o = ((psImage* )OUT)->data.TYPE[j]; \ - i1 = ((psImage* )IN)->data.TYPE[j]; \ - for(i = 0; i < numColsIn; i++, o++, i1++) { \ - *o = OP; \ - } \ - } \ -} - -// Preprocessor macro function to create arithmetic function based on input type -#define UNARY_TYPE(DIM,OUT,IN,OP) \ -switch (IN->type) { \ -case PS_TYPE_S32: \ - DIM(OUT,IN,OP,S32); \ - break; \ -case PS_TYPE_F32: \ - DIM(OUT,IN,OP,F32); \ - break; \ -case PS_TYPE_F64: \ - DIM(OUT,IN,OP,F64); \ - break; \ -case PS_TYPE_C32: \ - 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: { \ - char* strType; \ - PS_TYPE_NAME(strType, IN->type); \ - psError(PS_ERR_BAD_PARAMETER_TYPE, true, \ - PS_ERRORTEXT_psMatrix_TYPE_MISMATCH, \ - strType); \ - if (OUT != IN) { \ - psFree(OUT); \ - } \ - return NULL; \ - } \ -} - -// Preprocessor macro function to create arithmetic function operation name. Functions below that add -// FLT_EPSILON are done so to align results with a 64 bit computing architecture -#define UNARY_OP(DIM,OUT,IN,OP) \ -if(!strncmp(OP, "abs", 3)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,cabs(*i1)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,fabs(*i1)); \ - } \ -} else if(!strncmp(OP, "exp", 3)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,cexp(*i1)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,exp(*i1)); \ - } \ -} else if(!strncmp(OP, "ln", 2)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,clog(*i1)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,log(*i1)); \ - } \ -} else if(!strncmp(OP, "ten", 3)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,cpow(10.0,*i1)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,pow(10.0,*i1)); \ - } \ -} else if(!strncmp(OP, "log", 3)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,clog(*i1)/log(10.0)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,log10(*i1)); \ - } \ -} else if(!strncmp(OP, "sin", 3)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,csin(*i1)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,sin(*i1)); \ - } \ -} else if(!strncmp(OP, "dsin", 4)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,csin(*i1*D2R)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,sin(*i1*D2R)); \ - } \ -} else if(!strncmp(OP, "cos", 3)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,ccos(*i1)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,cos(*i1)); \ - } \ -} else if(!strncmp(OP, "dcos", 4)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,ccos(*i1*D2R)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,cos(*i1*D2R)); \ - } \ -} else if(!strncmp(OP, "tan", 3)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,ctan(*i1)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,tan(*i1)); \ - } \ -} else if(!strncmp(OP, "dtan", 4)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,ctan(*i1*D2R)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,tan(*i1*D2R)); \ - } \ -} else if(!strncmp(OP, "asin", 4)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,casin(*i1)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,asin(*i1)); \ - } \ -} else if(!strncmp(OP, "dasin", 5)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,R2D*casin(*i1)); \ - } else if(PS_IS_PSELEMTYPE_INT(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,(R2D*asin(*i1))); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,(R2D*asin(*i1))); \ - } \ -} else if(!strncmp(OP, "acos", 4)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,cacos(*i1)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,acos(*i1)); \ - } \ -} else if(!strncmp(OP, "dacos", 5)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,R2D*cacos(*i1)); \ - } else if(PS_IS_PSELEMTYPE_INT(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,(R2D*acos(*i1))); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,R2D*acos(*i1)); \ - } \ -} else if(!strncmp(OP, "atan", 4)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,catan(*i1)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,atan(*i1)); \ - } \ -} else if(!strncmp(OP, "datan", 5)) { \ - if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ - UNARY_TYPE(DIM,OUT,IN,R2D*catan(*i1)); \ - } else { \ - UNARY_TYPE(DIM,OUT,IN,R2D*atan(*i1)); \ - } \ -} else { \ - psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ - PS_ERRORTEXT_psMatrix_OPERATION_UNSUPPORTED, \ - OP); \ - if (OUT != IN) { \ - psFree(OUT); \ - } \ - return NULL; \ -} - -psPtr psUnaryOp(psPtr out, psPtr in, char *op) -{ - #define psUnaryOp_EXIT { \ - if (out != in) { \ - psFree(out); \ - } \ - return NULL; \ - } - - psType* psTypeIn = (psType* ) in; - - PS_PTR_CHECK_NULL_GENERAL(in, psUnaryOp_EXIT); - PS_PTR_CHECK_NULL_GENERAL(op, psUnaryOp_EXIT); - - psDimen dimIn = psTypeIn->dimen; - psElemType elTypeIn = psTypeIn->type; - - 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) { - psError(PS_ERR_BAD_PARAMETER_SIZE, true, - PS_ERRORTEXT_psMatrix_VECTOR_EMPTY); - psUnaryOp_EXIT; - } - - out = psVectorRecycle(out, - ((psVector*)in)->n, - elTypeIn); - if (out == NULL) { - psError(PS_ERR_UNKNOWN, false, - PS_ERRORTEXT_psMatrix_OUTPUT_VECTOR_NOT_CREATED); - psUnaryOp_EXIT; - } - - UNARY_OP(VECTOR, out, psTypeIn, op); // vector - break; - case PS_DIMEN_IMAGE: - if (((psImage* ) in)->numCols == 0 || ((psImage* ) in)->numRows == 0) { - psError(PS_ERR_BAD_PARAMETER_SIZE, true, - PS_ERRORTEXT_psMatrix_IMAGE_EMPTY); - psUnaryOp_EXIT; - } - - out = psImageRecycle(out, - ((psImage*)in)->numCols, - ((psImage*)in)->numRows, - elTypeIn); - if (out == NULL) { - psError(PS_ERR_UNKNOWN, false, - PS_ERRORTEXT_psMatrix_OUTPUT_IMAGE_NOT_CREATED); - psUnaryOp_EXIT; - } - - UNARY_OP(IMAGE, out, psTypeIn, op); // image - break; - default: - if (out != in) { - psFree(out); - } - psError(PS_ERR_BAD_PARAMETER_SIZE, true, - PS_ERRORTEXT_psMatrix_DIMEN_INVALID, - "in", dimIn); - psUnaryOp_EXIT; - } - - // Automtically free psScalar types, since they are usually allocated in the argument list when this - // function is called, provided that the input is not the output. - if(psTypeIn->dimen==PS_DIMEN_SCALAR && in!=out) { - psFree(in); - } - - return out; -} Index: trunk/psLib/src/dataManip/psMatrixVectorArithmetic.h =================================================================== --- trunk/psLib/src/dataManip/psMatrixVectorArithmetic.h (revision 3171) +++ (revision ) @@ -1,88 +1,0 @@ - -/** @file psMatrixVectorArithmetic.h - * - * @brief Provides unary and binary functions for simple matrix and vector element operations. Functions - * include: - * - * Addition (+) - * Subtraction (-) - * Multiplication (*) - * Division (/) - * Power (^) - * Minimum (min) - * Maximum (max) - * Absolute value (abs) - * Exponent (exp) - * Natural Log (ln) - * Power of 10 (ten) - * Log (log) - * Sine (sin or dsin) - * Cosine (cos or dcos) - * Tangent (tan or dtan) - * Arcsine (asin or dasin) - * Arccosine (acos or dacos) - * Arctan (atan or datan) - * - * Currently only vector-vector and image-image binary operations are supported. - * - * @ingroup MatrixVectorArithmetic - * - * @author Ross Harman, MHPCC - * - * @version $Revision: 1.7 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-01-12 22:17:01 $ - * - * Copyright 2004 Maui High Performance Computing Center, University of Hawaii - */ - -#ifndef PSMATRIX_VECTOR_ARITHMETIC_H -#define PSMATRIX_VECTOR_ARITHMETIC_H - -/// @addtogroup MatrixVectorArithmetic -/// @{ - -/** Perform simple binary arithmetic with images or vectors - * - * Performs addition, subtraction, multiplication, division, power, minumum, and maximum arithmetic - * operations with images and vectors. Uses the form: - * - * out = in1 op in2, - * - * Where op is: "=", "+", "-", "*", "/", "^", "min", or "max" - * - * This function only supports vector-vector or image-image operations. - * - * @return psType* : Pointer to either psImage or psVector. - */ -psType* psBinaryOp( - psPtr out, ///< Output type, either psImage or psVector. - psPtr in1, ///< First input, either psImage or psVector. - char *op, ///< Operator. - psPtr in2 ///< Second input, either psImage or psVector. -); - -/** Perform simple unary arithmetic with images or vectors - * - * Performs absolute value, exponent, natural log, power of 10, log, sine, cosine, tangent, arcsine, - * arccosine, or arctan. operations with images and vectors. Uses the form: - * - * out = op(in), - * - * Where op is: "abs", "exp", "ln", "ten", "log", "sin", "cos", "tan" "asin", "acos", "atan", "dsin", - * "dcos", dtan", "dasin", "dacos", or "datan". - * - * Trigometric Operations with "d" prefix use units of degrees. Those without are in radians. - * - * This function only supports vector-vector or image-image opertions. - * - * @return psType* : Pointer to either psImage or psVector. - */ -psType* psUnaryOp( - psPtr out, ///< Output type, either psImage or psVector. - psPtr in, ///< Input, either psImage or psVector. - char *op ///< Operator. -); - -/// @} - -#endif Index: trunk/psLib/src/dataManip/psMinimize.c =================================================================== --- trunk/psLib/src/dataManip/psMinimize.c (revision 3171) +++ trunk/psLib/src/dataManip/psMinimize.c (revision 3182) @@ -9,6 +9,6 @@ * @author GLG, MHPCC * - * @version $Revision: 1.104 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-02-03 00:54:10 $ + * @version $Revision: 1.105 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -538,5 +538,5 @@ psF64 flux = normalization * exp(-( u*u/(2.0 * sigmaX * sigmaX) + v*v/(2.0 * sigmaY * sigmaY)))/ - (2.0 * M_PI * sigmaX * sigmaY); + (2.0 * PS_PI * sigmaX * sigmaY); out->data.F32[i] = flux; @@ -865,5 +865,5 @@ for (psS32 i=0;idata.F64[k] * - cos(M_PI * ((psF32) j) * (0.5 + ((psF32) k)) / ((psF32) n)); + cos(PS_PI * ((psF32) j) * (0.5 + ((psF32) k)) / ((psF32) n)); } Index: trunk/psLib/src/dataManip/psUnaryOp.c =================================================================== --- trunk/psLib/src/dataManip/psUnaryOp.c (revision 3182) +++ trunk/psLib/src/dataManip/psUnaryOp.c (revision 3182) @@ -0,0 +1,387 @@ +/** @file psUnary.c + * + * @brief Provides unary functions for simple matrix and vector element operations. Functions + * include: + * + * Addition (+) + * Subtraction (-) + * Multiplication (*) + * Division (/) + * Power (^) + * Minimum (min) + * Maximum (max) + * Absolute value (abs) + * Exponent (exp) + * Natural Log (ln) + * Power of 10 (ten) + * Log (log) + * Sine (sin or dsin) + * Cosine (cos or dcos) + * Tangent (tan or dtan) + * Arcsine (asin or dasin) + * Arccosine (acos or dacos) + * Arctan (atan or datan) + * + * Currently only vector-vector and image-image binary operations are supported. + * + * @ingroup MatrixArithmetic + * + * @author Ross Harman, MHPCC + * @author Robert DeSonia, MHPCC + * + * @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ + * + * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii + */ + +/****************************************************************************** + * INCLUDE FILES * + ******************************************************************************/ +#include +#include +#include + +#include "psMemory.h" +#include "psError.h" +#include "psImage.h" +#include "psVector.h" +#include "psScalar.h" +#include "psLogMsg.h" +#include "psConstants.h" +#include "psDataManipErrors.h" + +/***************************************************************************** + * FUNCTION IMPLEMENTATION - LOCAL * + *****************************************************************************/ + +// Conversion for degrees to radians +#define D2R 0.01745329252111111 /* PI/180 */ + +// Conversion for radians to degrees +#define R2D 57.29577950924861 /* 180.0/PI */ + + +// Unary SCALAR operations +#define SCALAR(OUT,IN,OP,TYPE) \ +{ \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + o = &((psScalar* )OUT)->data.TYPE; \ + i1 = &((psScalar* )IN)->data.TYPE; \ + *o = OP; \ +} + +// Unary IMAGE operations +#define VECTOR(OUT,IN,OP,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 nIn = 0; \ + psS32 nOut = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + nIn = ((psVector* )IN)->n; \ + nOut = ((psVector* )OUT)->n; \ + if(nIn != nOut) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ + nIn, nOut); \ + if (OUT != IN) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + o = ((psVector* )OUT)->data.TYPE; \ + i1 = ((psVector* )IN)->data.TYPE; \ + for(i = 0; i < nIn; i++, o++, i1++) { \ + *o = OP; \ + } \ +} + +// Unary IMAGE operations +#define IMAGE(OUT,IN,OP,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 j = 0; \ + psS32 numRowsIn = 0; \ + psS32 numColsIn = 0; \ + psS32 numRowsOut = 0; \ + psS32 numColsOut = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + numRowsIn = ((psImage* )IN)->numRows; \ + numColsIn = ((psImage* )IN)->numCols; \ + numRowsOut = ((psImage* )OUT)->numRows; \ + numColsOut = ((psImage* )OUT)->numCols; \ + if(numRowsIn!=numRowsOut || numColsIn!=numColsOut) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_IMAGE_SIZE_DIFFERS, \ + numColsIn, numRowsIn, numColsOut, numRowsOut); \ + if (OUT != IN) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + for(j = 0; j < numRowsIn; j++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = ((psImage* )IN)->data.TYPE[j]; \ + for(i = 0; i < numColsIn; i++, o++, i1++) { \ + *o = OP; \ + } \ + } \ +} + +// Preprocessor macro function to create arithmetic function based on input type +#define UNARY_TYPE(DIM,OUT,IN,OP) \ +switch (IN->type) { \ +case PS_TYPE_S32: \ + DIM(OUT,IN,OP,S32); \ + break; \ +case PS_TYPE_F32: \ + DIM(OUT,IN,OP,F32); \ + break; \ +case PS_TYPE_F64: \ + DIM(OUT,IN,OP,F64); \ + break; \ +case PS_TYPE_C32: \ + 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: { \ + char* strType; \ + PS_TYPE_NAME(strType, IN->type); \ + psError(PS_ERR_BAD_PARAMETER_TYPE, true, \ + PS_ERRORTEXT_psMatrix_TYPE_MISMATCH, \ + strType); \ + if (OUT != IN) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ +} + +// Preprocessor macro function to create arithmetic function operation name. Functions below that add +// FLT_EPSILON are done so to align results with a 64 bit computing architecture +#define UNARY_OP(DIM,OUT,IN,OP) \ +if(!strncmp(OP, "abs", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,cabs(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,fabs(*i1)); \ + } \ +} else if(!strncmp(OP, "exp", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,cexp(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,exp(*i1)); \ + } \ +} else if(!strncmp(OP, "ln", 2)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,clog(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,log(*i1)); \ + } \ +} else if(!strncmp(OP, "ten", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,cpow(10.0,*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,pow(10.0,*i1)); \ + } \ +} else if(!strncmp(OP, "log", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,clog(*i1)/log(10.0)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,log10(*i1)); \ + } \ +} else if(!strncmp(OP, "sin", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,csin(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,sin(*i1)); \ + } \ +} else if(!strncmp(OP, "dsin", 4)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,csin(*i1*D2R)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,sin(*i1*D2R)); \ + } \ +} else if(!strncmp(OP, "cos", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,ccos(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,cos(*i1)); \ + } \ +} else if(!strncmp(OP, "dcos", 4)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,ccos(*i1*D2R)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,cos(*i1*D2R)); \ + } \ +} else if(!strncmp(OP, "tan", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,ctan(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,tan(*i1)); \ + } \ +} else if(!strncmp(OP, "dtan", 4)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,ctan(*i1*D2R)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,tan(*i1*D2R)); \ + } \ +} else if(!strncmp(OP, "asin", 4)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,casin(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,asin(*i1)); \ + } \ +} else if(!strncmp(OP, "dasin", 5)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,R2D*casin(*i1)); \ + } else if(PS_IS_PSELEMTYPE_INT(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,(R2D*asin(*i1))); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,(R2D*asin(*i1))); \ + } \ +} else if(!strncmp(OP, "acos", 4)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,cacos(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,acos(*i1)); \ + } \ +} else if(!strncmp(OP, "dacos", 5)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,R2D*cacos(*i1)); \ + } else if(PS_IS_PSELEMTYPE_INT(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,(R2D*acos(*i1))); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,R2D*acos(*i1)); \ + } \ +} else if(!strncmp(OP, "atan", 4)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,catan(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,atan(*i1)); \ + } \ +} else if(!strncmp(OP, "datan", 5)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,R2D*catan(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,R2D*atan(*i1)); \ + } \ +} else { \ + psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ + PS_ERRORTEXT_psMatrix_OPERATION_UNSUPPORTED, \ + OP); \ + if (OUT != IN) { \ + psFree(OUT); \ + } \ + return NULL; \ +} + +psPtr psUnaryOp(psPtr out, psPtr in, char *op) +{ + #define psUnaryOp_EXIT { \ + if (out != in) { \ + psFree(out); \ + } \ + return NULL; \ + } + + psType* psTypeIn = (psType* ) in; + + PS_PTR_CHECK_NULL_GENERAL(in, psUnaryOp_EXIT); + PS_PTR_CHECK_NULL_GENERAL(op, psUnaryOp_EXIT); + + psDimen dimIn = psTypeIn->dimen; + psElemType elTypeIn = psTypeIn->type; + + 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) { + psError(PS_ERR_BAD_PARAMETER_SIZE, true, + PS_ERRORTEXT_psMatrix_VECTOR_EMPTY); + psUnaryOp_EXIT; + } + + out = psVectorRecycle(out, + ((psVector*)in)->n, + elTypeIn); + if (out == NULL) { + psError(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_VECTOR_NOT_CREATED); + psUnaryOp_EXIT; + } + + UNARY_OP(VECTOR, out, psTypeIn, op); // vector + break; + case PS_DIMEN_IMAGE: + if (((psImage* ) in)->numCols == 0 || ((psImage* ) in)->numRows == 0) { + psError(PS_ERR_BAD_PARAMETER_SIZE, true, + PS_ERRORTEXT_psMatrix_IMAGE_EMPTY); + psUnaryOp_EXIT; + } + + out = psImageRecycle(out, + ((psImage*)in)->numCols, + ((psImage*)in)->numRows, + elTypeIn); + if (out == NULL) { + psError(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_IMAGE_NOT_CREATED); + psUnaryOp_EXIT; + } + + UNARY_OP(IMAGE, out, psTypeIn, op); // image + break; + default: + if (out != in) { + psFree(out); + } + psError(PS_ERR_BAD_PARAMETER_SIZE, true, + PS_ERRORTEXT_psMatrix_DIMEN_INVALID, + "in", dimIn); + psUnaryOp_EXIT; + } + + // Automtically free psScalar types, since they are usually allocated in the argument list when this + // function is called, provided that the input is not the output. + if(psTypeIn->dimen==PS_DIMEN_SCALAR && in!=out) { + psFree(in); + } + + return out; +} Index: trunk/psLib/src/dataManip/psUnaryOp.h =================================================================== --- trunk/psLib/src/dataManip/psUnaryOp.h (revision 3182) +++ trunk/psLib/src/dataManip/psUnaryOp.h (revision 3182) @@ -0,0 +1,69 @@ +/** @file psUnaryOp.h + * + * @brief Provides unary functions for simple matrix and vector element operations. Functions + * include: + * + * Addition (+) + * Subtraction (-) + * Multiplication (*) + * Division (/) + * Power (^) + * Minimum (min) + * Maximum (max) + * Absolute value (abs) + * Exponent (exp) + * Natural Log (ln) + * Power of 10 (ten) + * Log (log) + * Sine (sin or dsin) + * Cosine (cos or dcos) + * Tangent (tan or dtan) + * Arcsine (asin or dasin) + * Arccosine (acos or dacos) + * Arctan (atan or datan) + * + * Currently only vector-vector and image-image binary operations are supported. + * + * @ingroup MatrixArithmetic + * + * @author Ross Harman, MHPCC + * @author Robert DeSonia, MHPCC + * + * @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ + * + * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii + */ + +#ifndef PSUNARY_OP_H +#define PSUNARY_OP_H + +/// @addtogroup MatrixArithmetic +/// @{ + + +/** Perform simple unary arithmetic with images or vectors + * + * Performs absolute value, exponent, natural log, power of 10, log, sine, cosine, tangent, arcsine, + * arccosine, or arctan. operations with images and vectors. Uses the form: + * + * out = op(in), + * + * Where op is: "abs", "exp", "ln", "ten", "log", "sin", "cos", "tan" "asin", "acos", "atan", "dsin", + * "dcos", dtan", "dasin", "dacos", or "datan". + * + * Trigometric Operations with "d" prefix use units of degrees. Those without are in radians. + * + * This function only supports vector-vector or image-image opertions. + * + * @return psType* : Pointer to either psImage or psVector. + */ +psType* psUnaryOp( + psPtr out, ///< Output type, either psImage or psVector. + psPtr in, ///< Input, either psImage or psVector. + char *op ///< Operator. +); + +/// @} + +#endif Index: trunk/psLib/src/image/psImageStats.c =================================================================== --- trunk/psLib/src/image/psImageStats.c (revision 3171) +++ trunk/psLib/src/image/psImageStats.c (revision 3182) @@ -9,6 +9,6 @@ * @author GLG, MHPCC * - * @version $Revision: 1.64 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-02-03 00:54:11 $ + * @version $Revision: 1.65 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:42 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -139,5 +139,5 @@ for (i = 0; i < n; i++) { tmp = (double)(n - i); - tmp = (M_PI * (tmp - 0.5)) / ((double)n); + tmp = (PS_PI * (tmp - 0.5)) / ((double)n); scalingFactors[i] = cos(tmp); } @@ -348,10 +348,10 @@ // Chebyshev polynomials are 0. for (x = 0; x < input->numRows; x++) { - double xTmp = cos(M_PI * (0.5 + ((float) x)) / ((float) input->numRows)); + double xTmp = cos(PS_PI * (0.5 + ((float) x)) / ((float) input->numRows)); double xNode = - ((xTmp + bma + bpa) - 1.0); double xOrig = ((float) input->numRows) * (xNode - min) / (max - min); for (y = 0; y < input->numCols; y++) { - double yTmp = cos(M_PI * (0.5 + ((float) y)) / ((float) input->numCols)); + double yTmp = cos(PS_PI * (0.5 + ((float) y)) / ((float) input->numCols)); double yNode = - ((yTmp + bma + bpa) - 1.0); double yOrig = ((float) input->numCols) * (yNode - min) / (max - min); Index: trunk/psLib/src/imageops/psImageStats.c =================================================================== --- trunk/psLib/src/imageops/psImageStats.c (revision 3171) +++ trunk/psLib/src/imageops/psImageStats.c (revision 3182) @@ -9,6 +9,6 @@ * @author GLG, MHPCC * - * @version $Revision: 1.64 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-02-03 00:54:11 $ + * @version $Revision: 1.65 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:42 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -139,5 +139,5 @@ for (i = 0; i < n; i++) { tmp = (double)(n - i); - tmp = (M_PI * (tmp - 0.5)) / ((double)n); + tmp = (PS_PI * (tmp - 0.5)) / ((double)n); scalingFactors[i] = cos(tmp); } @@ -348,10 +348,10 @@ // Chebyshev polynomials are 0. for (x = 0; x < input->numRows; x++) { - double xTmp = cos(M_PI * (0.5 + ((float) x)) / ((float) input->numRows)); + double xTmp = cos(PS_PI * (0.5 + ((float) x)) / ((float) input->numRows)); double xNode = - ((xTmp + bma + bpa) - 1.0); double xOrig = ((float) input->numRows) * (xNode - min) / (max - min); for (y = 0; y < input->numCols; y++) { - double yTmp = cos(M_PI * (0.5 + ((float) y)) / ((float) input->numCols)); + double yTmp = cos(PS_PI * (0.5 + ((float) y)) / ((float) input->numCols)); double yNode = - ((yTmp + bma + bpa) - 1.0); double yOrig = ((float) input->numCols) * (yNode - min) / (max - min); Index: trunk/psLib/src/math/psBinaryOp.c =================================================================== --- trunk/psLib/src/math/psBinaryOp.c (revision 3182) +++ trunk/psLib/src/math/psBinaryOp.c (revision 3182) @@ -0,0 +1,582 @@ +/** @file psBinaryOp.c + * + * @brief Provides binary functions for simple matrix and vector element operations. Functions + * include: + * + * Addition (+) + * Subtraction (-) + * Multiplication (*) + * Division (/) + * Power (^) + * Minimum (min) + * Maximum (max) + * Absolute value (abs) + * Exponent (exp) + * Natural Log (ln) + * Power of 10 (ten) + * Log (log) + * Sine (sin or dsin) + * Cosine (cos or dcos) + * Tangent (tan or dtan) + * Arcsine (asin or dasin) + * Arccosine (acos or dacos) + * Arctan (atan or datan) + * + * Currently only vector-vector and image-image binary operations are supported. + * + * @ingroup MatrixArithmetic + * + * @author Ross Harman, MHPCC + * @author Robert DeSonia, MHPCC + * + * @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ + * + * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii + */ + +/****************************************************************************** + * INCLUDE FILES * + ******************************************************************************/ +#include +#include +#include + +#include "psMemory.h" +#include "psError.h" +#include "psImage.h" +#include "psVector.h" +#include "psScalar.h" +#include "psLogMsg.h" +#include "psConstants.h" +#include "psDataManipErrors.h" + +/***************************************************************************** + * FUNCTION IMPLEMENTATION - LOCAL * + *****************************************************************************/ + +// Conversion for degrees to radians +#define D2R 0.01745329252111111 /* PI/180 */ + +// Conversion for radians to degrees +#define R2D 57.29577950924861 /* 180.0/PI */ + +// Binary SCALAR_XXXX operations +#define SCALAR_SCALAR(OUT,IN1,OP,IN2,TYPE) \ +{ \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + o = &((psScalar* )OUT)->data.TYPE; \ + i1 = &((psScalar* )IN1)->data.TYPE; \ + i2 = &((psScalar* )IN2)->data.TYPE; \ + *o = OP; \ +} + +#define SCALAR_VECTOR(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 npt = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + npt = ((psVector* )IN2)->n; \ + o = ((psVector* )OUT)->data.TYPE; \ + i1 = &((psScalar* )IN1)->data.TYPE; \ + i2 = ((psVector* )IN2)->data.TYPE; \ + for (i=0; i < npt; i++, o++, i2++) { \ + *o = OP; \ + } \ +} + +#define SCALAR_IMAGE(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 j = 0; \ + psS32 numRows = 0; \ + psS32 numCols = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + numRows = ((psImage* )IN2)->numRows; \ + numCols = ((psImage* )IN2)->numCols; \ + for(j = 0; j < numCols; j++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = &((psScalar* )IN1)->data.TYPE; \ + i2 = ((psImage* )IN2)->data.TYPE[j]; \ + for(i = 0; i < numRows; i++, o++, i2++) { \ + *o = OP; \ + } \ + } \ +} + +// Binary VECTOR_XXXX operations +#define VECTOR_SCALAR(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 n1 = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + n1 = ((psVector* )IN1)->n; \ + o = ((psVector* )OUT)->data.TYPE; \ + i1 = ((psVector* )IN1)->data.TYPE; \ + i2 = &((psScalar* )IN2)->data.TYPE; \ + for (i=0; i < n1; i++, o++, i1++) { \ + *o = OP; \ + } \ +} + +#define VECTOR_VECTOR(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 n1 = 0; \ + psS32 n2 = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + n1 = ((psVector* )IN1)->n; \ + n2 = ((psVector* )IN2)->n; \ + if(n1 != n2) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ + n1, n2); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + o = ((psVector* )OUT)->data.TYPE; \ + i1 = ((psVector* )IN1)->data.TYPE; \ + i2 = ((psVector* )IN2)->data.TYPE; \ + for (i=0; i < n1; i++, o++, i1++, i2++) { \ + *o = OP; \ + } \ +} + +#define VECTOR_IMAGE(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 j = 0; \ + psS32 n1 = 0; \ + psS32 numRows2 = 0; \ + psS32 numCols2 = 0; \ + psDimen dim1 = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + n1 = ((psVector* )IN1)->n; \ + dim1 = ((psVector* )IN1)->type.dimen; \ + numRows2 = ((psImage* )IN2)->numRows; \ + numCols2 = ((psImage* )IN2)->numCols; \ + \ + if(dim1 == PS_DIMEN_VECTOR) { /* Regular vectors */ \ + if(n1!=numRows2) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ + n1, numRows2); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + \ + i1 = ((psVector* )IN1)->data.TYPE; \ + for(j = 0; j < numRows2; j++, i1++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i2 = ((psImage* )IN2)->data.TYPE[j]; \ + for(i = 0; i < numCols2; i++, o++, i2++) { \ + *o = OP; \ + } \ + } \ + } else { /* Transposed vectors */ \ + if(n1!=numCols2) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ + n1, numCols2); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + \ + for(j = 0; j < numRows2; j++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = ((psVector* )IN1)->data.TYPE; \ + i2 = ((psImage* )IN2)->data.TYPE[j]; \ + for(i = 0; i < numCols2; i++, o++, i1++, i2++) { \ + *o = OP; \ + } \ + } \ + } \ +} + +// Binary IMAGE_XXXX operations +#define IMAGE_SCALAR(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 j = 0; \ + psS32 numRows = 0; \ + psS32 numCols = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + numRows = ((psImage* )IN1)->numRows; \ + numCols = ((psImage* )IN1)->numCols; \ + for(j = 0; j < numRows; j++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = ((psImage* )IN1)->data.TYPE[j]; \ + i2 = &((psScalar* )IN2)->data.TYPE; \ + for(i = 0; i < numCols; i++, o++, i1++) { \ + *o = OP; \ + } \ + } \ +} + +#define IMAGE_VECTOR(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 j = 0; \ + psS32 n2 = 0; \ + psS32 numRows1 = 0; \ + psS32 numCols1 = 0; \ + psDimen dim2 = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + n2 = ((psVector* )IN2)->n; \ + dim2 = ((psVector* )IN2)->type.dimen; \ + numRows1 = ((psImage* )IN1)->numRows; \ + numCols1 = ((psImage* )IN1)->numCols; \ + \ + if(dim2 == PS_DIMEN_VECTOR) { /* Regular vectors */ \ + if(n2!=numRows1) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ + n2, numRows1); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + \ + i2 = ((psVector* )IN2)->data.TYPE; \ + for(j = 0; j < numRows1; j++, i1++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = ((psImage* )IN1)->data.TYPE[j]; \ + for(i = 0; i < numCols1; i++, o++, i1++) { \ + *o = OP; \ + } \ + } \ + } else { /* Transposed vectors */ \ + if(n2!=numCols1) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ + n2, numCols1); \ + if (OUT != IN1) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + \ + for(j = 0; j < numRows1; j++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = ((psVector* )IN2)->data.TYPE; \ + i2 = ((psImage* )IN1)->data.TYPE[j]; \ + for(i = 0; i < numCols1; i++, o++, i2++, i1++) { \ + *o = OP; \ + } \ + } \ + } \ +} + +#define IMAGE_IMAGE(OUT,IN1,OP,IN2,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 j = 0; \ + psS32 numRows1 = 0; \ + psS32 numCols1 = 0; \ + psS32 numRows2 = 0; \ + psS32 numCols2 = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + ps##TYPE *i2 = NULL; \ + numRows1 = ((psImage* )IN1)->numRows; \ + numCols1 = ((psImage* )IN1)->numCols; \ + numRows2 = ((psImage* )IN2)->numRows; \ + numCols2 = ((psImage* )IN2)->numCols; \ + if(numRows1!=numRows2 || numCols1!=numCols2) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_IMAGE_SIZE_DIFFERS, \ + numCols1, numRows1, numCols2, numRows2); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + for(j = 0; j < numRows1; j++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = ((psImage* )IN1)->data.TYPE[j]; \ + i2 = ((psImage* )IN2)->data.TYPE[j]; \ + for(i = 0; i < numCols1; i++, o++, i1++, i2++) { \ + *o = OP; \ + } \ + } \ +} + +// Preprocessor macro function to create arithmetic function based on input type +#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); \ + break; \ +case PS_TYPE_F64: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,F64); \ + break; \ +case PS_TYPE_C32: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,C32); \ + break; \ +case PS_TYPE_C64: \ + DIM1##_##DIM2(OUT,IN1,OP,IN2,C64); \ + break; \ +default: \ + /* char* strType; \ + PS_TYPE_NAME(strType,IN1->type); \ + psError(PS_ERR_BAD_PARAMETER_TYPE, true, \ + PS_ERRORTEXT_psMatrix_TYPE_MISMATCH, \ + strType); */ \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ +} + +// Preprocessor macro function to create arithmetic function operation name +#define BINARY_OP(DIM1,DIM2,OUT,IN1,OP,IN2) \ +if(!strncmp(OP, "=", 1)) { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i2,IN2); \ +} else if(!strncmp(OP, "+", 1)) { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 + *i2,IN2); \ +} else if(!strncmp(OP, "-", 1)) { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 - *i2,IN2); \ +} else if(!strncmp(OP, "*", 1)) { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 * *i2,IN2); \ +} else if(!strncmp(OP, "/", 1)) { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,*i1 / *i2,IN2); \ +} else if(!strncmp(OP, "^", 1)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN1->type)) { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,cpow(*i1,*i2),IN2); \ + } else { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,pow(*i1,*i2),IN2); \ + } \ +} else if(!strncmp(OP, "min", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN1->type)) { \ + psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ + PS_ERRORTEXT_psMatrix_MIN_COMPLEX_SUPPORT); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } else { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,fmin(*i1,*i2),IN2); \ + } \ +} else if(!strncmp(OP, "max", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN1->type)) { \ + psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ + PS_ERRORTEXT_psMatrix_MAX_COMPLEX_SUPPORT); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ + } else { \ + BINARY_TYPE(DIM1,DIM2,OUT,IN1,fmax(*i1,*i2),IN2); \ + } \ +} else { \ + psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ + PS_ERRORTEXT_psMatrix_OPERATION_UNSUPPORTED, \ + OP); \ + if (OUT != IN1 && OUT != IN2) { \ + psFree(OUT); \ + } \ + return NULL; \ +} + +psPtr psBinaryOp(psPtr out, psPtr in1, char *op, psPtr in2) +{ + + psVector* input1 = (psVector* ) in1; + psVector* input2 = (psVector* ) in2; + + #define psBinaryOp_EXIT { \ + if (out != in1 && out != in2) { \ + psFree(out); \ + } \ + return NULL; \ + } + + PS_PTR_CHECK_NULL_GENERAL(input1, psBinaryOp_EXIT); + PS_PTR_CHECK_NULL_GENERAL(input2, psBinaryOp_EXIT); + PS_PTR_CHECK_NULL_GENERAL(op, psBinaryOp_EXIT); + + PS_PTR_CHECK_TYPE_EQUAL_GENERAL(input1,input2, psBinaryOp_EXIT); + + PS_PTR_CHECK_DIMEN_GENERAL_NOT(input1, PS_DIMEN_OTHER, psBinaryOp_EXIT); + PS_PTR_CHECK_DIMEN_GENERAL_NOT(input2, PS_DIMEN_OTHER, psBinaryOp_EXIT); + + psType* psType1 = (psType*)in1; + psType* psType2 = (psType*)in2; + psDimen dim1 = psType1->dimen; + psDimen dim2 = psType2->dimen; + psElemType elType1 = psType1->type; + psElemType elType2 = psType2->type; + + if (dim1 == PS_DIMEN_VECTOR || dim1 == PS_DIMEN_TRANSV) { + if (((psVector* ) in1)->n == 0) { + psLogMsg(__func__, PS_LOG_WARN, "Vector contains zero elements"); + } + } else if (dim1 == PS_DIMEN_IMAGE) { + if (((psImage* ) in1)->numCols == 0 || ((psImage* ) in1)->numRows == 0) { + psLogMsg(__func__, PS_LOG_WARN, "Image contains zero length row or cols"); + } + } + + if (dim2 == PS_DIMEN_VECTOR || dim2 == PS_DIMEN_TRANSV) { + if (((psVector* ) in2)->n == 0) { + psLogMsg(__func__, PS_LOG_WARN, "Vector contains zero elements"); + } + } else if (dim2 == PS_DIMEN_IMAGE) { + if (((psImage* ) in2)->numCols == 0 || ((psImage* ) in2)->numRows == 0) { + psLogMsg(__func__, PS_LOG_WARN, "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(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_VECTOR_NOT_CREATED); + 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(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_IMAGE_NOT_CREATED); + return NULL; + } + BINARY_OP(SCALAR, IMAGE, out, psType1, op, psType2); // scalar op image + } else { + psError(PS_ERR_BAD_PARAMETER_TYPE, true, + PS_ERRORTEXT_psMatrix_DIMEN_INVALID, + "in2",dim2); + psBinaryOp_EXIT; + } + } 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(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_VECTOR_NOT_CREATED); + 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(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_VECTOR_NOT_CREATED); + 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(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_IMAGE_NOT_CREATED); + return NULL; + } + BINARY_OP(VECTOR, IMAGE, out, psType1, op, psType2); // vector op image + } else { + psError(PS_ERR_BAD_PARAMETER_TYPE, true, + PS_ERRORTEXT_psMatrix_DIMEN_INVALID, + "in2",dim2); + psBinaryOp_EXIT; + } + } else if (dim1 == PS_DIMEN_IMAGE) { + out = psImageRecycle(out, ((psImage*)in1)->numCols, ((psImage*)in1)->numRows, elType1); + if (out == NULL) { + psError(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_IMAGE_NOT_CREATED); + return NULL; + } + if (dim2 == PS_DIMEN_SCALAR) { + BINARY_OP(IMAGE, SCALAR, out, psType1, op, psType2); // image op scalar + } else if (dim2 == PS_DIMEN_VECTOR || dim2 == PS_DIMEN_TRANSV) { + BINARY_OP(IMAGE, VECTOR, out, psType1, op, psType2); // image op vector + } else if (dim2 == PS_DIMEN_IMAGE) { + BINARY_OP(IMAGE, IMAGE, out, psType1, op, psType2); // image op image + } else { + psError(PS_ERR_BAD_PARAMETER_TYPE, true, + PS_ERRORTEXT_psMatrix_DIMEN_INVALID, + "in2",dim2); + psBinaryOp_EXIT; + } + } else { + psError(PS_ERR_BAD_PARAMETER_TYPE, true, + PS_ERRORTEXT_psMatrix_DIMEN_INVALID, + "in1",dim1); + psBinaryOp_EXIT; + } + + // Automtically free psScalar types, since they are usually allocated in the argument list when this + // function is called, provided that the input is not the output. + if(psType1->dimen==PS_DIMEN_SCALAR && in1!=out) { + psFree(in1); + } + + if(psType2->dimen==PS_DIMEN_SCALAR && in2!=out) { + psFree(in2); + } + + return out; +} Index: trunk/psLib/src/math/psBinaryOp.h =================================================================== --- trunk/psLib/src/math/psBinaryOp.h (revision 3182) +++ trunk/psLib/src/math/psBinaryOp.h (revision 3182) @@ -0,0 +1,67 @@ + +/** @file psBinaryOp.h + * + * @brief Provides binary functions for simple matrix and vector element operations. Functions + * include: + * + * Addition (+) + * Subtraction (-) + * Multiplication (*) + * Division (/) + * Power (^) + * Minimum (min) + * Maximum (max) + * Absolute value (abs) + * Exponent (exp) + * Natural Log (ln) + * Power of 10 (ten) + * Log (log) + * Sine (sin or dsin) + * Cosine (cos or dcos) + * Tangent (tan or dtan) + * Arcsine (asin or dasin) + * Arccosine (acos or dacos) + * Arctan (atan or datan) + * + * Currently only vector-vector and image-image binary operations are supported. + * + * @ingroup MatrixArithmetic + * + * @author Ross Harman, MHPCC + * @author Robert DeSonia, MHPCC + * + * @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ + * + * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii + */ + +#ifndef PSBINARY_OP_H +#define PSBINARY_OP_H + +/// @addtogroup MatrixArithmetic +/// @{ + +/** Perform simple binary arithmetic with images or vectors + * + * Performs addition, subtraction, multiplication, division, power, minumum, and maximum arithmetic + * operations with images and vectors. Uses the form: + * + * out = in1 op in2, + * + * Where op is: "=", "+", "-", "*", "/", "^", "min", or "max" + * + * This function only supports vector-vector or image-image operations. + * + * @return psType* : Pointer to either psImage or psVector. + */ +psType* psBinaryOp( + psPtr out, ///< Output type, either psImage or psVector. + psPtr in1, ///< First input, either psImage or psVector. + char *op, ///< Operator. + psPtr in2 ///< Second input, either psImage or psVector. +); + +/// @} + +#endif Index: trunk/psLib/src/math/psConstants.h =================================================================== --- trunk/psLib/src/math/psConstants.h (revision 3171) +++ trunk/psLib/src/math/psConstants.h (revision 3182) @@ -6,6 +6,6 @@ * @author GLG, MHPCC * - * @version $Revision: 1.53 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-01-10 19:47:11 $ + * @version $Revision: 1.54 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -37,10 +37,10 @@ #define PS_COT(X) (1.0 / atan(X)) -#define DEG_TO_RAD(DEGREES) ((DEGREES) * M_PI / 180.0) -#define MIN_TO_RAD(MINUTES) ((MINUTES) * M_PI / (180.0 * 60.0)) -#define SEC_TO_RAD(SECONDS) ((SECONDS) * M_PI / (180.0 * 60.0 * 60.0)) -#define RAD_TO_DEG(RADIANS) ((RADIANS) * 180.0 / M_PI) -#define RAD_TO_MIN(RADIANS) ((RADIANS) * 180.0 * 60.0 / M_PI) -#define RAD_TO_SEC(RADIANS) ((RADIANS) * 180.0 * 60.0 * 60.0 / M_PI) +#define DEG_TO_RAD(DEGREES) ((DEGREES) * PS_PI / 180.0) +#define MIN_TO_RAD(MINUTES) ((MINUTES) * PS_PI / (180.0 * 60.0)) +#define SEC_TO_RAD(SECONDS) ((SECONDS) * PS_PI / (180.0 * 60.0 * 60.0)) +#define RAD_TO_DEG(RADIANS) ((RADIANS) * 180.0 / PS_PI) +#define RAD_TO_MIN(RADIANS) ((RADIANS) * 180.0 * 60.0 / PS_PI) +#define RAD_TO_SEC(RADIANS) ((RADIANS) * 180.0 * 60.0 * 60.0 / PS_PI) /***************************************************************************** Index: trunk/psLib/src/math/psMinimize.c =================================================================== --- trunk/psLib/src/math/psMinimize.c (revision 3171) +++ trunk/psLib/src/math/psMinimize.c (revision 3182) @@ -9,6 +9,6 @@ * @author GLG, MHPCC * - * @version $Revision: 1.104 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-02-03 00:54:10 $ + * @version $Revision: 1.105 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -538,5 +538,5 @@ psF64 flux = normalization * exp(-( u*u/(2.0 * sigmaX * sigmaX) + v*v/(2.0 * sigmaY * sigmaY)))/ - (2.0 * M_PI * sigmaX * sigmaY); + (2.0 * PS_PI * sigmaX * sigmaY); out->data.F32[i] = flux; @@ -865,5 +865,5 @@ for (psS32 i=0;idata.F64[k] * - cos(M_PI * ((psF32) j) * (0.5 + ((psF32) k)) / ((psF32) n)); + cos(PS_PI * ((psF32) j) * (0.5 + ((psF32) k)) / ((psF32) n)); } Index: trunk/psLib/src/math/psPolynomial.c =================================================================== --- trunk/psLib/src/math/psPolynomial.c (revision 3171) +++ trunk/psLib/src/math/psPolynomial.c (revision 3182) @@ -7,6 +7,6 @@ * polynomials. It also contains a Gaussian functions. * - * @version $Revision: 1.87 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-02-09 20:56:17 $ + * @version $Revision: 1.88 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -975,5 +975,5 @@ if (normal == true) { - tmp = 1.0 / PS_SQRT_F32(2.0 * M_PI * (sigma * sigma)); + tmp = 1.0 / PS_SQRT_F32(2.0 * PS_PI * (sigma * sigma)); } Index: trunk/psLib/src/math/psSpline.c =================================================================== --- trunk/psLib/src/math/psSpline.c (revision 3171) +++ trunk/psLib/src/math/psSpline.c (revision 3182) @@ -7,6 +7,6 @@ * polynomials. It also contains a Gaussian functions. * - * @version $Revision: 1.87 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-02-09 20:56:17 $ + * @version $Revision: 1.88 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -975,5 +975,5 @@ if (normal == true) { - tmp = 1.0 / PS_SQRT_F32(2.0 * M_PI * (sigma * sigma)); + tmp = 1.0 / PS_SQRT_F32(2.0 * PS_PI * (sigma * sigma)); } Index: trunk/psLib/src/math/psUnaryOp.c =================================================================== --- trunk/psLib/src/math/psUnaryOp.c (revision 3182) +++ trunk/psLib/src/math/psUnaryOp.c (revision 3182) @@ -0,0 +1,387 @@ +/** @file psUnary.c + * + * @brief Provides unary functions for simple matrix and vector element operations. Functions + * include: + * + * Addition (+) + * Subtraction (-) + * Multiplication (*) + * Division (/) + * Power (^) + * Minimum (min) + * Maximum (max) + * Absolute value (abs) + * Exponent (exp) + * Natural Log (ln) + * Power of 10 (ten) + * Log (log) + * Sine (sin or dsin) + * Cosine (cos or dcos) + * Tangent (tan or dtan) + * Arcsine (asin or dasin) + * Arccosine (acos or dacos) + * Arctan (atan or datan) + * + * Currently only vector-vector and image-image binary operations are supported. + * + * @ingroup MatrixArithmetic + * + * @author Ross Harman, MHPCC + * @author Robert DeSonia, MHPCC + * + * @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ + * + * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii + */ + +/****************************************************************************** + * INCLUDE FILES * + ******************************************************************************/ +#include +#include +#include + +#include "psMemory.h" +#include "psError.h" +#include "psImage.h" +#include "psVector.h" +#include "psScalar.h" +#include "psLogMsg.h" +#include "psConstants.h" +#include "psDataManipErrors.h" + +/***************************************************************************** + * FUNCTION IMPLEMENTATION - LOCAL * + *****************************************************************************/ + +// Conversion for degrees to radians +#define D2R 0.01745329252111111 /* PI/180 */ + +// Conversion for radians to degrees +#define R2D 57.29577950924861 /* 180.0/PI */ + + +// Unary SCALAR operations +#define SCALAR(OUT,IN,OP,TYPE) \ +{ \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + o = &((psScalar* )OUT)->data.TYPE; \ + i1 = &((psScalar* )IN)->data.TYPE; \ + *o = OP; \ +} + +// Unary IMAGE operations +#define VECTOR(OUT,IN,OP,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 nIn = 0; \ + psS32 nOut = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + nIn = ((psVector* )IN)->n; \ + nOut = ((psVector* )OUT)->n; \ + if(nIn != nOut) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_COUNT_DIFFERS, \ + nIn, nOut); \ + if (OUT != IN) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + o = ((psVector* )OUT)->data.TYPE; \ + i1 = ((psVector* )IN)->data.TYPE; \ + for(i = 0; i < nIn; i++, o++, i1++) { \ + *o = OP; \ + } \ +} + +// Unary IMAGE operations +#define IMAGE(OUT,IN,OP,TYPE) \ +{ \ + psS32 i = 0; \ + psS32 j = 0; \ + psS32 numRowsIn = 0; \ + psS32 numColsIn = 0; \ + psS32 numRowsOut = 0; \ + psS32 numColsOut = 0; \ + ps##TYPE *o = NULL; \ + ps##TYPE *i1 = NULL; \ + numRowsIn = ((psImage* )IN)->numRows; \ + numColsIn = ((psImage* )IN)->numCols; \ + numRowsOut = ((psImage* )OUT)->numRows; \ + numColsOut = ((psImage* )OUT)->numCols; \ + if(numRowsIn!=numRowsOut || numColsIn!=numColsOut) { \ + psError(PS_ERR_BAD_PARAMETER_SIZE, true, \ + PS_ERRORTEXT_psMatrix_IMAGE_SIZE_DIFFERS, \ + numColsIn, numRowsIn, numColsOut, numRowsOut); \ + if (OUT != IN) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ + for(j = 0; j < numRowsIn; j++) { \ + o = ((psImage* )OUT)->data.TYPE[j]; \ + i1 = ((psImage* )IN)->data.TYPE[j]; \ + for(i = 0; i < numColsIn; i++, o++, i1++) { \ + *o = OP; \ + } \ + } \ +} + +// Preprocessor macro function to create arithmetic function based on input type +#define UNARY_TYPE(DIM,OUT,IN,OP) \ +switch (IN->type) { \ +case PS_TYPE_S32: \ + DIM(OUT,IN,OP,S32); \ + break; \ +case PS_TYPE_F32: \ + DIM(OUT,IN,OP,F32); \ + break; \ +case PS_TYPE_F64: \ + DIM(OUT,IN,OP,F64); \ + break; \ +case PS_TYPE_C32: \ + 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: { \ + char* strType; \ + PS_TYPE_NAME(strType, IN->type); \ + psError(PS_ERR_BAD_PARAMETER_TYPE, true, \ + PS_ERRORTEXT_psMatrix_TYPE_MISMATCH, \ + strType); \ + if (OUT != IN) { \ + psFree(OUT); \ + } \ + return NULL; \ + } \ +} + +// Preprocessor macro function to create arithmetic function operation name. Functions below that add +// FLT_EPSILON are done so to align results with a 64 bit computing architecture +#define UNARY_OP(DIM,OUT,IN,OP) \ +if(!strncmp(OP, "abs", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,cabs(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,fabs(*i1)); \ + } \ +} else if(!strncmp(OP, "exp", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,cexp(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,exp(*i1)); \ + } \ +} else if(!strncmp(OP, "ln", 2)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,clog(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,log(*i1)); \ + } \ +} else if(!strncmp(OP, "ten", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,cpow(10.0,*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,pow(10.0,*i1)); \ + } \ +} else if(!strncmp(OP, "log", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,clog(*i1)/log(10.0)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,log10(*i1)); \ + } \ +} else if(!strncmp(OP, "sin", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,csin(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,sin(*i1)); \ + } \ +} else if(!strncmp(OP, "dsin", 4)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,csin(*i1*D2R)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,sin(*i1*D2R)); \ + } \ +} else if(!strncmp(OP, "cos", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,ccos(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,cos(*i1)); \ + } \ +} else if(!strncmp(OP, "dcos", 4)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,ccos(*i1*D2R)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,cos(*i1*D2R)); \ + } \ +} else if(!strncmp(OP, "tan", 3)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,ctan(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,tan(*i1)); \ + } \ +} else if(!strncmp(OP, "dtan", 4)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,ctan(*i1*D2R)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,tan(*i1*D2R)); \ + } \ +} else if(!strncmp(OP, "asin", 4)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,casin(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,asin(*i1)); \ + } \ +} else if(!strncmp(OP, "dasin", 5)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,R2D*casin(*i1)); \ + } else if(PS_IS_PSELEMTYPE_INT(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,(R2D*asin(*i1))); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,(R2D*asin(*i1))); \ + } \ +} else if(!strncmp(OP, "acos", 4)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,cacos(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,acos(*i1)); \ + } \ +} else if(!strncmp(OP, "dacos", 5)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,R2D*cacos(*i1)); \ + } else if(PS_IS_PSELEMTYPE_INT(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,(R2D*acos(*i1))); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,R2D*acos(*i1)); \ + } \ +} else if(!strncmp(OP, "atan", 4)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,catan(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,atan(*i1)); \ + } \ +} else if(!strncmp(OP, "datan", 5)) { \ + if(PS_IS_PSELEMTYPE_COMPLEX(IN->type)) { \ + UNARY_TYPE(DIM,OUT,IN,R2D*catan(*i1)); \ + } else { \ + UNARY_TYPE(DIM,OUT,IN,R2D*atan(*i1)); \ + } \ +} else { \ + psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ + PS_ERRORTEXT_psMatrix_OPERATION_UNSUPPORTED, \ + OP); \ + if (OUT != IN) { \ + psFree(OUT); \ + } \ + return NULL; \ +} + +psPtr psUnaryOp(psPtr out, psPtr in, char *op) +{ + #define psUnaryOp_EXIT { \ + if (out != in) { \ + psFree(out); \ + } \ + return NULL; \ + } + + psType* psTypeIn = (psType* ) in; + + PS_PTR_CHECK_NULL_GENERAL(in, psUnaryOp_EXIT); + PS_PTR_CHECK_NULL_GENERAL(op, psUnaryOp_EXIT); + + psDimen dimIn = psTypeIn->dimen; + psElemType elTypeIn = psTypeIn->type; + + 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) { + psError(PS_ERR_BAD_PARAMETER_SIZE, true, + PS_ERRORTEXT_psMatrix_VECTOR_EMPTY); + psUnaryOp_EXIT; + } + + out = psVectorRecycle(out, + ((psVector*)in)->n, + elTypeIn); + if (out == NULL) { + psError(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_VECTOR_NOT_CREATED); + psUnaryOp_EXIT; + } + + UNARY_OP(VECTOR, out, psTypeIn, op); // vector + break; + case PS_DIMEN_IMAGE: + if (((psImage* ) in)->numCols == 0 || ((psImage* ) in)->numRows == 0) { + psError(PS_ERR_BAD_PARAMETER_SIZE, true, + PS_ERRORTEXT_psMatrix_IMAGE_EMPTY); + psUnaryOp_EXIT; + } + + out = psImageRecycle(out, + ((psImage*)in)->numCols, + ((psImage*)in)->numRows, + elTypeIn); + if (out == NULL) { + psError(PS_ERR_UNKNOWN, false, + PS_ERRORTEXT_psMatrix_OUTPUT_IMAGE_NOT_CREATED); + psUnaryOp_EXIT; + } + + UNARY_OP(IMAGE, out, psTypeIn, op); // image + break; + default: + if (out != in) { + psFree(out); + } + psError(PS_ERR_BAD_PARAMETER_SIZE, true, + PS_ERRORTEXT_psMatrix_DIMEN_INVALID, + "in", dimIn); + psUnaryOp_EXIT; + } + + // Automtically free psScalar types, since they are usually allocated in the argument list when this + // function is called, provided that the input is not the output. + if(psTypeIn->dimen==PS_DIMEN_SCALAR && in!=out) { + psFree(in); + } + + return out; +} Index: trunk/psLib/src/math/psUnaryOp.h =================================================================== --- trunk/psLib/src/math/psUnaryOp.h (revision 3182) +++ trunk/psLib/src/math/psUnaryOp.h (revision 3182) @@ -0,0 +1,69 @@ +/** @file psUnaryOp.h + * + * @brief Provides unary functions for simple matrix and vector element operations. Functions + * include: + * + * Addition (+) + * Subtraction (-) + * Multiplication (*) + * Division (/) + * Power (^) + * Minimum (min) + * Maximum (max) + * Absolute value (abs) + * Exponent (exp) + * Natural Log (ln) + * Power of 10 (ten) + * Log (log) + * Sine (sin or dsin) + * Cosine (cos or dcos) + * Tangent (tan or dtan) + * Arcsine (asin or dasin) + * Arccosine (acos or dacos) + * Arctan (atan or datan) + * + * Currently only vector-vector and image-image binary operations are supported. + * + * @ingroup MatrixArithmetic + * + * @author Ross Harman, MHPCC + * @author Robert DeSonia, MHPCC + * + * @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:41 $ + * + * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii + */ + +#ifndef PSUNARY_OP_H +#define PSUNARY_OP_H + +/// @addtogroup MatrixArithmetic +/// @{ + + +/** Perform simple unary arithmetic with images or vectors + * + * Performs absolute value, exponent, natural log, power of 10, log, sine, cosine, tangent, arcsine, + * arccosine, or arctan. operations with images and vectors. Uses the form: + * + * out = op(in), + * + * Where op is: "abs", "exp", "ln", "ten", "log", "sin", "cos", "tan" "asin", "acos", "atan", "dsin", + * "dcos", dtan", "dasin", "dacos", or "datan". + * + * Trigometric Operations with "d" prefix use units of degrees. Those without are in radians. + * + * This function only supports vector-vector or image-image opertions. + * + * @return psType* : Pointer to either psImage or psVector. + */ +psType* psUnaryOp( + psPtr out, ///< Output type, either psImage or psVector. + psPtr in, ///< Input, either psImage or psVector. + char *op ///< Operator. +); + +/// @} + +#endif Index: trunk/psLib/src/psTest.c =================================================================== --- trunk/psLib/src/psTest.c (revision 3171) +++ trunk/psLib/src/psTest.c (revision 3182) @@ -11,7 +11,7 @@ // - +#include +#include #include -#include #include #include Index: trunk/psLib/src/pslib.h =================================================================== --- trunk/psLib/src/pslib.h (revision 3171) +++ trunk/psLib/src/pslib.h (revision 3182) @@ -9,6 +9,6 @@ * @author Eric Van Alst, MHPCC * -* @version $Revision: 1.30 $ $Name: not supported by cvs2svn $ -* @date $Date: 2005-02-03 00:54:10 $ +* @version $Revision: 1.31 $ $Name: not supported by cvs2svn $ +* @date $Date: 2005-02-10 02:36:41 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -126,7 +126,8 @@ #include "psMatrix.h" -/// @defgroup MatrixVectorArithmetic Matrix Vector Arithmetic Operations +/// @defgroup MatrixArithmetic Matrix Arithmetic Operations /// @ingroup DataManip -#include "psMatrixVectorArithmetic.h" +#include "psBinaryOp.h" +#include "psUnaryOp.h" #include "psRandom.h" Index: trunk/psLib/src/sys/psTrace.c =================================================================== --- trunk/psLib/src/sys/psTrace.c (revision 3171) +++ trunk/psLib/src/sys/psTrace.c (revision 3182) @@ -9,6 +9,6 @@ * @author GLG, MHPCC * - * @version $Revision: 1.41 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-02-03 00:54:11 $ + * @version $Revision: 1.42 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:42 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -40,4 +40,5 @@ #ifndef PS_NO_TRACE +#include #include #include @@ -146,7 +147,4 @@ /***************************************************************************** componentAdd(): Adds the component named "addNodeName" to the root tree. - -NOTE: replace the call to strsep() with a call to strtok(), which conforms -to ANSI-C. *****************************************************************************/ static psBool componentAdd(const char *addNodeName, psS32 level) @@ -180,5 +178,5 @@ strcpy(name, addNodeName); - pname = &name[1]; + pname = name+1; // Iterate through the components of addNodeName. Strip off the first // component of the name, find that in the root tree, or add it if it @@ -186,5 +184,10 @@ while (pname != NULL) { - firstComponent = strsep(&pname, "."); + firstComponent = pname; + pname = strtok(pname, "."); + if (pname != NULL) { + *pname = '\0'; + pname++; + } nodeExists = 0; for (i = 0; i < currentNode->n; i++) { @@ -324,5 +327,10 @@ pname = &name[1]; while (pname != NULL) { - firstComponent = strsep(&pname, "."); + firstComponent = pname; + pname = strtok(pname, "."); + if (pname != NULL) { + *pname = '\0'; + pname++; + } for (i = 0; i < currentNode->n; i++) { if (NULL == currentNode->subcomp[i]) { Index: trunk/psLib/src/sysUtils/psTrace.c =================================================================== --- trunk/psLib/src/sysUtils/psTrace.c (revision 3171) +++ trunk/psLib/src/sysUtils/psTrace.c (revision 3182) @@ -9,6 +9,6 @@ * @author GLG, MHPCC * - * @version $Revision: 1.41 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-02-03 00:54:11 $ + * @version $Revision: 1.42 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-02-10 02:36:42 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -40,4 +40,5 @@ #ifndef PS_NO_TRACE +#include #include #include @@ -146,7 +147,4 @@ /***************************************************************************** componentAdd(): Adds the component named "addNodeName" to the root tree. - -NOTE: replace the call to strsep() with a call to strtok(), which conforms -to ANSI-C. *****************************************************************************/ static psBool componentAdd(const char *addNodeName, psS32 level) @@ -180,5 +178,5 @@ strcpy(name, addNodeName); - pname = &name[1]; + pname = name+1; // Iterate through the components of addNodeName. Strip off the first // component of the name, find that in the root tree, or add it if it @@ -186,5 +184,10 @@ while (pname != NULL) { - firstComponent = strsep(&pname, "."); + firstComponent = pname; + pname = strtok(pname, "."); + if (pname != NULL) { + *pname = '\0'; + pname++; + } nodeExists = 0; for (i = 0; i < currentNode->n; i++) { @@ -324,5 +327,10 @@ pname = &name[1]; while (pname != NULL) { - firstComponent = strsep(&pname, "."); + firstComponent = pname; + pname = strtok(pname, "."); + if (pname != NULL) { + *pname = '\0'; + pname++; + } for (i = 0; i < currentNode->n; i++) { if (NULL == currentNode->subcomp[i]) {