Index: /branches/eam_branches/ipp-20110404/psphot/src/psphotKronMasked.c
===================================================================
--- /branches/eam_branches/ipp-20110404/psphot/src/psphotKronMasked.c	(revision 31314)
+++ /branches/eam_branches/ipp-20110404/psphot/src/psphotKronMasked.c	(revision 31314)
@@ -0,0 +1,111 @@
+# include "psphotInternal.h"
+
+bool psphotKronMasked (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);
+
+    // skip the chisq image (optionally?)
+    int chisqNum = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.CHISQ.NUM");
+    if (!status) chisqNum = -1;
+
+    // loop over the available readouts
+    for (int i = 0; i < num; i++) {
+        if (i == chisqNum) continue; // skip chisq image
+
+        // find the currently selected readout
+        pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, i); // File of interest
+        psAssert (file, "missing file?");
+
+        pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
+        psAssert (readout, "missing readout?");
+
+        pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS");
+        psAssert (detections, "missing detections?");
+
+        psArray *sources = detections->newSources;
+        psAssert (sources, "missing sources?");
+
+        pmPSF *psf = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.PSF");
+        // psAssert (psf, "missing psf?");
+
+        if (!psphotKronMaskedReadout (config, recipe, view, readout, sources, psf)) {
+            psError (PSPHOT_ERR_CONFIG, false, "failed to measure magnitudes for %s entry %d", filerule, i);
+            return false;
+        }
+    }
+    return true;
+}
+
+bool psphotKronMaskedReadout(pmConfig *config, psMetadata *recipe, const pmFPAview *view, pmReadout *readout, psArray *sources, pmPSF *psf) {
+
+    bool status = false;
+
+    if (!sources->n) {
+        psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping masked kron");
+        return true;
+    }
+
+    psTimerStart ("psphot.kron");
+
+    // determine the number of allowed threads
+    int nThreads = psMetadataLookupS32(&status, config->arguments, "NTHREADS"); // Number of threads
+    if (!status) {
+        nThreads = 0;
+    }
+
+    // bit-masks to test for good/bad pixels
+    psImageMaskType maskVal = psMetadataLookupImageMask(&status, recipe, "MASK.PSPHOT");
+    assert (maskVal);
+
+    // bit-mask to mark pixels not used in analysis
+    psImageMaskType markVal = psMetadataLookupImageMask(&status, recipe, "MARK.PSPHOT");
+    assert (markVal);
+
+    // maskVal is used to test for rejected pixels, and must include markVal
+    maskVal |= markVal;
+
+    // generate the mask image: increment counter for every source overlapping the pixel
+    psImage *kronMask = psImageAlloc (readout->image->numCols, readout->image->numRows, PS_TYPE_U8);
+    psImageInit (kronMask, 0);
+    int Nx = kronMask->numCols;
+    int Ny = kronMask->numRows;
+    for (int i = 0; i < sources->n; i++) {
+
+        pmSource *source = sources->data[i];
+        if (!source->peak) continue; // XXX how can we have a peak-less source?
+
+        // allocate space for moments
+        if (!source->moments) continue;
+
+	// generate the pixel masks
+	int Xo = source->moments->Mx;
+	int Yo = source->moments->My;
+	int Kr = 2.5*source->moments->Mrf;
+	int Kr2 = Kr*Kr;
+	
+	for (int iy = Yo - Kr; iy < Yo + Kr + 1; iy++) {
+	    if (iy < 0) continue;
+	    if (iy >= Ny) continue;
+	    for (int ix = Xo - Kr; ix < Xo + Kr + 1; ix++) {
+		if (ix < 0) continue;
+		if (ix >= Nx) continue;
+
+		if (PS_SQR(ix - Xo) + PS_SQR(iy - Yo) > Kr2) continue;
+		if (kronMask->data.U8[iy][ix] < 0xff) {
+		    kronMask->data.U8[iy][ix] ++;
+		}
+	    }
+	}
+    }		
+    psphotSaveImage (NULL, kronMask, "kronmask.fits");
+    psLogMsg ("psphot.kron", PS_LOG_DETAIL, "measure masked kron magnitudes : %f sec for %ld objects\n", psTimerMark ("psphot.kron"), sources->n);
+
+    psFree (kronMask);
+    return true;
+}