Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
probas.cxx
Go to the documentation of this file.
1#include "TCanvas.h"
2#include "TFile.h"
3#include "TH2F.h"
4#include "TIterator.h"
5#include "TKey.h"
6#include "TLegend.h"
7#include "TList.h"
8#include "TMVA/probas.h"
9#include "TMVA/tmvaglob.h"
10#include "TString.h"
11
12#include <iostream>
13
14using std::cout;
15using std::endl;
16
17// this macro plots the MVA probability distributions (Signal and
18// Background overlayed) of different MVA methods run in TMVA
19// (e.g. running TMVAnalysis.C).
20
21// input: - Input file (result from TMVA)
22// - use of TMVA plotting TStyle
23void TMVA::probas(TString dataset, TString fin , Bool_t useTMVAStyle )
24{
25 // set style and remove existing canvas'
26 TMVAGlob::Initialize( useTMVAStyle );
27
28 // switches
29 const Bool_t Draw_CFANN_Logy = kFALSE;
30 const Bool_t Save_Images = kTRUE;
31
32 // checks if file with name "fin" is already open, and if not opens one
33 TFile* file = TMVAGlob::OpenFile( fin );
34
35 const Int_t width = 600; // size of canvas
36
37 // this defines how many canvases we need
38 TCanvas *c = 0;
39
40 // counter variables
41 Int_t countCanvas = 0;
42
43 // list of existing MVAs
44 //const Int_t nveto = 1;
45 TString suffixSig = "_tr_S";
46 TString suffixBgd = "_tr_B";
47
48 // search for the right histograms in full list of keys
49 TList methods;
50 UInt_t nmethods = TMVAGlob::GetListOfMethods( methods,file->GetDirectory(dataset.Data()) );
51 if (nmethods==0) {
52 cout << "--- Probas.C: no methods found!" << endl;
53 return;
54 }
55 TIter next(&methods);
56 TKey *key, *hkey;
57 char fname[200];
58 TH1* sig(0);
59 TH1* bgd(0);
60
61
62 while ( (key = (TKey*)next()) ) {
63 TDirectory * mDir = (TDirectory*)key->ReadObj();
64 TList titles;
65 UInt_t ni = TMVAGlob::GetListOfTitles( mDir, titles );
66 TString methodName;
67 TMVAGlob::GetMethodName(methodName,key);
68 if (ni==0) {
69 cout << "+++ No titles found for classifier: " << methodName << endl;
70 return;
71 }
72 TIter nextTitle(&titles);
73 TKey *instkey;
74 TDirectory *instDir;
75
76 // iterate over all classifiers
77 while ( (instkey = (TKey *)nextTitle()) ) {
78 instDir = (TDirectory *)instkey->ReadObj();
79 TString instName = instkey->GetName();
80 TList h1hists;
81 UInt_t nhists = TMVAGlob::GetListOfKeys( h1hists, "TH1", instDir );
82 if (nhists==0) cout << "*** No histograms found!" << endl;
83 TIter nextInDir(&h1hists);
84 TString methodTitle;
85 TMVAGlob::GetMethodTitle(methodTitle,instDir);
86 Bool_t found = kFALSE;
87 while ( (hkey = (TKey*)nextInDir()) ) {
88 TH1 *th1 = (TH1*)hkey->ReadObj();
89 TString hname= th1->GetName();
90 if (hname.Contains( suffixSig ) && !hname.Contains( "Cut") &&
91 !hname.Contains("original") && !hname.Contains("smoothed")) {
92 // retrieve corresponding signal and background histograms
93 TString hnameS = hname;
94 TString hnameB = hname; hnameB.ReplaceAll("_S","_B");
95
96 sig = (TH1*)instDir->Get( hnameS );
97 bgd = (TH1*)instDir->Get( hnameB );
98
99 if (sig == 0 || bgd == 0) {
100 cout << "*** probas.C: big troubles in probas.... histogram: " << hname << " not found" << endl;
101 return;
102 }
103
104 TH1* sigF(0);
105 TH1* bkgF(0);
106
107 for (int i=0; i<= 5; i++) {
108 TString hspline = hnameS + Form("_smoothed_hist_from_spline%i",i);
109 sigF = (TH1*)instDir->Get( hspline );
110
111 if (sigF) {
112 bkgF = (TH1*)instDir->Get( hspline.ReplaceAll("_tr_S","_tr_B") );
113 break;
114 }
115 }
116 if (!sigF){
117 TString hspline = hnameS + TString("_smoothed_hist_from_KDE");
118 sigF = (TH1*)instDir->Get( hspline );
119
120 if (sigF) {
121 bkgF = (TH1*)instDir->Get( hspline.ReplaceAll("_tr_S","_tr_B") );
122 }
123 }
124
125 if ((sigF == NULL || bkgF == NULL) &&!hname.Contains("hist") ) {
126 cout << "*** probas.C: big troubles - did not find probability histograms" << endl;
127 return;
128 }
129 else {
130 // remove the signal suffix
131
132 // check that exist
133 if (NULL != sigF && NULL != bkgF && NULL!=sig && NULL!=bgd) {
134
135 found = kTRUE;
136 // chop off useless stuff
137 sig->SetTitle( TString("TMVA output for classifier: ") + methodTitle );
138
139 // create new canvas
140 cout << "--- Book canvas no: " << countCanvas << endl;
141 char cn[20];
142 sprintf( cn, "canvas%d", countCanvas+1 );
143 c = new TCanvas( cn, Form("TMVA Output Fit Variables %s",methodTitle.Data()),
144 countCanvas*50+200, countCanvas*20, width, width*0.78 );
145
146 // set the histogram style
147 TMVAGlob::SetSignalAndBackgroundStyle( sig, bgd );
148 TMVAGlob::SetSignalAndBackgroundStyle( sigF, bkgF );
149
150 // frame limits (choose judicuous x range)
151 Float_t nrms = 4;
152 Float_t xmin = TMath::Max( TMath::Min(sig->GetMean() - nrms*sig->GetRMS(),
153 bgd->GetMean() - nrms*bgd->GetRMS() ),
154 sig->GetXaxis()->GetXmin() );
155 Float_t xmax = TMath::Min( TMath::Max(sig->GetMean() + nrms*sig->GetRMS(),
156 bgd->GetMean() + nrms*bgd->GetRMS() ),
157 sig->GetXaxis()->GetXmax() );
158 Float_t ymin = 0;
159 Float_t ymax = TMath::Max( sig->GetMaximum(), bgd->GetMaximum() )*1.5;
160
161 if (Draw_CFANN_Logy && methodName == "CFANN") ymin = 0.01;
162
163 // build a frame
164 Int_t nb = 500;
165 TH2F* frame = new TH2F( TString("frame") + sig->GetName() + "_proba", sig->GetTitle(),
166 nb, xmin, xmax, nb, ymin, ymax );
167 frame->GetXaxis()->SetTitle(methodTitle);
168 frame->GetYaxis()->SetTitle("Normalized");
169 TMVAGlob::SetFrameStyle( frame );
170
171 // eventually: draw the frame
172 frame->Draw();
173
174 if (Draw_CFANN_Logy && methodName == "CFANN") c->SetLogy();
175
176 // overlay signal and background histograms
177 sig->SetMarkerColor( TMVAGlob::getSignalLine() );
178 sig->SetMarkerSize( 0.7 );
179 sig->SetMarkerStyle( 20 );
180 sig->SetLineWidth(1);
181
182 bgd->SetMarkerColor( TMVAGlob::getBackgroundLine() );
183 bgd->SetMarkerSize( 0.7 );
184 bgd->SetMarkerStyle( 24 );
185 bgd->SetLineWidth(1);
186
187 sig->Draw("samee");
188 bgd->Draw("samee");
189
190 sigF->SetFillStyle( 0 );
191 bkgF->SetFillStyle( 0 );
192 sigF->Draw("samehist");
193 bkgF->Draw("samehist");
194
195 // redraw axes
196 frame->Draw("sameaxis");
197
198 // Draw legend
199 TLegend *legend= new TLegend( c->GetLeftMargin(), 1 - c->GetTopMargin() - 0.2,
200 c->GetLeftMargin() + 0.4, 1 - c->GetTopMargin() );
201 legend->AddEntry(sig,"Signal data","P");
202 legend->AddEntry(sigF,"Signal PDF","L");
203 legend->AddEntry(bgd,"Background data","P");
204 legend->AddEntry(bkgF,"Background PDF","L");
205 legend->Draw("same");
206 legend->SetBorderSize(1);
207 legend->SetMargin( 0.3 );
208
209 // save canvas to file
210 c->Update();
211 TMVAGlob::plot_logo();
212 sprintf( fname, "%s/plots/mva_pdf_%s_c%i",dataset.Data(), methodTitle.Data(), countCanvas+1 );
213 if (Save_Images) TMVAGlob::imgconv( c, fname );
214 countCanvas++;
215 }
216 }
217 }
218
219 }
220 if(!found){
221 cout << "--- No PDFs found for method " << methodTitle << ". Did you request \"CreateMVAPdfs\" in the option string?" << endl;
222 }
223 }
224 }
225}
#define c(i)
Definition RSha256.hxx:101
const Bool_t kFALSE
Definition RtypesCore.h:92
float Float_t
Definition RtypesCore.h:57
const Bool_t kTRUE
Definition RtypesCore.h:91
include TDocParser_001 C image html pict1_TDocParser_001 png width
float xmin
float ymin
float xmax
float ymax
char * Form(const char *fmt,...)
The Canvas class.
Definition TCanvas.h:23
Describe directory structure in memory.
Definition TDirectory.h:45
virtual TObject * Get(const char *namecycle)
Return pointer to object identified by namecycle.
A ROOT file is a suite of consecutive data records (TKey instances) with a well defined format.
Definition TFile.h:54
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:58
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
Definition TH1.h:320
TAxis * GetYaxis()
Definition TH1.h:321
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition TH1.cxx:3073
2-D histogram with a float per channel (see TH1 documentation)}
Definition TH2.h:251
Book space in a file, create I/O buffers, to fill them, (un)compress them.
Definition TKey.h:28
virtual TObject * ReadObj()
To read a TObject* from the file.
Definition TKey.cxx:750
This class displays a legend box (TPaveText) containing several legend entries.
Definition TLegend.h:23
TLegendEntry * AddEntry(const TObject *obj, const char *label="", Option_t *option="lpf")
Add a new entry to this legend.
Definition TLegend.cxx:330
virtual void Draw(Option_t *option="")
Draw this legend with its current attributes.
Definition TLegend.cxx:423
void SetMargin(Float_t margin)
Definition TLegend.h:69
A doubly linked list.
Definition TList.h:44
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition TNamed.cxx:164
virtual const char * GetName() const
Returns name of object.
Definition TNamed.h:47
virtual void SetBorderSize(Int_t bordersize=4)
Definition TPave.h:73
Basic string class.
Definition TString.h:136
const char * Data() const
Definition TString.h:369
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition TString.h:692
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:624
void probas(TString dataset, TString fin="TMVA.root", Bool_t useTMVAStyle=kTRUE)
Short_t Max(Short_t a, Short_t b)
Definition TMathBase.h:212
Short_t Min(Short_t a, Short_t b)
Definition TMathBase.h:180
Definition file.py:1
auto * th1
Definition textalign.C:13