LSST Applications  21.0.0+04719a4bac,21.0.0-1-ga51b5d4+f5e6047307,21.0.0-11-g2b59f77+a9c1acf22d,21.0.0-11-ga42c5b2+86977b0b17,21.0.0-12-gf4ce030+76814010d2,21.0.0-13-g1721dae+760e7a6536,21.0.0-13-g3a573fe+768d78a30a,21.0.0-15-g5a7caf0+f21cbc5713,21.0.0-16-g0fb55c1+b60e2d390c,21.0.0-19-g4cded4ca+71a93a33c0,21.0.0-2-g103fe59+bb20972958,21.0.0-2-g45278ab+04719a4bac,21.0.0-2-g5242d73+3ad5d60fb1,21.0.0-2-g7f82c8f+8babb168e8,21.0.0-2-g8f08a60+06509c8b61,21.0.0-2-g8faa9b5+616205b9df,21.0.0-2-ga326454+8babb168e8,21.0.0-2-gde069b7+5e4aea9c2f,21.0.0-2-gecfae73+1d3a86e577,21.0.0-2-gfc62afb+3ad5d60fb1,21.0.0-25-g1d57be3cd+e73869a214,21.0.0-3-g357aad2+ed88757d29,21.0.0-3-g4a4ce7f+3ad5d60fb1,21.0.0-3-g4be5c26+3ad5d60fb1,21.0.0-3-g65f322c+e0b24896a3,21.0.0-3-g7d9da8d+616205b9df,21.0.0-3-ge02ed75+a9c1acf22d,21.0.0-4-g591bb35+a9c1acf22d,21.0.0-4-g65b4814+b60e2d390c,21.0.0-4-gccdca77+0de219a2bc,21.0.0-4-ge8a399c+6c55c39e83,21.0.0-5-gd00fb1e+05fce91b99,21.0.0-6-gc675373+3ad5d60fb1,21.0.0-64-g1122c245+4fb2b8f86e,21.0.0-7-g04766d7+cd19d05db2,21.0.0-7-gdf92d54+04719a4bac,21.0.0-8-g5674e7b+d1bd76f71f,master-gac4afde19b+a9c1acf22d,w.2021.13
LSST Data Management Base Package
CModel.h
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 
24 #ifndef LSST_MEAS_MODELFIT_CModelFit_h_INCLUDED
25 #define LSST_MEAS_MODELFIT_CModelFit_h_INCLUDED
26 
27 #include <bitset>
28 #include <string>
29 
30 #include "ndarray.h"
31 
32 #include "lsst/geom.h"
33 #include "lsst/pex/config.h"
35 #include "lsst/afw/table/Source.h"
45 
46 namespace lsst { namespace meas { namespace modelfit {
47 
143 
145  profileName("lux"),
146  priorSource("EMPIRICAL"),
147  priorName(),
148  nComponents(8),
149  maxRadius(0),
150  usePixelWeights(false),
151  weightsMultiplier(1.0),
152  doRecordHistory(true),
153  doRecordTime(true)
154  {}
155 
158  }
159 
160  PTR(Model) getModel() const;
161 
162  PTR(Prior) getPrior() const;
163 
166  "Name of the shapelet.RadialProfile that defines the model to fit"
167  );
168 
171  "One of 'FILE', 'LINEAR', 'EMPIRICAL', or 'NONE', indicating whether the prior should be loaded "
172  "from disk, created from one of the nested prior config/control objects, or None"
173  );
174 
177  "Name of the Prior that defines the model to fit (a filename in $MEAS_MODELFIT_DIR/data, "
178  "with no extension), if priorSource='FILE'. Ignored for forced fitting."
179  );
180 
183  "Configuration for a linear prior, used if priorSource='LINEAR'."
184  );
185 
188  "Configuration for an empirical prior, used if priorSource='EMPIRICAL'."
189  );
190 
191  LSST_CONTROL_FIELD(nComponents, int, "Number of Gaussian used to approximate the profile");
192 
194  maxRadius,
195  int,
196  "Maximum radius used in approximating profile with Gaussians (0=default for this profile)"
197  );
198 
201  bool,
202  "Use per-pixel variances as weights in the nonlinear fit (the final linear fit for"
203  " flux never uses per-pixel variances)"
204  );
205 
208  double,
209  "Scale the likelihood by this factor to artificially reweight it w.r.t. the prior."
210  );
211 
214  "Configuration for how the objective surface is explored. Ignored for forced fitting"
215  );
216 
218  doRecordHistory, bool,
219  "Whether to record the steps the optimizer takes (or just the number, if running as a plugin)"
220  );
221 
223  doRecordTime, bool,
224  "Whether to record the time spent in this stage"
225  );
226 
227 };
228 
234 
236  psfName("modelfit_DoubleShapeletPsfApprox"),
237  minInitialRadius(0.1),
239  {
240  initial.nComponents = 3; // use very rough model in initial fit
241  initial.optimizer.gradientThreshold = 1E-3; // with slightly coarser convergence criteria
243  initial.usePixelWeights = true;
244  dev.profileName = "luv";
245  exp.nComponents = 6;
247  }
248 
250  psfName,
251  std::string,
252  "Field name prefix of the Shapelet PSF approximation used to convolve the galaxy model; "
253  "must contain a set of fields matching the schema defined by shapelet.MultiShapeletFunctionKey."
254  );
255 
258  "Configuration parameters related to the determination of the pixels to include in the fit."
259  );
260 
263  "An initial fit (usually with a fast, approximate model) used to warm-start the exp and dev fits, "
264  "convolved with only the zeroth-order terms in the multi-shapelet PSF approximation."
265  );
266 
269  "Independent fit of the exponential component"
270  );
271 
274  "Independent fit of the de Vaucouleur component"
275  );
276 
278  minInitialRadius, double,
279  "Minimum initial radius in pixels (used to regularize initial moments-based PSF deconvolution)"
280  );
281 
284  "If the 2nd-moments shape used to initialize the fit failed, use the PSF moments multiplied by this."
285  " If <= 0.0, abort the fit early instead."
286  );
287 
288 };
289 
294 
296  enum FlagBit {
297  FAILED=0,
298  TR_SMALL,
302  NUMERIC_ERROR,
306  N_FLAGS
307  };
308 
310 
321 
322  ndarray::Array<Scalar const,1,1> nonlinear;
323  ndarray::Array<Scalar const,1,1> amplitudes;
324  ndarray::Array<Scalar const,1,1> fixed;
325 
328 };
329 
334 struct CModelResult {
335 
337  enum FlagBit {
338  FAILED=0,
352  NO_SHAPE,
354  SMALL_SHAPE,
360  N_FLAGS
361  };
362 
364 
371 
375 
378 
381 };
382 
392 public:
393 
396 
410  std::string const & name,
411  Control const & ctrl,
413  );
414 
431  std::string const & name,
432  Control const & ctrl,
433  afw::table::SchemaMapper & schemaMapper
434  );
435 
443  explicit CModelAlgorithm(Control const & ctrl);
444 
446  Control const & getControl() const { return _ctrl; }
447 
465  afw::image::Exposure<Pixel> const & exposure,
467  geom::Point2D const & center,
469  Scalar approxFlux=-1,
470  Scalar kronRadius=-1,
471  int footprintArea=-1
472  ) const;
473 
487  afw::image::Exposure<Pixel> const & exposure,
489  geom::Point2D const & center,
490  Result const & reference,
491  Scalar approxFlux=-1
492  ) const;
493 
506  void measure(
507  afw::table::SourceRecord & measRecord,
508  afw::image::Exposure<Pixel> const & exposure
509  ) const;
510 
526  void measure(
527  afw::table::SourceRecord & measRecord,
528  afw::image::Exposure<Pixel> const & exposure,
529  afw::table::SourceRecord const & refRecord
530  ) const;
531 
540  void fail(
541  afw::table::SourceRecord & measRecord,
543  ) const;
544 
547 
548 private:
549 
551 
552  // Actual implementations go here; we use an output argument for the result so we can get partial
553  // results to the plugin version when we throw.
554  void _applyImpl(
555  Result & result,
556  afw::image::Exposure<Pixel> const & exposure,
558  geom::Point2D const & center,
560  Scalar approxFlux,
561  Scalar kronRadius=-1,
562  int footprintArea=-1
563  ) const;
564 
565  // Actual implementations go here; we use an output argument for the result so we can get partial
566  // results to the SourceRecord version when we throw.
567  void _applyForcedImpl(
568  Result & result,
569  afw::image::Exposure<Pixel> const & exposure,
571  geom::Point2D const & center,
572  Result const & reference,
573  Scalar approxFlux
574  ) const;
575 
576  // gets/checks inputs from SourceRecord that are needed by both apply and applyForced
577  template <typename PixelT>
578  shapelet::MultiShapeletFunction _processInputs(
580  afw::image::Exposure<PixelT> const & exposure
581  ) const;
582 
583  class Impl;
584 
585  Control _ctrl;
586  PTR(Impl) _impl;
587 };
588 
589 }}} // namespace lsst::meas::modelfit
590 
591 #endif // !LSST_MEAS_MODELFIT_CModelFit_h_INCLUDED
py::object result
Definition: _schema.cc:430
table::Key< std::string > name
Definition: Amplifier.cc:116
table::Schema schema
Definition: python.h:134
#define PTR(...)
Definition: base.h:41
An ellipse core with quadrupole moments as parameters.
Definition: Quadrupole.h:47
Base class for all records.
Definition: BaseRecord.h:31
Defines the fields and offsets for a table.
Definition: Schema.h:50
A mapping between the keys of two Schemas, used to copy data between them.
Definition: SchemaMapper.h:21
Record class that contains measurements made on a single exposure.
Definition: Source.h:80
Exception to be thrown when a measurement algorithm experiences a known failure mode.
Definition: exceptions.h:48
Main public interface class for CModel algorithm.
Definition: CModel.h:391
void writeResultToRecord(Result const &result, afw::table::BaseRecord &record) const
Copy values from a Result struct to a BaseRecord object.
Result applyForced(afw::image::Exposure< Pixel > const &exposure, shapelet::MultiShapeletFunction const &psf, geom::Point2D const &center, Result const &reference, Scalar approxFlux=-1) const
Run the CModel algorithm in forced mode on an image, supplying inputs directly and returning outputs ...
void fail(afw::table::SourceRecord &measRecord, meas::base::MeasurementError *error) const
Handle an exception thrown by one of the measure() methods, setting the appropriate flag in the given...
CModelResult Result
Typedef to the master Result struct.
Definition: CModel.h:395
CModelAlgorithm(std::string const &name, Control const &ctrl, afw::table::SchemaMapper &schemaMapper)
Construct an algorithm instance suitable for forced photometry and add its fields to the Schema.
Control const & getControl() const
Return the control object the algorithm was constructed with.
Definition: CModel.h:446
CModelAlgorithm(std::string const &name, Control const &ctrl, afw::table::Schema &schema)
Construct an algorithm instance and add its fields to the Schema.
void measure(afw::table::SourceRecord &measRecord, afw::image::Exposure< Pixel > const &exposure, afw::table::SourceRecord const &refRecord) const
Run the CModel algorithm in forced mode on an image, using a SourceRecord for inputs and outputs.
void measure(afw::table::SourceRecord &measRecord, afw::image::Exposure< Pixel > const &exposure) const
Run the CModel algorithm on an image, using a SourceRecord for inputs and outputs.
CModelControl Control
Typedef to the master Control struct.
Definition: CModel.h:394
CModelAlgorithm(Control const &ctrl)
Construct an algorithm instance that cannot use SourceRecords for input/output.
Result apply(afw::image::Exposure< Pixel > const &exposure, shapelet::MultiShapeletFunction const &psf, geom::Point2D const &center, afw::geom::ellipses::Quadrupole const &moments, Scalar approxFlux=-1, Scalar kronRadius=-1, int footprintArea=-1) const
Run the CModel algorithm on an image, supplying inputs directly and returning outputs in a Result.
Abstract base class and concrete factories that define multi-shapelet galaxy models.
Definition: Model.h:56
Configuration object for Optimizer.
Definition: optimizer.h:187
double gradientThreshold
"If the maximum of the gradient falls below this threshold, consider the algorithm converged" ;
Definition: optimizer.h:207
double minTrustRadiusThreshold
"If the trust radius falls below this threshold, consider the algorithm converged" ;
Definition: optimizer.h:202
int maxOuterIterations
"maximum number of steps" ;
Definition: optimizer.h:272
Base class for objective functions for Optimizer.
Definition: optimizer.h:44
Base class for Bayesian priors.
Definition: Prior.h:36
A concrete Likelihood class that does not require its parameters and data to be in the same UnitSyste...
A multi-scale shapelet function.
Registry and utility class for multi-Gaussian approximations to radial profiles.
Definition: RadialProfile.h:60
static RadialProfile & get(std::string const &name)
Return a predefined radial profile given its name.
#define LSST_NESTED_CONTROL_FIELD(NAME, MODULE, TYPE, DOC)
A preprocessor macro used to define fields in C++ "control object" structs, for nested control object...
Definition: config.h:69
#define LSST_CONTROL_FIELD(NAME, TYPE, DOC)
A preprocessor macro used to define fields in C++ "control object" structs.
Definition: config.h:43
const char * source()
Source function that allows astChannel to source from a Stream.
Definition: Stream.h:224
double Scalar
Typedefs to be used for probability and parameter values.
Definition: common.h:44
A base class for image defects.
VectorQ moments
Definition: simpleShape.cc:151
The main control object for CModel, containing parameters for the final linear fit and aggregating th...
Definition: CModel.h:233
CModelStageControl exp
"Independent fit of the exponential component" ;
Definition: CModel.h:270
PixelFitRegionControl region
"Configuration parameters related to the determination of the pixels to include in the fit....
Definition: CModel.h:259
double fallbackInitialMomentsPsfFactor
"If the 2nd-moments shape used to initialize the fit failed, use the PSF moments multiplied by this....
Definition: CModel.h:286
double minInitialRadius
"Minimum initial radius in pixels (used to regularize initial moments-based PSF deconvolution)" ;
Definition: CModel.h:280
std::string psfName
"Field name prefix of the Shapelet PSF approximation used to convolve the galaxy model; " "must conta...
Definition: CModel.h:254
CModelStageControl initial
"An initial fit (usually with a fast, approximate model) used to warm-start the exp and dev fits,...
Definition: CModel.h:265
CModelStageControl dev
"Independent fit of the de Vaucouleur component" ;
Definition: CModel.h:275
Master result object for CModel, containing results for the final linear fit and three nested CModelS...
Definition: CModel.h:334
Scalar instFlux
Flux from the final linear fit.
Definition: CModel.h:365
Scalar objective
Objective value at the best-fit point (chisq/2)
Definition: CModel.h:370
Scalar fracDev
Fraction of flux from the final linear fit in the de Vaucouleur component (always between 0 and 1).
Definition: CModel.h:368
FlagBit
Flags that apply to all four CModel fits or just the last one.
Definition: CModel.h:337
@ REGION_USED_INITIAL_ELLIPSE_MAX
Fit region implied by the best-fit ellipse of the initial was too large, so we used the configuration...
Definition: CModel.h:350
@ BAD_REFERENCE
Reference fit failed, so forced fit will fail as well.
Definition: CModel.h:358
@ REGION_USED_PSF_AREA
Kron radius was unavailable or outside bounds, so the second-moment ellipse scaled to the PSF area wa...
Definition: CModel.h:346
@ NO_SHAPELET_PSF
Set if the Psf shapelet approximation failed.
Definition: CModel.h:356
@ FAILED
General failure flag for the linear fit flux; set if any other CModel flag is set,...
Definition: CModel.h:338
@ SMALL_SHAPE
Initial moments were sufficiently small that we used minInitialRadius to set the initial parameters.
Definition: CModel.h:354
@ REGION_USED_INITIAL_ELLIPSE_MIN
Fit region implied by the best-fit ellipse of the initial was too small, so we used the configuration...
Definition: CModel.h:348
@ N_FLAGS
Non-flag counter to indicate the number of flags.
Definition: CModel.h:360
@ REGION_USED_FOOTPRINT_AREA
Kron radius was unavailable or outside bounds, so the second-moment ellipse scaled to the footprint a...
Definition: CModel.h:343
@ NO_SHAPE
Set if the input SourceRecord had no valid shape slot with which to start the fit.
Definition: CModel.h:352
@ REGION_MAX_BAD_PIXEL_FRACTION
Set if we aborted early because the fit region had too many bad pixels.
Definition: CModel.h:341
@ BAD_CENTROID
Input centroid did not land within the fit region.
Definition: CModel.h:357
@ NO_FLUX
No flux was measured.
Definition: CModel.h:359
@ REGION_MAX_AREA
Set if we aborted early because the fit region was too large.
Definition: CModel.h:340
Scalar instFluxErr
Flux uncertainty from the final linear fit.
Definition: CModel.h:366
afw::geom::ellipses::Quadrupole finalFitRegion
Pixels used in the exp, dev, and linear fits.
Definition: CModel.h:377
CModelStageResult dev
Results from the de Vaucouleur (Sersic n=4) fit.
Definition: CModel.h:374
CModelStageResult initial
Results from the initial approximate nonlinear fit that feeds the others.
Definition: CModel.h:372
CModelStageResult exp
Results from the exponential (Sersic n=1) fit.
Definition: CModel.h:373
Scalar instFluxInner
Flux measured strictly within the fit region (no extrapolation).
Definition: CModel.h:367
std::bitset< N_FLAGS > flags
Array of flags.
Definition: CModel.h:380
afw::geom::ellipses::Quadrupole initialFitRegion
Pixels used in the initial fit.
Definition: CModel.h:376
LocalUnitTransform fitSysToMeasSys
Transforms to the coordinate system where parameters are defined.
Definition: CModel.h:379
Nested control object for CModel that configures one of the three ("initial", "exp",...
Definition: CModel.h:142
std::string profileName
"Name of the shapelet.RadialProfile that defines the model to fit" ;
Definition: CModel.h:167
int nComponents
"Number of Gaussian used to approximate the profile" ;
Definition: CModel.h:191
shapelet::RadialProfile const & getProfile() const
Definition: CModel.h:156
std::string priorName
"Name of the Prior that defines the model to fit (a filename in $MEAS_MODELFIT_DIR/data,...
Definition: CModel.h:179
OptimizerControl optimizer
"Configuration for how the objective surface is explored. Ignored for forced fitting" ;
Definition: CModel.h:215
bool doRecordTime
"Whether to record the time spent in this stage" ;
Definition: CModel.h:225
boost::shared_ptr< Prior > getPrior() const
bool usePixelWeights
"Use per-pixel variances as weights in the nonlinear fit (the final linear fit for" " flux never uses...
Definition: CModel.h:204
SemiEmpiricalPriorControl empiricalPriorConfig
"Configuration for an empirical prior, used if priorSource='EMPIRICAL'." ;
Definition: CModel.h:189
bool doRecordHistory
"Whether to record the steps the optimizer takes (or just the number, if running as a plugin)" ;
Definition: CModel.h:220
int maxRadius
"Maximum radius used in approximating profile with Gaussians (0=default for this profile)" ;
Definition: CModel.h:197
boost::shared_ptr< Model > getModel() const
double weightsMultiplier
"Scale the likelihood by this factor to artificially reweight it w.r.t. the prior....
Definition: CModel.h:210
std::string priorSource
"One of 'FILE', 'LINEAR', 'EMPIRICAL', or 'NONE', indicating whether the prior should be loaded " "fr...
Definition: CModel.h:173
SoftenedLinearPriorControl linearPriorConfig
"Configuration for a linear prior, used if priorSource='LINEAR'." ;
Definition: CModel.h:184
Result object for a single nonlinear fitting stage of the CModel algorithm.
Definition: CModel.h:293
ndarray::Array< Scalar const, 1, 1 > fixed
Opaque fixed parameters in specialized units.
Definition: CModel.h:324
Scalar instFluxErr
Flux uncertainty from just this stage fit.
Definition: CModel.h:316
Scalar time
Time spent in this fit in seconds.
Definition: CModel.h:319
boost::shared_ptr< UnitTransformedLikelihood > likelihood
Object used to evaluate models and compare to data.
Definition: CModel.h:314
boost::shared_ptr< Prior > prior
Bayesian priors on the parameters (defined fully by Control struct)
Definition: CModel.h:312
afw::table::BaseCatalog history
Trace of the optimizer's path, if enabled by diagnostic options.
Definition: CModel.h:326
boost::shared_ptr< Model > model
Model object that defines the parametrization (defined fully by Control struct)
Definition: CModel.h:311
FlagBit
Flags for a single CModel stage (note that there are additional flags for the full multi-stage fit)
Definition: CModel.h:296
@ TR_SMALL
Whether convergence was due to the optimizer trust region getting too small (not a failure!...
Definition: CModel.h:298
@ N_FLAGS
Non-flag counter to indicate the number of flags.
Definition: CModel.h:306
@ MAX_ITERATIONS
Whether the optimizer exceeded the maximum number of iterations.
Definition: CModel.h:300
@ NUMERIC_ERROR
Optimizer encountered a numerical error (something likely went to infinity).
Definition: CModel.h:302
@ BAD_REFERENCE
Reference fit failed, so forced fit will fail as well.
Definition: CModel.h:304
@ FAILED
General flag, indicating whether the flux for this stage can be trusted.
Definition: CModel.h:297
@ NO_FLUX
No flux was measured.
Definition: CModel.h:305
Scalar objective
Value of the objective function at the best fit point: chisq/2 - ln(prior)
Definition: CModel.h:318
ndarray::Array< Scalar const, 1, 1 > nonlinear
Opaque nonlinear parameters in specialized units.
Definition: CModel.h:322
std::bitset< N_FLAGS > flags
Array of flags.
Definition: CModel.h:327
Scalar instFluxInner
Flux measured strictly within the fit region (no extrapolation).
Definition: CModel.h:317
Scalar instFlux
Flux measured from just this stage fit.
Definition: CModel.h:315
boost::shared_ptr< OptimizerObjective > objfunc
Objective class used by the optimizer.
Definition: CModel.h:313
afw::geom::ellipses::Quadrupole ellipse
Best fit half-light ellipse in pixel coordinates.
Definition: CModel.h:320
ndarray::Array< Scalar const, 1, 1 > amplitudes
Opaque linear parameters in specialized units.
Definition: CModel.h:323
A local mapping between two UnitSystems.
Definition: UnitSystem.h:80
Key< int > psf
Definition: Exposure.cc:65