Index: /branches/eam_branches/ipp-20230313/psphot/src/psphotChoosePSF.c
===================================================================
--- /branches/eam_branches/ipp-20230313/psphot/src/psphotChoosePSF.c	(revision 42710)
+++ /branches/eam_branches/ipp-20230313/psphot/src/psphotChoosePSF.c	(revision 42711)
@@ -114,4 +114,10 @@
     options->psfTrendNy = psMetadataLookupS32 (&status, recipe, "PSF.TREND.NY");
     assert (status);
+
+    // use poissonian errors or local-sky errors
+    options->fitPSFstarCoords = psMetadataLookupBool (&status, recipe, "PSF_FIT_STAR_COORDS");
+    if (!status) {
+      options->fitPSFstarCoords = false;
+    }
 
     // get the fixed PSF fit radius
Index: /branches/eam_branches/ipp-20230313/psphot/src/psphotMakeResiduals.c
===================================================================
--- /branches/eam_branches/ipp-20230313/psphot/src/psphotMakeResiduals.c	(revision 42710)
+++ /branches/eam_branches/ipp-20230313/psphot/src/psphotMakeResiduals.c	(revision 42711)
@@ -369,9 +369,12 @@
 
                 if (fabs(resid->Ro->data.F32[oy][ox]) < pixelSN*dRo/sqrt(nKeep)) {
+		  // fprintf (stderr, "mask : %d, %d : %f +/- %f (%f for %d) : %d\n", ox, oy, resid->Ro->data.F32[oy][ox], dRo, dRo/sqrt(nKeep), nKeep, resid->mask->data.PM_TYPE_RESID_MASK_DATA[oy][ox]);
                   resid->mask->data.PM_TYPE_RESID_MASK_DATA[oy][ox] = 1;
 		  resid->Ro->data.F32[oy][ox] = 0.0;
 		  resid->Rx->data.F32[oy][ox] = 0.0;
 		  resid->Ry->data.F32[oy][ox] = 0.0;
-                }
+                } else {
+		  // fprintf (stderr, "keep : %d, %d : %f +/- %f (%f for %d) : %d\n", ox, oy, resid->Ro->data.F32[oy][ox], dRo, dRo/sqrt(nKeep), nKeep, resid->mask->data.PM_TYPE_RESID_MASK_DATA[oy][ox]);
+		}
             }
         }
Index: /branches/eam_branches/ipp-20230313/psphot/src/psphotMaskReadout.c
===================================================================
--- /branches/eam_branches/ipp-20230313/psphot/src/psphotMaskReadout.c	(revision 42710)
+++ /branches/eam_branches/ipp-20230313/psphot/src/psphotMaskReadout.c	(revision 42711)
@@ -58,4 +58,26 @@
             return false;
         }
+    }
+
+    // XXX test to repair variance
+    if (true) {
+      psImage *im = readout->image;
+      psImage *wt = readout->variance;
+      float maxVar = 0;
+      // find the max value
+      for (int j = 0; j < im->numRows; j++) {
+        for (int i = 0; i < im->numCols; i++) {
+	    if (!isfinite(wt->data.F32[j][i])) continue;
+	    maxVar = PS_MAX (wt->data.F32[j][i], maxVar);
+        }
+      }
+      // set pixels with nan variance and non-nan image to maxVar
+      for (int j = 0; j < im->numRows; j++) {
+        for (int i = 0; i < im->numCols; i++) {
+	    if (!isfinite(im->data.F32[j][i])) continue;
+	    if (isfinite(wt->data.F32[j][i])) continue;
+	    wt->data.F32[j][i] = maxVar;
+        }
+      }
     }
 
Index: /branches/eam_branches/ipp-20230313/psphot/src/psphotRoughClass.c
===================================================================
--- /branches/eam_branches/ipp-20230313/psphot/src/psphotRoughClass.c	(revision 42710)
+++ /branches/eam_branches/ipp-20230313/psphot/src/psphotRoughClass.c	(revision 42711)
@@ -198,5 +198,9 @@
         }
 
-        psfClump = pmSourcePSFClump (NULL, region, sources, PSF_SN_LIM, PSF_CLUMP_GRID_SCALE, MOMENTS_SX_MAX, MOMENTS_SY_MAX, MOMENTS_SX_MIN, MOMENTS_SY_MIN, MOMENTS_AR_MAX);
+	// use clipped stats or robust stats to define the clump
+	bool PSF_CLUMP_USE_CLIPPED_STATS = psMetadataLookupF32(&status, recipe, "PSF_CLUMP_USE_CLIPPED_STATS");
+	if (!status) { PSF_CLUMP_USE_CLIPPED_STATS = true; }
+
+        psfClump = pmSourcePSFClump (NULL, region, sources, PSF_SN_LIM, PSF_CLUMP_GRID_SCALE, MOMENTS_SX_MAX, MOMENTS_SY_MAX, MOMENTS_SX_MIN, MOMENTS_SY_MIN, MOMENTS_AR_MAX, PSF_CLUMP_USE_CLIPPED_STATS);
 
         psMetadataAddF32 (regionMD, PS_LIST_TAIL, "PSF.CLUMP.X",  PS_META_REPLACE, "psf clump center", psfClump.X);
Index: /branches/eam_branches/ipp-20230313/psphot/src/psphotSourceStats.c
===================================================================
--- /branches/eam_branches/ipp-20230313/psphot/src/psphotSourceStats.c	(revision 42710)
+++ /branches/eam_branches/ipp-20230313/psphot/src/psphotSourceStats.c	(revision 42711)
@@ -506,4 +506,8 @@
 	MOMENTS_SY_MIN = 0.5;
     }
+
+    // use clipped stats or robust stats to define the clump
+    bool PSF_CLUMP_USE_CLIPPED_STATS = psMetadataLookupF32(&status, recipe, "PSF_CLUMP_USE_CLIPPED_STATS");
+    if (!status) { PSF_CLUMP_USE_CLIPPED_STATS = true; }
 
     // when we set the window, we are not attempting to measure spatial variations; we can use a somewhat higher S/N limit
@@ -597,5 +601,5 @@
 
         // determine the PSF parameters from the source moment values
-        pmPSFClump psfClump = pmSourcePSFClump (NULL, NULL, sourceSubset, PSF_SN_LIM, PSF_CLUMP_GRID_SCALE, MOMENTS_SX_MAX, MOMENTS_SY_MAX, MOMENTS_SX_MIN, MOMENTS_SY_MIN, MOMENTS_AR_MAX);
+        pmPSFClump psfClump = pmSourcePSFClump (NULL, NULL, sourceSubset, PSF_SN_LIM, PSF_CLUMP_GRID_SCALE, MOMENTS_SX_MAX, MOMENTS_SY_MAX, MOMENTS_SX_MIN, MOMENTS_SY_MIN, MOMENTS_AR_MAX, PSF_CLUMP_USE_CLIPPED_STATS);
         psLogMsg ("psphot", 3, "sigma guess (pix) %.1f, nStars: %d of %d in clump, nSigma: %5.2f, X,  Y: %f, %f (%f, %f)\n", sigma[i], psfClump.nStars, psfClump.nTotal, psfClump.nSigma, psfClump.X, psfClump.Y, sqrt(psfClump.X) / sigma[i], sqrt(psfClump.Y) / sigma[i]);