lsst.jointcal
15.0-18-ga2cb10d+1
|
Class that handles the photometric least squares problem. More...
#include <PhotometryFit.h>
Public Member Functions | |
PhotometryFit (std::shared_ptr< Associations > associations, std::shared_ptr< PhotometryModel > photometryModel) | |
this is the only constructor More... | |
PhotometryFit (PhotometryFit const &)=delete | |
No copy or move: there is only ever one fitter of a given type. More... | |
PhotometryFit (PhotometryFit &&)=delete | |
PhotometryFit & | operator= (PhotometryFit const &)=delete |
PhotometryFit & | operator= (PhotometryFit &&)=delete |
void | assignIndices (std::string const &whatToFit) override |
Set parameters to fit and assign indices in the big matrix. More... | |
void | offsetParams (Eigen::VectorXd const &delta) override |
Offset the parameters by the requested quantities. More... | |
void | saveChi2MeasContributions (std::string const &baseName) const override |
Save a CSV file containing residuals of measurement terms. More... | |
void | saveChi2RefContributions (std::string const &baseName) const override |
Save a CSV file containing residuals of reference terms. More... | |
MinimizeResult | minimize (std::string const &whatToFit, double const nSigmaCut=0, bool const doRankUpdate=true, std::string const &dumpMatrixFile="") |
Does a 1 step minimization, assuming a linear model. More... | |
Chi2Statistic | computeChi2 () const |
Returns the chi2 for the current state. More... | |
void | leastSquareDerivatives (TripletList &tripletList, Eigen::VectorXd &grad) const |
Evaluates the chI^2 derivatives (Jacobian and gradient) for the current whatToFit setting. More... | |
virtual void | saveChi2Contributions (std::string const &baseName) const |
Save the full chi2 term per star that was used in the minimization, for debugging. More... | |
Protected Member Functions | |
unsigned | findOutliers (double nSigmaCut, MeasuredStarList &msOutliers, FittedStarList &fsOutliers) const |
Find Measurements and references contributing more than a cut, computed as \[ <chi2> + nSigmaCut + rms(chi2). \] The outliers are NOT removed, and no refit is done. More... | |
void | outliersContributions (MeasuredStarList &msOutliers, FittedStarList &fsOutliers, TripletList &tripletList, Eigen::VectorXd &grad) |
Contributions to derivatives from (presumably) outlier terms. More... | |
void | removeMeasOutliers (MeasuredStarList &outliers) |
Remove measuredStar outliers from the fit. No Refit done. More... | |
void | removeRefOutliers (FittedStarList &outliers) |
Remove refStar outliers from the fit. No Refit done. More... | |
Protected Attributes | |
std::shared_ptr< Associations > | _associations |
std::string | _whatToFit |
int | _lastNTrip |
unsigned int | _nParTot |
unsigned | _nMeasuredStars |
LOG_LOGGER | _log |
Class that handles the photometric least squares problem.
Definition at line 21 of file PhotometryFit.h.
|
inline |
this is the only constructor
Definition at line 24 of file PhotometryFit.h.
|
delete |
No copy or move: there is only ever one fitter of a given type.
|
delete |
|
overridevirtual |
Set parameters to fit and assign indices in the big matrix.
[in] | whatToFit | Valid strings : "Model", "Fluxes", which define which parameter sets are going to be fitted. whatToFit="Model Fluxes" will set both parameter sets variable when computing derivatives. Provided it contains "Model", whatToFit is passed over to the PhotometryModel, and can hence be used to control more finely which subsets of the photometric model are being fitted, if the the actual PhotometryModel implements such a possibility. |
Implements lsst::jointcal::FitterBase.
Definition at line 171 of file PhotometryFit.cc.
|
inherited |
Returns the chi2 for the current state.
Definition at line 16 of file FitterBase.cc.
|
protectedinherited |
Find Measurements and references contributing more than a cut, computed as
\[ <chi2> + nSigmaCut + rms(chi2). \]
The outliers are NOT removed, and no refit is done.
After returning from here, there are still measurements that contribute above the cut, but their contribution should be evaluated after a refit before discarding them.
[in] | nSigmaCut | Number of sigma to select on. |
[out] | msOutliers | list of MeasuredStar outliers to populate |
[out] | fsOutliers | list of FittedStar outliers to populate |
Definition at line 26 of file FitterBase.cc.
|
inherited |
Evaluates the chI^2 derivatives (Jacobian and gradient) for the current whatToFit setting.
The Jacobian is given as triplets in a sparse matrix, the gradient as a dense vector. The parameters which vary, and their indices, are to be set using assignIndices.
tripletList | tripletList of (row,col,value) representing the Jacobian of the chi2. |
grad | The gradient of the chi2. |
Definition at line 280 of file FitterBase.cc.
|
inherited |
Does a 1 step minimization, assuming a linear model.
This is a complete Newton Raphson step. Compute first and second derivatives, solve for the step and apply it, without a line search.
It calls assignIndices, leastSquareDerivatives, solves the linear system and calls offsetParams, then removes outliers in a loop if requested. Relies on sparse linear algebra.
[in] | whatToFit | See child method assignIndices for valid string values. |
[in] | nSigmaCut | How many sigma to reject outliers at. Outlier rejection ignored for nSigmaCut=0. |
[in] | doRankUpdate | Use CholmodSimplicialLDLT2.update() to do a fast rank update after outlier removal; otherwise do a slower full recomputation of the matrix. Only matters if nSigmaCut != 0. |
[in] | dumpMatrixFile | Write the pre-fit Hessian matrix and gradient to the files with "-mat.txt" and "-grad.txt". Be aware, this requires a large increase in memory usage to create a dense matrix before writing it; the output file may be large. Writing the matrix can be helpful for debugging bad fits. Read it and compute the real eigenvalues (recall that the Hessian is symmetric by construction) with numpy: hessian = np.matrix(np.loadtxt("dumpMatrixFile-mat.txt")) values, vectors = np.linalg.eigh(hessian) |
Definition at line 143 of file FitterBase.cc.
|
overridevirtual |
Offset the parameters by the requested quantities.
The used parameter layout is the one from the last call to assignIndices or minimize(). There is no easy way to check that the current setting of whatToFit and the provided Delta vector are compatible: we can only test the size.
[in] | delta | vector of offsets to apply |
Implements lsst::jointcal::FitterBase.
Definition at line 197 of file PhotometryFit.cc.
|
delete |
|
delete |
|
protectedinherited |
Contributions to derivatives from (presumably) outlier terms.
No discarding done.
Definition at line 255 of file FitterBase.cc.
|
protectedinherited |
Remove measuredStar outliers from the fit. No Refit done.
Definition at line 266 of file FitterBase.cc.
|
protectedinherited |
Remove refStar outliers from the fit. No Refit done.
Definition at line 274 of file FitterBase.cc.
|
virtualinherited |
Save the full chi2 term per star that was used in the minimization, for debugging.
Saves results to text files "baseName-meas.csv" and "baseName-ref.csv" for the MeasuredStar and RefStar contributions, respectively. This method is mostly useful for debugging: we will probably want to create a better persistence system for jointcal's internal representations in the future (see DM-12446).
Definition at line 288 of file FitterBase.cc.
|
overridevirtual |
Save a CSV file containing residuals of measurement terms.
Implements lsst::jointcal::FitterBase.
Definition at line 215 of file PhotometryFit.cc.
|
overridevirtual |
Save a CSV file containing residuals of reference terms.
Implements lsst::jointcal::FitterBase.
Definition at line 278 of file PhotometryFit.cc.
|
protectedinherited |
Definition at line 129 of file FitterBase.h.
|
protectedinherited |
Definition at line 132 of file FitterBase.h.
|
protectedinherited |
Definition at line 137 of file FitterBase.h.
|
protectedinherited |
Definition at line 134 of file FitterBase.h.
|
protectedinherited |
Definition at line 133 of file FitterBase.h.
|
protectedinherited |
Definition at line 130 of file FitterBase.h.