rivet is hosted by Hepforge, IPPP Durham

Rivet analyses reference

ZEUS_2008_I810112

Measurement of $D^{\pm}$ and $D^0$ production in deep inelastic scattering using a lifetime tag at HERA
Experiment: ZEUS (HERA)
Inspire ID: 810112
Status: VALIDATED
Authors:
  • Andrii Verbytskyi
References:
  • Eur.Phys.J.C 63 (2009) 171-188
  • DOI:10.1140/epjc/s10052-009-1088-x
  • arXiv: 0812.3775
Beams: e+ p+, p+ e+, e- p+, p+ e-
Beam energies: ANY
Run details:
  • Inclusive DIS

The production of $D^{\pm}$ and $D^{0}$ mesons has been measured with the ZEUS detector at HERA using an integrated luminosity of $133.6 pb^{-1}$. The measurements cover the kinematic range $5 < Q^2 < 1000 GeV^2$, $0.02 < y < 0.7$, $1.5 < p_T^D < 15 GeV$ and $\eta^D < 1.6$. Combinatorial background to the D meson signals is reduced by using the ZEUS microvertex detector to reconstruct displaced secondary vertices. Production cross sections are compared with the predictions of next-to-leading-order QCD which is found to describe the data well. Measurements are extrapolated to the full kinematic phase space in order to obtain the open-charm contribution, $F_2^c\bar{c}$, to the proton structure function, $F_2$.

Source code: ZEUS_2008_I810112.cc
  1// -*- C++ -*-
  2#include "Rivet/Analysis.hh"
  3#include "Rivet/Projections/UnstableParticles.hh"
  4#include "Rivet/Projections/FinalState.hh"
  5#include "Rivet/Projections/FastJets.hh"
  6#include "Rivet/Projections/DISFinalState.hh"
  7
  8namespace Rivet {
  9
 10
 11  /// @brief Measurement of $D^{\pm}$ and $D^0$ production in deep inelastic scattering using a lifetime tag at HERA
 12  ///
 13  /// @author Andrii Verbytskyi
 14  class ZEUS_2008_I810112 : public Analysis {
 15  public:
 16
 17    /// Constructor
 18    RIVET_DEFAULT_ANALYSIS_CTOR(ZEUS_2008_I810112);
 19
 20    /// @name Analysis methods
 21    /// @{
 22
 23    /// Book histograms and initialise projections before the run
 24    void init() {
 25      /// Initialise and register projections
 26      FinalState fs;
 27      const DISKinematics& diskin = DISKinematics();
 28      declare(diskin,"Kinematics");
 29      declare(UnstableParticles(), "UPS");
 30
 31      book(_h_Dp_q2, 3, 1, 1);
 32      book(_h_Dp_x, 4, 1, 1);
 33      book(_h_Dp_pt, 5, 1, 1);
 34      book(_h_Dp_eta, 6, 1, 1);
 35
 36      book(_h_Dp_yinq2[0], "TMP/Dp0", refData( 11, 1, 1));
 37      book(_h_Dp_yinq2[1], "TMP/Dp1", refData( 11, 1, 2));
 38      book(_h_Dp_yinq2[2], "TMP/Dp2", refData( 11, 1, 3));
 39
 40      book(_s_Dp_yinq2[0], 11, 1, 1);
 41      book(_s_Dp_yinq2[1], 11, 1, 2);
 42      book(_s_Dp_yinq2[2], 11, 1, 3);
 43
 44      book(_h_D0_q2, 7, 1, 1);
 45      book(_h_D0_x, 8, 1, 1);
 46      book(_h_D0_pt, 9, 1, 1);
 47      book(_h_D0_eta, 10, 1, 1);
 48
 49      book(_h_D0_yinq2[0], "TMP/D00", refData( 12, 1, 1));
 50      book(_h_D0_yinq2[1], "TMP/D01", refData( 12, 1, 2));
 51      book(_h_D0_yinq2[2], "TMP/D02", refData( 12, 1, 3));
 52
 53      book(_s_D0_yinq2[0], 12, 1, 1);
 54      book(_s_D0_yinq2[1], 12, 1, 2);
 55      book(_s_D0_yinq2[2], 12, 1, 3);
 56    }
 57
 58
 59    /// A routine that selects the bin in kinematic space
 60    int _getbinQ2_OK(const DISKinematics& dk) {
 61      if (inRange(dk.Q2()/GeV2, 5.0, 9.0)) return 0;
 62      if (inRange(dk.Q2()/GeV2, 9.0, 44.0)) return 1;
 63      if (inRange(dk.Q2()/GeV2, 44.0, 1000.0)) return 2;
 64      return -1;
 65    }
 66
 67
 68    /// Perform the per-event analysis
 69    void analyze(const Event& event) {
 70      /// DIS kinematics
 71      const DISKinematics& dk = apply<DISKinematics>(event, "Kinematics");
 72      double q2  = dk.Q2();
 73      double x   = dk.x();
 74      double y   = dk.y();
 75
 76      ///Assure the event falls into the kinematic range of the measurement.
 77      if (!inRange(q2/GeV2, 5, 1000)) vetoEvent;
 78      if (!inRange(y, 0.02, 0.7)) vetoEvent;
 79
 80      /// Find out in which bin the measurement falls.
 81      int bin = _getbinQ2_OK(dk);
 82      if ( bin < 0 ) vetoEvent;
 83
 84      const UnstableParticles& ufs = apply<UnstableParticles>(event, "UPS");
 85
 86      /// Get \f$D^0$\f particles
 87      for (const Particle& p : select(ufs.particles(), Cuts::abspid == PID::D0)) {
 88        ///But not not from \f$D^{*\pm}$\f decays
 89        if (p.hasAncestorWith(Cuts::pid == PID::DSTARPLUS)) continue;
 90        if (p.hasAncestorWith(Cuts::pid == PID::DSTARMINUS)) continue;
 91        /// Select particles only in the \f$\eta-p_{T}$\f region
 92        if (!inRange(p.eta(), -1.6, 1.6)) continue;
 93        if (!inRange(p.pt()/GeV, 1.5, 15.0)) continue;
 94        /// Fill the general histograms
 95        _h_D0_pt->fill(p.pt()/GeV);
 96        _h_D0_eta->fill(p.eta());
 97        _h_D0_q2->fill(q2/GeV2);
 98        _h_D0_x->fill(x);
 99        /// Fill the cross-sections in selected kinematic bins
100        _h_D0_yinq2[bin]->fill(y);
101      }
102
103      /// Get \f$D^{\pm}$\f particles
104      for (const Particle& p : select(ufs.particles(), Cuts::abspid == PID::DPLUS)) {
105        /// Select particles only in the \f$\eta-p_{T}$\f region
106        if (!inRange(p.eta(), -1.6, 1.6)) continue;
107        if (!inRange(p.pt()/GeV, 1.5, 15.0)) continue;
108        /// Fill the general histograms
109        _h_Dp_pt->fill(p.pt()/GeV);
110        _h_Dp_eta->fill(p.eta());
111        _h_Dp_q2->fill(q2/GeV2);
112        _h_Dp_x->fill(x);
113        /// Fill the cross-sections in selected kinematic bins
114        _h_Dp_yinq2[bin]->fill(y);
115      }
116    }
117
118
119    /// Normalise histograms etc., after the run
120    void finalize() {
121      const double sf = crossSection()/nanobarn/sumOfWeights();
122      /// For CS comparison the MC should be corrected to difference between BR used for data and in MC
123      scale(_h_Dp_pt, sf);
124      scale(_h_Dp_eta, sf);
125      scale(_h_Dp_q2, sf);
126      scale(_h_Dp_x, sf);
127
128      for (size_t i = 0; i < 3; i++) {
129        scale(_h_Dp_yinq2[i], sf);
130        barchart(_h_Dp_yinq2[i], _s_Dp_yinq2[i]);
131      }
132
133      scale(_h_D0_pt, sf);
134      scale(_h_D0_eta, sf);
135      scale(_h_D0_q2, sf);
136      scale(_h_D0_x, sf);
137
138      for (size_t i = 0; i < 3; i++) {
139        scale(_h_D0_yinq2[i], sf);
140        barchart(_h_D0_yinq2[i], _s_D0_yinq2[i]);
141      }
142    }
143
144    /// @}
145
146
147  private:
148
149    Histo1DPtr _h_D0_pt, _h_D0_eta, _h_D0_q2, _h_D0_x, _h_D0_yinq2[3];
150    Histo1DPtr _h_Dp_pt, _h_Dp_eta, _h_Dp_q2, _h_Dp_x, _h_Dp_yinq2[3];
151
152    Estimate1DPtr _s_Dp_yinq2[3];
153    Estimate1DPtr _s_D0_yinq2[3];
154
155  };
156
157
158  RIVET_DECLARE_PLUGIN(ZEUS_2008_I810112);
159
160}