Logo ROOT   6.08/07
Reference Guide
RootFinder.cxx
Go to the documentation of this file.
1 // @(#)root/tmva $Id$
2 // Author: Andreas Hoecker, Joerg Stelzer, Helge Voss
3 
4 /**********************************************************************************
5  * Project: TMVA - a Root-integrated toolkit for multivariate data analysis *
6  * Package: TMVA *
7  * Class : RootFinder *
8  * Web : http://tmva.sourceforge.net *
9  * *
10  * Description: *
11  * Implementation (see header file for description) *
12  * *
13  * Authors (alphabetical): *
14  * Andreas Hoecker <Andreas.Hocker@cern.ch> - CERN, Switzerland *
15  * Helge Voss <Helge.Voss@cern.ch> - MPI-K Heidelberg, Germany *
16  * Kai Voss <Kai.Voss@cern.ch> - U. of Victoria, Canada *
17  * *
18  * Copyright (c) 2005: *
19  * CERN, Switzerland *
20  * U. of Victoria, Canada *
21  * MPI-K Heidelberg, Germany *
22  * *
23  * Redistribution and use in source and binary forms, with or without *
24  * modification, are permitted according to the terms listed in LICENSE *
25  * (http://tmva.sourceforge.net/LICENSE) *
26  **********************************************************************************/
27 
28 #include "TMVA/RootFinder.h"
29 
30 #include "TMVA/MethodBase.h"
31 #include "TMVA/MsgLogger.h"
32 #include "TMVA/Types.h"
33 
34 #include "TMath.h"
35 #include "TObject.h"
36 
38 
39 ////////////////////////////////////////////////////////////////////////////////
40 /// constructor
41 
42 
44  Double_t rootMin,
45  Double_t rootMax,
46  Int_t maxIterations,
47  Double_t absTolerance )
48 : fRootMin( rootMin ),
49  fRootMax( rootMax ),
50  fMaxIter( maxIterations ),
51  fAbsTol ( absTolerance ),
52  fLogger ( new MsgLogger("RootFinder") )
53 {
54  fMethod=method;
55 }
56 
57 ////////////////////////////////////////////////////////////////////////////////
58 /// destructor
59 
61 {
62  delete fLogger;
63 }
64 
65 ////////////////////////////////////////////////////////////////////////////////
66 /// Root finding using Brents algorithm; taken from CERNLIB function RZERO
67 
69 {
71  Double_t fa = fMethod->GetValueForRoot( a ) - refValue;
72  Double_t fb = fMethod->GetValueForRoot( b ) - refValue;
73  if (fb*fa > 0) {
74  Log() << kWARNING << "<Root> initial interval w/o root: "
75  << "(a=" << a << ", b=" << b << "),"
76  << " (Eff_a=" << fMethod->GetValueForRoot( a )
77  << ", Eff_b=" << fMethod->GetValueForRoot( b ) << "), "
78  << "(fa=" << fa << ", fb=" << fb << "), "
79  << "refValue = " << refValue << Endl;
80  return 1;
81  }
82 
83  Bool_t ac_equal(kFALSE);
84  Double_t fc = fb;
85  Double_t c = 0, d = 0, e = 0;
86  for (Int_t iter= 0; iter <= fMaxIter; iter++) {
87  if ((fb < 0 && fc < 0) || (fb > 0 && fc > 0)) {
88 
89  // Rename a,b,c and adjust bounding interval d
90  ac_equal = kTRUE;
91  c = a; fc = fa;
92  d = b - a; e = b - a;
93  }
94 
95  if (TMath::Abs(fc) < TMath::Abs(fb)) {
96  ac_equal = kTRUE;
97  a = b; b = c; c = a;
98  fa = fb; fb = fc; fc = fa;
99  }
100 
101  Double_t tol = 0.5 * 2.2204460492503131e-16 * TMath::Abs(b);
102  Double_t m = 0.5 * (c - b);
103  if (fb == 0 || TMath::Abs(m) <= tol || TMath::Abs(fb) < fAbsTol) return b;
104 
105  // Bounds decreasing too slowly: use bisection
106  if (TMath::Abs (e) < tol || TMath::Abs (fa) <= TMath::Abs (fb)) { d = m; e = m; }
107  else {
108  // Attempt inverse cubic interpolation
109  Double_t p, q, r;
110  Double_t s = fb / fa;
111 
112  if (ac_equal) { p = 2 * m * s; q = 1 - s; }
113  else {
114  q = fa / fc; r = fb / fc;
115  p = s * (2 * m * q * (q - r) - (b - a) * (r - 1));
116  q = (q - 1) * (r - 1) * (s - 1);
117  }
118  // Check whether we are in bounds
119  if (p > 0) q = -q;
120  else p = -p;
121 
122  Double_t min1 = 3 * m * q - TMath::Abs (tol * q);
123  Double_t min2 = TMath::Abs (e * q);
124  if (2 * p < (min1 < min2 ? min1 : min2)) {
125  // Accept the interpolation
126  e = d; d = p / q;
127  }
128  else { d = m; e = m; } // Interpolation failed: use bisection.
129  }
130  // Move last best guess to a
131  a = b; fa = fb;
132  // Evaluate new trial root
133  if (TMath::Abs(d) > tol) b += d;
134  else b += (m > 0 ? +tol : -tol);
135 
136  fb = fMethod->GetValueForRoot( b ) - refValue;
137 
138  }
139 
140  // Return our best guess if we run out of iterations
141  Log() << kWARNING << "<Root> maximum iterations (" << fMaxIter
142  << ") reached before convergence" << Endl;
143 
144  return b;
145 }
146 
MsgLogger & Endl(MsgLogger &ml)
Definition: MsgLogger.h:162
return c
virtual Double_t GetValueForRoot(Double_t)
returns efficiency as function of cut
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
TArc * a
Definition: textangle.C:12
const Bool_t kFALSE
Definition: Rtypes.h:92
Double_t fAbsTol
Definition: RootFinder.h:68
Short_t Abs(Short_t d)
Definition: TMathBase.h:110
MethodBase * fMethod
Definition: RootFinder.h:71
static struct mg_connection * fc(struct mg_context *ctx)
Definition: civetweb.c:1956
Double_t fRootMax
Definition: RootFinder.h:66
MsgLogger * fLogger
Definition: RootFinder.h:73
Double_t Root(Double_t refValue)
Root finding using Brents algorithm; taken from CERNLIB function RZERO.
Definition: RootFinder.cxx:68
const double tol
Double_t fRootMin
Definition: RootFinder.h:65
TRandom2 r(17)
TMarker * m
Definition: textangle.C:8
MsgLogger & Log() const
message logger
Definition: RootFinder.h:74
#define ClassImp(name)
Definition: Rtypes.h:279
double Double_t
Definition: RtypesCore.h:55
virtual ~RootFinder(void)
destructor
Definition: RootFinder.cxx:60
you should not use this method at all Int_t Int_t Double_t Double_t Double_t e
Definition: TRolke.cxx:630
Abstract ClassifierFactory template that handles arbitrary types.
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t b
Definition: TRolke.cxx:630
const Bool_t kTRUE
Definition: Rtypes.h:91
float * q
Definition: THbookFile.cxx:87