Index: /branches/eam_branches/ipp-20140813/Ohana/src/addstar/Makefile =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/addstar/Makefile (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/addstar/Makefile (revision 37344) @@ -86,4 +86,5 @@ $(SRC)/ReadStarsTEXT.$(ARCH).o \ $(SRC)/ReadStarsSDSS.$(ARCH).o \ +$(SRC)/ReadXradFITS.$(ARCH).o \ $(SRC)/FilterStars.$(ARCH).o \ $(SRC)/ImageOptions.$(ARCH).o \ Index: /branches/eam_branches/ipp-20140813/Ohana/src/addstar/include/addstar.h =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/addstar/include/addstar.h (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/addstar/include/addstar.h (revision 37344) @@ -35,7 +35,9 @@ char *exthead; char *extdata; + char *extxrad; char *exttype; int extnum_head; int extnum_data; + int extnum_xrad; } HeaderSet; @@ -146,4 +148,5 @@ int OLD_RESORT; +int READ_XRAD_DATA; int PARALLEL; @@ -250,4 +253,7 @@ Stars *FilterStars PROTO((Stars *instars, Image *image, unsigned int imageID, const AddstarClientOptions *options)); Stars *MergeStars PROTO((Stars *stars, unsigned int *Nstars, Stars *instars, unsigned int Ninstars)); + +int ReadXradFITS PROTO((FILE *f, Header *theader, Stars *stars, unsigned int Nstars)); + double scat_subpix PROTO((double x, double y)); void update_coords PROTO((Average *average, Measure *measure, off_t *next)); Index: /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/FilterStars.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/FilterStars.c (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/FilterStars.c (revision 37344) @@ -9,5 +9,5 @@ // the imageID supplied here is the sequence **within this set** -// this value is updated based on the image table later +// this value is updated based on the image table later (in UpdateImageIDs) Stars *FilterStars (Stars *instars, Image *image, unsigned int imageID, const AddstarClientOptions *options) { @@ -127,4 +127,9 @@ stars[N].measure.imageID = imageID; // this value is updated in UpdateImageIDs + // add imageID to lensing entry, if it exists + if (stars[N].lensing) { + stars[N].lensing->imageID = imageID; + } + N ++; } Index: /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/LoadData.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/LoadData.c (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/LoadData.c (revision 37344) @@ -83,4 +83,19 @@ continue; } + + // XRAD : if we want to read the xrad table, skip to that table here: + if (headerSets[i].extnum_xrad != -1) { + int Nxrad = headerSets[i].extnum_xrad; + Nskip = 0; + for (j = 0; j < Nxrad; j++) { + Nskip += extsize[j]; + } + fseeko (f, Nskip, SEEK_SET); + + if (!ReadXradFITS (f, headers[Nxrad], inStars, images[0][Nvalid].nstar)) { + fprintf (stderr, "problem reading the radial flux data for %s\n", headerSets[i].extdata); + } + } + inStars = FilterStars (inStars, &images[0][Nvalid], Nvalid, options); *stars = MergeStars (*stars, Nstars, inStars, images[0][Nvalid].nstar); Index: /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/MatchHeaders.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/MatchHeaders.c (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/MatchHeaders.c (revision 37344) @@ -79,13 +79,33 @@ headerSets[Nimage].extdata = strcreate (extname); headerSets[Nimage].exthead = strcreate (exthead); + headerSets[Nimage].extxrad = NULL; headerSets[Nimage].extnum_data = i; headerSets[Nimage].extnum_head = -1; + headerSets[Nimage].extnum_xrad = -1; + + // XXX a special case for fforce xrad data + if (READ_XRAD_DATA) { + // extname is foobar.psf, convert to foobar.xrad + int Nchar = strlen (extname); // foobar.psf : Nchar = 10 + ALLOCATE (headerSets[Nimage].extxrad, char, Nchar + 2); // xrad is 1 longer than psf + memcpy (headerSets[Nimage].extxrad, extname, Nchar - 3); // Nchar - 3 = 7 + memcpy (&headerSets[Nimage].extxrad[Nchar-3], "xrad", 4); + headerSets[Nimage].extxrad[Nchar + 1] = 0; // put the 0 at element 11 for foobar.xrad\0 + } // find the matching exthead entry for (j = 0; j < Nheaders; j++) { if (!gfits_scan (headers[j], ExtnameKeyword, "%s", 1, extname)) continue; - if (strcmp (extname, headerSets[Nimage].exthead)) continue; - headerSets[Nimage].extnum_head = j; - break; + if (!strcmp (extname, headerSets[Nimage].exthead)) { + headerSets[Nimage].extnum_head = j; + if (!READ_XRAD_DATA) break; + } + if (READ_XRAD_DATA && !strcmp (extname, headerSets[Nimage].extxrad)) { + headerSets[Nimage].extnum_xrad = j; + } + if ((headerSets[Nimage].extnum_head > -1) && (headerSets[Nimage].extnum_xrad > -1)) { + // we can only get here if READ_XRAD_DATA is true + break; + } } Index: /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/ReadXradFITS.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/ReadXradFITS.c (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/ReadXradFITS.c (revision 37344) @@ -0,0 +1,88 @@ +# include "addstar.h" + +// given a file with the pointer at the start of the table block and the +// corresponding image header, load the xrad data from the table +int ReadXradFITS (FILE *f, Header *theader, Stars *stars, unsigned int Nstars) { + + FTable table; + + // we've already read in the header + table.header = theader; + off_t Nskip = theader[0].datasize; + fseeko (f, Nskip, SEEK_CUR); + + /* load the table data */ + if (!gfits_fread_ftable_data (f, &table, FALSE)) { + fprintf (stderr, "ERROR: can't read table header\n"); + exit (1); + } + + // I want to read the following columns from this table: + + off_t Nrow, NrowAlt; + int Ncol, NcolAlt; + char type[16], name[80]; + + strcpy (name, "APER_FLUX"); + float *AperFlux = gfits_get_bintable_column_data (theader, &table, name, type, &Nrow, &Ncol); + myAssert (!strcmp(type, "float"), "wrong column type for %s\n", name); + + strcpy (name, "APER_FLUX_ERR"); + float *AperFluxErr = gfits_get_bintable_column_data (theader, &table, name, type, &NrowAlt, &NcolAlt); + myAssert (!strcmp(type, "float"), "wrong column type for %s\n", name); + myAssert (Nrow == NrowAlt, "column mismatch?"); + myAssert (Ncol == NcolAlt, "column mismatch?"); + + strcpy (name, "APER_FLUX_STDEV"); + float *AperFluxStd = gfits_get_bintable_column_data (theader, &table, name, type, &NrowAlt, &NcolAlt); + myAssert (!strcmp(type, "float"), "wrong column type for %s\n", name); + myAssert (Nrow == NrowAlt, "column mismatch?"); + myAssert (Ncol == NcolAlt, "column mismatch?"); + + strcpy (name, "APER_FILL"); + float *AperFill = gfits_get_bintable_column_data (theader, &table, name, type, &NrowAlt, &NcolAlt); + myAssert (!strcmp(type, "float"), "wrong column type for %s\n", name); + myAssert (Nrow == NrowAlt, "column mismatch?"); + myAssert (Ncol == NcolAlt, "column mismatch?"); + + strcpy (name, "IPP_IDET"); + int *RadID = gfits_get_bintable_column_data (theader, &table, name, type, &NrowAlt, &NcolAlt); + myAssert (!strcmp(type, "int"), "wrong column type for %s\n", name); + myAssert (Nrow == NrowAlt, "column mismatch?"); + + if (Nrow > Nstars) { + myAbort("more radial measurements than stars? seems like a bug\n"); + } + + int i; + int Nap = 0; + for (i = 0; i < Nstars; i++) { + if (stars[i].measure.detID < RadID[Nap]) { + continue; // this star does not have an aperture + } + if (stars[i].measure.detID > RadID[Nap]) { + myAbort("radial apertures out of order? seems like a bug\n"); + } + + // we could allocate here, or just insist this is an error? + if (!stars[i].lensing) { + fprintf (stderr, "!"); + ALLOCATE (stars[i].lensing, Lensing, 1); + dvo_lensing_init (stars[i].lensing); + } + + stars[i].lensing-> F_ApR5 = AperFlux [Nap*Ncol + 3]; + stars[i].lensing->dF_ApR5 = AperFluxErr[Nap*Ncol + 3]; + stars[i].lensing->sF_ApR5 = AperFluxStd[Nap*Ncol + 3]; + stars[i].lensing->fF_ApR5 = AperFill [Nap*Ncol + 3]; + + stars[i].lensing-> F_ApR6 = AperFlux [Nap*Ncol + 4]; + stars[i].lensing->dF_ApR6 = AperFluxErr[Nap*Ncol + 4]; + stars[i].lensing->sF_ApR6 = AperFluxStd[Nap*Ncol + 4]; + stars[i].lensing->fF_ApR6 = AperFill [Nap*Ncol + 4]; + Nap ++; + } + + return TRUE; +} + Index: /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/UpdateImageIDs.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/UpdateImageIDs.c (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/UpdateImageIDs.c (revision 37344) @@ -55,4 +55,7 @@ for (i = 0; i < Nstars; i++) { stars[i].measure.imageID += imageID; + if (stars[i].lensing) { + stars[i].lensing->imageID += imageID; + } } Index: /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/args.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/args.c (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/addstar/src/args.c (revision 37344) @@ -81,4 +81,10 @@ PMM_CCD_TABLE = strcreate (argv[N]); remove_argument (N, &argc, argv); + } + + READ_XRAD_DATA = FALSE; + if ((N = get_argument (argc, argv, "-xrad"))) { + remove_argument (N, &argc, argv); + READ_XRAD_DATA = TRUE; } Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/Makefile =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/Makefile (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/Makefile (revision 37344) @@ -0,0 +1,52 @@ +default: dvolens dvolens_client +help: + @echo "make options: dvolens (default)" + +include ../../Makefile.System +HOME = $(ROOT)/src/dvolens +BIN = $(HOME)/bin +LIB = $(HOME)/lib +SRC = $(HOME)/src +MAN = $(HOME)/doc +INC = $(HOME)/include +include ../../Makefile.Common + +# programs may add their own internal requirements here +FULL_CFLAGS = $(BASE_CFLAGS) +FULL_CPPFLAGS = $(BASE_CPPFLAGS) +FULL_LDFLAGS = -lkapa -ldvo -lFITS -lohana $(BASE_LDFLAGS) + +dvolens: $(BIN)/dvolens.$(ARCH) +dvolens_client: $(BIN)/dvolens_client.$(ARCH) + +install: $(DESTBIN)/dvolens $(DESTBIN)/dvolens_client + +DVOLENS = \ +$(SRC)/ConfigInit.$(ARCH).o \ +$(SRC)/SetSignals.$(ARCH).o \ +$(SRC)/Shutdown.$(ARCH).o \ +$(SRC)/args.$(ARCH).o \ +$(SRC)/help.$(ARCH).o \ +$(SRC)/initialize.$(ARCH).o \ +$(SRC)/update_objects.$(ARCH).o \ +$(SRC)/update_objects_catalog.$(ARCH).o \ +$(SRC)/client_logger.$(ARCH).o \ +$(SRC)/dvolens.$(ARCH).o + +$(DVOLENS): $(INC)/dvolens.h +$(BIN)/dvolens.$(ARCH): $(DVOLENS) + +DVOLENS_CLIENT = \ +$(SRC)/ConfigInit.$(ARCH).o \ +$(SRC)/SetSignals.$(ARCH).o \ +$(SRC)/Shutdown.$(ARCH).o \ +$(SRC)/args.$(ARCH).o \ +$(SRC)/help.$(ARCH).o \ +$(SRC)/initialize.$(ARCH).o \ +$(SRC)/update_objects.$(ARCH).o \ +$(SRC)/update_objects_catalog.$(ARCH).o \ +$(SRC)/client_logger.$(ARCH).o \ +$(SRC)/dvolens_client.$(ARCH).o + +$(DVOLENS_CLIENT): $(INC)/dvolens.h +$(BIN)/dvolens_client.$(ARCH): $(DVOLENS_CLIENT) Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/include/dvolens.h =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/include/dvolens.h (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/include/dvolens.h (revision 37344) @@ -0,0 +1,64 @@ +# include +# include +# include +# include +# include + +typedef enum { + MODE_ERROR = 0, + MODE_UPDATE_OBJECTS, +} DvoLensMode; + +/* global variables set in parameter file */ +# define DVO_MAX_PATH 1024 +char *CATDIR; +char CATMODE[16]; /* raw, mef, split, mysql */ +char CATFORMAT[16]; /* internal, elixir, loneos, panstarrs */ + +int HOST_ID; +char *HOSTDIR; + +int PARALLEL; +int PARALLEL_MANUAL; +int PARALLEL_SERIAL; + +int VERBOSE; +int VERBOSE2; +int UPDATE; +int NTHREADS; + +DvoLensMode MODE; + +SkyRegion UserPatch; +char *UserCatalog; + +/*** dvolens prototypes ***/ +void ConfigInit PROTO((int *argc, char **argv)); +void GetConfig PROTO((char *config, char *field, char *format, int N, void *ptr)); + +DvoLensMode args PROTO((int argc, char **argv)); +int args_client PROTO((int argc, char **argv)); + +DvoLensMode initialize PROTO((int argc, char **argv)); +void initialize_client PROTO((int argc, char **argv)); + +void dvolens_usage PROTO((void)); +void dvolens_help PROTO((int argc, char **argv)); + +void dvolens_client_usage PROTO((void)); +void dvolens_client_help 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 main PROTO((int argc, char **argv)); + +void update_objects PROTO(()); +int update_objects_parallel PROTO((SkyList *sky)); +int update_objects_catalog PROTO((Catalog *catalog)); + +int client_logger_init PROTO((char *dirname)); +int client_logger_message PROTO((char *format,...)); + Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/ConfigInit.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/ConfigInit.c (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/ConfigInit.c (revision 37344) @@ -0,0 +1,50 @@ +# include "dvolens.h" + +void ConfigInit (int *argc, char **argv) { + + char CatdirPhotcodeFile[256]; + + /*** load configuration info ***/ + char *file = SelectConfigFile (argc, argv, "ptolemy"); + char *config = LoadConfigFile (file); + if (!config) { + fprintf (stderr, "ERROR: can't find configuration file %s\n", file); + if (file) 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); + + ScanConfig (config, "CATMODE", "%s", 0, CATMODE); + ScanConfig (config, "CATFORMAT", "%s", 0, CATFORMAT); + + if (*CATMODE == 0) strcpy (CATMODE, "RAW"); + if (*CATFORMAT == 0) strcpy (CATFORMAT, "ELIXIR"); + + sprintf (CatdirPhotcodeFile, "%s/Photcodes.dat", CATDIR); + if (!LoadPhotcodes (CatdirPhotcodeFile, NULL, FALSE)) { + fprintf (stderr, "error loading photcode table %s\n", CatdirPhotcodeFile); + exit (1); + } + + 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: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/SetSignals.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/SetSignals.c (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/SetSignals.c (revision 37344) @@ -0,0 +1,122 @@ +# include "dvolens.h" + +static int Protect = FALSE; +static int Trapped = FALSE; + +void TrapSignal (int sig) { + fprintf (stderr, "trapped signal %d\n", sig); + if (sig == 11) { + fprintf (stderr, "seg fault\n"); + exit (1); + } + if (Protect) { + Trapped = TRUE; + fprintf (stderr, "blocking until protected sections are clear\n"); + return; + } + Shutdown ("halted by signal (trapped)"); +} + +void SetProtect (int mode) { + Protect = mode; + if (Trapped && !Protect) Shutdown ("halted by signal (protect)"); +} + +int SetSignals () { + + int i; + + /* disable almost all signal interrupts */ + for (i = 0; i < 36; i++) { + switch (i) { + /* can't redirect these signals */ + case SIGKILL: /* kill -9: cannot be caught or ignored (POSIX.1-1990) */ + case SIGSTOP: /* SIGSTOP: cannot be caught or ignored (POSIX.1-1990) */ + /* ignore these signals */ + case SIGCHLD: /* child halted: ignore (POSIX.1-1990) */ + case SIGCONT: /* continue - maintain this action (POSIX.1-1990) */ + case SIGTSTP: /* stop signal sent from tty - why ignore? (POSIX.1-1990) */ + case SIGURG: /* socket signal, ignore this (POSIX.1-2001) */ +# ifdef SIGPWR + case SIGPWR: /* power failure - why ignore this? (Sys V) */ +# endif +# ifdef SIGWINCH + case SIGWINCH: /* window resized (4.3BSD) */ +# endif + break; + + default: + signal (i, TrapSignal); + } + } + return (TRUE); +} +/* + + Signal Value Action Comment + ------------------------------------------------------------------------- + SIGHUP 1 A Hangup detected on controlling terminal + or death of controlling process + SIGINT 2 A Interrupt from keyboard + SIGQUIT 3 A Quit from keyboard + SIGILL 4 A Illegal Instruction + SIGABRT 6 C Abort signal from abort(3) + SIGFPE 8 C Floating point exception + SIGKILL 9 AEF Kill signal + SIGSEGV 11 C Invalid memory reference + SIGPIPE 13 A Broken pipe: write to pipe with no readers + SIGALRM 14 A Timer signal from alarm(2) + SIGTERM 15 A Termination signal + SIGUSR1 30,10,16 A User-defined signal 1 + SIGUSR2 31,12,17 A User-defined signal 2 + SIGCHLD 20,17,18 B Child stopped or terminated + SIGCONT 19,18,25 Continue if stopped + SIGSTOP 17,19,23 DEF Stop process + SIGTSTP 18,20,24 D Stop typed at tty + SIGTTIN 21,21,26 D tty input for background process + SIGTTOU 22,22,27 D tty output for background process + + Next various other signals. + + Signal Value Action Comment + --------------------------------------------------------------------- + SIGTRAP 5 CG Trace/breakpoint trap + SIGIOT 6 CG IOT trap. A synonym for SIGABRT + SIGEMT 7,-,7 G + SIGBUS 10,7,10 AG Bus error + SIGSYS 12,-,12 G Bad argument to routine (SVID) + SIGSTKFLT -,16,- AG Stack fault on coprocessor + SIGURG 16,23,21 BG Urgent condition on socket (4.2 BSD) + SIGIO 23,29,22 AG I/O now possible (4.2 BSD) + SIGPOLL AG A synonym for SIGIO (System V) + SIGCLD -,-,18 G A synonym for SIGCHLD + SIGXCPU 24,24,30 AG CPU time limit exceeded (4.2 BSD) + SIGXFSZ 25,25,31 AG File size limit exceeded (4.2 BSD) + SIGVTALRM 26,26,28 AG Virtual alarm clock (4.2 BSD) + SIGPROF 27,27,29 AG Profile alarm clock + SIGPWR 29,30,19 AG Power failure (System V) + SIGINFO 29,-,- G A synonym for SIGPWR + SIGLOST -,-,- AG File lock lost + SIGWINCH 28,28,20 BG Window resize signal (4.3 BSD, Sun) + SIGUNUSED -,31,- AG Unused signal + (Here - denotes that a signal is absent; there where three values are given, the first one is usually valid for alpha and + sparc, the middle one for i386 and ppc, the last one for mips. Signal 29 is SIGINFO / SIGPWR on an alpha but SIGLOST on a + sparc.) + + The letters in the "Action" column have the following meanings: + + A Default action is to terminate the process. + + B Default action is to ignore the signal. + + C Default action is to dump core. + + D Default action is to stop the process. + + E Signal cannot be caught. + + F Signal cannot be ignored. + + G Not a POSIX.1 conformant signal. + +*/ Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/Shutdown.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/Shutdown.c (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/Shutdown.c (revision 37344) @@ -0,0 +1,19 @@ +# include "dvolens.h" + +int Shutdown (char *format, ...) { + va_list argp; + char *formatplus; + + ALLOCATE (formatplus, char, strlen(format) + 2); + strcpy (formatplus, format); + strcat (formatplus, "\n"); + + va_start (argp, format); + vfprintf (stderr, formatplus, argp); + free (formatplus); + va_end (argp); + + fprintf (stderr, "ERROR: dvolens halted\n"); + exit (1); +} + Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/args.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/args.c (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/args.c (revision 37344) @@ -0,0 +1,175 @@ +# include "dvolens.h" + +DvoLensMode args (int argc, char **argv) { + + int N; + + /* specify portion of the 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); + } + + /* specify region file by name (eg n0000/0000.00) */ + UserCatalog = NULL; + if ((N = get_argument (argc, argv, "-catalog"))) { + remove_argument (N, &argc, argv); + UserCatalog = strcreate (argv[N]); + remove_argument (N, &argc, argv); + } + + VERBOSE = VERBOSE2 = FALSE; + if ((N = get_argument (argc, argv, "-v"))) { + VERBOSE = TRUE; + remove_argument (N, &argc, argv); + } + if ((N = get_argument (argc, argv, "-vv"))) { + VERBOSE2 = VERBOSE = TRUE; + remove_argument (N, &argc, argv); + } + + NTHREADS = 0; + if ((N = get_argument (argc, argv, "-threads"))) { + remove_argument (N, &argc, argv); + NTHREADS = atof (argv[N]); + remove_argument (N, &argc, argv); + } + + // XXX for the moment, make this selection manual. it needs to be automatic + // based on the state of the SkyTable + HOST_ID = 0; + 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, + // dvolens 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 dvolens_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; + } + + DvoLensMode mode = MODE_ERROR; + if ((N = get_argument (argc, argv, "-update-objects"))) { + remove_argument (N, &argc, argv); + mode = MODE_UPDATE_OBJECTS; + } + + switch (mode) { + case MODE_UPDATE_OBJECTS: + if (argc != 1) dvolens_usage(); + break; + + default: + fprintf (stderr, "no valid mode selected\n"); + dvolens_usage(); + break; + } + if (argc != 1) dvolens_usage (); + + return mode; +} + +int args_client (int argc, char **argv) { + + int N; + + // by definition, the client is not parallel + PARALLEL = FALSE; + PARALLEL_MANUAL = FALSE; + PARALLEL_SERIAL = FALSE; + + 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) dvolens_client_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) dvolens_client_usage(); + + MODE = MODE_ERROR; + if ((N = get_argument (argc, argv, "-update-objects"))) { + MODE = MODE_UPDATE_OBJECTS; + remove_argument (N, &argc, argv); + } + if (!MODE) dvolens_client_usage(); + + /* specify portion of the 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); + } + + VERBOSE = VERBOSE2 = FALSE; + if ((N = get_argument (argc, argv, "-v"))) { + VERBOSE = TRUE; + remove_argument (N, &argc, argv); + } + if ((N = get_argument (argc, argv, "-vv"))) { + VERBOSE2 = VERBOSE = TRUE; + remove_argument (N, &argc, argv); + } + + UPDATE = FALSE; + if ((N = get_argument (argc, argv, "-update"))) { + remove_argument (N, &argc, argv); + UPDATE = TRUE; + } + + if ((MODE == MODE_UPDATE_OBJECTS) && (argc == 1)) return TRUE; + if (argc != 2) dvolens_client_usage (); + + return TRUE; +} + Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/client_logger.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/client_logger.c (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/client_logger.c (revision 37344) @@ -0,0 +1,37 @@ +# include "dvolens.h" + +// I'm getting unlogged errors and failures. I need a log ouput that the clients can +// write independent of the master + +static FILE *logfile = NULL; +int client_logger_init (char *dirname) { + + char filename[DVO_MAX_PATH]; + + snprintf (filename, DVO_MAX_PATH, "%s/log.dvolens.XXXXXX", dirname); + + int fd = mkstemp (filename); + if (fd == -1) { + fprintf (stderr, "failed to open client logger %s, exiting\n", filename); + exit (50); + } + + logfile = fdopen (fd, "w"); + if (!logfile) { + fprintf (stderr, "failed to fdopen client logger, exiting\n"); + exit (51); + } + return TRUE; +} + +int client_logger_message (char *format,...) { + + va_list argp; + + va_start (argp, format); + vfprintf (logfile, format, argp); + va_end (argp); + + fflush (logfile); + return TRUE; +} Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/dvolens.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/dvolens.c (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/dvolens.c (revision 37344) @@ -0,0 +1,20 @@ +# include "dvolens.h" + +int main (int argc, char **argv) { + + // get configuration info, args + DvoLensMode mode = initialize (argc, argv); + if (!mode) exit (2); + + switch (mode) { + case MODE_UPDATE_OBJECTS: + update_objects (); + exit (0); + + default: + fprintf (stderr, "ERROR: no valid dvolens mode chosen\n"); + exit (2); + } + fprintf (stderr, "IMPOSSIBLE: skipped out of switch?\n"); + exit (1); +} Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/dvolens_client.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/dvolens_client.c (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/dvolens_client.c (revision 37344) @@ -0,0 +1,28 @@ +# include "dvolens.h" + +// dvolens_client is run on a remote host and is responsible for updating the catalogs +// owned by that host. + +// dvolens_client -update-objects : calculate mean lensing parameters + +int main (int argc, char **argv) { + + // get configuration info, args, lockfile (set CATDIR, HOST_ID, HOSTDIR, etc) + initialize_client (argc, argv); + client_logger_init (HOSTDIR); + + switch (MODE) { + case MODE_UPDATE_OBJECTS: + client_logger_message ("start dvolens_client\n"); + update_objects (); + client_logger_message ("updated objects\n"); + break; + + default: + fprintf (stderr, "ERROR: no valid dvolens mode chosen\n"); + exit (2); + } + client_logger_message ("done with dvolens_client\n"); + + exit (0); +} Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/help.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/help.c (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/help.c (revision 37344) @@ -0,0 +1,57 @@ +# include "dvolens.h" + +void dvolens_usage (void) { + fprintf (stderr, "ERROR: USAGE: dvolens -update-objects\n"); + fprintf (stderr, " regions: -region RA RA DEC DEC)\n"); + fprintf (stderr, " or: -catalog (name)\n"); + fprintf (stderr, " use -h for more usage information\n"); + exit (2); +} + +void dvolens_help (int argc, char **argv) { + /* check for help request */ + if (get_argument (argc, argv, "-help")) goto show_help; + if (get_argument (argc, argv, "-h")) goto show_help; + if (argc == 1) dvolens_usage(); + return; + +show_help: + fprintf (stderr, "ERROR: USAGE: dvolens -update-objects\n"); + fprintf (stderr, " regions: -region RA RA DEC DEC)\n"); + fprintf (stderr, " or: -catalog (name)\n"); + fprintf (stderr, " options: \n"); + fprintf (stderr, " -time (start) (stop)\n"); + fprintf (stderr, " -v : verbose output\n"); + fprintf (stderr, " -vv : more verbose output\n"); + fprintf (stderr, " -statmode (mode)\n"); + fprintf (stderr, " -n (nloop)\n"); + fprintf (stderr, " -reset\n"); + fprintf (stderr, " -update\n"); + fprintf (stderr, " \n"); + exit (2); +} + +void dvolens_client_usage (void) { + fprintf (stderr, "ERROR: USAGE: dvolens -update-objects -hostID (hostID) -hostdir (hostdir) [options]\n"); + fprintf (stderr, " use -h for more usage information\n"); + exit (2); +} + +void dvolens_client_help (int argc, char **argv) { + + /* check for help request */ + if (get_argument (argc, argv, "-help")) goto show_help; + if (get_argument (argc, argv, "-h")) goto show_help; + if (argc == 1) dvolens_client_usage(); + return; + +show_help: + fprintf (stderr, "USAGE: dvolens_client -update-objects : determine average magnitudes for objects\n\n"); + fprintf (stderr, " db info : -hostID (hostID) -hostdir (hostdir) -catdir (catdir)\n"); + fprintf (stderr, "other options:\n"); + fprintf (stderr, " -v : verbose output\n"); + fprintf (stderr, " -vv : extra verbose output\n"); + fprintf (stderr, " \n"); + exit (2); +} + Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/initialize.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/initialize.c (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/initialize.c (revision 37344) @@ -0,0 +1,23 @@ +# include "dvolens.h" + +// dvolens has the following modes: +// * -update-objects + +DvoLensMode initialize (int argc, char **argv) { + + dvolens_help (argc, argv); + ConfigInit (&argc, argv); + DvoLensMode mode = args (argc, argv); + if (!mode) exit (2); + + return mode; +} + +void initialize_client (int argc, char **argv) { + + // XXX need to determine which globals can affect relphot_client in either mode and pass appropriately + + dvolens_client_help (argc, argv); + ConfigInit (&argc, argv); + args_client (argc, argv); +} Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/update_objects.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/update_objects.c (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/update_objects.c (revision 37344) @@ -0,0 +1,193 @@ +# include "dvolens.h" + +# define DEBUG 1 + +# define TIMESTAMP(TIME){ \ + gettimeofday (&stopTimer, (void *) NULL); \ + double dtime = DTIME (stopTimer, startTimer); \ + TIME += dtime; \ + gettimeofday (&startTimer, (void *) NULL); } + +void update_objects () { + + int i; + int status; + struct stat filestat; + Catalog catalog; + + double time1 = 0.0; + double time2 = 0.0; + double time3 = 0.0; + double time4 = 0.0; + + // load the current sky table (layout of all SkyRegions) + SkyTable *sky = SkyTableLoadOptimal (CATDIR, NULL, NULL, FALSE, 0, VERBOSE); + SkyTableSetFilenames (sky, CATDIR, "cpt"); + + // determine the populated SkyRegions overlapping the requested area (default depth) + SkyList *skylist = SkyListByPatch (sky, -1, &UserPatch); + + // XXX need to decide how to determine PARALLEL mode... + // in parallel mode, master: PARALLEL is true and HOST_ID == 0, client: PARALLEL is false, HOST_ID != 0 + // in non-parallel mode, master: PARALLEL is false and HOST_ID == 0 + if (PARALLEL && !HOST_ID) { + update_objects_parallel (skylist); + return; + } + + if (VERBOSE) fprintf (stderr, "re-loading catalog data\n"); + + /* load data from each region file */ + 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; + + INITTIME; + + // set up the basic catalog info + char hostfile[1024]; + snprintf (hostfile, 1024, "%s/%s.cpt", HOSTDIR, skylist[0].regions[i]->name); + catalog.filename = HOST_ID ? hostfile : skylist[0].filename[i]; + + // only update existing db tables + status = stat (catalog.filename, &filestat); + if ((status == -1) && (errno == ENOENT)) { + if (VERBOSE) fprintf (stderr, "no file %s, skipping\n", catalog.filename); + continue; + } + TIMESTAMP(time1); + + 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 | LOAD_LENSING | LOAD_LENSOBJ; + catalog.Nsecfilt = GetPhotcodeNsecfilt (); // set the desired number in case we need to create the catalog + + if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE, "w")) { + fprintf (stderr, "ERROR: failure reading catalog %s\n", catalog.filename); + exit (1); + } + if (VERBOSE && (catalog.Naves_disk == 0)) { + fprintf (stderr, "no data in %s, skipping\n", catalog.filename); + dvo_catalog_unlock (&catalog); + dvo_catalog_free (&catalog); + continue; + } + TIMESTAMP(time2); + + update_objects_catalog (&catalog); + TIMESTAMP(time3); + + struct timeval now; + gettimeofday (&now, (void *) NULL); + char *moddate = ohana_sec_to_date (now.tv_sec); + gfits_modify (&catalog.header, "DVOLENS", "%s", 1, moddate); + free (moddate); + + if (!UPDATE) { + dvo_catalog_unlock (&catalog); + dvo_catalog_free (&catalog); + continue; + } + + SetProtect (TRUE); + dvo_catalog_save (&catalog, VERBOSE); + dvo_catalog_unlock (&catalog); + SetProtect (FALSE); + + dvo_catalog_free (&catalog); + TIMESTAMP(time4); + + if (HOST_ID) { + client_logger_message ("updated catalog file %s\n", catalog.filename); + } + } + + fprintf (stderr, "time step 1 %10.3f sec : find catalog\n", time1); + fprintf (stderr, "time step 2 %10.3f sec : load catalog\n", time2); + fprintf (stderr, "time step 3 %10.3f sec : update catalog\n", time3); + fprintf (stderr, "time step 4 %10.3f sec : save catalog\n", time4); +} + +int update_objects_parallel (SkyList *sky) { + + // now launch the dvolens_client jobs to the parallel hosts + + // load the list of hosts + HostTable *table = HostTableLoad (CATDIR, sky->hosts); + if (!table) { + fprintf (stderr, "ERROR: problem with parallel host table\n"); + exit (2); + } + + int i, j; + for (i = 0; i < table->Nhosts; i++) { + + if (sky->Nregions < table->Nhosts) { + // do any of the regions want this host? + int wantThisHost = FALSE; + for (j = 0; j < sky->Nregions; j++) { + if (HostTableTestHost (sky->regions[j], table->hosts[i].hostID)) { + wantThisHost = TRUE; + break; + } + } + if (!wantThisHost) { + // fprintf (stderr, "skip host %s\n", table->hosts[i].hostname); + continue; + } + } + + // 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, "dvolens_client -update-objects -hostID %d -D CATDIR %s -hostdir %s -region %f %f %f %f", + table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, + UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax); + + // options & configs which affect dvolens_client -update-catalogs + // VERBOSE, VERBOSE2, UPDATE + + char tmpline[1024]; + if (VERBOSE) { snprintf (tmpline, 1024, "%s -v", command); strcpy (command, tmpline); } + if (VERBOSE2) { snprintf (tmpline, 1024, "%s -vv", command); strcpy (command, tmpline); } + if (UPDATE) { snprintf (tmpline, 1024, "%s -update", command); strcpy (command, tmpline); } + + fprintf (stderr, "command: %s\n", command); + + if (PARALLEL_MANUAL) continue; + + if (PARALLEL_SERIAL) { + int status = system (command); + if (status) { + fprintf (stderr, "ERROR running dvolens_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 + } + } + + if (PARALLEL_MANUAL) { + fprintf (stderr, "run the dvolens_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) { + fprintf (stderr, "at least one remote client job failed to load data, exiting\n"); + exit (3); + } + } + return (TRUE); +} Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/update_objects_catalog.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/update_objects_catalog.c (revision 37344) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvolens/src/update_objects_catalog.c (revision 37344) @@ -0,0 +1,167 @@ +# include "dvolens.h" + +int update_objects_catalog (Catalog *catalog) { + + off_t i, j, Lj, Mj; + + if (VERBOSE2) fprintf (stderr, "lensobj catalog %d : "OFF_T_FMT" ave, "OFF_T_FMT" meas, "OFF_T_FMT" lensing\n", + catalog->catID, catalog->Naverage, catalog->Nmeasure, catalog->Nlensing); + + // I am making this measurement per filter + int Nsecfilt = GetPhotcodeNsecfilt (); + + // I think we have already allocated lensobj + REALLOCATE (catalog->lensobj, Lensobj, catalog->Naverage * Nsecfilt); + + // Nlensobj tracks how many lensobj entries have I actually generated so far for this + // catalog This is not necessarily 1-to-1 with secfilt entries, but any average entry + // will either have 0 or Nsecfilt lensobj entries + int Nlensobj = 0; + + for (i = 0; i < catalog->Naverage; i++) { + Average *average = &catalog->average[i]; + if (average->Nlensing == 0) continue; + + // NOTE Nlensobj increments by Nsecfilt each time + Lensobj *lensobj = &catalog->lensobj[Nlensobj]; + for (j = 0; j < Nsecfilt; j++) { + dvo_lensobj_init (&lensobj[j], TRUE); // init accumulated values to 0 + } + + // reset the average values + + off_t Loff = average->lensingOffset; + off_t Moff = average->measureOffset; + + // loop over the lensing measurements. for each one, I need to find the corresponding measurement (make an index in lensing?) + for (Lj = 0; Lj < average->Nlensing; Lj++) { + + // find the corresponding measure + Lensing *lensing = &catalog->lensing[Loff + Lj]; + if (!isfinite(lensing->X11_sm_obj)) continue; + + Measure *measure = NULL; + int found = FALSE; + for (Mj = 0; !found && (Mj < average->Nmeasure); Mj++) { + measure = &catalog->measure[Moff + Mj]; + if (measure->imageID != lensing->imageID) continue; + found = TRUE; + } + if (!found) { + fprintf (stderr, "error: cannot match measurement with lensing parameter\n"); + fprintf (stderr, "objID: %d, catID: %d, detID: %d, imageID: %d, averef: %d\n", + lensing->objID, lensing->catID, lensing->detID, lensing->imageID, lensing->averef); + abort(); + } + + // skip measurements that do not match the current photcode + PhotCode *code = GetPhotcodebyCode (measure->photcode); + myAssert (code, "missing photcode?"); + myAssert (code->equiv > -1, "photcode not equivalent to secfilt?"); + + int Nsec = GetPhotcodeNsec (code->equiv); + myAssert (Nsec > -1, "cannot find Nsec?"); + + // if (measure->dbFlags & MEAS_BAD) SKIP_THIS_MEAS(Nbad); + if (measure->psfQF < 0.85) continue; + if (measure->psfQFperf < 0.85) continue; + + + lensobj[Nsec].X11_sm_obj += lensing->X11_sm_obj; + lensobj[Nsec].X12_sm_obj += lensing->X12_sm_obj; + lensobj[Nsec].X22_sm_obj += lensing->X22_sm_obj; + lensobj[Nsec]. E1_sm_obj += lensing-> E1_sm_obj; + lensobj[Nsec]. E2_sm_obj += lensing-> E2_sm_obj; + + lensobj[Nsec].X11_sh_obj += lensing->X11_sh_obj; + lensobj[Nsec].X12_sh_obj += lensing->X12_sh_obj; + lensobj[Nsec].X22_sh_obj += lensing->X22_sh_obj; + lensobj[Nsec]. E1_sh_obj += lensing-> E1_sh_obj; + lensobj[Nsec]. E2_sh_obj += lensing-> E2_sh_obj; + + lensobj[Nsec].X11_sm_psf += lensing->X11_sm_psf; + lensobj[Nsec].X12_sm_psf += lensing->X12_sm_psf; + lensobj[Nsec].X22_sm_psf += lensing->X22_sm_psf; + lensobj[Nsec]. E1_sm_psf += lensing-> E1_sm_psf; + lensobj[Nsec]. E2_sm_psf += lensing-> E2_sm_psf; + + lensobj[Nsec].X11_sh_psf += lensing->X11_sh_psf; + lensobj[Nsec].X12_sh_psf += lensing->X12_sh_psf; + lensobj[Nsec].X22_sh_psf += lensing->X22_sh_psf; + lensobj[Nsec]. E1_sh_psf += lensing-> E1_sh_psf; + lensobj[Nsec]. E2_sh_psf += lensing-> E2_sh_psf; + + lensobj[Nsec]. F_ApR5 += lensing-> F_ApR5; + lensobj[Nsec].dF_ApR5 += SQ(lensing->dF_ApR5); + lensobj[Nsec].sF_ApR5 += SQ(lensing->sF_ApR5); + lensobj[Nsec].fF_ApR5 += lensing->fF_ApR5; + + lensobj[Nsec]. F_ApR6 += lensing-> F_ApR6; + lensobj[Nsec].dF_ApR6 += SQ(lensing->dF_ApR6); + lensobj[Nsec].sF_ApR6 += SQ(lensing->sF_ApR6); + lensobj[Nsec].fF_ApR6 += lensing->fF_ApR6; + + lensobj[Nsec].Nmeas ++; + lensobj[Nsec].photcode = code->equiv; + lensobj[Nsec].objID = lensing->objID; + lensobj[Nsec].catID = lensing->catID; + } + + for (j = 0; j < Nsecfilt; j++) { + if (!lensobj[j].Nmeas) { + dvo_lensobj_init (&lensobj[j], FALSE); + continue; + } + float Nmeas = lensobj[j].Nmeas; + lensobj[j].X11_sm_obj /= Nmeas; + lensobj[j].X12_sm_obj /= Nmeas; + lensobj[j].X22_sm_obj /= Nmeas; + lensobj[j]. E1_sm_obj /= Nmeas; + lensobj[j]. E2_sm_obj /= Nmeas; + lensobj[j].X11_sh_obj /= Nmeas; + lensobj[j].X12_sh_obj /= Nmeas; + lensobj[j].X22_sh_obj /= Nmeas; + lensobj[j]. E1_sh_obj /= Nmeas; + lensobj[j]. E2_sh_obj /= Nmeas; + lensobj[j].X11_sm_psf /= Nmeas; + lensobj[j].X12_sm_psf /= Nmeas; + lensobj[j].X22_sm_psf /= Nmeas; + lensobj[j]. E1_sm_psf /= Nmeas; + lensobj[j]. E2_sm_psf /= Nmeas; + lensobj[j].X11_sh_psf /= Nmeas; + lensobj[j].X12_sh_psf /= Nmeas; + lensobj[j].X22_sh_psf /= Nmeas; + lensobj[j]. E1_sh_psf /= Nmeas; + lensobj[j]. E2_sh_psf /= Nmeas; + + lensobj[j]. F_ApR5 /= Nmeas; + lensobj[j].dF_ApR5 /= Nmeas; + lensobj[j].sF_ApR5 /= Nmeas; + lensobj[j].fF_ApR5 /= Nmeas; + lensobj[j]. F_ApR6 /= Nmeas; + lensobj[j].dF_ApR6 /= Nmeas; + lensobj[j].sF_ApR6 /= Nmeas; + lensobj[j].fF_ApR6 /= Nmeas; + + lensobj[j].dF_ApR5 = sqrt(lensobj[j].dF_ApR5); + lensobj[j].sF_ApR5 = sqrt(lensobj[j].sF_ApR5); + lensobj[j].dF_ApR6 = sqrt(lensobj[j].dF_ApR6); + lensobj[j].sF_ApR6 = sqrt(lensobj[j].sF_ApR6); + + // somewhere in here I should calculate the values of E1, E2, and gamma + } + + average->Nlensobj = Nsecfilt; + average->lensobjOffset = Nlensobj; + Nlensobj += Nsecfilt; + // XVERB |= (catalog->average[j].objID == OBJ_ID_SRC) && (catalog->average[j].catID == CAT_ID_SRC); + // XVERB |= (catalog->average[j].objID == OBJ_ID_DST) && (catalog->average[j].catID == CAT_ID_DST); + } + + catalog->Nlensobj = Nlensobj; + catalog->Nlensobj_disk = 0; + catalog->Nlensobj_off = 0; + REALLOCATE (catalog->lensobj, Lensobj, Nlensobj); + + return (TRUE); +} Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvopsps/src/insert_detections_dvopsps.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvopsps/src/insert_detections_dvopsps.c (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvopsps/src/insert_detections_dvopsps.c (revision 37344) @@ -68,5 +68,5 @@ catalog.Nsecfilt = GetPhotcodeNsecfilt (); - if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE, "w")) { + if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE, "r")) { fprintf (stderr, "ERROR: failure reading catalog %s\n", catalog.filename); exit (1); Index: /branches/eam_branches/ipp-20140813/Ohana/src/dvopsps/src/insert_objects_dvopsps.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/dvopsps/src/insert_objects_dvopsps.c (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/dvopsps/src/insert_objects_dvopsps.c (revision 37344) @@ -56,5 +56,5 @@ catalog.Nsecfilt = GetPhotcodeNsecfilt (); - if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE, "w")) { + if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE, "r")) { fprintf (stderr, "ERROR: failure reading catalog %s\n", catalog.filename); exit (1); Index: /branches/eam_branches/ipp-20140813/Ohana/src/getstar/src/getstar.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/getstar/src/getstar.c (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/getstar/src/getstar.c (revision 37344) @@ -14,5 +14,5 @@ set_db (&db); - sky = SkyTableLoadOptimal (CATDIR, SKY_TABLE, GSCFILE, TRUE, SKY_DEPTH, VERBOSE); + sky = SkyTableLoadOptimal (CATDIR, SKY_TABLE, GSCFILE, FALSE, SKY_DEPTH, VERBOSE); if (!sky) exit (1); Index: /branches/eam_branches/ipp-20140813/Ohana/src/libautocode/def/lensing-ps1-v5.d =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/libautocode/def/lensing-ps1-v5.d (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/libautocode/def/lensing-ps1-v5.d (revision 37344) @@ -2,5 +2,5 @@ EXTNAME DVO_LENSING_PS1_V5 TYPE BINTABLE -SIZE 128 +SIZE 136 DESCRIPTION DVO Lensing Table @@ -44,4 +44,7 @@ FIELD averef, AVE_REF, unsigned int, reference to average entry +FIELD imageID, IMAGE_ID, unsigned int, reference to DVO image ID +FIELD padding, PADDING, unsigned int, reference to DVO image ID + # 28 x float # 4 x int Index: /branches/eam_branches/ipp-20140813/Ohana/src/libautocode/def/lensing.d =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/libautocode/def/lensing.d (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/libautocode/def/lensing.d (revision 37344) @@ -2,5 +2,5 @@ EXTNAME DVO_LENSING TYPE BINTABLE -SIZE 128 +SIZE 136 DESCRIPTION DVO Lensing Table @@ -39,11 +39,14 @@ FIELD fF_ApR6, FLUX_AP_R6, float, Flux inside r = 6 -FIELD detID, DET_ID, unsigned int, detection ID -FIELD objID, OBJ_ID, unsigned int, unique ID for object in table -FIELD catID, CAT_ID, unsigned int, unique ID for table in which object was first realized -FIELD averef, AVE_REF, unsigned int, reference to average entry +FIELD detID, DET_ID, int, detection ID +FIELD objID, OBJ_ID, int, unique ID for object in table +FIELD catID, CAT_ID, int, unique ID for table in which object was first realized +FIELD averef, AVE_REF, int, reference to average entry + +FIELD imageID, IMAGE_ID, int, reference to DVO image ID +FIELD padding, PADDING, int, padding # 28 x float -# 4 x int -# = 128 bytes / detection +# 5 x int + padding +# = 136 bytes / detection # for 2.5G objects + 60 overlaps, expect 15TB for forced warp Index: /branches/eam_branches/ipp-20140813/Ohana/src/libautocode/def/lensobj-ps1-v5.d =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/libautocode/def/lensobj-ps1-v5.d (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/libautocode/def/lensobj-ps1-v5.d (revision 37344) @@ -43,7 +43,8 @@ FIELD E2, E2, float -FIELD objID, OBJ_ID, unsigned int, unique ID for object in table -FIELD catID, CAT_ID, unsigned int, unique ID for table in which object was first realized -FIELD pad, PAD, int +FIELD objID, OBJ_ID, unsigned int, unique ID for object in table +FIELD catID, CAT_ID, unsigned int, unique ID for table in which object was first realized +FIELD photcode, PHOTCODE, short +FIELD Nmeas, NMEAS, short # 31 x float Index: /branches/eam_branches/ipp-20140813/Ohana/src/libautocode/def/lensobj.d =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/libautocode/def/lensobj.d (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/libautocode/def/lensobj.d (revision 37344) @@ -29,21 +29,22 @@ FIELD E2_sh_psf, E2_SH_PSF, float, mean lensing shear psf stars E2 -FIELD F_ApR5, FLUX_AP_R5, float, Flux inside r = 5 -FIELD dF_ApR5, FLUX_ERR_AP_R5, float, Flux inside r = 5 -FIELD sF_ApR5, FLUX_STD_AP_R5, float, Flux inside r = 5 -FIELD fF_ApR5, FLUX_FIL_AP_R5, float, Flux inside r = 5 +FIELD F_ApR5, FLUX_AP_R5, float, Flux inside r = 5 +FIELD dF_ApR5, FLUX_ERR_AP_R5, float, Flux inside r = 5 +FIELD sF_ApR5, FLUX_STD_AP_R5, float, Flux inside r = 5 +FIELD fF_ApR5, FLUX_FIL_AP_R5, float, Flux inside r = 5 -FIELD F_ApR6, FLUX_AP_R6, float, Flux inside r = 5 -FIELD dF_ApR6, FLUX_ERR_AP_R6, float, Flux inside r = 5 -FIELD sF_ApR6, FLUX_AP_R6, float, Flux inside r = 5 -FIELD fF_ApR6, FLUX_AP_R6, float, Flux inside r = 5 +FIELD F_ApR6, FLUX_AP_R6, float, Flux inside r = 5 +FIELD dF_ApR6, FLUX_ERR_AP_R6, float, Flux inside r = 5 +FIELD sF_ApR6, FLUX_AP_R6, float, Flux inside r = 5 +FIELD fF_ApR6, FLUX_AP_R6, float, Flux inside r = 5 -FIELD gamma, GAMMA, float -FIELD E1, E1, float -FIELD E2, E2, float +FIELD gamma, GAMMA, float +FIELD E1, E1, float +FIELD E2, E2, float -FIELD objID, OBJ_ID, unsigned int, unique ID for object in table -FIELD catID, CAT_ID, unsigned int, unique ID for table in which object was first realized -FIELD pad, PAD, int +FIELD objID, OBJ_ID, unsigned int, unique ID for object in table +FIELD catID, CAT_ID, unsigned int, unique ID for table in which object was first realized +FIELD photcode, PHOTCODE, short +FIELD Nmeas, NMEAS, short # 31 x float Index: /branches/eam_branches/ipp-20140813/Ohana/src/libdvo/include/dvo.h =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/libdvo/include/dvo.h (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/libdvo/include/dvo.h (revision 37344) @@ -1039,5 +1039,5 @@ void dvo_lensing_init (Lensing *lensing); -void dvo_lensobj_init (Lensobj *lensobj); +void dvo_lensobj_init (Lensobj *lensobj, int toZero); void InitRegionHosts (RegionHostInfo *hosts, int Nhosts, int NHOSTS); Index: /branches/eam_branches/ipp-20140813/Ohana/src/libdvo/src/dvo_catalog.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/libdvo/src/dvo_catalog.c (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/libdvo/src/dvo_catalog.c (revision 37344) @@ -389,41 +389,39 @@ lensing->catID = -1; lensing->averef = 0; + + lensing->imageID = -1; + lensing->padding = 0; } // init all data, or just catalog data -void dvo_lensobj_init (Lensobj *lensobj) { - lensobj->X11_sm_obj = NAN; - lensobj->X12_sm_obj = NAN; - lensobj->X22_sm_obj = NAN; - lensobj->E1_sm_obj = NAN; - lensobj->E2_sm_obj = NAN; - - lensobj->X11_sh_obj = NAN; - lensobj->X12_sh_obj = NAN; - lensobj->X22_sh_obj = NAN; - lensobj->E1_sh_obj = NAN; - lensobj->E2_sh_obj = NAN; - - lensobj->X11_sm_psf = NAN; - lensobj->X12_sm_psf = NAN; - lensobj->X22_sm_psf = NAN; - lensobj->E1_sm_psf = NAN; - lensobj->E2_sm_psf = NAN; - - lensobj->X11_sh_psf = NAN; - lensobj->X12_sh_psf = NAN; - lensobj->X22_sh_psf = NAN; - lensobj->E1_sh_psf = NAN; - lensobj->E2_sh_psf = NAN; - - lensobj->F_ApR5 = NAN; - lensobj->dF_ApR5 = NAN; - lensobj->sF_ApR5 = NAN; - lensobj->fF_ApR5 = NAN; - - lensobj->F_ApR6 = NAN; - lensobj->dF_ApR6 = NAN; - lensobj->sF_ApR6 = NAN; - lensobj->fF_ApR6 = NAN; +void dvo_lensobj_init (Lensobj *lensobj, int toZero) { + lensobj->X11_sm_obj = toZero ? 0 : NAN; + lensobj->X12_sm_obj = toZero ? 0 : NAN; + lensobj->X22_sm_obj = toZero ? 0 : NAN; + lensobj->E1_sm_obj = toZero ? 0 : NAN; + lensobj->E2_sm_obj = toZero ? 0 : NAN; + lensobj->X11_sh_obj = toZero ? 0 : NAN; + lensobj->X12_sh_obj = toZero ? 0 : NAN; + lensobj->X22_sh_obj = toZero ? 0 : NAN; + lensobj->E1_sh_obj = toZero ? 0 : NAN; + lensobj->E2_sh_obj = toZero ? 0 : NAN; + lensobj->X11_sm_psf = toZero ? 0 : NAN; + lensobj->X12_sm_psf = toZero ? 0 : NAN; + lensobj->X22_sm_psf = toZero ? 0 : NAN; + lensobj->E1_sm_psf = toZero ? 0 : NAN; + lensobj->E2_sm_psf = toZero ? 0 : NAN; + lensobj->X11_sh_psf = toZero ? 0 : NAN; + lensobj->X12_sh_psf = toZero ? 0 : NAN; + lensobj->X22_sh_psf = toZero ? 0 : NAN; + lensobj->E1_sh_psf = toZero ? 0 : NAN; + lensobj->E2_sh_psf = toZero ? 0 : NAN; + lensobj->F_ApR5 = toZero ? 0 : NAN; + lensobj->dF_ApR5 = toZero ? 0 : NAN; + lensobj->sF_ApR5 = toZero ? 0 : NAN; + lensobj->fF_ApR5 = toZero ? 0 : NAN; + lensobj->F_ApR6 = toZero ? 0 : NAN; + lensobj->dF_ApR6 = toZero ? 0 : NAN; + lensobj->sF_ApR6 = toZero ? 0 : NAN; + lensobj->fF_ApR6 = toZero ? 0 : NAN; lensobj->gamma = NAN; @@ -434,5 +432,6 @@ lensobj->catID = -1; - lensobj->pad = 0; + lensobj->photcode = 0; + lensobj->Nmeas = 0; } Index: /branches/eam_branches/ipp-20140813/Ohana/src/libdvo/src/dvo_convert_PS1_V5.c =================================================================== --- /branches/eam_branches/ipp-20140813/Ohana/src/libdvo/src/dvo_convert_PS1_V5.c (revision 37343) +++ /branches/eam_branches/ipp-20140813/Ohana/src/libdvo/src/dvo_convert_PS1_V5.c (revision 37344) @@ -491,4 +491,6 @@ out[i].averef = in[i].averef; + + out[i].imageID = in[i].imageID; } return (out); @@ -543,4 +545,5 @@ out[i].averef = in[i].averef; + out[i].imageID = in[i].imageID; } return (out); @@ -555,5 +558,5 @@ for (i = 0; i < Nvalues; i++) { - dvo_lensobj_init (&out[i]); + dvo_lensobj_init (&out[i], FALSE); out[i].X11_sm_obj = in[i].X11_sm_obj; @@ -598,4 +601,6 @@ out[i].catID = in[i].catID; + out[i].photcode = in[i].photcode; + out[i].Nmeas = in[i].Nmeas; } return (out); @@ -651,4 +656,7 @@ out[i].objID = in[i].objID; out[i].catID = in[i].catID; + + out[i].photcode = in[i].photcode; + out[i].Nmeas = in[i].Nmeas; } return (out);