Index: trunk/Ohana/src/opihi/cmd.data/Makefile
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/Makefile	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/Makefile	(revision 38062)
@@ -61,4 +61,5 @@
 $(SRC)/ungridify.$(ARCH).o     \
 $(SRC)/histogram.$(ARCH).o	\
+$(SRC)/tdhistogram.$(ARCH).o	\
 $(SRC)/hermitian1d.$(ARCH).o	\
 $(SRC)/hermitian2d.$(ARCH).o	\
@@ -89,4 +90,5 @@
 $(SRC)/medacc.$(ARCH).o	\
 $(SRC)/mget.$(ARCH).o		\
+$(SRC)/mget3d.$(ARCH).o		\
 $(SRC)/minterpolate.$(ARCH).o	\
 $(SRC)/medimage.$(ARCH).o	\
@@ -130,4 +132,5 @@
 $(SRC)/imspline_apply.$(ARCH).o	\
 $(SRC)/imspline_construct.$(ARCH).o \
+$(SRC)/squash3d.$(ARCH).o	   \
 $(SRC)/imstats.$(ARCH).o	   \
 $(SRC)/style.$(ARCH).o		   \
Index: trunk/Ohana/src/opihi/cmd.data/init.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/init.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/init.c	(revision 38062)
@@ -51,4 +51,5 @@
 int ungridify        PROTO((int, char **));
 int histogram        PROTO((int, char **));
+int tdhistogram      PROTO((int, char **));
 int hermitian1d      PROTO((int, char **));
 int hermitian2d      PROTO((int, char **));
@@ -80,4 +81,5 @@
 int medacc           PROTO((int, char **));
 int mget             PROTO((int, char **));
+int mget3d           PROTO((int, char **));
 int minterp          PROTO((int, char **));
 int medimage_command PROTO((int, char **));
@@ -117,4 +119,5 @@
 int imspline_apply   PROTO((int, char **));
 int imspline_construct PROTO((int, char **));
+int squash3d         PROTO((int, char **));
 int stats            PROTO((int, char **));
 int imstats          PROTO((int, char **));
@@ -217,4 +220,5 @@
   {1, "header",       header,           "print image header"},
   {1, "histogram",    histogram,        "generate histogram from vector"},
+  {1, "tdhistogram",  tdhistogram,      "generate 2D histogram image from vector set"},
   {1, "hermitian1d",  hermitian1d,      "generate 1-D Hermitian Polynomial"},
   {1, "hermitian2d",  hermitian2d,      "generate 2-D Hermitian Polynomial"},
@@ -243,4 +247,5 @@
   {1, "medacc",       medacc,           "accumulate vector values in another vector"},
   {1, "mget",         mget,             "extract a vector from an image"},
+  {1, "mget3d",       mget3d,           "extract a vector from a 3D image"},
   {1, "imget",        mget,             "extract a vector from an image"},
   {1, "minterp",      minterp,          "interpolate image pixels"},
@@ -289,4 +294,5 @@
   {1, "imspline.apply", imspline_apply, "apply spline fit to generate an image"},
   {1, "imspline.const", imspline_construct, "create spline 2nd deriv. terms"},
+  {1, "squash3d",     squash3d,         "squash 3d buffer to 2d"},
   {1, "stats",        imstats,          "statistics on a portion of an image"},
   {1, "style",        style,            "set the style for graph plots"},
Index: trunk/Ohana/src/opihi/cmd.data/interpolate.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/interpolate.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/interpolate.c	(revision 38062)
@@ -1,4 +1,5 @@
 # include "data.h"
 
+// XXX use 'threshold' to interpolate to a value
 int interpolate (int argc, char **argv) {
 
@@ -11,6 +12,7 @@
     gprint (GP_ERR, "USAGE: interpolate Xi Yi Xo Yo\n");
     gprint (GP_ERR, "  Xi Yi - sorted reference vectors\n");
-    gprint (GP_ERR, "  Xo    - output positions\n");
-    gprint (GP_ERR, "  Yo    - output values\n");
+    gprint (GP_ERR, "  Xo    - output positions (vector)\n");
+    gprint (GP_ERR, "  Yo    - output values (vector)\n");
+    gprint (GP_ERR, "  (use 'threshold' to interpolate to a value)\n");
     return (FALSE);
   }
@@ -18,4 +20,6 @@
   if ((xin  = SelectVector (argv[1],  OLDVECTOR, TRUE)) == NULL) return (FALSE);
   if ((yin  = SelectVector (argv[2],  OLDVECTOR, TRUE)) == NULL) return (FALSE);
+
+  // target positions are a vector
   if ((xout = SelectVector (argv[3],  OLDVECTOR, TRUE)) == NULL) return (FALSE);
   if ((yout = SelectVector (argv[4],  ANYVECTOR, TRUE)) == NULL) return (FALSE);
Index: trunk/Ohana/src/opihi/cmd.data/mcreate.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/mcreate.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/mcreate.c	(revision 38062)
@@ -3,20 +3,32 @@
 int mcreate (int argc, char **argv) {
   
-  int Nx, Ny;
+  int N;
   Buffer *buf;
 
+  int Nz = 0;
+  if ((N = get_argument (argc, argv, "-nz"))) {
+    remove_argument (N, &argc, argv);
+    Nz = atoi (argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+
   if (argc != 4) {
-    gprint (GP_ERR, "USAGE: mcreate <buffer> Nx Ny\n");
+    gprint (GP_ERR, "USAGE: mcreate <buffer> Nx Ny [-nz Nz]\n");
     return (FALSE);
   }
 
   if ((buf = SelectBuffer (argv[1], ANYBUFFER, TRUE)) == NULL) return (FALSE);
-  Nx = atof (argv[2]);
-  Ny = atof (argv[3]);
+  int Nx = atof (argv[2]);
+  int Ny = atof (argv[3]);
 
   /* I should encapsulate this in a create_default_buffer */
   gfits_free_matrix (&buf[0].matrix);
   gfits_free_header (&buf[0].header);
-  CreateBuffer (buf, Nx, Ny, -32, 1.0, 0.0);
+
+  if (Nz) {
+    CreateBuffer3D (buf, Nx, Ny, Nz, -32, 1.0, 0.0);
+  } else {
+    CreateBuffer (buf, Nx, Ny, -32, 1.0, 0.0);
+  }
   return (TRUE);
 }
Index: trunk/Ohana/src/opihi/cmd.data/mget3d.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/mget3d.c	(revision 38062)
+++ trunk/Ohana/src/opihi/cmd.data/mget3d.c	(revision 38062)
@@ -0,0 +1,46 @@
+# include "data.h"
+
+int mget3d (int argc, char **argv) {
+  
+  int i;
+  Buffer *buf;
+  Vector *vec;
+
+  if (argc != 5) {
+    gprint (GP_ERR, "USAGE: mget <buffer> <vector> x y\n");
+    return (FALSE);
+  }
+
+  int x = atoi(argv[3]);
+  int y = atoi(argv[4]);
+
+  if ((buf = SelectBuffer (argv[1], OLDBUFFER, TRUE)) == NULL) return (FALSE);
+  if (buf[0].matrix.Naxes < 3) {
+    gprint (GP_ERR, "buffer is not 3D\n");
+    return FALSE;
+  }
+
+  int Nx = buf[0].matrix.Naxis[0];
+  int Ny = buf[0].matrix.Naxis[1];
+  int Nz = buf[0].matrix.Naxis[2];
+
+  int invalid = FALSE;
+  invalid = invalid || (x < 0);
+  invalid = invalid || (x >= Nx);
+  invalid = invalid || (y < 0);
+  invalid = invalid || (y >= Ny);
+  if (invalid) {
+    gprint (GP_ERR, "selection (%d,%d) out of range\n", x, y);
+    return (FALSE);
+  }
+
+  if ((vec = SelectVector (argv[2], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+
+  ResetVector (vec, OPIHI_FLT, Nz);
+  float *in  = (float *) buf[0].matrix.buffer + x + y*Nx;
+  opihi_flt *out = vec[0].elements.Flt;
+  for (i = 0; i < Nz; i++, in += Nx*Ny, out++) {
+    *out = *in;
+  }
+  return (TRUE);
+}
Index: trunk/Ohana/src/opihi/cmd.data/rd.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/rd.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/rd.c	(revision 38062)
@@ -4,11 +4,10 @@
 int rd (int argc, char **argv) {
   
-  int i, N, status, plane, Nplane, extend, Nextend, Nskip, JustHead, blank;
-  int ccdsel, done, Nword, IsCompressed;
-  char region[512], *ccdid, *filename;
-  FILE *f;
+  int i, N, Nskip, blank;
+  int done, Nword;
+  char region[512];
   Buffer *buf;
 
-  JustHead = FALSE;
+  int JustHead = FALSE;
   if ((N = get_argument (argc, argv, "-head"))) {
     remove_argument (N, &argc, argv);
@@ -16,5 +15,5 @@
   }
 
-  plane = 1;
+  int plane = -1;
   if ((N = get_argument (argc, argv, "-plane"))) {
     remove_argument (N, &argc, argv);
@@ -23,6 +22,6 @@
   }
 
-  extend = FALSE;
-  Nextend = -1;
+  int extend = FALSE;
+  int Nextend = -1;
   if ((N = get_argument (argc, argv, "-x"))) {
     remove_argument (N, &argc, argv);
@@ -32,6 +31,6 @@
   }
 
-  ccdsel = FALSE;
-  ccdid = (char *) NULL;
+  int ccdsel = FALSE;
+  char *ccdid = NULL;
   if ((N = get_argument (argc, argv, "-n"))) {
     remove_argument (N, &argc, argv);
@@ -52,5 +51,5 @@
 
   /* test if file exists */
-  f = fopen (argv[2], "r");
+  FILE *f = fopen (argv[2], "r");
   if (f == (FILE *) NULL) {
     gprint (GP_ERR, "file %s not found\n", argv[2]);
@@ -67,10 +66,10 @@
 
   /* save file name */
-  filename = filebasename (argv[2]);
+  char *filename = filebasename (argv[2]);
   strcpy (buf[0].file, filename);
   free (filename);
 
-  status = FALSE;
-  IsCompressed = FALSE;
+  int status = FALSE;
+  int IsCompressed = FALSE;
 
   /*** advance to the correct FITS extension ***/
@@ -168,15 +167,24 @@
 
   /* check for valid plane */
-  Nplane = buf[0].header.Naxis[2];
-  if (Nplane == 0) Nplane = 1;
-  if (plane > Nplane) {
-    gprint (GP_ERR, "-plane is too large: %d total planes\n", Nplane);
-    DeleteBuffer (buf);
-    fclose (f);
-    return (FALSE);
+  int Nz = buf[0].header.Naxis[2];
+  if (plane >= 0) {
+    // we are requesting a specific plane (-1 : all data)
+    int tooFar = Nz ? (plane >= Nz) : (plane > Nz);
+    if (tooFar) {
+      gprint (GP_ERR, "-plane is too large: %d total planes\n", Nz);
+      DeleteBuffer (buf);
+      fclose (f);
+      return (FALSE);
+    }
   }
 
   /* load matrix data */
   if (IsCompressed) {
+    if (plane > -1) {
+      gprint (GP_ERR, "-plane incompatible with compressed image\n");
+      DeleteBuffer (buf);
+      fclose (f);
+      return (FALSE);
+    }
     FTable ftable;
     Header theader;
@@ -191,6 +199,15 @@
     // XXX this currently does not work for a cube (we get a cube back, not a specific plane)
   } else {
-    sprintf (region, "-1 -1 -1 -1 %d %d", (plane - 1), plane);
-    status = gfits_fread_matrix_segment (f, &buf[0].matrix, &buf[0].header, region);
+    if (plane > -1) {
+      // read a single plane into a 2D image
+      sprintf (region, "-1 -1 -1 -1 %d %d", (plane - 1), plane);
+      status = gfits_fread_matrix_segment (f, &buf[0].matrix, &buf[0].header, region);
+      buf[0].header.Naxis[2] = 0;
+      buf[0].header.Naxes = 2;
+      gfits_modify (&buf[0].header, "NAXIS", "%d", 1, 2);
+      gfits_modify (&buf[0].header, "NAXIS3", "%d", 1, 0);
+    } else {
+      status = gfits_fread_matrix (f, &buf[0].matrix, &buf[0].header);
+    }
   }
   fclose (f);
@@ -200,12 +217,4 @@
     DeleteBuffer (buf);
     return (FALSE);
-  }
-
-  /* adjust buffer to represent 2D data */
-  if (Nplane > 1) {
-    buf[0].header.Naxis[2] = 0;
-    buf[0].header.Naxes = 2;
-    gfits_modify (&buf[0].header, "NAXIS", "%d", 1, 2);
-    gfits_delete (&buf[0].header, "NAXIS3", 1);
   }
 
@@ -215,4 +224,6 @@
     buf[0].header.Naxis[1] = 1;
     buf[0].matrix.Naxis[1] = 1;
+    gfits_modify (&buf[0].header, "NAXIS", "%d", 1, 2);
+    gfits_modify (&buf[0].header, "NAXIS2", "%d", 1, 1);
   }    
 
Index: trunk/Ohana/src/opihi/cmd.data/read_vectors.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/read_vectors.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/read_vectors.c	(revision 38062)
@@ -372,4 +372,10 @@
   }
 
+  int VERBOSE = FALSE;
+  if ((N = get_argument (argc, argv, "-v"))) {
+    remove_argument (N, &argc, argv);
+    VERBOSE = TRUE;
+  }
+
   int start = 0;
   int Nrows = -1; // -1 : read entire table
@@ -399,5 +405,6 @@
   // }
 
-  if (argc < 2) ESCAPE ("USAGE: read -fits extension [-extnum] [-keyword key] name name ...\n");
+  if ((argc < 2) && !getSizes) ESCAPE ("USAGE: read -fits extension [-extnum] [-keyword key] name name ...\n");
+  if ((argc != 1) && getSizes) ESCAPE ("USAGE: read -sizes -fits extension [-extnum] (does not read data values)\n");
 
   if (f == NULL) ESCAPE ("file not found\n");
@@ -436,5 +443,9 @@
     if (start >= header.Naxis[1]) ESCAPE ("invalid range: start >= Ny (%d)\n", header.Naxis[1]);
     if (Nrows < 0) ESCAPE ("invalid range: Nrows < 0\n");
-    if (start + Nrows > header.Naxis[1]) ESCAPE ("invalid range: start + Nrows > Ny (%d)\n", header.Naxis[1]);
+
+    // just a warning:
+    if (start + Nrows > header.Naxis[1]) {
+      if (VERBOSE) gprint (GP_ERR, "NOTE: reading last block will return only %d rows\n", header.Naxis[1] - start);
+    }
 
     // Ny = 100, start = 0, Nrows = -1 -> Nrows => 100
@@ -472,4 +483,11 @@
       }
 
+      if (getSizes) {
+	read_table_sizes (&header);
+	if (CCDKeyword != NULL) free (CCDKeyword); 
+	gfits_free_header (&header); 
+	return TRUE;
+      }
+
       if (Nrows == -1) {
 	Nrows = header.Naxis[1] - start;
@@ -478,5 +496,9 @@
       if (start >= header.Naxis[1]) ESCAPE ("invalid range: start >= Ny (%d)\n", header.Naxis[1]);
       if (Nrows < 0) ESCAPE ("invalid range: Nrows < 0\n");
-      if (start + Nrows > header.Naxis[1]) ESCAPE ("invalid range: start + Nrows > Ny (%d)\n", header.Naxis[1]);
+
+      // just a warning:
+      if (start + Nrows > header.Naxis[1]) {
+	if (VERBOSE) gprint (GP_ERR, "NOTE: reading last block will return only %d rows\n", header.Naxis[1] - start);
+      }
 
       if (!gfits_fread_ftable_range (f, padIfShort, &table, start, Nrows)) ESCAPE ("error reading table for extension %d\n", Nextend);
@@ -603,5 +625,5 @@
 
   set_int_variable ("table:Nx", header->Naxis[0]);
-  set_int_variable ("table:Nx", header->Naxis[0]);
+  set_int_variable ("table:Ny", header->Naxis[1]);
   set_int_variable ("table:Nfields", Nfields);
 
Index: trunk/Ohana/src/opihi/cmd.data/squash3d.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/squash3d.c	(revision 38062)
+++ trunk/Ohana/src/opihi/cmd.data/squash3d.c	(revision 38062)
@@ -0,0 +1,113 @@
+# include "data.h"
+
+enum {SQUASH_NONE, SQUASH_X, SQUASH_Y, SQUASH_Z};
+
+int squash3d (int argc, char **argv) {
+  
+  int ox, oy, oz;
+  Buffer *src;
+  Buffer *tgt;
+
+  if (argc != 4) {
+    gprint (GP_ERR, "USAGE: squash3d <3d> <2d> (dir)\n");
+    gprint (GP_ERR, " dir: x -> squash in x-dir, y -> squash in y-dir, z -> squash in z-dir\n");
+    return (FALSE);
+  }
+
+  int dir = SQUASH_NONE;
+  if (!strcasecmp (argv[3], "x")) dir = SQUASH_X;
+  if (!strcasecmp (argv[3], "y")) dir = SQUASH_Y;
+  if (!strcasecmp (argv[3], "z")) dir = SQUASH_Z;
+  if (!dir) {
+    gprint (GP_ERR, "invalid direction %s\n", argv[3]);
+    return FALSE;
+  }
+
+  if ((src = SelectBuffer (argv[1], OLDBUFFER, TRUE)) == NULL) return (FALSE);
+  if (src[0].matrix.Naxes < 3) {
+    gprint (GP_ERR, "buffer is not 3D\n");
+    return FALSE;
+  }
+
+  float *iBuf  = (float *) src[0].matrix.buffer;
+
+  if ((tgt = SelectBuffer (argv[2], ANYBUFFER, TRUE)) == NULL) return (FALSE);
+  gfits_free_matrix (&tgt[0].matrix);
+  gfits_free_header (&tgt[0].header);
+
+  switch (dir) {
+    case SQUASH_X:
+      {     
+	int iNx = src[0].matrix.Naxis[0];
+	int iNy = src[0].matrix.Naxis[1];
+	int iNz = src[0].matrix.Naxis[2];
+
+	// output is Nz,Ny
+	int oNx = iNz;
+	int oNy = iNy;
+	CreateBuffer (tgt, oNx, oNy, -32, 0.0, 1.0);
+	float *oBuf  = (float *) tgt[0].matrix.buffer;
+
+	for (ox = 0; ox < oNx; ox ++) {
+	  for (oy = 0; oy < oNy; oy ++) {
+	    float val = 0.0;
+	    for (oz = 0; oz < iNx; oz ++) { // src x is tgt z
+	      val += iBuf[oz + oy*iNx + ox*iNx*iNy];
+	    }
+	    oBuf[ox + oy*oNx] = val;
+	  }
+	}
+      }
+      break;
+      
+    case SQUASH_Y:
+      {     
+	int iNx = src[0].matrix.Naxis[0];
+	int iNy = src[0].matrix.Naxis[1];
+	int iNz = src[0].matrix.Naxis[2];
+
+	// output is Nx,Nz
+	int oNx = iNx;
+	int oNy = iNz;
+	CreateBuffer (tgt, oNx, oNy, -32, 0.0, 1.0);
+	float *oBuf  = (float *) tgt[0].matrix.buffer;
+
+	for (ox = 0; ox < oNx; ox ++) {
+	  for (oy = 0; oy < oNy; oy ++) { // src z is tgt y
+	    float val = 0.0;
+	    for (oz = 0; oz < iNy; oz ++) {
+	      val += iBuf[ox + oz*iNx + oy*iNx*iNy];
+	    }
+	    oBuf[ox + oy*oNx] = val;
+	  }
+	}
+      }
+      break;
+      
+    case SQUASH_Z:
+      {     
+	int iNx = src[0].matrix.Naxis[0];
+	int iNy = src[0].matrix.Naxis[1];
+	int iNz = src[0].matrix.Naxis[2];
+
+	// output is Nx,Ny
+	int oNx = iNx;
+	int oNy = iNy;
+	CreateBuffer (tgt, oNx, oNy, -32, 0.0, 1.0);
+	float *oBuf  = (float *) tgt[0].matrix.buffer;
+
+	for (ox = 0; ox < oNx; ox ++) {
+	  for (oy = 0; oy < oNy; oy ++) {
+	    float val = 0.0;
+	    for (oz = 0; oz < iNz; oz ++) {
+	      val += iBuf[ox + oy*iNx + oz*iNx*iNy];
+	    }
+	    oBuf[ox + oy*oNx] = val;
+	  }
+	}
+      }
+      break;
+  }
+
+  return (TRUE);
+}
Index: trunk/Ohana/src/opihi/cmd.data/tdhistogram.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/tdhistogram.c	(revision 38062)
+++ trunk/Ohana/src/opihi/cmd.data/tdhistogram.c	(revision 38062)
@@ -0,0 +1,176 @@
+# include "data.h"
+
+# define CHECKVAL(ARG) if (!isfinite(ARG)) { gprint (GP_ERR, "illegal value for %s: %f\n", #ARG, ARG); return (FALSE); }
+
+int tdhistogram (int argc, char **argv) {
+
+  int i, Nx, Ny, Nz, N;
+  float *val;
+  Buffer *bf;
+  Vector *vx, *vy, *vz, *range;
+  opihi_flt *x, *y, *z;
+
+  int reuse = FALSE;
+  if ((N = get_argument (argc, argv, "-reuse"))) {
+    remove_argument (N, &argc, argv);
+    reuse = TRUE;
+  }
+
+  range = NULL;
+  if ((N = get_argument (argc, argv, "-range"))) {
+    remove_argument (N, &argc, argv);
+    if ((range = SelectVector (argv[N], ANYVECTOR, TRUE)) == NULL) return (FALSE);    
+    remove_argument (N, &argc, argv);
+  }
+
+  double dx = NAN;
+  if ((N = get_argument (argc, argv, "-dx"))) {
+    remove_argument (N, &argc, argv);
+    dx = atof(argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+  double dy = NAN;
+  if ((N = get_argument (argc, argv, "-dy"))) {
+    remove_argument (N, &argc, argv);
+    dy = atof(argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+  double dz = NAN;
+  if ((N = get_argument (argc, argv, "-dz"))) {
+    remove_argument (N, &argc, argv);
+    dz = atof(argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+
+  int valid = (!reuse && (argc == 11)) || (reuse && (argc == 5));
+  if (!valid) {
+    gprint (GP_ERR, "USAGE: tdhistogram buffer x y z (Xmin) (Xmax) (Ymin) (Ymax) (Zmin) (Zmax) [-dx dx] [-dy dy] [-dz dz]\n");
+    gprint (GP_ERR, "   OR: tdhistogram buffer x y z -reuse\n");
+    gprint (GP_ERR, " output buffer is 3D\n");
+    return (FALSE);
+  }
+  
+  if ((bf = SelectBuffer (argv[1], ANYBUFFER, TRUE)) == NULL) return (FALSE);
+  if ((vx = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((vy = SelectVector (argv[3], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((vz = SelectVector (argv[4], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+
+  if (vx[0].Nelements != vy[0].Nelements) return (FALSE);
+  if (vx[0].Nelements != vz[0].Nelements) return (FALSE);
+
+  double Xmin, Xmax, Ymin, Ymax, Zmin, Zmax;
+
+  if (!reuse) {
+    Xmin = atof(argv[5]);
+    Xmax = atof(argv[6]);
+    if (!isfinite(dx)) {
+      Nx = 100;
+      dx = (Xmax - Xmin) / (Nx - 1);
+    } else {
+      Nx = (Xmax - Xmin) / dx + 1;
+    }
+    if (dx < 0) {
+      gprint (GP_ERR, "invalid value for delta: %f\n", dx);
+      return (FALSE);
+    }
+    Ymin = atof(argv[7]);
+    Ymax = atof(argv[8]);
+    if (!isfinite(dy)) {
+      Ny = 100;
+      dy = (Ymax - Ymin) / (Ny - 1);
+    } else {
+      Ny = (Ymax - Ymin) / dy + 1;
+    }
+    if (dy < 0) {
+      gprint (GP_ERR, "invalid value for delta: %f\n", dy);
+      return (FALSE);
+    }
+    Zmin = atof(argv[9]);
+    Zmax = atof(argv[10]);
+    if (!isfinite(dz)) {
+      Nz = 100;
+      dz = (Zmax - Zmin) / (Nz - 1);
+    } else {
+      Nz = (Zmax - Zmin) / dz + 1;
+    }
+    if (dz < 0) {
+      gprint (GP_ERR, "invalid value for delta: %f\n", dz);
+      return (FALSE);
+    }
+
+    if (Nz > 1000) {
+      gprint (GP_ERR, "warning: delta of %f will result in %d histogram bins\n", dz, Nz);
+      return (FALSE);
+    }
+  } else {
+    gfits_scan (&bf[0].header, "XMIN", "%lf", 1, &Xmin);
+    gfits_scan (&bf[0].header, "XMAX", "%lf", 1, &Xmax);
+    gfits_scan (&bf[0].header, "XDEL", "%lf", 1, &dx);
+    gfits_scan (&bf[0].header, "YMIN", "%lf", 1, &Ymin);
+    gfits_scan (&bf[0].header, "YMAX", "%lf", 1, &Ymax);
+    gfits_scan (&bf[0].header, "YDEL", "%lf", 1, &dy);
+    gfits_scan (&bf[0].header, "ZMIN", "%lf", 1, &Zmin);
+    gfits_scan (&bf[0].header, "ZMAX", "%lf", 1, &Zmax);
+    gfits_scan (&bf[0].header, "ZDEL", "%lf", 1, &dz);
+    Nx = bf[0].header.Naxis[0];
+    Ny = bf[0].header.Naxis[1];
+    Nz = bf[0].header.Naxis[2];
+  }
+
+  REQUIRE_VECTOR_FLT (vx, FALSE); 
+  REQUIRE_VECTOR_FLT (vy, FALSE); 
+  REQUIRE_VECTOR_FLT (vz, FALSE); 
+
+  CHECKVAL(Xmin);
+  CHECKVAL(Xmax);
+  CHECKVAL(dx);
+
+  CHECKVAL(Ymin);
+  CHECKVAL(Ymax);
+  CHECKVAL(dy);
+
+  if (range) {
+    ResetVector (range, OPIHI_FLT, Nz);
+    for (i = 0; i < range[0].Nelements; i++) {
+      range[0].elements.Flt[i] = Zmin + i*dz;
+    }
+  }
+
+  if (!reuse) {
+    gfits_free_matrix (&bf[0].matrix);
+    gfits_free_header (&bf[0].header);
+    CreateBuffer3D (bf, Nx, Ny, Nz, -32, 0.0, 1.0);
+    strcpy (bf[0].file, "(empty)");
+
+    gfits_modify (&bf[0].header, "XMIN", "%lf", 1, Xmin);
+    gfits_modify (&bf[0].header, "XMAX", "%lf", 1, Xmax);
+    gfits_modify (&bf[0].header, "XDEL", "%lf", 1, dx);
+    gfits_modify (&bf[0].header, "YMIN", "%lf", 1, Ymin);
+    gfits_modify (&bf[0].header, "YMAX", "%lf", 1, Ymax);
+    gfits_modify (&bf[0].header, "YDEL", "%lf", 1, dy);
+    gfits_modify (&bf[0].header, "ZMIN", "%lf", 1, Zmin);
+    gfits_modify (&bf[0].header, "ZMAX", "%lf", 1, Zmax);
+    gfits_modify (&bf[0].header, "ZDEL", "%lf", 1, dz);
+  }
+  
+  x = vx[0].elements.Flt;
+  y = vy[0].elements.Flt;
+  z = vz[0].elements.Flt;
+
+  val = (float *) bf[0].matrix.buffer;
+
+  for (i = 0; i < vx[0].Nelements; i++, x++, y++, z++) {
+    int Xb = (*x - Xmin) / dx;
+    int Yb = (*y - Ymin) / dy;
+    int Zb = (*z - Zmin) / dz;
+
+    if (Xb >= Nx) continue;
+    if (Yb >= Ny) continue;
+    if (Zb >= Nz) continue;
+    if (Xb < 0) continue;
+    if (Yb < 0) continue;
+    if (Zb < 0) continue;
+    val[Xb + Yb*Nx + Zb*Nx*Ny] ++;
+  }
+  return (TRUE);
+}
Index: trunk/Ohana/src/opihi/cmd.data/tv.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/tv.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/tv.c	(revision 38062)
@@ -28,4 +28,15 @@
   }
 
+  int plane = 0;
+  if ((N = get_argument (argc, argv, "-plane"))) {
+    remove_argument (N, &argc, argv);
+    plane = atoi (argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+  if (plane < 0) {
+    gprint (GP_ERR, " ERROR: -plane (plane) : cannot be negative\n");
+    return (FALSE);
+  }
+
   /* shell exits on pipe close, FIX */
   if ((N = get_argument (argc, argv, "-kill"))) {
@@ -55,7 +66,16 @@
   GetCoords (&coords, &buf[0].header);
   
-  image.data1d = (float *) buf[0].matrix.buffer;
   image.Nx = buf[0].matrix.Naxis[0];
   image.Ny = buf[0].matrix.Naxis[1];
+
+  int tooBig = buf[0].matrix.Naxis[2] ? (plane >= buf[0].matrix.Naxis[2]) : plane > 0;
+  if (tooBig) {
+    gprint (GP_ERR, " ERROR: -plane (plane) : out of bounds (%d vs %d)\n", plane, buf[0].matrix.Naxis[2]);
+    return (FALSE);
+  }
+  int Npix2D = image.Nx * image.Ny;
+
+  float *imdata = (float *) buf[0].matrix.buffer;
+  image.data1d = &imdata[plane*Npix2D];
 
   // send only the root of the file, not the full path
Index: trunk/Ohana/src/opihi/cmd.data/wd.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/wd.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/wd.c	(revision 38062)
@@ -81,14 +81,9 @@
   memcpy (temp_header.buffer, buf[0].header.buffer, temp_header.datasize);
 
-  if (temp_header.Naxes) {
-    // the inBlank value probably does not matter: temp_matrix is float, so nan is used
-    gfits_convert_format (&temp_header, &temp_matrix, outBitpix, outScale, outZero, 0xffff, outUnsign);
-  } else {
-    gfits_modify (&temp_header, "BITPIX", "%d", 1, outBitpix);
-    gfits_modify (&temp_header, "BSCALE", "%lf", 1, outScale);
-    gfits_modify (&temp_header, "BZERO",  "%lf", 1, outZero);
-    gfits_modify_alt (&temp_header, "UNSIGN", "%t", 1, outUnsign);
-  }
+  gfits_convert_format (&temp_header, &temp_matrix, outBitpix, outScale, outZero, 0xffff, outUnsign);
 
+  // Extend puts the output matrix in the first available non-PHU slot (ie, the last one)
+  // it updates NEXTEND and set EXTEND to TRUE, and modifies the PHU header (neither should happen)
+  // if those keywords do not exist...
   if (Extend) {
     Header Xhead;
Index: trunk/Ohana/src/opihi/cmd.data/zplot.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/zplot.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/zplot.c	(revision 38062)
@@ -4,11 +4,34 @@
   
   char *outname = NULL;
-  int i, kapa, valid, size;
+  int i, N, kapa, valid, size;
   opihi_flt *out;
   double min, range;
   Graphdata graphmode;
-  Vector *xvec, *yvec, *zvec, Zvec;
+  Vector *xvec, *yvec, *zvec, *dxmvec, *dxpvec, *dymvec, *dypvec, Zvec;
 
   if (!style_args (&graphmode, &argc, argv, &kapa)) return (FALSE);
+
+  /* decide on error bars */
+  dxmvec = dxpvec = dymvec = dypvec = NULL;
+  if ((N = get_argument (argc, argv, "-dx"))) {
+    remove_argument (N, &argc, argv);
+    if ((dxmvec = SelectVector (argv[N], OLDVECTOR, TRUE)) == NULL) return (FALSE);    
+    remove_argument (N, &argc, argv);
+  }
+  if ((N = get_argument (argc, argv, "+dx"))) {
+    remove_argument (N, &argc, argv);
+    if ((dxpvec = SelectVector (argv[N], OLDVECTOR, TRUE)) == NULL) return (FALSE);    
+    remove_argument (N, &argc, argv);
+  }
+  if ((N = get_argument (argc, argv, "-dy"))) {
+    remove_argument (N, &argc, argv);
+    if ((dymvec = SelectVector (argv[N], OLDVECTOR, TRUE)) == NULL) return (FALSE);    
+    remove_argument (N, &argc, argv);
+  }
+  if ((N = get_argument (argc, argv, "+dy"))) {
+    remove_argument (N, &argc, argv);
+    if ((dypvec = SelectVector (argv[N], OLDVECTOR, TRUE)) == NULL) return (FALSE);    
+    remove_argument (N, &argc, argv);
+  }
 
   valid  = (argc == 6);
@@ -40,4 +63,13 @@
   }
 
+  /* set errorbar mode (these are NOT sticky) */
+  graphmode.etype = 0;
+  if ((dymvec != NULL) && (dypvec == NULL)) dypvec = dymvec;
+  if ((dypvec != NULL) && (dymvec == NULL)) dymvec = dypvec;
+  if ((dypvec != NULL) || (dymvec != NULL)) graphmode.etype |= 0x01;
+  if ((dxmvec != NULL) && (dxpvec == NULL)) dxpvec = dxmvec;
+  if ((dxpvec != NULL) && (dxmvec == NULL)) dxmvec = dxpvec;
+  if ((dxpvec != NULL) || (dxmvec != NULL)) graphmode.etype |= 0x02;
+  
   /* find vectors */
   if ((xvec = SelectVector (argv[1], OLDVECTOR, TRUE)) == NULL) return (FALSE);
@@ -56,7 +88,13 @@
     return (FALSE);
   }
+  if (dypvec && (dypvec->Nelements != xvec->Nelements)) goto mismatch;
+  if (dymvec && (dymvec->Nelements != xvec->Nelements)) goto mismatch;
+  if (dxpvec && (dxpvec->Nelements != xvec->Nelements)) goto mismatch;
+  if (dxmvec && (dxmvec->Nelements != xvec->Nelements)) goto mismatch;
+
   SetVector (&Zvec, OPIHI_FLT, zvec[0].Nelements);
   out = Zvec.elements.Flt;
  
+  // note actual size is 3.3x plot -sz sizes. (DrawObjects.c:399)
   if (zvec[0].type == OPIHI_FLT) {
     opihi_flt *in = zvec[0].elements.Flt;
@@ -81,13 +119,31 @@
   graphmode.style = 2; /* plot points */
   graphmode.size = -1; /* point size determined by Zvec */
-  graphmode.etype = 0; /* no errorbars */
   PlotVectorTriplet (kapa, xvec, yvec, &Zvec, mask, &graphmode);
+  if (graphmode.etype & 0x01) {
+    PlotVectorSingle (kapa, dymvec, mask, "dym");
+    PlotVectorSingle (kapa, dypvec, mask, "dyp");
+  }
+  if (graphmode.etype & 0x02) {
+    PlotVectorSingle (kapa, dxmvec, mask, "dxm");
+    PlotVectorSingle (kapa, dxpvec, mask, "dxp");
+  }
 
   free (Zvec.elements.Ptr);
   if (mask) free (mask);
-  DeleteNamedVector (outname);
+
+  if (outname) {
+    DeleteNamedVector (outname);
+    free (outname);
+  }
 
   return (TRUE);
 
+mismatch:
+  gprint (GP_ERR, "error and data vector lengths are mismatched\n");
+  if (outname) {
+    DeleteNamedVector (outname);
+    free (outname);
+  }
+  return (FALSE);
 }