LSST Applications  21.0.0-172-gfb10e10a+18fedfabac,22.0.0+297cba6710,22.0.0+80564b0ff1,22.0.0+8d77f4f51a,22.0.0+a28f4c53b1,22.0.0+dcf3732eb2,22.0.1-1-g7d6de66+2a20fdde0d,22.0.1-1-g8e32f31+297cba6710,22.0.1-1-geca5380+7fa3b7d9b6,22.0.1-12-g44dc1dc+2a20fdde0d,22.0.1-15-g6a90155+515f58c32b,22.0.1-16-g9282f48+790f5f2caa,22.0.1-2-g92698f7+dcf3732eb2,22.0.1-2-ga9b0f51+7fa3b7d9b6,22.0.1-2-gd1925c9+bf4f0e694f,22.0.1-24-g1ad7a390+a9625a72a8,22.0.1-25-g5bf6245+3ad8ecd50b,22.0.1-25-gb120d7b+8b5510f75f,22.0.1-27-g97737f7+2a20fdde0d,22.0.1-32-gf62ce7b1+aa4237961e,22.0.1-4-g0b3f228+2a20fdde0d,22.0.1-4-g243d05b+871c1b8305,22.0.1-4-g3a563be+32dcf1063f,22.0.1-4-g44f2e3d+9e4ab0f4fa,22.0.1-42-gca6935d93+ba5e5ca3eb,22.0.1-5-g15c806e+85460ae5f3,22.0.1-5-g58711c4+611d128589,22.0.1-5-g75bb458+99c117b92f,22.0.1-6-g1c63a23+7fa3b7d9b6,22.0.1-6-g50866e6+84ff5a128b,22.0.1-6-g8d3140d+720564cf76,22.0.1-6-gd805d02+cc5644f571,22.0.1-8-ge5750ce+85460ae5f3,master-g6e05de7fdc+babf819c66,master-g99da0e417a+8d77f4f51a,w.2021.48
LSST Data Management Base Package
KernelPca.cc
Go to the documentation of this file.
1 // -*- lsst-c++ -*-
12 #include "lsst/afw/math.h"
13 #include "lsst/afw/image.h"
14 #include "lsst/log/Log.h"
16 
19 
20 namespace afwMath = lsst::afw::math;
21 namespace afwImage = lsst::afw::image;
23 
24 namespace lsst {
25 namespace ip {
26 namespace diffim {
27 namespace detail {
28 
55  template<typename PixelT>
58  ) :
59  afwMath::CandidateVisitor(),
60  _imagePca(imagePca),
61  _mean()
62  {};
63 
64  template<typename PixelT>
66  afwMath::KernelList kernelList;
67 
68  std::vector<std::shared_ptr<ImageT>> eigenImages = _imagePca->getEigenImages();
69  int ncomp = eigenImages.size();
70 
71  if (_mean) {
75  (*_mean, true))));
76  }
77  for (int i = 0; i < ncomp; i++) {
79  afwImage::Image<afwMath::Kernel::Pixel>(*eigenImages[i], true);
81  new afwMath::FixedKernel(img)
82  ));
83  }
84 
85  return kernelList;
86  }
87 
88  template<typename PixelT>
90 
91  KernelCandidate<PixelT> *kCandidate = dynamic_cast<KernelCandidate<PixelT> *>(candidate);
92  if (kCandidate == NULL) {
94  "Failed to cast SpatialCellCandidate to KernelCandidate");
95  }
96  LOGL_DEBUG("TRACE5.ip.diffim.SetPcaImageVisitor.processCandidate",
97  "Processing candidate %d", kCandidate->getId());
98 
99  try {
100  /* Normalize to unit sum */
101  std::shared_ptr<ImageT> kImage = kCandidate->getKernelSolution(
102  KernelCandidate<PixelT>::ORIG)->makeKernelImage();
103  *kImage /= kCandidate->getKernelSolution(
104  KernelCandidate<PixelT>::ORIG)->getKsum();
105  /* Tell imagePca they have the same weighting in the Pca */
106  _imagePca->addImage(kImage, 1.0);
107  } catch(pexExcept::Exception &e) {
108  return;
109  }
110  }
111 
112  template<typename PixelT>
114  /*
115  If we don't subtract off the mean before we do the Pca, the
116  subsequent terms carry less of the power than if you do subtract
117  off the mean. Explicit example:
118 
119  With mean subtraction:
120  DEBUG: Eigenvalue 0 : 0.010953 (0.373870 %)
121  DEBUG: Eigenvalue 1 : 0.007927 (0.270604 %)
122  DEBUG: Eigenvalue 2 : 0.001393 (0.047542 %)
123  DEBUG: Eigenvalue 3 : 0.001092 (0.037261 %)
124  DEBUG: Eigenvalue 4 : 0.000829 (0.028283 %)
125 
126  Without mean subtraction:
127  DEBUG: Eigenvalue 0 : 0.168627 (0.876046 %)
128  DEBUG: Eigenvalue 1 : 0.007935 (0.041223 %)
129  DEBUG: Eigenvalue 2 : 0.006049 (0.031424 %)
130  DEBUG: Eigenvalue 3 : 0.001188 (0.006173 %)
131  DEBUG: Eigenvalue 4 : 0.001050 (0.005452 %)
132 
133  After the first term above, which basically represents the mean,
134  the remaining terms carry less of the power than if you do
135  subtract off the mean. (0.041223/(1-0.876046) < 0.373870).
136  */
137  LOGL_DEBUG("TRACE5.ip.diffim.KernelPcaVisitor.subtractMean",
138  "Subtracting mean feature before Pca");
139 
140  _mean = _imagePca->getMean();
141  KernelPca<ImageT>::ImageList imageList = _imagePca->getImageList();
142  for (typename KernelPca<ImageT>::ImageList::const_iterator ptr = imageList.begin(),
143  end = imageList.end(); ptr != end; ++ptr) {
144  **ptr -= *_mean;
145  }
146  }
147 
162  template <typename ImageT>
164  {
165  Super::analyze();
166 
167  typename Super::ImageList const &eImageList = this->getEigenImages();
168  typename Super::ImageList::const_iterator iter = eImageList.begin(), end = eImageList.end();
169  for (size_t i = 0; iter != end; ++i, ++iter) {
170  std::shared_ptr<ImageT> eImage = *iter;
171 
172  /*
173  * Normalise eigenImages to have a maximum of 1.0. For n > 0 they
174  * (should) have mean == 0, so we can't use that to normalize
175  */
177  double const min = stats.getValue(afwMath::MIN);
178  double const max = stats.getValue(afwMath::MAX);
179 
180  double const extreme = (fabs(min) > max) ? min :max;
181  if (extreme != 0.0) {
182  *eImage /= extreme;
183  }
184  }
185  }
186 
187 
188  typedef float PixelT;
189  template class KernelPcaVisitor<PixelT>;
191 
192 }}}} // end of namespace lsst::ip::diffim::detail
int min
int end
int max
#define LSST_EXCEPT(type,...)
Create an exception with a given type.
Definition: Exception.h:48
afw::table::Key< afw::table::Array< ImagePixelT > > image
Class used by SpatialModelCell for spatial Kernel fitting.
Declaration of KernelPca and KernelPcaVisitor.
LSST DM logging module built on log4cxx.
#define LOGL_DEBUG(logger, message...)
Log a debug-level message using a varargs/printf style interface.
Definition: Log.h:515
uint64_t * ptr
Definition: RangeSet.cc:88
T begin(T... args)
A class to represent a 2-dimensional array of pixels.
Definition: Image.h:51
A kernel created from an Image.
Definition: Kernel.h:471
Base class for candidate objects in a SpatialCell.
Definition: SpatialCell.h:70
int getId() const
Return the candidate's unique ID.
Definition: SpatialCell.h:102
A class to evaluate image statistics.
Definition: Statistics.h:220
double getValue(Property const prop=NOTHING) const
Return the value of the desired property (if specified in the constructor)
Definition: Statistics.cc:1047
Class stored in SpatialCells for spatial Kernel fitting.
std::shared_ptr< StaticKernelSolution< PixelT > > getKernelSolution(CandidateSwitch cand) const
Overrides the analyze method of base class afwImage::ImagePca.
Definition: KernelPca.h:24
virtual void analyze()
Generate eigenimages that are normalised.
Definition: KernelPca.cc:163
A class to run a PCA on all candidate kernels (represented as Images).
Definition: KernelPca.h:40
void processCandidate(lsst::afw::math::SpatialCellCandidate *candidate)
Definition: KernelPca.cc:89
lsst::afw::math::KernelList getEigenKernels()
Definition: KernelPca.cc:65
KernelPcaVisitor(std::shared_ptr< KernelPca< ImageT > > imagePca)
Definition: KernelPca.cc:56
Provides consistent interface for LSST exceptions.
Definition: Exception.h:107
Reports errors in the logical structure of the program.
Definition: Runtime.h:46
T end(T... args)
Backwards-compatibility support for depersisting the old Calib (FluxMag0/FluxMag0Err) objects.
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)
Definition: Statistics.h:359
@ MIN
estimate sample minimum
Definition: Statistics.h:75
@ MAX
estimate sample maximum
Definition: Statistics.h:76
A base class for image defects.
T push_back(T... args)
T size(T... args)