rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

HERA_2015_I1377206

Measurement of sigma_red (F2) of H1 and ZEUS at different beam energies
Experiment: H1 (HERA)
Inspire ID: 1377206
Status: VALIDATED
Authors:
  • Hannes Jung (DESY)
  • Tymoteusz Strozniak
References: Beams: e+ p+, p+ e+, e- p+, p+ e-
Beam energies: ANY
Run details:
  • Cover full range in Bjorken-y

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
 31      Histo1DPtr dummy;
 32      string beamOpt = getOption<string>("BEAM","NONE");
 33
 34      if (beamOpt == "NONE") {
 35        const ParticlePair& beam = beams();
 36        if ( beam.first.pid() == PID::POSITRON || beam.second.pid() == PID::POSITRON ) {
 37          positron = true ;
 38        } else {
 39          positron = false ;
 40        }
 41      }
 42      else {
 43        if ( beamOpt == "EMINUS" )
 44          positron = false;
 45        else if ( beamOpt == "EPLUS" )
 46          positron = true;
 47        else {
 48          MSG_ERROR("Beam option error. You have specified an unsupported beam.");
 49          return;
 50        }
 51      }
 52
 53      double eps = 0.01 ;
 54      // NC e+ p at sqrts=318
 55      if (isCompatibleWithSqrtS(318., eps) && positron  ) {
 56        // cout << " NC e+ p sqrts = " << sqrtS() << endl ;
 57        _h_sigred.add( 0.1,     0.15, book(dummy,1,1,1)); // Q2=0.15
 58        _h_sigred.add( 0.15,     0.2, book(dummy,1,1,2)); // Q2=0.2
 59        _h_sigred.add( 0.2,      0.3, book(dummy,1,1,3)); // Q2=0.25
 60        _h_sigred.add( 0.3,      0.4, book(dummy,1,1,4)); // Q2=0.35
 61        _h_sigred.add( 0.4,      0.45, book(dummy,1,1,5)); // Q2=0.4
 62        _h_sigred.add( 0.45,     0.6, book(dummy,1,1,6)); // Q2=0.5
 63        _h_sigred.add( 0.6,      0.7, book(dummy,1,1,7)); // Q2=0.65
 64        _h_sigred.add( 0.7,      1.0, book(dummy,1,1,8)); // Q2=0.85
 65        _h_sigred.add( 1.1,      1.3, book(dummy,1,1,9)); // Q2=1.2
 66        _h_sigred.add( 1.3,      1.7, book(dummy,1,1,10)); // Q2=1.5
 67        _h_sigred.add( 1.7,      2.3, book(dummy,1,1,11)); // Q2=2
 68        _h_sigred.add( 2.3,      3.1, book(dummy,1,1,12)); // Q2=2.7
 69        _h_sigred.add( 3.1,      3.8, book(dummy,1,1,13)); // Q2=3.5
 70        _h_sigred.add( 3.8,      5.3, book(dummy,1,1,14)); // Q2=4.5
 71        _h_sigred.add( 5.3,      8.0, book(dummy,1,1,15)); // Q2=6.5
 72        _h_sigred.add( 8.0,      9.1, book(dummy,1,1,16)); // Q2=8.5
 73        _h_sigred.add( 9.1,      11., book(dummy,1,1,17)); // Q2=10
 74        _h_sigred.add( 11.,      13., book(dummy,1,1,18)); // Q2=12
 75        _h_sigred.add( 13.,     17.4, book(dummy,1,1,19)); // Q2=15
 76        _h_sigred.add( 17.4,    19.1, book(dummy,1,1,20)); // Q2=18
 77        _h_sigred.add( 19.1,    25.8, book(dummy,1,1,21)); // Q2=22
 78        _h_sigred.add( 25.8,     28., book(dummy,1,1,22)); // Q2=27
 79        _h_sigred.add( 30.,      42., book(dummy,1,1,23)); // Q2=35
 80        _h_sigred.add( 42.,      49., book(dummy,1,1,24)); // Q2=45
 81        _h_sigred.add( 54.,      65., book(dummy,1,1,25)); // Q2=60
 82        _h_sigred.add( 65.,      75., book(dummy,1,1,26)); // Q2=70
 83        _h_sigred.add( 75.,     108., book(dummy,1,1,27)); // Q2=90
 84        _h_sigred.add( 108.,    134., book(dummy,1,1,28)); // Q2=120
 85        _h_sigred.add( 134.,    180., book(dummy,1,1,29)); // Q2=150
 86        _h_sigred.add( 180.,    225., book(dummy,1,1,30)); // Q2=200
 87        _h_sigred.add( 225.,    280., book(dummy,1,1,31)); // Q2=250
 88        _h_sigred.add( 280.,    325., book(dummy,1,1,32)); // Q2=300
 89        _h_sigred.add( 355.,    455., book(dummy,1,1,33)); // Q2=400
 90        _h_sigred.add( 460.,    545., book(dummy,1,1,34)); // Q2=500
 91        _h_sigred.add( 560.,    765., book(dummy,1,1,35)); // Q2=650
 92        _h_sigred.add( 770.,    835., book(dummy,1,1,36)); // Q2=800
 93        _h_sigred.add( 900.,   1120., book(dummy,1,1,37)); // Q2=1000
 94        _h_sigred.add( 1120.,  1295., book(dummy,1,1,38)); // Q2=1200
 95        _h_sigred.add( 1300.,  1755., book(dummy,1,1,39)); // Q2=1500
 96        _h_sigred.add( 1800.,  2270., book(dummy,1,1,40)); // Q2=2000
 97        _h_sigred.add( 2500.,  3685., book(dummy,1,1,41)); // Q2=3000
 98        _h_sigred.add( 4000.,  6520., book(dummy,1,1,42)); // Q2=5000
 99        _h_sigred.add( 7000.,  9275., book(dummy,1,1,43)); // Q2=8000
100        _h_sigred.add( 10000.,15000., book(dummy,1,1,44)); // Q2=12000
101        _h_sigred.add( 17000.,24770., book(dummy,1,1,45)); // Q2=20000
102        _h_sigred.add( 25000.,42000., book(dummy,1,1,46)); // Q2=30000
103        // CC e+ p at sqrts=318
104        // cout << " CC e+ p sqrts = " << sqrtS() << endl ;
105        _h_sigred_cc.add( 280.,    325., book(dummy,6,1,1)); // Q2=300
106        _h_sigred_cc.add( 460.,    545., book(dummy,6,1,2)); // Q2=500
107        _h_sigred_cc.add( 900.,   1120., book(dummy,6,1,3)); // Q2=1000
108        _h_sigred_cc.add( 1300.,  1755., book(dummy,6,1,4)); // Q2=1500
109        _h_sigred_cc.add( 1800.,  2270., book(dummy,6,1,5)); // Q2=2000
110        _h_sigred_cc.add( 2500.,  3685., book(dummy,6,1,6)); // Q2=3000
111        _h_sigred_cc.add( 4000.,  6520., book(dummy,6,1,7)); // Q2=5000
112        _h_sigred_cc.add( 7000.,  9275., book(dummy,6,1,8)); // Q2=8000
113        _h_sigred_cc.add( 10000.,20000., book(dummy,6,1,9)); // Q2=15000
114        _h_sigred_cc.add( 20000.,42000., book(dummy,6,1,10)); // Q2=30000
115       }
116      // NC e+ p at sqrts=300
117      else if (isCompatibleWithSqrtS(300., eps) && positron  ) {
118        // cout << " NC e+ p sqrts = " << sqrtS() << endl ;
119         _h_sigred.add( 0.01,     0.05, book(dummy,2,1,1)); // Q2=0.045
120        _h_sigred.add( 0.05,     0.07, book(dummy,2,1,2)); // Q2=0.065
121        _h_sigred.add( 0.07,     0.09, book(dummy,2,1,3)); // Q2=0.085
122        _h_sigred.add( 0.09,     0.12, book(dummy,2,1,4)); // Q2=0.11
123        _h_sigred.add( 0.12,     0.18, book(dummy,2,1,5)); // Q2=0.15
124        _h_sigred.add( 0.18,     0.22, book(dummy,2,1,6)); // Q2=0.2
125        _h_sigred.add( 0.22,     0.32, book(dummy,2,1,7)); // Q2=0.25
126        _h_sigred.add( 0.32,     0.4,  book(dummy,2,1,8)); // Q2=0.35
127        _h_sigred.add( 0.4,      0.45, book(dummy,2,1,9)); // Q2=0.4
128        _h_sigred.add( 0.45,     0.6, book(dummy,2,1,10)); // Q2=0.5
129        _h_sigred.add( 0.6,      0.7, book(dummy,2,1,11)); // Q2=0.65
130        _h_sigred.add( 0.7,      1.0, book(dummy,2,1,12)); // Q2=0.85
131        _h_sigred.add( 1.1,      1.3, book(dummy,2,1,13)); // Q2=1.2
132        _h_sigred.add( 1.3,      1.7, book(dummy,2,1,14)); // Q2=1.5
133        _h_sigred.add( 1.7,      2.3, book(dummy,2,1,15)); // Q2=2
134        _h_sigred.add( 2.3,      3.1, book(dummy,2,1,16)); // Q2=2.7
135        _h_sigred.add( 3.1,      3.8, book(dummy,2,1,17)); // Q2=3.5
136        _h_sigred.add( 3.8,      5.3, book(dummy,2,1,18)); // Q2=4.5
137        _h_sigred.add( 5.3,      8.0, book(dummy,2,1,19)); // Q2=6.5
138        _h_sigred.add( 8.0,      9.1, book(dummy,2,1,20)); // Q2=8.5
139        _h_sigred.add( 9.1,      11., book(dummy,2,1,21)); // Q2=10
140        _h_sigred.add( 11.,      13., book(dummy,2,1,22)); // Q2=12
141        _h_sigred.add( 13.,     17.4, book(dummy,2,1,23)); // Q2=15
142        _h_sigred.add( 17.4,    19.1, book(dummy,2,1,24)); // Q2=18
143        _h_sigred.add( 19.1,    25.8, book(dummy,2,1,25)); // Q2=22
144        _h_sigred.add( 25.8,     28., book(dummy,2,1,26)); // Q2=27
145        _h_sigred.add( 30.,      42., book(dummy,2,1,27)); // Q2=35
146        _h_sigred.add( 42.,      49., book(dummy,2,1,28)); // Q2=45
147        _h_sigred.add( 54.,      65., book(dummy,2,1,29)); // Q2=60
148        _h_sigred.add( 65.,      75., book(dummy,2,1,30)); // Q2=70
149        _h_sigred.add( 75.,     108., book(dummy,2,1,31)); // Q2=90
150        _h_sigred.add( 108.,    134., book(dummy,2,1,32)); // Q2=120
151        _h_sigred.add( 134.,    180., book(dummy,2,1,33)); // Q2=150
152        _h_sigred.add( 180.,    225., book(dummy,2,1,34)); // Q2=200
153        _h_sigred.add( 225.,    280., book(dummy,2,1,35)); // Q2=250
154        _h_sigred.add( 280.,    325., book(dummy,2,1,36)); // Q2=300
155        _h_sigred.add( 355.,    455., book(dummy,2,1,37)); // Q2=400
156        _h_sigred.add( 460.,    545., book(dummy,2,1,38)); // Q2=500
157        _h_sigred.add( 560.,    765., book(dummy,2,1,39)); // Q2=650
158        _h_sigred.add( 770.,    835., book(dummy,2,1,40)); // Q2=800
159        _h_sigred.add( 900.,   1120., book(dummy,2,1,41)); // Q2=1000
160        _h_sigred.add( 1120.,  1295., book(dummy,2,1,42)); // Q2=1200
161        _h_sigred.add( 1300.,  1755., book(dummy,2,1,43)); // Q2=1500
162        _h_sigred.add( 1800.,  2270., book(dummy,2,1,44)); // Q2=2000
163        _h_sigred.add( 2500.,  3685., book(dummy,2,1,45)); // Q2=3000
164        _h_sigred.add( 4000.,  6520., book(dummy,2,1,46)); // Q2=5000
165        _h_sigred.add( 7000.,  9275., book(dummy,2,1,47)); // Q2=8000
166        _h_sigred.add( 10000.,15000., book(dummy,2,1,48)); // Q2=12000
167        _h_sigred.add( 17000.,24770., book(dummy,2,1,49)); // Q2=20000
168        _h_sigred.add( 25000.,42000., book(dummy,2,1,50)); // Q2=30000
169      }
170      else if (isCompatibleWithSqrtS(251., eps) && positron  ) {
171        // NC e+ p at sqrts=251
172        // cout << " NC e+ p sqrts = " << sqrtS() << endl ;
173        _h_sigred.add( 1.0,      1.7, book(dummy,3,1,1)); // Q2=1.5
174        _h_sigred.add( 1.7,      2.3, book(dummy,3,1,2)); // Q2=2
175        _h_sigred.add( 2.3,      3.1, book(dummy,3,1,3)); // Q2=2.5
176        _h_sigred.add( 3.1,      3.8, book(dummy,3,1,4)); // Q2=3.5
177        _h_sigred.add( 3.8,      5.3, book(dummy,3,1,5)); // Q2=5
178        _h_sigred.add( 5.3,      8.0, book(dummy,3,1,6)); // Q2=6.5
179        _h_sigred.add( 8.0,      9.1, book(dummy,3,1,7)); // Q2=8.5
180        _h_sigred.add( 11.,      13., book(dummy,3,1,8)); // Q2=12
181        _h_sigred.add( 13.,     17.4, book(dummy,3,1,9)); // Q2=15
182        _h_sigred.add( 17.4,    22.1, book(dummy,3,1,10)); // Q2=20
183        _h_sigred.add( 22.1,    28. , book(dummy,3,1,11)); // Q2=25
184        _h_sigred.add( 30.,      42., book(dummy,3,1,12)); // Q2=35
185        _h_sigred.add( 42.,      49., book(dummy,3,1,13)); // Q2=45
186        _h_sigred.add( 54.,      65., book(dummy,3,1,14)); // Q2=60
187        _h_sigred.add( 75.,     108., book(dummy,3,1,15)); // Q2=90
188        _h_sigred.add( 108.,    134., book(dummy,3,1,16)); // Q2=120
189        _h_sigred.add( 134.,    180., book(dummy,3,1,17)); // Q2=150
190        _h_sigred.add( 180.,    225., book(dummy,3,1,18)); // Q2=200
191        _h_sigred.add( 225.,    280., book(dummy,3,1,19)); // Q2=250
192        _h_sigred.add( 280.,    325., book(dummy,3,1,20)); // Q2=300
193        _h_sigred.add( 355.,    455., book(dummy,3,1,21)); // Q2=400
194        _h_sigred.add( 460.,    545., book(dummy,3,1,22)); // Q2=500
195        _h_sigred.add( 560.,    765., book(dummy,3,1,23)); // Q2=650
196        _h_sigred.add( 770.,    835., book(dummy,3,1,24)); // Q2=800
197      }
198      else if (isCompatibleWithSqrtS(225., eps) && positron  ) {
199        // NC e+ p at sqrts=225
200        // cout << " NC e+ p sqrts = " << sqrtS() << endl ;
201        _h_sigred.add( 1.0,      1.7, book(dummy,4,1,1)); // Q2=1.5
202        _h_sigred.add( 1.7,      2.3, book(dummy,4,1,2)); // Q2=2
203        _h_sigred.add( 2.3,      3.1, book(dummy,4,1,3)); // Q2=2.5
204        _h_sigred.add( 3.1,      3.8, book(dummy,4,1,4)); // Q2=3.5
205        _h_sigred.add( 3.8,      5.3, book(dummy,4,1,5)); // Q2=5
206        _h_sigred.add( 5.3,      8.0, book(dummy,4,1,6)); // Q2=6.5
207        _h_sigred.add( 8.0,      9.1, book(dummy,4,1,7)); // Q2=8.5
208        _h_sigred.add( 11.,      13., book(dummy,4,1,8)); // Q2=12
209        _h_sigred.add( 13.,     17.4, book(dummy,4,1,9)); // Q2=15
210        _h_sigred.add( 17.4,    22.1, book(dummy,4,1,10)); // Q2=20
211        _h_sigred.add( 22.1,    28. , book(dummy,4,1,11)); // Q2=25
212        _h_sigred.add( 30.,      42., book(dummy,4,1,12)); // Q2=35
213        _h_sigred.add( 42.,      49., book(dummy,4,1,13)); // Q2=45
214        _h_sigred.add( 54.,      65., book(dummy,4,1,14)); // Q2=60
215        _h_sigred.add( 75.,     108., book(dummy,4,1,15)); // Q2=90
216        _h_sigred.add( 108.,    134., book(dummy,4,1,16)); // Q2=120
217        _h_sigred.add( 134.,    180., book(dummy,4,1,17)); // Q2=150
218        _h_sigred.add( 180.,    225., book(dummy,4,1,18)); // Q2=200
219        _h_sigred.add( 225.,    280., book(dummy,4,1,19)); // Q2=250
220        _h_sigred.add( 280.,    325., book(dummy,4,1,20)); // Q2=300
221        _h_sigred.add( 355.,    455., book(dummy,4,1,21)); // Q2=400
222        _h_sigred.add( 460.,    545., book(dummy,4,1,22)); // Q2=500
223        _h_sigred.add( 560.,    765., book(dummy,4,1,23)); // Q2=650
224        _h_sigred.add( 770.,    835., book(dummy,4,1,24)); // Q2=800
225      }
226      else if (isCompatibleWithSqrtS(318., eps) && !positron  ) {
227        // NC e- p at sqrts=225
228        // cout << " NC e- p sqrts = " << sqrtS() << endl ;
229        _h_sigred.add( 54.,      65., book(dummy,5,1,1)); // Q2=60
230        _h_sigred.add( 75.,     108., book(dummy,5,1,2)); // Q2=90
231        _h_sigred.add( 108.,    134., book(dummy,5,1,3)); // Q2=120
232        _h_sigred.add( 134.,    180., book(dummy,5,1,4)); // Q2=150
233        _h_sigred.add( 180.,    225., book(dummy,5,1,5)); // Q2=200
234        _h_sigred.add( 225.,    280., book(dummy,5,1,6)); // Q2=250
235        _h_sigred.add( 280.,    325., book(dummy,5,1,7)); // Q2=300
236        _h_sigred.add( 355.,    455., book(dummy,5,1,8)); // Q2=400
237        _h_sigred.add( 460.,    545., book(dummy,5,1,9)); // Q2=500
238        _h_sigred.add( 560.,    765., book(dummy,5,1,10)); // Q2=650
239        _h_sigred.add( 770.,    835., book(dummy,5,1,11)); // Q2=800
240        _h_sigred.add( 900.,   1120., book(dummy,5,1,12)); // Q2=1000
241        _h_sigred.add( 1120.,  1295., book(dummy,5,1,13)); // Q2=1200
242        _h_sigred.add( 1300.,  1755., book(dummy,5,1,14)); // Q2=1500
243        _h_sigred.add( 1800.,  2270., book(dummy,5,1,15)); // Q2=2000
244        _h_sigred.add( 2500.,  3685., book(dummy,5,1,16)); // Q2=3000
245        _h_sigred.add( 4000.,  6520., book(dummy,5,1,17)); // Q2=5000
246        _h_sigred.add( 7000.,  9275., book(dummy,5,1,18)); // Q2=8000
247        _h_sigred.add( 10000.,15000., book(dummy,5,1,19)); // Q2=12000
248        _h_sigred.add( 17000.,24770., book(dummy,5,1,20)); // Q2=20000
249        _h_sigred.add( 25000.,42000., book(dummy,5,1,21)); // Q2=30000
250        _h_sigred.add( 42000.,70000., book(dummy,5,1,22)); // Q2=50000
251        // CC e- p at sqrts=318
252        // cout << " CC e- p sqrts = " << sqrtS() << endl ;
253        _h_sigred_cc.add( 280.,    325., book(dummy,7,1,1)); // Q2=300
254        _h_sigred_cc.add( 460.,    545., book(dummy,7,1,2)); // Q2=500
255        _h_sigred_cc.add( 900.,   1120., book(dummy,7,1,3)); // Q2=1000
256        _h_sigred_cc.add( 1300.,  1755., book(dummy,7,1,4)); // Q2=1500
257        _h_sigred_cc.add( 1800.,  2270., book(dummy,7,1,5)); // Q2=2000
258        _h_sigred_cc.add( 2500.,  3685., book(dummy,7,1,6)); // Q2=3000
259        _h_sigred_cc.add( 4000.,  6520., book(dummy,7,1,7)); // Q2=5000
260        _h_sigred_cc.add( 7000.,  9275., book(dummy,7,1,8)); // Q2=8000
261        _h_sigred_cc.add( 10000.,20000., book(dummy,7,1,9)); // Q2=15000
262        _h_sigred_cc.add( 20000.,42000., book(dummy,7,1,10)); // Q2=30000
263
264      }
265    }
266
267
268    void analyze(const Event& event) {
269
270      /// @todo Do the event by event analysis here
271      const DISKinematics& dk = applyProjection<DISKinematics>(event, "Kinematics");
272      const DISLepton& dl = applyProjection<DISLepton>(event,"Lepton");
273
274      // Get the DIS kinematics
275      double x  = dk.x();
276      double y = dk.y();
277      double Q2 = dk.Q2()/GeV;
278
279      // Flux factor
280      const double alpha = 7.29927e-3;
281      // GF = 1.16638e-5 Fermi constant
282      const double GF2 = 1.16638e-5*1.16638e-5;
283      // MW = 80.385 W-boson mass
284      const double MW2 = 80.385 * 80.385;
285
286      if (PID::isNeutrino(dl.out().abspid()) ) {
287        // fill histo for CC
288        double F = 2.0*M_PI*x/GF2 * pow((MW2 + Q2)/MW2,2);
289        _h_sigred_cc.fill(Q2,x,F); // fill histogram x,Q2
290      }
291      else {
292        // fill histo for NC
293        double F = x*pow(Q2,2.)/(2.0*M_PI*pow(alpha,2.)*(1.0+pow((1.-y),2.)));
294        _h_sigred.fill(Q2,x,F); // fill histogram x,Q2
295      }
296    }
297
298
299    /// Normalise histograms etc., after the run
300    void finalize() {
301      const double gev2nb =0.389e6;
302      const double scalefactor=crossSection()/nanobarn/sumOfWeights()/gev2nb ;
303      // with _h_sigred.scale also q2 bin width is scaled
304      _h_sigred.scale(scalefactor, this);
305      _h_sigred_cc.scale(scalefactor, this);
306    }
307
308    //@}
309
310
311    /// @name Histograms
312    //@{
313      BinnedHistogram _h_sigred, _h_sigred_cc;
314      Histo1DPtr _hist_Q2_10,_hist_Q2_100,_hist_Q2_1000,_hist_Q2_2000,_hist_Q2_3000;
315      bool positron ;
316    //@}
317
318
319  };
320
321
322  RIVET_DECLARE_PLUGIN(HERA_2015_I1377206);
323
324}