Index: trunk/Ohana/src/getstar/Makefile =================================================================== --- trunk/Ohana/src/getstar/Makefile (revision 24755) +++ trunk/Ohana/src/getstar/Makefile (revision 24914) @@ -23,9 +23,12 @@ dvoImageOverlaps.install: $(DESTBIN)/dvoImageOverlaps +dvoImagesAtCoords: $(BIN)/dvoImagesAtCoords.$(ARCH) +dvoImagesAtCoords.install: $(DESTBIN)/dvoImagesAtCoords + dvoImageExtract: $(BIN)/dvoImageExtract.$(ARCH) dvoImageExtract.install: $(DESTBIN)/dvoImageExtract -all: getstar dvoImageOverlaps dvoImageExtract -install: getstar.install dvoImageOverlaps.install dvoImageExtract.install +all: getstar dvoImageOverlaps dvoImageExtract dvoImagesAtCoords +install: getstar.install dvoImageOverlaps.install dvoImageExtract.install dvoImagesAtCoords.install GETSTAR = \ @@ -61,4 +64,18 @@ $(BIN)/dvoImageOverlaps.$(ARCH): $(OVERLAPS) +COVERLAPS = \ +$(SRC)/dvoImagesAtCoords.$(ARCH).o \ +$(SRC)/args_coords.$(ARCH).o \ +$(SRC)/ConfigInit_coords.$(ARCH).o \ +$(SRC)/Shutdown.$(ARCH).o \ +$(SRC)/SetSignals.$(ARCH).o \ +$(SRC)/ReadImageFiles.$(ARCH).o \ +$(SRC)/ReadImageHeader.$(ARCH).o \ +$(SRC)/GetFileMode.$(ARCH).o \ +$(SRC)/MatchCoords.$(ARCH).o + +$(COVERLAPS): $(INC)/dvoImagesAtCoords.h $(DVO_INCS) +$(BIN)/dvoImagesAtCoords.$(ARCH): $(COVERLAPS) + EXTRACT = \ $(SRC)/dvoImageExtract.$(ARCH).o \ Index: trunk/Ohana/src/getstar/include/dvoImagesAtCoords.h =================================================================== --- trunk/Ohana/src/getstar/include/dvoImagesAtCoords.h (revision 24914) +++ trunk/Ohana/src/getstar/include/dvoImagesAtCoords.h (revision 24914) @@ -0,0 +1,39 @@ +# include +# include +# include +# include + +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 ACCEPT_ASTROM; +char *astromFile; + +// These functions belong to dvoImageOverap +int args_coords PROTO((int argc, char **argv)); +int ConfigInit_coords PROTO((int *argc, char **argv)); +int Shutdown PROTO((char *format, ...)); +void TrapSignal PROTO((int sig)); +void SetProtect PROTO((int mode)); +int SetSignals PROTO((void)); + +Image *ReadImageFiles (char *filename, int *Nimages); +int ReadImageHeader (Header *header, Image *image); +int *MatchCoords(Image *, int, double, double, int *); + +int GetFileMode (Header *header); +// int edge_check (double *x1, double *y1, double *x2, double *y2); +// double opening_angle (double x1, double y1, double x2, double y2, double x3, double y3); Index: trunk/Ohana/src/getstar/src/ConfigInit_coords.c =================================================================== --- trunk/Ohana/src/getstar/src/ConfigInit_coords.c (revision 24914) +++ trunk/Ohana/src/getstar/src/ConfigInit_coords.c (revision 24914) @@ -0,0 +1,46 @@ +# include "dvoImagesAtCoords.h" + +int ConfigInit_coords (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/MatchCoords.c =================================================================== --- trunk/Ohana/src/getstar/src/MatchCoords.c (revision 24914) +++ trunk/Ohana/src/getstar/src/MatchCoords.c (revision 24914) @@ -0,0 +1,86 @@ +# include "dvoImagesAtCoords.h" +# ifndef FLT_MAX +# define FLT_MAX 1e32 +# endif + +/* given coordinate, find images in list that contain the point */ +int *MatchCoords (Image *dbImages, int NdbImages, double ra, double dec, int *Nmatch) { + + int i, j, N, addtolist, status; + int NMATCH, nmatch, *match; + Coords tcoords; + double r, d; + double Xi[4], Yi[4], Xo[4], Yo[4]; /* image and original corners */ + double Xmin, Xmax, Ymin, Ymax; + double xmin, xmax, ymin, ymax; + + *Nmatch = 0; + + /* 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]; + + + ymin = xmin = +FLT_MAX; + ymax = xmax = -FLT_MAX; + for (j = 0; j < 4; j++) { + xmin = MIN(xmin, Xi[j]); + xmax = MAX(xmax, Xi[j]); + ymin = MIN(ymin, Yi[j]); + ymax = MAX(ymax, Yi[j]); + } + + // transform input point to image coords + double x, y; + RD_to_XY(&x, &y, ra, dec, &dbImages[i].coords); + + if ((x >= xmin && x <= xmax) && (y >= ymin && y <= ymax)) { + + 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/args_coords.c =================================================================== --- trunk/Ohana/src/getstar/src/args_coords.c (revision 24914) +++ trunk/Ohana/src/getstar/src/args_coords.c (revision 24914) @@ -0,0 +1,86 @@ +# include "dvoImagesAtCoords.h" + +void help () { + fprintf (stderr, "USAGE: \n" + "dvoImagesAtCoords [-astrom astrom_file] (coords_file)\n" + ); + exit (2); +} + +int args_coords (int argc, char **argv) { + + int N; + + /* check for help request */ + if (get_argument (argc, argv, "-help") || + get_argument (argc, argv, "-h")) { + help (); + } + + astromFile = NULL; + if ((N = get_argument(argc, argv, "-astrom"))) { + remove_argument (N, &argc, argv); + if (argv[N] == NULL) { + fprintf (stderr, "ERROR: no file provided with -astrom\n"); + exit (1); + } + astromFile = strdup(argv[N]); + remove_argument (N, &argc, argv); + } + + /* 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); + } + /* accept bad header astrometry */ + ACCEPT_ASTROM = FALSE; + if ((N = get_argument (argc, argv, "-accept"))) { + ACCEPT_ASTROM = TRUE; + remove_argument (N, &argc, argv); + } + if ((N = get_argument (argc, argv, "-accept-astrom"))) { + ACCEPT_ASTROM = 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/dvoImagesAtCoords.c =================================================================== --- trunk/Ohana/src/getstar/src/dvoImagesAtCoords.c (revision 24914) +++ trunk/Ohana/src/getstar/src/dvoImagesAtCoords.c (revision 24914) @@ -0,0 +1,111 @@ +# include "dvoImagesAtCoords.h" + +typedef struct { + int id; + double ra; + double dec; +} Point; + +static int readPoints(char *filename, Point **pointsOut); +static int ListImagesAtCoords (Image *dbImages, int id, double ra, double dec, int *matches, int Nmatches); + +int main (int argc, char **argv) { + + int i, Nimages, NdbImages, status; + int Nmatches, *matches; + Image *dbImages; + int Npoints; + FITS_DB db; + + SetSignals (); + ConfigInit_coords (&argc, argv); + args_coords (argc, argv); + + Point *points; + Npoints = readPoints(argv[1], &points); + if (!Npoints) { + exit(0); + } + + if (astromFile) { + dbImages = ReadImageFiles(astromFile, &NdbImages); + } else { + /*** 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); + } + + Point *pt; + for (i = 0, pt = points; i < Npoints; i++, pt++) { + matches = MatchCoords (dbImages, NdbImages, pt->ra, pt->dec, &Nmatches); + ListImagesAtCoords(dbImages, pt->id, pt->ra, pt->dec, matches, Nmatches); + } + + exit (0); +} + +static int readPoints(char *filename, Point **pointsOut) +{ + FILE *f = fopen(filename, "r"); + if (!f) { + fprintf (stderr, "failed to open coordinate file %s", filename); + exit(1); + } + + int LEN = 20; + Point *pts; + ALLOCATE(pts, Point, LEN); + int Npoints = 0; + + int Nread; + // File Format is simple text file. + // Each Line has 3 values separated by spaces that represent the coordinates to match. + // + // 'unique integer id' 'ra in degrees' 'declination in degrees' + // we don't check the id for uniqueness (garbage in garbage out). + + while ((Nread = fscanf(f, "%d %lf %lf\n", &pts[Npoints].id, &pts[Npoints].ra, &pts[Npoints].dec)) == 3) { + Npoints++; + if (Npoints >= LEN) { + LEN += 20; + REALLOCATE(pts, Point, LEN); + } + } + if (Nread != -1) { + fprintf (stderr, "unexpected data on line %d of points file %s", Npoints, filename); + exit(1); + } + + *pointsOut = pts; + + return Npoints; +} + +static int ListImagesAtCoords (Image *dbImages, int id, double ra, double dec, int *matches, int Nmatches) +{ + int i; + for (i = 0; i < Nmatches; i++) { + int N = matches[i]; + fprintf (stdout, "%d (%lf, %lf) : %s\n", id, ra, dec, dbImages[N].name); + } + + return (TRUE); +} +