Logo ROOT   6.30.04
Reference Guide
 All Namespaces Files Pages
MnPrint.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 
10 #include "Minuit2/MnPrint.h"
11 #include "Minuit2/LAVector.h"
12 #include "Minuit2/LASymMatrix.h"
13 #include "Minuit2/FunctionMinimum.h"
14 #include "Minuit2/MnUserParameters.h"
15 #include "Minuit2/MnUserCovariance.h"
16 #include "Minuit2/MnGlobalCorrelationCoeff.h"
17 #include "Minuit2/MnUserParameterState.h"
18 #include "Minuit2/MinuitParameter.h"
19 #include "Minuit2/MnMachinePrecision.h"
20 #include "Minuit2/MinosError.h"
21 #include "Minuit2/ContoursError.h"
22 #include "Minuit2/MnPlot.h"
23 
24 #include <iomanip>
25 
26 #define PRECISION 13
27 #define WIDTH 20
28 
29 
30 namespace ROOT {
31 
32  namespace Minuit2 {
33 
34 #ifdef DEBUG
35 int gPrintLevel = 3;
36 #else
37 int gPrintLevel = 0;
38 #endif
39 
40 
41 int MnPrint::SetLevel(int level) {
42  int prevLevel = gPrintLevel;
43  gPrintLevel = level;
44  return prevLevel;
45 }
46 
47 int MnPrint::Level( ) {
48  return gPrintLevel;
49 }
50 
51 void MnPrint::PrintFcn(std::ostream & os, double value, bool endline) {
52  int pr = os.precision(PRECISION);
53  os << value;
54  if (endline) os << std::endl;
55  os.precision(pr);
56 }
57 
58 void MnPrint::PrintState(std::ostream & os, const MinimumState & state, const char * msg, int iter) {
59  // helper function to print function value, edm and ncalls from state in one single line
60  MnPrint::PrintState(os, state.Fval(), state.Edm(), state.NFcn(), msg, iter);
61 }
62 
63 void MnPrint::PrintState(std::ostream & os, double fval, double edm, int ncalls, const char * msg, int iter) {
64  // helper function to print function value, edm and ncalls and message in one single line
65  os << msg;
66  if (iter>=0) os << std::setw(3) << iter;
67  int pr = os.precision(PRECISION);
68  const int width = PRECISION+3;
69  os << " - FCN = " << std::setw(width) << fval;
70  os.precision(pr);
71  os << " Edm = " << std::setw(12) << edm << " NCalls = " << std::setw(6) << ncalls;
72  os << std::endl;
73 }
74 
75 
76 std::ostream& operator<<(std::ostream& os, const LAVector& vec) {
77  // print a vector
78  os << "LAVector parameters:" << std::endl;
79  int pr = os.precision(PRECISION);
80  {
81  int nrow = vec.size();
82  for (int i = 0; i < nrow; i++) {
83  os.width(WIDTH);
84  os << vec(i) << std::endl;
85  }
86  }
87  os.precision(pr);
88  return os;
89 }
90 
91 std::ostream& operator<<(std::ostream& os, const LASymMatrix& matrix) {
92  // print a matrix
93  os << "LASymMatrix parameters:" << std::endl;
94  int pr = os.precision(8);
95  {
96  //os << std::endl;
97  int n = matrix.Nrow();
98  for (int i = 0; i < n; i++) {
99  for (int j = 0; j < n; j++) {
100  os.width(15); os << matrix(i,j);
101  }
102  os << std::endl;
103  }
104  }
105  os.precision(pr);
106  return os;
107 }
108 
109 std::ostream& operator<<(std::ostream& os, const MnUserParameters& par) {
110  // print the MnUserParameter object
111  os << std::endl;
112 
113  os << "# ext. |" << "| Name |" << "| type |" << "| Value |" << "| Error +/- " << std::endl;
114 
115  os << std::endl;
116  int pr = os.precision();
117 
118  bool atLoLim = false;
119  bool atHiLim = false;
120  for(std::vector<MinuitParameter>::const_iterator ipar = par.Parameters().begin(); ipar != par.Parameters().end(); ++ipar) {
121  os << std::setw(4) << (*ipar).Number() << std::setw(5) << "||";
122  os << std::setw(10) << (*ipar).Name() << std::setw(3) << "||";
123  if((*ipar).IsConst()) {
124  os << " const ||" << std::setprecision(PRECISION) << std::setw(WIDTH) << (*ipar).Value() << " ||" << std::endl;
125  } else if((*ipar).IsFixed()) {
126  os << " fixed ||" << std::setprecision(PRECISION) << std::setw(WIDTH) << (*ipar).Value() << " ||" << std::endl;
127  } else if((*ipar).HasLimits()) {
128  if((*ipar).Error() > 0.) {
129  os << " limited ||" << std::setprecision(PRECISION) << std::setw(WIDTH) << (*ipar).Value();
130  if(fabs((*ipar).Value() - (*ipar).LowerLimit()) < par.Precision().Eps2()) {
131  os <<"*";
132  atLoLim = true;
133  }
134  if(fabs((*ipar).Value() - (*ipar).UpperLimit()) < par.Precision().Eps2()) {
135  os <<"**";
136  atHiLim = true;
137  }
138  os << " ||" << std::setw(12) << (*ipar).Error() << std::endl;
139  } else
140  os << " free ||" << std::setprecision(PRECISION) << std::setw(WIDTH) << (*ipar).Value() << " ||" << std::setw(12) << "no" << std::endl;
141  } else {
142  if((*ipar).Error() > 0.)
143  os << " free ||" << std::setprecision(PRECISION) << std::setw(WIDTH) << (*ipar).Value() << " ||" << std::setw(12) << (*ipar).Error() << std::endl;
144  else
145  os << " free ||" << std::setprecision(PRECISION) << std::setw(WIDTH) << (*ipar).Value() << " ||" << std::setw(12) << "no" << std::endl;
146 
147  }
148  }
149  os << std::endl;
150  if(atLoLim) os << "* Parameter is at Lower limit" << std::endl;
151  if(atHiLim) os << "** Parameter is at Upper limit" << std::endl;
152  os << std::endl;
153  os.precision(pr);
154 
155  return os;
156 }
157 
158 std::ostream& operator<<(std::ostream& os, const MnUserCovariance& matrix) {
159  // print the MnUserCovariance
160  os << std::endl;
161 
162  os << "MnUserCovariance: " << std::endl;
163  int pr = os.precision(6);
164  {
165  os << std::endl;
166  unsigned int n = matrix.Nrow();
167  for (unsigned int i = 0; i < n; i++) {
168  for (unsigned int j = 0; j < n; j++) {
169  os.width(13); os << matrix(i,j);
170  }
171  os << std::endl;
172  }
173  }
174 
175  os << std::endl;
176  os << "MnUserCovariance Parameter correlations: " << std::endl;
177 
178  {
179  os << std::endl;
180  unsigned int n = matrix.Nrow();
181  for (unsigned int i = 0; i < n; i++) {
182  double di = matrix(i,i);
183  for (unsigned int j = 0; j < n; j++) {
184  double dj = matrix(j,j);
185  os.width(13); os << matrix(i,j)/sqrt(fabs(di*dj));
186  }
187  os << std::endl;
188  }
189  }
190 
191  os.precision(pr);
192  return os;
193 }
194 
195 std::ostream& operator<<(std::ostream& os, const MnGlobalCorrelationCoeff& coeff) {
196  // print the global correlation coefficient
197  os << std::endl;
198 
199  os << "MnGlobalCorrelationCoeff: " << std::endl;
200  int pr = os.precision(6);
201  {
202  os << std::endl;
203  for (unsigned int i = 0; i < coeff.GlobalCC().size(); i++) {
204  os.width(13); os << coeff.GlobalCC()[i];
205  os << std::endl;
206  }
207  }
208 
209  os.precision(pr);
210  return os;
211 }
212 
213 std::ostream& operator<<(std::ostream& os, const MnUserParameterState& state) {
214  // print the MnUserParameterState
215  os << std::endl;
216 
217  if(!state.IsValid()) {
218  os << std::endl;
219  os <<"WARNING: MnUserParameterState is not valid."<<std::endl;
220  os << std::endl;
221  }
222  int pr = os.precision(PRECISION);
223 
224  os <<"# of function calls: "<<state.NFcn()<<std::endl;
225  os <<"function Value: " << state.Fval()<<std::endl;
226  os <<"expected distance to the Minimum (edm): " << state.Edm()<<std::endl;
227  os <<"external parameters: "<<state.Parameters()<<std::endl;
228  os <<"covariance matrix status: " << state.CovarianceStatus() << std::endl;
229  if(state.HasCovariance())
230  os <<"covariance matrix: "<<state.Covariance()<<std::endl;
231  else
232  os <<"covariance matrix is not present or not valid "<<std::endl;
233  if(state.HasGlobalCC())
234  os <<"global correlation coefficients : "<<state.GlobalCC()<<std::endl;
235 
236  if(!state.IsValid())
237  os <<"WARNING: MnUserParameterState is not valid."<<std::endl;
238 
239  os << std::endl;
240  os.precision(pr);
241 
242  return os;
243 }
244 
245 std::ostream& operator<<(std::ostream& os, const FunctionMinimum& min) {
246  // print the FunctionMinimum
247  os << std::endl;
248  if(!min.IsValid()) {
249  os <<"WARNING: Minuit did not converge."<<std::endl;
250  os << std::endl;
251  } else {
252  os <<"Minuit did successfully converge."<<std::endl;
253  }
254  int pr = os.precision(PRECISION);
255 
256  os <<"# of function calls: "<<min.NFcn()<<std::endl;
257  os <<"minimum function Value: " << min.Fval()<<std::endl;
258  os <<"minimum edm: " << min.Edm()<<std::endl;
259  os <<"minimum internal state vector: "<<min.Parameters().Vec()<<std::endl;
260  if(min.HasValidCovariance())
261  os <<"minimum internal covariance matrix: "<<min.Error().Matrix()<<std::endl;
262 
263  os << min.UserParameters() << std::endl;
264  //os << min.UserCovariance() << std::endl;
265  //os << min.UserState().GlobalCC() << std::endl;
266 
267  if(!min.IsValid()) {
268  os <<"WARNING: FunctionMinimum is invalid: " << std::endl;
269  if ( !min.State().IsValid() )
270  os << "\t State is invalid" << std::endl;
271  if ( min.IsAboveMaxEdm() )
272  os << "\t Edm is above max" << std::endl;
273  if ( min.HasReachedCallLimit() )
274  os << "\t Reached call limit" << std::endl;
275  }
276 
277 
278  os << std::endl;
279  os.precision(pr);
280 
281  return os;
282 }
283 
284 std::ostream& operator<<(std::ostream& os, const MinimumState& min) {
285 
286  os << std::endl;
287  int pr = os.precision(PRECISION);
288 
289  os <<"minimum function Value: " << min.Fval()<<std::endl;
290  os <<"minimum edm: " << min.Edm()<<std::endl;
291  os <<"minimum internal state vector: "<<min.Vec()<<std::endl;
292  os <<"minimum internal Gradient vector: "<<min.Gradient().Vec()<<std::endl;
293  if(min.HasCovariance())
294  os <<"minimum internal covariance matrix: "<<min.Error().Matrix()<<std::endl;
295 
296  os << std::endl;
297  os.precision(pr);
298 
299  return os;
300 }
301 
302 std::ostream& operator<<(std::ostream& os, const MnMachinePrecision& prec) {
303  // print the Precision
304  os << std::endl;
305 
306  int pr = os.precision(PRECISION);
307  os <<"current machine precision is set to "<<prec.Eps()<<std::endl;
308 
309  os << std::endl;
310  os.precision(pr);
311 
312  return os;
313 }
314 
315 std::ostream& operator<<(std::ostream& os, const MinosError& me) {
316  // print the Minos Error
317  os << std::endl;
318 
319  os <<"Minos # of function calls: "<<me.NFcn()<<std::endl;
320 
321  if(!me.IsValid())
322  os << "Minos Error is not valid." <<std::endl;
323  if(!me.LowerValid())
324  os << "lower Minos Error is not valid." <<std::endl;
325  if(!me.UpperValid())
326  os << "upper Minos Error is not valid." <<std::endl;
327  if(me.AtLowerLimit())
328  os << "Minos Error is Lower limit of Parameter "<<me.Parameter()<<"." <<std::endl;
329  if(me.AtUpperLimit())
330  os << "Minos Error is Upper limit of Parameter "<<me.Parameter()<<"." <<std::endl;
331  if(me.AtLowerMaxFcn())
332  os << "Minos number of function calls for Lower Error exhausted."<<std::endl;
333  if(me.AtUpperMaxFcn())
334  os << "Minos number of function calls for Upper Error exhausted."<<std::endl;
335  if(me.LowerNewMin()) {
336  os << "Minos found a new Minimum in negative direction."<<std::endl;
337  os << me.LowerState() <<std::endl;
338  }
339  if(me.UpperNewMin()) {
340  os << "Minos found a new Minimum in positive direction."<<std::endl;
341  os << me.UpperState() <<std::endl;
342  }
343 
344  int pr = os.precision();
345 
346  os << "# ext. |" << "| Name |" << "| Value@min |" << "| negative |" << "| positive " << std::endl;
347  os << std::setw(4) << me.Parameter() << std::setw(5) << "||";
348  os << std::setw(10) << me.LowerState().Name(me.Parameter()) << std::setw(3) << "||";
349  os << std::setprecision(PRECISION) << std::setw(WIDTH) << me.Min() << " ||" << std::setprecision(PRECISION) << std::setw(WIDTH) << me.Lower() << " ||" << std::setw(WIDTH) << me.Upper() << std::endl;
350 
351  os << std::endl;
352  os.precision(pr);
353 
354  return os;
355 }
356 
357 std::ostream& operator<<(std::ostream& os, const ContoursError& ce) {
358  // print the ContoursError
359  os << std::endl;
360  os <<"Contours # of function calls: "<<ce.NFcn()<<std::endl;
361  os << "MinosError in x: "<<std::endl;
362  os << ce.XMinosError() << std::endl;
363  os << "MinosError in y: "<<std::endl;
364  os << ce.YMinosError() << std::endl;
365  MnPlot plot;
366  plot(ce.XMin(), ce.YMin(), ce());
367  for(std::vector<std::pair<double,double> >::const_iterator ipar = ce().begin(); ipar != ce().end(); ++ipar) {
368  os << ipar - ce().begin() <<" "<< (*ipar).first <<" "<< (*ipar).second <<std::endl;
369  }
370  os << std::endl;
371 
372  return os;
373 }
374 
375  } // namespace Minuit2
376 
377 } // namespace ROOT