Changeset 5795

Timestamp:

Dec 16, 2005, 5:18:39 PM (21 years ago)

Author:

Paul Price

Message:

Importing PAP code from phase 2. NOTE: Reverted files in detrend/ and imsubtract/ to a previous version. DO NOT merge these files back!

Location:

branches/pap_branch_051214/psModules/src

Files:

• astrom/Makefile.am (modified) (1 diff)
• astrom/pmAstrometry.c (modified) (16 diffs)
• astrom/pmAstrometry.h (modified) (15 diffs)
• camera/Makefile.am (modified) (1 diff)
• detrend/pmFlatField.c (modified) (10 diffs)
• detrend/pmFlatField.h (modified) (3 diffs)
• detrend/pmMaskBadPixels.c (modified) (5 diffs)
• detrend/pmMaskBadPixels.h (modified) (4 diffs)
• detrend/pmMaskBadPixelsErrors.h (modified) (3 diffs)
• detrend/pmNonLinear.c (modified) (5 diffs)
• detrend/pmNonLinear.h (modified) (3 diffs)
• imsubtract/Makefile.am (modified) (1 diff)
• imsubtract/pmSubtractBias.c (modified) (11 diffs)
• imsubtract/pmSubtractBias.h (modified) (4 diffs)

branches/pap_branch_051214/psModules/src/astrom/Makefile.am

@@ -4 +4 @@
 libpsmoduleastrom_la_LDFLAGS  = -release $(PACKAGE_VERSION)
 libpsmoduleastrom_la_SOURCES  = \
-    pmAstrometry.c \
-    pmAstrometryObjects.c
+        pmAstrometry.c \
+        pmAstrometryObjects.c \
+        pmChipMosaic.c \
+        pmFPAConceptsGet.c \
+        pmFPAConceptsSet.c \
+        pmFPAConstruct.c \
+        pmFPARead.c \
+        pmFPAWrite.c \
+        pmReadout.c \
+        psAdditionals.c
 
 psmoduleincludedir = $(includedir)
 psmoduleinclude_HEADERS = \
-  pmAstrometry.h \
-  pmAstrometryObjects.h
+        pmAstrometry.h \
+        pmAstrometryObjects.h \
+        pmChipMosaic.h \
+        pmFPAConceptsGet.h \
+        pmFPAConceptsSet.h \
+        pmFPAConstruct.h \
+        pmFPARead.h \
+        pmFPAWrite.h \
+        pmReadout.h \
+        psAdditionals.h

branches/pap_branch_051214/psModules/src/astrom/pmAstrometry.c

@@ -10 +10 @@
 * XXX: We should review the extent of the warning messages on these functions
 * when the transformations are not successful.
-* 
+*
 * XXX: Should we implement non-linear cell->chip transforms?
-* 
-*  @version $Revision: 1.11 $ $Name: not supported by cvs2svn $
-*  @date $Date: 2005-12-08 00:00:57 $
+*
+*  @version $Revision: 1.11.2.1 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2005-12-17 03:18:38 $
 *
 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
@@ -60 +60 @@
         psFree(readout->mask);
         psFree(readout->weight);
-        //
-        // XXX: Not sure if this is the right way to do things.  Currently the psListAdd()
-        // increase the memory reference counter to the list data.  So, we
-        // iterate through the list, and decrement the reference counters.
-        //
-        if (1) {
-            if ((readout->bias != NULL) && (readout->bias->head != NULL)) {
-                psListElem *tmpElem = (psListElem *) readout->bias->head;
-                while (NULL != tmpElem) {
-                    psMemDecrRefCounter((psImage *) tmpElem->data);
-                    tmpElem = tmpElem->next;
-                }
-            }
-        }
-        psFree(readout->bias);
         psFree(readout->analysis);
         psFree(readout->parent);
@@ -88 +73 @@
         psFree(cell->concepts);
         psFree(cell->analysis);
+        psFree(cell->camera);
         //
         // Set the parent to NULL in all cell->readouts before psFree(cell->readouts)
@@ -101 +87 @@
         psFree(cell->readouts);
         psFree(cell->parent);
-        psFree(cell->private);
+        psFree(cell->hdu);
     }
 }
@@ -125 +111 @@
         psFree(chip->cells);
         psFree(chip->parent);
-        psFree(chip->private);
+        psFree(chip->hdu);
     }
 }
@@ -151 +137 @@
         }
         psFree(fpa->chips);
-        psFree(fpa->private);
+        psFree(fpa->hdu);
         psFree(fpa->phu);
     }
 }
 
+void p_pmHDUFree(p_pmHDU *hdu)
+{
+    if (hdu) {
+        psFree(hdu->extname);
+        psFree(hdu->header);
+        psFree(hdu->images);
+        psFree(hdu->masks);
+        psFree(hdu->weights);
+    }
+}
+
 // XXX: Verify these default values for row0, col0, rowBins, colBins
+// PAP: These values may disappear in the future in favour of values in parent->concepts?
 pmReadout *pmReadoutAlloc(pmCell *cell)
 {
     pmReadout *tmpReadout = (pmReadout *) psAlloc(sizeof(pmReadout));
 
-    tmpReadout->col0 = -1;
-    tmpReadout->row0 = -1;
-    tmpReadout->colBins = -1;
-    tmpReadout->rowBins = -1;
+    tmpReadout->col0 = 0;
+    tmpReadout->row0 = 0;
+    tmpReadout->colBins = 0;
+    tmpReadout->rowBins = 0;
     tmpReadout->image = NULL;
     tmpReadout->mask = NULL;
     tmpReadout->weight = NULL;
-    tmpReadout->bias = NULL;
     tmpReadout->analysis = psMetadataAlloc();
     tmpReadout->parent = cell;
@@ -179 +176 @@
 
 // XXX: Verify these default values for row0, col0.
+// PAP: These values may disappear in the future in favour of values in the "concepts"?
 pmCell *pmCellAlloc(
     pmChip *chip,
-    psMetadata *cameradata,
-    psString name)
+    psMetadata *cameraData,
+    const char *name)
 {
     pmCell *tmpCell = (pmCell *) psAlloc(sizeof(pmCell));
 
-    tmpCell->col0 = -1;
-    tmpCell->row0 = -1;
+    tmpCell->col0 = 0;
+    tmpCell->row0 = 0;
     tmpCell->toChip = NULL;
     tmpCell->toFPA = NULL;
@@ -196 +194 @@
         psLogMsg(__func__, PS_LOG_WARN, "WARNING: Could not add CELL.NAME to metadata.\n");
     }
-    tmpCell->camera = cameradata;
-    tmpCell->analysis = NULL;
+    tmpCell->camera = psMemIncrRefCounter(cameraData);
+    tmpCell->analysis = psMetadataAlloc();
     tmpCell->readouts = psArrayAlloc(0);
-    tmpCell->parent = chip;
+    tmpCell->parent = psMemIncrRefCounter(chip);
     if (chip != NULL) {
         chip->cells = psArrayAdd(chip->cells, 1, (psPtr) tmpCell);
     }
-    tmpCell->valid = false;
-    tmpCell->private = NULL;
+    tmpCell->valid = true;              // All cells are valid by default
+    tmpCell->hdu = NULL;
 
     psMemSetDeallocator(tmpCell, (psFreeFunc) cellFree);
@@ -211 +209 @@
 
 // XXX: Verify these default values for row0, col0.
+// PAP: row0, col0 may disappear in the future in favour of storing values in the "concepts".
 pmChip *pmChipAlloc(
     pmFPA *fpa,
-    psString name)
+    const char *name)
 {
     pmChip *tmpChip = (pmChip *) psAlloc(sizeof(pmChip));
 
-    tmpChip->col0 = -1;
-    tmpChip->row0 = -1;
+    tmpChip->col0 = 0;
+    tmpChip->row0 = 0;
     tmpChip->toFPA = NULL;
     tmpChip->fromFPA = NULL;
@@ -226 +225 @@
         psLogMsg(__func__, PS_LOG_WARN, "WARNING: Could not add CHIP.NAME to metadata.\n");
     }
-    tmpChip->analysis = NULL;
+    tmpChip->analysis = psMetadataAlloc();
     tmpChip->cells = psArrayAlloc(0);
-    tmpChip->parent = fpa;
+    tmpChip->parent = psMemIncrRefCounter(fpa);
     if (fpa != NULL) {
         fpa->chips = psArrayAdd(fpa->chips, 1, (psPtr) tmpChip);
     }
-    tmpChip->valid = false;
-    tmpChip->private = NULL;
+    tmpChip->valid = true;              // Work on all chips, by default
+    tmpChip->hdu = NULL;
 
     psMemSetDeallocator(tmpChip, (psFreeFunc) chipFree);
@@ -248 +247 @@
     tmpFPA->concepts = psMetadataAlloc();
     tmpFPA->analysis = NULL;
-    tmpFPA->camera = camera;
+    tmpFPA->camera = psMemIncrRefCounter((psPtr)camera);
     tmpFPA->chips = psArrayAlloc(0);
-    tmpFPA->private = NULL;
+    tmpFPA->hdu = NULL;
     tmpFPA->phu = NULL;
 
     psMemSetDeallocator(tmpFPA, (psFreeFunc) FPAFree);
     return(tmpFPA);
+}
+
+p_pmHDU *p_pmHDUAlloc(const char *extname)
+{
+    p_pmHDU *hdu = psAlloc(sizeof(p_pmHDU));
+    psMemSetDeallocator(hdu, (psFreeFunc)p_pmHDUFree);
+
+    hdu->extname = psStringCopy(extname);
+    hdu->header = NULL;
+    hdu->images = NULL;
+    hdu->masks = NULL;
+    hdu->weights = NULL;
+
+    return hdu;
 }
 
@@ -800 +813 @@
 }
 
+
+/*****************************************************************************
+ *****************************************************************************/
+
+// Set cells within a chip to a specified validity
+static bool setCellsValid(const pmChip *chip, // Chip of interest
+                          bool valid    // Valid?
+                         )
+{
+    PS_ASSERT_PTR_NON_NULL(chip, false);
+
+    psArray *cells = chip->cells;       // Component cells
+    if (! cells) {
+        return false;
+    }
+    for (int i = 0; i < cells->n; i++) {
+        pmCell *tmpCell = cells->data[i]; // Cell of interest
+        if (tmpCell) {
+            tmpCell->valid = valid;
+        }
+    }
+
+    return true;
+}
+
 /*****************************************************************************
  *****************************************************************************/
@@ -807 +845 @@
 {
     PS_ASSERT_PTR_NON_NULL(fpa, false);
-    if ((fpa->chips == NULL) || (chipNum >= fpa->chips->n)) {
+
+    psArray *chips = fpa->chips;        // Component chips
+    if ((chips == NULL) || (chipNum >= chips->n)) {
         return(false);
     }
-    psBool rc = true;
-
-    for (psS32 i = 0 ; i < fpa->chips->n ; i++) {
-        pmChip *tmpChip = (pmChip *) fpa->chips->data[i];
+
+    for (int i = 0 ; i < chips->n ; i++) {
+        pmChip *tmpChip = (pmChip *) chips->data[i];
         if (tmpChip == NULL) {
-            rc = false;
+            continue;
+        }
+        if (i == chipNum) {
+            tmpChip->valid = true;
+            setCellsValid(tmpChip, true);
         } else {
-            if (i == chipNum) {
-                tmpChip->valid = true;
-            } else {
-                tmpChip->valid = false;
-            }
-        }
-    }
-
-    return(rc);
+            tmpChip->valid = false;
+            setCellsValid(tmpChip, false);
+        }
+
+    }
+
+    return true;
 }
 
@@ -831 +872 @@
 /*****************************************************************************
 XXX: The SDRS is ambiguous on a few things:
-    Whether or not the other chips should be set valid=true.
-    Should we return the number of chip valid=true before or after they're set,
+    Whether or not the other chips should be set valid=true. [PAP: No]
+    Should we return the number of chip valid=true before or after they're set, [PAP: After]
  *****************************************************************************/
 /**
- * 
+ *
  * pmFPAExcludeChip shall set valid to false only for the specified chip
  * number (chipNum). In the event that the specified chip number does not exist
@@ -841 +882 @@
  * The function shall return the number of chips within the fpa that have valid
  * set to true.
- *  
+ *
  */
 int pmFPAExcludeChip(
@@ -849 +890 @@
     PS_ASSERT_PTR_NON_NULL(fpa, false);
 
-    if (fpa->chips == NULL) {
+    psArray *chips = fpa->chips;        // Component chips
+    if (chips == NULL) {
         psLogMsg(__func__, PS_LOG_WARN, "WARNING: fpa->chips == NULL\n");
         return(0);
     }
-    if ((chipNum >= fpa->chips->n) || (NULL == (pmChip *) fpa->chips->data[chipNum])) {
+    if ((chipNum >= chips->n) || (NULL == (pmChip *) chips->data[chipNum])) {
         psLogMsg(__func__, PS_LOG_WARN, "WARNING: the specified chip (%d) does not exist.\n", chipNum);
         return(0);
     }
 
-    psS32 numChips = 0;
-    for (psS32 i = 0 ; i < fpa->chips->n ; i++) {
-        pmChip *tmpChip = (pmChip *) fpa->chips->data[i];
+    int numChips = 0;                   // Number of valid chips
+    for (int i = 0 ; i < chips->n ; i++) {
+        pmChip *tmpChip = (pmChip *) chips->data[i]; // Chip of interest
         if (tmpChip != NULL) {
             if (i == chipNum) {
                 tmpChip->valid = false;
-            } else {
-                tmpChip->valid = true;
+                setCellsValid(tmpChip, false); // Wipe out the cell as well
+            } else if (tmpChip->valid) {
                 numChips++;
             }

branches/pap_branch_051214/psModules/src/astrom/pmAstrometry.h

@@ -8 +8 @@
 *  @author GLG, MHPCC
 *
-*  @version $Revision: 1.5 $ $Name: not supported by cvs2svn $
-*  @date $Date: 2005-12-08 00:00:57 $
+*  @version $Revision: 1.5.2.1 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2005-12-17 03:18:38 $
 *
 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
@@ -30 +30 @@
 {
     const char *extname;                // Extension name, if it corresponds to this level
-    psArray *pixels;                    // The pixel data, if it corresponds to this level
     psMetadata *header;                 // The FITS header, if it corresponds to this level
+    psArray *images;                    // The pixel data, if it corresponds to this level
+    psArray *masks;                     // The mask data, if it corresponds to this level
+    psArray *weights;                   // The weight data, if it corresponds to this level
 }
 p_pmHDU;
 
 /** Focal plane data structure
- * 
+ *
  *  A focal plane consists of one or more chips (according to the number of
  *  pieces of contiguous silicon). It contains metadata containers for the
@@ -49 +51 @@
  *  transformation which will be derived from numerically inverting the forward
  *  transformation.
- *  
+ *
  */
 typedef struct
@@ -62 +64 @@
     const psMetadata *camera;           ///< Camera configuration
     psArray *chips;                     ///< The chips
-    p_pmHDU *private;                   ///< FITS data
+    p_pmHDU *hdu;                       ///< FITS data
     psMetadata *phu;                    ///< Primary Header
 }
@@ -68 +70 @@
 
 /** Chip data structure
- *  
+ *
  *  A chip consists of one or more cells (according to the number of amplifiers
  *  on the device). The chip contains metadata containers for the concepts and
@@ -79 +81 @@
  *  inverting the forward transformation. A boolean indicates whether the chip is
  *  of interest, allowing it to be excluded from analysis.
- *  
+ *
  */
 typedef struct
@@ -95 +97 @@
     pmFPA *parent;                      ///< Parent FPA
     bool valid;                         ///< Do we bother about reading and working with this chip?
-    p_pmHDU *private;                   ///< FITS data
+    p_pmHDU *hdu;                       ///< FITS data
 }
 pmChip;
@@ -124 +126 @@
     pmChip *parent;                     ///< Parent chip
     bool valid;                         ///< Do we bother about reading and working with this cell?
-    p_pmHDU *private;                   ///< FITS data
+    p_pmHDU *hdu;                       ///< FITS data
 }
 pmCell;
@@ -141 +143 @@
 {
     // Position on the cell
+    // XXX: These may be removed in the future; use parent->concepts instead?
     int col0;                           ///< Offset from the left of chip.
     int row0;                           ///< Offset from the bottom of chip.
@@ -149 +152 @@
     psImage *mask;                      ///< Mask of input image
     psImage *weight;                    ///< Weight of input image
-    psList *bias;                       ///< List of bias section (sub-)images
     psMetadata *analysis;               ///< Readout-level analysis metadata
     pmCell *parent;                     ///< Parent cell
@@ -180 +182 @@
  */
 pmCell *pmCellAlloc(
-    pmChip *chip,                        ///< Parent chip
-    psMetadata *cameradata,
-    psString name
+    pmChip *chip,       ///< Parent chip
+    psMetadata *cameraData, ///< Camera data
+    const char *name    ///< Name of cell
 );
 
 /** Allocates a pmChip
- * 
+ *
  *  The constructor shall make an empty pmChip. If the parent fpa is not NULL,
  *  the parent link is made and the chip shall be placed in the parent's array
@@ -196 +198 @@
  */
 pmChip *pmChipAlloc(
-    pmFPA *fpa,
-    psString name
-
+    pmFPA *fpa,                         ///< FPA to which the chip belongs
+    const char *name                    ///< Name of chip
 );
 
 /** Allocates a pmFPA
- * 
+ *
  *  The constructor shall make an empty pmFPA. The chips array shall be
  *  allocated with a zero size, the camera and db pointers set to the values
  *  provided, and the concepts metadata constructed. All other pointers in the
  *  structure shall be initialized to NULL.
- * 
+ *
  */
 pmFPA *pmFPAAlloc(
@@ -213 +214 @@
 );
 
+/** Allocates a p_pmHDU
+ *
+ * XXX: More detailed description
+ *
+ * @return p_pmHDU*    newly allocated p_pmHDU
+ */
+p_pmHDU *p_pmHDUAlloc(const char *extname // Extension name
+                     );
+
 
 /** Verify parent links.
- * 
+ *
  *  This function checks the validity of the parent links in the FPA hierarchy.
  *  If a parent link is not set (or not set correctly), it is corrected, and the
  *  function shall return false. If all the parent pointers were correct, the
  *  function shall return true.
- * 
+ *
  */
 bool pmFPACheckParents(
@@ -229 +239 @@
 
 /** FUNC DESC
- * 
- * 
- * 
- * 
- */
-
+ *
+ *
+ *
+ *
+ */
+
+
+
+/**
+ *
+ * pmFPASelectChip shall set valid to true for the specified chip number
+ * (chipNum), and all other chips shall have valid set to false. In the event
+ * that the specified chip number does not exist within the fpa, the function
+ * shall return false.
+ *
+ */
+bool pmFPASelectChip(
+    pmFPA *fpa,
+    int chipNum
+);
+
+/**
+ *
+ * pmFPAExcludeChip shall set valid to false only for the specified chip
+ * number (chipNum). In the event that the specified chip number does not exist
+ * within the fpa, the function shall generate a warning, and perform no action.
+ * The function shall return the number of chips within the fpa that have valid
+ * set to true.
+ *
+ */
+int pmFPAExcludeChip(
+    pmFPA *fpa,
+    int chipNum
+);
 
 
@@ -423 +461 @@
 );
 
-/**
- * 
- * pmFPASelectChip shall set valid to true for the specified chip number
- * (chipNum), and all other chips shall have valid set to false. In the event
- * that the specified chip number does not exist within the fpa, the function
- * shall return false.
- *  
- */
-bool pmFPASelectChip(
-    pmFPA *fpa,
-    int chipNum
-);
-
-/**
- * 
- * pmFPAExcludeChip shall set valid to false only for the specified chip
- * number (chipNum). In the event that the specified chip number does not exist
- * within the fpa, the function shall generate a warning, and perform no action.
- * The function shall return the number of chips within the fpa that have valid
- * set to true.
- *  
- */
-int pmFPAExcludeChip(
-    pmFPA *fpa,
-    int chipNum
-);
-
 
 #endif // #ifndef PS_ASTROMETRY_H

branches/pap_branch_051214/psModules/src/camera/Makefile.am

@@ -3 +3 @@
 libpsmodulecamera_la_CPPFLAGS = $(SRCINC) $(PSMODULE_CFLAGS)
 libpsmodulecamera_la_LDFLAGS  = -release $(PACKAGE_VERSION)
-libpsmodulecamera_la_SOURCES  =
+libpsmodulecamera_la_SOURCES  = \
+        pmChipMosaic.c \
+        pmFPAConceptsGet.c \
+        pmFPAConceptsSet.c \
+        pmFPAConstruct.c \
+        pmFPAMorph.c \
+        pmFPARead.c \
+        pmFPAWrite.c
 
 psmoduleincludedir = $(includedir)
-psmoduleinclude_HEADERS =
-
+psmoduleinclude_HEADERS = \
+        pmChipMosaic.h \
+        pmFPAConceptsGet.h \
+        pmFPAConceptsSet.h \
+        pmFPAConstruct.h \
+        pmFPAMorph.h \
+        pmFPARead.h \
+        pmFPAWrite.h

branches/pap_branch_051214/psModules/src/detrend/pmFlatField.c

@@ +1 @@
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+
+
 /** @file  pmFlatField.c
  *
@@ -18 +24 @@
  *  @author Ross Harman, MHPCC
  *
- *  @version $Revision: 1.4 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-11-18 19:43:14 $
+ *  @version $Revision: 1.4.8.1 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-12-17 03:18:39 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -30 +36 @@
 #include<stdio.h>
 #include<math.h>
-#include <string.h>
+#include <strings.h>
 
 #include "pslib.h"
@@ -36 +42 @@
 #include "pmMaskBadPixels.h"
 #include "pmFlatFieldErrors.h"
-#include "pmSubtractBias.h"
 
-// XXX: This should be removed when the autoconf stuff handles psConstants.h correctly.
-#define PS_WARN_PTR_NON_NULL(NAME) \
-if ((NAME) == NULL) { \
-    psLogMsg(__func__, PS_LOG_WARN, "WARNING: %s is NULL.", #NAME); \
-} \
 
-bool pmFlatField(
-    pmReadout *in,
-    const pmReadout *flat)
+bool pmFlatField(pmReadout *in, pmReadout *mask, 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);
-    }
+    PS_ASSERT_PTR_NON_NULL(mask, false);
     int i = 0;
     int j = 0;
@@ -71 +55 @@
     psElemType flatType;
     psElemType maskType;
+    psImage *inImage = NULL;
     psImage *inMask = NULL;
     psImage *flatImage = NULL;
 
-    //
-    // Determine trimmed image from metadata.
-    //
-    psImage *trimmedImg = p_psDetermineTrimmedImage(in);
+
+    // Check for nulls
+    if (in == NULL) {
+        return true;       // Readout may not 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;
     flatImage = flat->image;
-    inMask = in->mask;
+    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 = mask->image;
 
     // Check input image and its mask are not larger than flat image
 
-    if (trimmedImg == NULL)
-        printf("XXX: 00\n");
-    if (flatImage == NULL)
-        printf("XXX 01\n");
-
-    if (trimmedImg->numRows>flatImage->numRows || trimmedImg->numCols>flatImage->numCols) {
+    if (inImage->numRows>flatImage->numRows || inImage->numCols>flatImage->numCols) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmFlatField_SIZE_INPUT_IMAGE,
-                 trimmedImg->numRows, trimmedImg->numCols, flatImage->numRows, flatImage->numCols);
+                 inImage->numRows, inImage->numCols, flatImage->numRows, flatImage->numCols);
         return false;
     }
@@ -102 +98 @@
 
     // Determine total offset based on image offset with chip offset
-    totOffCol = trimmedImg->col0 + in->col0;
-    totOffRow = trimmedImg->row0 + in->row0;
+    totOffCol = inImage->col0 + in->col0;
+    totOffRow = inImage->row0 + in->row0;
 
     // Check that offsets are within image limits
@@ -111 +107 @@
                  totOffRow, totOffCol, flatImage->numRows, flatImage->numCols);
         return false;
-    } else if(totOffRow>=trimmedImg->numRows || totOffCol>=trimmedImg->numCols) {
+    } else if(totOffRow>=inImage->numRows || totOffCol>=inImage->numCols) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmFlatField_OFFSET_INPUT_IMAGE,
-                 totOffRow, totOffCol, trimmedImg->numRows, trimmedImg->numCols);
+                 totOffRow, totOffCol, inImage->numRows, inImage->numCols);
         return false;
     } else if(totOffRow>=inMask->numRows || totOffCol>=inMask->numCols) {
@@ -124 +120 @@
 
     // Check for incorrect types
-    inType = trimmedImg->type.type;
+    inType = inImage->type.type;
     flatType = flatImage->type.type;
     maskType = inMask->type.type;
@@ -153 +149 @@
 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 < trimmedImg->numRows; j++) {                                                          \
-        for(i = totOffCol; i < trimmedImg->numCols; i++) {                                                      \
+    for(j = totOffRow; j < inImage->numRows; j++) {                                                          \
+        for(i = totOffCol; i < inImage->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 */           \
@@ -162 +158 @@
         }                                                                                                    \
     }                                                                                                        \
-    for(j = totOffRow; j < trimmedImg->numRows; j++) {                                                          \
-        for(i = totOffCol; i < trimmedImg->numCols; i++) {                                                      \
+    for(j = totOffRow; j < inImage->numRows; j++) {                                                          \
+        for(i = totOffCol; i < inImage->numCols; i++) {                                                      \
             if(!inMask->data.PS_TYPE_MASK_DATA[j][i]) {                                                      \
                 /* Module shall divide the input image by the flat-fielded image */                          \
-                trimmedImg->data.TYPE[j][i] /= flatImage->data.TYPE[j][i];                                      \
+                inImage->data.TYPE[j][i] /= flatImage->data.TYPE[j][i];                                      \
             }                                                                                                \
         }                                                                                                    \

branches/pap_branch_051214/psModules/src/detrend/pmFlatField.h

@@ +1 @@
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+
+
 /** @file  pmFlatField.h
  *
@@ -18 +24 @@
  *  @author Ross Harman, MHPCC
  *
- *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-10-20 23:06:24 $
+ *  @version $Revision: 1.2.8.1 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-12-17 03:18:39 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -37 +43 @@
 bool pmFlatField(
     pmReadout *in,          ///< Readout with input image
+    pmReadout *mask,        ///< Input image mask
     const pmReadout *flat   ///< Readout with flat image
 );

branches/pap_branch_051214/psModules/src/detrend/pmMaskBadPixels.c

@@ +1 @@
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+
 /** @file  pmMaskBadPixels.c
  *
@@ -19 +24 @@
  *  @author Ross Harman, MHPCC
  *
- *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-11-18 19:43:14 $
+ *  @version $Revision: 1.3.8.1 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-12-17 03:18:39 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -31 +36 @@
 #include<stdio.h>
 #include<math.h>
-#include<string.h>
+#include<strings.h>
 
 #include "pmMaskBadPixels.h"
 #include "pmMaskBadPixelsErrors.h"
-#include "pmSubtractBias.h"
 
-//XXX: REmove, autoconf is broken.
-#define PS_WARN_PTR_NON_NULL(NAME) \
-if ((NAME) == NULL) { \
-    psLogMsg(__func__, PS_LOG_WARN, "WARNING: %s is NULL.", #NAME); \
-} \
+bool pmMaskBadPixels(pmReadout *in, const psImage *mask, unsigned int maskVal, float sat,
+                     unsigned int growVal, int grow)
+{
+    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;
 
 
-/******************************************************************************
-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);
+    // 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;
+    }
 
-
-    for(psS32 i=rowMin; i<rowMax; i++) {
-        for(psS32 j=colMin; j<colMax; j++) {
-            // Old code:
-            //            psF32 rRound = 0.5 + sqrtf((psF32) (((col-j)*(col-j)) + ((row-i)*(row-i))));
-            //            if(rRound <= radius) {
-            psF32 squareDis = (psF32) (((col-j)*(col-j)) + ((row-i)*(row-i)));
-            if (squareDis <= squareRadius) {
-                inMask->data.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 pmReadout *mask,
-    unsigned int maskVal,
-    float sat,
-    unsigned int growVal,
-    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 totOffCol = 0;
-    int totOffRow = 0;
-    psElemType inType;
-    psElemType maskType;
-    psImage *inMask = NULL;
-
-    //
-    // 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);
+    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));
     }
     inMask = in->mask;
 
     // Check input image and its mask are not larger than mask
-    if(trimmedImg->numRows > mask->image->numRows || trimmedImg->numCols > mask->image->numCols) {
+    if(inImage->numRows > mask->numRows || inImage->numCols > mask->numCols) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmMaskBadPixels_SIZE_INPUT_IMAGE,
-                 trimmedImg->numRows, trimmedImg->numCols, mask->image->numRows, mask->image->numCols);
+                 inImage->numRows, inImage->numCols, mask->numRows, mask->numCols);
         return false;
-    } else if(inMask->numRows > mask->image->numRows || inMask->numCols > mask->image->numCols) {
+    } else if(inMask->numRows>mask->numRows || inMask->numCols>mask->numCols) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmMaskBadPixels_SIZE_MASK_IMAGE,
-                 inMask->numRows, inMask->numCols, mask->image->numRows, mask->image->numCols);
+                 inMask->numRows, inMask->numCols, mask->numRows, mask->numCols);
         return false;
     }
 
     // Determine total offset based on image offset with chip offset
-    totOffCol = trimmedImg->col0 + in->col0;
-    totOffRow = trimmedImg->row0 + in->row0;
+    totOffCol = inImage->col0 + in->col0;
+    totOffRow = inImage->row0 + in->row0;
 
     // Check that offsets are within image limits
-    if(totOffRow>=mask->image->numRows || totOffCol>=mask->image->numCols) {
+    if(totOffRow>=mask->numRows || totOffCol>=mask->numCols) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmMaskBadPixels_OFFSET_MASK_IMAGE,
-                 totOffRow, totOffCol, mask->image->numRows, mask->image->numCols);
+                 totOffRow, totOffCol, mask->numRows, mask->numCols);
         return false;
-    } else if(totOffRow>=trimmedImg->numRows || totOffCol>=trimmedImg->numCols) {
+    } else if(totOffRow>=inImage->numRows || totOffCol>=inImage->numCols) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmMaskBadPixels_OFFSET_INPUT_IMAGE,
-                 totOffRow, totOffCol, trimmedImg->numRows, trimmedImg->numCols);
+                 totOffRow, totOffCol, inImage->numRows, inImage->numCols);
         return false;
     } else if(totOffRow>=inMask->numRows || totOffCol>=inMask->numCols) {
@@ -193 +118 @@
 
     // Check for incorrect types
-    inType = trimmedImg->type.type;
-    maskType = mask->image->type.type;
+    inType = inImage->type.type;
+    maskType = mask->type.type;
     if(PS_IS_PSELEMTYPE_COMPLEX(inType)) {
         psError( PS_ERR_BAD_PARAMETER_TYPE, true,
@@ -207 +132 @@
     }
 
-    //
-    // 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; i<trimmedImg->numRows; i++) {
-        for(j=totOffCol; j<trimmedImg->numCols; j++) {
-            //
-            // A: Pixels which satisfy maskVal shall be masked
-            //
-            if (mask->image->data.U8[i][j] == maskVal) {
-                in->mask->data.U8[i][j] |= maskVal;
-            }
+    // Macro for all PS types
+    #define PM_BAD_PIXELS(TYPE)                                                                              \
+case PS_TYPE_##TYPE:                                                                                         \
+    for(j=totOffRow; j<inImage->numRows; j++) {                                                              \
+        for(i=totOffCol; i<inImage->numCols; 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; jj<rowMax; jj++) {                                                            \
+                    for(ii=colMin; ii<colMax; ii++) {                                                        \
+                        r = sqrtf((ii-i)*(ii-i)+(jj-j)*(jj-j));                                              \
+                        rRound = r + 0.5;                                                                    \
+                        if(rRound <= grow) {                                                                 \
+                            inMask->data.PS_TYPE_MASK_DATA[jj][ii] |=                                        \
+                                    (mask->data.PS_TYPE_MASK_DATA[j][i]&growVal);                            \
+                        }                                                                                    \
+                    }                                                                                        \
+                }                                                                                            \
+            }                                                                                                \
+        }                                                                                                    \
+    }                                                                                                        \
+    break;
 
-            //
-            // 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);
-            }
-        }
+    // 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 true;
+    return false;
 }

branches/pap_branch_051214/psModules/src/detrend/pmMaskBadPixels.h

@@ +1 @@
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+
 /** @file  pmMaskBadPixels.h
  *
@@ -19 +24 @@
  *  @author Ross Harman, MHPCC
  *
- *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-11-15 20:09:03 $
+ *  @version $Revision: 1.2.8.1 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-12-17 03:18:39 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -36 +41 @@
 } pmMaskValue;
 
-// XXX: Use PS_MIN, PS_MAX
 /** Macro to find maximum of two numbers */
 #define MAX(A,B)((A)>=(B)?(A):(B))
@@ -54 +58 @@
 bool pmMaskBadPixels(
     pmReadout *in,          ///< Readout containing input image data.
-    const pmReadout *mask,   ///< Mask data to be added to readout mask data.
+    const psImage *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.

branches/pap_branch_051214/psModules/src/detrend/pmMaskBadPixelsErrors.h

@@ +1 @@
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+
 /** @file  pmMaskBadPixelsErrors.h
  *
@@ -7 +12 @@
  *  @author Ross Harman, MHPCC
  *
- *  @version $Revision: 1.1 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-09-28 20:43:52 $
+ *  @version $Revision: 1.1.16.1 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-12-17 03:18:39 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -23 +28 @@
  *     $2  The error text (rest of the line in pmMaskBadPixelsErrors.h)
  *     $n  The order of the source line in pmMaskBadPixelsErrors.h (comments excluded)
- * 
+ *
  * DO NOT EDIT THE LINES BETWEEN //~Start and //~End!  ANY CHANGES WILL BE OVERWRITTEN.
  */

branches/pap_branch_051214/psModules/src/detrend/pmNonLinear.c

@@ +1 @@
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+
 /** @file  pmNonLinear.c
  *
@@ -5 +10 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.5 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-12-05 20:49:40 $
+ *  @version $Revision: 1.5.8.1 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-12-17 03:18:39 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -22 +27 @@
 
 #include "pmNonLinear.h"
-#include "pmSubtractBias.h"
 
-// XXX: Remove, autoconf must be
-#define PS_WARN_PTR_NON_NULL(NAME) \
-if ((NAME) == NULL) { \
-    psLogMsg(__func__, PS_LOG_WARN, "WARNING: %s is NULL.", #NAME); \
-} \
 /******************************************************************************
 pmNonLinearityLookup(): This routine will take an pmReadout image as input
 and a 1-D polynomial.  For each pixel in the input image, the polynomial will
-be evaluated at that pixels value, and the image pixel will then be set
-to
+be evaluated at that pixels value, and the image pixel will then be set to
 that value.
-*****************************************************************************/
+ *****************************************************************************/
 
-pmReadout *pmNonLinearityPolynomial(
-    pmReadout *inputReadout,
-    const psPolynomial1D *input1DPoly)
+pmReadout *pmNonLinearityPolynomial(pmReadout *inputReadout,
+                                    const psPolynomial1D *input1DPoly)
 {
     PS_ASSERT_PTR_NON_NULL(inputReadout, NULL);
@@ -45 +42 @@
     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 (psS32 i=0;i<trimmedImg->numRows;i++) {
-        for (psS32 j=0;j<trimmedImg->numCols;j++) {
-            trimmedImg->data.F32[i][j] = psPolynomial1DEval(input1DPoly,
-                                         trimmedImg->data.F32[i][j]);
+    for (i=0;i<inputReadout->image->numRows;i++) {
+        for (j=0;j<inputReadout->image->numCols;j++) {
+            inputReadout->image->data.F32[i][j] = psPolynomial1DEval(input1DPoly, inputReadout->image->data.F32[i][j]);
         }
     }
@@ -71 +61 @@
 inFluxe, and the corresponding value in outFlux.  The image pixel will then
 be set to the value from outFlux.
- 
-XXX: Must assert that filename exists.  This should probably happen in
-the lookup files.
  *****************************************************************************/
-pmReadout *pmNonLinearityLookup(
-    pmReadout *inputReadout,
-    const char *filename
-)
+pmReadout *pmNonLinearityLookup(pmReadout *inputReadout,
+                                const psVector *inFlux,
+                                const psVector *outFlux)
 {
     PS_ASSERT_PTR_NON_NULL(inputReadout,NULL);
     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);
+    PS_ASSERT_PTR_NON_NULL(inFlux,NULL);
+    if (inFlux->n < 2) {
+        psError(PS_ERR_UNKNOWN,true, "pmNonLinearityLookup(): input vector less than 2 elements.  Returning inputReadout image.");
+        return(inputReadout);
     }
-    //
-    // Determine trimmed image from metadata.
-    //
-    psImage *trimmedImg = p_psDetermineTrimmedImage(inputReadout);
+    PS_ASSERT_PTR_NON_NULL(outFlux,NULL);
+    psS32 tableSize = inFlux->n;
+    if (inFlux->n != outFlux->n) {
+        tableSize = PS_MIN(inFlux->n, outFlux->n);
+        psLogMsg(__func__, PS_LOG_WARN,
+                 "WARNING: pmNonLinear.c: pmNonLinearityLookup(): input vectors have different sizes (%d, %d)\n", inFlux->n, outFlux->n);
+    }
+    PS_ASSERT_VECTOR_TYPE(inFlux, PS_TYPE_F32, NULL);
+    PS_ASSERT_VECTOR_TYPE(outFlux, PS_TYPE_F32, NULL);
 
-    psLookupTable *tmpLT = psLookupTableAlloc(filename, "%f %f", 0);
-    psS32 numLines = psLookupTableRead(tmpLT);
+    psS32 i;
+    psS32 j;
+    psS32 binNum;
+    psScalar x;
     psS32 numPixels = 0;
-    if (numLines < 2) {
-        psLogMsg(__func__, PS_LOG_WARN,
-                 "WARNING: Lookup Table is too small.  Returning original pmReadout.\n");
-    } else {
-        for (psS32 i=0;i<trimmedImg->numRows;i++) {
-            for (psS32 j=0;j<trimmedImg->numCols;j++) {
-                psF64 tmpD = psLookupTableInterpolate(tmpLT, trimmedImg->data.F32[i][j], 1);
-                if (!isnan(tmpD)) {
-                    trimmedImg->data.F32[i][j] = tmpD;
-                } else {
-                    numPixels++;
-                }
+    psF32 slope;
+
+    x.type.type = PS_TYPE_F32;
+    for (i=0;i<inputReadout->image->numRows;i++) {
+        for (j=0;j<inputReadout->image->numCols;j++) {
+            x.data.F32 = inputReadout->image->data.F32[i][j];
+            binNum = p_psVectorBinDisect((psVector *)inFlux, &x);
+
+            if (binNum == -2) {
+                // We get here if x is below the table lookup range.
+                inputReadout->image->data.F32[i][j] = outFlux->data.F32[0];
+                numPixels++;
+
+            } else if (binNum == -1) {
+                // We get here if x is above the table lookup range.
+                inputReadout->image->data.F32[i][j] = outFlux->data.F32[tableSize-1];
+                numPixels++;
+
+            } else if (binNum < -2) {
+                // We get here if there was some other problem.
+                psError(PS_ERR_UNKNOWN,true, "pmNonLinearityLookup(): Could not perform p_psVectorBinDisect().  Returning inputReadout image.");
+                return(inputReadout);
+                numPixels++;
+            } else {
+                // Perform linear interpolation.
+                slope = (outFlux->data.F32[binNum+1] - outFlux->data.F32[binNum]) /
+                        (inFlux->data.F32[binNum+1]  - inFlux->data.F32[binNum]);
+                inputReadout->image->data.F32[i][j] = outFlux->data.F32[binNum] +
+                                                      ((x.data.F32 - inFlux->data.F32[binNum]) * slope);
             }
         }
-        if (numPixels > 0) {
-            psLogMsg(__func__, PS_LOG_WARN,
-                     "WARNING: pmNonLinear.c: pmNonLinearityLookup(): %d pixels outside table.", numPixels);
-        }
     }
-
+    if (numPixels > 0) {
+        psLogMsg(__func__, PS_LOG_WARN,
+                 "WARNING: pmNonLinear.c: pmNonLinearityLookup(): %d pixels outside table.", numPixels);
+    }
     return(inputReadout);
 }

branches/pap_branch_051214/psModules/src/detrend/pmNonLinear.h

@@ +1 @@
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+
 /** @file  pmNonLinear.h
  *
@@ -5 +10 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-10-20 23:06:24 $
+ *  @version $Revision: 1.2.8.1 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-12-17 03:18:39 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -18 +23 @@
 #include "pmAstrometry.h"
 
-pmReadout *pmNonLinearityPolynomial(
-    pmReadout *in,
-    const psPolynomial1D *coeff
-);
+pmReadout *pmNonLinearityPolynomial(pmReadout *in,
+                                    const psPolynomial1D *coeff);
 
-pmReadout *pmNonLinearityLookup(
-    pmReadout *in,
-    const char *filename
-);
+pmReadout *pmNonLinearityLookup(pmReadout *in,
+                                const psVector *inFlux,
+                                const psVector *outFlux);
 
 #endif

branches/pap_branch_051214/psModules/src/imsubtract/Makefile.am

@@ -4 +4 @@
 libpsmoduleimsubtract_la_LDFLAGS  = -release $(PACKAGE_VERSION)
 libpsmoduleimsubtract_la_SOURCES  = pmImageSubtract.c \
-    pmSubtractBias.c \
-    pmSubtractSky.c
+    pmSubtractBias.c
+#    pmSubtractSky.c
 
 psmoduleincludedir = $(includedir)
 psmoduleinclude_HEADERS = \
     pmImageSubtract.h \
-    pmSubtractBias.h \
-    pmSubtractSky.h
+    pmSubtractBias.h
+#    pmSubtractSky.h

branches/pap_branch_051214/psModules/src/imsubtract/pmSubtractBias.c

@@ +1 @@
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+
 /** @file  pmSubtractBias.c
  *
@@ -6 +11 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.6 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-11-19 00:55:18 $
+ *  @version $Revision: 1.6.8.1 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-12-17 03:18:39 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
  *
  */
-/*****************************************************************************/
-/* INCLUDE FILES                                                             */
-/*****************************************************************************/
-#include <stdio.h>
-#include <math.h>
-#include <string.h>
-#include "pslib.h"
+
 #if HAVE_CONFIG_H
 #include <config.h>
 #endif
+
 #include "pmSubtractBias.h"
 
-/*****************************************************************************/
-/* DEFINE STATEMENTS                                                         */
-/*****************************************************************************/
+#define PM_SUBTRACT_BIAS_POLYNOMIAL_ORDER 2
+#define PM_SUBTRACT_BIAS_SPLINE_ORDER 3
+
 // XXX: put these in psConstants.h
-void PS_POLY1D_PRINT(
-    psPolynomial1D *poly)
+void PS_POLY1D_PRINT(psPolynomial1D *poly)
 {
     printf("-------------- PS_POLY1D_PRINT() --------------\n");
@@ -57 +56 @@
 }\
 
-/*****************************************************************************/
-/* TYPE DEFINITIONS                                                          */
-/*****************************************************************************/
-
-/*****************************************************************************/
-/* GLOBAL VARIABLES                                                          */
-/*****************************************************************************/
-psS32 currentId = 0;                // XXX: remove
-psS32 memLeaks = 0;                 // XXX: remove
-//PRINT_MEMLEAKS(8); XXX
-/*****************************************************************************/
-/* FILE STATIC VARIABLES                                                     */
-/*****************************************************************************/
-
-/*****************************************************************************/
-/* FUNCTION IMPLEMENTATION - LOCAL                                           */
-/*****************************************************************************/
-
 /******************************************************************************
-psSubtractFrame(): this routine will take as input the pmReadout for the input
-image and a pmReadout for the bias image.  The bias image is subtracted in
-place from the input image.  We assume that sizes and types are checked
-elsewhere.
- 
-XXX: Verify that the image and readout offsets are being used the right way.
- 
-XXX: Ensure that it does the correct thing with image size.
+psSubtractFrame(): this routine will take as input a readout for the input
+image and a readout for the bias image.  The bias image is subtracted in
+place from the input image.
 *****************************************************************************/
-static pmReadout *SubtractFrame(
-    pmReadout *in,
-    const pmReadout *bias)
+static pmReadout *SubtractFrame(pmReadout *in,
+                                const pmReadout *bias)
 {
-    // XXX: When did the ->row0 and ->col0 offsets get coded?
-    for (psS32 i=0;i<in->image->numRows;i++) {
-        for (psS32 j=0;j<in->image->numCols;j++) {
+    psS32 i;
+    psS32 j;
+
+    if (bias == NULL) {
+        psLogMsg(__func__, PS_LOG_WARN,
+                 "WARNING: pmSubtractBias.c: SubtractFrame(): bias frame is NULL.  Returning original image.\n");
+        return(in);
+    }
+
+
+    if ((in->image->numRows + in->row0 - bias->row0) > bias->image->numRows) {
+        psError(PS_ERR_UNKNOWN,true, "bias image does not have enough rows.  Returning in image\n");
+        return(in);
+    }
+    if ((in->image->numCols + in->col0 - bias->col0) > bias->image->numCols) {
+        psError(PS_ERR_UNKNOWN,true, "bias image does not have enough columns.  Returning in image\n");
+        return(in);
+    }
+
+    for (i=0;i<in->image->numRows;i++) {
+        for (j=0;j<in->image->numCols;j++) {
             in->image->data.F32[i][j]-=
                 bias->image->data.F32[i+in->row0-bias->row0][j+in->col0-bias->col0];
-
             if ((in->mask != NULL) && (bias->mask != NULL)) {
                 (in->mask->data.U8[i][j])|=
@@ -105 +97 @@
 }
 
-
-/******************************************************************************
-psSubtractDarkFrame(): this routine will take as input the pmReadout for the
-input image and a pmReadout for the dark image.  The dark image is scaled and
-subtracted in place from the input image.
- 
-XXX: Verify that the image and readout offsets are being used the right way.
- 
-XXX: Ensure that it does the correct thing with image size.
-*****************************************************************************/
-static pmReadout *SubtractDarkFrame(
-    pmReadout *in,
-    const pmReadout *dark,
-    psF32 scale)
-{
-    // XXX: When did the ->row0 and ->col0 offsets get coded?
-    if (fabs(scale) > FLT_EPSILON) {
-        for (psS32 i=0;i<in->image->numRows;i++) {
-            for (psS32 j=0;j<in->image->numCols;j++) {
-                in->image->data.F32[i][j]-=
-                    (scale * dark->image->data.F32[i+in->row0-dark->row0][j+in->col0-dark->col0]);
-
-                if ((in->mask != NULL) && (dark->mask != NULL)) {
-                    (in->mask->data.U8[i][j])|=
-                        dark->mask->data.U8[i+in->row0-dark->row0][j+in->col0-dark->col0];
-                }
-            }
-        }
-    } else {
-        for (psS32 i=0;i<in->image->numRows;i++) {
-            for (psS32 j=0;j<in->image->numCols;j++) {
-                in->image->data.F32[i][j]-=
-                    dark->image->data.F32[i+in->row0-dark->row0][j+in->col0-dark->col0];
-
-                if ((in->mask != NULL) && (dark->mask != NULL)) {
-                    (in->mask->data.U8[i][j])|=
-                        dark->mask->data.U8[i+in->row0-dark->row0][j+in->col0-dark->col0];
-                }
-            }
-        }
-    }
-
-    return(in);
-}
-
 /******************************************************************************
 ImageSubtractScalar(): subtract a scalar from the input image.
  
-XXX: Is there a psLib function for this?
+XXX: Use a psLib function for this.
+ 
+XXX: This should
  *****************************************************************************/
-static psImage *ImageSubtractScalar(
-    psImage *image,
-    psF32 scalar)
+static psImage *ImageSubtractScalar(psImage *image,
+                                    psF32 scalar)
 {
     for (psS32 i=0;i<image->numRows;i++) {
@@ -221 +169 @@
 
     if (numOptions == 0) {
-        psError(PS_ERR_UNKNOWN,true, "No allowable statistics options have been specified.\n");
+        psError(PS_ERR_UNKNOWN,true, "No statistics options have been specified.\n");
     }
     if (numOptions != 1) {
@@ -230 +178 @@
 }
 
-/******************************************************************************
-Polynomial1DCopy(): This private function copies the members of the existing
-psPolynomial1D "in" into the existing psPolynomial1D "out".  The previous
-members of the existing psPolynomial1D "out" are psFree'ed.
- *****************************************************************************/
-static psBool Polynomial1DCopy(
-    psPolynomial1D *out,
-    psPolynomial1D *in)
-{
-    psFree(out->coeff);
-    psFree(out->coeffErr);
-    psFree(out->mask);
-
-    out->type = in->type;
-    out->nX = in->nX;
-
-    out->coeff = (psF64 *) psAlloc((in->nX + 1) * sizeof(psF64));
-    // XXX: use memcpy
-    for (psS32 i = 0 ; i < (in->nX + 1) ; i++) {
-        out->coeff[i] = in->coeff[i];
-    }
-
-    out->coeffErr = (psF64 *) psAlloc((in->nX + 1) * sizeof(psF64));
-    // XXX: use memcpy
-    for (psS32 i = 0 ; i < (in->nX + 1) ; i++) {
-        out->coeffErr[i] = in->coeffErr[i];
-    }
-
-    out->mask = (psMaskType *) psAlloc((in->nX + 1) * sizeof(psMaskType));
-    // XXX: use memcpy
-    for (psS32 i = 0 ; i < (in->nX + 1) ; i++) {
-        out->mask[i] = in->mask[i];
-    }
-
-    return(true);
-}
-
-/******************************************************************************
-Polynomial1DDup(): This private function duplicates and then returns the input
-psPolynomial1D "in".
- *****************************************************************************/
-static psPolynomial1D *Polynomial1DDup(
-    psPolynomial1D *in)
-{
-    psPolynomial1D *out = psPolynomial1DAlloc(in->nX, in->type);
-    Polynomial1DCopy(out, in);
-    return(out);
-}
-
-
-/******************************************************************************
-SplineCopy(): This private function copies the members of the existing
-psSpline in into the existing psSpline out.
- *****************************************************************************/
-static psBool SplineCopy(
-    psSpline1D *out,
-    psSpline1D *in)
-{
-    PS_ASSERT_PTR_NON_NULL(out, false);
-    PS_ASSERT_PTR_NON_NULL(in, false);
-
-    for (psS32 i = 0 ; i < out->n ; i++) {
-        psFree(out->spline[i]);
-    }
-    psFree(out->spline);
-    psFree(out->knots);
-    psFree(out->p_psDeriv2);
-
-    out->n = in->n;
-    out->spline = (psPolynomial1D **) psAlloc(in->n * sizeof(psPolynomial1D *));
-    for (psS32 i = 0 ; i < in->n ; i++) {
-        out->spline[i] = Polynomial1DDup(in->spline[i]);
-    }
-
-    // XXX: use psVectorCopy if they get it working.
-    out->knots = psVectorAlloc(in->knots->n, in->knots->type.type);
-    for (psS32 i = 0 ; i < in->knots->n ; i++) {
-        out->knots->data.F32[i] = in->knots->data.F32[i];
-    }
-    /*
-        out->knots = psVectorCopy(out->knots, in->knots, in->knots->type.type);
-    */
-
-    out->p_psDeriv2 = (psF32 *) psAlloc((in->n + 1) * sizeof(psF32));
-    // XXX: use memcpy
-    for (psS32 i = 0 ; i < (in->n + 1) ; i++) {
-        out->p_psDeriv2[i] = in->p_psDeriv2[i];
-    }
-
-    return(true);
-}
+
 
 /******************************************************************************
@@ -328 +186 @@
     PM_FIT_POLYNOMIAL: fit a polynomial to the entire input vector data.
     PM_FIT_SPLINE: fit splines to the input vector data.
-The resulting spline or polynomial is set in the fitSpec argument.
+XXX: Doesn't it make more sense to do polynomial interpolation on a few
+elements of the input vector, rather than fit a polynomial to the entire
+vector?
  *****************************************************************************/
-static psVector *ScaleOverscanVector(
-    psVector *overscanVector,
-    psS32 n,
-    void *fitSpec,
-    pmFit fit)
+static psVector *ScaleOverscanVector(psVector *overscanVector,
+                                     psS32 n,
+                                     void *fitSpec,
+                                     pmFit fit)
 {
     psTrace(".psModule.pmSubtracBias.ScaleOverscanVector", 4,
             "---- ScaleOverscanVector() begin (%d -> %d) ----\n", overscanVector->n, n);
+    //    PS_VECTOR_PRINT_F32(overscanVector);
 
     if (NULL == overscanVector) {
@@ -350 +210 @@
     //
     if (n == overscanVector->n) {
-        return(psVectorCopy(newVec, overscanVector, PS_TYPE_F32));
-    }
+        for (psS32 i = 0 ; i < n ; i++) {
+            newVec->data.F32[i] = overscanVector->data.F32[i];
+        }
+        return(newVec);
+    }
+    psPolynomial1D *myPoly;
+    psSpline1D *mySpline;
     psF32 x;
-
+    psS32 i;
     if (fit == PM_FIT_POLYNOMIAL) {
         // Fit a polynomial to the old overscan vector.
-        psPolynomial1D *myPoly = (psPolynomial1D *) fitSpec;
+        myPoly = (psPolynomial1D *) fitSpec;
         PS_ASSERT_POLY_NON_NULL(myPoly, NULL);
-        PS_ASSERT_POLY1D(myPoly, NULL);
         myPoly = psVectorFitPolynomial1D(myPoly, NULL, 0, overscanVector, NULL, NULL);
         if (myPoly == NULL) {
-            psError(PS_ERR_UNKNOWN, false, "Could not fit a polynomial to the psVector.\n");
+            psError(PS_ERR_UNKNOWN, false, "ScaleOverscanVector()(1): Could not fit a polynomial to the psVector.\n");
             return(NULL);
         }
@@ -368 +232 @@
         // of the old vector, use the fitted polynomial to determine the
         // interpolated value at that point, and set the new vector.
-        for (psS32 i=0;i<n;i++) {
+        for (i=0;i<n;i++) {
             x = ((psF32) i) * ((psF32) overscanVector->n) / ((psF32) n);
             newVec->data.F32[i] = psPolynomial1DEval(myPoly, x);
         }
     } else if (fit == PM_FIT_SPLINE) {
+        psS32 mustFreeSpline = 0;
+        // Fit a spline to the old overscan vector.
+        mySpline = (psSpline1D *) fitSpec;
+        // XXX: Does it make any sense to have a psSpline argument?
+        if (mySpline == NULL) {
+            mustFreeSpline = 1;
+        }
+
         //
         // NOTE: Since the X arg in the psVectorFitSpline1D() function is NULL,
@@ -378 +250 @@
         // properly when doing the spline eval.
         //
-        psSpline1D *mySpline = psVectorFitSpline1D(NULL, overscanVector);
+        //        mySpline = psVectorFitSpline1D(mySpline, NULL, overscanVector, NULL);
+        mySpline = psVectorFitSpline1D(NULL, overscanVector);
         if (mySpline == NULL) {
-            psError(PS_ERR_UNKNOWN, false, "Could not fit a spline to the psVector.\n");
+            psError(PS_ERR_UNKNOWN, false, "ScaleOverscanVector()(2): Could not fit a spline to the psVector.\n");
             return(NULL);
         }
+        //        PS_PRINT_SPLINE(mySpline);
 
         // For each element of the new vector, convert the x-ordinate to that
-        // of the old vector, use the fitted spline to determine the
+        // of the old vector, use the fitted polynomial to determine the
         // interpolated value at that point, and set the new vector.
-        for (psS32 i=0;i<n;i++) {
+        for (i=0;i<n;i++) {
             // Scale to [0 : overscanVector->n - 1]
             x = ((psF32) i) * ((psF32) (overscanVector->n-1)) / ((psF32) n);
             newVec->data.F32[i] = psSpline1DEval(mySpline, x);
         }
-
-        psSpline1D *ptrSpline = (psSpline1D *) fitSpec;
-        if (ptrSpline != NULL) {
-            // Copy the resulting spline fit into ptrSpline.
-            PS_ASSERT_SPLINE(ptrSpline, NULL);
-            SplineCopy(ptrSpline, mySpline);
-        }
-        psFree(mySpline);
+        if (mustFreeSpline ==1) {
+            psFree(mySpline);
+        }
+        //        PS_VECTOR_PRINT_F32(newVec);
+
+
     } else {
         psError(PS_ERR_UNKNOWN, true, "unknown fit type.  Returning NULL.\n");
@@ -412 +284 @@
 
 /******************************************************************************
+XXX: The SDRS does not specify type support.  F32 is implemented here.
  *****************************************************************************/
-static psS32 GetOverscanSize(
-    psImage *inImg,
-    pmOverscanAxis overScanAxis)
-{
-    if (overScanAxis == PM_OVERSCAN_ROWS) {
-        return(inImg->numCols);
-    } else if (overScanAxis == PM_OVERSCAN_COLUMNS) {
-        return(inImg->numRows);
-    } else if (overScanAxis == PM_OVERSCAN_ALL) {
-        return(1);
-    }
-    return(0);
-}
-
-/******************************************************************************
-GetOverscanAxis(in) this private routine determines the appropiate overscan
-axis from the parent cell metadata.
- 
-XXX: Verify the READDIR corresponds with my overscan axis.
- *****************************************************************************/
-static pmOverscanAxis GetOverscanAxis(pmReadout *in)
-{
-    psBool rc;
-    if ((in->parent != NULL) && (in->parent->concepts)) {
-        psS32 dir = psMetadataLookupS32(&rc, in->parent->concepts, "CELL.READDIR");
-        if (rc == true) {
-            if (dir == 1) {
-                return(PM_OVERSCAN_ROWS);
-            } else if (dir == 2) {
-                return(PM_OVERSCAN_COLUMNS);
-            } else if (dir == 3) {
-                return(PM_OVERSCAN_ALL);
-            }
-        }
-    }
-
-    psLogMsg(__func__, PS_LOG_WARN,
-             "WARNING: pmSubtractBias.(): could not determine CELL.READDIR from in->parent metadata.  Setting overscan axis to PM_OVERSCAN_NONE.\n");
-    return(PM_OVERSCAN_NONE);
-}
-
-/******************************************************************************
-psListLength(list): determine the length of a psList.
- 
-XXX: Put this elsewhere.
- *****************************************************************************/
-static psS32 psListLength(
-    psList *list)
-{
-    psS32 length = 0;
-    psListElem *tmpElem = (psListElem *) list->head;
-    while (NULL != tmpElem) {
-        tmpElem = tmpElem->next;
-        length++;
-    }
-    return(length);
-}
-
-/******************************************************************************
-Note: this isn't needed anymore as of psModule SDRS 12-09.
- *****************************************************************************/
-static psBool OverscanReducePixel(
-    psImage *in,
-    psList *bias,
-    psStats *myStats)
-{
-    PS_ASSERT_PTR_NON_NULL(in, NULL);
-    PS_ASSERT_PTR_NON_NULL(bias, NULL);
-    PS_ASSERT_PTR_NON_NULL(bias->head, NULL);
-    PS_ASSERT_PTR_NON_NULL(myStats, NULL);
-
-    // Allocate a psVector with one element per overscan image.
-    psS32 numOverscanImages = psListLength(bias);
-    psVector *statsAll = psVectorAlloc(numOverscanImages, PS_TYPE_F32);
-    psListElem *tmpOverscan = (psListElem *) bias->head;
-    psS32 i = 0;
-    psF64 statValue;
-    //
-    // We loop through each overscan image, calculating the specified
-    // statistic on that image.
-    //
-    while (NULL != tmpOverscan) {
-        psImage *myOverscanImage = (psImage *) tmpOverscan->data;
-
-        PS_ASSERT_IMAGE_TYPE(myOverscanImage, PS_TYPE_F32, NULL);
-        myStats = psImageStats(myStats, myOverscanImage, NULL, (psMaskType)0xffffffff);
-        if (myStats == NULL) {
-            psError(PS_ERR_UNKNOWN, false, "psImageStats(): could not perform requested statistical operation.  Returning in image.\n");
-            psFree(statsAll);
-            return(false);
-        }
-        if (false == p_psGetStatValue(myStats, &statValue)) {
-            psError(PS_ERR_UNKNOWN, false, "p_psGetStatValue(): could not determine result from requested statistical operation.  Returning in image.\n");
-            psFree(statsAll);
-            return(false);
-        }
-        statsAll->data.F32[i] = statValue;
-        i++;
-        tmpOverscan = tmpOverscan->next;
-    }
-
-    //
-    // We reduce the individual stats for each overscan image to
-    // a single psF32.
-    //
-    myStats = psVectorStats(myStats, statsAll, NULL, NULL, 0);
-    if (myStats == NULL) {
-        psError(PS_ERR_UNKNOWN, false, "psImageStats(): could not perform requested statistical operation.  Returning in image.\n");
-        psFree(statsAll);
-        return(false);
-    }
-    if (false == p_psGetStatValue(myStats, &statValue)) {
-        psError(PS_ERR_UNKNOWN, false, "p_psGetStatValue(): could not determine result from requested statistical operation.  Returning in image.\n");
-        psFree(statsAll);
-        return(false);
-    }
-
-    //
-    // Subtract the result and return.
-    //
-    ImageSubtractScalar(in, statValue);
-    psFree(statsAll);
-    return(in);
-}
-
-/******************************************************************************
-ReduceOverscanImageToCol(overscanImage, myStats): This private routine reduces
-a single psImage to a column by combining all pixels from each row into a
-single pixel via requested statistic in myStats.
- *****************************************************************************/
-static psVector *ReduceOverscanImageToCol(
-    psImage *overscanImage,
-    psStats *myStats)
-{
-    psF64 statValue;
-    psVector *tmpRow = psVectorAlloc(overscanImage->numCols, PS_TYPE_F32);
-    psVector *tmpCol = psVectorAlloc(overscanImage->numRows, PS_TYPE_F32);
-
-    //
-    // For each row, we store all pixels in that row into a temporary psVector,
-    // then we run psVectorStats() on that vector.
-    //
-    for (psS32 i=0;i<overscanImage->numRows;i++) {
-        for (psS32 j=0;j<overscanImage->numCols;j++) {
-            tmpRow->data.F32[j] = overscanImage->data.F32[i][j];
-        }
-
-        psStats *rc = psVectorStats(myStats, tmpRow, NULL, NULL, 0);
-        if (rc == NULL) {
-            psError(PS_ERR_UNKNOWN, true, "psVectorStats() could not perform requested statistical operation.  Returning in image.\n");
-            return(NULL);
-        }
-
-        if (false ==  p_psGetStatValue(rc, &statValue)) {
-            psError(PS_ERR_UNKNOWN, true, "p_psGetStatValue() could not determine result from requested statistical operation.  Returning in image.\n");
-            return(NULL);
-        }
-
-        tmpCol->data.F32[i] = (psF32) statValue;
-    }
-    psFree(tmpRow);
-
-    return(tmpCol);
-}
-
-/******************************************************************************
-ReduceOverscanImageToCol(overscanImage, myStats): This private routine reduces
-a single psImage to a row by combining all pixels from each column into a
-single pixel via requested statistic in myStats.
- *****************************************************************************/
-static psVector *ReduceOverscanImageToRow(
-    psImage *overscanImage,
-    psStats *myStats)
-{
-    psF64 statValue;
-    psVector *tmpRow = psVectorAlloc(overscanImage->numCols, PS_TYPE_F32);
-    psVector *tmpCol = psVectorAlloc(overscanImage->numRows, PS_TYPE_F32);
-
-    //
-    // For each column, we store all pixels in that column into a temporary psVector,
-    // then we run psVectorStats() on that vector.
-    //
-    for (psS32 i=0;i<overscanImage->numCols;i++) {
-        for (psS32 j=0;j<overscanImage->numRows;j++) {
-            tmpCol->data.F32[j] = overscanImage->data.F32[j][i];
-        }
-
-        psStats *rc = psVectorStats(myStats, tmpCol, NULL, NULL, 0);
-        if (rc == NULL) {
-            psError(PS_ERR_UNKNOWN, true, "psVectorStats() could not perform requested statistical operation.  Returning in image.\n");
-            return(NULL);
-        }
-
-        if (false ==  p_psGetStatValue(rc, &statValue)) {
-            psError(PS_ERR_UNKNOWN, true, "p_psGetStatValue() could not determine result from requested statistical operation.  Returning in image.\n");
-            return(NULL);
-        }
-
-        tmpRow->data.F32[i] = (psF32) statValue;
-    }
-    psFree(tmpCol);
-
-    return(tmpRow);
-}
-
-/******************************************************************************
-OverscanReduce(vecSize, bias, myStats): This private routine takes a psList of
-overscan images (in bias) and reduces them to a single psVector via the
-specified psStats struct.  The vector is then scaled to the length or the
-row/column in inImg.
- *****************************************************************************/
-static psVector* OverscanReduce(
-    psImage *inImg,
-    pmOverscanAxis overScanAxis,
-    psList *bias,
-    void *fitSpec,
-    pmFit fit,
-    psStats *myStats)
-{
-    if ((overScanAxis != PM_OVERSCAN_ROWS) && (overScanAxis != PM_OVERSCAN_COLUMNS)) {
-        psError(PS_ERR_UNKNOWN, true, "overScanAxis must be PM_OVERSCAN_ROWS or PM_OVERSCAN_COLUMNS\n");
-        return(NULL);
-    }
-    PS_ASSERT_PTR_NON_NULL(inImg, NULL);
-    PS_ASSERT_PTR_NON_NULL(bias, NULL);
-    PS_ASSERT_PTR_NON_NULL(bias->head, NULL);
-    PS_ASSERT_PTR_NON_NULL(myStats, NULL);
-    //
-    // Allocate a psVector for the output of this routine.
-    //
-    psS32 vecSize = GetOverscanSize(inImg, overScanAxis);
-    psVector *overscanVector = psVectorAlloc(vecSize, PS_TYPE_F32);
-
-    //
-    // Allocate an array of psVectors with one psVector per element of the
-    // final oversan column vector.  These psVectors will be used with
-    // psStats to reduce the multiple elements from each overscan column
-    // vector to a single final column vector.
-    //
-    psS32 numOverscanImages = psListLength(bias);
-    psVector **overscanVectors = (psVector **) psAlloc(numOverscanImages * sizeof(psVector *));
-    for (psS32 i = 0 ; i < numOverscanImages ; i++) {
-        overscanVectors[i] = NULL;
-    }
-
-    //
-    // We iterate through the list of overscan images.  For each image,
-    // we reduce it to a single column or row.  Save the overscan vector
-    // in overscanVectors[].
-    //
-    psListElem *tmpOverscan = (psListElem *) bias->head;
-    psS32 overscanID = 0;
-    while (tmpOverscan != NULL) {
-        psImage *tmpOverscanImage = (psImage *) tmpOverscan->data;
-        if (overScanAxis == PM_OVERSCAN_ROWS) {
-            overscanVectors[overscanID] = ReduceOverscanImageToRow(tmpOverscanImage, myStats);
-        } else if (overScanAxis == PM_OVERSCAN_COLUMNS) {
-            overscanVectors[overscanID] = ReduceOverscanImageToCol(tmpOverscanImage, myStats);
-        }
-
-        tmpOverscan = tmpOverscan->next;
-        overscanID++;
-    }
-
-    //
-    // For each overscan vector, if necessary, we scale that column or
-    // row to vecSize.  Note: we should have already ensured that the
-    // fit is poly or spline.
-    //
-    for (psS32 i = 0 ; i < numOverscanImages ; i++) {
-        psVector *tmpOverscanVector = overscanVectors[i];
-
-        if (tmpOverscanVector->n != vecSize) {
-            overscanVectors[i] = ScaleOverscanVector(tmpOverscanVector, vecSize, fitSpec, fit);
-            if (overscanVectors[i] == NULL) {
-                psError(PS_ERR_UNKNOWN, false, "ScaleOverscanVector(): could not scale the overscan vector.\n");
-                for (psS32 i = 0 ; i < numOverscanImages ; i++) {
-                    psFree(overscanVectors[i]);
-                }
-                psFree(overscanVectors);
-                psFree(tmpOverscanVector);
-                return(NULL);
-            }
-            psFree(tmpOverscanVector);
-        }
-    }
-
-    //
-    // We collect all elements in the overscan vectors for the various
-    // overscan images into a single psVector (tmpVec).  Then we call
-    // psStats on that vector to determine the final values for the
-    // overscan vector.
-    //
-    psVector *tmpVec = psVectorAlloc(numOverscanImages, PS_TYPE_F32);
-    psF64 statValue;
-    for (psS32 i = 0 ; i < vecSize ; i++) {
-        // Collect the i-th elements from each overscan vector into a single vector.
-        for (psS32 j = 0 ; j < numOverscanImages ; j++) {
-            tmpVec->data.F32[j] = overscanVectors[j]->data.F32[i];
-        }
-
-        if (NULL == psVectorStats(myStats, tmpVec, NULL, NULL, 0)) {
-            psError(PS_ERR_UNKNOWN, false, "psVectorStats(): could not perform requested statistical operation.  Returning in image.\n");
-            for (psS32 i = 0 ; i < numOverscanImages ; i++) {
-                psFree(overscanVectors[i]);
-            }
-            psFree(overscanVectors);
-            psFree(tmpVec);
-            return(NULL);
-        }
-        if (false == p_psGetStatValue(myStats, &statValue)) {
-            psError(PS_ERR_UNKNOWN, false, "p_psGetStatValue(): could not determine result from requested statistical operation.  Returning in image.\n");
-            for (psS32 i = 0 ; i < numOverscanImages ; i++) {
-                psFree(overscanVectors[i]);
-            }
-            psFree(overscanVectors);
-            psFree(tmpVec);
-            return(NULL);
-        }
-
-        overscanVector->data.F32[i] = (psF32) statValue;
-    }
-
-    //
-    // We're done.  Free the intermediate overscan vectors.
-    //
-    psFree(tmpVec);
-    for (psS32 i = 0 ; i < numOverscanImages ; i++) {
-        psFree(overscanVectors[i]);
-    }
-    psFree(overscanVectors);
-
-    //
-    // Return the computed overscanVector
-    //
-    return(overscanVector);
-}
-
-/******************************************************************************
-RebinOverscanVector(overscanVector, nBinOrig, myStats): this private routine
-takes groups of nBinOrig elements in the input vector, combines them into a
-single pixel via myStats and psVectorStats(), and then outputs a vector of
-those pixels.
- *****************************************************************************/
-static psS32 RebinOverscanVector(
-    psVector *overscanVector,
-    psS32 nBinOrig,
-    psStats *myStats)
-{
-    psF64 statValue;
-    psS32 nBin;
-    if ((nBinOrig > 1) && (nBinOrig < overscanVector->n)) {
-        psS32 numBins = 1+((overscanVector->n)/nBinOrig);
-        psVector *myBin = psVectorAlloc(numBins, PS_TYPE_F32);
-        psVector *binVec = psVectorAlloc(nBinOrig, PS_TYPE_F32);
-
-        for (psS32 i=0;i<numBins;i++) {
-            for(psS32 j=0;j<nBinOrig;j++) {
-                if (overscanVector->n > ((i*nBinOrig)+j)) {
-                    binVec->data.F32[j] = overscanVector->data.F32[(i*nBinOrig)+j];
-                } else {
-                    // XXX: we get here if nBinOrig does not evenly divide
-                    // the overscanVector vector.  This is the last bin.  Should
-                    // we change the binVec->n to acknowledge that?
-                    binVec->n = j;
-                }
-            }
-            psStats *rc = psVectorStats(myStats, binVec, NULL, NULL, 0);
-            if (rc == NULL) {
-                psError(PS_ERR_UNKNOWN, false, "psVectorStats(): could not perform requested statistical operation.  Returning in image.\n");
-                return(-1);
-            }
-            if (false ==  p_psGetStatValue(rc, &statValue)) {
-                psError(PS_ERR_UNKNOWN, false, "p_psGetStatValue(): could not determine result from requested statistical operation.  Returning in image.\n");
-                return(-1);
-            }
-            myBin->data.F32[i] = statValue;
-        }
-
-        // Change the effective size of overscanVector.
-        overscanVector->n = numBins;
-        for (psS32 i=0;i<numBins;i++) {
-            overscanVector->data.F32[i] = myBin->data.F32[i];
-        }
-        psFree(binVec);
-        psFree(myBin);
-        nBin = nBinOrig;
-    } else {
-        nBin = 1;
-    }
-
-    return(nBin);
-}
-
-/******************************************************************************
-FitOverscanVectorAndUnbin(inImg, overscanVector, overScanAxis, fitSpec, fit,
-nBin):  this private routine fits a psPolynomial or psSpline to the overscan
-vector.  It then creates a new vector, with a size determined by the input
-image, evaluates the psPolynomial or psSpline at each element in that vector,
-then returns that vector.
- *****************************************************************************/
-static psVector *FitOverscanVectorAndUnbin(
-    psImage *inImg,
-    psVector *overscanVector,
-    pmOverscanAxis overScanAxis,
-    void *fitSpec,
-    pmFit fit,
-    psS32 nBin)
-{
-    psPolynomial1D* myPoly = NULL;
-    psSpline1D *mySpline = NULL;
-    //
-    // Fit a polynomial or spline to the overscan vector.
-    //
-    if (fit == PM_FIT_POLYNOMIAL) {
-        myPoly = (psPolynomial1D *) fitSpec;
-        PS_ASSERT_POLY_NON_NULL(myPoly, NULL);
-        PS_ASSERT_POLY1D(myPoly, NULL);
-        myPoly = psVectorFitPolynomial1D(myPoly, NULL, 0, overscanVector, NULL, NULL);
-        if (myPoly == NULL) {
-            psError(PS_ERR_UNKNOWN, false, "Could not fit a polynomial to overscan vector.  Returning NULL.\n");
-            return(NULL);
-        }
-    } else if (fit == PM_FIT_SPLINE) {
-        mySpline = psVectorFitSpline1D(NULL, overscanVector);
-        if (mySpline == NULL) {
-            psError(PS_ERR_UNKNOWN, false, "Could not fit a spline to overscan vector.  Returning NULL.\n");
-            return(NULL);
-        }
-        if (fitSpec != NULL) {
-            // Copy the resulting spline fit into fitSpec.
-            psSpline1D *ptrSpline = (psSpline1D *) fitSpec;
-            PS_ASSERT_SPLINE(ptrSpline, NULL);
-            SplineCopy(ptrSpline, mySpline);
-        }
-    }
-
-    //
-    // Evaluate the poly/spline at each pixel in the overscan row/column.
-    //
-    psS32 vecSize = GetOverscanSize(inImg, overScanAxis);
-    psVector *newVec = psVectorAlloc(vecSize, PS_TYPE_F32);
-    if ((nBin > 1) && (nBin < overscanVector->n)) {
-        for (psS32 i = 0 ; i < vecSize ; i++) {
-            if (fit == PM_FIT_POLYNOMIAL) {
-                newVec->data.F32[i] = psPolynomial1DEval(myPoly, ((psF32) i) / ((psF32) nBin));
-            } else if (fit == PM_FIT_SPLINE) {
-                newVec->data.F32[i] = psSpline1DEval(mySpline, ((psF32) i) / ((psF32) nBin));
-            }
-        }
-    } else {
-        for (psS32 i = 0 ; i < vecSize ; i++) {
-            if (fit == PM_FIT_POLYNOMIAL) {
-                newVec->data.F32[i] = psPolynomial1DEval(myPoly, (psF32) i);
-            } else if (fit == PM_FIT_SPLINE) {
-                newVec->data.F32[i] = psSpline1DEval(mySpline, (psF32) i);
-            }
-        }
-    }
-
-    psFree(mySpline);
-    psFree(overscanVector);
-    return(newVec);
-}
-
-
-
-/******************************************************************************
-UnbinOverscanVector(inImg, overscanVector, overScanAxis, nBin):  this private
-routine takes a psVector overscanVector that was previously binned by a factor
-of nBin, and then expands it to its original size, duplicated elements nBin
-times for each element in the input vector overscanVector.
- *****************************************************************************/
-static psVector *UnbinOverscanVector(
-    psImage *inImg,
-    psVector *overscanVector,
-    pmOverscanAxis overScanAxis,
-    psS32 nBin)
-{
-    psS32 vecSize;
-
-    if (overScanAxis == PM_OVERSCAN_ROWS) {
-        vecSize = inImg->numCols;
-    } else if (overScanAxis == PM_OVERSCAN_COLUMNS) {
-        vecSize = inImg->numRows;
-    }
-
-    psVector *newVec = psVectorAlloc(vecSize, PS_TYPE_F32);
-    for (psS32 i = 0 ; i < vecSize ; i++) {
-        newVec->data.F32[i] = overscanVector->data.F32[i/nBin];
-    }
-
-    psFree(overscanVector);
-    return(newVec);
-}
-
-
-/******************************************************************************
-SubtractVectorFromImage(inImg, overscanVector, overScanAxis):  this private
-routine subtracts the overscanVector column-wise or row-wise from inImg.
- *****************************************************************************/
-static psImage *SubtractVectorFromImage(
-    psImage *inImg,
-    psVector *overscanVector,
-    pmOverscanAxis overScanAxis)
-{
-    //
-    // Subtract overscan vector row-wise from the image.
-    //
-    if (overScanAxis == PM_OVERSCAN_ROWS) {
-        for (psS32 i=0;i<inImg->numCols;i++) {
-            for (psS32 j=0;j<inImg->numRows;j++) {
-                inImg->data.F32[j][i]-= overscanVector->data.F32[i];
-            }
-        }
-    }
-
-    //
-    // Subtract overscan vector column-wise from the image.
-    //
-    if (overScanAxis == PM_OVERSCAN_COLUMNS) {
-        for (psS32 i=0;i<inImg->numRows;i++) {
-            for (psS32 j=0;j<inImg->numCols;j++) {
-                inImg->data.F32[i][j]-= overscanVector->data.F32[i];
-            }
-        }
-    }
-
-    return(inImg);
-}
-
-
-
-typedef enum {
-    PM_ERROR_NO_SUBTRACTION,
-    PM_WARNING_NO_SUBTRACTION,
-    PM_ERROR_NO_BIAS_SUBTRACT,
-    PM_WARNING_NO_BIAS_SUBTRACT,
-    PM_ERROR_NO_DARK_SUBTRACT,
-    PM_WARNING_NO_DARK_SUBTRACT,
-    PM_OKAY
-} pmSubtractBiasAssertStatus;
-/******************************************************************************
-AssertCodeOverscan(....) this private routine verifies that the various input
-parameters to pmSubtractBias() are correct for overscan subtraction.
- *****************************************************************************/
-pmSubtractBiasAssertStatus AssertCodeOverscan(
-    pmReadout *in,
-    void *fitSpec,
-    pmFit fit,
-    bool overscan,
-    psStats *stat,
-    int nBinOrig,
-    const pmReadout *bias,
-    const pmReadout *dark)
-{
-
-    PS_ASSERT_READOUT_NON_NULL(in, PM_ERROR_NO_SUBTRACTION);
-    PS_ASSERT_READOUT_NON_EMPTY(in, PM_ERROR_NO_SUBTRACTION);
-    PS_ASSERT_READOUT_TYPE(in, PS_TYPE_F32, PM_ERROR_NO_SUBTRACTION);
-    PS_WARN_PTR_NON_NULL(in->parent);
-    if (in->parent != NULL) {
-        PS_WARN_PTR_NON_NULL(in->parent->concepts);
-    }
-
-    if (overscan == true) {
-        pmOverscanAxis overScanAxis = GetOverscanAxis(in);
-        PS_ASSERT_PTR_NON_NULL(stat, PM_ERROR_NO_SUBTRACTION);
-        PS_ASSERT_PTR_NON_NULL(in->bias, PM_ERROR_NO_SUBTRACTION);
-        PS_ASSERT_PTR_NON_NULL(in->bias->head, PM_ERROR_NO_SUBTRACTION);
-        //
-        // Check the type, size of each bias image.
-        //
-        psListElem *tmpOverscan = (psListElem *) in->bias->head;
-        psS32 numOverscans = 0;
-        while (NULL != tmpOverscan) {
-            numOverscans++;
-            psImage *myOverscanImage = (psImage *) tmpOverscan->data;
-            PS_ASSERT_IMAGE_TYPE(myOverscanImage, PS_TYPE_F32, PM_ERROR_NO_SUBTRACTION);
-            // XXX: Get this right with the rows and columns.
-            if (overScanAxis == PM_OVERSCAN_ROWS) {
-                if (myOverscanImage->numRows != in->image->numRows) {
-                    psLogMsg(__func__, PS_LOG_WARN,
-                             "WARNING: pmSubtractBias.(): overscan image (# %d) has %d rows, input image has %d rows\n",
-                             numOverscans, myOverscanImage->numCols, in->image->numRows);
-                    if (fit == PM_FIT_NONE) {
-                        psError(PS_ERR_UNKNOWN, true, "Don't know how to scale the overscan vectors.  Set fit to PM_FIT_POLYNOMIAL or PM_FIT_SPLINE.\n");
-                        return(PM_ERROR_NO_SUBTRACTION);
-                    }
-                }
-            } else if (overScanAxis == PM_OVERSCAN_COLUMNS) {
-                if (myOverscanImage->numCols != in->image->numCols) {
-                    psLogMsg(__func__, PS_LOG_WARN,
-                             "WARNING: pmSubtractBias.(): overscan image (# %d) has %d columns, input image has %d columns\n",
-                             numOverscans, myOverscanImage->numCols, in->image->numCols);
-                    if (fit == PM_FIT_NONE) {
-                        psError(PS_ERR_UNKNOWN, true, "Don't know how to scale the overscan vectors.  Set fit to PM_FIT_POLYNOMIAL or PM_FIT_SPLINE.\n");
-                        return(PM_ERROR_NO_SUBTRACTION);
-                    }
-                }
-            } else if (overScanAxis != PM_OVERSCAN_ALL) {
-                psError(PS_ERR_UNKNOWN, true, "Must specify and overscan axis.\n");
-                return(PM_ERROR_NO_SUBTRACTION);
-            }
-            tmpOverscan = tmpOverscan->next;
-        }
-    } else {
-        if (fit != PM_FIT_NONE) {
-            psLogMsg(__func__, PS_LOG_WARN,
-                     "WARNING: pmSubtractBias.(): overscan is FALSE and fit is not PM_FIT_NONE.\n");
-            return(PM_WARNING_NO_SUBTRACTION);
-        }
-    }
-
-    // XXX: I do not like the following spec since it's useless to specify
-    // a psSpline as the fitSpec.
-    if (0) {
-        if ((fitSpec == NULL) &&
-                ((fit != PM_FIT_NONE) || (overscan == true))) {
-            psError(PS_ERR_UNKNOWN, true, "fitSpec is NULL and fit is not PM_FIT_NONE or overscan is TRUE.\n");
-            return(PM_ERROR_NO_SUBTRACTION);
-        }
-    }
-
-    return(PM_OKAY);
-}
-
-/******************************************************************************
-AssertCodeBias(....) this private routine verifies that the various input
-parameters to pmSubtractBias() are correct for bias subtraction.
- *****************************************************************************/
-static pmSubtractBiasAssertStatus AssertCodeBias(
-    pmReadout *in,
-    void *fitSpec,
-    pmFit fit,
-    bool overscan,
-    psStats *stat,
-    int nBinOrig,
-    const pmReadout *bias,
-    const pmReadout *dark)
-{
-    if ((in->image->numRows + in->row0 - bias->row0) > bias->image->numRows) {
-        psError(PS_ERR_UNKNOWN,true, "bias image does not have enough rows.  Returning in image\n");
-        return(PM_ERROR_NO_BIAS_SUBTRACT);
-    }
-    if ((in->image->numCols + in->col0 - bias->col0) > bias->image->numCols) {
-        psError(PS_ERR_UNKNOWN,true, "bias image does not have enough columns.  Returning in image\n");
-        return(PM_ERROR_NO_BIAS_SUBTRACT);
-    }
-
-    if (bias != NULL) {
-        PS_ASSERT_READOUT_NON_EMPTY(bias, PM_ERROR_NO_BIAS_SUBTRACT);
-        PS_ASSERT_READOUT_TYPE(bias, PS_TYPE_F32, PM_ERROR_NO_DARK_SUBTRACT);
-    }
-    return(PM_OKAY);
-}
-
-/******************************************************************************
-AssertCodeDark(....) this private routine verifies that the various input
-parameters to pmSubtractBias() are correct for dark subtraction.
- *****************************************************************************/
-pmSubtractBiasAssertStatus AssertCodeDark(
-    pmReadout *in,
-    void *fitSpec,
-    pmFit fit,
-    bool overscan,
-    psStats *stat,
-    int nBinOrig,
-    const pmReadout *bias,
-    const pmReadout *dark)
-{
-    if ((in->image->numRows + in->row0 - dark->row0) > dark->image->numRows) {
-        psError(PS_ERR_UNKNOWN, true, "dark image does not have enough rows.  Returning in image\n");
-        return(PM_ERROR_NO_DARK_SUBTRACT);
-    }
-    if ((in->image->numCols + in->col0 - dark->col0) > dark->image->numCols) {
-        psError(PS_ERR_UNKNOWN, true, "dark image does not have enough columns.  Returning in image\n");
-        return(PM_ERROR_NO_DARK_SUBTRACT);
-    }
-
-    if (dark != NULL) {
-        PS_ASSERT_READOUT_NON_EMPTY(dark, PM_ERROR_NO_DARK_SUBTRACT);
-        PS_ASSERT_READOUT_TYPE(dark, PS_TYPE_F32, PM_ERROR_NO_DARK_SUBTRACT);
-    }
-    return(PM_OKAY);
-}
-
-/******************************************************************************
-p_psDetermineTrimmedImage(): global routine: determines the region of the
-input pmReadout which will be operated on by the various detrend modules.  It
-does a metadata fetch on "CELL.TRIMSEC" for the parent cell of the pmReadout.
- 
-Use it this way:
-    PS_WARN_PTR_NON_NULL(in->parent);
-    if (in->parent != NULL) {
-        PS_WARN_PTR_NON_NULL(in->parent->concepts);
-    }
-    //
-    // Determine trimmed image from metadata.
-    //
-    psImage *trimmedImg = p_psDetermineTrimmedImage(in);
- 
-XXX: Create a pmUtils.c file and put this routine there.
- *****************************************************************************/
-psImage *p_psDetermineTrimmedImage(pmReadout *in)
-{
-    if ((in->parent == NULL) || (in->parent->concepts == NULL)) {
-        psLogMsg(__func__, PS_LOG_WARN,
-                 "WARNING: could not determine CELL.TRIMSEC from parent cell Metadata (NULL).\n");
-        return(in->image);
-    }
-
-    psBool rc = false;
-    psImage *trimmedImg = NULL;
-    psRegion *trimRegion = psMetadataLookupPtr(&rc, in->parent->concepts,
-                           "CELL.TRIMSEC");
-    if (rc == false) {
-        psLogMsg(__func__, PS_LOG_WARN,
-                 "WARNING: could not determine CELL.TRIMSEC from parent cell Metadata.\n");
-        trimmedImg = in->image;
-    } else {
-        trimmedImg = psImageSubset(in->image, *trimRegion);
-    }
-
-    return(trimmedImg);
-}
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-/******************************************************************************
-pmSubtractBias(....): see SDRS for complete specification.
- 
-XXX: Code and assert type support: U16, S32, F32.
-XXX: Add trace messages.
- *****************************************************************************/
-pmReadout *pmSubtractBias(
-    pmReadout *in,
-    void *fitSpec,
-    pmFit fit,
-    bool overscan,
-    psStats *stat,
-    int nBin,
-    const pmReadout *bias,
-    const pmReadout *dark)
+pmReadout *pmSubtractBias(pmReadout *in,
+                          void *fitSpec,
+                          const psList *overscans,
+                          pmOverscanAxis overScanAxis,
+                          psStats *stat,
+                          psS32 nBinOrig,
+                          pmFit fit,
+                          const pmReadout *bias)
 {
     psTrace(".psModule.pmSubtracBias.pmSubtractBias", 4,
             "---- pmSubtractBias() begin ----\n");
-    //
-    // Check input parameters, generate warnings and errors.
-    //
-    if (PM_OKAY != AssertCodeOverscan(in, fitSpec, fit, overscan, stat, nBin, bias, dark)) {
-        return(in);
-    }
-    //
-    // Determine trimmed image from metadata.
-    //
-    psImage *trimmedImg = p_psDetermineTrimmedImage(in);
-
-    //
-    // Subtract overscan frames if necessary.
-    //
-    if (overscan == true) {
-        pmOverscanAxis overScanAxis = GetOverscanAxis(in);
+    PS_ASSERT_READOUT_NON_NULL(in, NULL);
+    PS_ASSERT_READOUT_NON_EMPTY(in, NULL);
+    PS_ASSERT_READOUT_TYPE(in, PS_TYPE_F32, NULL);
+
+    //
+    // If the overscans != NULL, then check the type of each image.
+    //
+    if (overscans != NULL) {
+        psListElem *tmpOverscan = (psListElem *) overscans->head;
+        while (NULL != tmpOverscan) {
+            psImage *myOverscanImage = (psImage *) tmpOverscan->data;
+            PS_ASSERT_IMAGE_TYPE(myOverscanImage, PS_TYPE_F32, NULL);
+            tmpOverscan = tmpOverscan->next;
+        }
+    }
+
+    if ((overscans == NULL) && (overScanAxis != PM_OVERSCAN_NONE)) {
+        psError(PS_ERR_UNKNOWN,true, "(overscans == NULL) && (overScanAxis != PM_OVERSCAN_NONE).  Returning in image\n");
+        return(in);
+    }
+
+    // Check for an unallowable pmFit.
+    if ((fit != PM_OVERSCAN_NONE) &&
+            (fit != PM_OVERSCAN_ROWS) &&
+            (fit != PM_OVERSCAN_COLUMNS) &&
+            (fit != PM_OVERSCAN_ALL)) {
+        psError(PS_ERR_UNKNOWN, true, "fit is unallowable (%d).  Returning in image.\n", fit);
+        return(in);
+    }
+    // Check for an unallowable pmOverscanAxis.
+    if ((overScanAxis != PM_OVERSCAN_NONE) &&
+            (overScanAxis != PM_OVERSCAN_ROWS) &&
+            (overScanAxis != PM_OVERSCAN_COLUMNS) &&
+            (overScanAxis != PM_OVERSCAN_ALL)) {
+        psError(PS_ERR_UNKNOWN, true, "overScanAxis is unallowable (%d).  Returning in image.\n", overScanAxis);
+        return(in);
+    }
+    psS32 i;
+    psS32 j;
+    psS32 numBins = 0;
+    static psVector *overscanVector = NULL;
+    psVector *tmpRow = NULL;
+    psVector *tmpCol = NULL;
+    psVector *myBin = NULL;
+    psVector *binVec = NULL;
+    psListElem *tmpOverscan = NULL;
+    double statValue;
+    psImage *myOverscanImage = NULL;
+    psPolynomial1D *myPoly = NULL;
+    psSpline1D *mySpline = NULL;
+    psS32 nBin;
+
+    //
+    //  Create a static stats data structure and determine the highest
+    //  priority stats option.
+    //
+    static psStats *myStats = NULL;
+    if (myStats == NULL) {
+        myStats = psStatsAlloc(PS_STAT_SAMPLE_MEAN);
+        p_psMemSetPersistent(myStats, true);
+    }
+    if (stat != NULL) {
+        myStats->options = GenNewStatOptions(stat);
+    }
+
+
+    if (overScanAxis == PM_OVERSCAN_NONE) {
+        if (fit != PM_FIT_NONE) {
+            psLogMsg(__func__, PS_LOG_WARN,
+                     "WARNING: pmSubtractBias.(): overScanAxis equals NONE, and fit does not equal NONE.  Proceeding to full fram subtraction.\n");
+        }
+
+        if (overscans != NULL) {
+            psLogMsg(__func__, PS_LOG_WARN,
+                     "WARNING: pmSubtractBias.(): overScanAxis equals NONE and overscans does not equal NULL.  Proceeding to full fram subtraction.\n");
+        }
+        return(SubtractFrame(in, bias));
+    }
+
+    if ((overScanAxis == PM_OVERSCAN_ALL) && (fit != PM_FIT_NONE)) {
+        psLogMsg(__func__, PS_LOG_WARN,
+                 "WARNING: pmSubtractBias.(): overScanAxis equals ALL, and fit does not equal NONE.  Proceeding with the rest of the module.\n");
+    }
+
+
+    //
+    // We subtract each overscan region from the image data.
+    // If we get here we know that overscans != NULL.
+    //
+
+    if (overScanAxis == PM_OVERSCAN_ALL) {
+        tmpOverscan = (psListElem *) overscans->head;
+        while (NULL != tmpOverscan) {
+            myOverscanImage = (psImage *) tmpOverscan->data;
+
+            PS_ASSERT_IMAGE_TYPE(myOverscanImage, PS_TYPE_F32, NULL);
+            psStats *rc = psImageStats(myStats, myOverscanImage, NULL, (psMaskType)0xffffffff);
+            if (rc == NULL) {
+                psError(PS_ERR_UNKNOWN, false, "psImageStats(): could not perform requested statistical operation.  Returning in image.\n");
+                return(in);
+            }
+            if (false == p_psGetStatValue(myStats, &statValue)) {
+                psError(PS_ERR_UNKNOWN, false, "p_psGetStatValue(): could not determine result from requested statistical operation.  Returning in image.\n");
+                return(in);
+            }
+            ImageSubtractScalar(in->image, statValue);
+
+            tmpOverscan = tmpOverscan->next;
+        }
+        return(in);
+    }
+
+    // This check is redundant with above code.
+    if (!((overScanAxis == PM_OVERSCAN_ROWS) || (overScanAxis == PM_OVERSCAN_COLUMNS))) {
+        psError(PS_ERR_UNKNOWN, true, "overScanAxis is unallowable (%d).\nReturning in image.\n", overScanAxis);
+        return(in);
+    }
+
+    tmpOverscan = (psListElem *) overscans->head;
+    while (NULL != tmpOverscan) {
+        //        PS_IMAGE_PRINT_F32_HIDEF(in->image);
+        myOverscanImage = (psImage *) tmpOverscan->data;
+
+        if (overScanAxis == PM_OVERSCAN_ROWS) {
+            if (myOverscanImage->numCols != (in->image)->numCols) {
+                psLogMsg(__func__, PS_LOG_WARN,
+                         "WARNING: pmSubtractBias.(): overscan image has %d columns, input image has %d columns\n",
+                         myOverscanImage->numCols, in->image->numCols);
+            }
+
+            // We create a row vector and subtract this vector from image.
+            // XXX: Is there a better way to extract a psVector from a psImage without
+            // having to copy every element in that vector?
+            overscanVector = psVectorAlloc(myOverscanImage->numCols, PS_TYPE_F32);
+            for (i=0;i<overscanVector->n;i++) {
+                overscanVector->data.F32[i] = 0.0;
+            }
+            tmpRow = psVectorAlloc(myOverscanImage->numRows, PS_TYPE_F32);
+
+            // For each column of the input image, loop through every row,
+            // collect the pixel in that row, then performed the specified
+            // statistical op on those pixels.  Store this in overscanVector.
+            for (i=0;i<myOverscanImage->numCols;i++) {
+                for (j=0;j<myOverscanImage->numRows;j++) {
+                    tmpRow->data.F32[j] = myOverscanImage->data.F32[j][i];
+                }
+                psStats *rc = psVectorStats(myStats, tmpRow, NULL, NULL, 0);
+                if (rc == NULL) {
+                    psError(PS_ERR_UNKNOWN, false, "psVectorStats(): could not perform requested statistical operation.  Returning in image.\n");
+                    return(in);
+                }
+                if (false ==  p_psGetStatValue(rc, &statValue)) {
+                    psError(PS_ERR_UNKNOWN, false, "p_psGetStatValue(): could not determine result from requested statistical operation.  Returning in image.\n");
+                    return(in);
+                }
+                overscanVector->data.F32[i] = statValue;
+            }
+            psFree(tmpRow);
+
+            // Scale the overscan vector to the size of the input image.
+            if (overscanVector->n != in->image->numCols) {
+                if ((fit == PM_FIT_POLYNOMIAL) || (fit == PM_FIT_SPLINE)) {
+                    psVector *newVec = ScaleOverscanVector(overscanVector,
+                                                           in->image->numCols,
+                                                           fitSpec, fit);
+                    if (newVec == NULL) {
+                        psError(PS_ERR_UNKNOWN, false, "ScaleOverscanVector(): could not scale the overscan vector.  Returning in image.\n");
+                        return(in);
+                    }
+                    psFree(overscanVector);
+                    overscanVector = newVec;
+                } else {
+                    psError(PS_ERR_UNKNOWN, true, "Don't know how to scale the overscan vector.  Set fit to PM_FIT_SPLINE or PM_FIT_POLYNOMIAL.  Returning in image.\n");
+                    psFree(overscanVector);
+                    return(in);
+                }
+            }
+        }
+
+        if (overScanAxis == PM_OVERSCAN_COLUMNS) {
+            if (myOverscanImage->numRows != (in->image)->numRows) {
+                psLogMsg(__func__, PS_LOG_WARN,
+                         "WARNING: pmSubtractBias.(): overscan image has %d rows, input image has %d rows\n",
+                         myOverscanImage->numRows, in->image->numRows);
+            }
+
+            // We create a column vector and subtract this vector from image.
+            overscanVector = psVectorAlloc(myOverscanImage->numRows, PS_TYPE_F32);
+            for (i=0;i<overscanVector->n;i++) {
+                overscanVector->data.F32[i] = 0.0;
+            }
+            tmpCol = psVectorAlloc(myOverscanImage->numCols, PS_TYPE_F32);
+
+            // For each row of the input image, loop through every column,
+            // collect the pixel in that row, then performed the specified
+            // statistical op on those pixels.  Store this in overscanVector.
+            for (i=0;i<myOverscanImage->numRows;i++) {
+                for (j=0;j<myOverscanImage->numCols;j++) {
+                    tmpCol->data.F32[j] = myOverscanImage->data.F32[i][j];
+                }
+                psStats *rc = psVectorStats(myStats, tmpCol, NULL, NULL, 0);
+                if (rc == NULL) {
+                    psError(PS_ERR_UNKNOWN, false, "psVectorStats(): could not perform requested statistical operation.  Returning in image.\n");
+                    return(in);
+                }
+                if (false ==  p_psGetStatValue(rc, &statValue)) {
+                    psError(PS_ERR_UNKNOWN, false, "p_psGetStatValue(): could not determine result from requested statistical operation.  Returning in image.\n");
+                    return(in);
+                }
+                overscanVector->data.F32[i] = statValue;
+            }
+            psFree(tmpCol);
+
+            // Scale the overscan vector to the size of the input image.
+            if (overscanVector->n != in->image->numRows) {
+                if ((fit == PM_FIT_POLYNOMIAL) || (fit == PM_FIT_SPLINE)) {
+                    psVector *newVec = ScaleOverscanVector(overscanVector,
+                                                           in->image->numRows,
+                                                           fitSpec, fit);
+                    if (newVec == NULL) {
+                        psError(PS_ERR_UNKNOWN, false, "ScaleOverscanVector(): could not scale the overscan vector.  Returning in image.\n");
+                        return(in);
+                    }
+                    psFree(overscanVector);
+                    overscanVector = newVec;
+                } else {
+                    psError(PS_ERR_UNKNOWN, true, "Don't know how to scale the overscan vector.  Set fit to PM_FIT_SPLINE or PM_FIT_POLYNOMIAL.  Returning in image.\n");
+                    psFree(overscanVector);
+                    return(in);
+                }
+            }
+        }
+
         //
-        //  Create a psStats data structure and determine the highest
-        //  priority stats option.
+        // Re-bin the overscan vector (change its length).
         //
-        psStats *myStats = psStatsAlloc(PS_STAT_SAMPLE_MEAN);
-        if (stat != NULL) {
-            myStats->options = GenNewStatOptions(stat);
-        }
-
+        // Only if nBinOrig > 1.
+        if ((nBinOrig > 1) && (nBinOrig < overscanVector->n)) {
+            numBins = 1+((overscanVector->n)/nBinOrig);
+            myBin = psVectorAlloc(numBins, PS_TYPE_F32);
+            binVec = psVectorAlloc(nBinOrig, PS_TYPE_F32);
+
+            for (i=0;i<numBins;i++) {
+                for(j=0;j<nBinOrig;j++) {
+                    if (overscanVector->n > ((i*nBinOrig)+j)) {
+                        binVec->data.F32[j] = overscanVector->data.F32[(i*nBinOrig)+j];
+                    } else {
+                        // XXX: we get here if nBinOrig does not evenly divide
+                        // the overscanVector vector.  This is the last bin.  Should
+                        // we change the binVec->n to acknowledge that?
+                        binVec->n = j;
+                    }
+                }
+                psStats *rc = psVectorStats(myStats, binVec, NULL, NULL, 0);
+                if (rc == NULL) {
+                    psError(PS_ERR_UNKNOWN, false, "psVectorStats(): could not perform requested statistical operation.  Returning in image.\n");
+                    return(in);
+                }
+                if (false ==  p_psGetStatValue(rc, &statValue)) {
+                    psError(PS_ERR_UNKNOWN, false, "p_psGetStatValue(): could not determine result from requested statistical operation.  Returning in image.\n");
+                    return(in);
+                }
+                myBin->data.F32[i] = statValue;
+            }
+
+            // Change the effective size of overscanVector.
+            overscanVector->n = numBins;
+            for (i=0;i<numBins;i++) {
+                overscanVector->data.F32[i] = myBin->data.F32[i];
+            }
+            psFree(binVec);
+            psFree(myBin);
+            nBin = nBinOrig;
+        } else {
+            nBin = 1;
+        }
+
+        // At this point the number of data points in overscanVector should be
+        // equal to the number of rows/columns (whatever is appropriate) in the
+        // image divided by numBins.
         //
-        // Reduce overscan images to a single pixel, then subtract.
-        // This code is no longer required as of SDRS 12-09.
+
+
         //
-        if (overScanAxis == PM_OVERSCAN_ALL) {
-            if (false == OverscanReducePixel(trimmedImg, in->bias, myStats)) {
-                return(in);
-            }
-            psFree(myStats);
+        // This doesn't seem right.  The only way to do a spline fit is if,
+        // by SDRS requirements, fitSpec is not-NULL>  But in order for it
+        // to be non-NULL, someone must have called psSpline1DAlloc() with
+        // the min, max, and number of splines.
+        //
+        if (!((fitSpec == NULL) || (fit == PM_FIT_NONE))) {
+            //
+            // Fit a polynomial or spline to the overscan vector.
+            //
+            if (fit == PM_FIT_POLYNOMIAL) {
+                myPoly = (psPolynomial1D *) fitSpec;
+                myPoly = psVectorFitPolynomial1D(myPoly, NULL, 0, overscanVector, NULL, NULL);
+                if (myPoly == NULL) {
+                    psError(PS_ERR_UNKNOWN, false, "(3) Could not fit a polynomial to overscan vector.  Returning in image.\n");
+                    psFree(overscanVector);
+                    return(in);
+                }
+            } else if (fit == PM_FIT_SPLINE) {
+                // XXX: This makes no sense
+                // XXX: must free mySpline?
+                mySpline = (psSpline1D *) fitSpec;
+                mySpline = psVectorFitSpline1D(NULL, overscanVector);
+                if (mySpline == NULL) {
+                    psError(PS_ERR_UNKNOWN, false, "Could not fit a spline to overscan vector.  Returning in image.\n");
+                    psFree(overscanVector);
+                    return(in);
+                }
+            }
+
+            //
+            // Subtract fitted overscan vector row-wise from the image.
+            //
+            if (overScanAxis == PM_OVERSCAN_ROWS) {
+                for (i=0;i<(in->image)->numCols;i++) {
+                    psF32 tmpF32 = 0.0;
+                    if (fit == PM_FIT_POLYNOMIAL) {
+                        tmpF32 = psPolynomial1DEval(myPoly, ((psF32) i) / ((psF32) nBin));
+                    } else if (fit == PM_FIT_SPLINE) {
+                        tmpF32 = psSpline1DEval(mySpline, ((psF32) i) / ((psF32) nBin));
+                    }
+                    for (j=0;j<(in->image)->numRows;j++) {
+                        (in->image)->data.F32[j][i]-= tmpF32;
+                    }
+                }
+            }
+
+            //
+            // Subtract fitted overscan vector column-wise from the image.
+            //
+            if (overScanAxis == PM_OVERSCAN_COLUMNS) {
+                for (i=0;i<(in->image)->numRows;i++) {
+                    psF32 tmpF32 = 0.0;
+                    if (fit == PM_FIT_POLYNOMIAL) {
+                        tmpF32 = psPolynomial1DEval(myPoly, ((psF32) i) / ((psF32) nBin));
+                    } else if (fit == PM_FIT_SPLINE) {
+                        tmpF32 = psSpline1DEval(mySpline, ((psF32) i) / ((psF32) nBin));
+                    }
+
+                    for (j=0;j<(in->image)->numCols;j++) {
+                        (in->image)->data.F32[i][j]-= tmpF32;
+                    }
+                }
+            }
         } else {
             //
-            // Reduce the overscan images to a single overscan vector.
-            //
-            psVector *overscanVector = OverscanReduce(in->image, overScanAxis,
-                                       in->bias, fitSpec,
-                                       fit, myStats);
-            if (overscanVector == NULL) {
-                psError(PS_ERR_UNKNOWN, false, "Could not reduce overscan images to a single overscan vector.  Returning in image\n");
-                psFree(myStats);
-                return(in);
-            }
-
-            //
-            // Rebin the overscan vector if necessary.
-            //
-            psS32 newBin = RebinOverscanVector(overscanVector, nBin, myStats);
-            if (newBin < 0) {
-                psError(PS_ERR_UNKNOWN, false, "Could rebin the overscan vector.  Returning in image\n");
-                psFree(myStats);
-                return(in);
-            }
-
-            //
-            // If necessary, fit a psPolynomial or psSpline to the overscan vector.
-            // Then, unbin the overscan vector to appropriate length for the in image.
-            //
-            if ((fit == PM_FIT_POLYNOMIAL) || (fit == PM_FIT_SPLINE)) {
-                overscanVector = FitOverscanVectorAndUnbin(trimmedImg, overscanVector, overScanAxis, fitSpec, fit, newBin);
-                if (overscanVector == NULL) {
-                    psError(PS_ERR_UNKNOWN, false, "Could not fit the polynomial or spline to the overscan vector.  Returning in image\n");
-                    psFree(myStats);
-                    return(in);
-                }
-            } else {
-                overscanVector = UnbinOverscanVector(trimmedImg, overscanVector, overScanAxis, newBin);
-            }
-
-            //
-            // Subtract the overscan vector from the input image.
-            //
-            SubtractVectorFromImage(trimmedImg, overscanVector, overScanAxis);
-            psFree(myStats);
-            psFree(overscanVector);
-        }
-    }
-
-    //
-    // Perform bias subtraction if necessary.
-    //
-    if (bias != NULL) {
-        if (PM_OKAY == AssertCodeBias(in, fitSpec, fit, overscan, stat, nBin, bias, dark)) {
-            SubtractFrame(in, bias);
-        }
-    }
-
-    //
-    // Perform dark subtraction if necessary.
-    //
-    if (dark != NULL) {
-        if (PM_OKAY == AssertCodeDark(in, fitSpec, fit, overscan, stat, nBin, bias, dark)) {
-            psBool rc;
-            psF32 scale = 0.0;
-            if (in->parent != NULL) {
-                scale = psMetadataLookupS32(&rc, in->parent->concepts, "CELL.DARKTIME");
-                if (rc == false) {
-                    psLogMsg(__func__, PS_LOG_WARN,
-                             "WARNING: pmSubtractBias.(): could not determine CELL.FARKTIME from in->parent metadata.\n");
-                }
-            }
-            SubtractDarkFrame(in, dark, scale);
-        }
-    }
-
-    //
-    // All done.
-    //
+            // If we get here, then no polynomials were fit to the overscan
+            // vector.  We simply subtract it, taking into account binning,
+            // from the image.
+            //
+
+            //
+            // Subtract overscan vector row-wise from the image.
+            //
+            if (overScanAxis == PM_OVERSCAN_ROWS) {
+                for (i=0;i<(in->image)->numCols;i++) {
+                    for (j=0;j<(in->image)->numRows;j++) {
+                        (in->image)->data.F32[j][i]-= overscanVector->data.F32[i/nBin];
+                    }
+                }
+            }
+
+            //
+            // Subtract overscan vector column-wise from the image.
+            //
+            if (overScanAxis == PM_OVERSCAN_COLUMNS) {
+                for (i=0;i<(in->image)->numRows;i++) {
+                    for (j=0;j<(in->image)->numCols;j++) {
+                        (in->image)->data.F32[i][j]-= overscanVector->data.F32[i/nBin];
+                    }
+                }
+            }
+        }
+
+        psFree(overscanVector);
+
+        tmpOverscan = tmpOverscan->next;
+    }
+
     psTrace(".psModule.pmSubtracBias.pmSubtractBias", 4,
             "---- pmSubtractBias() exit ----\n");
+
+    if (bias != NULL) {
+        return(SubtractFrame(in, bias));
+    }
     return(in);
 }

branches/pap_branch_051214/psModules/src/imsubtract/pmSubtractBias.h

@@ +1 @@
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+
 /** @file  pmSubtractBias.h
  *
@@ -6 +11 @@
  *  @author GLG, MHPCC
  *
- *  @version $Revision: 1.4 $ $Name: not supported by cvs2svn $
- *  @date $Date: 2005-11-23 23:54:30 $
+ *  @version $Revision: 1.4.8.1 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2005-12-17 03:18:39 $
  *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
@@ -23 +28 @@
 #include<math.h>
 #include "pslib.h"
+
 #include "pmAstrometry.h"
 
@@ -38 +44 @@
 } pmFit;
 
-pmReadout *pmSubtractBias(
-    pmReadout *in,                      ///< The input pmReadout image
-    void *fitSpec,                      ///< A polynomial or spline, defining the fit type.
-    pmFit fit,                          ///< PM_FIT_SPLINE, PM_FIT_POLYNOMIAL, or PM_FIT_NONE
-    bool overscan,
-    psStats *stat,                      ///< The statistic to be used in combining overscan data
-    int nBin,                           ///< The amount of binning to be done image pixels.
-    const pmReadout *bias,              ///< A possibly NULL bias pmReadout which is to be subtracted
-    const pmReadout *dark               ///< A possibly NULL bias pmReadout which is to be subtracted
-);
-
-/******************************************************************************
-DetermineTrimmedImageg() This private routine determines the region of the
-input pmReadout which will be operated on by the various detrend modules.  It
-does a metadata fetch on "CELL.TRIMSEC" for the parent cell of the pmReadout.
- 
-XXX: Consider making a pmUtils.c file and put this routine there.
- *****************************************************************************/
-psImage *p_psDetermineTrimmedImage(
-    pmReadout *in
-);
+pmReadout *pmSubtractBias(pmReadout *in,                ///< The input pmReadout image
+                          void *fitSpec,                ///< A polynomial or spline, defining the fit type.
+                          const psList *overscans,      ///< A psList of overscan images
+                          pmOverscanAxis overScanAxis,  ///< Defines how overscans are applied
+                          psStats *stat,                ///< The statistic to be used in combining overscan data
+                          int nBin,                     ///< The amount of binning to be done image pixels.
+                          pmFit fit,                    ///< PM_FIT_SPLINE, PM_FIT_POLYNOMIAL, or PM_FIT_NONE
+                          const pmReadout *bias);       ///< A possibly NULL bias pmReadout which is to be subtracted
 
 #endif