Changeset 33274

Timestamp:

Feb 15, 2012, 2:30:58 PM (14 years ago)

Author:

eugene

Message:

deprecate the old relphot command of the form: relphot (catalog) (photcode); fix the BrightCatalog I/O; add parallel-manual option for testing

Location:

branches/eam_branches/ipp-20111122/Ohana/src/relphot

Files:

• include/relphot.h (modified) (7 diffs)
• src/BrightCatalog.c (modified) (5 diffs)
• src/args.c (modified) (8 diffs)
• src/bcatalog.c (modified) (1 diff)
• src/help.c (modified) (2 diffs)
• src/initialize.c (modified) (3 diffs)
• src/load_catalogs.c (modified) (3 diffs)
• src/load_images.c (modified) (3 diffs)
• src/reload_catalogs.c (modified) (2 diffs)
• src/relphot.c (modified) (1 diff)

branches/eam_branches/ipp-20111122/Ohana/src/relphot/include/relphot.h

@@ -41 +41 @@
 
 typedef struct {
-    AverageTiny *average;
-    MeasureTiny *measure;
-    SecFilt     *secfilt;
-    off_t Naverage;
-    off_t Nmeasure;
+  AverageTiny *average;       // array of (minimal) average data
+  MeasureTiny *measure;       // array of (minimal) measure data 
+  SecFilt     *secfilt;       // array of secfilt data (matched to average by Nsecfilt)
+  off_t       Naverage;
+  off_t       Nmeasure;
 } BrightCatalog;
 
@@ -57 +57 @@
   unsigned char trate;
 } ImageSubset;
+
+typedef struct {
+  Catalog *catalog;           // array of catalogs generated
+  int NCATALOG;               // number of catalogs allocated
+  int Ncatalog;               // number of catalogs generated
+  int Nsecfilt;               // number of catalogs generated
+  off_t *NAVERAGE;            // allocated Averages per catalog
+  off_t *NMEASURE;            // allocated Measures per catalog
+  int   *index;               // lookup table catID -> catalog[i]
+  int   *catIDs;              // lookup table catID <- catalog[i]
+  int    maxID;               // max catID value to date
+} CatalogSplitter;
 
 /* global variables set in parameter file */
@@ -85 +97 @@
 double IMFIT_SYS_SIGMA_LIM;
 double CLOUD_TOLERANCE;
+
+int    PARALLEL;
+int    PARALLEL_MANUAL;
 
 int    VERBOSE;
@@ -144 +159 @@
 
 SkyRegion UserPatch;
-int UserPatchSelect;
+char *UserCatalog;
 
 int USE_BASIC_CHECK;
@@ -214 +229 @@
 Catalog      *load_catalogs_parallel PROTO((SkyList *sky, int *Ncatalog));
 
-SkyList      *load_images         PROTO((FITS_DB *db, char *regionName, SkyRegion *region, int RegionSelect));
+SkyList      *load_images         PROTO((FITS_DB *db, char *regionName, SkyRegion *region));
 Image        *select_images       PROTO((SkyList *skylist, Image *timage, off_t Ntimage, off_t **LineNumber, off_t *Nimage));
 
@@ -285 +300 @@
 int BrightCatalogSave(char *filename, BrightCatalog *catalog);
 BrightCatalog *BrightCatalogMerge (Catalog *catalog, int Ncatalog);
+CatalogSplitter *BrightCatalogSplitInit (int Nsecfilt);
+int BrightCatalogSplitFree (CatalogSplitter *catalogs);
+int BrightCatalogSplit (CatalogSplitter *catalogs, BrightCatalog *bcatalog);
 
 int ImageSubsetSave(char *filename, ImageSubset *image, off_t Nimage);
@@ -291 +309 @@
 int CopyAverageTiny (AverageTiny *averageT, Average *average);
 int CopyMeasureTiny (MeasureTiny *measureT, Measure *measure);
+

branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/BrightCatalog.c

@@ -241 +241 @@
     gfits_create_table_header (&theader, "BINTABLE", "MEASURE_TINY");
 
-    gfits_define_bintable_column (&theader, "E", "RA_OFF",   "ra offset",              "arcsec", 1.0, 0.0);
-    gfits_define_bintable_column (&theader, "E", "DEC_OFF",  "dec offset",             "arcsec", 1.0, 0.0);
-    gfits_define_bintable_column (&theader, "E", "MAG_SYS",  "magnitude (sys)",                NULL,     1.0, 0.0);
-    gfits_define_bintable_column (&theader, "E", "MAG_CAL",  "magnitude (cal)",                NULL,     1.0, 0.0);
-    gfits_define_bintable_column (&theader, "E", "MAG_ERR",  "magnitude (err)",                NULL,     1.0, 0.0);
-    gfits_define_bintable_column (&theader, "E", "AIRMASS",  "airmass",                  NULL,   1.0, 0.0);
-    gfits_define_bintable_column (&theader, "E", "X_CCD",    "ccd x coord",            "pix",    1.0, 0.0);
-    gfits_define_bintable_column (&theader, "E", "Y_CCD",    "ccd y coord",            "pix",    1.0, 0.0);
-    gfits_define_bintable_column (&theader, "E", "EXPTIME",  "-2.5 * log (exposure time)", "sec",        1.0, 0.0);
-    gfits_define_bintable_column (&theader, "J", "TIME",     "time of exp",              "sec",          1.0, 1.0*0x8000);
-    gfits_define_bintable_column (&theader, "J", "AVE_REF",  "pointer to average table", NULL,   1.0, 1.0*0x8000);
-    gfits_define_bintable_column (&theader, "J", "IMAGE_ID", "image",                          NULL,     1.0, 1.0*0x8000);
-    gfits_define_bintable_column (&theader, "J", "DB_FLAGS", "flags",                          NULL,     1.0, 1.0*0x8000);
-    gfits_define_bintable_column (&theader, "J", "CAT_ID",   "catalog",                        NULL,     1.0, 0.0);
-    gfits_define_bintable_column (&theader, "K", "PHOTCODE", "photcode",                       NULL,     1.0, 0.0);
+    gfits_define_bintable_column (&theader, "E", "RA_OFF",   "ra offset",                  "arcsec", 1.0, 0.0);
+    gfits_define_bintable_column (&theader, "E", "DEC_OFF",  "dec offset",                 "arcsec", 1.0, 0.0);
+    gfits_define_bintable_column (&theader, "E", "MAG_SYS",  "magnitude (sys)",             NULL,    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "E", "MAG_CAL",  "magnitude (cal)",             NULL,    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "E", "MAG_ERR",  "magnitude (err)",             NULL,    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "E", "AIRMASS",  "airmass",                     NULL,    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "E", "X_CCD",    "ccd x coord",                "pix",    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "E", "Y_CCD",    "ccd y coord",                "pix",    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "E", "EXPTIME",  "-2.5 * log (exposure time)", "sec",    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "J", "TIME",     "time of exp",                "sec",    1.0, 1.0*0x8000);
+    gfits_define_bintable_column (&theader, "J", "AVE_REF",  "pointer to average table",    NULL,    1.0, 1.0*0x8000);
+    gfits_define_bintable_column (&theader, "J", "IMAGE_ID", "image",                       NULL,    1.0, 1.0*0x8000);
+    gfits_define_bintable_column (&theader, "J", "DB_FLAGS", "flags",                       NULL,    1.0, 1.0*0x8000);
+    gfits_define_bintable_column (&theader, "J", "CAT_ID",   "catalog",                     NULL,    1.0, 1.0*0x8000);
+    gfits_define_bintable_column (&theader, "I", "PHOTCODE", "photcode",                    NULL,    1.0, 1.0*0x80);
 
     // generate the output array that carries the data
@@ -340 +340 @@
     gfits_create_table_header (&theader, "BINTABLE", "AVERAGE_TINY");
 
-    gfits_define_bintable_column (&theader, "D", "RA",                "ra (J2000)",            "degree", 1.0, 0.0);
-    gfits_define_bintable_column (&theader, "D", "DEC",               "dec (J2000)",           "degree", 1.0, 0.0);
-    gfits_define_bintable_column (&theader, "J", "NMEAS",       "number of measures",     NULL,          1.0, 0.0);
-    gfits_define_bintable_column (&theader, "J", "MEASURE_OFF", "index to measurements",  NULL,          1.0, 0.0);
-    gfits_define_bintable_column (&theader, "J", "FLAGS",       "flags",                  NULL,          1.0, 0.0);
-    gfits_define_bintable_column (&theader, "J", "CAT_ID",      "catalog ref",            NULL,          1.0, 0.0);
+    gfits_define_bintable_column (&theader, "D", "RA",          "ra (J2000)",            "degree", 1.0, 0.0);
+    gfits_define_bintable_column (&theader, "D", "DEC",         "dec (J2000)",           "degree", 1.0, 0.0);
+    gfits_define_bintable_column (&theader, "J", "NMEAS",       "number of measures",     NULL,    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "J", "MEASURE_OFF", "index to measurements",  NULL,    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "J", "FLAGS",       "flags",                  NULL,    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "J", "CAT_ID",      "catalog ref",            NULL,    1.0, 0.0);
 
     // generate the output array that carries the data
@@ -396 +396 @@
     gfits_define_bintable_column (&theader, "E", "MAG",      "ra offset",                "arcsec", 1.0, 0.0);
     gfits_define_bintable_column (&theader, "E", "MAG_ERR",  "dec offset",               "arcsec", 1.0, 0.0);
-    gfits_define_bintable_column (&theader, "E", "MAG_CHI",  "magnitude (sys)",          NULL,     1.0, 0.0);
-    gfits_define_bintable_column (&theader, "J", "FLAGS",    "magnitude (cal)",          NULL,     1.0, 0.0);
-    gfits_define_bintable_column (&theader, "I", "NCODE",    "magnitude (err)",          NULL,     1.0, 0.0);
-    gfits_define_bintable_column (&theader, "I", "NUSED",    "airmass",                  NULL,     1.0, 0.0);
+    gfits_define_bintable_column (&theader, "E", "MAG_CHI",  "magnitude (sys)",           NULL,    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "J", "FLAGS",    "magnitude (cal)",           NULL,    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "I", "NCODE",    "magnitude (err)",           NULL,    1.0, 0.0);
+    gfits_define_bintable_column (&theader, "I", "NUSED",    "airmass",                   NULL,    1.0, 0.0);
     gfits_define_bintable_column (&theader, "I", "MAG_20",   "ccd x coord",              "pix",    1.0, 0.0);
     gfits_define_bintable_column (&theader, "I", "MAG_80",   "ccd y coord",              "pix",    1.0, 0.0);
@@ -460 +460 @@
 }
 
+// merge a list of catalogs into a single BrightCatalog array set
 BrightCatalog *BrightCatalogMerge (Catalog *catalog, int Ncatalog) {
 
@@ -507 +508 @@
   return bcatalog;
 }
+
+// distribute a bright catalog across separate catalogs
+CatalogSplitter *BrightCatalogSplitInit (int Nsecfilt) {
+
+  int i;
+
+  CatalogSplitter *catalogs = NULL;
+
+  ALLOCATE (catalogs, CatalogSplitter, 1);
+
+  // as we see new BrightCatalogs, we will update maxID and extend this array
+  catalogs->maxID = 0;
+  ALLOCATE (catalogs->index, int, catalogs->maxID + 1);
+  for (i = 0; i <= catalogs->maxID; i++) catalogs->index[i] = -1;
+
+  catalogs->Nsecfilt = Nsecfilt;
+
+  catalogs->Ncatalog =  0;
+  catalogs->NCATALOG = 16;
+  catalogs->catalog = NULL;
+  ALLOCATE (catalogs->catalog, Catalog, catalogs->NCATALOG);
+
+  ALLOCATE (catalogs->catIDs,   int,   catalogs->NCATALOG);
+  ALLOCATE (catalogs->NAVERAGE, off_t, catalogs->NCATALOG);
+  ALLOCATE (catalogs->NMEASURE, off_t, catalogs->NCATALOG);
+
+  for (i = 0; i < catalogs->NCATALOG; i++) {
+    catalogs->catIDs[i] = 0;
+    catalogs->NAVERAGE[i] = 100;
+    catalogs->NMEASURE[i] = 100;
+    catalogs->catalog[i].Naverage = 0;
+    catalogs->catalog[i].Nmeasure = 0;
+    ALLOCATE (catalogs->catalog[i].averageT, AverageTiny, catalogs->NAVERAGE[i]);
+    ALLOCATE (catalogs->catalog[i].measureT, MeasureTiny, catalogs->NMEASURE[i]);
+    ALLOCATE (catalogs->catalog[i].secfilt,  SecFilt,     catalogs->NAVERAGE[i]*Nsecfilt);
+  }
+  return catalogs;
+}
+
+// distribute a bright catalog across separate catalogs
+int BrightCatalogSplitFree (CatalogSplitter *catalogs) {
+
+  // XXX don't free the catalogs : 
+  // for (i = 0; i < NCATALOG; i++) {
+  //   free (catalogs->catalog[i].averageT);
+  //   free (catalogs->catalog[i].measureT);
+  //   free (catalogs->catalog[i].secfilt);
+  // }
+  // free (catalogs->catalog);
+
+  free (catalogs->catIDs);
+  free (catalogs->NAVERAGE);
+  free (catalogs->NMEASURE);
+  free (catalogs->index);
+  free (catalogs);
+  return TRUE;
+}
+
+// distribute a bright catalog across separate catalogs
+int BrightCatalogSplit (CatalogSplitter *catalogs, BrightCatalog *bcatalog) {
+
+  int i;
+
+  int Nsecfilt = catalogs->Nsecfilt;
+
+  // find the max value of catID in this BrightCatalog
+  int catIDmax = 0;
+  for (i = 0; i < bcatalog->Naverage; i++) {
+    catIDmax = MAX(catIDmax, bcatalog->average[i].catID);
+  }
+  // XXX validate the measure ID range here
+    
+  int maxIDold = catalogs->maxID;
+  catalogs->maxID = MAX (maxIDold, catIDmax);
+    
+  // extend the index array and init
+  REALLOCATE (catalogs->index, int, catalogs->maxID + 1);
+  for (i = maxIDold + 1; i <= catalogs->maxID; i++) catalogs->index[i] = -1;
+
+  // identify the new catID values
+  for (i = 0; i < bcatalog->Naverage; i++) {
+    int catID = bcatalog->average[i].catID;
+    assert (catID > 0);
+    assert (catID < catalogs->maxID + 1);
+    int idx = catalogs->index[catID];
+    if (idx != -1) continue; // already have seen this one
+
+    // a new catID value
+    int Ncat = catalogs->Ncatalog; // the next available slot
+    catalogs->catIDs[Ncat] = catID;
+    assert (Ncat >= 0);
+    assert (Ncat < catalogs->NCATALOG);
+
+    catalogs->index[catID] = Ncat;
+    assert (catID > 0);
+    assert (catID < catalogs->maxID + 1);
+
+    catalogs->Ncatalog ++;
+
+    if (catalogs->Ncatalog >= catalogs->NCATALOG) {
+      catalogs->NCATALOG += 16;
+      fprintf (stderr, "realloc catalogs->catalog: old: %llx  ", (long long) catalogs->catalog);
+      REALLOCATE (catalogs->catalog, Catalog, catalogs->NCATALOG);
+      fprintf (stderr, "new: %llx  -  %llx\n", (long long) catalogs->catalog, (long long) (catalogs->catalog + sizeof(Catalog)*catalogs->NCATALOG));
+
+      fprintf (stderr, "realloc catalogs->NAVERAGE: old: %llx  ", (long long) catalogs->NAVERAGE);
+      REALLOCATE (catalogs->NAVERAGE, off_t,  catalogs->NCATALOG);
+      fprintf (stderr, "new: %llx  -  %llx\n", (long long) catalogs->NAVERAGE, (long long) (catalogs->NAVERAGE + sizeof(off_t)*catalogs->NCATALOG));
+
+      fprintf (stderr, "realloc catalogs->NMEASURE: old: %llx  ", (long long) catalogs->NMEASURE);
+      REALLOCATE (catalogs->NMEASURE, off_t,  catalogs->NCATALOG);
+      fprintf (stderr, "new: %llx  -  %llx\n", (long long) catalogs->NMEASURE, (long long) (catalogs->NMEASURE + sizeof(off_t)*catalogs->NCATALOG));
+
+      REALLOCATE (catalogs->catIDs,   int,   catalogs->NCATALOG);
+
+      int j;
+      for (j = catalogs->NCATALOG - 16; j < catalogs->NCATALOG; j++) {
+        // XXX call the catalog init function in libdvo?
+        catalogs->catIDs[j] = 0;
+        catalogs->NAVERAGE[j] = 100;
+        catalogs->NMEASURE[j] = 100;
+        catalogs->catalog[j].Naverage = 0;
+        catalogs->catalog[j].Nmeasure = 0;
+        ALLOCATE (catalogs->catalog[j].averageT, AverageTiny, catalogs->NAVERAGE[j]);
+        ALLOCATE (catalogs->catalog[j].measureT, MeasureTiny, catalogs->NMEASURE[j]);
+        ALLOCATE (catalogs->catalog[j].secfilt,  SecFilt,     catalogs->NAVERAGE[j]*Nsecfilt);
+      }
+    }
+  }
+
+  // each bcatalog (from each host) has a different set of catID ranges
+  // they also (probably) have contiguous ranges.  But, it is a bit tricky to be sure I
+  // can use that info, so perhaps ignore it
+
+  // assign the averages to the corresponding catalog
+  for (i = 0; i < bcatalog->Naverage; i++) {
+    int ID = bcatalog->average[i].catID;
+    int Nc = catalogs->index[ID];
+    assert (Nc > -1);
+    assert (Nc < catalogs->NCATALOG);
+
+    int Na = catalogs->catalog[Nc].Naverage;
+    catalogs->catalog[Nc].averageT[Na] = bcatalog->average[i];
+
+    // secfilt entries are grouped in blocks for each average
+    int k;
+    for (k = 0; k < Nsecfilt; k++) {
+      catalogs->catalog[Nc].secfilt[Na*Nsecfilt + k] = bcatalog->secfilt[i*Nsecfilt + k];
+    }      
+
+    catalogs->catalog[Nc].Naverage ++;
+    if (catalogs->catalog[Nc].Naverage >= catalogs->NAVERAGE[Nc]) {
+      catalogs->NAVERAGE[Nc] += 100;
+      REALLOCATE (catalogs->catalog[Nc].averageT, AverageTiny, catalogs->NAVERAGE[Nc]);
+      REALLOCATE (catalogs->catalog[Nc].secfilt,  SecFilt,     catalogs->NAVERAGE[Nc]*Nsecfilt);
+    }          
+  }
+
+  // assign the measures to the corresponding catalog
+  // XXX what about averef and related links?  Do I need them?
+  for (i = 0; i < bcatalog->Nmeasure; i++) {
+    int ID = bcatalog->measure[i].catID;
+    int Nc = catalogs->index[ID];
+    assert (Nc > -1);
+    assert (Nc < catalogs->NCATALOG);
+    int Na = catalogs->catalog[Nc].Nmeasure;
+    catalogs->catalog[Nc].measureT[Na] = bcatalog->measure[i];
+    catalogs->catalog[Nc].Nmeasure ++;
+    if (catalogs->catalog[Nc].Nmeasure >= catalogs->NMEASURE[Nc]) {
+      catalogs->NMEASURE[Nc] += 100;
+      REALLOCATE (catalogs->catalog[Nc].measureT, MeasureTiny, catalogs->NMEASURE[Nc]);
+    }          
+  }
+
+  return TRUE;
+}

branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/args.c

@@ -38 +38 @@
   UserPatch.Dmin = -90;
   UserPatch.Dmax = +90;
-  UserPatchSelect = FALSE;
   if ((N = get_argument (argc, argv, "-region"))) {
     remove_argument (N, &argc, argv);
@@ -49 +48 @@
     UserPatch.Dmax = atof (argv[N]);
     remove_argument (N, &argc, argv);
-    UserPatchSelect = TRUE;
+  }
+
+  /* specify region file by name (eg n0000/0000.00) */
+  UserCatalog = NULL;
+  if ((N = get_argument (argc, argv, "-catalog"))) {
+    remove_argument (N, &argc, argv);
+    UserCatalog = strcreate (argv[N]);
+    remove_argument (N, &argc, argv);
   }
 
@@ -65 +71 @@
   if ((N = get_argument (argc, argv, "-vv"))) {
     VERBOSE2 = VERBOSE = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+
+  // XXX for the moment, make this selection manual.  it needs to be automatic 
+  // based on the state of the SkyTable
+  PARALLEL = FALSE;
+  if ((N = get_argument (argc, argv, "-parallel"))) {
+    PARALLEL = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+  // this is a test mode : rather than launching the remote jobs and waiting for completion,
+  // relphot will simply list the remote command and wait for the user to signal completion
+  PARALLEL_MANUAL = FALSE;
+  if ((N = get_argument (argc, argv, "-parallel-manual"))) {
+    PARALLEL = TRUE; // -parallel-manual implies -parallel
+    PARALLEL_MANUAL = TRUE;
     remove_argument (N, &argc, argv);
   }
@@ -260 +282 @@
 
   if (UpdateAverages && (argc == 1)) return TRUE;
-  if (UserPatchSelect && (argc == 2)) return TRUE;
-  if (argc != 3) relphot_usage ();
+  if (argc != 2) relphot_usage ();
 
   return TRUE;
@@ -269 +290 @@
 
   int N;
+
+  // by definition, the client is not parallel 
+  PARALLEL = FALSE;
 
   if ((N = get_argument (argc, argv, "-v"))) {
@@ -325 +349 @@
   UserPatch.Dmin = -90;
   UserPatch.Dmax = +90;
-  UserPatchSelect = FALSE;
   if ((N = get_argument (argc, argv, "-region"))) {
     remove_argument (N, &argc, argv);
@@ -336 +359 @@
     UserPatch.Dmax = atof (argv[N]);
     remove_argument (N, &argc, argv);
-    UserPatchSelect = TRUE;
   }
 
@@ -343 +365 @@
   return TRUE;
 }
+

branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/bcatalog.c

@@ -57 +57 @@
 
   Ncode = Ntime = Ndophot = Nmag = Nsigma = Nimag = Nfew = Npsfqf = Ngalaxy = 0;
+
+  // copy the following fields to the subcatalog:
+  subcatalog[0].catID = catalog[0].catID;
+  subcatalog[0].filename = catalog[0].filename;
 
   /* exclude stars not in range or with too few measurements */

branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/help.c

@@ -2 +2 @@
 
 void relphot_usage (void) {
-  fprintf (stderr, "ERROR: USAGE: relphot (region) (photcode)\n");
-  fprintf (stderr, "       or:    relphot (photcode) -region RA RA DEC DEC\n");
+  fprintf (stderr, "ERROR: USAGE: relphot (photcodes) -region RA RA DEC DEC\n");
+  fprintf (stderr, "       or:    relphot (photcodes) -catalog (name)\n");
   fprintf (stderr, "       or:    relphot -averages -region RA RA DEC DEC\n");
   fprintf (stderr, "  use -h for more usage information\n");
@@ -18 +18 @@
 
 show_help:
-  fprintf (stderr, "ERROR: USAGE: relphot (region) (photcode)\n");
-  fprintf (stderr, "       or:    relphot (photcode) -region RA RA DEC DEC\n");
+  fprintf (stderr, "ERROR: USAGE: relphot (photcode) -region RA RA DEC DEC\n");
+  fprintf (stderr, "       or:    relphot (photcode) -catalog (name)\n");
   fprintf (stderr, "       or:    relphot -averages -region RA RA DEC DEC\n");
   fprintf (stderr, "  options: \n");

branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/initialize.c

@@ -1 +1 @@
 # include "relphot.h"
 
-void ParsePhotcodeList (char *word) {
-
-  Nphotcodes = 0;
-  photcodes = NULL;
-  int NPHOTCODES = 10;
-  ALLOCATE (photcodes, PhotCode *, NPHOTCODES);
-
-  /* parse the comma-separated list of photcodesKeep */
-  char *myList = strcreate(word);
-  char *list = myList;
-  char *codename = NULL;
-  char *ptr = NULL;
-  while ((codename = strtok_r (list, ",", &ptr)) != NULL) {
-    list = NULL; // pass NULL on successive strtok_r calls
-    fprintf (stderr, "PHOTCODE LIST: %s\n", myList);
-    fprintf (stderr, "codename: %s\n", codename);
-    if ((photcodes[Nphotcodes] = GetPhotcodebyName (codename)) == NULL) {
-      fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", codename);
-      exit (1);
-    }
-    if (photcodes[Nphotcodes][0].type != PHOT_SEC) {
-      fprintf (stderr, "photcode %s is not an filter type (SEC)\n", codename);
-      exit (1);
-    }
-    Nphotcodes ++;
-    CHECK_REALLOCATE (photcodes, PhotCode *, NPHOTCODES, Nphotcodes, 10);
-  }
-}
+void ParsePhotcodeList (char *word);
 
 void initialize (int argc, char **argv) {
-
-  int N;
 
   relphot_help (argc, argv);
@@ -39 +10 @@
 
   if (!UpdateAverages) {
-
-    N = UserPatchSelect ? 1 : 2;
-
-# if (0)
-    // XXX DEP
-    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 is used to select the single average photcode being processed
-    PhotSec = GetPhotcodeNsec (photcode[0].code);
-# endif
-
     // load the list of photcodes into the globals (photcodes, Nphotcodes)
-    ParsePhotcodeList (argv[N]);
+    // only a single remaining argument in this mode (the list of photcodes, eg g,r,i)
+    ParsePhotcodeList (argv[1]);
   }
-  // XXX DEP PhotNsec = GetPhotcodeNsecfilt ();
   if (USE_GRID && (Nphotcodes > 1)) {
     fprintf (stderr, "grid correction analysis currently can only operate on a single photcode\n");
@@ -111 +65 @@
   ParsePhotcodeList (argv[1]);
 }
+
+void ParsePhotcodeList (char *word) {
+
+  Nphotcodes = 0;
+  photcodes = NULL;
+  int NPHOTCODES = 10;
+  ALLOCATE (photcodes, PhotCode *, NPHOTCODES);
+
+  /* parse the comma-separated list of photcodesKeep */
+  char *myList = strcreate(word);
+  char *list = myList;
+  char *codename = NULL;
+  char *ptr = NULL;
+  while ((codename = strtok_r (list, ",", &ptr)) != NULL) {
+    list = NULL; // pass NULL on successive strtok_r calls
+    fprintf (stderr, "PHOTCODE LIST: %s\n", myList);
+    fprintf (stderr, "codename: %s\n", codename);
+    if ((photcodes[Nphotcodes] = GetPhotcodebyName (codename)) == NULL) {
+      fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", codename);
+      exit (1);
+    }
+    if (photcodes[Nphotcodes][0].type != PHOT_SEC) {
+      fprintf (stderr, "photcode %s is not an filter type (SEC)\n", codename);
+      exit (1);
+    }
+    Nphotcodes ++;
+    CHECK_REALLOCATE (photcodes, PhotCode *, NPHOTCODES, Nphotcodes, 10);
+  }
+}
+

branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/load_catalogs.c

@@ -16 +16 @@
 
   // XXX need to decide how to determine PARALLEL mode...
-# define PARALLEL 1
   if (PARALLEL & !hostID) {
     catalog = load_catalogs_parallel (skylist, Ncatalog);
@@ -114 +113 @@
     fprintf (stderr, "command: %s\n", command);
 
-    // launch the job on the remote machine (no handshake)
-    int errorInfo = 0;
-    int pid = rconnect ("ssh", table->hosts[i].hostname, command, table->hosts[i].stdio, &errorInfo, FALSE);
-    if (!pid) {
-      if (DEBUG) fprintf (stderr, "failure to start %s (error %d)\n", table->hosts[i].hostname, errorInfo);
-      exit (1);
+    if (!PARALLEL_MANUAL) {
+      // launch the job on the remote machine (no handshake)
+      int errorInfo = 0;
+      int pid = rconnect ("ssh", table->hosts[i].hostname, command, table->hosts[i].stdio, &errorInfo, FALSE);
+      if (!pid) {
+        if (DEBUG) fprintf (stderr, "failure to start %s (error %d)\n", table->hosts[i].hostname, errorInfo);
+        exit (1);
+      }
+      table->hosts[i].pid = pid; // save for future reference
+      // check that all hosts started OK?
     }
-    table->hosts[i].pid = pid; // save for future reference
-
-    // check that all hosts started OK?
   }
 
-  // watch for stdout / stderr from those jobs...
-  // wait for all of them to complete...
+  if (PARALLEL_MANUAL) {
+    fprintf (stderr, "run the relphot_client commands above.  when these are done, hit return\n");
+    getchar();
+  } else {
+    // watch for stdout / stderr from those jobs...
+    // wait for all of them to complete...
+  }
 
-  // each host generates a BrightCatalog structure, with the measure, average, etc value loaded into a single 
-  // set of arrays (of MeasureTiny, AverageTiny, Secfilt).  I need to split out the
-  // per-catalog measurements into separate catalog entries.
+  // each host generates a BrightCatalog structure, with the measure, average, etc value
+  // loaded into a single set of arrays (of MeasureTiny, AverageTiny, Secfilt).  I need to
+  // split out the per-catalog measurements into separate catalog entries.
 
-  // we don't know how many catalogs we actually need.  allocate a single entry, then
-  // below we will allocate the max number needed after the BrightCatalog is loaded
-  int NCATALOG = 1;
-  Catalog *catalog = NULL;
-  ALLOCATE (catalog, Catalog, NCATALOG);
-
-  int *NAVERAGE;
-  int *NMEASURE;
-  ALLOCATE (NAVERAGE, int, NCATALOG);
-  ALLOCATE (NMEASURE, int, NCATALOG);
-  for (i = 0; i < NCATALOG; i++) {
-    NAVERAGE[i] = 100;
-    NMEASURE[i] = 100;
-    catalog[i].Naverage = 0;
-    catalog[i].Nmeasure = 0;
-    ALLOCATE (catalog[i].averageT, AverageTiny, NAVERAGE[i]);
-    ALLOCATE (catalog[i].measureT, MeasureTiny, NMEASURE[i]);
-    ALLOCATE (catalog[i].secfilt, SecFilt, Nsecfilt*NAVERAGE[i]);
-  }
+  // set up an initial array of catalogs
+  CatalogSplitter *catalogs = BrightCatalogSplitInit (Nsecfilt);
 
   for (i = 0; i < table->Nhosts; i++) {
@@ -162 +150 @@
     assert (bcatalog);
     
-    // XXX this is wrong : I need to distinguish between the array of catalogs
-    // and the range of catID values.
+    BrightCatalogSplit (catalogs, bcatalog);
 
-    // I need to determine the set of unique catID values and count them,
-    // I'll also need a lookup table from catID -> index.
+    free (bcatalog->average);
+    free (bcatalog->measure);
+    free (bcatalog->secfilt);
+    free (bcatalog);
+  }
 
-    // XXX put the code below into a function equivalent to BrightCatalogMerge?
+  Catalog *catalog = catalogs->catalog;
+  *Ncatalog = catalogs->Ncatalog;
+  BrightCatalogSplitFree (catalogs);
 
-    // find the max value of catID and realloc to match that (check measure & average or
-    // just average?)
-    int catIDmax = 0;
-    int j;
-    for (j = 0; j < bcatalog->Naverage; j++) {
-      catIDmax = MAX(catIDmax, bcatalog->average[j].catID);
-    }
-    // XXX validate the measure ID range here
-    
-    int NCATALOG_OLD = NCATALOG;
-    NCATALOG = MAX (NCATALOG, catIDmax);
-    REALLOCATE (catalog, Catalog, NCATALOG);
-    
-    REALLOCATE (NAVERAGE, int, NCATALOG);
-    REALLOCATE (NMEASURE, int, NCATALOG);
-    for (j = NCATALOG_OLD; j < NCATALOG; j++) {
-      NAVERAGE[j] = 100;
-      NMEASURE[j] = 100;
-      catalog[j].Naverage = 0;
-      catalog[j].Nmeasure = 0;
-      ALLOCATE (catalog[j].averageT, AverageTiny, NAVERAGE[j]);
-      ALLOCATE (catalog[j].measureT, MeasureTiny, NMEASURE[j]);
-      ALLOCATE (catalog[j].secfilt, SecFilt, Nsecfilt*NAVERAGE[j]);
-    }
-
-    // each bcatalog (from each host) has a different set of catID ranges
-    // they also (probably) have contiguous ranges.  But, it is a bit tricky to be sure I
-    // can use that info, so perhaps ignore it
-
-    // assign the averages to the corresponding catalog
-    for (j = 0; j < bcatalog->Naverage; j++) {
-      int Nc = bcatalog->average[j].catID;
-      assert (Nc < NCATALOG);
-      int Na = catalog[Nc].Naverage;
-      catalog[Nc].averageT[Na] = bcatalog->average[j];
-
-      // secfilt entries are grouped in blocks for each average
-      int k;
-      for (k = 0; k < Nsecfilt; k++) {
-        catalog[Nc].secfilt[Na*Nsecfilt + k] = bcatalog->secfilt[j*Nsecfilt + k];
-      }      
-
-      catalog[Nc].Naverage ++;
-      if (catalog[Nc].Naverage >= NAVERAGE[Nc]) {
-        NAVERAGE[Nc] += 100;
-        REALLOCATE (catalog[Nc].averageT, AverageTiny, NAVERAGE[Nc]);
-        REALLOCATE (catalog[Nc].secfilt, SecFilt, Nsecfilt*NAVERAGE[Nc]);
-      }        
-    }
-
-    // assign the measures to the corresponding catalog
-    // XXX what about averef and related links?  Do I need them?
-    for (j = 0; j < bcatalog->Nmeasure; j++) {
-      int Nc = bcatalog->measure[j].catID;
-      assert (Nc < NCATALOG);
-      int Na = catalog[Nc].Nmeasure;
-      catalog[Nc].measureT[Na] = bcatalog->measure[j];
-      catalog[Nc].Nmeasure ++;
-      if (catalog[Nc].Nmeasure >= NMEASURE[Nc]) {
-        NMEASURE[Nc] += 100;
-        REALLOCATE (catalog[Nc].measureT, MeasureTiny, NMEASURE[Nc]);
-      }        
-    }
-    // XXX free the input bcatalog entries
-  }
-  *Ncatalog = NCATALOG;
   return (catalog);
 }      

branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/load_images.c

@@ -9 +9 @@
 // This function generates a subset of the images based on selections.  Input db has already
 // been loaded with the raw fits table data
-SkyList *load_images (FITS_DB *db, char *regionName, SkyRegion *region, int RegionSelect) {
+SkyList *load_images (FITS_DB *db, char *regionName, SkyRegion *region) {
 
   Image     *image, *subset;
@@ -35 +35 @@
 
   // determine the populated SkyRegions overlapping the requested area
-  if (RegionSelect) {
+
+  // if the user selects a specific catalog, do that one; 
+  // otherwise limit to the selection region
+  if (regionName) {
+    Nchar = strlen(regionName);
+    if (!strcmp (&regionName[Nchar-4], ".cpt")) regionName[Nchar-4] = 0;
+    skylist = SkyListByName (sky, regionName);
+
+    // select the images which overlap the selected sky regions
+    // 'subset' points to a new copy of the data (different from 'image')
+    subset = select_images (skylist, image, Nimage, &LineNumber, &Nsubset);
+    MARKTIME("selected %d overlapping images: %f sec\n", (int) Nsubset, dtime);
+  } else {
     if (region[0].Rmin > region[0].Rmax) {
       SkyRegion subregion;
@@ -74 +86 @@
       MARKTIME("selected %d overlapping images: %f sec\n", (int) Nsubset, dtime);
     }
-  } else {
-    Nchar = strlen(regionName);
-    if (!strcmp (&regionName[Nchar-4], ".cpt")) regionName[Nchar-4] = 0;
-    skylist = SkyListByName (sky, regionName);
-
-    // select the images which overlap the selected sky regions
-    // 'subset' points to a new copy of the data (different from 'image')
-    subset = select_images (skylist, image, Nimage, &LineNumber, &Nsubset);
-    MARKTIME("selected %d overlapping images: %f sec\n", (int) Nsubset, dtime);
-  }
+  } 
 
   // save the subset of images in the static reference in ImageOps, set up indexes

branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/reload_catalogs.c

@@ -28 +28 @@
 
   // XXX need to decide how to determine PARALLEL mode...
-# define PARALLEL 1
   if (PARALLEL & !hostID) {
     reload_catalogs_parallel (skylist);
@@ -140 +139 @@
     fprintf (stderr, "command: %s\n", command);
 
-    // launch the job on the remote machine (no handshake)
-    int errorInfo = 0;
-    int pid = rconnect ("ssh", table->hosts[i].hostname, command, table->hosts[i].stdio, &errorInfo, FALSE);
-    if (!pid) {
-      if (DEBUG) fprintf (stderr, "failure to start %s (error %d)\n", table->hosts[i].hostname, errorInfo);
-      exit (1);
+    if (!PARALLEL_MANUAL) {
+      // launch the job on the remote machine (no handshake)
+      int errorInfo = 0;
+      int pid = rconnect ("ssh", table->hosts[i].hostname, command, table->hosts[i].stdio, &errorInfo, FALSE);
+      if (!pid) {
+        if (DEBUG) fprintf (stderr, "failure to start %s (error %d)\n", table->hosts[i].hostname, errorInfo);
+        exit (1);
+      }
+      table->hosts[i].pid = pid; // save for future reference
+      // check that all hosts started OK?
     }
-    table->hosts[i].pid = pid; // save for future reference
-
-    // check that all hosts started OK?
   }
 
-  // watch for stdout / stderr from those jobs...
-  // wait for all of them to complete...
+  if (PARALLEL_MANUAL) {
+    fprintf (stderr, "run the relphot_client commands above.  when these are done, hit return\n");
+    getchar();
+  } else {
+    // watch for stdout / stderr from those jobs...
+    // wait for all of them to complete...
+  }
 
   return (TRUE);

branches/eam_branches/ipp-20111122/Ohana/src/relphot/src/relphot.c

@@ -54 +54 @@
   MARKTIME("-- load image data: %f sec\n", dtime);
 
-  /* load regions and images based on specified sky patch */
-  // XXX need to mimic old-style load by passing patch name
-  // XXX need to reduce number of global variables in use.
-  // XXX this is fairly lame: argv[1] is photcode if UserPatchSelect is true
-  // XXX clean this up: deprecate region-name here (can be optional if desired)
-  skylist = load_images (&db, argv[1], &UserPatch, UserPatchSelect);
+  /* load regions and images based on specified sky patch and/or catalog */
+  skylist = load_images (&db, UserCatalog, &UserPatch);
   MARKTIME("-- load images: %f sec\n", dtime);