UvA-DARE (Digital Academic Repository) Test of the tau-model of Bose-Einstein correlations and reconstruction of the source function in hadronic Z-boson decay at LEP

Bose–Einstein correlations of pairs of identical charged pions produced in hadronic Z decays are analyzed in terms of various parametrizations. A good description is achieved using a Lévy stable distribution in conjunction with a model where a particle’s momentum is correlated with its space–time point of production, the τ -model. Using this description and the measured rapidity and transverse momentum distributions, the space–time evolution of particle emission in two-jet events is reconstructed. However, the elongation of the particle emission region previously observed is not accommodated in the τ -model, and this is investigated using an ad hoc modiﬁcation.

Abstract Bose-Einstein correlations of pairs of identical charged pions produced in hadronic Z decays are analyzed in terms of various parametrizations. A good description is achieved using a Lévy stable distribution in conjunction with a model where a particle's momentum is correlated with its space-time point of production, the τ -model. Using this description and the measured rapidity and transverse momentum distributions, the space-time evolution of particle emission in two-jet events is reconstructed. However, the elongation of the particle emission region previously observed is not accommodated in the τ -model, and this is investigated using an ad hoc modification.

Introduction
In particle and nuclear physics, intensity interferometry provides a direct experimental method for the determination of sizes, shapes and lifetimes of particle-emitting sources (for reviews see [1][2][3][4][5]). In particular, boson interferometry provides a powerful tool for the investigation of the space-time structure of particle production processes, since Bose-Einstein correlations (BEC) of two identical bosons, first observed in 1959 [6,7] and most recently at the Large Hadron Collider [8][9][10][11][12], reflect both geometrical and dynamical properties of the particle-radiating source.
In e + e − annihilation BEC have been observed [13] to be maximal when the invariant momentum difference of the bosons, Q = −(p 1 − p 2 ) 2 , is small, even when one of the relative momentum components is large. This is not the case either in hadron-hadron interactions [14] or in heavy-ion interactions [15,16], where BEC are found not to depend simply on Q, but to decrease when any of the relative momentum components is large, a behavior that can be described by hydrodynamical models of the source [5,17].
The size (radius) of the source in heavy-ion collisions has been found to decrease with increasing transverse momentum, p t , or transverse mass, m t = m 2 + p 2 t , of the bosons. This effect can also be explained by hydrodynamical models [17,18]. A similar effect has been seen in pp collisions [19], as well as in e + e − annihilation [20][21][22]. A model for BEC in e + e − annihilation, known as the τ -model, which results in BEC depending simply on Q, rather than on its components separately, has been proposed [23]. In this model the simple Q dependence is a consequence of a strong correlation between a particle's fourmomentum and its space-time point of production. This strong correlation also leads to a decrease of the observed source size with increasing m t [24].
On the other hand, studies of BEC in e + e − annihilation at LEP have found that the BEC correlation function does depend on components of Q, the shape of the source being elongated along the event axis [22,[25][26][27][28]. At HERA a similar elongation is observed in neutral current deep inelastic ep scattering [29]. Note that here "source" refers not to the entire volume in which particles are emitted, but rather the smaller "region of homogeneity" from which pions are emitted that have momenta similar enough to interfere and contribute to the correlation function. The size of this region of homogeneity is sometimes referred to as the correlation length.
The main purpose of the present paper is to test the τ -model, but the question of to what extent BEC in e + e − annihilation depend differently on different components of Q is also considered. Here we study BEC in hadronic Z decay. We investigate various static parametrizations in terms of Q. Fitting over a larger Q range than in previous studies, we find that none of these parametrizations gives an adequate description of the Bose-Einstein correlation function. Next we investigate the τ -model and find that it provides a good description of BEC, both for two-jet and three-jet events. The results for two-jet events, from fits in terms of Q and the transverse masses of the two pions with respect to the event axis, are used to reconstruct the complete spacetime picture of the particle emitting source of two-jet events in hadronic Z decay. Finally, the relevance of the elongation mentioned above is investigated using an ad hoc modification of the τ -model.

Data
The data used in the analysis were collected by the L3 detector [30][31][32][33][34] at an e + e − center-of-mass energy, √ s, of about 91.2 GeV. Calorimeter clusters having energy greater than 0.1 GeV are used to determine the event thrust [35] axis and to classify events as two-or three-jet events, which are analyzed separately, since differences in the space-time structure (and hence in the BEC parameters) could be expected, and indeed have been observed [26,36].
Backgrounds to hadronic Z decays such as leptonic Z decays, beam-wall and beam-gas interactions, and two-photon events are excluded by requiring that the visible energy, E vis , be within 0.5 √ s and 1.5 √ s, that the transverse and longitudinal energy imbalances be less than 0.6E vis and 0.4E vis , respectively, and that the number of calorimeter clusters be at least 15. In order to ensure well-measured charged tracks, events are required to have the thrust direction within the barrel region of the calorimeters and the acceptance region of the tracking chamber: | cos(Θ)| < 0.74, where Θ is the angle between the thrust axis and the beam direction. High precision charged tracks are selected by requiring • transverse momentum greater than 150 MeV, • at least one hit in the inner region of the tracking chamber (TEC), • more than 25 hits in the entire TEC spanning at least 40 wires of the possible 62, • a distance, in the plane transverse to the beam, of closest approach to the interaction vertex less than 10 mm.
Further, tracks lying in two small regions of azimuthal angle having less precise calibration are rejected: 45 • -52.5 • and 225 • -232.5 • . Since the resolution of the opening angle between pairs of tracks is crucial for the study of BEC, tracks are also rejected if there is no hit in the Z-chamber of the TEC. To further reject τ + τ − events, the second largest angle, φ 2 , between any two neighboring tracks in the transverse plane is required to be in the range 20 • -170 • . To further ensure that the event is well contained within the acceptance of the TEC, the thrust axis is determined using only the high precision tracks, and the event is rejected if | cos(Θ TEC )| > 0.7, where Θ TEC is the angle between this thrust axis and the beam. In total about 0.8 million events with an average number of about 12 high precision charged tracks are selected. This results in approximately 36 million like-sign pairs of charged tracks.
The number of jets in an event is determined using the Durham jet algorithm [37][38][39] with a jet resolution parameter y cut = 0.006, yielding about 0.5 million two-jet events and 0.3 million events having more than two jets. Since there are few events with more than three jets, the category of events with more than two jets is referred to as the threejet sample.

Bose-Einstein correlation function
The two-particle correlation function of two particles with four-momenta p 1 and p 2 is given by the ratio of the twoparticle number density, ρ 2 (p 1 , p 2 ), to the product of the two single-particle number densities, ρ 1 (p 1 )ρ 1 (p 2 ). Since we are here interested only in the correlation R 2 due to Bose-Einstein interference, the product of single-particle densities is replaced by ρ 0 (p 1 , p 2 ), the two-particle density that would occur in the absence of Bose-Einstein correlations: This ρ 2 is corrected for detector acceptance and efficiency using Monte Carlo events generated by the JETSET Monte Carlo generator [40] with the so-called BE 0 simulation of BEC [41] to which a full detector simulation has been applied, 1 by multiplying the measured ρ 2 , on a bin-by-bin basis, by the ratio of ρ 2 of the generated events to ρ 2 of the generated events after detector simulation. For the detectorsimulated events, as for the data, all charged tracks are used, whereas for the generator-level events only charged pions are used, since all measured tracks are assumed to be pions in the calculation of Q. Thus the Monte Carlo generator is used to extract the Q distribution for equally charged pion pairs from the observed Q distribution of all equally charged particle pairs. An event mixing technique is used to construct ρ 0 , whereby all tracks of each data event are replaced by tracks from different events having a multiplicity similar to that of the original event. This is accomplished by first rotating each event to a frame whose axes are the thrust, major, and minor [42] directions and then assigning the events to classes based on the track multiplicity. Each multiplicity defines a class except for very low and very high multiplicities, in which case several adjacent multiplicities are grouped into the same class in order to obtain sufficient statistics. The replacement track is randomly chosen from tracks of the same sign in a randomly chosen event either in the same class as the data event or a neighboring or next-to-neighboring class. Ideally one would define the classes based on the total particle (charged plus neutral) multiplicity, since the aim is to use events where the available phase space of the particles is the same. However, the ratio of neutral to charged multiplicities is not constant and both multiplicities are subject to detector efficiency losses. In an attempt to take this into account we therefore also use events of nearby classes.
A correction for detector acceptance and efficiency is applied to ρ 0 in the same way as to ρ 2 . The mixing technique removes all correlations, e.g., resonances and energymomentum conservation, not just Bose-Einstein correlations. Hence, ρ 0 is also corrected for this by a multiplicative factor which is the ratio of the densities of events to mixed events found using events generated by JETSET without BEC simulation.
Including all corrections, R 2 is measured by where data, gen, det, gen-noBE refer, respectively, to the data sample, a generator-level Monte Carlo sample, the same Monte Carlo sample passed through detector simulation and subjected to the same selection procedure as the data, and a generator-level sample of a Monte Carlo generated without BEC simulation. We note that R 2 is unaffected by single-particle acceptances and efficiencies. Being the same for ρ 2 and ρ 0 , these cancel in the ratio. Systematic uncertainties on R 2 are highly correlated point-to-point. Hence, in performing fits to R 2 , only statistical uncertainties will be used in calculating the χ 2 which is minimized.
In principle R 2 should also be corrected for final-state interactions (FSI), both Coulomb and strong, or alternatively, the functional form used to fit R 2 (Q) should be adapted to take account of FSI [43,44]. Coulomb interactions, being repulsive for like-sign charged pions, serve to increase the measured values of Q. This is often taken into account by weighting the pion pairs with the inverse of the socalled Gamow factor [1]. However, this factor strongly overcompensates, since it is derived assuming that the source is point-like and that all particles are pions [45]. Further the strong π -π FSI at short distances are attractive, further reducing the effect [46]. The net effect is small, amounting to about 3% at Q = 0 [44], comparable to the statistical uncertainty in the first bin of R 2 (Q) (cf. Figs. 2-3). Further, if R 2 data is weighted to account for FSI, R 2 gen , R 2 det and R 2 gen-noBE need to be similarly weighted. This is because the Monte Carlo programs are tuned to data, and this tuning presumably accounts, at least partially, for the effects of FSI, even though the Monte Carlo program does not explicitly do so. The result is that the various weights, to a large degree, cancel. We have used unlike-sign charged pions to check that applying Gamow factors is not needed in the present analysis. The uncertainty on FSI is included in the systematic uncertainties by varying the fit range (cf. Sect. 3.4).

Dependence on Q
Since the four-momenta of the two particles are on the massshell, the correlation function is defined in six-dimensional momentum space. Since BEC can be large only at small four-momentum difference Q, they are often parametrized in this one-dimensional distance measure. With a few assumptions [2,5,7], the two-particle correlation function, (1), is related to the Fourier transformed source distribution: where f (x) is the (space-time) density distribution of the source, andf (Q) is the Fourier transform of f (x). The parameter λ is introduced as a measure of the strength of the correlation. It may be less than unity for a variety of reasons such as inclusion of misidentified non-identical particles or the presence of long-lived resonance decays if the particle emission consists of a small, resolvable core and a halo having experimentally unresolvable large length scales [47,48]. The parameter γ and the (1 + Q) term parametrize possible long-range correlations not adequately accounted for in the reference sample. While there is no guarantee that (1 + Q) is the correct form, we will see that it does provide a good description of R 2 in the region Q > 1.5 GeV. In fact, if the reference sample indeed removes only BEC, we expect R 2 to contain no long-range correlations and hence to find = 0. There is no reason, however, to expect the hadron source to be spherically symmetric in jet fragmentation. Recent investigations have, in fact, found an elongation of the source along the jet axis [22,[25][26][27][28]. While this effect is well established, the elongation is actually only about 20%, which suggests that a parametrization in terms of the single variable Q, may be a good approximation. This is not the case in heavy-ion and hadron-hadron interactions, where BEC are found not to depend simply on Q, but on components of the four-momentum difference separately [5,[14][15][16][17]19]. However, in e + e − annihilation at lower energy [13] it has been observed that Q is the appropriate variable. We checked this [49] both for all and for two-jet events: We observe that R 2 does not decrease when both q 2 = ( p 1 − p 2 ) 2 and q 2 0 = (E 1 − E 2 ) 2 are large while Q 2 = q 2 − q 2 0 is small, but is maximal for Q 2 = q 2 − q 2 0 = 0, independent of the individual values of q and q 0 . In a different decomposition, 2 is the component transverse to the thrust axis and Q 2 LE = (p l1 − p l2 ) 2 − (E 1 − E 2 ) 2 combines the longitudinal momentum and energy differences, we find R 2 to be maximal along the line Q = 0, as is shown in Fig. 1, and fits, though of poor χ 2 (confidence level of about 1%), are consistent with equal radii, as has previously been observed by ALEPH [50].
We conclude that a parametrization in terms of Q can be considered a reasonable approximation for the purposes of this article. We shall return to the question of elongation in Sect. 5.

Symmetric parametrizations
The simplest assumption is that the source has a symmetric Gaussian distribution with mean μ = 0 and standard deviation R. In this casef (Q) = exp(iμQ − (RQ) 2 2 ) and However, this parametrization is often found to be inadequate in its description of data.
A model-independent way to study deviations from the Gaussian parametrization is to use [5,51,52] where κ is the third-order cumulant moment and 2RQ is the third-order Hermite polynomial. Note that the second-order cumulant corresponds to the radius R. Equation (5) reduces to (4) for κ = 0.
Another way to depart from the assumption of a Gaussian source is to replace the Gaussian by a symmetric Lévy stable distribution, which is characterized by three parameters: x 0 , R, and α. Its Fourier transform,f (Q), has the following form: where the index of stability, α, satisfies the inequality 0 < α ≤ 2. The case α = 2 corresponds to a Gaussian source distribution with mean x 0 and standard deviation R. For more details, see, e.g., [54]. Then R 2 has the following form [55]: These three parametrizations are fitted to the data for both two-and three-jet events. Fits of the Gaussian parametrization, (4), to the data (shown in Figs. 2a and 3a for twoand three-jet events, respectively) result in unacceptably low confidence levels: 10 −15 (χ 2 = 246 for 96 degrees of freedom) for two-jet and much worse (χ 2 = 456) for three-jet events. The fit is particularly bad at low values of Q, where both for two-and three-jet events R 2 is much steeper than a Gaussian.
A fit of the Edgeworth parametrization, (5), to the two-jet data, shown in Fig. 2b, finds κ = 0.74 ± 0.07, about 10 standard deviations from the Gaussian value of zero. The confidence level is indeed much better than that of the purely Gaussian fit, but is still poor, approximately 10 −5 . Close inspection of the figure shows that the fit curve is systematically above the data in the region 0.6-1.2 GeV and that the data for Q ≥ 1.5 GeV appear flatter than the curve, as is also the case for the purely Gaussian fit. Similar behavior, but with a worse χ 2 , is observed for three-jet events (Fig. 3b).
The fit of the Lévy parametrization, (7), to the two-jet data, shown in Fig. 2c, finds α = 1.44 ± 0.06, far from the Gaussian value of 2. The confidence level, although improved compared to the fit of (4) is still unacceptably low, approximately 10 −8 . Similar behavior, with a worse χ 2 , is seen for three-jet events (Fig. 3c).
Both the symmetric Lévy parametrization and the Edgeworth parametrization do a fair job of describing the region Q < 0.6 GeV, but fail at higher Q. R 2 in the region Q ≥ 1.5 GeV is nearly constant (≈1). However, in the region 0.6-1.5 GeV R 2 has a smaller value, dipping below Fig. 2 The Bose-Einstein correlation function R 2 for two-jet events with the results of fits of (a) the Gaussian, (4), (b) the Edgeworth, (5), and (c) the symmetric Lévy, (7), parametrizations. Also plotted is , the difference between the fit and the data. The dashed line represents the long-range part of the fit, i.e., γ (1 The dotted line represents a fit with λ = 0 to Q >

Fig. 3
The Bose-Einstein correlation function R 2 for three-jet events with the results of fits of (a) the Gaussian, (4), (b) the Edgeworth, (5), and (c) the symmetric Lévy, (7), parametrizations. Also plotted is , the difference between the fit and the data. The dashed line represents the long-range part of the fit, i.e., γ (1 The dotted line represents a fit with λ = 0 to Q > 1.52 GeV unity, 2 which is indicative of an anti-correlation. This is clearly seen in Figs. 2 and 3 by comparing the data in this region to an extrapolation of a linear fit, (4) with λ = 0, in the region Q ≥ 1.5 GeV. The inability to describe this dip in R 2 is the primary reason for the failure of both the Edgeworth and symmetric Lévy parametrizations. Failure to describe this anti-correlation region, i.e., the dip below unity, is an inherent feature of any parametrization of the form (3). We note that this dip is less apparent if one only plots (and fits) R 2 for Q < 2 GeV as has usually been done in the past, the dip being accommodated by a relatively large value of .

Time dependence of the source
The parametrizations discussed so far, which have proved insufficient to describe the BEC, all assume a static source. The parameter R is a constant. It has, however, been observed that R depends on the transverse mass, m t = [20][21][22]. It has been shown [56,57] that this dependence can be understood if the produced pions satisfy, approximately, the (generalized) Bjorken-Gottfried condition [58][59][60][61][62][63], whereby the fourmomentum of a produced particle and the space-time position at which it is produced are linearly related. Such a correlation between space-time and momentum-energy is also a feature of the Lund string model, which, as incorporated in JETSET, is very successful in describing detailed features of the hadronic final states of e + e − annihilation. Not only inclusive single-particle distributions and event shapes, but also correlations [64][65][66] are well described.

The τ -model
In Sect. 3.1 we have seen that BEC depend, at least approximately, only on Q and not on its components separately. This is a non-trivial result. For a hydrodynamical type of source, on the contrary, BEC decrease when any of the relative momentum components is large [5,17]. Further, we have seen that R 2 shows anti-correlations in the region 0.6-1.5 GeV, dipping below its values at higher Q.
A model which predicts such Q-dependence, as well as the absence of dependence on the components of Q separately, is the so-called τ -model [23]. Further it incorporates the Bjorken-Gottfried condition and predicts a specific transverse mass dependence of R 2 , which we subject to an experimental test here.
In this model, it is assumed that the average production point in the overall center-of-mass system, x = (t, r x , r y , r z ), of particles with a given four-momentum p = (E, p x , p y , p z ) is given by 2 More correctly, dipping below the value of γ (1 + ).
In the case of two-jet events, a = 1/m t where m t is the transverse mass and τ = t 2 − r 2 z is the longitudinal proper time. 3 For isotropically distributed particle production, the transverse mass is replaced by the mass in the definition of a and τ by the proper time, t 2 − r 2 x − r 2 y − r 2 z . In the case of three-jet events the relation is more complicated.
The second assumption is that the distribution of x μ (p μ ) about its average, δ (x μ (p μ ) − x μ (p μ )), is narrower than the proper-time distribution, H (τ ). Then the two-particle Bose-Einstein correlation function is indeed found to depend on the invariant relative momentum Q, rather than on its separate components, as well as on the values of a of the two particles [24]: where Since there is no particle production before the onset of the collision, H (τ ) should be a one-sided distribution. In the leading log approximation of QCD the parton shower is a fractal [67][68][69]. Further, a Lévy distribution arises naturally from a fractal [70]. We are thus led to choose a one-sided Lévy distribution for H (τ ) [24]. The characteristic function of H (τ ) can then be written [55] (for α = 1) as 4 where the parameter τ 0 is the proper time of the onset of particle production and τ is a measure of the width of the proper-time distribution. Using this characteristic function in (9), and incorporating the factor λ and the long-range parametrization, yields The terminology 'longitudinal' proper time and 'transverse' mass seems customary in the literature even though their definitions are analogous τ = t 2 − r 2 z and m t = E 2 − p 2 z . 4 Our notation, H (ω), corresponds to Nolan's [54] S(α, β, γ, δ; 1) parametrization with β = 1, γ =⇒ τ/2 1/α , and δ =⇒ τ 0 . The special case of α = 1 is also given by Nolan [54].
Note that the cosine factor generates oscillations corresponding to alternating correlated and anti-correlated regions, a feature clearly seen in the data (Figs. 2 and 3). Note also that since a = 1/m t for two-jet events, the τ -model predicts a decreasing effective source size with increasing m t . These features are subjected to a quantitative test in the following sections.

The τ -model for average a
Fits of (11) are difficult since R 2 depends on three variables: Q, a 1 and a 2 . Further, we have a simple expression for a only for two-jet events. Therefore, we first consider a simplification of (11) obtained by assuming (a) that particle production starts immediately, i.e., τ 0 = 0, and (b) an average a-dependence, which is implemented by introducing an effective radius defined by This results in where R a is related to R by Fits of (13) are first performed with R a as a free parameter. The fit results obtained for two-and three-jet events are listed in Table 1 and shown in Fig. 4 for two-jet events and in Fig. 5 for three-jet events. They have acceptable confidence levels, describing well the dip below unity in the 0.6-1.5 GeV region, as well as the peak at low values of Q.
As shown in Table 1, the estimates of some fit parameters are rather highly correlated. Taking these correlations into account, the fit parameters for the two-jet events satisfy (14), the difference between the left-and right-hand sides of the equation being less than a standard deviation. For three-jet events the difference amounts to about 1.5 standard deviations.
Note that no significant long-range correlation is observed: is zero within 1 standard deviation, and fits with fixed to zero find the same values, within 1 standard deviation, of the parameters as the fits of Table 1. Thus the method to remove non-Bose-Einstein correlations in the reference sample is apparently successful.
Fit results imposing (14) are given in Table 2. For two-jet events, the values of the parameters are comparable to those with R a free. For three-jet events, the imposition of (14) results in values of α and R closer to those for two-jet events. Unlike the fits with R a free, the values of differ somewhat from zero, which could indicate a slight deficiency in the description of BEC in the same way, though on a much smaller  . 4 The Bose-Einstein correlation function R 2 for two-jet events.
The curve corresponds to the fit of the one-sided Lévy parametrization, (13), with the parameter R a (a) free and (b) constrained by (14). The results of the fits are given in Tables 1 and 2, respectively. Also plotted is , the difference between the fit and the data. The dashed line represents the long-range part of the fit, i.e., γ (1 + Q) scale, as seen in the Edgeworth and symmetric Lévy fits of Figs. 2 and 3.

The τ -model with m t dependence
For two-jet events, a = 1/m t , while for three-jet events the situation is more complicated. We therefore limit fits of (11) to the two-jet data. For each bin in Q the average values of m t1 and m t2 are calculated, where m t1 and m t2 are the transverse masses of the two particles making up a pair, requiring m t1 > m t2 . Using these averages, (11) is fit to R 2 (Q). The fit results in τ 0 = 0.00 ± 0.02 fm, and we re-fit with τ 0 fixed to zero. The results are shown in Fig. 6 and Table 3. The parameters α, τ and λ are highly correlated. Equation (11) for two-jet events can be simplified by requiring m t1 ≈ m t2 . This leads to fit results [49] consistent with those without this simplification.
Since the τ -model describes the m t dependence of R 2 , an important test of the model is that its parameters, α, τ , and τ 0 , not depend on m t . Note that the parameter λ, which is not a parameter of the τ -model, can depend on m t , e.g., as a result of resonances [47,48]. The large correlations between the fit estimates of λ and those of the τ -model parameters complicate the testing of m t -independence. We perform fits in various regions of the m t1 -m t2 plane keeping α and τ fixed at the values obtained in the fit to the entire m t plane ( Table 3). The regions are chosen such that the numbers of pairs of particles in the regions are comparable. The m t regions and the confidence levels of the fits are shown in Table 4 along with the values of λ. The confidence levels are seen to be quite reasonably distributed in agreement with the hypothesis of m t -independence of the parameters of the τ -model (α, τ , and τ 0 ).

Systematic uncertainties
The following sources of systematic uncertainty are investigated:  The curve corresponds to the fit of the one-sided Lévy parametrization, (13), with the parameter R a (a) free and (b) constrained by (14). The results of the fits are given in Tables 1 and 2, respectively. Also plotted is , the difference between the fit and the data. The dashed line represents the long-range part of the fit, i.e., γ (1 + Q) Table 2 Results of fits of (13) imposing (14) for two-jet and three-jet events, which are shown in Figs. 4b and 5b . 6 The Bose-Einstein correlation function R 2 for two-jet events. The curve corresponds to the fit of (11). The results of this fit are given in Table 3. Also plotted is , the difference between the fit and the data. The dashed line represents the long-range part of the fit, i.e., γ (1 + Q) Table 3 Results of the fit of (11) for two-jet events, which is shown in Fig. 6. The parameter τ 0 is fixed to zero. MeV; (f) minimum of 21 hits spanning at least 48 wires to 31 hits spanning at least 32 wires; (g) maximum distance of closest approach from 5 mm to 15 mm. Fit range: The fits presented above were performed in the range 0 < Q < 4 GeV. The fits were repeated in the ranges 0.04 < Q < 4 GeV, 0 < Q < 3 GeV and 0 < Q < 5 GeV. Note that uncertainty on the contribution of FSI is included by considering the range 0.04 < Q < 4 GeV, FSI being most important at the smallest values of Q. Mixing: The event mixing method to construct ρ 0 uses tracks from events having similar multiplicity. The definition of "similar" was varied by changing the range of multiplicity classes from which the mixed tracks are chosen from 0 to 8, i.e., from demanding the same multiplicity class as the data event to including classes within ±4 of that of the data event.
Monte Carlo: Events generated by JETSET with no BEC simulation or by HERWIG [71] are used instead of JETSET with BEC for the correction of ρ 2 for detector acceptance, efficiency and non-pion pair background and PYTHIA [72]  For each source the root mean square of the deviations from the values measured in the standard selection is taken as the systematic uncertainty, positive and negative deviations being treated separately. The systematic uncertainties from the four sources are added in quadrature to obtain the total systematic uncertainty. For the physically interesting parameters (α, R, R a , τ ), the largest contributions come from Monte Carlo variation and/or mixing. The other parameters sometimes receive significant contributions from fit range and/or track and event selection.

The emission function of two-jet events
Within the framework of the τ -model, we now reconstruct the space-time picture of the emitting process for two-jet events. The emission function in space-time, S rτ ( r, τ ), normalized to unity, is given in the τ -model by [24]: 5 HERWIG was modified to use the decay and BEC routines of PYTHIA. Fragmentation parameters of both HERWIG and PYTHIA were tuned using L3 data with BEC simulated by the BE 32 Gaussian algorithm [74]. The events are generated using these parameters but without BEC simulation.
where n andn are, respectively, the number and average number of pions produced, and where ρ p ( p) is the singleparticle momentum distribution, ρ p ( p) = d 3 n dp x dp y dp z , which is normalized to the mean multiplicity,n. Note that, given the τ -model correlation between space-time and momentum space, Using τ 2 = t 2 − r 2 z and (17), we can rewrite (15) in terms of r and t: where is the Jacobian of the variable transformation. Given the symmetry of two-jet events, it is convenient to write S in cylindrical coordinates and average over the azimuthal angle. To simplify the reconstruction of S we assume that the momentum distribution ρ p can be factorized in the product of transverse and longitudinal distributions, e.g., ρ yp t = ρ y (y)ρ p t (p t ), where y is the rapidity and p t is the transverse momentum. This is found to be a reasonable approximation except at very high values of |y| or p t [49].
Using H (τ ) as obtained from the fit of (11) ( Table 3), which is shown in Fig. 7, together with the inclusive p t and rapidity distributions, 6 which are shown in Fig. 8, the full emission function is reconstructed.
Integrating (15) and (18) over the transverse coordinates results, respectively, in and These are plotted in Fig. 9. They exhibit a "boomerang" shape with maxima at low values of τ and η or t and r z , but with tails reaching out to very large values, a feature also observed in hadron-proton [76] and heavy-ion collisions [77]. Note however, that in the case of two-jet events, the emission function S ητ (η, τ ) has two maxima, as does S zt (r z , t), as is seen in Fig. 9. These two different maxima correspond to the two-jet structure of the events. This is in marked contrast with the reconstructed emission function in hadron-proton reactions at √ s = 22 GeV, where only a single maximum is observed [76].
The transverse part of the emission function is obtained by integrating (15) over r z and averaging over the azimuthal angle: where .
(23) Figure 10 shows the transverse part of the emission function, (22), for various proper times. Particle production starts immediately, increases rapidly and decreases slowly. A ringlike structure, similar to the expanding, ring-like wave created by a pebble in a pond, is observed. These pictures together form a movie representing this temporal process, which takes place, to our knowledge, at the shortest time scale ever measured. 7 In the movie two ring-like structures separate along the thrust axis, their radii increasing. They thus span the surface of two back-to-back cones whose tips meet at the origin. The rings are very faint at the beginning, then quickly brighten as the particle production probability increases to a sharp peak at τ ≈ 0.035 fm (≈0.12 × 10 −24 s). These separating and expanding rings then start to fade gradually with a characteristic duration of τ ≈ 1.56 fm (≈5.2 × 10 −24 s) following a power-law tail. Note however that the shape of these rings is in marked contrast to the more smoothly edged rings of fire found in hadron-proton reactions [76]. In our case, these rings are reconstructed from the Bose-Einstein correlation functions and from the measured p t and rapidity spectra of two-jet events without any reference to temperature or flow, indicating the non-thermal nature of particle production in e + e − annihilation.

Test of dependence of BEC on components of Q
In this section we return to the question of a possible elongation of the region of homogeneity and the possible dependence of the Bose-Einstein correlation function on components of Q.
The τ -model predicts that the two-particle BEC correlation function R 2 depends on the two-particle momentum   (r z , t), (21), normalized to unity difference only through Q, i.e., that a single radius parameter applies to all components of Q. However, previous studies [22,[25][26][27][28] have found that different radii are required for the different components of Q in order to fit the data, and indeed that the shape of the region of homogeneity is elongated along the event (thrust) axis. The question is whether this is an artifact of the Edgeworth or Gaussian parametrizations used in these studies or shows a limitation of the τ -model.

Test of elongation and the τ -model
In the previous studies of elongation, Q 2 is split into three components, where β = p 1out + p 2out E 1 + E 2 (25) in the Longitudinal Center of Mass System (LCMS) of the pair. The LCMS frame is defined as the frame, obtained by a Lorentz boost along the event axis, where the sum of the three-momenta of the two pions ( p 1 + p 2 ) is perpendicular to the event axis. Assuming azimuthal symmetry about the event axis suggests that the region of homogeneity have an ellipsoidal shape with one axis, the longitudinal axis, along the event axis. The longitudinal axis would then be longer, shorter, or equal to the other two (transverse) axes, which are of equal length. In the LCMS frame the event axis is referred to as the longitudinal direction; the direction of p 1 + p 2 as the out direction; and the direction perpendicular to these as the side direction. The components of | p 1 − p 2 | along these directions are denoted by Q L , Q out , and Q side . In the Gaussian and Edgeworth parametrizations R 2 Q 2 is then replaced in the fitting function by The longitudinal and transverse sizes of the source are measured by R L and R side , respectively, whereas the value of ρ out reflects both the transverse and temporal sizes. 8 Note that the LCMS and the rest frame of the pair differ only by a Lorentz boost of velocity β along the out direction. Hence R L and R side also measure the longitudinal and transverse size in the pair rest frame, while ρ out may be interpreted as 8 In the literature the coefficient of Q 2 out in (26) is usually denoted R 2 out . We prefer to use ρ 2 out to emphasize that, unlike R L and R side , ρ out contains a dependence on β and to differentiate it from R out in (31a) below. an average of the transverse size in the rest frame boosted to the LCMS: ρ 2 out = r 2 out 1 − β 2 . In our previous analysis [25] we found, using all events, R side /R L = 0.80 ± 0.02 +0.03 −0.18 and R side /R L = 0.81 ± 0.02 +0.03 −0.19 for the Gaussian and Edgeworth parametrizations, respectively. Since the present analysis uses a somewhat different event sample and mixing algorithm, we have repeated these fits using the same binning of 0.08 GeV for 0.00 < Q i < 1.04 GeV (i = L, side, out). We find consistent values [49] of R side /R L : 0.76 for both parametrizations. We also find the elongation to be greater for two-jet events than for three-jet events, as has previously been ob-served by OPAL [26]. Using the Edgeworth parametrization the elongation is (statistical uncertainties only) R side /R L = 0.65 ± 0.02 for two-jet events, 0.76 ± 0.02 for all events, and 0.84 ± 0.02 for three-jet events.
We next investigate whether the τ -model parametrizations could accommodate an elongation. This is done by incorporating (26) in (13). If (14) is imposed, this results in where and where we use the same long-range parametrization as in the above Edgeworth and Gaussian parametrizations. We have attempted the long-range parametrization (1+ Q), but that results in unacceptable χ 2 . To extend the range of Q included in the fit, the bin size is increased to 0.16 GeV for Q i > 0.88 GeV. Results for two-jet events of a fit of (27) for Q < 4 GeV are shown in Table 5. A clear preference for elongation (R side /R L < 1) is seen. The value of R found in the corresponding fit of R 2 (Q) ( Table 2) lies between the values of R side and R L found here. Furthermore, the value of the elongation agrees well with that found above for two-jet events using the Edgeworth parametrization.

Direct test of R 2 dependence on components of Q
To directly test the hypothesis of no separate dependence on components of Q, we investigate two decompositions of Q: where Q 2 LE = Q 2 L − ( E) 2 and q 2 out = Q 2 out − ( E) 2 . The decomposition of (29a) corresponds to the LCMS frame where the longitudinal and energy terms are combined; its three components of Q are invariant with respect to Lorentz boosts along the thrust axis. The decomposition of (29b) corresponds to the LCMS frame boosted to the rest frame of the pair; its three components of Q are invariant with respect to Lorentz boosts along the out direction.
The test then consists of replacing R 2 Q 2 by R 2 times one of the above equations and comparing a fit where the coefficients of all three terms are constrained to be equal to a fit where each of these coefficients is a free parameter.
Modifying (13) in this way, and imposing (14), results in where Note that, since E vanishes in the rest frame of the pair, the quantity r out is a direct measure of the spatial extent of the source in this frame whereas ρ out in (28) and R LE in (31a) are also sensitive to the temporal distribution of the source.
The results of the four fits of (30) are shown in Table 6 for Q < 4 GeV. For both parametrizations the fit allowing separate dependence on the components of Q has a higher confidence level than the fit which does not. The fit using (31) attains a confidence level of 2% when no separate dependence is allowed. While this is not poor enough to reject by itself the hypothesis of no elongation, the fit allowing separate dependence describes the data better with a confidence level of 38%. When (32) is used the fit not allowing separate dependence is rejected by a confidence level of 10 −7 , while allowing it results in the acceptable confidence level of 2%. For the fits of both (31) and (32) the differences in χ 2 between the fit allowing separate dependence and the fit not allowing it is huge, 296 and 464, respectively, while a difference of only 2 is expected if the hypothesis of no separate dependence is correct. This provides extremely strong evidence against the hypothesis of identical dependence of R 2 on the different components of Q.
The fits using (32) are compared to the data in Fig. 11 for small values of Q. As is apparent in Fig. 11 the shape of the dependence of R 2 on the different components of Q is different at small values of these components and can not be described by equal values of R L , R side and r out .  Fig. 11 Projections of the Bose-Einstein correlation function R 2 (Q L , Q side , q out ) for two-jet events using the region below 80 MeV for the non-projected components We have repeated the fits varying the upper limit of Q to as low as 1 GeV. While some parameters show a considerable dependence on the Q-range of the fit, R side /R L and r out /R L are found to be stable, varying by less than a standard deviation. The variation in the other parameters is thought to arise from correlations among the parameter values and the inability to determine well the values of the long-range parameters L , side and out when Q is confined to small values. The value of the elongation, R side /R L , is far from unity and agrees well with the value in the fit of (27) imposing (14) ( Table 5) and with that found using the Edgeworth parametrization. The value of r out /R L is also far from unity, but is greater than unity, in contrast to R side /R L . Even more striking is the inequality of R side and r out : r out /R side 2. This indicates the invalidity of the azimuthal symmetry assumed in the elongation analyses.
We also investigate elongation in the parametrization of (11) by replacing τ Q 2 by τ long Q 2 L + τ side Q 2 side + τ out q 2 out . Fits of the resulting parametrization are performed both with τ long , τ side and τ out as free parameters and with them constrained to be equal. In these fits only one m t bin is used, 0 < m t < 4 GeV. The results are shown in Table 7. The conclusions are the same as for the fits of (30). The value of τ side / τ L agrees with that of R side /R L , and τ out / τ long is greater than unity. Constraining the fit by τ L = τ side = τ out results in an unacceptable χ 2 .
We note that in all of the fits of ad hoc τ -model parametrizations the long-range correlation parameters, side and out , deviate substantially and significantly from zero. This could indicate that not all non-BEC correlations are removed from the reference sample, which in turn could influence the amount of elongation found by the fits. To investigate this we have performed fits of the region Q > 1.8 GeV fixing the BEC contribution to zero (λ = 0). The fits for Q < 4 GeV were then performed with the i fixed to the values obtained in the λ = 0 fits. Of course, the confidence lev- Table 7 Results for two-jet events of fits of (11) with τ Q 2 replaced by τ long Q 2 L + τ side Q 2 side + τ out q 2 out and long-range parametrization 1 + L Q L + side Q side + out q out , where L , side and out are free parameters. In the first fit all three τ parameters are free; in the second they are constrained to be equal. els were lower than those of the fits where the i could vary, but the elongation varied only slightly. For the fits of Table 6 R side /R LE = 0.58 ± 0.02 and R out /R LE = 0.984 ± 0.002 with a CL of 8%, while R side /R L = 0.64 ± 0.02 and R out /R L = 1.25 ± 0.03 with a CL of 0.09%. We also repeated the fits with the i fixed to zero. The confidence levels were much worse, but the elongation parameters proved quite robust: R side /R LE = 0.69 ± 0.02 and R out /R LE = 0.99 ± 0.04 with a CL of 10 −11 , and R side /R L = 0.72 ± 0.02 and R out /R L = 1.24 ± 0.03 with a CL of 10 −13 . Further, note in Fig. 11 that the data themselves show different dependencies on Q L , Q side , q out . We conclude that the elongation is real, i.e., not an artifact of the parametrizations previously used, and that the ad hoc modified τ -model provides a reasonable description of the elongation.

Discussion and conclusions
The usual parametrizations of BEC of pion pairs in terms of their four-momentum difference Q, such as the Gaussian, Edgeworth, or symmetric Lévy parametrizations, are found to be incapable of describing the Bose-Einstein correlation function R 2 , particularly the anti-correlation region, approximately 0.5 < Q < 1.5 GeV. This failure has not been previously so obvious, since previous analyses generally fit R 2 (Q) only up to 2 GeV or less and the anti-correlation region is then tacitly absorbed into the parametrization of long-range correlations. The existence of this anti-correlation region has the more general implication of ruling out the form usually proposed for BEC, Those parametrizations of BEC are based on a static view of pion emission. No time dependence is assumed, which means that either the pion emission volume is unchanging during pion emission or that the parameters describing the volume, which result from fitting R 2 , are some sort of time average. The τ -model assumes a very high degree of correlation between momentum space and space-time and introduces a time dependence explicitly. The Bose-Einstein correlation function R 2 measures the real part of products of Fouriertransformed proper-time distributions. Not only the positive correlations at low values Q but also negative correlations at somewhat higher values of Q are predicted. It is as though the Bose-Einstein symmetry pulls identical pions below some critical value of Q closer together creating an excess of pairs at lower Q and leaving a deficit of pairs around this value of Q.
In the τ -model R 2 is then a function of one two-particle variable, the invariant four-momentum difference of the two particles Q, and of two single-particle variables, the values of a of the two particles, where a is a parameter in the correlation between momentum space and space-time, (2). The time-dependence is assumed to be given by an asymmetric Lévy distribution, H (τ ), which imposes causality by allowing particle production only for τ > 0. The parameter a can be absorbed into an effective radius to obtain an expression, (13), for R 2 (Q) which successfully describes BEC, including the anti-correlation region, for both two-and three-jet events. We note that a similar anti-correlation region has recently been observed by the CMS Collaboration [10] in pp collisions at √ s = 0.9 and 7 TeV and successfully fit by the τ -model formula, (13).
For two-jet events a = 1/m t and the introduction of an effective radius is unnecessary. The BEC correlation function is then a function not only of Q but also of the transverse masses of the pions. This description of R 2 (Q, m t1 , m t2 ) also successfully describes the data.
Nevertheless, the τ -model description breaks down when confronted with data expressed in components of Q. The τ -model predicts that the only two-particle variable entering R 2 is Q, i.e., that there is a single radius parameter which is the same for each of the separate components of Q. This implies a spherical shape of the pion region of homogeneity, whereas previous analyses using static parametrizations have found a shape elongated along the event axis, as evidenced by R side being smaller than R L in the LCMS. Assuming azimuthal symmetry about the event axis, R side is the transverse radius of the region of homogeneity.
Accordingly, the τ -model equations for R 2 have been modified ad hoc to allow an elongation. Moreover, fits have been performed not only in the LCMS, but also in the rest frame of the pair. Fits not allowing elongation have much worse confidence levels than fits allowing elongation. The elongation found agrees with that found in the previous analyses. However, in the rest frame of the pair the radius parameter in the out direction is found to be approximately twice that in the side direction. Thus the assumption of azimuthal symmetry about the event axis is found not to hold.
We note that an elongation is expected [79] for a hadronizing string in the Lund model, where the transverse and longitudinal components of a particle's momentum arise from different mechanisms. Perhaps the absence of azimuthal symmetry arises from gluon radiation, which is represented by kinks in the Lund string. Therefore a possible improvement of the τ -model could be to allow the longitudinal and transverse components to have different degrees of correlation between a particle's momentum and its spacetime production point in (8).
The τ -model, as we have used it, incorporates a Lévy distribution. A Lévy distribution arises naturally from a fractal, or from a random walk or anomalous diffusion [70], and the parton shower of the leading log approximation of QCD is a fractal [67][68][69]. In this case, the Lévy index of stability is related to the strong coupling constant, α s , by [80,81] Assuming (generalized) local parton hadron duality [82][83][84], one can expect that the distribution of hadrons retains the features of the gluon distribution. For the value of α found in the fit of (11) for two-jet events (Table 3) we find α s = 0.46 ± 0.01 +0.04 −0.02 . This is a reasonable value for a scale of the order of 1 GeV, which is the scale at which the production of hadrons is thought to take place. For comparison, from τ decay, α s (m τ ≈ 1.8 GeV) = 0.34 ± 0.03 [85], and from the average value of α s (M Z ), α s (1 GeV) = 0.50 +0.06 −0.05 [86].
It is of particular interest to point out the m t dependence of the "width" of the source. In (11) the parameter associated with the width is τ . Note that it enters (11) as τ Q 2 /m t . In a Gaussian parametrization the radius R enters the parametrization as R 2 Q 2 . Our observation that τ is independent of m t thus corresponds to R ∝ 1/ √ m t and can be interpreted as a confirmation of the observation [20][21][22] of such a dependence of the Gaussian radii in 2and 3-dimensional analyses of Z decays. The lack of dependence of all the parameters of (11) on the transverse mass is in accordance with the τ -model.
Given these successes, the BEC fit results of the τ -model have been used to reconstruct the pion emission function of two-jet events. Note that all charged tracks are treated as pions. For BEC the main influence of this is a lowering of the value of the parameter λ. However, the influence of this as well as the influence of resonances on the rapidity and p t distributions may be more profound. Hence, the emission function is subject to some quantitative uncertainty. Nevertheless, the general features of the emission function should remain valid: Particle production begins immediately after collision, increases rapidly and then decreases slowly, occurring predominantly close to the light cone. In the transverse plane a ring-like structure expands outwards, which is similar to the picture in hadron-hadron interactions but unlike that of heavy-ion collisions [5]. Despite this similarity the physical process is much different. Reflecting the nonthermal nature of e + e − annihilation, the proper-time distribution and the space-time structure are reconstructed here without any reference to a temperature, in contrast to the results of earlier hadron-hadron and heavy-ion collisions.