File indexing completed on 2021-02-14 13:09:35
0001 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0002
0003 #include "DataFormats/Common/interface/Handle.h"
0004
0005 #include "DataFormats/JetReco/interface/CaloJet.h"
0006 #include "DataFormats/JetReco/interface/CaloJetCollection.h"
0007
0008 #include "FWCore/Framework/interface/MakerMacros.h"
0009 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0010
0011 #include "FWCore/Utilities/interface/InputTag.h"
0012
0013 #include "DQM/DataScouting/interface/DiJetVarProducer.h"
0014
0015 #include "TLorentzVector.h"
0016 #include "TVector3.h"
0017
0018 #include <memory>
0019 #include <vector>
0020
0021
0022
0023
0024 DiJetVarProducer::DiJetVarProducer(const edm::ParameterSet &iConfig)
0025 : inputJetTag_(iConfig.getParameter<edm::InputTag>("inputJetTag")),
0026 wideJetDeltaR_(iConfig.getParameter<double>("wideJetDeltaR")) {
0027
0028
0029 produces<std::vector<math::PtEtaPhiMLorentzVector>>("widejets");
0030
0031
0032 inputJetTagToken_ = consumes<reco::CaloJetCollection>(iConfig.getParameter<edm::InputTag>("inputJetTag"));
0033
0034 LogDebug("") << "Input Jet Tag: " << inputJetTag_.encode() << " ";
0035 LogDebug("") << "Radius Parameter Wide Jet: " << wideJetDeltaR_ << ".";
0036 }
0037
0038 DiJetVarProducer::~DiJetVarProducer() {}
0039
0040
0041 void DiJetVarProducer::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0042 using namespace std;
0043 using namespace edm;
0044 using namespace reco;
0045
0046
0047
0048
0049 std::unique_ptr<std::vector<math::PtEtaPhiMLorentzVector>> widejets(new std::vector<math::PtEtaPhiMLorentzVector>);
0050
0051
0052 edm::Handle<reco::CaloJetCollection> calojets_handle;
0053 iEvent.getByToken(inputJetTagToken_, calojets_handle);
0054
0055
0056
0057
0058 if (calojets_handle->size() >= 2) {
0059 TLorentzVector wj1_tmp;
0060 TLorentzVector wj2_tmp;
0061 TLorentzVector wj1;
0062 TLorentzVector wj2;
0063 TLorentzVector wdijet;
0064
0065
0066
0067
0068
0069
0070
0071
0072
0073
0074 TLorentzVector jet1, jet2;
0075
0076 reco::CaloJetCollection::const_iterator j1 = calojets_handle->begin();
0077 reco::CaloJetCollection::const_iterator j2 = j1;
0078 ++j2;
0079
0080 jet1.SetPtEtaPhiM(j1->pt(), j1->eta(), j1->phi(), j1->mass());
0081 jet2.SetPtEtaPhiM(j2->pt(), j2->eta(), j2->phi(), j2->mass());
0082
0083
0084
0085
0086
0087
0088 for (reco::CaloJetCollection::const_iterator it = calojets_handle->begin(); it != calojets_handle->end(); ++it) {
0089 TLorentzVector currentJet;
0090 currentJet.SetPtEtaPhiM(it->pt(), it->eta(), it->phi(), it->mass());
0091
0092 double DeltaR1 = currentJet.DeltaR(jet1);
0093 double DeltaR2 = currentJet.DeltaR(jet2);
0094
0095 if (DeltaR1 < DeltaR2 && DeltaR1 < wideJetDeltaR_) {
0096 wj1_tmp += currentJet;
0097 } else if (DeltaR2 < wideJetDeltaR_) {
0098 wj2_tmp += currentJet;
0099 }
0100 }
0101
0102
0103 if (wj1_tmp.Pt() > wj2_tmp.Pt()) {
0104 wj1 = wj1_tmp;
0105 wj2 = wj2_tmp;
0106 } else {
0107 wj1 = wj2_tmp;
0108 wj2 = wj1_tmp;
0109 }
0110
0111
0112 wdijet = wj1 + wj2;
0113
0114
0115
0116
0117
0118
0119
0120
0121
0122
0123
0124
0125
0126
0127
0128 math::PtEtaPhiMLorentzVector wj1math(wj1.Pt(), wj1.Eta(), wj1.Phi(), wj1.M());
0129 math::PtEtaPhiMLorentzVector wj2math(wj2.Pt(), wj2.Eta(), wj2.Phi(), wj2.M());
0130 widejets->push_back(wj1math);
0131 widejets->push_back(wj2math);
0132 }
0133
0134
0135
0136
0137
0138
0139
0140
0141
0142
0143 iEvent.put(std::move(widejets), "widejets");
0144 }
0145
0146 DEFINE_FWK_MODULE(DiJetVarProducer);