rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

CMS_2015_PAS_FSQ_15_007

Underlying Event Measurements with Leading Particles and Jets in pp collisions at 13 TeV
Experiment: CMS (LHC)
Status: VALIDATED
Authors:
  • Paolo Gunnellini
References:
  • CMS-PAS-FSQ-15-007
  • <http://cds.cern.ch/record/2104473?ln=en>
Beams: p+ p+
Beam energies: (6500.0, 6500.0) GeV
Run details:
  • Minimum Bias inelastic collisions

A measurement of the underlying event (UE) activity is performed in proton-proton collisions at the centre-of-mass energy of 13 TeV with the CMS experiment at the LHC. The measurement is performed using leading charged-particles as well as lead- ing charged-particle jets as reference objects. A leading charged-particle or charged- particle jet is required to be produced in the central pseudorapidity region ( |eta| < 2) jet and with transverse momentum pT greater than 0.5 (pT greater than 1) GeV for leading charged-particle (charge particle jet). The UE activity is measured in terms of the average multiplicity and scalar sum of pT of the charged-particles, in the azimuthal region transverse to the direction of the leading reference object.

Source code: CMS_2015_PAS_FSQ_15_007.cc
  1// -*- C++ -*-
  2
  3#include "Rivet/Analysis.hh"
  4#include "Rivet/Projections/FinalState.hh"
  5#include "Rivet/Projections/ChargedFinalState.hh"
  6#include "Rivet/Projections/FastJets.hh"
  7
  8using namespace std;
  9
 10namespace Rivet {
 11
 12  // UE charged particles vs. leading jet
 13  class CMS_2015_PAS_FSQ_15_007 : public Analysis {
 14  public:
 15    /// Constructor
 16    CMS_2015_PAS_FSQ_15_007() : Analysis("CMS_2015_PAS_FSQ_15_007") {}
 17
 18    void init() {
 19      const ChargedFinalState cfs(Cuts::abseta < 2.0 && Cuts::pt > 500 * MeV);
 20      declare(cfs, "CFS");
 21
 22      const ChargedFinalState cfsforjet(Cuts::abseta < 2.5 && Cuts::pt > 500 * MeV);
 23      const FastJets jetpro(cfsforjet, FastJets::SISCONE, 0.5);
 24      declare(jetpro, "Jets");
 25
 26      book(_h_PtSum_vs_leadTrackPt_transMin, 1, 1, 1);
 27      book(_h_PtSum_vs_leadTrackPt_transMax, 2, 1, 1);
 28      book(_h_PtSum_vs_leadTrackPt_transDiff, 3, 1, 1);
 29      book(_h_PtSum_vs_leadTrackPt_transAvg, 4, 1, 1);
 30
 31      book(_h_Nch_vs_leadTrackPt_transMin, 5, 1, 1);
 32      book(_h_Nch_vs_leadTrackPt_transMax, 6, 1, 1);
 33      book(_h_Nch_vs_leadTrackPt_transDiff, 7, 1, 1);
 34      book(_h_Nch_vs_leadTrackPt_transAvg, 8, 1, 1);
 35
 36      book(_h_PtSum_vs_leadJetPt_transMin, 9, 1, 1);
 37      book(_h_PtSum_vs_leadJetPt_transMax, 10, 1, 1);
 38      book(_h_PtSum_vs_leadJetPt_transDiff, 11, 1, 1);
 39      book(_h_PtSum_vs_leadJetPt_transAvg, 12, 1, 1);
 40
 41      book(_h_Nch_vs_leadJetPt_transMin, 13, 1, 1);
 42      book(_h_Nch_vs_leadJetPt_transMax, 14, 1, 1);
 43      book(_h_Nch_vs_leadJetPt_transDiff, 15, 1, 1);
 44      book(_h_Nch_vs_leadJetPt_transAvg, 16, 1, 1);
 45    }
 46
 47    /// Perform the per-event analysis
 48    void analyze(const Event& event) {
 49      // Find the lead jet, applying a restriction that the jets must be within |eta| < 2.
 50      FourMomentum p_leadjet, p_leadtrack;
 51      for (const Jet& j : apply<FastJets>(event, "Jets").jetsByPt(1.0 * GeV)) {
 52        if (j.abseta() < 2.0) {
 53          p_leadjet = j.momentum();
 54          break;
 55        }
 56      }
 57
 58      for (const Particle& j : apply<ChargedFinalState>(event, "CFS").particlesByPt(0.5 * GeV)) {
 59        if (j.abseta() < 2.0) {
 60          p_leadtrack = j.momentum();
 61          break;
 62        }
 63      }
 64
 65      if (p_leadjet.isZero() && p_leadtrack.isZero())
 66        vetoEvent;
 67      const double phileadjet = p_leadjet.phi();
 68      const double pTleadjet = p_leadjet.pT();
 69
 70      const double phileadtrack = p_leadtrack.phi();
 71      const double pTleadtrack = p_leadtrack.pT();
 72
 73      Particles particles = apply<ChargedFinalState>(event, "CFS").particlesByPt();
 74
 75      int nTransverse_leadjet = 0;
 76      double ptSumTransverse_leadjet = 0.;
 77      int nTransverse1_leadjet = 0;
 78      double ptSumTransverse1_leadjet = 0.;
 79      int nTransverse2_leadjet = 0;
 80      double ptSumTransverse2_leadjet = 0.;
 81      int nTransverseMin_leadjet = 0;
 82      double ptSumTransverseMin_leadjet = 0.;
 83      int nTransverseMax_leadjet = 0;
 84      double ptSumTransverseMax_leadjet = 0.;
 85      int nTowards_leadjet = 0;
 86      double ptSumTowards_leadjet = 0.;
 87      int nAway_leadjet = 0;
 88      double ptSumAway_leadjet = 0.;
 89
 90      int nTransverse_leadtrack = 0;
 91      double ptSumTransverse_leadtrack = 0.;
 92      int nTransverse1_leadtrack = 0;
 93      double ptSumTransverse1_leadtrack = 0.;
 94      int nTransverse2_leadtrack = 0;
 95      double ptSumTransverse2_leadtrack = 0.;
 96      int nTransverseMin_leadtrack = 0;
 97      double ptSumTransverseMin_leadtrack = 0.;
 98      int nTransverseMax_leadtrack = 0;
 99      double ptSumTransverseMax_leadtrack = 0.;
100      int nTowards_leadtrack = 0;
101      double ptSumTowards_leadtrack = 0.;
102      int nAway_leadtrack = 0;
103      double ptSumAway_leadtrack = 0.;
104
105      for (const Particle& p : particles) {
106        const double pT = p.pT() / GeV;
107
108        if (!p_leadjet.isZero()) {
109          double dphi_leadjet = p.phi() - phileadjet;
110          while (dphi_leadjet > PI) {
111            dphi_leadjet = dphi_leadjet - 2.0 * PI;
112          }
113          while (dphi_leadjet < -PI) {
114            dphi_leadjet = dphi_leadjet + 2. * PI;
115          }
116
117          if (dphi_leadjet > PI / 3. && dphi_leadjet < PI * 2. / 3.) {  // Transverse1 region
118            nTransverse_leadjet++;
119            ptSumTransverse_leadjet += pT;
120            nTransverse1_leadjet++;
121            ptSumTransverse1_leadjet += pT;
122          }
123
124          if (dphi_leadjet < -PI / 3. && dphi_leadjet > -PI * 2. / 3.) {  // Transverse2 region
125            nTransverse_leadjet++;
126            ptSumTransverse_leadjet += pT;
127            nTransverse2_leadjet++;
128            ptSumTransverse2_leadjet += pT;
129          }
130
131          if (fabs(dphi_leadjet) < PI / 3.) {  // Toward region
132            nTowards_leadjet++;
133            ptSumTowards_leadjet += pT;
134          }
135
136          if (fabs(dphi_leadjet) > 2. * PI / 3.) {  // Away region
137            nAway_leadjet++;
138            ptSumAway_leadjet += pT;
139          }
140
141        }  //jet found
142
143        if (!p_leadtrack.isZero()) {
144          double dphi_leadtrack = p.phi() - phileadtrack;
145          while (dphi_leadtrack > PI) {
146            dphi_leadtrack = dphi_leadtrack - 2.0 * PI;
147          }
148          while (dphi_leadtrack < -PI) {
149            dphi_leadtrack = dphi_leadtrack + 2. * PI;
150          }
151
152          if (dphi_leadtrack > PI / 3. && dphi_leadtrack < PI * 2. / 3.) {  // Transverse1 region
153            nTransverse_leadtrack++;
154            ptSumTransverse_leadtrack += pT;
155            nTransverse1_leadtrack++;
156            ptSumTransverse1_leadtrack += pT;
157          }
158
159          if (dphi_leadtrack < -PI / 3. && dphi_leadtrack > -PI * 2. / 3.) {  // Transverse2 region
160            nTransverse_leadtrack++;
161            ptSumTransverse_leadtrack += pT;
162            nTransverse2_leadtrack++;
163            ptSumTransverse2_leadtrack += pT;
164          }
165
166          if (fabs(dphi_leadtrack) < PI / 3.) {  // Toward region
167            nTowards_leadtrack++;
168            ptSumTowards_leadtrack += pT;
169          }
170
171          if (fabs(dphi_leadtrack) > 2. * PI / 3.) {  // Away region
172            nAway_leadtrack++;
173            ptSumAway_leadtrack += pT;
174          }
175
176        }  //track found
177
178      }  //Loop over particles
179
180      const double fullarea = 8. / 3. * PI;
181      const double halfarea = 4. / 3. * PI;
182
183      if (!p_leadjet.isZero()) {
184        if (nTransverse2_leadjet > nTransverse1_leadjet) {
185          nTransverseMax_leadjet = nTransverse2_leadjet;
186          nTransverseMin_leadjet = nTransverse1_leadjet;
187        }
188
189        else {
190          nTransverseMax_leadjet = nTransverse1_leadjet;
191          nTransverseMin_leadjet = nTransverse2_leadjet;
192        }
193
194        if (ptSumTransverse2_leadjet > ptSumTransverse1_leadjet) {
195          ptSumTransverseMax_leadjet = ptSumTransverse2_leadjet;
196          ptSumTransverseMin_leadjet = ptSumTransverse1_leadjet;
197        }
198
199        else {
200          ptSumTransverseMax_leadjet = ptSumTransverse1_leadjet;
201          ptSumTransverseMin_leadjet = ptSumTransverse2_leadjet;
202        }
203
204        _h_Nch_vs_leadJetPt_transDiff->fill(pTleadjet / GeV,
205                                            1. / halfarea * (nTransverseMax_leadjet - nTransverseMin_leadjet));
206        _h_PtSum_vs_leadJetPt_transDiff->fill(
207            pTleadjet / GeV, 1. / halfarea * (ptSumTransverseMax_leadjet - ptSumTransverseMin_leadjet));
208        _h_Nch_vs_leadJetPt_transAvg->fill(pTleadjet / GeV,
209                                           1. / fullarea * (nTransverseMax_leadjet + nTransverseMin_leadjet));
210        _h_PtSum_vs_leadJetPt_transAvg->fill(pTleadjet / GeV,
211                                             1. / fullarea * (ptSumTransverseMax_leadjet + ptSumTransverseMin_leadjet));
212        _h_Nch_vs_leadJetPt_transMax->fill(pTleadjet / GeV, 1. / halfarea * nTransverseMax_leadjet);
213        _h_PtSum_vs_leadJetPt_transMax->fill(pTleadjet / GeV, 1. / halfarea * ptSumTransverseMax_leadjet);
214        _h_Nch_vs_leadJetPt_transMin->fill(pTleadjet / GeV, 1. / halfarea * nTransverseMin_leadjet);
215        _h_PtSum_vs_leadJetPt_transMin->fill(pTleadjet / GeV, 1. / halfarea * ptSumTransverseMin_leadjet);
216
217      }  //for leading jet
218
219      if (!p_leadtrack.isZero()) {
220        if (nTransverse2_leadtrack > nTransverse1_leadtrack) {
221          nTransverseMax_leadtrack = nTransverse2_leadtrack;
222          nTransverseMin_leadtrack = nTransverse1_leadtrack;
223        }
224
225        else {
226          nTransverseMax_leadtrack = nTransverse1_leadtrack;
227          nTransverseMin_leadtrack = nTransverse2_leadtrack;
228        }
229
230        if (ptSumTransverse2_leadtrack > ptSumTransverse1_leadtrack) {
231          ptSumTransverseMax_leadtrack = ptSumTransverse2_leadtrack;
232          ptSumTransverseMin_leadtrack = ptSumTransverse1_leadtrack;
233        }
234
235        else {
236          ptSumTransverseMax_leadtrack = ptSumTransverse1_leadtrack;
237          ptSumTransverseMin_leadtrack = ptSumTransverse2_leadtrack;
238        }
239
240        _h_Nch_vs_leadTrackPt_transDiff->fill(pTleadtrack / GeV,
241                                              1. / halfarea * (nTransverseMax_leadtrack - nTransverseMin_leadtrack));
242        _h_PtSum_vs_leadTrackPt_transDiff->fill(
243            pTleadtrack / GeV, 1. / halfarea * (ptSumTransverseMax_leadtrack - ptSumTransverseMin_leadtrack));
244        _h_Nch_vs_leadTrackPt_transAvg->fill(pTleadtrack / GeV,
245                                             1. / fullarea * (nTransverseMax_leadtrack + nTransverseMin_leadtrack));
246        _h_PtSum_vs_leadTrackPt_transAvg->fill(
247            pTleadtrack / GeV, 1. / fullarea * (ptSumTransverseMax_leadtrack + ptSumTransverseMin_leadtrack));
248
249        _h_Nch_vs_leadTrackPt_transMax->fill(pTleadtrack / GeV, 1. / halfarea * nTransverseMax_leadtrack);
250        _h_PtSum_vs_leadTrackPt_transMax->fill(pTleadtrack / GeV, 1. / halfarea * ptSumTransverseMax_leadtrack);
251        _h_Nch_vs_leadTrackPt_transMin->fill(pTleadtrack / GeV, 1. / halfarea * nTransverseMin_leadtrack);
252        _h_PtSum_vs_leadTrackPt_transMin->fill(pTleadtrack / GeV, 1. / halfarea * ptSumTransverseMin_leadtrack);
253
254      }  //for leading track
255    }
256
257    /// Normalise histograms etc., after the run
258    void finalize() {}
259
260  private:
261    Profile1DPtr _h_Nch_vs_leadJetPt_transMax;
262    Profile1DPtr _h_PtSum_vs_leadJetPt_transMax;
263    Profile1DPtr _h_Nch_vs_leadJetPt_transMin;
264    Profile1DPtr _h_PtSum_vs_leadJetPt_transMin;
265    Profile1DPtr _h_Nch_vs_leadJetPt_transDiff;
266    Profile1DPtr _h_PtSum_vs_leadJetPt_transDiff;
267    Profile1DPtr _h_Nch_vs_leadJetPt_transAvg;
268    Profile1DPtr _h_PtSum_vs_leadJetPt_transAvg;
269
270    Profile1DPtr _h_Nch_vs_leadTrackPt_transMax;
271    Profile1DPtr _h_PtSum_vs_leadTrackPt_transMax;
272    Profile1DPtr _h_Nch_vs_leadTrackPt_transMin;
273    Profile1DPtr _h_PtSum_vs_leadTrackPt_transMin;
274    Profile1DPtr _h_Nch_vs_leadTrackPt_transDiff;
275    Profile1DPtr _h_PtSum_vs_leadTrackPt_transDiff;
276    Profile1DPtr _h_Nch_vs_leadTrackPt_transAvg;
277    Profile1DPtr _h_PtSum_vs_leadTrackPt_transAvg;
278  };
279
280  // This global object acts as a hook for the plugin system
281  RIVET_DECLARE_PLUGIN(CMS_2015_PAS_FSQ_15_007);
282}  // namespace Rivet