ROOT  6.07/01
Reference Guide
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
h1analysisTreeReader.C
Go to the documentation of this file.
1 #include "h1analysisTreeReader.h"
2 #include "TStyle.h"
3 #include "TCanvas.h"
4 #include "TPaveStats.h"
5 #include "TLine.h"
6 #include "TMath.h"
7 #include "TFile.h"
8 #include "TROOT.h"
9 
10 
11 const Double_t dxbin = (0.17-0.13)/40; // Bin-width
12 const Double_t sigma = 0.0012;
13 
14 //_____________________________________________________________________
16 {
17  Double_t x = xx[0];
18  if (x <= 0.13957) return 0;
19  Double_t xp3 = (x-par[3])*(x-par[3]);
20  Double_t res = dxbin*(par[0]*TMath::Power(x-0.13957, par[1])
21  + par[2] / 2.5066/par[4]*TMath::Exp(-xp3/2/par[4]/par[4]));
22  return res;
23 }
24 
25 //_____________________________________________________________________
27 {
28  Double_t x = xx[0];
29  if (x <= 0.13957) return 0;
30  Double_t xp3 = (x-0.1454)*(x-0.1454);
31  Double_t res = dxbin*(par[0]*TMath::Power(x-0.13957, 0.25)
32  + par[1] / 2.5066/sigma*TMath::Exp(-xp3/2/sigma/sigma));
33  return res;
34 }
35 //_____________________________________________________________________
37 // entry is the entry number in the current Tree
38 // Selection function to select D* and D0.
40  fProcessed++;
41  //in case one entry list is given in input, the selection has already been done.
42  if (!useList) {
43  // Return as soon as a bad entry is detected
44  if (TMath::Abs(*fMd0_d-1.8646) >= 0.04) return kFALSE;
45  if (*fPtds_d <= 2.5) return kFALSE;
46  if (TMath::Abs(*fEtads_d) >= 1.5) return kFALSE;
47  (*fIk)--; //original fIk used f77 convention starting at 1
48  (*fIpi)--;
49 
50 
51  if (fNhitrp.At(*fIk)* fNhitrp.At(*fIpi) <= 1) return kFALSE;
52 
53 
54  if (fRend.At(*fIk) -fRstart.At(*fIk) <= 22) return kFALSE;
55  if (fRend.At(*fIpi)-fRstart.At(*fIpi) <= 22) return kFALSE;
56  if (fNlhk.At(*fIk) <= 0.1) return kFALSE;
57  if (fNlhpi.At(*fIpi) <= 0.1) return kFALSE;
58  (*fIpis)--; if (fNlhpi.At(*fIpis) <= 0.1) return kFALSE;
59  if (*fNjets < 1) return kFALSE;
60  }
61  // if option fillList, fill the entry list
62  if (fillList) elist->Enter(entry);
63 
64  //fill some histograms
65  hdmd->Fill(*fDm_d);
66  h2->Fill(*fDm_d,*fRpd0_t/0.029979*1.8646/ *fPtd0_d);
67 
68  return kTRUE;
69 }
70 
72 // function called before starting the event loop
73 // -it performs some cleanup
74 // -it creates histograms
75 // -it sets some initialisation for the entry list
76 
77  Reset();
78 
79  //print the option specified in the Process function.
80  TString option = GetOption();
81  Info("Begin", "starting h1analysis with process option: %s", option.Data());
82 
83  delete gDirectory->GetList()->FindObject("elist");
84 
85  // case when one creates/fills the entry list
86  if (option.Contains("fillList")) {
87  fillList = kTRUE;
88  elist = new TEntryList("elist", "H1 selection from Cut");
89  // Add to the input list for processing in PROOF, if needed
90  if (fInput) {
91  fInput->Add(new TNamed("fillList",""));
92  // We send a clone to avoid double deletes when importing the result
93  fInput->Add(elist);
94  // This is needed to avoid warnings from output-to-members mapping
95  elist = 0;
96  }
97  }
98  if (fillList) Info("Begin", "creating an entry-list");
99  // case when one uses the entry list generated in a previous call
100  if (option.Contains("useList")) {
101  useList = kTRUE;
102  if (fInput) {
103  // Option "useList" not supported in PROOF directly
104  Warning("Begin", "option 'useList' not supported in PROOF - ignoring");
105  Warning("Begin", "the entry list must be set on the chain *before* calling Process");
106  } else {
107  TFile f("elist.root");
108  elist = (TEntryList*)f.Get("elist");
109  if (elist) elist->SetDirectory(0); //otherwise the file destructor will delete elist
110  }
111  }
112 }
113 
115 
116 // function called before starting the event loop
117 // -it performs some cleanup
118 // -it creates histograms
119 // -it sets some initialisation for the entry list
120 
121  Init(myTree);
122 
123  //print the option specified in the Process function.
124  TString option = GetOption();
125  Info("SlaveBegin",
126  "starting h1analysis with process option: %s (tree: %p)", option.Data(), myTree);
127 
128  //create histograms
129  hdmd = new TH1F("hdmd","Dm_d",40,0.13,0.17);
130  h2 = new TH2F("h2","ptD0 vs Dm_d",30,0.135,0.165,30,-3,6);
131 
132  fOutput->Add(hdmd);
133  fOutput->Add(h2);
134 
135  // Entry list stuff (re-parse option because on PROOF only SlaveBegin is called)
136  if (option.Contains("fillList")) {
137  fillList = kTRUE;
138  // Get the list
139  if (fInput) {
140  if ((elist = (TEntryList *) fInput->FindObject("elist")))
141  // Need to clone to avoid problems when destroying the selector
142  elist = (TEntryList *) elist->Clone();
143  if (elist)
144  fOutput->Add(elist);
145  else
146  fillList = kFALSE;
147  }
148  }
149  if (fillList) Info("SlaveBegin", "creating an entry-list");
150 }
151 
153  // function called at the end of the event loop
154 
155  hdmd = dynamic_cast<TH1F*>(fOutput->FindObject("hdmd"));
156  h2 = dynamic_cast<TH2F*>(fOutput->FindObject("h2"));
157 
158  if (hdmd == 0 || h2 == 0) {
159  Error("Terminate", "hdmd = %p , h2 = %p", hdmd, h2);
160  return;
161  }
162 
163  //create the canvas for the h1analysis fit
164  gStyle->SetOptFit();
165  TCanvas *c1 = new TCanvas("c1","h1analysis analysis",10,10,800,600);
166  c1->SetBottomMargin(0.15);
167  hdmd->GetXaxis()->SetTitle("m_{K#pi#pi} - m_{K#pi}[GeV/c^{2}]");
168  hdmd->GetXaxis()->SetTitleOffset(1.4);
169 
170  //fit histogram hdmd with function f5 using the loglfIkelihood option
171  if (gROOT->GetListOfFunctions()->FindObject("f5"))
172  delete gROOT->GetFunction("f5");
173  TF1 *f5 = new TF1("f5",fdm5,0.139,0.17,5);
174  f5->SetParameters(1000000, .25, 2000, .1454, .001);
175  hdmd->Fit("f5","lr");
176 
177  //create the canvas for tau d0
178  gStyle->SetOptFit(0);
179  gStyle->SetOptStat(1100);
180  TCanvas *c2 = new TCanvas("c2","tauD0",100,100,800,600);
181  c2->SetGrid();
182  c2->SetBottomMargin(0.15);
183 
184  // Project slices of 2-d histogram h2 along X , then fit each slice
185  // with function f2 and make a histogram for each fit parameter
186  // Note that the generated histograms are added to the list of objects
187  // in the current directory.
188  if (gROOT->GetListOfFunctions()->FindObject("f2"))
189  delete gROOT->GetFunction("f2");
190  TF1 *f2 = new TF1("f2",fdm2,0.139,0.17,2);
191  f2->SetParameters(10000, 10);
192  h2->FitSlicesX(f2,0,-1,1,"qln");
193  TH1D *h2_1 = (TH1D*)gDirectory->Get("h2_1");
194  h2_1->GetXaxis()->SetTitle("#tau[ps]");
195  h2_1->SetMarkerStyle(21);
196  h2_1->Draw();
197  c2->Update();
198  TLine *line = new TLine(0,0,0,c2->GetUymax());
199  line->Draw();
200 
201  // Have the number of entries on the first histogram (to cross check when running
202  // with entry lists)
203  TPaveStats *psdmd = (TPaveStats *)hdmd->GetListOfFunctions()->FindObject("stats");
204  psdmd->SetOptStat(1110);
205  c1->Modified();
206 
207  //save the entry list to a Root file if one was produced
208  if (fillList) {
209  if (!elist)
210  elist = dynamic_cast<TEntryList*>(fOutput->FindObject("elist"));
211  if (elist) {
212  Printf("Entry list 'elist' created:");
213  elist->Print();
214  TFile efile("elist.root","recreate");
215  elist->Write();
216  } else {
217  Error("Terminate", "entry list requested but not found in output");
218  }
219  }
220  // Notify the amount of processed events
221  if (!fInput) Info("Terminate", "processed %lld events", fProcessed);
222 }
223 
225 
226 }
227 
229 // called when loading a new file
230 // get branch pointers
231 
232  Info("Notify","processing file: %s",myTreeReader.GetTree()->GetCurrentFile()->GetName());
233 
234  if (elist && myTreeReader.GetTree()) {
235  if (fillList) {
237  } else if (useList) {
239  }
240  }
241  return kTRUE;
242 }
virtual const char * GetOption() const
Definition: TSelector.h:65
virtual void SetTitleOffset(Float_t offset=1)
Set distance between the axis and the axis title Offset is a correction factor with respect to the "s...
Definition: TAttAxis.cxx:245
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:823
double par[1]
Definition: unuranDistr.cxx:38
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3159
TList * GetListOfFunctions() const
Definition: TH1.h:244
virtual void SetParameters(const Double_t *params)
Definition: TF1.h:432
TSelectorList * fOutput
Definition: TSelector.h:50
long long Long64_t
Definition: RtypesCore.h:69
TTreeReaderValue< Int_t > fNjets
TLine * line
TCanvas * c1
Definition: legend1.C:2
TObject * FindObject(const char *name) const
Find object using its name.
Definition: THashList.cxx:213
TTree * GetTree() const
Definition: TTreeReader.h:166
R__EXTERN TStyle * gStyle
Definition: TStyle.h:423
tuple f2
Definition: surfaces.py:24
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition: TObject.cxx:892
TTreeReaderArray< Float_t > fNlhpi
A ROOT file is a suite of consecutive data records (TKey instances) with a well defined format...
Definition: TFile.h:45
Double_t GetUymax() const
Definition: TPad.h:229
Double_t fdm2(Double_t *xx, Double_t *par)
virtual TObject * Get(const char *namecycle)
Return pointer to object identified by namecycle.
T & At(size_t idx)
#define gROOT
Definition: TROOT.h:344
Basic string class.
Definition: TString.h:137
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:570
bool Bool_t
Definition: RtypesCore.h:59
TTreeReaderValue< Int_t > fIpis
TTreeReaderArray< Float_t > fRend
const Bool_t kFALSE
Definition: Rtypes.h:92
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition: TObject.cxx:254
virtual TObject * FindObject(const char *name) const
Find an object in this list using its name.
Definition: TList.cxx:497
The histogram statistics painter class.
Definition: TPaveStats.h:28
TFile * GetCurrentFile() const
Return pointer to the current file.
Definition: TTree.cxx:4987
Short_t Abs(Short_t d)
Definition: TMathBase.h:110
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...
TFile * f
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Definition: TMath.h:501
virtual void SetDirectory(TDirectory *dir)
Add reference to directory dir. dir can be 0.
const char * Data() const
Definition: TString.h:349
TTreeReaderValue< Float_t > fPtd0_d
Double_t x[n]
Definition: legend1.C:17
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:33
virtual void SetGrid(Int_t valuex=1, Int_t valuey=1)
Definition: TPad.h:326
TTreeReaderValue< Float_t > fPtds_d
virtual TObject * Clone(const char *newname="") const
Make a clone of an object using the Streamer facility.
Definition: TNamed.cxx:63
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:918
TTreeReaderValue< Float_t > fRpd0_t
Bool_t Process(Long64_t entry)
virtual void SetBottomMargin(Float_t bottommargin)
Set Pad bottom margin in fraction of the pad height.
Definition: TAttPad.cxx:97
TTreeReaderValue< Float_t > fDm_d
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2878
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:891
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:256
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:1204
void f5()
Definition: na49.C:91
A simple line.
Definition: TLine.h:41
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
virtual void SetMarkerStyle(Style_t mstyle=1)
Definition: TAttMarker.h:53
Double_t fdm5(Double_t *xx, Double_t *par)
1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:613
TTreeReaderArray< Int_t > fNhitrp
Long64_t entry
#define Printf
Definition: TGeoToOCC.h:18
const Double_t sigma
TTreeReaderArray< Float_t > fNlhk
Bool_t Notify()
This method must be overridden to handle object notification.
The Canvas class.
Definition: TCanvas.h:48
Double_t Exp(Double_t x)
Definition: TMath.h:495
const Double_t dxbin
return c2
Definition: legend2.C:14
EEntryStatus SetLocalEntry(Long64_t entry)
Definition: TTreeReader.h:160
virtual void SetEntryList(TEntryList *list, Option_t *opt="")
Set an EntryList.
Definition: TTree.cxx:8141
double Double_t
Definition: RtypesCore.h:55
TTreeReaderValue< Float_t > fMd0_d
TTreeReaderValue< Float_t > fEtads_d
TTreeReaderArray< Float_t > fRstart
virtual Bool_t Enter(Long64_t entry, TTree *tree=0)
Add entry entry to the list.
Definition: TEntryList.cxx:558
TTreeReaderValue< Int_t > fIk
void Init(TTree *myTree)
TList * fInput
Current object if processing object (vs. TTree)
Definition: TSelector.h:49
virtual void Add(TObject *obj)
Definition: TList.h:81
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:567
1-Dim function class
Definition: TF1.h:149
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:1252
A TTree object has a header with a name and a title.
Definition: TTree.h:98
#define gDirectory
Definition: TDirectory.h:221
virtual void Update()
Update canvas pad buffers.
Definition: TCanvas.cxx:2179
virtual void Print(const Option_t *option="") const
Print this list.
Definition: TEntryList.cxx:989
const Bool_t kTRUE
Definition: Rtypes.h:91
Int_t Fill(Double_t)
Invalid Fill method.
Definition: TH2.cxx:287
virtual void SetTitle(const char *title="")
Change (i.e. set) the title of the TNamed.
Definition: TNamed.cxx:152
A List of entry numbers in a TTree or TChain.
Definition: TEntryList.h:27
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:3607
void Modified(Bool_t flag=1)
Definition: TPad.h:407
TTreeReaderValue< Int_t > fIpi
TAxis * GetXaxis()
Definition: TH1.h:319
void SetOptStat(Int_t stat=1)
Set the stat option.
Definition: TPaveStats.cxx:303
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:904