File indexing completed on 2023-03-17 11:09:18
0001
0002
0003
0004
0005
0006
0007
0008 #include "HLTElectronEtFilter.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 #include "DataFormats/EgammaCandidates/interface/Electron.h"
0018 #include "DataFormats/EgammaCandidates/interface/ElectronIsolationAssociation.h"
0019 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
0020 #include "DataFormats/EgammaReco/interface/SuperClusterFwd.h"
0021
0022 #include "DataFormats/Common/interface/AssociationMap.h"
0023 #include "DataFormats/RecoCandidate/interface/RecoCandidate.h"
0024
0025
0026
0027
0028 HLTElectronEtFilter::HLTElectronEtFilter(const edm::ParameterSet& iConfig) : HLTFilter(iConfig) {
0029 candTag_ = iConfig.getParameter<edm::InputTag>("candTag");
0030 EtEB_ = iConfig.getParameter<double>("EtCutEB");
0031 EtEE_ = iConfig.getParameter<double>("EtCutEE");
0032
0033 ncandcut_ = iConfig.getParameter<int>("ncandcut");
0034
0035 l1EGTag_ = iConfig.getParameter<edm::InputTag>("l1EGCand");
0036
0037 candToken_ = consumes<trigger::TriggerFilterObjectWithRefs>(candTag_);
0038 }
0039
0040 HLTElectronEtFilter::~HLTElectronEtFilter() = default;
0041
0042 void HLTElectronEtFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0043 edm::ParameterSetDescription desc;
0044 makeHLTFilterDescription(desc);
0045 desc.add<edm::InputTag>("candTag", edm::InputTag("hltElectronPixelMatchFilter"));
0046 desc.add<double>("EtCutEB", 0.0);
0047 desc.add<double>("EtCutEE", 0.0);
0048 desc.add<int>("ncandcut", 1);
0049 desc.add<edm::InputTag>("l1EGCand", edm::InputTag("hltL1IsoRecoEcalCandidate"));
0050 descriptions.add("hltElectronEtFilter", desc);
0051 }
0052
0053
0054 bool HLTElectronEtFilter::hltFilter(edm::Event& iEvent,
0055 const edm::EventSetup& iSetup,
0056 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0057 using namespace trigger;
0058 if (saveTags()) {
0059 filterproduct.addCollectionTag(l1EGTag_);
0060 }
0061
0062
0063 reco::ElectronRef ref;
0064
0065 edm::Handle<trigger::TriggerFilterObjectWithRefs> PrevFilterOutput;
0066 iEvent.getByToken(candToken_, PrevFilterOutput);
0067
0068 std::vector<edm::Ref<reco::ElectronCollection> > elecands;
0069 PrevFilterOutput->getObjects(TriggerElectron, elecands);
0070
0071
0072 int n = 0;
0073
0074 for (auto& elecand : elecands) {
0075 ref = elecand;
0076 float Pt = ref->pt();
0077 float Eta = fabs(ref->eta());
0078
0079 if ((Eta < 1.479 && Pt > EtEB_) || (Eta >= 1.479 && Pt > EtEE_)) {
0080 n++;
0081 filterproduct.addObject(TriggerElectron, ref);
0082 }
0083 }
0084
0085 bool accept(n >= ncandcut_);
0086
0087 return accept;
0088 }
0089
0090
0091 #include "FWCore/Framework/interface/MakerMacros.h"
0092 DEFINE_FWK_MODULE(HLTElectronEtFilter);