File indexing completed on 2023-03-22 06:32:56
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022 #ifndef RecoTauTag_HLTProducers_L1TJetsMatching_h
0023 #define RecoTauTag_HLTProducers_L1TJetsMatching_h
0024
0025
0026 #include "FWCore/Framework/interface/Frameworkfwd.h"
0027 #include "FWCore/Framework/interface/global/EDProducer.h"
0028 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0029 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0030 #include "FWCore/Framework/interface/Event.h"
0031 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0032 #include "FWCore/Utilities/interface/InputTag.h"
0033 #include "DataFormats/Common/interface/Handle.h"
0034 #include "DataFormats/JetReco/interface/PFJetCollection.h"
0035 #include "DataFormats/TauReco/interface/PFTauFwd.h"
0036 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0037
0038 #include "Math/GenVector/VectorUtil.h"
0039 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0040 #include "FWCore/Utilities/interface/Exception.h"
0041 #include "DataFormats/JetReco/interface/PFJet.h"
0042
0043 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0044 #include "DataFormats/Math/interface/deltaR.h"
0045
0046 #include <vector>
0047 #include <utility>
0048 #include <tuple>
0049 #include <string>
0050
0051 template <typename T>
0052 class L1TJetsMatching : public edm::global::EDProducer<> {
0053 public:
0054 explicit L1TJetsMatching(const edm::ParameterSet&);
0055 ~L1TJetsMatching() override = default;
0056 void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0057 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0058
0059 private:
0060 std::pair<std::vector<T>, std::vector<T>> categorise(
0061 const std::vector<T>& pfMatchedJets, double pt1, double pt2, double pt3, double Mjj) const;
0062 std::tuple<std::vector<T>, std::vector<T>, std::vector<T>> categoriseVBFPlus2CentralJets(
0063 const std::vector<T>& pfMatchedJets, double pt1, double pt2, double pt3, double Mjj) const;
0064 const edm::EDGetTokenT<std::vector<T>> jetSrc_;
0065 const edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> jetTrigger_;
0066 const std::string matchingMode_;
0067 const double pt1Min_;
0068 const double pt2Min_;
0069 const double pt3Min_;
0070 const double mjjMin_;
0071 const double matchingR2_;
0072 };
0073
0074
0075
0076 template <typename T>
0077 std::pair<std::vector<T>, std::vector<T>> L1TJetsMatching<T>::categorise(
0078 const std::vector<T>& pfMatchedJets, double pt1, double pt2, double pt3, double Mjj) const {
0079 std::pair<std::vector<T>, std::vector<T>> output;
0080 unsigned int i1 = 0;
0081 unsigned int i2 = 0;
0082 double m2jj = 0;
0083 if (pfMatchedJets.size() > 1) {
0084 for (unsigned int i = 0; i < pfMatchedJets.size() - 1; i++) {
0085 const T& myJet1 = (pfMatchedJets)[i];
0086
0087 for (unsigned int j = i + 1; j < pfMatchedJets.size(); j++) {
0088 const T& myJet2 = (pfMatchedJets)[j];
0089
0090 const double m2jj_test = (myJet1.p4() + myJet2.p4()).M2();
0091
0092 if (m2jj_test > m2jj) {
0093 m2jj = m2jj_test;
0094 i1 = i;
0095 i2 = j;
0096 }
0097 }
0098 }
0099
0100 const T& myJet1 = (pfMatchedJets)[i1];
0101 const T& myJet2 = (pfMatchedJets)[i2];
0102 const double M2jj = (Mjj >= 0. ? Mjj * Mjj : -1.);
0103
0104 if ((m2jj > M2jj) && (myJet1.pt() >= pt1) && (myJet2.pt() > pt2)) {
0105 output.first.push_back(myJet1);
0106 output.first.push_back(myJet2);
0107 }
0108
0109 if ((m2jj > M2jj) && (myJet1.pt() < pt3) && (myJet1.pt() > pt2) && (myJet2.pt() > pt2)) {
0110 const T& myJetTest = (pfMatchedJets)[0];
0111 if (myJetTest.pt() > pt3) {
0112 output.second.push_back(myJet1);
0113 output.second.push_back(myJet2);
0114 output.second.push_back(myJetTest);
0115 }
0116 }
0117 }
0118
0119 return output;
0120 }
0121 template <typename T>
0122 std::tuple<std::vector<T>, std::vector<T>, std::vector<T>> L1TJetsMatching<T>::categoriseVBFPlus2CentralJets(
0123 const std::vector<T>& pfMatchedJets, double pt1, double pt2, double pt3, double Mjj) const {
0124 std::tuple<std::vector<T>, std::vector<T>, std::vector<T>> output;
0125 unsigned int i1 = 0;
0126 unsigned int i2 = 0;
0127
0128 double m2jj = 0;
0129 if (pfMatchedJets.size() > 1) {
0130 for (unsigned int i = 0; i < pfMatchedJets.size() - 1; i++) {
0131 const T& myJet1 = (pfMatchedJets)[i];
0132
0133 for (unsigned int j = i + 1; j < pfMatchedJets.size(); j++) {
0134 const T& myJet2 = (pfMatchedJets)[j];
0135
0136 const double m2jj_test = (myJet1.p4() + myJet2.p4()).M2();
0137
0138 if (m2jj_test > m2jj) {
0139 m2jj = m2jj_test;
0140 i1 = i;
0141 i2 = j;
0142 }
0143 }
0144 }
0145
0146 const T& myJet1 = (pfMatchedJets)[i1];
0147 const T& myJet2 = (pfMatchedJets)[i2];
0148 const double M2jj = (Mjj >= 0. ? Mjj * Mjj : -1.);
0149
0150 std::vector<T> vec4jets;
0151 vec4jets.reserve(4);
0152 std::vector<T> vec5jets;
0153 vec5jets.reserve(5);
0154 std::vector<T> vec6jets;
0155 vec6jets.reserve(6);
0156 if (pfMatchedJets.size() > 3) {
0157 if ((m2jj > M2jj) && (myJet1.pt() >= pt3) && (myJet2.pt() > pt2)) {
0158 vec4jets.push_back(myJet1);
0159 vec4jets.push_back(myJet2);
0160
0161 for (unsigned int i = 0; i < pfMatchedJets.size(); i++) {
0162 if (vec4jets.size() > 3)
0163 break;
0164 if (i == i1 or i == i2)
0165 continue;
0166 vec4jets.push_back(pfMatchedJets[i]);
0167 }
0168 }
0169
0170 if ((m2jj > M2jj) && (myJet1.pt() < pt1) && (myJet1.pt() < pt3) && (myJet1.pt() > pt2) &&
0171 (myJet2.pt() > pt2)) {
0172
0173 std::vector<unsigned int> idx_jets;
0174 idx_jets.reserve(pfMatchedJets.size() - 2);
0175
0176 for (unsigned int i = 0; i < pfMatchedJets.size(); i++) {
0177 if (i == i1 || i == i2)
0178 continue;
0179 if (pfMatchedJets[i].pt() > pt2) {
0180 idx_jets.push_back(i);
0181 }
0182 }
0183 if (idx_jets.size() == 3) {
0184 vec5jets.push_back(myJet1);
0185 vec5jets.push_back(myJet2);
0186 vec5jets.push_back(pfMatchedJets[idx_jets[0]]);
0187 vec5jets.push_back(pfMatchedJets[idx_jets[1]]);
0188 vec5jets.push_back(pfMatchedJets[idx_jets[2]]);
0189
0190 } else if (idx_jets.size() > 3) {
0191 vec6jets.push_back(myJet1);
0192 vec6jets.push_back(myJet2);
0193 vec6jets.push_back(pfMatchedJets[idx_jets[0]]);
0194 vec6jets.push_back(pfMatchedJets[idx_jets[1]]);
0195 vec6jets.push_back(pfMatchedJets[idx_jets[2]]);
0196 vec6jets.push_back(pfMatchedJets[idx_jets[3]]);
0197 }
0198 }
0199 }
0200
0201 output = std::make_tuple(vec4jets, vec5jets, vec6jets);
0202 }
0203
0204 return output;
0205 }
0206
0207 template <typename T>
0208 L1TJetsMatching<T>::L1TJetsMatching(const edm::ParameterSet& iConfig)
0209 : jetSrc_(consumes<std::vector<T>>(iConfig.getParameter<edm::InputTag>("JetSrc"))),
0210 jetTrigger_(consumes<trigger::TriggerFilterObjectWithRefs>(iConfig.getParameter<edm::InputTag>("L1JetTrigger"))),
0211 matchingMode_(iConfig.getParameter<std::string>("matchingMode")),
0212 pt1Min_(iConfig.getParameter<double>("pt1Min")),
0213 pt2Min_(iConfig.getParameter<double>("pt2Min")),
0214 pt3Min_(iConfig.getParameter<double>("pt3Min")),
0215 mjjMin_(iConfig.getParameter<double>("mjjMin")),
0216 matchingR2_(iConfig.getParameter<double>("matchingR") * iConfig.getParameter<double>("matchingR")) {
0217 if (matchingMode_ == "VBF") {
0218 produces<std::vector<T>>("TwoJets");
0219 produces<std::vector<T>>("ThreeJets");
0220 } else if (matchingMode_ == "VBFPlus2CentralJets") {
0221 produces<std::vector<T>>("FourJets");
0222 produces<std::vector<T>>("FiveJets");
0223 produces<std::vector<T>>("SixJets");
0224 } else {
0225 throw cms::Exception("InvalidConfiguration") << "invalid value for \"matchingMode\": " << matchingMode_
0226 << " (valid values are \"VBF\" and \"VBFPlus2CentralJets\")";
0227 }
0228 }
0229
0230 template <typename T>
0231 void L1TJetsMatching<T>::produce(edm::StreamID iSId, edm::Event& iEvent, const edm::EventSetup& iES) const {
0232 unique_ptr<std::vector<T>> pfMatchedJets(new std::vector<T>);
0233
0234
0235 edm::Handle<std::vector<T>> pfJets;
0236 iEvent.getByToken(jetSrc_, pfJets);
0237
0238 edm::Handle<trigger::TriggerFilterObjectWithRefs> l1TriggeredJets;
0239 iEvent.getByToken(jetTrigger_, l1TriggeredJets);
0240
0241
0242 l1t::JetVectorRef jetCandRefVec;
0243 l1TriggeredJets->getObjects(trigger::TriggerL1Jet, jetCandRefVec);
0244
0245 math::XYZPoint a(0., 0., 0.);
0246
0247
0248 for (unsigned int iJet = 0; iJet < pfJets->size(); iJet++) {
0249 const T& myJet = (*pfJets)[iJet];
0250 for (unsigned int iL1Jet = 0; iL1Jet < jetCandRefVec.size(); iL1Jet++) {
0251
0252
0253 if ((reco::deltaR2(myJet.p4(), jetCandRefVec[iL1Jet]->p4()) < matchingR2_) && (myJet.pt() > pt2Min_)) {
0254 pfMatchedJets->push_back(myJet);
0255 break;
0256 }
0257 }
0258 }
0259
0260 std::sort(pfMatchedJets->begin(), pfMatchedJets->end(), [](const T& j1, const T& j2) { return j1.pt() > j2.pt(); });
0261
0262 if (matchingMode_ == "VBF") {
0263 std::pair<std::vector<T>, std::vector<T>> output = categorise(*pfMatchedJets, pt1Min_, pt2Min_, pt3Min_, mjjMin_);
0264 auto output1 = std::make_unique<std::vector<T>>(output.first);
0265 auto output2 = std::make_unique<std::vector<T>>(output.second);
0266
0267 iEvent.put(std::move(output1), "TwoJets");
0268 iEvent.put(std::move(output2), "ThreeJets");
0269
0270 } else if (matchingMode_ == "VBFPlus2CentralJets") {
0271 std::tuple<std::vector<T>, std::vector<T>, std::vector<T>> output =
0272 categoriseVBFPlus2CentralJets(*pfMatchedJets, pt1Min_, pt2Min_, pt3Min_, mjjMin_);
0273 auto output1 = std::make_unique<std::vector<T>>(std::get<0>(output));
0274 auto output2 = std::make_unique<std::vector<T>>(std::get<1>(output));
0275 auto output3 = std::make_unique<std::vector<T>>(std::get<2>(output));
0276
0277 iEvent.put(std::move(output1), "FourJets");
0278 iEvent.put(std::move(output2), "FiveJets");
0279 iEvent.put(std::move(output3), "SixJets");
0280 }
0281 }
0282
0283 template <typename T>
0284 void L1TJetsMatching<T>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0285 edm::ParameterSetDescription desc;
0286 desc.add<edm::InputTag>("L1JetTrigger", edm::InputTag("hltL1DiJetVBF"))->setComment("Name of trigger filter");
0287 desc.add<edm::InputTag>("JetSrc", edm::InputTag("hltAK4PFJetsTightIDCorrected"))
0288 ->setComment("Input collection of PFJets");
0289 desc.add<std::string>("matchingMode", "VBF")
0290 ->setComment("Switch from Di/tri-jet (VBF) to Multi-jet (VBFPlus2CentralJets) matching");
0291 desc.add<double>("pt1Min", 110.0)->setComment("Minimal pT1 of PFJets to match");
0292 desc.add<double>("pt2Min", 35.0)->setComment("Minimal pT2 of PFJets to match");
0293 desc.add<double>("pt3Min", 110.0)->setComment("Minimum pT3 of PFJets to match");
0294 desc.add<double>("mjjMin", 650.0)->setComment("Minimal mjj of matched PFjets");
0295 desc.add<double>("matchingR", 0.5)->setComment("dR value used for matching");
0296 descriptions.setComment(
0297 "This module produces collection of PFJets matched to L1 Taus / Jets passing a HLT filter (Only p4 and vertex "
0298 "of returned PFJets are set).");
0299 descriptions.add(defaultModuleLabel<L1TJetsMatching<T>>(), desc);
0300 }
0301
0302 #endif