CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
HeavyFlavorValidation.cc
Go to the documentation of this file.
1 
8 // Original Author: Zoltan Gecse
9 
10 #include <memory>
20 
36 
38 
42 
44 
45 #include "TLorentzVector.h"
46 
47 using namespace std;
48 using namespace edm;
49 using namespace reco;
50 using namespace l1extra;
51 using namespace trigger;
52 
54  public:
57  protected:
58  void dqmBeginRun(const edm::Run&, const edm::EventSetup&);
59  void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override;
60  virtual void analyze(const edm::Event&, const edm::EventSetup&) override;
61  private:
62  int getMotherId( const Candidate * p );
63  void match( MonitorElement * me, vector<LeafCandidate> & from, vector<LeafCandidate> & to, double deltaRMatchingCut, vector<int> & map );
64  void myBook2D(DQMStore::IBooker & ibooker, TString name, vector<double> &xBins, TString xLabel, vector<double> &yBins, TString yLabel, TString title);
65  void myBook2D(DQMStore::IBooker & ibooker, TString name, vector<double> &xBins, TString xLabel, vector<double> &yBins, TString yLabel){
66  myBook2D(ibooker, name, xBins, xLabel, yBins, yLabel, name);
67  }
68  void myBookProfile2D(DQMStore::IBooker & ibooker, TString name, vector<double> &xBins, TString xLabel, vector<double> &yBins, TString yLabel, TString title);
69  void myBookProfile2D(DQMStore::IBooker & ibooker, TString name, vector<double> &xBins, TString xLabel, vector<double> &yBins, TString yLabel){
70  myBookProfile2D(ibooker, name, xBins, xLabel, yBins, yLabel, name);
71  }
72  void myBook1D(DQMStore::IBooker & ibooker, TString name, vector<double> &xBins, TString label, TString title );
73  void myBook1D(DQMStore::IBooker & ibooker, TString name, vector<double> &xBins, TString label ){
74  myBook1D(ibooker, name, xBins, label, name );
75  }
76  string dqmFolder;
79 
85 
86  vector<int> motherIDs;
91  vector<double> deltaEtaBins;
92  vector<double> deltaPhiBins;
93  vector<double> muonPtBins;
94  vector<double> muonEtaBins;
95  vector<double> muonPhiBins;
96  vector<double> dimuonPtBins;
97  vector<double> dimuonEtaBins;
98  vector<double> dimuonDRBins;
100  map<TString, MonitorElement *> ME;
101  vector<pair<string,int> > filterNamesLevels;
102  const double muonMass;
103 };
104 
106 //get parameters
107  dqmFolder(pset.getUntrackedParameter<string>("DQMFolder")),
108  triggerProcessName(pset.getUntrackedParameter<string>("TriggerProcessName")),
109  triggerPathName(pset.getUntrackedParameter<string>("TriggerPathName")),
110  motherIDs(pset.getUntrackedParameter<vector<int> >("MotherIDs")),
111  genGlobDeltaRMatchingCut(pset.getUntrackedParameter<double>("GenGlobDeltaRMatchingCut")),
112  globL1DeltaRMatchingCut(pset.getUntrackedParameter<double>("GlobL1DeltaRMatchingCut")),
113  globL2DeltaRMatchingCut(pset.getUntrackedParameter<double>("GlobL2DeltaRMatchingCut")),
114  globL3DeltaRMatchingCut(pset.getUntrackedParameter<double>("GlobL3DeltaRMatchingCut")),
115  deltaEtaBins(pset.getUntrackedParameter<vector<double> >("DeltaEtaBins")),
116  deltaPhiBins(pset.getUntrackedParameter<vector<double> >("DeltaPhiBins")),
117  muonPtBins(pset.getUntrackedParameter<vector<double> >("MuonPtBins")),
118  muonEtaBins(pset.getUntrackedParameter<vector<double> >("MuonEtaBins")),
119  muonPhiBins(pset.getUntrackedParameter<vector<double> >("MuonPhiBins")),
120  dimuonPtBins(pset.getUntrackedParameter<vector<double> >("DimuonPtBins")),
121  dimuonEtaBins(pset.getUntrackedParameter<vector<double> >("DimuonEtaBins")),
122  dimuonDRBins(pset.getUntrackedParameter<vector<double> >("DimuonDRBins")),
123  muonMass(0.106)
124 {
125  triggerSummaryRAWTag = consumes<TriggerEventWithRefs>(InputTag( pset.getUntrackedParameter<string>("TriggerSummaryRAW"), "", triggerProcessName));
126  triggerSummaryAODTag = consumes<TriggerEvent>(InputTag( pset.getUntrackedParameter<string>("TriggerSummaryAOD"), "", triggerProcessName));
127  triggerResultsTag = consumes<TriggerResults>(InputTag( pset.getUntrackedParameter<string>("TriggerResults"), "", triggerProcessName));
128  recoMuonsTag = consumes<MuonCollection>(pset.getParameter<InputTag>("RecoMuons"));
129  genParticlesTag = consumes<GenParticleCollection>(pset.getParameter<InputTag>("GenParticles"));
130 }
131 
132 void HeavyFlavorValidation::dqmBeginRun(const edm::Run& iRun, const edm::EventSetup & iSetup ) {
133  //discover HLT configuration
135  bool isChanged;
136  if(hltConfig.init(iRun, iSetup, triggerProcessName, isChanged)){
137  LogDebug("HLTriggerOfflineHeavyFlavor") << "Successfully initialized HLTConfigProvider with process name: "<<triggerProcessName<<endl;
138  }else{
139  LogWarning("HLTriggerOfflineHeavyFlavor") << "Could not initialize HLTConfigProvider with process name: "<<triggerProcessName<<endl;
140  return;
141  }
142  stringstream os;
143  vector<string> triggerNames = hltConfig.triggerNames();
144  for( size_t i = 0; i < triggerNames.size(); i++) {
145  TString triggerName = triggerNames[i];
146  if (triggerName.Contains(triggerPathName)){
147  vector<string> moduleNames = hltConfig.moduleLabels( triggerNames[i] );
148  for( size_t j = 0; j < moduleNames.size(); j++) {
149  TString name = moduleNames[j];
150  if(name.Contains("Filter")){
151  int level = 0;
152  if(name.Contains("L1"))
153  level = 1;
154  if(name.Contains("L2"))
155  level = 2;
156  if(name.Contains("L3"))
157  level = 3;
158  if(name.Contains("mumuFilter") || name.Contains("DiMuon") || name.Contains("MuonL3Filtered") || name.Contains("TrackMassFiltered"))
159  level = 4;
160  if(name.Contains("Vertex") || name.Contains("Dz"))
161  level = 5;
162  filterNamesLevels.push_back( pair<string,int>(moduleNames[j],level) );
163  os<<" "<<moduleNames[j];
164  }
165  }
166  break;
167  }
168  }
169  if(filterNamesLevels.size()==0){
170  LogDebug("HLTriggerOfflineHeavyFlavor")<<"Bad Trigger Path: "<<triggerPathName<<endl;
171  return;
172  }else{
173  LogDebug("HLTriggerOfflineHeavyFlavor")<<"Trigger Path: "<<triggerPathName<<" has filters:"<<os.str();
174  }
175 }
176 
177 
178 
180  edm::Run const & iRun,
181  edm::EventSetup const & iSetup) {
182  ibooker.cd();
184 
185  // create Monitor Elements
186  // Eta Pt Single
187  myBook2D(ibooker, "genMuon_genEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)");
188  myBook2D(ibooker, "globMuon_genEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)");
189  myBook2D(ibooker, "globMuon_recoEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)");
190 
191  for(size_t i=0; i<filterNamesLevels.size(); i++){
192  myBook2D(ibooker, TString::Format("filt%dMuon_recoEtaPt",int(i+1)), muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)", filterNamesLevels[i].first);
193  }
194  myBook2D(ibooker, "pathMuon_recoEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)", triggerPathName);
195  myBook2D(ibooker, "resultMuon_recoEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)");
196  // Eta Pt Single Resolution
197  myBookProfile2D(ibooker, "resGlobGen_genEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)");
198  for(size_t i=0; i<filterNamesLevels.size(); i++){
199  myBookProfile2D(ibooker, TString::Format("resFilt%dGlob_recoEtaPt",int(i+1)), muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)", filterNamesLevels[i].first);
200  }
201  myBookProfile2D(ibooker, "resPathGlob_recoEtaPt", muonEtaBins, "#mu eta", muonPtBins, " #mu pT (GeV)", triggerPathName);
202  // Eta Pt Double
203  myBook2D(ibooker, "genDimuon_genEtaPt", dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)");
204  myBook2D(ibooker, "globDimuon_genEtaPt", dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)");
205  myBook2D(ibooker, "globDimuon_recoEtaPt", dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)");
206  for(size_t i=0; i<filterNamesLevels.size(); i++){
207  myBook2D(ibooker, TString::Format("filt%dDimuon_recoEtaPt",int(i+1)), dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)", filterNamesLevels[i].first);
208  }
209  myBook2D(ibooker, "pathDimuon_recoEtaPt", dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)", triggerPathName);
210  myBook2D(ibooker, "resultDimuon_recoEtaPt", dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)");
211  for(size_t i=0; i<filterNamesLevels.size(); i++){
212  myBook2D(ibooker, TString::Format("diFilt%dDimuon_recoEtaPt",int(i+1)), dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)", filterNamesLevels[i].first);
213  }
214  myBook2D(ibooker, "diPathDimuon_recoEtaPt", dimuonEtaBins, "#mu#mu eta", dimuonPtBins, " #mu#mu pT (GeV)", triggerPathName);
215  // Eta Phi Single
216  myBook2D(ibooker, "genMuon_genEtaPhi", muonEtaBins, "#mu eta", muonPhiBins, "#mu phi");
217  myBook2D(ibooker, "globMuon_genEtaPhi", muonEtaBins, "#mu eta", muonPhiBins, "#mu phi");
218  myBook2D(ibooker, "globMuon_recoEtaPhi", muonEtaBins, "#mu eta", muonPhiBins, "#mu phi");
219  for(size_t i=0; i<filterNamesLevels.size(); i++){
220  myBook2D(ibooker, TString::Format("filt%dMuon_recoEtaPhi",int(i+1)), muonEtaBins, "#mu eta", muonPhiBins, "#mu phi", filterNamesLevels[i].first);
221  }
222  myBook2D(ibooker, "pathMuon_recoEtaPhi", muonEtaBins, "#mu eta", muonPhiBins, "#mu phi", triggerPathName);
223  myBook2D(ibooker, "resultMuon_recoEtaPhi", muonEtaBins, "#mu eta", muonPhiBins, "#mu phi");
224  // Rap Pt Double
225  myBook2D(ibooker, "genDimuon_genRapPt", dimuonEtaBins, "#mu#mu rapidity", dimuonPtBins, " #mu#mu pT (GeV)");
226  myBook2D(ibooker, "globDimuon_genRapPt", dimuonEtaBins, "#mu#mu rapidity", dimuonPtBins, " #mu#mu pT (GeV)");
227  myBook2D(ibooker, "globDimuon_recoRapPt", dimuonEtaBins, "#mu#mu rapidity", dimuonPtBins, " #mu#mu pT (GeV)");
228  for(size_t i=0; i<filterNamesLevels.size(); i++){
229  myBook2D(ibooker, TString::Format("filt%dDimuon_recoRapPt",int(i+1)), dimuonEtaBins, "#mu#mu rapidity", dimuonPtBins, " #mu#mu pT (GeV)", filterNamesLevels[i].first);
230  }
231  myBook2D(ibooker, "pathDimuon_recoRapPt", dimuonEtaBins, "#mu#mu rapidity", dimuonPtBins, " #mu#mu pT (GeV)", triggerPathName);
232  myBook2D(ibooker, "resultDimuon_recoRapPt", dimuonEtaBins, "#mu#mu rapidity", dimuonPtBins, " #mu#mu pT (GeV)");
233  for(size_t i=0; i<filterNamesLevels.size(); i++){
234  myBook2D(ibooker, TString::Format("diFilt%dDimuon_recoRapPt",int(i+1)), dimuonEtaBins, "#mu#mu rapidity", dimuonPtBins, " #mu#mu pT (GeV)", filterNamesLevels[i].first);
235  }
236  myBook2D(ibooker, "diPathDimuon_recoRapPt", dimuonEtaBins, "#mu#mu rapidity", dimuonPtBins, " #mu#mu pT (GeV)", triggerPathName);
237  // Pt DR Double
238  myBook2D(ibooker, "genDimuon_genPtDR", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R at IP");
239  myBook2D(ibooker, "globDimuon_genPtDR", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R at IP");
240  myBook2D(ibooker, "globDimuon_recoPtDR", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R at IP");
241  for(size_t i=0; i<filterNamesLevels.size(); i++){
242  myBook2D(ibooker, TString::Format("filt%dDimuon_recoPtDR",int(i+1)), dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R at IP", filterNamesLevels[i].first);
243  }
244  myBook2D(ibooker, "pathDimuon_recoPtDR", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R at IP", triggerPathName);
245  for(size_t i=0; i<filterNamesLevels.size(); i++){
246  myBook2D(ibooker, TString::Format("diFilt%dDimuon_recoPtDR",int(i+1)), dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R at IP", filterNamesLevels[i].first);
247  }
248  myBook2D(ibooker, "diPathDimuon_recoPtDR", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R at IP", triggerPathName);
249  myBook2D(ibooker, "resultDimuon_recoPtDR", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R at IP");
250  // Pt DRpos Double
251  myBook2D(ibooker, "globDimuon_recoPtDRpos", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R in MS");
252  for(size_t i=0; i<filterNamesLevels.size(); i++){
253  myBook2D(ibooker, TString::Format("filt%dDimuon_recoPtDRpos",int(i+1)), dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R in MS", filterNamesLevels[i].first);
254  }
255  myBook2D(ibooker, "pathDimuon_recoPtDRpos", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R in MS", triggerPathName);
256  for(size_t i=0; i<filterNamesLevels.size(); i++){
257  myBook2D(ibooker, TString::Format("diFilt%dDimuon_recoPtDRpos",int(i+1)), dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R in MS", filterNamesLevels[i].first);
258  }
259  myBook2D(ibooker, "diPathDimuon_recoPtDRpos", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R in MS", triggerPathName);
260  myBook2D(ibooker, "resultDimuon_recoPtDRpos", dimuonPtBins, " #mu#mu pT (GeV)", dimuonDRBins, "#mu#mu #Delta R in MS");
261 
262  // Matching
263  myBook2D(ibooker, "globGen_deltaEtaDeltaPhi", deltaEtaBins, "#Delta eta", deltaPhiBins, "#Delta phi");
264  for(size_t i=0; i<filterNamesLevels.size(); i++){
265  myBook2D(ibooker, TString::Format("filt%dGlob_deltaEtaDeltaPhi",int(i+1)), deltaEtaBins, "#Delta eta", deltaPhiBins, "#Delta phi", filterNamesLevels[i].first);
266  }
267  myBook2D(ibooker, "pathGlob_deltaEtaDeltaPhi", deltaEtaBins, "#Delta eta", deltaPhiBins, "#Delta phi", triggerPathName);
268  // Size of containers
269  vector<double> sizeBins; sizeBins.push_back(10); sizeBins.push_back(0); sizeBins.push_back(10);
270  myBook1D(ibooker, "genMuon_size", sizeBins, "container size" );
271  myBook1D(ibooker, "globMuon_size", sizeBins, "container size" );
272  for(size_t i=0; i<filterNamesLevels.size(); i++){
273  myBook1D(ibooker, TString::Format("filt%dMuon_size",int(i+1)), sizeBins, "container size", filterNamesLevels[i].first);
274  }
275  myBook1D(ibooker, "pathMuon_size", sizeBins, "container size", triggerPathName );
276 }
277 
279  if(filterNamesLevels.size()==0){
280  return;
281  }
282 //access the containers and create LeafCandidate copies
283  vector<LeafCandidate> genMuons;
285  iEvent.getByToken(genParticlesTag, genParticles);
286  if(genParticles.isValid()){
287  for(GenParticleCollection::const_iterator p=genParticles->begin(); p!= genParticles->end(); ++p){
288  if( p->status() == 1 && std::abs(p->pdgId())==13 &&
289  ( find( motherIDs.begin(), motherIDs.end(), -1 )!=motherIDs.end() || find( motherIDs.begin(), motherIDs.end(), getMotherId(&(*p)) )!=motherIDs.end() ) ){
290  genMuons.push_back( *p );
291  }
292  }
293  }else{
294  LogDebug("HLTriggerOfflineHeavyFlavor")<<"Could not access GenParticleCollection"<<endl;
295  }
296  sort(genMuons.begin(), genMuons.end(), GreaterByPt<LeafCandidate>());
297  ME["genMuon_size"]->Fill(genMuons.size());
298  LogDebug("HLTriggerOfflineHeavyFlavor")<<"GenParticleCollection from "<<genParticlesTag<<" has size: "<<genMuons.size()<<endl;
299 
300  vector<LeafCandidate> globMuons;
301  vector<LeafCandidate> globMuons_position;
302  Handle<MuonCollection> recoMuonsHandle;
303  iEvent.getByToken(recoMuonsTag, recoMuonsHandle);
304  if(recoMuonsHandle.isValid()){
305  for(MuonCollection::const_iterator p=recoMuonsHandle->begin(); p!= recoMuonsHandle->end(); ++p){
306  if(p->isGlobalMuon()){
307  globMuons.push_back(*p);
308  globMuons_position.push_back( LeafCandidate( p->charge(), math::XYZTLorentzVector(p->outerTrack()->innerPosition().x(), p->outerTrack()->innerPosition().y(), p->outerTrack()->innerPosition().z(), 0.) ) );
309  }
310  }
311  }else{
312  LogDebug("HLTriggerOfflineHeavyFlavor")<<"Could not access reco Muons"<<endl;
313  }
314  ME["globMuon_size"]->Fill(globMuons.size());
315  LogDebug("HLTriggerOfflineHeavyFlavor")<<"Global Muons from "<<recoMuonsTag<<" has size: "<<globMuons.size()<<endl;
316 
317 // access RAW trigger event
318  vector<vector<LeafCandidate> > muonsAtFilter;
319  vector<vector<LeafCandidate> > muonPositionsAtFilter;
320  for(size_t i=0; i<filterNamesLevels.size(); i++){
321  muonsAtFilter.push_back(vector<LeafCandidate>());
322  muonPositionsAtFilter.push_back(vector<LeafCandidate>());
323  }
324  Handle<TriggerEventWithRefs> rawTriggerEvent;
325  iEvent.getByToken( triggerSummaryRAWTag, rawTriggerEvent );
326  if( rawTriggerEvent.isValid() ){
327  for(size_t i=0; i<filterNamesLevels.size(); i++){
328  size_t index = rawTriggerEvent->filterIndex(InputTag( filterNamesLevels[i].first, "", triggerProcessName ));
329  if ( index < rawTriggerEvent->size() ){
330  if( filterNamesLevels[i].second==1 ){
331  vector<L1MuonParticleRef> l1Cands;
332  rawTriggerEvent->getObjects( index, TriggerL1Mu, l1Cands );
333  for(size_t j=0; j<l1Cands.size(); j++){
334  muonsAtFilter[i].push_back(*l1Cands[j]);
335  }
336  }else{
337  vector<RecoChargedCandidateRef> hltCands;
338  rawTriggerEvent->getObjects( index, TriggerMuon, hltCands );
339  for(size_t j=0; j<hltCands.size(); j++){
340  muonsAtFilter[i].push_back(*hltCands[j]);
341  if( filterNamesLevels[i].second==2 ){
342  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.) ) );
343  }
344  }
345  }
346  }
347  ME[TString::Format("filt%dMuon_size",int(i+1))]->Fill(muonsAtFilter[i].size());
348  LogDebug("HLTriggerOfflineHeavyFlavor")<<"Filter \""<<filterNamesLevels[i].first<<"\" has "<<muonsAtFilter[i].size()<<" muons"<<endl;
349  }
350  }else{
351  LogDebug("HLTriggerOfflineHeavyFlavor")<<"Could not access RAWTriggerEvent"<<endl;
352  }
353 
354 // access AOD trigger event
355  vector<LeafCandidate> pathMuons;
356  Handle<TriggerEvent> aodTriggerEvent;
357  iEvent.getByToken(triggerSummaryAODTag,aodTriggerEvent);
358  if(aodTriggerEvent.isValid()){
359  TriggerObjectCollection allObjects = aodTriggerEvent->getObjects();
360  for(int i=0; i<aodTriggerEvent->sizeFilters(); i++){
361  if(aodTriggerEvent->filterTag(i)==InputTag((filterNamesLevels.end()-1)->first,"",triggerProcessName)){
362  Keys keys = aodTriggerEvent->filterKeys(i);
363  for(size_t j=0; j<keys.size(); j++){
364  pathMuons.push_back( LeafCandidate( allObjects[keys[j]].id()>0 ? 1 : -1, math::PtEtaPhiMLorentzVector( allObjects[keys[j]].pt(), allObjects[keys[j]].eta(), allObjects[keys[j]].phi(), muonMass ) ) );
365  }
366  }
367  }
368  ME["pathMuon_size"]->Fill(pathMuons.size());
369  LogDebug("HLTriggerOfflineHeavyFlavor")<<"Path \""<<triggerPathName<<"\" has "<<pathMuons.size()<<" muons at last filter \""<<(filterNamesLevels.end()-1)->first<<"\""<<endl;
370  }else{
371  LogDebug("HLTriggerOfflineHeavyFlavor")<<"Could not access AODTriggerEvent"<<endl;
372  }
373 
374 // access Trigger Results
375  bool triggerFired = false;
377  iEvent.getByToken(triggerResultsTag,triggerResults);
378  if(triggerResults.isValid()){
379  LogDebug("HLTriggerOfflineHeavyFlavor")<<"Successfully initialized "<<triggerResultsTag<<endl;
380  const edm::TriggerNames & triggerNames = iEvent.triggerNames(*triggerResults);
381  bool hlt_exists = false;
382  for ( unsigned int i=0; i!=triggerNames.size(); i++) {
383  TString hlt_name = triggerNames.triggerName(i);
384  if (hlt_name.Contains(triggerPathName)) {
385  triggerFired = triggerResults->accept( i );
386  hlt_exists = true;
387  break;
388  }
389  }
390  if (!hlt_exists) {
391  LogDebug("HLTriggerOfflineHeavyFlavor")<<triggerResultsTag<<" has no trigger: "<<triggerPathName<<endl;
392  }
393  }else{
394  LogDebug("HLTriggerOfflineHeavyFlavor")<<"Could not initialize "<<triggerResultsTag<<endl;
395  }
396 
397 //create matching maps
398  vector<int> glob_gen(genMuons.size(),-1);
399  match( ME["globGen_deltaEtaDeltaPhi"], genMuons, globMuons, genGlobDeltaRMatchingCut, glob_gen );
400  vector<vector<int> > filt_glob;
401  for(size_t i=0; i<filterNamesLevels.size(); i++){
402  filt_glob.push_back( vector<int>(globMuons.size(),-1) );
403  if( filterNamesLevels[i].second == 1 ){
404  match( ME[TString::Format("filt%dGlob_deltaEtaDeltaPhi",int(i+1))], globMuons_position, muonsAtFilter[i] ,globL1DeltaRMatchingCut, filt_glob[i] );
405  }else if( filterNamesLevels[i].second == 2 ){
406  match( ME[TString::Format("filt%dGlob_deltaEtaDeltaPhi",int(i+1))], globMuons_position, muonPositionsAtFilter[i] ,globL2DeltaRMatchingCut, filt_glob[i] );
407  }else if( filterNamesLevels[i].second > 2){
408  match( ME[TString::Format("filt%dGlob_deltaEtaDeltaPhi",int(i+1))], globMuons, muonsAtFilter[i] ,globL3DeltaRMatchingCut, filt_glob[i] );
409  }
410  }
411  vector<int> path_glob(globMuons.size(),-1);
412  if( (filterNamesLevels.end()-1)->second == 1 ){
413  match( ME["pathGlob_deltaEtaDeltaPhi"], globMuons_position, pathMuons ,globL1DeltaRMatchingCut, path_glob );
414  }else if( (filterNamesLevels.end()-1)->second == 2 ){
415  match( ME["pathGlob_deltaEtaDeltaPhi"], globMuons, pathMuons ,globL2DeltaRMatchingCut, path_glob );
416  }else if( (filterNamesLevels.end()-1)->second > 2){
417  match( ME["pathGlob_deltaEtaDeltaPhi"], globMuons, pathMuons ,globL3DeltaRMatchingCut, path_glob );
418  }
419 
420 //fill histos
421  bool first = true;
422  for(size_t i=0; i<genMuons.size(); i++){
423  ME["genMuon_genEtaPt"]->Fill(genMuons[i].eta(), genMuons[i].pt());
424  ME["genMuon_genEtaPhi"]->Fill(genMuons[i].eta(), genMuons[i].phi());
425  if(glob_gen[i] != -1) {
426  ME["resGlobGen_genEtaPt"]->Fill(genMuons[i].eta(), genMuons[i].pt(), (globMuons[glob_gen[i]].pt()-genMuons[i].pt())/genMuons[i].pt() );
427  ME["globMuon_genEtaPt"]->Fill(genMuons[i].eta(), genMuons[i].pt());
428  ME["globMuon_genEtaPhi"]->Fill(genMuons[i].eta(), genMuons[i].phi());
429  ME["globMuon_recoEtaPt"]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].pt());
430  ME["globMuon_recoEtaPhi"]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].phi());
431  for(size_t f=0; f<filterNamesLevels.size(); f++) {
432  if(filt_glob[f][glob_gen[i]] != -1) {
433  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() );
434  ME[TString::Format("filt%dMuon_recoEtaPt",int(f+1))]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].pt());
435  ME[TString::Format("filt%dMuon_recoEtaPhi",int(f+1))]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].phi());
436  }else{
437  break;
438  }
439  }
440  if(path_glob[glob_gen[i]] != -1){
441  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() );
442  ME["pathMuon_recoEtaPt"]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].pt());
443  ME["pathMuon_recoEtaPhi"]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].phi());
444  }
445 //highest pt muon
446  if( first ){
447  first = false;
448  if( triggerFired ){
449  ME["resultMuon_recoEtaPt"]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].pt());
450  ME["resultMuon_recoEtaPhi"]->Fill(globMuons[glob_gen[i]].eta(), globMuons[glob_gen[i]].phi());
451  }
452  }
453  }
454  }
455 
456 //fill dimuon histograms (highest pT, opposite charge)
457  int secondMuon = 0;
458  for(size_t j=1; j<genMuons.size(); j++){
459  if(genMuons[0].charge()*genMuons[j].charge()==-1){
460  secondMuon = j;
461  break;
462  }
463  }
464  if(secondMuon > 0){
465 //two generated
466  double genDimuonPt = (genMuons[0].p4()+genMuons[secondMuon].p4()).pt();
467  double genDimuonEta = (genMuons[0].p4()+genMuons[secondMuon].p4()).eta();
468  double genDimuonRap = (genMuons[0].p4()+genMuons[secondMuon].p4()).Rapidity();
469  double genDimuonDR = deltaR<LeafCandidate,LeafCandidate>( genMuons[0], genMuons[secondMuon] );
470  bool highPt = genMuons[0].pt()>7. && genMuons[secondMuon].pt()>7;
471  ME["genDimuon_genEtaPt"]->Fill( genDimuonEta, genDimuonPt );
472  ME["genDimuon_genRapPt"]->Fill( genDimuonRap, genDimuonPt );
473  if(highPt) ME["genDimuon_genPtDR"]->Fill( genDimuonPt, genDimuonDR );
474 //two global
475  if(glob_gen[0]!=-1 && glob_gen[secondMuon]!=-1){
476  ME["globDimuon_genEtaPt"]->Fill( genDimuonEta, genDimuonPt );
477  ME["globDimuon_genRapPt"]->Fill( genDimuonRap, genDimuonPt );
478  if(highPt) ME["globDimuon_genPtDR"]->Fill( genDimuonPt, genDimuonDR );
479  double globDimuonPt = (globMuons[glob_gen[0]].p4()+globMuons[glob_gen[secondMuon]].p4()).pt();
480  double globDimuonEta = (globMuons[glob_gen[0]].p4()+globMuons[glob_gen[secondMuon]].p4()).eta();
481  double globDimuonRap = (globMuons[glob_gen[0]].p4()+globMuons[glob_gen[secondMuon]].p4()).Rapidity();
482  double globDimuonDR = deltaR<LeafCandidate,LeafCandidate>( globMuons[glob_gen[0]], globMuons[glob_gen[secondMuon]] );
483  double globDimuonDRpos = deltaR<LeafCandidate,LeafCandidate>( globMuons_position[glob_gen[0]], globMuons_position[glob_gen[secondMuon]] );
484  ME["globDimuon_recoEtaPt"]->Fill( globDimuonEta, globDimuonPt );
485  ME["globDimuon_recoRapPt"]->Fill( globDimuonRap, globDimuonPt );
486  if(highPt) ME["globDimuon_recoPtDR"]->Fill( globDimuonPt, globDimuonDR );
487  if(highPt) ME["globDimuon_recoPtDRpos"]->Fill( globDimuonPt, globDimuonDRpos );
488 //two filter objects
489  for(size_t f=0; f<filterNamesLevels.size(); f++){
490  if(filt_glob[f][glob_gen[0]] != -1 && filt_glob[f][glob_gen[secondMuon]] != -1){
491  ME[TString::Format("diFilt%dDimuon_recoEtaPt",int(f+1))]->Fill( globDimuonEta, globDimuonPt );
492  ME[TString::Format("diFilt%dDimuon_recoRapPt",int(f+1))]->Fill( globDimuonRap, globDimuonPt );
493  if(highPt) ME[TString::Format("diFilt%dDimuon_recoPtDR",int(f+1))]->Fill( globDimuonPt, globDimuonDR );
494  if(highPt) ME[TString::Format("diFilt%dDimuon_recoPtDRpos",int(f+1))]->Fill( globDimuonPt, globDimuonDRpos );
495  }else{
496  break;
497  }
498  }
499 //one filter object
500  for(size_t f=0; f<filterNamesLevels.size(); f++){
501  if(filt_glob[f][glob_gen[0]] != -1 || filt_glob[f][glob_gen[secondMuon]] != -1){
502  ME[TString::Format("filt%dDimuon_recoEtaPt",int(f+1))]->Fill( globDimuonEta, globDimuonPt );
503  ME[TString::Format("filt%dDimuon_recoRapPt",int(f+1))]->Fill( globDimuonRap, globDimuonPt );
504  if(highPt) ME[TString::Format("filt%dDimuon_recoPtDR",int(f+1))]->Fill( globDimuonPt, globDimuonDR );
505  if(highPt) ME[TString::Format("filt%dDimuon_recoPtDRpos",int(f+1))]->Fill( globDimuonPt, globDimuonDRpos );
506  }else{
507  break;
508  }
509  }
510 //two path objects
511  if(path_glob[glob_gen[0]] != -1 && path_glob[glob_gen[secondMuon]] != -1){
512  ME["diPathDimuon_recoEtaPt"]->Fill( globDimuonEta, globDimuonPt );
513  ME["diPathDimuon_recoRapPt"]->Fill( globDimuonRap, globDimuonPt );
514  if(highPt) ME["diPathDimuon_recoPtDR"]->Fill( globDimuonPt, globDimuonDR );
515  if(highPt) ME["diPathDimuon_recoPtDRpos"]->Fill( globDimuonPt, globDimuonDRpos );
516  }
517 //one path object
518  if(path_glob[glob_gen[0]] != -1 || path_glob[glob_gen[secondMuon]] != -1){
519  ME["pathDimuon_recoEtaPt"]->Fill( globDimuonEta, globDimuonPt );
520  ME["pathDimuon_recoRapPt"]->Fill( globDimuonRap, globDimuonPt );
521  if(highPt) ME["pathDimuon_recoPtDR"]->Fill( globDimuonPt, globDimuonDR );
522  if(highPt) ME["pathDimuon_recoPtDRpos"]->Fill( globDimuonPt, globDimuonDRpos );
523  }
524 //trigger result
525  if( triggerFired ){
526  ME["resultDimuon_recoEtaPt"]->Fill( globDimuonEta, globDimuonPt );
527  ME["resultDimuon_recoRapPt"]->Fill( globDimuonRap, globDimuonPt );
528  if(highPt) ME["resultDimuon_recoPtDR"]->Fill( globDimuonPt, globDimuonDR );
529  if(highPt) ME["resultDimuon_recoPtDRpos"]->Fill( globDimuonPt, globDimuonDRpos );
530  }
531  }
532  }
533 }
534 
536  const Candidate* mother = p->mother();
537  if( mother ){
538  if( mother->pdgId() == p->pdgId() ){
539  return getMotherId(mother);
540  }else{
541  return mother->pdgId();
542  }
543  }else{
544  return 0;
545  }
546 }
547 
548 void HeavyFlavorValidation::match( MonitorElement * me, vector<LeafCandidate> & from, vector<LeafCandidate> & to, double dRMatchingCut, vector<int> & map ){
549  vector<double> dR(from.size());
550  for(size_t i=0; i<from.size(); i++){
551  map[i] = -1;
552  dR[i] = 10.;
553  //find closest
554  for(size_t j=0; j<to.size(); j++){
555  double dRtmp = deltaR<double>(from[i].eta(), from[i].phi(), to[j].eta(), to[j].phi());
556  if( dRtmp < dR[i] ){
557  dR[i] = dRtmp;
558  map[i] = j;
559  }
560  }
561  //fill matching histo
562  if( map[i] != -1 ){
563  me->Fill( to[map[i]].eta()-from[i].eta(), deltaPhi<double>(to[map[i]].phi(), from[i].phi()) );
564  }
565  //apply matching cut
566  if( dR[i] > dRMatchingCut ){
567  map[i] = -1;
568  }
569  //remove duplication
570  if( map[i] != -1 ){
571  for(size_t k=0; k<i; k++){
572  if( map[k] != -1 && map[i] == map[k] ){
573  if( dR[i] < dR[k] ){
574  map[k] = -1;
575  }else{
576  map[i] = -1;
577  }
578  break;
579  }
580  }
581  }
582  }
583 }
584 
585 void HeavyFlavorValidation::myBook2D(DQMStore::IBooker & ibooker, TString name, vector<double> &ptBins, TString ptLabel, vector<double> &etaBins, TString etaLabel, TString title )
586 {
587 // dqmStore->setCurrentFolder(dqmFolder+"/"+folder);
588  int ptN = ptBins.size()==3 ? (int)ptBins[0]+1 : ptBins.size();
589  Double_t *pt = new Double_t[ ptN ];
590  for(int i=0; i<ptN; i++){
591  pt[i] = ptBins.size()==3 ? ptBins[1] + i*(ptBins[2]-ptBins[1])/ptBins[0] : ptBins[i] ;
592  }
593  int etaN = etaBins.size()==3 ? (int)etaBins[0]+1 : etaBins.size();
594  Double_t *eta = new Double_t[ etaN ];
595  for(int i=0; i<etaN; i++){
596  eta[i] = etaBins.size()==3 ? etaBins[1] + i*(etaBins[2]-etaBins[1])/etaBins[0] : etaBins[i] ;
597  }
598  TH2F *h = new TH2F( name, name, ptN-1, pt, etaN-1, eta );
599  h->SetXTitle(ptLabel);
600  h->SetYTitle(etaLabel);
601  h->SetTitle(title);
602  ME[name] = ibooker.book2D( name.Data(), h );
603  delete h;
604 }
605 
606 void HeavyFlavorValidation::myBookProfile2D(DQMStore::IBooker & ibooker, TString name, vector<double> &ptBins, TString ptLabel, vector<double> &etaBins, TString etaLabel, TString title )
607 {
608 // dqmStore->setCurrentFolder(dqmFolder+"/"+folder);
609  int ptN = ptBins.size()==3 ? (int)ptBins[0]+1 : ptBins.size();
610  Double_t *pt = new Double_t[ ptN ];
611  for(int i=0; i<ptN; i++){
612  pt[i] = ptBins.size()==3 ? ptBins[1] + i*(ptBins[2]-ptBins[1])/ptBins[0] : ptBins[i] ;
613  }
614  int etaN = etaBins.size()==3 ? (int)etaBins[0]+1 : etaBins.size();
615  Double_t *eta = new Double_t[ etaN ];
616  for(int i=0; i<etaN; i++){
617  eta[i] = etaBins.size()==3 ? etaBins[1] + i*(etaBins[2]-etaBins[1])/etaBins[0] : etaBins[i] ;
618  }
619  TProfile2D *h = new TProfile2D( name, name, ptN-1, pt, etaN-1, eta );
620  h->SetXTitle(ptLabel);
621  h->SetYTitle(etaLabel);
622  h->SetTitle(title);
623  ME[name] = ibooker.bookProfile2D( name.Data(), h );
624  delete h;
625 }
626 
627 void HeavyFlavorValidation::myBook1D(DQMStore::IBooker & ibooker, TString name, vector<double> &bins, TString label, TString title )
628 {
629 // dqmStore->setCurrentFolder(dqmFolder+"/"+folder);
630  int binsN = bins.size()==3 ? (int)bins[0]+1 : bins.size();
631  Double_t *myBins = new Double_t[ binsN ];
632  for(int i=0; i<binsN; i++){
633  myBins[i] = bins.size()==3 ? bins[1] + i*(bins[2]-bins[1])/bins[0] : bins[i] ;
634  }
635  TH1F *h = new TH1F( name, name, binsN-1, myBins );
636  h->SetXTitle(label);
637  h->SetTitle(title);
638  ME[name] = ibooker.book1D( name.Data(), h );
639  delete h;
640 }
641 
643 }
644 
645 //define this as a plug-in
#define LogDebug(id)
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
int i
Definition: DBlmapReader.cc:9
virtual edm::TriggerNames const & triggerNames(edm::TriggerResults const &triggerResults) const
Definition: Event.cc:208
void myBook2D(DQMStore::IBooker &ibooker, TString name, vector< double > &xBins, TString xLabel, vector< double > &yBins, TString yLabel, TString title)
virtual const Candidate * mother(size_type i=0) const =0
return pointer to mother
virtual void analyze(const edm::Event &, const edm::EventSetup &) override
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:449
enum start value shifted to 81 so as to avoid clashes with PDG codes
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
void cd(void)
Definition: DQMStore.cc:266
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
const std::vector< std::string > & triggerNames() const
names of trigger paths
vector< pair< string, int > > filterNamesLevels
const reco::GenParticle * mother(const reco::GenParticle &p, unsigned int imoth=0)
const double muonMass
Strings::size_type size() const
Definition: TriggerNames.cc:39
T eta() const
vector< double > dimuonEtaBins
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:7
float float float z
void bookHistograms(fwlite::EventContainer &eventCont)
EDGetTokenT< TriggerEvent > triggerSummaryAODTag
EDGetTokenT< GenParticleCollection > genParticlesTag
void Fill(long long x)
U second(std::pair< T, U > const &p)
PtEtaPhiMLorentzVectorD PtEtaPhiMLorentzVector
Lorentz vector with cartesian internal representation.
Definition: LorentzVector.h:25
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
Definition: ME.h:11
void myBookProfile2D(DQMStore::IBooker &ibooker, TString name, vector< double > &xBins, TString xLabel, vector< double > &yBins, TString yLabel, TString title)
int getMotherId(const Candidate *p)
int iEvent
Definition: GenABIO.cc:230
void match(MonitorElement *me, vector< LeafCandidate > &from, vector< LeafCandidate > &to, double deltaRMatchingCut, vector< int > &map)
MonitorElement * bookProfile2D(Args &&...args)
Definition: DQMStore.h:163
EDGetTokenT< TriggerEventWithRefs > triggerSummaryRAWTag
MonitorElement * book1D(Args &&...args)
Definition: DQMStore.h:115
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
int j
Definition: DBlmapReader.cc:9
void dqmBeginRun(const edm::Run &, const edm::EventSetup &)
double f[11][100]
void myBook1D(DQMStore::IBooker &ibooker, TString name, vector< double > &xBins, TString label, TString title)
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4
static std::string const triggerResults
Definition: EdmProvDump.cc:40
bool isValid() const
Definition: HandleBase.h:75
HeavyFlavorValidation(const edm::ParameterSet &)
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)
Definition: TriggerObject.h:81
EDGetTokenT< TriggerResults > triggerResultsTag
void setCurrentFolder(const std::string &fullpath)
Definition: DQMStore.cc:274
MonitorElement * book2D(Args &&...args)
Definition: DQMStore.h:133
void myBook2D(DQMStore::IBooker &ibooker, TString name, vector< double > &xBins, TString xLabel, vector< double > &yBins, TString yLabel)
std::string const & triggerName(unsigned int index) const
Definition: TriggerNames.cc:27
std::vector< size_type > Keys
bool init(const edm::Run &iRun, const edm::EventSetup &iSetup, const std::string &processName, bool &changed)
d&#39;tor
void myBookProfile2D(DQMStore::IBooker &ibooker, TString name, vector< double > &xBins, TString xLabel, vector< double > &yBins, TString yLabel)
std::pair< typename Association::data_type::first_type, double > match(Reference key, Association association, bool bestMatchByMaxValue)
Generic matching function.
Definition: Utils.h:10
tuple level
Definition: testEve_cfg.py:34
void myBook1D(DQMStore::IBooker &ibooker, TString name, vector< double > &xBins, TString label)
Definition: DDAxes.h:10
map< TString, MonitorElement * > ME
tuple size
Write out results.
EDGetTokenT< MuonCollection > recoMuonsTag
Definition: Run.h:41
Definition: DDAxes.h:10