LSSTApplications
20.0.0
LSSTDataManagementBasePackage
|
#include <KernelSolution.h>
template<typename InputT>
class lsst::ip::diffim::StaticKernelSolution< InputT >
Definition at line 83 of file KernelSolution.h.
◆ ImageT
◆ PixelT
◆ Ptr
template<typename InputT >
◆ ConditionNumberType
◆ KernelSolvedBy
Enumerator |
---|
NONE | |
CHOLESKY_LDLT | |
CHOLESKY_LLT | |
LU | |
EIGENVECTOR | |
Definition at line 37 of file KernelSolution.h.
◆ StaticKernelSolution()
template<typename InputT >
◆ ~StaticKernelSolution()
template<typename InputT >
◆ _setKernel()
template<typename InputT >
Set kernel after solution.
Definition at line 422 of file KernelSolution.cc.
427 unsigned int const nParameters =
_aVec.size();
429 unsigned int const nKernelParameters =
430 std::dynamic_pointer_cast<afwMath::LinearCombinationKernel>(
_kernel)->getKernelList().size();
431 if (nParameters != (nKernelParameters + nBackgroundParameters))
436 for (
unsigned int idx = 0; idx < nKernelParameters; idx++) {
439 str(
boost::format(
"Unable to determine kernel solution %d (nan)") % idx));
441 kValues[idx] =
_aVec(idx);
453 str(
boost::format(
"Unable to determine background solution %d (nan)") %
◆ _setKernelUncertainty()
template<typename InputT >
◆ build()
template<typename InputT >
Definition at line 261 of file KernelSolution.cc.
270 "Error: variance less than 0.0");
274 "Error: variance equals 0.0, cannot inverse variance weight");
278 std::dynamic_pointer_cast<afwMath::LinearCombinationKernel>(
_kernel)->getKernelList();
280 unsigned int const nKernelParameters = basisList.
size();
282 unsigned int const nParameters = nKernelParameters + nBackgroundParameters;
319 unsigned int const startCol = goodBBox.
getMinX();
320 unsigned int const startRow = goodBBox.
getMinY();
322 unsigned int endCol = goodBBox.
getMaxX() + 1;
323 unsigned int endRow = goodBBox.
getMaxY() + 1;
338 startRow, startCol, endRow-startRow, endCol-startCol
339 ).array().inverse().matrix();
342 eigenTemplate.resize(eigenTemplate.rows()*eigenTemplate.cols(), 1);
343 eigenScience.resize(eigenScience.rows()*eigenScience.cols(), 1);
344 eigeniVariance.resize(eigeniVariance.rows()*eigeniVariance.cols(), 1);
355 for (kiter = basisList.
begin(); kiter != basisList.
end(); ++kiter, ++eiter) {
362 cMat.resize(cMat.size(), 1);
367 double time = t.elapsed();
368 LOGL_DEBUG(
"TRACE3.ip.diffim.StaticKernelSolution.build",
369 "Total compute time to do basis convolutions : %.2f s", time);
376 Eigen::MatrixXd cMat(eigenTemplate.col(0).size(), nParameters);
379 for (
unsigned int kidxj = 0; eiterj != eiterE; eiterj++, kidxj++) {
380 cMat.col(kidxj) = eiterj->col(0);
384 cMat.col(nParameters-1).fill(1.);
387 _ivVec = eigeniVariance.col(0);
388 _iVec = eigenScience.col(0);
◆ getB()
Eigen::VectorXd const& lsst::ip::diffim::KernelSolution::getB |
( |
| ) |
|
|
inlineinherited |
◆ getBackground()
template<typename InputT >
◆ getConditionNumber() [1/2]
double lsst::ip::diffim::KernelSolution::getConditionNumber |
( |
ConditionNumberType |
conditionType | ) |
|
|
virtualinherited |
◆ getConditionNumber() [2/2]
double lsst::ip::diffim::KernelSolution::getConditionNumber |
( |
Eigen::MatrixXd const & |
mMat, |
|
|
ConditionNumberType |
conditionType |
|
) |
| |
|
virtualinherited |
Definition at line 98 of file KernelSolution.cc.
100 switch (conditionType) {
103 Eigen::SelfAdjointEigenSolver<Eigen::MatrixXd> eVecValues(mMat);
104 Eigen::VectorXd eValues = eVecValues.eigenvalues();
105 double eMax = eValues.maxCoeff();
106 double eMin = eValues.minCoeff();
107 LOGL_DEBUG(
"TRACE3.ip.diffim.KernelSolution.getConditionNumber",
108 "EIGENVALUE eMax / eMin = %.3e", eMax / eMin);
109 return (eMax / eMin);
114 Eigen::VectorXd sValues = mMat.jacobiSvd().singularValues();
115 double sMax = sValues.maxCoeff();
116 double sMin = sValues.minCoeff();
117 LOGL_DEBUG(
"TRACE3.ip.diffim.KernelSolution.getConditionNumber",
118 "SVD eMax / eMin = %.3e", sMax / sMin);
119 return (sMax / sMin);
125 "Undefined ConditionNumberType : only EIGENVALUE, SVD allowed.");
◆ getId()
int lsst::ip::diffim::KernelSolution::getId |
( |
| ) |
const |
|
inlineinherited |
◆ getKernel()
template<typename InputT >
◆ getKsum()
template<typename InputT >
◆ getM()
Eigen::MatrixXd const& lsst::ip::diffim::KernelSolution::getM |
( |
| ) |
|
|
inlineinherited |
◆ getSolutionPair()
template<typename InputT >
◆ getSolvedBy()
◆ makeKernelImage()
template<typename InputT >
◆ printA()
void lsst::ip::diffim::KernelSolution::printA |
( |
| ) |
|
|
inlineinherited |
◆ printB()
void lsst::ip::diffim::KernelSolution::printB |
( |
| ) |
|
|
inlineinherited |
◆ printM()
void lsst::ip::diffim::KernelSolution::printM |
( |
| ) |
|
|
inlineinherited |
◆ solve() [1/2]
template<typename InputT >
◆ solve() [2/2]
void lsst::ip::diffim::KernelSolution::solve |
( |
Eigen::MatrixXd const & |
mMat, |
|
|
Eigen::VectorXd const & |
bVec |
|
) |
| |
|
virtualinherited |
Definition at line 131 of file KernelSolution.cc.
141 Eigen::VectorXd aVec = Eigen::VectorXd::Zero(bVec.size());
146 LOGL_DEBUG(
"TRACE2.ip.diffim.KernelSolution.solve",
147 "Solving for kernel");
149 Eigen::FullPivLU<Eigen::MatrixXd> lu(mMat);
150 if (lu.isInvertible()) {
151 aVec = lu.solve(bVec);
153 LOGL_DEBUG(
"TRACE3.ip.diffim.KernelSolution.solve",
154 "Unable to determine kernel via LU");
159 Eigen::SelfAdjointEigenSolver<Eigen::MatrixXd> eVecValues(mMat);
160 Eigen::MatrixXd
const& rMat = eVecValues.eigenvectors();
161 Eigen::VectorXd eValues = eVecValues.eigenvalues();
163 for (
int i = 0; i != eValues.rows(); ++i) {
164 if (eValues(i) != 0.0) {
165 eValues(i) = 1.0/eValues(i);
169 aVec = rMat * eValues.asDiagonal() * rMat.transpose() * bVec;
173 LOGL_DEBUG(
"TRACE3.ip.diffim.KernelSolution.solve",
174 "Unable to determine kernel via eigen-values");
180 double time = t.elapsed();
181 LOGL_DEBUG(
"TRACE3.ip.diffim.KernelSolution.solve",
182 "Compute time for matrix math : %.2f s", time);
◆ _aVec
Eigen::VectorXd lsst::ip::diffim::KernelSolution::_aVec |
|
protectedinherited |
◆ _background
template<typename InputT >
◆ _bVec
Eigen::VectorXd lsst::ip::diffim::KernelSolution::_bVec |
|
protectedinherited |
◆ _cMat
template<typename InputT >
◆ _fitForBackground
bool lsst::ip::diffim::KernelSolution::_fitForBackground |
|
protectedinherited |
◆ _id
int lsst::ip::diffim::KernelSolution::_id |
|
protectedinherited |
◆ _iVec
template<typename InputT >
◆ _ivVec
template<typename InputT >
◆ _kernel
template<typename InputT >
◆ _kSum
template<typename InputT >
◆ _mMat
Eigen::MatrixXd lsst::ip::diffim::KernelSolution::_mMat |
|
protectedinherited |
◆ _SolutionId
int lsst::ip::diffim::KernelSolution::_SolutionId = 0 |
|
staticprotectedinherited |
◆ _solvedBy
The documentation for this class was generated from the following files:
- /j/snowflake/release/lsstsw/stack/1a1d771/Linux64/ip_diffim/20.0.0/include/lsst/ip/diffim/KernelSolution.h
- /j/snowflake/release/lsstsw/stack/1a1d771/Linux64/ip_diffim/20.0.0/src/KernelSolution.cc
Backwards-compatibility support for depersisting the old Calib (FluxMag0/FluxMag0Err) objects.
lsst::geom::Extent2I const getDimensions() const
Return the Kernel's dimensions (width, height)
std::shared_ptr< lsst::afw::math::Kernel > _kernel
Derived single-object convolution kernel.
lsst::afw::image::Image< lsst::afw::math::Kernel::Pixel > ImageT
#define LSST_EXCEPT_ADD(e, m)
Add the current location and a message to an existing exception before rethrowing it.
def format(config, name=None, writeSourceLine=True, prefix="", verbose=False)
double computeImage(lsst::afw::image::Image< Pixel > &image, bool doNormalize, double x=0.0, double y=0.0) const
Compute an image (pixellized representation of the kernel) in place.
Eigen::VectorXd _bVec
Derived least squares B vector.
double getValue(Property const prop=NOTHING) const
Return the value of the desired property (if specified in the constructor)
Statistics makeStatistics(lsst::afw::image::Image< Pixel > const &img, lsst::afw::image::Mask< image::MaskPixel > const &msk, int const flags, StatisticsControl const &sctrl=StatisticsControl())
Handle a watered-down front-end to the constructor (no variance)
bool _fitForBackground
Background terms included in fit.
double _background
Derived differential background estimate.
KernelSolvedBy _solvedBy
Type of algorithm used to make solution.
Eigen::MatrixXd _mMat
Derived least squares M matrix.
A kernel that is a linear combination of fixed basis kernels.
Eigen::MatrixXd imageToEigenMatrix(lsst::afw::image::Image< PixelT > const &img)
Turns a 2-d Image into a 2-d Eigen Matrix.
virtual double getConditionNumber(ConditionNumberType conditionType)
int _id
Unique ID for object.
Eigen::VectorXd _ivVec
Inverse variance.
#define LSST_EXCEPT(type,...)
Create an exception with a given type.
Eigen::VectorXd _iVec
Vectorized I.
lsst::geom::Extent2I getDimensions() const
Return the image's size; useful for passing to constructors.
void _setKernel()
Set kernel after solution.
Reports invalid arguments.
int getMaxY() const noexcept
#define LOGL_DEBUG(logger, message...)
@ MIN
estimate sample minimum
int getMaxX() const noexcept
Eigen::MatrixXd _cMat
K_i x R.
An integer coordinate rectangle.
int getMinX() const noexcept
Provides consistent interface for LSST exceptions.
lsst::geom::Box2I getBBox(ImageOrigin origin=PARENT) const
A class to represent a 2-dimensional array of pixels.
void setKernelParameters(std::vector< double > const ¶ms)
Set the kernel parameters of a spatially invariant kernel.
int getMinY() const noexcept
Eigen::VectorXd _aVec
Derived least squares solution matrix.
double _kSum
Derived kernel sum.