Index: trunk/psphot/src/setup.c
===================================================================
--- trunk/psphot/src/setup.c	(revision 4215)
+++ trunk/psphot/src/setup.c	(revision 4216)
@@ -1,19 +1,31 @@
 # include "psphot.h"
+
+// load the image
+// load or construct the mask
+//   apply X/Y border info
+//   apply the SATURATE mask?
+// load or construct the noise 
 
 psImageData *setup (psMetadata *config) 
 {
+    psMetadata *header = NULL;
+    psImage *image = NULL;
+    psImage *noise = NULL;
+    psImage *mask = NULL;
+
     psRegion  region = {0,0,0,0};	// a region representing the entire array
 
     psTimerStart ("psphot");
 
+    // setup header template, specify COMMENT and HISTORY as multis
+    header = psMetadataAlloc ();
+    psMetadataAdd (header, PS_LIST_HEAD, "COMMENT", PS_META_MULTI, "folder for comment", NULL);
+    psMetadataAdd (header, PS_LIST_HEAD, "HISTORY", PS_META_MULTI, "folder for history", NULL);
+
+    // read header and image data from INPUT 
     char *input = psMetadataLookupSTR (&status, config, "INPUT");
     psFits *file = psFitsAlloc (input);
-
-    psImageData *imdata = psAlloc(sizeof(psImageData));
-    psMetadata *header = psMetadataAlloc ();
-    psMetadataAdd (header, PS_LIST_HEAD, "COMMENT", PS_META_MULTI, "folder for comment", NULL);
-    psMetadataAdd (header, PS_LIST_HEAD, "HISTORY", PS_META_MULTI, "folder for history", NULL);
     header = psFitsReadHeader (header, file);
-    psImage *image = psFitsReadImage (NULL, file, region, 0);
+    image = psFitsReadImage (NULL, file, region, 0);
     psFree (file);
 
@@ -25,20 +37,17 @@
     psMetadataAdd (config, PS_LIST_TAIL, "GAIN",    PS_META_F32 | PS_META_REPLACE, "gain value used", GAIN);
 
-    imdata->image = image;
-    imdata->header = header;
-
     // load the noise image if it is supplied, otherwise build from input
     char *noiseName = psMetadataLookupSTR (&status, config, "NOISE");
     if (status == true) {
       file = psFitsAlloc (noiseName);
-      imdata->noise = psFitsReadImage  (NULL, file, region, 0);
+      noise = psFitsReadImage  (NULL, file, region, 0);
       psFree (file);
     } else {
       psScalar *value = psScalarAlloc (1.0 / GAIN, PS_TYPE_F64);
-      imdata->noise = BinaryOp (NULL, imdata->image, '/', value);
+      noise = BinaryOp (NULL, image, '/', value);
       psFree (value);
 
       psScalar *value = psScalarAlloc (PS_SQR(RDNOISE/GAIN), PS_TYPE_F64);
-      imdata->noise = BinaryOp (imdata->noise, imdata->noise, '+', value);
+      noise = BinaryOp (noise, noise, '+', value);
       psFree (value);
     }
@@ -47,9 +56,9 @@
     if (status == true) {
       file = psFitsAlloc (maskName);
-      imdata->mask  = psFitsReadImage  (NULL, file, region, 0);
-      // check if U8 & force?
+      mask  = psFitsReadImage  (NULL, file, region, 0);
+      // require U8
       psFree (file);
     } else {
-      imdata->mask = psImageAlloc (image->numCols, image->numRows, PS_TYPE_U8);
+      mask = psImageAlloc (image->numCols, image->numRows, PS_TYPE_U8);
     }
 
@@ -58,8 +67,24 @@
     psRegion *keep = psRegionAlloc (XBORDER, -XBORDER, YBORDER, -YBORDER);
     keep           = psRegionForImage (keep, image, keep);
-    psImageKeepRegion (imdata->mask, keep, 0x01);
+    psImageKeepRegion (mask, keep, 0x01);
+
+    // mask the saturated pixels
+    float SATURATE = psMetadataLookupF32 (&status, config, "SATURATE");
+    for (int i = 0; i < image->numRows; i++) {
+	for (int j = 0; j < image->numCols; j++) {
+	    if (image->data.F32[i][j] >= SATURATE) {
+		mask->data.U8[i][j] |= 0x02;
+	    }
+	}
+    }
 
     psLogMsg ("psphot", 3, "load data: %f sec\n", psTimerMark ("psphot"));
 
+    psImageData *imdata = psAlloc(sizeof(psImageData));
+    imdata->image = image;
+    imdata->header = header;
+    imdata->noise = noise;
+    imdata->mask = mask;
+
     return (true);
 }