File indexing completed on 2024-04-06 12:18:21
0001
0002
0003
0004
0005
0006
0007
0008 #include "HLTGenericFilter.h"
0009
0010 #include "DataFormats/Common/interface/Handle.h"
0011
0012 #include "FWCore/Framework/interface/MakerMacros.h"
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/EgammaReco/interface/SuperCluster.h"
0018 #include "DataFormats/EgammaReco/interface/SuperClusterFwd.h"
0019 #include "DataFormats/Common/interface/AssociationMap.h"
0020
0021 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0022
0023
0024
0025
0026 template <typename T1>
0027 HLTGenericFilter<T1>::HLTGenericFilter(const edm::ParameterSet& iConfig) : HLTFilter(iConfig) {
0028 candTag_ = iConfig.template getParameter<edm::InputTag>("candTag");
0029 varTag_ = iConfig.template getParameter<edm::InputTag>("varTag");
0030 l1EGTag_ = iConfig.template getParameter<edm::InputTag>("l1EGCand");
0031 rhoTag_ = iConfig.template getParameter<edm::InputTag>("rhoTag");
0032
0033 energyLowEdges_ = iConfig.template getParameter<std::vector<double>>("energyLowEdges");
0034 lessThan_ = iConfig.template getParameter<bool>("lessThan");
0035 useEt_ = iConfig.template getParameter<bool>("useEt");
0036 useAbs_ = iConfig.template getParameter<bool>("useAbs");
0037 thrRegularEB_ = iConfig.template getParameter<std::vector<double>>("thrRegularEB");
0038 thrRegularEE_ = iConfig.template getParameter<std::vector<double>>("thrRegularEE");
0039 thrOverEEB_ = iConfig.template getParameter<std::vector<double>>("thrOverEEB");
0040 thrOverEEE_ = iConfig.template getParameter<std::vector<double>>("thrOverEEE");
0041 thrOverE2EB_ = iConfig.template getParameter<std::vector<double>>("thrOverE2EB");
0042 thrOverE2EE_ = iConfig.template getParameter<std::vector<double>>("thrOverE2EE");
0043 ncandcut_ = iConfig.template getParameter<int>("ncandcut");
0044
0045 doRhoCorrection_ = iConfig.template getParameter<bool>("doRhoCorrection");
0046 rhoMax_ = iConfig.template getParameter<double>("rhoMax");
0047 rhoScale_ = iConfig.template getParameter<double>("rhoScale");
0048 effectiveAreas_ = iConfig.template getParameter<std::vector<double>>("effectiveAreas");
0049 absEtaLowEdges_ = iConfig.template getParameter<std::vector<double>>("absEtaLowEdges");
0050
0051 candToken_ = consumes<trigger::TriggerFilterObjectWithRefs>(candTag_);
0052 varToken_ = consumes<T1IsolationMap>(varTag_);
0053
0054 if (energyLowEdges_.size() != thrRegularEB_.size() or energyLowEdges_.size() != thrRegularEE_.size() or
0055 energyLowEdges_.size() != thrOverEEB_.size() or energyLowEdges_.size() != thrOverEEE_.size() or
0056 energyLowEdges_.size() != thrOverE2EB_.size() or energyLowEdges_.size() != thrOverE2EE_.size())
0057 throw cms::Exception("IncompatibleVects") << "energyLowEdges and threshold vectors should be of the same size. \n";
0058
0059 if (energyLowEdges_.at(0) != 0.0)
0060 throw cms::Exception("IncompleteCoverage") << "energyLowEdges should start from 0. \n";
0061
0062 for (unsigned int aIt = 0; aIt < energyLowEdges_.size() - 1; aIt++) {
0063 if (!(energyLowEdges_.at(aIt) < energyLowEdges_.at(aIt + 1)))
0064 throw cms::Exception("ImproperBinning") << "energyLowEdges entries should be in increasing order. \n";
0065 }
0066
0067 if (doRhoCorrection_) {
0068 rhoToken_ = consumes<double>(rhoTag_);
0069 if (absEtaLowEdges_.size() != effectiveAreas_.size())
0070 throw cms::Exception("IncompatibleVects") << "absEtaLowEdges and effectiveAreas should be of the same size. \n";
0071
0072 if (absEtaLowEdges_.at(0) != 0.0)
0073 throw cms::Exception("IncompleteCoverage") << "absEtaLowEdges should start from 0. \n";
0074
0075 for (unsigned int bIt = 0; bIt < absEtaLowEdges_.size() - 1; bIt++) {
0076 if (!(absEtaLowEdges_.at(bIt) < absEtaLowEdges_.at(bIt + 1)))
0077 throw cms::Exception("ImproperBinning") << "absEtaLowEdges entries should be in increasing order. \n";
0078 }
0079 }
0080 }
0081
0082 template <typename T1>
0083 void HLTGenericFilter<T1>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0084 edm::ParameterSetDescription desc;
0085 makeHLTFilterDescription(desc);
0086 desc.add<edm::InputTag>("candTag", edm::InputTag("hltSingleEgammaEtFilter"));
0087 desc.add<edm::InputTag>("varTag", edm::InputTag("hltSingleEgammaHcalIsol"));
0088 desc.add<edm::InputTag>("rhoTag", edm::InputTag(""));
0089 desc.add<std::vector<double>>("energyLowEdges", {0.0});
0090 desc.add<bool>("lessThan", true);
0091 desc.add<bool>("useEt", false);
0092 desc.add<bool>("useAbs", false);
0093 desc.add<std::vector<double>>("thrRegularEB", {0.0});
0094 desc.add<std::vector<double>>("thrRegularEE", {0.0});
0095 desc.add<std::vector<double>>("thrOverEEB", {-1.0});
0096 desc.add<std::vector<double>>("thrOverEEE", {-1.0});
0097 desc.add<std::vector<double>>("thrOverE2EB", {-1.0});
0098 desc.add<std::vector<double>>("thrOverE2EE", {-1.0});
0099 desc.add<int>("ncandcut", 1);
0100 desc.add<bool>("doRhoCorrection", false);
0101 desc.add<double>("rhoMax", 9.9999999E7);
0102 desc.add<double>("rhoScale", 1.0);
0103 desc.add<std::vector<double>>("effectiveAreas", {0.0, 0.0});
0104 desc.add<std::vector<double>>("absEtaLowEdges", {0.0, 1.479});
0105 desc.add<edm::InputTag>("l1EGCand", edm::InputTag("hltL1IsoRecoEcalCandidate"));
0106 descriptions.add(defaultModuleLabel<HLTGenericFilter<T1>>(), desc);
0107 }
0108
0109 template <typename T1>
0110 HLTGenericFilter<T1>::~HLTGenericFilter() {}
0111
0112 template <typename T1>
0113 float HLTGenericFilter<T1>::getEnergy(T1Ref candRef) const {
0114 return candRef->p();
0115 }
0116
0117 template <>
0118 float HLTGenericFilter<reco::RecoEcalCandidate>::getEnergy(T1Ref candRef) const {
0119 return candRef->superCluster()->energy();
0120 }
0121
0122 template <typename T1>
0123 float HLTGenericFilter<T1>::getEt(T1Ref candRef) const {
0124 return candRef->pt();
0125 }
0126
0127 template <>
0128 float HLTGenericFilter<reco::RecoEcalCandidate>::getEt(T1Ref candRef) const {
0129 return candRef->superCluster()->energy() * sin(2 * atan(exp(-candRef->eta())));
0130 }
0131
0132
0133 template <typename T1>
0134 bool HLTGenericFilter<T1>::hltFilter(edm::Event& iEvent,
0135 const edm::EventSetup& iSetup,
0136 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0137 using namespace trigger;
0138 if (saveTags()) {
0139 filterproduct.addCollectionTag(l1EGTag_);
0140 }
0141
0142
0143 int trigger_type = trigger::TriggerCluster;
0144 if (saveTags())
0145 trigger_type = trigger::TriggerPhoton;
0146
0147 edm::Handle<trigger::TriggerFilterObjectWithRefs> PrevFilterOutput;
0148 iEvent.getByToken(candToken_, PrevFilterOutput);
0149
0150 std::vector<T1Ref> recoCands;
0151 PrevFilterOutput->getObjects(TriggerCluster, recoCands);
0152 if (recoCands.empty())
0153 PrevFilterOutput->getObjects(TriggerPhoton,
0154 recoCands);
0155 if (recoCands.empty()) {
0156 PrevFilterOutput->getObjects(TriggerMuon, recoCands);
0157 trigger_type = trigger::TriggerMuon;
0158 }
0159
0160 edm::Handle<T1IsolationMap> depMap;
0161 iEvent.getByToken(varToken_, depMap);
0162
0163
0164 edm::Handle<double> rhoHandle;
0165 double rho = 0.0;
0166 if (doRhoCorrection_) {
0167 iEvent.getByToken(rhoToken_, rhoHandle);
0168 rho = *(rhoHandle.product());
0169 }
0170
0171 if (rho > rhoMax_)
0172 rho = rhoMax_;
0173 rho = rho * rhoScale_;
0174
0175
0176 int n = 0;
0177 for (unsigned int i = 0; i < recoCands.size(); i++) {
0178
0179 T1Ref ref = recoCands[i];
0180 typename T1IsolationMap::const_iterator mapi = (*depMap).find(ref);
0181
0182
0183
0184
0185 float vali = useAbs_ ? std::abs(mapi->val) : mapi->val;
0186
0187 float EtaSC = ref->eta();
0188
0189
0190 if (doRhoCorrection_) {
0191 auto cIt = std::lower_bound(absEtaLowEdges_.begin(), absEtaLowEdges_.end(), std::abs(EtaSC)) - 1;
0192 vali = vali - (rho * effectiveAreas_.at(std::distance(absEtaLowEdges_.begin(), cIt)));
0193 }
0194
0195 float energy;
0196 if (useEt_)
0197 energy = getEt(ref);
0198 else
0199 energy = getEnergy(ref);
0200
0201
0202 double cutRegularEB_ = 9999., cutRegularEE_ = 9999.;
0203 double cutOverEEB_ = 9999., cutOverEEE_ = 9999.;
0204 double cutOverE2EB_ = 9999., cutOverE2EE_ = 9999.;
0205
0206 auto dIt = std::lower_bound(energyLowEdges_.begin(), energyLowEdges_.end(), energy) - 1;
0207 unsigned iEn = std::distance(energyLowEdges_.begin(), dIt);
0208
0209 cutRegularEB_ = thrRegularEB_.at(iEn);
0210 cutRegularEE_ = thrRegularEE_.at(iEn);
0211 cutOverEEB_ = thrOverEEB_.at(iEn);
0212 cutOverEEE_ = thrOverEEE_.at(iEn);
0213 cutOverE2EB_ = thrOverE2EB_.at(iEn);
0214 cutOverE2EE_ = thrOverE2EE_.at(iEn);
0215
0216 if (lessThan_) {
0217 if ((std::abs(EtaSC) < 1.479 && vali <= cutRegularEB_) || (std::abs(EtaSC) >= 1.479 && vali <= cutRegularEE_)) {
0218 n++;
0219 filterproduct.addObject(trigger_type, ref);
0220 continue;
0221 }
0222 if (energy > 0. && (cutOverEEB_ > 0. || cutOverEEE_ > 0. || cutOverE2EB_ > 0. || cutOverE2EE_ > 0.)) {
0223 if ((std::abs(EtaSC) < 1.479 && vali / energy <= cutOverEEB_) ||
0224 (std::abs(EtaSC) >= 1.479 && vali / energy <= cutOverEEE_)) {
0225 n++;
0226 filterproduct.addObject(trigger_type, ref);
0227 continue;
0228 }
0229 if ((std::abs(EtaSC) < 1.479 && vali / (energy * energy) <= cutOverE2EB_) ||
0230 (std::abs(EtaSC) >= 1.479 && vali / (energy * energy) <= cutOverE2EE_)) {
0231 n++;
0232 filterproduct.addObject(trigger_type, ref);
0233 }
0234 }
0235 } else {
0236 if ((std::abs(EtaSC) < 1.479 && vali >= cutRegularEB_) || (std::abs(EtaSC) >= 1.479 && vali >= cutRegularEE_)) {
0237 n++;
0238 filterproduct.addObject(trigger_type, ref);
0239 continue;
0240 }
0241 if (energy > 0. && (cutOverEEB_ > 0. || cutOverEEE_ > 0. || cutOverE2EB_ > 0. || cutOverE2EE_ > 0.)) {
0242 if ((std::abs(EtaSC) < 1.479 && vali / energy >= cutOverEEB_) ||
0243 (std::abs(EtaSC) >= 1.479 && vali / energy >= cutOverEEE_)) {
0244 n++;
0245 filterproduct.addObject(trigger_type, ref);
0246 continue;
0247 }
0248 if ((std::abs(EtaSC) < 1.479 && vali / (energy * energy) >= cutOverE2EB_) ||
0249 (std::abs(EtaSC) >= 1.479 && vali / (energy * energy) >= cutOverE2EE_)) {
0250 n++;
0251 filterproduct.addObject(trigger_type, ref);
0252 }
0253 }
0254 }
0255 }
0256
0257
0258 bool accept(n >= ncandcut_);
0259
0260 return accept;
0261 }
0262
0263 typedef HLTGenericFilter<reco::RecoEcalCandidate> HLTEgammaGenericFilter;
0264 typedef HLTGenericFilter<reco::RecoChargedCandidate> HLTMuonGenericFilter;
0265 DEFINE_FWK_MODULE(HLTEgammaGenericFilter);
0266 DEFINE_FWK_MODULE(HLTMuonGenericFilter);