Changeset 18830

Timestamp:

Jul 31, 2008, 2:01:26 PM (18 years ago)

Author:

eugene

Message:

src/objects/Makefile.am

Location:

trunk/psModules/src

Files:

• camera/pmFPA.c (modified) (2 diffs)
• camera/pmFPA.h (modified) (2 diffs)
• camera/pmFPARead.c (modified) (14 diffs)
• camera/pmReadoutStack.c (modified) (7 diffs)
• camera/pmReadoutStack.h (modified) (2 diffs)
• config/pmConfig.c (modified) (1 diff)
• imcombine/pmReadoutCombine.c (modified) (9 diffs)
• imcombine/pmReadoutCombine.h (modified) (2 diffs)
• objects/Makefile.am (modified) (2 diffs)

trunk/psModules/src/camera/pmFPA.c

@@ -176 +176 @@
     readout->col0 = 0;
     readout->row0 = 0;
-    readout->imageScan = 0;
-    readout->maskScan = 0;
-    readout->weightScan = 0;
+
+    readout->thisImageScan = 0;
+    readout->thisMaskScan = 0;
+    readout->thisWeightScan = 0;
+
+    readout->lastImageScan = 0;
+    readout->lastMaskScan = 0;
+    readout->lastWeightScan = 0;
 }
 
@@ -269 +274 @@
     tmpReadout->col0 = 0;
 
-    tmpReadout->imageScan = 0;
-    tmpReadout->maskScan = 0;
-    tmpReadout->weightScan = 0;
+    tmpReadout->thisImageScan = 0;
+    tmpReadout->thisMaskScan = 0;
+    tmpReadout->thisWeightScan = 0;
+
+    tmpReadout->lastImageScan = 0;
+    tmpReadout->lastMaskScan = 0;
+    tmpReadout->lastWeightScan = 0;
 
     return(tmpReadout);

trunk/psModules/src/camera/pmFPA.h

@@ -6 +6 @@
  * @author Eugene Magnier, IfA
  *
- * @version $Revision: 1.22 $ $Name: not supported by cvs2svn $
- * @date $Date: 2008-06-05 01:31:33 $
+ * @version $Revision: 1.23 $ $Name: not supported by cvs2svn $
+ * @date $Date: 2008-08-01 00:01:25 $
  * Copyright 2005-2006 Institute for Astronomy, University of Hawaii
  */
@@ -127 +127 @@
     bool file_exists;                   ///< Does the file for this readout exist (read case only)?
     bool data_exists;                   ///< Does the data for this readout exist (read case only)?
-    int imageScan, maskScan, weightScan;///< Separate tracking numbers for reading images incrementally
+    int thisImageScan;                  ///< start scan for next/current read
+    int lastImageScan;                  ///< start scan of the last read
+    int thisMaskScan;                   ///< start scan for next/current read
+    int lastMaskScan;                   ///< start scan of the last read
+    int thisWeightScan;                 ///< start scan for next/current read
+    int lastWeightScan;                 ///< start scan of the last read
 } pmReadout;
 

trunk/psModules/src/camera/pmFPARead.c

@@ -45 +45 @@
 
 // Return pointer to appropriate value for tracking scans
-static int *readoutScansByType(pmReadout *readout, // Readout of interest
-                               fpaReadType type // Type of image
+static int readoutGetThisScan(pmReadout *readout, // Readout of interest
+                              fpaReadType type // Type of image
     )
 {
     switch (type) {
       case FPA_READ_TYPE_IMAGE:
-        return &readout->imageScan;
+        return readout->thisImageScan;
       case FPA_READ_TYPE_MASK:
-        return &readout->maskScan;
+        return readout->thisMaskScan;
       case FPA_READ_TYPE_WEIGHT:
-        return &readout->weightScan;
+        return readout->thisWeightScan;
       default:
         psAbort("Unknown read type: %x\n", type);
     }
+}
+
+// Return pointer to appropriate value for tracking scans
+static int readoutSetLastScan(pmReadout *readout, // Readout of interest
+                              fpaReadType type,   // Type of image
+                              int numScans        // requested number of scans
+    )
+{
+    switch (type) {
+      case FPA_READ_TYPE_IMAGE:
+        readout->lastImageScan = readout->thisImageScan + numScans;
+        return readout->lastImageScan;
+      case FPA_READ_TYPE_MASK:
+        readout->lastMaskScan = readout->thisMaskScan + numScans;
+        return readout->lastMaskScan;
+      case FPA_READ_TYPE_WEIGHT:
+        readout->lastWeightScan = readout->thisWeightScan + numScans;
+        return readout->lastWeightScan;
+      default:
+        psAbort("Unknown read type: %x\n", type);
+    }
+    return false;
 }
 
@@ -134 +156 @@
 }
 
-
-// Determine readout scan properties: the next and last scans
-// Requires that cellNumReadouts() has been called before (for header and concepts to have been read)
-// In the process, adjusts the TRIMSEC
-static bool readoutScanProperties(int *next, // Index of next scan
-                                  int *last, // Index of last scan
-                                  pmReadout *readout, // Readout of interest
-                                  int numScans, // Number of scans to read at a time
-                                  fpaReadType type, // Type of image
-                                  pmConfig *config // Configuration
+// Does the current readout, with scans set for a new read, represent any real data, or is it
+// beyond the end?  Requires that cellNumReadouts() has been called before (for header and
+// concepts to have been read) In the process, adjusts the TRIMSEC
+static bool readoutHaveMoreScans(bool *result,   // true : more data to read 
+                             pmReadout *readout, // Readout of interest
+                             int numScans, // Number of scans to read at a time
+                             fpaReadType type, // Type of image
+                             pmConfig *config // Configuration
     )
 {
-    assert(next);
-    assert(last);
+    assert(result);
     assert(readout);
 
-    psImage *image = *readoutImageByType(readout, type); // Appropriate image from readout
+    *result = false;
 
     if (!pmConceptsReadCell(readout->parent, PM_CONCEPT_SOURCE_DEFAULTS | PM_CONCEPT_SOURCE_DATABASE,
@@ -163 +182 @@
     PS_ASSERT_PTR_NON_NULL(cell, false);
     pmHDU *hdu = pmHDUFromCell(cell);   // HDU for data
+
     bool mdok = true;                   // Status of MD lookup
     psRegion *trimsec = psMetadataLookupPtr(&mdok, cell->concepts, "CELL.TRIMSEC"); // Trim sections
@@ -195 +215 @@
     // Calculate the segment offset and upper limit
     if (numScans == 0) {
-        *next = (readdir == 1) ? trimsec->y0 : trimsec->x0;
-    } else {
-        *next = image ? *readoutScansByType(readout, type) + numScans : 0;
-    }
-    *last = (readdir == 1) ? trimsec->y1 : trimsec->x1;
-
+        // in the case of numScans == 0, we never call this funtion unless the data has not yet 
+        // been read.  thus, only if the delta is  should we return false (ie, trimsec defines an empty region)
+        int thisScan = (readdir == 1) ? trimsec->y0 : trimsec->x0;
+        int lastScan = (readdir == 1) ? trimsec->y1 : trimsec->x1;
+        *result = (lastScan > thisScan);
+        return true;
+    } 
+
+    // allow multiple threads to read different segments into different readouts
+    // this thread may not have read a segment yet; but others have.  trust readout->imageScan
+
+    // is the start scan of the read less than the last possible scan?
+    int thisScan = readoutGetThisScan(readout, type);
+    int lastScan = (readdir == 1) ? trimsec->y1 : trimsec->x1;
+
+    *result = (lastScan > thisScan);
     return true;
+    // XXX fold this and the above case together
 }
 
@@ -217 +248 @@
     psImage *image = *readoutImageByType(readout, type);
 
+    // XXX this may not be the valid test in a multithread environment. consider a fileGroup of
+    // N readouts, but numScans set to 0.  only the first should report that it requires data, 
+    // even if all readouts lack the image pointer.
     if (!image) {
         return true;
-    } else if (numScans == 0) {
-        // Can only read the entire image once
+    } 
+
+    // If we have already read an image, this result implies we are done (no more to read)
+    if (numScans == 0) {
         return false;
     }
@@ -235 +271 @@
     }
 
-    int next;                           // Next position
-    int last;                           // Last position
-    if (!readoutScanProperties(&next, &last, readout, numScans, type, config)) {
+    bool haveData;
+    if (!readoutHaveMoreScans(&haveData, readout, numScans, type, config)) {
         psError(PS_ERR_UNKNOWN, false, "Unable to determine readout properties.");
         return false;
     }
 
-    return (next < last);
+    return haveData;
 }
 
@@ -423 +458 @@
     }
 
-    int next;                           // Next position
-    int last;                           // Last position
-    if (!readoutScanProperties(&next, &last, readout, numScans, type, config)) {
+    bool haveData;
+    if (!readoutHaveMoreScans (&haveData, readout, type, numScans, config)) {
         psError(PS_ERR_UNKNOWN, false, "Unable to determine readout properties.");
         return false;
     }
-    if (next >= last) {
+    if (!haveData) {
         psError(PS_ERR_IO, true, "No more of the readout to read.");
         return false;
     }
+
 
     pmHDU *hdu = pmHDUFromCell(cell);   // The HDU
@@ -459 +494 @@
     }
 
-
-
     // Check the third dimension
     int naxis = psMetadataLookupS32(&mdok, hdu->header, "NAXIS"); // The number of axes
@@ -478 +511 @@
     }
 
-    // Update number of scans read
-    int *scansTracker = readoutScansByType(readout, type); // Tracker for how many scans have been read
-    *scansTracker = next;
-    if (next == 0) {
+    // Determine the numbe of scans to read
+    int thisScan = readoutGetThisScan(readout, type);
+    int lastScan = readoutSetLastScan (readout, type, numScans);
+    if (thisScan == 0) {
         overlap = 0;
     }
-    next -= overlap;
+    thisScan -= overlap;
 
     // Calculate limits, adjust readout->row0,col0
@@ -491 +524 @@
     if (readdir == 1) {
         // Reading rows
-        readout->row0 = next;
+        readout->row0 = thisScan;
         readout->col0 = trimsec->x0;
         if (numScans == 0) {
@@ -498 +531 @@
     } else {
         // Reading cols
-        readout->col0 = next;
+        readout->col0 = thisScan;
         readout->row0 = trimsec->y0;
         if (numScans == 0) {
@@ -504 +537 @@
         }
     }
-    int upper = next + numScans + overlap;        // Upper limit to next section
 
     // Blow away existing data.
@@ -510 +542 @@
     psFree(*image);
     *image = NULL;
-    *image = readoutReadComponent(*image, fits, trimsec, readdir, next, upper, z, bad, pixelTypes[type]);
+    *image = readoutReadComponent(*image, fits, trimsec, readdir, thisScan, lastScan, z, bad, pixelTypes[type]);
 
     // Read overscans only for "image" type --- weights and masks shouldn't record overscans
@@ -528 +560 @@
         psRegion *biassec = NULL;           // Bias section from iteration
         while ((biassec = psListGetAndIncrement(biassecsIter))) {
-            psImage *bias = readoutReadComponent(NULL, fits, biassec, readdir, next, upper, z,
-                                                 bad, pixelTypes[type]); // The bias
+            psImage *bias = readoutReadComponent(NULL, fits, biassec, readdir, thisScan, lastScan, z, bad, pixelTypes[type]); // The bias
             psListAdd(readout->bias, PS_LIST_TAIL, bias);
             psFree(bias);                   // Drop reference

trunk/psModules/src/camera/pmReadoutStack.c

@@ -8 +8 @@
 
 #include "pmReadoutStack.h"
+
+// generate the specified image
+// XXX should it be an error for the image to exist?
+psImage *pmReadoutSetAnalysisImage(pmReadout *readout, // Readout containing image
+                                   const char *name, // Name of image in analysis metadata
+                                   int numCols, int numRows, // Expected size of image
+                                   psElemType type, // Expected type of image
+                                   double init // Initial value
+    )
+{
+    PS_ASSERT_PTR_NON_NULL(readout, false);
+    PS_ASSERT_STRING_NON_EMPTY(name, false);
+
+    psImage *image = psImageAlloc(numCols, numRows, type);
+
+    if (!psMetadataAddImage(readout->analysis, PS_LIST_TAIL, name, 0, "Analysis image from " __FILE__, image)) {
+        psAbort ("analysis image already exists");
+    }
+    psImageInit(image, init);
+
+    psFree (image); // we still have a view on readout->analysis
+    return image;
+}
+
+// retrieve the specified image
+// XXX not sure why this should call psMemIncrRefCounter
+psImage *pmReadoutGetAnalysisImage(pmReadout *readout, // Readout containing image
+                                   const char *name       // Name of image in analysis metadata
+    )
+{
+    PS_ASSERT_PTR_NON_NULL(readout, false);
+    PS_ASSERT_STRING_NON_EMPTY(name, false);
+
+    bool mdok;                          // Status of MD lookup
+    psImage *image = psMetadataLookupPtr(&mdok, readout->analysis, name);
+    return image;
+}
 
 psImage *pmReadoutAnalysisImage(pmReadout *readout, // Readout containing image
@@ -40 +77 @@
 }
 
-bool pmReadoutUpdateSize(pmReadout *readout, int minCols, int minRows,
-                         int numCols, int numRows, bool mask, bool weight,
-                         psMaskType blank)
-{
-    PS_ASSERT_PTR_NON_NULL(readout, false);
-
-    if (readout->image) {
-        readout->col0 = PS_MIN(minCols, readout->col0);
-        readout->row0 = PS_MIN(minRows, readout->row0);
-    } else {
-        readout->col0 = minCols;
-        readout->row0 = minRows;
-    }
-
-    // (reAllocate the images
-    if (!readout->image) {
-        readout->image = psImageAlloc(numCols, numRows, PS_TYPE_F32);
-        psImageInit(readout->image, NAN);
-    }
-    if (readout->image->numCols < numCols || readout->image->numRows < numRows) {
-        // Generate the new output image by extending the current one, or making a whole new one
-        psImage *newImage = psImageAlloc(numCols, numRows, PS_TYPE_F32);
-        psImageInit(newImage, NAN);
-        psImageOverlaySection(newImage, readout->image, readout->col0, readout->row0, "=");
-        psFree(readout->image);
-        readout->image = newImage;
-    }
+// XXX for the moment, use col0, row0, numCols, numRows supplied from the outside
+bool pmReadoutStackDefineOutput(pmReadout *readout, int col0, int row0, int numCols, int numRows, bool mask, bool weight, psMaskType blank)
+{
+    PS_ASSERT_PTR_NON_NULL(readout, false);
+
+    // XXX is this an error?
+    if (readout->image) return false;
+    readout->col0 = col0;
+    readout->row0 = row0;
+
+    // allocate the images
+    readout->image = psImageAlloc(numCols, numRows, PS_TYPE_F32);
+    psImageInit(readout->image, NAN);
 
     if (mask) {
-        if (!readout->mask) {
-            readout->mask = psImageAlloc(numCols, numRows, PS_TYPE_MASK);
-            psImageInit(readout->mask, blank);
-        }
-        if (readout->mask->numCols < numCols || readout->mask->numRows < numRows) {
-            psImage *newMask = psImageAlloc(numCols, numRows, PS_TYPE_MASK);
-            psImageInit(newMask, blank);
-            psImageOverlaySection(newMask, readout->mask, readout->col0, readout->row0, "=");
-            psFree(readout->mask);
-            readout->mask = newMask;
-        }
+        // XXX is this an error?
+        if (readout->mask) return false;
+        readout->mask = psImageAlloc(numCols, numRows, PS_TYPE_MASK);
+        psImageInit(readout->mask, blank);
     }
 
     if (weight) {
-        if (!readout->weight) {
-            readout->weight = psImageAlloc(numCols, numRows, PS_TYPE_F32);
-            psImageInit(readout->weight, NAN);
-        }
-        if (readout->weight->numCols < numCols || readout->weight->numRows < numRows) {
-            psImage *newWeight = psImageAlloc(numCols, numRows, PS_TYPE_F32);
-            psImageInit(newWeight, NAN);
-            psImageOverlaySection(newWeight, readout->weight, readout->col0, readout->row0, "=");
-            psFree(readout->weight);
-            readout->weight = newWeight;
-        }
+        // XXX is this an error?
+        if (readout->weight) return false;
+        readout->weight = psImageAlloc(numCols, numRows, PS_TYPE_F32);
+        psImageInit(readout->weight, NAN);
     }
 
@@ -99 +108 @@
 }
 
-bool pmReadoutStackValidate(int *minInputColsPtr, int *maxInputColsPtr, int *minInputRowsPtr,
-                            int *maxInputRowsPtr, int *numColsPtr, int *numRowsPtr,
-                            const psArray *inputs)
+bool pmReadoutStackSetOutputSize(int *col0, int *row0, int *numCols, int *numRows, const psArray *inputs)
 {
     PS_ASSERT_ARRAY_NON_NULL(inputs, false);
-    PS_ASSERT_PTR_NON_NULL(minInputColsPtr, false);
-    PS_ASSERT_PTR_NON_NULL(maxInputColsPtr, false);
-    PS_ASSERT_PTR_NON_NULL(minInputRowsPtr, false);
-    PS_ASSERT_PTR_NON_NULL(maxInputRowsPtr, false);
-    PS_ASSERT_PTR_NON_NULL(numColsPtr, false);
-    PS_ASSERT_PTR_NON_NULL(numRowsPtr, false);
-
-    // Step through each readout in the input image list to determine how big of an output image is needed to
-    // combine these input images.
-    int maxInputCols = 0;               // The largest input column value
-    int maxInputRows = 0;               // The largest input row value
-    int minInputCols = INT_MAX;         // The smallest input column value
-    int minInputRows = INT_MAX;         // The smallest input row value
-    int xSize = 0, ySize = 0;           // The size of the output image
+    PS_ASSERT_PTR_NON_NULL(col0, false);
+    PS_ASSERT_PTR_NON_NULL(row0, false);
+    PS_ASSERT_PTR_NON_NULL(numCols, false);
+    PS_ASSERT_PTR_NON_NULL(numRows, false);
+
+    // Step through each readout in the input image list to determine how big of an output
+    // image is needed to combine these input images.
 
     int xMin = INT_MAX;
@@ -123 +123 @@
     int xMax = 0;
     int yMax = 0;
+    int xSize = 0;
+    int ySize = 0;           // The size of the output image
 
     bool valid = false;                 // Do we have a single valid input?
@@ -128 +130 @@
         pmReadout *readout = inputs->data[i]; // Readout of interest
 
-        if (!readout) {
-            continue;
-        }
-        if (!readout->image) {
-            psError(PS_ERR_UNEXPECTED_NULL, true, "Input readout %ld has NULL image.\n", i);
-            return false;
-        }
+        if (!readout) continue;
 
         // use the trimsec to define the max full range of the output pixels
@@ -150 +146 @@
             yMax  = PS_MAX(yMax,  trimsec->y1);
         }
+        valid = true;
+        psTrace("psModules.camera", 7, "Readout %ld: trimsec: %f,%f - %f,%f\n", i, trimsec->x0, trimsec->y0, trimsec->x1, trimsec->y1);
+    }
+
+    *col0 = xMin;
+    *row0 = yMin;
+    *numCols = xSize;
+    *numRows = ySize;
+
+    if (!valid) {
+        psError(PS_ERR_UNKNOWN, false, "No valid input readouts.");
+    }
+    return valid;
+}
+
+bool pmReadoutUpdateSize(pmReadout *readout, int minCols, int minRows,
+                         int numCols, int numRows, bool mask, bool weight,
+                         psMaskType blank)
+{
+    PS_ASSERT_PTR_NON_NULL(readout, false);
+
+    if (readout->image) {
+        readout->col0 = PS_MIN(minCols, readout->col0);
+        readout->row0 = PS_MIN(minRows, readout->row0);
+    } else {
+        readout->col0 = minCols;
+        readout->row0 = minRows;
+    }
+
+    // (reAllocate the images
+    if (!readout->image) {
+        readout->image = psImageAlloc(numCols, numRows, PS_TYPE_F32);
+        psImageInit(readout->image, NAN);
+    }
+    if (readout->image->numCols < numCols || readout->image->numRows < numRows) {
+        // Generate the new output image by extending the current one, or making a whole new one
+        psImage *newImage = psImageAlloc(numCols, numRows, PS_TYPE_F32);
+        psImageInit(newImage, NAN);
+        psImageOverlaySection(newImage, readout->image, readout->col0, readout->row0, "=");
+        psFree(readout->image);
+        readout->image = newImage;
+    }
+
+    if (mask) {
+        if (!readout->mask) {
+            readout->mask = psImageAlloc(numCols, numRows, PS_TYPE_MASK);
+            psImageInit(readout->mask, blank);
+        }
+        if (readout->mask->numCols < numCols || readout->mask->numRows < numRows) {
+            psImage *newMask = psImageAlloc(numCols, numRows, PS_TYPE_MASK);
+            psImageInit(newMask, blank);
+            psImageOverlaySection(newMask, readout->mask, readout->col0, readout->row0, "=");
+            psFree(readout->mask);
+            readout->mask = newMask;
+        }
+    }
+
+    if (weight) {
+        if (!readout->weight) {
+            readout->weight = psImageAlloc(numCols, numRows, PS_TYPE_F32);
+            psImageInit(readout->weight, NAN);
+        }
+        if (readout->weight->numCols < numCols || readout->weight->numRows < numRows) {
+            psImage *newWeight = psImageAlloc(numCols, numRows, PS_TYPE_F32);
+            psImageInit(newWeight, NAN);
+            psImageOverlaySection(newWeight, readout->weight, readout->col0, readout->row0, "=");
+            psFree(readout->weight);
+            readout->weight = newWeight;
+        }
+    }
+
+    return true;
+}
+
+bool pmReadoutStackValidate(int *minInputColsPtr, int *maxInputColsPtr, int *minInputRowsPtr,
+                            int *maxInputRowsPtr, int *numColsPtr, int *numRowsPtr,
+                            const psArray *inputs)
+{
+    PS_ASSERT_ARRAY_NON_NULL(inputs, false);
+    PS_ASSERT_PTR_NON_NULL(minInputColsPtr, false);
+    PS_ASSERT_PTR_NON_NULL(maxInputColsPtr, false);
+    PS_ASSERT_PTR_NON_NULL(minInputRowsPtr, false);
+    PS_ASSERT_PTR_NON_NULL(maxInputRowsPtr, false);
+    PS_ASSERT_PTR_NON_NULL(numColsPtr, false);
+    PS_ASSERT_PTR_NON_NULL(numRowsPtr, false);
+
+    // Step through each readout in the input image list to determine how big of an output image is needed to
+    // combine these input images.
+    int maxInputCols = 0;               // The largest input column value
+    int maxInputRows = 0;               // The largest input row value
+    int minInputCols = INT_MAX;         // The smallest input column value
+    int minInputRows = INT_MAX;         // The smallest input row value
+    int xSize = 0, ySize = 0;           // The size of the output image
+
+    int xMin = INT_MAX;
+    int yMin = INT_MAX;
+    int xMax = 0;
+    int yMax = 0;
+
+    bool valid = false;                 // Do we have a single valid input?
+    for (long i = 0; i < inputs->n; i++) {
+        pmReadout *readout = inputs->data[i]; // Readout of interest
+
+        if (!readout) {
+            continue;
+        }
+        if (!readout->image) {
+            psError(PS_ERR_UNEXPECTED_NULL, true, "Input readout %ld has NULL image.\n", i);
+            return false;
+        }
+
+        // use the trimsec to define the max full range of the output pixels
+        pmCell *cell = readout->parent; // The parent cell
+        bool mdok = true;       // Status of MD lookup
+        psRegion *trimsec = psMetadataLookupPtr(&mdok, cell->concepts, "CELL.TRIMSEC"); // Trim section
+        if (!mdok || !trimsec || psRegionIsNaN(*trimsec)) {
+            psWarning("CELL.TRIMSEC is not set for readout %ld --- ignored.\n", i);
+        } else {
+            xSize = PS_MAX(xSize, trimsec->x1 - trimsec->x0);
+            ySize = PS_MAX(ySize, trimsec->y1 - trimsec->y0);
+            xMin  = PS_MIN(xMin,  trimsec->x0);
+            xMax  = PS_MAX(xMax,  trimsec->x1);
+            yMin  = PS_MIN(yMin,  trimsec->y0);
+            yMax  = PS_MAX(yMax,  trimsec->y1);
+        }
 
         valid = true;
@@ -158 +279 @@
         minInputRows = PS_MAX(yMin, PS_MIN(minInputRows, readout->row0));
         maxInputRows = PS_MIN(yMax, PS_MAX(maxInputRows, readout->row0 + readout->image->numRows));
-        psTrace("psModules.camera", 7, "Readout %ld: offset %d,%d; size %dx%d\n", i,
-                readout->col0, readout->row0, readout->image->numCols, readout->image->numRows);
+
+        psTrace("psModules.camera", 7, "Readout %ld: offset %d,%d; size %dx%d\n", i, readout->col0, readout->row0, readout->image->numCols, readout->image->numRows);
     }
 

trunk/psModules/src/camera/pmReadoutStack.h

@@ -24 +24 @@
     );
 
+psImage *pmReadoutSetAnalysisImage(pmReadout *readout, // Readout containing image
+                                   const char *name, // Name of image in analysis metadata
+                                   int numCols, int numRows, // Expected size of image
+                                   psElemType type, // Expected type of image
+                                   double init // Initial value
+    );
+
+// retrieve the specified image
+// XXX not sure why this should call psMemIncrRefCounter
+psImage *pmReadoutGetAnalysisImage(pmReadout *readout, // Readout containing image
+                                   const char *name       // Name of image in analysis metadata
+    );
+
+
 /// Return an image from analysis metadata, produced while stacking
 psImage *pmReadoutAnalysisImage(pmReadout *readout, // Readout containing image
@@ -32 +46 @@
     );
 
+// XXX for the moment, use col0, row0, numCols, numRows supplied from the outside
+bool pmReadoutStackDefineOutput(pmReadout *readout, int col0, int row0, int numCols, int numRows, bool mask, bool weight, psMaskType blank);
+
+bool pmReadoutStackSetOutputSize(int *col0, int *row0, int *numCols, int *numRows, const psArray *inputs);
+
 #endif

trunk/psModules/src/config/pmConfig.c

@@ -1627 +1627 @@
 
         if (trunc) {
-            if(truncate(filename, 0) != 0) {
-                psError(PS_ERR_IO, true, "Failed to truncate file, %s\n", filename);
-                nebServerFree(server);
+            if(truncate(path, 0) != 0) {
+                psError(PS_ERR_IO, true, "Failed to truncate Nebulous file %s (real name %s)\n", filename, path);
                 return NULL;
             }

trunk/psModules/src/imcombine/pmReadoutCombine.c

@@ -42 +42 @@
 }
 
+// check the input parameters and set up the output images
+bool pmReadoutCombinePrepare(pmReadout *output, const psArray *inputs, const pmCombineParams *params)
+{
+    // Check inputs
+    PS_ASSERT_PTR_NON_NULL(output, false);
+    PS_ASSERT_ARRAY_NON_NULL(inputs, false);
+    PS_ASSERT_PTR_NON_NULL(params, false);
+    PS_ASSERT_FLOAT_WITHIN_RANGE(params->fracLow, 0.0, 1.0, false);
+    PS_ASSERT_FLOAT_WITHIN_RANGE(params->fracHigh, 0.0, 1.0, false);
+
+    // valid combintion statistic?
+    bool valid = false;
+    valid |= (params->combine == PS_STAT_SAMPLE_MEAN);
+    valid |= (params->combine == PS_STAT_SAMPLE_MEDIAN);
+    valid |= (params->combine == PS_STAT_ROBUST_MEDIAN);
+    valid |= (params->combine == PS_STAT_FITTED_MEAN);
+    valid |= (params->combine == PS_STAT_CLIPPED_MEAN);
+    if (!valid) {
+        psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Combination method is not SAMPLE_MEAN, SAMPLE_MEDIAN, "
+                "ROBUST_MEDIAN, FITTED_MEAN or CLIPPED_MEAN.\n");
+        return false;
+    }
+
+    // weights exist if weights desired?
+    for (int i = 0; i < inputs->n; i++) {
+        pmReadout *readout = inputs->data[i]; // Readout of interest
+        if (params->weights && !readout->weight) {
+            psError(PS_ERR_UNEXPECTED_NULL, true,
+                    "Rejection based on weights requested, but no weights supplied for image %d.\n", i);
+            return false;
+        }
+    }
+
+    pmHDU *hdu = pmHDUFromReadout(output); // Output HDU
+    if (!hdu) {
+        psError(PS_ERR_UNEXPECTED_NULL, false, "Unable to find HDU for readout.\n");
+        return false;
+    }
+
+    //  set the output header metadata
+    psString comment = NULL;        // Comment to add to header
+    psStringAppend(&comment, "Combining using statistic: %x", params->combine);
+    if (!hdu->header) {
+        hdu->header = psMetadataAlloc();
+    }
+    psMetadataAddStr(hdu->header, PS_LIST_TAIL, "HISTORY", PS_META_DUPLICATE_OK, comment, "");
+    psFree(comment);
+
+    // note the clipping parameters, if used
+    if (params->combine == PS_STAT_CLIPPED_MEAN) {
+        psString comment = NULL;    // Comment to add to header
+        psStringAppend(&comment, "Combination clipping: %d iterations, rejection at %f sigma", params->iter, params->rej);
+        psMetadataAddStr(hdu->header, PS_LIST_TAIL, "HISTORY", PS_META_DUPLICATE_OK, comment, "");
+        psFree(comment);
+    }
+
+    // note the use of weights
+    if (params->weights) {
+        psMetadataAddStr(hdu->header, PS_LIST_TAIL, "HISTORY", PS_META_DUPLICATE_OK,
+                         "Using input weights to combine images", "");
+    }
+
+    // note the rejection fraction
+    float keepFrac = 1.0 - params->fracLow - params->fracHigh; // Fraction of pixels to keep
+    if (keepFrac != 1.0) {
+        psString comment = NULL;        // Comment to add to header
+        psStringAppend(&comment, "Min/max rejection: %f high, %f low, keep %d",
+                       params->fracHigh, params->fracLow, params->nKeep);
+        psMetadataAddStr(hdu->header, PS_LIST_TAIL, "HISTORY", PS_META_DUPLICATE_OK, comment, "");
+        psFree(comment);
+    }
+
+    // note the mask value actually used
+    psMaskType maskVal = params->maskVal; // The mask value
+    if (maskVal) {
+        psString comment = NULL;        // Comment to add to header
+        psStringAppend(&comment, "Mask for combination: %x", maskVal);
+        psMetadataAddStr(hdu->header, PS_LIST_TAIL, "HISTORY", PS_META_DUPLICATE_OK, comment, "");
+        psFree(comment);
+    }
+
+    // determine the output image size based on the input images
+    int row0, col0, numCols, numRows;
+    if (!pmReadoutStackSetOutputSize(&col0, &row0, &numCols, &numRows, inputs)) {
+        psError(PS_ERR_UNKNOWN, false, "problem setting output readout size.");
+        return false;
+    }
+
+    // generate the required output images based on the specified sizes
+    pmReadoutStackDefineOutput(output, col0, row0, numCols, numRows, true, params->weights, params->blank);
+    psTrace("psModules.imcombine", 7, "Output minimum: %d,%d\n", output->col0, output->row0);
+
+    // these calls allocate and save the requested images on the output analysis metadata
+    psImage *counts = pmReadoutSetAnalysisImage(output, PM_READOUT_STACK_ANALYSIS_COUNT, numCols, numRows, PS_TYPE_U16, 0);
+    if (!counts) {
+        return false;
+    }
+    psImage *sigma = pmReadoutSetAnalysisImage(output, PM_READOUT_STACK_ANALYSIS_SIGMA, numCols, numRows, PS_TYPE_F32, NAN);
+    if (!sigma) {
+        return false;
+    }
+
+    // Update the "concepts"
+    psList *inputCells = psListAlloc(NULL); // List of cells
+    for (long i = 0; i < inputs->n; i++) {
+        pmReadout *readout = inputs->data[i]; // Readout of interest
+        psListAdd(inputCells, PS_LIST_TAIL, readout->parent);
+    }
+    bool success = pmConceptsAverageCells(output->parent, inputCells, NULL, NULL, true);
+    psFree(inputCells);
+
+    // set these even though the values are not yet set
+    output->data_exists = true;
+    output->parent->data_exists = true;
+    output->parent->parent->data_exists = true;
+
+    return success;
+}
 
 // XXX: Maybe add support for S16 and S32 types.  Currently, only F32 supported.
@@ -61 +179 @@
     PS_ASSERT_FLOAT_WITHIN_RANGE(params->fracLow, 0.0, 1.0, false);
     PS_ASSERT_FLOAT_WITHIN_RANGE(params->fracHigh, 0.0, 1.0, false);
-    if (params->combine != PS_STAT_SAMPLE_MEAN && params->combine != PS_STAT_SAMPLE_MEDIAN &&
-            params->combine != PS_STAT_ROBUST_MEDIAN && params->combine != PS_STAT_FITTED_MEAN &&
-            params->combine != PS_STAT_CLIPPED_MEAN) {
-        psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Combination method is not SAMPLE_MEAN, SAMPLE_MEDIAN, "
-                "ROBUST_MEDIAN, FITTED_MEAN or CLIPPED_MEAN.\n");
-        return false;
-    }
-    for (int i = 0; i < inputs->n; i++) {
-        pmReadout *readout = inputs->data[i]; // Readout of interest
-        if (params->weights && !readout->weight) {
-            psError(PS_ERR_UNEXPECTED_NULL, true,
-                    "Rejection based on weights requested, but no weights supplied for image %d.\n", i);
-            return false;
-        }
-    }
-
-    bool first = !output->image;        // First pass through?
 
     pmHDU *hdu = pmHDUFromReadout(output); // Output HDU
@@ -85 +186 @@
     }
 
-    if (first) {
-        psString comment = NULL;        // Comment to add to header
-        psStringAppend(&comment, "Combining using statistic: %x", params->combine);
-        if (!hdu->header) {
-            hdu->header = psMetadataAlloc();
-        }
-        psMetadataAddStr(hdu->header, PS_LIST_TAIL, "HISTORY", PS_META_DUPLICATE_OK, comment, "");
-        psFree(comment);
-    }
+    pthread_t id = pthread_self();
+    char name[64];
+    sprintf (name, "%x", (unsigned int) id);
+    psTimerStart (name);
 
     psStatsOptions combineStdev = 0; // Statistics option for weights
@@ -118 +214 @@
         stats->clipSigma = params->rej;
         stats->clipIter = params->iter;
-
-        if (first) {
-            psString comment = NULL;    // Comment to add to header
-            psStringAppend(&comment, "Combination clipping: %d iterations, rejection at %f sigma",
-                           params->iter, params->rej);
-            psMetadataAddStr(hdu->header, PS_LIST_TAIL, "HISTORY", PS_META_DUPLICATE_OK, comment, "");
-            psFree(comment);
-        }
-    }
-    if (params->weights && first) {
-        psMetadataAddStr(hdu->header, PS_LIST_TAIL, "HISTORY", PS_META_DUPLICATE_OK,
-                         "Using input weights to combine images", "");
-    }
+    }
+
+    psImage *counts = pmReadoutGetAnalysisImage(output, PM_READOUT_STACK_ANALYSIS_COUNT);
+    if (!counts) {
+        return false;
+    }
+    psImage *sigma = pmReadoutGetAnalysisImage(output, PM_READOUT_STACK_ANALYSIS_SIGMA);
+    if (!sigma) {
+        return false;
+    }
+
+    stats->options |= combineStdev;
 
     int minInputCols, maxInputCols, minInputRows, maxInputRows; // Smallest and largest values to combine
     int xSize, ySize;                   // Size of the output image
-    if (!pmReadoutStackValidate(&minInputCols, &maxInputCols, &minInputRows, &maxInputRows, &xSize, &ySize,
-                                inputs)) {
+    if (!pmReadoutStackValidate(&minInputCols, &maxInputCols, &minInputRows, &maxInputRows, &xSize, &ySize, inputs)) {
         psError(PS_ERR_UNKNOWN, false, "No valid input readouts.");
         return false;
     }
-
-    pmReadoutUpdateSize(output, minInputCols, minInputRows, xSize, ySize, true, params->weights,
-                        params->blank);
-    psTrace("psModules.imcombine", 7, "Output minimum: %d,%d\n", output->col0, output->row0);
-
-    psImage *counts = pmReadoutAnalysisImage(output, PM_READOUT_STACK_ANALYSIS_COUNT, xSize, ySize,
-                                             PS_TYPE_U16, 0);
-    if (!counts) {
-        return false;
-    }
-    psImage *sigma = pmReadoutAnalysisImage(output, PM_READOUT_STACK_ANALYSIS_SIGMA, xSize, ySize,
-                                            PS_TYPE_F32, NAN);
-    if (!sigma) {
-        psFree(counts);
-        return false;
-    }
-
-    stats->options |= combineStdev;
 
     // We loop through each pixel in the output image.  We loop through each input readout.  We determine if
@@ -180 +256 @@
 
     float keepFrac = 1.0 - params->fracLow - params->fracHigh; // Fraction of pixels to keep
-    if (keepFrac != 1.0 && first) {
-        psString comment = NULL;        // Comment to add to header
-        psStringAppend(&comment, "Min/max rejection: %f high, %f low, keep %d",
-                       params->fracHigh, params->fracLow, params->nKeep);
-        psMetadataAddStr(hdu->header, PS_LIST_TAIL, "HISTORY", PS_META_DUPLICATE_OK, comment, "");
-        psFree(comment);
-    }
-
     psMaskType maskVal = params->maskVal; // The mask value
-    if (maskVal && first) {
-        psString comment = NULL;        // Comment to add to header
-        psStringAppend(&comment, "Mask for combination: %x", maskVal);
-        psMetadataAddStr(hdu->header, PS_LIST_TAIL, "HISTORY", PS_META_DUPLICATE_OK, comment, "");
-        psFree(comment);
-    }
 
     #ifndef PS_NO_TRACE
@@ -226 +288 @@
         int yOut = i - output->row0; // y position on output readout
         #ifdef SHOW_BUSY
-
+        
         if (psTraceGetLevel("psModules.imcombine") > 9) {
             printf("Processing row %d\r", i);
@@ -242 +304 @@
                 int xIn = j - readout->col0; // x position on input readout
                 psImage *image = readout->image; // The readout image
-
-                #if 0 // This should have been taken care of already:
-                // Check bounds
-                if (xIn < 0 || xIn >= image->numCols || yIn < 0 || yIn >= image->numRows) {
-                    continue;
-                }
-                #endif
 
                 pixelsData[r] = image->data.F32[yIn][xIn];
@@ -348 +403 @@
     }
     #endif
+
     psFree(index);
     psFree(pixels);
@@ -355 +411 @@
     psFree(stats);
     psFree(invScale);
-    psFree(counts);
-    psFree(sigma);
-
-    // Update the "concepts"
-    psList *inputCells = psListAlloc(NULL); // List of cells
-    for (long i = 0; i < inputs->n; i++) {
-        pmReadout *readout = inputs->data[i]; // Readout of interest
-        psListAdd(inputCells, PS_LIST_TAIL, readout->parent);
-    }
-    bool success = pmConceptsAverageCells(output->parent, inputCells, NULL, NULL, true);
-    psFree(inputCells);
-
-    output->data_exists = true;
-    output->parent->data_exists = true;
-    output->parent->parent->data_exists = true;
-
-    return success;
+
+    // fprintf (stderr, "done with combine %x : %f sec\n", (unsigned int) id, psTimerMark (name));
+    return true;
 }
 

trunk/psModules/src/imcombine/pmReadoutCombine.h

@@ -5 +5 @@
  * @author Paul Price, IfA
  *
- * @version $Revision: 1.13 $ $Name: not supported by cvs2svn $
- * @date $Date: 2008-03-29 03:10:17 $
+ * @version $Revision: 1.14 $ $Name: not supported by cvs2svn $
+ * @date $Date: 2008-08-01 00:01:26 $
  * Copyright 2004-2006 Institute for Astronomy, University of Hawaii
  */
@@ -36 +36 @@
                                      );
 
+// check the input parameters and set up the output images
+bool pmReadoutCombinePrepare(pmReadout *output, const psArray *inputs, const pmCombineParams *params);
+
 /// Combine multiple readouts, applying zero and scale, with optional minmax clipping
 bool pmReadoutCombine(pmReadout *output,///< Output readout; altered and returned

trunk/psModules/src/objects/Makefile.am

@@ -4 +4 @@
 libpsmodulesobjects_la_LDFLAGS  = -release $(PACKAGE_VERSION)
 libpsmodulesobjects_la_SOURCES  = \
+     pmDetections.c \
+     pmSpan.c \
+     pmFootprint.c \
+     pmFootprintArrayGrow.c \
+     pmFootprintArraysMerge.c \
+     pmFootprintAssignPeaks.c \
+     pmFootprintFind.c \
+     pmFootprintFindAtPoint.c \
+     pmFootprintIDs.c \
      pmPeaks.c \
      pmMoments.c \
@@ -45 +54 @@
 
 pkginclude_HEADERS = \
+     pmDetections.h \
+     pmSpan.h \
+     pmFootprint.h \
      pmPeaks.h \
      pmMoments.h \