First Measurement of the Double-inclusive B/B Hadron Energy Distribution in e + e − Annihilations, and of

We have made the first measurement of the double-inclusive B/Bbar energy distribution in e+e- annihilations, using a sample of 400,000 hadronic Z0 decay events recorded in the SLD experiment at SLAC between 1996 and 1998. The small and stable SLC beam spot and the CCD-based vertex detector were used to reconstruct B/Bbar decay vertices with high efficiency and purity, and to provide measurements of the kinematic quantities used to calculate the B energies in this novel technique. We measured the B/Bbar energies with good efficiency and resolution over the full kinematic range. We measured moments of the scaled energies of the B and Bbar hadrons vs. the opening angle between them. By comparing these results with perturbative QCD predictions we tested the ansatz of factorisation in heavy-quark production. A recent next-to-leading order calculation reproduces the data.


Introduction
The production of heavy hadrons (H ) in e + e − annihilation provides a laboratory for the study of heavyquark (Q) jet fragmentation. This is commonly char-Open access under CC BY license. acterized in terms of the observable x H ≡ 2E H / √ s, where E H is the energy of a B or D hadron containing a b or c quark, respectively, and √ s is the c.m. energy. The distribution of x H , D(x H ), is conventionally referred to as the heavy-quark 'fragmentation function'. 1 In recent publications we presented [1,2] the results of a new method for reconstructing B-hadron decays, and the B energy, inclusively, using only charged tracks, in the SLD experiment at SLAC. We used the upgraded charge-coupled device (CCD) vertex detector, installed in 1996, to reconstruct B-decay vertices in Z 0 decays with high efficiency and purity. Combined with the micron-sized SLC interaction point (IP), our precise vertexing allowed us to reconstruct the total transverse momentum of the tracks from Bdecays, and therefore the transverse momentum and mass associated with the neutral particles in the Bdecays. This allowed us to reconstruct accurately the energy of B hadrons. These studies yielded the most precise measurement of the b-quark fragmentation function, and allowed us to test models of heavy-quark fragmentation. Of the 9 models tested, only 4 were consistent with our precision data at better than the 1% level based on a χ 2 probability. This allowed us to reduce the model-dependent systematic uncertainty on the b-quark fragmentation function.
We have extended these studies and applied similar 'topological' vertexing techniques to tag events in which we reconstructed the energies of both leading B hadrons produced via e + e − → bb → BB + X. We measured the moments of the single-inclusive Bhadron scaled-energy distribution dN/dx B : as well as the moments of the double-inclusive scaledenergy distribution: where x B1 and x B2 are the scaled energies of the two B hadrons and the label is arbitrary, φ is the angle between their flight directions, and i and j are integers 1. We formed the normalised moments: Following the method proposed in [3] we used these quantities to test the ansatz of factorisation as applied to perturbative quantum chromodynamics (pQCD) calculations of e + e − → bb events.
Due to 'soft', or non-perturbative, effects, D i and D ij cannot be predicted absolutely in pQCD. Rather, the respective pQCD calculation must be folded with models of the non-perturbative ('hadronisation') process in order to derive predictions that can be compared with experimental data. However, provided that the ansatz of factorisation holds [3], namely that calculation of the perturbative and non-perturbative phases can be separated by (an arbitrary) factorisation scale µ F , the dependence on µ F cancels in Eq. (3) and hence G ij can be calculated absolutely in pQCD, up to possible 'higher twist' effects of order 1/ √ s, with no dependence on hadronisation models. Comparison of the measured G ij with pQCD predictions hence allows both a test of this ansatz and of the perturbative calculations. The G ij can be derived at leading order (LO) in pQCD using the numerical results in [3]. In addition, next-to-leading order (NLO) predictions for G ij have been calculated recently [4].
In Section 2 we describe the detector and the selection of e + e − → hadrons events used in this analysis. We present in Section 3 the first measurement of the double-inclusive B-hadron energy distribution and, in Section 4, of the normalised moments Eq. (3). In Section 5 we describe the estimation of the errors on our measurements. Finally, in Section 6, we test the ansatz of factorisation and compare our data with the pQCD predictions.

Apparatus and hadronic event selection
This analysis is based on roughly 400 000 hadronic Z 0 events produced in e + e − annihilations at a mean center-of-mass energy of √ s = 91.28 GeV at the SLAC Linear Collider (SLC), and recorded in the SLC Large Detector (SLD) between 1996 and 1998. A general description of the SLD can be found elsewhere [5]. The trigger and initial selection criteria for hadronic Z 0 decays are described in Ref. [6]. This analysis used charged tracks measured in the Central Drift Chamber (CDC) [7] and in the upgraded Vertex Detector (VXD3) [8]. Momentum measurement was enabled by a uniform axial magnetic field of 0.6 T. The CDC and VXD3 gave a momentum resolution of σ p ⊥ /p ⊥ = 0.01 ⊕ 0.0026 p ⊥ , where p ⊥ is the track momentum transverse to the beam axis in GeV/c. In the plane normal to the beamline the centroid of the micron-sized SLC IP was reconstructed from tracks in sets of approximately thirty sequential hadronic Z 0 decays to a precision of σ rφ IP 4 µm. The IP position along the beam axis was determined event by event, using charged tracks, with a resolution of σ z IP 20 µm. Including the uncertainty on the IP position, the resolution on the charged-track impact parameter (d) projected in the plane perpendicular to the beamline was σ rφ d = 8 ⊕ 33/(p sin 3/2 θ) µm, and the resolution in the plane containing the beam axis was σ z d = 10 ⊕33/(p sin 3/2 θ) µm, where θ is the track polar angle with respect to the beamline.
A set of cuts was applied to the data to select wellmeasured tracks and events well contained within the detector acceptance. Charged tracks were required to have a distance of closest approach transverse to the beam axis within 5 cm, and within 10 cm along the axis from the measured IP, as well as | cos θ | < 0.87, and p ⊥ > 0.15 GeV/c. Events were required to have a minimum of five such tracks, a thrust axis [9] polar angle w.r.t. the beamline, θ T , within | cos θ T | < 0.80, and a charged visible energy E vis of at least 20 GeV, which was calculated from the selected tracks assigned the charged pion mass. The efficiency for selecting a well-contained Z 0 → qq(g) event was estimated to be above 97% independent of quark flavor. The selected sample comprised 313 447 events, with an estimated 0.10 ± 0.05% background contribution dominated by For the purpose of estimating the efficiency and purity of the selection procedures we made use of a detailed Monte Carlo (MC) simulation of the detector. The JETSET 7.4 [10] event generator was used, with parameter values tuned to hadronic e + e − annihilation data [11], combined with a simulation of B-hadron decays tuned [12] to Υ (4S) data and a simulation of the SLD based on GEANT 3.21 [13]. Inclusive distributions of single-particle and event-topology observ-ables in hadronic events were found to be well described by the simulation [6]. Uncertainties in the simulation were taken into account in the systematic errors (Section 5).

B-hadron selection and energy measurement
The event sample for this analysis was selected using a 'topological' vertexing technique based on the detection and measurement of charged tracks, which is described in detail in Refs. [14][15][16]. We considered events in which we found decay vertices corresponding to both the leading B and B hadrons. First, the Durham algorithm [17] was applied to the selected hadronic events, with a y c parameter value of 0.015, in order to define a jet structure in each event. We found that this algorithm and y c value minimized the number of B (and D) decay tracks assigned to the wrong jet. This is an important feature for our analysis because we used vertex-related variables derived only from charged tracks. Events containing 2, 3, or 4 jets were retained for further analysis.
In each selected event, the vertexing algorithm was applied to the set of tracks in each jet. Vertices consistent with photon conversions or K 0 or Λ 0 decays were discarded. Events were retained in which a vertex was found in exactly two jets. 35 137 events were selected, of which 89.4% were estimated to be of bb origin. The efficiency for selecting true bb events was estimated to 36.3%.
The large masses of the B hadrons relative to light-flavor hadrons make it possible to distinguish B-hadron decay vertices from those vertices found in events of lighter flavors using the vertex invariant mass, M. However, due to those particles missed from the vertex, which are mainly neutrals, M cannot be fully determined. M can be written where M ch and M 0 are the total invariant masses of the set of vertex-associated tracks and the set of missing particles, respectively. P t is the total chargedtrack momentum transverse to the B flight direction, which, by momentum conservation, is identical to the transverse momentum of the set of missing particles. P chl and P 0l are the respective momenta along the B flight direction, which we take to be the vector joining the IP to the vertex. The lower bound for the mass of the decaying hadron, the 'P t -corrected vertex mass' [2], was used as the variable for selecting B hadrons. Fig. 1 shows the distribution of M Pt for vertices in the selected event sample, and the corresponding simulated distribution. Events were selected that contained at least one vertex with M Pt > 2.0 GeV/c 2 and M Pt 2 × M ch . The latter cut was found to reduce the contamination from fake vertices in light-quark events. In order to improve the bb purity of the sample, events were selected in which both vertices had a flight length, d vtx , such that 0.1 < d vtx < 2.3 cm; in which at least one vertex contained two 'significant' tracks, i.e., tracks with a normalised impact-parameter significance, d/σ d , of at least 2 units; and in which the angle between the vertex flight vectors, φ, satisfied cos φ < 0.99. The last cut was effective at removing events in which either a gluon had split into heavy quarks, or the jet-finder had artificially split a single heavy-quark jet into two jets.
The energy of each B hadron, E B , can be expressed as the sum of the reconstructed-vertex energy, E ch , and the energy of those particles not associated with the vertex, E 0 . We can write . The two unknowns, M 0 and P 0l , must be found in order to obtain E 0 . One kinematic constraint can be obtained by imposing the B-hadron mass on the vertex, where P B = P chl + P 0l is the total momentum of the B hadron. From Eq. (4) we derive the following inequality, where equality holds in the limit where both P 0l and P chl vanish in the B-hadron rest frame. Eq. (7) effectively sets an upper bound on M 0 , and a lower bound is given by zero Because M 0 peaks near M 0 max , [15] we set M 2 0 = M 2 0 max if M 2 0 max 0, and M 2 0 = 0 if M 2 0 max < 0. We calculated P 0l : and hence E 0 (Eq. (6)). We then reconstructed the B-hadron energy, E rec B = E 0 + E ch . Events were retained in which both reconstructed B energies satisfied E rec B < 60 GeV. A final sample of 19 809 events was obtained with an estimated bb selection efficiency of 21.7% and a background contribution of only 0.23%, which was almost entirely from Z 0 → cc events.
The energy resolution of the final B sample is shown in Fig. 2, where we plot the normalised residual on The distribution was fitted to a double Gaussian function for which the mean positions, widths and normalisations were allowed to vary. 79.5% of the population lies in the 'core' Gaussian, of width 21.3%; the remaining population is characterised by a 'tail' Gaussian of width 31.3%.   The angular resolution was similarly investigated and is shown in Fig. 3, where we plot the residual on φ: φ true − φ rec . The vast majority of angles are reconstructed to better than 5 mrad. The tail in the resolution function corresponds to B decays with shorter decay lengths, where the vertex position error had a larger relative effect on the determination of φ. The slight asymmetry in the resolution function is an artefact of the resolution folded with the steeply-falling distribution of φ, which causes a small bias towards larger reconstructed vaues; this bias is explicitly corrected in the analysis (see below).
The double-inclusive distribution of raw B-hadron energies is shown in Fig. 4. Since we do not distinguish between B and B hadrons, the reconstructed Bdecay vertices were labelled arbitrarily '1' and '2' on an event-by-event basis.
We divided E rec B by the beam energy, E beam = √ s/2, to obtain the reconstructed scaled B-hadron energy, x rec

Angle-dependent BB energy moments
In each event we quantified the correlations between the two B hadrons in terms of the angle de-pendent scaled-energy moments proposed in [3]. Using the raw measured scale B energies we first evaluated the moments (Eq. (1)) from the raw measured distribution: where the sum is over the set of reconstructed B hadrons and N is the number of events in the sample. Similarly, we evaluated in each cos φ bin the double moments (Eq. (2)) where the sum extends over the set of events in each cos φ bin. The normalised moments G ij (Eq. (3)) were evaluated: The first six moments, i = 1, 2, 3 and j = 1, . . ., i are shown in Fig. 5. The bin centers were defined by taking the average value of cos φ within each bin. Also shown in Fig. 5 is a comparison with the simulated normalised moments; the simulation reproduces the data. Given that we showed previously [1,2] that the Peterson function implemented in our simulation does not provide a good description of the b-quark fragmentation function, this agreement may naively appear to be surprising. However, if, as proposed in [3], the non-perturbative contributions to the normalised quantity G ij cancel, the agreement should be excellent, as observed.
We used our simulated event sample to correct for the effects of the detector acceptance, the efficiency of the technique for reconstructing B-hadron decays, the energy resolution, and bin migrations caused by the finite angular resolution. We defined a binwise correction factor: where G gen ij (φ)| MC is the normalised moment calculated using generated true e + e − → BB + X events, and G rec ij (φ)| MC is the corresponding moment calculated after simulation of the detector and application of the same analysis as applied to the data. For cos φ ∼ −1 the value of F ij is close to unity. As cos φ increases towards 1, F ij rises monotonically to approximately 8 (F 11 ) or 2 (F 33 ). The increase with cos φ reflects our decreasing efficiency to select BB events as the angle between the two B-decay vertices becomes smaller and the B energies decrease. For a given cos φ bin, F ij decreases as i and j increase due to the effective weighting of this efficiency by . We derived the true normalised moments: Fig 6. The uncertainties associated with this correction procedure are discussed in the next section.

Error estimation
Since it is not possible to calculate them analytically, we used our simulation to estimate the statistical error on the moments. The simulated event sample was divided into subsamples, each comprising the same number of events as in the data sample. The entire analysis was performed on each of these subsamples. Within each cos φ bin the r.m.s. deviation of the ensemble of results was evaluated and taken as the statistical error in that bin. For a given pair of i and j values the errors on the particular moment are not correlated between different bins of cos φ, but within a given cos φ bin the errors on the different moments are correlated, since the same set of events was used.
We considered sources of systematic uncertainty that potentially affect our measurements. These may be divided into uncertainties in modeling the detector and uncertainties on experimental measurements serving as input parameters to the underlying physics modeling. For these studies our simulation was used.
Since our energy reconstruction technique is strongly dependent on charged-track properties, four sources of uncertainty were investigated: our simulated tracking efficiency, transverse momentum resolution, and the resolutions on the track impact parameter and polar angles. In each case the simulation was corrected so as to reproduce the data. The full analysis was repeated on the data, and half the difference between the results obtained using the corrected and uncorrected simulations was taken as a symmetric systematic error.
A large number of measured quantities relating to the production and decay of charm and bottom hadrons are used as input to our simulation. In bb events we considered the uncertainties on: the branching fraction for Z 0 → bb; the rates of production of B ± , B 0 and B 0 s mesons, and B baryons; the lifetimes of B mesons and baryons; and the average Bhadron decay charged multiplicity. In cc events we considered the uncertainties on: the branching fraction for Z 0 → cc; the charmed hadron lifetimes, the charged multiplicity of charmed hadron decays, the production of K 0 from charmed hadron decays, and the fraction of charmed hadron decays containing no π 0 s. We also considered the rates of production of secondary bb and cc from gluon splitting. The uncertainty on the world-average value [2] of each quantity was used to rederive the corrected G ij (φ). In each bin of G ij (φ) the deviation between the rederived and standard values was taken as an estimate of the corresponding systematic error. Other relevant systematic effects such as variation of the event selection cuts and the assumed B-hadron mass were found to be very small.
The error associated with the choice of b-quark fragmentation function used in the simulation was estimated by rederiving the corrected data using in turn the UCLA, Kartvelishvili and Bowler fragmentation functions [2]. In each cos φ bin the r.m.s. deviation of the results w.r.t. the standard value, using the Peterson function, was taken as an estimate of the systematic error. This error was typically much smaller than that arising from the other error sources.
For each systematic error investigated, the shape of the angular dependence of the moments was not significantly changed, only the overall normalization was affected. For each moment, in each cos φ bin the dominant errors were typically those related to the uncertainties on the charged-track properties. The errors due to charm and bottom physics modeling were typically an order of magnitude smaller. In each cos φ bin all sources of systematic uncertainty were added in quadrature to obtain the total systematic error. The fully-corrected G ij , with statistical and systematic errors added in quadrature, are shown in Fig. 6.

Comparison with perturbative QCD calculations
The fully-corrected data [16] were compared (Fig. 6) with a recent calculation [4] of the normalised moments complete at NLO in pQCD. Both the LO and full NLO calculations are shown in Fig. 6; the calculations assume an α s (M 2 Z ) value of 0.120 and a pole bquark mass value of 5.0 GeV/c 2 . 2 It can be seen that the difference between the two is relatively small, and that the LO calculation lies systematically slightly below the NLO calculation. The small size of the NLO relative to the LO contributions is an indication that the normalised moments are perturbatively robust observables. For each moment shown, the LO calculation undershoots the data. The NLO calculation reproduces the data across the full range of cos φ.
This comparison does not rely on any convolution of the pQCD calculations with models of the nonperturbative hadronisation process. Hence the excellent agreement between the pQCD calculations and the data verifies the ansatz of factorization between the perturbative and non-perturbative phases that forms the basis for the pQCD calculation of heavy-hadron properties.