Index: trunk/psLib/src/image/.cvsignore =================================================================== --- trunk/psLib/src/image/.cvsignore (revision 4534) +++ (revision ) @@ -1,7 +1,0 @@ -Makefile.in -.deps -.libs -Makefile -*.lo -*.la - Index: trunk/psLib/src/image/Makefile.am =================================================================== --- trunk/psLib/src/image/Makefile.am (revision 4534) +++ (revision ) @@ -1,38 +1,0 @@ -#Makefile for image functions of psLib -# -INCLUDES = \ - -I$(top_srcdir)/src/astronomy \ - -I$(top_srcdir)/src/collections \ - -I$(top_srcdir)/src/dataManip \ - -I$(top_srcdir)/src/dataIO \ - -I$(top_srcdir)/src/sysUtils \ - $(all_includes) - -noinst_LTLIBRARIES = libpslibimage.la - -libpslibimage_la_SOURCES = \ - psImage.c \ - psImagePixelExtract.c \ - psImageStructManip.c \ - psImageGeomManip.c \ - psImagePixelManip.c \ - psImageStats.c \ - psImageFFT.c \ - psImageConvolve.c - -BUILT_SOURCES = psImageErrors.h -EXTRA_DIST = psImageErrors.dat psImageErrors.h image.i - -psImageErrors.h: psImageErrors.dat - $(top_srcdir)/src/psParseErrorCodes --data=$? $@ - -pslibincludedir = $(includedir) -pslibinclude_HEADERS = \ - psImage.h \ - psImagePixelExtract.h \ - psImageStructManip.h \ - psImageGeomManip.h \ - psImagePixelManip.h \ - psImageStats.h \ - psImageFFT.h \ - psImageConvolve.h Index: trunk/psLib/src/image/image.i =================================================================== --- trunk/psLib/src/image/image.i (revision 4534) +++ (revision ) @@ -1,11 +1,0 @@ -/* image headers */ -%include "psImageConvolve.h" -%include "psImageErrors.h" -%include "psImageStructManip.h" -%include "psImagePixelExtract.h" -%include "psImageFFT.h" -%include "psImage.h" -%include "psImageGeomManip.h" -%include "psImagePixelManip.h" -%include "psImageStats.h" - Index: trunk/psLib/src/image/psImage.c =================================================================== --- trunk/psLib/src/image/psImage.c (revision 4534) +++ (revision ) @@ -1,656 +1,0 @@ -/** @file psImage.c - * - * @brief Contains basic image definitions and operations. - * - * This file defines the basic type for an image struct and functions useful - * in manupulating images. - * - * @author Robert DeSonia, MHPCC - * @author Ross Harman, MHPCC - * - * @version $Revision: 1.75 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-07-11 21:38:19 $ - * - * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii - * - * That is the routine used to generate matrices. - */ - -#include -#include - -#include "psMemory.h" -#include "psError.h" -#include "psImage.h" -#include "psString.h" - -#include "psImageErrors.h" - -#define SQUARE(x) ((x)*(x)) -#define MIN(x,y) (((x) > (y)) ? (y) : (x)) -#define MAX(x,y) (((x) > (y)) ? (x) : (y)) - -static void imageFree(psImage* image) -{ - if (image == NULL) { - return; - } - - if (image->parent != NULL) { - psArrayRemove(image->parent->children,image); - image->parent = NULL; - } - - psImageFreeChildren(image); - - psFree(image->rawDataBuffer); - psFree(image->data.V); -} - -psImage* psImageAlloc(int numCols, - int numRows, - psElemType type) -{ - psS32 area = 0; - psS32 elementSize = PSELEMTYPE_SIZEOF(type); // element size in bytes - psS32 rowSize = numCols * elementSize; // row size in bytes. - - area = numCols * numRows; - - if (area < 1) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_AREA_NEGATIVE, - numRows, numCols); - return NULL; - } - - psImage* image = (psImage* ) psAlloc(sizeof(psImage)); - - psMemSetDeallocator(image, (psFreeFunc) imageFree); - - image->data.V = psAlloc(sizeof(psPtr ) * numRows); - - image->rawDataBuffer = psAlloc(area * elementSize); - - // set the row pointers. - image->data.V[0] = image->rawDataBuffer; - for (psS32 i = 1; i < numRows; i++) { - image->data.V[i] = (psPtr )((int8_t *) image->data.V[i - 1] + rowSize); - } - - *(psS32 *)&image->col0 = 0; - *(psS32 *)&image->row0 = 0; - *(psU32 *)&image->numCols = numCols; - *(psU32 *)&image->numRows = numRows; - *(psDimen* ) & image->type.dimen = PS_DIMEN_IMAGE; - *(psElemType* ) & image->type.type = type; - image->parent = NULL; - image->children = NULL; - - return image; -} - -psRegion psRegionSet(float x0, - float x1, - float y0, - float y1) -{ - psRegion out; - - out.x0 = x0; - out.y0 = y0; - out.x1 = x1; - out.y1 = y1; - - return out; -} - -psRegion psRegionFromString(const char* region) -{ - psS32 col0; - psS32 col1; - psS32 row0; - psS32 row1; - - // section should be of the form '[col0:col1,row0:row1]' - if (region == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_SUBSECTION_NULL); - return psRegionSet(NAN,NAN,NAN,NAN); - } - - if (sscanf(region,"[%d:%d,%d:%d]",&col0,&col1,&row0,&row1) < 4) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_SUBSECTION_INVALID, - region); - return psRegionSet(NAN,NAN,NAN,NAN); - } - - if (col0 > col1 || row0 > row1) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_SUBSET_RANGE_MALFORMED, - col0,col1,row0,row1); - return psRegionSet(NAN,NAN,NAN,NAN); - } - - return psRegionSet(col0,col1,row0,row1); -} - -psString psRegionToString(const psRegion region) -{ - char tmpText[256]; // big enough to store any region as text - - snprintf(tmpText,256,"[%g:%g,%g:%g]", - region.x0, region.x1, - region.y0, region.y1); - - return psStringCopy(tmpText); -} - -psImage* psImageRecycle(psImage* old, - int numCols, - int numRows, - const psElemType type) -{ - psS32 elementSize = PSELEMTYPE_SIZEOF(type); // element size in bytes - psS32 rowSize = numCols * elementSize; // row size in bytes. - - if (old == NULL) { - old = psImageAlloc(numCols, numRows, type); - return old; - } - - if (old->type.dimen != PS_DIMEN_IMAGE) { - psFree(old); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_NOT_AN_IMAGE); - return NULL; - } - - /* image already the right size/type? */ - if (numCols == old->numCols && numRows == old->numRows && - type == old->type.type) { - return old; - } - // Resize the image buffer - old->rawDataBuffer = psRealloc(old->data.V[0], - numCols * numRows * elementSize); - old->data.V = (psPtr *)psRealloc(old->data.V, numRows * sizeof(psPtr )); - - // recreate the row pointers - old->data.V[0] = old->rawDataBuffer; - for (psS32 i = 1; i < numRows; i++) { - old->data.V[i] = (psPtr )((int8_t *) old->data.V[i - 1] + rowSize); - } - - *(psU32 *)&old->numCols = numCols; - *(psU32 *)&old->numRows = numRows; - *(psElemType* ) & old->type.type = type; - - return old; -} - -bool p_psImageCopyToRawBuffer(void* buffer, - const psImage* input, - psElemType type) -{ - psElemType inDatatype; - psS32 numRows; - psS32 numCols; - - if (input == NULL || input->data.V == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - return false; - } - - if (input->type.dimen != PS_DIMEN_IMAGE) { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_NOT_AN_IMAGE); - return false; - } - - inDatatype = input->type.type; - numRows = input->numRows; - numCols = input->numCols; - - // cover the trival case of copy of the same - // datatype. - if (type == inDatatype) { - int rowSize = PSELEMTYPE_SIZEOF(inDatatype)*numCols; - for (psS32 row=0;rowdata.V[row], rowSize); - } - return true; - } - - #define PSIMAGE_BUFFER_COPY(INTYPE,OUTTYPE) { \ - ps##INTYPE *in; \ - ps##OUTTYPE *out = buffer; \ - for(psS32 row=0;rowdata.INTYPE[row]; \ - for (psS32 col=0;colchildren != NULL) { - psImage** children = (psImage**)image->children->data; - numFreed = image->children->n; - - // orphan the children first - // (so psFree doesn't try to modify the parent's children array while I'm using it) - for (psS32 i=0;iparent = NULL; - } - - psFree(image->children); - image->children = NULL; - } - - return numFreed; -} - -bool p_psImagePrint (FILE *f, psImage *a, char *name) -{ - - fprintf (f, "matrix: %s\n", name); - - for (int j = 0; j < a[0].numRows; j++) { - for (int i = 0; i < a[0].numCols; i++) { - fprintf (f, "%f ", p_psImageGetElementF64(a, i, j)); - } - fprintf (f, "\n"); - } - fprintf (f, "\n"); - return (true); -} - -complex psImagePixelInterpolate(const psImage* input, - float x, - float y, - const psImage* mask, - psMaskType maskVal, - complex unexposedValue, - psImageInterpolateMode mode) -{ - - if (input == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL,true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - return unexposedValue; - } - - #define PSIMAGE_PIXEL_INTERPOLATE_CASE(TYPE) \ -case PS_TYPE_##TYPE: \ - switch (mode) { \ - case PS_INTERPOLATE_FLAT: \ - return p_psImagePixelInterpolateFLAT_##TYPE( \ - input, \ - x, \ - y, \ - mask, \ - maskVal, \ - unexposedValue); \ - break; \ - case PS_INTERPOLATE_BILINEAR: \ - return p_psImagePixelInterpolateBILINEAR_##TYPE( \ - input, \ - x, \ - y, \ - mask, \ - maskVal, \ - unexposedValue); \ - break; \ - case PS_INTERPOLATE_BILINEAR_VARIANCE: \ - return p_psImagePixelInterpolateBILINEAR_VARIANCE_##TYPE( \ - input, \ - x, \ - y, \ - mask, \ - maskVal, \ - unexposedValue); \ - break; \ - default: \ - psError(PS_ERR_BAD_PARAMETER_VALUE,true, \ - PS_ERRORTEXT_psImage_INTERPOLATE_METHOD_INVALID, \ - mode); \ - } \ - break; - - switch (input->type.type) { - PSIMAGE_PIXEL_INTERPOLATE_CASE(U8); - PSIMAGE_PIXEL_INTERPOLATE_CASE(U16); - PSIMAGE_PIXEL_INTERPOLATE_CASE(U32); - PSIMAGE_PIXEL_INTERPOLATE_CASE(U64); - PSIMAGE_PIXEL_INTERPOLATE_CASE(S8); - PSIMAGE_PIXEL_INTERPOLATE_CASE(S16); - PSIMAGE_PIXEL_INTERPOLATE_CASE(S32); - PSIMAGE_PIXEL_INTERPOLATE_CASE(S64); - PSIMAGE_PIXEL_INTERPOLATE_CASE(F32); - PSIMAGE_PIXEL_INTERPOLATE_CASE(F64); - PSIMAGE_PIXEL_INTERPOLATE_CASE(C32); - PSIMAGE_PIXEL_INTERPOLATE_CASE(C64); - default: { - char* typeStr; - PS_TYPE_NAME(typeStr,input->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE,true, - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, - typeStr); - } - } - - return unexposedValue; -} - -psF64 p_psImageGetElementF64(psImage* image, - int col, - int row) -{ - if (image == NULL) { - return NAN; - } - if (col < 0 || col >= image->numCols) { - return NAN; - } - if (row < 0 || row >= image->numRows) { - return NAN; - } - - switch (image->type.type) { - case PS_TYPE_U8: - return image->data.U8[row][col]; - break; - case PS_TYPE_U16: - return image->data.U16[row][col]; - break; - case PS_TYPE_U32: - return image->data.U32[row][col]; - break; - case PS_TYPE_U64: - return image->data.U64[row][col]; - break; - case PS_TYPE_S8: - return image->data.S8[row][col]; - break; - case PS_TYPE_S16: - return image->data.S16[row][col]; - break; - case PS_TYPE_S32: - return image->data.S32[row][col]; - break; - case PS_TYPE_S64: - return image->data.S64[row][col]; - break; - case PS_TYPE_F32: - return image->data.F32[row][col]; - break; - case PS_TYPE_F64: - return image->data.F64[row][col]; - default: - return NAN; - } -} - -#define PSIMAGE_PIXEL_INTERPOLATE_FLAT(TYPE,RETURNTYPE) \ -inline RETURNTYPE p_psImagePixelInterpolateFLAT_##TYPE( \ - const psImage* input, \ - float x, \ - float y, \ - const psImage* mask, \ - psU32 maskVal, \ - RETURNTYPE unexposedValue) \ -{ \ - psS32 intX = (psS32) round((psF64)(x) - 0.5 + FLT_EPSILON); \ - psS32 intY = (psS32) round((psF64)(y) - 0.5 + FLT_EPSILON); \ - psS32 lastX = input->numCols - 1; \ - psS32 lastY = input->numRows - 1; \ - \ - if ((intX < 0) || \ - (intX > lastX) || \ - (intY < 0) || \ - (intY > lastY) || \ - ( (mask!=NULL) && \ - ((mask->data.PS_TYPE_MASK_DATA[intY][intX] & maskVal) != 0) ) ) { \ - return unexposedValue; \ - } \ - \ - return input->data.TYPE[intY][intX]; \ -} - -PSIMAGE_PIXEL_INTERPOLATE_FLAT(U8,psF64) -PSIMAGE_PIXEL_INTERPOLATE_FLAT(U16,psF64) -PSIMAGE_PIXEL_INTERPOLATE_FLAT(U32,psF64) -PSIMAGE_PIXEL_INTERPOLATE_FLAT(U64,psF64) -PSIMAGE_PIXEL_INTERPOLATE_FLAT(S8,psF64) -PSIMAGE_PIXEL_INTERPOLATE_FLAT(S16,psF64) -PSIMAGE_PIXEL_INTERPOLATE_FLAT(S32,psF64) -PSIMAGE_PIXEL_INTERPOLATE_FLAT(S64,psF64) -PSIMAGE_PIXEL_INTERPOLATE_FLAT(F32,psF64) -PSIMAGE_PIXEL_INTERPOLATE_FLAT(F64,psF64) -PSIMAGE_PIXEL_INTERPOLATE_FLAT(C32,psC64) -PSIMAGE_PIXEL_INTERPOLATE_FLAT(C64,psC64) - -#define PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(TYPE, RETURNTYPE, SUFFIX, FRACFUNC) \ -inline RETURNTYPE p_psImagePixelInterpolateBILINEAR_##SUFFIX( \ - const psImage* input, \ - float x, \ - float y, \ - const psImage* mask, \ - psU32 maskVal, \ - RETURNTYPE unexposedValue) \ -{ \ - double floorX = floor((psF64)(x) - 0.5); \ - double floorY = floor((psF64)(y) - 0.5); \ - psF64 fracX = x - 0.5 - floorX; \ - psF64 fracY = y - 0.5 - floorY; \ - psS32 intFloorX = (psS32) floorX; \ - psS32 intFloorY = (psS32) floorY; \ - psS32 lastX = input->numCols - 1; \ - psS32 lastY = input->numRows - 1; \ - ps##TYPE V00 = 0; \ - ps##TYPE V01 = 0; \ - ps##TYPE V10 = 0; \ - ps##TYPE V11 = 0; \ - psBool valid00 = false; \ - psBool valid01 = false; \ - psBool valid10 = false; \ - psBool 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 */ \ - return V00*FRACFUNC((1.0-fracX)*(1.0-fracY)) + V10*FRACFUNC(fracX*(1.0-fracY)) + \ - V01*FRACFUNC(fracY*(1.0-fracX)) + V11*FRACFUNC(fracX*fracY); \ - } \ - \ - /* OK, at least one pixel is not valid - need to do it piecemeal */ \ - \ - RETURNTYPE V0 = 0.0; \ - psBool valid0 = true; \ - if (valid00 && valid10) { \ - V0 = V00*FRACFUNC(1-fracX)+V10*FRACFUNC(fracX); \ - } else if (valid00) { \ - V0 = V00; \ - } else if (valid10) { \ - V0 = V10; \ - } else { \ - valid0 = false; \ - } \ - \ - RETURNTYPE V1 = 0.0; \ - psBool valid1 = true; \ - if (valid01 && valid11) { \ - V1 = V01*FRACFUNC(1-fracX)+V11*FRACFUNC(fracX); \ - } else if (valid01) { \ - V1 = V01; \ - } else if (valid11) { \ - V1 = V11; \ - } else { \ - valid1 = false; \ - } \ - \ - if (valid0 && valid1) { \ - return V0*FRACFUNC(1-fracY) + V1*FRACFUNC(fracY); \ - } else if (valid0) { \ - return V0; \ - } else if (valid1) { \ - return V1; \ - } \ - \ - return unexposedValue; \ -} - -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U8,psF64,U8,) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U16,psF64,U16,) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U32,psF64,U32,) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U64,psF64,U64,) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S8,psF64,S8,) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S16,psF64,S16,) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S32,psF64,S32,) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S64,psF64,S64,) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(F32,psF64,F32,) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(F64,psF64,F64,) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(C32,psC64,C32,) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(C64,psC64,C64,) - -// Variance Version -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U8,psF64,VARIANCE_U8,SQUARE) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U16,psF64,VARIANCE_U16,SQUARE) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U32,psF64,VARIANCE_U32,SQUARE) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(U64,psF64,VARIANCE_U64,SQUARE) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S8,psF64,VARIANCE_S8,SQUARE) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S16,psF64,VARIANCE_S16,SQUARE) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S32,psF64,VARIANCE_S32,SQUARE) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(S64,psF64,VARIANCE_S64,SQUARE) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(F32,psF64,VARIANCE_F32,SQUARE) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(F64,psF64,VARIANCE_F64,SQUARE) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(C32,psC64,VARIANCE_C32,SQUARE) -PSIMAGE_PIXEL_INTERPOLATE_BILINEAR(C64,psC64,VARIANCE_C64,SQUARE) Index: trunk/psLib/src/image/psImage.h =================================================================== --- trunk/psLib/src/image/psImage.h (revision 4534) +++ (revision ) @@ -1,244 +1,0 @@ -/** @file psImage.h - * - * @brief Contains basic image definitions and operations - * - * This file defines the basic type for an image struct and functions useful - * in manupulating images. - * - * @ingroup Image - * - * @author Robert DeSonia, MHPCC - * @author Ross Harman, MHPCC - * - * @version $Revision: 1.62 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-07-11 21:38:19 $ - * - * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii - */ -#ifndef PS_IMAGE_H -#define PS_IMAGE_H - -#include -#include -#include "psType.h" -#include "psArray.h" - -/// @addtogroup Image -/// @{ - -/** enumeration of options in interpolation - * - */ -typedef enum { - PS_INTERPOLATE_FLAT, ///< 'flat' interpolation (nearest pixel) - PS_INTERPOLATE_BILINEAR, ///< bi-linear interpolation - PS_INTERPOLATE_LANCZOS2, ///< Sinc interpolation with 4x4 pixel kernel - PS_INTERPOLATE_LANCZOS3, ///< Sinc interpolation with 6x6 pixel kernel - PS_INTERPOLATE_LANCZOS4, ///< Sinc interpolation with 8x8 pixel kernel - PS_INTERPOLATE_BILINEAR_VARIANCE, ///< Variance version of PS_INTERPOLATE_BILINEAR - PS_INTERPOLATE_LANCZOS2_VARIANCE, ///< Variance version of PS_INTERPOLATE_LANCZOS2 - PS_INTERPOLATE_LANCZOS3_VARIANCE, ///< Variance version of PS_INTERPOLATE_LANCZOS3 - PS_INTERPOLATE_LANCZOS4_VARIANCE, ///< Variance version of PS_INTERPOLATE_LANCZOS4 - // PS_INTERPOLATE_NUM_MODES ///< enum end-marker; does not coorespond to a interpolation mode -} psImageInterpolateMode; - -/** Basic image data structure. - * - * Struct for maintaining image data of varying types. It also contains - * information about image size, parent images and children images. - * - */ -typedef struct psImage -{ - const psMathType type; ///< Image data type and dimension. - const int numCols; ///< Number of columns in image - const int numRows; ///< Number of rows in image. - const int col0; ///< Column position relative to parent. - const int row0; ///< Row position relative to parent. - - union { - psU8** U8; ///< Unsigned 8-bit integer data. - psU16** U16; ///< Unsigned 16-bit integer data. - psU32** U32; ///< Unsigned 32-bit integer data. - psU64** U64; ///< Unsigned 64-bit integer data. - psS8** S8; ///< Signed 8-bit integer data. - psS16** S16; ///< Signed 16-bit integer data. - psS32** S32; ///< Signed 32-bit integer data. - psS64** S64; ///< Signed 64-bit integer data. - psF32** F32; ///< Single-precision float data. - psF64** F64; ///< Double-precision float data. - psC32** C32; ///< Single-precision complex data. - psC64** C64; ///< Double-precision complex data. - psPtr** PTR; ///< Void pointers. - psPtr* V; ///< Pointer to data. - } data; ///< Union for data types. - const struct psImage* parent; ///< Parent, if a subimage. - psArray* children; ///< Children of this region. - - psPtr rawDataBuffer; ///< Raw data buffer for Allocating/Freeing Images - void *lock; ///< Optional lock for thread safety -} -psImage; - -/** Basic image region structure. - * - * Struct for specifying a rectangular area in an image. - * - */ -typedef struct -{ - float x0; ///< the first column of the region. - float x1; ///< the last column of the region. - float y0; ///< the first row of the region. - float y1; ///< the last row of the region. -} -psRegion; - -/** Create an image of the specified size and type. - * - * Uses psLib memory allocation functions to create an image struct of the - * specified size and type. - * - * @return psImage* : Pointer to psImage. - * - */ -psImage* psImageAlloc( - int numCols, ///< Number of rows in image. - int numRows, ///< Number of columns in image. - psElemType type ///< Type of data for image. -) -; - -/** Create a psRegion with the specified attributes. - * - * @return psRegion : a cooresponding psRegion. - */ -psRegion psRegionSet( - float x0, ///< the first column of the region. - float x1, ///< the last column of the region + 1. - float y0, ///< the first row of the region. - float y1 ///< the last row of the region + 1. -); - -/** Create a psRegion with the attribute values given as a string. - * - * Create a psRegion with the attribute values given as a string. The format - * shall be of the standard IRAF form '[x0:x1,y0:y1]' - * - * @return psRegion: A new psRegion struct, or NULL is not successful. - */ -psRegion psRegionFromString( - const char* region ///< image rectangular region in the form '[x0:x1,y0:y1]' -); - -/** Create a string of the standard IRAF form '[x0:x1,y0:y1]' from a psRegion. - * - * @return psString: A new string representing the psRegion as text, or NULL - * is not successful. - */ -psString psRegionToString( - const psRegion region ///< the psRegion to convert to a string -); - -/** Resize a given image to the given size/type. - * - * @return psImage* Resized psImage. - * - */ -psImage* psImageRecycle( - psImage* old, ///< the psImage to recycle by resizing image buffer - int numCols, ///< the desired number of columns in image - int numRows, ///< the desired number of rows in image - const psElemType type ///< the desired datatype of the image -); - -/** Copy an image to a new buffer - * - * @return True if image copied or false if error - */ -bool p_psImageCopyToRawBuffer( - void* buffer, ///< the buffer used to copy the image - const psImage* input, ///< the input image to be copied - psElemType type ///< the datatype of the image to be copied -); - -/** Frees all children of a psImage. - * - * @return int Number of children freed. - * - */ -int psImageFreeChildren( - psImage* image ///< psImage in which all children shall be deallocated -); - -/** get an element of an image as a psF64. - * - * @return psF64 pixel value at specified location - */ -psF64 p_psImageGetElementF64( - psImage* image, ///< input image - int col, ///< pixel column - int row ///< pixel row -); - -/** print image pixel values. - * - * @return bool TRUE is successful, otherwise FALSE. - */ -bool p_psImagePrint( - FILE *f, ///< Destination stream - psImage *a, ///< image to print - char *name ///< name of the image (for title) -); - -/** Interpolate image pixel value given floating point coordinates. - * - * @return psF32 Pixel value interpolated from image or unexposedValue if - * given x,y doesn't coorespond to a valid image location - */ -complex psImagePixelInterpolate( - 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 - psMaskType maskVal, ///< the mask value - complex unexposedValue, ///< return value if x,y location is not in image. - psImageInterpolateMode mode ///< interpolation mode -); - -#define PIXEL_INTERPOLATE_FCN_PROTOTYPE(SUFFIX, RETURNTYPE) \ -inline RETURNTYPE p_psImagePixelInterpolate##SUFFIX( \ - 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 */ \ - psU32 maskVal, /**< the mask value */ \ - RETURNTYPE unexposedValue /**< return value if x,y location is not in image. */ \ - ); - -#define PIXEL_INTERPOLATE_FCNS(MODE) \ -PIXEL_INTERPOLATE_FCN_PROTOTYPE(MODE##_U8,psF64) \ -PIXEL_INTERPOLATE_FCN_PROTOTYPE(MODE##_U16,psF64) \ -PIXEL_INTERPOLATE_FCN_PROTOTYPE(MODE##_U32,psF64) \ -PIXEL_INTERPOLATE_FCN_PROTOTYPE(MODE##_U64,psF64) \ -PIXEL_INTERPOLATE_FCN_PROTOTYPE(MODE##_S8,psF64) \ -PIXEL_INTERPOLATE_FCN_PROTOTYPE(MODE##_S16,psF64) \ -PIXEL_INTERPOLATE_FCN_PROTOTYPE(MODE##_S32,psF64) \ -PIXEL_INTERPOLATE_FCN_PROTOTYPE(MODE##_S64,psF64) \ -PIXEL_INTERPOLATE_FCN_PROTOTYPE(MODE##_F32,psF64) \ -PIXEL_INTERPOLATE_FCN_PROTOTYPE(MODE##_F64,psF64) \ -PIXEL_INTERPOLATE_FCN_PROTOTYPE(MODE##_C32,psC64) \ -PIXEL_INTERPOLATE_FCN_PROTOTYPE(MODE##_C64,psC64) - -#ifndef SWIG -PIXEL_INTERPOLATE_FCNS(FLAT) -PIXEL_INTERPOLATE_FCNS(BILINEAR) -PIXEL_INTERPOLATE_FCNS(BILINEAR_VARIANCE) -#endif // ! SWIG - -#undef PIXEL_INTERPOLATE_FCN_PROTOTYPE -#undef PIXEL_INTERPOLATE_FCNS - -/// @} - -#endif // PS_IMAGE_H Index: trunk/psLib/src/image/psImageConvolve.c =================================================================== --- trunk/psLib/src/image/psImageConvolve.c (revision 4534) +++ (revision ) @@ -1,477 +1,0 @@ -/* @file psImageConvolve.c - * - * @brief Contains FFT transform related functions for psImage. - * - * @author Robert DeSonia, MHPCC - * - * @version $Revision: 1.21 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-07-11 21:38:19 $ - * - * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii - */ - -#include - -#include "psImageConvolve.h" -#include "psImageFFT.h" -#include "psImageStructManip.h" -#include "psBinaryOp.h" -#include "psMemory.h" -#include "psLogMsg.h" -#include "psError.h" - -#include "psImageErrors.h" - -#define FOURIER_PADDING 32 /* padding amount in every side of the image for fourier convolution */ - -static void freeKernel(psKernel* ptr); - -psKernel* psKernelAlloc(int xMin, int xMax, int yMin, int yMax) -{ - psKernel* result; - psS32 numRows; - psS32 numCols; - - // following is explicitly spelled out in the SDRS as a requirement - if (yMin > yMax) { - psLogMsg(__func__, PS_LOG_WARN, - "Specified yMin, %d, was greater than yMax, %d. Values swapped.", - yMin, yMax); - - psS32 temp = yMin; - yMin = yMax; - yMax = temp; - } - - // following is explicitly spelled out in the SDRS as a requirement - if (xMin > xMax) { - psLogMsg(__func__, PS_LOG_WARN, - "Specified xMin, %d, was greater than xMax, %d. Values swapped.", - xMin, xMax); - - psS32 temp = xMin; - xMin = xMax; - xMax = temp; - } - - numRows = yMax - yMin + 1; - numCols = xMax - xMin + 1; - - result = psAlloc(sizeof(psKernel)); - result->xMin = xMin; - result->xMax = xMax; - result->yMin = yMin; - result->yMax = yMax; - result->image = psImageAlloc(numCols,numRows,PS_TYPE_KERNEL); - memset(result->image->rawDataBuffer,0,numCols*numRows*PSELEMTYPE_SIZEOF(PS_TYPE_KERNEL)); - result->p_kernelRows = psAlloc(sizeof(float*)*numRows); - - float** kernelRows = result->p_kernelRows; - float** imageRows = result->image->data.PS_TYPE_KERNEL_DATA; - for (psS32 i = 0; i < numRows; i++) { - kernelRows[i] = imageRows[i] - xMin; - } - result->kernel = kernelRows - yMin; - - psMemSetDeallocator(result,(psFreeFunc)freeKernel); - - return result; -} - -void freeKernel(psKernel* ptr) -{ - if (ptr != NULL) { - psFree(ptr->image); - psFree(ptr->p_kernelRows); - } -} - -psKernel* psKernelGenerate(const psVector* tShifts, - const psVector* xShifts, - const psVector* yShifts, - bool relative) -{ - psS32 lastX; - psS32 lastY; - psS32 lastT; - psS32 x; - psS32 y; - psS32 t; - psS32 xMin = 0; - psS32 xMax = 0; - psS32 yMin = 0; - psS32 yMax = 0; - psS32 length = 0; - float normalizeTime = 1.0; // fraction of total time for each shift clock - psKernel* result = NULL; - float** kernel = NULL; - - // got non-NULL vectors? - if (tShifts == NULL || xShifts == NULL || yShifts == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImageConvolve_SHIFT_NULL); - return NULL; - } - - // types match? - if (xShifts->type.type != yShifts->type.type || - tShifts->type.type != xShifts->type.type) { - char* typeXStr; - char* typeYStr; - char* typeTStr; - PS_TYPE_NAME(typeXStr,xShifts->type.type); - PS_TYPE_NAME(typeYStr,yShifts->type.type); - PS_TYPE_NAME(typeTStr,tShifts->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImageConvolve_SHIFT_TYPE_MISMATCH, - typeTStr, typeXStr, typeYStr); - return NULL; - } - - // sizes match? - length = xShifts->n; - if (length != yShifts->n || - length != tShifts->n) { - psError(PS_ERR_BAD_PARAMETER_SIZE, true, - "Shift vectors can not be of different sizes."); - return NULL; - } - - // if no shifts, the kernel is just a 1 at 0,0 - if (length < 1) { - result = psKernelAlloc(0,0,0,0); - result->kernel[0][0] = 1; - return result; - } - - #define KERNEL_GENERATE_CASE(TYPE) \ -case PS_TYPE_##TYPE: { \ - ps##TYPE *tShiftData = tShifts->data.TYPE; \ - ps##TYPE *xShiftData = xShifts->data.TYPE; \ - ps##TYPE *yShiftData = yShifts->data.TYPE; \ - lastX = xShiftData[length-1]; \ - lastY = yShiftData[length-1]; \ - lastT = tShiftData[length-1]; \ - \ - for (int lcv = 0; lcv < length; lcv++) { \ - x = lastX - xShiftData[lcv]; \ - y = lastY - yShiftData[lcv]; \ - \ - if (x < xMin) { \ - xMin = x; \ - } else if (x > xMax) { \ - xMax = x; \ - } \ - if (y < yMin) { \ - yMin = y; \ - } else if (y > yMax) { \ - yMax = y; \ - } \ - } \ - \ - normalizeTime = 1.0 / (float)(tShiftData[length-1]); \ - result = psKernelAlloc(xMin,xMax,yMin,yMax); \ - kernel = result->kernel; \ - \ - psS32 prevT = 0; \ - for (int i = 0; i < length; i++) { \ - t = tShiftData[i] - prevT; \ - x = lastX - xShiftData[i]; \ - y = lastY - yShiftData[i]; \ - \ - kernel[y][x] += (float)t / (float)lastT; \ - prevT = tShiftData[i]; \ - } \ - break; \ - } - - #define RELATIVE_KERNEL_GENERATE_CASE(TYPE) \ -case PS_TYPE_##TYPE: { \ - ps##TYPE *tShiftData = tShifts->data.TYPE; \ - ps##TYPE *xShiftData = xShifts->data.TYPE; \ - ps##TYPE *yShiftData = yShifts->data.TYPE; \ - \ - x = 0; \ - y = 0; \ - t = 0; \ - \ - for (int lcv = length-1; lcv >= 0; lcv--) { \ - t += tShiftData[lcv]; \ - \ - if (x < xMin) { \ - xMin = x; \ - } else if (x > xMax) { \ - xMax = x; \ - } \ - if (y < yMin) { \ - yMin = y; \ - } else if (y > yMax) { \ - yMax = y; \ - } \ - x -= xShiftData[lcv]; \ - y -= yShiftData[lcv]; \ - \ - } \ - result = psKernelAlloc(xMin,xMax,yMin,yMax); \ - kernel = result->kernel; \ - \ - normalizeTime = 1.0 / (float)t; \ - x = 0; \ - y = 0; \ - for (psS32 i = length-1; i >= 0; i--) { \ - kernel[y][x] += (float)(tShiftData[i]) * normalizeTime; \ - x -= xShiftData[i]; \ - y -= yShiftData[i]; \ - \ - } \ - break; \ - } - - if (relative) { - switch (xShifts->type.type) { - RELATIVE_KERNEL_GENERATE_CASE(U8); - RELATIVE_KERNEL_GENERATE_CASE(U16); - RELATIVE_KERNEL_GENERATE_CASE(U32); - RELATIVE_KERNEL_GENERATE_CASE(U64); - RELATIVE_KERNEL_GENERATE_CASE(S8); - RELATIVE_KERNEL_GENERATE_CASE(S16); - RELATIVE_KERNEL_GENERATE_CASE(S32); - RELATIVE_KERNEL_GENERATE_CASE(S64); - RELATIVE_KERNEL_GENERATE_CASE(F32); - RELATIVE_KERNEL_GENERATE_CASE(F64); - RELATIVE_KERNEL_GENERATE_CASE(C32); - RELATIVE_KERNEL_GENERATE_CASE(C64); - - default: { - char* typeStr; - PS_TYPE_NAME(typeStr,xShifts->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, - typeStr); - } - } - } else { - switch (xShifts->type.type) { - KERNEL_GENERATE_CASE(U8); - KERNEL_GENERATE_CASE(U16); - KERNEL_GENERATE_CASE(U32); - KERNEL_GENERATE_CASE(U64); - KERNEL_GENERATE_CASE(S8); - KERNEL_GENERATE_CASE(S16); - KERNEL_GENERATE_CASE(S32); - KERNEL_GENERATE_CASE(S64); - KERNEL_GENERATE_CASE(F32); - KERNEL_GENERATE_CASE(F64); - KERNEL_GENERATE_CASE(C32); - KERNEL_GENERATE_CASE(C64); - - default: { - char* typeStr; - PS_TYPE_NAME(typeStr,xShifts->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, - typeStr); - } - } - } - - return result; -} - -psImage* psImageConvolve(psImage* out, const psImage* in, const psKernel* kernel, bool direct) -{ - if (in == NULL) { - psFree(out); - return NULL; - } - - if (kernel == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImageConvolve_KERNEL_NULL); - psFree(out); - return NULL; - } - psS32 xMin = kernel->xMin; - psS32 xMax = kernel->xMax; - psS32 yMin = kernel->yMin; - psS32 yMax = kernel->yMax; - float** kData = kernel->kernel; - - // make the output image to the proper size and type - psS32 numRows = in->numRows; - psS32 numCols = in->numCols; - - - - if (direct) { - // spatial convolution - - #define SPATIAL_CONVOLVE_CASE(TYPE) \ - case PS_TYPE_##TYPE: { \ - ps##TYPE** inData = in->data.TYPE; \ - out = psImageRecycle(out, numCols, numRows, PS_TYPE_##TYPE); \ - for (psS32 row=0;rowdata.TYPE[row]; \ - for (psS32 col=0;col= 0 && row-kRow < numRows) { \ - for (psS32 kCol = xMin; kCol < xMax; kCol++) { \ - if (col-kCol >= 0 && col-kCol < numCols) { \ - pixel += kData[kRow][kCol] * inData[row-kRow][col-kCol]; \ - } \ - } \ - } \ - } \ - outRow[col] = pixel; \ - } \ - } \ - } \ - break; - - switch (in->type.type) { - SPATIAL_CONVOLVE_CASE(U8) - SPATIAL_CONVOLVE_CASE(U16) - SPATIAL_CONVOLVE_CASE(U32) - SPATIAL_CONVOLVE_CASE(U64) - SPATIAL_CONVOLVE_CASE(S8) - SPATIAL_CONVOLVE_CASE(S16) - SPATIAL_CONVOLVE_CASE(S32) - SPATIAL_CONVOLVE_CASE(S64) - SPATIAL_CONVOLVE_CASE(F32) - SPATIAL_CONVOLVE_CASE(F64) - SPATIAL_CONVOLVE_CASE(C32) - SPATIAL_CONVOLVE_CASE(C64) - - default: { - char* typeStr; - PS_TYPE_NAME(typeStr,in->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, - typeStr); - psFree(out); - return NULL; - - } - } - - - } else { - // fourier convolution - psS32 paddedCols = numCols+2*FOURIER_PADDING; - psS32 paddedRows = numRows+2*FOURIER_PADDING; - - // check to see if kernel is smaller, otherwise padding it up will fail. - psS32 kRows = kernel->image->numRows; - psS32 kCols = kernel->image->numCols; - if (kRows >= numRows || kCols >= numCols) { - psError(PS_ERR_BAD_PARAMETER_SIZE, true, - PS_ERRORTEXT_psImageConvolve_KERNEL_TOO_LARGE, - kCols,kRows, - numCols, numRows); - psFree(out); - return NULL; - } - - // pad the image - psImage* paddedImage = psImageAlloc(paddedCols,paddedRows,in->type.type); - psS32 elementSize = PSELEMTYPE_SIZEOF(in->type.type); - - // zero out padded area on top and bottom - memset(paddedImage->data.U8[0],0,FOURIER_PADDING*paddedCols*elementSize); - memset(paddedImage->data.U8[FOURIER_PADDING+numRows-1],0,FOURIER_PADDING*paddedCols*elementSize); - - // fill in the image-containing rows. - psS32 sidePaddingSize = FOURIER_PADDING*elementSize; - psS32 imageRowSize = numCols*elementSize; - psU8* paddedData = paddedImage->data.U8[FOURIER_PADDING]; - for (psS32 row=0;rowdata.U8[row],imageRowSize); - paddedData += imageRowSize; - // zero out padded area on right edge. - memset(paddedData,0,sidePaddingSize); - paddedData += sidePaddingSize; - } - - // pad the kernel to the same size of paddedImage - psImage* paddedKernel = psImageAlloc(paddedCols,paddedRows,PS_TYPE_KERNEL); - memset(paddedKernel->data.U8[0],0,sizeof(float)*numCols*numRows); // zero-out image - psS32 yMax = kernel->yMax; - psS32 xMax = kernel->xMax; - for (psS32 row = kernel->yMin; row <= yMax;row++) { - psS32 padRow = row; - if (padRow < 0) { - padRow += paddedRows; - } - float* padData = paddedKernel->data.PS_TYPE_KERNEL_DATA[padRow]; - float* kernelRow = kernel->kernel[row]; - for (psS32 col = kernel->xMin; col <= xMax; col++) { - if (col < 0) { - padData[col+paddedCols] = kernelRow[col]; - } else { - padData[col] = kernelRow[col]; - } - } - } - - psImage* kernelFourier = psImageFFT(NULL, paddedKernel, PS_FFT_FORWARD); - if (kernelFourier == NULL) { - psError(PS_ERR_UNKNOWN, false, - PS_ERRORTEXT_psImageConvolve_KERNEL_FFT_FAILED); - psFree(out); - return NULL; - } - - psImage* inFourier = psImageFFT(NULL, paddedImage, PS_FFT_FORWARD); - if (inFourier == NULL) { - psError(PS_ERR_UNKNOWN, false, - PS_ERRORTEXT_psImageConvolve_FFT_FAILED); - psFree(out); - return NULL; - } - - // convolution in fourier domain is just a pixel-wise multiplication - for (int row = 0; row < paddedRows; row++) { - psC32* inRow = inFourier->data.C32[row]; - psC32* kRow = kernelFourier->data.C32[row]; - for (int col = 0; col < paddedCols; col++) { - inRow[col] *= kRow[col]; - } - } - - psImage* complexOut = psImageFFT(NULL, inFourier, - PS_FFT_REVERSE); - if (complexOut == NULL) { - psError(PS_ERR_UNKNOWN, false, - PS_ERRORTEXT_psImageConvolve_FFT_FAILED); - psFree(out); - return NULL; - } - - // subset out the padded area now. - psImage* complexOutSansPad = psImageSubset(complexOut, - psRegionSet(FOURIER_PADDING, FOURIER_PADDING+numCols,FOURIER_PADDING,FOURIER_PADDING+numRows)); - - out = psImageRecycle(out,numCols,numRows,PS_TYPE_F32); - float factor = 1.0f/(float)paddedCols/(float)paddedRows; - for (psS32 row = 0; row < numRows; row++) { - psF32* outRow = out->data.F32[row]; - psC32* resultRow = complexOutSansPad->data.C32[row]; - for (psS32 col = 0; col < numCols; col++) { - outRow[col] = crealf(resultRow[col])*factor; - } - } - - psFree(complexOut); // frees complexOutSansPad, as it is a child. - psFree(kernelFourier); - psFree(inFourier); - psFree(paddedImage); - psFree(paddedKernel); - - } - - return out; -} Index: trunk/psLib/src/image/psImageConvolve.h =================================================================== --- trunk/psLib/src/image/psImageConvolve.h (revision 4534) +++ (revision ) @@ -1,128 +1,0 @@ -/** @file psImageConvolve.h - * - * @brief image convolution functionality - * - * @ingroup Transform - * - * @author Robert DeSonia, MHPCC - * - * @version $Revision: 1.9 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-07-11 21:38:19 $ - * - * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii - */ - -#ifndef PS_IMAGE_CONVOLVE_H -#define PS_IMAGE_CONVOLVE_H - -#include "psImage.h" -#include "psVector.h" -#include "psType.h" - -#define PS_TYPE_KERNEL PS_TYPE_F32 /**< the data member to use for kernel image */ -#define PS_TYPE_KERNEL_DATA F32 /**< the data member to use for kernel image */ -#define PS_TYPE_KERNEL_NAME "psF32" /**< the data type for kernel as a string */ - -/** Kernel Type - * - * A floating-point data type used for storing kernel data. - * - */ -//typedef float psKernelType; - -/** A convolution kernel */ -typedef struct -{ - psImage* image; ///< Kernel data, in the form of an image - int xMin; ///< Most negative x index - int yMin; ///< Most negative y index - int xMax; ///< Most positive x index - int yMax; ///< Most positive y index - float** kernel; ///< Pointer to the kernel data - float** p_kernelRows; ///< Pointer to the rows of the kernel data; not intended for user use. -} -psKernel; - -/** Allocates a convolution kernel of the given range - * - * In order to perform a convolution, we need to define the convolution - * kernel. We need a more general object than a psImage so that we can - * incorporate the offset from the (0, 0) pixel to the (0, 0) value of the - * kernel. It might be convenient to allow both positive and negative - * indices to convey the positive and negative shifts. One might consider - * setting the x0 and y0 members of a psImage to the appropriate offsets, - * but this is not the purpose of these members, and doing so may affect the - * behavior of other psImage operations. - * - * This construction allows the kernel member to use negative indices, while - * preserving the location of psMemBlocks relative to allocated memory. - * - * The maximum extent of the kernel shifts shall be defined by the xMin, - * xMax, yMin and yMax members. Note that xMin and yMin, under normal - * circumstances, should be negative numbers. That is, - * myKernel->kernel[-3][-2] may be defined if yMin and xMin are equal to or - * more negative than -3 and -2, respectively. - * - * In the event that one of the minimum values is greater than the - * corresponding maximum value, the function shall generate a warning, and - * the offending values shall be exchanged. - * - * @return psKernel* A new kernel object - */ -psKernel* psKernelAlloc( - int xMin, ///< Most negative x index - int xMax, ///< Most positive x index - int yMin, ///< Most negative y index - int yMax ///< Most positive y index -); - -/** Generates a kernel given a list of shift values - * - * Given a list of values (e.g., shifts made in the course of OT guiding), - * psKernelGenerate shall return the appropriate kernel. The vectors xShifts - * and yShifts, which are a list of shifts relative to some starting point, - * will be supplied by the user. The elements of the vectors should be of an - * integer type; otherwise the values shall be truncated to integers. The - * output kernel shall be normalized such that the sum over the kernel is - * unity. - * - * If the vectors are not of the same number of elements, then the function - * shall generate a warning shall be generated, following which, the longer - * vector trimmed to the length of the shorter, and the function shall continue. - * - * @return psKernel* new Kernel object - */ -psKernel* psKernelGenerate( - const psVector* tShifts, ///< list of time shifts - const psVector* xShifts, ///< list of x-axis shifts - const psVector* yShifts, ///< list of y-axis shifts - bool relative - /**< specifies the starting point for the shifts; true=relative to previous shift - * false = relative to some other starting point. */ -); - -/** convolve an image with a kernel - * - * Given an input image and the convolution kernel, psImageConvolve shall - * convolve the input image, in, with the kernel, kernel and return the - * convolved image, out. - * - * Two methods shall be available for the convolution: if direct is true, - * then the convolution shall be performed in real space (appropriate for - * small kernels); otherwise, the convolution shall be performed using Fast - * Fourier Transforms (FFTs; appropriate for larger kernels). The latter - * option involves padding the input image, copying the kernel into an image - * of the same size as the padded input image, performing an FFT on each, - * multiplying the FFTs, and performing an inverse FFT before trimming the - * image back to the original size. - * - * @return psImage* resulting image - */ -psImage* psImageConvolve( - psImage* out, ///< a psImage to recycle. If NULL, a new psImage is made. - const psImage* in, ///< the psImage to convolve - const psKernel* kernel, ///< kernel to colvolve with - bool direct ///< specifies method, true=direct convolution, false=fourier -); - -#endif // #ifndef PS_IMAGE_CONVOLVE_H Index: trunk/psLib/src/image/psImageErrors.dat =================================================================== --- trunk/psLib/src/image/psImageErrors.dat (revision 4534) +++ (revision ) @@ -1,81 +1,0 @@ -# -# This file is used to generate psImageErrors.h content -# -# Format is: -# ERRORNAME(one word) ERRORTEXT -# -# N.B. in code, the ERRORNAME appears as PS_ERRORTEXT_ERRORNAME -#################################################################### -# psImage -psImage_AREA_NEGATIVE Specified number of rows (%d) or columns (%d) is invalid. -psImage_NOT_AN_IMAGE The input psImage must have a PS_DIMEN_IMAGE dimension type. -psImage_IMAGE_NULL Can not operate on a NULL psImage. -psImage_IMAGE_TYPE_UNSUPPORTED Specified psImage type, %s, is not supported. -psImage_INTERPOLATE_METHOD_INVALID Specified interpolation method (%d) is not supported. -psImage_REGION_NULL Specified psRegion is NULL. Operation could not be performed. -psImage_SUBSET_RANGE_INVALID Specified subset range, [%d:%d,%d:%d], is invalid or outside input psImage's boundaries, [0:%d,0:%d]. -psImage_SUBSET_RANGE_MALFORMED Specified subset range, [%d:%d,%d:%d], is invalid. Ranges must be incremental. -psImage_SUBSECTION_NULL Specified subsection string can not be NULL. -psImage_SUBSECTION_INVALID Specified subsection string, '%s', can not be parsed. Must be in the form '[x1:x2,y1:y2]'. -psImage_NOT_PARENT Specified psImage can not be a child of another psImage. -psImage_INPLACE_NOTSUPPORTED Specified input and output psImage can not reference the same psImage. -psImage_SUBSET_ZERO_SIZE Specified subset, [%d:%d,%d:%d], contains no pixel data. -psImage_IMAGE_MASK_SIZE Input psImage mask size, %dx%d, does not match psImage input size, %dx%d. -psImage_IMAGE_MASK_TYPE Input psImage mask type, %s, is not the supported mask datatype of %s. -psImage_BAD_STAT Specified statistic option, %d, is not valid. Must specify one and only one statistic type. -psImage_NO_STAT_OPTIONS Specified statistic option did not indicate any operation to perform. -psImage_STAT_NULL Specified statistic can not be NULL. -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, (%g,%g), 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. -psImageFFT_REVERSE_NOT_COMPLEX Input psImage (%s) is not complex. Reverse FFT operation requires a complex psImage input. -psImageFFT_FORWARD_NOT_REAL Input psImage (%s) is not real. Forward FFT operation requires a real psImage input. -psImageFFT_FFTW_PLAN_NULL Could not create a valid FFT plan to perform the transform. -psImageFFT_REAL_IMAG_TYPE_MISMATCH Real psImage type (%s) and imaginary psImage type (%s) must be the same. -psImageFFT_REAL_IMAG_SIZE_MISMATCH Real psImage size (%dx%d) and imaginary psImage size (%dx%d) must be the same. -psImageFFT_NONREAL_NOTSUPPORTED Input psImage type, %s, is required to be either psF32 or psF64. -psImageFFT_NONCOMPLEX_NOTSUPPORTED Input psImage type, %s, is required to be either psC32 or psC64. -psImageFFT_REAL_FORWARD_NOTSUPPORTED The PS_FFT_FORWARD and PS_FFT_REAL_RESULT combinition is not supported. -psImageFFT_FORWARD_REVERSE Can not specify both PS_FFT_FORWARD and PS_FFT_REVERSE options. -psImageFFT_NO_DIRECTION_OPTION Must specify either PS_FFT_FORWARD or PS_FFT_REVERSE option. -# -psImageIO_FILENAME_NULL Specified filename can not be NULL. -psImageIO_FILENAME_INVALID Could not open file,'%s'.\nCFITSIO Error: %s -psImageIO_EXTNAME_INVALID Could not find HDU with extension name '%s' in file %s.\nCFITSIO Error: %s -psImageIO_EXTNUM_INVALID Could not find HDU #%d in file %s.\nCFITSIO Error: %s -psImageIO_DATATYPE_UNKNOWN Could not determine image data type for file %s.\nCFITSIO Error: %s -psImageIO_IMAGE_DIM_UNKNOWN Could not determine image dimensions for file %s.\nCFITSIO Error: %s -psImageIO_IMAGE_SIZE_UNKNOWN Could not determine image size for file %s.\nCFITSIO Error: %s -psImageIO_IMAGE_DIMENSION_UNSUPPORTED Image number of dimensions, %d, is not valid. Only two or three dimensions supported for FITS I/O. -psImageIO_FITS_TYPE_UNSUPPORTED FITS image type, BITPIX=%d, in file %s is not supported. -psImageIO_READ_FAILED Reading from FITS file %s failed.\nCFITSIO Error: %s -psImageIO_TYPE_UNSUPPORTED Input psImage type, %s, is not supported. -psImageIO_WRITE_EXTNUM_INVALID Specified extension number, %d, must not exceed number of HDUs, %d, by more than one. -psImageIO_FILENAME_CREATE_FAILED Could not create file,'%s'.\nCFITSIO Error: %s -psImageIO_CREATE_EXTENSION_FAILED Could not create EXTNAME keyword for file,'%s'.\nCFITSIO Error: %s -psImageIO_CREATE_HDU_FAILED Could not create HDU for writing psImage in file,'%s'.\nCFITSIO Error: %s -psImageIO_WRITE_FAILED Could not write psImage data to file,'%s'.\nCFITSIO Error: %s -# -psImageManip_MAXMIN Specified min value, %g, can not be greater than the specified max value, %g. -psImageManip_MAXMIN_REAL Specified real-portion of min value, %g, can not be greater than the real-portion of max value, %g. -psImageManip_MAXMIN_IMAG Specified imaginary-portion of min value, %g, can not be greater than the imaginary-portion of max value, %g. -psImageManip_OPERATION_NULL Operation can not be NULL. -psImageManip_OVERLAY_TYPE_MISMATCH Input overlay psImage type, %s, must match input psImage type, %s. -psImageManip_CLIP_VALUE_INVALID Specified %s value, %g%+gi, is not the the range of input psImage's valid pixel values (%s), i.e. [%g:%g]. -psImageManip_OVERLAY_OPERATOR_INVALID Specified operation, '%s', is not supported. -psImageManip_SCALE_NOT_POSITIVE Specified scale value, %d, must be a positive value. -psImageManip_INTERPOLATION_MODE_UNSUPPORTED Specified interpolation mode, %d, is unsupported. -psImageConvolve_SHIFT_NULL Specified shift vectors can not be NULL. -psImageConvolve_KERNEL_NULL Specified psKernel can not be NULL. -psImageConvolve_SHIFT_TYPE_MISMATCH Input t-, x-, and y-shift vector types (%s/%s/%s) must match. -psImageConvolve_FFT_FAILED Failed to perform a fourier transform of input image. -psImageConvolve_KERNEL_FFT_FAILED Failed to perform a fourier transform of kernel. -psImageConvolve_KERNEL_TOO_LARGE Specified psKernel size, %dx%d, can not be larger than input psImage size, %dx%d. -psImage_COEFF_NULL Polynomial coefficients cannot be NULL. -psImageManip_TRANSFORM_NULL Specified input transform can not be NULL. Index: trunk/psLib/src/image/psImageErrors.h =================================================================== --- trunk/psLib/src/image/psImageErrors.h (revision 4534) +++ (revision ) @@ -1,101 +1,0 @@ -/** @file psImageErrors.h - * - * @brief Contains the error text for the image functions - * - * @ingroup ErrorHandling - * - * @author Robert DeSonia, MHPCC - * - * @version $Revision: 1.16 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-06-08 23:40:45 $ - * - * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii - */ - -#ifndef PS_IMAGE_ERRORS_H -#define PS_IMAGE_ERRORS_H - -/* N.B., lines between '//~Start' and '//~End' are automatic generated from - * the template following the '//~Start'. The template is used to generate - * the other lines by, for each error text in psImageErrors.dat, the following - * substitutions are made: - * $1 The error text macro name (first word in the psImageErrors.dat lines) - * $2 The error text (rest of the line in psImageErrors.dat) - * $n The order of the source line in psImageErrors.dat (comments excluded) - * - * DO NOT EDIT THE LINES BETWEEN //~Start and //~End! ANY CHANGES WILL BE OVERWRITTEN. - */ - -//~Start #define PS_ERRORTEXT_$1 "$2" -#define PS_ERRORTEXT_psImage_AREA_NEGATIVE "Specified number of rows (%d) or columns (%d) is invalid." -#define PS_ERRORTEXT_psImage_NOT_AN_IMAGE "The input psImage must have a PS_DIMEN_IMAGE dimension type." -#define PS_ERRORTEXT_psImage_IMAGE_NULL "Can not operate on a NULL psImage." -#define PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED "Specified psImage type, %s, is not supported." -#define PS_ERRORTEXT_psImage_INTERPOLATE_METHOD_INVALID "Specified interpolation method (%d) is not supported." -#define PS_ERRORTEXT_psImage_REGION_NULL "Specified psRegion is NULL. Operation could not be performed." -#define PS_ERRORTEXT_psImage_SUBSET_RANGE_INVALID "Specified subset range, [%d:%d,%d:%d], is invalid or outside input psImage's boundaries, [0:%d,0:%d]." -#define PS_ERRORTEXT_psImage_SUBSET_RANGE_MALFORMED "Specified subset range, [%d:%d,%d:%d], is invalid. Ranges must be incremental." -#define PS_ERRORTEXT_psImage_SUBSECTION_NULL "Specified subsection string can not be NULL." -#define PS_ERRORTEXT_psImage_SUBSECTION_INVALID "Specified subsection string, '%s', can not be parsed. Must be in the form '[x1:x2,y1:y2]'." -#define PS_ERRORTEXT_psImage_NOT_PARENT "Specified psImage can not be a child of another psImage." -#define PS_ERRORTEXT_psImage_INPLACE_NOTSUPPORTED "Specified input and output psImage can not reference the same psImage." -#define PS_ERRORTEXT_psImage_SUBSET_ZERO_SIZE "Specified subset, [%d:%d,%d:%d], contains no pixel data." -#define PS_ERRORTEXT_psImage_IMAGE_MASK_SIZE "Input psImage mask size, %dx%d, does not match psImage input size, %dx%d." -#define PS_ERRORTEXT_psImage_IMAGE_MASK_TYPE "Input psImage mask type, %s, is not the supported mask datatype of %s." -#define PS_ERRORTEXT_psImage_BAD_STAT "Specified statistic option, %d, is not valid. Must specify one and only one statistic type." -#define PS_ERRORTEXT_psImage_NO_STAT_OPTIONS "Specified statistic option did not indicate any operation to perform." -#define PS_ERRORTEXT_psImage_STAT_NULL "Specified statistic can not be NULL." -#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, (%g,%g), 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." -#define PS_ERRORTEXT_psImageFFT_FORWARD_NOT_REAL "Input psImage (%s) is not real. Forward FFT operation requires a real psImage input." -#define PS_ERRORTEXT_psImageFFT_FFTW_PLAN_NULL "Could not create a valid FFT plan to perform the transform." -#define PS_ERRORTEXT_psImageFFT_REAL_IMAG_TYPE_MISMATCH "Real psImage type (%s) and imaginary psImage type (%s) must be the same." -#define PS_ERRORTEXT_psImageFFT_REAL_IMAG_SIZE_MISMATCH "Real psImage size (%dx%d) and imaginary psImage size (%dx%d) must be the same." -#define PS_ERRORTEXT_psImageFFT_NONREAL_NOTSUPPORTED "Input psImage type, %s, is required to be either psF32 or psF64." -#define PS_ERRORTEXT_psImageFFT_NONCOMPLEX_NOTSUPPORTED "Input psImage type, %s, is required to be either psC32 or psC64." -#define PS_ERRORTEXT_psImageFFT_REAL_FORWARD_NOTSUPPORTED "The PS_FFT_FORWARD and PS_FFT_REAL_RESULT combinition is not supported." -#define PS_ERRORTEXT_psImageFFT_FORWARD_REVERSE "Can not specify both PS_FFT_FORWARD and PS_FFT_REVERSE options." -#define PS_ERRORTEXT_psImageFFT_NO_DIRECTION_OPTION "Must specify either PS_FFT_FORWARD or PS_FFT_REVERSE option." -#define PS_ERRORTEXT_psImageIO_FILENAME_NULL "Specified filename can not be NULL." -#define PS_ERRORTEXT_psImageIO_FILENAME_INVALID "Could not open file,'%s'.\nCFITSIO Error: %s" -#define PS_ERRORTEXT_psImageIO_EXTNAME_INVALID "Could not find HDU with extension name '%s' in file %s.\nCFITSIO Error: %s" -#define PS_ERRORTEXT_psImageIO_EXTNUM_INVALID "Could not find HDU #%d in file %s.\nCFITSIO Error: %s" -#define PS_ERRORTEXT_psImageIO_DATATYPE_UNKNOWN "Could not determine image data type for file %s.\nCFITSIO Error: %s" -#define PS_ERRORTEXT_psImageIO_IMAGE_DIM_UNKNOWN "Could not determine image dimensions for file %s.\nCFITSIO Error: %s" -#define PS_ERRORTEXT_psImageIO_IMAGE_SIZE_UNKNOWN "Could not determine image size for file %s.\nCFITSIO Error: %s" -#define PS_ERRORTEXT_psImageIO_IMAGE_DIMENSION_UNSUPPORTED "Image number of dimensions, %d, is not valid. Only two or three dimensions supported for FITS I/O." -#define PS_ERRORTEXT_psImageIO_FITS_TYPE_UNSUPPORTED "FITS image type, BITPIX=%d, in file %s is not supported." -#define PS_ERRORTEXT_psImageIO_READ_FAILED "Reading from FITS file %s failed.\nCFITSIO Error: %s" -#define PS_ERRORTEXT_psImageIO_TYPE_UNSUPPORTED "Input psImage type, %s, is not supported." -#define PS_ERRORTEXT_psImageIO_WRITE_EXTNUM_INVALID "Specified extension number, %d, must not exceed number of HDUs, %d, by more than one." -#define PS_ERRORTEXT_psImageIO_FILENAME_CREATE_FAILED "Could not create file,'%s'.\nCFITSIO Error: %s" -#define PS_ERRORTEXT_psImageIO_CREATE_EXTENSION_FAILED "Could not create EXTNAME keyword for file,'%s'.\nCFITSIO Error: %s" -#define PS_ERRORTEXT_psImageIO_CREATE_HDU_FAILED "Could not create HDU for writing psImage in file,'%s'.\nCFITSIO Error: %s" -#define PS_ERRORTEXT_psImageIO_WRITE_FAILED "Could not write psImage data to file,'%s'.\nCFITSIO Error: %s" -#define PS_ERRORTEXT_psImageManip_MAXMIN "Specified min value, %g, can not be greater than the specified max value, %g." -#define PS_ERRORTEXT_psImageManip_MAXMIN_REAL "Specified real-portion of min value, %g, can not be greater than the real-portion of max value, %g." -#define PS_ERRORTEXT_psImageManip_MAXMIN_IMAG "Specified imaginary-portion of min value, %g, can not be greater than the imaginary-portion of max value, %g." -#define PS_ERRORTEXT_psImageManip_OPERATION_NULL "Operation can not be NULL." -#define PS_ERRORTEXT_psImageManip_OVERLAY_TYPE_MISMATCH "Input overlay psImage type, %s, must match input psImage type, %s." -#define PS_ERRORTEXT_psImageManip_CLIP_VALUE_INVALID "Specified %s value, %g%+gi, is not the the range of input psImage's valid pixel values (%s), i.e. [%g:%g]." -#define PS_ERRORTEXT_psImageManip_OVERLAY_OPERATOR_INVALID "Specified operation, '%s', is not supported." -#define PS_ERRORTEXT_psImageManip_SCALE_NOT_POSITIVE "Specified scale value, %d, must be a positive value." -#define PS_ERRORTEXT_psImageManip_INTERPOLATION_MODE_UNSUPPORTED "Specified interpolation mode, %d, is unsupported." -#define PS_ERRORTEXT_psImageConvolve_SHIFT_NULL "Specified shift vectors can not be NULL." -#define PS_ERRORTEXT_psImageConvolve_KERNEL_NULL "Specified psKernel can not be NULL." -#define PS_ERRORTEXT_psImageConvolve_SHIFT_TYPE_MISMATCH "Input t-, x-, and y-shift vector types (%s/%s/%s) must match." -#define PS_ERRORTEXT_psImageConvolve_FFT_FAILED "Failed to perform a fourier transform of input image." -#define PS_ERRORTEXT_psImageConvolve_KERNEL_FFT_FAILED "Failed to perform a fourier transform of kernel." -#define PS_ERRORTEXT_psImageConvolve_KERNEL_TOO_LARGE "Specified psKernel size, %dx%d, can not be larger than input psImage size, %dx%d." -#define PS_ERRORTEXT_psImage_COEFF_NULL "Polynomial coefficients cannot be NULL." -#define PS_ERRORTEXT_psImageManip_TRANSFORM_NULL "Specified input transform can not be NULL" -//~End - -#endif // #ifndef PS_IMAGE_ERRORS_H Index: trunk/psLib/src/image/psImageFFT.c =================================================================== --- trunk/psLib/src/image/psImageFFT.c (revision 4534) +++ (revision ) @@ -1,456 +1,0 @@ -/** @file psImageFFT.c - * - * @brief Contains FFT transform related functions for psImage. - * - * @author Robert DeSonia, MHPCC - * - * @version $Revision: 1.16 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-06-18 03:13:02 $ - * - * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii - */ -#include -#include -#include -#include - -#include "psImageFFT.h" -#include "psError.h" -#include "psMemory.h" -#include "psLogMsg.h" -#include "psImageStructManip.h" - -#include "psImageErrors.h" - -#define PS_FFTW_PLAN_RIGOR FFTW_ESTIMATE - -static psBool p_fftwWisdomImported = false; - -psImage* psImageFFT(psImage* out, const psImage* image, psFFTFlags direction) -{ - psU32 numCols; - psU32 numRows; - psElemType type; - fftwf_plan plan; - - /* got good image data? */ - if (image == NULL) { - psFree(out); - return NULL; - } - - if ( ((direction & PS_FFT_FORWARD) != 0) ) { - if ((direction & PS_FFT_REVERSE) != 0) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImageFFT_FORWARD_REVERSE); - psFree(out); - return NULL; - } - if ((direction & PS_FFT_REAL_RESULT) != 0) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImageFFT_REAL_FORWARD_NOTSUPPORTED); - psFree(out); - return NULL; - } - } else if ((direction & PS_FFT_REVERSE) == 0) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImageFFT_NO_DIRECTION_OPTION); - psFree(out); - return NULL; - } - - type = image->type.type; - - /* make sure the system-level wisdom information is imported. */ - if (!p_fftwWisdomImported) { - fftwf_import_system_wisdom(); - p_fftwWisdomImported = true; - } - - numRows = image->numRows; - numCols = image->numCols; - - // n.b. FFTW can perform a in-place transform at the same rate or faster than out-of-place. - psS32 sign = ((direction & PS_FFT_FORWARD) != 0) ? FFTW_FORWARD : FFTW_BACKWARD; - - fftwf_complex* outBuffer = fftwf_malloc(numRows*numCols*sizeof(fftwf_complex)); - p_psImageCopyToRawBuffer(outBuffer, image, PS_TYPE_C32); - - plan = fftwf_plan_dft_2d(numRows, numCols, - outBuffer, - outBuffer, - sign, - PS_FFTW_PLAN_RIGOR); - - /* check if a plan exists now -- if not, it is a real problem at this point */ - if (plan == NULL) { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImageFFT_FFTW_PLAN_NULL); - psFree(out); - return NULL; - } - - /* finally, call FFTW with the plan made above */ - fftwf_execute(plan); - - fftwf_destroy_plan(plan); - - if (direction & PS_FFT_REAL_RESULT) { - // n.b., we do this instead of using fftwf_plan_dft_c2r because that - // plan requires a half-image, which would require a image reordering - // that is not as simple as performing a normal complex transform and - // then taking the real part of the result. If performance here - // becomes an issue, the use of fftwf_plan_dft_r2c should be considered - // as well as fftwf_plan_dft_c2r. - out = psImageRecycle(out,numCols,numRows,PS_TYPE_F32); - int index = 0; - for (int row=0; row < numRows; row++) { - psF32* outRow = out->data.F32[row]; - for (int col=0; col < numCols; col++) { - outRow[col] = crealf(outBuffer[index++]); // take just the real part - } - } - } else { - out = psImageRecycle(out,numCols,numRows,PS_TYPE_C32); - int index = 0; - for (int row=0; row < numRows; row++) { - psC32* outRow = out->data.C32[row]; - for (int col=0; col < numCols; col++) { - outRow[col] = outBuffer[index++]; // take just the real part - } - } - // memcpy(out->rawDataBuffer, outBuffer, numRows*numCols*sizeof(fftwf_complex)); - } - - fftwf_free(outBuffer); - - return out; - -} - -psImage* psImageReal(psImage* out, const psImage* in) -{ - psElemType type; - psU32 numCols; - psU32 numRows; - - if (in == NULL) { - psFree(out); - return NULL; - } - - type = in->type.type; - numCols = in->numCols; - numRows = in->numRows; - - /* if not a complex number, this is logically just a copy then */ - if (!PS_IS_PSELEMTYPE_COMPLEX(type)) { - return psImageCopy(out, in, type); - } - - if (type == PS_TYPE_C32) { - psF32* outRow; - psC32* inRow; - - out = psImageRecycle(out, numCols, numRows, PS_TYPE_F32); - for (psU32 row = 0; row < numRows; row++) { - outRow = out->data.F32[row]; - inRow = in->data.C32[row]; - - for (psU32 col = 0; col < numCols; col++) { - outRow[col] = crealf(inRow[col]); - } - } - } else if (type == PS_TYPE_C64) { - psF64* outRow; - psC64* inRow; - - out = psImageRecycle(out, numCols, numRows, PS_TYPE_F64); - for (psU32 row = 0; row < numRows; row++) { - outRow = out->data.F64[row]; - inRow = in->data.C64[row]; - - for (psU32 col = 0; col < numCols; col++) { - outRow[col] = creal(inRow[col]); - } - } - } else { - char* typeStr; - PS_TYPE_NAME(typeStr,type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, - typeStr); - - psFree(out); - return NULL; - } - - return out; -} - -psImage* psImageImaginary(psImage* out, const psImage* in) -{ - psElemType type; - psU32 numCols; - psU32 numRows; - - if (in == NULL) { - psFree(out); - return NULL; - } - - type = in->type.type; - numCols = in->numCols; - numRows = in->numRows; - - /* if not a complex image type, this is logically just zeroed image of same size */ - if (!PS_IS_PSELEMTYPE_COMPLEX(type)) { - out = psImageRecycle(out, numCols, numRows, type); - memset(out->data.V[0], 0, PSELEMTYPE_SIZEOF(type) * numCols * numRows); - return out; - } - - if (type == PS_TYPE_C32) { - psF32* outRow; - psC32* inRow; - - out = psImageRecycle(out, numCols, numRows, PS_TYPE_F32); - for (psU32 row = 0; row < numRows; row++) { - outRow = out->data.F32[row]; - inRow = in->data.C32[row]; - - for (psU32 col = 0; col < numCols; col++) { - outRow[col] = cimagf(inRow[col]); - } - } - } else if (type == PS_TYPE_C64) { - psF64* outRow; - psC64* inRow; - - out = psImageRecycle(out, numCols, numRows, PS_TYPE_F64); - for (psU32 row = 0; row < numRows; row++) { - outRow = out->data.F64[row]; - inRow = in->data.C64[row]; - - for (psU32 col = 0; col < numCols; col++) { - outRow[col] = cimag(inRow[col]); - } - } - } else { - char* typeStr; - PS_TYPE_NAME(typeStr,type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, - typeStr); - psFree(out); - return NULL; - } - - return out; -} - -psImage* psImageComplex(psImage* out, const psImage* real, const psImage* imag) -{ - psElemType type; - psU32 numCols; - psU32 numRows; - - if (real == NULL || imag == NULL) { - psFree(out); - return NULL; - } - - type = real->type.type; - numCols = real->numCols; - numRows = real->numRows; - - if (imag->type.type != type) { - char* typeStrReal; - char* typeStrImag; - PS_TYPE_NAME(typeStrReal,type); - PS_TYPE_NAME(typeStrImag,imag->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImageFFT_REAL_IMAG_TYPE_MISMATCH, - typeStrReal,typeStrImag); - psFree(out); - return NULL; - } - - if (imag->numCols != numCols || imag->numRows != numRows) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImageFFT_REAL_IMAG_SIZE_MISMATCH, - numCols, numRows, imag->numCols, imag->numRows); - psFree(out); - return NULL; - } - - if (type == PS_TYPE_F32) { - psC32* outRow; - psF32* realRow; - psF32* imagRow; - - out = psImageRecycle(out, numCols, numRows, PS_TYPE_C32); - - for (psU32 row = 0; row < numRows; row++) { - outRow = out->data.C32[row]; - realRow = real->data.F32[row]; - imagRow = imag->data.F32[row]; - - for (psU32 col = 0; col < numCols; col++) { - outRow[col] = realRow[col] + I * imagRow[col]; - } - } - } else if (type == PS_TYPE_F64) { - psC64* outRow; - psF64* realRow; - psF64* imagRow; - - out = psImageRecycle(out, numCols, numRows, PS_TYPE_C64); - for (psU32 row = 0; row < numRows; row++) { - outRow = out->data.C64[row]; - realRow = real->data.F64[row]; - imagRow = imag->data.F64[row]; - - for (psU32 col = 0; col < numCols; col++) { - outRow[col] = realRow[col] + I * imagRow[col]; - } - } - } else { - char* typeStr; - PS_TYPE_NAME(typeStr,type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImageFFT_NONREAL_NOTSUPPORTED, - typeStr); - psFree(out); - return NULL; - } - - - return out; -} - -psImage* psImageConjugate(psImage* out, const psImage* in) -{ - psElemType type; - psU32 numCols; - psU32 numRows; - - if (in == NULL) { - psFree(out); - return NULL; - } - - type = in->type.type; - numCols = in->numCols; - numRows = in->numRows; - - /* if not a complex image, this is logically just a image copy */ - if (!PS_IS_PSELEMTYPE_COMPLEX(type)) { - return psImageCopy(out, in, type); - } - - if (type == PS_TYPE_C32) { - psC32* outRow; - psC32* inRow; - - out = psImageRecycle(out, numCols, numRows, PS_TYPE_C32); - for (psU32 row = 0; row < numRows; row++) { - outRow = out->data.C32[row]; - inRow = in->data.C32[row]; - - for (psU32 col = 0; col < numCols; col++) { - outRow[col] = crealf(inRow[col]) - I * cimagf(inRow[col]); - } - } - } else if (type == PS_TYPE_C64) { - psC64* outRow; - psC64* inRow; - - out = psImageRecycle(out, numCols, numRows, PS_TYPE_C64); - for (psU32 row = 0; row < numRows; row++) { - outRow = out->data.C64[row]; - inRow = in->data.C64[row]; - - for (psU32 col = 0; col < numCols; col++) { - outRow[col] = creal(inRow[col]) - I * cimag(inRow[col]); - } - } - } else { - char* typeStr; - PS_TYPE_NAME(typeStr,type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImageFFT_NONCOMPLEX_NOTSUPPORTED, - typeStr); - psFree(out); - return NULL; - } - - return out; -} - -psImage* psImagePowerSpectrum(psImage* out, const psImage* in) -{ - psElemType type; - psU32 numCols; - psU32 numRows; - - if (in == NULL) { - psFree(out); - return NULL; - } - - type = in->type.type; - numCols = in->numCols; - numRows = in->numRows; - - if (type == PS_TYPE_C32) { - psF32* outRow; - psC32* inRow; - psF32 real; - psF32 imag; - psF32 fNumCols = numCols; - psF32 fNumRows = numRows; - psF32 numElementsSquared = fNumCols * fNumCols * fNumRows * fNumRows; - - out = psImageRecycle(out, numCols, numRows, PS_TYPE_F32); - for (psU32 row = 0; row < numRows; row++) { - outRow = out->data.F32[row]; - inRow = in->data.C32[row]; - - for (psU32 col = 0; col < numCols; col++) { - real = crealf(inRow[col]); - imag = cimagf(inRow[col]); - outRow[col] = (real * real + imag * imag) / numElementsSquared; - } - } - } else if (type == PS_TYPE_C64) { - psF64* outRow; - psC64* inRow; - psF64 real; - psF64 imag; - psF64 numElementsSquared = numCols * numCols * numRows * numRows; - - out = psImageRecycle(out, numCols, numRows, PS_TYPE_F64); - for (psU32 row = 0; row < numRows; row++) { - outRow = out->data.F64[row]; - inRow = in->data.C64[row]; - - for (psU32 col = 0; col < numCols; col++) { - real = creal(inRow[col]); - imag = cimag(inRow[col]); - outRow[col] = (real * real + imag * imag) / numElementsSquared; - } - } - } else { - char* typeStr; - PS_TYPE_NAME(typeStr,type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImageFFT_NONCOMPLEX_NOTSUPPORTED, - typeStr); - psFree(out); - return NULL; - } - - return out; - -} Index: trunk/psLib/src/image/psImageFFT.h =================================================================== --- trunk/psLib/src/image/psImageFFT.h (revision 4534) +++ (revision ) @@ -1,87 +1,0 @@ -/** @file psImageFFT.h - * - * @brief Contains FFT transform related functions for psImage - * - * @ingroup Transform - * - * @author Robert DeSonia, MHPCC - * - * @version $Revision: 1.6 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-07-11 21:38:19 $ - * - * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii - */ - -#ifndef PS_IMAGE_FFT_H -#define PS_IMAGE_FFT_H - -#include "psImage.h" -#include "psVectorFFT.h" // for psFFTFlags - -/// @addtogroup Transform -/// @{ - -/** Forward and reverse FFT calculations. - * - * This takes as input the image of interest (in) and the direction - * (direction), which is specified by an enumerated type psFftDirection. - * The input image may be of type psF32 or psC32, the result is always - * psC32. If the input vector is psF32, the direction must be forward. - * - * @return psImage* the FFT transformation result - */ -psImage* psImageFFT( - psImage* out, ///< a psImage to recycle. If NULL, a new psImage is made. - const psImage* image, ///< the psImage to apply transform to - psFFTFlags direction ///< the direction of the transform -); - -/** extract the real portion of a complex image - * - * @return psImage* real portion of the input image. - */ -psImage* psImageReal( - psImage* out, ///< a psImage to recycle. If NULL, a new psImage is made. - const psImage* in ///< the psImage to extract real portion from -); - -/** extract the imaginary portion of a complex image - * - * @return psImage* imaginary portion of the input image. - */ -psImage* psImageImaginary( - psImage* out, ///< a psImage to recycle. If NULL, a new psImage is made. - const psImage* in ///< the psImage to extract imaginary portion from -); - -/** creates a complex image from separate real and imaginary plane images - * - * @return psImage* resulting complex image - */ -psImage* psImageComplex( - psImage* out, ///< a psImage to recycle. If NULL, a new psImage is made. - const psImage* real, ///< the real plane image - const psImage* imag ///< the imaginary plane image -); - -/** computes the complex conjugate of an image - * - * @return psImage* the complex conjugate of the 'in' image - */ -psImage* psImageConjugate( - psImage* out, ///< a psImage to recycle. If NULL, a new psImage is made. - const psImage* in ///< the psImage to compute conjugate of -); - -/** computes the power spectrum of an image - * - * @return psImage* the power spectrum of the 'in' image - */ -psImage* psImagePowerSpectrum( - psImage* out, ///< a psImage to recycle. If NULL, a new psImage is made. - const psImage* in ///< the psImage to power spectrum of -); - -/// @} - -#endif // #ifndef PS_IMAGE_FFT_H Index: trunk/psLib/src/image/psImageGeomManip.c =================================================================== --- trunk/psLib/src/image/psImageGeomManip.c (revision 4534) +++ (revision ) @@ -1,871 +1,0 @@ -/** @file psImageGeomManip.c - * - * @brief Contains basic image pixel and geometry manipulation operations, as - * specified in the PSLIB SDRS sections "Image Pixel Manipulations" and - * "Image Geometry Manipulations". - * - * @ingroup Image - * - * @author Robert DeSonia, MHPCC - * @author Ross Harman, MHPCC - * - * @version $Revision: 1.12 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-07-11 21:38:19 $ - * - * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii - */ - -#include -#include // for isfinite(), etc. -#include -#include // for memcpy, etc. - -#include "psImageGeomManip.h" - -#include "psError.h" -#include "psImage.h" -#include "psImageStructManip.h" -#include "psStats.h" -#include "psMemory.h" -#include "psConstants.h" -#include "psImageErrors.h" -#include "psCoord.h" - -psImage* psImageRebin(psImage* out, - const psImage* in, - const psImage* restrict mask, - psMaskType maskVal, - int scale, - const psStats* stats) -{ - psS32 inRows; - psS32 inCols; - psS32 outRows; - psS32 outCols; - psVector* vec; // vector to hold the values of a single bin. - psVector* maskVec = NULL; // vector to hold the mask of a single bin. - psMaskType* maskData = NULL; - psStats* myStats; - double statVal; - - if (in == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - psFree(out); - return NULL; - } - - if (scale < 1) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImageManip_SCALE_NOT_POSITIVE, - scale); - psFree(out); - return NULL; - } - - if (stats == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_STAT_NULL); - psFree(out); - return NULL; - } - - if (p_psGetStatValue(stats, &statVal) == false) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_BAD_STAT, - stats->options); - psFree(out); - return NULL; - } - - vec = psVectorAlloc(scale * scale, in->type.type); - - if (mask != NULL) { - if (mask->type.type != PS_TYPE_MASK) { - char* typeStr; - PS_TYPE_NAME(typeStr,mask->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_MASK_TYPE, - typeStr, PS_TYPE_MASK_NAME); - psFree(out); - psFree(vec); - return NULL; - } - maskVec = psVectorAlloc(scale * scale, PS_TYPE_MASK); - maskData = maskVec->data.PS_TYPE_MASK_DATA; - } - - myStats = psAlloc(sizeof(psStats)); - *myStats = *stats; - - // create output image. - inRows = in->numRows; - inCols = in->numCols; - outRows = (inRows + scale - 1) / scale; // round-up for remainders - outCols = (inCols + scale - 1) / scale; // round-up for remainders - out = psImageRecycle(out, outCols, outRows, in->type.type); - - #define PS_IMAGE_REBIN_CASE(type) \ -case PS_TYPE_##type: { \ - ps##type* outRowData; \ - ps##type* vecData = vec->data.type; \ - psMaskType* inRowMask = NULL; \ - for (psS32 row = 0; row < outRows; row++) { \ - outRowData = out->data.type[row]; \ - psS32 inCurrentRow = row*scale; \ - psS32 inNextRow = (row+1)*scale; \ - for (psS32 col = 0; col < outCols; col++) { \ - psS32 inCurrentCol = col*scale; \ - psS32 inNextCol = (col+1)*scale; \ - psS32 n = 0; \ - for (psS32 inRow = inCurrentRow; inRow < inNextRow && inRow < inRows; inRow++) { \ - ps##type* inRowData = in->data.type[inRow]; \ - if (mask != NULL) { \ - inRowMask = mask->data.PS_TYPE_MASK_DATA[inRow]; \ - } \ - for (psS32 inCol = inCurrentCol; inCol < inNextCol && inCol < inCols; inCol++) { \ - if (maskData != NULL) { \ - maskData[n] = inRowMask[inCol]; \ - } \ - vecData[n++] = inRowData[inCol]; \ - } \ - } \ - vec->n = n; \ - myStats = psVectorStats(myStats,vec,NULL,maskVec,maskVal); \ - p_psGetStatValue(myStats,&statVal); \ - outRowData[col] = (ps##type)statVal; \ - } \ - } \ - } \ - break; - - switch (in->type.type) { - // PS_IMAGE_REBIN_CASE(U8); Not valid since psVectorStats doesn't allow - PS_IMAGE_REBIN_CASE(U16); - PS_IMAGE_REBIN_CASE(U32); // Not a requirement - PS_IMAGE_REBIN_CASE(U64); // Not a requirement - PS_IMAGE_REBIN_CASE(S8); - // PS_IMAGE_REBIN_CASE(S16); Not valid since psVectorStats doesn't allow - PS_IMAGE_REBIN_CASE(S32); // Not a requirement - PS_IMAGE_REBIN_CASE(S64); // Not a requirement - PS_IMAGE_REBIN_CASE(F32); - PS_IMAGE_REBIN_CASE(F64); - // PS_IMAGE_REBIN_CASE(C32); Not valid since psVectorStats doesn't allow - // PS_IMAGE_REBIN_CASE(C64); Not valid since psVectorStats doesn't allow - - default: { - char* typeStr; - PS_TYPE_NAME(typeStr,in->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, - typeStr); - psFree(out); - out = NULL; - } - } - - psFree(vec); - psFree(maskVec); - psFree(myStats); - - return out; -} - -psImage* psImageResample(psImage* out, - const psImage* in, - int scale, - psImageInterpolateMode mode) -{ - psS32 outRows; - psS32 outCols; - float invScale; - - if (in == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - psFree(out); - return NULL; - } - - if (scale < 1) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImageManip_SCALE_NOT_POSITIVE, - scale); - psFree(out); - return NULL; - } - - if (mode > PS_INTERPOLATE_LANCZOS4_VARIANCE ) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImageManip_INTERPOLATION_MODE_UNSUPPORTED, - mode); - psFree(out); - return NULL; - } - - // create an output image of the same size - // and type - outRows = in->numRows * scale; - outCols = in->numCols * scale; - invScale = 1.0f / (float)scale; - - #define PSIMAGE_RESAMPLE_CASE(TYPE) \ -case PS_TYPE_##TYPE: { \ - out = psImageRecycle(out,outCols, outRows, PS_TYPE_##TYPE); \ - for (psS32 row=0;rowdata.TYPE[row]; \ - float inRow = (float)row * invScale; \ - for (psS32 col=0;coltype.type) { - PSIMAGE_RESAMPLE_CASE(U8) - PSIMAGE_RESAMPLE_CASE(U16) - PSIMAGE_RESAMPLE_CASE(U32) - PSIMAGE_RESAMPLE_CASE(U64) - PSIMAGE_RESAMPLE_CASE(S8) - PSIMAGE_RESAMPLE_CASE(S16) - PSIMAGE_RESAMPLE_CASE(S32) - PSIMAGE_RESAMPLE_CASE(S64) - PSIMAGE_RESAMPLE_CASE(F32) - PSIMAGE_RESAMPLE_CASE(F64) - PSIMAGE_RESAMPLE_CASE(C32) - PSIMAGE_RESAMPLE_CASE(C64) - default: { - char* typeStr; - PS_TYPE_NAME(typeStr,in->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, - typeStr); - psFree(out); - out = NULL; - } - } - - return out; -} - -psImage* psImageRoll(psImage* out, - const psImage* input, - int dx, - int dy) -{ - psS32 outRows; - psS32 outCols; - psS32 elementSize; - - if (input == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - psFree(out); - return NULL; - } - // create an output image of the same size - // and type - outRows = input->numRows; - outCols = input->numCols; - elementSize = PSELEMTYPE_SIZEOF(input->type.type); - out = psImageRecycle(out, outCols, outRows, input->type.type); - - // make dx and dy between 0 and outCols or - // outRows, respectively - dx = dx % outCols; - dy = dy % outRows; - if (dx < 0) { - dx += outCols; - } - if (dy < 0) { - dy += outRows; - } - - psS32 segment1Size = elementSize * (outCols - dx); - psS32 segment2Size = elementSize * dx; - - for (psS32 row = 0; row < outRows; row++) { - psS32 inRowNumber = row + dy; - - if (inRowNumber >= outRows) { - inRowNumber -= outRows; - } - psU8* inRow = input->data.U8[inRowNumber]; // use byte arithmetic for all types - psU8* outRow = out->data.U8[row]; - - memcpy(outRow, inRow + segment2Size, segment1Size); - memcpy(outRow + segment1Size, inRow, segment2Size); - } - - return out; -} - -psImage* psImageRotate(psImage* out, - const psImage* input, - float angle, - complex exposed, - psImageInterpolateMode mode) -{ - if (input == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - psFree(out); - return NULL; - } - // put the angle in the range of 0...2PI. - angle = (float)((double)angle - (2.0*M_PI) * floor(angle / (2.0*M_PI))); - - if (fabsf(angle - M_PI_2) < FLT_EPSILON) { - // perform 1/4 rotate counter-clockwise - psS32 numRows = input->numCols; - psS32 numCols = input->numRows; - psS32 lastCol = numCols - 1; - psElemType type = input->type.type; - - out = psImageRecycle(out, numCols, numRows, type); - - #define PSIMAGE_ROTATE_LEFT_90(TYPE) \ - case PS_TYPE_##TYPE: { \ - ps##TYPE** inData = input->data.TYPE; \ - for (psS32 row=0;rowdata.TYPE[row]; \ - for (psS32 col=0;colnumRows; - psS32 lastRow = numRows - 1; - psS32 numCols = input->numCols; - psS32 lastCol = numCols - 1; - psElemType type = input->type.type; - - out = psImageRecycle(out, numCols, numRows, type); - - #define PSIMAGE_ROTATE_180_CASE(TYPE) \ - case PS_TYPE_##TYPE: { \ - for (psS32 row=0;rowdata.TYPE[row]; \ - ps##TYPE* inRow = input->data.TYPE[lastRow-row]; \ - for (psS32 col=0;colnumCols; - psS32 lastRow = numRows - 1; - psS32 numCols = input->numRows; - psElemType type = input->type.type; - - out = psImageRecycle(out, numCols, numRows, type); - - #define PSIMAGE_ROTATE_RIGHT_90(TYPE) \ - case PS_TYPE_##TYPE: { \ - ps##TYPE** inData = input->data.TYPE; \ - for (psS32 row=0;rowdata.TYPE[row]; \ - for (psS32 col=0;coltype.type); - } else { - psElemType type = input->type.type; - psS32 numRows = input->numRows; - psS32 numCols = input->numCols; - float centerX = (float)(numCols) / 2.0f; - float centerY = (float)(numRows) / 2.0f; - double cosT = cosf(angle); - double sinT = sinf(angle); - - // calculate the corners of the rotated image so we know the proper output image size. - // x' = x cos(t) + y sin(t); i.e, x' = (x-centerX)*cosT + (y-centerY)*sinT; - // y' = y cos(t) - x sin(t); i.e. y' = (y-centerY)*cosT - (x-centerX)*sinT; - - psS32 outCols = ceil(abs(numCols * cosT) + abs(numRows * sinT)) + 1; - psS32 outRows = ceil(abs(numCols * sinT) + abs(numRows * cosT)) + 1; - float minX = (float)outCols / -2.0f; - psS32 intMinY = outRows / -2; - - out = psImageRecycle(out, outCols, outRows, type); - - /* optimized public domain rotation routine by Karl Lager - * - * float cosT,sinT; - * cosT = cos(t); - * sinT = sin(t); - * for (y = min_y; y <= max_y; y++) { - * x' = min_x * cosT + y * sinT + x1'; - * y' = y * cosT - min_x * sinT + y1'; - * for (x = min_x; x <= max_x; x++) { - * if (x', y') is in the bounds of the bitmap, get pixel - * (x', y') and plot the pixel to (x, y) on screen. - * x' += cosT; - * y' -= sinT; - * } - * } - */ - - // precalculate some figures that are used within loop - float minXTimesCosTPlusCenterX = minX * cosT + centerX; - float CenterYMinusminXTimesSinT = centerY - minX * sinT; - - #define PSIMAGE_ROTATE_ARBITRARY_LOOP(TYPE,MODE) { \ - if (creal(exposed) < PS_MIN_##TYPE || \ - creal(exposed) > PS_MAX_##TYPE || \ - cimag(exposed) < PS_MIN_##TYPE || \ - cimag(exposed) > PS_MAX_##TYPE) { \ - psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ - PS_ERRORTEXT_psImageManip_CLIP_VALUE_INVALID, \ - "exposed", \ - creal(exposed),cimag(exposed), \ - PS_TYPE_##TYPE##_NAME, \ - (double)PS_MIN_##TYPE,(double)PS_MAX_##TYPE); \ - psFree(out); \ - out = NULL; \ - break; \ - } \ - float inX; \ - float inY; \ - ps##TYPE* outRow; \ - for (psS32 y = 0; y < outRows; y++) { \ - inX = minXTimesCosTPlusCenterX + (y+intMinY) * sinT; \ - inY = CenterYMinusminXTimesSinT + (y+intMinY) * cosT; \ - outRow = out->data.TYPE[y]; \ - for (psS32 x = 0; x < outCols; x++) { \ - outRow[x] = p_psImagePixelInterpolate##MODE##_##TYPE(input,inX,inY,NULL,0,exposed); \ - inX += cosT; \ - inY -= sinT; \ - } \ - } \ - } - - #define PSIMAGE_ROTATE_ARBITRARY_CASE(MODE) \ - case PS_INTERPOLATE_##MODE: \ - switch (type) { \ - case PS_TYPE_U8: \ - PSIMAGE_ROTATE_ARBITRARY_LOOP(U8,MODE); \ - break; \ - case PS_TYPE_U16: \ - PSIMAGE_ROTATE_ARBITRARY_LOOP(U16,MODE); \ - break; \ - case PS_TYPE_U32: /* Not a requirement */ \ - PSIMAGE_ROTATE_ARBITRARY_LOOP(U32,MODE); \ - break; \ - case PS_TYPE_U64: /* Not a requirement */ \ - PSIMAGE_ROTATE_ARBITRARY_LOOP(U64,MODE); \ - break; \ - case PS_TYPE_S8: \ - PSIMAGE_ROTATE_ARBITRARY_LOOP(S8,MODE); \ - break; \ - case PS_TYPE_S16: \ - PSIMAGE_ROTATE_ARBITRARY_LOOP(S16,MODE); \ - break; \ - case PS_TYPE_S32: /* Not a requirement */ \ - PSIMAGE_ROTATE_ARBITRARY_LOOP(S32,MODE); \ - break; \ - case PS_TYPE_S64: /* Not a requirement */ \ - PSIMAGE_ROTATE_ARBITRARY_LOOP(S64,MODE); \ - break; \ - case PS_TYPE_F32: \ - PSIMAGE_ROTATE_ARBITRARY_LOOP(F32,MODE); \ - break; \ - case PS_TYPE_F64: \ - PSIMAGE_ROTATE_ARBITRARY_LOOP(F64,MODE); \ - break; \ - case PS_TYPE_C32: \ - PSIMAGE_ROTATE_ARBITRARY_LOOP(C32,MODE); \ - break; \ - case PS_TYPE_C64: \ - PSIMAGE_ROTATE_ARBITRARY_LOOP(C64,MODE); \ - break; \ - default: { \ - char* typeStr; \ - PS_TYPE_NAME(typeStr,type); \ - psError(PS_ERR_BAD_PARAMETER_TYPE, true, \ - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, \ - typeStr); \ - psFree(out); \ - out = NULL; \ - } \ - } \ - break; - - switch (mode) { - PSIMAGE_ROTATE_ARBITRARY_CASE(FLAT); - PSIMAGE_ROTATE_ARBITRARY_CASE(BILINEAR); - PSIMAGE_ROTATE_ARBITRARY_CASE(BILINEAR_VARIANCE); - default: - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImageManip_INTERPOLATION_MODE_UNSUPPORTED, - mode); - psFree(out); - out = NULL; - } - } - - return out; -} - -psImage* psImageShift(psImage* out, - const psImage* input, - float dx, - float dy, - complex exposed, - psImageInterpolateMode mode) -{ - psS32 outRows; - psS32 outCols; - psS32 elementSize; - psElemType type; - - if (input == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - psFree(out); - return NULL; - } - // create an output image of the same size - // and type - outRows = input->numRows; - outCols = input->numCols; - type = input->type.type; - elementSize = PSELEMTYPE_SIZEOF(type); - out = psImageRecycle(out, outCols, outRows, type); - - #define PSIMAGE_SHIFT_CASE(MODE,TYPE) \ -case PS_TYPE_##TYPE: \ - if (creal(exposed) < PS_MIN_##TYPE || \ - creal(exposed) > PS_MAX_##TYPE || \ - cimag(exposed) < PS_MIN_##TYPE || \ - cimag(exposed) > PS_MAX_##TYPE) { \ - psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ - PS_ERRORTEXT_psImageManip_CLIP_VALUE_INVALID, \ - "exposed", \ - creal(exposed),cimag(exposed), \ - PS_TYPE_##TYPE##_NAME, \ - (double)PS_MIN_##TYPE,(double)PS_MAX_##TYPE); \ - psFree(out); \ - out = NULL; \ - break; \ - } \ - for (psS32 row=0;rowdata.TYPE[row]; \ - float y = dy+(float)row; \ - for (psS32 col=0;coltype.type) { \ - PSIMAGE_SHIFT_CASE(MODE,U8); \ - PSIMAGE_SHIFT_CASE(MODE,U16); \ - PSIMAGE_SHIFT_CASE(MODE,U32); /* Not a requirement */ \ - PSIMAGE_SHIFT_CASE(MODE,U64); /* Not a requirement */ \ - PSIMAGE_SHIFT_CASE(MODE,S8); \ - PSIMAGE_SHIFT_CASE(MODE,S16); \ - PSIMAGE_SHIFT_CASE(MODE,S32); /* Not a requirement */ \ - PSIMAGE_SHIFT_CASE(MODE,S64); /* Not a requirement */ \ - PSIMAGE_SHIFT_CASE(MODE,F32); \ - PSIMAGE_SHIFT_CASE(MODE,F64); \ - PSIMAGE_SHIFT_CASE(MODE,C32); \ - PSIMAGE_SHIFT_CASE(MODE,C64); \ - \ - default: { \ - char* typeStr; \ - PS_TYPE_NAME(typeStr,type); \ - psError(PS_ERR_BAD_PARAMETER_TYPE, true, \ - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, \ - typeStr); \ - psFree(out); \ - out = NULL; \ - } \ - } \ - break; - - switch (mode) { - PSIMAGE_SHIFT_ARBITRARY_CASE(FLAT); - PSIMAGE_SHIFT_ARBITRARY_CASE(BILINEAR); - PSIMAGE_SHIFT_ARBITRARY_CASE(BILINEAR_VARIANCE); - default: - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImageManip_INTERPOLATION_MODE_UNSUPPORTED, - mode); - psFree(out); - out = NULL; - } - - return out; -} - -psImage* psImageTransform(psImage *output, - psPixels* blankPixels, - const psImage *input, - const psImage *inputMask, - psMaskType inputMaskVal, - const psPlaneTransform *outToIn, - psRegion region, - const psPixels* pixels, - psImageInterpolateMode mode, - double exposedValue) -{ - if (input == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - psFree(output); - return NULL; - } - psElemType type = input->type.type; - - if (inputMask != NULL) { - if (input->numRows != inputMask->numRows || input->numCols != inputMask->numCols) { - psError(PS_ERR_BAD_PARAMETER_SIZE, true, - PS_ERRORTEXT_psImage_IMAGE_MASK_SIZE, - input->numCols, input->numRows, - inputMask->numCols, inputMask->numRows ); - psFree(output); - return NULL; - } - if (inputMask->type.type != PS_TYPE_MASK) { - char* typeStr; - PS_TYPE_NAME(typeStr,inputMask->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_MASK_TYPE, - typeStr, PS_TYPE_MASK_NAME); - psFree(output); - return NULL; - } - } - - if (outToIn == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImageManip_TRANSFORM_NULL); - return NULL; - } - - int row0; - int row1; - int col0; - int col1; - int numRows; - int numCols; - if (output == NULL) { // output image size is determined by psRegion - row0 = region.y0; - row1 = region.y1; - col0 = region.x0; - col1 = region.x1; - if (col1 < 1) { - col1 += input->numCols; - } - - if (row1 < 1) { - row1 += input->numRows; - } - - numRows = row1 - row0; - numCols = col1 - col0; - - if (numRows < 1 || numCols < 1) { - psError(PS_ERR_BAD_PARAMETER_SIZE, true, - "The specified region is invalid."); - psFree(output); - return NULL; - } - // create the output image. - output = psImageRecycle(output, numCols, numRows, input->type.type); - *(psS32*)&output->col0 = region.x0; - *(psS32*)&output->row0 = region.y0; - } else { // size of output is determined by output parameter - numRows = output->numRows; - numCols = output->numCols; - row0 = output->row0; - col0 = output->col0; - row1 = row0+numRows; - col1 = col0+numCols; - } - - // loop through the output image using the domain above and transform - // each output pixel to input coordinates and use psImagePixelInterpolate - // to determine the pixel value. - psPlane outPosition; - psPlane* inPosition = NULL; - - - - - - #define PSIMAGE_TRANSFORM_DOTRANSFORM(TYPE,MODE) \ - /* apply the transform to get the position in the input image */ \ - inPosition = psPlaneTransformApply(inPosition, outToIn, &outPosition); \ - \ - if (inPosition == NULL) { \ - psError(PS_ERR_UNKNOWN, false, \ - "Failed to apply the transform"); \ - psFree(output); \ - return NULL; \ - } \ - /* interpolate the cooresponding input pixel to get the output pixel value. */ \ - ps##TYPE value = p_psImagePixelInterpolate##MODE##_##TYPE(input, \ - inPosition->x, inPosition->y, \ - inputMask, inputMaskVal, NAN); \ - /* psFree(inPosition); */\ - if (isnan(value)) { \ - if (blankPixels != NULL) { \ - p_psPixelsAppend(blankPixels, blankPixels->nalloc, outPosition.x, outPosition.y); \ - } \ - value = exposedValue; \ - } \ - - #define PSIMAGE_TRANSFORM_LOOP(TYPE, MODE) { \ - for (int row = 0; row < numRows; row++) { \ - outPosition.y = row+row0; \ - ps##TYPE* outputData=output->data.TYPE[row]; \ - for (int col = 0; col < numCols; col++) { \ - outPosition.x = col+col0; \ - PSIMAGE_TRANSFORM_DOTRANSFORM(TYPE,MODE) \ - outputData[col] = value; \ - } \ - } \ - } - - #define PSIMAGE_TRANSFORM_FROMLIST(TYPE, MODE) { \ - int n = pixels->n; \ - for (int i= 0; i < n; i++) { \ - int x = pixels->data[i].x; \ - int y = pixels->data[i].y; \ - if (x >= col0 && x < col1 && y >= row0 && y < row1) { \ - outPosition.x = x; \ - outPosition.y = y; \ - PSIMAGE_TRANSFORM_DOTRANSFORM(TYPE,MODE) \ - output->data.TYPE[y][x] = value; \ - } \ - } \ - } - - #define PSIMAGE_TRANSFORM_CASE(MODE) \ -case PS_INTERPOLATE_##MODE: \ - switch (type) { \ - case PS_TYPE_F32: \ - PSIMAGE_TRANSFORM_LOOP(F32,MODE); \ - break; \ - case PS_TYPE_F64: \ - PSIMAGE_TRANSFORM_LOOP(F64,MODE); \ - break; \ - default: { \ - char* typeStr; \ - PS_TYPE_NAME(typeStr,type); \ - psError(PS_ERR_BAD_PARAMETER_TYPE, true, \ - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, \ - typeStr); \ - psFree(output); \ - return NULL; \ - } \ - } \ - break; - - switch (mode) { - PSIMAGE_TRANSFORM_CASE(FLAT); - PSIMAGE_TRANSFORM_CASE(BILINEAR); - PSIMAGE_TRANSFORM_CASE(BILINEAR_VARIANCE); - default: - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImageManip_INTERPOLATION_MODE_UNSUPPORTED, - mode); - psFree(output); - return NULL; - } - - psFree(inPosition); - - return output; -} - Index: trunk/psLib/src/image/psImageGeomManip.h =================================================================== --- trunk/psLib/src/image/psImageGeomManip.h (revision 4534) +++ (revision ) @@ -1,158 +1,0 @@ -/** @file psImageGeomManip.h - * - * @brief Contains basic image geometry manipulation operations, as - * specified in the PSLIB SDRS sections "Image Geometry Manipulations". - * - * @ingroup Image - * - * @author Robert DeSonia, MHPCC - * - * @version $Revision: 1.9 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-06-25 00:51:28 $ - * - * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii - */ -#ifndef PS_IMAGE_GEOM_MANIP_H -#define PS_IMAGE_GEOM_MANIP_H - -#include "psImage.h" -#include "psCoord.h" -#include "psStats.h" -#include "psPixels.h" - -/// @addtogroup Image -/// @{ - -/** Rebin image to new scale. - * - * A new image is constructed in which the dimensions are reduced by a factor of - * 1/scale. The scale, always a positive number, is equal in each dimension and - * specified the number of pixels used to define a new pixel in the output image. - * The output image is generated from all input image pixels. This function is - * defined for psU8, psS8, psS16, psF32, psF64, psC32, and psC64. - * - * @return psImage new image formed by rebinning input image. - */ -psImage* psImageRebin( - psImage* out, ///< an psImage to recycle. If NULL, a new image is created - const psImage* in, ///< input image - const psImage* mask, ///< mask for input image. If NULL, no masking is done. - psMaskType maskVal, ///< the bits to check in mask. - int scale, ///< the scale to rebin for each dimension - const psStats* stats - ///< the statistic to perform when rebinning. Only one method should be set. -); - -/** Resample image to new scale. - * - * A new image is constructed in which the dimensions are increased by a - * factor of scale. The scale, always a positive number, is equal in each - * dimension. The output image is generated from all input image pixels. - * Each pixel in the output image is derived by interpolating between - * neighboring pixels using the specified interpolation method (mode). - * - * @return psImage* resampled image result - */ -psImage* psImageResample( - psImage* out, ///< an psImage to recycle. If NULL, a new image is created - const psImage* in, ///< input image - int scale, ///< resample scaling factor - psImageInterpolateMode mode ///< the interpolation mode used in resampling -); - -/** Rotate the input image by given angle, specified in degrees. - * - * The output image must contain all of the pixels from the input image in - * their new frame. Pixels in the output image which do not map to input - * pixels should be set to exposed. The center of rotation is always the - * center pixel of the image. The rotation is specified in the sense that a - * positive angle is an anti-clockwise rotation. This function must be - * defined for the following types: psU8, psU16, psS8, psS16, psF32, psF64, - * psC32, psC64. - * - * @return psImage* the rotated image result. - */ -psImage* psImageRotate( - psImage* out, ///< an psImage to recycle. If NULL, a new image is created - const psImage* input, ///< input image - float angle, ///< the rotation angle in radians. - complex exposed, ///< the output image pixel values for non-imagery areas - psImageInterpolateMode mode ///< the interpolation mode used -); - -/** Shift image by an arbitrary number of pixels (dx,dy) in either direction. - * - * If the shift values are fractional, the output pixel values should - * interpolate between the input pixel values. The output image has the same - * dimensions as the input image. Pixels which fall off the edge of the - * output image are lost. Newly exposed pixels are set to the value given by - * exposed. This function must be defined for the following types: psU8, - * psU16, psS8, psS16, psF32, psF64, psC32, psC64. - * - * @return psImage* the shifted image result. - */ -psImage* psImageShift( - psImage* out, ///< an psImage to recycle. If NULL, a new image is created - const psImage* input, ///< input image - float dx, ///< the shift in x direction. - float dy, ///< the shift in y direction. - complex exposed, ///< the output image pixel values for non-imagery areas - psImageInterpolateMode mode ///< the interpolation mode to use -); - -/** Roll image by an integer number of pixels in either direction. - * - * The output image is the same dimensions as the input image. Edge pixels - * wrap to the other side (no values are lost). This function is - * defined for psU8, psS8, psS16, psF32, psF64, psC32, and psC64. - * - * @return psImage* the rolled version of the input image. - */ -psImage* psImageRoll( - psImage* out, ///< an psImage to recycle. If NULL, a new image is created - const psImage* input, ///< input image - int dx, ///< number of pixels to roll in the x-dimension - int dy ///< number of pixels to roll in the y-dimension -); - -/** Transform the input image according the supplied transformation. - * - * Transform the input image according the supplied transformation. The size - * of the transformed image is defined by the supplied output image, if - * non-NULL, or the region otherwise (size region.x1 - region.x0 by region.y1 - * region.y0, with out->x0 = region.x0 and out->y0 = region.y0). If the - * inputMask is non-NULL, those pixels in the inputMask matching inputMaskVal - * are to be ignored in the transformation. The inputMask must be of type - * psU8, and of the same size as the input, otherwise the function shall - * generate an error and return NULL. The transformation outToIn specifies the - * coordinates in the input image of a pixel in the output image — note that - * this is the reverse of what might be naively expected, but it is what is - * required in order to use psImagePixelInterpolate. If the pixels array is - * non-NULL, it shall consist of psPixelCoords, and only those pixels in the - * output image shall be transformed; otherwise, the entire image is - * generated. The interpolation is performed using the specified interpolation - * mode. Where a pixel in the output image does not correspond to a pixel in - * the input image (or all appropriate pixels in the input image are - * masked), the value shall be set to exposed, and the pixel added to the - * appropriate list of pixels (psPixels) in the array of blankPixels for - * return to the user. This function must be capable of handling the following - * types for the input (with corresponding types for the output): psF32, psF64. - - * - * @return psImage* The transformed image. - */ -psImage* psImageTransform( - psImage *output, ///< psImage to recycle, or NULL - psPixels* blankPixels, ///< list of pixels in output image not set, or NULL if no list is desired. - const psImage *input, ///< psImage to apply transform to - const psImage *inputMask, ///< if not NULL, mask of input psImage - psMaskType inputMaskVal, ///< masking value for inputMask - const psPlaneTransform *outToIn, ///< the transform to apply - psRegion region, ///< the size of the transformed image - const psPixels* pixels, /**< if not NULL, consists of psPixelCoords and specifies which pixels in - * output image shall be transformed; otherwise, entire image generated*/ - psImageInterpolateMode mode, ///< the interpolation scheme to be used - double exposedValue ///< Exposed value to which non-corresponding pixels are set -); - -#endif // #ifndef PS_IMAGE_GEOM_MANIP_H Index: trunk/psLib/src/image/psImagePixelExtract.c =================================================================== --- trunk/psLib/src/image/psImagePixelExtract.c (revision 4534) +++ (revision ) @@ -1,600 +1,0 @@ -/** @file psImagePixelExtract.c - * - * @brief Contains basic image extraction operations, as specified in the - * PSLIB SDRS sections "Image Pixel Extractions". - * - * @ingroup Image - * - * @author Robert DeSonia, MHPCC - * - * @version $Revision: 1.5 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-07-11 21:38:19 $ - * - * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii - * - */ - -#include - -#include "psMemory.h" -#include "psImagePixelExtract.h" -#include "psError.h" - -#include "psImageErrors.h" - -psVector* psImageSlice(psVector* out, - psVector* coords, - const psImage* restrict input, - const psImage* restrict mask, - psMaskType maskVal, - psRegion region, - psImageCutDirection direction, - const psStats* stats) -{ - double statVal; - psStats* myStats; - psElemType type; - psS32 inRows; - psS32 inCols; - psS32 delta = 1; - psF64* outData; - psS32 row0 = region.y0; - psS32 row1 = region.y1; - psS32 col0 = region.x0; - psS32 col1 = region.x1; - - if (input == NULL || input->data.V == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - psFree(out); - return NULL; - } - - if (col1 < 1) { - col1 += input->numCols; - } - - if (row1 < 1) { - row1 += input->numRows; - } - - if ( col0 < 0 || - row0 < 0 || - col1 > input->numCols || - row1 > input->numRows || - col0 >= col1 || - row0 >= row1) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_SUBSET_RANGE_INVALID, - col0, col1, row0, row1, - input->numCols, input->numRows); - psFree(out); - return NULL; - } - - type = input->type.type; - inRows = input->numRows; - inCols = input->numCols; - - if (direction == PS_CUT_X_NEG || direction == PS_CUT_Y_NEG) { - delta = -1; - } - - if (mask != NULL) { - if (inRows != mask->numRows || inCols != mask->numCols) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_IMAGE_MASK_SIZE, - mask->numCols,mask->numRows, - inCols, inRows); - psFree(out); - return NULL; - } - if (mask->type.type != PS_TYPE_MASK) { - char* typeStr; - PS_TYPE_NAME(typeStr,mask->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_MASK_TYPE, - typeStr, PS_TYPE_MASK_NAME); - psFree(out); - return NULL; - } - } - - if (stats == NULL) { - psError(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) { - psError(PS_ERR_BAD_PARAMETER_VALUE, false, - PS_ERRORTEXT_psImage_BAD_STAT,stats->options); - psFree(out); - return NULL; - } - // since stats input is const, I need to - // create a 'scratch' stats struct - myStats = psAlloc(sizeof(psStats)); - *myStats = *stats; - - psS32 numCols = col1-col0; - psS32 numRows = row1-row0; - - if (direction == PS_CUT_X_POS || direction == PS_CUT_X_NEG) { - psVector* imgVec = psVectorAlloc(numRows, type); - psVector* maskVec = NULL; - psMaskType* maskData = NULL; - psU32* outPosition = NULL; - - // recycle output to make a proper sized/type output structure - // n.b. type is double as that is the type given for all stats is - // psStats. - out = psVectorRecycle(out, numCols, PS_TYPE_F64); - if (coords != NULL) { - coords = psVectorRecycle(coords, numCols, PS_TYPE_U32); - outPosition = coords->data.U32; - } - outData = out->data.F64; - if (delta < 0) { - outData += numCols - 1; - if (outPosition != NULL) { - outPosition += numCols - 1; - } - } - - if (mask != NULL) { - maskVec = psVectorAlloc(numRows, mask->type.type); - } - #define PSIMAGE_CUT_VERTICAL(TYPE) \ - case PS_TYPE_##TYPE: { \ - psMaskType* maskVecData = NULL; \ - for (psS32 c=col0;cdata.TYPE[row0] + c; \ - ps##TYPE *imgVecData = imgVec->data.TYPE; \ - if (maskVec != NULL) { \ - maskVecData = maskVec->data.U8; \ - maskData = (psMaskType* )(mask->data.U8[row0]) + c; \ - } \ - for (psS32 r=row0;rtype = input->type; - imgVec->n = *(int*)&imgVec->nalloc = numCols; - if (mask != NULL) { - maskVec = psAlloc(sizeof(psVector)); - maskVec->type = mask->type; - maskVec->n = *(int*)&maskVec->nalloc = numCols; - } - // recycle output to make a proper sized/type output structure - // n.b. type is double as that is the type given for all stats in - // psStats. - out = psVectorRecycle(out, numRows, PS_TYPE_F64); - if (coords != NULL) { - coords = psVectorRecycle(coords, numRows, PS_TYPE_U32); - outPosition = coords->data.U32; - } - outData = out->data.F64; - if (delta < 0) { - outData += numRows-1; - if (outPosition != NULL) { - outPosition += numRows-1; - } - } - - for (psS32 r = row0; r < row1; r++) { - // point the vector struct to the - // data to calculate the stats - imgVec->data.U8 = (psPtr )(input->data.U8[r] + col0 * elementSize); - if (maskVec != NULL) { - maskVec->data.U8 = (psPtr )(mask->data.U8[r] + col0 * sizeof(psMaskType)); - } - myStats = psVectorStats(myStats, imgVec, NULL, maskVec, maskVal); - (void)p_psGetStatValue(myStats, &statVal); // we know it works cause we tested it above - *outData = statVal; - if (outPosition != NULL) { - *outPosition = r; - outPosition += delta; - - } - outData += delta; - } - psFree(imgVec); - psFree(maskVec); - } else { // don't know what the direction flag is - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_SLICE_DIRECTION_INVALID, - direction); - psFree(out); - out = NULL; - } - - psFree(myStats); - - return out; -} - -psVector* psImageCut(psVector* out, - psVector* cutCols, - psVector* cutRows, - const psImage* input, - const psImage* mask, - psMaskType maskVal, - psRegion region, - unsigned int nSamples, - psImageInterpolateMode mode) -{ - - if (input == NULL || input->data.V == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - psFree(out); - return NULL; - } - psS32 numCols = input->numCols; - psS32 numRows = input->numRows; - - if (nSamples < 2) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_nSamples_TOOSMALL, - nSamples); - psFree(out); - return NULL; - } - - float startCol = region.x0; - float startRow = region.y0; - float endCol = region.x1; - float endRow = region.y1; - if (startCol < 0 || startCol >= numCols || - startRow < 0 || startRow >= numRows || - endCol < 0 || endCol >= numCols || - endRow < 0 || endRow >= numRows) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_LINE_NOT_IN_IMAGE, - startCol,startRow,endCol,endRow, - numCols-1,numRows-1); - psFree(out); - return NULL; - } - - if (mode < PS_INTERPOLATE_FLAT ) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImageManip_INTERPOLATION_MODE_UNSUPPORTED, - mode); - psFree(out); - return NULL; - } - - if (mask != NULL) { - if (numRows != mask->numRows || numCols != mask->numCols) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_IMAGE_MASK_SIZE, - mask->numCols,mask->numRows, - numCols-1, numRows); - psFree(out); - return NULL; - } - if (mask->type.type != PS_TYPE_MASK) { - char* typeStr; - PS_TYPE_NAME(typeStr,mask->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_MASK_TYPE, - typeStr, PS_TYPE_MASK_NAME); - psFree(out); - return NULL; - } - } - - //resize the vectors for the coordinate output - psF32* cutColsData = NULL; - psF32* cutRowsData = NULL; - if (cutCols != NULL) { - cutCols = psVectorRecycle(cutCols, nSamples, PS_TYPE_F32); - cutColsData = cutCols->data.F32; - } - if (cutRows != NULL) { - cutRows = psVectorRecycle(cutRows, nSamples, PS_TYPE_F32); - cutRowsData = cutRows->data.F32; - } - - out = psVectorRecycle(out, nSamples, input->type.type); - - float dX = (endCol - startCol) / (float)(nSamples-1); - float dY = (endRow - startRow) / (float)(nSamples-1); - - #define LINEAR_CUT_CASE(TYPE) \ -case PS_TYPE_##TYPE: { \ - ps##TYPE* outData = out->data.TYPE; \ - for (psS32 i = 0; i < nSamples; i++) { \ - float x = startCol + (float)i*dX; \ - float y = startRow + (float)i*dY; \ - /* store off the location of the sample. */ \ - if (cutColsData != NULL) { \ - cutColsData[i] = x; \ - } \ - if (cutRowsData != NULL) { \ - cutRowsData[i] = y; \ - } \ - outData[i] = psImagePixelInterpolate(input,x,y,mask,maskVal,0,mode); \ - } \ - } \ - break; - - - switch (input->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,input->type.type); - psError(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* restrict mask, - psMaskType maskVal, - float x, - float y, - const psVector* radii, - const psStats* stats) -{ - double statVal; - - /* check the parameters */ - - if (input == NULL || input->data.V == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - psFree(out); - return NULL; - } - psS32 numCols = input->numCols; - psS32 numRows = input->numRows; - - if (mask != NULL) { - if (numRows != mask->numRows || numCols != mask->numCols) { - psError(PS_ERR_BAD_PARAMETER_SIZE, true, - PS_ERRORTEXT_psImage_IMAGE_MASK_SIZE, - mask->numCols,mask->numRows, - numCols, numRows); - psFree(out); - return NULL; - } - if (mask->type.type != PS_TYPE_MASK) { - char* typeStr; - PS_TYPE_NAME(typeStr,mask->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_MASK_TYPE, - typeStr, PS_TYPE_MASK_NAME); - psFree(out); - return NULL; - } - } - - if (x < 0 || x >= numCols || - y < 0 || y >= numRows) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_CENTER_NOT_IN_IMAGE, - x, y, - numCols-1, numRows-1); - psFree(out); - return NULL; - } - - if (radii == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_RADII_VECTOR_NULL); - psFree(out); - return NULL; - } - - if (radii->n < 2) { - psError(PS_ERR_BAD_PARAMETER_SIZE, true, - PS_ERRORTEXT_psImage_RADII_VECTOR_TOOSMALL, - radii->n); - psFree(out); - return NULL; - } - - if (stats == NULL) { - psError(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) { - psError(PS_ERR_BAD_PARAMETER_VALUE, false, - PS_ERRORTEXT_psImage_BAD_STAT, - stats->options); - psFree(out); - return NULL; - } - - /* completed checking the parameters */ - - // size the output vector to proper size. - psS32 numOut = radii->n - 1; - out = psVectorRecycle(out, numOut, PS_TYPE_F64); - psF64* outData = out->data.F64; - - psVector* rSqVec = psVectorCopy(NULL, radii, PS_TYPE_F32); - psF32* rSq = rSqVec->data.F32; - - psS32 startRow = y - rSq[numOut]; - psS32 endRow = y + rSq[numOut]; - psS32 startCol = x - rSq[numOut]; - psS32 endCol = x + 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 radii data - for (psS32 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 (psS32 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]), - input->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 (psS32 row=startRow; row <= endRow; row++) { - psF32* inRow = input->data.F32[row]; - psMaskType* maskRow = NULL; - if (mask != NULL) { - maskRow = mask->data.PS_TYPE_MASK_DATA[row]; - } - for (psS32 col=startCol; col <= endCol; col++) { - dX = x - (float)col - 0.5f; - dY = y - (float)row - 0.5f; - dist = dX*dX+dY*dY; - for (psS32 r = 0; r < numOut; r++) { - if (rSq[r] < dist && dist < rSq[r+1]) { - psS32 n = buffer[r]->n; - if (n == buffer[r]->nalloc) { // in case buffers already full, expand - buffer[r] = psVectorRealloc(buffer[r], n*2); - if (bufferMask[r] != NULL) { - bufferMask[r] = psVectorRealloc(bufferMask[r], n*2); - } - } - - 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 (psS32 r = 0; r < numOut; r++) { - myStats = psVectorStats(myStats,buffer[r], NULL, bufferMask[r],maskVal); - (void)p_psGetStatValue(myStats,&statVal); - outData[r] = statVal; - } - - psFree(myStats); - - for (psS32 lcv = 0; lcv < numOut; lcv++) { - psFree(buffer[lcv]); - psFree(bufferMask[lcv]); - } - psFree(buffer); - psFree(bufferMask); - psFree(rSqVec); - return out; -} Index: trunk/psLib/src/image/psImagePixelExtract.h =================================================================== --- trunk/psLib/src/image/psImagePixelExtract.h (revision 4534) +++ (revision ) @@ -1,126 +1,0 @@ -/** @file psImagePixelExtract.h -* -* @brief Contains basic image extraction operations, as specified in the -* PSLIB SDRS sections "Image Pixel Extractions". -* -* @ingroup Image -* -* @author Robert DeSonia, MHPCC -* -* @version $Revision: 1.5 $ $Name: not supported by cvs2svn $ -* @date $Date: 2005-06-25 00:51:28 $ -* -* Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii -*/ - -#ifndef PSIMAGE_PIXEL_EXTRACT_H -#define PSIMAGE_PIXEL_EXTRACT_H - -#include "psImage.h" -#include "psVector.h" -#include "psStats.h" - -/// @addtogroup Image -/// @{ - -/* Cut direction flag. Used with psImageCut function. - */ -typedef enum { - PS_CUT_X_POS, ///< Cut in the x dimension from left to right - PS_CUT_X_NEG, ///< Cut in the x dimension from rigth to left - PS_CUT_Y_POS, ///< Cut in the y dimension from bottom up - PS_CUT_Y_NEG, ///< Cut in the y dimension from top down. -} psImageCutDirection; - -/** Extract pixels from rectlinear region to a vector (array of floats). - * - * The output vector contains either col1-col0 or row1-row0 elements, based - * on the value of the direction: e.g., if direction is PS_CUT_X_POS, there - * are col1-col0 elements. The region to be sliced is defined by the - * lower-left corner, (col0,row0), and the upper-right corner, (col1,row1). - * Note that the row and column of the upper right-hand corner are NOT - * included in the region. In the event that col1 or row1 are negative, they - * shall be interpreted as being relative to the size of the parent image in - * that dimension. The input region is collapsed in the direction perpendicular - * to that specified by direction, and each element of the output vectors is - * derived from the statistics of the pixels at that direction coordinate. The - * statistic used to derive the output vector value is specified by stats. - * If mask is non-NULL, pixels for which the corresponding mask pixel - * matches maskVal are excluded from operations. If coords is not NULL, the - * calculated coordinates along the slice are returned in this vector. Only - * one of the statistics choices may be specified, otherwise the function - * must return an error. - * - * This function is defined for the following types: psS8, psU16, psF32, psF64. - * - * @return psVector the resulting vector - */ -psVector* psImageSlice( - psVector* out, ///< psVector to recycle, or NULL. - psVector* coords, - ///< If not NULL, it is populated with the coordinate in the slice dimension - ///< coorsponding to the output vector's value of the same position in the - ///< vector. This vector maybe resized and retyped as appropriate. - const psImage* input, ///< the input image in which to perform the slice - const psImage* mask, ///< the mask for the input image. - psMaskType maskVal, ///< the mask value to apply to the mask - psRegion region, ///< the slice region - psImageCutDirection direction, ///< the slice dimension and direction - const psStats* stats ///< the statistic to perform in slice operation -); - -/** Extract pixels from an image along a line to a vector (array of floats). - * - * The vector (xs,ys) - (xe,ye) forms the basis of the output vector. Pixels - * are considered in a rectangular region of width dw about this vector. The - * input region is collapsed in the perpendicular direction, and each element - * of the output vector represents pixel-sized boxes, where the value is - * derived from the statistics of the pixels interpolated along the - * perpendicular direction. The specific algorithm which must be used is - * described in the PSLib ADD (PSDC-430-006). The statistic used to derive - * the output vector value is specified by stats. Only one of the statistics - * choices may be specified, otherwise the function must return an error. - * This function must be defined for the following types: psS8, psU16, psF32, - * psF64. - * - * @return psVector* resulting vector - */ -psVector* psImageCut( - psVector* out, ///< psVector to recycle, or NULL. - 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* mask, ///< the mask for the input image. - psMaskType maskVal, ///< the mask value to apply to the mask - psRegion region, ///< the start and end points to cut along - unsigned int nSamples, ///< the number of samples along the cut - psImageInterpolateMode mode ///< the interpolation method to use -); - -/** Extract radial region data to a vector. A vector is constructed where each - * vector elements is derived from the statistics of the pixels which land - * within one of a sequence of radii. The radii are centered on the image - * pixel coordinate x,y, and are defined by the sequence of values in the - * vector radii. The specific algorithm which must be used is described in - * the PSLib ADD (PSDC-430-006). The statistic used to derive the output - * vector value is specified by stats. Only one of the statistics choices - * may be specified, otherwise the function must return an error. This - * function must be defined for the following types: psS8, psU16, psF32, - * psF64. - * - * @return psVector resulting vector - */ -psVector* psImageRadialCut( - psVector* out, ///< psVector to recycle, or NULL. - const psImage* input, ///< the input image in which to perform the cut - const psImage* mask, ///< the mask for the input image. - psMaskType maskVal, ///< the mask value to apply to the mask - float x, ///< the column of the center of the cut circle - float y, ///< the row of the center of the cut circle - const psVector* radii, ///< the radii of the cut circle - const psStats* stats ///< the statistic to perform in operation -); - -/// @} - -#endif // #ifndef PSIMAGE_PIXEL_EXTRACT_H Index: trunk/psLib/src/image/psImagePixelManip.c =================================================================== --- trunk/psLib/src/image/psImagePixelManip.c (revision 4534) +++ (revision ) @@ -1,397 +1,0 @@ -/** @file psImagePixelManip.c - * - * @brief Contains basic image pixel and geometry manipulation operations, as - * specified in the PSLIB SDRS sections "Image Pixel Manipulations" and - * "Image Geometry Manipulations". - * - * @ingroup Image - * - * @author Robert DeSonia, MHPCC - * @author Ross Harman, MHPCC - * - * @version $Revision: 1.5 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-07-11 21:38:19 $ - * - * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii - */ - -#include -#include // for isfinite(), etc. -#include -#include // for memcpy, etc. - -#include "psImagePixelManip.h" - -#include "psError.h" -#include "psImage.h" -#include "psStats.h" -#include "psMemory.h" -#include "psConstants.h" -#include "psImageErrors.h" -#include "psCoord.h" - -int psImageClip(psImage* input, - double min, - double vmin, - double max, - double vmax) -{ - psS32 numClipped = 0; - psU32 numRows; - psU32 numCols; - - if (input == NULL) { - return 0; - } - - if (max < min) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImageManip_MAXMIN, - (double)min,(double)max); - return 0; - } - - numRows = input->numRows; - numCols = input->numCols; - - switch (input->type.type) { - - #define psImageClipCase(type) \ - case PS_TYPE_##type: { \ - if (vmin < PS_MIN_##type || vmin > PS_MAX_##type) { \ - psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ - PS_ERRORTEXT_psImage_PARAMETER_OUTOF_TYPERANGE, \ - "vmin",vmin, PS_TYPE_##type##_NAME, \ - (psF64)PS_MIN_##type,(psF64)PS_MAX_##type); \ - } \ - if (vmax > PS_MAX_##type || vmax < PS_MIN_##type) { \ - psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ - PS_ERRORTEXT_psImage_PARAMETER_OUTOF_TYPERANGE, \ - "vmax",vmax, PS_TYPE_##type##_NAME, \ - (psF64)PS_MIN_##type,(psF64)PS_MAX_##type); \ - } \ - for (psU32 row = 0;rowdata.type[row]; \ - for (psU32 col = 0; col < numCols; col++) { \ - if ((psF64)inputRow[col] < min) { \ - inputRow[col] = (ps##type)vmin; \ - numClipped++; \ - } else if ((psF64)inputRow[col] > max) { \ - inputRow[col] = (ps##type)vmax; \ - numClipped++; \ - } \ - } \ - } \ - } \ - break; - - #define psImageClipCaseComplex(type,absfcn)\ - case PS_TYPE_##type: { \ - if (vmin < PS_MIN_##type || vmin > PS_MAX_##type) { \ - psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ - PS_ERRORTEXT_psImage_PARAMETER_OUTOF_TYPERANGE, \ - "vmin",vmin, PS_TYPE_##type##_NAME, \ - (psF64)PS_MIN_##type,(psF64)PS_MAX_##type); \ - } \ - if (vmax > PS_MAX_##type || vmax < PS_MIN_##type) { \ - psError(PS_ERR_BAD_PARAMETER_VALUE, true, \ - PS_ERRORTEXT_psImage_PARAMETER_OUTOF_TYPERANGE, \ - "vmax",vmax, PS_TYPE_##type##_NAME, \ - (psF64)PS_MIN_##type,(psF64)PS_MAX_##type); \ - } \ - for (psU32 row = 0;rowdata.type[row]; \ - for (psU32 col = 0; col < numCols; col++) { \ - if (absfcn(inputRow[col]) < min) { \ - inputRow[col] = (ps##type)vmin; \ - numClipped++; \ - } else if (absfcn(inputRow[col]) > max) { \ - inputRow[col] = (ps##type)vmax; \ - numClipped++; \ - } \ - } \ - } \ - } \ - break; - - psImageClipCase(S8) - psImageClipCase(S16) - psImageClipCase(S32) // Not a requirement - psImageClipCase(S64) // Not a requirement - psImageClipCase(U8) - psImageClipCase(U16) - psImageClipCase(U32) // Not a requirement - psImageClipCase(U64) // Not a requirement - psImageClipCase(F32) - psImageClipCase(F64) - psImageClipCaseComplex(C32, cabsf) - psImageClipCaseComplex(C64, cabs) - - default: { - char* typeStr; - PS_TYPE_NAME(typeStr,input->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, - typeStr); - } - } - - return numClipped; -} - -int psImageClipNaN(psImage* input, - float value) -{ - psS32 numClipped = 0; - psU32 numRows; - psU32 numCols; - - if (input == NULL) { - return 0; - } - numRows = input->numRows; - numCols = input->numCols; - - switch (input->type.type) { - - #define psImageClipNaNCase(type) \ - case PS_TYPE_##type: \ - for (psU32 row = 0;rowdata.type[row]; \ - for (psU32 col = 0; col < numCols; col++) { \ - if (! isfinite(inputRow[col])) { \ - inputRow[col] = (ps##type)value; \ - numClipped++; \ - } \ - } \ - } \ - break; - - psImageClipNaNCase(F32) - psImageClipNaNCase(F64) - psImageClipNaNCase(C32) - psImageClipNaNCase(C64) - - default: { - char* typeStr; - PS_TYPE_NAME(typeStr,input->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_IMAGE_TYPE_UNSUPPORTED, - typeStr); - } - } - - return numClipped; -} - -int psImageOverlaySection(psImage* image, - const psImage* overlay, - int col0, - int row0, - const char *op) -{ - psU32 imageNumRows; - psU32 imageNumCols; - psU32 overlayNumRows; - psU32 overlayNumCols; - psU32 imageRowLimit; - psU32 imageColLimit; - psElemType type; - psU32 pixelsOverlaid = 0; - - if (image == NULL || overlay == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - return pixelsOverlaid; - } - - if (op == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImageManip_OPERATION_NULL); - return pixelsOverlaid; - } - - type = image->type.type; - - if (type != overlay->type.type) { - char* typeStr; - char* typeStrOverlay; - PS_TYPE_NAME(typeStr,type); - PS_TYPE_NAME(typeStrOverlay,overlay->type.type); - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImageManip_OVERLAY_TYPE_MISMATCH, - typeStrOverlay, typeStr); - return pixelsOverlaid; - } - - imageNumRows = image->numRows; - imageNumCols = image->numCols; - overlayNumRows = overlay->numRows; - overlayNumCols = overlay->numCols; - imageRowLimit = row0 + overlayNumRows; - imageColLimit = col0 + overlayNumCols; - - /* check to see if overlay is within the input image */ - if ( row0 < 0 || - col0 < 0 || - imageRowLimit > imageNumRows || - imageColLimit > imageNumCols) { - - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_SUBSET_RANGE_INVALID, - col0, imageColLimit, row0, imageRowLimit, - imageNumCols, imageNumRows); - return pixelsOverlaid; - } - - - #define psImageOverlayLoop(DATATYPE,OP) { \ - for (int row=row0;rowdata.DATATYPE[row]; \ - ps##DATATYPE* overlayRow = overlay->data.DATATYPE[row-row0]; \ - for (int col=col0;coloptions, NULL); - PS_ASSERT_IMAGE_NON_NULL(in, NULL) - if (mask != NULL) { - PS_ASSERT_IMAGE_TYPE(mask, PS_TYPE_U8, NULL); - PS_ASSERT_IMAGES_SIZE_EQUAL(in, mask, NULL); - } - - if (in->parent == NULL) { - // stuff the image data into a psVector struct. - junkData = (psVector *) psAlloc(sizeof(psVector)); - junkData->type = in->type; - *(int*)&junkData->nalloc = in->numRows * in->numCols; - junkData->n = junkData->nalloc; - junkData->data.U8 = in->data.V[0]; // since psImage data is contiguous... - } else { - // image not necessarily contiguous - int numRows = in->numRows; - int numCols = in->numCols; - int rowSize = numCols * (PSELEMTYPE_SIZEOF(in->type.type)); - - junkData = psVectorAlloc(numRows*numCols, in->type.type); - junkData->n = junkData->nalloc; - - psU8* data = junkData->data.U8; - for (int row = 0; row < numRows; row++) { - memcpy(data, in->data.V[row], rowSize); - data += rowSize; - } - } - - if (mask != NULL) { - if (mask->parent == NULL) { - // stuff the mask data into a psVector struct. - junkMask = psAlloc(sizeof(psVector)); - junkMask->type = mask->type; - *(int*)&junkMask->nalloc = mask->numRows * mask->numCols; - junkMask->n = junkMask->nalloc; - junkMask->data.U8 = mask->data.V[0]; - } else { - // image not necessarily contiguous - int numRows = mask->numRows; - int numCols = mask->numCols; - int rowSize = numCols * (PSELEMTYPE_SIZEOF(mask->type.type)); - - junkMask = psVectorAlloc(numRows*numCols, mask->type.type); - junkMask->n = junkMask->nalloc; - - psU8* data = junkMask->data.U8; - for (int row = 0; row < numRows; row++) { - memcpy(data, mask->data.V[row], rowSize); - data += rowSize; - } - } - } - - stats = psVectorStats(stats, junkData, NULL, junkMask, maskVal); - - psFree(junkMask); - psFree(junkData); - return (stats); -} - -/***************************************************************************** -NOTE: We assume that the psHistogram structure out has already been allocated -and initialized. - *****************************************************************************/ -psHistogram* psImageHistogram(psHistogram* out, - const psImage* in, - const psImage* mask, - psMaskType maskVal) -{ - PS_ASSERT_PTR_NON_NULL(out, NULL); - PS_ASSERT_PTR_NON_NULL(in, NULL); - if (mask != NULL) { - PS_ASSERT_IMAGE_TYPE(mask, PS_TYPE_U8, NULL); - PS_ASSERT_IMAGES_SIZE_EQUAL(in, mask, NULL); - } - psVector* junkData = NULL; - psVector* junkMask = NULL; - - if (in->parent == NULL) { - // stuff the image data into a psVector struct. - junkData = (psVector *) psAlloc(sizeof(psVector)); - junkData->type = in->type; - *(int*)&junkData->nalloc = in->numRows * in->numCols; - junkData->n = junkData->nalloc; - junkData->data.U8 = in->data.V[0]; // since psImage data is contiguous... - } else { - // image not necessarily contiguous - int numRows = in->numRows; - int numCols = in->numCols; - int rowSize = numCols * (PSELEMTYPE_SIZEOF(in->type.type)); - - junkData = psVectorAlloc(numRows*numCols, in->type.type); - junkData->n = junkData->nalloc; - - psU8* data = junkData->data.U8; - for (int row = 0; row < numRows; row++) { - memcpy(data, in->data.V[row], rowSize); - data += rowSize; - } - } - - if (mask != NULL) { - if (mask->parent == NULL) { - // stuff the mask data into a psVector struct. - junkMask = psAlloc(sizeof(psVector)); - junkMask->type = mask->type; - *(int*)&junkMask->nalloc = mask->numRows * mask->numCols; - junkMask->n = junkMask->nalloc; - junkMask->data.U8 = mask->data.V[0]; - } else { - // image not necessarily contiguous - int numRows = mask->numRows; - int numCols = mask->numCols; - int rowSize = numCols * (PSELEMTYPE_SIZEOF(mask->type.type)); - - junkMask = psVectorAlloc(numRows*numCols, mask->type.type); - junkMask->n = junkMask->nalloc; - - psU8* data = junkMask->data.U8; - for (int row = 0; row < numRows; row++) { - memcpy(data, mask->data.V[row], rowSize); - data += rowSize; - } - } - } - - out = psVectorHistogram(out, junkData, NULL, junkMask, maskVal); - - psFree(junkMask); - psFree(junkData); - - return (out); -} - -/***************************************************************************** -calcScaleFactorsEval(n): The Chebyshev polynomials are defined over the -interval [-1.0 : 1.0]. Images typically have sizes of 512x512 or more. In -order to use Chebyshev polynomials, we must scale the coordinates from -0:512 to -1:1. This routine takes as input an integer N and produces as -output a vector of evenly spaced floating point values between -1.0:1.0. - -XXX: Use the p_psNormalizeVector here? - *****************************************************************************/ -double* calcScaleFactors(psS32 n) -{ - PS_ASSERT_INT_NONNEGATIVE(n, NULL); - psS32 i = 0; - double tmp = 0.0; - double *scalingFactors = (double *)psAlloc(n * sizeof(double)); - - for (i = 0; i < n; i++) { - tmp = (double)(n - i); - tmp = (M_PI * (tmp - 0.5)) / ((double)n); - scalingFactors[i] = cos(tmp); - } - - return (scalingFactors); -} - -// XXX: Use a static array of Chebyshev polynomials. -psPolynomial1D **p_psCreateChebyshevPolys(psS32 maxChebyPoly) -{ - PS_ASSERT_INT_POSITIVE(maxChebyPoly, NULL); - psPolynomial1D **chebPolys = NULL; - psS32 i = 0; - psS32 j = 0; - - chebPolys = (psPolynomial1D **) psAlloc(maxChebyPoly * sizeof(psPolynomial1D *)); - for (i = 0; i < maxChebyPoly; i++) { - chebPolys[i] = psPolynomial1DAlloc(i + 1, PS_POLYNOMIAL_ORD); - } - - // Create the Chebyshev polynomials. - // Polynomial i has i-th order. - chebPolys[0]->coeff[0] = 1; - chebPolys[1]->coeff[1] = 1; - for (i = 2; i < maxChebyPoly; i++) { - for (j = 0; j < chebPolys[i - 1]->n; j++) { - chebPolys[i]->coeff[j + 1] = 2 * chebPolys[i - 1]->coeff[j]; - } - for (j = 0; j < chebPolys[i - 2]->n; j++) { - chebPolys[i]->coeff[j] -= chebPolys[i - 2]->coeff[j]; - } - } - - return (chebPolys); -} - -/***************************************************************************** -psImageFitPolynomial(): This routine takes as input a 2-D image and produces -as output the coefficients of the Chebyshev polynomials which match that -input image. - Input: - Output: - Internal Data Structures: - chebPolys[i][j] - sums[i][j]: This will contain the sum of - input->data.F32[x][y] * - psPolynomial1DEval( -chebPolys[i], -(float) x) * - psPolynomial1DEval( -chebPolys[j], -(float) y, -); - over all pixels (x,y) in the image. - *****************************************************************************/ -psPolynomial2D* psImageFitPolynomial(psPolynomial2D* coeffs, - const psImage* input) -{ - PS_ASSERT_IMAGE_NON_NULL(input, NULL); - PS_ASSERT_IMAGE_NON_EMPTY(input, NULL); - if ((input->type.type != PS_TYPE_S8) && - (input->type.type != PS_TYPE_U16) && - (input->type.type != PS_TYPE_F32) && - (input->type.type != PS_TYPE_F64)) { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, "Unallowable image type.\n"); - } - PS_ASSERT_POLY_NON_NULL(coeffs, NULL); - PS_ASSERT_POLY_TYPE(coeffs, PS_POLYNOMIAL_CHEB, NULL); - psS32 x = 0; - psS32 y = 0; - psS32 i = 0; - psS32 j = 0; - double **sums = NULL; - psPolynomial1D* *chebPolys = NULL; - psS32 maxChebyPoly = 0; - double *cScalingFactors = NULL; - double *rScalingFactors = NULL; - - // Create the sums[][] data structure. This - // will hold the LHS of - // equation - // 29 in the ADD: sums[k][l] = SUM { - // image(x,y) * Tk(x) * Tl(y) } - sums = (double **)psAlloc(coeffs->nX * sizeof(double *)); - for (i = 0; i < coeffs->nX; i++) { - sums[i] = (double *)psAlloc(coeffs->nY * sizeof(double)); - } - // We scale the pixel positions to values - // between -1.0 and 1.0 - rScalingFactors = calcScaleFactors(input->numRows); - cScalingFactors = calcScaleFactors(input->numCols); - - // Determine how many Chebyshev polynomials - // are needed, then create them. - maxChebyPoly = coeffs->nX; - if (coeffs->nY > coeffs->nX) { - maxChebyPoly = coeffs->nY; - } - chebPolys = p_psCreateChebyshevPolys(maxChebyPoly); - - // Compute the sums[][] data structure. - for (i = 0; i < coeffs->nX; i++) { - for (j = 0; j < coeffs->nY; j++) { - sums[i][j] = 0.0; - for (x = 0; x < input->numRows; x++) { - for (y = 0; y < input->numCols; y++) { - double pixel = 0.0; - if (input->type.type == PS_TYPE_S8) { - pixel = (double) input->data.S8[x][y]; - } else if (input->type.type == PS_TYPE_U16) { - pixel = (double) input->data.U16[x][y]; - } else if (input->type.type == PS_TYPE_F32) { - pixel = (double) input->data.F32[x][y]; - } else if (input->type.type == PS_TYPE_F64) { - pixel = input->data.F64[x][y]; - } - sums[i][j] += pixel * psPolynomial1DEval(chebPolys[i],rScalingFactors[x]) * - psPolynomial1DEval(chebPolys[j], cScalingFactors[y]); - } - } - } - } - - for (i = 0; i < coeffs->nX; i++) { - for (j = 0; j < coeffs->nY; j++) { - coeffs->coeff[i][j] = sums[i][j]; - coeffs->coeff[i][j] /= (double)(input->numRows * input->numCols); - - if ((i != 0) && (j != 0)) { - coeffs->coeff[i][j] *= 4.0; - } else if ((i == 0) && (j == 0)) { - coeffs->coeff[i][j] *= 1.0; - } else { - coeffs->coeff[i][j] *= 2.0; - } - } - } - - // Free the Chebyshev polynomials that were - // created in this routine. - for (i = 0; i < maxChebyPoly; i++) { - psFree(chebPolys[i]); - } - psFree(chebPolys); - - // Free some data - for (i = 0; i < coeffs->nX; i++) { - psFree(sums[i]); - } - psFree(sums); - psFree(cScalingFactors); - psFree(rScalingFactors); - - return (coeffs); -} - - - - -/***************************************************************************** -psImageFitPolynomial(): This routine takes as input a 2-D image and produces -as output the coefficients of the Chebyshev polynomials which match that input -image. This is a TEST version of the code. It is not used by anything. - Input: - Output: - Internal Data Structures: - chebPolys[i][j] - sums[i][j]: This will contain the sum of - input->data.F32[x][y] * - psPolynomial1DEval( -chebPolys[i], -(float) x) * - psPolynomial1DEval( -chebPolys[j], -(float) y, -); - over all pixels (x,y) in the image. - *****************************************************************************/ -psPolynomial2D* psImageFitPolynomialTest(psPolynomial2D* coeffs, - const psImage* input) -{ - PS_ASSERT_IMAGE_NON_NULL(input, NULL); - PS_ASSERT_IMAGE_NON_EMPTY(input, NULL); - if ((input->type.type != PS_TYPE_S8) && - (input->type.type != PS_TYPE_U16) && - (input->type.type != PS_TYPE_F32) && - (input->type.type != PS_TYPE_F64)) { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, "Unallowable image type.\n"); - } - PS_ASSERT_POLY_NON_NULL(coeffs, NULL); - PS_ASSERT_POLY_TYPE(coeffs, PS_POLYNOMIAL_CHEB, NULL); - psS32 x = 0; - psS32 y = 0; - psS32 i = 0; - psS32 j = 0; - double **sums = NULL; - psPolynomial1D* *chebPolys = NULL; - psS32 maxChebyPoly = 0; - double *cScalingFactors = NULL; - double *rScalingFactors = NULL; - psImage *nodes = psImageAlloc(input->numCols, input->numRows, PS_TYPE_F64); - - double min = -1.0; - double max = 1.0; - double bma = 0.5 * (max-min); // 1 - double bpa = 0.5 * (max+min); // 0 - // We must calculate the value of the image at the nodes where the - // Chebyshev polynomials are 0. - for (x = 0; x < input->numRows; x++) { - double xTmp = cos(M_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 yNode = - ((yTmp + bma + bpa) - 1.0); - double yOrig = ((float) input->numCols) * (yNode - min) / (max - min); - - // nodes->data.F64[x][y] = psImagePixelInterpolate(input, yNode, xNode, NULL, 0, 0.0, PS_INTERPOLATE_BILINEAR); - // nodes->data.F64[x][y] = psImagePixelInterpolate(input, yTmp, xTmp, NULL, 0, 0.0, PS_INTERPOLATE_BILINEAR); - nodes->data.F64[x][y] = psImagePixelInterpolate(input, yOrig, xOrig, NULL, 0, 0.0, PS_INTERPOLATE_BILINEAR); - } - } - - // Create the sums[][] data structure. This - // will hold the LHS of - // equation - // 29 in the ADD: sums[k][l] = SUM { - // image(x,y) * Tk(x) * Tl(y) } - sums = (double **)psAlloc(coeffs->nX * sizeof(double *)); - for (i = 0; i < coeffs->nX; i++) { - sums[i] = (double *)psAlloc(coeffs->nY * sizeof(double)); - } - // We scale the pixel positions to values - // between -1.0 and 1.0 - rScalingFactors = calcScaleFactors(input->numRows); - cScalingFactors = calcScaleFactors(input->numCols); - - // Determine how many Chebyshev polynomials - // are needed, then create them. - maxChebyPoly = coeffs->nX; - if (coeffs->nY > coeffs->nX) { - maxChebyPoly = coeffs->nY; - } - chebPolys = p_psCreateChebyshevPolys(maxChebyPoly); - - // Compute the sums[][] data structure. - for (i = 0; i < coeffs->nX; i++) { - for (j = 0; j < coeffs->nY; j++) { - sums[i][j] = 0.0; - for (x = 0; x < input->numRows; x++) { - for (y = 0; y < input->numCols; y++) { - double pixel; - /* - if (input->type.type == PS_TYPE_S8) { - pixel = (double) input->data.S8[x][y]; - } else if (input->type.type == PS_TYPE_U16) { - pixel = (double) input->data.U16[x][y]; - } else if (input->type.type == PS_TYPE_F32) { - pixel = (double) input->data.F32[x][y]; - } else if (input->type.type == PS_TYPE_F64) { - pixel = input->data.F64[x][y]; - } - */ - pixel = nodes->data.F64[x][y]; - sums[i][j] += pixel * psPolynomial1DEval(chebPolys[i], rScalingFactors[x]) * - psPolynomial1DEval(chebPolys[j], cScalingFactors[y]); - } - } - } - } - - for (i = 0; i < coeffs->nX; i++) { - for (j = 0; j < coeffs->nY; j++) { - coeffs->coeff[i][j] = sums[i][j]; - coeffs->coeff[i][j] /= (double)(input->numRows * input->numCols); - - if ((i != 0) && (j != 0)) { - coeffs->coeff[i][j] *= 4.0; - } else if ((i == 0) && (j == 0)) { - coeffs->coeff[i][j] *= 1.0; - } else { - coeffs->coeff[i][j] *= 2.0; - } - } - } - - // Free the Chebyshev polynomials that were - // created in this routine. - for (i = 0; i < maxChebyPoly; i++) { - psFree(chebPolys[i]); - } - psFree(chebPolys); - - // Free some data - for (i = 0; i < coeffs->nX; i++) { - psFree(sums[i]); - } - psFree(sums); - psFree(cScalingFactors); - psFree(rScalingFactors); - psFree(nodes); - - return (coeffs); -} - -/***************************************************************************** -XXX: Use static variables for Chebyshev polynomials and scaling factors. - *****************************************************************************/ -psImage* p_psImageEvalPolynomialCheb(psImage* input, - const psPolynomial2D* coeffs) -{ - PS_ASSERT_POLY_TYPE(coeffs, PS_POLYNOMIAL_CHEB, NULL); - - psS32 x = 0; - psS32 y = 0; - psS32 i = 0; - psS32 j = 0; - psPolynomial1D* *chebPolys = NULL; - psS32 maxChebyPoly = 0; - double *cScalingFactors = NULL; - double *rScalingFactors = NULL; - float polySum = 0.0; - - // We scale the pixel positions to values between -1.0 and 1.0 - // Use static data structures here. - rScalingFactors = calcScaleFactors(input->numRows); - cScalingFactors = calcScaleFactors(input->numCols); - - // Determine how many Chebyshev polynomials - // are needed, then create them. - maxChebyPoly = coeffs->nX; - if (coeffs->nY > coeffs->nX) { - maxChebyPoly = coeffs->nY; - } - - chebPolys = p_psCreateChebyshevPolys(maxChebyPoly); - - for (x = 0; x < input->numRows; x++) { - for (y = 0; y < input->numCols; y++) { - polySum = 0.0; - for (i = 0; i < coeffs->nX; i++) { - for (j = 0; j < coeffs->nY; j++) { - polySum += - psPolynomial1DEval(chebPolys[i], rScalingFactors[x]) * - psPolynomial1DEval(chebPolys[j], cScalingFactors[y]) * - coeffs->coeff[i][j]; - } - } - - if (input->type.type == PS_TYPE_S8) { - input->data.S8[x][y] = (char) polySum; - } else if (input->type.type == PS_TYPE_U16) { - input->data.U16[x][y] = (short int) polySum; - } else if (input->type.type == PS_TYPE_F32) { - input->data.F32[x][y] = (float) polySum; - } else if (input->type.type == PS_TYPE_F64) { - input->data.F64[x][y] = polySum; - } - } - } - - // Free the Chebyshev polynomials that were - // created in this routine. - // XXX: Use static data structures here. - for (i = 0; i < maxChebyPoly; i++) { - psFree(chebPolys[i]); - } - psFree(chebPolys); - - psFree(cScalingFactors); - psFree(rScalingFactors); - - return input; -} - -psImage* p_psImageEvalPolynomialOrd(psImage* input, - const psPolynomial2D* coeffs) -{ - PS_ASSERT_POLY_TYPE(coeffs, PS_POLYNOMIAL_ORD, NULL); - - for (int row = 0; row < input->numRows ; row++) { - for (int col = 0; col < input->numCols ; col++) { - if (input->type.type == PS_TYPE_S8) { - input->data.S8[row][col] = (psS8) psPolynomial2DEval(coeffs, (psF32) row, (psF32) col); - } else if (input->type.type == PS_TYPE_U16) { - input->data.U16[row][col] = (psS16) psPolynomial2DEval(coeffs, (psF32) row, (psF32) col); - } else if (input->type.type == PS_TYPE_F32) { - input->data.F32[row][col] = psPolynomial2DEval(coeffs, (psF32) row, (psF32) col); - } else if (input->type.type == PS_TYPE_F64) { - input->data.F64[row][col] = (psF64) psPolynomial2DEval(coeffs, (psF32) row, (psF32) col); - } - } - } - - return(input); -} - - -/***************************************************************************** -XXX: I added normal polynomials to this routine. Let IfA know, put it in the -psLib SDR. - *****************************************************************************/ -psImage* psImageEvalPolynomial(psImage* input, - const psPolynomial2D* coeffs) -{ - PS_ASSERT_IMAGE_NON_NULL(input, NULL); - PS_ASSERT_IMAGE_NON_EMPTY(input, NULL); - if ((input->type.type != PS_TYPE_S8) && - (input->type.type != PS_TYPE_U16) && - (input->type.type != PS_TYPE_F32) && - (input->type.type != PS_TYPE_F64)) { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - "Unallowable image type.\n"); - } - PS_ASSERT_POLY_NON_NULL(coeffs, NULL); - - if (coeffs->type == PS_POLYNOMIAL_ORD) { - return(p_psImageEvalPolynomialOrd(input, coeffs)); - } else if (coeffs->type == PS_POLYNOMIAL_CHEB) { - return(p_psImageEvalPolynomialCheb(input, coeffs)); - } - printf("XXX: Error: wrong polynomial type\n"); - // XXX: psError() - return(NULL); -} - Index: trunk/psLib/src/image/psImageStats.h =================================================================== --- trunk/psLib/src/image/psImageStats.h (revision 4534) +++ (revision ) @@ -1,89 +1,0 @@ -/** @file psImageStats.h -* \brief Routines for calculating statistics on images. -* @ingroup ImageStats -* -* This file will hold the prototypes for procedures which calculate -* statistic on images, histograms on images, and fit/evaluate Chebyshev -* polynomials to images. -* -* @author GLG, MHPCC -* -* @version $Revision: 1.24 $ $Name: not supported by cvs2svn $ -* @date $Date: 2005-06-25 00:51:28 $ -* -* Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii -*/ -#ifndef PS_IMAGE_STATS_H -#define PS_IMAGE_STATS_H - -#include "psType.h" -#include "psVector.h" -#include "psImage.h" -#include "psStats.h" -#include "psFunctions.h" - -/// @addtogroup ImageStats -/// @{ - -/** This routine must determine the various statistics for the image. - * - * Determine statistics for image (or subimage). The statistics to be - * determined are specified by stats. The mask allows pixels to be excluded - * if their corresponding mask pixel value matches the value of maskVal. - * This function must be defined for the following types: psS8, psU16, psF32, - * psF64. - * - * @return psStats* the resulting statistics result(s) - */ -psStats* psImageStats( - psStats* stats, ///< defines statistics to be calculated - const psImage* in, ///< image (or subimage) to calculate stats - const psImage* mask, ///< mask data for image (NULL ok) - psMaskType maskVal ///< mask value for mask -); - -/** Construct a histogram from an image (or subimage). - * - * The histogram to generate is specified by psHistogram hist (see section - * 4.3.2 in SDRS). This function must be defined for the following types: - * psS8, psU16, psF32, psF64. - * - * @return psHistogram* the resulting histogram - */ -psHistogram* psImageHistogram( - psHistogram* out, ///< input histogram description & target - const psImage* in, ///< Image data to be histogramed. - const psImage* mask, ///< mask data for image (NULL ok) - psMaskType maskVal ///< mask Mask for mask -); - -/** Fit a 2-D polynomial surface to an image. - * - * The input structure coeffs contains the desired order and terms of - * interest. This function must be defined for the following types: psS8, - * psU16, psF32, psF64. - * - * @return psPolynomial2D* fitted polynomial result - * - */ -psPolynomial2D* psImageFitPolynomial( - psPolynomial2D* coeffs, ///< coefficient structure carries in desired terms & target - const psImage* input ///< input image -); - -/** Evaluate a 2-D polynomial surface for the image pixels. - * - * Given the input polynomial coefficients, set the image pixel values on the - * basis of the polynomial function. This function must be defined for the - * following types: psS8, psU16, psF32, psF64. - * - * @return psImage* the resulting image - */ -psImage* psImageEvalPolynomial( - psImage* input, ///< input image - const psPolynomial2D* coeffs ///< coefficient structure carries in desired terms -); - -/// @} - -#endif // #ifndef PS_IMAGE_STATS_H Index: trunk/psLib/src/image/psImageStructManip.c =================================================================== --- trunk/psLib/src/image/psImageStructManip.c (revision 4534) +++ (revision ) @@ -1,379 +1,0 @@ -/** @file psImageStructManip.c - * - * @brief Contains basic image structure manipulations, as specified in the - * PSLIB SDRS sections "Image Structure Manipulation". - * - * @ingroup Image - * - * @author Robert DeSonia, MHPCC - * - * @version $Revision: 1.4 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-07-09 01:24:30 $ - * - * Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii - * - */ - -#include - -#include "psMemory.h" -#include "psImageStructManip.h" -#include "psError.h" - -#include "psImageErrors.h" - -static psImage* imageSubset(psImage* out, - psImage* image, - psS32 col0, - psS32 row0, - psS32 col1, - psS32 row1) -{ - psU32 elementSize; // size of image element in bytes - psS32 inputColOffset; // offset in bytes to first subset pixel in input row - - if ( col0 < 0 || row0 < 0 ) { - // psError(PS_ERR_BAD_PARAMETER_VALUE, true, - // PS_ERRORTEXT_psImage_SUBSET_RANGE_INVALID); - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_SUBSET_RANGE_INVALID, - col0, col1-1, row0, row1-1, - image->numCols-1, image->numRows-1); - return NULL; - } - - if (image == NULL || image->data.V == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - return NULL; - } - - if (image->type.dimen != PS_DIMEN_IMAGE) { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_NOT_AN_IMAGE); - return NULL; - } - - if (col1 < 1) { - col1 = image->numCols + col1; - } - if (row1 < 1) { - row1 = image->numRows + row1; - } - - if ( col1 <= col0 || - row1 <= row0 || - col0 >= image->numCols || - row0 >= image->numRows || - col1 > image->numCols || - row1 > image->numRows ) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_SUBSET_RANGE_INVALID, - col0, col1-1, row0, row1-1, - image->numCols-1, image->numRows-1); - return NULL; - } - - - - psS32 numRows = row1-row0; - psS32 numCols = col1-col0; - - elementSize = PSELEMTYPE_SIZEOF(image->type.type); - - if (image->parent != NULL) { // if this is a child, we need to start working with parent. - col0 += image->col0; - col1 += image->col0; - row0 += image->row0; - row1 += image->row0; - image = (psImage*)image->parent; - } - - // increment the raw data buffer before freeing anything in the 'out' - psPtr rawData = psMemIncrRefCounter(image->rawDataBuffer); - - if (out != NULL) { - // if a child, need to orphan (disassociate from parent) first - if (out->parent != NULL) { - psArrayRemove(out->parent->children,out); // remove from parent's knowledge - out->parent = NULL; // break link to parent - } - - psFree(out->rawDataBuffer); // free the previous data reference - } else { - out = psAlloc(sizeof(psImage)); - out->data.V = NULL; - } - - out->data.V = psRealloc(out->data.V,sizeof(psPtr)*numRows); // resize row pointer array - *(psMathType*)&out->type = image->type; - // *(psU32*)&out->numCols = numCols; - // *(psU32*)&out->numRows = numRows; - *(psS32*)&out->numCols = numCols; - *(psS32*)&out->numRows = numRows; - *(psS32*)&out->row0 = row0; - *(psS32*)&out->col0 = col0; - out->parent = image; - out->children = NULL; - out->rawDataBuffer = rawData; - - // set the new psImage's deallocator to the same as the input image - psMemSetDeallocator(out,psMemGetDeallocator(image)); - - inputColOffset = elementSize * col0; - for (psS32 row = 0; row < numRows; row++) { - out->data.V[row] = image->data.U8[row0 + row] + inputColOffset; - } - - // add output image as a child of the input image. - psS32 n = 0; - psArray* children = image->children; - if (children == NULL) { - children = psArrayAlloc(16); // start with a reasonable size for growth - } else if (children->nalloc == children->n) { // full? - n = children->n; - children = psArrayRealloc(children,n*2); // double the array size - } else { - n = children->n; - } - children->data[n] = out; - children->n = n+1; - image->children = children; // push back any change (esp. if children==NULL before) - - return (out); -} - -psImage* psImageSubset(psImage* image, - psRegion region) -{ - return imageSubset(NULL,image,region.x0, region.y0, - region.x1, region.y1); -} - -psImage* psImageCopy(psImage* output, - const psImage* input, - psElemType type) -{ - psElemType inDatatype; - psS32 elementSize; - psS32 elements; - psS32 numRows; - psS32 numCols; - - if (input == NULL || input->data.V == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - psFree(output); - return NULL; - } - - if (input == output) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_INPLACE_NOTSUPPORTED); - psFree(output); - return NULL; - } - - if (input->type.dimen != PS_DIMEN_IMAGE) { - psError(PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_psImage_NOT_AN_IMAGE); - psFree(output); - return NULL; - } - - inDatatype = input->type.type; - numRows = input->numRows; - numCols = input->numCols; - elements = numRows * numCols; - elementSize = PSELEMTYPE_SIZEOF(inDatatype); - - output = psImageRecycle(output, numCols, numRows, type); - - // cover the trival case of copy of the same - // datatype. - if (type == inDatatype) { - for (psS32 row=0;rowdata.V[row], input->data.V[row], elementSize * numCols); - } - return output; - } - - #define PSIMAGE_ELEMENT_COPY(IN,INTYPE,OUT,OUTTYPE,ELEMENTS) { \ - ps##INTYPE *in; \ - ps##OUTTYPE *out; \ - for(psS32 row=0;rowdata.INTYPE[row]; \ - out = OUT->data.OUTTYPE[row]; \ - for (psS32 col=0;coldata.V == NULL) { - psError(PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_psImage_IMAGE_NULL); - return NULL; - } - - if (image->parent != NULL) { - return imageSubset(image, - (psImage*)image->parent, - region.x0+image->col0, - region.y0+image->row0, - region.x1+image->col0, - region.y1+image->row0); - } - - int col0 = region.x0; - int row0 = region.y0; - int col1 = region.x1; - int row1 = region.y1; - - if (col1 < 1) { - col1 += image->numCols; - } - - if (row1 < 1) { - row1 += image->numRows; - } - - if ( col0 < 0 || - row0 < 0 || - col1 > image->numCols || - row1 > image->numRows || - col0 >= col1 || - row0 >= row1 ) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - PS_ERRORTEXT_psImage_SUBSET_RANGE_INVALID, - col0, col1-1, row0, row1-1, - image->numCols-1, image->numRows-1); - psFree(image); - return NULL; - } - - psImageFreeChildren(image); - - psU32 elementSize = PSELEMTYPE_SIZEOF(image->type.type); - psU32 numCols = col1-col0; - psU32 numRows = row1-row0; - psU32 rowSize = elementSize*numCols; - psU32 colOffset = elementSize * col0; - psU8* imageData = image->rawDataBuffer; - for (psS32 row = row0; row < row1; row++) { - memmove(imageData,image->data.U8[row] + colOffset,rowSize); - imageData += rowSize; - } - - *(psU32*)&image->numRows = numRows; - *(psU32*)&image->numCols = numCols; - - // XXX: should I really resize the buffers? - image->data.V = psRealloc(image->data.V,sizeof(psPtr)*numRows); - image->rawDataBuffer = psRealloc(image->rawDataBuffer,rowSize*numRows); - - image->data.V[0] = image->rawDataBuffer; - for (psS32 r = 1; r < numRows; r++) { - image->data.U8[r] = image->data.U8[r-1] + rowSize; - } - - return (image); -} - Index: trunk/psLib/src/image/psImageStructManip.h =================================================================== --- trunk/psLib/src/image/psImageStructManip.h (revision 4534) +++ (revision ) @@ -1,84 +1,0 @@ -/** @file psImageStructManip.h -* -* @brief Contains basic image structure manipulation operations, as specified -* in the PSLIB SDRS sections "Image Structure Manipulation". -* -* @ingroup Image -* -* @author Robert DeSonia, MHPCC -* -* @version $Revision: 1.3 $ $Name: not supported by cvs2svn $ -* @date $Date: 2005-06-08 23:40:45 $ -* -* Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii -*/ - -#ifndef PSIMAGE_STRUCT_MANIP_H -#define PSIMAGE_STRUCT_MANIP_H - -#include "psImage.h" - -/// @addtogroup Image -/// @{ - -/** Create a subimage of the specified area. - * - * Define a subimage of the specified area of the given image. This function - * must raise an error if the requested subset area lies outside of the - * parent image and return NULL. The argument image is the parent image, - * region.x0, region.y0 specify the starting pixel of the subraster, and - * region.x1,region.y1 specify the extent of the desired subraster. Note - * that the row and column of this “upper right-hand corner” are NOT included - * in the region. In the event that x1 or y1 are negative, they shall be - * interpreted as being relative to the size of the parent image in that - * dimension. The entire subraster must be contained within the raster of the - * parent image. Note that the refCounter for the parent should be - * incremented. This function must be defined for the following types: psU8, - * psU16, psS8, psS16, psF32, psF64, psC32, psC64. - * - * @return psImage* : Pointer to psImage. - * - */ -psImage* psImageSubset( - psImage* image, ///< Parent image. - psRegion region ///< region of subimage -); - -/** Makes a copy of a psImage - * - * @return psImage* Copy of the input psImage. This may not be equal to the - * output parameter - * - */ -psImage* psImageCopy( - psImage* output, ///< if not NULL, a psImage that could be recycled. - const psImage* input, ///< the psImage to copy - psElemType type ///< the desired datatype of the returned copy -); - -/** Trim an image - * - * Trim the specified image in-place, which involves shuffling the pixels - * around in memory. The pixels in the region [col0:col1,row0:row1] shall consist - * the output image. The column col1 and row row1 are NOT included in the range. - * In the event that x1 or y1 are non-positive, they shall be interpreted as - * being relative to the size of the parent image in that dimension. - * - * If the entire specified subimage is not contained within the parent - * image, an error results and the return value will be NULL. - * - * N.B. If the input psImage is a child of another psImage, no pixel data - * will be trimmed, rather it equivalent to calling psImageSubset. If the input - * psImage is, however, a parent psImage, any children will be obliterated, - * i.e., freed from memory. - * - * @return psImage* trimmed image result - */ -psImage* psImageTrim( - psImage* image, ///< image to trim - psRegion region ///< trim region -); - -/// @} - -#endif // #ifndef PSIMAGE_STRUCT_MANIP_H