File indexing completed on 2023-03-17 11:09:17
0001
0002
0003
0004
0005
0006
0007
0008 #include "HLTEgammaEtFilterPairs.h"
0009
0010 #include "DataFormats/Common/interface/Handle.h"
0011
0012 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015
0016 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
0017
0018
0019
0020
0021 HLTEgammaEtFilterPairs::HLTEgammaEtFilterPairs(const edm::ParameterSet& iConfig) : HLTFilter(iConfig) {
0022 inputTag_ = iConfig.getParameter<edm::InputTag>("inputTag");
0023 etcutEB1_ = iConfig.getParameter<double>("etcut1EB");
0024 etcutEE1_ = iConfig.getParameter<double>("etcut1EE");
0025 etcutEB2_ = iConfig.getParameter<double>("etcut2EB");
0026 etcutEE2_ = iConfig.getParameter<double>("etcut2EE");
0027 l1EGTag_ = iConfig.getParameter<edm::InputTag>("l1EGCand");
0028 inputToken_ = consumes<trigger::TriggerFilterObjectWithRefs>(inputTag_);
0029 }
0030
0031 void HLTEgammaEtFilterPairs::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0032 edm::ParameterSetDescription desc;
0033 makeHLTFilterDescription(desc);
0034 desc.add<edm::InputTag>("inputTag", edm::InputTag("HLTEgammaL1MatchFilter"));
0035 desc.add<edm::InputTag>("l1EGCand", edm::InputTag("hltL1IsoRecoEcalCandidate"));
0036 desc.add<double>("etcut1EB", 1.0);
0037 desc.add<double>("etcut1EE", 1.0);
0038 desc.add<double>("etcut2EB", 1.0);
0039 desc.add<double>("etcut2EE", 1.0);
0040 descriptions.add("hltEgammaEtFilterPairs", desc);
0041 }
0042
0043 HLTEgammaEtFilterPairs::~HLTEgammaEtFilterPairs() = default;
0044
0045
0046 bool HLTEgammaEtFilterPairs::hltFilter(edm::Event& iEvent,
0047 const edm::EventSetup& iSetup,
0048 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0049 using namespace trigger;
0050
0051 if (saveTags()) {
0052 filterproduct.addCollectionTag(l1EGTag_);
0053 }
0054
0055 edm::Handle<trigger::TriggerFilterObjectWithRefs> PrevFilterOutput;
0056 iEvent.getByToken(inputToken_, PrevFilterOutput);
0057
0058 std::vector<edm::Ref<reco::RecoEcalCandidateCollection> > recoecalcands;
0059 PrevFilterOutput->getObjects(TriggerCluster, recoecalcands);
0060 if (recoecalcands.empty())
0061 PrevFilterOutput->getObjects(TriggerPhoton, recoecalcands);
0062
0063
0064
0065
0066
0067
0068
0069 int n(0);
0070
0071 for (unsigned int i = 0; i < recoecalcands.size(); i = i + 2) {
0072 edm::Ref<reco::RecoEcalCandidateCollection> r1 = recoecalcands[i];
0073 edm::Ref<reco::RecoEcalCandidateCollection> r2 = recoecalcands[i + 1];
0074
0075 bool first =
0076 (fabs(r1->eta()) < 1.479 && r1->et() >= etcutEB1_) || (fabs(r1->eta()) >= 1.479 && r1->et() >= etcutEE1_);
0077 bool second =
0078 (fabs(r2->eta()) < 1.479 && r2->et() >= etcutEB2_) || (fabs(r2->eta()) >= 1.479 && r2->et() >= etcutEE2_);
0079
0080 if (first && second) {
0081 n++;
0082 filterproduct.addObject(TriggerCluster, r1);
0083 filterproduct.addObject(TriggerCluster, r2);
0084 }
0085 }
0086
0087
0088 bool accept(n >= 1);
0089
0090 return accept;
0091 }
0092
0093
0094 #include "FWCore/Framework/interface/MakerMacros.h"
0095 DEFINE_FWK_MODULE(HLTEgammaEtFilterPairs);