Changeset 33094 for branches/eam_branches/ipp-20111122/psphot/src/psphotFitSourcesLinear.c

Timestamp:

Jan 11, 2012, 2:23:07 PM (15 years ago)

Author:

eugene

Message:

generate a background model of the variance image for use in the model-based variance; finish the model-based variance analysis

File:

• branches/eam_branches/ipp-20111122/psphot/src/psphotFitSourcesLinear.c (modified) (8 diffs)

branches/eam_branches/ipp-20111122/psphot/src/psphotFitSourcesLinear.c

@@ -12 +12 @@
 static bool SetBorderMatrixElements (psSparseBorder *border, pmReadout *readout, psArray *sources, bool constant_weights, int SKY_FIT_ORDER, psImageMaskType markVal);
 
-bool psphotGenerateModelVariance (pmConfig *config, const pmFPAview *view, int index, psMetadata *recipe, pmReadout *readout, psArray *sources);
+bool psphotGenerateModelVariance (pmConfig *config, const pmFPAview *view, pmFPAfile *file, int index, psMetadata *recipe, pmReadout *readout, psArray *sources);
 pmSourceFitVarMode psphotGetFitVarMode (psMetadata *recipe);
 bool psphotFreeModelVariance (pmReadout *readout, psArray *sources);
@@ -74 +74 @@
         if (fitVarMode == PM_SOURCE_PHOTFIT_MODEL_VAR) {
             // generate the model variance image & source pointers
-            if (!psphotGenerateModelVariance (config, view, i, recipe, readout, sources)) {
+            if (!psphotGenerateModelVariance (config, view, file, i, recipe, readout, sources)) {
                 psError (PSPHOT_ERR_CONFIG, false, "failed to fit sources (linear) for %s entry %d", filerule, i);
                 return false;
@@ -554 +554 @@
 }
 
-bool psphotGenerateModelVariance (pmConfig *config, const pmFPAview *view, int index, psMetadata *recipe, pmReadout *readout, psArray *sources) {
+bool psphotModelBackgroundReadout(psImage *model,  // Model image
+                                  psImage *modelStdev, // Model stdev image
+                                  psMetadata *analysis, // Analysis metadata for outputs
+                                  pmReadout *readout, // Readout for which to generate a background model
+                                  psImageBinning *binning, // Binning parameters
+                                  const pmConfig *config,// Configuration
+                                  bool useVarianceImage
+    );
+
+bool psphotGenerateModelVariance (pmConfig *config, const pmFPAview *view, pmFPAfile *file, int index, psMetadata *recipe, pmReadout *readout, psArray *sources) {
 
     bool status = false;
@@ -573 +582 @@
     psImage *modelVar = psImageCopy (NULL, readout->variance, PS_TYPE_F32);
 
-    // find the pmFPAfile for the readout (background model is saved on PSPHOT.BACKMDL regardless of 'filerule')
-    pmFPAfile *backModelFile = pmFPAfileSelectSingle(config->files, "PSPHOT.BACKMDL", index); // File of interest
-    assert (backModelFile);
-
-    // find the background model readout from the correct location
-    pmReadout *backModel = READOUT_OR_INTERNAL(view, backModelFile);
-    psAssert (backModel, "this must exist");
-
     // find the binning information
-    psImageBinning *backBinning = psMetadataLookupPtr(&status, backModel->analysis, "PSPHOT.BACKGROUND.BINNING");
+    psImageBinning *backBinning = psphotBackgroundBinning (modelVar, config);
     assert (backBinning);
+    
+    psImage *varModel = psImageAlloc(backBinning->nXruff, backBinning->nYruff, PS_TYPE_F32); // Background model
+    psImage *varModelStdev = psImageAlloc(backBinning->nXruff, backBinning->nYruff, PS_TYPE_F32); // Background model
+
+    if (!psphotModelBackgroundReadout(varModel, varModelStdev, NULL, readout, backBinning, config, true)) {
+        psError(PS_ERR_UNKNOWN, false, "Unable to generate background model");
+        psFree (varModel);
+        psFree (varModelStdev);
+        return false;
+    }
 
     // linear interpolation to full-scale
-    if (!psImageUnbin (modelVar, backModel->image, backBinning)) {
-        psError (PSPHOT_ERR_PROG, true, "inconsistent sizes for unbinning");
-        return false;
-    }
+    if (!psImageUnbin (modelVar, varModel, backBinning)) {
+        psError (PSPHOT_ERR_PROG, true, "inconsistent sizes for unbinning");
+        psFree (varModel);
+        psFree (varModelStdev);
+        return false;
+    }
+
+    psFree (varModel);
+    psFree (varModelStdev);
 
     // XXX for a test:
@@ -609 +625 @@
         psAssert (!source->modelVar, "programming error : modelVar should be NULL here");
         psAssert (source->modelFlux, "programming error : modelFlux should not be NULL here");
-        source->modelVar = psImageSubset(readout->image, region);
+        psAssert (source->modelFlux->data.F32, "programming error : modelFlux should not be NULL here");
+        source->modelVar = psImageSubset(modelVar, region);
 
         // add the source model to the model variance image
@@ -615 +632 @@
     }
 
-    // XXX add the readnoise to the image?  This is a bit problematic because I don't have the
-    // per-cell information.  Alternatively, I could use the median of the variance image or
-    // something equivalent to a background model to define the base level of the variance
-
     // XXX for a test:
     psphotSaveImage (NULL, modelVar, "model.var.fits");
-    psphotSaveImage (NULL, readout->image, "image.var.fits");
+    psphotSaveImage (NULL, readout->variance, "image.var.fits");
 
     // we save the model variance for future reference
@@ -638 +651 @@
     assert (modelVar);
 
+    psMetadataRemoveKey (readout->analysis, "PSPHOT.MODEL.VAR"); 
+
     // clear modelVar pointers for all of the source models
     for (int i = 0; i < sources->n; i++) {
@@ -645 +660 @@
         psFree (source->modelVar);
     }
-    psFree (modelVar);
 
     return true;