ROOT   Reference Guide
Loading...
Searching...
No Matches
MnCovarianceSqueeze.cxx
Go to the documentation of this file.
1// @(#)root/minuit2:$Id$
2// Authors: M. Winkler, F. James, L. Moneta, A. Zsenei 2003-2005
3
4/**********************************************************************
5 * *
6 * Copyright (c) 2005 LCG ROOT Math team, CERN/PH-SFT *
7 * *
8 **********************************************************************/
9
13#include "Minuit2/MnPrint.h"
14
15namespace ROOT {
16
17namespace Minuit2 {
18
20{
21 // squeeze MnUserCovariance class
22 // MnUserCovariance contasins the error matrix. Need to invert first to get the hessian, then
23 // after having squuezed the hessian, need to invert again to get the new error matrix
24 assert(cov.Nrow() > 0);
25 assert(n < cov.Nrow());
26
27 MnPrint print("MnCovarianceSqueeze");
28
29 MnAlgebraicSymMatrix hess(cov.Nrow());
30 for (unsigned int i = 0; i < cov.Nrow(); i++) {
31 for (unsigned int j = i; j < cov.Nrow(); j++) {
32 hess(i, j) = cov(i, j);
33 }
34 }
35
36 int ifail = Invert(hess);
37
38 if (ifail != 0) {
39 print.Warn("inversion failed; return diagonal matrix;");
40 MnUserCovariance result(cov.Nrow() - 1);
41 for (unsigned int i = 0, j = 0; i < cov.Nrow(); i++) {
42 if (i == n)
43 continue;
44 result(j, j) = cov(i, i);
45 j++;
46 }
47 return result;
48 }
49
50 MnAlgebraicSymMatrix squeezed = (*this)(hess, n);
51
52 ifail = Invert(squeezed);
53 if (ifail != 0) {
54 print.Warn("back-inversion failed; return diagonal matrix;");
55 MnUserCovariance result(squeezed.Nrow());
56 for (unsigned int i = 0; i < squeezed.Nrow(); i++) {
57 result(i, i) = 1. / squeezed(i, i);
58 }
59 return result;
60 }
61
62 return MnUserCovariance(std::vector<double>(squeezed.Data(), squeezed.Data() + squeezed.size()), squeezed.Nrow());
63}
64
66{
67
68 MnPrint print("MnCovarianceSqueeze");
69
70 // squueze the minimum error class
71 // Remove index-row on the Hessian matrix and the get the new correct error matrix
72 // (inverse of new Hessian)
73 MnAlgebraicSymMatrix hess = err.Hessian();
74 MnAlgebraicSymMatrix squeezed = (*this)(hess, n);
75 int ifail = Invert(squeezed);
76 if (ifail != 0) {
77 print.Warn("MinimumError inversion fails; return diagonal matrix.");
78
79 MnAlgebraicSymMatrix tmp(squeezed.Nrow());
80 for (unsigned int i = 0; i < squeezed.Nrow(); i++) {
81 tmp(i, i) = 1. / squeezed(i, i);
82 }
84 }
85
86 return MinimumError(squeezed, err.Dcovar());
87}
88
90{
91 // squueze a symmetrix matrix (remove entire row and column n)
92 assert(hess.Nrow() > 0);
93 assert(n < hess.Nrow());
94
95 MnAlgebraicSymMatrix hs(hess.Nrow() - 1);
96 for (unsigned int i = 0, j = 0; i < hess.Nrow(); i++) {
97 if (i == n)
98 continue;
99 for (unsigned int k = i, l = j; k < hess.Nrow(); k++) {
100 if (k == n)
101 continue;
102 hs(j, l) = hess(i, k);
103 l++;
104 }
105 j++;
106 }
107
108 return hs;
109}
110
111} // namespace Minuit2
112
113} // namespace ROOT
Class describing a symmetric matrix of size n.
Definition LASymMatrix.h:45
const double * Data() const
unsigned int Nrow() const
unsigned int size() const
MinimumError keeps the inv.
MnAlgebraicSymMatrix Hessian() const
MnUserCovariance operator()(const MnUserCovariance &, unsigned int) const
void Warn(const Ts &... args)
Definition MnPrint.h:126
Class containing the covariance matrix data represented as a vector of size n*(n+1)/2 Used to hide in...
const Int_t n
Definition legend1.C:16
int Invert(LASymMatrix &)
Definition LaInverse.cxx:21
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...
auto * l
Definition textangle.C:4