Next-to-leading order QCD corrections to fully charm tetraquark radiative decays

Introduction. Fully charm tetraquarks provide a rare opportunity to study how QCD binds four quarks within a nonrelativistic framework. Their heavy constituents permit an expansion in the relative velocity, while the absence of valence light quarks limits the role of light-meson exchange [1, 2, 3]. These features make their spectrum a test of short-distance color interactions and nonperturbative four-quark binding.

The LHCb Collaboration first observed a pronounced structure near 6.9 GeV in the di-$J/\psi$ spectrum [4]. ATLAS and CMS later confirmed the $X(6900)$ and reported additional structures near 6.6 and 7.1 GeV [5, 6]. An NRQCD study of hadroproduction had previously suggested spin 2 for the observed di-$J/\psi$ resonances [7]. Recent CMS analyses favor a common $J^{PC}=2^{++}$ assignment for the $X(6600)$, $X(6900)$, and $X(7100)$ [8, 9]. This assignment supports their interpretation as radial excitations [10, 11], but leaves their microscopic origin unsettled. Compact tetraquarks, hadronic molecules, and threshold effects can produce similar spectra [12, 13, 14, 15, 16, 17, 18, 19]. Discriminating among these possibilities requires complementary decay observables.

For the observed structures above the di-$J/\psi$ threshold, radiative decays into a photon and charmonium offer complementary information on their nature. The radiative process proceeds through short-distance $c\bar c$ annihilation and probes the quark contact density [20, 21], while the di-$J/\psi$ spectrum can also be affected by threshold effects and rescattering [19]. In NRQCD, the contact densities appear as long-distance matrix elements (LDMEs) of definite color-spin configurations, factorized from perturbative short-distance coefficients [22]. Together with the di-$J/\psi$ channel, radiative decays can also help constrain the spin and nature of the observed structures. Radiative width ratios have already been used to probe the nature of $\chi_{c1}(3872)$ [21, 23, 24]. The radiative-to-di-$J/\psi$ yield ratio provides a spin-sensitive constraint, while compatible pole positions extracted from both amplitudes would test a common resonance interpretation against a threshold enhancement [19].

Radiative decays also provide search channels for states below their lowest double-charmonium strong-decay thresholds. Beyond the observed $2^{++}$ states, the $0^{++}$ and $1^{+-}$ sectors remain unestablished. A recent lattice calculation [25] puts the $0^+$ and $1^+$ levels 33.2(13) and 10.9(14) MeV below the $\eta_c\eta_c$ and $\eta_cJ/\psi$ thresholds, respectively. If these levels lie below their corresponding thresholds, the strong decays are kinematically closed. The $0^+$ state could then decay only electromagnetically or weakly, with $\gamma J/\psi$ expected to dominate over the diphoton channel by orders of magnitude [26, 27]. The corresponding radiative search channel for the $1^{+-}$ state is $\eta_c$; charge-conjugation parity forbids interchanging the two charmonia in these channels.

The $\gamma J/\psi$ and $\gamma\eta_c$ widths have so far been calculated only at LO [28]. A complete one-loop calculation tests NRQCD factorization and calibrates the relation between radiative widths and tetraquark LDMEs. This is important at the charm scale, where NLO corrections to quarkonium annihilation and fully heavy tetraquark electromagnetic decays can reach tens of percent or order unity [27, 29, 30, 31]. The $0^{++}$ channel involves two $S$-wave color-spin operators and their interference, so its NLO hard coefficient is a $2\times2$ matrix, whereas the $2^{++}$ and $1^{+-}$ channels each have a single leading-$v$ operator. A common NLO normalization would preserve the LO color-spin hierarchy, whereas a nonuniform correction would change the relative weights and make the radiative width sensitive to the initial color-spin composition.

In this Letter, we present the first complete NLO QCD calculation of $T_{4c}(0^{++},2^{++})\to\gamma J/\psi$ and $T_{4c}(1^{+-})\to\gamma\eta_c$ at leading order in the relative velocity. We determine the short-distance coefficients for all relevant color-spin operators, including scalar interference, and predict the tensor $J/\psi$ polarization and dilepton angular distribution without an overall LDME normalization. For a spin-resolved tensor state, we show how the ratio of its $\gamma J/\psi$ and di-$J/\psi$ yields, combined with an inclusive-production fit and the total width, can constrain its LDME without a potential-model input. This extraction requires compatible production and decay normalizations. If scalar and tensor signals are unresolved, the same yield ratio can instead constrain their relative contributions, provided their respective radiative-to-hadronic decay ratios are known. We discuss these inputs and assumptions below.

NRQCD factorization and matching

The decays $T_{4c}(0^{++},2^{++}) \to\gamma J/\psi$ and $T_{4c}(1^{+-}) \to\gamma\eta_c$ are treated in NRQCD at leading order in the heavy-quark relative velocity $v$ [22, 28]. The $J/\psi$ and $\eta_c$ are projected onto the color-singlet ${}^{3}S_{1}$ and ${}^{1}S_{0}$ states, respectively. Their decay widths can be written as

Here $\rho_H = \langle O_1^{J/\psi}({}^{3}S_1)\rangle$ for $J=0,2$ and $\rho_H = \langle O_1^{\eta_c}({}^{1}S_0)\rangle$ for $J=1$. The matrices $h_J^{(i)}$ and $\rho_{T,J}$ are the hard coefficients and tetraquark LDMEs. The NLO matrices $h_J^{(1)}$ follow from Born–one-loop interference in each color-spin channel.

For identical charm quarks, Pauli statistics restrict the allowed $S$-wave color-spin bases to

For each operator $Q_n$, we define $L_n = \langle0|Q_n^{\mathrm{ren}}|T\rangle$ and $(\rho_T)_{nm} = L_n L_m^*$. The scalar hard coefficient and LDME are therefore $2 \times2$ matrices. Their off-diagonal entries describe interference between the two color-spin amplitudes. The transpose in Eq. (1) follows from $\left|\sum_n F_nL_n\right|^2$ and is retained for complex LDMEs. The axial and tensor channels each contain one leading-$v$ operator. The tensor amplitude nevertheless has three physical polarization components, whose spin-averaged contraction gives a single hard coefficient. Reduced nonrelativistic state normalization is used throughout. State-conversion factors and initial-spin averages are included once in Eq. (1).

The hard coefficients are determined by matching renormalized QCD amplitudes onto NRQCD with on-shell partonic states. The matching is performed at $P_T=P_H+k$, with $P_T^2=(4m_c)^2$, $P_H^2=(2m_c)^2$, and $k^2=0$. The potential-region contribution is isolated at finite relative velocity and subtracted before the zero-relative-momentum limit is taken [32]. Dimensional regularization with $D=4-2\epsilon$ is used in both theories. The charm-quark mass and field are renormalized in the on-shell scheme, while $\alpha_s$ is renormalized in the $\overline{\mathrm{MS}}$ scheme with $n_f=4$. The nonpotential NRQCD loop contribution vanishes in pure dimensional regularization. Since

$T_{4c} \to\gamma H_g$ vanishes by color conservation, the NLO correction is purely virtual. The BMHV prescription is adopted for $\gamma_5$ in the axial channel [33, 34].

We calculate the on-shell partonic amplitude $c c\bar{c}\bar{c} \to c\bar{c}+\gamma$ and project its external heavy-quark states onto the required color and spin configurations. Representative Born and one-loop diagrams are shown in Figure 1. At Born level, our implementation also reproduces the short-distance coefficients for $T_{4c} \to\gamma\psi$ using the same operator conventions [26]. The one-loop integrals reduce to 52 common master-integral classes. Their analytic finite combinations agree with independent numerical integral evaluations. The hard coefficients are finite after renormalization and matching. Photon Ward identities are verified before the transverse projection, and the $D$-dimensional spin factors are retained through the Laurent expansion. The logarithmic terms agree with the running of $\alpha_s^4$. The Supplemental Material provides operator conventions, finite amplitudes, matching checks, and numerical prescriptions [35].

NLO short-distance coefficients and radiative widths

For $N_c=3$ and $n_f=4$, the hard coefficients are evaluated at the leading-order matching point $M_T=4m_c$ and $M_H=2m_c$. The finite one-loop terms below correspond to $\mu_R=m_c$. In the ordered scalar basis $\{\bar{3}\otimes3,6\otimes\bar{6}\}$, the hard matrices are