Triangle singularities in ${J/\psi\rightarrow\eta\pi^0\phi}$ and ${\pi^0\pi^0\phi}$

The BESIII Collaboration recently reported the observation of the $a_0(980)^0-f_0(980)$ mixing in the isospin breaking decay $J/\psi\to \eta\pi^0\phi$. In the Dalitz plot for that decay with the $\eta$ reconstructed from two photons, there is a band around $1.4$~GeV on the $\pi^0 \phi$ distribution. In general, this peak can be due to a resonance or a kinematic effect. In this paper, we study the effects of a set of $K^*K\bar K$ triangle diagrams, and show that due to triangle singularities such diagrams can lead to a peak around 1.4~GeV in the $\pi^0\phi$ invariant mass distribution. The Dalitz plot induced by such a mechanism has a feature consistent with the BESIII observation, namely events along the band accumulate at both ends close to the Dalitz plot boundary. The effect of the same mechanism on the $J/\psi\to \pi^0\pi^0\phi$ and $J/\psi \to \eta\pi^0 K^+K^-$ decays are also investigated. We suggest to take more data for the $J/\psi\to \eta\pi^0\phi\to \eta\pi^0 K^+ K^-$ and check whether the structure around $1.4$~GeV persists for the $K^+K^-$ invariant mass away from the $\phi$ mass region. This is crucial for understanding whether the band is due to triangle singularities or due to a resonance. Were it the latter, the band should remain while it would not if it is due to the former.


I. INTRODUCTION
The BESIII Collaboration reported observations of the a 0 (980) 0 − f 0 (980) mixing in the processes J/ψ → φf 0 (980) → φa 0 (980) 0 → φ π 0 η and χ c1 → a 0 (980) 0 π 0 → f 0 (980)π 0 → π + π − π 0 [1,2]. The isospin of the a 0 (980) 0 is 1 and that of the f 0 (980) is 0, and thus their mixing breaks isospin. Both of the two resonances couple strongly to the KK meson pair [3], so that the mixing can happen through their coupling to the KK intermediate states and the mass difference between the charged and neutral kaons gives the isospin breaking. Then the mixing probability depends crucially on the coupling strengths of the a 0 (980) and f 0 (980) to the KK. Such a mixing mechanism was first proposed in the late 1970s [4], and was suggested to contain important information to clarify the nature of these two mesons. That is because the effective couplings of these mesons to the kaons can be related to their probabilities to be KK molecules (see Refs. [5][6][7][8] for discussions based on extensions of Weinberg's compositeness relations [9]). The J/ψ → ηπ 0 φ reaction as a probe of the a 0 (980) 0 − f 0 (980) mixing was suggested in Ref. [10], and was further analyzed using unitarized chiral approaches in Refs. [11,12].
In the Dalitz plot of the J/ψ → ηπ 0 φ process measured by the BESIII Collaboration (see Fig. 1 of Ref. [2]), it is clear that there is a peak near the KK threshold on the ηπ 0 distribution which can be interpreted as the a 0 (980) 0 − f 0 (980) mixing. In addition, there is also a clear accumulation of events for the π 0 φ invariant mass being around 1.4 GeV, which would be a peak if the Dalitz plot is projected to the π 0 φ invariant mass distribution. In general, this peak can be due to a resonance with a mass about 1.4 GeV (isospin I = 1 or I = 0) or a kinematic effect. From the first point of view, because the S-wave π 0 φ has J P C = 1 +− and this decay breaks isospin symmetry, that peak may be related to the isovector resonance claimed in Ref. [13] or the isoscalar h 1 (1380) whose mass is about 1.41 GeV [3,14,15]. 1 From the second point of view, because this is an isospin breaking process and that peak position about 1.4 GeV is near the threshold of K * K , similar to the enhancement of the isospin violation by triangle singularity as seen in the η(1405/1475) → π 0 f 0 (980) [19][20][21][22][23], triangle diagrams involving K * K K intermediate states can also lead to a peak and a significant isospin breaking enhancement here (see also Sec. VI.A.4 of Ref. [24] and references therein for the discussion of the possible role of triangle singularity on the h 1 (1380)). Triangle singularity is the leading Landau singularity [25] of a triangle diagram, and it depends crucially on the kinematics. To be more specific, when all of the three intermediate particles in a triangle diagram are on shell, moving collinearly, and all of the interaction vertices satisfy the energy-momentum conservation [26], the physical amplitude has a logarithmic triangle singularity leading to a peak in the invariant mass distributions (see also Ref. [27] for an intuitive picture of triangle singularity with a reformulation based on a diagrammatic approach). The production of π 0 φ in the reaction π − p → π 0 φn from the K * K K triangle singularity mechanism has been studied in Ref. [28]. In the J/ψ → ηπ 0 φ reaction, we focus on the diagrams that can lead to a peak around the K * K threshold, about 1.4 GeV, in the π 0 φ distribution, and thus consider the triangle diagrams shown in Fig. 1. The triangle singularity for each of these diagrams is at the physical boundary when the kinematics is such that the following processes happen: first, the J/ψ decays into ηK * K (the charge-conjugatedK * K diagrams are also included, see Fig. 1), and the K * decays into πK subsequently; then the K moves in the same direction with theK and catches up with it, and the KK pair finally forms the φ. Triangle singularity effects of the K * KK loops have been studied for processes with the KK pair forming an f 0 (980) or a 0 (980) instead of a φ in Refs. [19][20][21][22][23][29][30][31][32][33][34].
This paper is organized as follows. In Section II, we set up the formalism for calculating the decays J/ψ → ηπ 0 φ, J/ψ → π 0 π 0 φ, and J/ψ → ηπ 0 KK. For the last process, the KK can be either isovector or isoscalar and the isoscalar KK final state interaction (FSI) is taken into account using the inverse amplitude method. As will be shown in Section III, by considering the triangle diagrams shown in Fig. 1 for the J/ψ → ηπ 0 φ, a peak appears in the π 0 φ invariant mass distribution around 1.4 GeV by virtue of triangle singularity, which might be the physics behind the band at m π 0 φ 1.4 GeV in the Dalitz plot measured by the BESIII Collaboration as they have similar gross features. The results for the π 0 φ distribution of the J/ψ → π 0 π 0 φ reaction will also be shown, and the peak induced by the triangle diagrams is much broader since this is an isospin symmetry conserving reaction. A brief summary is given in Section IV. Some decay amplitudes and a derivation of the generic n-body phase space are relegated in Appendix A and Appendix B, respectively.
The diagrams we consider for the J/ψ → ηπ 0 φ are shown in Fig. 1. To obtain the decay amplitude, we use the following effective Lagrangian for the J/ψ → ηK * K , K * → πK, and KK → φ vertices: L int = L JV P P + L V P P , where g 1 and g 2 are coupling constants, and J µ is the field operator of the J/ψ, which as a charmonium is a light-flavor SU(3) singlet state. The light pseudoscalar meson octet and vector meson nonet matrices are denoted by P and V µ , respectively, where the ideal mixing between the ω and φ is assumed such that the valence quarks for them are (uū + dd)/ √ 2 and ss, respectively. The V P P Lagrangian 2 in Eq. (1) may be obtained The factor 2 is to take account of the charge conjugated contributions, i.e., the amplitudes for Fig. 1 (a, c) give the same contribution as those for Fig. 1 (b, d). As one can see from Eq. (4), the amplitudes M C and M N differ only by the masses of the intermediate mesons K and K * , and in the isospin limit, which is realized by using identical masses for the mesons in the same isospin multiplet, M C and M N in Eq. (5) would exactly cancel each other. We calculate this amplitude in the center-of-mass (CM) frame of the π 0 and φ pair, and choose the π 0 momentum direction as the z-direction.
The differential width of this process can be written as (see Appendix B for details of the phase space factor) where the quantities marked with * are evaluated in the CM frame of particle 2 and particle 3 in the final states, and the momentum (|p 1 |, Ω 1 ) is the momentum of particle 1 in the rest frame of the decay particle. Since we are interested in the invariant mass distribution of π 0 φ, we choose η as particle 1, π 0 as particle 2, and φ as particle 3. Then the differential width is given by with spin summing over polarizations of J/ψ and φ and the amplitude M J/ψ→ηπ 0 φ given in Eq. (5).
In the experimental measurement, the φ is reconstructed from its decays into the K + K − final states with the K + K − invariant mass in the region [m φ − 10 MeV, m φ + 10 MeV]. For a direct comparison with the experimental data, we need to take into account the subsequential decay of the φ into K + K − which means the ηπ 0 K + K − in the final state. In this case, since the K * decay into Kπ with an almost 100% branching fraction, tree diagrams should be included in addition to the triangle loop ones. The diagrams for this process considered here are shown in Fig. 2. The triangle diagrams involve the KK → KK final state interaction.
We consider the P -wave KK scattering amplitude with quantum numbers (I, J) = (0, 1) and (1, 1). For simplicity, only the contribution of the φ meson is considered in the (0, 1) channel. For the KK (I, J) = (1, 1) amplitude, we employ the inverse amplitude method as developed in Ref. [41], which is a unitary extension of chiral perturbation theory and can describe the meson-meson scattering data up to 1.2 GeV. In this approach, the unitarized T -matrix for a given partial wave can be written as where T 2 is the partial wave scattering amplitude from the leading order chiral Lagrangian, T 4 is the polynomial tree-level amplitude from the next-to-leading order chiral Lagrangian, and G is a diagonal matrix given by the loop integrals with two meson propagators. One can find the specific form of the T -matrix in Ref. [41]. Here, we use the isospin symmetric KK → KK amplitude to concentrate on the effect of triangle singularities. The amplitudes for the KK → KK FSI with quantum numbers (I, J) = (0, 1) and (1, 1) are written as: where m φ = 1019.5 MeV and Γ φ = 4.2 MeV are the mass and width of the φ meson, respectively [3], p φ = p out K + p out K , and θ is the angle between the 3-momenta of K mesons in the loop and final state in the KK CM frame. The momenta p in K and p out K (p in K and p out K ) denote the K (K) momenta in the loop and final state, respectively, and q in , q out are the magnitudes of the corresponding CM momenta. The T -matrix T 11 is a two-channel matrix with the quantum numbers (I = 1, J = 1), whose elements (T 11 ) 11 , (T 11 ) 12 and (T 11 ) 22 correspond to the P -wave scattering amplitudes for ππ → ππ, ππ → KK and KK → KK, respectively. The purpose of our calculation is to study the influence of triangle singularity on the π 0 KK invariant mass distribution, so we focus on the kinematic region where triangle singularity occurs, i.e. for all the intermediate particles being on-shell. We approximate the momentum q in in the denominator of Eq. (10) for the KK scattering vertex with the on- Before we apply this formula, we need to transform the particle basis from the isospin space to the charged space: where the phase convention |K − = − |I = 1/2, I z = −1/2 , the same as that in Ref. [41], is used. Then, the J/ψ → ηπ 0 K + K − amplitudes of the tree and triangle loop diagrams shown in Fig. 2 can be written as follows: Finally, as shown in Appendix B, the differential width of this 4-body decay process can be written as where the momenta (|p 1 |, Ω 1 ), (|p 2 |, Ω 2 ) and (|p 3 |, Ω 3 ) are evaluated in the rest frame of the decaying particle, the particles 2 + 3 + 4, and the particles 3 + 4, respectively. We choose η as particle 1, π 0 as particle 2, K as particle 3 andK as particle 4. Then the double differential distribution is given by C. J/ψ → π 0 π 0 φ Furthermore, we can consider the J/ψ → π 0 π 0 φ process with little additional effort because we only need to replace the η by the π 0 in Fig. 1. The difference is that the J/ψ → ηπ 0 φ breaks isospin symmetry while this process is isospin conserving. The J/ψ → ππφ reaction has been studied in Refs. [38,39,42,43] using the S-wave ππ final state interactions to explain the structures in the ππ invariant mass distributions of the J/ψ → φππ and J/ψ → ωππ [44][45][46].
From the Lagrangian in Eq. (1), we can obtain the following J/ψ → π 0 K * K amplitudes: Due to the difference of the signs from those in Eq. (3) for the J/ψ → ηK * K case, the total amplitude of the J/ψ → π 0 π 0 φ is given by One sees that the charged and neutral intermediate states give the same contributions in the isospin limit as expected since the transition preserves isospin symmetry.

A. General discussion
Triangle singularity at the physical boundary can be easily obtained by solving the equation derived in Ref. [27] which reads (see Fig. 3 for the particle assignment) where 23 ) are the energy and the magnitude of the 3-momentum of particle m 2 in the CM frame of (m 2 , m 3 ), β is the magnitude of the velocity of (m 2 , m 3 ) system in the CM frame of (m 1 , m 2 ), and γ = 1/ 1 − β 2 is the Lorentz boost factor from the (m 2 , m 3 ) CM frame to the (m 1 , m 2 ) CM frame.
For the triangle singularity to be on the physical boundary, the physical conditions of which are given by the Coleman-Norton theorem [26], the above equation needs to have a real solution, and all the arguments of the involved square root functions need to be positive.
Physically, this means that all the intermediate particles go on shell, moving collinearly, and the m 3 particle from the decay of m 1 should be fast enough to catch up with the m 2 from the decay of M , and then particles m 2 and m 3 interact like a classical process producing external particle(s) with an invariant mass m 23 . For more discussions, we refer to, e.g., Refs. [47,48].
We can compute where the logarithmic singularities of the charged and neutral K * KK Ref. [3] for the masses), and the latter is 2 MeV above the K * 0K 0 threshold at 1393.6 MeV.
Therefore, we will see two sharp peaks in the π 0 φ invariant mass distribution, 5 , and each of the peak has a cusp at the K * K threshold on its left shoulder. When the K * width is taken into account, the sharp peaks will be smeared to a smooth and much broader peak with a width dictated by the K * width. For a detailed study of the width effects on triangle singularity, we refer to Ref. [23,49].
In Fig. 4, we show the π 0 φ invariant mass distribution given by Eq. (7) with and without considering the K * width. The behavior of the curves is exactly as expected in the above. Since the difference between the two real parts is much smaller than twice the imaginary part, the two singularities lead to only a single broad peak shown as the dashed curve. Thus, the two peaks from the charged and neutral K * K K loops with isospin mass splitting turn into a single peak in the π 0 φ distribution with the inclusion of the K * width.
The BESIII Collaboration reported the Dalitz plot distribution for the decay J/ψ → ηπ 0 φ in the (m 2 ηπ 0 , m 2 π 0 φ ) plane [2]. Here, let us discuss the triangle singularity contribution to the Dalitz plot. Note that the a 0 (980)/f 0 (980) resonances seen in the Dalitz plot (along the π 0 η direction) in Ref. [2], as well as other resonances listed in the Review of Particle Physics [3], are not considered here, because they cannot produce any nontrivial narrow structure in 6 The structure at 1393.6 MeV is the K * 0K 0 threshold cusp, and the K * + K − threshold cusp is invisible in the plot because it is too close to (only 0.4 MeV below) the triangle singularity. the π 0 φ invariant mass distribution (as can be anticipated from the Dalitz plot projection), which is the focus here. The mechanism discussed here is only important in a small energy region for m π 0 φ not far from the K * K threshold. This means that the extension of our calculation to the whole phase space would not be adequate. However, as we discuss below the main feature of the triangle singularity contribution to the Dalitz plot is consistent with what was observed by the BESIII Collaboration.
With Eq. (7), the Dalitz plot distribution is given by with the amplitude M J/ψ→ηπ 0 φ given in Eq. (5). Considering the K * width, the resulting Dalitz plot is shown in Fig. 5. One can see that the peak shown as the dashed line in Except that in the low m π 0 η region there is a large contribution from the f 0 (980) and a 0 (980) The m π 0 φ distribution for the J/ψ → ηπ 0 φ including the K * width is shown in the left panel of Fig. 6, which is the same as Fig. 4 with the K * width up to an overall normalization factor. One can see a clear peak around 1.39 GeV with a width about 50 MeV. Because of the isospin breaking, the charged and neutral loops largely cancel each other in the region outside the peak, making the peak a prominent structure. The dependence on the scale µ which appears in the regularization of the UV divergence is also checked: the results with µ = 0.2, 1.2, and 2.2 GeV for J/ψ → ηπ 0 φ are shown as dashed, solid and dot-dashed lines, respectively, in the left panel of Fig. 6. As one can see, the peak in the case of J/ψ → ηπ 0 φ is affected little by changing µ. This is because the decay breaks isospin, and the UV divergence part gets largely canceled in the difference between the contributions of the charged and neutral meson loops, see Eq. (5). Such a cancellation has been discussed in Ref. [50] studying the charmed-meson loop contribution to the isospin breaking decay ψ → J/ψπ 0 : it is shown that at the leading order of nonrelativistic effective field theory, the UV divergence in the neutral and charge meson loops cancel with each other, leaving a finite piece depending on the masses of intermediate particles.
C. J/ψ → π 0 π 0 φ Differently, the J/ψ → π 0 π 0 φ reaction is isospin-symmetry allowed, and the charged and neutral loops add up to give the final result. Thus, one would expect the peak in the π 0 φ invariant mass distribution due to triangle singularities to be much more modest than that in the J/ψ → π 0 ηφ case, and the result should have a large scale dependence. This is indeed the case as can be seen from the right panel of Fig. 6. Nonetheless, the triangle singularities in this process have the same origin as the prominent one for the J/ψ → π 0 ηφ, and deserve to be studied in more detail to reveal another aspect of the J/ψ → π 0 π 0 φ process in addition to the ππ distribution that has been extensively studied both experimentally [46] and theoretically [38,39,42,43]. the singularities will damp quickly.
Therefore, in order to reveal the origin of the band in the Dalitz plot at m π 0 φ ∼ 1.4 GeV in the BESIII data [2], i.e., whether it is indeed due to triangle singularities, we can study the J/ψ → ηπ 0 KK reaction and investigate the correlation between the KK invariant mass and the triangle singularity effects. In Fig. 7   a clear peak at m π 0 K + K − ∼ 1.4 GeV due to isospin breaking, which is most prominent at m K + K − 1.02 GeV due to the presence of the φ resonance. When the KK invariant mass is away from the φ resonant region, the peak along m π 0 K + K − also becomes much less evident in the lower left plot. Were the band at m π 0 φ ∼ 1.4 GeV in the Dalitz plot observed by the BESIII Collaboration due to a resonance decaying into π 0 φ, the same resonance would also be able to decay into π 0 K + K − , and its signal should always be there no matter what value the K + K − invariant mass takes (as long as there are enough data).
In the left panel of Fig. 8 diagram is approximately given by the tree-level amplitude multiplied by a phase factor with the phase given by the scattering phase shift of the two particles that rescatter. As a result, the invariant mass distribution does not have a sharp peak, but behaves as that of the tree-level t-channel exchange diagram. One sees that the Schmid theorem is not applicable here for two reasons: (1) the K * meson has a width of about 50 MeV; (2) for the reaction of interest J/ψ → ηπ 0 K + K − , only the K + K − pair appears in the final state while both the K + K − and K 0K 0 appear in the loops.

IV. SUMMARY
In summary, by considering the triangle diagrams with K * KK intermediate particles for the isospin breaking process J/ψ → ηπ 0 φ shown in Fig. 1, we show that there appears a band in the Dalitz plot which is located at the π 0 φ invariant mass of around 1.4 GeV by virtue of triangle singularities. Events along the band accumulate around the two ends close to the Dalitz plot boundary. Such a feature is consistent with what was observed by the BESIII Collaboration with the η reconstructed from two photons in Ref. [2].
In order to check whether this is the genuine mechanism behind the experimental observation, we further study the processes J/ψ → π 0 π 0 φ and J/ψ → ηπ 0 K + K − . These processes are dominated by the isospin conserving contribution. It is shown that the triangle singular-ity mechanism has a much more modest effect for the J/ψ → π 0 π 0 φ and leads to a broader bump. As for the J/ψ → ηπ 0 K + K − , our results show that the triangle singularity effect is most evident with the K + K − invariant mass is in the φ resonance mass region where isospin breaking plays an important role. Thus, by taking events with m K + K − away from that region, we expect the band at m π 0 K + K − ∼ 1.4 GeV in the Dalitz plot would become weaker and eventually invisible if K + K − is sufficiently away from the φ mass region. On the other hand, if the band in the Dalitz plot is due to a resonance, which should also be able to decay into π 0 K + K − , one expects the band remains independent of the K + K − invariant mass.
Therefore, in order to reveal the physics behind the band at m π 0 φ ∼ 1.4 GeV observed by the BESIII Collaboration, we suggest experimentalists to take more data for the J/ψ → ηπ 0 φ to make the Dalitz plot more clear, and to check whether the structure persists in other With the technique of tensor decomposition of loop integrals, I µ , I µν , I µνρ , and I µνρσ can be rewritten as given in Ref. [56], I µνρσ =(g µν g ρσ + g µρ g νσ + g µσ g νρ )C 0000 Using these I id and I µ 1 µ 2 ...µn id , the loop amplitude M µν id is written as Since both the K * Kπ 0 and φKK vertices are P -wave and the K * propagator has a momentum dependent numerator, the amplitude of the triangle diagram is UV divergent. In all calculations of this paper, we use dimensional regularization and with the MS subtraction scheme for the subtraction of the UV divergence, and unless specified we take µ = 1.2 GeV for the scale introduced in dimensional regularization.
In Eq. (A8), for example, a term g ρσ I ρσ should be understood as The terms g ρσ I µρσ and g ρσ I µνρσ also need similar manipulations.
Appendix B: General phase space formula for 1 to n-body decay The partial decay rate of a particle of mass m into n bodies in its rest frame is given in terms of the Lorentz-invariant matrix element M by where dΦ n is an element of the n-body phase space given by To get the decay formula, we need to integrate the order-4 Dirac δ-function in Eq. (B2) which represents the four-momentum conservation. For instance, we can integrate the mo-menta of the last particle first: where we have defined an order-n ϑ-function as follows: Substituting Eq. (B3) into Eq. (B2), we can get So far, we still have an order-1 Dirac δ-function which represents the mass-shell condition of the last particle. Just like what we did in Eq. (B3), to eliminate the remaining order-1 Dirac δ-function, we can integrate over the momentum of the (n − 1)-th particle: Before we continue to simplify the expression, let us make some definitions to facilitate subsequent calculations: m (k) = p 2 (k) , E (k) = p 0 (k) = p 2 (k) + m 2 (k) .
In addition, we use X * (l) to represent the specific form of physical quantity X(X can be momentum p, energy E or mass m etc.) in the CM frame of the particle system l, l+1, · · · , n.
Furthermore, for simplicity, we redefine X * (l) l and X * (l) (l) as X ( * ) l and X * (l) , respectively. Since the element of final state phase space of each particle is Lorentz invariant, one can calculate each part in any reference frame. Then, let us come back to the integration of the phase space for the (n − 1)-th particle, and we calculate this part in the CM frame of (n − 1)-th and n-th particle system: (B7) Here we note that in the CM frame of the (n − 1)-th and n-th particle system we have where we have integrated out the last order-1 Dirac δ-function. Then substituting Eq. (B10) into Eq. (B6), we can get: dΦ n (p 1 , · · · , p n ) = 1 2(2π) 6 It is worth noting that the effect of the θ and ϑ functions is to limit the integrating range to the physical region, and it is easy to check that in the physical region,

(B13)
This expression is sufficient to get the result, but it is better to express it in terms of invariant masses. Then, let us consider the relation between energy and invariant masses for each particle: k ) is the momentum of particle k in the CM frame of the particle system k, k + 1, · · · , n, and m (k) is the invariant mass of the particle system k, k + 1, · · · , n.