# include "addstar.h"

// XXX where and when should I call RegisterMosaic??

Stars *gstars (char *filename, int *Nstars, Image **images, int *Nimages, int photcode) {

  int i, j, Nfile, Nheader, NHEADER, Nimage, NIMAGE;
  int Nskip, Nhead, Ndata, done, status, mode, NinStars;
  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;

  // parse the filename as a glob
  globList.gl_offs = 0;
  glob (filename, 0, NULL, &globList);

  // if the glob does not match, save the literal word:
  // otherwise save all glob matches
  if (globList.gl_pathc == 0) {
    Nfile = 1;
    ALLOCATE (file, char *, Nfile);
    file[0] = strcreate (filename);
  } else {
    Nfile = globList.gl_pathc;
    ALLOCATE (file, char *, Nfile);
    for (i = 0; i < Nfile; i++) {
      file[i] = strcreate (globList.gl_pathv[i]);
    }
  }

  // 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);

  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 = 0; i < Nfile; i++) {
      if (i > 0) {
	f = fopen (file[i], "r");
	if (f == NULL) {
	  fprintf (stderr, "ERROR: can't read header for %s\n", file[i]);
	  exit (1);
	}
	gfits_fread_header (f, header);
      }

      ReadImageHeader (header, &image[i], photcode);
    
      switch (mode) {
	case SIMPLE_CMP:
	case MOSAIC_CMP:
	  inStars = ReadStarsTEXT (f, &image[i].nstar);
	  inStars = FilterStars (inStars, &image[i]);
	  stars = MergeStars (stars, Nstars, inStars, image[i].nstar);
	  break;

	case SIMPLE_CMF:
	case MOSAIC_CMF:
	  inStars = ReadStarsFITS (f, header, NULL, &image[i].nstar);
	  inStars = FilterStars (inStars, &image[i]);
	  stars = MergeStars (stars, Nstars, inStars, image[i].nstar);
	  break;

	case MOSAIC_PHU:
	  inStars = NULL;
	  NinStars = 0;
	  break;
      }
      fclose (f);
      gfits_free_header (header);
    }
    *Nimages = Nimage;
    *images = image;
    
    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
  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);

  if (mode == MOSAIC_MEF) {
      exthead[Nimage] = strcreate ("PHU");
      extdata[Nimage] = strcreate ("NONE");
      extnum_data[Nimage] = -1;
      extnum_head[Nimage] = 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, "PS1DATA")) {
	  exttype[Nimage] = strcreate (tmpword);
	  gfits_scan (headers[i], "EXTNAME", "%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], "EXTNAME", "%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 ++;
      }
  }
  if (VERBOSE) fprintf (stderr, "file %s has %d headers, including %d images\n", file[0], Nheader, Nimage);

  // now run through the images, interpret the headers and read the stars
  for (i = 0; i < Nimage; i++) {
      Nhead = extnum_head[i];

      if (VERBOSE) fprintf (stderr, "reading header for %s (%s)\n", exthead[i], extdata[i]);
      ReadImageHeader (headers[Nhead], &image[i], 0);

      if (!strcmp(exthead[i], "PHU")) 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); 
	 
      inStars = ReadStarsFITS (f, headers[Nhead], headers[Ndata], &image[i].nstar);
      inStars = FilterStars (inStars, &image[i]);
      stars = MergeStars (stars, Nstars, inStars, image[i].nstar);
  }
  *Nimages = Nimage;
  *images = image;

  return stars;
}
