|
Rivet analyses reference
ZEUS_1997_I450085
Dijet photoproduction analysis
Experiment: ZEUS (HERA Run I)
Inspire ID: 450085
Status: VALIDATED
Authors:
- Juan Jose Juan Castella
- Jon Butterworth
- Matthew Wing
References:
- EEur. Phys. J C 1, (1998) 109–122
- DESY-97-196
- hep-ex/9710018
Beams: p+ e+
Beam energies: (820.0, 27.5) GeV
Run details:
- 820 GeV protons colliding with 27.5 GeV positrons; Direct and resolved photoproduction of dijet events; Both jets $E_T > 6$ GeV; Jet pseudorapidity $-1.375 < |\eta| < 1.875$; Pseudorapidity difference between jets $|\eta| < 0.5$;
Dijet cross sections are presented using photoproduction data obtained with the ZEUS detector during 1994. These measurements represent an extension of previous results, as the higher statistics allow cross sections to be measured at higher jet transverse energy (ETJ). Jets are identified in the hadronic final state, and the cross sections compared to complete next-to-leading order QCD calculations. Agreement with these calculations is seen for the pseudorapidity dependence of the direct photon events with ETJ $>$ 6 GeV and of the resolved photon events with ETJ $>$ 11 GeV. Calculated cross sections for resolved photon processes with 6 GeV $<$ ETJ $<$ 11 GeV lie below the data. The paper also studied reconstructing the jets with two different types of legacy cone algorithm (EUCELL and PUCELL). Histograms for these are not generated by this rivet routine as focus was on the the KTCLUS data.
Source code:
ZEUS_1997_I450085.cc
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119 | // -*- C++ -*-
#include "Rivet/Analysis.hh"
#include "Rivet/Projections/Beam.hh"
#include "Rivet/Projections/DISKinematics.hh"
#include "Rivet/Projections/FastJets.hh"
namespace Rivet {
// @brief ZEUS dijet photoproduction in different xgamma regions
class ZEUS_1997_I450085 : public Analysis {
public:
// Constructor
RIVET_DEFAULT_ANALYSIS_CTOR(ZEUS_1997_I450085);
// Initialization, books projections and histograms
void init() {
// Projections
// checking recombination scheme and radius checked with original code from M.Wing
FinalState fs;
declare(FastJets(fs, fastjet::JetAlgorithm::kt_algorithm, fastjet::RecombinationScheme::Et_scheme, 1.0), "Jets");
declare(DISKinematics(), "Kinematics");
//Histograms
// Table 1
book(_h_etabar_all[0], 1, 1, 1);
book(_h_etabar_all[1], 2, 1, 1);
book(_h_etabar_all[2], 3, 1, 1);
book(_h_etabar_all[3], 4, 1, 1);
// Table 6
book(_h_etabar[1][0], 21, 1, 1);
book(_h_etabar[1][1], 22, 1, 1);
book(_h_etabar[1][2], 23, 1, 1);
book(_h_etabar[1][3], 24, 1, 1);
// Table 7
book(_h_etabar[0][0], 25, 1, 1);
book(_h_etabar[0][1], 26, 1, 1);
book(_h_etabar[0][2], 27, 1, 1);
book(_h_etabar[0][3], 28, 1, 1);
}
// Analysis
void analyze(const Event & event) {
// Determine kinematics, including event orientation since ZEUS coord system is for +z = proton direction
const DISKinematics & kin = apply<DISKinematics>(event, "Kinematics");
if (kin.failed()) vetoEvent;
const int orientation = kin.orientation();
// Q2 and inelasticity cuts
if (kin.Q2() > 4 * GeV2) vetoEvent;
if (!inRange(kin.y(), 0.2, 0.8)) vetoEvent;
// Jet calculation
const Jets jets = apply<FastJets>(event, "Jets") \
.jets(Cuts::Et > 6 * GeV && Cuts::etaIn(-1.375 * orientation, 1.875 * orientation), cmpMomByEt);
MSG_DEBUG("Jet multiplicity = " << jets.size());
//Dijet event selection
if (jets.size() < 2) vetoEvent;
const Jet & j1 = jets[0];
const Jet & j2 = jets[1];
//Jet eta, average eta and eta difference calculation
const double eta1 = orientation * j1.eta(), eta2 = orientation * j2.eta();
const double etabar = (eta1 + eta2) / 2;
const double etadiff = eta1 - eta2;
//Cut in pseudorapidity difference
if (abs(etadiff) > 0.5) vetoEvent;
// Calculation of x_gamma^obs
/// note Assuming Ee is the lab frame positron momentum, not in proton rest frame cf.
const double xyobs = (j1.Et() * exp(-eta1) + j2.Et() * exp(-eta2)) / (2 * kin.y() * kin.beamLepton().E());
const size_t i_xyobs = (xyobs < 0.75) ? 0 : 1;
//Classify events according to minimum jet energies
size_t iE_min = 0;
if (j1.Et() > 8 * GeV && j2.Et() > 8 * GeV)
iE_min = 1;
if (j1.Et() > 11 * GeV && j2.Et() > 11 * GeV)
iE_min = 2;
if (j1.Et() > 15 * GeV && j2.Et() > 15 * GeV)
iE_min = 3;
//Fill histograms
for (size_t isel = 0; isel <= iE_min; ++isel) {
//T1
_h_etabar_all[isel]->fill(etabar);
//T6, T7
if (xyobs < 0.3) vetoEvent;
_h_etabar[i_xyobs][isel]->fill(etabar);
}
}
// Finalize
void finalize() {
const double sf = crossSection() / nanobarn / sumOfWeights();
for (size_t ix = 0; ix < 2; ++ix) {
for (auto& h : _h_etabar[ix]) scale(h, sf);
}
for (auto& h : _h_etabar_all) scale(h, sf);
}
private:
//name Histograms
Histo1DPtr _h_etabar_all[4], _h_etabar[2][4];
};
RIVET_DECLARE_PLUGIN(ZEUS_1997_I450085);
}
|
|