Changeset 40401 for trunk/Ohana/src/opihi/dvo/avperiodogram.c

Timestamp:

Apr 19, 2018, 9:41:13 PM (8 years ago)

Author:

eugene

Message:

updates to avperiodogram

File:

• trunk/Ohana/src/opihi/dvo/avperiodogram.c (modified) (12 diffs)

trunk/Ohana/src/opihi/dvo/avperiodogram.c

@@ -20 +20 @@
 
 void PeriodogramResultFree (PeriodogramResult *result);
+void PeriodogramResultSave (PeriodogramResult *result, char *extname, FILE *foutput, Average *average);
 PeriodogramResult *periodogram_fm_raw (opihi_flt *time, opihi_flt *flux, opihi_flt *dflux, int Npts);
 
@@ -27 +28 @@
 int avperiodogram (int argc, char **argv) {
   
-# if (0)
-
   int N, next, Nfields;
 
@@ -73 +72 @@
   // this does all the work of re-packaging the command, calling it on the remote machines, then loading in the results
   if (PARALLEL && !HOST_ID) {
-    // XXX this is not going to work : the measure args will have been stripped off
-    int status = HostTableParallelOps (skylist, argc, argv, RESULT_FILE, TRUE, 0, VERBOSE);
-
-    dbFreeFields (fields, Nfields);
-    dbFreeStack (stack, Nstack);
-    free (stack);
+    int status = HostTableParallelOps (skylist, argc, argv, NULL, TRUE, 0, VERBOSE);
+    SkyListFree (skylist);
     FreeSkyRegionSelection (selection);
-    dvo_catalog_free (&catalog);
-
     return status;
   }
@@ -157 +150 @@
   }
 
+  // command-line is of the form: avperiodogram (output) where (field op value)...
+  if (argc < 2) goto help;
+  char *output = strcreate (argv[1]);
+
+  // **** generate output file ****
+  // we save the results in a single FITS file, one extension per object
+  FILE *foutput = NULL;
+  if (RESULT_FILE) {
+    foutput = fopen (RESULT_FILE, "w");
+  } else {
+    foutput = fopen (output, "w");
+  }
+  // generate the PHU and write to disk
+  // XXX add some metadata here?
+  Header header;
+  Matrix matrix;
+  gfits_init_header (&header);
+  gfits_init_matrix (&matrix);
+  header.extend = TRUE;
+  gfits_create_header (&header);
+  gfits_create_matrix (&header, &matrix);
+  gfits_fwrite_header (foutput, &header);
+  gfits_fwrite_matrix (foutput, &matrix);
+  gfits_free_header (&header);
+  gfits_free_matrix (&matrix);
+  
   Catalog catalog;
   dvo_catalog_init (&catalog, TRUE);
@@ -166 +185 @@
   // init locally static variables (time refs)
   dbExtractAveragesInit (); 
-
-  // command-line is of the form: avperiodogram (output) where (field op value)...
-  if (argc < 2) goto help;
-  char *output = strcreate (argv[1]);
 
   // examine line for 'where' or 'match to'.  'match to' is forbidden
@@ -197 +212 @@
   ALLOCATE (values, dbValue, Nfields);
 
+  int Nobject = 0;
+  char extname[80];
+
   // grab data from all selected sky regions
   struct sigaction *old_sigaction = SetInterrupt();
@@ -221 +239 @@
     }
     dvo_catalog_unlock (&catalog);
-
-    // we save the results in a single FITS file, one extension per object
-    FILE *foutput = NULL;
-    if (RESULT_FILE) {
-      foutput = fopen (RESULT_FILE, "w");
-    } else {
-      foutput = fopen (output, "w");
-    }
-
-    // generate the PHU and write to disk
-    Header header;
-    Matrix matrix;
-    gfits_init_header (&header);
-    gfits_init_matrix (&matrix);
-    header.extend = TRUE;
-    gfits_create_header (&header);
-    gfits_create_matrix (&header, &matrix);
-    gfits_fwrite_header (foutput, &header);
-    gfits_fwrite_matrix (foutput, &matrix);
-    gfits_free_header (&header);
-    gfits_free_matrix (&matrix);
 
     // Scan the catalog for objects which match the WHERE clause
@@ -299 +296 @@
       PeriodogramResult *result = periodogram_fm_raw (time, mag, dmag, Npts);
 
-      PeriodogramResultSave (result, extname, foutput);
+      // XXX if we have 1e6 objects, we will have output files with sizes of ~80GB
+      // XXX warn or exit if Nboject > 1e6?
+      snprintf (extname, 80, "OBJ_%06d", Nobject);
+      PeriodogramResultSave (result, extname, foutput, average);
       PeriodogramResultFree (result);
+      Nobject ++;
     }
 
@@ -330 +331 @@
 
  help:
-
-# endif
   gprint (GP_ERR, "USAGE: avperiodogram (output) where (expression) -where-measure (expression)\n");
   gprint (GP_ERR, "  the where (expression) is a boolean to limit the objects analysed based on average properties\n");
   gprint (GP_ERR, "  NOTE: the optional -where-measure (expression) MUST come after the average where expression\n");
   return (FALSE);
-
-
 }
 
@@ -470 +467 @@
 }
 
-/*
-
 // write the period and power vectors to an extension in this output file
-void PeriodogramResultFree (PeriodogramResult *result, char *extname, FILE *foutput) {
+void PeriodogramResultSave (PeriodogramResult *result, char *extname, FILE *foutput, Average *average) {
 
   // generate a binary table 
@@ -482 +477 @@
   outtable.header = &outheader;
 
-  // XXX add some metadata to the header to identify the object
-  // extname should be related to the object IDs
   gfits_create_table_header (&outheader, "BINTABLE", extname);
+
+  // add some metadata derived from average:
+  gfits_modify (&outheader, "RA_OBJ",  "%lf", 1, average->R);
+  gfits_modify (&outheader, "DEC_OBJ", "%lf", 1, average->D);
+  gfits_modify (&outheader, "OBJ_ID",  "%d", 1, average->objID);
+  gfits_modify (&outheader, "CAT_ID",  "%d", 1, average->catID);
 
   gfits_define_bintable_column (&outheader, "D", "PERIOD", NULL, NULL, 1.0, 0.0);
   gfits_define_bintable_column (&outheader, "D", "POWER", NULL, NULL, 1.0, 0.0);
 
-  gfits_create_table (theader, ftable);
+  gfits_create_table (&outheader, &outtable);
 
   gfits_set_bintable_column (&outheader, &outtable, "PERIOD", result->period, result->Nperiods);
@@ -501 +500 @@
   gfits_free_table  (&outtable);
 }
-
-*/