NLL${'}$ resummation of jet mass

TL;DR

This paper develops a resummation framework for jet mass and energy distribution observables using effective field theory, including NLO corrections, improving data agreement and reducing uncertainties.

Contribution

It introduces a novel resummation method for non-global observables with full NLO corrections within an effective field theory framework.

Findings

Enhanced agreement with experimental data

Reduced theoretical uncertainties

First step towards higher-logarithmic resummation

Abstract

Starting from a factorization theorem in effective field theory, we present resummed results for two non-global observables: the invariant-mass distribution of jets and the energy distribution outside jets. Our results include the full next-to-leading-order corrections to the hard, jet and soft functions and are implemented in a parton-shower framework which generates the renormalization-group running in the effective theory. The inclusion of these matching corrections leads to an improved description of the data and reduced theoretical uncertainties. They will have to be combined with two-loop running in the future, but our results are an important first step towards the higher-logarithmic resummation of non-global observables.

Figures

Equations

$$\sigma (Q,Q_0)$$

$$\frac{d\sigma }{d\rho }=\frac{1}{2}\left(\frac{d\sigma }{d{\rho }_{\ell }}+\frac{d\sigma }{d{\rho }_h}\right).$$

$$\rho Q^2=M^2=(p_{\bar{c}}+p_s{)}^2=p_{\bar{c}}^2+Q\bar{n}\cdot p_s+\mathcal{O}(p_s^2).$$

$$\tilde{\sigma}(\tau)=\sum_{i=q,\bar{q},g}\tilde{j}_{i}(\tau Q,\mu)\,\sum_{m=1}^{\infty}\big{\langle}\bm{\mathcal{H}}^{i}_{m}(\{\underline{n}\},Q,\mu)\otimes\bm{\mathcal{\widetilde{S}}}_{m}(\{\underline{n}\},\tau,\mu)\big{\rangle}\,,$$

$$\displaystyle\sigma(Q,Q_{0})=\sum_{l=2}^{\infty}\big{\langle}\bm{\mathcal{H}}_{l}(\{\underline{n}^{\prime}\},Q,\mu_{h})\otimes\sum_{m\geq l}^{\infty}\bm{U}_{lm}(\{\underline{n}\},\mu_{s},\mu_{h})\,\hat{\otimes}\,\bm{\mathcal{S}}_{m}(\{\underline{n}\},Q_{0},\mu_{s})\big{\rangle}\,,$$

$$\mathbf{U}(\{\underline{n}\},{\mu }_s,{\mu }_h)=\mathbf{P}\exp \left[{\int }_{{\mu }_s}^{{\mu }_h}\frac{d\mu }{\mu }{\mathbf{\Gamma }}^H(\{\underline{n}\},\mu )\right].$$

$${\int }_{{\mu }_s}^{{\mu }_h}\frac{d\mu }{\mu }{\mathbf{\Gamma }}^H={\int }_{{\alpha }_s({\mu }_s)}^{{\alpha }_s({\mu }_h)}\frac{d\alpha }{\beta (\alpha )}\frac{\alpha }{4\pi }{\mathbf{\Gamma }}^{(1)}=\frac{1}{2{\beta }_0}\ln \frac{{\alpha }_s({\mu }_s)}{{\alpha }_s({\mu }_h)}{\mathbf{\Gamma }}^{(1)}\equiv t{\mathbf{\Gamma }}^{(1)}.$$

$$\bm{\Gamma}^{(1)}=\left(\begin{array}[]{ccccc}\,\bm{V}_{2}&\bm{R}_{2}&0&0&\ldots\\ 0&\bm{V}_{3}&\bm{R}_{3}&0&\ldots\\ 0&0&\bm{V}_{4}&\bm{R}_{4}&\ldots\\ 0&0&0&\bm{V}_{5}&\ldots\\ \vdots&\vdots&\vdots&\vdots&\ddots\\ \end{array}\right).$$

$${\mathbf{V}}_m$$

$${\mathbf{R}}_m$$

$$W_{ij}^k=\frac{n_i\cdot n_j}{(n_i\cdot n_k)(n_j\cdot n_k)},$$

$$\displaystyle\sigma^{\rm LL}(Q,Q_{0})=\sum_{m=2}^{\infty}\big{\langle}\bm{\mathcal{H}}_{2}(\{n_{1},n_{2}\},Q,\mu_{h})\otimes\bm{U}_{2m}(\{\underline{n}\},\mu_{s},\mu_{h})\,\hat{\otimes}\,\bm{1}\big{\rangle}\,.$$

$${\mathcal{H}}_2$$

$${\mathcal{S}}_m$$

$$\frac{{\sigma }^{\mathrm{L}{\mathrm{L}}^{\prime }}(Q,Q_0)}{{\sigma }_0}=$$

$$\displaystyle+\frac{\alpha_{s}(\mu_{h})}{4\pi}\sum_{m=2}^{\infty}\big{\langle}\bm{\mathcal{H}}_{2}^{(1)}(\{n_{1},n_{2}\},Q,\mu_{h})\,\otimes\,\bm{U}_{2m}(\{\underline{n}\},\mu_{s},\mu_{h})\,\hat{\otimes}\,\bm{1}\big{\rangle}$$

$$\displaystyle+\frac{\alpha_{s}(\mu_{h})}{4\pi}\sum_{m=3}^{\infty}\big{\langle}\bm{\mathcal{H}}_{3}^{(1)}(\{n_{1},n_{2},n_{3}\},Q,\mu_{h})\,\otimes\,\bm{U}_{3m}(\{\underline{n}\},\mu_{s},\mu_{h})\,\hat{\otimes}\,\bm{1}\big{\rangle}$$

$$\displaystyle+\frac{\alpha_{s}(\mu_{s})}{4\pi}\sum_{m=2}^{\infty}\big{\langle}\bm{\mathcal{H}}_{2}^{(0)}(\{n_{1},n_{2}\},Q,\mu_{h})\,\otimes\,\bm{U}_{2m}(\{\underline{n}\},\mu_{s},\mu_{h})\,\hat{\otimes}\,\bm{\mathcal{S}}_{m}^{(1)}(\{\underline{n}\},Q_{0},\mu_{s})\big{\rangle}.$$

$$⟨{\mathcal{H}}_2(\{n_1,n_2\},Q,\mu )\otimes {\mathcal{S}}_2(\{n_1,n_2\},Q_0,\mu )⟩={\sigma }_0{H}_2(Q^2,\mu )⟨{\mathcal{S}}_2(\{\bar{n},n\},Q_0,\mu )⟩,$$

$$H_2(Q^2,\mu )=1+\frac{{\alpha }_s}{4\pi }{C}_F\left[-8{\ln }^2\frac{\mu }{Q}-12\ln \frac{\mu }{Q}-16+\frac{7}{3}{\pi }^2\right],$$

$${\hat{\theta }}_2\equiv \tan \frac{{\theta }_2}{2}=uv,{\hat{\theta }}_3\equiv \tan \frac{{\theta }_3}{2}=v,$$

$${\mathcal{H}}_{3,\mathrm{I}}^{(1)}(u,v,Q,\mu )=$$

$$+8\delta (u){\left(\frac{\ln v}{v}\right)}_{+}+\left[-\ln \frac{\mu }{Q}F(u,0)+\frac{2u^2}{(1+u{)}^3}-F(u,0)\ln (1+u)\right]\delta (v){\left(\frac{1}{u}\right)}_{+}$$

$$\displaystyle+F(u,0)\delta(v)\left(\frac{\ln u}{u}\right)_{+}+F(u,v)\left(\frac{1}{u}\right)_{+}\left(\frac{1}{v}\right)_{+}\Bigg{\}}\Theta_{\rm in}(v).$$

$$F(u,v)=\frac{4\left[u\left[-2\left(u^2+u+1\right)v^2+u(2u(u+1)+1)v^4+u+2\right]+2\right]}{(u+1{)}^3}.$$

$${\mathcal{H}}_{3,\mathrm{II}}^{(1)}(u,v,Q,\mu )=$$

$$\displaystyle\hskip 34.14322pt+G(u,v)\left(\frac{1}{v}\right)_{+}\Bigg{\}}\Theta_{\rm in}(v),$$

$$G(u,v)=\frac{4\left[u\left[-2\left(u^2+u+1\right)v^2+u(u(u+2)+2)v^4+2(u+1)\right]+1\right]}{(u+1{)}^3}.$$

$${\mathcal{H}}_{3,\mathrm{III}}^{(1)}(u,v,Q,\mu )=$$

$$H(u,v)=\frac{4v\left(u^4{v}^4+u^2{v}^4+4u^2{v}^2+u^2+1\right)}{(u+1{)}^2\left(1-u{v}^2\right)}.$$

$$\displaystyle\big{\langle}\bm{\mathcal{H}}_{3}^{(1)}(\{\underline{n}\},Q,\mu_{h})\,\otimes\,\widehat{\bm{\mathcal{S}}}_{3}(\{\underline{n}\},\mu_{h})\big{\rangle}=\int_{0}^{1}\!du\int_{0}^{1}\!dv\,\big{\langle}\bm{\mathcal{H}}_{3}^{(1)}(u,v,Q,\mu_{h})\widehat{\bm{\mathcal{S}}}_{3}(u,v,\mu_{h})\big{\rangle}\,,$$

$${\hat{\mathcal{S}}}_3(u,v=0,{\mu }_h)={\hat{\mathcal{S}}}_2({\mu }_h)$$

$${\mathcal{H}}_3^{(1)}\otimes {\hat{\mathcal{S}}}_3={\int }_0^{1}dv\left[A\delta (v)+B(v)+\sum _{i=0}^1{C}_i(v){\left(\frac{{\ln }^{i}v}{v}\right)}_{+}\right]{\hat{\mathcal{S}}}_3(v),$$

$${\int }_0^1\frac{dv}{v}\left[{\hat{\mathcal{S}}}_3(v)-{\hat{\mathcal{S}}}_2\right]={\int }_{v_0}^1\frac{dv}{v}{\hat{\mathcal{S}}}_3(v)+\ln v_0{\hat{\mathcal{S}}}_2+\mathcal{O}(v_0),$$

Full text

††preprint:

CERN-TH-2019-006

1 Introduction

Up to now, higher-logarithmic resummations of collider observables have only been performed for the narrow class of global observables which constrain radiation uniformly over the entire phase space. This category includes very inclusive observables such as selected event shapes, but it excludes all observables with hard phase-space cuts or a fixed number of jets. In recent years, a lot of progress was made in the theoretical analysis of non-global observables Hatta:2008st ; Avsar:2009yb ; Hatta:2013iba ; Caron-Huot:2015bja ; Larkoski:2015zka ; Becher:2015hka ; Neill:2015nya ; Becher:2016mmh ; Larkoski:2016zzc ; Becher:2016omr ; Neill:2016stq ; Becher:2017nof ; Hatta:2017fwr ; Martinez:2018ffw ; Balsiger:2018ezi ; Neill:2018yet . This includes work on the structure of higher logarithms as well as studies of leading logarithms beyond the large-$N_{c}$ limit.

In this paper we start the computation of higher-logarithmic terms for non-global observables by analyzing two simple observables, the jet mass and the interjet energy flow, and presenting resummed predictions which include the full one-loop corrections to the relevant hard scattering processes, as well as the associated jet and soft functions. In the effective-theory framework we use for resummation Becher:2015hka ; Becher:2016mmh , these correspond to matching corrections and they will need to be supplemented by corrections to the renormalization-group (RG) running in the future to arrive at a complete higher-logarithmic treatment of the non-global part.

Our main goal in the present work is to develop the Monte Carlo methods to include these corrections as a step towards full higher-logarithmic resummation, but it is also interesting to study their numerical size, since they have never been computed for non-global observables and often dominate numerically in the global case. It is customary to add a prime to the logarithmic accuracy to indicate the presence of higher-order matching corrections. In this notation our next-to-leading-logarithmic results for the jet mass have NLL′ accuracy.

In Refs. Becher:2016mmh ; Becher:2016omr we have derived a factorization formula for interjet energy flow and light-jet mass. The key element is the presence of multi-Wilson-line operators which generate the intricate pattern of Non-Global Logarithms (NGLs). Explicitly, the result for interjet energy flow at a lepton collider has the form

[TABLE]

where $Q$ is the center-of-mass energy, and $Q_{0}=\beta Q$ is the energy scale above which we veto energy in the gap outside the jet cones. For simplicity, we choose the jet axis along the thrust axis. The above factorization formula neglects power corrections from $\mathcal{O}(\beta)$ terms. The hard functions $\bm{\mathcal{H}}_{m}$ describe hard radiation inside the jet cone, and their characteristic scale is $Q$ since radiation inside the cones is unrestricted. The index $m$ represents the number of hard partons inside the jet, which propagate along the directions $\{\underline{n}\}=\{n_{1},n_{2},\dots,n_{m}\}$. Each of these sources soft radiation, which we describe by a Wilson line along the direction of the hard parton. The matrix elements of these Wilson lines define the soft functions $\bm{\mathcal{S}}_{m}(\{\underline{n}\},Q_{0},\mu)$. To obtain the cross section, one integrates over the directions $\{\underline{n}\}$, which is indicated by the symbol $\otimes$. The hard and soft functions are matrices in the color space of the $m$ partons and one takes the color trace $\langle\dots\rangle$ after multiplying them. The operator definition for these functions and further explanations can be found in Becher:2016mmh .

The second observable we consider is the jet mass distribution at a lepton collider. To define the jet mass, we use the thrust axis to split every event into two hemispheres. One can then (randomly) select one of the two jets and compute its invariant mass $M$, which is usually discussed in terms of the dimensionless variable $\rho=M^{2}/Q^{2}$. Alternatively, one computes the mass in both hemispheres and chooses the heavier mass $\rho_{h}$ or lighter one $\rho_{\ell}$. Obviously, there is a relation among the these observables: the jet mass distribution is simply the average of heavy-jet mass and light-jet mass one

[TABLE]

We will call the hemisphere we select to measure the mass the left one, which means that the radiation in the right hemisphere is unconstrained.111In our previous paper, we called $\rho$ the left-jet mass and denoted it by $\rho_{L}$ Becher:2016omr . We introduce a light-like reference four-vector $n^{\mu}=(1,0,0,1)$ pointing to the right along the thrust axis and an opposite vector $\bar{n}^{\mu}=(1,0,0,-1)$ pointing to the left. The hard partons in the right hemisphere then generate the complicated pattern of soft radiation and associated NGLs. The main difference to formula (1) is that one also needs the standard inclusive jet functions to describe collinear radiation in the left hemisphere. Resummation effects in the jet mass distribution have been discussed in Refs. Dasgupta:2012hg ; Chien:2012ur ; Liu:2014oog ; Idilbi:2016hoa ; Kang:2018jwa , however only in Dasgupta:2012hg the leading NGLs were resummed. Our work is based on the factorization theorem for jet mass derived in Becher:2016omr . The invariant mass of the left jet is obtained from the momentum $p_{\bar{c}}$ of the energetic particles collinear to $\bar{n}$ and the soft partons in the left hemisphere,

[TABLE]

In the factorization theorem, the sum results in a convolution of the soft and jet functions. To avoid this, one can work in Laplace space, where the factorization formula has the product form

[TABLE]

where $\tau$ is the Laplace conjugate variable of $\rho$, and $\tilde{j}_{i}$ is the inclusive jet function Becher:2006qw ; Becher:2010pd , which by now is known to three loops Bruser:2018rad ; Banerjee:2018ozf . In (4) the index $m$ indicates the number of partons in the inclusive (right) hemisphere, so that $m=1$ at leading order (LO).

As long as we consider large jet cone sizes of $\mathcal{O}(1)$, the leading-logarithms (LLs) in interjet energy flow at a lepton collider are of the form $\alpha_{s}^{n}\ln^{n}\beta$. The interjet energy flow is a single logarithmic observable, because collinear logarithms cancel inside the large cone region and only soft logarithms remain. These logarithms arise from the multi-Wilson-line operators $\bm{\mathcal{S}}_{m}$ in (1) and one needs to use parton shower methods to resum the enhanced logarithms already at the LL level. In Balsiger:2018ezi we have written a dedicated parton-shower code to perform the resummation for such observables and have interfaced it with the MadGraph5$\_$aMC@NLO event generator Alwall:2014hca . This provides an automated framework to perform the LL resummation for single-logarithmic observables. However, collider observables are typically double logarithmic. The leading logarithms in the jet mass distribution, for example, are $\alpha_{s}^{n}\ln^{2n}\rho$. Even for non-global observables, these double logarithmic terms have a simple structure, and they can be factored out and treated separately. In the parton shower framework, we therefore subtract these “global” contributions and exponentiate them manually, as Dasgupta and Salam did in their original paper on NGLs Dasgupta:2001sh . Given their different nature, it is interesting to analyze both the interjet energy flow and the jet mass as examples and we will present LL′ and NLL′ improved results for single logarithmic and double logarithmic observables, separately. A second motivation to also analyze the jet mass, is that there are LEP measurements to which we can compare to, in contrast to the interjet energy flow. Unfortunately, the typical jet mass at LEP jet is quite low $M\lesssim 10\,{\rm GeV}$, which translates to a scale of the soft radiation of $Q_{0}\sim M^{2}/Q\lesssim 1\,{\rm GeV}$ so that non-perturbative effects are very important in the peak region of the distribution.

Our paper is organized as follows. In the next section, we will discuss ${\rm LL^{\prime}}$ resummation for interjet energy flow and show how one implements the one-loop corrections to the hard and soft functions. We then move to the jet mass distribution in Section 3, focussing on the differences to the single-logarithmic case. We will in particular show how to subtract global logarithms in the parton shower and in the soft function. After presenting numerical results in Section 4 and comparing to LEP data and PYTHIA results, we conclude in Section 5.

2 Interjet energy flow at ${\rm LL}^{\prime}$ accuracy

The perturbative expansion of the interjet energy flow in (1) suffers from large logarithms of the ratio of the hard scale $Q$ and the soft scale $Q_{0}$. To resum these, one solves the RG equation of the hard function and evolves it from its characteristic scale $\mu_{h}\sim Q$ down to a soft scale $\mu_{s}\sim Q_{0}$. This yields the RG-improved expression Becher:2016mmh

[TABLE]

where the evolution factor is defined as a path-ordered exponential of the anomalous dimension

[TABLE]

The RG-evolution generates additional partons and maps the $l$-parton configuration along the directions $\{\underline{n}^{\prime}\}=\{n_{1},\dots,n_{l}\}$ into an $m$-parton final state along the directions $\{\underline{n}\}=\{n_{1},\dots,n_{l},n_{l+1},\dots,n_{m}\}$. The symbol $\hat{\otimes}$ in (5) indicates the integral over the directions of the additional $m-l$ partons generated in the evolution.

At the leading logarithmic level, we only need the one-loop anomalous dimension and can rewrite the exponent as

[TABLE]

In the last step, we have introduced the evolution time $t\equiv t(\mu_{h},\mu_{s})$. For a given $\mu_{h}$, there is a one-to-one correspondence of the evolution time to the low scale $\mu_{s}$. Obviously, for $\mu_{h}=\mu_{s}$, we have $t=0$. During the evolution, $t$ grows and goes to infinity as $\mu_{s}$ hits the Landau pole. For $\mu_{h}=M_{Z}$ and two-loop running with a Landau pole at $\Lambda=0.230\,{\rm GeV}$, the choice $\mu_{s}=1\,{\rm GeV}$ corresponds to $t=0.08$. A plot connecting $t$ and $\mu_{s}$ for different values of $\mu_{h}$ can be found in Figure 1 of our previous paper Balsiger:2018ezi .

In Balsiger:2018ezi we implemented the RG evolution factor $\bm{U}(\{\underline{n}\},\mu_{s},\mu_{h})$ in the large-$N_{c}$ limit using the parton shower method proposed by Dasgupta and Salam in Dasgupta:2001sh . We don’t want to repeat the entire discussion here, but we give the algorithm in Appendix B, since we need to extend it to compute the soft functions, as discussed below. Let us also list the one-loop anomalous dimension, since its form will be relevant in the discussion of the jet mass below. It is given by Becher:2016mmh

[TABLE]

The entries $\bm{R}_{m}$ and $\bm{V}_{m}$ are angular functions associated with the emission of a real or virtual soft gluon and take the form

[TABLE]

where the color matrices $\bm{T}_{i,L}$ act on the hard function from the left, i.e. on the amplitude, while $\bm{T}_{i,R}$ acts on the conjugate amplitude. The sum runs over all unequal pairs $(ij)$ of the $m$ hard partons. The anomalous dimension involves the dipole radiator

[TABLE]

which is given by the product of the associated eikonal factors. In the virtual corrections, one integrates over the direction $n_{k}$ of the emission. We note that individually $\bm{R}_{m}$ and $\bm{V}_{m}$ suffer from collinear divergences, which cancel in the cross section. In the Monte Carlo implementation, one works with a collinear cutoff to regularize the divergences.

As long as we choose the $\mu_{h}$ and $\mu_{s}$ properly, the hard and soft functions will be free of large logarithms and the large logarithmic terms are resummed in the evolution factor. Because they are free of large logarithms, the higher-multiplicity hard functions are suppressed by $\alpha_{s}$ as $\bm{\mathcal{H}}_{l}\sim\alpha_{s}^{l-2}\bm{\mathcal{H}}_{2}$. At LL level, we thus only need to include the hard function $\bm{\mathcal{H}}_{2}$ and the soft function is given as the unit matrix in the color space $\bm{\mathcal{S}}_{m}\sim\bm{1}$. At LL accuracy, the RG-improved result (5) simplifies to

[TABLE]

To extend these results to NLL, one needs two ingredients: the one-loop matching corrections and the corrections to the RG running due to the two-loop anomalous dimensions. The present paper focuses on the first set of corrections, i.e. ${\rm LL^{\prime}}$ accuracy. Specifically, we need one-loop corrections to $\bm{\mathcal{H}}_{2}$, the tree-level result for $\bm{\mathcal{H}}_{3}$ and the one-loop soft functions $\bm{\mathcal{S}}_{m}$. We write their perturbative expansions in the form

[TABLE]

In this notation, the full ${\rm LL^{\prime}}$ resummed cross section takes the form

[TABLE]

We used here that the leading-order soft function $\bm{\mathcal{S}}_{m}^{(0)}$ is the unit matrix $\bm{1}$ in color space. The first line contains the LL result (11), and the remaining three lines show the different NLO corrections, which are depicted in Figure 2.

The hard functions $\bm{\mathcal{H}}_{m}$ include the momentum conservation and phase-space constraints on the hard partons. For two partons, these constraints render the integrals over the parton directions trivial. The momentum and jet direction constraints impose that the vectors $n_{1}$ and $n_{2}$ must point along the thrust axis and in opposite directions so that

[TABLE]

where we have used that also the color structure is trivial for two hard partons. The function $H_{2}(Q^{2},\mu)$ is the standard dijet hard function

[TABLE]

which arises also for global observables such as the event shape thrust. In the large-$N_{c}$ limit, we should replace $C_{F}\to N_{c}/2$.

In Becher:2016mmh we have derived an expression for the hard function $\bm{\mathcal{H}}_{3}^{(1)}$, which corresponds to the QCD process $\gamma^{\ast}\to q(p_{1})\bar{q}(p_{2})g(p_{3})$. By definition $\bm{\mathcal{H}}_{3}^{(1)}$ only depends on angular information of the three partons, since their energies have already been integrated over. For convenience we split the phase space integration into different regions according to the direction of the thrust axis, which for three-parton final states points in the opposite direction of the most energetic parton. Due to momentum conservation, the three partons must be in a plane. Using invariance of the cross section under rotation around the thrust axis, in Region I only the angles $\theta_{2}$ and $\theta_{3}$, between the partons and the thrust axis, are not fixed.

For convenience we parameterize these angles in terms of two variables $u$ and $v$ each going from [math] to $1$ and defined as

[TABLE]

where the variable $v$ is directly related to the larger angle $\theta_{3}$, while $u$ characterises the relative size of the angles. Please note that the variables $u$ and $v$ differ from the quantities of the same name used in Becher:2016mmh , where we defined the variables such that $v=1$ corresponded to the angle of the jet cone, rather than a $90^{\circ}$ angle as in (16). Because the same hard function $\bm{\mathcal{H}}_{3}^{(1)}$ also arises for the jet mass studied below, we prefer to not incorporate the specific phase-space constraint into its parameterization.

The bare hard function $\bm{\mathcal{H}}_{3}^{(1)}$ in terms of the angles $\hat{\theta}_{2}$ and $\hat{\theta}_{3}$ was given in (4.4) of Becher:2016omr . The corresponding representation includes a $\theta$-function constraint imposed to prevent the thrust axis from flipping. For simplicity, we choose the jet opening half-angle $\alpha\leq\frac{\pi}{3}$ so that the axis constraint is automatically fulfilled. The hard function suffers from divergences when $u$ and $v$ go to zero. In dimensional regularization after performing $\overline{\rm MS}$ subtraction, the contribution of Region I to the renormalized hard function $\bm{\mathcal{H}}_{3}^{(1)}$ is given by

[TABLE]

The function $\Theta_{\rm in}(v)$ ensures that all hard emissions are inside the jet. For the interjet energy flow it is given by $\Theta_{\rm in}(v)=\theta(\delta-v)$, with $\delta=\tan\frac{\alpha}{2}$, where $\alpha$ is the jet opening half-angle. In the large-$N_{c}$ limit, the color structure of the hard functions becomes trivial and we use non-bold symbols such as $\mathcal{H}_{3,{\rm I}}^{(1)}$ to indicate the scalar quantities which are relevant in this limit. The expression for the auxiliary function $F(u,v)$ is given by

[TABLE]

Similarly, in Region II we have

[TABLE]

with the parametrization $\hat{\theta}_{3}=uv$ and $\hat{\theta}_{2}=v$. The function $G(u,v)$ is defined as

[TABLE]

Region III describes the situation, where the gluon is the most energetic particle and we parameterize $\hat{\theta}_{1}=uv$, $\hat{\theta}_{2}=v$. The hard function reads

[TABLE]

with

[TABLE]

Next, we will discuss how to implement the above expressions into the parton shower code. We first rewrite the angular integral in the $\bm{\mathcal{H}}_{3}^{(1)}$ contribution as

[TABLE]

where we have defined $\widehat{\bm{\mathcal{S}}}_{3}(u,v,\mu_{h})=\sum_{m=3}^{\infty}\bm{U}_{3m}(\{\underline{n}\},\mu_{s},\mu_{h})\,\hat{\otimes}\,\bm{1}$, which is the LL RG evolution or parton shower soft function. To implement this formula into a Monte Carlo framework, we will randomly generate $u$ and $v$ and then run the shower $\widehat{\bm{\mathcal{S}}}_{3}(u,v,\mu_{h})$ for the given configuration. There is, however, one complication, namely that the hard function is a distribution and can therefore not be integrated point by point. One way to solve this problem is to evaluate $\widehat{\bm{\mathcal{S}}}_{3}(u,v,\mu_{h})$ on a grid, interpolate and then perform the integrations over $u$ and $v$. This works well because $\widehat{\bm{\mathcal{S}}}_{3}(u,v,\mu_{h})$ is a smooth function of the angles as can be seen from Figure 4. Note in particular that the limit $v\to 0$, in which both angles go to zero and the two Wilson lines become collinear, is completely smooth. In this limit the quark and gluon Wilson lines combine and produce the same radiation as a single quark Wilson line, encoded in the function $\widehat{\bm{\mathcal{S}}}_{2}$. The relation

[TABLE]

will lead to important simplifications below. In the right plot, we show the evolution time dependence of the soft function $\widehat{\bm{\mathcal{S}}}_{3}$ for fixed angles. One observes that the function falls off much faster when the hard partons approach the jet cone. In this configuration, more soft radiation exits the cone, explaining this suppression.

Interpolating the soft function $\widehat{\bm{\mathcal{S}}}_{3}$ gives accurate results, but is not efficient since the function depends on the phase-space constraints and thus needs to be recomputed when one changes the cone angle. It is much more natural to compute the convolution (23) directly in the Monte Carlo code. The simplest way to implement the plus distributions in the hard function into the Monte Carlo is to use a slicing method. To explain it in a simple setting, let us for the moment only consider the $v$ dependence and forget about the variable $u$. Then the convolution (23) takes the form

[TABLE]

where $B(v)$ represents a regular function. Thanks to relation (24) the $A$ term can be combined with the LL parton shower result involving $\widehat{\bm{\mathcal{S}}}_{2}$ and the contribution from $B(v)$ can be computed by randomly generating $v$-values and running the shower for each chosen configuration. The slicing method introduces a lower cutoff $v_{0}$ into the plus distribution integrals $C_{i}(v)$ to ensure that $v$ can not go to zero. With the cutoff in place, we can integrate the subtraction term, e.g.

[TABLE]

where one can use the same Monte Carlo method as for the $B(v)$ terms to simulate the first term with the collinear cutoff $v_{0}$, and then adds back the second term which is given by the LL parton shower result, multiplied by a logarithm of the cutoff parameter. The $v_{0}$ dependence will cancel out between the two terms up to power corrections. The power corrections in the artificial parameter $v_{0}$ can be neglected as long as one chooses it small enough. The slicing method involves large cancellations between the two terms on the right-hand side of (26), so for numerical stability reasons one should not choose $v_{0}$ too small. These two opposing requirements make slicing methods delicate, but we compared to the result using the interpolated soft function $\widehat{\bm{\mathcal{S}}}_{3}$ and found good consistence. The cutoff independence is demonstrated in Figure 12 in Appendix A.

Up to now we have disregarded the $u$ dependence, but the Monte Carlo implementation of the full equations (2), (2) and (21) involves nothing beyond the above discussion, except that we have to consider both integrations. As (24) shows, the soft function becomes trivial for $v\to 0$ and we can combine all $\delta(v)$ dependent terms with the parton shower for $\widehat{\bm{\mathcal{S}}}_{2}$. We thus only need to apply the slicing method to the $\delta(u)\left(\ln^{i}v/v\right)_{+}$ and $(1/u)_{+}(1/v)_{+}$ terms. The corresponding cutoff dependent compensation terms are collected in Appendix A.

The final ingredient we need to implement is the one-loop soft function, which is defined as a sum over all dipoles

[TABLE]

where the sum runs over all unordered pairs $(ij)$. In the large-$N_{c}$ limit only neighbouring legs give a contribution

[TABLE]

We evaluate the one-loop soft function numerically within our Monte Carlo code. It is well suited for this task since it generates emissions between neighbouring dipoles in an efficient way, by randomly choosing the rapidity $\hat{y}$ and azimuthal angle $\hat{\phi}$ of the emission in the COM (center-of-mass) frame of the emitting dipole $(n_{i},n_{j})$. Here and in the following, we will use hats to indicate kinematic quantities in the COM frame. Our hard function shower keeps emitting additional hard partons until one of them enters the veto region at which point it terminates. In our implementation, we use this last parton in the veto region to obtain the NLO correction to the soft function. At NLO, the renormalized soft function can be expressed as

[TABLE]

with $\Theta_{\rm out}^{\rm lab}(\hat{y},\hat{\phi})$ constraining soft radiation to be outside of the jet cone in the lab frame. In the Monte Carlo implementation, the factor in square brackets is a weight factor for the corresponding emission. The auxiliary function $f_{ij}(\hat{\phi},\hat{y})$ connects the transverse momentum $\hat{k}_{T}$ in the COM frame to the energy $Q_{0}$ in the lab frame, $\hat{k}_{T}f_{ij}(\hat{\phi},\hat{y})\leq Q_{0}$, and is given by $f_{ij}(\hat{\phi},\hat{y})=\frac{2}{M}\left(-\beta\cos\hat{\phi}+\cosh\hat{y}\right)$, where $M^{2}=2\,n_{i}\cdot n_{j}$ is the invariant “mass” of the dipole pair, and $\beta=\sqrt{1-M^{2}/4}$. The logarithm of $|\sin\hat{\phi}|$ arises from expanding the azimuthal angular integration in $\epsilon$, which is related to the space-time dimension through $d=4-2\epsilon$. A detailed derivation of expression (29) can be found in Appendix A.

While our slicing implementation of the hard function is simple but specific to the dijet processes and certainly not optimal, the above procedure to obtain the NLO soft function is simple, efficient and general. Compared to the LL parton shower code, including the one-loop soft function correction (29) yields

[TABLE]

where one evolves the hard function from hard scale to soft scale and multiplies it with the soft function $\bm{\mathcal{S}}_{m}^{(1)}$ of the corresponding multiplicity. When running our Monte Carlo code we fill three histograms, one for the LL shower, one for the logarithmic part of (29) and one for the non-logarithmic part. Further details of the Monte Carlo algorithm, including the implementation of the one-loop soft function are given in Appendix B.

The computer time needed to run the shower including the one-loop corrections depends on the maximum evolution time needed in the computation. For the interjet energy flow, we run the shower until $t=0.08$, corresponding to $\mu_{s}\approx 1\,{\rm GeV}$. For a collinear cutoff at $\eta_{\rm cut}=4$ ($\eta_{\rm cut}=5$) in the parton shower we then end up with about $15$ ($30$) hard partons per event on average. To resolve the peak region of the jet mass, discussed in the next section, we have to run to extremely low scales $\mu_{s}=0.275\,{\rm GeV}$, corresponding to $t=0.3$, near the Landau pole at $\Lambda=0.230\,{\rm GeV}$. At this scale, hundreds of partons are generated in each event and we need a few days of computer time on a cluster to obtain our numerical results, which will be presented in Section 4 below.

3 NLL′ resummation for jet mass

Our second task is to perform the resummation for the jet mass distribution at electron-positron colliders. In contrast to the interjet energy flow, this observable suffers from soft-collinear double logarithms. These then constitute the LL results, while the non-global structure only arises at NLL. The resummation of jet mass including the leading non-global logarithms has been discussed in Dasgupta:2001sh ; Banfi:2010pa ; Dasgupta:2012hg ; Becher:2016omr . At NLL level, the non-global logarithms yield a simple overall factor which multiplies the cross section. Beyond NLL this simple factorization does not hold anymore, and one needs to include the corrections piece by piece.222The recent paper Cal:2019hjc on the jet shape includes one-loop corrections only for the global part, which corresponds to $m=1$ in (4), and does therefore not reach full ${\rm NLL}^{\prime}$ accuracy. Including the non-global structure would result in a factorization formula similar to (4.18) in Balsiger:2018ezi . The basic structure of the corrections is of course the same as for the interjet energy flow, see (13) and Figure 2, and we therefore mainly focus on the differences to this case. In addition to the double logarithms, the most important new element is that the factorization arises in Laplace space. We use the same notation as Becher:2016omr , where we presented NLL resummation results. For $\rm NLL^{\prime}$ accuracy we need to keep one-loop matching corrections in the factorization formula (4) and the theorem then reads

[TABLE]

In the first line we must include one-loop corrections for the quark jet function $\tilde{j}_{q}$, the hard function $\bm{\mathcal{H}}_{1}$ and soft functions $\bm{\mathcal{\widetilde{S}}}_{m}$. We do not include the $\mathcal{O}(\alpha_{s}^{2})$ cross terms so that the first line turns into a sum of terms with the individual corrections. The hard function $\bm{\mathcal{H}}_{2}^{i}$ in the second line includes two hard partons in the right jet. Since it involves a power of $\alpha_{s}$ due to the hard emission, the remaining ingredients are only needed at LO. The second line also includes a gluon-jet contribution, for the case where the $q\bar{q}$ pair is in the right hemisphere. The one-loop hard functions are the same as for interjet energy flow, up to the different phase-space constraints. They are given in Appendix C.

In Laplace space, RG-evolution is multiplicative and we can factor out and exponentiate the double logarithms. Removing the double logarithmic part is important since our shower evolution, which also takes place in Laplace space, is purely soft. The subtraction of collinear contributions will also be needed for our numerical computation of the one-loop soft function. Using standard techniques introduced in Becher:2006nr , we can perform the inversion to momentum space analytically at the end and write a momentum space result directly in terms of Laplace-space ingredients.

The anomalous dimension $\bm{\Gamma}^{H}$ in (8) which drives the resummation of the logarithms in interjet energy flow (1) can be viewed in two ways: As the hard anomalous dimension, used to evolve the hard functions to the soft scale, or as the soft anomalous dimension which evolves the soft functions to a higher scale. RG invariance of the cross section implies that the two evolutions must agree. The situation is more interesting for the light-jet mass (4) which involves three ingredients. In this case RG invariance translates into the statement

[TABLE]

where

[TABLE]

The Casimir $C_{i}$ for the quark-jet channel is $C_{q}=C_{F}$, while the gluon configuration has $C_{g}=C_{A}$. In our paper Becher:2016omr , we have analyzed the one-loop soft anomalous dimension and found that it has the form

[TABLE]

where $\bm{\hat{\Gamma}}_{lm}$ is a regular non-logarithmic anomalous dimension, which takes the same form as (8), except for a subtraction to remove the collinear singularities, which give rise to the cusp piece in (34). The subtraction is achieved by replacing the diagonal elements in (8) by $\bm{V}_{m}\to\overline{\bm{V}}_{m}=\bm{V}_{m}-\bm{V}_{0}$, with

[TABLE]

where $\Theta_{\rm L}(n_{k})$ ensures that the emission is in the left hemisphere with the light jet. The trivial color structure arises from color conservation

[TABLE]

Note that $V_{0}$ is equal to the one-loop result (real plus virtual) for the case where there is only one hard parton on the right, which then, by momentum conservation, flies along $n$. The subtraction therefore removes the “global” one-loop part of the soft anomalous dimension. After this, the Monte Carlo result no longer involves collinear singularities. As before we regularize the collinear singularities in the individual entries of $\bm{\hat{\Gamma}}$ using a cutoff. The parton shower algorithm of Dasgupta and Salam Dasgupta:2001sh instead uses a veto algorithm to remove global logarithmic terms. Our subtraction of the global piece has the advantage that our Monte Carlo weights are always positive. Let us also note that the role of the subtraction is to separate out the collinear singularities, so that the same subtraction can be used for any process with the same double logarithmic structure, i.e. also in cases with more complicated geometry, where we cannot analytically compute the one-loop function.

To make use of the separation of the anomalous dimension into two pieces, we now factor the soft function as

[TABLE]

with

[TABLE]

The splitting of the soft function into single and double logarithmic pieces is of course not unique. We have chosen the double-logarithmic “global” part $\widetilde{S}_{G}^{i}$ such that it includes the full one-loop result, so that the “non-global” remainder function $\bm{\hat{S}}_{m}^{i}$ starts at two loops for $m=1$ partons in the right hemisphere. For the gluon case, we only need the tree-level result $\widetilde{S}_{G}^{g}=1$ since the hard function for this channel is suppressed by $\alpha_{s}$.

The global piece fulfills a standard RG-evolution equation driven by the cusp piece of (34) which can be immediately solved in Laplace space. Using the technique introduced in Becher:2006nr , the associated momentum-space solution takes the form

[TABLE]

with $\eta_{S}=2C_{i}A_{\gamma_{\rm cusp}}(\mu_{s},\mu)$, where the logarithm $L_{s}$ has been replaced by a derivative operator with respect to $\eta_{S}$.

With the global function at hand, the Monte Carlo simulation only needs to provide the remainder $\bm{\hat{S}}_{m}^{i}$. Its single logarithmic RG-evolution is obtained by the subtracted parton shower described above and the one-loop correction for an $m$-parton configuration is given by

[TABLE]

which, by construction, is free from collinear logarithms. We compute this difference in the large-$N_{c}$ limit by running the shower until it produces a parton in the left hemisphere, which is the veto region for the present case. The outside parton is the soft emission and we then compute the relevant one-loop weight factor precisely as in (29). The form of the Laplace space soft function can be found in the appendix in (66). When the emission arises from the first dipole, which involves the left parton along $n_{0}=\bar{n}$, we subtract the global part. For the quark-jet channel the subtraction is given by

[TABLE]

with a re-weighting factor

[TABLE]

The factor $X$ is simply the ratio of the radiator (10) associated with the original $(\bar{n},n)$ dipole and the one of the dipole $(\bar{n},n_{j})$ which emits the gluon and defines the frame in which $\hat{y}$ and $\hat{\phi}$ are generated. The subtraction removes the collinear divergence in the $(\bar{n},n_{j})$ dipole and yields $\bm{\hat{S}}_{m}^{q}$. The function $g_{ij}$ in (41) relates the momentum component $\bar{n}\cdot k$ in the lab frame to the transverse momentum $\hat{k}_{T}$ in the COM frame of the dipole $(n_{i},n_{j})$, analogously to the function $f_{ij}$ in (29). Its explicit form is given in the appendix in (C).

The final ingredients in (3) are the one-loop jet functions, which are well known. In Laplace space, the one-loop jet function is given by

[TABLE]

which translates to the momentum-space result Becher:2006nr

[TABLE]

with $\eta_{J}=2C_{i}A_{\gamma_{\rm cusp}}(\mu_{j},\mu)$. The relevant expressions for the ingredients are listed in Appendix E. Combining the global soft function with the jet function, we obtain

[TABLE]

where we define $\eta=\eta_{J}+\eta_{S}$. The full result is obtained after combining this with the subtracted shower evolution, the hard functions and the one-loop soft correction (40). To implement this expression in practice, we run the shower, tabulate the results for the individual contributions to (3) and then replace the global function $\widetilde{S}_{G}^{i}(\omega,\mu_{h})$ in (3) by the full result which includes the hard functions, evolution and one-loop corrections.

Up to NNLL, the integrated heavy-jet mass distribution is obtained as

[TABLE]

Using this result and relation (2) one obtains the light-jet mass.

4 Numerical results

In this section we will present numerical results, first for the interjet energy flow, then for the jet mass. For our plots, we work with $Q=M_{Z}$ and $\alpha_{s}(M_{Z})=0.1181$, and use two-loop $\alpha_{s}(\mu)$ running with $n_{f}=5$ quark flavors. To our knowledge, no measurements are available for the interjet energy flow, but we will compare our results for the jet mass to LEP measurements by ALEPH Buskulic:1992hq .

4.1 Interjet energy flow

For our numerical discussion we choose jet cone size parameter as $\alpha=\pi/3$. This is equivalent to $\delta=\tan\frac{\alpha}{2}=1/\sqrt{3}$, or rapidity gap size $\Delta y=-\ln\delta^{2}\approx 1.1$. We want to avoid small cone angles, or equivalently large rapidity gaps, in order not to have to deal with large collinear logarithms. In our plots we show the gap fraction

[TABLE]

which is the fraction of events in which the soft radiation outside the jets has an energy $E_{s}$ below the cutoff $Q_{0}$. By definition, the amount of energy in the gap must be below $Q/2$, otherwise the thrust axis, which defines our jet axis, would flip. The fixed order result is therefore $R(Q_{0}=Q/2)=1$ at any order in perturbation theory. The $\mathcal{O}(\alpha_{s}^{0})$ result with just two back-to-back partons is of course $R(Q_{0})=1$, a nontrivial $Q_{0}$ dependence only arises at $\mathcal{O}(\alpha_{s})$ when the third parton is inside the gap. We will refer to the $\mathcal{O}(\alpha_{s})$ result as LO.

As a first step, let us check the size of the individual corrections and investigate whether the scale dependence is reduced after including them. In Figure 5 we show the hard and soft corrections separately and then plot the scale bands from varying the associated scales by a factor two around their default values $\mu_{h}=Q$ and $\mu_{s}=Q_{0}$. Compared to the LL scale bands shown in red, the scale dependence is reduced in both cases after including the corrections. We observe that the hard corrections are quite significant and positive, while the soft corrections are moderate and negative. The hard corrections have two sources, virtual corrections to $\bm{\mathcal{H}}_{2}$ and real emission contributions encoded in $\bm{\mathcal{H}}_{3}$. The first of these is just a constant factor multiplying the LL result, while the second one comes together with the higher soft function $\bm{\mathcal{S}}_{3}$. Both corrections are positive. At high values of $Q_{0}$ the three parton contribution from $\bm{\mathcal{H}}_{3}$ is about twice as large as the one from the one-loop correction to $\bm{\mathcal{H}}_{2}$ and it becomes more dominant at smaller values.

It is clear that the large hard function corrections at $Q_{0}\lesssim Q/2$ must be compensated by terms which are power suppressed in $Q_{0}/Q$ and are not captured by the resummation based on the factorization formula (1), which arises in the limit $Q_{0}\to 0$. One can obtain these power suppressed terms by matching to the fixed-order result. More precisely, one adds to the resummed result the fixed-order prediction minus its expansion around $Q_{0}$. The subtraction removes the terms which are already included in the resummation. These power suppressed matching terms can be obtained as

[TABLE]

To evaluate this integral, one computes the cross section to find a parton inside the gap and subtracts from it its soft limit. The subtraction eliminates the virtual contributions and leads to a finite integral, which one can evaluate numerically. However, even after the matching to the fixed order result, the resummed result does not yet tend to $R(Q_{0})=1$ for $Q_{0}\to Q/2$ because we resum logarithms of $\mu_{s}/\mu_{h}\to 1/2$ for $\mu_{s}\approx Q_{0}$ and $\mu_{h}=Q$. To switch off the resummation, one can choose the soft scale in such a way that it approaches the hard scale $\mu_{h}$ as $Q_{0}\to Q_{\rm max}=Q/2$.

This can be achieved, for example with a profile function Abbate:2010xh of the form

[TABLE]

where $x_{s}=1$ corresponds to the default choice and the scale bands can be obtained by varying the parameter $x_{s}$ by a factor two. For low values of $Q_{0}$, this reduces to the standard choice $\mu_{s}(Q_{0})=Q_{0}x_{s}$. The power suppressed term in the denominator are chosen to switch off the resummation at the endpoint $Q_{0}=Q_{\rm max}$, similarly to what is usually achieved through a modification of the logarithms in traditional resummation. The simplest choice for (49) is $n=1$ and $c_{1}=-1$, but we observe that the approach to fixed order is relatively slow. To make it faster, we choose $n=4$ and impose that the first three derivatives at the end-point vanish, explicitly $c_{1}=-4$, $c_{2}=6$, $c_{3}=-4$, $c_{4}=1$. We plot the two different profile functions in Figure 6 and will use $n=4$ as the default in our numerical implementation. The choice of the profile function affects the resummation of power-suppressed contributions. If the shape is important, one should of course compute, or even resum, the power corrections to resolve the difference. The first step would be to include the matching up to NNLO, which would in principle be possible since the fixed-order results are available GehrmannDeRidder:2007hr ; Weinzierl:2009ms ; DelDuca:2016csb . In practice it would require some effort since we would need to compute the fixed-order expansion of our results (including the shower).

In Figure 7, we show an improved numerical result which includes the matching correction $\Delta R(Q_{0})$, shown as a black dotted line, and uses the scale choice (49) to switch off the resummation at the end-point. The matching correction is negative and compensates the large hard corrections near the end-point. The ${\rm LL}^{\prime}$ corrections lead to a larger gap fraction $R(Q_{0})$. As mentioned earlier, there is unfortunately no experimental data to which we can compare our results, but we compare to PYTHIA Sjostrand:2007gs . While the two results are similar at very low $Q_{0}$, PYTHIA is higher at intermediate values. We remind the reader, that the intermediate values heavily depend on the profile function used to switch off the resummation.

4.2 Jet mass

Let us now turn to the jet mass $\rho$. For interjet energy flow, we considered the integrated cross section, i.e. all events with energy in the gap below the veto, while we will look at the differential spectrum in the present case, since this is what was measured by the LEP experiments. We will however compute the spectrum by taking the derivative of the integrated cross section, which has the advantage that the spectrum is correctly normalized if the resummed prediction for the integrated cross section matches the fixed-order result at large $\rho$.

As a first step, we again separately plot the different ingredients and their scale dependence in Figure 8. In the first three plots we compare NLL to $\rm{NLL^{\prime}}$ with corrections from the jet, hard and soft functions. The red bands are the NLL result with scale variation, where we vary either the jet, hard or soft scale by a factor of two around the default values $\mu_{h}\sim Q$, $\mu_{j}\sim\sqrt{\rho}\,Q$ and $\mu_{s}\sim\rho\,Q$. The blue curves show contributions at $\rm NLL^{\prime}$ accuracy from one of the three ingredients with its associated scale variation. Obviously, the scale dependence is strongly reduced from NLL to $\rm NLL^{\prime}$ for jet and hard corrections. The soft scale dependence, on the other hand, is only modestly reduced after including one-loop soft function corrections. The scale bands mostly overlap with each other, which indicates that perturbative convergence is reasonably good in all the three cases.

In the last plot of Figure 8 we show the effect of adding the $\mathcal{O}(\alpha_{s})$ power corrections to the $\rm NLL^{\prime}$ results. The LO power corrections for the heavy-jet mass are known analytically and given in Appendix E. They are the same as for thrust, because the three-parton results for jet mass and thrust agree. Since the light-jet mass vanishes at $\mathcal{O}(\alpha_{s})$, we can immediately also obtain the LO power corrections for the jet mass distribution. From the plot, we observe that the difference between $\rm NLL^{\prime}$ and $\rm NLL^{\prime}+LO$ is very small, and that the contributions from power corrections will reduce the resummed result in the large jet mass region. In order to reproduce the full fixed order result, we use $C_{F}=4/3$ instead of the strict large-$N_{c}$ value $C_{F}=3/2$ for the hard, jet and soft one-loop corrections in the resummed results. We also use the exact color factors in the evolution factors of the global part (3).

The end-point of the jet mass distribution is at $\rho_{\rm max}=1/3$ at $\mathcal{O}(\alpha_{s})$, corresponding to a symmetrical configuration of the three partons. We will work with the same profile function (49) to switch off the higher-order terms at the end point. To adapt it to the present case, we set $Q_{0}=\rho\,Q$ and $Q_{\rm max}=Q/3$. For simplicity, we will adopt the canonical value $\mu_{j}=\sqrt{\mu_{s}\,\mu_{h}}$ in the following and only indirectly vary the jet scale through the variations of $\mu_{s}$ and $\mu_{h}$, which we vary independently by a factor of two around their default values.

At very low values of $\rho$, the scale $\mu_{s}(Q_{0})$ hits the Landau pole at $\Lambda=0.23\,{\rm GeV}$. Near the pole the soft corrections become large and negative, resulting in a negative cross section. To avoid this unphysical behaviour, we replace $\mu_{s}(Q_{0})\to\mu_{s}(Q_{0})+\Lambda$ so that the pole occurs at $\rho=0$. We also exponentiate the hard, jet and soft corrections to avoid the negative cross section. In the left plot of Figure 9 we show our result for the jet mass distribution after these modifications. In the right plot, we show the result with $\mu_{s}(Q_{0})=\rho\,Q$ and without exponentiation. We observe that the soft scale dependence changes sign at a point to the right of the peak. In this region the soft scale dependence becomes very small. With the modifications in $\mu_{s}$, we end up with quite small scale bands to the right of the peak, which are likely not an accurate characterization of the true uncertainties. The ${\rm NLL}^{\prime}$ peak in the right-hand plot is quite a bit higher because the cross section becomes negative below $\rho=0.004$ and our distributions are by construction normalized. An important feature of our result is that peak occurs at a very low value $\rho\approx 0.006$, which corresponds to $\mu_{s}\approx 0.5\,{\rm GeV}$ so that the peak region is strongly affected by nonperturbative effects. In Figure 9 we also show the PYTHIA Sjostrand:2007gs results, both on the parton level (dashed lines) and including hadronisation. The hadronisation effects shift the peak to the right by about $\Delta\rho\approx 0.006$, in accordance to what one expects from non-perturbative effects in the soft functions Korchemsky:1998ev ; Korchemsky:1999kt . The parton-level PYTHIA result is quite close to the ${\rm NLL}^{\prime}$ result.

In Figure 10 we compare the $\rm NLL^{\prime}+LO$ jet mass distribution with ALEPH results Buskulic:1992hq , obtained by combining their measurements for the light-jet and the heavy-jet mass using (2) and adding the uncertainties on the individual measurements in quadrature. One immediately sees that the experimental peak shifted to the right from non-perturbative effects and the shift is compatible with the PYTHIA hadronization result. We also observe that the jet mass distribution falls off quite rapidly and to make the region of larger $\rho$ visible, we include also a logarithmic plot in Figure 10. The plot also illustrates what motivated the profile function (49) with $n=4$. The choices ensures that we start switching off the resummation fairly quickly about half-way to the endpoint and go over to the fixed-order result. The plots show that, compared the LO fixed-order result, resummation greatly improved the description of the experimental data. On the other hand there is — if at all — only a relatively narrow region in $\rho$ in which both higher-order power corrections and non-perturbative corrections are small.

For completeness, we show in Figure 11 numerical results for the heavy-jet mass $\rho_{h}$ and the light-jet mass $\rho_{\ell}$. The heavy-jet mass is global and provides a reference variable at the same accuracy, but free from all the complications which arise for the jet mass. From the difference of the heavy-jet mass and the jet mass we obtain the light-jet mass. This is more sensitive to the non-global structure and also only has a nontrivial distribution at $\mathcal{O}(\alpha_{s}^{2})$ so that there is no matching at the accuracy we work. The end-point for the NLO light-jet mass is at $\rho_{\rm max}=1/6$, which is achieved when the four parton momenta form a tetrahedron, and we use this as the endpoint in our profile function (49). From the plot, one observes that also the heavy-jet distribution is affected by nonperturbative effects in the peak region, however, the peak is at a larger $\rho$ value than for the jet mass itself. Not surprisingly, the worst description of the data arises for the light-jet mass distribution. At larger $\rho$ values the description is worse because the fixed-order result starts at $\mathcal{O}(\alpha_{s}^{2})$ so that the matching corrections are beyond the accuracy of our computation. The peak region is not well described because it is in the nonperturbative regime and very narrow.

5 Conclusion and Outlook

In this paper we analyzed non-global observables and, for the first time, went beyond a resummation of only the leading non-global logarithms. Specifically, we analyzed the single-logarithmic interjet energy flow at ${\rm LL}^{\prime}$ and the double-logarithmic jet mass at ${\rm NLL}^{\prime}$. The prime indicates that we included the full next-to-leading-order corrections to the hard and soft functions, as well as the jet function in the case of jet mass. The practical implementation of these corrections is the main result of the present paper. To achieve full NLL resummation for the interjet energy flow, and NNLL accuracy for the jet mass, we will need to also include the two-loop corrections the RG running, but we observe that the inclusion of the one-loop matching corrections already leads to an improved description of these observables. Since the jet mass peaks at a low value corresponding to a soft scale of $M_{J}^{2}/Q\approx 0.5\,{\rm GeV}$ for LEP energies, the peak region is strongly affected by non-perturbative effects, similar to what is observed for other event shapes.

Due to the intricate structure of the soft emissions, factorization theorems for non-global observables and the associated RG evolution are much more complicated than in the global case. Instead of analytical computations, one needs to resort to a numerical Monte Carlo framework to perform the resummation. While the global heavy-jet mass only involves a soft function with two Wilson lines, the shower evolution for jet mass produces additional legs, and for low jet masses we can end up with soft emissions from hundreds of hard partons. However, concerning the NLO soft function, this is a minor complication, since we only connect pairs of legs at this accuracy. Indeed, the inclusion of the NLO corrections to the soft function is a minor modification of the leading-logarithmic shower framework. Using the shower emissions which end up in the veto region, we are able to compute the next-to-leading-order correction to the soft function in a general way, with almost no additional computer time.

The more involved part is the implementation of the NLO hard functions. These are in essence the usual real and virtual fixed-order corrections to the Born-level process, but individually suffer from collinear divergences. Computing them in dimensional regularization and renormalizing, one ends up with distributions in the angles of the hard partons which must be implemented into the Monte Carlo framework. We do this with a simple slicing scheme, which works well for two-jet production in $e^{+}e^{-}$ but is certainly not the most efficient method. The problem of combining a parton shower with fixed-order results arises of course also for general purpose showers and elegant solutions such as MC@NLO Frixione:2002ik and POHWEG Nason:2004rx are available and have by now been fully automated. A complication in our case is that our shower systematically neglects small soft momenta and therefore does not conserve momentum. As a result, its kinematics is different from the one in the hard functions. While more work is needed on the NLO hard functions, let us note that we have achieved full automation for the leading-order hard functions in our previous paper Balsiger:2018ezi by working with Les Houches event files generated by the tree-level generator in MadGraph5$\_$aMC@NLO. The same code also provides NLO shower matching and it would be very interesting to adapt it to our shower.

An important next step is of course the inclusion of second-order corrections into the RG-running to achieve the full resummation of subleading non-global logarithms. The corresponding anomalous dimension matrix involves three types of corrections: Double real emissions, real-virtual terms and fully virtual two-loop corrections. The relevant anomalous dimension matrix has been presented in a related framework by Caron-Huot Caron-Huot:2015bja . We are working on determining the anomalous dimension also in our formalism. The implementation into a Monte Carlo framework will be nontrivial, because one needs to numerically handle the collinear singularities of the individual entries. There are a number of recent papers addressing the issue of double emissions in general parton showers Li:2016yez ; Hoche:2017iem ; Hoche:2017hno ; Dulat:2018vuy .

A second interesting challenge is the inclusion of finite-$N_{c}$ effects, especially for non-global observables at hadron colliders. Our RG-evolution framework is in the general class of showers characterized in Nagy:2007ty and valid at finite $N_{c}$, but implementing the interference effects and complex phases which arise beyond $N_{c}\to\infty$ is challenging. Interesting progress towards the computation of such corrections has been made in Platzer:2012np ; Martinez:2018ffw .

We have analyzed two simple non-global observables in the present paper. This is a first step, but our ultimate goal is of course to use the same methods to understand jet structure at the LHC. For narrow jets, the non-global structure actually factorizes into a structure for each separate jet Dasgupta:2012hg ; Becher:2015hka ; Becher:2016mmh . Boosting our hemisphere jet mass result such that the left hemisphere transforms into a cone of radius $R$, one immediately obtains the non-global structure of the jet mass for an LHC jet of this radius. It will be interesting to analyze such observables in the future.

Acknowledgements.

The research of T.B. is supported by the Swiss National Science Foundation (SNF) under grant 200020_182038. T.B. and D.Y.S. would like to express a special thanks to the Mainz Institute for Theoretical Physics (MITP) for its hospitality and support. D.Y.S. acknowledges helpful discussions with Pier Monni on numerical Monte Carlo resummation methods. The authors thank Rudi Rahn and Thomas Rauh for comments on the manuscript.

Appendix A Hard and soft functions for interjet energy flow

A.1 Hard functions in the slicing scheme

We discussed in Section 2 that one can use a simple slicing method to implement the plus distribution terms inside hard function $\mathcal{H}_{3}^{(1)}(u,v)$ into the shower. In the main text, we have explained the procedure using the toy example (25) in which we disregarded the $u$ dependence. In this appendix we now provide the full expression for the hard function. As explained in Section 2, we can directly integrate over $u$ for the $\delta(v)$ terms and combine them with $\widehat{\mathcal{S}}_{2}$. Since there are no singularities inside Region III, we only give expressions for Regions I and II,

[TABLE]

where the cutoffs on $u$ and $v$ are chosen as $u_{0}=v_{0}=e^{-\eta_{\rm cut}-1}$ in the parton shower code. The cutoff $\eta_{\rm cut}$ is imposed in the parton shower on the rapidity of the emitted hard partons. It can be imposed in the lab frame or in the COM frame of the emitting dipole, see Balsiger:2018ezi for more discussions. We have checked that the cutoff dependence can be neglected, as can be observed in Figure 12 in which we show a numerical comparison between the results based on interpolating the soft function and the slicing method for different cutoffs.

A.2 One-loop soft functions

At the one-loop level, virtual corrections from soft gluons are scaleless (and therefore vanish in dimensional regularization), and we only need to include real-emission contributions. The soft function consists of a $d$-dimensional integral with phase-space cuts which ensure that the real emission is outside the jets (the inside part is again scaleless). The relevant soft integral is given by

[TABLE]

with $\tilde{\mu}=e^{\gamma_{E}}\mu^{2}/(4\pi)$ with $v^{\mu}=(1,0,0,0)$ and $v\cdot k=k^{0}$. To evaluate the contribution of the $(n_{i},n_{j})$ dipole, we Lorentz transform into a frame where the vectors $n_{i}$ and $n_{j}$ are back-to-back and the reference vectors take the form

[TABLE]

where $M^{2}=2\,n_{i}\cdot n_{j}$ is the invariant mass of the dipole pair, and $\beta=\sqrt{1-M^{2}/4}$. In this frame, we parameterize the integration momentum as

[TABLE]

With $d=4-2\epsilon$, the integral then reads

[TABLE]

where $\Omega_{d}$ is the surface of the $d$-dimensional unit sphere and $\Omega_{1}=2$. Introducing the auxiliary function $f_{ij}$ via

[TABLE]

we can perform the integral over $\hat{k}_{T}$. This integration yields a soft divergence, which is renormalized away in the ${\rm\overline{MS}}$ scheme. After expanding in $\epsilon$ we then immediately arrive at expression (29) which only involves a finite angular integration which we perform with the parton shower, which generates its emissions using the variables $\hat{y}$ and $\hat{\phi}$.

Appendix B Monte Carlo algorithm for the interjet energy flow

The inclusion of the NLO soft function is only a minor modification of the algorithm for LL resummation. In fact, the first three steps are identical to what was shown in Appendix B of Balsiger:2018ezi . The only difference arises in the last step, where we also compute the soft function. To record the results of the shower, we fill three histograms: $h_{U}$ contains the LL evolution, $h_{L}$ the coefficient of the logarithm of the soft function (29) and $h_{c}$ its non-logarithmic part.

The shower algorithm for the evolution of the function $\bm{\mathcal{H}}_{2}^{(0)}(\{n_{1},n_{2}\},Q,\mu_{h})$ to lower scales involves the following steps:

Start at evolution time $t=0$ from an initial event $E$ with vectors $\{n_{1},n_{2}\}$ and weight $w=1$. 2. 2.

Generate a random time step $\Delta t$ according to the probability distribution $\mathcal{P}_{E}(t)=V_{E}\exp(-V_{E}\Delta t)$, and insert the event weight $w$ into the histogram $h_{U}$ at time $t+\Delta t$. 3. 3.

Choose a dipole associated with a pair of neighbouring vectors $n_{i}$ and $n_{j}$ in $E$ with probability $V_{ij}/V_{E}$. Generate a new random vector $n_{k}$ and multiply the weight by the factor $R_{ij}^{k}/V_{ij}$, expressed in the random variables chosen to generate the direction of the new vector $n_{k}$. 4. 4.

If $n_{k}$ is outside the veto region, add this new vector to the event which then becomes $E^{\prime}=\{n_{1},\cdots,n_{i},n_{k},n_{j},\cdots,n_{2}\}$, multiply the weight by a factor $V_{E}/V_{E^{\prime}}$ and return to Step 2. Otherwise, add the weight factors

[TABLE]

to $h_{L}$ and $h_{c}$ at time $t$, go to Step 1 and start a new event.

In terms of these histograms, the soft function correction reads

[TABLE]

while the LL evolution factor is

[TABLE]

We discussed the implementation of $\bm{\mathcal{H}}_{3}^{(1)}$ in the main text. The shower algorithm is the same as the one described above, up to the fact that one starts the shower with a three parton configuration and does not need to compute the one-loop soft function.

Appendix C Hard and soft functions for the jet mass

In this appendix we list one-loop ingredients for the jet mass. The ingredients are closely related to the ones relevant for the interjet energy flow, but the notation is somewhat different. For the jet mass, the hard function $\mathcal{H}_{m}^{q(1)}$ denotes the configuration with a quark on the left and $m$ partons in the right hemisphere, while $m$ simply counts the total number of hard partons for the interjet energy flow. In the large $N_{c}$ limit the renormalized one-loop hard function $\mathcal{H}_{1}^{i}$ are thus given by

[TABLE]

The factor of one half is present because the LO total cross section is a sum of two identical contributions with the quark and anti-quark in the left hemisphere, respectively. The $\delta$-function of $\hat{\theta}_{1}=\tan(\theta_{1}/2)$ with $\theta_{i}$ ensures that the right parton flies along $n$-direction, opposite to the left parton along $\bar{n}$.

Since the thrust axis points along the opposite direction of the most energetic parton for a three-jet configuration, also the hard functions $\mathcal{H}_{1}^{i(1)}(\hat{\theta}_{1},Q,\mu)$ are the same as for the interjet energy flow. We use the same variables $u$ and $v$ introduced for the interjet energy flow to parameterize the angular variables in order to resolve the overlapping divergences inside the angular integration:

[TABLE]

where the regions are depicted in Figure 3. For the jet mass case, we no longer impose a cone constraint (i.e. we can set $\delta=1$), but we need to add the constraint

[TABLE]

to ensure that the thrust axis does not flip. On the level of the bare function, this constraint was given in (4.4) of Becher:2016omr , but was trivially fulfilled for our choice of the cone angle. Due to this constraint, the angle of any parton to the thrust axis cannot be larger than $\frac{\pi}{3}$.

Performing the variable transformation and writing the angular convolution as integrals over $u$ and $v$ as in (23) we have

[TABLE]

where the factor $\frac{1}{2}$ has the same source as in (60) and the interjet functions were given in (2), (2) and (21). The anti-quark hard function $\mathcal{H}_{2}^{\bar{q}(1)}$ is equal to the quark function. For the gluon function, there is also a region $\theta_{\bar{q}}>\theta_{q}$ which is parameterized analogously and gives an identical contribution.

As explained in Becher:2016omr , the soft function for the light-jet mass is directly related to the coft function in Sterman-Weinberg dijet cross section defined in Becher:2016mmh . In Laplace space, we have

[TABLE]

The evaluation of this expression proceeds along the same lines as for the interjet energy flow case derived in detail in Appendix A. If both emitting partons are in the right hemisphere, the renormalized one-loop result is given by

[TABLE]

with the measurement function $\Theta_{\rm L}^{\rm lab}(\hat{y},\hat{\phi})$ constraining the soft radiation to the left hemisphere, and a function

[TABLE]

If one of the two partons is on the left, the function has a collinear divergence, which can be subtracted, as detailed in Section 3. The subtraction was given in (41).

Appendix D Monte Carlo algorithm for the jet mass distribution

In this appendix we provide the Monte Carlo algorithm used for jet mass resummation, which is also applicable for other non-global observables with soft-collinear double logarithms. Compared to interjet energy flow, we need to subtract the global anomalous dimension and the one-loop global soft function. As for the interjet energy case, we fill three histograms: $h_{U}$ contains the LL evolution, $h_{L}$ the coefficient of the logarithm of the soft function (66) and $h_{c}$ its non-logarithmic part.

The algorithm for evolving $\mathcal{H}_{1}^{q}$ to lower scales involves the following steps:

Start at evolution time $t=0$ from an initial event $E$ with vectors $\{\bar{n},n_{1}\}$ and weight $w=1$. 2. 2.

Generate a random time step $\Delta t$ according to the probability distribution $\mathcal{P}_{E}(t)=\overline{V}_{\!E}\exp(-\overline{V}_{\!E}\,\Delta t)$, and insert the event weight $w$ into the histogram $h_{U}$ at time $t+\Delta t$. 3. 3.

Choose a dipole associated with a pair of neighbouring vectors $n_{i}$ and $n_{j}$ in $E$ with probability $V_{ij}/V_{E}$. Generate a new random vector $n_{k}$ and multiply the weight by the factor $R_{ij}^{k}/V_{ij}$, expressed in the random variables chosen to generate the direction of the new vector $n_{k}$. 4. 4.

If $n_{k}$ is in the right hemisphere, add this new vector to the event so that $E^{\prime}=\{\bar{n},\cdots,n_{i},n_{k},n_{j},\cdots,n_{1}\}$, multiply the weight by a factor $V_{E}/\overline{V}_{E^{\prime}}$ and return to Step 2. If $n_{k}$ is in the left hemisphere and was emitted from dipole $(\bar{n},n_{j})$, we need to subtract the global one-loop soft function $\widetilde{S}^{q\,(1)}_{G}$ in equation (41). This is achieved with the weight factors

[TABLE]

which are added to the histograms $h_{L}$ and $h_{C}$ at time $t$. After filling the histograms go to Step 1 and start a new event. Otherwise, add the unsubtracted weight factor

[TABLE]

to the respective histograms, go to Step 1 and start a new event.

The quantity $\overline{V}_{\!E}$ denotes the subtracted global anomalous dimension $\overline{V}_{\!E}=V_{E}-V_{0}$, where $V_{0}$ is the large-$N_{c}$ result for the subtraction (35) obtained by replacing the Casmir operator $C_{i}$ in this equation by $N_{c}/2$ for a quark jet, or $N_{c}$ for a gluon jet, respectively.

Appendix E Ingredients for jet mass resummation

For convenience, we collect here the perturbative results for ingredients used in the resummation formula for jet mass distribution. The evolution factors at NLL accuracy are given by

[TABLE]

with $r=\alpha_{s}(\mu)/\alpha_{s}\left(\nu\right)$. The expressions of the anomalous dimensions used in our paper are

[TABLE]

The LO integrated jet mass distribution is written as

[TABLE]

The integrated light-jet mass distribution is trivial at this order

[TABLE]

because the light jet has zero mass for three partons.