Index: trunk/psModules/src/detrend/pmFlatField.c =================================================================== --- trunk/psModules/src/detrend/pmFlatField.c (revision 5435) +++ trunk/psModules/src/detrend/pmFlatField.c (revision 5516) @@ -18,6 +18,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.2 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-10-20 23:06:24 $ + * @version $Revision: 1.3 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-11-15 20:09:03 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -36,4 +36,5 @@ #include "pmMaskBadPixels.h" #include "pmFlatFieldErrors.h" +#include "pmSubtractBias.h" @@ -42,6 +43,20 @@ const pmReadout *flat) { + // XXX: Use the proper image and readout asserts. + PS_ASSERT_PTR_NON_NULL(in, true); + PS_ASSERT_PTR_NON_NULL(in->image, false); + PS_ASSERT_PTR_NON_NULL(in->mask, false); + PS_ASSERT_PTR_NON_NULL(flat, false); + PS_ASSERT_PTR_NON_NULL(flat->image, false); + if (in == NULL) + printf("XXX: NULL\n"); + // XXX: Not sure if this is correct. Must consult with IfA. PS_ASSERT_PTR_NON_NULL(in->mask, false); + + PS_WARN_PTR_NON_NULL(in->parent); + if (in->parent != NULL) { + PS_WARN_PTR_NON_NULL(in->parent->concepts); + } int i = 0; int j = 0; @@ -51,37 +66,25 @@ psElemType flatType; psElemType maskType; - psImage *inImage = NULL; psImage *inMask = NULL; psImage *flatImage = NULL; - - // Check for nulls - if (in == NULL) { - return true; // Readout might have data in it - } else if(flat==NULL) { - psError( PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_pmFlatField_NULL_FLAT_READOUT); - return false; - } - - inImage = in->image; + // + // Determine trimmed image from metadata. + // + psImage *trimmedImg = p_psDetermineTrimmedImage(in); flatImage = flat->image; - if (inImage == NULL) { - psError( PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_pmFlatField_NULL_INPUT_IMAGE); - return false; - } else if(flatImage == NULL) { - psError( PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_pmFlatField_NULL_FLAT_IMAGE); - return false; - } inMask = in->mask; // Check input image and its mask are not larger than flat image - if (inImage->numRows>flatImage->numRows || inImage->numCols>flatImage->numCols) { + if (trimmedImg == NULL) + printf("XXX: 00\n"); + if (flatImage == NULL) + printf("XXX 01\n"); + + if (trimmedImg->numRows>flatImage->numRows || trimmedImg->numCols>flatImage->numCols) { psError( PS_ERR_BAD_PARAMETER_SIZE, true, PS_ERRORTEXT_pmFlatField_SIZE_INPUT_IMAGE, - inImage->numRows, inImage->numCols, flatImage->numRows, flatImage->numCols); + trimmedImg->numRows, trimmedImg->numCols, flatImage->numRows, flatImage->numCols); return false; } @@ -94,6 +97,6 @@ // Determine total offset based on image offset with chip offset - totOffCol = inImage->col0 + in->col0; - totOffRow = inImage->row0 + in->row0; + totOffCol = trimmedImg->col0 + in->col0; + totOffRow = trimmedImg->row0 + in->row0; // Check that offsets are within image limits @@ -103,8 +106,8 @@ totOffRow, totOffCol, flatImage->numRows, flatImage->numCols); return false; - } else if(totOffRow>=inImage->numRows || totOffCol>=inImage->numCols) { + } else if(totOffRow>=trimmedImg->numRows || totOffCol>=trimmedImg->numCols) { psError( PS_ERR_BAD_PARAMETER_SIZE, true, PS_ERRORTEXT_pmFlatField_OFFSET_INPUT_IMAGE, - totOffRow, totOffCol, inImage->numRows, inImage->numCols); + totOffRow, totOffCol, trimmedImg->numRows, trimmedImg->numCols); return false; } else if(totOffRow>=inMask->numRows || totOffCol>=inMask->numCols) { @@ -116,5 +119,5 @@ // Check for incorrect types - inType = inImage->type.type; + inType = trimmedImg->type.type; flatType = flatImage->type.type; maskType = inMask->type.type; @@ -145,6 +148,6 @@ case PS_TYPE_##TYPE: \ /* Per Eugene's request, use two sets of loops: first to fill mask, second to avoid div with bad pix */ \ - for(j = totOffRow; j < inImage->numRows; j++) { \ - for(i = totOffCol; i < inImage->numCols; i++) { \ + for(j = totOffRow; j < trimmedImg->numRows; j++) { \ + for(i = totOffCol; i < trimmedImg->numCols; i++) { \ if(flatImage->data.TYPE[j][i] <= 0.0) { \ /* Negative or zero flat pixels shall be masked in input image as PM_MASK_FLAT */ \ @@ -154,9 +157,9 @@ } \ } \ - for(j = totOffRow; j < inImage->numRows; j++) { \ - for(i = totOffCol; i < inImage->numCols; i++) { \ + for(j = totOffRow; j < trimmedImg->numRows; j++) { \ + for(i = totOffCol; i < trimmedImg->numCols; i++) { \ if(!inMask->data.PS_TYPE_MASK_DATA[j][i]) { \ /* Module shall divide the input image by the flat-fielded image */ \ - inImage->data.TYPE[j][i] /= flatImage->data.TYPE[j][i]; \ + trimmedImg->data.TYPE[j][i] /= flatImage->data.TYPE[j][i]; \ } \ } \ Index: trunk/psModules/src/detrend/pmMaskBadPixels.c =================================================================== --- trunk/psModules/src/detrend/pmMaskBadPixels.c (revision 5435) +++ trunk/psModules/src/detrend/pmMaskBadPixels.c (revision 5516) @@ -19,6 +19,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-09-28 20:43:52 $ + * @version $Revision: 1.2 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-11-15 20:09:03 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -35,8 +35,88 @@ #include "pmMaskBadPixels.h" #include "pmMaskBadPixelsErrors.h" +#include "pmSubtractBias.h" + +/****************************************************************************** +GrowPixel(inMask, row, col, radius, maskVal): This private routine takes an +input image mask and a pixel location, then sets (logical or) all pixels with +parameter radius if that pixel to maskVal parameter. + *****************************************************************************/ +psBool GrowPixel( + psImage *inMask, + psS32 col, + psS32 row, + psS32 radius, + psU32 maskVal) +{ + psS32 rowMin = PS_MAX(row-radius, 0); + psS32 rowMax = PS_MIN(row+radius+1, inMask->numRows); + psS32 colMin = PS_MAX(col-radius, 0); + psS32 colMax = PS_MIN(col+radius+1, inMask->numCols); + psF32 squareRadius = PS_SQR(radius); + + + for(psS32 i=rowMin; idata.U8[i][j] |= maskVal; + } + } + } + return(true); +} + +/****************************************************************************** +GetRadius(inImg, col, row, sat, growVal, grow): This private routine takes an +input image and pixel location and determines what radius that pixel must grow +by. + +//XXX: Inline this or macro it for speed. + *****************************************************************************/ +static psS32 GetRadius( + psImage *inImg, + psS32 col, + psS32 row, + psF32 sat, + psU32 growVal, + psS32 grow) +{ + psS32 growRadius = 0; + if (inImg->type.type == PS_TYPE_F32) { + if(inImg->data.F32[row][col] == growVal) { + growRadius = grow; + } + if (inImg->data.F32[row][col] > sat) { + growRadius = PS_MAX(growRadius+1, 2); + } + } else if (inImg->type.type == PS_TYPE_S32) { + if(inImg->data.S32[row][col] == growVal) { + growRadius = grow; + } + if (inImg->data.S32[row][col] > sat) { + growRadius = PS_MAX(growRadius+1, 2); + } + } else if (inImg->type.type == PS_TYPE_U16) { + if(inImg->data.U16[row][col] == growVal) { + growRadius = grow; + } + if (inImg->data.U16[row][col] > sat) { + growRadius = PS_MAX(growRadius+1, 2); + } + } else { + psError( PS_ERR_BAD_PARAMETER_TYPE, true, + PS_ERRORTEXT_pmMaskBadPixels_TYPE_MASK_IMAGE, + inImg->type.type); + } + return(growRadius); +} + bool pmMaskBadPixels( pmReadout *in, - const psImage *mask, + const pmReadout *mask, unsigned int maskVal, float sat, @@ -44,69 +124,57 @@ int grow) { + // XXX: Review this code then put all asserts at the top. + PS_ASSERT_PTR_NON_NULL(in, true); + PS_ASSERT_PTR_NON_NULL(in->image, false); + PS_WARN_PTR_NON_NULL(in->parent); + if (in->parent != NULL) { + PS_WARN_PTR_NON_NULL(in->parent->concepts); + } + PS_ASSERT_PTR_NON_NULL(mask, false); int i = 0; int j = 0; - int jj = 0; - int ii = 0; - int rRound = 0; - int rowMin = 0; - int rowMax = 0; - int colMin = 0; - int colMax = 0; int totOffCol = 0; int totOffRow = 0; - float r = 0.0f; psElemType inType; psElemType maskType; - psImage *inImage = NULL; psImage *inMask = NULL; - - // Check for nulls - if (in == NULL) { - return true; // Readout may not have data in it - } else if(mask==NULL) { - psError( PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_pmMaskBadPixels_NULL_MASK_IMAGE); - return false; - } - - inImage = in->image; - if (inImage == NULL) { - psError( PS_ERR_BAD_PARAMETER_NULL, true, - PS_ERRORTEXT_pmMaskBadPixels_NULL_INPUT_IMAGE); - return false; - } else if(in->mask == NULL) { - in->mask = psImageAlloc(inImage->numCols, inImage->numRows, PS_TYPE_MASK); - memset(in->mask->data.V[0], 0, inImage->numCols*inImage->numRows*PSELEMTYPE_SIZEOF(PS_TYPE_MASK)); + // + // Determine trimmed image from metadata. + // + psImage *trimmedImg = p_psDetermineTrimmedImage(in); + if(in->mask == NULL) { + in->mask = psImageAlloc(trimmedImg->numCols, trimmedImg->numRows, PS_TYPE_MASK); + PS_IMAGE_SET_U8(in->mask, 0); } inMask = in->mask; // Check input image and its mask are not larger than mask - if(inImage->numRows > mask->numRows || inImage->numCols > mask->numCols) { + if(trimmedImg->numRows > mask->image->numRows || trimmedImg->numCols > mask->image->numCols) { psError( PS_ERR_BAD_PARAMETER_SIZE, true, PS_ERRORTEXT_pmMaskBadPixels_SIZE_INPUT_IMAGE, - inImage->numRows, inImage->numCols, mask->numRows, mask->numCols); - return false; - } else if(inMask->numRows>mask->numRows || inMask->numCols>mask->numCols) { + trimmedImg->numRows, trimmedImg->numCols, mask->image->numRows, mask->image->numCols); + return false; + } else if(inMask->numRows > mask->image->numRows || inMask->numCols > mask->image->numCols) { psError( PS_ERR_BAD_PARAMETER_SIZE, true, PS_ERRORTEXT_pmMaskBadPixels_SIZE_MASK_IMAGE, - inMask->numRows, inMask->numCols, mask->numRows, mask->numCols); + inMask->numRows, inMask->numCols, mask->image->numRows, mask->image->numCols); return false; } // Determine total offset based on image offset with chip offset - totOffCol = inImage->col0 + in->col0; - totOffRow = inImage->row0 + in->row0; + totOffCol = trimmedImg->col0 + in->col0; + totOffRow = trimmedImg->row0 + in->row0; // Check that offsets are within image limits - if(totOffRow>=mask->numRows || totOffCol>=mask->numCols) { + if(totOffRow>=mask->image->numRows || totOffCol>=mask->image->numCols) { psError( PS_ERR_BAD_PARAMETER_SIZE, true, PS_ERRORTEXT_pmMaskBadPixels_OFFSET_MASK_IMAGE, - totOffRow, totOffCol, mask->numRows, mask->numCols); - return false; - } else if(totOffRow>=inImage->numRows || totOffCol>=inImage->numCols) { + totOffRow, totOffCol, mask->image->numRows, mask->image->numCols); + return false; + } else if(totOffRow>=trimmedImg->numRows || totOffCol>=trimmedImg->numCols) { psError( PS_ERR_BAD_PARAMETER_SIZE, true, PS_ERRORTEXT_pmMaskBadPixels_OFFSET_INPUT_IMAGE, - totOffRow, totOffCol, inImage->numRows, inImage->numCols); + totOffRow, totOffCol, trimmedImg->numRows, trimmedImg->numCols); return false; } else if(totOffRow>=inMask->numRows || totOffCol>=inMask->numCols) { @@ -118,6 +186,6 @@ // Check for incorrect types - inType = inImage->type.type; - maskType = mask->type.type; + inType = trimmedImg->type.type; + maskType = mask->image->type.type; if(PS_IS_PSELEMTYPE_COMPLEX(inType)) { psError( PS_ERR_BAD_PARAMETER_TYPE, true, @@ -132,57 +200,34 @@ } - // Macro for all PS types - #define PM_BAD_PIXELS(TYPE) \ -case PS_TYPE_##TYPE: \ - for(j=totOffRow; jnumRows; j++) { \ - for(i=totOffCol; inumCols; i++) { \ - \ - /* Pixels with flux greater than sat shall be masked */ \ - if(inImage->data.TYPE[j][i] > sat) { \ - inMask->data.PS_TYPE_MASK_DATA[j][i] |= PM_MASK_SAT; \ - } \ - \ - /* Pixels which satisfy maskVal shall be masked */ \ - inMask->data.PS_TYPE_MASK_DATA[j][i] |= (mask->data.PS_TYPE_MASK_DATA[j][i]&maskVal); \ - \ - /* Pixels which satisfy growVal and within the grow radius shall be masked */ \ - if(mask->data.PS_TYPE_MASK_DATA[j][i] & growVal) { \ - rowMin = MAX(j-grow, 0); \ - rowMax = MIN(j+grow+1, inImage->numRows); \ - colMin = MAX(i-grow, 0); \ - colMax = MIN(i+grow+1, inImage->numCols); \ - for(jj=rowMin; jjdata.PS_TYPE_MASK_DATA[jj][ii] |= \ - (mask->data.PS_TYPE_MASK_DATA[j][i]&growVal); \ - } \ - } \ - } \ - } \ - } \ - } \ - break; - - // Switch to call bad pixel masking macro defined above - switch(inType) { - PM_BAD_PIXELS(U8); - PM_BAD_PIXELS(U16); - PM_BAD_PIXELS(U32); - PM_BAD_PIXELS(U64); - PM_BAD_PIXELS(S8); - PM_BAD_PIXELS(S16); - PM_BAD_PIXELS(S32); - PM_BAD_PIXELS(S64); - PM_BAD_PIXELS(F32); - PM_BAD_PIXELS(F64); - default: - psError( PS_ERR_BAD_PARAMETER_TYPE, true, - PS_ERRORTEXT_pmMaskBadPixels_TYPE_UNSUPPORTED, - inType); - } - - return false; + // + // This is the main loop which examines each pixel and masks pixels if + // A: The mask matches the maskValue + // B: The pixel is larger than the saturation value + // C: The pixel equals the grow value (in which case a circle is masked) + // + for(i=totOffRow; inumRows; i++) { + for(j=totOffCol; jnumCols; j++) { + // + // A: Pixels which satisfy maskVal shall be masked + // + if (mask->image->data.U8[i][j] == maskVal) { + in->mask->data.U8[i][j] |= maskVal; + } + + // + // We first determine how much we need to grow by and store this in + // growRadius. + // + psS32 growRadius = GetRadius(trimmedImg, j, i, sat, growVal, grow); + + // + // Grow the pixel mask if necessary. + // + if (growRadius != 0) { + GrowPixel(in->mask, j, i, growRadius, maskVal); + } + } + } + + return true; } Index: trunk/psModules/src/detrend/pmMaskBadPixels.h =================================================================== --- trunk/psModules/src/detrend/pmMaskBadPixels.h (revision 5435) +++ trunk/psModules/src/detrend/pmMaskBadPixels.h (revision 5516) @@ -19,6 +19,6 @@ * @author Ross Harman, MHPCC * - * @version $Revision: 1.1 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-09-28 20:43:52 $ + * @version $Revision: 1.2 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-11-15 20:09:03 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -36,4 +36,5 @@ } pmMaskValue; +// XXX: Use PS_MIN, PS_MAX /** Macro to find maximum of two numbers */ #define MAX(A,B)((A)>=(B)?(A):(B)) @@ -53,5 +54,5 @@ bool pmMaskBadPixels( pmReadout *in, ///< Readout containing input image data. - const psImage *mask, ///< Mask data to be added to readout mask data. + const pmReadout *mask, ///< Mask data to be added to readout mask data. unsigned int maskVal, ///< Mask value to determine what to add to input mask. float sat, ///< Saturation limit to mask bad pixels. Index: trunk/psModules/src/detrend/pmNonLinear.c =================================================================== --- trunk/psModules/src/detrend/pmNonLinear.c (revision 5435) +++ trunk/psModules/src/detrend/pmNonLinear.c (revision 5516) @@ -5,6 +5,6 @@ * @author GLG, MHPCC * - * @version $Revision: 1.2 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-10-20 23:06:24 $ + * @version $Revision: 1.3 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-11-15 20:09:03 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -22,4 +22,5 @@ #include "pmNonLinear.h" +#include "pmSubtractBias.h" /****************************************************************************** @@ -37,11 +38,17 @@ PS_ASSERT_IMAGE_TYPE(inputReadout->image, PS_TYPE_F32, NULL); PS_ASSERT_PTR_NON_NULL(input1DPoly, NULL); + PS_WARN_PTR_NON_NULL(inputReadout->parent); + if (inputReadout->parent != NULL) { + PS_WARN_PTR_NON_NULL(inputReadout->parent->concepts); + } + // + // Determine trimmed image from metadata. + // + psImage *trimmedImg = p_psDetermineTrimmedImage(inputReadout); - psS32 i; - psS32 j; - - for (i=0;iimage->numRows;i++) { - for (j=0;jimage->numCols;j++) { - inputReadout->image->data.F32[i][j] = psPolynomial1DEval(input1DPoly, inputReadout->image->data.F32[i][j]); + for (psS32 i=0;inumRows;i++) { + for (psS32 j=0;jnumCols;j++) { + trimmedImg->data.F32[i][j] = psPolynomial1DEval(input1DPoly, + trimmedImg->data.F32[i][j]); } } @@ -68,4 +75,13 @@ PS_ASSERT_PTR_NON_NULL(inputReadout->image,NULL); PS_ASSERT_IMAGE_TYPE(inputReadout->image, PS_TYPE_F32, NULL); + PS_WARN_PTR_NON_NULL(inputReadout->parent); + if (inputReadout->parent != NULL) { + PS_WARN_PTR_NON_NULL(inputReadout->parent->concepts); + } + // + // Determine trimmed image from metadata. + // + psImage *trimmedImg = p_psDetermineTrimmedImage(inputReadout); + psLookupTable *tmpLT = psLookupTableAlloc(filename, "%f %f", 0); psS32 numLines = psLookupTableRead(tmpLT); @@ -75,9 +91,9 @@ "WARNING: Lookup Table is too small. Returning original pmReadout.\n"); } else { - for (psS32 i=0;iimage->numRows;i++) { - for (psS32 j=0;jimage->numCols;j++) { - psF64 tmpD = psLookupTableInterpolate(tmpLT, inputReadout->image->data.F32[i][j], 1); + for (psS32 i=0;inumRows;i++) { + for (psS32 j=0;jnumCols;j++) { + psF64 tmpD = psLookupTableInterpolate(tmpLT, trimmedImg->data.F32[i][j], 1); if (!isnan(tmpD)) { - inputReadout->image->data.F32[i][j] = tmpD; + trimmedImg->data.F32[i][j] = tmpD; } else { numPixels++;