Changeset 37360

Timestamp:

Sep 4, 2014, 7:57:02 AM (12 years ago)

Author:

eugene

Message:

dvoverify checks lensing indices, check image IDs in parallel context

Location:

trunk/Ohana/src/dvomerge

Files:

• include/dvoverify.h (modified) (1 diff)
• src/dvo_image_merge_dbs.c (modified) (2 diffs)
• src/dvoverify_args.c (modified) (2 diffs)
• src/dvoverify_catalogs.c (modified) (2 diffs)
• src/dvoverify_client.c (modified) (1 diff)
• src/dvoverify_utils.c (modified) (7 diffs)

trunk/Ohana/src/dvomerge/include/dvoverify.h

@@ -56 +56 @@
 int LoadImageIDs (char *catdir);
 int CheckImageID (Catalog *catalog);
+
+int SaveImageIDsSmall(char *filename);
+int LoadImageIDsSmall (char *filename);

trunk/Ohana/src/dvomerge/src/dvo_image_merge_dbs.c

@@ -179 +179 @@
   off_t lastID = IDmap->old[IDmap->Nmap-1];
 
+  // update measure.imageID
   for (i = 0; i < catalog[0].Nmeasure; i++) {
     oldID = catalog[0].measure[i].imageID;
@@ -200 +201 @@
     catalog[0].measure[i].imageID = newID;
   }
+
+  // also update lensing.imageID if lensing exists
+  for (i = 0; i < catalog[0].Nlensing; i++) {
+    oldID = catalog[0].lensing[i].imageID;
+    if (oldID == 0) continue;
+
+    newID = dvo_map_image_ID (IDmap, oldID);
+    if (newID == 0) {
+      if (oldID > lastID) {
+        fprintf (stderr, "problem with image IDs : input out of range\n");
+        fprintf (stderr, "old ID: "OFF_T_FMT", last ID: "OFF_T_FMT"\n", oldID, lastID);
+        exit (2);
+      }
+      if (!IDmap->notFound[oldID]) {
+        fprintf (stderr, "cannot find image ID "OFF_T_FMT"\n",  oldID);
+        IDmap->notFound[oldID] = TRUE;
+        // once we discover an imageID is not found, record that fact so we do not complain for every detection
+      }
+      // optionally exit here? or wait until end to report an error?
+      // exit (2);
+    }
+    catalog[0].lensing[i].imageID = newID;
+  }
   return TRUE;
 }

trunk/Ohana/src/dvomerge/src/dvoverify_args.c

@@ -56 +56 @@
   }
 
-  /*** dvoverify -parallel (dvoverify_client, actually) ignores this argument
-       dvoverify should save a subset table of just the IDs, and dvoverify_client should load it
-       dvoverify -cpt should have CHECK_IMAGE_ID = F as default ***/
+  /*** dvoverify -cpt should have CHECK_IMAGE_ID = F as default ***/
        
   CHECK_IMAGE_ID = TRUE;
@@ -150 +148 @@
   }
 
+  CHECK_IMAGE_ID = TRUE;
+  if ((N = get_argument (*argc, argv, "-skip-image-ids"))) {
+    CHECK_IMAGE_ID = FALSE;
+    remove_argument (N, argc, argv);
+  }
+
   // restrict to a portion of the sky
   UserPatch.Rmin = 0;

trunk/Ohana/src/dvomerge/src/dvoverify_catalogs.c

@@ -47 +47 @@
     }
 
+    sprintf (filename, "%s/%s.cpx", mycatdir, skylist[0].regions[i][0].name);
+    if (!VerifyTableFile (filename)) {
+      Nbad ++;
+      skipIndexCheck = TRUE;
+      AddFailures (filename);
+    }
+
+    sprintf (filename, "%s/%s.cpy", mycatdir, skylist[0].regions[i][0].name);
+    if (!VerifyTableFile (filename)) {
+      Nbad ++;
+      skipIndexCheck = TRUE;
+      AddFailures (filename);
+    }
+
     sprintf (filename, "%s/%s.cpt", mycatdir, skylist[0].regions[i][0].name);
     if (!skipIndexCheck) {
@@ -80 +94 @@
     exit (1);
   }    
+
+  // save the IDlist as a smaller FITS table
+  if (CHECK_IMAGE_ID) {
+    char filename[DVO_MAX_PATH];
+    snprintf (filename, DVO_MAX_PATH, "%s/ImageIDs.fits", CATDIR);
+    if (!SaveImageIDsSmall (filename)) {
+      fprintf (stderr, "ERROR: failure to save image IDs\n");
+      exit (2);
+    }
+  }
 
   int i;

trunk/Ohana/src/dvomerge/src/dvoverify_client.c

@@ -22 +22 @@
     // XXX in client mode, we should be reading a reduced table generated by the calling
     // serial program
-    LoadImageIDs (CATDIR);
+    char filename[DVO_MAX_PATH];
+    snprintf (filename, DVO_MAX_PATH, "%s/ImageIDs.fits", CATDIR);
+    if (!LoadImageIDsSmall (filename)) {
+      fprintf (stderr, "ERROR: failure to load image IDs\n");
+      exit (2);
+    }
   }
 

trunk/Ohana/src/dvomerge/src/dvoverify_utils.c

@@ -178 +178 @@
   catalog.catformat = DVO_FORMAT_UNDEF; // read-only,do not set the catformat
   catalog.catmode   = DVO_MODE_UNDEF; // read-only, do not set the catmode
-  catalog.catflags  = LOAD_AVES | LOAD_MEAS;
+  catalog.catflags  = LOAD_AVES | LOAD_MEAS | LOAD_LENSING | LOAD_LENSOBJ;
   catalog.Nsecfilt  = 0;
   
@@ -214 +214 @@
   int NmeasureTotal = 0;
   int measureOffsetOK = TRUE;
+  int NlensingTotal = 0;
+  int lensingOffsetOK = TRUE;
+  int NlensobjTotal = 0;
+  int lensobjOffsetOK = TRUE;
   for (i = 0; i < catalog.Naverage; i++) {
     NmeasureTotal += catalog.average[i].Nmeasure;
@@ -227 +231 @@
       fprintf (stderr, "measureOffset + Nmeasure > catalog.Nmeasure : %d %d %d\n", i, catalog.average[i].Nmeasure, (int) catalog.Nmeasure);
     }
+
+    // sum of Nlensing for each object < Nlensing for catalog
+    NlensingTotal += catalog.average[i].Nlensing;
+    if (VERBOSE && !(NlensingTotal <= catalog.Nlensing)) {
+      fprintf (stderr, "NlensingTotal > catalog.Nlensing: %d %d %d\n", i, catalog.average[i].Nlensing, (int) catalog.Nlensing);
+    }
+    // lensingOffset needs to be in range for each object
+    lensingOffsetOK &= (catalog.average[i].lensingOffset < catalog.Nlensing);
+    if (VERBOSE && !(catalog.average[i].lensingOffset < catalog.Nlensing)) {
+      fprintf (stderr, "lensingOffset >= catalog.Nlensing: %d %d %d\n", i, catalog.average[i].lensingOffset, (int) catalog.Nlensing);
+    }
+    lensingOffsetOK &= (catalog.average[i].lensingOffset + catalog.average[i].Nlensing <= catalog.Nlensing);
+    if (VERBOSE && !(catalog.average[i].lensingOffset + catalog.average[i].Nlensing <= catalog.Nlensing)) {
+      fprintf (stderr, "lensingOffset + Nlensing > catalog.Nlensing : %d %d %d\n", i, catalog.average[i].Nlensing, (int) catalog.Nlensing);
+    }
+
+    NlensobjTotal += catalog.average[i].Nlensobj;
+    if (VERBOSE && !(NlensobjTotal <= catalog.Nlensobj)) {
+      fprintf (stderr, "NlensobjTotal > catalog.Nlensobj: %d %d %d\n", i, catalog.average[i].Nlensobj, (int) catalog.Nlensobj);
+    }
+    lensobjOffsetOK &= (catalog.average[i].lensobjOffset < catalog.Nlensobj);
+    if (VERBOSE && !(catalog.average[i].lensobjOffset < catalog.Nlensobj)) {
+      fprintf (stderr, "lensobjOffset >= catalog.Nlensobj: %d %d %d\n", i, catalog.average[i].lensobjOffset, (int) catalog.Nlensobj);
+    }
+    lensobjOffsetOK &= (catalog.average[i].lensobjOffset + catalog.average[i].Nlensobj <= catalog.Nlensobj);
+    if (VERBOSE && !(catalog.average[i].lensobjOffset + catalog.average[i].Nlensobj <= catalog.Nlensobj)) {
+      fprintf (stderr, "lensobjOffset + Nlensobj > catalog.Nlensobj : %d %d %d\n", i, catalog.average[i].Nlensobj, (int) catalog.Nlensobj);
+    }
   }
 
@@ -236 +268 @@
   if (NmeasureTotal != catalog.Nmeasure) {
     fprintf (stderr, "ERROR: catalog %s has an invalid Nmeasure\n", catalog.filename);
+    status = FALSE;
+  }
+
+  if (!lensingOffsetOK) {
+    fprintf (stderr, "ERROR: catalog %s has an invalid lensingOffset\n", catalog.filename);
+    status = FALSE;
+  }
+
+  if (NlensingTotal != catalog.Nlensing) {
+    fprintf (stderr, "ERROR: catalog %s has an invalid Nlensing\n", catalog.filename);
+    status = FALSE;
+  }
+
+  if (!lensobjOffsetOK) {
+    fprintf (stderr, "ERROR: catalog %s has an invalid lensobjOffset\n", catalog.filename);
+    status = FALSE;
+  }
+
+  if (NlensobjTotal != catalog.Nlensobj) {
+    fprintf (stderr, "ERROR: catalog %s has an invalid Nlensobj\n", catalog.filename);
     status = FALSE;
   }
@@ -247 +299 @@
   }
 
-  int objIDsOK = TRUE;
-  int catIDsOK = TRUE;
-  int averefOK = TRUE;
-
-  for (i = 0; i < catalog.Naverage; i++) {
-    m = catalog.average[i].measureOffset;
-    for (j = 0; j < catalog.average[i].Nmeasure; j++) {
-      objIDsOK &= (catalog.average[i].objID == catalog.measure[m+j].objID);
-      catIDsOK &= (catalog.average[i].catID == catalog.measure[m+j].catID);
-      averefOK &= (catalog.measure[m+j].averef == i);
-    }
-  }
+  // check measure <-> average links
+  {
+    int objIDsOK = TRUE;
+    int catIDsOK = TRUE;
+    int averefOK = TRUE;
+
+    for (i = 0; i < catalog.Naverage; i++) {
+      m = catalog.average[i].measureOffset;
+      for (j = 0; j < catalog.average[i].Nmeasure; j++) {
+        objIDsOK &= (catalog.average[i].objID == catalog.measure[m+j].objID);
+        catIDsOK &= (catalog.average[i].catID == catalog.measure[m+j].catID);
+        averefOK &= (catalog.measure[m+j].averef == i);
+      }
+    }
     
-  if (!objIDsOK) {
-    fprintf (stderr, "ERROR: catalog %s has invalid obj IDs\n", catalog.filename);
-    status = FALSE;
-  }
-  if (!catIDsOK) {
-    fprintf (stderr, "ERROR: catalog %s has invalid cat IDs\n", catalog.filename);
-    status = FALSE;
-  }
-  if (!averefOK) {
-    fprintf (stderr, "ERROR: catalog %s has invalid averef values\n", catalog.filename);
-    status = FALSE;
+    if (!objIDsOK) {
+      fprintf (stderr, "ERROR: catalog %s has invalid obj IDs\n", catalog.filename);
+      status = FALSE;
+    }
+    if (!catIDsOK) {
+      fprintf (stderr, "ERROR: catalog %s has invalid cat IDs\n", catalog.filename);
+      status = FALSE;
+    }
+    if (!averefOK) {
+      fprintf (stderr, "ERROR: catalog %s has invalid averef values\n", catalog.filename);
+      status = FALSE;
+    }
+  }
+
+  // check lensing <-> average links
+  {
+    int objIDsOK = TRUE;
+    int catIDsOK = TRUE;
+    int averefOK = TRUE;
+
+    for (i = 0; i < catalog.Naverage; i++) {
+      m = catalog.average[i].lensingOffset;
+      for (j = 0; j < catalog.average[i].Nlensing; j++) {
+        objIDsOK &= (catalog.average[i].objID == catalog.lensing[m+j].objID);
+        catIDsOK &= (catalog.average[i].catID == catalog.lensing[m+j].catID);
+        averefOK &= (catalog.lensing[m+j].averef == i);
+      }
+    }
+    
+    if (!objIDsOK) {
+      fprintf (stderr, "ERROR: catalog %s has invalid lensing obj IDs\n", catalog.filename);
+      status = FALSE;
+    }
+    if (!catIDsOK) {
+      fprintf (stderr, "ERROR: catalog %s has invalid lensing cat IDs\n", catalog.filename);
+      status = FALSE;
+    }
+    if (!averefOK) {
+      fprintf (stderr, "ERROR: catalog %s has invalid lensing averef values\n", catalog.filename);
+      status = FALSE;
+    }
+  }
+
+  // check lensobj <-> average links
+  {
+    int objIDsOK = TRUE;
+    int catIDsOK = TRUE;
+
+    for (i = 0; i < catalog.Naverage; i++) {
+      m = catalog.average[i].lensobjOffset;
+      for (j = 0; j < catalog.average[i].Nlensobj; j++) {
+        objIDsOK &= (catalog.average[i].objID == catalog.lensobj[m+j].objID);
+        catIDsOK &= (catalog.average[i].catID == catalog.lensobj[m+j].catID);
+      }
+    }
+    
+    if (!objIDsOK) {
+      fprintf (stderr, "ERROR: catalog %s has invalid lensobj obj IDs\n", catalog.filename);
+      status = FALSE;
+    }
+    if (!catIDsOK) {
+      fprintf (stderr, "ERROR: catalog %s has invalid lensobj cat IDs\n", catalog.filename);
+      status = FALSE;
+    }
   }
 
@@ -321 +428 @@
     }
   }
+
+  for (i = 0; i < catalog[0].Naverage; i++) {
+    m = catalog[0].average[i].lensingOffset;
+    for (j = 0; j < catalog[0].average[i].Nlensing; j++) {
+      id = catalog[0].lensing[m+j].imageID;
+      if (id == 0) continue; // detections from ref photcodes can (should) have unset image IDs
+      if (id > maxID) {
+        Nfail ++;
+        continue;
+        // is this sufficient to catch IDs set without an image table?
+      }
+      if (IDlist) {
+        if (IDlist[id] < 0) {
+          Nfail ++;
+          continue;
+        }
+      } else {
+        if (id > 0) {
+          Nfail ++;
+          continue;
+        }
+      }
+    }
+  }
+
   return Nfail;
 }
@@ -382 +514 @@
   return TRUE;
 }
+
+# define GET_COLUMN(OUT,NAME,TYPE) \
+  OUT = gfits_get_bintable_column_data (&theader, &ftable, NAME, type, &Nrow, &Ncol); \
+  myAssert (!strcmp(type, #TYPE), "wrong column type");
+
+// STATUS is value expected for success
+# define CHECK_STATUS(STATUS,MSG,...)                                   \
+  if (!(STATUS)) {                                                      \
+    fprintf (stderr, MSG, __VA_ARGS__);                                 \
+    return FALSE;                                                       \
+  }
+
+// write out the IDmap data for clients to read
+int SaveImageIDsSmall(char *filename) {
+
+  Header header;
+  Header theader;
+  Matrix matrix;
+  FTable ftable;
+
+  gfits_init_header (&header);
+  header.extend = TRUE;
+  gfits_create_header (&header);
+  gfits_create_matrix (&header, &matrix);
+
+  gfits_create_table_header (&theader, "BINTABLE", "IMAGE_IDS");
+
+  gfits_define_bintable_column (&theader, "J", "IMAGE_IDS", "image IDs", NULL, 1.0, FT_BZERO_INT32);
+
+  // generate the output array that carries the data
+  gfits_create_table (&theader, &ftable);
+
+  // add the columns to the output array
+  gfits_set_bintable_column (&theader, &ftable, "IMAGE_IDS", IDlist, maxID + 1);
+
+  FILE *f = fopen (filename, "w");
+  if (!f) {
+    fprintf (stderr, "ERROR: cannot open image ID file for output %s\n", filename);
+    return FALSE;
+  }
+
+  int status;
+  status = gfits_fwrite_header  (f, &header);
+  CHECK_STATUS (status, "ERROR: cannot write header for image ID file %s\n", filename);
+
+  status = gfits_fwrite_matrix  (f, &matrix);
+  CHECK_STATUS (status, "ERROR: cannot write matrix for image ID file %s\n", filename);
+
+  status = gfits_fwrite_Theader (f, &theader);
+  CHECK_STATUS (status, "ERROR: cannot write table header for image ID file %s\n", filename);
+
+  status = gfits_fwrite_table  (f, &ftable);
+  CHECK_STATUS (status, "ERROR: cannot write table data for image ID file %s\n", filename);
+
+  gfits_free_header (&header);
+  gfits_free_matrix (&matrix);
+  gfits_free_header (&theader);
+  gfits_free_table (&ftable);
+
+  int fd = fileno (f);
+
+  status = fflush (f);
+  CHECK_STATUS (!status, "ERROR: cannot flush file image ID file %s\n", filename);
+
+  status = fsync (fd);
+  CHECK_STATUS (!status, "ERROR: cannot flush file image ID file %s\n", filename);
+
+  status = fclose (f);
+  CHECK_STATUS (!status, "ERROR: problem closing image ID file file %s\n", filename);
+
+  return TRUE;
+}
+
+int LoadImageIDsSmall (char *filename) {
+
+  int Ncol;
+  off_t Nrow;
+  Header header;
+  Header theader;
+  Matrix matrix;
+  FTable ftable;
+
+  FILE *f = fopen (filename, "r");
+  if (!f) {
+    fprintf (stderr, "ERROR: cannot open image subset file %s\n", filename);
+    return FALSE;
+  }
+
+  /* load in PHU segment (ignore) */
+  if (!gfits_fread_header (f, &header)) {
+    if (VERBOSE) fprintf (stderr, "can't read image subset header\n");
+    fclose (f);
+    return FALSE;
+  }
+  if (!gfits_fread_matrix (f, &matrix, &header)) {
+    if (VERBOSE) fprintf (stderr, "can't read image subset matrix\n");
+    gfits_free_header (&header);
+    fclose (f);
+    return FALSE;
+  }
+
+  ftable.header = &theader;
+
+  // load data for this header 
+  if (!gfits_load_header (f, &theader)) {
+    fclose (f);
+    return FALSE;
+  }
+
+  // read the fits table bytes
+  if (!gfits_fread_ftable_data (f, &ftable, FALSE)) {
+    fclose (f);
+    return FALSE;
+  }
+  fclose (f);
+
+  char type[16];
+
+  GET_COLUMN (IDlist, "IMAGE_IDS", int);
+  maxID = Nrow - 1;
+  fprintf (stderr, "loaded data for %lld images\n", (long long) Nrow);
+
+  return TRUE;
+}
+