File indexing completed on 2024-04-06 12:18:29
0001 #ifndef HLTrigger_HLTfilters_L1TJetFilterT_h
0002 #define HLTrigger_HLTfilters_L1TJetFilterT_h
0003
0004 #include <vector>
0005 #include <cmath>
0006 #include <iterator>
0007
0008 #include "DataFormats/Common/interface/Ref.h"
0009 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0010 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0014 #include "HLTrigger/HLTcore/interface/HLTFilter.h"
0015 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0016
0017 template <class T>
0018 class L1TJetFilterT : public HLTFilter {
0019 public:
0020 explicit L1TJetFilterT(const edm::ParameterSet&);
0021 ~L1TJetFilterT() override;
0022 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0023 bool hltFilter(edm::Event&,
0024 const edm::EventSetup&,
0025 trigger::TriggerFilterObjectWithRefs& filterproduct) const override;
0026
0027 private:
0028 edm::InputTag const l1tJetTag_;
0029 edm::EDGetTokenT<std::vector<T>> const l1tJetToken_;
0030
0031 double const minPt_;
0032 double const minEta_;
0033 double const maxEta_;
0034 int const minN_;
0035
0036 edm::ParameterSet scalings_;
0037 std::vector<double> barrelScalings_;
0038 std::vector<double> overlapScalings_;
0039 std::vector<double> endcapScalings_;
0040
0041 double offlineJetPt(double const pt, double const eta) const;
0042 };
0043
0044 template <class T>
0045 L1TJetFilterT<T>::L1TJetFilterT(const edm::ParameterSet& iConfig)
0046 : HLTFilter(iConfig),
0047 l1tJetTag_(iConfig.getParameter<edm::InputTag>("inputTag")),
0048 l1tJetToken_(consumes<std::vector<T>>(l1tJetTag_)),
0049 minPt_(iConfig.getParameter<double>("MinPt")),
0050 minEta_(iConfig.getParameter<double>("MinEta")),
0051 maxEta_(iConfig.getParameter<double>("MaxEta")),
0052 minN_(iConfig.getParameter<int>("MinN")) {
0053 scalings_ = iConfig.getParameter<edm::ParameterSet>("Scalings");
0054 barrelScalings_ = scalings_.getParameter<std::vector<double>>("barrel");
0055 overlapScalings_ = scalings_.getParameter<std::vector<double>>("overlap");
0056 endcapScalings_ = scalings_.getParameter<std::vector<double>>("endcap");
0057 }
0058
0059 template <class T>
0060 L1TJetFilterT<T>::~L1TJetFilterT() = default;
0061
0062 template <class T>
0063 void L1TJetFilterT<T>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0064 edm::ParameterSetDescription desc;
0065 makeHLTFilterDescription(desc);
0066 desc.add<edm::InputTag>("inputTag", edm::InputTag("ak4PFL1PuppiCorrected"));
0067 desc.add<double>("MinPt", -1.0);
0068 desc.add<double>("MinEta", -5.0);
0069 desc.add<double>("MaxEta", 5.0);
0070 desc.add<int>("MinN", 1);
0071
0072 edm::ParameterSetDescription descScalings;
0073 descScalings.add<std::vector<double>>("barrel", {0.0, 1.0, 0.0});
0074 descScalings.add<std::vector<double>>("overlap", {0.0, 1.0, 0.0});
0075 descScalings.add<std::vector<double>>("endcap", {0.0, 1.0, 0.0});
0076 desc.add<edm::ParameterSetDescription>("Scalings", descScalings);
0077
0078 descriptions.add(defaultModuleLabel<L1TJetFilterT<T>>(), desc);
0079 }
0080
0081 template <class T>
0082 bool L1TJetFilterT<T>::hltFilter(edm::Event& iEvent,
0083 const edm::EventSetup& iSetup,
0084 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0085
0086
0087
0088
0089
0090 if (saveTags()) {
0091 filterproduct.addCollectionTag(l1tJetTag_);
0092 }
0093
0094 auto const& l1tJets = iEvent.getHandle(l1tJetToken_);
0095
0096 int nJet(0);
0097 for (auto iJet = l1tJets->begin(); iJet != l1tJets->end(); ++iJet) {
0098 if (offlineJetPt(iJet->pt(), iJet->eta()) >= minPt_ && iJet->eta() <= maxEta_ && iJet->eta() >= minEta_) {
0099 ++nJet;
0100 edm::Ref<std::vector<T>> ref(l1tJets, std::distance(l1tJets->begin(), iJet));
0101 filterproduct.addObject(trigger::TriggerObjectType::TriggerL1PFJet, ref);
0102 }
0103 }
0104
0105
0106 return nJet >= minN_;
0107 }
0108
0109 template <class T>
0110 double L1TJetFilterT<T>::offlineJetPt(double const pt, double const eta) const {
0111 if (std::abs(eta) < 1.5)
0112 return (barrelScalings_.at(0) + pt * barrelScalings_.at(1) + pt * pt * barrelScalings_.at(2));
0113 else if (std::abs(eta) < 2.4)
0114 return (overlapScalings_.at(0) + pt * overlapScalings_.at(1) + pt * pt * overlapScalings_.at(2));
0115 else
0116 return (endcapScalings_.at(0) + pt * endcapScalings_.at(1) + pt * pt * endcapScalings_.at(2));
0117
0118
0119
0120
0121
0122
0123
0124
0125
0126
0127
0128
0129
0130 }
0131
0132 #endif