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 book(_hpt, {2.5, 3.5, 5.0, 10.0}, {"d11-x01-y01", "d11-x01-y02", "d11-x01-y03"});
52 book(_h["1211"], 12,1,1);
53 book(_h["1212"], 12,1,2);
54 book(_h["1311"], 13,1,1);
55 }
56
57
58 /// Perform the per-event analysis
59 void analyze(const Event& event) {
60
61
62 const DISKinematics& dk = apply<DISKinematics>(event, "Kinematics");
63
64
65 bool isDIS = false;
66 bool ETAG44 = false;
67 bool ETAG33 = false;
68
69 const double y = dk.y();
70 const double Q2 = dk.Q2();
71
72 if (Q2 > 2 && Q2 <100 && y > 0.05 && y < 0.7) isDIS = true;
73 if (Q2 < 0.009 && y > 0.02 && y <0.32 ) ETAG44 = true;
74 if (Q2 < 0.01 && y > 0.29 && y <0.62 ) ETAG33 = true;
75 if (isDIS == false && ETAG44 == false && ETAG33 == false) vetoEvent;
76
77
78 //Creating array of D*
79 Cut cuts = isDIS? (Cuts::pT > 1.5*GeV && Cuts::abseta < 1.5) : (Cuts::pT > 2*GeV && Cuts::absrap < 1.5);
80 Particles unstables = apply<ParticleFinder>(event, "DStars").particles(cuts);
81 const Particles dstars = select(unstables, [](const Particle& p){
82 return p.abspid() == PID::DSPLUS;
83 });
84
85 if(dstars.empty() ) vetoEvent;
86 MSG_DEBUG("D*" << dstars.size());
87
88 const Particle& dstar = dstars.front();
89 // boosting the system
90 const LorentzTransform hcmboost = dk.boostHCM();
91 const FourMomentum hcmMom = hcmboost.transform(dstar.momentum());
92
93 //discriminate between dis and photoprod, and between ETA33 and ETA 44
94
95 //kinematics quantities
96 const double E = dstar.E();
97 const double p_z = dstar.pz();
98 //std::cout<<"y: "<<y<<endl;
99 if (y<0.02) vetoEvent;
100 const double m2 = 2.25; // charm mass^2
101 const double E_e = dk.beamLepton().E();
102 const double z = (E - p_z)/(2*y*E_e);
103 if (z > 1) {
104 MSG_DEBUG("Momentum fraction greater than unity! This should not happen. Vetoing event.");
105 vetoEvent;
106 }
107
108 const double M2 = (1.44*hcmMom.pT2() + m2)/(z*(1-z));
109 const double x_g = (M2 + Q2)/(y*dk.s());
110
111 const double y_capp = dstar.rapidity();
112 const double W = sqrt(dk.W2());
113
114 //perform the cuts
115 if (isDIS == true){
116 _h["211"]->fill(dstar.pT());
117 _h["411"]->fill(dstar.eta());
118 _h["511"]->fill(Q2);
119 _h["611"]->fill(log10(x_g));
120 //boosting to the hcm frame
121 _h["311"]->fill(hcmMom.pT());
122 }
123 if (ETAG33 == true && abs(y_capp) < 1.5 && dstar.pT() > 2.5*GeV ) {
124
125 _h["rap194"]->fill(y_capp);
126 _h["pt194"] ->fill(dstar.pT());
127 _hpt->fill(dstar.pT(), y_capp);
128
129 if ( W > 173 && W<273) {
130 _h["1211"]->fill(log10(x_g));
131 }
132 if ( W > 130 && W<230) {
133 _h["1212"]->fill(log10(x_g));
134 }
135 }
136 if (ETAG44 == true && abs(y_capp) < 1.5 && dstar.pT()> 2) {
137
138 _h["pt88"]->fill(dstar.pT());
139 _h["rap88"]->fill(y_capp);
140 _h["1311"]->fill(log10(x_g));
141
142 }
143 }
144
145 /// Normalise histograms etc., after the run
146 void finalize() {
147 // conversion factors from ep to gamma p xsections (as given in publication)
148 const double F_etag33 = 0.0128;
149 const double F_etag44 = 0.0838;
150
151 double norm = crossSection()/nanobarn/sumW() ;
152 scale( _h["211"], norm);
153
154 scale(_h["311"], norm);
155 scale(_h["411"], norm);
156 scale(_h["511"], norm);
157 scale(_h["611"], norm);
158 double norm_mub = crossSection()/microbarn/sumW();
159
160 scale(_h["rap194"],norm_mub/F_etag33 );
161 scale(_h["pt194"], norm_mub/F_etag33 );
162 scale(_h["rap88"], norm_mub/F_etag44 );
163 scale(_h["pt88"], norm_mub/F_etag44 );
164
165 scale(_hpt, norm/F_etag33);
166 scale(_h["1211"], norm_mub/F_etag33);
167 scale(_h["1212"], norm_mub/F_etag33);
168 scale(_h["1311"], norm_mub/F_etag44);
169 }
170
171 ///@}
172
173
174 /// @name Histograms
175 ///@{
176 map<string, Histo1DPtr> _h;
177 Histo1DGroupPtr _hpt;
178 ///@}
179
180
181 };
182
183
184 RIVET_DECLARE_PLUGIN(H1_1999_I481112);
185
186}
|