Index: branches/pap/psModules/src/objects/pmSourcePhotometry.c
===================================================================
--- branches/pap/psModules/src/objects/pmSourcePhotometry.c	(revision 28179)
+++ branches/pap/psModules/src/objects/pmSourcePhotometry.c	(revision 28484)
@@ -107,12 +107,12 @@
         // the source peak pixel is guaranteed to be on the image, and only minimally different from the source center
         double fluxScale = pmTrend2DEval (psf->FluxScale, (float)source->peak->x, (float)source->peak->y);
-	psAssert (isfinite(fluxScale), "how can the flux scale be invalid? source at %d, %d\n", source->peak->x, source->peak->y);
-	psAssert (fluxScale > 0.0, "how can the flux scale be negative? source at %d, %d\n", source->peak->x, source->peak->y);
-	source->psfFlux = fluxScale * source->modelPSF->params->data.F32[PM_PAR_I0];
-	source->psfFluxErr = fluxScale * source->modelPSF->dparams->data.F32[PM_PAR_I0];
-	source->psfMag = -2.5*log10(source->psfFlux);
+        psAssert (isfinite(fluxScale), "how can the flux scale be invalid? source at %d, %d\n", source->peak->x, source->peak->y);
+        psAssert (fluxScale > 0.0, "how can the flux scale be negative? source at %d, %d\n", source->peak->x, source->peak->y);
+        source->psfFlux = fluxScale * source->modelPSF->params->data.F32[PM_PAR_I0];
+        source->psfFluxErr = fluxScale * source->modelPSF->dparams->data.F32[PM_PAR_I0];
+        source->psfMag = -2.5*log10(source->psfFlux);
     } else {
         status = pmSourcePhotometryModel (&source->psfMag, &source->psfFlux, source->modelPSF);
-	source->psfFluxErr = source->psfFlux * (source->modelPSF->dparams->data.F32[PM_PAR_I0] / source->modelPSF->params->data.F32[PM_PAR_I0]);
+        source->psfFluxErr = source->psfFlux * (source->modelPSF->dparams->data.F32[PM_PAR_I0] / source->modelPSF->params->data.F32[PM_PAR_I0]);
     }
 
@@ -120,18 +120,18 @@
     if (source->modelFits) {
         bool foundEXT = false;
-	for (int i = 0; i < source->modelFits->n; i++) {
-	    pmModel *model = source->modelFits->data[i];
-	    status = pmSourcePhotometryModel (&model->mag, NULL, model);
-	    if (model == source->modelEXT) foundEXT = true;
-	}
-	if (foundEXT) {
-	    source->extMag = source->modelEXT->mag;
-	} else {
-	    status = pmSourcePhotometryModel (&source->extMag, NULL, source->modelEXT);
-	}
+        for (int i = 0; i < source->modelFits->n; i++) {
+            pmModel *model = source->modelFits->data[i];
+            status = pmSourcePhotometryModel (&model->mag, NULL, model);
+            if (model == source->modelEXT) foundEXT = true;
+        }
+        if (foundEXT) {
+            source->extMag = source->modelEXT->mag;
+        } else {
+            status = pmSourcePhotometryModel (&source->extMag, NULL, source->modelEXT);
+        }
     } else {
-	if (source->modelEXT) {
-	    status = pmSourcePhotometryModel (&source->extMag, NULL, source->modelEXT);
-	}
+        if (source->modelEXT) {
+            status = pmSourcePhotometryModel (&source->extMag, NULL, source->modelEXT);
+        }
     }
 
@@ -204,5 +204,5 @@
         }
         if (mode & PM_SOURCE_PHOT_APCORR) {
-	    // XXX this should be removed -- we no longer fit for the 'sky bias'
+            // XXX this should be removed -- we no longer fit for the 'sky bias'
             rflux   = pow (10.0, 0.4*source->psfMag);
             source->apMag -= PS_SQR(source->apRadius)*rflux * psf->skyBias + psf->skySat / rflux;
@@ -236,11 +236,11 @@
     flux = model->modelFlux (model->params);
     if (flux > 0) {
-	mag = -2.5*log10(flux);
+        mag = -2.5*log10(flux);
     }
     if (fitMag) {
-	*fitMag = mag;
+        *fitMag = mag;
     }
     if (fitFlux) {
-	*fitFlux = flux;
+        *fitFlux = flux;
     }
 
@@ -380,5 +380,5 @@
 
     if (source->diffStats == NULL) {
-	source->diffStats = pmSourceDiffStatsAlloc();
+        source->diffStats = pmSourceDiffStatsAlloc();
     }
 
@@ -388,5 +388,5 @@
     int   nMask = 0;
     int   nBad  = 0;
-    
+
     psImage *flux     = source->pixels;
     psImage *variance = source->variance;
@@ -394,28 +394,28 @@
 
     for (int iy = 0; iy < flux->numRows; iy++) {
-	for (int ix = 0; ix < flux->numCols; ix++) {
+        for (int ix = 0; ix < flux->numCols; ix++) {
             if (mask->data.PS_TYPE_IMAGE_MASK_DATA[iy][ix] & maskVal) {
-		nMask ++;
-                continue;
-	    }
-
-	    float SN = flux->data.F32[iy][ix] / sqrt(variance->data.F32[iy][ix]);
-
-	    if (SN > +FLUX_LIMIT) { 
-		nGood ++;
-		fGood += fabs(flux->data.F32[iy][ix]);
-	    }
-
-	    if (SN < -FLUX_LIMIT) { 
-		nBad ++;
-		fBad += fabs(flux->data.F32[iy][ix]);
-	    }
-	}
+                nMask ++;
+                continue;
+            }
+
+            float SN = flux->data.F32[iy][ix] / sqrt(variance->data.F32[iy][ix]);
+
+            if (SN > +FLUX_LIMIT) {
+                nGood ++;
+                fGood += fabs(flux->data.F32[iy][ix]);
+            }
+
+            if (SN < -FLUX_LIMIT) {
+                nBad ++;
+                fBad += fabs(flux->data.F32[iy][ix]);
+            }
+        }
     }
 
     source->diffStats->nGood      = nGood;
-    source->diffStats->fRatio     = (fGood + fBad         == 0.0) ? NAN : fGood / (fGood + fBad);	   
-    source->diffStats->nRatioBad  = (nGood + nBad         == 0)   ? NAN : nGood / (float) (nGood + nBad);	   
-    source->diffStats->nRatioMask = (nGood + nMask        == 0)   ? NAN : nGood / (float) (nGood + nMask);	   
+    source->diffStats->fRatio     = (fGood + fBad         == 0.0) ? NAN : fGood / (fGood + fBad);
+    source->diffStats->nRatioBad  = (nGood + nBad         == 0)   ? NAN : nGood / (float) (nGood + nBad);
+    source->diffStats->nRatioMask = (nGood + nMask        == 0)   ? NAN : nGood / (float) (nGood + nMask);
     source->diffStats->nRatioAll  = (nGood + nMask + nBad == 0)   ? NAN : nGood / (float) (nGood + nMask + nBad);
 
@@ -628,5 +628,5 @@
         }
     }
-    model->nPix = Npix; 
+    model->nPix = Npix;
     model->nDOF = Npix - 1;
     model->chisq = dC;
@@ -636,5 +636,5 @@
 
 
-double pmSourceModelWeight(const pmSource *Mi, int term, const bool unweighted_sum, const float covarFactor)
+double pmSourceModelWeight(const pmSource *Mi, int term, const bool unweighted_sum, const float covarFactor, psImageMaskType maskVal)
 {
     PS_ASSERT_PTR_NON_NULL(Mi, NAN);
@@ -652,5 +652,5 @@
     for (int yi = 0; yi < Pi->numRows; yi++) {
         for (int xi = 0; xi < Pi->numCols; xi++) {
-            if (Ti->data.PS_TYPE_IMAGE_MASK_DATA[yi][xi])
+            if (Ti->data.PS_TYPE_IMAGE_MASK_DATA[yi][xi] & maskVal)
                 continue;
             if (!unweighted_sum) {
@@ -684,5 +684,5 @@
 }
 
-double pmSourceModelDotModel (const pmSource *Mi, const pmSource *Mj, const bool unweighted_sum, const float covarFactor)
+double pmSourceModelDotModel (const pmSource *Mi, const pmSource *Mj, const bool unweighted_sum, const float covarFactor, psImageMaskType maskVal)
 {
     PS_ASSERT_PTR_NON_NULL(Mi, NAN);
@@ -727,7 +727,7 @@
     for (yi = yIs, yj = yJs; yi < yIe; yi++, yj++) {
         for (xi = xIs, xj = xJs; xi < xIe; xi++, xj++) {
-            if (Ti->data.PS_TYPE_IMAGE_MASK_DATA[yi][xi])
-                continue;
-            if (Tj->data.PS_TYPE_IMAGE_MASK_DATA[yj][xj])
+            if (Ti->data.PS_TYPE_IMAGE_MASK_DATA[yi][xi] & maskVal)
+                continue;
+            if (Tj->data.PS_TYPE_IMAGE_MASK_DATA[yj][xj] & maskVal)
                 continue;
 
@@ -746,5 +746,5 @@
 }
 
-double pmSourceDataDotModel (const pmSource *Mi, const pmSource *Mj, const bool unweighted_sum, const float covarFactor)
+double pmSourceDataDotModel (const pmSource *Mi, const pmSource *Mj, const bool unweighted_sum, const float covarFactor, psImageMaskType maskVal)
 {
     PS_ASSERT_PTR_NON_NULL(Mi, NAN);
@@ -789,7 +789,7 @@
     for (yi = yIs, yj = yJs; yi < yIe; yi++, yj++) {
         for (xi = xIs, xj = xJs; xi < xIe; xi++, xj++) {
-            if (Ti->data.PS_TYPE_IMAGE_MASK_DATA[yi][xi])
-                continue;
-            if (Tj->data.PS_TYPE_IMAGE_MASK_DATA[yj][xj])
+            if (Ti->data.PS_TYPE_IMAGE_MASK_DATA[yi][xi] & maskVal)
+                continue;
+            if (Tj->data.PS_TYPE_IMAGE_MASK_DATA[yj][xj] & maskVal)
                 continue;