Changeset 24244 for branches/cnb_branches/cnb_branch_20090301/ppMerge/src

Timestamp:

May 26, 2009, 1:59:32 PM (17 years ago)

Author:

beaumont

Message:

merged with trunk

Location:

branches/cnb_branches/cnb_branch_20090301

Files:

• . (modified) (1 prop)
• ppMerge/src (modified) (1 prop)
• ppMerge/src/Makefile.am (modified) (1 diff)
• ppMerge/src/ppMerge.c (modified) (1 diff)
• ppMerge/src/ppMergeArguments.c (modified) (4 diffs)
• ppMerge/src/ppMergeErrorCodes.dat (modified) (1 diff)
• ppMerge/src/ppMergeErrorCodes.h.in (modified) (1 diff)
• ppMerge/src/ppMergeLoop.c (modified) (6 diffs)
• ppMerge/src/ppMergeMask.c (modified) (5 diffs)
• ppMerge/src/ppMergeVersion.c (modified) (3 diffs)
• ppMerge/src/ppMergeVersionDefinitions.h.in (copied) (copied from trunk/ppMerge/src/ppMergeVersionDefinitions.h.in )

branches/cnb_branches/cnb_branch_20090301/ppMerge/src

@@ -7 +7 @@
 stamp-h1
 .libs
+ppMergeVersionDefinitions.h

@@ -7 +7 @@
 stamp-h1
 .libs
+ppMergeVersionDefinitions.h

branches/cnb_branches/cnb_branch_20090301/ppMerge/src/Makefile.am

@@ -1 +1 @@
 bin_PROGRAMS = ppMerge
 
-# PPMERGE_VERSION=`if [ -e ../../VERSION ]; then cat ../../VERSION; else svnversion; fi`
-# PPMERGE_BRANCH=`if [ -e ../../BRANCH ]; then cat ../../BRANCH; else svn info | sed -n -e '/URL:/ h' -e '/Repository Root:/ { x; H; x; s|Repository Root: \(.*\)\nURL: \1\(.*\)|\2| ; s|^/|| ; s|/[a-zA-Z]*/src.*|| ; p }'; fi`
-# PPMERGE_SOURCE=`if [ -e ../../SOURCE ]; then cat ../../SOURCE; else svn info | sed -n -e 's/Repository UUID: // p'; fi`
-# 
-# # Force recompilation of ppMergeVersion.c, since it gets the version information
-# ppMergeVersion.c: FORCE
-#       touch ppMergeVersion.c
-# FORCE: ;
+if HAVE_SVNVERSION
+PPMERGE_VERSION=`$(SVNVERSION) ..`
+else
+PPMERGE_VERSION="UNKNOWN"
+endif
 
-ppMerge_CFLAGS = $(PPMERGE_CFLAGS) $(PPSTATS_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS) -DPPMERGE_VERSION=$(SVN_VERSION) -DPPMERGE_BRANCH=$(SVN_BRANCH) -DPPMERGE_SOURCE=$(SVN_SOURCE)
+if HAVE_SVN
+PPMERGE_BRANCH=`$(SVN) info .. | $(SED) -n -e '/URL:/ h' -e '/Repository Root:/ { x; H; x; s|Repository Root: \(.*\)\nURL: \1\(.*\)|\2| ; s|^/|| ; s|/[a-zA-Z]*/src.*|| ; p }'`
+PPMERGE_SOURCE=`$(SVN) info | $(SED) -n -e 's/Repository UUID: // p'`
+else
+PPMERGE_BRANCH="UNKNOWN"
+PPMERGE_SOURCE="UNKNOWN"
+endif
+
+# Force recompilation of ppMergeVersion.c, since it gets the version information
+ppMergeVersion.c: ppMergeVersionDefinitions.h
+ppMergeVersionDefinitions.h: ppMergeVersionDefinitions.h.in FORCE
+        -$(RM) ppMergeVersionDefinitions.h
+        $(SED) -e "s|@PPMERGE_VERSION@|\"$(PPMERGE_VERSION)\"|" -e "s|@PPMERGE_BRANCH@|\"$(PPMERGE_BRANCH)\"|" -e "s|@PPMERGE_SOURCE@|\"$(PPMERGE_SOURCE)\"|" ppMergeVersionDefinitions.h.in > ppMergeVersionDefinitions.h
+FORCE: ;
+
+
+ppMerge_CFLAGS = $(PPMERGE_CFLAGS) $(PPSTATS_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
 ppMerge_LDFLAGS = $(PPMERGE_LIBS) $(PPSTATS_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
 

branches/cnb_branches/cnb_branch_20090301/ppMerge/src/ppMerge.c

@@ -95 +95 @@
         psFree(statsOut);
         fclose(statsFile);
+        pmConfigRunFilenameAddWrite(config, "STATS", statsName);
     }
 

branches/cnb_branches/cnb_branch_20090301/ppMerge/src/ppMergeArguments.c

@@ -156 +156 @@
     psMetadataAddS32(arguments, PS_LIST_TAIL, "-nkeep",    0, "Minimum number of pixels in stack to keep", 0);
     psMetadataAddBool(arguments, PS_LIST_TAIL, "-variances", 0, "Use image variances in combination?", false);
+    psMetadataAddBool(arguments, PS_LIST_TAIL, "-use-masks", 0, "Use image masks in combination?", false);
 
     // XXX EAM : not clear this should be allowed on the command line.
@@ -179 +180 @@
 
     /** Mask construction parameters */
-    psMetadataAddF32(arguments, PS_LIST_TAIL, "-mask-suspect",  0, "Threshold for suspect pixels (sigma)", NAN);
+    psMetadataAddF32(arguments, PS_LIST_TAIL, "-mask-suspect-sigma",  0, "Threshold for suspect pixels (sigma)", NAN);
     psMetadataAddF32(arguments, PS_LIST_TAIL, "-mask-bad",      0, "Threshold for bad pixels (sigma)", NAN);
     psMetadataAddStr(arguments, PS_LIST_TAIL, "-mask-mode",     0, "Mode to identify bad pixels", NULL);
@@ -308 +309 @@
     VALUE_ARG_RECIPE_INT("-nkeep",      "NKEEP",    S32, 0);
     VALUE_ARG_RECIPE_BOOL("-variances", "VARIANCES");
+    VALUE_ARG_RECIPE_BOOL("-use-masks", "INPUTS.MASKS.USE");
 
     // XXX we do not supply this on the command line
@@ -329 +331 @@
 
     /** Mask construction parameters */
-    VALUE_ARG_RECIPE_FLOAT("-mask-suspect", "MASK.SUSPECT", F32);
+    VALUE_ARG_RECIPE_FLOAT("-mask-suspect-sigma", "MASK.SUSPECT.SIGMA", F32);
+    VALUE_ARG_RECIPE_FLOAT("-mask-suspect-min",   "MASK.SUSPECT.MIN", F32);
+    VALUE_ARG_RECIPE_FLOAT("-mask-suspect-max",   "MASK.SUSPECT.MAX", F32);
+    VALUE_ARG_RECIPE_STR("-mask-mode", "MASK.SUSPECT.MODE");
     VALUE_ARG_RECIPE_FLOAT("-mask-bad",     "MASK.BAD",     F32);
     VALUE_ARG_RECIPE_INT("-mask-grow",      "MASK.GROW",    S32, 0);

branches/cnb_branches/cnb_branch_20090301/ppMerge/src/ppMergeErrorCodes.dat

@@ -2 +2 @@
 # This file is used to generate ppMergeErrorClasses.h
 #
-BASE = 700              First value we use; lower values belong to psLib
+BASE = 6000             First value we use; lower values belong to psLib
 # these errors correspond to standard exit conditions
 ARGUMENTS               Incorrect arguments

branches/cnb_branches/cnb_branch_20090301/ppMerge/src/ppMergeErrorCodes.h.in

@@ -21 +21 @@
  */
 typedef enum {
-    PPMERGE_ERR_BASE = 512,
+    PPMERGE_ERR_BASE = 6000,
     PPMERGE_ERR_${ErrorCode},
     PPMERGE_ERR_NERROR

branches/cnb_branches/cnb_branch_20090301/ppMerge/src/ppMergeLoop.c

@@ -36 +36 @@
     bool mdok;                          ///< Status of MD lookup
     bool haveMasks = psMetadataLookupBool(&mdok, arguments, "INPUTS.MASKS"); // Do we have masks?
+    bool useMasks = psMetadataLookupBool(&mdok, arguments, "INPUTS.MASKS.USE"); // Do we have masks?
     bool haveVariances = psMetadataLookupBool(&mdok, arguments, "INPUTS.VARIANCES"); // Do we have variances?
 
@@ -73 +74 @@
 
     pmCombineParams *combination = pmCombineParamsAlloc(combineStat); ///< Combination parameters
-    combination->maskVal = maskVal;
+    combination->maskVal = useMasks ? maskVal : 0;
     combination->blank = pmConfigMaskGet("BLANK", config);
     combination->nKeep = nKeep;
@@ -149 +150 @@
         }
         pmCell *outCell;                ///< Cell of interest
+
+        // XXX TEST : force a single loop
+        // outCell = pmFPAviewNextCell(view, outFPA, 1); {
+
         while ((outCell = pmFPAviewNextCell(view, outFPA, 1))) {
             if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) {
@@ -164 +169 @@
                 pmHDU *hdu = pmHDUGetHighest(outFPA, outChip, outCell); // HDU for file
                 if (hdu && hdu != lastHDU) {
+                    if (!hdu->header) {
+                        hdu->header = psMetadataAlloc();
+                    }
                     ppMergeVersionHeader(hdu->header);
                     lastHDU = hdu;
@@ -229 +237 @@
                 psAssert (fileGroups->n > 0, "no valid file groups defined");
                 fileGroup = fileGroups->data[0];
-                if (!pmShutterCorrectionGeneratePrepare(outRO, pattern, fileGroup->readouts, maskVal)) {
+                if (!pmShutterCorrectionGeneratePrepare(outRO, pattern, fileGroup->readouts, combination->maskVal)) {
                     goto ERROR;
                 }
@@ -411 +419 @@
             }
 
-            // calculate CTEMASK after stats so stats reflect median image
+            // calculate CTEMASK after stats so stats reflect median image
             if (type == PPMERGE_TYPE_CTEMASK && outRO) {
-                // need to apply range cuts on the output image
-                psAssert (outRO->mask, "mask is not defined");
-                psAssert (outRO->mask->numCols == outRO->image->numCols, "mismatch between image and mask");
-                psAssert (outRO->mask->numRows == outRO->image->numRows, "mismatch between image and mask");
-
-                // CTEMASK parameters
-                float cteMin = psMetadataLookupF32(NULL, arguments, "CTE.MIN"); // Number of fringe points
-
-                char *cteMaskName = psMetadataLookupStr (&mdok, config->arguments, "MASK.SET.VALUE");
-                psImageMaskType cteMaskValue = pmConfigMaskGet(cteMaskName, config);
-
-                psF32 **outputImage = outRO->image->data.F32; 
-                psImageMaskType **outputMask = outRO->mask->data.PS_TYPE_IMAGE_MASK_DATA; 
-                for (int iy = 0; iy < outRO->image->numRows; iy++) {
-                    for (int ix = 0; ix < outRO->image->numCols; ix++) {
-                        if (outputImage[iy][ix] < cteMin) {
-                            outputMask[iy][ix] |= cteMaskValue;
-                        }
-                    }
+                // need to apply range cuts on the output image
+                psAssert (outRO->mask, "mask is not defined");
+                psAssert (outRO->mask->numCols == outRO->image->numCols, "mismatch between image and mask");
+                psAssert (outRO->mask->numRows == outRO->image->numRows, "mismatch between image and mask");
+
+                // CTEMASK parameters
+                float cteMin = psMetadataLookupF32(NULL, arguments, "CTE.MIN"); // Number of fringe points
+
+                char *cteMaskName = psMetadataLookupStr (&mdok, config->arguments, "MASK.SET.VALUE");
+                psImageMaskType cteMaskValue = pmConfigMaskGet(cteMaskName, config);
+
+                if (0) {
+                  psFits *fits = NULL;
+                  fits = psFitsOpen ("combine.fits", "w");
+                  psFitsWriteImage (fits, NULL, outRO->image, 0, NULL);
+                  psFitsClose (fits);
+
+                  fits = psFitsOpen ("inmask.fits", "w");
+                  psFitsWriteImage (fits, NULL, outRO->mask, 0, NULL);
+                  psFitsClose (fits);
                 }
+
+                psF32 **outputImage = outRO->image->data.F32;
+                psImageMaskType **outputMask = outRO->mask->data.PS_TYPE_IMAGE_MASK_DATA;
+                for (int iy = 0; iy < outRO->image->numRows; iy++) {
+                    for (int ix = 0; ix < outRO->image->numCols; ix++) {
+                        if (outputImage[iy][ix] < cteMin) {
+                            outputMask[iy][ix] |= cteMaskValue;
+                        }
+                    }
+                }
+
+                if (0) {
+                  psFits *fits = NULL;
+                  fits = psFitsOpen ("otmask.fits", "w");
+                  psFitsWriteImage (fits, NULL, outRO->mask, 0, NULL);
+                  psFitsClose (fits);
+                }
+
             }
 

branches/cnb_branches/cnb_branch_20090301/ppMerge/src/ppMergeMask.c

@@ -31 +31 @@
     int sample = psMetadataLookupS32(NULL, config->arguments, "SAMPLE"); ///< Size of sample for statistics
     bool chipStats = psMetadataLookupBool(&mdok, config->arguments, "MASK.CHIPSTATS"); ///< Statistics on chip?
-    float maskSuspect = psMetadataLookupF32(NULL, config->arguments, "MASK.SUSPECT"); ///< Threshold for suspect pixels
+
+    char *maskSuspectMode  = psMetadataLookupStr(NULL, config->arguments, "MASK.SUSPECT.MODE"); ///< Threshold for suspect pixels
+    if (strcasecmp(maskSuspectMode, "SIGMA") && strcasecmp(maskSuspectMode, "VALUE")) {
+        psError (PS_ERR_UNKNOWN, true, "Invalid choice for MASK.SUSPECT.MODE: %s\n", maskSuspectMode);
+        return false;
+    }
+    float maskSuspectSigma = psMetadataLookupF32(NULL, config->arguments, "MASK.SUSPECT.SIGMA"); ///< Threshold for suspect pixels
+    float maskSuspectMin   = psMetadataLookupF32(NULL, config->arguments, "MASK.SUSPECT.MIN"); ///< Threshold for suspect pixels
+    float maskSuspectMax   = psMetadataLookupF32(NULL, config->arguments, "MASK.SUSPECT.MAX"); ///< Threshold for suspect pixels
+
     float maskBad = psMetadataLookupF32(NULL, config->arguments, "MASK.BAD"); ///< Threshold for bad pixels
     pmMaskIdentifyMode maskMode = psMetadataLookupS32(NULL, config->arguments, "MASK.MODE"); ///< Mode for identifying bad pixels
@@ -111 +120 @@
             }
 
-            if (!ppMergeFileOpenInput(config, view, i)) {
+            if (!ppMergeFileOpenInput(config, inView, i)) {
                 psError(PS_ERR_UNKNOWN, false, "Unable to open file %d", i);
                 psFree(inView);
@@ -154 +163 @@
                 valueIndex = 0;
             }
-            for (int i = 0; i < num; i++) {
-                int pixel = numPix * psRandomUniform(rng);
-                int x = pixel % numCols;
-                int y = pixel / numCols;
-                if (mask && (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskValRaw)) continue;
-                if (outMask && (outMask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskValOut)) continue;
-                if (!isfinite(image->data.F32[y][x])) continue;
-
-                values->data.F32[valueIndex++] = image->data.F32[y][x];
-            }
-
-            if (!chipStats) {
-                values->n = valueIndex;
-                if (!psVectorStats(statistics, values, NULL, NULL, 0)) {
-                    psError(PS_ERR_UNKNOWN, false, "Unable to do statistics on readout.");
-                    psFree(inView);
-                    psFree(readout);
-                    goto MERGE_MASK_ERROR;
-                }
-
-                float mean = psStatsGetValue(statistics, meanStat);
-                float stdev = psStatsGetValue(statistics, stdevStat);
-
-                // this function increments the count for each suspect pixel in each input plane
-                // maskValRaw is used to test for valid input pixels
-                if (!pmMaskFlagSuspectPixels(outRO, readout, mean, stdev, maskSuspect, maskValRaw)) {
-                    psError(PS_ERR_UNKNOWN, false, "Unable to find suspect values in file %d", i);
-                    psFree(inView);
-                    psFree(readout);
-                    goto MERGE_MASK_ERROR;
-                }
-                pmCellFreeData(inCell);
-
-                if (!pmFPAfileIOChecks(config, inView, PM_FPA_AFTER)) {
-                    psFree(inView);
-                    psFree(readout);
-                    goto MERGE_MASK_ERROR;
-                }
-            }
+
+            if (!strcasecmp(maskSuspectMode, "VALUE")) {
+                // this function increments the count for each suspect pixel in each input plane
+                // maskValRaw is used to test for valid input pixels
+                if (!pmMaskFlagSuspectPixelsByValue(outRO, readout, maskSuspectMin, maskSuspectMax, maskValRaw)) {
+                    psError(PS_ERR_UNKNOWN, false, "Unable to find suspect values in file %d", i);
+                    psFree(inView);
+                    psFree(readout);
+                    goto MERGE_MASK_ERROR;
+                }
+                pmCellFreeData(inCell);
+
+                if (!pmFPAfileIOChecks(config, inView, PM_FPA_AFTER)) {
+                    psFree(inView);
+                    psFree(readout);
+                    goto MERGE_MASK_ERROR;
+                }
+            } else {
+                // extract a subset of pixels for stats measurement -- don't use pixels which are masked for this calculation
+                for (int i = 0; i < num; i++) {
+                    int pixel = numPix * psRandomUniform(rng);
+                    int x = pixel % numCols;
+                    int y = pixel / numCols;
+                    if (mask && (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskValRaw)) continue;
+                    if (outMask && (outMask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskValOut)) continue;
+                    if (!isfinite(image->data.F32[y][x])) continue;
+
+                    values->data.F32[valueIndex] = image->data.F32[y][x];
+                    valueIndex++;
+                }
+
+                // for per-readout stats, measure the stats and find the suspect pixels
+                if (!chipStats) {
+                    values->n = valueIndex;
+                    if (!psVectorStats(statistics, values, NULL, NULL, 0)) {
+                        psError(PS_ERR_UNKNOWN, false, "Unable to do statistics on readout.");
+                        psFree(inView);
+                        psFree(readout);
+                        goto MERGE_MASK_ERROR;
+                    }
+                    float mean = psStatsGetValue(statistics, meanStat);
+                    float stdev = psStatsGetValue(statistics, stdevStat);
+
+                    // this function increments the count for each suspect pixel in each input plane
+                    // maskValRaw is used to test for valid input pixels
+                    if (!pmMaskFlagSuspectPixelsBySigma(outRO, readout, mean, stdev, maskSuspectSigma, maskValRaw)) {
+                        psError(PS_ERR_UNKNOWN, false, "Unable to find suspect values in file %d", i);
+                        psFree(inView);
+                        psFree(readout);
+                        goto MERGE_MASK_ERROR;
+                    }
+                    pmCellFreeData(inCell);
+
+                    if (!pmFPAfileIOChecks(config, inView, PM_FPA_AFTER)) {
+                        psFree(inView);
+                        psFree(readout);
+                        goto MERGE_MASK_ERROR;
+                    }
+                }
+            }
             psFree(readout);
             psFree(outRO);
@@ -198 +228 @@
 
         // Additional run through cells if we want chip-level statistics
-        if (chipStats && valueIndex > 0) {
+        // only used for MASK.SUSPECT.MODE == SIGMA
+        if (!strcasecmp(maskSuspectMode, "SIGMA") && chipStats && valueIndex > 0) {
             values->n = valueIndex;
             if (!psVectorStats(statistics, values, NULL, NULL, 0)) {
@@ -221 +252 @@
                 float stdev = psStatsGetValue(statistics, stdevStat);
 
-                if (!pmMaskFlagSuspectPixels(outRO, readout, mean, stdev, maskSuspect, maskValRaw)) {
+                if (!pmMaskFlagSuspectPixelsBySigma(outRO, readout, mean, stdev, maskSuspectSigma, maskValRaw)) {
                     psError(PS_ERR_UNKNOWN, false, "Unable to find suspect values in file %d", i);
                     goto MERGE_MASK_ERROR;

branches/cnb_branches/cnb_branch_20090301/ppMerge/src/ppMergeVersion.c

@@ -21 +21 @@
 
 #include "ppMergeVersion.h"
+#include "ppMergeVersionDefinitions.h"
 
 #ifndef PPMERGE_VERSION
@@ -32 +33 @@
 #endif
 
-#define xstr(s) str(s)
-#define str(s) #s
-
 psString ppMergeVersion(void)
 {
     char *value = NULL;
-    psStringAppend(&value, "%s@%s", xstr(PPMERGE_BRANCH), xstr(PPMERGE_VERSION));
+    psStringAppend(&value, "%s@%s", PPMERGE_BRANCH, PPMERGE_VERSION);
     return value;
 }
@@ -44 +42 @@
 psString ppMergeSource(void)
 {
-    return psStringCopy(xstr(PPMERGE_SOURCE));
+    return psStringCopy(PPMERGE_SOURCE);
 }
 

@@ -7 +7 @@
 stamp-h1
 .libs
+ppMergeVersionDefinitions.h