Observation of $\psi(3686) \to p \bar{p} \eta^{\prime}$ and improved measurement of $J/\psi \to p \bar{p} \eta^{\prime}$

We observe the process $\psi(3686) \to p \bar{p} \eta^{\prime}$ for the first time, with a statistical significance higher than 10$\sigma$, and measure the branching fraction of $J/\psi \to p \bar{p} \eta^{\prime}$ with an improved accuracy compared to earlier studies. The measurements are based on $4.48 \times 10^8$ $\psi(3686)$ and $1.31 \times 10^{9}$ $J/\psi$ events collected by the BESIII detector operating at the BEPCII. The branching fractions are determined to be $B(\psi(3686) \to p \bar{p} \eta^{\prime}) = (1.10\pm0.10\pm0.08)\times10^{-5}$ and $B(J/\psi \to p \bar{p} \eta^{\prime})=(1.26\pm0.02\pm 0.07)\times10^{-4}$, where the first uncertainties are statistical and the second ones systematic. Additionally, the $\eta-\eta^{\prime}$ mixing angle is determined to be $-24^{\circ} \pm 11^{\circ}$ based on $\psi(3686) \to p \bar{p} \eta^{\prime}$, and $-24^{\circ} \pm 9^{\circ}$ based on $J/\psi \to p \bar{p} \eta^{\prime}$, respectively.


I. INTRODUCTION
Quantum Chromodynamics (QCD), the theory describing the strong interaction, has been tested thoroughly at high energy. However, in the medium energy region, theoretical calculations based on first principles are still unreliable since the non-perturbative contribution is significant and calculations have to rely on models. Experimental measurements in this energy region are helpful to validate or falsify models, constrain parameters, and inspire new calculations. Charmonium states are on the boundary between perturbative and non-perturbative regimes in QCD, therefore, their decays, especially the hadronic decays, provide ideal inputs to study the QCD. The availability of very large samples of vector charmonia, produced via electron-positron annihilation, such as J/ψ and ψ(3686), makes experimental studies of rare processes and decay channels with complicated intermediate structures possible.
Among these hadronic decays, scenarios of ψ (in the following, ψ denotes either J/ψ or ψ(3686)) decaying into baryon pairs have been understood via cc annihilation into three gluons or a virtual photon [1]. But its natural extension, the three-body decays, ψ → ppP , where P represents a pseudoscalar meson such as π 0 , η, or η ′ , still need more studies since intermediate states contribute significantly here. Specific models based on nucleon and N * pole diagrams have been proposed to deal with these problems [2][3][4]. However, recent studies have focused on the final states ppπ 0 and ppη, and not so much on ppη ′ , partially due to the limited experimental measurements.
In addition, using the process ψ → ppη ′ , we are able to test the "12% rule". The ratio Q of the branching fractions of J/ψ and ψ(3686) can be written in terms of their total and leptonic widths under the assumption that the charmonium systems are non-relativistic and decay to hadrons predominantly via point-like annihilation into three gluons [6,7]: Q = B(ψ(3686)→ggg) B(J/ψ→ggg) = Γ(ψ(3686)→e + e − )·Γ(J/ψ) Γ(J/ψ→e + e − )·Γ(ψ(3686)) = (12.2 ± 2.4)% [8]. This relation was extended to exclusive processes later, ignoring other factors associated with each exclusive mode such as multiplicity and phase space factors. Although the "12% rule" has been confirmed experimentally for many decay modes, severe violation has been found in several channels [8]. Many theoretical explanations [9] have been proposed to explain the violation of the "12% rule", but none is satisfactory.

II. BESIII DETECTOR AND DATA SAMPLES
The BESIII detector, described in detail in Ref. [15], has a geometrical acceptance of 93% of 4π solid angle. It can be subdivided into four main sub-detectors. A helium-based multi-layer drift chamber (MDC) determines the momentum of charged particles, traveling through a 1 T (for J/ψ sample 0.9 T in 2012) magnetic field, with a resolution 0.5% at 1 GeV/c, as well as the ionization energy loss (dE/dx) with a resolution better than 6.0% for electrons from Bhabha scattering. A time-of-flight system (TOF) made of plastic scintillators with a time resolution of 80 ps (110 ps) in the barrel (end caps) is used for particle identification (PID). An electromagnetic calorimeter (EMC) consisting of 6240 CsI(Tl) crystals measures the energies of photons with a resolution of 2.5% (5.0%) in the barrel (end caps) at 1 GeV, and their positions with a resolution of 6 mm (9 mm) in the barrel (end caps). A muon counter (MUC) based on resistive plate chambers with 2 cm position resolution provides information for muon identification.

III. EVENT SELECTIONS AND BACKGROUND ANALYSIS
Charged tracks are reconstructed from hits in the MDC. For each track, the polar angle must satisfy | cos θ| < 0.93 and the point of closest approach to the interaction point must be within ±1 cm in the plane perpendicular to the beam and ±10 cm along the beam direction. The TOF and dE/dx information is combined to calculate PID likelihoods for the pion and proton hypotheses, and the PID hypothesis with the largest likelihood is assigned to the track.
Photons are reconstructed from isolated electromagnetic showers in the EMC. The angle between the directions of any cluster and its nearest charged track must be larger than 10 or 30 degrees to pion or (anti-)proton tracks, respectively. The efficiency and energy resolution are improved by including the energy deposited in nearby TOF counters. A photon candidate must deposit at least 25 MeV (50 MeV) in the barrel (end caps) region, corresponding to an angular coverage of | cos θ| < 0.80 (0.86 < | cos θ| < 0.92). The timing information obtained from the EMC is required to be 0 ≤ t EMC ≤ 700 ns to suppress electronic noise and beam backgrounds unrelated to the event.
Signal candidates must have four charged tracks identified as p,p, π + and π − , as well as at least one (two) photon(s) for the η ′ → γπ + π − (η ′ → ηπ + π − ) mode. The events with 920 < M γπ + π − (ηπ + π − ) < 1000 MeV/c 2 are accepted for further analysis, where M γπ + π − (ηπ + π − ) is the invariant mass of γπ + π − or ηπ + π − , respectively. To improve the resolution and to suppress backgrounds, a kinematic fit to all final state particles with a constraint on the initial e + e − fourmomentum is performed. In addition, for the η ′ → π + π − η mode, the invariant mass of the two photons is constrained to the nominal mass of η. The χ 2 of the kinematic fit for each decay mode is required to be less than an optimized value obtained by maximizing the figure of merit S/ √ S + B, where S is the number of signal events from a signal MC sample normalized to the preliminary measurements and B is the number of background events in the η ′ signal region obtained from inclusive MC samples. When there are more photon candidates than required in an event, we loop over all possible combinations and keep the one with the smallest kinematic fit χ 2 . After the kinematic fit, the η ′ signal region is defined as 948.2 < M γπ + π − (ηπ + π − ) < 967.4 MeV/c 2 , while the sideband regions are defined as 932.4 < M γπ + π − (ηπ + π − ) < 942.0 MeV/c 2 and 974.0 < M γπ + π − (ηπ + π − ) < 983.6 MeV/c 2 .
To remove background events, we apply the following requirements: Here, M rec γ , M rec π + π − , and M rec γγ are the recoil mass of γ, π + π − , and γγ, while m χcJ and m J/ψ are the nominal χ cJ and J/ψ masses [21], respectively. The mass window for each requirement is determined based on the exclusive MC simulation.
The backgrounds from ψ(3686) and J/ψ decays are studied with inclusive MC samples. For ψ(3686) → ppη ′ with η ′ → γπ + π − , even after the χ cJ mass window requirements, the main remaining background is the decay ψ(3686) → γχ cJ with χ cJ → ppπ + π − , which has the same final state as the signal process. A study of MC simulated events shows that the M γπ + π − distribution from ψ(3686) → γχ cJ → γppπ + π − is smooth. Therefore, its contribution can be easily determined in a fit. For the other three decay modes, ψ(3686) → ppη ′ with η ′ → ηπ + π − , J/ψ → ppη ′ with η ′ → γπ + π − and η ′ → ηπ + π − , there are no dominant background processes, but many decay channels with a small contribution each. The backgrounds from the continuum process e + e − → qq are studied with data samples taken at √ s = 3.080 and 3.650 GeV. The background level is found to be very low, and the background events do not peak in the signal region.
The M γπ + π − and M ηπ + π − distributions of the events that pass all selection criteria are shown in Fig. 1. Peaks originating from η ′ decays are observed. Figure 2 shows the Dalitz plots of the events in the η ′ signal region, and Figs. 3 and 4 show the invariant mass projections, where the side band backgrounds have been subtracted. Based on these plots, no obvious intermediate structures in invariant mass of pη ′ ,pη ′ , or pp are observed.

IV. SIGNAL YIELDS AND BRANCHING FRACTIONS
To determine the branching fractions, simultaneous unbinned maximum likelihood fits to the γπ + π − and ηπ + π − invariant mass spectra are performed for the ψ(3686) data and  for the J/ψ data. The signal shape is represented by the MCsimulated η ′ mass distribution, convolved with a Gaussian function with free mean and width to account for the mass and resolution difference between data and MC simulation. The background is parameterized as a second-order Chebyshev polynomial with free parameters. In the simultaneous fit, the ratio of the number of η ′ → γπ + π − events to that of η ′ → ηπ + π − events is fixed to where ǫ η ′ →γπ + π − and ǫ η ′ →π + π − η are the global efficiencies for each η ′ decay mode. Due to differences in tracking and PID efficiencies between data and MC simulation for protons and anti-protons, the MC-determined global efficiencies are corrected by multiplying factors 1.030 (1.038) and 0.980 (0.984) for tracking and PID, respectively, for ψ(3686) → ppη ′ (J/ψ → ppη ′ ). These correction factors are ratios of efficiencies between data and MC simulation obtained by studying the control samples ψ → ppπ + π − , where the efficiencies are weighted according to the distributions of transverse momentum (for tracking) or momentum (for PID) of protons and anti-protons.

V. SYSTEMATIC UNCERTAINTIES
The systematic uncertainties mainly come from the MDC tracking, photon and η reconstruction, PID, kinematic fit, mass windows, branching fractions of the decay modes used to reconstruct the η ′ , the number of ψ decays, fitting procedure, and the physics model used to determine the efficiency. All the contributions are given in Table II. The overall systematic uncertainties are obtained by adding all systematic uncertainties, taking the correlations into account.
The uncertainty in the MDC tracking efficiency for each pion is estimated with the control sample ψ(3686) → π + π − J/ψ, and a 1.0% systematic uncertainty per pion is obtained [24]. This gives a total of 2.0% for each decay mode. The tracking efficiencies of protons and anti-protons are studied with the control sample ψ → ppπ + π − . The MC efficiencies for the signal processes are corrected using the results from the control samples, and the uncertainties  , and pp (c) for J/ψ → ppη ′ with η ′ → γπ + π − , and those of pη ′ (d),pη ′ (e), and pp (f) for J/ψ → ppη ′ with η ′ → π + π − η. The dots with error bars show background subtracted data, and the red lines are the corresponding distributions from signal MC.
The uncertainty in the PID efficiency for pions is estimated to be 1.0% per pion [27], and the total PID uncertainty for two pions is 2.0% for each decay mode. The efficiencies of the proton and the anti-proton identification are studied with the control samples ψ → ppπ + π − . The MC efficiencies for the signal processes are corrected using the results from the control samples, and the uncertainties of the corrections are taken as the systematic uncertainties, which are 0.5% (0.6%) per proton (anti-proton) for the ψ(3686), and 0.6% (0.8%) per proton (anti-proton) for the J/ψ samples. Assuming the uncertainties of proton and anti-proton are totally correlated, the PID uncertainty of a proton and anti-proton pair is calculated to be 1.1%(1.4%) for the ψ(3686) (J/ψ) data samples. The PID efficiencies of pion and proton are independent, and the total PID uncertainty is determined to be 2.3% (2.5%) for the ψ(3686) (J/ψ) samples.
The uncertainty due to the mass windows used to veto background events originates from the differences in the mass resolutions between data and MC simulation. We repeat the analysis by enlarging or reducing the mass window. The largest difference is used as an estimate of the corresponding systematic uncertainty. The η ′ signal region and sideband regions are not used to veto the background events, so they have no effect on the branching fraction determination. The uncertainties due to different mass windows are considered to be independent, so we add them in quadrature. For the decay ψ(3686) → ppη ′ with η ′ → γπ + π − (η ′ → ηπ + π − ), the uncertainty is 4.9% (1.8%). The uncertainty for the masswindow selection for J/ψ → ppη ′ is found negligible.
The fit range, signal shape, and background shape are considered as the sources of the systematic uncertainty related with the fit procedure. In the nominal fit, the mass range is [0.90, 1.04] GeV/c 2 , and we repeat the fit by changing the range by ±10 MeV/c 2 . The largest change in the final result is taken as the uncertainty due to the fit range, which is 1.4% and 0.6% for ψ(3686) → ppη ′ and J/ψ → ppη ′ , respectively. For the signal shape, we change the nominal shape to a double-Gaussian or a Breit-Wigner function convolved with a Gaussian function, and the largest difference from the nominal result is taken as the uncertainty of the signal shape, which is 3.9% and 1.2% for ψ(3686) → ppη ′ and J/ψ → ppη ′ , respectively. We replace the background shape with a firstorder Chebyshev or second-order polynomial function, the largest differences from the nominal result, 2.7% and 1.3% for ψ(3686) → ppη ′ and J/ψ → ppη ′ , respectively, are taken as systematic uncertainties.
The signal MC sample is generated assuming pure phase space distribution, in which possible intermediate states and non-flat angular distributions are ignored. Although no strong structure is visible in the Dalitz plots shown in Fig. 2, the phase space MC does not provide a very good description of the data, as shown in Fig. 3 and Fig. 4, which results in a large systematic uncertainty, especially for the J/ψ decay modes. For ψ(3686) → ppη ′ , since the statistics are limited, we weight each MC-generated phase space event by the M pp distribution, and a difference of 1.5% in the efficiency between the nominal and weighted MC samples is taken as the uncertainty. For J/ψ → ppη ′ , we re-generate signal MC events based on BODY3 [30], a data-driven MC generator, and a difference of 3.4% in the efficiencies is taken as the uncertainty.
Our results for the branching fractions of ψ(3686) → ppη ′ and J/ψ → ppη ′ result in the ratio B(ψ(3686)→ppη ′ ) B(J/ψ→ppη ′ ) = (8.7 ± 1.0)%, where the common uncertainties have been canceled. Even though the ratio is in reasonabe agreement with 12%, we note that the kinematics of the two processes are very different, and the "12% rule" may be too naive in this case. The phase space ratio is Ω ψ(3686)→ppη ′ /Ω J/ψ→ppη ′ = 8.13, if any possible intermediate structure is ignored. Furthermore, the Dalitz plots of J/ψ and ψ(3686) decays, shown in Fig. 2, indicate that many events in ψ(3686) decays, possibly via N * N + c.c. intermediate states with pη ′ orpη ′ mass greater than 2.13 GeV/c 2 , are not kinematically possible in J/ψ decays. Taken these factors into account, the Q value is suppressed a lot, implying that the "12% rule" is violated significantly.