Index: branches/eam_branches/ipp-20130419/psModules/src/astrom/pmAstrometryUtils.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/astrom/pmAstrometryUtils.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/astrom/pmAstrometryUtils.c	(revision 35421)
@@ -23,4 +23,7 @@
 #include "pmAstrometryUtils.h"
 
+// this is used by the test output block
+static int Nout = 0;
+
 // given a 2D transformation -- L(x,y),M(x,y) -- find the coordinates x,y
 // for which L,M = 0,0. tol is the allowed error on x,y.
@@ -33,4 +36,20 @@
     double Xo = (trans->y->coeff[0][0]*trans->x->coeff[0][1] - trans->x->coeff[0][0]*trans->y->coeff[0][1])/R;
     double Yo = (trans->x->coeff[0][0]*trans->y->coeff[1][0] - trans->y->coeff[0][0]*trans->x->coeff[1][0])/R;
+
+    if (0) {
+      // this is a test output block, not used in normal operations
+      char filename[64];
+      snprintf (filename, 64, "trans.%03d.md", Nout);
+      FILE *f = fopen (filename, "w");
+
+      psMetadata *md = psMetadataAlloc();
+
+      psPolynomial2DtoMetadata (md, trans->x, "X");
+      psPolynomial2DtoMetadata (md, trans->y, "Y");
+      psMetadataConfigPrint (f, md);
+      psFree (md);
+      fclose (f);
+      Nout ++;
+    }
 
     // iterate to actual solution: requires small non-linear terms
@@ -52,11 +71,13 @@
         for (int i = 0; (dPos > tol) && (i < 20); i++) {
             // NOTE: order for Alpha is: [y][x]
+            // original Alpha->data.F32[1][0] = psPolynomial2DEval (XdY, Xo, Yo);
+            // original Alpha->data.F32[0][1] = psPolynomial2DEval (YdX, Xo, Yo);
             Alpha->data.F32[0][0] = psPolynomial2DEval (XdX, Xo, Yo);
-            Alpha->data.F32[1][0] = psPolynomial2DEval (XdY, Xo, Yo);
-            Alpha->data.F32[0][1] = psPolynomial2DEval (YdX, Xo, Yo);
+            Alpha->data.F32[0][1] = psPolynomial2DEval (XdY, Xo, Yo);
+            Alpha->data.F32[1][0] = psPolynomial2DEval (YdX, Xo, Yo);
             Alpha->data.F32[1][1] = psPolynomial2DEval (YdY, Xo, Yo);
 
-            Beta->data.F32[0] = psPolynomial2DEval (trans->x, Xo, Yo);
-            Beta->data.F32[1] = psPolynomial2DEval (trans->y, Xo, Yo);
+            Beta->data.F32[0] = -1.0 * psPolynomial2DEval (trans->x, Xo, Yo);
+            Beta->data.F32[1] = -1.0 * psPolynomial2DEval (trans->y, Xo, Yo);
 
             if (!psMatrixGJSolve (Alpha, Beta)) {
@@ -71,8 +92,8 @@
 	    }
 
-            Xo -= Beta->data.F32[0];
-            Yo -= Beta->data.F32[1];
+            Xo += Beta->data.F32[0];
+            Yo += Beta->data.F32[1];
             dPos = hypot(Beta->data.F32[0], Beta->data.F32[1]);
-
+	    psTrace ("psastro", 5, "i: %d, Xo,Yo : %f %f  dX, dY: %f %f dPos: %f\n", i, Xo, Yo, Beta->data.F32[0], Beta->data.F32[1], dPos);
         }
         psFree (Alpha);
Index: branches/eam_branches/ipp-20130419/psModules/src/camera/pmFPAMaskWeight.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/camera/pmFPAMaskWeight.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/camera/pmFPAMaskWeight.c	(revision 35421)
@@ -301,6 +301,8 @@
     psRegion *trimsec = psMetadataLookupPtr(&mdok, cell->concepts, "CELL.TRIMSEC"); // Trim section
     if (!mdok || psRegionIsNaN(*trimsec)) {
-        psError(PS_ERR_IO, true, "CELL.TRIMSEC is not set --- unable to set variance.\n");
-        return false;
+      // if trimsec is not defined, use the full image
+      trimsec = psRegionAlloc(0,0,0,0);
+      // psError(PS_ERR_IO, true, "CELL.TRIMSEC is not set --- unable to set variance.\n");
+      // return false;
     }
 
Index: branches/eam_branches/ipp-20130419/psModules/src/camera/pmFPAfileDefine.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/camera/pmFPAfileDefine.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/camera/pmFPAfileDefine.c	(revision 35421)
@@ -446,7 +446,131 @@
 }
 
+// given a filename, convert to UNIX namespace and read the PHU 
+psMetadata *readPHUfromFilename (char *filename, pmConfig *config) {
+
+    // Need to generate an FPA
+    psString realName = pmConfigConvertFilename(filename, config, false, false);
+    if (!realName) {
+	psError(psErrorCodeLast(), false, "Failed to convert file name %s", filename);
+	return NULL;
+    }
+
+    // load the header of the first image
+    // EXTWORD (fits->extword) is not relevant to the PHU
+    psFits *fits = psFitsOpen(realName, "r"); // FITS file
+    if (!fits) {
+	psError(psErrorCodeLast(), false, "Failed to open file %s", realName);
+	psFree(realName);
+	return NULL;
+    }
+
+    psMetadata *phu = psFitsReadHeader (NULL, fits); // Primary header
+    if (!phu) {
+	psError(psErrorCodeLast(), false, "Failed to read file header %s", realName);
+	psFree(realName);
+	return NULL;
+    }
+
+    if (!psFitsClose(fits)) {
+	psError(psErrorCodeLast(), false, "Failed to close file %s", realName);
+	psFree(realName);
+	psFree(phu);
+	return NULL;
+    }
+
+    return phu;
+}
+
+// this this function wants to return:
+// pmFPA, PHU, fileLevel, outConfig
+// camera, cameraName, formatName
+typedef struct {
+    pmFPA *fpa;
+    psMetadata *phu;
+    psMetadata *format;
+    pmFPALevel fileLevel;
+    psString cameraName;
+    psString formatName;
+} pmFPAfromFilenameOutput;
+
+// for the given filename, read PHU and determine camera format; build an FPA for the file
+bool pmFPAfromFilename (pmFPAfromFilenameOutput *output, pmConfig **outConfig, pmConfig *sysConfig, char *filename){
+
+    // Need to generate an FPA
+    psMetadata *phu = readPHUfromFilename (filename, sysConfig);
+    if (!phu) {
+	psError(psErrorCodeLast(), false, "Failed to read PHU for %s", filename);
+	return false;
+    }
+
+    // if we expect the loaded FPA to differ in configuration from the current system configuration
+    // generate an output config for this FPA
+    pmConfig *config = NULL;
+    if (outConfig) {
+	config = pmConfigAlloc();
+	config->user = psMemIncrRefCounter(sysConfig->user);
+	config->system = psMemIncrRefCounter(sysConfig->system);
+
+	psFree (config->files);
+	config->files = psMemIncrRefCounter(sysConfig->files);
+	psFree (config->arguments);
+	config->arguments = psMemIncrRefCounter(sysConfig->arguments);
+
+	*outConfig = config;
+    } else {
+	config = sysConfig;
+    }
+
+    // values which are returned to calling function
+    psString formatName = NULL;	// Name of camera format
+    psString cameraName = NULL;	// Name of camera
+    psMetadata *camera = NULL;	// Camera configuration
+
+    // Determine the current format from the header; determine camera if not specified already.
+    psMetadata *format = pmConfigCameraFormatFromHeader(&camera, &cameraName, &formatName, config, phu, true);
+    if (!format) {
+	psError(psErrorCodeLast(), false, "Failed to determine camera format for %s", filename);
+	psFree(camera);
+	psFree(formatName);
+	psFree(phu);
+	return false;
+    }
+
+    pmFPALevel fileLevel = pmFPAPHULevel(format);
+    if (fileLevel == PM_FPA_LEVEL_NONE) {
+	psError(PM_ERR_CONFIG, true, "Unable to determine file level for %s", filename);
+	psFree(camera);
+	psFree(formatName);
+	psFree(phu);
+	return false;
+    }
+
+    // build the template fpa, set up the basic view
+    // we supply the metaCamera name (if NULL, baseCamera name is used)
+    pmFPA *fpa = pmFPAConstruct(camera, cameraName);
+    psFree(camera);
+
+    if (!fpa) {
+	psError(psErrorCodeLast(), false, "Failed to construct FPA from %s", filename);
+	psFree(formatName);
+	psFree(format);
+	psFree(phu);
+	return NULL;
+    }
+
+    output->fpa = fpa;
+    output->phu = phu;
+    output->format = format;
+    output->fileLevel = fileLevel;
+    output->cameraName = cameraName;
+    output->formatName = formatName;
+
+    return true;
+
+}
 
 /// Define a file from an array of filenames
-static pmFPAfile *fpaFileDefineFromArray(pmConfig *config, // Configuration
+static pmFPAfile *fpaFileDefineFromArray(pmConfig **outConfig, // output configuration
+					 pmConfig *sysConfig, // global configuration
                                          pmFPAfile *bind, // File to bind to, or NULL
                                          const char *name, // Name of file
@@ -454,5 +578,5 @@
     )
 {
-    PS_ASSERT_PTR_NON_NULL(config, NULL);
+    PS_ASSERT_PTR_NON_NULL(sysConfig, NULL);
     PS_ASSERT_STRING_NON_EMPTY(name, NULL);
 
@@ -463,4 +587,5 @@
     pmFPALevel fileLevel = PM_FPA_LEVEL_NONE; // Level for files
     psMetadata *phu = NULL;             // Primary header
+
     if (bind) {
         // Use the FPA we're binding to
@@ -468,67 +593,17 @@
         fileLevel = bind->fileLevel;
     } else {
-        // Need to generate an FPA
-        psString realName = pmConfigConvertFilename(filenames->data[0], config, false, false);
-        if (!realName) {
-            psError(psErrorCodeLast(), false, "Failed to convert file name %s", (char *)filenames->data[0]);
-            return NULL;
-        }
-
-        // load the header of the first image
-        // EXTWORD (fits->extword) is not relevant to the PHU
-        psFits *fits = psFitsOpen(realName, "r"); // FITS file
-        if (!fits) {
-            psError(psErrorCodeLast(), false, "Failed to open file %s", realName);
-            psFree(realName);
-            return NULL;
-        }
-        phu = psFitsReadHeader (NULL, fits); // Primary header
-        if (!phu) {
-            psError(psErrorCodeLast(), false, "Failed to read file header %s", realName);
-            psFree(realName);
-            return NULL;
-        }
-        if (!psFitsClose(fits)) {
-            psError(psErrorCodeLast(), false, "Failed to close file %s", realName);
-            psFree(realName);
-            return NULL;
-        }
-
-        // Determine the current format from the header; determine camera if not specified already.
-        psMetadata *camera = NULL;      // Camera configuration
-        format = pmConfigCameraFormatFromHeader(&camera, &cameraName, &formatName, config, phu, true);
-        if (!format) {
-            psError(psErrorCodeLast(), false, "Failed to determine camera format for %s", realName);
-            psFree(camera);
-            psFree(formatName);
-            psFree(realName);
-            psFree(phu);
-            return NULL;
-        }
-
-        fileLevel = pmFPAPHULevel(format);
-        if (fileLevel == PM_FPA_LEVEL_NONE) {
-            psError(PM_ERR_CONFIG, true, "Unable to determine file level for %s", realName);
-            psFree(camera);
-            psFree(formatName);
-            psFree(realName);
-            psFree(phu);
-            return NULL;
-        }
-
-        // build the template fpa, set up the basic view
-        // XXX do we want this to be the baseCamera name or the metaCamera name?
-        fpa = pmFPAConstruct(camera, cameraName);
-        psFree(camera);
-        if (!fpa) {
-            psError(psErrorCodeLast(), false, "Failed to construct FPA from %s", realName);
-            psFree(formatName);
-            psFree(realName);
-            psFree(format);
-            psFree(phu);
-            return NULL;
-        }
-        psFree(realName);
-    }
+	pmFPAfromFilenameOutput output;
+	if (!pmFPAfromFilename (&output, outConfig, sysConfig, filenames->data[0])) {
+	    return NULL;
+	}
+	fpa = output.fpa;
+	phu = output.phu;
+	format = output.format;
+	fileLevel = output.fileLevel;
+	cameraName = output.cameraName;
+	formatName = output.formatName;
+    }
+
+    pmConfig *config = outConfig ? *outConfig : sysConfig;
 
     // load the given filerule (from config->camera) and bind it to the fpa
@@ -561,27 +636,9 @@
         // Check that the file corresponds to the same camera and format
         if (!phu) {
-            psString realName = pmConfigConvertFilename(filenames->data[i], config, false, false);
-            if (!realName) {
-                psError(psErrorCodeLast(), false, "Failed to convert file name %s", (char*)filenames->data[i]);
+	    phu = readPHUfromFilename (filenames->data[i], config);
+            if (!phu) {
+                psError(psErrorCodeLast(), false, "Failed to read PHU for %s", (char *)filenames->data[i]);
                 return NULL;
             }
-            psFits *fits = psFitsOpen(realName, "r"); // FITS file
-            if (!fits) {
-                psError(psErrorCodeLast(), false, "Failed to open file %s", realName);
-                psFree(realName);
-                return NULL;
-            }
-            phu = psFitsReadHeader(NULL, fits);
-            if (!psFitsClose(fits)) {
-                psError(psErrorCodeLast(), false, "Failed to close file %s", realName);
-                psFree(realName);
-                return NULL;
-            }
-            if (!phu) {
-                psError(psErrorCodeLast(), false, "Failed to read file header %s", realName);
-                psFree(realName);
-                return NULL;
-            }
-            psFree(realName);
         }
 
@@ -651,7 +708,6 @@
 }
 
-
-pmFPAfile *pmFPAfileDefineFromArgs(bool *success, pmConfig *config,
-                                   const char *filename, const char *argname)
+// find the file associated with the argname & generate a pmFPAfile for it based on the filerule
+pmFPAfile *pmFPAfileDefineFromArgs(bool *success, pmConfig *config, const char *filename, const char *argname)
 {
     PS_ASSERT_PTR_NON_NULL(config, NULL);
@@ -676,5 +732,5 @@
     }
 
-    pmFPAfile *file = fpaFileDefineFromArray(config, NULL, filename, filenames); // File of interest
+    pmFPAfile *file = fpaFileDefineFromArray(NULL, config, NULL, filename, filenames); // File of interest
 
     if (success) {
@@ -685,6 +741,6 @@
 }
 
-pmFPAfile *pmFPAfileBindFromArgs(bool *success, pmFPAfile *input, pmConfig *config,
-                                 const char *filename, const char *argname)
+// find the file associated with the argname & bind it to the given pmFPAfile for it based on the filerule
+pmFPAfile *pmFPAfileBindFromArgs(bool *success, pmFPAfile *input, pmConfig *config, const char *filename, const char *argname)
 {
     PS_ASSERT_PTR_NON_NULL(input, NULL);
@@ -710,5 +766,5 @@
     }
 
-    pmFPAfile *file = fpaFileDefineFromArray(config, input, filename, filenames); // File of interest
+    pmFPAfile *file = fpaFileDefineFromArray(NULL, config, input, filename, filenames); // File of interest
 
     if (success) {
@@ -719,6 +775,6 @@
 }
 
-pmFPAfile *pmFPAfileDefineSingleFromArgs(bool *success, pmConfig *config, const char *filename,
-                                         const char *argname, int entry)
+// find the specific file associated with the argname & generate a pmFPAfile for it based on the filerule
+pmFPAfile *pmFPAfileDefineSingleFromArgs(bool *success, pmConfig *config, const char *filename, const char *argname, int entry)
 {
     PS_ASSERT_PTR_NON_NULL(config, NULL);
@@ -746,5 +802,5 @@
     psArray *single = psArrayAlloc(1);  // Array of single filename of interest
     single->data[0] = psMemIncrRefCounter(filenames->data[entry]);
-    pmFPAfile *file = fpaFileDefineFromArray(config, NULL, filename, single); // File of interest
+    pmFPAfile *file = fpaFileDefineFromArray(NULL, config, NULL, filename, single); // File of interest
     psFree(single);
 
@@ -756,4 +812,5 @@
 }
 
+// find the file in the config list & generate a pmFPAfile for it based on the filerule
 pmFPAfile *pmFPAfileDefineFromRun(bool *success, pmFPAfile *bind, pmConfig *config, const char *filename)
 {
@@ -769,5 +826,5 @@
     }
 
-    pmFPAfile *file = fpaFileDefineFromArray(config, bind, filename, filenames); // File of interest
+    pmFPAfile *file = fpaFileDefineFromArray(NULL, config, bind, filename, filenames); // File of interest
     psFree(filenames);
 
@@ -779,4 +836,5 @@
 }
 
+// find the files in the config list & generate an array of pmFPAfiles for them based on the filerule
 psArray *pmFPAfileDefineMultipleFromRun(bool *success, psArray *bind, pmConfig *config, const char *filename)
 {
@@ -808,5 +866,5 @@
         dummy->data[0] = files->data[i];
         pmFPAfile *bindFile = bind ? bind->data[i] : NULL; // File to which to bind
-        files->data[i] = psMemIncrRefCounter(fpaFileDefineFromArray(config, bindFile, filename, dummy));
+        files->data[i] = psMemIncrRefCounter(fpaFileDefineFromArray(NULL, config, bindFile, filename, dummy));
         if (!files->data[i]) {
             psError(psErrorCodeLast(), false, "Unable to define file %s %d", filename, i);
@@ -823,4 +881,38 @@
 
     return files;
+}
+
+// find the file associated with the argname & generate a pmFPAfile for it based on the filerule
+pmFPAfile *pmFPAfileDefineNewConfig(bool *success, pmConfig **outConfig, pmConfig *sysConfig, const char *filename, const char *argname)
+{
+    PS_ASSERT_PTR_NON_NULL(outConfig, NULL);
+    PS_ASSERT_PTR_NON_NULL(sysConfig, NULL);
+    PS_ASSERT_STRING_NON_EMPTY(filename, NULL);
+    PS_ASSERT_STRING_NON_EMPTY(argname, NULL);
+
+    // Search the argument data for the named fileset (argname)
+    bool status;                        // Status of MD lookup
+    psArray *filenames = psMetadataLookupPtr(&status, sysConfig->arguments, argname); // Filenames for file
+    if (!status) {
+        if (success) {
+            *success = true;
+        }
+        return NULL;
+    }
+    if (filenames->n == 0) {
+        psError(PM_ERR_CONFIG, true, "No files in array in %s in arguments", argname);
+        if (success) {
+            *success = false;
+        }
+        return NULL;
+    }
+
+    pmFPAfile *file = fpaFileDefineFromArray(outConfig, sysConfig, NULL, filename, filenames); // File of interest
+
+    if (success) {
+        *success = file ? true : false;
+    }
+
+    return file;
 }
 
Index: branches/eam_branches/ipp-20130419/psModules/src/camera/pmFPAfileDefine.h
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/camera/pmFPAfileDefine.h	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/camera/pmFPAfileDefine.h	(revision 35421)
@@ -81,4 +81,14 @@
     const char *filename                ///< Name of file
     );
+
+// find the file associated with the argname & generate a pmFPAfile for it based on the filerule
+pmFPAfile *pmFPAfileDefineNewConfig(
+    bool *success,                      ///< Found files?
+    pmConfig **outConfig, 		///< output configuration for this file
+    pmConfig *sysConfig, 		///< existing system config info
+    const char *filename, 		///< name of filerule
+    const char *argname			///< argument entry
+    );
+
 
 // look for the given argname on the argument list.  find the give filename from the file rules
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_DEV.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_DEV.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_DEV.c	(revision 35421)
@@ -295,5 +295,5 @@
 // make an initial guess for parameters
 // 0.5 PIX: moments and peaks are in pixel coords, thus so are model parameters
-bool PM_MODEL_GUESS (pmModel *model, pmSource *source)
+bool PM_MODEL_GUESS (pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal)
 {
     psF32 *PAR  = model->params->data.F32;
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_DEV.h
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_DEV.h	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_DEV.h	(revision 35421)
@@ -5,5 +5,5 @@
 psF32 pmModelFunc_DEV(psVector *deriv, const psVector *params, const psVector *pixcoord);
 bool pmModelLimits_DEV(psMinConstraintMode mode, int nParam, float *params, float *beta);
-bool pmModelGuess_DEV(pmModel *model, pmSource *source);
+bool pmModelGuess_DEV(pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal);
 psF64 pmModelFlux_DEV(const psVector *params);
 psF64 pmModelRadius_DEV(const psVector *params, psF64 flux);
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_EXP.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_EXP.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_EXP.c	(revision 35421)
@@ -287,5 +287,5 @@
 // make an initial guess for parameters
 // 0.5 PIX: moments and peaks are in pixel coords, thus so are model parameters
-bool PM_MODEL_GUESS (pmModel *model, pmSource *source)
+bool PM_MODEL_GUESS (pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal)
 {
     psF32 *PAR  = model->params->data.F32;
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_EXP.h
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_EXP.h	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_EXP.h	(revision 35421)
@@ -5,5 +5,5 @@
 psF32 pmModelFunc_EXP(psVector *deriv, const psVector *params, const psVector *pixcoord);
 bool pmModelLimits_EXP(psMinConstraintMode mode, int nParam, float *params, float *beta);
-bool pmModelGuess_EXP(pmModel *model, pmSource *source);
+bool pmModelGuess_EXP(pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal);
 psF64 pmModelFlux_EXP(const psVector *params);
 psF64 pmModelRadius_EXP(const psVector *params, psF64 flux);
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_GAUSS.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_GAUSS.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_GAUSS.c	(revision 35421)
@@ -192,5 +192,5 @@
 // make an initial guess for parameters
 // 0.5 PIX: moments and peaks are in pixel coords, thus so are model parameters
-bool PM_MODEL_GUESS (pmModel *model, pmSource *source)
+bool PM_MODEL_GUESS (pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal)
 {
     psF32 *PAR  = model->params->data.F32;
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_GAUSS.h
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_GAUSS.h	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_GAUSS.h	(revision 35421)
@@ -5,5 +5,5 @@
 psF32 pmModelFunc_GAUSS(psVector *deriv, const psVector *params, const psVector *pixcoord);
 bool pmModelLimits_GAUSS(psMinConstraintMode mode, int nParam, float *params, float *beta);
-bool pmModelGuess_GAUSS(pmModel *model, pmSource *source);
+bool pmModelGuess_GAUSS(pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal);
 psF64 pmModelFlux_GAUSS(const psVector *params);
 psF64 pmModelRadius_GAUSS(const psVector *params, psF64 flux);
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_PGAUSS.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_PGAUSS.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_PGAUSS.c	(revision 35421)
@@ -193,5 +193,5 @@
 // make an initial guess for parameters
 // 0.5 PIX: moments and peaks are in pixel coords, thus so are model parameters
-bool PM_MODEL_GUESS (pmModel *model, pmSource *source)
+bool PM_MODEL_GUESS (pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal)
 {
     psF32 *PAR  = model->params->data.F32;
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_PGAUSS.h
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_PGAUSS.h	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_PGAUSS.h	(revision 35421)
@@ -5,5 +5,5 @@
 psF32 pmModelFunc_PGAUSS(psVector *deriv, const psVector *params, const psVector *pixcoord);
 bool pmModelLimits_PGAUSS(psMinConstraintMode mode, int nParam, float *params, float *beta);
-bool pmModelGuess_PGAUSS(pmModel *model, pmSource *source);
+bool pmModelGuess_PGAUSS(pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal);
 psF64 pmModelFlux_PGAUSS(const psVector *params);
 psF64 pmModelRadius_PGAUSS(const psVector *params, psF64 flux);
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_PS1_V1.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_PS1_V1.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_PS1_V1.c	(revision 35421)
@@ -212,5 +212,5 @@
 // make an initial guess for parameters
 // 0.5 PIX: moments and peaks are in pixel coords, thus so are model parameters
-bool PM_MODEL_GUESS (pmModel *model, pmSource *source)
+bool PM_MODEL_GUESS (pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal)
 {
     psF32 *PAR  = model->params->data.F32;
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_PS1_V1.h
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_PS1_V1.h	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_PS1_V1.h	(revision 35421)
@@ -5,5 +5,5 @@
 psF32 pmModelFunc_PS1_V1(psVector *deriv, const psVector *params, const psVector *pixcoord);
 bool pmModelLimits_PS1_V1(psMinConstraintMode mode, int nParam, float *params, float *beta);
-bool pmModelGuess_PS1_V1(pmModel *model, pmSource *source);
+bool pmModelGuess_PS1_V1(pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal);
 psF64 pmModelFlux_PS1_V1(const psVector *params);
 psF64 pmModelRadius_PS1_V1(const psVector *params, psF64 flux);
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_QGAUSS.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_QGAUSS.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_QGAUSS.c	(revision 35421)
@@ -213,5 +213,5 @@
 // make an initial guess for parameters
 // 0.5 PIX: moments and peaks are in pixel coords, thus so are model parameters
-bool PM_MODEL_GUESS (pmModel *model, pmSource *source)
+bool PM_MODEL_GUESS (pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal)
 {
     psF32 *PAR  = model->params->data.F32;
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_QGAUSS.h
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_QGAUSS.h	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_QGAUSS.h	(revision 35421)
@@ -5,5 +5,5 @@
 psF32 pmModelFunc_QGAUSS(psVector *deriv, const psVector *params, const psVector *pixcoord);
 bool pmModelLimits_QGAUSS(psMinConstraintMode mode, int nParam, float *params, float *beta);
-bool pmModelGuess_QGAUSS(pmModel *model, pmSource *source);
+bool pmModelGuess_QGAUSS(pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal);
 psF64 pmModelFlux_QGAUSS(const psVector *params);
 psF64 pmModelRadius_QGAUSS(const psVector *params, psF64 flux);
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_RGAUSS.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_RGAUSS.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_RGAUSS.c	(revision 35421)
@@ -202,5 +202,5 @@
 // make an initial guess for parameters
 // 0.5 PIX: moments and peaks are in pixel coords, thus so are model parameters
-bool PM_MODEL_GUESS (pmModel *model, pmSource *source)
+bool PM_MODEL_GUESS (pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal)
 {
     psF32 *PAR  = model->params->data.F32;
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_RGAUSS.h
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_RGAUSS.h	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_RGAUSS.h	(revision 35421)
@@ -5,5 +5,5 @@
 psF32 pmModelFunc_RGAUSS(psVector *deriv, const psVector *params, const psVector *pixcoord);
 bool pmModelLimits_RGAUSS(psMinConstraintMode mode, int nParam, float *params, float *beta);
-bool pmModelGuess_RGAUSS(pmModel *model, pmSource *source);
+bool pmModelGuess_RGAUSS(pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal);
 psF64 pmModelFlux_RGAUSS(const psVector *params);
 psF64 pmModelRadius_RGAUSS(const psVector *params, psF64 flux);
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_SERSIC.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_SERSIC.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_SERSIC.c	(revision 35421)
@@ -307,5 +307,5 @@
 // make an initial guess for parameters
 // 0.5 PIX: moments and peaks are in pixel coords, thus so are model parameters
-bool PM_MODEL_GUESS (pmModel *model, pmSource *source)
+bool PM_MODEL_GUESS (pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal)
 {
     pmMoments *moments = source->moments;
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_SERSIC.h
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_SERSIC.h	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_SERSIC.h	(revision 35421)
@@ -5,5 +5,5 @@
 psF32 pmModelFunc_SERSIC(psVector *deriv, const psVector *params, const psVector *pixcoord);
 bool pmModelLimits_SERSIC(psMinConstraintMode mode, int nParam, float *params, float *beta);
-bool pmModelGuess_SERSIC(pmModel *model, pmSource *source);
+bool pmModelGuess_SERSIC(pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal);
 psF64 pmModelFlux_SERSIC(const psVector *params);
 psF64 pmModelRadius_SERSIC(const psVector *params, psF64 flux);
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_TRAIL.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_TRAIL.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_TRAIL.c	(revision 35421)
@@ -168,4 +168,13 @@
         dPAR[PM_PAR_THETA]  = PAR[PM_PAR_I0] * dPdT;
         dPAR[PM_PAR_SIGMA]  = 0;	// we don't actually allow this to vary, so we do not need to calculate it
+
+	for (int i = 0; i < 7; i++) {
+	  if (isnan(dPAR[i])) {
+	    fprintf (stderr, "*");
+	  }
+	}
+    }
+    if (isnan(f)) {
+      fprintf (stderr, "!");
     }
     return(f);
@@ -228,7 +237,111 @@
 }
 
+# define NA 21
+# define NR 21
+static float flux[NA][NR];
+static float npix[NA][NR];
+
+bool pmTrailGetAngle (float *To, float *Io, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal, float sigma) {
+
+  float Xo, Yo;
+  if (!pmModelSetPosition(&Xo, &Yo, source)) return false;
+
+  psImage *image = source->pixels;
+  psImage *mask = source->maskObj;
+  psF32 **imData = image->data.F32;
+  psImageMaskType **mkData = mask->data.PS_TYPE_IMAGE_MASK_DATA;
+
+  // do a loop over the pixels, generating (dX,dY) dot (cos(theta),sin(theta))
+  int dP;
+  int dX = Xo - mask->col0;
+  dP = mask->numCols - dX;
+  int DX = PS_MAX(dP, dX);
+  int NX = mask->numCols;
+
+  int dY = Yo - mask->row0;
+  dP = mask->numRows - dY;
+  int DY = PS_MAX(dP, dY);
+  int NY = mask->numRows;
+
+  // just hard wire this for now...
+  float radius = 10.0;
+  float radius2 = PS_SQR(radius);
+
+  // we have an array of Angles x Radii
+  float dT = M_PI / NA;
+  for (int na = 0; na < NA; na++) {
+    memset (flux[na], 0, NR*sizeof(float));
+    memset (npix[na], 0, NR*sizeof(float));
+  }
+
+  // we skip any pixels [real or virtual] outside of the specified radius (nominally the aperture radius)
+  // ix and iy track pixels relative to the centroid
+  for (int ix = -DX; ix < DX + 1; ix++) {
+    if (ix > radius) continue;
+    int mx = ix + dX;
+    for (int iy = -DY; iy < DY + 1; iy++) {
+      if (iy > radius) continue;
+      if (ix*ix + iy*iy > radius2) continue;
+      int my = iy + dY;
+      
+      // include count only the unmasked pixels within the image area
+      if (mx < 0) continue;
+      if (my < 0) continue;
+      if (mx >= NX) continue;
+      if (my >= NY) continue;
+      
+      // count pixels which are masked only with bad pixels
+      if (mkData[my][mx] & maskVal)continue;
+
+      // we have defined NA to be 21
+      int na = 0;
+      for (float angle = 0.0; na < NA; angle += dT, na ++) {
+
+	// XXX optimization : pre-compute the angle sines and cosines
+	float Rad = (ix * cos(angle)) + (iy * sin(angle));
+	int nr = PS_MAX (PS_MIN (NR, Rad + 0.5*NR), 0);
+
+	flux[na][nr] += imData[my][mx];
+	npix[na][nr] ++;
+      }
+    }
+  }
+
+  // generate a psf model (with integral of 1.0)
+  float Ao = 1.0 / sqrt(2*M_PI) / sigma;
+  float psf[NR];
+  for (int nr = 0; nr < NR; nr++) {
+    psf[nr] = Ao*exp(-0.5*PS_SQR((nr - 0.5*NR)/sigma));
+  }
+
+  float fangle[NA];
+  for (int na = 0; na < NA; na++) {
+    float fpsf = 0.0;
+    for (int nr = 0; nr < NR; nr++) {
+      fpsf += psf[nr]*flux[na][nr];
+    }
+    fangle[na] = fpsf;
+    // fprintf (stderr, "fpsf: %f, theta = %f\n", fpsf, PS_DEG_RAD*dT*na);
+  }
+
+  float peak = fangle[0];
+  int pbin = 0;
+  for (int na = 0; na < NA; na++) {
+    if (fangle[na] > peak) {
+      peak = fangle[na];
+      pbin = na;
+    }
+  }
+
+  // result is peak, pbin -> angle
+  *To = dT * pbin;
+  *Io = peak * Ao;
+  
+  return true;
+}
+
 // make an initial guess for parameters
 // 0.5 PIX: moments and peaks are in pixel coords, thus so are model parameters
-bool PM_MODEL_GUESS (pmModel *model, pmSource *source)
+bool PM_MODEL_GUESS (pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal)
 {
     psF32 *PAR  = model->params->data.F32;
@@ -260,4 +373,7 @@
     //else { size = psfAxes.major; }
 
+    float theta, peak;
+    pmTrailGetAngle (&theta, &peak, source, maskVal, markVal, psfAxes.major);
+
     // axes.major is a sigma in the major direction; scale to 
     PAR[PM_PAR_LENGTH] = 1.5*2.35*size; // a tophat of length L has L = 1.5 * 2.35 * sigma
@@ -266,7 +382,8 @@
 
     // set the model normalization
-    if (!pmModelSetNorm(&PAR[PM_PAR_I0], source)) {
-      return false;
-    }
+    // if (!pmModelSetNorm(&PAR[PM_PAR_I0], source)) {
+    //   return false;
+    // }
+    PAR[PM_PAR_I0] = peak;
 
     // set the model position
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_TRAIL.h
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_TRAIL.h	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/models/pmModel_TRAIL.h	(revision 35421)
@@ -5,5 +5,5 @@
 psF32 pmModelFunc_TRAIL(psVector *deriv, const psVector *params, const psVector *pixcoord);
 bool pmModelLimits_TRAIL(psMinConstraintMode mode, int nParam, float *params, float *beta);
-bool pmModelGuess_TRAIL(pmModel *model, pmSource *source);
+bool pmModelGuess_TRAIL(pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal);
 psF64 pmModelFlux_TRAIL(const psVector *params);
 psF64 pmModelRadius_TRAIL(const psVector *params, psF64 flux);
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/pmModelFuncs.h
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/pmModelFuncs.h	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/pmModelFuncs.h	(revision 35421)
@@ -107,5 +107,5 @@
 //  This function provides the model guess parameters based on the details of
 //  the given source.
-typedef bool (*pmModelGuessFunc)(pmModel *model, pmSource *source);
+typedef bool (*pmModelGuessFunc)(pmModel *model, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal);
 
 //  This function constructs the PSF model for the given source based on the
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/pmPSF.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/pmPSF.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/pmPSF.c	(revision 35421)
@@ -350,5 +350,7 @@
 	shape.sx  = modelPar[PM_PAR_SXX];
 	shape.sy  = modelPar[PM_PAR_SYY];
-	shape.sxy = modelPar[PM_PAR_SXY];
+	shape.sxy = modelPar[PM_PAR_SXY] / 2.0;
+	// XXX I *think* dividing by 2.0 is the right direction, but this 
+	// needs to be checked with a real test
     } else {
 	shape.sx  = modelPar[PM_PAR_SXX] / M_SQRT2;
@@ -392,5 +394,5 @@
 	modelPar[PM_PAR_SXX] = shape.sx;
 	modelPar[PM_PAR_SYY] = shape.sy;
-	modelPar[PM_PAR_SXY] = shape.sxy;
+	modelPar[PM_PAR_SXY] = shape.sxy; // XXX NEED factor of 2 here for correct angle conversion
     } else {
 	modelPar[PM_PAR_SXX] = shape.sx * M_SQRT2;
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/pmPSFtryFitEXT.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/pmPSFtryFitEXT.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/pmPSFtryFitEXT.c	(revision 35421)
@@ -83,5 +83,5 @@
         }
 
-        source->modelEXT = pmSourceModelGuess (source, options->type);
+        source->modelEXT = pmSourceModelGuess (source, options->type, maskVal, markVal);
         if (source->modelEXT == NULL) {
             psTrace ("psModules.objects", 4, "masking %d (%d,%d) : failed to generate model guess\n", i, source->peak->x, source->peak->y);
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceFitModel.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceFitModel.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceFitModel.c	(revision 35421)
@@ -75,5 +75,7 @@
     PS_ASSERT_PTR_NON_NULL(source->pixels, false);
     PS_ASSERT_PTR_NON_NULL(source->maskObj, false);
-    PS_ASSERT_PTR_NON_NULL(source->variance, false);
+
+    // XXX if variance is NULL, use pixels instead (wrong, but not badly wrong)
+    // PS_ASSERT_PTR_NON_NULL(source->variance, false);
 
     psBool fitStatus = true;
@@ -93,4 +95,7 @@
     float Xo = model->params->data.F32[PM_PAR_XPOS];
     float Yo = model->params->data.F32[PM_PAR_YPOS];
+
+    // if variance is NULL, we pretend pixels == variance
+    float **vWgt = source->variance ? source->variance->data.F32 : source->pixels->data.F32;
 
     // fill in the coordinate and value entries
@@ -103,15 +108,16 @@
             }
             // skip zero-variance points
-            if (source->variance->data.F32[i][j] == 0) {
+            if (vWgt[i][j] == 0) {
                 continue;
             }
             // skip nan values in image
             if (!isfinite(source->pixels->data.F32[i][j])) {
+		fprintf (stderr, "WARNING: unmasked nan in image : %x vs %x\n", source->maskObj->data.PS_TYPE_IMAGE_MASK_DATA[i][j], maskVal);
                 continue;
             }
 
             // skip nan values in image
-            if (!isfinite(source->variance->data.F32[i][j])) {
-		fprintf (stderr, "impossible! %x vs %x\n", source->maskObj->data.PS_TYPE_IMAGE_MASK_DATA[i][j], maskVal);
+            if (!isfinite(vWgt[i][j])) {
+		fprintf (stderr, "WARNING: unmasked nan in variance : %x vs %x\n", source->maskObj->data.PS_TYPE_IMAGE_MASK_DATA[i][j], maskVal);
 		continue;
             }
@@ -140,5 +146,5 @@
             // source sky variance
             if (options->poissonErrors) {
-                yErr->data.F32[nPix] = 1.0 / source->variance->data.F32[i][j];
+                yErr->data.F32[nPix] = 1.0 / vWgt[i][j];
             } else {
                 yErr->data.F32[nPix] = 1.0 / options->weight;
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceFitPCM.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceFitPCM.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceFitPCM.c	(revision 35421)
@@ -145,5 +145,5 @@
 bool pmSourceModelGuessPCM (pmPCMdata *pcm, pmSource *source, psImageMaskType maskVal, psImageMaskType markVal) {
 
-    if (!pcm->modelConv->modelGuess(pcm->modelConv, source)) {
+  if (!pcm->modelConv->modelGuess(pcm->modelConv, source, maskVal, markVal)) {
 	return false;
     }
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceIO.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceIO.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceIO.c	(revision 35421)
@@ -980,4 +980,5 @@
 	    return false;
         }
+
         // if this is not TRUE, the output files only contain the psf measurements.
         bool XSRC_OUTPUT = psMetadataLookupBool(&status, recipe, "EXTENDED_SOURCE_ANALYSIS");
@@ -1005,4 +1006,15 @@
             hdu->header = psFitsReadHeader (NULL, file->fits);
         }
+
+	// "WCS" is CMF without detected objects
+	if (file->type == PM_FPA_FILE_WCS) {
+	  psTrace("psModules.objects", 6, "read CMF table from %s : %s : %s", file->filename, headname, dataname);
+	  psFree (headname);
+	  psFree (dataname);
+	  psFree (deteffname);
+	  break;
+	}
+
+	// EXTDATA is the PSF data associated with this image header
 
         // we need to find the corresponding table EXTNAME.
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceMoments.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceMoments.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceMoments.c	(revision 35421)
@@ -151,5 +151,5 @@
 
 	float *vPix = source->pixels->data.F32[row];
-	float *vWgt = source->variance->data.F32[row];
+	float *vWgt = source->variance ? source->variance->data.F32[row] : source->pixels->data.F32[row];
 
 	psImageMaskType  *vMsk = (source->maskObj == NULL) ? NULL : source->maskObj->data.PS_TYPE_IMAGE_MASK_DATA[row];
@@ -266,6 +266,6 @@
 
     float **vPix = source->pixels->data.F32;
-    float **vWgt = source->variance->data.F32;
-    psImageMaskType  **vMsk = (source->maskObj == NULL) ? NULL : source->maskObj->data.PS_TYPE_IMAGE_MASK_DATA;
+    float **vWgt = source->variance ? source->variance->data.F32 : source->pixels->data.F32;
+    psImageMaskType **vMsk = (source->maskObj == NULL) ? NULL : source->maskObj->data.PS_TYPE_IMAGE_MASK_DATA;
 
     float RF = 0.0;
@@ -467,5 +467,5 @@
 
 	float *vPix = source->pixels->data.F32[row];
-	float *vWgt = source->variance->data.F32[row];
+	float *vWgt = source->variance ? source->variance->data.F32[row] : source->pixels->data.F32[row];
 
 	psImageMaskType *vMsk = (source->maskObj == NULL) ? NULL : source->maskObj->data.PS_TYPE_IMAGE_MASK_DATA[row];
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceSky.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceSky.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceSky.c	(revision 35421)
@@ -122,8 +122,13 @@
     psTrace("psModules.objects", 10, "---- %s() begin ----\n", __func__);
     PS_ASSERT_PTR_NON_NULL(source, false);
-    PS_ASSERT_IMAGE_NON_NULL(source->variance, false);
     PS_ASSERT_IMAGE_NON_NULL(source->maskObj, false);
     PS_ASSERT_PTR_NON_NULL(source->peak, false);
     PS_ASSERT_INT_POSITIVE(Radius, false);
+
+    if (source->variance == NULL) {
+      // XXX this is needed so psphotCheckRadiusPSF has a real value
+      source->moments->dSky = 1.0;
+      return true;
+    }
 
     // maskVal is used to test for rejected pixels, and must include markVal
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceUtils.c
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceUtils.c	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceUtils.c	(revision 35421)
@@ -46,8 +46,8 @@
 
 /******************************************************************************
-    pmSourceModelGuess(source, model): This function allocates a new
-    pmModel structure based on the given modelType specified in the argument list.
-    The corresponding pmModelGuess function is returned, and used to
-    supply the values of the params array in the pmModel structure.
+    pmSourceModelGuess(source, model, maskVal, markVal): This function allocates a new
+    pmModel structure based on the given modelType specified in the argument list.  The
+    corresponding pmModelGuess function is returned, and used to supply the values of the
+    params array in the pmModel structure.
  
     XXX: Many parameters are based on the src->moments structure, which is in
@@ -55,6 +55,5 @@
     functions will be in image, not subImage coords.  Remember this.
 *****************************************************************************/
-pmModel *pmSourceModelGuess(pmSource *source,
-                            pmModelType modelType)
+pmModel *pmSourceModelGuess(pmSource *source, pmModelType modelType, psImageMaskType maskVal, psImageMaskType markVal)
 {
     psTrace("psModules.objects", 10, "---- %s() begin ----\n", __func__);
@@ -65,5 +64,5 @@
     pmModel *model = pmModelAlloc(modelType);
 
-    if (!model->modelGuess(model, source)) {
+    if (!model->modelGuess(model, source, maskVal, markVal)) {
 	psFree (model);
 	return NULL;
Index: branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceUtils.h
===================================================================
--- branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceUtils.h	(revision 35420)
+++ branches/eam_branches/ipp-20130419/psModules/src/objects/pmSourceUtils.h	(revision 35421)
@@ -27,5 +27,7 @@
 pmModel *pmSourceModelGuess(
     pmSource *source,   ///< The input pmSource
-    pmModelType model   ///< The type of model to be created.
+    pmModelType model,   ///< The type of model to be created.
+    psImageMaskType maskVal, 
+    psImageMaskType markVal
 );