Index: trunk/Ohana/src/relphot/include/relphot.h =================================================================== --- trunk/Ohana/src/relphot/include/relphot.h (revision 42807) +++ trunk/Ohana/src/relphot/include/relphot.h (revision 42821) @@ -417,6 +417,11 @@ int IS_DIFF_DB; +double MAG_MIN; double MAG_LIM; double SIGMA_LIM; + +double STAR_MAX_PSF_KRON; +double PSF_QF_MIN; + double IMAGE_SCATTER; double IMAGE_OFFSET; @@ -471,4 +476,5 @@ int RESET_FLATCORR; int REPAIR_WARPS; +int USE_REF_EQUIV; int PRESERVE_PS1; int REQUIRE_PSFFIT; @@ -486,4 +492,6 @@ double IMAGE_GOOD_FRACTION; int CALIBRATE_STACKS_AND_WARPS; +int CALIBRATE_WARPS_LIKE_CHIPS; +int CALIBRATE_STACKS_LIKE_CHIPS; int KEEP_UBERCAL; Index: trunk/Ohana/src/relphot/src/ConfigInit.c =================================================================== --- trunk/Ohana/src/relphot/src/ConfigInit.c (revision 42807) +++ trunk/Ohana/src/relphot/src/ConfigInit.c (revision 42821) @@ -23,6 +23,12 @@ if (VERBOSE) fprintf (stderr, "loaded config file: %s\n", file); - GetConfig (config, "MAG_LIM", "%lf", 0, &MAG_LIM); - GetConfig (config, "SIGMA_LIM", "%lf", 0, &SIGMA_LIM); + // the following cuts are limit the objects used to calculate the image zero points. + // These are used in bcatalog.c and apply to Mcat (nominal measurement magnitude) + GetConfig (config, "SIGMA_LIM", "%lf", 0, &SIGMA_LIM); + GetConfig (config, "MAG_LIM", "%lf", 0, &MAG_LIM); + + DefConfig ("MAG_MIN", "%lf", 15.0, MAG_MIN); + DefConfig ("STAR_MAX_PSF_KRON", "%lf", 0.2, STAR_MAX_PSF_KRON); + DefConfig ("PSF_QF_MIN", "%lf", 0.95, PSF_QF_MIN); DefConfig ("STAR_SCATTER", "%lf", 0.1, STAR_SCATTER); Index: trunk/Ohana/src/relphot/src/ImageOps.c =================================================================== --- trunk/Ohana/src/relphot/src/ImageOps.c (revision 42807) +++ trunk/Ohana/src/relphot/src/ImageOps.c (revision 42821) @@ -531,6 +531,9 @@ // until the analysis has converged a bit, do not use the IRLS analysis - // default is MaxIterations = 10 - if (UseStandardOLS(ZPT_IMAGES)) { + // default is MaxIterations = 10. + + // However, for stack and warp calibration, we are not iteratively solving for the + // zero points. In this case, use IRLS in a single pass + if (!CALIBRATE_STACKS_AND_WARPS && UseStandardOLS(ZPT_IMAGES)) { brightStars.MaxIterations = 0; kronStars.MaxIterations = 0; @@ -672,5 +675,5 @@ if ((image[i].imageID == TEST_IMAGE1) || (image[i].imageID == TEST_IMAGE2)) { - fprintf (stderr, "%1d, %3d : %3d, %3d : %10.6f %10.6f : %6.3f %6.3f %6.3f %6.3f %6.3f : %6.3f\n", (int) i, (int) j, (int) c, (int) m, catalog[c].averageT[n].R, catalog[c].averageT[n].D, MsysKron, MrelKron, Mmos, Mgrid, Mflat, kronStars.alldata->yVector[Nkron]); + fprintf (stderr, "%1d, %3d : %3d, %3d : %10.6f %10.6f : %6.3f %6.3f | %6.3f %6.3f | %6.3f %6.3f %6.3f %6.3f | %5.3f\n", (int) i, (int) j, (int) c, (int) m, catalog[c].averageT[n].R, catalog[c].averageT[n].D, MsysPSF, MrelPSF, MsysKron, MrelKron, Mmos, Mgrid, Mflat, Moff, psfStars.alldata->dyVector[Nref]); } @@ -702,17 +705,5 @@ if (mark) continue; - // use liststats to find the 20-pct, median, 80-pct points - StatType stats; - liststats_setmode (&stats, "MEDIAN"); - liststats (psfStars.alldata->yVector, NULL, NULL, Nref, &stats); - double altSigma = (stats.Upper80 - stats.Lower20) / 1.6; // 20% to 80% encompasses 60% of the values, corresponds to the range (-0.85 sigma : +0.85 sigma) - - // soften the individual errors with 10% of the scatter above - for (j = 0; j < Nref; j++) { - double newSigma = hypot(psfStars.alldata->dyVector[j], 0.1*altSigma); - psfStars.alldata->dyVector[j] = newSigma; - } - - // soften the errors based on the scatter + // soften the errors based on the core scatter FitDataSetSoften (&psfStars, Nref); Index: trunk/Ohana/src/relphot/src/StarOps.c =================================================================== --- trunk/Ohana/src/relphot/src/StarOps.c (revision 42807) +++ trunk/Ohana/src/relphot/src/StarOps.c (revision 42821) @@ -34,5 +34,8 @@ } -// +// Nsecfilt is a property of the photcode table. we only need to retrieve it once. +static int Nsecfilt = -1; + +// return the average photometry for this measurement float getMrel (Catalog *catalog, off_t meas, int cat, dvoMagClassType class, dvoMagSourceType source) { @@ -41,6 +44,16 @@ int ecode = GetPhotcodeEquivCodebyCode (photcode); + if (USE_REF_EQUIV) { + // tie this photometry to an alternative refernce + // (e.g., tie UDR4.i to PV3.i) + int refcode = GetPhotcodeEquivCodebyCode (ecode); + if (refcode < 0) return NAN; + ecode = refcode; + } + int Nsec = GetPhotcodeNsec(ecode); - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsec < 0) return NAN; + + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } int entry = Nsecfilt*ave+Nsec; @@ -312,5 +325,5 @@ } - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } SetMrelInfo summary, results; @@ -334,5 +347,5 @@ int i; - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } SetMrelInfo summary, results; @@ -424,5 +437,5 @@ ThreadInfo *threadinfo = data; - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } SetMrelInfo results; @@ -451,5 +464,5 @@ off_t k; - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } off_t m = average[0].measureOffset; @@ -559,5 +572,5 @@ off_t j; - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } for (i = 0; i < Ncatalog; i++) { @@ -601,5 +614,5 @@ ALLOCATE (slist, double, Ntot); - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } // eliminate bad stars using the stats for a single secfilt at a time @@ -666,5 +679,5 @@ // N1 = N2 = N3 = N4 = N0 = 0; - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } Ntot = 0; @@ -719,5 +732,5 @@ StatType stats; - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } Ntot = 0; @@ -760,5 +773,5 @@ StatType stats; - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } Ntot = 0; @@ -805,5 +818,5 @@ ALLOCATE (Mlist, double, NBIN); - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } int Ns; @@ -840,5 +853,5 @@ Graphdata graphdata; - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } Ntotal = 0; @@ -911,5 +924,5 @@ return; - int Nsecfilt = GetPhotcodeNsecfilt (); + if (Nsecfilt < 0) { Nsecfilt = GetPhotcodeNsecfilt (); } FILE *fout = NULL; Index: trunk/Ohana/src/relphot/src/args.c =================================================================== --- trunk/Ohana/src/relphot/src/args.c (revision 42807) +++ trunk/Ohana/src/relphot/src/args.c (revision 42821) @@ -237,4 +237,9 @@ } + USE_REF_EQUIV = FALSE; + if ((N = get_argument (argc, argv, "-use-ref-equiv"))) { + remove_argument (N, &argc, argv); + USE_REF_EQUIV = TRUE; + } PRESERVE_PS1 = FALSE; if ((N = get_argument (argc, argv, "-preserve-ps1"))) { @@ -521,4 +526,19 @@ remove_argument (N, &argc, argv); USE_REFERENCE_WEIGHT = TRUE; + } + + /* Normally, stacks and warps are treated differently from chips. In this case, + WHAT IS DIFFERENT? + */ + + CALIBRATE_WARPS_LIKE_CHIPS = FALSE; + if ((N = get_argument (argc, argv, "-calibrate-warps-like-chips"))) { + remove_argument (N, &argc, argv); + CALIBRATE_WARPS_LIKE_CHIPS = TRUE; + } + CALIBRATE_STACKS_LIKE_CHIPS = FALSE; + if ((N = get_argument (argc, argv, "-calibrate-stacks-like-chips"))) { + remove_argument (N, &argc, argv); + CALIBRATE_STACKS_LIKE_CHIPS = TRUE; } @@ -836,4 +856,9 @@ } + USE_REF_EQUIV = FALSE; + if ((N = get_argument (argc, argv, "-use-ref-equiv"))) { + remove_argument (N, &argc, argv); + USE_REF_EQUIV = TRUE; + } PRESERVE_PS1 = FALSE; if ((N = get_argument (argc, argv, "-preserve-ps1"))) { @@ -978,4 +1003,15 @@ } + CALIBRATE_WARPS_LIKE_CHIPS = FALSE; + if ((N = get_argument (argc, argv, "-calibrate-warps-like-chips"))) { + remove_argument (N, &argc, argv); + CALIBRATE_WARPS_LIKE_CHIPS = TRUE; + } + CALIBRATE_STACKS_LIKE_CHIPS = FALSE; + if ((N = get_argument (argc, argv, "-calibrate-stacks-like-chips"))) { + remove_argument (N, &argc, argv); + CALIBRATE_STACKS_LIKE_CHIPS = TRUE; + } + STAGES = STAGE_CHIP | STAGE_WARP | STAGE_STACK; if ((N = get_argument (argc, argv, "-skip-chip"))) { Index: trunk/Ohana/src/relphot/src/bcatalog.c =================================================================== --- trunk/Ohana/src/relphot/src/bcatalog.c (revision 42807) +++ trunk/Ohana/src/relphot/src/bcatalog.c (revision 42821) @@ -1,3 +1,23 @@ # include "relphot.h" + +# define BAD_COUNT(FIELD) \ + badCount.FIELD ++; \ + if (doTest1 && (catalog[0].measure[offset].imageID == TEST_IMAGE1)) { badTest1.FIELD ++; } \ + if (doTest2 && (catalog[0].measure[offset].imageID == TEST_IMAGE2)) { badTest2.FIELD ++; } + +typedef struct { + int Ncode; + int Ntime; + int Ndophot; + int Nmag; + int Nsigma; + int Nimag; + int Nfew; + int Ngalaxy; + int Npsfqf; + int Nnan; + int Nbad; + int Npoor; +} MyBadCountType; int LimitDensityCatalog_ByNmeasure (Catalog *subcatalog, Catalog *oldcatalog); @@ -9,5 +29,6 @@ off_t NAVERAGE, NMEASURE, Naverage, Nmeasure, Nm; float mag; - int Ncode, Ntime, Ndophot, Nmag, Nsigma, Nimag, Nfew, Ngalaxy, Npsfqf, Nnan, Nbad, Npoor; + + // int Ncode, Ntime, Ndophot, Nmag, Nsigma, Nimag, Nfew, Ngalaxy, Npsfqf, Nnan, Nbad, Npoor; int Nsecfilt = GetPhotcodeNsecfilt (); @@ -26,5 +47,15 @@ int NmMin = catalog[0].Nmeasure; - Ncode = Ntime = Ndophot = Nmag = Nsigma = Nimag = Nfew = Npsfqf = Ngalaxy = Nnan = Nbad = Npoor = 0; + MyBadCountType badCount = {0,0,0, 0,0,0, 0,0,0, 0,0,0}; + MyBadCountType badTest1 = {0,0,0, 0,0,0, 0,0,0, 0,0,0}; + MyBadCountType badTest2 = {0,0,0, 0,0,0, 0,0,0, 0,0,0}; + + int doTest1 = (TEST_IMAGE1 >= 0); + int doTest2 = (TEST_IMAGE2 >= 0); + + FILE *ftest1 = (doTest1) ? fopen ("test1.meas.dat", "a") : NULL; + FILE *ftest2 = (doTest2) ? fopen ("test2.meas.dat", "a") : NULL; + + // Ncode = Ntime = Ndophot = Nmag = Nsigma = Nimag = Nfew = Npsfqf = Ngalaxy = Nnan = Nbad = Npoor = 0; // copy the following fields to the subcatalog: @@ -57,4 +88,17 @@ offset = catalog[0].average[i].measureOffset + j; + + if (doTest1) { + Measure *m = &catalog[0].measure[offset]; + if (m->imageID == TEST_IMAGE1) { + fprintf (ftest1, "%5d %5d %5d | %5d %12d %5.3f %5.3f %5.3f %5.3f %4.3f %4.3f 0x%08x\n", m->detID, m->objID, m->catID, m->photcode, m->t, m->M, m->Map, m->Mkron, m->dM, m->psfQF, m->psfQFperf, m->photFlags); + } + } + if (doTest2) { + Measure *m = &catalog[0].measure[offset]; + if (m->imageID == TEST_IMAGE2) { + fprintf (ftest2, "%5d %5d %5d | %5d %12d %5.3f %5.3f %5.3f %5.3f %4.3f %4.3f 0x%08x\n", m->detID, m->objID, m->catID, m->photcode, m->t, m->M, m->Map, m->Mkron, m->dM, m->psfQF, m->psfQFperf, m->photFlags); + } + } /* select measurements by photcode */ @@ -69,5 +113,5 @@ } if (!found) { - Ncode ++; + BAD_COUNT(Ncode); continue; } @@ -75,34 +119,34 @@ /* select measurements by time */ if (TimeSelect) { - if (catalog[0].measure[offset].t < TSTART) { Ntime ++; continue; } - if (catalog[0].measure[offset].t > TSTOP) { Ntime ++; continue; } + if (catalog[0].measure[offset].t < TSTART) { BAD_COUNT(Ntime); continue; } + if (catalog[0].measure[offset].t > TSTOP) { BAD_COUNT(Ntime); continue; } } /* select measurements by quality */ - if (DophotSelect && ((catalog[0].measure[offset].photFlags >> 16) != DophotValue)) { Ndophot ++; continue; } + if (DophotSelect && ((catalog[0].measure[offset].photFlags >> 16) != DophotValue)) { BAD_COUNT(Ndophot); continue; } float Maper = USE_APER_FOR_STARGAL ? catalog[0].measure[offset].Map : catalog[0].measure[offset].Mkron; // skip garbage measurements - if (isnan(catalog[0].measure[offset].psfQF) || (catalog[0].measure[offset].psfQF < 0.95)) { Npsfqf ++; continue; } - if (isnan(catalog[0].measure[offset].psfQFperf) || (catalog[0].measure[offset].psfQFperf < 0.95)) { Npsfqf ++; continue; } - if (isnan(catalog[0].measure[offset].M)) { Nnan ++; continue; } - if (isnan(Maper)) { Nnan ++; continue; } + if (isnan(catalog[0].measure[offset].psfQF) || (catalog[0].measure[offset].psfQF < PSF_QF_MIN)) { BAD_COUNT(Npsfqf); continue; } + if (isnan(catalog[0].measure[offset].psfQFperf) || (catalog[0].measure[offset].psfQFperf < PSF_QF_MIN)) { BAD_COUNT(Npsfqf); continue; } + if (isnan(catalog[0].measure[offset].M)) { BAD_COUNT(Nnan); continue; } + if (isnan(Maper)) { BAD_COUNT(Nnan); continue; } // require 0x01 in photFlags (fitted with a PSF) -- add to the test data - if (REQUIRE_PSFFIT && (catalog[0].measure[offset].photFlags & 0x01) == 0) { Nbad ++; continue; } + if (REQUIRE_PSFFIT && (catalog[0].measure[offset].photFlags & 0x01) == 0) { BAD_COUNT(Nbad); continue; } // very loose cut on PSF - APER mag (Map or Mkron) -- this lightly rejects galaxies & other oddities float Mkp = catalog[0].measure[offset].M - Maper; - if (fabs(Mkp) > 1.0) { Nbad ++; continue; } + if (fabs(Mkp) > 1.0) { BAD_COUNT(Nbad); continue; } if (catalog[0].measure[offset].photFlags & code->photomBadMask) { - Nbad++; + BAD_COUNT(Nbad); continue; } - if (catalog[0].measure[offset].photFlags & code->photomPoorMask) { Npoor++; continue;} + if (catalog[0].measure[offset].photFlags & code->photomPoorMask) { BAD_COUNT(Npoor); continue;} // weak test for galaxies - if (Mkp > 0.5) { + if (Mkp > STAR_MAX_PSF_KRON) { nEXT ++; } else { @@ -112,15 +156,23 @@ /* select measurements by mag limit */ mag = PhotCat (&catalog[0].measure[offset], MAG_CLASS_PSF); - if (mag > MAG_LIM) { Nmag ++; continue; } + if (mag > MAG_LIM) { BAD_COUNT(Nmag); continue; } + if (mag < MAG_MIN) { BAD_COUNT(Nmag); continue; } + // note that these two cuts are applied to all photcodes. these are also applied to + // the nominal apparent mags : if the expected sensitivity is very different between + // different filters, or if exposure times have a large range, these cuts may not be + // effective (or may exclude too much). /* select measurements by measurement error */ - if ((SIGMA_LIM > 0) && (catalog[0].measure[offset].dM > SIGMA_LIM)) { Nsigma ++; continue; } + if ((SIGMA_LIM > 0) && (catalog[0].measure[offset].dM > SIGMA_LIM)) { BAD_COUNT(Nsigma); continue; } /* select measurements by mag limit */ if (ImagSelect) { mag = PhotInst (&catalog[0].measure[offset], MAG_CLASS_PSF); - if (mag < ImagMin) { Nimag ++; continue; } - if (mag > ImagMax) { Nimag ++; continue; } - } + if (mag < ImagMin) { BAD_COUNT(Nimag); continue; } + if (mag > ImagMax) { BAD_COUNT(Nimag); continue; } + } + // note that instrumental mag cuts are reasonable for individual exposures using the same + // detector, but are probably not appropriate for stacks with a large range of the number + // of inputs. // count this measurement as valid for this secfilt entry @@ -147,5 +199,5 @@ if (nEXT >= nPSF) { Nmeasure -= Nm; - Ngalaxy ++; + badCount.Ngalaxy++; continue; } @@ -157,5 +209,5 @@ if (Nm <= STAR_TOOFEW) { /* enough measurements total? */ Nmeasure -= Nm; - Nfew ++; + badCount.Nfew++; continue; } @@ -167,5 +219,5 @@ if (!anySecfiltGood) { Nmeasure -= Nm; - Nfew ++; + badCount.Nfew++; continue; } @@ -204,5 +256,14 @@ subcatalog[0].Naverage, subcatalog[0].Nmeasure, catalog[0].Naverage, catalog[0].Nmeasure, catalog[0].filename, NmMin, NmMax); fprintf (stderr, "rejections: %d stars have too few measures:\n %d code, %d time, %d dophot, %d mag, %d sigma, %d imag, %d psfqf, %d Nnan, %d galaxies, %d bad, %d poor\n", - Nfew, Ncode, Ntime, Ndophot, Nmag, Nsigma, Nimag, Npsfqf, Nnan, Ngalaxy, Nbad, Npoor); + badCount.Nfew, badCount.Ncode, badCount.Ntime, badCount.Ndophot, badCount.Nmag, badCount.Nsigma, badCount.Nimag, badCount.Npsfqf, badCount.Nnan, badCount.Ngalaxy, badCount.Nbad, badCount.Npoor); + } + + if (doTest1) { + fprintf (stderr, "TEST1 rejections: %d stars have too few measures: %d code, %d time, %d dophot, %d mag, %d sigma, %d imag, %d psfqf, %d Nnan, %d galaxies, %d bad, %d poor\n", + badTest1.Nfew, badTest1.Ncode, badTest1.Ntime, badTest1.Ndophot, badTest1.Nmag, badTest1.Nsigma, badTest1.Nimag, badTest1.Npsfqf, badTest1.Nnan, badTest1.Ngalaxy, badTest1.Nbad, badTest1.Npoor); + } + if (doTest2) { + fprintf (stderr, "TEST2 rejections: %d stars have too few measures: %d code, %d time, %d dophot, %d mag, %d sigma, %d imag, %d psfqf, %d Nnan, %d galaxies, %d bad, %d poor\n", + badTest2.Nfew, badTest2.Ncode, badTest2.Ntime, badTest2.Ndophot, badTest2.Nmag, badTest2.Nsigma, badTest2.Nimag, badTest2.Npsfqf, badTest2.Nnan, badTest2.Ngalaxy, badTest2.Nbad, badTest2.Npoor); } @@ -212,5 +273,8 @@ } - free (Nvalid) + free (Nvalid); + + if (ftest1) { fclose (ftest1); } + if (ftest2) { fclose (ftest2); } return (TRUE); Index: trunk/Ohana/src/relphot/src/extra.c =================================================================== --- trunk/Ohana/src/relphot/src/extra.c (revision 42807) +++ trunk/Ohana/src/relphot/src/extra.c (revision 42821) @@ -70,13 +70,20 @@ } -// for now (20140710) I need to identify gpc1 stacks explicitly. generalize in the future +// for now (20140710) I need to identify gpc1 & gpc2 stacks explicitly. generalize in the future int isGPC1stack (int photcode) { - if ((photcode == 11000) || (photcode == 14100)) return TRUE; // g-band - if ((photcode == 11100) || (photcode == 14200)) return TRUE; // r-band - if ((photcode == 11200) || (photcode == 14300)) return TRUE; // i-band - if ((photcode == 11300) || (photcode == 14400)) return TRUE; // z-band - if ((photcode == 11400) || (photcode == 14500)) return TRUE; // y-band - if ((photcode == 11500) || (photcode == 14600)) return TRUE; // w-band + // 11000 : GPC1.g.SkyChip + // 11001 : VARIOUS.g.SkyChip (output of xcamera stacks) + // 14100 : SIMTEST.g.SkyChip + // 31000 : GPC2.g.SkyChip + + if (CALIBRATE_STACKS_LIKE_CHIPS) return FALSE; + + if ((photcode == 11000) || (photcode == 11001) || (photcode == 14100) || (photcode == 31000)) return TRUE; // g-band + if ((photcode == 11100) || (photcode == 11101) || (photcode == 14200) || (photcode == 31100)) return TRUE; // r-band + if ((photcode == 11200) || (photcode == 11201) || (photcode == 14300) || (photcode == 31200)) return TRUE; // i-band + if ((photcode == 11300) || (photcode == 11301) || (photcode == 14400) || (photcode == 31300)) return TRUE; // z-band + if ((photcode == 11400) || (photcode == 11401) || (photcode == 14500) || (photcode == 31400)) return TRUE; // y-band + if ((photcode == 11500) || (photcode == 11501) || (photcode == 14600) || (photcode == 31500)) return TRUE; // w-band return FALSE; @@ -85,4 +92,13 @@ // for now (20140710) I need to identify gpc1 stacks explicitly. generalize in the future int isGPC1warp (int photcode) { + + if (CALIBRATE_WARPS_LIKE_CHIPS) return FALSE; + + // 11000 : GPC1.g.SkyChip + // 14100 : SIMTEST.g.SkyChip + + // 12000 : GPC1.g.ForcedWarp + // 15100 : SIMTEST.g.ForcedWarp + // 32000 : GPC2.g.ForcedWarp // diff warps get stack-like photcodes (kind of lame) @@ -97,10 +113,10 @@ } - if ((photcode == 12000) || (photcode == 15100)) return TRUE; // g-band - if ((photcode == 12100) || (photcode == 15200)) return TRUE; // r-band - if ((photcode == 12200) || (photcode == 15300)) return TRUE; // i-band - if ((photcode == 12300) || (photcode == 15400)) return TRUE; // z-band - if ((photcode == 12400) || (photcode == 15500)) return TRUE; // y-band - if ((photcode == 12500) || (photcode == 15600)) return TRUE; // w-band + if ((photcode == 12000) || (photcode == 15100) || (photcode == 32000)) return TRUE; // g-band + if ((photcode == 12100) || (photcode == 15200) || (photcode == 32100)) return TRUE; // r-band + if ((photcode == 12200) || (photcode == 15300) || (photcode == 32200)) return TRUE; // i-band + if ((photcode == 12300) || (photcode == 15400) || (photcode == 32300)) return TRUE; // z-band + if ((photcode == 12400) || (photcode == 15500) || (photcode == 32400)) return TRUE; // y-band + if ((photcode == 12500) || (photcode == 15600) || (photcode == 32500)) return TRUE; // w-band return FALSE; Index: trunk/Ohana/src/relphot/src/fit1d_irls.c =================================================================== --- trunk/Ohana/src/relphot/src/fit1d_irls.c (revision 42807) +++ trunk/Ohana/src/relphot/src/fit1d_irls.c (revision 42821) @@ -374,5 +374,6 @@ StatType stats; liststats_setmode (&stats, "MEDIAN"); - liststats (dataset->alldata->yVector, NULL, NULL, Nvalues, &stats); + liststats (dataset->alldata->yVector, dataset->alldata->dyVector, NULL, Nvalues, &stats); + // liststats sorts yVector and dyVector together double altSigma = (stats.Upper80 - stats.Lower20) / 1.6; // 20% to 80% encompasses 60% of the values, corresponds to the range (-0.85 sigma : +0.85 sigma) @@ -380,5 +381,5 @@ // soften the individual errors with 10% of the scatter above for (int j = 0; j < Nvalues; j++) { - double newSigma = hypot(dataset->alldata->dyVector[j], 0.1*altSigma); + double newSigma = hypot(dataset->alldata->dyVector[j], 0.5*altSigma); dataset->alldata->dyVector[j] = newSigma; } Index: trunk/Ohana/src/relphot/src/initialize.c =================================================================== --- trunk/Ohana/src/relphot/src/initialize.c (revision 42807) +++ trunk/Ohana/src/relphot/src/initialize.c (revision 42821) @@ -50,11 +50,12 @@ fprintf (stderr, "DophotSelect: FALSE\n"); } - fprintf (stderr, "PSF_QF limit: 0.85 (hardwired)\n"); + fprintf (stderr, "PSF_QF limit: %f\n", PSF_QF_MIN); // fprintf (stderr, "Photom Bad Mask: 0x%08x, Photom Poor Mask: 0x%08x\n"); - fprintf (stderr, "MAG_LIM: %f, SIGMA_LIM: %f\n", MAG_LIM, SIGMA_LIM); + fprintf (stderr, "MAG_MIN: %f, MAG_LIM: %f, SIGMA_LIM: %f\n", MAG_MIN, MAG_LIM, SIGMA_LIM); fprintf (stderr, "INST_MAG_MIN: %f, INST_MAG_MAX: %f\n", ImagMin, ImagMax); + fprintf (stderr, "STAR_MAX_PSF_KRON: %f\n", STAR_MAX_PSF_KRON); fprintf (stderr, "STAR_TOOFEW: %d\n", STAR_TOOFEW); Index: trunk/Ohana/src/relphot/src/load_catalogs.c =================================================================== --- trunk/Ohana/src/relphot/src/load_catalogs.c (revision 42807) +++ trunk/Ohana/src/relphot/src/load_catalogs.c (revision 42821) @@ -151,5 +151,7 @@ // TimeSelect -time // DophotSelect - // (note that psfQF is applied rigidly at 0.85, as is the galaxy test) + // PSF_QF_MIN + // STAR_MAX_PSF_KRON + // MAG_MIN // MAG_LIM // SIGMA_LIM @@ -158,6 +160,8 @@ char *command = NULL; - strextend (&command, "relphot_client %s -load %s -hostID %d -D CATDIR %s -hostdir %s -region %f %f %f %f -D CAMERA %s -D MAG_LIM %f -D SIGMA_LIM %f", - PhotcodeList, table->hosts[i].results, table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax, CAMERA, MAG_LIM, SIGMA_LIM); + strextend (&command, "relphot_client %s -load %s -hostID %d -D CATDIR %s -hostdir %s", PhotcodeList, table->hosts[i].results, table->hosts[i].hostID, CATDIR, table->hosts[i].pathname); + strextend (&command, "-region %f %f %f %f -D CAMERA %s", UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax, CAMERA); + strextend (&command, "-D MAG_MIN %f -D MAG_LIM %f -D SIGMA_LIM %f", MAG_MIN, MAG_LIM, SIGMA_LIM); + strextend (&command, "-D STAR_MAX_PSF_KRON %f -D PSF_QF_MIN %f", STAR_MAX_PSF_KRON, PSF_QF_MIN); if (VERBOSE) { strextend (&command, "-v"); } Index: trunk/Ohana/src/relphot/src/reload_catalogs.c =================================================================== --- trunk/Ohana/src/relphot/src/reload_catalogs.c (revision 42807) +++ trunk/Ohana/src/relphot/src/reload_catalogs.c (revision 42821) @@ -251,4 +251,6 @@ if (!(STAGES & STAGE_WARP)) { strextend (&command, "-skip-warp"); } if (!(STAGES & STAGE_STACK)) { strextend (&command, "-skip-stack"); } + if (CALIBRATE_WARPS_LIKE_CHIPS) { strextend (&command, "-calibrate-warps-like-chips"); } + if (CALIBRATE_STACKS_LIKE_CHIPS) { strextend (&command, "-calibrate-stacks-like-chips"); } // deprecate Index: trunk/Ohana/src/relphot/src/relphot.c =================================================================== --- trunk/Ohana/src/relphot/src/relphot.c (revision 42807) +++ trunk/Ohana/src/relphot/src/relphot.c (revision 42821) @@ -33,7 +33,8 @@ switch (mode) { - case UPDATE_IMAGES: + case UPDATE_IMAGES: // -images // calculate zero points for images (may group by exposure [mosaic], night [tgroups]; may measure flat-correction) - // IF CALLED WITH NLOOP == 0, DOES NOT LOAD bcatalog, just loads image table and generates ImageSubset.dat + // IF CALLED WITH NLOOP == 0 (and NOT -only-stacks-and-warps), this step DOES NOT LOAD bcatalog. Instead it just loads + // the image table and applies the Mcal values to the measurements. In parallel context, it generates ImageSubset.dat. // If called with -update, calls reload_catalogs() just like apply_offsets, UNLESS -only-stacks-and-warps is selected relphot_images (skylist);