Rivet analyses referenceCMS_2011_I916908Traditional leading jet UE measurement at $\sqrt{s} = 0.9$ and 7 TeVExperiment: CMS (LHC) Inspire ID: 916908 Status: VALIDATED Authors:
Beam energies: (450.0, 450.0); (3500.0, 3500.0) GeV Run details:
A measurement of the underlying activity in scattering processes with a hard scale in the several-GeV region is performed in proton-proton collisions at Energies of 0.9 and 7 TeV, using data collected by the CMS experiment at the LHC. The production of charged particles with pseudorapidity |eta| < 2 and transverse momentum $p_\perp > 0.5$ GeV/$c$ is studied in the azimuthal region transverse to that of the leading set of charged particles forming a track-jet. Various comparisons are made between the two different energies and also beteen two sets of cuts on p_\perp for leading track jet p_\perp-leading $> 3$ GeV and pT-leading $> 20$ GeV. The activity is studied using 5 types of plots. Two profile plots for the multiplicity of charged particles and the scalar sum of p_\perp. and three distributions for the two previous quantities as well as p_\perp for all the particles in the transverse region. Beam energy must be specified (in GeV) as analysis option "ENERGY" when rivet-merging samples. Source code: CMS_2011_I916908.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
8namespace Rivet {
9
10
11 /// UE charged particles vs. leading jet
12 class CMS_2011_I916908 : public Analysis {
13 public:
14
15 RIVET_DEFAULT_ANALYSIS_CTOR(CMS_2011_I916908);
16
17
18 void init() {
19 const ChargedFinalState cfs((Cuts::etaIn(-2.0, 2.0) && Cuts::pT >= 500*MeV));
20 declare(cfs, "CFS");
21
22 const ChargedFinalState cfsforjet((Cuts::etaIn(-2.5, 2.5) && Cuts::pT >= 500*MeV));
23 const FastJets jetpro(cfsforjet, JetAlg::SISCONE, 0.5);
24 declare(jetpro, "Jets");
25
26 if (isCompatibleWithSqrtS(7000*GeV)) {
27 book(_h_Nch_vs_pT ,1, 1, 1); // Nch vs. pT_max
28 book(_h_Sum_vs_pT ,2, 1, 1); // sum(pT) vs. pT_max
29 book(_h_pT3_Nch ,5, 1, 1); // transverse Nch, pT_max > 3GeV
30 book(_h_pT3_Sum ,6, 1, 1); // transverse sum(pT), pT_max > 3GeV
31 book(_h_pT3_pT ,7, 1, 1); // transverse pT, pT_max > 3GeV
32 book(_h_pT20_Nch ,8, 1, 1); // transverse Nch, pT_max > 20GeV
33 book(_h_pT20_Sum ,9, 1, 1); // transverse sum(pT), pT_max > 20GeV
34 book(_h_pT20_pT ,10, 1, 1); // transverse pT, pT_max > 20GeV
35 }
36
37 if (isCompatibleWithSqrtS(900*GeV)) {
38 book(_h_Nch_vs_pT ,3, 1, 1); // Nch vs. pT_max
39 book(_h_Sum_vs_pT ,4, 1, 1); // sum(pT) vs. pT_max
40 book(_h_pT3_Nch ,11, 1, 1); // transverse Nch, pT_max > 3GeV
41 book(_h_pT3_Sum ,12, 1, 1); // transverse sum(pT), pT_max > 3GeV
42 book(_h_pT3_pT ,13, 1, 1); // transverse pT, pT_max > 3GeV
43 }
44
45 book(sumOfWeights3, "TMP/sumOfWeights3");
46 book(sumOfWeights20, "TMP/sumOfWeights20");
47 book(_nch_tot_pT3, "TMP/nch_tot_pT3");
48 book(_nch_tot_pT20, "TMP/nch_tot_pT20");
49 }
50
51
52 /// Perform the per-event analysis
53 void analyze(const Event& event) {
54
55 // Find the lead jet, applying a restriction that the jets must be within |eta| < 2.
56 FourMomentum p_lead;
57 for (const Jet& j : apply<FastJets>(event, "Jets").jetsByPt(Cuts::pT > 1.0*GeV && Cuts::abseta < 2.0)) {
58 p_lead = j.momentum();
59 break;
60 }
61 if (p_lead.isZero()) vetoEvent;
62 const double philead = p_lead.phi();
63 const double pTlead = p_lead.pT();
64
65 Particles particles = apply<ChargedFinalState>(event, "CFS").particlesByPt();
66
67 int nTransverse = 0;
68 double ptSumTransverse = 0.;
69 for (const Particle& p : particles) {
70 double dphi = fabs(deltaPhi(philead, p.phi()));
71 if (dphi>PI/3. && dphi<PI*2./3.) { // Transverse region
72 nTransverse++;
73
74 const double pT = p.pT()/GeV;
75 ptSumTransverse += pT;
76
77 if (pTlead > 3.0*GeV) _h_pT3_pT->fill(pT);
78 if (isCompatibleWithSqrtS(7000*GeV) && pTlead > 20.0*GeV) _h_pT20_pT->fill(pT);
79 }
80 }
81
82 const double area = 8./3. * PI;
83 _h_Nch_vs_pT->fill(pTlead/GeV, 1./area*nTransverse);
84 _h_Sum_vs_pT->fill(pTlead/GeV, 1./area*ptSumTransverse);
85 if(pTlead > 3.0*GeV) {
86 _h_pT3_Nch->fill(nTransverse);
87 _h_pT3_Sum->fill(ptSumTransverse);
88 sumOfWeights3->fill();
89 _nch_tot_pT3->fill(nTransverse);
90 }
91 if (isCompatibleWithSqrtS(7000.) && pTlead > 20.0*GeV) {
92 _h_pT20_Nch->fill(nTransverse);
93 _h_pT20_Sum->fill(ptSumTransverse);
94 sumOfWeights20->fill();
95 _nch_tot_pT20->fill(nTransverse);
96 }
97 }
98
99
100 /// Normalise histograms etc., after the run
101 void finalize() {
102 normalize(_h_pT3_Nch);
103 normalize(_h_pT3_Sum);
104 if (sumOfWeights3->val() != 0.0) normalize(_h_pT3_pT, *_nch_tot_pT3 / *sumOfWeights3);
105
106 if (isCompatibleWithSqrtS(7000*GeV)) {
107 normalize(_h_pT20_Nch);
108 normalize(_h_pT20_Sum);
109 if (sumOfWeights20->val() != 0.0) normalize(_h_pT20_pT, *_nch_tot_pT20 / *sumOfWeights20);
110 }
111 }
112
113
114
115 private:
116
117 /// @{
118 CounterPtr sumOfWeights3;
119 CounterPtr sumOfWeights20;
120
121 CounterPtr _nch_tot_pT3;
122 CounterPtr _nch_tot_pT20;
123
124 Profile1DPtr _h_Nch_vs_pT;
125 Profile1DPtr _h_Sum_vs_pT;
126 Histo1DPtr _h_pT3_Nch;
127 Histo1DPtr _h_pT3_Sum;
128 Histo1DPtr _h_pT3_pT;
129 Histo1DPtr _h_pT20_Nch;
130 Histo1DPtr _h_pT20_Sum;
131 Histo1DPtr _h_pT20_pT;
132 /// @}
133
134 };
135
136
137
138 RIVET_DECLARE_ALIASED_PLUGIN(CMS_2011_I916908, CMS_2011_S9120041);
139
140}
|