lsst.meas.base g72c761de85+f36129b12f
SdssCentroid.cc
Go to the documentation of this file.
1// -*- lsst-c++ -*-
2/*
3 * LSST Data Management System
4 * Copyright 2008-2013 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#include <iostream>
24#include <cmath>
25#include <numeric>
26#include "ndarray/eigen.h"
27#include "lsst/geom/Angle.h"
29#include "lsst/pex/exceptions.h"
34
35namespace lsst {
36namespace meas {
37namespace base {
38namespace {
39FlagDefinitionList flagDefinitions;
40} // namespace
41
42FlagDefinition const SdssCentroidAlgorithm::FAILURE = flagDefinitions.addFailureFlag();
43FlagDefinition const SdssCentroidAlgorithm::EDGE =
44 flagDefinitions.add("flag_edge", "Object too close to edge");
46 flagDefinitions.add("flag_noSecondDerivative", "Vanishing second derivative");
48 flagDefinitions.add("flag_almostNoSecondDerivative", "Almost vanishing second derivative");
49FlagDefinition const SdssCentroidAlgorithm::NOT_AT_MAXIMUM =
50 flagDefinitions.add("flag_notAtMaximum", "Object is not at a maximum");
51
53
54namespace {
55
56/************************************************************************************************************/
57
58float const AMPAST4 = 1.33; // amplitude of `4th order' corr compared to theory
59
60/*
61 * Do the Gaussian quartic interpolation for the position
62 * of the maximum for three equally spaced values vm,v0,vp, assumed to be at
63 * abscissae -1,0,1; the answer is returned as *cen
64 *
65 * Return 0 is all is well, otherwise 1
66 */
67static int inter4(float vm, float v0, float vp, float *cen, bool negative = false) {
68 float const sp = v0 - vp;
69 float const sm = v0 - vm;
70 float const d2 = sp + sm;
71 float const s = 0.5 * (vp - vm);
72
73 if ((!negative && (d2 <= 0.0f || v0 <= 0.0f)) || (negative && (d2 >= 0.0f || v0 >= 0.0f))) {
74 return (1);
75 }
76
77 *cen = s / d2 * (1.0 + AMPAST4 * sp * sm / (d2 * v0));
78
79 return fabs(*cen) < 1 ? 0 : 1;
80}
81
82/*****************************************************************************/
83/*
84 * Calculate error in centroid
85 */
86float astrom_errors(float skyVar, // variance of pixels at the sky level
87 float sourceVar, // variance in peak due to excess counts over sky
88 float A, // abs(peak value in raw image)
89 float tau2, // Object is N(0,tau2)
90 float As, // abs(peak value in smoothed image)
91 float s, // slope across central pixel
92 float d, // curvature at central pixel
93 float sigma, // width of smoothing filter
94 int quarticBad) { // was quartic estimate bad?
95
96 double const k = quarticBad ? 0 : AMPAST4; /* quartic correction coeff */
97 double const sigma2 = sigma * sigma; /* == sigma^2 */
98 double sVar, dVar; /* variances of s and d */
99 double xVar; /* variance of centroid, x */
100
102 return (1e3);
103 }
104
105 if (sigma <= 0) { /* no smoothing; no covariance */
106 sVar = 0.5 * skyVar; /* due to sky */
107 dVar = 6 * skyVar;
108
109 sVar += 0.5 * sourceVar * exp(-1 / (2 * tau2));
110 dVar += sourceVar * (4 * exp(-1 / (2 * tau2)) + 2 * exp(-1 / (2 * tau2)));
111 } else { /* smoothed */
112 sVar = skyVar / (8 * geom::PI * sigma2) * (1 - exp(-1 / sigma2));
113 dVar = skyVar / (2 * geom::PI * sigma2) * (3 - 4 * exp(-1 / (4 * sigma2)) + exp(-1 / sigma2));
114
115 sVar += sourceVar / (12 * geom::PI * sigma2) * (exp(-1 / (3 * sigma2)) - exp(-1 / sigma2));
116 dVar += sourceVar / (3 * geom::PI * sigma2) * (2 - 3 * exp(-1 / (3 * sigma2)) + exp(-1 / sigma2));
117 }
118
119 xVar = sVar * pow(1 / d + k / (4 * As) * (1 - 12 * s * s / (d * d)), 2) +
120 dVar * pow(s / (d * d) - k / (4 * As) * 8 * s * s / (d * d * d), 2);
121
122 return (xVar >= 0 ? sqrt(xVar) : NAN);
123}
124
125/************************************************************************************************************/
126/*
127 * Estimate the position of an object, assuming we know that it's approximately the size of the PSF
128 */
129
130template <typename ImageXy_locatorT, typename VarImageXy_locatorT>
131void doMeasureCentroidImpl(double *xCenter, // output; x-position of object
132 double *dxc, // output; error in xCenter
133 double *yCenter, // output; y-position of object
134 double *dyc, // output; error in yCenter
135 double *sizeX2, double *sizeY2, // output; object widths^2 in x and y directions
136 ImageXy_locatorT im, // Locator for the pixel values
137 VarImageXy_locatorT vim, // Locator for the image containing the variance
138 double smoothingSigma, // Gaussian sigma of already-applied smoothing filter
139 FlagHandler flagHandler) {
140 /*
141 * find a first quadratic estimate
142 */
143 double const d2x = 2 * im(0, 0) - im(-1, 0) - im(1, 0);
144 double const d2y = 2 * im(0, 0) - im(0, -1) - im(0, 1);
145 double const sx = 0.5 * (im(1, 0) - im(-1, 0));
146 double const sy = 0.5 * (im(0, 1) - im(0, -1));
147
148 if (d2x == 0.0 || d2y == 0.0) {
151 }
152 if (d2x < 0.0 || d2y < 0.0) {
153 throw LSST_EXCEPT(MeasurementError,
155 (boost::format(": d2I/dx2, d2I/dy2 = %g %g") % d2x % d2y).str(),
157 }
158
159 double const dx0 = sx / d2x;
160 double const dy0 = sy / d2y; // first guess
161
162 if (fabs(dx0) > 10.0 || fabs(dy0) > 10.0) {
163 throw LSST_EXCEPT(
164 MeasurementError,
166 (boost::format(": sx, d2x, sy, d2y = %f %f %f %f") % sx % d2x % sy % d2y).str(),
168 }
169
170 double vpk = im(0, 0) + 0.5 * (sx * dx0 + sy * dy0); // height of peak in image
171 if (vpk < 0) {
172 vpk = -vpk;
173 }
174 /*
175 * now evaluate maxima on stripes
176 */
177 float m0x = 0, m1x = 0, m2x = 0;
178 float m0y = 0, m1y = 0, m2y = 0;
179
180 int quarticBad = 0;
181 quarticBad += inter4(im(-1, -1), im(0, -1), im(1, -1), &m0x);
182 quarticBad += inter4(im(-1, 0), im(0, 0), im(1, 0), &m1x);
183 quarticBad += inter4(im(-1, 1), im(0, 1), im(1, 1), &m2x);
184
185 quarticBad += inter4(im(-1, -1), im(-1, 0), im(-1, 1), &m0y);
186 quarticBad += inter4(im(0, -1), im(0, 0), im(0, 1), &m1y);
187 quarticBad += inter4(im(1, -1), im(1, 0), im(1, 1), &m2y);
188
189 double xc, yc; // position of maximum
190 double sigmaX2, sigmaY2; // widths^2 in x and y of smoothed object
191
192 if (quarticBad) { // >= 1 quartic interpolator is bad
193 xc = dx0;
194 yc = dy0;
195 sigmaX2 = vpk / d2x; // widths^2 in x
196 sigmaY2 = vpk / d2y; // and y
197 } else {
198 double const smx = 0.5 * (m2x - m0x);
199 double const smy = 0.5 * (m2y - m0y);
200 double const dm2x = m1x - 0.5 * (m0x + m2x);
201 double const dm2y = m1y - 0.5 * (m0y + m2y);
202 double const dx = m1x + dy0 * (smx - dy0 * dm2x); // first quartic approx
203 double const dy = m1y + dx0 * (smy - dx0 * dm2y);
204 double const dx4 = m1x + dy * (smx - dy * dm2x); // second quartic approx
205 double const dy4 = m1y + dx * (smy - dx * dm2y);
206
207 xc = dx4;
208 yc = dy4;
209 sigmaX2 = vpk / d2x - (1 + 6 * dx0 * dx0) / 4; // widths^2 in x
210 sigmaY2 = vpk / d2y - (1 + 6 * dy0 * dy0) / 4; // and y
211 }
212 /*
213 * Now for the errors.
214 */
215 float tauX2 = sigmaX2; // width^2 of _un_ smoothed object
216 float tauY2 = sigmaY2;
217 tauX2 -= smoothingSigma * smoothingSigma; // correct for smoothing
218 tauY2 -= smoothingSigma * smoothingSigma;
219
220 if (tauX2 <= smoothingSigma * smoothingSigma) { // problem; sigmaX2 must be bad
221 tauX2 = smoothingSigma * smoothingSigma;
222 }
223 if (tauY2 <= smoothingSigma * smoothingSigma) { // sigmaY2 must be bad
224 tauY2 = smoothingSigma * smoothingSigma;
225 }
226
227 float const skyVar = (vim(-1, -1) + vim(0, -1) + vim(1, -1) + vim(-1, 0) + vim(1, 0) + vim(-1, 1) +
228 vim(0, 1) + vim(1, 1)) /
229 8.0; // Variance in sky
230 float const sourceVar = vim(0, 0); // extra variance of peak due to its photons
231 float const A = vpk * sqrt((sigmaX2 / tauX2) * (sigmaY2 / tauY2)); // peak of Unsmoothed object
232
233 *xCenter = xc;
234 *yCenter = yc;
235
236 *dxc = astrom_errors(skyVar, sourceVar, A, tauX2, vpk, sx, d2x, fabs(smoothingSigma), quarticBad);
237 *dyc = astrom_errors(skyVar, sourceVar, A, tauY2, vpk, sy, d2y, fabs(smoothingSigma), quarticBad);
238
239 *sizeX2 = tauX2; // return the estimates of the (object size)^2
240 *sizeY2 = tauY2;
241}
242
243template <typename MaskedImageXy_locatorT>
244void doMeasureCentroidImpl(double *xCenter, // output; x-position of object
245 double *dxc, // output; error in xCenter
246 double *yCenter, // output; y-position of object
247 double *dyc, // output; error in yCenter
248 double *sizeX2, double *sizeY2, // output; object widths^2 in x and y directions
249 double *peakVal, // output; peak of object
250 MaskedImageXy_locatorT mim, // Locator for the pixel values
251 double smoothingSigma, // Gaussian sigma of already-applied smoothing filter
252 bool negative, FlagHandler flagHandler) {
253 /*
254 * find a first quadratic estimate
255 */
256 double const d2x = 2 * mim.image(0, 0) - mim.image(-1, 0) - mim.image(1, 0);
257 double const d2y = 2 * mim.image(0, 0) - mim.image(0, -1) - mim.image(0, 1);
258 double const sx = 0.5 * (mim.image(1, 0) - mim.image(-1, 0));
259 double const sy = 0.5 * (mim.image(0, 1) - mim.image(0, -1));
260
261 if (d2x == 0.0 || d2y == 0.0) {
264 }
265 if ((!negative && (d2x < 0.0 || d2y < 0.0)) || (negative && (d2x > 0.0 || d2y > 0.0))) {
266 throw LSST_EXCEPT(MeasurementError,
268 (boost::format(": d2I/dx2, d2I/dy2 = %g %g") % d2x % d2y).str(),
270 }
271
272 double const dx0 = sx / d2x;
273 double const dy0 = sy / d2y; // first guess
274
275 if (fabs(dx0) > 10.0 || fabs(dy0) > 10.0) {
276 throw LSST_EXCEPT(
277 MeasurementError,
279 (boost::format(": sx, d2x, sy, d2y = %f %f %f %f") % sx % d2x % sy % d2y).str(),
281 }
282
283 double vpk = mim.image(0, 0) + 0.5 * (sx * dx0 + sy * dy0); // height of peak in image
284 // if (vpk < 0) {
285 // vpk = -vpk;
286 //}
287 /*
288 * now evaluate maxima on stripes
289 */
290 float m0x = 0, m1x = 0, m2x = 0;
291 float m0y = 0, m1y = 0, m2y = 0;
292
293 int quarticBad = 0;
294 quarticBad += inter4(mim.image(-1, -1), mim.image(0, -1), mim.image(1, -1), &m0x, negative);
295 quarticBad += inter4(mim.image(-1, 0), mim.image(0, 0), mim.image(1, 0), &m1x, negative);
296 quarticBad += inter4(mim.image(-1, 1), mim.image(0, 1), mim.image(1, 1), &m2x, negative);
297
298 quarticBad += inter4(mim.image(-1, -1), mim.image(-1, 0), mim.image(-1, 1), &m0y, negative);
299 quarticBad += inter4(mim.image(0, -1), mim.image(0, 0), mim.image(0, 1), &m1y, negative);
300 quarticBad += inter4(mim.image(1, -1), mim.image(1, 0), mim.image(1, 1), &m2y, negative);
301
302 double xc, yc; // position of maximum
303 double sigmaX2, sigmaY2; // widths^2 in x and y of smoothed object
304
305 if (quarticBad) { // >= 1 quartic interpolator is bad
306 xc = dx0;
307 yc = dy0;
308 sigmaX2 = vpk / d2x; // widths^2 in x
309 sigmaY2 = vpk / d2y; // and y
310 } else {
311 double const smx = 0.5 * (m2x - m0x);
312 double const smy = 0.5 * (m2y - m0y);
313 double const dm2x = m1x - 0.5 * (m0x + m2x);
314 double const dm2y = m1y - 0.5 * (m0y + m2y);
315 double const dx = m1x + dy0 * (smx - dy0 * dm2x); // first quartic approx
316 double const dy = m1y + dx0 * (smy - dx0 * dm2y);
317 double const dx4 = m1x + dy * (smx - dy * dm2x); // second quartic approx
318 double const dy4 = m1y + dx * (smy - dx * dm2y);
319
320 xc = dx4;
321 yc = dy4;
322 sigmaX2 = vpk / d2x - (1 + 6 * dx0 * dx0) / 4; // widths^2 in x
323 sigmaY2 = vpk / d2y - (1 + 6 * dy0 * dy0) / 4; // and y
324 }
325 /*
326 * Now for the errors.
327 */
328 float tauX2 = sigmaX2; // width^2 of _un_ smoothed object
329 float tauY2 = sigmaY2;
330 tauX2 -= smoothingSigma * smoothingSigma; // correct for smoothing
331 tauY2 -= smoothingSigma * smoothingSigma;
332
333 if (tauX2 <= smoothingSigma * smoothingSigma) { // problem; sigmaX2 must be bad
334 tauX2 = smoothingSigma * smoothingSigma;
335 }
336 if (tauY2 <= smoothingSigma * smoothingSigma) { // sigmaY2 must be bad
337 tauY2 = smoothingSigma * smoothingSigma;
338 }
339
340 float const skyVar =
341 (mim.variance(-1, -1) + mim.variance(0, -1) + mim.variance(1, -1) + mim.variance(-1, 0) +
342 mim.variance(1, 0) + mim.variance(-1, 1) + mim.variance(0, 1) + mim.variance(1, 1)) /
343 8.0; // Variance in sky
344 float const sourceVar = mim.variance(0, 0); // extra variance of peak due to its photons
345 float const A = vpk * sqrt((sigmaX2 / tauX2) * (sigmaY2 / tauY2)); // peak of Unsmoothed object
346
347 *xCenter = xc;
348 *yCenter = yc;
349
350 *dxc = astrom_errors(skyVar, sourceVar, A, tauX2, vpk, sx, d2x, fabs(smoothingSigma), quarticBad);
351 *dyc = astrom_errors(skyVar, sourceVar, A, tauY2, vpk, sy, d2y, fabs(smoothingSigma), quarticBad);
352
353 *sizeX2 = tauX2; // return the estimates of the (object size)^2
354 *sizeY2 = tauY2;
355
356 *peakVal = vpk;
357}
358
359template <typename MaskedImageT>
361 const int y, MaskedImageT const &mimage, int binX, int binY,
362 FlagHandler _flagHandler) {
363 geom::Point2D const center(x + mimage.getX0(), y + mimage.getY0());
364 afw::geom::ellipses::Quadrupole const &shape = psf->computeShape(center);
365 double const smoothingSigma = shape.getDeterminantRadius();
366#if 0
367 double const nEffective = psf->computeEffectiveArea(); // not implemented yet (#2821)
368#else
369 double const nEffective = 4 * M_PI * smoothingSigma * smoothingSigma; // correct for a Gaussian
370#endif
371
372 std::shared_ptr<afw::math::Kernel const> kernel = psf->getLocalKernel(center);
373 int const kWidth = kernel->getWidth();
374 int const kHeight = kernel->getHeight();
375
376 geom::BoxI bbox(geom::Point2I(x - binX * (2 + kWidth / 2), y - binY * (2 + kHeight / 2)),
377 geom::ExtentI(binX * (3 + kWidth + 1), binY * (3 + kHeight + 1)));
378
379 // image to smooth, a shallow copy
381 try {
382 subImage.reset(new MaskedImageT(mimage, bbox, afw::image::LOCAL));
383 } catch (pex::exceptions::LengthError &err) {
384 throw LSST_EXCEPT(MeasurementError, SdssCentroidAlgorithm::EDGE.doc,
386 }
387 std::shared_ptr<MaskedImageT> binnedImage = afw::math::binImage(*subImage, binX, binY, afw::math::MEAN);
388 binnedImage->setXY0(subImage->getXY0());
389 // image to smooth into, a deep copy.
390 MaskedImageT smoothedImage = MaskedImageT(*binnedImage, true);
391 assert(smoothedImage.getWidth() / 2 == kWidth / 2 + 2); // assumed by the code that uses smoothedImage
392 assert(smoothedImage.getHeight() / 2 == kHeight / 2 + 2);
393
394 afw::math::convolve(smoothedImage, *binnedImage, *kernel, afw::math::ConvolutionControl());
395 *smoothedImage.getVariance() *= binX * binY * nEffective; // We want the per-pixel variance, so undo the
396 // effects of binning and smoothing
397
398 return std::make_pair(smoothedImage, smoothingSigma);
399}
400
401} // end anonymous namespace
402
404 afw::table::Schema &schema)
405 : _ctrl(ctrl),
406 _centroidKey(CentroidResultKey::addFields(schema, name, "centroid from Sdss Centroid algorithm",
407 SIGMA_ONLY)),
408 _flagHandler(FlagHandler::addFields(schema, name, getFlagDefinitions())),
409 _centroidExtractor(schema, name, true),
410 _centroidChecker(schema, name, ctrl.doFootprintCheck, ctrl.maxDistToPeak) {}
412 afw::image::Exposure<float> const &exposure) const {
413 // get our current best guess about the centroid: either a centroider measurement or peak.
414 geom::Point2D center = _centroidExtractor(measRecord, _flagHandler);
415 CentroidResult result;
416 result.x = center.getX();
417 result.y = center.getY();
418 measRecord.set(_centroidKey, result); // better than NaN
419
421 typedef MaskedImageT::Image ImageT;
422 typedef MaskedImageT::Variance VarianceT;
423 bool negative = false;
424 try {
425 negative = measRecord.get(measRecord.getSchema().find<afw::table::Flag>("flags_negative").key);
426 } catch (pexExcept::Exception &e) {
427 }
428
429 MaskedImageT const &mimage = exposure.getMaskedImage();
430 ImageT const &image = *mimage.getImage();
432
433 int const x = image.positionToIndex(center.getX(), afw::image::X).first;
434 int const y = image.positionToIndex(center.getY(), afw::image::Y).first;
435
436 if (!image.getBBox().contains(geom::Extent2I(x, y) + image.getXY0())) {
438 }
439
440 // Algorithm uses a least-squares fit (implemented via a convolution) to a symmetrized PSF model.
441 // If you don't have a Psf, you need to use another centroider, such as GaussianCentroider.
442 if (!psf) {
443 throw LSST_EXCEPT(FatalAlgorithmError, "SdssCentroid algorithm requires a Psf with every exposure");
444 }
445
446 int binX = 1;
447 int binY = 1;
448 double xc = 0., yc = 0., dxc = 0., dyc = 0.; // estimated centre and error therein
449 for (int binsize = 1; binsize <= _ctrl.binmax; binsize *= 2) {
451 smoothAndBinImage(psf, x, y, mimage, binX, binY, _flagHandler);
452 MaskedImageT const smoothedImage = result.first;
453 double const smoothingSigma = result.second;
454
455 MaskedImageT::xy_locator mim =
456 smoothedImage.xy_at(smoothedImage.getWidth() / 2, smoothedImage.getHeight() / 2);
457
458 double sizeX2, sizeY2; // object widths^2 in x and y directions
459 double peakVal; // peak intensity in image
460
461 doMeasureCentroidImpl(&xc, &dxc, &yc, &dyc, &sizeX2, &sizeY2, &peakVal, mim, smoothingSigma, negative,
462 _flagHandler);
463
464 if (binsize > 1) {
465 // dilate from the lower left corner of central pixel
466 xc = (xc + 0.5) * binX - 0.5;
467 dxc *= binX;
468 sizeX2 *= binX * binX;
469
470 yc = (yc + 0.5) * binY - 0.5;
471 dyc *= binY;
472 sizeY2 *= binY * binY;
473 }
474
475 xc += x; // xc, yc are measured relative to pixel (x, y)
476 yc += y;
477
478 double const fac = _ctrl.wfac * (1 + smoothingSigma * smoothingSigma);
479 double const facX2 = fac * binX * binX;
480 double const facY2 = fac * binY * binY;
481
482 if (sizeX2 < facX2 && ::pow(xc - x, 2) < facX2 && sizeY2 < facY2 && ::pow(yc - y, 2) < facY2) {
483 if (binsize > 1 || _ctrl.peakMin < 0.0 || peakVal > _ctrl.peakMin) {
484 break;
485 }
486 }
487
488 if (sizeX2 >= facX2 || ::pow(xc - x, 2) >= facX2) {
489 binX *= 2;
490 }
491 if (sizeY2 >= facY2 || ::pow(yc - y, 2) >= facY2) {
492 binY *= 2;
493 }
494 }
495 result.x = afw::image::indexToPosition(xc + image.getX0());
496 result.y = afw::image::indexToPosition(yc + image.getY0());
497
498 result.xErr = sqrt(dxc * dxc);
499 result.yErr = sqrt(dyc * dyc);
500 measRecord.set(_centroidKey, result);
501 _centroidChecker(measRecord);
502}
503
505 _flagHandler.handleFailure(measRecord, error);
506}
507
513 if (flag == SdssCentroidAlgorithm::FAILURE) continue;
514 if (mapper.getInputSchema().getNames().count(mapper.getInputSchema().join(name, flag.name)) == 0)
515 continue;
517 mapper.getInputSchema().find<afw::table::Flag>(name + "_" + flag.name).key;
518 mapper.addMapping(key);
519 }
520}
521
522} // namespace base
523} // namespace meas
524} // namespace lsst
table::Key< std::string > name
AmpInfoBoxKey bbox
double x
#define LSST_EXCEPT(type,...)
table::Key< double > sigma2
afw::table::Key< double > sigma
afw::table::Key< afw::table::Array< ImagePixelT > > image
SchemaMapper * mapper
This implements the SdssCentroid algorithm within the meas_base measurement framework.
int y
MaskedImageT getMaskedImage()
std::shared_ptr< lsst::afw::detection::Psf const > getPsf() const
Field< T >::Value get(Key< T > const &key) const
void set(Key< T > const &key, U const &value)
SchemaItem< T > find(std::string const &name) const
A FunctorKey for CentroidResult.
Base for centroid measurement transformations.
Exception to be thrown when a measurement algorithm experiences a fatal error.
Definition: exceptions.h:76
vector-type utility class to build a collection of FlagDefinitions
Definition: FlagHandler.h:60
std::size_t size() const
return the current size (number of defined elements) of the collection
Definition: FlagHandler.h:125
Utility class for handling flag fields that indicate the failure modes of an algorithm.
Definition: FlagHandler.h:148
Exception to be thrown when a measurement algorithm experiences a known failure mode.
Definition: exceptions.h:48
static FlagDefinition const NOT_AT_MAXIMUM
Definition: SdssCentroid.h:77
virtual void fail(afw::table::SourceRecord &measRecord, MeasurementError *error=nullptr) const
Handle an exception thrown by the current algorithm by setting flags in the given record.
SdssCentroidAlgorithm(Control const &ctrl, std::string const &name, afw::table::Schema &schema)
static FlagDefinition const NO_SECOND_DERIVATIVE
Definition: SdssCentroid.h:75
static FlagDefinition const FAILURE
Definition: SdssCentroid.h:73
virtual void measure(afw::table::SourceRecord &measRecord, afw::image::Exposure< float > const &exposure) const
Called to measure a single child source in an image.
static FlagDefinitionList const & getFlagDefinitions()
Definition: SdssCentroid.cc:52
static FlagDefinition const ALMOST_NO_SECOND_DERIVATIVE
Definition: SdssCentroid.h:76
static FlagDefinition const EDGE
Definition: SdssCentroid.h:74
A C++ control class to handle SdssCentroidAlgorithm's configuration.
Definition: SdssCentroid.h:48
double wfac
"fiddle factor for adjusting the binning" ;
Definition: SdssCentroid.h:52
int binmax
"maximum allowed binning" ;
Definition: SdssCentroid.h:50
double peakMin
"if the peak's less than this insist on binning at least once" ;
Definition: SdssCentroid.h:51
SdssCentroidTransform(Control const &ctrl, std::string const &name, afw::table::SchemaMapper &mapper)
T exp(T... args)
T fabs(T... args)
T make_pair(T... args)
double indexToPosition(double ind)
void convolve(OutImageT &convolvedImage, InImageT const &inImage, KernelT const &kernel, ConvolutionControl const &convolutionControl=ConvolutionControl())
std::shared_ptr< ImageT > binImage(ImageT const &inImage, int const binX, int const binY, lsst::afw::math::Property const flags=lsst::afw::math::MEAN)
constexpr double PI
Extent< int, 2 > ExtentI
Point< double, 2 > Point2D
Point< int, 2 > Point2I
@ SIGMA_ONLY
Only the diagonal elements of the covariance matrix are provided.
Definition: constants.h:45
T pow(T... args)
T reset(T... args)
T sqrt(T... args)
A reusable struct for centroid measurements.
CentroidElement y
y (row) coordinate of the measured position
CentroidElement x
x (column) coordinate of the measured position
ErrElement yErr
standard deviation of y
ErrElement xErr
standard deviation of x
Simple class used to define and document flags The name and doc constitute the identity of the FlagDe...
Definition: FlagHandler.h:40
Key< int > psf
table::Schema schema