CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_1_8_patch9/src/DQM/HLTEvF/plugins/HLTMonHcalIsoTrack.cc

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 //
00003 // Package:    HLTriggerOffline/special
00004 // Class:      DQMHcalIsoTrackAlCaRaw
00005 // 
00013 //
00014 // Original Author:  Grigory SAFRONOV
00015 //         Created:  Mon Oct  6 10:10:22 CEST 2008
00016 // $Id: HLTMonHcalIsoTrack.cc,v 1.5 2010/08/07 14:55:56 wmtan Exp $
00017 //
00018 //
00019 
00020 
00021 // user include files
00022 
00023 #include "DQM/HLTEvF/interface/HLTMonHcalIsoTrack.h"
00024 #include "DataFormats/HLTReco/interface/TriggerEventWithRefs.h"
00025 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
00026 
00027 #include "DataFormats/L1Trigger/interface/L1JetParticle.h"
00028 #include "DataFormats/L1Trigger/interface/L1JetParticleFwd.h"
00029 
00030 #include "DataFormats/TrackReco/interface/TrackFwd.h"
00031 #include "DataFormats/TrackReco/interface/Track.h"
00032 
00033 #include "DataFormats/HcalIsolatedTrack/interface/IsolatedPixelTrackCandidate.h"
00034 #include "DataFormats/HcalIsolatedTrack/interface/IsolatedPixelTrackCandidateFwd.h"
00035 
00036 #include "DataFormats/Math/interface/deltaR.h"
00037 
00038 double HLTMonHcalIsoTrack::getDist(double eta1, double phi1, double eta2, double phi2)
00039 {
00040   double dphi = fabs(phi1 - phi2); 
00041   if(dphi>acos(-1)) dphi = 2*acos(-1)-dphi;
00042   double dr = sqrt(dphi*dphi + pow(eta1-eta2,2));
00043   return dr;
00044 }
00045 
00046 std::pair<int,int> HLTMonHcalIsoTrack::towerIndex(double eta, double phi) 
00047 {
00048   int ieta = 0;
00049   int iphi = 0;
00050   if (eta!=0)
00051     {
00052       for (int i=1; i<21; i++)
00053         {
00054           if (fabs(eta)<(i*0.087)&&fabs(eta)>(i-1)*0.087) ieta=int(fabs(eta)/eta)*i;
00055         }
00056       if (fabs(eta)>1.740&&fabs(eta)<=1.830) ieta=int(fabs(eta)/eta)*21;
00057       if (fabs(eta)>1.830&&fabs(eta)<=1.930) ieta=int(fabs(eta)/eta)*22;
00058       if (fabs(eta)>1.930&&fabs(eta)<=2.043) ieta=int(fabs(eta)/eta)*23;
00059       if (fabs(eta)>2.043&&fabs(eta)<=2.172) ieta=int(fabs(eta)/eta)*24;
00060     }
00061   
00062   double delta=phi+0.174532925;
00063   if (delta<0) delta=delta+2*acos(-1);
00064   if (fabs(eta)<1.740) 
00065     {
00066       for (int i=0; i<72; i++)
00067         {
00068           if (delta<(i+1)*0.087266462&&delta>i*0.087266462) iphi=i;
00069         }
00070     }
00071   else 
00072     {
00073       for (int i=0; i<36; i++)
00074         {
00075           if (delta<2*(i+1)*0.087266462&&delta>2*i*0.087266462) iphi=2*i;
00076         }
00077     }
00078   
00079   return std::pair<int,int>(ieta,iphi);
00080 }
00081 
00082 
00083 HLTMonHcalIsoTrack::HLTMonHcalIsoTrack(const edm::ParameterSet& iConfig)
00084 {
00085   folderName_ = iConfig.getParameter<std::string>("folderName");
00086   outRootFileName_=iConfig.getParameter<std::string>("outputRootFileName");
00087   
00088   useProducerCollections_=iConfig.getParameter<bool>("useProducerCollections");
00089   hltRAWEventTag_=iConfig.getParameter<std::string>("hltRAWTriggerEventLabel");
00090   hltAODEventTag_=iConfig.getParameter<std::string>("hltAODTriggerEventLabel");
00091   
00092   hltProcess_=iConfig.getParameter<std::string>("hltProcessName");
00093   
00094   saveToRootFile_=iConfig.getParameter<bool>("SaveToRootFile");
00095   
00096   triggers = iConfig.getParameter<std::vector<edm::ParameterSet> >("triggers");
00097   for (std::vector<edm::ParameterSet>::iterator inTrig = triggers.begin(); inTrig != triggers.end(); inTrig++)
00098     {
00099       trigNames_.push_back(inTrig->getParameter<std::string>("triggerName"));
00100       l1filterLabels_.push_back(inTrig->getParameter<std::string>("hltL1filterLabel"));
00101       l2filterLabels_.push_back(inTrig->getParameter<std::string>("hltL2filterLabel"));
00102       l3filterLabels_.push_back(inTrig->getParameter<std::string>("hltL3filterLabel"));
00103       l2collectionLabels_.push_back(inTrig->getParameter<std::string>("l2collectionLabel"));
00104       l3collectionLabels_.push_back(inTrig->getParameter<std::string>("l3collectionLabel"));
00105     }
00106 }
00107 
00108 
00109 HLTMonHcalIsoTrack::~HLTMonHcalIsoTrack()
00110 {}
00111 
00112 void HLTMonHcalIsoTrack::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
00113 {
00114   edm::Handle<trigger::TriggerEventWithRefs> triggerObj;
00115   edm::InputTag toLab=edm::InputTag(hltRAWEventTag_,"",hltProcess_);
00116   iEvent.getByLabel(toLab,triggerObj); 
00117   if(!triggerObj.isValid()) return;
00118   
00119   for (unsigned int trInd=0; trInd<triggers.size(); trInd++)
00120     {
00121       bool l1pass=false;
00122       std::vector<l1extra::L1JetParticleRef> l1CenJets;
00123       std::vector<l1extra::L1JetParticleRef> l1ForJets;
00124       std::vector<l1extra::L1JetParticleRef> l1TauJets;
00125       edm::InputTag l1Tag = edm::InputTag(l1filterLabels_[trInd], "",hltProcess_);
00126       trigger::size_type l1filterIndex=triggerObj->filterIndex(l1Tag);
00127 
00128       if (l1filterIndex<triggerObj->size())
00129         {
00130           triggerObj->getObjects(l1filterIndex, trigger::TriggerL1CenJet, l1CenJets);
00131           triggerObj->getObjects(l1filterIndex, trigger::TriggerL1ForJet, l1ForJets);
00132           triggerObj->getObjects(l1filterIndex, trigger::TriggerL1TauJet, l1TauJets);
00133         }
00134       
00135       if (l1CenJets.size()>0||l1ForJets.size()>0||l1TauJets.size()>0) 
00136         {
00137           hL2L3acc[trInd]->Fill(1+0.001,1);
00138           l1pass=true;
00139         }
00140       
00141       if (!l1pass) continue;
00142       std::vector<reco::IsolatedPixelTrackCandidateRef> l2tracks;
00143       edm::InputTag l2Tag = edm::InputTag(l2filterLabels_[trInd],"",hltProcess_);
00144       trigger::size_type l2filterIndex=triggerObj->filterIndex(l2Tag);
00145       if (l2filterIndex<triggerObj->size()) triggerObj->getObjects(l2filterIndex, trigger::TriggerTrack, l2tracks); 
00146       
00147       std::vector<reco::IsolatedPixelTrackCandidateRef> l3tracks;
00148       edm::InputTag l3Tag = edm::InputTag(l3filterLabels_[trInd], "",hltProcess_);
00149       trigger::size_type l3filterIndex=triggerObj->filterIndex(l3Tag);
00150       if (l3filterIndex<triggerObj->size()) triggerObj->getObjects(l3filterIndex, trigger::TriggerTrack, l3tracks);
00151       
00152       edm::Handle<reco::IsolatedPixelTrackCandidateCollection> l3col;
00153       edm::InputTag l3colTag=edm::InputTag(l3collectionLabels_[trInd],"",hltProcess_);
00154       
00155       if (l2tracks.size()>0) 
00156         {
00157           hL2L3acc[trInd]->Fill(2+0.0001,1);
00158           if (useProducerCollections_)
00159             {
00160               iEvent.getByLabel(l3collectionLabels_[trInd],l3col);
00161               if(!l3col.isValid()) continue;
00162               
00163               for (reco::IsolatedPixelTrackCandidateCollection::const_iterator l3it=l3col->begin(); l3it!=l3col->end(); ++l3it)
00164                 {
00165                   double drmin=100;
00166                   int selL2tr=-1;
00167                   for (unsigned int il2=0; il2<l2tracks.size(); il2++)
00168                     {
00169                       double drl2l3=reco::deltaR(l3it->eta(),l3it->phi(),l2tracks[il2]->eta(),l2tracks[il2]->phi());
00170                       if (drl2l3<drmin)
00171                         {
00172                           drmin=drl2l3;
00173                           selL2tr=il2;
00174                         }
00175                     }
00176                   if (selL2tr!=-1&&drmin<0.03&&l2tracks[selL2tr]->p()!=0) hL3candL2rat[trInd]->Fill(l3it->p()/l2tracks[selL2tr]->p(),1);
00177                   if (selL2tr!=-1) hL3L2trackMatch[trInd]->Fill(drmin,1);
00178                 }
00179             }
00180         }
00181       
00182       if (l3tracks.size()>0) hL2L3acc[trInd]->Fill(3+0.0001,1);
00183       
00184       l1extra::L1JetParticleRef maxPtJet;
00185       
00186       double l1maxPt=-1;
00187       for (unsigned int i=0; i<l1CenJets.size(); i++)
00188         {
00189           if (l1CenJets[i]->pt()>l1maxPt) 
00190             {
00191               l1maxPt=l1CenJets[i]->pt();
00192               maxPtJet=l1CenJets[i];
00193             }
00194         }
00195       for (unsigned int i=0; i<l1ForJets.size(); i++)
00196         {
00197           if (l1ForJets[i]->pt()>l1maxPt) 
00198             {
00199               l1maxPt=l1ForJets[i]->pt();
00200               maxPtJet=l1ForJets[i];
00201             }
00202         }
00203       for (unsigned int i=0; i<l1TauJets.size(); i++)
00204         {
00205           if (l1TauJets[i]->pt()>l1maxPt) 
00206             {
00207               l1maxPt=l1TauJets[i]->pt();
00208               maxPtJet=l1TauJets[i];
00209             }
00210         }
00211       
00212       if (maxPtJet.isNonnull()) hL1eta[trInd]->Fill(maxPtJet->eta(),1);
00213       if (maxPtJet.isNonnull()) hL1phi[trInd]->Fill(maxPtJet->phi(),1);
00214       
00215       edm::Handle<reco::IsolatedPixelTrackCandidateCollection> l2col;
00216       edm::InputTag l2colTag=edm::InputTag(l2collectionLabels_[trInd],"",hltProcess_);
00217       if (useProducerCollections_)
00218         {
00219           iEvent.getByLabel(l2collectionLabels_[trInd],l2col);
00220           if(!l2col.isValid()) continue;
00221           for (reco::IsolatedPixelTrackCandidateCollection::const_iterator l2it=l2col->begin(); l2it!=l2col->end(); l2it++)
00222             {
00223               hL2isolationP[trInd]->Fill(l2it->maxPtPxl(),1);
00224             }
00225         }
00226       
00227       for (unsigned int i=0; i<l2tracks.size(); i++)
00228         {
00229           std::pair<int, int> tower=towerIndex(l2tracks[i]->eta(), l2tracks[i]->phi());
00230           hL2TowerOccupancy[trInd]->Fill(tower.first,tower.second,1);
00231         }
00232       for (unsigned int i=0; i<l3tracks.size(); i++)
00233         {
00234           std::pair<int, int> tower=towerIndex(l3tracks[i]->eta(), l3tracks[i]->phi());
00235           hL3TowerOccupancy[trInd]->Fill(tower.first,tower.second,1);
00236         }
00237     }
00238 }
00239 
00240 void HLTMonHcalIsoTrack::beginJob()
00241 {
00242   dbe_ = edm::Service<DQMStore>().operator->();
00243   dbe_->setCurrentFolder(folderName_);
00244 
00245   char tmp1[100];
00246   char tmp2[100];
00247   for (unsigned int i=0; i<triggers.size(); i++)
00248     {
00249       std::sprintf(tmp1,"hL2L3acc_%s",trigNames_[i].c_str());
00250       std::sprintf(tmp2,"number of L1, L2 and L3 accepts; [%s]",trigNames_[i].c_str());
00251       MonitorElement* hL2L3accBuf=dbe_->book1D(tmp1,tmp2,3,1,4);
00252       hL2L3acc.push_back(hL2L3accBuf);
00253       hL2L3acc[i]->setTitle(tmp1);
00254       hL2L3acc[i]->setAxisTitle("trigger level",1);
00255       
00256       std::sprintf(tmp1,"hL2L3trackMatch_%s",trigNames_[i].c_str());
00257       std::sprintf(tmp2,"R from L3 object to L2 object; [%s]",trigNames_[i].c_str());
00258       MonitorElement* hL3L2trackMatchBuf=dbe_->book1D(tmp1,tmp2,1000,0,1);
00259       hL3L2trackMatch.push_back(hL3L2trackMatchBuf);
00260       hL3L2trackMatch[i]->setAxisTitle("R(eta,phi)",1);
00261       
00262       std::sprintf(tmp1,"hL2L3rat_%s",trigNames_[i].c_str());
00263       std::sprintf(tmp2,"ratio of L3 to L2 momentum measurement; [%s]",trigNames_[i].c_str());
00264       MonitorElement* hL3L2ratBuf=dbe_->book1D(tmp1,tmp2,1000,0,10);
00265       hL3candL2rat.push_back(hL3L2ratBuf);
00266       hL3candL2rat[i]->setAxisTitle("P_L3/P_L2",1);
00267     
00268       std::sprintf(tmp1,"hL1eta_%s",trigNames_[i].c_str());
00269       std::sprintf(tmp2,"eta distribution of L1 triggers; [%s]",trigNames_[i].c_str());
00270       MonitorElement* hL1etaBuf=dbe_->book1D(tmp1,tmp2,1000,-7,7);
00271       hL1eta.push_back(hL1etaBuf);
00272       hL1eta[i]->setAxisTitle("eta",1);
00273       
00274       std::sprintf(tmp1,"hL1phi_%s",trigNames_[i].c_str());
00275       std::sprintf(tmp2,"phi distribution of L1 triggers; [%s]",trigNames_[i].c_str());
00276       MonitorElement* hL1phiBuf=dbe_->book1D(tmp1,tmp2,1000,-4,4);
00277       hL1phi.push_back(hL1phiBuf);
00278       hL1phi[i]->setAxisTitle("phi",1);
00279       
00280       std::sprintf(tmp1,"hL2isolation_%s",trigNames_[i].c_str());
00281       std::sprintf(tmp2,"isolation momentum at L2; [%s]",trigNames_[i].c_str());
00282       MonitorElement* hL2isolationBuf=dbe_->book1D(tmp1,tmp2,1000,0,10);
00283       hL2isolationP.push_back(hL2isolationBuf);
00284       hL2isolationP[i]->setAxisTitle("max P (GeV)",1);
00285       
00286       std::sprintf(tmp1,"hL2occupancy_%s",trigNames_[i].c_str());
00287       std::sprintf(tmp2,"tower occupancy at L2; [%s]",trigNames_[i].c_str());
00288       MonitorElement* hL2TowerOccupancyBuf=dbe_->book2D(tmp1,tmp2,48,-25,25,73,0,73);
00289       hL2TowerOccupancy.push_back(hL2TowerOccupancyBuf);
00290       hL2TowerOccupancy[i]->setAxisTitle("ieta",1);
00291       hL2TowerOccupancy[i]->setAxisTitle("iphi",2);
00292       hL2TowerOccupancy[i]->getTH2F()->SetOption("colz");
00293       hL2TowerOccupancy[i]->getTH2F()->SetStats(kFALSE);
00294   
00295       
00296       std::sprintf(tmp1,"hL3occupancy_%s",trigNames_[i].c_str());
00297       std::sprintf(tmp2,"tower occupancy at L3; [%s]",trigNames_[i].c_str());
00298       MonitorElement* hL3TowerOccupancyBuf=dbe_->book2D(tmp1,tmp2,48,-25,25,73,0,73);
00299       hL3TowerOccupancy.push_back(hL3TowerOccupancyBuf);
00300       hL3TowerOccupancy[i]->setAxisTitle("ieta",1);
00301       hL3TowerOccupancy[i]->setAxisTitle("iphi",2);
00302       hL3TowerOccupancy[i]->getTH2F()->SetOption("colz");
00303       hL3TowerOccupancy[i]->getTH2F()->SetStats(kFALSE);
00304     }
00305 }
00306 
00307 void HLTMonHcalIsoTrack::endJob() {
00308 
00309 if(dbe_&&saveToRootFile_) 
00310   {  
00311     dbe_->save(outRootFileName_);
00312   }
00313 }
00314 
00315 DEFINE_FWK_MODULE(HLTMonHcalIsoTrack);