File indexing completed on 2024-04-06 12:18:21
0001
0002
0003
0004
0005
0006
0007 #include "HLTElectronMissingHitsFilter.h"
0008 #include "DataFormats/Common/interface/Handle.h"
0009 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0012 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
0013 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidateFwd.h"
0014
0015 #include "DataFormats/EgammaCandidates/interface/Electron.h"
0016 #include "DataFormats/TrackReco/interface/Track.h"
0017 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0018 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
0019
0020 HLTElectronMissingHitsFilter::HLTElectronMissingHitsFilter(const edm::ParameterSet& iConfig) : HLTFilter(iConfig) {
0021 candTag_ = iConfig.getParameter<edm::InputTag>("candTag");
0022 candToken_ = consumes<trigger::TriggerFilterObjectWithRefs>(candTag_);
0023 electronTag_ = iConfig.getParameter<edm::InputTag>("electronProducer");
0024 electronToken_ = consumes<reco::ElectronCollection>(electronTag_);
0025 barrelcut_ = iConfig.getParameter<int>("barrelcut");
0026 endcapcut_ = iConfig.getParameter<int>("endcapcut");
0027 ncandcut_ = iConfig.getParameter<int>("ncandcut");
0028 }
0029
0030 HLTElectronMissingHitsFilter::~HLTElectronMissingHitsFilter() = default;
0031
0032 void HLTElectronMissingHitsFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0033 edm::ParameterSetDescription desc;
0034 makeHLTFilterDescription(desc);
0035 desc.add<edm::InputTag>("candTag", edm::InputTag("hltSingleElectronOneOEMinusOneOPFilter"));
0036 desc.add<edm::InputTag>("electronProducer",
0037 edm::InputTag("hltL1NonIsoHLTNonIsoSingleElectronEt15LTIPixelMatchFilte"));
0038 desc.add<int>("barrelcut", 0);
0039 desc.add<int>("endcapcut", 0);
0040 desc.add<int>("ncandcut", 1);
0041 descriptions.add("hltElectronMissingHitsFilter", desc);
0042 }
0043
0044 bool HLTElectronMissingHitsFilter::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(electronTag_);
0051
0052 edm::Handle<trigger::TriggerFilterObjectWithRefs> PrevFilterOutput;
0053 iEvent.getByToken(candToken_, PrevFilterOutput);
0054
0055 std::vector<edm::Ref<reco::RecoEcalCandidateCollection> > recoecalcands;
0056 PrevFilterOutput->getObjects(TriggerCluster, recoecalcands);
0057 if (recoecalcands.empty())
0058 PrevFilterOutput->getObjects(TriggerPhoton, recoecalcands);
0059
0060 edm::Handle<reco::ElectronCollection> electronHandle;
0061 iEvent.getByToken(electronToken_, electronHandle);
0062
0063 int n(0);
0064
0065 edm::RefToBase<reco::Candidate> candref;
0066 for (auto& recoecalcand : recoecalcands) {
0067 reco::SuperClusterRef scCand = recoecalcand->superCluster();
0068 for (auto iElectron = electronHandle->begin(); iElectron != electronHandle->end(); iElectron++) {
0069 reco::ElectronRef electronref(reco::ElectronRef(electronHandle, iElectron - electronHandle->begin()));
0070 const reco::SuperClusterRef scEle = electronref->superCluster();
0071 if (scCand == scEle) {
0072 int missinghits = 0;
0073 if (electronref->gsfTrack().isNonnull()) {
0074 missinghits = electronref->gsfTrack()->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS);
0075 } else if (electronref->track().isNonnull()) {
0076 missinghits = electronref->track()->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS);
0077 } else {
0078 std::cerr << "Electron without track..." << std::endl;
0079 }
0080
0081 if (fabs(electronref->eta()) < 1.479) {
0082 if (missinghits < barrelcut_) {
0083 n++;
0084 filterproduct.addObject(TriggerElectron, electronref);
0085 }
0086 }
0087
0088 if (fabs(electronref->eta()) > 1.479) {
0089 if (missinghits < endcapcut_) {
0090 n++;
0091 filterproduct.addObject(TriggerElectron, electronref);
0092 }
0093 }
0094 }
0095 }
0096 }
0097
0098 bool accept(n >= ncandcut_);
0099
0100 return accept;
0101 }
0102
0103
0104 #include "FWCore/Framework/interface/MakerMacros.h"
0105 DEFINE_FWK_MODULE(HLTElectronMissingHitsFilter);