LSSTApplications  19.0.0-14-gb0260a2+72efe9b372,20.0.0+7927753e06,20.0.0+8829bf0056,20.0.0+995114c5d2,20.0.0+b6f4b2abd1,20.0.0+bddc4f4cbe,20.0.0-1-g253301a+8829bf0056,20.0.0-1-g2b7511a+0d71a2d77f,20.0.0-1-g5b95a8c+7461dd0434,20.0.0-12-g321c96ea+23efe4bbff,20.0.0-16-gfab17e72e+fdf35455f6,20.0.0-2-g0070d88+ba3ffc8f0b,20.0.0-2-g4dae9ad+ee58a624b3,20.0.0-2-g61b8584+5d3db074ba,20.0.0-2-gb780d76+d529cf1a41,20.0.0-2-ged6426c+226a441f5f,20.0.0-2-gf072044+8829bf0056,20.0.0-2-gf1f7952+ee58a624b3,20.0.0-20-geae50cf+e37fec0aee,20.0.0-25-g3dcad98+544a109665,20.0.0-25-g5eafb0f+ee58a624b3,20.0.0-27-g64178ef+f1f297b00a,20.0.0-3-g4cc78c6+e0676b0dc8,20.0.0-3-g8f21e14+4fd2c12c9a,20.0.0-3-gbd60e8c+187b78b4b8,20.0.0-3-gbecbe05+48431fa087,20.0.0-38-ge4adf513+a12e1f8e37,20.0.0-4-g97dc21a+544a109665,20.0.0-4-gb4befbc+087873070b,20.0.0-4-gf910f65+5d3db074ba,20.0.0-5-gdfe0fee+199202a608,20.0.0-5-gfbfe500+d529cf1a41,20.0.0-6-g64f541c+d529cf1a41,20.0.0-6-g9a5b7a1+a1cd37312e,20.0.0-68-ga3f3dda+5fca18c6a4,20.0.0-9-g4aef684+e18322736b,w.2020.45
LSSTDataManagementBasePackage
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"
32 #include "lsst/jointcal/CcdImage.h"
35 #include "lsst/jointcal/Frame.h"
36 #include "lsst/jointcal/FatPoint.h"
39 
40 #include "lsst/afw/image/Image.h"
43 
44 #include "lsst/pex/exceptions.h"
45 #include "lsst/geom/Box.h"
46 #include "lsst/geom/Point.h"
47 #include "lsst/afw/image/Calib.h"
48 
49 #include "lsst/sphgeom/LonLat.h"
50 #include "lsst/sphgeom/Circle.h"
52 
53 namespace jointcal = lsst::jointcal;
54 
55 namespace {
56 LOG_LOGGER _log = LOG_GET("jointcal.Associations");
57 }
58 
59 namespace lsst {
60 namespace jointcal {
61 
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  LOGLS_DEBUG(_log, "Catalog " << ccdImage->getName() << " has " << ccdImage->getWholeCatalog().size()
73  << " objects.");
74 }
75 
78  centers.reserve(ccdImageList.size());
79  for (auto const &ccdImage : ccdImageList) {
80  centers.push_back(ccdImage->getBoresightRaDec());
81  }
82  auto commonTangentPoint = geom::averageSpherePoint(centers);
83  LOGLS_DEBUG(_log, "Using common tangent point: " << commonTangentPoint.getPosition(geom::degrees));
84  setCommonTangentPoint(commonTangentPoint.getPosition(geom::degrees));
85 }
86 
88  _commonTangentPoint = Point(commonTangentPoint.getX(), commonTangentPoint.getY()); // a jointcal::Point
89  for (auto &ccdImage : ccdImageList) ccdImage->setCommonTangentPoint(_commonTangentPoint);
90 }
91 
93  // Compute the frame on the common tangent plane that contains all input images.
94  Frame tangentPlaneFrame;
95 
96  for (auto const &ccdImage : ccdImageList) {
97  Frame CTPFrame = ccdImage->getPixelToCommonTangentPlane()->apply(ccdImage->getImageFrame(), false);
98  if (tangentPlaneFrame.getArea() == 0)
99  tangentPlaneFrame = CTPFrame;
100  else
101  tangentPlaneFrame += CTPFrame;
102  }
103 
104  // Convert tangent plane coordinates to RaDec.
105  AstrometryTransformLinear identity;
106  TanPixelToRaDec commonTangentPlaneToRaDec(identity, _commonTangentPoint);
107  Frame raDecFrame = commonTangentPlaneToRaDec.apply(tangentPlaneFrame, false);
108 
110  points.reserve(4);
111  // raDecFrame is in on-sky (RA,Dec) degrees stored as an x/y box:
112  // the on-sky bounding box it represents is given by the corners of that box.
113  points.emplace_back(sphgeom::LonLat::fromDegrees(raDecFrame.xMin, raDecFrame.yMin));
114  points.emplace_back(sphgeom::LonLat::fromDegrees(raDecFrame.xMax, raDecFrame.yMin));
115  points.emplace_back(sphgeom::LonLat::fromDegrees(raDecFrame.xMin, raDecFrame.yMax));
116  points.emplace_back(sphgeom::LonLat::fromDegrees(raDecFrame.xMax, raDecFrame.yMax));
117 
119 }
120 
121 void Associations::associateCatalogs(const double matchCutInArcSec, const bool useFittedList,
122  const bool enlargeFittedList) {
123  // clear reference stars
124  refStarList.clear();
125 
126  // clear measurement counts and associations to refstars, but keep fittedStars themselves.
127  for (auto &item : fittedStarList) {
128  item->clearBeforeAssoc();
129  }
130  // clear fitted stars
131  if (!useFittedList) fittedStarList.clear();
132 
133  for (auto &ccdImage : ccdImageList) {
134  std::shared_ptr<AstrometryTransform> toCommonTangentPlane = ccdImage->getPixelToCommonTangentPlane();
135 
136  // Clear the catalog to fit and copy the whole catalog into it.
137  // This allows reassociating from scratch after a fit.
138  ccdImage->resetCatalogForFit();
139  MeasuredStarList &catalog = ccdImage->getCatalogForFit();
140 
141  // Associate with previous lists.
142  /* To speed up the match (more precisely the contruction of the FastFinder), select in the
143  fittedStarList the objects that are within reach of the current ccdImage */
144  Frame ccdImageFrameCPT = toCommonTangentPlane->apply(ccdImage->getImageFrame(), false);
145  ccdImageFrameCPT = ccdImageFrameCPT.rescale(1.10); // add 10 % margin.
146  // We cannot use FittedStarList::ExtractInFrame, because it does an actual copy, which we don't want
147  // here: we want the pointers in the StarMatch to refer to fittedStarList elements.
148  FittedStarList toMatch;
149 
150  for (auto const &fittedStar : fittedStarList) {
151  if (ccdImageFrameCPT.inFrame(*fittedStar)) {
152  toMatch.push_back(fittedStar);
153  }
154  }
155 
156  // divide by 3600 because coordinates in CTP are in degrees.
157  auto starMatchList = listMatchCollect(Measured2Base(catalog), Fitted2Base(toMatch),
158  toCommonTangentPlane.get(), matchCutInArcSec / 3600.);
159 
160  /* should check what this removeAmbiguities does... */
161  LOGLS_DEBUG(_log, "Measured-to-Fitted matches before removing ambiguities " << starMatchList->size());
162  starMatchList->removeAmbiguities(*toCommonTangentPlane);
163  LOGLS_DEBUG(_log, "Measured-to-Fitted matches after removing ambiguities " << starMatchList->size());
164 
165  // Associate MeasuredStar -> FittedStar using the surviving matches.
166 
167  int matchedCount = 0;
168  for (auto const &starMatch : *starMatchList) {
169  auto bs = starMatch.s1;
170  auto ms_const = std::dynamic_pointer_cast<const MeasuredStar>(bs);
171  auto ms = std::const_pointer_cast<MeasuredStar>(ms_const);
172  auto bs2 = starMatch.s2;
173  auto fs_const = std::dynamic_pointer_cast<const FittedStar>(bs2);
174  auto fs = std::const_pointer_cast<FittedStar>(fs_const);
175  ms->setFittedStar(fs);
176  matchedCount++;
177  }
178  LOGLS_INFO(_log, "Matched " << matchedCount << " objects in " << ccdImage->getName());
179 
180  // add unmatched objets to FittedStarList
181  int unMatchedCount = 0;
182  for (auto const &mstar : catalog) {
183  // to check if it was matched, just check if it has a fittedStar Pointer assigned
184  if (mstar->getFittedStar()) continue;
185  if (enlargeFittedList) {
186  auto fs = std::make_shared<FittedStar>(*mstar);
187  // transform coordinates to CommonTangentPlane
188  toCommonTangentPlane->transformPosAndErrors(*fs, *fs);
190  mstar->setFittedStar(fs);
191  }
192  unMatchedCount++;
193  }
194  LOGLS_INFO(_log, "Unmatched objects: " << unMatchedCount);
195  } // end of loop on CcdImages
196 
197  // !!!!!!!!!!!!!!!!!
198  // TODO: DO WE REALLY NEED THIS???
199  // Why do we need to do this, instead of directly computing them in normalizeFittedStars?
200  // What makes the magnitudes special here?
201  // !!!!!!!!!!!!!!!!!
202  // assignMags();
203 }
204 
206  std::string const &fluxField, float refCoordinateErr,
207  bool rejectBadFluxes) {
208  if (refCat.size() == 0) {
210  " reference catalog is empty : stop here "));
211  }
212 
213  afw::table::CoordKey coordKey = refCat.getSchema()["coord"];
214  // Handle reference catalogs that don't have position errors.
215  afw::table::Key<float> raErrKey;
216  afw::table::Key<float> decErrKey;
217  if (std::isnan(refCoordinateErr)) {
218  raErrKey = refCat.getSchema()["coord_raErr"];
219  decErrKey = refCat.getSchema()["coord_decErr"];
220  }
221 
222  auto fluxKey = refCat.getSchema().find<double>(fluxField).key;
223  // Handle reference catalogs that don't have flux errors.
224  afw::table::Key<double> fluxErrKey;
225  try {
226  fluxErrKey = refCat.getSchema().find<double>(fluxField + "Err").key;
227  } catch (pex::exceptions::NotFoundError &) {
228  LOGLS_WARN(_log, "Flux error field ("
229  << fluxField << "Err"
230  << ") not found in reference catalog. Not using ref flux errors.");
231  }
232 
233  refStarList.clear();
234  for (size_t i = 0; i < refCat.size(); i++) {
235  auto const &record = refCat.get(i);
236 
237  auto coord = record->get(coordKey);
238  double flux = record->get(fluxKey);
239  double fluxErr;
240  if (fluxErrKey.isValid()) {
241  fluxErr = record->get(fluxErrKey);
242  } else {
244  }
245  double ra = lsst::geom::radToDeg(coord.getLongitude());
246  double dec = lsst::geom::radToDeg(coord.getLatitude());
247  auto star = std::make_shared<RefStar>(ra, dec, flux, fluxErr);
248 
249  if (std::isnan(refCoordinateErr)) {
250  star->vx = record->get(raErrKey);
251  star->vy = record->get(decErrKey);
252  } else {
253  // Compute and use the fake errors
254  star->vx = std::pow(refCoordinateErr / 1000. / 3600. / std::cos(coord.getLatitude()), 2);
255  star->vy = std::pow(refCoordinateErr / 1000. / 3600., 2);
256  }
257  // TODO: cook up a covariance as none of our current refcats have it
258  star->vxy = 0.;
259 
260  // Reject sources with non-finite fluxes and flux errors, and fluxErr=0 (which gives chi2=inf).
261  if (rejectBadFluxes && (!std::isfinite(flux) || !std::isfinite(fluxErr) || fluxErr <= 0)) continue;
262  refStarList.push_back(star);
263  }
264 
265  // project on CTP (i.e. RaDec2CTP), in degrees
266  AstrometryTransformLinear identity;
267  TanRaDecToPixel raDecToCommonTangentPlane(identity, _commonTangentPoint);
268 
269  associateRefStars(matchCut.asArcseconds(), &raDecToCommonTangentPlane);
270 }
271 
272 void Associations::associateRefStars(double matchCutInArcSec, const AstrometryTransform *transform) {
273  // associate with FittedStars
274  // 3600 because coordinates are in degrees (in CTP).
276  matchCutInArcSec / 3600.);
277 
278  LOGLS_DEBUG(_log, "Refcat matches before removing ambiguities " << starMatchList->size());
279  starMatchList->removeAmbiguities(*transform);
280  LOGLS_DEBUG(_log, "Refcat matches after removing ambiguities " << starMatchList->size());
281 
282  // actually associate things
283  for (auto const &starMatch : *starMatchList) {
284  const BaseStar &bs = *starMatch.s1;
285  const RefStar &rs_const = dynamic_cast<const RefStar &>(bs);
286  RefStar &rs = const_cast<RefStar &>(rs_const);
287  const BaseStar &bs2 = *starMatch.s2;
288  const FittedStar &fs_const = dynamic_cast<const FittedStar &>(bs2);
289  FittedStar &fs = const_cast<FittedStar &>(fs_const);
290  // rs->setFittedStar(*fs);
291  fs.setRefStar(&rs);
292  }
293 
294  LOGLS_INFO(_log,
295  "Associated " << starMatchList->size() << " reference stars among " << refStarList.size());
296 }
297 
298 void Associations::prepareFittedStars(int minMeasurements) {
299  selectFittedStars(minMeasurements);
300  normalizeFittedStars();
301 }
302 
303 void Associations::selectFittedStars(int minMeasurements) {
304  LOGLS_INFO(_log, "Fitted stars before measurement # cut: " << fittedStarList.size());
305 
306  int totalMeasured = 0, validMeasured = 0;
307 
308  // first pass: remove objects that have less than a certain number of measurements.
309  for (auto const &ccdImage : ccdImageList) {
310  MeasuredStarList &catalog = ccdImage->getCatalogForFit();
311  // Iteration happens internal to the loop, as we may delete measuredStars from catalog.
312  for (MeasuredStarIterator mi = catalog.begin(); mi != catalog.end();) {
313  MeasuredStar &mstar = **mi;
314  ++totalMeasured;
315 
316  auto fittedStar = mstar.getFittedStar();
317  // measuredStar has no fittedStar: move on.
318  if (fittedStar == nullptr) {
319  ++mi;
320  continue;
321  }
322 
323  // keep FittedStars which either have a minimum number of
324  // measurements, or are matched to a RefStar
325  if (!fittedStar->getRefStar() && fittedStar->getMeasurementCount() < minMeasurements) {
326  fittedStar->getMeasurementCount()--;
327  mi = catalog.erase(mi); // mi now points to the next measuredStar.
328  } else {
329  ++validMeasured;
330  ++mi;
331  }
332  } // end loop on objects in catalog
333  } // end loop on catalogs
334 
335  // now FittedStars with less than minMeasurements should have zero measurementCount.
336  for (FittedStarIterator fi = fittedStarList.begin(); fi != fittedStarList.end();) {
337  if ((*fi)->getMeasurementCount() == 0) {
338  fi = fittedStarList.erase(fi);
339  } else {
340  ++fi;
341  }
342  }
343 
344  LOGLS_INFO(_log, "Fitted stars after measurement # cut: " << fittedStarList.size());
345  LOGLS_INFO(_log, "Total, valid number of Measured stars: " << totalMeasured << ", " << validMeasured);
346 }
347 
348 void Associations::normalizeFittedStars() const {
349  // Clear positions in order to take the average of the measuredStars.
350  for (auto &fittedStar : fittedStarList) {
351  fittedStar->x = 0.0;
352  fittedStar->y = 0.0;
353  fittedStar->setFlux(0.0);
354  fittedStar->getMag() = 0.0;
355  }
356 
357  // Iterate over measuredStars to add their values into their fittedStars
358  for (auto const &ccdImage : ccdImageList) {
359  std::shared_ptr<AstrometryTransform> toCommonTangentPlane = ccdImage->getPixelToCommonTangentPlane();
360  MeasuredStarList &catalog = ccdImage->getCatalogForFit();
361  for (auto &mi : catalog) {
362  auto fittedStar = mi->getFittedStar();
363  if (fittedStar == nullptr)
364  throw(LSST_EXCEPT(
365  pex::exceptions::RuntimeError,
366  "All measuredStars must have a fittedStar: did you call selectFittedStars()?"));
367  auto point = toCommonTangentPlane->apply(*mi);
368  fittedStar->x += point.x;
369  fittedStar->y += point.y;
370  fittedStar->getFlux() += mi->getFlux();
371  }
372  }
373 
374  for (auto &fi : fittedStarList) {
375  auto measurementCount = fi->getMeasurementCount();
376  fi->x /= measurementCount;
377  fi->y /= measurementCount;
378  fi->getFlux() /= measurementCount;
379  fi->getMag() = utils::nanojanskyToABMagnitude(fi->getFlux());
380  }
381 }
382 
383 void Associations::assignMags() {
384  for (auto const &ccdImage : ccdImageList) {
385  MeasuredStarList &catalog = ccdImage->getCatalogForFit();
386  for (auto const &mstar : catalog) {
387  auto fstar = mstar->getFittedStar();
388  if (!fstar) continue;
389  fstar->addMagMeasurement(mstar->getMag(), mstar->getMagWeight());
390  }
391  }
392 }
393 
395  // By default, Associations::fittedStarList is expressed on the Associations::commonTangentPlane.
396  // For AstrometryFit, we need it on the sky.
398  LOGLS_WARN(_log,
399  "DeprojectFittedStars: Fitted stars are already in sidereal coordinates, nothing done ");
400  return;
401  }
402 
406 }
407 
410  return item->getCatalogForFit().size() > 0;
411  });
412 }
413 
415  size_t count = 0;
416  for (auto const &fittedStar : fittedStarList) {
417  if ((fittedStar != nullptr) & (fittedStar->getRefStar() != nullptr)) count++;
418  }
419  return count;
420 }
421 
422 #ifdef TODO
423 void Associations::collectMCStars(int realization) {
424  CcdImageIterator I;
425  StarMatchIterator smI;
426 
427  for (I = ccdImageList.begin(); I != ccdImageList.end(); I++) {
428  CcdImage &ccdImage = **I;
429  string dbimdir = ccdImage.Dir();
430  string mctruth = dbimdir + "/mc/mctruth.list";
431 
432  if (realization >= 0) {
433  stringstream sstrm;
434  sstrm << dbimdir << "/mc/mctruth_" << realization << ".list";
435  mctruth = sstrm.str();
436  }
437 
438  AstrometryTransformIdentity gti;
439  MeasuredStarList &catalog = ccdImage.getCatalogForFit();
440 
441  // BaseStarWithErrorList mctruthlist(mctruth);
442  DicStarList mctruthlist(mctruth);
443  auto starMatchList =
444  listMatchCollect(Measured2Base(catalog), Dic2Base(mctruthlist), &gti, 1. /* pixel ? */);
445  if (starMatchList)
446  for (smI = starMatchList->begin(); smI != starMatchList->end(); smI++) {
447  StarMatch &sm = *smI;
448  BaseStar *bs = sm.s1;
449  MeasuredStar *mstar = dynamic_cast<MeasuredStar *>(bs);
450  bs = sm.s2;
451  DicStar *dstar = dynamic_cast<DicStar *>(bs);
452  std::unique_ptr<BaseStarWithError> mcstar(new BaseStarWithError(*bs));
453  mcstar->GetMCInfo().iflux = dstar->getval("iflux");
454  mcstar->GetMCInfo().tflux = dstar->getval("sflux");
455  /*
456  mstar->SetMCTruth(mcstar);
457  mstar->SetMCMeas(mcstar);
458  */
459  }
460  else
461  LOGLS_FATAL(_log, "CollectMCStars Unable to match MCTruth w/ catalog!");
462  }
463 }
464 
465 void Associations::setFittedStarColors(std::string dicStarListName, std::string color,
466  double matchCutArcSec) {
467  // decode color string in case it is x-y
468  size_t pos_minus = color.find('-');
469  bool compute_diff = (pos_minus != string::npos);
470  std::string c1, c2;
471  c1 = color.substr(0, pos_minus); // if pos_minus == npos, means "up to the end"
472  if (compute_diff) c2 = color.substr(pos_minus + 1, string::npos);
473  DicStarList cList(dicStarListName);
474  if (!cList.HasKey(c1))
475  throw(GastroException("Associations::SetFittedstarColors : " + dicStarListName +
476  " misses a key named \"" + c1 + "\""));
477  if (compute_diff && !cList.HasKey(c2))
478  throw(GastroException("Associations::SetFittedstarColors : " + dicStarListName +
479  " misses a key named \"" + c2 + "\""));
480  // we associate in some tangent plane. The reference catalog is expressed on the sky,
481  // but FittedStar's may be still in this tangent plane.
483  AstrometryTransformIdentity id;
484  TanRaDecToPixel proj(AstrometryTransformLinear(), getCommonTangentPoint());
485  // project or not ?
486  AstrometryTransform *id_or_proj = &proj;
487  if (fittedStarList.inTangentPlaneCoordinates) id_or_proj = &id;
488  // The color List is to be projected:
489  TStarList projected_cList((BaseStarList &)cList, proj);
490  // Associate
491  auto starMatchList = listMatchCollect(Fitted2Base(fittedStarList), (const BaseStarList &)projected_cList,
492  id_or_proj, matchCutArcSec / 3600);
493 
494  LOGLS_INFO(_log, "Matched " << starMatchList->size() << '/' << fittedStarList.size()
495  << " FittedStars to color catalog");
496  // Evaluate and assign colors.
497  for (auto i = starMatchList->begin(); i != starMatchList->end(); ++i) {
498  BaseStar *s1 = i->s1;
499  FittedStar *fs = dynamic_cast<FittedStar *>(s1);
500  BaseStar *s2 = i->s2;
501  const TStar *ts = dynamic_cast<const TStar *>(s2);
502  const DicStar *ds = dynamic_cast<const DicStar *>(ts->get_original());
503  fs->color = ds->getval(c1);
504  if (compute_diff) fs->color -= ds->getval(c2);
505  }
506 }
507 
508 #endif /* TODO */
509 } // namespace jointcal
510 } // namespace lsst
lsst::geom::degrees
constexpr AngleUnit degrees
constant with units of degrees
Definition: Angle.h:109
lsst::jointcal::MeasuredStar
objects measured on actual images.
Definition: MeasuredStar.h:46
lsst::geom::radToDeg
constexpr double radToDeg(double x) noexcept
Definition: Angle.h:52
LOG_LOGGER
#define LOG_LOGGER
Definition: Log.h:703
lsst::geom::averageSpherePoint
SpherePoint averageSpherePoint(std::vector< SpherePoint > const &coords)
Return the average of a list of coordinates.
Definition: SpherePoint.cc:235
lsst::jointcal::Associations::associateCatalogs
void associateCatalogs(const double matchCutInArcsec=0, const bool useFittedList=false, const bool enlargeFittedList=true)
incrementaly builds a merged catalog of all image catalogs
Definition: Associations.cc:121
Frame.h
std::string
STL class.
std::shared_ptr< lsst::afw::geom::SkyWcs >
Associations.h
lsst::jointcal::Associations::deprojectFittedStars
void deprojectFittedStars()
Sends back the fitted stars coordinates on the sky FittedStarsList::inTangentPlaneCoordinates keeps t...
Definition: Associations.cc:394
lsst::jointcal::Frame::xMax
double xMax
Definition: Frame.h:41
lsst::jointcal::MeasuredStar::getMagWeight
double getMagWeight() const
the inverse of the mag variance
Definition: MeasuredStar.h:106
lsst.pipe.tasks.assembleCoadd.filter
filter
Definition: assembleCoadd.py:342
CcdImage.h
lsst::jointcal::Frame::rescale
Frame rescale(const double factor) const
rescale it. The center does not move.
Definition: Frame.cc:110
lsst::jointcal::Associations::refStarList
RefStarList refStarList
Definition: Associations.h:57
AstrometryTransform.h
std::vector::reserve
T reserve(T... args)
lsst::jointcal::StarList::applyTransform
void applyTransform(const Operator &op)
enables to apply a geometrical transform if Star is Basestar or derives from it.
Definition: StarList.h:98
lsst::jointcal::AstrometryTransformLinear
implements the linear transformations (6 real coefficients).
Definition: AstrometryTransform.h:426
Box.h
wcs
table::Key< table::Array< std::uint8_t > > wcs
Definition: SkyWcs.cc:71
std::cos
T cos(T... args)
std::numeric_limits::quiet_NaN
T quiet_NaN(T... args)
Calib.h
lsst::afw::table::CoordKey
A FunctorKey used to get or set celestial coordinates from a pair of lsst::geom::Angle keys.
Definition: aggregates.h:210
lsst::sphgeom::LonLat::fromDegrees
static LonLat fromDegrees(double lon, double lat)
Definition: LonLat.h:50
std::vector
STL class.
std::string::find
T find(T... args)
MeasuredStar.h
std::list::size
T size(T... args)
lsst.pex::exceptions::NotFoundError
Reports attempts to access elements using an invalid key.
Definition: Runtime.h:151
lsst::jointcal::Measured2Base
BaseStarList & Measured2Base(MeasuredStarList &This)
Definition: MeasuredStar.cc:58
lsst::jointcal::Associations::collectRefStars
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.
Definition: Associations.cc:205
lsst::jointcal::TanRaDecToPixel
This one is the Tangent Plane (called gnomonic) projection (from celestial sphere to tangent plane)
Definition: AstrometryTransform.h:707
dec
double dec
Definition: Match.cc:41
ListMatch.h
Combinatorial searches for linear transformations to go from list1 to list2.
std::stringstream
STL class.
LOG_GET
#define LOG_GET(logger)
Definition: Log.h:75
std::shared_ptr::get
T get(T... args)
lsst::jointcal::MeasuredStarIterator
MeasuredStarList::iterator MeasuredStarIterator
Definition: MeasuredStar.h:154
lsst::jointcal::StarMatchIterator
::std::list< StarMatch >::iterator StarMatchIterator
Definition: StarMatch.h:133
LOGLS_INFO
#define LOGLS_INFO(logger, message)
Definition: Log.h:628
lsst::jointcal::FittedStar
The objects which have been measured several times.
Definition: FittedStar.h:61
lsst::jointcal::Associations::getCommonTangentPoint
Point getCommonTangentPoint() const
can be used to project sidereal coordinates related to the image set on a plane.
Definition: Associations.h:105
lsst::jointcal::RefStar
Objects used as position anchors, typically USNO stars.
Definition: RefStar.h:39
lsst::jointcal::Associations::fittedStarList
FittedStarList fittedStarList
Definition: Associations.h:58
FatPoint.h
Image.h
lsst::jointcal::Ref2Base
BaseStarList & Ref2Base(RefStarList &This)
Definition: RefStar.cc:34
lsst::jointcal::TanPixelToRaDec
The transformation that handles pixels to sideral transformations (Gnomonic, possibly with polynomial...
Definition: AstrometryTransform.h:627
lsst::jointcal::Associations::ccdImageList
CcdImageList ccdImageList
Definition: Associations.h:56
lsst::jointcal::Frame::xMin
double xMin
coordinate of boundary.
Definition: Frame.h:41
Circle.h
This file declares a class for representing circular regions on the unit sphere.
lsst::jointcal::Fitted2Base
BaseStarList & Fitted2Base(FittedStarList &This)
Definition: FittedStar.cc:64
std::list::clear
T clear(T... args)
std::list::push_back
T push_back(T... args)
Point.h
LOGLS_WARN
#define LOGLS_WARN(logger, message)
Definition: Log.h:648
std::isnan
T isnan(T... args)
lsst::afw::table::SortedCatalogT::find
iterator find(RecordId id)
Return an iterator to the record with the given ID.
Definition: SortedCatalog.h:77
lsst::jointcal::listMatchCollect
std::unique_ptr< StarMatchList > listMatchCollect(const BaseStarList &list1, const BaseStarList &list2, const AstrometryTransform *guess, const double maxDist)
assembles star matches.
Definition: ListMatch.cc:569
std::isfinite
T isfinite(T... args)
lsst::jointcal::FittedStarList::inTangentPlaneCoordinates
bool inTangentPlaneCoordinates
Definition: FittedStar.h:125
id
table::Key< int > id
Definition: Detector.cc:162
lsst::jointcal::CcdImage
Handler of an actual image from a single CCD.
Definition: CcdImage.h:64
lsst::jointcal::FittedStar::setRefStar
void setRefStar(const RefStar *_refStar)
Set the astrometric reference star associated with this star.
Definition: FittedStar.cc:46
LOGLS_FATAL
#define LOGLS_FATAL(logger, message)
Definition: Log.h:688
lsst::afw::table._source.SourceCatalog
Definition: _source.py:32
PropertySet.h
lsst::jointcal::Frame::yMin
double yMin
Definition: Frame.h:41
lsst::afw::table::Key< float >
std::list::erase
T erase(T... args)
lsst::jointcal::TanPixelToRaDec::apply
virtual void apply(const double xIn, const double yIn, double &xOut, double &yOut) const=0
lsst::jointcal::BaseStar::getMag
double getMag() const
Definition: BaseStar.h:103
lsst::jointcal::MeasuredStar::getFittedStar
std::shared_ptr< FittedStar > getFittedStar() const
Definition: MeasuredStar.h:113
lsst::jointcal::MeasuredStarList
A list of MeasuredStar. They are usually filled in Associations::createCcdImage.
Definition: MeasuredStar.h:146
lsst::jointcal
Definition: Associations.h:49
lsst::jointcal::Associations::prepareFittedStars
void prepareFittedStars(int minMeasurements)
Set the color field of FittedStar 's from a colored catalog.
Definition: Associations.cc:298
lsst.pipe.tasks.assembleCoadd.ccd
ccd
Definition: assembleCoadd.py:343
lsst
A base class for image defects.
Definition: imageAlgorithm.dox:1
lsst::jointcal::Frame::yMax
double yMax
Definition: Frame.h:41
LSST_EXCEPT
#define LSST_EXCEPT(type,...)
Create an exception with a given type.
Definition: Exception.h:48
lsst::utils::nanojanskyToABMagnitude
double nanojanskyToABMagnitude(double flux)
Convert a flux in nanojansky to AB magnitude.
Definition: Magnitude.cc:30
LOGLS_DEBUG
#define LOGLS_DEBUG(logger, message)
Definition: Log.h:608
std::string::substr
T substr(T... args)
lsst.pipe.tasks.insertFakes.flux
flux
Definition: insertFakes.py:455
lsst::jointcal::Associations::setCommonTangentPoint
void setCommonTangentPoint(lsst::geom::Point2D const &commonTangentPoint)
Sets a shared tangent point for all ccdImages.
Definition: Associations.cc:87
photoCalib
Key< int > photoCalib
Definition: Exposure.cc:67
std::vector::emplace_back
T emplace_back(T... args)
lsst::sphgeom::ConvexPolygon::getBoundingCircle
Circle getBoundingCircle() const override
getBoundingCircle returns a bounding-circle for this region.
Definition: ConvexPolygon.cc:300
lsst::jointcal::Frame::inFrame
bool inFrame(double x, double y) const
inside?
Definition: Frame.cc:120
lsst.pex::exceptions::InvalidParameterError
Reports invalid arguments.
Definition: Runtime.h:66
lsst::jointcal::Frame::getArea
double getArea() const
Definition: Frame.cc:118
std::list::begin
T begin(T... args)
lsst::jointcal::JointcalControl::sourceFluxField
std::string sourceFluxField
"name of flux field in source catalog" ;
Definition: JointcalControl.h:40
detector
table::Key< int > detector
Definition: DetectorCollection.cc:172
key
Key< U > key
Definition: Schema.cc:281
lsst::geom::Point< double, 2 >
lsst::jointcal::JointcalControl
Definition: JointcalControl.h:39
lsst::geom::Angle
A class representing an angle.
Definition: Angle.h:127
std::count_if
T count_if(T... args)
lsst::geom::Box2I
An integer coordinate rectangle.
Definition: Box.h:55
LonLat.h
This file contains a class representing spherical coordinates.
lsst::jointcal::Point
A point in a plane.
Definition: Point.h:36
transform
table::Key< int > transform
Definition: TransformMap.cc:299
lsst::jointcal::FittedStarList
A list of FittedStar s. Such a list is typically constructed by Associations.
Definition: FittedStar.h:123
std::stringstream::str
T str(T... args)
lsst::jointcal::Associations::nCcdImagesValidForFit
int nCcdImagesValidForFit() const
return the number of CcdImages with non-empty catalogs to-be-fit.
Definition: Associations.cc:408
lsst::jointcal::CcdImage::getCatalogForFit
MeasuredStarList const & getCatalogForFit() const
Gets the catalog to be used for fitting, which may have been cleaned-up.
Definition: CcdImage.h:94
std::list::end
T end(T... args)
lsst::jointcal::Associations::createCcdImage
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.
Definition: Associations.cc:62
lsst::jointcal::BaseStar
The base class for handling stars. Used by all matching routines.
Definition: BaseStar.h:50
lsst::jointcal::AstrometryTransform
a virtual (interface) class for geometric transformations.
Definition: AstrometryTransform.h:65
lsst.pipe.tasks.assembleCoadd.visit
visit
Definition: assembleCoadd.py:343
lsst::sphgeom::ConvexPolygon::convexHull
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...
Definition: ConvexPolygon.h:65
lsst::jointcal::Associations::nFittedStarsWithAssociatedRefStar
size_t nFittedStarsWithAssociatedRefStar() const
Return the number of fittedStars that have an associated refStar.
Definition: Associations.cc:414
lsst::sphgeom::Circle
Circle is a circular region on the unit sphere that contains its boundary.
Definition: Circle.h:46
std::unique_ptr
STL class.
lsst::geom::Angle::asArcseconds
constexpr double asArcseconds() const noexcept
Return an Angle's value in arcseconds.
Definition: Angle.h:178
lsst::jointcal::Associations::computeBoundingCircle
lsst::sphgeom::Circle computeBoundingCircle() const
Return the bounding circle in on-sky (RA, Dec) coordinates containing all CcdImages.
Definition: Associations.cc:92
StarMatch.h
pairs of points
Log.h
LSST DM logging module built on log4cxx.
lsst::afw::table::SortedCatalogT
Custom catalog class for record/table subclasses that are guaranteed to have an ID,...
Definition: fwd.h:63
lsst::jointcal::Associations::computeCommonTangentPoint
void computeCommonTangentPoint()
Sets a shared tangent point for all ccdImages, using the mean of the centers of all ccdImages.
Definition: Associations.cc:76
visitInfo
Key< int > visitInfo
Definition: Exposure.cc:70
lsst::jointcal::Frame
rectangle with sides parallel to axes.
Definition: Frame.h:38
bbox
AmpInfoBoxKey bbox
Definition: Amplifier.cc:117
lsst::jointcal::FittedStarIterator
FittedStarList::iterator FittedStarIterator
Definition: FittedStar.h:132
ConvexPolygon.h
This file declares a class for representing convex polygons with great circle edges on the unit spher...
exceptions.h
VisitInfo.h
lsst::jointcal::BaseStarList
StarList< BaseStar > BaseStarList
Definition: BaseStar.h:119
std::pow
T pow(T... args)