rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ZEUS_1996_I420332

Measurement of the F2 structure function in deep inelastic e+ p scattering using 1994 data from the ZEUS detector at HERA
Experiment: ZEUS (HERA)
Inspire ID: 420332
Status: VALIDATED
Authors:
  • Hannes Jung
References:
  • Z.Phys.C 72 (1996) 399, DOI:10.1007/s002880050260,10.1007/BF02909169, - hep-ex/9607002
Beams: e+ p+, p+ e+
Beam energies: (27.5, 820.0); (820.0, 27.5) GeV
    No run details listed

We present measurements of the structure function $F_2$ in ep scattering at HERA in the range $ 3.5 < Q^2 < 5000 $ GeV$^2$. A new reconstruction method has allowed a significant improvement in the resolution of the kinematic variables and an extension of the kinematic region covered by the experiment. At $Q^2 < 35 $ GeV$^2$ the range in $x$ now spans $ 6.3 \cdot 10^{-5} < x < 0.08$ providing overlap with measurements from fixed target experiments. At values of $Q^2$ above 1000 GeV$^2$ the $x$ range extends to 0.5.

Source code: ZEUS_1996_I420332.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
  7namespace Rivet {
  8
  9
 10  /// @brief   F2 structure function in DIS e+ p (ZEUS)
 11  
 12  class ZEUS_1996_I420332 : public Analysis {
 13  public:
 14
 15    /// Constructor
 16    RIVET_DEFAULT_ANALYSIS_CTOR(ZEUS_1996_I420332);
 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      // Initialise and register projections
 27      declare(DISLepton(), "Lepton");
 28      declare(DISKinematics(), "Kinematics");
 29		
 30	
 31       Histo1DPtr dummy;
 32	 _h_f2.add( 3.2,    4.0,book(dummy,1,1,1)); 
 33	 _h_f2.add( 4.0,    5.0,book(dummy,2,1,1));
 34	 _h_f2.add( 5.0,    7.0,book(dummy,3,1,1));
 35	 _h_f2.add( 7.0,    9.0,book(dummy,5,1,1));
 36	 _h_f2.add( 9.0,    11.0,book(dummy,6,1,1));
 37	 _h_f2.add( 11.0,   13.0,book(dummy,8,1,1));
 38	 _h_f2.add( 13.0,   16.0,book(dummy,9,1,1));
 39	 _h_f2.add( 16.0,   20.0,book(dummy,11,1,1));
 40	 _h_f2.add( 20.,    32,  book(dummy,13,1,1));
 41	 _h_f2.add( 32.,    40., book(dummy,15,1,1));
 42	 _h_f2.add( 40.,    50., book(dummy,17,1,1));
 43	 _h_f2.add( 50.,    65., book(dummy,18,1,1));
 44	 _h_f2.add( 65.,    85., book(dummy,19,1,1));
 45	 _h_f2.add( 85.,    110.,book(dummy,20,1,1));
 46	 _h_f2.add( 110.,   140.,book(dummy,21,1,1));
 47	 _h_f2.add( 140.,   185.,book(dummy,22,1,1));
 48	 _h_f2.add( 185.,   240.,book(dummy,23,1,1));
 49	 _h_f2.add( 240.,   310.,book(dummy,24,1,1));
 50	 _h_f2.add( 310.,   410.,book(dummy,25,1,1));
 51	 _h_f2.add( 410.,   530.,book(dummy,26,1,1));
 52	 _h_f2.add( 530.,   710.,book(dummy,27,1,1));
 53	 _h_f2.add( 710.,   900.,book(dummy,28,1,1));
 54	 _h_f2.add( 900.,  1300.,book(dummy,29,1,1));
 55	 _h_f2.add( 1300., 1800.,book(dummy,30,1,1));
 56	 _h_f2.add( 1800., 2500.,book(dummy,31,1,1));
 57	 _h_f2.add( 2500., 3500.,book(dummy,32,1,1));
 58	 _h_f2.add( 3500., 15000.,book(dummy,33,1,1));
 59
 60/*
 61        book(_hist_Q2_10, "Q2_10",100,1., 11.0);
 62        book(_hist_Q2_100, "Q2_100",100,10., 100.0);
 63        book(_hist_Q2_1000, "Q2_1000",100,100., 1000.0);
 64        book(_hist_Q2_2000, "Q2_2000",100,800., 5000.0);
 65        book(_hist_Q2_3000, "Q2_3000",100,3000., 10000.0);
 66*/
 67
 68    }
 69
 70
 71    /// Perform the per-event analysis
 72    void analyze(const Event& event) {
 73
 74      const DISKinematics& dk = applyProjection<DISKinematics>(event, "Kinematics");
 75      //const DISLepton& dl = applyProjection<DISLepton>(event,"Lepton");
 76
 77      // Get the DIS kinematics
 78      double x  = dk.x();
 79      double y = dk.y();
 80      double Q2 = dk.Q2()/GeV;
 81	
 82	// Flux factor
 83	const double alpha = 7.29927e-3;
 84	double F = x*pow(Q2,2.)/(2.0*M_PI*pow(alpha,2.)*(1.0+pow((1.-y),2.)));
 85/*
 86	_hist_Q2_10-> fill(Q2) ;
 87	_hist_Q2_100-> fill(Q2) ;
 88	_hist_Q2_1000-> fill(Q2) ;
 89	_hist_Q2_2000-> fill(Q2) ;
 90	_hist_Q2_3000-> fill(Q2) ;
 91*/
 92	_h_f2.fill(Q2,x,F); // wypelniamy histogram x w skali Q2
 93
 94
 95    }
 96
 97
 98    /// Normalise histograms etc., after the run
 99    void finalize() {
100	double gev2nb =0.389e6;
101      double scalefactor=crossSection()/nanobarn/sumOfWeights()/gev2nb ;
102      // with _h_f2.scale also q2 bin width is scaled
103      _h_f2.scale(scalefactor, this);
104
105    }
106
107    ///@}
108
109
110    /// @name Histograms
111    ///@{
112	BinnedHistogram _h_f2;
113      Histo1DPtr _hist_Q2_10,_hist_Q2_100,_hist_Q2_1000,_hist_Q2_2000,_hist_Q2_3000;
114    ///@}
115
116
117  };
118
119
120  RIVET_DECLARE_PLUGIN(ZEUS_1996_I420332);
121
122}