Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
h1analysisTreeReader.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_tree
3/// \notebook
4/// H1 analysis example expressed in terms of TTreeReader
5///
6/// \macro_code
7///
8/// \author Anders Eie, 2013
9
11#include "TStyle.h"
12#include "TCanvas.h"
13#include "TPaveStats.h"
14#include "TLine.h"
15#include "TMath.h"
16#include "TFile.h"
17#include "TROOT.h"
18
19
20const Double_t dxbin = (0.17-0.13)/40; // Bin-width
21const Double_t sigma = 0.0012;
22
23//_____________________________________________________________________
25{
26 Double_t x = xx[0];
27 if (x <= 0.13957) return 0;
28 Double_t xp3 = (x-par[3])*(x-par[3]);
29 Double_t res = dxbin*(par[0]*TMath::Power(x-0.13957, par[1])
30 + par[2] / 2.5066/par[4]*TMath::Exp(-xp3/2/par[4]/par[4]));
31 return res;
32}
33
34//_____________________________________________________________________
36{
37 Double_t x = xx[0];
38 if (x <= 0.13957) return 0;
39 Double_t xp3 = (x-0.1454)*(x-0.1454);
40 Double_t res = dxbin*(par[0]*TMath::Power(x-0.13957, 0.25)
41 + par[1] / 2.5066/sigma*TMath::Exp(-xp3/2/sigma/sigma));
42 return res;
43}
44//_____________________________________________________________________
46// entry is the entry number in the current Tree
47// Selection function to select D* and D0.
49 fProcessed++;
50 //in case one entry list is given in input, the selection has already been done.
51 if (!useList) {
52 // Return as soon as a bad entry is detected
53 if (TMath::Abs(*fMd0_d-1.8646) >= 0.04) return kFALSE;
54 if (*fPtds_d <= 2.5) return kFALSE;
55 if (TMath::Abs(*fEtads_d) >= 1.5) return kFALSE;
56 (*fIk)--; //original fIk used f77 convention starting at 1
57 (*fIpi)--;
58
59
60 if (fNhitrp.At(*fIk)* fNhitrp.At(*fIpi) <= 1) return kFALSE;
61
62
63 if (fRend.At(*fIk) -fRstart.At(*fIk) <= 22) return kFALSE;
64 if (fRend.At(*fIpi)-fRstart.At(*fIpi) <= 22) return kFALSE;
65 if (fNlhk.At(*fIk) <= 0.1) return kFALSE;
66 if (fNlhpi.At(*fIpi) <= 0.1) return kFALSE;
67 (*fIpis)--; if (fNlhpi.At(*fIpis) <= 0.1) return kFALSE;
68 if (*fNjets < 1) return kFALSE;
69 }
70 // if option fillList, fill the entry list
71 if (fillList) elist->Enter(entry);
72
73 //fill some histograms
74 hdmd->Fill(*fDm_d);
75 h2->Fill(*fDm_d,*fRpd0_t/0.029979*1.8646/ *fPtd0_d);
76
77 return kTRUE;
78}
79
80void h1analysisTreeReader::Begin(TTree* /*myTree*/) {
81// function called before starting the event loop
82// -it performs some cleanup
83// -it creates histograms
84// -it sets some initialisation for the entry list
85
86 Reset();
87
88 //print the option specified in the Process function.
89 TString option = GetOption();
90 Info("Begin", "starting h1analysis with process option: %s", option.Data());
91
92 delete gDirectory->GetList()->FindObject("elist");
93
94 // case when one creates/fills the entry list
95 if (option.Contains("fillList")) {
97 elist = new TEntryList("elist", "H1 selection from Cut");
98 // Add to the input list for processing in PROOF, if needed
99 if (fInput) {
100 fInput->Add(new TNamed("fillList",""));
101 // We send a clone to avoid double deletes when importing the result
102 fInput->Add(elist);
103 // This is needed to avoid warnings from output-to-members mapping
104 elist = 0;
105 }
106 Info("Begin", "creating an entry-list");
107 }
108 // case when one uses the entry list generated in a previous call
109 if (option.Contains("useList")) {
110 useList = kTRUE;
111 if (fInput) {
112 // In PROOF option "useList" is processed in SlaveBegin and we do not need
113 // to do anything here
114 } else {
115 TFile f("elist.root");
116 elist = (TEntryList*)f.Get("elist");
117 if (elist) elist->SetDirectory(0); //otherwise the file destructor will delete elist
118 }
119 }
120}
121
123
124// function called before starting the event loop
125// -it performs some cleanup
126// -it creates histograms
127// -it sets some initialisation for the entry list
128
129 Init(myTree);
130
131 //print the option specified in the Process function.
132 TString option = GetOption();
133 Info("SlaveBegin",
134 "starting h1analysis with process option: %s (tree: %p)", option.Data(), myTree);
135
136 //create histograms
137 hdmd = new TH1F("hdmd","Dm_d",40,0.13,0.17);
138 h2 = new TH2F("h2","ptD0 vs Dm_d",30,0.135,0.165,30,-3,6);
139
140 fOutput->Add(hdmd);
141 fOutput->Add(h2);
142
143 // Entry list stuff (re-parse option because on PROOF only SlaveBegin is called)
144 if (option.Contains("fillList")) {
145 fillList = kTRUE;
146 // Get the list
147 if (fInput) {
148 if ((elist = (TEntryList *) fInput->FindObject("elist")))
149 // Need to clone to avoid problems when destroying the selector
150 elist = (TEntryList *) elist->Clone();
151 if (elist)
152 fOutput->Add(elist);
153 else
155 }
156 }
157 if (fillList) Info("SlaveBegin", "creating an entry-list");
158 if (option.Contains("useList")) useList = kTRUE;
159}
160
162 // function called at the end of the event loop
163
164 hdmd = dynamic_cast<TH1F*>(fOutput->FindObject("hdmd"));
165 h2 = dynamic_cast<TH2F*>(fOutput->FindObject("h2"));
166
167 if (hdmd == 0 || h2 == 0) {
168 Error("Terminate", "hdmd = %p , h2 = %p", hdmd, h2);
169 return;
170 }
171
172 //create the canvas for the h1analysis fit
173 gStyle->SetOptFit();
174 TCanvas *c1 = new TCanvas("c1","h1analysis analysis",10,10,800,600);
175 c1->SetBottomMargin(0.15);
176 hdmd->GetXaxis()->SetTitle("m_{K#pi#pi} - m_{K#pi}[GeV/c^{2}]");
178
179 //fit histogram hdmd with function f5 using the loglfIkelihood option
180 if (gROOT->GetListOfFunctions()->FindObject("f5"))
181 delete gROOT->GetFunction("f5");
182 TF1 *f5 = new TF1("f5",fdm5,0.139,0.17,5);
183 f5->SetParameters(1000000, .25, 2000, .1454, .001);
184 hdmd->Fit("f5","lr");
185
186 //create the canvas for tau d0
187 gStyle->SetOptFit(0);
188 gStyle->SetOptStat(1100);
189 TCanvas *c2 = new TCanvas("c2","tauD0",100,100,800,600);
190 c2->SetGrid();
191 c2->SetBottomMargin(0.15);
192
193 // Project slices of 2-d histogram h2 along X , then fit each slice
194 // with function f2 and make a histogram for each fit parameter
195 // Note that the generated histograms are added to the list of objects
196 // in the current directory.
197 if (gROOT->GetListOfFunctions()->FindObject("f2"))
198 delete gROOT->GetFunction("f2");
199 TF1 *f2 = new TF1("f2",fdm2,0.139,0.17,2);
200 f2->SetParameters(10000, 10);
201 h2->FitSlicesX(f2,0,-1,1,"qln");
202 TH1D *h2_1 = (TH1D*)gDirectory->Get("h2_1");
203 h2_1->GetXaxis()->SetTitle("#tau[ps]");
204 h2_1->SetMarkerStyle(21);
205 h2_1->Draw();
206 c2->Update();
207 TLine *line = new TLine(0,0,0,c2->GetUymax());
208 line->Draw();
209
210 // Have the number of entries on the first histogram (to cross check when running
211 // with entry lists)
212 TPaveStats *psdmd = (TPaveStats *)hdmd->GetListOfFunctions()->FindObject("stats");
213 psdmd->SetOptStat(1110);
214 c1->Modified();
215
216 //save the entry list to a Root file if one was produced
217 if (fillList) {
218 if (!elist)
219 elist = dynamic_cast<TEntryList*>(fOutput->FindObject("elist"));
220 if (elist) {
221 Printf("Entry list 'elist' created:");
222 elist->Print();
223 TFile efile("elist.root","recreate");
224 elist->Write();
225 } else {
226 Error("Terminate", "entry list requested but not found in output");
227 }
228 }
229 // Notify the amount of processed events
230 if (!fInput) Info("Terminate", "processed %lld events", fProcessed);
231}
232
234
235}
236
238// called when loading a new file
239// get branch pointers
240
241 Info("Notify","processing file: %s",myTreeReader.GetTree()->GetCurrentFile()->GetName());
242
243 if (elist && myTreeReader.GetTree()) {
244 if (fillList) {
246 } else if (useList) {
248 }
249 }
250 return kTRUE;
251}
#define f(i)
Definition RSha256.hxx:104
const Bool_t kFALSE
Definition RtypesCore.h:92
bool Bool_t
Definition RtypesCore.h:63
double Double_t
Definition RtypesCore.h:59
long long Long64_t
Definition RtypesCore.h:73
const Bool_t kTRUE
Definition RtypesCore.h:91
#define gDirectory
Definition TDirectory.h:290
#define gROOT
Definition TROOT.h:406
void Printf(const char *fmt,...)
R__EXTERN TStyle * gStyle
Definition TStyle.h:412
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title.
Definition TAttAxis.cxx:293
virtual void SetMarkerStyle(Style_t mstyle=1)
Set the marker style.
Definition TAttMarker.h:40
The Canvas class.
Definition TCanvas.h:23
A List of entry numbers in a TTree or TChain.
Definition TEntryList.h:26
virtual void SetTree(const TTree *tree)
If a list for a tree with such name and filename exists, sets it as the current sublist If not,...
virtual void SetDirectory(TDirectory *dir)
Add reference to directory dir. dir can be 0.
virtual Bool_t Enter(Long64_t entry, TTree *tree=0)
Add entry #entry to the list.
virtual void Print(const Option_t *option="") const
Print this list.
1-Dim function class
Definition TF1.h:213
virtual void SetParameters(const Double_t *params)
Definition TF1.h:644
A ROOT file is a suite of consecutive data records (TKey instances) with a well defined format.
Definition TFile.h:54
1-D histogram with a double per channel (see TH1 documentation)}
Definition TH1.h:618
1-D histogram with a float per channel (see TH1 documentation)}
Definition TH1.h:575
TAxis * GetXaxis()
Get the behaviour adopted by the object about the statoverflows. See EStatOverflows for more informat...
Definition TH1.h:320
virtual TFitResultPtr Fit(const char *formula, Option_t *option="", Option_t *goption="", Double_t xmin=0, Double_t xmax=0)
Fit histogram with function fname.
Definition TH1.cxx:3892
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition TH1.cxx:3350
TList * GetListOfFunctions() const
Definition TH1.h:243
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
Int_t Fill(Double_t)
Invalid Fill method.
Definition TH2.cxx:294
virtual void FitSlicesX(TF1 *f1=0, Int_t firstybin=0, Int_t lastybin=-1, Int_t cut=0, Option_t *option="QNR", TObjArray *arr=0)
Project slices along X in case of a 2-D histogram, then fit each slice with function f1 and make a hi...
Definition TH2.cxx:915
TObject * FindObject(const char *name) const
Find object using its name.
A simple line.
Definition TLine.h:22
virtual void Add(TObject *obj)
Definition TList.h:87
virtual TObject * FindObject(const char *name) const
Find an object in this list using its name.
Definition TList.cxx:578
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition TNamed.cxx:164
virtual TObject * Clone(const char *newname="") const
Make a clone of an object using the Streamer facility.
Definition TNamed.cxx:74
virtual const char * GetName() const
Returns name of object.
Definition TNamed.h:47
virtual Int_t Write(const char *name=0, Int_t option=0, Int_t bufsize=0)
Write this object to the current directory.
Definition TObject.cxx:798
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:893
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition TObject.cxx:197
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition TObject.cxx:867
The histogram statistics painter class.
Definition TPaveStats.h:18
void SetOptStat(Int_t stat=1)
Set the stat option.
TList * fInput
List of objects available during processing.
Definition TSelector.h:41
TSelectorList * fOutput
! List of objects created during processing
Definition TSelector.h:42
virtual const char * GetOption() const
Definition TSelector.h:57
Basic string class.
Definition TString.h:136
const char * Data() const
Definition TString.h:369
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition TString.h:624
void SetOptStat(Int_t stat=1)
The type of information printed in the histogram statistics box can be selected via the parameter mod...
Definition TStyle.cxx:1589
void SetOptFit(Int_t fit=1)
The type of information about fit parameters printed in the histogram statistics box can be selected ...
Definition TStyle.cxx:1541
T & At(std::size_t idx)
TTree * GetTree() const
EEntryStatus SetLocalEntry(Long64_t entry)
Set the next local tree entry.
A TTree represents a columnar dataset.
Definition TTree.h:79
TFile * GetCurrentFile() const
Return pointer to the current file.
Definition TTree.cxx:5459
virtual void SetEntryList(TEntryList *list, Option_t *opt="")
Set an EntryList.
Definition TTree.cxx:9005
TTreeReaderValue< Int_t > fIpi
void Init(TTree *myTree)
TTreeReaderArray< Int_t > fNhitrp
Bool_t Process(Long64_t entry)
TTreeReaderValue< Float_t > fRpd0_t
TTreeReaderArray< Float_t > fNlhpi
TTreeReaderValue< Float_t > fMd0_d
void SlaveBegin(TTree *)
TTreeReaderValue< Float_t > fPtd0_d
TTreeReaderValue< Int_t > fIk
Bool_t Notify()
This method must be overridden to handle object notification.
TTreeReaderValue< Int_t > fIpis
TTreeReaderValue< Int_t > fNjets
TTreeReaderArray< Float_t > fRend
TTreeReaderValue< Float_t > fEtads_d
void Begin(TTree *)
TTreeReaderArray< Float_t > fRstart
TTreeReaderValue< Float_t > fDm_d
TTreeReaderArray< Float_t > fNlhk
TTreeReaderValue< Float_t > fPtds_d
TLine * line
Double_t fdm5(Double_t *xx, Double_t *par)
const Double_t sigma
const Double_t dxbin
Double_t fdm2(Double_t *xx, Double_t *par)
return c1
Definition legend1.C:41
Double_t x[n]
Definition legend1.C:17
return c2
Definition legend2.C:14
Double_t Exp(Double_t x)
Definition TMath.h:727
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition TMath.h:735
Short_t Abs(Short_t d)
Definition TMathBase.h:120