 ROOT   Reference Guide 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
12 #include "Minuit2/MinimumError.h"
13 #include "Minuit2/MnPrint.h"
14
15 namespace ROOT {
16
17 namespace 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
l
auto * l
Definition: textangle.C:4
n
const Int_t n
Definition: legend1.C:16
ROOT::Minuit2::MinimumError::Hessian
MnAlgebraicSymMatrix Hessian() const
Definition: MinimumError.h:56
ROOT::Minuit2::Invert
int Invert(LASymMatrix &)
Definition: LaInverse.cxx:21
ROOT::Minuit2::LASymMatrix::Data
const double * Data() const
Definition: LASymMatrix.h:268
ROOT::Minuit2::LASymMatrix
Class describing a symmetric matrix of size n.
Definition: LASymMatrix.h:45
ROOT::Minuit2::LASymMatrix::size
unsigned int size() const
Definition: LASymMatrix.h:272
ROOT::Minuit2::MnUserCovariance
Class containing the covariance matrix data represented as a vector of size n*(n+1)/2 Used to hide in...
Definition: MnUserCovariance.h:26
ROOT::Minuit2::MnUserCovariance::Nrow
unsigned int Nrow() const
Definition: MnUserCovariance.h:84
ROOT::Minuit2::MnPrint::Warn
void Warn(const Ts &... args)
Definition: MnPrint.h:126
ROOT::Minuit2::MinimumError::Dcovar
double Dcovar() const
Definition: MinimumError.h:70
MnUserCovariance.h
ROOT::Minuit2::LASymMatrix::Nrow
unsigned int Nrow() const
Definition: LASymMatrix.h:274
ROOT::Minuit2::MinimumError::MnInvertFailed
@ MnInvertFailed
Definition: MinimumError.h:35
MnCovarianceSqueeze.h
ROOT::Minuit2::MinimumError
MinimumError keeps the inv.
Definition: MinimumError.h:28
MinimumError.h
ROOT::Minuit2::MnCovarianceSqueeze::operator()
MnUserCovariance operator()(const MnUserCovariance &, unsigned int) const
Definition: MnCovarianceSqueeze.cxx:19
MnPrint.h
ROOT
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...
Definition: EExecutionPolicy.hxx:4
ROOT::Minuit2::MnPrint
Definition: MnPrint.h:73