CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
TTbar_Kinematics.cc
Go to the documentation of this file.
5 
6 
7 using namespace edm;
9  genEventInfoProductTag_(iConfig.getParameter<edm::InputTag>("genEventInfoProductTag"))
10 {
11  dbe = 0;
13 }
14 
15 
17 {
18 
19  // do anything here that needs to be done at desctruction time
20  // (e.g. close files, deallocate resources etc.)
21 
22 }
23 
24 
25 //
26 // member functions
27 //
28 
29 // ------------ method called for each event ------------
30 void
32 {
33 
34  // --- the MC weights ---
36  iEvent.getByLabel(genEventInfoProductTag_, evt_info);
37  if(!evt_info.isValid()) return;
38  weight = evt_info->weight() ;
39 
40  // --- get TopQuarkAnalysis TtGenEvent
42  iEvent.getByLabel("genEvt", genEvt);
43 
44  if(!genEvt.isValid())return;
45 
46  const reco::GenParticle* top = 0;
47  const reco::GenParticle* antitop = 0;
48  const reco::GenParticle* bottom = 0;
49  const reco::GenParticle* antibottom = 0;
50  const reco::GenParticle* Wplus = 0;
51  const reco::GenParticle* Wmin = 0;
52 
53  top = genEvt->top();
54  antitop = genEvt->topBar();
55  bottom = genEvt->b();
56  antibottom = genEvt->bBar();
57  Wplus = genEvt->wPlus();
58  Wmin = genEvt->wMinus();
59 
60  tlv_Top = TLorentzVector(0,0,0,0) ;
61  tlv_TopBar = TLorentzVector(0,0,0,0) ;
62  tlv_Bottom = TLorentzVector(0,0,0,0) ;
63  tlv_BottomBar = TLorentzVector(0,0,0,0) ;
64  tlv_Wplus = TLorentzVector(0,0,0,0) ;
65  tlv_Wmin = TLorentzVector(0,0,0,0) ;
66  tlv_TTbar = TLorentzVector(0,0,0,0);
67 
68  if(top) tlv_Top.SetPxPyPzE(top->p4().px(),top->p4().py(),top->p4().pz(),top->p4().e());
69  if(antitop) tlv_TopBar.SetPxPyPzE(antitop->p4().px(),antitop->p4().py(),antitop->p4().pz(),antitop->p4().e());
70  if(bottom) tlv_Bottom.SetPxPyPzE(bottom->p4().px(),bottom->p4().py(),bottom->p4().pz(),bottom->p4().e());
71  if(antibottom) tlv_BottomBar.SetPxPyPzE(antibottom->p4().px(),antibottom->p4().py(),antibottom->p4().pz(),antibottom->p4().e());
72  if(Wplus) tlv_Wplus.SetPxPyPzE(Wplus->p4().px(),Wplus->p4().py(),Wplus->p4().pz(),Wplus->p4().e());
73  if(Wmin) tlv_Wmin.SetPxPyPzE(Wmin->p4().px(),Wmin->p4().py(),Wmin->p4().pz(),Wmin->p4().e());
75 
76  //---topquarkquantities---
77  nEvt->Fill(0.5,weight);
78  if(top && antitop){
79  hTopPt->Fill(tlv_Top.Pt(),weight);
80  hTopPt->Fill(tlv_TopBar.Pt(),weight);
81 
82  hTopY->Fill(tlv_Top.Rapidity(),weight);
83  hTopY->Fill(tlv_TopBar.Rapidity(),weight);
84 
86  hTopMass->Fill(tlv_TopBar.M(),weight);
87 
88  //---ttbarpairquantities---
91 
92  hTTbarY->Fill(tlv_TTbar.Rapidity(),weight);
93  hTTbarY->Fill(tlv_TTbar.Rapidity(),weight);
94 
97  }
98  if(bottom && antibottom){
101 
104 
105  //hBottomY->Fill(math::XYZTLorentzVector(bottom->momentum()).Rapidity(),weight);
106  //hBottomY->Fill(math::XYZTLorentzVector(antibottom->momentum()).Rapidity(),weight);
107 
108  hBottomY->Fill(tlv_Bottom.Rapidity(),weight);
109  hBottomY->Fill(tlv_BottomBar.Rapidity(),weight);
110 
113 
116 
119 
122 
125 
128 
129  hBottomYPz->Fill(tlv_Bottom.Rapidity(),tlv_Bottom.Pz(),weight);
131 
134 
137 
138  hBottomMassY->Fill(tlv_Bottom.M(),tlv_Bottom.Rapidity(),weight);
140 
143  }
144  if(Wplus && Wmin){
145  hWplusPz->Fill(tlv_Wplus.Pz(),weight);
146  hWminPz->Fill(tlv_Wmin.Pz(),weight);
147  }
148 }
149 
150 
151 // ------------ method called once each job just before starting event loop ------------
152 void
154 {
155  if(!dbe) return;
156  dbe->setCurrentFolder("Generator/TTbar");
157 
158  nEvt = dbe->book1D("nEvt", "n analyzed Events", 1, 0., 1.);
159 
160  hTopPt = dbe->book1D("TTbar_TopPt","t quark transverse momentum",1000,0.,1000.); hTopPt->setAxisTitle("t quark transverse momentum");
161  hTopY = dbe->book1D("TTbar_TopY","t quark rapidity",200,-5.,5.); hTopY->setAxisTitle("t quark rapidity");
162  hTopMass = dbe->book1D("TTbar_TopMass","t quark mass",500,0.,500.); hTopMass->setAxisTitle("t quark mass");
163 
164  hTTbarPt = dbe->book1D("TTbar_TTbarPt","tt pair transverse momentum",1000,0.,1000.); hTTbarPt->setAxisTitle("tt pair transverse momentum");
165  hTTbarY = dbe->book1D("TTbar_TTbarY","tt pair rapidity",200,-5.,5.); hTTbarY->setAxisTitle("tt pair rapidity");
166  hTTbarMass = dbe->book1D("TTbar_TTbarMass","tt pair mass",1000,0.,1000.); hTTbarMass->setAxisTitle("tt pair mass");
167 
168  hBottomPt = dbe->book1D("TTbar_BottomPt","b quark transverse momentum",1000,0.,1000.); hBottomPt->setAxisTitle("b quark transverse momentum");
169  hBottomEta = dbe->book1D("TTbar_BottomEta","b quark pseudo-rapidity",200,-5.,5.); hBottomEta->setAxisTitle("b quark pseudo-rapidity");
170  hBottomY = dbe->book1D("TTbar_BottomY","b quark rapidity",200,-5.,5.); hBottomY->setAxisTitle("b quark rapidity");
171  hBottomPz = dbe->book1D("TTbar_BottomPz","b quark longitudinal momentum",200,-100.,100.); hBottomPz->setAxisTitle("b quark longitudinal momentum");
172  hBottomE = dbe->book1D("TTbar_BottomE","b quark energy",1000,0.,1000.); hBottomE->setAxisTitle("b quark energy");
173  hBottomMass = dbe->book1D("TTbar_BottomMass","b quark mass",50,0.,5.); hBottomMass->setAxisTitle("b quark mass");
174 
175  hBottomPtPz = dbe->book2D("TTbar_BottomPtPz","b quark longitudinal vs transverse momentum",1000,0.,1000.,200,-100.,100.); hBottomPtPz->setAxisTitle("P_{z} (GeV)",1);hBottomPtPz->setAxisTitle("P_{t} (GeV)",2);
176  hBottomEtaPz = dbe->book2D("TTbar_BottomEtaPz","b quark longitudinal momentum vs pseudorapidity",200,-5.,5.,200,-100.,100.); hBottomEtaPz->setAxisTitle("#eta",1);hBottomEtaPz->setAxisTitle("P_{z} (GeV)",1);
177  hBottomEtaPt = dbe->book2D("TTbar_BottomEtaPt"," quark transveral momentum vs pseudorapidity",200,-5.,5.,1000,0.,1000.); hBottomEtaPt->setAxisTitle("#eta");hBottomEtaPt->setAxisTitle("P_{t} (GeV)");
178  hBottomYPz = dbe->book2D("TTbar_BottomYPz","b quark longitudinal momentum vs rapidity",200,-5.,5.,200,-100.,100.); hBottomYPz->setAxisTitle("Y");hBottomYPz->setAxisTitle("P_{z} (GeV)");
179  hBottomMassPz = dbe->book2D("TTbar_BottomMassPz","b quark longitudinal momentum vs mass",50,0.,5.,200,-100.,100.); hBottomMassPz->setAxisTitle("M (GeV)");hBottomMassPz->setAxisTitle("P_{z} (GeV)");
180  hBottomMassEta = dbe->book2D("TTbar_BottomMassEta","b quark pseudorapidity vs mass",50,0.,5.,200,-5.,5.); hBottomMassEta->setAxisTitle("M (GeV)");hBottomMassEta->setAxisTitle("#eta");
181  hBottomMassY = dbe->book2D("TTbar_BottomMassY","b quark rapidity vs mass",50,0.,5.,200,-5.,5.); hBottomMassY->setAxisTitle("M (GeV)"); hBottomMassY->setAxisTitle("Y");
182  hBottomMassDeltaY = dbe->book2D("TTbar_BottomMassDeltaY","b quark pseudorapidity - rapidity vs mass",50,0.,50.,2000,-5.,5.); hBottomMassDeltaY->setAxisTitle("M (GeV)");hBottomMassDeltaY->setAxisTitle("Y");
183 
184  hWplusPz = dbe->book1D("TTbar_WplusPz","W+ boson longitudinal momentum",200,-100.,100.); hWplusPz->setAxisTitle("W+ boson longitudinal momentum");
185  hWminPz = dbe->book1D("TTbar_WminPz","W- boson longitudinal momentum",200,-100.,100.); hWminPz->setAxisTitle("W- boson longitudinal momentum");
186 
187 }
188 
189 // ------------ method called once each job just after ending the event loop ------------
190 void
192 {
193 }
194 
195 // ------------ method called when starting to processes a run ------------
196 void
198 {
199 }
200 
201 // ------------ method called when ending the processing of a run ------------
202 void
204 {
205 }
206 
207 // ------------ method called when starting to processes a luminosity block ------------
208 void
210 {
211 }
212 
213 // ------------ method called when ending the processing of a luminosity block ------------
214 void
216 {
217 }
218 
219 // ------------ method fills 'descriptions' with the allowed parameters for the module ------------
220 void
222  //The following says we do not know what parameters are allowed so do no validation
223  // Please change this to state exactly what you do use, even if it is no parameters
225  desc.setUnknown();
226  descriptions.addDefault(desc);
227 }
MonitorElement * hBottomYPz
MonitorElement * nEvt
TLorentzVector tlv_Bottom
MonitorElement * book1D(const char *name, const char *title, int nchX, double lowX, double highX)
Book 1D histogram.
Definition: DQMStore.cc:722
MonitorElement * hBottomMass
virtual const LorentzVector & p4() const GCC11_FINAL
four-momentum Lorentz vector
MonitorElement * hBottomPtPz
virtual void endRun(edm::Run const &, edm::EventSetup const &)
virtual void endLuminosityBlock(edm::LuminosityBlock const &, edm::EventSetup const &)
MonitorElement * hTTbarMass
virtual void beginLuminosityBlock(edm::LuminosityBlock const &, edm::EventSetup const &)
virtual void analyze(const edm::Event &, const edm::EventSetup &)
MonitorElement * hBottomMassY
TLorentzVector tlv_Wplus
MonitorElement * hBottomEtaPt
virtual void endJob()
void Fill(long long x)
virtual void beginJob()
MonitorElement * hBottomEta
TLorentzVector tlv_TopBar
int iEvent
Definition: GenABIO.cc:243
void addDefault(ParameterSetDescription const &psetDescription)
MonitorElement * hTopMass
MonitorElement * hTTbarY
TLorentzVector tlv_TTbar
MonitorElement * hWminPz
bool isValid() const
Definition: HandleBase.h:76
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:361
TLorentzVector tlv_BottomBar
MonitorElement * hBottomMassPz
MonitorElement * hBottomMassEta
edm::InputTag genEventInfoProductTag_
MonitorElement * hBottomPz
MonitorElement * hTopY
TLorentzVector tlv_Wmin
MonitorElement * hTTbarPt
MonitorElement * hBottomMassDeltaY
MonitorElement * hBottomPt
MonitorElement * hTopPt
TLorentzVector tlv_Top
MonitorElement * hBottomE
MonitorElement * hWplusPz
DQMStore * dbe
ME&#39;s &quot;container&quot;.
virtual void beginRun(edm::Run const &, edm::EventSetup const &)
MonitorElement * book2D(const char *name, const char *title, int nchX, double lowX, double highX, int nchY, double lowY, double highY)
Book 2D histogram.
Definition: DQMStore.cc:850
void setAxisTitle(const std::string &title, int axis=1)
set x-, y- or z-axis title (axis=1, 2, 3 respectively)
MonitorElement * hBottomEtaPz
TTbar_Kinematics(const edm::ParameterSet &)
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:434
Definition: Run.h:36
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
MonitorElement * hBottomY