Observation of $\psi(3686)\to\Xi(1530)^{0}\bar{\Xi}(1530)^{0}$ and $\Xi(1530)^{0}\bar{\Xi}^0$

Using $(448.1\pm2.9)\times 10^{6}$ $\psi(3686)$ events collected with the BESIII detector and a single-baryon tagging technique, we present the first observation of the decays $\psi(3686)\to\Xi(1530)^{0}\bar{\Xi}(1530)^{0}$ and $\Xi(1530)^{0}\bar{\Xi}^0$. The branching fractions are measured to be ${\cal{B}}(\psi(3686)\to\Xi(1530)^{0}\bar{\Xi}(1530)^{0}) = (6.77\pm0.14\pm0.39)\times10^{-5}$ and ${\cal{B}}(\psi(3686)\to\Xi(1530)^{0}\bar{\Xi}^{0}) = (0.53\pm0.04\pm0.03)\times10^{-5}$. Here, the first and second uncertainties are statistical and systematic, respectively. In addition, the parameter associated with the angular distribution for the decay $\psi(3686)\to\Xi(1530)^{0}\bar{\Xi}(1530)^0$ is determined to be $\alpha = 0.32\pm0.19\pm0.07$, in agreement with theoretical predictions within one standard deviation.


I. INTRODUCTION
The production of the J/ψ and ψ(3686) resonances (here both denoted by ψ) in e + e − annihilation and their subsequent two-body hadronic decays can be used to test perturbative Quantum Chromodynamics (QCD) and QCD-based calculations [1], which concern topics such as gluon spin, quark-distribution amplitudes in baryonantibaryon (BB) pairs, and total hadron-helicity conservation. In particular, the decay of charmonium to BB pairs is a suitable avenue for testing these calculations due to its simple topology. However, only a few ψ → BB decay modes have been precisely measured until now. Many decay modes, such as ψ(3686) → Ξ(1530) 0Ξ (1530) 0 and Ξ(1530) 0Ξ0 decays, are still either unknown or less precisely measured [2].
A precise measurement of the angular distributions of ψ → BB is desirable to test helicity conservation rules and theoretical models [3,4], especially the color-octet contributions in these decays. The angular distribution for the process ψ → BB can be expressed as where θ B is the angle of the momentum vector of one of the baryons in the center-of-mass (CM) system with respect to the direction of the positron beam. Theoretically, the values of the angular-distribution parameter α have been predicted by first-order QCD calculations [3,4]. In the prediction of Claudson et al. [3], the mass of the final baryon is taken into consideration unitarily, while the constituent quarks inside the baryon are taken as massless when computing the decay amplitude. Yet Carimalo et al. [4] deemed that quark mass effects are more sensitive than electromagnetic contributions to the α value , and the mass effects at the quark level are taken into account. Experimentally, the values of α have been measured for several baryon-pair final states [5][6][7][8][9]. Positive values of α were found for the processes ψ(3686) → Ξ(1530) −Ξ (1530) + [10] and ψ(3686) → Ξ 0Ξ0 , while for other known processes, such as J/ψ → Σ 0Σ0 and Σ(1385)Σ(1385), the values were found to be negative. Chen and Ping [11] noted that the angular distribution parameter for ψ → BB could be negative when rescattering effects of BB in heavy quarkonium decays are taken into account. The BESIII experiment has collected a large data sample at the ψ(3686) peak, which can be used to investigate several baryon decay channels e.g., ψ(3686) → Ξ(1530) 0Ξ (1530) 0 , where the α value is predicted to be 0.32 by C. Carimalo [4].
An additional motivation for studying the baryonic decays of ψ mesons is that perturbative QCD predicts that the ratio of the branching fractions of the ψ(3686) and the J/ψ decaying into the same final state should follow the so-called '12% rule' [12]: Although the ratio has been measured with a wide variety of final states [13], the review in Ref. [14], which includes a comparison of theory and experiment, concludes that the current theoretical explanations are unsatisfactory. More experimental results, such as the decay of ψ(3686) → Ξ(1530) 0Ξ0 , are required to test the '12% rule' [12] and theoretical models [14].

II. BESIII DETECTOR AND MONTE CARLO SIMULATION
The BESIII detector [20] records symmetric e + e − collisions provided by the BEPCII storage ring [21], which operates with a peak luminosity of 1 × 10 33 cm −2 s −1 in the center-of-mass energy range from 2.0 to 4.9 GeV. The cylindrical core of the BESIII 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 with resistive plate counter muon identification modules interleaved with steel. The charged-particle momentum resolution at 1 GeV/c is 0.5%, and the specific energy loss (dE/dx) resolution 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 [22] Monte Carlo (MC) package, which includes the geometric description of the BESIII detector and the detector response, are used to determine detection efficiencies and to estimate background contributions. The simulation models the beam-energy spread and the initial state radiation (ISR) in the e + e − annihilations with the generator kkmc [23]. The inclusive MC sample includes the production of the ψ(3686) resonance, the ISR production of the J/ψ, and the continuum processes incorporated in kkmc [23]. The known decay modes are modelled with evtgen [24] using branching fractions taken from the Particle Data Group (PDG) [13], and the remaining unknown charmonium decays are modelled with lundcharm [25]. Final-state radiation (FSR) from charged final state particles is incorporated using the photos package [26]. To determine the detection efficiencies for the signal modes ψ(3686) → Ξ(1530) 0Ξ (1530) 0 and ψ(3686) → Ξ(1530) 0Ξ0 , 1 million simulated signal MC samples are generated for each reconstructed mode, taking into account the angular distributions measured in this analysis. The processes Ξ(1530) 0 → Ξ − π + with Ξ − → Λπ − and Λ → pπ − are generated uniformly in phase space. The data set collected at the CM energy of 3.65 (3.572) GeV with an integrated luminosity of 43.88 (23.14) pb −1 is used to estimate the contamination from the continuum processes e + e − → Ξ(1530) 0Ξ (1530) 0 and e + e − → Ξ(1530) 0Ξ0 [18].

III. DATA ANALYSIS
Charged tracks detected in the MDC are required to be within a polar angle (θ) range of | cos θ| < 0.93, where θ is defined with respect to the axis of the MDC. Particle identification (PID) for charged tracks combines measurements of the specific energy loss dE/dx in the MDC and the flight time measured by the TOF to form likelihoods L h (h = p, K, π) for each hadron h hypothesis. Tracks are identified as protons (pions) when the proton (pion) hypothesis has the greatest likelihood, i.e. L p(π) > L K and L p(π) > L π(p) .
The Λ candidates are reconstructed from pπ − pairs with an invariant mass lying within ±5 MeV/c 2 of the nominal Λ mass [13]. This mass window is chosen by optimizing the figure-of-merit (FOM) S √ S+B , where S is the number of signal events and B is the number of background events. In this case, S and B are taken from signal MC and the inclusive MC samples, respectively. A secondary-vertex fit [27] is performed with all the pπ − combinations; those with χ 2 < 500 are kept for further analysis. To further suppress background events, the decay length of the Λ candidate is required to be positive. In the case of multiple candidates, the one with an unconstrained mass closest to the nominal mass is retained. The Ξ − candidates are reconstructed from Λπ − pairs with an invariant mass within ±8 MeV/c 2 of the nominal Ξ − mass [13]. A secondary-vertex fit [27] is performed on all Λπ − combinations. Only the candidate with an invariant mass closest to the nominal mass is retained if there is more than one candidate in an event. The decay length of the Ξ − is required to be positive to further suppress background. The Ξ(1530) 0 candidates are reconstructed by combining Ξ − candidates with an additional π + , minimizing the variable |M Ξπ − M Ξ(1530) 0 |, where M Ξπ is the invariant mass of Ξ − π + and M Ξ(1530) 0 is the nominal mass of the Ξ(1530) 0 [13]. On average, there are 1.07 (1.02) candidates per event for the Ξ(1530) 0Ξ (1530) 0 (Ξ(1530) 0Ξ0 ) final state before the best candidate selection.
The presence of the anti-baryon candidatesΞ 0 andΞ(1530) 0 is inferred from the mass recoiling against the tagged Ξ − π + system, M recoil where E CM is CM energy, E Ξ − π + and P Ξ − π + are the energy and momentum of the selected Ξ − π + system in the e + e − rest frame. Analysis of the ψ(3686) inclusive MC sample with a generic event-type investigation tool, TopoAna [28], indicates that the main background events come from ψ(3686) → π + π − (π 0 π 0 )J/ψ with J/ψ → Ξ −Ξ+ . The background shape from the π 0 π 0 J/ψ channel is flat, while the background shape from the π + π − J/ψ channel is not smooth. We also investigate the contribution from continuum processes with the off-peak data at √ s = 3.65 (3.572) GeV. It is found that few events survive the selection criteria, and no Ξ(1530) 0 (Ξ 0 ) peaking contribution is seen due to the rather low sample size. Taking into account the normalization of the integrated luminosity and CM energy dependence of the cross section, the contribution from continuum processes is expected to be small and is neglected in the remainder of the analysis. In addition, there is the possibility of miscombinations of negative pions with similar momenta in both the baryon and anti-baryon decay chains in the signal events. Including these in the Ξ(1530) 0 or Ξ 0 reconstruction leads to a wrong-combination background (WCB).

A. Branching fractions
The signal yields for the two decays ψ(3686) → Ξ(1530) 0Ξ (1530) 0 and Ξ(1530) 0Ξ0 are determined by performing a simultaneous extended unbinned maximum likelihood fit to the M recoil Ξ − π + and M recoil Ξ + π − spectra, respectively. In the fit, the signal shapes are represented by the MC-simulated shape sampled from a multidimensional histogram using the signal MC simulation, convolved with a Gaussian function to take into account the detector-resolution difference between data and MC simulation, where the parameters of the Gaussian function are left free but constrained to be the same for the two charge conjugation modes. The background events from ψ(3686) → π + π − J/ψ with J/ψ → Ξ −Ξ+ are described by their respective MC-simulated shapes and the corresponding number of background events, computed according to their respective branching fractions, are fixed. The WCB is also described by the MC-simulated shape, and its contribution is fixed according to the MC simulation. The shape of the other sources of the background (Other-Bkg) looks smooth and is described by a third-order Chebychev polynomial function. The fit results are shown in Fig. 2. The statistical significance for both ψ(3686) → Ξ(1530) 0Ξ (1530) 0 and ψ(3686) → Ξ(1530) 0Ξ0 is found to be greater than 10σ.
The signal yields (N sig ), the detection efficiencies (ǫ) and the branching fractions (B(ψ(3686) → X)) for the analyzed decay modes are summarized in Table I. The systematic uncertainties on the measurements for the branching fractions will be discussed in Section V A.

B. Angular distribution
The angular-distribution parameter α for the decay ψ(3686) → Ξ(1530) 0Ξ (1530) 0 is determined by performing a least-squares fit based on Eq. (1) to the efficiencycorrected distribution of N sig versus cos θ B . The signal yields are extracted in eight equal bins of cos θ B in the range from −0.8 to 0.8 using the method described in Section IV A. The corresponding signal yield in each bin of cos θ B is corrected with an efficiency obtained from a signal MC sample generated uniformly in phase space.
The distribution of the efficiency-corrected signal yield from data together with the results of the fit are shown in Fig. 3, obtained by separate fits to the corresponding distributions. The goodness of the fit is calculated to be χ 2 /NDF = 3.58/6 for the Ξ(1530) 0 tag, and 4.82/6 for theΞ(1530) 0 tag, where NDF is the number of degrees of freedom. The values of α obtained from the fit are summarized in Table I. The systematic uncertainty on the measurement of α will be discussed in Sect. V B. Due to the limited sample size, α for the decay ψ(3686) → Ξ(1530) 0Ξ0 is not measured in this study.

A. Branching fractions
The systematic uncertainties in the branching fraction measurements are associated with knowledge of the tracking efficiency, PID efficiency, trigger efficiency, Ξ − reconstruction efficiency, the choice of mass windows and decay lengths for the Λ and Ξ − candidates, fit range, signal and background shapes, the knowledge of the branching fraction of the intermediate decay, and the uncertainty on the total number of ψ(3686) events. The uncertainties due to the tracking and PID efficiency (except for those associated with the π + from the Ξ(1530) 0 decay) and the Λ reconstruction efficiency have been incorporated in the Ξ − reconstruction efficiency, and the uncertainty due to the trigger efficiency is negligible [29].
Below is a summary of all the sources of the systematic uncertainties identified in this analysis.
• Tracking and PID: The uncertainties due to the π + tracking (coming from the Ξ(1530) 0 mother particle) efficiency and the PID efficiency are investigated with a control sample of ψ(3686) → π + π − pp decays. The observed difference in efficiency between the data and MC simulation, 1.0%, for both quantities, is assigned as the uncertainty for the tracking and PID efficiency separately [10].
• Ξ − reconstruction: The uncertainty due to the Ξ − reconstruction efficiency is estimated from a study of a control sample of ψ(3686) → Ξ −Ξ+ decays [10]. The efficiency differences between the data and MC simulation are found to be 3.0% for the Ξ − reconstruction and 6.0% for theΞ + reconstruction, which are taken as the systematic uncertainties.
• Decay length of Λ (Ξ − ): The uncertainty associated with the efficiency for the decay-length requirement on Λ (Ξ − ) candidates is estimated as the difference of selection efficiency on the data and MC simulation with a control sample J/ψ → Ξ −Ξ+ decays.
• Mass window of Λ (Ξ − ): The uncertainty associated with the mass-window requirement on the Λ (Ξ − ) candidates is estimated as the difference of the selection efficiency on the data and MC simulation, evaluated with the J/ψ → Ξ −Ξ+ control sample.
• Mass window of Ξ(1530) 0 : The uncertainty associated with the mass-window requirement on the Ξ(1530) 0 candidates is estimated by varying the requirement by ±1 MeV/c 2 . The largest difference in efficiency between data and MC simulation is taken as the systematic uncertainty.
• Fit range: Any potential bias associated with the choice of fit range is estimated by varying this range by ±10 MeV/c 2 . The maximum observed difference on the branching fraction is taken as the systematic uncertainty.
• Signal shape: The uncertainty due to the choice of signal shape is estimated by changing the MCsimulated shape with a kernel density estimation [30] of the unbinned signal MC derived shape.  The difference in the measured branching fraction is taken as the systematic uncertainty.
• Background shape: The uncertainty associated with the use of a Chebychev polynomial function to describe the remaining background events from miscellaneous sources is estimated by performing alternative fits with a second or fourthorder Chebychev polynomial function. The differ-ence in the measured branching fraction is taken as the systematic uncertainty.
• π + π − J/ψ background: The uncertainty associated with the π + π − J/ψ background has two components: the impact of the mass-window requirement on M recoil π + π − and the knowledged of the size of background from ψ(3686) → π + π − J/ψ with J/ψ → Ξ −Ξ+ . The uncertainty related to the M recoil π + π − mass window is estimated by changing the requirement by ±1 MeV/c 2 . The largest difference on branching fraction is taken as the systematic uncertainty. The uncertainty related to the amount of ψ(3686) → π + π − J/ψ with J/ψ → Ξ −Ξ+ background is estimated by changing the branching fraction of J/ψ → Ξ −Ξ+ within one standard deviation. The resulting difference on branching fraction yield is taken as the uncertainty. The total uncertainty associated with the π + π − J/ψ background is calculated by taking the quadratic sum of the two individual contributions.
• WCB: Any possible bias related to WCB is assessed by comparing the signal yields with and without the corresponding component included in the fit. The difference on the branching fraction is found to be negligible.
• Angular distribution: The MC simulation for the signal channels requires as input the angulardistribution parameters of the the baryon pairs. For the ψ(3686) → Ξ(1530) 0Ξ (1530) 0 mode, this uncertainty is estimated by varying the measured α values by ±1σ uncertainty in the MC simulation, and the resulting largest difference compared with the nominal efficiency is taken as the uncertainty. For the mode Ξ(1530) 0Ξ0 , we set α = 0 and generate a new MC sample to obtain the detection efficiency, whose difference relative to the nominal efficiency is taken as the uncertainty. The uncertainty of the MC model for simulating the decays of Ξ(1530), Ξ and Λ baryon is estimated by comparing the efficiency with and without considering the angular distribution of baryon pair via evtgen [24]. The total uncertainty related to the angular distribution is assigned to be the quadratic sum of all individual terms.
The total systematic uncertainty in the branchingfraction measurement is obtained by summing the individual contributions in quadrature, which is summarized in Table II.

B. Angular distribution
The sources of systematic uncertainty for the measurement of the α value include those from the signal yields in different cos θ B intervals and those associated with the α-value fit. The uncertainties for the signal yield arise from the choice of fit range, the signal and background shapes, the MC model, the level of the WCB and that of the background size from ψ(3686) → π + π − J/ψ with J/ψ → Ξ −Ξ+ . These uncertainties are estimated with the same method as described in Sec. V A, where the shifts in the value of α with respect to that from the default fit are taken as the systematic uncertainties. The uncertainties associated with the α-value fit are assigned as follows.
• The cos θ B binning: This uncertainty is estimated by comparing the different binning schemes.
In the nominal binning scheme, the cos θ B ∈ [−0.8, 0.8] is divided into 8 bins, while an alternative binning scheme divides the |cosθ B | ∈ [0, 0.8] into 8 bins. The resulting differences on the α value are taken as the systematic uncertainties.
• The cos θ B range: This uncertainty is estimated by repeating the fit after changing the cos θ B range to [−0.7, 0.7], with the same bin size as the nominal fit. The differences on the α value are found to be negligible for both tagging modes.
All the systematic uncertainties for the α measurement are summarized in Table III, where the total systematic uncertainty is the quadratic sum of the contributions.

VI. RESULTS
Combined branching fractions and α values are calculated with a weighted least-squares method [31]. The single-baryon recoil mass method leads to some double counting of the Ξ(1530) 0Ξ (1530) 0 final state; MC studies indicate this occurs at a rate of about 10%, which is taken as a common statistical uncertainty when calculating the combined branching fractions and α value. The systematic uncertainties are weighted to properly account for common and uncommon systematic uncertainties using the method in the Ref. [10]. The results are summarized in Table I. The measured value of α is in agreement with the predictions of Refs. [4].
To test isospin symmetry, the ratios of the combined branching fractions for the decay ψ(3686) → X between neutral and charged-modes, R ψ(3686)→X = B com neutral B com charged , are calculated by assuming that the common systematic uncertainties between the neutral and charged decay modes cancel in the measurement. The ratio is determined to be 0.59 ± 0.03 ± 0.06 for the decay ψ(3686) → Ξ(1530)Ξ(1530), in significant disagreement with the expectation from isospin symmetry. The corresponding result for ψ(3686) → Ξ(1530)Ξ is found to be 0.76±0.13±0.06, which is compatible with that expected from isospin symmetry. In both cases, the first uncertainty is statistical and the second is systematic.

ACKNOWLEDGMENTS
The BESIII collaboration thanks the staff of BEPCII and the IHEP computing center for their strong support. This work is supported in part by National Key Basic Research Program of China under Contracts