238 std::cout <<
"Usage : " <<
argv[0] <<
" [parameters.py]" << std::endl;
242 cout <<
"Hello from " <<
argv[0] <<
"!" << endl;
246 gSystem->Load(
"libFWCoreFWLite" );
250 cout <<
"Getting parameters" << endl;
253 std::shared_ptr<edm::ProcessDesc>
b = builder.processDesc();
254 std::shared_ptr<edm::ParameterSet>
parameters = b->getProcessPSet();
255 parameters->registerIt();
265 cout <<
"Making RunLumiSelector" << endl;
268 cout <<
"setting up TFileService" << endl;
273 cout <<
"Setting up chain event" << endl;
279 cout <<
"Booking histograms" << endl;
282 std::map<std::string, TH1*>
hists;
284 hists[
"hist_nJet" ] = theDir.
make<TH1D>(
"hist_nJet" ,
"Number of Calo Jets", 10, 0, 10 ) ;
285 hists[
"hist_nPFJet" ] = theDir.
make<TH1D>(
"hist_nPFJet" ,
"Number of PF Jets", 10, 0, 10 ) ;
287 hists[
"hist_jetPt" ] = theDir.
make<TH1D>(
"hist_jetPt" ,
"Jet p_{T}", 400, 0, 400 ) ;
288 hists[
"hist_jetEtaVsPhi" ] = theDir.
make<TH2D>(
"hist_jetEtaVsPhi" ,
"Jet #phi versus #eta;#eta;#phi", 50, -5.0, 5.0, 50, -
TMath::Pi(),
TMath::Pi() ) ;
289 hists[
"hist_jetNTracks" ] = theDir.
make<TH1D>(
"hist_jetNTracks" ,
"Jet N_{TRACKS}", 20, 0, 20 ) ;
290 hists[
"hist_jetNTracksVsPt" ] = theDir.
make<TH2D>(
"hist_jetNTracksVsPt" ,
"Number of Tracks versus Jet p_{T};Jet p_{T}(GeV/c) ;N_{Tracks}",20, 0, 200, 20, 0, 20 ) ;
291 hists[
"hist_jetEMF" ] = theDir.
make<TH1D>(
"hist_jetEMF" ,
"Jet EMF", 200, 0, 1) ;
292 hists[
"hist_jetPdgID" ] = theDir.
make<TH1D>(
"hist_jetPdgID" ,
"PDG Id of Jet Constituents", 10000, 0, 10000 ) ;
293 hists[
"hist_jetGenEmE" ] = theDir.
make<TH1D>(
"hist_jetGenEmE" ,
"Gen Jet EM Energy", 200, 0, 200 ) ;
294 hists[
"hist_jetGenHadE" ] = theDir.
make<TH1D>(
"hist_jetGenHadE" ,
"Gen Jet HAD Energy", 200, 0, 200 ) ;
295 hists[
"hist_jetGenEMF" ] = theDir.
make<TH1D>(
"hist_jetGenEMF" ,
"Gen Jet EMF", 200, 0, 1) ;
296 hists[
"hist_jetEoverGenE" ] = theDir.
make<TH1D>(
"hist_jetEoverGenE" ,
"Energy of reco Jet / Energy of gen Jet", 200, 0, 2.0) ;
297 hists[
"hist_jetCorr" ] = theDir.
make<TH1D>(
"hist_jetCorr" ,
"Jet Correction Factor", 200, 0, 1.0 ) ;
298 hists[
"hist_n90Hits" ] = theDir.
make<TH1D>(
"hist_n90Hits" ,
"Jet n90Hits", 20, 0, 20) ;
299 hists[
"hist_fHPD" ] = theDir.
make<TH1D>(
"hist_fHPD" ,
"Jet fHPD", 200, 0, 1) ;
300 hists[
"hist_nConstituents" ] = theDir.
make<TH1D>(
"hist_nConstituents" ,
"Jet nConstituents", 20, 0, 20 ) ;
301 hists[
"hist_jetCHF" ] = theDir.
make<TH1D>(
"hist_jetCHF" ,
"Jet Charged Hadron Fraction", 200, 0, 1.0) ;
303 hists[
"hist_good_jetPt" ] = theDir.make<TH1D>(
"hist_good_jetPt" ,
"Jet p_{T}", 400, 0, 400 ) ;
304 hists[
"hist_good_jetEtaVsPhi" ] = theDir.make<TH2D>(
"hist_good_jetEtaVsPhi" ,
"Jet #phi versus #eta;#eta;#phi", 50, -5.0, 5.0, 50, -
TMath::Pi(),
TMath::Pi() ) ;
305 hists[
"hist_good_jetNTracks" ] = theDir.make<TH1D>(
"hist_good_jetNTracks" ,
"Jet N_{TRACKS}", 20, 0, 20 ) ;
306 hists[
"hist_good_jetNTracksVsPt" ] = theDir.make<TH2D>(
"hist_good_jetNTracksVsPt" ,
"Number of Tracks versus Jet p_{T};Jet p_{T}(GeV/c) ;N_{Tracks}",20, 0, 200, 20, 0, 20 ) ;
307 hists[
"hist_good_jetEMF" ] = theDir.make<TH1D>(
"hist_good_jetEMF" ,
"Jet EMF", 200, 0, 1) ;
308 hists[
"hist_good_jetPdgID" ] = theDir.make<TH1D>(
"hist_good_jetPdgID" ,
"PDG Id of Jet Constituents", 10000, 0, 10000 ) ;
309 hists[
"hist_good_jetGenEmE" ] = theDir.make<TH1D>(
"hist_good_jetGenEmE" ,
"Gen Jet EM Energy", 200, 0, 200 ) ;
310 hists[
"hist_good_jetGenHadE" ] = theDir.make<TH1D>(
"hist_good_jetGenHadE" ,
"Gen Jet HAD Energy", 200, 0, 200 ) ;
311 hists[
"hist_good_jetGenEMF" ] = theDir.make<TH1D>(
"hist_good_jetGenEMF" ,
"Gen Jet EMF", 200, 0, 1) ;
312 hists[
"hist_good_jetEoverGenE" ] = theDir.make<TH1D>(
"hist_good_jetEoverGenE" ,
"Energy of reco Jet / Energy of gen Jet", 200, 0, 2.0) ;
313 hists[
"hist_good_jetCorr" ] = theDir.make<TH1D>(
"hist_good_jetCorr" ,
"Jet Correction Factor", 200, 0, 1.0 ) ;
314 hists[
"hist_good_n90Hits" ] = theDir.make<TH1D>(
"hist_good_n90Hits" ,
"Jet n90Hits", 20, 0, 20) ;
315 hists[
"hist_good_fHPD" ] = theDir.make<TH1D>(
"hist_good_fHPD" ,
"Jet fHPD", 200, 0, 1) ;
316 hists[
"hist_good_nConstituents" ] = theDir.make<TH1D>(
"hist_good_nConstituents" ,
"Jet nConstituents", 20, 0, 20 ) ;
317 hists[
"hist_good_jetCHF" ] = theDir.make<TH1D>(
"hist_good_jetCHF" ,
"Jet Charged Hadron Fraction", 200, 0, 1.0) ;
320 hists[
"hist_pf_jetPt" ] = theDir.make<TH1D>(
"hist_pf_jetPt" ,
"PFJet p_{T}", 400, 0, 400 ) ;
321 hists[
"hist_pf_jetEtaVsPhi" ] = theDir.make<TH2D>(
"hist_pf_jetEtaVsPhi" ,
"PFJet #phi versus #eta;#eta;#phi", 50, -5.0, 5.0, 50, -
TMath::Pi(),
TMath::Pi() ) ;
322 hists[
"hist_pf_jetNTracks" ] = theDir.make<TH1D>(
"hist_pf_jetNTracks" ,
"PFJet N_{TRACKS}", 20, 0, 20 ) ;
323 hists[
"hist_pf_jetNTracksVsPt" ] = theDir.make<TH2D>(
"hist_pf_jetNTracksVsPt" ,
"Number of Tracks versus Jet p_{T};Jet p_{T}(GeV/c) ;N_{Tracks}",20, 0, 200, 20, 0, 20 ) ;
324 hists[
"hist_pf_jetEMF" ] = theDir.make<TH1D>(
"hist_pf_jetEMF" ,
"PFJet EMF", 200, 0, 1) ;
325 hists[
"hist_pf_jetCHF" ] = theDir.make<TH1D>(
"hist_pf_jetCHF" ,
"PFJet CHF", 200, 0, 1) ;
326 hists[
"hist_pf_jetNHF" ] = theDir.make<TH1D>(
"hist_pf_jetNHF" ,
"PFJet NHF", 200, 0, 1) ;
327 hists[
"hist_pf_jetCEF" ] = theDir.make<TH1D>(
"hist_pf_jetCEF" ,
"PFJet CEF", 200, 0, 1) ;
328 hists[
"hist_pf_jetNEF" ] = theDir.make<TH1D>(
"hist_pf_jetNEF" ,
"PFJet NEF", 200, 0, 1) ;
329 hists[
"hist_pf_jetPdgID" ] = theDir.make<TH1D>(
"hist_pf_jetPdgID" ,
"PDG Id of Jet Constituents", 10000, 0, 10000 ) ;
330 hists[
"hist_pf_jetGenEmE" ] = theDir.make<TH1D>(
"hist_pf_jetGenEmE" ,
"Gen Jet EM Energy", 200, 0, 200 ) ;
331 hists[
"hist_pf_jetGenHadE" ] = theDir.make<TH1D>(
"hist_pf_jetGenHadE" ,
"Gen Jet HAD Energy", 200, 0, 200 ) ;
332 hists[
"hist_pf_jetGenEMF" ] = theDir.make<TH1D>(
"hist_pf_jetGenEMF" ,
"Gen Jet EMF", 200, 0, 1) ;
333 hists[
"hist_pf_jetEoverGenE" ] = theDir.make<TH1D>(
"hist_pf_jetEoverGenE" ,
"Energy of reco Jet / Energy of gen Jet", 200, 0, 2.0) ;
334 hists[
"hist_pf_jetCorr" ] = theDir.make<TH1D>(
"hist_pf_jetCorr" ,
"PFJet Correction Factor", 200, 0, 1.0 ) ;
335 hists[
"hist_pf_nConstituents" ] = theDir.make<TH1D>(
"hist_pf_nConstituents" ,
"PFJet nConstituents", 20, 0, 20 ) ;
337 hists[
"hist_pf_good_jetPt" ] = theDir.make<TH1D>(
"hist_pf_good_jetPt" ,
"PFJet p_{T}", 400, 0, 400 ) ;
338 hists[
"hist_pf_good_jetEtaVsPhi" ] = theDir.make<TH2D>(
"hist_pf_good_jetEtaVsPhi" ,
"PFJet #phi versus #eta;#eta;#phi", 50, -5.0, 5.0, 50, -
TMath::Pi(),
TMath::Pi() ) ;
339 hists[
"hist_pf_good_jetNTracks" ] = theDir.make<TH1D>(
"hist_pf_good_jetNTracks" ,
"PFJet N_{TRACKS}", 20, 0, 20 ) ;
340 hists[
"hist_pf_good_jetNTracksVsPt" ] = theDir.make<TH2D>(
"hist_pf_good_jetNTracksVsPt",
"Number of Tracks versus Jet p_{T};Jet p_{T}(GeV/c) ;N_{Tracks}",20, 0, 200, 20, 0, 20 ) ;
341 hists[
"hist_pf_good_jetEMF" ] = theDir.make<TH1D>(
"hist_pf_good_jetEMF" ,
"PFJet EMF", 200, 0, 1) ;
342 hists[
"hist_pf_good_jetCHF" ] = theDir.make<TH1D>(
"hist_pf_good_jetCHF" ,
"PFJet CHF", 200, 0, 1) ;
343 hists[
"hist_pf_good_jetNHF" ] = theDir.make<TH1D>(
"hist_pf_good_jetNHF" ,
"PFJet NHF", 200, 0, 1) ;
344 hists[
"hist_pf_good_jetCEF" ] = theDir.make<TH1D>(
"hist_pf_good_jetCEF" ,
"PFJet CEF", 200, 0, 1) ;
345 hists[
"hist_pf_good_jetNEF" ] = theDir.make<TH1D>(
"hist_pf_good_jetNEF" ,
"PFJet NEF", 200, 0, 1) ;
346 hists[
"hist_pf_good_jetPdgID" ] = theDir.make<TH1D>(
"hist_pf_good_jetPdgID" ,
"PDG Id of Jet Constituents", 10000, 0, 10000 ) ;
347 hists[
"hist_pf_good_jetGenEmE" ] = theDir.make<TH1D>(
"hist_pf_good_jetGenEmE" ,
"Gen Jet EM Energy", 200, 0, 200 ) ;
348 hists[
"hist_pf_good_jetGenHadE" ] = theDir.make<TH1D>(
"hist_pf_good_jetGenHadE" ,
"Gen Jet HAD Energy", 200, 0, 200 ) ;
349 hists[
"hist_pf_good_jetGenEMF" ] = theDir.make<TH1D>(
"hist_pf_good_jetGenEMF" ,
"Gen Jet EMF", 200, 0, 1) ;
350 hists[
"hist_pf_good_jetEoverGenE" ] = theDir.make<TH1D>(
"hist_pf_good_jetEoverGenE" ,
"Energy of reco Jet / Energy of gen Jet", 200, 0, 2.0) ;
351 hists[
"hist_pf_good_jetCorr" ] = theDir.make<TH1D>(
"hist_pf_good_jetCorr" ,
"PFJet Correction Factor", 200, 0, 1.0 ) ;
352 hists[
"hist_pf_good_nConstituents" ] = theDir.make<TH1D>(
"hist_pf_good_nConstituents" ,
"PFJet nConstituents", 20, 0, 20 ) ;
354 hists[
"hist_mjj" ] = theDir.make<TH1D>(
"hist_mjj" ,
"Dijet mass", 300, 0, 300 ) ;
355 hists[
"hist_dR_jj" ] = theDir.make<TH1D>(
"hist_dR_jj" ,
"#Delta R_{JJ}", 200, 0,
TMath::TwoPi() ) ;
356 hists[
"hist_imbalance_jj" ] = theDir.make<TH1D>(
"hist_imbalance_jj" ,
"Dijet imbalance", 200, -1.0, 1.0 ) ;
358 hists[
"hist_pf_mjj" ] = theDir.make<TH1D>(
"hist_pf_mjj" ,
"Dijet mass", 300, 0, 300 ) ;
359 hists[
"hist_pf_dR_jj" ] = theDir.make<TH1D>(
"hist_pf_dR_jj" ,
"#Delta R_{JJ}", 200, 0,
TMath::TwoPi() ) ;
360 hists[
"hist_pf_imbalance_jj" ] = theDir.make<TH1D>(
"hist_pf_imbalance_jj" ,
"Dijet imbalance", 200, -1.0, 1.0 ) ;
364 cout <<
"Making functors" << endl;
371 pfJetStudiesParams );
376 cout <<
"About to begin looping" << endl;
380 for (
ev.toBegin(); !
ev.atEnd(); ++
ev, ++nev) {
385 if ( runLumiSel(
ev) ==
false )
continue;
387 if ( nev % 100 == 0 )
cout <<
"Processing run " <<
event.id().run() <<
", lumi " <<
event.id().luminosityBlock() <<
", event " <<
event.id().event() << endl;
393 caloSelector(
event, retCalo );
397 pfSelector(
event, retPF );
402 if ( retCalo.
test( caloSelector.caloKin() ) ) {
403 if ( retCalo.
test( caloSelector.caloDeltaPhi() ) ) {
404 vector<pat::Jet>
const & allCaloJets = caloSelector.allCaloJets();
406 for ( std::vector<pat::Jet>::const_iterator jetBegin = allCaloJets.begin(),
407 jetEnd = jetBegin + 2, ijet = jetBegin;
408 ijet != jetEnd; ++ijet ) {
412 double pt = jet.
pt();
416 hists[
"hist_jetPt"]->Fill( pt );
417 hists[
"hist_jetEtaVsPhi"]->Fill( jet.
eta(), jet.
phi() );
418 hists[
"hist_jetNTracks"]->Fill( jetTracks.
size() );
419 hists[
"hist_jetNTracksVsPt"]->Fill( pt, jetTracks.
size() );
421 hists[
"hist_jetCorr"]->Fill( jet.
jecFactor(
"Uncorrected") );
422 hists[
"hist_n90Hits"]->Fill( static_cast<int>(jet.
jetID().
n90Hits) );
423 hists[
"hist_fHPD"]->Fill( jet.
jetID().
fHPD );
426 if ( useMC && jet.
genJet() != 0 ) {
433 TLorentzVector p4_tracks(0,0,0,0);
435 itrkEnd = jetTracks.
end();
436 itrk != itrkEnd; ++itrk ) {
437 TLorentzVector p4_trk;
438 double M_PION = 0.140;
439 p4_trk.SetPtEtaPhiM( (*itrk)->pt(), (*itrk)->eta(), (*itrk)->phi(), M_PION );
442 hists[
"hist_jetCHF"]->Fill( p4_tracks.Energy() / jet.
energy() );
450 if ( retCalo.
test( caloSelector.caloJetID() ) ) {
451 pat::Jet const & jet0 = caloSelector.caloJet0();
452 pat::Jet const & jet1 = caloSelector.caloJet1();
454 TLorentzVector p4_j0( jet0.
px(), jet0.
py(), jet0.
pz(), jet0.
energy() );
455 TLorentzVector p4_j1( jet1.
px(), jet1.
py(), jet1.
pz(), jet1.
energy() );
457 TLorentzVector p4_jj = p4_j0 + p4_j1;
459 hists[
"hist_mjj"]->Fill( p4_jj.M() );
460 hists[
"hist_dR_jj"]->Fill( p4_j0.DeltaR( p4_j1 ) );
461 hists[
"hist_imbalance_jj"]->Fill( (p4_j0.Perp() - p4_j1.Perp() ) /
462 (p4_j0.Perp() + p4_j1.Perp() ) );
464 hists[
"hist_good_jetPt"]->Fill( jet0.
pt() );
465 hists[
"hist_good_jetEtaVsPhi"]->Fill( jet0.
eta(), jet0.
phi() );
469 hists[
"hist_good_jetCorr"]->Fill( jet0.
jecFactor(
"Uncorrected") );
470 hists[
"hist_good_n90Hits"]->Fill( static_cast<int>(jet0.
jetID().
n90Hits) );
471 hists[
"hist_good_fHPD"]->Fill( jet0.
jetID().
fHPD );
472 hists[
"hist_good_nConstituents"]->Fill( jet0.
nConstituents() );
475 hists[
"hist_good_jetPt"]->Fill( jet1.
pt() );
476 hists[
"hist_good_jetEtaVsPhi"]->Fill( jet1.
eta(), jet1.
phi() );
480 hists[
"hist_good_jetCorr"]->Fill( jet1.
jecFactor(
"Uncorrected") );
481 hists[
"hist_good_n90Hits"]->Fill( static_cast<int>(jet1.
jetID().
n90Hits) );
482 hists[
"hist_good_fHPD"]->Fill( jet1.
jetID().
fHPD );
483 hists[
"hist_good_nConstituents"]->Fill( jet1.
nConstituents() );
493 if ( retPF.
test( pfSelector.pfDeltaPhi() ) ) {
495 vector<pat::Jet>
const & allPFJets = pfSelector.allPFJets();
497 for ( std::vector<pat::Jet>::const_iterator jetBegin = allPFJets.begin(),
498 jetEnd = jetBegin + 2, ijet = jetBegin;
499 ijet != jetEnd; ++ijet ) {
503 double pt = jet.
pt();
505 hists[
"hist_pf_jetPt"]->Fill( pt );
506 hists[
"hist_pf_jetEtaVsPhi"]->Fill( jet.
eta(), jet.
phi() );
514 if ( useMC && jet.
genJet() != 0 ) {
526 if ( retPF.
test( pfSelector.pfJetID() ) ) {
527 pat::Jet const & jet0 = pfSelector.pfJet0();
528 pat::Jet const & jet1 = pfSelector.pfJet1();
530 TLorentzVector p4_j0( jet0.
px(), jet0.
py(), jet0.
pz(), jet0.
energy() );
531 TLorentzVector p4_j1( jet1.
px(), jet1.
py(), jet1.
pz(), jet1.
energy() );
533 TLorentzVector p4_jj = p4_j0 + p4_j1;
535 hists[
"hist_pf_mjj"]->Fill( p4_jj.M() );
536 hists[
"hist_pf_dR_jj"]->Fill( p4_j0.DeltaR( p4_j1 ) );
537 hists[
"hist_pf_imbalance_jj"]->Fill( (p4_j0.Perp() - p4_j1.Perp() ) /
538 (p4_j0.Perp() + p4_j1.Perp() ) );
540 hists[
"hist_pf_good_jetPt"]->Fill( jet0.
pt() );
541 hists[
"hist_pf_good_jetEtaVsPhi"]->Fill( jet0.
eta(), jet0.
phi() );
542 hists[
"hist_pf_good_nConstituents"]->Fill( jet0.
nConstituents() );
549 hists[
"hist_pf_good_jetPt"]->Fill( jet1.
pt() );
550 hists[
"hist_pf_good_jetEtaVsPhi"]->Fill( jet1.
eta(), jet1.
phi() );
551 hists[
"hist_pf_good_nConstituents"]->Fill( jet1.
nConstituents() );
563 cout <<
"Calo jet selection" << endl;
565 cout <<
"PF jet selection" << endl;
T getParameter(std::string const &) const
virtual double pt() const final
transverse momentum
float neutralHadronEnergyFraction() const
neutralHadronEnergyFraction (relative to uncorrected jet energy)
virtual double eta() const final
momentum pseudorapidity
float chargedHadronEnergyFraction() const
chargedHadronEnergyFraction (relative to uncorrected jet energy)
const reco::TrackRefVector & associatedTracks() const
method to return a vector of refs to the tracks associated to this jet
float chargedEmEnergyFraction() const
chargedEmEnergyFraction (relative to uncorrected jet energy)
const_iterator end() const
Termination of iteration.
const_iterator begin() const
Initialize an iterator over the RefVector.
virtual double phi() const final
momentum azimuthal angle
virtual double px() const final
x coordinate of momentum vector
virtual double energy() const final
energy
static void enable()
enable automatic library loading
T * make(const Args &...args) const
make new ROOT object
float emEnergyFraction() const
returns the jet electromagnetic energy fraction
virtual double pz() const final
z coordinate of momentum vector
reco::JetID const & jetID() const
accessing Jet ID information
TFileDirectory mkdir(const std::string &dir, const std::string &descr="")
create a new subdirectory
Analysis-level calorimeter jet class.
bool test(std::string s) const
test
virtual int nConstituents() const
of constituents
const reco::GenJet * genJet() const
return the matched generated jet
size_type size() const
Size of the RefVector.
float jecFactor(const std::string &level, const std::string &flavor="none", const std::string &set="") const
float neutralEmEnergyFraction() const
neutralEmEnergyFraction (relative to uncorrected jet energy)
virtual double py() const final
y coordinate of momentum vector