File indexing completed on 2023-05-06 02:46:27
0001
0002
0003
0004
0005
0006
0007 #include "HLTEgammaL1TMatchFilterRegional.h"
0008 #include "DataFormats/Common/interface/Handle.h"
0009 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0010 #include "DataFormats/L1Trigger/interface/Jet.h"
0011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0012 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0014 #include "FWCore/Framework/interface/EventSetup.h"
0015
0016 #define TWOPI 6.283185308
0017
0018
0019
0020 HLTEgammaL1TMatchFilterRegional::HLTEgammaL1TMatchFilterRegional(const edm::ParameterSet& iConfig)
0021 : HLTFilter(iConfig) {
0022 candIsolatedTag_ = iConfig.getParameter<edm::InputTag>("candIsolatedTag");
0023 l1EGTag_ = iConfig.getParameter<edm::InputTag>(
0024 "l1IsolatedTag");
0025 candNonIsolatedTag_ = iConfig.getParameter<edm::InputTag>("candNonIsolatedTag");
0026 L1SeedFilterTag_ = iConfig.getParameter<edm::InputTag>("L1SeedFilterTag");
0027 l1JetsTag_ = iConfig.getParameter<edm::InputTag>(
0028 "l1CenJetsTag");
0029 l1TausTag_ =
0030 iConfig.getParameter<edm::InputTag>("l1TausTag");
0031 ncandcut_ = iConfig.getParameter<int>("ncandcut");
0032 doIsolated_ = iConfig.getParameter<bool>("doIsolated");
0033 region_eta_size_ = iConfig.getParameter<double>("region_eta_size");
0034 region_eta_size_ecap_ = iConfig.getParameter<double>("region_eta_size_ecap");
0035 region_phi_size_ = iConfig.getParameter<double>("region_phi_size");
0036 barrel_end_ = iConfig.getParameter<double>("barrel_end");
0037 endcap_end_ = iConfig.getParameter<double>("endcap_end");
0038
0039 candIsolatedToken_ = consumes<reco::RecoEcalCandidateCollection>(candIsolatedTag_);
0040 if (!doIsolated_)
0041 candNonIsolatedToken_ = consumes<reco::RecoEcalCandidateCollection>(candNonIsolatedTag_);
0042 L1SeedFilterToken_ = consumes<trigger::TriggerFilterObjectWithRefs>(L1SeedFilterTag_);
0043 }
0044
0045 HLTEgammaL1TMatchFilterRegional::~HLTEgammaL1TMatchFilterRegional() = default;
0046
0047 void HLTEgammaL1TMatchFilterRegional::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0048 edm::ParameterSetDescription desc;
0049 makeHLTFilterDescription(desc);
0050 desc.add<edm::InputTag>("candIsolatedTag", edm::InputTag("hltRecoIsolatedEcalCandidate"));
0051 desc.add<edm::InputTag>("l1IsolatedTag",
0052 edm::InputTag("hltCaloStage2Digis"));
0053 desc.add<edm::InputTag>("candNonIsolatedTag", edm::InputTag("hltRecoNonIsolatedEcalCandidate"));
0054 desc.add<edm::InputTag>("l1NonIsolatedTag",
0055 edm::InputTag("l1extraParticles", "NonIsolated"));
0056 desc.add<edm::InputTag>("L1SeedFilterTag", edm::InputTag("theL1SeedFilter"));
0057 desc.add<edm::InputTag>("l1CenJetsTag", edm::InputTag("hltCaloStage2Digis"));
0058 desc.add<edm::InputTag>("l1TausTag", edm::InputTag("hltCaloStage2Digis", "Tau"));
0059 desc.add<int>("ncandcut", 1);
0060 desc.add<bool>("doIsolated", true);
0061 desc.add<double>("region_eta_size", 0.522);
0062 desc.add<double>("region_eta_size_ecap", 1.0);
0063 desc.add<double>("region_phi_size", 1.044);
0064 desc.add<double>("barrel_end", 1.4791);
0065 desc.add<double>("endcap_end", 2.65);
0066 descriptions.add("HLTEgammaL1TMatchFilterRegional", desc);
0067 }
0068
0069
0070
0071
0072
0073 bool HLTEgammaL1TMatchFilterRegional::hltFilter(edm::Event& iEvent,
0074 const edm::EventSetup& iSetup,
0075 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0076
0077 using namespace trigger;
0078
0079
0080 if (saveTags()) {
0081 filterproduct.addCollectionTag(l1EGTag_);
0082 filterproduct.addCollectionTag(l1JetsTag_);
0083 filterproduct.addCollectionTag(l1TausTag_);
0084 filterproduct.addCollectionTag(candIsolatedTag_);
0085 if (!doIsolated_ && !candNonIsolatedTag_.label().empty()) {
0086 filterproduct.addCollectionTag(candNonIsolatedTag_);
0087 }
0088 }
0089
0090 edm::Ref<reco::RecoEcalCandidateCollection> ref;
0091
0092
0093 int n(0);
0094
0095
0096 edm::Handle<reco::RecoEcalCandidateCollection> recoIsolecalcands;
0097 iEvent.getByToken(candIsolatedToken_, recoIsolecalcands);
0098
0099 edm::Handle<trigger::TriggerFilterObjectWithRefs> L1SeedOutput;
0100 iEvent.getByToken(L1SeedFilterToken_, L1SeedOutput);
0101
0102 std::vector<l1t::EGammaRef> l1EGs;
0103 L1SeedOutput->getObjects(TriggerL1EG, l1EGs);
0104
0105 std::vector<l1t::JetRef> l1Jets;
0106 L1SeedOutput->getObjects(TriggerL1Jet, l1Jets);
0107
0108 std::vector<l1t::TauRef> l1Taus;
0109 L1SeedOutput->getObjects(TriggerL1Tau, l1Taus);
0110
0111 for (auto recoecalcand = recoIsolecalcands->begin(); recoecalcand != recoIsolecalcands->end(); recoecalcand++) {
0112 if (std::abs(recoecalcand->eta()) < endcap_end_) {
0113
0114
0115
0116
0117 bool matchedSCEG = matchedToL1Cand(l1EGs, recoecalcand->eta(), recoecalcand->phi());
0118 bool matchedSCJet = matchedToL1Cand(l1Jets, recoecalcand->eta(), recoecalcand->phi());
0119 bool matchedSCTau = matchedToL1Cand(l1Taus, recoecalcand->eta(), recoecalcand->phi());
0120
0121 if (matchedSCEG || matchedSCJet || matchedSCTau) {
0122 n++;
0123 ref = edm::Ref<reco::RecoEcalCandidateCollection>(recoIsolecalcands,
0124 distance(recoIsolecalcands->begin(), recoecalcand));
0125 filterproduct.addObject(TriggerCluster, ref);
0126 }
0127
0128 }
0129
0130 }
0131
0132
0133
0134
0135 if (!doIsolated_ && !candNonIsolatedTag_.label().empty()) {
0136 edm::Handle<reco::RecoEcalCandidateCollection> recoNonIsolecalcands;
0137 iEvent.getByToken(candNonIsolatedToken_, recoNonIsolecalcands);
0138
0139 for (auto recoecalcand = recoNonIsolecalcands->begin(); recoecalcand != recoNonIsolecalcands->end();
0140 recoecalcand++) {
0141 if (std::abs(recoecalcand->eta()) < endcap_end_) {
0142 bool matchedSCEG = matchedToL1Cand(l1EGs, recoecalcand->eta(), recoecalcand->phi());
0143 bool matchedSCJet = matchedToL1Cand(l1Jets, recoecalcand->eta(), recoecalcand->phi());
0144 bool matchedSCTau = matchedToL1Cand(l1Taus, recoecalcand->eta(), recoecalcand->phi());
0145
0146 if (matchedSCEG || matchedSCJet || matchedSCTau) {
0147 n++;
0148 ref = edm::Ref<reco::RecoEcalCandidateCollection>(recoNonIsolecalcands,
0149 distance(recoNonIsolecalcands->begin(), recoecalcand));
0150 filterproduct.addObject(TriggerCluster, ref);
0151 }
0152
0153 }
0154
0155 }
0156 }
0157
0158
0159 bool accept(n >= ncandcut_);
0160
0161 return accept;
0162 }
0163
0164 bool HLTEgammaL1TMatchFilterRegional::matchedToL1Cand(const std::vector<l1t::EGammaRef>& l1Cands,
0165 const float scEta,
0166 const float scPhi) const {
0167 for (auto const& l1Cand : l1Cands) {
0168
0169 double etaBinLow = 0.;
0170 double etaBinHigh = 0.;
0171 if (std::abs(scEta) < barrel_end_) {
0172 etaBinLow = l1Cand->eta() - region_eta_size_ / 2.;
0173 etaBinHigh = etaBinLow + region_eta_size_;
0174 } else {
0175 etaBinLow = l1Cand->eta() - region_eta_size_ecap_ / 2.;
0176 etaBinHigh = etaBinLow + region_eta_size_ecap_;
0177 }
0178
0179 float deltaphi = std::abs(scPhi - l1Cand->phi());
0180 if (deltaphi > TWOPI)
0181 deltaphi -= TWOPI;
0182 if (deltaphi > TWOPI / 2.)
0183 deltaphi = TWOPI - deltaphi;
0184
0185 if (scEta < etaBinHigh && scEta > etaBinLow && deltaphi < region_phi_size_ / 2.) {
0186 return true;
0187 }
0188 }
0189 return false;
0190 }
0191
0192 bool
0193
0194 HLTEgammaL1TMatchFilterRegional::matchedToL1Cand(const std::vector<l1t::JetRef>& l1Cands,
0195 const float scEta,
0196 const float scPhi) const {
0197 for (auto const& l1Cand : l1Cands) {
0198
0199 double etaBinLow = 0.;
0200 double etaBinHigh = 0.;
0201 if (std::abs(scEta) < barrel_end_) {
0202 etaBinLow = l1Cand->eta() - region_eta_size_ / 2.;
0203 etaBinHigh = etaBinLow + region_eta_size_;
0204 } else {
0205 etaBinLow = l1Cand->eta() - region_eta_size_ecap_ / 2.;
0206 etaBinHigh = etaBinLow + region_eta_size_ecap_;
0207 }
0208
0209 float deltaphi = std::abs(scPhi - l1Cand->phi());
0210 if (deltaphi > TWOPI)
0211 deltaphi -= TWOPI;
0212 if (deltaphi > TWOPI / 2.)
0213 deltaphi = TWOPI - deltaphi;
0214
0215 if (scEta < etaBinHigh && scEta > etaBinLow && deltaphi < region_phi_size_ / 2.) {
0216 return true;
0217 }
0218 }
0219 return false;
0220 }
0221 bool
0222
0223 HLTEgammaL1TMatchFilterRegional::matchedToL1Cand(const std::vector<l1t::TauRef>& l1Cands,
0224 const float scEta,
0225 const float scPhi) const {
0226 for (auto const& l1Cand : l1Cands) {
0227
0228 double etaBinLow = 0.;
0229 double etaBinHigh = 0.;
0230 if (std::abs(scEta) < barrel_end_) {
0231 etaBinLow = l1Cand->eta() - region_eta_size_ / 2.;
0232 etaBinHigh = etaBinLow + region_eta_size_;
0233 } else {
0234 etaBinLow = l1Cand->eta() - region_eta_size_ecap_ / 2.;
0235 etaBinHigh = etaBinLow + region_eta_size_ecap_;
0236 }
0237
0238 float deltaphi = std::abs(scPhi - l1Cand->phi());
0239 if (deltaphi > TWOPI)
0240 deltaphi -= TWOPI;
0241 if (deltaphi > TWOPI / 2.)
0242 deltaphi = TWOPI - deltaphi;
0243
0244 if (scEta < etaBinHigh && scEta > etaBinLow && deltaphi < region_phi_size_ / 2.) {
0245 return true;
0246 }
0247 }
0248 return false;
0249 }
0250
0251
0252 #include "FWCore/Framework/interface/MakerMacros.h"
0253 DEFINE_FWK_MODULE(HLTEgammaL1TMatchFilterRegional);