Index: trunk/Ohana/src/uniphot/Makefile =================================================================== --- trunk/Ohana/src/uniphot/Makefile (revision 33653) +++ trunk/Ohana/src/uniphot/Makefile (revision 33654) @@ -1,3 +1,3 @@ -default: uniphot setphot setfwhm +default: uniphot setphot setphot_client setfwhm help: @echo "make options: uniphot setphot setfwhm default help install default (uniphot setphot setfwhm)" @@ -18,7 +18,9 @@ uniphot: $(BIN)/uniphot.$(ARCH) +setfwhm: $(BIN)/setfwhm.$(ARCH) setphot: $(BIN)/setphot.$(ARCH) -setfwhm: $(BIN)/setfwhm.$(ARCH) -install: $(DESTBIN)/uniphot $(DESTBIN)/setphot $(DESTBIN)/setfwhm +setphot_client: $(BIN)/setphot_client.$(ARCH) + +install: $(DESTBIN)/uniphot $(DESTBIN)/setfwhm $(DESTBIN)/setphot $(DESTBIN)/setphot_client UNIPHOT = \ @@ -33,6 +35,6 @@ $(SRC)/find_image_sgroups.$(ARCH).o \ $(SRC)/fit_groups.$(ARCH).o \ -$(SRC)/update.$(ARCH).o \ -$(SRC)/update_catalog.$(ARCH).o \ +$(SRC)/update_dvo_uniphot.$(ARCH).o \ +$(SRC)/update_catalog_uniphot.$(ARCH).o \ $(SRC)/convert.$(ARCH).o \ $(SRC)/SetSignals.$(ARCH).o \ @@ -45,17 +47,27 @@ SETPHOT = \ $(SRC)/setphot.$(ARCH).o \ -$(SRC)/initialize.$(ARCH).o \ -$(SRC)/ConfigInit.$(ARCH).o \ -$(SRC)/args.$(ARCH).o \ -$(SRC)/liststats.$(ARCH).o \ +$(SRC)/initialize_setphot.$(ARCH).o \ $(SRC)/load_zpt_table.$(ARCH).o \ -$(SRC)/load_images.$(ARCH).o \ +$(SRC)/load_images_setphot.$(ARCH).o \ +$(SRC)/ImageSubset.$(ARCH).o \ $(SRC)/match_zpts_to_images.$(ARCH).o \ +$(SRC)/match_flatcorr_to_images.$(ARCH).o \ +$(SRC)/update_dvo_setphot.$(ARCH).o \ $(SRC)/update_catalog_setphot.$(ARCH).o \ $(SRC)/SetSignals.$(ARCH).o \ $(SRC)/Shutdown.$(ARCH).o -$(SETPHOT): $(INC)/uniphot.h +$(SETPHOT): $(INC)/setphot.h $(BIN)/setphot.$(ARCH): $(SETPHOT) + +SETPHOT_CLIENT = \ +$(SRC)/setphot_client.$(ARCH).o \ +$(SRC)/update_dvo_setphot.$(ARCH).o \ +$(SRC)/update_catalog_setphot.$(ARCH).o \ +$(SRC)/initialize_setphot_client.$(ARCH).o \ +$(SRC)/ImageSubset.$(ARCH).o + +$(SETPHOT_CLIENT): $(INC)/setphot.h +$(BIN)/setphot_client.$(ARCH): $(SETPHOT_CLIENT) SETFWHM = \ Index: trunk/Ohana/src/uniphot/include/setphot.h =================================================================== --- trunk/Ohana/src/uniphot/include/setphot.h (revision 33654) +++ trunk/Ohana/src/uniphot/include/setphot.h (revision 33654) @@ -0,0 +1,74 @@ +# include +# include +# include + +typedef struct { + float zpt; + float zpt_err; + e_time time; + int found; +} ZptTable; + +typedef struct { + float Mcal; + float dMcal; + unsigned int imageID; + unsigned int photom_map_id; + unsigned int flags; +} ImageSubset; + +/* global variables set in parameter file */ +# define MAX_PATH_LENGTH 1024 +char ImageCat[MAX_PATH_LENGTH]; +char *CATDIR; +int HOST_ID; +char *HOSTDIR; +char *IMAGES; +int VERBOSE; +int RESET; +int UBERCAL; // load the supplied ubercal zero point fits table (with flat-field corrections) +int NO_METADATA; // the supplied ubercal data has no descriptive metadata +int UPDATE; +int PARALLEL; +int PARALLEL_MANUAL; +int PARALLEL_SERIAL; +int IMAGES_ONLY; + +/***** prototypes ****/ +int main PROTO((int argc, char **argv)); + +void ConfigInit PROTO((int *argc, char **argv)); +void initialize_setphot PROTO((int argc, char **argv)); +int args_setphot PROTO((int argc, char **argv)); + +void initialize_setphot_client PROTO((int argc, char **argv)); +int args_setphot_client PROTO((int argc, char **argv)); +int update_dvo_setphot_client PROTO((ImageSubset *image, off_t Nimage, FlatCorrectionTable *flatcorr)); + +ImageSubset *ImageSubsetLoad PROTO((char *filename, off_t *nimage)); +int ImageSubsetSave PROTO((char *filename, Image *image, off_t Nimage)); +Image *ImagesFromSubset PROTO((ImageSubset *subset, off_t N)); + +void lock_image_db PROTO((FITS_DB *db, char *filename)); +void unlock_image_db PROTO((FITS_DB *db)); +void create_image_db PROTO((FITS_DB *db)); +void set_db PROTO((FITS_DB *in)); +int Shutdown PROTO((char *format, ...) OHANA_FORMAT(printf, 1, 2) ); +void TrapSignal PROTO((int sig)); +void SetProtect PROTO((int mode)); +int SetSignals PROTO((void)); + +// setphot-specific prototypes +ZptTable *load_zpt_table PROTO((char *filename, int *nzpts)); +Image *load_images_setphot PROTO((FITS_DB *db, off_t *Nimage)); +int match_zpts_to_images PROTO((Image *image, off_t Nimage, ZptTable *zpts, int Nzpts)); +int update_dvo_setphot PROTO((Image *image, off_t Nimage, FlatCorrectionTable *flatcorr)); +void update_catalog_setphot PROTO((Catalog *catalog, Image *image, off_t *index, off_t Nimage, FlatCorrectionTable *flatcorr)); +void update_catalog_setphot_client PROTO((Catalog *catalog, ImageSubset *image, off_t *index, off_t Nimage, FlatCorrectionTable *flatcorr)); + +ZptTable *load_zpt_ubercal PROTO((char *filename, int *nzpts, FlatCorrectionTable *flatcorrTable)); +ZptTable *load_zpt_ubercal_nometadata PROTO((char *filename, int *nzpts, FlatCorrectionTable *flatcorrTable)); +int match_flatcorr_to_images PROTO((Image *image, off_t Nimage, FlatCorrectionTable *flatcorrTable)); + +int parse_zpt_offsets PROTO((char *ZPT_OFFSET_FILTERS, char *ZPT_OFFSET_VALUES)); +float apply_zpt_offset PROTO((short code)); Index: trunk/Ohana/src/uniphot/include/setphot_client.h =================================================================== --- trunk/Ohana/src/uniphot/include/setphot_client.h (revision 33654) +++ trunk/Ohana/src/uniphot/include/setphot_client.h (revision 33654) @@ -0,0 +1,168 @@ +# include +# include +# include + +# ifdef ANSI +# define F_SETFL 4 +# define O_NONBLOCK 0200000 +# define AF_UNIX 1 +# define SOCK_STREAM 1 +# define ENOENT 2 /* No such file or directory */ +# endif /* ANSI */ + +typedef struct { + double median; + double mean; + double sigma; + double error; + double chisq; + double min; + double max; + double total; + int Nmeas; +} StatType; + +typedef struct { + void *tgroup; + void *sgroup; +} ImageLink; + +typedef struct { + char tstart[64]; + char tstop[64]; + char label[64]; + float M; + float dM; + float dMsub; + double v1, v2; + Image **image; + ImageLink **imlink; + int Nimage, Ngood; +} Group; + +typedef struct { + float zpt; + float zpt_err; + e_time time; + int found; +} ZptTable; + +typedef struct { + float fwhm_major; + float fwhm_minor; + e_time time; + int found; + unsigned short photcode; +} FWHMTable; + +/* global variables set in parameter file */ +char ImageCat[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; /** XXX EAM : depth of catalog tables, fix usage */ +char GSCFILE[256]; +char STATMODE[64]; +int VERBOSE; +int UBERCAL; // load the supplied ubercal zero point fits table (with flat-field corrections) +int NO_METADATA; // the supplied ubercal data has no descriptive metadata +int NLOOP; +int TimeSelect; +int VERBOSE; +int UPDATE; +int IMAGE_BAD; +double RADIUS; +double TRANGE; +time_t TSTART; +time_t TSTOP; +PhotCode *photcode; + +enum {black, white, red, orange, yellow, green, blue, indigo, violet}; + +typedef struct { + double xmin, xmax, ymin, ymax; + int style, ptype, ltype, etype, color; + double lweight, size; +} UPGraphdata; + +/***** prototypes ****/ +void ConfigInit PROTO((int *argc, char **argv)); +void DonePlotting PROTO((UPGraphdata *graphmode, int N)); +void JpegPlot PROTO((UPGraphdata *graphmode, int N, char *filename)); +void PSPlot PROTO((UPGraphdata *graphmode, int N, char *filename)); +void PlotLabel PROTO((char *string, int N)); +void PlotVector PROTO((int Npts, double *vect, int mode, int N)); +void PrepPlotting PROTO((int Npts, UPGraphdata *graphmode, int N)); +void XDead PROTO((void)); +int args_uniphot PROTO((int argc, char **argv)); +int args_setphot PROTO((int argc, char **argv)); +void dumpresult PROTO((void)); +Group *find_image_sgroups PROTO((FITS_DB *db, ImageLink **imlink, int *Nsgroup)); +Group *find_image_tgroups PROTO((FITS_DB *db, ImageLink **imlink, int *Ntgroup)); +void fit_sgroup PROTO((Group *sgroup, int Nsgroup)); +void fit_tgroup PROTO((Group *tgroup, int Ntgroup)); +int gcatalog PROTO((Catalog *catalog)); +void initialize_uniphot PROTO((int argc, char **argv)); +void initialize_setphot PROTO((int argc, char **argv)); +void initstats PROTO((char *mode)); +int liststats PROTO((double *value, double *dvalue, int N, StatType *stats)); +int load_images_uniphot PROTO((FITS_DB *db)); +int main PROTO((int argc, char **argv)); +int open_graph PROTO((int N)); +void sort PROTO((unsigned int *X, int N)); +void sortB PROTO((double *X, double *Y, int N)); +void sortD PROTO((double *X, double *Y, double *Z, int N)); +void wcatalog PROTO((Catalog *catalog)); +void wimages PROTO((Image *image, int Nimage)); +void check_permissions PROTO((char *basefile)); +void lock_image_db PROTO((FITS_DB *db, char *filename)); +void unlock_image_db PROTO((FITS_DB *db)); +void create_image_db PROTO((FITS_DB *db)); +void set_db PROTO((FITS_DB *in)); +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 subset_images PROTO((FITS_DB *db)); +void update_dvo_uniphot PROTO((FITS_DB *db, Group *sgroup, int Nsgroup)); +void update_catalog_uniphot PROTO((Catalog *catalog, Group *sgroup, int warn)); +void sort_time PROTO((unsigned int *value, int N)); + +// setphot-specific prototypes +ZptTable *load_zpt_table PROTO((char *filename, int *nzpts)); +Image *load_images_setphot PROTO((FITS_DB *db, off_t *Nimage)); +int match_zpts_to_images PROTO((Image *image, off_t Nimage, ZptTable *zpts, int Nzpts)); +int update_dvo_setphot PROTO((Image *image, off_t Nimage, FlatCorrectionTable *flatcorr)); +void update_catalog_setphot PROTO((Catalog *catalog, Image *image, off_t *index, off_t Nimage, FlatCorrectionTable *flatcorr)); + +ZptTable *load_zpt_ubercal (char *filename, int *nzpts, FlatCorrectionTable *flatcorrTable); +ZptTable *load_zpt_ubercal_nometadata (char *filename, int *nzpts, FlatCorrectionTable *flatcorrTable); +int match_flatcorr_to_images (Image *image, off_t Nimage, FlatCorrectionTable *flatcorrTable); + +/*** time/coord conversion functions not supplied by libohana ***/ +time_t TimeRef PROTO((double time, time_t TimeReference, int TimeFormat)); +double TimeValue PROTO((time_t time, time_t TimeReference, int TimeFormat)); + +int hh_hms PROTO((double hh, int *hr, int *mn, double *sc)); +int dd_dms PROTO((double dd, int *dg, int *mn, double *sc)); +int hms_format PROTO((char *line, double value)); +int dms_format PROTO((char *line, double value)); +int hh_hm PROTO((double hh, int *hr, double *mn)); +int day_to_sec PROTO((char *string, time_t *second)); +int hms_to_sec PROTO((char *string, time_t *second)); +char *ohana_sec_to_hms PROTO((time_t second)); +char *ohana_sec_to_day PROTO((time_t second)); + +char *meade_deg_to_str PROTO((double deg)); +char *meade_ra_to_str PROTO((double deg)); +char *meade_dec_to_str PROTO((double deg)); +char *strptime PROTO((const char *s, const char *format, struct tm *tm)); +time_t GetTimeReference PROTO((char *reference)); +int GetTimeUnits PROTO((char *name)); + +void initialize_setfwhm PROTO((int argc, char **argv)); +int args_setfwhm PROTO((int argc, char **argv)); +FWHMTable *load_fwhm_table PROTO((char *filename, int *nfwhm)); +int match_fwhm_to_images PROTO((Image *image, off_t Nimage, FWHMTable *fwhm, int Nfwhm)); + Index: trunk/Ohana/src/uniphot/include/uniphot.h =================================================================== --- trunk/Ohana/src/uniphot/include/uniphot.h (revision 33653) +++ trunk/Ohana/src/uniphot/include/uniphot.h (revision 33654) @@ -66,4 +66,6 @@ char STATMODE[64]; int VERBOSE; +int UBERCAL; // load the supplied ubercal zero point fits table (with flat-field corrections) +int NO_METADATA; // the supplied ubercal data has no descriptive metadata int NLOOP; int TimeSelect; @@ -95,5 +97,4 @@ void XDead PROTO((void)); int args_uniphot PROTO((int argc, char **argv)); -int args_setphot PROTO((int argc, char **argv)); void dumpresult PROTO((void)); Group *find_image_sgroups PROTO((FITS_DB *db, ImageLink **imlink, int *Nsgroup)); @@ -103,8 +104,8 @@ int gcatalog PROTO((Catalog *catalog)); void initialize_uniphot PROTO((int argc, char **argv)); -void initialize_setphot PROTO((int argc, char **argv)); void initstats PROTO((char *mode)); int liststats PROTO((double *value, double *dvalue, int N, StatType *stats)); int load_images_uniphot PROTO((FITS_DB *db)); +Image *load_images_setfwhm PROTO((FITS_DB *db, off_t *Nimage)); int main PROTO((int argc, char **argv)); int open_graph PROTO((int N)); @@ -112,5 +113,4 @@ void sortB PROTO((double *X, double *Y, int N)); void sortD PROTO((double *X, double *Y, double *Z, int N)); -void update_dvo_catalog PROTO((Catalog *catalog, Group *sgroup, int warn)); void wcatalog PROTO((Catalog *catalog)); void wimages PROTO((Image *image, int Nimage)); @@ -125,32 +125,8 @@ int SetSignals PROTO((void)); int subset_images PROTO((FITS_DB *db)); -void update PROTO((FITS_DB *db, Group *sgroup, int Nsgroup)); +void update_dvo_uniphot PROTO((FITS_DB *db, Group *sgroup, int Nsgroup)); +void update_catalog_uniphot PROTO((Catalog *catalog, Group *sgroup, int warn)); void sort_time PROTO((unsigned int *value, int N)); -// setphot-specific prototypes -ZptTable *load_zpt_table PROTO((char *filename, int *nzpts)); -Image *load_images_setphot PROTO((FITS_DB *db, off_t *Nimage)); -int match_zpts_to_images PROTO((Image *image, off_t Nimage, ZptTable *zpts, int Nzpts)); -int update_setphot PROTO((Image *image, off_t Nimage)); -void update_catalog_setphot PROTO((Catalog *catalog, Image *image, off_t *index, off_t Nimage)); - -/*** time/coord conversion functions not supplied by libohana ***/ -time_t TimeRef PROTO((double time, time_t TimeReference, int TimeFormat)); -double TimeValue PROTO((time_t time, time_t TimeReference, int TimeFormat)); - -int hh_hms PROTO((double hh, int *hr, int *mn, double *sc)); -int dd_dms PROTO((double dd, int *dg, int *mn, double *sc)); -int hms_format PROTO((char *line, double value)); -int dms_format PROTO((char *line, double value)); -int hh_hm PROTO((double hh, int *hr, double *mn)); -int day_to_sec PROTO((char *string, time_t *second)); -int hms_to_sec PROTO((char *string, time_t *second)); -char *ohana_sec_to_hms PROTO((time_t second)); -char *ohana_sec_to_day PROTO((time_t second)); - -char *meade_deg_to_str PROTO((double deg)); -char *meade_ra_to_str PROTO((double deg)); -char *meade_dec_to_str PROTO((double deg)); -char *strptime PROTO((const char *s, const char *format, struct tm *tm)); time_t GetTimeReference PROTO((char *reference)); int GetTimeUnits PROTO((char *name)); Index: trunk/Ohana/src/uniphot/src/ConfigInit.c =================================================================== --- trunk/Ohana/src/uniphot/src/ConfigInit.c (revision 33653) +++ trunk/Ohana/src/uniphot/src/ConfigInit.c (revision 33654) @@ -7,6 +7,4 @@ double ZERO_POINT; char *config, *file; - char CatdirPhotcodeFile[256]; - char MasterPhotcodeFile[256]; /*** load configuration info ***/ @@ -23,5 +21,4 @@ ScanConfig (config, "CATMODE", "%s", 0, CATMODE); ScanConfig (config, "CATFORMAT", "%s", 0, CATFORMAT); - ScanConfig (config, "PHOTCODE_FILE", "%s", 0, MasterPhotcodeFile); sprintf (ImageCat, "%s/Images.dat", CATDIR); @@ -43,12 +40,4 @@ 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, TRUE)) { - fprintf (stderr, "error loading photcode table %s or master file %s\n", CatdirPhotcodeFile, MasterPhotcodeFile); - exit (1); - } - } Index: trunk/Ohana/src/uniphot/src/ImageSubset.c =================================================================== --- trunk/Ohana/src/uniphot/src/ImageSubset.c (revision 33654) +++ trunk/Ohana/src/uniphot/src/ImageSubset.c (revision 33654) @@ -0,0 +1,183 @@ +# include "setphot.h" + +# define GET_COLUMN(OUT,NAME,TYPE) \ + TYPE *OUT = gfits_get_bintable_column_data (&theader, &ftable, NAME, type, &Nrow, &Ncol); \ + myAssert (!strcmp(type, #TYPE), "wrong column type"); + +ImageSubset *ImageSubsetLoad(char *filename, off_t *nimage) { + + int i, Ncol; + off_t Nrow; + Header header; + Header theader; + Matrix matrix; + FTable ftable; + + *nimage = 0; + ImageSubset *image = NULL; + + FILE *f = fopen (filename, "r"); + if (!f) { + fprintf (stderr, "ERROR: cannot open image subset file %s\n", filename); + return NULL; + } + + /* load in PHU segment (ignore) */ + if (!gfits_fread_header (f, &header)) { + if (VERBOSE) fprintf (stderr, "can't read image subset header\n"); + fclose (f); + return NULL; + } + if (!gfits_fread_matrix (f, &matrix, &header)) { + if (VERBOSE) fprintf (stderr, "can't read image subset matrix\n"); + gfits_free_header (&header); + fclose (f); + return NULL; + } + + ftable.header = &theader; + + // load data for this header + if (!gfits_load_header (f, &theader)) { + fclose (f); + return NULL; + } + + // read the fits table bytes + if (!gfits_fread_ftable_data (f, &ftable, FALSE)) { + fclose (f); + return NULL; + } + fclose (f); + + // a bit annoying : we read the entire block of data (180M), then extract the columns, then set the image structure values. + // this means I need 3 copies in memory at some point. ugh. + + char type[16]; + + GET_COLUMN (Mcal, "MCAL", float); + GET_COLUMN (dMcal, "MCAL_ERR", float); + GET_COLUMN (imageID, "IMAGE_ID", int); + GET_COLUMN (map, "PHOTOM_MAP", int); + GET_COLUMN (flags, "FLAGS", int); + + // XXX free the fits table data here + + ALLOCATE (image, ImageSubset, Nrow); + for (i = 0; i < Nrow; i++) { + image[i].imageID = imageID[i]; + image[i].photom_map_id = map[i]; + image[i].flags = flags[i]; + image[i].Mcal = Mcal[i]; + image[i].dMcal = dMcal[i]; + } + fprintf (stderr, "loaded data for %lld images\n", (long long) Nrow); + + free (Mcal); + free (dMcal); + free (imageID); + free (map); + free (flags); + + *nimage = Nrow; + return image; +} + +int ImageSubsetSave(char *filename, Image *image, off_t Nimage) { + + int i; + Header header; + Header theader; + Matrix matrix; + FTable ftable; + + gfits_init_header (&header); + header.extend = TRUE; + gfits_create_header (&header); + gfits_create_matrix (&header, &matrix); + + gfits_create_table_header (&theader, "BINTABLE", "IMAGE_SUBSET"); + + gfits_define_bintable_column (&theader, "E", "MCAL", "zero point offset", "magnitudes", 1.0, 0.0); + gfits_define_bintable_column (&theader, "E", "MCAL_ERR", "zero point error", "magnitudes", 1.0, 0.0); + + // an unsigned int needs to have bzero of 0x8000 + gfits_define_bintable_column (&theader, "J", "IMAGE_ID", "image ID", NULL, 1.0, 1.0*0x8000); + gfits_define_bintable_column (&theader, "J", "PHOTOM_MAP", "map", NULL, 1.0, 1.0*0x8000); + gfits_define_bintable_column (&theader, "J", "FLAGS", "flags", NULL, 1.0, 1.0*0x8000); + + // generate the output array that carries the data + gfits_create_table (&theader, &ftable); + + float *Mcal, *dMcal; + unsigned int *imageID, *map, *flags; + + // create intermediate storage arrays + ALLOCATE (Mcal, float, Nimage); + ALLOCATE (dMcal, float, Nimage); + ALLOCATE (imageID, unsigned int, Nimage); + ALLOCATE (map, unsigned int, Nimage); + ALLOCATE (flags, unsigned int, Nimage); + + // assign the storage arrays + for (i = 0; i < Nimage; i++) { + imageID[i] = image[i].imageID; + map[i] = image[i].photom_map_id; + flags[i] = image[i].flags; + Mcal[i] = image[i].Mcal; + dMcal[i] = image[i].dMcal; + } + + // add the columns to the output array + gfits_set_bintable_column (&theader, &ftable, "MCAL", Mcal, Nimage); + gfits_set_bintable_column (&theader, &ftable, "MCAL_ERR", dMcal, Nimage); + gfits_set_bintable_column (&theader, &ftable, "IMAGE_ID", imageID, Nimage); + gfits_set_bintable_column (&theader, &ftable, "PHOTOM_MAP", map, Nimage); + gfits_set_bintable_column (&theader, &ftable, "FLAGS", flags, Nimage); + + free (Mcal); + free (dMcal); + free (imageID); + free (map); + free (flags); + + FILE *f = fopen (filename, "w"); + if (!f) { + fprintf (stderr, "ERROR: cannot open image subset file for output %s\n", filename); + return FALSE; + } + + gfits_fwrite_header (f, &header); + gfits_fwrite_matrix (f, &matrix); + gfits_fwrite_Theader (f, &theader); + gfits_fwrite_table (f, &ftable); + + gfits_free_header (&header); + gfits_free_matrix (&matrix); + gfits_free_header (&theader); + gfits_free_table (&ftable); + + fclose (f); + fflush (f); + + return TRUE; +} + +Image *ImagesFromSubset (ImageSubset *subset, off_t N) { + + off_t i; + + // we have been given an ImageSubset array, containing a reduced set of image fields + // create full a Image array and save the needed values + Image *image = NULL; + ALLOCATE (image, Image, N); + for (i = 0; i < N; i++) { + image[i].imageID = subset[i].imageID ; + image[i].photom_map_id = subset[i].photom_map_id; + image[i].flags = subset[i].flags ; + image[i].Mcal = subset[i].Mcal ; + image[i].dMcal = subset[i].dMcal ; + } + return image; +} + Index: trunk/Ohana/src/uniphot/src/args.c =================================================================== --- trunk/Ohana/src/uniphot/src/args.c (revision 33653) +++ trunk/Ohana/src/uniphot/src/args.c (revision 33654) @@ -72,28 +72,4 @@ } -int args_setphot (int argc, char **argv) { - - int N; - - 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; - } - - if (argc != 2) { - fprintf (stderr, "ERROR: USAGE: setphot (zptfile) [options]\n"); - exit (2); - } - - return (TRUE); -} - int args_setfwhm (int argc, char **argv) { Index: trunk/Ohana/src/uniphot/src/convert.c =================================================================== --- trunk/Ohana/src/uniphot/src/convert.c (revision 33653) +++ trunk/Ohana/src/uniphot/src/convert.c (revision 33654) @@ -1,7 +1,3 @@ # include "uniphot.h" -# define _XOPEN_SOURCE /* glibc2 (strptime) needs this */ -# include - -/** additional time / coordinate conversions not supplied by libohana **/ time_t GetTimeReference (char *reference) { @@ -36,302 +32,2 @@ return (Units); } - -int hh_hms (double hh, int *hr, int *mn, double *sc) { - - int N, flag; - - flag = SIGN(hh); - hh = fabs(hh); - - // rationalize hh to range -24.0 < hh < 24.0 - if (hh >= 24.0) { - N = (int)(hh/24.0); - hh -= 24.0*N; - } - - *hr = (int) hh; - *mn = (int) 60*(hh - *hr); - *sc = 3600.0*(hh - *hr - *mn / 60.0); - if (*sc > 59.99) { - *sc = 0.0; - *mn += 1.0; - } - *hr *= flag; - return (TRUE); -} - -int dd_dms (double dd, int *dg, int *mn, double *sc) { - - int flag; - - flag = SIGN(dd); - dd = fabs (dd); - *dg = (int) dd; - *mn = (int) 60*(dd - *dg); - *sc = 3600.0*(dd - *dg - *mn/60.0); - if (*sc > 59.99) { - *sc = 0; - *mn += 1.0; - } - *dg *= flag; - return (TRUE); -} - -int hms_format (char *line, double value) { - - int hr, mn; - double sc; - - hh_hms (value, &hr, &mn, &sc); - hr = (int) value; - if (isnan (value)) - sprintf (line, "xx:xx:xx.xx"); - else { - if (value < 0) { - sprintf (line, "-%02d:%02d:%05.2f", abs(hr), mn, sc); - } else { - sprintf (line, "+%02d:%02d:%05.2f", hr, mn, sc); - } - } - return (TRUE); -} - -int dms_format (char *line, double value) { - - int dg, mn; - double sc; - - dd_dms (value, &dg, &mn, &sc); - if (value < 0) { - sprintf (line, "-%02d:%02d:%05.2f", abs(dg), mn, sc); - } else { - sprintf (line, "+%02d:%02d:%05.2f", dg, mn, sc); - } - return (TRUE); -} - -/***** convert 00:00:00 or 00:00 to 0 - 86400 ****/ -int hms_to_sec (char *string, time_t *second) { - - char *p; - struct tm time; - - p = strptime (string, "%H:%M:%S", &time); - if (p != NULL) goto valid; - - p = strptime (string, "%H:%M", &time); - if (p != NULL) goto valid; - - return (FALSE); - -valid: - if (*p) return (FALSE); - *second = time.tm_hour*3600 + time.tm_min*60 + time.tm_sec; - return (TRUE); -} - -/***** convert Mon[@00:00:00] or 00:00 to 0 - 86400*7 ****/ -int day_to_sec (char *string, time_t *second) { - - char *p; - struct tm time; - - bzero (&time, sizeof(time)); - p = strptime (string, "%A@%H:%M:%S", &time); - if (p != NULL) goto valid; - - p = strptime (string, "%A@%H:%M", &time); - if (p != NULL) goto valid; - - p = strptime (string, "%A@%H", &time); - if (p != NULL) goto valid; - - p = strptime (string, "%A", &time); - if (p != NULL) goto valid; - - return (FALSE); - -valid: - if (*p) return (FALSE); - *second = time.tm_wday*86400 + time.tm_hour*3600 + time.tm_min*60 + time.tm_sec; - return (TRUE); -} - -/***** convert seconds to HH:MM:SS ****/ -char *ohana_sec_to_hms (time_t second) { - - struct tm *gmt; - char *line; - - ALLOCATE (line, char, 64); - gmt = gmtime (&second); - sprintf (line, "%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); - return (line); -} - -/***** convert seconds to Day@HH:MM:SS ****/ -char *ohana_sec_to_day (time_t second) { - - struct tm *gmt; - char *line; - - ALLOCATE (line, char, 64); - gmt = gmtime (&second); - switch (gmt[0].tm_wday) { - case 0: - sprintf (line, "Sun@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); - break; - case 1: - sprintf (line, "Mon@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); - break; - case 2: - sprintf (line, "Tue@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); - break; - case 3: - sprintf (line, "Wed@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); - break; - case 4: - sprintf (line, "Thu@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); - break; - case 5: - sprintf (line, "Fri@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); - break; - case 6: - sprintf (line, "Sat@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); - break; - } - return (line); -} - -int hh_hm (double hh, int *hr, double *mn) { - - int flag; - - flag = SIGN(hh); - hh = fabs (hh); - - *mn = 60.0*(hh - (int)hh); - *hr = (int) hh; - *hr *= flag; - return (TRUE); -} - -char *meade_deg_to_str (double deg) { - - int hr; - double mn; - char *line; - - ALLOCATE (line, char, 16); - - hh_hm (deg, &hr, &mn); - - sprintf (line, "%03d:%04.1f", abs(hr), mn); - return (line); -} - -char *meade_ra_to_str (double deg) { - - int hr; - double mn; - char *line; - - ALLOCATE (line, char, 16); - - hh_hm (deg/15.0, &hr, &mn); - - sprintf (line, "%02d:%04.1f", abs(hr), mn); - return (line); -} - -char *meade_dec_to_str (double deg) { - - int hr; - double mn; - char *line; - - ALLOCATE (line, char, 16); - - hh_hm (deg, &hr, &mn); - - if (deg < 0) { - sprintf (line, "-%02d:%04.1f", abs(hr), mn); - } else { - sprintf (line, "+%02d:%04.1f", hr, mn); - } - return (line); -} - -/* convert UNIX time to a value referenced to the TimeReference in the given unit */ -double TimeValue (time_t time, time_t TimeReference, int TimeFormat) { - - double value, dt; - - dt = (time > TimeReference) ? (time - TimeReference) : -1 * (double)(TimeReference - time); - switch (TimeFormat) { - case TIME_JD: - value = time / 86400.0 + 2440587.5; - break; - case TIME_MJD: - value = time / 86400.0 + 40587.0; - break; - case TIME_DAYS: - value = dt / 86400.0; - break; - case TIME_HOURS: - value = dt / 3600.0; - break; - case TIME_MINUTES: - value = dt / 60.0; - break; - case TIME_SECONDS: - default: - value = dt; - break; - } - return (value); -} - -/* convert time value referenced to the TimeReference in the given unit to UNIX time */ -time_t TimeRef (double value, time_t TimeReference, int TimeFormat) { - - int dt; - time_t time; - - switch (TimeFormat) { - case TIME_JD: - time = (value - 2440587.5) * 86400.0; - return (time); - break; - case TIME_MJD: - time = (value - 40587.0) * 86400.0; - return (time); - break; - case TIME_DAYS: - dt = value * 86400.0; - break; - case TIME_HOURS: - dt = value * 3600.0; - break; - case TIME_MINUTES: - dt = value * 60.0; - break; - case TIME_SECONDS: - default: - dt = value; - break; - } - - time = TimeReference + dt; - return (time); -} - -/* times may be in forms as: - * 20040200450s (N seconds since 1970.0) - * 2440900.232j (julian date) - * 99/02/23,03:22:18 (date string) - * (separators may be anything except space, +, -) - * 99:02:15:12:23:30 - * 99:02:15:12h23m30s - */ - Index: trunk/Ohana/src/uniphot/src/initialize.c =================================================================== --- trunk/Ohana/src/uniphot/src/initialize.c (revision 33653) +++ trunk/Ohana/src/uniphot/src/initialize.c (revision 33654) @@ -21,11 +21,4 @@ } -void initialize_setphot (int argc, char **argv) { - - /* are these set correctly? */ - ConfigInit (&argc, argv); - args_setphot (argc, argv); -} - void initialize_setfwhm (int argc, char **argv) { Index: trunk/Ohana/src/uniphot/src/initialize_setphot.c =================================================================== --- trunk/Ohana/src/uniphot/src/initialize_setphot.c (revision 33654) +++ trunk/Ohana/src/uniphot/src/initialize_setphot.c (revision 33654) @@ -0,0 +1,166 @@ +# include "setphot.h" + +static char *ZPT_OFFSETS_FILTERS = NULL; +static char *ZPT_OFFSETS_VALUES = NULL; + +void initialize_setphot (int argc, char **argv) { + + struct stat statbuffer; + + char CatdirPhotcodeFile[256]; + + /* are these set correctly? */ + ConfigInit (&argc, argv); + args_setphot (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, TRUE)) { + fprintf (stderr, "error loading photcode table %s\n", CatdirPhotcodeFile); + exit (1); + } + + // if we have a zpt offset list, parse it here + parse_zpt_offsets (ZPT_OFFSETS_FILTERS, ZPT_OFFSETS_VALUES); +} + +int args_setphot (int argc, char **argv) { + + int N; + + VERBOSE = FALSE; + if ((N = get_argument (argc, argv, "-v"))) { + VERBOSE = TRUE; + remove_argument (N, &argc, argv); + } + + IMAGES_ONLY = FALSE; + if ((N = get_argument (argc, argv, "-images-only"))) { + IMAGES_ONLY = TRUE; + remove_argument (N, &argc, argv); + } + + RESET = FALSE; + if ((N = get_argument (argc, argv, "-reset"))) { + RESET = 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, + // setphot 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 setphot_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); + } + + UBERCAL = FALSE; + if ((N = get_argument (argc, argv, "-ubercal"))) { + remove_argument (N, &argc, argv); + UBERCAL = TRUE; + } + + NO_METADATA = FALSE; + if ((N = get_argument (argc, argv, "-no-metadata"))) { + remove_argument (N, &argc, argv); + NO_METADATA = TRUE; + } + + UPDATE = FALSE; + if ((N = get_argument (argc, argv, "-update"))) { + remove_argument (N, &argc, argv); + UPDATE = TRUE; + } + + // -zpt-offsets code[,code,etc] value[,value,etc] + // eg: -zpt-offsets g,r,i 0.01,0.02,-0.02 + if ((N = get_argument (argc, argv, "-zpt-offsets"))) { + remove_argument (N, &argc, argv); + ZPT_OFFSETS_FILTERS = strcreate (argv[N]); + remove_argument (N, &argc, argv); + ZPT_OFFSETS_VALUES = strcreate (argv[N]); + remove_argument (N, &argc, argv); + } + + + if (argc != 2) { + fprintf (stderr, "USAGE: setphot (zptfile) [options]\n"); + fprintf (stderr, " options:\n"); + fprintf (stderr, " -v : verbose mode\n"); + fprintf (stderr, " -ubercal : interpret zpt file as a GPC1 ubercal FITS table\n"); + fprintf (stderr, " -no-metadata : assume the ubercal FITS table has the layout defined by Eddie Schlafly in Dec 2011\n"); + fprintf (stderr, " -update : actually write results to detections tables\n"); + fprintf (stderr, " Note that the dvo db can be specified by -D CATDIR (directory)\n"); + exit (2); + } + + return (TRUE); +} + +void GetConfig (char *config, char *field, char *format, int N, void *ptr); + +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, MAX_PATH_LENGTH); + GetConfig (config, "CATDIR", "%s", 0, tmpcatdir); + CATDIR = abspath (tmpcatdir, MAX_PATH_LENGTH); + 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_client.c =================================================================== --- trunk/Ohana/src/uniphot/src/initialize_setphot_client.c (revision 33654) +++ trunk/Ohana/src/uniphot/src/initialize_setphot_client.c (revision 33654) @@ -0,0 +1,95 @@ +# include "setphot.h" + +void initialize_setphot_client (int argc, char **argv) { + + struct stat statbuffer; + + char CatdirPhotcodeFile[256]; + + SetZeroPoint (25.0); // XXX is this needed? + args_setphot_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, TRUE)) { + fprintf (stderr, "error loading photcode table %s\n", CatdirPhotcodeFile); + exit (1); + } +} + +int args_setphot_client (int argc, char **argv) { + + int N; + + VERBOSE = FALSE; + if ((N = get_argument (argc, argv, "-v"))) { + VERBOSE = TRUE; + remove_argument (N, &argc, argv); + } + + RESET = FALSE; + if ((N = get_argument (argc, argv, "-reset"))) { + RESET = TRUE; + remove_argument (N, &argc, argv); + } + + UBERCAL = FALSE; + if ((N = get_argument (argc, argv, "-ubercal"))) { + remove_argument (N, &argc, argv); + UBERCAL = TRUE; + } + + 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) goto usage; + + 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) goto usage; + + if ((N = get_argument (argc, argv, "-catdir"))) { + remove_argument (N, &argc, argv); + CATDIR = strcreate (argv[N]); + remove_argument (N, &argc, argv); + } + if (!CATDIR) goto usage; + + IMAGES = NULL; + if ((N = get_argument (argc, argv, "-images"))) { + remove_argument (N, &argc, argv); + IMAGES = strcreate (argv[N]); + remove_argument (N, &argc, argv); + } + if (!IMAGES) goto usage; + + if (argc != 1) goto usage; + return (TRUE); + +usage: + fprintf (stderr, "USAGE: setphot_client -hostID (hostID) -catdir (catdir) -hostdir (hostdir) -images (images) [options]\n"); + fprintf (stderr, " options:\n"); + fprintf (stderr, " -v : verbose mode\n"); + fprintf (stderr, " -ubercal : interpret zpt file as a GPC1 ubercal FITS table\n"); + fprintf (stderr, " -update : actually write results to detections tables\n"); + exit (2); +} + Index: trunk/Ohana/src/uniphot/src/load_images.c =================================================================== --- trunk/Ohana/src/uniphot/src/load_images.c (revision 33653) +++ trunk/Ohana/src/uniphot/src/load_images.c (revision 33654) @@ -11,5 +11,5 @@ } -Image *load_images_setphot (FITS_DB *db, off_t *Nimage) { +Image *load_images_setfwhm (FITS_DB *db, off_t *Nimage) { Image *image; Index: trunk/Ohana/src/uniphot/src/load_images_setphot.c =================================================================== --- trunk/Ohana/src/uniphot/src/load_images_setphot.c (revision 33654) +++ trunk/Ohana/src/uniphot/src/load_images_setphot.c (revision 33654) @@ -0,0 +1,22 @@ +# include "setphot.h" + +Image *load_images_setphot (FITS_DB *db, off_t *Nimage) { + + Image *image; + + if (VERBOSE) fprintf (stderr, "finding images\n"); + + /* read entire db table */ + if (!dvo_image_load (db, VERBOSE, FALSE)) Shutdown ("can't read image catalog %s", db[0].filename); + + /* use a vtable to keep the images to be calibrated */ + image = gfits_table_get_Image (&db[0].ftable, Nimage, &db[0].swapped); + if (!image) { + fprintf (stderr, "ERROR: failed to read images\n"); + exit (2); + } + + fprintf (stderr, "loaded "OFF_T_FMT" images\n", *Nimage); + + return (image); +} Index: trunk/Ohana/src/uniphot/src/load_zpt_table.c =================================================================== --- trunk/Ohana/src/uniphot/src/load_zpt_table.c (revision 33653) +++ trunk/Ohana/src/uniphot/src/load_zpt_table.c (revision 33654) @@ -1,3 +1,3 @@ -# include "uniphot.h" +# include "setphot.h" ZptTable *load_zpt_table (char *filename, int *nzpts) { @@ -41,2 +41,411 @@ return zpts; } + +/* Slightly more generic loader than the ubercal version. Still assumes Nfilters x Nseasons + the input file must contain the following: + PHU Header : NFILTER, NSEASON, NCHIP_X, NCHIP_Y, NCELL_X, NCELL_Y, TS0_nnnn (season nnnn start mjd), TS1_nnnn (season nnnn end mjd) + + NSEASON * 2 extensions with + TABLE : mjd, zpt; Header: FILTER + IMAGE : 1D array of flat offsets +*/ + +ZptTable *load_zpt_ubercal(char *filename, int *nzpts, FlatCorrectionTable *flatcorrTable) { + + int i, nfilter, nseason, ix, iy, ixc, iyc, Ncol; + off_t Nrow; + char type[16], filter[80]; + int Nzpts, NZPTS; + ZptTable *zpts; + Header header; + Header theader; + Matrix matrix; + FTable ftable; + + // parameters describing the flat-field correction + int NSEASON; + int NFILTER; + int NCHIP_X; + int NCHIP_Y; + int NCELL_X; + int NCELL_Y; + int CHIP_DX; + int CHIP_DY; + + *nzpts = 0; + + FILE *f = fopen (filename, "r"); + if (!f) { + fprintf (stderr, "ERROR: cannot open zpt table file %s\n", filename); + exit (1); + } + + /* load in PHU segment (ignore) */ + if (!gfits_fread_header (f, &header)) { + if (VERBOSE) fprintf (stderr, "can't read Flat Correction header\n"); + fclose (f); + return (NULL); + } + if (!gfits_fread_matrix (f, &matrix, &header)) { + if (VERBOSE) fprintf (stderr, "can't read Flat Correction matrix\n"); + gfits_free_header (&header); + fclose (f); + return (NULL); + } + + Nzpts = 0; + NZPTS = 0; + ALLOCATE (zpts, ZptTable, NZPTS); + + // this function would be better if we read the list of filters, seasons, and the dimensions from the header + // for current testing, make fake smfs that correspond to specific chips, filter, and mjd ranges? + + // the simple files from Eddie have no internal metadata describing the corrections, + // so they must be manually encoded + + // hard-wired values which describe the ubercal analysis + if (NO_METADATA) { + NFILTER = 5; + NSEASON = 4; + NCHIP_X = 8; + NCHIP_Y = 8; + NCELL_X = 2; + NCELL_Y = 2; + CHIP_DX = 2*2424; // != 4880 + CHIP_DY = 2*2430; // != 4864 + // Note that Eddie has identified the center of the boundary between cells xy3n and + // xy4n and used that to split the chip. I am setting CHIP_DX,DY to force the same split + } else { + if (!gfits_scan (&header, "NSEASON", "%d", 1, &NSEASON)) { + fprintf (stderr, "cannot find NSEASON in header of %s\n", filename); + fclose (f); + return NULL; + } + if (!gfits_scan (&header, "NFILTER", "%d", 1, &NFILTER)) { + fprintf (stderr, "cannot find NFILTER in header of %s\n", filename); + fclose (f); + return NULL; + } + if (!gfits_scan (&header, "NCHIP_X", "%d", 1, &NCHIP_X)) { + fprintf (stderr, "cannot find NCHIP_X in header of %s\n", filename); + fclose (f); + return NULL; + } + if (!gfits_scan (&header, "NCHIP_Y", "%d", 1, &NCHIP_Y)) { + fprintf (stderr, "cannot find NCHIP_Y in header of %s\n", filename); + fclose (f); + return NULL; + } + if (!gfits_scan (&header, "NCELL_X", "%d", 1, &NCELL_X)) { + fprintf (stderr, "cannot find NCELL_X in header of %s\n", filename); + fclose (f); + return NULL; + } + if (!gfits_scan (&header, "NCELL_Y", "%d", 1, &NCELL_Y)) { + fprintf (stderr, "cannot find NCELL_Y in header of %s\n", filename); + fclose (f); + return NULL; + } + if (!gfits_scan (&header, "CHIP_DX", "%d", 1, &CHIP_DX)) { + fprintf (stderr, "cannot find CHIP_DX in header of %s\n", filename); + fclose (f); + return NULL; + } + if (!gfits_scan (&header, "CHIP_DY", "%d", 1, &CHIP_DY)) { + fprintf (stderr, "cannot find CHIP_DY in header of %s\n", filename); + fclose (f); + return NULL; + } + } + + flatcorrTable->Nseason = NSEASON; + ALLOCATE (flatcorrTable->tstart, e_time, NSEASON); + ALLOCATE (flatcorrTable->tstop, e_time, NSEASON); + + // ubercal hard-coded values + char filters_uc[5][3] = {"g", "r", "i", "z", "y"}; + double tstart_uc[] = {55000.0, 55296.0, 55327.0, 55662.0}; + double tstop_uc[] = {55296.0, 55327.0, 55662.0, 60000.0}; + + double mjdstart, mjdstop; + for (i = 0; i < NSEASON; i++) { + if (NO_METADATA) { + mjdstart = tstart_uc[i]; + mjdstop = tstop_uc[i]; + } else { + char name[9]; + snprintf (name, 9, "TS0_%04d", i); + if (!gfits_scan (&header, name, "%lf", 1, &mjdstart)) { + fprintf (stderr, "cannot find %s in header of %s\n", name, filename); + fclose (f); + return NULL; + } + + snprintf (name, 9, "TS1_%04d", i); + if (!gfits_scan (&header, name, "%lf", 1, &mjdstop)) { + fprintf (stderr, "cannot find %s in header of %s\n", name, filename); + fclose (f); + return NULL; + } + } + flatcorrTable->tstart[i] = ohana_mjd_to_sec(mjdstart); + flatcorrTable->tstop[i] = ohana_mjd_to_sec(mjdstop); + } + + // we have 5 filters, and 4 flat-field correction sets for each + flatcorrTable->Ncorr = NFILTER*NSEASON*NCHIP_X*NCHIP_Y*NCELL_X*NCELL_Y; + flatcorrTable->Nimage = NFILTER*NSEASON*NCHIP_X*NCHIP_Y; + + ALLOCATE (flatcorrTable->corr, FlatCorrection, flatcorrTable->Ncorr); + ALLOCATE (flatcorrTable->image, FlatCorrectionImage, flatcorrTable->Nimage); + memset (flatcorrTable->corr, 0, flatcorrTable->Ncorr*sizeof(FlatCorrection)); + + int corrID = 1; + int Nimage = 0; + ftable.header = &theader; + for (nfilter = 0; nfilter < NFILTER; nfilter++) { + // *** load the ZERO POINT table *** + + // load data for this header + if (!gfits_load_header (f, &theader)) return (NULL); + + // read the fits table bytes + if (!gfits_fread_ftable_data (f, &ftable, FALSE)) return (NULL); + + // skip over remaining bytes in data segment + fseeko (f, ftable.datasize - ftable.validsize, SEEK_CUR); + + // need to create and assign to flat-field correction + double *mjd = gfits_get_bintable_column_data (&theader, &ftable, "mjd_obs", type, &Nrow, &Ncol); + assert (!strcmp(type, "double")); + + double *zp = gfits_get_bintable_column_data (&theader, &ftable, "zp", type, &Nrow, &Ncol); + assert (!strcmp(type, "double")); + + // float *zperr = gfits_get_bintable_column_data (&theader, &ftable, "resid", type, &Nrow, &Ncol); + // assert (!strcmp(type, "float")); + + NZPTS += Nrow; + REALLOCATE (zpts, ZptTable, NZPTS); + for (i = 0; i < Nrow; i++) { + zpts[i+Nzpts].time = ohana_mjd_to_sec(mjd[i]); + zpts[i+Nzpts].zpt = zp[i]; + zpts[i+Nzpts].zpt_err = 0.0; + zpts[i+Nzpts].found = FALSE; + } + Nzpts += Nrow; + + // *** load the flat-field correction image *** + + // the image contains the flat-field corrections for a specific filter + + // load data for this header + if (!gfits_load_header (f, &header)) return (NULL); + + if (NO_METADATA) { + strcpy (filter, filters_uc[nfilter]); + } else { + if (!gfits_scan (&header, "FILTER", "%s", 1, filter)) { + fprintf (stderr, "cannot find FILTER in header of %s\n", filename); + fclose (f); + return NULL; + } + } + + // read the fits table bytes + double *offset64 = NULL; + float *offset32 = NULL; + int use32 = TRUE; + + if (!gfits_fread_matrix (f, &matrix, &header)) return (NULL); + switch (header.bitpix) { + case -32: // float + offset32 = (float *) matrix.buffer; + use32 = TRUE; + break; + case -64: // double + offset64 = (double *) matrix.buffer; + use32 = FALSE; + break; + default: + fprintf (stderr, "invalid bitpix for flat-field correction image: %d\n", header.bitpix); + exit (2); + } + + // XXX the initial hacked-together table from Eddie is missing the last 2 elements. they should be zero + if (NO_METADATA) { + matrix.Naxis[0] += 2; + REALLOCATE (offset64, double, matrix.Naxis[0]); + matrix.buffer = (char *) offset64; + } + + for (nseason = 0; nseason < NSEASON; nseason++) { // seasons + for (iy = 0; iy < NCHIP_Y; iy++) { // y-chip + for (ix = 0; ix < NCHIP_X; ix++) { // x-chip + + // photcode name + char photname[64]; + snprintf (photname, 64, "GPC1.%s.XY%d%d", filter, ix, iy); + // note that the XY00, XY07, etc, chips will have photcode values of 0 + + flatcorrTable->image[Nimage].photcode = GetPhotcodeCodebyName(photname); + flatcorrTable->image[Nimage].Nx = NCELL_X; + flatcorrTable->image[Nimage].Ny = NCELL_Y; + flatcorrTable->image[Nimage].ID = corrID; + flatcorrTable->image[Nimage].DX = CHIP_DX; + flatcorrTable->image[Nimage].DY = CHIP_DY; + flatcorrTable->image[Nimage].tstart = flatcorrTable->tstart[nseason]; + flatcorrTable->image[Nimage].tstop = flatcorrTable->tstop[nseason]; + + int seq_full = -1; // sequence number within the full table (all filters concatenated together) + int seq_filt = -1; // sequence number for just this filter + + // XXX we should have a NCHIP_X * NCHIP_Y array of x and y parity values + + // This enforces a 180 chip rotation for XY3n - XY7n & is only known to be valid for GPC1 (the XYnn names as well) + for (iyc = 0; iyc < NCELL_Y; iyc++) { + for (ixc = 0; ixc < NCELL_X; ixc++) { + if (ix < 4) { + // chips ix < 4 should be flipped + seq_full = nfilter*NSEASON*NCHIP_X*NCHIP_Y*NCELL_X*NCELL_Y + nseason*NCHIP_X*NCHIP_Y*NCELL_X*NCELL_Y + (iy*NCELL_Y + NCELL_Y - 1 - iyc)*NCHIP_X*NCELL_X + (ix*NCELL_X + NCELL_X - 1 - ixc); + seq_filt = nseason*NCHIP_X*NCHIP_Y*NCELL_X*NCELL_Y + (iy*NCELL_Y + NCELL_Y - 1 - iyc)*NCHIP_X*NCELL_X + (ix*NCELL_X + NCELL_X - 1 - ixc); + // = .... + y_parity*iyc - (y_parity - 1)*(NCELL_Y - 1) / 2 + } else { + seq_full = nfilter*NSEASON*NCHIP_X*NCHIP_Y*NCELL_X*NCELL_Y + nseason*NCHIP_X*NCHIP_Y*NCELL_X*NCELL_Y + (iy*NCELL_Y + iyc)*NCHIP_X*NCELL_X + (ix*NCELL_X + ixc); + seq_filt = nseason*NCHIP_X*NCHIP_Y*NCELL_X*NCELL_Y + (iy*NCELL_Y + iyc)*NCHIP_X*NCELL_X + (ix*NCELL_X + ixc); + } + assert (seq_full > -1); + assert (seq_filt > -1); + assert (seq_full < flatcorrTable->Ncorr); + if (matrix.Naxes == 1) { + assert (seq_filt < matrix.Naxis[0]); + } + if (matrix.Naxes == 2) { + assert (seq_filt < matrix.Naxis[0]*matrix.Naxis[1]); + } + assert (!flatcorrTable->corr[seq_full].ID); + flatcorrTable->corr[seq_full].x = ixc; + flatcorrTable->corr[seq_full].y = iyc; + if (use32) { + flatcorrTable->corr[seq_full].offset = offset32[seq_filt]; + } else { + flatcorrTable->corr[seq_full].offset = offset64[seq_filt]; + } + flatcorrTable->corr[seq_full].ID = corrID; + } + } + corrID ++; + Nimage ++; + } + } + } + } + + /*** convert from corr,image format to offsets ***/ + FlatCorrectionInternal (flatcorrTable); + + fprintf (stderr, "loaded %d zero points\n", Nzpts); + + *nzpts = Nzpts; + return zpts; +} + +static short maxCode = 0; +static short *zpt_codeval = NULL; +static float *zpt_offsets = NULL; +static short *zpt_index = NULL; + +int parse_zpt_offsets (char *ZPT_OFFSET_FILTERS, char *ZPT_OFFSET_VALUES) { + + if (!ZPT_OFFSET_FILTERS) return TRUE; + assert (ZPT_OFFSET_FILTERS); + assert (ZPT_OFFSET_VALUES); + + char *filters = strcreate(ZPT_OFFSET_FILTERS); + char *values = strcreate(ZPT_OFFSET_VALUES); + + char *p1 = filters; + char *p2 = values; + + char *filter = NULL; + char *value = NULL; + + char *s1 = NULL; + char *s2 = NULL; + + // save an array of the zero point offsets and their photcodes + int Ncode = 0; + int NCODE = 8; + ALLOCATE (zpt_codeval, short, NCODE); + ALLOCATE (zpt_offsets, float, NCODE); + + while (1) { + filter = strtok_r (p1, ",", &s1); + if (!filter) break; + + value = strtok_r (p2, ",", &s2); + if (!value) { + fprintf (stderr, "ERROR: mismatch between list of photcodes and list of offsets\n"); + exit (1); + } + + // do something here + PhotCode *code = GetPhotcodebyName (filter); + if (!code) { + fprintf (stderr, "ERROR: unknown photcode %s\n", filter); + exit (1); + } + if (code->type != PHOT_SEC) { + fprintf (stderr, "ERROR: photcode %s is not an average (SEC) type\n", filter); + exit (1); + } + + zpt_codeval[Ncode] = code->code; + zpt_offsets[Ncode] = atof(value); + Ncode ++; + + if (Ncode >= NCODE) { + NCODE += 8; + REALLOCATE (zpt_codeval, short, NCODE); + REALLOCATE (zpt_offsets, float, NCODE); + } + + maxCode = MAX(code->code, maxCode); + + p1 = NULL; + p2 = NULL; + } + + // generate an index to the zpt offsets by the photcode + ALLOCATE (zpt_index, short, maxCode + 1); + + int i; + for (i = 0; i < maxCode + 1; i++) { + zpt_index[i] = -1; + } + + for (i = 0; i < Ncode; i++) { + short seq = zpt_codeval[i]; + if (zpt_index[seq] != -1) { + fprintf (stderr, "duplicate photcode in -zpt-offset list\n"); + exit (2); + } + zpt_index[seq] = i; + } + return TRUE; +} + +float apply_zpt_offset (short code) { + + if (!zpt_offsets) return 0.0; + + // code is a primary photcode (but may be out of range, in which case no correction is supplied) + if (code <= 0) return 0.0; + if (code > maxCode) return 0.0; + + short index = zpt_index[code]; + if (index == -1) return 0.0; + + float offset = zpt_offsets[index]; + return offset; +} Index: trunk/Ohana/src/uniphot/src/match_flatcorr_to_images.c =================================================================== --- trunk/Ohana/src/uniphot/src/match_flatcorr_to_images.c (revision 33654) +++ trunk/Ohana/src/uniphot/src/match_flatcorr_to_images.c (revision 33654) @@ -0,0 +1,62 @@ +# include "setphot.h" + +// the date/time of the image is used to find the 'season' +// the photcode is used to find the actual flat correction + +int match_flatcorr_to_images (Image *image, off_t Nimage, FlatCorrectionTable *flatcorrTable) { + + int i, j; + + // make 4 lookup tables for the photcodes (one for each season) + + short **index; + + short maxCode = 0; + for (i = 0; i < flatcorrTable->Nimage; i++) { + maxCode = MAX (maxCode, flatcorrTable->image[i].photcode); + } + maxCode ++; // we want the outer bound, not the last value + + ALLOCATE (index, short *, flatcorrTable->Nseason); + for (i = 0; i < flatcorrTable->Nseason; i++) { + ALLOCATE (index[i], short, maxCode); + for (j = 0; j < maxCode; j++) { + index[i][j] = -1; + } + } + + for (i = 0; i < flatcorrTable->Nimage; i++) { + // the flat-field correction is defined for non-existent images (eg, XY00) + // these have photcode == 0, so skip them + if (!flatcorrTable->image[i].photcode) continue; + + // which season? + for (j = 0; j < flatcorrTable->Nseason; j++) { + if (flatcorrTable->image[i].tstart == flatcorrTable->tstart[j]) { + assert (index[j][flatcorrTable->image[i].photcode] == -1); + index[j][flatcorrTable->image[i].photcode] = i; + break; + } + } + } + + for (i = 0; i < Nimage; i++) { + if (!image[i].photcode) continue; // skip PHU images + + int found = FALSE; + for (j = 0; !found && (j < flatcorrTable->Nseason); j++) { + if (image[i].tzero < flatcorrTable->tstart[j]) continue; + if (image[i].tzero > flatcorrTable->tstop[j]) continue; + + int seq = index[j][image[i].photcode]; + if (seq == -1) break; + + image[i].photom_map_id = flatcorrTable->image[seq].ID; + found = TRUE; + } + if (!found) { + fprintf (stderr, "image does not match flatcorr Table ranges: %s\n", image[i].name); + } + } + return TRUE; +} Index: trunk/Ohana/src/uniphot/src/match_zpts_to_images.c =================================================================== --- trunk/Ohana/src/uniphot/src/match_zpts_to_images.c (revision 33653) +++ trunk/Ohana/src/uniphot/src/match_zpts_to_images.c (revision 33654) @@ -1,3 +1,3 @@ -# include "uniphot.h" +# include "setphot.h" /* sort a coordinate pair (X,Y) and the associated index (S) */ @@ -39,9 +39,25 @@ sort_zpts_by_time (zpts, Nzpts); + if (RESET) { + for (i = 0; i < Nimage; i++) { + image[i].Mcal = 0.0; + image[i].dMcal = NAN; + image[i].flags &= ~ID_IMAGE_PHOTOM_NOCAL; // clear the NOCAL flag + if (UBERCAL) { + image[i].flags &= ~ID_IMAGE_PHOTOM_UBERCAL; // clear the UBERCAL flag + } + } + } + NImatch = 0; // matched images for (i = Nz = 0; (i < Nimage) && (Nz < Nzpts); ) { - if (i % 1000 == 0) fprintf (stderr, "."); - if (Nz % 100 == 0) fprintf (stderr, "!"); + // if (i % 1000 == 0) fprintf (stderr, "."); + // if (Nz % 100 == 0) fprintf (stderr, "!"); + + if (!isfinite(zpts[Nz].zpt)) { + Nz++; + continue; + } Ni = index[i]; @@ -69,7 +85,24 @@ // we have a match: set zpt and record the match - image[Ni].Mcal = SCALE*code[0].C - zpts[Nz].zpt; + // is the zero point supplied nominally corrected or not? + // UBERCAL includes 2.5log(exptime) + K*airmass in the zero point + if (UBERCAL) { + image[Ni].Mcal = SCALE*code[0].C - zpts[Nz].zpt + 2.5*log10(image[Ni].exptime) + code[0].K*(image[Ni].secz - 1.000); + myAssert (isfinite(image[Ni].Mcal), "oops, ubercal made a nan image"); + } else { + image[Ni].Mcal = SCALE*code[0].C - zpts[Nz].zpt; + } + + // if we have defined zero point offsets, then apply them here + float offset = apply_zpt_offset (code[0].equiv); + assert (isfinite(offset)); + image[Ni].Mcal += offset; + image[Ni].dMcal = zpts[Nz].zpt_err; image[Ni].flags &= ~ID_IMAGE_PHOTOM_NOCAL; // clear the NOCAL flag + // image[Ni].flags |= ID_IMAGE_PHOTOM_EXTERN; XXX do we want some flag like this? + if (UBERCAL) { + image[Ni].flags |= ID_IMAGE_PHOTOM_UBERCAL; + } zpts[Nz].found = TRUE; NImatch ++; Index: trunk/Ohana/src/uniphot/src/setfwhm.c =================================================================== --- trunk/Ohana/src/uniphot/src/setfwhm.c (revision 33653) +++ trunk/Ohana/src/uniphot/src/setfwhm.c (revision 33654) @@ -20,5 +20,5 @@ // load images - image = load_images_setphot (&db, &Nimage); + image = load_images_setfwhm (&db, &Nimage); if (!UPDATE) dvo_image_unlock (&db); Index: trunk/Ohana/src/uniphot/src/setphot.c =================================================================== --- trunk/Ohana/src/uniphot/src/setphot.c (revision 33653) +++ trunk/Ohana/src/uniphot/src/setphot.c (revision 33654) @@ -1,3 +1,3 @@ -# include "uniphot.h" +# include "setphot.h" int main (int argc, char **argv) { @@ -8,4 +8,5 @@ ZptTable *zpts; Image *image; + FlatCorrectionTable flatcorrTable; /* get configuration info, args, lockfile */ @@ -17,5 +18,17 @@ if (db.dbstate == LCK_EMPTY) Shutdown ("ERROR: No images in catalog %s (1)", db.filename); - zpts = load_zpt_table (argv[1], &Nzpts); + if (UBERCAL) { + zpts = load_zpt_ubercal (argv[1], &Nzpts, &flatcorrTable); + + char flatcorrfile[MAX_PATH_LENGTH]; + int size = snprintf (flatcorrfile, MAX_PATH_LENGTH, "%s/flatcorr.fits", CATDIR); + assert (size < MAX_PATH_LENGTH); + FlatCorrectionSave(&flatcorrTable, flatcorrfile); + // XXX should this program save any old copy of this file? + } else { + zpts = load_zpt_table (argv[1], &Nzpts); + } + + if (!zpts) Shutdown ("failed to load zero points, or empty table"); // load images @@ -25,11 +38,18 @@ match_zpts_to_images (image, Nimage, zpts, Nzpts); - update_setphot(image, Nimage); + if (UBERCAL) { + match_flatcorr_to_images (image, Nimage, &flatcorrTable); + } - // write image table + if (!IMAGES_ONLY) { + status = update_dvo_setphot (image, Nimage, &flatcorrTable); + } + + // write image table (even if some remote clients failed) if (UPDATE) { dvo_image_save (&db, VERBOSE); } + if (!status) exit (1); exit (0); } @@ -39,10 +59,23 @@ setphot (zpt_table) - * load text table of zpts, time (photcode?) - * load images - * match images to zpts - * set zpts - * load catalogs - * update detection (& averages?) + setphot has 2 modes : basic & ubercal + + * outline of basic mode: + + ** load text table of zpts, time (photcode?) + ** load images + ** match images to zpts + ** set zpts + ** load catalogs + ** update detection (& averages?) + + * outline of ubercal mode: + + ** load fits table of zpts & time + table of flat-field corrections + ** load images + ** match images to zpts + ** set zpts + ** load catalogs + ** update detection (& averages?) */ Index: trunk/Ohana/src/uniphot/src/setphot_client.c =================================================================== --- trunk/Ohana/src/uniphot/src/setphot_client.c (revision 33654) +++ trunk/Ohana/src/uniphot/src/setphot_client.c (revision 33654) @@ -0,0 +1,30 @@ +# include "setphot.h" + +// setphot_client is run on a remote host and is responsible for updating the catalogs +// owned by that host. + +// setphot_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_setphot_client (argc, argv); + + // load the image subset table from the specified location + off_t Nimage; + ImageSubset *subset = ImageSubsetLoad (IMAGES, &Nimage); + Image *image = ImagesFromSubset (subset, Nimage); + + // load the flat-field correction table from CATDIR + char flatcorrFile[1024]; + snprintf (flatcorrFile, 1024, "%s/flatcorr.fits", CATDIR); + FlatCorrectionTable *flatcorrTable = FlatCorrectionLoad (flatcorrFile, VERBOSE); + + update_dvo_setphot (image, Nimage, flatcorrTable); + + exit (0); +} Index: trunk/Ohana/src/uniphot/src/uniphot.c =================================================================== --- trunk/Ohana/src/uniphot/src/uniphot.c (revision 33653) +++ trunk/Ohana/src/uniphot/src/uniphot.c (revision 33654) @@ -74,5 +74,5 @@ if (!UPDATE) exit (0); - update (&db, sgroup, Nsgroup); + update_dvo_uniphot (&db, sgroup, Nsgroup); dvo_image_unlock (&db); Index: trunk/Ohana/src/uniphot/src/update.c =================================================================== --- trunk/Ohana/src/uniphot/src/update.c (revision 33653) +++ (revision ) @@ -1,141 +1,0 @@ -# include "uniphot.h" -# include - -void update (FITS_DB *db, Group *sgroup, int Nsgroup) { - - off_t i, Nimage, Nkeep, *keep; - int j, status, Nmin; - char line[256]; - glob_t pglob; - double Rmin, Rmax, Dmin, Dmax, x, y, radius; - Image *image; - Catalog catalog; - Coords coords; - - image = gfits_table_get_Image (&db[0].ftable, &Nimage, &db[0].swapped); - if (!image) { - fprintf (stderr, "ERROR: failed to read images\n"); - exit (2); - } - - // create a subset list so we can make a vtable - Nkeep = 0; - ALLOCATE (keep, off_t, Nimage); - - // identify the images used and clear the NOCAL flags on the rest - for (i = 0; i < Nimage; i++) { - if (image[i].flags & ID_IMAGE_PHOTOM_NOCAL) { - image[i].flags &= ~ID_IMAGE_PHOTOM_NOCAL; - continue; - } - keep[Nkeep] = i; - Nkeep ++; - } - - /* apply calculated space-group offset to image Mcal values */ - for (i = 0; i < Nsgroup; i++) { - for (j = 0; j < sgroup[i].Nimage; j++) { - sgroup[i].image[j][0].Mcal -= sgroup[i].M; - } - } - - // save the rows in the image table which were used in this analysis - gfits_vtable_from_ftable (&db[0].ftable, &db[0].vtable, keep, Nkeep); - - // write image table - dvo_image_update (db, VERBOSE); - - // XXX need to fix the update for the catalog (or make it optional) - return; - - // XXX this process uses the existence of the file to perform the update - // XXX convert this to an examination of the SkyTable - /** update catalog tables **/ - pglob.gl_offs = 0; - sprintf (line, "%s/*/*.cpt", CATDIR); - status = glob (line, 0, NULL, &pglob); - - coords.crpix1 = coords.crpix2 = 0.0; - coords.cdelt1 = coords.cdelt2 = 1.0; /* scale is degrees, radius in degrees */ - coords.pc1_1 = coords.pc2_2 = 1.0; - coords.pc1_2 = coords.pc2_1 = 0.0; - strcpy (coords.ctype, "RA---TAN"); - - for (i = 0; i < pglob.gl_pathc; i++) { - catalog.filename = pglob.gl_pathv[i]; - catalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF; - catalog.Nsecfilt = GetPhotcodeNsecfilt (); - - // an error exit status here is a significant error - if (!dvo_catalog_open (&catalog, NULL, VERBOSE, "a")) { - fprintf (stderr, "ERROR: failure to open/create catalog file %s\n", catalog.filename); - exit (2); - } - if (!catalog.Naves_disk) { - dvo_catalog_unlock (&catalog); - dvo_catalog_free (&catalog); - continue; - } - - gfits_scan (&catalog.header, "RA0", "%lf", 1, &Rmin); - gfits_scan (&catalog.header, "RA1", "%lf", 1, &Rmax); - gfits_scan (&catalog.header, "DEC0", "%lf", 1, &Dmin); - gfits_scan (&catalog.header, "DEC1", "%lf", 1, &Dmax); - - Rmin = ohana_normalize_angle (Rmin); - Rmax = ohana_normalize_angle (Rmax); - - coords.crval1 = 0.5*(Rmin + Rmax); - coords.crval2 = 0.5*(Dmin + Dmax); - - Nmin = 0; - Rmin = 1000; - /* primitive version: match catalog with closest sgroup */ - for (j = 0; j < Nsgroup; j++) { - if (!RD_to_XY (&x, &y, sgroup[j].v1, sgroup[j].v2, &coords)) continue; - radius = hypot (x, y); - if ((j == 0) || (radius < Rmin)) { - Rmin = radius; - Nmin = j; - } - } - - fprintf (stderr, "catalog: %s sgroup: %d %s %f\n", catalog.filename, Nmin, sgroup[Nmin].label, Rmin); - update_dvo_catalog (&catalog, &sgroup[Nmin], (Rmin > 2*RADIUS)); - dvo_catalog_save (&catalog, VERBOSE); - dvo_catalog_unlock (&catalog); - dvo_catalog_free (&catalog); - } -} - -/* loop over all average - - if (source equiv photcode) - - loop over sgroups - - if (radius < RADIUS) - - average.M (+/-) source.M -*/ - -/* loop over all measure - - if (source equiv photcode) - - loop over sgroups - - if (radius < RADIUS) - - measure.Mcal (+/-) source.M -*/ - -/* alternative (slower, more robust) - - loop over all average - loop over all average.Nm - select measure with source equiv photcode - match measure to image, get Mcal & offset - if any offset for average[i] is different, give error - apply offset to average.M -*/ Index: trunk/Ohana/src/uniphot/src/update_catalog.c =================================================================== --- trunk/Ohana/src/uniphot/src/update_catalog.c (revision 33653) +++ (revision ) @@ -1,34 +1,0 @@ -# include "uniphot.h" - -void update_dvo_catalog (Catalog *catalog, Group *sgroup, int warn) { - - int i, j, m, found; - int Nsec, Nsecfilt; - PhotCode *code; - - Nsec = GetPhotcodeNsec (photcode[0].code); - Nsecfilt = GetPhotcodeNsecfilt (); - - found = 0; - for (i = 0; i < catalog[0].Naverage; i++) { - - if (!isnan(catalog[0].secfilt[i*Nsecfilt+Nsec].M)) { - catalog[0].secfilt[i*Nsecfilt+Nsec].M += sgroup[0].M; - } - - m = catalog[0].average[i].measureOffset; - for (j = 0; j < catalog[0].average[i].Nmeasure; j++, m++) { - code = GetPhotcodebyCode (catalog[0].measure[m].photcode); - if (code == NULL) continue; - if (code[0].type != PHOT_DEP) continue; - if (code[0].equiv != photcode[0].code) continue; - catalog[0].measure[m].Mcal -= sgroup[0].M; - found ++; - } - } - - if (found) { - fprintf (stderr, "found %d matches\n", found); - if (warn) fprintf (stderr, "warning: updated values, perhaps out of range?\n"); - } -} Index: trunk/Ohana/src/uniphot/src/update_catalog_setphot.c =================================================================== --- trunk/Ohana/src/uniphot/src/update_catalog_setphot.c (revision 33653) +++ trunk/Ohana/src/uniphot/src/update_catalog_setphot.c (revision 33654) @@ -1,93 +1,44 @@ -# include "uniphot.h" -# include +# include "setphot.h" -int update_setphot (Image *image, off_t Nimage) { +void update_catalog_setphot (Catalog *catalog, Image *image, off_t *index, off_t Nimage, FlatCorrectionTable *flatcorr) { - SkyRegion UserPatch; - SkyTable *sky = NULL; - SkyList *skylist = NULL; - Catalog catalog; - off_t i, maxID, *index; + off_t i, j, found; - // create an index for the image IDs - maxID = 0; - for (i = 0; i < Nimage; i++) { - maxID = MAX(maxID, image[i].imageID); - } - ALLOCATE (index, off_t, maxID + 1); - for (i = 0; i < Nimage; i++) { - if (!image[i].imageID) continue; // images with ID == 0 are virtual - index[image[i].imageID] = i; - } - - - // load the current sky table (layout of all SkyRegions) - sky = SkyTableLoadOptimal (CATDIR, SKY_TABLE, GSCFILE, TRUE, SKY_DEPTH, VERBOSE); - SkyTableSetFilenames (sky, CATDIR, "cpt"); - - // 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); - - // update measurements for each populated catalog - for (i = 0; i < skylist[0].Nregions; i++) { - - // set up the basic catalog info - catalog.filename = skylist[0].filename[i]; - catalog.catformat = dvo_catalog_catformat (CATFORMAT); // set the default catformat from config data - catalog.catmode = dvo_catalog_catmode (CATMODE); // set the default catmode from config data - 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_setphot (&catalog, image, index, Nimage); - - if (!UPDATE) { - dvo_catalog_unlock (&catalog); - dvo_catalog_free (&catalog); - continue; - } - - 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); -} - -void update_catalog_setphot (Catalog *catalog, Image *image, off_t *index, off_t Nimage) { - - off_t i, j, m, id, idx, found; - float Mcal, dMcal; + // if we are resetting, reset all flags + DVOMeasureFlags photomFlags = ID_MEAS_POOR_PHOTOM | ID_MEAS_SKIP_PHOTOM | ID_MEAS_PHOTOM_UBERCAL; found = 0; for (i = 0; i < catalog[0].Naverage; i++) { - m = catalog[0].average[i].measureOffset; + off_t m = catalog[0].average[i].measureOffset; for (j = 0; j < catalog[0].average[i].Nmeasure; j++, m++) { - idx = catalog[0].measure[m].imageID; + off_t idx = catalog[0].measure[m].imageID; if (idx <= 0) continue; // detections with imageID == 0 do not have a valid image (eg, ref photcode) - id = index[idx]; - Mcal = image[id].Mcal; - dMcal = image[id].dMcal; + off_t id = index[idx]; + float Mcal = image[id].Mcal; + float dMcal = image[id].dMcal; + float Mcal_offset = 0.0; - catalog[0].measure[m].Mcal = Mcal; + // if we know about a flat-field correction, then we need to apply the sub-chip correction + int flat_id = image[id].photom_map_id; + if (flat_id > 0) { + Mcal_offset = FlatCorrectionOffset (flatcorr, flat_id, catalog[0].measure[m].Xccd, catalog[0].measure[m].Yccd); + } + + catalog[0].measure[m].Mcal = Mcal - Mcal_offset; catalog[0].measure[m].dMcal = dMcal; + myAssert(isfinite(catalog[0].measure[m].Mcal), "oops: ubercal made a nan"); + + if (RESET) { + catalog[0].measure[m].dbFlags &= ~photomFlags; + } + + // if we are setting the zero points from an UBERCAL database, and this detection is from one of those images, + // then tag the measurement as well + if (UBERCAL && (image[id].flags & ID_IMAGE_PHOTOM_UBERCAL)) { + catalog[0].measure[m].dbFlags |= ID_MEAS_PHOTOM_UBERCAL; + } + found ++; } Index: trunk/Ohana/src/uniphot/src/update_catalog_uniphot.c =================================================================== --- trunk/Ohana/src/uniphot/src/update_catalog_uniphot.c (revision 33654) +++ trunk/Ohana/src/uniphot/src/update_catalog_uniphot.c (revision 33654) @@ -0,0 +1,35 @@ +# include "uniphot.h" + +// update a single catalog for uniphot +void update_catalog_uniphot (Catalog *catalog, Group *sgroup, int warn) { + + int i, j, m, found; + int Nsec, Nsecfilt; + PhotCode *code; + + Nsec = GetPhotcodeNsec (photcode[0].code); + Nsecfilt = GetPhotcodeNsecfilt (); + + found = 0; + for (i = 0; i < catalog[0].Naverage; i++) { + + if (!isnan(catalog[0].secfilt[i*Nsecfilt+Nsec].M)) { + catalog[0].secfilt[i*Nsecfilt+Nsec].M += sgroup[0].M; + } + + m = catalog[0].average[i].measureOffset; + for (j = 0; j < catalog[0].average[i].Nmeasure; j++, m++) { + code = GetPhotcodebyCode (catalog[0].measure[m].photcode); + if (code == NULL) continue; + if (code[0].type != PHOT_DEP) continue; + if (code[0].equiv != photcode[0].code) continue; + catalog[0].measure[m].Mcal -= sgroup[0].M; + found ++; + } + } + + if (found) { + fprintf (stderr, "found %d matches\n", found); + if (warn) fprintf (stderr, "warning: updated values, perhaps out of range?\n"); + } +} Index: trunk/Ohana/src/uniphot/src/update_dvo_setphot.c =================================================================== --- trunk/Ohana/src/uniphot/src/update_dvo_setphot.c (revision 33654) +++ trunk/Ohana/src/uniphot/src/update_dvo_setphot.c (revision 33654) @@ -0,0 +1,155 @@ +# include "setphot.h" +int update_dvo_setphot_parallel (SkyTable *sky, Image *image, off_t Nimage); + +int update_dvo_setphot (Image *image, off_t Nimage, FlatCorrectionTable *flatcorr) { + + SkyRegion UserPatch; + SkyTable *sky = NULL; + SkyList *skylist = NULL; + Catalog catalog; + off_t i, maxID, *index; + + // 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) { + update_dvo_setphot_parallel (sky, image, Nimage); + return TRUE; + } + + // create an index for the image IDs + maxID = 0; + for (i = 0; i < Nimage; i++) { + maxID = MAX(maxID, image[i].imageID); + } + ALLOCATE (index, off_t, maxID + 1); + for (i = 0; i < maxID + 1; i++) { + index[i] = -1; + } + for (i = 0; i < Nimage; i++) { + if (!image[i].imageID) continue; // images with ID == 0 are virtual + index[image[i].imageID] = i; + } + + // 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); + + // 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_setphot (&catalog, image, index, Nimage, flatcorr); + + if (!UPDATE) { + dvo_catalog_unlock (&catalog); + dvo_catalog_free (&catalog); + continue; + } + + 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_setphot_parallel (SkyTable *sky, Image *image, off_t Nimage) { + + // write out the subset table of image information + char imageFile[512]; + snprintf (imageFile, 512, "%s/Images.subset.dat", CATDIR); + if (!ImageSubsetSave (imageFile, image, Nimage)) { + fprintf (stderr, "failed to write image subset\n"); + exit (1); + } + + // now launch the setphot_client jobs to the parallel hosts + + // 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, MAX_PATH_LENGTH); + free (table->hosts[i].pathname); + table->hosts[i].pathname = tmppath; + + 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); + + fprintf (stderr, "command: %s\n", command); + + char tmpline[1024]; + if (VERBOSE) { snprintf (tmpline, 1024, "%s -v", command); strcpy (command, tmpline); } + if (RESET) { snprintf (tmpline, 1024, "%s -reset", command); strcpy (command, tmpline); } + if (UPDATE) { snprintf (tmpline, 1024, "%s -update", command); strcpy (command, tmpline); } + if (UBERCAL) { snprintf (tmpline, 1024, "%s -ubercal", command); strcpy (command, tmpline); } + + if (PARALLEL_MANUAL) continue; + + if (PARALLEL_SERIAL) { + int status = system (command); + if (status) { + fprintf (stderr, "ERROR running setphot_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 + } + } + + // wait for the remote jobs to be completed + if (PARALLEL_MANUAL) { + fprintf (stderr, "run the setphot_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_uniphot.c =================================================================== --- trunk/Ohana/src/uniphot/src/update_dvo_uniphot.c (revision 33654) +++ trunk/Ohana/src/uniphot/src/update_dvo_uniphot.c (revision 33654) @@ -0,0 +1,142 @@ +# include "uniphot.h" +# include + +// update images & catalogs for uniphot +void update_dvo_uniphot (FITS_DB *db, Group *sgroup, int Nsgroup) { + + off_t i, Nimage, Nkeep, *keep; + int j, status, Nmin; + char line[256]; + glob_t pglob; + double Rmin, Rmax, Dmin, Dmax, x, y, radius; + Image *image; + Catalog catalog; + Coords coords; + + image = gfits_table_get_Image (&db[0].ftable, &Nimage, &db[0].swapped); + if (!image) { + fprintf (stderr, "ERROR: failed to read images\n"); + exit (2); + } + + // create a subset list so we can make a vtable + Nkeep = 0; + ALLOCATE (keep, off_t, Nimage); + + // identify the images used and clear the NOCAL flags on the rest + for (i = 0; i < Nimage; i++) { + if (image[i].flags & ID_IMAGE_PHOTOM_NOCAL) { + image[i].flags &= ~ID_IMAGE_PHOTOM_NOCAL; + continue; + } + keep[Nkeep] = i; + Nkeep ++; + } + + /* apply calculated space-group offset to image Mcal values */ + for (i = 0; i < Nsgroup; i++) { + for (j = 0; j < sgroup[i].Nimage; j++) { + sgroup[i].image[j][0].Mcal -= sgroup[i].M; + } + } + + // save the rows in the image table which were used in this analysis + gfits_vtable_from_ftable (&db[0].ftable, &db[0].vtable, keep, Nkeep); + + // write image table + dvo_image_update (db, VERBOSE); + + // XXX need to fix the update for the catalog (or make it optional) + return; + + // XXX this process uses the existence of the file to perform the update + // XXX convert this to an examination of the SkyTable + /** update catalog tables **/ + pglob.gl_offs = 0; + sprintf (line, "%s/*/*.cpt", CATDIR); + status = glob (line, 0, NULL, &pglob); + + coords.crpix1 = coords.crpix2 = 0.0; + coords.cdelt1 = coords.cdelt2 = 1.0; /* scale is degrees, radius in degrees */ + coords.pc1_1 = coords.pc2_2 = 1.0; + coords.pc1_2 = coords.pc2_1 = 0.0; + strcpy (coords.ctype, "RA---TAN"); + + for (i = 0; i < pglob.gl_pathc; i++) { + catalog.filename = pglob.gl_pathv[i]; + catalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF; + catalog.Nsecfilt = GetPhotcodeNsecfilt (); + + // an error exit status here is a significant error + if (!dvo_catalog_open (&catalog, NULL, VERBOSE, "a")) { + fprintf (stderr, "ERROR: failure to open/create catalog file %s\n", catalog.filename); + exit (2); + } + if (!catalog.Naves_disk) { + dvo_catalog_unlock (&catalog); + dvo_catalog_free (&catalog); + continue; + } + + gfits_scan (&catalog.header, "RA0", "%lf", 1, &Rmin); + gfits_scan (&catalog.header, "RA1", "%lf", 1, &Rmax); + gfits_scan (&catalog.header, "DEC0", "%lf", 1, &Dmin); + gfits_scan (&catalog.header, "DEC1", "%lf", 1, &Dmax); + + Rmin = ohana_normalize_angle (Rmin); + Rmax = ohana_normalize_angle (Rmax); + + coords.crval1 = 0.5*(Rmin + Rmax); + coords.crval2 = 0.5*(Dmin + Dmax); + + Nmin = 0; + Rmin = 1000; + /* primitive version: match catalog with closest sgroup */ + for (j = 0; j < Nsgroup; j++) { + if (!RD_to_XY (&x, &y, sgroup[j].v1, sgroup[j].v2, &coords)) continue; + radius = hypot (x, y); + if ((j == 0) || (radius < Rmin)) { + Rmin = radius; + Nmin = j; + } + } + + fprintf (stderr, "catalog: %s sgroup: %d %s %f\n", catalog.filename, Nmin, sgroup[Nmin].label, Rmin); + update_catalog_uniphot (&catalog, &sgroup[Nmin], (Rmin > 2*RADIUS)); + dvo_catalog_save (&catalog, VERBOSE); + dvo_catalog_unlock (&catalog); + dvo_catalog_free (&catalog); + } +} + +/* loop over all average + + if (source equiv photcode) + + loop over sgroups + + if (radius < RADIUS) + + average.M (+/-) source.M +*/ + +/* loop over all measure + + if (source equiv photcode) + + loop over sgroups + + if (radius < RADIUS) + + measure.Mcal (+/-) source.M +*/ + +/* alternative (slower, more robust) + + loop over all average + loop over all average.Nm + select measure with source equiv photcode + match measure to image, get Mcal & offset + if any offset for average[i] is different, give error + apply offset to average.M +*/