Changeset 20083

Timestamp:

Oct 12, 2008, 4:05:13 PM (18 years ago)

Author:

eugene

Message:

move psf star dump to own function

File:

• trunk/psphot/src/psphotOutput.c (modified) (1 diff)

trunk/psphot/src/psphotOutput.c

@@ -226 +226 @@
     return header;
 }
+
+// XXX add args as needed
+bool psphotDumpPSFStars (pmReadout *readout, pmPSFtry *try, float radius, psMaskType maskVal, psMaskType markVal) {
+
+    psphotSaveImage (NULL, readout->image,  "rawstars.fits");
+
+    for (int i = 0; i < try->sources->n; i++) {
+        // masked for: bad model fit, outlier in parameters
+        if (try->mask->data.U8[i] & PSFTRY_MASK_ALL)
+            continue;
+
+        pmSource *source = try->sources->data[i];
+        float x = source->modelPSF->params->data.F32[PM_PAR_XPOS];
+        float y = source->modelPSF->params->data.F32[PM_PAR_YPOS];
+
+        // set the mask and subtract the PSF model
+        // XXX should we be using maskObj? should we be unsetting the mask?
+        // use pmModelSub because modelFlux has not been generated
+        assert (source->maskObj);
+        psImageKeepCircle (source->maskObj, x, y, radius, "OR", markVal);
+        pmModelSub (source->pixels, source->maskObj, source->modelPSF, PM_MODEL_OP_FULL, maskVal);
+        psImageKeepCircle (source->maskObj, x, y, radius, "AND", PS_NOT_U8(markVal));
+    }
+
+    FILE *f = fopen ("shapes.dat", "w");
+    for (int i = 0; i < try->sources->n; i++) {
+        psF32 inPar[10];  // must be psF32 to pmPSF_FitToModel
+
+        // masked for: bad model fit, outlier in parameters
+        if (try->mask->data.U8[i] & PSFTRY_MASK_ALL) continue;
+
+        pmSource *source = try->sources->data[i];
+        psF32 *outPar = source->modelEXT->params->data.F32;
+
+        psEllipseShape shape;
+
+        shape.sx  = outPar[PM_PAR_SXX] / M_SQRT2;
+        shape.sy  = outPar[PM_PAR_SYY] / M_SQRT2;
+        shape.sxy = outPar[PM_PAR_SXY];
+
+        psEllipsePol pol = pmPSF_ModelToFit (outPar);
+        inPar[PM_PAR_E0] = pol.e0;
+        inPar[PM_PAR_E1] = pol.e1;
+        inPar[PM_PAR_E2] = pol.e2;
+        pmPSF_FitToModel (inPar, 0.1);
+
+        psEllipseAxes axes1 = psEllipseShapeToAxes (shape, 20.0);
+        psEllipseAxes axes2 = psEllipsePolToAxes(pol, 0.1);
+
+        psEllipsePol pol2 = psEllipseAxesToPol (axes1);
+
+        fprintf (f, "%3d  %7.2f %7.2f  %7.4f %7.4f %7.4f  --  %7.4f %7.4f %7.4f  :  %7.4f %7.4f %7.4f  --  %7.4f %7.4f %7.4f : %7.4f %7.4f %6.1f : %7.4f %7.4f %6.1f\n",
+                 i, outPar[PM_PAR_XPOS], outPar[PM_PAR_YPOS],
+                 outPar[PM_PAR_SXX], outPar[PM_PAR_SXY], outPar[PM_PAR_SYY],
+                 pol.e0, pol.e1, pol.e2,
+                 pol2.e0, pol2.e1, pol2.e2,
+                 inPar[PM_PAR_SXX], inPar[PM_PAR_SXY], inPar[PM_PAR_SYY],
+                 axes1.major, axes1.minor, axes1.theta*PM_DEG_RAD,
+                 axes2.major, axes2.minor, axes2.theta*PM_DEG_RAD
+            );
+    }
+    fclose (f);
+
+    psphotSaveImage (NULL, readout->image,  "psfstars.fits");
+    pmSourcesWritePSFs (try->sources, "psfstars.dat");
+    pmSourcesWriteEXTs (try->sources, "extstars.dat", false);
+    // XXX need alternative output function
+    // psMetadata *psfData = pmPSFtoMetadata (NULL, try->psf);
+    // psMetadataConfigWrite (psfData, "psfmodel.dat");
+    psLogMsg ("psphot.choosePSF", PS_LOG_INFO, "wrote out psf-subtracted image, psf data, exiting\n");
+
+    return true;
+}