gtsam 4.2.0
gtsam
Loading...
Searching...
No Matches
gtsam::NoiseModelFactorN< ValueTypes > Class Template Referenceabstract

Detailed Description

template<class... ValueTypes>
class gtsam::NoiseModelFactorN< ValueTypes >

A convenient base class for creating your own NoiseModelFactor with n variables.

To derive from this class, implement evaluateError().

For example, a 2-way factor that computes the difference in x-translation between a Pose3 and Point3 could be implemented like so:

class MyFactor : public NoiseModelFactorN<Pose3, Point3> {
public:
MyFactor(Key pose_key, Key point_key, const SharedNoiseModel& noiseModel)
: Base(noiseModel, pose_key, point_key) {}
Vector evaluateError(
const Pose3& T, const Point3& p,
boost::optional<Matrix&> H_T = boost::none,
boost::optional<Matrix&> H_p = boost::none) const override {
Matrix36 t_H_T; // partial derivative of translation w.r.t. pose T
// Only compute t_H_T if needed:
Point3 t = T.translation(H_T ? &t_H_T : 0);
double a = t(0); // a_H_t = [1, 0, 0]
double b = p(0); // b_H_p = [1, 0, 0]
double error = a - b; // H_a = 1, H_b = -1
// H_T = H_a * a_H_t * t_H_T = the first row of t_H_T
if (H_T) *H_T = (Matrix(1, 6) << t_H_T.row(0)).finished();
// H_p = H_b * b_H_p = -1 * [1, 0, 0]
if (H_p) *H_p = (Matrix(1, 3) << -1., 0., 0.).finished();
return Vector1(error);
}
};
// Unit Test
TEST(NonlinearFactor, MyFactor) {
MyFactor f(X(1), X(2), noiseModel::Unit::Create(1));
EXPECT_DOUBLES_EQUAL(-8., f.evaluateError(Pose3(), Point3(8., 7., 6.))(0),
1e-9);
Values values;
values.insert(X(1), Pose3(Rot3::RzRyRx(0.1, 0.2, 0.3), Point3(1, 2, 3)));
values.insert(X(2), Point3(1, 2, 3));
EXPECT_CORRECT_FACTOR_JACOBIANS(f, values, 1e-5, 1e-5);
}
#define EXPECT_CORRECT_FACTOR_JACOBIANS(factor, values, numerical_derivative_step, tolerance)
Check the Jacobians produced by a factor against finite differences.
Definition factorTesting.h:114
Vector3 Point3
As of GTSAM 4, in order to make GTSAM more lean, it is now possible to just typedef Point3 to Vector3...
Definition Point3.h:36
noiseModel::Base::shared_ptr SharedNoiseModel
Aliases.
Definition NoiseModel.h:724
std::uint64_t Key
Integer nonlinear key type.
Definition types.h:100
A 3D pose (R,t) : (Rot3,Point3)
Definition Pose3.h:37
static Rot3 RzRyRx(double x, double y, double z, OptionalJacobian< 3, 1 > Hx=boost::none, OptionalJacobian< 3, 1 > Hy=boost::none, OptionalJacobian< 3, 1 > Hz=boost::none)
Rotations around Z, Y, then X axes as in http://en.wikipedia.org/wiki/Rotation_matrix,...
Definition Rot3M.cpp:85
static shared_ptr Create(size_t dim)
Create a unit covariance noise model.
Definition NoiseModel.h:597
A nonlinear sum-of-squares factor with a zero-mean noise model implementing the density Templated on...
Definition NonlinearFactor.h:174
A convenient base class for creating your own NoiseModelFactor with n variables.
Definition NonlinearFactor.h:400
In nonlinear factors, the error function returns the negative log-likelihood as a non-linear function...

These factors are templated on a values structure type. The values structures are typically more general than just vectors, e.g., Rot3 or Pose3, which are objects in non-linear manifolds (Lie groups).

+ Inheritance diagram for gtsam::NoiseModelFactorN< ValueTypes >:

Public Member Functions

template<int I = 1>
Key key () const
 Returns a key.
 
Constructors
 NoiseModelFactorN ()
 Default Constructor for I/O.
 
 NoiseModelFactorN (const SharedNoiseModel &noiseModel, KeyType< ValueTypes >... keys)
 Constructor.
 
template<typename CONTAINER = std::initializer_list<Key>, typename = IsContainerOfKeys<CONTAINER>>
 NoiseModelFactorN (const SharedNoiseModel &noiseModel, CONTAINER keys)
 Constructor.
 
NoiseModelFactor methods
Vector unwhitenedError (const Values &x, boost::optional< std::vector< Matrix > & > H=boost::none) const override
 This implements the unwhitenedError virtual function by calling the n-key specific version of evaluateError, which is pure virtual so must be implemented in the derived class.
 
Virtual methods
virtual Vector evaluateError (const ValueTypes &... x, OptionalMatrix< ValueTypes >... H) const =0
 Override evaluateError to finish implementing an n-way factor.
 
Convenience method overloads
Vector evaluateError (const ValueTypes &... x) const
 No-Jacobians requested function overload.
 
template<typename... OptionalJacArgs, typename = IndexIsValid<sizeof...(OptionalJacArgs) + 1>>
Vector evaluateError (const ValueTypes &... x, OptionalJacArgs &&... H) const
 Some (but not all) optional Jacobians are omitted (function overload)
 
Shortcut functions <tt>key1()</tt> -> <tt>key\<1\>()</tt>
Key key1 () const
 
template<int I = 2>
Key key2 () const
 
template<int I = 3>
Key key3 () const
 
template<int I = 4>
Key key4 () const
 
template<int I = 5>
Key key5 () const
 
template<int I = 6>
Key key6 () const
 
- Public Member Functions inherited from gtsam::NoiseModelFactor
 NoiseModelFactor ()
 Default constructor for I/O only.
 
 ~NoiseModelFactor () override
 Destructor.
 
template<typename CONTAINER >
 NoiseModelFactor (const SharedNoiseModel &noiseModel, const CONTAINER &keys)
 Constructor.
 
void print (const std::string &s="", const KeyFormatter &keyFormatter=DefaultKeyFormatter) const override
 Print.
 
bool equals (const NonlinearFactor &f, double tol=1e-9) const override
 Check if two factors are equal.
 
size_t dim () const override
 get the dimension of the factor (number of rows on linearization)
 
const SharedNoiseModelnoiseModel () const
 access to the noise model
 
Vector whitenedError (const Values &c) const
 Vector of errors, whitened This is the raw error, i.e., i.e.
 
Vector unweightedWhitenedError (const Values &c) const
 Vector of errors, whitened, but unweighted by any loss function.
 
double weight (const Values &c) const
 Compute the effective weight of the factor from the noise model.
 
double error (const Values &c) const override
 Calculate the error of the factor.
 
boost::shared_ptr< GaussianFactorlinearize (const Values &x) const override
 Linearize a non-linearFactorN to get a GaussianFactor, \( Ax-b \approx h(x+\delta x)-z = h(x) + A \delta x - z \) Hence \( b = z - h(x) = - \mathtt{error\_vector}(x) \).
 
shared_ptr cloneWithNewNoiseModel (const SharedNoiseModel newNoise) const
 Creates a shared_ptr clone of the factor with a new noise model.
 
- Public Member Functions inherited from gtsam::NonlinearFactor
 NonlinearFactor ()
 Default constructor for I/O only.
 
template<typename CONTAINER >
 NonlinearFactor (const CONTAINER &keys)
 Constructor from a collection of the keys involved in this factor.
 
void print (const std::string &s="", const KeyFormatter &keyFormatter=DefaultKeyFormatter) const override
 print
 
virtual ~NonlinearFactor ()
 Destructor.
 
double error (const HybridValues &c) const override
 All factor types need to implement an error function.
 
virtual bool active (const Values &) const
 Checks whether a factor should be used based on a set of values.
 
virtual shared_ptr clone () const
 Creates a shared_ptr clone of the factor - needs to be specialized to allow for subclasses.
 
virtual shared_ptr rekey (const std::map< Key, Key > &rekey_mapping) const
 Creates a shared_ptr clone of the factor with different keys using a map from old->new keys.
 
virtual shared_ptr rekey (const KeyVector &new_keys) const
 Clones a factor and fully replaces its keys.
 
virtual bool sendable () const
 Should the factor be evaluated in the same thread as the caller This is to enable factors that has shared states (like the Python GIL lock)
 
- Public Member Functions inherited from gtsam::Factor
virtual ~Factor ()=default
 Default destructor.
 
bool empty () const
 Whether the factor is empty (involves zero variables).
 
Key front () const
 First key.
 
Key back () const
 Last key.
 
const_iterator find (Key key) const
 find
 
const KeyVectorkeys () const
 Access the factor's involved variable keys.
 
const_iterator begin () const
 Iterator at beginning of involved variable keys.
 
const_iterator end () const
 Iterator at end of involved variable keys.
 
size_t size () const
 
virtual void printKeys (const std::string &s="Factor", const KeyFormatter &formatter=DefaultKeyFormatter) const
 print only keys
 
bool equals (const This &other, double tol=1e-9) const
 check equality
 
KeyVectorkeys ()
 
iterator begin ()
 Iterator at beginning of involved variable keys.
 
iterator end ()
 Iterator at end of involved variable keys.
 

Public Types

enum  { N = sizeof...(ValueTypes) }
 N is the number of variables (N-way factor)
 
template<int I, typename = IndexIsValid<I>>
using ValueType = typename std::tuple_element< I - 1, std::tuple< ValueTypes... > >::type
 The type of the I'th template param can be obtained as ValueType.
 
- Public Types inherited from gtsam::NoiseModelFactor
typedef boost::shared_ptr< Thisshared_ptr
 Noise model.
 
- Public Types inherited from gtsam::NonlinearFactor
typedef boost::shared_ptr< Thisshared_ptr
 
- Public Types inherited from gtsam::Factor
typedef KeyVector::iterator iterator
 Iterator over keys.
 
typedef KeyVector::const_iterator const_iterator
 Const iterator over keys.
 

Protected Types

using Base = NoiseModelFactor
 
using This = NoiseModelFactorN< ValueTypes... >
 
template<typename T >
using OptionalMatrix = boost::optional< Matrix & >
 
template<typename T >
using KeyType = Key
 
SFINAE aliases
template<typename From , typename To >
using IsConvertible = typename std::enable_if< std::is_convertible< From, To >::value, void >::type
 
template<int I>
using IndexIsValid = typename std::enable_if<(I >=1) &&(I<=N), void >::type
 
template<typename Container >
using ContainerElementType = typename std::decay< decltype(*std::declval< Container >().begin())>::type
 
template<typename Container >
using IsContainerOfKeys = IsConvertible< ContainerElementType< Container >, Key >
 
- Protected Types inherited from gtsam::NoiseModelFactor
typedef NonlinearFactor Base
 
typedef NoiseModelFactor This
 
- Protected Types inherited from gtsam::NonlinearFactor
typedef Factor Base
 
typedef NonlinearFactor This
 

Friends

class boost::serialization::access
 Serialization function.
 

Additional Inherited Members

- Protected Member Functions inherited from gtsam::NoiseModelFactor
 NoiseModelFactor (const SharedNoiseModel &noiseModel)
 Constructor - only for subclasses, as this does not set keys.
 
- Protected Member Functions inherited from gtsam::Factor
 Factor ()
 Default constructor for I/O.
 
template<typename CONTAINER >
 Factor (const CONTAINER &keys)
 Construct factor from container of keys.
 
template<typename ITERATOR >
 Factor (ITERATOR first, ITERATOR last)
 Construct factor from iterator keys.
 
- Static Protected Member Functions inherited from gtsam::Factor
template<typename CONTAINER >
static Factor FromKeys (const CONTAINER &keys)
 Construct factor from container of keys.
 
template<typename ITERATOR >
static Factor FromIterators (ITERATOR first, ITERATOR last)
 Construct factor from iterator keys.
 
- Protected Attributes inherited from gtsam::NoiseModelFactor
SharedNoiseModel noiseModel_
 
- Protected Attributes inherited from gtsam::Factor
KeyVector keys_
 The keys involved in this factor.
 

Member Typedef Documentation

◆ ValueType

template<class... ValueTypes>
template<int I, typename = IndexIsValid<I>>
using gtsam::NoiseModelFactorN< ValueTypes >::ValueType = typename std::tuple_element<I - 1, std::tuple<ValueTypes...> >::type

The type of the I'th template param can be obtained as ValueType.

I is 1-indexed for backwards compatibility/consistency! So for example,

Factor::ValueType<1> // Pose3
Factor::ValueType<2> // Point3
// Factor::ValueType<0> // ERROR! Will not compile.
// Factor::ValueType<3> // ERROR! Will not compile.
Definition Factor.h:68

You can also use the shortcuts X1, ..., X6 which are the same as ValueType<1>, ..., ValueType<6> respectively (see detail::NoiseModelFactorAliases).

Note that, if your class is templated AND you want to use ValueType<1> inside your class, due to dependent types you need the template keyword: typename MyFactor<T>::template ValueType<1>.

Constructor & Destructor Documentation

◆ NoiseModelFactorN() [1/2]

template<class... ValueTypes>
gtsam::NoiseModelFactorN< ValueTypes >::NoiseModelFactorN ( const SharedNoiseModel noiseModel,
KeyType< ValueTypes >...  keys 
)
inline

Constructor.

Example usage: NoiseModelFactorN(noise, key1, key2, ..., keyN)

Parameters
noiseModelShared pointer to noise model.
keysKeys for the variables in this factor, passed in as separate arguments.

◆ NoiseModelFactorN() [2/2]

template<class... ValueTypes>
template<typename CONTAINER = std::initializer_list<Key>, typename = IsContainerOfKeys<CONTAINER>>
gtsam::NoiseModelFactorN< ValueTypes >::NoiseModelFactorN ( const SharedNoiseModel noiseModel,
CONTAINER  keys 
)
inline

Constructor.

Example usage: NoiseModelFactorN(noise, {key1, key2, ..., keyN}) Example usage: NoiseModelFactorN(noise, keys) where keys is a vector<Key>

Parameters
noiseModelShared pointer to noise model.
keysA container of keys for the variables in this factor.

Member Function Documentation

◆ evaluateError() [1/3]

template<class... ValueTypes>
Vector gtsam::NoiseModelFactorN< ValueTypes >::evaluateError ( const ValueTypes &...  x) const
inline

No-Jacobians requested function overload.

This specializes the version below to avoid recursive calls since this is commonly used.

e.g. const Vector error = factor.evaluateError(pose, point);

◆ evaluateError() [2/3]

template<class... ValueTypes>
template<typename... OptionalJacArgs, typename = IndexIsValid<sizeof...(OptionalJacArgs) + 1>>
Vector gtsam::NoiseModelFactorN< ValueTypes >::evaluateError ( const ValueTypes &...  x,
OptionalJacArgs &&...  H 
) const
inline

Some (but not all) optional Jacobians are omitted (function overload)

e.g. const Vector error = factor.evaluateError(pose, point, Hpose);

◆ evaluateError() [3/3]

template<class... ValueTypes>
virtual Vector gtsam::NoiseModelFactorN< ValueTypes >::evaluateError ( const ValueTypes &...  x,
OptionalMatrix< ValueTypes >...  H 
) const
pure virtual

Override evaluateError to finish implementing an n-way factor.

Both the x and H arguments are written here as parameter packs, but when overriding this method, you probably want to explicitly write them out. For example, for a 2-way factor with variable types Pose3 and Point3, you should implement:

const Pose3& x1, const Point3& x2,
boost::optional<Matrix&> H1 = boost::none,
boost::optional<Matrix&> H2 = boost::none) const override { ... }
virtual Vector evaluateError(const ValueTypes &... x, OptionalMatrix< ValueTypes >... H) const =0
Override evaluateError to finish implementing an n-way factor.

If any of the optional Matrix reference arguments are specified, it should compute both the function evaluation and its derivative(s) in the requested variables.

Parameters
xThe values of the variables to evaluate the error for. Passed in as separate arguments.
[out]HThe Jacobian with respect to each variable (optional).

◆ key()

template<class... ValueTypes>
template<int I = 1>
Key gtsam::NoiseModelFactorN< ValueTypes >::key ( ) const
inline

Returns a key.

Usage: key<I>() returns the I'th key. I is 1-indexed for backwards compatibility/consistency! So for example,

NoiseModelFactorN<Pose3, Point3> factor(noise, key1, key2);
key<1>() // = key1
key<2>() // = key2
// key<0>() // ERROR! Will not compile
// key<3>() // ERROR! Will not compile
Key key() const
Returns a key.
Definition NonlinearFactor.h:518

Note that, if your class is templated AND you are trying to call key<1> inside your class, due to dependent types you need the template keyword: this->key1().

◆ unwhitenedError()

template<class... ValueTypes>
Vector gtsam::NoiseModelFactorN< ValueTypes >::unwhitenedError ( const Values x,
boost::optional< std::vector< Matrix > & >  H = boost::none 
) const
inlineoverridevirtual

This implements the unwhitenedError virtual function by calling the n-key specific version of evaluateError, which is pure virtual so must be implemented in the derived class.

Example usage:

values.insert(...) // populate values
std::vector<Matrix> Hs(2); // this will be an optional output argument
const Vector error = factor.unwhitenedError(values, Hs);
double error(const Values &c) const override
Calculate the error of the factor.
Definition NonlinearFactor.cpp:138
A non-templated config holding any types of Manifold-group elements.
Definition Values.h:65
void insert(Key j, const Value &val)
Add a variable with the given j, throws KeyAlreadyExists<J> if j is already present.
Definition Values.cpp:157
Parameters
[in]xA Values object containing the values of all the variables used in this factor
[out]HA vector of (dynamic) matrices whose size should be equal to n. The Jacobians w.r.t. each variable will be output in this parameter.

Implements gtsam::NoiseModelFactor.


The documentation for this class was generated from the following file: