Index: /trunk/psphot/src/psphotApResid.c
===================================================================
--- /trunk/psphot/src/psphotApResid.c	(revision 9733)
+++ /trunk/psphot/src/psphotApResid.c	(revision 9734)
@@ -4,5 +4,5 @@
 
 // measure the aperture residual statistics
-bool psphotApResid (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf) { 
+bool psphotApResid (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf) {
 
     int Nfail = 0;
@@ -18,5 +18,5 @@
     float MAX_AP_OFFSET = psMetadataLookupF32 (&status, recipe, "MAX_AP_OFFSET");
 
-    // set limits on the aperture magnitudes 
+    // set limits on the aperture magnitudes
     pmSourceMagnitudesInit (recipe);
 
@@ -25,74 +25,73 @@
     psf->growth = pmGrowthCurveAlloc (3.0, 100.0, REF_RADIUS);
     psphotGrowthCurve (readout, psf);
-    
-    psVector *mask    = psVectorAlloc (300, PS_TYPE_U8);
-    psVector *xPos    = psVectorAlloc (300, PS_TYPE_F64);
-    psVector *yPos    = psVectorAlloc (300, PS_TYPE_F64);
-    psVector *flux    = psVectorAlloc (300, PS_TYPE_F64);
-    psVector *r2rflux = psVectorAlloc (300, PS_TYPE_F64);
-    psVector *apResid = psVectorAlloc (300, PS_TYPE_F64);
-    psVector *dMag    = psVectorAlloc (300, PS_TYPE_F64);
-    mask->n = xPos->n = yPos->n = flux->n = r2rflux->n = apResid->n = dMag->n = 0;
+
+    psVector *mask    = psVectorAllocEmpty (300, PS_TYPE_U8);
+    psVector *xPos    = psVectorAllocEmpty (300, PS_TYPE_F64);
+    psVector *yPos    = psVectorAllocEmpty (300, PS_TYPE_F64);
+    psVector *flux    = psVectorAllocEmpty (300, PS_TYPE_F64);
+    psVector *r2rflux = psVectorAllocEmpty (300, PS_TYPE_F64);
+    psVector *apResid = psVectorAllocEmpty (300, PS_TYPE_F64);
+    psVector *dMag    = psVectorAllocEmpty (300, PS_TYPE_F64);
     Npsf = 0;
 
     // select all good PM_SOURCE_TYPE_STAR entries
     for (int i = 0; i < sources->n; i++) {
-	source = sources->data[i];
-	model = source->modelPSF;
-
-	if (source->type != PM_SOURCE_TYPE_STAR) continue; 
-	if (source->mode &  PM_SOURCE_MODE_SATSTAR) continue;
-	if (source->mode &  PM_SOURCE_MODE_BLEND) continue;
-	if (source->mode &  PM_SOURCE_MODE_FAIL) continue;
-	if (source->mode &  PM_SOURCE_MODE_POOR) continue;
-
-	// get growth-corrected, apTrend-uncorrected magnitudes in scaled apertures
-	// will fail if below S/N threshold or model is missing
-	if (!pmSourceMagnitudes (source, psf, PM_SOURCE_PHOT_GROWTH)) {
-	    Nskip ++;
-	    continue;
-	}
-
-	apResid->data.F64[Npsf] = source->apMag - source->psfMag;
-
-	xPos->data.F64[Npsf]    = model->params->data.F32[PM_PAR_XPOS];
-	yPos->data.F64[Npsf]    = model->params->data.F32[PM_PAR_YPOS];
-
-	flux->data.F64[Npsf]    = pow(10.0, -0.4*source->psfMag);
-	r2rflux->data.F64[Npsf] = PS_SQR(model->radiusTMP) / flux->data.F64[Npsf];
-
-	mask->data.U8[Npsf] = 0;
-
-	// XXX sanity clip?
-	// XXX need to see if all data were tossed?
-	// XXX need to subtract median?
-	// XXX need to put this in the config data...
-	// if (fabs(apResid->data.F64[Npsf]) > 0.2) continue;
-	if ((MAX_AP_OFFSET > 0) && (fabs(apResid->data.F64[Npsf]) > MAX_AP_OFFSET)) {
-	    Nfail ++;
-	    continue;
-	}
-
-	dMag->data.F64[Npsf] = model->dparams->data.F32[PM_PAR_I0] / model->params->data.F32[PM_PAR_I0];
-
-	psVectorExtend (mask, 	 100, 1);
-	psVectorExtend (xPos, 	 100, 1);
-	psVectorExtend (yPos, 	 100, 1);
-	psVectorExtend (flux,    100, 1);
-	psVectorExtend (r2rflux, 100, 1);
-	psVectorExtend (dMag,    100, 1);
-	psVectorExtend (apResid, 100, 1);
-	Npsf ++;
-    }
-    psLogMsg ("psphot.apresid", 4, "measure aperture residuals : %f sec for %d objects (%d skipped, %d failed, %ld invalid)\n", 
-	      psTimerMark ("psphot"), Npsf, Nskip, Nfail, sources->n - Npsf - Nskip - Nfail);
+        source = sources->data[i];
+        model = source->modelPSF;
+
+        if (source->type != PM_SOURCE_TYPE_STAR) continue;
+        if (source->mode &  PM_SOURCE_MODE_SATSTAR) continue;
+        if (source->mode &  PM_SOURCE_MODE_BLEND) continue;
+        if (source->mode &  PM_SOURCE_MODE_FAIL) continue;
+        if (source->mode &  PM_SOURCE_MODE_POOR) continue;
+
+        // get growth-corrected, apTrend-uncorrected magnitudes in scaled apertures
+        // will fail if below S/N threshold or model is missing
+        if (!pmSourceMagnitudes (source, psf, PM_SOURCE_PHOT_GROWTH)) {
+            Nskip ++;
+            continue;
+        }
+
+        apResid->data.F64[Npsf] = source->apMag - source->psfMag;
+
+        xPos->data.F64[Npsf]    = model->params->data.F32[PM_PAR_XPOS];
+        yPos->data.F64[Npsf]    = model->params->data.F32[PM_PAR_YPOS];
+
+        flux->data.F64[Npsf]    = pow(10.0, -0.4*source->psfMag);
+        r2rflux->data.F64[Npsf] = PS_SQR(model->radiusTMP) / flux->data.F64[Npsf];
+
+        mask->data.U8[Npsf] = 0;
+
+        // XXX sanity clip?
+        // XXX need to see if all data were tossed?
+        // XXX need to subtract median?
+        // XXX need to put this in the config data...
+        // if (fabs(apResid->data.F64[Npsf]) > 0.2) continue;
+        if ((MAX_AP_OFFSET > 0) && (fabs(apResid->data.F64[Npsf]) > MAX_AP_OFFSET)) {
+            Nfail ++;
+            continue;
+        }
+
+        dMag->data.F64[Npsf] = model->dparams->data.F32[PM_PAR_I0] / model->params->data.F32[PM_PAR_I0];
+
+        psVectorExtend (mask,    100, 1);
+        psVectorExtend (xPos,    100, 1);
+        psVectorExtend (yPos,    100, 1);
+        psVectorExtend (flux,    100, 1);
+        psVectorExtend (r2rflux, 100, 1);
+        psVectorExtend (dMag,    100, 1);
+        psVectorExtend (apResid, 100, 1);
+        Npsf ++;
+    }
+    psLogMsg ("psphot.apresid", 4, "measure aperture residuals : %f sec for %d objects (%d skipped, %d failed, %ld invalid)\n",
+              psTimerMark ("psphot"), Npsf, Nskip, Nfail, sources->n - Npsf - Nskip - Nfail);
 
     // XXX choose a better value here?
     if (Npsf < 5) {
-	psError(PSPHOT_ERR_APERTURE, true, "no valid aperture residual sources, giving up");
-	return false;
-    }
-
-    // APTREND options : NONE SKYBIAS XY_LIN XY_QUAD SKY_XY_LIN FULL 
+        psError(PSPHOT_ERR_APERTURE, true, "no valid aperture residual sources, giving up");
+        return false;
+    }
+
+    // APTREND options : NONE SKYBIAS XY_LIN XY_QUAD SKY_XY_LIN FULL
     // APTREND options are used in the switch block below
     pmPSFApTrendOptions ApTrendOption = DEFAULT_OPTION;
@@ -100,6 +99,6 @@
     if (status) ApTrendOption = pmPSFApTrendOptionFromName (optionName);
     if (ApTrendOption == PM_PSF_APTREND_ERROR) {
-	psError(PSPHOT_ERR_APERTURE, true, "invalid aperture residual trend %s", optionName);
-	return false;
+        psError(PSPHOT_ERR_APERTURE, true, "invalid aperture residual trend %s", optionName);
+        return false;
     }
 
@@ -112,177 +111,177 @@
     switch (ApTrendOption) {
       case PM_PSF_APTREND_NONE:
-	// remove ApTrend fit from pmPSFtry
-	psf->ApTrend->coeff[0][0][0][0] = 0;
-	break;
+        // remove ApTrend fit from pmPSFtry
+        psf->ApTrend->coeff[0][0][0][0] = 0;
+        break;
       case PM_PSF_APTREND_CONSTANT:
-	stats->clipIter = 2;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
-	psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
-	    return false;
-	}
-	// apply the fit
-	stats->clipIter = 3;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
-	psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "Fitting aperture correction");
-	    return false;
-	}
-	break;
+        stats->clipIter = 2;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
+        psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
+            return false;
+        }
+        // apply the fit
+        stats->clipIter = 3;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
+        psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "Fitting aperture correction");
+            return false;
+        }
+        break;
       case PM_PSF_APTREND_SKYBIAS:
-	// first clip out objects which are too far from the median 
-	stats->clipIter = 2;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
-	psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
-	    return false;
-	}
-	// apply the fit
-	stats->clipIter = 3;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKYBIAS);
-	psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting sky bias");
-	    return false;
-	}
-	break;
+        // first clip out objects which are too far from the median
+        stats->clipIter = 2;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
+        psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
+            return false;
+        }
+        // apply the fit
+        stats->clipIter = 3;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKYBIAS);
+        psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting sky bias");
+            return false;
+        }
+        break;
       case PM_PSF_APTREND_SKYSAT:
-	// first clip out objects which are too far from the median 
-	stats->clipIter = 2;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
-	psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
-	    return false;
-	}
-	// apply the fit
-	stats->clipIter = 2;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKYBIAS);
-	psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting sky bias");
-	    return false;
-	}
-	// apply the fit
-	stats->clipIter = 3;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKYSAT);
-	psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting skysat");
-	    return false;
-	}
-	break;
+        // first clip out objects which are too far from the median
+        stats->clipIter = 2;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
+        psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
+            return false;
+        }
+        // apply the fit
+        stats->clipIter = 2;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKYBIAS);
+        psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting sky bias");
+            return false;
+        }
+        // apply the fit
+        stats->clipIter = 3;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKYSAT);
+        psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting skysat");
+            return false;
+        }
+        break;
       case PM_PSF_APTREND_XY_LIN:
-	// first clip out objects which are too far from the median 
-	stats->clipIter = 2;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
-	psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
-	    return false;
-	}
-	// apply the fit
-	stats->clipIter = 3;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_XY_LIN);
-	psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "fitting, XY_LIN");
-	    return false;
-	}
-	break;
+        // first clip out objects which are too far from the median
+        stats->clipIter = 2;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
+        psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
+            return false;
+        }
+        // apply the fit
+        stats->clipIter = 3;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_XY_LIN);
+        psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "fitting, XY_LIN");
+            return false;
+        }
+        break;
       case PM_PSF_APTREND_XY_QUAD:
-	// first clip out objects which are too far from the median 
-	stats->clipIter = 2;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
-	psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
-	    return false;
-	}
-	// apply the fit
-	stats->clipIter = 3;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_XY_QUAD);
-	psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "Fitting XY_QUAD");
-	    return false;
-	}
-	break;
+        // first clip out objects which are too far from the median
+        stats->clipIter = 2;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
+        psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
+            return false;
+        }
+        // apply the fit
+        stats->clipIter = 3;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_XY_QUAD);
+        psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "Fitting XY_QUAD");
+            return false;
+        }
+        break;
       case PM_PSF_APTREND_SKY_XY_LIN:
-	// first clip out objects which are too far from the median 
-	stats->clipIter = 2;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
-	psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
-	    return false;
-	}
-	// apply the fit
-	stats->clipIter = 3;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKY_XY_LIN);
-	psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "Fitting sky xy_lin");
-	    return false;
-	}
-	break;
+        // first clip out objects which are too far from the median
+        stats->clipIter = 2;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
+        psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
+            return false;
+        }
+        // apply the fit
+        stats->clipIter = 3;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKY_XY_LIN);
+        psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "Fitting sky xy_lin");
+            return false;
+        }
+        break;
       case PM_PSF_APTREND_SKYSAT_XY_LIN:
-	// first clip out objects which are too far from the median 
-	stats->clipIter = 2;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
-	psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
-	    return false;
-	}
-	// apply the fit
-	stats->clipIter = 3;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKYBIAS);
-	psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting sky bias");
-	    return false;
-	}
-	// apply the fit
-	stats->clipIter = 3;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKYSAT_XY_LIN);
-	psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "Fitting skyset xy_lin");
-	    return false;
-	}
-	break;
+        // first clip out objects which are too far from the median
+        stats->clipIter = 2;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
+        psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting nothing");
+            return false;
+        }
+        // apply the fit
+        stats->clipIter = 3;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKYBIAS);
+        psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "clipping, fitting sky bias");
+            return false;
+        }
+        // apply the fit
+        stats->clipIter = 3;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKYSAT_XY_LIN);
+        psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "Fitting skyset xy_lin");
+            return false;
+        }
+        break;
       case PM_PSF_APTREND_ALL:
-	// first clip out objects which are too far from the median 
-	stats->clipIter = 2;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
-	psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "Failed to measure apTrend");
-	    return false;
-	}
-	// fit just SkyBias and clip out objects which are too far from the median 
-	stats->clipIter = 2;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKYBIAS);
-	psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "fitting skyBias");
-	    return false;
-	}
-	// finally, fit x, y, SkyBias and clip out objects which are too far from the median 
-	stats->clipIter = 3;
-	pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_ALL);
-	psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
-	if (psf->ApTrend == NULL) {
-	    psError(PSPHOT_ERR_PHOTOM, false, "fitting all");
-	    return false;
-	}
-	break;
+        // first clip out objects which are too far from the median
+        stats->clipIter = 2;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_CONSTANT);
+        psf->ApTrend  = psVectorClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "Failed to measure apTrend");
+            return false;
+        }
+        // fit just SkyBias and clip out objects which are too far from the median
+        stats->clipIter = 2;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_SKYBIAS);
+        psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "fitting skyBias");
+            return false;
+        }
+        // finally, fit x, y, SkyBias and clip out objects which are too far from the median
+        stats->clipIter = 3;
+        pmPSFMaskApTrend (psf->ApTrend, PM_PSF_APTREND_ALL);
+        psf->ApTrend  = psVectorChiClipFitPolynomial4D (psf->ApTrend, stats, mask, PSFTRY_MASK_ALL, apResid, dMag, xPos, yPos, r2rflux, flux);
+        if (psf->ApTrend == NULL) {
+            psError(PSPHOT_ERR_PHOTOM, false, "fitting all");
+            return false;
+        }
+        break;
       default:
-	psError(PSPHOT_ERR_PHOTOM, true, "Unknown APTREND value: %s", optionName);
-	return false;
+        psError(PSPHOT_ERR_PHOTOM, true, "Unknown APTREND value: %s", optionName);
+        return false;
     }
 
@@ -295,8 +294,8 @@
     psStats *residStats = psStatsAlloc (PS_STAT_SAMPLE_STDEV);
     for (int i = 0; i < dMag->n; i++) {
-	if (dMag->data.F64[i] > 0.01) {
-	    mask->data.U8[i] |= 0x02;
-	} 
-	if (! mask->data.U8[i]) Nkeep ++;
+        if (dMag->data.F64[i] > 0.01) {
+            mask->data.U8[i] |= 0x02;
+        }
+        if (! mask->data.U8[i]) Nkeep ++;
     }
     residStats  = psVectorStats (residStats, resid, NULL, mask, 0x03);
@@ -320,12 +319,12 @@
 
     psLogMsg ("psphot.apresid", 3, "measure full-frame aperture residual: %f sec\n", psTimerMark ("psphot"));
-    psLogMsg ("psphot.apresid", 4, "aperture residual: %f +/- %f : %f bias, %f skysat (%d of %d used)\n", 
-	      psf->ApResid, psf->dApResid, psf->skyBias, psf->skySat, Nkeep, Npsf);
-    psLogMsg ("psphot.apresid", 4, "apresid trends: %f %f %f %f %f\n", 
-	      1e3*psf->ApTrend->coeff[1][0][0][0],
-	      1e6*psf->ApTrend->coeff[2][0][0][0],
-	      1e6*psf->ApTrend->coeff[1][1][0][0],
-	      1e3*psf->ApTrend->coeff[0][1][0][0],
-	      1e6*psf->ApTrend->coeff[0][2][0][0]);
+    psLogMsg ("psphot.apresid", 4, "aperture residual: %f +/- %f : %f bias, %f skysat (%d of %d used)\n",
+              psf->ApResid, psf->dApResid, psf->skyBias, psf->skySat, Nkeep, Npsf);
+    psLogMsg ("psphot.apresid", 4, "apresid trends: %f %f %f %f %f\n",
+              1e3*psf->ApTrend->coeff[1][0][0][0],
+              1e6*psf->ApTrend->coeff[2][0][0][0],
+              1e6*psf->ApTrend->coeff[1][1][0][0],
+              1e3*psf->ApTrend->coeff[0][1][0][0],
+              1e6*psf->ApTrend->coeff[0][2][0][0]);
 
     psFree (mask);
Index: /trunk/psphot/src/psphotBasicDeblend.c
===================================================================
--- /trunk/psphot/src/psphotBasicDeblend.c	(revision 9733)
+++ /trunk/psphot/src/psphotBasicDeblend.c	(revision 9734)
@@ -2,5 +2,5 @@
 
 // 2006.02.07 : no leaks
-bool psphotBasicDeblend (psArray *sources, psMetadata *recipe) { 
+bool psphotBasicDeblend (psArray *sources, psMetadata *recipe) {
 
     int N;
@@ -25,8 +25,7 @@
     // we use an index for this so the spatial sorting is kept
     psVector *SN = psVectorAlloc (sources->n, PS_DATA_F32);
-    SN->n = SN->nalloc;
     for (int i = 0; i < SN->n; i++) {
-	source = sources->data[i];
-	SN->data.F32[i] = source->peak->counts;
+        source = sources->data[i];
+        SN->data.F32[i] = source->peak->counts;
     }
     psVector *index = psVectorSortIndex (NULL, SN);
@@ -35,91 +34,89 @@
     // examine sources in decreasing SN order
     for (int i = sources->n - 1; i >= 0; i--) {
-	N = index->data.U32[i];
-	source = sources->data[N];
+        N = index->data.U32[i];
+        source = sources->data[N];
 
-	if (source->mode & PM_SOURCE_MODE_BLEND) continue;
+        if (source->mode & PM_SOURCE_MODE_BLEND) continue;
 
-	// temporary array for overlapping objects we find
-	psArray *overlap = psArrayAlloc (100);
-	// DROP overlap->n = 0;
+        // temporary array for overlapping objects we find
+        psArray *overlap = psArrayAllocEmpty (100);
 
-	// search backwards for overlapping sources
-	for (int j = N - 1; j >= 0; j--) {
-	    testSource = sources->data[j];
-	    if (testSource->peak->x <  source->pixels->col0) continue;
-	    if (testSource->peak->x >= source->pixels->col0 + source->pixels->numCols) continue;
-	    if (testSource->peak->y <  source->pixels->row0) break;
-	    if (testSource->peak->y >= source->pixels->row0 + source->pixels->numRows) {
-		fprintf (stderr, "warning: invalid condition\n");
-		continue;
-	    }
-	    psArrayAdd (overlap, 100, testSource);
-	}
+        // search backwards for overlapping sources
+        for (int j = N - 1; j >= 0; j--) {
+            testSource = sources->data[j];
+            if (testSource->peak->x <  source->pixels->col0) continue;
+            if (testSource->peak->x >= source->pixels->col0 + source->pixels->numCols) continue;
+            if (testSource->peak->y <  source->pixels->row0) break;
+            if (testSource->peak->y >= source->pixels->row0 + source->pixels->numRows) {
+                fprintf (stderr, "warning: invalid condition\n");
+                continue;
+            }
+            psArrayAdd (overlap, 100, testSource);
+        }
 
-	// search forwards for overlapping sources
-	for (int j = N + 1; j < sources->n; j++) {
-	    testSource = sources->data[j];
-	    if (testSource->peak->x <  source->pixels->col0) continue;
-	    if (testSource->peak->x >= source->pixels->col0 + source->pixels->numCols) continue;
-	    if (testSource->peak->y <  source->pixels->row0) {
-		fprintf (stderr, "warning: invalid condition\n");
-		continue;
-	    }
-	    if (testSource->peak->y >= source->pixels->row0 + source->pixels->numRows) break;
-	    psArrayAdd (overlap, 100, testSource);
-	}
+        // search forwards for overlapping sources
+        for (int j = N + 1; j < sources->n; j++) {
+            testSource = sources->data[j];
+            if (testSource->peak->x <  source->pixels->col0) continue;
+            if (testSource->peak->x >= source->pixels->col0 + source->pixels->numCols) continue;
+            if (testSource->peak->y <  source->pixels->row0) {
+                fprintf (stderr, "warning: invalid condition\n");
+                continue;
+            }
+            if (testSource->peak->y >= source->pixels->row0 + source->pixels->numRows) break;
+            psArrayAdd (overlap, 100, testSource);
+        }
 
-	if (overlap->n == 0) {
-	    psFree (overlap);
-	    continue;
-	}
+        if (overlap->n == 0) {
+            psFree (overlap);
+            continue;
+        }
 
-	// this source has overlapping neighbors, check for actual blends
-	// generate source contour (1/4 peak counts)
-	// set the threshold based on user inputs
-	
-	// threshold is fraction of the source peak flux
-	// image is background subtracted; source->moments->Sky should always be 0.0
-	threshold = FRACTION * source->moments->Peak;
-	// threshold is no less than NSIGMA above the local median sigma?
-	threshold = PS_MAX (threshold, NSIGMA*sqrt(source->moments->dSky));
+        // this source has overlapping neighbors, check for actual blends
+        // generate source contour (1/4 peak counts)
+        // set the threshold based on user inputs
 
-	// generate a basic contour (set of x,y coordinates at-or-below flux level)
-	psArray *contour = pmSourceContour (source->pixels, source->peak->x, source->peak->y, threshold);
-	if (contour == NULL) {
-	    psFree (overlap);
-	    continue;
-	}
+        // threshold is fraction of the source peak flux
+        // image is background subtracted; source->moments->Sky should always be 0.0
+        threshold = FRACTION * source->moments->Peak;
+        // threshold is no less than NSIGMA above the local median sigma?
+        threshold = PS_MAX (threshold, NSIGMA*sqrt(source->moments->dSky));
 
-	// the source contour consists of two vectors, xv and yv.  the contour is 
-	// a series of coordinate pairs, (xv[i],yv[i]) & (xv[i+1],yv[i+1]).  both
-	// coordinate pairs have the same yv[] value, with xv[i] corresponding to
-	// the left boundary and xv[i+1] corresponding to the right boundary
+        // generate a basic contour (set of x,y coordinates at-or-below flux level)
+        psArray *contour = pmSourceContour (source->pixels, source->peak->x, source->peak->y, threshold);
+        if (contour == NULL) {
+            psFree (overlap);
+            continue;
+        }
 
-	psVector *xv = contour->data[0];
-	psVector *yv = contour->data[1];
-	for (int k = 0; k < overlap->n; k++) {
-	    testSource = overlap->data[k];
-	    if (testSource->peak->counts > source->peak->counts) continue;
-	    for (int j = 0; j < xv->n; j+=2) {
-		if (fabs(yv->data.F32[j] - testSource->peak->y) > 0.5) continue;
-		if (xv->data.F32[j+0] > testSource->peak->x) break;
-		if (xv->data.F32[j+1] < testSource->peak->x) break;
+        // the source contour consists of two vectors, xv and yv.  the contour is
+        // a series of coordinate pairs, (xv[i],yv[i]) & (xv[i+1],yv[i+1]).  both
+        // coordinate pairs have the same yv[] value, with xv[i] corresponding to
+        // the left boundary and xv[i+1] corresponding to the right boundary
 
-		testSource->mode |= PM_SOURCE_MODE_BLEND;
+        psVector *xv = contour->data[0];
+        psVector *yv = contour->data[1];
+        for (int k = 0; k < overlap->n; k++) {
+            testSource = overlap->data[k];
+            if (testSource->peak->counts > source->peak->counts) continue;
+            for (int j = 0; j < xv->n; j+=2) {
+                if (fabs(yv->data.F32[j] - testSource->peak->y) > 0.5) continue;
+                if (xv->data.F32[j+0] > testSource->peak->x) break;
+                if (xv->data.F32[j+1] < testSource->peak->x) break;
 
-		// add this to the list of source->blends
-		if (source->blends == NULL) {
-		    source->blends = psArrayAlloc (16);
-		    // DROP source->blends->n = 0;
-		}
-		psArrayAdd (source->blends, 16, testSource);
+                testSource->mode |= PM_SOURCE_MODE_BLEND;
 
-		Nblend ++;
-		j = xv->n;
-	    }
-	}  
-	psFree (overlap);
-	psFree (contour);
+                // add this to the list of source->blends
+                if (source->blends == NULL) {
+                    source->blends = psArrayAllocEmpty (16);
+                }
+                psArrayAdd (source->blends, 16, testSource);
+
+                Nblend ++;
+                j = xv->n;
+            }
+        }
+        psFree (overlap);
+        psFree (contour);
     }
     psLogMsg ("psphot.deblend", 3, "identified %d blended objects (%f sec)\n", Nblend, psTimerMark ("psphot"));
Index: /trunk/psphot/src/psphotChoosePSF.c
===================================================================
--- /trunk/psphot/src/psphotChoosePSF.c	(revision 9733)
+++ /trunk/psphot/src/psphotChoosePSF.c	(revision 9734)
@@ -50,5 +50,5 @@
     pmSourceFitModelInit (15, 0.1, POISSON_ERRORS);
 
-    stars = psArrayAlloc (sources->n);
+    stars = psArrayAllocEmpty (sources->n);
 
     // select the candidate PSF stars (pointers to original sources)
@@ -87,5 +87,4 @@
     // set up an array to store the results
     psArray *models = psArrayAlloc (list->n);
-    models->n = list->n;
 
     // try each model option listed in config
Index: /trunk/psphot/src/psphotEnsemblePSF.c
===================================================================
--- /trunk/psphot/src/psphotEnsemblePSF.c	(revision 9733)
+++ /trunk/psphot/src/psphotEnsemblePSF.c	(revision 9734)
@@ -27,6 +27,6 @@
 
     // pre-calculate all model pixels
-    psArray  *models = psArrayAlloc (sources->n);
-    psVector *index  = psVectorAlloc (sources->n, PS_TYPE_U32);
+    psArray  *models = psArrayAllocEmpty (sources->n);
+    psVector *index  = psVectorAllocEmpty (sources->n, PS_TYPE_U32);
     // DROP models->n = index->n = 0;
 
@@ -141,5 +141,4 @@
     psSparse *sparse = psSparseAlloc (models->n, 100);
     psVector *weight = psVectorAlloc (models->n, PS_TYPE_F32);
-    weight->n = models->n;
 
     for (int i = 0; i < models->n; i++) {
Index: /trunk/psphot/src/psphotFakeSources.c
===================================================================
--- /trunk/psphot/src/psphotFakeSources.c	(revision 9733)
+++ /trunk/psphot/src/psphotFakeSources.c	(revision 9734)
@@ -6,25 +6,24 @@
 
     psArray *sources = psArrayAlloc (50);
-    sources->n = 50;
 
     for (int i = 0; i < sources->n; i++) {
-	pmSource *source = pmSourceAlloc ();
-	source->moments = pmMomentsAlloc ();
-	source->moments->x = 10;
-	source->moments->y = 10;
-	source->moments->Sx = 1;
-	source->moments->Sy = 1;
-	source->moments->Sxy = 0;
-	source->moments->Sum = 1000;
-	source->moments->Peak = 100;
-	source->moments->Sky = 10;
-	source->moments->nPixels = 10;
+        pmSource *source = pmSourceAlloc ();
+        source->moments = pmMomentsAlloc ();
+        source->moments->x = 10;
+        source->moments->y = 10;
+        source->moments->Sx = 1;
+        source->moments->Sy = 1;
+        source->moments->Sxy = 0;
+        source->moments->Sum = 1000;
+        source->moments->Peak = 100;
+        source->moments->Sky = 10;
+        source->moments->nPixels = 10;
 
-	source->peak = pmPeakAlloc (10, 10, 0, 0);
-	source->type = PM_SOURCE_TYPE_STAR;
+        source->peak = pmPeakAlloc (10, 10, 0, 0);
+        source->type = PM_SOURCE_TYPE_STAR;
 
-	pmModelType modelType = pmModelSetType ("PS_MODEL_QGAUSS");
-	source->modelPSF = pmSourceModelGuess (source, modelType);
-	sources->data[i] = source;
+        pmModelType modelType = pmModelSetType ("PS_MODEL_QGAUSS");
+        source->modelPSF = pmSourceModelGuess (source, modelType);
+        sources->data[i] = source;
     }
     return sources;
Index: /trunk/psphot/src/psphotFitSet.c
===================================================================
--- /trunk/psphot/src/psphotFitSet.c	(revision 9733)
+++ /trunk/psphot/src/psphotFitSet.c	(revision 9734)
@@ -8,18 +8,17 @@
     if (f == NULL) return false;
 
-    psArray *modelSet = psArrayAlloc (16);
-    // DROP modelSet->n = 0;
+    psArray *modelSet = psArrayAllocEmpty (16);
 
     while (fscanf (f, "%lf %lf %lf", &x, &y, &Io) == 3) {
-	pmModel *model = pmModelAlloc (oneModel->type);
-	
-	for (psS32 i = 0; i < model->params->n; i++) {
-	    model->params->data.F32[i] = oneModel->params->data.F32[i];
-	    model->dparams->data.F32[i] = oneModel->dparams->data.F32[i];
-	}
-	model->params->data.F32[1] = Io;
-	model->params->data.F32[2] = x;
-	model->params->data.F32[3] = y;
-	psArrayAdd (modelSet, 16, model);
+        pmModel *model = pmModelAlloc (oneModel->type);
+
+        for (psS32 i = 0; i < model->params->n; i++) {
+            model->params->data.F32[i] = oneModel->params->data.F32[i];
+            model->dparams->data.F32[i] = oneModel->dparams->data.F32[i];
+        }
+        model->params->data.F32[1] = Io;
+        model->params->data.F32[2] = x;
+        model->params->data.F32[3] = y;
+        psArrayAdd (modelSet, 16, model);
     }
 
@@ -31,14 +30,14 @@
     // subtract object, leave local sky
     for (int i = 0; i < modelSet->n; i++) {
-	pmModel *model = modelSet->data[i];
-	pmModelSub (source->pixels, source->mask, model, false, false);
-    
-	fprintf (stderr, "output parameters (obj %d):\n", i);
-	for (int n = 0; n < model->params->n; n++) {
-	    fprintf (stderr, "%d : %f\n", n, model->params->data.F32[n]);
-	}
+        pmModel *model = modelSet->data[i];
+        pmModelSub (source->pixels, source->mask, model, false, false);
+
+        fprintf (stderr, "output parameters (obj %d):\n", i);
+        for (int n = 0; n < model->params->n; n++) {
+            fprintf (stderr, "%d : %f\n", n, model->params->data.F32[n]);
+        }
     }
 
-    // write out 
+    // write out
     psphotSaveImage (NULL, source->pixels, "resid.fits");
     psphotSaveImage (NULL, source->mask, "mask.fits");
Index: /trunk/psphot/src/psphotSourceFits.c
===================================================================
--- /trunk/psphot/src/psphotSourceFits.c	(revision 9733)
+++ /trunk/psphot/src/psphotSourceFits.c	(revision 9734)
@@ -21,8 +21,8 @@
     y = PSF->params->data.F32[3];
 
-    psArray *modelSet = psArrayAlloc (source->blends->n + 1);
+    psArray *modelSet = psArrayAllocEmpty (source->blends->n + 1);
     psArrayAdd (modelSet, 16, PSF);
 
-    psArray *sourceSet = psArrayAlloc (source->blends->n + 1);
+    psArray *sourceSet = psArrayAllocEmpty (source->blends->n + 1);
     psArrayAdd (sourceSet, 16, source);
 
@@ -300,5 +300,4 @@
 
     modelSet = psArrayAlloc (2);
-    modelSet->n = 2;
 
     DBL = pmModelCopy (PSF);
Index: /trunk/psphot/src/psphotSourceStats.c
===================================================================
--- /trunk/psphot/src/psphotSourceStats.c	(revision 9733)
+++ /trunk/psphot/src/psphotSourceStats.c	(revision 9734)
@@ -2,5 +2,5 @@
 
 // 2006.02.02 : no leaks
-psArray *psphotSourceStats (pmReadout *readout, psMetadata *recipe, psArray *peaks) 
+psArray *psphotSourceStats (pmReadout *readout, psMetadata *recipe, psArray *peaks)
 {
     bool     status  = false;
@@ -16,63 +16,63 @@
     char *breakPt  = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
 
-    sources = psArrayAlloc (peaks->n);
+    sources = psArrayAllocEmpty (peaks->n);
 
     for (int i = 0; i < peaks->n; i++) {
 
-	// create a new source, add peak
-	pmSource *source = pmSourceAlloc();
-	source->peak = (pmPeak *)psMemIncrRefCounter(peaks->data[i]);
+        // create a new source, add peak
+        pmSource *source = pmSourceAlloc();
+        source->peak = (pmPeak *)psMemIncrRefCounter(peaks->data[i]);
 
-	// allocate image, weight, mask arrays for each peak (square of radius OUTER)
-	pmSourceDefinePixels (source, readout, source->peak->x, source->peak->y, OUTER);
-	if (!strcasecmp (breakPt, "PEAKS")) { 
-	    psArrayAdd (sources, 100, source);
-	    psFree (source);
-	    continue;
-	}
+        // allocate image, weight, mask arrays for each peak (square of radius OUTER)
+        pmSourceDefinePixels (source, readout, source->peak->x, source->peak->y, OUTER);
+        if (!strcasecmp (breakPt, "PEAKS")) {
+            psArrayAdd (sources, 100, source);
+            psFree (source);
+            continue;
+        }
 
-	// XXX skip faint sources?
+        // XXX skip faint sources?
 
-	// measure a local sky value
-	// XXX EAM : this should use ROBUST not SAMPLE median, but it is broken
-	status = pmSourceLocalSky (source, PS_STAT_SAMPLE_MEDIAN, INNER);
-	if (!status) {
-	  psFree (source);
-	  continue;
-	}
+        // measure a local sky value
+        // XXX EAM : this should use ROBUST not SAMPLE median, but it is broken
+        status = pmSourceLocalSky (source, PS_STAT_SAMPLE_MEDIAN, INNER);
+        if (!status) {
+          psFree (source);
+          continue;
+        }
 
-	// measure the local sky variance (needed if noise is not sqrt(signal))
-	// XXX EAM : this should use ROBUST not SAMPLE median, but it is broken
-	status = pmSourceLocalSkyVariance (source, PS_STAT_SAMPLE_MEDIAN, INNER);
-	if (!status) {
-	  psFree (source);
-	  continue;
-	}
+        // measure the local sky variance (needed if noise is not sqrt(signal))
+        // XXX EAM : this should use ROBUST not SAMPLE median, but it is broken
+        status = pmSourceLocalSkyVariance (source, PS_STAT_SAMPLE_MEDIAN, INNER);
+        if (!status) {
+          psFree (source);
+          continue;
+        }
 
-	// measure basic source moments
-	status = pmSourceMoments (source, RADIUS);
-	if (status) {
-	    // add to the source array
-	    psArrayAdd (sources, 100, source);
-	    psFree (source);
-	    continue;
-	}
+        // measure basic source moments
+        status = pmSourceMoments (source, RADIUS);
+        if (status) {
+            // add to the source array
+            psArrayAdd (sources, 100, source);
+            psFree (source);
+            continue;
+        }
 
-	// if no valid pixels, or massive swing, likely saturated source,
-	// try a much larger box
-	BIG_RADIUS = PS_MIN (INNER, 3*RADIUS);
-	psTrace ("psphot", 4, "retrying moments for %d, %d\n", source->peak->x, source->peak->y);
-	status = pmSourceMoments (source, BIG_RADIUS);
-	if (status) {
-	    // add to the source array
-	    psArrayAdd (sources, 100, source);
-	    psFree (source);
-	    continue;
-	}
+        // if no valid pixels, or massive swing, likely saturated source,
+        // try a much larger box
+        BIG_RADIUS = PS_MIN (INNER, 3*RADIUS);
+        psTrace ("psphot", 4, "retrying moments for %d, %d\n", source->peak->x, source->peak->y);
+        status = pmSourceMoments (source, BIG_RADIUS);
+        if (status) {
+            // add to the source array
+            psArrayAdd (sources, 100, source);
+            psFree (source);
+            continue;
+        }
 
-	psFree (source);
-	continue;
+        psFree (source);
+        continue;
     }
-     
+
     psLogMsg ("psphot", 3, "%ld moments: %f sec\n", sources->n, psTimerMark ("psphot"));
 
Index: /trunk/psphot/src/psphotTest.c
===================================================================
--- /trunk/psphot/src/psphotTest.c	(revision 9733)
+++ /trunk/psphot/src/psphotTest.c	(revision 9734)
@@ -15,5 +15,5 @@
 int main (int argc, char **argv) {
 
-    psRegion region = {0,0,0,0};	// a region representing the entire array
+    psRegion region = {0,0,0,0};        // a region representing the entire array
     psphotTestArguments (&argc, argv);
 
@@ -56,22 +56,21 @@
     mdi->data.S32 = 0;
     mdi->type = PS_DATA_S32;
-    
+
     // create a test image
     // psImage *tmpimage = psImageAlloc (10, 10, PS_DATA_F32);
 
     // create a test table
-    table = psArrayAlloc (10);
-table->n = 0;
+    table = psArrayAllocEmpty (10);
 
     for (int i = 0; i < 10; i++) {
-	row = psMetadataAlloc ();
-	psMetadataAdd (row, PS_LIST_TAIL, "ROW",   PS_DATA_S32,    "", i);
-	psMetadataAdd (row, PS_LIST_TAIL, "FROW",  PS_TYPE_F32,    "", 0.1*i);
-	psMetadataAdd (row, PS_LIST_TAIL, "DUMMY", PS_DATA_STRING, "", "test line");
-    
-	table->data[i] = row;
+        row = psMetadataAlloc ();
+        psMetadataAdd (row, PS_LIST_TAIL, "ROW",   PS_DATA_S32,    "", i);
+        psMetadataAdd (row, PS_LIST_TAIL, "FROW",  PS_TYPE_F32,    "", 0.1*i);
+        psMetadataAdd (row, PS_LIST_TAIL, "DUMMY", PS_DATA_STRING, "", "test line");
+
+        table->data[i] = row;
     }
     table->n = 10;
-    
+
     psMetadata *theader = psMetadataAlloc ();
     psMetadataAdd (theader, PS_LIST_HEAD, "EXTNAME", PS_DATA_STRING, "extension name", "SMPFILE");