Index: /branches/eam_branches/ipp-20111122/Ohana/src/addstar/src/mkcmf.c
===================================================================
--- /branches/eam_branches/ipp-20111122/Ohana/src/addstar/src/mkcmf.c	(revision 33134)
+++ /branches/eam_branches/ipp-20111122/Ohana/src/addstar/src/mkcmf.c	(revision 33135)
@@ -19,4 +19,6 @@
 void writeStars_PS1_DEV_0 (FTable *ftable, double *X, double *Y, double *M, int Nstars);
 
+int ADDNOISE = TRUE;
+
 int main (int argc, char **argv) {
 
@@ -34,4 +36,18 @@
   Coords coords;
 
+  static char *photcode = "SIMTEST.r.Chip";
+  if ((N = get_argument (argc, argv, "-photcode"))) {
+    remove_argument (N, &argc, argv);
+    photcode = strcreate (argv[N]);
+    remove_argument (N, &argc, argv);
+  }
+
+  static double mjd = NAN;
+  if ((N = get_argument (argc, argv, "-mjd"))) {
+    remove_argument (N, &argc, argv);
+    mjd = strtod (argv[N], NULL);
+    remove_argument (N, &argc, argv);
+  }
+
   static char *date = "2001-01-01";
   if ((N = get_argument (argc, argv, "-date"))) {
@@ -66,58 +82,22 @@
   }
 
+  // add support for all cmf types
+  int FROM_COORDS = FALSE;
+  if ((N = get_argument (argc, argv, "-coords"))) {
+    remove_argument (N, &argc, argv);
+    FROM_COORDS = TRUE;
+  }
+
+  // add support for all cmf types
+  ADDNOISE = TRUE;
+  if ((N = get_argument (argc, argv, "-no-noise"))) {
+    remove_argument (N, &argc, argv);
+    ADDNOISE = FALSE;
+  }
+
   if (argc != 3) {
     fprintf (stderr, "USAGE mkcmf (input) (output) [-date date] [-time time] [-radec ra dec] [-cmftype type]\n");
     exit (2);
   }
-
-  f = fopen (argv[1], "r");
-  if (f == NULL) {
-    fprintf (stderr, "unable to open input file %s\n", argv[1]);
-    exit (1);
-  }
-    
-  // load test stars from a file:
-  Nstars = 0;
-  NSTARS = 100;
-  ALLOCATE (X, double, NSTARS);
-  ALLOCATE (Y, double, NSTARS);
-  ALLOCATE (M, double, NSTARS);
-
-  Xmax = Ymax = 0;
-  while (fscanf (f, "%lf %lf %lf", &X[Nstars], &Y[Nstars], &M[Nstars]) != EOF) {
-    Xmax  = MAX(Xmax, X[Nstars]);
-    Ymax  = MAX(Ymax, Y[Nstars]);
-
-    if (Nstars == NSTARS - 1) {
-      NSTARS += 100;
-      REALLOCATE (X, double, NSTARS);
-      REALLOCATE (Y, double, NSTARS);
-      REALLOCATE (M, double, NSTARS);
-    }
-    Nstars ++;
-  }
-
-  // create primary header
-  gfits_init_header (&header);    
-  header.extend = TRUE;
-  gfits_create_header (&header);
-  gfits_create_matrix (&header, &matrix);
-  gfits_print (&header, "NEXTEND", "%d", 1, 1);
-
-  // XXX add minimum needed header fields 
-  gfits_print (&header, "IMNAXIS1", "%d", 1, (int)(Xmax + 50));
-  gfits_print (&header, "IMNAXIS2", "%d", 1, (int)(Ymax + 50));
-
-  gfits_modify (&header, "NSTARS",   "%d", 1, Nstars);
-  gfits_modify (&header, "PHOTCODE", "%s", 1, "SIMTEST.r.Chip");
-  gfits_modify (&header, "DATE-OBS", "%s", 1, date);
-  gfits_modify (&header, "UTC-OBS",  "%s", 1, time);
-  gfits_modify (&header, "ZERO_PT", "%lf", 1, 25.0);
-  gfits_modify (&header, "EXPTIME", "%lf", 1, 1.0);
-  gfits_modify (&header, "AIRMASS", "%lf", 1, 1.0);
-  gfits_modify (&header, "NASTRO",   "%d", 1, 10); 
-
-  gfits_modify (&header, "IMAGEID",  "%d", 1, (int)(1000*drand48())); 
-  gfits_modify (&header, "SOURCEID", "%d", 1, (int)(100*drand48())); 
 
   /* bore site center guess */
@@ -137,4 +117,72 @@
   coords.Npolyterms = 1;
 
+  // load stars and generate complete output fields
+  f = fopen (argv[1], "r");
+  if (f == NULL) {
+    fprintf (stderr, "unable to open input file %s\n", argv[1]);
+    exit (1);
+  }
+    
+  // load test stars from a file:
+  Nstars = 0;
+  NSTARS = 100;
+  ALLOCATE (X, double, NSTARS);
+  ALLOCATE (Y, double, NSTARS);
+  ALLOCATE (M, double, NSTARS);
+
+  Xmax = Ymax = 0;
+  while (1) {
+    int status;
+    if (FROM_COORDS) {
+      double ra, dec, mag, xraw, yraw, mraw;
+      status = fscanf (f, "%lf %lf %lf %lf %lf %lf", &ra, &dec, &mag, &xraw, &yraw, &mraw);
+      if (status != EOF) {
+	RD_to_XY (&X[Nstars], &Y[Nstars], ra, dec, &coords);
+      }
+      M[Nstars] = mraw;
+    } else {
+      status = fscanf (f, "%lf %lf %lf", &X[Nstars], &Y[Nstars], &M[Nstars]);
+    }
+    if (status == EOF) break;
+    
+    Xmax  = MAX(Xmax, X[Nstars]);
+    Ymax  = MAX(Ymax, Y[Nstars]);
+
+    if (Nstars == NSTARS - 1) {
+      NSTARS += 100;
+      REALLOCATE (X, double, NSTARS);
+      REALLOCATE (Y, double, NSTARS);
+      REALLOCATE (M, double, NSTARS);
+    }
+    Nstars ++;
+  }
+
+  // create primary header
+  gfits_init_header (&header);    
+  header.extend = TRUE;
+  gfits_create_header (&header);
+  gfits_create_matrix (&header, &matrix);
+  gfits_print (&header, "NEXTEND", "%d", 1, 1);
+
+  // XXX add minimum needed header fields 
+  gfits_print (&header, "IMNAXIS1", "%d", 1, (int)(Xmax + 50));
+  gfits_print (&header, "IMNAXIS2", "%d", 1, (int)(Ymax + 50));
+
+  gfits_modify (&header, "NSTARS",   "%d", 1, Nstars);
+  gfits_modify (&header, "PHOTCODE", "%s", 1, photcode);
+  if (isfinite(mjd)) {
+    gfits_modify (&header, "MJD-OBS", "%lf", 1, mjd);
+  } else {
+    gfits_modify (&header, "DATE-OBS", "%s", 1, date);
+    gfits_modify (&header, "UTC-OBS",  "%s", 1, time);
+  }
+  gfits_modify (&header, "ZERO_PT", "%lf", 1, 25.0);
+  gfits_modify (&header, "EXPTIME", "%lf", 1, 1.0);
+  gfits_modify (&header, "AIRMASS", "%lf", 1, 1.0);
+  gfits_modify (&header, "NASTRO",   "%d", 1, 10); 
+
+  gfits_modify (&header, "IMAGEID",  "%d", 1, (int)(1000*drand48())); 
+  gfits_modify (&header, "SOURCEID", "%d", 1, (int)(100*drand48())); 
+
   PutCoords (&coords, &header);
   gfits_modify (&header, "EXTNAME",   "%s", 1, "Chip.hdr");
@@ -273,7 +321,13 @@
     fSN = 1.0 / sqrt(flux);
 
-    stars[i].X = X[i] + FX * fSN * rnd_gauss(0.0, 1.0);
-    stars[i].Y = Y[i] + FY * fSN * rnd_gauss(0.0, 1.0);
-    stars[i].M = M[i] + fSN*rnd_gauss(0.0, 1.0);
+    stars[i].X = X[i];
+    stars[i].Y = Y[i];
+    stars[i].M = M[i];
+
+    if (ADDNOISE) {
+      stars[i].X += FX * fSN * rnd_gauss(0.0, 1.0);
+      stars[i].Y += FY * fSN * rnd_gauss(0.0, 1.0);
+      stars[i].M += fSN*rnd_gauss(0.0, 1.0);
+    }
 
     stars[i].dX = FX * fSN;
@@ -312,7 +366,13 @@
     fSN = 1.0 / sqrt(flux);
 
-    stars[i].X = X[i] + FX * fSN * rnd_gauss(0.0, 1.0);
-    stars[i].Y = Y[i] + FY * fSN * rnd_gauss(0.0, 1.0);
-    stars[i].M = M[i] + fSN*rnd_gauss(0.0, 1.0);
+    stars[i].X = X[i];
+    stars[i].Y = Y[i];
+    stars[i].M = M[i];
+
+    if (ADDNOISE) {
+      stars[i].X += FX * fSN * rnd_gauss(0.0, 1.0);
+      stars[i].Y += FY * fSN * rnd_gauss(0.0, 1.0);
+      stars[i].M += fSN*rnd_gauss(0.0, 1.0);
+    }
 
     stars[i].dX = FX * fSN;
@@ -354,7 +414,13 @@
     fSN = 1.0 / sqrt(flux);
 
-    stars[i].X = X[i] + FX * fSN * rnd_gauss(0.0, 1.0);
-    stars[i].Y = Y[i] + FY * fSN * rnd_gauss(0.0, 1.0);
-    stars[i].M = M[i] + fSN*rnd_gauss(0.0, 1.0);
+    stars[i].X = X[i];
+    stars[i].Y = Y[i];
+    stars[i].M = M[i];
+
+    if (ADDNOISE) {
+      stars[i].X += FX * fSN * rnd_gauss(0.0, 1.0);
+      stars[i].Y += FY * fSN * rnd_gauss(0.0, 1.0);
+      stars[i].M += fSN*rnd_gauss(0.0, 1.0);
+    }
 
     stars[i].dX = FX * fSN;
@@ -396,7 +462,13 @@
     fSN = 1.0 / sqrt(flux);
 
-    stars[i].X = X[i] + FX * fSN * rnd_gauss(0.0, 1.0);
-    stars[i].Y = Y[i] + FY * fSN * rnd_gauss(0.0, 1.0);
-    stars[i].M = M[i] + fSN*rnd_gauss(0.0, 1.0);
+    stars[i].X = X[i];
+    stars[i].Y = Y[i];
+    stars[i].M = M[i];
+
+    if (ADDNOISE) {
+      stars[i].X += FX * fSN * rnd_gauss(0.0, 1.0);
+      stars[i].Y += FY * fSN * rnd_gauss(0.0, 1.0);
+      stars[i].M += fSN*rnd_gauss(0.0, 1.0);
+    }
 
     stars[i].dX = FX * fSN;