ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/algomez/GenAnalyzer/src/GenAnalyzer.cc
Revision: 1.1
Committed: Mon Jun 3 20:36:36 2013 UTC (11 years, 11 months ago) by algomez
Content type: text/plain
Branch: MAIN
Log Message:
*** empty log message ***

File Contents

# Content
1 // -*- C++ -*-
2 //
3 // Package: GenAnalyzer
4 // Class: GenAnalyzer
5 //
6 /**\class GenAnalyzer GenAnalyzer.cc UserCode/GenAnalyzer/src/GenAnalyzer.cc
7
8 Description: [one line class summary]
9
10 Implementation:
11 [Notes on implementation]
12 */
13 //
14 // Original Author: Alejandro Gomez
15 // Created: Mon May 20 14:38:03 CDT 2013
16 // $Id$
17 //
18 //
19
20
21 // system include files
22 #include <memory>
23 #include <vector>
24
25 // user include files
26 #include "FWCore/Framework/interface/Frameworkfwd.h"
27 #include "FWCore/Framework/interface/EDAnalyzer.h"
28
29 #include "FWCore/Framework/interface/Event.h"
30 #include "FWCore/Framework/interface/MakerMacros.h"
31 #include "DataFormats/FWLite/interface/Handle.h"
32
33 #include "DataFormats/HepMCCandidate/interface/GenParticleFwd.h"
34 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
35
36 #include "DataFormats/Candidate/interface/Candidate.h"
37 #include "DataFormats/Candidate/interface/CandidateFwd.h"
38 #include "DataFormats/Common/interface/Ref.h"
39
40 #include "FWCore/ServiceRegistry/interface/Service.h"
41 #include "CommonTools/UtilAlgos/interface/TFileService.h"
42
43 #include "FWCore/Utilities/interface/InputTag.h"
44 #include "FWCore/ParameterSet/interface/ParameterSet.h"
45 #include "DataFormats/Math/interface/LorentzVector.h"
46 #include "DataFormats/Math/interface/Vector3D.h"
47
48 #include "TLorentzVector.h"
49 #include "TVector3.h"
50 #include "TH1D.h"
51 //
52 // class declaration
53 //
54
55 class GenAnalyzer : public edm::EDAnalyzer {
56 public:
57 explicit GenAnalyzer(const edm::ParameterSet&);
58 ~GenAnalyzer();
59
60
61 private:
62 virtual void beginJob() ;
63 virtual void analyze(const edm::Event&, const edm::EventSetup&);
64 virtual void endJob() ;
65
66 // ----------member data ---------------------------
67 edm::InputTag src_;
68 double stop1Mass_;
69 double stop2Mass_;
70
71 TH1D * histo;
72 // Higgs
73 TH1D * h_higgs_mass;
74 TH1D * h_higgsTrue_mass;
75 TH1D * h_higgs_pt;
76 TH1D * h_higgsBB_deltaR;
77 TH1D * h_higgs_num;
78
79 // Stop1
80 TH1D * h_stop1_mass;
81 TH1D * h_stop1True_mass;
82 TH1D * h_stop1_pt;
83 TH1D * h_stop1Bj_deltaR;
84
85 // Stop2
86 TH1D * h_stop2_mass;
87 TH1D * h_stop2True_mass;
88
89 // Jets + B
90 TH1D * h_jetsB_ht;
91 TH1D * h_jetsB_pt;
92 TH1D * h_jetsB1_pt;
93 TH1D * h_jetsB2_pt;
94 TH1D * h_jetsB3_pt;
95 TH1D * h_jetsB4_pt;
96 TH1D * h_jetsB_num;
97 TH1D * h_jetsB1_eta;
98 TH1D * h_jetsB2_eta;
99 TH1D * h_jetsB3_eta;
100 TH1D * h_jetsB4_eta;
101 TH1D * h_jetsB1_phi;
102 TH1D * h_jetsB2_phi;
103 TH1D * h_jetsB3_phi;
104 TH1D * h_jetsB4_phi;
105 TH1D * h_jetsB1jetsB2_deltaR;
106 TH1D * h_jetsB1jetsB3_deltaR;
107
108 // Bjets
109 TH1D * h_Bjet1_pt;
110 TH1D * h_Bjet2_pt;
111 TH1D * h_Bjet_num;
112 TH1D * h_Bjet1Bjet2_deltaEta;
113 TH1D * h_Bjet1Bjet2_deltaPhi;
114 TH1D * h_Bjet1Bjet2_cosDeltaPhi;
115 TH1D * h_Bjet1Bjet2_deltaR;
116
117 // Jets (without b)
118 TH1D * h_jet_pt;
119 TH1D * h_jet1_pt;
120 TH1D * h_jet2_pt;
121 TH1D * h_jet_num;
122
123 // Fancy stuff
124 TH1D * h_jet1Bjet1_eta;
125 TH1D * h_jetfromStop1;
126 TH1D * h_BjetfromStop1;
127 TH1D * h_Bjet1fromStop1;
128 TH1D * h_Bjet2fromStop1;
129 TH1D * h_jetsBfromStop1;
130 TH1D * h_jetsB1fromStop1;
131 TH1D * h_jetsB2fromStop1;
132 TH1D * h_jetsB3fromStop1;
133 TH1D * h_jetsB4fromStop1;
134
135 };
136
137 //
138 // constants, enums and typedefs
139 //
140
141 //
142 // static data member definitions
143 //
144
145 //
146 // constructors and destructor
147 //
148 GenAnalyzer::GenAnalyzer(const edm::ParameterSet& iConfig) :
149 src_( iConfig.getParameter<edm::InputTag>( "src" ) ), // Obtain inputs
150 stop1Mass_( iConfig.getParameter<double>( "stop1Mass" ) ),
151 stop2Mass_( iConfig.getParameter<double>( "stop2Mass" ) )
152 {
153 //now do what ever initialization is needed
154 edm::Service<TFileService> fs; // Output File
155
156 // Initialize Histograms
157 histo = fs->make<TH1D>("jetpt" , "Jet p_{T}" , 30 ,0., 300. ); // test histo
158
159 // Higgs
160 h_higgs_mass = fs->make<TH1D>("higgs_mass" , "Higgs Mass" , 12.5 , 99.0, 150. );
161 h_higgsTrue_mass = fs->make<TH1D>("higgsTrue_mass" , "True Higgs Mass" , 12.5, 99., 150. );
162 h_higgs_pt = fs->make<TH1D>("higgs_pt" , "Higgs p_{T}" , 98 , 0., 500. );
163 //h_higgsBB_deltaR = fs->make<TH1D>("higgs_deltaR" , "#Delta R (b,b) [from H]" , 50 , 0., 5. );
164 h_higgs_num = fs->make<TH1D>("higgs_num" , "Number of Higgs" , 10 , -0.5, 9.5 );
165
166 // Stop1
167 h_stop1_mass = fs->make<TH1D>("stop1_mass" , "Stop1 Mass" , 41 , 70., 530. );
168 h_stop1True_mass = fs->make<TH1D>("stop1True_mass" , "True Stop1 Mass" , 41 , 70, 530. );
169 h_stop1_pt = fs->make<TH1D>("stop1_pt" , "Stop1 p_{T}" , 98 , 0., 500. );
170 //h_stop1Bj_deltaR = fs->make<TH1D>("stop1Bj_deltaR" , "#Delta R (b,j) [from Stop1]" , 50 , 0., 5. );
171
172 // Stop2
173 h_stop2_mass = fs->make<TH1D>("stop2_mass" , "Stop2 Mass" , 31 , 220., 780. );
174 h_stop2True_mass = fs->make<TH1D>("stop2True_mass" , "True Stop2 Mass" , 31 , 220., 780. );
175
176 // Jets + B
177 h_jetsB_ht = fs->make<TH1D>("jetsB_ht" , "H_{T} (Jets + B)" , 40 , 200., 1000. );
178 h_jetsB_pt = fs->make<TH1D>("jetsB_pt" , "Jets + B p_{T}" , 98 , 0., 500. );
179 h_jetsB1_pt = fs->make<TH1D>("jetsB1_pt" , "Leading Jets + B p_{T}" , 98 , 0., 500. );
180 h_jetsB2_pt = fs->make<TH1D>("jetsB2_pt" , "2nd Leading Jets + B p_{T}" , 98 , 0., 500. );
181 h_jetsB3_pt = fs->make<TH1D>("jetsB3_pt" , "3rd Leading Jets + B p_{T}" , 98 , 0., 500. );
182 h_jetsB4_pt = fs->make<TH1D>("jetsB4_pt" , "4th Leading Jets + B p_{T}" , 98 , 0., 500. );
183 h_jetsB_num = fs->make<TH1D>("jetsB_num" , "Number of Jets + B" , 16 , -0.5, 15.5 );
184 h_jetsB1_eta = fs->make<TH1D>("jetsB1_eta" , "Leading Jets + B #eta" , 28 , -3.1, 3.1 );
185 h_jetsB2_eta = fs->make<TH1D>("jetsB2_eta" , "2nd Leading Jets + B #eta" , 28 , -3.1, 3.1 );
186 h_jetsB3_eta = fs->make<TH1D>("jetsB3_eta" , "3rd Leading Jets + B #eta" , 28 , -3.1, 3.1 );
187 h_jetsB4_eta = fs->make<TH1D>("jetsB4_eta" , "4th Leading Jets + B #eta" , 28 , -3.1, 3.1 );
188 h_jetsB1_phi = fs->make<TH1D>("jetsB1_phi" , "Leading Jets + B #phi" , 58 , -5.1, 5.1 );
189 h_jetsB2_phi = fs->make<TH1D>("jetsB2_phi" , "2nd Leading Jets + B #phi" , 58 , -5.1, 5.1 );
190 h_jetsB3_phi = fs->make<TH1D>("jetsB3_phi" , "3rd Leading Jets + B #phi" , 58 , -5.1, 5.1 );
191 h_jetsB4_phi = fs->make<TH1D>("jetsB4_phi" , "4th Leading Jets + B #phi" , 58 , -5.1, 5.1 );
192 h_jetsB1jetsB2_deltaR = fs->make<TH1D>("jetsB1jetsB2_deltaR" , "#Delta R (jb_{1},jb_{2})" , 50 , 0., 5. );
193 h_jetsB1jetsB3_deltaR = fs->make<TH1D>("jetsB1jetsB3_deltaR" , "#Delta R (jb_{1},jb_{3})" , 50 , 0., 5. );
194
195 // Bjets
196 h_Bjet1_pt = fs->make<TH1D>("Bjet1_pt" , "Leading B-Jet p_{T}" , 98 , 0., 500. );
197 h_Bjet2_pt = fs->make<TH1D>("Bjet2_pt" , "2nd Leading B-Jet p_{T}" , 98 , 0., 500. );
198 h_Bjet_num = fs->make<TH1D>("Bjet_num" , "Number of B-Jets" , 10 , -0.5, 9.5 );
199 //h_Bjet1Bjet2_deltaEta = fs->make<TH1D>("jet1Bjet1_deltaEta" , "#Delta #eta (bj_{1},bj_{2})" , 28 , -3.1, 3.1 );
200 h_Bjet1Bjet2_deltaPhi = fs->make<TH1D>("Bjet1Bjet2_deltaPhi" , "#Delta #phi (bj_{1}, bj_{2})" , 50 , 0., 5. );
201 h_Bjet1Bjet2_cosDeltaPhi = fs->make<TH1D>("Bjet1Bjet2_cosDeltaPhi" , "cos [#Delta #phi (bj_{1}, bj_{2})]" , 26 , -1.2, 1.2 );
202 h_Bjet1Bjet2_deltaR = fs->make<TH1D>("Bjet1Bjet2_deltaR" , "#Delta R (b_{1},b_{2})" , 50 , 0., 5. );
203
204 // Jets (without b)
205 h_jet_pt = fs->make<TH1D>("jet_pt" , "Jets p_{T}" , 98 , 0., 500. );
206 h_jet1_pt = fs->make<TH1D>("jet1_pt" , "Leading Jet p_{T}" , 98 , 0., 500. );
207 h_jet2_pt = fs->make<TH1D>("jet2_pt" , "2nd Leading Jet p_{T}" , 98 , 0., 500. );
208 h_jet_num = fs->make<TH1D>("jet_num" , "Number of Jets" , 10 , -0.5, 9.5 );
209
210 // Fancy stuff
211 //h_jet1Bjet1_eta = fs->make<TH1D>("jet1Bjet1_eta" , "#eta (j_{1},bj_{1})" , 28 , -3.1, 3.1 );
212 h_jetfromStop1 = fs->make<TH1D>("jetfromStop1" , "Position of Jets (w/o b) from Stop1" , 16, -0.5, 15.5 );
213 h_BjetfromStop1 = fs->make<TH1D>("BjetfromStop1" , "Position of B-Jets from Stop1" , 16, -0.5, 15.5 );
214 h_Bjet1fromStop1 = fs->make<TH1D>("Bjet1fromStop1" , "Position of Leading B-Jet from Stop1" , 16, -0.5, 15.5 );
215 h_Bjet2fromStop1 = fs->make<TH1D>("Bjet2fromStop1" , "Position of 2nd Leading B-Jet from Stop1" , 16, -0.5, 15.5 );
216 h_jetsBfromStop1 = fs->make<TH1D>("jetsBfromStop1" , "Position of Jets from Stop1" , 16, -0.5, 15.5 );
217 h_jetsB1fromStop1 = fs->make<TH1D>("jetsB1fromStop1" , "Position of Leading Jet from Stop1" , 16, -0.5, 15.5 );
218 h_jetsB2fromStop1 = fs->make<TH1D>("jetsB2fromStop1" , "Position of 2nd Leading Jet from Stop1" , 16, -0.5, 15.5 );
219 h_jetsB3fromStop1 = fs->make<TH1D>("jetsB3fromStop1" , "Position of 3rd Leading Jet from Stop1" , 16, -0.5, 15.5 );
220 h_jetsB4fromStop1 = fs->make<TH1D>("jetsB4fromStop1" , "Position of 4th Leading Jet from Stop1" , 16, -0.5, 15.5 );
221
222 }
223
224 // Bool to Sort by Pt
225 bool ComparePt(TLorentzVector a, TLorentzVector b) { return a.Pt() > b.Pt(); }
226
227 GenAnalyzer::~GenAnalyzer()
228 {
229
230 // do anything here that needs to be done at desctruction time
231 // (e.g. close files, deallocate resources etc.)
232
233 }
234
235
236 //
237 // member functions
238 //
239
240 // ------------ method called to for each event ------------
241 void GenAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup){
242
243 using namespace std;
244 using namespace edm;
245 using namespace reco;
246
247 // Way to call GenParticle from root
248 edm::Handle<std::vector<reco::GenParticle>> particles;
249 iEvent.getByLabel( src_ , particles );
250 const std::vector<reco::GenParticle> & p = *particles;
251
252 // Inizialize
253 size_t higgs_num = 0;
254
255 vector< TLorentzVector > p4jet;
256 vector< TLorentzVector > p4Bjet;
257 vector< TLorentzVector > p4jetsB;
258 vector< TLorentzVector > p4HiggsB;
259 vector< TLorentzVector > p4Higgs;
260 vector< TLorentzVector > p4Stop1B;
261 vector< TLorentzVector > p4Stop1jet;
262 vector< TLorentzVector > p4Stop1jetsB;
263 vector< TLorentzVector > p4Stop1;
264
265 // Begin Loop for GenParticles
266 for(unsigned int i = 0; i < particles->size(); ++ i) {
267
268 const Candidate * mom = p[i].mother(); // call mother particle
269 int idAbs = abs( p[i].pdgId() ); // call pdg Id
270 if( idAbs == 5 ) histo->Fill( p[i].pt() ); // test
271
272
273 //////////// True Mass plots
274 if( idAbs == 25 ){
275 higgs_num++; // Check number of Higgs
276 h_higgsTrue_mass->Fill(p[i].mass()); // Higgs true mass
277 }
278 if( idAbs == 1000002 ) h_stop1True_mass->Fill(p[i].mass()); // Stop1 true mass
279 if( idAbs == 1000004 ) h_stop2True_mass->Fill(p[i].mass()); // Stop2 true mass
280 ////////////
281
282 if( p[i].status()!= 3 ) continue; // Make sure only "final" particles are present
283
284 //////////// Storing TLorentz Vectors
285 TLorentzVector tmpjet, tmpBjet, tmpjetsB, tmpHiggs, tmpStop1B, tmpStop1jet; // tmp TLorentzVectors
286 if( mom ){
287 // Jets without b
288 if( ( idAbs == 1 ) | ( idAbs == 2 ) | ( idAbs == 3 ) | ( idAbs == 4) ){
289 tmpjet.SetPxPyPzE(p[i].px(),p[i].py(),p[i].pz(),p[i].energy());
290 p4jet.push_back( tmpjet );
291 }
292 std::sort(p4jet.begin(), p4jet.end(), ComparePt);
293 // Bjets
294 if( idAbs == 5 ){
295 tmpBjet.SetPxPyPzE(p[i].px(),p[i].py(),p[i].pz(),p[i].energy());
296 p4Bjet.push_back( tmpBjet );
297 }
298 std::sort(p4Bjet.begin(), p4Bjet.end(), ComparePt);
299 // Jets with b
300 if( ( idAbs == 1 ) | ( idAbs == 2 ) | ( idAbs == 3 ) | ( idAbs == 4) | ( idAbs == 5 )){
301 tmpBjet.SetPxPyPzE(p[i].px(),p[i].py(),p[i].pz(),p[i].energy());
302 p4jetsB.push_back( tmpBjet );
303 }
304 std::sort(p4jetsB.begin(), p4jetsB.end(), ComparePt);
305
306 // Higgs
307 if( ( idAbs == 5 ) & ( mom->pdgId() == 25 ) ){
308 tmpHiggs.SetPxPyPzE(p[i].px(),p[i].py(),p[i].pz(),p[i].energy());
309 p4HiggsB.push_back( tmpHiggs );
310 }
311 std::sort(p4HiggsB.begin(), p4HiggsB.end(), ComparePt);
312 // Stop1
313 if( abs( mom->pdgId() )== 1000002 ){
314 if( ( idAbs == 1 ) | ( idAbs == 2 ) | ( idAbs == 3 ) | ( idAbs == 4) ){
315 tmpStop1B.SetPxPyPzE(p[i].px(),p[i].py(),p[i].pz(),p[i].energy());
316 p4Stop1B.push_back( tmpStop1B );
317 p4Stop1jetsB.push_back( tmpStop1B );
318 }
319 else {
320 tmpStop1jet.SetPxPyPzE(p[i].px(),p[i].py(),p[i].pz(),p[i].energy());
321 p4Stop1jet.push_back( tmpStop1jet );
322 p4Stop1jetsB.push_back( tmpStop1jet );
323 }
324 }
325 std::sort(p4Stop1B.begin(), p4Stop1B.end(), ComparePt);
326 std::sort(p4Stop1jet.begin(), p4Stop1jet.end(), ComparePt);
327 std::sort(p4Stop1jetsB.begin(), p4Stop1jetsB.end(), ComparePt);
328 }
329 } // end GenParticles loop
330
331 // Test number of true Higgs
332 h_higgs_num->Fill( higgs_num );
333
334
335 ///// Higgs reconstructed mass
336 for(unsigned int j = 0; j < p4HiggsB.size()-1; ++j) {
337 for(unsigned int k=j+1; k < p4HiggsB.size(); ++k) {
338 if(j==k) continue; // avoid repetition
339 TLorentzVector p4CandidateHiggs = p4HiggsB[j] + p4HiggsB[k]; // Higgs TLorentzVector
340 //cout << p4HiggsB[j].E() << " " << p4HiggsB[k].E() << " " << p4CandidateHiggs.M() << endl;
341 Double_t massHiggs = 125.0; // nominal mass value
342 Double_t diffmassHiggs = p4CandidateHiggs.M() - massHiggs; // mass difference
343 //cout << diffmassHiggs << endl;
344 if( abs(diffmassHiggs) < 0.1 ){ // diffmass presicion
345 h_higgs_mass->Fill( p4CandidateHiggs.M() ); // Higgs mass
346 h_higgs_pt->Fill( p4CandidateHiggs.Pt() ); // Higgs pT
347 p4Higgs.push_back( p4CandidateHiggs ); // store final TLorentzVector
348 //cout << p4CandidateHiggs.M() << endl;
349 }
350 }
351 }
352 //////////////////////////////
353
354
355 /////// Stop1 reconstructed mass
356 for(unsigned int j = 0; j < p4Stop1B.size(); ++j) {
357 for(unsigned int k= 0; k < p4Stop1jet.size(); ++k) {
358 TLorentzVector p4CandidateStop1 = p4Stop1B[j] + p4Stop1jet[k]; // Higgs TLorentzVector
359 Double_t massStop1 = stop1Mass_ ; // nominal mass value
360 Double_t diffmassStop1 = p4CandidateStop1.M() - massStop1; // mass difference
361 if( abs(diffmassStop1) < 1.0 ){ // diffmass presicion
362 h_stop1_mass->Fill( p4CandidateStop1.M() ); // Stop1 mass
363 h_stop1_pt->Fill( p4CandidateStop1.Pt() ); // Stop1 pt
364 p4Stop1.push_back( p4CandidateStop1 );
365 }
366 }
367 }
368 //////////////////////////////////
369
370 // Stop2 reconstructed mass
371 for(unsigned int j = 0; j < p4Higgs.size(); ++j) {
372 for(unsigned int k= 0; k < p4Stop1.size(); ++k) {
373 TLorentzVector p4CandidateStop2 = p4Higgs[j] + p4Stop1[k]; // Higgs TLorentzVector
374 Double_t massStop2 = stop2Mass_; // nominal mass value
375 Double_t diffmassStop2 = p4CandidateStop2.M() - massStop2; // mass difference
376 if( abs(diffmassStop2) < 1.0 ){ // diffmass presicion
377 h_stop2_mass->Fill( p4CandidateStop2.M() ); // Stop1 mass
378 }
379 }
380 }
381
382 // Jets + B
383 Double_t Ht = 0;
384 for(unsigned int j = 0; j < p4jetsB.size(); ++j) {
385
386 Ht += p4jetsB[j].Pt();
387
388 h_jetsB_pt->Fill( p4jetsB[j].Pt() );
389 h_jetsB1_pt->Fill( p4jetsB[0].Pt() );
390 h_jetsB2_pt->Fill( p4jetsB[1].Pt() );
391 h_jetsB3_pt->Fill( p4jetsB[2].Pt() );
392 h_jetsB4_pt->Fill( p4jetsB[3].Pt() );
393 h_jetsB_num->Fill( p4jetsB.size() );
394 h_jetsB1_eta->Fill( p4jetsB[0].Eta() );
395 h_jetsB2_eta->Fill( p4jetsB[1].Eta() );
396 h_jetsB3_eta->Fill( p4jetsB[2].Eta() );
397 h_jetsB4_eta->Fill( p4jetsB[3].Eta() );
398 h_jetsB1_phi->Fill( p4jetsB[0].Phi() );
399 h_jetsB2_phi->Fill( p4jetsB[1].Phi() );
400 h_jetsB3_phi->Fill( p4jetsB[2].Phi() );
401 h_jetsB4_phi->Fill( p4jetsB[3].Phi() );
402 h_jetsB1jetsB2_deltaR->Fill( p4jetsB[0].DeltaR( p4jetsB[1] ) );
403 h_jetsB1jetsB3_deltaR->Fill( p4jetsB[0].DeltaR( p4jetsB[2] ) );
404
405 for(unsigned int k= 0; k < p4Stop1jet.size(); ++k) {
406 if( p4Stop1jet[k].Pt() == p4jetsB[j].Pt() ) h_jetfromStop1->Fill( j+1 );
407 //cout << k << " "<< j << " " <<p4Stop1jet[k].Pt() << " " << p4jetsB[j].Pt() << endl;
408 }
409 if( p4Stop1B[0].Pt() == p4jetsB[j].Pt() ) h_Bjet1fromStop1->Fill( j+1 );
410 if( p4Stop1B[1].Pt() == p4jetsB[j].Pt() ) h_Bjet2fromStop1->Fill( j+1 );
411 for(unsigned int kk= 0; kk < p4Stop1B.size(); ++kk) {
412 if( p4Stop1B[kk].Pt() == p4jetsB[j].Pt() ) h_BjetfromStop1->Fill( j+1 );
413 //cout << kk << " "<< j << " " <<p4Stop1B[kk].Pt() << " " << p4jetsB[j].Pt() << endl;
414 }
415 if( p4Stop1jetsB[0].Pt() == p4jetsB[j].Pt() ) h_jetsB1fromStop1->Fill( j+1 );
416 if( p4Stop1jetsB[1].Pt() == p4jetsB[j].Pt() ) h_jetsB2fromStop1->Fill( j+1 );
417 if( p4Stop1jetsB[2].Pt() == p4jetsB[j].Pt() ) h_jetsB3fromStop1->Fill( j+1 );
418 if( p4Stop1jetsB[3].Pt() == p4jetsB[j].Pt() ) h_jetsB4fromStop1->Fill( j+1 );
419 for(unsigned int kkk= 0; kkk < p4Stop1jetsB.size(); ++kkk) {
420 if( p4Stop1jetsB[kkk].Pt() == p4jetsB[j].Pt() ) h_jetsBfromStop1->Fill( j+1 );
421 //cout << kkk << " "<< j << " " <<p4Stop1jetsB[kkk].Pt() << " " << p4jetsB[j].Pt() << endl;
422 }
423 }
424 h_jetsB_ht->Fill( Ht );
425 /////////////////////////////
426
427 // Bjets
428 for(unsigned int j = 0; j < p4Bjet.size(); ++j) {
429 h_Bjet1_pt->Fill( p4Bjet[0].Pt() );
430 h_Bjet2_pt->Fill( p4Bjet[1].Pt() );
431 h_Bjet_num->Fill( p4Bjet.size() );
432 h_Bjet1Bjet2_deltaPhi->Fill( p4Bjet[0].DeltaPhi( p4Bjet[1] ) );
433 h_Bjet1Bjet2_cosDeltaPhi->Fill( cos( p4Bjet[0].DeltaPhi( p4Bjet[1] ) ) );
434 h_Bjet1Bjet2_deltaR->Fill( p4Bjet[0].DeltaR( p4Bjet[1] ) );
435 }
436 /////////////////////////////
437
438 // Jets
439 for(unsigned int j = 0; j < p4jet.size(); ++j) {
440 h_jet_pt->Fill( p4jet[j].Pt() );
441 h_jet1_pt->Fill( p4jet[0].Pt() );
442 h_jet2_pt->Fill( p4jet[1].Pt() );
443 h_jet_num->Fill( p4jet.size() );
444 }
445 /////////////////////////////
446
447 }
448
449 // ------------ method called once each job just after ending the event loop ------------
450 void
451 GenAnalyzer::beginJob() {
452 }
453 void
454 GenAnalyzer::endJob() {
455 }
456
457 //define this as a plug-in
458 DEFINE_FWK_MODULE(GenAnalyzer);