Changeset 30141 for branches/eam_branches/ipp-20101205/psphot/src/psphotReplaceUnfit.c

Timestamp:

Dec 21, 2010, 3:01:12 PM (16 years ago)

Author:

eugene

Message:

for psf-matched images: re-determine the psf; re-generate the models with the new psf; re-fit the fluxes; subtract all sources; as radial aperture fluxes are measured only replace the current source of interest; replace all sources before convolving to the next psf

File:

• branches/eam_branches/ipp-20101205/psphot/src/psphotReplaceUnfit.c (modified) (2 diffs)

branches/eam_branches/ipp-20101205/psphot/src/psphotReplaceUnfit.c

@@ -75 +75 @@
       pmSourceAdd (source, PM_MODEL_OP_FULL, maskVal);
     }
+
+    psphotVisualShowImage(readout);
     psLogMsg ("psphot.replace", PS_LOG_INFO, "replaced models for %ld objects: %f sec\n", sources->n, psTimerMark ("psphot.replace"));
     return true;
@@ -101 +103 @@
     return true;
 }
+
+// modify the sources to point at the corresponding pixels for the given filerule
+bool psphotRedefinePixels (pmConfig *config, const pmFPAview *view, const char *filerule)
+{
+    bool status = true;
+
+    // select the appropriate recipe information
+    psMetadata *recipe  = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE);
+    psAssert (recipe, "missing recipe?");
+
+    int num = psphotFileruleCount(config, filerule);
+
+    // loop over the available readouts
+    for (int i = 0; i < num; i++) {
+        if (!psphotRedefinePixelsReadout (config, view, filerule, i, recipe)) {
+            psError (PSPHOT_ERR_CONFIG, false, "failed to replace all sources for %s entry %d", filerule, i);
+            return false;
+        }
+    }
+    return true;
+}
+
+bool psphotRedefinePixelsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) {
+
+    bool status;
+    pmSource *source;
+
+    psTimerStart ("psphot.replace");
+
+    // find the currently selected readout
+    pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest
+    psAssert (file, "missing file?");
+
+    pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
+    psAssert (readout, "missing readout?");
+
+    // XXX the sources have already been copied (merge into here?)
+    pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS");
+    psAssert (detections, "missing detections?");
+
+    psArray *sources = detections->allSources;
+    psAssert (sources, "missing sources?");
+
+    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
+    psImageMaskType maskVal = psMetadataLookupImageMask(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
+    psAssert (maskVal, "missing mask value?");
+
+    for (int i = 0; i < sources->n; i++) {
+      source = sources->data[i];
+
+      // sources have not yet been subtracted in this image (but this flag may be raised)
+      source->tmpFlags &= ~PM_SOURCE_TMPF_SUBTRACTED;
+
+      float Xo = source->modelPSF->params->data.F32[PM_PAR_XPOS];
+      float Yo = source->modelPSF->params->data.F32[PM_PAR_YPOS];
+      float radius = source->modelPSF->fitRadius;
+      if (radius < 1) {
+          fprintf (stderr, "!");
+      }
+
+      pmSourceRedefinePixels (source, readout, Xo, Yo, radius, true);
+    }
+    return true;
+}
+
+// for now, let's store the detections on the readout->analysis for each readout
+bool psphotResetModels (pmConfig *config, const pmFPAview *view, const char *filerule)
+{
+    bool status = true;
+
+    // select the appropriate recipe information
+    psMetadata *recipe  = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE);
+    psAssert (recipe, "missing recipe?");
+
+    int num = psphotFileruleCount(config, filerule);
+
+    // loop over the available readouts
+    for (int i = 0; i < num; i++) {
+        if (!psphotResetModelsReadout (config, view, filerule, i, recipe)) {
+            psError (PSPHOT_ERR_CONFIG, false, "failed to replace all sources for %s entry %d", filerule, i);
+            return false;
+        }
+    }
+    return true;
+}
+
+bool psphotResetModelsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) {
+
+    bool status;
+    pmSource *source;
+
+    psTimerStart ("psphot.replace");
+
+    // find the currently selected readout
+    pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest
+    psAssert (file, "missing file?");
+
+    pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
+    psAssert (readout, "missing readout?");
+
+    // XXX the sources have already been copied (merge into here?)
+    pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS");
+    psAssert (detections, "missing detections?");
+
+    psArray *sources = detections->allSources;
+    psAssert (sources, "missing sources?");
+
+    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
+    psImageMaskType maskVal = psMetadataLookupImageMask(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
+    psAssert (maskVal, "missing mask value?");
+
+    pmPSF *psf = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.PSF");
+    psAssert (psf, "missing psf?");
+
+    for (int i = 0; i < sources->n; i++) {
+      source = sources->data[i];
+
+      // sources have not yet been subtracted in this image (but this flag may be raised)
+      source->tmpFlags &= ~PM_SOURCE_TMPF_SUBTRACTED;
+
+      // the guess central intensity comes from the peak:
+      float Io = source->peak->flux;
+      float Xo = source->modelPSF->params->data.F32[PM_PAR_XPOS];
+      float Yo = source->modelPSF->params->data.F32[PM_PAR_YPOS];
+      float radius = source->modelPSF->fitRadius;
+
+      // generate a model for this object with Io = 1.0
+      pmModel *modelPSF = pmModelFromPSFforXY(psf, Xo, Yo, Io);
+      if (modelPSF == NULL) {
+          psWarning ("Failed to determine PSF model for source at (%f,%f), skipping", Xo, Yo);
+          continue;
+      }
+
+      // set the source PSF model
+      psFree (source->modelPSF);
+      source->modelPSF = modelPSF;
+      source->modelPSF->fitRadius = radius;
+
+      pmSourceCacheModel (source, maskVal);  // ALLOC x14 (!)
+    }
+
+    // psphotSaveImage(NULL, readoutIn->image, "image.in.fits");
+    // psphotSaveImage(NULL, readoutOut->image, "image.out.sub.fits");
+    // psphotVisualShowImage (readout);
+
+    psLogMsg ("psphot.replace", PS_LOG_INFO, "subtracted models for %ld objects: %f sec\n", sources->n, psTimerMark ("psphot.replace"));
+    return true;
+}
+