Changeset 29937 for trunk

Timestamp:

Dec 5, 2010, 9:45:46 PM (15 years ago)

Author:

eugene

Message:

add forced photometry of positive sources; generate a background model for PPSUB.OUTPUT even if PSPHOT.BACKMDL exists (is certain to be for the wrong image)

Location:

trunk/ppSub

Files:

• doc/psphot.txt (copied) (copied from branches/eam_branches/ipp-20101103/ppSub/doc/psphot.txt )
• src/Makefile.am (modified) (2 diffs)
• src/ppSub.h (modified) (2 diffs)
• src/ppSubArguments.c (modified) (1 diff)
• src/ppSubBackground.c (modified) (2 diffs)
• src/ppSubCamera.c (modified) (3 diffs)
• src/ppSubFlagNeighbors.c (modified) (1 diff)
• src/ppSubInputDetections.c (copied) (copied from branches/eam_branches/ipp-20101103/ppSub/src/ppSubInputDetections.c )
• src/ppSubLoop.c (modified) (4 diffs)
• src/ppSubMakePSF.c (modified) (2 diffs)
• src/ppSubMatchPSFs.c (modified) (1 diff)
• src/ppSubReadoutForcedPhot.c (copied) (copied from branches/eam_branches/ipp-20101103/ppSub/src/ppSubReadoutForcedPhot.c )
• src/ppSubReadoutPhotometry.c (modified) (4 diffs)

trunk/ppSub/src/Makefile.am

@@ -38 +38 @@
         ppSubLoop.c                     \
         ppSubDefineOutput.c             \
+        ppSubInputDetections.c          \
         ppSubExtras.c                   \
         ppSubFlagNeighbors.c            \
@@ -46 +47 @@
         ppSubReadoutJpeg.c              \
         ppSubReadoutPhotometry.c        \
+        ppSubReadoutForcedPhot.c        \
         ppSubReadoutStats.c             \
         ppSubReadoutSubtract.c          \

trunk/ppSub/src/ppSub.h

@@ -45 +45 @@
     bool photometry;                    // Perform photometry?
     bool inverse;                       // Output inverse subtraction as well?
+    bool forcedPhot1;                   // perform forced photometry?
+    bool forcedPhot2;                   // perform forced photometry?
     bool saveInConv;                    // Save convolved input?
     bool saveRefConv;                   // Save convolved reference?
@@ -105 +107 @@
                             ppSubData *data ///< Processing data
     );
+
+bool ppSubInputDetections (bool *foundDetections, const char *sourcesName, const char *imageName, ppSubData *data);
+bool ppSubReadoutForcedPhot(const char *outputName, const char *targetName, const char *sourceName, ppSubData *data);
+bool psphotCopyResults (bool *foundDetections, pmFPAfile *target, pmFPAfile *source, pmFPAview *view);
 
 /// Higher-order background subtraction

trunk/ppSub/src/ppSubArguments.c

@@ -86 +86 @@
     psMetadataAddS32(arguments, PS_LIST_TAIL, "-convolve", 0, "Image to convolve [1 or 2]", 0);
     psMetadataAddBool(arguments, PS_LIST_TAIL, "-photometry", 0, "Perform photometry?", NULL);
+    psMetadataAddBool(arguments, PS_LIST_TAIL, "-forced-phot", 0, "Perform forced photometry?", NULL);
+    psMetadataAddBool(arguments, PS_LIST_TAIL, "-forced-input1", 0, "Perform forced photometry?", NULL);
+    psMetadataAddBool(arguments, PS_LIST_TAIL, "-forced-input2", 0, "Perform forced photometry?", NULL);
     psMetadataAddF32(arguments, PS_LIST_TAIL, "-zp", 0, "Zero point for photometry", NAN);
     psMetadataAddBool(arguments, PS_LIST_TAIL, "-inverse", 0, "Generate inverse subtractions?", false);

trunk/ppSub/src/ppSubBackground.c

@@ -37 +37 @@
     pmFPAview *view = ppSubViewReadout(); // View to readout
     pmReadout *outRO = pmFPAfileThisReadout(config->files, view, "PPSUB.OUTPUT"); // Output image
-    pmReadout *modelRO = pmFPAfileThisReadout(config->files, view, "PSPHOT.BACKMDL"); // Background model
 
-    // Generate the background model, if required
+    // Generate the background model
+    if (!psphotModelBackground(config, view, "PPSUB.OUTPUT")) {
+      psError(psErrorCodeLast(), false, "Unable to model background");
+      psFree(view);
+      return false;
+    }
+
+    // select the model readout (should now exist)
+    pmReadout *modelRO = pmFPAfileThisReadout(config->files, view, "PSPHOT.BACKMDL");
     if (!modelRO) {
-        // Create the background model
-        if (!psphotModelBackgroundReadoutFileIndex(config, view, "PPSUB.OUTPUT", 0)) {
-            psError(psErrorCodeLast(), false, "Unable to model background");
-            psFree(view);
-            return false;
-        }
-        // select the model readout (should now exist)
-        modelRO = pmFPAfileThisReadout(config->files, view, "PSPHOT.BACKMDL");
-        if (!modelRO) {
-            psError(psErrorCodeLast(), false, "Unable to find background model");
-            psFree(view);
-            return false;
-        }
+      psError(psErrorCodeLast(), false, "Unable to find background model");
+      psFree(view);
+      return false;
     }
     psFree(view);
 
-    psImageBinning *binning = psMetadataLookupPtr(&mdok, modelRO->analysis,
-                                                  "PSPHOT.BACKGROUND.BINNING"); // Binning for model
+    psImageBinning *binning = psMetadataLookupPtr(&mdok, modelRO->analysis, "PSPHOT.BACKGROUND.BINNING"); // Binning for model
     psImage *modelImage = modelRO->image; // Background model
     psImage *image = outRO->image; // Image of interest
@@ -83 +79 @@
     }
 
-    pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL");
-    pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL.STDEV");
-
     return true;
 }

trunk/ppSub/src/ppSubCamera.c

@@ -221 +221 @@
                           "Generate inverse subtractions?", true);
     }
-
-    data->inverse = psMetadataLookupBool(NULL, recipe, "INVERSE");
-    data->photometry = psMetadataLookupBool(NULL, recipe, "PHOTOMETRY");
+    if (psMetadataLookupBool(NULL, config->arguments, "-forced-phot")) {
+        psMetadataAddBool(recipe, PS_LIST_TAIL, "FORCED.PHOTOMETRY.BOTH", PS_META_REPLACE, "Perform forced photometry?", true);
+    }
+    if (psMetadataLookupBool(NULL, config->arguments, "-forced-input1")) {
+        psMetadataAddBool(recipe, PS_LIST_TAIL, "FORCED.PHOTOMETRY.INPUT1", PS_META_REPLACE, "Perform forced photometry?", true);
+    }
+    if (psMetadataLookupBool(NULL, config->arguments, "-forced-input2")) {
+        psMetadataAddBool(recipe, PS_LIST_TAIL, "FORCED.PHOTOMETRY.INPUT2", PS_META_REPLACE, "Perform forced photometry?", true);
+    }
+
+    data->inverse     = psMetadataLookupBool(NULL, recipe, "INVERSE");
+    data->photometry  = psMetadataLookupBool(NULL, recipe, "PHOTOMETRY");
+    data->forcedPhot1 = psMetadataLookupBool(NULL, recipe, "FORCED.PHOTOMETRY.BOTH") || psMetadataLookupBool(NULL, recipe, "FORCED.PHOTOMETRY.INPUT1");
+    data->forcedPhot2 = psMetadataLookupBool(NULL, recipe, "FORCED.PHOTOMETRY.BOTH") || psMetadataLookupBool(NULL, recipe, "FORCED.PHOTOMETRY.INPUT2");
 
     // Convolved input image
@@ -363 +374 @@
             return false;
         }
-        // specify the number of psphot input images
-        psMetadataAddS32 (config->arguments, PS_LIST_TAIL, "PSPHOT.INPUT.NUM", PS_META_REPLACE, "number of inputs", 1);
         pmFPAfileActivate(config->files, false, "PSPHOT.INPUT");
 
@@ -405 +414 @@
             invSources->save = true;
         }
+
+        // files need to do the forced photometry on the positions of sources in the positive images
+        if (data->forcedPhot1) {
+            // this pmFPAfile is used to carry sources detected in the positive image #1
+            pmFPAfile *posSources1 = defineOutputFile(config, output, true, "PPSUB.POS1.SOURCES", PM_FPA_FILE_CMF);
+            if (!posSources1) {
+                psError(psErrorCodeLast(), false, "Unable to set up forced source file.");
+                return false;
+            }
+            posSources1->save = true;
+
+            // this pmFPAfile is used to carry sources detected in the diff image @ the positions from positive image #1
+            pmFPAfile *frcSources1 = defineOutputFile(config, output, true, "PPSUB.FORCED1.SOURCES", PM_FPA_FILE_CMF);
+            if (!frcSources1) {
+                psError(psErrorCodeLast(), false, "Unable to set up forced source file.");
+                return false;
+            }
+            frcSources1->save = true;
+        }
+
+        if (data->forcedPhot2) {
+            // this pmFPAfile is used to carry sources detected in the positive image #2
+            pmFPAfile *posSources2 = defineOutputFile(config, inverse, true, "PPSUB.POS2.SOURCES", PM_FPA_FILE_CMF);
+            if (!posSources2) {
+                psError(psErrorCodeLast(), false, "Unable to set up forced source file.");
+                return false;
+            }
+            posSources2->save = true;
+
+            // this pmFPAfile is used to carry sources detected in the diff image @ the positions from positive image #2
+            pmFPAfile *frcSources2 = defineOutputFile(config, inverse, true, "PPSUB.FORCED2.SOURCES", PM_FPA_FILE_CMF);
+            if (!frcSources2) {
+                psError(psErrorCodeLast(), false, "Unable to set up forced source file.");
+                return false;
+            }
+            frcSources2->save = true;
+        }
     }
 

trunk/ppSub/src/ppSubFlagNeighbors.c

@@ -235 +235 @@
         src = sources->data[Imin]; 
 
-        fprintf (stderr, "j: %d, Imin: %d, obj x,y: %f, %f  src x,y: %f, %f, SN: %f, ID: %d\n", j, Imin, obj->x, obj->y, src->peak->xf, src->peak->yf, src->peak->SN, src->id);
+        // fprintf (stderr, "j: %d, Imin: %d, obj x,y: %f, %f  src x,y: %f, %f, SN: %f, ID: %d\n", j, Imin, obj->x, obj->y, src->peak->xf, src->peak->yf, src->peak->SN, src->id);
 
         // add source to object

trunk/ppSub/src/ppSubLoop.c

@@ -57 +57 @@
     }
 
+    if (data->forcedPhot1) {
+        bool foundDetections = false;
+        if (!ppSubInputDetections(&foundDetections, "PPSUB.POS1.SOURCES", "PPSUB.INPUT", data)) {
+            psError(psErrorCodeLast(), false, "Unable to measure positive detections (1)");
+            return false;
+        }
+        // if nothing was found, don't bother doing the forced photometry below
+        if (!foundDetections) {
+            psWarning ("no sources found in positive image 1, skipping forced photometry");
+            data->forcedPhot1 = false;
+        }
+    }
+    if (data->forcedPhot2) {
+        bool foundDetections = false;
+        if (!ppSubInputDetections(&foundDetections, "PPSUB.POS2.SOURCES", "PPSUB.REF", data)) {
+            psError(psErrorCodeLast(), false, "Unable to measure positive detections (2)");
+            return false;
+        }
+        // if nothing was found, don't bother doing the forced photometry below
+        if (!foundDetections) {
+            psWarning ("no sources found in positive image 1, skipping forced photometry");
+            data->forcedPhot1 = false;
+        }
+    }
+
+    // XXX if it exists, use the POS1, POS2 results for the FWHMs
     if (!ppSubMatchPSFs(data)) {
         psError(psErrorCodeLast(), false, "Unable to match PSFs.");
@@ -70 +96 @@
         return false;
     }
+
+    // XXX add in a positive image detection step here (if needed)
+    
 
     psMetadataAddF32(data->stats, PS_LIST_TAIL, "TIME_MATCH", 0, "Time to match PSFs",
@@ -153 +182 @@
     }
     // dumpout(config, "diff.3.fits");
+
+    // forced photometry for positive image 1
+    if (data->forcedPhot1 && !data->quality) {
+        if (!ppSubReadoutForcedPhot("PPSUB.FORCED1.SOURCES", "PPSUB.OUTPUT", "PPSUB.POS1.SOURCES", data)) {
+            psError(psErrorCodeLast(), false, "Unable to perform photometry.");
+            return false;
+        }
+    }
+
+    // forced photometry for positive image 2
+    if (data->forcedPhot2 && !data->quality) {
+        if (!ppSubReadoutForcedPhot("PPSUB.FORCED2.SOURCES", "PPSUB.OUTPUT", "PPSUB.POS2.SOURCES", data)) {
+            psError(psErrorCodeLast(), false, "Unable to perform photometry.");
+            return false;
+        }
+    }
 
     if (!ppSubFilesIterateUp(config, PPSUB_FILES_PHOT_SUB)) {
@@ -227 +272 @@
     }
 
+    pmFPAfileDropInternal(config->files, "PSPHOT.BACKGND");
+    pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL");
+    pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL.STDEV");
+
     return true;
 }

trunk/ppSub/src/ppSubMakePSF.c

@@ -35 +35 @@
 
     psTimerStart("PPSUB_PHOT");
-
-    psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSUB_RECIPE); // Recipe for ppSub
-    psAssert(recipe, "We checked this earlier, so it should be here.");
-
-    if (!psMetadataLookupBool(NULL, recipe, "PHOTOMETRY")) {
-        return true;
-    }
 
     bool reverse = psMetadataLookupBool(NULL, config->arguments, "REVERSE"); // Reverse sense of subtraction?
@@ -93 +86 @@
     // use flags to toss totally bogus entries?
     psArray *goodSources = ppSubSelectPSFSources (sources);
-    if (!psphotReadoutFindPSF(config, view, goodSources)) {
+    if (!psphotReadoutFindPSF(config, view, "PSPHOT.INPUT", goodSources)) {
         // This is likely a data quality issue
         // XXX Split into multiple cases using error codes?

trunk/ppSub/src/ppSubMatchPSFs.c

@@ -72 +72 @@
     // Extract the loaded sources from the associated readout, and generate PSF
     // Here, we assume the image is background-subtracted
-    if (!psphotReadoutFindPSF(config, view, sources)) {
+    if (!psphotReadoutFindPSF(config, view, "PSPHOT.INPUT", sources)) {
         psErrorStackPrint(stderr, "Unable to determine PSF");
         psWarning("Unable to determine PSF.");

trunk/ppSub/src/ppSubReadoutPhotometry.c

@@ -31 +31 @@
 
     if (!data->photometry) {
-        return true;
-    }
-
-    // Look up recipe values
-    psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSUB_RECIPE); // Recipe for ppSim
-    psAssert(recipe, "We checked this earlier, so it should be here.");
-
-    if (!psMetadataLookupBool(NULL, recipe, "PHOTOMETRY")) {
-        pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL");
-        pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL.STDEV");
         return true;
     }
@@ -81 +71 @@
     }
 
-    if (!psphotReadoutMinimal(config, view)) {
+    if (!psphotReadoutMinimal(config, view, "PSPHOT.INPUT")) {
         // This is likely a data quality issue
         // XXX Split into multiple cases using error codes?
@@ -90 +80 @@
 
     // If no sources were found, there's no error,  but we want to trigger 'bad quality'
+    psWarning("no sources found: why is this being set to bad quality??");
     pmDetections *detections = psMetadataLookupPtr(NULL, photRO->analysis, "PSPHOT.DETECTIONS"); // Sources
     if (!detections) {
@@ -136 +127 @@
         }
     }
-
-    pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL");
-    pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL.STDEV");
 
     return true;