Index: trunk/stac/src/stacRead.c
===================================================================
--- trunk/stac/src/stacRead.c	(revision 3375)
+++ trunk/stac/src/stacRead.c	(revision 3610)
@@ -3,4 +3,6 @@
 #include "pslib.h"
 #include "stac.h"
+
+#define BUFFER 100			// Size of buffer for incrementally reading coordinates
 
 psArray *stacReadImages(stacConfig *config)
@@ -40,4 +42,123 @@
 }
 
+psArray *stacReadCoords(const char *filename)
+{
+    FILE *file = fopen(filename, "r");
+    if (file == NULL) {
+	psLogMsg("stac.read.coords", PS_LOG_ERROR, "Cannot open coordinate file, %s\n", filename);
+	return NULL;
+    }
+
+    psTrace("stac.read.coords", 5, "Reading coordinate file, %s\n", filename);
+
+    psArray *coords = psArrayAlloc(BUFFER); // The array of coordinates to be returned
+    float x, y;				// Coordinates to read
+    int num = 0;			// Number of coordinates read
+    while (fscanf(file, "%f %f\n", &x, &y) == 2) {
+	psPlane *coord = psPlaneAlloc();// A coordinate
+	coord->x = x;
+	coord->y = y;
+	coords->data[num] = coord;
+	num++;
+	if (num % BUFFER) {
+	    coords = psArrayRealloc(coords, num + BUFFER);
+	}
+    }
+    coords->n = num;
+
+    psTrace("stac.read.coords", 5, "%d coordinates read.\n", num);
+
+    fclose(file);
+    return coords;
+}
+
+
+psPlaneTransform *stacReadMap(const char *filename)
+{
+    FILE *mapfp = fopen(filename, "r");
+    if (mapfp == NULL) {
+	psLogMsg("stac.read.map", PS_LOG_ERROR, "Cannot open map file, %s\n", filename);
+	return NULL;
+    }
+    // Read the file
+    psTrace("stac.read.map", 5, "Reading map file %s....\n", filename);
+    
+    // Format is now:
+    // order
+    // x coefficients
+    // y coefficients
+    // 
+    // where the order is 1 for linear, 2 for quadratic, 3 for cubic.
+    // and the coefficients are read by the following pseudo-code:
+    // 
+    // 	for (int k = 0; k < order + 1; k++)
+    // 	    for (int j = 0; j < k; j++)
+    // 		int i = k - j;
+    //              read coefficient of x^i y^j
+    // 
+    // This produces the ordering:
+    // x^0 y^0, x^1 y^0, x^0 y^1, x^2 y^0, x^1 y^1, x^0 y^2, x^3 y^0, x^2 y^1, x^1 y^2, x^0 y^3
+    // 
+    // This is, of course, for third order polynomials.
+    // For lower orders, the list is truncated at the appropriate level.
+
+    int order = 0;			// Polynomial order
+    if (fscanf(mapfp, "%d", &order) != 1) {
+	psLogMsg("stac.read.map", PS_LOG_ERROR, "Unable to read map file %s\n", filename);
+	fclose(mapfp);
+	return NULL;
+    }
+    
+    psTrace("stac.read.map", 5, "Polynomial order: %d\n", order);
+    
+    psPlaneTransform *map = psPlaneTransformAlloc(order + 1, order + 1); // The transformation
+    // Set coefficient masks
+    for (int i = 0; i < order + 1; i++) {
+	for (int j = 0; j < order + 1; j++) {
+	    if (i + j > order + 1) {
+		map->x->mask[i][j] = 1;
+		map->y->mask[i][j] = 1;
+	    } else {
+		map->x->mask[i][j] = 0;
+		map->y->mask[i][j] = 0;
+	    }
+	}
+    }
+    
+    // Read x coefficients
+    psTrace("stac.read.map", 7, "x' = \n");
+    for (int k = 0; k < order + 1; k++) {
+	for (int j = 0; j < k + 1; j++) {
+	    int i = k - j;
+	    if (fscanf(mapfp, "%lf", &map->x->coeff[i][j]) != 1) {
+		psLogMsg("stac.read.map", PS_LOG_ERROR, "Unable to read map file %s\n", filename);
+		fclose(mapfp);
+		psFree(map);
+		return NULL;
+	    }
+	    psTrace("stac.read.map", 7, "      %f x^%d y^%d\n", map->x->coeff[i][j], i, j);
+	}
+    }
+    // Read y coefficients
+    psTrace("stac.read.maps", 7, "y' = \n");
+    for (int k = 0; k < order + 1; k++) {
+	for (int j = 0; j < k + 1; j++) {
+	    int i = k - j;
+	    if (fscanf(mapfp, "%lf", &map->y->coeff[i][j]) != 1) {
+		psLogMsg("stac.read.map", PS_LOG_ERROR, "Unable to read map file %s\n", filename);
+		fclose(mapfp);
+		psFree(map);
+		return NULL;
+	    }
+	    psTrace("stac.read.map", 7, "      %f x^%d y^%d\n", map->y->coeff[i][j], i, j);
+	}
+    }
+    
+    fclose(mapfp);
+
+    return map;
+}
+    
+
 
 psArray *stacReadMaps(stacConfig *config)
@@ -54,84 +175,10 @@
 	    exit(EXIT_FAILURE);
 	}
-	// Open the file
+	// Read the file
 	sprintf(mapfile,"%s.map",filenames->data[i]);
-	FILE *mapfp = fopen(mapfile,"r");
-	if (mapfp == NULL) {
-	    psLogMsg("stac.read.maps", PS_LOG_ERROR, "Cannot open map file, %s\n",mapfile);
-	    exit(EXIT_FAILURE);
+	maps->data[i] = stacReadMap(mapfile);
+	if (maps->data[i] == NULL) {
+	    psLogMsg("stac.read.maps", PS_LOG_ERROR, "Unable to read map: %s\n", mapfile);
 	}
-	// Read the file
-	psTrace("stac.read.maps",5,"Reading map file %s....\n",mapfile);
-
-	// Format is now:
-	// order
-	// x coefficients
-	// y coefficients
-	// 
-	// where the order is 1 for linear, 2 for quadratic, 3 for cubic.
-	// and the coefficients are read by the following pseudo-code:
-	// 
-	// 	for (int k = 0; k < order + 1; k++)
-	// 	    for (int j = 0; j < k; j++)
-	// 		int i = k - j;
-	//              read coefficient of x^i y^j
-	// 
-	// This produces the ordering:
-	// x^0 y^0, x^1 y^0, x^0 y^1, x^2 y^0, x^1 y^1, x^0 y^2, x^3 y^0, x^2 y^1, x^1 y^2, x^0 y^3
-	// 
-	// This is, of course, for third order polynomials.
-	// For lower orders, the list is truncated at the appropriate level.
-
-	int order = 0;			// Polynomial order
-	if (fscanf(mapfp, "%d", &order) != 1) {
-	    psLogMsg("stac.read.maps", PS_LOG_ERROR, "Unable to read map file %s\n", mapfile);
-	    exit(EXIT_FAILURE);
-	}
-
-	psTrace("stac.read.maps", 5, "Polynomial order: %d\n", order);
-
-	psPlaneTransform *map = psPlaneTransformAlloc(order + 1, order + 1); // The transformation
-	// Set coefficient masks
-	for (int i = 0; i < order + 1; i++) {
-	    for (int j = 0; j < order + 1; j++) {
-		if (i + j > order + 1) {
-		    map->x->mask[i][j] = 1;
-		    map->y->mask[i][j] = 1;
-		} else {
-		    map->x->mask[i][j] = 0;
-		    map->y->mask[i][j] = 0;
-		}
-	    }
-	}
-
-	// Read x coefficients
-	psTrace("stac.read.maps", 7, "x' = \n");
-	for (int k = 0; k < order + 1; k++) {
-	    for (int j = 0; j < k + 1; j++) {
-		int i = k - j;
-		if (fscanf(mapfp, "%lf", &map->x->coeff[i][j]) != 1) {
-		    psLogMsg("stac.read.maps", PS_LOG_ERROR, "Unable to read map file %s\n", mapfile);
-		    exit(EXIT_FAILURE);
-		}
-		psTrace("stac.read.maps", 7, "      %f x^%d y^%d\n", map->x->coeff[i][j], i, j);
-	    }
-	}
-	// Read y coefficients
-	psTrace("stac.read.maps", 7, "y' = \n");
-	for (int k = 0; k < order + 1; k++) {
-	    for (int j = 0; j < k + 1; j++) {
-		int i = k - j;
-		if (fscanf(mapfp, "%lf", &map->y->coeff[i][j]) != 1) {
-		    psLogMsg("stac.read.maps", PS_LOG_ERROR, "Unable to read map file %s\n", mapfile);
-		    exit(EXIT_FAILURE);
-		}
-		psTrace("stac.read.maps", 7, "      %f x^%d y^%d\n", map->y->coeff[i][j], i, j);
-	    }
-	}
-
-   	fclose(mapfp);
-
-	// Plug it into the array
-	maps->data[i] = map;
 
     }