Changeset 25986

Timestamp:

Nov 1, 2009, 4:00:09 PM (17 years ago)

Author:

eugene

Message:

add psphotVisualPlotChisq

File:

• trunk/psphot/src/psphotVisual.c (modified) (10 diffs)

trunk/psphot/src/psphotVisual.c

@@ -106 +106 @@
 }
 
-bool psphotVisualScaleImage (int kapaFD, psImage *inImage, const char *name, int channel) {
+bool psphotVisualScaleImage (int kapaFD, psImage *inImage, psImage *inMask, const char *name, int channel) {
 
     KiiImage image;
@@ -116 +116 @@
     psStats *stats = psStatsAlloc (PS_STAT_ROBUST_MEDIAN | PS_STAT_ROBUST_STDEV);
     psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS);
-    if (!psImageBackground(stats, NULL, inImage, NULL, 0, rng)) {
+    if (!psImageBackground(stats, NULL, inImage, inMask, 0xffff, rng)) {
         fprintf (stderr, "failed to get background values\n");
         return false;
@@ -172 +172 @@
 
     // psphotVisualShowMask (kapa, readout->mask, "mask", 2);
-    psphotVisualScaleImage (kapa, readout->variance, "variance", 1);
-    psphotVisualScaleImage (kapa, readout->image, "image", 0);
+    psphotVisualScaleImage (kapa, readout->variance, readout->mask, "variance", 1);
+    psphotVisualScaleImage (kapa, readout->image, readout->mask, "image", 0);
 
     pmVisualAskUser(NULL);
@@ -197 +197 @@
     }
 
-    psphotVisualScaleImage (kapa, backgnd->image, "backgnd", 2);
-    psphotVisualScaleImage (kapa, readout->image, "backsub", 0);
+    psphotVisualScaleImage (kapa, backgnd->image, readout->mask, "backgnd", 2);
+    psphotVisualScaleImage (kapa, readout->image, readout->mask, "backsub", 0);
 
     pmVisualAskUser(NULL);
@@ -495 +495 @@
     Ymin = PS_MAX(Ymin, -0.1);
 
+    // XXX test: hardwire plot limits
+    // Xmin = -0.1; Ymin = -0.1;
+    // Xmax = 20.1; Ymax = 20.1;
+
     // storage vectors for data to be plotted
     psVector *xBright = psVectorAllocEmpty (sources->n, PS_TYPE_F32);
@@ -654 +658 @@
 
         int nRegions = psMetadataLookupS32 (&status, analysis, "PSF.CLUMP.NREGIONS");
-        float PSF_CLUMP_NSIGMA = psMetadataLookupF32 (&status, analysis, "PSF_CLUMP_NSIGMA");
+        float PSF_CLUMP_NSIGMA = psMetadataLookupF32 (&status, recipe, "PSF_CLUMP_NSIGMA");
 
         graphdata.color = KapaColorByName ("blue");
@@ -1107 +1111 @@
     }
 
-    psphotVisualScaleImage (myKapa, outsat, "satstar", 2);
+    psphotVisualScaleImage (myKapa, outsat, NULL, "satstar", 2);
 
     pmVisualAskUser(NULL);
@@ -1862 +1866 @@
 
         int nRegions = psMetadataLookupS32 (&status, analysis, "PSF.CLUMP.NREGIONS");
-        float PSF_CLUMP_NSIGMA = psMetadataLookupF32 (&status, analysis, "PSF_CLUMP_NSIGMA");
+        float PSF_CLUMP_NSIGMA = psMetadataLookupF32 (&status, recipe, "PSF_CLUMP_NSIGMA");
 
         graphdata.color = KapaColorByName ("blue");
@@ -1929 +1933 @@
     if (myKapa == -1) return false;
 
-    psphotVisualScaleImage (myKapa, readout->image, "resid", 1);
+    psphotVisualScaleImage (myKapa, readout->image, readout->mask, "resid", 1);
 
     pmVisualAskUser(NULL);
@@ -2049 +2053 @@
     psFree (y);
     psFree (dy);
+
+    pmVisualAskUser(NULL);
+    return true;
+}
+
+bool psphotVisualPlotChisq (psArray *sources) {
+
+    Graphdata graphdata;
+
+    // XXX TEST if (!pmVisualIsVisual()) return true;
+
+    int myKapa = psphotKapaChannel (2);
+    if (myKapa == -1) return false;
+
+    KapaClearPlots (myKapa);
+    KapaInitGraph (&graphdata);
+
+    psVector *x = psVectorAllocEmpty (sources->n, PS_TYPE_F32);
+    psVector *y = psVectorAllocEmpty (sources->n, PS_TYPE_F32);
+
+    graphdata.xmin = +32.0;
+    graphdata.xmax = -32.0;
+    graphdata.ymin = +32.0;
+    graphdata.ymax = -32.0;
+
+    FILE *f = fopen ("chisq.dat", "w");
+
+    // construct the plot vectors
+    int n = 0;
+    for (int i = 0; i < sources->n; i++) {
+        pmSource *source = sources->data[i];
+        if (!source) continue;
+        if (source->type != PM_SOURCE_TYPE_STAR) continue;
+        if (!source->moments) continue;
+        if (!isfinite(source->moments->Sum)) continue;
+        if (!source->modelPSF) continue;
+        if (!isfinite(source->modelPSF->chisq)) continue;
+
+        x->data.F32[n] = -2.5*log10(source->moments->Sum);
+        y->data.F32[n] = source->modelPSF->chisq / source->modelPSF->nDOF;
+        graphdata.xmin = PS_MIN(graphdata.xmin, x->data.F32[n]);
+        graphdata.xmax = PS_MAX(graphdata.xmax, x->data.F32[n]);
+        graphdata.ymin = PS_MIN(graphdata.ymin, y->data.F32[n]);
+        graphdata.ymax = PS_MAX(graphdata.ymax, y->data.F32[n]);
+
+        fprintf (f, "%d %d %f %f\n", i, n, x->data.F32[n], y->data.F32[n]);
+
+        n++;
+    }
+    x->n = y->n = n;
+    fclose (f);
+
+    float range;
+    range = graphdata.xmax - graphdata.xmin;
+    graphdata.xmax += 0.05*range;
+    graphdata.xmin -= 0.05*range;
+    range = graphdata.ymax - graphdata.ymin;
+    graphdata.ymax += 0.05*range;
+    graphdata.ymin -= 0.05*range;
+
+    // XXX test
+    graphdata.xmin = -17.0;
+    graphdata.xmax =  -3.0;
+    graphdata.ymin =  -0.1;
+    graphdata.ymax = +10.1;
+
+    KapaSetLimits (myKapa, &graphdata);
+
+    KapaSetFont (myKapa, "helvetica", 14);
+    KapaBox (myKapa, &graphdata);
+    KapaSendLabel (myKapa, "PSF Mag", KAPA_LABEL_XM);
+    KapaSendLabel (myKapa, "ChiSq", KAPA_LABEL_YM);
+
+    graphdata.color = KapaColorByName ("black");
+    graphdata.ptype = 2;
+    graphdata.size = 0.5;
+    graphdata.style = 2;
+    KapaPrepPlot (myKapa, n, &graphdata);
+    KapaPlotVector (myKapa, n, x->data.F32, "x");
+    KapaPlotVector (myKapa, n, y->data.F32, "y");
+
+    psFree (x);
+    psFree (y);
 
     pmVisualAskUser(NULL);