Index: trunk/Ohana/src/opihi/cmd.data/Makefile
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/Makefile	(revision 41715)
+++ trunk/Ohana/src/opihi/cmd.data/Makefile	(revision 41891)
@@ -28,4 +28,6 @@
 $(SRC)/bisection.$(ARCH).o		\
 $(SRC)/cast.$(ARCH).o		\
+$(SRC)/chebyshev.$(ARCH).o		\
+$(SRC)/chebyshev_commands.$(ARCH).o	\
 $(SRC)/center.$(ARCH).o	\
 $(SRC)/clear.$(ARCH).o		\
Index: trunk/Ohana/src/opihi/cmd.data/applyfit2d.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/applyfit2d.c	(revision 41715)
+++ trunk/Ohana/src/opihi/cmd.data/applyfit2d.c	(revision 41891)
@@ -116,2 +116,117 @@
 
 
+int applyfit2d_full (int argc, char **argv) {
+  
+  int i, j, n, order;
+  char *c, name[64];
+  double **C, X, Y;
+  opihi_flt *z;
+  Vector *xvec, *yvec, *zvec;
+
+  if (argc != 4) {
+    gprint (GP_ERR, "USAGE: applyfit2d_full x y z\n");
+    return (FALSE);
+  }
+
+  c = get_variable ("Cnn");
+  if (c == NULL) {
+    gprint (GP_ERR, "no fit available\n");
+    return (FALSE);
+  }
+  order = atof (c);
+  free (c);
+
+  if ((xvec = SelectVector (argv[1], OLDVECTOR, TRUE)) == NULL) return (FALSE);    
+  if ((yvec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);    
+  if ((zvec = SelectVector (argv[3], ANYVECTOR, TRUE)) == NULL) return (FALSE);    
+
+  ALLOCATE (C, double *, order+1);
+  for (i = 0; i <= order; i++) {
+    ALLOCATE (C[i], double, order+1);
+    for (j = 0; j <= order; j++) {
+      sprintf (name, "CX%dY%d", i, j);
+      c = get_variable (name);
+      if (c == NULL) {
+	gprint (GP_ERR, "missing fit term %d,%d\n", i, j);
+	for (j = 0; j < i; j++) free (C[j]);
+	free (C);
+	return (FALSE);
+      }
+      C[i][j] = atof (c);
+      free (c);
+    }
+  }
+
+  ResetVector (zvec, OPIHI_FLT, xvec[0].Nelements);
+  bzero (zvec[0].elements.Flt, sizeof(opihi_flt)*zvec[0].Nelements);
+  z = zvec[0].elements.Flt;
+
+  // we have four cases (x == flt or int) and (y == flt or int)
+  if ((xvec[0].type == OPIHI_FLT) && (yvec[0].type == OPIHI_FLT)) {
+    opihi_flt *x = xvec[0].elements.Flt;
+    opihi_flt *y = yvec[0].elements.Flt;
+    for (n = 0; n < xvec[0].Nelements; n++, x++, y++, z++) {
+      Y = X = 1;
+      for (j = 0; j <= order; j++) {
+	X = Y;
+	for (i = 0; i <= order; i++) {
+	  *z += C[i][j]*X;
+	  X = X * (*x);
+	}
+	Y = Y * (*y);
+      }
+    }
+  }
+  if ((xvec[0].type == OPIHI_FLT) && (yvec[0].type == OPIHI_INT)) {
+    opihi_flt *x = xvec[0].elements.Flt;
+    opihi_int *y = yvec[0].elements.Int;
+    for (n = 0; n < xvec[0].Nelements; n++, x++, y++, z++) {
+      Y = X = 1;
+      for (j = 0; j <= order; j++) {
+	X = Y;
+	for (i = 0; i <= order; i++) {
+	  *z += C[i][j]*X;
+	  X = X * (*x);
+	}
+	Y = Y * (*y);
+      }
+    }
+  }
+  if ((xvec[0].type == OPIHI_INT) && (yvec[0].type == OPIHI_FLT)) {
+    opihi_int *x = xvec[0].elements.Int;
+    opihi_flt *y = yvec[0].elements.Flt;
+    for (n = 0; n < xvec[0].Nelements; n++, x++, y++, z++) {
+      Y = X = 1;
+      for (j = 0; j <= order; j++) {
+	X = Y;
+	for (i = 0; i <= order; i++) {
+	  *z += C[i][j]*X;
+	  X = X * (*x);
+	}
+	Y = Y * (*y);
+      }
+    }
+  }
+  if ((xvec[0].type == OPIHI_INT) && (yvec[0].type == OPIHI_INT)) {
+    opihi_int *x = xvec[0].elements.Int;
+    opihi_int *y = yvec[0].elements.Int;
+    for (n = 0; n < xvec[0].Nelements; n++, x++, y++, z++) {
+      Y = X = 1;
+      for (j = 0; j <= order; j++) {
+	X = Y;
+	for (i = 0; i <= order; i++) {
+	  *z += C[i][j]*X;
+	  X = X * (*x);
+	}
+	Y = Y * (*y);
+      }
+    }
+  }
+
+  for (i = 0; i < order + 1; i++) free (C[i]);
+  free (C);
+
+  return (TRUE);
+}
+
+
Index: trunk/Ohana/src/opihi/cmd.data/chebyshev.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/chebyshev.c	(revision 41891)
+++ trunk/Ohana/src/opihi/cmd.data/chebyshev.c	(revision 41891)
@@ -0,0 +1,167 @@
+# include "data.h"
+
+int chebyshev_help (int argc, char **argv);
+int chebyshev_list (int argc, char **argv);
+int chebyshev_delete (int argc, char **argv);
+int chebyshev_rename (int argc, char **argv);
+
+int chebyshev_poly (int argc, char **argv);
+int chebyshev_scale (int argc, char **argv);
+int chebyshev_applyscale (int argc, char **argv);
+
+int chebyshev_fit1d (int argc, char **argv);
+int chebyshev_applyfit1d (int argc, char **argv);
+
+int chebyshev_fit2d (int argc, char **argv);
+int chebyshev_applyfit2d (int argc, char **argv);
+
+// collection of chebyshev commands:
+static Command chebyshev_commands[] = {
+  {1, "help",       chebyshev_help,        "chebyshev help info"},
+  {1, "list",       chebyshev_list,        "list define chebyshev polynomial sets"},
+  {1, "delete",     chebyshev_delete,      "delete a chebyshev"},
+  {1, "rename",     chebyshev_rename,      "rename a chebyshev"},
+
+  {1, "poly",       chebyshev_poly,        "return evaluated chebyshev polynomial of given order for supplied coordinate vector"},
+
+  {1, "scale",      chebyshev_scale,       "set domain scale based on supplied coordinate vector"},
+  {1, "applyscale", chebyshev_applyscale,  "apply domain scale to supplied coordinate vector"},
+
+  {1, "fit1d",      chebyshev_fit1d,       "fit chebyshev polynomial to supplied data vector"},
+  {1, "applyfit1d", chebyshev_applyfit1d,  "apply chebyshev polynomial fit to supplied coordinate vector"},
+
+  {1, "fit2d",      chebyshev_fit2d,       "fit chebyshev polynomial to supplied data vector"},
+  {1, "applyfit2d", chebyshev_applyfit2d,  "apply chebyshev polynomial fit to supplied coordinate vector"},
+};
+
+int chebyshev_command (int argc, char **argv) {
+
+  int N;
+
+  if (argc < 2) {
+    chebyshev_help(0,NULL);
+    return (FALSE);
+  }
+
+  N = sizeof (chebyshev_commands) / sizeof (Command);
+
+  /* find the chebyshev sub-command which matches */
+  for (int i = 0; i < N; i++) {
+    if (!strcmp (chebyshev_commands[i].name, argv[1])) {
+      int status = (*chebyshev_commands[i].func) (argc - 1, argv + 1);
+      return (status);
+    }
+  }
+
+  gprint (GP_ERR, "unknown chebyshev command %s\n", argv[1]);
+  return (FALSE);
+}
+
+/*** below are the utility commands which just manage the collection of chebyshevs ***/
+
+int chebyshev_help (int argc, char **argv) {
+  OHANA_UNUSED_PARAM(argc);
+  OHANA_UNUSED_PARAM(argv);
+
+  gprint (GP_ERR, "USAGE: chebyshev (command)\n");
+  gprint (GP_ERR, "    chebyshev help                         : this listing\n");
+  gprint (GP_ERR, "    chebyshev list                         : list chebyshevs\n");
+  gprint (GP_ERR, "    chebyshev delete     (chebyshev)       : delete named chebyshev\n");
+  gprint (GP_ERR, "    chebyshev rename     (chebyshev) (new) : change chebyshev name to new name\n");
+
+  gprint (GP_ERR, "    chebyshev poly       (chebyshev) (new)    : change chebyshev name to new name\n");
+  gprint (GP_ERR, "    chebyshev scale      (chebyshev) (new)    : change chebyshev name to new name\n");
+  gprint (GP_ERR, "    chebyshev applyscale (chebyshev) (new)    : change chebyshev name to new name\n");
+
+  gprint (GP_ERR, "    chebyshev fit1d      (chebyshev) (new)    : change chebyshev name to new name\n");
+  gprint (GP_ERR, "    chebyshev applyfit1d (chebyshev) (new)    : change chebyshev name to new name\n");
+
+  gprint (GP_ERR, "    chebyshev fit2d      (chebyshev) (new)    : change chebyshev name to new name\n");
+  gprint (GP_ERR, "    chebyshev applyfit2d (chebyshev) (new)    : change chebyshev name to new name\n");
+
+  return FALSE;
+}
+
+int chebyshev_list (int argc, char **argv) {
+  OHANA_UNUSED_PARAM(argv);
+
+  if (argc != 1) {
+    gprint (GP_ERR, "USAGE: chebyshev list\n");
+    return FALSE;
+  }
+
+  ListChebyshevs();
+  return TRUE;
+}
+
+int chebyshev_delete (int argc, char **argv) {
+
+  int N;
+
+  int QUIET = FALSE;
+  if ((N = get_argument (argc, argv, "-q"))) {
+    QUIET = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+
+  if (argc != 2) {
+    gprint (GP_ERR, "USAGE: chebyshev delete [-q] (chebyshev)\n");
+    return FALSE;
+  }
+
+  ChebyshevType *chebyshev = FindChebyshev (argv[1]);
+  if (chebyshev == NULL) {
+    if (QUIET) return TRUE;
+    gprint (GP_ERR, "chebyshev %s not found\n", argv[1]);
+    return FALSE;
+  }
+
+  int status = DeleteChebyshev (chebyshev); 
+  if (!status) abort (); // status = FALSE means the cheb selected above was not found
+  return TRUE;
+}
+
+int chebyshev_rename (int argc, char **argv) {
+
+  ChebyshevType *chebyshev;
+
+  if (argc != 3) {
+    gprint (GP_ERR, "USAGE: chebyshev rename (chebyshev) (newname)\n");
+    return FALSE;
+  }
+
+  chebyshev = FindChebyshev (argv[1]);
+  if (chebyshev == NULL) {
+    gprint (GP_ERR, "chebyshev %s not found\n", argv[1]);
+    return FALSE;
+  }
+
+  free (chebyshev->name);
+  chebyshev->name = strcreate (argv[2]);
+  return TRUE;
+}
+
+/* 
+
+   typedef struct {
+     char *name;
+     opihi_flt scale;
+     opihi_flt zero;
+     int   order;
+     opihi_flt *A;
+   } ChebyshevType;
+
+
+   examples of using the 
+
+   cheb fit1d t1 x y 4     -- fit the vector y = f(x) using 4th order
+    (saves the scaling for x, order, coeffs)
+   cheb applyfit1d t1 x yfit  -- apply the fit defined in t1 to x
+   cheb applyfit1d t1 X Yfit  -- apply the fit defined in t1 to a different vector X
+
+   cheb fit t2 x y 2
+
+
+
+ */
+
Index: trunk/Ohana/src/opihi/cmd.data/chebyshev_commands.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/chebyshev_commands.c	(revision 41891)
+++ trunk/Ohana/src/opihi/cmd.data/chebyshev_commands.c	(revision 41891)
@@ -0,0 +1,344 @@
+# include "data.h"
+
+int chebyshev_poly (int argc, char **argv) {
+
+  if (argc != 4) {
+    gprint (GP_ERR, "USAGE: chebyshev poly (xvec) (yvec) (order)\n");
+    gprint (GP_ERR, "       return requested chebyshev polynomial for given order\n");
+    gprint (GP_ERR, "       xvec should be normalized (range of -1 to +1)\n");
+    return FALSE;
+  }
+
+  Vector *xvec = NULL, *yvec = NULL;
+  if ((xvec = SelectVector (argv[1], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((yvec = SelectVector (argv[2], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  int order = atoi(argv[3]);
+
+  ResetVector (yvec, OPIHI_FLT, xvec->Nelements);
+
+  ChebyshevPolyVector (yvec, xvec, order);
+
+  return TRUE;
+}
+
+int chebyshev_fit1d (int argc, char **argv) {
+
+  int N;
+
+  int VERBOSE = FALSE;
+  if ((N = get_argument (argc, argv, "-v"))) {
+    VERBOSE = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+
+  if (argc != 5) {
+    gprint (GP_ERR, "USAGE: chebyshev fit1d (name) (xv) (yv) (order)\n");
+    gprint (GP_ERR, "       1D fit of xv to yv using given order\n");
+    return FALSE;
+  }
+
+  ChebyshevType *cheb = FindChebyshev (argv[1]);
+  if (!cheb) {
+    if (VERBOSE) gprint (GP_ERR, "chebyshev %s not found, creating it\n", argv[1]);
+    cheb = CreateChebyshev (argv[1]);
+  }
+
+  // select the x,y vectors
+  Vector *xvec = NULL, *yvec = NULL;
+  if ((xvec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((yvec = SelectVector (argv[3], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+
+  cheb->dimen = 1;
+  cheb->order = atoi (argv[4]);
+  int nterm = cheb->order + 1;
+
+  // check type and size:
+  if (xvec->Nelements < nterm) {
+    gprint (GP_ERR, "insufficient data to support the requested order\n");
+    return FALSE;
+  }
+  if (xvec->Nelements != yvec->Nelements) {
+    gprint (GP_ERR, "vector lengths %s (%d) and %s (%d) do not match\n", xvec->name, xvec->Nelements, yvec->name, yvec->Nelements);
+    return FALSE;
+  }
+  CastVector (xvec, OPIHI_FLT);
+  CastVector (yvec, OPIHI_FLT);
+
+  ChebyshevSetScale (cheb, xvec, 0);
+
+  // generate the normalized vector xNorm
+  // use a local vector?
+  Vector *xNorm = ChebyshevNormVector (cheb, xvec, 0);
+
+  // generate the N cheb polynomials based on xNorm
+  ALLOCATE_PTR (poly, Vector *, nterm);
+  for (int i = 0; i < nterm; i++) {
+    poly[i] = ChebyshevPolyVector (NULL, xNorm, i);
+  }
+
+  // fit the data to the cheb polynomials
+  ChebyshevPolyFit1D (cheb, yvec, poly);
+
+  for (int i = 0; i < nterm; i++) {
+    fprintf (stderr, "%d : %e\n", i, cheb->A[i]);
+  }
+
+  // free local temporary variables
+  for (int i = 0; i < nterm; i++) {
+    FREE (poly[i]);
+  }
+  FREE (poly);
+  FREE (xNorm);
+
+  return TRUE;
+}
+
+int chebyshev_applyfit1d (int argc, char **argv) {
+
+  if (argc != 4) {
+    gprint (GP_ERR, "USAGE: chebyshev applyfit1d (name) (xvec) (yfit)\n");
+    gprint (GP_ERR, "       apply 1D fit to xvec position to yield yfit\n");
+    return FALSE;
+  }
+
+  ChebyshevType *cheb = FindChebyshev (argv[1]);
+  if (!cheb) {
+    gprint (GP_ERR, "chebyshev %s not found\n", argv[1]);
+    return FALSE;
+  }
+  int nterm = cheb->order + 1;
+
+  Vector *xvec = NULL, *yvec = NULL;
+  if ((xvec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((yvec = SelectVector (argv[3], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+
+  ResetVector (yvec, OPIHI_FLT, xvec->Nelements);
+
+  // generate the normalized vector xNorm
+  Vector *xNorm = ChebyshevNormVector (cheb, xvec, 0);
+
+  // generate the N cheb polynomials based on xNorm
+  ALLOCATE_PTR (poly, Vector *, nterm);
+  for (int i = 0; i < nterm; i++) {
+    poly[i] = ChebyshevPolyVector (NULL, xNorm, i);
+  }
+
+  ChebyshevPolyApplyFit1D (cheb, yvec, poly);
+
+  // free local temporary variables
+  for (int i = 0; i < nterm; i++) {
+    FREE (poly[i]);
+  }
+  FREE (poly);
+  FREE (xNorm);
+
+  return TRUE;
+}
+
+int chebyshev_fit2d (int argc, char **argv) {
+
+  int N;
+
+  int VERBOSE = FALSE;
+  if ((N = get_argument (argc, argv, "-v"))) {
+    VERBOSE = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+
+  if (argc != 6) {
+    gprint (GP_ERR, "USAGE: chebyshev fit1d (name) (xv) (yv) (zv) (order)\n");
+    gprint (GP_ERR, "       2D fit of (xv,yv) to zv using given order\n");
+    return FALSE;
+  }
+
+  ChebyshevType *cheb = FindChebyshev (argv[1]);
+  if (!cheb) {
+    if (VERBOSE) gprint (GP_ERR, "chebyshev %s not found, creating it\n", argv[1]);
+    cheb = CreateChebyshev (argv[1]);
+  }
+
+  // select the x,y vectors
+  Vector *xvec = NULL, *yvec = NULL, *zvec = NULL;
+  if ((xvec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((yvec = SelectVector (argv[3], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((zvec = SelectVector (argv[4], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+
+  cheb->dimen = 2;
+  cheb->order = atoi (argv[5]);
+  int nterm = cheb->order + 1;
+
+  // check type and size:
+  if (xvec->Nelements < nterm) {
+    gprint (GP_ERR, "insufficient data to support the requested order\n");
+    return FALSE;
+  }
+  if (xvec->Nelements != yvec->Nelements) {
+    gprint (GP_ERR, "vector lengths %s (%d) and %s (%d) do not match\n", xvec->name, xvec->Nelements, yvec->name, yvec->Nelements);
+    return FALSE;
+  }
+  CastVector (xvec, OPIHI_FLT);
+  CastVector (yvec, OPIHI_FLT);
+  CastVector (zvec, OPIHI_FLT);
+
+  ChebyshevSetScale (cheb, xvec, 0);
+  ChebyshevSetScale (cheb, yvec, 1);
+
+  // generate the normalized vector xNorm
+  // use a local vector?
+  Vector *xNorm = ChebyshevNormVector (cheb, xvec, 0);
+  Vector *yNorm = ChebyshevNormVector (cheb, yvec, 1);
+
+  // generate the N cheb polynomials based on xNorm, yNorm
+  ALLOCATE_PTR (xPoly, Vector *, nterm);
+  ALLOCATE_PTR (yPoly, Vector *, nterm);
+  for (int i = 0; i < nterm; i++) {
+    xPoly[i] = ChebyshevPolyVector (NULL, xNorm, i);
+    yPoly[i] = ChebyshevPolyVector (NULL, yNorm, i);
+  }
+
+  // fit the data to the cheb polynomials
+  ChebyshevPolyFit2D (cheb, zvec, xPoly, yPoly);
+
+  int n = 0;
+  for (int ix = 0; ix < nterm; ix++) {
+    for (int iy = 0; iy < nterm; iy++) {
+      fprintf (stderr, "%d %d : %e\n", ix, iy, cheb->A[n]);
+      n ++;
+    }
+  }
+
+  // free local temporary variables
+  for (int i = 0; i < nterm; i++) {
+    FREE (xPoly[i]);
+    FREE (yPoly[i]);
+  }
+  FREE (xPoly);
+  FREE (yPoly);
+  FREE (xNorm);
+  FREE (yNorm);
+
+  return TRUE;
+}
+
+int chebyshev_applyfit2d (int argc, char **argv) {
+
+  if (argc != 5) {
+    gprint (GP_ERR, "USAGE: chebyshev applyfit1d (name) (xvec) (yvec) (xfit)\n");
+    gprint (GP_ERR, "       apply 1D fit to xvec,yvec positions to yield zfit\n");
+    return FALSE;
+  }
+
+  ChebyshevType *cheb = FindChebyshev (argv[1]);
+  if (!cheb) {
+    gprint (GP_ERR, "chebyshev %s not found\n", argv[1]);
+    return FALSE;
+  }
+  int nterm = cheb->order + 1;
+
+  Vector *xvec = NULL, *yvec = NULL, *zvec = NULL;
+  if ((xvec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((yvec = SelectVector (argv[3], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((zvec = SelectVector (argv[4], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+
+  if (xvec->Nelements != yvec->Nelements) {
+    gprint (GP_ERR, "vector lengths %s (%d) and %s (%d) do not match\n", xvec->name, xvec->Nelements, yvec->name, yvec->Nelements);
+    return FALSE;
+  }
+
+  ResetVector (zvec, OPIHI_FLT, xvec->Nelements);
+
+  // generate the normalized vector xNorm
+  // use a local vector?
+  Vector *xNorm = ChebyshevNormVector (cheb, xvec, 0);
+  Vector *yNorm = ChebyshevNormVector (cheb, yvec, 1);
+
+  // generate the N cheb polynomials based on xNorm, yNorm
+  ALLOCATE_PTR (xPoly, Vector *, nterm);
+  ALLOCATE_PTR (yPoly, Vector *, nterm);
+  for (int i = 0; i < nterm; i++) {
+    xPoly[i] = ChebyshevPolyVector (NULL, xNorm, i);
+    yPoly[i] = ChebyshevPolyVector (NULL, yNorm, i);
+  }
+
+
+  ChebyshevPolyApplyFit2D (cheb, zvec, xPoly, yPoly);
+
+  // free local temporary variables
+  for (int i = 0; i < nterm; i++) {
+    FREE (xPoly[i]);
+    FREE (yPoly[i]);
+  }
+  FREE (xPoly);
+  FREE (yPoly);
+  FREE (xNorm);
+  FREE (yNorm);
+
+  return TRUE;
+}
+
+int chebyshev_scale (int argc, char **argv) {
+
+  if (argc != 3) {
+    gprint (GP_ERR, "USAGE: chebyshev scale (name) (vector)\n");
+    gprint (GP_ERR, "       determine scale factors for the names chebyshev based on the vector\n");
+    return FALSE;
+  }
+
+  ChebyshevType *cheb = FindChebyshev (argv[1]);
+  if (!cheb) {
+    gprint (GP_ERR, "chebyshev %s not found, creating it\n", argv[1]);
+    cheb = CreateChebyshev (argv[1]);
+  } else {
+    gprint (GP_ERR, "resetting scale factors for %s\n", argv[1]);
+  }
+
+  Vector *vec = NULL;
+  if ((vec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+
+  // check type and size:
+  if (vec->type != OPIHI_FLT) {
+    gprint (GP_ERR, "cannot generate chebyshev scale factors for non-float vector\n");
+    return FALSE;
+  }
+  if (!vec->Nelements) {
+    gprint (GP_ERR, "cannot generate chebyshev scale factors based on empty vector\n");
+    return FALSE;
+  }
+
+  int dir = 0;
+  ChebyshevSetScale (cheb, vec, dir);
+
+  return TRUE;
+}
+
+int chebyshev_applyscale (int argc, char **argv) {
+
+  if (argc != 4) {
+    gprint (GP_ERR, "USAGE: chebyshev applyscale (name) (vector) (output)\n");
+    gprint (GP_ERR, "       re-scale vector for chebyshev\n");
+    return FALSE;
+  }
+
+  ChebyshevType *cheb = FindChebyshev (argv[1]);
+  if (!cheb) {
+    gprint (GP_ERR, "chebyshev %s not found\n", argv[1]);
+    return FALSE;
+  }
+
+  Vector *vec = NULL, *out = NULL;
+  if ((vec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((out = SelectVector (argv[3], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  CastVector (vec, OPIHI_FLT);
+
+  ResetVector (out, OPIHI_FLT, vec->Nelements);
+
+  int dir = 0;
+
+  for (int i = 0; i < vec->Nelements; i++) {
+    if (!isfinite(vec->elements.Flt[i])) { out->elements.Flt[i] = NAN; continue; }
+    out->elements.Flt[i] = cheb->scale[dir]*vec->elements.Flt[i] + cheb->zero[dir];
+  }
+
+  return TRUE;
+}
+
Index: trunk/Ohana/src/opihi/cmd.data/fit2d.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/fit2d.c	(revision 41715)
+++ trunk/Ohana/src/opihi/cmd.data/fit2d.c	(revision 41891)
@@ -294,2 +294,272 @@
   return (TRUE);
 }
+
+int fit2d_full (int argc, char **argv) {
+  
+  int N;
+  char name[64];
+
+  int VERBOSE = FALSE;
+  if ((N = get_argument (argc, argv, "-v"))) {
+    VERBOSE = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+
+  int Quiet = FALSE;
+  if ((N = get_argument (argc, argv, "-q"))) {
+    Quiet = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+  if ((N = get_argument (argc, argv, "-quiet"))) {
+    Quiet = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+
+  float ClipNSigma = 0;
+  int ClipNiter  = 1;
+  if ((N = get_argument (argc, argv, "-clip"))) {
+    remove_argument (N, &argc, argv);
+    ClipNSigma = atof(argv[N]);
+    remove_argument (N, &argc, argv);
+    ClipNiter  = atof(argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+
+  opihi_flt *dz = NULL;
+  Vector *dzvec = NULL;
+  int Weight = FALSE;
+  if ((N = get_argument (argc, argv, "-dz"))) {
+    remove_argument (N, &argc, argv);
+    if ((dzvec = SelectVector (argv[N], OLDVECTOR, TRUE)) == NULL) return (FALSE);    
+    remove_argument (N, &argc, argv);
+    Weight = TRUE;
+  }
+
+  if (argc != 5) {
+    gprint (GP_ERR, "USAGE: fit2d_full x y z order [-dz wt] [-quiet/-q] [-clip Nsigma Niter]\n");
+    return (FALSE);
+  }
+
+  Vector *xvec, *yvec, *zvec;
+  if ((xvec = SelectVector (argv[1], OLDVECTOR, TRUE)) == NULL) return (FALSE);    
+  if ((yvec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);    
+  if ((zvec = SelectVector (argv[3], OLDVECTOR, TRUE)) == NULL) return (FALSE);    
+
+  if (xvec[0].Nelements != yvec[0].Nelements) {
+    gprint (GP_ERR, "vectors must have same length\n");
+    return (FALSE);
+  }
+  if (xvec[0].Nelements != zvec[0].Nelements) {
+    gprint (GP_ERR, "vectors must have same length\n");
+    return (FALSE);
+  }
+  CastVector (xvec, OPIHI_FLT);
+  CastVector (yvec, OPIHI_FLT);
+  CastVector (zvec, OPIHI_FLT);
+
+  if (Weight) {
+    CastVector (dzvec, OPIHI_FLT);
+    if (xvec[0].Nelements != dzvec[0].Nelements) {
+      gprint (GP_ERR, "vectors must have same length\n");
+      return (FALSE);
+    }
+  }
+  
+  int order = atof (argv[4]);
+  int nterm = SQ(order + 1);
+
+  ALLOCATE_PTR (zfit, opihi_flt, xvec[0].Nelements);
+  ALLOCATE_PTR (mask, char, xvec[0].Nelements);
+  memset (mask, 0, xvec[0].Nelements);
+
+  /* allocate the summation matrices */
+  ALLOCATE_PTR (b, double *, nterm);
+  ALLOCATE_PTR (c, double *, nterm);
+  for (int i = 0; i < nterm; i++) {
+    ALLOCATE (c[i], double, nterm);
+    ALLOCATE (b[i], double, 1);
+  }
+
+  for (N = 0; N < ClipNiter; N++) {
+
+    /* init registers for current pass */
+    for (int i = 0; i < nterm; i++) {
+      memset (c[i], 0, nterm*sizeof(double));
+      memset (b[i], 0, sizeof(double));
+    }
+
+    opihi_flt *x = xvec[0].elements.Flt;
+    opihi_flt *y = yvec[0].elements.Flt;
+    opihi_flt *z = zvec[0].elements.Flt;
+    if (Weight) dz = dzvec[0].elements.Flt;
+
+    /* add up the x,y values */
+    for (int i = 0; i < xvec[0].Nelements; i++, x++, y++, z++) {
+      if (mask[i]) continue;
+      if (!finite(*x) || !finite(*y)) continue;
+      double dZ = 1.0;
+      if (Weight) { 
+	dZ = 1.0 / SQ(*dz);
+	dz ++;
+      }
+
+      int n = 0;
+      double iY = 1*dZ;
+      double iX = iY;
+
+      for (int iy = 0; iy <= order; iy++) {
+	iX = iY;
+	for (int ix = 0; ix <= order; ix++) {
+	  b[n][0] += iX * (*z);
+
+	  int m = 0;
+	  double jY = 1;
+	  double jX = jY;
+
+	  for (int jy = 0; jy <= order; jy++) {
+	    jX = jY;
+	    for (int jx = 0; jx <= order; jx++) {
+	      c[n][m] += iX*jX;
+	      jX = jX * (*x);
+	      m ++;
+	    }
+	    jY = jY * (*y);
+	  }
+	  iX = iX * (*x);
+	  n ++;
+	}
+	iY = iY * (*y);
+      }
+    }
+
+    /** test print **/
+    if (VERBOSE) {
+      for (int i = 0; i < nterm; i++) {
+	for (int j = 0; j < nterm; j++) {
+	  gprint (GP_ERR, "%g  ", c[i][j]);
+	}
+	gprint (GP_ERR, "\n");
+      }
+      gprint (GP_ERR, "-----\n");
+    }
+
+    dgaussjordan (c, b, nterm, 1);
+
+    /** test print **/
+    if (VERBOSE) {
+      int n = 0;
+      for (int ny = 0; ny <= order; ny++) {
+	for (int nx = 0; nx <= order; nx++) {
+	  gprint (GP_ERR, "x^%d y^%d: %g\n", nx, ny, b[n][0]);
+	  n ++;
+	}
+      }
+    }
+
+    /** test print **/
+    if (VERBOSE) {
+      for (int i = 0; i < nterm; i++) {
+	for (int j = 0; j < nterm; j++) {
+	  gprint (GP_ERR, "%g  ", c[i][j]);
+	}
+	gprint (GP_ERR, "\n");
+      }
+      gprint (GP_ERR, "-----\n");
+    }
+
+    /* the b[][0] terms are in the following order:
+       y^0 x^0, y^0 x^1, ... y^0 x^N
+       y^1 x^0, y^1 x^1, ... y^1 x^N
+       ...
+       y^N x^0, y^N x^1, ... y^N x^N
+    */
+
+    /* generate fitted values */
+    x = xvec[0].elements.Flt;
+    y = yvec[0].elements.Flt;
+    opihi_flt *zf = zfit;
+    for (int n = 0; n < xvec[0].Nelements; n++, x++, y++, zf++) {
+      if (!finite(*x) || !finite(*y)) continue;
+      *zf = 0;
+      double Y = 1;
+      double X = Y;
+
+      int n = 0;
+      for (int ny = 0; ny <= order; ny++) {
+	X = Y;
+	for (int nx = 0; nx <= order; nx++) {
+	  *zf += b[n][0]*X;
+	  X = X * (*x);
+	  n ++;
+	}
+	Y = Y * (*y);
+      }
+    }
+
+    /* measure fit residual scatter */
+    x  = xvec[0].elements.Flt;
+    y  = yvec[0].elements.Flt;
+    z  = zvec[0].elements.Flt;
+    zf = zfit;
+    double dZ = 0.0;
+    double dZ2 = 0.0;
+    int Npt = 0;
+    for (int i = 0, Npt = 0; i < xvec[0].Nelements; i++, x++, y++, z++, zf++) {
+      if (mask[i]) continue;
+      if (!finite(*x) || !finite(*y)) continue;
+      dZ  += (*z - *zf);
+      dZ2 += SQ(*z - *zf);
+      Npt ++;
+    }
+    double mean  = dZ / Npt;
+    double sigma = sqrt (fabs(dZ2/Npt - SQ(mean)));
+    double maxsigma = ClipNSigma * sigma;
+
+    if (VERBOSE) gprint (GP_ERR, "mean: %g, sigma: %g, maxsigma: %g\n", mean, sigma, maxsigma);
+
+    /* mask outlier points */
+    x  = xvec[0].elements.Flt;
+    y  = yvec[0].elements.Flt;
+    z  = zvec[0].elements.Flt;
+    zf = zfit;
+    int Nmask = 0;
+    for (int i = 0; ClipNSigma && (i < xvec[0].Nelements); i++, x++, y++, z++, zf++) {
+      dZ = (*z - *zf);
+      if (fabs(dZ) > maxsigma) {
+	mask[i] = TRUE;
+	Nmask ++;
+      } else {
+	mask[i] = FALSE;
+      }	
+    }
+    if (VERBOSE) gprint (GP_ERR, "pass: %d, Nmask: %d\n", N, Nmask);
+  }
+
+  if (!Quiet) gprint (GP_ERR, "z = ");
+  int n = 0;
+  for (int ny = 0; ny <= order; ny++) {
+    for (int nx = 0; nx <= order; nx++) {
+      sprintf (name, "CX%dY%d", nx, ny);
+      set_variable (name, b[n][0]);
+      if (!Quiet) gprint (GP_ERR, "%f x^%d y^%d  ", b[n][0], nx, ny);
+      n++;
+    }
+  }
+  sprintf (name, "Cnn");
+  set_variable (name, (double) order);
+  
+  if (!Quiet) gprint (GP_ERR, "\n");
+
+  /* free internal data */
+  free (zfit);
+  free (mask);
+
+  for (int i = 0; i < nterm; i++) {
+    free (c[i]);
+    free (b[i]);
+  }
+  free (b);
+  free (c);
+
+  return (TRUE);
+}
Index: trunk/Ohana/src/opihi/cmd.data/init.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/init.c	(revision 41715)
+++ trunk/Ohana/src/opihi/cmd.data/init.c	(revision 41891)
@@ -7,4 +7,5 @@
 int applyfit1d       PROTO((int, char **));
 int applyfit2d       PROTO((int, char **));
+int applyfit2d_full  PROTO((int, char **));
 int applyfit3d       PROTO((int, char **));
 int box              PROTO((int, char **));
@@ -13,4 +14,5 @@
 int center           PROTO((int, char **));
 int cast             PROTO((int, char **));
+int chebyshev_command PROTO((int, char **));
 int circstats        PROTO((int, char **));
 int clear            PROTO((int, char **));
@@ -45,4 +47,5 @@
 int fit1d_irls       PROTO((int, char **));
 int fit2d            PROTO((int, char **));
+int fit2d_full       PROTO((int, char **));
 int fit3d            PROTO((int, char **));
 int gaussjordan      PROTO((int, char **));
@@ -204,4 +207,5 @@
   {1, "applyfit1d",   applyfit1d,       "apply 1-d fit to new vector"},
   {1, "applyfit2d",   applyfit2d,       "apply 2-d fit to new vector"},
+  {1, "applyfit2d_full", applyfit2d_full, "apply 2-d fit to new vector"},
   {1, "applyfit3d",   applyfit3d,       "apply 3-d fit to new vector"},
   {1, "bisection",    bisection,        "use bisection to find threshold in vector"},
@@ -211,4 +215,5 @@
   {1, "center",       center,           "center image on coords"},
   {1, "cast",         cast,             "cast input vector to specified type"},
+  {1, "chebyshev",    chebyshev_command, "Chebyshev polynomial commands"},
   {1, "circstats",    circstats,        "circular statistics"},
   {1, "clear",        clear,            "erase plot"},
@@ -244,4 +249,5 @@
   {1, "fit1d_irls",   fit1d_irls,       "fit polynomial to vector pair"},
   {1, "fit2d",        fit2d,            "fit 2-d polynomial to vector triplet"},
+  {1, "fit2d_full",   fit2d_full,       "fit 2-d polynomial to vector triplet (all X^n Y^m coeffs)"},
   {1, "fit3d",        fit3d,            "fit 3-d polynomial to vector quad"},
   {1, "gaussdev",     gaussdeviate,     "generate a gaussian deviate vector"},
Index: trunk/Ohana/src/opihi/cmd.data/medimage_calc.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/medimage_calc.c	(revision 41715)
+++ trunk/Ohana/src/opihi/cmd.data/medimage_calc.c	(revision 41891)
@@ -19,5 +19,5 @@
 static float *irls_testval = NULL;
 
-enum {CALC_MEDIAN, CALC_MEAN, CALC_IRLS, CALC_WTMEAN};
+enum {CALC_MEDIAN, CALC_MEAN, CALC_IRLS, CALC_WTMEAN, CALC_INNER_FRACTION};
 
 int medimage_calc (int argc, char **argv) {
@@ -39,4 +39,7 @@
   }
 
+  float minRange = NAN;
+  float maxRange = NAN;
+
   int mode = CALC_MEDIAN;
   if ((N = get_argument (argc, argv, "-mean"))) {
@@ -45,15 +48,28 @@
   }
   if ((N = get_argument (argc, argv, "-irls"))) {
-    if (mode != CALC_MEDIAN) { gprint (GP_ERR, "supply only one of -mean, -median, -irls, -wtmean\n"); return FALSE; }
+    if (mode != CALC_MEDIAN) goto choose_one;
     mode = CALC_IRLS;
     remove_argument (N, &argc, argv);
   }
   if ((N = get_argument (argc, argv, "-wtmean"))) {
-    if (mode != CALC_MEDIAN) { gprint (GP_ERR, "supply only one of -mean, -median, -irls, -wtmean\n"); return FALSE; }
+    if (mode != CALC_MEDIAN) goto choose_one;
     mode = CALC_WTMEAN;
     remove_argument (N, &argc, argv);
   }
+  if ((N = get_argument (argc, argv, "-inner-fraction"))) {
+    if (mode != CALC_MEDIAN) goto choose_one;
+    if (argc < N + 3) goto fraction_options;
+    mode = CALC_INNER_FRACTION;
+    remove_argument (N, &argc, argv);
+    minRange = atof(argv[N]);
+    if ((minRange < 0.0) || (minRange > 1.0)) goto fraction_options;
+    remove_argument (N, &argc, argv);
+    maxRange = atof(argv[N]);
+    if ((maxRange < 0.0) || (maxRange > 1.0)) goto fraction_options;
+    if (maxRange <= minRange) goto fraction_options;
+    remove_argument (N, &argc, argv);
+  }
   if ((N = get_argument (argc, argv, "-median"))) {
-    if (mode != CALC_MEDIAN) { gprint (GP_ERR, "supply only one of -mean, -median, -irls, -wtmean\n"); return FALSE; }
+    if (mode != CALC_MEDIAN) goto choose_one;
     mode = CALC_MEDIAN;
     remove_argument (N, &argc, argv);
@@ -157,4 +173,25 @@
 	  }
 	  break;
+	case CALC_INNER_FRACTION:
+	  fsort (val, N);
+	  int Ns = MIN(MAX(0, minRange * N),N); // e.g., 0.1 * 50 = 5
+	  int Ne = MIN(MAX(0, maxRange * N),N); // e.g., 0.9 * 50 = 45 
+	  int Npt = Ne - Ns;
+
+	  float sum = 0.0;
+	  for (n = Ns; n < Ne; n++) {
+	    sum += val[n];
+	  }
+	  outvalue[Npix] = sum / (float) Npt;
+
+	  if (varvalue) {
+	    float sum = 0.0;
+	    for (n = Ns; n < Ne; n++) {
+	      sum += SQ(val[n] - outvalue[Npix]);
+	    }
+	    // variance on the mean (stdev / sqrt(N))^2
+	    varvalue[Npix] = sum / (Npt - 1) / Npt;
+	  }
+	  break;
 	case CALC_MEAN: {
 	  float sum = 0.0;
@@ -195,4 +232,12 @@
   irls_free();
   return TRUE;
+
+ choose_one:
+  gprint (GP_ERR, "supply only one of -mean, -median, -irls, -wtmean, -inner-fraction\n");
+  return FALSE; 
+
+ fraction_options:
+  gprint (GP_ERR, "USE: -inner-fraction (min) (max) : min & max in range 0.0 to 1.0\n");
+  return FALSE; 
 }
 
Index: trunk/Ohana/src/opihi/cmd.data/tdhistogram.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/tdhistogram.c	(revision 41715)
+++ trunk/Ohana/src/opihi/cmd.data/tdhistogram.c	(revision 41891)
@@ -45,6 +45,8 @@
   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, "USAGE: tdhistogram buffer x y z (Xmin) (Xmax) (Ymin) (Ymax) (Zmin) (Zmax) [-dx dx] [-dy dy] [-dz dz] [-range range]\n");
     gprint (GP_ERR, "   OR: tdhistogram buffer x y z -reuse\n");
+    gprint (GP_ERR, "   -dx, -dy, -dy specify the bin size in these directions\n");
+    gprint (GP_ERR, "   -range : generate an output vector corresponding to the elements in the z-direction\n");
     gprint (GP_ERR, " output buffer is 3D\n");
     return (FALSE);
Index: trunk/Ohana/src/opihi/cmd.data/vstats.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/vstats.c	(revision 41715)
+++ trunk/Ohana/src/opihi/cmd.data/vstats.c	(revision 41891)
@@ -48,6 +48,8 @@
   }
 
-  if (argc != 2) {
-    gprint (GP_ERR, "USAGE: vstat (vector) [-ignore value] [-q] [-quiet] [-iter Niter] [-sigma-clip Nsigma]\n");
+  int valid = (argc == 2);
+  valid |= (argc > 3) && !strcmp (argv[2], "where");
+  if (!valid) {
+    gprint (GP_ERR, "USAGE: vstat (vector) [-ignore value] [-q] [-quiet] [-iter Niter] [-sigma-clip Nsigma] [where (logical expression)]\n");
     gprint (GP_ERR, "  default is 1 iteration without sigma clipping or 3 with sigma clipping\n");
     return (FALSE);
@@ -58,21 +60,48 @@
   /* we need two passes, one for max, min, mean, sum, one for median, stdev, etc */
 
-  // set a good / bad mask
+  // tvec is used for logical test (truth vector)
+  Vector *tvec = NULL;
+  if (argc > 3) {
+    int size;
+    char *out = dvomath (argc - 3, &argv[3], &size, 1);
+    if (out == NULL) {
+      print_error ();
+      return FALSE;
+    }
+    if ((tvec = SelectVector (out, OLDVECTOR, TRUE)) == NULL) {
+      gprint (GP_ERR, " invalid logic result\n");
+      DeleteNamedVector (out);
+      free (out);
+      return (FALSE);
+    }
+    if (tvec->Nelements != vec->Nelements) {
+      gprint (GP_ERR, "logical vector does not match in length\n");
+      return FALSE;
+    }
+  }
+
+  // set a good / bad mask (mask == 1 means 'bad')
   ALLOCATE (mask, char, vec[0].Nelements);
   if (vec[0].type == OPIHI_FLT) {
     opihi_flt *X = vec[0].elements.Flt;
     for (i = 0; i < vec[0].Nelements; i++, X++) {
-      mask[i] = 1;
-      if (!finite (*X)) continue;
-      if (Ignore && (*X == IgnoreValue)) continue;
-      mask[i] = 0;
+      mask[i] = 0; // do not mask unless we have a reason below
+      if (tvec) {
+	// note the logical vector is 0 == bad (ignore) while mask[i] has the opposite sense (0 is good)
+	mask[i] = (tvec->type == OPIHI_FLT) ? (tvec->elements.Flt[i] == 0.0) : (tvec->elements.Int[i] == 0);
+      }
+      if (!finite (*X)) { mask[i] = 1; }
+      if (Ignore && (*X == IgnoreValue)) { mask[i] = 1; }
     }      
   } else {
     opihi_int *X = vec[0].elements.Int;
     for (i = 0; i < vec[0].Nelements; i++, X++) {
-      mask[i] = 1;
-      if (!finite (*X)) continue;
-      if (Ignore && (*X == IgnoreValue)) continue;
-      mask[i] = 0;
+      mask[i] = 0; // do not mask unless we have a reason below
+      if (tvec) {
+	// note the logical vector is 0 == bad (ignore) while mask[i] has the opposite sense (0 is good)
+	mask[i] = (tvec->type == OPIHI_FLT) ? (tvec->elements.Flt[i] == 0.0) : (tvec->elements.Int[i] == 0);
+      }
+      if (!finite (*X)) { mask[i] = 1; }
+      if (Ignore && (*X == IgnoreValue)) { mask[i] = 1; }
     }      
   }
