Threshold Factorization Redux

We reanalyze the factorization theorems for Drell-Yan process and for deep inelastic scattering near threshold, as constructed in the framework of the soft-collinear effective theory (SCET), from a new, consistent perspective. In order to formulate the factorization near threshold in SCET, we should include an additional degree of freedom with small energy, collinear to the beam direction. The corresponding collinear-soft mode is included to describe the parton distribution function (PDF) near threshold. The soft function is modified by subtracting the contribution of the collinear-soft modes in order to avoid double counting on the overlap region. As a result, the proper soft function becomes infrared finite, and all the factorized parts are free of rapidity divergence. Furthermore, the separation of the relevant scales in each factorized part becomes manifest. We apply the same idea to the dihadron production in $e^+ e^-$ annihilation near threshold, and show that the resultant soft function is also free of infrared and rapidity divergences.


I. INTRODUCTION
Factorization theorems in which high-energy processes are divided into hard, collinear, and soft parts are essential in providing precise theoretical predictions. Though it is difficult to probe the threshold region experimentally, it is theoretically both interesting and tantalizing to exploit the factorization near threshold in Drell-Yan (DY) process and in deep inelastic scattering (DIS). A prominent distinction near threshold is that there exists nonvanishing soft interaction since the real contribution does not cancel the virtual contribution completely due to the kinematic constraint near threshold [1,2]. This distinctive feature near threshold was also discussed in Refs. [3][4][5][6] in the framework of the soft-collinear effective theory (SCET) [7][8][9][10].
It is well established in full quantum chromodynamics (QCD) that the structure functions F DY for DY process and F 1 for DIS can be schematically written in a factorized form as [1,2] where the kinematical variables τ = q 2 /s = Q 2 /s, z = τ /(x 1 x 2 ), and x = −q 2 /2P · q = Q 2 /2P · q are all close to 1. Here q µ is the hard momentum carried by a photon, s is the momentum squared of incoming hadrons in DY, and P µ is the momentum of an incoming hadron. The hard functions H DY (Q) and H DIS (Q) describe the physics at the large scale Q. The PDF f i/N (x, µ) is the collinear part from the incoming hadrons, which represents the probability that a specific parton of type i in a hadron N has a longitudinal momentum fraction x. The jet function J(Q √ 1 − x, µ) describes the energetic collinear particles in the final state in DIS. Finally the soft functions S DY (Q(1 − z), µ) and S DIS (Q(1 − x), µ) describe the emission of soft particles. (Actually S DIS (Q(1 − x), µ) = 1 to all orders in α s , as will be discussed later.) And '⊗' implies an appropriate convolution.
If we consider the factorization in SCET, the surmised form of the factorized structure functions in full QCD near threshold in Eqs. (1) and (2) looks reasonable at first glance, but there are delicate and discomfiting aspects. In fulll QCD, the PDF has the collinear divergences yielding the DGLAP equation through special definitions of the PDF or with other techniques [1]. And the soft part with the eikonal cross sections is IR finite. In SCET, if we naively separate the collinear, soft modes, there appear IR and rapidity divergences in each factorized part. In this paper, we address this problem by considering additional modes required in SCET near threshold.
The issues on the factorization near threshold in SCET are summarized as follows: First, since the incoming active partons take almost all the hadron momenta, the emission of additional collinear particles is prohibited. It means that only the virtual correction contributes to the collinear part [11,12]. Therefore, if we consider the collinear interaction alone, the PDF near threshold should be definitely different from the PDF away from threshold, where the latter includes the effect of real gluon emissions.
Second, even though we accept Eqs. (1) and (2) and compute the factorized parts in SCET perturbatively, we encounter infrared (IR) divergences not only in the PDFs but also in the soft functions. The IR divergence in the PDF can be safely absorbed in the nonperturbative part, but the IR divergence in the soft function is a serious problem since it destroys the factorization and prevents a legitimate resummation of large logarithms of 1 − x or 1 − z.
The existence of the IR divergence has been casually disregarded in the belief that the final physical result should be free of it. But it was pointed out in Ref. [13,14] that the soft functions indeed contain IR divergence by carefully separating the IR and ultraviolet (UV) divergences. In Ref. [13,14], we have suggested how some of the divergences can be transferred from the soft part to the collinear part to make the soft function IR finite, while the collinear part reproduces the PDF near threshold. Basically this amounts to reshuffling divergences based on physics near threshold, but it was difficult to explain how the scale dependence in each part can be established consistently to resum large logarithms. For example, the typical scale for the PDF is µ ∼ Λ QCD or larger, and the relevant scale to the soft function is Third, as we will see later, the soft parts in Eqs. (1) and (2) include the rapidity divergence. The rapidity divergence arises when the product of the lightcone momenta remains constant, while each component goes to zero or infinity [15,16]. Though the rapidity divergence exists in each factorized part, the scattering cross section, which is a convolution of the factorized parts, is free of the rapidity divergence. However, the cancellation of the rapidity divergence occurs only when the invariant masses of the different modes are of the same order. Near threshold, the invariant masses of the collinear particles and the soft particles are different and there is no reason for the rapidity divergence to cancel in the sum of the collinear and the soft parts.
The naive extension of the factorized form in Eqs. (1) and (2) to SCET contains all these problems. And the questions are how we can obtain a consistent factorization formula near threshold, and how it can be connected to the factorization away from threshold. The predicament can be resolved by noticing that, near threshold, there is an additional degree of freedom with small energy, collinear to the beam direction. The small energy scale is given by ω = Q(1 − z) in DY process or ω = Q(1 − x) in DIS with z ∼ x ∼ 1, and lies between the large scale Q and the low scale Λ QCD . The main points of our paper are to identify the new degrees of freedom, to incorporate them in the definition of the PDF and the soft functions, to calculate the perturbative corrections at order α s , and to show that we obtain the proper factorization near threshold.
The high-energy processes including the threshold region can be efficiently described by SCET. In SCET, the n-collinear momentum scales as (n · p, n · p, p µ ⊥ ) ∼ Q(1, λ 2 , λ), where Q is the large energy scale and λ is the small parameter for power counting in SCET. The n-collinear momentum scales as Q(λ 2 , 1, λ). Here n µ and n µ denote the lightcone vectors satisfying n 2 = 0, n 2 = 0 and n · n = 2. In order to describe the threshold region, from the n-and n-collinear interactions we decouple the n-collinear-soft (csoft) and the n-csoft modes respectively, which scale as p µ n,cs ∼ ω(1, α 2 , α), p μ n,cs ∼ ω(α 2 , 1, α).
Here the power-counting parameter α for the csoft modes satisfies the relation ωα ∼ Qλ, such that the collinear and csoft particles have p 2 ∼ ω 2 α 2 ∼ Q 2 λ 2 . The additional partition of the collinear modes depends on the new scale introduced near threshold. The csoft modes are soft since the overall scale is governed by the small scale ω, but the momentum components scale like collinear momenta. The nomenclature for the collinear-soft modes varies, as in the collinear-soft (csoft) modes [17,18], the coft modes [19], and the soft-collinear modes [20], referring to the modes with similar momentum scaling in different situations. Here we will simply call these modes the csoft modes. The important feature near threshold is that the incoming active parton cannot emit real collinear particles, but the particles in the csoft modes can be emitted. And we define the PDF near threshold including the csoft modes. The new definition of the PDF covers the threshold region as well as the region away from threshold since the effect of the csoft modes away from threshold is cancelled to all orders, while it correctly describes the PDF near threshold. To avoid double counting on the overlap region, the contribution of the csoft modes should be subtracted from the soft part to obtain the soft function. Here the soft modes near threshold scale as p s ∼ ω(1, 1, 1), hence the csoft mode can be also considered to be a subset of the soft mode.
The effect of the csoft modes in the collinear and soft modes is more interesting when we consider the details of the higher order calculations. Without the csoft modes, both the soft part and the PDF contain the IR and rapidity divergences. But when the csoft modes are subtracted from the soft part, the resultant soft functions are free of the IR and rapidity divergences, and the PDF is free of rapidity divergence when the contribution of the csoft modes is included.
The structure of the paper is as follows: In Sec. II, the main idea of incorporating the csoft modes is presented in SCET. The PDF and the soft functions are defined near threshold in DY and in DIS processes. In Sec. III, the soft functions and the PDF are computed at order α s with the csoft modes. In Sec. IV, we consider the factorization of the dihadron production in e + e − annihilation near threshold, in which the effect of the csoft modes is included in the fragmentation functions. Finally we conclude in Sec. V.

A. Extension of the PDF to threshold
The main issue in constructing the PDF near threshold is how to implement the tight kinematic constraint, and how to relate it to the PDF away from threshold. Near threshold, the incoming partons cannot emit real collinear particles, hence only the virtual corrections contribute. On the other hand, the emission of the csoft modes is allowed. Therefore we start from defining the PDF near threshold by subdividing the collinear field into the collinear and the csoft modes. In Ref. [21], the decomposition of the collinear and the csoft modes has been performed in order to describe a fragmenting process to a jet with a large momentum fraction z. It can be adopted in defining the PDF near threshold.
We first decouple the soft mode ∼ Qζ(1, 1, 1) near threshold from the collinear mode, where ζ is a small parameter to characterize 1 − z or 1 − x. This is obtained by redefining the collinear fields in terms of the soft Wilson line [9] as where the soft Wilson line Y n is given by Then we extract the csoft mode ∼ Qζ(1, α 2 , α) in the collinear sector. That is, the collinear gluon A µ n is further decomposed into A µ n → A µ n +A µ n,cs . The resultant collinear mode scales with the large energy Q, while the csoft mode scales with Qζ. After the decomposition, the covariant derivative can be written as iD µ = iD µ c + iD µ cs = P µ + gA µ n + i∂ µ + gA µ cs . Here P µ (i∂ µ ) is the operator extracting the collinear (csoft) momentum, which applies only to the collinear (csoft) operator.
The collinear quark distribution function, which is the PDF away from threshold, is defined as where P + ≡ n · P is the operator extracting the largest momentum component from the collinear field. The average over spin and color is included in the matrix element. The combination χ n = W † n ξ n is the collinear gauge-invariant block with the collinear Wilson line W n and the collinear quark field ξ n .
For the proper treatment of the PDF near threshold, we need to include the csoft mode, which describes the emission along the beam direction. It is implemented by replacing n·iD c with n · iD = n · iD c + n · iD cs . Then we define the PDF as which covers all the regions, near and away from threshold. Note that the expression for φ q/N is invariant under the collinear and csoft gauge transformations. In order to show the gauge invariance order by order in power counting in a manifest way [22], we redefine the collinear gluon as where the collinear Wilson lineŴ n is expressed in terms of the newly defined collinear gluon fieldÂ µ . Then the covariant derivative can be written as where the hats in W n and A n are removed for simplicity. Then the delta function in Eq. (7) is written as We can decouple the csoft interaction from the collinear part by redefining the collinear field as ξ n → Y n,cs ξ n , A µ n → Y n,cs A µ n Y † n,cs , which is similar to the decoupling of the soft interaction in Eq. (4). Here the csoft Wilson line Y n,cs is defined as and Yn ,cs is obtained by switching n ↔ n. Using the relation n · iD cs = Yn ,cs n · i∂Y † n,cs , the final expression for the PDF is given by Note that this new definition of the PDF is also valid in the region away from threshold. In this case the term i∂ + in the delta function is much smaller than xP + − P + , and it can be safely neglected. Then, due to the unitarity of the csoft Wilson line, the effect of the csoft modes cancels to all orders. Hence we can recover Eq. (6) and describe the PDF away from threshold with the collinear interactions only. Near threshold, we can put the label momentum in Eq. (13) as P + = P + , and obtain the PDF as The same result near threshold has been also derived in SCET + approach [23]. Therefore the expression in Eq. (13) can be regarded as the definition of the PDF over all kinematic regions. Using the expression in Eq. (14), we can calculate the PDF near threshold at the parton level, i.e., φ q/q . As we will see later, the calculation exactly reproduces the standard PDF in the limit x → 1 and it satisfies the Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) evolution. From Eq. (14), the fluctuation of the csoft mode is estimated as p 2 cs ∼ Λ 2 QCD . Therefore this mode scales as p µ cs ∼ Qζ(1, α 2 , α) with α = Λ QCD /(Qζ). When we discuss the decomposition of the collinear and the csoft modes below Eq. (5), the csoft mode can be regarded as a subset of the collinear mode. Hence the collinear mode would scale as Q(1, λ 2 , λ ) with λ ∼ ζ 1/2 α. And the offshellness is given by p 2 c ∼ Λ 2 QCD /ζ, which is much larger than the typical hadronic scale squared. However, near threshold this collinear mode contributes to the PDF only through the virtual corrections without any specific scale. So there is no impact on integrating out the mode ∼ Q(1, λ 2 , λ ) and we can scale it down to p c ∼ Q(1, λ 2 , λ) with λ ∼ Λ QCD /Q. Then the collinear mode at the lower scale has the offshellness p 2 c ∼ Λ 2 QCD .

B. Prescription of the factorization near threshold
With the new definition of the PDF, the naive factorization formulae which are schematically given in Eqs. (1) and (2) can be cast into an appropriate form near threshold. In order to construct correct factorization theorems near threshold, we take the following steps: • After integrating out the hard interactions, we construct the naive factorization formalism by decomposing the collinear and soft interactions.
• Next we decouple the csoft mode from the collinear mode, and express the PDF in terms of the collinear and the csoft fields. 1 • We define the soft functions by subtracting the csoft contributions to avoid double counting.
The derivation of the naive factorization for the structure functions in SCET is not repeated here. Instead we refer to Ref. [14], where the details of the derivation are presented and the naively factorized results are shown in Eqs. (2.13) and (2.33) for DY and DIS respectively.
The naive factorization for the structure function in DY process is written as where N c is the number of colors, J µ is the electromagnetic current, and Q 2 = q 2 is the invariant mass squared of the lepton pair. τ = Q 2 /s and z = Q 2 /ŝ, whereŝ is the centerof-mass energy squared for the incoming partons. The cross section is given by is the Born cross section for the quark flavor f with the electric charge Q f . The naive soft functionS DY is defined as We call this naive soft function since the subtraction of the csoft modes is not included yet. The naive factorization for the structure function F 1 in the Breit frame is given by where P µ = P + n µ /2 is the momentum of the hadron N along the beam and the final-state jet function in the n direction is defined as [9] Xn χn|Xn Xn|χn = / n 2 Here the naive soft functionS DIS is given as Based on the naive factorization in Eq. (15) and (17), we construct the proper factorization theorem as follows: As shown in section. II A, we decompose the collinear and the csoft modes in the collinear sector and express the PDF in Eq. (14). This procedure can be achieved by replacing the collinear PDF f i=q,q/N with the appropriate PDF φ i/N . Then we note that the csoft modes scaling as p n,cs ∼ Qζ(1, α 2 , α) and pn ,cs ∼ Qζ(α 2 , 1, α) are also the subsets of the soft mode ∼ Qζ(1, 1, 1), Therefore, in order to avoid double counting for the overlap region, the csoft contribution should be subtracted from the naive soft function.
Finally the correct factorization theorem for DY and DIS processes are given as Note that there is no correlation between φ q/N 1 and φq /N 2 since the n-collinear(-csoft) and ncollinear(-csoft) fields do not interact with each other at leading order in the power counting of the collinear(-csoft) limit. 2 Here S DY is the soft function for DY after subtracting the csoft contribution from the naive soft function. The subtracting procedure follows the basic idea of the zero-bin subtraction [28], which can be legitimately applied to the csoft sector [29]. As a result the soft function contains neither IR nor rapidity divergence contrary to the naive soft function. In sec. III we will see the details of the computation for the soft function at order α s . When we compare Eq. (21) with Eq. (17), we see thatS DIS is not present in the final factorization theorem since the csoft contribution cancelsS DIS to all orders in α s . At tree level, both S DY in Eq. (20) and the naive soft function are normalized to δ(1 − z). At order α s , S DY is given as S where the superscript denotes the order in α s . And S cs and S cs are the n-csoft and n-csoft contributions respectively, which are given as Here Y n,cs and Y n,cs are the csoft Wilson lines in terms of the n-csoft and n-csoft gluons respectively. Note that 2i∂ 0 in the argument of the naive soft function in Eq. (16) is replaced by n · i∂ (n · i∂) in S cs (S cs ) according to the power counting. In order to specify the soft region completely in Eq. (22), we may introduce and add the contribution of the so-called 'soft-soft (ssoft)' mode scaling as p ss ∼ Qζ(α 2 , α 2 , α 2 ). But the ssoft mode does not contribute to Eq. (22). In general we can divide the full soft region into the 'hard-soft (hsoft)', the csoft, and the ssoft regions. This division and the partition of the full soft region are similar to the procedure in constructing the hard, collinear, and soft modes in SCET from QCD. The only difference is that the large energy Q in the full theory is replaced with Qζ here. Therefore, we can systematically factorize the full soft modes into the hsoft, csoft and ssoft degrees of freedom. In this respect, S DY in Eq. (22) can be regarded as the one-loop correction to the hsoft function obtained from the matching between the full soft and the csoft contributions. Since the csoft modes can reproduce the low energy behavior of the full soft function as the collinear modes do from the full theory, we argue that the hsoft function remains IR finite at higher orders.
For the soft function in DIS, the nonzero csoft contribution only comes from S cs in Eq. (23). And it is the same as the naive soft function, as well as the csoft contribution to the PDF in Eq. (14). At order α s , S

III. THE SOFT FUNCTIONS AND THE PDF NEAR THRESHOLD
We compute the soft functions and the PDF explicitly at one loop in order to verify the statements in the previous section. We regulate the UV divergence using the dimensional regularization with D = 4 − 2 and the MS scheme. We introduce the fictitious gluon mass m g to regulate the IR divergence. We also consider the rapidity divergence [15,16], which appear as the loop momentum k + or k − goes to infinity while k + k − remains finite. We employ the Wilson lines as [16] W n = perms.
exp − g n · P and the rapidity divergence appears as poles in η. The Wilson lines in the n direction can be obtained by switching n and n. For the csoft mode, we use the form of the soft Wilson line, but the soft field is replaced by the csoft field.

A. The soft functions near threshold
We first consider the soft function for DY process near threshold. The naive soft function is defined in Eq. (16), and the correct soft function can be obtained through the csoft subtractions, given by Eq. (22). In order to see the scale dependence clearly, we introduce the dimensionful soft function, which is given as where ω = Q(1 − z). The soft virtual contribution at one loop is given as = α s C F π δ(ω) where k + = n · k and k − = n · k. The real gluon emission at order α s is given as Here Θ is the step function, and we put = η = 0 since the integral has neither the UV nor the rapidity divergence. The final result is expressed in terms of the Λ-distribution. It is defined as where f is a smooth function at ω = 0. In defining the Λ-distribution, Λ is an arbitrary scale larger than ω, but it suffices that Λ is slightly larger than ω.
In obtaining the final result in Eq. (28) with the Λ-distribution, we write M R S,DY as Note that this expression is independent of Λ. But the scale Λ is chosen such that the integration over the delta function δ(ω − k + − k − ) should yield a nonzero value. That means Λ can be any value larger than ω, but from physics Λ is slightly larger than ω, but of the same order, i.e., Λ ∼ Q(1 − z). Fig. 1 shows the phase space for the real gluon emission, where the shaded green region denotes the integration region for the part proportional to δ(ω) in Eq. (28) or (30). The dashed line represents the constraint for nonzero ω. Combining Eqs. (27) and (28), we obtain the naive DY soft function at order α s as We can clearly see that the naive one-loop result contains the IR divergence as the logarithm of m g . This is due to the incomplete cancellation of the virtual and real corrections. If the phase space for the real gluon emission spanned all the region with no constraint, the virtual and the real corrections would cancel, M R S + M V S = 0. In this case, the soft function would become zero at order α s , which holds true to all orders due to the fact that Y † n Y n = Y † n Yn = 1. However, as can be seen in Fig. 1, the phase space for the real gluon emission does not cover all the IR region (near the red line) available to the virtual corrections near threshold. Therefore the incomplete cancellation yields IR divergence in the naive DY soft function. The existence of the IR divergence, which was pointed out in Ref. [13,14], could invalidate the factorization near threshold. Furthermore, as shown in Fig. 1, there is no rapidity divergence in the real gluon emission since the phase space responsible for the rapidity divergence is not included in the phase space for the real gluon emission. The rapidity divergence in the naive soft function comes solely from the virtual contribution.
These problems posed by the naive soft contribution can be resolved by introducing the csoft modes. The csoft contribution is included in the definition of the PDF, but the csoft momentum is also a subset of the soft momentum. In order to avoid double counting, the csoft contribution is subtracted from the naive soft function to define the true soft function near threshold, given by Eq. (22). The subtraction removes both the IR and the rapidity divergences in the soft function.
Let us first consider the contribution of the n-csoft mode at order α s . In DY process, the contribution of the n-csoft mode is the same as the n-csoft case due to the symmetry under n ↔ n. We calculate the dimensionful csoft function S cs (ω) = S cs /Q at order α s . Here ω = Q(1 − z) and the dimensionless csoft function S cs is defined in Eq. (23). The virtual contribution of the csoft mode is the same as that of the soft mode M V cs = M V S , which is presented in Eq. (27).
The real gluon emission at order α s is written as where the result is written by invoking the Λ-distribution, in which the IR divergence as ω → 0 is extracted in the term proportional to δ(ω). In Fig. 2 we show the structure of the phase space for the real emission both in the k + -k − and in the k + -k 2 L planes, where k 2 L ≡ k + k − . The shaded green region is the integration region for the part proportional to δ(ω) and the dashed line is the phase space for M R cs (ω = 0). When the real and the virtual contributions are combined, the contribution of the part proportional to δ(ω) comes from the region B in the k + -k 2 L plane of Fig. 2 since the virtual contribution covers whole region above the line k + k − = m 2 g with the minus sign relative to the real contribution. Therefore the csoft contribution at order α s can be written as Here I B is the integral over the region B in Fig. 2, and it is given as Here the original rapidity regulator |k + − k − | −η is replaced with k −η + . In the region B of Fig. 2, k − in the regulator can be safely ignored since the rapidity divergence occurs only when k + goes to infinity and k − goes to zero while k 2 L = k + k − remains finite. Actually keeping k − in the regulator has no effect on the calculation of the region B in the limit η → 0. Therefore there is no difference between the original regulator and k −η + as far as we integrate over the region B. In the k + -k 2 L plane, I B includes the UV (IR) divergence as k 2 FIG. 2. The structure of the phase space for the real n-csoft gluon emission. The green region is the integration region for the part proportional to δ(ω) with the Λ-distribution. The dotted line k + = ω is the constraint from the delta function for nonzero ω. In the k + -k 2 L plane, the integration region looks simple, and the regions where UV and the rapidity divergences arise are shown.
goes to infinity (m 2 g ). The computation of M R cs (ω = 0) is straightforward. Finally the bare csoft contribution at order α s is given as The n-csoft contribution S cs is the same as S cs due to the symmetry under n ↔ n. We finally obtain the proper soft function for DY process at order α s as Note that this soft function contains only the UV divergence with neither rapidity nor IR divergences. This function is governed by a single scale ω ∼ Λ ∼ Q(1 − z), hence the scale µ to minimize the large logarithms is given by Q (1 − z).
The dimensionful soft function in Eq. (36) can be easily converted to the dimensionless soft function. From the definition of the Λ-distribution in Eq. (29), we have the following relations to the standard plus distribution: Then the renormalized dimensionless soft function to next-to-leading order (NLO) in α s is given as which is the same as W DY in Ref. [13,14]. We now consider the soft function in DIS. The naive soft function has been defined in Eq. (19). Note that the naive soft function is the same as the n-csoft function in Eq. (23) except that the Wilson lines involved are the soft fields for the naive soft function, and the csoft fields for the csoft function. The soft and n-csoft momenta scale as p µ s ∼ Qζ(1, 1, 1) and p µ cs ∼ Qζ(1, α 2 , α). However, since both functions involve the same scale Q(1 − z) ∼ Qζ in the delta functions, they are the same to all orders in α s . 3 Unlike DY process, there is no n-csoft contribution in DIS. When we consider n-csoft contribution from Eq. (19), the delta function becomes δ(1 − z) since Q(1 − z) is much larger than n · p cs in power counting. The n-csoft real emission is the same as the virtual n-csoft contribution but with the opposite sign. Therefore the n-csoft contribution cancels at one loop, and to all orders in α s . As a result, the proper soft function in DIS remains as δ(1 − z) to all orders in α s .

B. The PDF at NLO near threshold
From the definition of the PDF near threshold in Eq. (14), we compute the correction at order α s . If we take the perturbative limit, the PDF at the parton level can be additionally factorized as [5,11,23,30] where Q = p + is the large quark momentum, C q is the collinear part, and S sc is the csoft function defined in Eq. (23). At one loop, the virtual contribution for C q is given as The wave function renormalization and the residue are given by And we obtain the collinear part to NLO as It contains the rapidity, the UV and the IR divergences. From Eqs. (35) and (37), we also obtain the dimensionless csoft function to NLO as Therefore the bare PDF to order α s is given as The PDF near threshold is free of rapidity divergence, and it is the same as the PDF away from threshold when we take the limit x → 1. Obviously, the renormalization-group behavior satisfies the DGLAP evolution in the limit x → 1, and we do not repeat solving the renormalization group equations and refer to Refs. [5,6].

IV. DIHADRON PRODUCTION IN e + e − ANNIHILATION NEAR THRESHOLD
We consider the dihadron production in e + e − annihilation near threshold: e − e + → h 1 + h 2 + X, where X denotes soft particles in the final state. Here the final hadrons h 1 and h 2 in the n and n direction take almost all the energies of the mother partons p + h1 = x 1 Q, where Q is the center-of-mass energy and x 1,2 are close to 1. The scattering cross section is factorized into the two fragmentation functions and the soft function. If we naively compute the soft function without the csoft subtraction, here we also have the IR divergence as well as the rapidity divergence, invalidating the factorization theorem. Therefore we define the fragmentation functions in terms of the collinear and the csoft fields describing the radiations in the directions of the observed hadrons. Then we can properly subtract the csoft interactions from the naive soft function, and as a result the factorization theorem can be written as dσ dp + h1 dp − (45) Here σ 0 is the Born scattering cross section and the threshold region corresponds to x 1,2 → 1. The hard function H DH is given by the same as H DY in Eq. (20). S DH is the soft function for the dihadron production and it can be obtained after the csoft subtraction from the naive soft function.
Following the analysis of the fragmentation function to a jet (FFJ) in the large z limit [21], we define the hadron fragmentation function as where we set the transverse momentum of the hadron p ⊥ h1 as zero. Similarly, the fragmentation function from the antiquark D h 2 /q (z 2 , µ) can be obtained in terms of the n-collinear and n-csoft fields. As in the case of the PDF, this definition is valid near, and away from threshold. Near threshold, putting P † + = p + h 1 we can further simplify Eq. (46) as We can easily check that the result at order α s for Eq. (47) at the parton level is given by the same as the result on the PDF and its renormalization behavior follows the DGLAP evolution in the large z limit. For the soft function, we start with the naive soft function, which is defined as The csoft functions can be obtained by taking the n-and n-csoft limits onS DH , which are subtracted from the naive soft function. As a result the NLO correction to the soft function for dihadron production is given by where the csoft functions S cs and S cs have been defined in Eqs. (23) and (24) respectively. In order to see the scale dependence transparently, we consider the NLO calculation with the dimensionful soft function. It is defined as where ω + = Q(1 − z 1 ) and ω − = Q(1 − z 2 ). Then the naive dimensionful soft function is expressed asS The virtual one-loop contribution toS DH is given by M V S,DH = δ(ω + )δ(ω − )M V S , where the one loop resultM V S is given in Eq. (27). And the real contribution is given as Here we put = η = 0 since the integral has no UV, rapidity divergences.
M R S,DH in Eq. (51) is IR divergent as ω ± → 0. In order to extract the IR divergence, we employ the Λ-distribution defined in Eq. (29), where the upper values for ω ± are set as Λ ± . Then Eq. (51) can be rewritten as Combining the virtual and the real contributions, we obtain the naive soft function at order α s as This result contains the IR and the rapidity divergences as expected. From Eq. (49) the csoft contributions to be subtracted from the naive dimensionful soft function are given by δ(ω − )S (1) cs (ω + ) + δ(ω + )S (1) cs (ω − ). And using Eq. (35) we write the NLO csoft contributions as Finally, subtracting Eq. (54) from Eq. (53) we obtain the NLO result for the bare dimensionful soft function as We clearly see that the problematic IR and rapidity divergences are removed by the csoft subtraction as in DY process. Using Eq. (37), we can convert Eq. (55) to the dimensionless soft function. And the scales for the logarithms are determined by µ ∼ Λ ± . As a result, the renormalized soft function at NLO is given as In the framework of SCET, we have scrutinized the factorization theorems near threshold in DY, DIS processes and in the dihadron production in e + e − annihilation by introducing the csoft modes in SCET. The important point in analyzing these processes near threshold is that there appears a csoft mode governed by the scale ω = Q(1 − z) Q. Near threshold, real collinear particles for the PDF cannot be emitted due to the kinematical constraint, while the csoft modes can. The effect of the csoft modes can be implemented by decoupling the csoft modes from the collinear fields. The resultant PDF consists of the csoft Wilson lines along with the collinear fields. Note that the definitions of the PDF and the fragmentation function in Eqs. (13) and (46) with the csoft modes are valid not only near threshold but also away from threshold. Away from threshold, the effect of the csoft mode is cancelled, and the collinear mode alone describes the whole process.
The naive collinear contribution to the PDF contains the rapidity divergence, but it is cancelled when the contribution of the csoft modes is included. And the resultant PDFs are the PDF obtained from the full QCD. In Ref. [14], the same conclusion has been obtained by reshuffling suitable divergences in the collinear and the soft parts, but here the physics becomes more elaborate. It is also true for the fragmentation function that the definition including the csoft modes yields the fragmentation function in the full QCD, which can be extended to the threshold region.
The csoft mode is also a subset of the soft mode, hence the contribution of the csoft modes should be subtracted in the soft part to avoid double counting. And the soft functions after the subtraction are free of IR and rapidity divergences, and can be handled perturbatively. It is because the phase space for the soft modes responsible for the IR and rapidity divergences coincides with the phase space for the csoft modes. Therefore the IR and rapidity divergences are cancelled when the csoft modes are subtracted from the soft part.
The renormalization group evolution of the soft functions and the PDF can be obtained, for example, as in Refs. [5,6], so we will not repeat it here. But the important point is that the introduction of the csoft modes justifies the use of the renormalization group equation since we have explicitly verified that the higher-order correction to the soft function in DY process (it is zero in DIS) is indeed IR finite and free of rapidity divergences. So far, it has been assumed that the soft functions should be IR finite. But, in fact, the removal of the IR and rapidity divergences is sophisticated and it is accomplished by including the csoft modes. In this paper, we give a firm basis to the use of the renormalization group equation in resumming large logarithms near threshold.
In Refs. [11,12], the authors have considered the same processes, i.e., DY and DIS processes near threshold. However, their results are different from the ones presented here. The main difference is that the PDFs in Refs. [11,12] are formulated by combining the collinear and the soft mode. As a result there arises a correlation between the two collinear sectors in DY process. The prescription in Refs. [11,12] might hold only when Qζ ∼ ω is close to Λ QCD , where the csoft modes become identical to the soft mode. In this paper we set Qζ ∼ ω as a free small energy scale. As far as Qζ ∼ ω Λ QCD , with the help of the csoft modes, we can formulate the factorization theorem where no correlation arises between two collinear and soft sectors.
The final results presented here are the same as those in previous literature [13,14], so we may wonder what can be learned from this analysis. In spite of the same result, there are some illuminating points in our factorization procedure, which are worth commenting. First, the new scale Q(1 − z) is introduced, and it does not have to be related to Λ QCD . Previously it has been considered to accommodate the new scale with the power counting of λ ∼ Λ QCD /Q.
Second, the introduction of the csoft modes yields the soft function free of IR and rapidity divergences. Without the csoft modes, careful analysis shows that the real emission in the soft function contains IR divergence which is not cancelled by the virtual correction. Furthermore, there exists rapidity divergence. We have confirmed that the soft functions are indeed free of IR and rapidity divergences, and the PDF turns out to be the same as the PDF in full QCD.
Third, it becomes apparent which scale governs in each factorized parts. The naive soft function includes the contributions from different scales, but by separating the csoft modes, it now becomes apparent that each part receives the appropriate scale dependence. For the PDF, the scale in the logarithm is of order µ ∼ Λ QCD or larger and for the soft function, µ ∼ ω = Q(1 − z), or Q(1 − x). From the study of the dijet production [19,20], we now know that there should be additional degrees of freedom to account for the resummation with respect to the different scales in different factorized parts. And it is also true in DY and DIS processes, as well as the dihadron production, near threshold.
In conclusion, we have analyzed the factorization near threshold in SCET by including the csoft modes in defining the PDF, and subtracting the csoft contributions for the soft functions. The newly defined PDF can be properly extended to threshold and the resultant factorization theorem takes care of the problem of the IR and rapidity divergences in the soft function, which enables the resummation of the large logarithms through the renormalization group equation. Only after the inclusion of the csoft modes, the factorized result in SCET, Eqs. (20) and (21), is consistent with the result in full QCD.