![]() |
![]() |
00001 /* 00002 * =========================================================================== 00003 * 00004 * Filename: RecoTauTwoProngFilter 00005 * 00006 * Description: Modify taus remove low-pt second prongs 00007 * 00008 * Author: Evan K. Friis (UC Davis) 00009 * 00010 * =========================================================================== 00011 */ 00012 00013 #include "RecoTauTag/RecoTau/interface/RecoTauBuilderPlugins.h" 00014 #include "RecoTauTag/RecoTau/interface/RecoTauCommonUtilities.h" 00015 00016 namespace reco { namespace tau { 00017 00018 namespace { 00019 // Delete an element from a ref vector 00020 void deleteFrom(const PFCandidateRef ref, PFCandidateRefVector* collection) { 00021 PFCandidateRefVector::iterator todelete = collection->end(); 00022 for (PFCandidateRefVector::iterator cand = collection->begin(); 00023 cand != collection->end(); ++cand) { 00024 if (*cand == ref) { 00025 todelete = cand; 00026 break; 00027 } 00028 } 00029 if (todelete != collection->end()) 00030 collection->erase(todelete); 00031 } 00032 } 00033 00034 class RecoTauTwoProngFilter : public RecoTauModifierPlugin { 00035 public: 00036 explicit RecoTauTwoProngFilter(const edm::ParameterSet& pset); 00037 virtual ~RecoTauTwoProngFilter() {} 00038 void operator()(PFTau&) const; 00039 private: 00040 double minPtFractionForSecondProng_; 00041 }; 00042 00043 RecoTauTwoProngFilter::RecoTauTwoProngFilter(const edm::ParameterSet& pset):RecoTauModifierPlugin(pset) { 00044 minPtFractionForSecondProng_ = pset.getParameter<double>("minPtFractionForSecondProng"); 00045 } 00046 00047 void RecoTauTwoProngFilter::operator()(PFTau& tau) const { 00048 if (tau.signalPFChargedHadrCands().size() == 2) { 00049 const PFCandidateRefVector &signalCharged = tau.signalPFChargedHadrCands(); 00050 size_t indexOfHighestPt = 00051 (signalCharged[0]->pt() > signalCharged[1]->pt()) ? 0 : 1; 00052 size_t indexOfLowerPt = ( indexOfHighestPt ) ? 0 : 1; 00053 double ratio = signalCharged[indexOfLowerPt]->pt()/ 00054 signalCharged[indexOfHighestPt]->pt(); 00055 00056 if (ratio < minPtFractionForSecondProng_) { 00057 PFCandidateRef keep = signalCharged[indexOfHighestPt]; 00058 PFCandidateRef filter = signalCharged[indexOfLowerPt]; 00059 // Make our new signal charged candidate collection 00060 PFCandidateRefVector newSignalCharged; 00061 newSignalCharged.push_back(keep); 00062 PFCandidateRefVector newSignal = tau.signalPFCands(); 00063 deleteFrom(filter, &newSignal); 00064 00065 // Copy our filtered cand to isolation 00066 PFCandidateRefVector newIsolationCharged = 00067 tau.isolationPFChargedHadrCands(); 00068 newIsolationCharged.push_back(filter); 00069 PFCandidateRefVector newIsolation = tau.isolationPFCands(); 00070 newIsolation.push_back(filter); 00071 00072 // Update tau members 00073 tau.setP4(tau.p4() - filter->p4()); 00074 tau.setisolationPFChargedHadrCandsPtSum( 00075 tau.isolationPFChargedHadrCandsPtSum() - filter->pt()); 00076 tau.setCharge(tau.charge() - filter->charge()); 00077 // Update tau constituents 00078 tau.setsignalPFChargedHadrCands(newSignalCharged); 00079 tau.setsignalPFCands(newSignal); 00080 tau.setisolationPFChargedHadrCands(newIsolationCharged); 00081 tau.setisolationPFCands(newIsolation); 00082 } 00083 } 00084 } 00085 }} // end namespace reco::tau 00086 #include "FWCore/Framework/interface/MakerMacros.h" 00087 DEFINE_EDM_PLUGIN(RecoTauModifierPluginFactory, 00088 reco::tau::RecoTauTwoProngFilter, 00089 "RecoTauTwoProngFilter");