Rivet analyses referenceATLAS_2011_S8994773Calo-based underlying event at 900 GeV and 7 TeV in ATLASExperiment: ATLAS (LHC) Inspire ID: 891834 Status: VALIDATED Authors:
Beam energies: (450.0, 450.0); (3500.0, 3500.0) GeV Run details:
Underlying event measurements with the ATLAS detector at the LHC at center-of-mass energies of 900 GeV and 7 TeV, using calorimeter clusters rather than charged tracks. Beam energy must be specified (in GeV) as analysis option "ENERGY" when rivet-merging samples. Source code: ATLAS_2011_S8994773.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4
5namespace Rivet {
6
7
8 /// @brief Calo-based underlying event at 900 GeV and 7 TeV in ATLAS
9 ///
10 /// @author Jinlong Zhang
11 class ATLAS_2011_S8994773 : public Analysis {
12 public:
13
14 RIVET_DEFAULT_ANALYSIS_CTOR(ATLAS_2011_S8994773);
15
16
17 void init() {
18 const FinalState fs500((Cuts::etaIn(-2.5, 2.5) && Cuts::pT >= 500*MeV));
19 declare(fs500, "FS500");
20 const FinalState fslead((Cuts::etaIn(-2.5, 2.5) && Cuts::pT >= 1.0*GeV));
21 declare(fslead, "FSlead");
22
23 // Get an index for the beam energy
24 isqrts = -1;
25 if (isCompatibleWithSqrtS(900)) isqrts = 0;
26 else if (isCompatibleWithSqrtS( 7000)) isqrts = 1;
27 assert(isqrts >= 0);
28
29 // N profiles, 500 MeV pT cut
30 book(_hist_N_transverse_500 ,1+isqrts, 1, 1);
31 // pTsum profiles, 500 MeV pT cut
32 book(_hist_ptsum_transverse_500 ,3+isqrts, 1, 1);
33 // N vs. Delta(phi) profiles, 500 MeV pT cut
34 book(_hist_N_vs_dPhi_1_500 ,13+isqrts, 1, 1);
35 book(_hist_N_vs_dPhi_2_500 ,13+isqrts, 1, 2);
36 book(_hist_N_vs_dPhi_3_500 ,13+isqrts, 1, 3);
37 }
38
39
40 void analyze(const Event& event) {
41 // Require at least one cluster in the event with pT >= 1 GeV
42 const FinalState& fslead = apply<FinalState>(event, "FSlead");
43 if (fslead.size() < 1) {
44 vetoEvent;
45 }
46
47 // These are the particles with pT > 500 MeV
48 const FinalState& chargedNeutral500 = apply<FinalState>(event, "FS500");
49
50 // Identify leading object and its phi and pT
51 Particles particles500 = chargedNeutral500.particlesByPt();
52 Particle p_lead = particles500[0];
53 const double philead = p_lead.phi();
54 const double etalead = p_lead.eta();
55 const double pTlead = p_lead.pT();
56 MSG_DEBUG("Leading object: pT = " << pTlead << ", eta = " << etalead << ", phi = " << philead);
57
58 // Iterate over all > 500 MeV particles and count particles and scalar pTsum in the three regions
59 vector<double> num500(3, 0), ptSum500(3, 0.0);
60 // Temporary histos that bin N in dPhi.
61 // NB. Only one of each needed since binnings are the same for the energies and pT cuts
62 Histo1D hist_num_dphi_500(refData(13+isqrts,1,1));
63 for (const Particle& p : particles500) {
64 const double pT = p.pT();
65 const double dPhi = deltaPhi(philead, p.phi());
66 const int ir = region_index(dPhi);
67 num500[ir] += 1;
68 ptSum500[ir] += pT;
69
70 // Fill temp histos to bin N in dPhi
71 if (p.genParticle() != p_lead.genParticle()) { // We don't want to fill all those zeros from the leading track...
72 hist_num_dphi_500.fill(dPhi, 1);
73 }
74 }
75
76
77 // Now fill underlying event histograms
78 // The densities are calculated by dividing the UE properties by dEta*dPhi
79 // -- each region has a dPhi of 2*PI/3 and dEta is two times 2.5
80 const double dEtadPhi = (2*2.5 * 2*PI/3.0);
81 _hist_N_transverse_500->fill(pTlead/GeV, num500[1]/dEtadPhi);
82 _hist_ptsum_transverse_500->fill(pTlead/GeV, ptSum500[1]/GeV/dEtadPhi);
83
84 // Update the "proper" dphi profile histograms
85 // Note that we fill dN/dEtadPhi: dEta = 2*2.5, dPhi = 2*PI/nBins
86 // The values tabulated in the note are for an (undefined) signed Delta(phi) rather than
87 // |Delta(phi)| and so differ by a factor of 2: we have to actually norm for angular range = 2pi
88 const size_t nbins = refData(13+isqrts,1,1).numPoints();
89 for (size_t i = 0; i < nbins; ++i) {
90 double mean = hist_num_dphi_500.bin(i).xMid();
91 double value = 0.;
92 if (hist_num_dphi_500.bin(i).numEntries() > 0) {
93 mean = hist_num_dphi_500.bin(i).xMean();
94 value = hist_num_dphi_500.bin(i).area()/hist_num_dphi_500.bin(i).xWidth()/10.0;
95 }
96 if (pTlead/GeV >= 1.0) _hist_N_vs_dPhi_1_500->fill(mean, value);
97 if (pTlead/GeV >= 2.0) _hist_N_vs_dPhi_2_500->fill(mean, value);
98 if (pTlead/GeV >= 3.0) _hist_N_vs_dPhi_3_500->fill(mean, value);
99 }
100
101 }
102
103
104 private:
105
106 // Little helper function to identify Delta(phi) regions
107 inline int region_index(double dphi) {
108 assert(inRange(dphi, 0.0, PI, CLOSED, CLOSED));
109 if (dphi < PI/3.0) return 0;
110 if (dphi < 2*PI/3.0) return 1;
111 return 2;
112 }
113
114
115 int isqrts;
116
117 Profile1DPtr _hist_N_transverse_500;
118
119 Profile1DPtr _hist_ptsum_transverse_500;
120
121 Profile1DPtr _hist_N_vs_dPhi_1_500;
122 Profile1DPtr _hist_N_vs_dPhi_2_500;
123 Profile1DPtr _hist_N_vs_dPhi_3_500;
124
125 };
126
127
128
129 RIVET_DECLARE_ALIASED_PLUGIN(ATLAS_2011_S8994773, ATLAS_2011_I891834);
130
131}
|