LOGO

RestFrames  v1.0.1
RestFrames HEP Event Analysis Software Library
example_Zll.C
Go to the documentation of this file.
1 // RestFrames: particle physics event analysis library
3 // --------------------------------------------------------------------
4 // Copyright (c) 2014-2016, Christopher Rogan
13 //
14 // This file is part of RestFrames.
15 //
16 // RestFrames is free software; you can redistribute it and/or modify
17 // it under the terms of the GNU General Public License as published by
18 // the Free Software Foundation; either version 2 of the License, or
19 // (at your option) any later version.
20 //
21 // RestFrames is distributed in the hope that it will be useful,
22 // but WITHOUT ANY WARRANTY; without even the implied warranty of
23 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24 // GNU General Public License for more details.
25 //
26 // You should have received a copy of the GNU General Public License
27 // along with RestFrames. If not, see <http://www.gnu.org/licenses/>.
29 
30 #define COMPILER (!defined(__CINT__) && !defined(__CLING__))
31 #if defined(__MAKECINT__) || defined(__ROOTCLING__) || COMPILER
32 #include "RestFrames/RestFrames.hh"
33 #else
34 RestFrames::RFKey ensure_autoload(1);
35 #endif
36 
37 using namespace RestFrames;
38 
39 void example_Zll(const std::string& output_name = "output_Zll.root"){
40 
41  double mZ = 91.188; // GeV, PDG 2016
42  double wZ = 2.495;
43 
44  // Number of events to generate
45  int Ngen = 100000;
46 
48  g_Log << LogInfo << "Initializing generator frames and tree..." << LogEnd;
50  LabGenFrame LAB_Gen("LAB_Gen","LAB");
51  ResonanceGenFrame Z_Gen("Z_Gen","Z");
52  VisibleGenFrame Lp_Gen("Lp_Gen","#it{l}^{+}");
53  VisibleGenFrame Lm_Gen("Lm_Gen","#it{l}^{-}");
54 
55  //-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//
56 
57  LAB_Gen.SetChildFrame(Z_Gen);
58  Z_Gen.AddChildFrame(Lp_Gen);
59  Z_Gen.AddChildFrame(Lm_Gen);
60 
61  if(LAB_Gen.InitializeTree())
62  g_Log << LogInfo << "...Successfully initialized generator tree" << LogEnd;
63  else
64  g_Log << LogError << "...Failed initializing generator tree" << LogEnd;
65 
66  //-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//
67 
68  // Set Z pole mass and width
69  Z_Gen.SetMass(mZ); Z_Gen.SetWidth(wZ);
70 
71  // set lepton pT and eta cuts
72  Lp_Gen.SetPtCut(15.); Lp_Gen.SetEtaCut(2.5);
73  Lm_Gen.SetPtCut(15.); Lm_Gen.SetEtaCut(2.5);
74 
75  if(LAB_Gen.InitializeAnalysis())
76  g_Log << LogInfo << "...Successfully initialized generator analysis" << std::endl << LogEnd;
77  else
78  g_Log << LogError << "...Failed initializing generator analysis" << LogEnd;
81 
83  g_Log << LogInfo << "Initializing reconstruction frames and trees..." << LogEnd;
85  LabRecoFrame LAB("LAB","LAB");
86  DecayRecoFrame Z("Z","Z");
87  VisibleRecoFrame Lp("Lp","#it{l}^{+}");
88  VisibleRecoFrame Lm("Lm","#it{l}^{-}");
89 
90  //-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//
91 
92  LAB.SetChildFrame(Z);
93  Z.AddChildFrame(Lp);
94  Z.AddChildFrame(Lm);
95 
96  if(LAB.InitializeTree())
97  g_Log << LogInfo << "...Successfully initialized reconstruction trees" << LogEnd;
98  else
99  g_Log << LogError << "...Failed initializing reconstruction trees" << LogEnd;
100 
101  if(LAB.InitializeAnalysis())
102  g_Log << LogInfo << "...Successfully initialized analyses" << LogEnd;
103  else
104  g_Log << LogError << "...Failed initializing analyses" << LogEnd;
105 
108 
109  TreePlot* tree_plot = new TreePlot("TreePlot","TreePlot");
110 
111  // generator tree
112  tree_plot->SetTree(LAB_Gen);
113  tree_plot->Draw("GenTree", "Generator Tree", true);
114 
115  // reco tree
116  tree_plot->SetTree(LAB);
117  tree_plot->Draw("RecoTree", "Reconstruction Tree");
118 
119  //-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//-//
120 
121  // Declare observables for histogram booking
122  HistPlot* hist_plot = new HistPlot("HistPlot","Z #rightarrow #it{l}^{+} #it{l}^{-}");
123 
124  const HistPlotVar& MZ = hist_plot->GetNewVar("MZ", "M_{Z}", 70., 110., "[GeV]");
125  const HistPlotVar& cosZ = hist_plot->GetNewVar("cosZ","cos #theta_{Z}", -1., 1.);
126  const HistPlotVar& dphiZ = hist_plot->GetNewVar("dphiZ", "#Delta #phi_{Z}", 0., 2.*acos(-1.));
127 
128  hist_plot->AddPlot(MZ);
129  hist_plot->AddPlot(cosZ);
130  hist_plot->AddPlot(dphiZ);
131  hist_plot->AddPlot(cosZ, MZ);
132 
133  for(int igen = 0; igen < Ngen; igen++){
134  if(igen%((std::max(Ngen,10))/10) == 0)
135  g_Log << LogInfo << "Generating event " << igen << " of " << Ngen << LogEnd;
136 
137  // generate event
138  LAB_Gen.ClearEvent(); // clear the gen tree
139 
140  double PTZ = mZ*gRandom->Rndm();
141  LAB_Gen.SetTransverseMomentum(PTZ); // give the Z some Pt
142  double PzZ = mZ*(2.*gRandom->Rndm()-1.);
143  LAB_Gen.SetLongitudinalMomentum(PzZ); // give the Z some Pz
144 
145  LAB_Gen.AnalyzeEvent(); // generate a new event
146 
147  // analyze event
148  LAB.ClearEvent(); // clear the reco tree
149 
150  Lp.SetLabFrameFourVector(Lp_Gen.GetFourVector(), 1); // Set lepton 4-vec and charge
151  Lm.SetLabFrameFourVector(Lm_Gen.GetFourVector(),-1); // Set lepton 4-vec and charge
152 
153  LAB.AnalyzeEvent(); // analyze the event
154 
155  // calculate observables
156  MZ = Z.GetMass();
157  cosZ = Z.GetCosDecayAngle();
158  dphiZ = LAB.GetDeltaPhiDecayPlanes(Z);
159 
160  hist_plot->Fill();
161  }
162 
163  hist_plot->Draw();
164 
165  TFile fout(output_name.c_str(),"RECREATE");
166  fout.Close();
167  hist_plot->WriteOutput(output_name);
168  hist_plot->WriteHist(output_name);
169  tree_plot->WriteOutput(output_name);
170 
171 }
172 
173 # ifndef __CINT__ // main function for stand-alone compilation
174 int main(){
175  example_Zll();
176  return 0;
177 }
178 #endif