gtsam 4.2.0
gtsam
GaussianConditional.h
Go to the documentation of this file.
1/* ----------------------------------------------------------------------------
2
3 * GTSAM Copyright 2010, Georgia Tech Research Corporation,
4 * Atlanta, Georgia 30332-0415
5 * All Rights Reserved
6 * Authors: Frank Dellaert, et al. (see THANKS for the full author list)
7
8 * See LICENSE for the license information
9
10 * -------------------------------------------------------------------------- */
11
18// \callgraph
19
20#pragma once
21
22#include <boost/utility.hpp>
23
27#include <gtsam/inference/Conditional-inst.h>
29
30#include <random> // for std::mt19937_64
31
32namespace gtsam {
33
40 class GTSAM_EXPORT GaussianConditional :
41 public JacobianFactor,
42 public Conditional<JacobianFactor, GaussianConditional>
43 {
44 public:
46 typedef boost::shared_ptr<This> shared_ptr;
49
52
55
57 GaussianConditional(Key key, const Vector& d, const Matrix& R,
58 const SharedDiagonal& sigmas = SharedDiagonal());
59
61 GaussianConditional(Key key, const Vector& d, const Matrix& R, Key parent1,
62 const Matrix& S,
63 const SharedDiagonal& sigmas = SharedDiagonal());
64
66 GaussianConditional(Key key, const Vector& d, const Matrix& R, Key parent1,
67 const Matrix& S, Key parent2, const Matrix& T,
68 const SharedDiagonal& sigmas = SharedDiagonal());
69
73 template<typename TERMS>
74 GaussianConditional(const TERMS& terms,
75 size_t nrFrontals, const Vector& d,
76 const SharedDiagonal& sigmas = SharedDiagonal());
77
82 template<typename KEYS>
84 const KEYS& keys, size_t nrFrontals, const VerticalBlockMatrix& augmentedMatrix,
85 const SharedDiagonal& sigmas = SharedDiagonal());
86
88 static GaussianConditional FromMeanAndStddev(Key key, const Vector& mu,
89 double sigma);
90
92 static GaussianConditional FromMeanAndStddev(Key key, const Matrix& A,
93 Key parent, const Vector& b,
94 double sigma);
95
98 static GaussianConditional FromMeanAndStddev(Key key, //
99 const Matrix& A1, Key parent1,
100 const Matrix& A2, Key parent2,
101 const Vector& b, double sigma);
102
104 template<typename... Args>
105 static shared_ptr sharedMeanAndStddev(Args&&... args) {
106 return boost::make_shared<This>(FromMeanAndStddev(std::forward<Args>(args)...));
107 }
108
116 template<typename ITERATOR>
117 static shared_ptr Combine(ITERATOR firstConditional, ITERATOR lastConditional);
118
122
124 void print(
125 const std::string& = "GaussianConditional",
126 const KeyFormatter& formatter = DefaultKeyFormatter) const override;
127
129 bool equals(const GaussianFactor&cg, double tol = 1e-9) const override;
130
134
139 double logNormalizationConstant() const override;
140
148 double logProbability(const VectorValues& x) const;
149
155 double evaluate(const VectorValues& x) const;
156
158 double operator()(const VectorValues& x) const {
159 return evaluate(x);
160 }
161
175 VectorValues solve(const VectorValues& parents) const;
176
177 VectorValues solveOtherRHS(const VectorValues& parents, const VectorValues& rhs) const;
178
180 void solveTransposeInPlace(VectorValues& gy) const;
181
184 const VectorValues& frontalValues) const;
185
187 JacobianFactor::shared_ptr likelihood(const Vector& frontal) const;
188
195 VectorValues sample(std::mt19937_64* rng) const;
196
204 VectorValues sample(const VectorValues& parentsValues,
205 std::mt19937_64* rng) const;
206
208 VectorValues sample() const;
209
211 VectorValues sample(const VectorValues& parentsValues) const;
212
216
218 constABlock R() const { return Ab_.range(0, nrFrontals()); }
219
221 constABlock S() const { return Ab_.range(nrFrontals(), size()); }
222
224 constABlock S(const_iterator it) const { return BaseFactor::getA(it); }
225
227 const constBVector d() const { return BaseFactor::getb(); }
228
240 inline double determinant() const { return exp(logDeterminant()); }
241
253 double logDeterminant() const;
254
258
263 double logProbability(const HybridValues& x) const override;
264
269 double evaluate(const HybridValues& x) const override;
270
271 using Conditional::operator(); // Expose evaluate(const HybridValues&) method..
272 using JacobianFactor::error; // Expose error(const HybridValues&) method..
273
275
276
277#ifdef GTSAM_ALLOW_DEPRECATED_SINCE_V42
280
282 void GTSAM_DEPRECATED scaleFrontalsBySigma(VectorValues& gy) const;
284#endif
285
286 private:
288 friend class boost::serialization::access;
289 template<class Archive>
290 void serialize(Archive & ar, const unsigned int /*version*/) {
291 ar & BOOST_SERIALIZATION_BASE_OBJECT_NVP(BaseFactor);
292 ar & BOOST_SERIALIZATION_BASE_OBJECT_NVP(BaseConditional);
293 }
294 }; // GaussianConditional
295
297template<>
298struct traits<GaussianConditional> : public Testable<GaussianConditional> {};
299
300} // \ namespace gtsam
301
303
Base class for conditional densities.
Conditional Gaussian Base class.
Factor Graph Values.
Included from all GTSAM files.
Global functions in a separate testing namespace.
Definition: chartTesting.h:28
std::string serialize(const T &input)
serializes to a string
Definition: serialization.h:113
void print(const Matrix &A, const string &s, ostream &stream)
print without optional string, must specify cout yourself
Definition: Matrix.cpp:156
std::uint64_t Key
Integer nonlinear key type.
Definition: types.h:100
std::function< std::string(Key)> KeyFormatter
Typedef for a function to format a key, i.e. to convert it to a string.
Definition: Key.h:35
A manifold defines a space in which there is a notion of a linear tangent space that can be centered ...
Definition: concepts.h:30
Template to create a binary predicate.
Definition: Testable.h:111
A helper that implements the traits interface for GTSAM types.
Definition: Testable.h:151
This class stores a dense matrix and allows it to be accessed as a collection of vertical blocks.
Definition: VerticalBlockMatrix.h:43
HybridValues represents a collection of DiscreteValues and VectorValues.
Definition: HybridValues.h:38
Definition: Conditional.h:64
KeyVector::const_iterator const_iterator
Const iterator over keys.
Definition: Factor.h:80
A GaussianConditional functions as the node in a Bayes network.
Definition: GaussianConditional.h:43
GaussianConditional This
Typedef to this class.
Definition: GaussianConditional.h:45
constABlock S(const_iterator it) const
Get a view of the S matrix for the variable pointed to by the given key iterator.
Definition: GaussianConditional.h:224
constABlock R() const
Return a view of the upper-triangular R block of the conditional.
Definition: GaussianConditional.h:218
JacobianFactor BaseFactor
Typedef to our factor base class.
Definition: GaussianConditional.h:47
GaussianConditional()
default constructor needed for serialization
Definition: GaussianConditional.h:54
Conditional< BaseFactor, This > BaseConditional
Typedef to our conditional base class.
Definition: GaussianConditional.h:48
static shared_ptr Combine(ITERATOR firstConditional, ITERATOR lastConditional)
Combine several GaussianConditional into a single dense GC.
double determinant() const
Compute the determinant of the R matrix.
Definition: GaussianConditional.h:240
double operator()(const VectorValues &x) const
Evaluate probability density, sugar.
Definition: GaussianConditional.h:158
static shared_ptr sharedMeanAndStddev(Args &&... args)
Create shared pointer by forwarding arguments to fromMeanAndStddev.
Definition: GaussianConditional.h:105
constABlock S() const
Get a view of the parent blocks.
Definition: GaussianConditional.h:221
boost::shared_ptr< This > shared_ptr
shared_ptr to this class
Definition: GaussianConditional.h:46
const constBVector d() const
Get a view of the r.h.s.
Definition: GaussianConditional.h:227
An abstract virtual base class for JacobianFactor and HessianFactor.
Definition: GaussianFactor.h:39
A Gaussian factor in the squared-error form.
Definition: JacobianFactor.h:91
boost::shared_ptr< This > shared_ptr
shared_ptr to this class
Definition: JacobianFactor.h:96
VectorValues represents a collection of vector-valued variables associated each with a unique integer...
Definition: VectorValues.h:74