Changeset 15210

Timestamp:

Oct 3, 2007, 4:56:44 PM (19 years ago)

Author:

eugene

Message:

cleanup code, allow multi-file, multi-extension addstar

Location:

trunk/Ohana/src/addstar

Files:

• include/addstar.h (modified) (2 diffs)
• src/LoadStars.c (modified) (2 diffs)

trunk/Ohana/src/addstar/include/addstar.h

@@ -39 +39 @@
   char *refcat;
 } DVO_DATA;
+
+typedef struct {
+  char *exthead;
+  char *extdata;
+  char *exttype;
+  int extnum_head;
+  int extnum_data;
+} HeaderSet;
 
 typedef struct sockaddr_in SockAddress;
@@ -154 +162 @@
 Stars     *grefcat                PROTO((char *Refcat, SkyRegion *catstats, int photcode, int *nstars));
 Stars     *grefstars              PROTO((char *file, int photcode, int *Nstars));
+
 Stars     *LoadStars              PROTO((char *file, int *Nstars, Image **images, int *Nimages, int photcode));
+Header   **LoadHeaders            PROTO((FILE *f, int *mode, int *Nheader));
+HeaderSet *MatchHeaders           PROTO((int **extsize, int *nimage, int mode, Header **headers, int Nheaders));
+int        LoadData               PROTO((FILE *f, char *file, Image **images, int *nvalid, Stars **stars, int *Nstars, Header **headers, int *extsize, HeaderSet *headerSets, int NheaderSets));
+
 int        in_image               PROTO((double r, double d, Image *image));
 int        load_pt_catalog        PROTO((Catalog *catalog, SkyRegion *region));  /*** choose new name ***/

trunk/Ohana/src/addstar/src/LoadStars.c

@@ -3 +3 @@
 Stars *LoadStars (char *filename, int *Nstars, Image **images, int *Nimages, int photcode) {
 
-  int i, j, N, Nfile, Nheader, NHEADER, Nimage, NIMAGE;
-  int Nskip, Nhead, Ndata, done, status, mode, NinStars;
-  char **file, *name;
+  int i, Nfile, Nheaders, NheaderSets, mode, *extsize;
+  char **file;
   FILE *f;
   glob_t globList;
-  char **exthead, **extdata, **exttype, tmpword[80];
-  int *extnum_head, *extnum_data, *extsize;
-  Header *header, **headers;
-  Image *image;
-  Stars *inStars, *stars;
+  Header **headers;
+  Stars *stars;
+  HeaderSet *headerSets;
 
   // parse the filename as a glob
@@ -32 +29 @@
   }
 
-  // open the first file, read the PHU header
-  f = fopen (file[0], "r");
-  if (f == NULL) {
-    fprintf (stderr, "ERROR: can't read header for %s\n", file[0]);
-    exit (1);
-  }
-  ALLOCATE (header, Header, 1);
-  gfits_fread_header (f, header);
-
-  mode = GetFileMode (header);
-
+  *Nimages = 0;
+  *Nstars = 0;
+  *images = NULL;
   stars = NULL;
-  *Nstars = 0;
-
-  /*** load data from a single PHU or a collection of PHU files ***/
-  if ((mode != SIMPLE_MEF) && (mode != MOSAIC_MEF)) {
-    Nimage = Nfile;
-    ALLOCATE (image, Image, Nimage);
-    for (i = N = 0; i < Nfile; i++) {
-      // XXX we should not drop all input images if one fails
-      if (i > 0) {
-        f = fopen (file[i], "r");
-        if (f == NULL) {
-          fprintf (stderr, "can't read header for %s, skipping\n", file[i]);
-          continue;
-        }
-        gfits_fread_header (f, header);
+
+  for (i = 0; i < Nfile; i++) {
+    f = fopen (file[i], "r");
+    if (f == NULL) {
+      fprintf (stderr, "can't read file %s, skipping\n", file[i]);
+      continue;
+    }
+
+    headers = LoadHeaders (f, &mode, &Nheaders);
+    headerSets = MatchHeaders (&extsize, &NheaderSets, mode, headers, Nheaders);
+    if (headerSets == NULL) {
+      fprintf (stderr, "ERROR: can't read headers for %s\n", file[i]);
+      continue;
+    }
+    if (NheaderSets == 0) {
+      fprintf (stderr, "no object data in file %s, skipping\n", file[i]);
+      continue;
+    }
+    if (VERBOSE) fprintf (stderr, "file %s has %d headers, including %d images\n", file[i], Nheaders, NheaderSets);
+
+    /* supplied photcode is incompatible with multi-chip images */
+    if ((NheaderSets > 1) && photcode) {
+      fprintf (stderr, "ERROR: photcode cannot be supplied to multi-chip images -- manually adjust the headers\n");
+      exit (1);
+    }
+
+    LoadData (f, file[i], images, Nimages, &stars, Nstars, headers, extsize, headerSets, NheaderSets);
+  }
+
+  if (*Nimages == 0) {
+    if (Nfile == 1) 
+      fprintf (stderr, "no valid image data in any of these files, giving up\n");
+    else 
+      fprintf (stderr, "no valid image data in this file, giving up\n");
+    exit (0);
+  }
+
+  return stars;
+}
+
+// load all of the headers, jump in file to skip data segments
+Header **LoadHeaders (FILE *f, int *mode, int *Nheaders) {
+
+  int i, status, Nskip, NHEADERS;
+  Header **headers;
+
+  /* we need to examine the extensions to determine the headers and the data */
+  NHEADERS = 10;
+  ALLOCATE (headers, Header *, NHEADERS);
+
+  // load all headers into memory
+  for (i = 0;; i++) {
+    ALLOCATE (headers[i], Header, 1);
+    status = gfits_fread_header (f, headers[i]);
+    if (!status) { 
+      *Nheaders = i;
+      return (headers);
+    }
+
+    // check the mode for this file
+    if (i == 0) {
+      *mode = GetFileMode (headers[0]);
+      if ((*mode == SIMPLE_CMP) || (*mode == MOSAIC_CMP)) {
+        *Nheaders = i;
+        return (headers);
       }
-
-      if (!ReadImageHeader (header, &image[N], photcode)) {
-        fprintf (stderr, "skipping %s\n", file[i]);
-        continue;
-      }
-      
-      /* find image rootname */
-      name = filebasename (file[i]);
-      snprintf (image[N].name, 64, name);
-      free (name);
-    
-      switch (mode) {
-        case SIMPLE_CMP:
-        case MOSAIC_CMP:
-          inStars = ReadStarsTEXT (f, &image[N].nstar);
-          inStars = FilterStars (inStars, &image[N]);
-          stars = MergeStars (stars, Nstars, inStars, image[N].nstar);
-          break;
-
-        case SIMPLE_CMF:
-        case MOSAIC_CMF:
-          inStars = ReadStarsFITS (f, header, NULL, &image[N].nstar);
-          inStars = FilterStars (inStars, &image[N]);
-          stars = MergeStars (stars, Nstars, inStars, image[N].nstar);
-          break;
-
-        case MOSAIC_PHU:
-          inStars = NULL;
-          NinStars = 0;
-          break;
-      }
-      fclose (f);
-      gfits_free_header (header);
-      N++;
-    }
-    *Nimages = N;
-    *images = image;
-    
-    if (N == 0) {
-      // XXX how do we track this error?
-      fprintf (stderr, "no valid image data in %s, giving up\n", filename);
-      exit (0);
-    }
-
-    return stars;
-  }
-    
-  /* we have a multi-chip image */
-
-  /* supplied photcode is incompatible with multi-chip images */
-  if (photcode) {
-    fprintf (stderr, "ERROR: photcode cannot be supplied to multi-chip images -- manually adjust the headers\n");
-    exit (1);
-  }
-
-  /* we need to examine the extensions to determine the headers and the data */
-  NHEADER = 10;
-  ALLOCATE (headers, Header *, NHEADER);
-
-  // the first header is already loaded
-  headers[0] = header;
-  Nskip = gfits_matrix_size (header);
-  fseek (f, Nskip, SEEK_CUR); 
-
-  // load all headers into memory
-  done = FALSE;
-  for (i = 1; !done; i++) {
-      ALLOCATE (headers[i], Header, 1);
-      status = gfits_fread_header (f, headers[i]);
-      if (!status) { 
-          done = TRUE;
-      } else {
-          Nskip = gfits_matrix_size (headers[i]);
-          fseek (f, Nskip, SEEK_CUR); 
-      }
-      if (i == NHEADER - 1) {
-          NHEADER += 10;
-          REALLOCATE (headers, Header *, NHEADER);
-      }
-  }
-  Nheader = i - 1; /* we failed on the last loop */
-    
-  // space to store the images, indexes to the matching headers
+    }
+
+    // advance to the next header
+    Nskip = gfits_matrix_size (headers[i]);
+    fseek (f, Nskip, SEEK_CUR); 
+    if (i == NHEADERS - 1) {
+      NHEADERS += 10;
+      REALLOCATE (headers, Header *, NHEADERS);
+    }
+  }
+}
+
+HeaderSet *MatchHeaders (int **extsize, int *nimage, int mode, Header **headers, int Nheaders) {
+
+  int i, j, Nimage, NIMAGE;
+  char extname[80], exttype[80], exthead[80];
+  HeaderSet *headerSets;
+
+  ALLOCATE (extsize[0], int, Nheaders);
+
   Nimage = 0;
-  NIMAGE = Nheader;
-  ALLOCATE (image, Image, NIMAGE);
-  ALLOCATE (exthead, char *, NIMAGE);
-  ALLOCATE (extdata, char *, NIMAGE);
-  ALLOCATE (exttype, char *, NIMAGE);
-  ALLOCATE (extnum_head, int, NIMAGE);
-  ALLOCATE (extnum_data, int, NIMAGE);
-  ALLOCATE (extsize, int, Nheader);
+  NIMAGE = 10;
+  ALLOCATE (headerSets, HeaderSet, NIMAGE);
+
+  // what is the mode of the first header (ie, do we have a PHU DIS image?)
+  mode = GetFileMode (headers[0]);
 
   if (mode == MOSAIC_MEF) {
-      exthead[Nimage] = strcreate ("PHU");
-      extdata[Nimage] = strcreate ("NONE");
-      extnum_data[Nimage] = -1;
-      extnum_head[Nimage] = 0;
-      Nimage ++;
+    headerSets[Nimage].exthead     = strcreate ("PHU");
+    headerSets[Nimage].extdata     = strcreate ("NONE");
+    headerSets[Nimage].extnum_data = -1;
+    headerSets[Nimage].extnum_head =  0;
+    Nimage ++;
   }
 
   // now examine the headers, count the table entries, find corresponding headers
-  for (i = 0; i < Nheader; i++) {
-      extsize[i] = headers[i][0].size + gfits_matrix_size (headers[i]);
-      gfits_scan (headers[i], "EXTTYPE", "%s", 1, tmpword);
-
-      if (!strcmp (tmpword, "SMPDATA")   ||  
-          !strcmp (tmpword, "PS1_DEV_0") ||  
-          !strcmp (tmpword, "PS1_DEV_1")) {
-
-          exttype[Nimage] = strcreate (tmpword);
-          gfits_scan (headers[i], ExtnameKeyword, "%s", 1, tmpword);
-          extdata[Nimage] = strcreate (tmpword);
-          gfits_scan (headers[i], "EXTHEAD", "%s", 1, tmpword);
-          exthead[Nimage] = strcreate (tmpword);
-          extnum_data[Nimage] = i;
-          extnum_head[Nimage] = -1;
-          // find the matching exthead entry
-          for (j = 0; j < Nheader; j++) {
-            if (gfits_scan (headers[j], ExtnameKeyword, "%s", 1, tmpword)) {
-              if (!strcmp (tmpword, exthead[Nimage])) {
-                extnum_head[Nimage] = j;
-              }
-            }
-          }
-          // skip or crash on table with missing matching header?
-          if (extnum_head[Nimage] == -1) {
-              fprintf (stderr, "ERROR: can't read header for %s\n", file[0]);
-              exit (1);
-          }
-          Nimage ++;
-      }
-  }
+  for (i = 0; i < Nheaders; i++) {
+    if (mode == SIMPLE_CMP) {
+      extsize[0][i] = headers[i][0].size;
+    } else {
+      extsize[0][i] = headers[i][0].size + gfits_matrix_size (headers[i]);
+    }
+    gfits_scan (headers[i], "EXTTYPE", "%s", 1, exttype);
+
+    if (!strcmp (exttype, "SMPDATA")) goto keep;
+    if (!strcmp (exttype, "PS1_DEV_0")) goto keep;
+    if (!strcmp (exttype, "PS1_DEV_1")) goto keep;
+    continue;
+
+  keep:
+    headerSets[Nimage].exttype = strcreate (exttype);
+
+    gfits_scan (headers[i], ExtnameKeyword, "%s", 1, extname);
+    gfits_scan (headers[i], "EXTHEAD", "%s", 1, exthead);
+
+    headerSets[Nimage].extdata     = strcreate (extname);
+    headerSets[Nimage].exthead     = strcreate (exthead);
+    headerSets[Nimage].extnum_data = i;
+    headerSets[Nimage].extnum_head = -1;
+
+    // find the matching exthead entry
+    for (j = 0; j < Nheaders; j++) {
+      if (!gfits_scan (headers[j], ExtnameKeyword, "%s", 1, extname)) continue;
+      if (strcmp (extname, headerSets[Nimage].exthead)) continue;
+      headerSets[Nimage].extnum_head = j;
+      break;
+    }
+
+    // skip or crash on table with missing matching header?
+    if (headerSets[Nimage].extnum_head == -1) {
+      return NULL;
+    }
+    Nimage ++;
+    if (Nimage == NIMAGE) {
+      NIMAGE += 10;
+      REALLOCATE (headerSets, HeaderSet, NIMAGE);
+    }
+  }
+
   // some old format files did not write EXTTYPE.  they have a single table in the first
   // extension matched to the header in the PHU
   if (Nimage == 0) {
-      extsize[0] = headers[0][0].size + gfits_matrix_size (headers[0]);
-      extsize[1] = headers[1][0].size + gfits_matrix_size (headers[1]);
-      gfits_scan (headers[1], ExtnameKeyword, "%s", 1, tmpword);
-      if (!strcmp (tmpword, "SMPFILE")) {
-          extdata[Nimage] = strcreate (tmpword);
-          exttype[Nimage] = strcreate ("SMPDATA");
-          exthead[Nimage] = strcreate ("PHU");
-          extnum_head[Nimage] = 0;
-          extnum_data[Nimage] = 1;
-          Nimage = 1;
+    extsize[0][0] = headers[0][0].size + gfits_matrix_size (headers[0]);
+    extsize[0][1] = headers[1][0].size + gfits_matrix_size (headers[1]);
+    gfits_scan (headers[1], ExtnameKeyword, "%s", 1, extname);
+    if (!strcmp (extname, "SMPFILE")) {
+      headerSets[Nimage].extdata     = strcreate (extname);
+      headerSets[Nimage].exttype     = strcreate ("SMPDATA");
+      headerSets[Nimage].exthead     = strcreate ("PHU");
+      headerSets[Nimage].extnum_head = 0;
+      headerSets[Nimage].extnum_data = 1;
+      Nimage = 1;
+    }
+  }
+  
+  *nimage = Nimage;
+  return (headerSets);
+}
+
+// examine the header sets and set the Image entries for the the valid images
+int LoadData (FILE *f, char *file, Image **images, int *nvalid, Stars **stars, int *Nstars, Header **headers, int *extsize, HeaderSet *headerSets, int Nimages) {
+
+  char *name;
+  int i, j, Nvalid, Nhead, Ndata, Nskip;
+  Stars *inStars;
+
+  if (images[0] == NULL) {
+    Nvalid = 0;
+    NVALID = 10;
+    ALLOCATE (images[0], Image, NVALID);
+  } else {
+    Nvalid = *nvalid;
+    NVALID = Nvalid + 10;
+    REALLOCATE (images[0], Image, NVALID);
+  }    
+
+  // find image rootname
+  name = filebasename (file);
+
+  // now run through the images, interpret the headers and read the stars
+  for (i = 0; i < Nimages; i++) {
+    Nhead = headerSets[i].extnum_head;
+
+    if (VERBOSE) fprintf (stderr, "reading header for %s (%s)\n", headerSets[i].exthead, headerSets[i].extdata);
+    if (!ReadImageHeader (headers[Nhead], &images[0][Nvalid], 0)) {
+      fprintf (stderr, "skipping %s\n", headerSets[i].exthead);
+      continue;
+    }
+
+    // XXX use something to set the chip name? EXTNAME?
+    if (!strcmp(headerSets[i].exthead, "PHU") && (Nimages == 1)) {
+      snprintf (images[0][Nvalid].name, 64, "%s", name);
+    } else {
+      snprintf (images[0][Nvalid].name, 64, "%s[%s]", name, headerSets[i].exthead);
+    }
+
+    // skip the table if there is no data segment (eg, mosaic WRP image)
+    if (!strcmp(headerSets[i].extdata, "NONE")) {
+      Nvalid++;
+      if (Nvalid == NVALID) {
+        NVALID += 10;
+        REALLOCATE (images[0], Image, NVALID);
       }
-  }
-  if (Nimage == 0) Shutdown ("no object data in file");
-    
-  if (VERBOSE) fprintf (stderr, "file %s has %d headers, including %d images\n", file[0], Nheader, Nimage);
-
-  /* find image rootname */
-  name = filebasename (file[0]);
-
-  // now run through the images, interpret the headers and read the stars
-  for (i = N = 0; i < Nimage; i++) {
-      Nhead = extnum_head[i];
-
-      if (VERBOSE) fprintf (stderr, "reading header for %s (%s)\n", exthead[i], extdata[i]);
-      if (!ReadImageHeader (headers[Nhead], &image[N], 0)) {
-          fprintf (stderr, "skipping %s\n", exthead[i]);
-          continue;
-      }
-
-      // XXX use something to set the chip name? EXTNAME?
-      if (!strcmp(exthead[i], "PHU") && (Nimage == 1)) {
-        snprintf (image[N].name, 64, "%s", name);
-      } else {
-        snprintf (image[N].name, 64, "%s[%s]", name, exthead[i]);
-      }
-
-      // skip the table if there is not data segment (eg, mosaic WRP image)
-      if (!strcmp(extdata[i], "NONE")) {
-          N++;
-          continue;
-      }
-
-      // advance the pointer to the start of the corresponding table block
-      Ndata = extnum_data[i];
-      Nskip = 0;
-      for (j = 0; j < Ndata; j++) {
-          Nskip += extsize[j];
-      }
-      fseek (f, Nskip, SEEK_SET); 
+      continue;
+    }
+
+    // advance the pointer to the start of the corresponding table block
+    Ndata = headerSets[i].extnum_data;
+    Nskip = 0;
+    for (j = 0; j < Ndata; j++) {
+      Nskip += extsize[j];
+    }
+    fseek (f, Nskip, SEEK_SET); 
          
-      inStars = ReadStarsFITS (f, headers[Nhead], headers[Ndata], &image[N].nstar);
-      inStars = FilterStars (inStars, &image[N]);
-      stars = MergeStars (stars, Nstars, inStars, image[N].nstar);
-      N++;
+    inStars = ReadStarsFITS (f, headers[Nhead], headers[Ndata], &images[0][Nvalid].nstar);
+    inStars = FilterStars (inStars, &images[0][Nvalid]);
+    *stars = MergeStars (*stars, Nstars, inStars, images[0][Nvalid].nstar);
+    Nvalid++;
   }
   free (name);
-  *Nimages = N;
-  *images = image;
-
-  if (N == 0) {
-    fprintf (stderr, "no valid image data in %s, giving up\n", filename);
-    exit (0);
-  }
-
-  return stars;
-}
+  *nvalid = Nvalid;
+  return (TRUE);
+}
+