Index: /trunk/Ohana/src/getstar/Makefile
===================================================================
--- /trunk/Ohana/src/getstar/Makefile	(revision 12773)
+++ /trunk/Ohana/src/getstar/Makefile	(revision 12774)
@@ -18,5 +18,10 @@
 
 getstar: $(BIN)/getstar.$(ARCH)
-install: $(DESTBIN)/getstar
+getstar.install: $(DESTBIN)/getstar
+
+dvoImageOverlaps: $(BIN)/dvoImageOverlaps.$(ARCH)
+dvoImageOverlaps.install: $(DESTBIN)/dvoImageOverlaps
+
+install: getstar.install dvoImageOverlaps.install
 
 GETSTAR = \
@@ -30,2 +35,19 @@
 $(GETSTAR): $(INC)/getstar.h
 $(BIN)/getstar.$(ARCH): $(GETSTAR)
+
+OVERLAPS = \
+$(SRC)/dvoImageOverlaps.$(ARCH).o 	 \
+$(SRC)/args_overlaps.$(ARCH).o		 \
+$(SRC)/ConfigInit_overlaps.$(ARCH).o 	 \
+$(SRC)/Shutdown.$(ARCH).o	         \
+$(SRC)/SetSignals.$(ARCH).o	         \
+$(SRC)/ReadImageFiles.$(ARCH).o          \
+$(SRC)/ReadImageHeader.$(ARCH).o         \
+$(SRC)/ListImageOverlaps.$(ARCH).o       \
+$(SRC)/GetFileMode.$(ARCH).o             \
+$(SRC)/edge_check.$(ARCH).o              \
+$(SRC)/opening_angle.$(ARCH).o           \
+$(SRC)/MatchImages.$(ARCH).o
+
+$(OVERLAPS): $(INC)/dvoImageOverlaps.h
+$(BIN)/dvoImageOverlaps.$(ARCH): $(OVERLAPS)
Index: /trunk/Ohana/src/getstar/include/dvoImageOverlaps.h
===================================================================
--- /trunk/Ohana/src/getstar/include/dvoImageOverlaps.h	(revision 12774)
+++ /trunk/Ohana/src/getstar/include/dvoImageOverlaps.h	(revision 12774)
@@ -0,0 +1,36 @@
+# include <ohana.h>
+# include <dvo.h>
+# include <signal.h>
+# include <glob.h>
+
+enum {NONE, SIMPLE_CMP, SIMPLE_CMF, SIMPLE_MEF, MOSAIC_CMP, MOSAIC_CMF, MOSAIC_MEF, MOSAIC_PHU};
+
+int       VERBOSE;
+
+char OUTPUT[256];
+char GSCFILE[256];
+char CATDIR[256];
+char CATMODE[16];    /* raw, mef, split, mysql */
+char CATFORMAT[16];  /* internal, elixir, loneos, panstarrs */
+char SKY_TABLE[256];
+int  SKY_DEPTH;
+char   ImageCat[256];
+Coords *MOSAIC;         // carries the mosaic into ReadImageHeader
+
+int WITH_PHU;
+int SOLO_PHU;
+
+int  args_overlaps    	 PROTO((int argc, char **argv));
+int  ConfigInit_overlaps PROTO((int *argc, char **argv));
+int  Shutdown         PROTO((char *format, ...));
+void TrapSignal       PROTO((int sig));
+void SetProtect       PROTO((int mode));
+int  SetSignals       PROTO(());
+
+Image *ReadImageFiles (char *filename, int *Nimages);
+int ReadImageHeader (Header *header, Image *image);
+int *MatchImage (Image *dbImages, int NdbImages, Image *image, int *Nmatch);
+int ListImageOverlaps (Image *dbImages, Image *image, int *matches, int Nmatches);
+
+int GetFileMode (Header *header);
+int edge_check (double *x1, double *y1, double *x2, double *y2);
Index: /trunk/Ohana/src/getstar/src/ConfigInit_overlaps.c
===================================================================
--- /trunk/Ohana/src/getstar/src/ConfigInit_overlaps.c	(revision 12774)
+++ /trunk/Ohana/src/getstar/src/ConfigInit_overlaps.c	(revision 12774)
@@ -0,0 +1,46 @@
+# include "dvoImageOverlaps.h"
+
+int ConfigInit_overlaps (int *argc, char **argv) {
+
+  char *config, *file;
+  char CatdirPhotcodeFile[256];
+  char MasterPhotcodeFile[256];
+
+  /*** load configuration info ***/
+  file = SelectConfigFile (argc, argv, "ptolemy");
+  config = LoadConfigFile (file);
+  if (config == (char *) NULL) {
+    fprintf (stderr, "ERROR: can't find configuration file %s\n", file);
+    if (file != (char *) NULL) free (file);
+    exit (1);
+  }
+  if (VERBOSE) fprintf (stderr, "loaded config file: %s\n", file);
+
+  ScanConfig (config, "GSCFILE",                "%s", 0, GSCFILE);
+  ScanConfig (config, "CATDIR",                 "%s", 0, CATDIR);
+  ScanConfig (config, "CATMODE",                "%s",  0, CATMODE);
+  ScanConfig (config, "CATFORMAT",              "%s",  0, CATFORMAT);
+  ScanConfig (config, "PHOTCODE_FILE",         	"%s",  0, MasterPhotcodeFile);
+  if (!ScanConfig (config, "SKY_DEPTH",         "%d",  0, &SKY_DEPTH)) {
+    SKY_DEPTH = 2;
+  }
+  if (!ScanConfig (config, "SKY_TABLE",         "%s",  0, SKY_TABLE)) {
+    SKY_TABLE[0] = 0;
+  }
+  sprintf (ImageCat, "%s/Images.dat", CATDIR);
+
+  if (*CATMODE == 0) strcpy (CATMODE, "RAW");
+  if (*CATFORMAT == 0) strcpy (CATFORMAT, "ELIXIR");
+
+  /* XXX this does not yet write out the master photcode table */
+  sprintf (CatdirPhotcodeFile, "%s/Photcodes.dat", CATDIR);
+  if (!LoadPhotcodes (CatdirPhotcodeFile, MasterPhotcodeFile)) {
+    fprintf (stderr, "error loading photcode table %s or master file %s\n", CatdirPhotcodeFile, MasterPhotcodeFile);
+    exit (1);
+  }
+
+  free (config);
+  free (file);
+
+  return (TRUE);
+}
Index: /trunk/Ohana/src/getstar/src/GetFileMode.c
===================================================================
--- /trunk/Ohana/src/getstar/src/GetFileMode.c	(revision 12774)
+++ /trunk/Ohana/src/getstar/src/GetFileMode.c	(revision 12774)
@@ -0,0 +1,45 @@
+# include "dvoImageOverlaps.h"
+
+// examine the PHU of this file and determine the file mode
+int GetFileMode (Header *header) {
+
+  char ctype[80];
+  int Naxis, Nx, Ny;
+  int simple, extend, haveNaxis, haveNx, haveNy, haveCTYPE;
+
+  gfits_scan (header, "SIMPLE", "%t", 1, &simple);
+  haveNaxis = gfits_scan (header, "NAXIS",  "%d", 1, &Naxis);
+  haveNx    = gfits_scan (header, "NAXIS1", "%d", 1, &Nx);
+  haveNy    = gfits_scan (header, "NAXIS2", "%d", 1, &Ny);
+  haveCTYPE = gfits_scan (header, "CTYPE1", "%s", 1, ctype);
+
+  gfits_scan (header, "EXTEND", "%t", 1, &extend);
+    
+  if ((Naxis == 2) || !simple) {
+    if (!strcmp (&ctype[4], "-WRP")) {
+      return MOSAIC_CMP;
+    }
+    return SIMPLE_CMP;
+  }
+
+  if (!extend && strcmp (&ctype[4], "-DIS")) {
+    if (!strcmp (&ctype[4], "-WRP")) {
+      return MOSAIC_CMF;
+    }
+    return SIMPLE_CMF;
+  }
+
+  if (!extend && !strcmp (&ctype[4], "-DIS")) {
+    return MOSAIC_PHU;
+  }
+
+  if (extend && strcmp (&ctype[4], "-DIS")) {
+    return SIMPLE_MEF;
+  }
+
+  if (extend && !strcmp (&ctype[4], "-DIS")) {
+    return MOSAIC_MEF;
+  }
+
+  return (NONE);
+}
Index: /trunk/Ohana/src/getstar/src/ListImageOverlaps.c
===================================================================
--- /trunk/Ohana/src/getstar/src/ListImageOverlaps.c	(revision 12774)
+++ /trunk/Ohana/src/getstar/src/ListImageOverlaps.c	(revision 12774)
@@ -0,0 +1,14 @@
+# include "dvoImageOverlaps.h"
+
+/* given image, find catalog images which overlap it */
+int ListImageOverlaps (Image *dbImages, Image *image, int *matches, int Nmatches) {
+  
+  int i, N;
+
+  for (i = 0; i < Nmatches; i++) {
+    N = matches[i];
+    fprintf (stdout, "%s  :  %s\n", image[0].name, dbImages[N].name);
+  }
+  
+  return (TRUE);
+}
Index: /trunk/Ohana/src/getstar/src/MatchImages.c
===================================================================
--- /trunk/Ohana/src/getstar/src/MatchImages.c	(revision 12774)
+++ /trunk/Ohana/src/getstar/src/MatchImages.c	(revision 12774)
@@ -0,0 +1,151 @@
+# include "dvoImageOverlaps.h"
+# ifndef FLT_MAX
+# define FLT_MAX 1e32
+# endif
+
+/* given image, find catalog images which overlap it */
+int *MatchImage (Image *dbImages, int NdbImages, Image *image, int *Nmatch) {
+  
+  int i, j, k, N, addtolist, status;
+  int NMATCH, nmatch, *match;
+  Coords tcoords;
+  double r, d, x, y;
+  double Xi[4], Yi[4], Xo[4], Yo[4];  /* image and original corners */
+  double Xmin, Xmax, Ymin, Ymax;
+  double xmin, xmax, ymin, ymax;
+
+  if (!WITH_PHU && !strcmp (&image[0].coords.ctype[4], "-DIS")) return FALSE;
+  if ( SOLO_PHU &&  strcmp (&image[0].coords.ctype[4], "-DIS")) return FALSE;
+
+  /* project onto rectilinear grid with 1 arcsec pixels */
+  /* we keep the original crpix1,2 and crref1,2 */
+  /* for mosaic astrometry, the grid should be w.r.t. the tangent-plane, not chip coords */
+  if (!strcmp (&image[0].coords.ctype[4], "-WRP")) {
+    // if the input image is a mosaic, we need to register that mosaic for the calculations below
+    if (MOSAIC == NULL) {
+      fprintf (stderr, "no mosaic for WRP image (use -mosaic)\n");
+      return (FALSE);
+    }
+    RegisterMosaic (MOSAIC);
+
+    tcoords = MOSAIC[0];
+    tcoords.cdelt1 = tcoords.cdelt2 = 1.0 / 3600.0;
+    tcoords.pc1_1 = tcoords.pc2_2 = 1.0;
+    tcoords.pc1_2 = tcoords.pc2_1 = 0.0;
+    tcoords.Npolyterms = 1;
+    strcpy (tcoords.ctype, "RA---TAN");
+    /* register so image->sky conversions below have correct mosaic */
+  } else {
+    tcoords = image[0].coords;
+    tcoords.cdelt1 = tcoords.cdelt2 = 1.0 / 3600.0;
+    tcoords.pc1_1 = tcoords.pc2_2 = 1.0;
+    tcoords.pc1_2 = tcoords.pc2_1 = 0.0;
+    strcpy (tcoords.ctype, "RA---TAN");
+  }
+
+  /* define original corners */
+  SetImageCorners (Xo, Yo, &image[0]);
+  
+  Ymin = Xmin = +FLT_MAX;
+  Ymax = Xmax = -FLT_MAX;
+  for (j = 0; j < 4; j++) {
+    /* XY-to_RD is two-level if ctype == WRP */
+    XY_to_RD (&r, &d, Xo[j], Yo[j], &image[0].coords);
+    RD_to_XY (&Xo[j], &Yo[j], r, d, &tcoords);
+    Xmin = MIN (Xmin, Xo[j]);
+    Xmax = MAX (Xmax, Xo[j]);
+    Ymin = MIN (Ymin, Yo[j]);
+    Ymax = MAX (Ymax, Yo[j]);
+  }
+
+  /* match represents the subset of overlapping images */
+  nmatch = 0;
+  NMATCH = 20;
+  ALLOCATE (match, int, NMATCH);
+
+  /* setup links for mosaic WRP and DIS entries */
+  BuildChipMatch (dbImages, NdbImages);
+
+  /* run through image table and search for overlaps
+     also define the vtable entries for the images we keep  */
+
+  for (i = 0; i < NdbImages; i++) {
+
+    if (!WITH_PHU && !strcmp (&dbImages[i].coords.ctype[4], "-DIS")) continue;
+    if ( SOLO_PHU &&  strcmp (&dbImages[i].coords.ctype[4], "-DIS")) continue;
+
+    /* if any of these images are WRP images, need to find matching DIS */
+    if (!FindMosaicForImage (dbImages, NdbImages, i)) continue;
+
+    /* define image corners */
+    SetImageCorners (Xi, Yi, &dbImages[i]);
+    // Xi[4] = Xi[0]; Yi[4] = Yi[0];
+
+    /* transform to tcoords, skip corners off image */
+    /*** XXX this will fail for very large images which extend beyond 180deg */
+    ymin = xmin = +FLT_MAX;
+    ymax = xmax = -FLT_MAX;
+    for (j = N = 0; j < 4; j++) {
+      status = XY_to_RD (&r, &d, Xi[j], Yi[j], &dbImages[i].coords);
+      if (!status) continue;
+      status = RD_to_XY (&Xi[N], &Yi[N], r, d, &tcoords);
+      if (!status) continue;
+      xmin = MIN (xmin, Xi[N]);
+      xmax = MAX (xmax, Xi[N]);
+      ymin = MIN (ymin, Yi[N]);
+      ymax = MAX (ymax, Yi[N]);
+      N++;
+    }
+
+    /* check if one corner of dbImages[i] is inside image[0] */
+    for (j = 0; j < N; j++) {
+      addtolist = TRUE;
+      addtolist &= (Xi[j] >= Xmin);
+      addtolist &= (Xi[j] <= Xmax);
+      addtolist &= (Yi[j] >= Ymin);
+      addtolist &= (Yi[j] <= Ymax);
+      if (addtolist) goto addtolist;
+    }
+
+    /* or else, check if one corner of image[0] is inside dbImages[i] */
+    for (j = 0; j < 4; j++) {
+      addtolist = TRUE;
+      addtolist &= (Xo[j] >= xmin);
+      addtolist &= (Xo[j] <= xmax);
+      addtolist &= (Yo[j] >= ymin);
+      addtolist &= (Yo[j] <= ymax);
+      if (addtolist) goto addtolist;
+    }
+
+    // no match, skip this dbImage
+    continue;
+
+  addtolist:
+    match[nmatch] = i;
+    nmatch ++;
+    if (nmatch == NMATCH) {
+      NMATCH += 20;
+      REALLOCATE (match, int, NMATCH);
+    }
+  }
+  
+  if (VERBOSE) fprintf (stderr, "found %d overlapping images\n", nmatch);
+
+  *Nmatch = nmatch;
+  return (match);
+}
+  
+void SetImageCorners (double *X, double *Y, Image *image) {
+
+  if (!strcmp(&image[0].coords.ctype[4], "-DIS")) {
+    X[0] = -0.5*image[0].NX; Y[0] = -0.5*image[0].NY;
+    X[1] = +0.5*image[0].NX; Y[1] = -0.5*image[0].NY;
+    X[2] = +0.5*image[0].NX; Y[2] = +0.5*image[0].NY;
+    X[3] = -0.5*image[0].NX; Y[3] = +0.5*image[0].NY;
+  } else {
+    X[0] = 0;           Y[0] = 0;
+    X[1] = image[0].NX; Y[1] = 0;
+    X[2] = image[0].NX; Y[2] = image[0].NY;
+    X[3] = 0;           Y[3] = image[0].NY;
+  }
+}
Index: /trunk/Ohana/src/getstar/src/ReadImageFiles.c
===================================================================
--- /trunk/Ohana/src/getstar/src/ReadImageFiles.c	(revision 12774)
+++ /trunk/Ohana/src/getstar/src/ReadImageFiles.c	(revision 12774)
@@ -0,0 +1,225 @@
+# include "dvoImageOverlaps.h"
+
+Image *ReadImageFiles (char *filename, int *Nimages) {
+
+  int i, j, N, Nfile, Nheader, NHEADER, Nimage, NIMAGE;
+  int Nskip, Nhead, Ndata, done, status, mode, NinStars;
+  char **file, *name;
+  FILE *f;
+  glob_t globList;
+  char **exthead, **extdata, **exttype, tmpword[80];
+  int *extnum_head, *extnum_data, *extsize;
+  Header *header, **headers;
+  Image *image;
+
+  // 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);
+
+  /*** 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++) {
+      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);
+      }
+
+      if (!ReadImageHeader (header, &image[N])) {
+	fprintf (stderr, "skipping %s\n", file[i]);
+	continue;
+      }
+      
+      /* find image rootname */
+      name = filebasename (file[i]);
+      snprintf (image[N].name, 32, name);
+      free (name);
+
+      fclose (f);
+      gfits_free_header (header);
+      N++;
+    }
+    *Nimages = N;
+    
+    if (N == 0) {
+      fprintf (stderr, "ERROR: no valid image data in %s, giving up\n", filename);
+      exit (1);
+    }
+    return image;
+  }
+    
+  /* we have a multi-chip image */
+
+  /* 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, "PS1_DEV_0")) {
+
+	  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 ++;
+      }
+  }
+
+  // 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], "EXTNAME", "%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;
+      }
+  }
+  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];
+
+      // XXX do I need to advance the file pointer, or does ReadImageHeader do this?
+      if (VERBOSE) fprintf (stderr, "reading header for %s (%s)\n", exthead[i], extdata[i]);
+      if (!ReadImageHeader (headers[Nhead], &image[N])) {
+	  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, 32, "%s", name);
+      } else {
+	snprintf (image[N].name, 32, "%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); 
+      N++;
+  }
+  free (name);
+  *Nimages = N;
+
+  if (N == 0) {
+    fprintf (stderr, "ERROR: no valid image data in %s, giving up\n", filename);
+    exit (1);
+  }
+
+  return image;
+}
Index: /trunk/Ohana/src/getstar/src/ReadImageHeader.c
===================================================================
--- /trunk/Ohana/src/getstar/src/ReadImageHeader.c	(revision 12774)
+++ /trunk/Ohana/src/getstar/src/ReadImageHeader.c	(revision 12774)
@@ -0,0 +1,157 @@
+# include "dvoImageOverlaps.h"
+
+/* read an image header corresponding to a CMF / CMP data block */
+int ReadImageHeader (Header *header, Image *image) {
+
+  int Nastro, ccdnum, hour, min, Nx, Ny, photcode;
+  double tmp, sec, Cerror;
+  char *c, photname[64], line[80];
+
+  /* get astrometry information */
+  if (!GetCoords (&image[0].coords, header)) {
+    fprintf (stderr, "no astrometric solution in header\n");
+    return (FALSE);
+  }
+  // XXX currently, image uses an unsigned short for NX,XY. this is rather restrictive
+  // and needs to be at least checked.
+  gfits_scan (header, "NAXIS1",   "%d", 1, &Nx);
+  gfits_scan (header, "NAXIS2",   "%d", 1, &Ny);
+  if ((Nx < 0) || (Nx > 0xffff)) {
+    fprintf (stderr, "WARNING: NX, NY out of range : image boundary will be wrong\n");
+  }
+  image[0].NX = Nx;
+  image[0].NY = Ny;
+
+  // if (!gfits_scan (header, "TZERO",   "%d",  1, &image[0].tzero)) {
+  // image[0].tzero = parse_time (header);
+  // }
+
+  /* only load astrometry, NAXIS1,2, and time if this is a MOSAIC_PHU (ctype is ....-DIS) */
+  if (!strcmp (&image[0].coords.ctype[4], "-DIS")) {
+    MOSAIC = &image[0].coords;
+    return (TRUE);
+  }
+
+  /* require Nastro > 0 unless or ACCEPT_ASTROM */
+  Nastro = 0;
+  gfits_scan (header, "NASTRO", "%d", 1, &Nastro);
+  if (Nastro == 0) {
+    fprintf (stderr, "bad astrometric solution in header\n");
+    return (FALSE);
+  }
+  if (!strcmp (&image[0].coords.ctype[4], "-WRP")) {
+    if (MOSAIC == NULL) {
+      fprintf (stderr, "no mosaic for WRP image (use -mosaic)\n");
+      return (FALSE);
+    }
+    RegisterMosaic (MOSAIC);
+  } else {
+    /* force image to lie in 0-360 range */
+    while (image[0].coords.crval1 < 0) image[0].coords.crval1 += 360.0;
+    while (image[0].coords.crval1 > 360.0) image[0].coords.crval1 -= 360.0;
+  }
+
+  { 
+    double R, D;
+    /* sanity check on the image coordinates */
+    XY_to_RD (&R, &D, 0.5*Nx, 0.5*Ny, &image[0].coords);
+    if (!finite(R) || !finite(D)) {
+      fprintf (stderr, "corrupted header coordinates, skipping\n");
+      return (FALSE);
+    }
+  }
+    
+  /* CERROR in data file is in arcsec */
+  if (!gfits_scan (header, "CERROR",   "%lf", 1, &tmp)) tmp = 1.0;
+  image[0].cerror = tmp * 50.0;
+ 
+  /* get photcode from header */
+  /*** XXX is photcode needed?? ***/
+  if (!gfits_scan (header, "PHOTCODE", "%s", 1, photname)) {
+    fprintf (stderr, "photcode not supplied in header\n");
+    return (FALSE);
+  }
+  photcode = GetPhotcodeCodebyName (photname);
+  if (photcode == 0) { 
+    fprintf (stderr, "no valid photcode is supplied\n");
+    return (FALSE);
+  }
+  image[0].photcode = photcode;
+
+  // XXX we don't really need these
+  // gfits_scan (header, ExptimeKeyword,  "%lf", 1, &tmp);
+  // image[0].exptime = tmp;
+  
+  /*** why are we no longer using APMIFIT?? ***/
+  tmp = 0;
+  /* gfits_scan (header, "APMIFIT",  "%lf", 1, &tmp); */
+  image[0].apmifit = tmp;
+
+  tmp = 0;
+  /* gfits_scan (header, "dAPMIFIT", "%lf", 1, &tmp); */
+  image[0].dapmifit = tmp;
+
+  tmp = 0;
+  gfits_scan (header, "FLIMIT",   "%lf", 1, &tmp);
+  image[0].detection_limit = tmp * 10.0;
+
+  tmp = 0;
+  gfits_scan (header, "FSATUR",   "%lf", 1, &tmp);
+  image[0].saturation_limit = tmp * 10.0;
+
+  tmp = 0;
+  gfits_scan (header, "FWHM_X",   "%lf", 1, &tmp);
+  image[0].fwhm_x = tmp * 25.0 * image[0].coords.cdelt1 * 3600.0;
+
+  tmp = 0;
+  gfits_scan (header, "FWHM_Y",   "%lf", 1, &tmp);
+  image[0].fwhm_y = tmp * 25.0 * image[0].coords.cdelt1 * 3600.0;
+
+  // XXX we don't really need the following values
+# if (0)  
+  if (STKeyword[0]) {
+    /* get ST (used for airmass calculation) */
+    gfits_scan (header, STKeyword, "%s", 1, line);
+    /* remove ':' characters */
+    for (c = strchr (line, ':'); c != (char *) NULL; c = strchr (line, ':')) { *c = ' '; }
+    sscanf (line, "%d %d %lf", &hour, &min, &sec);
+    image[0].sidtime = hour + min/60.0 + sec/3600.0;
+  } else {
+    double jd;
+    jd = sec_to_jd (image[0].tzero);
+    image[0].sidtime  = ohana_lst (jd, Longitude);
+  }
+  image[0].latitude = Latitude;
+
+  tmp = 0;
+  gfits_scan (header, "TRATE",   "%lf", 1, &tmp);
+  image[0].trate = 10000 * tmp;
+
+  tmp = 0;
+  gfits_scan (header, AirmassKeyword, "%lf", 1, &tmp);
+  image[0].secz = MIN (NO_MAG, tmp);
+
+  if (!gfits_scan (header, CCDNumKeyword, "%d", 1, &ccdnum)) {
+    image[0].ccdnum = 0xff;
+  } else {
+    image[0].ccdnum = ccdnum;
+  }
+
+  gfits_scan (header, "ZERO_PT", "%lf", 1, &ZeroPt);
+  SetZeroPoint (ZeroPt);
+# endif
+
+  /* secz is in units milli-airmass */
+  image[0].Mcal = 0.0;
+  image[0].Xm   = NO_MAG;
+  image[0].code = 0;
+  bzero (image[0].dummy, sizeof(image[0].dummy));
+
+  /* find expected number of stars */
+  if (!gfits_scan (header, "NSTARS", "%d", 1, &image[0].nstar)) {
+    fprintf (stderr, "can't get NSTARS from header\n");
+    return (FALSE);
+  }
+
+  return (TRUE);
+}
Index: /trunk/Ohana/src/getstar/src/args_overlaps.c
===================================================================
--- /trunk/Ohana/src/getstar/src/args_overlaps.c	(revision 12774)
+++ /trunk/Ohana/src/getstar/src/args_overlaps.c	(revision 12774)
@@ -0,0 +1,65 @@
+# include "dvoImageOverlaps.h"
+
+void help () {
+  fprintf (stderr, "USAGE: \n"
+	   "dvoImageOverlaps (image)\n"
+    );
+  exit (2);
+}
+
+int args_overlaps (int argc, char **argv) {
+  
+  int N;
+
+  /* check for help request */
+  if (get_argument (argc, argv, "-help") ||
+      get_argument (argc, argv, "-h")) {
+    help ();
+  }
+
+  /* check for command line options */
+  WITH_PHU = FALSE;
+  if ((N = get_argument (argc, argv, "+phu"))) {
+    WITH_PHU = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+  SOLO_PHU = FALSE;
+  if ((N = get_argument (argc, argv, "-phu"))) {
+    WITH_PHU = TRUE;
+    SOLO_PHU = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+
+  /* check for command line options */
+  VERBOSE = FALSE;
+  if ((N = get_argument (argc, argv, "-v"))) {
+    VERBOSE = TRUE;
+    remove_argument (N, &argc, argv);
+  }
+  /* provide a mosaic for distortion */
+  MOSAIC = NULL;
+  if ((N = get_argument (argc, argv, "-mosaic"))) {
+    Header header;
+    ALLOCATE (MOSAIC, Coords, 1);
+
+    remove_argument (N, &argc, argv);
+    if (!gfits_read_header (argv[N], &header)) {
+      fprintf (stderr, "ERROR: can't read header for mosaic %s\n", argv[N]);
+      exit (1);
+    }
+    if (!GetCoords (MOSAIC, &header)) {
+      fprintf (stderr, "ERROR: no astrometric solution in header\n");
+      exit (1);
+    }
+    if (strcmp(&MOSAIC[0].ctype[4], "-DIS")) {
+      fprintf (stderr, "ERROR: not a mosaic distortion header\n");
+      exit (1);
+    }
+    remove_argument (N, &argc, argv);
+    gfits_free_header (&header);
+  }
+
+  if (argc != 2) help();
+
+  return (TRUE);
+}
Index: /trunk/Ohana/src/getstar/src/dvoImageOverlaps.c
===================================================================
--- /trunk/Ohana/src/getstar/src/dvoImageOverlaps.c	(revision 12774)
+++ /trunk/Ohana/src/getstar/src/dvoImageOverlaps.c	(revision 12774)
@@ -0,0 +1,45 @@
+# include "dvoImageOverlaps.h"
+
+int main (int argc, char **argv) {
+
+  int i, Nimages, NdbImages, status;
+  int Nmatches, *matches;
+  Image *images, *dbImages;
+  Catalog catalog;
+  FITS_DB db;
+
+  SetSignals ();
+  ConfigInit_overlaps (&argc, argv);
+  args_overlaps (argc, argv);
+  
+  images = ReadImageFiles (argv[1], &Nimages); 
+
+  /*** update the image table ***/
+  /* setup image table format and lock */
+  db.mode   = dvo_catalog_catmode (CATMODE);
+  db.format = dvo_catalog_catformat (CATFORMAT);
+  status    = dvo_image_lock (&db, ImageCat, 3600.0, LCK_SOFT);  // shorter timeout?
+  if (!status) Shutdown ("ERROR: failure to lock image catalog %s", db.filename);
+
+  /* load or create the image table */
+  if (db.dbstate == LCK_EMPTY) {
+    fprintf (stderr, "no images in database (%s)\n", ImageCat);
+    exit (1);
+  } else {
+    if (!dvo_image_load (&db, VERBOSE, FALSE)) {
+      Shutdown ("can't read image catalog %s", db.filename);
+    }
+  }
+  dvo_image_unlock (&db);
+
+  // convert database table to internal structure
+  dbImages = gfits_table_get_Image (&db.ftable, &NdbImages, &db.swapped);
+  
+  // for (i = 1; i < 2; i++) {
+  for (i = 0; i < Nimages; i++) {
+    matches = MatchImage (dbImages, NdbImages, &images[i], &Nmatches);
+    ListImageOverlaps (dbImages, &images[i], matches, Nmatches);
+  }
+
+  exit (0);
+}
Index: /trunk/Ohana/src/getstar/src/edge_check.c
===================================================================
--- /trunk/Ohana/src/getstar/src/edge_check.c	(revision 12774)
+++ /trunk/Ohana/src/getstar/src/edge_check.c	(revision 12774)
@@ -0,0 +1,36 @@
+# include "dvoImageOverlaps.h"
+
+int edge_check (double *x1, double *y1, double *x2, double *y2) {
+
+  double theta1, theta2;
+  double Theta1, Theta2;
+
+  theta1 = opening_angle (x1[0], y1[0], x2[0], y2[0], x1[1], y1[1]); 
+  theta2 = opening_angle (x1[0], y1[0], x2[0], y2[0], x2[1], y2[1]); 
+
+  if (theta1*theta2 < 0.0) {
+    return (FALSE);
+  }
+
+  if (fabs(theta1) < fabs(theta2)) {
+    return (FALSE);
+  }
+
+  Theta1 = theta1;
+  Theta2 = theta2;
+  theta1 = opening_angle (x2[0], y2[0], x1[1], y1[1], x2[1], y2[1]); 
+  theta2 = opening_angle (x2[0], y2[0], x1[1], y1[1], x1[0], y1[0]); 
+  
+ 
+  if (theta1*theta2 < 0.0) {
+    return (FALSE);
+  }
+
+  if (fabs(theta1) < fabs(theta2)) {
+    return (FALSE);
+  }
+
+  return (TRUE);
+
+}
+
Index: /trunk/Ohana/src/getstar/src/opening_angle.c
===================================================================
--- /trunk/Ohana/src/getstar/src/opening_angle.c	(revision 12774)
+++ /trunk/Ohana/src/getstar/src/opening_angle.c	(revision 12774)
@@ -0,0 +1,20 @@
+# include "dvoImageOverlaps.h"
+
+double opening_angle (double x1, double y1, double x2, double y2, double x3, double y3) {
+
+  double dx1, dy1, dx2, dy2, ct, st, theta;
+
+  dx1 = x1 - x2;
+  dy1 = y1 - y2;
+  
+  dx2 = x3 - x2;
+  dy2 = y3 - y2;
+  
+  ct = (dx1*dx2 + dy1*dy2);
+  st = (dx1*dy2 - dx2*dy1);
+
+  theta = atan2 (st, ct);
+
+  return (theta);
+
+}
Index: /trunk/Ohana/src/getstar/src/parse_time.c
===================================================================
--- /trunk/Ohana/src/getstar/src/parse_time.c	(revision 12774)
+++ /trunk/Ohana/src/getstar/src/parse_time.c	(revision 12774)
@@ -0,0 +1,146 @@
+# include "addstar.h"
+
+int parse_time (Header *header) {
+
+  double jd;
+  int Ny, Nf, mode;
+  int Nsec, hour, min, sec, year, month, day;
+  char *py, *pm, *pd, *c;
+  char line[256];
+
+  /* we want to find JD or MJD to get Nsec (seconds since 01/01/1970) */
+
+  /* try JD first */
+  if (strcasecmp (JDKeyword, "NONE")) {
+    uppercase (JDKeyword);
+    if (!gfits_scan (header, JDKeyword, "%lf", 1, &jd)) {
+      fprintf (stderr, "ERROR: missing JD Keyword %s\n", JDKeyword);
+      exit (1);
+    }
+    Nsec = (jd - 2440587.5)*86400;
+    return (Nsec);
+  }
+
+  /* try MJD next */
+  if (strcasecmp (MJDKeyword, "NONE")) {
+    uppercase (MJDKeyword);
+    if (!gfits_scan (header, MJDKeyword, "%lf", 1, &jd)) {
+      fprintf (stderr, "ERROR: missing MJD Keyword %s\n", MJDKeyword);
+      exit (1);
+    }
+    Nsec = (jd - 40587.0)*86400;
+    return (Nsec);
+  }
+    
+  if (!strcasecmp (UTKeyword, "NONE")) {
+      fprintf (stderr, "ERROR: no valid Date/Time keywords\n");
+      exit (1);
+  }
+  if (!strcasecmp (DateKeyword, "NONE")) {
+      fprintf (stderr, "ERROR: no valid Date/Time keywords\n");
+      exit (1);
+  }
+  if (!strcasecmp (DateMode, "NONE")) {
+      fprintf (stderr, "ERROR: no valid Date/Time keywords\n");
+      exit (1);
+  }
+
+  /* get UT and DATE */
+  uppercase (UTKeyword);
+  if (!gfits_scan (header, UTKeyword, "%s", 1, line)) {
+      fprintf (stderr, "ERROR: missing UT Keyword %s\n", UTKeyword);
+      exit (1);
+    }
+  /* remove ':' characters */
+  for (c = strchr (line, 0x3a); c != (char *) NULL; c = strchr (line, 0x3a)) { *c = ' '; }
+  sscanf (line, "%d %d %d", &hour, &min, &sec);
+
+  /* parse mode line */
+  uppercase (DateMode);
+  for (Ny = 0, c = strchr (DateMode, 'Y'); c != (char ) NULL; c = strchr (c + 1, 'Y'), Ny++);
+  if ((Ny != 2) && (Ny != 4)) {
+    fprintf (stderr, "error in DATE-MODE format: %s\n", DateMode);
+    exit (1);
+  }
+  py = strchr (DateMode, 'Y');
+  pm = strchr (DateMode, 'M');
+  pd = strchr (DateMode, 'D');
+  if ((py == (char *) NULL) || (pm == (char *) NULL) || (pd == (char *) NULL)) {
+    fprintf (stderr, "error in DATE-MODE format: %s\n", DateMode);
+    exit (1);
+  }
+  if ((py > pm) && (py < pd)) {
+    fprintf (stderr, "error in DATE-MODE format: %s\n", DateMode);
+    exit (1);
+  }
+  if ((py > pd) && (py < pm)) {
+    fprintf (stderr, "error in DATE-MODE format: %s\n", DateMode);
+    exit (1);
+  }
+  mode = 0;
+  if ((py < pm) && (pm < pd)) { mode = 1; }  /* yyyy-mm-dd */
+  if ((py < pm) && (pm > pd)) { mode = 2; }  /* yyyy-dd-mm */
+  if ((py > pm) && (pm < pd)) { mode = 3; }  /* mm-dd-yyyy */
+  if ((py > pm) && (pm > pd)) { mode = 4; }  /* dd-mm-yyyy */
+  if (!mode) {
+    fprintf (stderr, "error in DATE-MODE format: %s\n", DateMode);
+    exit (1);
+  }
+
+  /* parse date entry */
+  uppercase (DateKeyword);
+  if (!gfits_scan (header, DateKeyword, "%s",  1, line)) {
+    fprintf (stderr, "ERROR: missing DATE Keyword %s\n", DateKeyword);
+    exit (1);
+  }
+  /* remove possible separators: ':', '/' '.', '-' */
+  for (c = strchr (line, 0x3a); c != (char *) NULL; c = strchr (line, 0x3a)) { *c = ' '; }
+  for (c = strchr (line, 0x2f); c != (char *) NULL; c = strchr (line, 0x2f)) { *c = ' '; }
+  for (c = strchr (line, 0x2e); c != (char *) NULL; c = strchr (line, 0x2e)) { *c = ' '; }
+  for (c = strchr (line, 0x2d); c != (char *) NULL; c = strchr (line, 0x2d)) { *c = ' '; }
+
+  Nf = 0;
+  switch (mode) {
+  case 1:
+    Nf = sscanf (line, "%d %d %d", &year, &month, &day);
+    break;
+  case 2:
+    Nf = sscanf (line, "%d %d %d", &year, &day, &month);
+    break;
+  case 3:
+    Nf = sscanf (line, "%d %d %d", &month, &day, &year);
+    break;
+  case 4:
+    Nf = sscanf (line, "%d %d %d", &day, &month, &year);
+    break;
+  }
+  if (Nf != 3) {
+    fprintf (stderr, "error in date entry (%s) or DATE-MODE format (%s)\n", line, DateMode);
+    exit (1);
+  }
+
+  if (year > 1000) {
+    if (Ny == 2) {
+      fprintf (stderr, "warning: mode line claims 2 digit year, but 4 digit year found\n");
+    }
+  } else {
+    if (Ny == 4) {
+      fprintf (stderr, "warning: mode line claims 4 digit year, but 2 digit year found\n");
+    }
+    if (year < 50) year += 100;
+    year += 1900;
+  }    
+
+  /* this should probably use localtime */
+
+  /* convert yy.mm.dd hh.mm.ss to Nsec since 1970 (jd = 2440587.5) */
+  /* note that in this section, tm_mon has range 1-12, unlike for gmtime () */
+  jd = day - 32075 + (int)(1461*(year + 4800 + (int)(((month)-14)/12))/4)
+    + (int)(367*((month) - 2 - (int)(((month) - 14)/12)*12)/12)
+    - (int)(3*(int)((year + 4900 + (int)(((month) - 14)/12))/100)/4) - 0.5;
+  /* jd is the julian day of the whole day only not the time */
+  Nsec = (jd - 2440587.5)*86400 + 3600.0*hour + min*60.0 + sec;
+  
+  return (Nsec);
+
+}