Rivet analyses referenceATLAS_2012_I1091481Azimuthal ordering of charged hadronsExperiment: ATLAS (LHC) Inspire ID: 1091481 Status: VALIDATED Authors:
Beam energies: (450.0, 450.0); (3500.0, 3500.0) GeV Run details:
Measurement of the ordering of charged hadrons in the azimuthal angle relative to the beam axis at the Large Hadron Collider (LHC). A spectral analysis of correlations between longitudinal and transverse components of the momentum of the charged hadrons is performed. Data were recorded with the ATLAS detector at centre-of-mass energies of $\sqrt{s} = 900 \text{GeV}$ and $\sqrt{s} = 7 \text{TeV}$. The correlations measured in a phase space region dominated by low-$p_\perp$ particles are not well described by conventional models of hadron production. The measured spectra show features consistent with the fragmentation of a QCD string represented by a helix-like ordered gluon chain. Source code: ATLAS_2012_I1091481.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/ChargedFinalState.hh"
4
5namespace Rivet {
6
7
8 /// @brief Azimuthal ordering of charged hadrons
9 class ATLAS_2012_I1091481 : public Analysis {
10 public:
11
12 /// Constructor
13 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2012_I1091481);
14
15 /// Book histograms and initialise projections before the run
16 void init() {
17
18 ChargedFinalState cfs100(Cuts::abseta < 2.5 && Cuts::pT > 0.1*GeV);
19 declare(cfs100,"CFS100");
20 ChargedFinalState cfs500(Cuts::abseta < 2.5 && Cuts::pT > 0.5*GeV);
21 declare(cfs500,"CFS500");
22
23 for (double eVal : allowedEnergies()) {
24 const string en = toString(int(eVal));
25 if (isCompatibleWithSqrtS(eVal)) _sqs = en;
26 size_t ih = bool(en == "900") + 1;
27
28 book(_h[en+"E_10_100"], ih, 1, 1);
29 book(_h[en+"E_1_100"], ih, 1, 2);
30 book(_h[en+"E_10_500"], ih, 1, 3);
31
32 book(_h[en+"eta_10_100"], ih, 2, 1);
33 book(_h[en+"eta_1_100"], ih, 2, 2);
34 book(_h[en+"eta_10_500"], ih, 2, 3);
35
36 book(_c[en+"inclusive"], "_norm_inclusive_"+en);
37 book(_c[en+"lowPt"], "_norm_lowPt_"+en);
38 book(_c[en+"pt500"], "_norm_pt500_"+en);
39 }
40 if (_sqs == "" && !merging()) {
41 throw BeamError("Invalid beam energy for " + name() + "\n");
42 }
43 }
44
45
46 // Recalculate particle energy assuming pion mass
47 double getPionEnergy(const Particle& p) {
48 double m_pi = 0.1396*GeV;
49 double p2 = p.p3().mod2()/(GeV*GeV);
50 return sqrt(sqr(m_pi) + p2);
51 }
52
53
54 // S_eta core for one event
55 //
56 // -1 + 1/Nch * |sum_j^Nch exp[i*(xi eta_j - Phi_j)]|^2
57 //
58 double getSeta(const Particles& part, double xi) {
59 std::complex<double> c_eta (0.0, 0.0);
60 for (const Particle& p : part) {
61 double eta = p.eta();
62 double phi = p.phi();
63 double arg = xi*eta-phi;
64 std::complex<double> temp(cos(arg), sin(arg));
65 c_eta += temp;
66 }
67 return std::norm(c_eta)/part.size() - 1.0;
68 }
69
70
71 // S_E core for one event
72 //
73 // -1 + 1/Nch * |sum_j^Nch exp[i*(omega X_j - Phi_j)]|^2
74 //
75 double getSE(const Particles& part, double omega) {
76 double Xj = 0.0;
77 std::complex<double> c_E (0.0, 0.0);
78 for (unsigned int i=0; i < part.size(); ++i) {
79 Xj += 0.5*getPionEnergy(part[i]);
80 double phi = part[i].phi();
81 double arg = omega*Xj - phi;
82 std::complex<double> temp(cos(arg), sin(arg));
83 c_E += temp;
84 Xj += 0.5*getPionEnergy(part[i]);
85 }
86 return std::norm(c_E)/part.size() - 1.0;
87 }
88
89
90 // Convenient fill function
91 void fillS(Histo1DPtr h, const Particles& part, bool SE=true) {
92 // Loop over bins, take bin centers as parameter values
93 // @todo use barchart()
94 for (auto& b : h->bins()) {
95 double x = b.xMid();
96 double width = b.xWidth();
97 double y;
98 if (SE) y = getSE(part, x);
99 else y = getSeta(part, x);
100 h->fill(x, y * width);
101 }
102 }
103
104
105 /// Perform the per-event analysis
106 void analyze(const Event& event) {
107 // Charged fs
108 const ChargedFinalState& cfs100 = apply<ChargedFinalState>(event, "CFS100");
109 const Particles part100 = cfs100.particles(cmpMomByEta);
110 const ChargedFinalState& cfs500 = apply<ChargedFinalState>(event, "CFS500");
111 const Particles& part500 = cfs500.particles(cmpMomByEta);
112
113 // Veto event if the most inclusive phase space has less than 10 particles and the max pT is > 10 GeV
114 if (part100.size() < 11) vetoEvent;
115 double ptmax = cfs100.particlesByPt()[0].pT()/GeV;
116 if (ptmax > 10.0) vetoEvent;
117
118 // Fill the pt>100, pTmax<10 GeV histos
119 fillS(_h[_sqs+"E_10_100"], part100, true);
120 fillS(_h[_sqs+"eta_10_100"], part100, false);
121 _c[_sqs+"inclusive"]->fill();
122
123 // Fill the pt>100, pTmax<1 GeV histos
124 if (ptmax < 1.0) {
125 fillS(_h[_sqs+"E_1_100"], part100, true);
126 fillS(_h[_sqs+"eta_1_100"], part100, false);
127 _c[_sqs+"lowPt"]->fill();
128 }
129
130 // Fill the pt>500, pTmax<10 GeV histos
131 if (part500.size() > 10) {
132 fillS(_h[_sqs+"E_10_500"], part500, true );
133 fillS(_h[_sqs+"eta_10_500"], part500, false);
134 _c[_sqs+"pt500"]->fill();
135 }
136 }
137
138 /// Normalise histograms etc., after the run
139 void finalize() {
140 // The scaling takes the multiple fills per event into account
141 for (double eVal : allowedEnergies()) {
142 const string en = toString(int(eVal));
143
144 scale(_h[en+"E_10_100"], 1.0/ *_c[en+"inclusive"]);
145 scale(_h[en+ "E_1_100"], 1.0/ *_c[en+"lowPt"]);
146 scale(_h[en+"E_10_500"], 1.0/ *_c[en+"pt500"]);
147
148 scale(_h[en+"eta_10_100"], 1.0/ *_c[en+"inclusive"]);
149 scale(_h[en+ "eta_1_100"], 1.0/ *_c[en+"lowPt"]);
150 scale(_h[en+"eta_10_500"], 1.0/ *_c[en+"pt500"]);
151 }
152 }
153
154
155 private:
156
157 map<string,Histo1DPtr> _h;
158 map<string,CounterPtr> _c;
159
160 string _sqs = "";
161 };
162
163
164 RIVET_DECLARE_PLUGIN(ATLAS_2012_I1091481);
165
166}
|