CMS 3D CMS Logo

/data/doxygen/doxygen-1.7.3/gen/CMSSW_4_2_8/src/DQMOffline/Trigger/src/EgHLTOffHelper.cc

Go to the documentation of this file.
00001 #include "DQMOffline/Trigger/interface/EgHLTOffHelper.h"
00002 
00003 
00004 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
00005 #include "DataFormats/DetId/interface/DetId.h"
00006 #include "DataFormats/JetReco/interface/CaloJetCollection.h"
00007 #include "FWCore/Common/interface/TriggerNames.h"
00008 
00009 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterTools.h"
00010 
00011 #include "DataFormats/TrackReco/interface/Track.h"
00012 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
00013 
00014 #include "RecoEgamma/EgammaHLTAlgos/interface/EgammaHLTTrackIsolation.h"
00015 #include "RecoEgamma/EgammaIsolationAlgos/interface/EgammaTowerIsolation.h"
00016 #include "RecoEgamma/EgammaIsolationAlgos/interface/EgammaRecHitIsolation.h"
00017 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgo.h"
00018 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgoRcd.h"
00019 
00020 #include "DQMOffline/Trigger/interface/EgHLTTrigCodes.h"
00021 #include "DQMOffline/Trigger/interface/EgHLTTrigTools.h"
00022 #include "DQMOffline/Trigger/interface/EgHLTErrCodes.h"
00023 
00024 #include <iostream>
00025 
00026 using namespace egHLT;
00027 
00028 
00029 
00030 OffHelper::~OffHelper()
00031 {
00032   if(hltEleTrkIsolAlgo_) delete hltEleTrkIsolAlgo_;
00033   if(hltPhoTrkIsolAlgo_) delete hltPhoTrkIsolAlgo_;
00034 }
00035 
00036 void OffHelper::setup(const edm::ParameterSet& conf)
00037 {
00038 
00039   ecalRecHitsEBTag_ = conf.getParameter<edm::InputTag>("BarrelRecHitCollection");
00040   ecalRecHitsEETag_ = conf.getParameter<edm::InputTag>("EndcapRecHitCollection");
00041   caloJetsTag_ = conf.getParameter<edm::InputTag>("CaloJetCollection");
00042   isolTrkTag_ = conf.getParameter<edm::InputTag>("IsolTrackCollection");
00043   hbheHitsTag_ = conf.getParameter<edm::InputTag>("HBHERecHitCollection");
00044   hfHitsTag_ = conf.getParameter<edm::InputTag>("HFRecHitCollection");
00045   electronsTag_ = conf.getParameter<edm::InputTag>("ElectronCollection");
00046   photonsTag_ = conf.getParameter<edm::InputTag>("PhotonCollection");
00047   triggerSummaryLabel_ = conf.getParameter<edm::InputTag>("triggerSummaryLabel");
00048   hltTag_ = conf.getParameter<std::string>("hltTag");
00049   beamSpotTag_ = conf.getParameter<edm::InputTag>("BeamSpotProducer");
00050   caloTowersTag_ = conf.getParameter<edm::InputTag>("CaloTowers");
00051   trigResultsTag_ = conf.getParameter<edm::InputTag>("TrigResults");
00052 
00053   eleCuts_.setup(conf.getParameter<edm::ParameterSet>("eleCuts"));
00054   eleLooseCuts_.setup(conf.getParameter<edm::ParameterSet>("eleLooseCuts"));
00055   phoCuts_.setup(conf.getParameter<edm::ParameterSet>("phoCuts"));
00056   phoLooseCuts_.setup(conf.getParameter<edm::ParameterSet>("phoLooseCuts"));
00057  
00058   //now we have the isolations completely configurable via python
00059   hltEMIsolOuterCone_ = conf.getParameter<double>("hltEMIsolOuterCone");
00060   hltEMIsolInnerConeEB_ = conf.getParameter<double>("hltEMIsolInnerConeEB");
00061   hltEMIsolEtaSliceEB_ = conf.getParameter<double>("hltEMIsolEtaSliceEB");
00062   hltEMIsolEtMinEB_ = conf.getParameter<double>("hltEMIsolEtMinEB");
00063   hltEMIsolEMinEB_ = conf.getParameter<double>("hltEMIsolEMinEB");
00064   hltEMIsolInnerConeEE_ = conf.getParameter<double>("hltEMIsolInnerConeEE");
00065   hltEMIsolEtaSliceEE_ = conf.getParameter<double>("hltEMIsolEtaSliceEE");
00066   hltEMIsolEtMinEE_ = conf.getParameter<double>("hltEMIsolEtMinEE");
00067   hltEMIsolEMinEE_ = conf.getParameter<double>("hltEMIsolEMinEE");
00068 
00069   hltPhoTrkIsolPtMin_ = conf.getParameter<double>("hltPhoTrkIsolPtMin");
00070   hltPhoTrkIsolOuterCone_ = conf.getParameter<double>("hltPhoTrkIsolOuterCone");
00071   hltPhoTrkIsolInnerCone_ = conf.getParameter<double>("hltPhoTrkIsolInnerCone"); 
00072   hltPhoTrkIsolZSpan_ = conf.getParameter<double>("hltPhoTrkIsolZSpan");
00073   hltPhoTrkIsolRSpan_ = conf.getParameter<double>("hltPhoTrkIsolZSpan");
00074   hltPhoTrkIsolCountTrks_ = conf.getParameter<bool>("hltPhoTrkIsolCountTrks");
00075 
00076   hltEleTrkIsolPtMin_ = conf.getParameter<double>("hltEleTrkIsolPtMin");
00077   hltEleTrkIsolOuterCone_ = conf.getParameter<double>("hltEleTrkIsolOuterCone");
00078   hltEleTrkIsolInnerCone_ = conf.getParameter<double>("hltEleTrkIsolInnerCone"); 
00079   hltEleTrkIsolZSpan_ = conf.getParameter<double>("hltEleTrkIsolZSpan");
00080   hltEleTrkIsolRSpan_ = conf.getParameter<double>("hltEleTrkIsolZSpan");
00081 
00082   hltHadIsolOuterCone_ = conf.getParameter<double>("hltHadIsolOuterCone");
00083   hltHadIsolInnerCone_ = conf.getParameter<double>("hltHadIsolInnerCone");
00084   hltHadIsolEtMin_ = conf.getParameter<double>("hltHadIsolEtMin");
00085   hltHadIsolDepth_ = conf.getParameter<int>("hltHadIsolDepth"); 
00086 
00087   calHLTHcalIsol_ = conf.getParameter<bool>("calHLTHcalIsol");
00088   calHLTEmIsol_ = conf.getParameter<bool>("calHLTEmIsol");
00089   calHLTEleTrkIsol_ = conf.getParameter<bool>("calHLTEleTrkIsol");
00090   calHLTPhoTrkIsol_ = conf.getParameter<bool>("calHLTPhoTrkIsol");
00091  
00092   trigCutParams_ = conf.getParameter<std::vector<edm::ParameterSet> >("triggerCuts"); //setupTriggers used to be in this function but had to be moved due to HLTConfigChanges (has to be called beginRun) so we have to save this for later.
00093 
00094   hltEleTrkIsolAlgo_ = new EgammaHLTTrackIsolation(hltEleTrkIsolPtMin_,hltEleTrkIsolOuterCone_,hltEleTrkIsolZSpan_,hltEleTrkIsolRSpan_,hltEleTrkIsolInnerCone_);
00095   hltPhoTrkIsolAlgo_ = new EgammaHLTTrackIsolation(hltPhoTrkIsolPtMin_,hltPhoTrkIsolOuterCone_,hltPhoTrkIsolZSpan_,hltPhoTrkIsolRSpan_,hltPhoTrkIsolInnerCone_);
00096                                                
00097 
00098 }
00099 
00100 //this code was taken out of OffHelper::setup due to HLTConfigProvider changes
00101 //it still assumes that this is called only once
00102 void OffHelper::setupTriggers(const HLTConfigProvider& hltConfig,const std::vector<std::string>& hltFiltersUsed)
00103 {
00104   hltFiltersUsed_ = hltFiltersUsed; //expensive but only do this once and faster ways could make things less clear
00105   //now work out how many objects are requires to pass filter for it to accept
00106   hltFiltersUsedWithNrCandsCut_.clear();
00107   for(size_t filterNr=0;filterNr<hltFiltersUsed_.size();filterNr++){
00108     hltFiltersUsedWithNrCandsCut_.push_back(std::make_pair(hltFiltersUsed_[filterNr],egHLT::trigTools::getMinNrObjsRequiredByFilter(hltFiltersUsed_[filterNr])));
00109   }
00110 
00111   //now loading the cuts for every trigger into our vector which stores them
00112   //only load cuts for triggers that are in hltFiltersUsed
00113   
00114   for(size_t trigNr=0;trigNr<trigCutParams_.size();trigNr++) {
00115     std::string trigName = trigCutParams_[trigNr].getParameter<std::string>("trigName");
00116     if(std::find(hltFiltersUsed_.begin(),hltFiltersUsed_.end(),trigName)!=hltFiltersUsed_.end()){ //perhaps I should sort hltFiltersUsed_....
00117       trigCuts_.push_back(std::make_pair(TrigCodes::getCode(trigName),OffEgSel(trigCutParams_[trigNr])));
00118       //   std::cout<<trigName<<std::endl<<"between"<<std::endl<<trigCutParams_[trigNr]<<std::endl<<"after"<<std::endl;
00119     }
00120   }
00121   trigCutParams_.clear();//dont need it any more, get rid of it
00122 
00123   //to make my life difficult, the scaled l1 paths are special
00124   //and arent stored in trigger event
00125   //to I have to figure out the path, see if it passes
00126   //and then hunt down the l1 seed filter and use that to match to the pho/ele
00127   //matching on l1 seed filter is not enough as that will be passed for normal 
00128   //electron triggers even if pre-scale hasnt fired
00129   l1PreScaledFilters_.clear();
00130   l1PreScaledPaths_.clear();
00131   l1PreAndSeedFilters_.clear();
00132   for(size_t filterNr=0;filterNr<hltFiltersUsed_.size();filterNr++){
00133     if(hltFiltersUsed_[filterNr].find("hltPreL1")==0){ //l1 prescaled path
00134       l1PreScaledFilters_.push_back(hltFiltersUsed_[filterNr]);
00135     }
00136   }
00137 
00138   egHLT::trigTools::translateFiltersToPathNames(hltConfig,l1PreScaledFilters_,l1PreScaledPaths_);
00139   if(l1PreScaledPaths_.size()==l1PreScaledFilters_.size()){
00140     for(size_t pathNr=0;pathNr<l1PreScaledPaths_.size();pathNr++){
00141      
00142       std::string l1SeedFilter =egHLT::trigTools::getL1SeedFilterOfPath(hltConfig,l1PreScaledPaths_[pathNr]);
00143       //---Morse====
00144       //std::cout<<l1PreScaledFilters_[pathNr]<<"  "<<l1PreScaledPaths_[pathNr]<<"  "<<l1SeedFilter<<std::endl;
00145       //------------
00146       l1PreAndSeedFilters_.push_back(std::make_pair(l1PreScaledFilters_[pathNr],l1SeedFilter));
00147     }
00148   }
00149 }
00150 
00151 int OffHelper::makeOffEvt(const edm::Event& edmEvent,const edm::EventSetup& setup,egHLT::OffEvt& offEvent)
00152 {
00153   offEvent.clear();
00154   int errCode=0; //excution stops as soon as an error is flagged
00155   if(errCode==0) errCode = getHandles(edmEvent,setup);
00156   if(errCode==0) errCode = fillOffEleVec(offEvent.eles());
00157   if(errCode==0) errCode = fillOffPhoVec(offEvent.phos());
00158   if(errCode==0) errCode =  setTrigInfo(edmEvent, offEvent);
00159   if(errCode==0) offEvent.setJets(recoJets_);
00160   return errCode;
00161 }
00162 
00163 int OffHelper::getHandles(const edm::Event& event,const edm::EventSetup& setup)
00164 {
00165   try { 
00166     setup.get<CaloGeometryRecord>().get(caloGeom_);
00167     setup.get<CaloTopologyRecord>().get(caloTopology_);
00168     //setup.get<EcalSeverityLevelAlgoRcd>().get(ecalSeverityLevel_);
00169   }catch(...){
00170     return errCodes::Geom;
00171   }
00172   try {
00173     setup.get<IdealMagneticFieldRecord>().get(magField_);
00174   }catch(...){
00175     return errCodes::MagField;
00176   }
00177 
00178   //get objects
00179   if(!getHandle(event,triggerSummaryLabel_,trigEvt_)) return errCodes::TrigEvent; //must have this, otherwise skip event
00180   if(!getHandle(event,trigResultsTag_,trigResults_)) return errCodes::TrigEvent; //re using bit to minimise bug fix code changes
00181   if(!getHandle(event,electronsTag_,recoEles_)) return errCodes::OffEle; //need for electrons
00182   if(!getHandle(event,photonsTag_, recoPhos_)) return errCodes::OffPho; //need for photons
00183   if(!getHandle(event,caloJetsTag_,recoJets_)) return errCodes::OffJet; //need for electrons and photons
00184 
00185   //need for HLT isolations (rec hits also need for sigmaIPhiIPhi (ele/pho) and r9 pho)
00186   if(!getHandle(event,ecalRecHitsEBTag_,ebRecHits_)) return errCodes::EBRecHits;
00187   if(!getHandle(event,ecalRecHitsEETag_,eeRecHits_)) return errCodes::EERecHits;
00188   if(!getHandle(event,isolTrkTag_,isolTrks_)) return errCodes::IsolTrks;
00189   if(!getHandle(event,hbheHitsTag_, hbheHits_)) return errCodes::HBHERecHits; //I dont think we need hbhe rec-hits any more
00190   if(!getHandle(event,hfHitsTag_, hfHits_)) return errCodes::HFRecHits;//I dont think we need hf rec-hits any more
00191   if(!getHandle(event,beamSpotTag_,beamSpot_)) return errCodes::BeamSpot;
00192   if(!getHandle(event,caloTowersTag_,caloTowers_)) return errCodes::CaloTowers;
00193 
00194   
00195   return 0;
00196 }
00197 
00198 //this function coverts GsfElectrons to a format which is actually useful to me
00199 int OffHelper::fillOffEleVec(std::vector<OffEle>& egHLTOffEles)
00200 {
00201   egHLTOffEles.clear();
00202   egHLTOffEles.reserve(recoEles_->size());
00203   for(reco::GsfElectronCollection::const_iterator gsfIter=recoEles_->begin(); gsfIter!=recoEles_->end();++gsfIter){
00204     if(!gsfIter->ecalDrivenSeed()) continue; //avoid PF electrons (this is Eg HLT validation and HLT is ecal driven)
00205 
00206     OffEle::IsolData isolData;   
00207     fillIsolData(*gsfIter,isolData);
00208     
00209     OffEle::ClusShapeData clusShapeData;
00210     fillClusShapeData(*gsfIter,clusShapeData);
00211 
00212     OffEle::HLTData hltData;
00213     fillHLTData(*gsfIter,hltData);
00214 
00215     egHLTOffEles.push_back(OffEle(*gsfIter,clusShapeData,isolData,hltData));
00216     
00217     //now we would like to set the cut results
00218     OffEle& ele =  egHLTOffEles.back();
00219     ele.setCutCode(eleCuts_.getCutCode(ele));
00220     ele.setLooseCutCode(eleLooseCuts_.getCutCode(ele));
00221     
00222     std::vector<std::pair<TrigCodes::TrigBitSet,int> >trigCutsCutCodes;
00223     for(size_t i=0;i<trigCuts_.size();i++) trigCutsCutCodes.push_back(std::make_pair(trigCuts_[i].first,trigCuts_[i].second.getCutCode(ele)));
00224     ele.setTrigCutsCutCodes(trigCutsCutCodes);
00225   }//end loop over gsf electron collection
00226   return 0;
00227 }
00228 
00229 void OffHelper::fillIsolData(const reco::GsfElectron& ele,OffEle::IsolData& isolData)
00230 {
00231   EgammaTowerIsolation hcalIsolAlgo(hltHadIsolOuterCone_,hltHadIsolInnerCone_,hltHadIsolEtMin_,hltHadIsolDepth_,caloTowers_.product());
00232   EcalRecHitMetaCollection ebHits(*ebRecHits_);
00233   EcalRecHitMetaCollection eeHits(*eeRecHits_);
00234   EgammaRecHitIsolation ecalIsolAlgoEB(hltEMIsolOuterCone_,hltEMIsolInnerConeEB_,hltEMIsolEtaSliceEB_,
00235                                        hltEMIsolEtMinEB_,hltEMIsolEMinEB_,caloGeom_,&ebHits,ecalSeverityLevel_.product(),DetId::Ecal);
00236   EgammaRecHitIsolation ecalIsolAlgoEE(hltEMIsolOuterCone_,hltEMIsolInnerConeEE_,hltEMIsolEtaSliceEE_,
00237                                        hltEMIsolEtMinEE_,hltEMIsolEMinEE_,caloGeom_,&eeHits,ecalSeverityLevel_.product(),DetId::Ecal);
00238   
00239   isolData.ptTrks=ele.dr03TkSumPt();
00240   isolData.nrTrks=999; //no longer supported
00241   isolData.em= ele.dr03EcalRecHitSumEt();
00242   isolData.hadDepth1 =  ele.dr03HcalDepth1TowerSumEt();
00243   isolData.hadDepth2 =  ele.dr03HcalDepth2TowerSumEt();   
00244 
00245   //now time to do the HLT algos
00246   if(calHLTHcalIsol_) isolData.hltHad=hcalIsolAlgo.getTowerESum(&ele);
00247   else isolData.hltHad = 0.;
00248   if(calHLTEleTrkIsol_) isolData.hltTrksEle=hltEleTrkIsolAlgo_->electronPtSum(&(*(ele.gsfTrack())),isolTrks_.product());
00249   else isolData.hltTrksEle = 0.;
00250   if(calHLTPhoTrkIsol_){
00251     if(hltPhoTrkIsolCountTrks_) isolData.hltTrksPho=hltPhoTrkIsolAlgo_->photonTrackCount(&ele,isolTrks_.product(),false);
00252     else isolData.hltTrksPho=hltPhoTrkIsolAlgo_->photonPtSum(&ele,isolTrks_.product(),false);
00253   }
00254   else isolData.hltTrksPho = 0.;
00255   if(calHLTEmIsol_) isolData.hltEm = ecalIsolAlgoEB.getEtSum(&ele) + 
00256                                      ecalIsolAlgoEE.getEtSum(&ele);
00257   else isolData.hltEm = 0.;
00258   
00259 }
00260 
00261 
00262 void OffHelper::fillClusShapeData(const reco::GsfElectron& ele,OffEle::ClusShapeData& clusShapeData)
00263 {
00264   clusShapeData.sigmaEtaEta = ele.sigmaEtaEta();
00265   clusShapeData.sigmaIEtaIEta = ele.sigmaIetaIeta();
00266   double e5x5 = ele.e5x5();
00267   if(e5x5!=0.){
00268     clusShapeData.e1x5Over5x5 = ele.e1x5()/e5x5;
00269     clusShapeData.e2x5MaxOver5x5 = ele.e2x5Max()/e5x5;
00270   }else{
00271     clusShapeData.e1x5Over5x5 = -1;
00272     clusShapeData.e2x5MaxOver5x5 = -1;
00273   }
00274   
00275   //want to calculate r9, sigmaPhiPhi and sigmaIPhiIPhi, have to do old fashioned way
00276   const reco::BasicCluster& seedClus = *(ele.superCluster()->seed());
00277   const DetId seedDetId = seedClus.hitsAndFractions()[0].first; //note this may not actually be the seed hit but it doesnt matter because all hits will be in the barrel OR endcap
00278   if(seedDetId.subdetId()==EcalBarrel){
00279     std::vector<float> stdCov = EcalClusterTools::covariances(seedClus,ebRecHits_.product(),caloTopology_.product(),caloGeom_.product());
00280     std::vector<float> crysCov = EcalClusterTools::localCovariances(seedClus,ebRecHits_.product(),caloTopology_.product());
00281     clusShapeData.sigmaPhiPhi = sqrt(stdCov[2]);
00282       clusShapeData.sigmaIPhiIPhi =  sqrt(crysCov[2]);
00283       if(ele.superCluster()->rawEnergy()!=0.){
00284         clusShapeData.r9 = EcalClusterTools::e3x3(seedClus,ebRecHits_.product(),caloTopology_.product()) / ele.superCluster()->rawEnergy();
00285       }else clusShapeData.r9 = -1.;
00286        
00287   }else{
00288     std::vector<float> stdCov = EcalClusterTools::covariances(seedClus,eeRecHits_.product(),caloTopology_.product(),caloGeom_.product()); 
00289     std::vector<float> crysCov = EcalClusterTools::localCovariances(seedClus,eeRecHits_.product(),caloTopology_.product());
00290     clusShapeData.sigmaPhiPhi = sqrt(stdCov[2]);
00291     clusShapeData.sigmaIPhiIPhi = sqrt(crysCov[2]);
00292     if(ele.superCluster()->rawEnergy()!=0.){
00293         clusShapeData.r9 = EcalClusterTools::e3x3(seedClus,eeRecHits_.product(),caloTopology_.product()) / ele.superCluster()->rawEnergy();
00294     }else clusShapeData.r9 = -1.;
00295   } 
00296 }
00297 
00298 //reco approximations of hlt quantities
00299 void OffHelper::fillHLTData(const reco::GsfElectron& ele,OffEle::HLTData& hltData)
00300 {
00301   if(ele.closestCtfTrackRef().isNonnull() && 
00302      ele.closestCtfTrackRef()->extra().isNonnull()){
00303     reco::TrackRef ctfTrack = ele.closestCtfTrackRef();
00304     reco::SuperClusterRef scClus = ele.superCluster();
00305 
00306     //dEta
00307     const reco::BeamSpot::Point& bsPos = beamSpot_->position();     
00308     math::XYZPoint scPosWRTVtx(scClus->x()-bsPos.x(), scClus->y()-bsPos.y() , scClus->z()-ctfTrack->vz());
00309     hltData.dEtaIn = fabs(scPosWRTVtx.eta()-ctfTrack->eta());
00310 
00311     //dPhi: lifted straight from hlt code
00312     float deltaPhi=fabs(ctfTrack->outerPosition().phi()-scClus->phi());
00313     if(deltaPhi>6.283185308) deltaPhi -= 6.283185308;
00314     if(deltaPhi>3.141592654) deltaPhi = 6.283185308-deltaPhi;
00315     hltData.dPhiIn = deltaPhi;
00316     
00317     //invEInvP
00318     if(ele.ecalEnergy()!=0 && ctfTrack->p()!=0) hltData.invEInvP= 1/ele.ecalEnergy() - 1/ctfTrack->p();
00319     else hltData.invEInvP = 0;
00320   }else{
00321     hltData.dEtaIn =999;
00322     hltData.dPhiIn =999;
00323     hltData.invEInvP = 999;
00324 
00325   }
00326 }
00327 
00328 //this function coverts Photons to a format which more useful to me
00329 int OffHelper::fillOffPhoVec(std::vector<OffPho>& egHLTOffPhos)
00330 {
00331   egHLTOffPhos.clear();
00332   egHLTOffPhos.reserve(recoPhos_->size());
00333   for(reco::PhotonCollection::const_iterator phoIter=recoPhos_->begin(); phoIter!=recoPhos_->end();++phoIter){
00334 
00335     OffPho::IsolData isolData;  
00336     OffPho::ClusShapeData clusShapeData;
00337   
00338     fillIsolData(*phoIter,isolData);
00339     fillClusShapeData(*phoIter,clusShapeData);
00340  
00341     egHLTOffPhos.push_back(OffPho(*phoIter,clusShapeData,isolData));
00342     OffPho& pho =  egHLTOffPhos.back();
00343     pho.setCutCode(phoCuts_.getCutCode(pho));
00344     pho.setLooseCutCode(phoLooseCuts_.getCutCode(pho));
00345 
00346     std::vector<std::pair<TrigCodes::TrigBitSet,int> >trigCutsCutCodes;
00347     for(size_t i=0;i<trigCuts_.size();i++) trigCutsCutCodes.push_back(std::make_pair(trigCuts_[i].first,trigCuts_[i].second.getCutCode(pho)));
00348     pho.setTrigCutsCutCodes(trigCutsCutCodes); 
00349 
00350 
00351   }//end loop over photon collection
00352   return 0;
00353 }
00354 
00355 
00356 void OffHelper::fillIsolData(const reco::Photon& pho,OffPho::IsolData& isolData)
00357 {
00358   EgammaTowerIsolation hcalIsolAlgo(hltHadIsolOuterCone_,hltHadIsolInnerCone_,hltHadIsolEtMin_,hltHadIsolDepth_,caloTowers_.product());
00359   EcalRecHitMetaCollection ebHits(*ebRecHits_);
00360   EcalRecHitMetaCollection eeHits(*ebRecHits_);
00361   EgammaRecHitIsolation ecalIsolAlgoEB(hltEMIsolOuterCone_,hltEMIsolInnerConeEB_,hltEMIsolEtaSliceEB_,
00362                                        hltEMIsolEtMinEB_,hltEMIsolEMinEB_,caloGeom_,&ebHits,ecalSeverityLevel_.product(),DetId::Ecal);
00363   EgammaRecHitIsolation ecalIsolAlgoEE(hltEMIsolOuterCone_,hltEMIsolInnerConeEE_,hltEMIsolEtaSliceEE_,
00364                                        hltEMIsolEtMinEE_,hltEMIsolEMinEE_,caloGeom_,&eeHits,ecalSeverityLevel_.product(),DetId::Ecal);
00365   
00366   isolData.nrTrks = pho.nTrkHollowConeDR03();
00367   isolData.ptTrks = pho.trkSumPtHollowConeDR03();
00368   isolData.em = pho.ecalRecHitSumEtConeDR03();
00369   isolData.had = pho.hcalTowerSumEtConeDR03();  
00370 
00371   //now calculate hlt algos
00372   if(calHLTHcalIsol_) isolData.hltHad=hcalIsolAlgo.getTowerESum(&pho);
00373   else isolData.hltHad = 0.;
00374   if(calHLTPhoTrkIsol_){
00375     if(hltPhoTrkIsolCountTrks_) isolData.hltTrks=hltPhoTrkIsolAlgo_->photonTrackCount(&pho,isolTrks_.product(),false);
00376     else isolData.hltTrks=hltPhoTrkIsolAlgo_->photonPtSum(&pho,isolTrks_.product(),false);
00377   }
00378   else isolData.hltTrks = 0.;
00379   if(calHLTEmIsol_) isolData.hltEm = ecalIsolAlgoEB.getEtSum(&pho) + 
00380                                      ecalIsolAlgoEE.getEtSum(&pho);
00381   else isolData.hltEm = 0.;
00382   
00383 }
00384 
00385 void OffHelper::fillClusShapeData(const reco::Photon& pho,OffPho::ClusShapeData& clusShapeData)
00386 {
00387   clusShapeData.sigmaEtaEta = pho.sigmaEtaEta();
00388   clusShapeData.sigmaIEtaIEta =  pho.sigmaIetaIeta();
00389   double e5x5 =  pho.e5x5();
00390   if(e5x5!=0.){ //even though it is almost impossible for this to be 0., this code can never ever crash under any situation
00391     clusShapeData.e1x5Over5x5 =  pho.e1x5()/e5x5;
00392     clusShapeData.e2x5MaxOver5x5 = pho.e2x5()/e5x5;
00393   }else{
00394     clusShapeData.e1x5Over5x5 = -1;
00395     clusShapeData.e2x5MaxOver5x5 = -1;
00396   }
00397   clusShapeData.r9 = pho.r9();
00398   
00399   //sigmaPhiPhi and sigmaIPhiIPhi are not in object (and nor should they be) so have to get them old fashioned way
00400   //need to figure out if its in the barrel or endcap
00401   //get the first hit of the cluster and figure out if its barrel or endcap 
00402   const reco::BasicCluster& seedClus = *(pho.superCluster()->seed());
00403   const DetId seedDetId = seedClus.hitsAndFractions()[0].first; //note this may not actually be the seed hit but it doesnt matter because all hits will be in the barrel OR endcap (it is also incredably inefficient as it getHitsByDetId passes the vector by value not reference
00404   if(seedDetId.subdetId()==EcalBarrel){
00405     std::vector<float> stdCov = EcalClusterTools::covariances(seedClus,ebRecHits_.product(),caloTopology_.product(),caloGeom_.product());
00406     std::vector<float> crysCov = EcalClusterTools::localCovariances(seedClus,ebRecHits_.product(),caloTopology_.product());
00407     clusShapeData.sigmaPhiPhi = sqrt(stdCov[2]);
00408     clusShapeData.sigmaIPhiIPhi =  sqrt(crysCov[2]);
00409     }else{
00410     std::vector<float> stdCov = EcalClusterTools::covariances(seedClus,eeRecHits_.product(),caloTopology_.product(),caloGeom_.product()); 
00411     std::vector<float> crysCov = EcalClusterTools::localCovariances(seedClus,eeRecHits_.product(),caloTopology_.product());
00412     
00413     clusShapeData.sigmaPhiPhi = sqrt(stdCov[2]);
00414     clusShapeData.sigmaIPhiIPhi = sqrt(crysCov[2]); 
00415   }
00416 }  
00417 
00418 int OffHelper::setTrigInfo(const edm::Event & edmEvent, egHLT::OffEvt& offEvent)
00419 {
00420   TrigCodes::TrigBitSet evtTrigBits = trigTools::getFiltersPassed(hltFiltersUsedWithNrCandsCut_,trigEvt_.product(),hltTag_);
00421   //the l1 prescale paths dont have a filter with I can figure out if it passed or failed with so have to use TriggerResults
00422   if(l1PreScaledPaths_.size()==l1PreScaledFilters_.size()){ //check to ensure both vectors have same number of events incase of screw ups     
00423     const edm::TriggerNames & triggerNames = edmEvent.triggerNames(*trigResults_);
00424     for(size_t pathNr=0;pathNr<l1PreScaledPaths_.size();pathNr++){ //now we have to check the prescaled l1 trigger paths
00425       unsigned int pathIndex = triggerNames.triggerIndex(l1PreScaledPaths_[pathNr]);
00426       if(pathIndex<trigResults_->size() && trigResults_->accept(pathIndex)){
00427         evtTrigBits |=TrigCodes::getCode(l1PreScaledFilters_[pathNr]);
00428       }    
00429     }
00430   }
00431 
00432   offEvent.setEvtTrigBits(evtTrigBits);
00433 
00434   trigTools::setFiltersObjPasses(offEvent.eles(),hltFiltersUsed_,l1PreAndSeedFilters_,evtTrigBits,trigEvt_.product(),hltTag_);
00435   trigTools::setFiltersObjPasses(offEvent.phos(),hltFiltersUsed_,l1PreAndSeedFilters_,evtTrigBits,trigEvt_.product(),hltTag_); 
00436   return 0;
00437 }