Index: /trunk/Ohana/src/addstar/include/addstar.h =================================================================== --- /trunk/Ohana/src/addstar/include/addstar.h (revision 37035) +++ /trunk/Ohana/src/addstar/include/addstar.h (revision 37036) @@ -41,6 +41,7 @@ typedef struct { - Average average; - Measure measure; + Average average; + Measure measure; + Lensing *lensing; // optionally carry out the lensing measurements int found; } Stars; @@ -250,10 +251,14 @@ void update_coords PROTO((Average *average, Measure *measure, off_t *next)); off_t *init_measure_links PROTO((Average *average, off_t Naverage, Measure *measure, off_t Nmeasure)); +off_t *init_lensing_links PROTO((Average *average, off_t Naverage, Lensing *lensing, off_t Nlensing)); off_t *init_missing_links PROTO((Average *average, off_t Naverage, Missing *missing, off_t Nmissing)); off_t add_meas_link PROTO((Average *average, off_t *next_meas, off_t Nmeasure, off_t NMEASURE)); off_t add_miss_link PROTO((Average *average, off_t *next_miss, off_t Nmissing)); +off_t add_lens_link PROTO((Average *average, off_t *next_lens, off_t Nlensing, off_t NLENSING)); off_t *build_measure_links PROTO((Average *average, off_t Naverage, Measure *measure, off_t Nmeasure)); +off_t *build_lensing_links PROTO((Average *average, off_t Naverage, Lensing *lensing, off_t Nlensing)); Measure *sort_measure PROTO((Average *average, off_t Naverage, Measure *measure, off_t Nmeasure, off_t *next_meas)); Missing *sort_missing PROTO((Average *average, off_t Naverage, Missing *missing, off_t Nmissing, off_t *next_miss)); +Lensing *sort_lensing PROTO((Average *average, off_t Naverage, Lensing *lensing, off_t Nlensing, off_t *next_lens)); int ImageOptions PROTO((AddstarClientOptions *options, Image *images, off_t Nimages)); int GetFileMode PROTO((Header *header)); Index: /trunk/Ohana/src/addstar/src/FilterStars.c =================================================================== --- /trunk/Ohana/src/addstar/src/FilterStars.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/FilterStars.c (revision 37036) @@ -99,4 +99,10 @@ stars[N].measure.dFluxKron /= image[0].exptime; } + if (!isnan(stars[N].measure.FluxAp)) { + stars[N].measure.FluxAp /= image[0].exptime; + } + if (!isnan(stars[N].measure.dFluxAp)) { + stars[N].measure.dFluxAp /= image[0].exptime; + } // the external ID is supplied, but do we trust it? Index: /trunk/Ohana/src/addstar/src/MatchHeaders.c =================================================================== --- /trunk/Ohana/src/addstar/src/MatchHeaders.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/MatchHeaders.c (revision 37036) @@ -59,7 +59,10 @@ if (!strcmp (exttype, "PS1_V3")) goto keep; if (!strcmp (exttype, "PS1_V4")) goto keep; + if (!strcmp (exttype, "PS1_V5")) goto keep; if (!strcmp (exttype, "PS1_SV1")) goto keep; if (!strcmp (exttype, "PS1_SV2")) goto keep; - if (!strcmp (exttype, "PS1_DV3")) { + if (!strcmp (exttype, "PS1_SV3")) goto keep; + if (!strcmp (exttype, "PS1_DV3")) goto keep; + if (!strcmp (exttype, "PS1_DV4")) { goto keep; } Index: /trunk/Ohana/src/addstar/src/ReadImageHeader.c =================================================================== --- /trunk/Ohana/src/addstar/src/ReadImageHeader.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/ReadImageHeader.c (revision 37036) @@ -45,5 +45,6 @@ image[0].NY = Ny; - image[0].refColor = NAN; + image[0].refColorBlue = NAN; + image[0].refColorRed = NAN; if (!gfits_scan (header, "TZERO", "%d", 1, &image[0].tzero) && !ACCEPT_TIME) { @@ -135,5 +136,6 @@ tmp = 0; - gfits_scan (header, "FLIMIT", "%lf", 1, &tmp); + // gfits_scan (header, "FLIMIT", "%lf", 1, &tmp); + gfits_scan (header, "DETEFF.MAGREF", "%lf", 1, &tmp); image[0].detection_limit = tmp * 10.0; Index: /trunk/Ohana/src/addstar/src/ReadStarsFITS.c =================================================================== --- /trunk/Ohana/src/addstar/src/ReadStarsFITS.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/ReadStarsFITS.c (revision 37036) @@ -12,8 +12,12 @@ Stars *Convert_PS1_V3 PROTO((FTable *table, unsigned int *nstars)); Stars *Convert_PS1_V4 PROTO((FTable *table, unsigned int *nstars)); +Stars *Convert_PS1_V5 PROTO((FTable *table, unsigned int *nstars)); +Stars *Convert_PS1_V5_Lensing 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)); Stars *Convert_PS1_SV2 PROTO((FTable *table, unsigned int *nstars)); +Stars *Convert_PS1_SV3 PROTO((FTable *table, unsigned int *nstars)); Stars *Convert_PS1_DV3 PROTO((FTable *table, unsigned int *nstars)); +Stars *Convert_PS1_DV4 PROTO((FTable *table, unsigned int *nstars)); // given a file with the pointer at the start of the table block and the @@ -73,4 +77,11 @@ stars = Convert_PS1_V4 (&table, &Nstars); } + if (!strcmp (type, "PS1_V5")) { + if (table.header[0].Naxis[0] == 312) { + stars = Convert_PS1_V5_Lensing (&table, &Nstars); + } else { + stars = Convert_PS1_V5 (&table, &Nstars); + } + } if (!strcmp (type, "PS1_SV1")) { stars = Convert_PS1_SV1 (&table, &Nstars); @@ -79,7 +90,13 @@ stars = Convert_PS1_SV2 (&table, &Nstars); } + if (!strcmp (type, "PS1_SV3")) { + stars = Convert_PS1_SV3 (&table, &Nstars); + } if (!strcmp (type, "PS1_DV3")) { stars = Convert_PS1_DV3 (&table, &Nstars); } + if (!strcmp (type, "PS1_DV4")) { + stars = Convert_PS1_DV4 (&table, &Nstars); + } if (stars == NULL) { fprintf (stderr, "invalid table type %s\n", type); @@ -91,4 +108,21 @@ return stars; +} + +float GetFluxFromFluxOrMag (float flux, float mag) { + + if (isnan(mag) && isnan(flux)) return NAN; + + if (isnan(flux)) return pow(10.0,-0.4*mag); + return flux; +} + +float GetFluxErrFromFluxOrMag (float dFlux, float flux, float dMag) { + + if (isnan(dMag) && isnan(dFlux)) return NAN; + if (isnan(dFlux) && isnan(flux)) return NAN; + + if (isnan(dFlux)) return (fabs(dMag * flux)); + return dFlux; } @@ -360,4 +394,5 @@ stars[i].measure.dMcal = ps1data[i].dMcal; stars[i].measure.Map = ps1data[i].Map + ZeroPt; + stars[i].measure.dMap = ps1data[i].dM; // a proxy measure stars[i].measure.Mkron = NAN; // not provided by PS1_V1: @@ -445,4 +480,5 @@ stars[i].measure.dMcal = ps1data[i].dMcal; stars[i].measure.Map = ps1data[i].Map + ZeroPt; + stars[i].measure.dMap = ps1data[i].dM; // a proxy measure stars[i].measure.Mkron = NAN; // not provided by PS1_V1_Alt: @@ -522,4 +558,5 @@ stars[i].measure.dMcal = ps1data[i].dMcal; stars[i].measure.Map = ps1data[i].Map + ZeroPt; + stars[i].measure.dMap = ps1data[i].dM; // a proxy measure stars[i].measure.Mkron = NAN; // not provided by PS1_V2: @@ -595,4 +632,5 @@ stars[i].measure.dMcal = ps1data[i].dMcal; stars[i].measure.Map = ps1data[i].Map + ZeroPt; + stars[i].measure.dMap = ps1data[i].dM; // a proxy measure stars[i].measure.Mkron = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN; @@ -600,8 +638,10 @@ // these fluxes are converted from counts to counts/sec in FilterStars.c - stars[i].measure.FluxPSF = ps1data[i].Flux; - stars[i].measure.dFluxPSF = ps1data[i].dFlux; + stars[i].measure.FluxPSF = GetFluxFromFluxOrMag (ps1data[i].Flux, ps1data[i].M); + stars[i].measure.dFluxPSF = GetFluxErrFromFluxOrMag (ps1data[i].dFlux, stars[i].measure.FluxPSF, ps1data[i].dM); stars[i].measure.FluxKron = ps1data[i].kronFlux; stars[i].measure.dFluxKron = ps1data[i].kronFluxErr; + stars[i].measure.FluxAp = GetFluxFromFluxOrMag (NAN, ps1data[i].Map); + stars[i].measure.dFluxAp = GetFluxErrFromFluxOrMag (NAN, stars[i].measure.FluxAp, stars[i].measure.dMap); stars[i].measure.Sky = ps1data[i].sky; @@ -676,4 +716,5 @@ stars[i].measure.dMcal = ps1data[i].dMcal; stars[i].measure.Map = ps1data[i].Map + ZeroPt; + stars[i].measure.dMap = (ps1data[i].apFlux > 0.0) ? fabs(ps1data[i].apFluxErr / ps1data[i].apFlux) : NAN; stars[i].measure.Mkron = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN; @@ -681,8 +722,10 @@ // these fluxes are converted from counts to counts/sec in FilterStars.c - stars[i].measure.FluxPSF = ps1data[i].Flux; - stars[i].measure.dFluxPSF = ps1data[i].dFlux; + stars[i].measure.FluxPSF = GetFluxFromFluxOrMag (ps1data[i].Flux, ps1data[i].M); + stars[i].measure.dFluxPSF = GetFluxErrFromFluxOrMag (ps1data[i].dFlux, stars[i].measure.FluxPSF, ps1data[i].dM); stars[i].measure.FluxKron = ps1data[i].kronFlux; stars[i].measure.dFluxKron = ps1data[i].kronFluxErr; + stars[i].measure.FluxAp = GetFluxFromFluxOrMag (ps1data[i].apFlux, ps1data[i].Map); + stars[i].measure.dFluxAp = GetFluxErrFromFluxOrMag (ps1data[i].apFluxErr, stars[i].measure.FluxAp, stars[i].measure.dMap); stars[i].measure.Sky = ps1data[i].sky; @@ -724,5 +767,5 @@ } -Stars *Convert_PS1_SV1 (FTable *table, unsigned int *nstars) { +Stars *Convert_PS1_V5 (FTable *table, unsigned int *nstars) { off_t Nstars; @@ -730,12 +773,7 @@ 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); + CMF_PS1_V5 *ps1data; + + ps1data = gfits_table_get_CMF_PS1_V5 (table, &Nstars, NULL); if (!ps1data) { fprintf (stderr, "skipping inconsistent entry\n"); @@ -763,13 +801,16 @@ stars[i].measure.dMcal = ps1data[i].dMcal; stars[i].measure.Map = ps1data[i].Map + ZeroPt; + stars[i].measure.dMap = (ps1data[i].apFlux > 0.0) ? fabs(ps1data[i].apFluxErr / ps1data[i].apFlux) : NAN; stars[i].measure.Mkron = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN; stars[i].measure.dMkron = (ps1data[i].kronFlux > 0.0) ? ps1data[i].kronFluxErr / ps1data[i].kronFlux : NAN; - + // these fluxes are converted from counts to counts/sec in FilterStars.c - stars[i].measure.FluxPSF = ps1data[i].Flux; - stars[i].measure.dFluxPSF = ps1data[i].dFlux; + stars[i].measure.FluxPSF = GetFluxFromFluxOrMag (ps1data[i].Flux, ps1data[i].M); + stars[i].measure.dFluxPSF = GetFluxErrFromFluxOrMag (ps1data[i].dFlux, stars[i].measure.FluxPSF, ps1data[i].dM); stars[i].measure.FluxKron = ps1data[i].kronFlux; stars[i].measure.dFluxKron = ps1data[i].kronFluxErr; + stars[i].measure.FluxAp = GetFluxFromFluxOrMag (ps1data[i].apFlux, ps1data[i].Map); + stars[i].measure.dFluxAp = GetFluxErrFromFluxOrMag (ps1data[i].apFluxErr, stars[i].measure.FluxAp, stars[i].measure.dMap); stars[i].measure.Sky = ps1data[i].sky; @@ -799,6 +840,9 @@ // 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 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 @@ -808,5 +852,5 @@ } -Stars *Convert_PS1_SV1_Alt (FTable *table, unsigned int *nstars) { +Stars *Convert_PS1_V5_Lensing (FTable *table, unsigned int *nstars) { off_t Nstars; @@ -814,9 +858,7 @@ 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); + CMF_PS1_V5_Lensing *ps1data; + + ps1data = gfits_table_get_CMF_PS1_V5_Lensing (table, &Nstars, NULL); if (!ps1data) { fprintf (stderr, "skipping inconsistent entry\n"); @@ -844,13 +886,230 @@ stars[i].measure.dMcal = ps1data[i].dMcal; stars[i].measure.Map = ps1data[i].Map + ZeroPt; + stars[i].measure.dMap = (ps1data[i].apFlux > 0.0) ? fabs(ps1data[i].apFluxErr / ps1data[i].apFlux) : NAN; stars[i].measure.Mkron = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN; stars[i].measure.dMkron = (ps1data[i].kronFlux > 0.0) ? ps1data[i].kronFluxErr / ps1data[i].kronFlux : NAN; - + // these fluxes are converted from counts to counts/sec in FilterStars.c - stars[i].measure.FluxPSF = ps1data[i].Flux; - stars[i].measure.dFluxPSF = ps1data[i].dFlux; + stars[i].measure.FluxPSF = GetFluxFromFluxOrMag (ps1data[i].Flux, ps1data[i].M); + stars[i].measure.dFluxPSF = GetFluxErrFromFluxOrMag (ps1data[i].dFlux, stars[i].measure.FluxPSF, ps1data[i].dM); stars[i].measure.FluxKron = ps1data[i].kronFlux; stars[i].measure.dFluxKron = ps1data[i].kronFluxErr; + stars[i].measure.FluxAp = GetFluxFromFluxOrMag (ps1data[i].apFlux, ps1data[i].Map); + stars[i].measure.dFluxAp = GetFluxErrFromFluxOrMag (ps1data[i].apFluxErr, stars[i].measure.FluxAp, stars[i].measure.dMap); + + stars[i].measure.Sky = ps1data[i].sky; + stars[i].measure.dSky = ps1data[i].dSky; + + stars[i].measure.psfChisq = ps1data[i].psfChisq; + stars[i].measure.psfQF = ps1data[i].psfQF; + stars[i].measure.psfQFperf = ps1data[i].psfQFperf; + + 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; + + ALLOCATE (stars[i].lensing, Lensing, 1); + dvo_lensing_init (stars[i].lensing); + + stars[i].lensing->X11_sm_obj = ps1data[i].X11_sm_obj; + stars[i].lensing->X12_sm_obj = ps1data[i].X12_sm_obj; + stars[i].lensing->X22_sm_obj = ps1data[i].X22_sm_obj; + stars[i].lensing->E1_sm_obj = ps1data[i].E1_sm_obj; + stars[i].lensing->E2_sm_obj = ps1data[i].E2_sm_obj; + + stars[i].lensing->X11_sh_obj = ps1data[i].X11_sh_obj; + stars[i].lensing->X12_sh_obj = ps1data[i].X12_sh_obj; + stars[i].lensing->X22_sh_obj = ps1data[i].X22_sh_obj; + stars[i].lensing->E1_sh_obj = ps1data[i].E1_sh_obj; + stars[i].lensing->E2_sh_obj = ps1data[i].E2_sh_obj; + + stars[i].lensing->X11_sm_psf = ps1data[i].X11_sm_psf; + stars[i].lensing->X12_sm_psf = ps1data[i].X12_sm_psf; + stars[i].lensing->X22_sm_psf = ps1data[i].X22_sm_psf; + stars[i].lensing->E1_sm_psf = ps1data[i].E1_sm_psf; + stars[i].lensing->E2_sm_psf = ps1data[i].E2_sm_psf; + + stars[i].lensing->X11_sh_psf = ps1data[i].X11_sh_psf; + stars[i].lensing->X12_sh_psf = ps1data[i].X12_sh_psf; + stars[i].lensing->X22_sh_psf = ps1data[i].X22_sh_psf; + stars[i].lensing->E1_sh_psf = ps1data[i].E1_sh_psf; + stars[i].lensing->E2_sh_psf = ps1data[i].E2_sh_psf; + + // stars[i].lensing->F_ApR5 = ps1data[i].F_ApR5; + // stars[i].lensing->dF_ApR5 = ps1data[i].dF_ApR5; + // stars[i].lensing->sF_ApR5 = ps1data[i].sF_ApR5; + // stars[i].lensing->fF_ApR5 = ps1data[i].fF_ApR5; + // + // stars[i].lensing->F_ApR6 = ps1data[i].F_ApR6; + // stars[i].lensing->dF_ApR6 = ps1data[i].dF_ApR6; + // stars[i].lensing->sF_ApR6 = ps1data[i].sF_ApR6; + // stars[i].lensing->fF_ApR6 = ps1data[i].fF_ApR6; + + // this is may optionally be replaced by the internal sequence (see FilterStars.c) + stars[i].lensing->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 (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.dMap = ps1data[i].dM; // a proxy measure + + stars[i].measure.Mkron = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN; + stars[i].measure.dMkron = (ps1data[i].kronFlux > 0.0) ? ps1data[i].kronFluxErr / ps1data[i].kronFlux : NAN; + + // these fluxes are converted from counts to counts/sec in FilterStars.c + stars[i].measure.FluxPSF = GetFluxFromFluxOrMag (ps1data[i].Flux, ps1data[i].M); + stars[i].measure.dFluxPSF = GetFluxErrFromFluxOrMag (ps1data[i].dFlux, stars[i].measure.FluxPSF, ps1data[i].dM); + stars[i].measure.FluxKron = ps1data[i].kronFlux; + stars[i].measure.dFluxKron = ps1data[i].kronFluxErr; + stars[i].measure.FluxAp = GetFluxFromFluxOrMag (NAN, ps1data[i].Map); + stars[i].measure.dFluxAp = GetFluxErrFromFluxOrMag (NAN, stars[i].measure.FluxAp, stars[i].measure.dMap); + + stars[i].measure.Sky = ps1data[i].sky; + stars[i].measure.dSky = ps1data[i].dSky; + + stars[i].measure.psfChisq = ps1data[i].psfChisq; + stars[i].measure.psfQF = ps1data[i].psfQF; + stars[i].measure.psfQFperf = ps1data[i].psfQFperf; + + 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.dMap = ps1data[i].dM; // a proxy measure + + stars[i].measure.Mkron = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN; + stars[i].measure.dMkron = (ps1data[i].kronFlux > 0.0) ? ps1data[i].kronFluxErr / ps1data[i].kronFlux : NAN; + + // these fluxes are converted from counts to counts/sec in FilterStars.c + stars[i].measure.FluxPSF = GetFluxFromFluxOrMag (ps1data[i].Flux, ps1data[i].M); + stars[i].measure.dFluxPSF = GetFluxErrFromFluxOrMag (ps1data[i].dFlux, stars[i].measure.FluxPSF, ps1data[i].dM); + stars[i].measure.FluxKron = ps1data[i].kronFlux; + stars[i].measure.dFluxKron = ps1data[i].kronFluxErr; + stars[i].measure.FluxAp = GetFluxFromFluxOrMag (NAN, ps1data[i].Map); + stars[i].measure.dFluxAp = GetFluxErrFromFluxOrMag (NAN, stars[i].measure.FluxAp, stars[i].measure.dMap); stars[i].measure.Sky = ps1data[i].sky; @@ -921,4 +1180,5 @@ stars[i].measure.dMcal = ps1data[i].dMcal; stars[i].measure.Map = ps1data[i].Map + ZeroPt; + stars[i].measure.dMap = (ps1data[i].apFlux > 0.0) ? fabs(ps1data[i].apFluxErr / ps1data[i].apFlux) : NAN; stars[i].measure.Mkron = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN; @@ -926,8 +1186,10 @@ // these fluxes are converted from counts to counts/sec in FilterStars.c - stars[i].measure.FluxPSF = ps1data[i].Flux; - stars[i].measure.dFluxPSF = ps1data[i].dFlux; + stars[i].measure.FluxPSF = GetFluxFromFluxOrMag (ps1data[i].Flux, ps1data[i].M); + stars[i].measure.dFluxPSF = GetFluxErrFromFluxOrMag (ps1data[i].dFlux, stars[i].measure.FluxPSF, ps1data[i].dM); stars[i].measure.FluxKron = ps1data[i].kronFlux; stars[i].measure.dFluxKron = ps1data[i].kronFluxErr; + stars[i].measure.FluxAp = GetFluxFromFluxOrMag (ps1data[i].apFlux, ps1data[i].Map); + stars[i].measure.dFluxAp = GetFluxErrFromFluxOrMag (ps1data[i].apFluxErr, stars[i].measure.FluxAp, stars[i].measure.dMap); stars[i].measure.Sky = ps1data[i].sky; @@ -966,5 +1228,5 @@ } -Stars *Convert_PS1_DV3 (FTable *table, unsigned int *nstars) { +Stars *Convert_PS1_SV3 (FTable *table, unsigned int *nstars) { off_t Nstars; @@ -972,7 +1234,7 @@ double ZeroPt; Stars *stars; - CMF_PS1_DV3 *ps1data; - - ps1data = gfits_table_get_CMF_PS1_DV3 (table, &Nstars, NULL); + CMF_PS1_SV3 *ps1data; + + ps1data = gfits_table_get_CMF_PS1_SV3 (table, &Nstars, NULL); if (!ps1data) { fprintf (stderr, "skipping inconsistent entry\n"); @@ -1000,13 +1262,98 @@ stars[i].measure.dMcal = ps1data[i].dMcal; stars[i].measure.Map = ps1data[i].Map + ZeroPt; + stars[i].measure.dMap = (ps1data[i].apFlux > 0.0) ? fabs(ps1data[i].apFluxErr / ps1data[i].apFlux) : NAN; stars[i].measure.Mkron = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN; stars[i].measure.dMkron = (ps1data[i].kronFlux > 0.0) ? ps1data[i].kronFluxErr / ps1data[i].kronFlux : NAN; - + // these fluxes are converted from counts to counts/sec in FilterStars.c - stars[i].measure.FluxPSF = ps1data[i].Flux; - stars[i].measure.dFluxPSF = ps1data[i].dFlux; + stars[i].measure.FluxPSF = GetFluxFromFluxOrMag (ps1data[i].Flux, ps1data[i].M); + stars[i].measure.dFluxPSF = GetFluxErrFromFluxOrMag (ps1data[i].dFlux, stars[i].measure.FluxPSF, ps1data[i].dM); stars[i].measure.FluxKron = ps1data[i].kronFlux; stars[i].measure.dFluxKron = ps1data[i].kronFluxErr; + stars[i].measure.FluxAp = GetFluxFromFluxOrMag (ps1data[i].apFlux, ps1data[i].Map); + stars[i].measure.dFluxAp = GetFluxErrFromFluxOrMag (ps1data[i].apFluxErr, stars[i].measure.FluxAp, stars[i].measure.dMap); + + stars[i].measure.Sky = ps1data[i].sky; + stars[i].measure.dSky = ps1data[i].dSky; + + stars[i].measure.psfChisq = ps1data[i].psfChisq; + stars[i].measure.psfQF = ps1data[i].psfQF; + stars[i].measure.psfQFperf = ps1data[i].psfQFperf; + + 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_DV3 (FTable *table, unsigned int *nstars) { + + off_t Nstars; + unsigned int i; + double ZeroPt; + Stars *stars; + CMF_PS1_DV3 *ps1data; + + ps1data = gfits_table_get_CMF_PS1_DV3 (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.dMap = (ps1data[i].apFlux > 0.0) ? fabs(ps1data[i].apFluxErr / ps1data[i].apFlux) : NAN; + + stars[i].measure.Mkron = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN; + stars[i].measure.dMkron = (ps1data[i].kronFlux > 0.0) ? ps1data[i].kronFluxErr / ps1data[i].kronFlux : NAN; + + // these fluxes are converted from counts to counts/sec in FilterStars.c + stars[i].measure.FluxPSF = GetFluxFromFluxOrMag (ps1data[i].Flux, ps1data[i].M); + stars[i].measure.dFluxPSF = GetFluxErrFromFluxOrMag (ps1data[i].dFlux, stars[i].measure.FluxPSF, ps1data[i].dM); + stars[i].measure.FluxKron = ps1data[i].kronFlux; + stars[i].measure.dFluxKron = ps1data[i].kronFluxErr; + stars[i].measure.FluxAp = GetFluxFromFluxOrMag (ps1data[i].apFlux, ps1data[i].Map); + stars[i].measure.dFluxAp = GetFluxErrFromFluxOrMag (ps1data[i].apFluxErr, stars[i].measure.FluxAp, stars[i].measure.dMap); stars[i].measure.Sky = ps1data[i].sky; @@ -1047,4 +1394,90 @@ } - - +Stars *Convert_PS1_DV4 (FTable *table, unsigned int *nstars) { + + off_t Nstars; + unsigned int i; + double ZeroPt; + Stars *stars; + CMF_PS1_DV4 *ps1data; + + ps1data = gfits_table_get_CMF_PS1_DV4 (table, &Nstars, NULL); + if (!ps1data) { + fprintf (stderr, "skipping inconsistent entry\n"); + return (NULL); + } + ZeroPt = GetZeroPoint(); + + fprintf (stderr, "WARNING: Convert_PS1_DV4 not yet updated to match real format\n"); + + 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.dMap = (ps1data[i].apFlux > 0.0) ? fabs(ps1data[i].apFluxErr / ps1data[i].apFlux) : NAN; + + stars[i].measure.Mkron = (ps1data[i].kronFlux > 0.0) ? -2.5*log10(ps1data[i].kronFlux) + ZeroPt : NAN; + stars[i].measure.dMkron = (ps1data[i].kronFlux > 0.0) ? ps1data[i].kronFluxErr / ps1data[i].kronFlux : NAN; + + // these fluxes are converted from counts to counts/sec in FilterStars.c + stars[i].measure.FluxPSF = GetFluxFromFluxOrMag (ps1data[i].Flux, ps1data[i].M); + stars[i].measure.dFluxPSF = GetFluxErrFromFluxOrMag (ps1data[i].dFlux, stars[i].measure.FluxPSF, ps1data[i].dM); + stars[i].measure.FluxKron = ps1data[i].kronFlux; + stars[i].measure.dFluxKron = ps1data[i].kronFluxErr; + stars[i].measure.FluxAp = GetFluxFromFluxOrMag (ps1data[i].apFlux, ps1data[i].Map); + stars[i].measure.dFluxAp = GetFluxErrFromFluxOrMag (ps1data[i].apFluxErr, stars[i].measure.FluxAp, stars[i].measure.dMap); + + stars[i].measure.Sky = ps1data[i].sky; + stars[i].measure.dSky = ps1data[i].dSky; + + stars[i].measure.psfChisq = ps1data[i].psfChisq; + stars[i].measure.psfQF = ps1data[i].psfQF; + stars[i].measure.psfQFperf = ps1data[i].psfQFperf; + 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/SEDfit.c =================================================================== --- /trunk/Ohana/src/addstar/src/SEDfit.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/SEDfit.c (revision 37036) @@ -162,9 +162,12 @@ } + double R = incat[0].average[i].R; + double D = incat[0].average[i].D; + for (j = 0; valid && (j < Nmodel); j++) { n = modelRow[j]; dvo_measure_init (&outcat[0].measure[Nmeas]); - outcat[0].measure[Nmeas].dR = 0.0; - outcat[0].measure[Nmeas].dD = 0.0; + outcat[0].measure[Nmeas].R = R; + outcat[0].measure[Nmeas].D = D; outcat[0].measure[Nmeas].M = table[0].row[minFit.row][0].mags[n] + minFit.Md; outcat[0].measure[Nmeas].dM = 0.0; Index: /trunk/Ohana/src/addstar/src/StarOps.c =================================================================== --- /trunk/Ohana/src/addstar/src/StarOps.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/StarOps.c (revision 37036) @@ -8,4 +8,5 @@ star[0].found = -1; // found == -1 -> not yet found (use enums?) + star[0].lensing = NULL; // we only populate this if needed return TRUE; } Index: /trunk/Ohana/src/addstar/src/addstar.c =================================================================== --- /trunk/Ohana/src/addstar/src/addstar.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/addstar.c (revision 37036) @@ -7,5 +7,5 @@ int Nmatch, status, loadObjects; off_t i, Nimages; - off_t Naverage, Nmeasure; + off_t Naverage, Nmeasure, Nlensing; Stars *stars, **subset; Image *images; @@ -100,5 +100,5 @@ /* match stars to existing catalog data (or otherwise manipulate catalog data) */ - Nmatch = Naverage = Nmeasure = 0; + Nmatch = Naverage = Nmeasure = Nlensing = 0; for (i = 0; loadObjects && (i < skylist[0].Nregions); i++) { @@ -107,5 +107,5 @@ catalog.catformat = dvo_catalog_catformat (CATFORMAT); // set the default catformat from config data catalog.catmode = dvo_catalog_catmode (CATMODE); // set the default catmode from config data - catalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF; + catalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF | LOAD_LENSING; catalog.Nsecfilt = GetPhotcodeNsecfilt (); @@ -182,4 +182,5 @@ Naverage += catalog.Naverage; Nmeasure += catalog.Nmeasure; + Nlensing += catalog.Nlensing; // write out catalog, if appropriate @@ -187,5 +188,5 @@ SetProtect (TRUE); if (options.update) { - catalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF; + catalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF | LOAD_LENSING; dvo_catalog_update (&catalog, VERBOSE); } else { @@ -231,5 +232,5 @@ gettimeofday (&stop, NULL); dtime = DTIME (stop, start); - fprintf (stderr, "SUCCESS: elapsed time %9.4f sec for %5d stars (%5d matches), "OFF_T_FMT" average, "OFF_T_FMT" measure\n", dtime, Nstars, Nmatch, Naverage, Nmeasure); + fprintf (stderr, "SUCCESS: elapsed time %9.4f sec for %5d stars (%5d matches), "OFF_T_FMT" average, "OFF_T_FMT" measure, "OFF_T_FMT" lensing\n", dtime, Nstars, Nmatch, Naverage, Nmeasure, Nlensing); exit (0); Index: /trunk/Ohana/src/addstar/src/build_links.c =================================================================== --- /trunk/Ohana/src/addstar/src/build_links.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/build_links.c (revision 37036) @@ -225,2 +225,141 @@ } +/*******************************************************************************************/ + +/* build the initial links assuming the table is sorted, + not partial, and has a correct set of average[].lensingOffset,Nlensing values */ +off_t *init_lensing_links (Average *average, off_t Naverage, Lensing *lensing, off_t Nlensing) { + + off_t i, j, N; + off_t *next_lens; + + N = 0; + + ALLOCATE (next_lens, off_t, Nlensing); + for (i = 0; i < Naverage; i++) { + if (!average[i].Nlensing) continue; + for (j = 0; j < average[i].Nlensing - 1; j++, N++) { + next_lens[N] = N + 1; + if (N >= Nlensing) { + fprintf (stderr, "WARNING: N out of bounds (1)\n"); + } + } + next_lens[N] = -1; + if (N >= Nlensing) { + fprintf (stderr, "WARNING: N out of bounds (2)\n"); + } + + if (N >= Nlensing) { + fprintf (stderr, "overflow in init_lensing_links\n"); + abort (); + } + N++; + } + return (next_lens); +} + +/* construct lensing links which are valid FOR THIS LOAD + * - if we have a full load, we will get links which can + * be used by other programs (eg, relphot, etc) + * - if we have a partial load, the links are only valid + * for that partial load + */ + +off_t *build_lensing_links (Average *average, off_t Naverage, Lensing *lensing, off_t Nlensing) { + + off_t i, m, k, Nm, averef; + off_t *next_lens; + + ALLOCATE (next_lens, off_t, Nlensing); + + /* reset the Nm, offset values for average */ + for (i = 0; i < Naverage; i++) { + average[i].lensingOffset = -1; + average[i].Nlensing = 0; + } + + for (Nm = 0; Nm < Nlensing; Nm++) { + averef = lensing[Nm].averef; + m = average[averef].lensingOffset; + next_lens[Nm] = -1; + + if (m == -1) { /* no links yet for source */ + average[averef].lensingOffset = Nm; + average[averef].Nlensing = 1; + continue; + } + + for (k = 0; next_lens[m] != -1; k++) { + m = next_lens[m]; + if (m >= Nlensing) { + fprintf (stderr, "WARNING: m out of bounds (1)\n"); + } + } + + average[averef].Nlensing = k + 2; + next_lens[m] = Nm; + if (m >= Nlensing) { + fprintf (stderr, "WARNING: m out of bounds (2)\n"); + } + } + return (next_lens); +} + +/* average[].lensingOffset, average[].Nlensing are valid within an addstar run */ +off_t add_lens_link (Average *average, off_t *next_lens, off_t Nlensing, off_t NLENSING) { + + off_t k, m; + + /* if we have trouble, check validity of next_lens[m] : m < Nlensing */ + m = average[0].lensingOffset; + + for (k = 0; k < average[0].Nlensing - 1; k++) { + m = next_lens[m]; + if (m >= NLENSING) { + fprintf (stderr, "WARNING: m out of bounds (3)\n"); + } + } + + /* set up references */ + next_lens[Nlensing] = -1; + if (Nlensing >= NLENSING) { + fprintf (stderr, "WARNING: Nlensing out of bounds (1)\n"); + } + + if (m == -1) { + average[0].lensingOffset = Nlensing; + } else { + next_lens[m] = Nlensing; + if (m >= NLENSING) { + fprintf (stderr, "WARNING: m out of bounds (4)\n"); + } + } + + return (TRUE); +} + +Lensing *sort_lensing (Average *average, off_t Naverage, Lensing *lensing, off_t Nlensing, off_t *next_lens) { + + off_t i, k, n, N; + Lensing *tmplensing; + + /* fix order of Lensing (memory intensive, but fast) */ + N = 0; + ALLOCATE (tmplensing, Lensing, Nlensing); + for (i = 0; i < Naverage; i++) { + if (!average[i].Nlensing) continue; + n = average[i].lensingOffset; + average[i].lensingOffset = N; + for (k = 0; k < average[i].Nlensing; k++, N++) { + if (n == -1) abort(); + tmplensing[N] = lensing[n]; + if (lensing[n].averef != i) abort(); + tmplensing[N].averef = i; + n = next_lens[n]; + } + } + free (lensing); + return (tmplensing); +} + + Index: /trunk/Ohana/src/addstar/src/fakeimage.c =================================================================== --- /trunk/Ohana/src/addstar/src/fakeimage.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/fakeimage.c (revision 37036) @@ -6,5 +6,4 @@ double pltscale, pixscale; double Rmin, Rmax, Dmin, Dmax; - // double Xmin, Xmax, Ymin, Ymax; double dX, dY, r, d; char chipname[80], chipdata[256], name[80]; @@ -168,6 +167,4 @@ image[0].cerror = 0.0; - // RD_to_XY (&Xmax, &Ymax, Rmax, Dmax, MOSAIC); - // RD_to_XY (&Xmin, &Ymin, Rmin, Dmin, MOSAIC); image[0].NX = Rmax - Rmin; image[0].NY = Dmax - Dmin; Index: /trunk/Ohana/src/addstar/src/find_matches.c =================================================================== --- /trunk/Ohana/src/addstar/src/find_matches.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/find_matches.c (revision 37036) @@ -14,6 +14,6 @@ if (NSTAR_GROUP <= 0) { - fprintf (stderr, "ERROR: NSTAR_GROUP NOT SET!\n"); - exit (1); + fprintf (stderr, "ERROR: NSTAR_GROUP NOT SET!\n"); + exit (1); } @@ -36,5 +36,5 @@ ALLOCATE (Y2, double, NAVE); ALLOCATE (N2, off_t, NAVE); - ALLOCATE (catalog[0].found, off_t, NAVE); + ALLOCATE (catalog[0].found_t, off_t, NAVE); /* for secfilt j and star i, secfilt[i*Nsecfilt+j] */ @@ -77,5 +77,4 @@ if (Nstars < 1) { if (VERBOSE) fprintf (stderr, "skipping %s, no overlapping stars\n", catalog[0].filename); - free (catalog[0].found); free (X1); free (Y1); @@ -92,5 +91,5 @@ RD_to_XY (&X2[i], &Y2[i], catalog[0].average[i].R, catalog[0].average[i].D, &tcoords); N2[i] = i; - catalog[0].found[N2[i]] = -1; + catalog[0].found_t[N2[i]] = -1; } if (Nave > 1) sort_coords_index (X2, Y2, N2, Nave); @@ -160,6 +159,6 @@ // the following measure elements cannot be set until here: - catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N].average.R); - catalog[0].measure[Nmeas].dD = 3600.0*(catalog[0].average[n].D - stars[N].average.D); + catalog[0].measure[Nmeas].R = stars[N].average.R; + catalog[0].measure[Nmeas].D = stars[N].average.D; catalog[0].measure[Nmeas].dbFlags = 0; catalog[0].measure[Nmeas].averef = n; // this must be an absolute sequence number, if partial average is loaded @@ -167,16 +166,21 @@ catalog[0].measure[Nmeas].catID = catalog[0].catID; - // rationalize dR: - if (catalog[0].measure[Nmeas].dR > +180.0*3600.0) { - // average on high end of boundary, move star up - stars[N].average.R += 360.0; - catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N].average.R); - } - if (catalog[0].measure[Nmeas].dR < -180.0*3600.0) { - // average on low end of boundary, move star down - stars[N].average.R -= 360.0; - catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N].average.R); - } - if (fabs(catalog[0].measure[Nmeas].dR) > 10*RADIUS) { + float dRoff = dvoOffsetR(&catalog[0].measure[Nmeas], &catalog[0].average[n]); + + // rationalize R: + if (dRoff > +180.0*3600.0) { + // average on high end of boundary, move star up + catalog[0].measure[Nmeas].R += 360.0; + dRoff -= 360.0*3600.0; + } + if (dRoff < -180.0*3600.0) { + // average on low end of boundary, move star down + catalog[0].measure[Nmeas].R -= 360.0; + dRoff += 360.0*3600.0; + } + if (fabs(dRoff) > 10*RADIUS) { + // take declination into account and check again. + double cosD = cos(RAD_DEG*catalog[0].average[n].D); + if (fabs(dRoff*cosD) > 10*RADIUS) { fprintf (stderr, "error: %10.6f,%10.6f vs %10.6f,%10.6f (%f,%f vs %f,%f)\n", catalog[0].average[n].R, catalog[0].average[n].D, @@ -184,4 +188,5 @@ X1[i], X2[J], Y1[i], Y2[J]); + } } @@ -196,5 +201,5 @@ if (Nsec > -1) { if (isnan(catalog[0].secfilt[n*Nsecfilt+Nsec].M)) { - catalog[0].secfilt[n*Nsecfilt+Nsec].M = PhotCat (&catalog[0].measure[Nmeas]); + catalog[0].secfilt[n*Nsecfilt+Nsec].M = PhotCat (&catalog[0].measure[Nmeas], MAG_CLASS_PSF); } } @@ -213,9 +218,9 @@ } /* this catalog star matches more than one image star */ - if (catalog[0].found[n] > -1) { - catalog[0].measure[catalog[0].found[n]].dbFlags |= ID_MEAS_BLEND_OBJ; + if (catalog[0].found_t[n] > -1) { + catalog[0].measure[catalog[0].found_t[n]].dbFlags |= ID_MEAS_BLEND_OBJ; catalog[0].measure[Nmeas].dbFlags |= ID_MEAS_BLEND_OBJ; } else { - catalog[0].found[n] = Nmeas; + catalog[0].found_t[n] = Nmeas; } /* Nm is updated, but not written out in -update mode (for existing entries) @@ -232,6 +237,6 @@ } - /* incorporate unmatched image stars, if this star is in field of this catalog */ - /* these new entries are all written out in UPDATE mode */ +/* incorporate unmatched image stars, if this star is in field of this catalog */ +/* these new entries are all written out in UPDATE mode */ for (i = 0; (i < Nstars) && !options.only_match; i += NSTAR_GROUP) { /* make sure there is space for next entry */ @@ -260,5 +265,5 @@ if (PSPS_ID) { - catalog[0].average[Nave].extID = CreatePSPSObjectID(catalog[0].average[Nave].R, catalog[0].average[Nave].D); + catalog[0].average[Nave].extID = CreatePSPSObjectID(catalog[0].average[Nave].R, catalog[0].average[Nave].D); } @@ -274,6 +279,6 @@ // the following measure elements cannot be set until here: - catalog[0].measure[Nmeas].dR = 0.0; - catalog[0].measure[Nmeas].dD = 0.0; + catalog[0].measure[Nmeas].R = stars[i].average.R; + catalog[0].measure[Nmeas].D = stars[i].average.D; catalog[0].measure[Nmeas].dbFlags = 0; catalog[0].measure[Nmeas].averef = Nave; // XXX EAM : must be absolute Nave if partial read @@ -284,9 +289,9 @@ /* in UPDATE mode, this value is not saved; use relphot to recalculate */ if (Nsec > -1) { - catalog[0].secfilt[Nave*Nsecfilt+Nsec].M = PhotCat (&catalog[0].measure[Nmeas]); + catalog[0].secfilt[Nave*Nsecfilt+Nsec].M = PhotCat (&catalog[0].measure[Nmeas], MAG_CLASS_PSF); } /* next[Nmeas] should always be -1 in this context (it is always the only - measurement for the star) */ + measurement for the star) */ stars[i].found = Nmeas; next_meas[Nmeas] = -1; // initial value here update below @@ -309,5 +314,5 @@ } - /* note stars which have been found in this catalog */ +/* note stars which have been found in this catalog */ for (i = 0; i < NstarsIn; i++) { if (stars[i].found > -1) { @@ -316,5 +321,5 @@ } - /* check if the catalog has changed? if no change, no need to write */ +/* check if the catalog has changed? if no change, no need to write */ catalog[0].objID = objID; // new max value, save on catalog close catalog[0].Naverage = Nave; @@ -323,5 +328,4 @@ if (VERBOSE) fprintf (stderr, "Nstars, Nave, Nmeas: "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT", ("OFF_T_FMT" matches)\n", Nstars, Nave, Nmeas, Nmatch); - free (catalog[0].found); free (X1); free (Y1); Index: /trunk/Ohana/src/addstar/src/find_matches_closest.c =================================================================== --- /trunk/Ohana/src/addstar/src/find_matches_closest.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/find_matches_closest.c (revision 37036) @@ -7,6 +7,6 @@ double *X1, *Y1, *X2, *Y2; double dX, dY, dR; - off_t *N1, *N2, *next_meas; - off_t Nave, NAVE, Nmeas, NMEAS, Nmatch; + off_t *N1, *N2, *next_meas, *next_lens; + off_t Nave, NAVE, Nmeas, NMEAS, Nmatch, Nlens, NLENS; int Nsecfilt, Nsec; unsigned int objID, catID; @@ -36,5 +36,5 @@ ALLOCATE (Y2, double, NAVE); ALLOCATE (N2, off_t, NAVE); - ALLOCATE (catalog[0].found, off_t, NAVE); + ALLOCATE (catalog[0].found_t, off_t, NAVE); /* for secfilt j and star i, secfilt[i*Nsecfilt+j] */ @@ -42,4 +42,5 @@ Nmatch = 0; NMEAS = Nmeas = catalog[0].Nmeasure; + NLENS = Nlens = catalog[0].Nlensing; // current max obj ID for this catalog @@ -76,5 +77,4 @@ if (Nstars < 1) { if (VERBOSE) fprintf (stderr, "skipping %s, no overlapping stars\n", catalog[0].filename); - free (catalog[0].found); free (X1); free (Y1); @@ -91,5 +91,5 @@ RD_to_XY (&X2[i], &Y2[i], catalog[0].average[i].R, catalog[0].average[i].D, &tcoords); N2[i] = i; - catalog[0].found[N2[i]] = -1; + catalog[0].found_t[N2[i]] = -1; } if (Nave > 1) sort_coords_index (X2, Y2, N2, Nave); @@ -100,6 +100,8 @@ // is sorted while processed next_meas = init_measure_links (catalog[0].average, Nave, catalog[0].measure, Nmeas); + next_lens = init_lensing_links (catalog[0].average, Nave, catalog[0].lensing, Nlens); } else { next_meas = build_measure_links (catalog[0].average, Nave, catalog[0].measure, Nmeas); + next_lens = build_lensing_links (catalog[0].average, Nave, catalog[0].lensing, Nlens); } @@ -177,4 +179,9 @@ REALLOCATE (catalog[0].measure, Measure, NMEAS); } + if (Nlens >= NLENS) { + NLENS = Nlens + 1000; + REALLOCATE (next_lens, off_t, NLENS); + REALLOCATE (catalog[0].lensing, Lensing, NLENS); + } /* add to end of measurement list */ @@ -186,7 +193,8 @@ catalog[0].measure[Nmeas] = stars[N].measure; - /** dR,dD now represent arcsec **/ - catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N].average.R); - catalog[0].measure[Nmeas].dD = 3600.0*(catalog[0].average[n].D - stars[N].average.D); + // measure now carries R,D (not dR,dD) + // note that ReadStarsFITS does not set measure.R,D and average.R,D + catalog[0].measure[Nmeas].R = stars[N].average.R; + catalog[0].measure[Nmeas].D = stars[N].average.D; catalog[0].measure[Nmeas].dbFlags = 0; catalog[0].measure[Nmeas].averef = n; @@ -194,22 +202,40 @@ catalog[0].measure[Nmeas].catID = catalog[0].catID; + float dRoff = dvoOffsetR(&catalog[0].measure[Nmeas], &catalog[0].average[n]); + // rationalize dR - if (catalog[0].measure[Nmeas].dR > +180.0*3600.0) { + if (dRoff > +180.0*3600.0) { // average on high end of boundary, move star up - stars[N].average.R += 360.0; - catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N].average.R); - } - if (catalog[0].measure[Nmeas].dR < -180.0*3600.0) { + catalog[0].measure[Nmeas].R += 360.0; + dRoff -= 360.0*3600.0; + } + if (dRoff < -180.0*3600.0) { // average on low end of boundary, move star down - stars[N].average.R -= 360.0; - catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N].average.R); - } - if (fabs(catalog[0].measure[Nmeas].dR) > 10*RADIUS) { - fprintf (stderr, "error: %10.6f,%10.6f vs %10.6f,%10.6f (%f,%f vs %f,%f)\n", - catalog[0].average[n].R, catalog[0].average[n].D, - stars[N].average.R, stars[N].average.D, - X1[i], X2[Jmin], - Y1[i], Y2[Jmin]); - // XXX abort on this? -- this is a bad failure... + catalog[0].measure[Nmeas].R -= 360.0; + dRoff += 360.0*3600.0; + } + if (fabs(dRoff) > 10*RADIUS) { + // take declination into account and check again. + double cosD = cos(RAD_DEG*catalog[0].average[n].D); + if (fabs(dRoff*cosD) > 10*RADIUS) { + fprintf (stderr, "error: %10.6f,%10.6f vs %10.6f,%10.6f (%f,%f vs %f,%f)\n", + catalog[0].average[n].R, catalog[0].average[n].D, + stars[N].average.R, stars[N].average.D, + X1[i], X2[Jmin], + Y1[i], Y2[Jmin]); + // XXX abort on this? -- this is a bad failure... + } + } + + // add the lensing values if they exist + if (stars[N].lensing) { + add_lens_link (&catalog[0].average[n], next_lens, Nlens, NLENS); // ? + catalog[0].lensing[Nlens] = stars[N].lensing[0]; + + catalog[0].lensing[Nlens].averef = n; + catalog[0].lensing[Nlens].objID = catalog[0].average[n].objID; + catalog[0].lensing[Nlens].catID = catalog[0].catID; + catalog[0].average[n].Nlensing ++; + Nlens ++; } @@ -224,5 +250,5 @@ if (Nsec > -1) { if (isnan(catalog[0].secfilt[n*Nsecfilt+Nsec].M)) { - catalog[0].secfilt[n*Nsecfilt+Nsec].M = PhotCat (&catalog[0].measure[Nmeas]); + catalog[0].secfilt[n*Nsecfilt+Nsec].M = PhotCat (&catalog[0].measure[Nmeas], MAG_CLASS_PSF); } } @@ -234,5 +260,5 @@ Nm is recalculated in build_meas_links if loaded table is not sorted */ stars[N].found = Nmeas; - catalog[0].found[n] = Nmeas; + catalog[0].found_t[n] = Nmeas; catalog[0].average[n].Nmeasure ++; Nmeas ++; @@ -249,4 +275,9 @@ REALLOCATE (catalog[0].measure, Measure, NMEAS); } + if (Nlens >= NLENS - NSTAR_GROUP) { + NLENS = Nlens + 1000; + REALLOCATE (next_lens, off_t, NLENS); + REALLOCATE (catalog[0].lensing, Lensing, NLENS); + } if (Nave >= NAVE) { NAVE = Nave + 1000; @@ -278,30 +309,54 @@ for (j = 0; j < NSTAR_GROUP; j++) { - // supply the measurments from this detection - catalog[0].measure[Nmeas] = stars[i + j].measure; - - // the following measure elements cannot be set until here: - catalog[0].measure[Nmeas].dR = 0.0; // astrometric offset, not error - catalog[0].measure[Nmeas].dD = 0.0; // astrometric offset, not error - catalog[0].measure[Nmeas].dbFlags = 0; - catalog[0].measure[Nmeas].averef = Nave; - catalog[0].measure[Nmeas].objID = catalog[0].average[Nave].objID; - catalog[0].measure[Nmeas].catID = catalog[0].catID; - - /* set the average magnitude if not already set and the photcode.equiv is not 0 */ - /* in UPDATE mode, this value is not saved; use relphot to recalculate */ - if (Nsec > -1) { - catalog[0].secfilt[Nave*Nsecfilt+Nsec].M = PhotCat (&catalog[0].measure[Nmeas]); - } - - /* next[Nmeas] should always be -1 in this context (it is always the only - measurement for the star) */ - stars[i+j].found = Nmeas; - next_meas[Nmeas] = -1; // inital value here update below - Nmeas ++; + // supply the measurments from this detection + dvo_measure_init (&catalog[0].measure[Nmeas]); + catalog[0].measure[Nmeas] = stars[i + j].measure; + + // the following measure elements cannot be set until here: + catalog[0].measure[Nmeas].R = stars[i].average.R; + catalog[0].measure[Nmeas].D = stars[i].average.D; + catalog[0].measure[Nmeas].dbFlags = 0; + catalog[0].measure[Nmeas].averef = Nave; + catalog[0].measure[Nmeas].objID = catalog[0].average[Nave].objID; + catalog[0].measure[Nmeas].catID = catalog[0].catID; + + /* set the average magnitude if not already set and the photcode.equiv is not 0 */ + /* in UPDATE mode, this value is not saved; use relphot to recalculate */ + if (Nsec > -1) { + catalog[0].secfilt[Nave*Nsecfilt+Nsec].M = PhotCat (&catalog[0].measure[Nmeas], MAG_CLASS_PSF); + } + + /* next[Nmeas] should always be -1 in this context (it is always the only + measurement for the star) */ + stars[i+j].found = Nmeas; + next_meas[Nmeas] = -1; // inital value here update below + Nmeas ++; } for (j = 0; j < NSTAR_GROUP - 1; j++) { - next_meas[Nmeas - NSTAR_GROUP + j] = Nmeas - NSTAR_GROUP + j + 1; - } + next_meas[Nmeas - NSTAR_GROUP + j] = Nmeas - NSTAR_GROUP + j + 1; + } + + // if we have lensing data, insert that as well + if (stars[i].lensing) { + catalog[0].average[Nave].Nlensing = NSTAR_GROUP; + catalog[0].average[Nave].lensingOffset = Nlens; + for (j = 0; j < NSTAR_GROUP; j++) { + // add the lensing values if they exist + if (stars[i + j].lensing) { + dvo_lensing_init (&catalog[0].lensing[Nlens]); + catalog[0].lensing[Nlens] = stars[i + j].lensing[0]; + + catalog[0].lensing[Nlens].averef = Nave; + catalog[0].lensing[Nlens].objID = catalog[0].average[Nave].objID; + catalog[0].lensing[Nlens].catID = catalog[0].catID; + next_lens[Nlens] = -1; + Nlens ++; + } + } + for (j = 0; j < NSTAR_GROUP - 1; j++) { + next_lens[Nlens - NSTAR_GROUP + j] = Nlens - NSTAR_GROUP + j + 1; + } + } + Nave ++; } @@ -309,4 +364,5 @@ REALLOCATE (catalog[0].average, Average, Nave); REALLOCATE (catalog[0].measure, Measure, Nmeas); + REALLOCATE (catalog[0].lensing, Lensing, Nlens); if (options.nosort) { @@ -315,4 +371,5 @@ catalog[0].sorted = TRUE; catalog[0].measure = sort_measure (catalog[0].average, Nave, catalog[0].measure, Nmeas, next_meas); + catalog[0].lensing = sort_lensing (catalog[0].average, Nave, catalog[0].lensing, Nlens, next_lens); } @@ -328,8 +385,8 @@ catalog[0].Naverage = Nave; catalog[0].Nmeasure = Nmeas; + catalog[0].Nlensing = Nlens; catalog[0].Nsecf_mem = Nave*Nsecfilt; - if (VERBOSE) fprintf (stderr, "Nstars, Nave, Nmeas: "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT", ("OFF_T_FMT" matches)\n", Nstars, Nave, Nmeas, Nmatch); - - free (catalog[0].found); + if (VERBOSE) fprintf (stderr, "Nstars, Nave, Nmeas, Nlens: "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT", ("OFF_T_FMT" matches)\n", Nstars, Nave, Nmeas, Nlens, Nmatch); + free (X1); free (Y1); @@ -339,4 +396,5 @@ free (Y2); free (next_meas); + free (next_lens); return (Nmatch); Index: /trunk/Ohana/src/addstar/src/find_matches_closest_refstars.c =================================================================== --- /trunk/Ohana/src/addstar/src/find_matches_closest_refstars.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/find_matches_closest_refstars.c (revision 37036) @@ -31,5 +31,5 @@ ALLOCATE (Y2, double, NAVE); ALLOCATE (N2, off_t, NAVE); - ALLOCATE (catalog[0].found, off_t, NAVE); + ALLOCATE (catalog[0].found_t, off_t, NAVE); /* for secfilt j and star i, secfilt[i*Nsecfilt+j] */ @@ -71,5 +71,4 @@ if (Nstars < 1) { if (VERBOSE) fprintf (stderr, "skipping %s, no overlapping stars\n", catalog[0].filename); - free (catalog[0].found); free (X1); free (Y1); @@ -86,5 +85,5 @@ RD_to_XY (&X2[i], &Y2[i], catalog[0].average[i].R, catalog[0].average[i].D, &tcoords); N2[i] = i; - catalog[0].found[N2[i]] = -1; + catalog[0].found_t[N2[i]] = -1; } if (Nave > 1) sort_coords_index (X2, Y2, N2, Nave); @@ -180,6 +179,6 @@ /** *** dR,dD now in arcsec *** **/ - catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N][0].average.R); - catalog[0].measure[Nmeas].dD = 3600.0*(catalog[0].average[n].D - stars[N][0].average.D); + catalog[0].measure[Nmeas].R = stars[N][0].average.R; + catalog[0].measure[Nmeas].D = stars[N][0].average.D; catalog[0].measure[Nmeas].dbFlags = 0; catalog[0].measure[Nmeas].averef = n; @@ -187,14 +186,16 @@ catalog[0].measure[Nmeas].catID = catalog[0].catID; + float dRoff = dvoOffsetR(&catalog[0].measure[Nmeas], &catalog[0].average[n]); + // rationalize dR: - if (catalog[0].measure[Nmeas].dR > +180.0*3600.0) { + if (dRoff > +180.0*3600.0) { // average on high end of boundary, move star up - stars[N][0].average.R += 360.0; - catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N][0].average.R); - } - if (catalog[0].measure[Nmeas].dR < -180.0*3600.0) { + catalog[0].measure[Nmeas].R += 360.0; + dRoff -= 360.0*3600.0; + } + if (dRoff < -180.0*3600.0) { // average on low end of boundary, move star down - stars[N][0].average.R -= 360.0; - catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N][0].average.R); + catalog[0].measure[Nmeas].R -= 360.0; + dRoff += 360.0*3600.0; } @@ -218,9 +219,9 @@ /* this catalog star matches more than one image star */ // XXX should this be an average flag? - if (catalog[0].found[n] > -1) { - catalog[0].measure[catalog[0].found[n]].dbFlags |= ID_MEAS_BLEND_OBJ; + if (catalog[0].found_t[n] > -1) { + catalog[0].measure[catalog[0].found_t[n]].dbFlags |= ID_MEAS_BLEND_OBJ; catalog[0].measure[Nmeas].dbFlags |= ID_MEAS_BLEND_OBJ; } else { - catalog[0].found[n] = Nmeas; + catalog[0].found_t[n] = Nmeas; } @@ -290,6 +291,6 @@ catalog[0].measure[Nmeas] = stars[N][0].measure; - catalog[0].measure[Nmeas].dR = 0.0; // astrometric offset, not error - catalog[0].measure[Nmeas].dD = 0.0; // astrometric offset, not error + catalog[0].measure[Nmeas].R = catalog[0].average[Nave].R; + catalog[0].measure[Nmeas].D = catalog[0].average[Nave].D; catalog[0].measure[Nmeas].dbFlags = 0; catalog[0].measure[Nmeas].averef = Nave; @@ -329,5 +330,4 @@ if (VERBOSE) fprintf (stderr, "Nstars, Nave, Nmeas: "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT", ("OFF_T_FMT" matches)\n", Nstars, Nave, Nmeas, Nmatch); - free (catalog[0].found); free (X1); free (Y1); Index: /trunk/Ohana/src/addstar/src/find_matches_refstars.c =================================================================== --- /trunk/Ohana/src/addstar/src/find_matches_refstars.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/find_matches_refstars.c (revision 37036) @@ -37,5 +37,5 @@ ALLOCATE (Y2, double, NAVE); ALLOCATE (N2, off_t, NAVE); - ALLOCATE (catalog[0].found, off_t, NAVE); + ALLOCATE (catalog[0].found_t, off_t, NAVE); REALLOCATE (catalog[0].average, Average, NAVE); REALLOCATE (catalog[0].secfilt, SecFilt, NAVE*catalog[0].Nsecfilt); @@ -69,5 +69,5 @@ RD_to_XY (&X2[i], &Y2[i], catalog[0].average[i].R, catalog[0].average[i].D, &tcoords); N2[i] = i; - catalog[0].found[N2[i]] = -1; + catalog[0].found_t[N2[i]] = -1; } if (Nave > 1) sort_coords_index (X2, Y2, N2, Nave); @@ -141,7 +141,7 @@ catalog[0].measure[Nmeas] = stars[N][0].measure; - /** *** dR,dD now in arcsec *** **/ - catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N][0].average.R); - catalog[0].measure[Nmeas].dD = 3600.0*(catalog[0].average[n].D - stars[N][0].average.D); + /** measure now stores R,D **/ + catalog[0].measure[Nmeas].R = stars[N][0].average.R; + catalog[0].measure[Nmeas].D = stars[N][0].average.D; catalog[0].measure[Nmeas].dbFlags = 0; catalog[0].measure[Nmeas].averef = n; @@ -149,14 +149,16 @@ catalog[0].measure[Nmeas].catID = catalog[0].catID; + float dRoff = dvoOffsetR(&catalog[0].measure[Nmeas], &catalog[0].average[n]); + // rationalize dR: - if (catalog[0].measure[Nmeas].dR > +180.0*3600.0) { + if (dRoff > +180.0*3600.0) { // average on high end of boundary, move star up - stars[N][0].average.R += 360.0; - catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N][0].average.R); - } - if (catalog[0].measure[Nmeas].dR < -180.0*3600.0) { + catalog[0].measure[Nmeas].R += 360.0; + dRoff -= 360.0*3600.0; + } + if (dRoff < -180.0*3600.0) { // average on low end of boundary, move star down - stars[N][0].average.R -= 360.0; - catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N][0].average.R); + catalog[0].measure[Nmeas].R -= 360.0; + dRoff += 360.0*3600.0; } @@ -188,9 +190,9 @@ } /* this catalog star matches more than one image star */ - if (catalog[0].found[n] > -1) { - catalog[0].measure[catalog[0].found[n]].dbFlags |= ID_MEAS_BLEND_OBJ; + if (catalog[0].found_t[n] > -1) { + catalog[0].measure[catalog[0].found_t[n]].dbFlags |= ID_MEAS_BLEND_OBJ; catalog[0].measure[Nmeas].dbFlags |= ID_MEAS_BLEND_OBJ; } else { - catalog[0].found[n] = Nmeas; + catalog[0].found_t[n] = Nmeas; } @@ -262,6 +264,6 @@ catalog[0].measure[Nmeas] = stars[N][0].measure; - catalog[0].measure[Nmeas].dR = 0.0; - catalog[0].measure[Nmeas].dD = 0.0; + catalog[0].measure[Nmeas].R = stars[N][0].average.R; + catalog[0].measure[Nmeas].D = stars[N][0].average.D; catalog[0].measure[Nmeas].t = (stars[N][0].measure.t == 0) ? TIMEREF : stars[N][0].measure.t; /** careful : time_t vs e_time **/ @@ -303,5 +305,4 @@ if (VERBOSE) fprintf (stderr, "Nstars, Nave, Nmeas: %d "OFF_T_FMT" "OFF_T_FMT", ("OFF_T_FMT" matches)\n", Nstars, Nave, Nmeas, Nmatch); - free (catalog[0].found); free (X1); free (Y1); Index: /trunk/Ohana/src/addstar/src/findskycell.c =================================================================== --- /trunk/Ohana/src/addstar/src/findskycell.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/findskycell.c (revision 37036) @@ -161,10 +161,4 @@ exit (0); } - -# define MARKTIME(MSG,...) { \ - float dtime; \ - gettimeofday (&stop, (void *) NULL); \ - dtime = DTIME (stop, start); \ - fprintf (stderr, MSG, __VA_ARGS__); } int mktree (char *treefile, char *catdir) { @@ -402,6 +396,5 @@ } - struct timeval start, stop; - gettimeofday (&start, (void *) NULL); + INITTIME; int Npts = 10000000; Index: /trunk/Ohana/src/addstar/src/load2mass_catalog.c =================================================================== --- /trunk/Ohana/src/addstar/src/load2mass_catalog.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/load2mass_catalog.c (revision 37036) @@ -17,10 +17,13 @@ for (i = 0; i < Nstars; i+=3) { + double R = stars[i].average.R; + double D = stars[i].average.D; + // construct an average object for this object // XXX for now, the output objects will have limited astrometric interpretation... // XXX every 3 stars represents 3 measurements and 1 average dvo_average_init (&catalog[0].average[Nave]); - catalog[0].average[Nave].R = stars[i].average.R; - catalog[0].average[Nave].D = stars[i].average.D; + catalog[0].average[Nave].R = R; + catalog[0].average[Nave].D = D; catalog[0].average[Nave].measureOffset = Nmeas; @@ -33,6 +36,6 @@ catalog[0].measure[Nmeas] = stars[i+j].measure; - catalog[0].measure[Nmeas].dR = 0.0; - catalog[0].measure[Nmeas].dD = 0.0; + catalog[0].measure[Nmeas].R = R; + catalog[0].measure[Nmeas].D = D; catalog[0].measure[Nmeas].dt = NAN_S_SHORT; Index: /trunk/Ohana/src/addstar/src/loadsupercos_ops.c =================================================================== --- /trunk/Ohana/src/addstar/src/loadsupercos_ops.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/loadsupercos_ops.c (revision 37036) @@ -1,7 +1,4 @@ # include "addstar.h" # include "supercos.h" - -// this is defined in libohana/src/string.c but not generally exposed -char *_parse_nextword_csv (char *string); int loadsupercos_getFilterInfo (char *line, char *emulsion, char *filterID) { @@ -13,5 +10,5 @@ // emulsion is field 11, filterID is field 12 for (i = 1; i < 11; i++) { - p1 = _parse_nextword_csv (p1); + p1 = parse_nextword_csv (p1); if (!p1) { fprintf (stderr, "error parsing filter info for line %s\n", line); @@ -24,5 +21,5 @@ } - char *p2 = _parse_nextword_csv (p1); + char *p2 = parse_nextword_csv (p1); if (!p2) { fprintf (stderr, "error parsing filter info for line %s\n", line); @@ -34,5 +31,5 @@ } - char *p3 = _parse_nextword_csv (p2); + char *p3 = parse_nextword_csv (p2); if (!p3) { fprintf (stderr, "error parsing filter info for line %s\n", line); @@ -60,5 +57,5 @@ // lstObs is field 20 for (i = 1; i < 20; i++) { - p1 = _parse_nextword_csv (p1); + p1 = parse_nextword_csv (p1); if (!p1) { fprintf (stderr, "error parsing filter info for line %s\n", line); Index: /trunk/Ohana/src/addstar/src/mkcmf.c =================================================================== --- /trunk/Ohana/src/addstar/src/mkcmf.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/mkcmf.c (revision 37036) @@ -1,4 +1,5 @@ # include "mkcmf.h" +# define ZERO_POINT 25.0 # define SKY 100.0 # define DSKY 2.0 @@ -14,8 +15,11 @@ void gauss_init (int Nbin); double rnd_gauss (double mean, double sigma); +void writeStars_PS1_V5_Lensing (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars); +void writeStars_PS1_V5 (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars); void writeStars_PS1_V4 (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars); void writeStars_PS1_V3 (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars); void writeStars_PS1_V2 (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars); void writeStars_PS1_V1 (FTable *ftable, double *X, double *Y, double *M, int Nstars); +void writeStars_PS1_SV3 (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars); void writeStars_PS1_DEV_1 (FTable *ftable, double *X, double *Y, double *M, int Nstars); void writeStars_PS1_DEV_0 (FTable *ftable, double *X, double *Y, double *M, int Nstars); @@ -23,4 +27,7 @@ int ADDNOISE = TRUE; float BAD_PSFQF_FRAC = 0.0; + +static float exptime = 1.0; +static Coords coords; int main (int argc, char **argv) { @@ -38,5 +45,4 @@ Matrix matrix; FTable ftable; - Coords coords; int APPEND = FALSE; @@ -105,5 +111,4 @@ } - float exptime = 1.0; if ((N = get_argument (argc, argv, "-exptime"))) { remove_argument (N, &argc, argv); @@ -288,5 +293,5 @@ gfits_modify (&header, "UTC-OBS", "%s", 1, time); } - gfits_modify (&header, "ZERO_PT", "%lf", 1, 25.0); + gfits_modify (&header, "ZERO_PT", "%lf", 1, ZERO_POINT); gfits_modify (&header, "EXPTIME", "%lf", 1, exptime); gfits_modify (&header, "AIRMASS", "%lf", 1, airmass); @@ -327,4 +332,16 @@ if (!strcmp(type, "PS1_V4")) { writeStars_PS1_V4 (&ftable, X, Y, M, Flag, Nstars); + found = TRUE; + } + if (!strcmp(type, "PS1_V5")) { + writeStars_PS1_V5 (&ftable, X, Y, M, Flag, Nstars); + found = TRUE; + } + if (!strcmp(type, "PS1_V5_Lensing")) { + writeStars_PS1_V5_Lensing (&ftable, X, Y, M, Flag, Nstars); + found = TRUE; + } + if (!strcmp(type, "PS1_SV3")) { + writeStars_PS1_SV3 (&ftable, X, Y, M, Flag, Nstars); found = TRUE; } @@ -764,2 +781,353 @@ } +void writeStars_PS1_V5 (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars) { + + int i; + CMF_PS1_V5 *stars; + float flux, fSN; + + // XXX add gaussian-distributed noise based on counts + // this needs to make different output 'stars' entries depending on the desired type + ALLOCATE (stars, CMF_PS1_V5, Nstars); + for (i = 0; i < Nstars; i++) { + + flux = pow (10.0, -0.4*M[i]); + fSN = 1.0 / sqrt(flux); + + stars[i].detID = i; + stars[i].X = X[i]; + stars[i].Y = Y[i]; + stars[i].dX = FX * fSN; + stars[i].dY = FY * fSN; + + stars[i].posangle = 10.0; + stars[i].pltscale = 0.25; + + stars[i].M = M[i]; + stars[i].dM = fSN; + + stars[i].Flux = flux; + stars[i].dFlux = flux * fSN; + + stars[i].Map = M[i] - 0.05; + stars[i].MapRaw = M[i] - 0.10; + + stars[i].apRadius = 8.0; + + stars[i].apFlux = pow(10.0, -0.4*stars[i].Map); + stars[i].apFluxErr = stars[i].apFlux * fSN; + + stars[i].Mcalib = M[i] + ZERO_POINT + 2.5*log10(exptime); + stars[i].dMcal = 0.05; + + XY_to_RD (&stars[i].RA, &stars[i].DEC, X[i], Y[i], &coords); + stars[i].apNpix = 3.14*8.0*8.0; + + stars[i].Mpeak = M[i] + 1.0; + stars[i].sky = SKY; + stars[i].dSky = DSKY; + stars[i].psfChisq = PSFCHI; + stars[i].crNsigma = CRN; + stars[i].extNsigma = EXTN; + stars[i].fx = FX; + stars[i].fy = FY; + stars[i].df = DF; + + stars[i].k = 1.0; + stars[i].fwhmMaj = FX*2.8; + stars[i].fwhmMin = FY*2.8; + + // randomly give poor PSFQF values + if ((BAD_PSFQF_FRAC > 0.0) && (drand48() < BAD_PSFQF_FRAC)) { + stars[i].psfQF = 0.25; + stars[i].psfQFperf = 0.24; + } else { + stars[i].psfQF = PSFQUAL; + stars[i].psfQFperf = MAX(PSFQUAL - 0.01, 0.0); + } + + stars[i].psfNdof = 1; + stars[i].psfNpix = 2; + + stars[i].Mxx = FX; + stars[i].Mxy = 0.01; + stars[i].Myy = FX; + stars[i].M3c = FX; + stars[i].M3s = FX; + stars[i].M4c = FX; + stars[i].M4s = FX; + stars[i].Mr1 = FX; + stars[i].Mrh = FX; + + stars[i].kronFlux = flux * 1.25; + stars[i].kronFluxErr = fSN * flux * 1.25; + + stars[i].kronInner = fSN * flux * 0.9; + stars[i].kronOuter = fSN * flux * 1.5; + + stars[i].skyLimitRad = 1; + stars[i].skyLimitFlux = 2; + stars[i].skyLimitSlope = 0.2; + + stars[i].flags = Flag[i]; + stars[i].flags2 = 0x80; + + stars[i].nFrames = 1; + + if (ADDNOISE) { + stars[i].X += FX * fSN * rnd_gauss(0.0, 1.0); + stars[i].Y += FY * fSN * rnd_gauss(0.0, 1.0); + float Moff = fSN*rnd_gauss(0.0, 1.0); + stars[i].M += Moff; + stars[i].Map += Moff; + stars[i].MapRaw += Moff; + stars[i].Mcalib += Moff; + } + } + + gfits_table_set_CMF_PS1_V5 (ftable, stars, Nstars); + gfits_modify (ftable->header, "EXTTYPE", "%s", 1, "PS1_V5"); +} + +void writeStars_PS1_SV3 (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars) { + + int i; + CMF_PS1_SV3 *stars; + float flux, fSN; + + // XXX add gaussian-distributed noise based on counts + // this needs to make different output 'stars' entries depending on the desired type + ALLOCATE (stars, CMF_PS1_SV3, Nstars); + for (i = 0; i < Nstars; i++) { + + flux = pow (10.0, -0.4*M[i]); + fSN = 1.0 / sqrt(flux); + + stars[i].detID = i; + stars[i].X = X[i]; + stars[i].Y = Y[i]; + stars[i].dX = FX * fSN; + stars[i].dY = FY * fSN; + + stars[i].posangle = 10.0; + stars[i].pltscale = 0.25; + + stars[i].M = M[i]; + stars[i].dM = fSN; + + stars[i].Flux = flux; + stars[i].dFlux = flux * fSN; + + stars[i].Map = M[i] - 0.05; + stars[i].MapRaw = M[i] - 0.10; + + stars[i].apRadius = 8.0; + + stars[i].apFlux = pow(10.0, -0.4*stars[i].Map); + stars[i].apFluxErr = stars[i].apFlux * fSN; + + stars[i].Mcalib = M[i] + ZERO_POINT + 2.5*log10(exptime); + stars[i].dMcal = 0.05; + + XY_to_RD (&stars[i].RA, &stars[i].DEC, X[i], Y[i], &coords); + stars[i].apNpix = 3.14*8.0*8.0; + + stars[i].Mpeak = M[i] + 1.0; + stars[i].sky = SKY; + stars[i].dSky = DSKY; + stars[i].psfChisq = PSFCHI; + stars[i].crNsigma = CRN; + stars[i].extNsigma = EXTN; + stars[i].fx = FX; + stars[i].fy = FY; + stars[i].df = DF; + + stars[i].k = 1.0; + stars[i].fwhmMaj = FX*2.8; + stars[i].fwhmMin = FY*2.8; + + // randomly give poor PSFQF values + if ((BAD_PSFQF_FRAC > 0.0) && (drand48() < BAD_PSFQF_FRAC)) { + stars[i].psfQF = 0.25; + stars[i].psfQFperf = 0.24; + } else { + stars[i].psfQF = PSFQUAL; + stars[i].psfQFperf = MAX(PSFQUAL - 0.01, 0.0); + } + + stars[i].psfNdof = 1; + stars[i].psfNpix = 2; + + stars[i].Mxx = FX; + stars[i].Mxy = 0.01; + stars[i].Myy = FX; + stars[i].M3c = FX; + stars[i].M3s = FX; + stars[i].M4c = FX; + stars[i].M4s = FX; + stars[i].Mr1 = FX; + stars[i].Mrh = FX; + + stars[i].kronFlux = flux * 1.25; + stars[i].kronFluxErr = fSN * flux * 1.25; + + stars[i].kronInner = fSN * flux * 0.9; + stars[i].kronOuter = fSN * flux * 1.5; + + // stars[i].skyLimitRad = 1; + // stars[i].skyLimitFlux = 2; + // stars[i].skyLimitSlope = 0.2; + + stars[i].flags = Flag[i]; + stars[i].flags2 = 0x80; + + stars[i].nFrames = 1; + + if (ADDNOISE) { + stars[i].X += FX * fSN * rnd_gauss(0.0, 1.0); + stars[i].Y += FY * fSN * rnd_gauss(0.0, 1.0); + float Moff = fSN*rnd_gauss(0.0, 1.0); + stars[i].M += Moff; + stars[i].Map += Moff; + stars[i].MapRaw += Moff; + stars[i].Mcalib += Moff; + } + } + + gfits_table_set_CMF_PS1_SV3 (ftable, stars, Nstars); + gfits_modify (ftable->header, "EXTTYPE", "%s", 1, "PS1_SV3"); +} + +void writeStars_PS1_V5_Lensing (FTable *ftable, double *X, double *Y, double *M, unsigned int *Flag, int Nstars) { + + int i; + CMF_PS1_V5_Lensing *stars; + float flux, fSN; + + // XXX add gaussian-distributed noise based on counts + // this needs to make different output 'stars' entries depending on the desired type + ALLOCATE (stars, CMF_PS1_V5_Lensing, Nstars); + for (i = 0; i < Nstars; i++) { + + flux = pow (10.0, -0.4*M[i]); + fSN = 1.0 / sqrt(flux); + + stars[i].detID = i; + stars[i].X = X[i]; + stars[i].Y = Y[i]; + stars[i].dX = FX * fSN; + stars[i].dY = FY * fSN; + + stars[i].posangle = 10.0; + stars[i].pltscale = 0.25; + + stars[i].M = M[i]; + stars[i].dM = fSN; + + stars[i].Flux = flux; + stars[i].dFlux = flux * fSN; + + stars[i].Map = M[i] - 0.05; + stars[i].MapRaw = M[i] - 0.10; + + stars[i].apRadius = 8.0; + + stars[i].apFlux = pow(10.0, -0.4*stars[i].Map); + stars[i].apFluxErr = stars[i].apFlux * fSN; + + stars[i].Mcalib = M[i] + ZERO_POINT + 2.5*log10(exptime); + stars[i].dMcal = 0.05; + + XY_to_RD (&stars[i].RA, &stars[i].DEC, X[i], Y[i], &coords); + stars[i].apNpix = 3.14*8.0*8.0; + + stars[i].Mpeak = M[i] + 1.0; + stars[i].sky = SKY; + stars[i].dSky = DSKY; + stars[i].psfChisq = PSFCHI; + stars[i].crNsigma = CRN; + stars[i].extNsigma = EXTN; + stars[i].fx = FX; + stars[i].fy = FY; + stars[i].df = DF; + + stars[i].k = 1.0; + stars[i].fwhmMaj = FX*2.8; + stars[i].fwhmMin = FY*2.8; + + // randomly give poor PSFQF values + if ((BAD_PSFQF_FRAC > 0.0) && (drand48() < BAD_PSFQF_FRAC)) { + stars[i].psfQF = 0.25; + stars[i].psfQFperf = 0.24; + } else { + stars[i].psfQF = PSFQUAL; + stars[i].psfQFperf = MAX(PSFQUAL - 0.01, 0.0); + } + + stars[i].psfNdof = 1; + stars[i].psfNpix = 2; + + stars[i].Mxx = FX; + stars[i].Mxy = 0.01; + stars[i].Myy = FX; + stars[i].M3c = FX; + stars[i].M3s = FX; + stars[i].M4c = FX; + stars[i].M4s = FX; + stars[i].Mr1 = FX; + stars[i].Mrh = FX; + + stars[i].X11_sm_obj = FX; + stars[i].X12_sm_obj = FX; + stars[i].X22_sm_obj = FY; + stars[i].E1_sm_obj = FX; + stars[i].E2_sm_obj = FX; + + stars[i].X11_sh_obj = FX; + stars[i].X12_sh_obj = 0.2; + stars[i].X22_sh_obj = FY; + stars[i].E1_sh_obj = FX*0.9; + stars[i].E2_sh_obj = FX*0.8; + + stars[i].X11_sm_psf = FX; + stars[i].X12_sm_psf = FX; + stars[i].X22_sm_psf = FX; + stars[i].E1_sm_psf = FX; + stars[i].E2_sm_psf = FX; + + stars[i].X11_sh_psf = FX; + stars[i].X12_sh_psf = FX; + stars[i].X22_sh_psf = FX; + stars[i].E1_sh_psf = FX; + stars[i].E2_sh_psf = FX; + + stars[i].kronFlux = flux * 1.25; + stars[i].kronFluxErr = fSN * flux * 1.25; + + stars[i].kronInner = fSN * flux * 0.9; + stars[i].kronOuter = fSN * flux * 1.5; + + stars[i].skyLimitRad = 1; + stars[i].skyLimitFlux = 2; + stars[i].skyLimitSlope = 0.2; + + stars[i].flags = Flag[i]; + stars[i].flags2 = 0x80; + + stars[i].nFrames = 1; + + if (ADDNOISE) { + stars[i].X += FX * fSN * rnd_gauss(0.0, 1.0); + stars[i].Y += FY * fSN * rnd_gauss(0.0, 1.0); + float Moff = fSN*rnd_gauss(0.0, 1.0); + stars[i].M += Moff; + stars[i].Map += Moff; + stars[i].MapRaw += Moff; + stars[i].Mcalib += Moff; + } + } + + gfits_table_set_CMF_PS1_V5_Lensing (ftable, stars, Nstars); + gfits_modify (ftable->header, "EXTTYPE", "%s", 1, "PS1_V5"); +} + Index: /trunk/Ohana/src/addstar/src/replace_match.c =================================================================== --- /trunk/Ohana/src/addstar/src/replace_match.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/replace_match.c (revision 37036) @@ -12,6 +12,6 @@ j = i + m; if (measure[j].photcode != star[0].measure.photcode) continue; - measure[j].dR = 3600.0*(average[0].R - star[0].average.R); - measure[j].dD = 3600.0*(average[0].D - star[0].average.D); + measure[j].R = star[0].average.R; + measure[j].D = star[0].average.D; measure[j].M = star[0].measure.M; measure[j].dM = star[0].measure.dM; Index: /trunk/Ohana/src/addstar/src/resort_catalog.c =================================================================== --- /trunk/Ohana/src/addstar/src/resort_catalog.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/resort_catalog.c (revision 37036) @@ -1,3 +1,7 @@ # include "addstar.h" + +void SortAveMatch (off_t *MEAS, off_t *AVE, off_t N); +void resort_catalog_measure (Catalog *catalog); +void resort_catalog_lensing (Catalog *catalog); void resort_catalog_old (Catalog *catalog) { @@ -5,10 +9,8 @@ off_t *next_meas; off_t Naves, Nmeas; - double dtime; - struct timeval start, stop; if (catalog[0].sorted == TRUE) return; - gettimeofday (&start, NULL); + INITTIME; /* internal counters */ @@ -22,34 +24,18 @@ catalog[0].measure = sort_measure (catalog[0].average, Naves, catalog[0].measure, Nmeas, next_meas); - gettimeofday (&stop, NULL); - dtime = DTIME (stop, start); - fprintf (stderr, " match time %9.4f sec for %7lld measures, %6lld average\n", dtime, (long long) Nmeas, (long long) Naves); + MARKTIME (" match time %9.4f sec for %7lld measures, %6lld average\n", dtime, (long long) Nmeas, (long long) Naves); return; } -// sort the measure Sequence based on the average Sequence entries -void SortAveMeasMatch (off_t *MEAS, off_t *AVE, off_t N) { - -# define SWAPFUNC(A,B){ off_t tmp_meas; off_t tmp_ave; \ - tmp_meas = MEAS[A]; MEAS[A] = MEAS[B]; MEAS[B] = tmp_meas; \ - tmp_ave = AVE[A]; AVE[A] = AVE[B]; AVE[B] = tmp_ave; \ - } -# define COMPARE(A,B)(AVE[A] < AVE[B]) - OHANA_SORT (N, COMPARE, SWAPFUNC); -# undef SWAPFUNC -# undef COMPARE -} - -# define MARKTIME(MSG,...) { \ - float dtime; \ - gettimeofday (&stop, (void *) NULL); \ - dtime = DTIME (stop, start); start = stop; \ - fprintf (stderr, MSG, __VA_ARGS__); } - -// XXX : where is the time going? perhaps the ALLOCATE? -// XXX : I don't thnk his is getting the right answer yet. - void resort_catalog (Catalog *catalog) { + + if (catalog[0].sorted == TRUE) return; + + resort_catalog_measure (catalog); + resort_catalog_lensing (catalog); +} + +void resort_catalog_measure (Catalog *catalog) { off_t Naverage, Nmeasure; @@ -62,6 +48,4 @@ Measure *measureTMP = NULL; - if (catalog[0].sorted == TRUE) return; - // struct timeval start, stop; // gettimeofday (&start, NULL); @@ -70,4 +54,6 @@ Nmeasure = catalog[0].Nmeasure; Naverage = catalog[0].Naverage; + + if (!Nmeasure) return; measure = catalog[0].measure; @@ -113,5 +99,5 @@ // fprintf (stderr, "\n"); - SortAveMeasMatch(measureSeq, averageSeq, Nmeasure); + SortAveMatch(measureSeq, averageSeq, Nmeasure); // MARKTIME("sort : %f sec\n", dtime); @@ -188,2 +174,152 @@ } +void resort_catalog_lensing (Catalog *catalog) { + + off_t Naverage, Nlensing; + Lensing *lensing; + Average *average; + off_t i, j, N, currentAve; + + off_t *lensingSeq = NULL; + off_t *averageSeq = NULL; + Lensing *lensingTMP = NULL; + + // struct timeval start, stop; + // gettimeofday (&start, NULL); + + /* internal counters */ + Nlensing = catalog[0].Nlensing; + Naverage = catalog[0].Naverage; + + if (!Nlensing) return; + + lensing = catalog[0].lensing; + average = catalog[0].average; + + // we have a table of average objects and an unsorted table of measurements. each measurement + // has a reference to the average object sequence (as well as an ID) + // lensing[i].averef -> average[averef] + // lensing[i].objID = average[averef].objID + // lensing[i].catID = average[averef].catID + + // we want a sorted lensing array with all averef entries in sequence + + ALLOCATE (lensingSeq, off_t, Nlensing); + ALLOCATE (averageSeq, off_t, Nlensing); + + for (i = 0; i < Nlensing; i++) { + lensingSeq[i] = i; + averageSeq[i] = lensing[i].averef; + + 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]].catID == lensing[lensingSeq[i]].catID, "object / detection mismatch"); +# if (1) + myAssert(average[averageSeq[i]].objID == lensing[lensingSeq[i]].objID, "object / detection mismatch"); +# else + // for reasons I do not understand, the mini dvodbs generated on stsci1X had a handful of detections with an inconsistency between averef and objID. + // this happened for 28 detections in the dbs on /data/stsci1?.0/eugene/dvo3pi.20130616, but not at all (as far as I know) in the rest of LAP DVO + if (average[averageSeq[i]].objID != lensing[lensingSeq[i]].objID) { + fprintf (stderr, "R"); + lensing[lensingSeq[i]].objID = average[averageSeq[i]].objID; // XXX I don't really like this... + } +# endif + } + } + + // check that averageSeq is now in order + // for (i = 1; i < Nlensing; i++) { + // if (averageSeq[i] < averageSeq[i-1]) { + // fprintf (stderr, "%d ", (int) i); + // } + // } + // fprintf (stderr, "\n"); + + SortAveMatch(lensingSeq, averageSeq, Nlensing); + // MARKTIME("sort : %f sec\n", dtime); + + // check that averageSeq is now in order + // for (i = 1; i < Nlensing; i++) { + // if (averageSeq[i] < averageSeq[i-1]) { + // fprintf (stderr, "%d ", (int) i); + // } + // } + // fprintf (stderr, "\n"); + + // copy the lensing entries in the sorted order + ALLOCATE (lensingTMP, Lensing, Nlensing); + for (i = 0; i < Nlensing; i++) { + j = lensingSeq[i]; + lensingTMP[i] = lensing[j]; + } + // MARKTIME("assign lensing : %f sec\n", dtime); + + // update the values of average.lensingOffset and average.Nlensing + FREE(lensing); + catalog[0].lensing = lensingTMP; + + N = 0; + currentAve = averageSeq[0]; + average[currentAve].lensingOffset = 0; + for (i = 0; i < Nlensing; i++) { + if (averageSeq[i] != currentAve) { + // we have hit the next entry in the list + average[currentAve].Nlensing = N; + N = 0; + currentAve = averageSeq[i]; + average[currentAve].lensingOffset = i; + } + N++; + } + // N++; + average[currentAve].Nlensing = N; + // MARKTIME("update Nlensing : %f sec\n", dtime); + + int NlensingTotal = 0; + int lensingOffsetOK = TRUE; + for (i = 0; i < Naverage; i++) { + NlensingTotal += catalog[0].average[i].Nlensing; + if (VERBOSE && !(NlensingTotal <= catalog[0].Nlensing)) { + fprintf (stderr, "too few lensing: %d %d %d\n", (int) i, NlensingTotal, (int) catalog[0].Nlensing); + } + lensingOffsetOK &= (catalog[0].average[i].lensingOffset < catalog[0].Nlensing); + if (VERBOSE && !(catalog[0].average[i].lensingOffset < catalog[0].Nlensing)) { + fprintf (stderr, "offset too large: %d %d %d\n", (int) i, catalog[0].average[i].Nlensing, (int) catalog[0].Nlensing); + } + lensingOffsetOK &= (catalog[0].average[i].lensingOffset + catalog[0].average[i].Nlensing <= catalog[0].Nlensing); + if (VERBOSE && !(catalog[0].average[i].lensingOffset + catalog[0].average[i].Nlensing <= catalog[0].Nlensing)) { + fprintf (stderr, "orrset + Nlensing too large: %d + %d > %d %d\n", (int) i, catalog[0].average[i].lensingOffset, catalog[0].average[i].Nlensing, (int) catalog[0].Nlensing); + } + } + + if (!lensingOffsetOK) { + fprintf (stderr, "ERROR: catalog %s has an invalid lensingOffset\n", catalog[0].filename); + } + + if (NlensingTotal != catalog[0].Nlensing) { + fprintf (stderr, "ERROR: catalog %s has an invalid Nlensing\n", catalog[0].filename); + } + + // MARKTIME(" match time %9.4f sec for %7lld lensings, %6lld average\n", dtime, (long long) Nlensing, (long long) Naverage); + + catalog[0].sorted = TRUE; + + FREE (lensingSeq); + FREE (averageSeq); + + return; +} + +// sort the measure or lensing Sequence based on the average Sequence entries +void SortAveMatch (off_t *MEAS, off_t *AVE, off_t N) { + +# define SWAPFUNC(A,B){ off_t tmp_meas; off_t tmp_ave; \ + tmp_meas = MEAS[A]; MEAS[A] = MEAS[B]; MEAS[B] = tmp_meas; \ + tmp_ave = AVE[A]; AVE[A] = AVE[B]; AVE[B] = tmp_ave; \ + } +# define COMPARE(A,B)(AVE[A] < AVE[B]) + OHANA_SORT (N, COMPARE, SWAPFUNC); +# undef SWAPFUNC +# undef COMPARE +} + Index: /trunk/Ohana/src/addstar/src/resort_threaded.c =================================================================== --- /trunk/Ohana/src/addstar/src/resort_threaded.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/resort_threaded.c (revision 37036) @@ -41,5 +41,5 @@ catalog.catformat = dvo_catalog_catformat (CATFORMAT); // set the default catformat from config data catalog.catmode = dvo_catalog_catmode (CATMODE); // set the default catmode from config data - catalog.catflags = LOAD_AVES | LOAD_MEAS; + catalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_LENSING; catalog.Nsecfilt = GetPhotcodeNsecfilt (); Index: /trunk/Ohana/src/addstar/src/resort_unthreaded_catalogs.c =================================================================== --- /trunk/Ohana/src/addstar/src/resort_unthreaded_catalogs.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/resort_unthreaded_catalogs.c (revision 37036) @@ -30,5 +30,5 @@ catalog.catformat = dvo_catalog_catformat (CATFORMAT); // set the default catformat from config data catalog.catmode = dvo_catalog_catmode (CATMODE); // set the default catmode from config data - catalog.catflags = LOAD_AVES | LOAD_MEAS; + catalog.catflags = LOAD_AVES | LOAD_MEAS | LOAD_LENSING; catalog.Nsecfilt = GetPhotcodeNsecfilt (); Index: /trunk/Ohana/src/addstar/src/update_coords.c =================================================================== --- /trunk/Ohana/src/addstar/src/update_coords.c (revision 37035) +++ /trunk/Ohana/src/addstar/src/update_coords.c (revision 37036) @@ -19,6 +19,6 @@ continue; } - R = measure[m].dR; - D = measure[m].dD; + R = measure[m].R; + D = measure[m].D; r += R; d += D; @@ -33,13 +33,8 @@ r = r / Npt; /* these are corrections in 1/100 arcsec to RA and DEC */ d = d / Npt; - average[0].R -= r / 3600.0; - average[0].D -= d / 3600.0; + average[0].R = r; + average[0].D = d; m = average[0].measureOffset; /* first measurement of this star */ - for (i = 0; i < average[0].Nmeasure; i++) { - measure[m].dR -= r; - measure[m].dD -= d; - m = next[m]; - } - + /* measure scatter, if possible */ if (Npt < 2) return; @@ -47,5 +42,5 @@ dR2 = r2 / Npt - r*r; dD2 = d2 / Npt - d*d; - average[0].ChiSqAve = sqrt (dD2 + dR2 / SQ(cos(d*RAD_DEG))); + average[0].ChiSqAve = 3600.0*sqrt (dD2 + dR2 / SQ(cos(d*RAD_DEG))); /* ChiSqAve is supposed to be a chisq */ Index: /trunk/Ohana/src/addstar/test/simple.dvo =================================================================== --- /trunk/Ohana/src/addstar/test/simple.dvo (revision 37035) +++ /trunk/Ohana/src/addstar/test/simple.dvo (revision 37036) @@ -7,30 +7,34 @@ test.fields PS1_DEV_0 PS1_V1 test.fields PS1_DEV_1 PS1_V1 - test.fields PS1_V1 PS1_V1 - test.fields PS1_V2 PS1_V1 - test.fields PS1_V3 PS1_V1 - test.fields PS1_V4 PS1_V1 + test.fields PS1_V1 PS1_V1 + test.fields PS1_V2 PS1_V1 + test.fields PS1_V3 PS1_V1 + test.fields PS1_V4 PS1_V1 test.fields PS1_DEV_0 PS1_V2 test.fields PS1_DEV_1 PS1_V2 - test.fields PS1_V1 PS1_V2 - test.fields PS1_V2 PS1_V2 - test.fields PS1_V3 PS1_V2 - test.fields PS1_V4 PS1_V4 + test.fields PS1_V1 PS1_V2 + test.fields PS1_V2 PS1_V2 + test.fields PS1_V3 PS1_V2 + test.fields PS1_V4 PS1_V4 test.fields PS1_DEV_0 PS1_V3 test.fields PS1_DEV_1 PS1_V3 - test.fields PS1_V1 PS1_V3 - test.fields PS1_V2 PS1_V3 - test.fields PS1_V3 PS1_V3 - test.fields PS1_V4 PS1_V4 + test.fields PS1_V1 PS1_V3 + test.fields PS1_V2 PS1_V3 + test.fields PS1_V3 PS1_V3 + test.fields PS1_V4 PS1_V4 test.fields PS1_DEV_0 PS1_V4 test.fields PS1_DEV_1 PS1_V4 - test.fields PS1_V1 PS1_V4 - test.fields PS1_V2 PS1_V4 - test.fields PS1_V3 PS1_V4 - test.fields PS1_V4 PS1_V4 + test.fields PS1_V1 PS1_V4 + test.fields PS1_V2 PS1_V4 + test.fields PS1_V3 PS1_V4 + test.fields PS1_V4 PS1_V4 end + +if (not($?NO_NOISE)) + $NO_NOISE = "" +end # create a populated catdir with a single cmf -- test each field @@ -50,9 +54,11 @@ mkinput - exec mkcmf test.in.txt test.cmf -date 2008/1/1 -time 01:00:00 -radec $RA $DEC -type $1 + echo mkcmf test.in.txt test.cmf -date 2008/1/1 -time 01:00:00 -radec $RA $DEC -type $1 $NO_NOISE + exec mkcmf test.in.txt test.cmf -date 2008/1/1 -time 01:00:00 -radec $RA $DEC -type $1 $NO_NOISE if ($TAP_VERBOSE) echo exec addstar -D CATDIR catdir.test -D CAMERA simtest test.cmf -D CATFORMAT $2 -quick-airmass exec addstar -D CATDIR catdir.test -D CAMERA simtest test.cmf -D CATFORMAT $2 -quick-airmass else + echo addstar -D CATDIR catdir.test -D CAMERA simtest test.cmf -D CATFORMAT $2 -quick-airmass exec addstar -D CATDIR catdir.test -D CAMERA simtest test.cmf -D CATFORMAT $2 -quick-airmass >& tmp.log end @@ -112,9 +118,17 @@ # some fields require arithmetic manipulations - if ("$name:0" == "KRON_FLUX") - set v1 = -2.5*log(v1) - end - if ("$name:0" == "KRON_FLUX_ERR") - set v1 = KRON_FLUX_ERR / KRON_FLUX + ## if ("$name:0" == "KRON_FLUX") + ## set v1 = -2.5*log(v1) + ## end + ## if ("$name:0" == "KRON_FLUX_ERR") + ## set v1 = KRON_FLUX_ERR / KRON_FLUX + ## end + if (("$name:0" == "X_PSF_SIG") || ("$name:0" == "Y_PSF_SIG")) + set v1 = int(100*(v1)) + set v2 = int(100*(v2 + 0.0001)) + end + if ("$name:0" == "POSANGLE") + set v1 = int((0xffff/360.0)*(v1)) + set v2 = int((0xffff/360.0)*(v2 + 0.0028)) end @@ -206,7 +220,7 @@ # this list is good for PS1_V1 list testfields_PS1_V1 - IPP_IDET : detid - X_PSF : xccd - Y_PSF : yccd + IPP_IDET : detid + X_PSF : xccd + Y_PSF : yccd X_PSF_SIG : xccd:err # FAIL Y_PSF_SIG : yccd:err # FAIL @@ -215,6 +229,6 @@ POSANGLE : SKIP # astrometry is not calibrated in the cmf PLTSCALE : SKIP # astrometry is not calibrated in the cmf - PSF_INST_MAG : mag:inst - PSF_INST_MAG_SIG : mag:err + PSF_INST_MAG : mag:inst + PSF_INST_MAG_SIG : mag:err AP_MAG_STANDARD : mag:ap # FAIL AP_MAG_RADIUS : SKIP # no accessor @@ -222,13 +236,13 @@ CAL_PSF_MAG : SKIP # photometry is not calibrated in the cmf CAL_PSF_MAG_SIG : SKIP # photometry is not calibrated in the cmf - SKY : sky - SKY_SIG : sky_err - PSF_CHISQ : psf_chisq - CR_NSIGMA : cr_nsigma + SKY : sky + SKY_SIG : sky_err + PSF_CHISQ : psf_chisq + CR_NSIGMA : cr_nsigma EXT_NSIGMA : ext_nsigma - PSF_MAJOR : FWHM_MAJ - PSF_MINOR : FWHM_MIN + PSF_MAJOR : FWHM_MAJ + PSF_MINOR : FWHM_MIN PSF_THETA : THETA # FAIL - PSF_QF : PSF_QF + PSF_QF : PSF_QF PSF_NDOF : SKIP # no accessor PSF_NPIX : SKIP # no accessor @@ -237,12 +251,12 @@ MOMENTS_YY : SKIP # no accessor FLAGS : phot_flags - N_FRAMES : SKIP # no accessor + N_FRAMES : SKIP # no accessor end # this list is good for PS1_V2 list testfields_PS1_V2 - IPP_IDET : detid - X_PSF : xccd - Y_PSF : yccd + IPP_IDET : detid + X_PSF : xccd + Y_PSF : yccd X_PSF_SIG : xccd:err # FAIL Y_PSF_SIG : yccd:err # FAIL @@ -251,6 +265,6 @@ POSANGLE : SKIP # astrometry is not calibrated in the cmf PLTSCALE : SKIP # astrometry is not calibrated in the cmf - PSF_INST_MAG : mag:inst - PSF_INST_MAG_SIG : mag:err + PSF_INST_MAG : mag:inst + PSF_INST_MAG_SIG : mag:err AP_MAG_STANDARD : mag:ap # FAIL AP_MAG_RADIUS : SKIP # no accessor @@ -258,13 +272,13 @@ CAL_PSF_MAG : SKIP # photometry is not calibrated in the cmf CAL_PSF_MAG_SIG : SKIP # photometry is not calibrated in the cmf - SKY : sky - SKY_SIG : sky_err - PSF_CHISQ : psf_chisq - CR_NSIGMA : cr_nsigma + SKY : sky + SKY_SIG : sky_err + PSF_CHISQ : psf_chisq + CR_NSIGMA : cr_nsigma EXT_NSIGMA : ext_nsigma - PSF_MAJOR : FWHM_MAJ - PSF_MINOR : FWHM_MIN + PSF_MAJOR : FWHM_MAJ + PSF_MINOR : FWHM_MIN PSF_THETA : THETA # FAIL - PSF_QF : PSF_QF + PSF_QF : PSF_QF PSF_NDOF : SKIP # no accessor PSF_NPIX : SKIP # no accessor @@ -273,18 +287,18 @@ MOMENTS_YY : SKIP # no accessor FLAGS : phot_flags - N_FRAMES : SKIP # no accessor + N_FRAMES : SKIP # no accessor end # this list is good for PS1_V3 list testfields_PS1_V3 - IPP_IDET : detid - X_PSF : xccd - Y_PSF : yccd + IPP_IDET : detid + X_PSF : xccd + Y_PSF : yccd X_PSF_SIG : xccd:err # FAIL Y_PSF_SIG : yccd:err # FAIL POSANGLE : SKIP # astrometry is not calibrated in the cmf PLTSCALE : SKIP # astrometry is not calibrated in the cmf - PSF_INST_MAG : mag:inst - PSF_INST_MAG_SIG : mag:err + PSF_INST_MAG : mag:inst + PSF_INST_MAG_SIG : mag:err PSF_INST_FLUX : SKIP # not ingested into DVO PSF_INST_FLUX_SIG : SKIP # not ingested into DVO @@ -297,13 +311,13 @@ DEC_PSF : SKIP # astrometry is not calibrated in the cmf PEAK_FLUX_AS_MAG : SKIP # not ingested into DVO - SKY : sky - SKY_SIG : sky_err - PSF_CHISQ : psf_chisq - CR_NSIGMA : cr_nsigma + SKY : sky + SKY_SIG : sky_err + PSF_CHISQ : psf_chisq + CR_NSIGMA : cr_nsigma EXT_NSIGMA : ext_nsigma - PSF_MAJOR : FWHM_MAJ - PSF_MINOR : FWHM_MIN + PSF_MAJOR : FWHM_MAJ + PSF_MINOR : FWHM_MIN PSF_THETA : THETA # FAIL - PSF_QF : PSF_QF + PSF_QF : PSF_QF PSF_QF_PERFECT : SKIP # not ingested into DVO PSF_NDOF : PSF_NDOF @@ -326,51 +340,118 @@ end -# this list is good for PS1_V4 -list testfields_PS1_V4 - IPP_IDET : detid - X_PSF : xccd - Y_PSF : yccd - X_PSF_SIG : xccd:err # FAIL - Y_PSF_SIG : yccd:err # FAIL - POSANGLE : SKIP # astrometry is not calibrated in the cmf - PLTSCALE : SKIP # astrometry is not calibrated in the cmf - PSF_INST_MAG : mag:inst - PSF_INST_MAG_SIG : mag:err - PSF_INST_FLUX : SKIP # not ingested into DVO - PSF_INST_FLUX_SIG : SKIP # not ingested into DVO - AP_MAG : mag:aperinst # FAIL +# this list is good for PS1_V5 +list testfields_PS1_V5 + IPP_IDET : detid + X_PSF : xccd + Y_PSF : yccd + X_PSF_SIG : xccd:err + Y_PSF_SIG : yccd:err + POSANGLE : posangle + PLTSCALE : platescale + PSF_INST_MAG : mag:inst + PSF_INST_MAG_SIG : mag:err + PSF_INST_FLUX : flux:inst + PSF_INST_FLUX_SIG : flux:inst:err + AP_MAG : mag:aper:inst AP_MAG_RAW : SKIP # not ingested into DVO AP_MAG_RADIUS : SKIP # not ingested into DVO - CAL_PSF_MAG : SKIP # photometry is not calibrated in the cmf - CAL_PSF_MAG_SIG : SKIP # photometry is not calibrated in the cmf - RA_PSF : SKIP # astrometry is not calibrated in the cmf - DEC_PSF : SKIP # astrometry is not calibrated in the cmf - PEAK_FLUX_AS_MAG : SKIP # not ingested into DVO - SKY : sky - SKY_SIGMA : sky_err - PSF_CHISQ : psf_chisq - CR_NSIGMA : cr_nsigma - EXT_NSIGMA : ext_nsigma - PSF_MAJOR : FWHM_MAJ - PSF_MINOR : FWHM_MIN - PSF_THETA : THETA # FAIL - PSF_QF : PSF_QF - PSF_QF_PERFECT : SKIP # not ingested into DVO - PSF_NDOF : PSF_NDOF - PSF_NPIX : PSF_NPIX - MOMENTS_XX : MXX - MOMENTS_XY : MXY - MOMENTS_YY : MYY - MOMENTS_M3C : SKIP # not ingested into DVO - MOMENTS_M3S : SKIP # not ingested into DVO - MOMENTS_M4C : SKIP # not ingested into DVO - MOMENTS_M4S : SKIP # not ingested into DVO - MOMENTS_R1 : SKIP # not ingested into DVO - MOMENTS_RH : SKIP # not ingested into DVO - KRON_FLUX : mag:kroninst - KRON_FLUX_ERR : mag:kronerr - KRON_FLUX_INNER : SKIP # not ingested into DVO - KRON_FLUX_OUTER : SKIP # not ingested into DVO - FLAGS : phot_flags - N_FRAMES : SKIP # not ingested into DVO -end + AP_FLUX : flux:aper:inst + AP_FLUX_SIG : flux:aper:inst:err + CAL_PSF_MAG : mag:cal + CAL_PSF_MAG_SIG : mag:cal:err + RA_PSF : ra + DEC_PSF : dec + AP_NPIX : SKIP # not ingested into DVO + PEAK_FLUX_AS_MAG : SKIP # not ingested into DVO + SKY : sky + SKY_SIGMA : sky_err + PSF_CHISQ : psf_chisq + CR_NSIGMA : cr_nsigma + EXT_NSIGMA : ext_nsigma + PSF_MAJOR : fwhm_maj # XXX technically, this is not FWHM but Sigma + PSF_MINOR : fwhm_min # XXX technically, this is not FWHM but Sigma + PSF_THETA : theta # FAIL + PSF_CORE : SKIP # not ingested into DVO + PSF_FWHM_MAJ : SKIP # not ingested into DVO + PSF_FWHM_MIN : SKIP # not ingested into DVO + PSF_QF : psf_qf + PSF_QF_PERFECT : SKIP # not ingested into DVO + PSF_NDOF : psf_ndof + PSF_NPIX : psf_npix + MOMENTS_XX : Mxx + MOMENTS_XY : Mxy + MOMENTS_YY : Myy + MOMENTS_M3C : SKIP # not ingested into DVO + MOMENTS_M3S : SKIP # not ingested into DVO + MOMENTS_M4C : SKIP # not ingested into DVO + MOMENTS_M4S : SKIP # not ingested into DVO + MOMENTS_R1 : SKIP # not ingested into DVO + MOMENTS_RH : SKIP # not ingested into DVO + KRON_FLUX : flux:kron:inst + KRON_FLUX_ERR : flux:kron:inst:err + KRON_FLUX_INNER : SKIP # not ingested into DVO + KRON_FLUX_OUTER : SKIP # not ingested into DVO + SKY_LIMIT_RAD : SKIP # not ingested into DVO + SKY_LIMIT_FLUX : SKIP # not ingested into DVO + SKY_LIMIT_SLOPE : SKIP # not ingested into DVO + FLAGS : phot_flags + FLAGS2 : SKIP # not ingested into DVO + N_FRAMES : SKIP # not ingested into DVO +end + +# this list is good for PS1_SV3 +list testfields_PS1_SV3 + IPP_IDET : detid + X_PSF : xccd + Y_PSF : yccd + X_PSF_SIG : xccd:err + Y_PSF_SIG : yccd:err + POSANGLE : posangle + PLTSCALE : platescale + PSF_INST_MAG : mag:inst + PSF_INST_MAG_SIG : mag:err + PSF_INST_FLUX : flux:inst + PSF_INST_FLUX_SIG : flux:inst:err + AP_MAG : mag:aper:inst + AP_MAG_RAW : SKIP # not ingested into DVO + AP_MAG_RADIUS : SKIP # not ingested into DVO + AP_FLUX : flux:aper:inst + AP_FLUX_SIG : flux:aper:inst:err + CAL_PSF_MAG : mag:cal + CAL_PSF_MAG_SIG : mag:cal:err + RA_PSF : ra + DEC_PSF : dec + AP_NPIX : SKIP # not ingested into DVO + PEAK_FLUX_AS_MAG : SKIP # not ingested into DVO + SKY : sky + SKY_SIGMA : sky_err + PSF_CHISQ : psf_chisq + CR_NSIGMA : cr_nsigma + EXT_NSIGMA : ext_nsigma + PSF_MAJOR : fwhm_maj # XXX technically, this is not FWHM but Sigma + PSF_MINOR : fwhm_min # XXX technically, this is not FWHM but Sigma + PSF_THETA : theta # FAIL + PSF_CORE : SKIP # not ingested into DVO + PSF_FWHM_MAJ : SKIP # not ingested into DVO + PSF_FWHM_MIN : SKIP # not ingested into DVO + PSF_QF : psf_qf + PSF_QF_PERFECT : SKIP # not ingested into DVO + PSF_NDOF : psf_ndof + PSF_NPIX : psf_npix + MOMENTS_XX : Mxx + MOMENTS_XY : Mxy + MOMENTS_YY : Myy + MOMENTS_M3C : SKIP # not ingested into DVO + MOMENTS_M3S : SKIP # not ingested into DVO + MOMENTS_M4C : SKIP # not ingested into DVO + MOMENTS_M4S : SKIP # not ingested into DVO + MOMENTS_R1 : SKIP # not ingested into DVO + MOMENTS_RH : SKIP # not ingested into DVO + KRON_FLUX : flux:kron:inst + KRON_FLUX_ERR : flux:kron:inst:err + KRON_FLUX_INNER : SKIP # not ingested into DVO + KRON_FLUX_OUTER : SKIP # not ingested into DVO + FLAGS : phot_flags + FLAGS2 : SKIP # not ingested into DVO + N_FRAMES : SKIP # not ingested into DVO +end