File indexing completed on 2024-04-06 12:18:31
0001 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0002
0003 #include "DataFormats/Common/interface/Handle.h"
0004
0005 #include "DataFormats/Common/interface/Ref.h"
0006 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0007 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0008
0009 #include "DataFormats/JetReco/interface/CaloJet.h"
0010 #include "DataFormats/JetReco/interface/CaloJetCollection.h"
0011
0012 #include "DataFormats/MuonReco/interface/Muon.h"
0013
0014 #include "FWCore/Framework/interface/MakerMacros.h"
0015 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0016
0017 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0018 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0019 #include "FWCore/Utilities/interface/InputTag.h"
0020
0021 #include "HLTrigger/JetMET/interface/HLTRFilter.h"
0022 #include "DataFormats/Common/interface/Ref.h"
0023
0024
0025
0026
0027 HLTRFilter::HLTRFilter(const edm::ParameterSet& iConfig)
0028 : HLTFilter(iConfig),
0029 inputTag_(iConfig.getParameter<edm::InputTag>("inputTag")),
0030 inputMetTag_(iConfig.getParameter<edm::InputTag>("inputMetTag")),
0031 doMuonCorrection_(iConfig.getParameter<bool>("doMuonCorrection")),
0032 min_R_(iConfig.getParameter<double>("minR")),
0033 min_MR_(iConfig.getParameter<double>("minMR")),
0034 DoRPrime_(iConfig.getParameter<bool>("doRPrime")),
0035 accept_NJ_(iConfig.getParameter<bool>("acceptNJ")),
0036 R_offset_((iConfig.existsAs<double>("R2Offset") ? iConfig.getParameter<double>("R2Offset") : 0)),
0037 MR_offset_((iConfig.existsAs<double>("MROffset") ? iConfig.getParameter<double>("MROffset") : 0)),
0038 R_MR_cut_((iConfig.existsAs<double>("RMRCut") ? iConfig.getParameter<double>("RMRCut") : -999999.))
0039
0040 {
0041 m_theInputToken = consumes<std::vector<math::XYZTLorentzVector>>(inputTag_);
0042 m_theMETToken = consumes<edm::View<reco::MET>>(inputMetTag_);
0043 LogDebug("") << "Inputs/minR/minMR/doRPrime/acceptNJ/R2Offset/MROffset/RMRCut : " << inputTag_.encode() << " "
0044 << inputMetTag_.encode() << " " << min_R_ << " " << min_MR_ << " " << accept_NJ_ << R_offset_ << " "
0045 << MR_offset_ << " " << R_MR_cut_ << ".";
0046
0047
0048 produces<reco::METCollection>();
0049 }
0050
0051 HLTRFilter::~HLTRFilter() = default;
0052
0053 void HLTRFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0054 edm::ParameterSetDescription desc;
0055 makeHLTFilterDescription(desc);
0056 desc.add<edm::InputTag>("inputTag", edm::InputTag("hltRHemisphere"));
0057 desc.add<edm::InputTag>("inputMetTag", edm::InputTag("hltMet"));
0058 desc.add<bool>("doMuonCorrection", false);
0059 desc.add<double>("minR", 0.3);
0060 desc.add<double>("minMR", 100.0);
0061 desc.add<bool>("doRPrime", false);
0062 desc.add<bool>("acceptNJ", true);
0063 desc.add<double>("R2Offset", 0.0);
0064 desc.add<double>("MROffset", 0.0);
0065 desc.add<double>("RMRCut", -999999.0);
0066 descriptions.add("hltRFilter", desc);
0067 }
0068
0069
0070
0071
0072
0073
0074 bool HLTRFilter::hltFilter(edm::Event& iEvent,
0075 const edm::EventSetup& iSetup,
0076 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0077 using namespace std;
0078 using namespace edm;
0079 using namespace reco;
0080
0081
0082 Handle<vector<math::XYZTLorentzVector>> hemispheres;
0083 iEvent.getByToken(m_theInputToken, hemispheres);
0084
0085
0086 edm::Handle<edm::View<reco::MET>> inputMet;
0087 iEvent.getByToken(m_theMETToken, inputMet);
0088
0089
0090 if (not hemispheres.isValid() or not inputMet.isValid())
0091 return false;
0092
0093 if (hemispheres
0094 ->empty()) {
0095 return accept_NJ_;
0096 }
0097
0098
0099
0100
0101 int nMuons;
0102 switch (hemispheres->size()) {
0103 case 2:
0104 nMuons = 0;
0105 break;
0106 case 5:
0107 nMuons = 1;
0108 break;
0109 case 10:
0110 nMuons = 2;
0111 break;
0112 default:
0113 return false;
0114 }
0115
0116
0117 TLorentzVector ja(hemispheres->at(0).x(), hemispheres->at(0).y(), hemispheres->at(0).z(), hemispheres->at(0).t());
0118 TLorentzVector jb(hemispheres->at(1).x(), hemispheres->at(1).y(), hemispheres->at(1).z(), hemispheres->at(1).t());
0119
0120 std::vector<math::XYZTLorentzVector> muonVec;
0121
0122 double MR = CalcMR(ja, jb);
0123 double R = CalcR(MR, ja, jb, inputMet, muonVec);
0124
0125 if (MR >= min_MR_ && R >= min_R_ && ((R * R - R_offset_) * (MR - MR_offset_)) >= R_MR_cut_) {
0126 addObjects(iEvent, filterproduct, MR, R * R);
0127 return true;
0128 }
0129 if (nMuons == 0) {
0130 addObjects(iEvent, filterproduct, MR, R * R);
0131 return false;
0132 }
0133
0134 ja.SetXYZT(hemispheres->at(3).x(), hemispheres->at(3).y(), hemispheres->at(3).z(), hemispheres->at(3).t());
0135 jb.SetXYZT(hemispheres->at(4).x(), hemispheres->at(4).y(), hemispheres->at(4).z(), hemispheres->at(4).t());
0136 muonVec.push_back(hemispheres->at(2));
0137
0138 MR = CalcMR(ja, jb);
0139 R = CalcR(MR, ja, jb, inputMet, muonVec);
0140 if (MR >= min_MR_ && R >= min_R_ && ((R * R - R_offset_) * (MR - MR_offset_)) >= R_MR_cut_) {
0141 addObjects(iEvent, filterproduct, MR, R * R);
0142 return true;
0143 }
0144
0145 if (nMuons == 1) {
0146 addObjects(iEvent, filterproduct, MR, R * R);
0147 return false;
0148 }
0149
0150 muonVec.pop_back();
0151
0152 ja.SetXYZT(hemispheres->at(6).x(), hemispheres->at(6).y(), hemispheres->at(6).z(), hemispheres->at(6).t());
0153 jb.SetXYZT(hemispheres->at(7).x(), hemispheres->at(7).y(), hemispheres->at(7).z(), hemispheres->at(7).t());
0154 muonVec.push_back(hemispheres->at(5));
0155
0156 MR = CalcMR(ja, jb);
0157 R = CalcR(MR, ja, jb, inputMet, muonVec);
0158 if (MR >= min_MR_ && R >= min_R_ && ((R * R - R_offset_) * (MR - MR_offset_)) >= R_MR_cut_) {
0159 addObjects(iEvent, filterproduct, MR, R * R);
0160 return true;
0161 }
0162
0163 ja.SetXYZT(hemispheres->at(8).x(), hemispheres->at(8).y(), hemispheres->at(8).z(), hemispheres->at(8).t());
0164 jb.SetXYZT(hemispheres->at(9).x(), hemispheres->at(9).y(), hemispheres->at(9).z(), hemispheres->at(9).t());
0165 muonVec.push_back(hemispheres->at(2));
0166
0167 MR = CalcMR(ja, jb);
0168 R = CalcR(MR, ja, jb, inputMet, muonVec);
0169 if (MR >= min_MR_ && R >= min_R_ && ((R * R - R_offset_) * (MR - MR_offset_)) >= R_MR_cut_) {
0170 addObjects(iEvent, filterproduct, MR, R * R);
0171 return true;
0172 }
0173
0174
0175 addObjects(iEvent, filterproduct, MR, R * R);
0176 return false;
0177 }
0178
0179 void HLTRFilter::addObjects(edm::Event& iEvent,
0180 trigger::TriggerFilterObjectWithRefs& filterproduct,
0181 double MR,
0182 double Rsq) const {
0183
0184 std::unique_ptr<reco::METCollection> razorObject(new reco::METCollection());
0185
0186 reco::MET::LorentzVector mrRsqP4(MR, Rsq, 0, 0);
0187 reco::MET::Point vtx(0, 0, 0);
0188
0189 reco::MET mrRsq(mrRsqP4, vtx);
0190 razorObject->push_back(mrRsq);
0191
0192 edm::RefProd<reco::METCollection> ref_before_put = iEvent.getRefBeforePut<reco::METCollection>();
0193
0194
0195 iEvent.put(std::move(razorObject));
0196 edm::Ref<reco::METCollection> mrRsqRef(ref_before_put, 0);
0197
0198
0199 if (saveTags())
0200 filterproduct.addCollectionTag(edm::InputTag(*moduleLabel()));
0201 filterproduct.addObject(trigger::TriggerMET, mrRsqRef);
0202 }
0203
0204 double HLTRFilter::CalcMR(TLorentzVector ja, TLorentzVector jb) {
0205 if (ja.Pt() <= 0.1)
0206 return -1;
0207
0208 ja.SetPtEtaPhiM(ja.Pt(), ja.Eta(), ja.Phi(), 0.0);
0209 jb.SetPtEtaPhiM(jb.Pt(), jb.Eta(), jb.Phi(), 0.0);
0210
0211 if (ja.Pt() > jb.Pt()) {
0212 TLorentzVector temp = ja;
0213 ja = jb;
0214 jb = temp;
0215 }
0216
0217 double A = ja.P();
0218 double B = jb.P();
0219 double az = ja.Pz();
0220 double bz = jb.Pz();
0221 TVector3 jaT, jbT;
0222 jaT.SetXYZ(ja.Px(), ja.Py(), 0.0);
0223 jbT.SetXYZ(jb.Px(), jb.Py(), 0.0);
0224 double ATBT = (jaT + jbT).Mag2();
0225
0226 double MR = sqrt((A + B) * (A + B) - (az + bz) * (az + bz) -
0227 (jbT.Dot(jbT) - jaT.Dot(jaT)) * (jbT.Dot(jbT) - jaT.Dot(jaT)) / (jaT + jbT).Mag2());
0228
0229 double mybeta = (jbT.Dot(jbT) - jaT.Dot(jaT)) / sqrt(ATBT * ((A + B) * (A + B) - (az + bz) * (az + bz)));
0230
0231 double mygamma = 1. / sqrt(1. - mybeta * mybeta);
0232
0233
0234 return MR * mygamma;
0235 }
0236
0237 double HLTRFilter::CalcR(double MR,
0238 TLorentzVector ja,
0239 TLorentzVector jb,
0240 edm::Handle<edm::View<reco::MET>> inputMet,
0241 const std::vector<math::XYZTLorentzVector>& muons) {
0242
0243 TVector3 met;
0244 met.SetPtEtaPhi((inputMet->front()).pt(), 0.0, (inputMet->front()).phi());
0245
0246 std::vector<math::XYZTLorentzVector>::const_iterator muonIt;
0247 for (muonIt = muons.begin(); muonIt != muons.end(); muonIt++) {
0248 TVector3 tmp;
0249 tmp.SetPtEtaPhi(muonIt->pt(), 0, muonIt->phi());
0250 met -= tmp;
0251 }
0252
0253 double MTR = sqrt(0.5 * (met.Mag() * (ja.Pt() + jb.Pt()) - met.Dot(ja.Vect() + jb.Vect())));
0254
0255
0256 return float(MTR) / float(MR);
0257 }
0258 DEFINE_FWK_MODULE(HLTRFilter);