Changeset 17242

Timestamp:

Mar 31, 2008, 12:34:12 PM (18 years ago)

Author:

eugene

Message:

adding option to calculate averages

Location:

trunk/Ohana/src/relphot

Files:

• Makefile (modified) (1 diff)
• include/relphot.h (modified) (2 diffs)
• src/StarOps.c (modified) (1 diff)
• src/args.c (modified) (2 diffs)
• src/initialize.c (modified) (1 diff)
• src/relphot.c (modified) (1 diff)
• src/relphot_objects.c (added)

trunk/Ohana/src/relphot/Makefile

@@ -39 +39 @@
 $(SRC)/reload_catalogs.$(ARCH).o \
 $(SRC)/relphot.$(ARCH).o         \
+$(SRC)/relphot_objects.$(ARCH).o         \
 $(SRC)/select_images.$(ARCH).o   \
 $(SRC)/setExclusions.$(ARCH).o   \

trunk/Ohana/src/relphot/include/relphot.h

@@ -78 +78 @@
 int    USE_GRID;
 int    PLOTDELAY;
+int    UpdateAverages;
 
 int    RELPHOT_GRID_X;
@@ -191 +192 @@
 void          setMrelFinal        PROTO((Catalog *catalog));
 int           setMrelOutput       PROTO((Catalog *catalog, int Ncatalog, int mark));
+int           setMave             PROTO((Catalog *catalog, int Ncatalog));
 void          set_ZP              PROTO((double ZERO));
 int           setrefcode          PROTO((Image *image, int Nimage)); 

trunk/Ohana/src/relphot/src/StarOps.c

@@ -146 +146 @@
 }
 
+// for each average object, set the average mags based on existing equiv photometry
+int setMave (Catalog *catalog, int Ncatalog) {
+
+  int i, j, k, m, Ns, Nsecfilt, N, Nmax, Nc;
+  float Msys;
+  double *list, *dlist;
+  StatType stats;
+  PhotCode *code;
+
+  // pre-allocate a list for stats purposes
+  Nmax = 0;
+  for (i = 0; i < Ncatalog; i++) {
+    for (j = 0; j < catalog[i].Naverage; j++) {
+      Nmax = MAX (Nmax, catalog[i].average[j].Nmeasure);
+    }
+  }
+  ALLOCATE (list, double, MAX (1, Nmax));
+  ALLOCATE (dlist, double, MAX (1, Nmax));
+
+  Nsecfilt = catalog[0].Nsecfilt;
+
+  for (i = 0; i < Ncatalog; i++) {
+    for (j = 0; j < catalog[i].Naverage; j++) {
+      for (Ns = 0; Ns < Nsecfilt; Ns++) {
+
+        code = GetPhotcodebyNsec (Ns);
+        Nc = code[0].code;
+        
+        N = 0;
+        m = catalog[i].average[j].measureOffset;
+        for (k = 0; k < catalog[i].average[j].Nmeasure; k++, m++) {
+          if (catalog[i].measure[m].dbFlags & MEAS_BAD) continue;
+          if (GetPhotcodeEquivCodebyCode (catalog[i].measure[m].photcode) != Nc) continue;
+
+          Msys = PhotSys (&catalog[i].measure[m], &catalog[i].average[j], &catalog[i].secfilt[j*PhotNsec]);
+          list[N] = Msys - catalog[i].measure[m].Mcal;
+          dlist[N] = MAX (catalog[i].measure[m].dM, MIN_ERROR);
+          N++;
+        }
+        if (N < 1) continue;
+
+        liststats (list, dlist, N, &stats);
+
+        /* use sigma or error in dM for output? */
+        catalog[i].secfilt[Nsecfilt*j+Ns].M  = stats.mean;
+        catalog[i].secfilt[Nsecfilt*j+Ns].dM = MAX (stats.error, stats.sigma);
+        catalog[i].secfilt[Nsecfilt*j+Ns].Xm = (stats.Nmeas > 1) ? 100.0*log10(stats.chisq) : NAN_S_SHORT;
+        catalog[i].secfilt[Nsecfilt*j+Ns].Ncode = N;
+        catalog[i].secfilt[Nsecfilt*j+Ns].Nused = stats.Nmeas;
+      }
+    }
+  }
+
+  free (list);
+  free (dlist);
+  return (TRUE);
+}
+
 /* set measure.Mcal for all measures except ID_MEAS_NOCAL and ID_IMAGE_NOCAL */
 int setMcalOutput (Catalog *catalog, int Ncatalog) {

trunk/Ohana/src/relphot/src/args.c

@@ -181 +181 @@
   }
 
-  if ( UserPatchSelect && (argc != 2)) usage ();
-  if (!UserPatchSelect && (argc != 3)) usage ();
+  UpdateAverages = FALSE;
+  if ((N = get_argument (argc, argv, "-averages"))) {
+    remove_argument (N, &argc, argv);
+    UpdateAverages = TRUE;
+  }
+
+  if (UpdateAverages && (argc == 1)) return TRUE;
+  if (UserPatchSelect && (argc == 2)) return TRUE;
+  if (argc != 3) usage ();
 
   return TRUE;
@@ -190 +197 @@
   fprintf (stderr, "ERROR: USAGE: relphot (region) (photcode)\n");
   fprintf (stderr, "       or:    relphot (photcode) -region RA RA DEC DEC\n");
+  fprintf (stderr, "       or:    relphot -averages -region RA RA DEC DEC\n");
   fprintf (stderr, "  options: \n");
   fprintf (stderr, "  -time (start) (stop)\n");

trunk/Ohana/src/relphot/src/initialize.c

@@ -8 +8 @@
   args (argc, argv);
 
-  N = UserPatchSelect ? 1 : 2;
-  if ((photcode = GetPhotcodebyName (argv[N])) == NULL) {
-    fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", argv[N]);
-    exit (1);
-  }
-  if (photcode[0].type != PHOT_SEC) {
-    fprintf (stderr, "photcode %s is not a primary or secondary filter\n", argv[N]);
-    exit (1);
+  if (!UpdateAverages) {
+
+    N = UserPatchSelect ? 1 : 2;
+
+    if ((photcode = GetPhotcodebyName (argv[N])) == NULL) {
+      fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", argv[N]);
+      exit (1);
+    }
+    if (photcode[0].type != PHOT_SEC) {
+      fprintf (stderr, "photcode %s is not a primary or secondary filter\n", argv[N]);
+      exit (1);
+    }
+    PhotSec = GetPhotcodeNsec (photcode[0].code);
   }
 
   PhotNsec = GetPhotcodeNsecfilt ();
-  PhotSec = GetPhotcodeNsec (photcode[0].code);
 
   initstats (STATMODE);

trunk/Ohana/src/relphot/src/relphot.c

@@ -11 +11 @@
   /* get configuration info, args */
   initialize (argc, argv);
+
+  /* the object analysis is a separate process iterating over catalogs */
+  if (UpdateAverages) {
+    relphot_objects ();
+    exit (0);
+  }
 
   /* register database handle with shutdown procedure */