Rivet analyses referenceCMS_2014_I1303894Differential cross-section of $W$ bosons + jets in $pp$ collisions at $\sqrt{s}=7$ TeVExperiment: CMS (LHC) Inspire ID: 1303894 Status: VALIDATED Authors:
Beam energies: (3500.0, 3500.0) GeV Run details:
A study of jet production in association with $W$ bosons has been performed, in events with the $W$ decaying to a muon. Jets are required to have $pT > 30$ GeV and $|eta| < 2.4$. Muons are required to have $pT > 25$ and $|eta| < 2.1$. Jets are only considered if they are separated from the muon by $\Delta{R} > 0.5$. Muons are dressed with photons in a cone of $0.1$ around the muon. Source code: CMS_2014_I1303894.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Projections/VetoedFinalState.hh"
6#include "Rivet/Projections/PromptFinalState.hh"
7#include "Rivet/Projections/InvMassFinalState.hh"
8#include "Rivet/Projections/MissingMomentum.hh"
9#include "Rivet/Projections/WFinder.hh"
10#include "Rivet/Projections/DressedLeptons.hh"
11
12namespace Rivet {
13
14
15 /// @brief Differential cross-section of W bosons + jets in pp collisions at sqrt(s)=7 TeV
16 /// @author Darin Baumgartel (darinb@cern.ch)
17 ///
18 /// Based on Rivet analysis originally created by Anil Singh (anil@cern.ch), Lovedeep Saini (lovedeep@cern.ch)
19 class CMS_2014_I1303894 : public Analysis {
20 public:
21
22 /// Constructor
23 CMS_2014_I1303894()
24 : Analysis("CMS_2014_I1303894")
25 { }
26
27
28 // Book histograms and initialise projections before the run
29 void init() {
30 // Prompt leptons only, no test on nu flavour.
31 // Projections
32 const FinalState fs;
33 declare(fs, "FS");
34
35 MissingMomentum missing(fs);
36 declare(missing, "MET");
37
38 PromptFinalState pfs(fs);
39 IdentifiedFinalState bareMuons(pfs);
40 bareMuons.acceptIdPair(PID::MUON);
41 DressedLeptons muonClusters(fs, bareMuons, -1); //, Cuts::open(), false, false);
42 declare(muonClusters, "muonClusters");
43
44 IdentifiedFinalState neutrinos(pfs);
45 neutrinos.acceptIdPair(PID::NU_MU);
46 declare(neutrinos, "neutrinos");
47
48 VetoedFinalState jetFS(fs);
49 jetFS.addVetoOnThisFinalState(muonClusters);
50 jetFS.addVetoOnThisFinalState(neutrinos);
51 jetFS.vetoNeutrinos();
52 FastJets jetprojection(jetFS, FastJets::ANTIKT, 0.5);
53 declare(jetprojection, "Jets");
54
55 // Histograms
56 book(_histDPhiMuJet1 ,1,1,1);
57 book(_histDPhiMuJet2 ,2,1,1);
58 book(_histDPhiMuJet3 ,3,1,1);
59 book(_histDPhiMuJet4 ,4,1,1);
60
61 book(_histEtaJet1 ,5,1,1);
62 book(_histEtaJet2 ,6,1,1);
63 book(_histEtaJet3 ,7,1,1);
64 book(_histEtaJet4 ,8,1,1);
65
66 book(_histHT1JetInc ,9,1,1);
67 book(_histHT2JetInc ,10,1,1);
68 book(_histHT3JetInc ,11,1,1);
69 book(_histHT4JetInc ,12,1,1);
70
71 book(_histJet30MultExc ,13,1,1);
72 book(_histJet30MultInc ,14,1,1);
73
74 book(_histPtJet1 ,15,1,1);
75 book(_histPtJet2 ,16,1,1);
76 book(_histPtJet3 ,17,1,1);
77 book(_histPtJet4 ,18,1,1);
78
79 // Counters
80 book(_n_1jet, "n_1jet");
81 book(_n_2jet, "n_2jet");
82 book(_n_3jet, "n_3jet");
83 book(_n_4jet, "n_4jet");
84 book(_n_inclusivebinsummation, "n_inclusivebinsummation");
85 }
86
87
88 void analyze(const Event& event) {
89 // Get the dressed muon
90 const DressedLeptons& muonClusters = apply<DressedLeptons>(event, "muonClusters");
91 int nmu = muonClusters.dressedLeptons().size();
92 if (nmu < 1) vetoEvent;
93 DressedLepton dressedmuon = muonClusters.dressedLeptons()[0];
94 if (dressedmuon.momentum().abseta() > 2.1) vetoEvent;
95 if (dressedmuon.momentum().pT() < 25.0*GeV) vetoEvent;
96
97 // Get the muon neutrino
98 //const Particles& neutrinos = apply<FinalState>(event, "neutrinos").particlesByPt();
99
100 // Check that the muon and neutrino are not decay products of tau
101 if (dressedmuon.constituentLepton().hasAncestor( PID::TAU)) vetoEvent;
102 if (dressedmuon.constituentLepton().hasAncestor(-PID::TAU)) vetoEvent;
103
104 // Get the missing momentum
105 const MissingMomentum& met = apply<MissingMomentum>(event, "MET");
106 const double ptmet = met.visibleMomentum().pT();
107 const double phimet = (-met.visibleMomentum()).phi();
108
109 // Calculate MET and MT(mu,MET), and remove events with MT < 50 GeV
110 const double ptmuon = dressedmuon.pT();
111 const double phimuon = dressedmuon.phi();
112 const double mt_mumet = sqrt(2*ptmuon*ptmet*(1.0 - cos(phimet-phimuon)));
113
114 // Remove events in MT < 50 region
115 if (mt_mumet < 50*GeV) vetoEvent;
116
117 // Loop over jets and fill pt/eta/phi quantities in vectors
118 const Jets& jets_filtered = apply<FastJets>(event, "Jets").jetsByPt(0.0*GeV);
119 vector<float> finaljet_pT_list, finaljet_eta_list, finaljet_phi_list;
120 double htjets = 0.0;
121 for (size_t ii = 0; ii < jets_filtered.size(); ++ii) {
122 // Jet pT/eta/phi
123 double jet_pt = jets_filtered[ii].pT();
124 double jet_eta = jets_filtered[ii].eta();
125 double jet_phi = jets_filtered[ii].phi();
126
127 // Kinemetic cuts for jet acceptance
128 if (fabs(jet_eta) > 2.4) continue;
129 if (jet_pt < 30.0*GeV) continue;
130 if (deltaR(dressedmuon, jets_filtered[ii]) < 0.5) continue;
131
132 // Add jet to jet list and increases the HT variable
133 finaljet_pT_list.push_back(jet_pt);
134 finaljet_eta_list.push_back(jet_eta);
135 finaljet_phi_list.push_back(jet_phi);
136 htjets += fabs(jet_pt);
137 }
138
139
140 // Filling of histograms:
141 // Fill as many jets as there are into the exclusive jet multiplicity
142 if (!finaljet_pT_list.empty())
143 _histJet30MultExc->fill(finaljet_pT_list.size());
144
145 for (size_t ij = 0; ij < finaljet_pT_list.size(); ++ij) {
146 _histJet30MultInc->fill(ij+1);
147 _n_inclusivebinsummation->fill();
148 }
149
150 if (finaljet_pT_list.size() >= 1) {
151 _histPtJet1->fill(finaljet_pT_list[0]);
152 _histEtaJet1->fill(fabs(finaljet_eta_list[0]));
153 _histDPhiMuJet1->fill(deltaPhi(finaljet_phi_list[0], phimuon));
154 _histHT1JetInc->fill(htjets);
155 _n_1jet->fill();
156 }
157
158 if (finaljet_pT_list.size() >= 2) {
159 _histPtJet2->fill(finaljet_pT_list[1]);
160 _histEtaJet2->fill(fabs(finaljet_eta_list[1]));
161 _histDPhiMuJet2->fill(deltaPhi(finaljet_phi_list[1], phimuon));
162 _histHT2JetInc->fill(htjets);
163 _n_2jet->fill();
164 }
165
166 if (finaljet_pT_list.size() >= 3) {
167 _histPtJet3->fill(finaljet_pT_list[2]);
168 _histEtaJet3->fill(fabs(finaljet_eta_list[2]));
169 _histDPhiMuJet3->fill(deltaPhi(finaljet_phi_list[2], phimuon));
170 _histHT3JetInc->fill(htjets);
171 _n_3jet->fill();
172 }
173
174 if (finaljet_pT_list.size() >=4 ) {
175 _histPtJet4->fill(finaljet_pT_list[3]);
176 _histEtaJet4->fill(fabs(finaljet_eta_list[3]));
177 _histDPhiMuJet4->fill(deltaPhi(finaljet_phi_list[3], phimuon));
178 _histHT4JetInc-> fill(htjets);
179 _n_4jet->fill();
180 }
181
182 }
183
184
185 // Finalize the histograms.
186 void finalize() {
187
188 const double inclusive_cross_section = crossSection();
189 const double norm_1jet_histo = inclusive_cross_section*dbl(*_n_1jet)/sumOfWeights();
190 const double norm_2jet_histo = inclusive_cross_section*dbl(*_n_2jet)/sumOfWeights();
191 const double norm_3jet_histo = inclusive_cross_section*dbl(*_n_3jet)/sumOfWeights();
192 const double norm_4jet_histo = inclusive_cross_section*dbl(*_n_4jet)/sumOfWeights();
193 const double norm_incmultiplicity = inclusive_cross_section*dbl(*_n_inclusivebinsummation)/sumOfWeights();
194
195 normalize(_histJet30MultExc, norm_1jet_histo);
196 normalize(_histJet30MultInc, norm_incmultiplicity);
197
198 normalize(_histPtJet1, norm_1jet_histo);
199 normalize(_histHT1JetInc, norm_1jet_histo);
200 normalize(_histEtaJet1, norm_1jet_histo);
201 normalize(_histDPhiMuJet1, norm_1jet_histo);
202
203 normalize(_histPtJet2, norm_2jet_histo);
204 normalize(_histHT2JetInc, norm_2jet_histo);
205 normalize(_histEtaJet2, norm_2jet_histo);
206 normalize(_histDPhiMuJet2, norm_2jet_histo);
207
208 normalize(_histPtJet3, norm_3jet_histo);
209 normalize(_histHT3JetInc, norm_3jet_histo);
210 normalize(_histEtaJet3, norm_3jet_histo);
211 normalize(_histDPhiMuJet3, norm_3jet_histo);
212
213 normalize(_histPtJet4, norm_4jet_histo);
214 normalize(_histHT4JetInc, norm_4jet_histo);
215 normalize(_histEtaJet4, norm_4jet_histo);
216 normalize(_histDPhiMuJet4, norm_4jet_histo);
217 }
218
219
220 private:
221
222 Histo1DPtr _histJet30MultExc, _histJet30MultInc;
223 Histo1DPtr _histPtJet1, _histPtJet2, _histPtJet3, _histPtJet4;
224 Histo1DPtr _histEtaJet1, _histEtaJet2, _histEtaJet3, _histEtaJet4;
225 Histo1DPtr _histDPhiMuJet1, _histDPhiMuJet2, _histDPhiMuJet3, _histDPhiMuJet4;
226 Histo1DPtr _histHT1JetInc, _histHT2JetInc, _histHT3JetInc, _histHT4JetInc;
227
228 CounterPtr _n_1jet, _n_2jet, _n_3jet, _n_4jet, _n_inclusivebinsummation;
229
230 };
231
232
233 RIVET_DECLARE_PLUGIN(CMS_2014_I1303894);
234
235}
|