File indexing completed on 2024-04-06 12:18:20
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "HLTEgammaDoubleEtDeltaPhiFilter.h"
0010
0011 #include "DataFormats/Common/interface/Handle.h"
0012
0013 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0016
0017 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
0018
0019
0020
0021
0022 HLTEgammaDoubleEtDeltaPhiFilter::HLTEgammaDoubleEtDeltaPhiFilter(const edm::ParameterSet& iConfig)
0023 : HLTFilter(iConfig) {
0024 inputTag_ = iConfig.getParameter<edm::InputTag>("inputTag");
0025 etcut_ = iConfig.getParameter<double>("etcut");
0026 minDeltaPhi_ = iConfig.getParameter<double>("minDeltaPhi");
0027 l1EGTag_ = iConfig.getParameter<edm::InputTag>("l1EGCand");
0028 inputToken_ = consumes<trigger::TriggerFilterObjectWithRefs>(inputTag_);
0029 }
0030
0031 HLTEgammaDoubleEtDeltaPhiFilter::~HLTEgammaDoubleEtDeltaPhiFilter() = default;
0032
0033 void HLTEgammaDoubleEtDeltaPhiFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0034 edm::ParameterSetDescription desc;
0035 makeHLTFilterDescription(desc);
0036 desc.add<edm::InputTag>("inputTag", edm::InputTag("hltDoublePhotonEt5L1MatchFilterRegional"));
0037 desc.add<edm::InputTag>("l1EGCand", edm::InputTag("hltL1IsoRecoEcalCandidate"));
0038 desc.add<double>("etcut", 5.0);
0039 desc.add<double>("minDeltaPhi", 2.5);
0040 descriptions.add("hltEgammaDoubleEtDeltaPhiFilter", desc);
0041 }
0042
0043
0044 bool HLTEgammaDoubleEtDeltaPhiFilter::hltFilter(edm::Event& iEvent,
0045 const edm::EventSetup& iSetup,
0046 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0047 using namespace trigger;
0048
0049 if (saveTags()) {
0050 filterproduct.addCollectionTag(l1EGTag_);
0051 }
0052
0053
0054 edm::Handle<trigger::TriggerFilterObjectWithRefs> PrevFilterOutput;
0055 iEvent.getByToken(inputToken_, PrevFilterOutput);
0056
0057 std::vector<edm::Ref<reco::RecoEcalCandidateCollection> > recoecalcands;
0058 PrevFilterOutput->getObjects(TriggerCluster, recoecalcands);
0059 if (recoecalcands.empty())
0060 PrevFilterOutput->getObjects(TriggerPhoton,
0061 recoecalcands);
0062
0063
0064 edm::Ref<reco::RecoEcalCandidateCollection> ref1, ref2;
0065
0066
0067 int n(0);
0068 for (auto& ref : recoecalcands) {
0069 if (ref->et() >= etcut_) {
0070 ++n;
0071 if (n == 1)
0072 ref1 = ref;
0073 if (n == 2)
0074 ref2 = ref;
0075 }
0076 }
0077
0078
0079 double deltaPhi(0.0);
0080 if (n == 2) {
0081 deltaPhi = fabs(ref1->phi() - ref2->phi());
0082 if (deltaPhi > M_PI)
0083 deltaPhi = 2 * M_PI - deltaPhi;
0084
0085 filterproduct.addObject(TriggerCluster, ref1);
0086 filterproduct.addObject(TriggerCluster, ref2);
0087 }
0088
0089
0090 bool accept(n == 2 && deltaPhi > minDeltaPhi_);
0091
0092 return accept;
0093 }
0094
0095
0096
0097
0098
0099
0100 #include "FWCore/Framework/interface/MakerMacros.h"
0101 DEFINE_FWK_MODULE(HLTEgammaDoubleEtDeltaPhiFilter);