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 = applyProjection<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 Particles unstables ;
82 if (isDIS ) { unstables = apply<ParticleFinder>(event, "DStars").particles(Cuts::pT > 1.5*GeV && Cuts::absetaIn(0,1.5));}
83 else { unstables = apply<ParticleFinder>(event, "DStars").particles(Cuts::pT > 2.0*GeV && Cuts::absrapIn(0,1.5));}
84 const Particles dstars = filter_select(unstables, [](const Particle& p){return p.abspid() == PID::DSPLUS;});
85
86 if(dstars.empty() ) vetoEvent;
87 MSG_DEBUG("D*" << dstars.size());
88
89 const Particle& dstar = dstars.front();
90 // boosting the system
91 const LorentzTransform hcmboost = dk.boostHCM();
92 const FourMomentum hcmMom = hcmboost.transform(dstar.momentum());
93
94 //discriminate between dis and photoprod, and between ETA33 and ETA 44
95
96 //kinematics quantities
97 const double E = dstar.E();
98 const double p_z = dstar.pz();
99 //std::cout<<"y: "<<y<<endl;
100 if(y<0.02) vetoEvent;
101 const double m2 = 2.25; // charm mass^2
102 const double E_e = dk.beamLepton().E();
103 const double z = (E - p_z)/(2*y*E_e);
104
105
106 const double M2 = (1.44*hcmMom.pT2() + m2)/(z*(1-z));
107 const double x_g = (M2 + Q2)/(y*dk.s());
108 // std::cout<<"s: "<<dk.s()<<endl;
109 // std::cout<<"M2: "<<M2<<endl;
110 // std::cout<<"z: "<<z<<endl;
111
112 const double y_capp = dstar.rapidity();
113 const double W = sqrt(dk.W2());
114 //std::cout<<"x_g: "<<x_g<<endl;
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 }
135 if( W > 130 && W<230){
136 _h["1212"]->fill(log10(x_g));
137
138 }
139
140
141 }
142 if(ETAG44 == true && abs(y_capp) < 1.5 && dstar.pT()> 2) {
143
144 _h["pt88"]->fill(dstar.pT());
145 _h["rap88"]->fill(y_capp);
146 _h["1311"]->fill(log10(x_g));
147
148 }
149 }
150
151 /// Normalise histograms etc., after the run
152 void finalize() {
153 // conversion factors from ep to gamma p xsections (as given in publication)
154 const double F_etag33 = 0.0128;
155 const double F_etag44 = 0.0838;
156
157 double norm = crossSection()/nanobarn/sumW() ;
158 scale( _h["211"], norm);
159
160 scale(_h["311"], norm);
161 scale(_h["411"], norm);
162 scale(_h["511"], norm);
163 scale(_h["611"], norm);
164 double norm_mub = crossSection()/microbarn/sumW();
165
166 scale(_h["rap194"],norm_mub/F_etag33 );
167 scale(_h["pt194"], norm_mub/F_etag33 );
168 scale(_h["rap88"], norm_mub/F_etag44 );
169 scale(_h["pt88"], norm_mub/F_etag44 );
170
171 _hpt.scale( norm/F_etag33, this);
172 scale(_h["1211"], norm_mub/F_etag33);
173 scale(_h["1212"], norm_mub/F_etag33);
174 scale(_h["1311"], norm_mub/F_etag44);
175 }
176
177 ///@}
178
179
180 /// @name Histograms
181 ///@{
182 map<string, Histo1DPtr> _h;
183 map<string, Profile1DPtr> _p;
184 map<string, CounterPtr> _c;
185 BinnedHistogram _hpt;
186 ///@}
187
188
189 };
190
191
192 RIVET_DECLARE_PLUGIN(H1_1999_I481112);
193
194}
|