CMS 3D CMS Logo

HLTTauCaloDQMOfflineSource.cc

Go to the documentation of this file.
00001 #include "DQMOffline/Trigger/interface/HLTTauCaloDQMOfflineSource.h"
00002 #include "Math/GenVector/VectorUtil.h"
00003 #include <iostream>
00004 #include <iomanip>
00005 #include <fstream>
00006 
00007 HLTTauCaloDQMOfflineSource::HLTTauCaloDQMOfflineSource(const edm::ParameterSet& iConfig):
00008  l2TauInfoAssoc_(iConfig.getParameter<edm::InputTag>("L2InfoAssociationInput")),
00009  mcColl_(iConfig.getParameter<edm::InputTag>("refCollection")),
00010  met_(iConfig.getParameter<edm::InputTag >("MET")),
00011  doRef_(iConfig.getParameter<bool>("doReference")),
00012  matchDeltaRMC_(iConfig.getParameter<double>("MatchDeltaR")),
00013  triggerTag_((iConfig.getParameter<std::string>("DQMFolder"))),
00014  l2Isolated_(iConfig.getParameter<edm::InputTag>("L2IsolatedJets")),
00015  outFile_(iConfig.getParameter<std::string>("OutputFileName")),
00016  EtMin_(iConfig.getParameter<double>("EtMin")),
00017  EtMax_(iConfig.getParameter<double>("EtMax")),
00018  NBins_(iConfig.getParameter<int>("NBins"))
00019 
00020 {
00021 
00022    store = &*edm::Service<DQMStore>();
00023   
00024   if(store)
00025     {
00026       //Create the histograms
00027       store->setCurrentFolder(triggerTag_);
00028       jetEt= store->book1D("L2tauCandEt","tauCandEt",NBins_,EtMin_,EtMax_);
00029       jetEta= store->book1D("L2tauCandEta","tauCandEta",50,-2.5,2.5);
00030       jetPhi= store->book1D("L2tauCandPhi","tauCandPhi",63,-3.14,3.14);
00031       ecalIsolEt=store->book1D("L2ecalIsolEt","ecalIsolEt",40,0,20);
00032       towerIsolEt=store->book1D("L2towerIsolEt","towerIsolEt",40,0,20);
00033       seedTowerEt=store->book1D("L2seedTowerEt","seedTowerEt",40,0,80);
00034       nClusters=store->book1D("L2nClusters","nClusters",20,0,20);
00035       clusterEtaRMS=store->book1D("L2clusterEtaRMS","clusterEtaRMS",25,0,0.5);
00036       clusterPhiRMS=store->book1D("L2clusterPhiRMS","clusterPhiRMS",25,0,0.5);
00037       clusterDeltaRRMS=store->book1D("L2clusterDeltaRRMS","clusterDeltaRRMS",25,0,0.5);
00038       EtEffNum=store->book1D("L2EtEffNum","Efficiency vs E_{t}(Numerator)",NBins_,EtMin_,EtMax_);
00039       EtEffDenom=store->book1D("L2EtEffDenom","Efficiency vs E_{t}(Denominator)",NBins_,EtMin_,EtMax_);
00040       MET=store->book1D("MET","Missing E_{t}",NBins_,EtMin_,EtMax_);
00041       
00042     }
00043  
00044 }
00045 
00046 
00047 HLTTauCaloDQMOfflineSource::~HLTTauCaloDQMOfflineSource()
00048 {
00049 }
00050 
00051 
00052 
00053 void
00054 HLTTauCaloDQMOfflineSource::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
00055 {
00056    using namespace edm;
00057    using namespace reco;
00058 
00059    Handle<L2TauInfoAssociation> l2TauInfoAssoc; //Handle to the input (L2 Tau Info Association)
00060    iEvent.getByLabel(l2TauInfoAssoc_,l2TauInfoAssoc);
00061 
00062    Handle<LVColl> McInfo; //Handle To The Truth!!!!
00063    iEvent.getByLabel(mcColl_,McInfo);
00064 
00065    Handle<reco::CaloJetCollection> l2Isolated;
00066    iEvent.getByLabel(l2Isolated_,l2Isolated);
00067 
00068    if (!l2TauInfoAssoc.isValid() || !l2Isolated.isValid())
00069      {
00070       edm::LogInfo("HLTTauCaloDQMOfflineSource") << "l2TauInfoAssoc object not found, "
00071       "skipping event"; 
00072     return;
00073      }
00074 
00075 
00076    std::vector<l1extra::L1JetParticleRef> tauCandRefVec;
00077 
00078    if(l2TauInfoAssoc.isValid())//get the Association class
00079      {
00080 
00081        //Lets see if we have MC w matching or real data
00082        if(McInfo.isValid())
00083        //If the Collection exists do work
00084        if(l2TauInfoAssoc->size()>0)
00085          for(L2TauInfoAssociation::const_iterator p = l2TauInfoAssoc->begin();p!=l2TauInfoAssoc->end();++p)
00086            {
00087              //Retrieve The L2TauIsolationInfo Class from the AssociationMap
00088              const L2TauIsolationInfo l2info = p->val;
00089        
00090              //Retrieve the Jet From the AssociationMap
00091              const Jet& jet =*(p->key);
00092   
00093              if((doRef_&&match(jet,*McInfo))||(!doRef_))
00094                  {
00095                    ecalIsolEt->Fill(l2info.ECALIsolConeCut);
00096                    towerIsolEt->Fill(l2info.TowerIsolConeCut);
00097                    nClusters->Fill(l2info.ECALClusterNClusters);
00098                    seedTowerEt->Fill(l2info.SeedTowerEt);
00099                    clusterEtaRMS->Fill(l2info.ECALClusterEtaRMS);
00100                    clusterPhiRMS->Fill(l2info.ECALClusterPhiRMS);
00101                    clusterDeltaRRMS->Fill(l2info.ECALClusterDRRMS);
00102                    jetEt->Fill(jet.et());
00103                    jetEta->Fill(jet.eta());
00104                    jetPhi->Fill(jet.phi());
00105               
00106                    EtEffDenom->Fill(jet.et());
00107 
00108                    if(l2Isolated.isValid())
00109                      {
00110                        if(matchJet(jet,*l2Isolated)) 
00111                             EtEffNum->Fill(jet.et());
00112 
00113                      }
00114 
00115                  }
00116 
00117            
00118            } 
00119                
00120      }
00121 
00122    //Plot the missing Et. To be used in SingleTau mainly
00123    Handle<CaloMETCollection> met; iEvent.getByLabel(met_,met);
00124    if(met.isValid())//get the Association class
00125      {
00126        MET->Fill((*met)[0].pt());
00127      }
00128 
00129 
00130 }
00131 
00132 
00133 
00134 void 
00135 HLTTauCaloDQMOfflineSource::beginJob(const edm::EventSetup&)
00136 {
00137 
00138 }
00139 
00140 
00141 void 
00142 HLTTauCaloDQMOfflineSource::endJob() {
00143  
00144   //Get Efficiency
00145 
00146   EtEffNum->getTH1F()->Sumw2();
00147   EtEffDenom->getTH1F()->Sumw2();
00148 
00149   //Write file
00150   if(outFile_.size()>0)
00151     if (store) store->save (outFile_);
00152 
00153 }
00154 
00155 bool 
00156 HLTTauCaloDQMOfflineSource::match(const reco::Jet& jet,const LVColl& McInfo)
00157 {
00158 
00159   //Loop On the Collection and see if your tau jet is matched to one there
00160  //Also find the nearest Matched MC Particle to your Jet (to be complete)
00161  
00162  bool matched=false;
00163 
00164  if(McInfo.size()>0)
00165   for(std::vector<LV>::const_iterator it = McInfo.begin();it!=McInfo.end();++it)
00166    {
00167           double delta = ROOT::Math::VectorUtil::DeltaR(jet.p4().Vect(),*it);
00168           if(delta<matchDeltaRMC_)
00169             {
00170               matched=true;
00171              
00172             }
00173    }
00174 
00175 
00176 
00177  return matched;
00178 }
00179 
00180 bool 
00181 HLTTauCaloDQMOfflineSource::matchJet(const reco::Jet& jet,const reco::CaloJetCollection& McInfo)
00182 {
00183 
00184   //Loop On the Collection and see if your tau jet is matched to one there
00185  //Also find the nearest Matched MC Particle to your Jet (to be complete)
00186  
00187  bool matched=false;
00188 
00189  if(McInfo.size()>0)
00190   for(reco::CaloJetCollection::const_iterator it = McInfo.begin();it!=McInfo.end();++it)
00191    {
00192           double delta = ROOT::Math::VectorUtil::DeltaR(jet.p4().Vect(),it->p4().Vect());
00193           if(delta<matchDeltaRMC_)
00194             {
00195               matched=true;          
00196             }
00197    }
00198 
00199 
00200 
00201  return matched;
00202 }
00203 
00204 

Generated on Tue Jun 9 17:34:10 2009 for CMSSW by  doxygen 1.5.4