File indexing completed on 2024-04-06 12:18:34
0001
0002
0003
0004
0005
0006
0007
0008 #include "FWCore/Framework/interface/ESHandle.h"
0009 #include "FWCore/Framework/interface/EventSetup.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0012 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014 #include "FWCore/Utilities/interface/InputTag.h"
0015 #include "DataFormats/Common/interface/Ref.h"
0016 #include "DataFormats/Common/interface/Handle.h"
0017 #include "DataFormats/Math/interface/deltaPhi.h"
0018 #include "HLTrigger/JetMET/interface/HLTMonoJetFilter.h"
0019 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0020
0021
0022
0023
0024 template <typename T>
0025 HLTMonoJetFilter<T>::HLTMonoJetFilter(const edm::ParameterSet& iConfig)
0026 : HLTFilter(iConfig),
0027 inputJetTag_(iConfig.template getParameter<edm::InputTag>("inputJetTag")),
0028 maxPtSecondJet_(iConfig.template getParameter<double>("maxPtSecondJet")),
0029 maxDeltaPhi_(iConfig.template getParameter<double>("maxDeltaPhi")),
0030 triggerType_(iConfig.template getParameter<int>("triggerType")) {
0031 m_theObjectToken = consumes<std::vector<T>>(inputJetTag_);
0032 LogDebug("") << "HLTMonoJetFilter: Input/maxPtSecondJet/maxDeltaPhi/triggerType : " << inputJetTag_.encode() << " "
0033 << maxPtSecondJet_ << " " << maxDeltaPhi_ << " " << triggerType_;
0034 }
0035
0036 template <typename T>
0037 HLTMonoJetFilter<T>::~HLTMonoJetFilter() = default;
0038
0039 template <typename T>
0040 void HLTMonoJetFilter<T>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0041 edm::ParameterSetDescription desc;
0042 makeHLTFilterDescription(desc);
0043 desc.add<edm::InputTag>("inputJetTag", edm::InputTag("hltAntiKT5ConvPFJets"));
0044 desc.add<double>("maxPtSecondJet", 9999.);
0045 desc.add<double>("maxDeltaPhi", 99.);
0046 desc.add<int>("triggerType", trigger::TriggerJet);
0047 descriptions.add(defaultModuleLabel<HLTMonoJetFilter<T>>(), desc);
0048 }
0049
0050
0051
0052
0053 template <typename T>
0054 bool HLTMonoJetFilter<T>::hltFilter(edm::Event& iEvent,
0055 const edm::EventSetup& iSetup,
0056 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0057 using namespace std;
0058 using namespace edm;
0059 using namespace reco;
0060 using namespace trigger;
0061
0062 typedef vector<T> TCollection;
0063 typedef Ref<TCollection> TRef;
0064
0065
0066 if (saveTags())
0067 filterproduct.addCollectionTag(inputJetTag_);
0068
0069
0070 TRef ref1, ref2;
0071
0072
0073 Handle<TCollection> objects;
0074 iEvent.getByToken(m_theObjectToken, objects);
0075
0076
0077 int n(0);
0078
0079 if (!objects->empty()) {
0080 int countJet(0);
0081 double jet1Phi = 0.;
0082 double jet2Phi = 0.;
0083 double jet2Pt = 0.;
0084
0085 typename TCollection::const_iterator i(objects->begin());
0086 for (; i != objects->end(); i++) {
0087 if (countJet == 0) {
0088 ref1 = TRef(objects, distance(objects->begin(), i));
0089 jet1Phi = i->phi();
0090 }
0091 if (countJet == 1) {
0092 ref2 = TRef(objects, distance(objects->begin(), i));
0093 jet2Pt = i->pt();
0094 jet2Phi = i->phi();
0095 }
0096 countJet++;
0097 if (countJet >= 2)
0098 break;
0099 }
0100
0101 if (countJet == 1) {
0102 n = 1;
0103 } else if (countJet > 1 && jet2Pt < maxPtSecondJet_) {
0104 n = 1;
0105 } else if (countJet > 1 && jet2Pt >= maxPtSecondJet_) {
0106 double Dphi = std::abs(deltaPhi(jet1Phi, jet2Phi));
0107 if (Dphi >= maxDeltaPhi_)
0108 n = -1;
0109 else
0110 n = 1;
0111 } else {
0112 n = -1;
0113 }
0114
0115 if (n == 1) {
0116 filterproduct.addObject(triggerType_, ref1);
0117 if (countJet > 1)
0118 filterproduct.addObject(triggerType_, ref2);
0119 }
0120 }
0121
0122 bool accept(n == 1);
0123
0124 return accept;
0125 }