The Belle Collaboration

# Measurement of branching fractions of
$\Lambda_{c}^{+}\\!\rightarrow\\!{}pK_{S}^{0}K_{S}^{0}$ and
$\Lambda_{c}^{+}\\!\rightarrow\\!{}pK_{S}^{0}\eta$ at Belle

###### Abstract

We present a study of a singly Cabibbo-suppressed decay
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and a Cabibbo-
favored decay ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ based on
$980~{}\mbox{\,fb}^{-1}$ of data collected by the Belle detector, operating at
the KEKB energy-asymmetric $e^{+}e^{-}$ collider. We measure their branching
fractions relative to ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$:
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}})/\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})={(1.47\pm
0.08\pm 0.04)\times 10^{-2}}$ and
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta})/\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})={(2.72\pm
0.06\pm 0.13)\times 10^{-1}}$. Combining with the world average
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})$, we have the
absolute branching fractions,
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}})={(2.34\pm
0.12\pm 0.06\pm 0.12)\times 10^{-4}}$ and
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta})={(4.32\pm
0.10\pm 0.20\pm 0.22)\times 10^{-3}}$. The first and second uncertainties are
statistical and systematic, respectively, while the third ones arise from the
uncertainty on
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})$. The
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ is observed for
the first time and has a statistical significance more than $10\sigma$, and
the branching fraction of
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ is consistent with the
world average but has a threefold improvement of precision.

##  I Introduction

The weak decays of charmed baryons provide an excellent platform for
understanding Quantum Chromodynamics with transitions involving the charm
quark. Their decay amplitudes consist of factorizable and non-factorizable
contributions. The latter may play a non-trivial or essential role in the
decay amplitude and are described by various approaches, including the pole
model [1, 2], the covariant confined quark model [3, 4], current algebra [5,
6, 7] and SU(3)F symmetry [8, 9, 10]. To date, there is no established
phenomenological model to describe the complicated physics of baryon decays.
Precise measurements of branching fractions of charmed baryon weak decays are
useful to study the dynamics of charmed baryons and test the predictions of
theoretical models. In addition, the singly Cabibbo-suppressed (SCS) charm
decays are essential probes to search for $C\\!P$ violation in the charm
sector [11, 12, 13] and hunt for new physics beyond the standard model [14,
15, 16].

Experimentally, the investigation of weak decays of charmed baryons is more
challenging than that of charmed mesons, mainly due to lower production rates.
For the lightest state, $\Lambda_{c}^{+}$, hadronic modes have been studied at
several experiments, but some of them have yet to be observed or the precision
can be much improved [17]. The Cabibbo-favored (CF)
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ and the SCS
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ decays are two
examples 111Throughout this paper, the charge-conjugate modes are implied..
For ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, the world average
branching fraction,
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta})\\!=\\!{(4.15\pm
0.90)\times 10^{-3}}$ [17], still has a large uncertainty ($22\%$). The
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$, of which the
predicted branching fraction is
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}})\\!=\\!{(1.9\pm
0.4)\times 10^{-3}}$ based on SU(3)F symmetry [19], has not previously been
observed.

In this paper, we present a precise measurement of
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}})$ and
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta})$ based on
the full Belle dataset. These three-body decay modes also allow the study of
intermediate states, such as the $N^{*}(1535)$, the nature of which is still
unknown. Understanding the nature of $N^{*}(1535)$ is very challenging and
important in hadron physics. The $N^{*}(1535)$, with spin parity
$J^{P}\\!=\\!1/2^{-}$, has a larger mass than the radial excitation
$N^{*}(1440)$, which is opposite to the prediction in classical constituent
quark models [20]. It also couples strongly to channels with strangeness, such
as $\eta{}N$ and $K\Lambda$, which is difficult to explain within the naive
constituent quark models [21, 22]. The inclusion of five-quark components
gives a natural explanation for these properties [23]. The
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ decay, in which the
final-state $p\eta$ is in a pure isospin $I\\!=\\!1/2$ state, is an ideal
process to study the $N^{*}(1535)$ resonance, since $N^{*}(1535)$ has a large
branching ratio to $p\eta$, in $S$-wave. The other interesting intermediate
resonances are $a_{0}(980)$ and $f_{0}(980)$ in
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ since they both
couple to $K\overline{K}$. The nature of the light scalars $f_{0}(980)$ and
$a_{0}(980)$ remains not fully understood and continues to cause controversy
[24, 25, 26], but they are often interpreted as compact tetraquark states [27,
28, 29] or $K\overline{K}$ bound states [30, 31]. Therefore, we reconstruct
the Dalitz plots of
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ decays to check such
interesting intermediate resonances.

##  II Detector and data set

This analysis uses the full dataset recorded by the Belle detector [32]
operating at the KEKB energy-asymmetric $e^{+}e^{-}$ collider [33]. This data
sample corresponds to a total integrated luminosity of
$980~{}\mbox{\,fb}^{-1}$ collected at or near the $\Upsilon(nS)$
($n\\!=\\!1,\,2,\,3,\,4,\,5$) resonances. The Belle detector is a large-solid-
angle magnetic spectrometer consisting of a silicon vertex detector, a central
drift chamber (CDC), an array of aerogel threshold Cherenkov counters (ACC), a
barrel-like arrangement of time-of-flight scintillation counters (TOF), and an
electromagnetic calorimeter (ECL) consisting of CsI(Tl) crystals. These
components are all located inside a superconducting solenoid coil that
provides a 1.5 T magnetic field. The iron flux-return of the magnet is
instrumented to detect $K^{0}_{L}$ mesons and to identify muons. The detector
is described in detail elsewhere [32].

Monte Carlo (MC) simulated events are generated with evtgen [34] and pythia
[35], and are subsequently processed through the full detector simulation
based on geant3 [36]. Final-state radiation from charged particles is included
at the event generation stage using photos [37]. “Generic” MC samples include
$B\overline{B}$ events and continuum processes
$e^{+}e^{-}\\!\rightarrow\\!{}q\bar{q}$ ($q\\!=\\!u,\,d,\,s,\,c)$
corresponding to an integrated luminosity three times that of the data.
Samples of MC events of $\Lambda_{c}^{+}$ signal decay modes are produced in
the $e^{+}e^{-}\\!\rightarrow\\!{}c\overline{c}$ process, decayed uniformly in
three-body phase space, and used to study the efficiency.

##  III Event selection

We reconstruct the two signal modes
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ and their reference
mode ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$. The event selections
are optimized based on a figure of merit (FOM), defined as ${\rm
FOM}\\!=\\!N_{S}/\sqrt{N_{B}}$ for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ due to its
branching fraction having not yet been measured, and ${\rm
FOM}\\!=\\!N_{S}/\sqrt{N_{S}+N_{B}}$ for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ assuming its current
world average branching fraction [17]. Here $N_{S}$ and $N_{B}$ are the
expected yields of signal and background, respectively, based on numbers of
candidates in the $M(\Lambda_{c}^{+})$ signal regions, where
$M(\Lambda_{c}^{+})$ is the invariant mass of reconstructed $\Lambda_{c}^{+}$
candidates. These signal regions are defined to be within 10, 22, and 18
MeV/$c^{2}$ of the nominal $\Lambda_{c}^{+}$ mass [17] for the
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$,
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$ channels, respectively;
each signal band includes $\approx$98% of the signal. For the expected
background, $N_{B}$, the number found in MC is multiplied by the data/MC yield
ratio in the $M(\Lambda_{c}^{+})$ sideband region
($30<|M(\Lambda_{c}^{+})-m_{\Lambda_{c}^{+}}|\\!<\\!{50~{}{\rm MeV}/c^{2}}$),
where $m_{\Lambda_{c}^{+}}$ is the nominal $\Lambda_{c}^{+}$ mass [17].

The particle identification (PID) likelihood for a given particle hypothesis,
$\mathcal{L}_{i}$ ($i\\!=\\!\pi,\,K,\,p$), is calculated from the photon yield
in the ACC, energy-loss measurements in the CDC, and time-of-flight
information from the TOF [38]. Charged tracks satisfying
$\mathcal{R}(p|K)\\!=\\!\mathcal{L}_{p}/(\mathcal{L}_{p}+\mathcal{L}_{K})\\!>\\!0.6$
and
$\mathcal{R}(p|\pi)\\!=\\!\mathcal{L}_{p}/(\mathcal{L}_{p}+\mathcal{L}_{\pi})\\!>\\!0.6$,
are identified as protons. These PID requirements have signal efficiencies of
94% for ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and 97%
for ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$.

For proton candidates, the point on the track nearest to the axis defined by
the positron beam and in the direction opposite to it (“$z$-axis”) is required
to be within 3.0 cm of the interaction point in the $z$-direction and within
1.0 cm on the transverse ($x$-$y$) plane. This requirement rejects tracks not
originating at the interaction point (IP) and introduces a negligible signal
efficiency loss ($<0.01\%$).

Candidate $K{}^{0}_{S}$’s are reconstructed from pairs of oppositely-charged
tracks, treated as pions, using an artificial neural network (NN) [39]. The NN
utilizes the following 13 input variables: the $K{}^{0}_{S}$ momentum in the
laboratory frame; the separation in $z$ between the two $\pi^{\pm}$ tracks at
their intersection in the $x$-$y$ plane; for each track, the nearest distance
to the IP in the $x$-$y$ plane; the $K{}^{0}_{S}$ flight length in the $x$-$y$
plane; the angle between the $K{}^{0}_{S}$ momentum and the vector joining the
IP to the $K{}^{0}_{S}$ decay vertex; in the $K{}^{0}_{S}$ rest frame, the
angle between the $\pi^{+}$ momentum and the laboratory-frame boost direction;
and, for each $\pi^{\pm}$ track, the number of CDC hits in both stereo and
axial views, and the presence or absence of SVD hits. Detailed information is
provided elsewhere [40]. The invariant mass of the reconstructed
$K{}^{0}_{S}\\!\rightarrow\\!\pi^{+}\pi^{-}$ candidate is required to lie
within $10~{}{\rm MeV}/c^{2}$ of the nominal $K{}^{0}_{S}$ mass [17]; this
includes 99.9% of the $K{}^{0}_{S}$ signal. The two pion tracks from each
$K{}^{0}_{S}$ candidate are refitted to originate from a common vertex and
constrained to have invariant mass equal to the nominal $K{}^{0}_{S}$ mass
[17]. The corresponding fit quality $\chi_{\rm mv}^{2}(K{}^{0}_{S})$ is
required to be smaller than 100. The selected $K{}^{0}_{S}$ sample has a
purity of greater than 98%.

Photon candidates are identified as energy clusters in the ECL that are not
associated with any charged track. The ratio of the energy deposited in the
3$\times$3 array of crystals centered on the crystal with the highest energy,
to the energy deposited in the corresponding 5$\times$5 array of crystals, is
required to be greater than 0.8. The photon energy is required to be greater
than 50 MeV in the barrel region (covering the polar angle
$32^{\circ}\\!<\\!\theta\\!<\\!129^{\circ}$), and greater than 100 MeV in the
endcap region ($12^{\circ}\\!<\\!\theta\\!<\\!31^{\circ}$ or
$132^{\circ}\\!<\\!\theta\\!<\\!157^{\circ}$).

Candidate $\eta\\!\rightarrow\\!\gamma\gamma$ decays are reconstructed from
photon pairs having an invariant mass satisfying $500~{}{\rm
MeV}/c^{2}\\!<\\!M(\gamma\gamma)\\!<\\!580~{}{\rm MeV}/c^{2}$ ($3\sigma$ in
$M_{\eta}(\gamma\gamma)$ resolution). The invariant mass of each $\eta$
candidate is constrained to the nominal $\eta$ mass [17] at the
$\Lambda_{c}^{+}$ decay vertex (described below). The fit quality of this mass
constraint is required to satisfy $\chi_{m}^{2}(\eta)\\!<\\!8$, and the
resulting $\eta$ momentum in the laboratory frame is required to be greater
than 0.4 GeV/$c$. To further suppress the background, $\eta$ candidates are
vetoed if either of daughters can be paired with another photon such that the
$\gamma\gamma$ pair has an invariant mass within $2.5\sigma$ of the nominal
$\pi^{0}$ mass ($\sigma\\!=\\!5$ MeV/$c^{2}$). This $\pi^{0}$-veto results in
a signal loss of 28% and removes 72% of background.

The $\Lambda_{c}^{+}$ candidates are assembled by forming combinations of the
final-state particles for each mode. The $p$ and $K{}^{0}_{S}$ are required to
originate from a common vertex (denoted the $\Lambda_{c}^{+}$ decay vertex and
the $K{}^{0}_{S}$ production vertex) with a fit quality $\chi_{\rm
vtx}^{2}\\!<\\!24$. To reduce combinatorial background, the scaled momentum of
the $\Lambda_{c}^{+}$ candidate, defined as
$x_{p}\\!=\\!p^{*}c/\sqrt{s/4-M^{2}(\Lambda_{c}^{+})\cdot c^{4}}$, is required
to be greater than 0.48, where $s$ is the square of the center-of-mass energy
and $p^{*}$ is the momentum of reconstructed $\Lambda_{c}^{+}$ candidates in
the $e^{+}e^{-}$ center-of-mass frame.

For the SCS decay
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$, a
non-$K{}^{0}_{S}$ peaking background from the CF decay
$\Lambda_{c}^{+}\\!\rightarrow\\!pK{}^{0}_{S}\pi^{+}\pi^{-}$ exists, even
though it is suppressed by the vertex fit and $K{}^{0}_{S}$ selection. The
$K{}^{0}_{S}$ decay length $L$ is determined by the projection of the vector
joining the $K{}^{0}_{S}$ production and decay vertices onto the $K{}^{0}_{S}$
momentum direction, and its corresponding uncertainty $\sigma_{L}$ is
calculated by propagating uncertainties in the vertices and the $K{}^{0}_{S}$
momentum, including their correlations. To suppress the non-$K{}^{0}_{S}$
peaking CF background, we require the significance of the $K{}^{0}_{S}$ decay
length $L/\sigma_{L}(K{}^{0}_{S})\\!>\\!10$ for the slower of the two
$K{}^{0}_{S}$’s in
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$. This
requirement reduces the signal efficiency by 3%, and rejects 80% of
non-$K{}^{0}_{S}$ peaking background. The remaining non-$K{}^{0}_{S}$ peaking
background is ignored in the $M(\Lambda_{c}^{+})$ fits because it has a tiny
ratio 0.4% to signal based on the MC studies with the branching fraction
$(1.6\pm 0.12)\%$ [17], but considered in the systematic uncertainty.

After applying all selection criteria to the data, we find 1.03, 1.06, and
1.01 candidates per event for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$,
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$, respectively, in
candidates selected from the entire $M(\Lambda_{c}^{+})$ fit region
($|M(\Lambda_{c}^{+})-m_{\Lambda_{c}^{+}}|\\!<\\!0.05~{}{\rm GeV}/c^{2}$).
Correspondingly, about 3.1%, 5.7% and 1.2% of events have multiple signal
candidates, which do not introduce any peaking background. We retain all
candidates for this branching fraction measurement.

##  IV Yield extraction

The signal yield is extracted by an unbinned extended maximum likelihood fit
to the $M(\Lambda_{c}^{+})$ distribution. The signal probability density
function (PDF) is a sum of three symmetric Gaussian functions for the
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}K{}^{0}_{S}$ mode, a sum of
one symmetric Gaussian and two asymmetric Gaussians for the
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ mode, and a sum of one
symmetric Gaussian and three asymmetric Gaussians for the
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$ mode. The Gaussian
functions share a common mean parameter but have different width parameters.
The fit is first performed on truth-matched signal MC events.

In fitting data, the mean is allowed a common shift ($\delta_{\mu}$) from the
value found in MC, and the widths are those found in MC, multiplied by a
common scaling factor ($k_{\sigma}$). The background PDF is a first-order
polynomial function for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and a second-
order polynomial function for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$. The background parameters
are floated to account for differences between the experimental data and MC
simulated samples. The results are shown in Fig. 1, along with the pulls
$(N_{\rm data}-N_{\rm fit})/\sigma_{\rm data}$ where $\sigma_{\rm data}$ is
the error on $N_{\rm data}$. The pull distributions demonstrate that the data
are statistically consistent with the fitted shapes. The signal and background
yields are listed in Table 1.

For the ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ mode, we
obtain the difference in the log likelihoods obtained from fits performed with
and without a signal PDF, ${\Delta\ln\mathcal{L}\\!=\\!524}$; as the number of
degrees of freedom without a signal component is three less than that with a
signal component (parameters $N_{\rm sig}$, $\delta_{\mu}$ and $k_{\sigma}$
are dropped), and this value of ${\Delta\ln\mathcal{L}}$ corresponds to a
statistical significance greater than $10\sigma$. This measurement constitutes
the first observation of this SCS $\Lambda_{c}^{+}$ decay.

Table 1: The fitted yields of signal and background in the overall fit region (FR) and the signal region (SR) for the ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$, ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, and ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$ modes. For the definition of these regions, see the text. Yields | ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ | ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ | ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$
---|---|---|---
$N_{\rm sig}^{\rm FR}$ | $\,\,2442\pm 103$ | $12877\pm 317$ | $515296\pm 1129$
$N_{\rm bkg}^{\rm FR}$ | $41138\pm 222$ | $75144\pm 403$ | $627427\pm 1177$
$N_{\rm sig}^{\rm SR}$ | $2391\pm 101$ | $12641\pm 311$ | $500457\pm 1096$
$N_{\rm bkg}^{\rm SR}$ | $8228\pm 44$ | $32935\pm 177$ | $226055\pm 424\,\,$

\begin{overpic}[width=144.39223pt]{MassLc_LcTopKsKs_data.eps}\put(78.0,65.0){\large(a)}\end{overpic}\begin{overpic}[width=144.39223pt]{MassLc_LcTopKsEta_data.eps}\put(78.0,65.0){\large(b)}\end{overpic}\begin{overpic}[width=144.39223pt]{MassLc_LcTopKs_data.eps}\put(78.0,65.0){\large(c)}\end{overpic}

Figure 1: The distributions of invariant mass of $\Lambda_{c}^{+}$ candidates
(points with error bars) and corresponding fit results (red curves) for (a)
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$, (b)
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, and (c)
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$, respectively. The red
(blue) dashed histograms show the signal (background).

##  V Branching fraction

For the three-body decay modes, the Dalitz plots for candidates in the
$M(\Lambda_{c}^{+})$ signal region and sideband region are shown in Figs. 2(a,
b) for ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and Figs.
2(d, e) for ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$. For
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$, Bose symmetry
requires invariance under the exchange of the two $K{}^{0}_{S}$’s, hence the
Dalitz plot for two $pK{}^{0}_{S}$ masses is symmetric. We plot
$M^{2}(pK{}^{0}_{S})_{\rm max}$ versus $M^{2}(pK{}^{0}_{S})_{\rm min}$ in half
of the Dalitz plot, as shown in Figs. 2(a–c), and use it to measure the
branching fraction.

For each mode, a large MC sample of signal events, generated uniformly across
the decay phase space, is used to determine the reconstruction efficiency. For
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, the efficiencies are
calculated in bins across the phase space, based on truth-matched signal yield
in the $M(\Lambda_{c}^{+})$ signal region. The results are shown in Fig. 2(c)
for ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and Fig.
2(f) for ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$. For the
reference mode ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$, the
efficiency is $\varepsilon_{0}\\!=\\!(33.85\pm 0.05)\%$ in the
$M(\Lambda_{c}^{+})$ fit region.

\begin{overpic}[width=144.39223pt]{DPinMSR_LcTopKsKs_unblind_halfDP.eps}\put(72.0,62.0){\large(a)}\end{overpic}\begin{overpic}[width=144.39223pt]{DPinMSB_LcTopKsKs_unblind_halfDP.eps}\put(72.0,62.0){\large(b)}\end{overpic}\begin{overpic}[width=144.39223pt]{DPeff_pKsVspKs_LcTopKsKs.eps}\put(72.0,62.0){\large(c)}\end{overpic}

\begin{overpic}[width=144.39223pt]{DPinMSR_LcTopKsEta_unblind_final.eps}\put(72.0,62.0){\large(d)}\end{overpic}\begin{overpic}[width=144.39223pt]{DPinMSB_LcTopKsEta_unblind_final.eps}\put(72.0,62.0){\large(e)}\end{overpic}\begin{overpic}[width=144.39223pt]{DPeff_LcTopKsEta.eps}\put(72.0,62.0){\large(f)}\end{overpic}

Figure 2: Plots (a, d) show the Dalitz plots in the $M(\Lambda_{c}^{+})$
signal region, and (b, e) show the Dalitz plots in the $M(\Lambda_{c}^{+})$
sideband region for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ (top) and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ (bottom). Plots (c, f)
show the average signal efficiency in bins across the Dalitz plane. The red
curves show the edges of kinematic phase-space region of the decays.

For ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, we correct for
efficiency bin-by-bin across the phase space. The Dalitz plots are divided
uniformly into 7$\times$7 bins for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and 5$\times$5
bins for ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, as shown in
Figs. 2(c, f) respectively. The efficiency-corrected yields are

$\displaystyle N_{\rm{corr}}$ $\displaystyle=$
$\displaystyle\sum_{i}(N_{i}^{\rm tot}-N_{\rm bkg}^{\rm
SR}f_{i}^{\rm{bkg}})/\varepsilon_{i}\,,$ (1)

where $N_{i}^{\rm{tot}}$ is the raw yield in the $i^{\rm th}$ bin of the
Dalitz plot in $M(\Lambda_{c}^{+})$ signal region, $N_{\rm bkg}^{\rm SR}$ is
the fitted background yield as listed in Table 1, $f_{i}^{\rm{bkg}}$ is the
fraction of background in the $i^{\rm th}$-bin, with
$\sum_{i}{}f_{i}\\!=\\!1$. These fractions are obtained from the Dalitz plot
distribution of events in the $M(\Lambda_{c}^{+})$ sideband region, shown in
Fig. 2(b) for ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$
and Fig. 2(e) for ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$.
Using the generic MC sample, we find that the Dalitz plot in the chosen
$M(\Lambda_{c}^{+})$ sideband region is consistent with the generic background
in the $M(\Lambda_{c}^{+})$ signal region. The uncertainties on each variable
in Eq. (1) have been considered and are propagated into the efficiency-
corrected yields, $N_{\rm corr}$. We obtain

$\displaystyle
N_{\rm{corr}}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}})$
$\displaystyle=$ $\displaystyle(1.55\pm 0.08)\times 10^{4}\,,$ (2)
$\displaystyle
N_{\rm{corr}}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta})$
$\displaystyle=$ $\displaystyle(1.63\pm 0.04)\times 10^{5}\,.$ (3)

The relative branching fractions of signal modes to reference mode are
determined by Eqs. (4, 5).

$\displaystyle\hskip-20.0pt\frac{\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}})}{\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})}$
$\displaystyle=$
$\displaystyle\frac{N_{\rm{corr}}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}})}{\mathcal{B}(K{}^{0}_{S}\\!\rightarrow\\!\pi^{+}\pi^{-})N^{\rm
FR}_{\rm{sig}}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})/\varepsilon_{0}}\,,$
(4)
$\displaystyle\hskip-20.0pt\frac{\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta})}{\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})}$
$\displaystyle=$
$\displaystyle\frac{N_{\rm{corr}}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta})}{\mathcal{B}(\eta\\!\rightarrow\\!\gamma\gamma)N^{\rm
FR}_{\rm{sig}}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})/\varepsilon_{0}}\,.$
(5)

Inserting the efficiency-corrected yields in Eqs. (2, 3), $N^{\rm
FR}_{\rm{sig}}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})$ in Table 1
and $\varepsilon_{0}$ for ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$,
and the world averages
$\mathcal{B}(K{}^{0}_{S}\\!\rightarrow\\!\pi^{+}\pi^{-})\\!=\\!(69.20\pm
0.05)\%$ and $\mathcal{B}(\eta\\!\rightarrow\\!\gamma\gamma)\\!=\\!(39.41\pm
0.20)\%$ [17], we find

$\displaystyle\frac{\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}})}{\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})}$
$\displaystyle=$ $\displaystyle(1.47\pm 0.08)\times 10^{-2}\,,$ (6)
$\displaystyle\frac{\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta})}{\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})}$
$\displaystyle=$ $\displaystyle(2.72\pm 0.06)\times 10^{-1}\,,$ (7)

Combining with the world average branching fraction of reference mode
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})\\!=\\!(1.59\pm
0.08)\%$ [17], we have the absolute branching fractions:

$\displaystyle\hskip-10.0pt\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}})$
$\displaystyle=$ $\displaystyle(2.34\pm 0.12\pm 0.12)\times 10^{-4}\,,$ (8)
$\displaystyle\hskip-10.0pt\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta})$
$\displaystyle=$ $\displaystyle(4.32\pm 0.10\pm 0.22)\times 10^{-3}\,,$ (9)

where the uncertainties are statistical and from the uncertainty on
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})$.

We examine the Dalitz plots for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, after background
subtraction and efficiency correction, for intermediate resonances. In
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$, clear evidence
for $f_{0}(980)$ or $a_{0}(980)^{0}$ (labeled as $S_{0}(980)$) near the
$K{}^{0}_{S}K{}^{0}_{S}$ threshold is seen, as shown in Fig. 3. In
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, a significant
enhancement consistent with $N^{*}(1535)$ is found near the $p\eta$ threshold,
as shown in Fig. 4. In the future, amplitude analyses of these decays can be
expected to improve our understanding of the nature of $S_{0}(980)$ and
$N^{*}(1535)$.

\begin{overpic}[width=108.405pt]{DPinMSR_LcTopKsKs_unblind_sig.eps}\put(71.0,61.0){\large(a)}\end{overpic}\begin{overpic}[width=108.405pt]{DPproj12_LcTopKsKs_unblind_sig.eps}\put(75.0,58.0){\large(b)}\end{overpic}

\begin{overpic}[width=108.405pt]{DPproj23_LcTopKsKs_unblind_sig.eps}\put(75.0,58.0){\large(c)}\end{overpic}

Figure 3: For ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$,
the Dalitz plot after background subtraction and efficiency correction bin-by-
bin and its projections superimposing with signal MC produced by phase space
mode (blue histograms). This symmetric Dalitz plot and its projections show
two entries per candidate, one for each possible $pK{}^{0}_{S}$ combination. A
dominant structure near the $K{}^{0}_{S}K{}^{0}_{S}$ threshold, which we
identify with $f_{0}(980)$ or $a_{0}(980)^{0}$, is clearly seen.

\begin{overpic}[width=108.405pt]{DPinMSR_LcTopKsEta_unblind_sig.eps}\put(71.0,61.0){\large(a)}\end{overpic}\begin{overpic}[width=108.405pt]{DPproj23_LcTopKsEta_unblind_sig.eps}\put(75.0,58.0){\large(b)}\end{overpic}

\begin{overpic}[width=108.405pt]{DPproj12_LcTopKsEta_unblind_sig.eps}\put(75.0,58.0){\large(c)}\end{overpic}\begin{overpic}[width=108.405pt]{DPproj13_LcTopKsEta_unblind_sig.eps}\put(75.0,58.0){\large(d)}\end{overpic}

Figure 4: For ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, the
Dalitz plot after background subtraction and efficiency correction bin-by-bin
and its projections superimposing with signal MC produced by phase space mode
(blue histograms). A significant structure of $N^{*}(1535)$ near the $p\eta$
threshold is found.

##  VI Systematic uncertainty

In measuring the ratio of branching fractions, many systematic uncertainties
cancel, as they affect both the signal and reference modes. The remaining
systematic uncertainties are summarized in Table 2 and introduced in detail
below.

Table 2: Relative systematic uncertainties of the branching fractions of ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, and the uncertainty from the branching fraction of the reference mode. sources | $\mathcal{B}_{{\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}}$ | $\mathcal{B}_{{\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}}$
---|---|---
$K{}^{0}_{S}$ reconstruction | 1.4% | 0.4%
proton PID efficiency | 0.9% | 0.5%
$\eta$ reconstruction | – | 4.0%
$M(\Lambda_{c}^{+})$ fit procedure | 1.9% | 2.3%
efficiency-correction procedure | 0.6% | 0.4%
non-$K{}^{0}_{S}$ peaking background | 0.8% | –
$\delta\mathcal{B}/\mathcal{B}(K{}^{0}_{S}\\!\rightarrow\\!\pi^{+}\pi^{-},\,\eta\\!\rightarrow\\!\gamma\gamma)$ | 0.1% | 0.5%
total syst. uncertainty | 2.7% | 4.7%
$\delta\mathcal{B}/\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})$ | 5.0% | 5.0%

The systematic uncertainty associated with the $K{}^{0}_{S}$ reconstruction is
considered as follows. A table of $K{}^{0}_{S}$ efficiency ratios of data to
MC in eight bins of the $K{}^{0}_{S}$ momentum distribution,
$R{}_{\varepsilon}^{K{}^{0}_{S}}$, is determined based on a control sample
${D^{*\pm}\\!\rightarrow\\!(D^{0}\\!\rightarrow\\!K{}^{0}_{S}\pi^{0})\pi^{\pm}}$.
The unfolded momentum distribution in data of $K{}^{0}_{S}$ from signal is
obtained using the ${}_{s}\mathcal{P}lot$ technique [41]. From one
$R{}_{\varepsilon}^{K{}^{0}_{S}}$ table, we can determine the average ratios:
(1) for $\Lambda_{c}^{+}\\!\rightarrow\\!{}pK^{0}_{S,{\rm fast}}K^{0}_{S,{\rm
slow}}$ where the subscript ‘fast’ (‘slow’) indicates the faster (slower) of
two $K{}^{0}_{S}$’s in the final state,
$\overline{R}{}_{\varepsilon}^{K{}^{0}_{S}}\\!=\\!\sum_{i}^{8}\sum_{j}^{8}N_{ij}(R{}_{\varepsilon,i}^{K^{0}_{S,{\rm
fast}}}R{}_{\varepsilon,j}^{K^{0}_{S,{\rm
slow}}})/\sum_{i}^{8}\sum_{j}^{8}N_{ij}$ calculated on the two-dimensional
$(p_{\scriptscriptstyle{K^{0}_{S,{\rm
fast}}}},\,p_{\scriptscriptstyle{K^{0}_{S,{\rm slow}}}})$ distribution due to
the correlations between the momenta of two $K{}^{0}_{S}$’s. Here $N_{ij}$ and
$(R{}_{\varepsilon,i}^{K^{0}_{S,{\rm fast}}}R{}_{\varepsilon,j}^{K^{0}_{S,{\rm
slow}}})$ are the yield and the averaged $R{}_{\varepsilon}^{K{}^{0}_{S}}$,
respectively, in the bin of $i^{\rm th}$ raw and $j^{\rm th}$ column of such
two-dimensional momenta distribution; (2) for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}(\eta)$,
${\overline{R}{}_{\varepsilon}^{K{}^{0}_{S}}\\!=\\!\sum_{i}^{8}N_{i}R{}_{\varepsilon,i}^{K{}^{0}_{S}}/\sum_{i}^{8}N_{i}}$
calculated on the one-dimensional $p_{\scriptscriptstyle K{}^{0}_{S}}$
distribution. Here $N_{i}$ and $R{}_{\varepsilon,i}^{K{}^{0}_{S}}$ are the
yield and the averaged $R{}_{\varepsilon}^{K{}^{0}_{S}}$, respectively, in the
$i^{\rm th}$ bin of such one-dimensional distribution. We build 10000
$R{}_{\varepsilon}^{K{}^{0}_{S}}$ tables by randomly fluctuating
$R{}_{\varepsilon,i}^{K{}^{0}_{S}}$ in each bin according to its uncertainty
and calculate $\overline{R}{}_{\varepsilon}^{K{}^{0}_{S}}$ for each. We take
the mean and root-mean-square (RMS) values from the distribution of
$\overline{R}{}_{\varepsilon,\,{\rm
sig.}}^{K{}^{0}_{S}}/\overline{R}{}_{\varepsilon,\,{\rm
ref.}}^{K{}^{0}_{S}}-1$, where the subscripts ‘sig.’ and ‘ref.’ refer to the
signal and reference modes, respectively, and add in quadrature as the
estimate of the systematic uncertainty.

Since the protons in the signal and reference modes have different kinematic
distributions, the systematic effects due to PID do not cancel completely. The
data/MC ratio of proton PID efficiency depends on the proton momentum and
polar angle: $R{}_{\varepsilon}^{p}(p,\cos\theta)$. Such a
$R{}_{\varepsilon}^{p}$ map is determined based on an inclusive sample of
$\Lambda\\!\rightarrow\\!{}p\pi^{-}$. Following steps similar to those used
above for $K{}^{0}_{S}$ efficiency, we obtain the unfolded ($p,\,\cos\theta$)
two-dimensional distribution for protons using the ${}_{s}\mathcal{P}lot$
technique [41], and plot the $\overline{R}{}_{\varepsilon,\,{\rm
sig.}}^{p}/\overline{R}{}_{\varepsilon,\,{\rm ref.}}^{p}\\!-\\!1$ values based
on 10000 maps of $R_{\varepsilon}^{p}(p,\,\cos\theta)$. The systematic
uncertainty due to PID is obtained by adding in quadrature the mean and RMS
values of the $\overline{R}{}_{\varepsilon,\,{\rm
sig.}}^{p}/\overline{R}{}_{\varepsilon,\,{\rm ref.}}^{p}-1$ distribution.

The uncertainty due to $\eta\\!\rightarrow\\!\gamma\gamma$ reconstruction is
estimated to be 4%, considering 2% per photon according to a study of
radiative Bhabha events.

The systematic uncertainties from the $M(\Lambda_{c}^{+})$ fits for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ channels are evaluated
to be 1.8% and 2.3%, respectively, after considering two sources below. (a)
The uncertainty due to fixing the signal parameters in the fits is estimated
by randomly varying them via a multiple-dimensional Gaussian function
(including these parameters’ uncertainties and their correlation matrix from
the $M(\Lambda_{c}^{+})$ fit of truth-matched signals). We produce 1000 sets
of such signal parameters and repeat the $M(\Lambda_{c}^{+})$ fits. We take
the ratio of RMS to mean value of the distribution of fitted yield as the
relative systematic uncertainty: 0.2% for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$, 0.4% for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, and 0.2% for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$. (b) To evaluate the
potential fit bias, we perform a bias check for the fitted signal yield based
on 1000 sets of MC samples, of which the signals are randomly sampled from a
large signal MC sample and the backgrounds from the generic $B\overline{B}$
and continuum MC samples. Their sampled yields are equal to the fitted yields
in Table 1. We perform $M(\Lambda_{c}^{+})$ fits for these samples. The fitted
signal yields are plotted and fitted with a Gaussian function. The shifts of
the fitted mean values of the Gaussian functions from the corresponding input
values are assigned as systematic uncertainties: 1.9% for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$, 2.3% for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, and 0.1% for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$. The uncertainties for
signal modes and reference mode are added in quadrature, as listed in Table 2.

The systematic effects from the efficiency corrections for the
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ channels are evaluated
to be 0.6% and 0.4%, respectively, which are obtained by taking the quadratic
sum of the following sources: (a) Varying bin size: the 7$\times$7 bins are
changed to 6$\times$6 and 8$\times$8 bins for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and the
5$\times$5 bins are changed to 4$\times$4 and 6$\times$6 bins for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$. The changes of
efficiency-corrected yields, $0.2\%$ for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and $0.1\%$ for
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$, are assigned as the
systematic uncertainties. (b) To estimate the uncertainties due to the
background Dalitz plot, we shift the $M(\Lambda_{c}^{+})$ sideband region by
$\pm 5$ MeV, and repeat the efficiency correction. The resulting changes of
efficiency-corrected yields, 0.1% for both channels, are assigned as
systematic uncertainty. (c) The signal efficiency effects due to the
additional requirements in the signal mode with respect to the reference mode,
such as $p(\eta)$, $\chi^{2}_{m}(\eta)$, and $L/\sigma_{L}(K{}^{0}_{S})$, are
neglected, as the signal distributions unfolded from data using the
${}_{s}\mathcal{P}lot$ technique [41] and truth-matched signal distributions
from MC are consistent. (d) Systematic effects from the $\chi_{\rm vtx}^{2}$
requirement are considered, since the signal and reference modes have
different $\chi_{\rm vtx}^{2}$ distributions. We change the requirement to
$\chi_{\rm vtx}^{2}\\!<\\!21$ and repeat our measurement. The resulting
changes to the nominal results, 0.6% and 0.3%, are small as expected and
assigned as the corresponding systematic uncertainties. (e) The uncertainty
due to the $\pi^{0}$ veto for $\eta$ candidates in
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ is estimated by
enlarging the veto region from $\pm 12.5$ MeV/$c^{2}$ to be $\pm 15$
MeV/$c^{2}$. The resulting change on the branching fraction is 0.2%, and is
assigned as a systematic uncertainty. (f) The uncertainty due to limited MC
statistics for efficiency value is 0.1%.

The uncertainty due to the non-$K{}^{0}_{S}$ peaking background is estimated
based on the generic MC sample aforementioned. As the rate of this background
may depend on intermediate processes, we double its size, and take the
resulting ratio with the signal yield, 0.8%, as the associated systematic
uncertainty. The uncertainties on
$\mathcal{B}(K{}^{0}_{S}\\!\rightarrow\\!\pi^{+}\pi^{-})\\!=\\!{(69.20\pm
0.05)\%}$ ($\delta\mathcal{B}/\mathcal{B}=0.1\%$) and
$\mathcal{B}(\eta\\!\rightarrow\\!\gamma\gamma)\\!=\\!{(39.41\pm 0.20)\%}$
($\delta\mathcal{B}/\mathcal{B}=0.5\%$) are also considered. All uncertainties
above are added in quadrature to give an overall systematic uncertainty, as
listed in Table 2. Additionally, the uncertainty from the world average
branching fraction of the reference mode ($5.0\%$) is considered.

##  VII Summary

In summary, based on the entire dataset with an integrated luminosity of
$980~{}\mbox{\,fb}^{-1}$ collected by the Belle detector at the KEKB energy-
asymmetric $e^{+}e^{-}$ collider, we present the first observation of SCS
decay ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ with a
large statistical significance of ${>\\!10\sigma}$, and measure the branching
fractions of ${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ relative to
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}}$:

$\displaystyle\hskip-20.0pt\frac{\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}})}{\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})}$
$\displaystyle=$ $\displaystyle(1.47\pm 0.08\pm 0.04)\times 10^{-2}\,,$ (10)
$\displaystyle\hskip-20.0pt\frac{\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta})}{\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})}$
$\displaystyle=$ $\displaystyle(2.72\pm 0.06\pm 0.13)\times 10^{-1}\,,$ (11)

where the uncertainties are statistical and systematic, respectively, Using
the world average
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})=(1.59\pm
0.08)\%$ [17], we obtain the absolute branching fractions

$\displaystyle\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}})=\hskip
120.0pt$ $\displaystyle\hskip 40.0pt(2.34\pm 0.12\pm 0.06\pm 0.12)\times
10^{-4}\,,$ (12)
$\displaystyle\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta})=\hskip
129.0pt$ $\displaystyle\hskip 40.0pt(4.32\pm 0.10\pm 0.20\pm 0.22)\times
10^{-3}\,,$ (13)

where the first uncertainties are statistical, the second systematic, and the
third from the uncertainty on
$\mathcal{B}({\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}})$. The first of
these branching fractions is measured for the first time and found to be much
smaller than the theoretical prediction of $(1.9\pm 0.4)\\!\times\\!10^{-3}$
[19]. The latter one is consistent with the world average, $(4.15\pm
0.90)\\!\times\\!10^{-3}$ [17], and has a precision improved by threefold.

We also show the Dalitz plots of
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}K{}^{0}_{S}}$ and
${\Lambda_{c}^{+}\\!\rightarrow\\!{}pK{}^{0}_{S}\eta}$ after background
subtraction and efficiency correction, and find two clear structures:
$f_{0}(980)$ or $a_{0}(980)$ decaying into $K{}^{0}_{S}K{}^{0}_{S}$ and
$N^{*}(1535)$ decaying into $p\eta$. This indicates that amplitude analyses of
these decays, which will be possible using the large datasets that will become
available at BESIII [42] or Belle II [43], may help understand the nature of
such intermediate resonances.

## Acknowledgments

L. K. Li thanks Li-Sheng Geng and Ju-Jun Xie for helpful discussions on the
$N^{*}(1535)$. We thank the KEKB group for the excellent operation of the
accelerator, and the KEK cryogenics group for the efficient operation of the
solenoid.

## References

  * Xu and Kamal [1992] Q. P. Xu and A. N. Kamal, Phys. Rev. D 46, 270 (1992).
  * Cheng and Tseng [1993] H.-Y. Cheng and B. Tseng, Phys. Rev. D 48, 4188 (1993).
  * Korner and Kramer [1992] J. G. Korner and M. Kramer, Z. Phys. C 55, 659 (1992).
  * Ivanov  _et al._ [1998] M. A. Ivanov, J. G. Korner, V. E. Lyubovitskij, and A. G. Rusetsky, Phys. Rev. D 57, 5632 (1998).
  * Sharma and Verma [1999] K. K. Sharma and R. C. Verma, Eur. Phys. J. C 7, 217 (1999).
  * Zou  _et al._ [2020] J. Zou, F. Xu, G. Meng, and H.-Y. Cheng, Phys. Rev. D 101, 014011 (2020).
  * Cheng  _et al._ [2018] H.-Y. Cheng, X.-W. Kang, and F. Xu, Phys. Rev. D 97, 074028 (2018).
  * Lü  _et al._ [2016] C.-D. Lü, W. Wang, and F.-S. Yu, Phys. Rev. D 93, 056008 (2016).
  * Geng  _et al._ [2019a] C. Q. Geng, C.-W. Liu, and T.-H. Tsai, Phys. Lett. B 794, 19 (2019a).
  * Geng  _et al._ [2019b] C.-Q. Geng, C.-W. Liu, T.-H. Tsai, and Y. Yu, Phys. Rev. D 99, 114022 (2019b).
  * Brod  _et al._ [2012] J. Brod, A. L. Kagan, and J. Zupan, Phys. Rev. D 86, 014023 (2012).
  * Cheng and Chiang [2012] H.-Y. Cheng and C.-W. Chiang, Phys. Rev. D 85, 034036 (2012).
  * Li  _et al._ [2012] H.-n. Li, C.-D. Lu, and F.-S. Yu, Phys. Rev. D 86, 036012 (2012).
  * Grossman  _et al._ [2007] Y. Grossman, A. L. Kagan, and Y. Nir, Phys. Rev. D 75, 036008 (2007).
  * Grossman  _et al._ [2012] Y. Grossman, A. L. Kagan, and J. Zupan, Phys. Rev. D 85, 114036 (2012).
  * Altmannshofer  _et al._ [2012] W. Altmannshofer, R. Primulando, C.-T. Yu, and F. Yu, JHEP 04, 049 (2012).
  * Workman  _et al._ [2022] R. L. Workman _et al._ (Particle Data Group), Prog. Theor. Exp. Phys. 2022, 083C01 (2022).
  * Note [1] Throughout this paper, the charge-conjugate modes are implied.
  * Cen  _et al._ [2019] J.-Y. Cen, C.-Q. Geng, C.-W. Liu, and T.-H. Tsai, Eur. Phys. J. C 79, 946 (2019).
  * Capstick and Roberts [2000] S. Capstick and W. Roberts, Prog. Part. Nucl. Phys. 45, S241 (2000).
  * Xie and Geng [2017] J.-J. Xie and L.-S. Geng, Phys. Rev. D 96, 054009 (2017).
  * Pavao  _et al._ [2018] R. Pavao, S. Sakai, and E. Oset, Phys. Rev. C 98, 015201 (2018).
  * Zou [2010] B.-S. Zou, Nucl. Phys. A 835, 199 (2010).
  * Guo  _et al._ [2018] F.-K. Guo, C. Hanhart, U.-G. Meißner, Q. Wang, Q. Zhao, and B.-S. Zou, Rev. Mod. Phys. 90, 015004 (2018).
  * Achasov  _et al._ [2021] N. N. Achasov, J. V. Bennett, A. V. Kiselev, E. A. Kozyrev, and G. N. Shestakov, Phys. Rev. D 103, 014010 (2021).
  * Wang  _et al._ [2022] Z.-Q. Wang, X.-W. Kang, J. A. Oller, and L. Zhang, Phys. Rev. D 105, 074016 (2022).
  * Alford and Jaffe [2000] M. G. Alford and R. L. Jaffe, Nucl. Phys. B 578, 367 (2000).
  * Maiani  _et al._ [2004] L. Maiani, F. Piccinini, A. D. Polosa, and V. Riquer, Phys. Rev. Lett. 93, 212002 (2004).
  * ’t Hooft  _et al._ [2008] G. ’t Hooft, G. Isidori, L. Maiani, A. D. Polosa, and V. Riquer, Phys. Lett. B 662, 424 (2008).
  * Weinstein and Isgur [1990] J. D. Weinstein and N. Isgur, Phys. Rev. D 41, 2236 (1990).
  * Baru  _et al._ [2004] V. Baru, J. Haidenbauer, C. Hanhart, Y. Kalashnikova, and A. E. Kudryavtsev, Phys. Lett. B 586, 53 (2004).
  * [32] A. Abashian et al. (Belle Collaboration), Nucl. Instrum. Meth. A 479, 117 (2002); also see Section 2 in J. Brodzicka et al., Prog. Theor. Exp. Phys. 2012, 04D001 (2012).
  * [33] S. Kurokawa and E. Kikutani, Nucl. Instrum. Meth. A 499, 1 (2003), and other papers included in this Volume; T. Abe et al., Prog. Theor. Exp. Phys. 2013, 03A001 (2013) and references therein.
  * Lange [2001] D. J. Lange, Nucl. Instrum. Meth. A 462, 152 (2001).
  * Sjöstrand  _et al._ [2001] T. Sjöstrand, P. Edén, C. Friberg, L. Lönnblad, G. Miu, S. Mrenna, and E. Norrbin, Comput. Phys. Commun. 135, 238 (2001).
  * Brun  _et al._ [1987] R. Brun _et al._ , CERN Report No. CERN-DD-EE-84-1 (1987).
  * Barberio and Wa̧s [1994] E. Barberio and Z. Wa̧s, Comput. Phys. Commun. 79, 291 (1994).
  * Nakano [2002] E. Nakano, Nucl. Instrum. Meth. A 494, 402 (2002).
  * Feindt and Kerzel [2006] M. Feindt and U. Kerzel, Nucl. Instrum. Meth. A 559, 190 (2006).
  * [40] H. Nakano, Ph.D Thesis, Tohoku University (2014), Chapter 4, http://hdl.handle.net/10097/58814.
  * Pivk and Le Diberder [2005] M. Pivk and F. R. Le Diberder, Nucl. Instrum. Meth. A 555, 356 (2005).
  * Ablikim  _et al._ [2020] M. Ablikim _et al._ (BESIII Collaboration), Chin. Phys. C 44, 040001 (2020).
  * Kou  _et al._ [2019] E. Kou _et al._ , Prog. Theor. Exp. Phys. 2019, 123C01 (2019).