Index: /branches/eam_branches/ipp-20110710/psphot/src/psphotKronIterate.c
===================================================================
--- /branches/eam_branches/ipp-20110710/psphot/src/psphotKronIterate.c	(revision 32235)
+++ /branches/eam_branches/ipp-20110710/psphot/src/psphotKronIterate.c	(revision 32236)
@@ -37,5 +37,5 @@
         psAssert (detections, "missing detections?");
 
-        psArray *sources = detections->allSources;
+        psArray *sources = detections->newSources ? detections->newSources : detections->allSources;
         psAssert (sources, "missing sources?");
 
@@ -127,5 +127,5 @@
     // psphotVisualRangeImage (kapa, kronWindow, "kronwin", 0, 0.0, 1.0);
 
-    psphotSaveImage (NULL, kronWindow, "kron.window.v0.fits");
+    // psphotSaveImage (NULL, kronWindow, "kron.window.v0.fits");
 
     for (int j = 0; j < 5; j++) {
@@ -139,6 +139,8 @@
 
 	    // replace object in image
+	    bool reSubtract = false;
 	    if (source->tmpFlags & PM_SOURCE_TMPF_SUBTRACTED) {
 		pmSourceAdd (source, PM_MODEL_OP_FULL, maskVal);
+		reSubtract = true;
 	    }
 
@@ -151,11 +153,5 @@
 	    // clear the window function for this source based on the moments
 	    psphotKronWindowSetSource (source, kronWindow, (j > 0), false);
-	    // psphotKronWindowSetSource (source, kronWindow, false, false);
 	    // psphotVisualRangeImage (kapa, kronWindow, "kronwin", 1, 0.0, 1.0);
-
-	    // 165, 539;
-	    if ((fabs(source->peak->xf - 165) < 3) && (fabs(source->peak->yf - 539) < 3)) {
-		fprintf (stderr, "test obj\n");
-	    }
 
 	    // this function populates moments->Mrf,KronFlux,KronFluxErr
@@ -165,26 +161,13 @@
 	    // set a window function for each source based on the moments
 	    psphotKronWindowSetSource (source, kronWindow, true, true);
-	    // psphotKronWindowSetSource (source, kronWindow, false, true);
-
-	    // test source fluxes
-	    pmSourceMagnitudes (source, psf, photMode, maskVal, markVal, source->apRadius);
-	    float kmag = -2.5*log10(source->moments->KronFlux);
-# define TEST_X1 167
-# define TEST_Y1 299
-# define TEST_X2 180
-# define TEST_Y2 300
-	    if ((fabs(source->peak->xf - TEST_X1) < 3) && (fabs(source->peak->yf - TEST_Y1) < 3)) {
-		fprintf (stderr, "R1: %f vs %f  (%f) (%f)\n", source->moments->KronRadiusPSF, source->moments->Mrf, kmag, windowRadius);
+
+	    // if we subtracted it above, re-subtract the object, leave local sky
+	    if (reSubtract) {
+		pmSourceSub (source, PM_MODEL_OP_FULL, maskVal);
 	    }
-	    if ((fabs(source->peak->xf - TEST_X2) < 3) && (fabs(source->peak->yf - TEST_Y2) < 3)) {
-		fprintf (stderr, "R2: %f vs %f  (%f) (%f)\n", source->moments->KronRadiusPSF, source->moments->Mrf, kmag, windowRadius);
-	    }
-
-	    // re-subtract the object, leave local sky
-	    pmSourceSub (source, PM_MODEL_OP_FULL, maskVal);
 	}
-	char name[64];
-	sprintf (name, "kron.window.v%d.fits", j+1);
-	psphotSaveImage (NULL, kronWindow, name);
+	// char name[64];
+	// sprintf (name, "kron.window.v%d.fits", j+1);
+	// psphotSaveImage (NULL, kronWindow, name);
     }
     psFree (kronWindow);
@@ -275,9 +258,7 @@
 	    assert (z >= 0.0);
 
-	    // float weight  = 1.0;
+	    // weight by window image and wide Gaussian
 	    float weight1  = vWin[row+Ywo][col+Xwo]*exp(-z);
 	    float weight2  = vWin[yFlip+Ywo][xFlip+Xwo]*exp(-z);
-	    // float weight1  = vWin[row+Ywo][col+Xwo];
-	    // float weight2  = vWin[yFlip+Ywo][xFlip+Xwo];
 
 	    float fDiff1 = vPix[row][col]*weight1;
@@ -315,9 +296,4 @@
 	if (fabs(yDiff) > radKron) continue;
 
-	// coordinate of mirror pixel
-	int yFlip = yCM - yDiff;
-	if (yFlip < 0) continue;
-	if (yFlip >= source->pixels->numRows) continue;
-
 	for (psS32 col = 0; col < source->pixels->numCols ; col++) {
 	    // check mask and value for this pixel
@@ -328,30 +304,12 @@
 	    if (fabs(xDiff) > radKron) continue;
 
-	    // coordinate of mirror pixel
-	    int xFlip = xCM - xDiff;
-	    if (xFlip < 0) continue;
-	    if (xFlip >= source->pixels->numCols) continue;
-
-	    // check mask and value for mirror pixel
-	    if (vMsk && (vMsk[yFlip][xFlip] & maskVal)) continue;
-	    if (isnan(vPix[yFlip][xFlip])) continue;
-
 	    // radKron is just a function of (xDiff, yDiff)
 	    psF32 r2  = PS_SQR(xDiff) + PS_SQR(yDiff);
 	    if (r2 > radKron2) continue;
 
-	    // float z = r2 * rsigma2;
-	    // assert (z >= 0.0);
-
-	    // float weight  = 1.0;
-	    // float weight1  = vWin[row+Ywo][col+Xwo]*exp(-z);
-	    // float weight2  = vWin[yFlip+Ywo][xFlip+Xwo]*exp(-z);
 	    float weight1  = vWin[row+Ywo][col+Xwo];
-	    float weight2  = vWin[yFlip+Ywo][xFlip+Xwo];
-
 	    float fDiff1 = vPix[row][col]*weight1;
-	    float fDiff2 = vPix[yFlip][xFlip]*weight2;
-
-	    float pDiff = (fDiff1 > 0.0) ? sqrt(fabs(fDiff1*fDiff2)) : -sqrt(fabs(fDiff1*fDiff2));
+
+	    float pDiff = fDiff1;
 	    psF32 wDiff = vWgt[row][col] * weight1;
 
@@ -365,6 +323,5 @@
     source->moments->Mrf = Mrf;
     source->moments->KronFlux    = Sum;
-    source->moments->KronFluxErr = Sum * (nKronPix / Win);
-    // source->moments->KronFluxErr = sqrt(Var / Win);
+    source->moments->KronFluxErr = sqrt(Var);
 
     return true;
@@ -393,9 +350,6 @@
     float Mminor = 0.5*(Mxx + Myy) - 0.5*sqrt(PS_SQR(Mxx - Myy) + 4.0*PS_SQR(Mxy));
 
-    // float kratio = source->moments->KronFinner / source->moments->KronFlux;
-
+    // Mxx, Mxy, Myy define the elliptical shape, but Mrf defines the width 
     float scale = PS_SQR(0.5 * source->moments->Mrf) / Mmajor;
-    // float scale = useKronRadius ? 2.0 * source->moments->Mrf / Mmajor : 2.0;
-    // float scale = (kratio > 0.4) ? 9.0 * source->moments->Mrf / Mmajor : 3.0 * source->moments->Mrf / Mmajor;
 
     float Sxx = scale * Mmajor * Mminor / Myy; // sigma_x^2