Index: trunk/Ohana/src/opihi/cmd.astro/Makefile
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/Makefile	(revision 37246)
+++ trunk/Ohana/src/opihi/cmd.astro/Makefile	(revision 37807)
@@ -71,4 +71,6 @@
 $(SRC)/star.$(ARCH).o		   \
 $(SRC)/transform.$(ARCH).o         \
+$(SRC)/vshimage.$(ARCH).o         \
+$(SRC)/shimage.$(ARCH).o         \
 $(SRC)/imsub.$(ARCH).o		   \
 $(SRC)/imfit.$(ARCH).o		   \
Index: trunk/Ohana/src/opihi/cmd.astro/cdensify.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/cdensify.c	(revision 37246)
+++ trunk/Ohana/src/opihi/cmd.astro/cdensify.c	(revision 37807)
@@ -97,11 +97,5 @@
   // generate the PSF in a local tangent plane
   Coords coords;
-  coords.crpix1 = coords.crpix2 = 0.0;
-  coords.crval1 = coords.crval2 = 0.0;
-  coords.cdelt1 = coords.cdelt2 = 1.0;
-  coords.pc1_1  = coords.pc2_2  = 1.0;
-  coords.pc1_2  = coords.pc2_1  = 0.0;
-  coords.Npolyterms = 0;
-  strcpy (coords.ctype, "DEC--TAN");
+  InitCoords (&coords, "DEC--TAN");
 
   r = vr[0].elements.Flt;
Index: trunk/Ohana/src/opihi/cmd.astro/coords.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/coords.c	(revision 37246)
+++ trunk/Ohana/src/opihi/cmd.astro/coords.c	(revision 37807)
@@ -71,5 +71,5 @@
     if ((mosbuffer = SelectBuffer (MOSAIC, OLDBUFFER, TRUE)) == NULL) goto escape;
     GetCoords (&moscoords, &mosbuffer[0].header);
-    RegisterMosaic (&moscoords);
+    coords.mosaic = &moscoords;
   }
   
Index: trunk/Ohana/src/opihi/cmd.astro/csystem.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/csystem.c	(revision 37246)
+++ trunk/Ohana/src/opihi/cmd.astro/csystem.c	(revision 37807)
@@ -13,16 +13,24 @@
   if (argc != 5) goto syntax;
 
-  switch (argv[1][0]) {
-    case 'C': input = COORD_CELESTIAL; break;
-    case 'G': input = COORD_GALACTIC; break;
-    case 'E': input = COORD_ECLIPTIC; break;
-    default: goto syntax;
+  if (!strcmp(argv[1], "G2004")) {
+    input = COORD_GALACTIC_REID_2004;
+  } else {
+    switch (argv[1][0]) {
+      case 'C': input = COORD_CELESTIAL; break;
+      case 'G': input = COORD_GALACTIC; break;
+      case 'E': input = COORD_ECLIPTIC; break;
+      default: goto syntax;
+    }
   }
 
-  switch (argv[2][0]) {
-    case 'C': output = COORD_CELESTIAL; break;
-    case 'G': output = COORD_GALACTIC; break;
-    case 'E': output = COORD_ECLIPTIC; break;
-    default: goto syntax;
+  if (!strcmp(argv[2], "G2004")) {
+    output = COORD_GALACTIC_REID_2004;
+  } else {
+    switch (argv[2][0]) {
+      case 'C': output = COORD_CELESTIAL; break;
+      case 'G': output = COORD_GALACTIC; break;
+      case 'E': output = COORD_ECLIPTIC; break;
+      default: goto syntax;
+    }
   }
 
Index: trunk/Ohana/src/opihi/cmd.astro/fitplx.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/fitplx.c	(revision 37246)
+++ trunk/Ohana/src/opihi/cmd.astro/fitplx.c	(revision 37807)
@@ -104,13 +104,8 @@
   /* project coordinates to a plane centered on the object with units of arcsec */
   Coords coords;
+  InitCoords (&coords, "DEC--SIN");
   coords.crval1 = Rmean;
   coords.crval2 = Dmean;
-  coords.crpix1 = 0;
-  coords.crpix2 = 0;
   coords.cdelt1 = coords.cdelt2 = 1.0 / 3600.0;
-  coords.pc1_1  = coords.pc2_2 = 1.0;
-  coords.pc1_2  = coords.pc2_1 = 0.0;
-  coords.Npolyterms = 1;
-  strcpy (coords.ctype, "DEC--SIN");
 
   double *X, *Y, *t, *pX, *pY, *dX, *dY;
Index: trunk/Ohana/src/opihi/cmd.astro/fitpm.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/fitpm.c	(revision 37246)
+++ trunk/Ohana/src/opihi/cmd.astro/fitpm.c	(revision 37807)
@@ -100,13 +100,8 @@
   /* project coordinates to a plane centered on the object with units of arcsec */
   Coords coords;
+  InitCoords (&coords, "DEC--SIN");
   coords.crval1 = Rmean;
   coords.crval2 = Dmean;
-  coords.crpix1 = 0;
-  coords.crpix2 = 0;
   coords.cdelt1 = coords.cdelt2 = 1.0 / 3600.0;
-  coords.pc1_1  = coords.pc2_2 = 1.0;
-  coords.pc1_2  = coords.pc2_1 = 0.0;
-  coords.Npolyterms = 1;
-  strcpy (coords.ctype, "DEC--SIN");
 
   double *X, *Y, *t, *dX, *dY;
Index: trunk/Ohana/src/opihi/cmd.astro/gauss.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/gauss.c	(revision 37246)
+++ trunk/Ohana/src/opihi/cmd.astro/gauss.c	(revision 37807)
@@ -30,4 +30,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);
+  }
 
   Nborder = 3;
Index: trunk/Ohana/src/opihi/cmd.astro/getcoords.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/getcoords.c	(revision 37246)
+++ trunk/Ohana/src/opihi/cmd.astro/getcoords.c	(revision 37807)
@@ -54,5 +54,5 @@
   if (!strcmp(&coords.ctype[4], "-WRP")) {
     if (MOSAIC == NULL) {
-      gprint (GP_ERR, "must supply mosaic for WRP coords\n");
+      gprint (GP_ERR, "must supply mosaic for WRP coords with -mosaic [buffer]\n");
       return (FALSE);
     }
@@ -65,5 +65,5 @@
     }
     GetCoords (&moscoords, &mosbuffer[0].header);
-    RegisterMosaic (&moscoords);
+    coords.mosaic = &moscoords;
   }
   
Index: trunk/Ohana/src/opihi/cmd.astro/init.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/init.c	(revision 37246)
+++ trunk/Ohana/src/opihi/cmd.astro/init.c	(revision 37807)
@@ -61,4 +61,6 @@
 int times                   PROTO((int, char **));
 int transform               PROTO((int, char **));
+int vshimage                PROTO((int, char **));
+int shimage                 PROTO((int, char **));
 
 static Command cmds[] = {  
@@ -121,4 +123,6 @@
   {1, "star",        star,         "star stats at rough coords"},
   {1, "transform",   transform,    "geometric transformation of image"},
+  {1, "vshimage",    vshimage,     "generate images for vector spherical harmonic terms"},
+  {1, "shimage",     shimage,      "generate images for spherical harmonic terms"},
 }; 
 
Index: trunk/Ohana/src/opihi/cmd.astro/objload.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/objload.c	(revision 37246)
+++ trunk/Ohana/src/opihi/cmd.astro/objload.c	(revision 37807)
@@ -19,4 +19,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);
+  }
 
   Objtype = 0;
Index: trunk/Ohana/src/opihi/cmd.astro/region.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/region.c	(revision 37246)
+++ trunk/Ohana/src/opihi/cmd.astro/region.c	(revision 37807)
@@ -17,4 +17,11 @@
   if (!GetGraph (&graphmode, &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);
+  }
 
   NoClear = FALSE;
@@ -80,5 +87,6 @@
   if (!ohana_str_to_radec (&Ra, &Dec, argv[1], argv[2])) return (FALSE);
   Radius = atof (argv[3]);
-  strcpy (graphmode.coords.ctype, "DEC--TAN");
+  InitCoords (&graphmode.coords, "DEC--TAN");
+
   if (argc == 5) {
     if (!strcasecmp (argv[4], "TAN")) 
@@ -99,6 +107,11 @@
       strcpy (graphmode.coords.ctype, "DEC--PAR");
   }
+  
+  graphmode.coords.crval1 = Ra;
+  graphmode.coords.crval2 = Dec;
 
-  
+  graphmode.coords.pc1_1 = (graphmode.flipeast) ? -1 : 1;
+  graphmode.coords.pc2_2 = (graphmode.flipnorth) ? -1 : 1;
+
   /* ask kapa for coordinate limits, to get the right aspect ratio */
   KapaGetLimits (kapa, &dx, &dy);
@@ -132,14 +145,4 @@
   set_int_variable ("NORTH_UP", !graphmode.flipnorth);
 
-  graphmode.coords.pc1_1 = (graphmode.flipeast) ? -1 : 1;
-  graphmode.coords.pc2_2 = (graphmode.flipnorth) ? -1 : 1;
-
-  graphmode.coords.pc1_2 = graphmode.coords.pc2_1 = 0.0;
-  graphmode.coords.crval1 = Ra;
-  graphmode.coords.crval2 = Dec;
-  graphmode.coords.crpix1 = 0.0;
-  graphmode.coords.crpix2 = 0.0;
-  graphmode.coords.cdelt1 = graphmode.coords.cdelt2 = 1.0;
-
   if (!NoClear) KapaClearSections (kapa);
   KapaSetLimits (kapa, &graphmode);
Index: trunk/Ohana/src/opihi/cmd.astro/shimage.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/shimage.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.astro/shimage.c	(revision 37807)
@@ -0,0 +1,64 @@
+# include "data.h"
+
+int shimage (int argc, char **argv) {
+  
+  int ix, iy;
+  Buffer *FrIm, *FiIm, *src;
+
+  if (argc != 6) {
+    gprint (GP_ERR, "USAGE: shimage (image) l m Fr Fi\n");
+    gprint (GP_ERR, "  set the Fr and Fi values for pixels in the given image for (l,m)\n");
+    return (FALSE);
+  }
+
+  if ((src = SelectBuffer (argv[1], OLDBUFFER, TRUE)) == NULL) return (FALSE);
+
+  Coords coords;
+  GetCoords (&coords, &src[0].header);
+  if (!strcmp(&coords.ctype[4], "-WRP")) {
+    fprintf (stderr, "mosaic astrometry not yet supported for vshimage\n");
+    return FALSE;
+  }
+
+  int l = atoi(argv[2]);
+  int m = atoi(argv[3]);
+
+  if ((FrIm = SelectBuffer (argv[4], ANYBUFFER, TRUE)) == NULL) return (FALSE);
+  if ((FiIm = SelectBuffer (argv[5], ANYBUFFER, TRUE)) == NULL) return (FALSE);
+
+  int Nx = src[0].header.Naxis[0];
+  int Ny = src[0].header.Naxis[1];
+
+  gfits_free_matrix (&FrIm[0].matrix); gfits_free_header (&FrIm[0].header); CreateBuffer (FrIm, Nx, Ny, -32, 1.0, 0.0);
+  gfits_free_matrix (&FiIm[0].matrix); gfits_free_header (&FiIm[0].header); CreateBuffer (FiIm, Nx, Ny, -32, 1.0, 0.0);
+
+  // l=0 allocates space for a single value
+  SHterms *terms = SHtermsInit (0);
+
+  float *Fr = (float *) FrIm[0].matrix.buffer;
+  float *Fi = (float *) FiIm[0].matrix.buffer;
+
+  for (ix = 0; ix < Nx; ix++) {
+    for (iy = 0; iy < Ny; iy++) {
+
+      int Npix = ix + iy*Nx;
+
+      double R, D;
+      int status = XY_to_RD (&R, &D, ix, iy, &coords);
+
+      if (!status) {
+	Fr[Npix] = NAN;
+	Fi[Npix] = NAN;
+	continue;
+      }
+      
+      SHtermsForLM (terms, R, D, l, m); 
+
+      Fr[Npix] = terms->Fr[0];
+      Fi[Npix] = terms->Fi[0];
+    }
+  }
+  SHtermsFree (terms);
+
+  return (TRUE);
+}
Index: trunk/Ohana/src/opihi/cmd.astro/vshimage.c
===================================================================
--- trunk/Ohana/src/opihi/cmd.astro/vshimage.c	(revision 37807)
+++ trunk/Ohana/src/opihi/cmd.astro/vshimage.c	(revision 37807)
@@ -0,0 +1,74 @@
+# include "data.h"
+
+int vshimage (int argc, char **argv) {
+  
+  int ix, iy;
+  Buffer *Rbim, *Reim, *Dbim, *Deim, *src;
+
+  if (argc != 8) {
+    gprint (GP_ERR, "USAGE: vsh image l m dRb dRe dDb dDe\n");
+    gprint (GP_ERR, "  set the dRb, dRe, dDb, dDe values for pixels in the given image for (l,m)\n");
+    return (FALSE);
+  }
+
+  if ((src = SelectBuffer (argv[1], OLDBUFFER, TRUE)) == NULL) return (FALSE);
+
+  Coords coords;
+  GetCoords (&coords, &src[0].header);
+  if (!strcmp(&coords.ctype[4], "-WRP")) {
+    fprintf (stderr, "mosaic astrometry not yet supported for vshimage\n");
+    return FALSE;
+  }
+
+  int l = atoi(argv[2]);
+  int m = atoi(argv[3]);
+
+  if ((Rbim = SelectBuffer (argv[4], ANYBUFFER, TRUE)) == NULL) return (FALSE);
+  if ((Reim = SelectBuffer (argv[5], ANYBUFFER, TRUE)) == NULL) return (FALSE);
+  if ((Dbim = SelectBuffer (argv[6], ANYBUFFER, TRUE)) == NULL) return (FALSE);
+  if ((Deim = SelectBuffer (argv[7], ANYBUFFER, TRUE)) == NULL) return (FALSE);
+
+  int Nx = src[0].header.Naxis[0];
+  int Ny = src[0].header.Naxis[1];
+
+  gfits_free_matrix (&Rbim[0].matrix); gfits_free_header (&Rbim[0].header); CreateBuffer (Rbim, Nx, Ny, -32, 1.0, 0.0);
+  gfits_free_matrix (&Reim[0].matrix); gfits_free_header (&Reim[0].header); CreateBuffer (Reim, Nx, Ny, -32, 1.0, 0.0);
+  gfits_free_matrix (&Dbim[0].matrix); gfits_free_header (&Dbim[0].header); CreateBuffer (Dbim, Nx, Ny, -32, 1.0, 0.0);
+  gfits_free_matrix (&Deim[0].matrix); gfits_free_header (&Deim[0].header); CreateBuffer (Deim, Nx, Ny, -32, 1.0, 0.0);
+
+  // l=0 allocates space for a single value
+  VSHterms *terms = VSHtermsInit (0);
+
+  float *Rb = (float *) Rbim[0].matrix.buffer;
+  float *Re = (float *) Reim[0].matrix.buffer;
+  float *Db = (float *) Dbim[0].matrix.buffer;
+  float *De = (float *) Deim[0].matrix.buffer;
+
+  for (ix = 0; ix < Nx; ix++) {
+    for (iy = 0; iy < Ny; iy++) {
+
+      int Npix = ix + iy*Nx;
+
+      double R, D;
+      int status = XY_to_RD (&R, &D, ix, iy, &coords);
+
+      if (!status) {
+	Rb[Npix] = NAN;
+	Re[Npix] = NAN;
+	Db[Npix] = NAN;
+	De[Npix] = NAN;
+	continue;
+      }
+      
+      VSHtermsForLM (terms, R, D, l, m); 
+
+      Rb[Npix] = terms->dR_B[0];
+      Re[Npix] = terms->dR_E[0];
+      Db[Npix] = terms->dD_B[0];
+      De[Npix] = terms->dD_E[0];
+    }
+  }
+  VSHtermsFree (terms);
+
+  return (TRUE);
+}