Index: trunk/Ohana/configure.tcsh =================================================================== --- trunk/Ohana/configure.tcsh (revision 31027) +++ trunk/Ohana/configure.tcsh (revision 31160) @@ -321,6 +321,6 @@ set sysincpath = "$sysincpath /sw/include" set dlltype = dylib - set CFLAGS = "$CFLAGS -D_DARWIN_C_SOURCE" - set defines = "-D_DARWIN_C_SOURCE" + # set CFLAGS = "$CFLAGS -D_DARWIN_C_SOURCE" + # set defines = "-D_DARWIN_C_SOURCE" breaksw; case HP-UX: Index: trunk/Ohana/src/addstar/include/addstar.h =================================================================== --- trunk/Ohana/src/addstar/include/addstar.h (revision 31027) +++ trunk/Ohana/src/addstar/include/addstar.h (revision 31160) @@ -263,4 +263,6 @@ Stars *Convert_PS1_V2 PROTO((FTable *table, unsigned int *nstars)); Stars *Convert_PS1_V3 PROTO((FTable *table, unsigned int *nstars)); +Stars *Convert_PS1_SV1 PROTO((FTable *table, unsigned int *nstars)); +Stars *Convert_PS1_SV1_Alt PROTO((FTable *table, unsigned int *nstars)); int InitStar PROTO((Stars *star)); Index: trunk/Ohana/src/addstar/src/MatchHeaders.c =================================================================== --- trunk/Ohana/src/addstar/src/MatchHeaders.c (revision 31027) +++ trunk/Ohana/src/addstar/src/MatchHeaders.c (revision 31160) @@ -58,4 +58,5 @@ if (!strcmp (exttype, "PS1_V2")) goto keep; if (!strcmp (exttype, "PS1_V3")) goto keep; + if (!strcmp (exttype, "PS1_SV1")) goto keep; continue; Index: trunk/Ohana/src/addstar/src/ReadStarsFITS.c =================================================================== --- trunk/Ohana/src/addstar/src/ReadStarsFITS.c (revision 31027) +++ trunk/Ohana/src/addstar/src/ReadStarsFITS.c (revision 31160) @@ -54,4 +54,7 @@ stars = Convert_PS1_V3 (&table, &Nstars); } + if (!strcmp (type, "PS1_SV1")) { + stars = Convert_PS1_SV1 (&table, &Nstars); + } if (stars == NULL) { fprintf (stderr, "invalid table type %s\n", type); @@ -501,2 +504,146 @@ return (stars); } + +Stars *Convert_PS1_SV1 (FTable *table, unsigned int *nstars) { + + off_t Nstars; + unsigned int i; + double ZeroPt; + Stars *stars; + CMF_PS1_SV1 *ps1data; + + if (table[0].header[0].Naxis[0] == 196) { + stars = Convert_PS1_SV1_Alt (table, nstars); + return (stars); + } + + ps1data = gfits_table_get_CMF_PS1_SV1 (table, &Nstars, NULL); + if (!ps1data) { + fprintf (stderr, "skipping inconsistent entry\n"); + return (NULL); + } + ZeroPt = GetZeroPoint(); + + ALLOCATE (stars, Stars, Nstars); + for (i = 0; i < Nstars; i++) { + InitStar (&stars[i]); + stars[i].measure.Xccd = ps1data[i].X; + stars[i].measure.Yccd = ps1data[i].Y; + stars[i].measure.dXccd = ToShortPixels(ps1data[i].dX); + stars[i].measure.dYccd = ToShortPixels(ps1data[i].dY); + + stars[i].measure.posangle = ToShortDegrees(ps1data[i].posangle); + stars[i].measure.pltscale = ps1data[i].pltscale; + + if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) { + stars[i].measure.M = NAN; + } else { + stars[i].measure.M = ps1data[i].M + ZeroPt; + } + stars[i].measure.dM = ps1data[i].dM; + stars[i].measure.dMcal = ps1data[i].dMcal; + stars[i].measure.Map = ps1data[i].Map + ZeroPt; + + stars[i].measure.Sky = ps1data[i].sky; + stars[i].measure.dSky = ps1data[i].dSky; + + stars[i].measure.psfChisq = ps1data[i].psfChisq; + stars[i].measure.psfQual = ps1data[i].psfQual; + stars[i].measure.psfNdof = ps1data[i].psfNdof; + stars[i].measure.psfNpix = ps1data[i].psfNpix; + stars[i].measure.crNsigma = ps1data[i].crNsigma; + stars[i].measure.extNsigma = ps1data[i].extNsigma; + + stars[i].measure.FWx = ToShortPixels(ps1data[i].fx); + stars[i].measure.FWy = ToShortPixels(ps1data[i].fy); + stars[i].measure.theta = ToShortDegrees(ps1data[i].df); + + stars[i].measure.Mxx = ToShortPixels(ps1data[i].Mxx); + stars[i].measure.Mxy = ToShortPixels(ps1data[i].Mxy); + stars[i].measure.Myy = ToShortPixels(ps1data[i].Myy); + + stars[i].measure.photFlags = ps1data[i].flags; + + // this is may optionally be replaced by the internal sequence (see FilterStars.c) + stars[i].measure.detID = ps1data[i].detID; + + // the Average fields and the following Measure fields are set in FilterStars after + // the image metadata is in hand: dR, dD, Mcal, dt, airmass, az, t, imageID, extID, + // averef is set in find_matches, dbFlags is zero on ingest. + + // the following fields are currently not being set anywhere: t_msec + } + *nstars = Nstars; + return (stars); +} + +Stars *Convert_PS1_SV1_Alt (FTable *table, unsigned int *nstars) { + + off_t Nstars; + unsigned int i; + double ZeroPt; + Stars *stars; + CMF_PS1_SV1 *ps1data; + + // some test output files were produced called CMF_PS1_SV1 but with mismatch byte boundaries + + ps1data = gfits_table_get_CMF_PS1_SV1_Alt (table, &Nstars, NULL); + if (!ps1data) { + fprintf (stderr, "skipping inconsistent entry\n"); + return (NULL); + } + ZeroPt = GetZeroPoint(); + + ALLOCATE (stars, Stars, Nstars); + for (i = 0; i < Nstars; i++) { + InitStar (&stars[i]); + stars[i].measure.Xccd = ps1data[i].X; + stars[i].measure.Yccd = ps1data[i].Y; + stars[i].measure.dXccd = ToShortPixels(ps1data[i].dX); + stars[i].measure.dYccd = ToShortPixels(ps1data[i].dY); + + stars[i].measure.posangle = ToShortDegrees(ps1data[i].posangle); + stars[i].measure.pltscale = ps1data[i].pltscale; + + if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) { + stars[i].measure.M = NAN; + } else { + stars[i].measure.M = ps1data[i].M + ZeroPt; + } + stars[i].measure.dM = ps1data[i].dM; + stars[i].measure.dMcal = ps1data[i].dMcal; + stars[i].measure.Map = ps1data[i].Map + ZeroPt; + + stars[i].measure.Sky = ps1data[i].sky; + stars[i].measure.dSky = ps1data[i].dSky; + + stars[i].measure.psfChisq = ps1data[i].psfChisq; + stars[i].measure.psfQual = ps1data[i].psfQual; + stars[i].measure.psfNdof = ps1data[i].psfNdof; + stars[i].measure.psfNpix = ps1data[i].psfNpix; + stars[i].measure.crNsigma = ps1data[i].crNsigma; + stars[i].measure.extNsigma = ps1data[i].extNsigma; + + stars[i].measure.FWx = ToShortPixels(ps1data[i].fx); + stars[i].measure.FWy = ToShortPixels(ps1data[i].fy); + stars[i].measure.theta = ToShortDegrees(ps1data[i].df); + + stars[i].measure.Mxx = ToShortPixels(ps1data[i].Mxx); + stars[i].measure.Mxy = ToShortPixels(ps1data[i].Mxy); + stars[i].measure.Myy = ToShortPixels(ps1data[i].Myy); + + stars[i].measure.photFlags = ps1data[i].flags; + + // this is may optionally be replaced by the internal sequence (see FilterStars.c) + stars[i].measure.detID = ps1data[i].detID; + + // the Average fields and the following Measure fields are set in FilterStars after + // the image metadata is in hand: dR, dD, Mcal, dt, airmass, az, t, imageID, extID, + // averef is set in find_matches, dbFlags is zero on ingest. + + // the following fields are currently not being set anywhere: t_msec + } + *nstars = Nstars; + return (stars); +} + Index: trunk/Ohana/src/addstar/src/resort_catalog.c =================================================================== --- trunk/Ohana/src/addstar/src/resort_catalog.c (revision 31027) +++ trunk/Ohana/src/addstar/src/resort_catalog.c (revision 31160) @@ -92,6 +92,9 @@ averageSeq[i] = measure[i].averef; - myAssert(average[averageSeq[i]].objID == measure[measureSeq[i]].objID, "object / detection mismatch"); - myAssert(average[averageSeq[i]].catID == measure[measureSeq[i]].catID, "object / detection mismatch"); + if (catalog[0].catformat >= DVO_FORMAT_PS1_V1) { + // earlier formats did not carry the objID or catID, so they are not available (we could assign on load, but we don't) + myAssert(average[averageSeq[i]].objID == measure[measureSeq[i]].objID, "object / detection mismatch"); + myAssert(average[averageSeq[i]].catID == measure[measureSeq[i]].catID, "object / detection mismatch"); + } } Index: trunk/Ohana/src/dvomerge/src/dvoverify.c =================================================================== --- trunk/Ohana/src/dvomerge/src/dvoverify.c (revision 31027) +++ trunk/Ohana/src/dvomerge/src/dvoverify.c (revision 31160) @@ -13,6 +13,5 @@ # define DEBUG 0 -int VERBOSE; -int CHECKSORTED; +int VERBOSE = FALSE; int NNotSorted = 0; @@ -28,14 +27,22 @@ // Catalog catalog; - VERBOSE = FALSE; + int CHECKSORTED; + if ((N = get_argument (argc, argv, "-v"))) { VERBOSE = TRUE; remove_argument (N, &argc, argv); } + if ((N = get_argument (argc, argv, "-verbose"))) { + VERBOSE = TRUE; + remove_argument (N, &argc, argv); + } if ((N = get_argument (argc, argv, "-s"))) { CHECKSORTED = TRUE; remove_argument (N, &argc, argv); } - + if ((N = get_argument (argc, argv, "-sorted"))) { + CHECKSORTED = TRUE; + remove_argument (N, &argc, argv); + } // restrict to a portion of the sky @@ -57,6 +64,8 @@ if (argc != 2) { - fprintf (stderr, "USAGE: dvoverify (catdir) [-region Rmin Rmax Dmin Dmax] [-v -s]\n\n -v = verbose \n -s = checks if sorted, return error if not\n"); - fprintf (stderr, " catdir : database of interest\n"); + fprintf (stderr, "USAGE: dvoverify (catdir) [-region Rmin Rmax Dmin Dmax] [-v] [-s]\n\n"); + fprintf (stderr, " -v : VERBOSE\n"); + fprintf (stderr, " -s : checks if sorted, return error if not\n"); + fprintf (stderr, " (catdir) : database of interest\n"); exit (2); } @@ -96,5 +105,4 @@ for (i = 0; i < inlist[0].Nregions; i++) { if (!inlist[0].regions[i][0].table) continue; - if ((NNotSorted > 0) && CHECKSORTED) continue; if (i % 1000 == 0) fprintf (stderr, "."); @@ -116,4 +124,10 @@ if (!CheckCatalogIndexes(catdir, inlist[0].filename[i], inlist[0].regions[i])){ Nbad ++; + } + + // exit immediately if any file are unsorted and we require sorted tables + if (CHECKSORTED && NNotSorted) { + fprintf (stderr, "ERROR: files are not sorted\n"); + exit (1); } } @@ -123,10 +137,9 @@ exit (1); } - if ((NNotSorted > 0) && CHECKSORTED) { - fprintf (stderr, "ERROR: files are not sorted\n"); - exit (1); - } fprintf (stderr, "SUCCESS: no files are bad\n"); + if (NNotSorted) { + fprintf (stderr, "NOTE: %d files are not sorted\n", NNotSorted); + } exit (0); } Index: trunk/Ohana/src/kapa2/src/PSObjects.c =================================================================== --- trunk/Ohana/src/kapa2/src/PSObjects.c (revision 31027) +++ trunk/Ohana/src/kapa2/src/PSObjects.c (revision 31160) @@ -270,5 +270,5 @@ float *x, *y, *z; double mxi, mxj, myi, myj, bxi, bxj, byi, byj, bx, by; - double sx, sy, d, sx1, sy1, sx2, sy2; + double sx, sy, d, sx1, sy1, sx2, sy2, X0, Y0, X1, Y1; mxi = graph[0].axis[0].dfx / (object[0].x1 - object[0].x0); @@ -432,4 +432,9 @@ } if (object[0].ptype == 100) { /* connect a pair of points */ + X0 = graph[0].axis[0].fx; + X1 = graph[0].axis[0].fx + graph[0].axis[0].dfx; + Y0 = graph[0].axis[1].fy; + Y1 = graph[0].axis[1].fy + graph[0].axis[1].dfy; + for (i = 0; i + 1 < object[0].Npts; i+=2) { if (!(finite(x[i]) && finite(y[i]))) continue; @@ -438,5 +443,5 @@ sx2 = x[i+1]*mxi + y[i+1]*mxj + bx; sy2 = x[i+1]*myi + y[i+1]*myj + by; - DrawLine (sx1, sy1, sx2, sy2); + ClipLinePS (sx1, sy1, sx2, sy2, X0, Y0, X1, Y1, f); } } @@ -585,4 +590,9 @@ } if (object[0].ptype == 100) { + X0 = graph[0].axis[0].fx; + X1 = graph[0].axis[0].fx + graph[0].axis[0].dfx; + Y0 = graph[0].axis[1].fy; + Y1 = graph[0].axis[1].fy + graph[0].axis[1].dfy; + for (i = 0; i + 1 < object[0].Npts; i+=2) { if (!(finite(x[i]) && finite(y[i]))) continue; @@ -591,5 +601,5 @@ sx2 = x[i+1]*mxi + y[i+1]*mxj + bx; sy2 = x[i+1]*myi + y[i+1]*myj + by; - DrawLine (sx1, sy1, sx2, sy2); + ClipLinePS (sx1, sy1, sx2, sy2, X0, Y0, X1, Y1, f); } } Index: trunk/Ohana/src/kapa2/src/SetImageData.c =================================================================== --- trunk/Ohana/src/kapa2/src/SetImageData.c (revision 31027) +++ trunk/Ohana/src/kapa2/src/SetImageData.c (revision 31160) @@ -128,5 +128,5 @@ Picture_to_Image (&Xmax, &Ymax, image[0].picture.dx, image[0].picture.dy, &image[0].picture); - KiiSendMessage (sock, "%g %g %g %g", Xmin, Xmax, Ymin, Ymax); + KiiSendMessage (sock, "%g %g %g %g %d %d", Xmin, Xmax, Ymin, Ymax, image[0].picture.dx, image[0].picture.dy); return (TRUE); Index: trunk/Ohana/src/kapa2/src/bDrawObjects.c =================================================================== --- trunk/Ohana/src/kapa2/src/bDrawObjects.c (revision 31027) +++ trunk/Ohana/src/kapa2/src/bDrawObjects.c (revision 31160) @@ -415,4 +415,10 @@ } if (object[0].ptype == 100) { /* connect a pair of points */ + + double X0 = graph[0].axis[0].fx; + double X1 = graph[0].axis[0].fx + graph[0].axis[0].dfx; + double Y0 = graph[0].axis[1].fy; + double Y1 = graph[0].axis[1].fy + graph[0].axis[1].dfy; + for (i = 0; i + 1 < object[0].Npts; i+=2) { if (!(finite(x[i]) && finite(y[i]))) continue; @@ -421,5 +427,5 @@ sx2 = x[i+1]*mxi + y[i+1]*mxj + bx; sy2 = x[i+1]*myi + y[i+1]*myj + by; - DrawLine (buffer, sx1, sy1, sx2, sy2); + bDrawClipLine (buffer, sx1, sy1, sx2, sy2, X0, Y0, X1, Y1); } } @@ -561,4 +567,10 @@ } if (object[0].ptype == 100) { + + double X0 = graph[0].axis[0].fx; + double X1 = graph[0].axis[0].fx + graph[0].axis[0].dfx; + double Y0 = graph[0].axis[1].fy; + double Y1 = graph[0].axis[1].fy + graph[0].axis[1].dfy; + for (i = 0; i + 1 < object[0].Npts; i+=2) { if (!(finite(x[i]) && finite(y[i]))) continue; @@ -567,5 +579,5 @@ sx2 = x[i+1]*mxi + y[i+1]*mxj + bx; sy2 = x[i+1]*myi + y[i+1]*myj + by; - DrawLine (buffer, sx1, sy1, sx2, sy2); + bDrawClipLine (buffer, sx1, sy1, sx2, sy2, X0, Y0, X1, Y1); } } Index: trunk/Ohana/src/libautocode/Makefile.Targets =================================================================== --- trunk/Ohana/src/libautocode/Makefile.Targets (revision 31027) +++ trunk/Ohana/src/libautocode/Makefile.Targets (revision 31160) @@ -60,4 +60,5 @@ $(ASRC)/cmf-ps1-v2.$(ARCH).o \ $(ASRC)/cmf-ps1-v3.$(ARCH).o \ +$(ASRC)/cmf-ps1-sv1.$(ARCH).o \ $(ASRC)/cmf-smpdata.$(ARCH).o \ $(ASRC)/getstar-ps1-dev-0.$(ARCH).o \ @@ -131,4 +132,5 @@ $(AINC)/cmf-ps1-v2.h \ $(AINC)/cmf-ps1-v3.h \ +$(AINC)/cmf-ps1-sv1.h \ $(AINC)/cmf-smpdata.h \ $(AINC)/getstar-ps1-dev-0.h \ Index: trunk/Ohana/src/libautocode/def/cmf-ps1-sv1.d =================================================================== --- trunk/Ohana/src/libautocode/def/cmf-ps1-sv1.d (revision 31160) +++ trunk/Ohana/src/libautocode/def/cmf-ps1-sv1.d (revision 31160) @@ -0,0 +1,127 @@ +# name of structure type +STRUCT CMF_PS1_SV1 +EXTNAME CMF_PS1_SV1 +TYPE BINTABLE +SIZE 200 + +# elements of data structure / FITS table +FIELD detID, IPP_IDET, unsigned int, detection ID +FIELD X, X_PSF, float, x coord, pixels +FIELD Y, Y_PSF, float, y coord, pixels +FIELD dX, X_PSF_SIG, float, x coord error, pixels + +FIELD dY, Y_PSF_SIG, float, y coord error, pixels +FIELD posangle, POSANGLE, float, Posangle at source, degrees +FIELD pltscale, PLTSCALE, float, Plate Scale at source, arcsec/pixel +FIELD M, PSF_INST_MAG, float, inst mags, mags + +FIELD dM, PSF_INST_MAG_SIG, float, inst mag error, mags +FIELD Flux, PSF_INST_FLUX, float, psf flux, counts +FIELD dFlux, PSF_INST_FLUX_SIG, float, psf flux error, counts +FIELD Map, AP_MAG_STANDARD, float, standard aperture mag, mags + +FIELD MapRaw, AP_MAG_RAW, float, raw aperture mag, mags +FIELD apRadius, AP_MAG_RADIUS, float, radius used for fit, pixels +FIELD Mcalib, CAL_PSF_MAG, float, calibrated psf mag, mags +FIELD dMcal, CAL_PSF_MAG_SIG, float, zero point scatter, mags + +# NOTE: RA & DEC (both double) need to be on an 8-byte boundary... +FIELD RA, RA_PSF, double, PSF RA coord, degrees +FIELD DEC, DEC_PSF, double, PSF DEC coord, degrees + +FIELD Mpeak, PEAK_FLUX_AS_MAG, float, peak flux as a mag, mags +FIELD sky, SKY, float, sky flux, cnts/sec +FIELD dSky, SKY_SIG, float, sky flux error, cnts/sec +FIELD psfChisq, PSF_CHISQ, float, psf fit chisq + +FIELD crNsigma, CR_NSIGMA, float, Nsigma deviations from PSF to CF +FIELD extNsigma, EXT_NSIGMA, float, Nsigma deviations from PSF to EXT +FIELD fx, PSF_MAJOR, float, psf fit major axis, pixels +FIELD fy, PSF_MINOR, float, psf fit minor axis, pixels + +FIELD df, PSF_THETA, float, ellipse angle, degrees +FIELD psfQual, PSF_QF, float, quality factor +FIELD psfQualPerfect, PSF_QF_PERFECT, float, quality factor perfect +FIELD psfNdof, PSF_NDOF, int, psf degrees of freedom + +FIELD psfNpix, PSF_NPIX, int, psf number of pixels +FIELD Mxx, MOMENTS_XX, float, second moment X, pixels^2 +FIELD Mxy, MOMENTS_XY, float, second moment Y, pixels^2 +FIELD Myy, MOMENTS_YY, float, second moment XY, pixels^2 + +FIELD M3c, MOMENTS_M3C, float, third moment cos(t), pixels^3 +FIELD M3s, MOMENTS_M3S, float, third moment sin(t), pixels^3 +FIELD M4c, MOMENTS_M4C, float, fourth moment cos(t), pixels^4 +FIELD M4s, MOMENTS_M4S, float, fourth moment sin(t), pixels^4 + +FIELD Mr1, MOMENTS_R1, float, first radial moment, pixels +FIELD Mrh, MOMENTS_RH, float, half radial moment, pixels^1/2 +FIELD kronFlux, KRON_FLUX, float, kron flux, counts +FIELD kronFluxErr, KRON_FLUX_ERR, float, kron flux error, counts + +FIELD kronInner, KRON_FLUX_INNER, float, kron flux 1 +//# include +//# include + +# include +# include +# include # include # include -# include -# include # include # include @@ -151,5 +154,5 @@ int KapaSetImageCoords (int fd, Coords *coords); int KapaGetImageCoords (int fd, Coords *coords); -int KapaGetImageRange (int fd, double *Xmin, double *Xmax, double *Ymin, double *Ymax); +int KapaGetImageRange (int fd, double *Xmin, double *Xmax, double *Ymin, double *Ymax, int *dX, int *dY); /* KiiOverlay.c */ Index: trunk/Ohana/src/libkapa/src/KiiPicture.c =================================================================== --- trunk/Ohana/src/libkapa/src/KiiPicture.c (revision 31027) +++ trunk/Ohana/src/libkapa/src/KiiPicture.c (revision 31160) @@ -203,10 +203,10 @@ } -int KapaGetImageRange (int fd, double *Xmin, double *Xmax, double *Ymin, double *Ymax) { +int KapaGetImageRange (int fd, double *Xmin, double *Xmax, double *Ymin, double *Ymax, int *dX, int *dY) { /* tell kapa to look for the incoming image */ KiiSendCommand (fd, 4, "GIMR"); - KiiScanMessage (fd, "%lf %lf %lf %lf", Xmin, Xmax, Ymin, Ymax); + KiiScanMessage (fd, "%lf %lf %lf %lf %d %d", Xmin, Xmax, Ymin, Ymax, dX, dY); KiiWaitAnswer (fd, "DONE"); Index: trunk/Ohana/src/libohana/include/ohana.h =================================================================== --- trunk/Ohana/src/libohana/include/ohana.h (revision 31027) +++ trunk/Ohana/src/libohana/include/ohana.h (revision 31160) @@ -155,7 +155,10 @@ # define OFF_T_FMT "%jd" # endif -# ifdef _DARWIN_C_SOURCE +// # ifdef _DARWIN_C_SOURCE +// # define OFF_T_FMT "%lld" +// # endif +# ifdef darwin_x86 # define OFF_T_FMT "%lld" -#endif +# endif # ifndef OFF_T_FMT # define OFF_T_FMT "%ld" Index: trunk/Ohana/src/libohana/test/typetest.c =================================================================== --- trunk/Ohana/src/libohana/test/typetest.c (revision 31027) +++ trunk/Ohana/src/libohana/test/typetest.c (revision 31160) @@ -68,6 +68,12 @@ fprintf (stderr, "STDC_VERSION is not set\n"); # endif + } - } + // 3) have we defined OFF_T_FMT correctly? + off_t big_value; + + big_value = 10; + fprintf (stderr, "this is a bit int: "OFF_T_FMT" n'est pas?\n", big_value); + exit (status); } Index: trunk/Ohana/src/opihi/cmd.astro/region.c =================================================================== --- trunk/Ohana/src/opihi/cmd.astro/region.c (revision 31027) +++ trunk/Ohana/src/opihi/cmd.astro/region.c (revision 31160) @@ -5,5 +5,5 @@ double Ra, Dec, Radius; float dx, dy; - int N, kapa, NoClear; + int N, kapa, NoClear, dXpix, dYpix; char *name; Graphdata graphmode; @@ -27,5 +27,5 @@ remove_argument (N, &argc, argv); KapaGetImageCoords (kapa, &graphmode.coords); - KapaGetImageRange (kapa, &graphmode.xmin, &graphmode.xmax, &graphmode.ymax, &graphmode.ymin); + KapaGetImageRange (kapa, &graphmode.xmin, &graphmode.xmax, &graphmode.ymax, &graphmode.ymin, &dXpix, &dYpix); set_variable ("XMIN", graphmode.xmin); Index: trunk/Ohana/src/opihi/cmd.data/densify.c =================================================================== --- trunk/Ohana/src/opihi/cmd.data/densify.c (revision 31027) +++ trunk/Ohana/src/opihi/cmd.data/densify.c (revision 31160) @@ -1,8 +1,10 @@ # include "data.h" + +# define CHECKVAL(ARG) if (!isfinite(ARG)) { gprint (GP_ERR, "illegal value for %s: %f\n", #ARG, ARG); return (FALSE); } int densify (int argc, char **argv) { - int i, Nx, Ny, Xb, Yb, Normalize, N; - float Xmin, Xmax, dX, Ymin, Ymax, dY; + int i, Nx, Ny, Xb, Yb, Normalize, N, Xpix, Ypix, good, UseGraph; + double Xmin, Xmax, dX, Ymin, Ymax, dY; float *val; Buffer *bf; @@ -16,6 +18,14 @@ } - if (argc != 10) { + UseGraph = FALSE; + if ((N = get_argument (argc, argv, "-graph"))) { + remove_argument (N, &argc, argv); + UseGraph = TRUE; + } + + good = UseGraph ? (argc == 4) : (argc == 10); + if (!good) { gprint (GP_ERR, "USAGE: densify buffer x y Xmin Xmax dX Ymin Ymax dY\n"); + gprint (GP_ERR, " OR: densify buffer x y -graph\n"); return (FALSE); } @@ -30,11 +40,32 @@ REQUIRE_VECTOR_FLT (vy, FALSE); - Xmin = atof (argv[4]); - Xmax = atof (argv[5]); - dX = atof (argv[6]); + if (UseGraph) { + int kapa; + Graphdata graphmode; + if (!GetGraph (&graphmode, &kapa, NULL)) return (FALSE); + KapaGetImageRange (kapa, &Xmin, &Xmax, &Ymax, &Ymin, &Xpix, &Ypix); + Xmax = graphmode.xmax; + Xmin = graphmode.xmin; + Ymax = graphmode.ymax; + Ymin = graphmode.ymin; + dX = (Xmax - Xmin) / (Xpix - 1); + dY = (Ymax - Ymin) / (Ypix - 1); + } else { + Xmin = atof (argv[4]); + Xmax = atof (argv[5]); + dX = atof (argv[6]); - Ymin = atof (argv[7]); - Ymax = atof (argv[8]); - dY = atof (argv[9]); + Ymin = atof (argv[7]); + Ymax = atof (argv[8]); + dY = atof (argv[9]); + } + + CHECKVAL(Xmin); + CHECKVAL(Xmax); + CHECKVAL(dX); + + CHECKVAL(Ymin); + CHECKVAL(Ymax); + CHECKVAL(dY); Nx = (Xmax - Xmin) / dX + 1; Index: trunk/Ohana/src/opihi/cmd.data/fit1d.c =================================================================== --- trunk/Ohana/src/opihi/cmd.data/fit1d.c (revision 31027) +++ trunk/Ohana/src/opihi/cmd.data/fit1d.c (revision 31160) @@ -42,5 +42,5 @@ if (argc != 4) { - gprint (GP_ERR, "USAGE: fit x y order [-dy wt] [-quiet/-q] [-clip Nsigma Niter]\n"); + gprint (GP_ERR, "USAGE: fit1d x y order [-dy wt] [-quiet/-q] [-clip Nsigma Niter]\n"); return (FALSE); } Index: trunk/Ohana/src/opihi/cmd.data/fit2d.c =================================================================== --- trunk/Ohana/src/opihi/cmd.data/fit2d.c (revision 31027) +++ trunk/Ohana/src/opihi/cmd.data/fit2d.c (revision 31160) @@ -48,5 +48,5 @@ if (argc != 5) { - gprint (GP_ERR, "USAGE: fit x y z order [-dz wt]\n"); + gprint (GP_ERR, "USAGE: fit2d x y z order [-dz wt] [-quiet/-q] [-clip Nsigma Niter]\n"); return (FALSE); } Index: trunk/Ohana/src/opihi/cmd.data/histogram.c =================================================================== --- trunk/Ohana/src/opihi/cmd.data/histogram.c (revision 31027) +++ trunk/Ohana/src/opihi/cmd.data/histogram.c (revision 31160) @@ -61,5 +61,4 @@ opihi_int *V = xvec[0].elements.Int; for (i = 0; i < xvec[0].Nelements; i++, V++) { - if (isnan(*V)) continue; bin = MIN (MAX (0, (*V - start) / delta), Nbins - 1); OUT[bin]++; Index: trunk/Ohana/src/opihi/cmd.data/limits.c =================================================================== --- trunk/Ohana/src/opihi/cmd.data/limits.c (revision 31027) +++ trunk/Ohana/src/opihi/cmd.data/limits.c (revision 31160) @@ -3,5 +3,5 @@ int limits (int argc, char **argv) { - int N, APPLY; + int N, APPLY, dX, dY; int kapa; char *name; @@ -28,5 +28,5 @@ if ((N = get_argument (argc, argv, "-image"))) { remove_argument (N, &argc, argv); - KapaGetImageRange (kapa, &graphmode.xmin, &graphmode.xmax, &graphmode.ymax, &graphmode.ymin); + KapaGetImageRange (kapa, &graphmode.xmin, &graphmode.xmax, &graphmode.ymax, &graphmode.ymin, &dX, &dY); set_variable ("XMIN", graphmode.xmin); @@ -34,4 +34,12 @@ set_variable ("YMIN", graphmode.ymin); set_variable ("YMAX", graphmode.ymax); + + set_variable ("KAPA_XMIN", graphmode.xmin); + set_variable ("KAPA_XMAX", graphmode.xmax); + set_variable ("KAPA_YMIN", graphmode.ymin); + set_variable ("KAPA_YMAX", graphmode.ymax); + + set_variable ("KAPA_XPIX", dX); + set_variable ("KAPA_YPIX", dY); // if (!NoClear) KapaClearSections (kapa); Index: trunk/Ohana/src/opihi/cmd.data/section.c =================================================================== --- trunk/Ohana/src/opihi/cmd.data/section.c (revision 31027) +++ trunk/Ohana/src/opihi/cmd.data/section.c (revision 31160) @@ -68,5 +68,5 @@ if ((argc == 1) && (action == NONE)) { KapaGetSection (kapa, "*"); - gprint (GP_ERR, "USAGE: section name [x y dx dy] [options]\n"); + gprint (GP_ERR, "USAGE: section name [x y dx dy] or [options]\n"); gprint (GP_ERR, "OPTIONS: -list : show properties of all sections\n"); gprint (GP_ERR, " -up : move section up in display stack\n"); @@ -76,7 +76,15 @@ gprint (GP_ERR, " -imtool (position) : set location of image zoom / status box\n"); gprint (GP_ERR, " (position may be: -x, +x, -y, +y, none)\n"); + gprint (GP_ERR, " -image : define section upper-right corner so image fills the section\n"); + gprint (GP_ERR, " -bg (color) : set background color (see style -c help for color choices)\n"); + return (TRUE); } + if ((argc != 2) && (action != NONE)) { + gprint (GP_ERR, "can only use dash options without numbers\n"); + return (FALSE); + } + if (argc == 2) { /* select / show section */ Index: trunk/Ohana/src/opihi/dvo/dbExtractImages.c =================================================================== --- trunk/Ohana/src/opihi/dvo/dbExtractImages.c (revision 31027) +++ trunk/Ohana/src/opihi/dvo/dbExtractImages.c (revision 31160) @@ -59,4 +59,5 @@ time_t t; dbValue value; + off_t Nmosaic; value.Flt = NAN; @@ -236,4 +237,21 @@ value.Flt = image[N].fwhm_y / 25.0; break; + + case IMAGE_FWHM_MEDIAN: + if (!(Nmosaic = FindMosaicForImage (image, Nimage, N))) return value; + Nmosaic --; // XXX kind of a hack: FindMosaicForImage returns 0 or the mosaic seq number + 1 + value.Flt = (image[Nmosaic].fwhm_x + image[Nmosaic].fwhm_y) / 50.0; + break; + case IMAGE_FWHM_MAJ_MEDIAN: + if (!(Nmosaic = FindMosaicForImage (image, Nimage, N))) return value; + Nmosaic --; // XXX kind of a hack: FindMosaicForImage returns 0 or the mosaic seq number + 1 + value.Flt = image[Nmosaic].fwhm_x / 25.0; + break; + case IMAGE_FWHM_MIN_MEDIAN: + if (!(Nmosaic = FindMosaicForImage (image, Nimage, N))) return value; + Nmosaic --; // XXX kind of a hack: FindMosaicForImage returns 0 or the mosaic seq number + 1 + value.Flt = image[Nmosaic].fwhm_y / 25.0; + break; + case IMAGE_TRATE: value.Flt = image[N].trate / 10000.0; Index: trunk/Ohana/src/opihi/dvo/dbFields.c =================================================================== --- trunk/Ohana/src/opihi/dvo/dbFields.c (revision 31027) +++ trunk/Ohana/src/opihi/dvo/dbFields.c (revision 31160) @@ -395,7 +395,16 @@ if (!strcasecmp (fieldName, "cerror" )) ESCAPE (IMAGE_CERROR, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "FWHM" )) ESCAPE (IMAGE_FWHM, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "FWHM_MAJ" )) ESCAPE (IMAGE_FWHM_MAJ, MAG_NONE, OPIHI_FLT); - if (!strcasecmp (fieldName, "FWHM_MIN" )) ESCAPE (IMAGE_FWHM_MIN, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "FWHM" )) ESCAPE (IMAGE_FWHM, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "FWHM_MAJ" )) ESCAPE (IMAGE_FWHM_MAJ, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "FWHM_MIN" )) ESCAPE (IMAGE_FWHM_MIN, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "FWHM_MAJOR" )) ESCAPE (IMAGE_FWHM_MAJ, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "FWHM_MINOR" )) ESCAPE (IMAGE_FWHM_MIN, MAG_NONE, OPIHI_FLT); + + if (!strcasecmp (fieldName, "FWHM_MEDIAN" )) ESCAPE (IMAGE_FWHM_MEDIAN, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "FWHM_MAJ_MEDIAN")) ESCAPE (IMAGE_FWHM_MAJ_MEDIAN, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "FWHM_MIN_MEDIAN")) ESCAPE (IMAGE_FWHM_MIN_MEDIAN, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "FWHM_MAJOR_MEDIAN")) ESCAPE (IMAGE_FWHM_MAJ_MEDIAN, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "FWHM_MINOR_MEDIAN")) ESCAPE (IMAGE_FWHM_MIN_MEDIAN, MAG_NONE, OPIHI_FLT); + if (!strcasecmp (fieldName, "trate" )) ESCAPE (IMAGE_TRATE, MAG_NONE, OPIHI_FLT); Index: trunk/Ohana/src/opihi/dvo/imextract.c =================================================================== --- trunk/Ohana/src/opihi/dvo/imextract.c (revision 31027) +++ trunk/Ohana/src/opihi/dvo/imextract.c (revision 31160) @@ -197,7 +197,17 @@ gprint (GP_ERR, " cerror : astrometric scatter\n"); - gprint (GP_ERR, " FWHM : mean fwhm of exposure\n"); - gprint (GP_ERR, " FWHM_MAJ : fwhm of major axis\n"); - gprint (GP_ERR, " FWHM_MIN : fwhm of minor axis\n"); + gprint (GP_ERR, " -- Note: the follow FWHM are from the PSF model --\n"); + gprint (GP_ERR, " FWHM : mean fwhm of chip\n"); + gprint (GP_ERR, " FWHM_MAJ : fwhm of chip (major axis)\n"); + gprint (GP_ERR, " FWHM_MIN : fwhm of chip (minor axis)\n"); + gprint (GP_ERR, " FWHM_MAJOR : fwhm of chip (major axis)\n"); + gprint (GP_ERR, " FWHM_MININ : fwhm of chip (minor axis)\n"); + + gprint (GP_ERR, " FWHM_MEDIAN : median fwhm of exposure\n"); + gprint (GP_ERR, " FWHM_MAJ_MEDIAN : median fwhm of major axis\n"); + gprint (GP_ERR, " FWHM_MIN_MEDIAN : median fwhm of minor axis\n"); + gprint (GP_ERR, " FWHM_MAJOR_MEDIAN : median fwhm of major axis\n"); + gprint (GP_ERR, " FWHM_MININ_MEDIAN : median fwhm of minor axis\n"); + gprint (GP_ERR, " trate : tracking rate for TDI images\n"); Index: trunk/Ohana/src/opihi/dvo/skycoverage.c =================================================================== --- trunk/Ohana/src/opihi/dvo/skycoverage.c (revision 31027) +++ trunk/Ohana/src/opihi/dvo/skycoverage.c (revision 31160) @@ -5,5 +5,5 @@ int skycoverage (int argc, char **argv) { - int WITH_MOSAIC, SOLO_MOSAIC; + int WITH_MOSAIC, SOLO_MOSAIC, ShowDensity; off_t i, Nimage; int N, status, TimeSelect, ByName, xs, ys; @@ -60,4 +60,10 @@ } + ShowDensity = FALSE; + if ((N = get_argument (argc, argv, "-density"))) { + remove_argument (N, &argc, argv); + ShowDensity = TRUE; + } + ByName = FALSE; if ((N = get_argument (argc, argv, "-name"))) { @@ -141,5 +147,5 @@ if (argc != 3) { gprint (GP_ERR, "USAGE: skycoverage (buffer) (Npts)\n"); - gprint (GP_ERR, " options: [-scale pixscale] [-center ra dec] [-size Nx Nx] [-proj projection] [-time start range] [-trange start stop] [-name name] [-photcode name] [+mosaic] [-mosaic]\n"); + gprint (GP_ERR, " options: [-scale pixscale] [-center ra dec] [-size Nx Nx] [-proj projection] [-time start range] [-trange start stop] [-name name] [-photcode name] [+mosaic] [-mosaic] [-density]\n"); gprint (GP_ERR, " (buffer) saves bitmapped image\n"); gprint (GP_ERR, " (Npts) gives the number of test points per image in each dimension\n"); @@ -148,4 +154,5 @@ gprint (GP_ERR, " -size (Nx) (Ny) : specifies the size of the image [360/scale, 180/scale]\n"); gprint (GP_ERR, " -proj (projection) : specifies the projection choice [AIT]\n"); + gprint (GP_ERR, " -density : create image with relative density (else binary on/off)\n"); gprint (GP_ERR, " note: we need 64800 / (pixscale)^2 pixels to represent the sky\n"); return (FALSE); @@ -186,5 +193,4 @@ V = (float *)buf[0].matrix.buffer; - bzero (V, Nx*Ny*sizeof(float)); for (ys = 0; ys < Ny; ys++) { @@ -194,5 +200,7 @@ status &= (r <= 360); if (status) { - V[ys*Nx + xs] = 2; + V[ys*Nx + xs] = ShowDensity ? 0 : 2; + } else { + V[ys*Nx + xs] = ShowDensity ? -1 : 0; } } @@ -249,5 +257,9 @@ xs = (int)Xs; ys = (int)Ys; - V[ys*Nx + xs] = 1; + if (ShowDensity) { + V[ys*Nx + xs] += 1; + } else { + V[ys*Nx + xs] = 1; + } } } Index: trunk/Ohana/src/opihi/include/dvoshell.h =================================================================== --- trunk/Ohana/src/opihi/include/dvoshell.h (revision 31027) +++ trunk/Ohana/src/opihi/include/dvoshell.h (revision 31160) @@ -179,4 +179,5 @@ IMAGE_TIME, IMAGE_FWHM, + IMAGE_FWHM_MEDIAN, IMAGE_EXPTIME, IMAGE_NSTAR, @@ -200,4 +201,6 @@ IMAGE_FWHM_MAJ, IMAGE_FWHM_MIN, + IMAGE_FWHM_MAJ_MEDIAN, + IMAGE_FWHM_MIN_MEDIAN, IMAGE_TRATE, IMAGE_IMAGE_ID, Index: trunk/Ohana/src/opihi/lib.data/graphtools.c =================================================================== --- trunk/Ohana/src/opihi/lib.data/graphtools.c (revision 31027) +++ trunk/Ohana/src/opihi/lib.data/graphtools.c (revision 31160) @@ -59,4 +59,9 @@ set_variable ("YMIN", graphmode[0].ymin); set_variable ("YMAX", graphmode[0].ymax); + + set_variable ("KAPA_XMIN", graphmode[0].xmin); + set_variable ("KAPA_XMAX", graphmode[0].xmax); + set_variable ("KAPA_YMIN", graphmode[0].ymin); + set_variable ("KAPA_YMAX", graphmode[0].ymax); } @@ -95,4 +100,9 @@ SetGraph (graphmode); + set_variable ("KAPA_XMIN", graphmode[0].xmin); + set_variable ("KAPA_XMAX", graphmode[0].xmax); + set_variable ("KAPA_YMIN", graphmode[0].ymin); + set_variable ("KAPA_YMAX", graphmode[0].ymax); + set_variable ("XMIN", graphmode[0].xmin); set_variable ("XMAX", graphmode[0].xmax); Index: trunk/Ohana/src/opihi/lib.shell/VectorOps.c =================================================================== --- trunk/Ohana/src/opihi/lib.shell/VectorOps.c (revision 31027) +++ trunk/Ohana/src/opihi/lib.shell/VectorOps.c (revision 31160) @@ -218,5 +218,5 @@ } free (vec[0].elements.Int); - vec[0].elements.Flt = vo; + vec[0].elements.Flt = temp; vec[0].type = OPIHI_FLT; } else { @@ -229,5 +229,5 @@ } free (vec[0].elements.Flt); - vec[0].elements.Int = vo; + vec[0].elements.Int = temp; vec[0].type = OPIHI_INT; } Index: trunk/Ohana/src/relastro/include/relastro.h =================================================================== --- trunk/Ohana/src/relastro/include/relastro.h (revision 31027) +++ trunk/Ohana/src/relastro/include/relastro.h (revision 31160) @@ -89,4 +89,5 @@ char GSCFILE[256]; char CATDIR[256]; +char *HIGH_SPEED_DIR; char CATMODE[16]; /* raw, mef, split, mysql */ char CATFORMAT[16]; /* internal, elixir, loneos, panstarrs */ @@ -95,5 +96,5 @@ double SIGMA_LIM; -int SRC_MEAS_TOOFEW; //catalog objects wich fewer detections then this are ignored +int SRC_MEAS_TOOFEW; //catalog objects wich fewer detections then this are ignored double MIN_ERROR; Index: trunk/Ohana/src/relastro/src/UpdateObjects.c =================================================================== --- trunk/Ohana/src/relastro/src/UpdateObjects.c (revision 31027) +++ trunk/Ohana/src/relastro/src/UpdateObjects.c (revision 31160) @@ -92,4 +92,8 @@ } # endif + + if (catalog[i].average[j].Nmeasure == 0) { + continue; + } N = 0; Index: trunk/Ohana/src/relastro/src/args.c =================================================================== --- trunk/Ohana/src/relastro/src/args.c (revision 31027) +++ trunk/Ohana/src/relastro/src/args.c (revision 31160) @@ -39,4 +39,6 @@ remove_argument (N, &argc, argv); RADIUS = atof(argv[N]); + remove_argument (N, &argc, argv); + HIGH_SPEED_DIR = strcreate(argv[N]); remove_argument (N, &argc, argv); } Index: trunk/Ohana/src/relastro/src/high_speed_objects.c =================================================================== --- trunk/Ohana/src/relastro/src/high_speed_objects.c (revision 31027) +++ trunk/Ohana/src/relastro/src/high_speed_objects.c (revision 31160) @@ -22,7 +22,7 @@ int zcode, zNsec, ycode, yNsec, jcode, jNsec, hcode, hNsec, kcode, kNsec, USNO_R, USNO_N, Nsecfilt; char filename[1024]; - char outdir[]="/data/ipp022.0/ndeacon/hispeedzy"; Noff = strlen(CATDIR); - sprintf (filename, "%s/%s", outdir, &catalog[0].filename[Noff]); + sprintf (filename, "%s/%s", HIGH_SPEED_DIR, &catalog[0].filename[Noff]); + printf("%s\n",filename); dvo_catalog_init(&catalog1, TRUE); /*initialise new catalogue*/ catalog1.filename = strcreate(filename); @@ -121,8 +121,4 @@ foundA = FALSE; for (j = 0; !foundA && (j < catalog[0].average[i].Nmeasure); j++, m++) { - if((catalog[0].average[i].R>204.1923)&&(catalog[0].average[i].R<204.1924)&&(catalog[0].average[i].D>11.376)&&(catalog[0].average[i].D<11.377)) - { - printf("Hello"); - } if (MeasMatchesPhotcode(&catalog[0].measure[m], photcodesGroupA, NphotcodesGroupA)) { @@ -135,9 +131,5 @@ foundB = FALSE; for (j = 0; !foundB && (j < catalog[0].average[i].Nmeasure); j++, m++) { - - if((catalog[0].average[i].R>204.192)&&(catalog[0].average[i].R<204.1925)&&(catalog[0].average[i].D>11.376)&&(catalog[0].average[i].D<11.377)) - { - printf("Hello"); - } + if (MeasMatchesPhotcode(&catalog[0].measure[m], photcodesGroupB, NphotcodesGroupB)) { @@ -158,9 +150,5 @@ if (foundA && !foundB) { // average-based tests: - - if((catalog[0].average[i].R>204.1923)&&(catalog[0].average[i].R<204.1924)&&(catalog[0].average[i].D>11.376)&&(catalog[0].average[i].D<11.377)) - { - printf("Hello"); - } + valid = TRUE; @@ -204,8 +192,4 @@ // average-based tests: - if((catalog[0].average[i].R>204.192)&&(catalog[0].average[i].R<204.193)&&(catalog[0].average[i].D>11.372)&&(catalog[0].average[i].D<11.373)) - { - printf("Hello"); - } valid = TRUE; valid &= ((catalog[0].average[i].flags & 0x01000000) == 0); @@ -214,6 +198,6 @@ valid &= ((catalog[0].secfilt[i*Nsecfilt + jNsec].M < 1.0)||(isnan(catalog[0].secfilt[i*Nsecfilt + jNsec].M))); - valid &= (catalog[0].secfilt[i*Nsecfilt + yNsec].Nused > 1); - valid &= (catalog[0].secfilt[i*Nsecfilt + yNsec].dM < 0.2); + valid &= ((catalog[0].secfilt[i*Nsecfilt + yNsec].Nused > 1)||(catalog[0].secfilt[i*Nsecfilt + zNsec].Nused > 1)); + valid &= ((catalog[0].secfilt[i*Nsecfilt + yNsec].dM < 0.2)||(catalog[0].secfilt[i*Nsecfilt + zNsec].dM < 0.2)); /*if ((catalog[0].secfilt[i*Nsecfilt + zNsec].M < 1.0) || (catalog[0].secfilt[i*Nsecfilt + zNsec].Nused == 1)) { Index: trunk/Ohana/src/relastro/src/select_images.c =================================================================== --- trunk/Ohana/src/relastro/src/select_images.c (revision 31027) +++ trunk/Ohana/src/relastro/src/select_images.c (revision 31160) @@ -144,6 +144,7 @@ } + // this adds 1.3 sec for 3M images if (!FindMosaicForImage (timage, Ntimage, i)) { - fprintf (stderr, "cannot find mosaic for "OFF_T_FMT"\n", i); + fprintf (stderr, "cannot find mosaic for "OFF_T_FMT"\n", i); continue; } @@ -165,10 +166,9 @@ found = FALSE; - /* transform corners to ra,dec */ + /* transform corners to ra,dec -- costs ~3sec for 3M images (pikake) */ double RminImage = 360.0; double RmaxImage = 0.0; double DminImage = +90.0; double DmaxImage = -90.0; - // int leftside = FALSE; for (j = 0; j < 5; j++) { XY_to_RD (&Ri[j], &Di[j], Xi[j], Yi[j], &timage[i].coords); @@ -197,5 +197,5 @@ if (USE_BASIC_CHECK) goto found_it; - // RA(nStart) is guaranteed to be < RminImage: + // RA(nStart) is guaranteed to be < RminImage: -- costs 0.5sec for 3M images nStart = getRegionStartByRA (RminImage, RmaxSky, skylist[0].Nregions); @@ -258,6 +258,6 @@ } MARKTIME("finish image selection: %f sec\n", dtime); - - if (VERBOSE) fprintf (stderr, "found "OFF_T_FMT" images\n", nimage); + + if (VERBOSE) fprintf (stderr, "found "OFF_T_FMT" images\n", nimage); REALLOCATE (image, Image, MAX (nimage, 1)); Index: trunk/Ohana/src/relphot/src/args.c =================================================================== --- trunk/Ohana/src/relphot/src/args.c (revision 31027) +++ trunk/Ohana/src/relphot/src/args.c (revision 31160) @@ -78,4 +78,5 @@ remove_argument (N, &argc, argv); PLOTDELAY = 1e6*atof(argv[N]); + PLOTSTUFF = TRUE; // always turn on plotting if i request a plot delay remove_argument (N, &argc, argv); } Index: trunk/Ohana/src/relphot/src/bcatalog.c =================================================================== --- trunk/Ohana/src/relphot/src/bcatalog.c (revision 31027) +++ trunk/Ohana/src/relphot/src/bcatalog.c (revision 31160) @@ -7,5 +7,5 @@ off_t NAVERAGE, NMEASURE, Naverage, Nmeasure, Nm; float mag; - int Ncode, Ntime, Ndophot, Nmag, Nsigma, Nimag, Nfew; + int Ncode, Ntime, Ndophot, Nmag, Nsigma, Nimag, Nfew, Ngalaxy, Npsfqf; // XXX PhotNsec as a global is a bad idea; either get it from catalog @@ -23,5 +23,5 @@ Nmeasure = Naverage = 0; - Ncode = Ntime = Ndophot = Nmag = Nsigma = Nimag = Nfew = 0; + Ncode = Ntime = Ndophot = Nmag = Nsigma = Nimag = Nfew = Npsfqf = Ngalaxy = 0; /* exclude stars not in range or with too few measurements */ @@ -44,4 +44,6 @@ Nm = 0; + int nEXT = 0; + int nPSF = 0; for (j = 0; j < catalog[0].average[i].Nmeasure; j++) { @@ -62,4 +64,16 @@ // XXX chnage this to select by bitflags if (DophotSelect && ((catalog[0].measure[offset].photFlags >> 16) != DophotValue)) { Ndophot ++; continue; } + + // skip garbage measurements + if (catalog[0].measure[offset].psfQual < 0.85) { Npsfqf ++; continue; } + + // check for galaxies + if (!isnan(catalog[0].measure[offset].Map)) { + if (catalog[0].measure[offset].M - catalog[0].measure[offset].Map > 0.15) { + nEXT ++; + } else { + nPSF ++; + } + } /* select measurements by mag limit */ @@ -95,4 +109,11 @@ } + // skip object if it is likely to be a galaxy + if (nEXT >= nPSF) { + Nmeasure -= Nm; + Ngalaxy ++; + continue; + } + // XXXX test : what checks do I need to make elsewhere to avoid problems here? if (Nm <= STAR_TOOFEW) { /* enough measurements in band? */ @@ -121,6 +142,6 @@ fprintf (stderr, "using "OFF_T_FMT" stars ("OFF_T_FMT" measures) of "OFF_T_FMT" for catalog %s\n", subcatalog[0].Naverage, subcatalog[0].Nmeasure, i, catalog[0].filename); - fprintf (stderr, "rejections: %d code, %d time, %d dophot, %d mag, %d sigma, %d imag, %d few\n", - Ncode, Ntime, Ndophot, Nmag, Nsigma, Nimag, Nfew); + fprintf (stderr, "rejections: %d code, %d time, %d dophot, %d mag, %d sigma, %d imag, %d few, %d psfqf, %d galaxies\n", + Ncode, Ntime, Ndophot, Nmag, Nsigma, Nimag, Nfew, Npsfqf, Ngalaxy); } return (TRUE); Index: trunk/Ohana/src/relphot/src/load_images.c =================================================================== --- trunk/Ohana/src/relphot/src/load_images.c (revision 31027) +++ trunk/Ohana/src/relphot/src/load_images.c (revision 31160) @@ -7,4 +7,6 @@ fprintf (stderr, MSG, __VA_ARGS__); } +// This function generates a subset of the images based on selections. Input db has already +// been loaded with the raw fits table data SkyList *load_images (FITS_DB *db, char *regionName, SkyRegion *region, int RegionSelect) { @@ -32,5 +34,5 @@ } - // convert database table to internal structure + // convert database table to internal structure (binary to Image) image = gfits_table_get_Image (&db[0].ftable, &Nimage, &db[0].swapped); if (!image) { @@ -44,10 +46,15 @@ MARKTIME("selected images: %f sec\n", dtime); + // generate db->vtable from db->ftable based on the selection + // XXX does this simply duplicate the memory needlessly? we recreate these lines + // in reload_images. If we had saved the line numbers, we could avoid this gfits_vtable_from_ftable (&db[0].ftable, &db[0].vtable, LineNumber, Nsubset); MARKTIME("converted ftable to vtable: %f sec\n", dtime); + // save the subset of images in the static reference in ImageOps, set up indexes initImages (subset, Nsubset); MARKTIME("init images: %f sec\n", dtime); + // match chips to mosaics (if applicable) initMosaics (subset, Nsubset); MARKTIME("init mosaics: %f sec\n", dtime); Index: trunk/Ohana/src/relphot/src/select_images.c =================================================================== --- trunk/Ohana/src/relphot/src/select_images.c (revision 31027) +++ trunk/Ohana/src/relphot/src/select_images.c (revision 31160) @@ -33,5 +33,5 @@ struct timeval start, stop; - double RmaxSkyRegion, RminSkyRegion, DminSkyRegion, DmaxSkyRegion, RmidSkyRegion; + double RmaxSkyRegion, RminSkyRegion, RmidSkyRegion, DminSkyRegion, DmaxSkyRegion; double *RmaxSky; @@ -131,13 +131,21 @@ } - /* 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; + /* define image corners - note the DIS images (mosaic phu) are special */ + if (!strcmp(&timage[i].coords.ctype[4], "-DIS")) { + Xi[0] = -0.5*timage[i].NX; Yi[0] = -0.5*timage[i].NY; + Xi[1] = +0.5*timage[i].NX; Yi[1] = -0.5*timage[i].NY; + Xi[2] = +0.5*timage[i].NX; Yi[2] = +0.5*timage[i].NY; + Xi[3] = -0.5*timage[i].NX; Yi[3] = +0.5*timage[i].NY; + Xi[4] = -0.5*timage[i].NX; Yi[4] = -0.5*timage[i].NY; + } else { + 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; + } found = FALSE; - /* transform corners to ra,dec -- costs ~3sec for 3M images */ + /* transform corners to ra,dec -- costs ~3sec for 3M images (pikake) */ double RminImage = 360.0; double RmaxImage = 0.0; @@ -147,5 +155,5 @@ XY_to_RD (&Ri[j], &Di[j], Xi[j], Yi[j], &timage[i].coords); Ri[j] = ohana_normalize_angle_to_midpoint (Ri[j], RmidSkyRegion); - + RminImage = MIN(RminImage, Ri[j]); RmaxImage = MAX(RmaxImage, Ri[j]); @@ -158,5 +166,5 @@ RminImage = tmp - 360.0; } - + // check that this image is even in range of the searched region if (DminImage > DmaxSkyRegion) continue; @@ -182,5 +190,5 @@ tcoords.crval2 = skycoords[m].Dc; - /* transform to ra,dec */ + /* transform corner coords to X,Y in this catalog system */ InRange = TRUE; for (j = 0; (j < 5) && InRange; j++) { @@ -352,21 +360,2 @@ return (Nlo); } - -off_t getRegionStopByRA (double R, double *Rref, off_t Nregions) { - - // use bisection to find the overlapping mosaic - - off_t Nlo, Nhi, N; - - // find the last mosaic before start - Nlo = 0; Nhi = Nregions; - while (Nhi - Nlo > 10) { - N = 0.5*(Nlo + Nhi); - if (Rref[N] < R) { - Nlo = MAX(N, 0); - } else { - Nhi = MIN(N, Nregions); - } - } - return (Nlo); -} Index: trunk/Ohana/src/tools/src/fhead.c =================================================================== --- trunk/Ohana/src/tools/src/fhead.c (revision 31027) +++ trunk/Ohana/src/tools/src/fhead.c (revision 31160) @@ -1,20 +1,42 @@ # include # include +# include int main (int argc, char **argv) { - int N, Extend, Nextend, status; + int N, Extnum, Nextend, status; int i, j; off_t nbytes; Header head; - char *p; + char *p, *CCDKeyword, *Extname, extname[80]; + FILE *f; + off_t Nbytes; + regex_t preg; - Extend = FALSE; + CCDKeyword = NULL; + if ((N = get_argument (argc, argv, "-keyword"))) { + remove_argument (N, &argc, argv); + CCDKeyword = strcreate (argv[N]); + remove_argument (N, &argc, argv); + } + if (CCDKeyword == NULL) { + CCDKeyword = strcreate ("EXTNAME"); + } + + Extnum = FALSE; Nextend = 0; if ((N = get_argument (argc, argv, "-x"))) { - Extend = TRUE; + Extnum = TRUE; remove_argument (N, &argc, argv); Nextend = atoi (argv[N]); remove_argument (N, &argc, argv); + } + + Extname = NULL; + if ((N = get_argument (argc, argv, "-n"))) { + remove_argument (N, &argc, argv); + Extname = strcreate (argv[N]); + remove_argument (N, &argc, argv); + regcomp (&preg, Extname, REG_EXTENDED); } @@ -23,22 +45,52 @@ fprintf (stdout, "------> %s <------\n", argv[i]); - if (Extend) { - status = gfits_read_Xheader (argv[i], &head, Nextend); - } else { - status = gfits_read_header (argv[i], &head); - } - - if (!status) continue; - - for (j = 79; j < head.datasize; j+= 80) { - head.buffer[j] = 10; + status = FALSE; + if (!Extnum && !Extname) { + if (!gfits_read_header (argv[i], &head)) { + continue; + } + } + if (Extnum) { + if (!gfits_read_Xheader (argv[i], &head, Nextend)) { + continue; + } + } + if (Extname) { + /* keep reading headers until we reach the one we want */ + Nextend = 0; + f = fopen (argv[i], "r"); + if (f == NULL) continue; + while (gfits_fread_header (f, &head)) { + /* extract the EXTNAME (or other CCDKeyword) for this component (set to PHU for 0th component) */ + if (!gfits_scan (&head, CCDKeyword, "%s", 1, extname)) { + if (Nextend == 0) { + strcpy (extname, "PHU"); + } else { + strcpy (extname, "UNKNOWN"); + } + } + if (!regexec (&preg, extname, 0, NULL, 0)) { + goto done; + } + + Nbytes = gfits_data_size (&head); + fseeko (f, Nbytes, SEEK_CUR); + Nextend ++; + } + // failed to find the desired header + continue; } - p = gfits_header_field (&head, "END", 1); - nbytes = p - head.buffer; - fwrite (head.buffer, nbytes, 1, stdout); - gfits_free_header (&head); + done: + for (j = 79; j < head.datasize; j+= 80) { + head.buffer[j] = 10; + } + p = gfits_header_field (&head, "END", 1); + nbytes = p - head.buffer; + fwrite (head.buffer, nbytes, 1, stdout); + gfits_free_header (&head); + + } + exit (0); } - exit (0); -} Index: trunk/Ohana/src/tools/src/ftable.c =================================================================== --- trunk/Ohana/src/tools/src/ftable.c (revision 31027) +++ trunk/Ohana/src/tools/src/ftable.c (revision 31160) @@ -3,5 +3,5 @@ # include "inttypes.h" -char *print_table_row (char *row, Header *header); +int print_table_rows (FTable *table, int start, int Nrows); FILE *load_extension (char *file, int Nextend, char *Extname, Header *header); void print_column (FTable *table, int Column, char *Colname); @@ -13,7 +13,7 @@ int main (int argc, char **argv) { - off_t i, Nx, Ny, Nbytes, Nread, Row; + off_t Nx, Ny, Nbytes, Nread, Row; int N, Nextend, Column, ListExtname, Layout; - char *Extname, *Colname, *line, ttype[80]; + char *Extname, *Colname, ttype[80]; FTable table; Header header; @@ -104,45 +104,150 @@ /* print a row */ if (Row) { - line = print_table_row (&table.buffer[Nx*Row], table.header); - fprintf (stdout, "%s\n", line); - free (line); + if (!print_table_rows (&table, Row, 1)) { + fprintf (stderr, "failed to print row\n"); + } exit (0); } /* print complete table */ - for (i = 0; i < Ny; i++) { - line = print_table_row (&table.buffer[Nx*i], table.header); - fprintf (stdout, "%s\n", line); - free (line); + if (!print_table_rows (&table, 0, Ny)) { + fprintf (stderr, "failed to print table\n"); } exit (0); } -/* print an ASCII table to a row with single spaces separating value */ -char *print_table_row (char *row, Header *header) { +# define SWAP_BYTE(BYTE) \ + tmp = BYTE[0]; BYTE[0] = BYTE[1]; BYTE[1] = tmp; +# define SWAP_WORD(BYTE) \ + tmp = BYTE[0]; BYTE[0] = BYTE[3]; BYTE[3] = tmp; \ + tmp = BYTE[1]; BYTE[1] = BYTE[2]; BYTE[2] = tmp; +# define SWAP_DBLE(BYTE) \ + tmp = BYTE[0]; BYTE[0] = BYTE[7]; BYTE[7] = tmp; \ + tmp = BYTE[1]; BYTE[1] = BYTE[6]; BYTE[6] = tmp; \ + tmp = BYTE[2]; BYTE[2] = BYTE[5]; BYTE[5] = tmp; \ + tmp = BYTE[3]; BYTE[3] = BYTE[4]; BYTE[4] = tmp; + +/* print Nrows of the given table starting at row 'start' */ +int print_table_rows (FTable *table, int start, int Nrows) { - off_t Nx; - int i, j, Nfields, Nbytes, Nvals, Oout, Oin; - char field[16], type[16], format[16], *line; - - gfits_scan (header, "NAXIS1", OFF_T_FMT, 1, &Nx); - gfits_scan (header, "TFIELDS", "%d", 1, &Nfields); + off_t Nx, Ny; + int n, i, j, Nfields, *Nbyte, *Nvals, Oout, Oin, Nv, Nb, byte, status; + char field[16], **types, format[16], type[80], *line, *row, tmp; + double *Tzero, *Tscal; + + gfits_scan (table->header, "NAXIS1", OFF_T_FMT, 1, &Nx); + gfits_scan (table->header, "NAXIS2", OFF_T_FMT, 1, &Ny); + gfits_scan (table->header, "TFIELDS", "%d", 1, &Nfields); + + if (start < 0) return FALSE; + if (start >= Ny) return FALSE; + if (Nrows < 0) return FALSE; + if (start + Nrows > Ny) return FALSE; /* assume we have one space per byte column */ - ALLOCATE (line, char, 2*Nx + 1); - - Oin = Oout = 0; - for (i = 1; i <= Nfields; i++) { - sprintf (field, "TFORM%d", i); - gfits_scan (header, field, "%s", 1, format); /* get field format */ - gfits_table_format (format, type, &Nvals, &Nbytes); /* convert to c-style */ - memcpy (&line[Oout], &row[Oin], Nvals*Nbytes); - for (j = 0; j < Nvals*Nbytes; j++) if (line[Oout+j] == 0) line[Oout+j] = ' '; - line[Oout+Nvals*Nbytes] = ' '; - Oout += Nvals*Nbytes + 1; - Oin += Nvals*Nbytes; - } - - return (line); + ALLOCATE (line, char, MAX(2*Nx+1,512)); + + ALLOCATE (types, char *, Nfields); + ALLOCATE (Nvals, int, Nfields); + ALLOCATE (Nbyte, int, Nfields); + ALLOCATE (Tzero, double, Nfields); + ALLOCATE (Tscal, double, Nfields); + + // determine the layout of the columns + for (i = 0; i < Nfields; i++) { + sprintf (field, "TFORM%d", i+1); + gfits_scan (table->header, field, "%s", 1, format); /* get field format */ + + if (Binary) { + gfits_bintable_format (format, type, &Nv, &Nb); /* convert to c-style */ + + sprintf (field, "TZERO%d", i+1); + status = gfits_scan (table[0].header, field, "%lf", 1, &Tzero[i]); /* get field format */ + if (!status) Tzero[i] = 0.0; + + sprintf (field, "TSCAL%d", i+1); + status = gfits_scan (table[0].header, field, "%lf", 1, &Tscal[i]); /* get field format */ + if (!status) Tscal[i] = 1.0; + } else { + gfits_table_format (format, type, &Nv, &Nb); /* convert to c-style */ + } + + types[i] = strcreate (type); + Nvals[i] = Nv; + Nbyte[i] = Nb; + } + + for (n = start; n < start + Nrows; n++) { + + row = &table->buffer[Nx*n]; + + if (Binary) { + byte = 0; // counter for byte element of this row + for (i = 0; i < Nfields; i++) { + int found = FALSE; + if (!strcmp (types[i], "char")) { + memcpy (line, &row[byte], Nvals[i]*Nbyte[i]); + fprintf (stdout, "%s ", line); + found = TRUE; + } else { + for (j = 0; j < Nvals[i]; j++) { + memcpy (line, &row[byte + Nbyte[i]*j], Nbyte[i]); + if (!strcmp (types[i], "int")) { +# ifdef BYTE_SWAP + SWAP_WORD (line); +# endif + fprintf (stdout, "%d ", (int)(*(int *)line * Tscal[i] + Tzero[i])); + found = TRUE; + } + if (!strcmp (types[i], "short")) { +# ifdef BYTE_SWAP + SWAP_BYTE (line); +# endif + fprintf (stdout, "%d ", (int)(*(short *)line * Tscal[i] + Tzero[i])); + found = TRUE; + } + if (!strcmp (types[i], "int64_t")) { +# ifdef BYTE_SWAP + SWAP_DBLE (line); +# endif + fprintf (stdout, "%" PRId64" ", (int64_t)(*(int64_t*)line * Tscal[i] + Tzero[i])); + found = TRUE; + } + if (!strcmp (types[i], "float")) { +# ifdef BYTE_SWAP + SWAP_WORD (line); +# endif + fprintf (stdout, "%e ", (*(float *)line * Tscal[i] + Tzero[i])); + found = TRUE; + } + if (!strcmp (types[i], "double")) { +# ifdef BYTE_SWAP + SWAP_DBLE (line); +# endif + fprintf (stdout, "%e ", (*(double *)line * Tscal[i] + Tzero[i])); + found = TRUE; + } + } + } + byte += Nvals[i]*Nbyte[i]; + if (!found) { + fprintf (stderr, "failed to find format for %d : %s\n", i, types[i]); + } + } + } else { + Oout = 0; + Oin = 0; + for (i = 0; i < Nfields; i++) { + memcpy (&line[Oout], &row[Oin], Nvals[i]*Nbyte[i]); + for (j = 0; j < Nvals[i]*Nbyte[i]; j++) if (line[Oout+j] == 0) line[Oout+j] = ' '; + line[Oout+Nvals[i]*Nbyte[i]] = ' '; + Oout += Nvals[i]*Nbyte[i] + 1; + Oin += Nvals[i]*Nbyte[i]; + } + fprintf (stdout, "%s ", line); + } + fprintf (stdout, "\n"); + } + return (TRUE); } @@ -325,8 +430,9 @@ sprintf (field, "TFORM%d", i); gfits_scan (header, field, "%s", 1, format); - if (Binary) + if (Binary) { + gfits_bintable_format (format, type, &Nv, &Nb); + } else { gfits_table_format (format, type, &Nv, &Nb); - else - gfits_bintable_format (format, type, &Nv, &Nb); + } Nstart += Nv*Nb; } @@ -334,9 +440,9 @@ sprintf (field, "TFORM%d", Column); gfits_scan (header, field, "%s", 1, format); - if (Binary) + if (Binary) { gfits_bintable_format (format, type, &Nv, &Nb); /* convert to c-style */ - else + } else { gfits_table_format (format, type, &Nv, &Nb); /* convert to c-style */ - + } ALLOCATE (line, char, Nv*Nb + 1); Index: trunk/Ohana/src/uniphot/Makefile =================================================================== --- trunk/Ohana/src/uniphot/Makefile (revision 31027) +++ trunk/Ohana/src/uniphot/Makefile (revision 31160) @@ -1,5 +1,5 @@ -default: uniphot setphot +default: uniphot setphot setfwhm help: -@echo "make options: uniphot setphot default help install default (uniphot setphot)" +@echo "make options: uniphot setphot setfwhm default help install default (uniphot setphot setfwhm)" include ../../Makefile.System @@ -19,5 +19,6 @@ uniphot: $(BIN)/uniphot.$(ARCH) setphot: $(BIN)/setphot.$(ARCH) -install: $(DESTBIN)/uniphot $(DESTBIN)/setphot +setfwhm: $(BIN)/setfwhm.$(ARCH) +install: $(DESTBIN)/uniphot $(DESTBIN)/setphot $(DESTBIN)/setfwhm UNIPHOT = \ @@ -34,4 +35,5 @@ $(SRC)/update.$(ARCH).o \ $(SRC)/update_catalog.$(ARCH).o \ +$(SRC)/convert.$(ARCH).o \ $(SRC)/SetSignals.$(ARCH).o \ $(SRC)/Shutdown.$(ARCH).o \ @@ -56,2 +58,17 @@ $(SETPHOT): $(INC)/uniphot.h $(BIN)/setphot.$(ARCH): $(SETPHOT) + +SETFWHM = \ +$(SRC)/setfwhm.$(ARCH).o \ +$(SRC)/initialize.$(ARCH).o \ +$(SRC)/ConfigInit.$(ARCH).o \ +$(SRC)/args.$(ARCH).o \ +$(SRC)/liststats.$(ARCH).o \ +$(SRC)/load_fwhm_table.$(ARCH).o \ +$(SRC)/load_images.$(ARCH).o \ +$(SRC)/match_fwhm_to_images.$(ARCH).o \ +$(SRC)/SetSignals.$(ARCH).o \ +$(SRC)/Shutdown.$(ARCH).o + +$(SETFWHM): $(INC)/uniphot.h +$(BIN)/setfwhm.$(ARCH): $(SETFWHM) Index: trunk/Ohana/src/uniphot/include/uniphot.h =================================================================== --- trunk/Ohana/src/uniphot/include/uniphot.h (revision 31027) +++ trunk/Ohana/src/uniphot/include/uniphot.h (revision 31160) @@ -29,4 +29,6 @@ typedef struct { + char tstart[64]; + char tstop[64]; char label[64]; float M; @@ -45,4 +47,12 @@ int found; } ZptTable; + +typedef struct { + float fwhm_major; + float fwhm_minor; + e_time time; + int found; + unsigned short photcode; +} FWHMTable; /* global variables set in parameter file */ @@ -124,2 +134,29 @@ int update_setphot PROTO((Image *image, off_t Nimage)); void update_catalog_setphot PROTO((Catalog *catalog, Image *image, off_t *index, off_t Nimage)); + +/*** time/coord conversion functions not supplied by libohana ***/ +time_t TimeRef PROTO((double time, time_t TimeReference, int TimeFormat)); +double TimeValue PROTO((time_t time, time_t TimeReference, int TimeFormat)); + +int hh_hms PROTO((double hh, int *hr, int *mn, double *sc)); +int dd_dms PROTO((double dd, int *dg, int *mn, double *sc)); +int hms_format PROTO((char *line, double value)); +int dms_format PROTO((char *line, double value)); +int hh_hm PROTO((double hh, int *hr, double *mn)); +int day_to_sec PROTO((char *string, time_t *second)); +int hms_to_sec PROTO((char *string, time_t *second)); +char *ohana_sec_to_hms PROTO((time_t second)); +char *ohana_sec_to_day PROTO((time_t second)); + +char *meade_deg_to_str PROTO((double deg)); +char *meade_ra_to_str PROTO((double deg)); +char *meade_dec_to_str PROTO((double deg)); +char *strptime PROTO((const char *s, const char *format, struct tm *tm)); +time_t GetTimeReference PROTO((char *reference)); +int GetTimeUnits PROTO((char *name)); + +void initialize_setfwhm PROTO((int argc, char **argv)); +int args_setfwhm PROTO((int argc, char **argv)); +FWHMTable *load_fwhm_table PROTO((char *filename, int *nfwhm)); +int match_fwhm_to_images PROTO((Image *image, off_t Nimage, FWHMTable *fwhm, int Nfwhm)); + Index: trunk/Ohana/src/uniphot/src/args.c =================================================================== --- trunk/Ohana/src/uniphot/src/args.c (revision 31027) +++ trunk/Ohana/src/uniphot/src/args.c (revision 31160) @@ -96,2 +96,26 @@ } +int args_setfwhm (int argc, char **argv) { + + int N; + + VERBOSE = FALSE; + if ((N = get_argument (argc, argv, "-v"))) { + VERBOSE = TRUE; + remove_argument (N, &argc, argv); + } + + UPDATE = FALSE; + if ((N = get_argument (argc, argv, "-update"))) { + remove_argument (N, &argc, argv); + UPDATE = TRUE; + } + + if (argc != 2) { + fprintf (stderr, "ERROR: USAGE: setfwhm (fwhmfile) [options]\n"); + exit (2); + } + + return (TRUE); +} + Index: trunk/Ohana/src/uniphot/src/convert.c =================================================================== --- trunk/Ohana/src/uniphot/src/convert.c (revision 31160) +++ trunk/Ohana/src/uniphot/src/convert.c (revision 31160) @@ -0,0 +1,337 @@ +# include "uniphot.h" +# define _XOPEN_SOURCE /* glibc2 (strptime) needs this */ +# include + +/** additional time / coordinate conversions not supplied by libohana **/ + +time_t GetTimeReference (char *reference) { + + time_t TimeReference; + + if (!ohana_str_to_time (reference, &TimeReference)) { + fprintf (stderr, "error in time reference %s\n", reference); + exit (2); + } + + return (TimeReference); +} + +int GetTimeUnits (char *name) { + + int Units; + + Units = FALSE; + if (!strcasecmp (name, "JD")) Units = TIME_JD; + if (!strcasecmp (name, "MJD")) Units = TIME_MJD; + if (!strcasecmp (name, "date")) Units = TIME_DATE; + if (!strcasecmp (name, "days")) Units = TIME_DAYS; + if (!strcasecmp (name, "hours")) Units = TIME_HOURS; + if (!strcasecmp (name, "min")) Units = TIME_MINUTES; + if (!strcasecmp (name, "sec")) Units = TIME_SECONDS; + if (!Units) { + fprintf (stderr, "error in time units %s\n", name); + exit (2); + } + + return (Units); +} + +int hh_hms (double hh, int *hr, int *mn, double *sc) { + + int N, flag; + + flag = SIGN(hh); + hh = fabs(hh); + + // rationalize hh to range -24.0 < hh < 24.0 + if (hh >= 24.0) { + N = (int)(hh/24.0); + hh -= 24.0*N; + } + + *hr = (int) hh; + *mn = (int) 60*(hh - *hr); + *sc = 3600.0*(hh - *hr - *mn / 60.0); + if (*sc > 59.99) { + *sc = 0.0; + *mn += 1.0; + } + *hr *= flag; + return (TRUE); +} + +int dd_dms (double dd, int *dg, int *mn, double *sc) { + + int flag; + + flag = SIGN(dd); + dd = fabs (dd); + *dg = (int) dd; + *mn = (int) 60*(dd - *dg); + *sc = 3600.0*(dd - *dg - *mn/60.0); + if (*sc > 59.99) { + *sc = 0; + *mn += 1.0; + } + *dg *= flag; + return (TRUE); +} + +int hms_format (char *line, double value) { + + int hr, mn; + double sc; + + hh_hms (value, &hr, &mn, &sc); + hr = (int) value; + if (isnan (value)) + sprintf (line, "xx:xx:xx.xx"); + else { + if (value < 0) { + sprintf (line, "-%02d:%02d:%05.2f", abs(hr), mn, sc); + } else { + sprintf (line, "+%02d:%02d:%05.2f", hr, mn, sc); + } + } + return (TRUE); +} + +int dms_format (char *line, double value) { + + int dg, mn; + double sc; + + dd_dms (value, &dg, &mn, &sc); + if (value < 0) { + sprintf (line, "-%02d:%02d:%05.2f", abs(dg), mn, sc); + } else { + sprintf (line, "+%02d:%02d:%05.2f", dg, mn, sc); + } + return (TRUE); +} + +/***** convert 00:00:00 or 00:00 to 0 - 86400 ****/ +int hms_to_sec (char *string, time_t *second) { + + char *p; + struct tm time; + + p = strptime (string, "%H:%M:%S", &time); + if (p != NULL) goto valid; + + p = strptime (string, "%H:%M", &time); + if (p != NULL) goto valid; + + return (FALSE); + +valid: + if (*p) return (FALSE); + *second = time.tm_hour*3600 + time.tm_min*60 + time.tm_sec; + return (TRUE); +} + +/***** convert Mon[@00:00:00] or 00:00 to 0 - 86400*7 ****/ +int day_to_sec (char *string, time_t *second) { + + char *p; + struct tm time; + + bzero (&time, sizeof(time)); + p = strptime (string, "%A@%H:%M:%S", &time); + if (p != NULL) goto valid; + + p = strptime (string, "%A@%H:%M", &time); + if (p != NULL) goto valid; + + p = strptime (string, "%A@%H", &time); + if (p != NULL) goto valid; + + p = strptime (string, "%A", &time); + if (p != NULL) goto valid; + + return (FALSE); + +valid: + if (*p) return (FALSE); + *second = time.tm_wday*86400 + time.tm_hour*3600 + time.tm_min*60 + time.tm_sec; + return (TRUE); +} + +/***** convert seconds to HH:MM:SS ****/ +char *ohana_sec_to_hms (time_t second) { + + struct tm *gmt; + char *line; + + ALLOCATE (line, char, 64); + gmt = gmtime (&second); + sprintf (line, "%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); + return (line); +} + +/***** convert seconds to Day@HH:MM:SS ****/ +char *ohana_sec_to_day (time_t second) { + + struct tm *gmt; + char *line; + + ALLOCATE (line, char, 64); + gmt = gmtime (&second); + switch (gmt[0].tm_wday) { + case 0: + sprintf (line, "Sun@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); + break; + case 1: + sprintf (line, "Mon@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); + break; + case 2: + sprintf (line, "Tue@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); + break; + case 3: + sprintf (line, "Wed@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); + break; + case 4: + sprintf (line, "Thu@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); + break; + case 5: + sprintf (line, "Fri@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); + break; + case 6: + sprintf (line, "Sat@%02d:%02d:%02d", gmt[0].tm_hour, gmt[0].tm_min, gmt[0].tm_sec); + break; + } + return (line); +} + +int hh_hm (double hh, int *hr, double *mn) { + + int flag; + + flag = SIGN(hh); + hh = fabs (hh); + + *mn = 60.0*(hh - (int)hh); + *hr = (int) hh; + *hr *= flag; + return (TRUE); +} + +char *meade_deg_to_str (double deg) { + + int hr; + double mn; + char *line; + + ALLOCATE (line, char, 16); + + hh_hm (deg, &hr, &mn); + + sprintf (line, "%03d:%04.1f", abs(hr), mn); + return (line); +} + +char *meade_ra_to_str (double deg) { + + int hr; + double mn; + char *line; + + ALLOCATE (line, char, 16); + + hh_hm (deg/15.0, &hr, &mn); + + sprintf (line, "%02d:%04.1f", abs(hr), mn); + return (line); +} + +char *meade_dec_to_str (double deg) { + + int hr; + double mn; + char *line; + + ALLOCATE (line, char, 16); + + hh_hm (deg, &hr, &mn); + + if (deg < 0) { + sprintf (line, "-%02d:%04.1f", abs(hr), mn); + } else { + sprintf (line, "+%02d:%04.1f", hr, mn); + } + return (line); +} + +/* convert UNIX time to a value referenced to the TimeReference in the given unit */ +double TimeValue (time_t time, time_t TimeReference, int TimeFormat) { + + double value, dt; + + dt = (time > TimeReference) ? (time - TimeReference) : -1 * (double)(TimeReference - time); + switch (TimeFormat) { + case TIME_JD: + value = time / 86400.0 + 2440587.5; + break; + case TIME_MJD: + value = time / 86400.0 + 40587.0; + break; + case TIME_DAYS: + value = dt / 86400.0; + break; + case TIME_HOURS: + value = dt / 3600.0; + break; + case TIME_MINUTES: + value = dt / 60.0; + break; + case TIME_SECONDS: + default: + value = dt; + break; + } + return (value); +} + +/* convert time value referenced to the TimeReference in the given unit to UNIX time */ +time_t TimeRef (double value, time_t TimeReference, int TimeFormat) { + + int dt; + time_t time; + + switch (TimeFormat) { + case TIME_JD: + time = (value - 2440587.5) * 86400.0; + return (time); + break; + case TIME_MJD: + time = (value - 40587.0) * 86400.0; + return (time); + break; + case TIME_DAYS: + dt = value * 86400.0; + break; + case TIME_HOURS: + dt = value * 3600.0; + break; + case TIME_MINUTES: + dt = value * 60.0; + break; + case TIME_SECONDS: + default: + dt = value; + break; + } + + time = TimeReference + dt; + return (time); +} + +/* times may be in forms as: + * 20040200450s (N seconds since 1970.0) + * 2440900.232j (julian date) + * 99/02/23,03:22:18 (date string) + * (separators may be anything except space, +, -) + * 99:02:15:12:23:30 + * 99:02:15:12h23m30s + */ + Index: trunk/Ohana/src/uniphot/src/dumpresult.c =================================================================== --- trunk/Ohana/src/uniphot/src/dumpresult.c (revision 31027) +++ trunk/Ohana/src/uniphot/src/dumpresult.c (revision 31160) @@ -17,5 +17,5 @@ Mgrp = tgrp[0].M; fprintf (f, "%7.4f %7.4f %7.4f %7.4f %10.6f %10.6f %f %s\n", - 0.001*Mcal, 0.001*Mgrp, 0.001*Mset, 0.001*sgroup[i].image[j][0].dMcal, + Mcal, Mgrp, Mset, sgroup[i].image[j][0].dMcal, sgroup[i].image[j][0].coords.crval1, sgroup[i].image[j][0].coords.crval2, (sgroup[i].image[j][0].tzero-915148800)/86400.0, tgrp[0].label); } Index: trunk/Ohana/src/uniphot/src/find_image_sgroups.c =================================================================== --- trunk/Ohana/src/uniphot/src/find_image_sgroups.c (revision 31027) +++ trunk/Ohana/src/uniphot/src/find_image_sgroups.c (revision 31160) @@ -3,6 +3,6 @@ Group *find_image_sgroups (FITS_DB *db, ImageLink **Imlink, int *Nsgroup) { - off_t i, Nimage; - int j, Ngroup, Nentry, NENTRY; + off_t i, j, Nimage; + int Ngroup, Nentry, NENTRY; double r, d, x, y, radius; Group *group; @@ -27,4 +27,8 @@ ALLOCATE (group, Group, Nimage); + if (VERBOSE) fprintf (stderr, "finding images\n"); + BuildChipMatch (image, Nimage); + // MARKTIME("build chip match: %f sec\n", dtime); + /* set imlink.sgroups = NULL as a marker */ for (i = 0; i < Nimage; i++) imlink[i].sgroup = NULL; @@ -34,6 +38,22 @@ if (image[i].flags & ID_IMAGE_PHOTOM_NOCAL) continue; + // XXX optionally, we should be able to use ONLY the DIS or NOT the DIS images + // NOCAL above is used to mark images which do not match the photcode (including the DIS) + // if (!strcmp(&image[i].coords.ctype[4], "-DIS")) continue; + + // this adds 1.3 sec for 3M images + if (!FindMosaicForImage (image, Nimage, i)) { + fprintf (stderr, "cannot find mosaic for "OFF_T_FMT"\n", i); + continue; + } + + /* define image center - note the DIS images (mosaic phu) are special */ + if (!strcmp(&image[i].coords.ctype[4], "-DIS")) { + XY_to_RD (&r, &d, 0.0, 0.0, &image[i].coords); + } else { + XY_to_RD (&r, &d, 0.5*image[i].NX, 0.5*image[i].NX, &image[i].coords); + } + /* new sgroup, set ref coords */ - XY_to_RD (&r, &d, 0.5*image[i].NX, 0.5*image[i].NX, &image[i].coords); coords.crval1 = r; coords.crval2 = d; @@ -59,8 +79,22 @@ if (image[j].flags & ID_IMAGE_PHOTOM_NOCAL) continue; if (imlink[j].sgroup != NULL) continue; + // XXX optionally, we should be able to use ONLY the DIS or NOT the DIS images + // NOCAL above is used to mark images which do not match the photcode (including the DIS) + // if (!strcmp(&image[j].coords.ctype[4], "-DIS")) continue; + + // this adds 1.3 sec for 3M images + if (!FindMosaicForImage (image, Nimage, j)) { + fprintf (stderr, "cannot find mosaic for "OFF_T_FMT"\n", j); + continue; + } /* project image center to local coords, check radius */ - XY_to_RD (&r, &d, 0.5*image[j].NX, 0.5*image[j].NX, &image[j].coords); + if (!strcmp(&image[j].coords.ctype[4], "-DIS")) { + XY_to_RD (&r, &d, 0.0, 0.0, &image[j].coords); + } else { + XY_to_RD (&r, &d, 0.5*image[j].NX, 0.5*image[j].NX, &image[j].coords); + } if (!RD_to_XY (&x, &y, r, d, &coords)) continue; + /* RD_to_XY returns FALSE if opposite hemispheres */ radius = hypot (x, y); Index: trunk/Ohana/src/uniphot/src/find_image_tgroups.c =================================================================== --- trunk/Ohana/src/uniphot/src/find_image_tgroups.c (revision 31027) +++ trunk/Ohana/src/uniphot/src/find_image_tgroups.c (revision 31160) @@ -5,5 +5,5 @@ char *start, *stop; int j, Ngroup, NGROUP, Nentry, NENTRY; - off_t i, Nimage; + off_t i, Nimage, Ntime; unsigned int *time, *tmin, *tmax; Group *group; @@ -21,8 +21,12 @@ /* sort time list (use only valid images?) */ ALLOCATE (time, unsigned int, Nimage); + Ntime = 0; for (i = 0; i < Nimage; i++) { - time[i] = image[i].tzero; + if (image[i].flags & ID_IMAGE_PHOTOM_NOCAL) continue; + if (!strcmp(&image[i].coords.ctype[4], "-DIS")) continue; + time[Ntime] = image[i].tzero; + Ntime ++; } - sort_time (time, Nimage); + sort_time (time, Ntime); /* find groups with dt < TRANGE */ @@ -34,5 +38,5 @@ /* generate tgroups */ - for (i = 0; i < Nimage - 1; i++) { + for (i = 0; i < Ntime - 1; i++) { if (time[i+1] - time[i] < TRANGE) continue; @@ -47,5 +51,5 @@ tmin[Ngroup] = time[i + 1]; } - tmax[Ngroup] = time[Nimage - 1]; + tmax[Ngroup] = time[Ntime - 1]; Ngroup ++; ALLOCATE (group, Group, Ngroup); @@ -63,4 +67,6 @@ stop = ohana_sec_to_date (tmax[i]); snprintf (group[i].label, 64, "%s - %s", start, stop); + strcpy(group[i].tstart, start); + strcpy(group[i].tstop, stop); free (start); free (stop); @@ -70,4 +76,5 @@ if (image[j].tzero > tmax[i]) continue; if (image[j].flags & ID_IMAGE_PHOTOM_NOCAL) continue; + if (!strcmp(&image[j].coords.ctype[4], "-DIS")) continue; group[i].image[Nentry] = &image[j]; Index: trunk/Ohana/src/uniphot/src/fit_groups.c =================================================================== --- trunk/Ohana/src/uniphot/src/fit_groups.c (revision 31027) +++ trunk/Ohana/src/uniphot/src/fit_groups.c (revision 31160) @@ -29,4 +29,6 @@ tgroup[i].Ngood = stats.Nmeas; + // fprintf (stderr, "tgroup %d : %f +/- %f : %d stars\n", i, stats.mean, stats.sigma, stats.Nmeas); + initstats ("MEAN"); liststats (mlist, dlist, Nlist, &stats); @@ -67,4 +69,6 @@ sgroup[i].Ngood = stats.Nmeas; + // fprintf (stderr, "sgroup %d : %f +/- %f : %d stars\n", i, stats.mean, stats.sigma, stats.Nmeas); + initstats ("MEAN"); liststats (mlist, dlist, Nlist, &stats); Index: trunk/Ohana/src/uniphot/src/initialize.c =================================================================== --- trunk/Ohana/src/uniphot/src/initialize.c (revision 31027) +++ trunk/Ohana/src/uniphot/src/initialize.c (revision 31160) @@ -28,2 +28,9 @@ } +void initialize_setfwhm (int argc, char **argv) { + + /* are these set correctly? */ + ConfigInit (&argc, argv); + args_setfwhm (argc, argv); +} + Index: trunk/Ohana/src/uniphot/src/load_fwhm_table.c =================================================================== --- trunk/Ohana/src/uniphot/src/load_fwhm_table.c (revision 31160) +++ trunk/Ohana/src/uniphot/src/load_fwhm_table.c (revision 31160) @@ -0,0 +1,49 @@ +# include "uniphot.h" + +FWHMTable *load_fwhm_table (char *filename, int *nfwhm) { + + char name[1024], photcode[256]; + int Nfwhm, NFWHM; + FWHMTable *fwhm; + double fwhm_major, fwhm_minor, mjd; + + FILE *f; + + f = fopen (filename, "r"); + if (!f) { + fprintf (stderr, "ERROR: cannot open fwhm table file %s\n", filename); + exit (1); + } + + Nfwhm = 0; + NFWHM = 100; + ALLOCATE (fwhm, FWHMTable, NFWHM); + + // format is fixed: (time in mjd) (fwhm_major) (fwhm_minor) + int status; + while ((status = fscanf (f, "%s %lf %s %lf %lf", name, &mjd, photcode, &fwhm_major, &fwhm_minor)) == 5) { + fwhm[Nfwhm].fwhm_major = fwhm_major; + fwhm[Nfwhm].fwhm_minor = fwhm_minor; + fwhm[Nfwhm].time = ohana_mjd_to_sec (mjd); + fwhm[Nfwhm].photcode = GetPhotcodeCodebyName (photcode); + fwhm[Nfwhm].found = FALSE; + + if (fwhm[Nfwhm].photcode == 0) { + fprintf (stderr, "error in photcode %s\n", photcode); + abort(); + } + + Nfwhm ++; + CHECK_REALLOCATE (fwhm, FWHMTable, NFWHM, Nfwhm, 100); + } + + if (status != EOF) { + fprintf (stderr, "unexpected formatting on line %d (status = %d)\n", Nfwhm, status); + exit (2); + } + + fprintf (stderr, "loaded %d fwhm values\n", Nfwhm); + + *nfwhm = Nfwhm; + return fwhm; +} Index: trunk/Ohana/src/uniphot/src/match_fwhm_to_images.c =================================================================== --- trunk/Ohana/src/uniphot/src/match_fwhm_to_images.c (revision 31160) +++ trunk/Ohana/src/uniphot/src/match_fwhm_to_images.c (revision 31160) @@ -0,0 +1,175 @@ +# include "uniphot.h" + +/* sort a coordinate pair (X,Y) and the associated index (S) */ +static void sort_fwhm_by_time (FWHMTable *fwhm, int Nfwhm) { + +# define SWAPFUNC(A,B){ FWHMTable tmp; \ + tmp = fwhm[A]; fwhm[A] = fwhm[B]; fwhm[B] = tmp; \ + } +# define COMPARE(A,B)(fwhm[A].time < fwhm[B].time) + + OHANA_SORT (Nfwhm, COMPARE, SWAPFUNC); + +# undef SWAPFUNC +# undef COMPARE +} + +# define SCALE 25.0 + +int match_fwhm_to_images (Image *image, off_t Nimage, FWHMTable *fwhm, int Nfwhm) { + + // I have two lists. I need to match images->tzero to fwhm->time. multiple images may match a single fwhm + + // sort both lists by time (or at least get sorted indices) + // sweep through both lists, advancing the one that is lagging + // apply the matches + + int dT, NImatch, Nmatch; + off_t i, Ni, Nf, *index; + float *fwhmMajorSet, *fwhmMinorSet; + + // create index and sort + ALLOCATE (index, off_t, Nimage); + for (i = 0; i < Nimage; i++) { + index[i] = i; + } + sort_image_subset(image, index, Nimage); // slightly misnamed sort function from libdvo + + // XXX GPC1 hack: accumulate the values within a single exposure (time) to provide the + // mosaic exposure with a mean or median + int Ncode = 0; + int NCODE = 100; + ALLOCATE (fwhmMajorSet, float, NCODE); + ALLOCATE (fwhmMinorSet, float, NCODE); + + // sort the fwhm + sort_fwhm_by_time (fwhm, Nfwhm); + + NImatch = 0; // matched images + for (i = Nf = 0; (i < Nimage) && (Nf < Nfwhm); ) { + + if (i % 1000 == 0) fprintf (stderr, "."); + if (Nf % 100 == 0) fprintf (stderr, "!"); + + Ni = index[i]; + dT = image[Ni].tzero - fwhm[Nf].time; + + // negative dT, i is too small (allow a 1sec overlap window) + if (dT < -1) { + i++; + continue; + } + + // XXX careful about definition of image->tzero and fwhm->time + // negative dT, i is too small (allow a 1sec overlap window) + if (dT > +1) { + Nf++; + continue; + } + + { + // loop over the exposure in this time block to get a median fwhm + float fwhmMajorMedian = 0; + float fwhmMinorMedian = 0; + + int NfS; + int dTf = 0; + Ncode = 0; + for (NfS = Nf; (dTf < +1) && (NfS < Nfwhm); NfS++) { + dTf = image[Ni].tzero - fwhm[NfS].time; + if (dTf < -1 ) { + break; // too far in fwhm sequence + } + + fwhmMajorSet[Ncode] = fwhm[NfS].fwhm_major; + fwhmMinorSet[Ncode] = fwhm[NfS].fwhm_minor; + Ncode ++; + if (Ncode == NCODE) { + NCODE += 100; + REALLOCATE(fwhmMajorSet, float, NCODE); + REALLOCATE(fwhmMinorSet, float, NCODE); + } + } + + // find the medians + fsort (fwhmMajorSet, Ncode); + fsort (fwhmMinorSet, Ncode); + if (Ncode) { + if (Ncode % 2) { + fwhmMajorMedian = fwhmMajorSet[(int)(Ncode/2)]; + fwhmMinorMedian = fwhmMinorSet[(int)(Ncode/2)]; + } else { + fwhmMajorMedian = 0.5*(fwhmMajorSet[(int)(Ncode/2)-1] + fwhmMajorSet[(int)(Ncode/2)]); + fwhmMinorMedian = 0.5*(fwhmMinorSet[(int)(Ncode/2)-1] + fwhmMinorSet[(int)(Ncode/2)]); + } + } + + // now find the mosaic image and set the median values + int iS; + int dTi = 0; + for (iS = i; (dTi < +1) && (iS < Nimage); iS++) { + int NiS = index[iS]; + dTi = image[NiS].tzero - fwhm[Nf].time; + if (dTi > 1) { + break; // too far in image sequence + } + // XXX possible hack: if the mosaics ever get a photcode, this will break + if (image[NiS].photcode) continue; + image[NiS].fwhm_x = fwhmMajorMedian * SCALE; + image[NiS].fwhm_y = fwhmMinorMedian * SCALE; + break; + } + } + + // we have two sets: images[Ni,Ni+N] and fwhm[Nf,Nf+M], now we need to match them by photcode + int iNext = i + 1; + int NfNext = Nf + 1; + int dTi = 0; + int iS; + for (iS = i; (dTi < +1) && (iS < Nimage); iS++) { + int NiS = index[iS]; + dTi = image[NiS].tzero - fwhm[Nf].time; + if (dTi > 1) { + iNext = iS; + break; // too far in image sequence + } + + int dTf = 0; + int NfS; + for (NfS = Nf; (dTf < +1) && (NfS < Nfwhm); NfS++) { + dTf = image[Ni].tzero - fwhm[NfS].time; + if (dTf < -1 ) { + NfNext = NfS; + break; // too far in fwhm sequence + } + + if (image[NiS].photcode == fwhm[NfS].photcode) { + // we have a match: set zpt and record the match + image[NiS].fwhm_x = fwhm[NfS].fwhm_major * SCALE; + image[NiS].fwhm_y = fwhm[NfS].fwhm_minor * SCALE; + fwhm[NfS].found = TRUE; + NImatch ++; + } + } + } + + // advance the image counter only -- a single zpt may match more than one image + Nf = NfNext; + i = iNext; + } + + // how many fwhm have we matched? + Nmatch = 0; + for (Nf = 0; Nf < Nfwhm; Nf++) { + if (fwhm[Nf].found) { + Nmatch ++; + } else { + // fprintf (stderr, "%d %d\n", fwhm[Nf].time, fwhm[Nf].photcode); + } + } + + fprintf (stderr, "found %d zpt matches, %d image matches\n", Nmatch, NImatch); + + return (TRUE); +} + Index: trunk/Ohana/src/uniphot/src/setfwhm.c =================================================================== --- trunk/Ohana/src/uniphot/src/setfwhm.c (revision 31160) +++ trunk/Ohana/src/uniphot/src/setfwhm.c (revision 31160) @@ -0,0 +1,48 @@ +# include "uniphot.h" + +int main (int argc, char **argv) { + + off_t Nimage; + int status, Nfwhm; + FITS_DB db; + FWHMTable *fwhm; + Image *image; + + /* get configuration info, args, lockfile */ + initialize_setfwhm (argc, argv); + + set_db (&db); + status = dvo_image_lock (&db, ImageCat, 60.0, (UPDATE ? LCK_XCLD : LCK_SOFT)); + if (!status) Shutdown ("ERROR: failure to lock image catalog %s", db.filename); + if (db.dbstate == LCK_EMPTY) Shutdown ("ERROR: No images in catalog %s (1)", db.filename); + + fwhm = load_fwhm_table (argv[1], &Nfwhm); + + // load images + image = load_images_setphot (&db, &Nimage); + if (!UPDATE) dvo_image_unlock (&db); + + // apply the newly loaded fwhm values to the corresponding images + match_fwhm_to_images (image, Nimage, fwhm, Nfwhm); + + // write image table + if (UPDATE) { + dvo_image_save (&db, VERBOSE); + } + + exit (0); +} + + +/* setphot : set the zero points for images in the db (perhaps based on external information) + setphot (zpt_table) + + * load text table of zpts, time (photcode?) + * load images + * match images to zpts + * set zpts + * load catalogs + * update detection (& averages?) + + */ + Index: trunk/Ohana/src/uniphot/src/subset_images.c =================================================================== --- trunk/Ohana/src/uniphot/src/subset_images.c (revision 31027) +++ trunk/Ohana/src/uniphot/src/subset_images.c (revision 31160) @@ -3,9 +3,9 @@ int subset_images (FITS_DB *db) { - off_t i, Nimage, Nkeep, *keep; + off_t i, Nimage; int equiv; Image *image; - /* use a vtable to keep the images to be calibrated */ + // convert from the binary I/O format to the internal structure (Image) image = gfits_table_get_Image (&db[0].ftable, &Nimage, &db[0].swapped); if (!image) { @@ -13,7 +13,4 @@ exit (2); } - - Nkeep = 0; - ALLOCATE (keep, off_t, Nimage); /* mark images to be calibrated */ @@ -32,9 +29,5 @@ } image[i].flags &= ~ID_IMAGE_PHOTOM_NOCAL; - keep[Nkeep] = i; - Nkeep ++; } - - gfits_vtable_from_ftable (&db[0].ftable, &db[0].vtable, keep, Nkeep); return (TRUE); } Index: trunk/Ohana/src/uniphot/src/uniphot.c =================================================================== --- trunk/Ohana/src/uniphot/src/uniphot.c (revision 31027) +++ trunk/Ohana/src/uniphot/src/uniphot.c (revision 31160) @@ -1,3 +1,6 @@ # include "uniphot.h" + +static char *timeref = "2000/01/01,00:00:00"; +static char *timeunit = "days"; int main (int argc, char **argv) { @@ -8,4 +11,8 @@ FITS_DB db; + // set up time format stuff + time_t TimeReference = GetTimeReference (timeref); + int TimeUnits = GetTimeUnits (timeunit); + /* get configuration info, args, lockfile */ initialize_uniphot (argc, argv); @@ -15,8 +22,8 @@ if (!status) Shutdown ("ERROR: failure to lock image catalog %s", db.filename); if (db.dbstate == LCK_EMPTY) Shutdown ("ERROR: No images in catalog %s (1)", db.filename); - if (!UPDATE) dvo_image_unlock (&db); /* load images */ load_images_uniphot (&db); + if (!UPDATE) dvo_image_unlock (&db); /* filter image list by selection */ @@ -38,7 +45,23 @@ fprintf (stdout, "# NLOOP: %d\n", NLOOP); fprintf (stdout, "# time groups : %d\n", Ntgroup); + fprintf (stdout, "# TIMEREF : %s\n", timeref); + fprintf (stdout, "# TIMEFORMAT : %s\n", timeunit); for (i = 0; i < Ntgroup; i++) { - fprintf (stdout, "%s %5d %5d %7.4f %7.4f %7.4f\n", tgroup[i].label, - tgroup[i].Nimage, tgroup[i].Ngood, 0.001*tgroup[i].M, 0.001*tgroup[i].dM, 0.001*tgroup[i].dMsub); + + double tstart = NAN; + double tstop = NAN; + + time_t time; + if (ohana_str_to_time (tgroup[i].tstart, &time)) { + tstart = TimeValue (time, TimeReference, TimeUnits); + } + if (ohana_str_to_time (tgroup[i].tstop, &time)) { + tstop = TimeValue (time, TimeReference, TimeUnits); + } + + fprintf (stdout, "%s : %12.6f %12.6f : %5d %5d %7.4f %7.4f %7.4f\n", + tgroup[i].label, tstart, tstop, + tgroup[i].Nimage, tgroup[i].Ngood, + tgroup[i].M, tgroup[i].dM, tgroup[i].dMsub); } fprintf (stdout, "\n"); @@ -47,5 +70,5 @@ for (i = 0; i < Nsgroup; i++) { fprintf (stdout, "%s %5d %5d %7.4f %7.4f %7.4f\n", sgroup[i].label, - sgroup[i].Nimage, sgroup[i].Ngood, 0.001*sgroup[i].M, 0.001*sgroup[i].dM, 0.001*sgroup[i].dMsub); + sgroup[i].Nimage, sgroup[i].Ngood, sgroup[i].M, sgroup[i].dM, sgroup[i].dMsub); } if (!UPDATE) exit (0); Index: trunk/Ohana/src/uniphot/src/update.c =================================================================== --- trunk/Ohana/src/uniphot/src/update.c (revision 31027) +++ trunk/Ohana/src/uniphot/src/update.c (revision 31160) @@ -4,5 +4,5 @@ void update (FITS_DB *db, Group *sgroup, int Nsgroup) { - off_t i, Nimage; + off_t i, Nimage, Nkeep, *keep; int j, status, Nmin; char line[256]; @@ -19,7 +19,16 @@ } - /* clear the NOCAL flags */ + // create a subset list so we can make a vtable + Nkeep = 0; + ALLOCATE (keep, off_t, Nimage); + + // identify the images used and clear the NOCAL flags on the rest for (i = 0; i < Nimage; i++) { - image[i].flags &= ~ID_IMAGE_PHOTOM_NOCAL; + if (image[i].flags & ID_IMAGE_PHOTOM_NOCAL) { + image[i].flags &= ~ID_IMAGE_PHOTOM_NOCAL; + continue; + } + keep[Nkeep] = i; + Nkeep ++; } @@ -31,6 +40,12 @@ } - /** write image table **/ + // save the rows in the image table which were used in this analysis + gfits_vtable_from_ftable (&db[0].ftable, &db[0].vtable, keep, Nkeep); + + // write image table dvo_image_update (db, VERBOSE); + + // XXX need to fix the update for the catalog (or make it optional) + return; // XXX this process uses the existence of the file to perform the update