Index: /trunk/Ohana/src/uniphot/Makefile =================================================================== --- /trunk/Ohana/src/uniphot/Makefile (revision 35102) +++ /trunk/Ohana/src/uniphot/Makefile (revision 35103) @@ -23,6 +23,8 @@ setposangle: $(BIN)/setposangle.$(ARCH) setposangle_client: $(BIN)/setposangle_client.$(ARCH) +setastrom: $(BIN)/setastrom.$(ARCH) +setastrom_client: $(BIN)/setastrom_client.$(ARCH) -install: $(DESTBIN)/uniphot $(DESTBIN)/setfwhm $(DESTBIN)/setphot $(DESTBIN)/setphot_client $(DESTBIN)/setposangle $(DESTBIN)/setposangle_client +install: $(DESTBIN)/uniphot $(DESTBIN)/setfwhm $(DESTBIN)/setphot $(DESTBIN)/setphot_client $(DESTBIN)/setposangle $(DESTBIN)/setposangle_client $(DESTBIN)/setastrom $(DESTBIN)/setastrom_client UNIPHOT = \ @@ -99,4 +101,26 @@ $(BIN)/setposangle_client.$(ARCH): $(SETPOSANGLE_CLIENT) +SETASTROM = \ +$(SRC)/setastrom.$(ARCH).o \ +$(SRC)/initialize_setastrom.$(ARCH).o \ +$(SRC)/astrom_correction.$(ARCH).o \ +$(SRC)/update_dvo_setastrom.$(ARCH).o \ +$(SRC)/update_catalog_setastrom.$(ARCH).o \ +$(SRC)/SetSignals.$(ARCH).o \ +$(SRC)/Shutdown.$(ARCH).o + +$(SETASTROM): $(INC)/setastrom.h +$(BIN)/setastrom.$(ARCH): $(SETASTROM) + +SETASTROM_CLIENT = \ +$(SRC)/setastrom_client.$(ARCH).o \ +$(SRC)/update_dvo_setastrom.$(ARCH).o \ +$(SRC)/update_catalog_setastrom.$(ARCH).o \ +$(SRC)/astrom_correction.$(ARCH).o \ +$(SRC)/initialize_setastrom.$(ARCH).o + +$(SETASTROM_CLIENT): $(INC)/setastrom.h +$(BIN)/setastrom_client.$(ARCH): $(SETASTROM_CLIENT) + SETFWHM = \ $(SRC)/setfwhm.$(ARCH).o \ Index: /trunk/Ohana/src/uniphot/include/setastrom.h =================================================================== --- /trunk/Ohana/src/uniphot/include/setastrom.h (revision 35103) +++ /trunk/Ohana/src/uniphot/include/setastrom.h (revision 35103) @@ -0,0 +1,56 @@ +# include +# include +# include + +// the interpolating spline information: this is based on the spline used by opihi/cmd.data/spline.c +typedef struct { + int Nknots; + double *xk; + double *yk; + double *y2; +} Spline; + +/* global variables set in parameter file */ +# define DVO_MAX_PATH 1024 +char ImageCat[DVO_MAX_PATH]; +char *CATDIR; +int HOST_ID; +char *HOSTDIR; +int VERBOSE; +int UPDATE; +int PARALLEL; +int PARALLEL_MANUAL; +int PARALLEL_SERIAL; +char *UPDATE_CATFORMAT; +char *SINGLE_CPT; + +SkyRegion UserPatch; + +/***** prototypes ****/ +int main PROTO((int argc, char **argv)); + +void GetConfig PROTO((char *config, char *field, char *format, int N, void *ptr)); +void ConfigInit PROTO((int *argc, char **argv)); + +void usage_setastrom PROTO(()); +void initialize_setastrom PROTO((int argc, char **argv)); +int args_setastrom PROTO((int argc, char **argv)); + +void usage_setastrom_client PROTO(()); +void initialize_setastrom_client PROTO((int argc, char **argv)); +int args_setastrom_client PROTO((int argc, char **argv)); + +// int Shutdown PROTO((char *format, ...) OHANA_FORMAT(printf, 1, 2) ); +// void TrapSignal PROTO((int sig)); +// void SetProtect PROTO((int mode)); +// int SetSignals PROTO((void)); + +int update_dvo_setastrom PROTO((char *corrfile)); +int update_dvo_setastrom_parallel PROTO((SkyTable *sky, char *corrfile)); + +int update_catalog_setastrom PROTO((Catalog *catalog)); + +// spline correction functions +int load_astrom_correction (char *filename); +int get_astrom_correction (int sub, int chip, double *dX, double *dY, float Minst); +double spline_apply_dbl (double *x, double *y, double *y2, int N, double X); Index: /trunk/Ohana/src/uniphot/include/setphot.h =================================================================== --- /trunk/Ohana/src/uniphot/include/setphot.h (revision 35102) +++ /trunk/Ohana/src/uniphot/include/setphot.h (revision 35103) @@ -25,4 +25,5 @@ char *HOSTDIR; char *IMAGES; +char *SINGLE_CPT; int VERBOSE; int RESET; @@ -34,4 +35,6 @@ int PARALLEL_SERIAL; int IMAGES_ONLY; + +SkyRegion UserPatch; /***** prototypes ****/ Index: /trunk/Ohana/src/uniphot/include/setphot_client.h =================================================================== --- /trunk/Ohana/src/uniphot/include/setphot_client.h (revision 35102) +++ /trunk/Ohana/src/uniphot/include/setphot_client.h (revision 35103) @@ -78,4 +78,6 @@ time_t TSTOP; PhotCode *photcode; + +SkyRegion UserPatch; enum {black, white, red, orange, yellow, green, blue, indigo, violet}; Index: /trunk/Ohana/src/uniphot/include/setposangle.h =================================================================== --- /trunk/Ohana/src/uniphot/include/setposangle.h (revision 35102) +++ /trunk/Ohana/src/uniphot/include/setposangle.h (revision 35103) @@ -25,4 +25,5 @@ char *HOSTDIR; char *IMAGES; +char *SINGLE_CPT; int VERBOSE; int RESET; @@ -31,4 +32,6 @@ int PARALLEL_MANUAL; int PARALLEL_SERIAL; + +SkyRegion UserPatch; /***** prototypes ****/ Index: /trunk/Ohana/src/uniphot/src/astrom_correction.c =================================================================== --- /trunk/Ohana/src/uniphot/src/astrom_correction.c (revision 35103) +++ /trunk/Ohana/src/uniphot/src/astrom_correction.c (revision 35103) @@ -0,0 +1,182 @@ +# include "setastrom.h" + +// astrometry correction functions: +// load the correction set from a file + +// the correction is saved as a set of splines, one per file extension +// extensions have names of the form XYnn.Dx.Td + +// where nn = 00 - 77, the chip number +// x is X or Y, the correction direction (dir) +// d is 0 or 1, pre or post camera correction (sub) + +// these are stored in an N-D array of (Spline *) splineset[sub][dir][chip] + +# define N_SUB 2 +# define N_DIR 2 +# define N_CHIP 78 +// note that we are allocating pointers for XY00 - XY77, but only the octal elements are set (and not the 00,07,70,77 ones) + +Spline ****splineset; + +int load_astrom_correction (char *filename) { + + int Ncol; + off_t Nrow; + + char type[16]; + char extname[80]; + + Header header; + Header theader; + FTable ftable; + + /* open file for input */ + FILE *f = fopen (filename, "r"); + if (f == (FILE *) NULL) { + gprint (GP_ERR, "can't open file for read : %s\n", filename); + return FALSE; + } + + // read the PHU header + if (!gfits_load_header (f, &header)) return (FALSE); + + int Nbytes = gfits_data_size (&header); + fseeko (f, Nbytes, SEEK_CUR); + // XXX gfits_free_header (&header); + + ftable.header = &theader; + + int sub, dir, chip; + + // we have a N-D array of splines: + // splineset[sub][dir][chip] + + ALLOCATE (splineset, Spline ***, N_SUB); + for (sub = 0; sub < N_SUB; sub++) { + ALLOCATE (splineset[sub], Spline **, N_DIR); + for (dir = 0; dir < N_DIR; dir++) { + ALLOCATE (splineset[sub][dir], Spline *, N_CHIP); + for (chip = 0; chip < N_CHIP; chip ++) { + splineset[sub][dir][chip] = NULL; + } + } + } + + while (TRUE) { + // load data for this header : if not found, assume we hit the end of the file + if (!gfits_load_header (f, &theader)) break; + + if (!gfits_fread_ftable_data (f, &ftable, FALSE)) return (FALSE); + + if (!gfits_scan (&theader, "EXTNAME", "%s", 1, extname)) return (FALSE); + + // figure out the chip ID, diretion, and sequence from the EXTNAME + + // EXTNAME is of the form XY%d%d.D%s.T%s where + char chipname[3]; + chipname[0] = extname[2]; + chipname[1] = extname[3]; + chipname[2] = 0; + int chip = atoi(chipname); + myAssert (chip > 0, "error in extname %s not(chip > 0)\n", extname); + myAssert (chip < 77, "error in extname %s not(chip < 77)\n", extname); + myAssert (chip != 7, "error in extname %s not(chip != 7)\n", extname); + myAssert (chip != 70, "error in extname %s not(chip != 70)\n", extname); + + int dir = -1; + if (extname[6] == 'X') { + dir = 0; + } + if (extname[6] == 'Y') { + dir = 1; + } + myAssert (dir != -1, "error in extname %s (dir)\n", extname); + + int sub = -1; + if (extname[9] == '0') { + sub = 0; + } + if (extname[9] == '1') { + sub = 1; + } + myAssert (sub != -1, "error in extname %s (sub)\n", extname); + + // need to create and assign to flat-field correction + double *xk = gfits_get_bintable_column_data (&theader, &ftable, "X_KNOT", type, &Nrow, &Ncol); + myAssert (!strcmp(type, "double"), "error reading X_KNOT"); + + // need to create and assign to flat-field correction + double *yk = gfits_get_bintable_column_data (&theader, &ftable, "Y_KNOT", type, &Nrow, &Ncol); + myAssert (!strcmp(type, "double"), "error reading X_KNOT"); + + // need to create and assign to flat-field correction + double *y2 = gfits_get_bintable_column_data (&theader, &ftable, "DY2_DX", type, &Nrow, &Ncol); + myAssert (!strcmp(type, "double"), "error reading X_KNOT"); + + ALLOCATE (splineset[sub][dir][chip], Spline, 1); + splineset[sub][dir][chip]->xk = xk; + splineset[sub][dir][chip]->yk = yk; + splineset[sub][dir][chip]->y2 = y2; + splineset[sub][dir][chip]->Nknots = Nrow; + } + return (TRUE); +} + +int get_astrom_correction (int sub, int chip, double *dX, double *dY, float Minst) { + + Spline *spline = NULL; + + *dX = 0.0; + *dY = 0.0; + + // XXX need to include some limits on Minst? + + spline = splineset[sub][0][chip]; + if (spline) { + *dX = spline_apply_dbl (spline->xk, spline->yk, spline->y2, spline->Nknots, Minst); + } + + spline = splineset[sub][1][chip]; + if (spline) { + *dY = spline_apply_dbl (spline->xk, spline->yk, spline->y2, spline->Nknots, Minst); + } + return (TRUE); +} + +/* evaluate spline for x, y, y2 at X */ +double spline_apply_dbl (double *x, double *y, double *y2, int N, double X) { + + int i, lo, hi; + double dx, a, b, value; + + // saturate correction at high and low ends + if (X < x[0]) return y[0]; + if (X > x[N-1]) return y[N-1]; + + /* find correct element in array (x must be sorted) */ + lo = 0; + hi = N-1; + while (hi - lo > 1) { + i = 0.5*(hi+lo); + if (x[i] > X) { + hi = i; + } else { + lo = i; + } + } + + /* error condition: duplicate abssisca */ + dx = x[hi] - x[lo]; + if (dx == 0.0) { + return (HUGE_VAL); + } + + /* evaluate spline */ + a = (x[hi] - X) / dx; + b = (X - x[lo]) / dx; + + value = a*y[lo] + b*y[hi] + ((a*a*a - a)*y2[lo] + (b*b*b - b)*y2[hi])*(dx*dx) / 6.0; + return (value); + +} Index: /trunk/Ohana/src/uniphot/src/initialize_setastrom.c =================================================================== --- /trunk/Ohana/src/uniphot/src/initialize_setastrom.c (revision 35103) +++ /trunk/Ohana/src/uniphot/src/initialize_setastrom.c (revision 35103) @@ -0,0 +1,270 @@ +# include "setastrom.h" + +void usage_setastrom () { + fprintf (stderr, "USAGE: setastrom (corrections) [options]\n"); + fprintf (stderr, " options:\n"); + fprintf (stderr, " -v : verbose mode\n"); + fprintf (stderr, " -region Rmin Rmax Dmin Dmax\n"); + fprintf (stderr, " -update-catformat (format) : change database schema on output\n"); + fprintf (stderr, " -update : actually write results to detections tables\n"); + fprintf (stderr, " -parallel : run in parallel mode\n"); + fprintf (stderr, " -h : this help list\n"); + fprintf (stderr, " -help : this help list\n"); + fprintf (stderr, " --h : this help list\n"); + fprintf (stderr, " --help : this help list\n"); + fprintf (stderr, " Note that the dvo db can be specified by -D CATDIR (directory)\n"); + exit (2); +} + +void initialize_setastrom (int argc, char **argv) { + + int N; + struct stat statbuffer; + char CatdirPhotcodeFile[256]; + + if ((N = get_argument (argc, argv, "-h"))) usage_setastrom(); + if ((N = get_argument (argc, argv, "-help"))) usage_setastrom(); + if ((N = get_argument (argc, argv, "--h"))) usage_setastrom(); + if ((N = get_argument (argc, argv, "--help"))) usage_setastrom(); + + /* are these set correctly? */ + ConfigInit (&argc, argv); + args_setastrom (argc, argv); + + if (stat (CATDIR, &statbuffer)) { + fprintf (stderr, "error accessing dvo database directory '%s'\n", CATDIR); + exit (1); + } + + // load the photcode table + sprintf (CatdirPhotcodeFile, "%s/Photcodes.dat", CATDIR); + if (!LoadPhotcodes (CatdirPhotcodeFile, NULL, FALSE)) { + fprintf (stderr, "error loading photcode table %s\n", CatdirPhotcodeFile); + exit (1); + } +} + +int args_setastrom (int argc, char **argv) { + + int N; + + SINGLE_CPT = NULL; + if ((N = get_argument (argc, argv, "-cpt"))) { + remove_argument (N, &argc, argv); + SINGLE_CPT = strcreate (argv[N]); + remove_argument (N, &argc, argv); + } + + /* specify portion of the sky : allow default of all sky? */ + UserPatch.Rmin = 0; + UserPatch.Rmax = 360; + UserPatch.Dmin = -90; + UserPatch.Dmax = +90; + if ((N = get_argument (argc, argv, "-region"))) { + remove_argument (N, &argc, argv); + UserPatch.Rmin = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Rmax = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Dmin = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Dmax = atof (argv[N]); + remove_argument (N, &argc, argv); + } + + UPDATE_CATFORMAT = NULL; + if ((N = get_argument (argc, argv, "-update-catformat"))) { + remove_argument (N, &argc, argv); + UPDATE_CATFORMAT = strcreate (argv[N]); + remove_argument (N, &argc, argv); + } + + VERBOSE = FALSE; + if ((N = get_argument (argc, argv, "-v"))) { + VERBOSE = TRUE; + remove_argument (N, &argc, argv); + } + + PARALLEL = FALSE; + if ((N = get_argument (argc, argv, "-parallel"))) { + PARALLEL = TRUE; + remove_argument (N, &argc, argv); + } + // this is a test mode : rather than launching the remote jobs and waiting for completion, + // setastrom will simply list the remote command and wait for the user to signal completion + PARALLEL_MANUAL = FALSE; + if ((N = get_argument (argc, argv, "-parallel-manual"))) { + PARALLEL = TRUE; // -parallel-manual implies -parallel + PARALLEL_MANUAL = TRUE; + remove_argument (N, &argc, argv); + } + // this is a test mode : rather than launching the setastrom_client jobs remotely, they are + // run in serial via 'system' + PARALLEL_SERIAL = FALSE; + if ((N = get_argument (argc, argv, "-parallel-serial"))) { + if (PARALLEL_MANUAL) { + fprintf (stderr, "ERROR: cannot mix -parallel-manual and -parallel-serial\n"); + exit (1); + } + PARALLEL = TRUE; // -parallel-serial implies -parallel + PARALLEL_SERIAL = TRUE; + remove_argument (N, &argc, argv); + } + + UPDATE = FALSE; + if ((N = get_argument (argc, argv, "-update"))) { + remove_argument (N, &argc, argv); + UPDATE = TRUE; + } + + if (argc != 2) usage_setastrom(); + + return (TRUE); +} + +void usage_setastrom_client () { + fprintf (stderr, "USAGE: setastrom_client (correction) -hostID (hostID) -catdir (catdir) -hostdir (hostdir) [options]\n"); + fprintf (stderr, " options:\n"); + fprintf (stderr, " -region Rmin Rmax Dmin Dmax\n"); + fprintf (stderr, " -update-catformat (format) : change database schema on output\n"); + fprintf (stderr, " -v : verbose mode\n"); + fprintf (stderr, " -update : actually write results to detections tables\n"); + exit (2); +} + +void initialize_setastrom_client (int argc, char **argv) { + + int N; + struct stat statbuffer; + char CatdirPhotcodeFile[256]; + + if ((N = get_argument (argc, argv, "-h"))) usage_setastrom_client(); + if ((N = get_argument (argc, argv, "-help"))) usage_setastrom_client(); + if ((N = get_argument (argc, argv, "--h"))) usage_setastrom_client(); + if ((N = get_argument (argc, argv, "--help"))) usage_setastrom_client(); + + SetZeroPoint (25.0); // XXX is this needed? + args_setastrom_client (argc, argv); + + if (stat (CATDIR, &statbuffer)) { + fprintf (stderr, "error accessing dvo database directory '%s'\n", CATDIR); + exit (1); + } + + // load the photcode table : XXX needed? + sprintf (CatdirPhotcodeFile, "%s/Photcodes.dat", CATDIR); + if (!LoadPhotcodes (CatdirPhotcodeFile, NULL, FALSE)) { + fprintf (stderr, "error loading photcode table %s\n", CatdirPhotcodeFile); + exit (1); + } +} + +int args_setastrom_client (int argc, char **argv) { + + int N; + + /* specify portion of the sky : allow default of all sky? */ + UserPatch.Rmin = 0; + UserPatch.Rmax = 360; + UserPatch.Dmin = -90; + UserPatch.Dmax = +90; + if ((N = get_argument (argc, argv, "-region"))) { + remove_argument (N, &argc, argv); + UserPatch.Rmin = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Rmax = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Dmin = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Dmax = atof (argv[N]); + remove_argument (N, &argc, argv); + } + + UPDATE_CATFORMAT = NULL; + if ((N = get_argument (argc, argv, "-update-catformat"))) { + remove_argument (N, &argc, argv); + UPDATE_CATFORMAT = strcreate (argv[N]); + remove_argument (N, &argc, argv); + } + + VERBOSE = FALSE; + if ((N = get_argument (argc, argv, "-v"))) { + VERBOSE = TRUE; + remove_argument (N, &argc, argv); + } + + UPDATE = FALSE; + if ((N = get_argument (argc, argv, "-update"))) { + remove_argument (N, &argc, argv); + UPDATE = TRUE; + } + + HOST_ID = 0; + if ((N = get_argument (argc, argv, "-hostID"))) { + remove_argument (N, &argc, argv); + HOST_ID = atoi (argv[N]); + remove_argument (N, &argc, argv); + } + if (!HOST_ID) usage_setastrom_client(); + + HOSTDIR = NULL; + if ((N = get_argument (argc, argv, "-hostdir"))) { + remove_argument (N, &argc, argv); + HOSTDIR = strcreate (argv[N]); + remove_argument (N, &argc, argv); + } + if (!HOSTDIR) usage_setastrom_client(); + + if ((N = get_argument (argc, argv, "-catdir"))) { + remove_argument (N, &argc, argv); + CATDIR = strcreate (argv[N]); + remove_argument (N, &argc, argv); + } + if (!CATDIR) usage_setastrom_client(); + + if (argc != 2) usage_setastrom_client(); + return (TRUE); +} + +void ConfigInit (int *argc, char **argv) { + + double ZERO_POINT; + char *config, *file; + + /*** 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 (0); + } + if (VERBOSE) fprintf (stderr, "loaded config file: %s\n", file); + + // force CATDIR to be absolute (so parallel mode will work) + char *tmpcatdir = NULL; + ALLOCATE (tmpcatdir, char, DVO_MAX_PATH); + GetConfig (config, "CATDIR", "%s", 0, tmpcatdir); + CATDIR = abspath (tmpcatdir, DVO_MAX_PATH); + free (tmpcatdir); + + sprintf (ImageCat, "%s/Images.dat", CATDIR); + + ScanConfig (config, "ZERO_PT", "%lf", 0, &ZERO_POINT); + SetZeroPoint (ZERO_POINT); + + free (config); + free (file); +} + +void GetConfig (char *config, char *field, char *format, int N, void *ptr) { + + char *status; + + status = ScanConfig (config, field, format, N, ptr); + if (status == NULL) { + fprintf (stderr, "error in config, cannot find %s\n", field); + exit (1); + } + return; +} Index: /trunk/Ohana/src/uniphot/src/initialize_setphot.c =================================================================== --- /trunk/Ohana/src/uniphot/src/initialize_setphot.c (revision 35102) +++ /trunk/Ohana/src/uniphot/src/initialize_setphot.c (revision 35103) @@ -49,4 +49,28 @@ if ((N = get_argument (argc, argv, "-reset"))) { RESET = TRUE; + remove_argument (N, &argc, argv); + } + + SINGLE_CPT = NULL; + if ((N = get_argument (argc, argv, "-cpt"))) { + remove_argument (N, &argc, argv); + SINGLE_CPT = strcreate (argv[N]); + remove_argument (N, &argc, argv); + } + + // region of interest + UserPatch.Rmin = 0; + UserPatch.Rmax = 360; + UserPatch.Dmin = -90; + UserPatch.Dmax = +90; + if ((N = get_argument (argc, argv, "-region"))) { + remove_argument (N, &argc, argv); + UserPatch.Rmin = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Rmax = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Dmin = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Dmax = atof (argv[N]); remove_argument (N, &argc, argv); } Index: /trunk/Ohana/src/uniphot/src/initialize_setphot_client.c =================================================================== --- /trunk/Ohana/src/uniphot/src/initialize_setphot_client.c (revision 35102) +++ /trunk/Ohana/src/uniphot/src/initialize_setphot_client.c (revision 35103) @@ -36,4 +36,21 @@ if ((N = get_argument (argc, argv, "-reset"))) { RESET = TRUE; + remove_argument (N, &argc, argv); + } + + // region of interest + UserPatch.Rmin = 0; + UserPatch.Rmax = 360; + UserPatch.Dmin = -90; + UserPatch.Dmax = +90; + if ((N = get_argument (argc, argv, "-region"))) { + remove_argument (N, &argc, argv); + UserPatch.Rmin = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Rmax = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Dmin = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Dmax = atof (argv[N]); remove_argument (N, &argc, argv); } Index: /trunk/Ohana/src/uniphot/src/initialize_setposangle.c =================================================================== --- /trunk/Ohana/src/uniphot/src/initialize_setposangle.c (revision 35102) +++ /trunk/Ohana/src/uniphot/src/initialize_setposangle.c (revision 35103) @@ -58,4 +58,28 @@ // XXX for now, hardwire gpc1 photcodes + SINGLE_CPT = NULL; + if ((N = get_argument (argc, argv, "-cpt"))) { + remove_argument (N, &argc, argv); + SINGLE_CPT = strcreate (argv[N]); + remove_argument (N, &argc, argv); + } + + // region of interest + UserPatch.Rmin = 0; + UserPatch.Rmax = 360; + UserPatch.Dmin = -90; + UserPatch.Dmax = +90; + if ((N = get_argument (argc, argv, "-region"))) { + remove_argument (N, &argc, argv); + UserPatch.Rmin = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Rmax = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Dmin = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Dmax = atof (argv[N]); + remove_argument (N, &argc, argv); + } + PARALLEL = FALSE; if ((N = get_argument (argc, argv, "-parallel"))) { @@ -146,4 +170,21 @@ } + // region of interest + UserPatch.Rmin = 0; + UserPatch.Rmax = 360; + UserPatch.Dmin = -90; + UserPatch.Dmax = +90; + if ((N = get_argument (argc, argv, "-region"))) { + remove_argument (N, &argc, argv); + UserPatch.Rmin = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Rmax = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Dmin = atof (argv[N]); + remove_argument (N, &argc, argv); + UserPatch.Dmax = atof (argv[N]); + remove_argument (N, &argc, argv); + } + HOST_ID = 0; if ((N = get_argument (argc, argv, "-hostID"))) { Index: /trunk/Ohana/src/uniphot/src/setastrom.c =================================================================== --- /trunk/Ohana/src/uniphot/src/setastrom.c (revision 35103) +++ /trunk/Ohana/src/uniphot/src/setastrom.c (revision 35103) @@ -0,0 +1,23 @@ +# include "setastrom.h" + +int main (int argc, char **argv) { + + int status; + + /* get configuration info, args, lockfile */ + initialize_setastrom (argc, argv); + + // argv[1] holds the correction file + status = update_dvo_setastrom (argv[1]); + + if (!status) exit (1); + exit (0); +} + + +/* setastrom : set the posangle & pltscale for measurements in the db + + ** load images + ** load catalogs + ** update detection (& averages?) + */ Index: /trunk/Ohana/src/uniphot/src/setastrom_client.c =================================================================== --- /trunk/Ohana/src/uniphot/src/setastrom_client.c (revision 35103) +++ /trunk/Ohana/src/uniphot/src/setastrom_client.c (revision 35103) @@ -0,0 +1,20 @@ +# include "setastrom.h" + +// setastrom_client is run on a remote host and is responsible for updating the catalogs +// owned by that host. + +// setastrom_client -hostname (hostname) -catdir (catdir) -hostdir (hostdir) -images (images) +// load the SkyTable and HostTable from (catdir) [contains the host responsibilities] +// load the image data from (images) [contains the zp infomation to apply] +// loop over tables from SkyTable with my host ID +// the calling program tells the client which host they are (or host ID?); we do not actually have to run this on the real host + +int main (int argc, char **argv) { + + // get configuration info, args, lockfile (set CATDIR, HOST_ID, HOSTDIR, IMAGES) + initialize_setastrom_client (argc, argv); + + update_dvo_setastrom (argv[1]); + + exit (0); +} Index: /trunk/Ohana/src/uniphot/src/update_catalog_setastrom.c =================================================================== --- /trunk/Ohana/src/uniphot/src/update_catalog_setastrom.c (revision 35103) +++ /trunk/Ohana/src/uniphot/src/update_catalog_setastrom.c (revision 35103) @@ -0,0 +1,116 @@ +# include "setastrom.h" + +// XXX this function is dependent in two ways on the specific numerical values of the photcodes +// 1) the examined measurements are limited to GPC1 values based on the photcode range 10000 - 10600 +// 2) the chip ID is determined from the numerical value of the photcode (chip = photcode % 100) + +int update_catalog_setastrom (Catalog *catalog) { + + off_t i, found; + + time_t timeRef = ohana_date_to_sec("2011/05/11,00:00:00"); + + found = 0; + + for (i = 0; i < catalog[0].Nmeasure; i++) { + + Measure *measure = &catalog[0].measure[i]; + + // XXX only do GPC1 data for now + if (measure[0].photcode < 10000) continue; + if (measure[0].photcode > 10600) continue; + + // XXX hardwire the list of chips to correct + // some notes: + // XY24 surpisingly does not need correction + // XY27 has poor astrometry with a weird correction for all mags + // XY36 surpisingly does not need correction + // XY67 has poor astrometry with a weird correction for all mags + + int chipID = measure[0].photcode % 100; + switch (chipID) { + case 4: + case 5: + case 6: + case 14: + case 15: + case 16: + case 17: + case 25: + case 26: + case 34: + case 35: + case 37: + case 40: + case 41: + case 42: + case 43: + case 50: + case 51: + case 52: + case 53: + case 60: + case 61: + case 62: + case 63: + case 71: + case 72: + case 73: + break; + default: + measure[0].Xfix = measure[0].Xccd; + measure[0].Yfix = measure[0].Yccd; + continue; + } + + double Xccd = measure[0].Xccd; + double Yccd = measure[0].Yccd; + double pltScale = fabs(measure[0].pltscale); + + // this is not needed, is it? + // double posAngle = FromShortDegrees(measure[0].posangle); + + // I need to know which chip we have. use photcode assuming range + int chip = measure[0].photcode % 100; + + // correction may be Minst or Minst + 5.0*log(fwhm + float fwhm_maj = FromShortPixels(measure[0].FWx); + float Minst = PhotInst (measure); + float MinstSB = Minst + 5.0*log10(fwhm_maj); + + // I should validate that the spline I'm using is SB corrected not Minst... + + // XXX this is very ad hoc. the splines are measured based on an extraction from dvo + // the dvo values of fwhm_maj are sigmas (NOT FWHM!) in pixels. + + // saturate fwhm_maj in the range 1.0 < fwhm_maj < 5.0 + // (images with really odd psf probably have poor astrometry anyway...) + fwhm_maj = MIN(fwhm_maj,5.0); + fwhm_maj = MAX(fwhm_maj,1.0); + + double dX; + double dY; + + // correction is in arcseconds + if (measure[0].t < timeRef) { + get_astrom_correction (0, chip, &dX, &dY, MinstSB); + } else { + // XXX do not correct new data + // get_astrom_correction (1, chip, &dX, &dY, MinstSB); + dX = 0.0; + dY = 0.0; + } + + // do not modify the original Xccd,Yccd values + measure[0].Xfix = Xccd + dX / pltScale; + measure[0].Yfix = Yccd + dY / pltScale; + + found ++; + } + + if (found) { + fprintf (stderr, "found "OFF_T_FMT" matches\n", found); + } + + return (TRUE); +} Index: /trunk/Ohana/src/uniphot/src/update_dvo_setastrom.c =================================================================== --- /trunk/Ohana/src/uniphot/src/update_dvo_setastrom.c (revision 35103) +++ /trunk/Ohana/src/uniphot/src/update_dvo_setastrom.c (revision 35103) @@ -0,0 +1,141 @@ +# include "setastrom.h" + +int update_dvo_setastrom (char *corrfile) { + + SkyTable *sky = NULL; + SkyList *skylist = NULL; + Catalog catalog; + off_t i; + + // load the current sky table (layout of all SkyRegions) + sky = SkyTableLoadOptimal (CATDIR, NULL, NULL, FALSE, -1, VERBOSE); + SkyTableSetFilenames (sky, CATDIR, "cpt"); + + if (PARALLEL && !HOST_ID) { + int status = update_dvo_setastrom_parallel (sky, corrfile); + return status; + } + + load_astrom_correction (corrfile); + + // determine the populated SkyRegions overlapping the requested area (default depth) + if (SINGLE_CPT) { + skylist = SkyRegionByCPT (sky, SINGLE_CPT); + } else { + skylist = SkyListByPatch (sky, -1, &UserPatch); + } + myAssert (skylist, "ooops!"); + + // update measurements for each populated catalog + for (i = 0; i < skylist[0].Nregions; i++) { + + // does this host ID match the desired location for the table? + if (!HostTableTestHost(skylist[0].regions[i], HOST_ID)) continue; + + char localFilename[1024]; + snprintf (localFilename, 1024, "%s/%s.cpt", HOSTDIR, skylist->regions[i]->name); + + // set up the basic catalog info + catalog.filename = HOST_ID ? localFilename : skylist[0].filename[i]; + catalog.catformat = DVO_FORMAT_UNDEF; // not needed since we skip empty catalogs + catalog.catmode = DVO_MODE_UNDEF; // not needed since we skip empty catalogs + catalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF; + catalog.Nsecfilt = GetPhotcodeNsecfilt (); + + if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE, "w")) { + fprintf (stderr, "ERROR: failure reading catalog %s\n", catalog.filename); + exit (1); + } + if (!catalog.Naves_disk) { + if (VERBOSE) fprintf (stderr, "no data in %s, skipping\n", catalog.filename); + dvo_catalog_unlock (&catalog); + dvo_catalog_free (&catalog); + continue; + } + + update_catalog_setastrom (&catalog); + + if (!UPDATE) { + dvo_catalog_unlock (&catalog); + dvo_catalog_free (&catalog); + continue; + } + + // modify the output format as desired (ignore current format on disk) + if (UPDATE_CATFORMAT) { + catalog.catformat = dvo_catalog_catformat (UPDATE_CATFORMAT); + } + + if (VERBOSE) fprintf (stderr, "saving catalog %s\n", catalog.filename); + + dvo_catalog_save (&catalog, VERBOSE); + dvo_catalog_unlock (&catalog); + dvo_catalog_free (&catalog); + } + return (TRUE); +} + +# define DEBUG 1 + +int update_dvo_setastrom_parallel (SkyTable *sky, char *corrfile) { + + // now launch the setastrom_client jobs to the parallel hosts + + char *abscorrfile = abspath (corrfile, DVO_MAX_PATH); + + // load the list of hosts + HostTable *table = HostTableLoad (CATDIR, sky->hosts); + + int i; + for (i = 0; i < table->Nhosts; i++) { + + // ensure that the paths are absolute path names + char *tmppath = abspath (table->hosts[i].pathname, DVO_MAX_PATH); + free (table->hosts[i].pathname); + table->hosts[i].pathname = tmppath; + + char command[1024]; + snprintf (command, 1024, "setastrom_client %s -hostID %d -catdir %s -hostdir %s -region %f %f %f %f", + abscorrfile, table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, + UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax); + + fprintf (stderr, "command: %s\n", command); + + char tmpline[1024]; + if (VERBOSE) { snprintf (tmpline, 1024, "%s -v", command); strcpy (command, tmpline); } + if (UPDATE) { snprintf (tmpline, 1024, "%s -update", command); strcpy (command, tmpline); } + if (UPDATE_CATFORMAT) { snprintf (tmpline, 1024, "%s -update-catformat %s", command, UPDATE_CATFORMAT); strcpy (command, tmpline); } + + if (PARALLEL_MANUAL) continue; + + if (PARALLEL_SERIAL) { + int status = system (command); + if (status) { + fprintf (stderr, "ERROR running setastrom_client\n"); + exit (2); + } + } else { + // launch the job on the remote machine (no handshake) + int errorInfo = 0; + int pid = rconnect ("ssh", table->hosts[i].hostname, command, table->hosts[i].stdio, &errorInfo, FALSE); + if (!pid) { + if (DEBUG) fprintf (stderr, "failure to start %s (error %d)\n", table->hosts[i].hostname, errorInfo); + exit (1); + } + table->hosts[i].pid = pid; // save for future reference + } + } + free (abscorrfile); + + // wait for the remote jobs to be completed + if (PARALLEL_MANUAL) { + fprintf (stderr, "run the setastrom_client commands above. when these are done, hit return\n"); + getchar(); + } + if (!PARALLEL_MANUAL && !PARALLEL_SERIAL) { + int status = HostTableWaitJobsGetIO (table, __FILE__, __LINE__, VERBOSE); + if (!status) return FALSE; + } + + return (TRUE); +} Index: /trunk/Ohana/src/uniphot/src/update_dvo_setphot.c =================================================================== --- /trunk/Ohana/src/uniphot/src/update_dvo_setphot.c (revision 35102) +++ /trunk/Ohana/src/uniphot/src/update_dvo_setphot.c (revision 35103) @@ -4,5 +4,4 @@ int update_dvo_setphot (Image *image, off_t Nimage, FlatCorrectionTable *flatcorr) { - SkyRegion UserPatch; SkyTable *sky = NULL; SkyList *skylist = NULL; @@ -34,9 +33,10 @@ // determine the populated SkyRegions overlapping the requested area (default depth) - UserPatch.Rmin = 0; - UserPatch.Rmax = 360; - UserPatch.Dmin = -90; - UserPatch.Dmax = +90; - skylist = SkyListByPatch (sky, -1, &UserPatch); + if (SINGLE_CPT) { + skylist = SkyRegionByCPT (sky, SINGLE_CPT); + } else { + skylist = SkyListByPatch (sky, -1, &UserPatch); + } + myAssert (skylist, "ooops!"); // update measurements for each populated catalog @@ -110,6 +110,7 @@ char command[1024]; - snprintf (command, 1024, "setphot_client -hostID %d -catdir %s -hostdir %s -images %s", - table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, imageFile); + snprintf (command, 1024, "setphot_client -hostID %d -catdir %s -hostdir %s -images %s -region %f %f %f %f", + table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, imageFile, + UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax); fprintf (stderr, "command: %s\n", command); Index: /trunk/Ohana/src/uniphot/src/update_dvo_setposangle.c =================================================================== --- /trunk/Ohana/src/uniphot/src/update_dvo_setposangle.c (revision 35102) +++ /trunk/Ohana/src/uniphot/src/update_dvo_setposangle.c (revision 35103) @@ -3,5 +3,4 @@ int update_dvo_setposangle (ImageSubset *image, off_t Nimage) { - SkyRegion UserPatch; SkyTable *sky = NULL; SkyList *skylist = NULL; @@ -36,9 +35,10 @@ // determine the populated SkyRegions overlapping the requested area (default depth) - UserPatch.Rmin = 0; - UserPatch.Rmax = 360; - UserPatch.Dmin = -90; - UserPatch.Dmax = +90; - skylist = SkyListByPatch (sky, -1, &UserPatch); + if (SINGLE_CPT) { + skylist = SkyRegionByCPT (sky, SINGLE_CPT); + } else { + skylist = SkyListByPatch (sky, -1, &UserPatch); + } + myAssert (skylist, "ooops!"); // update measurements for each populated catalog @@ -113,6 +113,7 @@ char command[1024]; - snprintf (command, 1024, "setposangle_client -hostID %d -catdir %s -hostdir %s -images %s", - table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, imageFile); + snprintf (command, 1024, "setposangle_client -hostID %d -catdir %s -hostdir %s -images %s -region %f %f %f %f", + table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, imageFile, + UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax); fprintf (stderr, "command: %s\n", command);