Index: /branches/eam_branches/ipp-20140904/Ohana/src/libdvo/src/AstromOffsetMapOps.c
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/libdvo/src/AstromOffsetMapOps.c	(revision 37673)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/libdvo/src/AstromOffsetMapOps.c	(revision 37674)
@@ -3,4 +3,6 @@
 /* AstromOffsetMap functions:
  */
+
+int dump_map_data (float *x, float *y, float *f, int Npts, char *filename);
 
 float AstromOffsetMapValue (AstromOffsetMap *map, float x, float y, int xdir) {
@@ -108,12 +110,12 @@
   int Ny = map->Ny;
 
-  float **A, **B;
-  ALLOCATE (A, float *, Nx*Ny);
-  ALLOCATE (B, float *, Nx*Ny);
+  double **A, **B;
+  ALLOCATE (A, double *, Nx*Ny);
+  ALLOCATE (B, double *, Nx*Ny);
   for (i = 0; i < Nx*Ny; i++) {
-    ALLOCATE (A[i], float, Nx*Ny);
-    memset (A[i], 0, sizeof(float)*Nx*Ny);
-    ALLOCATE (B[i], float, 1);
-    memset (B[i], 0, sizeof(float));
+    ALLOCATE (A[i], double, Nx*Ny);
+    memset (A[i], 0, sizeof(double)*Nx*Ny);
+    ALLOCATE (B[i], double, 1);
+    memset (B[i], 0, sizeof(double));
   }    
 
@@ -135,24 +137,24 @@
     for (iy = 0; iy < Ny; iy++) {
       // define & init summing variables
-      float rx_rx_ry_ry = 0;
-      float rx_rx_dy_ry = 0;
-      float dx_rx_ry_ry = 0;
-      float dx_rx_dy_ry = 0;
-      float fi_rx_ry    = 0;
-      float rx_rx_py_py = 0;
-      float rx_rx_qy_py = 0;
-      float dx_rx_py_py = 0;
-      float dx_rx_qy_py = 0;
-      float fi_rx_py    = 0;
-      float px_px_ry_ry = 0;
-      float px_px_dy_ry = 0;
-      float qx_px_ry_ry = 0;
-      float qx_px_dy_ry = 0;
-      float fi_px_ry    = 0;
-      float px_px_py_py = 0;
-      float px_px_qy_py = 0;
-      float qx_px_py_py = 0;
-      float qx_px_qy_py = 0;
-      float fi_px_py    = 0;
+      double rx_rx_ry_ry = 0;
+      double rx_rx_dy_ry = 0;
+      double dx_rx_ry_ry = 0;
+      double dx_rx_dy_ry = 0;
+      double fi_rx_ry    = 0;
+      double rx_rx_py_py = 0;
+      double rx_rx_qy_py = 0;
+      double dx_rx_py_py = 0;
+      double dx_rx_qy_py = 0;
+      double fi_rx_py    = 0;
+      double px_px_ry_ry = 0;
+      double px_px_dy_ry = 0;
+      double qx_px_ry_ry = 0;
+      double qx_px_dy_ry = 0;
+      double fi_px_ry    = 0;
+      double px_px_py_py = 0;
+      double px_px_qy_py = 0;
+      double qx_px_py_py = 0;
+      double qx_px_qy_py = 0;
+      double fi_px_py    = 0;
 
       // generate the sums for the fitting matrix element I,J
@@ -167,9 +169,9 @@
 
 	// base coordinate offset for this point (x,y) relative to this map element (n,m)
-	// float dx = psImageBinningGetRuffX (map->binning, x->data.F32[i]) - (n + 0.5);
-	// float dy = psImageBinningGetRuffY (map->binning, y->data.F32[i]) - (m + 0.5);
-
-	float dx = x[i] * map->dX - ix - 0.5;
-	float dy = y[i] * map->dY - iy - 0.5;
+	// double dx = psImageBinningGetRuffX (map->binning, x->data.F32[i]) - (n + 0.5);
+	// double dy = psImageBinningGetRuffY (map->binning, y->data.F32[i]) - (m + 0.5);
+
+	double dx = x[i] * map->dX - ix - 0.5;
+	double dy = y[i] * map->dY - iy - 0.5;
 
 	// XXX do I need to do something different at the edge or not?  I want to allow
@@ -195,10 +197,10 @@
 
 	// related offset values
-	float rx = 1.0 - dx;
-	float ry = 1.0 - dy;
-	float px = 1.0 + dx;
-	float py = 1.0 + dy;
-	float qx = 0.0 - dx;
-	float qy = 0.0 - dy;
+	double rx = 1.0 - dx;
+	double ry = 1.0 - dy;
+	double px = 1.0 + dx;
+	double py = 1.0 + dy;
+	double qx = 0.0 - dx;
+	double qy = 0.0 - dy;
 
 	// sum the appropriate elements for the different quadrants
@@ -306,5 +308,5 @@
   }
 
-  if (0) {
+  if (1) {
     FILE *fd = fopen ("matrix.dat", "w");
     for (i = 0; i < Nx*Ny; i++) {
@@ -317,5 +319,5 @@
   }
 
-  if (!fgaussjordan(A, B, Nx*Ny, 1)) {
+  if (!dgaussjordan(A, B, Nx*Ny, 1)) {
     fprintf (stderr, "FAIL\n");
     exit (1);
@@ -350,2 +352,14 @@
 }
 
+int dump_map_data (float *x, float *y, float *f, int Npts, char *filename) {
+
+  FILE *fout = fopen (filename, "w");
+  
+  int i;
+  for (i = 0; i < Npts; i++) {
+    fprintf (fout, "%d %f %f %f\n", i, x[i], y[i], f[i]);
+  }
+  fclose (fout);
+  return TRUE;
+} 
+
