Renormalization of the next-to-leading-power γγ → h and gg → h soft quark functions

,


Introduction
Soft functions -vacuum expectation values of soft fields -appear in factorization theorems for QCD high-energy scattering whenever soft effects do not simply cancel.In leading power in the high-energy expansion, the soft fields appear through semiinfinite light-like Wilson lines defined as where P is the path-ordering, g s denotes the strong coupling, and T a the SU (N c ) generator in the representation of the energetic particle interacting with the soft gluon.This is a consequence of the fact that at leading power soft gluons couple only through the eikonal vertex ig s T a n µ , where n µ is the direction of the four-momentum of the particle.A well-known example is the soft operator for the Drell-Yan process q q → γ * near the kinematic threshold [1].T denotes the time-ordering.A similar soft operator weighted by the final-state measurement function is relevant to the two-jet limit of event shapes in e + e − annihilation.
In soft-collinear effective theory (SCET), the soft Wilson lines arise from the field redefinition that removes the eikonal interaction from the leading-power softcollinear Lagrangian [2].This transformation does not remove the soft-gluon interactions from the sub-leading power Lagrangian [3,4], resulting in soft field-strength insertions Y † n F µν Y n in the soft operators.Generalized soft-operators of this type appeared first in the context of heavy-quark physics [5][6][7] and in the context of highenergy scattering in [8,9].In addition, there exist sub-leading power soft operators that contain soft light quarks, rather than gluon field-strength insertions.These soft functions are key ingredients in factorization formulas at next-to-leading power (NLP), but relatively little is known about them beyond their definition and leadingorder expressions from direct computation.An exception is the Drell-Yan threshold for which both the soft-gluon and soft-quark NLP soft function have been computed at O(α 2 s ) [10][11][12].Their anomalous dimension is however only known at the level of the leading double pole from direct computation, and only for the case of the soft gluon field insertion [9,13].Developing methods to compute the anomalous dimension (renormalization-group kernels) of generalized soft functions is a prerequisite for extending resummation at NLP to the next-to-leading logarithmic accuracy.
In this work, we focus on the NLP soft quark operators O γ and O g that appear in the factorization formula for the subleading-power γγ → h and gg → h form factors mediated by the bottom quark [14][15][16].In the context of the factorization formula, the leading-order renormalization group (RG) kernels for both operators have been inferred from RG consistency [16,17], that is, from the requirement that the full process is factorization scale independent and the prior knowledge of the RG kernels and anomalous dimensions of all other functions in the factorization formula.For the case of O γ , the consistency has been confirmed by a direct calculation of the kernel [18], but a corresponding result for the non-abelian generalization O g is not available.In the following we show that the RG kernels for both soft operators can be obtained with position-space renormalization methods by a relatively straightforward computation.The case of O γ relevant to γγ → h turns out to be very similar to the renormalization of the leading-twist heavy-light quark light-ray operator [19,20] that defines the so-called B-meson LCDA.The gg → h soft operator O g presents a complication, related to "colour charges at infinity".In this respect, it resembles the QED-generalization of the B-meson LCDA for decays into two electrically charged light particles [21,22].
Recent years have witnessed rapid developments in the application of conformal symmetry to perturbative QCD calculations (see [23] for an early pedagogical introduction and the recent treatment [24] with technical details in [25]).The basic idea is to consider QCD at the critical point in non-integer dimensions where the β function vanishes identically, allowing one to take advantage of the conformal symmetry to the fullest extent.In this work, we extend the use of conformal symmetry, for the first time, to the study of non-string operators O γ and O g , where the quark fields are located on different light rays.
The outline of this paper is as follows: We set up the definitions of the soft operators in Sec. 2, then calculate the RG kernel for O γ in Sec. 3. The result coincides with the one from [18], but our derivation is considerably simpler.The renormalization of the soft operator O g for the soft-quark contribution to the gg → h form factor is presented in Sec. 4. The naive definition of the operator leads to a rapidity-divergent anomalous dimension, which requires an extra subtraction.As a consequence, the factorization formula presented in [16] needs a rearrangement between its soft and (anti-) collinear parts that resembles the factorization of electromagnetic effects in two-body B-meson decays to electrically charged mesons [21,26].Regarding the renormalization, the soft-quark operator for γγ → h is essentially the "double-copy" of the leading-twist B-LCDA operator.We exploit this observation in Sec. 5 to infer its two-loop RG kernel from the known result for the B-LCDA [27] using the conformal symmetry technique.We construct an ansatz for the kernel of O g at higher orders and test it against the consistency requirement at two-loop.We conclude in Sec. 6.Three appendices provide further technical details.In particular, in App.B we demonstrate that the collinear matching coefficients, which enter the subleading-power factorization formulas, remain unaltered and gauge-independent in the presence of an off-shell infrared regulator in general covariant R ξ gauge.

Setup
The soft-quark functions appearing in the Higgs production through a light-quark loop in two photon (γγ → h) and two gluon (gg → h) fusion at the subleading power are defined as vacuum matrix elements of operators made of semi-infinite, light-like Wilson lines (1.1) and soft-quark fields.The bare operator definitions in position space read, for γγ → h and gg → h, respectively, where n µ ± are two light-like vectors, n 2 ± = 0, normalized to n − • n + = 2. q denotes the light quark field.The superscript "uns" on O g (s, t) indicates that it is an unsubtracted operator corresponding to the definition in [16] 1 .A necessary modification will be explained in Sec. 4. Notice that the two adjoint colour indices a, b in O uns g are kept open, and we keep them implicit for convenience.Here, the light-like Wilson lines Y n ± given in (1.1) are in the fundamental representation of SU (N c ).While the operator O γ stays the same for the h → γγ decay process due to (2.3) below, the Wilson lines Y n ± in O uns g for h → gg decay extend from 0 to +∞n ± instead, as required for the coupling of soft gluons to out-going (anti-) collinear gluon fields.The anomalous dimensions, however, remain the same under the exchange of initial and final state due to the time-reversal symmetry of QCD.
In the γγ → h case, the semi-infinite Wilson lines combine to form finite-length segments as in the case of light-cone string operators, with The case of O uns g (s, t) is fundamentally different due to two inserted colour matrices corresponding to colour charges at the infinity, which arise from the energetic back-to-back gluons in the initial state of the gg → h process.We can bring O uns g (s, t) into a form that resembles O γ (s, t) by applying the identity to convert the fundamental Wilson lines to adjoint ones, defined as where 0 + guarantees convergence of the integral and is related to the pole prescription for the in-going eikonal propagator.We further introduced the infrared (rapidity) regulator e −iλδ ± with δ i ∈ R in the definition of the semi-infinite adjoint Wilson lines [28] .In momentum space, this regulator is equivalent to modifying the eikonal propagator generated by the λ-integration as with ℓ denoting the soft momentum absorbed by the Wilson line.The need for this regulator will become clear in Sec.4.1.In the following, we often drop the +i0 + with the implicit understanding that δ ± means δ ± + i0 + .The δ ± regulator can be further related to the off-shellness of the energetic external gluons from which the semi-infinite Wilson lines are induced.More specifically, the eikonal propagators describing soft gluon absorption in the (anti-) collinear sector are constructed as follows: (2.8) Hence we identify the relations where p c (p c) denotes the (anti-) collinear external momentum.These relations will play an important role in Sec.4.3.It is then possible to rewrite O uns g (s, t) as The expression in the second line can be derived directly from the subleading power SCET Lagrangian insertion without applying the soft-decoupling transformation.In practice, we find it is more convenient to adopt this expression for calculating the renormalization factor of the operator O uns g (s, t).In phenomenological studies of Higgs decay or production through a bottomquark loop, the soft functions enter the factorization formulae as the vacuum matrix elements of the relevant operators in momentum space.The required Fourier transform reads,2 The fact that the momentum-space soft functions depend on the single variable w = ℓ − ℓ + is a consequence of reparameterization/Lorentz invariance."Disc" denotes the discontinuity of the Fourier transform, and ℓ ± ≡ n ± •ℓ > 0 are the soft quark lightcone momentum components conjugate to s and t, respectively.They are positive since their physical interpretation corresponds to twice the soft quark energy as the multipole expansion puts the soft quark position on the light cone.We anticipated that the Fourier transform for the subtracted soft gluon operator O g (s, t) must be used instead of the unsubtracted one, O uns g (s, t).The definition for O g (s, t) will be given later in (4.11).
To convert the anomalous dimension from position space to momentum space, it proves to be more convenient to perform the Fourier transform at the operator level using the expression [29] (2.12) Similar to the case of the B-meson LCDA, it is necessary to take s (t) as a variable in the lower (upper) half of the complex plane, hence the replacement s → s − i0 + (t → t + i0 + ) in the Fourier integral.This arises directly from ℓ ± > 0 dictating the inverse transformation to be3 which technically constitutes a Laplace transform.Then, the ∓i0 + prescription for s and t in (2.12) follows from the requirement of the existence of the inverse Laplace transform, which guarantees the invertibility between the soft operators in momentum space O γ,g (w) and position space O γ,g (s, t).In the following we implicitly assume that s (t) means s − i0 + (t + i0 + ).The renormalization of these soft operators corresponds to integrating out the quantum fluctuations above the soft scale µ (of order m b in phenomenological applications) on which the renormalized soft operators depend.This can be achieved in the background field method [31], which separates fields into (soft) background (with momenta p < µ) and quantum (with momenta p > µ) components, Integrating out the "quantum" modes in a gauge-invariant operator O leads to an effective operator O eff that only depends on the "soft background" degrees of freedom, where S is the QCD action.The advantage of the background field method is that it retains the gauge invariance of O eff in the evaluation of the path integral (2.15).Due background gauge (2.12) and the associated Feynman rules can be found in [32].Since (2.11) is invariant under local gauge transformations of soft background fields, the background technique allows us to choose gauges for the soft background and quantum fields independently.
3 Renormalization of O (s, t) The above figure is still not a single pdf that can be scaled.]Since Wilson lines are finite-length ones in O (s, t), only ultraviolaet (UV) singularities exist.The potential infrared singularities are naturally regularized by the non-vanishing o↵-shellness (or mass) of the soft quark.In this section, we provide details on renormalizing the soft operator directly in the position space, which turns out to be rather compact and easy to generalize to higher loops.The one-loop contributions to O (s, t) are depicted in Fig. 1.There are also mirror parts for diagrams (a) and (c).Diagram (a) and its mirror give rise to non-local parts. 2 , and diagram (b) gives rise to the local cusp terms, while (c) and (d) do not have UV poles.To simplify the notation, we also use q(p) and q(p) to denote the external quark and anti-quark fields in momentum space.We can put either the quark mass or the o↵shellness (but not both) to zero from the start since we focus on UV poles, and we choose to drop the quark mass.

Coordinate-space calculation
Inserting the QCD interaction Lagrangian once, expanding the n -Wilson line to first order and performing Wick contractions, diagram (a) reads (We drop the subscripts to the presence of the soft background fields, in the simplest case, the conventional gauge-fixing condition needs to be modified to where D µ [A s ] is the covariant derivative acting on the quantum component in the presence of soft background.To the lowest order in the soft background field A s , the gauge-fixing condition (2.16) coincides with the Feynman gauge.Eq. (2.16) serves as a simple demonstration of the background field method principle, but it is important to note that the gauge-fixing condition can be generalized to accommodate the general covariant (R ξ ) and axial (e.g., light-cone) gauges as well [31].This allows us, in practice, to adopt the R ξ gauge in our calculations for integrating out the quantum fields.The expression for the action S[q s + q q , A s + A q ] − S[q s , A s ] in the background gauge (2.16) and the associated Feynman rules can be found in [31].Since (2.15) is invariant under local gauge transformations of soft background fields, the background technique allows us to choose gauges for the soft background and quantum fields independently.
3 Renormalization of O γ (s, t) Since the Wilson lines in O γ (s, t) are finite in length, no rapidity divergences exist.While the potential infrared singularities are naturally regularized by the nonvanishing off-shellness (or mass) of the soft quark, the UV divergence is exclusively captured by dimensional regularization.In this section, we provide details on renormalizing the soft operator directly in the position space, which turns out to be rather compact and easy to generalize to higher loops.The one-loop contributions to O γ (s, t) are depicted in Fig. 1.There are also mirror ones for diagrams (a) and (c).Diagram (a) and its mirror give rise to non-local terms, 4 diagram (b) gives rise to a local cusp term, diagram (c) and its mirror contribute in background R ξ gauge only when ξ ̸ = 1, and diagram (d) does not have UV poles.To simplify the notation, we also use q(p) and q(p) to denote the external quark and anti-quark fields in momentum space.We can put either the quark mass or the off-shellness (but not both) to zero from the start since we focus on UV poles, and we choose to drop the quark mass.

Position-space calculation
Inserting the QCD interaction Lagrangian once, expanding the n − -Wilson line to first order and performing Wick contractions, diagram (a) reads (we drop the subscripts s and q for brevity, and abbreviate ū = 1 − u), (3.1)We employ dimensional regularization with D = 4 − 2ε for the UV divergences and the momentum-space MS scheme for renormalization, which amounts to subtracting poles in ε.The n + -Wilson line [0, sn + ] does not contribute to the integral I a hence is not written for brevity.The Wick contractions are between quantum fields and the integration over u comes from separating the Wilson line into soft background and quantum modes [32], The two external quark fields q(z) and q(sn + ) in the first line of (3.1) are soft, and we employed the light-cone gauge n − • A s (x) = 0 for the soft background fields.This sets the soft background Wilson line in the n µ − direction to 1 and drastically simplifies the calculation.As explained earlier, we effectively employ the Feynman gauge in detailing our calculations for integrating out the quantum fields (loop integrals) and comment on the results in general R ξ gauge.After some algebra we obtain where only 1/ε divergent terms are kept after the second equal sign.The pole arises from the integral over the loop momentum l, while integrating over p facilitates Fourier transform back into position space.The non-vanishing off-shellness p 2 of the external quark provides the necessary IR regulator.In the last step, we elevated the result to operator form without presuming any gauge fixing condition on the soft background fields.This is justified by the gauge invariance of the computation with the background field method.Hereafter, we adopt the standard plus-prescription defined as Hereafter, ln(stµ 2 e 2γ E ) ≡ ln(µ(s − i0 + ) e γ E ) + ln(µ(t + i0 + ) e γ E ) serves as a shorthand notation without ambiguities regarding the branch cuts of the logarithms.Unlike diagram (a), Wilson lines of different directions participate in the interaction, making it impossible to remove the soft background Wilson lines completely by a gauge fixing.Therefore, technically, one has to first expand each Wilson line using (2.14), then Wick contract between terms with a single quantum field from each direction n µ ± in the expansion using the background-field Feynman rules, and finally resum the infinite series.We can, however, bypass all these complications by focusing on computing the corrections to A s at the lowest order (namely, setting A s = 0 in the Wilson lines) knowing that the soft background Wilson lines ensuring gauge invariance will fall in place automatically from higher-order A s contributions in the background field technique due to its manifest gauge invariance.We apply this trick in calculating diagrams hereafter.The expression in the first line of (3.5) comes from Wick-contracting the two quantum gluon fields, n − •A q (utn − ) and n + •A q (vsn + ) from the Wilson line expansion (3.2) setting A µ s = 0.The result is then matched to the full gauge-invariant operator O γ yielding the last line of (3.5).In general R ξ gauge, the result for diagram (b) reads, By setting the background gluon field A s directly to 0, diagram (c) in Feynman gauge gives rise to the loop integral where K n (z) denotes the Bessel-K function.While the virtuality p 2 ̸ = 0 regulates the IR behaviour of the integral, the UV is tamed by the off-light-cone distance between q(tn − ) and [0, sn + ] entering the phase factor.We use the asymptotic behaviour to check that the remaining integrals over u and α do not produce singularities from the endpoints u → 0, α → 0, 1.This argument allows us to conclude that diagram (c) and its mirror do not contribute to the renormalization of O γ (s, t) in Feynman gauge.However, diagram (c) and its mirror diagram do contribute in the general R ξ gauge, which will become important in Sec. 4. The R ξ gauge result reads, and the mirror diagram contributes identically.Finally, the one-loop integral for diagram (d) in Feynman gauge reads, (3.10) Similar to diagram (c), the IR-finiteness of diagram (d) is a consequence of the nonvanishing off-shellness p 2 1 , p 2 2 ̸ = 0 of the two external quarks.The superficial (here) logarithmic UV divergence of the integrand in the second line is tamed by the offlight-cone phase factor e i[tn − •l−usn + •(p+l)] .Thus, diagram (d) does not contribute to the renormalization of O γ (s, t).In this case, this remains true in general R ξ gauge.
In addition, one must account for the MS quark field renormalization constant.In the general R ξ gauge, it reads Summing up all contributions in R ξ gauge, including the mirror diagrams and the quark field renormalization, and then multiplying an overall factor (−1), we find that the ξ-dependence cancels, resulting in the one-loop gauge-invariant counterterm/renormalization factor 12) The superscript "bare" indicates the bare operator (in terms of renormalized fields).
The position-space anomalous dimension / evolution kernel for O γ (s, t) at one-loop is then immediate: The operator O γ then satisfies the following renormalization group equation (RGE) in position space, We observe that the one-loop anomalous dimension in position space is rather simple.Most importantly, it essentially constitutes a double-copy of the evolution kernel of the leading-twist B-meson LCDA [19,20], and we will investigate this observation further in Sec.5.1.On the one hand, our position space calculational method is straightforwardly applicable to other operators and higher loops.On the other hand, position-space formalism allows us to take advantage of conformal symmetry techniques, which are useful not only for solving the renormalization-group equations (RGEs) but also for performing higher-order calculations.

Anomalous dimension in momentum space
In phenomenological studies, the momentum-space representation of the renormalization factor/anomalous dimension is often preferred.We therefore apply the Fourier transform to (3.12), which can be accomplished straightforwardly by employing where Õ(ω ′ ) represents a momentum space "test function".To apply the above formulae, it is necessary to follow the pole prescription s → s − i0 + , t → t + i0 + established in (2.12).For the latter it is helpful to take complex conjugate of (3.15) and (3.16) to accommodate the +i0 + prescription.We adopt the standard plusprescription in the momentum space defined as, where f (ω) is a test function.Following the procedure outlined in appendix A for performing the Fourier transform yields the momentum-space renormalization factor Z γ which reads, where L w = ln w µ 2 and is the same distribution that appears in the anomalous dimension of the leading-twist B-meson LCDA in momentum space [19].
The anomalous dimension γ γ that controls the scale dependence of the vacuum expectation value of the soft operator O γ (2.11), S γ , through the renormalizationgroup equation (RGE) is then immediate from (3.18), This concludes the calculation of the one-loop renormalization and anomalous dimension of the soft-quark photon operator O γ .Eq. (3.21) coincides with the result previously obtained in [18], but has been derived from standard position-space methods without going through transverse-momentum unintegrated functions.The position-space formalism proves to be even more advantageous for calculating the anomalous dimension of the more complicated soft-quark gluon operator, as we will detail in the subsequent section.

Applies to all figures!]
There are additional one-loop corrections (shown in the second row of the figure) due to gluon exchanges with the semi-infinite adjoint Wilson lines.Here, we adopt the second representation given in (2.8).Due to the two semi-infinite Wilson lines, O uns g (s, t) exhibits rapidity divergences.Consequently, one needs to separate UV poles from the singular behaviour in rapidity by introducing additional regulators.In this work, we introduce an o↵-shellness ± as in (2.6), (2.7) to regulate the rapidity divergence induced by semi-infinite Wilson lines.
[• Q7 Change to ± to conform with definition (2.6).] [• Q8 the following is no longer true and should be removed: "This regulator appears as an exponential suppression when the adjoint Wilson line approaches infinity."]4.1 Unsubtracted operator O uns g (s, t) The total one-loop corrections generated by diagrams in the first row of Fig. 2 is immediate from Sec. 3, 4 Renormalization of O g (s, t) In this section, we extend the calculation of anomalous dimension to the soft operator relevant to the gg → h case.The one-loop diagrams contributing to the renormalization of O uns g (s, t) in the form of (2.10) are shown in Fig. 2. The first row displays diagrams that also appeared for O γ (s, t) and their contributions can be obtained from the previous calculation as follows: The colour factor C F in diagram (a) of Fig. 1 and its mirror diagram is replaced by C F − C A /2, while the colour factor C F in diagram (b) remains the same.The colour factor of diagram (c) and its mirror changes to C F − C A /2, and the contribution of diagram (d) remains zero.This is a consequence of the colour generators placed between the soft quarks and the fundamental finite-length Wilson lines displayed in the second line of (2.10).
There are additional one-loop corrections (shown in the second row of the figure) due to gluon exchanges involving the semi-infinite adjoint Wilson lines dictated by the second line of (2.10).Due to the two semi-infinite Wilson lines, O uns g (s, t) exhibits rapidity divergences.Consequently, one needs to separate UV poles from the singular behaviour in rapidity by introducing additional regulators.In this work, we introduce an off-shellness δ ± as in (2.6), (2.7) to regulate the rapidity divergence induced by the semi-infinite Wilson lines.The presence of the semi-infinite Wilson lines and the consequent introduction of the δ ± regulators can potentially interfere with gauge invariance as will be seen below.Hence it is instructive to employ the general covariant R ξ background field gauge for the one-loop calculation of O g .

Unsubtracted operator O uns
g (s, t) As mentioned earlier, the total one-loop corrections generated by diagrams in the first row of Fig. 2 is immediate from Sec. 3.For the pole part of the sum of all diagrams, including the quark field renormalization (3.11), we find in R ξ gauge which is IR/rapidity finite.Note that the ξ-dependence proportional to C F cancels between diagrams (b), (c) and its mirror 5 in Fig. 1, and the quark field renormalization.However the ξ dependence proportional to C A , which originates from diagram (c) in Fig. 1 and its mirror diagram, still remains in the sum of the first-row diagrams.
We now proceed to compute the one-loop corrections induced by the semi-infinite Wilson lines.The contribution of diagram (a) in Fig. 2 takes the form One can now keep only the leading ε term in the integrand as in the photon case to simplify the calculations without changing the UV and rapidity behavior captured by ε and δ − , respectively.Thus we obtain where ∂ t = ∂/∂t.We will return to this expression later when we discuss how it relates to conformal symmetry.It is important to note that the ε, δ ± → 0 limit must be taken in the order of lim . The mirror diagram is obtained by changing −i∂ t to i∂ s and replacing −δ − by δ + in the square bracket above.The result in (4.3) remains the same in general R ξ gauge.
For diagram (b) in Fig. 2, we have It is worth pointing out that the λ-integral is well-defined despite the potential complications in the vicinity of λ = −t and the branch cut generated by λ + t < 0 thanks to the +i0 + prescription introduced for t necessary to make the Fourier transform of the anomalous dimension invertible (see remarks below (2.11)).This prescription moves the variable (λ + t) off the real axis into the upper complex plane, thereby circumventing the problematic region of λ ≤ −t.The imaginary unit in ln(−iµt e γ E ) of (4.4) is canceled out by the mirror of diagram (b), rendering the total sum proportional to ln(µ 2 st). 6Carrying out the calculation in general R ξ gauge, we find identical results.Diagram (c) translates into the loop integral Technically, λ 1 + t (λ 2 + s) is a variable in the upper (lower) half of the complex plane, as explained above, which bypasses the potential singularities at λ 1 = −t and λ 2 = −s.Following the procedure of taking the limit in the order of lim , we find I c to be finite.From a different perspective, the UV divergences are associated with fluctuations in the small or light-like distances, both of which are absent in diagram (c) as the two vertices connected by the gluon propagator never coincide, nor does their separation become light-like.We therefore immediately conclude that diagram (c) does not contribute to the anomalous dimension.This conclusion remains unchanged when R ξ gauge is adopted.
The same argument also applies to diagram (d), hence it does not contribute to the RG kernel either.To see this, we start from the loop integral Since (λ + t) is properly defined thanks to the +i0 + prescription for t, we can split the integration region as to simplify the calculation for t > 0. We can take δ − = 0 from the start in I 1 without changing the UV behaviour, which allows us to conclude that I 1 is UV-finite following (3.8).The contribution of I 2 can be calculated by a change of variable λ → −t + λ 2 /(−4p 2 α ᾱs) after which the λ-integral can be performed using the identity where E n (z) denotes the exponential-integral function.The resulting α-integral takes a relatively simple form and an analysis of the integrand, particularly in the endpoint region α → 0, 1, finally reveals the UV-finiteness of diagram (d) in the lim Furthermore, explicit calculations also show that diagram (d) does not contribute to the renormalization of O g in general R ξ gauge.Summing up all contributions together, we obtain the following one-loop result for the unsubtracted operator O g (s, t) in general R ξ gauge: We added the "tree-level" contribution O g (s, t) in the above for later convenience.Note that with the inclusion of the quark self-energy, the δ ± and ξ dependence is proportional only to the colour factor C A .The former is a result of interactions involving the adjoint semi-infinite Wilson lines, while the latter originates solely from the third diagram in the first row of Fig. 2 contributing with colour factor C F − C A /2 to the gluon case.
4.2 Subtracted IR/rapidity-finite O g (s, t) The dependence of the UV poles in (4.9) on the IR/rapidity regulators δ ± implies that one cannot consistently define an IR-finite anomalous dimension for O uns g , unlike for O γ .The problem arises from the two incoming energetic gluons, which give rise to the semi-infinite adjoint Wilson lines, representing the incoming colour charges "from infinity".The situation is analogous to the factorization of QED effects in the rare decay B s → ℓ + ℓ − [21], where the electric charges moving to infinity cause a similarly ill-defined soft-function anomalous dimension.The remedy is to define a subtracted soft function, which rearranges terms between the soft and (anti-) collinear functions of the factorization formula for the process, which cures a similar problem in the (anti-) collinear sector.
For the gg → h soft function, the subtraction operator is naturally defined as The position 0 is where the q qh interaction happens and the two free colour indices a and b reside at the infinity.The explicit definition of the adjoint Wilson lines in the subtraction operator follows (2.6) including the δ ± regulators.We are now ready to define the rapidity-finite soft operator entering the factorization formula for gg → h process mediated by a light-quark loop as where ⟨S g (0)⟩ denotes the vacuum expectation of S g (0) with the trivial colour tensor removed, with R ± defined here for later convenience.The superscript (ξ) emphasizes that R where all divergences in either " !0 or ± !0 limit (or both) are kept.
[• Q15 sound like "and" but should be "or"?]Combined with (4.5), we see that the rapidity regulators cancel entirely as expected and the correct double pole for the non-abelian part is acquired, Note that the subtraction factor R R + is in the same form of [21,31], in which case the QED counterpart of (4.6) was constructed.We will demonstrate in Sec.4.4 their relations to the two radiative jet functions calculated in [16].
The one-loop renormalization for the subtracted gluon operator O g (4.6) can then be obtained immediately from (4.10).In position space, it takes the following form, and like in the photon case, we have added the quark self-energy (3.9).The appearance of ln(@ s @ t /µ 2 ) in (4.10) in the momentum space is local, giving rise to ln ` `+ µ 2 = ln w µ 2 ⌘ L w in momentum space via (2.9).It is then straightforward to find the renormalization constant in momentum space Z g (w, w 0 ) by applying the Fourier transform identities (3.12) and (3.13), is quite simple.In R ξ gauge, we find where all divergences in either the ε → 0 or the δ ± → 0 limit (or both) are kept. 7 Combining with (4.9), we obtain We see clearly that both the ξ and δ ± dependence are completely removed by the subtraction.Note that after substituting takes the same form in Feynman gauge as in [21,26] where the QED counterpart of (4.10) was constructed.The one-loop renormalization for the subtracted gluon operator O g (4.10) can now be obtained immediately from (4.14).In position space, it takes the form We point out that the gauge dependence of ⟨S g (0)⟩ using the δ-regulator was already observed in [28], stemming from δ playing the role of an IR regulator when ξ ̸ = 1 at the one-loop order.While, in principle, the ξ dependence can be removed by a manipulation of the integrand [28], this would be contrary to the spirit of the subtraction procedure, where the gauge dependence of ⟨S g (0)⟩ is required to cancel the one already present in the unsubtracted soft function.It would be interesting to check whether the complete cancellation of δ ± and ξ dependence between + observed below persists to higher loop orders.
The expression ln(∂ s ∂ t /µ 2 ) is local in the momentum space, namely ln(ℓ − ℓ + /µ 2 ) = ln(w/µ 2 ) ≡ L w via (2.11).It is then straightforward to find the renormalization kernel in momentum space Z g (w, w ′ ) by applying the Fourier transform identities (3.15) and (3.16) and following the procedure in appendix A, resulting in where Γ(w, w ′ ) is defined in (3.19).This expression for the renormalization factor is in agreement with the conjecture in [16]. 8Taking C A → 0, we recover the renormalization factor (3.18) for the photon case.Finally, we obtain the anomalous dimensions: In position space, from (4.15), The momentum-space anomalous dimension follows readily from (4.16), We note that the analytical solution to the momentum-space one-loop RGE for O g , governed by γ g , has been found in [16].

Factorization formula for gg → h
From the perspective of factorization, the need for the additional regulators δ ± implies a failure of naive soft-collinear factorization.In the presence of the in-coming colour charges from −∞, representing the (anti-) collinear initial-state gluons, soft emissions require an additional rearrangement implemented by dividing the soft operator by the vacuum expectation value of the subtraction operator S g (0) (see (4.11)).
Only with this subtraction does the soft operator have a well-defined renormalizationgroup equation on its own.The rearrangement means that the modification of the soft operator must be compensated by a corresponding modification of the collinear functions in the factorization formula [21].
To implement this rearrangement for the process gg → h, we begin by quoting the factorization formula for the gg → h form factor with a light quark loop from [16], The ⊗ symbol refers to the convolution in light-cone momentum components ℓ ± .For the following, only the last term, T 3 , in the formula is relevant.The point is that the soft function S g is defined here as the matrix element of the unsubtracted operator (2.10), which cannot be consistently renormalized on its own.The dimensionally regulated, on-shell gg → h matrix element can be written (focusing on the T 3 term above) as where the factors A ±,(ξ) g (k 2 i ) arise from the (soft-decoupled) collinear matrix element of the initial-state gluon(s), For convenience, we strip the matrix element of the factors g s T a ε µ ⊥ (k).We also understand the factors H 3 , S g , J g , and A ±,(ξ) g as the regulated but unrenormalized hard, soft, and collinear functions.With this convention, (4.21) coincides with (4.19) for on-shell gluons, since A ±,(ξ) g (k 2 )| k 2 =0 = 1 to all orders, because the loop integrals are scaleless.Making the factors A ±,(ξ) g (k 2 ) explicit is crucial for the discussion.Like S g these factors cannot be consistently UV-renormalized on their own, because their anomalous dimension depends on the non-dimensional IR regulator employed to extract the UV poles.We further mention that the on-shell amplitude M gg→h is UV-finite but IR divergent.In dimensional regularization the IR divergence can be removed by multiplying the matrix element with Z −1 gg , where Z gg has an alternative interpretation as the UV renormalization constant of the SCET gluonic two-jet operator of the non-decoupled (anti-) collinear SCET gluon fields A µ ⊥,n ± = W † n ± iD µ ⊥ W n ± .The computation of this renormalization factor involves collinear, anti-collinear and soft loops and is well-defined.In the MS scheme [33], it is given by9 Before proceeding with the rearrangement, we discuss the (anti-) collinear (more precisely (anti-) hardcollinear) functions J g (±M h ℓ ∓ ) in the factorization formula.Quite generally, the collinear functions are defined as matching coefficients of timeordered products with the subleading power SCET Lagrangian [10].For the T 3 term, the relevant definition of the (position-space) collinear function reads where z − = (n + • z/2)n − , and the superscript "(0)" denotes the soft-decoupled fields, q s is the soft quark field which will eventually enter the definition of the soft operator, and ξ n − stands for the collinear quark field in the n − direction. 10The collinear function can be extracted by sandwiching the collinear part of the above matrix element between a collinear gluon and the vacuum state.Since the renormalization of the soft-quark function requires the introduction of a non-dimensional IR regulator δ ± , the collinear function must be computed with a consistent regulator as well.As pointed out in (2.9), the δ ± -regulators are related to the off-shellnesses of the (anti-) collinear gluons through The technical details of this computation in general R ξ gauge can be found in the App.B. The main result is that the collinear functions J g (±M h ℓ ∓ ), defined as above as the matching coefficient from SCET I onto SCET II , remain the same as the on-shell ones and are, in fact, gauge-invariant even with an external off-shell gluon.To arrive at this result, it is important that the matrix element of / A (0) ⊥,n − (sn + ), which coincides with (4.21) up to a phase factor, and the Dirac matrix, is not unity for an off-shell gluon, but involves (see (B.3) for the ε-exact expression) in background field R ξ gauge.The dependence of this expression on the off-shellness k 2 and gauge-fixing parameter ξ is exactly what is required to make the collinear functions gauge-independent and finite as k 2 → 0.
We are now in the position to return to (4.20) when the matrix element is computed with a small off-shellness k 2 i ≪ M 2 h of the external gluons.As said, J g (±M h ℓ ∓ ) remain unaltered, but S g turns into the unsubtracted soft function S uns,(ξ) g , which depends on the δ ± -regulators as well as on ξ.On the other hand, from (4.11), (4.12), we have and from (4.13) we can separate − in a symmetric way as11 R (ξ) This allows us to rewrite (4.20) in the form where now the soft function is independent of the δ-regulator and ξ, and can be renormalized on its own.The subtraction factors are multiplied back to the externalcollinear gluon factors.Combining (4.28) with (4.26) and making use of where the gauge and off-shell regulator dependence also cancels exactly.Thus the rearrangement not only allows to UV-renormalize the soft function on its own, but also the external collinear factors, implementing a restoration of soft-collinear factorization in very much the same way as in [21,26] for external fermions.Incidentally, the second line of (4.29) equals (R These poles are removed by the UV renormalization of the external-collinear factors.

Conformal symmetry and higher-order results
In this section, we demonstrate the power of conformal symmetry in promoting the anomalous dimensions of the soft operators O γ and O g to higher orders.In particular, we establish an all-order relation between γ γ , the anomalous dimension of O γ , and the evolution kernel of the leading-twist non-local heavy-light quark operator.In addition, we construct a general ansatz for the anomalous dimension of O g at higher orders in α s , test it against the RG consistency using known two-loop hard and gluon jet functions.Another major advantage of the conformal symmetry technique is that it offers a systematic approach to solving the RGEs and therefrom reveals the behaviours of the soft functions, which we discuss in App. C. We start by recalling that the conformal transformations in the light-like directions form an SL(2, R) group generated by the so-called collinear conformal generators.The n + -collinear conformal generators can be written explicitly as [23], where j = 1 is the conformal spin of the soft quark and θ s = s∂ s the Euler operator.
The generators of n − -collinear conformal transformations are written analogously as These operators are induced by projecting the D-dimensional conformal group generators onto the n ± light-like directions.The commutation relations read, S + , S − = 2 S 0 , S 0 , S ± = ± S ± ; T + , T − = 2 T 0 , T 0 , T ± = ± T ± , (5.3) while S i , T j = 0 for i, j ∈ {0, ±}.Clearly, these operators form a direct product of two SL(2, R) algebras.At the one-loop order, the anomalous dimension of the soft operators O γ and O g in (3.13) and (4.17) can be expressed equivalently in terms of these conformal operators in the compact form [36,37], ) The simplicity of the anomalous dimensions in the conformal representation shown above motivates us to explore the application of conformal symmetry techniques beyond the leading order in the subsequent sections.

Relating evolution kernels of soft operators and the B-LCDA
Extending the one-loop results (3.21), (4.18), the position-space anomalous dimensions for O γ (s, t) and O g (s, t), which we denote collectively as γ S , can be tentatively written in the position-space representation to all orders in perturbation theory as where the s ↔ t symmetry has been implemented.The function ω S (α s , α) (S = γ, g) has no power-like divergences at the end-points, and Γ S (α s ) is a constant.Γ g cusp ≡ Γ F cusp −Γ A cusp /2 and Γ γ cusp ≡ Γ F cusp where "F" ("A") stands for the fundamental (adjoint) representation. 12The appearance of the cusp anomalous dimension is a consequence of the universality of IR divergences in QCD.By analysis of the Feynman diagrams and requiring that γ S is compatible with total space-time translations, the anomalous dimension can in principle have an additional contribution of the form, which, at the one-loop order, could come from diagram (d) in Fig. 1.Such a contribution, however, necessarily moves the two soft quarks off the light cone, contradicting the fact that soft interactions do not alter the (anti-) collinear directions of the underlying fields.From another perspective, such contributions will lead to momentum-space anomalous dimensions depending on ℓ + and ℓ − separately instead of only their product w = ℓ + ℓ − , as can be clearly seen from the Fourier transform detailed in App.A, violating boost invariance of the gg → h decay process.Hence, we conclude that the anomalous dimension for both O γ and O g take the form (5.6) to all orders.Namely, we can write, which allows us to apply the collinear conformal symmetry technique to the collinear and anti-collinear sectors independently.We note that the form of the all-order anomalous dimension for O γ and O g in (5.6) preserves the location of the Higgs production vertex hqq conveniently set at 0 where the two light-like Wilson lines connect.This is similar to the leading-twist B-meson LCDA evolution, under which the location of the heavy quark is also fixed at the origin.Consider the anomalous dimension of O γ at the critical point of D-dimensional QCD with the critical coupling α * s such that β(α * s ) = −2α * s (ε + β(α * s )) = 0 [24] where β(α s ) is the QCD β-function in four dimensions.First, we note that the full special conformal generator K µ acting on O γ receives quantum corrections, where we have applied the sn − ↔ tn + symmetry.Clearly, the tree-level special collinear conformal generators S + and T + are induced by n ∓ • K, respectively.Similarly, the dilatation operator of the collinear conformal group acting on O γ takes the form with the α * s dependence of ∆D(α * s , s) fully determined by the anomalous dimension γ γ (s, µ) and the QCD β-function.This is in contrast to ∆K(α * s ; s), which receives additional contributions from the so-called "conformal anomaly" that can be calculated systematically through a diagrammatic approach [25].Since the conformal symmetry at the critical point is exact, we can consequently establish the following exact commutation relations [27,38] The commutation relations (5.11) are readily solved [27] yielding the solution, where we replaced α * s by an arbitrary α s and traded the ε-dependence with the ordinary QCD β-function to obtain γ γ (s, µ) in four-dimensions thanks to the properties of MS-like schemes [24].
We can further identify K(α s ; z) and D(α s ; z) as the special conformal and dilatation generator acting on the leading-twist heavy-light operator [27] whose matrix element defines the leading-twist B-meson LCDA.The reason is that K(α s ; s) and D(α s ; s) can be systematically computed by considering their UV divergences when acting on the gauge-invariant operator X(0)[0, sn + ]q(sn + ) from which the commutation relations (5.11) follow.Here X(0) must be built of fundamental QCD objects such as the quark field, gluon field strength tensor, or a Wilson line so that definite conformal properties apply.Since each local light field (q, G µν ) is associated with a conformal generator [23] which depends on the location of the light field, X(0) can only be a non-local object, i.e., a Wilson line.On the one hand, (5.11) implies that the special conformal symmetry is preserved only in one particular direction n + or n − , on the other hand, the dilatation symmetry is broken by the cusp term Γ cusp .This implies that the Wilson line, which we still denote as X(0), must lie in a direction v µ different from n µ + , hence forming a cusp with the Wilson line [0, sn + ], and extend to infinity to ensure gauge invariance.In principle, no further constraints can be imposed on v µ from the perspective of conformal symmetry.But as we have explicitly shown in Sec.4.1, light-like Wilson lines induce unnecessary complications in the form of rapidity divergences.Furthermore, space-and time-like Wilson lines are equivalent for our consideration.Therefore, we are free to choose v µ to be a time-like vector fixing X(0) = [−∞v, 0] with v 2 = 1.We can further identify X(0) as the heavy quark field hv (0) in the heavy-quark effective theory (HQET), i.e., X(0) = hv (0), because hv (0) = Q(−∞v)[−∞v, 0] with Q being a sterile scalar field at infinity decoupled from the QCD interactions [39].It is now evident that the collinear conformal generators (anomalous dimension) associated with O B = hv (0)[0, sn + ]q(sn + ) coincide with (5.1) ((5.4)) at the one-loop order up to a constant that cannot be determined from symmetry considerations.This serves as a consistency check.We can thus make full use of the known results for the anomalous dimension and conformal generators of O B [27,37].
Comparing to the evolution kernel of the operator hv (0)[0, sn + ]q(sn + ) [27], we see that γ γ (s, t, µ) and H B (s, µ) + H * B (t, µ), with " * " denoting the complex conjugate, may differ by the constant 2(Γ γ (α s ) − Γ + (α s )).It can be fixed by considering the local limit s, t → 0 of the identity This difference removes the locally divergent Γ cusp ln(μ 2 st) term in the positionspace representation, making the limit well-defined.While the remaining terms in γ γ (s, t, µ) approach 2γ q , the non-cusp part of the quark collinear anomalous dimension [40], H B (s, µ) + H * B (t, µ) becomes a double-copy of the anomalous dimension 2(γ q + γ Q ) that appears for the heavy-to-light operator.We therefore conclude that13 Putting everything together, we find, which is valid to all orders in perturbation theory.
As an application of (5.16), we obtain the two-loop O(α 2 s ) correction to the anomalous dimension γ γ (s, t) by first fixing the constant part.The explicit expression for Γ + in (5.13) is known to the two-loop order and reads [27] where T F = 1/2 is the trace of the SU (N c ) generators in the fundamental representation and n f is the number of active flavors.The heavy-quark anomalous dimension γ Q at two loops reads [41] Then from (5.15), the constant 2Γ γ entering the conformal representation of the anomalous dimension γ γ (s, t) is The full expression for the anomalous dimension in the conformal representation then follows from (5.12), The cusp anomalous dimension to the two-loop order reads [38], and the O(α s ) correction to the special conformal generator K is available from [27], with O(s) being a test function.
Similarly, the constant term Γ γ in the position-space representation (5.6) (with S = γ) of γ γ (s, t) can be computed using Γ γ = γ + − γ Q with γ + given explicitly in Eq. ( 27) of [27] yielding, The anomalous dimension in position-space then takes the form at two loops, where We confirmed that the Fourier transform of (5.24) reproduces the two-loop momentum-space anomalous dimension γ γ reported in [17], where it was extracted from RG consistency of the γγ → h factorization formula.We note that the twoloop expressions of γ γ , which participates in the light-quark induced γγ → h process at the three-loop order (i.e., to O(α em α 2 s )), has also been verified by comparing the predicted large logarithms with numerical calculations [15].
It is also interesting to point out that we expect the Casimir scaling for O γ to hold at least to the two-loop order.Namely, by having the soft fields and the two finite-distance Wilson lines in (2.1) in an arbitrary representation R of SU (N c ), 14the anomalous dimension γ R γ of the resulting operator O R γ can be obtained by simply replacing C F → C R in (5.20) and (5.24),where C R denotes the Casimir value of representation R.This follows from the well-known Casimir scaling of the cusp anomalous dimension, 15 the Casimir scaling of the heavy-quark anomalous dimension at this order [45], the trivial change of colour factors in the one-loop correction to the conformal operator K in (5.16), and the Casimir scaling for Γ + corresponding to C F → C R in (5.17).

Bootstrap the two-loop anomalous dimensions of O g
In this subsection, we investigate the possibility of bootstrapping the anomalous dimension of O g to higher orders without explicit calculations.It turns out that, despite the additional S − and T − contributions in comparison with the case of O γ (see (5.4) and (5.5)), it is still feasible to some extent.First we note that S − and T − are induced by the adjoint Wilson line extending to infinity and their appearance is a direct consequence of the rapidity divergence in the original O uns g .The subsequent introduction of the subtraction operator (4.10) brings extra constants to γ g as well.
To be more precise, the first commutation relation in (5.11) becomes invalid, and we can only expect a single commutator identity, which we conjecture to hold at least to two loops with the conformal operators unchanged from the O γ case: Here, Γ A cusp reflects the presence of adjoint Wilson lines.In this case, the anomalous dimension γ g (α * s ; s, µ) can no longer be solved for unambiguously, and we resort to the general ansatz which is motivated by the one-loop result, satisfies the commutation relation (5.26), and takes the form of (5.6) in the position-space representation.Here K is taken to be identical to the B-meson LCDA case, Γ g is a constant, and we used to construct the above ansatz.(Note that S − receives no quantum corrections and hence its canonical form holds to all orders.)This ansatz also complies with the fact that the ln( S − /µ) contribution is exclusively generated by the adjoint Wilson lines interacting with the soft quarks, giving rise to the adjoint cusp anomalous dimension.
In addition, Eq. (5.27) also satisfies the equation the significance of which calls for further investigation.We did not find an exact way to extract the O(α 2 s ) correction to Γ g from known results of anomalous dimensions of local operators.However, it is reasonable to expect that, 2Γ g = 2Γ γ − Γ ∧ + ∆ (5.30) with 2Γ g being the constant part of γ g (s, t) in the conformal representation (see (5.8) and (5.27)).Γ ∧ denotes the non-cusp anomalous dimension of the vacuum expectation value of the subtraction operator (4.10) which can be found in [46,47], Interestingly, it differs from the so-called eikonal function f eik [47] extractable from the Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) kernel, the anomalous dimension γ DY s for the Drell-Yan soft function [48], and the soft anomalous dimension near threshold [45] only by a ζ 3 term.It is also worth noting that all these anomalous dimensions exhibit Casimir scaling at the two-loop order.∆ in (5.30) accounts for the mismatch arising from contributions beyond Γ γ and Γ ∧ that start to appear from the two-loop order.By imposing the consistency relation required by the factorization theorem and using the two-loop hard and jet (J g ) function available at the two-loop order in [35], we can fix ∆ in (5.30) to be which again contains only a ζ 3 term.

Conclusions
We established a framework to obtain the renormalization of the soft-quark functions that appear in the context of Higgs form factors generated by light-quark loops by direct computation from their operator definitions.While the renormalization-group equations for these functions were already known under the assumption that the factorization formula for the physical form factors holds [16,17], the direct computation provides an independent confirmation of this assumption.For the γγ → h case, the position-space renormalization technique employed in the present work provides a significant simplification relative to the previous direct computation [18].In case of gg → h, it reveals that the operator definition that follows from the naive application of SCET factorization must be amended by a subtraction, related to the in-coming colour charges, which turns out to be conceptually identical to the subtraction necessary for the factorization of electromagnetic corrections to the form factor in B s → µ + µ − decay [21].Perhaps most importantly, the methods developed here open a path to compute the renormalization of other soft functions that have recently appeared in the context of hard-scattering factorization beyond the leading power.
For the case of the soft operator O γ relevant to γγ → h, we established an all-order relation of its anomalous dimension to the one that appears for the leadingtwist B-meson light-cone distribution amplitude.The present work also explores for the first time the consequences of conformal symmetry, which has been successfully applied to the computation of the two-loop renormalization of the leading-twist heavy-light string operator [27], for non-string light-cone operators such as O γ and O g with the latter containing additional complications due to semi-infinite Wilson lines.Thus, it would be highly beneficial from both the theoretical and phenomenological perspective to have a direct two-loop calculation for γ g to verify the conjecture (5.27) and track down the origin of the discrepancy ∆ in (5.32).Further investigations and applications of collinear conformal and position-space techniques may expose "hidden relations" that explain why various anomalous dimensions differ at the two-loop order only by ζ 3 terms.
where we applied (3.15).We next change variables to ℓ − = w/ℓ + , ℓ ′ − = w/ℓ ′ + , leading to We used the "initial condition" O γ (ℓ ′ + , ℓ ′ − ) = O γ (w ′ ) thanks to Lorentz invariance, namely, O γ in momentum space only depends on the scalar product The Fourier transform over t in (A.1) can be carried out straightforwardly as a complex conjugate to (A.3) yielding, Thus in total, the logarithmic contribution ln (stµ 2 e 2γ E ) in momentum space reads Following the same procedure and making use of the identity (3.16), the Fourier transform for the remaining position-space renormalization factor in (3.13) can be worked out.

B Collinear function with an off-shell external gluon
To extract the collinear function, we sandwich the matching equation (4.24) between a soft quark state with momentum p s and a slightly off-shell gluon.We then project out the soft quark and the soft Wilson line leaving a collinear vacuum to gluon matrix element of collinear fields.This turns (4.24) into The soft quark momentum enters the collinear function only with its n + -component p µ s + = 1 2 (n − • p s )n µ + due to multipole expansion.The off-shell momentum-space collinear function then depends on p The dashed line represents the (decoupled) soft quark field, the symbols "⊗" and "•" denote the L (1) ξq insertion and [ ξn − W n − ](0), respectively.Note the presence of a symmetry factor 1/2 of diagram (i).The contribution from the one-loop gluon vacuum polarization cancels between the left-and right-hand side of (B.1) and is therefore not shown.
Implementing the relevant Feynman rules listed in the appendix of [49] for the subleading power SCET Lagrangian, we find that the diagrams in Fig. 4 contribute to the off-shell jet function in MS-scheme as The dashed line represents the (decoupled) soft quark field, the symbol "⌦" ("•") denotes the L (1) insertion ([ ⇠(0) n W n ](0)).Note the presence of a symmetry factor 1/2 in diagram (i).
-38 -order of ↵ s in R ⇠ gauge using the original subleading power SCET Lagrangian L (1) n (z) insertion.The dashed line represents the (decoupled) soft quark field, the symbol "⌦" ("•") denotes the L (1) insertion ([ ⇠(0) n W n ](0)).Note the presence of a symmetry factor 1/2 in diagram (i).4, it also comes with a symmetry factor 1/2.The vacuum polarization diagram cancels between the left-and right-hand side of the matching equation and, hence, can be dropped.
The calculation is performed in the background field R ξ gauge [31].Alternatively, conducting the calculation in R ξ gauge without employing the background field method leads to identical result.The virtuality of the external gluon serves as an infrared (IR) regulator, implying k 2 ≪ p 2 .Consequently, in (B.2), we retain only the non-vanishing terms in the k 2 → 0 limit.The right-hand side of the matching equation (B.1) receives quantum corrections when the external gluon is off-shell.At the one-loop level, they are generated by the diagram depicted in Fig. 5. Its contribution to the off-shell collinear function at NLO reads, Here, we have neglected the contribution from the gluon vacuum polarization as it cancels between the left-and right-hand side of (B.1).From (B.2) and (B.3), we find the off-shell collinear function to O(α s ) to be, Feynman diagrams needed for computing the o↵-shell gluon jet function J g at the one-loop level using the gauge-invariant SCET building blocks A ?µ n and n .The symbols "⌦", " " and "•" denote the SCET collinear fields / A ? n (z), n (z), and n (0), respectively.Diagrams only involving A ?µ n (i.e.Fig. 5 and the one-loop gluon vacuum polarization correction to the first diagram in Fig. 4) are neglected since their contributions are canceled between the two sides of eq.(B.1).
Employing the alternative approach of computing the jet function in terms of the SCET building blocks provides a valuable cross-check to our result in (B.4).In particular, in this case, there is no need to calculate corrections to the right-hand side of the matching equation (B.1), as they are entirely o↵set by the corresponding one-loop diagrams originating from the left-hand side, where only / A ? n contributes to the interaction (not shown in Fig. 6).Performing the one-loop calculation in R ⇠ gauge using the background method and summing up all contributing diagrams confirms our result of J OfS,(1) g in (B.4) to all orders in ".Furthermore, we have also computed the o↵-shell jet function in the light-cone gauge where we find complete agreement to J OfS,(1) g in (B.4).The main advantage of adopting the light-cone gauge is that the number of diagrams is reduced to a minimum.

C Solving RGE using conformal symmetry
In this appendix, we investigate the application of conformal symmetry to solving the RGE for O and O g inspired by the study of B-LCDA [17,34,[50][51][52][53].First we note that in terms of the collinear conformal generators in (5.1) and (5.2), the position--35 -Figure 6. Feynman diagrams contributing to the off-shell collinear function J g at the one-loop level using the SCET building A ⊥µ n − and χ n − .The symbols "⊗", "•" and "•" denote the fields / A ⊥ n − (z), χ n − (z), and χ n − (0), respectively.Diagrams only involving A ⊥µ n − (i.e.Fig. 5 and the one-loop gluon vacuum polarization correction to the first diagram in Fig. 4) are neglected since their contributions are canceled between the two sides of (B.1).where we performed the ε-expansion to O(1/ε), since we are only interested in the UV renormalization of the function.Remarkably, the ε-divergence of the off-shell collinear function coincides with the on-shell one.It is also important to note that while both sides of the matching equation (B.1) exhibit gauge-dependence in the single-pole term 1/ε, the off-shell collinear function defined as a matching coefficient remains gauge-independent.
Equivalently, the collinear function can be computed in terms of the gaugeinvariant collinear building blocks in SCET, in which case the subleading SCET Lagrangian takes the form L n − ξ n − ](z).The relevant Feynman diagrams in this approach are collected in Fig. 6.Computing the jet function using the SCET building blocks provides a valuable cross-check of the result (B.4).In particular, in this case, there is no need to calculate corrections to the right-hand side of the matching equation (B.1), as they are entirely canceled by the corresponding one-loop diagrams on the left-hand side, where only / A ⊥ n− contributes to the interaction (not shown in Fig. 6).Performing the one-loop calculation in background-field R ξ gauge and summing up all contributing diagrams confirms the expression (B.4) for J OfS,(1) g to all orders in ε.
Furthermore, we also computed the off-shell collinear function in the light-cone gauge and again found complete agreement with (B.4).The main advantage of adopting the light-cone gauge is that the number of diagrams is smaller.

C Solving RGE using conformal symmetry
In this appendix, we investigate the application of conformal symmetry to solving the RGE for O γ and O g inspired by the study of B-LCDA evolution [29,37,[50][51][52].First, we note that in terms of the collinear conformal generators in (5.1) and (5.2), the position-space renormalization factor for O g in the MS-scheme can be rewritten compactly as The corresponding one-loop anomalous dimension reads, H (1) . (C.2) The position-space RGE is simply given by d d ln µ + γ g O g (s, t) = 0 . (C.3) The non-local information of γ g is hidden in the conformal operators.Since the two sets of generators { S + , S − , S 0 } and { T + , T − , T 0 } commute with each other, the one-loop kernel can be disentangled as H (1) = H s + H t , with similarly for H (1) t with S ± → T ± .Consequently, the properties of the full one-loop kernel H (1) are completely determined by H (1) s .Before delving into the full solution of (C.4), we take a look at three degenerate cases.

Abelian limit: C
In this case, ln( S − /µ) is absent, which aligns with the evolution kernel of the leading-twist B-LCDA and O γ .Therefore, the eigenfunctions of H (1) s are completely determined by the eigenequation of S + [17,37], since the underlying system has only one degree of freedom (one variable s for H where ψ 0 (z) denotes the polygamma function of the zeroth order.The eigenfunctions of H s are determined by solving the eigenequation of θ s with powerlaw solutions, expressed as θ s s λ = λ s λ .This scenario makes no distinction between the fundamental and adjoint Wilson lines because all fields live in the adjoint representation of the gauge group.
We now turn to the general case (C.4).We start from the representation The final form of the eigenfunction for H (1) is simply the product ϕ λ 1 (s, µ)ϕ λ 2 (t, µ).
Compared to other techniques, the conformal symmetry-inspired Mellin-space approach to solving the RGE is systematic and well-suited for accommodating higherorder corrections to the anomalous dimensions [53].The analytic solution to the one-loop RGE in momentum space can be derived by applying the Fourier and Mellin transform to (C.14) using (C.10).We refer [16] for the explicit expression.

Figure 1 .
Figure 1.One-loop corrections to O (s, t).Double lines denote the two finite Wilson line segments, the two solid lines the soft quark fields and the wiggle lines are gluon propagators.Mirror diagrams to (a) and (c) as well as the self-energy of the external quarks are not shown.The self-energy of the lightlike Wilson lines vanishes in Feynman gauge and therefore are neglected.

Figure 1 .
Figure 1.One-loop corrections to O γ (s, t).Double lines denote the two finite Wilson line segments, the two horizontal solid lines correspond to the soft quark fields, and the curly lines represent gluon propagators.Mirror diagrams to (a) and (c), the self-energy of the external quarks, and the self-energy of the light-like Wilson lines (which vanishes in R ξ gauge) are not shown.

Figure 2 .
Figure 2. Diagrams responsible for O g (s, t) renormalization at one loop.The double lines superposed with wiggly ones are semi-infinite adjoint Wilson lines.Mirror diagrams to (a), (b) and (d) are not shown.Contributions from diagrams in the first row are immediately obtainable from eq. (3.10) with appropriate changes in the color factor.

Figure 2 .
Figure 2. Diagrams responsible for O uns g (s, t) renormalization at one loop using the representation in the second line of (2.10).The double lines superposed with curly ones are semi-infinite adjoint Wilson lines.Mirror diagrams to (a), (b) and (d) are not shown.Contributions from diagrams in the first row are immediately obtainable from eq. (3.12) with appropriate changes in the colour factor.

0 Figure 3 .
Figure 3. Double lines superposed by wiggly lines represent light-like adjoint Wilson lines.Simple wiggly line denotes the gluon propagator.

Figure 3 .
Figure 3. Double lines superposed by curly ones represent light-like adjoint Wilson lines.The simple curly line denotes the gluon propagator.The self-energy of the Wilson lines vanishes in R ξ gauge hence not shown.

Figure 4 .
Figure 4. Feynman diagrams necessary for computing the LHS of eq.(B.1) up to the order of ↵ s in R ⇠ gauge (with and without the background field formalism) using the original subleading power SCET Lagrangian L (1) q⇠(z) = qs (z )W † n (z)i / D ?n ⇠ (0)n (z) insertion.The dashed line represents the (decoupled) soft quark field, the symbols "⌦" and "•" denote the L(1) q⇠ insertion and [ ⇠(0) n W n ](0), respectively.Note the presence of a symmetry factor 1/2 in diagram (i).The soft quark momentum enters the calculation of the jet function only with its n + -component p µ s + = 1 2 (n • p s )n µ + due to multi-pole expansion.The contribution from the one-loop gluon vacuum polarization cancels between the LHS and RHS of (B.1) hence not shown.

Figure 4 .
Figure 4. Feynman diagrams contributing to the left-hand side of (B.1) O(α s ) in R ξ gauge (with and without the background field formalism) using the original subleading power SCET Lagrangian L

Figure 6 .
Figure 6.Feynman diagram responsible for the RHS contribution of the matching equation (4.28).Similar to diagram (i) in Fig.5, it also comes with a symmetry factor of 1/2.

Figure 6 .
Figure 6.Feynman diagram responsible for the RHS contribution of the matching equation (4.28).Similar to diagram (i) in Fig.5, it also comes with a symmetry factor of 1/2.The vacuum polarization diagram cancels between the LHS and RHS of the matching equation hence can be skipped.

Figure 5 .
Figure 5. Feynman diagram responsible for the right-hand side (RHS) contribution of the matching equation (B.1), i.e., the matrix element of the gauge-invariant building block A (4.21).Similar to diagram (i) in Fig.4, it also comes with a symmetry factor 1/2.The vacuum polarization diagram cancels between the left-and right-hand side of the matching equation and, hence, can be dropped.