Rivet analyses referenceZEUS_1997_I450085Dijet photoproduction analysisExperiment: ZEUS (HERA Run I) Inspire ID: 450085 Status: VALIDATED Authors:
Beam energies: (820.0, 27.5) GeV Run details:
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// -*- C++ -*-
2#include "Rivet/Analysis.hh"
3#include "Rivet/Projections/Beam.hh"
4#include "Rivet/Projections/DISKinematics.hh"
5#include "Rivet/Projections/FastJets.hh"
6
7namespace Rivet {
8 // @brief ZEUS dijet photoproduction in different xgamma regions
9
10 class ZEUS_1997_I450085 : public Analysis {
11 public:
12
13 // Constructor
14 RIVET_DEFAULT_ANALYSIS_CTOR(ZEUS_1997_I450085);
15
16 // Initialization, books projections and histograms
17 void init() {
18
19 // Projections
20 // checking recombination scheme and radius checked with original code from M.Wing
21 FinalState fs;
22 declare(FastJets(fs, fastjet::JetAlgorithm::kt_algorithm, fastjet::RecombinationScheme::Et_scheme, 1.0), "Jets");
23 declare(DISKinematics(), "Kinematics");
24
25 //Histograms
26 // Table 1
27 book(_h_etabar_all[0], 1, 1, 1);
28 book(_h_etabar_all[1], 2, 1, 1);
29 book(_h_etabar_all[2], 3, 1, 1);
30 book(_h_etabar_all[3], 4, 1, 1);
31 // Table 6
32 book(_h_etabar[1][0], 21, 1, 1);
33 book(_h_etabar[1][1], 22, 1, 1);
34 book(_h_etabar[1][2], 23, 1, 1);
35 book(_h_etabar[1][3], 24, 1, 1);
36 // Table 7
37 book(_h_etabar[0][0], 25, 1, 1);
38 book(_h_etabar[0][1], 26, 1, 1);
39 book(_h_etabar[0][2], 27, 1, 1);
40 book(_h_etabar[0][3], 28, 1, 1);
41
42 }
43 // Analysis
44 void analyze(const Event & event) {
45
46 // Determine kinematics, including event orientation since ZEUS coord system is for +z = proton direction
47 const DISKinematics & kin = apply<DISKinematics>(event, "Kinematics");
48 if (kin.failed()) vetoEvent;
49 const int orientation = kin.orientation();
50
51 // Q2 and inelasticity cuts
52 if (kin.Q2() > 4 * GeV2) vetoEvent;
53 if (!inRange(kin.y(), 0.2, 0.8)) vetoEvent;
54
55 // Jet calculation
56 const Jets jets = apply<FastJets>(event, "Jets") \
57 .jets(Cuts::Et > 6 * GeV && Cuts::etaIn(-1.375 * orientation, 1.875 * orientation), cmpMomByEt);
58 MSG_DEBUG("Jet multiplicity = " << jets.size());
59 //Dijet event selection
60 if (jets.size() < 2) vetoEvent;
61 const Jet & j1 = jets[0];
62 const Jet & j2 = jets[1];
63
64 //Jet eta, average eta and eta difference calculation
65 const double eta1 = orientation * j1.eta(), eta2 = orientation * j2.eta();
66 const double etabar = (eta1 + eta2) / 2;
67 const double etadiff = eta1 - eta2;
68
69 //Cut in pseudorapidity difference
70 if (abs(etadiff) > 0.5) vetoEvent;
71
72 // Calculation of x_gamma^obs
73 /// note Assuming Ee is the lab frame positron momentum, not in proton rest frame cf.
74 const double xyobs = (j1.Et() * exp(-eta1) + j2.Et() * exp(-eta2)) / (2 * kin.y() * kin.beamLepton().E());
75 const size_t i_xyobs = (xyobs < 0.75) ? 0 : 1;
76
77 //Classify events according to minimum jet energies
78 size_t iE_min = 0;
79 if (j1.Et() > 8 * GeV && j2.Et() > 8 * GeV)
80 iE_min = 1;
81 if (j1.Et() > 11 * GeV && j2.Et() > 11 * GeV)
82 iE_min = 2;
83 if (j1.Et() > 15 * GeV && j2.Et() > 15 * GeV)
84 iE_min = 3;
85
86 //Fill histograms
87 for (size_t isel = 0; isel <= iE_min; ++isel) {
88
89 //T1
90 _h_etabar_all[isel]->fill(etabar);
91 //T6, T7
92 if (xyobs < 0.3) vetoEvent;
93 _h_etabar[i_xyobs][isel]->fill(etabar);
94 }
95
96
97 }
98
99 // Finalize
100 void finalize() {
101 const double sf = crossSection() / nanobarn / sumOfWeights();
102 for (size_t ix = 0; ix < 2; ++ix) {
103 for (auto& h : _h_etabar[ix]) scale(h, sf);
104
105 }
106 for (auto& h : _h_etabar_all) scale(h, sf);
107
108 }
109
110
111 private:
112 //name Histograms
113 Histo1DPtr _h_etabar_all[4], _h_etabar[2][4];
114
115 };
116
117 RIVET_DECLARE_PLUGIN(ZEUS_1997_I450085);
118
119}
|