LSSTApplications  19.0.0+10,19.0.0+80,19.0.0-1-g20d9b18+31,19.0.0-1-g49a97f9+4,19.0.0-1-g8c57eb9+31,19.0.0-1-g9a028c0+13,19.0.0-1-ga72da6b+3,19.0.0-1-gb77924a+15,19.0.0-1-gbfe0924+66,19.0.0-1-gc3516c3,19.0.0-1-gefe1d0d+49,19.0.0-1-gf8cb8b4+3,19.0.0-14-g7511ce4+6,19.0.0-16-g3dc1a33c+6,19.0.0-17-g59f0e8a+4,19.0.0-17-g9c22e3c+9,19.0.0-18-g35bb99870+2,19.0.0-19-g2772d4a+9,19.0.0-2-g260436e+53,19.0.0-2-g31cdcee,19.0.0-2-g9675b69+3,19.0.0-2-gaa2795f,19.0.0-2-gbcc4de1,19.0.0-2-gd6f004e+6,19.0.0-2-gde8e5e3+5,19.0.0-2-gff6972b+15,19.0.0-21-ge306cd8,19.0.0-21-gf122e69+2,19.0.0-3-g2d43a51+2,19.0.0-3-gf3b1435+6,19.0.0-4-g56feb96,19.0.0-4-gac56cce+17,19.0.0-49-gce872c1+1,19.0.0-5-g66946eb+6,19.0.0-5-gd8897ba+6,19.0.0-51-gfc4a647e,19.0.0-7-g686a884+5,19.0.0-7-g886f805+5,19.0.0-8-g305ff64,w.2020.17
LSSTDataManagementBasePackage
AssessSpatialKernelVisitor.cc
Go to the documentation of this file.
1 // -*- lsst-c++ -*-
12 #include "lsst/afw/math.h"
13 #include "lsst/afw/image.h"
14 #include "lsst/log/Log.h"
17 
21 
22 #define DEBUG_IMAGES 0
23 
24 namespace afwMath = lsst::afw::math;
25 namespace afwImage = lsst::afw::image;
26 namespace dafBase = lsst::daf::base;
28 
29 namespace lsst {
30 namespace ip {
31 namespace diffim {
32 namespace detail {
53  template<typename PixelT>
58  ) :
59  afwMath::CandidateVisitor(),
60  _spatialKernel(spatialKernel),
61  _spatialBackground(spatialBackground),
62  _ps(ps.deepCopy()),
63  _imstats(ImageStatistics<PixelT>(ps)),
64  _nGood(0),
65  _nRejected(0),
66  _nProcessed(0),
67  _useCoreStats(ps.getAsBool("useCoreStats")),
68  _coreRadius(ps.getAsInt("candidateCoreRadius"))
69  {};
70 
71  template<typename PixelT>
74  ) {
75 
76  KernelCandidate<PixelT> *kCandidate = dynamic_cast<KernelCandidate<PixelT> *>(candidate);
77  if (kCandidate == NULL) {
79  "Failed to cast SpatialCellCandidate to KernelCandidate");
80  }
81  if (!(kCandidate->isInitialized())) {
83  LOGL_DEBUG("TRACE2.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
84  "Cannot process candidate %d, continuing", kCandidate->getId());
85  return;
86  }
87 
88  LOGL_DEBUG("TRACE1.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
89  "Processing candidate %d", kCandidate->getId());
90 
91  /*
92  Note - this is a hack until the Kernel API is upgraded by the
93  Davis crew. I need a "local" version of the spatially varying
94  Kernel
95  */
96  afwImage::Image<double> kImage(_spatialKernel->getDimensions());
97  double kSum = _spatialKernel->computeImage(kImage, false,
98  kCandidate->getXCenter(), kCandidate->getYCenter());
100  kernelPtr(new afwMath::FixedKernel(kImage));
101  /* </hack> */
102 
103  double background = (*_spatialBackground)(kCandidate->getXCenter(), kCandidate->getYCenter());
104 
105  MaskedImageT diffim = kCandidate->getDifferenceImage(kernelPtr, background);
106 
107  if (DEBUG_IMAGES) {
108  kImage.writeFits(str(boost::format("askv_k%d.fits") % kCandidate->getId()));
109  diffim.writeFits(str(boost::format("askv_d%d.fits") % kCandidate->getId()));
110  }
111 
112  /* Official resids */
113  try {
114  if (_useCoreStats)
115  _imstats.apply(diffim, _coreRadius);
116  else
117  _imstats.apply(diffim);
118  } catch (pexExcept::Exception& e) {
119  LOGL_DEBUG("TRACE2.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
120  "Unable to calculate imstats for Candidate %d", kCandidate->getId());
122  return;
123  }
124 
125  _nProcessed += 1;
126 
127  LOGL_DEBUG("TRACE4.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
128  "Chi2 = %.3f", _imstats.getVariance());
129  LOGL_DEBUG("TRACE4.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
130  "X = %.2f Y = %.2f",
131  kCandidate->getXCenter(),
132  kCandidate->getYCenter());
133  LOGL_DEBUG("TRACE4.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
134  "Kernel Sum = %.3f", kSum);
135  LOGL_DEBUG("TRACE4.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
136  "Background = %.3f", background);
137  LOGL_DEBUG("TRACE2.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
138  "Candidate %d resids = %.3f +/- %.3f sigma (%d pix)",
139  kCandidate->getId(),
140  _imstats.getMean(),
141  _imstats.getRms(),
142  _imstats.getNpix());
143 
144  bool meanIsNan = std::isnan(_imstats.getMean());
145  bool rmsIsNan = std::isnan(_imstats.getRms());
146  if (meanIsNan || rmsIsNan) {
148  LOGL_DEBUG("TRACE3.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
149  "Rejecting candidate %d, encountered NaN",
150  kCandidate->getId());
151  _nRejected += 1;
152  return;
153  }
154 
155  if (_ps->getAsBool("spatialKernelClipping")) {
156  if (fabs(_imstats.getMean()) > _ps->getAsDouble("candidateResidualMeanMax")) {
158  LOGL_DEBUG("TRACE3.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
159  "Rejecting candidate %d; bad mean residual : |%.3f| > %.3f",
160  kCandidate->getId(),
161  _imstats.getMean(),
162  _ps->getAsDouble("candidateResidualMeanMax"));
163  _nRejected += 1;
164  }
165  else if (_imstats.getRms() > _ps->getAsDouble("candidateResidualStdMax")) {
167  LOGL_DEBUG("TRACE3.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
168  "Rejecting candidate %d; bad residual rms : %.3f > %.3f",
169  kCandidate->getId(),
170  _imstats.getRms(),
171  _ps->getAsDouble("candidateResidualStdMax"));
172  _nRejected += 1;
173  }
174  else {
176  LOGL_DEBUG("TRACE3.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
177  "Spatial kernel OK");
178  _nGood += 1;
179  }
180  }
181  else {
183  LOGL_DEBUG("TRACE5.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
184  "Sigma clipping not enabled");
185  _nGood += 1;
186  }
187 
188  /* Core resids for debugging */
189  if (!(_useCoreStats)) {
190  try {
191  _imstats.apply(diffim, _coreRadius);
192  } catch (pexExcept::Exception& e) {
193  LOGL_DEBUG("TRACE2.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
194  "Unable to calculate core imstats for Candidate %d",
195  kCandidate->getId());
197  return;
198  }
199  LOGL_DEBUG("TRACE3.ip.diffim.AssessSpatialKernelVisitor.processCandidate",
200  "Candidate %d core resids = %.3f +/- %.3f sigma (%d pix)",
201  kCandidate->getId(),
202  _imstats.getMean(),
203  _imstats.getRms(),
204  _imstats.getNpix());
205  }
206  }
207 
208  typedef float PixelT;
209  template class AssessSpatialKernelVisitor<PixelT>;
210 
211 }}}} // end of namespace lsst::ip::diffim::detail
def format(config, name=None, writeSourceLine=True, prefix="", verbose=False)
Definition: history.py:174
void setStatus(Status status)
Set the candidate&#39;s status.
Definition: SpatialCell.cc:54
afw::table::Key< afw::table::Array< ImagePixelT > > image
Asseses the quality of a candidate given a spatial kernel and background model.
Class stored in SpatialCells for spatial Kernel fitting.
float getXCenter() const
Return the object&#39;s column-centre.
Definition: SpatialCell.h:90
afw::image::MaskedImage< PixelT > getDifferenceImage(CandidateSwitch cand)
Calculate associated difference image using internal solutions.
Provides consistent interface for LSST exceptions.
Definition: Exception.h:107
int getId() const
Return the candidate&#39;s unique ID.
Definition: SpatialCell.h:104
#define LOGL_DEBUG(logger, message...)
Log a debug-level message using a varargs/printf style interface.
Definition: Log.h:504
LSST DM logging module built on log4cxx.
AssessSpatialKernelVisitor(std::shared_ptr< lsst::afw::math::LinearCombinationKernel > spatialKernel, lsst::afw::math::Kernel::SpatialFunctionPtr spatialBackground, lsst::daf::base::PropertySet const &ps)
A base class for image defects.
#define DEBUG_IMAGES
A class to manipulate images, masks, and variance as a single object.
Definition: MaskedImage.h:73
Reports errors in the logical structure of the program.
Definition: Runtime.h:46
Class used by SpatialModelCell for spatial Kernel fitting.
Class to calculate difference image statistics.
#define LSST_EXCEPT(type,...)
Create an exception with a given type.
Definition: Exception.h:48
float getYCenter() const
Return the object&#39;s row-centre.
Definition: SpatialCell.h:93
Class for storing generic metadata.
Definition: PropertySet.h:67
void processCandidate(lsst::afw::math::SpatialCellCandidate *candidate)
T isnan(T... args)
Declaration of AssessSpatialKernelVisitor.
Base class for candidate objects in a SpatialCell.
Definition: SpatialCell.h:70
Backwards-compatibility support for depersisting the old Calib (FluxMag0/FluxMag0Err) objects...
A class to represent a 2-dimensional array of pixels.
Definition: Image.h:58
Image Subtraction helper functions.
A kernel created from an Image.
Definition: Kernel.h:518