Index: trunk/Ohana/src/relastro/Makefile =================================================================== --- trunk/Ohana/src/relastro/Makefile (revision 36630) +++ trunk/Ohana/src/relastro/Makefile (revision 36833) @@ -127,4 +127,5 @@ $(SRC)/plotstuff.$(ARCH).o \ $(SRC)/relastroVisual.$(ARCH).o \ +$(SRC)/syncfile.$(ARCH).o \ $(SRC)/BrightCatalog.$(ARCH).o Index: trunk/Ohana/src/relastro/include/relastro.h =================================================================== --- trunk/Ohana/src/relastro/include/relastro.h (revision 36630) +++ trunk/Ohana/src/relastro/include/relastro.h (revision 36833) @@ -83,7 +83,8 @@ double L, M; /* Focal Plane - pixels */ double X, Y; /* Chip Coords - pixels */ - double Mag; - double dMag; - double dPos; + float Mag; + float Color; + float dMag; + float dPos; int mask; int Nmeas; @@ -207,7 +208,7 @@ double MaxDensityValue; -char *PHOTCODE_KEEP_LIST, *PHOTCODE_SKIP_LIST; -int NphotcodesKeep, NphotcodesSkip; -PhotCode **photcodesKeep, **photcodesSkip; +char *PHOTCODE_KEEP_LIST, *PHOTCODE_SKIP_LIST, *PHOTCODE_RESET_LIST; +int NphotcodesKeep, NphotcodesSkip, NphotcodesReset; +PhotCode **photcodesKeep, **photcodesSkip, **photcodesReset; char *PHOTCODE_A_LIST, *PHOTCODE_B_LIST; @@ -217,4 +218,8 @@ SkyRegionSelection SELECTION; +char *DCR_COLOR_POS, *DCR_COLOR_NEG; +PhotCode *DCR_PHOTCODE_POS, *DCR_PHOTCODE_NEG; +int DCR_NSEC_POS, DCR_NSEC_NEG; + int ImagSelect; double ImagMin, ImagMax; @@ -233,4 +238,7 @@ int CLIP_THRESH; int USE_BASIC_CHECK; + +int ExcludeBogus; +double ExcludeBogusRadius; FitMode FIT_MODE; @@ -307,5 +315,5 @@ int liststats PROTO((double *value, double *dvalue, int N, StatType *stats)); int liststats_pos PROTO((double *value, double *dvalue, int N, StatType *stats, int XVERB)); -Catalog *load_catalogs PROTO((SkyList *skylist, int *Ncatalog, int subselect, int hostID, char *hostpath)); +Catalog *load_catalogs PROTO((SkyList *skylist, int *Ncatalog, int subselect, int hostID, char *hostpath, char *syncfile)); int load_images PROTO((FITS_DB *db, SkyList *skylist, int UseFullOverlap)); Image *select_images PROTO((SkyList *skylist, Image *timage, off_t Ntimage, off_t **LineNumber, off_t *Nimage, int UseFullOverlap)); @@ -417,4 +425,6 @@ int setMeanR (double ra_fit, MeasureTiny *measure, Average *average, SecFilt *secfilt); int setMeanD (double dec_fit, MeasureTiny *measure, Average *average, SecFilt *secfilt); +double getMeanR_Big (Measure *measure, Average *average, SecFilt *secfilt); +double getMeanD_Big (Measure *measure, Average *average, SecFilt *secfilt); int setMeanR_Big (double ra_fit, Measure *measure, Average *average, SecFilt *secfilt); int setMeanD_Big (double dec_fit, Measure *measure, Average *average, SecFilt *secfilt); @@ -521,2 +531,9 @@ int ImageTableSave (char *filename, Image *images, off_t Nimages); int select_mosaics_hostregion (RegionHostTable *regionHosts, Image *image, off_t Nimage); + +float getColor (off_t meas, int cat); + +int strextend (char *input, char *format,...); + +int areImagesLoaded (); +int areImagesMatched (); Index: trunk/Ohana/src/relastro/src/BrightCatalog.c =================================================================== --- trunk/Ohana/src/relastro/src/BrightCatalog.c (revision 36630) +++ trunk/Ohana/src/relastro/src/BrightCatalog.c (revision 36833) @@ -167,5 +167,5 @@ GET_COLUMN(measureOffset, "OFF_MEASURE", int); GET_COLUMN(missingOffset, "OFF_MISSING", int); - GET_COLUMN(extendOffset, "OFF_EXTEND", int); + GET_COLUMN(refColor, "REF_COLOR", float); GET_COLUMN(flags, "FLAGS", int); GET_COLUMN(photFlagsUpper, "PHOTFLAGS_U", int); @@ -202,5 +202,5 @@ average[i].measureOffset = measureOffset[i] ; average[i].missingOffset = missingOffset[i] ; - average[i].extendOffset = extendOffset[i] ; + average[i].refColor = refColor[i] ; average[i].flags = flags[i] ; average[i].photFlagsUpper = photFlagsUpper[i] ; @@ -234,5 +234,5 @@ free (measureOffset); free (missingOffset); - free (extendOffset); + free (refColor); free (flags); free (photFlagsUpper); @@ -476,5 +476,5 @@ gfits_define_bintable_column (&theader, "J", "OFF_MEASURE", "offset to first psf measurement", "", 1.0, 1.0*0x8000); gfits_define_bintable_column (&theader, "J", "OFF_MISSING", "offset to first missing obs", "", 1.0, 1.0*0x8000); - gfits_define_bintable_column (&theader, "J", "OFF_EXTEND", "offset to first extended measurement", "", 1.0, 1.0*0x8000); + gfits_define_bintable_column (&theader, "E", "REF_COLOR", "reference color", "", 1.0, 0.0); gfits_define_bintable_column (&theader, "J", "FLAGS", "average object flags (star; ghost; etc)", "", 1.0, 1.0*0x8000); gfits_define_bintable_column (&theader, "J", "PHOTFLAGS_U", "upper bit of 2 bit summary of per-measure photflags", "", 1.0, 1.0*0x8000); @@ -510,5 +510,5 @@ int *measureOffset ; ALLOCATE (measureOffset , int , catalog->Naverage); int *missingOffset ; ALLOCATE (missingOffset , int , catalog->Naverage); - int *extendOffset ; ALLOCATE (extendOffset , int , catalog->Naverage); + float *refColor ; ALLOCATE (refColor , float , catalog->Naverage); int *flags ; ALLOCATE (flags , int , catalog->Naverage); int *photFlagsUpper; ALLOCATE (photFlagsUpper, int , catalog->Naverage); @@ -543,5 +543,5 @@ measureOffset[i] = average[i].measureOffset ; missingOffset[i] = average[i].missingOffset ; - extendOffset[i] = average[i].extendOffset ; + refColor[i] = average[i].refColor ; flags[i] = average[i].flags ; photFlagsUpper[i] = average[i].photFlagsUpper ; @@ -575,5 +575,5 @@ gfits_set_bintable_column (&theader, &ftable, "OFF_MEASURE", measureOffset, catalog->Naverage); gfits_set_bintable_column (&theader, &ftable, "OFF_MISSING", missingOffset, catalog->Naverage); - gfits_set_bintable_column (&theader, &ftable, "OFF_EXTEND", extendOffset, catalog->Naverage); + gfits_set_bintable_column (&theader, &ftable, "REF_COLOR", refColor, catalog->Naverage); gfits_set_bintable_column (&theader, &ftable, "FLAGS", flags, catalog->Naverage); gfits_set_bintable_column (&theader, &ftable, "PHOTFLAGS_U", photFlagsUpper, catalog->Naverage); @@ -605,5 +605,5 @@ free (measureOffset); free (missingOffset); - free (extendOffset); + free (refColor); free (flags); free (photFlagsUpper); Index: trunk/Ohana/src/relastro/src/FitChip.c =================================================================== --- trunk/Ohana/src/relastro/src/FitChip.c (revision 36630) +++ trunk/Ohana/src/relastro/src/FitChip.c (revision 36833) @@ -4,5 +4,5 @@ int FitChip (StarData *raw, StarData *ref, int Nmatch, Image *image) { - int i, NstatFull, Nstat, Niter; + int i, NstatFull, Nstat, Niter, Nkeep; float dLsig, dMsig, dRsig; float dLsigFull, dMsigFull, dRsigFull; @@ -24,4 +24,5 @@ // measure the scatter distribution (use only the bright end detections) + Nkeep = 0; for (i = 0; i < Nmatch; i++) { if (raw[i].mask) continue; @@ -42,5 +43,9 @@ continue; } - } + Nkeep ++; + } + + // I'm rejecting some points from the fit above; I count the remainders and actually + // use that count to set order_use below // figures to assess the fitting process: @@ -50,16 +55,15 @@ relastroVisualPlotChipFit(raw, ref, dRmax, Nmatch); - // fit the requested order polynomial - if (CHIPORDER > 0) { - int Ndof_min = 3; - int order_max = 0.5*(sqrt(4*Nmatch - 4*Ndof_min + 1) - 3); - int order_use = MIN (CHIPORDER, order_max); - if (order_use < 1) { - if (VERBOSE2) fprintf (stderr, "insufficient measurements (%d) for linear fit\n", Nmatch); - image[0].flags |= ID_IMAGE_ASTROM_FEW; - return FALSE; - } - image[0].coords.Npolyterms = order_use; - } + // set the maximum order for the polynomial (based on number of stars kept above) + int order_use = 0; + if (Nkeep > 5) order_use = 1; // 4 stars per polynomial term (per dimension) + if (Nkeep > 30) order_use = 2; // 5 stars per polynomial term (per dimension) + if (Nkeep > 60) order_use = 3; // 6 stars per polynomial term (per dimension) + if (order_use < 1) { + if (VERBOSE2) fprintf (stderr, "insufficient measurements (%d) for linear fit\n", Nkeep); + image[0].flags |= ID_IMAGE_ASTROM_FEW; + return FALSE; + } + image[0].coords.Npolyterms = order_use; if (fit) fit_free (fit); @@ -71,30 +75,4 @@ fit_add (fit, raw[i].X, raw[i].Y, ref[i].L, ref[i].M, raw[i].dPos); } - -# if (0) - // check if the fit has enough data points for the polynomial order - skip = FALSE; - switch (image[0].coords.Npolyterms) { - case 0: - case 1: - skip = (fit[0].Npts < 8); - break; - case 2: - skip = (fit[0].Npts < 11); - break; - case 3: - skip = (fit[0].Npts < 15); - break; - default: - fprintf (stderr, "invalid chip order %d\n", image[0].coords.Npolyterms); - skip = TRUE; - } - if (skip) { - if (VERBOSE2) fprintf (stderr, "insufficient measurements (%d) for requested order (%d)\n", fit[0].Npts, image[0].coords.Npolyterms); - fit_free (fit); - image[0].flags |= ID_IMAGE_ASTROM_FEW; - return FALSE; - } -# endif // measure the fit, update the coords & object coordinates @@ -134,4 +112,20 @@ } + int Ncolor = 0; + float *colorList = NULL; + ALLOCATE (colorList, float, Nmatch); + + for (i = 0; i < Nmatch; i++) { + if (!raw[i].mask) continue; + if (!isnan(ref[i].Color)) { + colorList[Ncolor] = ref[i].Color; + Ncolor ++; + } + } + fsort (colorList, Ncolor); + float colorMedian = (Ncolor > 0) ? colorList[(int)(0.5*Ncolor)] : NAN; + image[0].refColor = colorMedian; + free (colorList); + GetScatterRawRef(&dLsigFull, &dMsigFull, &dRsigFull, &NstatFull, raw, ref, Nmatch, SIGMA_LIM); GetScatterRawRef(&dLsig, &dMsig, &dRsig, &Nstat, raw, ref, Nmatch, IMFIT_SYS_SIGMA_LIM); @@ -151,4 +145,6 @@ image[0].dYpixSys = dMsig; image[0].nFitAstrom = fit[0].Npts; + + // fprintf (stderr, "%s %6.3f %4d %4d\n", image[0].name, image[0].refColor, Ncolor, image[0].nFitAstrom); if (fit) fit_free (fit); Index: trunk/Ohana/src/relastro/src/GetAstromError.c =================================================================== --- trunk/Ohana/src/relastro/src/GetAstromError.c (revision 36630) +++ trunk/Ohana/src/relastro/src/GetAstromError.c (revision 36833) @@ -30,4 +30,5 @@ code = GetPhotcodebyCode (measure[0].photcode); + if (!code) return NAN; // do not raise an exception, just send back the result @@ -72,4 +73,5 @@ code = GetPhotcodebyCode (measure[0].photcode); + if (!code) return NAN; // do not raise an exception, just send back the result Index: trunk/Ohana/src/relastro/src/ImageOps.c =================================================================== --- trunk/Ohana/src/relastro/src/ImageOps.c (revision 36630) +++ trunk/Ohana/src/relastro/src/ImageOps.c (revision 36833) @@ -39,4 +39,16 @@ // N_onImage was 'Nlist' // N_ONIMAGE was 'NLIST' + +int areImagesLoaded () { + + if (image) return TRUE; + return FALSE; +} + +int areImagesMatched () { + + if (MeasureToImage) return TRUE; + return FALSE; +} Image *getimages (off_t *N, off_t **line_number) { @@ -333,4 +345,13 @@ */ +float getColor (off_t meas, int cat) { + + off_t i; + + i = MeasureToImage[cat][meas]; + if (i == -1) return (NAN); + return (image[i].refColor); +} + void plot_images () { @@ -468,10 +489,15 @@ } + X = measure[0].Xccd; + Y = measure[0].Yccd; if (USE_FIXED_PIXCOORDS) { - X = isfinite(measure[0].Xfix) ? measure[0].Xfix : measure[0].Xccd; - Y = isfinite(measure[0].Yfix) ? measure[0].Yfix : measure[0].Yccd; - } else { - X = measure[0].Xccd; - Y = measure[0].Yccd; + if (isfinite(measure[0].Xfix) && isfinite(measure[0].Yfix)) { + float dX = measure[0].Xfix - measure[0].Xccd; + float dY = measure[0].Yfix - measure[0].Yccd; + if (hypot(dX,dY) < 2.0) { + X = measure[0].Xfix; + Y = measure[0].Yfix; + } + } } n = measure[0].averef; @@ -560,5 +586,5 @@ int Noff = NoffRAave + NoffDECave + NoffRAori + NoffDECori; - if (VERBOSE && (Noff > 0)) fprintf (stderr, "Noff ave RA %d, Noff ave DEC %d, Noff ori RA %d, Noff ori DEC %d\n", NoffRAave, NoffDECave, NoffRAori, NoffDECori); + if (VERBOSE2 && (Noff > 0)) fprintf (stderr, "Noff ave RA %d, Noff ave DEC %d, Noff ori RA %d, Noff ori DEC %d\n", NoffRAave, NoffDECave, NoffRAori, NoffDECori); saveOffsets (dPos, nPos, im); @@ -609,10 +635,15 @@ Measure *measure = &catalog[c].measure[m]; + X = measure[0].Xccd; + Y = measure[0].Yccd; if (USE_FIXED_PIXCOORDS) { - X = isfinite(measure[0].Xfix) ? measure[0].Xfix : measure[0].Xccd; - Y = isfinite(measure[0].Yfix) ? measure[0].Yfix : measure[0].Yccd; - } else { - X = measure[0].Xccd; - Y = measure[0].Yccd; + if (isfinite(measure[0].Xfix) && isfinite(measure[0].Yfix)) { + float dX = measure[0].Xfix - measure[0].Xccd; + float dY = measure[0].Yfix - measure[0].Yccd; + if (hypot(dX,dY) < 2.0) { + X = measure[0].Xfix; + Y = measure[0].Yfix; + } + } } n = measure[0].averef; @@ -679,12 +710,16 @@ /* apply the current image transformation or use the current value of R+dR, D+dD? */ + raw[i].X = measure[0].Xccd; + raw[i].Y = measure[0].Yccd; if (USE_FIXED_PIXCOORDS) { - raw[i].X = isfinite(measure[0].Xfix) ? measure[0].Xfix : measure[0].Xccd; - raw[i].Y = isfinite(measure[0].Yfix) ? measure[0].Yfix : measure[0].Yccd; - } else { - raw[i].X = measure[0].Xccd; - raw[i].Y = measure[0].Yccd; - } - + if (isfinite(measure[0].Xfix) && isfinite(measure[0].Yfix)) { + float dX = measure[0].Xfix - measure[0].Xccd; + float dY = measure[0].Yfix - measure[0].Yccd; + if (hypot(dX,dY) < 2.0) { + raw[i].X = measure[0].Xfix; + raw[i].Y = measure[0].Yfix; + } + } + } raw[i].Mag = measure[0].M; raw[i].dMag = measure[0].dM; @@ -738,4 +773,6 @@ StarData *ref; + int Nsecfilt = GetPhotcodeNsecfilt(); + mosaic = NULL; moscoords = NULL; @@ -765,4 +802,8 @@ ref[i].dMag = measure[0].dM; ref[i].dPos = GetAstromErrorTiny (&measure[0], ERROR_MODE_POS); + + if ((DCR_NSEC_POS >= 0) && (DCR_NSEC_NEG >= -1)) { + ref[i].Color = catalog[c].secfilt[n*Nsecfilt + DCR_NSEC_POS].M - catalog[c].secfilt[n*Nsecfilt + DCR_NSEC_NEG].M; + } ref[i].mask = FALSE; @@ -859,5 +900,5 @@ N = 0; for (k = 0; k < catalog[0].average[j].Nmeasure; k++) { - //reset flag on each invocation + // reset flag on each invocation measure[k].dbFlags &= ~ID_MEAS_POOR_ASTROM; @@ -898,4 +939,7 @@ double *R, *D, *dR, *dD, *d2; StatType statsR, statsD; + + // XXX we are not going to use this for now + return; Nsecfilt = GetPhotcodeNsecfilt(); @@ -1063,4 +1107,5 @@ } else { code = GetPhotcodebyCode (measure[0].photcode); + if (!code) return FALSE; mask = code[0].astromBadMask; } @@ -1129,4 +1174,5 @@ } else { code = GetPhotcodebyCode (measure[0].photcode); + if (!code) return FALSE; mask = code[0].astromBadMask; } Index: trunk/Ohana/src/relastro/src/UpdateObjectOffsets.c =================================================================== --- trunk/Ohana/src/relastro/src/UpdateObjectOffsets.c (revision 36630) +++ trunk/Ohana/src/relastro/src/UpdateObjectOffsets.c (revision 36833) @@ -126,41 +126,44 @@ table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax, STATMODE, MIN_ERROR); - char tmpline[1024]; - if (FIT_MODE == FIT_PM_ONLY) { snprintf (tmpline, 1024, "%s -pm", command); strcpy (command, tmpline); } - if (FIT_MODE == FIT_PAR_ONLY) { snprintf (tmpline, 1024, "%s -par", command); strcpy (command, tmpline); } - if (FIT_MODE == FIT_PM_AND_PAR) { snprintf (tmpline, 1024, "%s -pmpar", command); strcpy (command, tmpline); } - - if (VERBOSE) { snprintf (tmpline, 1024, "%s -v", command); strcpy (command, tmpline); } - if (VERBOSE2) { snprintf (tmpline, 1024, "%s -vv", command); strcpy (command, tmpline); } - if (RESET) { snprintf (tmpline, 1024, "%s -reset", command); strcpy (command, tmpline); } - if (UPDATE) { snprintf (tmpline, 1024, "%s -update", command); strcpy (command, tmpline); } - - if (RESET_BAD_IMAGES) { snprintf (tmpline, 1024, "%s -reset-bad-images", command); strcpy (command, tmpline); } - - if (ImagSelect) { snprintf (tmpline, 1024, "%s -instmag %f %f", command, ImagMin, ImagMax); strcpy (command, tmpline); } - if (MaxDensityUse) { snprintf (tmpline, 1024, "%s -max-density %f", command, MaxDensityValue); strcpy (command, tmpline); } + if (FIT_MODE == FIT_PM_ONLY) strextend (command, "-pm"); + if (FIT_MODE == FIT_PAR_ONLY) strextend (command, "-par"); + if (FIT_MODE == FIT_PM_AND_PAR) strextend (command, "-pmpar"); + + if (VERBOSE) strextend (command, "-v"); + if (VERBOSE2) strextend (command, "-vv"); + if (RESET) strextend (command, "-reset"); + if (UPDATE) strextend (command, "-update"); + + if (RESET_BAD_IMAGES) strextend (command, "-reset-bad-images"); + + if (ImagSelect) strextend (command, "-instmag %f %f", ImagMin, ImagMax); + if (MaxDensityUse) strextend (command, "-max-density %f", MaxDensityValue); - if (USE_BASIC_CHECK) { snprintf (tmpline, 1024, "%s -basic-image-search", command); strcpy (command, tmpline); } - if (FlagOutlier) { snprintf (tmpline, 1024, "%s -clip %d", command, CLIP_THRESH); strcpy (command, tmpline); } + if (USE_BASIC_CHECK) strextend (command, "-basic-image-search"); + if (FlagOutlier) strextend (command, "-clip %d", CLIP_THRESH); + if (ExcludeBogus) strextend (command, "-exclude-bogus %f", ExcludeBogusRadius); - if (USE_FIXED_PIXCOORDS) { snprintf (tmpline, 1024, "%s -D USE_FIXED_PIXCOORDS 1", command); strcpy (command, tmpline); } - - if (PHOTCODE_KEEP_LIST) { snprintf (tmpline, 1024, "%s +photcode %s", command, PHOTCODE_KEEP_LIST); strcpy (command, tmpline); } - if (PHOTCODE_SKIP_LIST) { snprintf (tmpline, 1024, "%s -photcode %s", command, PHOTCODE_SKIP_LIST); strcpy (command, tmpline); } - if (PhotFlagSelect) { snprintf (tmpline, 1024, "%s +photflags", command); strcpy (command, tmpline); } - if (PhotFlagBad) { snprintf (tmpline, 1024, "%s +photflagbad %d", command, PhotFlagBad); strcpy (command, tmpline); } - if (PhotFlagPoor) { snprintf (tmpline, 1024, "%s +photflagpoor %d", command, PhotFlagPoor); strcpy (command, tmpline); } + if (USE_FIXED_PIXCOORDS) strextend (command, "-D USE_FIXED_PIXCOORDS 1"); + + if (PHOTCODE_KEEP_LIST) strextend (command, "+photcode %s", PHOTCODE_KEEP_LIST); + if (PHOTCODE_SKIP_LIST) strextend (command, "-photcode %s", PHOTCODE_SKIP_LIST); + if (PhotFlagSelect) strextend (command, "+photflags"); + if (PhotFlagBad) strextend (command, "+photflagbad %d", PhotFlagBad); + if (PhotFlagPoor) strextend (command, "+photflagpoor %d", PhotFlagPoor); // XXX note that the above pass in the flag as decimal -- also note that args.c cannot handle 0xHEX values - if (MinBadQF > 0.0) { snprintf (tmpline, 1024, "%s -min-bad-psfqf %f", command, MinBadQF); strcpy (command, tmpline); } - if (MaxMeanOffset != 10.0) { snprintf (tmpline, 1024, "%s -max-mean-offset %f", command, MaxMeanOffset); strcpy (command, tmpline); } + if (DCR_COLOR_POS && DCR_COLOR_NEG) { + strextend (command, "-dcr-color %s %s", DCR_COLOR_POS, DCR_COLOR_NEG); + } + + if (MinBadQF > 0.0) strextend (command, "-min-bad-psfqf %f", MinBadQF); + if (MaxMeanOffset != 10.0) strextend (command, "-max-mean-offset %f", MaxMeanOffset); if (TimeSelect) { char *tstart = ohana_sec_to_date (TSTART); char *tstop = ohana_sec_to_date (TSTOP); - snprintf (tmpline, 1024, "%s -time %s %s", command, tstart, tstop); + strextend (command, "-time %s %s", tstart, tstop); free (tstart); free (tstop); - strcpy (command, tmpline); } fprintf (stderr, "command: %s\n", command); Index: trunk/Ohana/src/relastro/src/UpdateObjects.c =================================================================== --- trunk/Ohana/src/relastro/src/UpdateObjects.c (revision 36630) +++ trunk/Ohana/src/relastro/src/UpdateObjects.c (revision 36833) @@ -10,4 +10,5 @@ static double *T; static double *dT; +static double *C; void initObjectData (Catalog *catalog, int Ncatalog) { @@ -36,4 +37,6 @@ ALLOCATE (pX, double, MAX (1, Nmax)); ALLOCATE (pY, double, MAX (1, Nmax)); + + ALLOCATE (C, double, MAX (1, Nmax)); } @@ -54,4 +57,6 @@ free (pX); free (pY); + + free (C); } @@ -61,5 +66,5 @@ off_t j, k, m; - int i, N, Nsecfilt, mode, result, status, XVERB; + int i, N, Nc, Nsecfilt, mode, result, status, XVERB; StatType statsR, statsD; Coords coords; @@ -75,4 +80,6 @@ memset (&fitPAR, 0, sizeof(fitPAR)); initObjectData (catalog, Ncatalog); + + int setRefColor = areImagesMatched(); /* project coordinates to a plane centered on the object with units of arcsec */ @@ -112,8 +119,16 @@ fitPAR.chisq = NAN; + // if we fail to fit the astrometry for some reason, we need to set/reset these + catalog[i].average[j].flags |= ID_STAR_NO_ASTROM; + catalog[i].average[j].ChiSqAve = NAN; + catalog[i].average[j].ChiSqPM = NAN; + catalog[i].average[j].ChiSqPar = NAN; + catalog[i].average[j].Npos = 0; + if (catalog[i].average[j].Nmeasure == 0) { - continue; - } - + continue; + } + + Nc = 0; N = 0; m = catalog[i].average[j].measureOffset; @@ -138,22 +153,38 @@ } + // reset the bit to note that a detection was used (or not) + measure[k].dbFlags &= ~ID_MEAS_USED_OBJ; + if (measureBig) { measureBig[k].dbFlags &= ~ID_MEAS_USED_OBJ; } + // does the measurement pass the supplied filtering constraints? // MeasFilterTestTiny does not test psfQF + // exclude bad detections based on: photcodes, psfQF, time range, photflags & astromBadMask, mag_inst int keepMeasure = measureBig ? MeasFilterTest(&measureBig[k], FALSE) : MeasFilterTestTiny(&measure[k], FALSE); if (!keepMeasure) { - measure[k].dbFlags &= ~ID_MEAS_USED_OBJ; - if (measureBig) { measureBig[k].dbFlags &= ~ID_MEAS_USED_OBJ; } continue; } + // mark (as POOR) any measurements which are deviant from the mean by > ExcludeBogusRadius + if (ExcludeBogus) { + double Ri = getMeanR (&measure[k], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]); + double Di = getMeanD (&measure[k], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]); + coords.crval1 = catalog[i].average[j].R; + coords.crval2 = catalog[i].average[j].D; + double Xi, Yi; + RD_to_XY (&Xi, &Yi, Ri, Di, &coords); + double radius = hypot(Xi, Yi); + if (radius > ExcludeBogusRadius) { + measure[k].dbFlags |= ID_MEAS_POOR_ASTROM; + if (measureBig) { measureBig[k].dbFlags |= ID_MEAS_POOR_ASTROM; } + continue; + } + measure[k].dbFlags &= ~ID_MEAS_POOR_ASTROM; + if (measureBig) { measureBig[k].dbFlags &= ~ID_MEAS_POOR_ASTROM; } + } + // outlier rejection - if (FlagOutlier && (measure[k].dbFlags & ID_MEAS_POOR_ASTROM)) { - measure[k].dbFlags &= ~ID_MEAS_USED_OBJ; - if (measureBig) { measureBig[k].dbFlags &= ~ID_MEAS_USED_OBJ; } + if (FALSE && FlagOutlier && (measure[k].dbFlags & ID_MEAS_POOR_ASTROM)) { continue; } - - measure[k].dbFlags |= ID_MEAS_USED_OBJ; - if (measureBig) { measureBig[k].dbFlags |= ID_MEAS_USED_OBJ; } R[N] = getMeanR (&measure[k], &catalog[i].average[j], &catalog[i].secfilt[j*Nsecfilt]); @@ -170,14 +201,6 @@ // allow a given photcode or measurement to be // ignored if the error is NAN (for photcode, set astromErrSys to NaN) - if (isnan(dX[N])) { - measure[k].dbFlags &= ~ID_MEAS_USED_OBJ; - if (measureBig) { measureBig[k].dbFlags &= ~ID_MEAS_USED_OBJ; } - continue; - } - if (isnan(dY[N])) { - measure[k].dbFlags &= ~ID_MEAS_USED_OBJ; - if (measureBig) { measureBig[k].dbFlags &= ~ID_MEAS_USED_OBJ; } - continue; - } + if (isnan(dX[N])) continue; + if (isnan(dY[N])) continue; // add systematic error in quadrature, if desired @@ -201,12 +224,21 @@ dD[N] = dY[N] / 3600.0; + if (setRefColor) { + float color = getColor (m+k, i); + if (!isnan(color)) { + C[Nc] = color; + Nc++; + } + } + + measure[k].dbFlags |= ID_MEAS_USED_OBJ; + if (measureBig) { measureBig[k].dbFlags |= ID_MEAS_USED_OBJ; } + N++; - } + } // loop over measurements : catalog[i].average[j].Nmeasure // if we have too few good detections for the desired fit, or too limited a // baseline, use a fit with fewer parameters. XXX if we have too few measurements // for even the average position, consider including the lower-quality detections? - - catalog[i].average[j].flags &= ~ID_STAR_NO_ASTROM; // find Tmin & Tmax from the list of accepted measurements @@ -231,10 +263,4 @@ // too few measurements for average position (require 2 values) if (N < SRC_MEAS_TOOFEW) { - // XXX need to define PHOTOM and ASTROM object flags - // XXX reset the average value fields? - catalog[i].average[j].flags |= ID_STAR_NO_ASTROM; - catalog[i].average[j].ChiSqAve = NAN; - catalog[i].average[j].ChiSqPM = NAN; - catalog[i].average[j].ChiSqPar = NAN; if (N < 2) continue; } @@ -268,9 +294,13 @@ catalog[i].average[j].flags |= ID_STAR_FIT_PM; Npm ++; + + // XXX a hard-wired hack... + if ((fabs(fitPM.uR) > 2.0) || (fabs(fitPM.uD) > 2.0)) { + mode = FIT_AVERAGE; + catalog[i].average[j].flags |= ID_STAR_BAD_PM; + } } if (mode == FIT_PM_AND_PAR) { - // fprintf (stderr, "parallax fitting is still untested (%s, %d)\n", __FILE__, __LINE__); - float pXmin = +2.0; float pXmax = -2.0; @@ -296,4 +326,10 @@ catalog[i].average[j].flags |= ID_STAR_FIT_PAR; Npar ++; + + // XXX a hard-wired hack... + if ((fabs(fitPAR.uR) > 2.0) || (fabs(fitPAR.uD) > 2.0)) { + mode = FIT_AVERAGE; + catalog[i].average[j].flags |= ID_STAR_BAD_PM; + } } else { // need to set mode = FIT_PM_ONLY if we do not fit for parallax @@ -321,4 +357,10 @@ catalog[i].average[j].flags |= ID_STAR_FIT_AVE; Nave ++; + } + + if (setRefColor) { + dsort (C, Nc); + float colorMedian = (Nc > 0) ? C[(int)(0.5*Nc)] : NAN; + catalog[i].average[j].refColor = colorMedian; } @@ -411,4 +453,8 @@ catalog[i].average[j].Trange = (Trange * 86400 * 365.25); catalog[i].average[j].Npos = fit.Nfit; + + // unset the NO_ASTROM bit (not(NO_ASTROM) == HAVE_ASTROM) + catalog[i].average[j].flags &= ~ID_STAR_NO_ASTROM; + if (XVERB) fprintf (stderr, "%f %f -> %f %f (%f,%f) pm=(%f %f) chisq=(%f, %f, %f)\n", catalog[i].average[j].R, Index: trunk/Ohana/src/relastro/src/args.c =================================================================== --- trunk/Ohana/src/relastro/src/args.c (revision 36630) +++ trunk/Ohana/src/relastro/src/args.c (revision 36833) @@ -179,8 +179,18 @@ FlagOutlier = FALSE; if ((N = get_argument (argc, argv, "-clip"))) { + fprintf (stderr, "-clip is currently disabled\n"); remove_argument (N, &argc, argv); CLIP_THRESH = atof (argv[N]); remove_argument (N, &argc, argv); FlagOutlier = TRUE; + } + + ExcludeBogus = FALSE; + ExcludeBogusRadius = 0.0; + if ((N = get_argument (argc, argv, "-exclude-bogus"))) { + remove_argument (N, &argc, argv); + ExcludeBogusRadius = atof (argv[N]); + remove_argument (N, &argc, argv); + ExcludeBogus = TRUE; } @@ -284,4 +294,21 @@ } + DCR_COLOR_POS = NULL; + DCR_COLOR_NEG = NULL; + if ((N = get_argument (argc, argv, "-dcr-color"))) { + remove_argument (N, &argc, argv); + DCR_COLOR_POS = strcreate(argv[N]); + remove_argument (N, &argc, argv); + DCR_COLOR_NEG = strcreate(argv[N]); + remove_argument (N, &argc, argv); + } + + PHOTCODE_RESET_LIST = NULL; + if ((N = get_argument (argc, argv, "-reset-to-photcode"))) { + remove_argument (N, &argc, argv); + PHOTCODE_RESET_LIST = strcreate(argv[N]); + remove_argument (N, &argc, argv); + } + PHOTCODE_KEEP_LIST = NULL; if ((N = get_argument (argc, argv, "+photcode"))) { @@ -291,19 +318,21 @@ } - PHOTCODE_SKIP_LIST = strcreate("SCOS.103a.E,SCOS.4414.OG590,SCOS.4415.OG590,SCOS.IIIaF.OG590,SCOS.IIIaF.RG610,SCOS.IIIaF.RG630,SCOS.IIIaJ.GG385,SCOS.IIIaJ.GG395,SCOS.IVN.RG715,SCOS.IVN.RG9"); + char *SuperCOSMOS_SKIP = strcreate("SCOS.103a.E,SCOS.4414.OG590,SCOS.4415.OG590,SCOS.IIIaF.OG590,SCOS.IIIaF.RG610,SCOS.IIIaF.RG630,SCOS.IIIaJ.GG385,SCOS.IIIaJ.GG395,SCOS.IVN.RG715,SCOS.IVN.RG9"); + + PHOTCODE_SKIP_LIST = NULL; if ((N = get_argument (argc, argv, "-photcode"))) { remove_argument (N, &argc, argv); - char *tmp1 = strcreate(argv[N]); - - int Ntotal = strlen(tmp1) + strlen(PHOTCODE_SKIP_LIST) + 5; - - char *tmp2 = NULL; - ALLOCATE (tmp2, char, Ntotal); - snprintf (tmp2, Ntotal, "%s,%s", PHOTCODE_SKIP_LIST, tmp1); - - free (tmp1); - free (PHOTCODE_SKIP_LIST); - - PHOTCODE_SKIP_LIST = tmp2; + char *RawSkip = strcreate(argv[N]); + + char *GotSkip = strstr (RawSkip, SuperCOSMOS_SKIP); + if (!GotSkip) { + int Ntotal = strlen(RawSkip) + strlen(SuperCOSMOS_SKIP) + 5; + ALLOCATE (PHOTCODE_SKIP_LIST, char, Ntotal); + snprintf (PHOTCODE_SKIP_LIST, Ntotal, "%s,%s", SuperCOSMOS_SKIP, RawSkip); + free (RawSkip); + free (SuperCOSMOS_SKIP); + } else { + PHOTCODE_SKIP_LIST = RawSkip; + } remove_argument (N, &argc, argv); } @@ -582,8 +611,18 @@ FlagOutlier = FALSE; if ((N = get_argument (argc, argv, "-clip"))) { + fprintf (stderr, "-clip is currently disabled\n"); remove_argument (N, &argc, argv); CLIP_THRESH = atof (argv[N]); remove_argument (N, &argc, argv); FlagOutlier = TRUE; + } + + ExcludeBogus = FALSE; + ExcludeBogusRadius = 0.0; + if ((N = get_argument (argc, argv, "-exclude-bogus"))) { + remove_argument (N, &argc, argv); + ExcludeBogusRadius = atof (argv[N]); + remove_argument (N, &argc, argv); + ExcludeBogus = TRUE; } @@ -655,4 +694,21 @@ remove_argument (N, &argc, argv); TimeSelect = TRUE; + } + + DCR_COLOR_POS = NULL; + DCR_COLOR_NEG = NULL; + if ((N = get_argument (argc, argv, "-dcr-color"))) { + remove_argument (N, &argc, argv); + DCR_COLOR_POS = strcreate(argv[N]); + remove_argument (N, &argc, argv); + DCR_COLOR_NEG = strcreate(argv[N]); + remove_argument (N, &argc, argv); + } + + PHOTCODE_RESET_LIST = NULL; + if ((N = get_argument (argc, argv, "-reset-to-photcode"))) { + remove_argument (N, &argc, argv); + PHOTCODE_RESET_LIST = strcreate(argv[N]); + remove_argument (N, &argc, argv); } @@ -861,2 +917,15 @@ exit (2); } + +int strextend (char *input, char *format,...) { + + char tmpextra[1024], tmpline[1024]; + va_list argp; + + va_start (argp, format); + vsnprintf (tmpextra, 1024, format, argp); + snprintf (tmpline, 1024, "%s %s", input, tmpextra); + strcpy (input, tmpline); + + return TRUE; +} Index: trunk/Ohana/src/relastro/src/bcatalog.c =================================================================== --- trunk/Ohana/src/relastro/src/bcatalog.c (revision 36630) +++ trunk/Ohana/src/relastro/src/bcatalog.c (revision 36833) @@ -1,8 +1,21 @@ # include "relastro.h" +static int Nkeep1 = 0; +static int Nkeep2 = 0; static int Nskip1 = 0; static int Nskip2 = 0; -static unsigned int Tref = 1323353985; -static short Cref = 10355; + +FILE *fbogus = NULL; +static int NskipBogus = 0; + +// test image: 2013/06/15,13:25:51, GPC1.r.XY50 +static int CHECK_TEST_IMAGE = FALSE; +// static unsigned int Tref = 1378812312; +// static short Cref = 10001; + +static unsigned int Tref = 1379570672; +static short Cref = 10341; + +int LimitDensityCatalog_ByNmeasureGrid (Catalog *subcatalog, Catalog *oldcatalog); int bcatalog (Catalog *subcatalog, Catalog *catalog) { @@ -11,4 +24,16 @@ off_t NAVERAGE, NMEASURE, Naverage, Nmeasure, Nm; int Nsecfilt; + Coords coords; + + /* for outlier rejection, project coordinates to a plane centered on the object with units of arcsec */ + coords.crval1 = 0; + coords.crval2 = 0; + coords.crpix1 = 0; + coords.crpix2 = 0; + coords.cdelt1 = coords.cdelt2 = 1.0 / 3600.0; + coords.pc1_1 = coords.pc2_2 = 1.0; + coords.pc1_2 = coords.pc2_1 = 0.0; + coords.Npolyterms = 1; + strcpy (coords.ctype, "DEC--SIN"); // XXX in the future, use catalog[0].Nsecfilt only? allow catalogs to have variable Nsecfilt? @@ -34,4 +59,13 @@ ID_STAR_NO_ASTROM ; + if (VERBOSE2 && ExcludeBogus && (fbogus == NULL)) { + char name[1024]; + snprintf (name, 1024, "%s/bogus.%02d.dat", CATDIR, HOST_ID); + fbogus = fopen (name, "w"); + if (!fbogus) { + fprintf (stderr, "trouble opening bogus detection dump : %s\n", name); + } + } + /* exclude stars not in range or with too few measurements */ for (i = 0; i < catalog[0].Naverage; i++) { @@ -60,12 +94,16 @@ if (!MeasFilterTest(&catalog[0].measure[offset], TRUE)) { catalog[0].measure[offset].dbFlags &= ~ID_MEAS_USED_CHIP; - if (FALSE && (abs(catalog[0].measure[offset].t - Tref) < 10) && (catalog[0].measure[offset].photcode == Cref)) { + if (CHECK_TEST_IMAGE && (abs(catalog[0].measure[offset].t - Tref) < 10) && (catalog[0].measure[offset].photcode == Cref)) { Nskip1 ++; } continue; } - - // filter out outliers - if (FlagOutlier && (catalog[0].measure[offset].dbFlags & ID_MEAS_POOR_ASTROM)) { + if (CHECK_TEST_IMAGE && (abs(catalog[0].measure[offset].t - Tref) < 10) && (catalog[0].measure[offset].photcode == Cref)) { + Nkeep1 ++; + } + + // filter out outliers - these are detections inconsistent with the offset distribution + // XXX disable this for now + if (FALSE && FlagOutlier && (catalog[0].measure[offset].dbFlags & ID_MEAS_POOR_ASTROM)) { catalog[0].measure[offset].dbFlags &= ~ID_MEAS_USED_CHIP; if (FALSE && (abs(catalog[0].measure[offset].t - Tref) < 10) && (catalog[0].measure[offset].photcode == Cref)) { @@ -76,4 +114,25 @@ catalog[0].measure[offset].dbFlags |= ID_MEAS_USED_CHIP; + // exclude bogus + if (ExcludeBogus) { + double Ri = getMeanR_Big (&catalog[0].measure[offset], &catalog[0].average[i], &catalog[0].secfilt[i*Nsecfilt]); + double Di = getMeanD_Big (&catalog[0].measure[offset], &catalog[0].average[i], &catalog[0].secfilt[i*Nsecfilt]); + coords.crval1 = catalog[0].average[i].R; + coords.crval2 = catalog[0].average[i].D; + double Xi, Yi; + RD_to_XY (&Xi, &Yi, Ri, Di, &coords); + double radius = hypot(Xi, Yi); + if (radius > ExcludeBogusRadius) { + NskipBogus ++; + if (VERBOSE2) { + FILE *foutput = fbogus ? fbogus : stderr; + char *date = ohana_sec_to_date(catalog[0].measure[offset].t); + fprintf (foutput, "exclude bogus: %10.6f %10.6f : %10.6f %10.6f : %6.2f %6.2f : %6.2f : %5d %s\n", catalog[0].average[i].R, catalog[0].average[i].D, Ri, Di, Xi, Yi, radius, catalog[0].measure[offset].photcode, date); + free (date); + } + continue; + } + } + // re-assess on each run of relastro if a measurement should be used @@ -87,4 +146,8 @@ // but before the final average properties are calculated, these measurements may be // allowed. + + if (CHECK_TEST_IMAGE && (abs(catalog[0].measure[offset].t - Tref) < 10) && (catalog[0].measure[offset].photcode == Cref)) { + Nkeep2 ++; + } CopyMeasureToTiny (&subcatalog[0].measureT[Nmeasure], &catalog[0].measure[offset]); @@ -125,5 +188,5 @@ // limit the total number of stars in the catalog if (MaxDensityUse) { - LimitDensityCatalog_ByNmeasure (subcatalog, catalog); + LimitDensityCatalog_ByNmeasureGrid (subcatalog, catalog); } else { if (VERBOSE2) { @@ -134,9 +197,18 @@ } } + if (CHECK_TEST_IMAGE && (Nkeep1 + Nkeep2 + Nskip1 + Nskip2 > 0)) { + fprintf (stderr, "kept %d %d, skipped %d %d\n", Nkeep1, Nkeep2, Nskip1, Nskip2); + } + return (TRUE); } void bcatalog_show_skips () { - fprintf (stderr, "Nskip: %d, %d\n", Nskip1, Nskip2); + if (ExcludeBogus) { + if (fbogus) fclose (fbogus); + fprintf (stderr, "NskipBogus: %d\n", NskipBogus); + // fprintf (stderr, "Nskip: %d, %d\n", Nskip1, Nskip2); + // fprintf (stderr, "Nkeep: %d, %d\n", Nkeep1, Nkeep2); + } } @@ -155,4 +227,270 @@ # undef COMPARE +} + +/* generate a grid in a locally-projected coordinate system, try to select average entries + from each grid cell in decending Nmeasure order. + */ +int LimitDensityCatalog_ByNmeasureGrid (Catalog *subcatalog, Catalog *oldcatalog) { + + off_t i, j; + int ix, iy; + + Catalog tmpcatalog; + + double Rmin, Rmax, Dmin, Dmax; + + int Nsecfilt = GetPhotcodeNsecfilt (); + + gfits_scan (&oldcatalog[0].header, "RA0", "%lf", 1, &Rmin); + gfits_scan (&oldcatalog[0].header, "DEC0", "%lf", 1, &Dmin); + gfits_scan (&oldcatalog[0].header, "RA1", "%lf", 1, &Rmax); + gfits_scan (&oldcatalog[0].header, "DEC1", "%lf", 1, &Dmax); + + if (VERBOSE2) fprintf (stderr, "extracting from catalog covering region %f,%f to %f,%f\n", Rmin, Dmin, Rmax, Dmax); + + float AREA = fabs(Dmax - Dmin) * fabs(Rmax - Rmin) * cos (0.5*RAD_DEG*(Dmax + Dmin)); + assert (AREA > 0); + + off_t Nmax = MaxDensityValue * AREA; + if (subcatalog[0].Naverage <= Nmax) { + if (VERBOSE) { + fprintf (stderr, "subcatalog has less than the max density\n"); + } + return (TRUE); + } + + off_t Naverage = subcatalog[0].Naverage; + + // generate a grid in locally projected space + double Rc = 0.5*(Rmin + Rmax); + double Dc = 0.5*(Dmin + Dmax); + + /* project coordinates to a plane centered on the object with units of arcsec */ + Coords coords; + coords.crval1 = Rc; + coords.crval2 = Dc; + coords.crpix1 = 0; + coords.crpix2 = 0; + coords.cdelt1 = coords.cdelt2 = 1.0 / 3600.0; + coords.pc1_1 = coords.pc2_2 = 1.0; + coords.pc1_2 = coords.pc2_1 = 0.0; + coords.Npolyterms = 1; + strcpy (coords.ctype, "DEC--SIN"); + + // convert all average R,D values to X,Y: + double *X, *Y; + ALLOCATE (X, double, Naverage); + ALLOCATE (Y, double, Naverage); + float Xmin = +10000.0, Ymin = +10000.0; + float Xmax = -10000.0, Ymax = -10000.0; + for (i = 0; i < Naverage; i++) { + X[i] = NAN; + Y[i] = NAN; + // skip any stars which are outside of nominal catalog range + if (subcatalog[0].average[i].R < Rmin) continue; + if (subcatalog[0].average[i].R > Rmax) continue; + if (subcatalog[0].average[i].D < Dmin) continue; + if (subcatalog[0].average[i].D > Dmax) continue; + RD_to_XY (&X[i], &Y[i], subcatalog[0].average[i].R, subcatalog[0].average[i].D, &coords); + Xmin = MIN (Xmin, X[i]); + Xmax = MAX (Xmax, X[i]); + Ymin = MIN (Ymin, Y[i]); + Ymax = MAX (Ymax, Y[i]); + } + + // how many grid cells? what is the grid spacing? + float dX = Xmax - Xmin; + float dY = Ymax - Ymin; + + // *** XXX for the moment, I'm using a hard-wired cell size (200 arcsec ~ 3.3 arcmin) + int NX = (int)(dX / 200) + 1; + int NY = (int)(dY / 200) + 1; + // fprintf (stderr, "Density Grid: %d x %d\n", NX, NY); + // XXX check that NX,NY are sensible (5 degrees / 200 arcsec seems like the absolute max) + if (NX > 1000) { + fprintf (stderr, "serious problem with %s: NX = %d\n", subcatalog[0].filename, NX); + exit (3); + } + if (NY > 1000) { + fprintf (stderr, "serious problem with %s: NY = %d\n", subcatalog[0].filename, NY); + exit (3); + } + + // kind of ugly : generate a grid of index, Nmeasure arrays + // to be filled below (I also need NN and Nn to track the number of + // entries in each). + int **NN_grid; + int **Nn_grid; + int ***Nm_grid; + off_t ***idxgrid; + ALLOCATE (NN_grid, int *, NX); + ALLOCATE (Nn_grid, int *, NX); + ALLOCATE (Nm_grid, int **, NX); + ALLOCATE (idxgrid, off_t **, NX); + + for (ix = 0; ix < NX; ix++) { + ALLOCATE (NN_grid[ix], int, NY); + ALLOCATE (Nn_grid[ix], int, NY); + ALLOCATE (Nm_grid[ix], int *, NY); + ALLOCATE (idxgrid[ix], off_t *, NY); + for (iy = 0; iy < NY; iy++) { + Nn_grid[ix][iy] = 0; + NN_grid[ix][iy] = 100; + ALLOCATE (Nm_grid[ix][iy], int, NN_grid[ix][iy]); + ALLOCATE (idxgrid[ix][iy], off_t, NN_grid[ix][iy]); + } + } + + // assign all of the average entries to a grid cell + for (i = 0; i < Naverage; i++) { + if (isnan(X[i])) continue; + if (isnan(Y[i])) continue; + ix = (X[i] - Xmin) / 200.0; + iy = (Y[i] - Ymin) / 200.0; + int Nn = Nn_grid[ix][iy]; + Nm_grid[ix][iy][Nn] = subcatalog[0].average[i].Nmeasure; + + // if we are resetting to a given photcode, we need to have that photcode... + if (NphotcodesReset) { + int k; + int foundReset = FALSE; + int m = subcatalog[0].average[i].measureOffset; + MeasureTiny *measure = &subcatalog[0].measureT[m]; + for (j = 0; (j < subcatalog[0].average[i].Nmeasure) && !foundReset; j++) { + if (CHECK_TEST_IMAGE && (abs(measure[j].t - Tref) < 10) && (measure[j].photcode == Cref)) { + fprintf (stderr, "."); + } + for (k = 0; (k < NphotcodesReset) && !foundReset; k++) { + if (photcodesReset[k][0].code == measure[j].photcode) foundReset = TRUE; + } + } + if (!foundReset) { + Nm_grid[ix][iy][Nn] = 0; + } + } + + idxgrid[ix][iy][Nn] = i; + Nn_grid[ix][iy] ++; + if (Nn_grid[ix][iy] >= NN_grid[ix][iy]) { + NN_grid[ix][iy] += 100; + REALLOCATE (Nm_grid[ix][iy], int, NN_grid[ix][iy]); + REALLOCATE (idxgrid[ix][iy], off_t, NN_grid[ix][iy]); + } + } + + // sort all of the grid cells + for (ix = 0; ix < NX; ix++) { + for (iy = 0; iy < NY; iy++) { + sort_by_Nmeasure (Nm_grid[ix][iy], idxgrid[ix][iy], Nn_grid[ix][iy]); + NN_grid[ix][iy] = 0; // I'm going to use this array to track which element I've already selected + } + } + + // cycle over the grid until we ready Nmax + off_t *keepidx = NULL; + ALLOCATE (keepidx, off_t, Naverage); + memset (keepidx, 0, Naverage*sizeof(off_t)); + int Nkeep = 0; + + for (i = 0; (i < 20) && (Nkeep < Nmax); i++) { + for (ix = 0; (ix < NX) && (Nkeep < Nmax); ix++) { + for (iy = 0; (iy < NY) && (Nkeep < Nmax); iy++) { + if (NN_grid[ix][iy] >= Nn_grid[ix][iy]) continue; // all used up! + int Nn = NN_grid[ix][iy]; + keepidx[Nkeep] = idxgrid[ix][iy][Nn]; + Nkeep ++; + NN_grid[ix][iy] ++; + } + } + } + + // count the number of measurements this selection will yield + off_t ave, NMEASURE = 0; + for (i = 0; i < Nkeep; i++) { + ave = keepidx[i]; + NMEASURE += subcatalog[0].average[ave].Nmeasure; + } + + // test catID : 37262 37261 37257 37258 + int dumpit = FALSE; + dumpit |= (oldcatalog[0].catID == 37007); + // dumpit |= (oldcatalog[0].catID == 37261); + // dumpit |= (oldcatalog[0].catID == 37257); + // dumpit |= (oldcatalog[0].catID == 37258); + if (dumpit) { + char name[64]; + snprintf (name, 64, "cat.%05d.dump.dat", oldcatalog[0].catID); + FILE *fdump = fopen (name, "w"); + for (i = 0; i < Nkeep; i++) { + ave = keepidx[i]; + fprintf (fdump, "%10.6f %10.6f %d\n", subcatalog[0].average[ave].R, subcatalog[0].average[ave].D, subcatalog[0].average[ave].Nmeasure); + } + fclose (fdump); + } + + // allocate the output data + ALLOCATE (tmpcatalog.average, Average, Nkeep); + ALLOCATE (tmpcatalog.measureT, MeasureTiny, NMEASURE); + ALLOCATE (tmpcatalog.secfilt, SecFilt, Nkeep * Nsecfilt); + + off_t Nmeasure = 0; + + // copy the Nkeep selected entries from subcatalog to tmpcatalog (adjusting links) + for (i = 0; i < Nkeep; i++) { + ave = keepidx[i]; + tmpcatalog.average[i] = subcatalog[0].average[ave]; + tmpcatalog.average[i].measureOffset = Nmeasure; + for (j = 0; j < tmpcatalog.average[i].Nmeasure; j++) { + off_t offset = subcatalog[0].average[ave].measureOffset + j; + tmpcatalog.measureT[Nmeasure] = subcatalog[0].measureT[offset]; + tmpcatalog.measureT[Nmeasure].averef = i; + Nmeasure ++; + } + for (j = 0; j < Nsecfilt; j++) { + tmpcatalog.secfilt[i*Nsecfilt + j] = subcatalog[0].secfilt[ave*Nsecfilt + j]; + } + } + + if (VERBOSE2) { + char *basename = filebasename (oldcatalog[0].filename); + fprintf (stderr, "limited to %d ("OFF_T_FMT" subset, "OFF_T_FMT" total) stars, "OFF_T_FMT" ("OFF_T_FMT" subset, "OFF_T_FMT" total) measures for catalog %s\n", + Nkeep, subcatalog[0].Naverage, oldcatalog[0].Naverage, Nmeasure, subcatalog[0].Nmeasure, oldcatalog[0].Nmeasure, basename); + free (basename); + } + + free (X); + free (Y); + + for (ix = 0; ix < NX; ix++) { + for (iy = 0; iy < NY; iy++) { + free (Nm_grid[ix][iy]); + free (idxgrid[ix][iy]); + } + free (NN_grid[ix]); + free (Nn_grid[ix]); + free (Nm_grid[ix]); + free (idxgrid[ix]); + } + free (NN_grid); + free (Nn_grid); + free (Nm_grid); + free (idxgrid); + + free (keepidx); + + free (subcatalog[0].average); + free (subcatalog[0].measureT); + free (subcatalog[0].secfilt); + + subcatalog[0].average = tmpcatalog.average; + subcatalog[0].measureT = tmpcatalog.measureT; + subcatalog[0].secfilt = tmpcatalog.secfilt; + subcatalog[0].Naverage = Nkeep; + subcatalog[0].Nmeasure = Nmeasure; + subcatalog[0].Nsecfilt = oldcatalog[0].Nsecfilt; + subcatalog[0].Nsecf_mem = Naverage * oldcatalog[0].Nsecfilt; + + return (TRUE); } @@ -205,4 +543,15 @@ NMEASURE += subcatalog[0].average[ave].Nmeasure; } + +# if (0) + if (oldcatalog[0].catID == 59962) { + FILE *fdump = fopen ("cat.dump.dat", "w"); + for (i = 0; i < Nmax; i++) { + ave = index[i]; + fprintf (fdump, "%10.6f %10.6f %d\n", subcatalog[0].average[ave].R, subcatalog[0].average[ave].D, subcatalog[0].average[ave].Nmeasure); + } + fclose (fdump); + } +# endif // allocate the output data Index: trunk/Ohana/src/relastro/src/dvo_astrom_ops.c =================================================================== --- trunk/Ohana/src/relastro/src/dvo_astrom_ops.c (revision 36630) +++ trunk/Ohana/src/relastro/src/dvo_astrom_ops.c (revision 36833) @@ -91,4 +91,44 @@ } +double getMeanR_Big (Measure *measure, Average *average, SecFilt *secfilt) { + + double ra; + + /* the measure carries the instantaneous mean position at the epoch t */ + ra = average[0].R - measure[0].dR / 3600.0; + + /* possible corrections to mean ra: + + - proper-motion and parallax + - abberation + - precession and nutation, etc + - refraction + - DCR + + */ + + return (ra); +} + +double getMeanD_Big (Measure *measure, Average *average, SecFilt *secfilt) { + + double dec; + + /* the measure carries the instantaneous mean position at the epoch t */ + dec = average[0].D - measure[0].dD / 3600.0; + + /* possible corrections to mean ra: + + - proper-motion and parallax + - abberation + - precession and nutation, etc + - refraction + - DCR + + */ + + return (dec); +} + int setMeanR_Big (double ra_fit, Measure *measure, Average *average, SecFilt *secfilt) { @@ -144,2 +184,3 @@ return (TRUE); } + Index: trunk/Ohana/src/relastro/src/initialize.c =================================================================== --- trunk/Ohana/src/relastro/src/initialize.c (revision 36630) +++ trunk/Ohana/src/relastro/src/initialize.c (revision 36833) @@ -8,13 +8,36 @@ if (RELASTRO_OP == OP_MERGE_SOURCE) return; - fprintf (stderr, "PHOTCODE_KEEP_LIST: %s\n", PHOTCODE_KEEP_LIST); - fprintf (stderr, "PHOTCODE_SKIP_LIST: %s\n", PHOTCODE_SKIP_LIST); - fprintf (stderr, "PHOTCODE_A_LIST: %s\n", PHOTCODE_A_LIST); - fprintf (stderr, "PHOTCODE_B_LIST: %s\n", PHOTCODE_B_LIST); + if (DCR_COLOR_POS) fprintf (stderr, "DCR_COLOR_POS: %s\n", DCR_COLOR_POS); + if (DCR_COLOR_NEG) fprintf (stderr, "DCR_COLOR_NEG: %s\n", DCR_COLOR_NEG); - photcodesKeep = ParsePhotcodeList (PHOTCODE_KEEP_LIST, &NphotcodesKeep, FALSE); - photcodesSkip = ParsePhotcodeList (PHOTCODE_SKIP_LIST, &NphotcodesSkip, FALSE); - photcodesGroupA = ParsePhotcodeList (PHOTCODE_A_LIST, &NphotcodesGroupA, TRUE); - photcodesGroupB = ParsePhotcodeList (PHOTCODE_B_LIST, &NphotcodesGroupB, TRUE); + if (PHOTCODE_KEEP_LIST) fprintf (stderr, "PHOTCODE_KEEP_LIST: %s\n", PHOTCODE_KEEP_LIST); + if (PHOTCODE_SKIP_LIST) fprintf (stderr, "PHOTCODE_SKIP_LIST: %s\n", PHOTCODE_SKIP_LIST); + if (PHOTCODE_RESET_LIST) fprintf (stderr, "PHOTCODE_RESET_LIST: %s\n", PHOTCODE_RESET_LIST); + if (PHOTCODE_A_LIST) fprintf (stderr, "PHOTCODE_A_LIST: %s\n", PHOTCODE_A_LIST); + if (PHOTCODE_B_LIST) fprintf (stderr, "PHOTCODE_B_LIST: %s\n", PHOTCODE_B_LIST); + + photcodesKeep = ParsePhotcodeList (PHOTCODE_KEEP_LIST, &NphotcodesKeep, FALSE); + photcodesSkip = ParsePhotcodeList (PHOTCODE_SKIP_LIST, &NphotcodesSkip, FALSE); + photcodesReset = ParsePhotcodeList (PHOTCODE_RESET_LIST, &NphotcodesReset, FALSE); + photcodesGroupA = ParsePhotcodeList (PHOTCODE_A_LIST, &NphotcodesGroupA, TRUE); + photcodesGroupB = ParsePhotcodeList (PHOTCODE_B_LIST, &NphotcodesGroupB, TRUE); + + DCR_NSEC_POS = DCR_NSEC_NEG = -1; + if (DCR_COLOR_POS) { + DCR_PHOTCODE_POS = GetPhotcodebyName (DCR_COLOR_POS); + if (!DCR_PHOTCODE_POS) { + fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", DCR_COLOR_POS); + exit (1); + } + DCR_NSEC_POS = GetPhotcodeNsec (DCR_PHOTCODE_POS[0].code); + } + if (DCR_COLOR_NEG) { + DCR_PHOTCODE_NEG = GetPhotcodebyName (DCR_COLOR_NEG); + if (!DCR_PHOTCODE_NEG) { + fprintf (stderr, "ERROR: photcode %s not found in photcode table\n", DCR_COLOR_NEG); + exit (1); + } + DCR_NSEC_NEG = GetPhotcodeNsec (DCR_PHOTCODE_NEG[0].code); + } initstats (STATMODE); @@ -42,13 +65,15 @@ args_client (argc, argv); - fprintf (stderr, "PHOTCODE_KEEP_LIST: %s\n", PHOTCODE_KEEP_LIST); - fprintf (stderr, "PHOTCODE_SKIP_LIST: %s\n", PHOTCODE_SKIP_LIST); - fprintf (stderr, "PHOTCODE_A_LIST: %s\n", PHOTCODE_A_LIST); - fprintf (stderr, "PHOTCODE_B_LIST: %s\n", PHOTCODE_B_LIST); + if (PHOTCODE_KEEP_LIST) fprintf (stderr, "PHOTCODE_KEEP_LIST: %s\n", PHOTCODE_KEEP_LIST); + if (PHOTCODE_SKIP_LIST) fprintf (stderr, "PHOTCODE_SKIP_LIST: %s\n", PHOTCODE_SKIP_LIST); + if (PHOTCODE_RESET_LIST) fprintf (stderr, "PHOTCODE_RESET_LIST: %s\n", PHOTCODE_RESET_LIST); + if (PHOTCODE_A_LIST) fprintf (stderr, "PHOTCODE_A_LIST: %s\n", PHOTCODE_A_LIST); + if (PHOTCODE_B_LIST) fprintf (stderr, "PHOTCODE_B_LIST: %s\n", PHOTCODE_B_LIST); - photcodesKeep = ParsePhotcodeList (PHOTCODE_KEEP_LIST, &NphotcodesKeep, FALSE); - photcodesSkip = ParsePhotcodeList (PHOTCODE_SKIP_LIST, &NphotcodesSkip, FALSE); - photcodesGroupA = ParsePhotcodeList (PHOTCODE_A_LIST, &NphotcodesGroupA, TRUE); - photcodesGroupB = ParsePhotcodeList (PHOTCODE_B_LIST, &NphotcodesGroupB, TRUE); + photcodesKeep = ParsePhotcodeList (PHOTCODE_KEEP_LIST, &NphotcodesKeep, FALSE); + photcodesSkip = ParsePhotcodeList (PHOTCODE_SKIP_LIST, &NphotcodesSkip, FALSE); + photcodesReset = ParsePhotcodeList (PHOTCODE_RESET_LIST, &NphotcodesReset, FALSE); + photcodesGroupA = ParsePhotcodeList (PHOTCODE_A_LIST, &NphotcodesGroupA, TRUE); + photcodesGroupB = ParsePhotcodeList (PHOTCODE_B_LIST, &NphotcodesGroupB, TRUE); initstats (STATMODE); Index: trunk/Ohana/src/relastro/src/launch_region_hosts.c =================================================================== --- trunk/Ohana/src/relastro/src/launch_region_hosts.c (revision 36630) +++ trunk/Ohana/src/relastro/src/launch_region_hosts.c (revision 36833) @@ -1,17 +1,4 @@ # include "relastro.h" # define DEBUG 0 - -int strextend (char *input, char *format,...) { - - char tmpextra[1024], tmpline[1024]; - va_list argp; - - va_start (argp, format); - vsnprintf (tmpextra, 1024, format, argp); - snprintf (tmpline, 1024, "%s %s", input, tmpextra); - strcpy (input, tmpline); - - return TRUE; -} int launch_region_hosts (RegionHostTable *regionHosts) { @@ -88,7 +75,14 @@ if (PHOTCODE_KEEP_LIST) strextend (command, "+photcode %s", PHOTCODE_KEEP_LIST); if (PHOTCODE_SKIP_LIST) strextend (command, "-photcode %s", PHOTCODE_SKIP_LIST); + if (PHOTCODE_RESET_LIST) strextend (command, "-reset-to-photcode %s", PHOTCODE_RESET_LIST); if (MaxDensityUse) strextend (command, "-max-density %f", MaxDensityValue); if (ImagSelect) strextend (command, "-instmag %f %f", ImagMin, ImagMax); + if (ExcludeBogus) strextend (command, "-exclude-bogus %f", ExcludeBogusRadius); + + if (DCR_COLOR_POS && DCR_COLOR_NEG) { + strextend (command, "-dcr-color %s %s", DCR_COLOR_POS, DCR_COLOR_NEG); + } + if (PhotFlagSelect) strextend (command, "+photflags"); if (PhotFlagBad) strextend (command, "+photflagbad %d", PhotFlagBad); @@ -97,4 +91,11 @@ if (MinBadQF > 0.0) strextend (command, "-min-bad-psfqf %f", MinBadQF); if (MaxMeanOffset != 10.0) strextend (command, "-max-mean-offset %f", MaxMeanOffset); + + strextend (command, "-D RELASTRO_SRC_MEAS_TOOFEW %d", SRC_MEAS_TOOFEW); + strextend (command, "-D RELASTRO_SIGMA_LIM %f", SIGMA_LIM); + strextend (command, "-D RELASTRO_DPOS_MAX %f", DPOS_MAX); + strextend (command, "-D ADDSTAR_RADIUS %f", ADDSTAR_RADIUS); + + if (USE_FIXED_PIXCOORDS) strextend (command, "-D USE_FIXED_PIXCOORDS 1"); if (TimeSelect) { Index: trunk/Ohana/src/relastro/src/load_catalogs.c =================================================================== --- trunk/Ohana/src/relastro/src/load_catalogs.c (revision 36630) +++ trunk/Ohana/src/relastro/src/load_catalogs.c (revision 36833) @@ -1,15 +1,16 @@ # include "relastro.h" -Catalog *load_catalogs_parallel (SkyList *sky, int *Ncatalog); +Catalog *load_catalogs_parallel (SkyList *sky, int *Ncatalog, char *syncfile); void bcatalog_show_skips (); -Catalog *load_catalogs (SkyList *skylist, int *Ncatalog, int subselect, int hostID, char *hostpath) { - - int i, j, k, m, Nstar; +Catalog *load_catalogs (SkyList *skylist, int *Ncatalog, int subselect, int hostID, char *hostpath, char *syncfile) { + + int i, j, Nstar; + // int k, m; Catalog *catalog, *pcatalog, tcatalog; // XXX need to decide how to determine PARALLEL mode... if (PARALLEL && !hostID) { - catalog = load_catalogs_parallel (skylist, Ncatalog); + catalog = load_catalogs_parallel (skylist, Ncatalog, syncfile); return catalog; } @@ -73,4 +74,5 @@ if (RESET) { for (j = 0; j < catalog[Ncat].Naverage; j++) { +# if (0) catalog[Ncat].average[j].flags = 0; m = catalog[Ncat].average[j].measureOffset; @@ -78,4 +80,5 @@ catalog[Ncat].measure[m+k].dbFlags = 0; } +# endif } } @@ -84,5 +87,5 @@ } - // XXX TEST : bcatalog_show_skips(); + bcatalog_show_skips(); Nstar = 0; @@ -113,5 +116,5 @@ // CATDIR is supplied globally # define DEBUG 1 -Catalog *load_catalogs_parallel (SkyList *sky, int *Ncatalog) { +Catalog *load_catalogs_parallel (SkyList *sky, int *Ncatalog, char *syncfile) { char uniquer[12]; @@ -150,24 +153,26 @@ table->hosts[i].results, table->hosts[i].hostID, CATDIR, table->hosts[i].pathname, UserPatch.Rmin, UserPatch.Rmax, UserPatch.Dmin, UserPatch.Dmax, STATMODE, MIN_ERROR, SIGMA_LIM); - char tmpline[1024]; - if (FIT_MODE == FIT_PM_ONLY) { snprintf (tmpline, 1024, "%s -pm", command); strcpy (command, tmpline); } - if (FIT_MODE == FIT_PAR_ONLY) { snprintf (tmpline, 1024, "%s -par", command); strcpy (command, tmpline); } - if (FIT_MODE == FIT_PM_AND_PAR) { snprintf (tmpline, 1024, "%s -pmpar", command); strcpy (command, tmpline); } - - if (VERBOSE) { snprintf (tmpline, 1024, "%s -v", command); strcpy (command, tmpline); } - if (VERBOSE2) { snprintf (tmpline, 1024, "%s -vv", command); strcpy (command, tmpline); } - if (RESET) { snprintf (tmpline, 1024, "%s -reset", command); strcpy (command, tmpline); } - if (ImagSelect) { snprintf (tmpline, 1024, "%s -instmag %f %f", command, ImagMin, ImagMax); strcpy (command, tmpline); } - if (MaxDensityUse) { snprintf (tmpline, 1024, "%s -max-density %f", command, MaxDensityValue); strcpy (command, tmpline); } + if (FIT_MODE == FIT_PM_ONLY) strextend (command, "-pm"); + if (FIT_MODE == FIT_PAR_ONLY) strextend (command, "-par"); + if (FIT_MODE == FIT_PM_AND_PAR) strextend (command, "-pmpar"); + + if (VERBOSE) strextend (command, "-v"); + if (VERBOSE2) strextend (command, "-vv"); + if (RESET) strextend (command, "-reset"); + if (ImagSelect) strextend (command, "-instmag %f %f", ImagMin, ImagMax); + if (MaxDensityUse) strextend (command, "-max-density %f", MaxDensityValue); + if (FlagOutlier) strextend (command, "-clip %d", CLIP_THRESH); + if (ExcludeBogus) strextend (command, "-exclude-bogus %f", ExcludeBogusRadius); - if (FlagOutlier) { snprintf (tmpline, 1024, "%s -clip %d", command, CLIP_THRESH); strcpy (command, tmpline); } - - if (USE_FIXED_PIXCOORDS) { snprintf (tmpline, 1024, "%s -D USE_FIXED_PIXCOORDS 1", command); strcpy (command, tmpline); } - - if (PHOTCODE_KEEP_LIST) { snprintf (tmpline, 1024, "%s +photcode %s", command, PHOTCODE_KEEP_LIST); strcpy (command, tmpline); } - if (PHOTCODE_SKIP_LIST) { snprintf (tmpline, 1024, "%s -photcode %s", command, PHOTCODE_SKIP_LIST); strcpy (command, tmpline); } - if (PhotFlagSelect) { snprintf (tmpline, 1024, "%s +photflags", command); strcpy (command, tmpline); } - if (PhotFlagBad) { snprintf (tmpline, 1024, "%s +photflagbad %d", command, PhotFlagBad); strcpy (command, tmpline); } - if (PhotFlagPoor) { snprintf (tmpline, 1024, "%s +photflagpoor %d", command, PhotFlagPoor); strcpy (command, tmpline); } + if (DCR_COLOR_POS && DCR_COLOR_NEG) { + strextend (command, "-dcr-color %s %s", DCR_COLOR_POS, DCR_COLOR_NEG); + } + + if (USE_FIXED_PIXCOORDS) strextend (command, "-D USE_FIXED_PIXCOORDS 1"); + if (PHOTCODE_KEEP_LIST) strextend (command, "+photcode %s", PHOTCODE_KEEP_LIST); + if (PHOTCODE_SKIP_LIST) strextend (command, "-photcode %s", PHOTCODE_SKIP_LIST); + if (PhotFlagSelect) strextend (command, "+photflags"); + if (PhotFlagBad) strextend (command, "+photflagbad %d", PhotFlagBad); + if (PhotFlagPoor) strextend (command, "+photflagpoor %d", PhotFlagPoor); // XXX note that the above pass in the flag as decimal -- also note that args.c cannot handle 0xHEX values @@ -175,8 +180,7 @@ char *tstart = ohana_sec_to_date (TSTART); char *tstop = ohana_sec_to_date (TSTOP); - snprintf (tmpline, 1024, "%s -time %s %s", command, tstart, tstop); + strextend (command, "-time %s %s", tstart, tstop); free (tstart); free (tstop); - strcpy (command, tmpline); } @@ -211,4 +215,13 @@ } + // update syncfile here (save lots of I/O time) + + // NOTE: if I let all hosts load blindly, I saturate the data clients with too many + // relastro_client requests. I need to have the master mediate this. the master + // will not launch the next remote job until this one says it is done + if (syncfile) { + update_sync_file (syncfile, 1); + } + // each host generates a BrightCatalog structure, with the measure, average, etc value // loaded into a single set of arrays (of MeasureTiny, AverageTiny, Secfilt). I need to Index: trunk/Ohana/src/relastro/src/relastro_client.c =================================================================== --- trunk/Ohana/src/relastro/src/relastro_client.c (revision 36630) +++ trunk/Ohana/src/relastro/src/relastro_client.c (revision 36833) @@ -33,5 +33,5 @@ // USAGE: relastro_client -load-objects int Ncatalog; - Catalog *catalog = load_catalogs (skylist, &Ncatalog, TRUE, HOST_ID, HOSTDIR); + Catalog *catalog = load_catalogs (skylist, &Ncatalog, TRUE, HOST_ID, HOSTDIR, NULL); if (!catalog) { fprintf (stderr, "ERROR loading catalogs from %s\n", CATDIR); Index: trunk/Ohana/src/relastro/src/relastro_images.c =================================================================== --- trunk/Ohana/src/relastro/src/relastro_images.c (revision 36630) +++ trunk/Ohana/src/relastro/src/relastro_images.c (revision 36833) @@ -28,6 +28,14 @@ // XXX pass in the image table // XXX who carries the image grid? - catalog = load_catalogs (skylist, &Ncatalog, TRUE, 0, NULL); + + // photcodesKeep is used here to allow measurements from the images being calibrated + // note if -reset-to-photcode is selected, photocodesKeep is replaced with below with photcodesReset + catalog = load_catalogs (skylist, &Ncatalog, TRUE, 0, NULL, NULL); MARKTIME("load catalog data: %f sec\n", dtime); + + if (photcodesReset) { + photcodesKeep = photcodesReset; + NphotcodesKeep = NphotcodesReset; + } if (Ncatalog == 0) { @@ -50,4 +58,6 @@ // set test points based on the starmap createStarMap (catalog, Ncatalog); + + // XXX NOTE : for 2mass reset, photcodesKeep should now limit to 2MASS measurements /* major modes */ Index: trunk/Ohana/src/relastro/src/relastro_objects.c =================================================================== --- trunk/Ohana/src/relastro/src/relastro_objects.c (revision 36630) +++ trunk/Ohana/src/relastro/src/relastro_objects.c (revision 36833) @@ -40,4 +40,6 @@ snprintf (hostfile, 1024, "%s/%s.cpt", hostpath, skylist[0].regions[i]->name); catalog.filename = hostID ? hostfile : skylist[0].filename[i]; + + // set up the basic catalog info 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 Index: trunk/Ohana/src/relastro/src/relastro_parallel_images.c =================================================================== --- trunk/Ohana/src/relastro/src/relastro_parallel_images.c (revision 36630) +++ trunk/Ohana/src/relastro/src/relastro_parallel_images.c (revision 36833) @@ -45,13 +45,14 @@ /* load catalog data from region files (hostID is 0 since we are not a client */ - catalog = load_catalogs (skylist, &Ncatalog, TRUE, 0, NULL); + char *syncfile = make_filename (CATDIR, regionHosts->hosts[myHost].hostname, REGION_HOST_ID, "loadcat.sync"); + catalog = load_catalogs (skylist, &Ncatalog, TRUE, 0, NULL, syncfile); MARKTIME("-- load catalog data: %f sec\n", dtime); + free (syncfile); - // NOTE: if I let all hosts load blindly, I saturate the data clients with too many - // relastro_client requests. I need to have the master mediate this. the master - // will not launch the next remote job until this one says it is done - char *syncfile = make_filename (CATDIR, regionHosts->hosts[myHost].hostname, REGION_HOST_ID, "loadcat.sync"); - update_sync_file (syncfile, 1); - + if (photcodesReset) { + photcodesKeep = photcodesReset; + NphotcodesKeep = NphotcodesReset; + } + // generate tables go from catID,objID -> catSeq,objSeq indexCatalogs (catalog, Ncatalog);