Rivet analyses referenceHERA_2015_I1377206Measurement of sigma_red (F2) of H1 and ZEUS at different beam energiesExperiment: H1 (HERA) Inspire ID: 1377206 Status: VALIDATED Authors:
Beam energies: ANY Run details:
A combination is presented of all inclusive deep inelastic cross sections previously published by the H1 and ZEUS collaborations at HERA for neutral and charged current ep scattering for zero beam polarisation. The data were taken at proton beam energies of 920, 820, 575 and 460 GeV and an electron beam energy of 27.5GeV. The data correspond to an integrated luminosity of about 1 fb^-1 and span six orders of magnitude in negative four-momentum-transfer squared, Q2, and Bjorken x . The correlations of the systematic uncertainties were evaluated and taken into account for the combination. The combined cross sections were input to QCD analyses at leading order, next-to-leading order and at next-to-next-to-leading order, providing a new set of parton distribution functions, called HERAPDF2.0.. The energy and beam options need to be specified when using rivet-merge. Source code: HERA_2015_I1377206.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/Beam.hh"
7
8namespace Rivet {
9
10
11 /// @brief Measurement of sigma_red (F2) of H1 and ZEUS at different beam energies
12 class HERA_2015_I1377206 : public Analysis {
13 public:
14
15 /// Constructor
16 RIVET_DEFAULT_ANALYSIS_CTOR(HERA_2015_I1377206);
17
18
19 /// @name Analysis methods
20 /// @{
21
22 /// Book histograms and initialise projections before the run
23 void init() {
24
25 // Initialise and register projections
26 declare(FinalState(Cuts::abseta < 5 && Cuts::pT > 100*MeV), "FS");
27 declare(DISLepton(), "Lepton");
28 declare(DISKinematics(), "Kinematics");
29
30 Histo1DPtr dummy;
31 string beamOpt = getOption<string>("BEAM","NONE");
32
33 if (beamOpt == "NONE") {
34 const ParticlePair& beam = beams();
35 _positron = (beam.first.pid() == PID::POSITRON || beam.second.pid() == PID::POSITRON );
36 } else {
37 if (beamOpt == "EMINUS") _positron = false;
38 else if (beamOpt == "EPLUS") _positron = true;
39 else {
40 MSG_ERROR("Beam option error. You have specified an unsupported beam.");
41 return;
42 }
43 }
44
45 // Book beams-dependent histograms
46 const double eps = 0.01;
47 if (isCompatibleWithSqrtS(318*GeV, eps) && _positron ) {
48 // NC e+ p at sqrts=318
49 const vector<double> Q2edges = {
50 0.1, 0.15, 0.2, 0.3, 0.4, 0.45, 0.6, 0.7, 1.0, 1.1, 1.3, 1.7, 2.3,
51 3.1, 3.8, 5.3, 8.0, 9.1, 11., 13., 17.4, 19.1, 25.8, 28., 30., 42.,
52 49., 54., 65., 75., 108., 134., 180., 225., 280., 325., 355., 455., 460.,
53 545., 560., 765., 770., 835., 900., 1120., 1295., 1300., 1755., 1800.,
54 2270., 2500., 3685., 4000., 6520., 7000., 9275., 10000., 15000., 17000.,
55 24770., 25000., 42000.
56 };
57 book(_h_sigred, Q2edges);
58 _h_sigred->maskBins({9, 24, 27, 36, 38, 40, 42, 44, 47, 49, 51, 53, 55, 57, 59, 61});
59 size_t idx = 0;
60 for (auto& b : _h_sigred->bins()) {
61 book(b, 1, 1, ++idx);
62 }
63 // CC e+ p at sqrts=318
64 book(_h_sigred_cc, {280., 325., 460., 545., 900., 1120., 1300., 1755., 1800., 2270.,
65 2500., 3685., 4000., 6520., 7000., 9275., 10000., 20000., 42000.});
66 _h_sigred_cc->maskBins({2, 4, 6, 8, 10, 12, 14, 16});
67 idx = 0;
68 for (auto& b : _h_sigred_cc->bins()) {
69 book(b, 6, 1, ++idx);
70 }
71 }
72 else if (isCompatibleWithSqrtS(300*GeV, eps) && _positron) {
73 // NC e+ p at sqrts=300
74 const vector<double> Q2edges = {
75 0.01, 0.05, 0.07, 0.09, 0.12, 0.18, 0.22, 0.32, 0.4, 0.45, 0.6, 0.7, 1.0, 1.1, 1.3,
76 1.7, 2.3, 3.1, 3.8, 5.3, 8., 9.1, 11., 13., 17.4, 19.1, 25.8, 28., 30., 42., 49.,
77 54., 65., 75., 108., 134., 180., 225., 280., 325., 355., 455., 460., 545., 560.,
78 765., 770., 835., 900., 1120., 1295., 1300., 1755., 1800., 2270., 2500., 3685.,
79 4000., 6520., 7000., 9275., 10000., 15000., 17000., 24770., 25000., 42000.
80 };
81 book(_h_sigred, Q2edges);
82 _h_sigred->maskBins({13, 28, 31, 40, 42, 44, 46, 48, 51, 53, 55, 57, 59, 61, 63, 65});
83 size_t idx = 0;
84 for (auto& b : _h_sigred->bins()) {
85 book(b, 2, 1, ++idx);
86 }
87 }
88 else if (isCompatibleWithSqrtS(251*GeV, eps) && _positron ) {
89 // NC e+ p at sqrts=251
90 const vector<double> Q2edges = {
91 1., 1.7, 2.3, 3.1, 3.8, 5.3, 8., 9.1, 11., 13., 17.4, 22.1, 28.,
92 30., 42., 49., 54., 65., 75., 108., 134., 180., 225., 280., 325.,
93 355., 455., 460., 545., 560., 765., 770., 835.
94 };
95 book(_h_sigred, Q2edges);
96 _h_sigred->maskBins({8, 13, 16, 18, 25, 29, 31});
97 size_t idx = 0;
98 for (auto& b : _h_sigred->bins()) {
99 book(b, 3, 1, ++idx);
100 }
101 }
102
103 // e+- p at sqrts=225
104 else if (isCompatibleWithSqrtS(225*GeV, eps)) {
105 if (_positron) {
106 // NC e+ p at sqrts=225
107 const vector<double> Q2edges = {
108 1., 1.7, 2.3, 3.1, 3.8, 5.3, 8., 9.1, 11., 13., 17.4, 22.1, 28.,
109 30., 42., 49., 54., 65., 75., 108., 134., 180., 225., 280., 325.,
110 355., 455., 460., 545., 560., 765., 770., 835.
111 };
112 book(_h_sigred, Q2edges);
113 _h_sigred->maskBins({8, 13, 16, 18, 25, 29, 31});
114 size_t idx = 0;
115 for (auto& b : _h_sigred->bins()) {
116 book(b, 4, 1, ++idx);
117 }
118 }
119 else {
120 // NC e- p at sqrts=225
121 const vector<double> Q2edges = {
122 54., 65., 75., 108., 134., 180., 225., 280., 325., 355., 455.,
123 460., 545., 560., 765., 770., 835., 900., 1120., 1295., 1300.,
124 1755., 1800., 2270., 2500., 3685., 4000., 6520., 7000., 9275.,
125 10000., 15000., 17000., 24770., 25000., 42000., 70000.
126 };
127 book(_h_sigred, Q2edges);
128 _h_sigred->maskBins({2, 9, 11, 13, 15, 17, 20, 22, 24, 26, 28, 30, 32, 34});
129 size_t idx = 0;
130 for (auto& b : _h_sigred->bins()) {
131 book(b, 5, 1, ++idx);
132 }
133 // CC e- p at sqrts=225
134 book(_h_sigred_cc, {280., 325., 460., 545., 900., 1120., 1300., 1755., 1800., 2270.,
135 2500., 3685., 4000., 6520., 7000., 9275., 10000., 20000., 42000.});
136 _h_sigred_cc->maskBins({2, 4, 6, 8, 10, 12, 14, 16});
137 idx = 0;
138 for (auto& b : _h_sigred_cc->bins()) {
139 book(b, 7, 1, ++idx);
140 }
141 }
142 }
143 }
144
145
146 void analyze(const Event& event) {
147
148 /// @todo Do the event by event analysis here
149 const DISKinematics& dk = apply<DISKinematics>(event, "Kinematics");
150 const DISLepton& dl = apply<DISLepton>(event,"Lepton");
151
152 // Get the DIS kinematics
153 double x = dk.x();
154 double y = dk.y();
155 double Q2 = dk.Q2()/GeV;
156
157 // Flux factor
158 const double alpha = 7.29927e-3;
159 // GF = 1.16638e-5 Fermi constant
160 const double GF2 = 1.16638e-5*1.16638e-5;
161 // MW = 80.385 W-boson mass
162 const double MW2 = 80.385 * 80.385;
163
164 if (PID::isNeutrino(dl.out().abspid()) ) {
165 // fill histo for CC
166 double F = 2.0*M_PI*x/GF2 * sqr((MW2 + Q2)/MW2);
167 _h_sigred_cc->fill(Q2,x,F); // fill histogram x,Q2
168 }
169 else {
170 // fill histo for NC
171 double F = x*sqr(Q2)/(2.0*M_PI*sqr(alpha)*(1.0+sqr(1-y)));
172 _h_sigred->fill(Q2,x,F); // fill histogram x,Q2
173 }
174 }
175
176
177 /// Normalise histograms etc., after the run
178 void finalize() {
179 const double gev2nb = 0.389e6;
180 const double scalefactor=crossSection()/nanobarn/sumOfWeights()/gev2nb ;
181 // with _h_sigred.scale also q2 bin width is scaled
182 scale(_h_sigred, scalefactor);
183 scale(_h_sigred_cc, scalefactor);
184 divByGroupWidth({_h_sigred, _h_sigred_cc});
185 }
186
187 /// @}
188
189
190 /// @name Histograms
191 /// @{
192 Histo1DGroupPtr _h_sigred, _h_sigred_cc;
193 Histo1DPtr _hist_Q2_10,_hist_Q2_100,_hist_Q2_1000,_hist_Q2_2000,_hist_Q2_3000;
194 bool _positron;
195 /// @}
196
197 };
198
199
200 RIVET_DECLARE_PLUGIN(HERA_2015_I1377206);
201
202}
|