On the interference of $ggH$ and $c\bar{c}H$ Higgs production mechanisms and the determination of charm Yukawa coupling at the LHC

Higgs boson production in association with a charm-quark jet proceeds through two different mechanisms - one that involves the charm Yukawa coupling and the other that involves direct Higgs coupling to gluons. The interference of the two contributions requires a helicity flip and, therefore, cannot be computed with massless charm quarks. In this paper, we consider QCD corrections to the interference contribution starting from charm-gluon collisions with massive charm quarks and taking the massless limit, $m_c \to 0$. The behavior of QCD cross sections in that limit differs from expectations based on the canonical QCD factorization. This implies that QCD corrections to the interference term necessarily involve logarithms of the ratio $M_H/m_c$ whose resummation is currently unknown. Although the explicit next-to-leading order QCD computation does confirm the presence of up to two powers of $\ln(M_H/m_c)$ in the interference contribution, their overall impact on the magnitude of QCD corrections to the interference turns out to be moderate due to a cancellation between double and single logarithmic terms.


Introduction
Studies of Yukawa couplings play an important role in the verification of the mechanism of electroweak symmetry breaking as described by the Standard Model. By now, Higgs couplings to bottom and top quarks, as well as to tau leptons and muons, have been measured to a precision of about twenty percent [1][2][3][4][5][6]. Within the error bars, the measured values for all four Yukawa couplings are consistent with the Standard Model predictions.
However, the Yukawa couplings to lighter fermions have not been studied experimentally. Although it is generally agreed that the Yukawa couplings of electrons and up, down and strange quarks can be observed if and only if they enormously deviate from their Standard Model values, the situation with the charm Yukawa is not so hopeless. In fact, it appears that with the full LHC luminosity, the charm Yukawa coupling can be measured if its value deviates from the Standard Model expectation by an order one factor [7]. Different observables to measure the charm Yukawa coupling at the LHC have been proposed; they include inclusive (H → cc) and exclusive (H → J/ψ + γ and similar) decays of the Higgs boson [8][9][10], the modifications of the Higgs transverse momentum distribution [11] in the gg → H + X process and, finally, Higgs boson production cross section in association with a charm jet [12].
In this paper we focus on the latter process, pp → H + jet c . At leading order in perturbative QCD, Higgs bosons are produced in association with charm jets in the partonic process cg → Hc. The amplitude of this process receives contributions proportional to the charm Yukawa coupling and to an effective ggH coupling M ∼ g Yuk M 1 + g ggH M 2 , (1.1) see Figure 1. As a result, the pp → H + jet c cross section contains the interference term σ Hc ∼ g 2 Yukσ 1 + g 2 ggHσ 2 + g Yuk g ggHσInt . (1. 2) It can be expected that a reliable description of Higgs boson production in association with a charm jet can be obtained by systematically computing the different terms in Eq. (1.2) to higher orders in perturbative QCD. In fact, it is emphasized in Ref. [12] that the largest theoretical uncertainty in using H + jet c production cross section to constrain charm Yukawa coupling is related to perturbative QCD uncertainties so that it seems natural to compute higher order QCD corrections to σ Hc in Eq. (1.2). However, pursuing this program for the interference term in Eq. (1.2) is quite subtle as we now discuss. Indeed, perturbative computations in QCD are performed with massless incoming partons. In case of the massless charm quark that, however, has non-vanishing Yukawa coupling to the Higgs boson, the interference term in Eq. (1.2) vanishes and we obtain lim mc→0 σ Hc ∼ g 2 Yukσ 1 + g 2 ggHσ 2 .  Figure 1: Leading-order Feynman diagrams contributing to the pp → Hc process. We distinguish two separate production mechanisms: one that is driven by the Yukawa coupling (left) and the other one that requires direct coupling of Higgs to gluons (right).
For the massive charm quark the interference does not vanish and is proportional to the charm mass in the first power. Whether or not the interference contribution is negligible depends on the relative magnitude of the two amplitudes in Eq. (1.1). Leading-order computations with massless quarks show that the charm-Yukawa independent amplitude g ggH M 2 in Eq. (1.1) is larger than the charm-Yukawa dependent one g Yuk M 1 suggesting that the interference may be non-negligible. It is straightforward to calculate the interference at leading order in perturbative QCD. Indeed, the interference requires one helicity flip on a charm line that connects initial and final states; this flip is accomplished by a single mass insertion. This implies that one can compute the interference of the two amplitudes using massive charm quarks, take the m c → 0 limit and account for the first non-vanishing term proportional to m c . Since we require a charm jet in the final state, none of the kinematic invariants of the cg → Hc process can be small. Hence, once one power of m c is extracted, the rest of the leading-order calculation of the interference contribution can be performed using the standard approximation of massless (charm) quarks. Such calculation, that we describe in Section 5, shows that the leading-order interference amounts to about ten percent of the contribution to the H + jet c cross section that is proportional to the Yukawa coupling squared.
Although the interference contribution is not large, it is worth thinking about it at nextto-leading order (NLO) in perturbative QCD since there are reasons to believe that the interference contribution is perturbatively unstable, at variance with the two other contributions to pp → H + jet c cross section. Indeed, even if we require an energetic charm jet in the final state, soft and collinear kinematic configurations lead to logarithmic sensitivity of the interference to the charm mass m c . Hence, before the m c → 0 approximation can be taken, the quasi-singular contributions proportional to logarithms of m c have to be extracted from both real and virtual corrections to the interference part of the production cross section.
One may argue that, since the finite charm mass provides yet another way to regulate collinear divergences, it is to be expected that the procedure described above will lead to a familiar picture of (quasi)-collinear factorization of QCD amplitudes. If so, all ln(m c )dependent terms should disappear once infrared safe cross sections and distributions are computed using short-distance quantities, including conventional parton distribution functions (PDFs). However, we will show that for the interference contribution this expectation is invalid and that well-known formulas that describe collinear factorization of mass singularities are not applicable in that case. We will also show that the helicity flip leads to an appearance of soft-quark singularities that, interestingly, make jet algorithms logarithmicallysensitive to m c .
There are two consequences of the above discussion. First, the problem of estimating the magnitude of the interference contribution to the production of a Higgs boson in association with a charm jet turns into an interesting problem in perturbative QCD that borders on such important issues as soft and collinear QCD factorization for mass power corrections [13][14][15][16][17]. Second, a more complex pattern of this factorization, as compared to the canonical collinear and soft cases [18], implies that NLO QCD corrections to leading-order interference are enhanced by up to two powers of a large logarithm ln Q/m c where Q is a typical hard scale in the process pp → H + jet c . For this reason NLO QCD corrections to the interference may be expected to be significant and it becomes essential to explicitly compute them. This is what we set out to do in this paper.
The rest of the paper is organized as follows. In the next section we derive a relation between MS-regulated and mass-regulated parton distribution functions at O(α s ) using the process of Higgs boson production in cc annihilation. We use the established relation to remove "conventional" collinear logarithms from NLO QCD corrections to the interference contribution to the production of Higgs boson in association with charm jet. In Section 3 we discuss factorization of mass singularities in the interference contribution to cg → Hc process and show that it works differently as compared to the standard case [18]. In Section 4 we briefly describe the technical details of the calculation of NLO QCD corrections to the interference contribution. In Section 5 we present phenomenological results and discuss the relative importance of logarithmically-enhanced terms. We conclude in Section 6. Additional discussion of soft and collinear limits of the interference contributions as well as some relevant soft integrals can be found in several appendices.

Matching parton distribution functions
It is well-known that quark masses screen collinear singularities. For this reason we can think about small quark masses as a particular choice of a collinear regulator. Since, when describing "leading-twist" inclusive partonic processes, collinear sensitivity either cancels out or is absorbed into parton distribution functions, it is possible to derive relations between parton distribution functions that are used for computations with nearly massive and strictly massless quarks by requiring that predictions for physical processes are independent of a collinear regulator. 1 To derive a relation between "massive" and "massless" PDFs, we start with the production of a Higgs boson in an annihilation of two massive charm quarks and write the differential cross section as are parton distribution functions and the superscript m implies that all relevant quantities should be computed using quark masses as collinear regulators. Also,σ (m) ij→H+X is the partonic differential cross-section. At leading order i(j) = c,c; at higher orders other channels also contribute.
Calculation at leading order in α s is straightforward since the leading-order cross section σ (m) cc→H has a regular m c → 0 limit. It follows that at leading order in α s there is no difference between f (m) i and conventional MS parton distribution functions, i.e. f The situation becomes more complicated at next-to-leading order where the charm quark mass screens collinear singularities; hence, our goal is to re-write the NLO QCD contributions to the cross section cc → H + X in such a way that logarithms of m c are extracted explicitly.
We begin by considering the process c(p 1 ) +c(p 2 ) → H + g(p 3 ) and treating charm quarks as massive. Kinematic regions that lead to soft and (quasi-)collinear singularities are well understood. The behavior of matrix elements in these limits is described by conventional factorization formulas [18]. We can define a hard m c -independent cross section by subtracting the singular limits. When the subtracted terms are added back and integrated over unresolved parts of the Hg phase space, logarithms of the charm mass appear. This procedure is identical to methods developed to extract infrared and collinear singularities from real emission contributions to partonic cross sections. Its application in the present context allows us to explicitly extract logarithms of the charm mass.
To organize the calculation, we follow the nested soft-collinear subtraction scheme [20][21][22][23] which, at next-to-leading order, is equivalent to the FKS scheme [24,25]. We use dimensional regularization 2 to regularize soft singularities and the charm mass to regularize the collinear ones. Using notations from Ref. [20], we write the partonic cross section for c(p 1 ) +c(p 2 ) → H + g(p 3 ) as The key observation is that since the fully-regulated (last) term in Eq. (2.2) is free of both soft and quasi-collinear singularities, the limit m c → 0 can be safely taken there. On the contrary, both soft and collinear subtraction terms exhibit mass singularities; these mass singularities need to be extracted.
Consider the soft limit defined as p 3 ·p 1 ∼ p 3 ·p 2 → 0. It reads [18] where g s is the unrenormalized strong coupling constant. Since the soft gluon decouples from the function F LM (1 c , 2c) we can integrate Eq. (2.3) over the gluon phase space. We work in the center-of-mass frame of the colliding charm partons and parametrize their energies as E 1 = E 2 = E. The center of mass energy squared in the massless approximation is then s = 4E 2 . We also cut integrals over gluon energy at E 3 = E max , cf. Ref [20]. Integrating Eq. (2.3) over gluon phase space [dg 3 ] and taking the m c → 0 limit, we find and Ω (d−2) is the solid angle of the (d − 2)-dimensional space. 3 The notation1 c (2c) implies that the corresponding four-momenta should be taken in the massless approximation. The two integrals I 1m(2m) in Eq. (2.4) read where β = 1 − m 2 c /E 2 and we neglected all power-suppressed terms when writing the results.
Collinear subtraction terms contain quasi-collinear singularities. The two collinear limits correspond to two distinct cases, p 1 ·p 3 ∼ m 2 c → 0 and p 2 ·p 3 ∼ m 2 c → 0. They read and is the collinear splitting function. The variable z is defined as z = (E i − E 3 )/E i with i = 1, 2, as appropriate.
Since, when computing the collinear limits, we do not change the gluon phase space [20], integrated collinear subtraction terms are still described by angular integrals I 1m(2m) shown in Eq. (2.5). Performing the soft subtraction of the collinear-subtracted cross section, we find andP and L 1 = ln(E max /E). A similar expression can be written for It is straightforward to combine soft and collinear contributions and to expand them in . At this point, it is convenient to switch to a strong coupling constant renormalized at the scale µ. We obtain where L c = ln(M 2 H /m 2 c ) − 1. To determine full NLO QCD correction to cc → H cross section, we need to include virtual corrections. They are computed in a standard way (see e.g. Ref. [26] where such a computation is reported); the result is then expanded around m c = 0. We renormalize the Yukawa coupling in the MS scheme at the scale µ. The result reads Upon combining virtual, soft, collinear and fully-regulated terms, we obtain the following NLO QCD contribution to cc → H + X cross section where The first term on the right hand side of Eq. (2.14) is the hard inelastic contribution; it can be computed directly in the massless limit, m c = 0. The second term on the right hand side of Eq. (2.14) is the soft-virtual piece; it describes kinematic configuration that is equivalent to the leading-order one. The third term in Eq. (2.14) describes kinematic configurations that are boosted relative to the leading-order ones; we note that the residual logarithmic dependence on m c is present in these contributions only.
A similar computation for massless charm partons requires, in addition, a collinear PDF renormalization to make the partonic cross section collinear-finite and independent of the regularization parameter . The result reads The partonic cross sections in Eqs. (2.14) and (2.16) should be convoluted with different parton distribution functions to obtain hadronic cross sections: in case of Eq. (2.16) we must use conventional MS PDFs whereas in case when the incoming charm quarks are massive a special set of PDFs is required. Nevertheless, since m c is just a collinear regulator, results for short-distance hadronic cross sections should be the same, independent of whether one starts with nearly massive or massless charm quarks. This requirement allows us to derive a relation between the "massive" and the MS PDFs. It reads The computation that we just described allows us to determine the coefficient G cc (z). We find (2.20) We can also compute "off-diagonal"coefficients G ab that involve charm quark and gluon PDFs; they are important for removing mass singularities that arise in g → c and c → g transitions. Computations proceed along the same lines as described above except that we employ other partonic processes for the analysis. Namely, we derive a relation for g → c transition by considering a cg → Hc process in a theory where only Yukawa coupling is present and no direct ggH coupling is allowed. To derive a relation for c → g transition, we again consider cg → Hc process but now only allow for the ggH coupling. In both cases only (quasi)-collinear singularities are present; this simplifies the required computations significantly. We find The results for the functions G ab (z) reported in Eqs. (2.20,2.21) are important for the calculation of NLO QCD corrections to the interference contribution to Higgs boson production in association with a charm jet. Indeed, as explained in the Introduction, to access the interference, we need to start with the massive incoming charm quarks and carefully study the massless limit. Since the charm mass serves as a collinear regulator, we are forced to use parton distribution functions f (m) i . We then use the relation Eq. (2.18) to express these functions through the conventional MS PDFs and, in doing so, remove logarithms of m c that are associated with the radiation by the incoming charm quarks. Because the interference contribution to pp → H + jet c involves a helicity flip, standard collinear logarithms associated with initial state emissions are not the only logarithms of the charm mass that appear in the cross section. We elaborate on this statement in the next section.

Interference contributions and factorization in the quasi-collinear limit
The discussion in the previous section shows that the dependence on m c disappears from hard cross sections provided that conventional factorization formulas for soft and quasi-collinear singularities, Eqs. (2.3,2.6), hold true. However, since the interference contribution requires a helicity flip, its soft and quasi-collinear limits are different from the conventional ones. As we explain below, such limits can still be described by simpler matrix elements but these matrix elements do not always correspond to processes with reduced multiplicities of final state particles.
To discuss and illustrate these subtleties in more detail, consider the process The first point that needs to be emphasized is that, if we work with a finite charm mass, true soft and collinear limits of the process in Eq. (3.1) are, in fact, conventional. These contributions can be extracted and combined with the virtual corrections to cg → Hc and renormalized gluon parton distribution function giving a finite result for the partonic cross section. Such a procedure is identical to what is usually done in NLO QCD computations [24,25,27] that are traditionally performed using dimensional regularization for soft and collinear divergences. However, cancellation of "true" infrared and collinear divergences does not tell us anything about non-analytic dependence of partonic cross sections on the charm mass that we need to extract before taking the m c → 0 limit. In this paper, we adopt a pragmatic approach and extract all the terms that are singular in the m c → 0 limit by studying interference contributions to squares of scattering amplitudes, computed explicitly with massive charm quarks, for all relevant partonic processes, cg → Hcg , gg → Hcc , cq → Hcq , cc → Hcc , cc → Hcc . In that sense, we do not attempt to develop an understanding of infrared and quasi-collinear factorization for generic processes that involve a helicity flip. However, to illustrate main differences with the conventional collinear factorization we discuss an emission of a collinear gluon off an incoming charm quark in case of the interference contribution in some detail.
To this end, we consider the process in Eq. (3.1) in the quasi-collinear limit p 1 · p 4 ∼ m 2 c . To describe this limit, we divide the amplitude for the full process cg → Hcg into two parts where M sing refers to diagrams where the gluon is emitted off the incoming quark with the momentum p 1 and M fin refers to the remaining diagrams. The first contribution (M sing ) is singular in the p 1 ·p 4 ∼ m 2 c → 0 limit whereas the second one (M fin ) is not. Upon squaring the amplitude, Eq. (3.2), and summing over polarizations of initial and final state particles, we obtain where the ellipsis stands for contributions that are finite in the p 1 ·p 4 ∼ m 2 c → 0 limit. We note that the product of singular and non-singular contributions to the amplitude squared that we retain in Eq. (3.3) is known to be non-singular in the conventional quasi-collinear limits provided that physical polarizations are used to describe the emitted gluon. As we will show below, this is not the case for the interference contributions considered in this paper.
To extract the quasi-collinear singularities from the square of the amplitude in Eq. (3.3) we need to analyze the quasi-collinear kinematics; for this analysis there is no difference between helicity-conserving and helicity-flipping contributions. Indeed, following the standard approach, we re-write the four-momentum of the incoming charm quark and the fourmomentum of the emitted gluon through massless momentap 1 and p 2 and find In the above equation, s = 2p 1 · p 2 and p 4,⊥ ·p 1 = p 4,⊥ · p 2 = 0. We use the on-shell condition p 2 4 = 0 and obtain It follows that We conclude that the kinematic region where p 2 4,⊥ ∼ m 2 c provides unsuppressed contributions to the cross section in the m c → 0 limit.
We write the singular contribution as follows and use the decomposition of the four-momenta p 1,4 given in Eq. (3.4) to obtain We note that in deriving Eq. (3.8) we have used p 4 · 4 = 0 and the gauge fixing condition The result for the singular contribution, Eq. (3.8), is generic; it does not distinguish between helicity-conserving and helicity-flipping contributions. However, it is easy to see that there is a difference between the two. For example, since the helicity-flipping contribution requires one additional power of m c , one can convince oneself that a combination of terms labeled as κ in Eq. (3.9) may contribute to the collinear limit of the interference but it cannot contribute to the collinear limit of the helicity-conserving amplitudes.
Hence, performing standard manipulations and paying attention to subtleties indicated above, we obtain the contribution to the interference that is non-analytic in the m c → 0 limit. It reads They have to be computed in the m c = 0 limit. 4 It is instructive to discuss the origin of the different terms in Eq. (3.10). The first term on the right-hand side of Eq. (3.10) contains the leading-order interference multiplied with the standard massive collinear splitting function; this is the conventional quasi-collinear limit applied to the interference. If only this term were present in Eq. (3.10), there would be no differences in the collinear factorization between helicity-changing and helicity-conserving contributions.
The second and the third terms on the right-hand side of Eq. (3.12) The interference requires a helicity flip that is facilitated by a single mass insertion. The above equation shows that this mass insertion can occur either in the (p 1 + m c ) density matrix of the external quark or "inside" theÂÂ † term. The structure that appears in the second term on the right hand side of Eq. (3.10) originates from the mass term in the density matrix. Once the mass term is extracted, the rest can be computed in the massless approximation. We find (3.13) The last term on the right hand side in Eq. (3.10) describes a quasi-collinear singularity that originates from the interference of singular and regular contributions in Eq. (3.3); it is particular to the helicity-flipping case and does not appear in the conventional collinear limits. As a consequence, this contribution still depends on the part of the reduced matrix element of the original 2 → 3 process calculated in the strict collinear limit of the incoming massless charm quark and the emitted gluon.
We emphasize once again that Eq. (3.10) shows clear differences between conventional factorization of quasi-collinear singularities and the factorization in case of the helicity flip. These differences lead to a peculiar structure of the logarithms of the charm quark mass in the interference contribution as they do not follow canonical pattern and cannot be removed by a transition to MS parton distribution functions. In addition, we also find that the interference contributions exhibit quasi-soft quark singularities that also lead to logarithms of the charm mass. Although it would be interesting to understand factorization of mass singularities in processes with the helicity flip from a more general perspective, our strategy for now is to explicitly compute all relevant contributions within fixed-order perturbation theory extracting all non-analytic m c -dependent terms along the way. Additional details of our approach can be found in several appendices.

Technical details of the calculation
In this section we briefly describe calculation of one-loop and real emission contributions to the interference part of the cg → Hc process. We begin with the discussion of the virtual corrections.
We compute the relevant one-loop diagrams keeping charm-quark masses finite. We employ the standard Passarino-Veltman reduction [28] to express the cg → Hc amplitude in terms of one-loop scalar integrals. 5 After computing the one-loop contribution to the interference, we expand the expression around m c = 0 and keep the leading O(m c ) term in this expansion. 6 The one-loop amplitudes contain ultraviolet and infrared singularities. The ultraviolet singularities are removed by the renormalization. We closely follow the discussion in Appendix A of Ref. [26] where many of the required renormalization constants are presented. Similar to Ref. [26], we renormalize the charm-quark mass in the on-shell scheme but employ the MS renormalization for the Yukawa coupling constant. In addition to the discussion in that reference, we require the one-loop renormalization constant of the effective ggH vertex that we take from Ref. [29]. After the ultraviolet renormalization is performed, the cg → Hc amplitude still contains 1/ poles of infrared origin. These poles satisfy the Catani's formula [30] and cancel with similar poles in real emission contributions to the partonic cross section.
According to the discussion in Section 3, factorization of quasi-collinear and quasi-soft singularities in the interference contribution is not canonical. This implies that even if we take the m c → 0 limit and switch to MS parton distribution functions, the NLO QCD corrections to the interference still contain logarithms of the charm mass. Since it is currently unknown how these logarithms can be resummed, we follow a pragmatic approach. Namely, we compute relevant virtual and real emission contributions, extract from them logarithms of the charm mass and take the m c → 0 limit once the mass logarithms have been extracted. To accomplish this, we construct subtraction terms for the real emission contributions for soft, collinear, quasi-collinear and quasi-soft singularities by direct inspection of the relevant matrix elements. The subtraction terms are then integrated over unresolved real emission phase space and combined with the PDF renormalization, including the transition from "massive" to MS PDFs, and the virtual corrections. The only difference with respect to the canonical procedure for NLO QCD computations is that in our case the subtraction terms are directly obtained from the squared amplitude and are not written in terms of easily recognizable universal functions; see Section 3 and Appendix A for further details. In particular, even contributions that are enhanced by two powers of a logarithm of the charm mass, O(α s ln 2 (m c )), do not appear to be proportional to the leading-order interference contribution to the cross section.

Numerical results
To present numerical results we consider proton-proton collisions at 13 TeV. We take M H = 125 GeV for the Higgs-boson mass and m c = 1.3 GeV for the pole mass of the charm quark. The charm Yukawa coupling is calculated using the MS charm mass, m c (M H ) = 0.81 GeV. 7 We use NNPDF31 lo as 0118 and NNPDF31 nlo as 0118 parton distribution functions [34,35] for leading and next-to-leading order computations, respectively. The value of the strong coupling constant α s is calculated using dedicated routines provided with NNPDF sets.
To define jets we use standard anti-k ⊥ algorithm with ∆R = 0.4; charm jets are required to contain at least one c orc quark. For numerical computations we require at least one charm jet with p t,j > 20 GeV and |η j | < 2.5. Moreover, we demand that the charm parton inside the charm jet carries at least 75% of the jet's transverse momentum. 8 The latter requirement removes kinematic cases where a soft charm is clustered together with a hard gluon into a charm jet in spite of a large angular separation between the two. Since, as we explained earlier, our calculation is logarithmically sensitive to soft emissions of charm quarks, defining charm jets with an additional cut on the charm quark transverse momentum allows us to avoid jet-algorithm dependent logarithms of m c that may appear otherwise. We note that we apply all these requirements even in the subtraction terms where c andc momenta are computed in the collinear and/or soft approximations.
We start by presenting fiducial cross sections for the three terms in Eq. for the ggH-dependent cross section, the Yukawa-dependent cross section and the interference, respectively. In Figure 2 we show the comparison between µ = M H cross sections and the interference in dependence on the cut of the charm jet transverse momentum. We observe that the ratio of the interference to the Yukawa-dependent contribution is about ten percent for all values of the p t,j -cut.     It is instructive to separate the NLO contributions to the interference into parts that are independent of m c and parts that are logarithmically enhanced for all the partonic channels. The relevant results are shown in Table 1. We find that the largest contribution at NLO comes from the gluon-gluon channel which is enhanced by the large gluon luminosity. Also, the charm-gluon (cg) and the charm-quark channels (cq) provide relatively large contributions. 9 Note that the (cq) channel is free of logarithmic contributions since there are no singular limits that involve charm quarks. Contributions related to the PDF transformation do not feature the double-logarithmic part since the O(ln 2 m c ) terms originate exclusively from soft-collinear limits that involve c-quarks.
It follows from Table 1 that double-logarithmic terms and single-logarithmic terms provide nearly equal, but opposite in sign, contributions to the NLO QCD interference. This cancellation between terms with different parametric dependence on m c should be considered as an artifact but it does emphasize that studying only the leading logarithmic O(ln 2 m c ) contribution in this case is insufficient for phenomenology. We also note that the O(ln 2 m c ) term in the cg channel is quite small reflecting the fact that there is a very strong -but incomplete -cancellation between double-logarithmic contributions to real and virtual corrections in this case. Finally, we emphasize that it is unclear to what extent these various cancellations persist in higher orders; for this reason, a resummation of charm-mass logarithms for the interference contribution is desirable.
We continue with the discussion of kinematic distributions. We focus on the transverse momentum and the rapidity distributions of Higgs bosons in the interference contribution to pp → Hc cross section. They are shown in Figure 3. We first discuss the transverse momentum distribution, Figure 3a; when interpreting this figure it is important to recall that the absolute value of both LO and NLO distributions is plotted there and that the LO distribution is always negative. We observe in Figure 3a that the leading-order distribution (blue) is large and negative at small p t,H ; as p t,H increases, the distribution goes to zero. The NLO QCD corrections affect the shape of the p t,H distribution. Indeed, a sharp edge at p t,H = 20 GeV, present at leading order, gets smeared at NLO. At moderate values of transverse momenta p t,H ∼ 60 GeV the K-factor is equal to one, while there is a large O(+50%) reduction 10 at p t,H ∼ 100 GeV. Second, at p t,H ∼ 150 GeV the NLO distribution goes through zero and becomes positive for larger values of p t,H . Asymptotically, at even higher p t,H the LO and NLO distributions appear to be equal in absolute values but opposite in sign. Of course, all this happens at such high values of p t,H that are irrelevant for phenomenology, but it is quite a peculiar feature nevertheless.  We only consider the interference contribution. We note that the absolute value of dσ Int /dp t,H is displayed in the left panel. This implies that this distribution actually changes sign at around p t,H ∼ 150 GeV. The lower panels show ratios to the LO interference contribution.
Compared to Higgs transverse momentum distribution, the rapidity distribution of the Higgs boson in the interference contribution is much less volatile. Indeed, it follows from Figure 3b that the difference between leading and next-to-leading-order distributions is welldescribed by a constant K-factor all the way up to |y H | ∼ 2. Beyond this value of the rapidity, the NLO distribution goes to zero faster than the LO one.

Conclusions
Production of Higgs bosons in association with charm jets at the LHC is mediated by two distinct mechanisms, one that involves the charm Yukawa coupling and the other one that involves an effective ggH vertex. Their interference involves a helicity flip and, for this reason, vanishes in the limit of massless charm quarks.
Since partonic cross sections are routinely computed for massless incoming partons and since the charm quark appears in the initial state in the main process cg → Hc, it is interesting to understand how to circumvent the problem of having to deal with a massive parton in the initial state and to provide reliable estimate of the interference contribution.
We have addressed this problem by studying the m c → 0 limit of the helicity-flipping interference contribution including NLO QCD corrections. We have shown that the factorization of quasi-collinear and quasi-soft singularities in this case differs from the canonical pattern. We used explicit expressions for real and virtual matrix elements to extract logarithms of the charm quark mass and, having accomplished this, took the m c → 0 limit in the remaining parts of the computation. We removed parts of the O(ln m c ) contributions by expressing results through conventional MS parton distribution functions valid for massless partons. Nevertheless, given an unconventional behavior of the interference in quasi-soft and quasi-collinear limits, logarithms of the charm quark mass survive in the final result for the NLO QCD corrections.
We have found that the absolute value of the leading-order interference is reduced by about fifty percent once NLO QCD corrections are accounted for. This significant but still "perturbatively acceptable" reduction is the result of a very strong cancellation between terms that involve double and single logarithms of the charm quark mass. We have observed that the NLO QCD corrections to the interference are kinematics-dependent and may change shapes of certain kinematic distributions in a significant way.
Higgs boson production in association with a charm jet is a promising way to study charm Yukawa coupling at the LHC [12]. The interference contribution, that is estimated to be about 10 percent of the Yukawa contribution at leading order, could have been perturbatively unstable given the required helicity flip and an unconventional pattern of quasi-soft and quasi-collinear limits. We addressed this question by performing a dedicated NLO QCD computation for the interference term and did not find a strong indication that this might be the case. Nevertheless, the moderate size of the NLO QCD corrections is the consequence of a very strong cancellations between double and single logarithms of the charm mass. It is unclear if this cancellation persists in higher orders. Hence, resummation of O(ln m c )enhanced terms for this process is quite desirable. In this appendix we describe a procedure to extract O(ln m c ) contributions to real emission corrections. They arise because of the quasi-singular behavior of real emission amplitudes in the soft or collinear limits involving charm quarks. The potential singularities in these limits are regulated by the charm mass leading to an appearance of O(ln m c )-enhanced terms when integrated over relevant phase spaces. To extract logarithms of m c , we subtract approximate expressions from exact matrix elements that make the difference integrable in the m c → 0 limit and integrate the subtraction terms over unresolved phase space to explicitly extract logarithms of m c .
As an example, we consider the gluon-gluon partonic channel, i.e.
and discuss the extraction of O(ln m c ) terms in detail. This channel is suitable for such a discussion since, if the charm-quark mass is kept finite, it is free of soft and collinear divergences. Hence, all relevant contributions can be computed numerically for small but finite m c , and used to validate formulas where logarithms of m c have been extracted and m c → 0 limit has been taken where appropriate. As we already mentioned in the main text, we use the nested soft-collinear subtraction scheme which, at this order, is equivalent to the FKS subtraction scheme [24,25]. The details of the subtraction scheme can be found in the literature [20][21][22][23] and we do not repeat this discussion here. Nevertheless, the treatment of quasi-collinear and quasi-soft singularities related to the emission of massive charm quarks is new and requires an explanation.
We focus on the interference contribution between the Yukawa-like and the ggH-like production mechanisms in the process Eq. (A.1). The interference term is non-zero only if helicity flip on the charm line occurs. Furthermore, the presence of such a helicity flip causes the usual factorization formulas to break down and the singular limits need to be explicitly analyzed. We note that, thanks to the symmetry of the squared amplitude for the process in Eq. (A.1) under the exchange of c andc, we can consider only the case wherec quark becomes soft or quasi-collinear to one of the other partons. The case when both c andc become unresolved is impossible since we require a charm jet in the final state.
The quasi-singular limits which appear in this channel are related to the soft-quark limit S 4 with E 4 ∼ m c and the three collinear limits C 4i with i = 1, 2, 3 where (p 4 ·p i ) ∼ m 2 c . Performing an iterative subtraction of these singular limits, we find where the first term on the right-hand side denotes the fully-regulated contribution and the second and third terms are the collinear and the soft integrated subtraction terms. The factors ω

A.1 Integration of the soft-quark subtraction terms
Consider the soft-quark subtraction term To compute it, we need to know the behavior of the amplitude in the limit p 4 ∼ m c → 0 and then integrate it over the phase space of the charm anti-quark with momentum p 4 . Although, normally, soft (gluon) emissions factorize into a product of an eikonal factor and a tree-level matrix element squared, a similar formula for soft-quark emission, relevant for helicity-flipping processes, does not exist. Hence, we determine the soft-quark limit of the interference by studying an explicit expression of the amplitude for the process in Eq. (A.1) in the limit p 4 ∼ m c → 0. We find where functions F ij depend on the momenta p 1 , p 2 and p 3 only. We emphasize that these functions are different from the leading-order interference contribution. The massless limit, m c → 0, can be now taken everywhere except for the eikonal factors and the phase-space measure of the unresolved parton [dp 4 ]. We write the integrated soft subtraction term as follows where . . . denotes the phase space integration and the relevant soft integrals can be found in Appendix B. We stress that soft integrals are finite in four dimensions since they are naturally regulated by the charm-quark mass m c .

A.2 Integration of the quasi-collinear subtraction terms
In this subsection, we discuss how to define and compute the soft-subtracted quasi-collinear limits of the interference contribution using the process in Eq. (A.1) as an example. We focus on the sector 43 where c andc become collinear to each other. The relevant quantity reads 11 To proceed further, we split the above formula into collinear and soft-collinear terms We first discuss the collinear subtraction term C 43 F LM (1 g , 2 g ; 3 c , 4c) defined as follows In the above equation, the functions C 1,2 depend on the momenta p 1 , p 2 , p 3 andp 4 . The bar over momentum p 4 indicates that the relevant collinear limit has been taken, i.e.
Note that in Eq. (A.10) β 4 = 1 − m 2 c /E 2 4 is the velocity ofc and n 3 is a unit vector pointing in the direction of momentum p 3 . We note that in Eq. (A.9) the massless limit m c → 0 has not been taken. We also note that the functions C 1,2 are regular in the soft-quark limit, Our goal is to extract all O(ln m c ) terms arising from Eq. (A.9) and take the massless limit after that. To do so, we add and subtract the soft limits of the functions C i where C i,soft (1 c , 2 g , 3 c ) = C i (1 g , 2 g , 3 c , 0). The above procedure splits the integral in Eq. (A.9) into two parts: the regulated integral that contains the expression in the square bracket in Eq. (A.11) and the soft part. In the regulated part, the soft divergence at E 4 = 0 has been regularized. This implies that, after integrating 1/(p 3 + p 4 ) 2 over the relative angle between p 3 and p 4 and extracting logarithms of m c from this angular integral, we can set m c to zero everywhere else right away. We obtain where we have used the fact that in m c → 0 limit we can write p 3 = zp 34 andp 4 = (1 − z)p 34 for p 2 34 = 0.
We will now discuss the soft part of the collinear subtraction term. It reads [dp H ][dp 3 ][dp 4 ](2π) 4 δ (p 12 − p H − p 3 −p 4 ) (A. 13) We emphasize that this term still contains soft singularity and, for this reason, the m c → 0 limit cannot be taken. However, it is convenient to combine this integral with the soft-collinear subtraction term C 43 S 4 F LM (1 g , 2 g ; 3 c , 4c) , c.f. Eq. (A.8); if this is done, the required computations simplify significantly. The soft-collinear integrated subtraction term in sector 43 reads (A.14) To derive this result we used the soft-limit of the interference amplitude reported in Eq. (A.5).
We emphasize that, since the soft operator is present on the left hand side in the above equation, the soft-quark momentum p 4 is removed from the energy-momentum conserving delta-function. Moreover, since the two integrals in Eqs. (A.13,A.14) appear to be the same up to the argument of the delta-functions. We combine the two integrals and find [dp H ][dp 3 ][dp 4 ](2π) 4 δ ( To proceed further, we note that it is straightforward to integrate over directions of the quark with momentum p 4 but integration over its energy is more involved. It is convenient to split the E 4 integration into two regions by introducing an auxiliary parameter σ We choose σ to satisfy the following inequality m c σ E 3 . For small energies, E 4 < σ E 3 , we can drop the momentump 4 from the energy momentum conserving delta-function which leads to This relation implies that the integrand in Eq. (A.16) is non-vanishing only in the high-energy domain where E 4 > σ m c and, therefore, the limit m c → 0 can be taken. This leads to the following expression To arrive at Eq. (A.19) we introduced the four-momentum p 34 = p 3 +p 4 and a variable z such that p 3 = zp 34 in terms that contain the delta-function δ(p 12 − p H − p 3 −p 4 ). In terms that contain the delta-function δ(p 12 − p H − p 3 ), we set (1 − z) = E 4 /E 3 , rename p 3 into p 34 and set σ → 0. The lower integration boundary z min is given by z min = 1 − E max /E 34 < 0. In total, the integrated collinear term C 43 (1 − S 4 )F LM (1 g , 2 g ; 3 c , 4c) is given by a sum of expressions in Eqs. (A.12,A.19). We describe a numerical check of validity of this result in the following section.

A.3 Numerical checks
Since the cross section of the gluon-gluon channel, Eq. (A.1), is finite as long as we keep the non-zero charm mass, analytical results derived in the previous section can be checked numerically by computing σ gg→Hcc explicitly for small values of the charm mass without any approximation.
The comparison is shown in Figure 4. We use fiducial cuts described in the main text and compare hadronic contributions to the interference for gg partonic channel computed in two different ways. Green points (rectangles) show the results of the computation without any approximation, i.e. by directly integrating the matrix element squared. Blue points (circles) show the results of the computation that relies on the expansion around m c → 0 limit, as described in previous subsections. The two results should agree for small values of m c . The upper panel of Figure 4 shows the absolute values of the interference cross section in the gg partonic channel obtained with the two methods, while their difference is shown in the lower panel. We see a better and better agreement between the two results as we mover to smaller and smaller values of the charm-quark mass. This indicates that the m c -dependence of the interference contribution is properly reconstructed.

Cross-section [fb]
σ rec σ real Figure 4: The cross section of the gg → Hcc process calculated by a direct integration of the matrix element with non-zero charm-quark mass, σ real (green rectangles), and reconstructed using procedure described in previous subsections, σ rec (blue circles). We employ the same parameters and kinematic constraints as in the main text.

B Soft-quark integrals
In this section we list integrals that are required for the integrated soft-quark subtraction terms. We need a number of integrals depending on the type and configuration of the emitters p a and p b as well as the propagator appearing in the eikonal factor.
We note that we are interested only in the terms that contain logarithms of the charmquark mass and constant terms, but we drop all power-suppressed terms which vanish in the m c → 0 limit. All integrals are computed in d = 4 dimensions since all singularities are naturally regulated by the charm-quark mass.
The phase-space measure for a massive-quark emission, p 2 4 = m 2 c , reads where k 4 is the length of p 4 momentum, dΩ (3) denotes angular integration and E max is the usual energy cutoff of the nested soft-collinear subtraction scheme. In the remaining part of this section, we list soft-quark integrals that are needed to obtain integrated soft-quark subtraction terms, see Section A.1 for details. 12 Two massless emitters: Two emitters a, b have four-momenta p a = E a (1, n a ) and p b = E b (1, n b ), respectively. Both four-momenta are light-like p 2 a = p 2 b = 0. Vectors n a and n b describe direction of flight of the emitters; we refer to the opening angle between n a and n b as θ ab .
One massive and one massless emitters: Two emitters a, b have four-momenta p a = E a (1, n a ) and p b = E b (1, β b n b ), respectively. They satisfy p 2 a = 0 and p 2 b = m 2 c . We refer to the opening angle between n a and n b as θ ab . We require three soft integrals of this type where we used s ab = sin(θ ab /2) and c ab = cos(θ ab /2).
Two massive emitters: Two emitters a, b have four-momenta p a = E a (1, β a n a ) and p b = E b (1, β b n b ), respectively. They satisfy p 2 a = p 2 b = m 2 c . We refer to the opening angle between n a and n b as θ ab . In this case, we use E max = E a . We find  where we used s ab = sin(θ ab /2) and c ab = cos(θ ab /2).