Charged Particle Momentum Spectra in e+e- annihilation at sqrt(s) = 192-209 GeV

Charged particle momentum distributions are studied in the reaction e+e- ->hadrons, using data collected with the OPAL detector at centre-of-mass energies from 192 GeV to 209 GeV. The data correspond to an average centre-of- mass energy of 201.7 GeV and a total integrated luminosity of 433 pb-1. The measured distributions and derived quantities, in combination with corresponding results obtained at lower centre-of-mass energies, are compared to QCD predictions in various theoretical approaches to study the energy dependence of the strong interaction and to test QCD as the theory describing it. In general, a good agreement is found between the measurements and the corresponding QCD predictions.


Introduction
Charged particle momentum distributions are studied in the hadronic decays of virtual photons or Z bosons (hereafter referred to as (Z/γ) * decays) at centre-of-mass (c.m.) energies, √ s, between 192 and 209 GeV, the highest available e + e − energies.These distributions can be measured with relatively high precision.Combined with lower energy data they provide a powerful test of quantum chromodynamics (QCD) as the theory describing the strong interaction.
Previous studies using e + e − annihilation data at c.m. energies from 130 up to 189 GeV have shown that QCD based models and calculations give a good description of most of the measured distributions [1,2,3,4,5,6].With the yet higher energy data used in this paper, we obtain even more stringent tests of the models and theory.In addition we use these highest energy data in combination with earlier measurements at lower c.m. energies to study the energy dependence of the strong interaction.We measure the charged particle momentum distribution, 1/σ • dσ ch /dx p , and the ξ p distribution, 1/σ • dσ ch /dξ p , where x p = 2p/ √ s, ξ p = ln(1/x p ), p is the particle momentum, σ is the cross-section for non-radiative hadronic (Z/γ) * decays and σ ch is the charged particle crosssection in these events.We also measure the distribution of the rapidity, y = 1 2 ln E+p E−p , where p is the momentum component parallel to the thrust axis and E is the energy of the particle, and the distributions of the 3-momentum components in, p in ⊥ , and perpendicular to, p out ⊥ , the event plane.This plane is defined by the eigenvectors associated with the two largest eigenvalues of the momentum tensor, as in reference [1].
We present comparisons of the fractional momentum distributions 1/σ • dσ ch /dx p and 1/σ • dσ ch /dξ p in the range √ s = 14−202 GeV to QCD predictions in different theoretical approaches.
As in earlier publications [1,2,3] we also study the evolution with the c.m. energy of the peak position ξ 0 in the ξ p distribution.
In Section 2 a brief description of the OPAL detector is given.The samples of data and simulated events used in the analysis are described in Section 3. In Section 4 the event selection and analysis procedure are described.In Section 5 we present the observed distributions and, together with those obtained at lower energies, compare them to the predictions of different Monte Carlo programs and analytic QCD predictions.A summary and conclusions are given in Section 6.

2
The OPAL detector The OPAL detector was operated at the LEP e + e − collider at CERN.A detailed description can be found in [7].The analysis presented here relies mainly on the measurement of momenta and directions of charged particles in the tracking chambers and of energy deposited in the electromagnetic calorimeters of the detector.
All tracking systems were located inside a solenoidal magnet which provided a uniform axial magnetic field of 0.435 T along the beam axis 1 .The magnet was surrounded by a lead glass electromagnetic calorimeter and a hadron calorimeter of the sampling type.Outside the hadron calorimeter, the detector was surrounded by a system of muon chambers.There were similar layers of detectors in the forward and backward endcaps.
The tracking system consisted of a silicon microvertex detector, an inner vertex chamber, a large volume jet chamber, and specialised chambers at the outer radius of the jet chamber to improve the measurements in the z direction.The main tracking detector was the central jet chamber.This device was approximately 4 m long and had an outer radius of about 1.85 m.It had 24 sectors with radial planes of 159 sense wires spaced by 1 cm.The momenta p of tracks in the r-φ plane were measured with a precision σ p /p = 0.02 2 + (0.0015 • p[GeV/c]) 2 .
The electromagnetic calorimeters in the barrel and the endcap sections of the detector consisted of 11704 lead glass blocks with a depth of 24.6 radiation lengths in the barrel and typically 22 radiation lengths in the endcaps.The barrel section covered the angular region | cos θ| < 0.82 and the endcap section the region 0.82 < | cos θ| < 0.98.

Data and Monte Carlo samples
The data used in this measurement were recorded in 1999 and 2000, the final two years of the LEP 2 program.During this time e + e − interactions were collected at c.m. energies in the range 192-209 GeV.The total integrated luminosity of the data sample, as evaluated using small angle Bhabha collisions, is 433 pb −1 .The luminosity collected at the different c.m. energies is detailed in Figure 1 and in Table 1.
Monte Carlo event samples were generated at c.m. energies of 192,196,200,202,206 and 208 GeV and were processed using a full simulation of the OPAL detector [8].Distributions at intermediate energies were obtained by linear interpolation between the distributions at the nearest neighbouring available c.m. energies.In general, the Monte Carlo samples contain at least 10 times the statistics of the data.Signal e + e − → (Z/γ) * → qq(g) events were generated using the PYTHIA 6.125 [9,10] Monte Carlo program for the parton shower and hadronisation stage, which was interfaced with the KK2F program [11] to obtain a more accurate description of initial state radiation (ISR).The simulation parameters were tuned to OPAL data taken at the Z peak [12].The PYTHIA Monte Carlo model has previously been shown to provide a good description of e + e − annihilation data for c.m. energies up to 189 GeV (see e.g.[1,2,3,4,5,6,12]).
As an alternative to the string fragmentation model implemented in PYTHIA, events were generated using the HERWIG 5.9 [13] Monte Carlo program, also tuned to OPAL data as described in [2,12].The HERWIG Monte Carlo program implements the cluster fragmentation model.
In addition, we generated events of the type e + e − → 4 fermions (diagrams without intermediate gluons).These 4-fermion events, in particular those with four quarks in the final state, constitute the major background in this analysis.Simulated 4-fermion events with hadronic and leptonic final states were generated using the GRC4F 2.1 [14] Monte Carlo program.The final states were produced via s-channel or t-channel diagrams, including W + W − and ZZ production.This generator was interfaced to JETSET 7.4 [9] using the same parameter set [12] for the fragmentation and decays as used for (Z/γ) * events.The JETSET program implements the string fragmentation model.Two other possible sources of background events were simulated.Hadronic two-photon processes were evaluated using the VERMASEREN [15], HERWIG and PHOJET [16] programs.Production of e + e − → (Z/γ) * → τ + τ − was evaluated using the KK2F event generator in combination with the TAUOLA [17] decay library to simulate the decay of the tau leptons.
In addition to PYTHIA and HERWIG we also used the ARIADNE 4.08 [18] event generator to compare to our corrected data distributions.The ARIADNE program implements the colour dipole model to describe the parton shower process.For the fragmentation stage the generator was interfaced to the JETSET 7.4 program.The parameter set used for ARIADNE is documented in [19].This model provides a good description of global e + e − event properties at √ s = M Z , as do PYTHIA and HERWIG.
4 Data analysis

Selection of events
The majority of hadronic events produced at c.m. energies above the Z resonance are radiative events in which initial state photon radiation reduces the invariant mass of the hadronic system to about M Z .An experimental separation between radiative and non-radiative events is therefore required.In addition, at the present energies above 190 GeV, ZZ and W + W − production are kinematically possible and form a significant background to the (Z/γ) * → qq events.In this paper we use similar techniques to those of our previous analyses of e + e − annihilation data at lower c.m. energies [1,2,3], to select non-radiative (Z/γ) * → qq events.

Preselection
Hadronic events are identified using criteria as described in [20].The efficiency of selecting non-radiative hadronic events is greater than 98% at all c.m. energies, as can be seen in Table 1.
We define as particles tracks recorded in the tracking chambers and clusters recorded in the electromagnetic calorimeter.The tracks are required to have transverse momentum relative to the beam axis, p T > 150 MeV/c, a number of hits in the jet chamber, N hits ≥ 40, a distance of the point of closest approach to the collision point in the r-φ plane, d 0 ≤ 2 cm, and along the z axis, z 0 ≤ 25 cm.The clusters in the electromagnetic calorimeter are required to have a minimum energy of 100 MeV in the barrel and 250 MeV in the endcap sections.

ISR-fit selection
To reject radiative events, we determine the effective c.m. energy √ s ′ of the observed hadronic system as follows [21].
First isolated photons are identified by looking for energy deposits greater than 3 GeV in the electromagnetic calorimeter, having less than 1 GeV of additional energy deposited in a cone of 0.2 radians around its direction.The remaining particles are formed into jets using the Durham [22] algorithm with a value for the resolution parameter y cut = 0.02.A matching algorithm [23] is employed to reduce double counting of energy in cases where charged tracks point towards electromagnetic clusters.The energy of additional photons emitted close to the beam direction is estimated by performing three separate kinematic fits assuming zero, one or two such photons.Of the acceptable fits, the one with the lowest number of photons is selected.The value of √ s ′ is computed from the fitted momenta of the jets, excluding photons identified in the detector or close to the beam directions.The value of √ s ′ is set to √ s if the fit assuming zero initial state photons was selected.The 4-momenta of all measured particles are boosted into the rest frame of the observed hadronic system.
To reject events with large initial-state radiation (ISR), we require √ s ′ > √ s − 10 GeV.This is referred to as the "ISR-fit" selection.In Figure 2a we compare the √ s ′ distribution of the data set at √ s=200 GeV, after the preselection was applied, to simulated (Z/γ) * events and to 4-fermion and other background events.Simulated (Z/γ) * events are classified into radiative events, s ′ true < √ s − 1 GeV, where s ′ true is the true effective c.m. energy, and non-radiative events, which is the complement.Though about 27% of the selected (Z/γ) * events are by this definition radiative events, the fraction of selected events with s ′ true < √ s − 10 GeV is only 5%.The background from 4-fermion events and the efficiency of selecting non-radiative events are given in Table 1.

Final selection
The estimated background from e + e − → τ + τ − and two-photon events of the type γγ → qq is small at high √ s ′ .To reject this background further and to ensure events are well contained in the OPAL detector we require at least seven accepted tracks, and the cosine of the polar angle of the thrust axis, | cos θ T | < 0.9.The background from the τ + τ − and two-photon events in the final selected sample is estimated from Monte Carlo samples to be about 0.1% and is neglected.
To reduce the background of 4-fermion events in the remaining sample, we test the compatibility of the events with QCD-like production processes.A QCD event weight W QCD is computed as follows.We force each event into a four-jet configuration in the Durham jet scheme and use the EVENT2 [24] program to calculate the O(α 2 s ) matrix element |M(p 1 , p 2 , p 3 , p 4 )| 2 for the processes e + e − → qqqq, qqgg [25], where p 1 , p 2 , p 3 and p 4 are the momenta of the reconstructed jets.Since neither quark nor gluon identification is performed on the jets, we calculate the matrix element for each permutation of the jet momenta and use the permutation with the largest value for the matrix element to define the event weight, Note that the definition of the event weight contains kinematic information only and is independent of the value of α s .
The weight W QCD is expected to have large values for processes described by the QCD matrix element, originating from (Z/γ) * → qq, and smaller values for W + W − or ZZ events.In Figure 2b we compare the data distribution of W QCD , after the ISR-fit selection, to the expectations of our simulation.A good separation between the (Z/γ) * and 4-fermion events is achieved by requiring W QCD ≥ −0.5.
In addition to the cut on W QCD , used also in our earlier papers, we apply a likelihood based rejection of hadronic and semi-leptonic ZZ and W + W − decays (see [26] and references therein).The distributions for the semi-leptonic and hadronic likelihoods, after all selection cuts described above, are shown in Figure 3a and 3b, for the data and for the signal and background Monte Carlo expectations.Events are rejected if the semi-leptonic likelihood is greater than 0.5 or the hadronic likelihood is greater than 0.25, thus reducing the remaining 4-fermion background by more than a factor two while reducing the expected signal by only about 2.2%.After the final selection criteria the estimated efficiency for selecting signal events is around 76% while the remaining 4-fermion background in estimated to be around 5%.The individual numbers for each of the six energy ranges are given Table 1.

Correction procedure
The remaining 4-fermion background in each bin of each observable was estimated by Monte Carlo simulation and is subtracted from the observed bin content.A bin-by-bin multiplication procedure is then used to correct the observed distributions for the effects of detector resolution and acceptance as well as for the presence of remaining radiative (Z/γ) * events.A bin-by-bin correction procedure is suitable for the measured distributions as the effects of finite resolution and acceptance cause only limited migration (and therefore correlation) between bins and the data are well described by the Monte Carlo programs used to determine the corrections.
For the multiplicative correction each bin, after background subtraction, is corrected from the "detector level" to the "hadron level" using two samples of Monte Carlo (Z/γ) * events at each c.m. energy.The hadron level sample does not include initial state radiation or detector effects and allows all particles with lifetimes shorter than 3× 10 −10 s to decay.The detector level sample includes full simulation of the OPAL detector and initial state radiation and contains only those events which pass the same cuts as are applied to the data.The bin-by-bin correction factors are derived from the ratio of the distributions at the hadron level to those at the detector level.The correction factors tend to be around 1.15 in most bins, with values further away from 1 at the edges of the distributions.

Systematic uncertainties
The experimental systematic uncertainties are estimated by repeating the analysis with varied experimental conditions.To reduce statistical fluctuations in the magnitudes of the systematic uncertainties the average over three neighbouring bins is taken.In addition, all results and errors obtained for the different energy ranges are combined in a single set of results, as described in Section 4.4, thus further reducing any statistical fluctuations in the systematic uncertainties.
Possible inadequacies in the simulation of the response of the detector in the endcap regions are accounted for by restricting the analysis to the barrel region of the detector, requiring the thrust axis of accepted events to lie within the range | cos θ T | < 0.7.This reduces the event sample by approximately 26%.The corresponding systematic error is the deviation of the results from those of the standard analysis.
To evaluate uncertainties due to inadequacies in the track modelling the selection criteria for tracks were modified.The maximum allowed distance of the point of closest approach of a track to the collision point in the r-φ plane, d 0 , was changed from 2 to 5 cm, the maximum distance in the z direction, z 0 , from 25 to 10 cm and the minimum number of hits from 40 to 80.These modifications change the number of selected tracks by up to approximately 12%.In addition, track momenta in the Monte Carlo samples were smeared to degrade the resolution by 10%.The quadratic sum over the deviations from the standard result, obtained from each of these variations, is included in the total systematic uncertainty.
Uncertainties arising from the selection of non-radiative events are estimated by repeating the analysis using a different technique [2] to determine the value for √ s ′ .This technique differs from our standard √ s ′ algorithm in that in this case the kinematic fit always assumes one photon, either unobserved close to the beam direction, or detected in the ECAL.The final event sample with this √ s ′ algorithm has an overlap of approximately 97% with the standard sample and is approximately the same size.The difference relative to the standard result is included in the total systematic error.Systematic uncertainties associated with the subtraction of the 4-fermion background events are estimated by varying the cut on W QCD and on the 4-fermion likelihoods.We use W QCD ≥ −0.8 which increases the event sample by approximately 2%, and W QCD ≥ 0 which reduces the event sample by approximately 8%.We also vary the cuts on the likelihood values between 0.1 and 0.4 for hadronic 4-fermion events and between 0.25 and 0.75 for semi-leptonic events.The maximum deviation from the standard result obtained from the variations of the cuts on W QCD and the hadronic 4-fermion likelihood is included in the total systematic uncertainty, as is the maximum deviation obtained from the variations of the cut on the semi-leptonic likelihood.
In addition, we vary the predicted background to be subtracted, up and down by 5%, slightly more than its measured uncertainty at √ s= 189 GeV of 4% [26], and include the largest difference from the standard result in the systematic uncertainty.
The difference in the results when we determine the bin-to-bin correction factors from simulated (Z/γ) * events generated using HERWIG instead of PYTHIA is taken as the uncertainty in the modelling of the (Z/γ) * events.
All systematic uncertainties are added in quadrature to obtain the total systematic errors.The largest contribution to the total uncertainty is that obtained using HERWIG as an alternative hadronisation model in the correction procedure.Other significant contributions to the systematic uncertainties arise from varying the quality criteria on tracks and from varying the 4-fermion rejection cuts.

Combining results
The event selection, the correction procedure and the estimate of systematic uncertainties are performed independently for each of the six energy regions detailed in Table 1.Afterwards the results are combined to obtain a single set of results with the best possible statistical precision.The combination procedure assumes all systematic uncertainties to be fully correlated between the energy regions.This means data points obtained at different energies are combined into a weighted average, using the statistical errors to determine the weights.
To obtain the systematic uncertainty on the combined results all individual contributions to the systematic uncertainty, described in Section 4.3, are first determined separately in each energy region.For the combined results the individual contributions to the systematic uncertainty are then taken as the average over the corresponding values obtained in each of the six energy regions.Finally, the combination of the thus averaged individual contributions into the total systematic uncertainty on the combined results is done, as in Section 4.3, by adding all contributions in quadrature.

5
Results and comparison to theory There are different theoretical approaches in which fractional energy2 spectra of final state hadrons can be calculated.At high x, conventional perturbation theory, utilising expansions in powers of α s , is applicable.Non-perturbative effects occurring at the final stage of the hadronisation process are not calculable but they can be factorised out and left to be determined experimentally.At low x, perturbation theory can also be applied, provided that logarithmically enhanced terms (e.g.∼ ln(1/x)) are resummed to all orders.Remaining singular effects are not accounted for and as a consequence the predictions obtained are formally only valid for asymptotically high energies.Corrections to these asymptotic predictions can be calculated as an expansion in powers of √ α s .An important precondition to the validity of these asymptotic predictions is that the non-perturbative hadron formation process is local, i.e. the properties of the hadronic final state closely follow the properties of the partons before hadronisation.This presumed equivalence is referred to as Local Parton-Hadron Duality (LPHD) [27].For a recent review of particle production in the low x region we refer to [28].
In the following, we compare x p and ξ p distributions over a large range of c.m. energies to theoretical predictions applicable to either the low x or high x regions, in order to test the validity of these QCD approaches and the LPHD hypothesis.Since the ranges in x p or ξ p over which the different predictions are valid depend not only on the terms included in these predictions, but also on those that are missing, we do not quantitatively know these ranges a priori.Our approach, therefore, is to fit predictions to the data in a range where they can be observed to describe the data well.We subsequently vary these ranges and extrapolate the fitted predictions outside them to establish the limits of the validity of the predictions.

Fong and Webber prediction
At low x, QCD predicts destructive interference effects in soft gluon emissions, causing the asymptotic shape of the ξ distribution to be Gaussian [29].The position of the peak in the distribution is predicted to increase linearly with τ = ln ( √ s/2Λ), where Λ is the QCD scale, defined using the one loop expression: α s (τ ) = 2π/βτ , with β = 11 − 2N f /N c , N c the number of colours and N f the number of active flavours.Corrections to this asymptotic prediction were calculated up to O(α s ) by Fong and Webber [30], yielding a skewed Gaussian shape for the ξ spectrum: where δ = (ξ − ξ)/σ and, for quark initiated jets: with ρ = 11 + 2N f /N 3 c and ω = 1 + N f /N 3 c .The term ξ0 is a non-perturbative offset to the ξ distribution.
We have fitted Equation (2) simultaneously to the ξ p distribution measured in the present study and in previous studies by OPAL [31,1,2,3] and lower energy experiments [32,33], covering the range √ s = 14 − 202 GeV.As the energy scale involved in coherent parton emissions, as described by the formalism, is expected to be low, we apply the formalism assuming three active flavours.The free parameters in the prediction are Λ and ξ0 .The overall normalisation N (τ ) is also not predicted in the formalism of Equation (2).We therefore vary the normalisation separately at each c.m. energy in the fit.
When comparing the measured ξ p distributions to the predicted ξ distributions, we neglect the effects of the finite masses of the hadrons produced in the final state.For this reason the highest ξ p values, corresponding to momenta lower than 300 MeV/c, are excluded from the fit, avoiding the region where such mass effects could play a significant role.Moreover, as the formalism of Equation ( 2) is not expected to be valid at low ξ (high x), where logarithms of the form ln(1 − x) are important, data below ξ p = 2.0 are also not considered in the fit.
Figure 6 shows the ξ p data compared to the fitted predictions of Equation (2).For clarity, not all OPAL data included in the fit are shown in the figure, specifically our data at 161, 172, 183 and 189 GeV.The data in the range included in the fit (full lines) are well described by the prediction.At low ξ p the prediction clearly lies above the data, while at high ξ p it appears to give a reasonable description of the measurements beyond the fitted range.The total χ 2 of the fit is 121 for 174 degrees of freedom when including the full experimental errors (statistical and systematic errors added in quadrature) in the fit and assuming them to be uncorrelated.This is the only way to treat all data consistently since for the lower energy results only the combined statistical and systematic uncertainties are available.As a consequence, the quoted χ 2 value presumably underestimates the true χ 2 of the fit.
The obtained values for Λ and ξ0 are 153 ± 7 ± 55 MeV and −0.70 ± 0.03 ± 0.23, respectively, where the first errors are the fit uncertainties, defined as the square root of the diagonal elements of the covariance matrix, and the second errors were obtained by varying the fit range.The lower limit on ξ p was varied between 1.5 and 2.5 and the upper limit between the ξ p values corresponding to particle momenta of 200 and 400 MeV/c.The errors on the parameters do not include an estimate of how they may be affected by missing terms in the predictions.Therefore these parameters are only meaningful within the formalism.A more detailed test of the formalism of Equation ( 2) than the one presented here would need to account for more subtle effects such as e.g. the changing flavour composition with c.m. energy.For this, additional information would need to be included in the fit, which goes beyond the scope of our analysis.

MLLA prediction
A prediction for the ξ distribution can also be formulated in the so-called modified leading-logapproximation (MLLA, for a review see e.g.[28]).This approximation constitutes a complete resummation of single and double logarithmic terms [29,34].In the MLLA approach a so-called limiting spectrum prediction for the ξ distribution can be formulated in which the QCD scale, Λ, and the scale at which the perturbative process terminates, Q 0 , are taken to have the same value.In this study we use the equations given in [28]: where C F = (N 2 c − 1)/2N c and K h is a hadronisation constant which accounts for the number of hadrons of type h produced per final state parton.In some earlier publications the factor of 2 in Equation ( 7), accounting for the two quark jets present in e + e − annihilation, was incorporated into K h .From the LPHD hypothesis K h is expected to be independent of the underlying process, including its scale.For the expression for the limiting spectrum distribution, D lim , we refer to [28].In this formalism τ = ln ( √ s/2Λ eff ), where Λ eff is an effective scale factor representing both the QCD scale, Λ, and the cut-off scale, Q 0 .Contrary to the Fong and Webber prediction, the MLLA limiting spectrum formalism predicts not only the shape of the ξ distribution but, apart from the hadronisation correction, also its normalisation.
In Figure 7a, the same ξ p data shown in Figure 6 are compared to the fitted predictions of Equation (7), where the free parameters are the effective scale, Λ eff , and a separate hadronisation constant, K h , for each c.m. energy.We choose to fit a separate hadronisation constant at every energy in order subsequently to examine the obtained K h values to establish how well the formalism and the LPHD hypothesis work.Higher order effects, not accounted for in the formalism, are expected to give rise to some energy dependence in the K h factors, in particular at the lowest c.m. energies (see e.g.[35]).As in Section 5.1.1,we apply the formalism assuming three active flavours.The fit is restricted to the region around the peaks defined by 0.75 + 0.33 log( √ s ) < ξ p < 0.9 + 0.8 log( √ s ).
We find the fit describes the data well within the fitted range (full line) as illustrated by the total χ 2 of the fit of 122 for 132 degrees of freedom.Outside the fit range, both at low and high ξ p , deviations of the predictions from the data are observed.The MLLA fit is seen to describe the low ξ p data better, but the high ξ p data worse, than the fit of the Fong and Webber prediction discussed in Section 5.1.1 (see Figure 6).The value of Λ eff obtained from the fit is 254 ± 3 ± 31 MeV, where the first error is the fit uncertainty and the second error is obtained from varying the fit range.Both the lower and upper limits on ξ p were lowered and raised by 0.5 units in ξ p .
In Figure 7b the obtained K h factors are shown as a function of the c.m. energy.The errors shown are the fit error and the combined error from the fit and from varying the fit range, as above.Above 130 GeV the K h factors obtained from our fit are consistent with being constant.Below 50 GeV a modest rise in the values of the K h factors is observed.

The energy dependence of ξ 0
The MLLA also provides a prediction for the energy evolution of the peak position, ξ 0 , of the ξ distribution [36,35]: with τ defined as in Section 5.1.1 and C = a 2 /16N c b.
To determine the peak position ξ 0 of the measured ξ p distributions in each of the energy regions of the present study, we fitted Equation ( 2) to these distributions in the range 2.0 < ξ p < 6.2, taking the position of the maximum of the fitted function as the result for ξ 0 .The systematic uncertainties are determined by repeating the fit using the systematic variations described in Section 4.3.In addition, the uncertainty due to the choice of the fit range is estimated by repeating the fit using two alternative fit ranges: 3.2 < ξ p < 4.8 and 1.6 < ξ p < 6.6, taking the larger deviation from the nominal result as the systematic uncertainty.The uncertainties are added in quadrature to obtain the total systematic uncertainty.As for the momentum distributions, we combine the ξ 0 results from the different energies to obtain: ξ 0 (201.7 GeV) = 4.158 ± 0.007(stat.)± 0.036(syst.).At c.m. energies below 91 GeV, with the exception of the TOPAZ result [33], no ξ 0 determinations were published.Therefore, we determined ξ 0 values for all energies below 91 GeV by fitting Equation (2) to the published ξ p (or x p ) results in a narrow region around the peak.This procedure ensures that the ξ 0 results are obtained in a uniform manner over the entire range of energies studied.The error on ξ 0 for the lower energies is taken from the variation in the peak position when modifying the fit range.For energies above 91 GeV this was found to be the dominant source of uncertainty.
In Figure 8 our ξ 0 result is shown together with earlier ξ 0 results from OPAL [31, 1, 2, 3] and other LEP experiments [5,37,38] and the ξ 0 values we determined for the lower energy experiments [32,33,39,40,41].The results are compared to the predictions of PYTHIA, HERWIG and ARIADNE and to the fitted prediction of Equation (8).For the fit three active flavours were assumed and the scale, Λ, was the only free parameter.The fit yields Λ = 203 ± 2 MeV, where the error is the fit uncertainty.The fitted MLLA prediction is in good agreement with the data.The PYTHIA and ARIADNE predictions are slightly below the data, while the HERWIG prediction, in particular at the highest c.m. energies, is much lower than the data.

Comparisons to theory at high x
Another theoretical framework for momentum distributions is formulated in terms of the scale evolution of fragmentation functions.In this framework, a factorisation between the perturbative scattering process on one hand, and the parton emission cascade and non-perturbative hadronisation stage on the other hand, is formulated as follows (see e.g.[42,43,44]): The C f (x, α s , √ s) are coefficient functions, describing the probability to obtain a parton f (gluon, g, or quark, u, d, s, c, b) with energy fraction x in an e + e − collision, and D h f (x, µ) are fragmentation functions, describing the expected distribution of final state hadrons h from an initial parton f at scale µ.The x dependence of fragmentation functions is not known a priori, but their dependence on the energy scale is predicted by theory.In QCD this scale dependence, also known as scaling violation, is described by parton evolution equations [45], or "DGLAP" equations.These equations, as well as the coefficient functions of Equation ( 9), are currently known to next-to-leading order (NLO) accuracy in α s [46,42,43,44].
Using a computer program [47] written for scaling violation analyses of proton structure functions in NLO QCD, to which we added NLO evolution equations for fragmentation functions and NLO coefficient functions for e + e − annihilation, we performed a fit of Equation ( 9) to OPAL and lower energy data.These data include the inclusive charged particle x p distributions measured here and at lower energies [50,1,2,3,32,39,48,49], flavour tagged x p distributions [50], and results on the helicity components in charged particle production [51].
The fit involves defining fragmentation functions D h f (x, µ 0 ) at an arbitrarily chosen input scale, µ 0 , evolving these over the entire relevant range of scales using the parton evolution equations, and using them to calculate the theory expectation for the measured distributions (Equation ( 9)).The a priori unknown input consists of the input fragmentation functions and the strong coupling constant α s .The fragmentation functions are parameterised, similarly to reference [51], as where the parameters N i , α i , β i and c i are determined independently for i = g, u, (ds), c, b.
The separation between the fragmentation functions of the different quark flavours relies on the (uds), c and b flavour tagged data included in the fit and on the fact that as a consequence of the different electric and weak charges the relative contributions from d, s and b quarks compared those from u and c quarks, vary with the c.m. energy.No data included in the fit provide information to separate the contributions from d and s quarks.Therefore only one fragmentation function is defined for these two flavours.
In earlier studies [42,44,52,53,54], similar fits were performed, including data up to c.m. energies of 91 GeV, to obtain measurements of α s .The higher energy LEP 2 data do not add significantly to the achievable precision on α s , which is mostly constrained by the very precise low energy and LEP 1 results.We therefore choose, in our default fit, to fix α s at its current best estimate of α s (M Z ) = 0.1181 ± 0.002 [55], to test how well the fitted predictions describe the data over the now extended range of scales.The data included in the fit are restricted to x p values between 0.07 and 0.80.In addition, the particle momenta were required to be greater than 1 GeV/c, effectively raising the lower x p limit for the lowest energy experiments.This cut is motivated by the fact that the present formalism does not account for destructive interference effects in soft gluon emissions (see Section 5.1).
In Figure 9, x p distributions measured at c.m. energies from 14 to 202 GeV are compared to the fitted predictions.The theory predictions are in good agreement with the data in the fitted range (full lines).The total χ 2 of the fit is 198 for 188 degrees of freedom.As for the fits described in Section 5.1, this χ 2 result is based on the full (statistical and systematic) uncertainties on the data and neglects any correlations between data points.Outside the fitted range, at low x p and low c.m. energies, charged particle production is suppressed.This suppression is not described by the predictions, which do not include soft gluon interference effects, as stated above.
Scaling violations in charged particle production are shown in Figure 10, where the 1/σ • dσ ch /dx p results are shown as a function of √ s for a number of x p bins.Negative scaling violations are visible for the higher x p values.At low x p the charged particle production rate demonstrates little dependence on the c.m. energy.Changes in the flavour composition with c.m. energy also influence the observed level of scaling violations, especially at high x p and around √ s = 91 GeV where a relatively large fraction of events is expected to contain b quarks.
The data in Figure 10 are also compared to the fitted NLO predictions.Overall, a good agreement between data and theory is found in the region of the fit (full lines).We thus conclude that inclusive charged particle production data can be described by the QCD based parton evolution formalism, using the world-average value for α s .When we leave α s as a free parameter in the fit, we obtain α s (M Z ) = 0.113 ± 0.005 ± 0.007, where the first error is the fit uncertainty and the second the uncertainty obtained by independently varying the renormalisation and factorisation scales up and down by a factor 2. For the dependence of Equation ( 9) and the parton evolution equations on these scales we refer to [45,46,42,43,44].

Summary and conclusion
We have presented a measurement of charged particle momentum distributions in hadronic e + e − annihilation events produced at LEP 2 at c.m. energies between 192 and 209 GeV.The results are combined into a single set of distributions, corresponding to an average c.m. energy of √ s = 201.7 GeV.
The results, in combination with those obtained at lower c.m. energies, were compared to QCD Monte Carlo models and to analytical QCD predictions calculated in various approaches.In general, good agreement is observed between the data and these theory predictions in the regions where the latter are expected to be valid.
The p in ⊥ , p out ⊥ , y, x p and ξ p distributions, averaged over the full c.m. energy range √ s = 192 − 209 GeV, are tabulated in Tables 2 to 6 and shown in Figures 4 and 5.These distributions correspond to an average c.m. energy of √ s = 201.7 GeV with a standard deviation of 4.8 GeV.They are compared to the corresponding predictions of the PYTHIA, HERWIG and ARIADNE Monte Carlo programs.In general, a reasonable agreement is observed between the predictions and the data for the different distributions.The HERWIG Monte Carlo program agrees better with the data than the PYTHIA and ARIADNE models.The p out ⊥ distribution is harder in the data than predicted by any of the Monte Carlo programs.

Figure 1 :
Figure 1: Integrated luminosity collected by OPAL in 1999 and 2000 at different c.m. energies.

Figure 2 :
Figure 2: Distributions of (a) the effective c.m. energy √ s ′ and (b) the QCD event weight W QCD in the data at √ s = 200 GeV (full points).The data are shown with statistical errors only.Estimated signal and background (BG) contributions are shown as detailed in the upper figure.The vertical lines indicate where the selection cuts are applied.

Figure 3 :
Figure 3: Distributions of the likelihood for events to be (a) a semi-leptonic 4-fermion event or (b) a hadronic 4-fermion event, in the data at √ s = 200 GeV (full points).The data are shown with statistical errors only.Estimated signal and background (BG) contributions are shown as detailed in the upper figure.The vertical lines indicate where the selection cuts are applied.

Figure 4 :
Figure 4: Distributions of the charged particle momenta (a) p in ⊥ and (b) p out ⊥ and of the rapidity (c) y at √ s = 201.7 GeV, compared to the corresponding PYTHIA, HERWIG and ARIADNE predictions.The statistical errors and the combined statistical and systematic errors are shown as the inner and outer error bars, respectively.

Figure 5 :Figure 6 :
Figure 5: Distributions of the charged particle fractional momenta (a) x p and (b) ξ p at √ s = 201.7 GeV, compared to the corresponding PYTHIA, HERWIG and ARIADNE predictions.The statistical errors and the combined statistical and systematic errors are shown as the inner and outer error bars, respectively.

Figure 7 :
Figure 7: a) ξ p spectra measured in the range √ s = 14 − 202 GeV, compared to MLLA predictions in the limiting spectrum approach ( [28], Equation (7) in this paper), fitted to the data.The full lines indicate the region of the fit.For clarity, not all OPAL data included in the fit are shown in the figure.The error bars on the data represent the combined statistical and systematic uncertainties.b) The hadronisation constants K h at different √ s as obtained from the fit, described in Section 5.1.2.The inner error bar represents the fit error and the outer error bar the combined fit error and the error obtained from varying the fit range.

Figure 8 :
Figure 8: Evolution of the position of the peak in the ξ p distribution, ξ 0 , with c.m. energy, in the range √ s = 14 − 202 GeV, compared to the fitted MLLA prediction (Equation (8) in this paper) and to predictions of PYTHIA, HERWIG and ARIADNE.The error bars of the data represent the combined statistical and systematic uncertainties.

Figure 10 :
Figure 10: Inclusive charged particle rate, 1/σ • dσ ch /dx p , as a function of √ s in bins of x p , compared to fitted NLO predictions (Equation (9) in this paper) using the world average value of α s (M Z ) = 0.1181.The full lines indicate the region of the fit.The uncertainties of the data represent the combined statistical and systematic errors.

Table 1 :
Data samples corresponding to six ranges in c.m. energy.The luminosities and the luminosity weighted mean c.m. energies are given in the first two rows.The rows labelled as "Preselection", "ISR-fit selection" and "Final selection" correspond to the number of events passing these criteria.For each of these selection stages the expected number of events, the remaining fraction of 4-fermion background and the signal efficiency are also given.The errors on the integrated luminosity include statistical and systematic uncertainties added in quadrature.All other errors are the statistical errors only.