Index: /tags/delstar-1-2/Ohana/src/delstar/Makefile =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/Makefile (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/Makefile (revision 4683) @@ -0,0 +1,75 @@ +include ../../Configure +CC = gcc -g # override this for now for minimal test +HOME = $(ROOT)/src/delstar +CONFIG = $(ROOT)/config +PROGRAM = delstar + +BIN = $(HOME)/bin +INC = $(HOME)/include +SRC = $(HOME)/src +MAN = $(HOME)/doc +DESTBIN = $(LBIN) +DESTLIB = $(LLIB) +DESTINC = $(LINC) +DESTMAN = $(LMAN) + +# +INCS = -I$(INC) -I$(LINC) -I$(XINC) +LIBS = -L$(LLIB) -lFITS -lohana -lm +CFLAGS = -o $*.$(ARCH).o $(INCS) -O2 +CCFLAGS = $(INCS) $(LIBS) + +DELSTAR = \ +$(SRC)/delstar.$(ARCH).o \ +$(SRC)/SetSignals.$(ARCH).o \ +$(SRC)/ConfigInit.$(ARCH).o \ +$(SRC)/args.$(ARCH).o \ +$(SRC)/image-db.$(ARCH).o \ +$(SRC)/delete_imagename.$(ARCH).o \ +$(SRC)/delete_imagefile.$(ARCH).o \ +$(SRC)/delete_times.$(ARCH).o \ +$(SRC)/gimages.$(ARCH).o \ +$(SRC)/gregion_image.$(ARCH).o \ +$(SRC)/gregion_patch.$(ARCH).o \ +$(SRC)/RegionOps.$(ARCH).o \ +$(SRC)/gcatalog.$(ARCH).o \ +$(SRC)/wcatalog.$(ARCH).o \ +$(SRC)/find_image_db.$(ARCH).o \ +$(SRC)/find_matches.$(ARCH).o \ +$(SRC)/parse_time.$(ARCH).o \ +$(SRC)/convert.$(ARCH).o \ +$(SRC)/sort_lists.$(ARCH).o \ +$(SRC)/check_permissions.$(ARCH).o \ +$(SRC)/dimages.$(ARCH).o + +OBJ = $(DELSTAR) + +default: $(PROGRAM) + +$(DELSTAR): $(INC)/delstar.h + +# dependancy rules for binary code ########################## +$(PROGRAM): $(BIN)/$(PROGRAM).$(ARCH) + +$(BIN)/$(PROGRAM).$(ARCH): $(OBJ) + @if [ ! -d $(BIN) ]; then mkdir -p $(BIN); fi + $(CC) $(OBJ) -o $(BIN)/$(PROGRAM).$(ARCH) $(CCFLAGS) + +install: $(DESTBIN)/$(PROGRAM) + +$(DESTBIN)/$(PROGRAM): $(BIN)/$(PROGRAM).$(ARCH) + @if [ ! -d $(DESTBIN) ]; then mkdir -p $(DESTBIN); fi + rm -f $(DESTBIN)/$(PROGRAM) + cp $(BIN)/$(PROGRAM).$(ARCH) $(DESTBIN)/$(PROGRAM) + +# utilities ################################################# +clean: + rm -f $(BIN)/*.$(ARCH) + rm -f `find . -name "*.o"` + rm -f `find . -name "*~"` + rm -f `find . -name "#*"` + +.SUFFIXES: .$(ARCH).o + +.c.$(ARCH).o: + $(CC) $(CFLAGS) -c $< Index: /tags/delstar-1-2/Ohana/src/delstar/bin/.cvsignore =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/bin/.cvsignore (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/bin/.cvsignore (revision 4683) @@ -0,0 +1,1 @@ +*.linux *.lin64 *.sol *.sun *.sid *.hp *.irix Index: /tags/delstar-1-2/Ohana/src/delstar/include/addstar.h =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/include/addstar.h (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/include/addstar.h (revision 4683) @@ -0,0 +1,13 @@ +# include +# include + +/* global variables set in parameter file */ +char ImageCat[256]; +char ImageTemplate[256]; +char CatTemplate[256]; +char GSCFILE[256]; +char CATDIR[256]; +double NSIGMA; +double ALPHA; +int VERBOSE; + Index: /tags/delstar-1-2/Ohana/src/delstar/include/delstar.h =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/include/delstar.h (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/include/delstar.h (revision 4683) @@ -0,0 +1,90 @@ +# include +# include +# include + +typedef struct { + Coords coords; + float *X, *Y; + int *N; + double RA[2], DEC[2]; + double Area, density, spacing; +} CatStats; + +/* global variables set in parameter file */ +char ImageCat[256]; +char ImageTemplate[256]; +char CatTemplate[256]; +char GSCFILE[256]; +char CATDIR[256]; +double NSIGMA; +double ALPHA; +int VERBOSE; +int ORPHAN; +int MISSED; + +time_t START; +time_t END; +PhotCode *PHOTCODE; + +int MODE; +enum {MODE_IMAGENAME, MODE_IMAGEFILE, MODE_TIME, MODE_ORPHAN, MODE_MISSED}; + +enum {TIME_NONE, TIME_JD, TIME_MJD, TIME_DATE, TIME_DAYS, TIME_HOURS, TIME_MINUTES, TIME_SECONDS}; + +char DateKeyword[64], DateMode[64], UTKeyword[64], MJDKeyword[64], JDKeyword[64]; + +/*** conversion functions ***/ +FILE *GetDB PROTO((int *state)); +Image *GetImages PROTO((int *nimage)); +GSCRegion *LoadRegions PROTO((int *nregions)); +Image *find_images PROTO((FILE *f, CatStats *catstats, int *Nimages)); +int *find_images_data PROTO((Image *timage, int *nlist)); +int *find_images_name PROTO((char *filename, int *nlist)); +int *find_images_time PROTO((time_t start, time_t end, PhotCode *code, int *nlist)); +Image *gimages PROTO((char *filename)); +GSCRegion *gregion_image PROTO((Image *image, int *Nregions)); +GSCRegion *gregion_patch PROTO((GSCRegion *patch, int *nregions)); +Image *gtimes PROTO((int *NIMAGE)); +char *sec_to_date PROTO((time_t second)); +void ConfigInit PROTO((int *argc, char **argv)); +int FindDecBand PROTO((double dec, double *DEC0, double *DEC1)); +int FindRegionByPoint PROTO((GSCRegion *region, double ra, double dec, GSCRegion *gsc, int Ngsc)); +int FindRegionDecBandStart PROTO((GSCRegion *region, int Nregion, double dec)); +int FindRegionDecBandStop PROTO((GSCRegion *region, int Nregion, double dec)); +int SetImages PROTO((Image *new, int Nnew)); +void SetProtect PROTO((int mode)); +int SetSignals PROTO(()); +int Shutdown PROTO((char *format, ...)); +void TrapSignal PROTO((int sig)); +int args PROTO((int *argc, char **argv)); +void check_permissions PROTO((char *basefile)); +int chk_time PROTO((char *line)); +time_t date_to_sec PROTO((char *date)); +void delete_imagefile PROTO((char *filename)); +void delete_imagename PROTO((char *name)); +void delete_times PROTO(()); +int dimages PROTO((int *imlist, int Nimlist)); +int dms_to_ddd PROTO((double *Value, char *string)); +int edge_check PROTO((double *x1, double *y1, double *x2, double *y2)); +void find_matches PROTO((Catalog *catalog, int photcode, int start, int end)); +int gcatalog PROTO((Catalog *catalog)); +void help PROTO(()); +int hms_to_deg PROTO((double *h0, double *h1, double *d0, double *d1, char *string)); +time_t jd_to_sec PROTO((double jd)); +int load_image_db PROTO(()); +void lock_image_db PROTO(()); +void match_images PROTO((Catalog *catalog, Image *image, int Nimage)); +time_t mjd_to_sec PROTO((double mjd)); +double opening_angle PROTO((double x1, double y1, double x2, double y2, double x3, double y3)); +int parse_time PROTO((Header *header)); +int save_image_db PROTO(()); +double sec_to_jd PROTO((time_t second)); +double sec_to_mjd PROTO((time_t second)); +void sort_lists PROTO((float *X, float *Y, int *S, int N)); +void sort_regions PROTO((GSCRegion *region, int N)); +int str_to_dtime PROTO((char *line, double *second)); +int str_to_radec PROTO((double *ra, double *dec, char *str1, char *str2)); +int str_to_time PROTO((char *line, time_t *second)); +void unlock_image_db PROTO(()); +void usage PROTO(()); +int wcatalog PROTO((Catalog *catalog)); Index: /tags/delstar-1-2/Ohana/src/delstar/src/ConfigInit.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/ConfigInit.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/ConfigInit.c (revision 4683) @@ -0,0 +1,40 @@ +# include "delstar.h" + +void ConfigInit (int *argc, char **argv) { + + char *config, *file; + char PhotCodeFile[256]; + + /*** load configuration info ***/ + file = SelectConfigFile (argc, argv, "ptolemy"); + config = LoadConfigFile (file); + if (config == (char *) NULL) { + fprintf (stderr, "ERROR: can't find configuration file %s\n", file); + if (file != (char *) NULL) free (file); + exit (1); + } + if (VERBOSE) fprintf (stderr, "loaded config file: %s\n", file); + + ScanConfig (config, "NSIGMA", "%lf", 0, &NSIGMA); + ScanConfig (config, "ALPHA", "%lf", 0, &ALPHA); + ScanConfig (config, "IMAGE_CATALOG", "%s", 0, ImageCat); + ScanConfig (config, "IMAGE_CATALOG_TEMPLATE", "%s", 0, ImageTemplate); + ScanConfig (config, "CATALOG_TEMPLATE", "%s", 0, CatTemplate); + ScanConfig (config, "GSCFILE", "%s", 0, GSCFILE); + ScanConfig (config, "CATDIR", "%s", 0, CATDIR); + ScanConfig (config, "PHOTCODE_FILE", "%s", 0, PhotCodeFile); + + ScanConfig (config, "DATE-KEYWORD", "%s", 0, DateKeyword); + ScanConfig (config, "DATE-MODE", "%s", 0, DateMode); + ScanConfig (config, "UT-KEYWORD", "%s", 0, UTKeyword); + ScanConfig (config, "MJD-KEYWORD", "%s", 0, MJDKeyword); + ScanConfig (config, "JD-KEYWORD", "%s", 0, JDKeyword); + + if (!LoadPhotcodes (PhotCodeFile)) { + fprintf (stderr, "error loading photcodes\n"); + exit (1); + } + + free (config); + free (file); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/RegionOps.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/RegionOps.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/RegionOps.c (revision 4683) @@ -0,0 +1,235 @@ +# include "delstar.h" + +# define NBANDS 24 + +double DecBands[] = { + -90.0, -82.5, -75.0, -67.5, -60.0, + -52.5, -45.0, -37.5, -30.0, -22.5, + -15.0, -07.5, +00.0, +07.5, +15.0, + +22.5, +30.0, +37.5, +45.0, +52.5, + +60.0, +67.5, +75.0, +82.5, +90.0 +}; + +char *DecNames[] = { + "s8230", "s7500", "s6730", "s6000", "s5230", + "s4500", "s3730", "s3000", "s2230", "s1500", + "s0730", "s0000", "n0000", "n0730", "n1500", + "n2230", "n3000", "n3730", "n4500", "n5230", + "n6000", "n6730", "n7500", "n8230", +}; + +GSCRegion *LoadRegions (int *nregions) { + + int Nregions, i, NDec, Nx, Ny; + char temp[50], file[256]; + double RA0, RA1, DEC0, DEC1, D0, D1; + FILE *f; + GSCRegion *regions; + Header header; + Header theader; + Matrix matrix; + FTable table; + + f = fopen (GSCFILE, "r"); + if (f == NULL) { + fprintf (stderr, "ERROR: can't find GSC region file %s\n", GSCFILE); + exit (1); + } + + /* init & load in table data */ + table.header = &theader; + if (!fits_fread_header (f, &header)) Shutdown ("can't read primary header for GSC Region table"); + if (!fits_fread_matrix (f, &matrix, &header)) Shutdown ("can't read primary matrix for GSC Region table"); + if (!fits_fread_ftable (f, &table, "REGIONS")) Shutdown ("can't read GSC Region table"); + + fits_scan (table.header, "NAXIS1", "%d", 1, &Nx); + fits_scan (table.header, "NAXIS2", "%d", 1, &Ny); + + Nregions = 0; + ALLOCATE (regions, GSCRegion, Ny + 2); /* allow room for N & S poles */ + + /* pole region is artificial - not in table */ + sprintf (file, "%s/n8230/pole.cpt", CATDIR); + regions[Nregions].DEC[0] = 86.25; + regions[Nregions].DEC[1] = 93.75; + regions[Nregions].RA[0] = 0.0; + regions[Nregions].RA[1] = 360.0; + strcpy (regions[Nregions].filename, file); + Nregions ++; + + /* fill out GSCRegion structure */ + for (i = 0; i < Ny; i++) { + strncpy (temp, &table.buffer[i*48], 48); + temp[49] = 0; + hms_to_deg (&RA0, &RA1, &DEC0, &DEC1, &temp[7]); + if (RA1 < RA0) RA1 += 360.0; + + if (DEC0 < DEC1) { + regions[Nregions].DEC[0] = DEC0; + regions[Nregions].DEC[1] = DEC1; + } else { + regions[Nregions].DEC[0] = DEC1; + regions[Nregions].DEC[1] = DEC0; + } + /* skip the pole entries */ + if (regions[Nregions].DEC[1] > 86.25) continue; + + regions[Nregions].RA[0] = RA0; + regions[Nregions].RA[1] = RA1; + + /** convert DEC0, DEC1 to NBigDec **/ + NDec = FindDecBand (0.5*(DEC0 + DEC1), &D0, &D1); + if (NDec == -1) Shutdown ("programming / table error"); + + temp[5] = 0; + sprintf (file, "%s/%s/%s.cpt", CATDIR, DecNames[NDec], &temp[1]); + + strcpy (regions[Nregions].filename, file); + Nregions ++; + } + if (Nregions > Ny + 2) { + fprintf (stderr, "ERROR: too many regions\n"); + exit (2); + } + + sort_regions (regions, Nregions); + *nregions = Nregions; + return (regions); +} + +/* given dec value, find upper and lower bounds of enclosing dec band */ +int FindDecBand (double dec, double *DEC0, double *DEC1) { + + int i; + + for (i = 0; i < NBANDS; i++) { + if (DecBands[i + 0] > dec) continue; + if (DecBands[i + 1] <= dec) continue; + *DEC0 = DecBands[i+0]; + *DEC1 = DecBands[i+1]; + return (i); + } + return (-1); +} + +/* note that region.DEC[0] is sorted, but not necessarily DEC[1] */ +int FindRegionDecBandStart (GSCRegion *region, int Nregion, double dec) { + + int N, Nlo, Nhi; + double D0, D1; + + if (dec <= 90) return (0); + + /* D0, D1 are upper lower band boundaries */ + /* all regions in this band have DEC[0] >= D0 */ + N = FindDecBand (dec, &D0, &D1); + if (N == -1) Shutdown ("programming / table error"); + + /* bracket value of interest */ + Nlo = 0; Nhi = Nregion - 1; + while (Nhi - Nlo > 10) { + N = 0.5*(Nlo + Nhi); + if (region[N].DEC[0] < D0) { + Nlo = N; + } else { + Nhi = N; + } + } + /* Nlo - Nhi bracket the transition across dec */ + /* v[Nlo] < v <= v[Nhi] */ + + /* find the last entry < dec */ + for (N = Nhi; N >= Nlo; N--) { + if (region[N].DEC[0] < D0) return (N); + } + return (-1); + /* this will fail if transition value is not in the data range */ +} + +/* note that region.DEC[0] is sorted, but not necessarily DEC[1] */ +int FindRegionDecBandStop (GSCRegion *region, int Nregion, double dec) { + + int N, Nlo, Nhi; + double D0, D1; + + if (dec >= 90) return (Nregion - 1); + + /* D0, D1 are upper lower band boundaries */ + /* all regions in this band have DEC[0] < D1 */ + N = FindDecBand (dec, &D0, &D1); + if (N == -1) Shutdown ("programming / table error"); + + /* bracket value of interest */ + Nlo = 0; Nhi = Nregion - 1; + while (Nhi - Nlo > 10) { + N = 0.5*(Nlo + Nhi); + if (region[N].DEC[0] > D1) { + Nhi = N; + } else { + Nlo = N; + } + } + /* Nlo - Nhi bracket the transition across dec */ + /* v[Nlo] <= v < v[Nhi] */ + + /* find the last entry < dec */ + for (N = Nlo; N <= Nhi; N++) { + if (region[N].DEC[0] > D1) return (N); + } + fprintf (stderr, "fail dec band stop: %f: %f %f, %d - %d, %f %f\n", + dec, D0, D1, Nlo, Nhi, region[Nlo].DEC[0], region[Nlo].DEC[1]); + return (-1); + /* this will fail if transition value is not in the data range */ +} + +int FindRegionByPoint (GSCRegion *region, double ra, double dec, GSCRegion *gsc, int Ngsc) { + + int i, Ns; + + Ns = FindRegionDecBandStart (gsc, Ngsc, dec); + if (Ns < 0) return (FALSE); + + for (i = Ns; i < Ngsc; i++) { + if (gsc[i].RA[0] >= ra) continue; + if (gsc[i].RA[1] < ra) continue; + if (gsc[i].DEC[0] >= dec) continue; + if (gsc[i].DEC[1] < dec) continue; + region[0] = gsc[i]; + return (TRUE); + } + return (FALSE); +} + +/* bracket search: + + Nlo T N Nhi + X X X X X X X V X X X X X X X X X X X X X X + Nhi + + Nlo N T Nhi + X X X X X X X X X X X X X X X V X X X X X X + Nlo + +*/ + +/* + if (buffer[i*48 + 19] == ' ') continue; + strncpy (temp, &buffer[i*48 + 19], 20); temp[20] = 0; + hms_to_deg (&RA, &DEC, temp, ' ', 3); + + + dBigRA = 360.0 / (int)(0.5 + 48*cos(3.1415927*0.5*(BigDecBounds[NBigDec] + BigDecBounds[NBigDec + 1])/180.0)); + + NBigRA = ra / dBigRA; + + NBig = NBigRA; + for (i = 0; (i < 12) && (i < NBigDec); i++) { + NBig += NBigRASections[i]; + } + for (i = 13; (i < 24) && (i < NBigDec); i++) { + NBig += NBigRASections[i]; + } + + + fprintf (stderr, "%d %d %d %d %f %f -> %f %f %f\n", NBigDec, NBigRA, NBigRASections[NBigDec], NBig, ra, dec, BigDecBounds[NBigDec], BigDecBounds[NBigDec + 1], dBigRA); +*/ Index: /tags/delstar-1-2/Ohana/src/delstar/src/SetSignals.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/SetSignals.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/SetSignals.c (revision 4683) @@ -0,0 +1,118 @@ +# include "delstar.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 this signals */ + case SIGKILL: /* kill -9: cannot be caught or ignored */ + case SIGSTOP: /* SIGSTOP: cannot be caught or ignored */ + /* ignore these signals */ + case SIGCHLD: /* child halted: ignore */ + case SIGPWR: /* power failure - why ignore this? */ + case SIGWINCH: /* window resized */ + case SIGCONT: /* continue - maintain this action */ + case SIGTSTP: /* stop signal sent from tty - why ignore? */ + case SIGURG: /* socket signal, ignore this */ + 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: /tags/delstar-1-2/Ohana/src/delstar/src/aregion.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/aregion.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/aregion.c (revision 4683) @@ -0,0 +1,186 @@ +# include "delstar.h" + +double BigDecBounds[] = {0.0, 7.5, 15.0, 22.5, 30.0, 37.5, 45.0, + 52.5, 60.0, 67.5, 75.0, 82.5, 90.0, + 0.0, -7.5, -15.0, -22.5, -30.0, -37.5, -45.0, + -52.5, -60.0, -67.5, -75.0, -82.5, -90.0}; +char *DecSections[] = {"N0000", "N0730", "N1500", "N2230", "N3000", "N3730", "N4500", + "N5230", "N6000", "N6730", "N7500", "N8230", "weirdness", + "S0000", "S0730", "S1500", "S2230", "S3000", "S3730", "S4500", + "S5230", "S6000", "S6730", "S7500", "S8230", "weirdness"}; + +char *Dec2Sections[] = {"n0000", "n0730", "n1500", "n2230", "n3000", "n3730", "n4500", + "n5230", "n6000", "n6730", "n7500", "n8230", "weirdness", + "s0000", "s0730", "s1500", "s2230", "s3000", "s3730", "s4500", + "s5230", "s6000", "s6730", "s7500", "s8230", "weirdness"}; + +char *disk[] = {"disk 1", "disk 1", "disk 1", "disk 1", "disk 1", "disk 1", "disk 1", + "disk 1", "disk 1", "disk 1", "disk 1", "disk 1", "weirdness", + "disk 1", "disk 2", "disk 2", "disk 2", "disk 2", "disk 2", "disk 2", + "disk 2", "disk 2", "disk 2", "disk 2", "disk 2", "weirdness"}; + +int NBigRASections [] = {48, 47, 45, 43, 40, 36, 32, 28, 21, 15, 9, 3, 3, 48, 47, 45, 43, 40, 36, 32, 28, 21, 15, 9, 3, 3}; + +int NDecLines[] = {593, 584, 551, 530, 522, 465, 406, 362, 280, 198, 123, 24, + 0, 597, 578, 574, 577, 534, 499, 442, 376, 294, 212, 144, 48}; + +/* find region file which contains ra, dec */ +aregion (region, f, ra, dec) +GSCRegion region[]; +FILE *f; +double ra, dec; +{ + + + char buffer[28800], temp[50], file[50]; + double dr, dd; + double RA0, RA1, DEC0, DEC1; + int i, NBigDec, NLINES, done, nregion; + + while (ra < 0) { ra += 360.0; } + while (ra >= 360.0) { ra -= 360.0; } + + if (dec >= 86.25) { + sprintf (file, "%s/n8230/pole.cpt\0", CATDIR); + region[0].DEC[0] = 86.25; + region[0].DEC[1] = 93.75; + region[0].RA[0] = -180.0; + region[0].RA[1] = 540.0; + strcpy (region[0].filename, file); + return; + } + + NBigDec = -1; + for (i = 0; i < 12; i++) { +# ifdef DEBUG + fprintf (stderr, "%d %f %f %f\n", i, dec, BigDecBounds[i], BigDecBounds[i+1]); +# endif + if ((dec >= BigDecBounds[i]) && (dec < BigDecBounds[i+1])) { + NBigDec = i; + break; + } + } + if (NBigDec < 0) { + for (i = 13; i < 24; i++) { +# ifdef DEBUG + fprintf (stderr, "%d %f %f %f\n", i, dec, BigDecBounds[i], BigDecBounds[i+1]); +# endif + if ((dec < BigDecBounds[i]) && (dec >= BigDecBounds[i+1])) { + NBigDec = i; + break; + } + } + } + if (NBigDec < 0) { + fprintf (stderr, "dec out of range: %f\n", dec); + } + + NLINES = 0; + for (i = 0; i < NBigDec; i++) { + NLINES += NDecLines[i]; + } + fseek (f, 5*2880 + 48*NLINES, SEEK_SET); + + done = FALSE; + Fread (buffer, 1, 48*NDecLines[NBigDec], f, "char"); + for (i = 0; !done && (i < NDecLines[NBigDec]); i++) { + strncpy (temp, &buffer[i*48], 48); + temp[49] = 0; + hms_to_deg (&RA0, &RA1, &DEC0, &DEC1, &temp[7]); + if (RA1 < RA0) RA1 += 360.0; +# ifdef DEBUG + fprintf (stderr, "%f %f %f %f %f %f %s\n", DEC0, dec, DEC1, RA0, ra, RA1, temp); +# endif + if ((dec >= 0) && (dec >= DEC0) && (dec < DEC1) && (ra >= RA0) && (ra < RA1)) { + done = TRUE; + } + if ((dec < 0) && (dec < DEC0) && (dec >= DEC1) && (ra >= RA0) && (ra < RA1)) { + done = TRUE; + } + } + if (!done) { + fprintf (stderr, "error in search: %f %f\n", ra, dec); + exit (1); + } + temp[5] = 0; + sprintf (file, "%s/%s/%s.cpt\0", CATDIR, Dec2Sections[NBigDec],&temp[1]); + if (DEC0 < DEC1) { + region[0].DEC[0] = DEC0; + region[0].DEC[1] = DEC1; + } else { + region[0].DEC[0] = DEC1; + region[0].DEC[1] = DEC0; + } + region[0].RA[0] = RA0; + region[0].RA[1] = RA1; + strcpy (region[0].filename, file); + return; +} + + +/**********/ +int hms_to_deg (h0, h1, d0, d1, string) + char *string; + double *h0, *h1, *d0, *d1; +{ + + int flag_d0, flag_d1, flag_h0, flag_h1; + double tmp; + + *d0 = *h0 = *d1 = *h1 = 0; + + flag_h0 = dparse (h0, 1, string); + flag_h1 = dparse (h1, 4, string); + flag_d0 = dparse (d0, 7, string); + flag_d1 = dparse (d1, 9, string); + *h0 *= flag_h0; + *h1 *= flag_h1; + *d0 *= flag_d0; + *d1 *= flag_d1; + + dparse (&tmp, 2, string); + *h0 += tmp/60.0; + dparse (&tmp, 3, string); + *h0 += tmp/3600.0; + + dparse (&tmp, 5, string); + *h1 += tmp/60.0; + dparse (&tmp, 6, string); + *h1 += tmp/3600.0; + + dparse (&tmp, 8, string); + *d0 += tmp/60.0; + + dparse (&tmp, 10, string); + *d1 += tmp/60.0; + + *h0 *= 15*flag_h0; + *h1 *= 15*flag_h1; + *d0 *= flag_d0; + *d1 *= flag_d1; + + return (TRUE); +} + + +/* + if (buffer[i*48 + 19] == ' ') continue; + strncpy (temp, &buffer[i*48 + 19], 20); temp[20] = 0; + hms_to_deg (&RA, &DEC, temp, ' ', 3); + + + dBigRA = 360.0 / (int)(0.5 + 48*cos(3.1415927*0.5*(BigDecBounds[NBigDec] + BigDecBounds[NBigDec + 1])/180.0)); + + NBigRA = ra / dBigRA; + + NBig = NBigRA; + for (i = 0; (i < 12) && (i < NBigDec); i++) { + NBig += NBigRASections[i]; + } + for (i = 13; (i < 24) && (i < NBigDec); i++) { + NBig += NBigRASections[i]; + } + + + fprintf (stderr, "%d %d %d %d %f %f -> %f %f %f\n", NBigDec, NBigRA, NBigRASections[NBigDec], NBig, ra, dec, BigDecBounds[NBigDec], BigDecBounds[NBigDec + 1], dBigRA); +*/ Index: /tags/delstar-1-2/Ohana/src/delstar/src/args.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/args.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/args.c (revision 4683) @@ -0,0 +1,100 @@ +# include "delstar.h" + +void help () { + + fprintf (stderr, "USAGE:\n"); + fprintf (stderr, " delstar (filename) or (start range) or -orphan (region)\n"); + fprintf (stderr, " delstar -time (start) (stop/range)\n"); + fprintf (stderr, " delstar -orphan (region)\n"); + fprintf (stderr, " delstar -missed (region)\n\n"); + fprintf (stderr, " optional flags:\n"); + fprintf (stderr, " -v : verbose mode\n"); + fprintf (stderr, " -photcode (code) : restrict by photcode\n"); + fprintf (stderr, "\n"); + exit (2); + +} + +void usage () { + fprintf (stderr, "USAGE: delstar (filename) / [optional mode] : -h for help\n"); + exit (2); +} + +int args (int *argc, char **argv) { + + int N; + double trange; + time_t tmp; + + /* check for help request */ + if (get_argument (*argc, argv, "-help") || + get_argument (*argc, argv, "-h")) { + help (); + } + + VERBOSE = FALSE; + if ((N = get_argument (*argc, argv, "-v"))) { + VERBOSE = TRUE; + remove_argument (N, argc, argv); + } + + MODE = MODE_IMAGEFILE; + if ((N = get_argument (*argc, argv, "-name"))) { + if (MODE != MODE_IMAGEFILE) usage(); + MODE = MODE_IMAGENAME; + remove_argument (N, argc, argv); + } + if ((N = get_argument (*argc, argv, "-orphan"))) { + if (MODE != MODE_IMAGEFILE) usage(); + MODE = MODE_ORPHAN; + remove_argument (N, argc, argv); + } + if ((N = get_argument (*argc, argv, "-missed"))) { + if (MODE != MODE_IMAGEFILE) usage(); + MODE = MODE_MISSED; + remove_argument (N, argc, argv); + } + if ((N = get_argument (*argc, argv, "-time"))) { + if (MODE != MODE_IMAGEFILE) usage(); + MODE = MODE_TIME; + remove_argument (N, argc, argv); + + if (!str_to_time (argv[N], &START)) usage (); + remove_argument (N, argc, argv); + + /* interpret second value */ + if (str_to_dtime (argv[N], &trange)) { + if (trange < 0) { + END = START; + START = END + trange; + } else { + END = START + trange; + } + remove_argument (N, argc, argv); + goto goodtime; + } + if (str_to_time (argv[N], &END)) { + if (START > END) { + tmp = START; + START = END; + END = START; + } + remove_argument (N, argc, argv); + goto goodtime; + } + usage (); + } +goodtime: + + /* restrict to a single photcode (not compatible with -image) */ + PHOTCODE = NULL; + if ((N = get_argument (*argc, argv, "-photcode"))) { + remove_argument (N, argc, argv); + PHOTCODE = GetPhotcodebyName (argv[N]); + remove_argument (N, argc, argv); + } + + if ((MODE == MODE_TIME) && (*argc != 1)) usage (); + if ((MODE != MODE_TIME) && (*argc != 2)) usage (); + return (TRUE); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/check_permissions.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/check_permissions.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/check_permissions.c (revision 4683) @@ -0,0 +1,67 @@ +# include "addstar.h" + +void check_permissions (char *basefile) { + + char *c, dir[256], filename[256]; + struct stat filestat; + uid_t uid; + gid_t gid; + int status, cmode; + + uid = getuid(); + gid = getgid(); + + /* check permission to write to directory */ + sprintf (filename, "%s", basefile); + c = strrchr (filename, '/'); + if (c == (char *) NULL) { + strcpy (dir, "."); + } else { + *c = 0; + strcpy (dir, filename); + } + status = stat (dir, &filestat); + if (status == -1) { + fprintf (stderr, "directory %s does not exist, creating...\n", dir); + cmode = S_IRWXU | S_IRWXG | S_IRWXO; + status = mkdir (dir, cmode); + if (status == -1) { + fprintf (stderr, "ERROR: can't create %s\n", dir); + exit (1); + } + } + status = stat (dir, &filestat); + if (((uid == filestat.st_uid) && (filestat.st_mode & S_IRWXU)) || + ((gid == filestat.st_gid) && (filestat.st_mode & S_IRWXG)) || + (filestat.st_mode & S_IRWXO)) { + } else { + fprintf (stderr, "ERROR: can't write to %s\n", dir); + exit (1); + } + + /* check permission to write to file */ + sprintf (filename, "%s", basefile); + status = stat (filename, &filestat); + if (status == 0) { /* file exists, are permissions OK? */ + if (((uid == filestat.st_uid) && (filestat.st_mode & S_IRUSR) && (filestat.st_mode & S_IWUSR)) || + ((gid == filestat.st_gid) && (filestat.st_mode & S_IRGRP) && (filestat.st_mode & S_IWGRP)) || + ((filestat.st_mode & S_IROTH) && (filestat.st_mode & S_IWOTH))) { + } else { + fprintf (stderr, "ERROR: can't write to %s\n", filename); + exit (1); + } + } + + /* check permission to write to backup file */ + sprintf (filename, "%s~", basefile); + status = stat (filename, &filestat); + if (status == 0) { /* file exists, are permissions OK? */ + if (((uid == filestat.st_uid) && (filestat.st_mode & S_IRUSR) && (filestat.st_mode & S_IWUSR)) || + ((gid == filestat.st_gid) && (filestat.st_mode & S_IRGRP) && (filestat.st_mode & S_IWGRP)) || + ((filestat.st_mode & S_IROTH) && (filestat.st_mode & S_IWOTH))) { + } else { + fprintf (stderr, "ERROR: can't write to %s\n", filename); + exit (1); + } + } +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/convert.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/convert.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/convert.c (revision 4683) @@ -0,0 +1,360 @@ +# include "delstar.h" + +/***** convert [-]00:00:00 to 0.0000 ****/ +int dms_to_ddd (double *Value, char *string) { + + int valid, neg, status; + double tmp, value; + char *p1, *p2, *px; + + valid = FALSE; + neg = FALSE; + stripwhite (string); + p1 = string; + px = string + strlen(string); + + if (string[0] == '-') { + valid = TRUE; + neg = TRUE; + p1 = &string[1]; + } + if (string[0] == '+') { + valid = TRUE; + neg = FALSE; + p1 = &string[1]; + } + if (isdigit(string[0])) { + valid = TRUE; + p1 = &string[0]; + } + if (!valid) { return (FALSE); } + + status = 1; + tmp = strtod (p1, &p2); + if (p2 == p1) return (FALSE); /* entry not a number: +fred */ + value = tmp; + if (p2 == px) goto escape; /* entry only number: +1.0 */ + p1 = p2 + 1; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; /* entry not a number: +1:fred */ + status = 2; + value += tmp / 60.0; + if (p2 == px) goto escape; /* entry only number: +1:1 */ + p1 = p2 + 1; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; /* entry not a number: +1:1:fred */ + value += tmp / 3600.0; + + escape: + if (neg) { + value *= -1; + } + *Value = value; + + return (status); +} + +/**********/ +int str_to_radec (double *ra, double *dec, char *str1, char *str2) { + + double Ra, Dec; + + *ra = *dec = 0; + switch (dms_to_ddd (&Ra, str1)) { + case 0: + fprintf (stderr, "syntax error in RA\n"); + return (FALSE); + case 1: + break; + case 2: + Ra = Ra * 15; + break; + } + switch (dms_to_ddd (&Dec, str2)) { + case 0: + fprintf (stderr, "syntax error in DEC\n"); + return (FALSE); + case 1: + case 2: + break; + } + *ra = Ra; + *dec = Dec; + return (TRUE); +} + +/**********/ +int chk_time (char *line) { + + char *p1, *p2; + double tmp; + int mode; + + p1 = line; + tmp = strtod (p1, &p2); + if (p2 == p1 + strlen (p1) - 1) { + if (*p2 == 'j') { + mode = TIME_JD; + } + if (*p2 == 'J') { + mode = TIME_MJD; + } + if (*p2 == 'd') { + mode = TIME_DAYS; + } + if (*p2 == 'h') { + mode = TIME_HOURS; + } + if (*p2 == 'm') { + mode = TIME_MINUTES; + } + if (*p2 == 's') { + mode = TIME_SECONDS; + } + } else { + mode = TIME_DATE; + } + return (mode); +} + +/**********/ +int str_to_time (char *line, time_t *second) { + + char *tmpline; + struct tm *gmt; + struct timeval now; + double jd; + time_t tsec; + + if (!strcasecmp (line, "NOW")) { + gettimeofday (&now, (struct timezone *) NULL); + *second = now.tv_sec; + return (TRUE); + } + + if (!strncasecmp (line, "TODAY", 5)) { + gettimeofday (&now, (struct timezone *) NULL); + if (line[5]) { /* line has extra data (ie, hh:mm:ss) */ + tsec = now.tv_sec; + ALLOCATE (tmpline, char, 64); + gmt = gmtime (&tsec); + sprintf (tmpline, "%04d/%02d/%02d,%s", + 1900 + gmt[0].tm_year, gmt[0].tm_mon+1, gmt[0].tm_mday, &line[6]); + *second = date_to_sec (tmpline); + free (tmpline); + return (TRUE); + } else { + *second = 86400 * ((int)(now.tv_sec / 86400)); + return (TRUE); + } + } + + switch (chk_time (line)) { + case 0: + return (FALSE); + case TIME_JD: + jd = strtod (line, 0); + *second = jd_to_sec (jd); + return (TRUE); + case TIME_MJD: + jd = strtod (line, 0); + *second = mjd_to_sec (jd); + return (TRUE); + case TIME_DATE: + *second = date_to_sec (line); + return (TRUE); + case TIME_DAYS: + *second = strtod (line, 0) * 86400.0; + return (TRUE); + case TIME_HOURS: + *second = strtod (line, 0) * 3600.0; + return (TRUE); + case TIME_MINUTES: + *second = strtod (line, 0) * 60.0; + return (TRUE); + case TIME_SECONDS: + *second = strtod (line, 0); + return (TRUE); + } + return (FALSE); +} + + +/**********/ +int str_to_dtime (char *line, double *second) { + + switch (chk_time (line)) { + case 0: + case TIME_JD: + case TIME_MJD: + case TIME_DATE: + return (FALSE); + case TIME_DAYS: + *second = strtod (line, 0) * 86400.0; + return (TRUE); + case TIME_HOURS: + *second = strtod (line, 0) * 3600.0; + return (TRUE); + case TIME_MINUTES: + *second = strtod (line, 0) * 60.0; + return (TRUE); + case TIME_SECONDS: + *second = strtod (line, 0); + return (TRUE); + } + return (FALSE); +} + +/**********/ +double sec_to_jd (time_t second) { + + double jd; + + jd = second/86400.0 + 2440587.5; + return (jd); +} + +/**********/ +time_t jd_to_sec (double jd) { + + time_t second; + + second = (jd - 2440587.5)*86400; + return (second); +} + +/**********/ +double sec_to_mjd (time_t second) { + + double mjd; + + mjd = second/86400.0 + 40587.0; + return (mjd); +} + +/**********/ +time_t mjd_to_sec (double mjd) { + + time_t second; + + second = (mjd - 40587.0)*86400; + return (second); +} + +/**********/ +char *sec_to_date (time_t second) { + + struct tm *gmt; + char *line; + + ALLOCATE (line, char, 64); + gmt = gmtime (&second); + sprintf (line, "%04d/%02d/%02d,%02d:%02d:%02d", + 1900 + gmt[0].tm_year, gmt[0].tm_mon+1, gmt[0].tm_mday, + gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); + return (line); + +} + +/***** date in format yyyy/mm/dd,hh:mm:ss *****/ +time_t date_to_sec (char *date) { + + time_t second; + double tmp, jd; + struct tm now; + char *p1, *p2, *px; + + p1 = date; + px = date + strlen(date); + bzero (&now, sizeof(now)); + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; + now.tm_year = tmp; + if (now.tm_year > 1000) now.tm_year -= 1900; + if (now.tm_year < 50) now.tm_year += 100; + if (p2 == px) goto escape; + p1 = p2 + 1; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; + now.tm_mon = tmp - 1; /* mon runs from 0 - 11 */ + if (p2 == px) goto escape; + p1 = p2 + 1; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; + now.tm_mday = tmp; + if (p2 == px) goto escape; + p1 = p2 + 1; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; + p1 = p2 + 1; + now.tm_hour = tmp; + if (p2 == px) goto escape; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; + now.tm_min = tmp; + if (p2 == px) goto escape; + p1 = p2 + 1; + + tmp = strtod (p1, &p2); + if (p2 == p1) goto escape; + now.tm_sec = tmp; + if (p2 == px) goto escape; + p1 = p2 + 1; + + escape: + jd = now.tm_mday - 32075 + (int)(1461*(1900 + now.tm_year + 4800 + (int)(((now.tm_mon+1)-14)/12))/4) + + (int)(367*((now.tm_mon+1) - 2 - (int)(((now.tm_mon+1) - 14)/12)*12)/12) + - (int)(3*(int)((1900 + now.tm_year + 4900 + (int)(((now.tm_mon+1) - 14)/12))/100)/4) - 0.5; + + second = (jd - 2440587.5)*86400 + 3600.0*now.tm_hour + now.tm_min*60.0 + now.tm_sec; + + return (second); +} + +/**********/ +int hms_to_deg (double *h0, double *h1, double *d0, double *d1, char *string) { + + int flag_d0, flag_d1, flag_h0, flag_h1; + double tmp; + + *d0 = *h0 = *d1 = *h1 = 0; + + flag_h0 = dparse (h0, 1, string); + flag_h1 = dparse (h1, 4, string); + flag_d0 = dparse (d0, 7, string); + flag_d1 = dparse (d1, 9, string); + *h0 *= flag_h0; + *h1 *= flag_h1; + *d0 *= flag_d0; + *d1 *= flag_d1; + + dparse (&tmp, 2, string); + *h0 += tmp/60.0; + dparse (&tmp, 3, string); + *h0 += tmp/3600.0; + + dparse (&tmp, 5, string); + *h1 += tmp/60.0; + dparse (&tmp, 6, string); + *h1 += tmp/3600.0; + + dparse (&tmp, 8, string); + *d0 += tmp/60.0; + + dparse (&tmp, 10, string); + *d1 += tmp/60.0; + + *h0 *= 15*flag_h0; + *h1 *= 15*flag_h1; + *d0 *= flag_d0; + *d1 *= flag_d1; + + return (TRUE); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/delete_image.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/delete_image.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/delete_image.c (revision 4683) @@ -0,0 +1,43 @@ + +delete_imagefile (char *filename) { + + int i, Nregions; + Image *image; + GSCRegion *region; + Catalog catalog; + + image = gimages (filename); + + fprintf (stderr, "deleting %s\n", image[0].name); + region = gregions (image, &Nregions); + + for (i = 0; i < Nregions; i++) { + fprintf (stderr, "deleting from %s\n", region[i].filename); + catalog.filename = region[i].filename; /* don't free region before catalog! */ + switch (lock_catalog (&catalog, LCK_SOFT)) { + case 0: + fprintf (stderr, "ERROR: can't lock file\n"); + exit (1); + case 1: + gcatalog (&catalog); + break; + case 2: + fprintf (stderr, "no data in %s, skipping\n", catalog.filename); + catalog.Naverage = 0; + catalog.Nmeasure = 0; + break; + default: + fprintf (stderr, "weird lock_catalog exit state\n"); + exit (1); + } + find_matches (&catalog, image); + wcatalog (&catalog); + } + + dimages (f, image, 1); + + /* this probably gets moved inside dimages */ + unlock_image_db (ImageCat, f, LCK_HARD, &dbstate); + fprintf (stderr, "SUCCESS\n"); + exit (0); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/delete_imagefile.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/delete_imagefile.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/delete_imagefile.c (revision 4683) @@ -0,0 +1,58 @@ +# include "delstar.h" + +void delete_imagefile (char *filename) { + + int i, j; + int Nregions, Nimlist; + int *imlist; + double trange; + time_t start, stop; + Image *image; + GSCRegion *region; + Catalog catalog; + + /* load information about file - time/photcode */ + image = gimages (filename); + + if (VERBOSE) fprintf (stderr, "deleting %s\n", image[0].name); + region = gregion_image (image, &Nregions); + + for (i = 0; i < Nregions; i++) { + if (VERBOSE) fprintf (stderr, "deleting from %s\n", region[i].filename); + catalog.filename = region[i].filename; /* don't free region before catalog! */ + switch (lock_catalog (&catalog, LCK_XCLD)) { + case 0: + fprintf (stderr, "ERROR: can't lock file %s\n", catalog.filename); + exit (1); + case 1: + gcatalog (&catalog); + break; + case 2: + if (VERBOSE) fprintf (stderr, "no data in %s, skipping\n", catalog.filename); + unlock_catalog (&catalog); + continue; + default: + fprintf (stderr, "weird lock_catalog exit state\n"); + exit (1); + } + /* pad exposure time by 1 sec to require a valid time */ + /* trate is in 0.1 msec / row - stop is the latest exposure end time */ + trange = 1e-4*image[0].NY*image[0].trate + image[0].exptime + 1; + start = image[0].tzero; + stop = image[0].tzero + trange; + find_matches (&catalog, image[j].source, start, stop); + wcatalog (&catalog); + unlock_catalog (&catalog); + } + + /* find and delete matching images */ + imlist = find_images_data (image, &Nimlist); + if (!Nimlist) Shutdown ("image %s not found in db", filename); + + dimages (imlist, Nimlist); + + save_image_db (); + unlock_image_db (); + fprintf (stderr, "SUCCESS\n"); + exit (0); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/delete_imagename.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/delete_imagename.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/delete_imagename.c (revision 4683) @@ -0,0 +1,63 @@ +# include "delstar.h" + +void delete_imagename (char *name) { + + int i, j, k; + int Nregions, Nimlist, Nimage; + int *imlist; + double trange; + time_t start, stop; + Image *image; + GSCRegion *region; + Catalog catalog; + + image = GetImages (&Nimage); + + /* find image in db by name */ + imlist = find_images_name (name, &Nimlist); + if (!Nimlist) Shutdown ("image %s not found in db", name); + + for (k = 0; k < Nimlist; k++) { + + j = imlist[k]; + if (VERBOSE) fprintf (stderr, "deleting %s\n", image[j].name); + region = gregion_image (&image[j], &Nregions); + + for (i = 0; i < Nregions; i++) { + if (VERBOSE) fprintf (stderr, "deleting from %s\n", region[i].filename); + catalog.filename = region[i].filename; /* don't free region before catalog! */ + switch (lock_catalog (&catalog, LCK_XCLD)) { + case 0: + fprintf (stderr, "ERROR: can't lock file %s\n", catalog.filename); + exit (1); + case 1: + gcatalog (&catalog); + break; + case 2: + if (VERBOSE) fprintf (stderr, "no data in %s, skipping\n", catalog.filename); + unlock_catalog (&catalog); + continue; + default: + fprintf (stderr, "weird lock_catalog exit state\n"); + exit (1); + } + + /* trate is in 0.1 msec / row - stop is the latest exposure end time */ + /* pad exposure time by 1 sec to require a valid time */ + trange = 1e-4*image[j].NY*image[j].trate + image[j].exptime + 1; + start = image[j].tzero; + stop = image[j].tzero + trange; + find_matches (&catalog, image[j].source, start, stop); + wcatalog (&catalog); + unlock_catalog (&catalog); + } + } + + /* delete the identified images */ + dimages (imlist, Nimlist); + + save_image_db (); + unlock_image_db (); + fprintf (stderr, "SUCCESS\n"); + exit (0); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/delete_missed.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/delete_missed.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/delete_missed.c (revision 4683) @@ -0,0 +1,24 @@ +# include "delstar.h" + +/* drop all MISSED values for the given catalog */ + +delete_missed (Catalog *catalog) { + + int i; + int Nave, Nmeas, Nmiss; + + Nave = catalog[0].Naverage; + Nmeas = catalog[0].Nmeasure; + Nmiss = catalog[0].Nmissing; + + if (VERBOSE) fprintf (stderr, "starting with Nave, Nmeas, Nmiss: %d %d %d\n", Nave, Nmeas, Nmiss); + + /* set up references for missing to average */ + for (i = 0; i < Nave; i++) { + catalog[0].average[i].Nn = 0; + } + REALLOCATE (catalog[0].missing, Missing, 1); + catalog[0].Nmissing = 0; + if (VERBOSE) fprintf (stderr, " ending with Nave, Nmeas, Nmiss: %d %d %d\n", Nave, Nmeas, Nmiss); +} + Index: /tags/delstar-1-2/Ohana/src/delstar/src/delete_orphans.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/delete_orphans.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/delete_orphans.c (revision 4683) @@ -0,0 +1,256 @@ +# include "delstar.h" + +delete_orphans (char *name) { + + int i, j, k, n, m, N, M, found; + int *N1, *N2, *next, *next_miss, *ave_miss, last, last_miss; + int Nave, NAVE, Nmeas, NMEAS, Nmiss, NMISS, Nmatch; + int start, end, Nmeasfound, Nsecfilt; + unsigned int flags; + Measure *tmpmeasure; + Missing *tmpmissing; + Coords tcoords; + Catalog catalog; + + /* find and load catalog file */ + catalog.filename = name; + gcatalog (&catalog); + gcatstats (&catalog, &catstats); + + /* find images overlapping catalog */ + image = find_images_region (&catstats, &Nimage); + match_images (&catalog, image, Nimage); + + /** allocate local arrays **/ + Nave = catalog.Naverage; + + Nmeas = catalog.Nmeasure; + ALLOCATE (next, int, Nmeas); + + Nmiss = catalog.Nmissing; + ALLOCATE (next_miss, int, Nmiss); + ALLOCATE (ave_miss, int, Nmiss); + + if (VERBOSE) fprintf (stderr, "starting with Nave, Nmeas, Nmiss: %d %d %d\n", Nave, Nmeas, Nmiss); + + /* set up pointers for linked list of measure */ + for (i = 0; i < Nmeas - 1; i++) { + next[i] = i+1; + } + next[i] = -1; + last = i; + /* set up pointers for linked list of missing */ + for (i = 0; i < Nmiss - 1; i++) { + next_miss[i] = i+1; + } + next_miss[i] = -1; + last_miss = i; + /* set up references for missing to average */ + for (i = 0; i < Nave; i++) { + for (j = 0; j < catalog.average[i].Nn; j++) { + ave_miss[catalog.average[i].missing + j] = i; + } + } + Nmeasfound = 0; + Nsecfilt = catalog.Nsecfilt; + + /* fprintf (stderr, "fixing the measures...\n"); */ + for (i = 0; (i < Nmeas); i++) { + if ((catalog.measure[i].t != 0) && (catalog.image[i] == -1)) { + /* this star is an orphan */ + Nmeasfound ++; + next[i] = -2; /* we delete this one */ + /* fix the list links: connect the previous valid link to the next valid link */ + for (j = i; (j >= 0) && (next[j] == -2); j--); /* find previous entry to fix link */ + if (j >= 0) { /* if j < 0, there is no previous valid link, ignore this step */ + if (next[j] != i) { + fprintf (stderr, "error? this link seems to have been lost\n"); + exit (1); + } + /* find next valid entry to fix link */ + for (k = i; (k < Nmeas) && (next[k] == -2); k++); + if (k < Nmeas) + next[j] = k; + else + next[j] = -1; /* last link in list gets a -1 */ + } + + /*** fix the corresponding average entry ***/ + n = catalog.measure[i].averef; + if (catalog.average[n].Nm == 0) { /* this should never happen */ + fprintf (stderr, "error? we deleted one too many objects?\n"); + exit (1); + } + catalog.average[n].Nm --; + /* this was only entry in list: will be deleted below. meanwhile, delete all missing entries*/ + if ((catalog.average[n].Nm < 1) && (catalog.average[n].Nn > 0)) { + m = catalog.average[n].missing; + for (j = 0; j < catalog.average[n].Nn; j++) { + M = next_miss[m]; + next_miss[m] = -2; + m = M; + } + m = catalog.average[n].missing; + /* fix the list links: connect the previous valid link to the next valid link */ + for (j = m; (j >= 0) && (next_miss[j] == -2); j--); /* find previous entry to fix link */ + if (j >= 0) { /* if j < 0, there is no previous valid link, ignore this step */ + if (next_miss[j] != m) { + fprintf (stderr, "error? this link seems to have been lost\n"); + exit (1); + } + /* find next valid entry to fix link */ + for (k = m; (k < Nmiss) && (next_miss[k] == -2); k++); + if (k < Nmiss) + next_miss[j] = k; + else + next_miss[j] = -1; /* last link in list gets a -1 */ + } + } + /* this was first entry in list */ + if ((catalog.average[n].offset == i) && (catalog.average[n].Nm > 0)) { + m = catalog.average[n].offset; + /* find next valid entry -- notice lack of error checking... */ + for (j = 0; (j < Nmeas) && (next[m+j] == -2); j++); + catalog.average[n].offset = m + j; + } + + } + } + fprintf (stderr, "found %d meas to remove\n", Nmeasfound); + + /* fprintf (stderr, "fixing the missing...\n"); */ + /** find missing in time range of image **/ + for (i = 0; (i < Nmiss); i++) { + if ((next_miss[i] != -2) && (catalog.missing[i].t >= start) && (catalog.missing[i].t <= end)) { + /* this star is in this image */ + + next_miss[i] = -2; /* we delete this one */ + /* fix the list links: connect the previous valid link to the next valid link */ + for (j = i; (j >= 0) && (next_miss[j] == -2); j--); /* find previous entry to fix link */ + if (j >= 0) { /* if j < 0, there is no previous valid link, ignore this step */ + if (next_miss[j] != i) { + fprintf (stderr, "error? this link seems to have been lost\n"); + exit (1); + } + /* find next valid entry to fix link */ + for (k = i; (k < Nmiss) && (next_miss[k] == -2); k++); + if (k < Nmiss) + next_miss[j] = k; + else + next_miss[j] = -1; /* last link in list gets a -1 */ + } + + /* find the corresponding avearge entry + found = FALSE; + for (n = 0; !found && (n < Nave); n++) { + if ((catalog.average[n].missing > 0) && (catalog.average[n].missing <= i) && (catalog.average[n].missing + catalog.average[n].Nn > i)) + found = TRUE; + } + n--; */ + /*** fix the corresponding average entry ***/ + n = ave_miss[i]; + if (catalog.average[n].Nn == 0) { /* this should never happen */ + fprintf (stderr, "error? we deleted one too many missing?\n"); + exit (1); + } + catalog.average[n].Nn --; + /* this was first entry in list */ + if ((catalog.average[n].missing == i) && (catalog.average[n].Nn > 0)) { + m = catalog.average[n].missing; + for (j = 0; (j < Nmiss) && (next_miss[m+j] == -2); j++); + catalog.average[n].missing = m + j; + } + + } + } + + /* currently not worked out, but we will need to delete the references to blended image and cat stars */ +# if 0 + /*** handle multiple stars */ + /* this image star matches more than one catalog star */ + if (stars[N].found > -1) { + catalog.measure[stars[N].found].flags |= BLEND_IMAGE; + catalog.measure[Nmeas].flags |= BLEND_IMAGE; + } + if (stars[N].found == -2) { /* this image star matches a catalog star on a neighboring catalog */ + catalog.measure[Nmeas].flags |= BLEND_IMAGE_NEIGHBOR; + } + if (stars[N].found == -1) { /* this image star matches only this star */ + stars[N].found = Nmeas; /* save first match, in case coincidences are found */ + } + /* this catalog star matches more than one image star */ + if (catalog.found[n] > -1) { + catalog.measure[catalog.found[n]].flags |= BLEND_CATALOG; + catalog.measure[Nmeas].flags |= BLEND_CATALOG; + } else { + catalog.found[n] = Nmeas; + } +# endif + + /* fprintf (stderr, "fixing the averages...\n"); */ + /* fix Average list: delete entries with Nm == 0 */ + for (i = j = 0; (i < Nave) && (j < Nave); i++, j++) { + for (; (j < Nave) && (catalog.average[j].Nm == 0); j++); + if ((i != j) && (j < Nave)) { + catalog.average[i] = catalog.average[j]; + for (k = 0; k < catalog.Nsecfilt; k++) { + catalog.secfilt[i*Nsecfilt + k] = catalog.secfilt[j*Nsecfilt + k]; + } + } + if (j == Nave) i--; + } + Nave = i; + REALLOCATE (catalog.average, Average, Nave); + REALLOCATE (catalog.secfilt, SecFilt, MAX (1, Nave*Nsecfilt)); + + /* fprintf (stderr, "fixing the measure order...\n"); */ + /* fix order of Measure (memory intensive, but fast) */ + N = 0; + ALLOCATE (tmpmeasure, Measure, Nmeas); + for (i = 0; i < Nave; i++) { + n = catalog.average[i].offset; + catalog.average[i].offset = N; + for (k = 0; k < catalog.average[i].Nm; k++, N++) { + if ((n == -1) || (n == -2)) { + fprintf (stderr, "error: linked list is confused\n"); + exit (1); + } + tmpmeasure[N] = catalog.measure[n]; + tmpmeasure[N].averef = i; + n = next[n]; + } + } + Nmeas = N; + free (catalog.measure); + catalog.measure = tmpmeasure; + REALLOCATE (catalog.measure, Measure, Nmeas); + + /* fprintf (stderr, "fixing the mising order...\n"); */ + /* fix order of Missing (memory intensive, but fast) */ + N = 0; + ALLOCATE (tmpmissing, Missing, Nmiss); + for (i = 0; i < Nave; i++) { + if (catalog.average[i].Nn > 0) { + n = catalog.average[i].missing; + catalog.average[i].missing = N; + for (k = 0; k < catalog.average[i].Nn; k++, N++) { + if ((n == -1) || (n == -2)) { + fprintf (stderr, "error: linked list is confused\n"); + exit (1); + } + tmpmissing[N] = catalog.missing[n]; + n = next_miss[n]; + } + } + } + Nmiss = N; + free (catalog.missing); + catalog.missing = tmpmissing; + REALLOCATE (catalog.missing, Missing, Nmiss); + + catalog.Naverage = Nave; + catalog.Nmeasure = Nmeas; + catalog.Nmissing = Nmiss; + if (VERBOSE) fprintf (stderr, " ending with Nave, Nmeas, Nmiss: %d %d %d\n", Nave, Nmeas, Nmiss); + +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/delete_times.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/delete_times.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/delete_times.c (revision 4683) @@ -0,0 +1,83 @@ +# include "delstar.h" + +void delete_times () { + + int i, j, k, n; + int Nimage, Nimlist, found, code; + int Ntregions, Nregions, NREGIONS; + int *imlist; + GSCRegion *region, *tregion; + Image *image; + Catalog catalog; + + code = (PHOTCODE == NULL) ? -1 : PHOTCODE[0].code; + + Nregions = 0; + NREGIONS = 10; + ALLOCATE (region, GSCRegion, NREGIONS); + + image = GetImages (&Nimage); + + /* find images for time range, delete each image */ + imlist = find_images_time (START, END, PHOTCODE, &Nimlist); + + /* find all overlapping regions */ + for (n = 0; n < Nimlist; n++) { + j = imlist[n]; + if (VERBOSE) fprintf (stderr, "finding regions for %s\n", image[j].name); + + tregion = gregion_image (&image[j], &Ntregions); + for (i = 0; i < Ntregions; i++) { + found = FALSE; + for (k = 0; (k < Nregions) && !found; k++) { + found = !strcmp (region[k].filename, tregion[i].filename); + } + if (!found) { + strcpy (region[Nregions].filename, tregion[i].filename); + region[Nregions].RA[0] = tregion[i].RA[0]; + region[Nregions].RA[1] = tregion[i].RA[1]; + region[Nregions].DEC[0] = tregion[i].DEC[0]; + region[Nregions].DEC[1] = tregion[i].DEC[1]; + Nregions ++; + if (Nregions >= NREGIONS) { + NREGIONS += 10; + REALLOCATE (region, GSCRegion, NREGIONS); + } + } + } + free (tregion); + } + + /* delete from all identified regions */ + for (i = 0; i < Nregions; i++) { + if (VERBOSE) fprintf (stderr, "deleting from %s\n", region[i].filename); + catalog.filename = region[i].filename; /* don't free region before catalog! */ + switch (lock_catalog (&catalog, LCK_XCLD)) { + case 0: + fprintf (stderr, "ERROR: can't lock file\n"); + exit (1); + case 1: + gcatalog (&catalog); + break; + case 2: + fprintf (stderr, "no data in %s, skipping\n", catalog.filename); + unlock_catalog (&catalog); + continue; + default: + fprintf (stderr, "weird lock_catalog exit state\n"); + exit (1); + } + + find_matches (&catalog, code, START, END); + wcatalog (&catalog); + unlock_catalog (&catalog); + } + + /* delete the identified images */ + dimages (imlist, Nimlist); + + save_image_db (); + unlock_image_db (); + fprintf (stderr, "SUCCESS\n"); + exit (0); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/delstar.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/delstar.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/delstar.c (revision 4683) @@ -0,0 +1,35 @@ +# include "delstar.h" + +int main (int argc, char **argv) { + + SetSignals (); + ConfigInit (&argc, argv); + args (&argc, argv); + + lock_image_db (); + load_image_db (); + + /* delete orphaned measurements (no image) from catalog */ + switch (MODE) { + case MODE_IMAGEFILE: + delete_imagefile (argv[1]); + break; + case MODE_IMAGENAME: + delete_imagename (argv[1]); + break; + case MODE_TIME: + delete_times (); + break; + case MODE_ORPHAN: + fprintf (stderr, "delete orphans not available\n"); + // delete_orphans (argv[1]); + break; + case MODE_MISSED: + fprintf (stderr, "delete missed not available\n"); + // delete_missed (argv[1]); + break; + default: + usage (); + } + exit (1); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/dimages.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/dimages.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/dimages.c (revision 4683) @@ -0,0 +1,32 @@ +# include "delstar.h" + +/* delete images by sequence number */ + +int dimages (int *imlist, int Nimlist) { + + int i, j, Nimage; + Image *image; + + image = GetImages (&Nimage); + + /* first, mark the deleted images with a flag (name = 0) */ + for (i = 0; i < Nimlist; i++) { + j = imlist[i]; + image[j].name[0] = 0; + /* this puts a restriction on the DB: image name must not be NULL */ + } + + /* delete the marked images */ + for (i = 0, j = 0; j < Nimage; j++) { + if (!image[j].name[0]) continue; + if (i != j) { + image[i] = image[j]; + } + i++; + } + Nimage = i; + REALLOCATE (image, Image, Nimage); + + SetImages (image, Nimage); + return (TRUE); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/find_image_db.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/find_image_db.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/find_image_db.c (revision 4683) @@ -0,0 +1,94 @@ +# include "delstar.h" + +int *find_images_name (char *filename, int *nlist) { + + int i, Nimage, Nlist, NLIST; + int *list; + char *p, *name; + Image *image; + + /* strip off all but the filename */ + p = strrchr (filename, '/'); + if (p == NULL) { + name = filename; + } else { + name = p + 1; + } + + image = GetImages (&Nimage); + + Nlist = 0; + NLIST = 100; + ALLOCATE (list, int, NLIST); + + for (i = 0; i < Nimage; i++) { + if (strcmp (image[i].name, name)) continue; + list[Nlist] = i; + Nlist ++; + CHECK_REALLOCATE (list, int, NLIST, Nlist, 100); + } + + *nlist = Nlist; + return (list); +} + +/* find images in db by image data (time/photcode) */ +int *find_images_data (Image *timage, int *nlist) { + + int i, Nimage, Nlist, NLIST; + int *list; + Image *image; + time_t start, stop; + int code; + + start = timage[0].tzero - MAX(0.05*timage[0].trate*timage[0].NY, 1); + stop = timage[0].tzero + MAX(1.05*timage[0].trate*timage[0].NY, 1); + code = timage[0].source; + + image = GetImages (&Nimage); + + Nlist = 0; + NLIST = 100; + ALLOCATE (list, int, NLIST); + + for (i = 0; i < Nimage; i++) { + if (image[i].tzero < start) continue; + if (image[i].tzero > stop) continue; + if (image[i].source != code) continue; + list[Nlist] = i; + Nlist ++; + CHECK_REALLOCATE (list, int, NLIST, Nlist, 100); + } + + *nlist = Nlist; + return (list); +} + +/* find images in db by image data (time/photcode) */ +int *find_images_time (time_t start, time_t end, PhotCode *code, int *nlist) { + + int i, Nimage, Nlist, NLIST; + int *list; + Image *image; + + image = GetImages (&Nimage); + + Nlist = 0; + NLIST = 100; + ALLOCATE (list, int, NLIST); + + for (i = 0; i < Nimage; i++) { + if (image[i].tzero < START) continue; + if (image[i].tzero > END) continue; + if (code != NULL) { + if (image[i].source != code[0].code) continue; + } + list[Nlist] = i; + Nlist ++; + CHECK_REALLOCATE (list, int, NLIST, Nlist, 100); + } + + *nlist = Nlist; + return (list); +} + Index: /tags/delstar-1-2/Ohana/src/delstar/src/find_images.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/find_images.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/find_images.c (revision 4683) @@ -0,0 +1,185 @@ +# include "delstar.h" + +Image *find_images (FILE *f, CatStats *catstats, int *Nimages) { + + Header header; + Image *timage, *image; + int i, j, k, found, nimage, Nimage, NIMAGE, NTIMAGE, Nloop, Nlast; + int n, Nim, status; + double Xc[6], Yc[6], Xi[6], Yi[6], r, d, x, y, dx, dy; + Coords *tcoords; + + /* we make positional comparisons in the projection of catalog */ + tcoords = &catstats[0].coords; + /* define catalog corners */ + Xc[0] = catstats[0].RA[0]; Yc[0] = catstats[0].DEC[0]; + Xc[1] = catstats[0].RA[1]; Yc[1] = catstats[0].DEC[0]; + Xc[2] = catstats[0].RA[1]; Yc[2] = catstats[0].DEC[1]; + Xc[3] = catstats[0].RA[0]; Yc[3] = catstats[0].DEC[1]; + Xc[4] = catstats[0].RA[0]; Yc[4] = catstats[0].DEC[0]; + Xc[5] = catstats[0].RA[1]; Yc[5] = catstats[0].DEC[1]; + for (j = 0; j < 6; j++) { + r = Xc[j]; d = Yc[j]; + RD_to_XY (&Xc[j], &Yc[j], r, d, tcoords); + } + /* find Y positions of RA center */ + r = 0.5*(catstats[0].RA[0] + catstats[0].RA[1]); + d = catstats[0].DEC[0]; + RD_to_XY (&x, &y, r, d, tcoords); + Yc[0] = MIN (y, Yc[0]); + Yc[1] = MIN (y, Yc[1]); + Yc[4] = MIN (y, Yc[4]); + /* find Y positions of RA center */ + r = 0.5*(catstats[0].RA[0] + catstats[0].RA[1]); + d = catstats[0].DEC[1]; + RD_to_XY (&x, &y, r, d, tcoords); + Yc[2] = MAX (y, Yc[2]); + Yc[3] = MAX (y, Yc[3]); + Yc[5] = MAX (y, Yc[5]); + + dx = 0.02*(Xc[2] - Xc[0]); + dy = 0.02*(Yc[2] - Yc[0]); + Xc[0] -= dx; Yc[0] -= dy; + Xc[1] += dx; Yc[1] -= dy; + Xc[2] += dx; Yc[2] += dy; + Xc[3] -= dx; Yc[3] += dy; + Xc[4] -= dx; Yc[4] -= dy; + Xc[5] -= dx; Yc[5] -= dy; + + /* check if image datafile exists, get header, number of images */ + if (!fits_fread_header (f, &header)) { + fprintf (stderr, "ERROR: error reading image catalog header\n"); + exit (1); + } + Nimage = 0; + fits_scan (&header, "NIMAGES", "%d", 1, &Nimage); + if (Nimage == 0) { + fprintf (stderr, "ERROR: No images in catalog %s (1)\n", ImageCat); + Shutdown (); + } + fseek (f, header.size, SEEK_SET); + + /* set up buffers for images, temporary storage */ + NTIMAGE = 100; + ALLOCATE (timage, Image, NTIMAGE); + NIMAGE = 100; + ALLOCATE (image, Image, NIMAGE); + nimage = 0; + Nloop = Nimage / NTIMAGE + 1; + Nlast = Nimage % NTIMAGE; + + /* read in images in groups of NTIMAGE (100) */ + for (n = 0; n < Nloop; n++) { + Nim = (n == Nloop - 1) ? Nlast : NTIMAGE; + status = Fread (timage, sizeof(Image), Nim, f, "image"); + if (status != Nim) { + fprintf (stderr, "ERROR: couldn't read images from image catalog: %s\n", ImageCat); + exit (1); + } + /* test each image in block */ + for (i = 0; i < Nim; i++) { + /* define image corners */ + Xi[0] = 0; Yi[0] = 0; + Xi[1] = timage[i].NX; Yi[1] = 0; + Xi[2] = timage[i].NX; Yi[2] = timage[i].NY; + Xi[3] = 0; Yi[3] = timage[i].NY; + Xi[4] = 0; Yi[4] = 0; + Xi[5] = timage[i].NX; Yi[5] = timage[i].NY; + found = FALSE; + /* transform to tcoords */ + if (catstats[0].DEC[1] > 86.25) { /* pole */ + for (j = 0; j < 6; j++) { + XY_to_RD (&r, &d, Xi[j], Yi[j], &timage[i].coords); + if (d > catstats[0].DEC[0] - 0.5) found = TRUE; + } + } else { + for (j = 0; j < 6; j++) { + XY_to_RD (&r, &d, Xi[j], Yi[j], &timage[i].coords); + RD_to_XY (&Xi[j], &Yi[j], r, d, tcoords); + } + /* check if edges cross */ + for (j = 0; (j < 5) && !found; j++) { + for (k = 0; (k < 5) && !found; k++) { + found |= edge_check (&Xi[j], &Yi[j], &Xc[k], &Yc[k]); + } + } + } + if (found) { + image[nimage] = timage[i]; + image[nimage].code = 0; + nimage ++; + if (nimage == NIMAGE) { + NIMAGE += 100; + REALLOCATE (image, Image, NIMAGE); + } + } + } + } + + if (VERBOSE) { + for (i = 0; i < nimage; i++) { + XY_to_RD (&r, &d, 0.5*image[i].NX, 0.5*image[i].NY, &image[i].coords); + fprintf (stderr, "associated images: %d %8.4f %8.4f %10d %6d %5.3f %6.3f %6.3f\n", + i, r, d, image[i].tzero, image[i].nstar, 0.001*image[i].secz, + 0.001*image[i].Mcal, 0.001*image[i].dMcal); + } + } + + REALLOCATE (image, Image, MAX (nimage, 1)); + free (timage); + *Nimages = nimage; + return (image); +} + +int edge_check (double *x1, double *y1, double *x2, double *y2) { + + double theta1, theta2; + double Theta1, Theta2; + + theta1 = opening_angle (x1[0], y1[0], x2[0], y2[0], x1[1], y1[1]); + theta2 = opening_angle (x1[0], y1[0], x2[0], y2[0], x2[1], y2[1]); + + if (theta1*theta2 < 0.0) { + return (FALSE); + } + + if (fabs(theta1) < fabs(theta2)) { + return (FALSE); + } + + Theta1 = theta1; + Theta2 = theta2; + theta1 = opening_angle (x2[0], y2[0], x1[1], y1[1], x2[1], y2[1]); + theta2 = opening_angle (x2[0], y2[0], x1[1], y1[1], x1[0], y1[0]); + + + if (theta1*theta2 < 0.0) { + return (FALSE); + } + + if (fabs(theta1) < fabs(theta2)) { + return (FALSE); + } + + return (TRUE); + +} + +double opening_angle (double x1, double y1, double x2, double y2, double x3, double y3) { + + double dx1, dy1, dx2, dy2, ct, st, theta; + + dx1 = x1 - x2; + dy1 = y1 - y2; + + dx2 = x3 - x2; + dy2 = y3 - y2; + + ct = (dx1*dx2 + dy1*dy2); + st = (dx1*dy2 - dx2*dy1); + + theta = atan2 (st, ct); + + return (theta); + +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/find_matches.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/find_matches.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/find_matches.c (revision 4683) @@ -0,0 +1,243 @@ +# include "delstar.h" + +void find_matches (Catalog *catalog, int photcode, int start, int end) { + + int i, j, k, n, m, N, M, drop; + int *next, *next_miss, *ave_miss, last, last_miss; + int Nave, Nmeas, Nmiss; + int Nmeasfound, Nsecfilt; + int this, prev; + Measure *tmpmeasure; + Missing *tmpmissing; + + /** allocate local arrays **/ + Nave = catalog[0].Naverage; + + Nmeas = catalog[0].Nmeasure; + ALLOCATE (next, int, MAX(Nmeas,1)); + + Nmiss = catalog[0].Nmissing; + ALLOCATE (next_miss, int, MAX(Nmiss,1)); + ALLOCATE (ave_miss, int, MAX(Nmiss,1)); + + if (VERBOSE) fprintf (stderr, "starting with Nave, Nmeas, Nmiss: %d %d %d\n", Nave, Nmeas, Nmiss); + + /* set up pointers for linked list of measure */ + for (i = 0; i < Nmeas - 1; i++) { + next[i] = i+1; + } + next[i] = -1; + last = i; + /* set up pointers for linked list of missing */ + for (i = 0; i < Nmiss - 1; i++) { + next_miss[i] = i+1; + } + next_miss[i] = -1; + last_miss = i; + /* set up references for missing to average */ + for (i = 0; i < Nave; i++) { + for (j = 0; j < catalog[0].average[i].Nn; j++) { + ave_miss[catalog[0].average[i].missing + j] = i; + } + } + + if (VERBOSE) fprintf (stderr, "deleting for range %d to %d, photcode %d\n", start, end, photcode); + Nmeasfound = 0; + Nsecfilt = catalog[0].Nsecfilt; + + if (VERBOSE) fprintf (stderr, "fixing measure...\n"); + + /** find measure in time range **/ + this = prev = -1; + for (i = 0; (i < Nmeas); i++) { + if (VERBOSE && !(i % 10000)) fprintf (stderr, ". "); + drop = TRUE; + drop &= (catalog[0].measure[i].t >= start); + drop &= (catalog[0].measure[i].t <= end); + drop &= ((photcode == -1) || (photcode == catalog[0].measure[i].source)); + if (!drop) { + prev = i; + continue; + } + Nmeasfound ++; + + /* this star is in this image */ + this = next[i]; + next[i] = -2; /* we delete this one */ + if (prev != -1) { next[prev] = this; } + +# if (0) + /* why is this section disabled? */ + /* fix the list links: connect the previous valid link to the next valid link */ + for (j = i; (j >= 0) && (next[j] == -2); j--); /* find previous entry to fix link */ + if (j >= 0) { /* if j < 0, there is no previous valid link, ignore this step */ + if (next[j] != i) { + fprintf (stderr, "error? (1) this link seems to have been lost\n"); + fprintf (stderr, "j: %d, next[j]: %d, i: %d\n", j, next[j], i); + exit (1); + } + /* find next valid entry to fix link */ + for (k = i; (k < Nmeas) && (next[k] == -2); k++); + if (k < Nmeas) + next[j] = k; + else + next[j] = -1; /* last link in list gets a -1 */ + } +# endif + + /*** fix the corresponding average entry ***/ + n = catalog[0].measure[i].averef; + if (catalog[0].average[n].Nm == 0) { /* this should never happen */ + fprintf (stderr, "error? we deleted one too many objects?\n"); + exit (1); + } + catalog[0].average[n].Nm --; + /* this was only entry in list: will be deleted below. meanwhile, delete all missing entries*/ + if ((catalog[0].average[n].Nm < 1) && (catalog[0].average[n].Nn > 0)) { + m = catalog[0].average[n].missing; + for (j = 0; j < catalog[0].average[n].Nn; j++) { + M = next_miss[m]; + next_miss[m] = -2; + m = M; + } + m = catalog[0].average[n].missing; + /* fix the list links: connect the previous valid link to the next valid link */ + for (j = m; (j >= 0) && (next_miss[j] == -2); j--); /* find previous entry to fix link */ + if (j >= 0) { /* if j < 0, there is no previous valid link, ignore this step */ + if (next_miss[j] != m) { + fprintf (stderr, "error? (2) this link seems to have been lost\n"); + fprintf (stderr, "j: %d, next_miss[j]: %d, i: %d\n", j, next_miss[j], i); + exit (1); + } + /* find next valid entry to fix link */ + for (k = m; (k < Nmiss) && (next_miss[k] == -2); k++); + if (k < Nmiss) + next_miss[j] = k; + else + next_miss[j] = -1; /* last link in list gets a -1 */ + } + } + /* this was first entry in list */ + if ((catalog[0].average[n].offset == i) && (catalog[0].average[n].Nm > 0)) { + m = catalog[0].average[n].offset; + /* find next valid entry -- notice lack of error checking... */ + for (j = 0; (j < Nmeas) && (next[m+j] == -2); j++); + catalog[0].average[n].offset = m + j; + } + } + fprintf (stderr, "found %d meas to remove\n", Nmeasfound); + + if (VERBOSE) fprintf (stderr, "fixing missing..."); + /** find missing in time range of image **/ + for (i = 0; (i < Nmiss); i++) { + if (next_miss[i] == -2) continue; + if (catalog[0].missing[i].t < start) continue; + if (catalog[0].missing[i].t > end) continue; + + next_miss[i] = -2; /* we delete this one */ + /* fix the list links: connect the previous valid link to the next valid link */ + for (j = i; (j >= 0) && (next_miss[j] == -2); j--); /* find previous entry to fix link */ + if (j >= 0) { /* if j < 0, there is no previous valid link, ignore this step */ + if (next_miss[j] != i) { + fprintf (stderr, "error? (3) this link seems to have been lost\n"); + fprintf (stderr, "j: %d, next_miss[j]: %d, i: %d\n", j, next_miss[j], i); + exit (1); + } + /* find next valid entry to fix link */ + for (k = i; (k < Nmiss) && (next_miss[k] == -2); k++); + if (k < Nmiss) + next_miss[j] = k; + else + next_miss[j] = -1; /* last link in list gets a -1 */ + } + + /*** fix the corresponding average entry ***/ + n = ave_miss[i]; + if (catalog[0].average[n].Nn == 0) { /* this should never happen */ + fprintf (stderr, "error? we deleted one too many missing?\n"); + exit (1); + } + catalog[0].average[n].Nn --; + /* this was first entry in list */ + if ((catalog[0].average[n].missing == i) && (catalog[0].average[n].Nn > 0)) { + m = catalog[0].average[n].missing; + for (j = 0; (j < Nmiss) && (next_miss[m+j] == -2); j++); + catalog[0].average[n].missing = m + j; + } + } + + /* we should delete the references to blended image and cat stars ?? */ + /* or drop since we are changing this concept ?? */ + + /* fix Average list: delete entries with Nm == 0 */ + for (i = j = 0; (i < Nave) && (j < Nave); i++, j++) { + for (; (j < Nave) && (catalog[0].average[j].Nm == 0); j++); + if ((i != j) && (j < Nave)) { + catalog[0].average[i] = catalog[0].average[j]; + for (k = 0; k < catalog[0].Nsecfilt; k++) { + catalog[0].secfilt[i*Nsecfilt + k] = catalog[0].secfilt[j*Nsecfilt + k]; + } + } + if (j == Nave) i--; + } + Nave = i; + REALLOCATE (catalog[0].average, Average, Nave); + REALLOCATE (catalog[0].secfilt, SecFilt, MAX (1, Nave*Nsecfilt)); + + /* fix order of Measure (memory intensive, but fast) */ + N = 0; + ALLOCATE (tmpmeasure, Measure, Nmeas); + for (i = 0; i < Nave; i++) { + n = catalog[0].average[i].offset; + catalog[0].average[i].offset = N; + for (k = 0; k < catalog[0].average[i].Nm; k++, N++) { + if ((n == -1) || (n == -2)) { + fprintf (stderr, "error: linked list is confused\n"); + exit (1); + } + tmpmeasure[N] = catalog[0].measure[n]; + tmpmeasure[N].averef = i; + n = next[n]; + } + } + Nmeas = N; + free (catalog[0].measure); + catalog[0].measure = tmpmeasure; + REALLOCATE (catalog[0].measure, Measure, Nmeas); + + /* fprintf (stderr, "fixing the mising order...\n"); */ + /* fix order of Missing (memory intensive, but fast) */ + N = 0; + ALLOCATE (tmpmissing, Missing, Nmiss); + for (i = 0; i < Nave; i++) { + if (catalog[0].average[i].Nn > 0) { + n = catalog[0].average[i].missing; + catalog[0].average[i].missing = N; + for (k = 0; k < catalog[0].average[i].Nn; k++, N++) { + if ((n == -1) || (n == -2)) { + fprintf (stderr, "error: linked list is confused\n"); + exit (1); + } + tmpmissing[N] = catalog[0].missing[n]; + n = next_miss[n]; + } + } + } + Nmiss = N; + free (catalog[0].missing); + catalog[0].missing = tmpmissing; + REALLOCATE (catalog[0].missing, Missing, Nmiss); + + fprintf (stderr, "\n"); + + catalog[0].Naverage = Nave; + catalog[0].Nmeasure = Nmeas; + catalog[0].Nmissing = Nmiss; + if (VERBOSE) fprintf (stderr, " ending with Nave, Nmeas, Nmiss: %d %d %d\n", Nave, Nmeas, Nmiss); + + free (next); + free (next_miss); + free (ave_miss); + return; +} + Index: /tags/delstar-1-2/Ohana/src/delstar/src/gcatalog.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/gcatalog.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/gcatalog.c (revision 4683) @@ -0,0 +1,45 @@ +# include "delstar.h" + +int gcatalog (Catalog *catalog) { + + int Nsecfilt; + char mode; + + /* read catalog header */ + mode = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF; + if (!load_catalog (catalog, mode, VERBOSE)) { + fprintf (stderr, "ERROR: failure loading catalog\n"); + exit (1); + } + + /* check Nsecfile value, update if needed */ + Nsecfilt = GetPhotcodeNsecfilt (); + if (catalog[0].Nsecfilt < Nsecfilt) { + + int i, j, Nextra, in, out; + SecFilt *insec, *outsec; + + Nextra = Nsecfilt - catalog[0].Nsecfilt; + insec = catalog[0].secfilt; + ALLOCATE (outsec, SecFilt, catalog[0].Naverage * Nsecfilt); + for (in = out = i = 0; i < catalog[0].Naverage; i++) { + for (j = 0; j < catalog[0].Nsecfilt; j++, in++, out++) { + outsec[out].M = insec[in].M; + outsec[out].Xm = insec[in].Xm; + } + for (j = 0; j < Nextra; j++, out++) { + outsec[out].M = NO_MAG; + outsec[out].Xm = NO_MAG; + } + } + free (catalog[0].secfilt); + catalog[0].secfilt = outsec; + catalog[0].Nsecfilt = Nsecfilt; + } + + if (catalog[0].Nsecfilt > Nsecfilt) { + fprintf (stderr, "ERROR: can't reduce number of secondary filters\n"); + exit (1); + } + return (TRUE); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/gcatstats.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/gcatstats.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/gcatstats.c (revision 4683) @@ -0,0 +1,57 @@ +# include "delstar.h" + +gcatstats (catalog, catstats) +Catalog catalog[]; +CatStats catstats[]; +{ + + int i; + double RaCenter, DecCenter; + double MinRA, MaxRA, MinDEC, MaxDEC; + float *X1, *Y1; + int *N1; + Coords tcoords; + + fits_scan (&catalog[0].header, "RA0", "%lf", 1, &MinRA); + fits_scan (&catalog[0].header, "RA1", "%lf", 1, &MaxRA); + fits_scan (&catalog[0].header, "DEC0", "%lf", 1, &MinDEC); + fits_scan (&catalog[0].header, "DEC1", "%lf", 1, &MaxDEC); + + /* double check on region RA and DEC ranges */ + DecCenter = 0.5*(MinDEC + MaxDEC); + RaCenter = 0.5*(MinRA + MaxRA); + if (MaxDEC > 86.25) { /* we are on the pole */ + DecCenter = 90.0; + RaCenter = 0.0; + } + + catstats[0].RA[0] = MinRA; + catstats[0].RA[1] = MaxRA; + catstats[0].DEC[0] = MinDEC; + catstats[0].DEC[1] = MaxDEC; + + /** allocate local arrays **/ + ALLOCATE (catstats[0].X, float, catalog[0].Naverage); + ALLOCATE (catstats[0].Y, float, catalog[0].Naverage); + ALLOCATE (catstats[0].N, int, catalog[0].Naverage); + + /* project onto rectilinear grid with 1 arcsec pixels, sort by X */ + /* reference for coords is center of field */ + catstats[0].coords.crval1 = RaCenter; + catstats[0].coords.crval2 = DecCenter; + catstats[0].coords.crpix1 = catstats[0].coords.crpix2 = 0.0; + catstats[0].coords.cdelt1 = catstats[0].coords.cdelt2 = 1.0 / 3600.0; + catstats[0].coords.pc1_1 = catstats[0].coords.pc2_2 = 1.0; + catstats[0].coords.pc1_2 = catstats[0].coords.pc2_1 = 0.0; + strcpy (catstats[0].coords.ctype, "RA---TAN"); + + X1 = catstats[0].X; + Y1 = catstats[0].Y; + for (i = 0; i < catalog[0].Naverage; i++, X1++, Y1++) { + fRD_to_XY (X1, Y1, catalog[0].average[i].R, catalog[0].average[i].D, &catstats[0].coords); + catstats[0].N[i] = i; + } + if (catalog[0].Naverage > 1) sort_lists (catstats[0].X, catstats[0].Y, catstats[0].N, catalog[0].Naverage); + +} + Index: /tags/delstar-1-2/Ohana/src/delstar/src/gimages.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/gimages.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/gimages.c (revision 4683) @@ -0,0 +1,87 @@ +# include "delstar.h" + +/* load information about image from image header + * this should be the same as addstar/gstars, but it is not... + */ + +Image *gimages (char *filename) { + + FILE *f; + Header header; + char photcode[64], *c; + double tmp; + Image *image; + int Nc; + + ALLOCATE (image, Image, 1); + /* load header */ + if (!fits_read_header (filename, &header)) { + Shutdown ("ERROR: can't find image file %s", filename); + } + + /* open file */ + f = fopen (filename, "r"); + if (f == NULL) { + Shutdown ("ERROR: can't find data file %s", filename); + } + fseek (f, header.size, SEEK_SET); + + /* add file name to image structure */ + c = strrchr (filename, 0x2f); + if (c == (char *) NULL) { + strcpy (image[0].name, filename); + } else { + strcpy (image[0].name, (c+1)); + } + + /* get astrometry information */ + if (!GetCoords (&image[0].coords, &header)) { + Shutdown ("ERROR: no astrometric solution in header"); + } + while (image[0].coords.crval1 < 0) image[0].coords.crval1 += 360.0; + while (image[0].coords.crval1 > 360.0) image[0].coords.crval1 -= 360.0; + + /* CERROR in data file is in pixels, convert to 20*arcsec */ + image[0].cerror = tmp * 50.0 * image[0].coords.cdelt1 * 3600.0; + while (image[0].coords.crval1 < 0) image[0].coords.crval1 += 360.0; + while (image[0].coords.crval1 > 360.0) image[0].coords.crval1 -= 360.0; + + /* get other header info */ + fits_scan (&header, "NAXIS1", "%hd", 1, &image[0].NX); + fits_scan (&header, "NAXIS2", "%hd", 1, &image[0].NY); + + fits_scan (&header, "PHOTCODE", "%s", 1, photcode); + Nc = GetPhotcodeCodebyName (photcode); + if (!Nc) { + Shutdown ("ERROR: photcode %s not found in photcode table", photcode); + } + image[0].source = Nc; + + tmp = 0; + fits_scan (&header, "FLIMIT", "%lf", 1, &tmp); + image[0].detection_limit = tmp * 10.0; + + tmp = 0; + fits_scan (&header, "FSATUR", "%lf", 1, &tmp); + image[0].saturation_limit = tmp * 10.0; + + if (!fits_scan (&header, "TZERO", "%d", 1, &image[0].tzero)) { + image[0].tzero = parse_time (&header); + } + + tmp = 0; + fits_scan (&header, "TRATE", "%lf", 1, &tmp); + image[0].trate = 10000 * tmp; + + tmp = 0; + fits_scan (&header, "AIRMASS", "%lf", 1, &tmp); + image[0].secz = 1000*tmp; + + /* secz is in units milli-airmass */ + image[0].Mcal = ALPHA*(image[0].secz - 1000); + image[0].Xm = NO_MAG; + + free (header.buffer); + + return (image); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/gregion_image.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/gregion_image.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/gregion_image.c (revision 4683) @@ -0,0 +1,33 @@ +# include "delstar.h" + +/* given image with coords, find regions which overlap image */ +GSCRegion *gregion_image (Image *image, int *Nregions) { + + GSCRegion *regions, impatch; + double x, y, X, Y, ra, dec; + + /* this has trouble for images which overlap 0,360 */ + impatch.RA[0] = 360; + impatch.RA[1] = 0; + impatch.DEC[0] = +90; + impatch.DEC[1] = -90; + + if (!strcmp (&image[0].coords.ctype[4], "-WRP")) { + Shutdown ("method for deleting mosaic data not defined"); + // RegisterMosaic (MOSAIC); + } + + for (x = 0; x <= 1; x++) { + for (y = 0; y <= 1; y++) { + X = image[0].NX*(1.1*x - 0.05); + Y = image[0].NY*(1.1*y - 0.05); + XY_to_RD (&ra, &dec, X, Y, &image[0].coords); + impatch.RA[0] = MIN (ra, impatch.RA[0]); + impatch.RA[1] = MAX (ra, impatch.RA[1]); + impatch.DEC[0] = MIN (dec, impatch.DEC[0]); + impatch.DEC[1] = MAX (dec, impatch.DEC[1]); + } + } + regions = gregion_patch (&impatch, Nregions); + return (regions); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/gregion_patch.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/gregion_patch.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/gregion_patch.c (revision 4683) @@ -0,0 +1,43 @@ +# include "delstar.h" + +GSCRegion *gregion_patch (GSCRegion *patch, int *nregions) { + + int i, Ns, Ne, Ngsc, Nregions; + GSCRegion *gsc, *regions; + + gsc = LoadRegions (&Ngsc); + + /* find lower and upper boundaries */ + Ns = FindRegionDecBandStart (gsc, Ngsc, patch[0].DEC[0]); + Ne = FindRegionDecBandStop (gsc, Ngsc, patch[0].DEC[1]); + if (Ns < 0) Shutdown ("invalid patch lower bound %f", patch[0].DEC[0]); + if (Ne < 0) Shutdown ("invalid patch upper bound %f", patch[0].DEC[1]); + + Nregions = 0; + ALLOCATE (regions, GSCRegion, Ne - Ns); + + /* find all regions in range within ra range */ + for (i = Ns; i < Ne; i++) { + if (gsc[i].RA[0] >= patch[0].RA[1]) continue; + if (gsc[i].RA[1] < patch[0].RA[0]) continue; + if (gsc[i].DEC[0] >= patch[0].DEC[1]) continue; + if (gsc[i].DEC[1] < patch[0].DEC[0]) continue; + /* add to the list */ + regions[Nregions] = gsc[i]; + Nregions ++; + } + + if (VERBOSE) fprintf (stderr, "Nregions in patch: %d, Ne: %d, Ns: %d\n", Nregions, Ne, Ns); + if (Nregions > (Ne - Ns)) { + fprintf (stderr, "ERROR: too many regions found\n"); + exit (1); + } + + *nregions = Nregions; + return (regions); +} + +/* given a ra,dec bounded region, find all region files which overlap it + - split region into two if overlapping 0,360 boundary +*/ + Index: /tags/delstar-1-2/Ohana/src/delstar/src/gregions.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/gregions.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/gregions.c (revision 4683) @@ -0,0 +1,53 @@ +# include "delstar.h" + +GSCRegion *gregions (Image *image, int *Nregions) { + + GSCRegion *region; + FILE *f; + double x, y; + double dr, dd, dec, ra; + int i, j, done, nregion, NREGION; + + f = fopen (GSCFILE, "r"); + if (f == NULL) { + fprintf (stderr, "ERROR: can't find GSC region file %s\n", GSCFILE); + exit (1); + } + + /* find regions at image corners */ + NREGION = 10; + ALLOCATE (region, GSCRegion, NREGION); + nregion = 0; + + /* look for new regions on grid across image */ + for (x = 0; x <= 1; x+=0.25) { + for (y = 0; y <= 1; y+=0.25) { + XY_to_RD (&ra, &dec, image[0].NX*(1.1*x - 0.05), image[0].NY*(1.1*y - 0.05), &image[0].coords); + aregion (®ion[nregion], f, ra, dec); + done = FALSE; + for (j = 0; (j < nregion) && !done; j++) { + if (!strcmp (region[nregion].filename, region[j].filename)) { + nregion --; + done = TRUE; + } + } + nregion ++; + if (nregion == NREGION) { + NREGION += 10; + REALLOCATE (region, GSCRegion, NREGION); + } + } + } + + if (VERBOSE) { + fprintf (stderr, "found %d region files:\n", nregion); + for (i = 0; i < nregion; i++) { + fprintf (stderr, " %d %s\n", i, region[i].filename); + } + } + *Nregions = nregion; + + fclose (f); + return (region); + +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/gtimes.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/gtimes.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/gtimes.c (revision 4683) @@ -0,0 +1,68 @@ +# include "delstar.h" + +Image *gtimes (int *NIMAGE) { + + int i, npimage, NPIMAGE, ntimage, Ntimage, NTIMAGE, Nimage; + time_t tstart, tend; + double range; + Image *pimage, *timage; + Header header; + FILE *f; + + /* globals: START, END, PHOTCODE */ + + /* get image db file pointer */ + f = GetDB (&dbstate); + if (dbstate == LCK_EMPTY) { + Shutdown ("database is empty"); + } + fseek (f, 0, SEEK_SET); + + /* read header */ + if (!fits_fread_header (f, &header)) { + Shutdown ("can't read image catalog %s", ImageCat); + } + fseek (f, header.size, SEEK_SET); + + Nimage = 0; + fits_scan (&header, "NIMAGES", "%d", 1, &Nimage); + + NTIMAGE = 100; + ALLOCATE (timage, Image, NTIMAGE); + + npimage = 0; + NPIMAGE = 20; + ALLOCATE (pimage, Image, NPIMAGE); + + /* add match on photcode */ + for (Ntimage = 0; Ntimage < Nimage; Ntimage += ntimage) { + ntimage = Fread (timage, sizeof(Image), NTIMAGE, f, "image"); + for (i = 0; i < ntimage; i++) { + tstart = timage[i].tzero; + tend = timage[i].tzero + 1e-4*timage[i].trate*timage[i].NY; + if ((tstart <= END) && (tend >= START)) { + pimage[npimage] = timage[i]; + npimage ++; + if (npimage == NPIMAGE) { + NPIMAGE += 20; + REALLOCATE (pimage, Image, NPIMAGE); + } + } + } + } + free (timage); + + { + double Ra, Dec; + if (VERBOSE) fprintf (stderr, "found the following images:\n"); + for (i = 0; VERBOSE && (i < npimage); i++) { + XY_to_RD (&Ra, &Dec, 0.5*pimage[i].NX, 0.5*pimage[i].NY, &pimage[i].coords); + fprintf (stderr, "%3d %s %8.4f %8.4f %d %d\n", i, pimage[i].name, Ra, Dec, + pimage[i].tzero, pimage[i].nstar); + } + } + + *NIMAGE = npimage; + return (pimage); + +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/image-db.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/image-db.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/image-db.c (revision 4683) @@ -0,0 +1,119 @@ +# include "delstar.h" + +# define LOCK LCK_XCLD + +/* file-pointer to Image.db table. if image db does not yet exist */ +static FILE *f = (FILE *) NULL; +static int dbstate = LCK_UNLOCK; + +static Header header; +static Image *image = NULL; +static int Nimage; + +FILE *GetDB (int *state) { + *state = dbstate; + return (f); +} + +Image *GetImages (int *nimage) { + *nimage = Nimage; + return (image); +} + +int SetImages (Image *new, int Nnew) { + Nimage = Nnew; + image = new; + fits_modify (&header, "NIMAGES", "%d", 1, Nimage); + return (TRUE); +} + +/* clean up open / locked ImageCat before shutting down */ +int Shutdown (char *format, ...) { + va_list argp; + char *formatplus; + + ALLOCATE (formatplus, char, strlen(format)); + strcpy (formatplus, format); + strcat (formatplus, "\n"); + + va_start (argp, format); + vfprintf (stderr, formatplus, argp); + va_end (argp); + + SetProtect (TRUE); + fclearlockfile (ImageCat, f, LOCK, &dbstate); + fprintf (stderr, "ERROR: addstar halted\n"); + exit (1); +} + +void lock_image_db () { + /* lock the image catalog */ + check_permissions (ImageCat); + f = fsetlockfile (ImageCat, 3600.0, LOCK, &dbstate); + if (f == NULL) { + Shutdown ("ERROR: can't lock image catalog"); + } + if (dbstate == LCK_EMPTY) { + Shutdown ("ERROR: database is empty"); + } + fseek (f, 0, SEEK_SET); +} + +int load_image_db () { + + int nimage; + + /* read header */ + fseek (f, 0, SEEK_SET); + if (!fits_fread_header (f, &header)) { + Shutdown ("can't read image catalog %s", ImageCat); + } + fseek (f, header.size, SEEK_SET); + + Nimage = 0; + fits_scan (&header, "NIMAGES", "%d", 1, &Nimage); + ALLOCATE (image, Image, Nimage); + + nimage = Fread (image, sizeof(Image), Nimage, f, "image"); + if (nimage != Nimage) { + Shutdown ("error reading images from db file"); + } + + return (TRUE); +} + +int save_image_db () { + + int status; + char line[1024]; + mode_t mode; + + /* use cp to make backup copy */ + sprintf (line, "cp %s %s~", ImageCat, ImageCat); + status = system (line); + if (status) { + Shutdown ("ERROR: unable to create %s~, exiting", ImageCat); + } + mode = S_IRUSR | S_IWUSR | S_IRGRP | S_IWGRP | S_IROTH | S_IWOTH; + sprintf (line, "%s~", ImageCat); + chmod (line, mode); + + /* update header */ + fits_modify (&header, "NIMAGES", "%d", 1, Nimage); + + /* write out data to already opened file */ + fseek (f, 0, SEEK_SET); + Fwrite (header.buffer, 1, header.size, f, "char"); + Fwrite (image, sizeof(Image), Nimage, f, "image"); + + return (TRUE); +} + +void unlock_image_db () { + + mode_t mode; + + fclearlockfile (ImageCat, f, LOCK, &dbstate); + mode = S_IRUSR | S_IWUSR | S_IRGRP | S_IWGRP | S_IROTH | S_IWOTH; + chmod (ImageCat, mode); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/match_images.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/match_images.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/match_images.c (revision 4683) @@ -0,0 +1,52 @@ +# include "delstar.h" + +void match_images (Catalog *catalog, Image *image, int Nimage) { + + int j, k, found; + unsigned int *start, *stop; + + /* this must be allocated so future free will not fail */ + ALLOCATE (catalog[0].image, int, MAX (catalog[0].Nmeasure, 1)); + if (catalog[0].Naverage == 0) { + if (VERBOSE) fprintf (stderr, "no stars in catalog, skipping\n"); + return (FALSE); + } + + ALLOCATE (start, unsigned int, Nimage); + ALLOCATE (stop, unsigned int, Nimage); + for (j = 0; j < Nimage; j++) { + start[j] = image[j].tzero - MAX(0.05*image[j].trate*image[j].NY, 1); + stop[j] = image[j].tzero + MAX(1.05*image[j].trate*image[j].NY, 1); + } + + for (j = 0; j < catalog[0].Nmeasure; j++) { + found = FALSE; + if (catalog[0].measure[j].t == 0) { + catalog[0].image[j] = -1; + found = TRUE; + } + for (k = 0; (k < Nimage) && !found; k++) { + if ((catalog[0].measure[j].t >= start[k]) && + (catalog[0].measure[j].t <= stop[k])) { + catalog[0].image[j] = k; + found = TRUE; + } + } + if (!found) { + catalog[0].image[j] = -1; + /* fprintf (stderr, "missing: %d %d\n", catalog[0].image[j], catalog[0].measure[j].t); */ + } + } + free (start); + free (stop); +} + + /* this routine uses the time of each measurement to match the +measurement with an image. Since the measurement is only store to 1 +sec accuracy, which corresponds to roughly 30 rows at nominal speed, +we can't tell exactly which image the star come from. However, this +doesn't matter, and in fact this helps a bit: a measurement from the +top of one image is the same as from the bottom of the next. +Therefore, we intentionally blur the edges of the images by 5%, which +will help to tie together neighboring images... */ + Index: /tags/delstar-1-2/Ohana/src/delstar/src/missed.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/missed.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/missed.c (revision 4683) @@ -0,0 +1,24 @@ +delete_missed () { + catalog.filename = argv[1]; + switch (lock_catalog (&catalog, LCK_XCLD)) { + case 0: + fprintf (stderr, "ERROR: can't lock file\n"); + exit (1); + case 1: + gcatalog (&catalog); + break; + case 2: + fprintf (stderr, "no data in %s, skipping\n", catalog.filename); + catalog.Naverage = 0; + catalog.Nmeasure = 0; + break; + default: + fprintf (stderr, "weird lock_catalog exit state\n"); + exit (1); + } + delete_missed (&catalog); + wcatalog (&catalog); + fclearlockfile (ImageCat, f, LCK_HARD, &dbstate); + fprintf (stderr, "SUCCESS\n"); + exit (0); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/orphans.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/orphans.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/orphans.c (revision 4683) @@ -0,0 +1,14 @@ + +delete_orphans () { + ALLOCATE (region, GSCRegion, 1); + catalog.filename = argv[1]; + gcatalog (&catalog); + gcatstats (&catalog, &catstats); + image = find_images_region (&catstats, &Nimage); + match_images (&catalog, image, Nimage); + delete_orphans (&catalog); + wcatalog (&catalog); + fclearlockfile (ImageCat, f, LCK_HARD, &dbstate); + fprintf (stderr, "SUCCESS\n"); + exit (0); +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/parse_time.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/parse_time.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/parse_time.c (revision 4683) @@ -0,0 +1,112 @@ +# include "delstar.h" + +int parse_time (Header *header) { + + double jd; + int Ny, Nf, mode; + int Nsec, hour, min, sec, year, month, day; + char *py, *pm, *pd, *c; + char line[256]; + + /* we want to find JD or MJD to get Nsec (seconds since 01/01/1970) */ + + /* try JD first */ + if (strcasecmp (JDKeyword, "NONE")) { + uppercase (JDKeyword); + fits_scan (header, JDKeyword, "%lf", 1, &jd); + Nsec = (jd - 2440587.5)*86400; + return (Nsec); + } + + /* try MJD next */ + if (strcasecmp (MJDKeyword, "NONE")) { + uppercase (MJDKeyword); + fits_scan (header, MJDKeyword, "%lf", 1, &jd); + Nsec = (jd - 40587.0)*86400; + return (Nsec); + } + + /* get UT and DATE */ + uppercase (UTKeyword); + fits_scan (header, UTKeyword, "%s", 1, line); + /* remove ':' characters */ + for (c = strchr (line, 0x3a); c != (char *) NULL; c = strchr (line, 0x3a)) { *c = ' '; } + sscanf (line, "%d %d %d", &hour, &min, &sec); + + /* parse mode line */ + uppercase (DateMode); + for (Ny = 0, c = strchr (DateMode, 'Y'); c != (char ) NULL; c = strchr (c + 1, 'Y'), Ny++); + if ((Ny != 2) && (Ny != 4)) { + Shutdown ("error in DATE-MODE format: %s", DateMode); + } + py = strchr (DateMode, 'Y'); + pm = strchr (DateMode, 'M'); + pd = strchr (DateMode, 'D'); + if ((py == (char *) NULL) || (pm == (char *) NULL) || (pd == (char *) NULL)) { + Shutdown ("error in DATE-MODE format: %s", DateMode); + } + if ((py > pm) && (py < pd)) { + Shutdown ("error in DATE-MODE format: %s", DateMode); + } + if ((py > pd) && (py < pm)) { + Shutdown ("error in DATE-MODE format: %s", DateMode); + } + mode = 0; + if ((py < pm) && (pm < pd)) { mode = 1; } /* yyyy-mm-dd */ + if ((py < pm) && (pm > pd)) { mode = 2; } /* yyyy-dd-mm */ + if ((py > pm) && (pm < pd)) { mode = 3; } /* mm-dd-yyyy */ + if ((py > pm) && (pm > pd)) { mode = 4; } /* dd-mm-yyyy */ + if (!mode) { + Shutdown ("error in DATE-MODE format: %s", DateMode); + } + + /* parse date entry */ + uppercase (DateKeyword); + fits_scan (header, DateKeyword, "%s", 1, line); + /* remove possible separators: ':', '/' '.', '-' */ + for (c = strchr (line, 0x3a); c != (char *) NULL; c = strchr (line, 0x3a)) { *c = ' '; } + for (c = strchr (line, 0x2f); c != (char *) NULL; c = strchr (line, 0x2f)) { *c = ' '; } + for (c = strchr (line, 0x2e); c != (char *) NULL; c = strchr (line, 0x2e)) { *c = ' '; } + for (c = strchr (line, 0x2d); c != (char *) NULL; c = strchr (line, 0x2d)) { *c = ' '; } + + switch (mode) { + case 1: + Nf = sscanf (line, "%d %d %d", &year, &month, &day); + break; + case 2: + Nf = sscanf (line, "%d %d %d", &year, &day, &month); + break; + case 3: + Nf = sscanf (line, "%d %d %d", &month, &day, &year); + break; + case 4: + Nf = sscanf (line, "%d %d %d", &day, &month, &year); + break; + } + if (Nf != 3) { + Shutdown ("error in date entry (%s) or DATE-MODE format (%s)", line, DateMode); + } + + if (year > 1000) { + if (Ny == 2) { + fprintf (stderr, "warning: mode line claims 2 digit year, but 4 digit year found\n"); + } + } else { + if (Ny == 4) { + fprintf (stderr, "warning: mode line claims 4 digit year, but 2 digit year found\n"); + } + if (year < 50) year += 100; + year += 1900; + } + + /* convert yy.mm.dd hh.mm.ss to Nsec since 1970 (jd = 2440587.5) */ + /* note that in this section, tm_mon has range 1-12, unlike for gmtime () */ + jd = day - 32075 + (int)(1461*(year + 4800 + (int)(((month)-14)/12))/4) + + (int)(367*((month) - 2 - (int)(((month) - 14)/12)*12)/12) + - (int)(3*(int)((year + 4900 + (int)(((month) - 14)/12))/100)/4) - 0.5; + /* jd is the julian day of the whole day only not the time */ + Nsec = (jd - 2440587.5)*86400 + 3600.0*hour + min*60.0 + sec; + + return (Nsec); + +} Index: /tags/delstar-1-2/Ohana/src/delstar/src/sort_lists.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/sort_lists.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/sort_lists.c (revision 4683) @@ -0,0 +1,80 @@ +# include "delstar.h" + +void sort_lists (float *X, float *Y, int *S, int N) { + + int l,j,ir,i; + double tX, tY, tS; + + l = N >> 1; + ir = N - 1; + for (;;) { + if (l > 0) { + l--; + tX = X[l]; + tY = Y[l]; + tS = S[l]; + } else { + tX = X[ir]; + X[ir] = X[0]; + tY = Y[ir]; + Y[ir] = Y[0]; + tS = S[ir]; + S[ir] = S[0]; + if (--ir == 0) { + X[0] = tX; + Y[0] = tY; + S[0] = tS; + return; + } + } + i = l; + j = (l << 1) + 1; + while (j <= ir) { + if (j < ir && X[j] < X[j+1]) j++; + if (tX < X[j]) { + X[i] = X[j]; + Y[i] = Y[j]; + S[i] = S[j]; + j += (i=j) + 1; + } + else j = ir + 1; + } + X[i] = tX; + Y[i] = tY; + S[i] = tS; + } +} + +void sort_regions (GSCRegion *region, int N) { + + int l,j,ir,i; + GSCRegion tmp; + + l = N >> 1; + ir = N - 1; + for (;;) { + if (l > 0) { + l--; + tmp = region[l]; + } else { + tmp = region[ir]; + region[ir] = region[0]; + if (--ir == 0) { + region[0] = tmp; + return; + } + } + i = l; + j = (l << 1) + 1; + while (j <= ir) { + if (j < ir && region[j].DEC[0] < region[j+1].DEC[0]) j++; + if (tmp.DEC[0] < region[j].DEC[0]) { + region[i] = region[j]; + j += (i=j) + 1; + } + else j = ir + 1; + } + region[i] = tmp; + } +} + Index: /tags/delstar-1-2/Ohana/src/delstar/src/wcatalog.c =================================================================== --- /tags/delstar-1-2/Ohana/src/delstar/src/wcatalog.c (revision 4683) +++ /tags/delstar-1-2/Ohana/src/delstar/src/wcatalog.c (revision 4683) @@ -0,0 +1,37 @@ +# include "delstar.h" + +int wcatalog (Catalog *catalog) { + + /* we'll leave these for now, but they are LONEOS specific */ + fits_modify (&catalog[0].header, "MARKSTAR", "%t", 1, FALSE); + fits_modify (&catalog[0].header, "ADDUSNO", "%t", 1, FALSE); + fits_modify (&catalog[0].header, "MARKROCK", "%t", 1, FALSE); + + if (!save_catalog (catalog, VERBOSE)) { + fprintf (stderr, "ERROR: failure saving catalog\n"); + exit (1); + } + + /* free, initialize data structures */ + + if (catalog[0].Naverage) { + free (catalog[0].average); + catalog[0].Naverage = 0; + } + if (catalog[0].Nmeasure) { + free (catalog[0].measure); + catalog[0].Nmeasure = 0; + } + if (catalog[0].Nmissing) { + free (catalog[0].missing); + catalog[0].Nmissing = 0; + } + if (catalog[0].Nsecfilt) { + free (catalog[0].secfilt); + catalog[0].Nsecfilt = 0; + } + fits_free_header (&catalog[0].header); + + return (TRUE); +} +