# Quark mass dependence of the low-lying charmed mesons at one loop in
HH$\chi$PT

F. Gil-Domínguez<EMAIL_ADDRESS>Departamento de Física Teórica and
IFIC, Centro Mixto Universidad de Valencia-CSIC, Parc Científic UV, C/
Catedrático José Beltrán, 2, 46980 Paterna, Spain  R. Molina
<EMAIL_ADDRESS>Departamento de Física Teórica and IFIC, Centro Mixto
Universidad de Valencia-CSIC, Parc Científic UV, C/ Catedrático José Beltrán,
2, 46980 Paterna, Spain

###### Abstract

We study the light and heavy quark mass dependence of the low-lying charmed
mesons in the framework of one-loop HH$\chi$PT. The low energy constants are
determined by analyzing the available lattice data from different LQCD
simulations. Model selection tools are implemented to determine the relevant
parameters as required by data with a higher precision. Discretization and
other effects due to the charm quark mass setting are discussed.

##  I Introduction

The discovery of exotic hadrons in the heavy quark sector, which cannot be
accomodated in terms of $q\bar{q}$ mesons or $qqq$ baryons, like those with a
tetraquark or pentaquark structure, has manifested the relevance of hadronic
loops in order to explain the masses and other properties of many states in
the hadron spectrum Guo et al. (2018a).

Even there has been claims since long time of such effects for some low-lying
states. Such is the case, for example, of the scalar mesons. An approach for
(pseudoscalar) meson - (pseudoscalar) meson scattering based on chiral
symmetry and unitarity, the so-called Unitarized Chiral Perturbation Theory
(UChPT) Truong (1988); Dobado and Pelaez (1993, 1997); Nieves and Ruiz Arriola
(2000); Garcia-Recio et al. (2004); Oller et al. (1999); Gomez Nicola and
Pelaez (2002) has predicted the properties of the $\sigma$, $a_{0}$ and
$f_{0}$, being in reasonably good agreement with experiment Oller and Oset
(1997); Oller et al. (1998); Kaiser (1998); Locher et al. (1998); Nieves and
Ruiz Arriola (2000). The study of the large $N_{c}$ behaviour concludes that,
while the low-lying scalars are non-ordinary mesons, the vector mesons are
dominated by the $q\bar{q}$ interaction Pelaez (2016). From the Lattice QCD
(LQCD) side, however, both, $q\bar{q}$ and $\pi\pi$ operators are needed in
order to obtain the properties of the $\rho$ meson Aoki et al. (2007);
Gockeler et al. (2008); Feng et al. (2011); Lang et al. (2011); Pelissier and
Alexandru (2013); Bali et al. (2016); Guo et al. (2016); Erben et al. (2020);
Wilson et al. (2015); Dudek et al. (2013); Bulava et al. (2016); Feng et al.
(2015); Alexandrou et al. (2017); Fu and Wang (2016); Metivet (2015); Fischer
et al. (2020). The role of the coupling of the $K\bar{K}$ channel to $\pi\pi$
in the $\rho$-meson mass has been studied in Guo et al. (2016); Hu et al.
(2016, 2017); Molina and Ruiz de Elvira (2020).

In the charm sector, it could be easier to dissentangle ordinary from non-
ordinary mesons, since many of the observed particles are narrow states. In
fact, a plethora of exotic states, for which hadron loops might play a
relevant role, especially when these states are close to thresholds, have been
discovered in the recent years Olsen et al. (2018); Chen et al. (2022).
Unexpectedly, it has been found that, when hadron loops are incorporated in
the quark model, mass shifts can be quite large even for the low-lying
charmonium states Barnes and Swanson (2008).

In the heavy-light meson sector, heavy quark symmetry arises in the heavy-
quark limit ($m_{Q}\to\infty$). This involves both Heavy Quark Spin Symmetry
(HQSS) and Heavy Quark Flavor Symmetry (HQFS) Isgur and Wise (1989, 1990,
1991); Manohar and Wise (2000). Heavy quark spin symmetry emerges in this
limit because the color magnetic moment of a heavy quark is proportional to
the inverse of its mass, $\mu^{c}\propto 1/m_{Q}$, and thus the interaction
becomes invariant under the SU(2) spin symmetry, i. e., for each $Q\bar{q}$
state there will be another one degenerated in mass which can be obtained by a
flip of the spin of the heavy quark. The spin-parity of the heavy and light
degrees of freedom are then separately conserved. Beyond that, the exact value
of the heavy quark mass in the interactions with the light degrees of freedom
plays no role and this gives rise to heavy flavor symmetry. Obviously, these
symmetries are broken for finite quark masses and terms that go like $1/m_{Q}$
could be relevant for the charm quark.

Thus, in the heavy quark limit, one can use the light degrees of freedom
(d.o.f) to classify the heavy-light meson states. The total angular momentum
of the light d.o.f are $j_{l}=L\pm\frac{1}{2}$, being $L$ the angular
momentum. In the ground state, two degenerate states with
$J^{P}=\\{0,1\\}^{-}$ appear for $L=0$ ($j_{l}^{P}=\frac{1}{2}^{-}$), while
for $L=1$ there are a couple of doublets, one with $J^{P}=\\{0,1\\}^{+}$
($j_{l}^{P}=\frac{1}{2}^{+}$), and another for $J^{P}=\\{1,2\\}^{+}$
($j_{l}^{P}=3/2^{+}$). Actually, one can find in the PDG the $D$ and $D^{*}$
mesons which correspond to the ground state doublet, while the
$D^{*}_{0}(2300)$ and $D_{1}(2420)$ could be related to the
$J^{P}=\\{0,1\\}^{+}$ doublet, and the $D_{1}(2430)$ and $D^{*}_{2}(2460)$
would correspond to the doublet with $J^{P}=\\{1,2\\}^{+}$.

Similarly, we find in the PDG their strange counterparts, the
$(D_{s},D^{*}_{s})$, $(D_{s0}(2317)$, $D_{s1}(2460))$, and $(D_{s1}(2536)$,
$D_{s2}(2573))$, which should correspond to the doublets with
$J^{P}=\\{0,1\\}^{-}$, $\\{0,1\\}^{+}$ and $\\{1,2\\}^{+}$ respectively.
Regarding these states, the constituent quark model predicts broad states
decaying to $D^{(*)}K$ for the $J^{P}=\\{0,1\\}^{+}$ doublet with masses
higher in about $100$ MeV than the experimental ones Godfrey and Isgur (1985);
Godfrey and Kokoski (1991); Di Pierro and Eichten (2001). In contrast, the
observed states, the $D_{s0}(2317)$ and $D_{s1}(2460)$, are very narrow
resonances which are placed relatively close to the $DK$ and $DK^{*}$
thresholds respectively in the spectrum Aubert et al. (2003); Besson et al.
(2003). For these reasons, molecular Barnes et al. (2003); Kolomeitsev and
Lutz (2004); Faessler et al. (2007a, b); Lutz and Soyeur (2008); Albaladejo et
al. (2018) and also tetraquark explanations Cheng and Hou (2003); Terasaki
(2003); Maiani et al. (2005); Dmitrasinovic (2005) have been proposed. On the
other hand, the predictions of the constituent quark model for the
$J^{P}=\\{1,2\\}^{+}$ ($j_{l}^{P}=3/2^{+}$) doublet are in agreement with
experiment Godfrey and Isgur (1985); Godfrey and Kokoski (1991); Di Pierro and
Eichten (2001). Still, molecular components for the same quantum numbers
cannot be discarded. In fact, in Molina et al. (2010), assuming that the
$D^{*}_{2}(2573)$ is a $D^{*}K^{*}$ molecular state, the recenlty observed
$T_{cs}$ and $T_{c\bar{s}}$ states Aaij et al. (2020a, b); Collaboration
(2022) are predicted, see also Molina and Oset (2020, 2022).

In general, states with an additional strange quark are in principle expected
to be $100$ MeV higher in mass, since $m_{s}/m_{d}\simeq 20$ MeV with
$m_{d}\approx 5$ MeV, however, this is not the case for the charmed mesons. As
reported by the PDG, the masses of the $D^{*}_{0}(2300)$ and $D_{1}(2420)$ are
similar to those of their strange partners, $D_{s0}^{*}(2317)$ and
$D_{s1}(2460)$. The situation for the $L=1$ charmed mesons is more complex.
The measurements of the masses of the $0^{+}$ state lie between $2300$ and
$2400$ MeV Abe et al. (2004); Aubert et al. (2009); Link et al. (2004).
Theoretically, a two pole structure has been proposed in coupled-channel
approaches based on chiral lagrangians, unitarity and heavy quark symmetry,
resulting in two states with masses of around $2100$ MeV and $2400-2450$ MeV
Kolomeitsev and Lutz (2004); Guo et al. (2006, 2009); Albaladejo et al.
(2017); Du et al. (2018). The existence of a lower mass pole is supported by
lattice QCD, which consistently finds the mass of the $D^{*}_{0}(2300)$ well
below that of the $D^{*}_{s0}(2317)$, in contrast to the experiment Moir et
al. (2016); Gayer et al. (2021).

A strongly interacting system at low energies composed of heavy and light
mesons with the heavy meson being of $Q\bar{q}$ type can be studied by means
of Heavy Hadron Chiral Perturbation Theory (HH$\chi$PT) Wise (1992); Burdman
and Donoghue (1992); Yan et al. (1992). This is an effective field theory
based on both chiral and heavy quark symmetries, where a systematic two double
expansions in the parameters $Q/\Lambda_{\chi}$ and $\Lambda_{QCD}/m_{Q}$ can
be made, being $Q$ the momentum transfer (soft scale), $Q\sim m_{\pi}\sim
p_{\pi}$, $m_{Q}$, is the mass of the heavy quark, and $\Lambda_{\chi}=4\pi
f\simeq 1$ GeV (hard scale). HH$\chi$PT has been applied to evaluate the
ground state heavy meson masses in the heavy quark limit Isgur and Wise (1989,
1990). Leading order corrections to these masses due to the finite heavy quark
mass have been studied in Rosner and Wise (1993); Randall and Sather (1993);
Di Bartolomeo et al. (1995); Blok et al. (1997); Jenkins (1994); Yeh and Lee
(1995) for non-zero light quark masses. The $j^{P}_{l}=\frac{1}{2}^{+}$
doublet has also been introduced as an explicit degree of freedom in Mehen and
Springer (2005); Ananthanarayan et al. (2007); Alhakami (2016) where the
masses of both doublets were studied including one-loop corrections.

Given the large amount of parameters in one-loop HH$\chi$PT it has not been
possible yet to pin down precisely the parameters of the theory. Because of
the recent progress in lattice QCD it is worth coming back to this issue. In
this article we study the quark mass dependence of the ground state charmed
mesons, $D$, $D^{*}$, $D_{s}$ and $D^{*}_{s}$ within one-loop HH$\chi$PT
Jenkins (1994) by analyzing the available lattice data. There are several
reasons to do this: (1) to test the predictive power of HH$\chi$PT with the
available lattice data, (2) to better understand the significance of one-loop
corrections involving light mesons for the ground state charmed mesons, (3) to
provide an overview of the trends of charmed mesons in different lattice QCD
simulations, and (4), since many exotic states have been observed in the
recent years Guo et al. (2018a); Aubert et al. (2003); Besson et al. (2003);
Aaij et al. (2020a, b, 2022), the analysis done here can be used in the future
to investigate the quark mass dependence of these new states. This could help
to discriminate between diferent pictures (tetraquark, molecule …).

We investigate both, the light and heavy quark (near the physical point) mass
dependence of the ground state charmed mesons, $D_{(s)}$ and $D_{(s)}^{*}$. We
investigate which parameters of one-loop HH$\chi$PT are more relevant as
required by the LQCD data implementing model selection tools, as the Least
Absolute Shrinkage and Selection Operator (LASSO) method in combination with
cross validation and information criteria James and Tibshirani (2009);
Tibshirani (1996); Landay et al. (2017). This can be possible now since there
are sufficient lattice simulations to better constraint the parameters of
HH$\chi$PT. For this analysis, we have taken as input the lattice data on the
ground state charmed mesons of Ref. Kalinowski and Wagner (2015); Carrasco et
al. (2014) (ETMC), the PACS-CS ensembles Mohler and Woloshyn (2011); Aoki et
al. (2009), the ensembles from Hadron Spectrum Collaboration (HSC) Cheung et
al. (2016, 2021), the data of Table 1 of Prelovsek et al. (2021), the data
from Ref. Bruno et al. (2017) (CLS), and that of Refs. Bali et al. (2017,
2015) (RQCD). Most of these data are collected in Tables tables VIII-XIV of
Guo et al. (2018b), where a previous lattice data analysis using the chiral
lagrangians developed in Lutz and Soyeur (2008); Kolomeitsev and Lutz (2004)
was performed. In addition, we incorporate the latest studies of charmed meson
masses from MILC in $N_{f}=2+1+1$ Bazavov et al. (2018a, b), which have been
done for different quark masses and lattice spacings.

##  II Formalism

### Heavy meson masses to one-loop order in HH$\chi$PT

In HH$\chi$PT Wise (1992); Burdman and Donoghue (1992); Yan et al. (1992) it
is customary to introduce a single field, $H_{l}$, with $l=1,2,3$, to refer to
$u,d,s$ respectively, including both the pseudoscalar and vector meson fields
defined as, $P_{l}^{(*)}=(D^{(*)0},D^{(*)+},D_{s}^{(*)})$, for the charmed
mesons (and similarly for the bottom mesons),

$H^{(Q)}_{l}=\frac{(1+\not{v})}{2}(P^{*(Q)}_{l\mu}\gamma^{\mu}-P^{(Q)}_{l}\gamma_{5})\
.$ (1)

with $Q=c,b$. Similarly, the conjugated meson field is defined as,

$\bar{H}^{(Q)}_{l}=\gamma^{0}H_{l}^{(Q)\dagger}\gamma^{0}=(P^{*(Q)\dagger}_{l\mu}\gamma^{\mu}+P^{(Q)\dagger}_{l}\gamma_{5})\frac{(1+\not{v})}{2}\
.$ (2)

These heavy meson fields transforms as a doublet under heavy quark spin
symmetry and like an anti-triplet under $SU(3)_{V}$ flavour symmetry, i.e.,

$H_{l}^{(Q)}\to S_{Q}(H^{(Q)}U^{\dagger})_{l}$ (3)

and

$\bar{H}^{(Q)}_{l}\to(U\bar{H}^{(Q)})_{l}S^{\dagger}_{Q}\ ,$ (4)

for the heavy meson field and its conjugate, respectively. A velocity
dependent field can also be defined as Georgi (1990),

$H^{(Q)}_{v}=\sqrt{m_{H}}e^{i\,m_{H}}\not{v}v_{\mu}x^{\mu}H^{(Q)}(x)\ ,$ (5)

which satisfies $\not{v}H_{v}=H_{v}$ and $H_{v}\not{v}=-H_{v}$. The mass
$m_{H}$ in Eq. (5) can be expanded in terms of powers of $1/m_{Q}$ starting in
$m_{H}=m_{Q}$. The chiral lagrangian involving heavy meson fields defined as
in Eq. (5) and pseudoscalar fields, has a systematic expansion in $1/m_{H}$,
derivative and loop expansion. At leading order, the chiral lagrangian reads
Jenkins (1994),

$\displaystyle{\cal L}^{0}=$
$\displaystyle-i\,\mathrm{Tr}\bar{H}_{v}(v\cdot\partial)H_{v}+i\,\mathrm{Tr}\bar{H}_{v}H_{v}(v\cdot
V)$ $\displaystyle+2g\,\mathrm{Tr}\bar{H}_{v}H_{v}(S_{lv}\cdot A)$ (6)

where $S_{lv}$ stands for the velocity dependent spin operator $S^{\mu}_{lv}$
Jenkins and Manohar (1991), which satisfies $v\cdot S_{lv}=0$, and in the rest
frame of the heavy meson, $v^{\mu}=(1,0)$, this is,
$S^{\mu}_{l}=(0,\vec{\sigma}/2)$. The coupling $g$ can be fixed to obtain the
experimental radiative and pion decay widths of the charmed vector mesons, we
set $g^{2}=0.55$. The pseudoscalar meson fields appear in the combinations,

$V^{\mu}=\frac{1}{2}(\xi\partial^{\mu}\xi^{\dagger}+\xi^{\dagger}\partial^{\mu}\xi)\
,\quad
A^{\mu}=\frac{i}{2}(\xi\partial^{\mu}\xi^{\dagger}-\xi^{\dagger}\partial^{\mu}\xi)\
,$ (7)

where $\xi=\sqrt{U}$, $U=e^{i\phi/\sqrt{2}f}$, and

$\displaystyle\phi=$
$\displaystyle\left({\begin{array}[]{ccc}\tfrac{\pi^{0}}{\sqrt{2}}+\tfrac{1}{\sqrt{6}}\eta&\pi^{+}&K^{+}\\\
\pi^{-}&-\tfrac{\pi^{0}}{\sqrt{2}}+\tfrac{1}{\sqrt{6}}\eta&K^{0}\\\
K^{-}&\bar{K}^{0}&-\tfrac{2}{\sqrt{6}}\eta\\\ \end{array}}\right).$ (11)

For the one-loop calculation of the heavy meson masses, terms with up to two
insertions of the light quark matrix $M=\mathrm{diag}(m_{u},m_{d},m_{s})$, are
needed because counterterms at one-loop order are proportional to two powers
of the light quark masses. Thus, the relevant lagrangian needed for the one-
loop calculation is given by Jenkins (1994),

$\displaystyle{\cal L}^{\not{s}}_{v}=$
$\displaystyle\,\sigma\,\mathrm{Tr}M_{\xi}\mathrm{Tr}\bar{H}_{v}H_{v}+a\,\mathrm{Tr}\bar{H}_{v}H_{v}M_{\xi}+b\,\mathrm{Tr}\bar{H}_{v}H_{v}M_{\xi}M_{\xi}+c\,\mathrm{Tr}M_{\xi}\mathrm{Tr}\bar{H}_{v}H_{v}M_{\xi}+d\,\mathrm{Tr}M_{\xi}M_{\xi}\mathrm{Tr}\bar{H}_{v}H_{v}$
$\displaystyle-\frac{\Delta}{8}\mathrm{Tr}\bar{H}_{v}\sigma^{\mu\nu}H_{v}\sigma_{\mu\nu}-\frac{\Delta^{(a)}}{8}\mathrm{Tr}\bar{H}_{v}\sigma^{\mu\nu}H_{v}\sigma_{\mu\nu}M_{\xi}-\frac{\Delta^{(\sigma)}}{8}\mathrm{Tr}M_{\xi}\mathrm{Tr}\bar{H}_{v}\sigma^{\mu\nu}H_{v}\sigma_{\mu\nu}\
,$ (12)

where $M_{\xi}=\frac{1}{2}(\xi M\xi+\xi^{\dagger}M\xi^{\dagger})$ and $M$ is
the (diagonal) light quark mass matrix. Here it is assumed $m=m_{u}=m_{d}$. If
considering isosping breaking, extra terms should be added, see Eq. (2.14) of
Jenkins (1994). In general, operators which respect HQSS have an expansion
which starts at ${\cal O}(1)$, while the ones violating HQSS, start at ${\cal
O}(\frac{1}{m_{Q}})$. The coefficients $a$ and $\sigma$ are dimensionless.
These are functions of $\frac{1}{m_{Q}}$ starting at ${\cal O}(1)$ and
therefore terms with these coefficients respect HQSS. These terms with $a$ and
$\sigma$ give rise to the $SU(3)_{V}$ violating mass splittings among the
$P^{(*)}$ mesons and to a singlet contribution to the masses, respectively.
Terms accompanying the coefficients $b,c$ and $d$, correspond to chiral
lagrangian terms with two insertions of the light quark mass matrix,
preserving also heavy quark spin symmetry. The hiperfine splittings of the
heavy mesons depend on the coefficients
$(\Delta,\Delta^{(\sigma)},\Delta^{(a)})$, which start at ${\cal
O}(\frac{1}{m_{Q}})$, and therefore, terms accompanying these coefficients
violate HQSS. If one considers only terms with no insertions of the light
quark matrix, one has $M_{P^{*}_{a}}-M_{P_{a}}=\Delta$, which is the hyperfine
mass splitting at tree level. Terms with $\Delta^{(\sigma)}$ and
$\Delta^{(a)}$ in Eq. (12) correspond to one insertion of the light quark mass
matrix. While the former yields a light quark flavor dependent term, the
latter leads to a flavor singlet contribution. Both of them contribute to
counterterm coefficients in the $M_{P^{*}}-M_{P}$ splitting.

The one-loop masses of heavy mesons can be written in terms of two linear
combinations, the spin average term, $\frac{1}{4}(M_{P_{a}}+3M_{P^{*}_{a}})$
and the hyperfine splitting, $M_{P^{*}}-M_{P}$, which respects and violate
heavy quark spin symmetry, respectively. These, which have been calculated to
one-loop in Ref. Jenkins (1994), are

$\displaystyle\frac{1}{4}(M_{P_{l}}+3M_{P^{*}_{l}})$ $\displaystyle=$
$\displaystyle\,m_{H}+\alpha_{l}-\sum_{X=\pi,K,\eta}\beta_{l}^{(X)}\frac{M_{X}^{3}}{16\pi
f^{2}}+\sum_{X=\pi,K,\eta}\left(\gamma_{l}^{(X)}-\lambda_{l}^{(X)}\alpha_{l}\right)\frac{M_{X}^{2}}{16\pi^{2}f^{2}}\log\left(M_{X}^{2}/\mu^{2}\right)+c_{l}$
(13) $\displaystyle M_{P^{*}_{l}}-M_{P_{l}}$ $\displaystyle=$
$\displaystyle\,\Delta+\sum_{X=\pi,K,\eta}\left(\gamma_{l}^{(X)}-\lambda_{l}^{(X)}\Delta\right)\frac{M_{X}^{2}}{16\pi^{2}f^{2}}\log\left(M_{X}^{2}/\mu^{2}\right)+\delta
c_{l}\ $ (14)

where the coefficients $\alpha_{l}$, $\beta_{l}^{(X)}$,
$(\gamma_{l}^{(X)}-\lambda_{l}^{(X)}\alpha_{l})$, $c_{l}$,
$(\gamma_{l}^{(X)}-\lambda_{l}^{(X)}\Delta)$ and $\delta c_{l}$ are given in
the appendix A. The coefficients $\alpha_{l}$,
$(\gamma_{l}^{(X)}-\lambda_{l}^{(X)}\alpha_{l})$, $c_{l}$, and $\delta c_{l}$
are proportional to powers of the light quark masses while $\beta_{l}^{(X)}$
and $(\gamma_{l}^{(X)}-\lambda_{l}^{(X)}\Delta)$ accompany terms proportional
to $M_{X}$, the mass of the $X$ pseudoscalar meson, $\pi$, $K$ or $\eta$.
Thus, $m_{H}$ and $\Delta$ can be interpreted as the spin-average mass and
hyperfine splitting of the heavy mesons in the SU(3) chiral limit. The scale
$\mu$ in Eqs. (13) and (14) is fixed to $770~{}$MeV.

##  III Fitting procedure with the LASSO method

We collect the available lattice data on low-lying charmed mesons included
here in Tables 9-15, of the Appendix C. These include,

  * •

Data from Ref. Kalinowski and Wagner (2015), which uses Wilson twisted mass
lattice QCD with 2+1+1 dynamical quark flavors for pion masses in the range of
$225-445$ MeV, $L\sim 2-3$ fm and three lattice spacings $a=0.0619,0.0815$ and
$0.0885$ fm Carrasco et al. (2014) [ETMC] (Tables 9 and 10).

  * •

The PACS ensembles with $N_{f}=2+1$ flavor Clover-Wilson configurations Mohler
and Woloshyn (2011) for pion masses in the range $m_{\pi}\sim 150-410$ MeV,
$L\simeq 2.9$ fm and a lattice spacing $a=0.0907$ MeV Aoki et al. (2009)
[PACS-CS] (Table 11).

  * •

The ensembles from Hadron Spectrum Collaboration for $N_{f}=2+1$ flavours of
dynamical quarks with an anisotropic (clover) action Cheung et al. (2016),
with $a_{s}\sim 0.12$ fm, $L\sim 1.9-3.8$ fm and $m_{\pi}\simeq 240,390$ MeV.
See also Table 2 of Cheung et al. (2021) [HSC] (Table 12).

  * •

The data of Table 1 of Prelovsek et al. (2021) for a $N_{f}=2+1$ simulation
with Wilson dynamical fermions provided by the Coordinated Lattice Simulations
(CLS) with a pion mass of $m_{\pi}=280$ MeV, $L\simeq 2.1$ fm and a lattice
spacing of $a=0.08636$ fm Bruno et al. (2017) [CLS] (Table 13).

  * •

Data collected in Lang et al. (2014) for $N_{F}=2$ dynamical light quarks with
improved Wilson fermions, a lattice spacing $a=0.1239$, $L\simeq 2$ fm and
$m_{\pi}\simeq 266$ MeV [Prelovsek et al.] (Table 14).

  * •

The data of the simulation with $N_{f}=2$ improved clover Wilson sea quarks
collected in Table 1 of Bali et al. (2017) for $m_{\pi}=150,290$ MeV, $L\sim
1.7-4.5$ fm and $a=0.07$ fm Bali et al. (2015) [RQCD] (Table 15).

  * •

The latest studies of charmed meson masses from MILC with highly improved
staggered-quark ensembles with four flavors of dynamical quarks Bazavov et al.
(2018a, b), which have been done for different quark masses, $m_{\pi}\sim
130-330$, volumes $L\sim 2.5-6$ fm and lattice spacings ranging from
$0.03-0.15$ fm (Table I of Bazavov et al. (2018b)) allowing for a continuous
extrapolation Brambilla et al. (2018) [MILC]111We use as input the function in
the continuum limit from MILC Brambilla et al. (2018) for the $D$ and $D_{s}$
masses..

Some of these data were collected in Tables VIII-XIV of Guo et al. (2018b).
Indeed, the data of PACS-CS and ETMC used here correspond to Tables VIII, X
and XI of Guo et al. (2018b), where a previous lattice data analysis using the
chiral lagrangians developed in Lutz and Soyeur (2008); Kolomeitsev and Lutz
(2004) was performed. The data of Tables XII and XIII of Guo et al. (2018b)
based on MILC ensembles of the HPQCD and LHPC collaborations are not included
here, because instead we include the latest MILC studies of Bazavov et al.
(2018a, b) which superseed the older data in finests lattice spacing and
bigger volumes that yielded a continuous extrapolation Brambilla et al.
(2018).

In this work we will analyze the above data with one-loop order HH$\chi$PT
given by Eqs. (13) and (14), looking for consistency between the different
sets of data and trying to find a unique solution of the charmed meson masses
that provides both, the heavy and light quark mass dependence of those. At
present, this can only be achieved considering data from different
collaborations. Regarding this, several aspects need to be discussed. First,
we observe that, how the charm quark is taken to the physical point and which
value is taken the observable to do so, can explain largely the deviations
observed by the different collaborations, concluding that these actually come
mostly from different charm quark masses at the so-called “physical point”. To
compare the charm quark mass settings of the LQCD data used here, in Table 1
we compare the value of $m_{D_{s}}$ in physical units for the lightest pion
mass and finest lattice spacing of every collaboration.

Col. | $a$ | $L$ | $m_{\pi}$ | $m_{D_{s}}$
---|---|---|---|---
ETMC1 | 0.0619 | 3.0 | 224 | 1988
ETMC2 | 0.0619 | 3.0 | 224 | 2001
PACS-CS | 0.0907 | 2.9 | 156 | 1809
HSC | 0.12 | 3.8 | 239 | 1967
CLS | 0.08636 | 2.1 | 280 | 1981
Prelovsek et al. | 0.1239 | 2.0 | 266 | 1657
RQCD | 0.071 | 4.5 | 150 | 1977
Table 1: Comparison of the charm quark mass settings of the different
collaborations through their value of $m_{D_{s}}$. Note that the physical
value is $m_{D_{s}}=1968.35\pm 0.07$ MeV Zyla et al. (2020). $L$ and $a$ are
given in fm, and the masses in units of MeV.

For the value of the lattice spacing given in Cheung et al. (2021) determinded
through the $\Omega$ baryon mass, the mass of the $D_{s}$ for HSC is physical
for its light pion mass, the one for RQCD and CLS is around $10$ MeV above the
physical one, $20-30$ MeV above for ETMC, while being much lower for the PACS-
CS, and [Prelovsek et al.] data sets in $\sim 160$ and $\sim 300$ MeV,
respectively. Second, there might be discretization effects of different size.
In general, we consider that the systematic effects due to discretization are
larger for the charm quark mass and that these can be absorbed in the
parameters $m_{H}$ and $\Delta$, as follows,

$\displaystyle m_{H}=m_{H}^{\prime}+r_{H}{\cal O}(a^{k})$ (15)
$\displaystyle\Delta=\Delta^{\prime}+r_{\Delta}{\cal O}(a^{k})\ ,$ (16)

where the term dependent on the lattice spacing $a$ will be different for
every collaboration depending on the action used, and $m_{H}^{\prime}$,
$\Delta^{\prime}$, are just the limit in the continuum of the tree level
parameters. Here, we focus on fitting the $m_{H}$ and $\Delta$ parameters,
introducing a new couple of parameters for every different lattice spacing,
collaboration, and charm quark mass determination. There are two
collaborations which perform simulations at different lattice spacings, MILC
and ETMC. For MILC we perform a refit of the continuous extrapolation of
Brambilla et al. (2018) that already takes into account the possible
systematic effects of the staggered fermions, while for ETMC we analyze the
data first and then study the discretization effects a posteriori making
${\cal O}(a^{k})\to a^{2}$ in Eqs. (15) and (16).

In the next section, first we perform individual fits for every collaboration.
Then, we make an attempt to perform a global fit of data combining the
different data sets. When doing so, we realize that the inclusion of the PACS-
CS data set increases tremendously the $\chi^{2}$. For this reason we do not
include these data in the global fit. In the case of MILC, we have a
continuous function which is possible to fit in combination with the rest of
data (we take $10$ data points for the $m_{D}$ mass of this function spread
equally in the $m_{\pi}$ range of the lattice data, Table I of Bazavov et al.
(2018b)). Since the discretization effects are embedded in the different
$m_{H}$ and $\Delta$ parameters, there is no reason of concern to mix these
various data sets, and indeed, as we will see, we obtain a reasonably good
value of the $\chi^{2}_{\mathrm{d.o.f}}$. Thus, for the global fit, we have,
$7+2n$ parameters, with $n=12$, equals to $31$ parameters, $40$ data points
for $m_{D_{(s)}}$ and $30$ data for $m_{D_{(s)}^{*}}$222In the case of RQCD
(Table 15) there are no data for $m_{D_{s}^{(*)}}$., and the two equations
Eqs. (13) and (14). Apart from the various $(m_{H},\Delta)$ parameters
introduced, the other parameters to be fitted in Eqs. (13) and (14) are,
($\frac{\sigma m_{\pi}}{B_{0}}$, $\frac{am_{\pi}}{B_{0}}$,
$\frac{bm^{3}_{\pi}}{B_{0}^{2}}$, $\frac{cm^{3}_{\pi}}{B_{0}^{2}}$,
$\frac{dm^{3}_{\pi}}{B_{0}^{2}}$, $\frac{\Delta^{(\sigma)}m_{\pi}}{B_{0}}$,
$\frac{\Delta^{(a)}m_{\pi}}{B_{0}}$, where $m_{\pi}$ is the physical pion
mass, and it is introduced to redefine the parameters in terms of adimensional
quantities.

Given that the size of the data sample is not very large compared to the
number of parameters, one could easily have problems with overfitting the
data. In order to avoid that, we follow the procedure of Ref. Landay et al.
(2017), where the Least Absolute Shrinkage and Selection Operator (LASSO)
method in combination with cross validation James and Tibshirani (2009);
Tibshirani (1996) and information criteria were applied for the analysis of
photoproduction data with a set of models with different number of parameters.
With this method, only those models which contain an optimal number of
parameters able to describe the data with enough accuracy are selected. To
implement the LASSO, a penalty term is added to the $\chi^{2}$

$\displaystyle\chi^{2}(\vec{p},\vec{q})=\sum_{i}^{N}\frac{\left(y_{i}-f_{i}(\vec{p},\vec{q})\right)^{2}}{\sigma_{y_{i}}^{2}}+P,$
$\displaystyle P=\frac{\lambda}{10}\sum_{i=1}^{7}|p_{i}|\ .$ (17)

In order to select the best “$\lambda$” parameter, we proceed as follows: in
every iteration we divide randomly the original data set into three subsets.
We keep two subsets as the training data set and one subset as the validation
data set. These sets will provide the $\chi^{2}$ of the training and
validation sets, $\chi^{2}_{T}$, $\chi^{2}_{V}$, respectively. As explained in
Refs. James and Tibshirani (2009); Tibshirani (1996); Landay et al. (2017),
while $\chi^{2}_{T}$ is a monotonously increasing function of $\lambda$,
$\chi^{2}_{V}$ has a sweet spot in which the data are not underfitted (large
$\lambda$) and not overfitted (small $\lambda$). The fluctuations given by the
random selection of the validation set provide the error of this minimum, and
also the one of the $\lambda$ value, which is chosen as the
$\lambda_{\mathrm{opt}}>\lambda_{\mathrm{min}}$, compatible within errors with
the minimum of $\chi^{2}_{V}(\lambda_{\mathrm{min}})$.

The result for the $\lambda_{\mathrm{opt}}$ from cross validation is compared
with the one obtained from three information criteria, which are, AIC, AICc,
and BIC, defined as,

$\displaystyle\text{AIC}=2k-2\log(L)=2k+\chi^{2}\ ,$
$\displaystyle\text{AICc}=\text{AIC}+\frac{2k(k+1)}{n-k-1}\ ,$
$\displaystyle\text{BIC}=k\log(n)-2\log(L)=k\log(n)+\chi^{2},$

In the above equations, $k$ is the number of parameters that changes with
$\lambda$, $n$ is the length of the data, and $L$ is the likelihood. The
optimal value of $\lambda$ is given by the respective minimum of the criteria.

##  IV Results

Fit | $\frac{\sigma m_{\pi}}{{B_{0}}^{2}}\cdot 10^{-4}$ | $\frac{am_{\pi}}{{B_{0}}^{2}}\cdot 10^{-4}$ | $\frac{bm^{3}_{\pi}}{{B_{0}}^{2}}\cdot 10^{-10}$ | $\frac{cm^{3}_{\pi}}{{B_{0}}^{2}}\cdot 10^{-10}$ | $\frac{dm^{3}_{\pi}}{{B_{0}}^{2}}\cdot 10^{-10}$ | $\frac{\Delta^{(\sigma)}m_{\pi}}{{B_{0}}}\cdot 10^{-5}$ | $\frac{\Delta^{(a)}m_{\pi}}{{B_{0}}}\cdot 10^{-5}$
---|---|---|---|---|---|---|---
ETMC${}_{\textrm{(I)}}$ | $0.89\pm 0.08$ | $2.08\pm 0.13$ | $-0.55\pm 0.21$ | $6.51\pm 0.20$ | $3.39\pm 0.18$ | $5.78\pm 3.43$ | $-2.16\pm 2.27$
ETMC${}_{\textrm{(II)}}$ | $1.02\pm 0.01$ | $2.02\pm 1.18$ | - | $6.16\pm 0.24$ | $3.24\pm 0.27$ | $5.88\pm 3.85$ | $-2.17\pm 2.24$
MILC${}_{\textrm{(I)}}$ | $3\pm 22$ | $3\pm 18$ | $7\pm 108$ | $2\pm 55$ | $1.56\pm 7.08$ | $4\pm 255$ | $36\pm 478$
MILC${}_{\textrm{(II)}}$ | $1.74\pm 0.40$ | $2.40\pm 0.73$ | - | $5.61\pm 1.66$ | - | - | -
PACS${}_{\textrm{(I)}}$ | $0.90\pm 0.24$ | $2.62\pm 0.20$ | $-2.38\pm 1.81$ | $7.03\pm 1.24$ | $3.53\pm 0.74$ | $3.60\pm 0.16$ | $-0.76\pm 0.44$
PACS${}_{\textrm{(II)}}$ | $1.22\pm 0.02$ | $2.36\pm 0.03$ | - | $5.40\pm 0.03$ | $4.50\pm 0.03$ | $3.59\pm 0.16$ | $-0.75\pm 0.44$
GLOBAL${}_{\textrm{(I)}}$ | $1.84\pm 0.01$ | $2.54\pm 0.02$ | $0.57\pm 0.05$ | $4.47\pm 0.04$ | $2.34\pm 0.03$ | $2.70\pm 0.21$ | $-1.43\pm 0.47$
GLOBAL${}_{\textrm{(II)}}$ | $1.17\pm 0.01$ | $2.58\pm 0.02$ | - | $4.89\pm 0.05$ | $2.27\pm 0.03$ | $2.56\pm 0.20$ | $-1.47\pm 0.46$
Table 2: Collection of parameters from different fits. The subindex (I) means
”before LASSO” and (II) means ”after LASSO”.

In this section, we present, first, the results obtained from indiviudal fits
for data sets larger enough to be fitted by means Eqs. (13), (14) (see also
Eqs. (29), (30), (31)). These data sets are: ETMC, PACS and MILC. Second, we
show the result of a global fit including also the rest of the data sets, HQS,
CLS, [Prelovsek et al.] Lang et al. (2014), with few data points for heavy
meson masses.

All the fits are carried out with quantities in units of the lattice spacing
($am_{\pi}$, $am_{D}$…). As commented in the previous section, we assume that
discretization effects are much larger for the charm quark mass and absorb the
overall dependence of the charmed meson masses with the lattice spacing in the
$m_{H}$ and $\Delta$ parameters (we introduce one pair of ($m_{H}$, $\Delta$)
for every lattice spacing). For the data sets with more than one lattice
spacing, ETMC data, we also analyze this discretization effect after the fit.
We will compare this result at $a=0$ with the one from MILC (extrapolated data
to the continuum) and the result of the global fit.

###  IV.1 ETMC

The results after performing the LASSO are collected in Fig. 1.

Figure 1: LASSO method plots for ETMC data analysis.

As can be seen in Fig. 1, in the first figure, from top to bottom, the
$\chi^{2}_{\mathrm{dof}}$ is a monotonously increasing function with
$\lambda$. Fig. 1 (second) shows the relevance of the parameters according to
the lattice data. In this figure it is clear that the data do not constrain
well the parameter $b$, which according to the analysis with LASSO, it is a
superfluous degree of freedom. Next, $d$ and $\Delta_{\sigma}$ are also less
relevant for the data. The $b$ and $d$ parameters correspond to chiral
lagrangian terms with two insertions of the light quark mass matrix preserving
HQSS, while $\Delta_{\sigma}$ violates HQSS and provides higher order
corrections to the hyperfine splitting at tree level. Notice that parameters
of these both types, like $c$ and $\Delta_{a}$, are, however, required by the
data, i. e., not all the parameters of the same order are equally relevant.

One should have a criteria in order to decide whether one or the three
parameters can be removed. We follow the criteria used in Landay et al.
(2017), explained in the previous section, cross validation and information
criteria. These are given in the third and fourth panels of Fig. 1,
respectively. Information criteria are evaluated with all data, while for
cross validation, we divide the data sets in three sets for every lattice
spacing, $a$, and select $2$ sets for the training set and $1$ for the
validation (test) set. As one can see, the information criteria have three
minima corresponding to leaving out, one parameter, $b$, which is also the
global minima, two parameters, $b$ and $d$, and three parameters, $b$, $d$ and
$\Delta_{\sigma}$, respectively. The global minima coincides with the one from
cross validation, in Fig. 1 (second). Thus, we select this parametrization
dropping out parameter $b$ as the most predictive one.

The parameters obtained before and after the LASSO, that we call Fits I and
II, respectively, are given in the two first rows of Tables 2 and 3.

Fit |  $m_{H}$ (MeV) |  $\Delta$ (MeV)
---|---|---
Col. |  $a$ (fm) | (I) | (II) | (I) | (II)
ETMC1 | $0.0619$ | $1995\pm 12$ | $1988\pm 16$ | $122\pm 21$ | $121\pm 24$
$0.0815$ | $2014\pm 13$ | $2007\pm 17$ | $129\pm 25$ | $128\pm 27$
$0.0885$ | $2015\pm 14$ | $2008\pm 18$ | $132\pm 28$ | $131\pm 29$
ETMC2 | $0.0619$ | $2008\pm 12$ | $2001\pm 16$ | $123\pm 21$ | $123\pm 23$
$0.0815$ | $2038\pm 13$ | $2031\pm 16$ | $126\pm 26$ | $126\pm 27$
$0.0885$ | $2047\pm 14$ | $2041\pm 18$ | $131\pm 25$ | $130\pm 27$
ETMC | 0 | $1972$ | $1965$ | $114$ | $114$
Table 3: Tree level parameters of the ETMC fits. The remaining parameters are
shown in the first two rows of Table 2.

We see that the result between the two different fits is very similar and
that, in fact, the inclusion of the $b$ parameters, do not affect much the
values of the rest of the parameters in this case.

For these data sets, one can also study the discretization effects. Assuming,

$\displaystyle m_{H}=m_{H}^{\prime}+r_{H}a^{2}$
$\displaystyle\Delta=\Delta^{\prime}+r_{\Delta}a^{2}\ ,$ (18)

with the value of $m_{H}$ and $\Delta$ for different $a$ obtained in the fit
and its errors plotted in Fig. 2.

Figure 2: Continuum extrapolation of $m_{H}$ and $\Delta$ parameters from the
ETMC analysis (Fit II).

We fit the two different data from fixing the charm quark mass with the
$\eta_{c}$ or $J/\psi$ mass with Eq. (18). We observe that, in these data
sets, the difference between fixing the charm quark mass with the $\eta_{c}$
or $J/\psi$ mass is basically due to discretization effects, since this
difference is close to zero for $a\sim 0$, see Fig. 2 (top). For the
extrapolation to the continuum of $m_{H}$, we take the average of the two
lines at $a\to 0$. For the hyperfine splitting, we obtain the line in Fig. 2
(bottom). Since this parameter starts its expansion at ${\cal O}(1/m_{Q})$,
there is no difference between setting the charm quark mass with $\eta_{c}$ or
$J/\psi$, and the error due to discretization effects of this parameter is
much smaller, being much more stable when varying the lattice spacing. For
$a=0$, we give the values of $m_{H}$ and $\Delta$ in the last row of Table 3.

We observe that the parameter $m_{H}$ is indeed very close to the physical
value of $(m_{D}+3m_{D^{*}})/4=1973$ MeV, however, $\Delta$ differs $\sim
20\%$ from the physical hyperfine splitting, $m_{D^{*}}-m_{D}=141~{}$MeV. This
motivates us to perform also fits with the expressions of the charmed meson
masses at tree level, keeping terms with one insertion of the light quark mass
matrix in Eq. (12). The results are given in the Appendix B. While for the
one-loop expressions we obtain $\Delta=114$ MeV for $a=0$, for the tree level
expressions, Eqs. (33) and (34), the continuous extrapolation of this
parameter that we obtain is $\Delta=129$ MeV. The values of $m_{H}$ are quite
similar for these different fits, $m_{H}\simeq 1970$. As expected, the one-
loop expressions reduce the value of the hyperfine splitting in the chiral
limit. We show the results of Fit II in Fig. 3. Eqs. (13) and (14) provide a
very good description of the data. The $\chi^{2}_{\mathrm{dof}}$ obtained is
0.44. To evaluate the errors, here and in the following, we use automatic
differentiation Ramos (2019) with the Julia ADerrors package.

Finally, we depict in Fig. 4 the extrapolation to the continuum limit of the
spin average mass and the hyperfine splitting (red band) in comparison with
the ETMC data set. As can be inferred, discretization effects appear to be of
$\sim 2-4\%$ in the spin average mass for these data sets. However, for the
hyperfine splitting, which suffers from large errors of $\sim 20-30\%$, most
of the LQCD data fall inside the error band. In the continuum limit, the LQCD
data agree well with the experimental value, shown with a star.

Figure 3: Results of the $D$, $D^{*}$, $D_{s}$ and $D^{*}_{s}$ meson masses
from the ETMC analysis (Fit II).

Figure 4: Comparison between ETMC data analysis and the continuum
extrapolation obtained.

###  IV.2 MILC

These lattice data are based on the simulations of Bazavov et al. (2018a, b)
for different pion masses and lattice spacings, which made it possible to
extrapolate to the continuum the masses of the charmed mesons using EFT’s,
with the theoretical framework being based on a merger of one-loop HQET and
HMrAS$\chi$PT, where a generic lattice spacing dependence, as well as higher
order terms, where incorporated Brambilla et al. (2018). This function has
been provided to us and here we conduct a refit of it with Eqs. (13) and (14).
In this case, the $\chi^{2}_{\mathrm{test}}$ is a monotonously increasing
function and does not show a minimum except for $\lambda=0$. However, the
information criteria, shown in Fig. 5 (top), exhibit a minimum around
$\lambda\simeq 2$. Actually, for this value of $\lambda$, four parameters have
been dropped out, $b$, $d$, $\Delta_{\sigma}$ and $\Delta_{a}$.

Figure 5: LASSO method plots for MILC data analysis. Collaboration |  $m_{H}$ (MeV) |  $\Delta$ (MeV)
---|---|---
(I) | (II) | (I) | (II)
MILC | $1863\pm 1271$ | $1959\pm 40$ | $110$ | $110$
Table 4: Tree level parameters of the MILC fit. The remaining parameters are
shown in the first two rows of Table 2.

Figure 6: Continuum extrapolation of charmed meson masses from ETMC (red),
MILC data (purple) and physical trajectory extracted from the global analysis
(blue).

Note that the first three parameters are also the less relevant parameters for
the ETMC data. The values of the parameters, before and after performing the
LASSO method, Fits I and II, are given in Tables 2 and 4. In this case the
data and continuous extrapolation only contain the values for the
$m_{D_{(s)}}$ masses. Since the masses of the vector counterparts are absent
we have to fix the value of $\Delta$. We have tried several values in the
range of the quantities obtained in the other fits, and here, fixing $\Delta$
to $110$ MeV seems to provide the best result. Since four parameters are
dropped after the LASSO, the values of the parameters in Fits I and II are
quite different. The result of the parameters after the LASSO method has much
smaller errors, telling that the method is working well. Concretely, in Fit
II, the value and error obtained for $m_{H}$ is quite reasonable and also
compatible with the one from the ETMC fit. The $\chi^{2}_{\mathrm{dof}}$
obtained in this case is small, 0.04, because of the refit of the function. In
Fig. 6 we depict the charmed meson masses in comparison with the continuous
extrapolation of the ETMC data that we did previously and the result of the
physical extrapolation for the global fit presented in Sec. IV.4. We observe
that the ETMC data show less dependence of the charmed meson masses with the
pion mass than MILC data and the result from the global fit, but the error
bands overlap each other and with the experimental charmed meson masses, shown
with stars. As could be expected, $m_{D_{s}}$ does not depend much with the
pion mass, since the strange quark mass is close to physical.

We can also compare our result for the charmed meson mass in the SU(3) chiral
limit with the one obtained in Bazavov et al. (2018a); Brambilla et al.
(2018). In Bazavov et al. (2018a) the value of $m_{D}$ in the SU(3) chiral
limit obtained is $1842.7\pm 5.4$ MeV, while our fit results in
$m_{D}^{SU(3)}=1877.0\pm 39.7$ MeV. These two results are consistent.

###  IV.3 PACS-CS

In this case the information criteria did not show a clear minimum, except for
$\lambda=0$. However, cross validation gives a clear minimum in
$\chi^{2}_{\mathrm{test}}$ around $\lambda\simeq 2$, as shown in Fig. 7 (top),
which corresponds to dropping out the $b$ parameter, see Fig. 7 (bottom).

Figure 7: LASSO method plots for PACS data analysis.

Nevertheless, the $\chi^{2}_{\mathrm{dof}}$ obtained is very high (much larger
than $1$). This is clearly not consistent with the previous data sets, where
very good values of the $\chi^{2}_{\mathrm{dof}}$ were obtained. We conclude
that this lattice data set for the charmed meson masses is in tension with the
rest of the lattice data sets. We notice that the high value of the $\chi^{2}$
for this collaboration is caused by one light pion mass. If we eliminate this
point and we refit the data we obtain a $\chi^{2}_{\mathrm{dof}}$ of $4.4$. We
show the plots of the charmed meson masses as a function of the pion mass in
lattice units in Fig. 8,

Figure 8: Results of the $D$, $D^{*}$, $D_{s}$ and $D^{*}_{s}$ meson masses
from the PACS analysis (Fit II) after eliminating the fourth pion mass of
Table 11.

were the conflicting data point is shown in red color. The values of the
parameters obtained are given in Tables 2 and 5. Still, the parameters
obtained in Table 2 are not so different from the other data sets. However,
$m_{H}$ in Table 5 is significantly lower, suggesting that the charm quark
mass in these simulations is much lower than the experimental, as also
suggested in Table 1.

Collaboration |  $m_{H}$ (MeV) |  $\Delta$ (MeV)
---|---|---
(I) | (II) | (I) | (II)
PACS | $1767\pm 47$ | $1706\pm 3$ | $99\pm 2$ | $99\pm 2$
Table 5: Tree level parameters of the PACS fit. The remaining parameters are
shown in the first two rows of Table 2.

###  IV.4 Global fit

This fit includes the data of ETMC, HSC, CLS, [Prelovsek et al.] and RQCD
given in Tables 9, 10, 12, 14 and 15. The data of PACS-CS are not included
because of the tensions observed with the other data sets. In this case the
information criteria, shown in Fig. 9 (top), have two minima, for
$\lambda\simeq 0.5$ (global minima) and $\lambda\simeq 2.5$. This corresponds
again to leaving out one parameter, $b$, or three, $b$, $d$ and
$\Delta_{\sigma}$, see Fig. 9 (bottom). The results of the parameters before
and after the LASSO are given in Tables 2 and 6.

Figure 9: LASSO method plots for the global analysis. Collaboration |  $m_{H}$ (MeV) |  $\Delta$ (MeV)
---|---|---
Label |  $a$ (fm) | (I) | (II) | (I) | (II)
ETMC1 | $0.0619$ | $1937\pm 4$ | $1945\pm 4$ | $140\pm 7$ | $141\pm 7$
$0.0815$ | $1949\pm 5$ | $1956\pm 5$ | $148\pm 14$ | $148\pm 14$
$0.0885$ | $1955\pm 5$ | $1963\pm 5$ | $152\pm 12$ | $153\pm 12$
ETMC2 | $0.0619$ | $1950\pm 3$ | $1958\pm 3$ | $139\pm 7$ | $140\pm 7$
$0.0815$ | $1973\pm 5$ | $1980\pm 5$ | $145\pm 13$ | $146\pm 13$
$0.0885$ | $1988\pm 4$ | $1996\pm 4$ | $148\pm 13$ | $149\pm 13$
ETMC | 0 | $1916$ | $1925$ | $131$ | $132$
HSC | ${1/6079}^{*}$ | $1903\pm 5$ | $1910\pm 5$ | $99\pm 7$ | $100\pm 7$
${1/5667}^{*}$ | $1874\pm 5$ | $1881\pm 5$ | $91\pm 8$ | $92\pm 8$
CLS | $0.08636$ | $1940\pm 1$ | $1945\pm 1$ | $95\pm 2$ | $96\pm 2$
Prelovsek | $0.1239$ | $1574\pm 2$ | $1582\pm 2$ | $97\pm 2$ | $98\pm 2$
RQCD | $0.0714$ | $1910\pm 1$ | $1917\pm 1$ | $98\pm 1$ | $98\pm 1$
MILC | $0$ | $1913\pm 4$ | $1920\pm 4$ | $111\pm 4$ | $112\pm 4$
Table 6: Tree level parameters of the global fit. The lattice spacing with a
∗, corresponding to HSC, are the temporal lattice spacing and they are shown
in units of GeV-1. The remaining parameters are shown in the first two rows of
Table 2.

In Figs. 10 and 11 we plot the results for the charmed meson masses in
comparison with the data. In Fig. 10 we plot the spin averaged masses and the
hyperfine splittings of the charmed mesons as a function of the pion mass in
physical units. We leave out ETMC and [Prelovsek et al.] data from this plot
for a better view. In general, the error bands for ETMC data overlap with
those from HSC, and the band for [Prelovsek et al.] data results in much
smaller charmed masses. We show the results for these two data sets separately
in Fig. 11 in units of the lattice spacing. As can be seen all data sets
included in the global fit are very well described. The
$\chi^{2}_{\mathrm{dof}}$ obtained is $0.85$. In Fig. 10 we also show the
result for the global fit fixing the $m_{H}$ and $\Delta$ parameters according
to the experimental data. As can be seen, the physical trajectories for the
spin averaged charmed meson masses are closer to the ones of RQCD data, which
includes an almost physical pion mass. Still, the hyperfine splittings of the
charmed mesons in the lattice data are smaller than the physical ones in $\sim
8-15\%$.

Figure 10: Results of the $D$, $D^{*}$, $D_{s}$ and $D^{*}_{s}$ meson masses
for the global analysis of some collaborations.

Figure 11: As Fig. 10 but with the remaining collaborations.

We also performed an analysis with the tree level expressions, Eqs. (33) and
(34), and the results are given in Appendix B. The $\chi^{2}_{\mathrm{dof}}$
obtained is a bit larger than with the one-loop expressions, $1.11$, but it
also admits a good fit. In general, the values of $m_{H}$ and $\Delta$ are
larger than the ones obtained with the one-loop expressions, as expected. The
values obtained for $\Delta$ are more reasonable for the one-loop expressions
than for the tree level expressions, as can be seen comparing Table 6 and
Table 8, where for ETMC, larger values of $\Delta$ than the physical hyperfine
splitting, are obtained using the tree level expressions.

It is interesting to see how the $m_{H}/M_{\mathrm{av}}^{c\bar{c}}$ and
$\Delta$ parameters, with
$M_{\mathrm{av}}^{c\bar{c}}=(m_{\eta_{c}}+3m_{J/\psi})/4$, $m_{H}$ and
$\Delta$ being the spin average mass and hyperfine splitting in the chiral
limit, depend on the lattice spacing. This is shown in Fig. 12.

Figure 12: $m_{H}$ and $\Delta$ parameters as a function of the lattice
spacing.

From Fig. 12, one obtains that $m_{H}\simeq 0.63M_{\mathrm{av}}^{c\bar{c}}$,
where we expect this relation to be valid for charm quark masses beyond the
physical one, but not very far from the range of the charmed quark masses of
the LQCD simulations studied here333At the physical point, this relation
deviates $\sim 10\%$ for the bottom counterparts assuming that the physical
masses behave as leading order.. Since $\Delta$ is a heavy quark spin
violating parameter, one does not expect it to depend heavily on $m_{H}$. In
Fig. 12, taking into account the dependence of this parameter with the lattice
spacing for the ETMC data, one sees that when $a\to 0$ this would lead to a
range, $\Delta\sim 90-130$ MeV. In the future, more precise LQCD data can help
to determine more precisely this parameter.

##  V Conclusions

We have investigated the light and heavy quark mass dependence of the
$D^{(*)}_{(s)}$ mesons by means of an analysis of the available lattice data
from different LQCD Collaborations withing one-loop HH$\chi$PT. Model
selection tools were implemented to determine the low energy constants with a
higher precision.

We have found that, concerning the higher order parameters,
$b,c,d,\Delta_{\sigma}$ and $\Delta_{a}$, while $c$ and $\Delta_{a}$ are
required by the precision of the LQCD data, $b,d$ and $\Delta_{\sigma}$ are
less relevant, and, in general, LQCD data do not restrict the $b$ parameter.
We have found more precise values of the one-loop parameters, once we have
applied these model selection techniques.

We also observe that, the largest difference in the various simulations comes
mainly from the setting of the charm quark mass to the physical point, and
that this difference becomes also larger due to discretization effects.
Improvement on the methods to set the physical charm quark mass in the LQCD
simulations are required at present, according to the present analysis.
However, once these differences are taken into account, most of the LQCD data
(all except for PACS) can be very well described together by one-loop
HH$\chi$PT, considering different values of the tree level parameters, $m_{H}$
and $\Delta$. The data can be described well also with the tree level
expressions taking into account one insertion of the light quark matrix
instead of two, however, too high values of $\Delta$ are obtained in some
cases.

The description with one loop involving pseudoscalar Goldstone bosons
describes better the LQCD data and provides more reasonable values of the tree
level parameters after the LASSO method is implemented. Once the charm quark
mass is fixed within the theory as required by the experimental data, the
charmed meson masses obtained from the global fit of LQCD data are consistent
with the experimental data.

##  VI Acknowledgments

R. M. and F. Gil thank J. Komijani for providing the MILC continuous
extrapolation to the physical point of the $D_{(s)}$ meson. We also thank A.
Ramos and F. Romero for useful discussions. We acknowledge support from the
CIDEGENT program with Ref. CIDEGENT/2019/015 and from the spanish national
grants PID2019-106080GB-C21 and PID2020-112777GB-I00. This project has
received funding from the European Union’s Horizon 2020 programme No. 824093
for the STRONG-2020 project.

##  Appendix A Coefficients in Eqs. (13) and (14)

The coefficients $\alpha_{l}$, $\beta_{l}^{(X)}$,
$(\gamma_{l}^{(X)}-\lambda_{l}^{(X)}\alpha_{l})$, $c_{l}$,
$(\gamma_{l}^{(X)}-\lambda_{l}^{(X)}\Delta)$ and $\delta c_{l}$, are given by,

$\displaystyle\alpha_{l}=2\sigma(2m+m_{s})+2a(m\delta_{1,2l}+m_{s}\delta_{3l})$
(19) $\displaystyle\beta^{(\pi)}=g^{2}\left({\begin{array}[]{c}\frac{3}{2}\\\
\frac{3}{2}\\\
0\end{array}}\right);\,\beta^{(K)}=g^{2}\left({\begin{array}[]{c}1\\\ 1\\\
2\end{array}}\right);\beta^{(\eta)}=g^{2}\left({\begin{array}[]{c}\frac{1}{6}\\\
\frac{1}{6}\\\ \frac{2}{3}\end{array}}\right)$ (29)
$\displaystyle(\gamma_{l}^{(\pi)}-\lambda_{l}^{(\pi)}\alpha_{l})=$
$\displaystyle-a(3m\delta_{1l}+\delta_{2l})-6\sigma m$
$\displaystyle(\gamma_{l}^{(K)}-\lambda_{l}^{(K)}\alpha_{l})=$
$\displaystyle-a(m+m_{s})(\delta_{1l}+\delta_{2l})-2\sigma(m+m_{s})$
$\displaystyle+3g^{2}a(m_{s}-m)(\delta_{1l}-\delta_{3l})$
$\displaystyle(\gamma_{l}^{(\eta)}-\lambda_{l}^{(\eta)}\alpha_{l})=$
$\displaystyle-a\left\\{\frac{1}{3}m(\delta_{1l}+\delta_{2l})+\frac{4}{3}m_{s}\delta_{3l}\right\\}$
$\displaystyle-\frac{\sigma}{3}(2m+4m_{s})$
$\displaystyle(\gamma_{l}^{(X)}-\lambda_{l}^{(X)}\Delta)=$
$\displaystyle-\Delta\frac{2}{3}\lambda^{(X)}_{l}=-\Delta\beta^{(X)}_{l}$ (30)
$\displaystyle c_{l}=$
$\displaystyle\,4b(m^{2}\delta_{1,2l}+m_{s}^{2}\delta_{3l})+4c(2m+m_{s})(m\delta_{1,2l}+m_{s}\delta_{3l})$
$\displaystyle+4d(2m^{2}+m_{s}^{2})$ $\displaystyle\delta c_{l}=$
$\displaystyle\,2\Delta^{(\sigma)}(2m+m_{s})+2\Delta^{(a)}(m\delta_{1,2l}+m_{s}\delta_{3l})\,$
(31)

where we have assumed the isospin limit $m=m_{u}=m_{d}$ and
$\delta_{1,2l}=\delta_{1,l}+\delta_{2,l}$. For the general formulas taking
into account isospin breaking see Ref. Jenkins (1994). The light quark masses
can be related to the leading order pseudoscalar masses, $M_{X}^{0}$, like

$\displaystyle m=\frac{M^{02}_{\pi}}{2B_{0}},\qquad
m_{s}=\frac{M^{02}_{K}}{B_{0}}-\frac{M^{02}_{\pi}}{2B_{0}}.$ (32)

For simplicity, we identify the physical masses of the light pseudoscalar
meson masses with the leading order pseudoscalar masses, $M_{\pi}=M^{0}_{\pi}$
and $M_{K}=M^{0}_{K}$,
$M_{\eta}^{2}=M_{\eta}^{02}=(4M_{K}^{02}-M_{\pi}^{02})/3$.

##  Appendix B Results for the tree level fits of the charmed meson masses

In this section we present the results of analyzing the charmed meson masses
from LQCD with the tree level expressions including one insertion of the light
quark mass matrix. The relations in these case are,

$\displaystyle\frac{1}{4}(M_{P_{l}}$ $\displaystyle+3M_{P_{l}^{*}})=$
$\displaystyle m_{H}+2\sigma(2m+m_{s})+2a(m\delta_{1,2l}+m_{s}\delta_{3l})\ ,$
(33) $\displaystyle M_{P_{l}}^{*}$ $\displaystyle-M_{P_{l}}=$
$\displaystyle\Delta+2\Delta^{(\sigma)}(2m+m_{s})+2\Delta^{(a)}(m\delta_{1,2l}+m_{s}\delta_{3l})\
.$ (34)

The parameters obtained in the various analyses using these relations are
given in Tables 8 and 7. The $\chi^{2}_{dof}$ obtained is also good, $0.42$
for the ETMC fit and $1.11$ for the Global fit. For MILC we obtain a value
close to zero because of the refitting of the continuously extrapolated
function.

Fit |  $\frac{\sigma m_{\pi}}{{B_{0}}^{2}}$$\cdot 10^{-4}$ |  $\frac{am_{\pi}}{{B_{0}}^{2}}$$\cdot 10^{-4}$ |  $\frac{\Delta^{(\sigma)}m_{\pi}}{{B_{0}}}$$\cdot 10^{-5}$ |  $\frac{\Delta^{(a)}m_{\pi}}{{B_{0}}}$$\cdot 10^{-5}$
---|---|---|---|---
ETMC | $0.16\pm 0.17$ | $2.03\pm 0.07$ | $6.22\pm 4.78$ | $-2.02\pm 2.24$
MILC | $0.83\pm 0.26$ | $2.45\pm 0.37$ | $6.07\pm 1.75$ | $3.75\pm 4.47$
GLOBAL | $0.69\pm 0.01$ | $1.98\pm 0.03$ | $0.48\pm 0.23$ | $-1.88\pm 0.47$
Table 7: Results of the analyses with the tree level formulas, Eqs. (33) and (34). The remaining parameters are shown in Table 8 Fit | Label |  $a$ (fm) |  $m_{H}$ (MeV) |  $\Delta$ (MeV)
---|---|---|---|---
ETMC | ETMC1 | $0.0619$ | $1994\pm 12$ | $137\pm 34$
$0.0815$ | $2010\pm 14$ | $145\pm 39$
$0.0885$ | $2014\pm 14$ | $147\pm 42$
ETMC2 | $0.0619$ | $2006\pm 12$ | $138\pm 33$
$0.0815$ | $2034\pm 13$ | $143\pm 39$
$0.0885$ | $2047\pm 14$ | $146\pm 39$
ETMC | $0$ | $1971$ | $129$
MILC | MILC | $0$ | $1952\pm 21$ | $140$
GLOBAL | ETMC1 | $0.0619$ | $1959\pm 3$ | $178\pm 8$
$0.0815$ | $1972\pm 5$ | $188\pm 16$
$0.0885$ | $1974\pm 4$ | $194\pm 14$
ETMC2 | $0.0619$ | $1971\pm 3$ | $177\pm 8$
$0.0815$ | $1996\pm 5$ | $185\pm 14$
$0.0885$ | $2007\pm 4$ | $190\pm 15$
HSC | $1/6079$ | $1929\pm 5$ | $129\pm 8$
$1/5667$ | $1899\pm 5$ | $124\pm 9$
CLS | $0.08636$ | $1966\pm 1$ | $124\pm 2$
Prelovsek | $0.1239$ | $1596\pm 1$ | $128\pm 3$
RQCD | $0.0714$ | $1934\pm 1$ | $127\pm 2$
MILC | $0$ | $1938\pm 3$ | $144\pm 5$
Table 8: Results of the analyses with the tree level formulas, Eqs. (33) and
(34). The remaining parameters are shown in Table 7.

##  Appendix C Data included in the fits of IV

$a$ (fm) | $(L/a)^{3}\times T/a$ | $am_{\pi}$ | $am_{K}$ | $am_{D}$ | $am_{D_{s}}$ | $am_{D^{*}}$ | $am_{D_{s}^{*}}$ | $am_{J/\Psi}$
---|---|---|---|---|---|---|---|---
$0.0619(18)$ | $48^{3}\times 96$ | 0.0703 | 0.1697 | 0.5905 | 0.6236 | 0.6466 | 0.6770 | 0.9715
$48^{3}\times 96$ | 0.0806 | 0.1738 | 0.5906 | 0.6234 | 0.6506 | 0.6763 | 0.9697
$48^{3}\times 96$ | 0.0975 | 0.1768 | 0.5913 | 0.6229 | 0.6486 | 0.6764 | 0.9703
$0.0815(30)$ | $32^{3}\times 64$ | 0.1074 | 0.2133 | 0.7840 | 0.8159 | 0.8568 | 0.8905 | 1.2791
$32^{3}\times 64$ | 0.1549 | 0.2279 | 0.7895 | 0.8183 | 0.8678 | 0.8950 | 1.2828
$24^{3}\times 48$ | 0.1935 | 0.2430 | 0.7934 | 0.8175 | 0.8745 | 0.8965 | 1.2818
$0.0885(36)$ | $32^{3}\times 64$ | 0.1240 | 0.2512 | 0.8514 | 0.8953 | 0.9356 | 0.9806 | 1.3890
$32^{3}\times 64$ | 0.1412 | 0.2569 | 0.8544 | 0.8972 | 0.9363 | 0.9802 | 1.3895
$24^{3}\times 48$ | 0.1440 | 0.2589 | 0.8552 | 0.8978 | 0.9403 | 0.9844 | 1.3906
| $24^{3}\times 48$ | 0.1988 | 0.2764 | 0.8599 | 0.8950 | 0.9487 | 0.9841 | 1.3882
Table 9: $D^{(*)}_{(s)}$ and $J/\Psi$ meson masses. Data collected from Guo et al. (2018b) provided by the authors of Kalinowski and Wagner (2015). The charm quark mass is fixed by reproducing the physical $J/\psi$ meson mass. Lattice spacing uncertainties are discused on Carrasco et al. (2014). $N_{F}=2+1+1$. $a$ (fm) | $(L/a)^{3}\times T/a$ | $am_{\pi}$ | $am_{K}$ | $am_{D}$ | $am_{D_{s}}$ | $am_{D^{*}}$ | $am_{D_{s}^{*}}$ | $am_{\eta_{c}}$
---|---|---|---|---|---|---|---|---
$0.0619(18)$ | $48^{3}\times 96$ | 0.0703 | 0.1697 | 0.5947 | 0.6279 | 0.6506 | 0.6809 | 0.9351
$48^{3}\times 96$ | 0.0806 | 0.1738 | 0.5949 | 0.6277 | 0.6546 | 0.6803 | 0.9332
$48^{3}\times 96$ | 0.0975 | 0.1768 | 0.5955 | 0.6271 | 0.6526 | 0.6804 | 0.9335
$0.0815(30)$ | $32^{3}\times 64$ | 0.1074 | 0.2133 | 0.7946 | 0.8263 | 0.8664 | 0.9001 | 1.2312
$32^{3}\times 64$ | 0.1549 | 0.2279 | 0.8004 | 0.8291 | 0.8777 | 0.9049 | 1.2342
$24^{3}\times 48$ | 0.1935 | 0.2430 | 0.8039 | 0.8278 | 0.8840 | 0.9059 | 1.2314
$0.0885(36)$ | $32^{3}\times 64$ | 0.1240 | 0.2512 | 0.8677 | 0.9114 | 0.9506 | 0.9953 | 1.3370
$32^{3}\times 64$ | 0.1412 | 0.2569 | 0.8708 | 0.9132 | 0.9511 | 0.9949 | 1.3379
$24^{3}\times 48$ | 0.1440 | 0.2589 | 0.8714 | 0.9137 | 0.9545 | 0.9990 | 1.3382
| $24^{3}\times 48$ | 0.1988 | 0.2764 | 0.8753 | 0.9102 | 0.9627 | 0.9980 | 1.3325
Table 10: The same as Table 9 but with $D^{(*)}_{(s)}$ and $\eta_{c}$ meson masses. The charm quark mass is fixed by reproducing the physical $\eta_{c}$ meson mass. $a$ (fm) | $am_{\pi}$ | $am_{K}$ | $am_{D}$ | $am_{D_{s}}$ | $am_{D^{*}}$ | $am_{D_{s}^{*}}$ | $am_{\eta_{c}}$ | $am_{J/\Psi}$
---|---|---|---|---|---|---|---|---
$a=0.0907(13)$ | 0.32242 | 0.36269 | 0.84022 | 0.86273 | 0.91237 | 0.93297 | 1.22924 | 1.27877
0.26191 | 0.32785 | - | 0.85158 | - | 0.91864 | 1.22098 | 1.26860
0.18903 | 0.29190 | 0.79580 | 0.84000 | 0.86327 | 0.90429 | 1.21369 | 1.25988
| 0.13593 | 0.27282 | 0.78798 | 0.83929 | 0.85776 | 0.90429 | 1.21268 | 1.25904
| 0.07162 | 0.25454 | 0.77646 | 0.83149 | 0.83656 | 0.89268 | 1.20765 | 1.25258
Table 11: Data for heavy meson masses from Mohler and Woloshyn (2011) and
light meson masses collected from Aoki et al. (2009) [PACS-CS]. The charm
quark hopping pa rameter $\kappa_{c}$ has been tuned to the value where the
spin-averaged kinetic mass $(M_{D_{s}}+3M_{D_{s}^{*}})/4$ takes its physical
value. The lowest pion mass data is also named as Ensemble 2 in Lang et al.
(2014). Also, $(L/a)^{3}\times T/a=32^{3}\times 64$ and $N_{F}=2+1$.

$a_{t}^{-1}$ (MeV) | $(L/a_{s})^{3}\times(T/a_{t})$ | $m_{\pi}$ | $m_{K}$ | $m_{D}$ | $m_{D_{s}}$ | $m_{D^{*}}$ | $m_{D_{s}^{*}}$
---|---|---|---|---|---|---|---
6079 | $32^{3}\times 256$ | 0.03928(18) | 0.08344(7) | 0.30923(11) | 0.32356(12) | 0.33058(24) | 0.34448(15)
5667 | $(16,20,24)^{3}\times 128$ | 0.06906(13) | 0.09698(9) | 0.33303(31) | 0.34441(29) | 0.35494(46) | 0.36587(35)
Table 12: Data collected from Cheung et al. (2021) for $N_{F}=2+1$ [HSC]. The charm quark mass parameters have been determined using the $\eta_{c}$ physical mass, obtaining the values, $(a_{t}m_{\eta_{c}})^{2}=0.2412$ and $0.2735$ for the light and heavy pion mass respectively (see Figs. 1 of Cheung et al. (2016); Liu et al. (2012)).  $a$ (fm) | $m_{\pi}$ | $m_{K}$ | $m_{D}$ | $m_{D_{s}}$ | $m_{D^{*}}$ | $m_{D_{s}^{*}}$ | $M_{avg}$
---|---|---|---|---|---|---|---
$a=0.08636(138)$ | 280 | 467 | 1927 | 1981 | 2048 | 2102 | 3103
Table 13: Data collected from Prelovsek et al. (2021) where the spin average $M_{avg}=(m_{\eta_{c}}+3m_{J/\Psi})/4$ has been tune to the physical point, for $(L/a)^{3}\times T/a=24^{3}\times 128$ and $N_{F}=2+1$ [CLS]. Masses given in MeV. $a$ (fm) | $am_{\pi}$ | $am_{K}$ | $am_{D}$ | $am_{D_{s}}$ | $am_{D^{*}}$ | $am_{D_{s}^{*}}$ | $aM_{avg}$
---|---|---|---|---|---|---|---
$a=0.1239(13)$ | 0.1673 | 0.3467 | 0.9801 | 1.04075433 | 1.0629 | 1.11635248 | 1.52499
Table 14: Data collected from Ensemble 1 in Lang et al. (2014) for $(L/a)^{3}\times T/a=16^{3}\times 32$ and $N_{F}=2$ [Prelovsek et al.]. The spin-average $(m_{\eta_{c}}+3m_{J/\psi})/4$ is tuned to reproduce its physical value.  $a$ (fm) | $(L/a)^{3}\times T/a$ | $m_{\pi}$ | $m_{K}$ | $m_{D}$ | $m_{D_{s}}$ | $m_{D^{*}}$ | $m_{D_{s}^{*}}$
---|---|---|---|---|---|---|---
$a=0.071$ | $24^{3}\times 48$ | 306.9 | 540 | 1907 | - | 2038 | -
$32^{3}\times 64$ | 294.6 | 528 | 1902 | - | 2030 | -
$40^{3}\times 64$ | 288.8 | 527 | 1901 | - | 2030 | -
$64^{3}\times 64$ | 289.5 | 526 | 1898 | - | 2030 | -
$48^{3}\times 64$ | 159.7 | 500 | 1880 | - | 2007 | -
$64^{3}\times 64$ | 149.7 | 497 | 1877 | 1976.9 | 1996 | 2094.9
Table 15: Data collected from table I in Bali et al. (2017) for $N_{F}=2$
[RQCD]. The charm quark mass is fixed through the experimental value of the
spin-averaged $1S$ charmonium mass, $m_{1S}=3068.5$ MeV. Masses given in MeV.

## References

  * Guo et al. (2018a) F.-K. Guo, C. Hanhart, U.-G. Meißner, Q. Wang, Q. Zhao, and B.-S. Zou, Rev. Mod. Phys. 90, 015004 (2018a), [Erratum: Rev.Mod.Phys. 94, 029901 (2022)], eprint 1705.00141.
  * Truong (1988) T. N. Truong, Phys. Rev. Lett. 61, 2526 (1988).
  * Dobado and Pelaez (1993) A. Dobado and J. R. Pelaez, Phys. Rev. D 47, 4883 (1993), eprint hep-ph/9301276.
  * Dobado and Pelaez (1997) A. Dobado and J. R. Pelaez, Phys. Rev. D 56, 3057 (1997), eprint hep-ph/9604416.
  * Nieves and Ruiz Arriola (2000) J. Nieves and E. Ruiz Arriola, Nucl. Phys. A 679, 57 (2000), eprint hep-ph/9907469.
  * Garcia-Recio et al. (2004) C. Garcia-Recio, M. F. M. Lutz, and J. Nieves, Phys. Lett. B 582, 49 (2004), eprint nucl-th/0305100.
  * Oller et al. (1999) J. A. Oller, E. Oset, and J. R. Pelaez, Phys. Rev. D 59, 074001 (1999), [Erratum: Phys.Rev.D 60, 099906 (1999), Erratum: Phys.Rev.D 75, 099903 (2007)], eprint hep-ph/9804209.
  * Gomez Nicola and Pelaez (2002) A. Gomez Nicola and J. R. Pelaez, Phys. Rev. D 65, 054009 (2002), eprint hep-ph/0109056.
  * Oller and Oset (1997) J. A. Oller and E. Oset, Nucl. Phys. A 620, 438 (1997), [Erratum: Nucl.Phys.A 652, 407–409 (1999)], eprint hep-ph/9702314.
  * Oller et al. (1998) J. A. Oller, E. Oset, and J. R. Pelaez, Phys. Rev. Lett. 80, 3452 (1998), eprint hep-ph/9803242.
  * Kaiser (1998) N. Kaiser, Eur. Phys. J. A 3, 307 (1998).
  * Locher et al. (1998) M. P. Locher, V. E. Markushin, and H. Q. Zheng, Eur. Phys. J. C 4, 317 (1998), eprint hep-ph/9705230.
  * Pelaez (2016) J. R. Pelaez, Phys. Rept. 658, 1 (2016), eprint 1510.00653.
  * Aoki et al. (2007) S. Aoki et al. (CP-PACS), Phys. Rev. D76, 094506 (2007), eprint 0708.3705.
  * Gockeler et al. (2008) M. Gockeler, R. Horsley, Y. Nakamura, D. Pleiter, P. E. L. Rakow, G. Schierholz, and J. Zanotti (QCDSF), PoS LATTICE2008, 136 (2008), eprint 0810.5337.
  * Feng et al. (2011) X. Feng, K. Jansen, and D. B. Renner, Phys. Rev. D83, 094505 (2011), eprint 1011.5288.
  * Lang et al. (2011) C. B. Lang, D. Mohler, S. Prelovsek, and M. Vidmar, Phys. Rev. D84, 054503 (2011), [Erratum: Phys. Rev.D89,no.5,059903(2014)], eprint 1105.5636.
  * Pelissier and Alexandru (2013) C. Pelissier and A. Alexandru, Phys. Rev. D87, 014503 (2013), eprint 1211.0092.
  * Bali et al. (2016) G. S. Bali, S. Collins, A. Cox, G. Donald, M. Göckeler, C. B. Lang, and A. Schäfer (RQCD), Phys. Rev. D93, 054509 (2016), eprint 1512.08678.
  * Guo et al. (2016) D. Guo, A. Alexandru, R. Molina, and M. Döring, Phys. Rev. D 94, 034501 (2016), eprint 1605.03993.
  * Erben et al. (2020) F. Erben, J. R. Green, D. Mohler, and H. Wittig, Phys. Rev. D 101, 054504 (2020), eprint 1910.01083.
  * Wilson et al. (2015) D. J. Wilson, R. A. Briceño, J. J. Dudek, R. G. Edwards, and C. E. Thomas, Phys. Rev. D92, 094502 (2015), eprint 1507.02599.
  * Dudek et al. (2013) J. J. Dudek, R. G. Edwards, and C. E. Thomas (Hadron Spectrum), Phys. Rev. D87, 034505 (2013), [Erratum: Phys. Rev.D90,no.9,099902(2014)], eprint 1212.0830.
  * Bulava et al. (2016) J. Bulava, B. Fahy, B. Horz, K. J. Juge, C. Morningstar, and C. H. Wong, Nucl. Phys. B910, 842 (2016), eprint 1604.05593.
  * Feng et al. (2015) X. Feng, S. Aoki, S. Hashimoto, and T. Kaneko, Phys. Rev. D91, 054504 (2015), eprint 1412.6319.
  * Alexandrou et al. (2017) C. Alexandrou, L. Leskovec, S. Meinel, J. Negele, S. Paul, M. Petschlies, A. Pochinsky, G. Rendon, and S. Syritsyn, Phys. Rev. D96, 034525 (2017), eprint 1704.05439.
  * Fu and Wang (2016) Z. Fu and L. Wang, Phys. Rev. D94, 034505 (2016), eprint 1608.07478.
  * Metivet (2015) T. Metivet (Budapest-Marseille-Wuppertal), PoS LATTICE2014, 079 (2015), eprint 1410.8447.
  * Fischer et al. (2020) M. Fischer, B. Kostrzewa, M. Mai, M. Petschlies, F. Pittler, M. Ueding, C. Urbach, and M. Werner (ETM) (2020), eprint 2006.13805.
  * Hu et al. (2016) B. Hu, R. Molina, M. Döring, and A. Alexandru, Phys. Rev. Lett. 117, 122001 (2016), eprint 1605.04823.
  * Hu et al. (2017) B. Hu, R. Molina, M. Döring, M. Mai, and A. Alexandru, Phys. Rev. D 96, 034520 (2017), eprint 1704.06248.
  * Molina and Ruiz de Elvira (2020) R. Molina and J. Ruiz de Elvira, JHEP 11, 017 (2020), eprint 2005.13584.
  * Olsen et al. (2018) S. L. Olsen, T. Skwarnicki, and D. Zieminska, Rev. Mod. Phys. 90, 015003 (2018), eprint 1708.04012.
  * Chen et al. (2022) H.-X. Chen, W. Chen, X. Liu, Y.-R. Liu, and S.-L. Zhu (2022), eprint 2204.02649.
  * Barnes and Swanson (2008) T. Barnes and E. S. Swanson, Phys. Rev. C 77, 055206 (2008), eprint 0711.2080.
  * Isgur and Wise (1989) N. Isgur and M. B. Wise, Phys. Lett. B 232, 113 (1989).
  * Isgur and Wise (1990) N. Isgur and M. B. Wise, Phys. Lett. B 237, 527 (1990).
  * Isgur and Wise (1991) N. Isgur and M. B. Wise, Phys. Rev. Lett. 66, 1130 (1991).
  * Manohar and Wise (2000) A. V. Manohar and M. B. Wise, _Heavy quark physics_ , vol. 10 (2000), ISBN 978-0-521-03757-0.
  * Godfrey and Isgur (1985) S. Godfrey and N. Isgur, Phys. Rev. D 32, 189 (1985).
  * Godfrey and Kokoski (1991) S. Godfrey and R. Kokoski, Phys. Rev. D 43, 1679 (1991).
  * Di Pierro and Eichten (2001) M. Di Pierro and E. Eichten, Phys. Rev. D 64, 114004 (2001), eprint hep-ph/0104208.
  * Aubert et al. (2003) B. Aubert et al. (BaBar), Phys. Rev. Lett. 90, 242001 (2003), eprint hep-ex/0304021.
  * Besson et al. (2003) D. Besson et al. (CLEO), Phys. Rev. D 68, 032002 (2003), [Erratum: Phys.Rev.D 75, 119908 (2007)], eprint hep-ex/0305100.
  * Barnes et al. (2003) T. Barnes, F. E. Close, and H. J. Lipkin, Phys. Rev. D 68, 054006 (2003), eprint hep-ph/0305025.
  * Kolomeitsev and Lutz (2004) E. E. Kolomeitsev and M. F. M. Lutz, Phys. Lett. B 582, 39 (2004), eprint hep-ph/0307133.
  * Faessler et al. (2007a) A. Faessler, T. Gutsche, V. E. Lyubovitskij, and Y.-L. Ma, Phys. Rev. D 76, 014005 (2007a), eprint 0705.0254.
  * Faessler et al. (2007b) A. Faessler, T. Gutsche, V. E. Lyubovitskij, and Y.-L. Ma, Phys. Rev. D 76, 114008 (2007b), eprint 0709.3946.
  * Lutz and Soyeur (2008) M. F. M. Lutz and M. Soyeur, Nucl. Phys. A 813, 14 (2008), eprint 0710.1545.
  * Albaladejo et al. (2018) M. Albaladejo, P. Fernandez-Soler, J. Nieves, and P. G. Ortega, Eur. Phys. J. C 78, 722 (2018), eprint 1805.07104.
  * Cheng and Hou (2003) H.-Y. Cheng and W.-S. Hou, Phys. Lett. B 566, 193 (2003), eprint hep-ph/0305038.
  * Terasaki (2003) K. Terasaki, Phys. Rev. D 68, 011501 (2003), eprint hep-ph/0305213.
  * Maiani et al. (2005) L. Maiani, F. Piccinini, A. D. Polosa, and V. Riquer, Phys. Rev. D 71, 014028 (2005), eprint hep-ph/0412098.
  * Dmitrasinovic (2005) V. Dmitrasinovic, Phys. Rev. Lett. 94, 162002 (2005).
  * Molina et al. (2010) R. Molina, T. Branz, and E. Oset, Phys. Rev. D 82, 014010 (2010), eprint 1005.0335.
  * Aaij et al. (2020a) R. Aaij et al. (LHCb), Phys. Rev. Lett. 125, 242001 (2020a), eprint 2009.00025.
  * Aaij et al. (2020b) R. Aaij et al. (LHCb), Phys. Rev. D 102, 112003 (2020b), eprint 2009.00026.
  * Collaboration (2022) L. Collaboration (2022), eprint 2212.02717.
  * Molina and Oset (2020) R. Molina and E. Oset, Phys. Lett. B 811, 135870 (2020), eprint 2008.11171.
  * Molina and Oset (2022) R. Molina and E. Oset (2022), eprint 2211.01302.
  * Abe et al. (2004) K. Abe et al. (Belle), Phys. Rev. D 69, 112002 (2004), eprint hep-ex/0307021.
  * Aubert et al. (2009) B. Aubert et al. (BaBar), Phys. Rev. D 79, 112004 (2009), eprint 0901.1291.
  * Link et al. (2004) J. M. Link et al. (FOCUS), Phys. Lett. B 586, 11 (2004), eprint hep-ex/0312060.
  * Guo et al. (2006) F.-K. Guo, P.-N. Shen, H.-C. Chiang, R.-G. Ping, and B.-S. Zou, Phys. Lett. B 641, 278 (2006), eprint hep-ph/0603072.
  * Guo et al. (2009) F.-K. Guo, C. Hanhart, and U.-G. Meissner, Eur. Phys. J. A 40, 171 (2009), eprint 0901.1597.
  * Albaladejo et al. (2017) M. Albaladejo, P. Fernandez-Soler, F.-K. Guo, and J. Nieves, Phys. Lett. B 767, 465 (2017), eprint 1610.06727.
  * Du et al. (2018) M.-L. Du, M. Albaladejo, P. Fernández-Soler, F.-K. Guo, C. Hanhart, U.-G. Meißner, J. Nieves, and D.-L. Yao, Phys. Rev. D 98, 094018 (2018), eprint 1712.07957.
  * Moir et al. (2016) G. Moir, M. Peardon, S. M. Ryan, C. E. Thomas, and D. J. Wilson, JHEP 10, 011 (2016), eprint 1607.07093.
  * Gayer et al. (2021) L. Gayer, N. Lang, S. M. Ryan, D. Tims, C. E. Thomas, and D. J. Wilson (Hadron Spectrum), JHEP 07, 123 (2021), eprint 2102.04973.
  * Wise (1992) M. B. Wise, Phys. Rev. D 45, R2188 (1992).
  * Burdman and Donoghue (1992) G. Burdman and J. F. Donoghue, Phys. Lett. B 280, 287 (1992).
  * Yan et al. (1992) T.-M. Yan, H.-Y. Cheng, C.-Y. Cheung, G.-L. Lin, Y. C. Lin, and H.-L. Yu, Phys. Rev. D 46, 1148 (1992), [Erratum: Phys.Rev.D 55, 5851 (1997)].
  * Rosner and Wise (1993) J. L. Rosner and M. B. Wise, Phys. Rev. D 47, 343 (1993).
  * Randall and Sather (1993) L. Randall and E. Sather, Phys. Lett. B 303, 345 (1993), eprint hep-ph/9211267.
  * Di Bartolomeo et al. (1995) N. Di Bartolomeo, R. Gatto, F. Feruglio, and G. Nardulli, Phys. Lett. B 347, 405 (1995), eprint hep-ph/9411210.
  * Blok et al. (1997) B. Blok, J. G. Korner, D. Pirjol, and J. C. Rojas, Nucl. Phys. B 496, 358 (1997), eprint hep-ph/9607233.
  * Jenkins (1994) E. E. Jenkins, Nucl. Phys. B 412, 181 (1994), eprint hep-ph/9212295.
  * Yeh and Lee (1995) T.-W. Yeh and C.-e. Lee (1995), eprint hep-ph/9510241.
  * Mehen and Springer (2005) T. Mehen and R. P. Springer, Phys. Rev. D 72, 034006 (2005), eprint hep-ph/0503134.
  * Ananthanarayan et al. (2007) B. Ananthanarayan, S. Banerjee, K. Shivaraj, and A. Upadhyay, Phys. Lett. B 651, 124 (2007), eprint 0706.0942.
  * Alhakami (2016) M. H. Alhakami, Phys. Rev. D 93, 094007 (2016), eprint 1603.08848.
  * Aaij et al. (2022) R. Aaij et al. (LHCb), Nature Phys. 18, 751 (2022), eprint 2109.01038.
  * James and Tibshirani (2009) D. W. T. H. James, G. and R. Tibshirani, _The Elements of Statistical Learning_ , vol. XXII (Springer-Verlag, New York, USA, 2009), ISBN 978-0-387-84857-0.
  * Tibshirani (1996) R. Tibshirani, Journal of the Royal Statistical Society. Series B (Methodological) 58, 267 (1996), ISSN 00359246, URL http://www.jstor.org/stable/2346178.
  * Landay et al. (2017) J. Landay, M. Döring, C. Fernández-Ramírez, B. Hu, and R. Molina, Phys. Rev. C 95, 015203 (2017), eprint 1610.07547.
  * Kalinowski and Wagner (2015) M. Kalinowski and M. Wagner, Phys. Rev. D 92, 094508 (2015), eprint 1509.02396.
  * Carrasco et al. (2014) N. Carrasco et al. (European Twisted Mass), Nucl. Phys. B 887, 19 (2014), eprint 1403.4504.
  * Mohler and Woloshyn (2011) D. Mohler and R. M. Woloshyn, Phys. Rev. D 84, 054505 (2011), eprint 1103.5506.
  * Aoki et al. (2009) S. Aoki et al. (PACS-CS), Phys. Rev. D 79, 034503 (2009), eprint 0807.1661.
  * Cheung et al. (2016) G. K. C. Cheung, C. O’Hara, G. Moir, M. Peardon, S. M. Ryan, C. E. Thomas, and D. Tims (Hadron Spectrum), JHEP 12, 089 (2016), eprint 1610.01073.
  * Cheung et al. (2021) G. K. C. Cheung, C. E. Thomas, D. J. Wilson, G. Moir, M. Peardon, and S. M. Ryan (Hadron Spectrum), JHEP 02, 100 (2021), eprint 2008.06432.
  * Prelovsek et al. (2021) S. Prelovsek, S. Collins, D. Mohler, M. Padmanath, and S. Piemonte, JHEP 06, 035 (2021), eprint 2011.02542.
  * Bruno et al. (2017) M. Bruno, T. Korzec, and S. Schaefer, Phys. Rev. D 95, 074504 (2017), eprint 1608.08900.
  * Bali et al. (2017) G. S. Bali, S. Collins, A. Cox, and A. Schäfer, Phys. Rev. D 96, 074501 (2017), eprint 1706.01247.
  * Bali et al. (2015) G. S. Bali, S. Collins, B. Glässle, M. Göckeler, J. Najjar, R. H. Rödl, A. Schäfer, R. W. Schiel, W. Söldner, and A. Sternbeck, Phys. Rev. D 91, 054501 (2015), eprint 1412.7336.
  * Guo et al. (2018b) X.-Y. Guo, Y. Heo, and M. F. M. Lutz, Phys. Rev. D 98, 014510 (2018b), eprint 1801.10122.
  * Bazavov et al. (2018a) A. Bazavov et al. (Fermilab Lattice, MILC, TUMQCD), Phys. Rev. D 98, 054517 (2018a), eprint 1802.04248.
  * Bazavov et al. (2018b) A. Bazavov et al., Phys. Rev. D 98, 074512 (2018b), eprint 1712.09262.
  * Georgi (1990) H. Georgi, Phys. Lett. B 240, 447 (1990).
  * Jenkins and Manohar (1991) E. E. Jenkins and A. V. Manohar, Phys. Lett. B 255, 558 (1991).
  * Lang et al. (2014) C. B. Lang, L. Leskovec, D. Mohler, S. Prelovsek, and R. M. Woloshyn, Phys. Rev. D 90, 034510 (2014), eprint 1403.8103.
  * Brambilla et al. (2018) N. Brambilla, J. Komijani, A. S. Kronfeld, and A. Vairo (TUMQCD), Phys. Rev. D 97, 034503 (2018), eprint 1712.04983.
  * Zyla et al. (2020) P. A. Zyla et al. (Particle Data Group), PTEP 2020, 083C01 (2020).
  * Ramos (2019) A. Ramos, Comput. Phys. Commun. 238, 19 (2019), eprint 1809.01289.
  * Liu et al. (2012) L. Liu, G. Moir, M. Peardon, S. M. Ryan, C. E. Thomas, P. Vilaseca, J. J. Dudek, R. G. Edwards, B. Joo, and D. G. Richards (Hadron Spectrum), JHEP 07, 126 (2012), eprint 1204.5425.