LSST Applications 27.0.0,g0265f82a02+469cd937ee,g02d81e74bb+21ad69e7e1,g1470d8bcf6+cbe83ee85a,g2079a07aa2+e67c6346a6,g212a7c68fe+04a9158687,g2305ad1205+94392ce272,g295015adf3+81dd352a9d,g2bbee38e9b+469cd937ee,g337abbeb29+469cd937ee,g3939d97d7f+72a9f7b576,g487adcacf7+71499e7cba,g50ff169b8f+5929b3527e,g52b1c1532d+a6fc98d2e7,g591dd9f2cf+df404f777f,g5a732f18d5+be83d3ecdb,g64a986408d+21ad69e7e1,g858d7b2824+21ad69e7e1,g8a8a8dda67+a6fc98d2e7,g99cad8db69+f62e5b0af5,g9ddcbc5298+d4bad12328,ga1e77700b3+9c366c4306,ga8c6da7877+71e4819109,gb0e22166c9+25ba2f69a1,gb6a65358fc+469cd937ee,gbb8dafda3b+69d3c0e320,gc07e1c2157+a98bf949bb,gc120e1dc64+615ec43309,gc28159a63d+469cd937ee,gcf0d15dbbd+72a9f7b576,gdaeeff99f8+a38ce5ea23,ge6526c86ff+3a7c1ac5f1,ge79ae78c31+469cd937ee,gee10cc3b42+a6fc98d2e7,gf1cff7945b+21ad69e7e1,gfbcc870c63+9a11dc8c8f
LSST Data Management Base Package
Loading...
Searching...
No Matches
Associations.cc
Go to the documentation of this file.
1// -*- LSST-C++ -*-
2/*
3 * This file is part of jointcal.
4 *
5 * Developed for the LSST Data Management System.
6 * This product includes software developed by the LSST Project
7 * (https://www.lsst.org).
8 * See the COPYRIGHT file at the top-level directory of this distribution
9 * for details of code ownership.
10 *
11 * This program is free software: you can redistribute it and/or modify
12 * it under the terms of the GNU General Public License as published by
13 * the Free Software Foundation, either version 3 of the License, or
14 * (at your option) any later version.
15 *
16 * This program is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 * GNU General Public License for more details.
20 *
21 * You should have received a copy of the GNU General Public License
22 * along with this program. If not, see <https://www.gnu.org/licenses/>.
23 */
24
25#include <cmath>
26#include <iostream>
27#include <limits>
28#include <sstream>
29
30#include "lsst/log/Log.h"
35#include "lsst/jointcal/Frame.h"
39
43
44#include "lsst/pex/exceptions.h"
45#include "lsst/geom/Box.h"
46#include "lsst/geom/Point.h"
48
49#include "lsst/sphgeom/LonLat.h"
50#include "lsst/sphgeom/Circle.h"
52
53namespace jointcal = lsst::jointcal;
54
55namespace {
56LOG_LOGGER _log = LOG_GET("lsst.jointcal.Associations");
57}
58
59namespace lsst {
60namespace jointcal {
61
65 lsst::geom::Box2I const &bbox, std::string const &filter,
67 const std::shared_ptr<afw::cameraGeom::Detector>& detector, int visit, int ccd,
68 lsst::jointcal::JointcalControl const &control) {
69 auto ccdImage = std::make_shared<CcdImage>(catalog, wcs, visitInfo, bbox, filter, photoCalib, detector,
70 visit, ccd, control.sourceFluxField);
71 ccdImageList.push_back(ccdImage);
72}
73
76 centers.reserve(ccdImageList.size());
77 for (auto const &ccdImage : ccdImageList) {
78 centers.push_back(ccdImage->getBoresightRaDec());
79 }
80 auto commonTangentPoint = geom::averageSpherePoint(centers);
81 LOGLS_DEBUG(_log, "Using common tangent point: " << commonTangentPoint.getPosition(geom::degrees));
82 setCommonTangentPoint(commonTangentPoint.getPosition(geom::degrees));
83}
84
86 _commonTangentPoint = Point(commonTangentPoint.getX(), commonTangentPoint.getY()); // a jointcal::Point
87 for (auto &ccdImage : ccdImageList) ccdImage->setCommonTangentPoint(_commonTangentPoint);
88}
89
91 // Compute the frame on the common tangent plane that contains all input images.
92 Frame tangentPlaneFrame;
93
94 for (auto const &ccdImage : ccdImageList) {
95 Frame CTPFrame = ccdImage->getPixelToCommonTangentPlane()->apply(ccdImage->getImageFrame(), false);
96 if (tangentPlaneFrame.getArea() == 0)
97 tangentPlaneFrame = CTPFrame;
98 else
99 tangentPlaneFrame += CTPFrame;
100 }
101
102 // Convert tangent plane coordinates to RaDec.
104 TanPixelToRaDec commonTangentPlaneToRaDec(identity, _commonTangentPoint);
105 Frame raDecFrame = commonTangentPlaneToRaDec.apply(tangentPlaneFrame, false);
106
108 points.reserve(4);
109 // raDecFrame is in on-sky (RA,Dec) degrees stored as an x/y box:
110 // the on-sky bounding box it represents is given by the corners of that box.
111 LOGLS_INFO(_log, "Computed tangent plane box for data (degrees): " << raDecFrame);
112 points.emplace_back(sphgeom::LonLat::fromDegrees(raDecFrame.xMin, raDecFrame.yMin));
113 points.emplace_back(sphgeom::LonLat::fromDegrees(raDecFrame.xMax, raDecFrame.yMin));
114 points.emplace_back(sphgeom::LonLat::fromDegrees(raDecFrame.xMin, raDecFrame.yMax));
115 points.emplace_back(sphgeom::LonLat::fromDegrees(raDecFrame.xMax, raDecFrame.yMax));
116
118}
119
120void Associations::associateCatalogs(const double matchCutInArcSec, const bool useFittedList,
121 const bool enlargeFittedList) {
122 // clear reference stars
124
125 // clear measurement counts and associations to refstars, but keep fittedStars themselves.
126 for (auto &item : fittedStarList) {
127 item->clearBeforeAssoc();
128 }
129 // clear fitted stars
130 if (!useFittedList) fittedStarList.clear();
131
132 for (auto &ccdImage : ccdImageList) {
133 std::shared_ptr<AstrometryTransform> toCommonTangentPlane = ccdImage->getPixelToCommonTangentPlane();
134
135 // Clear the catalog to fit and copy the whole catalog into it.
136 // This allows reassociating from scratch after a fit.
137 ccdImage->resetCatalogForFit();
138 MeasuredStarList &catalog = ccdImage->getCatalogForFit();
139
140 // Associate with previous lists.
141 /* To speed up the match (more precisely the contruction of the FastFinder), select in the
142 fittedStarList the objects that are within reach of the current ccdImage */
143 Frame ccdImageFrameCPT = toCommonTangentPlane->apply(ccdImage->getImageFrame(), false);
144 ccdImageFrameCPT = ccdImageFrameCPT.rescale(1.10); // add 10 % margin.
145 // We cannot use FittedStarList::ExtractInFrame, because it does an actual copy, which we don't want
146 // here: we want the pointers in the StarMatch to refer to fittedStarList elements.
147 FittedStarList toMatch;
148
149 for (auto const &fittedStar : fittedStarList) {
150 if (ccdImageFrameCPT.inFrame(*fittedStar)) {
151 toMatch.push_back(fittedStar);
152 }
153 }
154
155 // divide by 3600 because coordinates in CTP are in degrees.
156 auto starMatchList = listMatchCollect(Measured2Base(catalog), Fitted2Base(toMatch),
157 toCommonTangentPlane.get(), matchCutInArcSec / 3600.);
158
159 /* should check what this removeAmbiguities does... */
160 LOGLS_DEBUG(_log, "Measured-to-Fitted matches before removing ambiguities " << starMatchList->size());
161 starMatchList->removeAmbiguities(*toCommonTangentPlane);
162 LOGLS_DEBUG(_log, "Measured-to-Fitted matches after removing ambiguities " << starMatchList->size());
163
164 // Associate MeasuredStar -> FittedStar using the surviving matches.
165
166 int matchedCount = 0;
167 for (auto const &starMatch : *starMatchList) {
168 auto bs = starMatch.s1;
169 auto ms_const = std::dynamic_pointer_cast<const MeasuredStar>(bs);
170 auto ms = std::const_pointer_cast<MeasuredStar>(ms_const);
171 auto bs2 = starMatch.s2;
172 auto fs_const = std::dynamic_pointer_cast<const FittedStar>(bs2);
173 auto fs = std::const_pointer_cast<FittedStar>(fs_const);
174 ms->setFittedStar(fs);
175 matchedCount++;
176 }
177 LOGLS_DEBUG(_log, "Matched " << matchedCount << " objects in " << ccdImage->getName());
178
179 // add unmatched objets to FittedStarList
180 int unMatchedCount = 0;
181 for (auto const &mstar : catalog) {
182 // to check if it was matched, just check if it has a fittedStar Pointer assigned
183 if (mstar->getFittedStar()) continue;
184 if (enlargeFittedList) {
185 auto fs = std::make_shared<FittedStar>(*mstar);
186 // transform coordinates to CommonTangentPlane
187 toCommonTangentPlane->transformPosAndErrors(*fs, *fs);
189 mstar->setFittedStar(fs);
190 }
191 unMatchedCount++;
192 }
193 LOGLS_DEBUG(_log, "Unmatched objects: " << unMatchedCount);
194 } // end of loop on CcdImages
195
196 // !!!!!!!!!!!!!!!!!
197 // TODO: DO WE REALLY NEED THIS???
198 // Why do we need to do this, instead of directly computing them in normalizeFittedStars?
199 // What makes the magnitudes special here?
200 // !!!!!!!!!!!!!!!!!
201 // assignMags();
202}
203
205 std::string const &fluxField, float refCoordinateErr,
206 bool rejectBadFluxes) {
207 if (refCat.size() == 0) {
209 " reference catalog is empty : stop here "));
210 }
211
212 afw::table::CoordKey coordKey = refCat.getSchema()["coord"];
213 // Handle reference catalogs that don't have position errors.
214 afw::table::Key<float> raErrKey;
215 afw::table::Key<float> decErrKey;
216 if (std::isnan(refCoordinateErr)) {
217 raErrKey = refCat.getSchema()["coord_raErr"];
218 decErrKey = refCat.getSchema()["coord_decErr"];
219 }
220
221 auto fluxKey = refCat.getSchema().find<double>(fluxField).key;
222 // Handle reference catalogs that don't have flux errors.
223 afw::table::Key<double> fluxErrKey;
224 try {
225 fluxErrKey = refCat.getSchema().find<double>(fluxField + "Err").key;
227 LOGLS_WARN(_log, "Flux error field ("
228 << fluxField << "Err"
229 << ") not found in reference catalog. Not using ref flux errors.");
230 }
231
232 // Handle reference catalogs that don't have proper motion & error
233 afw::table::Key<geom::Angle> pmRaKey, pmDecKey;
234 afw::table::Key<float> pmRaErrKey, pmDecErrKey, pmRaDecCovKey;
235 try {
236 pmRaKey = refCat.getSchema().find<geom::Angle>("pm_ra").key;
237 pmDecKey = refCat.getSchema().find<geom::Angle>("pm_dec").key;
239 LOGLS_WARN(_log, "Not loading proper motions: (pm_ra,pm_dec) fields not found in reference catalog.");
240 } catch (pex::exceptions::TypeError &ex) {
241 LOGLS_WARN(_log, "Not loading proper motions: RA/Dec proper motion values must be `geom:Angle`: "
242 << ex.what());
243 }
244 try {
245 pmRaErrKey = refCat.getSchema().find<float>("pm_raErr").key;
246 pmDecErrKey = refCat.getSchema().find<float>("pm_decErr").key;
248 LOGLS_WARN(_log, "Not loading proper motions: error fields not available: " << ex.what());
249 }
250
251 // TODO: we aren't getting covariances from Gaia yet, so maybe ignore this for now?
252 try {
253 pmRaDecCovKey = refCat.getSchema().find<float>("pm_ra_pm_dec_Cov").key;
255 LOGLS_WARN(_log, "No ra/dec proper motion covariances in refcat: " << ex.what());
256 }
257
259 for (size_t i = 0; i < refCat.size(); i++) {
260 auto const &record = refCat.get(i);
261
262 auto coord = record->get(coordKey);
263 double flux = record->get(fluxKey);
264 double fluxErr;
265 if (fluxErrKey.isValid()) {
266 fluxErr = record->get(fluxErrKey);
267 } else {
269 }
270 double ra = lsst::geom::radToDeg(coord.getLongitude());
271 double dec = lsst::geom::radToDeg(coord.getLatitude());
272 auto star = std::make_shared<RefStar>(ra, dec, flux, fluxErr);
273
274 if (std::isnan(refCoordinateErr)) {
275 // refcat errors are unitless but stored as radians: convert to deg**2
276 star->vx = std::pow(lsst::geom::radToDeg(record->get(raErrKey)), 2);
277 star->vy = std::pow(lsst::geom::radToDeg(record->get(decErrKey)), 2);
278 } else {
279 // Convert the fake errors from mas to deg**2
280 star->vx = std::pow(refCoordinateErr / 1000. / 3600. / std::cos(coord.getLatitude()), 2);
281 star->vy = std::pow(refCoordinateErr / 1000. / 3600., 2);
282 }
283
284 if (pmRaKey.isValid()) {
285 if (pmRaDecCovKey.isValid()) {
286 star->setProperMotion(std::make_unique<ProperMotion const>(
287 record->get(pmRaKey).asRadians(), record->get(pmDecKey).asRadians(),
288 record->get(pmRaErrKey), record->get(pmDecErrKey), record->get(pmRaDecCovKey)));
289 } else {
290 star->setProperMotion(std::make_unique<ProperMotion const>(
291 record->get(pmRaKey).asRadians(), record->get(pmDecKey).asRadians(),
292 record->get(pmRaErrKey), record->get(pmDecErrKey)));
293 }
294 }
295
296 // TODO: cook up a covariance as none of our current refcats have it
297 star->vxy = 0.;
298
299 // Reject sources with non-finite fluxes and flux errors, and fluxErr=0 (which gives chi2=inf).
300 if (rejectBadFluxes && (!std::isfinite(flux) || !std::isfinite(fluxErr) || fluxErr <= 0)) continue;
302 }
303
304 // project on CTP (i.e. RaDec2CTP), in degrees
306 TanRaDecToPixel raDecToCommonTangentPlane(identity, _commonTangentPoint);
307
308 associateRefStars(matchCut.asArcseconds(), &raDecToCommonTangentPlane);
309}
310
311void Associations::associateRefStars(double matchCutInArcSec, const AstrometryTransform *transform) {
312 // associate with FittedStars
313 // 3600 because coordinates are in degrees (in CTP).
314 auto starMatchList = listMatchCollect(Ref2Base(refStarList), Fitted2Base(fittedStarList), transform,
315 matchCutInArcSec / 3600.);
316
317 LOGLS_DEBUG(_log, "Refcat matches before removing ambiguities " << starMatchList->size());
318 starMatchList->removeAmbiguities(*transform);
319 LOGLS_DEBUG(_log, "Refcat matches after removing ambiguities " << starMatchList->size());
320
321 // actually associate things
322 for (auto const &starMatch : *starMatchList) {
323 const BaseStar &bs = *starMatch.s1;
324 const auto &rs_const = dynamic_cast<const RefStar &>(bs);
325 auto &rs = const_cast<RefStar &>(rs_const);
326 const BaseStar &bs2 = *starMatch.s2;
327 const auto &fs_const = dynamic_cast<const FittedStar &>(bs2);
328 auto &fs = const_cast<FittedStar &>(fs_const);
329 // rs->setFittedStar(*fs);
330 fs.setRefStar(&rs);
331 }
332
333 LOGLS_INFO(_log,
334 "Associated " << starMatchList->size() << " reference stars among " << refStarList.size());
335}
336
337void Associations::prepareFittedStars(int minMeasurements) {
338 selectFittedStars(minMeasurements);
339 normalizeFittedStars();
340}
341
342void Associations::selectFittedStars(int minMeasurements) {
343 LOGLS_INFO(_log, "Fitted stars before measurement # cut: " << fittedStarList.size());
344
345 int totalMeasured = 0, validMeasured = 0;
346
347 // first pass: remove objects that have less than a certain number of measurements.
348 for (auto const &ccdImage : ccdImageList) {
349 MeasuredStarList &catalog = ccdImage->getCatalogForFit();
350 // Iteration happens internal to the loop, as we may delete measuredStars from catalog.
351 for (auto mi = catalog.begin(); mi != catalog.end();) {
352 MeasuredStar &mstar = **mi;
353 ++totalMeasured;
354
355 auto fittedStar = mstar.getFittedStar();
356 // measuredStar has no fittedStar: move on.
357 if (fittedStar == nullptr) {
358 ++mi;
359 continue;
360 }
361
362 // keep FittedStars which either have a minimum number of
363 // measurements, or are matched to a RefStar
364 if (!fittedStar->getRefStar() && fittedStar->getMeasurementCount() < minMeasurements) {
365 fittedStar->getMeasurementCount()--;
366 mi = catalog.erase(mi); // mi now points to the next measuredStar.
367 } else {
368 ++validMeasured;
369 ++mi;
370 }
371 } // end loop on objects in catalog
372 } // end loop on catalogs
373
374 // now FittedStars with less than minMeasurements should have zero measurementCount.
375 for (auto fi = fittedStarList.begin(); fi != fittedStarList.end();) {
376 if ((*fi)->getMeasurementCount() == 0) {
377 fi = fittedStarList.erase(fi);
378 } else {
379 ++fi;
380 }
381 }
382
383 LOGLS_INFO(_log, "Fitted stars after measurement # cut: " << fittedStarList.size());
384 LOGLS_INFO(_log, "Total, valid number of Measured stars: " << totalMeasured << ", " << validMeasured);
385}
386
387void Associations::normalizeFittedStars() {
388 // Clear positions in order to take the average of the measuredStars.
389 for (auto &fittedStar : fittedStarList) {
390 fittedStar->x = 0.0;
391 fittedStar->y = 0.0;
392 fittedStar->setFlux(0.0);
393 fittedStar->getMag() = 0.0;
394 }
395
396 // Iterate over measuredStars to add their values into their fittedStars
397 for (auto const &ccdImage : ccdImageList) {
398 std::shared_ptr<AstrometryTransform> toCommonTangentPlane = ccdImage->getPixelToCommonTangentPlane();
399 MeasuredStarList &catalog = ccdImage->getCatalogForFit();
400 for (auto &mi : catalog) {
401 auto fittedStar = mi->getFittedStar();
402 if (fittedStar == nullptr)
403 throw(LSST_EXCEPT(
404 pex::exceptions::RuntimeError,
405 "All measuredStars must have a fittedStar: did you call selectFittedStars()?"));
406 auto point = toCommonTangentPlane->apply(*mi);
407 fittedStar->x += point.x;
408 fittedStar->y += point.y;
409 fittedStar->getFlux() += mi->getFlux();
410 }
411 }
412
413 _maxMeasuredStars = 0;
414 for (auto &fi : fittedStarList) {
415 auto measurementCount = fi->getMeasurementCount();
416 _maxMeasuredStars += measurementCount;
417 fi->x /= measurementCount;
418 fi->y /= measurementCount;
419 fi->getFlux() /= measurementCount;
420 fi->getMag() = utils::nanojanskyToABMagnitude(fi->getFlux());
421 }
422}
423
425 auto iter = fittedStarList.begin();
426 auto end = fittedStarList.end();
427 size_t count = 0;
428 while (iter != end) {
429 auto fittedStar = *iter;
430 if (fittedStar->getMeasurementCount() == 0) {
431 LOGLS_TRACE(_log, "Deleting FittedStar (has no measuredStars): " << *fittedStar);
432 iter = fittedStarList.erase(iter);
433 count++;
434 } else {
435 ++iter;
436 }
437 }
438 if (count > 0) {
439 LOGLS_INFO(_log, "Removed " << count << " fittedStars that had no associated measuredStar.");
440 }
441}
442
443void Associations::assignMags() {
444 for (auto const &ccdImage : ccdImageList) {
445 MeasuredStarList &catalog = ccdImage->getCatalogForFit();
446 for (auto const &mstar : catalog) {
447 auto fstar = mstar->getFittedStar();
448 if (!fstar) continue;
449 fstar->addMagMeasurement(mstar->getMag(), mstar->getMagWeight());
450 }
451 }
452}
453
455 // By default, Associations::fittedStarList is expressed on the Associations::commonTangentPlane.
456 // For AstrometryFit, we need it on the sky.
458 LOGLS_WARN(_log,
459 "DeprojectFittedStars: Fitted stars are already in sidereal coordinates, nothing done ");
460 return;
461 }
462
466}
467
470 return item->getCatalogForFit().size() > 0;
471 });
472}
473
475 size_t count = 0;
476 for (auto const &fittedStar : fittedStarList) {
477 if ((fittedStar != nullptr) && (fittedStar->getRefStar() != nullptr)) count++;
478 }
479 return count;
480}
481
482#ifdef TODO
483void Associations::collectMCStars(int realization) {
484 CcdImageIterator I;
486
487 for (I = ccdImageList.begin(); I != ccdImageList.end(); I++) {
488 CcdImage &ccdImage = **I;
489 string dbimdir = ccdImage.Dir();
490 string mctruth = dbimdir + "/mc/mctruth.list";
491
492 if (realization >= 0) {
493 stringstream sstrm;
494 sstrm << dbimdir << "/mc/mctruth_" << realization << ".list";
495 mctruth = sstrm.str();
496 }
497
498 AstrometryTransformIdentity gti;
499 MeasuredStarList &catalog = ccdImage.getCatalogForFit();
500
501 // BaseStarWithErrorList mctruthlist(mctruth);
502 DicStarList mctruthlist(mctruth);
503 auto starMatchList =
504 listMatchCollect(Measured2Base(catalog), Dic2Base(mctruthlist), &gti, 1. /* pixel ? */);
505 if (starMatchList)
506 for (smI = starMatchList->begin(); smI != starMatchList->end(); smI++) {
507 StarMatch &sm = *smI;
508 BaseStar *bs = sm.s1;
509 MeasuredStar *mstar = dynamic_cast<MeasuredStar *>(bs);
510 bs = sm.s2;
511 DicStar *dstar = dynamic_cast<DicStar *>(bs);
512 std::unique_ptr<BaseStarWithError> mcstar(new BaseStarWithError(*bs));
513 mcstar->GetMCInfo().iflux = dstar->getval("iflux");
514 mcstar->GetMCInfo().tflux = dstar->getval("sflux");
515 /*
516 mstar->SetMCTruth(mcstar);
517 mstar->SetMCMeas(mcstar);
518 */
519 }
520 else
521 LOGLS_FATAL(_log, "CollectMCStars Unable to match MCTruth w/ catalog!");
522 }
523}
524
525#endif /* TODO */
526} // namespace jointcal
527} // namespace lsst
AmpInfoBoxKey bbox
Definition Amplifier.cc:117
int end
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...
#define LSST_EXCEPT(type,...)
Create an exception with a given type.
Definition Exception.h:48
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.
Definition Log.h:659
#define LOGLS_FATAL(logger, message)
Log a fatal-level message using an iostream-based interface.
Definition Log.h:699
#define LOG_GET(logger)
Returns a Log object associated with logger.
Definition Log.h:75
#define LOGLS_INFO(logger, message)
Log a info-level message using an iostream-based interface.
Definition Log.h:639
#define LOG_LOGGER
Definition Log.h:714
#define LOGLS_DEBUG(logger, message)
Log a debug-level message using an iostream-based interface.
Definition Log.h:619
#define LOGLS_TRACE(logger, message)
Log a trace-level message using an iostream-based interface.
Definition Log.h:599
This file contains a class representing spherical coordinates.
double dec
Definition Match.cc:41
pairs of points
T begin(T... args)
Tag types used to declare specialized field types.
Definition misc.h:31
A FunctorKey used to get or set celestial coordinates from a pair of lsst::geom::Angle keys.
Definition aggregates.h:292
Custom catalog class for record/table subclasses that are guaranteed to have an ID,...
A class representing an angle.
Definition Angle.h:128
constexpr double asArcseconds() const noexcept
Return an Angle's value in arcseconds.
Definition Angle.h:185
An integer coordinate rectangle.
Definition Box.h:55
size_t nFittedStarsWithAssociatedRefStar() const
Return the number of fittedStars that have an associated refStar.
void cleanFittedStars()
Remove FittedStars that have no measured stars; this can happen after outlier rejection.
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, const std::shared_ptr< lsst::afw::geom::SkyWcs > &wcs, const std::shared_ptr< lsst::afw::image::VisitInfo > &visitInfo, const lsst::geom::Box2I &bbox, const std::string &filter, const std::shared_ptr< afw::image::PhotoCalib > &photoCalib, const std::shared_ptr< afw::cameraGeom::Detector > &detector, int visit, int ccd, const lsst::jointcal::JointcalControl &control)
Create a ccdImage from an exposure catalog and metadata, and add it to the list.
void associateCatalogs(double matchCutInArcsec=0, bool useFittedList=false, bool enlargeFittedList=true)
incrementaly builds a merged catalog of all image catalogs
void setCommonTangentPoint(lsst::geom::Point2D const &commonTangentPoint)
Shared tangent point for all ccdImages (decimal degrees).
Point getCommonTangentPoint() const
Shared tangent point for all ccdImages (decimal degrees).
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.
void prepareFittedStars(int minMeasurements)
Prepare the fittedStar list by making quality cuts and normalizing measurements.
void deprojectFittedStars()
Sends back the fitted stars coordinates on the sky FittedStarsList::inTangentPlaneCoordinates keeps t...
a virtual (interface) class for geometric transformations.
implements the linear transformations (6 real coefficients).
The base class for handling stars. Used by all matching routines.
Definition BaseStar.h:51
double getMag() const
Definition BaseStar.h:105
Handler of an actual image from a single CCD.
Definition CcdImage.h:64
MeasuredStarList const & getCatalogForFit() const
Gets the catalog to be used for fitting, which may have been cleaned-up.
Definition CcdImage.h:94
FittedStars are objects whose position or flux is going to be fitted, and which come from the associa...
Definition FittedStar.h:50
void setRefStar(const RefStar *_refStar)
Set the astrometric reference star associated with this star.
Definition FittedStar.cc:46
A list of FittedStar s. Such a list is typically constructed by Associations.
Definition FittedStar.h:116
rectangle with sides parallel to axes.
Definition Frame.h:38
bool inFrame(double x, double y) const
inside?
Definition Frame.cc:120
double xMin
coordinate of boundary.
Definition Frame.h:41
Frame rescale(double factor) const
rescale it. The center does not move.
Definition Frame.cc:110
double getArea() const
Definition Frame.cc:118
Sources measured on 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.
A point in a plane.
Definition Point.h:37
Objects used as position/flux anchors (e.g.
Definition RefStar.h:45
void applyTransform(const Operator &op)
enables to apply a geometrical transform if Star is Basestar or derives from it.
Definition StarList.h:98
The transformation that handles pixels to sideral transformations (Gnomonic, possibly with polynomial...
virtual void apply(double xIn, double yIn, double &xOut, double &yOut) const=0
Transform pixels to ICRS RA, Dec in degrees.
This one is the Tangent Plane (called gnomonic) projection (from celestial sphere to tangent plane)
Reports invalid arguments.
Definition Runtime.h:66
Reports attempts to access elements using an invalid key.
Definition Runtime.h:151
Reports errors from accepting an object of an unexpected or inappropriate type.
Definition Runtime.h:167
Circle is a circular region on the unit sphere that contains its boundary.
Definition Circle.h:53
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)
Definition LonLat.h:57
T clear(T... args)
T cos(T... args)
T count(T... args)
T emplace_back(T... args)
T end(T... args)
T erase(T... args)
T get(T... args)
T isfinite(T... args)
T isnan(T... args)
double nanojanskyToABMagnitude(double flux)
Convert a flux in nanojansky to AB magnitude.
Definition Magnitude.cc:30
SpherePoint averageSpherePoint(std::vector< SpherePoint > const &coords)
Return the average of a list of coordinates.
AngleUnit constexpr degrees
constant with units of degrees
Definition Angle.h:110
constexpr double radToDeg(double x) noexcept
Definition Angle.h:53
BaseStarList & Measured2Base(MeasuredStarList &This)
BaseStarList & Fitted2Base(FittedStarList &This)
Definition FittedStar.cc:64
std::unique_ptr< StarMatchList > listMatchCollect(const BaseStarList &list1, const BaseStarList &list2, const AstrometryTransform *guess, double maxDist)
assembles star matches.
Definition ListMatch.cc:568
::std::list< StarMatch >::iterator StarMatchIterator
Definition StarMatch.h:134
BaseStarList & Ref2Base(RefStarList &This)
Definition RefStar.cc:42
T pow(T... args)
T push_back(T... args)
T quiet_NaN(T... args)
T reserve(T... args)
T size(T... args)
T str(T... args)
std::string sourceFluxField
"name of flux field in source catalog" ;
Key< int > visitInfo
Definition Exposure.cc:70