Index: trunk/psLib/src/image/psImage.c =================================================================== --- trunk/psLib/src/image/psImage.c (revision 1894) +++ trunk/psLib/src/image/psImage.c (revision 1897) @@ -10,6 +10,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.45 $ $Name: not supported by cvs2svn $ - * @date $Date: 2004-09-21 23:15:04 $ + * @version $Revision: 1.46 $ $Name: not supported by cvs2svn $ + * @date $Date: 2004-09-25 02:06:12 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -199,4 +199,6 @@ float x, float y, + const psImage* mask, + unsigned int maskVal, psC64 unexposedValue, psImageInterpolateMode mode) @@ -214,9 +216,19 @@ switch (mode) { \ case PS_INTERPOLATE_FLAT: \ - return p_psImagePixelInterpolateFLAT_##TYPE(input, x, y, \ + return p_psImagePixelInterpolateFLAT_##TYPE( \ + input, \ + x, \ + y, \ + mask, \ + maskVal, \ unexposedValue); \ break; \ case PS_INTERPOLATE_BILINEAR: \ - return p_psImagePixelInterpolateBILINEAR_##TYPE(input,x,y, \ + return p_psImagePixelInterpolateBILINEAR_##TYPE( \ + input, \ + x, \ + y, \ + mask, \ + maskVal, \ unexposedValue); \ break; \ @@ -259,4 +271,6 @@ float x, \ float y, \ + const psImage* mask, \ + unsigned int maskVal, \ psF64 unexposedValue) \ { \ @@ -269,5 +283,7 @@ (intX > lastX) || \ (intY < 0) || \ - (intY > lastY)) { \ + (intY > lastY) || \ + ( (mask!=NULL) && \ + ((mask->data.PS_TYPE_MASK_DATA[intY][intX] & maskVal) != 0) ) ) { \ return unexposedValue; \ } \ @@ -280,4 +296,6 @@ float x, \ float y, \ + const psImage* mask, \ + unsigned int maskVal, \ psC64 unexposedValue) \ { \ @@ -290,5 +308,7 @@ (intX > lastX) || \ (intY < 0) || \ - (intY > lastY)) { \ + (intY > lastY) || \ + ( (mask!=NULL) && \ + ((mask->data.PS_TYPE_MASK_DATA[intY][intX] & maskVal) != 0) ) ) { \ return unexposedValue; \ } \ @@ -314,109 +334,186 @@ float x, \ float y, \ + const psImage* mask, \ + unsigned int maskVal, \ psF64 unexposedValue) \ { \ double floorX = floor((psF64)(x) - 0.5); \ double floorY = floor((psF64)(y) - 0.5); \ - double fracX = x - 0.5 - floorX; \ - double fracY = y - 0.5 - floorY; \ + psF64 fracX = x - 0.5 - floorX; \ + psF64 fracY = y - 0.5 - floorY; \ int intFloorX = (int) floorX; \ int intFloorY = (int) floorY; \ int lastX = input->numCols - 1; \ int lastY = input->numRows - 1; \ - double rx = 0.0; \ - psF64 pixel = 0.0; \ - ps##TYPE* currentRow; \ - ps##TYPE* nextRow; \ - \ - if ((intFloorX < 0) || \ - (intFloorX > lastX) || \ - (intFloorY < 0) || \ - (intFloorY > lastY)) { \ - return unexposedValue; \ - } \ - \ - currentRow = input->data.TYPE[intFloorY]; \ - if (intFloorY == lastY) { \ - pixel = currentRow[intFloorX]; \ - if (intFloorX < lastX) { \ - pixel+= fracY * ((psF64)currentRow[intFloorX+1] - \ - (psF64)currentRow[intFloorX]); \ - } \ - return(pixel); \ - } \ - nextRow = input->data.TYPE[intFloorY+1]; \ - if (intFloorX == lastX) { \ - pixel = currentRow[intFloorX]; \ - if (intFloorY < lastY) { \ - pixel+= fracX * ((psF64)nextRow[intFloorX] - \ - (psF64)currentRow[intFloorX]); \ - } \ - return(pixel); \ - } \ - \ - rx = currentRow[intFloorX] + \ - fracX * ((psF64)currentRow[intFloorX+1] - \ - (psF64)currentRow[intFloorX]); \ - \ - pixel = rx + fracY * ((psF64)nextRow[intFloorX] + \ - fracX * ((psF64)nextRow[intFloorX+1] - \ - (psF64)nextRow[intFloorX]) - rx); \ - \ - return(pixel); \ -} + ps##TYPE V00; \ + ps##TYPE V01; \ + ps##TYPE V10; \ + ps##TYPE V11; \ + bool valid00 = false; \ + bool valid01 = false; \ + bool valid10 = false; \ + bool valid11 = false; \ + \ + if (intFloorY >= 0 && intFloorY <= lastY) { \ + if (intFloorX >= 0 && intFloorX <= lastX) { \ + V00 = input->data.TYPE[intFloorY][intFloorX]; \ + valid00 = (mask == NULL) || \ + ((mask->data.PS_TYPE_MASK_DATA[intFloorY][intFloorX] & maskVal) == 0); \ + } \ + if (intFloorX >= -1 && intFloorX < lastX) { \ + V10 = input->data.TYPE[intFloorY][intFloorX+1]; \ + valid10 = (mask == NULL) || \ + ((mask->data.PS_TYPE_MASK_DATA[intFloorY][intFloorX+1] & maskVal) == 0); \ + } \ + } \ + if (intFloorY >= -1 && intFloorY < lastY) { \ + if (intFloorX >= 0 && intFloorX <= lastX) { \ + V01 = input->data.TYPE[intFloorY+1][intFloorX]; \ + valid01 = (mask == NULL) || \ + ((mask->data.PS_TYPE_MASK_DATA[intFloorY+1][intFloorX] & maskVal) == 0); \ + } \ + if (intFloorX >= -1 && intFloorX < lastX) { \ + V11 = input->data.TYPE[intFloorY+1][intFloorX+1]; \ + valid11 = (mask == NULL) || \ + ((mask->data.PS_TYPE_MASK_DATA[intFloorY+1][intFloorX+1] & maskVal) == 0); \ + } \ + } \ + \ + /* cover likely case of all pixels being valid more efficiently */ \ + if (valid00 && valid10 && valid01 && valid11) { \ + /* formula from the ADD */ \ + psF64 rx = V00 + (V10-V00)*fracX; \ + return (rx + (V01 + (V11 - V01) * fracX - rx) * fracY); \ + } \ + \ + /* OK, at least one pixel is not valid - need to do it piecemeal */ \ + \ + psF64 V0; \ + bool valid0 = true; \ + if (valid00 && valid10) { \ + V0 = V00*(1-fracX)+V10*fracX; \ + } else if (valid00) { \ + V0 = V00; \ + } else if (valid10) { \ + V0 = V10; \ + } else { \ + valid0 = false; \ + } \ + \ + psF64 V1; \ + bool valid1 = true; \ + if (valid01 && valid11) { \ + V0 = V01*(1-fracX)+V11*fracX; \ + } else if (valid01) { \ + V0 = V01; \ + } else if (valid11) { \ + V0 = V11; \ + } else { \ + valid1 = false; \ + } \ + \ + if (valid0 && valid1) { \ + return ( V0*(1-fracY) + V1*fracY ); \ + } else if (valid0) { \ + return V0; \ + } else if (valid1) { \ + return V1; \ + } \ + \ + return unexposedValue; \ +} + #define PSIMAGE_PIXEL_INTERPOLATE_BILINEAR_COMPLEX(TYPE) \ inline psC64 p_psImagePixelInterpolateBILINEAR_##TYPE(const psImage* input, \ float x, \ float y, \ + const psImage* mask, \ + unsigned int maskVal, \ psC64 unexposedValue) \ { \ double floorX = floor((psF64)(x) - 0.5); \ double floorY = floor((psF64)(y) - 0.5); \ - double fracX = x - 0.5 - floorX; \ - double fracY = y - 0.5 - floorY; \ + psF64 fracX = x - 0.5 - floorX; \ + psF64 fracY = y - 0.5 - floorY; \ int intFloorX = (int) floorX; \ int intFloorY = (int) floorY; \ int lastX = input->numCols - 1; \ int lastY = input->numRows - 1; \ - double rx = 0.0; \ - psC64 pixel = 0.0; \ - ps##TYPE* currentRow; \ - ps##TYPE* nextRow; \ - \ - if ((intFloorX < 0) || \ - (intFloorX > lastX) || \ - (intFloorY < 0) || \ - (intFloorY > lastY)) { \ - return unexposedValue; \ - } \ - \ - currentRow = input->data.TYPE[intFloorY]; \ - if (intFloorY == lastY) { \ - pixel = currentRow[intFloorX]; \ - if (intFloorX < lastX) { \ - pixel+= fracY * (currentRow[intFloorX+1] - \ - currentRow[intFloorX]); \ - } \ - return(pixel); \ - } \ - nextRow = input->data.TYPE[intFloorY+1]; \ - if (intFloorX == lastX) { \ - pixel = currentRow[intFloorX]; \ - if (intFloorY < lastY) { \ - pixel+= fracX * (nextRow[intFloorX] - \ - currentRow[intFloorX]); \ - } \ - return(pixel); \ - } \ - \ - rx = currentRow[intFloorX] + \ - fracX * (currentRow[intFloorX+1] - \ - currentRow[intFloorX]); \ - \ - pixel = rx + fracY * ((psF64)nextRow[intFloorX] + \ - fracX * (nextRow[intFloorX+1] - \ - nextRow[intFloorX]) - rx); \ - \ - return(pixel); \ + ps##TYPE V00; \ + ps##TYPE V01; \ + ps##TYPE V10; \ + ps##TYPE V11; \ + bool valid00 = false; \ + bool valid01 = false; \ + bool valid10 = false; \ + bool valid11 = false; \ + \ + if (intFloorY >= 0 && intFloorY <= lastY) { \ + if (intFloorX >= 0 && intFloorX <= lastX) { \ + V00 = input->data.TYPE[intFloorY][intFloorX]; \ + valid00 = (mask == NULL) || \ + ((mask->data.PS_TYPE_MASK_DATA[intFloorY][intFloorX] & maskVal) == 0); \ + } \ + if (intFloorX >= -1 && intFloorX < lastX) { \ + V10 = input->data.TYPE[intFloorY][intFloorX+1]; \ + valid10 = (mask == NULL) || \ + ((mask->data.PS_TYPE_MASK_DATA[intFloorY][intFloorX+1] & maskVal) == 0); \ + } \ + } \ + if (intFloorY >= -1 && intFloorY < lastY) { \ + if (intFloorX >= 0 && intFloorX <= lastX) { \ + V01 = input->data.TYPE[intFloorY+1][intFloorX]; \ + valid01 = (mask == NULL) || \ + ((mask->data.PS_TYPE_MASK_DATA[intFloorY+1][intFloorX] & maskVal) == 0); \ + } \ + if (intFloorX >= -1 && intFloorX < lastX) { \ + V11 = input->data.TYPE[intFloorY+1][intFloorX+1]; \ + valid11 = (mask == NULL) || \ + ((mask->data.PS_TYPE_MASK_DATA[intFloorY+1][intFloorX+1] & maskVal) == 0); \ + } \ + } \ + \ + /* cover likely case of all pixels being valid */ \ + if (valid00 && valid10 && valid01 && valid11) { \ + /* formula from the ADD */ \ + psC64 rx = V00 + (V10-V00)*fracX; \ + return (rx + (V01 + (V11 - V01) * fracX - rx) * fracY); \ + } \ + \ + /* OK, at least one pixel is not valid - need to do it piecemeal */ \ + \ + psC64 V0; \ + bool valid0 = true; \ + if (valid00 && valid10) { \ + V0 = V00*(1-fracX)+V10*fracX; \ + } else if (valid00) { \ + V0 = V00; \ + } else if (valid10) { \ + V0 = V10; \ + } else { \ + valid0 = false; \ + } \ + \ + psC64 V1; \ + bool valid1 = true; \ + if (valid01 && valid11) { \ + V0 = V01*(1-fracX)+V11*fracX; \ + } else if (valid01) { \ + V0 = V01; \ + } else if (valid11) { \ + V0 = V11; \ + } else { \ + valid1 = false; \ + } \ + \ + if (valid0 && valid1) { \ + return ( V0*(1-fracY) + V1*fracY ); \ + } else if (valid0) { \ + return V0; \ + } else if (valid1) { \ + return V1; \ + } \ + \ + return unexposedValue; \ } Index: trunk/psLib/src/image/psImage.h =================================================================== --- trunk/psLib/src/image/psImage.h (revision 1894) +++ trunk/psLib/src/image/psImage.h (revision 1897) @@ -11,6 +11,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.36 $ $Name: not supported by cvs2svn $ - * @date $Date: 2004-09-21 19:51:41 $ + * @version $Revision: 1.37 $ $Name: not supported by cvs2svn $ + * @date $Date: 2004-09-25 02:06:12 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -118,4 +118,6 @@ float x, ///< column location to derive value of float y, ///< row location ot derive value of + const psImage* mask, ///< if not NULL, the mask of the input image + unsigned int maskVal, ///< the mask value psC64 unexposedValue, ///< return value if x,y location is not in image. psImageInterpolateMode mode ///< interpolation mode @@ -124,28 +126,36 @@ #define p_psImagePixelInterpolateFcns(TYPE) \ inline psF64 p_psImagePixelInterpolateFLAT_##TYPE( \ - const psImage* input, \ - float x, \ - float y, \ - psF64 unexposedValue \ + const psImage* input, /**< input image for interpolation */ \ + float x, /**< column location to derive value of */ \ + float y, /**< row location ot derive value of */ \ + const psImage* mask, /**< if not NULL, the mask of the input image */ \ + unsigned int maskVal, /**< the mask value */ \ + psF64 unexposedValue /**< return value if x,y location is not in image. */ \ ); \ inline psF64 p_psImagePixelInterpolateBILINEAR_##TYPE( \ - const psImage* input, \ - float x, \ - float y, \ - psF64 unexposedValue \ + const psImage* input, /**< input image for interpolation */ \ + float x, /**< column location to derive value of */ \ + float y, /**< row location ot derive value of */ \ + const psImage* mask, /**< if not NULL, the mask of the input image */ \ + unsigned int maskVal, /**< the mask value */ \ + psF64 unexposedValue /**< return value if x,y location is not in image. */ \ ); #define p_psImagePixelInterpolateComplexFcns(TYPE) \ inline psC64 p_psImagePixelInterpolateFLAT_##TYPE( \ - const psImage* input, \ - float x, \ - float y, \ - psC64 unexposedValue \ + const psImage* input, /**< input image for interpolation */ \ + float x, /**< column location to derive value of */ \ + float y, /**< row location ot derive value of */ \ + const psImage* mask, /**< if not NULL, the mask of the input image */ \ + unsigned int maskVal, /**< the mask value */ \ + psC64 unexposedValue /**< return value if x,y location is not in image. */ \ ); \ inline psC64 p_psImagePixelInterpolateBILINEAR_##TYPE( \ - const psImage* input, \ - float x, \ - float y, \ - psC64 unexposedValue \ + const psImage* input, /**< input image for interpolation */ \ + float x, /**< column location to derive value of */ \ + float y, /**< row location ot derive value of */ \ + const psImage* mask, /**< if not NULL, the mask of the input image */ \ + unsigned int maskVal, /**< the mask value */ \ + psC64 unexposedValue /**< return value if x,y location is not in image. */ \ ); @@ -162,4 +172,6 @@ p_psImagePixelInterpolateComplexFcns(C32) p_psImagePixelInterpolateComplexFcns(C64) + /// @} + #endif Index: trunk/psLib/src/image/psImageErrors.dat =================================================================== --- trunk/psLib/src/image/psImageErrors.dat (revision 1894) +++ trunk/psLib/src/image/psImageErrors.dat (revision 1897) @@ -26,4 +26,9 @@ psImage_SLICE_DIRECTION_INVALID Specified slice direction, %d, is invalid. psImage_PARAMETER_OUTOF_TYPERANGE Specified %s value, %g, is outside of psImage type's range (%s: %g to %g). +psImage_nSamples_TOOSMALL Specified number of samples, %d, must be greater than 1 to make a line. +psImage_LINE_NOT_IN_IMAGE Specified line, (%f,%f)->(%f,%f), does not entirely lie in psImage's boundaries, [0:%d,0:%d]. +psImage_RADII_VECTOR_NULL Specified radii vector can not be NULL. +psImage_CENTER_NOT_IN_IMAGE Specified center, (%d,%d), is outside of the psImage boundaries, [0:%d,0:%d]. +psImage_RADII_VECTOR_TOOSMALL Input radii vector size, %d, can not be less than 2. # psImageFFT_IMAGE_TYPE_UNSUPPORTED Input psImage type (%s) is not supported. Valid image types are psF32 and psC32. Index: trunk/psLib/src/image/psImageErrors.h =================================================================== --- trunk/psLib/src/image/psImageErrors.h (revision 1894) +++ trunk/psLib/src/image/psImageErrors.h (revision 1897) @@ -7,6 +7,6 @@ * @author Robert DeSonia, MHPCC * - * @version $Revision: 1.4 $ $Name: not supported by cvs2svn $ - * @date $Date: 2004-09-23 18:31:49 $ + * @version $Revision: 1.5 $ $Name: not supported by cvs2svn $ + * @date $Date: 2004-09-25 02:06:12 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -48,4 +48,9 @@ #define PS_ERRORTEXT_psImage_SLICE_DIRECTION_INVALID "Specified slice direction, %d, is invalid." #define PS_ERRORTEXT_psImage_PARAMETER_OUTOF_TYPERANGE "Specified %s value, %g, is outside of psImage type's range (%s: %g to %g)." +#define PS_ERRORTEXT_psImage_nSamples_TOOSMALL "Specified number of samples, %d, must be greater than 1 to make a line." +#define PS_ERRORTEXT_psImage_LINE_NOT_IN_IMAGE "Specified line, (%f,%f)->(%f,%f), does not entirely lie in psImage's boundaries, [0:%d,0:%d]." +#define PS_ERRORTEXT_psImage_RADII_VECTOR_NULL "Specified radii vector can not be NULL." +#define PS_ERRORTEXT_psImage_CENTER_NOT_IN_IMAGE "Specified center, (%d,%d), is outside of the psImage boundaries, [0:%d,0:%d]." +#define PS_ERRORTEXT_psImage_RADII_VECTOR_TOOSMALL "Input radii vector size, %d, can not be less than 2." #define PS_ERRORTEXT_psImageFFT_IMAGE_TYPE_UNSUPPORTED "Input psImage type (%s) is not supported. Valid image types are psF32 and psC32." #define PS_ERRORTEXT_psImageFFT_REVERSE_NOT_COMPLEX "Input psImage (%s) is not complex. Reverse FFT operation requires a complex psImage input." Index: trunk/psLib/src/image/psImageExtraction.c =================================================================== --- trunk/psLib/src/image/psImageExtraction.c (revision 1894) +++ trunk/psLib/src/image/psImageExtraction.c (revision 1897) @@ -10,6 +10,6 @@ * @author Robert DeSonia, MHPCC * -* @version $Revision: 1.14 $ $Name: not supported by cvs2svn $ -* @date $Date: 2004-09-23 18:31:49 $ +* @version $Revision: 1.15 $ $Name: not supported by cvs2svn $ +* @date $Date: 2004-09-25 02:06:12 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -609,6 +609,7 @@ psVector* psImageCut(psVector* out, - psVector* coords, - const psImage* input, + psVector* cutCols, + psVector* cutRows, + const psImage* in, const psImage* restrict mask, unsigned int maskVal, @@ -620,10 +621,125 @@ psImageInterpolateMode mode) { - - return NULL; + if (in == NULL || in->data.V == NULL) { + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageCut", + PS_ERR_BAD_PARAMETER_NULL, true, + PS_ERRORTEXT_psImage_IMAGE_NULL); + psFree(out); + return NULL; + } + int numCols = in->numCols; + int numRows = in->numRows; + + if (nSamples < 2) { + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageCut", + PS_ERR_BAD_PARAMETER_VALUE, true, + PS_ERRORTEXT_psImage_nSamples_TOOSMALL, + nSamples); + psFree(out); + return NULL; + } + + if (startCol < 0 || startCol >= numCols || + startRow < 0 || startRow >= numRows || + endCol < 0 || endCol >= numCols || + endRow < 0 || endRow >= numRows) { + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageCut", + PS_ERR_BAD_PARAMETER_VALUE, true, + PS_ERRORTEXT_psImage_LINE_NOT_IN_IMAGE, + startCol,startRow,endCol,endRow, + numCols,numRows); + psFree(out); + } + + if (mask != NULL) { + if (numRows != mask->numRows || numCols != mask->numCols) { + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageCut", + PS_ERR_BAD_PARAMETER_VALUE, true, + PS_ERRORTEXT_psImage_IMAGE_MASK_SIZE, + mask->numCols,mask->numRows, + numCols, numRows); + psFree(out); + } + if (mask->type.type != PS_TYPE_MASK) { + char* typeStr; + PS_TYPE_NAME(typeStr,mask->type.type); + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageCut", + PS_ERR_BAD_PARAMETER_TYPE, true, + PS_ERRORTEXT_psImage_IMAGE_MASK_TYPE, + typeStr, PS_TYPE_MASK_NAME); + psFree(out); + } + } + + //resize the vectors for the coordinate output + psF32* cutColsData = NULL; + psF32* cutRowsData = NULL; + if (cutCols != NULL) { + (void)psVectorRecycle(cutCols, nSamples, PS_TYPE_F32); + cutColsData = cutCols->data.F32; + } + if (cutRows != NULL) { + (void)psVectorRecycle(cutRows, nSamples, PS_TYPE_F32); + cutRowsData = cutRows->data.F32; + } + + out = psVectorRecycle(out,nSamples,in->type.type); + + float dX = (endCol - startCol) / (nSamples-1); + float dY = (endRow - startRow) / (nSamples-1); + + float x = startCol; + float y = startRow; + + #define LINEAR_CUT_CASE(TYPE) \ +case PS_TYPE_##TYPE: { \ + ps##TYPE* outData = out->data.TYPE; \ + for (int i = 0; i < nSamples; i++) { \ + /* store off the location of the sample. */ \ + if (cutColsData != NULL) { \ + cutColsData[i] = x; \ + } \ + if (cutRowsData != NULL) { \ + cutRowsData[i] = y; \ + } \ + outData[i] = psImagePixelInterpolate(in,x,y,mask,maskVal,0,mode); \ + x += dX; \ + y += dY; \ + } \ + } \ + break; + + + switch (in->type.type) { + LINEAR_CUT_CASE(U8); + LINEAR_CUT_CASE(U16); + LINEAR_CUT_CASE(U32); + LINEAR_CUT_CASE(U64); + LINEAR_CUT_CASE(S8); + LINEAR_CUT_CASE(S16); + LINEAR_CUT_CASE(S32); + LINEAR_CUT_CASE(S64); + LINEAR_CUT_CASE(F32); + LINEAR_CUT_CASE(F64); + LINEAR_CUT_CASE(C32); + LINEAR_CUT_CASE(C64); + + default: { + char* typeStr; + PS_TYPE_NAME(typeStr,in->type.type); + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageCut", + PS_ERR_BAD_PARAMETER_TYPE, true, + PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, + typeStr); + psFree(out); + out = NULL; + } + } + + return out; } psVector* psImageRadialCut(psVector* out, - const psImage* input, + const psImage* in, const psImage* restrict mask, unsigned int maskVal, @@ -633,5 +749,194 @@ const psStats* stats) { - - return NULL; + double statVal; + + /* check the parameters */ + + if (in == NULL || in->data.V == NULL) { + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageRadialCut", + PS_ERR_BAD_PARAMETER_NULL, true, + PS_ERRORTEXT_psImage_IMAGE_NULL); + psFree(out); + return NULL; + } + int numCols = in->numCols; + int numRows = in->numRows; + + if (mask != NULL) { + if (numRows != mask->numRows || numCols != mask->numCols) { + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageRadialCut", + PS_ERR_BAD_PARAMETER_VALUE, true, + PS_ERRORTEXT_psImage_IMAGE_MASK_SIZE, + mask->numCols,mask->numRows, + numCols, numRows); + psFree(out); + } + if (mask->type.type != PS_TYPE_MASK) { + char* typeStr; + PS_TYPE_NAME(typeStr,mask->type.type); + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageRadialCut", + PS_ERR_BAD_PARAMETER_TYPE, true, + PS_ERRORTEXT_psImage_IMAGE_MASK_TYPE, + typeStr, PS_TYPE_MASK_NAME); + psFree(out); + } + } + + if (centerCol < 0 || centerCol >= numCols || + centerRow < 0 || centerRow >= numRows) { + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageRadialCut", + PS_ERR_BAD_PARAMETER_VALUE, true, + PS_ERRORTEXT_psImage_CENTER_NOT_IN_IMAGE, + centerCol, centerRow, + numCols, numRows); + psFree(out); + } + + if (radii == NULL) { + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageRadialCut", + PS_ERR_BAD_PARAMETER_NULL, true, + PS_ERRORTEXT_psImage_RADII_VECTOR_NULL); + psFree(out); + return NULL; + } + + if (radii->n < 2) { + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageRadialCut", + PS_ERR_BAD_PARAMETER_VALUE, true, + PS_ERRORTEXT_psImage_RADII_VECTOR_TOOSMALL, + radii->n); + psFree(out); + return NULL; + } + + if (stats == NULL) { + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageRadialCut", + PS_ERR_BAD_PARAMETER_NULL, true, + PS_ERRORTEXT_psImage_STAT_NULL); + psFree(out); + return NULL; + } + + // verify that the stats struct specifies a + // single stats operation + if (p_psGetStatValue(stats, &statVal) == false) { + psErrorMsg(PS_ERRORNAME_DOMAIN "psImageRadialCut", + PS_ERR_BAD_PARAMETER_VALUE, false, + PS_ERRORTEXT_psImage_BAD_STAT); + psFree(out); + return NULL; + } + + /* completed checking the parameters */ + + // size the output vector to proper size. + int numOut = radii->n - 1; + out = psVectorRecycle(out, numOut, PS_TYPE_F64); + psF64* outData = out->data.F64; + + // sort the radii by value + psVector* rSqVec = psVectorCopy(NULL, radii, PS_TYPE_F32); + psVectorSort(rSqVec,rSqVec); + psF32* rSq = rSqVec->data.F32; + + int startRow = centerRow - rSq[numOut]; + int endRow = centerRow + rSq[numOut]; + int startCol = centerCol - rSq[numOut]; + int endCol = centerCol + rSq[numOut]; + + if (startRow < 0) { + startRow = 0; + } + + if (startCol < 0) { + startCol = 0; + } + + if (endRow >= numRows) { + endRow = numRows - 1; + } + + if (endCol >= numCols) { + endCol = numCols - 1; + } + + // Square the data + for (int d = 0; d <= numOut; d++) { + rSq[d] *= rSq[d]; + } + + // create temporary vectors for the data binning step + psVector** buffer = psAlloc(sizeof(psVector*)*numOut); + psVector** bufferMask = psAlloc(sizeof(psVector*)*numOut); + for (int lcv = 0; lcv < numOut; lcv++) { + // n.b. alloc enough for the data by making the vectors slightly larger + // than the area of the region of interest. + buffer[lcv] = psVectorAlloc(1+4*(rSq[lcv+1]-rSq[lcv]), + in->type.type); + buffer[lcv]->n = 0; + + bufferMask[lcv] = NULL; + if (mask != NULL) { + bufferMask[lcv] = psVectorAlloc(1+4*(rSq[lcv+1]-rSq[lcv]), + PS_TYPE_MASK); + bufferMask[lcv]->n = 0; + } + } + + float dX; + float dY; + float dist; + for (int row=startRow; row <= endRow; row++) { + psF32* inRow = in->data.F32[row]; + psMaskType* maskRow = NULL; + if (mask != NULL) { + maskRow = mask->data.PS_TYPE_MASK_DATA[row]; + } + for (int col=startCol; col <= endCol; col++) { + dX = centerCol - (float)col - 0.5f; + dY = centerRow - (float)row - 0.5f; + dist = dX*dX+dY*dY; + for (int r = 0; r < numOut;) { + if (rSq[r] < dist && dist < rSq[++r]) { + int n = buffer[r]->n; + if (n == buffer[r]->nalloc) { // in case buffers already full, expand + buffer[r] = psVectorRealloc(n*2, buffer[r]); + if (bufferMask[r] != NULL) { + bufferMask[r] = psVectorRealloc(n*2, bufferMask[r]); + } + } + + buffer[r]->data.F32[n] = inRow[col]; + buffer[r]->n = n+1; + + if (maskRow != NULL) { + bufferMask[r]->data.PS_TYPE_MASK_DATA[n] = maskRow[col]; + bufferMask[r]->n = n+1; + } + + break; + } + } + } + } + + psStats* myStats = psAlloc(sizeof(psStats)); + *myStats = *stats; + + for (int r = 0; r < numOut; r++) { + myStats = psVectorStats(myStats,buffer[r], bufferMask[r],maskVal); + (void)p_psGetStatValue(myStats,&statVal); + outData[r] = statVal; + } + + psFree(myStats); + + for (int lcv = 0; lcv < numOut; lcv++) { + psFree(buffer[lcv]); + psFree(bufferMask[lcv]); + } + psFree(buffer); + psFree(bufferMask); + + return out; } Index: trunk/psLib/src/image/psImageExtraction.h =================================================================== --- trunk/psLib/src/image/psImageExtraction.h (revision 1894) +++ trunk/psLib/src/image/psImageExtraction.h (revision 1897) @@ -10,6 +10,6 @@ * @author Robert DeSonia, MHPCC * -* @version $Revision: 1.12 $ $Name: not supported by cvs2svn $ -* @date $Date: 2004-09-23 18:31:49 $ +* @version $Revision: 1.13 $ $Name: not supported by cvs2svn $ +* @date $Date: 2004-09-25 02:06:12 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -145,5 +145,6 @@ psVector* psImageCut( psVector* out, ///< psVector to recycle, or NULL. - psVector* coords, ///< if not NULL, the calculated coordinates along the slice (output) + psVector* cutCols, ///< if not NULL, the calculated column values along the slice (output) + psVector* cutRows, ///< if not NULL, the calculated row values along the slice (output) const psImage* input, ///< the input image in which to perform the cut const psImage* restrict mask, ///< the mask for the input image. Index: trunk/psLib/src/image/psImageManip.c =================================================================== --- trunk/psLib/src/image/psImageManip.c (revision 1894) +++ trunk/psLib/src/image/psImageManip.c (revision 1897) @@ -10,6 +10,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.22 $ $Name: not supported by cvs2svn $ - * @date $Date: 2004-09-25 00:47:54 $ + * @version $Revision: 1.23 $ $Name: not supported by cvs2svn $ + * @date $Date: 2004-09-25 02:06:12 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -580,5 +580,5 @@ float inRow = (float)row * invScale; \ for (int col=0;coldata.TYPE[y]; \ for (int x = 0; x < outCols; x++) { \ - outRow[x] = p_psImagePixelInterpolate##MODE##_##TYPE(in,inX,inY,unexposedValue); \ + outRow[x] = p_psImagePixelInterpolate##MODE##_##TYPE(in,inX,inY,NULL,0,unexposedValue); \ inX += cosT; \ inY -= sinT; \ @@ -1013,5 +1013,5 @@ for (int col=0;col