Non-Global Logarithms, Factorization, and the Soft Substructure of Jets

An outstanding problem in QCD and jet physics is the factorization and resummation of logarithms that arise due to phase space constraints, so-called non-global logarithms (NGLs). In this paper, we show that NGLs can be factorized and resummed down to an unresolved infrared scale by making sufficiently many measurements on a jet or other restricted phase space region. Resummation is accomplished by renormalization group evolution of the objects in the factorization theorem and anomalous dimensions can be calculated to any perturbative accuracy and with any number of colors. To connect with the NGLs of more inclusive measurements, we present a novel perturbative expansion which is controlled by the volume of the allowed phase space for unresolved emissions. Arbitrary accuracy can be obtained by making more and more measurements so to resolve lower and lower scales. We find that even a minimal number of measurements produces agreement with Monte Carlo methods for leading-logarithmic resummation of NGLs at the sub-percent level over the full dynamical range relevant for the Large Hadron Collider. We also discuss other applications of our factorization theorem to soft jet dynamics and how to extend to higher-order accuracy.

C Factorization For Two Strongly Ordered Soft Jets 63

Introduction
A fundamental problem in QCD and collider physics is the identification of hierarchical scales in a system defined by some number of observations made on that system. Generically, ratios of these scales appear in logarithms at every order of the perturbative expansion of the cross section, and can become large in the soft or collinear regions of phase space. To tame these large logarithms and so to improve the convergence of the perturbative expansion, resummation of the logarithms to all orders is required. This resummation requires the factorization of the physics that dominates at each scale from one another, so that one can guarantee that an all-orders description is possible in a particular region of phase space. This program has seen enormous success with applications to predictions for deep inelastic scattering (e.g. [1][2][3][4][5][6]), e + e − (e.g. [7][8][9][10][11][12]), pp and pp collision experiments (e.g. [13][14][15][16][17][18][19][20]), and weak decays (e.g. [21][22][23][24][25]), amongst others.
Strictly speaking the picture outlined above for the resummation of large logarithms is only known to completely capture all logarithms to a given accuracy if all radiation in an event contributes to the observables. Such observables are referred to as "global" if all final state particles contribute to their value. Global observables at particle colliders include thrust [26], angularities [27], or weak boson p T distributions in Drell-Yan production [28], for example. However, global observables are only a subset of observables that are interesting for studying QCD or new physics. Jets, and observables defined on their constituents, rapidity gaps [29][30][31][32], or any observable that is only defined by radiation in a limited region of the full phase space are referred to as "non-global" [33]. Especially with the phase space available for high p T jets at the Large Hadron Collider, there has been substantial effort in defining jet substructure observables [34][35][36], and with the discovery of the Higgs boson [37,38], identifying vector boson fusion events with forward jets is essential for determining properties of the Higgs. Thus, non-global observables are widely used and therefore require detailed theoretical understanding.
Unlike global observables, non-global observables are sensitive to both the relevant scales within the jet or identified phase space region as well as the scale outside the jet. 1 While the out-of-jet region is not directly measured, emissions originating from the outside can contaminate and affect the in-jet region of phase space on which measurements are performed. Therefore, for an accurate description of the cross section to a given logarithmic accuracy of a non-global observable, we must resum not only logarithms of ratios of in-jet scales (global logarithms), but also logarithms of ratios of in-jet to out-of-jet scales. The latter logarithms are referred to as non-global logarithms (NGLs).
While the existence of NGLs has been known for some time, and NGLs have been well studied in the literature [3,33,[39][40][41][42][43][44][45][46][47][48][49][50][51][52][53], they have proved challenging to understand. The leading NGLs in the large N c limit can be resummed by Monte Carlo methods [33,39] or the Banfi-Marchesini-Smye (BMS) equation [40], but a systematic understanding of NGLs to all logarithmic orders is lacking. 2 Other effects, such as finite N c at leading logarithmic order [42,50] and at fixed-order [46-49, 51, 52] have been studied in detail. A systematic understanding of NGLs to all orders has so far not been possible because it has not been shown how to factorize the ratio of in-jet to out-of-jet scales from one another. Therefore, to understand and resum NGLs to arbitrary accuracy requires proving factorization of injet and out-of-jet scales from one another. In this paper, we will do this and present a procedure for systematic improvement of the accuracy to which the NGLs are computed.
Rather than understanding NGLs directly, recently it has been emphasized that the effects of NGLs can be reduced or power-suppressed in some cases [43,[62][63][64]. To do this, one can use jet grooming techniques [64][65][66][67][68] that remove soft, wide-angle radiation in the jet which is (potentially) likely sensitive to out-of-jet scales. While the removal of NGLs using these grooming techniques is indeed one potential solution to the problem (although factorization and resummation for groomed observables is currently also not understood to all orders) here we will attack them directly and work toward their all-orders description.
To understand the problem of NGLs more precisely, consider measuring the mass m and energy E jet of a jet and an inclusive observable on the region outside the jet, such as the out-of-jet energy, E out . Importantly, note that only soft radiation can be sensitive to out-of-jet scales by the collinear safety of the jet finding algorithm, as only soft radiation can cross phase space boundaries. Assuming that soft and collinear physics factorizes, we introduce a soft function S which encodes the pattern of soft radiation from the dipoles (pairs of eikonal Wilson lines) in the event. For a global observable, the soft function only depends on scales set directly by the measurements on the event. However, in this case, the soft function is sensitive to in-jet and out-of-jet scales [46][47][48]: S ≡ S (m, E jet , E out ) . (1.1) No measurement has been done on the jet to determine if the emissions setting the mass come from in the jet or outside the jet, and so it would seem like this soft function cannot be further factorized. Therefore, without isolating the in-jet and out-of-jet scales in some way, the NGLs of this system cannot be resummed. In the case just discussed, the soft function could not be further factorized because no measurement was done on the jet to isolate the region of phase space where the NGLs are important. For the jet mass, this region of phase space is the emission of a single soft gluon from outside the jet into the jet, near the jet boundary. A single gluon is not an infrared and collinear (IRC) safe quantity, so we should really think of this as a soft subjet located near the jet boundary. 3 Nevertheless, just measuring the jet mass is not sufficient to uniquely specify this region of phase space. The resolution to the problem of NGLs is therefore clear: we must measure sufficiently many observables on the jet so as to isolate the region of phase space in which the NGLs are important. This is the goal of this paper. By measuring several observables on a jet we are able to identify the region of phase space where the NGLs live, refactorize the soft function, and resum the NGLs by renormalization group evolution of the now-factorized soft function.

Overview of the paper
Here, we present a detailed overview of the content and reasoning of this paper so that our logic is not lost to the details of calculation. For simplicity, in this paper we restrict ourselves to jet production in e + e − collisions, although the approach can be extended to other situations. Our approach in this paper to accomplishing the resummation of NGLs is the following. First, we find a jet in an e + e − → hadrons collision event by identifying the broadening axis [71][72][73] of the event and including those particles that lie within a cone of fixed radius R of the broadening axis. For much of the jet's phase space relevant for NGLs, this is identical to finding jets with the anti-k T jet algorithm [74] with radius R and the Winner-Take-All (WTA) recombination scheme [73,75,76]. The broadening axis is insensitive to recoil effects on the jet axis, and so the jet axis aligns with the direction of the hardest radiation in the jet. This is necessary to eliminate back-reaction on the jet direction from wide-angle, soft radiation. For the region of phase space outside the jet, we measure some quantity, which we refer to as B. This measurement sets the scale of the out-of-jet radiation and we require B 1 which enforces soft and collinear dynamics to dominate the out-of-jet region. Additionally, we will assume that the out-of-jet scale is much lower than the in-jet scale which is the phase space region in which the NGLs are large and must be resummed.
Within the jet, we want to guarantee that the jet contains a soft subjet approaching the jet boundary which is sensitive to the out-of-jet scale B. We do this by measuring several IRC safe n-point energy correlation functions on the jet [77,78]. In particular, to uniquely identify a single soft subjet and determine its energy fraction and angle from the hard jet core, we measure the two-and three-point energy correlation functions e (α) 2 , e (β) 2 and e (α) 3 , for angular exponents α, β. The measurement of three energy correlation functions is required to enforce that the soft subjet is not collinear with the hard jet core. The soft subjet region of this three-dimensional phase space is parametrically defined by [78] e where we have assumed that α > β. This soft subjet region is an essential component of the description of the full phase space defined by the energy correlation functions and is required to make predictions of distributions of jet discrimination observables such as D 2 , defined in Ref. [78]. The phase space formed from the simultaneous measurement of the two-and three-point energy correlation functions is described in Ref. [78], and we will present analytic predictions for the full double-differential cross section in a future publication [79]. Identification of the soft subjet region of phase space enables a factorization of the cross section in this region of phase space by identification of the dominant modes that contribute to the particular values of the measured energy correlation functions. Assuming that B is an additive observable and B e , and S n sjnsj (e (α) 3 ; R), each of which will be discussed in more detail in Sec. 3. In Fig. 1a, we present an illustration of the modes that contribute in the soft subjet region of phase space. Note that our factorization theorem in the soft subjet region factorizes the in-jet scales defined by e from the out-of-jet scale B, and therefore the NGLs of ratios of the soft subjet energy to the out-of-jet scale B can be resummed. In the effective field theory language, the additional measurement has converted one of the soft scales to a hard scale, allowing for the resummation of the NGLs by standard renormalization group techniques.
By studying the dynamics of the soft subjet in this region of phase space, we are led to introduce what we term the dressed gluon approximation, 4 which captures the resummation of NGLs due to unresolved emissions associated with the soft subjet. 5 In particular, in the region of phase space with a single soft subjet, we have one-dressed gluon. We demonstrate that the dressed gluon approximation can be used to calculate the NGLs of a more inclusive observable by marginalizing over the factorization theorem. Importantly, the one-dressed gluon approximation can be calculated to arbitrary perturbative accuracy and for any number of colors, N c .
While the measurement of e has allowed us to successfully factorize the jet scales set by e (α) 2 and e (β) 2 from B, because we have not resolved any further emissions 4 The term "dressed gluon" is also used in approaches to renormalon resummation (see e.g. [80]). These approaches, although similar in the spirit of associating additional dynamics with a single gluon, attempt to describe completely distinct physical phenomena. 5 A similar approximation was used in calculation of jets with rapidity gaps [81][82][83][84], termed the "out-ofgap" expansion. The connection of these dressed gluon expansions to physically measurable subprocesses, as well as their realization as an expansion of the BMS equation and the role of the buffer region was not addressed in these works.  Figure 1: (a) Illustration of the phase space configuration and dominant modes for a jet containing a hard core and a soft subjet. Here the gray radiation denotes global soft radiation S nnn sj , and the green radiation denotes collinear radiation along the direction of the energetic jet axes, J n and Jn. The soft subjet dynamics is described by soft jet modes, J n sj shown in blue, and boundary soft modes shown in red, S n sjnsj . (b) Schematic of the ladder of factorization theorems defined by increasingly differential measurements made on the jet. With each additional measurement, the NGLs are pushed to the soft function at a lower unresolved scale. The S , and S are schematic, typically being a product of multiple functions, but depend only on a single scale.
in the jet, there still exist NGLs at the scale of e (α) 3 . To resum those NGLs requires resolving two soft subjets located near the jet boundary. To isolate this region of phase space and resum the NGLs including e (α) 3 , we can measure the two-, three-and four-point energy correlation functions, factorize the cross section and renormalize, which produces the two-dressed gluon approximation. This then pushes the NGLs to the unresolved scale set by the four-point energy correlation function. The procedure can then be repeated, by measuring higher-point energy correlation functions to resolve more and more soft subjets in the jet, resumming the NGLs down to some unresolved scale below which no soft subjets are identified. We illustrate this increasingly differential factorization theorem ladder for resumming NGLs in Fig. 1b. We discuss the convergence of the reorganization of the traditional perturbative expansion in terms of the number of dressed gluons, and show that the contribution from higher numbers of dressed gluons is highly suppressed by the available phase space volume. We also relate the dressed gluon expansion to an expansion of the BMS equation. We stress that these factorization theorems can be calculated to any perturbative accuracy and for arbitrary numbers of colors N c , allowing for the extension to the resummation of subleading logarithmic corrections.
To justify that this step-by-step resummation procedure of the NGLs accurately captures those NGLs known to exist, we compare our dressed gluon approximation of hemisphere jet masses in the large-N c limit to Monte Carlo resummation, and the fixed-order expansion of the BMS equation. By only including the one-and two-dressed gluon approximations, we find agreement with solutions of the BMS equation at the sub-percent level for phenomenological values of the NGLs. This demonstrates that the one-and twodressed gluon approximations capture the dominant contributions to the leading NGLs in the large-N c limit, with small corrections due to the presence of NGLs at lower resolution scales. The dressed gluon approximation is easily incorporated analytically into existing factorization theorems for multi-jet processes. 6 The dressed gluon approximation also provides analytical understanding of many features of jet physics and NGLs; for example, dressed gluons manifest the "buffer region" [39] near the jet boundary in which emissions are forbidden.
The outline of this paper is as follows. In Sec. 2 we define the soft subjet region of phase space via measurements of e 2 and e 3 . In Sec. 3, we present the factorization of the cross section in the soft subjet phase space region within the context of soft-collinear effective theory, calculating anomalous dimensions and renormalizing the functions appearing in the cross section. In Sec. 4, we introduce the dressed gluon approximation, which follows from our factorization theorem in the soft subjet region of phase space. We demonstrate how the dressed gluon approximation can be used in the calculation of NGLs by computing the NGLs for the hemisphere invariant mass with both one and two dressed gluons, and provide a numerical comparison to the BMS equation, and various approximations found in the literature. We also discuss analytic insights into the dynamics of NGLs which are realized in the dressed gluon approximation. In Sec. 5, we discuss how our approach can be extended beyond leading logarithmic accuracy and discuss the necessity of capturing NGLs arising from collinear splitting along the boundary of phase space. We conclude in Sec. 6 and discuss directions for further understanding of NGLs and soft subjet dynamics to all-orders. Calculational details are presented in appendices.

Observables and Phase Space
As discussed in the introduction, our strategy for resumming NGLs is to isolate the region of jet phase space which is sensitive to both in-jet and out-of-jet scales, using IRC safe measurements for which we can prove a factorization theorem. We will consider the process e + e − → hadrons on which we apply a broadening axis cone algorithm. 7 For our purposes, the implementation of this algorithm is to identify the broadening axes in the event, which can be done in an inclusive fashion by minimizing 2-jettiness [71,72,85] with angular exponent β = 1, draw fixed cones of radius R around the axes, and study the largest energy jet, defining the out-of-jet region to be the complement of the cone of radius R. We will often take R = π/2, in which case this divides the event into hemispheres. The choice of jet algorithm is vital for simplifying the analysis of NGLs; the use of the broadening axes and geometric cones guarantees the jet axis lies along the direction of the hardest radiation in the jet with a circular shape. This remains true for the complete set of possible soft 6 This assumes that the jet definitions in those factorization theorems are robust to soft subjets approaching their boundaries. 7 We thank Jesse Thaler and Daniele Bertolini for discussions about such algorithms.
subjets, even when the soft subjet lies near or on the jet boundary. This scheme is largely equivalent to finding the jets with the anti-k T jet algorithm [74] with the Winner-Take-All (WTA) recombination scheme [73,75,76]. For jets containing soft subjets, the two algorithms will give identical factorizations at leading power for much of the soft subjet phase space. However, because it is a sequential recombination algorithm, the boundary of anti-k T jets will distort in the presence of soft subjets located sufficiently near the boundary. The power counting of the factorization theorems give a precise definition of the phase space boundaries between these regimes. At lowest order, the region of phase space sensitive to both in-jet and out-of-jet scales consists of a single soft gluon in the jet, located near the jet boundary. A single gluon is not an IRC safe object, and so the natural IRC safe generalization of a soft gluon is a soft subjet. Therefore, we wish to isolate jets which have a hard core of radiation and a soft subjet at large angle from the jet core. Because we use a recoil-free jet algorithm, this soft, wide angle subjet does not displace the jet axis from the hard jet core.
To identify the soft subjet region of the jet's phase space, we will measure a number of the n-point energy correlation functions [77,78] on the jet. The n-point energy correlation function is an IRC safe observable and is sensitive to n-prong structure in a jet. To identify the soft subjet region of phase space we will need the two-and three-point energy correlation functions which we define for e + e − collisions as [77] e (β) where J represents the jet, E i and p i are the energy and four momentum of particle i in the jet J and β is an angular exponent that is required to be greater than 0 for IRC safety. The four-point and higher energy correlation functions are defined as the natural generalization. The n-point energy correlation function vanishes in all soft and collinear limits of an n particle configuration. In the soft or collinear limit, the e (β) 2 are equivalent to the (recoil-free) angularities [27,73,86,87] and when measured event-wide in e + e − collisions, e 3 is equivalent to the C-parameter [88,89] computed to O(α s ). The measurements that we perform on these events are as follows. First, we measure some observable B on the region outside the identified jet J. On the jet J, we measure the observables e for angular exponents α, β and we will assume that α > β. To demand that the dynamics of the jet are dominated by soft and collinear radiation, we require that e (α) 2 1. A jet with a hard core and a single soft subjet has 2-prong substructure, and so to identify 2-prong jets we require [78] which follows straightforwardly from power counting. The measurement of the two-and three-point energy correlation functions on a jet resolve either 1-or 2-prong substructure, which is described in detail in Ref. [78] and is displayed in Fig. 2a. Importantly, the energy correlation functions provide a parametric separation of the 1-and 2-prong regions of phase space, defined by the precise scaling relation of Eq. (2.2). Due to this parametric separation, well defined factorization theorems exist in both regions of phase space. To enforce that the subjet is both soft and at a wide angle from the jet core, we therefore require e (α) Combined with the condition e (α) 2 1, this forces z sj 1 and θ sj ∼ R, where R is the jet radius, which we assume to be an order 1 number. In particular, this additional measurement allows us to distinguish the case of a soft subjet from the case of two energetic collinear subjets, which has been studied in Ref. [90]. The e (α) 2 , e (β) 2 phase space is described in detail in Ref. [91] (see also [92]) and displayed in Fig. 2b. The two different subjet configurations, which exist on the boundaries of the allowed phase space defined by e  and e (β) 2 , the energy fraction z sj and angle from the jet core θ sj of the soft subjet are well-defined and IRC safe quantities. Once the soft subjet region of phase space has been identified using the energy correlation functions, the observables we will consider as measured on the jet are z sj , θ sj and e (α) 3 . Fig. 1a illustrates the structure of the event we are studying and the observables that we measure on the in-jet and out-of-jet regions.  Figure 3: The two distinct subjet configurations which exist in the two prong region of phase space. (a) Two energetic collinear subjets, which has been studied in Ref. [90], and populates the region of phase space e

Effective Field Theory Description and Factorization
In this section we present a factorization theorem in the soft subjet region of phase space described in Sec. 2. We use the formalism of soft-collinear effective theory (SCET) [93][94][95][96], an effective theory of QCD in the soft and collinear limits. Because the factorization theorem involves many novel features, we will discuss its structure in detail. We begin with a power counting analysis in Sec. 3.1 to determine the modes required in the low energy effective theory. The mode structure dictates the functions appearing in the factorization theorem, which is presented in Sec. 3.2. The novel feature of the factorization theorem is the presence of modes whose virtuality is set by their angle to the boundary of the jet and whose resummation is directly tied to the resummation of NGLs.

Modes of the Factorization
The simultaneous measurement of e defines a multi-scale structure in the low energy effective theory required for the complete factorized description of the soft subjet region of phase space. A proper understanding of the relative scalings of the modes is essential to specify the structure of the factorization theorem we will present in Eq. (3.10), and to understand the structure of zero bin subtractions [97] needed to remove overlaps. The scaling of the modes in the low energy effective theory can be determined by power counting arguments, which have been considered for multi-differential observables resolving multiple subjets in Refs. [78,98]. We will follow the conventions used in those papers, writing the scalings of all modes in terms of the physical observables, instead of the traditional λ, as we find this to be more transparent. Typically, the modes of a factorization theorem are set by the observables measured. We will find in addition that the existence of the jet boundary itself will play a vital role in determining all the modes required for the complete factorized description.
Because we only measure the observable B in the out-of-jet region, its mode structure is simple and so we will discuss it first. Importantly, we require that B is an additive observable, like the energy correlation functions, angularities, mass, etc., so that the factorization theorem is of a universal form. We assume that the out-of-jet region has order-1 angular size and B 1 and so the measured value of B is dominated by collinear and global soft modes. Given a momentum p, we will adopt the following notation for its components expressed in light-cone coordinates defined by the directions a andā: As an example, if we assume that B measures the out-of-jet broadening, the scaling of the collinear and soft modes is where Q is a proxy for the energy in the out-of-jet region and is of the same order as the total scattering energy in the event. n is the direction of the jet of interest, whilen is the direction of the jet that contributes to the measured value of B. Now, we turn to the in-jet modes that contribute to the measurement of the energy correlation functions. As shown schematically in Fig. 3, in the soft subjet region of phase space, we assign the power counting θ sj ∼ R ∼ 1, where θ sj is the angle between the jet axis and the soft subjet direction. In this region of phase space, the two-point energy correlation functions are dominated by soft, wide angle radiation and so we have where z sj is the energy fraction of the soft subjet. As there is no parametric difference between e (α) 2 and e (β) 2 , these observables are redundant from a power-counting perspective and in the following, we will express the scaling of all modes in terms of e on the jet resolves the hard collinear structure of the jet core, the structure of the soft subjet, and global soft radiation from the hard dipoles present in the event. 8 We will denote the angular size of the hard core by θ c , the angular size of the soft subjet by θ cs and the energy fraction of the global soft radiation by z s . Therefore, the parametric scaling of e (α) 3 in the soft subjet region of phase space is set by these three contributions: e (α) 3 ∼ z sj (θ α c + z sj θ α cs + z s ) .
where E J is the energy of the jet and n andn are the light-like directions of the jet of interest and the other jet in the event, respectively. The soft subjet mode's momentum scales like in the light-cone coordinates defined by the direction of the soft subjet, n sj . These are the complete set of modes defined by the scales set by the measurements of e alone. If these measurements were global, that is, sensitive to all radiation in the event, this would be the complete enumeration of the modes that contribute to the measured observables. However, because the energy correlation functions are only measured on the jet, the boundary of the jet plays an important role in defining the relevant modes as well.
In particular, since we are considering the case where the out-of-jet scale is much less than the in-jet scale, namely B e (α) 2 , 9 for the modes in the soft subjet, the angle between the soft subjet axis and the jet boundary ∆θ sj ≡ R − θ sj places additional constraints on the soft subjet dynamics, much like an additional measurement would. In the region of phase space in which the NGLs are parametrically large, and should be resummed, this must be taken into account. In Sec. 5 we will briefly discuss the case when the NGLs are not parametrically large, and how this factorization theorem is modified. We therefore expect that ∆θ sj defines a relevant scale in the effective theory, and should be included in the power counting analysis. We must consider the possibility of a mode whose angular scale with respect to the soft subjet axis is not set by the measurement of the two-and three-point energy correlation functions, but rather by the jet boundary itself. This mode does not contribute to the two-point energy correlation functions e  3 . This therefore defines an additional soft mode which is localized around the soft subjet's direction and constrained by the jet boundary. We therefore refer to this new mode as a boundary soft mode. 9 Formally we take the scaling B ∼ The presence of this mode is absolutely necessary for understanding NGLs. Importantly, within the fat jet, it is the only mode that contributes to the cross section singular logarithmic terms of the form The necessity of this mode, justified here by power counting, also appears from explicit calculation of the functions appearing in the factorization theorem in Eq. (3.10). Logarithms of ∆θ sj can also arise from global soft radiation in the out of jet region (see App. B.5), and so the boundary soft mode will be required for renormalization group consistency of the factorization theorem. It is critical for the NGL resummation that the two different factorized functions at different energy scales are both sensitive to the jet boundary. The scaling of the momentum of the boundary soft mode is determined by considering its contribution to e (α) 3 , given that its angular scale is set by ∆θ sj . The dominant contribution to e (α) 3 from the boundary soft modes is where the energy fraction of the boundary soft radiation is z bs . Therefore the boundary soft mode's momentum components scale like written in the light-cone coordinates defined by the soft subjet axis. For consistency of the factorization, we must enforce that the soft subjet modes cannot resolve the jet boundary and that the boundary soft modes are localized near the jet boundary. That is, the angular size of the soft subjet modes must be parametrically smaller than that of the boundary soft modes: Therefore, the factorization theorem applies in a region of the phase space where the soft subjet is becoming pinched against the boundary of the jet, but lies far enough away that the modes of the soft subjet do not touch the boundary. 10

The Factorization Theorem for a Soft Subjet
We will prove the following factorization theorem for the production of a soft subjet:

10)
10 For the anti-kT jet algorithm, if the modes of the soft subjet resolve the jet boundary (so ∆θsj ∼ θcs), this will in general result in distortions of the jet boundary due to clustering effects. In a strict cone algorithm, as we use here, this boundary-collinear regime is still factorizable, and is relevant for the resummation of subleading NGLs, as will be discussed in Sec. 5. valid under the assumptions on the phase space described in Sec. 3.1. Here convolutions are implicit in any variable that is twice repeated, and we have explicitly indicated the dependence on the jet boundaries with the jet radius R. The operator definitions of the functions are given in App. A, but the physical origin of each function is straightforward to understand, with each function describing the dynamics of one of the modes described in Sec. 3.1. A brief description of the functions appearing in Eq. (3.10) is as follows: • H(Q 2 ) is the hard function for the production of a dijet event at an e + e − collider.
• H sj nn e • J n e • Jn(B) is a jet function describing the collinear modes of the out-of-jet region of the event.
• S nnn sj e (α) 3 ; B; R is the global soft function, involving three Wilson line directions, n,n, n sj .
is a jet function describing the dynamics of the soft subjet modes, which carry the bulk of the energy in the soft subjet.
• S n sjnsj (e  The calculation to one-loop of various objects in this factorization theorem is presented in App. B. We now describe how the factorization theorem of Eq. (3.10) arises in a multi-stage matching onto the effective theory involving the modes of Sec. 3.1. 11 This multi-stage matching procedure is shown schematically in Fig. 4. In the first step, QCD is matched onto an SCET theory at the hard scale of the e + e − event by matching the electroweak currents 11 In the case that zsj B, such a multi-stage matching is not necessary and one can construct the factorization theorem along the lines of Refs. [85,87]. However, here we pursue the multi-stage matching since it will be necessary in a forthcoming publication [99] dealing with the factorization of soft subjets when zsj B.  Figure 4: A schematic of the multi-stage matching procedure used to prove the factorization theorem of Eq. (3.10) for the soft subjet region of phase space. As discussed in the text, the matching proceeds in three stages: matching to SCET, refactorizing the soft function to describe the soft jet production, and factorizing the boundary soft mode. The canonical scales of the modes in the final factorization theorem are shown on the right, ordered in virtuality. Here we have chosen an angular exponent α = 2 for concreteness. of both theories [100][101][102]. This theory is then evolved down to the soft scales defined by e (α) 2 ∼ e (β) 2 , where these soft modes can be decoupled via the BPS field redefinition [96] from the collinear modes. The soft scale associated with this measurement is e (α) 2 ∼ z sj , so the soft subjet is not resolved, and is simply described as part of the soft function. At this stage we have the usual SCET factorization formula for dijet production: where we have decoupled the collinear modes in the jet functions from the soft radiation at this scale using a BPS field redefinition. We have not yet performed the full multipole expansion required to separate all infrared scales [103,104], since from the power-counting arguments given above we have not included enough modes to separate all scales. Hence we let all modes resolved at this scale contribute to the jet measurements.
As we lower the scale below e (α) 2 , we resolve the soft subjet, and therefore we must further expand the soft function. However, since the collinear dynamics of the jets along the axes n andn are at the scale e (α) 3 and B, respectively, and they have been decoupled from the soft radiation, their scales can simply be lowered without any matching. 12 This should be contrasted with the factorization theorem of Ref. [90] for the case of two collinear subjets, where additional matching must be performed in the jet function. In our case, an additional matching step must instead be performed in the soft function, corresponding to the known fact that the NGLs appear in the soft function. At the scale e (α) 2 , we expand the soft function and match onto a hard function H sj nn describing the production of a soft subjet, a jet functionJ describing the dynamics of the soft subjet, and a global soft function involving three Wilson lines. This refactorization is shown schematically as the transition from the first to the second column in Fig. 4. This is analogous to the construction of the PDFs in SCET, see Ref. [96]. The factorization theorem then becomes In the final step of the matching, we lower to the scale set by the measurement of e (α) 3 . At this scale we probe the dynamics of the soft subjet, and must perform a final factorization ofJ into a function describing the boundary soft modes and a function describing the soft subjet modes:J (3.14) The boundary soft modes can be BPS decoupled from the soft subjet modes, thus resulting in the final form of the factorization theorem in Eq. (3.10). The virtualities of all the modes present in the final factorization theorem are listed in Fig. 4 for the specific case of α = 2. The final factorization of the soft subjet function in Eq. (3.14) is essential to resum all logarithms. In particular, the soft subjet modes, described by J n sj e (α) 3 , are not sensitive to the jet boundary, as should be expected for a collinear mode, whereas the boundary soft modes are sensitive to the jet boundary. Therefore, S n sjnsj (e (α) 3 ; R) will have a logarithmic dependence on the angular distance of the soft subjet to the boundary ∆θ sj , diverging as the soft subjet approaches the boundary. The same logarithmic dependence of the boundary angle is found in the naïve global soft function S nnn sj e (α) 3 ; B; R both in the injet and out-of-jet regions of phase space. To avoid double counting, one must subtract the contribution of the boundary soft region of phase space from the global soft function (this subtraction is implemented in the effective theory via a zero bin subtraction [97]), which in turn removes the logarithmic dependence on ∆θ sj from the global soft contribution to the identified jet. However, this dependence on ∆θ sj will exist in the boundary softs and the global soft radiation in the out-of-jet region. This will be important for the resummation of the NGLs associated with the soft subjet production by running from the boundary soft scale down to the out-of-jet scale.
In summary, we have presented a factorization theorem describing the region of phase space in which a soft subjet is identified within a jet. By performing multiple measurements 12 At this stage, the collinear modes cannot contribute to e to isolate a specific region of phase space, we were able to refactorize the multi-scale soft function. Effectively, the additional measurement converted soft scales to hard scales, so that logarithms associated with these ratios of scales can be resummed by standard renormalization group evolution, hence achieving the resummation of NGLs in this particular region of phase space. While this factorization theorem is interesting in its own right for studying the dynamics of the soft subjet, as relevant for example for the factorized description of jet substructure variables, we will not pursue this direction further in this paper, leaving it to a future publication [79]. Instead, in this paper we will focus on using the factorization theorem in this region of phase space to understand properties of NGLs, and applying this understanding to the resummation of NGLs for more inclusive observables.

Dressed Gluon Approximation
In the previous section, we have shown how multiple measurements can be used to isolate a region of phase space involving a single soft subjet, and we presented a factorization theorem, Eqs. (3.10) and (3.11), in the framework of SCET describing this region of phase space. By making multiple measurements, we are able to refactorize the soft function allowing for the resummation NGLs of e (α) 2 /B in this particular region of phase space. In this section, we discuss how the resummation of the NGLs in the soft subjet region of phase space allows us to understand more general properties of NGLs. We propose a reorganization of the perturbative expansion for calculating NGLs of more inclusive observables. We call this the dressed gluon approximation, which we define in Sec. 4.1. We then explicitly demonstrate how the dressed gluon approximation can be used in a calculation, by calculating the one-and two-dressed gluon contribution to the NGLs for the hemisphere jet mass in Sec. 4.2. In this section, we also develop an understanding for the emission of a soft subjet off of N eikonal lines, and present a conjectured factorization theorem. In Sec. 4.3 we compare our dressed gluon approximations to both resummed and fixed-order solutions of the BMS equation, which captures the leading NGLs in the large N c limit. Using the properties of the dressed gluon approximation, in Sec. 4.4 we discuss some insights into features of NGLs, and compare the expansion in the number of dressed gluons to other expansions of the BMS equation.

4.1
Putting the Pants on a Gluon One Leg at a Time The dressed gluon approximation originates from the observation that the factorization theorem for the soft subjet resums a set of higher order corrections to the matrix element for the production of a single soft gluon from the nn dipole. The matrix element for gluon emission from the dipole is given by the hard function, H sj nn (z sj , θ sj ), in the factorization theorem of Eq. (3.11). These corrections include, however, more than just the threshold virtual corrections of the soft gluon current, as they also contain an arbitrary number of soft gluon emissions in the out-of-jet region. Therefore, these are precisely the corrections associated with the NGLs.
However, our factorization theorem is multi-differential, and so to determine the NGLs for a more inclusive measurement requires marginalizing or integrating over observables that are not included in the measurement. This integration can potentially be problematic from the point of view of our factorization theorem. We have only resummed NGLs down to the scale set by e (α) 3 , but for an inclusive measurement we must integrate over all regions of phase space, including scales lower than e (α) 3 , beyond which no jet structure is resolved. This would seem to indicate a loss of formal accuracy, as there may be large logarithms the factorization theorem is ignorant to that should be resummed. Indeed, traditional strict logarithm counting is breaking down in this case, as there would be terms scaling like α s ln ∼ 1 that are not resummed by our factorization theorem.
The resolution of this is the realization that resummation of NGLs in our factorization theorem down to the scale e (α) 3 exponentially suppresses the phase space for emissions at scales lower than e (α) 3 . This feature of the expansion in the number of dressed gluons will be discussed in detail in Sec. 4.4.3, after we have shown that the dressed gluon approximation gives rise to the buffer region [39] in Sec. 4.4.1. Therefore, while a formal logarithmic accuracy of the resummation of NGLs may not necessarily be guaranteed, effects from emissions at unresolved scales are suppressed by their allowed phase space. We refer to the resummation of a soft subjet according to the factorization theorem of Eq. (3.11) as the dressed gluon approximation. 13 Resummation in the factorization theorem dresses the soft subjet at a scale defined by e by an arbitrary number of soft emissions down to a scale set by e (α) 3 . By the structure of the factorization theorem, the one-dressed gluon approximation is guaranteed to include the correct NGL at O(α 2 s ), and resums a tower of NGLs at higher orders in α s . The accuracy of this approximation is controlled by the volume of allowed phase space for emissions at scales lower than e (α) 3 , so the one-dressed gluon approximation does not fully include, for example, the NGL at O(α 3 s ). To fully describe this NGL, we must include the two-dressed gluon approximation, by resolving emissions down to a lower scale in the factorization theorem with further measurements. This then correctly describes the NGL at O(α 3 s ), and resums a tower of NGLs to higher orders in α s . One can continue the procedure, adding more and more dressed gluons to obtain an arbitrarily accurate description of the NGLs. However, as we will demonstrate, the phase space suppression accompanied by an increased number of dressed gluons causes the dressed gluon expansion to converge rapidly, so that for phenomenologically relevant values of the NGLs, only one or two dressed gluons are required for an accurate description. We will explicitly consider up to two dressed gluons in this paper.
We will first present a detailed explanation of the one-dressed gluon approximation and its construction from the factorization theorem of Eq. (3.11), before discussing the extension to multiple dressed gluons. Because we are interested in the NGLs for more inclusive jet measurements, we must integrate over the unresolved scale set by e (β) 3 . We will do this by Laplace transforming the multi-differential cross section for the soft subjet phase space region as: .
(4.1) 13 We hope that this nomenclature explains the unusual title of this section.
The cross section fully inclusive over e (α) 3 is then found by the limitẽ (α) 3 → 0. This limit gives a prediction for a single soft gluon matrix element, with all possible low energy unresolved configurations produced by its subsequent splittings.
The limitẽ → 0 of the cross section defined in Eq. (4.1) is formally singular at any fixed-order and is only regulated when resummed to all-orders. In particular, the fixed-order anomalous dimensions present in the factorized form of Eq. (4.1) are singular asẽ (α) 3 → 0, which makes their interpretation challenging. However, by reorganizing the functions present in the factorization theorem into in-jet and out-of-jet contributions, all dependence on the observableẽ can be controlled and finite anomalous dimensions can be identified in the limitẽ (α) 3 → 0. Specifically, we rewrite the factorization theorem in the suggestive form: where we have refactorized the global soft function along the lines of Refs. [46,47], separating out the global logarithms ofẽ and B that are resummable within the factorization theorem. All convolutions are implicit. The (in) and (out) labels denote the in-jet and out-of-jet phase space regions, and (NG) denotes the non-global contributions. For compactness, we have used the notation of Eq. (3.14), whereJ n sj is the unfactorized soft subjet function, which contains both the boundary soft and jet modes of the soft subjet.
Using the renormalization group equation that resums the global logarithms ofẽ nnn sj . We also use the shorthand notation S (out+NG) nnn sj for the product of the out-of-jet and non-global soft functions. In Eq. (4.2), we have removed any global contribution to B from the initial nn dipole by the appropriate global part of the soft function S nn (ẽ (α) 3 ; B; R). This soft function has a similar factorization: where we have explicitly included dependence on the renormalization scale µ. The first factor in Eq. (4.2), This fact has deep consequences. We now introduce the two functions that define the dressed gluon's factorization theorem: , (4.5) which are the second and third factors, respectively, in square brackets in Eq. (4.2). By the renormalization group invariance of the total cross section and the first factor of Eq. (4.2), the product W nn G nnn sj must also be renormalization group invariant. That is, these functions have the renormalization group equations: where γ D is the anomalous dimension, which is given to one-loop in App. B.7. The resummed dressed gluon with one-loop anomalous dimensions is where the scale at which G nnn sj is evaluated has been set to µ i . Taking the tree-level matrix-elements then gives: Note that the dressed gluon's matrix element vanishes as it approaches the jet boundary, where θ sj → R. Therefore, emissions are suppressed near the jet boundary corresponding to the buffer region identified in Monte Carlo simulations of NGLs [39].

Calculating with a Dressed Gluon
From the suggestive form of Eq. (4.2), we are able to define a generic procedure for incorporating non-global effects into the resummation of an arbitrary additive observable measured on a jet or other restricted phase space region. 14 As a concrete example, we will 14 Additivity of the energy correlation functions and the out-of-jet observable B was necessary for the original form of the factorization theorem and its rewriting in Eq. (4.2). However, we strongly suspect that the resummation of NGLs for non-additive observables, such as the fractional jet multiplicity [105] defined with the jets-without-jets algorithm [75], can be accomplished by extending the methods discussed here. We thank Jesse Thaler for discussions on this point.
use the dressed gluon approximation, with one and two dressed gluons, to include nonglobal effects in the factorization theorem for the hemisphere mass observables in e + e − collisions. This will be sufficient to clearly illustrate how the procedure can be extended to an arbitrary number of dressed gluons, and for an arbitrary additive observable. We will denote the mass of the left (right) hemisphere as m L (m R ), and consider the cumulative cross section defined as where σ 0 is the Born cross section and is the double differential cross section of the hemisphere masses. Non-global effects only appear in the soft function, so, in this section, we will demonstrate how to include the oneand two-dressed gluon approximations into the soft function for hemisphere mass measurements, with the generalization to other non-global measurements being straightforward.
We emphasize that the use of the dressed gluon will depend on the relation between m L and m R . It is therefore important to clarify the notation used in describing the soft subjet factorization theorem compared with that used for the left-right hemisphere mass distribution. When we presented the soft subjet factorization theorem, we referred to an in-jet region, in which we were multi-differential, and an out-of-jet region, in which a single observable B was measured. Importantly, as discussed in Sec. 3.1, the factorization theorem presented is valid in the case that the out-of-jet scale is lower than the in-jet scale. 15 To make the correspondence with the left-right hemisphere mass distribution, we must make an assumption for the relation between m L and m R . For the case of the hemisphere masses considered here, this is somewhat of a trivial point since the distinction is just a matter of relabeling; however, for general geometries it is important. In this section, we will therefore write explicit expressions for the case m L > m R , and simply indicate with R ↔ L the opposite case. For m L > m R , the measurement m R corresponds to the out-ofjet measurement B in the soft subjet factorization theorem, as it sets the lower scale. The multi-differential measurement is then made in the left hemisphere, which in the dressed gluon approximation, will give rise to a dressed gluon which will be integrated over the phase space in the left hemisphere. We will use somewhat interchangeably the notation left (right) and in (out) depending on the context, with the hope that this does not cause confusion.

A Single Dressed Gluon
We begin with the simplest case of one-dressed gluon. It is important to note that one cannot just simply replace the standard fixed-order calculation for the dressed gluon, since the presence of the resummation of the dressed gluon would change the renormalization of the global divergences. These we do not want to change, since in the hemisphere mass factorization theorem these divergences are tied to the resummation of the large global logarithms. The domain of consistency of our soft subjet factorization automatically imposes this constraint: we must only dress the gluon when it is energetic enough. Therefore, we only dress the gluon in the left hemisphere if m L > m R ; otherwise, we dress the gluon in the right hemisphere.
At one-loop, the cumulative soft function is therefore dressed as: The Θ-functions enforce the phase space constraints and, in particular, only turn the dressing factor on if the mass in a hemisphere is not set by the identified gluon emission in that hemisphere. The purely global contributions to the mass can be separated out by simple rearrangement: The global contribution, corresponding to the "1" in Eq. (4.12), has been studied in great detail in the literature, so for our purposes here, we will ignore it. Further, because there are no divergences associated with the dressed gluon, we can work in strict d = 4. Then, the one-dressed gluon contribution is The "NG" notation in the superscript denotes that we are only considering the non-global contribution to the soft function as captured by the one-dressed gluon.
The integrals can be evaluated and one finds where γ E is the Euler-Mascheroni constant and Γ[x] is the Euler Gamma function. Expanding to the first few orders, we find where ζ(3) = 1.202... is the Riemann ζ-function and we have used The α 2 s term is correct, while higher order terms are in general not, but this is expected for reasons discussed earlier as the dressed gluon approximation is not an expansion with a fixed logarithmic counting. In particular, the term at α 3 s differs by a factor of 2 from the true result in the large-N c limit, where C F → N c /2 [51]. Nevertheless, the two-dressed gluon approximation will fully capture this term, and produce a more accurate approximation for terms at even higher orders. We will shortly discuss in more detail the organization of the perturbative expansion in terms of dressed gluons.
In this section we have focused on extracting only the NGLs for the hemisphere mass distribution; however, it should be clear from the presentation that the dressed gluon approximation can be used to perform a complete calculation also including global logarithms. In particular, the expansion in dressed gluons reorganizes the perturbative expansion in terms of ordinary gluons in the low scale matrix element. The global resummation factor U (µ f , µ i ) multiplies the perturbative expansion in dressed gluons. Explicitly, if we have calculated to the th loop order, then formally we have We have subtracted all the fixed order NGLs c i L i that are included in the i dressed gluons, denoted by DG i . This is apparent from the factor "1" appearing in Eq. (4.12), which was ignored in this section, and will be further seen in subtractions necessary to extract the non-global contribution from the two-dressed gluon approximation in Sec. 4.2.4.

Generalization to N Eikonal Lines
In this section we generalize the construction of the dressed gluon approximation from the emission of a soft subjet from two eikonal lines in the n,n directions, to the case of the emission of a soft subjet from N eikonal lines. This construction is necessary to go beyond the one-dressed gluon approximation. Since the leading NGLs arise in the stronglyordered limit, we are interested in studying multiple strongly-ordered soft subjets, which will become the multiple dressed gluons. In the strongly-ordered limit, the factorization theorem can be obtained by performing a sequence of matchings, where at each stage, all more energetic subjets can be treated as eikonal lines, and all less energetic subjets are unresolved. This is shown schematically in Fig. 6 for the case of two strongly-ordered subjets, which will be discussed in detail in Sec. 4.2.4 where we consider the calculation of the two-dressed gluon contribution to the NGLs for the hemisphere jet mass. Because we can perform this sequence of matchings, to generalize the dressed gluon approximation it suffices to understand how to add a soft subjet to N eikonal lines. When we considered the factorization for a single soft subjet in Sec. 3, we emphasized that the soft subjet region of phase space can be isolated in an IRC safe manner by a multi-differential measurement of the energy correlation functions. This construction can be generalized to isolating m soft subjets by measuring the energy correlation functions e m+2 . The explicit condition on the phase space in terms of the energy correlation functions is not of particular interest, but has been discussed for isolating three prong structure in [98]. We therefore simply assume that sufficiently many IRC safe measurements have been made to isolate the desired region of phase space.
We now consider the addition of a soft subjet to N eikonal lines. We first discuss the organization of color before describing the anomalous dimension of the generalized dressed gluon. All color matrices encoding the color entanglement of the eikonal lines resides in the hard matching coefficient describing the soft subjet production from the initial N eikonal lines. To understand this organization of color, we start with an N (sub)jet factorization theorem in the region of phase space that has been isolated by a sufficiently differential measurement. This defines the initial hard matching coefficient H N and soft function S N .
To understand how to add the soft subjet to this factorization theorem, we consider the diagrams that can contribute to the matching coefficient describing the soft subjet production off of the N eikonal lines. Only diagrams which are fully color connected 16 that contain the real soft gluon that will form the soft subjet will contribute [106][107][108]. The disconnected diagrams will not be color entangled with the soft subjet evolution, and are reproduced by a soft function containing only the original N eikonal lines with the appropriate measurement constraints. Thus we are lead to conjecture the following factorization theorem that describes the addition of a soft subjet: can always be given a topological condition on the diagram of two-eikonal particle irreducible. These are the diagrams that appear naturally in the logarithm of the soft function. Here e res is the observable that sets the unresolved infrared scale and the B N are all the possible out-of-jet measurements. Explicit dependence on the jet radius and the arguments of the jet functions have been suppressed. The S i...k denotes a soft function with eikonal lines i...k. Note that for compactness, we have written the soft subjet jet function following the notation of Eq. (3.14), and have not explicitly written it as factorized into jet and boundary soft modes, although this factorization must be performed for a completely factorized description of the soft subjet dynamics. The notation tr sj denotes that we are to trace over all the color indicies of the eikonal lines i, ..., k. The adjoint indices AB are tied to the soft subjet's Wilson line, and contracts with the T color matrices of the matching coefficient H sj AB . The form of the factorization at the level of cut diagrams is given in Fig. 5.
The key feature of Eq. This matrix element selects all connected diagrams with a single soft jet state crossing the cut, while all disconnected diagrams cancel between the time-ordered and anti-time ordered products. Thus it is directly related to the logarithm of the diagrammatic expansion of a soft function. In App. C, we will use this conjecture for the factorization theorem to explicitly write the factorization theorem for two soft subjet production. There is substantial evidence for this factorization theorem. First, it preserves the renormalization group structure of the parent factorization theorem: there is a precise cancellation in the anomalous dimensions of the N -jet hard function and the N -jet soft function. Were we to introduce a soft eikonal line into the global soft function, this would violate the cancellation. In effect, we would not be able to factorize the soft jet production matching coefficient from the underlying hard process. The underlying N -hard jets have been fixed by the measurements imposed, and cannot be modified by subsequent infrared evolution. Second, the factorization formula is consistent with the known expression for the matrix elements for two soft partons, as given in Ref. [109], which will be discussed in detail below. Third, the generation of the hard matching coefficient for soft jet production can be considered as a differential operator acting on the soft function of parent eikonal lines that generate the soft subjet. Since the soft function is known to exponentiate only a strict subset of the diagrams involved in its calculation, this naturally leads to the hard matching given in Eq. (4.20).
From this factorization theorem, we can then define generalized dressed gluon factors W and G as As with the one-dressed gluon, we Laplace transform e res so that the limitẽ res → 0 of its Laplace conjugate results in being fully inclusive over e res . These W functions depend on the global color structure of the N eikonal lines in the system, but importantly is independent of the out-of-jet scales B N . Similar to the factorization of the global soft function in Eq. (4.3), in-jet and out-of-jet scales are split into different soft functions in the refactorization of the color multipole function. The G function depends on the directions of the lines participating in the soft subjet product. As was the case with the factorization theorem for one soft subjet, for hemisphere jets, the product of each W and G are renormalization group invariant. The structure of W can be expressed in a color multipole expansion as where W ij (z sj , Ω sj ) describes the emission of the dressed gluon from a color dipole, W ijk (z sj , Ω sj ) describes the connected emission of the dressed gluon from three Wilson lines, etc.
A possible modification to the conjecture is that soft subjet soft functions do not themselves form a color singlet as postulated above in the formulation of the color multipole function, so that is a generic color tensor over the color indices of the eikonal lines i...k. This color tensor must then be inserted between the contractions of the soft and hard functions H N and S N , following the T matrix conventions of Ref. [109] (since effectively it can be expressed as a sum over such color matrices). For soft functions with 2 or 3 Wilson lines, as considered in this paper, the function can always be written as a color singlet, consistent with the conjectured form discussed above.

Anomalous Dimension for Dressed Gluons with N Eikonal Lines
To dress the resolved gluon requires renormalization of the generalized W and G functions, as defined in Eq. (4.21). At one-loop, the anomalous dimension of these objects can be determined by the one and two soft gluon emission matrix elements from an arbitrary Npoint squared amplitude. This can be determined directly from the factorization equation (4.18), or from Ref. [109]. 17 In the notation of Ref. [109], for soft gluons with momenta q 1 , q 2 these squared amplitudes are

24)
A(q 1 , q 2 , p a 1 1 , ..., p a N N ) respectively, where the indices in parentheses denote the dipole from which the soft gluon has been emitted and p i and a i are the momenta and color of particle i. While the factorization above is valid for arbitrary soft gluon emissions, in the strongly-ordered limit, the soft gluon emission factors are where (s-o) denotes strongly-ordered.
To continue, we assume that at least one of the eikonal lines lives in the jet whose substructure we wish to probe. Then, Eq. (4.24) has the interpretation as the tree-level hard matching coefficient for the soft subjet while Eq. (4.25) describes the real emission contribution from the soft subjet. This emission can be either in the jet, where it contributes at the resolution scale of the jet e res , or it is outside the jet, where it contributes to the out-of-jet measurements B N . The soft real emission in the jet at the resolution scale could be either global soft, boundary soft, or radiation from other soft modes formed by having multiple eikonal lines in the jet. To determine the anomalous dimension of the W and G functions, and therefore to dress the gluon, we are only interested in emissions from the soft subjet that leave the jet.
The anomalous dimensions of W and G are determined solely by the term in Eq. (4.25) that is explicitly proportional to C A . Contributions from this term can entangle in-and outof-jet scales, corresponding to sensitivity to non-global structure. By contrast, the emissions from first term of Eq. (4.25) are color disconnected and therefore can only contribute to the global structure of either the in-or out-of-jet regions. These observations allow us to rewrite the C A term for emission of a strongly-ordered soft gluon with momentum q from the soft subjet in the direction n sj as We immediately recognize the factor in square brackets as the hard function for the soft gluon created from the eikonal line i, j. Therefore, to one-loop, the anomalous dimension of the generalized dressed gluon is where q = (1,q) andq is a unit vector in the direction of the emission from the soft subjet and "out" means that the angular integral is only evaluated in the out-of-jet region.
The coefficient of the dressed gluon anomalous dimension is given the by the gluon color Casimir, C A , multiplied by the O(α s ) cusp anomalous dimension [110][111][112], Γ 0 cusp (α s ) = α s /π. We conjecture that to all orders in perturbation theory the coefficient of the dressed gluon anomalous dimension is given by C A Γ cusp . While this may seem unmotivated from the above calculation, it appears naturally from the factorization theorem for the soft subjet in Eq. (3.10). Recall that the dynamics of the soft subjet which are sensitive to the non-global structure is described by the boundary soft mode. Because it is soft, the boundary soft mode sees the soft subjet as eikonalized and because it is collinear-soft, only sees the eikonal likes n sj andn sj . Thus it is natural that the anomalous dimension associated with the dressed gluon should be given by the cusp anomalous dimension. 18 For hemisphere jets, the renormalization group evolution of the W and G functions is straightforward to determine. Recall that we can expressed the W function in a basis of color multipoles as where the . . . denotes higher multipoles. The tree-level matrix element for the dipole factor is W where z sj is the energy fraction of the soft subjet and n sj is a lightlike vector along the direction of the soft subjet. The renormalization group evolution of W (ij) is which, to one-loop, the anomalous dimension is given by Eq. (4.28). At one-and two-loop order, only dipoles contribute to the evolution of the soft subjet [109,[117][118][119], and so only the evolution of W (ij) (z sj , Ω sj ) is necessary to resum logarithms associated with dressing a gluon through next-to-leading logarithmic accuracy. For reasons which will become clear when we discuss the relation of the dressed gluon to the BMS equation, we will use the following notation for the exponentiated anomalous dimension with fixed coupling where the normalization is chosen to correspond to the definition of the logarithm appearing in Eq. (4.16).

Two Dressed Gluons
We now continue and present the two-dressed gluon approximation as applied to the hemisphere jet mass. This discussion will be limited to leading logarithmic accuracy for the two-dressed gluon, but the extension to higher logarithmic orders can be accomplished straightforwardly by calculating the objects in the appropriate factorization theorem to higher perturbative orders. Unlike the one-dressed gluon, for which we presented explicit calculations, in this section we will only present the form of the integrand for the nonglobal component of the soft function for hemisphere jet masses which comes from the two-dressed gluon approximation. Numerical comparison of the two-dressed gluon to the BMS equation will be presented in Sec. 4.3.
As with the identification of the one soft subjet region of phase space, we perform a series of measurements on the jet to identify two soft subjets and a hard core of radiation. As discussed in Ref. [98], this three-prong region of phase space can be isolated by measuring e on the jet and demanding parametric relations between them. As the precise relationships are not vital to the two-dressed gluon calculation, we do not present them. To ensure that the soft subjets are well-separated, we can measure an additional three-point energy correlation function, e (β) 3 , and demand that e 3 . These measurements fully isolate the two soft subjet region, and further measurements can be performed to determine the subjet energy fractions and relative angles. The two soft subjet region of the cross section can then be factorized into appropriate hard, jet and soft functions describing the various modes that contribute to the various observables. As in the case of a single soft subjet, we must assume that we are in a region of phase space where B e (α) 3 , so that the NGLs are parametrically large. The precise structure of this factorization theorem, while potentially interesting for particular applications, is not relevant for the dressed gluon approximation and will not be presented here. As with the one-dressed gluon, to obtain the description of the two-dressed gluon from the factorization theorem, we associate functions in the factorization theorem and integrate over unresolved scales. Fig. 6a illustrates the two soft subjet region of phase space, where we assume that the energy of the subjets is strongly-ordered.
Once we have isolated two strongly-ordered soft subjets, we can then connect to the two-dressed gluon approximation by a sequence of matchings between different effective theories. This is illustrated in Fig. 6b, where we assume the subjets have momenta p and q with p q, and the procedure is similar in spirit to that proposed for describing the parton shower in SCET [120][121][122]. At the lowest scale, the softest subjet with momenta q is produced from the np ornp dipoles (this assumes that the soft subjets are color connected, or else they are produced from the nn dipole, a point we will discuss in more detail shortly). This emission is dressed by the anomalous dimension of the dressed gluon, Eq. (4.28), integrated up to the scale of the harder subjet, L (where the strongly-ordered limit breaks down). Once integrated to the scale of the harder subjet with momentum p, the procedure is repeated, with this subjet dressed by the anomalous dimension of the dressed gluon, and integrated up to the scale L, corresponding to the highest soft subjet scale.
To illustrate this procedure, we will first apply it to the one-dressed gluon, and show that it reproduces the calculation from Sec. 4.2.1. Using the exponentiated anomalous dimension from Eq.  We now immediately recognize the first term in Eq. (4.34) as the global hemisphere soft function, and therefore, the second term is the non-global soft function for the one-dressed gluon: One can verify that Eq. (4.35) agrees with the analytic expression for the non-global soft function for the one-dressed gluon in Eq. (4.14). Now, we apply the same procedure to determine the two-dressed gluon contribution, in the strongly-ordered limit. Unlike with a single dressed gluon, with two dressed gluons there are two possible color structures, depending on whether the dressed gluons are color connected or not. Following the discussion of Sec. 4.2.3, we know that the dressing in the two different cases is distinct, as it depends on the number of Wilson lines which are color connected to the soft subjet. In the case that the two dressed gluons are not color connected, they are only connected to the original nn dipole, and are therefore each dressed in the same way as a single dressed gluon. The contribution from two color disconnected dressed gluons is then given by the square of the one dressed gluon result, appropriately symmetrized. If instead the two dressed gluons are color connected, the lowest scale soft subjet, q, sees both the np, and pn dipoles, which must be included in the dressing. It is this case of two color connected strongly ordered dressed gluons, which we will now discuss in detail, however, we must remember that to obtain the complete result, we must add the color disconnected component.
The two-dressed gluon soft function is where we assume that p q, and the superscript "cc" indicates that this is the color connected component. Again, to isolate the non-global contribution, we re-associate the Θ-functions which yields The color connected, strongly-ordered global hemisphere soft function for two gluons is immediately identified as the first line of Eq. (4.37). The second line is the stronglyordered non-global soft function for the one-dressed gluon, with either gluon p or gluon q dressed. Because of the subtractions in this term, when dressed gluon p is collinear with n, its contribution vanishes, and similar for the contribution when dressed gluon q is collinear to n. The last line in Eq. (4.37) is the color connected, strongly-ordered non-global soft function for the two-dressed gluon: The subtractions that appear in the expression above from re-associating the phase space constraints from dressing are necessary from the effective theory perspective to force the non-global soft function to be restricted to the two-dressed gluon phase space region. These remove appropriate collinear limits of the two resolved gluons in the soft function, acting as zero bin subtractions in the effective theory [97].
To obtain the description of the full two-dressed gluon, including both color connected and disconnected contributions, we add Eq. (4.38) to the square of the one-dressed gluon, appropriately symmetrized. Then, the full two-dressed gluon is (4.39) The form of the two-dressed gluon calculation suggests that the all-orders non-global soft function calculated in terms of dressed gluons takes a form motivated by non-abelian exponentiation of the soft function [106,123,124]. We conjecture that this exponentiation of the color connected dressed gluons holds to all orders, and that the full non-global soft function can be written schematically as whereC i , denotes an exponentiated color factor [123]. A detailed study of this conjecture is, however, beyond the scope of this paper.

Numerical Comparison to the BMS Equation
In this section we compare our dressed gluon approximations with different calculations of the leading, large-N c NGLs for the hemisphere jet mass. We will compare to both fixed-order and resummed distributions for the NGLs. In the large-N c limit, the leading fixed-order NGLs were calculated in Ref. [51] to 5 loops by explicit iteration of the BMS equation. The expansion is where we recall that Leading logarithmic resummation of NGLs in the large-N c limit can be achieved by a numerical or Monte Carlo solution to the BMS equation. Dasgupta and Salam (DS) [33] proposed the following fit to their Monte Carlo result: We have written our own implementation of the DS Monte Carlo algorithm and will include it in the comparisons. We begin in Fig. 7 where we compare the one-dressed gluon approximation to the fixed-order and resummed results listed above. In the comparison, we include the 2-and 3loop fixed-order expansions, the DS fit, and our Monte Carlo output. As expected, because  it reproduces the full 2-loop result, the one-dressed gluon agrees with the 2-loop expansion at small values of L. As discussed earlier, the one-dressed gluon does not include the full 3-loop result and so the 3-loop expansion is a better approximation of the full resummed result out to about L = 1. Nevertheless, the one-dressed gluon is accurate to better than 5% of the resummed distribution out to L = 1.5, corresponding to a ratio of several hundred between the hemisphere jet masses. A distinction between the resummed distribution and the one-dressed gluon is that the one-dressed gluon diverges at sufficiently large L. This behavior will be present at any fixed order in the dressed gluon expansion, and therefore does not produce the physical large L distribution. However, such large L values where the dressed gluon diverges are well beyond the range of phenomenological applications. Fig. 7 also shows the output of our Monte Carlo, which agrees to within 1% with the DS fit out to L = 2.
Improved accuracy can be obtained by including the two-dressed gluon approximation which we show in Fig. 8. As with the one-dressed gluon, we compare with the DS fit and our implementation of the Monte Carlo, but now we compare with the 3-and 4-loop fixed-order expansions. The inclusion of the two-dressed gluon approximation fully accounts for the 3-loop fixed order result, and so nicely agrees out to about L = 1. Perhaps more interesting is that the two-dressed gluon approximation agrees with the Monte Carlo to better than 2% out to L = 2, and near there is actually a better approximation of the all-orders result than the 4-loop expansion. The accuracy of the two-dressed gluon is impressive and may suggest that the dressed gluon approximation is a convergent expansion. We also compare to the output of our Monte Carlo, which begins to diverge from the DS fit near L = 2.5, where finite cutoff effects or finite statistics could be important.
Starting at the two-dressed gluon approximation, there are contributions from both color connected and color disconnected soft subjets. Fig. 8 Figure 8: Comparison of the two-dressed gluon approximation to the 3-and 4-loop fixed order results for the non-global piece of the hemisphere mass soft function, the Dasgupta-Salam fit for the leading logarithmic resummation of the NGLs, and our implementation of the DS Monte Carlo. The two-dressed gluon contains both color connected and disconnected contributions, and we include in the comparison the color connected as well as the total two-dressed gluon approximation.
two-dressed gluon as well as the full two-dressed gluon approximation. The disconnected term of the two-dressed gluon approximation is the square of the one-dressed gluon approximation, appropriately symmetrized, and so first contributes at O(α 4 s ). The connected piece of the two-dressed gluon approximation contributes at O(α 3 s ), and so one may think that the disconnected contribution can be formally ignored until higher orders. However, in the dressed gluon expansion, both appear at the same order. Indeed, as shown in Fig. 8, the accuracy is enormously improved by including both the connected and disconnected contributions to the two-dressed gluon approximation. This re-emphasizes that the dressed gluon approximation is not an expansion in coupling or logarithmic counting, but rather in the number of resolved subjets as defined by the measurements performed on the jet. This is the picture that arises from the effective field theory framework, and the excellent convergence of the dressed gluon approximation provides strong support that this is the correct way to organize the expansion.

Insights into features of NGLs and the BMS Equation
In addition to the highly accurate description of NGLs, the dressed gluon approximation also provides a physical picture for the effects of a phase space boundary. In this section we discuss some features of NGLs and methods of expansion of the BMS equation, for which the dressed gluon approximation provides insight. From these examples, we are able to more precisely define the formal expansion of the dressed gluon approximation as an expansion in the unresolved phase space volume, which we discuss in Sec   Figure 9: (a) A schematic depiction of the buffer region of Ref. [39] which arises analytically in the dressed gluon approximation. As the out-of-jet scale is lowered (or correspondingly the in-jet scale is raised), the buffer region grows. (b) The width of the buffer region η buffer in pseudorapidity as a function of L, which exhibits linear growth. The width of the buffer region is only plotted for values of L for which the one-dressed gluon approximation provides a valid approximation to the NGL dynamics.

Buffer Region
In Ref. [39] a scenario for the possible underlying dynamics for NGLs was proposed, which identified a "buffer region" near the phase space boundaries where emissions are suppressed. Consider some phase space region Ω in which an energy veto is applied. Ref. [39] proposed that the mechanism for suppressing radiation emitted into Ω was due to a buffer region around the boundary of Ω which itself contained little radiation. A particularly interesting consequence of this proposal is the approximate geometry independence of the NGLs, as the buffer region smoothes the detailed shape of the boundary of Ω.
In Ref. [39] an evolution equation for the width of the buffer region as a function of the in-Ω and out-of-Ω scales was proposed based on some simple assumptions. The solution they found was where η buffer is the width of the buffer region in pseudorapidity, L and L are logarithms of two scales in Ω, and δη δL acts as an average speed of the evolution of the border of the buffer region in η, assumed to be independent of L and postulated to be proportional to C A . A Monte Carlo study was performed, which provided some qualitative support for the buffer mechanism; nevertheless, the exact linear relation of Eq. (4.43) was not observed. However, it was not clear if this was due to L dependence of δη δL for some unknown dynamics or simply because the system had not reached its asymptotic behavior.
The existence and properties of such a buffer region can be addressed by studying the soft jet region of phase space, as it describes the dynamics of a single gluon which itself radiates more gluons. The boundary soft function, its anomalous dimension, and its manifestation in the dressed gluon approximation provide strong support for the buffer mechanism. The buffer region is manifest in the one-dressed gluon approximation through the tree-level matrix-element which vanishes as the dressed gluon approaches the boundary of the jet. A schematic depiction of the buffer region for our geometrical setup is shown in Fig. 9a. As a proxy for the width of the buffer region as defined by the one-dressed gluon approximation, we consider the half maximum of the suppression factor in Eq. (4.44). That is, we define the width of the buffer region η buffer for hemisphere jet masses via The width of the buffer region η buffer is plotted in Fig. 9b for the range of L over which the one-dressed gluon approximation is valid. It appears to rapidly asymptote to linear growth as a function of L, which continues throughout the range of validity of the onedressed gluon approximation. While these values of L are smaller than those considered in Ref. [39], and so it is possible that different dynamics are involved, we believe that the dressed gluon approximation provides support for the buffer region. At higher values of L, additional dressed gluons are required to accurately describe NGLs, which is beyond the scope of this paper.

Expansions of the BMS Equation
In this section we relate our dressed gluon approximation with the BMS equation [40], and comment on various expansions proposed in the literature. The BMS equation for the purely non-global piece of the hemisphere mass distribution, g ab (L), is given by with the boundary condition g ab (0) = 1. For simplicity, in this section we will follow closely the notation in the BMS literature, in particular Ref. [51]. We will use the left/right instead of in/out labeling, and we define , (4.47) and Note that Eq. (4.48) is the exponentiated dressed gluon anomalous dimension from Eq. (4.32). The evolution equation resums all leading NGLs in the large N c limit. As it is a non-linear equation, solving it exactly analytically is challenging, and so several expansions have been proposed to systematically approximate it. One expansion that we presented in Eq. (4.41) is the fixed-order expansion to 5-loops from Ref. [51]. While this expansion provides insight at small values of NGLs, it exhibits poor convergence and may only be an asymptotic series for L > 1. Another approximation used in the literature [40,51] is to rewrite the BMS equation in the form (4.49) The second term does not contribute to three loops, and the first term gives a linear evolution equation which is straightforward to solve. Specializing to the case where a, b = n,n and performing the integrals, this linear equation becomes While the linear equation is easily solved in closed form, it has no particular region of validity and a numerical comparison to the solution to the BMS equation demonstrates that it is a poor approximation above L ∼ 0.5. The dressed gluon approximation, by contrast to either of these other expansions, is not an expansion in the coupling, but rather in the number of resolved soft subjets. It is therefore interesting to ask how the dressed gluon approximation manifests as an expansion of the BMS equation. Motivated by the assumed non-abelian exponentiation of the dressed gluon soft function discussed in Sec. 4.2.4, we expand the non-global function g ab (L) as g ab (L) = exp g (1) ab (L) +g (2) ab (L) + · · · . (4.51) We will relateg (1) ab (L) to the one-dressed gluon,g (2) ab (L) to the color connected piece of the two-dressed gluon, and so forth. While we have written g ab (L) in an exponentiated form, when working to a particular order, we expand g ab (L) to a fixed number of dressed gluons. For example, the expansion to second order is g ab (L) = 1 +g  Inserting this expansion into the BMS equation, the differential equation for the onedressed gluon functiong (1) ab (L) is then which is equivalent to the one-dressed gluon result calculated in Eq. (4.14). This truncation is also very similar to the linearized equation, Eq. (4.50). However, unlike Eq. (4.50), the one-dressed gluon is not itself exponentiated. It is interesting that this simple approximation to the BMS equation, which is a valid expansion independent of its interpretation in terms of dressed gluons, to our knowledge has not been well-studied in the literature. The approximation of Eq. (4.53) is also significantly more accurate than the linearized BMS approximation, Eq. (4.50).
While we have shown that this expansion reproduces the one-dressed gluon, it is interesting to compare the interpretation of the subtraction term in the BMS equation and the dressed gluon approximation. Note that the −1 subtraction appears explicitly in both the dressed gluon integrand of Eq. (4.13), and in the above expansion of the BMS equation Eq. (4.53), where it arises from the subtraction term −g ab (L) in the full BMS equation, Eq. (4.46). In the derivation of the BMS equation presented in Ref. [40], this term was included as a virtual subtraction by unitarity. While this is indeed the unique subtraction which renders the BMS evolution equation infrared finite, it does not necessarily correspond to virtual corrections from our point of view since jets can genuinely be collinear. On the other hand, in the effective field theory approach the origin of this subtraction is clear: it arises from restricting the effective field theory description of the soft subjet region to its regime of validity, similar to a zero bin subtraction [97]. 19 Indeed, in the effective field theory approach, because we have completely factorized the dynamics describing the soft jet, the real and virtual infrared divergences cancel separately within each function in the factorization theorem of Eq. (3.10).
The two-dressed gluon is found by expanding to higher order. Naïvely, the BMS equation for the two-dressed gluon would include derivatives of both one-and two-dressed gluon functions as:g    However, to obtain the dressed gluon approximation, we drop the term involving derivatives ofg (1) ab (L) as it first contributes at higher order in α s thang (2) ab (L). Then, the differential equation forg (2) ab (L) is , which upon integrating, is exactly the equation for the color connected part of the twodressed gluon in Eq. (4.38) with a, b = n,n. One can work to higher orders in the dressed gluon expansion, only keeping those terms in the differential equation that contribute to the same order in α s , and build up the full solution to the BMS equation.

Perturbative Expansion vs. Dressed Gluon Expansion
We conclude this section by discussing the distinction between a fixed-order expansion of the BMS equation and the dressed gluon expansion, defined in Eq. (4.51). We argue that the dressed gluon approximation is the correct way to organize the perturbative expansion of NGLs, as supported by the convergence and high accuracy of the dressed gluon approximation demonstrated in Sec. 4.3.
The dressed gluon approximation is not an expansion in the coupling or traditional logarithmic counting, but in the number of resolved soft subjets in a jet as defined by the soft subjet factorization theorem. An arbitrary number of unresolved gluons are emitted from the soft subjets so that the dressed gluon approximation includes terms to all orders in the coupling, yet for a fixed number of dressed gluons, it does not fully capture the complete logarithmic series to any formal accuracy using a traditional logarithmic counting. There are always contributions at the same logarithmic accuracy which arise from a higher number of dressed gluons that are ignored. As discussed in the study of the buffer region, the anomalous dimension of the soft subjet/dressed gluon suppresses emissions near the boundary of the jet, as well as suppressing the region of phase space when multiple resolved soft subjets approach one another. The interpretation of this is that higher-order dressed gluons have a significantly reduced phase space volume in which they can live. This has important consequences for how the formal expansion of the dressed gluon should be addressed.
In a fixed logarithmic counting, one implicitly counts the phase space for soft gluons as O(1), and so it does not affect the parametric scaling of the logarithms or the manner in which they should be counted. 20 This is appropriate when resolved emissions do not affect the allowed phase space volume for further emissions in the system; however, this is not the case with NGLs, as is nicely illustrated with the one-dressed gluon and its buffer region shown in Fig. 9. To take the phase space suppression into account in the logarithmic counting, we write the expansion in the schematic form where the ∆η (i) correspond to the phase space suppression from the buffer region with i dressed gluons, with ∆η (i+1) > ∆η (i) . The ∆η (i) are treated as negligible in a traditional logarithmic counting, but will be important for NGLs. For small L, the buffer region is small, and therefore the ∆η (i) can be formally neglected. In this case, higher powers of L are suppressed and we expect that the resummation of subleading terms, like α s L n , will be as important as the leading logarithms. For 20 Historically, the vast majority of global observables that were studied, like thrust [26], C-parameter [88,89,126], heavy jet mass [127], broadening [128][129][130], etc., define scales for which soft emissions are at lower, or at least the same, virtuality as collinear emissions. Thus the traditional logarithmic counting αs ln ∼ 1 accurately captures the dominant singular structure. However, there are phenomenologicallyrelevant examples of observables for which collinear emissions have lower virtuality than soft emissions; for example, recoil-free angularities with angular exponent β < 1 [73]. For these observables the naïve logarithmic counting must also be modified, instead taking (αs/β) ln ∼ 1, otherwise perturbative predictions will be extended well beyond their range of applicability. L ∼ 1, the buffer region is also O(1), and so the available phase space for resolved gluons is suppressed. With traditional logarithmic counting, L ∼ 1, and the entire first column of Eq. (4.56) must be resummed, as done by the BMS equation. However, including the power counting of the volume of phase space, higher order logarithmic terms are increasingly suppressed by the allowed phase space, a fact which should be taken into account in the organization of the perturbative expansion. Consistent incorporation of logarithmic and phase space counting is accomplished by the dressed gluon approximation through the resummation of the unresolved emissions associated with the dressed gluon. This explains why even the one-dressed gluon approximation exhibits such good convergence even to L ∼ 1. As L increases beyond 1, there is a competition between the higher powers of L and the phase space suppression, and so convergence of the expansion at large L requires including greater numbers of dressed gluons. However, since the total phase space volume is finite and the buffer region becomes large at large L, we expect that the expansion in the number of dressed gluons will rapidly converge.
These dynamics are not incorporated in a fixed-order expansion in which there is no distinction between the counting of resolved and unresolved gluons. Indeed, the fixed-order expansion does not seem to uniformly converge to the leading logarithmic resummation, suggesting it is an asymptotic series. A possible cause of the asymptotic nature of the series could be the perturbative expansion of the factor In the dressed gluon approximation, this arises from renormalization group evolution of the generalized dressed gluon anomalous dimension of Eq. (4.28). Traditional logarithmic counting assumes that the integrand appearing in the exponential is O(1) in all regions of phase space, which is not true. Maintaining this factor, as the dressed gluon does, appears to be vital for convergence of the expansion to arbitrary L values.

Resummation of NGLs to Higher Accuracy
The resummation of NGLs has thus far been restricted to leading logarithmic accuracy, where it is described by the BMS equation. We have described how the leading NGLs can be captured, and systematically calculated, by a sequence of factorization theorems producing the dressed gluon approximation. In this section we discuss how our approach can be extended beyond leading logarithmic accuracy.

Subleading Soft Corrections
We begin by briefly discussing the NGLs that can be calculated and understood using the soft subjet factorization theorem and the corresponding dressed gluon approximation. Working within the soft approximation, there are three ways in which the accuracy of the dressed gluon approximation can be improved: by including a greater number of dressed gluons, by calculating the anomalous dimensions and matching for the dressing to higher orders, or including soft dijets, whose energies are not strongly ordered. Because we have established a factorization theorem for soft subjets, and a set of observables that factorize in all regions of phase space, all improvements are well defined and in principle are straightforward to implement. However, the various types of improvements will affect the distributions in different ways, correspond to different expansions, and be relevant in different regions of phase space. As discussed in Sec. 4.4.3 the expansion to higher orders in the dressed gluon approximation is an expansion in the allowed phase space volume, while the expansion to higher orders in the dressed gluon anomalous dimensions and matching corresponds to a more familiar logarithmic expansion of the factorization theorem.
At large values of L, the most important corrections will be those from including an increasing number of dressed gluons at leading-logarithmic accuracy. However, calculating higher dressed gluons is of limited phenomenological interest, as this simply resums a more complete set of leading NGLs, and the two-dressed gluon approximation is already accurate at the percent-level for a wide dynamic range. Further, we argued that the dynamics of the buffer region and its phase space suppression modifies the naïve logarithmic counting, and therefore we expect convergence and high accuracy of the dressed gluon approximation to large values of L even for a limited number of dressed gluons.
On the other hand, for L 1, we expect the most important corrections from the resummation of subleading NGLs are captured by calculating the dressed gluon anomalous dimensions to higher order, and capturing effects that are simply beyond the soft approximation. This is evident from the numerical comparison of the one-dressed gluon approximation to the leading logarithmic resummation of the NGLs, which agree at the percent level for L 1.

Going Beyond the Soft Approximation
Subleading corrections to the leading logarithmic resummation are most important at small L. As L becomes increasingly small, the out-of-jet scale B approaches the soft jet scale, e effectively implemented a measurement on the soft jet emissions. This required a refactorization of the dynamics of the soft jet into soft jet modes, which do not resolve the boundary of the jet, and boundary soft modes, which are lower energy but do resolve the boundary. In the case that B e (α) 2 this argument is no longer valid, the dynamics of the soft jet should no longer be refactorized in this manner, and instead one should have a soft jet function which itself is sensitive to the jet boundary. Therefore, in this region of phase space, collinear splittings from the soft subjet are sensitive to the boundary of the jet. The contributions to subleading NGLs from such splittings is apparent from a detailed study of the fixed order calculation at two loops [46,48]. This region of phase space is schematically illustrated in Fig. 10.
While the NGLs arising in this region of phase space are formally subleading by traditional logarithmic counting, they may nevertheless be as important or even more important than the multiple dressed gluon effects of the leading NGLs. Furthermore, in the soft subjet factorization theorem, logarithms of the ratios of the in-jet to out-of-jet scale are resummed. 2 , but the relation is not parametric. In this region of phase space a soft subjet is pressed up against the boundary of the jet, and subleading NGLs arise from collinear splittings which cross the jet boundary.
In the collinear region of phase space, however, this ratio is order 1, and instead one must resum logarithms of the angle between the subjet axis and the jet boundary. Contributions from this region of phase space could therefore be important phenomenologically, and are also important for achieving a greater understanding of the dynamics of the buffer region. The approach to the resummation of NGLs presented in this paper, namely isolating a region of phase space, providing a factorized description of the dynamics, and resuming the NGLs by renormalization group evolution, can again be applied to understand this region of phase space. A factorization theorem describing this resummation, and allowing for the resummation of subleading logarithms arising from collinear splittings, along with a numerical study of their relevance will be presented in a forthcoming paper [131].

Ingredients for Subleading NGLs
Given these considerations, we summarize the needed ingredients for resumming the subleading NGLs: • One real emission with a virtual correction, and two real emissions in the color multipole function of a single dressed gluon.
• The matching of a single dressed gluon to two-loops. This requires the two loop soft gluon current, calculated in Refs. [132,133].
• Soft dijet matching in the not strongly-ordered regime, with a single real emission in their color multipole function.
• Incorporation of the subleading collinear effects at the jet boundary matched to the wide angle soft emissions.
Each of these contributions form indepedent physical processes, and due to the dynamics of the buffer region modifying the logarithmic counting, may not all contribute equally to subleading NGLs. Since the resummation of subleading NGLs has not yet been accomplished, it would be extremely interesting to determine if such a resummation has phenomenological implications.

Conclusions
In this paper we have presented a novel approach to the resummation of NGLs. By performing a multi-differential measurement on a jet, we are able to identify a phase space region involving a soft subjet, in which the NGLs can be factorized and resummed down to an unresolved infrared scale. Resummation is accomplished by renormalization group evolution in the factorization theorem and anomalous dimensions can be calculated to any perturbative accuracy. An understanding of the dynamics of the NGLs from the factorization theorem in the soft subjet region of phase space led us to introduce the dressed gluon approximation. We demonstrated how the dressed gluon approximation can be used to calculate NGLs by explicitly calculating the NGLs for the hemisphere jet mass distribution in the one-and two-dressed gluon approximations. These computations were compared with numerical solutions to the BMS equation, and expansions of it. The dressed gluon approximation exhibited excellent convergence over phenomenological values of L, and provided considerably better convergence than other expansions, including the fixed-order perturbative result. Indeed, the dressed gluon approximation does not correspond to a fixed L counting, but rather a novel expansion in the number of resolved gluons, which can be thought of as a perturbative expansion in distinct factorization theorems. We showed how this expansion at leading logarithmic accuracy can be obtained from the BMS equation. The dressed gluon approximation also gives an analytic realization of the buffer region, a proposed underlying description of the dynamics involved in the physics of NGLs. Futhermore, we have also discussed how our approach of isolating regions of phase space with multi-differential measurements to resum NGLs can be extended systematically beyond leading logarithmic accuracy.
We have also realized the BMS equation through a sequence of effective theories that produce the dressed gluon approximation. Part of this was due to the organization of the infrared degrees of freedom into the modes of soft-collinear effective theory. Repeating the soft subjet factorization decorates the original factorization theorem with new terms whose renormalization group invariance is independent of the parent factorization theorem's renormalization group structure. 21 Thus each term corresponding to a different number of dressed gluons in the dressed gluon approximation is realized as a renormalization group invariant object, using the SCET modes. This helps to illuminate why the derivation of the BMS equation has resisted a renormalization group derivation within the SCET framework. It would be interesting to see if one could formulate a different effective theory that gave the BMS equation directly as a renormalization group equation. A potential starting point may be reorganizing the infrared degrees of freedom using the so-called group space variables of [134] used to derive and solve the B-JIMWLK hierarchy [59][60][61][135][136][137] (see also [138]). The relationship between the the B-JIMWLK hierarchy and the BMS equation has been widely studied [42,50,53,70,139,140] and has been used to calculate leading NGLs with full color dependence. Nevertheless, it is difficult to see how one would incorporate the collinear degrees of freedom in such an approach, as well as the observable dependence of the global renormalization group that is the hallmark of soft-collinear factorization.
Understanding the factorization and all orders resummation properties of non-global observables is essential for connecting with many phenomenologically relevant jet observables. In this paper we have presented a first step towards this goal by presenting an effective field theory understanding of NGLs, and their relation to the soft substructure of jets. Since this paper presents the first step in understanding the soft substructure of jets, we conclude by discussing several important applications where we believe that our factorization theorem and understanding of the soft subjet region can be fruitfully applied.

→ 1 Jet Bin Transition for Electroweak Boson Production
Jet binning plays an important role in many LHC analyses, for example H → W W . In this example, the experimental sensitivity is highest in the exclusive zero-and one-jet bins due to the large tt background. There has been considerable study of the resummation for the exclusive zero-jet bin [19,[141][142][143][144][145][146], as well as for the exclusive one- [147,148] , and even two-jet bin [149]. However, in these cases a factorization theorem only exists in the case that the jets are at the hard scale. An important open problem is how to describe the transition from the zero-jet to one-jet region, where the jet has small p T . Attempts at understanding this region by combining information from different jet bins has been discussed in Ref. [150].
To understand the zero-jet to one-jet transition requires understanding the factorization theorem in the regime that a soft (sub)jet is formed. Our soft subjet factorization theorem provides a description of the dynamics in this region of phase space, and therefore can be used to study the transition.

Soft Subjet Region for Jet Substructure Observables
Another important application of our factorization theorem for the soft subjet region of phase space is towards the analytic understanding of jet substructure observables which resolve a two-prong structure, as required for boosted W/Z/H tagging. A complete description of the relevant phase space requires factorization theorems for one-prong jets, jets with hard, collinear subjets (described by the SCET + effective theory [90]), and jets with a hard core and a soft subjet, as presented here. While the focus of this paper has been of those aspects of the soft subjet factorization theorem as relevant for understanding NGLs, our factorization theorem provides a complete description of the dynamics in the soft subjet region of phase space. It can therefore be incorporated into a complete study of the phase space for two-prong jets. An analytic calculation for the substructure observable D 2 [78] combining the factorization theorems in each relevant region of phase space will be presented in a forthcoming publication [79].

Improving Monte Carlo Generators
Monte Carlo generators play a vital role in the accurate and realistic description of QCD processes at colliders. The soft subjet factorization theorem may have consequences for developing Monte Carlos that are accurate to beyond leading logarithmic accuracy. With the one-prong and collinear subjet regions, the soft subjet region completes the description of the e  n+1 } which completely characterize a 1 → n splitting. By randomly choosing a point in the multi-dimensional phase space and implementing the appropriate resummation for that region of phase space, one could envision a fully differential parton shower, accurate to the logarithmic accuracy of the factorization theorems in all regions of phase space.
results for the case that the soft subjet is produced by a gluon, off of the initial qq pair in e + e − annihilation. Other partonic configurations are straightforward, and obey the same type of factorization, but their hard production coefficient is not enhanced by the soft singularity 1/z sj .
The operator definitions in this section are given in terms of the collinear gauge invariant quark and gluon SCET fields [93,94], which we denote B µ ⊥n sj , χ n , as well as (lightlike) Wilson lines, S q . The Wilson lines extend from the origin to infinity along the direction of their specifying vector, q. Explicitly where P denotes path ordering, and A is the appropriate gauge field, and the color representation has been suppressed. Since we only consider the case of e + e − , all Wilson lines are outgoing. The soft Wilson lines carry the color representation of their parent collinear sectors, that is, adjoint representation for gluons and fundamental representation for quarks. Since we have no more than three Wilson lines in a soft function, the soft functions can always be written as color-singlet traces. In the more general case, the soft function is a color matrix, which must be traced against the hard functions, H(Q 2 ) and H sj appearing in the factorization theorem ( see e.g. Refs. [87,151] for more details). We will also use the large label momentum operator P µ [93] in the function definitions, which extracts the large component of the momentum for a particle in a given sector. We denote by Q the center of mass energy of the e + e − collisions, so that Q/2 is the energy deposited in a hemisphere, and Q SJ Q is the large component of the soft jet momentum. The functions appearing in the soft subjet factorization theorem of Eq. (3.10) have the following SCET operator definitions: • Soft Subjet Jet Function: 3 SJ B ⊥n sj µ (0)|0 • Jet Function: • Boundary Soft Function: The definitions of these functions include measurement operators, which when acting on the final state, return the value of a given observable. The operator e (β) 3 measures the contribution to e (β) 3 from final states, and must be appropriately expanded following the power counting of the sector on which it acts. The operators Θ F J , and Θ O constrain the measured radiation to be in the jet or out of the jet, respectively, and will be defined shortly.
The action of the measurement function e (β) 3 on a arbitrary state for each of the factorized sectors contributing to the three-point energy correlation function measurement is given by where, for simplicity, we have extracted the normalization factors Q hj =n · p hj , Q sj =n sj · p sj . (A.13) Q hj and Q sj are the large light-cone momentum components for the hard jet and the soft subjet, respectively. These expressions follow from properly expanding the definition of the energy correlation function measurements in the power counting of each of the sectors. Note that on the jet sectors, the three-point correlation measurement becomes an effective two-point correlation measurement, since the two-point energy correlation function is set by the initial splitting of the subjet. The in-jet restriction, Θ F J , is given by (A.14) The jet restriction must also be expanded following the power counting of the given sector.
We will see that this is actually quite subtle for the soft subjet modes, since the angle between the soft subjet axis and the boundary of the jet has a non-trivial power counting.
In particular, the expansion of Θ F J (k) is different for the soft subjet jet and boundary soft modes, and will demonstrate the necessity of performing the complete factorization of the soft subjet dynamics into jet and boundary soft modes. The explicit expansions in each sector's power counting will be given in App. B, when we consider the one-loop calculation of the functions appearing in the factorization theorem. Finally, since we are considering the case where the out-of-jet scale B is much less than the in-jet scale, the operator must also be included in the definition of the soft subjet functions. This operators vetoes out-of-jet radiation above the scale B. The explicit expression for Θ O (B) expanded in the power counting of each of the factorized sectors will be given in the one-loop calculations of App. B.

B One-Loop Calculations of Soft Subjet Functions
In this appendix we present the one-loop calculation of all the functions appearing in the soft subjet factorization theorem of Sec. 3.2, whose operator definitions are given in App. A. As discussed when defining the soft subjet functions App. A, we will only give results for the case that the soft subjet is produced by a gluon, although it is straightforward to extend the calculation to other partonic configurations. Throughout this section, we will make use of the convenient shorthand notation for the integration measure of an on-shell, massless, final-state parton. For the jet and soft functions, we only give the final expressions in the Laplace space of e (α) 3 , where they satisfy a multiplicative renormalization group evolution. This allows for a straightforward comparison of the anomalous dimensions.

B.1 Hard Matching for Dijet Production
The hard matching coefficient, H(Q 2 ), is the well known hard function for the production of a qq pair in e + e − annihilation. It is defined by where C(Q 2 , µ) is the Wilson coefficient obtained from matching the full theory QCD currentψγ µ ψ onto the SCET dijet operatorχ n γ µ ⊥ χn. This Wilson coefficient is well known (see e.g. [4,87,90,152] ), and is given at one-loop by The branch cut in the logarithms must be taken as −Q 2 → −Q 2 − i .

B.2 Hard Matching for Soft Jet Production
The hard matching coefficient H sj (z sj , θ sj ) is determined by the finite parts of the logarithm of the soft matrix element for a single soft state The virtual corrections of the effective theory cancel the IR divergences of this matrix element, giving a finite matching coefficient. This matrix element can be calculated from the square of the soft gluon current [153,154], which is known to two loop order [132,133].
Here, for simplicity, we restrict ourselves to one-loop accuracy. The tree level and one-loop hard matching coefficients for the soft subjet production are given by H sj(tree) nn (z sj , n sj ) = α s C F 4π 2 z sj n ·n n · n sj n sj ·n , sj n · n sj n sj ·n The results of [154] can be used to determine the soft-jet production matching from an arbitrary number of hard jets at one loop.

B.3 Jet Function
In this section we calculate the jet function for the energetic subjet along the n direction. The one-loop expression for the naïve (before zero bin subtraction) jet function is Here we have chosen to calculate the jet function by integrating against the splitting function [155]. Since we have assumed the partonic configuration in which the soft subjet is a gluon jet, the jet in the n direction is assumed to be described by a collinear quark field. For the splitting functions we use the (slightly unconventional) notation where the denote that the splitting functions are spin averaged.
The jet algorithm and measurement constraint are given by where we have used the expression for the action of the e (β) 3 measurement on a hard jet state from Eq. (A.7).
In the power counting of the n collinear sector, the second two terms vanish upon performing the multipole expansion on the jet function constraint. The first term simplifies since for any particle i in the n collinear sector. The phase space for the two partons in the jet with these constraints imposed is then given by where Z defines the large momentum fractions of the partons as n · k 1 = Q J Z ,n · k 2 = Q J (1 − Z) , (B.14) and the angle φ is defined by Substituting this into the expression for the jet function gives (B. 16) Performing the integrals as an expansion in and transforming to Laplace space, we find 3 is the Laplace conjugate of e (β) 3 , and we have explicitly separated the finite and divergent pieces. The argument of the logarithms is given by Note that all zero bins for the jet function vanish.

B.4 Soft Subjet Jet Function
In this section we calculate the jet function for the soft subjet itself. Since the soft subjet is near the boundary of the jet, we will see that we must carefully treat the jet boundary constraint, emphasizing the role of the boundary soft mode. The one-loop expression for the naïve (before zero bin subtraction) jet function is where n f denotes the number of light flavors, and T F = 1/2 specifies our normalization convention for the SU(3) algebra. Here we have taken the soft subjet to be a gluon jet, and have again chosen to calculate the jet function by integrating over the splitting functions, where the arguments of the splitting function denote the energy fraction of the two partons, as defined in Eq. (B.8). The jet algorithm and measurement constraint are given by where we have used the expression for the action of the e (β) 3 measurement on a soft jet state from Eq. (A.7). Since we are considering the case where the out-of-jet scale B is lower than the in-jet scale, we can multipole expand the constraint in the out-of-jet region as which eliminates the second two terms in Eq. (B.11). This implies that the jet boundary effectively acts as a hard wall for radiation in the soft subjet jet function. For the jet modes of the soft subjet, we can also multipole expand the jet function constraints The fact that this constraint can be multipole expanded follows from the power counting in Sec. 3 2 . This can also be seen from expanding the jet constraints in the local soft subjet coordinates, where we find n · n sj n · n sj + 4 k ⊥sj · n ⊥sj (n · n sj ) 2n sj · k + ... , where we have used that n ⊥ sj = −n ⊥ sj . The ability to perform this multipole expansion relies crucially on the fact that we have fully factorized the dynamics of the soft subjet into jet modes and boundary soft modes. For the boundary soft modes, we cannot perform the above multipole expansion. This implies that the soft subjet jet function will not depend on the factor tan 2 R 2 − tan 2 θ sj 2 , as should be the case for a collinear function, while the boundary soft function carries the entire dependence of the soft subjet dynamics on the difference tan 2 R 2 − tan 2 θ sj 2 . Since both particles are constrained to lie within the jet, and the jet boundary constraint is multipole expanded, the phase space constraints for the soft subjet jet function are identical as for the standard jet function for the variable e (β) 3 , but for a gluon jet. The explicit expression will be given shortly. Alternatively, it is possible to calculate the soft subjet jet function without performing the multipole expansion on the jet constraint. In this case one find that the phase space for the jet function is corrected by a term which depends on tan 2 R 2 − tan 2 θ sj 2 , leading to a correction to the jet function depending on tan 2 R 2 − tan 2 θ sj 2 . However, we have explicitly checked that performing the appropriate boundary soft zero bin subtraction entirely removes this correction, again emphasizing the importance of this mode. We therefore stress the importance of a proper power counting analysis when analyzing the effective theories for more complicated jet configurations. A similar feature was also noted in Ref. [87] for the calculation of different individual jet functions with a jet algorithm constraint.
We now give explicit expressions for the gluon jet function. Since the phase space is identical to that given in Eq. (B.13), but with the integration performed against the splitting functions as indicated in Eq. (B.20), we simply give the final result. Performing the integrals as an expansion in and transforming to Laplace space, we find where β 0 is defined with the normalization and whereẽ (β) 3 is the Laplace conjugate of e (β) 3 . We have explicitly separated the finite and divergent pieces. The argument of the logarithms is given by (B. 28) Note that all zero bins for the soft subjet jet function vanish.

B.5 Global Soft Function
In this section we calculate the one-loop global soft function. The soft function involves three eikonal lines in the n,n, and n sj directions, since the angle between the soft subjet axis and the n andn axes is O(1), and is therefore resolved by the soft radiation. This is distinct from the situation in the SCET + factorization theorem of Ref. [90]. We will see the importance of the performing the appropriate zero bin subtractions, and the role of the boundary soft mode. Indeed, the fact that the soft function has a non-trivial zero bin is itself unusual. The general form of the one-loop soft function is (see e.g. Ref. [87]) where T i is the color generator of leg i, and the sum runs over all pairs of legs. The global soft radiation is at a scale such that it can contribute to both the in-jet and out-of-jet observables. Since we work only to one-loop in this appendix, the integral in the soft function is over the phase space for a single parton. We can therefore straightforwardly separate the in and out-of-jet contributions through the measurement functions in: where k denotes the momentum of the soft parton. In this section we will split the calculation into two pieces, considering first the in-jet contribution, and then the out-of-jet contribution. This is important to emphasize that contributions to the soft function which depend on large logarithms of tan 2 R 2 − tan 2 θ sj 2 arise only from the out-of-jet region of integration. Although such logs naively appear in the in-jet contribution to the soft function, they are removed by the boundary soft zero-bin subtraction.
To one-loop, the soft function for the exchange between the eikonal lines n a and n b is given by for the in-jet region, and for the out-of-jet region. Following the decomposition in Eq. (B.29), we have explicitly extracted the color factor, so that it does not appear in these expressions. The dressed gluon approximation holds for an arbitrary additive observable, B, for example, in Sec. 4.2 we used the mass as an example. Here, for simplicity we have chosen to measure the energy in the out-of-jet region. Since the soft subjet soft function contains in its definition the three eikonal lines n,n, n sj , we must sum over contributions from exchanges between all possible pairs.

Naïve In-Jet Soft Function
We begin by calculating the naïve (i.e. without zero-bin subtraction) contributions to the in-jet soft function. For simplicity, we give only the finite pieces, dropping -divergences.
The anomalous dimensions will be given in App. B.7. The contributions from the exchange between the three possible pairs of eikonal lines are given bỹ Here we have extracted the common factor 36) as well as the functions R (1) , B (1) , and constants C (1) . The functions R (1) depend only on θ sj and R, and are given by Where the integrals I are given as: The functions B (1) contain singular dependence on the difference between θ sj and R, that is, the angle of the soft jet to the jet boundary, and are given as:

Boundary Soft Zero-Bin of In-Jet Soft Function
We now calculate the boundary soft zero bin of the in-jet soft function. This is the only nonvanishing zero bin. Both constraints in the soft measurement function can be expanded in the zero bin. The jet boundary constraint can be expanded as θ tan 2 R 2 − n · k n · k → θ tan 2 R 2 − tan 2 θ sj 2 + 4 k ⊥sj · n ⊥sj (n · n sj ) 2n sj · k , (B.49) where we have used the expression given in Eq. (B.24) for the expansion of the jet constraint. Note importantly that for the boundary soft modes, this cannot be multipole expanded, unlike for the jet modes of the soft subjet, as was discussed in Sec. B.4. For the measurement function, we have the expansion δ e (β) 3 − 2 −1+ β 2 N Sn sj · k Q n sj · k n sj · k β 2 (n · n sj ) β 2 (n · n sj ) 1−β . (B.50) Furthermore, in the integrand we can make the following expansions in the zero bin n · n sj n · k k · n sj → 1 n sj · k k · n sj ,n · n sj n · k k · n sj → 1 n sj · k k · n sj . (B.51) Performing the integration, we find the the zero bin contribution arising from the exchange between the n andn Wilson lines vanishesS  Here the superscript "b.s.b." indicates that this is the contribution from the boundary soft zero bin.

Zero Bin Subtracted In-Jet Soft Function
We now give the expression for the in-jet soft subjet soft function after performing the zero bin subtraction. TheS nn terms are unaffected by the zero bin subtraction, however we include them so that we can rewrite all contributions in a similar form. After zero bin subtraction, the contributions from the three different exchanges are given bỹ S (1, fin) n nsj (ẽ  We see that the potentially large logarithm of tan 2 R 2 − tan 2 θsj 2 , which was present in the in-jet soft function before zero bin subtraction has been removed by the boundary soft zero bin subtraction from the terms that contribute to the anomalous dimension, again emphasizing its crucial role in the factorization theorem. We also emphasize that the presence of a non-trivial zero bin for the soft function is an interesting feature of this factorization theorem.

Out-of-Jet Contribution to Soft Function
In this section we calculate the contribution to the soft function from out-of-jet radiation. While we have seen that for the in-jet contribution the large logarithm of tan 2 R 2 −tan 2 θsj 2 was removed by the zero bin, this will not be the case for the out-of-jet radiation. When performing the calculation we will integrate over the entire out-of-jet region, except a region of radius R B around the axis of the jet in the right hemisphere. This acts as a regulator, allowing us to calculate each of the contributions to the soft function, S (1, out) ij (B) separately. For reference, we take the out of jet measurement to be the cumulative energy deposited, however, we have explicitly checked that using other out-of-jet measurements lead to the same dressed gluon anomolous dimension as given in Sec. 4.1.
We begin by calculating the naïve (non zero bin subtracted) soft function, whose integrand was given in Eq. Again, we see the explicit appearance of tan 2 R 2 − tan 2 θsj 2 in the out-of-jet contribution to the soft function. Although R B is required as a regulator in each of the S (1, out) ij (B), the sum of Eqs. B.60, B.61, and B.62 are non-singular as R B → 0, and all anomalous dimensions should be expanded in this limit. For a further discussion, see [87]. In other words, B ∼ Q z 2 s , with z s the energy fraction of the global soft radiation. All other modes appearing in the factorization, in particular, the boundary soft, and soft subjet modes, are parametrically more energetic. For all possible zero bins to the soft function in the out-of-jet region, we can therefore multipole expand the measurement function appearing in Eq. (B.32). Therefore, all such zero bins give a vanishing contribution. This implies that the dependence on tan 2 R 2 − tan 2 θsj 2 is not zero bin subtracted in the out-of-jet contribution to the soft function, unlike for the in-jet contribution. Therefore, while the global soft function does depend on tan 2 R 2 − tan 2 θsj 2 , it comes entirely from the out-of-jet region of integration. As we will see in App. B.6, the boundary soft function also depends on tan 2 R 2 − tan 2 θsj 2 . Indeed, tan 2 R 2 − tan 2 θsj 2 appears in the anomalous dimensions for both these functions, but with opposite signs, as is required by the renormalization group consistency of the factorization theorem.

B.6 Boundary Soft Function
In this section we calculate the boundary soft function. In the multi-stage matching of Sec. 3 which gave rise to the boundary soft mode, the boundary soft modes were decoupled from the soft subjet collinear modes via a BPS field redefinition. The boundary soft function therefore has the form of a global soft function, in particular it is calculated with eikonal Feynman rules, but it only has Wilson lines in the n sj andn sj directions. This is important, as it implies that the boundary soft function has the same color structure as the soft subjet jet function. In this appendix we have assumed that the soft subjet is a gluon jet. This can be understood intuitively since the boundary softs are a collinear soft mode, and hence are genuinely boosted in the n sj direction, so all other Wilson lines collapse to then sj . Thus the color structure of the boundary soft modes is simply that of the dipole formed by the soft subjet, and all other eikonal lines merged into one.
The one loop expression for the boundary soft function is given by 3 − N BSn sj · k Q n sj · k n sj · k β 2 Θ tan 2 R 2 − n · k 1 n · k 1 n sj ·n sj n sj · k k ·n sj . (B.69) Here we have already multipole expanded away any possible out-of-jet contributions, since the boundary soft scale is higher than the out-of-jet scale. We must again take care in expanding the jet radius constraint. From Eq. (B.24), we have Θ tan 2 R 2 − n · k 1 n · k 1 → Θ tan 2 R 2 − tan 2 θ sj 2 + 4 k ⊥sj · n ⊥sj (n · n sj ) 2n sj · k , (B.70) which, unlike for the soft subjet jet function, cannot be further expanded. The one loop expression for the boundary soft function is then given by 3 − N BSn sj · k Q n sj · k n sj · k β 2 Θ tan 2 R 2 − tan 2 θ sj 2 + 4 k ⊥sj · n ⊥sj (n · n sj ) 2n sj · k n sj ·n sj n sj · k k ·n sj . (B.71) We therefore see that the boundary soft contribution is identical to the soft subjet collinear-bin of the global softs, given in (B.53), with the substitutionn →n sj , and changing the normalization of the measurement function. We can therefore immediately write down the one-loop boundary soft functionS (1) nsjnsj (ẽ (β) 3 ) = α s C A 2π   π 2 6 + π 2 8(1 − β) − π 2 β 12 (B.72) where for simplicity we have given only the finite pieces, dropping -divergences. We see that the boundary soft mode carries the dependence of the soft subjet dynamics on the difference, tan 2 R 2 − tan 2 θsj 2 , which is completely factorized from the collinear dynamics of the soft subjet. However, importantly, the color structure of the boundary soft is determined by the color structure of the soft subjet, showing that it is indeed describing its dynamics. The difference tan 2 R 2 −tan 2 θsj 2 therefore appears in both the boundary soft function, and in the out-of-jet contribution to the soft function, as seen in App. B.5. The fact that it appears in both these functions is required for the renormalization group consistency of the factorization theorem.

B.7 Anomalous Dimensions
In this section we collect the one-loop anomalous dimensions for all the functions calculated in this appendix. The two hard functions satisfy multiplicative renormalization group equations. For the dijet production hard function, we have For the soft subjet production hard function, we have µ d dµ ln H sj nn (z sj , n sj , µ) = − α s C A π ln 2µ 2n · n Q 2 sj n · n sj n sj ·n + α s π β 0 . (B.75) the soft subjet production. Again, both the q, and p soft subjet jet functions must be refactorized into boundary soft and jet functions to achieve a complete factorization of the soft subjet dynamics. Finally, the ... terms in Eq. (C.3) denote terms which involve all three eikonal lines n,n, and p in the production of the second soft subjet. These terms do not appear in the tree-level matching.