Changeset 12689

Timestamp:

Mar 29, 2007, 4:48:20 PM (19 years ago)

Author:

eugene

Message:

added deblend for saturated stars, fixed AddNoise

Location:

trunk/psphot/src

Files:

• Makefile.am (modified) (1 diff)
• psphot.h (modified) (1 diff)
• psphotAddNoise.c (modified) (4 diffs)
• psphotBasicDeblend.c (modified) (1 diff)
• psphotDeblendSatstars.c (added)
• psphotEvalPSF.c (modified) (4 diffs)
• psphotFindPeaks.c (modified) (2 diffs)
• psphotReadout.c (modified) (6 diffs)

trunk/psphot/src/Makefile.am

@@ -51 +51 @@
         psphotSourcePlots.c      \
         psphotRadialPlot.c       \
+        psphotDeblendSatstars.c  \
         psphotMosaicSubimage.c   \
         psphotAddNoise.c         \

trunk/psphot/src/psphot.h

@@ -111 +111 @@
 bool psphotSubWithTest (pmSource *source, bool useState);
 bool psphotSetState (pmSource *source, bool curState);
+bool psphotDeblendSatstars (psArray *sources, psMetadata *recipe);

trunk/psphot/src/psphotAddNoise.c

@@ -4 +4 @@
 
     bool status = false;
-    psTimerStart ("psphot");
+    psEllipseShape oldshape;
+    psEllipseShape newshape;
+    psEllipseAxes axes;
 
     PS_ASSERT (readout, false);
     PS_ASSERT (readout->parent, false);
     PS_ASSERT (readout->parent->concepts, false);
+
+    psTimerStart ("psphot");
 
     // increase variance by factor*(object noise):
@@ -15 +19 @@
 
     float FACTOR = psMetadataLookupF32 (&status, recipe, "NOISE.FACTOR");
+    PS_ASSERT (status, false);
+    float SIZE = psMetadataLookupF32 (&status, recipe, "NOISE.SIZE");
     PS_ASSERT (status, false);
 
@@ -39 +45 @@
         }
 
-        model->params->data.F32[PM_PAR_I0] *= FACTOR;
+        psF32 *PAR = model->params->data.F32;
+
+        // save original values
+        float oldI0  = PAR[PM_PAR_I0];
+        oldshape.sx  = PAR[PM_PAR_SXX];
+        oldshape.sy  = PAR[PM_PAR_SYY];
+        oldshape.sxy = PAR[PM_PAR_SXY];
+
+        // increase size and height of source
+        axes = psEllipseShapeToAxes (oldshape, 20.0);
+        axes.major *= SIZE;
+        axes.minor *= SIZE;
+        newshape = psEllipseAxesToShape (axes);
+        PAR[PM_PAR_I0]  = FACTOR*PS_SQR(oldI0);
+        PAR[PM_PAR_SXX] = newshape.sx;
+        PAR[PM_PAR_SYY] = newshape.sy;
+        PAR[PM_PAR_SXY] = newshape.sxy;
+
         if (add) {
             pmModelAdd (source->weight, source->mask, model, false, false);
@@ -45 +68 @@
             pmModelSub (source->weight, source->mask, model, false, false);
         }
+        
+        // restore original values
+        PAR[PM_PAR_I0]  = oldI0;        
+        PAR[PM_PAR_SXX] = oldshape.sx;
+        PAR[PM_PAR_SYY] = oldshape.sy;
+        PAR[PM_PAR_SXY] = oldshape.sxy;
     }
     if (add) {

trunk/psphot/src/psphotBasicDeblend.c

@@ -1 +1 @@
 # include "psphot.h"
 
-// 2006.02.07 : no leaks
 bool psphotBasicDeblend (psArray *sources, psMetadata *recipe) {
 

trunk/psphot/src/psphotEvalPSF.c

@@ -73 +73 @@
       case PM_MODEL_UNTRIED:
         source->mode &= ~PM_SOURCE_MODE_FITTED; 
+        if (source->mode & PM_SOURCE_MODE_SATSTAR) {
+            psTrace ("psphot", 5, "satstar was not fitted");
+        }
         return false;
       case PM_MODEL_BADARGS:
@@ -79 +82 @@
       default:
         source->mode |= PM_SOURCE_MODE_FAIL;
+        if (source->mode & PM_SOURCE_MODE_SATSTAR) {
+            psTrace ("psphot", 5, "satstar failed fit");
+        }
         return false;
     }
@@ -103 +109 @@
     if (model->params->data.F32[PM_PAR_I0] <= 0) {
         source->mode |= PM_SOURCE_MODE_FAIL;
+        if (source->mode & PM_SOURCE_MODE_SATSTAR) {
+            psTrace ("psphot", 5, "satstar failed fit (peak below 0)");
+        }
         return false;
     } 
@@ -130 +139 @@
     keep &= (SN > PSF_MIN_SN);
     keep &= (Chi < PSF_MAX_CHI);
+
+    if (source->mode & PM_SOURCE_MODE_SATSTAR) {
+        psTrace ("psphot", 5, "satstar fit results: %f, %f  %d :  %f %f : %f %f : %f %f\n", 
+                     model->params->data.F32[PM_PAR_XPOS], model->params->data.F32[PM_PAR_YPOS], keep,
+                 dSX, nSx, dSY, nSy, SN, Chi);
+    }
+
     if (keep) {
         if (source->mode & PM_SOURCE_MODE_PSFSTAR) {

trunk/psphot/src/psphotFindPeaks.c

@@ -39 +39 @@
     // optionally save example images under trace 
     if (psTraceGetLevel("psphot") > 5) {
-        psphotSaveImage (NULL, smooth_im, "imsmooth.fits");
-        psphotSaveImage (NULL, smooth_wt, "wtsmooth.fits");
+        char name[64];
+        sprintf (name, "imsmooth.v%d.fits", pass);
+        psphotSaveImage (NULL, smooth_im, name);
+        sprintf (name, "wtsmooth.v%d.fits", pass);
+        psphotSaveImage (NULL, smooth_wt, name);
     }
 
@@ -58 +61 @@
     // optionally save example images under trace 
     if (psTraceGetLevel("psphot") > 5) {
-        psphotSaveImage (NULL, smooth_im, "snsmooth.fits");
+        char name[64];
+        sprintf (name, "snsmooth.v%d.fits", pass);
+        psphotSaveImage (NULL, smooth_im, name);
     }
 

trunk/psphot/src/psphotReadout.c

@@ -61 +61 @@
         return psphotReadoutCleanup (config, readout, recipe, NULL, NULL);
     }
+    pmPeaksWriteText (peaks, "oldpeaks.dat");
 
     // construct sources and measure basic stats
@@ -70 +71 @@
         return psphotReadoutCleanup(config, readout, recipe, NULL, sources);
     }
+
+    psphotDeblendSatstars (sources, recipe);
 
     // mark blended peaks PS_SOURCE_BLEND
@@ -102 +105 @@
     psphotGuessModels (readout, sources, recipe, psf);
 
+    psphotSaveImage (NULL, readout->image,  "image.v0.fits");
+
     // linear PSF fit to peaks
     psphotEnsemblePSF (readout, sources, recipe, psf, FALSE);
+    psphotSaveImage (NULL, readout->image,  "image.v1.fits");
     if (!strcasecmp (breakPt, "ENSEMBLE")) {
         goto finish;
     }
 
-    // plot positive + negative sources (replace, then remove)
-    psphotSourcePlots (readout, sources, recipe);
-
     // non-linear PSF and EXT fit to brighter sources
     psphotBlendFit (readout, sources, recipe, psf);
+    psphotSaveImage (NULL, readout->image,  "image.v2.fits");
 
     // replace all sources
     psphotReplaceAll (sources);
-
-    // plot positive sources 
-    // psphotSourcePlots (readout, sources, recipe);
+    psphotSaveImage (NULL, readout->image,  "image.v3.fits");
 
     // linear PSF fit to remaining peaks
     psphotEnsemblePSF (readout, sources, recipe, psf, TRUE);
+    psphotSaveImage (NULL, readout->image,  "image.v4.fits");
     if (!strcasecmp (breakPt, "PASS1")) {
         goto finish;
@@ -137 +140 @@
     // find the peaks in the image
     psArray *newPeaks = psphotFindPeaks (readout, recipe, 2);
+    pmPeaksWriteText (newPeaks, "newpeaks.dat");
 
     // remove noise for subtracted objects
@@ -153 +157 @@
     // replace all sources
     psphotReplaceAll (sources);
+    psphotSaveImage (NULL, readout->image,  "image.v5.fits");
 
     // merge the newly selected peaks into the existing list
@@ -160 +165 @@
     // linear PSF fit to remaining peaks
     psphotEnsemblePSF (readout, sources, recipe, psf, TRUE);
+    psphotSaveImage (NULL, readout->image,  "image.v6.fits");
 
 finish:
+
+    // plot positive sources 
+    psphotSourcePlots (readout, sources, recipe);
 
     // measure aperture photometry corrections