44 #include "TLorentzVector.h"
49 using namespace l1extra;
50 using namespace trigger;
57 virtual void beginRun(
const Run & iRun,
const EventSetup & iSetup);
59 virtual void endJob();
61 void match(
MonitorElement * me, vector<LeafCandidate> & from, vector<LeafCandidate> &
to,
double deltaRMatchingCut, vector<int> &
map );
62 void myBook2D( TString
name, vector<double> &xBins, TString xLabel, vector<double> &yBins, TString yLabel, TString
title);
63 void myBook2D( TString name, vector<double> &xBins, TString xLabel, vector<double> &yBins, TString yLabel){
64 myBook2D( name, xBins, xLabel, yBins, yLabel, name);
66 void myBookProfile2D( TString name, vector<double> &xBins, TString xLabel, vector<double> &yBins, TString yLabel, TString title);
67 void myBookProfile2D( TString name, vector<double> &xBins, TString xLabel, vector<double> &yBins, TString yLabel){
68 myBookProfile2D( name, xBins, xLabel, yBins, yLabel, name);
70 void myBook1D( TString name, vector<double> &xBins, TString
label, TString title );
72 myBook1D( name, xBins, label, name );
96 map<TString, MonitorElement *>
ME;
103 dqmFolder(pset.getUntrackedParameter<
string>(
"DQMFolder")),
104 triggerProcessName(pset.getUntrackedParameter<
string>(
"TriggerProcessName")),
105 triggerPathName(pset.getUntrackedParameter<
string>(
"TriggerPathName")),
106 triggerSummaryRAWTag(
InputTag( pset.getUntrackedParameter<
string>(
"TriggerSummaryRAW"),
"", triggerProcessName)),
107 triggerSummaryAODTag(
InputTag( pset.getUntrackedParameter<
string>(
"TriggerSummaryAOD"),
"", triggerProcessName)),
108 triggerResultsTag(
InputTag( pset.getUntrackedParameter<
string>(
"TriggerResults"),
"", triggerProcessName)),
109 recoMuonsTag(pset.getParameter<
InputTag>(
"RecoMuons")),
110 genParticlesTag(pset.getParameter<
InputTag>(
"GenParticles")),
111 motherIDs(pset.getUntrackedParameter<vector<int> >(
"MotherIDs")),
112 genGlobDeltaRMatchingCut(pset.getUntrackedParameter<double>(
"GenGlobDeltaRMatchingCut")),
113 globL1DeltaRMatchingCut(pset.getUntrackedParameter<double>(
"GlobL1DeltaRMatchingCut")),
114 globL2DeltaRMatchingCut(pset.getUntrackedParameter<double>(
"GlobL2DeltaRMatchingCut")),
115 globL3DeltaRMatchingCut(pset.getUntrackedParameter<double>(
"GlobL3DeltaRMatchingCut")),
116 deltaEtaBins(pset.getUntrackedParameter<vector<double> >(
"DeltaEtaBins")),
117 deltaPhiBins(pset.getUntrackedParameter<vector<double> >(
"DeltaPhiBins")),
118 muonPtBins(pset.getUntrackedParameter<vector<double> >(
"MuonPtBins")),
119 muonEtaBins(pset.getUntrackedParameter<vector<double> >(
"MuonEtaBins")),
120 muonPhiBins(pset.getUntrackedParameter<vector<double> >(
"MuonPhiBins")),
121 dimuonPtBins(pset.getUntrackedParameter<vector<double> >(
"DimuonPtBins")),
122 dimuonEtaBins(pset.getUntrackedParameter<vector<double> >(
"DimuonEtaBins")),
123 dimuonDRBins(pset.getUntrackedParameter<vector<double> >(
"DimuonDRBins")),
132 LogDebug(
"HLTriggerOfflineHeavyFlavor") <<
"Successfully initialized HLTConfigProvider with process name: "<<
triggerProcessName<<endl;
139 for(
size_t i = 0;
i < triggerNames.size();
i++) {
140 TString triggerName = triggerNames[
i];
142 vector<string> moduleNames = hltConfig.
moduleLabels( triggerNames[
i] );
143 for(
size_t j = 0;
j < moduleNames.size();
j++) {
144 TString
name = moduleNames[
j];
145 if(name.Contains(
"Filter")){
147 if(name.Contains(
"L1"))
149 else if(name.Contains(
"L2"))
151 else if(name.Contains(
"L3"))
153 else if(name.Contains(
"mumuFilter") || name.Contains(
"JpsiTrackMass"))
156 os<<
" "<<moduleNames[
j];
172 LogError(
"HLTriggerOfflineHeavyFlavor") <<
"Could not find DQMStore service\n";
268 vector<double> sizeBins; sizeBins.push_back(10); sizeBins.push_back(0); sizeBins.push_back(10);
269 myBook1D(
"genMuon_size", sizeBins,
"container size" );
270 myBook1D(
"globMuon_size", sizeBins,
"container size" );
279 LogDebug(
"HLTriggerOfflineHeavyFlavor")<<
"Could not access DQM Store service"<<endl;
286 vector<LeafCandidate> genMuons;
290 for(GenParticleCollection::const_iterator
p=genParticles->begin();
p!= genParticles->end(); ++
p){
291 if(
p->status() == 1 &&
std::abs(
p->pdgId())==13 &&
293 genMuons.push_back( *
p );
297 LogDebug(
"HLTriggerOfflineHeavyFlavor")<<
"Could not access GenParticleCollection"<<endl;
300 ME[
"genMuon_size"]->Fill(genMuons.size());
301 LogDebug(
"HLTriggerOfflineHeavyFlavor")<<
"GenParticleCollection from "<<
genParticlesTag<<
" has size: "<<genMuons.size()<<endl;
303 vector<LeafCandidate> globMuons;
304 vector<LeafCandidate> globMuons_position;
308 for(MuonCollection::const_iterator
p=recoMuonsHandle->begin();
p!= recoMuonsHandle->end(); ++
p){
309 if(
p->isGlobalMuon()){
310 globMuons.push_back(*
p);
311 globMuons_position.push_back(
LeafCandidate(
p->charge(),
math::XYZTLorentzVector(
p->outerTrack()->innerPosition().x(),
p->outerTrack()->innerPosition().y(),
p->outerTrack()->innerPosition().z(), 0.) ) );
315 LogDebug(
"HLTriggerOfflineHeavyFlavor")<<
"Could not access reco Muons"<<endl;
317 ME[
"globMuon_size"]->Fill(globMuons.size());
318 LogDebug(
"HLTriggerOfflineHeavyFlavor")<<
"Global Muons from "<<
recoMuonsTag<<
" has size: "<<globMuons.size()<<endl;
321 vector<vector<LeafCandidate> > muonsAtFilter;
322 vector<vector<LeafCandidate> > muonPositionsAtFilter;
324 muonsAtFilter.push_back(vector<LeafCandidate>());
325 muonPositionsAtFilter.push_back(vector<LeafCandidate>());
329 if( rawTriggerEvent.
isValid() ){
332 if ( index < rawTriggerEvent->
size() ){
334 vector<L1MuonParticleRef> l1Cands;
335 rawTriggerEvent->getObjects( index,
TriggerL1Mu, l1Cands );
336 for(
size_t j=0;
j<l1Cands.size();
j++){
337 muonsAtFilter[
i].push_back(*l1Cands[
j]);
340 vector<RecoChargedCandidateRef> hltCands;
341 rawTriggerEvent->getObjects( index,
TriggerMuon, hltCands );
342 for(
size_t j=0;
j<hltCands.size();
j++){
343 muonsAtFilter[
i].push_back(*hltCands[
j]);
345 muonPositionsAtFilter[
i].push_back(
LeafCandidate( hltCands[j]->
charge(),
math::XYZTLorentzVector(hltCands[j]->track()->innerPosition().
x(), hltCands[j]->track()->innerPosition().
y(), hltCands[j]->track()->innerPosition().
z(), 0.) ) );
350 ME[TString::Format(
"filt%dMuon_size",
int(
i+1))]->Fill(muonsAtFilter[
i].
size());
351 LogDebug(
"HLTriggerOfflineHeavyFlavor")<<
"Filter \""<<
filterNamesLevels[
i].first<<
"\" has "<<muonsAtFilter[
i].size()<<
" muons"<<endl;
354 LogDebug(
"HLTriggerOfflineHeavyFlavor")<<
"Could not access RAWTriggerEvent"<<endl;
358 vector<LeafCandidate> pathMuons;
363 for(
int i=0;
i<aodTriggerEvent->sizeFilters();
i++){
365 Keys keys = aodTriggerEvent->filterKeys(
i);
366 for(
size_t j=0;
j<keys.size();
j++){
371 ME[
"pathMuon_size"]->Fill(pathMuons.size());
374 LogDebug(
"HLTriggerOfflineHeavyFlavor")<<
"Could not access AODTriggerEvent"<<endl;
378 bool triggerFired =
false;
384 bool hlt_exists =
false;
385 for (
unsigned int i=0;
i!=triggerNames.
size();
i++) {
388 triggerFired = triggerResults->accept(
i );
401 vector<int> glob_gen(genMuons.size(),-1);
403 vector<vector<int> > filt_glob;
405 filt_glob.push_back( vector<int>(globMuons.size(),-1) );
414 vector<int> path_glob(globMuons.size(),-1);
425 for(
size_t i=0;
i<genMuons.size();
i++){
426 ME[
"genMuon_genEtaPt"]->Fill(genMuons[
i].
eta(), genMuons[
i].pt());
427 ME[
"genMuon_genEtaPhi"]->Fill(genMuons[
i].
eta(), genMuons[
i].
phi());
428 if(glob_gen[
i] != -1){
429 ME[
"resGlobGen_genEtaPt"]->Fill(genMuons[
i].
eta(), genMuons[
i].pt(), (globMuons[glob_gen[
i]].pt()-genMuons[
i].pt())/genMuons[
i].pt() );
430 ME[
"globMuon_genEtaPt"]->Fill(genMuons[
i].
eta(), genMuons[
i].pt());
431 ME[
"globMuon_genEtaPhi"]->Fill(genMuons[
i].
eta(), genMuons[
i].
phi());
432 ME[
"globMuon_recoEtaPt"]->Fill(globMuons[glob_gen[
i]].
eta(), globMuons[glob_gen[
i]].pt());
433 ME[
"globMuon_recoEtaPhi"]->Fill(globMuons[glob_gen[
i]].
eta(), globMuons[glob_gen[
i]].
phi());
435 if(filt_glob[
f][glob_gen[
i]] != -1){
436 ME[TString::Format(
"resFilt%dGlob_recoEtaPt",
int(
f+1))]->Fill(globMuons[glob_gen[
i]].
eta(), globMuons[glob_gen[
i]].pt(), (muonsAtFilter[
f][filt_glob[
f][glob_gen[
i]]].pt()-globMuons[glob_gen[
i]].pt())/globMuons[glob_gen[
i]].pt() );
437 ME[TString::Format(
"filt%dMuon_recoEtaPt",
int(
f+1))]->Fill(globMuons[glob_gen[
i]].
eta(), globMuons[glob_gen[
i]].pt());
438 ME[TString::Format(
"filt%dMuon_recoEtaPhi",
int(
f+1))]->Fill(globMuons[glob_gen[
i]].
eta(), globMuons[glob_gen[
i]].
phi());
443 if(path_glob[glob_gen[
i]] != -1){
444 ME[
"resPathGlob_recoEtaPt"]->Fill(globMuons[glob_gen[
i]].
eta(), globMuons[glob_gen[
i]].pt(), (pathMuons[path_glob[glob_gen[
i]]].pt()-globMuons[glob_gen[
i]].pt())/globMuons[glob_gen[
i]].pt() );
445 ME[
"pathMuon_recoEtaPt"]->Fill(globMuons[glob_gen[
i]].
eta(), globMuons[glob_gen[
i]].pt());
446 ME[
"pathMuon_recoEtaPhi"]->Fill(globMuons[glob_gen[
i]].
eta(), globMuons[glob_gen[
i]].
phi());
452 ME[
"resultMuon_recoEtaPt"]->Fill(globMuons[glob_gen[
i]].
eta(), globMuons[glob_gen[
i]].pt());
453 ME[
"resultMuon_recoEtaPhi"]->Fill(globMuons[glob_gen[
i]].
eta(), globMuons[glob_gen[
i]].
phi());
461 for(
size_t j=1;
j<genMuons.size();
j++){
462 if(genMuons[0].
charge()*genMuons[
j].charge()==-1){
471 double genDimuonPt = (genMuons[0].p4()+genMuons[secondMuon].p4()).pt();
472 double genDimuonEta = (genMuons[0].p4()+genMuons[secondMuon].p4()).
eta();
473 double genDimuonRap = (genMuons[0].p4()+genMuons[secondMuon].p4()).Rapidity();
474 double genDimuonDR = deltaR<LeafCandidate,LeafCandidate>( genMuons[0], genMuons[secondMuon] );
475 bool highPt = genMuons[0].pt()>7. && genMuons[secondMuon].pt()>7;
476 ME[
"genDimuon_genEtaPt"]->Fill( genDimuonEta, genDimuonPt );
477 ME[
"genDimuon_genRapPt"]->Fill( genDimuonRap, genDimuonPt );
478 if(highPt)
ME[
"genDimuon_genPtDR"]->Fill( genDimuonPt, genDimuonDR );
480 if(glob_gen[0]!=-1 && glob_gen[secondMuon]!=-1){
481 ME[
"globDimuon_genEtaPt"]->Fill( genDimuonEta, genDimuonPt );
482 ME[
"globDimuon_genRapPt"]->Fill( genDimuonRap, genDimuonPt );
483 if(highPt)
ME[
"globDimuon_genPtDR"]->Fill( genDimuonPt, genDimuonDR );
484 double globDimuonPt = (globMuons[glob_gen[0]].p4()+globMuons[glob_gen[secondMuon]].p4()).pt();
485 double globDimuonEta = (globMuons[glob_gen[0]].p4()+globMuons[glob_gen[secondMuon]].p4()).
eta();
486 double globDimuonRap = (globMuons[glob_gen[0]].p4()+globMuons[glob_gen[secondMuon]].p4()).Rapidity();
487 double globDimuonDR = deltaR<LeafCandidate,LeafCandidate>( globMuons[glob_gen[0]], globMuons[glob_gen[secondMuon]] );
488 double globDimuonDRpos = deltaR<LeafCandidate,LeafCandidate>( globMuons_position[glob_gen[0]], globMuons_position[glob_gen[secondMuon]] );
489 ME[
"globDimuon_recoEtaPt"]->Fill( globDimuonEta, globDimuonPt );
490 ME[
"globDimuon_recoRapPt"]->Fill( globDimuonRap, globDimuonPt );
491 if(highPt)
ME[
"globDimuon_recoPtDR"]->Fill( globDimuonPt, globDimuonDR );
492 if(highPt)
ME[
"globDimuon_recoPtDRpos"]->Fill( globDimuonPt, globDimuonDRpos );
495 if(filt_glob[
f][glob_gen[0]] != -1 && filt_glob[
f][glob_gen[secondMuon]] != -1){
496 ME[TString::Format(
"diFilt%dDimuon_recoEtaPt",
int(
f+1))]->Fill( globDimuonEta, globDimuonPt );
497 ME[TString::Format(
"diFilt%dDimuon_recoRapPt",
int(
f+1))]->Fill( globDimuonRap, globDimuonPt );
498 if(highPt)
ME[TString::Format(
"diFilt%dDimuon_recoPtDR",
int(
f+1))]->Fill( globDimuonPt, globDimuonDR );
499 if(highPt)
ME[TString::Format(
"diFilt%dDimuon_recoPtDRpos",
int(
f+1))]->Fill( globDimuonPt, globDimuonDRpos );
506 if(filt_glob[
f][glob_gen[0]] != -1 || filt_glob[
f][glob_gen[secondMuon]] != -1){
507 ME[TString::Format(
"filt%dDimuon_recoEtaPt",
int(
f+1))]->Fill( globDimuonEta, globDimuonPt );
508 ME[TString::Format(
"filt%dDimuon_recoRapPt",
int(
f+1))]->Fill( globDimuonRap, globDimuonPt );
509 if(highPt)
ME[TString::Format(
"filt%dDimuon_recoPtDR",
int(
f+1))]->Fill( globDimuonPt, globDimuonDR );
510 if(highPt)
ME[TString::Format(
"filt%dDimuon_recoPtDRpos",
int(
f+1))]->Fill( globDimuonPt, globDimuonDRpos );
516 if(path_glob[glob_gen[0]] != -1 && path_glob[glob_gen[secondMuon]] != -1){
517 ME[
"diPathDimuon_recoEtaPt"]->Fill( globDimuonEta, globDimuonPt );
518 ME[
"diPathDimuon_recoRapPt"]->Fill( globDimuonRap, globDimuonPt );
519 if(highPt)
ME[
"diPathDimuon_recoPtDR"]->Fill( globDimuonPt, globDimuonDR );
520 if(highPt)
ME[
"diPathDimuon_recoPtDRpos"]->Fill( globDimuonPt, globDimuonDRpos );
523 if(path_glob[glob_gen[0]] != -1 || path_glob[glob_gen[secondMuon]] != -1){
524 ME[
"pathDimuon_recoEtaPt"]->Fill( globDimuonEta, globDimuonPt );
525 ME[
"pathDimuon_recoRapPt"]->Fill( globDimuonRap, globDimuonPt );
526 if(highPt)
ME[
"pathDimuon_recoPtDR"]->Fill( globDimuonPt, globDimuonDR );
527 if(highPt)
ME[
"pathDimuon_recoPtDRpos"]->Fill( globDimuonPt, globDimuonDRpos );
531 ME[
"resultDimuon_recoEtaPt"]->Fill( globDimuonEta, globDimuonPt );
532 ME[
"resultDimuon_recoRapPt"]->Fill( globDimuonRap, globDimuonPt );
533 if(highPt)
ME[
"resultDimuon_recoPtDR"]->Fill( globDimuonPt, globDimuonDR );
534 if(highPt)
ME[
"resultDimuon_recoPtDRpos"]->Fill( globDimuonPt, globDimuonDRpos );
549 return mother->
pdgId();
557 vector<double>
dR(from.size());
558 for(
size_t i=0;
i<from.size();
i++){
562 for(
size_t j=0;
j<to.size();
j++){
563 double dRtmp = deltaR<double>(from[
i].eta(), from[
i].phi(), to[
j].eta(), to[
j].phi());
574 if(
dR[
i] > dRMatchingCut ){
579 for(
size_t k=0;
k<
i;
k++){
580 if( map[
k] != -1 && map[i] == map[
k] ){
596 int ptN = ptBins.size()==3 ? (int)ptBins[0]+1 : ptBins.size();
597 Double_t *pt =
new Double_t[ ptN ];
598 for(
int i=0;
i<ptN;
i++){
599 pt[
i] = ptBins.size()==3 ? ptBins[1] +
i*(ptBins[2]-ptBins[1])/ptBins[0] : ptBins[
i] ;
601 int etaN = etaBins.size()==3 ? (int)etaBins[0]+1 : etaBins.size();
602 Double_t *
eta =
new Double_t[ etaN ];
603 for(
int i=0;
i<etaN;
i++){
604 eta[
i] = etaBins.size()==3 ? etaBins[1] +
i*(etaBins[2]-etaBins[1])/etaBins[0] : etaBins[
i] ;
606 TH2F *
h =
new TH2F( name, name, ptN-1, pt, etaN-1, eta );
607 h->SetXTitle(ptLabel);
608 h->SetYTitle(etaLabel);
617 int ptN = ptBins.size()==3 ? (int)ptBins[0]+1 : ptBins.size();
618 Double_t *pt =
new Double_t[ ptN ];
619 for(
int i=0;
i<ptN;
i++){
620 pt[
i] = ptBins.size()==3 ? ptBins[1] +
i*(ptBins[2]-ptBins[1])/ptBins[0] : ptBins[
i] ;
622 int etaN = etaBins.size()==3 ? (int)etaBins[0]+1 : etaBins.size();
623 Double_t *
eta =
new Double_t[ etaN ];
624 for(
int i=0;
i<etaN;
i++){
625 eta[
i] = etaBins.size()==3 ? etaBins[1] +
i*(etaBins[2]-etaBins[1])/etaBins[0] : etaBins[
i] ;
627 TProfile2D *
h =
new TProfile2D( name, name, ptN-1, pt, etaN-1, eta );
628 h->SetXTitle(ptLabel);
629 h->SetYTitle(etaLabel);
638 int binsN = bins.size()==3 ? (int)bins[0]+1 : bins.size();
639 Double_t *myBins =
new Double_t[ binsN ];
640 for(
int i=0;
i<binsN;
i++){
641 myBins[
i] = bins.size()==3 ? bins[1] +
i*(bins[2]-bins[1])/bins[0] : bins[
i] ;
643 TH1F *
h =
new TH1F( name, name, binsN-1, myBins );
double globL1DeltaRMatchingCut
vector< double > dimuonPtBins
virtual edm::TriggerNames const & triggerNames(edm::TriggerResults const &triggerResults) const
InputTag triggerSummaryAODTag
virtual const Candidate * mother(size_type i=0) const =0
return pointer to mother
MonitorElement * book1D(const char *name, const char *title, int nchX, double lowX, double highX)
Book 1D histogram.
enum start value shifted to 81 so as to avoid clashes with PDG codes
vector< double > deltaPhiBins
#define DEFINE_FWK_MODULE(type)
const std::vector< std::string > & triggerNames() const
names of trigger paths
vector< pair< string, int > > filterNamesLevels
void myBook1D(TString name, vector< double > &xBins, TString label, TString title)
Strings::size_type size() const
double globL3DeltaRMatchingCut
string triggerProcessName
vector< double > dimuonEtaBins
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
void myBook2D(TString name, vector< double > &xBins, TString xLabel, vector< double > &yBins, TString yLabel)
U second(std::pair< T, U > const &p)
PtEtaPhiMLorentzVectorD PtEtaPhiMLorentzVector
Lorentz vector with cartesian internal representation.
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
vector< double > dimuonDRBins
int getMotherId(const Candidate *p)
InputTag triggerSummaryRAWTag
void match(MonitorElement *me, vector< LeafCandidate > &from, vector< LeafCandidate > &to, double deltaRMatchingCut, vector< int > &map)
void setVerbose(unsigned level)
double genGlobDeltaRMatchingCut
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
vector< double > muonPtBins
HeavyFlavorValidation(const edm::ParameterSet &)
double globL2DeltaRMatchingCut
const std::vector< std::string > & moduleLabels(unsigned int trigger) const
label(s) of module(s) on a trigger path
virtual int pdgId() const =0
PDG identifier.
std::vector< TriggerObject > TriggerObjectCollection
collection of trigger physics objects (e.g., all isolated muons)
vector< double > deltaEtaBins
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
std::string const & triggerName(unsigned int index) const
std::vector< size_type > Keys
bool init(const edm::Run &iRun, const edm::EventSetup &iSetup, const std::string &processName, bool &changed)
d'tor
vector< double > muonPhiBins
virtual void beginRun(const Run &iRun, const EventSetup &iSetup)
void myBookProfile2D(TString name, vector< double > &xBins, TString xLabel, vector< double > &yBins, TString yLabel, TString title)
void myBook1D(TString name, vector< double > &xBins, TString label)
void myBook2D(TString name, vector< double > &xBins, TString xLabel, vector< double > &yBins, TString yLabel, TString title)
std::pair< typename Association::data_type::first_type, double > match(Reference key, Association association, bool bestMatchByMaxValue)
Generic matching function.
InputTag triggerResultsTag
map< TString, MonitorElement * > ME
MonitorElement * book2D(const char *name, const char *title, int nchX, double lowX, double highX, int nchY, double lowY, double highY)
Book 2D histogram.
tuple size
Write out results.
void myBookProfile2D(TString name, vector< double > &xBins, TString xLabel, vector< double > &yBins, TString yLabel)
vector< double > muonEtaBins
void setCurrentFolder(const std::string &fullpath)
virtual void analyze(const edm::Event &, const edm::EventSetup &)
MonitorElement * bookProfile2D(const char *name, const char *title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, int nchZ, double lowZ, double highZ, const char *option="s")