Rivet analyses referenceZEUS_2010_I875006ZEUS dijet cross-sections in neutral-current DISExperiment: ZEUS (HERA) Inspire ID: 875006 Status: VALIDATED Authors:
Beam energies: ANY Run details:
Single- and double-differential inclusive dijet cross sections in neutral current deep inelastic ep scattering have been measured with the ZEUS detector using an integrated luminosity of 374 pb$^{-1}$. The measurement was performed at large values of the photon virtuality, Q$^2$, between 125 and 20000 GeV$^2$. The jets were reconstructed with the k$_T$ cluster algorithm in the Breit reference frame and selected by requiring their transverse energies in the Breit frame, E$_T^{B,jet}$, to be larger than 8 GeV. In addition, the invariant mass of the dijet system, M$_{jj}$, was required to be greater than 20 GeV. The cross sections are described by the predictions of next-to-leading-order QCD. Source code: ZEUS_2010_I875006.cc 1// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/FinalState.hh"
4#include "Rivet/Projections/FastJets.hh"
5#include "Rivet/Projections/DISKinematics.hh"
6#include "Rivet/Projections/DISFinalState.hh"
7#include <cmath>
8
9namespace Rivet {
10
11
12 /// @brief DIS dijets in the Breit frame
13 class ZEUS_2010_I875006 : public Analysis {
14 public:
15
16 /// Constructor
17 RIVET_DEFAULT_ANALYSIS_CTOR(ZEUS_2010_I875006);
18
19
20 /// @name Analysis methods
21 /// @{
22
23 /// Book histograms and initialise projections before the run
24 void init() {
25
26 // Initialise and register projections
27 declare(DISKinematics(), "Kinematics");
28
29 // All final state particles boosted to Breit frame then clustered
30 //using FastJet KT algorithm with jet radius parameter 1
31 const DISFinalState DISfs(DISFrame::BREIT);
32 FastJets DISjetfs(DISfs, JetAlg::KT, 1.0);
33 declare(DISjetfs, "DISjets");
34
35 // Book histograms
36 // specify custom binning
37
38 //Non-Grouped histgrams
39 book(_h_Q2, 1,1,1);
40 book(_h_XBj, 2,1,1);
41 book(_h_Et, 3,1,1);
42 book(_h_Mjj, 4,1,1);
43 book(_h_Eta, 5,1,1);
44 book(_h_Zeta, 6,1,1);
45
46 //Zeta values seperated into Q2 ranges
47 book(_h_ZetaQ2[0], 7,1,1);
48 book(_h_ZetaQ2[1], 8,1,1);
49 book(_h_ZetaQ2[2], 9,1,1);
50 book(_h_ZetaQ2[3], 10,1,1);
51 book(_h_ZetaQ2[4], 11,1,1);
52 book(_h_ZetaQ2[5], 12,1,1);
53
54 //Transverse jet energy seperated into Q2 ranges
55 book(_h_EtQ2[0], 13,1,1);
56 book(_h_EtQ2[1], 14,1,1);
57 book(_h_EtQ2[2], 15,1,1);
58 book(_h_EtQ2[3], 16,1,1);
59 book(_h_EtQ2[4], 17,1,1);
60 book(_h_EtQ2[5], 18,1,1);
61 }
62
63 /// Perform the per-event analysis
64 void analyze(const Event& event) {
65
66 //First Lorentz invariant quantities in Lab frame
67 DISKinematics dis = apply<DISKinematics>(event, "Kinematics");
68 double Q2 = dis.Q2();
69 double xbj = dis.x();
70 double y = dis.y();
71
72 //Perform required cut on Q2 and y
73 if (!inRange(Q2, 125*GeV2, 20000*GeV2)) vetoEvent;
74 if (!inRange(y, 0.2, 0.6)) vetoEvent;
75
76 //Get Lorentz transforms for Breit Boost and Lab Boost
77 const LorentzTransform breitboost = dis.boostBreit();
78 const LorentzTransform labboost = breitboost.inverse();
79
80 //Get jets clustered in Breit frame
81 Jets jets = apply<FastJets>(event, "DISjets").jetsByPt();
82
83 //Boost jets to lab frame
84 for(std::vector<int>::size_type i=0; i<jets.size(); i++){
85 jets[i].transformBy(labboost);
86 }
87
88 //Cut on Pseudorapidity in lab frame
89 const int orientation = dis.orientation();
90 vector<Jet> cutJets;
91 for(std::vector<int>::size_type i=0; i<jets.size(); i++){
92 double etaJet = jets[i].eta()*orientation;
93 if(etaJet < 2.5 && etaJet > -1){
94 cutJets.push_back(jets[i]);
95 }
96 }
97
98 //veto event if only single jet
99 if(cutJets.size()<2){
100 vetoEvent;
101 }
102
103 //Boost jets to Breit frame
104 for(std::vector<int>::size_type i=0; i<cutJets.size(); i++){
105 cutJets[i].transformBy(breitboost);
106 }
107
108 //Sort jets by et in descending order
109 std::sort(cutJets.begin(), cutJets.end(),[](const Jet& j1, const Jet& j2){
110 return j1.Et()>j2.Et();
111 });
112
113 //Ensure two hardest jets have Et>8GeV in Breit frame
114 const Jet& jet1 = cutJets[0];
115 const Jet& jet2 = cutJets[1];
116
117 if(jet1.Et()<8*GeV || jet2.Et()<8*GeV){
118 vetoEvent;
119 }
120
121 //Extract required quantities in Breit frame
122 //Dijet mean transverse energy
123 const double dijetEt = (jet1.Et() + jet2.Et())/2;
124
125 //Invariant dijet mass of hardest transverse jets > 20GeV
126 const double Mjj = FourMomentum(jet1.momentum() + jet2.momentum()).mass();
127 if(Mjj<20*GeV){
128 vetoEvent;
129 }
130
131 const double eta1 = orientation*jet1.eta();
132 const double eta2 = orientation*jet2.eta();
133 const double etastar = abs(eta1 - eta2)/2;
134
135 const double logZeta = log10(xbj*(1 + pow(Mjj,2)/Q2));
136
137 //Fill histograms
138 _h_Q2->fill(Q2);
139 _h_XBj->fill(xbj);
140 _h_Et->fill(dijetEt);
141 _h_Mjj->fill(Mjj);
142 _h_Eta->fill(etastar);
143 _h_Zeta->fill(logZeta);
144
145 //Fill histograms for different Q2 ranges
146 if(Q2>125*GeV2 && Q2<=250*GeV2){
147 _h_ZetaQ2[0]->fill(logZeta);
148 _h_EtQ2[0]->fill(dijetEt);
149 }else if (Q2>250*GeV2 && Q2<=500*GeV2){
150 _h_ZetaQ2[1]->fill(logZeta);
151 _h_EtQ2[1]->fill(dijetEt);
152 }else if (Q2>500*GeV2 && Q2<=1000*GeV2){
153 _h_ZetaQ2[2]->fill(logZeta);
154 _h_EtQ2[2]->fill(dijetEt);
155 }else if (Q2>1000*GeV2 && Q2<=2000*GeV2){
156 _h_ZetaQ2[3]->fill(logZeta);
157 _h_EtQ2[3]->fill(dijetEt);
158 }else if (Q2>2000*GeV2 && Q2<=5000*GeV2){
159 _h_ZetaQ2[4]->fill(logZeta);
160 _h_EtQ2[4]->fill(dijetEt);
161 }else if (Q2>5000*GeV2 && Q2<=20000*GeV2){
162 _h_ZetaQ2[5]->fill(logZeta);
163 _h_EtQ2[5]->fill(dijetEt);
164 }
165 }
166
167 /// Normalise histograms etc., after the run
168 void finalize() {
169 //Calculate scaling factor from cross section
170 const double sf = crossSection()/picobarn/sumW(); //Scale factor with cuts
171
172 scale(_h_Q2, sf);
173 scale(_h_XBj, sf);
174 scale(_h_Et, sf);
175 scale(_h_Mjj, sf);
176 scale(_h_Eta, sf);
177 scale(_h_Zeta, sf);
178
179 for(int i = 0; i<6;i++){
180 scale(_h_ZetaQ2[i], sf);
181 scale(_h_EtQ2[i],sf);
182 }
183 }
184
185 /// @}
186
187
188 /// @name Histograms
189 /// @{
190 Histo1DPtr _h_Q2;
191 Histo1DPtr _h_XBj;
192 Histo1DPtr _h_Et;
193 Histo1DPtr _h_Mjj;
194 Histo1DPtr _h_Eta;
195 Histo1DPtr _h_Zeta;
196 Histo1DPtr _h_ZetaQ2[6];
197 Histo1DPtr _h_EtQ2[6];
198 /// @}
199
200 };
201
202
203 RIVET_DECLARE_PLUGIN(ZEUS_2010_I875006);
204
205}
|