Rivet analyses referenceH1_1999_I481112Measurement of D* meson cross-sections at HERA and determination of the gluon density in the proton using NLO QCDExperiment: H1 (HERA) Inspire ID: 481112 Status: VALIDATED Authors:
Beam energies: (27.5, 820.0); (820.0, 27.5) GeV
With the H1 detector at the $ep$ collider HERA, $D^*$ meson production cross sections have been measured in deep inelastic scattering with four-momentum transfers $Q^2>2$ GeV$^2$ and in photoproduction at energies around $W_{\gamma p} \sim 88$ GeV and 194 GeV. Next-to-Leading Order QCD calculations are found to describe the differential cross sections within theoretical and experimental uncertainties. Using these calculations, the NLO gluon momentum distribution in the proton, $x_g g(x_g)$, has been extracted in the momentum fraction range $7.5 \times 10^{-4}< x_g <4 \times 10^{-2}$ at average scales $mu^2 =25$ to 50 GeV$^2$. The gluon momentum fraction $x_g$ has been obtained from the measured kinematics of the scattered electron and the $D^*$ meson in the final state. Source code: H1_1999_I481112.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/UnstableParticles.hh"
7
8namespace Rivet {
9
10
11 /// @brief Measurement of D* meson cross-sections at HERA
12 class H1_1999_I481112 : public Analysis {
13 public:
14
15 /// Constructor
16 RIVET_DEFAULT_ANALYSIS_CTOR(H1_1999_I481112);
17
18
19 /// @name Analysis methods
20 ///@{
21
22
23 /// Book histograms and initialise projections before the run
24 void init() {
25
26 // Initialise and register projections
27
28 // The basic final-state projection:
29 // all final-state particles within
30 // the given eta acceptance
31 const FinalState fs(Cuts::abseta < 1.5);
32 declare(fs, "fs");
33 // The final-state particles declared above are clustered using FastJet with
34
35 //Initialize quantities needed for cuts
36 declare(DISKinematics(), "Kinematics");
37 declare(UnstableParticles(), "DStars");
38
39
40 Histo1DPtr dummy;
41
42 book(_h["211"],2,1,1);
43 book(_h["311"],3,1,1);
44 book(_h["411"],4,1,1);
45 book(_h["511"],5,1,1);
46 book(_h["611"],6,1,1);
47 book(_h["rap194"], 7,1,1);
48 book(_h["pt194"], 8,1,1);
49 book(_h["rap88"], 9,1,1);
50 book(_h["pt88"], 10,1,1);
51 _hpt.add(2.5, 3.5, book(dummy, 11,1,1));
52 _hpt.add(3.5, 5.0, book(dummy, 11,1,2));
53 _hpt.add(5.0, 10.0, book(dummy, 11,1,3));
54 book(_h["1211"],12,1,1);
55 book(_h["1212"],12,1,2);
56 book(_h["1311"],13,1,1);
57 }
58
59
60 /// Perform the per-event analysis
61 void analyze(const Event& event) {
62
63
64 const DISKinematics& dk = apply<DISKinematics>(event, "Kinematics");
65
66
67 bool isDIS = false;
68 bool ETAG44 = false;
69 bool ETAG33 = false;
70
71 const double y = dk.y();
72 const double Q2 = dk.Q2();
73
74 if (Q2 > 2 && Q2 <100 && y > 0.05 && y < 0.7) isDIS = true;
75 if (Q2 < 0.009 && y > 0.02 && y <0.32 ) ETAG44 = true;
76 if (Q2 < 0.01 && y > 0.29 && y <0.62 ) ETAG33 = true;
77 if (isDIS == false && ETAG44 == false && ETAG33 == false) vetoEvent;
78
79
80 //Creating array of D*
81 Cut cuts = isDIS? (Cuts::pT > 1.5*GeV && Cuts::abseta < 1.5) : (Cuts::pT > 2*GeV && Cuts::absrap < 1.5);
82 Particles unstables = apply<ParticleFinder>(event, "DStars").particles(cuts);
83 const Particles dstars = filter_select(unstables, [](const Particle& p){
84 return p.abspid() == PID::DSPLUS;
85 });
86
87 if(dstars.empty() ) vetoEvent;
88 MSG_DEBUG("D*" << dstars.size());
89
90 const Particle& dstar = dstars.front();
91 // boosting the system
92 const LorentzTransform hcmboost = dk.boostHCM();
93 const FourMomentum hcmMom = hcmboost.transform(dstar.momentum());
94
95 //discriminate between dis and photoprod, and between ETA33 and ETA 44
96
97 //kinematics quantities
98 const double E = dstar.E();
99 const double p_z = dstar.pz();
100 //std::cout<<"y: "<<y<<endl;
101 if (y<0.02) vetoEvent;
102 const double m2 = 2.25; // charm mass^2
103 const double E_e = dk.beamLepton().E();
104 const double z = (E - p_z)/(2*y*E_e);
105 if (z > 1) {
106 MSG_DEBUG("Momentum fraction greater than unity! This should not happen. Vetoing event.");
107 vetoEvent;
108 }
109
110 const double M2 = (1.44*hcmMom.pT2() + m2)/(z*(1-z));
111 const double x_g = (M2 + Q2)/(y*dk.s());
112
113 const double y_capp = dstar.rapidity();
114 const double W = sqrt(dk.W2());
115
116 //perform the cuts
117 if (isDIS == true){
118 _h["211"]->fill(dstar.pT());
119 _h["411"]->fill(dstar.eta());
120 _h["511"]->fill(Q2);
121 _h["611"]->fill(log10(x_g));
122 //boosting to the hcm frame
123 _h["311"]->fill(hcmMom.pT());
124 }
125 if (ETAG33 == true && abs(y_capp) < 1.5 && dstar.pT() > 2.5*GeV ) {
126
127 _h["rap194"]->fill(y_capp);
128 _h["pt194"] ->fill(dstar.pT());
129 _hpt.fill(dstar.pT(), y_capp);
130
131 if ( W > 173 && W<273) {
132 _h["1211"]->fill(log10(x_g));
133 }
134 if ( W > 130 && W<230) {
135 _h["1212"]->fill(log10(x_g));
136 }
137 }
138 if (ETAG44 == true && abs(y_capp) < 1.5 && dstar.pT()> 2) {
139
140 _h["pt88"]->fill(dstar.pT());
141 _h["rap88"]->fill(y_capp);
142 _h["1311"]->fill(log10(x_g));
143
144 }
145 }
146
147 /// Normalise histograms etc., after the run
148 void finalize() {
149 // conversion factors from ep to gamma p xsections (as given in publication)
150 const double F_etag33 = 0.0128;
151 const double F_etag44 = 0.0838;
152
153 double norm = crossSection()/nanobarn/sumW() ;
154 scale( _h["211"], norm);
155
156 scale(_h["311"], norm);
157 scale(_h["411"], norm);
158 scale(_h["511"], norm);
159 scale(_h["611"], norm);
160 double norm_mub = crossSection()/microbarn/sumW();
161
162 scale(_h["rap194"],norm_mub/F_etag33 );
163 scale(_h["pt194"], norm_mub/F_etag33 );
164 scale(_h["rap88"], norm_mub/F_etag44 );
165 scale(_h["pt88"], norm_mub/F_etag44 );
166
167 _hpt.scale( norm/F_etag33, this);
168 scale(_h["1211"], norm_mub/F_etag33);
169 scale(_h["1212"], norm_mub/F_etag33);
170 scale(_h["1311"], norm_mub/F_etag44);
171 }
172
173 ///@}
174
175
176 /// @name Histograms
177 ///@{
178 map<string, Histo1DPtr> _h;
179 BinnedHistogram _hpt;
180 ///@}
181
182
183 };
184
185
186 RIVET_DECLARE_PLUGIN(H1_1999_I481112);
187
188}
|