Index: trunk/pswarp/src/Makefile.am
===================================================================
--- trunk/pswarp/src/Makefile.am	(revision 10952)
+++ trunk/pswarp/src/Makefile.am	(revision 10954)
@@ -19,4 +19,5 @@
 pswarpErrorCodes.c         \
 pswarpParseCamera.c   	   \
+pswarpTransformReadout.c   \
 pswarpVersion.c            
 
Index: trunk/pswarp/src/pswarp.h
===================================================================
--- trunk/pswarp/src/pswarp.h	(revision 10952)
+++ trunk/pswarp/src/pswarp.h	(revision 10954)
@@ -13,4 +13,15 @@
 # define PSWARP_RECIPE "PSWARP" // Name of the recipe to use
 
+typedef struct {
+  double Xo, Xx, Xy;
+  double Yo, Yx, Yy;
+} pswarpMap;
+
+typedef struct {
+  pswarpMap **grid;
+  int Nx, Ny;
+  int nXpix, nYpix;
+} pswarpMapGrid;
+
 pmConfig *pswarpArguments (int argc, char **argv);
 bool pswarpParseCamera (pmConfig *config);
@@ -19,2 +30,4 @@
 bool pswarpDataSave (pmConfig *config);
 void pswarpCleanup (pmConfig *config);
+bool pswarpTransformReadout (pmReadout *output, pmReadout *input, pmConfig *config);
+
Index: trunk/pswarp/src/pswarpArguments.c
===================================================================
--- trunk/pswarp/src/pswarpArguments.c	(revision 10952)
+++ trunk/pswarp/src/pswarpArguments.c	(revision 10954)
@@ -50,14 +50,19 @@
     pmConfigFileSetsMD (config->arguments, config, "WEIGHT", "-weight", "-weightlist");
     
-    if (argc != 2) {
+    if (argc != 3) {
 	psError(PSWARP_ERR_ARGUMENTS, true, "Incorrect arguments supplied");
 	return NULL;
     }
     
+    psArray *array;
+
     // output position is fixed
     psMetadataAddStr (config->arguments, PS_LIST_TAIL, "OUTPUT", 0, "", argv[1]);
 
     // skycell position is fixed
-    psMetadataAddStr (config->arguments, PS_LIST_TAIL, "SKYCELL", 0, "", argv[2]);
+    array = psArrayAlloc(1);
+    array->data[0] = psStringCopy (argv[2]);
+    status = psMetadataAddPtr (config->arguments, PS_LIST_TAIL, "SKYCELL", PS_DATA_ARRAY, "", array);
+    psFree (array);
 
     psTrace("pswarp", 1, "Done with pswarpArguments...\n");
Index: trunk/pswarp/src/pswarpDataLoad.c
===================================================================
--- trunk/pswarp/src/pswarpDataLoad.c	(revision 10952)
+++ trunk/pswarp/src/pswarpDataLoad.c	(revision 10954)
@@ -32,9 +32,24 @@
 
     // de-activate PSWARP.SKYCELL and PSWARP.OUTPUT
-    // 
+    pmFPAfileActivate (config->files, false, "PSWARP.SKYCELL");
     pmFPAfileActivate (config->files, false, "PSWARP.OUTPUT");
-    pmFPAfileActivate (config->files, false, "PSWARP.SKYCELL");
 
     pmFPAview *view = pmFPAviewAlloc (0);
+
+    // find the FPA phu
+    // XXX wrap the test below into a function
+    // XXX need to read only the headers for the skycell
+    // XXX need to optionally load the astrometry datafile
+    bool bilevelAstrometry = false;
+    pmHDU *phu = pmFPAviewThisPHU (view, input->fpa);
+    if (phu) {
+      char *ctype = psMetadataLookupStr (NULL, phu->header, "CTYPE1");
+      if (ctype) {
+	bilevelAstrometry = !strcmp (&ctype[4], "-DIS");
+      }
+    }
+    if (bilevelAstrometry) {
+      pmAstromReadBilevelMosaic (input->fpa, phu->header);
+    } 
 
     // files associated with the science image
@@ -45,4 +60,14 @@
         if (!chip->process || !chip->file_exists) { continue; }
 	pmFPAfileIOChecks (config, view, PM_FPA_BEFORE);
+
+        // read WCS data from the corresponding header
+	pmHDU *hdu = pmFPAviewThisHDU (view, input->fpa);
+	if (bilevelAstrometry) {
+	  pmAstromReadBilevelChip (chip, hdu->header); 
+	} else {
+	  // XXX get pixelScale from recipes.  does it matter?
+	  float pixelScale = 13.5;
+	  pmAstromReadWCS (input->fpa, chip, hdu->header, pixelScale);
+	}
 
 	while ((cell = pmFPAviewNextCell (view, input->fpa, 1)) != NULL) {
@@ -56,5 +81,5 @@
 		if (! readout->data_exists) { continue; }
 
-		pswarpConvertReadout (output, readout, config);
+		pswarpTransformReadout (output, readout, config);
 
 		pmFPAfileIOChecks (config, view, PM_FPA_AFTER);
Index: trunk/pswarp/src/pswarpDefine.c
===================================================================
--- trunk/pswarp/src/pswarpDefine.c	(revision 10952)
+++ trunk/pswarp/src/pswarpDefine.c	(revision 10954)
@@ -23,4 +23,25 @@
     pmFPAfileRead (skycell, view, config);
 
+    // XXX this is not a sufficient test
+    view->chip = 0;
+    pmChip *chip = pmFPAviewThisChip (view, skycell->fpa);
+    pmHDU *phu = pmFPAviewThisPHU (view, skycell->fpa);
+    pmHDU *hdu = pmFPAviewThisHDU (view, skycell->fpa);
+    bool bilevelAstrometry = false;
+    if (phu) {
+      char *ctype = psMetadataLookupStr (NULL, phu->header, "CTYPE1");
+      if (ctype) {
+	bilevelAstrometry = !strcmp (&ctype[4], "-DIS");
+      }
+    }
+    if (bilevelAstrometry) {
+      pmAstromReadBilevelMosaic (skycell->fpa, phu->header);
+      pmAstromReadBilevelChip (chip, hdu->header); 
+    } else {
+      // XXX get pixelScale from recipes.  does it matter?
+      float pixelScale = 13.5;
+      pmAstromReadWCS (skycell->fpa, chip, hdu->header, pixelScale);
+    }
+
     // construct the output fpa here
     pmFPAfile *output = psMetadataLookupPtr (NULL, config->files, "PSWARP.OUTPUT");
Index: trunk/pswarp/src/pswarpMapGrid.c
===================================================================
--- trunk/pswarp/src/pswarpMapGrid.c	(revision 10954)
+++ trunk/pswarp/src/pswarpMapGrid.c	(revision 10954)
@@ -0,0 +1,15 @@
+# include "pswarp.h"
+
+static pswarpMapFree (pswarpMap *map) {
+  return;
+}
+
+pswarpMap *pswarpMapAlloc () {
+
+  pswarpMap *map = (pswarpMap *) psAlloc (sizeof(pswarpMap));
+  psMemSetDeallocator(map, (psFreeFunc) pswarpMapFree);
+
+  return map;
+}
+
+  
Index: trunk/pswarp/src/pswarpParseCamera.c
===================================================================
--- trunk/pswarp/src/pswarpParseCamera.c	(revision 10952)
+++ trunk/pswarp/src/pswarpParseCamera.c	(revision 10954)
@@ -3,9 +3,10 @@
 bool pswarpParseCamera (pmConfig *config) {
 
-    bool status = false;
+    bool status;
 
     // the input image(s) are required arguments; they define the camera
+    status = false;
     pmFPAfile *input = pmFPAfileDefineFromArgs (&status, config, "PSWARP.INPUT", "INPUT");
-    if (!status) {
+    if (!input) {
 	psError(PSWARP_ERR_CONFIG, false, "Failed to build FPA from PSWARP.INPUT");
 	return false;
@@ -13,11 +14,13 @@
 
     // the input skycell is a required argument: it defines the output image
+    status = false;
     pmFPAfile *skycell = pmFPAfileDefineFromArgs (&status, config, "PSWARP.SKYCELL", "SKYCELL");
-    if (!status) {
-	psError(PSWARP_ERR_CONFIG, false, "Failed to build FPA from PSWARP.INPUT");
+    if (!skycell) {
+	psError(PSWARP_ERR_CONFIG, false, "Failed to build FPA from PSWARP.SKYCELL");
 	return false;
     }
 
     // the mask is not required; the  skycell is a required argument: it defines the output image
+    status = false;
     pmFPAfileDefineFromArgs (&status, config, "PSWARP.MASK", "MASK");
     if (!status) {
@@ -26,4 +29,5 @@
 
     // the weight is not required; the  skycell is a required argument: it defines the output image
+    status = false;
     pmFPAfileDefineFromArgs (&status, config, "PSWARP.WEIGHT", "WEIGHT");
     if (!status) {
@@ -32,5 +36,8 @@
 
     // these calls bind the I/O handle to the specified fpa
-    pmFPAfileDefineOutput (config, skycell->fpa, "PSWARP.OUTPUT");
+    if (!pmFPAfileDefineOutput (config, skycell->fpa, "PSWARP.OUTPUT")) {
+	psError(PSWARP_ERR_CONFIG, false, "Failed to build FPA from PSWARP.OUTPUT");
+	return false;
+    }
 
     // Chip selection: turn on only the chips specified
Index: trunk/pswarp/src/pswarpTransformReadout.c
===================================================================
--- trunk/pswarp/src/pswarpTransformReadout.c	(revision 10952)
+++ trunk/pswarp/src/pswarpTransformReadout.c	(revision 10954)
@@ -8,5 +8,5 @@
 
     // select the current recipe
-    psMetadata *recipe  = psMetadataLookupPtr (NULL, config->recipes, PSPHOT_RECIPE);
+    // psMetadata *recipe  = psMetadataLookupPtr (NULL, config->recipes, PSWARP_RECIPE);
 
     int outNx = output->image->numCols;
@@ -27,4 +27,7 @@
     pmChip *chipOutput = cell->parent;
     pmFPA *fpaOutput = chipOutput->parent;
+
+    psF32 **outData = output->image->data.F32;
+    psImage *inImage = input->image;
 
     // Iterate over the output image pixels
@@ -48,6 +51,9 @@
 
 		// XXX get interpolation method from the recipe
-		outImage->data.F32[y][x] = (psF32)psImagePixelInterpolate(image, inPix->x, inPix->y, mask, 1, NAN, PS_INTERPOLATE_BILINEAR);
+		outData[y][x] = (psF32)psImagePixelInterpolate(inImage, inPix->x, inPix->y, NULL, 1, NAN, PS_INTERPOLATE_BILINEAR);
+		// modify zero and scale?
 
+		// outData[y][x] = (psF32)psImagePixelInterpolate(inImage, inPix->x, inPix->y, mask, 1, NAN, PS_INTERPOLATE_BILINEAR);
+		
 # if (0)
 		if (error) {
@@ -58,8 +64,4 @@
 											    mask, 1, NAN);
 		}
-# endif
-
-		outImage->data.F32[y][x] = (outImage->data.F32[y][x] - offset) / scale;
-# if (0)
 		if (error) {
 		    outError->data.F32[y][x] = outError->data.F32[y][x] / SQUARE(scale);
Index: trunk/pswarp/src/pswarpTransformReadout_Opt.c
===================================================================
--- trunk/pswarp/src/pswarpTransformReadout_Opt.c	(revision 10954)
+++ trunk/pswarp/src/pswarpTransformReadout_Opt.c	(revision 10954)
@@ -0,0 +1,86 @@
+# include "pswarp.h"
+
+bool pswarpTransformReadout (pmReadout *output, pmReadout *input, pmConfig *config) {
+
+    // XXX this implementation currently ignores the use of the region 
+    psImage *region = NULL;
+    pmCell *cell = NULL;
+
+    // select the current recipe
+    psMetadata *recipe  = psMetadataLookupPtr (NULL, config->recipes, PSPHOT_RECIPE);
+
+    int outNx = output->image->numCols;
+    int outNy = output->image->numRows;
+
+    psPlane *inPix = psPlaneAlloc();	// Coordinates on the input detector
+    psPlane *outPix = psPlaneAlloc();	// Coordinates on the output detector
+
+    psPlane *FP = psPlaneAlloc();	// Coordinates on the focal plane
+    psPlane *TP = psPlaneAlloc();	// Coordinates on the tangent plane
+    psSphere *sky = psSphereAlloc();	// Coordinates on the sky
+
+    cell = input->parent;
+    pmChip *chipInput = cell->parent;
+    pmFPA *fpaInput = chipInput->parent;
+
+    cell = output->parent;
+    pmChip *chipOutput = cell->parent;
+    pmFPA *fpaOutput = chipOutput->parent;
+
+    psF32 **outData = output->image->data.F32;
+    psImage *inImage = input->image;
+
+    // we might want to do the rectangular regions outside of the selection independently
+    psImageInit (output->image, NAN);
+
+    // find the output pixel range
+    pswarpMatchRange (&minX, &minY, &maxX, &maxY, input, output);
+
+    // Iterate over the output image pixels
+    for (int y = 0; y < outNy; y++) {
+	for (int x = 0; x < outNx; x++) {
+	    // Only transform those pixels requested
+	    if (!region || (region && region->data.U8[y][x])) {
+		// Transform!
+
+		// XXX double check this 1/2 pixel offset
+		outPix->x = (double)x + 0.5;
+		outPix->y = (double)y + 0.5;
+
+		psPlaneTransformApply(FP, chipOutput->toFPA, outPix);
+		psPlaneTransformApply (TP, fpaOutput->toTPA, FP);
+		psDeproject (sky, TP, fpaOutput->toSky);
+		
+		psProject (TP, sky, fpaInput->toSky);
+		psPlaneTransformApply (FP, fpaInput->fromTPA, TP);
+		psPlaneTransformApply (inPix, chipInput->fromFPA, FP);
+
+		// XXX get interpolation method from the recipe
+		outImage->data.F32[y][x] = (psF32)psImagePixelInterpolate(image, inPix->x, inPix->y, mask, 1, NAN, PS_INTERPOLATE_BILINEAR);
+
+# if (0)
+		if (error) {
+		    // Error is actually the variance
+		    outError->data.F32[y][x] = (psF32)p_psImageErrorInterpolateBILINEAR_F32(error,
+											    detector->x,
+											    detector->y,
+											    mask, 1, NAN);
+		}
+# endif
+
+		outImage->data.F32[y][x] = (outImage->data.F32[y][x] - offset) / scale;
+# if (0)
+		if (error) {
+		    outError->data.F32[y][x] = outError->data.F32[y][x] / SQUARE(scale);
+		}
+# endif
+
+	    } // Pixels of interest
+
+	}
+    } // Iterating over output pixels
+
+    return true;
+}
+
+    
