File indexing completed on 2023-03-17 11:09:36
0001
0002
0003
0004
0005
0006
0007
0008 #include "HLTrigger/JetMET/interface/HLTNVFilter.h"
0009
0010 #include "DataFormats/Common/interface/Handle.h"
0011
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013
0014 #include "DataFormats/METReco/interface/CaloMET.h"
0015
0016 #include "FWCore/Framework/interface/ESHandle.h"
0017 #include "FWCore/Framework/interface/EventSetup.h"
0018 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0019 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0020
0021
0022
0023
0024 HLTNVFilter::HLTNVFilter(const edm::ParameterSet& iConfig) : HLTFilter(iConfig) {
0025 inputJetTag_ = iConfig.getParameter<edm::InputTag>("inputJetTag");
0026 inputMETTag_ = iConfig.getParameter<edm::InputTag>("inputMETTag");
0027 minNV_ = iConfig.getParameter<double>("minNV");
0028 minEtjet1_ = iConfig.getParameter<double>("minEtJet1");
0029 minEtjet2_ = iConfig.getParameter<double>("minEtJet2");
0030 m_theJetToken = consumes<reco::CaloJetCollection>(inputJetTag_);
0031 m_theMETToken = consumes<trigger::TriggerFilterObjectWithRefs>(inputMETTag_);
0032 }
0033
0034 HLTNVFilter::~HLTNVFilter() = default;
0035
0036 void HLTNVFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0037 edm::ParameterSetDescription desc;
0038 makeHLTFilterDescription(desc);
0039 desc.add<edm::InputTag>("inputJetTag", edm::InputTag("iterativeCone5CaloJets"));
0040 desc.add<edm::InputTag>("inputMETTag", edm::InputTag("hlt1MET60"));
0041 desc.add<double>("minEtJet2", 20.);
0042 desc.add<double>("minEtJet1", 80.);
0043 desc.add<double>("minNV", 0.1);
0044 descriptions.add("hltNVFilter", desc);
0045 }
0046
0047
0048 bool HLTNVFilter::hltFilter(edm::Event& iEvent,
0049 const edm::EventSetup& iSetup,
0050 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0051 using namespace std;
0052 using namespace edm;
0053 using namespace reco;
0054 using namespace trigger;
0055
0056
0057 if (saveTags()) {
0058 filterproduct.addCollectionTag(inputJetTag_);
0059 filterproduct.addCollectionTag(inputMETTag_);
0060 }
0061
0062 Handle<CaloJetCollection> recocalojets;
0063 iEvent.getByToken(m_theJetToken, recocalojets);
0064 Handle<trigger::TriggerFilterObjectWithRefs> metcal;
0065 iEvent.getByToken(m_theMETToken, metcal);
0066
0067
0068 int n(0);
0069
0070 if (recocalojets->size() > 1) {
0071
0072
0073 double etjet1 = 0.;
0074 double etjet2 = 0.;
0075 double etmiss = 0.;
0076 int countjets = 0;
0077
0078 VRcalomet vrefMET;
0079 metcal->getObjects(TriggerMET, vrefMET);
0080 CaloMETRef metRef = vrefMET.at(0);
0081 etmiss = vrefMET.at(0)->et();
0082
0083 CaloJetRef ref1, ref2;
0084 for (auto recocalojet = recocalojets->begin(); recocalojet <= (recocalojets->begin() + 1); recocalojet++) {
0085 if (countjets == 0) {
0086 etjet1 = recocalojet->et();
0087 ref1 = CaloJetRef(recocalojets, distance(recocalojets->begin(), recocalojet));
0088 }
0089 if (countjets == 1) {
0090 etjet2 = recocalojet->et();
0091 ref2 = CaloJetRef(recocalojets, distance(recocalojets->begin(), recocalojet));
0092 }
0093 countjets++;
0094 }
0095
0096 double NV = (etmiss * etmiss - (etjet1 - etjet2) * (etjet1 - etjet2)) / (etjet2 * etjet2);
0097 if (etjet1 > minEtjet1_ && etjet2 > minEtjet2_ && NV > minNV_) {
0098 filterproduct.addObject(TriggerMET, metRef);
0099 filterproduct.addObject(TriggerJet, ref1);
0100 filterproduct.addObject(TriggerJet, ref2);
0101 n++;
0102 }
0103
0104 }
0105
0106
0107 bool accept(n >= 1);
0108
0109 return accept;
0110 }