Measurement of $\psi(3686)\to\Lambda\bar{\Lambda}\eta$ and $\psi(3686)\to\Lambda\bar{\Lambda}\pi^0$ decays

Based on a sample of $448.1\times10^6\ \psi(3686)$ events collected with the BESIII detector, a study of $\psi(3686)\to\Lambda\bar{\Lambda}\pi^0$ and $\psi(3686)\to\Lambda\bar{\Lambda}\eta$ is performed. Evidence of the isospin-violating decay $\psi(3686)\to\Lambda\bar{\Lambda}\pi^0$ is found for the first time with a statistical significance of $3.7\sigma$, the branching fraction $\mathcal{B}(\psi(3686)\to\Lambda\bar{\Lambda}\pi^0)$ is measured to be $(1.42\pm0.39\pm0.59)\times10^{-6}$, and its corresponding upper limit is determined to be $2.47\times10^{-6}$ at 90\% confidence level. A partial wave analysis of $\psi(3686)\to\Lambda\bar{\Lambda}\eta$ shows that the peak around $\Lambda\eta$ invariant mass threshold favors a $\Lambda^*$ resonance with mass and width in agreement with the $\Lambda(1670)$. The branching fraction of the $\psi(3686)\to\Lambda\bar{\Lambda}\eta$ is measured to be $(2.34\pm0.18\pm0.52)\times10^{-5}$. The first uncertainties are statistical and the second are systematic.

In 2012, another sample of ψ(3686) events was collected at the BESIII detector [2]. The total data set of 448 million ψ(3686) events, corresponding to a four-fold increase of 2009 data, allows for an in-depth investigation on the decays ψ(3686) → ΛΛπ 0 and ψ(3686) → ΛΛη, and searches for intermediate Λ states in the Λπ 0 and Λη mass spectra. In addition, these branching fractions may also be used to test the "12%" rule [3][4][5], which predicts that the ratio of branching fractions of ψ(3686) and J/ψ decays into the same light hadron final states is around 12%. In this paper, the charge-conjugate process is always implied unless explicitly mentioned.

II. DETECTOR AND MONTE CARLO SAMPLES
The BESIII detector [6] records symmetric e + e − collisions provided by the BEPCII storage ring [7], in the center-of-mass energy range from 2.0 GeV to 4.95 GeV, with a peak luminosity of 1 × 10 33 cm −2 s −1 achieved at √ s = 3.77 GeV. BESIII has collected large data samples in this energy region [8]. The cylindrical core of the BE-SIII detector covers 93% of the full solid angle and consists of a helium-based multilayer drift chamber (MDC), a plastic scintillator time-of-flight system (TOF), and a CsI(Tl) electromagnetic calorimeter (EMC), which are all enclosed in a superconducting solenoidal magnet providing a 1.0 T magnetic field. The solenoid is supported by an octagonal flux-return yoke equipped with resistive plate counter muon identification modules interleaved with steel. The charged-particle momentum resolution at 1 GeV/c is 0.5%, and the resolution of the specific energy loss dE/dx in the MDC is 6% for electrons from Bhabha scattering. The EMC measures photon energies with a resolution of 2.5% (5%) at 1 GeV in the barrel (end cap) region. The time resolution in the TOF barrel region is 68 ps, while that in the end cap region is 110 ps.
Simulated data samples produced with a GEANT4based [9] Monte Carlo (MC) package, which includes the geometric description of the BESIII detector and the detector response, are used to optimize the event selection criteria, determine detection efficiencies and estimate backgrounds. The simulation models the beam energy spread and initial state radiation (ISR) in e + e − annihilations with the generator kkmc [10,11]. The inclusive MC sample aims to include all possible processes involving the production of the J/ψ and ψ(3686) resonances, and the continuum processes incorporated in kkmc [10,11]. The known decay modes are modeled with evtgen [12,13] using branching fractions taken from the Particle Data Group (PDG) [14], and the remaining unknown charmonium decays are modeled with lundcharm [15,16].
Final state radiation (FSR) from charged particles is incorporated using photos [17]. Signal MC samples of ψ(3686) → ΛΛπ 0 decays are generated with uniform phase space (PHSP), while ψ(3686) → ΛΛη decays are generated according to the results of the partial wave analysis reported later in this paper.

III. EVENT SELECTION
The processes ψ(3686) → ΛΛπ 0 and ψ(3686) → ΛΛη are reconstructed with Λ → pπ − ,Λ →pπ + , π 0 → γγ and η → γγ. Since the final state for both channels is ppπ + π − γγ, the number of charged tracks is required to be four with net charge zero. Each track must satisfy |cos θ| < 0.93, where θ is the polar angle of the track measured by the MDC with respect to the direction of the positron beam.
Each of the photon candidates is required to have an energy deposit in the EMC of at least 25 MeV in the barrel (|cos θ| < 0.80) or 50 MeV in the end caps (0.86 < |cos θ| < 0.92). To eliminate showers from charged tracks, the angle between the position of each shower in the EMC and any charged track must be greater than 10 degrees. To suppress electronic noise and showers unrelated to the event, the EMC time difference from the event start time is required to be within [0, 700] ns. At least two photon candidates are required.
The Λ andΛ candidates are reconstructed by combining pairs of oppositely charged tracks with pion and proton mass hypotheses, fulfilling a secondary vertex constraint [18]. Events with at least one pπ − (Λ) and onē pπ + (Λ) candidate are selected. In the case of multiple ΛΛ pair candidates, the one with the minimum value of χ 2 svtx (Λ) + χ 2 svtx (Λ) is chosen, where χ 2 svtx (Λ) and χ 2 svtx (Λ) are the fit qualities of the secondary vertex fits for Λ andΛ, respectively. To improve the momentum and energy resolution and to reduce background contributions, a four-constraint (4C) energy-momentum conservation kinematic fit is applied to the event candidates under the hypothesis of ΛΛγγ (i.e., not considering the γγ mass), and the corresponding χ 2 4C is required to be less than 40. For events with more than two photon candidates, the combination with the best fit quality is selected from all possible combinations. To reject possible background contributions from ψ(3686) → ΛΛγ and ψ(3686) → ΛΛγγγ, we further require that the χ 2 of the 4C fit for the ψ(3686) → ΛΛγγ assignment is smaller than those of ΛΛγ and ΛΛγγγ. The final pπ − andpπ + mass distributions in two progresses are shown in Figs. 1 (a) and (b), respectively, where clear Λ andΛ signals are visible.
In the case of ψ(3686) → ΛΛπ 0 , additional requirements are applied to further reduce the contamination from background. The χ 2 4C is required to be less than 15, to further suppress background events with one or more than two photons in the final states. The veto cut M(ΛΛ) < 3.4 GeV/c 2 is applied in order to suppress background from ψ(3686) → Σ 0Σ0 . Two other veto cuts, M(pπ 0 ,pπ 0 ) < 1.17 GeV/c 2 and M(pπ 0 ,pπ 0 ) > 1.2 GeV/c 2 , are applied in order to suppress background from ψ(3686) → ΛΣ − π + . The invariant masses M(γ low Λ) and M(γ lowΛ ) are both required to be outside of (1.183, 1.203) GeV/c 2 to suppress the ψ(3686) → ΛΣ 0 π 0 background, where γ low represents the less energetic candidate photon.

IV. BACKGROUND STUDY
To investigate the possible background contributions, the same selection criteria are applied to an inclusive MC sample of 506 million ψ(3686) events. A topological analysis of the surviving events is performed with the generic tool TopoAna [19], and the results indicate that the background peaking at the π 0 invariant mass mainly comes from ψ(3686) → ΛΣ 0 π 0 , while the other background sources present a flat distribution. Thus, a PHSP MC sample of ψ(3686) → ΛΣ 0 π 0 + c.c. is generated, giving a background estimate of 20.4 ± 1.9 events, by using a branching fraction of (1.54±0.04±0.13)×10 −4 obtained from B(ψ(3686) → ΛΣ − π + ) with isospin symmetry considerations.
To estimate the background from e + e − continuum processes, the same procedure is performed on data taken at √ s = 3.773 GeV, with an integrated luminosity of 2.92 fb −1 [20]. The background events are extracted by fitting the M γγ mass distribution, normalized to the ψ(3686) data taking into account the luminosity and energy-dependent cross section of the quantum electrodynamics (QED) processes. The normalization factor f is calculated as where N , L, σ, and ǫ refer to the number of observed events, integrated luminosity of data, cross section, and detection efficiency at the two center of mass energies, respectively. The details on the cross section values can be found in Ref. [21]. The detection efficiency ratio ǫ ψ(3686) /ǫ ψ(3770) has been determined by Monte Carlo simulations. After normalization, the background contributions from e + e − → ΛΛπ 0 and e + e − → ΛΛη at 3.686 GeV are determined to be 13.2 ± 1.7 and 19.1 ± 2.0 events, respectively.
Due to the identical event topology, these background events are indistinguishable from signal events and are subtracted directly by fixing their magnitudes in the fit when extracting signal yields. Here we assume that the interference between ψ(3686) decay and continuum process is negligible.
The ψ(3686) → ΛΛπ 0 signal yield is obtained from an extended unbinned maximum likelihood fit to the γγ invariant mass distribution. The total probability density function consists of a signal and various background contributions. The signal component is modeled with the MC simulated signal shape convolved with a Gaussian function to account for a possible difference in the mass resolution between data and MC simulation. The background events from e + e − annihilations are described by the shape obtained from the data taken at √ s = 3.773 GeV, while the peaking background from ψ(3686) → ΛΣ 0 π 0 is modeled with the MC simulation shape. In the fit, the contributions of these two background sources are fixed to the values discussed above. In addition, the nonpeaking background is parameterized by a first order Chebychev function. From the fit, shown in Fig. 2, we estimate 23.0 ± 6.3 ΛΛπ 0 events with a statistical significance of 3.7σ which is evaluated by comparing the likelihood values with and without the π 0 signal included in the fit. The detection efficiency obtained from MC simulation events is 9.0% and these results are summarized in Table III.

VI. ANALYSIS OF ψ(3686) → ΛΛη
The distribution of M (γγ) in the η mass region is shown in Fig. 3. A fit to the η signal with the MC simulated signal shape convolved with a Gaussian function is performed, and the background contribution is described by the shape obtained from the continuum data plus a first order Chebychev function. The fitting results are shown in Fig. 3, with a total of 218 ± 17 ΛΛη signal  Using the Feynman diagram calculation package [22], a partial wave analysis (PWA) is performed based on an unbinned maximum likelihood fit. In the global fit, resonances are described by a relativistic Breit-Wigner propagator, with the mass and width as free parameters, where s is the squared invariant mass.
To describe the Λη andΛη mass spectra, all kinematically-allowed resonances of Λ * and Σ * listed in the PDG [14] are considered. Only components with a statistical significance larger than 5σ are kept in the baseline solution. PWA results indicate that the Λ(1670) plus the nonresonant contribution could provide a good description of data, as illustrated in Figs. 5 and 6. The fitted mass and width of Λ(1670), (1672±5) MeV/c 2 and (38 ± 10) MeV, are also in agreement with the world average values; a total of 116 ± 28 ψ(3686) → Λ(1670)Λ candidate events are measured (based on the PWA amplitude "fit fraction"), and the detection efficiency is determined to be 12.5% by using a PWA-weighted MC sample. The measured yield and detection efficiency are summarized in Table III. The hypothesis of a Λ(1690) state instead of Λ(1670) in the model has been tested, leading to a reasonable description of the data, but with a sightly worse fit quality and with resonance parameters not consistent with the PDG values; it has thus been rejected.
To obtain the detection efficiency of ψ(3686) → ΛΛη, a MC sample is generated in accordance with the above PWA results. The corrected detection efficiency, 12.9%, and the number of signal events, 218 ± 17 are presented in Table III.

VII. SYSTEMATIC UNCERTAINTIES
In this analysis, the systematic uncertainties on the branching fractions mainly come from the following sources: The efficiency of Λ(Λ) reconstruction is studied using the control sample of ψ(3686) → ΛΛ decays, and a correction factor of 0.980 ± 0.011 [23] is applied to the efficiencies obtained from MC simulation. The uncertainty of the correction factor, 1.1%, which includes the uncertainties of MDC tracking and Λ(Λ) reconstruction, is considered as the uncertainty of the efficiency of Λ(Λ) reconstruction.
• Photon detection The photon detection efficiency has been studied using a high-purity control sample of J/ψ → ρ 0 π 0 [24]. The difference between the detection efficiencies of data and MC is around 1% per photon. Thus, 2% is assigned as the total systematic uncertainty for the detection of the two photons.
• Kinematic fit The uncertainty associated with the 4C kinematic fit comes from the inconsistency between data and MC simulation in the fit. This difference is reduced by correcting the track helix parameters of the MC simulation, with parameters from [25] and [26]. Following the method described in Ref. [27], we obtain the systematic uncertainties of the 4C kinematic fit as 3.8% and 1.8% for ψ(3686) → ΛΛπ 0 and ψ(3686) → ΛΛη, respectively.
• Mass window requirements The systematic uncertainties related to each individual mass window requirement are estimated by varying the size of the mass window by one standard deviation of the corresponding mass resolution. For the mass window of M(ΛΛ) < 3.4 GeV/c 2 , the uncertainty is estimated by decreasing the required mass threshold by 10 MeV/c 2 . The largest variation of branching fraction for each mass requirement is considered as the related systematic uncertainty.
• Signal shape In order to estimate the systematic uncertainty due to the signal shape, alternative fits are performed to determine the yields of signal events; the MC shape is replaced with a Breit-Wigner function convolved with a Gaussian function or a single Gaussian function, by varying the fits of the invariant mass distributions by either contracting, expanding or shifting the fit range by ±10 MeV. The maximum differences with the nominal results are assigned as the corresponding systematic uncertainties.
• Background uncertainty To estimate the uncertainty of the nonpeaking background shape in the fit to M(γγ), we performed alternative fits by replacing the first-order Chebychev function with a second-order Chebychev function for ψ(3686) data. The maximum changes of 2.0% and 3.5% are considered as systematic uncertainties. The uncertainties of background from continuum events and the decay ψ(3686) →ΛΣ 0 π 0 are propagated from the statistical uncertainties quoted in Sec. IV.
• Interference between ψ(3686) and continuum amplitudes To estimate the effect from interference of the continuum amplitude with the resonance amplitude, we use the method from Ref. [28]. The maximum impact from interference term with respect to the resonance term is defined as r max where c is the conversion constant, σ f c (s) is the cross section of the continuum process measured from data, B f is the branching fraction of ψ(3686) → ΛΛπ 0 and ψ(3686) → ΛΛη that we measured in this paper and the factor B is constant depending on the resonance parameters quoted from Ref. [28]. The r max R , 40.3% and 20.6% of ψ(3686) → ΛΛπ 0 and ψ(3686) → ΛΛη, are taken as the uncertainty of interference, respectively. Since the Λ(1670)Λ cannot be studied well in continuum with our current statistics, no systematic is provided for this final state.

• Physics model
To have a good description of data from ψ(3686) → ΛΛη, an event generator based on the PWA results is developed to determine the detection efficiency. We vary the default configuration to a setup either with only the Λ(1690) or with a combination of Λ(1670) and Λ(1690), and consider the largest change in the detection efficiency as systematic uncertainty of the physical model. For ψ(3686) → ΛΛπ 0 , we use PHSP as the nominal event generator. The change in detection efficiency using an alternative model within the allowed phase space of the Λπ 0 system resonances is assigned as systematic uncertainty for the ΛΛπ 0 model.

• Intermediate decays
The uncertainties of the quoted decay branching fractions for the intermediate particles from PDG [14] are taken as systematic uncertainties.
• Number of ψ(3686) events The number of ψ(3686) events is determined from an analysis of inclusive hadronic ψ(3686) decays. The uncertainty of the number of ψ(3686) events, 0.6% [2], is taken as systematic uncertainty.
In addition, the systematic uncertainties associated with the PWA, which contribute to the measurement of the Λ(1670) mass and width and of the corresponding production branching fraction, are described below.
• Additional resonances To investigate the impact on the PWA results from other possible components, the analysis has been performed including additional possible states (e.g. Λ(1690)). The changes of the mass, width, and fitted fraction of Λ(1670) are considered as systematic uncertainties, and the largest one was chosen.
• Background uncertainty In the ψ(3686) → ΛΛη decays, the background level is quite low, and the events from the η sidebands are considered in the PWA. To estimate the uncertainty, the scale factor of background events from η sidebands has been varied by ±50%, and the largest variation of the results is assigned as systematic uncertainty.
• PHSP parametrization In the partial wave analysis, PHSP is parametrized as a resonance with an extremely large width, and fixed values of spin and parity. The contribution to the systematic uncertainty is estimated by replacing the spin parity of All the systematic uncertainty sources and values are summarized in Tables I and II, respectively, for the ψ(3686) → ΛΛπ 0 , ψ(3686) → ΛΛη decays, where the total uncertainties are given by the quadratic sum, assuming statistical independence of all the contributions. The distinction between additive and multiplicative sources of systematic uncertainties are indicated in the table.

VIII. RESULTS
The branching fractions of the decays of interest are calculated as where X is π 0 or η, N obs X is the number of signal candidates, N ψ(3686) is the number of ψ(3686) events determined with inclusive hadronic events, ǫ is the detection efficiency obtained from the MC simulation. B(Λ → pπ − ), B(π 0 → γγ) and B(η → γγ) are the corresponding branching fractions from PDG [14]. Using the numbers given in Table III,  Based on the PWA results, it was found that the ev-ident structure around the Λη mass threshold could be described by the Λ(1670). The mass and width are determined to be M = (1672 ± 5 ± 6) MeV/c 2 and Γ = (38 ± 10 ± 19) MeV, which are consistent with those in PDG [14]. The corresponding product branching fraction is calculated to be B(ψ(3686) → Λ(1670)Λ) × B(Λ(1670) → Λη) = (1.29 ± 0.31 ± 0.62) × 10 −5 , where the first uncertainty is statistical and the second is systematic. Due to the limited statistical significance of the ψ(3686) → ΛΛπ 0 signal (3.7σ), the upper limit of this branching fraction has been determined. We repeat the maximum-likelihood fits by varying the signal shape, nonpeaking background, peaking background as well as interference between ψ(3686) and continuum amplitudes, and take the most conservative upper limit among different choices. To incorporate the multiplicative systematic uncertainties in the calculation of the upper limit, the likelihood distribution is smeared by a Gaussian function with a mean of zero and a width equal to σ ǫ as shown below [29,30] L ′ (n ′ ) ∝ where L(n) is the likelihood distribution as a function of the yield n, ǫ 0 is the detection efficiency and σ ǫ is the multiplicative systematic uncertainty. Figure 7 shows the likelihood function without and with incorporating the systematic uncertainties. The upper limit on the number of ψ(3686) → ΛΛπ 0 events, N ′ UL , is determined to be 40, and the corresponding upper limit of the branching fraction is obtained to be ψ(3686) → ΛΛπ 0 < 2.47 × 10 −6 at the 90% confidence level (C.L.).
Compared with the branching fraction of J/ψ → ΛΛπ 0 and J/ψ → ΛΛη [1], the ratio between the branching fractions of ψ(3686) and J/ψ decaying to the same hadronic final state is defined as Q h . Taking PHSP factors into account, the ratio of ΛΛπ 0 and ΛΛη, Q h (π 0 ) and Q h (η), are calculated to be (1.4 ± 0.7)% and (2.3 ± 0.6)%, respectively, both contradicting the 12% rule significantly.