Index: /trunk/psLib/src/imageops/psImageUnbin.c
===================================================================
--- /trunk/psLib/src/imageops/psImageUnbin.c	(revision 11794)
+++ /trunk/psLib/src/imageops/psImageUnbin.c	(revision 11795)
@@ -11,6 +11,8 @@
 
 // interpolate from model to background (~bresenham linear interpolation)
-// XXX this code skips the initial pixels?
-// XXX this code may be off by 0.5,0.5 pixels
+// DX, DY are the binning factor
+// dx, dy is the distance is high-res pixels to the 0,0 corner of the first
+// binned pixel. 
+// XXX check that this is still consistent with psphotImageMedian...
 psImage *psImageUnbin(psImage *out, const psImage *in, int DX, int DY, int dx, int dy)
 {
@@ -30,5 +32,5 @@
     psF32 **vOut = out->data.F32;
 
-    // loop over all pixels with
+    // loop over all input pixels excluding the last
     for (int Iy = 0; Iy < ny-1; Iy ++) {
         for (int Ix = 0; Ix < nx-1; Ix ++) {
@@ -126,13 +128,20 @@
         }
     }
+    // return out;
 
     // the four corners
     {
         float V;
+	// center of last pixel 
+	int Xs = PS_MAX (0, PS_MIN (Nx, ( 0 + 0.5)*DX - dx));
+	int Xe = PS_MAX (0, PS_MIN (Nx, (nx - 0.5)*DX - dx));
+	int Ys = PS_MAX (0, PS_MIN (Ny, ( 0 + 0.5)*DY - dy));
+	int Ye = PS_MAX (0, PS_MIN (Ny, (ny - 0.5)*DY - dy));
+
         // 0,0
         V = vIn[0][0];
-        for (int iy = 0; iy < DY - dy; iy++)
-        {
-            for (int ix = 0; ix < DX - dx; ix++) {
+        for (int iy = 0; iy < Ys; iy++)
+        {
+            for (int ix = 0; ix < Xs; ix++) {
                 vOut[iy][ix] = V;
             }
@@ -140,7 +149,7 @@
         // Nx,0
         V = vIn[0][nx-1];
-        for (int iy = 0; iy < DY - dy; iy++)
-        {
-            for (int ix = nx*DX - dx; ix < Nx; ix++) {
+        for (int iy = 0; iy < Ys; iy++)
+        {
+            for (int ix = Xe; ix < Nx; ix++) {
                 vOut[iy][ix] = V;
             }
@@ -148,7 +157,7 @@
         // 0,Ny
         V = vIn[ny-1][0];
-        for (int iy = ny*DY - dy; iy < Ny; iy++)
-        {
-            for (int ix = 0; ix < DX - dx; ix++) {
+        for (int iy = Ye; iy < Ny; iy++)
+        {
+            for (int ix = 0; ix < Xs; ix++) {
                 vOut[iy][ix] = V;
             }
@@ -156,7 +165,7 @@
         // Nx,Ny
         V = vIn[ny-1][nx-1];
-        for (int iy = ny*DY - dy; iy < Ny; iy++)
-        {
-            for (int ix = nx*DX - dx; ix < Nx; ix++) {
+        for (int iy = Ye; iy < Ny; iy++)
+        {
+            for (int ix = Xe; ix < Nx; ix++) {
                 vOut[iy][ix] = V;
             }
