lsst.meas.extensions.photometryKron  15.0
KronPhotometry.cc
Go to the documentation of this file.
1 // -*- LSST-C++ -*-
2 /*
3  * LSST Data Management System
4  * Copyright 2008-2015 LSST Corporation.
5  *
6  * This product includes software developed by the
7  * LSST Project (http://www.lsst.org/).
8  *
9  * This program is free software: you can redistribute it and/or modify
10  * it under the terms of the GNU General Public License as published by
11  * the Free Software Foundation, either version 3 of the License, or
12  * (at your option) any later version.
13  *
14  * This program is distributed in the hope that it will be useful,
15  * but WITHOUT ANY WARRANTY; without even the implied warranty of
16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17  * GNU General Public License for more details.
18  *
19  * You should have received a copy of the LSST License Statement and
20  * the GNU General Public License along with this program. If not,
21  * see <http://www.lsstcorp.org/LegalNotices/>.
22  */
23 
24 #include <numeric>
25 #include <cmath>
26 #include <functional>
27 #include "boost/math/constants/constants.hpp"
28 #include "lsst/pex/exceptions.h"
29 #include "lsst/afw/geom/Point.h"
30 #include "lsst/afw/geom/Box.h"
31 #include "lsst/afw/geom/SpanSet.h"
33 #include "lsst/afw/table/Source.h"
37 #include "lsst/afw/detection/Psf.h"
38 #include "lsst/afw/coord/Coord.h"
40 #include "lsst/afw/geom/ellipses.h"
41 #include "lsst/meas/base.h"
43 
45 
46 namespace lsst {
47 namespace meas {
48 namespace extensions {
49 namespace photometryKron {
50 
51 namespace {
52 base::FlagDefinitionList flagDefinitions;
53 } // end anonymous
54 
55 base::FlagDefinition const KronFluxAlgorithm::FAILURE = flagDefinitions.addFailureFlag( "general failure flag, set if anything went wrong");
56 base::FlagDefinition const KronFluxAlgorithm::EDGE = flagDefinitions.add("flag_edge", "bad measurement due to image edge");
57 base::FlagDefinition const KronFluxAlgorithm::BAD_SHAPE_NO_PSF = flagDefinitions.add("flag_bad_shape_no_psf", "bad shape and no PSF");
58 base::FlagDefinition const KronFluxAlgorithm::NO_MINIMUM_RADIUS = flagDefinitions.add("flag_no_minimum_radius", "minimum radius could not enforced: no minimum value or PSF");
59 base::FlagDefinition const KronFluxAlgorithm::NO_FALLBACK_RADIUS = flagDefinitions.add("flag_no_fallback_radius", "no minimum radius and no PSF provided");
60 base::FlagDefinition const KronFluxAlgorithm::BAD_RADIUS = flagDefinitions.add("flag_bad_radius", "bad Kron radius");
61 base::FlagDefinition const KronFluxAlgorithm::USED_MINIMUM_RADIUS = flagDefinitions.add("flag_used_minimum_radius", "used the minimum radius for the Kron aperture");
62 base::FlagDefinition const KronFluxAlgorithm::USED_PSF_RADIUS = flagDefinitions.add("flag_used_psf_radius", "used the PSF Kron radius for the Kron aperture");
63 base::FlagDefinition const KronFluxAlgorithm::SMALL_RADIUS = flagDefinitions.add("flag_small_radius", "measured Kron radius was smaller than that of the PSF");
64 base::FlagDefinition const KronFluxAlgorithm::BAD_SHAPE = flagDefinitions.add("flag_bad_shape", "shape for measuring Kron radius is bad; used PSF shape");
65 
67  return flagDefinitions;
68 }
69 
72 
73 namespace {
74 
75 template <typename MaskedImageT>
76  class FootprintFlux {
77 public:
78  explicit FootprintFlux() : _sum(0.0), _sumVar(0.0) {}
79 
81  void reset() {
82  _sum = _sumVar = 0.0;
83  }
84  void reset(afw::detection::Footprint const&) {}
85 
87  void operator()(afw::geom::Point2I const & pos,
88  typename MaskedImageT::Image::Pixel const & ival,
89  typename MaskedImageT::Variance::Pixel const & vval) {
90  _sum += ival;
91  _sumVar += vval;
92  }
93 
95  double getSum() const { return _sum; }
96 
98  double getSumVar() const { return _sumVar; }
99 
100 private:
101  double _sum;
102  double _sumVar;
103 };
104 
105 /************************************************************************************************************/
116 template <typename MaskedImageT, typename WeightImageT>
117 class FootprintFindMoment {
118 public:
119  FootprintFindMoment(MaskedImageT const& mimage,
120  afw::geom::Point2D const& center, // center of the object
121  double const ab, // axis ratio
122  double const theta // rotation of ellipse +ve from x axis
123  ) : _xcen(center.getX()), _ycen(center.getY()),
124  _ab(ab),
125  _cosTheta(::cos(theta)),
126  _sinTheta(::sin(theta)),
127  _sum(0.0), _sumR(0.0),
128 #if 0
129  _sumVar(0.0), _sumRVar(0.0),
130 #endif
131  _imageX0(mimage.getX0()), _imageY0(mimage.getY0())
132  {}
133 
135  void reset() {}
136  void reset(afw::detection::Footprint const& foot) {
137  _sum = _sumR = 0.0;
138 #if 0
139  _sumVar = _sumRVar = 0.0;
140 #endif
141 
142  MaskedImageT const& mimage = this->getImage();
143  afw::geom::Box2I const& bbox(foot.getBBox());
144  int const x0 = bbox.getMinX(), y0 = bbox.getMinY(), x1 = bbox.getMaxX(), y1 = bbox.getMaxY();
145 
146  if (x0 < _imageX0 || y0 < _imageY0 ||
147  x1 >= _imageX0 + mimage.getWidth() || y1 >= _imageY0 + mimage.getHeight()) {
149  (boost::format("Footprint %d,%d--%d,%d doesn't fit in image %d,%d--%d,%d")
150  % x0 % y0 % x1 % y1
151  % _imageX0 % _imageY0
152  % (_imageX0 + mimage.getWidth() - 1) % (_imageY0 + mimage.getHeight() - 1)
153  ).str());
154  }
155  }
156 
158  void operator()(afw::geom::Point2I const & pos, typename MaskedImageT::Image::Pixel const & ival) {
159  double x = static_cast<double>(pos.getX());
160  double y = static_cast<double>(pos.getY());
161  double const dx = x - _xcen;
162  double const dy = y - _ycen;
163  double const du = dx*_cosTheta + dy*_sinTheta;
164  double const dv = -dx*_sinTheta + dy*_cosTheta;
165 
166  double r = ::hypot(du, dv*_ab); // ellipsoidal radius
167 #if 1
168  if (::hypot(dx, dy) < 0.5) { // within a pixel of the centre
169  /*
170  * We gain significant precision for flattened Gaussians by treating the central pixel specially
171  *
172  * If the object's centered in the pixel (and has constant surface brightness) we have <r> == eR;
173  * if it's at the corner <r> = 2*eR; we interpolate between these exact results linearily in the
174  * displacement. And then add in quadrature which is also a bit dubious
175  *
176  * We could avoid all these issues by estimating <r> using the same trick as we use for
177  * the sinc fluxes; it's not clear that it's worth it.
178  */
179 
180  double const eR = 0.38259771140356325; // <r> for a single square pixel, about the centre
181  r = ::hypot(r, eR*(1 + ::hypot(dx, dy)/afw::geom::ROOT2));
182  }
183 #endif
184 
185  _sum += ival;
186  _sumR += r*ival;
187 #if 0
188  typename MaskedImageT::Variance::Pixel vval = iloc.variance(0, 0);
189  _sumVar += vval;
190  _sumRVar += r*r*vval;
191 #endif
192  }
193 
195  double getIr() const { return _sumR/_sum; }
196 
197 #if 0
198 // double getIrVar() const { return _sumRVar/_sum - getIr()*getIr(); } // Wrong?
200  double getIrVar() const { return _sumRVar/(_sum*_sum) + _sumVar*_sumR*_sumR/::pow(_sum, 4); }
201 #endif
202 
204  bool getGood() const { return _sum > 0 && _sumR > 0; }
205 
206 private:
207  double const _xcen; // center of object
208  double const _ycen; // center of object
209  double const _ab; // axis ratio
210  double const _cosTheta, _sinTheta; // {cos,sin}(angle from x-axis)
211  double _sum; // sum of I
212  double _sumR; // sum of R*I
213 #if 0
214  double _sumVar; // sum of Var(I)
215  double _sumRVar; // sum of R*R*Var(I)
216 #endif
217  int const _imageX0, _imageY0; // origin of image we're measuring
218 
219 };
220 } // end anonymous namespace
221 
223  afw::geom::ellipses::Axes const& shape,
224  afw::geom::LinearTransform const& transformation,
225  double const radius
226  )
227 {
228  afw::geom::ellipses::Axes axes(shape);
229  axes.scale(radius/axes.getDeterminantRadius());
230  return axes.transform(transformation);
231 }
232 
233 template<typename ImageT>
235  ImageT const& image,
237  afw::geom::Point2D const& center,
238  KronFluxControl const& ctrl
239  )
240 {
241  //
242  // We might smooth the image because this is what SExtractor and Pan-STARRS do. But I don't see much gain
243  //
244  double const sigma = ctrl.smoothingSigma; // Gaussian width of smoothing sigma to apply
245  bool const smoothImage = sigma > 0;
246  int kSize = smoothImage ? 2*int(2*sigma) + 1 : 1;
247  afw::math::GaussianFunction1<afw::math::Kernel::Pixel> gaussFunc(smoothImage ? sigma : 100);
248  afw::math::SeparableKernel kernel(kSize, kSize, gaussFunc, gaussFunc);
249  bool const doNormalize = true, doCopyEdge = false;
250  afw::math::ConvolutionControl convCtrl(doNormalize, doCopyEdge);
251  double radius0 = axes.getDeterminantRadius();
252  double radius = std::numeric_limits<double>::quiet_NaN();
253  float radiusForRadius = std::nanf("");
254  for (int i = 0; i < ctrl.nIterForRadius; ++i) {
255  axes.scale(ctrl.nSigmaForRadius);
256  radiusForRadius = axes.getDeterminantRadius(); // radius we used to estimate R_K
257  //
258  // Build an elliptical Footprint of the proper size
259  //
261  afw::geom::ellipses::Ellipse(axes, center)));
262  afw::geom::Box2I bbox = !smoothImage ?
263  foot.getBBox() :
264  kernel.growBBox(foot.getBBox()); // the smallest bbox needed to convolve with Kernel
265  bbox.clip(image.getBBox());
266  ImageT subImage(image, bbox, afw::image::PARENT, smoothImage);
267  if (smoothImage) {
268  afw::math::convolve(subImage, ImageT(image, bbox, afw::image::PARENT, false), kernel, convCtrl);
269  }
270  //
271  // Find the desired first moment of the elliptical radius, which corresponds to the major axis.
272  //
273  FootprintFindMoment<ImageT, afw::detection::Psf::Image> iRFunctor(
274  subImage, center, axes.getA()/axes.getB(), axes.getTheta()
275  );
276 
277  try {
278  foot.getSpans()->applyFunctor(
279  iRFunctor, *(subImage.getImage()));
281  if (i == 0) {
282  LSST_EXCEPT_ADD(e, "Determining Kron aperture");
283  }
284  break; // use the radius we have
285  }
286 
287  if (!iRFunctor.getGood()) {
288  throw LSST_EXCEPT(BadKronException, "Bad integral defining Kron radius");
289  }
290 
291  radius = iRFunctor.getIr()*sqrt(axes.getB()/axes.getA());
292  if (radius <= radius0) {
293  break;
294  }
295  radius0 = radius;
296 
297  axes.scale(radius/axes.getDeterminantRadius()); // set axes to our current estimate of R_K
298  iRFunctor.reset();
299  }
300 
301  return std::make_shared<KronAperture>(center, axes, radiusForRadius);
302 }
303 
304 // Photometer an image with a particular aperture
305 template<typename ImageT>
307  ImageT const& image, // Image to measure
308  afw::geom::ellipses::Ellipse const& aperture, // Aperture in which to measure
309  double const maxSincRadius // largest radius that we use sinc apertures to measure
310  )
311 {
312  afw::geom::ellipses::Axes const& axes = aperture.getCore();
313  if (axes.getB() > maxSincRadius) {
314  FootprintFlux<ImageT> fluxFunctor;
315  auto spans = afw::geom::SpanSet::fromShape(aperture);
316  spans->applyFunctor(
317  fluxFunctor, *(image.getImage()), *(image.getVariance()));
318  return std::make_pair(fluxFunctor.getSum(), ::sqrt(fluxFunctor.getSumVar()));
319  }
320  try {
321  base::ApertureFluxResult fluxResult = base::ApertureFluxAlgorithm::computeSincFlux<float>(image, aperture);
322  return std::make_pair(fluxResult.flux, fluxResult.fluxSigma);
323  } catch(pex::exceptions::LengthError &e) {
324  LSST_EXCEPT_ADD(e, (boost::format("Measuring Kron flux for object at (%.3f, %.3f);"
325  " aperture radius %g,%g theta %g")
326  % aperture.getCenter().getX() % aperture.getCenter().getY()
327  % axes.getA() % axes.getB() % afw::geom::radToDeg(axes.getTheta())).str());
328  throw e;
329  }
330 }
331 
332 
334  CONST_PTR(afw::detection::Psf) const& psf, // PSF to measure
335  afw::geom::Point2D const& center, // Centroid of source on parent image
336  double smoothingSigma=0.0 // Gaussian sigma of smoothing applied
337  )
338 {
339  assert(psf);
340  double const radius = psf->computeShape(center).getDeterminantRadius();
341  // For a Gaussian N(0, sigma^2), the Kron radius is sqrt(pi/2)*sigma
342  return ::sqrt(afw::geom::PI/2)*::hypot(radius, std::max(0.0, smoothingSigma));
343 }
344 
345 template<typename ImageT>
347  ImageT const& image,
348  double const nRadiusForFlux,
349  double const maxSincRadius
350  ) const
351 {
352  afw::geom::ellipses::Axes axes(getAxes()); // Copy of ellipse core, so we can scale
353  axes.scale(nRadiusForFlux);
354  afw::geom::ellipses::Ellipse const ellip(axes, getCenter());
355 
356  return photometer(image, ellip, maxSincRadius);
357 }
358 
359 /************************************************************************************************************/
360 
367  KronFluxControl const & ctrl,
368  std::string const & name,
369  afw::table::Schema & schema,
370  daf::base::PropertySet & metadata
371 ) : _name(name),
372  _ctrl(ctrl),
373  _fluxResultKey(
374  meas::base::FluxResultKey::addFields(schema, name, "flux from Kron Flux algorithm")
375  ),
376  _radiusKey(schema.addField<float>(name + "_radius", "Kron radius (sqrt(a*b))")),
377  _radiusForRadiusKey(schema.addField<float>(name + "_radius_for_radius",
378  "radius used to estimate <radius> (sqrt(a*b))")),
379  _psfRadiusKey(schema.addField<float>(name + "_psf_radius", "Radius of PSF")),
380  _centroidExtractor(schema, name, true)
381 {
382  _flagHandler = meas::base::FlagHandler::addFields(schema, name, getFlagDefinitions());
383  metadata.add(name + "_nRadiusForFlux", ctrl.nRadiusForFlux);
384 }
385 
387  afw::table::SourceRecord & measRecord,
389 ) const {
390  _flagHandler.handleFailure(measRecord, error);
391 }
392 
393 void KronFluxAlgorithm::_applyAperture(
394  afw::table::SourceRecord & source,
395  afw::image::Exposure<float> const& exposure,
396  KronAperture const& aperture
397  ) const
398 {
399  double const rad = aperture.getAxes().getDeterminantRadius();
401  throw LSST_EXCEPT(
403  BAD_RADIUS.doc,
405  );
406  }
407 
409  try {
410  result = aperture.measureFlux(exposure.getMaskedImage(), _ctrl.nRadiusForFlux, _ctrl.maxSincRadius);
411  } catch (pex::exceptions::LengthError const& e) {
412  // We hit the edge of the image; there's no reasonable fallback or recovery
413  throw LSST_EXCEPT(
415  EDGE.doc,
416  EDGE.number
417  );
419  throw LSST_EXCEPT(
421  EDGE.doc,
422  EDGE.number
423  );
424  }
425 
426  // set the results in the source object
427  meas::base::FluxResult fluxResult;
428  fluxResult.flux = result.first;
429  fluxResult.fluxSigma = result.second;
430  source.set(_fluxResultKey, fluxResult);
431  source.set(_radiusKey, aperture.getAxes().getDeterminantRadius());
432  //
433  // REMINDER: In the old code, the psfFactor is calculated using getPsfFactor,
434  // and the values set for _fluxCorrectionKeys. See old meas_algorithms version.
435 }
436 
437 void KronFluxAlgorithm::_applyForced(
438  afw::table::SourceRecord & source,
439  afw::image::Exposure<float> const & exposure,
440  afw::geom::Point2D const & center,
441  afw::table::SourceRecord const & reference,
442  afw::geom::AffineTransform const & refToMeas
443  ) const
444 {
445  float const radius = reference.get(reference.getSchema().find<float>(_ctrl.refRadiusName).key);
446  KronAperture const aperture(reference, refToMeas, radius);
447  _applyAperture(source, exposure, aperture);
448  if (exposure.getPsf()) {
449  source.set(_psfRadiusKey, calculatePsfKronRadius(exposure.getPsf(), center, _ctrl.smoothingSigma));
450  }
451 }
452 
454  afw::table::SourceRecord & source,
455  afw::image::Exposure<float> const& exposure
456  ) const {
457  afw::geom::Point2D center = _centroidExtractor(source, _flagHandler);
458 
459  // Did we hit a condition that fundamentally prevented measuring the Kron flux?
460  // Such conditions include hitting the edge of the image and bad input shape, but not low signal-to-noise.
461  bool bad = false;
462 
463  afw::image::MaskedImage<float> const& mimage = exposure.getMaskedImage();
464 
465  double R_K_psf = -1;
466  if (exposure.getPsf()) {
467  R_K_psf = calculatePsfKronRadius(exposure.getPsf(), center, _ctrl.smoothingSigma);
468  }
469 
470  //
471  // Get the shape of the desired aperture
472  //
474  if (!source.getShapeFlag()) {
475  axes = source.getShape();
476  } else {
477  bad = true;
478  if (!exposure.getPsf()) {
479  throw LSST_EXCEPT(
483  );
484  }
485  axes = exposure.getPsf()->computeShape();
486  _flagHandler.setValue(source, BAD_SHAPE.number, true);
487  }
488  if (_ctrl.useFootprintRadius) {
489  afw::geom::ellipses::Axes footprintAxes(source.getFootprint()->getShape());
490  // if the Footprint's a disk of radius R we want footRadius == R.
491  // As <r^2> = R^2/2 for a disk, we need to scale up by sqrt(2)
492  footprintAxes.scale(::sqrt(2));
493 
494  double radius0 = axes.getDeterminantRadius();
495  double const footRadius = footprintAxes.getDeterminantRadius();
496 
497  if (footRadius > radius0*_ctrl.nSigmaForRadius) {
498  radius0 = footRadius/_ctrl.nSigmaForRadius; // we'll scale it up by nSigmaForRadius
499  axes.scale(radius0/axes.getDeterminantRadius());
500  }
501  }
502 
503  PTR(KronAperture) aperture;
504  if (_ctrl.fixed) {
505  aperture.reset(new KronAperture(source));
506  } else {
507  try {
508  aperture = KronAperture::determineRadius(mimage, axes, center, _ctrl);
509  } catch (pex::exceptions::OutOfRangeError& e) {
510  // We hit the edge of the image: no reasonable fallback or recovery possible
511  throw LSST_EXCEPT(
513  EDGE.doc,
514  EDGE.number
515  );
516  } catch (BadKronException& e) {
517  // Not setting bad=true because we only failed due to low S/N
518  aperture = _fallbackRadius(source, R_K_psf, e);
519  } catch(pex::exceptions::Exception& e) {
520  bad = true; // There's something fundamental keeping us from measuring the Kron aperture
521  aperture = _fallbackRadius(source, R_K_psf, e);
522  }
523  }
524 
525  /*
526  * Estimate the minimum acceptable Kron radius as the Kron radius of the PSF or the
527  * provided minimum radius
528  */
529 
530  // Enforce constraints on minimum radius
531  double rad = aperture->getAxes().getDeterminantRadius();
532  if (_ctrl.enforceMinimumRadius) {
533  double newRadius = rad;
534  if (_ctrl.minimumRadius > 0.0) {
535  if (rad < _ctrl.minimumRadius) {
536  newRadius = _ctrl.minimumRadius;
537  _flagHandler.setValue(source, USED_MINIMUM_RADIUS.number, true);
538  }
539  } else if (!exposure.getPsf()) {
540  throw LSST_EXCEPT(
544  );
545  } else if (rad < R_K_psf) {
546  newRadius = R_K_psf;
547  _flagHandler.setValue(source, USED_PSF_RADIUS.number, true);
548  }
549  if (newRadius != rad) {
550  aperture->getAxes().scale(newRadius/rad);
551  _flagHandler.setValue(source, SMALL_RADIUS.number, true); // guilty after all
552  }
553  }
554 
555  _applyAperture(source, exposure, *aperture);
556  source.set(_radiusForRadiusKey, aperture->getRadiusForRadius());
557  source.set(_psfRadiusKey, R_K_psf);
558  if (bad) _flagHandler.setValue(source, FAILURE.number, true);
559 }
560 
562  afw::table::SourceRecord & measRecord,
563  afw::image::Exposure<float> const & exposure,
564  afw::table::SourceRecord const & refRecord,
565  afw::geom::SkyWcs const & refWcs
566  ) const {
567  afw::geom::Point2D center = _centroidExtractor(measRecord, _flagHandler);
568  auto xytransform = afw::geom::makeWcsPairTransform(refWcs, *exposure.getWcs());
569  _applyForced(measRecord, exposure, center, refRecord,
570  linearizeTransform(*xytransform, refRecord.getCentroid())
571  );
572 
573 }
574 
575 
576 PTR(KronAperture) KronFluxAlgorithm::_fallbackRadius(afw::table::SourceRecord& source, double const R_K_psf,
577  pex::exceptions::Exception& exc) const
578 {
579  _flagHandler.setValue(source, BAD_RADIUS.number, true);
580  double newRadius;
581  if (_ctrl.minimumRadius > 0) {
582  newRadius = _ctrl.minimumRadius;
583  _flagHandler.setValue(source, USED_MINIMUM_RADIUS.number, true);
584  } else if (R_K_psf > 0) {
585  newRadius = R_K_psf;
586  _flagHandler.setValue(source, USED_PSF_RADIUS.number, true);
587  } else {
588  throw LSST_EXCEPT(
592  );
593  }
594  PTR(KronAperture) aperture(new KronAperture(source));
595  aperture->getAxes().scale(newRadius/aperture->getAxes().getDeterminantRadius());
596  return aperture;
597 }
598 
599 
600 #define INSTANTIATE(TYPE) \
601 template PTR(KronAperture) KronAperture::determineRadius<afw::image::MaskedImage<TYPE> >( \
602  afw::image::MaskedImage<TYPE> const&, \
603  afw::geom::ellipses::Axes, \
604  afw::geom::Point2D const&, \
605  KronFluxControl const& \
606  ); \
607 template std::pair<double, double> KronAperture::measureFlux<afw::image::MaskedImage<TYPE> >( \
608  afw::image::MaskedImage<TYPE> const&, \
609  double const, \
610  double const \
611  ) const;
612 
613 INSTANTIATE(float);
614 
615 }}}} // namespace lsst::meas::extensions::photometryKron
double const getB() const
static meas::base::FlagDefinition const SMALL_RADIUS
static meas::base::FlagDefinition const NO_MINIMUM_RADIUS
std::shared_ptr< TransformPoint2ToPoint2 > makeWcsPairTransform(SkyWcs const &src, SkyWcs const &dst)
static std::shared_ptr< geom::SpanSet > fromShape(int r, Stencil s=Stencil::CIRCLE, Point2I offset=Point2I())
virtual void measureForced(afw::table::SourceRecord &measRecord, afw::image::Exposure< float > const &exposure, afw::table::SourceRecord const &refRecord, afw::geom::SkyWcs const &refWcs) const
ShapeSlotDefinition::MeasValue getShape() const
std::string const & _name
std::pair< double, double > photometer(ImageT const &image, afw::geom::ellipses::Ellipse const &aperture, double const maxSincRadius)
double smoothingSigma
"Smooth image with N(0, smoothingSigma^2) Gaussian while estimating R_K" ;
double const getA() const
double calculatePsfKronRadius(boost::shared_ptr< afw::detection::Psf const > const &psf, afw::geom::Point2D const &center, double smoothingSigma=0.0)
double minimumRadius
"Minimum Kron radius (if == 0.0 use PSF&#39;s Kron radius) if enforceMinimumRadius. " "Also functions as ...
std::shared_ptr< Footprint > getFootprint() const
BaseCore const & getCore() const
#define PTR(...)
SchemaItem< T > find(std::string const &name) const
Transformer transform(LinearTransform const &transform)
KronFluxAlgorithm(Control const &ctrl, std::string const &name, afw::table::Schema &schema, daf::base::PropertySet &metadata)
A class that knows how to calculate fluxes using the KRON photometry algorithm.
bool fixed
"if true, use existing shape and centroid measurements instead of fitting" ;
void setValue(afw::table::BaseRecord &record, std::size_t i, bool value) const
static boost::shared_ptr< KronAperture > determineRadius(ImageT const &image, afw::geom::ellipses::Axes axes, afw::geom::Point2D const &center, KronFluxControl const &ctrl)
Determine the Kron Aperture from an image.
Field< T >::Value get(Key< T > const &key) const
#define CONST_PTR(...)
afw::table::Key< double > sigma
STL class.
AffineTransform linearizeTransform(TransformPoint2ToPoint2 const &original, Point2D const &inPoint)
static meas::base::FlagDefinition const USED_MINIMUM_RADIUS
double constexpr PI
static meas::base::FlagDefinition const USED_PSF_RADIUS
MaskedImageT getMaskedImage()
static afw::geom::ellipses::Axes getKronAxes(afw::geom::ellipses::Axes const &shape, afw::geom::LinearTransform const &transformation, double const radius)
Determine Kron axes from a reference image.
static meas::base::FlagDefinition const NO_FALLBACK_RADIUS
static meas::base::FlagDefinition const BAD_SHAPE_NO_PSF
double nSigmaForRadius
"Multiplier of rms size for aperture used to initially estimate the Kron radius" ; ...
bool enforceMinimumRadius
"If true check that the Kron radius exceeds some minimum" ;
T make_pair(T... args)
std::pair< double, double > measureFlux(ImageT const &image, double const nRadiusForFlux, double const maxSincRadius) const
Photometer within the Kron Aperture on an image.
double constexpr ROOT2
std::string refRadiusName
"Name of field specifying reference Kron radius for forced measurement" ;
double x
void convolve(OutImageT &convolvedImage, InImageT const &inImage, KernelT const &kernel, ConvolutionControl const &convolutionControl=ConvolutionControl())
T max(T... args)
static FlagHandler addFields(afw::table::Schema &schema, std::string const &prefix, FlagDefinitionList const &flagDefs, FlagDefinitionList const &exclDefs=FlagDefinitionList::getEmptyList())
static meas::base::FlagDefinition const FAILURE
#define LSST_EXCEPT(type,...)
geom::ellipses::Quadrupole computeShape(geom::Point2D position=makeNullPoint(), image::Color color=image::Color()) const
double maxSincRadius
"Largest aperture for which to use the slow, accurate, sinc aperture code" ;
static meas::base::FlagDefinition const BAD_RADIUS
Point2D const & getCenter() const
std::shared_ptr< lsst::afw::detection::Psf > getPsf()
int nIterForRadius
"Number of times to iterate when setting the Kron radius" ;
bool useFootprintRadius
"Use the Footprint size as part of initial estimate of Kron radius" ;
static meas::base::FlagDefinition const EDGE
T nanf(T... args)
constexpr double radToDeg(double x) noexcept
void handleFailure(afw::table::BaseRecord &record, MeasurementError const *error=nullptr) const
static meas::base::FlagDefinitionList const & getFlagDefinitions()
std::shared_ptr< geom::SkyWcs const > getWcs() const
void set(Key< T > const &key, U const &value)
table::Box2IKey bbox
std::shared_ptr< geom::SpanSet > getSpans() const
virtual void fail(afw::table::SourceRecord &measRecord, meas::base::MeasurementError *error=NULL) const
geom::Box2I getBBox() const
void clip(Box2I const &other)
T quiet_NaN(T... args)
double const getTheta() const
void add(std::string const &name, T const &value)
lsst::afw::geom::Box2I growBBox(lsst::afw::geom::Box2I const &bbox) const
virtual void measure(afw::table::SourceRecord &measRecord, afw::image::Exposure< float > const &exposure) const
#define LSST_EXCEPTION_TYPE(t, b, c)
#define LSST_EXCEPT_ADD(e, m)
CentroidSlotDefinition::MeasValue getCentroid() const
double nRadiusForFlux
"Number of Kron radii for Kron flux" ;
#define INSTANTIATE(TYPE)
static meas::base::FlagDefinition const BAD_SHAPE