00001 #include "DQMOffline/Trigger/interface/HLTTauRefProducer.h"
00002 #include "TLorentzVector.h"
00003
00004 #include "DataFormats/TauReco/interface/PFTau.h"
00005 #include "DataFormats/TauReco/interface/PFTauDiscriminator.h"
00006 #include "DataFormats/TauReco/interface/CaloTau.h"
00007 #include "DataFormats/TauReco/interface/CaloTauDiscriminatorByIsolation.h"
00008
00009 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
00010 #include "DataFormats/EgammaCandidates/interface/GsfElectronFwd.h"
00011 #include "DataFormats/EgammaCandidates/interface/Electron.h"
00012 #include "DataFormats/EgammaCandidates/interface/ElectronFwd.h"
00013 #include "AnalysisDataFormats/Egamma/interface/ElectronIDAssociation.h"
00014 #include "AnalysisDataFormats/Egamma/interface/ElectronID.h"
00015 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
00016 #include "DataFormats/TrackReco/interface/Track.h"
00017
00018 #include "DataFormats/MuonReco/interface/Muon.h"
00019 #include "DataFormats/JetReco/interface/CaloJet.h"
00020 #include "TLorentzVector.h"
00021 #include "DataFormats/Math/interface/deltaR.h"
00022
00023 #include "DataFormats/CaloTowers/interface/CaloTower.h"
00024 #include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
00025 #include "DataFormats/CaloTowers/interface/CaloTowerFwd.h"
00026 #include "Math/GenVector/VectorUtil.h"
00027
00028
00029 using namespace edm;
00030 using namespace reco;
00031 using namespace std;
00032
00033 HLTTauRefProducer::HLTTauRefProducer(const edm::ParameterSet& iConfig)
00034 {
00035
00036
00037
00038
00039 ParameterSet pfTau = iConfig.getUntrackedParameter<edm::ParameterSet>("PFTaus");
00040 PFTaus_ = pfTau.getUntrackedParameter<InputTag>("PFTauProducer");
00041 PFTauDis_ = pfTau.getUntrackedParameter<std::vector<InputTag> >("PFTauDiscriminators");
00042 doPFTaus_ = pfTau.getUntrackedParameter<bool>("doPFTaus",false);
00043 ptMinPFTau_= pfTau.getUntrackedParameter<double>("ptMin",15.);
00044
00045 ParameterSet electrons = iConfig.getUntrackedParameter<edm::ParameterSet>("Electrons");
00046 Electrons_ = electrons.getUntrackedParameter<InputTag>("ElectronCollection");
00047 doElectrons_ = electrons.getUntrackedParameter<bool>("doElectrons",false);
00048 e_idAssocProd_ = electrons.getUntrackedParameter<InputTag>("IdCollection");
00049 e_ctfTrackCollection_= electrons.getUntrackedParameter<InputTag>("TrackCollection");
00050 ptMinElectron_= electrons.getUntrackedParameter<double>("ptMin",15.);
00051 e_doID_ = electrons.getUntrackedParameter<bool>("doID",false);
00052 e_doTrackIso_ = electrons.getUntrackedParameter<bool>("doTrackIso",false);
00053 e_trackMinPt_= electrons.getUntrackedParameter<double>("ptMinTrack",1.5);
00054 e_lipCut_= electrons.getUntrackedParameter<double>("lipMinTrack",1.5);
00055 e_minIsoDR_= electrons.getUntrackedParameter<double>("InnerConeDR",0.02);
00056 e_maxIsoDR_= electrons.getUntrackedParameter<double>("OuterConeDR",0.6);
00057 e_isoMaxSumPt_= electrons.getUntrackedParameter<double>("MaxIsoVar",0.02);
00058
00059 ParameterSet muons = iConfig.getUntrackedParameter<edm::ParameterSet>("Muons");
00060 Muons_ = muons.getUntrackedParameter<InputTag>("MuonCollection");
00061 doMuons_ = muons.getUntrackedParameter<bool>("doMuons",false);
00062 ptMinMuon_= muons.getUntrackedParameter<double>("ptMin",15.);
00063
00064 ParameterSet jets = iConfig.getUntrackedParameter<edm::ParameterSet>("Jets");
00065 Jets_ = jets.getUntrackedParameter<InputTag>("JetCollection");
00066 doJets_ = jets.getUntrackedParameter<bool>("doJets");
00067 ptMinJet_= jets.getUntrackedParameter<double>("etMin");
00068
00069 ParameterSet towers = iConfig.getUntrackedParameter<edm::ParameterSet>("Towers");
00070 Towers_ = towers.getUntrackedParameter<InputTag>("TowerCollection");
00071 doTowers_ = towers.getUntrackedParameter<bool>("doTowers");
00072 ptMinTower_= towers.getUntrackedParameter<double>("etMin");
00073 towerIsol_= towers.getUntrackedParameter<double>("towerIsolation");
00074
00075 ParameterSet photons = iConfig.getUntrackedParameter<edm::ParameterSet>("Photons");
00076 Photons_ = photons.getUntrackedParameter<InputTag>("PhotonCollection");
00077 doPhotons_ = photons.getUntrackedParameter<bool>("doPhotons");
00078 ptMinPhoton_= photons.getUntrackedParameter<double>("etMin");
00079 photonEcalIso_= photons.getUntrackedParameter<double>("ECALIso");
00080
00081
00082 etaMax = iConfig.getUntrackedParameter<double>("EtaMax",2.5);
00083
00084
00085
00086 produces<LorentzVectorCollection>("PFTaus");
00087 produces<LorentzVectorCollection>("Electrons");
00088 produces<LorentzVectorCollection>("Muons");
00089 produces<LorentzVectorCollection>("Jets");
00090 produces<LorentzVectorCollection>("Photons");
00091 produces<LorentzVectorCollection>("Towers");
00092
00093 }
00094
00095 HLTTauRefProducer::~HLTTauRefProducer(){ }
00096
00097 void HLTTauRefProducer::produce(edm::Event& iEvent, const edm::EventSetup& iES)
00098 {
00099 if(doPFTaus_)
00100 doPFTaus(iEvent,iES);
00101 if(doElectrons_)
00102 doElectrons(iEvent,iES);
00103 if(doMuons_)
00104 doMuons(iEvent,iES);
00105 if(doJets_)
00106 doJets(iEvent,iES);
00107 if(doPhotons_)
00108 doPhotons(iEvent,iES);
00109 if(doTowers_)
00110 doTowers(iEvent,iES);
00111
00112 }
00113
00114 void
00115 HLTTauRefProducer::doPFTaus(edm::Event& iEvent,const edm::EventSetup& iES)
00116 {
00117 auto_ptr<LorentzVectorCollection> product_PFTaus(new LorentzVectorCollection);
00118
00119 edm::Handle<PFTauCollection> pftaus;
00120 if(iEvent.getByLabel(PFTaus_,pftaus))
00121 {
00122 for(unsigned int i=0;i<pftaus->size();++i)
00123 if((*pftaus)[i].pt()>ptMinPFTau_&&fabs((*pftaus)[i].eta())<etaMax)
00124 {
00125 reco::PFTauRef thePFTau(pftaus,i);
00126 for(unsigned int j=0;j<PFTauDis_.size();++j)
00127 {
00128 edm::Handle<PFTauDiscriminator> pftaudis;
00129 if(iEvent.getByLabel(PFTauDis_[j],pftaudis))
00130 {
00131 if((*pftaudis)[thePFTau]>0.5)
00132 {
00133 LorentzVector vec((*pftaus)[i].px(),(*pftaus)[i].py(),(*pftaus)[i].pz(),(*pftaus)[i].energy());
00134 product_PFTaus->push_back(vec);
00135 }
00136 }
00137 }
00138 }
00139 }
00140
00141 iEvent.put(product_PFTaus,"PFTaus");
00142
00143 }
00144
00145
00146 void
00147 HLTTauRefProducer::doElectrons(edm::Event& iEvent,const edm::EventSetup& iES)
00148 {
00149 auto_ptr<LorentzVectorCollection> product_Electrons(new LorentzVectorCollection);
00150
00151
00152 edm::Handle<reco::ElectronIDAssociationCollection> pEleID;
00153 if(e_doID_){
00154
00155 iEvent.getByLabel(e_idAssocProd_,pEleID);
00156
00157 if (!pEleID.isValid()){
00158 edm::LogInfo("")<< "Error! Can't get electronIDAssocProducer by label. ";
00159 e_doID_ = false;
00160 }
00161 }
00162 edm::Handle<reco::TrackCollection> pCtfTracks;
00163 iEvent.getByLabel(e_ctfTrackCollection_, pCtfTracks);
00164 if (!pCtfTracks.isValid()) {
00165 edm::LogInfo("")<< "Error! Can't get " << e_ctfTrackCollection_.label() << " by label. ";
00166 iEvent.put(product_Electrons,"Electrons");
00167 return;
00168 }
00169 const reco::TrackCollection * ctfTracks = pCtfTracks.product();
00170 edm::Handle<GsfElectronCollection> electrons;
00171 if(iEvent.getByLabel(Electrons_,electrons))
00172 for(size_t i=0;i<electrons->size();++i)
00173 {
00174 edm::Ref<reco::GsfElectronCollection> electronRef(electrons,i);
00175 bool idDec=false;
00176 if(e_doID_){
00177 reco::ElectronIDAssociationCollection::const_iterator tagIDAssocItr;
00178 tagIDAssocItr = pEleID->find(electronRef);
00179 const reco::ElectronIDRef& id_tag = tagIDAssocItr->val;
00180 idDec=id_tag->cutBasedDecision();
00181 }else idDec=true;
00182 if((*electrons)[i].pt()>ptMinElectron_&&fabs((*electrons)[i].eta())<etaMax&&idDec)
00183 {
00184 if(e_doTrackIso_){
00185 reco::TrackCollection::const_iterator tr = ctfTracks->begin();
00186 double sum_of_pt_ele=0;
00187 for(;tr != ctfTracks->end();++tr)
00188 {
00189 double lip = (*electrons)[i].gsfTrack()->dz() - tr->dz();
00190 if(tr->pt() > e_trackMinPt_ && fabs(lip) < e_lipCut_){
00191 double dphi=fabs(tr->phi()-(*electrons)[i].trackMomentumAtVtx().phi());
00192 if(dphi>acos(-1.))dphi=2*acos(-1.)-dphi;
00193 double deta=fabs(tr->eta()-(*electrons)[i].trackMomentumAtVtx().eta());
00194 double dr_ctf_ele = sqrt(deta*deta+dphi*dphi);
00195 if((dr_ctf_ele>e_minIsoDR_) && (dr_ctf_ele<e_maxIsoDR_)){
00196 double cft_pt_2 = (tr->pt())*(tr->pt());
00197 sum_of_pt_ele += cft_pt_2;
00198 }
00199 }
00200 }
00201 double isolation_value_ele = sum_of_pt_ele/((*electrons)[i].trackMomentumAtVtx().Rho()*(*electrons)[i].trackMomentumAtVtx().Rho());
00202 if(isolation_value_ele<e_isoMaxSumPt_){
00203 LorentzVector vec((*electrons)[i].px(),(*electrons)[i].py(),(*electrons)[i].pz(),(*electrons)[i].energy());
00204 product_Electrons->push_back(vec);
00205 }
00206
00207 }
00208 else{
00209 LorentzVector vec((*electrons)[i].px(),(*electrons)[i].py(),(*electrons)[i].pz(),(*electrons)[i].energy());
00210 product_Electrons->push_back(vec);
00211 }
00212 }
00213 }
00214
00215 iEvent.put(product_Electrons,"Electrons");
00216 }
00217
00218 void
00219 HLTTauRefProducer::doMuons(edm::Event& iEvent,const edm::EventSetup& iES)
00220 {
00221 auto_ptr<LorentzVectorCollection> product_Muons(new LorentzVectorCollection);
00222
00223 edm::Handle<MuonCollection> muons;
00224 if(iEvent.getByLabel(Muons_,muons))
00225
00226 for(size_t i = 0 ;i<muons->size();++i)
00227 {
00228
00229 if((*muons)[i].pt()>ptMinMuon_&&fabs((*muons)[i].eta())<etaMax)
00230 {
00231 LorentzVector vec((*muons)[i].px(),(*muons)[i].py(),(*muons)[i].pz(),(*muons)[i].energy());
00232 product_Muons->push_back(vec);
00233 }
00234 }
00235
00236
00237 iEvent.put(product_Muons,"Muons");
00238
00239 }
00240
00241
00242 void
00243 HLTTauRefProducer::doJets(edm::Event& iEvent,const edm::EventSetup& iES)
00244 {
00245 auto_ptr<LorentzVectorCollection> product_Jets(new LorentzVectorCollection);
00246
00247 edm::Handle<CaloJetCollection> jets;
00248 if(iEvent.getByLabel(Jets_,jets))
00249 for(size_t i = 0 ;i<jets->size();++i)
00250 {
00251 if((*jets)[i].et()>ptMinJet_&&fabs((*jets)[i].eta())<etaMax)
00252 {
00253 LorentzVector vec((*jets)[i].px(),(*jets)[i].py(),(*jets)[i].pz(),(*jets)[i].energy());
00254 product_Jets->push_back(vec);
00255 }
00256 }
00257 iEvent.put(product_Jets,"Jets");
00258 }
00259
00260 void
00261 HLTTauRefProducer::doTowers(edm::Event& iEvent,const edm::EventSetup& iES)
00262 {
00263 auto_ptr<LorentzVectorCollection> product_Towers(new LorentzVectorCollection);
00264
00265 edm::Handle<CaloTowerCollection> towers;
00266 if(iEvent.getByLabel(Towers_,towers))
00267 for(size_t i = 0 ;i<towers->size();++i)
00268 {
00269 if((*towers)[i].pt()>ptMinTower_&&fabs((*towers)[i].eta())<etaMax)
00270 {
00271
00272 double isolET=0;
00273 for(unsigned int j=0;j<towers->size();++j)
00274 {
00275 if(ROOT::Math::VectorUtil::DeltaR((*towers)[i].p4(),(*towers)[j].p4())<0.5)
00276 isolET+=(*towers)[j].pt();
00277 }
00278 isolET-=(*towers)[i].pt();
00279 if(isolET<towerIsol_)
00280 {
00281 LorentzVector vec((*towers)[i].px(),(*towers)[i].py(),(*towers)[i].pz(),(*towers)[i].energy());
00282 product_Towers->push_back(vec);
00283 }
00284 }
00285 }
00286 iEvent.put(product_Towers,"Towers");
00287 }
00288
00289
00290 void
00291 HLTTauRefProducer::doPhotons(edm::Event& iEvent,const edm::EventSetup& iES)
00292 {
00293 auto_ptr<LorentzVectorCollection> product_Gammas(new LorentzVectorCollection);
00294
00295 edm::Handle<PhotonCollection> photons;
00296 if(iEvent.getByLabel(Photons_,photons))
00297 for(size_t i = 0 ;i<photons->size();++i)
00298 if((*photons)[i].ecalRecHitSumEtConeDR04()<photonEcalIso_)
00299 {
00300 if((*photons)[i].et()>ptMinPhoton_&&fabs((*photons)[i].eta())<etaMax)
00301 {
00302 LorentzVector vec((*photons)[i].px(),(*photons)[i].py(),(*photons)[i].pz(),(*photons)[i].energy());
00303 product_Gammas->push_back(vec);
00304 }
00305 }
00306 iEvent.put(product_Gammas,"Photons");
00307 }
00308