Belle II Software development
Public Member Functions | Private Attributes | List of all members
MultivariateNormalGenerator Class Reference
Packages » framework

Class to generate normal distributed, correlated random numbers given the mean values and the covariance matrix of all dimensions. More...

#include <MultivariateNormalGenerator.h>

Public Member Functions

 MultivariateNormalGenerator ()
 default constructor to allow later initialization
 
 MultivariateNormalGenerator (int n, const double *mean, const double *cov)
 constructor with array interface: mean and covariance are passed as double arrays where the covariance is expected to be an NxN matrix in row major layout
 
 MultivariateNormalGenerator (const Eigen::VectorXd &mean, const Eigen::MatrixXd &cov)
 constructor with Eigen matrix interface.
 
Eigen::VectorXd generate () const
 Generate a set of correlated random numbers with the previouly set mean and covariance.
 
void reset ()
 reset the generator setting the size to 0.
 
size_t size () const
 Return the number of elements to be generated on generate()
 
void generate (double *output) const
 Generate a set of correlated random numbers with the previouly set mean and covariance and store them in buffer output.
 
TVector3 generateVec3 () const
 Generate a set of correlated random numbers with the previouly set mean and covariance and return a TVector3.
 
TVectorD generateVecT () const
 Generate a set of correlated random numbers with the previouly set mean and covariance and return a TVectorT<double>
 
bool setMeanCov (int n, const double *mean, const double *cov)
 set the mean and covariance for the distribution with array interface: mean and covariance are passed as double arrays where the covariance is expected to be an NxN matrix in row major layout
 
bool setMeanCov (const Eigen::VectorXd &mean, const Eigen::MatrixXd &cov)
 set the mean and covariance for the distribution.
 
template<class value_type >
bool setMeanCov (const TVectorT< value_type > &mean, const TMatrixTBase< value_type > &cov)
 set mean and covariance matrix from ROOT vector/matrix objects, e.g.
 
template<class value_type >
bool setMeanCov (const ROOT::Math::XYZVector &mean, const TMatrixTBase< value_type > &cov)
 set the mean and covariance for the distribution.
 

Private Attributes

Eigen::VectorXd m_mean
 Member to store the mean values of the distribution.
 
Eigen::MatrixXd m_transform
 Member to store the transformation matrix for standard normal distributed random values.
 

Detailed Description

Class to generate normal distributed, correlated random numbers given the mean values and the covariance matrix of all dimensions.

This class can be used to generate normal distributed random values according to a given covariance matrix (assuming the covariance matrix is positive semi-definite).

To use it first set the desired mean values and covariance matrix using setMeanCov() and then call generate() to generate one set of values.

Warning
: setMeanCov() will not work for all matrices. Please check the return value when setting the covariance matrix.

To generate normal distributed random values according to a covariance matrix we need to decompose the covariance matrix $M$ into $M= A A^T$. Given the vector of mean values as $\mu$ and a vector of standard normal distributed random values $(\mu=0, \sigma=1)$ as n we can obtain a set of correlated random values $x = \mu + A * n$.

Usually the Cholesky decomposition is chosen as it computes $M = L L^T$. However the Cholesky composition only works for positive definite matrices so it does not work if for example one of the values is fixed and has no error.

To ease this restriction a little we use the LDLT decomposition given as $M = P^T L D L^T P$ where $P$ is a permutation matrix and D a diagonal matrix. We then can use $A = P^T L \sqrt(D)$ to caluclate the correlated values also for positive semi-definite covariance matrices if the elements of D are positive.

Definition at line 53 of file MultivariateNormalGenerator.h.

Constructor & Destructor Documentation

◆ MultivariateNormalGenerator() [1/3]

MultivariateNormalGenerator ( )
inline

default constructor to allow later initialization

Definition at line 56 of file MultivariateNormalGenerator.h.

56{}

◆ MultivariateNormalGenerator() [2/3]

MultivariateNormalGenerator ( int  n,
const double *  mean,
const double *  cov 
)
inline

constructor with array interface: mean and covariance are passed as double arrays where the covariance is expected to be an NxN matrix in row major layout

Parameters
ndimensionality
meanpointer to the n mean values of the distribution
covpointer to the n*n covariance values in row major layout

Definition at line 64 of file MultivariateNormalGenerator.h.

65 {
66 setMeanCov(n, mean, cov);
67 }
Belle2::MultivariateNormalGenerator::setMeanCov
bool setMeanCov(int n, const double *mean, const double *cov)
set the mean and covariance for the distribution with array interface: mean and covariance are passed...
Definition: MultivariateNormalGenerator.cc:20

◆ MultivariateNormalGenerator() [3/3]

MultivariateNormalGenerator ( const Eigen::VectorXd &  mean,
const Eigen::MatrixXd &  cov 
)
inline

constructor with Eigen matrix interface.

Parameters
meanVector of mean values
covMatrix containing the covariance values

Definition at line 72 of file MultivariateNormalGenerator.h.

73 {
74 setMeanCov(mean, cov);
75 }

Member Function Documentation

◆ generate() [1/2]

Eigen::VectorXd generate ( ) const
inline

Generate a set of correlated random numbers with the previouly set mean and covariance.

Returns
Vector containing the generated random numbers

Definition at line 80 of file MultivariateNormalGenerator.h.

81 {
82 //To get the correlated multivariate normal distribution, we multiply
83 //standard normal distributed values (mean=0, sigma=1) with a
84 //transformation matrix we obtained from an LDLT decomposition and add
85 //the mean values.
86 Eigen::VectorXd x(m_mean.rows());
87 for (int i = 0; i < m_mean.rows(); ++i) x(i) = gRandom->Gaus();
88 return m_mean + (m_transform * x);
89 }
Belle2::MultivariateNormalGenerator::m_mean
Eigen::VectorXd m_mean
Member to store the mean values of the distribution.
Definition: MultivariateNormalGenerator.h:170
Belle2::MultivariateNormalGenerator::m_transform
Eigen::MatrixXd m_transform
Member to store the transformation matrix for standard normal distributed random values.
Definition: MultivariateNormalGenerator.h:173

◆ generate() [2/2]

void generate ( double *  output) const
inline

Generate a set of correlated random numbers with the previouly set mean and covariance and store them in buffer output.

Parameters
outputpointer to array where generated values will be stored.

Definition at line 100 of file MultivariateNormalGenerator.h.

101 {
102 Eigen::VectorXd x = generate();
103 for (int i = 0; i < x.rows(); ++i) { output[i] = x(i); }
104 }
Belle2::MultivariateNormalGenerator::generate
Eigen::VectorXd generate() const
Generate a set of correlated random numbers with the previouly set mean and covariance.
Definition: MultivariateNormalGenerator.h:80

◆ generateVec3()

TVector3 generateVec3 ( ) const
inline

Generate a set of correlated random numbers with the previouly set mean and covariance and return a TVector3.

Optimally, the set mean and covariance matrix should be of dimension three, otherwise just the first size() elements of the TVector3 are set and the remaining elements are zero. If size() is bigger than 3 the remaining values will be discarded.

Definition at line 112 of file MultivariateNormalGenerator.h.

113 {
114 Eigen::VectorXd x = generate();
115 TVector3 output(0, 0, 0);
116 for (unsigned int i = 0; i < std::min(3u, (unsigned int)size()); ++i) {
117 output[i] = x(i);
118 }
119 return output;
120 }
Belle2::MultivariateNormalGenerator::size
size_t size() const
Return the number of elements to be generated on generate()
Definition: MultivariateNormalGenerator.h:95

◆ generateVecT()

TVectorD generateVecT ( ) const
inline

Generate a set of correlated random numbers with the previouly set mean and covariance and return a TVectorT<double>

Definition at line 125 of file MultivariateNormalGenerator.h.

126 {
127 Eigen::VectorXd x = generate();
128 TVectorD output(x.rows());
129 output.SetElements(x.data());
130 return output;
131 }

◆ reset()

void reset ( )

reset the generator setting the size to 0.

Subsequent calls to generate will return 0-sized results until the generator is reinitialized using setMeanCov()

Definition at line 14 of file MultivariateNormalGenerator.cc.

15{
16 m_mean.resize(0);
17 m_transform.resize(0, 0);
18}

◆ setMeanCov() [1/2]

bool setMeanCov ( const Eigen::VectorXd &  mean,
const Eigen::MatrixXd &  cov 
)

set the mean and covariance for the distribution.

Parameters
meanVector of mean values
covMatrix containing the covariance values
Returns
true if covariance could be decomposited, false otherwise

Definition at line 36 of file MultivariateNormalGenerator.cc.

37{
38 reset();
39 if (mean.rows() != cov.rows()) {
40 B2ERROR("Mean values and covariance matrix need to be of the same dimension");
41 return false;
42 }
43 if (cov.rows() != cov.cols()) {
44 B2ERROR("Covariance matrix needs to be a square matrix");
45 return false;
46 }
47
48 // The usual way to calulate multivariate normal distributed values is to use
49 // the cholesky decomposition C = ll' and then calculate y = mean + l * x
50 // where x is a vector of standard normal distributed values. However this
51 // only works for positive definite matrices so we use the LDLT composition
52 // which gives us C = P'LDL'P and we compute L=P'Lsqrt(D) to get the same
53 // result but in a more robust way which works also for semi-definite
54 // matrices
55 auto ldlt = cov.ldlt();
56 if (ldlt.info() != Eigen::Success) {
57 B2ERROR("Cannot compute LDLT decomposition of covariance "
58 "matrix, maybe not positive semi-definite?");
59 return false;
60 }
61 Eigen::MatrixXd L = ldlt.matrixL();
62 Eigen::MatrixXd D = ldlt.vectorD().asDiagonal();
63 if (D.minCoeff() < 0) {
64 B2ERROR("MultivariateNormalGenerator: Negative values when computing LDL^T "
65 "decomposition, cannot compute M=AA^T, resulting random numbers "
66 "will not be correct");
67 return false;
68 }
69 auto P = ldlt.transpositionsP().transpose();
70 m_transform = P * L * D.cwiseSqrt();
71 //And save the mean values
72 m_mean = mean;
73 return true;
74}
Belle2::MultivariateNormalGenerator::reset
void reset()
reset the generator setting the size to 0.
Definition: MultivariateNormalGenerator.cc:14

◆ setMeanCov() [2/2]

bool setMeanCov ( int  n,
const double *  mean,
const double *  cov 
)

set the mean and covariance for the distribution with array interface: mean and covariance are passed as double arrays where the covariance is expected to be an NxN matrix in row major layout

Parameters
ndimensionality
meanpointer to the n mean values of the distribution
covpointer to the n*n covariance values in row major layout
Returns
true if covariance could be decomposited, false otherwise

Definition at line 20 of file MultivariateNormalGenerator.cc.

21{
22 Eigen::VectorXd emean(n);
23 Eigen::MatrixXd ecov(n, n);
24 //Fill the vector and matrix from the buffers in row major layout
25 for (int i = 0; i < n; ++i) {
26 emean[i] = mean[i];
27 for (int j = 0; j < n; ++j) {
28 ecov(i, j) = cov[i * n + j];
29 }
30 }
31 //And delegate to the correct function
32 return setMeanCov(emean, ecov);
33}

◆ size()

size_t size ( ) const
inline

Return the number of elements to be generated on generate()

Definition at line 95 of file MultivariateNormalGenerator.h.

95{ return m_mean.rows(); }

Member Data Documentation

◆ m_mean

Eigen::VectorXd m_mean
private

Member to store the mean values of the distribution.

Definition at line 170 of file MultivariateNormalGenerator.h.

◆ m_transform

Eigen::MatrixXd m_transform
private

Member to store the transformation matrix for standard normal distributed random values.

Definition at line 173 of file MultivariateNormalGenerator.h.

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