CMSSW Reference Manual

00001 
00002 //===========Original Author: Matthias Mozer
00003 //===========Modified for purpose of Tau Group: Konstantinos A. Petridis
00004 
00005 
00006 #include "DQMOffline/Trigger/interface/HLTTauElDQMOfflineSource.h"
00007 
00008 // Collaborating Class Header
00009 #include "FWCore/Framework/interface/Frameworkfwd.h"
00010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00011 #include "FWCore/Framework/interface/MakerMacros.h"
00012 #include "DataFormats/HLTReco/interface/TriggerEventWithRefs.h"
00013 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
00014 #include "DataFormats/EgammaCandidates/interface/Electron.h"
00015 #include "DataFormats/EgammaCandidates/interface/PhotonFwd.h"
00016 #include "DataFormats/EgammaCandidates/interface/Photon.h"
00017 #include "DataFormats/L1Trigger/interface/L1EmParticle.h"
00018 #include "DataFormats/L1Trigger/interface/L1EmParticleFwd.h"
00019 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00020 #include "DataFormats/Common/interface/AssociationMap.h"
00021 
00022 #include "DataFormats/Common/interface/Handle.h"
00023 #include "DataFormats/Common/interface/RefToBase.h"
00024 #include "FWCore/ServiceRegistry/interface/Service.h"
00025 #include "PhysicsTools/UtilAlgos/interface/TFileService.h"
00026 #include "FWCore/Utilities/interface/Exception.h"
00027 
00028 #include "TFile.h"
00029 #include "TDirectory.h"
00030 #include "TH1F.h"
00031 #include <iostream>
00032 #include <string>
00033 #include <Math/VectorUtil.h>
00034 using namespace ROOT::Math::VectorUtil ;
00035 
00037 HLTTauElDQMOfflineSource::HLTTauElDQMOfflineSource(const edm::ParameterSet& pset)  
00038 {
00039 
00040   //paramters for generator study
00041 
00042   refCollection_=pset.getUntrackedParameter<edm::InputTag>("refCollection");
00043   reqNum_ = pset.getParameter<unsigned int>("reqNum");
00044   pdgGen_ =  pset.getParameter<int>("pdgGen");
00045   genEtaAcc_ = pset.getParameter<double>("genEtaAcc");
00046   genEtAcc_ = pset.getParameter<double>("genEtAcc");
00047   outputFile_= pset.getParameter<std::string>("outputFile");
00048   triggerName_=pset.getParameter<std::string>("triggerName");
00049   //plotting paramters
00050   thePtMin_ = pset.getUntrackedParameter<double>("PtMin",0.);
00051   thePtMax_ = pset.getUntrackedParameter<double>("PtMax",1000.);
00052   theNbins_ = pset.getUntrackedParameter<unsigned int>("Nbins",40);  
00053   //info for each filter-step
00054   std::vector<edm::ParameterSet> filters_ = pset.getParameter<std::vector<edm::ParameterSet> >("filters");
00055 
00056   for(std::vector<edm::ParameterSet>::iterator filterconf = filters_.begin() ; filterconf != filters_.end() ; filterconf++){
00057     m_theHLTCollectionLabels.push_back(filterconf->getParameter<edm::InputTag>("HLTCollectionLabels"));
00058     m_theHLTOutputTypes.push_back(filterconf->getParameter<unsigned int>("theHLTOutputTypes"));
00059     std::vector<double> bounds = filterconf->getParameter<std::vector<double> >("PlotBounds");
00060     assert(bounds.size() == 2);
00061     m_plotBounds.push_back(std::pair<double,double>(bounds[0],bounds[1]));
00062     m_isoNames.push_back(filterconf->getParameter<std::vector<edm::InputTag> >("IsoCollections"));
00063     assert(m_isoNames.back().size()>0);
00064     if (m_isoNames.back().at(0).label()=="none")
00065       m_plotiso.push_back(false);
00066     else{
00067       m_plotiso.push_back(true);
00068        //std::cout << "plotting isolation for: " <<  m_isoNames.back().at(0).label() << std::endl;
00069     }
00070 
00071   }
00072 }
00073 
00074 
00075 void HLTTauElDQMOfflineSource::beginJob(const edm::EventSetup&){
00076   std::string histoname="total eff";
00077   DQMStore* store = &*edm::Service<DQMStore>();
00078 
00079   if(store){
00080 
00081 
00082     store->setCurrentFolder(triggerName_);
00083     m_total = store->book1D(histoname.c_str(),histoname.c_str(),m_theHLTCollectionLabels.size()+2,0,m_theHLTCollectionLabels.size()+2);
00084     m_total->setBinLabel(m_theHLTCollectionLabels.size()+1,"Total");
00085     m_total->setBinLabel(m_theHLTCollectionLabels.size()+2,"Gen");
00086     for (unsigned int u=0; u<m_theHLTCollectionLabels.size(); u++){m_total->setBinLabel(u+1,m_theHLTCollectionLabels[u].label().c_str());}
00087     
00088     MonitorElement* tmphisto;
00089     
00090     histoname = "gen et";
00091     m_etgen =  store->book1D(histoname.c_str(),histoname.c_str(),theNbins_,thePtMin_,thePtMax_);
00092     histoname = "gen eta";
00093     m_etagen = store->book1D(histoname.c_str(),histoname.c_str(),theNbins_,-2.7,2.7);
00094     
00095     for(unsigned int i = 0; i< m_theHLTCollectionLabels.size() ; i++){
00096       histoname = m_theHLTCollectionLabels[i].label()+"et";
00097       tmphisto = store->book1D(histoname.c_str(),histoname.c_str(),theNbins_,thePtMin_,thePtMax_);
00098       m_ethist.push_back(tmphisto);
00099       
00100       histoname = m_theHLTCollectionLabels[i].label()+"eta";
00101       tmphisto =  store->book1D(histoname.c_str(),histoname.c_str(),theNbins_,-2.7,2.7);
00102       m_etahist.push_back(tmphisto);          
00103       
00104       histoname = m_theHLTCollectionLabels[i].label()+"etMatched";
00105       tmphisto = store->book1D(histoname.c_str(),histoname.c_str(),theNbins_,thePtMin_,thePtMax_);
00106       m_ethistmatch.push_back(tmphisto);
00107       
00108       histoname = m_theHLTCollectionLabels[i].label()+"etaMatched";
00109       tmphisto =  store->book1D(histoname.c_str(),histoname.c_str(),theNbins_,-2.7,2.7);
00110       m_etahistmatch.push_back(tmphisto);          
00111     }
00112   }
00113 }
00114 
00116 HLTTauElDQMOfflineSource::~HLTTauElDQMOfflineSource(){
00117 }
00118 
00119 void HLTTauElDQMOfflineSource::analyze(const edm::Event & event , const edm::EventSetup& setup){
00120 
00121 
00122 
00123   // fill L1 and HLT info
00124   // get objects possed by each filter
00125   edm::Handle<trigger::TriggerEventWithRefs> triggerObj;
00126   event.getByLabel("hltTriggerSummaryRAW",triggerObj); 
00127   if(!triggerObj.isValid()) { 
00128     edm::LogWarning("HLTTauElDQMOfflineSource") << "RAW-type HLT results not found, skipping event";
00129     return;
00130   }
00131 // throw(cms::Exception("Release Validation Error")<< "RAW-type HLT results not found" );
00132 
00133 
00134 
00135 //  for(int i = 0; i<triggerObj->size() ;i++ ){
00136 //    std::cout << triggerObj->filterTag(i) << std::endl;
00137 //  }  
00138 
00139   // total event number
00140   m_total->Fill(m_theHLTCollectionLabels.size()+0.5);
00141   // get matching objects
00142   edm::Handle<LVColl> refC;
00143   bool ref_ok=false;
00144   LVColl mcparts;
00145   if(event.getByLabel(refCollection_,refC)){
00146     ref_ok=true;
00147     if(refC->size())m_total->Fill(m_theHLTCollectionLabels.size()+1.5);
00148     for (size_t i=0;i<refC->size();++i)
00149       {
00150         double eta=refC->at(i).Eta();double et=refC->at(i).Et();
00151         m_etgen->Fill(et);
00152         m_etagen->Fill(eta);
00153         mcparts.push_back(refC->at(i));
00154       }
00155   }
00156   
00157   for(unsigned int n=0; n < m_theHLTCollectionLabels.size() ; n++) { //loop over filter modules
00158 
00159     switch(m_theHLTOutputTypes[n]){
00160     case 82: // non-iso L1
00161       fillHistos<l1extra::L1EmParticleCollection>(triggerObj,event,n,mcparts);break;
00162     case 83: // iso L1
00163       fillHistos<l1extra::L1EmParticleCollection>(triggerObj,event,n,mcparts);break;
00164     case 91: //photon 
00165       fillHistos<reco::RecoEcalCandidateCollection>(triggerObj,event,n,mcparts);break;
00166     case 92: //electron 
00167       fillHistos<reco::ElectronCollection>(triggerObj,event,n,mcparts);break;
00168     case 100: // TriggerCluster
00169       fillHistos<reco::RecoEcalCandidateCollection>(triggerObj,event,n,mcparts);break;
00170     default: throw(cms::Exception("Release Validation Error")<< "HLT output type not implemented: theHLTOutputTypes[n]" );
00171     }
00172   }
00173 
00174 }
00175 
00176 template <class T> void HLTTauElDQMOfflineSource::fillHistos(edm::Handle<trigger::TriggerEventWithRefs>& triggerObj,const edm::Event& iEvent ,unsigned int n,LVColl& mcparts){
00177   
00178   std::vector<edm::Ref<T> > recoecalcands;
00179   if (!( triggerObj->filterIndex(m_theHLTCollectionLabels[n])>=triggerObj->size() )){ // only process if availabel
00180     
00181     // retrieve saved filter objects
00182     triggerObj->getObjects(triggerObj->filterIndex(m_theHLTCollectionLabels[n]),m_theHLTOutputTypes[n],recoecalcands);
00183     //Danger: special case, L1 non-isolated
00184     // needs to be merged with L1 iso
00185     if(m_theHLTOutputTypes[n]==82){
00186       std::vector<edm::Ref<T> > isocands;
00187       triggerObj->getObjects(triggerObj->filterIndex(m_theHLTCollectionLabels[n]),83,isocands);
00188       if(isocands.size()>0)
00189         for(unsigned int i=0; i < isocands.size(); i++)
00190           recoecalcands.push_back(isocands[i]);
00191     }
00192     
00193     //fill filter objects into histos
00194     if (recoecalcands.size()&&mcparts.size()){
00195 
00196       if(recoecalcands.size() >= reqNum_) 
00197         m_total->Fill(n+0.5);
00198       for (unsigned int i=0; i<recoecalcands.size(); i++) {
00199         //unmatched
00200         m_ethist[n]->Fill(recoecalcands[i]->et() );
00201         m_etahist[n]->Fill(recoecalcands[i]->eta() );
00202       }
00203     
00204       
00205       for(size_t i=0;i<mcparts.size();++i)
00206         {
00207           double closestDr=1000.;
00208           LV closest;
00209           for (size_t ii=0; ii<recoecalcands.size(); ++ii) {
00210             //matched
00211             math::XYZVector candDir=recoecalcands[ii]->momentum();
00212             double dr = DeltaR(mcparts[i],candDir);
00213             if(dr < closestDr){
00214                 closestDr = dr;
00215                 //closest = mcparts[i];
00216             }
00217           }
00218           if( (closestDr<0.5&&(m_theHLTOutputTypes[n]==82 || m_theHLTOutputTypes[n]==83)) ||
00219               (closestDr<0.1&&m_theHLTOutputTypes[n]>90) ){
00220             float eta   =mcparts[i].Eta();//closest.Eta();
00221             float et    =mcparts[i].Et();//closest.Et();
00222             m_ethistmatch[n]->Fill( et );
00223             m_etahistmatch[n]->Fill( eta );
00224           }
00225         }
00226     }
00227   }
00228 }
00229   
00230 
00231 
00232 void HLTTauElDQMOfflineSource::endJob(){
00233 
00234   if(outputFile_.size()>0)
00235   if (&*edm::Service<DQMStore>()) edm::Service<DQMStore>()->save (outputFile_);
00236   //  total->Scale(1./total->GetBinContent(1));
00237   //for(unsigned int n= m_theHLTCollectionLabels.size()-1 ; n>0;n--){
00238   //  ethist[n]->Divide(ethist[n-1]);
00239   //  etahist[n]->Divide(etahist[n-1]);
00240   //}
00241 }
00242
HLTTauElDQMOfflineSource.cc
