CMS 3D CMS Logo

DiJetVarProducer.cc
Go to the documentation of this file.
2 
4 
7 
10 
12 
14 
15 #include "TLorentzVector.h"
16 #include "TVector3.h"
17 
18 #include <memory>
19 #include <vector>
20 
21 //
22 // constructors and destructor
23 //
25  : inputJetTag_(iConfig.getParameter<edm::InputTag>("inputJetTag")),
26  wideJetDeltaR_(iConfig.getParameter<double>("wideJetDeltaR")) {
27  // register your products
28  // produces<std::vector<double> >("dijetvariables");
29  produces<std::vector<math::PtEtaPhiMLorentzVector>>("widejets");
30 
31  // set Token(-s)
32  inputJetTagToken_ = consumes<reco::CaloJetCollection>(iConfig.getParameter<edm::InputTag>("inputJetTag"));
33 
34  LogDebug("") << "Input Jet Tag: " << inputJetTag_.encode() << " ";
35  LogDebug("") << "Radius Parameter Wide Jet: " << wideJetDeltaR_ << ".";
36 }
37 
39 
40 // ------------ method called to produce the data ------------
42  using namespace std;
43  using namespace edm;
44  using namespace reco;
45 
46  // ## The output collections
47  // std::unique_ptr<std::vector<double> > dijetvariables(new
48  // std::vector<double>);
49  std::unique_ptr<std::vector<math::PtEtaPhiMLorentzVector>> widejets(new std::vector<math::PtEtaPhiMLorentzVector>);
50 
51  // ## Get jet collection
53  iEvent.getByToken(inputJetTagToken_, calojets_handle);
54  // cout << "size: " << calojets_handle->size() << endl;
55 
56  // ## Wide Jet Algorithm
57  // At least two jets
58  if (calojets_handle->size() >= 2) {
59  TLorentzVector wj1_tmp;
60  TLorentzVector wj2_tmp;
61  TLorentzVector wj1;
62  TLorentzVector wj2;
63  TLorentzVector wdijet;
64 
65  // // Loop over all the input jets
66  // for(reco::CaloJetCollection::const_iterator it =
67  // calojets_handle->begin(); it != calojets_handle->end(); ++it)
68  // {
69  // cout << "jet: " << it->pt() << " " << it->eta() << " " << it->phi()
70  // << endl;
71  // }
72 
73  // Find two leading jets
74  TLorentzVector jet1, jet2;
75 
76  reco::CaloJetCollection::const_iterator j1 = calojets_handle->begin();
77  reco::CaloJetCollection::const_iterator j2 = j1;
78  ++j2;
79 
80  jet1.SetPtEtaPhiM(j1->pt(), j1->eta(), j1->phi(), j1->mass());
81  jet2.SetPtEtaPhiM(j2->pt(), j2->eta(), j2->phi(), j2->mass());
82 
83  // cout << "j1: " << jet1.Pt() << " " << jet1.Eta() << " " << jet1.Phi() <<
84  // endl; cout << "j2: " << jet2.Pt() << " " << jet2.Eta() << " " <<
85  // jet2.Phi() << endl;
86 
87  // Create wide jets (radiation recovery algorithm)
88  for (reco::CaloJetCollection::const_iterator it = calojets_handle->begin(); it != calojets_handle->end(); ++it) {
89  TLorentzVector currentJet;
90  currentJet.SetPtEtaPhiM(it->pt(), it->eta(), it->phi(), it->mass());
91 
92  double DeltaR1 = currentJet.DeltaR(jet1);
93  double DeltaR2 = currentJet.DeltaR(jet2);
94 
95  if (DeltaR1 < DeltaR2 && DeltaR1 < wideJetDeltaR_) {
96  wj1_tmp += currentJet;
97  } else if (DeltaR2 < wideJetDeltaR_) {
98  wj2_tmp += currentJet;
99  }
100  }
101 
102  // Re-order the wide jets in pT
103  if (wj1_tmp.Pt() > wj2_tmp.Pt()) {
104  wj1 = wj1_tmp;
105  wj2 = wj2_tmp;
106  } else {
107  wj1 = wj2_tmp;
108  wj2 = wj1_tmp;
109  }
110 
111  // Create dijet system
112  wdijet = wj1 + wj2;
113 
114  // cout << "j1 wide: " << wj1.Pt() << " " << wj1.Eta() << " " <<
115  // wj1.Phi() << " " << wj1.M() << endl; cout << "j2 wide: " <<
116  // wj2.Pt() << " " << wj2.Eta() << " " << wj2.Phi() << " " << wj2.M()
117  // << endl; cout << "MJJWide: " << wdijet.M() << endl; cout <<
118  // "DeltaEtaJJWide: " << fabs(wj1.Eta()-wj2.Eta()) << endl; cout <<
119  // "DeltaPhiJJWide: " << fabs(wj1.DeltaPhi(wj2)) << endl;
120 
121  // // Put variables in the container
122  // dijetvariables->push_back( wdijet.M() ); //0 =
123  // MJJWide dijetvariables->push_back( fabs(wj1.Eta()-wj2.Eta()) );
124  // //1 = DeltaEtaJJWide dijetvariables->push_back(
125  // fabs(wj1.DeltaPhi(wj2)) ); //2 = DeltaPhiJJWide
126 
127  // Put widejets in the container
128  math::PtEtaPhiMLorentzVector wj1math(wj1.Pt(), wj1.Eta(), wj1.Phi(), wj1.M());
129  math::PtEtaPhiMLorentzVector wj2math(wj2.Pt(), wj2.Eta(), wj2.Phi(), wj2.M());
130  widejets->push_back(wj1math);
131  widejets->push_back(wj2math);
132  }
133  // else
134  // {
135  // // Put variables in the container
136  // dijetvariables->push_back( -1 ); //0 = MJJWide
137  // dijetvariables->push_back( -1 ); //1 = DeltaEtaJJWide
138  // dijetvariables->push_back( -1 ); //2 = DeltaPhiJJWide
139  // }
140 
141  // ## Put objects in the Event
142  // iEvent.put(std::move(dijetvariables), "dijetvariables");
143  iEvent.put(std::move(widejets), "widejets");
144 }
145 
~DiJetVarProducer() override
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
std::string encode() const
Definition: InputTag.cc:159
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
edm::EDGetTokenT< reco::CaloJetCollection > inputJetTagToken_
PtEtaPhiMLorentzVectorD PtEtaPhiMLorentzVector
Lorentz vector with cartesian internal representation.
Definition: LorentzVector.h:25
int iEvent
Definition: GenABIO.cc:224
edm::InputTag inputJetTag_
void produce(edm::Event &, const edm::EventSetup &) override
fixed size matrix
HLT enums.
DiJetVarProducer(const edm::ParameterSet &)
def move(src, dest)
Definition: eostools.py:511
#define LogDebug(id)