LSST Applications g0b6bd0c080+a72a5dd7e6,g1182afd7b4+2a019aa3bb,g17e5ecfddb+2b8207f7de,g1d67935e3f+06cf436103,g38293774b4+ac198e9f13,g396055baef+6a2097e274,g3b44f30a73+6611e0205b,g480783c3b1+98f8679e14,g48ccf36440+89c08d0516,g4b93dc025c+98f8679e14,g5c4744a4d9+a302e8c7f0,g613e996a0d+e1c447f2e0,g6c8d09e9e7+25247a063c,g7271f0639c+98f8679e14,g7a9cd813b8+124095ede6,g9d27549199+a302e8c7f0,ga1cf026fa3+ac198e9f13,ga32aa97882+7403ac30ac,ga786bb30fb+7a139211af,gaa63f70f4e+9994eb9896,gabf319e997+ade567573c,gba47b54d5d+94dc90c3ea,gbec6a3398f+06cf436103,gc6308e37c7+07dd123edb,gc655b1545f+ade567573c,gcc9029db3c+ab229f5caf,gd01420fc67+06cf436103,gd877ba84e5+06cf436103,gdb4cecd868+6f279b5b48,ge2d134c3d5+cc4dbb2e3f,ge448b5faa6+86d1ceac1d,gecc7e12556+98f8679e14,gf3ee170dca+25247a063c,gf4ac96e456+ade567573c,gf9f5ea5b4d+ac198e9f13,gff490e6085+8c2580be5c,w.2022.27
LSST Data Management Base Package
mixture.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
24#include "pybind11/pybind11.h"
25#include "pybind11/eigen.h"
26#include "pybind11/stl.h"
27
28#include <sstream> // Python.h must come before even system headers
29
30#include "ndarray/pybind11.h"
31#include "ndarray/eigen.h"
32
33#include "lsst/utils/python.h"
35
36namespace py = pybind11;
37using namespace pybind11::literals;
38
39namespace lsst {
40namespace meas {
41namespace modelfit {
42namespace {
43
44using PyMixtureComponent = py::class_<MixtureComponent>;
45using PyMixtureUpdateRestriction =
46 py::class_<MixtureUpdateRestriction, std::shared_ptr<MixtureUpdateRestriction>>;
47using PyMixture = py::class_<Mixture, std::shared_ptr<Mixture>>;
48
49static PyMixtureComponent declareMixtureComponent(py::module &mod) {
50 PyMixtureComponent cls(mod, "MixtureComponent");
51 cls.def("getDimension", &MixtureComponent::getDimension);
52 cls.def_readwrite("weight", &MixtureComponent::weight);
53 cls.def("getMu", &MixtureComponent::getMu);
54 cls.def("setMu", &MixtureComponent::setMu);
55 cls.def("getSigma", &MixtureComponent::getSigma);
56 cls.def("setSigma", &MixtureComponent::setSigma);
57 cls.def("project", (MixtureComponent (MixtureComponent::*)(int) const) & MixtureComponent::project,
58 "dim"_a);
59 cls.def("project", (MixtureComponent (MixtureComponent::*)(int, int) const) & MixtureComponent::project,
60 "dim1"_a, "dim2"_a);
61 cls.def(py::init<int>(), "dim"_a);
62 cls.def(py::init<Scalar, Vector const &, Matrix const &>(), "weight"_a, "mu"_a, "sigma"_a);
63 utils::python::addOutputOp(cls, "__str__");
64 utils::python::addOutputOp(cls, "__repr__");
65 return cls;
66}
67
68static PyMixtureUpdateRestriction declareMixtureUpdateRestriction(py::module &mod) {
69 PyMixtureUpdateRestriction cls(mod, "MixtureUpdateRestriction");
70 cls.def("getDimension", &MixtureUpdateRestriction::getDimension);
71 cls.def(py::init<int>(), "dim"_a);
72 // The rest of this interface isn't usable in Python, and doesn't need to be.
73 return cls;
74}
75
76static PyMixture declareMixture(py::module &mod) {
77 PyMixture cls(mod, "Mixture");
78 afw::table::io::python::addPersistableMethods<Mixture>(cls);
79 cls.def("__iter__", [](Mixture &self) { return py::make_iterator(self.begin(), self.end()); },
80 py::keep_alive<0, 1>());
81 cls.def("__getitem__",
82 [](Mixture &self, std::ptrdiff_t i) { return self[utils::python::cppIndex(self.size(), i)]; },
83 py::return_value_policy::reference_internal);
84 cls.def("__len__", &Mixture::size);
85 cls.def("getComponentCount", &Mixture::getComponentCount);
86 cls.def("project", (std::shared_ptr<Mixture> (Mixture::*)(int) const) & Mixture::project, "dim"_a);
87 cls.def("project", (std::shared_ptr<Mixture> (Mixture::*)(int, int) const) & Mixture::project, "dim1"_a,
88 "dim2"_a);
89 cls.def("getDimension", &Mixture::getDimension);
90 cls.def("normalize", &Mixture::normalize);
91 cls.def("shift", &Mixture::shift, "dim"_a, "offset"_a);
92 cls.def("clip", &Mixture::clip, "threshold"_a = 0.0);
93 cls.def("getDegreesOfFreedom", &Mixture::getDegreesOfFreedom);
94 cls.def("setDegreesOfFreedom", &Mixture::setDegreesOfFreedom,
96 cls.def("evaluate",
97 [](Mixture const &self, MixtureComponent const &component,
98 ndarray::Array<Scalar, 1, 0> const &array) -> Scalar {
99 return self.evaluate(component, ndarray::asEigenMatrix(array));
100 },
101 "component"_a, "x"_a);
102 cls.def("evaluate",
103 [](Mixture const &self, ndarray::Array<Scalar, 1, 0> const &array) -> Scalar {
104 return self.evaluate(ndarray::asEigenMatrix(array));
105 },
106 "x"_a);
107 cls.def("evaluate", (void (Mixture::*)(ndarray::Array<Scalar const, 2, 1> const &,
108 ndarray::Array<Scalar, 1, 0> const &) const) &
110 "x"_a, "p"_a);
111 cls.def("evaluateComponents", &Mixture::evaluateComponents, "x"_a, "p"_a);
112 cls.def("evaluateDerivatives",
113 py::overload_cast<ndarray::Array<Scalar const, 1, 1> const &,
114 ndarray::Array<Scalar,1,1> const &,
115 ndarray::Array<Scalar,2,1> const &>(&Mixture::evaluateDerivatives, py::const_),
116 "x"_a, "gradient"_a, "hessian"_a);
117 cls.def("draw", &Mixture::draw, "rng"_a, "x"_a);
118 cls.def("updateEM", (void (Mixture::*)(ndarray::Array<Scalar const, 2, 1> const &,
119 ndarray::Array<Scalar const, 1, 0> const &, Scalar, Scalar)) &
121 "x"_a, "w"_a, "tau1"_a = 0.0, "tau2"_a = 0.5);
122 cls.def("updateEM", (void (Mixture::*)(ndarray::Array<Scalar const, 2, 1> const &,
123 ndarray::Array<Scalar const, 1, 0> const &,
124 MixtureUpdateRestriction const &restriction, Scalar, Scalar)) &
126 "x"_a, "w"_a, "restriction"_a, "tau1"_a = 0.0, "tau2"_a = 0.5);
127 cls.def("updateEM", (void (Mixture::*)(ndarray::Array<Scalar const, 2, 1> const &,
128 MixtureUpdateRestriction const &restriction, Scalar, Scalar)) &
130 "x"_a, "restriction"_a, "tau1"_a = 0.0, "tau2"_a = 0.5);
131 cls.def("clone", &Mixture::clone);
132 cls.def(py::init<int, Mixture::ComponentList &, Scalar>(), "dim"_a, "components"_a,
134 utils::python::addOutputOp(cls, "__str__");
135 utils::python::addOutputOp(cls, "__repr__");
136 return cls;
137}
138
139PYBIND11_MODULE(mixture, mod) {
140 py::module::import("lsst.afw.math");
141
142 auto clsMixtureComponent = declareMixtureComponent(mod);
143 auto clsMixtureUpdateRestriction = declareMixtureUpdateRestriction(mod);
144 auto clsMixture = declareMixture(mod);
145 clsMixture.attr("Component") = clsMixtureComponent;
146 clsMixture.attr("UpdateRestriction") = clsMixtureUpdateRestriction;
147}
148
149}
150}
151}
152} // namespace lsst::meas::modelfit::anonymous
int end
T begin(T... args)
MixtureComponent project(int dim) const
Project the distribution onto the given dimension (marginalize over all others)
Vector getMu() const
Get/set the location parameter (mean/median/mode) of this component.
Definition: Mixture.h:58
Scalar weight
Weight of this distribution in the mixture.
Definition: Mixture.h:54
void setSigma(Matrix const &sigma)
void setMu(Vector const &mu)
Definition: Mixture.h:59
Matrix getSigma() const
Get/set the shape/size parameter.
Definition: Mixture.h:71
int getDimension() const
Return the number of dimensions.
Definition: Mixture.h:51
void updateEM(ndarray::Array< Scalar const, 2, 1 > const &x, ndarray::Array< Scalar const, 1, 0 > const &w, Scalar tau1=0.0, Scalar tau2=0.5)
Perform an Expectation-Maximization step, updating the component parameters to match the given weight...
Scalar getDegreesOfFreedom() const
Get the number of degrees of freedom in the component Student's T distributions (inf=Gaussian)
Definition: Mixture.h:187
std::size_t size() const
Return the number of components.
Definition: Mixture.h:157
std::shared_ptr< Mixture > project(int dim) const
Project the distribution onto the given dimensions (marginalize over all others)
Scalar evaluate(Component const &component, Eigen::MatrixBase< Derived > const &x) const
Evaluate the probability density at the given point for the given component distribution.
Definition: Mixture.h:198
void normalize()
Iterate over all components, rescaling their weights so they sum to one.
void evaluateComponents(ndarray::Array< Scalar const, 2, 1 > const &x, ndarray::Array< Scalar, 2, 1 > const &p) const
Evaluate the contributions of each component to the full probability at the given points.
void evaluateDerivatives(ndarray::Array< Scalar const, 1, 1 > const &x, ndarray::Array< Scalar, 1, 1 > const &gradient, ndarray::Array< Scalar, 2, 1 > const &hessian) const
Evaluate the derivative of the distribution at the given point.
void shift(int dim, Scalar offset)
Shift the mixture in the given dimension, adding the given offset to all mu vectors.
virtual std::shared_ptr< Mixture > clone() const
Polymorphic deep copy.
virtual int getComponentCount() const
Return the number of components.
Definition: Mixture.h:160
void draw(afw::math::Random &rng, ndarray::Array< Scalar, 2, 1 > const &x) const
Draw random variates from the distribution.
int getDimension() const
Return the number of dimensions.
Definition: Mixture.h:169
void setDegreesOfFreedom(Scalar df=std::numeric_limits< Scalar >::infinity())
Set the number of degrees of freedom in the component Student's T distributions (inf=Gaussian)
std::size_t clip(Scalar threshold=0.0)
Iterate over all components, removing those with weight less than or equal to threshold.
T infinity(T... args)
PYBIND11_MODULE(_cameraGeom, mod)
Definition: _cameraGeom.cc:38
void addOutputOp(PyClass &cls, std::string const &method)
Add __str__ or __repr__ method implemented by operator<<.
Definition: python.h:87
std::size_t cppIndex(std::ptrdiff_t size, std::ptrdiff_t i)
Compute a C++ index from a Python index (negative values count from the end) and range-check.
Definition: python.h:124
double Scalar
Typedefs to be used for probability and parameter values.
Definition: common.h:44
A base class for image defects.