Non-global logarithms at finite Nc beyond leading order

We analytically compute non-global logarithms at finite Nc fully up to 4 loops and partially at 5 loops, for the hemisphere mass distribution in e+e- to di-jets to leading logarithmic accuracy. Our method of calculation relies solely on integrating the eikonal squared-amplitudes for the emission of soft energy-ordered real-virtual gluons over the appropriate phase space. We show that the series of non-global logarithms in the said distribution exhibits a pattern of exponentiation thus confirming - by means of brute force - previous findings. In the large-Nc limit, our results coincide with those recently reported in literature. A comparison of our proposed exponential form with all-orders numerical solutions is performed and the phenomenological impact of the finite-Nc corrections is discussed.


Introduction
As the LHC started colliding protons with unprecedented beam energy, interest has risen in topics that had not received much attention previously, with the aim of uncovering new physics signals. Of these topics, that have recently seen substantial development, is the substructure of "fat" jets originating from the almost-collinear decay products of heavy resonances that are highly boosted (see for example refs. [1][2][3][4][5][6][7][8][9][10][11][12][13]). Many substructure techniques, such as filtering [1], pruning [4] and trimming [14], have been developed for the purpose of improving the discrimination of signals from QCD background. The latter substructure techniques aid in providing cleaner and more accurate measurements of the properties of these resonances through: first, identifying the origin of the jet (decayed massive particle -signal -or plain QCD radiation -background-). Second, mitigating away the jet constituent particles that have most likely originated from initial-state radiation, underlying-event and pile-up.
These techniques require in many situations calculations of QCD observables (e.g., jet mass, jet shapes, etc) which need special attention particularly in the vicinity of the threshold limit where they become highly contaminated with perturbative large logarithms as well as non-perturbative corrections. There has recently been some analytical work on a handful of the said substructure techniques with the aim to pave the way to a better understanding of their analytical properties (see [15] and references therein). Nonetheless and for the majority of QCD observables and substructure techniques the only other option available is resorting to numerical simulations which are based on Monte Carlo (MC) integration methods, and which use several approximations, e.g., Herwig [16,17], Pythia [18,19] and Sherpa [20]. These MC event generators have been very successful in describing collider data and are commonly used in the extraction of crucial information to boost the search for new physics.
An important issue that needs addressing is the accuracy of the said MC algorithms and the range of validity of the approximations used therein. For instance, amongst the widely adopted approximations in the said MC generators is that of large-N c limit (with N c the dimension of SU(N c ) group). The latter limit, which corresponds to neglecting nonplanar Feynman diagrams, greatly simplifies the otherwise tremendously complex colour structure, especially at high multiplicities. However, MC generators are generally tuned with data from collider experiments for parameters that account for non-perturbative effects such as hadronisation, underlying event, etc. The process of tuning is itself vulnerable to erroneously ascribing neglected perturbative (observable-dependent) components, which might be originating from finite-N c corrections, to universal non-perturbative parameters. This could then potentially be a source of major discrepancy between the data and the predictions by the MC generators. It is thus of great importance to assess the validity of these approximations and make sure that neglected terms would not affect precision measurements.
Amongst the issues that MC generators are meant to tackle is that of the resummation of large logarithms typically inherent in the distributions of most observables. These large logarithms are a manifestation of the miscancellation of infrared/collinear singularities at the matrix-element level, due to the exclusion of real-emission events in certain regions of phase space. For several observables of sufficiently inclusive nature, 1 i.e., global observables, the resummation of these logarithms is relatively straightforward and has even been achieved analytically to NNNLL accuracy [21]. In fact semi-numerical programs have been developed with the power of resumming a wide range of global observables up to NLL (CAESAR [22]) and even to NNLL recently (ARES [23]). However the extension of the resummation programme, up NNLL or even to just NLL accuracy in some cases, has seen slow progress for another class of observables, namely non-global observables [24,25].
Non-global observables are observables that are sensitive to emissions in restricted angular regions of the phase space. The distributions of such observables contain logarithms (named non-global logarithms (NGLs)) of the scales present in the process. For instance, the hemisphere mass distribution contains logarithms of the ratio Q/(Qρ) , where Q is the hard scale of the process and ρ is the normalised hemisphere mass squared. In the region where ρ ≪ 1, these NGLs can form large contributions to the said distributions, and should thus be resummed to all-orders. Up to very recently, their resummation was only possible numerically at large N c by means of: an MC program [24,25] or solutions to the non-linear integro-differential Banfi-Marchesini-Smye (BMS) evolution equation [26]. The large-N c approximation significantly simplifies the colour flow in multiple gluon branchings enabling the possibility of the resummation of NGLs at least numerically. Much effort has recently been advocated to achieving numerical (analytical) resummation of NGLs at finite (large) N c . The work of Hatta and Ueda [27] exploits the suggestion of Weigert [28] to use an analogy between the resummation of small-Bjorken-x (BFKL) logarithms and that of NGLs at finite N c in a numerical fashion. They have noticed that the neglected finite-N c corrections are indeed negligible in the context of e + e − → di-jets. They have however speculated that the situation may be drastically different for hadronic collisions. Furthermore, Rubin [29] numerically computed the NGLs series for both filtered Higgs-jet mass as well as interjet energy flow observables up to five-and six-loops, respectively, at large N c . In the same limit, Schwartz and Zhu [30] worked on the analytical solution to the BMS equation by means of an iterative series-solution up to five-loops.
The major hindrance that one inevitably faces when attempting to compute NGLs analytically at finite N c is twofold. Firstly, the colour topology of a multi-gluon event requires evaluations of non-trivial traces of colour matrices in SU(N c ), which become increasingly cumbersome starting from four-loops. Secondly and not less important, the non-Abelian gluon branchings increase the number of Feynman diagrams factorially at each escalating order to the extent that an automated way of accounting for all possible branchings becomes inescapable. 2 Besides, there is also the issue of the various possible real, virtual and real-virtual gluon configurations that are eventually responsible for the miscancellation of soft singularities, thus leading to the appearance of large logarithms. These difficulties may have been the main reason for the slow progress in the resummation of NGLs at finite N c .
In this paper we overcome the above-mentioned difficulties and present the first analytical calculation of NGLs at finite N c beyond leading order. Working in the eikonal approximation [31][32][33][34][35][36] for soft (strongly) energy-ordered partons, the first problem, i.e, that of colour structure, is resolved via the use of the Mathematica package ColorMath developed by Sjödahl [37,38]. The latter program performs the summation of SU(N c ) colour matrices in an automated way at any loop order. For the second obstacle we developed a Mathematica code that accounts for all possible gluon branchings (and thus for all possible antenna functions) in an automated way. 3 Consequently we have been able to analytically calculate all squared amplitudes for the emission of soft energy-ordered gluons (for all possible real, virtual and real-virtual configurations) in the eikonal approximation fully (at finite N c ) up to five-loops. We leave the presentation of these squared-amplitudes and the method of calculation to a forthcoming paper [39].
With these squared-amplitudes at hand, we provide in this paper a calculation of NGLs at finite N c to single logarithmic accuracy for single-hemisphere mass distribution in e + e − → di-jets up to five-loops. While our calculation is full at four-loops it is incomplete at five-loops due to missing terms for which the squared amplitudes are not so simple to simplify and/or integrate. We find that the aforementioned distribution exhibits a pattern of exponentiation both for global and non-global logarithms. We consequently write the all-orders distribution as a product of two exponentials; the first being the usual Sudakov form factor and the second represents the "resummed form factor" for NGLs. For the sake of cross-checking we take the large-N c limit of our result and compare it with previous calculations obtained by Schwartz and Zhu [30]. We find complete agreement up to our accuracy, which is four-loops. Furthermore, we compare our analytical resummed factor to the available all-orders numerical results [24] and discuss the phenomenological implications of our findings, particularly the issue of the accuracy of the large-N c limit, by assessing the importance of neglected finite-N c corrections up to four-loops.
This paper is organised as follows. In the next section we outline the usual procedure of calculating NGLs by defining the observable, kinematics and the general relation for the hemisphere mass distribution in terms of the squared amplitudes and a "measurement operator". 4 We present, in the same section, the calculation of NGLs at leading order (two-loops) to warm up for higher loops. In section 3 we explicitly calculate NGLs beyond leading order at three, four and five-loops. The difficulties associated with calculations at five-loops will be addressed therein. We compare our findings, in section 4, to those obtained at large N c in ref. [30] as well as to the all-orders parametrised form reported by Dasgupta and Salam, which they obtained by fitting to the output of their MC program [24]. We also assess the relative size of the corrections due to finite N c up to four-loops and discuss our findings in the same section. Finally, we conclude our work in section 5.

Hemisphere mass distribution at one and two-loops
Our aim in this paper is to calculate NGLs at finite N c to single logarithmic accuracy up to the fifth order in the strong coupling α s (or equivalently up to five-loops). Our calculation is performed using QCD squared-amplitudes for the emission of energy-ordered gluons in the eikonal approximation. The latter is sufficient to capture all single logarithms α n s L n , with L being the large NGL. As stated in the introduction, we do not show herein explicit formulae for the said squared-amplitudes and refer the reader to our coming paper [39]. Moreover, for the purpose of this paper, we do not consider the role of any jet algorithm, and postpone such work to future publications.

Observable definition and kinematics
For the sake of illustration and to avoid unnecessary complications from a hadronic environment, we choose to work with the same observable that was used in the original paper on NGLs by Dasgupta and Salam [24] within the framework of QCD, that is, the hemisphere mass distribution in e + e − → di-jets. This very observable was also considered in refs. [40][41][42] in the context of soft-collinear effective theory (SCET). In both refs. [24,[40][41][42] NGLs were only computed up to two-loops. In the eikonal approximation, sufficient for our purpose, we consider energy-ordered soft gluon emissions: 5 Q ≫ k t1 ≫ k t2 ≫ · · · ≫ k tn , with Q the centre of mass energy and k ti the transverse momenta of emitted gluons k i . We note that gluon decay into quarks has a sub-leading contribution to NGLs as was found at two-loops in refs. [40,43].
The four-momenta of the outgoing quark, anti-quark and gluons are expressed in rapidity parametrisation as: where recoil effects are negligible to single logarithmic accuracy. Here η i and φ i are the rapidity and azimuthal angle of the i th emission and k ti its transverse momentum with respect to the z-axis, which we choose to coincide with the outgoing quark direction. We have k ti = ω i sin θ i , with ω i the energy of gluon k i and θ i its polar angle. The rapidity is related to the polar angle θ i through the relation η i = − ln tan(θ i /2). We define the following "antenna" functions relevant to the squared amplitudes that we use in this paper: The quark and anti-quark directions define two coaxial back-to-back hemispheres (H L and H R ) whose axis coincides with the thrust axis at single logarithmic accuracy (see figure 1). We pick for measurement the hemisphere pointing in the positive z-axis (quark direction). The normalised hemisphere mass (squared) ρ is then defined by: where we introduced the transverse momenta fractions x i = k ti /Q , and the sum over the index i extends over all emitted real gluons in the measured hemisphere H R .  We compute the integrated hemisphere mass distribution (cross-section) normalised to the Born cross-section, defined by: where W X 12···m = W X (k 1 , k 2 , · · · , k m ) is the eikonal matrix-element squared for the emission of m energy-ordered soft gluons of configuration X off the primary qq pair at m th order, normalised to the Born cross-section. The sum over X extends over all possible real (R) and/or virtual (V) configurations of all the gluons {k j }. For instance, at 2 loops (m = 2) the eikonal squared-amplitudes W X 12 over which the sum is taken are: W RR 12 , W RV 12 , W VR 12 , and W VV 12 . The quantity W RV 12 , for example, is read as: the squared amplitude for the emission of two energy-ordered gluons, k 1 and k 2 , with gluon k 1 real and gluon k 2 virtual. Notice that, in the eikonal approximation, the squared amplitude for the softest gluon being virtual is simply minus the squared amplitude for it being real. In other words: where x could either be R or V. At one-and two-loops, for example, one has: (2.7) The phase space factor for the emission of m gluons is: where g s = √ 4πα s andᾱ s = α s /π. The factor m i=1 k 2 ti /2g 2 s multiplies the squared amplitude W X 12···m to produce the purely angular squared-amplitude W X 12···m (i.e., W X 12···m depends only on η and φ variables and the corresponding colour factor). In other words, one may write: (2.10) The non-linear "measurement operator"Û m acts on the squared amplitudes W X 12···m and plays the role of excluding gluon emission events for which the hemisphere mass is greater than ρ. It is not, however, equivalent to a simple heaviside step function Θ(ρ − i x i e −η i ), since it requires non-numerical input (information about the real-virtual nature of the various gluons). Due to strong ordering, the measurement operatorÛ m factorises into a product of individual measurement operators;Û m = m i=1û i . The squared amplitudes W X 12···m are eigenfunctions of the measurement operatorsû i with eigenvalues 0 or 1 such that: • if gluon k i is real and outside H R thenû i W X 12···m = W X 12···m , This means that events with real emissions inside the measured hemisphere and which contribute more than ρ to the hemisphere mass are excluded (i.e., not integrated over). That is, Σ(ρ) represents the probability that the measured hemisphere mass be less than ρ , as is expressed in eq. (2.4). We therefore write the measurement operator as: , and vice versa. In the above we used the relations Θ in i + Θ out

One-loop calculation and the Sudakov exponentiation
Having properly set up the form of the integrated distribution, we can now proceed with the calculation at each loop level. To warm up for higher loops we start with the one-loop case. At one-loop the squared amplitude for the emission of a single real gluon, multiplied by the corresponding phase space, is given by: The corresponding virtual contribution is The measurement operator reads: which when acting on the squared eikonal amplitudes yields: means that the real and virtual contributions completely cancel out in sufficiently inclusive cross-sections). Substituting into the expression of Σ 1 (e.q. (2.5) with m = 1), we are left with the uncancelled integration: where the step function Θ ρ 1 restricts x 1 to be greater than ρ e η 1 , and since x 1 < 1 then one also has the restriction on the rapidity such that η 1 < L, with L = ln(1/ρ). Performing the integration to single logarithmic accuracy we find: We note that the leading logarithms in the hemisphere mass distribution are double logarithms, which originate from soft and collinear (to the direction of the outgoing quark) singularities of the squared amplitudes for primary gluon emissions off the initiating hard qq pair. It has long been known that the resummed distribution accounting for these primary emissions (or global logarithms) to all-orders is entirely generated from the leading-order result by simple exponentiation (Sudakov form factor). However, Dasgupta and Salam [24] showed that a new class of large single logarithms appears starting at two gluons emission, and which they termed NGLs. We may thus express the resummed hemisphere mass distribution as follows: where Σ P is the primary Sudakov form factor, and Σ NG is the resummed non-global factor, Our aim in this paper is to compute the non-global functions Σ NG m (ρ) for m = 2, 3, 4 and 5.

Two-loops calculation and non-global logarithms
The various squared-amplitudes for the emission of two energy-ordered (real or virtual) gluons with x 2 ≪ x 1 ≪ 1 are expressed as (see for instance [44]): where the irreducible term W RR 12 is: and the measurement operator in this case is given by: Acting on the squared amplitudes yields: The sum of these terms then gives: where we used the squared amplitudes from eqs. (2.20a) and (2.20b). Substituting into the expression of Σ 2 (ρ) we obtain: where we introduced the shorthand notation dΠ 12···m = m i=1 dΦ i Θ ρ i . In the first integral in the right-hand-side of eq. (2.25), the integrand is symmetric under the exchange k 1 ↔ k 2 , which means that we can relax the condition x 1 > x 2 and divide the integral by a factor 2!. Hence this integral factors out into the product of two separate identical contributions from gluons k 1 and k 2 , which both have exactly the same form as the one-loop result (eq. (2.15)). Thus we find for this term: which is just the expansion of the Sudakov Σ P (ρ) at second order. Hence: The latter expression is the pure non-global contribution at this order. It is given by: , since no collinear (double logarithms) are present for the pure non-global contribution. Hence the x integration easily factors from the rapidity integration and we just set the lower limit on x 2 to ρ (to single logarithmic accuracy). Performing the trivial integration over x i we obtain the result L 2 /2!. We note that at n th order we have: Using the result of integration over φ 2 and η 2 from eq. (A.6a) of appendix A (with {j, m} → {1, 2}), we obtain: whereL =ᾱ s L and ζ is the Riemann-Zeta function. This is exactly the result obtained by Dasgupta and Salam [24] for NGLs at two-loops. To the best of our knowledge the analytical calculation of NGLs at finite N c beyond this order has not been performed before, and it is this very task that we do in the next section for the first time in the literature.
3 Non-global logarithms beyond leading order

Three-loops calculation
Having reproduced the well-known result for NGLs at leading order (two-loops), we proceed to compute NGLs at finite N c at next-to-leading order, namely triple gluons emission. As usual we begin by the simplification of the measurement operator which will help us identify both angular and real-virtual configurations giving rise to large logarithms: whereˆ U 3 is the collection of terms which when acting on the squared amplitudes W X 123 , with X summed over, yields a zero. The action of the measurement operator on the various squared-amplitudes summed over X gives: As stated in the introduction, the explicit expressions for the various squared amplitudes above (together with those at higher loops) will be presented in our forthcoming work [39].
Here we restrict ourselves to showing the simplification of the above squared-amplitudes in terms of the antenna functions defined previously in eq. (2.2). We have: In eq.
The first contribution is: which is simply the expansion of the Sudakov Σ P (ρ) at this order. The factor 1/3! accounts for the fact that the integrand in eq. (3.4) is completely symmetric under the exchange of gluons, which means that the condition x 1 > x 2 > x 3 can be relaxed and the result multiplied by 1/3!. The integral is then factored out into the product of three identical integrals, each of them resembling the one-loop result eq. (2.15). The second contribution to Σ 3 (ρ) is: By swapping k 1 ↔ k 2 in the second integral of eq. (3.5), and performing the successive permutations: k 1 ↔ k 2 then k 1 ↔ k 3 in the third integral of the same equation, Σ B 3 (ρ) becomes: The three integrands in eq. (3.6) are identical except for the region of integration over transverse momenta fractions. Thus we unify them into a single integral with the region of integration expressed by: Hence we write: factors out into a product of the one-loop primary cross-section and the twoloop NGLs cross-section (eqs. (2.15) and (2.27b) respectively). This result is expected from the expansion of the Sudakov form factor at one-loop times the leading NGLs. Said differently, Σ B 3 (ρ) is just an "interference" term related to previous orders. The remaining term Σ C 3 (ρ), which is the pure irreducible NGLs contribution at this order, Σ NG 3 (ρ), is proportional to C F C 2 A and given by: (3.9) From the above equation one sees that the new irreducible NGLs contribution at three-loop order is generated by two mechanisms: (a) the energy-ordered real gluons k 1 and k 2 outside H R coherently emit the softest gluon on the other hand. Both of these contributions are not related to previous orders. It seems, at first inspection, that the second mechanism mentioned above (particularly the case where gluon k 2 is inside H R , as clearly shown in eq. (3.10) below) is in contradiction with the common picture about the origin of NGLs. The latter picture dictates that, to all-orders, NGLs are entirely generated from a soft emission into H R that is coherently radiated by arbitrary ensembles of soft, but harder, large-angle energy-ordered gluons outside H R [24,25]. Nonetheless, and though NGLs contribution from the said mechanism comes from both gluons k 2 and k 3 inside H R , gluon k 2 is actually virtual. We shall see later that this mechanism persists at higher loops too. Hence whenever a contribution to NGLs comes from configurations whereby gluons other than the softest are inside H R , then these "other" gluons must be virtual. 6 In fact the contributions of the two terms W RVR 123 and W RRR 123 in eq. (3.9) are separately divergent but their sum is finite. The integral (3.9) can be expressed as a sum of two finite terms in the following way: where we performed the trivial integration over transverse momenta fractions to obtain L 3 /3! and used the results of rapidity and azimuthal integrations shown in appendix A. The terms A 23 q1 (η 1 ) and A 23 1q (η 1 ) are given in eqs. (A.7b) and (A.7c) with {i, j, m} → {1, 2, 3}. The integration over η 1 in the first line of eq. (3.11) yields the result 8ζ 3 and that in the second line gives 4ζ 3 . Thus the pure non-global contribution at this order reads: Hence, up to this order we have: It is intriguing to note that the coefficients of NGLs at finite N c for two and three-loops (i.e., ζ 2 and ζ 3 ) are identical to those found at large N c [26,30]. 7 This is due to the fact that the combination of real/virtual squared amplitudes strangely produces identical integrands in the expressions of the NGLs contributions Σ NG 2 and Σ NG 3 (eqs. (2.27b) and (3.10) respectively). It would have been tremendously easy to resum NGLs for the hemisphere mass distribution to all-orders if the pattern in (3.13) persisted at higher loops. Although, we shall encounter Zeta functions at higher loop orders, the pattern itself unfortunately breaks down starting at four-loops, as we shall see in the next subsection.

Four-loops calculation
For four gluons emission, the computation of NGLs for the hemisphere mass distribution proceeds in an analogous manner to that of two and three gluons emission. At this (fourloops) order the measurement operator reads: whereˆ U 4 is the sum of all terms which when operate on the squared amplitudes W X 1234 and X is summed over give zero. Acting by the measurement operator on the various squared amplitudes and summing over configurations we obtain: • the softest gluon is always inside H R , • each real gluon k i is associated with a step function Θ out i , • virtual gluons are associated with neither Θ in nor Θ out (i.e., they can either be in or out of H R ).
The hemisphere mass distribution at fourth order may then be cast in the form (2.5). Substituting the various matrix-elements squared we can split the hemisphere mass distribution, Σ 4 (ρ), into five parts: where terms in the last line are the four-loop irreducible components of the squared amplitudes for the corresponding gluon configurations. The parts Σ A 4 , Σ B 4 , Σ C 4 , and Σ D 4 completely reduce to integrals we calculated at previous orders, while the remaining Σ E 4 part is the new NGLs contribution. Let us evaluate each of these integrals separately starting with the reducible parts.

Reducible parts
For the first part Σ A 4 we can, as usual, relax the condition x 1 > x 2 > x 3 > x 4 and multiply the result by a factor of 1/4!. This part then factors out into the product of four identical integrals of the form we met at O(α s ) (eqs. (2.15) and (2.16)), thus we obtain for this term: which is just the expansion of the Sudakov at the fourth order. The second part Σ B 4 is carried out in a fashion analogous to that of Σ B 3 in eq. (3.5). The five integrands are transformed into the first integrand, W RR 12 W R 3 W R 4 , with the appropriate changes in the integration limits. Thus we have six integrals having identical integrands but different regions of integrations. Writing these integration regions as step functions and simplifying we obtain Θ(x 1 − x 2 )Θ(x 3 − x 4 ). This means that Σ B 4 factors out into the product of two integrals, one over k 1 and k 2 and the other over k 3 and k 4 , as follows: where we have used eqs. (2.26) and (2.27b) to arrive at the second line of the above equation. Eq. (3.18) is in fact the interference of the expansion of the Sudakov form factor Σ P with NGLs at two-loops Σ NG 2 . Performing the integrations in the third part Σ C 4 along the same lines outlined above for Σ B 4 (i.e., by making appropriate changes of variables) we obtain: where we used eqs. (2.15) and (3.9). This result is the interference between the expansion of the Sudakov and Σ NG 3 . The part Σ D 4 can be written as follows: which indicates a possible pattern of exponentiation of NGLs since this term resembles the structure of the expansion of exp{Σ NG 2 } at this (fourth) order.

Irreducible part
The irreducible part at four-loops, Σ E 4 (ρ), is given by: ) even breaking the colour pattern. Therefore the loss of the pattern of NGLs is in fact a manifestation of the break in the structure of the eikonal amplitudes. This might be related to the failure of the "probabilistic scheme" discussed by Dokshitzer et al. in ref. [32], where such (irreducible) contributions to the eikonal squared-amplitude were dubbed "monster" terms. They were traced back to be originating from the "colour polarisability" of jets [32]. 8 Moreover, we note here, as can be seen from the form of Σ E 4 , that contributions to NGLs at this order are generated when the softest gluon k 4 is emitted inside the measured hemisphere H R , whilst the hardest gluon k 1 is always outside. The other two gluons, k 2 and k 3 , may be emitted inside H R provided they are virtual (in accordance with the observation made in the previous subsection).
The contribution Σ E 4 is not related to the expansion of the Sudakov nor to NGLs at previous orders, and represents the new non-global contribution at four-loops. Together with Σ D 4 they form the total non-global contribution, Σ NG 4 , to the hemisphere mass distribution at this order. To evaluate the part Σ E 4 we first integrate over transverse momenta fractions, which as usual yields L 4 /4!, and then perform the azimuthal and rapidity integrations according to the angular configurations indicated in eq. (3.21). The azimuthal integrations are not as straightforward as at three-loops and we find the method used in ref. [30] of contour integration very useful in reducing the number of integrals to be performed. The reduced integrals are then carried out analytically whenever possible, otherwise numerically. Since all analytical integrations yield results that are explicitly proportional to ζ 4 , the resultant values from numerical integrations were interpreted in terms of ζ 4 . In addition to this semi-numerical approach we also verify our results by numerically integrating each of the finite terms in the following form of Σ E 4 (ρ): over the full (7-dimensional) phase space using the multi-dimensional numerical-integration library Cuba [45]. The final result reads: which may also be rewritten in the following two alternative forms: The expression (3.24a) explicitly shows the finite-N c correction to the large-N c result, while (3.24b) emphasises the pattern C F C A n ζ n+1 seen at two and three-loops. It also reveals that even though Σ E 4 is a new irreducible contribution at four-loops, it still contains factors related to lower-order NGL contributions (the term (C F C A ζ 2 ) 2 ). Observe that the size of the finite-N c correction in (3.24a) is about ∼ 1.5% that of the large-N c result. This is in agreement with the conclusion arrived at in [46] for the impact of finite-N c corrections at all-orders for e + e − processes. The total non-global contribution at this order is then given by: and thus the hemisphere mass distribution up to this order is expressed as: In the next subsection we discuss the five-loops case and the possibility of resummation of NGLs.

Five-loops and beyond
Following the same steps as before we write the measurement operator at five-loops as follows: where, as usual,ˆ U 5 is the sum of all terms that yield vanishing contributions to the hemisphere mass distribution. Acting by the measurement operator on the various squared amplitudes and summing over configurations we obtain: . (3.28) The hemisphere mass distribution at five-loops is then given in eq. (2.5) with m = 5.
Substituting the various matrix-elements squared and following the same procedure outlined at four-loops we again obtain two types of contributions; reducible, Σ r 5 , and irreducible, Σ irr 5 . The former contains all the interference terms between the Sudakov factor Σ P and NGLs at previous orders as well as interference terms between two and three-loops NGLs. Explicitly written it reads: Note that the penultimate term in the above equation contains the contribution Σ P 1 × (Σ NG 2 ) 2 /2! (eqs. (3.25) and (3.20)). The irreducible contribution Σ irr 5 is expressed as: , which all contain "monster" terms.
The calculation of Σ irr 5 turns out to be trickier and more involved than anticipated. In particular we have not yet been able to simplify (like we did at four-loops) the monster parts of the irreducible amplitudes of the second type above, to forms that can readily be integrated. Other than the monster parts, all remaining terms (either of the first or second type above) are in fact integrable.
Till the full expression of the irreducible contribution is simplified and integrated, we express the result of Σ irr 5 in the following form (based on the pattern seen at four-loops (3.24b) and the pieces found in the integrable amplitudes of eq. (3.30)): where a = ζ 5 /ζ 2 ζ 3 ≃ 0.5244 and the constant coefficients α and β are yet to be determined. We discuss the possible values of these constants when we compare our results with those at large N c in the next section. The form (3.31a) explicitly shows the finite-N c correction. The total NGLs contribution at five-loops then reads: The results we obtained up to five-loops, particularly eq. (3.29), in fact suggest a possible resummation of NGLs into an exponential function of the form: Eq. (3.33) may actually be rewritten in a form analogous to that found in ref. [47] (eqs. (5.10) and (5.11)) for clustering logarithms. To this end we write: where The above-mentioned similarity between clustering logarithms and NGLs emphasises the common (non-global) origin of the two types of logarithms. Moreover, following the pattern in (3.31b), eq. (3.33) may also be recast into the form: The expansion of the above exponential exactly reproduces the terms we have calculated up to five-loops including all interference terms in the distribution. At each higher order one simply adds a new irreducible NGLs term in the exponent. In fact, if the pattern deduced in eqs. (3.24a) and (3.24b) persists at higher-loop orders, then one can put forth the following ansatz for the general form of the n th order contribution to the exponent of the resummed NGLs factor: where ⌊n⌋ represents the floor function of n, and γ n , σ n ,γ n , andσ k are constant coefficients to be determined from integrations. The two formulae presented above for the ansatz are equivalent up to the constant coefficients. The first form stresses the finite-N c correction while the second preserves the pattern seen at two-, three-and four-loops (eq. (3.24b)).
The above formulae may only be verified once higher-loop orders are carried out explicitly. We hope to perform such calculations in the near future. We note that in the exponent of (3.33) (and (3.36)) the series inL has alternating signs at each escalating order. To assess the relative size of the three and four-loops corrections to the leading two-loops result, we plot in figure 2 the ratio Σ NG / exp(Σ NG 2 ) for various truncations of the series in the exponent in eq. (3.33). The leading NGLs coefficient seems to dominate for only relatively small values ofL (L 0.15). For larger values the series seems to depart from the leading term in an alternating way (towards larger (smaller) values for odd (even) loop orders). These significant variations mean that the terms computed thus far are insufficient to capture the full behaviour of the all-orders resummed distribution. We expect, however, that adding few more terms in the exponent may lead to a convergent and more stable behaviour, since one could argue that higher-order terms are suppressed byL n /n!, while ζ n saturates at 1 as n becomes larger.
To cross-check our results, eqs. (3.33) and (3.36), we compare them, in the next section, to previous calculations at large N c both at fixed order and to all-orders.

Comparison with analytical results at large N c
Having calculated the coefficients of NGLs at finite N c fully up to four-loops and partially at five-loops, we can now compare our findings to those of Schwartz and Zhu [30] obtained through the analytical solution to the BMS equation [26] in the large-N c limit. To go from the finite-N c case to the large-N c approximation we simply invoke the replacement . This is equivalent to expanding C F to first order in colour: Then the expansion of our result (3.36) to leading order in colour, i.e., at large N c , up to five-loops is: where we have written the explicit values of ζ 2 = π 2 /6 and ζ 4 = π 4 /90. The full result reported by Schwartz and Zhu (SZ) at large N c , Σ NG SZ , is [30]: where L is simply N cL . The two results are thus identical up to four-loops. Recalling that they were arrived at using different approaches, their equality provides a solid cross-check of the correctness of the computed NGLs coefficients (at least up to four-loops). Even though we are unable to fully compare our result to that of Schwartz and Zhu at five-loops, due to the missing values of α and β, it is ironic to note that the pattern spotted at four-loops, eq. (3.24b), seems to hold true at five-loops. Whilst the term ζ 5 is apparent in (4.3), the product ζ 2 ζ 3 is disguising in the factor π 2 ζ 3 /360. A quick comparison reveals the values: Given the above values we expect the finite-N c result of NGLs up to five-loops to be expressed as: We are not, however, claiming to have fully accounted for NGLs at this order since we have not explicitly calculated the coefficients of NGLs with colour factors C 2 F C 3 A and C F C 4 A . We hope that further research on this can help verify the above equation (and eq. (3.37) in general) in which case it may actually be possible to find a key to the analytical resummation of NGLs both at large and finite N c .
Our approach additionally has the benefit that it sheds light on the possibility of assessing the validity of the large-N c approximation, by means of judging the impact of neglected finite-N c corrections. An important note in this regard is that at two-and threeloops, as can be seen by comparing eqs. (3.26) and (4.2) at O(L 2 ) and O(L 3 ), there are no hidden terms buried by the large-N c approximation, and the finite-N c result can simply be obtained from the solution of the BMS equation by just restoring the full colour factors through: N 2 c → 2C F C A (at two-loops) and N 3 c → 2 C F C 2 A (at three-loops). At four-loops this is not true and in fact there is a hidden correction that is given plainly in (3.24a). We regard this as the first-order proper finite-N c correction which introduces new terms that are entirely absent at large N c . The second-order proper finite-N c correction occurs at five-loops and is shown in (3.31a).

Comparison with all-orders numerical results
In order to verify our resummed formula (3.33), and even (4.5) which includes the fiveloops term, it is instructive to compare it to the all-orders numerical solution of the finite-N c Weigert equation [28]. We have not, unfortunately, been able to obtain the output of the MC program, written by Hatta and Ueda [27], for the hemisphere mass distribution. 9 We thus postpone this discussion till the said numerical distribution becomes available. Furthermore, to assess the importance of the missing higher-loop terms in (3.33) we compare it to either the results obtained by the numerical solution of the BMS equation [26] or to the output of the numerical MC program of Dasgupta and Salam (DS) [24]. The latter two numerical solutions are in fact identical within a percent accuracy [26,30], and we thus restrict ourselves to the DS MC program.
Let us introduce the standard evolution parameter t [24], which accounts for the running of the coupling: The hemisphere mass distribution has not yet been coded into the MC program [46]. where β 0 is the one-loop coefficient of the QCD β function. At fixed order one has t = α s L/2π =L/2. Hence substitutingL by 2 t into eq. (3.33) up to four-loops we find: We compare the result (4.7) with the parametrisation for NGLs to all-orders obtained in ref. [24] by fitting to the output of the aforementioned DS MC program [24]: (4.8) In figure 3 we plot our approximate resummed result (4.7) for various truncations along with the DS resummed factor (4.8) for a range of t ∈ [0, 0.5] both at large (left) and finite (right) N c . Recall that a value of t = 0.3 corresponds to a value of L = 19 and ρ ∼ 10 −8 for α s ∼ 0.1, which is sufficient for phenomenological purposes. Few points to note from the plots. Firstly, as expected for finite N c , all curves are shifted up due to the fact that one is using C F = 4/3 ≃ 1.33 instead of C F = C A /2 = 3/2 = 1.5 (recall that C F is in the exponent). Secondly, it is striking to observe that the best approximation to the all-orders result for quite a large range of t is the leading two-loops result, 10 exp Σ NG 2 , for both large and finite-N c cases. This suggests that the alternating, positive and negative, higher-loop contributions to Σ NG somehow balance out.
Moreover, the main feature of the plots and which has direct link to the purpose of this paper is actually seen at small values of t. We see that the interval of t over which the four-loops result and the all-orders resummed factor overlap is 0 ≤ t 0.12. This interval of overlapping is smaller for three-loops, 0 ≤ t 0.08, and even smallest for two-loops 0 ≤ t 0.05. The latter feature may be seen more clearly in figures 4 and 5. One would therefore expect that adding more terms in the exponent of (3.33) leads to increasingly larger intervals of overlapping.  A similar observation was also made in ref. [29] for "filtering analysis" in the case of the filtering parameters n filt = 2 and η filt = 0.1, 0.3 (figure 18 of ref. [29]), as well as n filt = 3 and η filt = 0.3 ( figure 21). However the author of ref. [29], having plotted the expansion of the full filtered Higgs-jet mass distribution including both primary and non-global logarithms in the Cambridge-Aachen jet algorithm [48,49], ascribed the convergence of the series, as one adds higher-loop terms, to the dominance of the primary series. The author then verified this explanation by plotting the same distribution for higher values of η filt ( figure  19) where collinear logarithms are expected to be absent and NGLs become of the same order as primary logarithms. Two issues to point out regarding our work compared to that of ref. [29]: firstly, we are plotting purely the NGLs resummed exponential factor and hence the convergence seen in figures 4 and 5 has nothing to do with primary logarithms. Secondly, it is well known [43,[50][51][52][53] that employing the Cambridge-Aachen jet algorithm not only reduces the size of NGLs but also introduces clustering logarithms that are as important as NGLs. Thus plotting the full distribution, which includes primary, non-global and clustering logarithms, would not tell much about the convergence of the NGLs series.
Moreover, the author of ref. [29] also plotted (figure 24) the pure NGLs resummed factor for the interjet energy flow distribution and concluded that, up to six-loops, the NGLs series seems to be divergent. Recalling that the coefficient of the two-loops NGLs depends on the rapidity gap ∆η [25,43], it is likely that higher-loop NGLs coefficients depend on ∆η too. The divergence may thus be due to the presence of the ∆η terms. For the hemisphere mass distribution that we have treated in this paper there is no such rapidity gap dependence. A proper answer, however, may only be given once the formerinterjet energy flow -distribution is carefully considered, a task which we hope to perform in coming publications.
Notice finally that the so far discussed NGLs behaviour is in contrast to that seen for clustering logarithms [47] where the whole structure of the all-orders result is mostly captured by the first few terms in the exponent.

Conclusions
In this paper we have considered the calculation of the leading NGLs at finite N c up to five-loops for the hemisphere mass distribution in e + e − → di-jet events. We performed the calculation by means of integrating the squared amplitudes for the emission of energyordered soft gluons in the eikonal approximation, valid at single logarithmic accuracy, over a suitable phase space achieved through a measurement operator. The two and three-loops results were shown to be relatively straightforward to obtain, and were found to be directly related to the Riemann-Zeta function. We noticed that, up to this loop order, finite-N c corrections are absent. This was a direct consequence of the relatively simple structure of the eikonal amplitudes (up to this order) as well as the combination of real/virtual amplitudes induced by the phase space measurement operator. Within the same eikonal framework we computed the four-loops contribution to NGLs distribution. The latter turns out to be much harder than the previous two orders and the simple result in terms of a product of a single colour factor and a Zeta function breaks down. This failure originates from a break in the simple structure of the corresponding real-virtual eikonal amplitudes at four-loops order, a phenomenon that was noticed more than two decays ago [32]. Nevertheless, we were able to overcome this complexity, compute NGLs and even spot a new pattern for NGLs at and beyond this order. This pattern helped us to successfully write down the five-loops contribution to NGLs up to constant coefficients which we extracted from comparisons to previous large-N c results [30]. We hope to be able to fully compute these constants in the near future. The five-loops calculation reveals that the NGLs distribution seems to exhibit a pattern of exponentiation. To this end the said distribution was cast in an exponential form with full NGLs coefficients and colour factors up to four-loops in the exponent.
Comparisons to large-N c results obtained by other authors [30] confirmed our findings, at least in the latter limit. We then took the step forth and compared our exponential function to the all-orders numerical resummed result reported by Dasgupta and Salam [24]. To our surprise, the shape of the all-orders result was best represented by the twoloops approximation for a wide range of the evolution parameter t. In the region of small t, however, adding more terms in the exponent of our resummed result yielded better agreement, then the two-loops result, with the all-orders numerical result. This suggests that more higher-loop contributions are needed for our result to be of any phenomenological significance (i.e., till the agreement extends to values of t up to ∼ 0.2 -0.3). The task of computing these higher-loop terms might not be impossible after all given that we have developed, over the course of preparing this paper, the machinery for: computing eikonal amplitudes at finite N c to theoretically any loop order, reducing the dimension of the phase space over which to integrate, and spotting a pattern for NGLs at each order.
The φ i -azimuthal averaging over the inverse of the angular function (ij) is given by: In the case where a gluon k m is constrained within the measured hemisphere region and gluons k i and k j are constrained outside, the k m angular integration over w m ij yields [30]: With the same conditions (η i < 0, η j < 0 and η m > 0) we can perform the azimuthal and rapidity integrations for the following antenna functions: We can also perform further integrations over φ j and η j : with Li 2 the polylogarithm function of order 2. Similarly integrating over the angles of the softest particle k n in the antenna A mn ij yields: A mn ij =w m ij (ln(1 + im ) + ln(1 + jm ) − ln(1 + ij )) .
Then we have: × sinh |η 1 − η 3 |csch|η 2 − η 3 | − cosh(η 1 − η 2 ) sinh(η 1 − η 2 ) ln coth η 1 − 1 coth η 2 − 1 + B A note on NGLs-ζ n relation As a byproduct, we notice from eqs. (2.30) and (3.11) that one may define the following possibly "new" logarithmic-integral representation for the Riemann-Zeta function: where Γ(s) = (s − 1)! is the Gamma function and the variable s is greater than 1. The above formula is valid if s is an integer. In the case of non-integer real values one has to take the modulus of the right-hand-side in eq. (B.1). In terms of the polar variables (θ, φ) the Zeta function admits the integral formula: where c ≡ cos θ. Notice that the form (B.2) seems to fail (in Mathematica 9) for s > 10. Moreover, it is interesting to note that the first non-divergent value of the Zeta function is for s = 2, and so is the first non-vanishing coefficient of NGLs. If we let S s denote the NGLs coefficient at the s th loop order then we can write S s as the Mellin transform of the function (e η − 1) −1 [54]. That is: which is at least true for two and three-loops NGLs coefficients. Recall that the Mellin transform techniques were employed in ref. [55] to compute the first resummed result for event-shape distributions. Whether there exists a more profound relation between NGLs and the Zeta function (and its related functions such as the polylogarithms) is a subject that requires further investigations.