Improved Limits on the Lepton-Flavor Violating Decays tau- -->l-l+l-

A search for the neutrinoless, lepton-flavor violating decay of the tau lepton into three charged leptons has been performed using 376 fb^{-1} of data collected at an e^+e^- center-of-mass energy around 10.58 GeV with the BaBar detector at the PEP-II storage rings. In all six decay modes considered, the numbers of events found in data are compatible with the background expectations. Upper limits on the branching fractions are set in the range (4-8) x 10^{-8} at 90% confidence level.

considered, the numbers of events found in data are compatible with the background expectations. Upper limits on the branching fractions are set in the range (4 − 8) × 10 −8 at 90% confidence level.
In tau decays, the most stringent limit on LFV is B(τ → µγ) < 4.5 × 10 −8 at 90% CL [5]. Many descriptions of physics beyond the Standard Model (SM), particularly models seeking to describe neutrino mixing, predict enhanced LFV in tau decays over muon decays with branching fractions from 10 −10 up to the current experimental limits [6,7,8]. An observation of LFV in tau decays would be a clear signature of non-SM physics, while improved limits will provide further constraints on theoretical models. This paper presents a search for LFV in the neutrinoless decay τ − → ℓ − ℓ + ℓ − , where ℓ is an electron or muon. All possible lepton combinations consistent with charge conservation are considered, leading to six distinct de- . The analysis is based on data recorded by the BABAR detector at the PEP-II asymmetric-energy e + e − storage rings operated at the Stanford Linear Accelerator Center. The data sample consists of 339 fb −1 recorded at √ s = 10.58 GeV, and 37 fb −1 recorded at √ s = 10.54 GeV. With an expected cross section for tau pairs at the luminosity-weighted √ s of σ τ τ = 0.919 ± 0.003 nb [10], this data sample contains about 690 million tau decays. The BABAR detector is described in detail in Ref. [11]. Charged-particle (track) momenta are measured with a 5-layer double-sided silicon vertex tracker and a 40layer helium-isobutane drift chamber inside a 1.5-T superconducting solenoidal magnet. The transverse momentum resolution is parameterized as σ p T /p T = (0.13 · p T /[ GeV/c] + 0.45)%. An electromagnetic calorimeter consisting of 6580 CsI(Tl) crystals is used to identify electrons and photons, a ring-imaging Cherenkov detector is used to identify charged hadrons, and the instrumented magnetic flux return (IFR), embedded with limited streamer tubes and resistive plate chambers, is used to identify muons.
A Monte Carlo (MC) simulation of lepton-flavor violating tau decays is used to optimize the parameter space for the search. Simulated tau-pair events including higherorder radiative corrections are generated using KK2f [12] with one tau decaying to three leptons with a 3-body phase space distribution, while the other tau decays according to measured rates [13] simulated with Tauola [14]. Final state radiative effects are simulated for all decays using Photos [15]. The detector response is simulated with GEANT4 [16], and the simulated events are then reconstructed in the same manner as data.
The signature of the decay τ − → ℓ − ℓ + ℓ − is a set of three charged particles, each identified as either an electron or muon, with an invariant mass and energy equal to that of the parent tau lepton. Candidate signal events in this analysis are required to have a "1-3 topology," where one tau decay yields three charged particles, while the second tau decay yields one charged particle. Events with four well-reconstructed tracks and zero net charge are selected, and the tracks are required to point toward a common region consistent with τ + τ − production and decay. The polar angle of all four tracks in the laboratory frame is required to be within the calorimeter acceptance range. Pairs of oppositely-charged tracks are ignored if their invariant mass, assuming electron mass hypotheses, is less than 30 MeV/c 2 , as these tracks are likely to be from photon conversions in the traversed material. The event is divided into hemispheres in the e + e − centerof-mass (c.m.) frame using the plane perpendicular to the thrust axis, as calculated from the observed tracks and neutral energy deposits. The signal hemisphere must contain exactly three tracks while the other hemisphere must contain exactly one.
Each of the charged particles found in the signal hemisphere must be identified as either an electron or muon candidate. Electrons are identified using the ratio of calorimeter energy to track momentum (E/p), the ionization loss in the tracking system (dE/dx), and the shape of the shower in the calorimeter. Muon identification makes use of a neural net, inputs to which include the number of hits in the IFR, the number of interaction lengths traversed, and the energy deposition in the calorimeter. Muons with momentum less than 500 MeV/c do not penetrate far enough into the IFR to be identified. For the lepton momentum spectrum predicted by the signal MC, the electron and muon identification requirements are found to have an average efficiency per lepton of 91% and 65%, respectively. The probability for a pion to be misidentified as an electron in 3-prong tau decays is 2.7%, while the probability to be misidentified as a muon is 2.9%.
The particle identification (PID) requirements are not sufficient to suppress certain backgrounds, particularly those from light quark pair production and higher-order radiative Bhabha and µ + µ − events that can have four leptons in the final state. To reduce these backgrounds, additional selection criteria are applied to the six differ-ent decay modes. For all decay modes, the momentum of the 1-prong track is required to be less than 4.8 GeV/c in the c.m. frame. Additionally, the track in the 1-prong hemisphere is assigned the most-likely mass hypothesis, and the mass of the 1-prong hemisphere is calculated from the four-momentum of that track and the missing momentum in the event. This mass is required to be in the range 0.3 − 3.0 GeV/c 2 for all channels except e − e + e − and µ − e + e − , for which the mass is required to be in the range 0.5 − 2.5 GeV/c 2 . For the e − e + e − and µ − e + e − decay modes, radiative Bhabha events are further suppressed by rejecting events with pairs of oppositelycharged electron tracks in the 3-prong hemisphere with invariant mass less than 250 MeV/c 2 . For the e − e + e − and e − µ + µ − decay modes, the charged particle in the 1-prong hemisphere is required to deposit energy in the calorimeter, and must not be identified as an electron, while for the µ − e + e − and µ − µ + µ − decay modes this track must not be identified as a muon. For the e − e + e − and e − µ + µ − decay modes, the net transverse momentum of the four tracks must be greater than 400 MeV/c, while for the µ − e + e − mode it must be greater than 200 MeV/c. Events in all six decay modes are required to have no track in the 3-prong hemisphere that is consistent with being a kaon.
To reduce backgrounds further, candidate signal events are required to have an invariant mass and total energy in the 3-prong hemisphere consistent with a parent tau lepton. These quantities are calculated from the observed track momenta assuming lepton masses that correspond to the specific decay mode. The energy difference is defined as ∆E ≡ E ⋆ rec − E ⋆ beam , where E ⋆ rec is the total energy of the tracks observed in the 3-prong hemisphere and E ⋆ beam is the beam energy, with both quantities measured in the c.m. frame. The mass difference is defined as ∆M ≡ M rec − m τ where M rec is the reconstructed invariant mass of the three tracks and m τ = 1.777 GeV/c 2 is the tau mass [13].
The signal distributions in the (∆M, ∆E) plane (see Fig. 1) are broadened by detector resolution and radiative effects. In all decay modes, the radiation of photons from the incoming e + e − particles and from the outgoing tau particles leads to a tail at low values of ∆E. Radiation from the final-state leptons, which is more likely for electrons than muons, produces a tail at low values of ∆M as well.  Fig. 1 shows the observed data in the (∆M, ∆E) plane, along with the signal region boundaries and the expected signal distributions. To avoid bias, a blinded analysis procedure was followed with the number of data events in the signal region remaining unknown until the selection criteria were finalized and all cross checks were performed.
There are three main classes of background remaining after the selection criteria are applied: low multiplicity qq events (mainly continuum light-quark production); QED events (Bhabha and µ + µ − ); and SM τ + τ − events. These three background classes have distinctive distributions in the (∆E, ∆M ) plane. The qq events tend to populate the plane uniformly, while QED backgrounds are restricted to a narrow band at positive values of ∆E, and τ + τ − backgrounds are restricted to negative values of both ∆E and ∆M . A negligible two-photon background remains.
The expected background rates for each decay mode are determined by fitting a set of probability density functions (PDFs) to the observed data in the grand sideband (GS) region of the (∆E, ∆M ) plane. The GS region, shown in Fig. 1, lies between −600 and 400 MeV/c 2 in ∆M and −700 and 400 MeV in ∆E, excluding the signal region. For the qq background, a PDF is constructed from the product of two PDFs P M ′ and P E ′ , where The (∆M ′ , ∆E ′ ) axes have been slightly rotated from (∆M, ∆E) to take into account the observed correlation between ∆E and ∆M for the distribution. The resulting PDF has a total of eight fit parameters, including the rotation angle, all of which are determined by fits to MC qq background samples for each decay mode. For the τ + τ − background PDF, the function P M ′′ (∆M ′′ ) is the sum of two Gaussians with common mean, while the functional form of P E ′′ (∆E ′′ ) is the same as that for the qq PDF. To properly model the wedge-shaped distribution due to the kinematic limit in tau decays, a coordinate transformation of the form ∆M ′′ = cosβ 1 ∆M + sinβ 1 ∆E and ∆E ′′ = cosβ 2 ∆E −sinβ 2 ∆M is performed. In total there are 11 free parameters describing this PDF, and all are determined by fits to the MC τ + τ − samples.
For the three decay channels in which there is a significant QED background, an analytic PDF is constructed from the product of a Crystal Ball function [17] in ∆E ′ and a third-order polynomial in ∆M ′ , where again the (∆M ′ , ∆E ′ ) axes have been rotated slightly from (∆E, ∆M ) to fit the observed distribution. The six parameters of this PDF, including the rotation angle, are obtained by fitting data control samples that are enhanced in Bhabha or µ + µ − events.
With the shapes of the three background PDFs determined, an unbinned maximum likelihood fit to the data in the GS region is used to find the expected background rate in the signal region, shown in Table I  e − e + e − 8.9 ± 0.2 1.33 ± 0.25 4.9 1 4.3 µ − e + e − 8.3 ± 0.6 0.89 ± 0.27 5.0 2 8.0 µ + e − e − 12.4 ± 0.8 0.30 ± 0.55 2.7 2 5.8 e + µ − µ − 8.8 ± 0.8 0.54 ± 0.21 4.6 1 5.6 e − µ + µ − 6.2 ± 0.5 0.81 ± 0.31 6.6 0 3.7 µ − µ + µ − 5.5 ± 0.7 0.33 ± 0.19 6.7 0 5.3 shape determinations and background fits are performed separately for each of the six decay modes. The efficiency of the selection for signal events is estimated with a MC simulation of lepton-flavor violating tau decays. About 40% of the MC signal events pass the 1-3 topology requirement. The total efficiency for signal events to be found in the signal region is shown in Table I for each decay mode and ranges from 5.5% to 12.4%. This efficiency includes the 85% branching fraction for 1prong tau decays.
The PID efficiencies and misidentification probabilities have been measured with control samples both for data and for MC events, as a function of particle momentum, polar angle, and azimuthal angle in the laboratory frame. The systematic uncertainties related to PID have been estimated from the statistical uncertainties of the efficiency measurements and from their discrepancies between data and Monte Carlo, and range from 2.3% for e − e + e − to 12.5% for µ − µ + µ − [18]. The modeling of the tracking efficiency contributes an additional 1% uncertainty. All other sources of uncertainty in the signal efficiency are found to be small, including the statistical limitation of the MC signal samples, modeling in the generator of radiative effects, track momentum resolution, trigger performance, observables used in the selection criteria, and knowledge of the tau 1-prong branching fractions. The signal efficiency has been estimated using a 3-body phase space model and no additional uncertainty is assigned for possible model dependence. Despite the effect of slow muons on the PID efficiency, the total selection efficiency is found to be uniform within 20% across 95% of the phase space for the three leptons.
Since the background levels are extracted directly from the data, systematic uncertainties on the background estimation are directly related to the background parameterization and the fit technique used. The finite data available in the GS region to determine the background rates contributes a significant uncertainty in all decay channels. Uncertainties related to the background PDFs are estimated by varying the background shape parameters within their errors and repeating the fits, and from changing the functional form of the PDFs. The total uncertainties on the background estimates are shown in Table I. Cross checks of the background estimation are performed by considering the number of events expected and observed in sideband regions immediately neighboring the signal region for each decay mode.
The numbers of events observed (N obs ) and the background expectations (N bgd ) are shown in Table I, with no significant excess found in any decay mode. Upper limits on the branching fractions are calculated according to B 90 UL = N 90 UL /(2ε Lσ τ τ ), where N 90 UL is the 90% CL upper limit for the number of signal events when N obs events are observed with N bgd background events expected. The values ε, L, and σ τ τ are the selection efficiency, luminosity, and τ + τ − cross section, respectively. The uncertainty on the product Lσ τ τ is 1.0%. The branching fraction upper limits are calculated including all uncertainties using the technique of Cousins and Highland [19] following the implementation of Barlow [20]. The sensitivity or expected upper limit UL exp 90 , defined as the mean upper limit expected in the background-only hypothesis, is included in Table I. The 90% CL upper limits on the τ − → ℓ − ℓ + ℓ − branching fractions are in the range (4 − 8) × 10 −8 . These limits represent up to an order of magnitude improvement over the previous experimental bounds [21,22].