78 centers.
push_back(ccdImage->getBoresightRaDec());
86 _commonTangentPoint =
Point(commonTangentPoint.getX(), commonTangentPoint.getY());
87 for (
auto &ccdImage :
ccdImageList) ccdImage->setCommonTangentPoint(_commonTangentPoint);
92 Frame tangentPlaneFrame;
95 Frame CTPFrame = ccdImage->getPixelToCommonTangentPlane()->apply(ccdImage->getImageFrame(),
false);
96 if (tangentPlaneFrame.
getArea() == 0)
97 tangentPlaneFrame = CTPFrame;
99 tangentPlaneFrame += CTPFrame;
104 TanPixelToRaDec commonTangentPlaneToRaDec(identity, _commonTangentPoint);
105 Frame raDecFrame = commonTangentPlaneToRaDec.
apply(tangentPlaneFrame,
false);
120 const bool enlargeFittedList) {
126 item->clearBeforeAssoc();
136 ccdImage->resetCatalogForFit();
142 Frame ccdImageFrameCPT = toCommonTangentPlane->apply(ccdImage->getImageFrame(),
false);
143 ccdImageFrameCPT = ccdImageFrameCPT.
rescale(1.10);
149 if (ccdImageFrameCPT.
inFrame(*fittedStar)) {
156 toCommonTangentPlane.
get(), matchCutInArcSec / 3600.);
159 LOGLS_DEBUG(_log,
"Measured-to-Fitted matches before removing ambiguities " << starMatchList->size());
160 starMatchList->removeAmbiguities(*toCommonTangentPlane);
161 LOGLS_DEBUG(_log,
"Measured-to-Fitted matches after removing ambiguities " << starMatchList->size());
165 int matchedCount = 0;
166 for (
auto const &starMatch : *starMatchList) {
167 auto bs = starMatch.s1;
168 auto ms_const = std::dynamic_pointer_cast<const MeasuredStar>(bs);
169 auto ms = std::const_pointer_cast<MeasuredStar>(ms_const);
170 auto bs2 = starMatch.s2;
171 auto fs_const = std::dynamic_pointer_cast<const FittedStar>(bs2);
172 auto fs = std::const_pointer_cast<FittedStar>(fs_const);
173 ms->setFittedStar(fs);
176 LOGLS_INFO(_log,
"Matched " << matchedCount <<
" objects in " << ccdImage->getName());
179 int unMatchedCount = 0;
180 for (
auto const &mstar : catalog) {
182 if (mstar->getFittedStar())
continue;
183 if (enlargeFittedList) {
184 auto fs = std::make_shared<FittedStar>(*mstar);
186 toCommonTangentPlane->transformPosAndErrors(*fs, *fs);
188 mstar->setFittedStar(fs);
192 LOGLS_INFO(_log,
"Unmatched objects: " << unMatchedCount);
204 std::string const &fluxField,
float refCoordinateErr,
205 bool rejectBadFluxes) {
206 if (refCat.
size() == 0) {
208 " reference catalog is empty : stop here "));
216 raErrKey = refCat.
getSchema()[
"coord_raErr"];
217 decErrKey = refCat.
getSchema()[
"coord_decErr"];
227 << fluxField <<
"Err"
228 <<
") not found in reference catalog. Not using ref flux errors.");
232 for (
size_t i = 0; i < refCat.
size(); i++) {
233 auto const &record = refCat.
get(i);
235 auto coord = record->
get(coordKey);
236 double flux = record->get(fluxKey);
238 if (fluxErrKey.isValid()) {
239 fluxErr = record->get(fluxErrKey);
245 auto star = std::make_shared<RefStar>(ra,
dec,
flux, fluxErr);
253 star->vx =
std::pow(refCoordinateErr / 1000. / 3600. /
std::cos(coord.getLatitude()), 2);
254 star->vy =
std::pow(refCoordinateErr / 1000. / 3600., 2);
266 TanRaDecToPixel raDecToCommonTangentPlane(identity, _commonTangentPoint);
268 associateRefStars(matchCut.
asArcseconds(), &raDecToCommonTangentPlane);
275 matchCutInArcSec / 3600.);
277 LOGLS_DEBUG(_log,
"Refcat matches before removing ambiguities " << starMatchList->size());
278 starMatchList->removeAmbiguities(*
transform);
279 LOGLS_DEBUG(_log,
"Refcat matches after removing ambiguities " << starMatchList->size());
282 for (
auto const &starMatch : *starMatchList) {
286 const BaseStar &bs2 = *starMatch.s2;
294 "Associated " << starMatchList->size() <<
" reference stars among " <<
refStarList.
size());
298 selectFittedStars(minMeasurements);
299 normalizeFittedStars();
302 void Associations::selectFittedStars(
int minMeasurements) {
305 int totalMeasured = 0, validMeasured = 0;
317 if (fittedStar ==
nullptr) {
324 if (!fittedStar->getRefStar() && fittedStar->getMeasurementCount() < minMeasurements) {
325 fittedStar->getMeasurementCount()--;
326 mi = catalog.
erase(mi);
336 if ((*fi)->getMeasurementCount() == 0) {
344 LOGLS_INFO(_log,
"Total, valid number of Measured stars: " << totalMeasured <<
", " << validMeasured);
347 void Associations::normalizeFittedStars()
const {
352 fittedStar->setFlux(0.0);
353 fittedStar->getMag() = 0.0;
359 MeasuredStarList &catalog = ccdImage->getCatalogForFit();
360 for (
auto &mi : catalog) {
361 auto fittedStar = mi->getFittedStar();
362 if (fittedStar ==
nullptr)
364 pex::exceptions::RuntimeError,
365 "All measuredStars must have a fittedStar: did you call selectFittedStars()?"));
366 auto point = toCommonTangentPlane->apply(*mi);
367 fittedStar->x += point.x;
368 fittedStar->y += point.y;
369 fittedStar->getFlux() += mi->getFlux();
374 auto measurementCount = fi->getMeasurementCount();
375 fi->x /= measurementCount;
376 fi->y /= measurementCount;
377 fi->getFlux() /= measurementCount;
382 void Associations::assignMags() {
384 MeasuredStarList &catalog = ccdImage->getCatalogForFit();
385 for (
auto const &mstar : catalog) {
387 if (!fstar)
continue;
398 "DeprojectFittedStars: Fitted stars are already in sidereal coordinates, nothing done ");
409 return item->getCatalogForFit().size() > 0;
416 if ((fittedStar !=
nullptr) & (fittedStar->getRefStar() !=
nullptr))
count++;
422 void Associations::collectMCStars(
int realization) {
428 string dbimdir = ccdImage.Dir();
429 string mctruth = dbimdir +
"/mc/mctruth.list";
431 if (realization >= 0) {
433 sstrm << dbimdir <<
"/mc/mctruth_" << realization <<
".list";
434 mctruth = sstrm.
str();
437 AstrometryTransformIdentity gti;
441 DicStarList mctruthlist(mctruth);
445 for (smI = starMatchList->begin(); smI != starMatchList->end(); smI++) {
446 StarMatch &sm = *smI;
447 BaseStar *bs = sm.s1;
448 MeasuredStar *mstar =
dynamic_cast<MeasuredStar *
>(bs);
450 DicStar *dstar =
dynamic_cast<DicStar *
>(bs);
452 mcstar->GetMCInfo().iflux = dstar->getval(
"iflux");
453 mcstar->GetMCInfo().tflux = dstar->getval(
"sflux");
460 LOGLS_FATAL(_log,
"CollectMCStars Unable to match MCTruth w/ catalog!");
465 double matchCutArcSec) {
467 size_t pos_minus = color.
find(
'-');
468 bool compute_diff = (pos_minus != string::npos);
470 c1 = color.
substr(0, pos_minus);
471 if (compute_diff) c2 = color.
substr(pos_minus + 1, string::npos);
472 DicStarList cList(dicStarListName);
473 if (!cList.HasKey(c1))
474 throw(GastroException(
"Associations::SetFittedstarColors : " + dicStarListName +
475 " misses a key named \"" + c1 +
"\""));
476 if (compute_diff && !cList.HasKey(c2))
477 throw(GastroException(
"Associations::SetFittedstarColors : " + dicStarListName +
478 " misses a key named \"" + c2 +
"\""));
482 AstrometryTransformIdentity
id;
485 AstrometryTransform *id_or_proj = &proj;
491 id_or_proj, matchCutArcSec / 3600);
494 <<
" FittedStars to color catalog");
496 for (
auto i = starMatchList->begin(); i != starMatchList->end(); ++i) {
497 BaseStar *s1 = i->s1;
498 FittedStar *fs =
dynamic_cast<FittedStar *
>(s1);
499 BaseStar *s2 = i->s2;
500 const TStar *ts =
dynamic_cast<const TStar *
>(s2);
501 const DicStar *ds =
dynamic_cast<const DicStar *
>(ts->get_original());
502 fs->color = ds->getval(c1);
503 if (compute_diff) fs->color -= ds->getval(c2);
This file declares a class for representing circular regions on the unit sphere.
This file declares a class for representing convex polygons with great circle edges on the unit spher...
table::Key< int > detector
#define LSST_EXCEPT(type,...)
Create an exception with a given type.
Combinatorial searches for linear transformations to go from list1 to list2.
LSST DM logging module built on log4cxx.
#define LOGLS_WARN(logger, message)
Log a warn-level message using an iostream-based interface.
#define LOGLS_FATAL(logger, message)
Log a fatal-level message using an iostream-based interface.
#define LOG_GET(logger)
Returns a Log object associated with logger.
#define LOGLS_INFO(logger, message)
Log a info-level message using an iostream-based interface.
#define LOGLS_DEBUG(logger, message)
Log a debug-level message using an iostream-based interface.
This file contains a class representing spherical coordinates.
table::Key< table::Array< std::uint8_t > > wcs
size_type size() const
Return the number of elements in the catalog.
Schema getSchema() const
Return the schema associated with the catalog's table.
std::shared_ptr< RecordT > const get(size_type i) const
Return a pointer to the record at index i.
A FunctorKey used to get or set celestial coordinates from a pair of lsst::geom::Angle keys.
SchemaItem< T > find(std::string const &name) const
Find a SchemaItem in the Schema by name.
Custom catalog class for record/table subclasses that are guaranteed to have an ID,...
A class representing an angle.
constexpr double asArcseconds() const noexcept
Return an Angle's value in arcseconds.
An integer coordinate rectangle.
size_t nFittedStarsWithAssociatedRefStar() const
Return the number of fittedStars that have an associated refStar.
void computeCommonTangentPoint()
Sets a shared tangent point for all ccdImages, using the mean of the centers of all ccdImages.
lsst::sphgeom::Circle computeBoundingCircle() const
Return the bounding circle in on-sky (RA, Dec) coordinates containing all CcdImages.
int nCcdImagesValidForFit() const
return the number of CcdImages with non-empty catalogs to-be-fit.
void createCcdImage(afw::table::SourceCatalog &catalog, std::shared_ptr< lsst::afw::geom::SkyWcs > wcs, std::shared_ptr< lsst::afw::image::VisitInfo > visitInfo, lsst::geom::Box2I const &bbox, std::string const &filter, std::shared_ptr< afw::image::PhotoCalib > photoCalib, std::shared_ptr< afw::cameraGeom::Detector > detector, int visit, int ccd, lsst::jointcal::JointcalControl const &control)
Create a ccdImage from an exposure catalog and metadata, and add it to the list.
void setCommonTangentPoint(lsst::geom::Point2D const &commonTangentPoint)
Sets a shared tangent point for all ccdImages.
void associateCatalogs(const double matchCutInArcsec=0, const bool useFittedList=false, const bool enlargeFittedList=true)
incrementaly builds a merged catalog of all image catalogs
Point getCommonTangentPoint() const
can be used to project sidereal coordinates related to the image set on a plane.
void collectRefStars(afw::table::SimpleCatalog &refCat, geom::Angle matchCut, std::string const &fluxField, float refCoordinateErr, bool rejectBadFluxes=false)
Collect stars from an external reference catalog and associate them with fittedStars.
CcdImageList ccdImageList
FittedStarList fittedStarList
void prepareFittedStars(int minMeasurements)
Set the color field of FittedStar 's from a colored catalog.
void deprojectFittedStars()
Sends back the fitted stars coordinates on the sky FittedStarsList::inTangentPlaneCoordinates keeps t...
The base class for handling stars. Used by all matching routines.
Handler of an actual image from a single CCD.
MeasuredStarList const & getCatalogForFit() const
Gets the catalog to be used for fitting, which may have been cleaned-up.
The objects which have been measured several times.
void setRefStar(const RefStar *_refStar)
Set the astrometric reference star associated with this star.
A list of FittedStar s. Such a list is typically constructed by Associations.
bool inTangentPlaneCoordinates
rectangle with sides parallel to axes.
Frame rescale(const double factor) const
rescale it. The center does not move.
bool inFrame(double x, double y) const
inside?
double xMin
coordinate of boundary.
objects measured on actual images.
std::shared_ptr< FittedStar > getFittedStar() const
double getMagWeight() const
the inverse of the mag variance
A list of MeasuredStar. They are usually filled in Associations::createCcdImage.
Objects used as position anchors, typically USNO stars.
void applyTransform(const Operator &op)
enables to apply a geometrical transform if Star is Basestar or derives from it.
The transformation that handles pixels to sideral transformations (Gnomonic, possibly with polynomial...
virtual void apply(const double xIn, const double yIn, double &xOut, double &yOut) const=0
This one is the Tangent Plane (called gnomonic) projection (from celestial sphere to tangent plane)
Reports invalid arguments.
Reports attempts to access elements using an invalid key.
Circle is a circular region on the unit sphere that contains its boundary.
static ConvexPolygon convexHull(std::vector< UnitVector3d > const &points)
convexHull returns the convex hull of the given set of points if it exists and throws an exception ot...
Circle getBoundingCircle() const override
getBoundingCircle returns a bounding-circle for this region.
static LonLat fromDegrees(double lon, double lat)
T emplace_back(T... args)
SpherePoint averageSpherePoint(std::vector< SpherePoint > const &coords)
Return the average of a list of coordinates.
constexpr AngleUnit degrees
constant with units of degrees
constexpr double radToDeg(double x) noexcept
StarList< BaseStar > BaseStarList
BaseStarList & Measured2Base(MeasuredStarList &This)
std::unique_ptr< StarMatchList > listMatchCollect(const BaseStarList &list1, const BaseStarList &list2, const AstrometryTransform *guess, const double maxDist)
assembles star matches.
BaseStarList & Fitted2Base(FittedStarList &This)
MeasuredStarList::iterator MeasuredStarIterator
::std::list< StarMatch >::iterator StarMatchIterator
BaseStarList & Ref2Base(RefStarList &This)
FittedStarList::iterator FittedStarIterator
double nanojanskyToABMagnitude(double flux)
Convert a flux in nanojansky to AB magnitude.
A base class for image defects.
std::string sourceFluxField
"name of flux field in source catalog" ;