Index: trunk/Ohana/src/opihi/cmd.data/Makefile
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/Makefile	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/Makefile	(revision 37807)
@@ -39,4 +39,5 @@
 $(SRC)/densify.$(ARCH).o	\
 $(SRC)/device.$(ARCH).o	\
+$(SRC)/dft2d.$(ARCH).o	\
 $(SRC)/dimendown.$(ARCH).o	\
 $(SRC)/dimenup.$(ARCH).o	\
@@ -160,4 +161,10 @@
 $(SRC)/vstats.$(ARCH).o		   \
 $(SRC)/xsection.$(ARCH).o          \
+$(SRC)/vsh.$(ARCH).o	   \
+$(SRC)/vshfit.$(ARCH).o	   \
+$(SRC)/shterms.$(ARCH).o	   \
+$(SRC)/shfit.$(ARCH).o	   \
+$(SRC)/shdot.$(ARCH).o	   \
+$(SRC)/shapply.$(ARCH).o	   \
 $(SRC)/wd.$(ARCH).o		   \
 $(SRC)/write_vectors.$(ARCH).o	   \
Index: trunk/Ohana/src/opihi/cmd.data/center.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/center.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/center.c	(revision 37807)
@@ -16,4 +16,11 @@
   if (!GetImage (NULL, &kapa, name)) return (FALSE);
   FREE (name);
+
+  int channel = 0;
+  if ((N = get_argument (argc, argv, "-ch"))) {
+    channel = GetKapaChannelFromString (argv[N]);
+    if (!channel) return FALSE;
+    KiiSetChannel (kapa, channel - 1);
+  }
 
   // XXX need an option to center the image based on the current plot limits
Index: trunk/Ohana/src/opihi/cmd.data/dft2d.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/dft2d.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/dft2d.c	(revision 37807)
@@ -0,0 +1,142 @@
+# include "data.h"
+
+// perform a 2D fourier transform.  
+
+static int NxLast = 0;
+static int NyLast = 0;
+static float *Fterms = NULL;
+
+int dft2d (int argc, char **argv) {
+  
+  int ix, iy, wx, wy;
+  Buffer *src = NULL;
+  Buffer *tgt = NULL;
+
+  if (argc != 4) goto usage;
+  
+  if (strcasecmp(argv[2], "to")) goto usage;
+  if ((src = SelectBuffer (argv[1], OLDBUFFER, TRUE)) == NULL) goto usage;
+  if ((tgt = SelectBuffer (argv[3], ANYBUFFER, TRUE)) == NULL) goto usage;
+
+  int Nx = src[0].matrix.Naxis[0];
+  int Ny = src[0].matrix.Naxis[1];
+  ResetBuffer (tgt, Nx, Ny, -32, 0.0, 1.0);
+
+  if (Nx % 2) goto usage;
+  if (Ny % 2) goto usage;
+
+  // generate and save the fourier terms, if needed
+  if (Fterms && (NxLast == Nx) && (NyLast == Ny)) goto reuse;
+
+  if (!Fterms) {
+    ALLOCATE (Fterms, float, SQ(Nx*Ny));
+    NxLast = Nx;
+    NyLast = Ny;
+  } else {
+    REALLOCATE (Fterms, float, SQ(Nx*Ny));
+    NxLast = Nx;
+    NyLast = Ny;
+  }
+
+  // x-dir cosine terms
+  float Wx = 2*M_PI/(float)Nx;
+  float Wy = 2*M_PI/(float)Ny;
+  for (wx = 0; wx <= Nx / 2; wx++) {
+    for (wy = 0; wy <= Ny / 2; wy++) {
+      float rx = ((wx == 0) || (wx == Nx/2)) ? (1.0/Nx) : (2.0/Nx);
+      float ry = ((wy == 0) || (wy == Ny/2)) ? (1.0/Ny) : (2.0/Ny);
+      float rn = rx*ry;
+      for (ix = 0; ix < Nx; ix++) {
+	float fxc = cos(wx*Wx*ix);
+	for (iy = 0; iy < Ny; iy++) {
+	  float fyc = cos(wy*Wy*iy);
+	  int Nf = wx + wy*Ny;
+	  int out = ix + iy*Nx + Nf*Nx*Ny;
+	  Fterms[out] = fxc*fyc*rn;
+	  if (wy == 0) continue;
+	  if (wy == Ny/2) continue;
+	  float fys = sin(wy*Wy*iy);
+	  Nf = wx + (Ny - wy)*Ny;
+	  out = ix + iy*Nx + Nf*Nx*Ny;
+	  Fterms[out] = fxc*fys*rn;
+	}
+	if (wx == 0) continue;
+	if (wx == Nx/2) continue;
+	float fxs = sin(wx*Wx*ix);
+	for (iy = 0; iy < Ny; iy++) {
+	  float fyc = cos(wy*Wy*iy);
+	  int Nf = (Nx - wx) + wy*Ny;
+	  int out = ix + iy*Nx + Nf*Nx*Ny;
+	  Fterms[out] = fxs*fyc*rn;
+	  if (wy == 0) continue;
+	  if (wy == Ny/2) continue;
+	  float fys = sin(wy*Wy*iy);
+	  Nf = (Nx - wx) + (Ny - wy)*Ny;
+	  out = ix + iy*Nx + Nf*Nx*Ny;
+	  Fterms[out] = fxs*fys*rn;
+	}
+      }
+    }
+  }
+
+ reuse:
+  {
+    // generate and save the dot products
+    float *sv = (float *)src->matrix.buffer;
+    float *tv = (float *)tgt->matrix.buffer;
+    for (wx = 0; wx <= Nx/2; wx++) {
+      for (wy = 0; wy <= Ny/2; wy++) {
+	float fsum = 0.0;
+	int Nf = (wx + wy*Nx)*Nx*Ny;
+	for (ix = 0; ix < Nx; ix++) {
+	  for (iy = 0; iy < Ny; iy++) {
+	    float Ft = Fterms[Nf + ix + iy*Nx];
+	    float Fi = sv[ix + iy*Nx];
+	    fsum += Ft*Fi;
+	    // fprintf (stderr, "%d,%d : %d,%d : %f %f %f\n", ix, iy, wx, wy, Ft, Fi, fsum);
+	  }
+	}
+	tv[wx + wy*Nx] = fsum;
+      }
+    }
+  }
+  return (TRUE);
+
+ usage:
+  gprint (GP_ERR, "USAGE: dft2d (input) to (output)\n");
+  gprint (GP_ERR, "  NOTE: Nx & Ny must both be even\n");
+  return (FALSE);
+}
+
+
+/*** 
+
+     2D direct fourier transform:
+
+     we have an image src of size [Nx,Ny].  in each direction,
+     we can generate the following fourier components:
+
+     cos(0*2pi*x/Nx), cos(1*2pi*x/Nx),... cos((Nx/2)  *2pi*x/Nx), 
+     sin((Nx/2-1)*2pi*x/Nx), sin((Nx/2-2)*2pi*x/Nx),.. sin(1*2pi*x/Nx)
+
+     Note that there are 2 fewer sin elements than there are cos
+     elements:
+       cos(0) = 1 -> DC term
+       sin(0) = 0
+       sin(pi*x) = 0 for integer values of x
+		     
+     we represent the output with the highest frequency in the middle and the DC element
+     in the 0,0 corner
+
+     to generate the output elements, we need to take the dot product 
+     of the input image with an image for the given x,y frequency
+ 
+     thus, for an Nx*Ny input image we need a cube of fourier terms of dimension
+     Nx*Ny*(Nx*Ny) (only one of which is trivial)
+
+     In the assumption that we will likely re-do the same size image multiple times, I
+     store the fourier factor array once it is generated.
+
+     NOTE: do I need to require Nx & Ny even?
+
+ ***/
Index: trunk/Ohana/src/opihi/cmd.data/dimenup.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/dimenup.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/dimenup.c	(revision 37807)
@@ -23,8 +23,5 @@
     return (FALSE);
   }
-
-  gfits_free_matrix (&buf[0].matrix);
-  gfits_free_header (&buf[0].header);
-  CreateBuffer (buf, Nx, Ny, -32, 0.0, 1.0);
+  ResetBuffer (buf, Nx, Ny, -32, 0.0, 1.0);
 
   out = (float *) buf[0].matrix.buffer;
Index: trunk/Ohana/src/opihi/cmd.data/erase.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/erase.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/erase.c	(revision 37807)
@@ -15,4 +15,11 @@
   if (!GetImage (NULL, &kapa, name)) return (FALSE);
   FREE (name);
+
+  int channel = 0;
+  if ((N = get_argument (argc, argv, "-ch"))) {
+    channel = GetKapaChannelFromString (argv[N]);
+    if (!channel) return FALSE;
+    KiiSetChannel (kapa, channel - 1);
+  }
 
   if (argc < 2) {
Index: trunk/Ohana/src/opihi/cmd.data/fit2d.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/fit2d.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/fit2d.c	(revision 37807)
@@ -138,4 +138,6 @@
     y = yvec[0].elements.Flt;
     z = zvec[0].elements.Flt;
+    if (Weight) dz = dzvec[0].elements.Flt;
+
     for (i = 0; i < xvec[0].Nelements; i++, x++, y++, z++) {
       if (mask[i]) continue;
Index: trunk/Ohana/src/opihi/cmd.data/init.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/init.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/init.c	(revision 37807)
@@ -1,3 +1,4 @@
 # include "data.h"
+// XXX adding a comment
 
 int accum            PROTO((int, char **));
@@ -31,4 +32,5 @@
 int densify          PROTO((int, char **));
 int device           PROTO((int, char **));
+int dft2d            PROTO((int, char **));
 int dimendown        PROTO((int, char **));
 int dimenup          PROTO((int, char **));
@@ -147,4 +149,10 @@
 int vpop             PROTO((int, char **));
 int vsmooth          PROTO((int, char **));
+int vsh              PROTO((int, char **));
+int vshfit           PROTO((int, char **));
+int shterms          PROTO((int, char **));
+int shfit            PROTO((int, char **));
+int shdot            PROTO((int, char **));
+int shapply          PROTO((int, char **));
 int wd               PROTO((int, char **));
 int write_vectors    PROTO((int, char **));
@@ -189,4 +197,5 @@
   {1, "densify",      densify,          "create an image histogram from a set of vectors"},
   {1, "device",       device,           "set / get current graphics device"},
+  {1, "dft2d",        dft2d,            "2D discrete fourier transform"},
   {1, "dimendown",    dimendown,        "convert image to vector"},
   {1, "dimenup",      dimenup,          "convert vector to image"},
@@ -314,4 +323,10 @@
   {1, "vstats",       vstats,           "statistics on a vector"},
   {1, "vzload",       vzload,           "load vectors as overlay on image display (scaled points)"},
+  {1, "vsh",          vsh,              "Vector Spherical Harmonics"},
+  {1, "vshfit",       vshfit,           "Vector Spherical Harmonics fits"},
+  {1, "shterms",      shterms,          "Spherical Harmonics terms"},
+  {1, "shfit",        shfit,            "Spherical Harmonics fits"},
+  {1, "shdot",        shdot,            "Spherical Harmonics dot product"},
+  {1, "shapply",      shapply,          "Spherical Harmonics fit application"},
   {1, "wd",           wd,               "write an image to a file"},
   {1, "write",        write_vectors,    "write vectors to datafile"},
Index: trunk/Ohana/src/opihi/cmd.data/load.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/load.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/load.c	(revision 37807)
@@ -23,4 +23,11 @@
   if (!GetImageData (&data, &kapa, name)) return (FALSE);
   FREE (name);
+
+  int channel = 0;
+  if ((N = get_argument (argc, argv, "-ch"))) {
+    channel = GetKapaChannelFromString (argv[N]);
+    if (!channel) return FALSE;
+    KiiSetChannel (kapa, channel - 1);
+  }
 
   ISCEL = FALSE;
Index: trunk/Ohana/src/opihi/cmd.data/point.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/point.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/point.c	(revision 37807)
@@ -21,4 +21,11 @@
   if (!GetImageData (&data, &kapa, name)) return (FALSE);
   FREE (name);
+
+  int channel = 0;
+  if ((N = get_argument (argc, argv, "-ch"))) {
+    channel = GetKapaChannelFromString (argv[N]);
+    if (!channel) return FALSE;
+    KiiSetChannel (kapa, channel - 1);
+  }
 
   celestial = FALSE;
Index: trunk/Ohana/src/opihi/cmd.data/read_vectors.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/read_vectors.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/read_vectors.c	(revision 37807)
@@ -5,4 +5,5 @@
 
 void read_vectors_cleanup ();
+int read_table_sizes (Header *header);
 
 int datafile (int argc, char **argv) {
@@ -331,15 +332,14 @@
 
   off_t Nbytes;
-  int i, j, N, Nextend, Ny, Binary, vecType, padIfShort;
-  char type[16], ID[80], *CCDKeyword;
+  int i, j, N, Ny, Binary, vecType;
+  char type[16], ID[80];
   FTable table;
   Header header;
   Vector **vec;
-  int FITS_TRANSPOSE;
 
   table.buffer = NULL;
   header.buffer = NULL;
 
-  FITS_TRANSPOSE = FALSE;
+  int FITS_TRANSPOSE = FALSE;
   if ((N = get_argument (argc, argv, "-transpose"))) {
     remove_argument (N, &argc, argv);
@@ -347,5 +347,5 @@
   }
 
-  CCDKeyword = NULL;
+  char *CCDKeyword = NULL;
   if ((N = get_argument (argc, argv, "-keyword"))) {
     remove_argument (N, &argc, argv);
@@ -354,5 +354,5 @@
   }
 
-  padIfShort = FALSE;
+  int padIfShort = FALSE;
   if ((N = get_argument (argc, argv, "-pad-if-short"))) {
     remove_argument (N, &argc, argv);
@@ -360,8 +360,33 @@
   }
 
-  Nextend = -1;
+  int getSizes = FALSE;
+  if ((N = get_argument (argc, argv, "-sizes"))) {
+    remove_argument (N, &argc, argv);
+    getSizes = TRUE;
+  }
+
+  int Nextend = -1;
   if ((N = get_argument (argc, argv, "-extnum"))) {
     remove_argument (N, &argc, argv);
     Nextend = atoi (extname);
+  }
+
+  int start = 0;
+  int Nrows = -1; // -1 : read entire table
+  if ((N = get_argument (argc, argv, "-range"))) {
+    remove_argument (N, &argc, argv);
+    start = atoi (argv[N]);
+    remove_argument (N, &argc, argv);
+    Nrows = atoi (argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+
+  // if CharAsVectors, char fields will be saved as vectors NAME:0 -- NAME:n for n characters
+  // else char fields will be saved as $NAME:0 - $NAME:m for m rows
+  // if (Ny > 10000), force CharAsVectors
+  int CharAsVectors = FALSE;
+  if ((N = get_argument (argc, argv, "-char-vectors"))) {
+    remove_argument (N, &argc, argv);
+    CharAsVectors = TRUE;
   }
 
@@ -397,6 +422,24 @@
     }
     if (!gfits_load_header (f, &header)) ESCAPE ("error reading header for extension %d\n", Nextend);
-    if (!gfits_fread_ftable_data (f, &table, padIfShort)) ESCAPE ("error reading table for extension %d\n", Nextend);
-
+
+    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;
+    }
+    if (start < 0) ESCAPE ("invalid range: start < 0\n");
+    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]);
+
+    // Ny = 100, start = 0, Nrows = -1 -> Nrows => 100
+    // Ny = 100, start = 10, Nrows = 90
+
+    if (!gfits_fread_ftable_range (f, padIfShort, &table, start, Nrows)) ESCAPE ("error reading table for extension %d\n", Nextend);
   } else {
     if (CCDKeyword == NULL) {
@@ -428,5 +471,17 @@
 	continue;
       }
-      if (!gfits_fread_ftable_data (f, &table, padIfShort)) ESCAPE ("error reading table for extension\n");
+
+      if (Nrows == -1) {
+	Nrows = header.Naxis[1] - start;
+      }
+      if (start < 0) ESCAPE ("invalid range: start < 0\n");
+      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]);
+
+      if (!gfits_fread_ftable_range (f, padIfShort, &table, start, Nrows)) ESCAPE ("error reading table for extension %d\n", Nextend);
+
+      // if (!gfits_fread_ftable_data (f, &table, padIfShort)) ESCAPE ("error reading table for extension\n");
+
       break;
     }
@@ -465,5 +520,5 @@
     if (!FITS_TRANSPOSE) {
       // read string column into a list rather than a vector
-      if (!strcmp (type, "char")) {
+      if (!strcmp (type, "char") && !CharAsVectors && (Ny < 3000)) {
 	char *fieldName = argv[i];
 	char *Ptr = data;
@@ -537,2 +592,19 @@
   if (vec) free (vec);
 }
+
+// read -fits foo -sizes -- Nx, Ny, Nfields -> $table:Nx, $table:Ny, $table:$Nfields
+// read -fits foo -fields 
+
+int read_table_sizes (Header *header) {
+  
+  int Nfields;
+
+  gfits_scan (header, "TFIELDS", "%d", 1, &Nfields);
+
+  set_int_variable ("table:Nx", header->Naxis[0]);
+  set_int_variable ("table:Nx", header->Naxis[0]);
+  set_int_variable ("table:Nfields", Nfields);
+
+  return TRUE;
+}
+
Index: trunk/Ohana/src/opihi/cmd.data/relocate.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/relocate.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/relocate.c	(revision 37807)
@@ -17,4 +17,11 @@
   FREE (name);
 
+  int channel = 0;
+  if ((N = get_argument (argc, argv, "-ch"))) {
+    channel = GetKapaChannelFromString (argv[N]);
+    if (!channel) return FALSE;
+    KiiSetChannel (kapa, channel - 1);
+  }
+
   if (argc != 3) {
     gprint (GP_ERR, "USAGE: relocate x y [-n]\n");
Index: trunk/Ohana/src/opihi/cmd.data/resize.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/resize.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/resize.c	(revision 37807)
@@ -17,4 +17,11 @@
   if (!GetImage (NULL, &kapa, name)) return (FALSE);
   FREE (name);
+
+  int channel = 0;
+  if ((N = get_argument (argc, argv, "-ch"))) {
+    channel = GetKapaChannelFromString (argv[N]);
+    if (!channel) return FALSE;
+    KiiSetChannel (kapa, channel - 1);
+  }
 
   if ((N = get_argument (argc, argv, "-by-image"))) {
Index: trunk/Ohana/src/opihi/cmd.data/save.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/save.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/save.c	(revision 37807)
@@ -16,4 +16,11 @@
   FREE (name);
 
+  int channel = 0;
+  if ((N = get_argument (argc, argv, "-ch"))) {
+    channel = GetKapaChannelFromString (argv[N]);
+    if (!channel) return FALSE;
+    KiiSetChannel (kapa, channel - 1);
+  }
+
   celestial = FALSE;
   if ((N = get_argument (argc, argv, "-c"))) {
Index: trunk/Ohana/src/opihi/cmd.data/shapply.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/shapply.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/shapply.c	(revision 37807)
@@ -0,0 +1,63 @@
+# include "data.h"
+
+int shapply (int argc, char **argv) {
+  
+  int i, j;
+  Vector *Rvec, *Dvec, *Frvec, *Fivec, *Vrvec, *Vivec;
+
+  if (argc != 8) {
+    gprint (GP_ERR, "USAGE: shapply R D lmax Fr Fi Vr Vi\n");
+    gprint (GP_ERR, "  apply the Ylm values (Fr,Fi) at the given R,D values to get the fit at those positions\n");
+    return (FALSE);
+  }
+
+  if ((Rvec = SelectVector (argv[1], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Dvec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if (Rvec->Nelements != Dvec->Nelements) {
+    gprint (GP_ERR, "R and D sizes do not match\n");
+    return FALSE;
+  }
+
+  int lmax = atoi(argv[3]);
+
+  if ((Frvec = SelectVector (argv[4], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Fivec = SelectVector (argv[5], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if (Frvec->Nelements != Fivec->Nelements) { 
+    gprint (GP_ERR, "Fr and Fi sizes do not match\n");
+    return FALSE;
+  }
+
+  if ((Vrvec = SelectVector (argv[6], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Vivec = SelectVector (argv[7], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+
+  SHterms *terms = SHtermsInit (lmax);
+  if (Frvec->Nelements != terms->Nterms) {
+    gprint (GP_ERR, "Lmax (%d) terms (%d) does not match size of Fr,Fi vectors (%d) \n", lmax, terms->Nterms, Frvec->Nelements);
+    return FALSE;
+  }
+
+  ResetVector (Vrvec, OPIHI_FLT, Rvec->Nelements);
+  ResetVector (Vivec, OPIHI_FLT, Rvec->Nelements);
+
+  // measure the dot product \sum(F_i * Ylm_i)
+  for (i = 0; i < Rvec->Nelements; i++) {
+
+    SHtermsForRD (terms, Rvec->elements.Flt[i], Dvec->elements.Flt[i]); 
+
+    double *Fr = Frvec->elements.Flt;
+    double *Fi = Fivec->elements.Flt;
+
+    double Vr = 0.0;
+    double Vi = 0.0;
+    for (j = 0; j < terms->Nterms; j++) {
+      Vr += Fr[j] * terms->Fr[j];
+      Vi += Fi[j] * terms->Fi[j];
+    }
+    Vrvec->elements.Flt[i] = Vr;
+    Vivec->elements.Flt[i] = Vi;
+  }
+
+  SHtermsFree (terms);
+
+  return (TRUE);
+}
Index: trunk/Ohana/src/opihi/cmd.data/shdot.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/shdot.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/shdot.c	(revision 37807)
@@ -0,0 +1,69 @@
+# include "data.h"
+
+int shdot (int argc, char **argv) {
+  
+  int i, j;
+  Vector *Rvec, *Dvec, *Fvec, *Lvec, *Mvec, *Frvec, *Fivec;
+  double *Fr, *Fi;
+
+  if (argc != 9) {
+    gprint (GP_ERR, "USAGE: shdot R D value lmax Lout Mout Fr Fi\n");
+    gprint (GP_ERR, "  find the dot product of the scalar field value (at points R,D) to spherical harmonics up to the given lmax\n");
+    return (FALSE);
+  }
+
+  if ((Rvec = SelectVector (argv[1], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Dvec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Fvec = SelectVector (argv[3], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+
+  int lmax = atoi(argv[4]);
+
+  if ((Lvec  = SelectVector (argv[5], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Mvec  = SelectVector (argv[6], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Frvec = SelectVector (argv[7], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Fivec = SelectVector (argv[8], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+
+  SHterms *terms = SHtermsInit (lmax);
+
+  ALLOCATE_ZERO (Fr, double, terms->Nterms);
+  ALLOCATE_ZERO (Fi, double, terms->Nterms);
+
+  ResetVector (Frvec, OPIHI_FLT, terms->Nterms);
+  ResetVector (Fivec, OPIHI_FLT, terms->Nterms);
+  ResetVector (Lvec,  OPIHI_FLT, terms->Nterms);
+  ResetVector (Mvec,  OPIHI_FLT, terms->Nterms);
+
+  // measure the dot product \sum(F_i * Ylm_i)
+  for (i = 0; i < Rvec->Nelements; i++) {
+
+    SHtermsForRD (terms, Rvec->elements.Flt[i], Dvec->elements.Flt[i]); 
+
+    double Fv = Fvec->elements.Flt[i];
+
+    for (j = 0; j < terms->Nterms; j++) {
+      Fr[j] += Fv * terms->Fr[j];
+      Fi[j] += Fv * terms->Fi[j];
+    }
+  }
+
+  // for (j = 0; j < terms->Nterms; j++) {
+  //   fprintf (stderr, "%d : %d %d : %f %f\n", j, terms->l[j], terms->m[j], Fr[j], Fi[j]);
+  // }
+
+  for (j = 0; j < terms->Nterms; j++) {
+    Fr[j] /= Fvec->Nelements;
+    Fi[j] /= Fvec->Nelements;
+  }
+
+  for (j = 0; j < terms->Nterms; j++) {
+    // fprintf (stderr, "%d : %d %d : %f %f\n", j, terms->l[j], terms->m[j], 4*M_PI*Fr[j], 4*M_PI*Fi[j]);
+    Lvec[0].elements.Flt[j]  = terms->l[j];
+    Mvec[0].elements.Flt[j]  = terms->m[j];
+    Frvec[0].elements.Flt[j] = 4*M_PI*Fr[j];
+    Fivec[0].elements.Flt[j] = 4*M_PI*Fi[j];
+  }
+
+  SHtermsFree (terms);
+
+  return (TRUE);
+}
Index: trunk/Ohana/src/opihi/cmd.data/shfit.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/shfit.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/shfit.c	(revision 37807)
@@ -0,0 +1,163 @@
+# include "data.h"
+
+int shfit (int argc, char **argv) {
+  
+  int i, j, k;
+  Vector *Rvec, *Dvec, *Fvec, *Lvec, *Mvec, *Frvec, *Fivec;
+  double *Fr, *Fi;
+
+  if (argc != 9) {
+    gprint (GP_ERR, "USAGE: shfit R D value lmax Lout Mout Fr Fi\n");
+    gprint (GP_ERR, "  fit the scalar field value (at points R,D) to spherical harmonics up to the given lmax\n");
+    return (FALSE);
+  }
+
+  if ((Rvec = SelectVector (argv[1], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Dvec = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Fvec = SelectVector (argv[3], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+
+  int lmax = atoi(argv[4]);
+
+  if ((Lvec  = SelectVector (argv[5], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Mvec  = SelectVector (argv[6], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Frvec = SelectVector (argv[7], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Fivec = SelectVector (argv[8], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+
+  SHterms *terms = SHtermsInit (lmax);
+
+  ALLOCATE_ZERO (Fr, double, terms->Nterms);
+  ALLOCATE_ZERO (Fi, double, terms->Nterms);
+
+  ResetVector (Frvec, OPIHI_FLT, terms->Nterms);
+  ResetVector (Fivec, OPIHI_FLT, terms->Nterms);
+  ResetVector (Lvec,  OPIHI_FLT, terms->Nterms);
+  ResetVector (Mvec,  OPIHI_FLT, terms->Nterms);
+
+  // we only fit the linearlly independent terms: Re(m >= 0), Im(m > 0)
+  int Nre = 0;
+  int Nim = 0;
+  for (i = 0; i < terms->Nterms; i++) {
+    if (terms->m[i] >= 0) Nre ++;
+    if (terms->m[i] >  0) Nim ++;
+  }
+  int *Jre = NULL;
+  int *Jim = NULL;
+  ALLOCATE (Jre, int, Nre);
+  ALLOCATE (Jim, int, Nim);
+
+  Nre = Nim = 0;
+  for (i = 0; i < terms->Nterms; i++) {
+    if (terms->m[i] >= 0) {
+      Jre[Nre] = i;
+      Nre ++;
+    }
+    if (terms->m[i] >  0) {
+      Jim[Nim] = i;
+      Nim ++;
+    }
+  }
+
+  double **Are, **bre, **Aim, **bim;
+  ALLOCATE (Are, double *, Nre);
+  ALLOCATE (bre, double *, Nre);
+  ALLOCATE (Aim, double *, Nim);
+  ALLOCATE (bim, double *, Nim);
+  for (i = 0; i < Nre; i++) {
+    ALLOCATE_ZERO (Are[i], double, Nre);
+    ALLOCATE_ZERO (bre[i], double, 1);
+  }
+  for (i = 0; i < Nim; i++) {
+    ALLOCATE_ZERO (Aim[i], double, Nim);
+    ALLOCATE_ZERO (bim[i], double, 1);
+  }
+
+  // measure the dot product \sum(F_i * Ylm_i) and the cross terms (\sum(Y_lm * Y_jk))
+  for (i = 0; i < Rvec->Nelements; i++) {
+
+    SHtermsForRD (terms, Rvec->elements.Flt[i], Dvec->elements.Flt[i]); 
+
+    double Fv = Fvec->elements.Flt[i];
+
+    for (j = 0; j < Nre; j++) {
+      int jre = Jre[j];
+      bre[j][0] += Fv * terms->Fr[jre];
+    }
+    for (j = 0; j < Nim; j++) {
+      int jim = Jim[j];
+      bim[j][0] += Fv * terms->Fi[jim];
+    }
+
+    for (j = 0; j < Nre; j++) {
+      int jre = Jre[j];
+      for (k = j; k < Nre; k++) {
+	int kre = Jre[k];
+	Are[j][k] += terms->Fr[jre] * terms->Fr[kre];
+      }
+    }
+
+    for (j = 0; j < Nim; j++) {
+      int jim = Jim[j];
+      for (k = j; k < Nim; k++) {
+	int kim = Jim[k];
+	Aim[j][k] += terms->Fi[jim] * terms->Fi[kim];
+      }
+    }
+  }
+
+  for (j = 1; j < Nre; j++) {
+    for (k = 0; k < j; k++) {
+      Are[j][k] = Are[k][j];
+    }	
+  }
+  for (j = 1; j < Nim; j++) {
+    for (k = 0; k < j; k++) {
+      Aim[j][k] = Aim[k][j];
+    }	
+  }
+
+  fprintf (stderr, "--- Are --- : bre \n");
+  for (j = 0; j < Nre; j++) {
+    fprintf (stderr, "%10.6f ", Are[j][j]);
+    fprintf (stderr, " : %10.6f\n", bre[j][0]);
+  }
+
+  fprintf (stderr, "--- Aim --- : bim \n");
+  for (i = 0; i < Nim; i++) {
+    fprintf (stderr, "%10.6f : ", Aim[i][i]);
+    fprintf (stderr, " : %10.6f\n", bim[i][0]);
+  }
+
+  if (!dgaussjordan (Are, bre, Nre, 1)) {
+    gprint (GP_ERR, "failed to fit data : ill-conditioned matrix\n");
+    return FALSE;
+  }
+  if (!dgaussjordan (Aim, bim, Nim, 1)) {
+    gprint (GP_ERR, "failed to fit data : ill-conditioned matrix\n");
+    return FALSE;
+  }
+
+  for (j = 0; j < terms->Nterms; j++) { 
+    Fr[j] = 0.0;
+    Fi[i] = 0.0;
+  }
+  for (j = 0; j < Nre; j++) {
+    int jre = Jre[j];
+    Fr[jre] = bre[j][0];
+  }
+  for (j = 0; j < Nim; j++) {
+    int jim = Jim[j];
+    Fi[jim] = bim[j][0];
+  }
+
+  for (j = 0; j < terms->Nterms; j++) {
+    // fprintf (stderr, "%d : %d %d : %f %f\n", j, terms->l[j], terms->m[j], 4*M_PI*Fr[j], 4*M_PI*Fi[j]);
+    Lvec[0].elements.Flt[j]  = terms->l[j];
+    Mvec[0].elements.Flt[j]  = terms->m[j];
+    Frvec[0].elements.Flt[j] = Fr[j];
+    Fivec[0].elements.Flt[j] = Fi[j];
+  }
+
+  SHtermsFree (terms);
+
+  return (TRUE);
+}
Index: trunk/Ohana/src/opihi/cmd.data/shterms.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/shterms.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/shterms.c	(revision 37807)
@@ -0,0 +1,43 @@
+# include "data.h"
+
+int shterms (int argc, char **argv) {
+  
+  int i;
+  Vector *Frvec, *Fivec, *lvec, *mvec;
+
+  if (argc != 8) {
+    gprint (GP_ERR, "USAGE: shterms Fr Fi l m lmax R D\n");
+    gprint (GP_ERR, "  set the Fr and Fi values for sh terms up to lmax at R,D\n");
+    return (FALSE);
+  }
+
+  if ((Frvec = SelectVector (argv[1], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Fivec = SelectVector (argv[2], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((lvec  = SelectVector (argv[3], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((mvec  = SelectVector (argv[4], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+
+  int lmax = atoi(argv[5]);
+
+  SHterms *terms = SHtermsInit (lmax);
+
+  double R = atof(argv[6]);
+  double D = atof(argv[7]);
+
+  ResetVector (Frvec, OPIHI_FLT, terms->Nterms);
+  ResetVector (Fivec, OPIHI_FLT, terms->Nterms);
+  ResetVector (lvec,  OPIHI_FLT, terms->Nterms);
+  ResetVector (mvec,  OPIHI_FLT, terms->Nterms);
+
+  SHtermsForRD (terms, R, D); 
+
+  for (i = 0; i < terms->Nterms; i++) {
+    Frvec[0].elements.Flt[i] = terms->Fr[i];
+    Fivec[0].elements.Flt[i] = terms->Fi[i];
+    lvec[0].elements.Flt[i]  = terms->l[i];
+    mvec[0].elements.Flt[i]  = terms->m[i];
+  }
+
+  SHtermsFree (terms);
+
+  return (TRUE);
+}
Index: trunk/Ohana/src/opihi/cmd.data/tv.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/tv.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/tv.c	(revision 37807)
@@ -18,4 +18,13 @@
   if (!GetImage (&data, &kapa, name)) return (FALSE);
   FREE (name);
+
+  int channel = 0;
+  if ((N = get_argument (argc, argv, "-ch"))) {
+    remove_argument (N, &argc, argv);
+    channel = GetKapaChannelFromString (argv[N]);
+    remove_argument (N, &argc, argv);
+    if (!channel) return FALSE;
+    KiiSetChannel (kapa, channel - 1);
+  }
 
   /* shell exits on pipe close, FIX */
Index: trunk/Ohana/src/opihi/cmd.data/tvchannel.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/tvchannel.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/tvchannel.c	(revision 37807)
@@ -22,23 +22,31 @@
   }
 
-  Nchannel = atoi (argv[1]);
+  Nchannel = GetKapaChannelFromString (argv[1]);
+  if (!Nchannel) return FALSE;
+
+  KiiSetChannel (kapa, Nchannel - 1);
+  return (TRUE);
+}
+
+int GetKapaChannelFromString (char *string) {
+
+  int Nchannel = atoi (string);
   if (Nchannel == 0) {
     // try the string values R/Red, G/Green, B/Blue
-    if (!strcasecmp (argv[1], "R") || !strcasecmp (argv[1], "RED")) {
+    if (!strcasecmp (string, "R") || !strcasecmp (string, "RED")) {
       Nchannel = 1;
     }
-    if (!strcasecmp (argv[1], "G") || !strcasecmp (argv[1], "GREEN")) {
+    if (!strcasecmp (string, "G") || !strcasecmp (string, "GREEN")) {
       Nchannel = 2;
     }
-    if (!strcasecmp (argv[1], "B") || !strcasecmp (argv[1], "BLUE")) {
+    if (!strcasecmp (string, "B") || !strcasecmp (string, "BLUE")) {
       Nchannel = 3;
     }
   }
-  if ((Nchannel < 1) || (Nchannel > 3)) {
-    gprint (GP_ERR, "invalid channel : use 1 - 3 or (R)ed, (G)reen, (B)lue\n");
-    return (FALSE);
+  if ((Nchannel < 1) || (Nchannel > 10)) {
+    gprint (GP_ERR, "invalid channel : use 1 - 10 or (R)ed, (G)reen, (B)lue\n");
+    gprint (GP_ERR, "   (R)ed, (G)reen, (B)lue == (1,2,3)\n");
+    return (0);
   }
-    
-  KiiSetChannel (kapa, Nchannel - 1);
-  return (TRUE);
+  return Nchannel;
 }
Index: trunk/Ohana/src/opihi/cmd.data/tvcontour.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/tvcontour.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/tvcontour.c	(revision 37807)
@@ -71,4 +71,11 @@
   if (!GetImage (NULL, &kapa, name)) return (FALSE);
   FREE (name);
+
+  int channel = 0;
+  if ((N = get_argument (argc, argv, "-ch"))) {
+    channel = GetKapaChannelFromString (argv[N]);
+    if (!channel) return FALSE;
+    KiiSetChannel (kapa, channel - 1);
+  }
 
   if ((argc != 4) && (argc != 5)) {
Index: trunk/Ohana/src/opihi/cmd.data/tvgrid.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/tvgrid.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/tvgrid.c	(revision 37807)
@@ -24,4 +24,11 @@
   if (!GetImage (NULL, &kapa, name)) return (FALSE);
   FREE (name);
+
+  int channel = 0;
+  if ((N = get_argument (argc, argv, "-ch"))) {
+    channel = GetKapaChannelFromString (argv[N]);
+    if (!channel) return FALSE;
+    KiiSetChannel (kapa, channel - 1);
+  }
 
   if (argc != 3) {
Index: trunk/Ohana/src/opihi/cmd.data/vload.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/vload.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/vload.c	(revision 37807)
@@ -18,4 +18,11 @@
   if (!GetImage (NULL, &kapa, name)) return (FALSE);
   FREE (name);
+
+  int channel = 0;
+  if ((N = get_argument (argc, argv, "-ch"))) {
+    channel = GetKapaChannelFromString (argv[N]);
+    if (!channel) return FALSE;
+    KiiSetChannel (kapa, channel - 1);
+  }
 
   type = KII_OVERLAY_BOX;
Index: trunk/Ohana/src/opihi/cmd.data/vsh.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/vsh.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/vsh.c	(revision 37807)
@@ -0,0 +1,43 @@
+# include "data.h"
+
+int vsh (int argc, char **argv) {
+  
+  int i;
+  Vector *Rbvec, *Revec, *Dbvec, *Devec;
+
+  if (argc != 8) {
+    gprint (GP_ERR, "USAGE: vsh dRb dRe dDb dDe lmax R D\n");
+    gprint (GP_ERR, "  set the dRb, dRe, dDb, dDe values for vsh terms up to lmax at R,D\n");
+    return (FALSE);
+  }
+
+  if ((Rbvec = SelectVector (argv[1], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Revec = SelectVector (argv[2], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Dbvec = SelectVector (argv[3], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Devec = SelectVector (argv[4], ANYVECTOR, TRUE)) == NULL) return (FALSE);
+
+  int lmax = atoi(argv[5]);
+
+  VSHterms *terms = VSHtermsInit (lmax);
+
+  double R = atof(argv[6]);
+  double D = atof(argv[7]);
+
+  ResetVector (Rbvec, OPIHI_FLT, terms->Nterms);
+  ResetVector (Revec, OPIHI_FLT, terms->Nterms);
+  ResetVector (Dbvec, OPIHI_FLT, terms->Nterms);
+  ResetVector (Devec, OPIHI_FLT, terms->Nterms);
+
+  VSHtermsForRD (terms, R, D); 
+
+  for (i = 0; i < terms->Nterms; i++) {
+    Rbvec[0].elements.Flt[i] = terms->dR_B[i];
+    Revec[0].elements.Flt[i] = terms->dR_E[i];
+    Dbvec[0].elements.Flt[i] = terms->dD_B[i];
+    Devec[0].elements.Flt[i] = terms->dD_E[i];
+  }
+
+  VSHtermsFree (terms);
+
+  return (TRUE);
+}
Index: trunk/Ohana/src/opihi/cmd.data/vshfit.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/vshfit.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.data/vshfit.c	(revision 37807)
@@ -0,0 +1,62 @@
+# include "data.h"
+
+int vshfit (int argc, char **argv) {
+  
+  int i, j;
+  Vector *Rvec, *Dvec, *dRvec, *dDvec;
+  double *Re, *Rb, *De, *Db;
+
+  if (argc != 6) {
+    gprint (GP_ERR, "USAGE: vshfit Rvec Dvec dRvec dDvec lmax\n");
+    gprint (GP_ERR, "  fit the vector field dRvec, dDvec (at points Rvec,Dvec) to vector spherical harmonics up to the given lmax\n");
+    return (FALSE);
+  }
+
+  if ((Rvec  = SelectVector (argv[1], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((Dvec  = SelectVector (argv[2], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((dRvec = SelectVector (argv[3], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+  if ((dDvec = SelectVector (argv[4], OLDVECTOR, TRUE)) == NULL) return (FALSE);
+
+  int lmax = atoi(argv[5]);
+
+  VSHterms *terms = VSHtermsInit (lmax);
+
+  ALLOCATE_ZERO (Re, double, terms->Nterms);
+  ALLOCATE_ZERO (Rb, double, terms->Nterms);
+  ALLOCATE_ZERO (De, double, terms->Nterms);
+  ALLOCATE_ZERO (Db, double, terms->Nterms);
+
+  for (i = 0; i < Rvec->Nelements; i++) {
+
+    VSHtermsForRD (terms, Rvec->elements.Flt[i], Dvec->elements.Flt[i]); 
+
+    double dR = dRvec->elements.Flt[i];
+    double dD = dDvec->elements.Flt[i];
+
+    for (j = 0; j < terms->Nterms; j++) {
+      Rb[j] += dR * terms->dR_B[j];
+      Re[j] += dR * terms->dR_E[j];
+      Db[j] += dD * terms->dD_B[j];
+      De[j] += dD * terms->dD_E[j];
+    }
+  }
+
+  for (j = 0; j < terms->Nterms; j++) {
+    fprintf (stderr, "%d : %d %d : %f %f %f %f\n", j, terms->l[j], terms->m[j], Rb[j], Re[j], Db[j], De[j]);
+  }
+
+  for (j = 0; j < terms->Nterms; j++) {
+    Rb[j] /= terms->Nterms;
+    Re[j] /= terms->Nterms;
+    Db[j] /= terms->Nterms;
+    De[j] /= terms->Nterms;
+  }
+
+  for (j = 0; j < terms->Nterms; j++) {
+    fprintf (stderr, "%d : %d %d : %f %f %f %f\n", j, terms->l[j], terms->m[j], Rb[j], Re[j], Db[j], De[j]);
+  }
+
+  VSHtermsFree (terms);
+
+  return (TRUE);
+}
Index: trunk/Ohana/src/opihi/cmd.data/vzload.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.data/vzload.c	(revision 37049)
+++ trunk/Ohana/src/opihi/cmd.data/vzload.c	(revision 37807)
@@ -18,4 +18,11 @@
   if (!GetImage (NULL, &kapa, name)) return (FALSE);
   FREE (name);
+
+  int channel = 0;
+  if ((N = get_argument (argc, argv, "-ch"))) {
+    channel = GetKapaChannelFromString (argv[N]);
+    if (!channel) return FALSE;
+    KiiSetChannel (kapa, channel - 1);
+  }
 
   MAX_OUTPUT_SIZE = 10.0;
