Logarithmic Accuracy of Angular-Ordered Parton Showers

We study the logarithmic accuracy of angular-ordered parton showers by considering the singular limits of multiple emission matrix elements. This allows us to consider different choices for the evolution variable and propose a new choice which has both the correct logarithmic behaviour and improved performance away from the singular regions. In particular the description of $e^+e^-$ event shapes in the non-logarithmic region is significantly improved.


Introduction
Monte Carlo event generators [1][2][3][4], which provide a complete description of the complicated hadronic final state observed in high-energy particle collisions, are essential tools as their results can be directly compared with experimental measurements. These simulations combine a calculation of the hard scattering process, usually at next-to-leading order accuracy, with parton shower (PS) evolution from the scale of the hard process to a low energy scale where non-perturbative hadronization models describe the formation of hadrons from the quarks and gluons of the perturbative calculation. Together with a non-perturbative model of multiple parton scattering and decay of the primary hadrons, these generators simulate the final hadronic state. 1 Most of the progress made in this field over the last decade came from matching the parton shower approximation of QCD radiation with fixed-order matrix elements. This increased the accuracy of the cross-section calculation and improved the description of hard radiation, which is not adequately described by the soft and collinear approximations used in parton shower algorithms. In the last few years however there has been a revival of work [6][7][8][9] to improve the accuracy of the parton shower algorithm in antenna [10][11][12] and dipole [13][14][15] showers, as well as work on amplitude-based evolution to treat subleading colour effects [16,17].
A recent work [18] showed that two popular dipole shower algorithms, used in PYTHIA 8 [19] and Dire [20], have issues even at leading-logarithmic accuracy due to the way the singular emissions are split between different dipole contributions and how recoils are handled. The authors considered an initial qq dipole and the emission of two gluons g 1 and g 2 that are both soft and collinear to either of the hard partons and widely separated in rapidity from each other. Given these requirements, the two emissions must be independent and the double-emission probability is dP (2) where y i is the rapidity of the gluon i and p T i is its transverse momentum, all computed in the original qq dipole frame, where the z axis is aligned with the q direction. The second gluon, g 2 , can be emitted either from theq − g 1 or from the q − g 1 dipole. However, although g 2 may be further from g 1 than g 1 is from q or q, when the event is looked at in the emitting-dipole frame, g 2 may be closer in angle to g 1 , which will thus play the role of the emitter. This results in an incorrect colour factor, since C A /2 is assigned instead of C F . This mistake has no effect at leading colour, since C F → C A /2 in the large number of colours limit, though it does correspond to an error in the subleading colour contribution. Furthermore, if g 1 is identified as the emitting particle in the emitting dipole, it has to balance the transverse momentum of g 2 and where the bold symbol indicates it is a two-momentum. This implies that p T 1 can receive a substantial modification if the transverse momentum of the second gluon is only marginally smaller than that of the first emission, thus violating Eqn. (1.1).
In this paper we will use a similar approach to that of Ref. [18] in order to analyse the behaviour of the improved angular-ordered shower of Ref. [21]. The subleading colour issue does not affect an angular-ordered parton shower, which implements colour coherence by construction, so that in the above example g 2 can only be emitted, with the correct colour factor, in a cone around q or g 1 that is smaller than the angle that separates q and g 1 . However, the effect of the recoil must be carefully taken into account. The angular-ordered parton shower, which uses a "global" recoil (the momenta of all partons in the shower are changed to ensure momentum conservation) and 1 → 2 splittings, is significantly different from the dipole showers, which implement "local" recoil (where only the momenta of colour-connected partons change to ensure momentum conservation), as considered in Ref. [18]. While some of the issues considered in Ref. [18] are irrelevant for parton showers using 1 → 2 kinematics and global recoil, some of the underlying physics issues addressed can occur in the angular-ordered parton shower, although they manifest themselves in different ways.
In the next section we briefly introduce the relevant features of a massless parton shower algorithm, including a definition of logarithmic accuracy which will guide our analysis. In section 3 we present the definitions of the parton momenta and kinematics used in the angular-ordered parton shower. These are then used to construct three different interpretations of the evolution variable and consider the logarithmic accuracy of each. We then discuss the tuning procedure used for the Herwig 7 angular-ordered parton shower to ensure a like-for-like comparison between new and old evolution variables. Finally we present our conclusions. In Appendix A we discuss a technical detail related to the splitting g → qq and in Appendix B we explicitly show that the current default recoil scheme implemented in Herwig 7 only correctly describes the double-logarithmically enhanced terms, thus justifying the proposal of a new recoil prescription.

Definition of Logarithmic Accuracy
Fixed-order calculations quickly become cumbersome when we increase the particle multiplicity to take into account the emission of extra jets. However, the leading contribution from such emissions arises in the soft and collinear regions of the phase space, i.e. when we consider the emission of a gluon with vanishing energy or of a parton whose momentum is parallel to the momentum of the emitter. In this latter limit the cross section for the emission of an extra parton is fully factorised, so that we can easily derive the emission probability where z is the light-cone momentum fraction (see eq. (3.3)), P (z) are the collinear splitting functions and t is a scale that approaches 0 in the collinear limit. We see that if we try to integrate the collinear emission probability in (2.1) over the available phase space, there is a logarithmic divergence for t → 0. When we consider the emission of a soft gluon, i.e. when z → 1, we have another source of logarithmic singularities as the splitting kernels all behave like where C = C A in case of gluon splitting and C = C F if the gluon is emitted from a quark line. This simple approximation allows us to correctly take into account double logarithms associated with soft-collinear gluon emission and single logarithms associated with a collinear branching. When we consider n branchings, we can have at most 2n large logarithms, L, of widely disparate scales of the problem, which arise if all the emissions are simultaneously soft and collinear: this means that the emission probability is proportional to α n s L 2n , and we call such contributions leading logarithms (LL). The use of quasicollinear splitting functions [22] gives the first next-to-leading (i.e. single) collinear logarithms (NLL), i.e. α n s L 2n−1 , and together with the choice of the two loop running coupling evaluated using the Catani-Marchesini-Webber scheme [23] at the transverse momentum of the radiated partons [24], includes all leading (double) and next-toleading (single) logarithmic contributions, except for those due to soft wide angle gluon emissions.
In general defining a strict logarithmic accuracy for a parton shower algorithm is difficult. Formally a parton shower algorithm has only leading logarithmic accuracy, although it is able to capture many next-to-leading contributions. There are some classes of infrared-safe observables where an improved coherent branching formalism leads to full next-to-leading log accuracy (e.g. in semi-inclusive hard processes such as deep inelastic scattering and Drell-Yan at large x [23]). In Ref. [18] it was shown that in some regions of the phase space the double-soft-gluon emission probability is not correctly described by dipole showers. In practice, neglecting subleading colour contributions 2 , the parton shower approximation of eq. (1.1) fails only when the transverse momenta of the two emitted gluons are commensurate and thus the recoil procedure quite significantly changes the transverse momentum of the first emission. Since logarithms of commensurate scales are small, it was also found that, for a wide range of event-shape observables, the leading terms are correct but the nextto-leading logarithmic terms are wrong.
Based on this observation, a necessary (but not sufficient) condition for an algorithm to be next-to-leading log accurate is that the singularity structure of the spectrum in Eqn. (1.1) is reproduced in all the regions of the Lund plane [26], which describes the available phase space in terms of the transverse momenta and rapidities of the emitted gluons relative to a suitably-defined frame/axis. As was first pointed out in Ref. [26], and exploited in Ref. [18] to understand the logarithmic accuracy of parton showers, the leading-logarithmic gluon emission is uniform in the plane defined by the logarithm of this transverse momentum and rapidity. Specific corrections to the uniform distribution can be made in specific phase space regions, to promote this description to next-to-leading logarithmic. In more detail, as the cut-off of a parton shower, or value of an event shape observable, is made logarithmically smaller (O < −L ), the area of the Lund plane increases as the square of this logarithm, ∼ L 2 . If a parton shower algorithm makes an order 1 error over an area of the Lund plane, i.e. a region that grows at rate proportional to L 2 , we say that it is not leading-logarithmically accurate. Conversely, if it does not make such an error, we say that it has the potential to be leading-logarithmically accurate. If a parton shower algorithm makes an order 1 error only along a line in the Lund plane, i.e. a region that grows at rate proportional to L, we say that it is leading-logarithmically accurate but not next-to-leading-logarithmically accurate. Our aim is to construct an algorithm that makes order 1 errors only at isolated points in the Lund plane, i.e. regions that do not grow with L, and therefore give rise only to errors in event shape distributions of either next-to-next-to-leading logarithmic or power-suppressed order. Emission of two gluons of similar transverse momenta corresponds to a line in the Lund plane and therefore careful consideration of this configuration is required to reach next-to-leading logarithmic accuracy. The importance of recoil effects for correct description of this region was first pointed out in Ref. [27].
In the following we will consider three recoil scheme prescriptions, one of which leads to an incorrect kinematic mapping in the soft limit. In Appendix. B we explicitly show how this leads to incorrect NLL contributions in the thrust distribution as an example event shape observable.

Kinematics
We will define all momenta in terms of the Sudakov basis such that the 4-momentum of particle l is q l = α l p + β l n + k ⊥l , where the reference vectors p and n are the momentum of the parent parton with onshell mass m 0 and a lightlike vector that points in the direction of its colour partner. They obey so that the transverse momenta are defined with reference to the direction of p and n and the transverse momentum 4-vector k ⊥l is spacelike. If we consider a particle ij that splits into a pair of particles i and j, the light-cone momentum fractions of particles i and j are defined as The relative transverse momentum of the branching is given by 4) and the magnitude of the spatial component is therefore given by The parton shower evolution terminates when where p 2 T min is an infrared cutoff tuned to data of the order of 1 GeV. For many results we will not need a specific representation of the reference vectors. If we do need a representation we will use the choice made in Ref. [21] for final-state radiation with a final-state colour partner, i.e.
where Q is the invariant mass of the radiating particle and its colour partner, and m 0 , m s are the masses of the radiating particle and its colour partner, respectively. Figure 1. The kinematics of two branchings in the angular-ordered parton shower. The off-shell momenta (q i ), on-shell masses (m i ) and light-cone momentum fractions of the partons are shown together with the evolution variable (q i ), transverse momentum (p T i ) and azimuthal angle (φ i ) of each branching.

Single Emission
For the branching 0 → 12, with no further emission we have: (3.9a) where, q ⊥ is the transverse momentum 4-vector, m 0,1,2 are the on-shell masses of the particles, z is the light-cone momentum defined in Eqn. 3.3, β 1,2 are determined by the on-shell condition q 2 1,2 = m 2 1,2 and β 0 by momentum conservation. The virtuality of the parton initiating the branching is therefore where q 2 ⊥ = −p 2 T .

Second emission
We now consider two emissions, the first with z 1 ,q 1 , φ 1 and the second from the first outgoing parton of the first branching with z 2 ,q 2 , φ 2 , as shown in Fig. 1. We define the off-shell momenta of the four partons after the branchings as: where p 2 = m 2 0 , the β i coefficients are fixed by the on-shell condition and momentum conservation and the space-like transverse momentum (3.12) such that q 2 ⊥i = −p 2 T i = −p 2 T i . The virtualities of the branching partons are: In all the cases we will consider parton 4 will be a gluon, m 4 = 0, so that partons 1 and 3 must have the same mass, i.e. m 1 = m 3 . It will also prove useful to define a unit vector in the direction of the transverse momentum, i.e. (3.14)

Interpretation of the Evolution Variable
In Ref. [21] the extension of the original angular-ordered parton shower [28] to include mass effects and longitudinal boost invariance along the jet axis was presented. In this algorithm the evolution variable is in order to include mass effects, in particular the correct mass in the propagator, retain angular-ordering and have a simple single emission probability where P i→jk (z,q) is the quasi-collinear splitting function [22], z is the light-cone momentum fraction and φ is the azimuthal angle of the transverse momentum generated in the splitting. The strong coupling α S is evaluated at the scale from Eqns. (4.1) and (3.10) we can see that µ coincides with the transverse momentum of the splitting [23,29], which we label p T , if m 1 = m 2 = 0. For a single emission (or the last emission in an extended shower) where the children are on their mass-shell, the kinematics are unambiguously defined by Eqn. 4.1 and the ordering variable can be expressed equivalently in terms of q 2 and p 2 T : However, when the children of a branching go on to branch further so that they are off-shell, it is clear from Eqn. (3.13) that we cannot preserve simultaneously q 2 0 and p 2 T . The choice of the preserved quantity will determine the interpretation ofq 2 . The procedure used by Herwig is to first generate a value ofq 2 , z and φ for a branching and calculate the preserved kinematic variable from them. Then the upper limit ofq 2 is calculated for each of the children and the shower proceeds to the next branching. Only at the end of the whole shower evolution, is the generation of each branching completed by constructing its kinematics from its (now off-shell) children's momenta, using the kinematic variable that had been constructed fromq 2 . Thus any other kinematic variables are shifted slightly, to accommodate the change from on-shell to off-shell kinematics. The interested reader can find further details concerning the kinematic reconstruction in Sec. 6.1 of Ref. [29]. As the virtuality acquired from the new partons does not depend upon the azimuthal angle, as can be seen from Eqn. (3.13), we can already anticipate that the shift in the other kinematic variables is not affected by the value of φ.
We will investigate three different choices for the kinematic variable that is preserved.

p T preserving scheme
The original choice of Ref. [21] was to use Eqn. 4.1 together with the expression of the virtuality in Eqn. 3.10, to define the transverse momentum of the branching 0 → 12, where on-shell masses m 1,2 3 are used for the particles produced in the branching. As observed in Ref. [30] this choice tends to give too much hard radiation in the parton shower, as the virtuality of the parent parton can arbitrarily grow after multiple emissions.

q 2 preserving scheme
Ref. [30] suggested that the virtuality of the branching should be determined using the virtualities the particles produced in the branching develop after subsequent evolution, such that Clearly this is the same as Eqn. 4.5 if there is no further emission, i.e. q 2 1,2 = m 2 1,2 .
This choice, however, has the problem that the subsequent evolution of the partons is not guaranteed to result in a physical, i.e. a p 2 T ≥ 0, solution of Eqn. 4.6. In Ref. [30] it was noted that the vetoing of emissions that give a non-physical solution affected the logarithmic evolution of the total number of particles, i.e. the leadinglogarithmic evolution was not correct. Hence, if there was no physical solution the transverse momentum was set to zero such that the virtuality of the branching particle is We remark that, even if the transverse momentum p T of the previous emission changes, the strong coupling of that splitting remains evaluated at z 1 (1 − z 1 )q 1 , i.e. the original transverse momentum in case of massless splitting. Analogously, each emission can be vetoed only when it is generated, so subsequent emissions will not affect this veto.

Dot-product preserving scheme
Motivated by the original massless angular-ordered parton shower of Ref. [28], where the evolution variable was related to the dot product of the outgoing momenta, we investigate the choiceq where the inclusion of the masses is required to give the correct propagator in the general case. However, it is not needed for gluon emission, m 0 = m 1 and m 2 = 0, and only becomes relevant in g → qq branching. In this case As before this reduces to the same result in the case of no further emission. The major advantage of the original massless algorithm [28] was that the subsequent evolution would always have a physical solution for the transverse momentum. If we consider gluon emission the conditioñ is sufficient, but not necessary, for there to be a solution for the transverse momentum in Eqn. 4.9.
If this inequality is satisfied, the virtuality of the branching parton is Assuming that the branching parton was produced in a previous branching, with light-cone momentum fraction z i and evolution scaleq i , the angular-ordering condition ensures thatq < z iqi . Hence so that if Eqn. 4.10 is satisfied for one branching it will also be satisfied for previous branchings. So provided that we requirẽ where m 1,2 are either the physical, or cut-off masses of the partons, the subsequent evolution will be guaranteed to have a physical solution for the transverse momentum.
There are two issues with this choice. The first is that if we impose Eqn. 4.13 on radiation from a heavy quark with mass m, the transverse momentum of the branching must satisfy which, since p T ∼ (1−z)Eθ corresponds to θ ≥ m/E, i.e. the hard dead-cone [31,32] the new algorithm was designed to avoid [21]. In practice we use a cut-off on the transverse momentum of the emission which is fine for radiation from gluons and light quarks, and also for the charm quark since the cut-off is close to the charm mass. For the 3 rd generation quarks we get a small fraction of events where the kinematics cannot be reconstructed ( 0.2 per mille and 0.5% of q → qg branchings for bottom and top quarks, respectively, hardly varying with centre-of-mass energy). However this region is subleading, i.e. does not give rise to either soft or collinear logarithms, and therefore we adopt the approach of setting the transverse momentum of the emission to zero as above in this case.
The second, although less important, issue is the g → qq branching. The limit in this case is presented in Appendix A. For massive quarks, in particular the bottom quark, this limit is stricter than the cut-off on the transverse momentum we use. We therefore have some g → bb branchings where we are forced to set the transverse momentum to zero. Again this region is subleading ( 0.5% of g → bb branchings, again hardly varying with centre-of-mass energy) and therefore does not affect the logarithmic accuracy. In this case the g → qq only gives logarithms of the quark mass, and the neglected region does not contribute to these logarithms.
A full study of these mass effects is beyond the scope of this work, although very important and we hope to return to it in the future.

Phase-space corrections
The angular ordering of the parton shower, which allows a consistent treatment of colour coherence effects, leads to regions of phase space without any gluon emissions. This is the so-called dead zone.  The choice of the preserved quantity in the presence of multiple emissions can significantly affect the phase-space region that is filled by the shower. Fig. 2 illustrates the Dalitz plot for e + e − → qq. We have clustered the partons using the FastJet [33] implementation of the k T jet algorithm [34] and we have switched off g → qq splittings in order to unambiguously define the q andq jets. We can appreciate how little the q 2 -preserving scheme populates the dead zone, coloured in yellow, in opposition to the p T -preserving scheme. This feature is essential when matching to higher order computations, like matrix element corrections, since they will take care to fill this hard region of the phase space. We notice that the dotproduct-preserving scheme (bottom-right pane) displays an intermediate behaviour between the two older schemes, with the number of points in the dead zone for the dot-product-preserving scheme about half of that in p T -preserving scheme.
In order to enforce the similarities between the dot-product preserving scheme and the q 2 one, that is the current Herwig default, we implemented a rejection veto to avoid generating too large virtualities. Indeed the virtuality of the shower progenitor, i.e. the emitter particle that was present prior to the shower, increases when multiple emissions are generated, only in the q 2 -preserving scheme is it kept fixed. To this end, let us consider the two-body phase space for the process e + e − → qq, which reads where Ω is the solid angle that describes the direction of the quark and λ is the Källén function introduced in Eqn. (3.8). When n emissions are generated the phase space becomes where k 2 l is the virtuality developed by the shower progenitor l = q,q. Thus, if we want to factorize the phase space over the original two-body one, we need to include the Jacobian factor . (4.17) Since J < 1, we can simply implement a reweighting procedure: at the end of the showering phase we generate a random number r smaller than 1 and we accept the event only if r < J, otherwise we shower the event anew. Looking at the Dalitz plots (bottom panel of Fig. 2), we see that while this has only a modest effect, it does somewhat suppress, about a 10% reduction, the events in the dead zone. Note that these plots are all made with the same set of parameters.

Assessing the Logarithmic Accuracy
The angular-ordered parton shower has the correct single-emission probability by construction. However it is still instructive to calculate the Lund variables, i.e. the transverse momentum k ⊥ and rapidity y, to see how the Herwig variables relate to the physical ones. For a single gluon emission, m 0 = m 1 = m and m 2 = 0, all three choices for the interpretation of the evolution variable are identical, giving where λ = λ (1, b, c). The first approximation is that both the radiating particle and the spectator are massless, i.e. m → 0, and the second approximation is that the emitted gluon is soft, i.e. z = 1 − with → 0. The Herwig soft collinear gluon emission probability from a massless quark line is given by if we rearrange the above expression in terms of the Lund variables k T and y we reproduce the correct form of the soft collinear emission probability We now need to investigate the accuracy for two successive gluon emissions, i.e. m 0,1,3 = m, m 2,4 = 0. In particular, in angular-ordered parton showers, one can obtain strongly disordered regions in which a second emission is much harder (in energy, contribution to jet virtuality or transverse momentum) than the first. We therefore have to check that the kinematics of the softer first gluon are not disturbed by the second harder one.
The different schemes only affect the relationship between the transverse momenta and the evolution variable, this means that the kinematics are the same in all three schemes when expressed in terms of the transverse momenta. The Lund variables for the two emissions are therefore: (5.4d) All three choices of evolution variable are identical for one emission, therefore 5) and the virtuality of the branching parton is For the first branching the relationships depend on our choice of reconstruction scheme.

p T preserving scheme
If we use the p T preserving scheme the final virtual mass of the original parton is and where we recall thatn i is a unit vector parallel to p T i , see Eqn. (3.14).
In the massless and soft limits, z 1,2 → 1 such that z 1,2 = 1 − 1,2 and 1,2 1, the Lund variables are (5.10d) In the soft limit As the limit from angular-ordering isq 1 ≥q 2 we see that for there is a disordered region where the contribution of a second harder gluon to the virtuality of the original parton is dominant. In this disordered region, k ⊥2 k ⊥1 so that we can neglect 1q1 relative toq 2 and the kinematics are effectively independent. However, there is a region in which the transverse momentum of the first emission overwhelms that of the second, ifq 2 < 1q1 = k ⊥1 . This is the region in which the emission angle of the second gluon is smaller than the recoil angle of the quark from the first gluon (Fig. 3). It is an issue because we have measured the transverse momentum and rapidity relative to the fixed jet axis, not the local axis of emission 4 . fixed jet axis If we calculate the Lund variables using q 3 as the axis: The second gluon variables are now the same as the single emission case, Eqn. 5.1, thus retaining the correct behaviour in the soft limit. The first gluon variables are correct this time, becauseq 2 is always smaller thanq 1 and the factor of 2 makes it arbitrarily smaller. Thus, this scheme is accurate to leading logarithmic order as it reproduces the correct behaviour of the soft, collinear splitting function.

q 2 preserving scheme
For the q 2 preserving scheme so that the transverse momentum is non-zero if In the limit that both z 1,2 → 1 then which is effectively the requirement that the generated virtualities are ordered, which is clearly violated in the disordered region we are concerned about.
In the ordered region in which a solution is possible, the Lund variables, calculated relative to the q 3 axis are: In the bulk of the region, theq 2 2 terms are negligible. However, along the "line" 2q 2 2 ∼ 1q 2 1 the generated k 2 ⊥1 value is wrong by a factor of order 1. Moreover, for most reasonable event shapes, e.g. thrust, the first gluon is the dominant one. Therefore this is a next-to-leading-logarithmic (NLL) error, i.e. the double logarithmic behaviour is correct, while the single soft logarithm is incorrect. An explicit derivation for the case of the thrust is given in Appendix B.
In the disordered region, p T 1 = 0, therefore the Lund variables are: with k 2 ⊥2 and y 2 given by Eqn. 5.18. While the kinematics of the second gluon are correct, kinematics of the first gluon are completely wrong in this region in the Lund plane. This could, in principle, be a leading-log effect. However, for the example of the thrust distribution, in this region the second gluon is the hardest one and the first gluon gives a sub-leading contribution to the observable. Therefore, again, it is only along the line at the edge of this region that one gets a significant effect and it is a NLL error. We conclude that the q 2 preserving looks undesirable, in reconstructing incorrect kinematics over a finite area of the Lund plane. In practice this leads to a NLL error in the thrust distribution (see Appendix B). Related problems with the q 2 -preserving scheme were also noted in Ref. [36].

Dot-product preserving scheme
In the dot-product preserving scheme the transverse momentum of the second branching is unchanged but for the first it becomes The difference relative to Eqn. 5.14 looks minor, but now we have to compareq 2 1 with 2q 2 2 ,q 2 2 has to be smaller thanq 2 1 and the factor of 2 makes it parametrically smaller. The second term can therefore never be as large as the first.
The virtuality of the first parton is which for soft emissions can be dominated by the second emission for 2 > 1 . In this case the transverse momentum of the second branching is In the massless and soft limits the Lund variables, with respect to the direction of p, are while with respect to the direction of q 3 they become k 2 ⊥1 ≈ 2 1 (q 2 1 + 2q 2 2 ); (5.24a)

Global recoil
We also need to consider the impact of the implementation of the global recoil in Herwig 7. For simplicity we will consider the case of two final-state particles, the generic case can be found in Ref. [29]. We have a particle a with momentum 25) which splits into particles b and c, whose momenta are given by (5.26) where λ is the Källén function defined in Eqn. (3.8) and b = m 2 b /s, c = m 2 c /s. During the shower evolution the particles acquire a virtuality q 2 b = b s and q 2 c = c s and their momenta are modified 29) and However, if we want to have two particles with invariant mass q 2 b and q 2 c , whose threemomentum is parallel to the direction of p b and p c respectively, the two particles must have four-momentum equal to (5.31) As q b + q c = q b + q c , they can be simply related by a Lorentz transform along the p b (p c ) direction. The boost parameter for b is where we have used the shorthand notation λ = λ(1, b, c) and λ = λ (1, b , c ). The expression may look complicated, but if we consider that b, c, b and c are all much smaller than 1, we get Also the partons which have q b (q c ) as shower progenitor need to be boosted along the direction of the progenitor. This boost will leave the transverse momentum, the light-cone momentum z and the ordering variableq (since it is expressed in terms of scalar products and z) invariant, but not the rapidity of the particles. Indeed the rapidities of partons having the b as shower progenitor are slightly shift towards smaller values 34) and the rapidities of those coming from the c cascade are slightly pulled in the opposite direction where we expand the result because the boost parameter is generally much smaller than 1, being of the order of (q 2 − m 2 )/s, where q 2 is the virtuality developed by the colour partner of the shower progenitor and m 2 its mass. Let us now discuss the impact of global recoil for soft emission in the massless limit, i.e. for b = c = 0. Let us assume for simplicity that b is a quark q and c is an anti-quarkq. If we use the default Herwig 7 settings, partons originated from b will all have positive rapidity and the single emission probability in the soft limit is while the probability of a soft-emission originated from c is given by 37) and the sum of the two contributions yields However, after we apply our global recoil, the rapidity of the partons gets shifted, to the left for partons coming from b and to the right for those coming from c, causing a double counting of the central-rapidity region. If we callβ the average boost-parameter that is applied after the global recoil, Eqn. (5.38) will be modified to Nevertheless, given the fact thatβ is of the order q 2 /s and for soft emission typically q 2 s, this is a power-suppressed effect, i.e. non-logarithmic, and therefore does not alter the logarithmic accuracy of the parton shower.

Tuning
The new interpretation of the evolution variable means that the hadronization parameters (which are highly sensitive to the PS algorithm) have to be retuned. In order to do so, we follow the same strategy as in Ref. [30]: simulated events are analysed with Rivet [37], which also enables a comparison with experimental results. The dependence on the hadronization and parton shower parameters [38] is interpolated by the Professor program [39], which also finds the set of values which best fit the experimental measurements. In our case, where observables were measured by multiple experiments, only the most recent set of data is used. We have not included LHC data in the tuning due to the high CPU-time requirement. We consider only the transverse momentum (pTmin) and not the virtuality as a cutoff parameter.
Professor offers the possibility to weight each observable differently: we adopted the same weights as in Ref. [30]. Furthermore, as in [30], to prevent the fit being dominated by a few observables with very small experimental uncertainty, we impose a minimum relative error of 5% in the computation of the chi-squared χ 2 .
The following procedure is adopted to tune Herwig 7.
1. First the strong coupling computed in the CMW scheme [23] α CMW s , the minimum transverse momentum allowed in the showering phase p min T , and the light quark hadronization parameters are tuned to event shapes, charged-particle multiplicity and identified-particle spectra and rates which only involve light quark hadrons. This class of observables is labelled as "general" in Tab. 2. 4. We then vary one parameter at a time to see if our tune corresponds to the minimum of the χ 2 . In case any of the parameters are significantly displaced from the minimum, we retune them all (this time considering all the experimental distributions for light, bottom and charm quarks together).
5. We repeat the previous step except that now if any parameters are too far from the minimum of the χ 2 , their values are adjusted by hand. In particular, this is needed for bottom quark hadronization parameters like ClMaxBottom which Professor is not able to tune: this behaviour was also found in Ref. [30].
The values of the default parameters and the new ones we find with our tuning procedure are shown in Tab. 1. The χ 2 per degree of freedom computed with the observables used for the tune, together with some recent data from the ATLAS experiment [88] which is sensitive to both quark and gluon jet properties, are shown in Tab. 2.
From Tab. 1 we can notice that the four reconstruction choices correspond to four significantly different values of the strong coupling, where smaller values correspond to the schemes that give a poorer description of the non-logarithmically enhanced region of the spectrum. The introduction of the veto procedure in the dot preserving scheme indeed induces a 4% enhancement in α s .

Results
In this section we present the results of our simulations, in order to compare the predictions obtained with the several implementations of the recoil discussed above. We first discuss the LEP results, for which Herwig provides matrix-element corrections (MEC), and then LHC ones for which Herwig does not.

LEP results
The first event-shape distribution we consider is thrust, Fig. 4. We find the wellknown behaviour of the p T -preserving scheme, which overpopulates the non-logarithmically-enhanced region of phase space that is already filled by MEC and corresponds to the tail of the distribution. Although the dot-product scheme performs better than the p T one it still overpopulates the dead zone, however the description of the tail of the spectrum improves if we include the rejection veto described in Sec. 4.3.1. In the right panel of Fig. 4 an expanded view of the small 1 − T region is displayed, where we notice that the new choice of the recoil yields a better agreement with data.
Very similar conclusions can be drawn from the thrust major and minor (Fig. 5) distributions, and from the plots of the C-and D-parameters (Fig. 6). For all the event shape distributions except for D, all the options over-populate the first bin, Preserved p T in [30] q 2 in [30] p T q 2 q i · q j q i ·q j +veto Light-quark hadronization and shower parameters   but the q 2 and dot-product-plus-veto are similar to each other and closest to the data.
Looking at the behaviour of the jet resolution parameter in Fig. 7 we observe that the p T -scheme most closely matches the data in the large − log(y 23 ) (small y 23 )  Table 2. The χ 2 per degree of freedom for different choices of the preserved quantity in the angular-ordered shower, obtained with the distributions we used to tune the light, bottom and charm parameters respectively. The χ 2 corresponding to ATLAS jets, particle multiplicities (mult), event shapes (event), identified-particle spectra (ident), quark jets (jet), gluon jets (gluon) and charged particle distributions (charged) are also shown.
Data p T q 2 q 1 · q 2 q 1 · q 2 +veto tail of the distribution. However, in the small − log(y 23 ) region the q 2 scheme yields a better description of the data. The dot-product scheme with the veto behaves very similar to the q 2 scheme, while the scheme without the veto is similar to the p T scheme in the tail of the distribution and to the q 2 one in the opposite limit, thus retaining the best description of the data over the whole range.
In Fig. 8 we show the multiplicity distribution of charged particles in gluon jets Data p T q 2 q 1 · q 2 q 1 · q 2 +veto  for two different gluon energies. We see that the differences between all of the recoil schemes are much smaller than the experimental error and in general they all give a good agreement with the data.
The schemes all fail to describe the peak region of the b fragmentation function, with the different options making little difference, see Fig. 9. Nevertheless, the dotproduct-plus-veto scheme gives the best overall description of b data, as can be seen from Tab. 2.
While all the data shown for e + e − collisions was used as part of the tuning, this is true for all the tunes and therefore the differences are due to the improvements in the parton shower.  Figure 8. Multiplicity distribution of charged particles in gluons jets for two different gluon energies compared with data from OPAL [78].

LHC results
Data from jets at the LHC seem to prefer the p T scheme as shown in Fig. 10. However, this behaviour is due to the absence of MEC in Herwig for the events we are simulating. This implies that the dead zone remains unpopulated and the migration of events in this region partially solves the lack of hard emission generation. Nevertheless we do expect that matching with higher order computations will lead to the same behaviour that we find in LEP observables, i.e. that the p T scheme yields too much hard radiation, while for the q 2 scheme, for which the kinematics of subsequent soft emissions are not guaranteed to be independent, we expect worse behaviour in the opposite region of the spectrum, and the dot-product-preserving scheme features intermediate properties. This data was not included in the tuning.

Conclusions
The pioneering work in Ref. [18] investigated the logarithmic accuracy of dipole showers by focusing on the pattern of multiple emissions. Driven by this work, we have studied how different choices of the recoil scheme in Herwig can impact the logarithmic accuracy of the distributions. We investigated the original choice of Ref. [21], where the transverse momentum of the emission is preserved during the shower evolution, and the alternative proposal to preserve the virtuality of the splitting, introduced in Ref. [30]. We observed that although the latter prescription retains in general a good description of the experimental data, it breaks the formal logarithmic accuracy of the parton shower, as multiple soft emissions well separated in rapidity are not independent. On the other hand, the older recoil scheme overpopulates the non-logarithmically-enhanced region of the phase space, which should not be filled by the parton shower, but instead by higher order computations.
Due to the undesirable features of these recoil schemes, we proposed an alternative interpretation of the angular-ordering variable that well describes the process of multiple independent soft emission and retains a good agreement with data while also considering the hard tail of the distributions. In order to enforce the correct behaviour in the hard region of the spectrum, we implemented a veto that suppresses large virtualities at the end of the parton shower. This veto applies only to final state radiation and in the future we plan to propose an extension which also includes initial state radiation. In the present work we mainly focused on the case of a massless emitter. The study of mass effects is crucial in assessing the accuracy of the parton shower and will be considered in future works.
A g → qq branching in the dot-product preserving scheme In the case of g → qq branching the transverse momentum of the splitting, Eqn. 4.9, becomes where m is the quark mass. So requiring is sufficient, but not necessary, for there to be a physical solution in this case. In this case the virtuality of the branching is which again will allow a solution but give a stricter limit.
B Impact of the recoil scheme on the logarithmic accuracy of the thrust distribution In this appendix we prove that the thrust is described only to LL accuracy in the q 2 -preserving scheme, as this recoil scheme prescription introduces incorrect NLL terms at order α 2 s . To do so, we make use of the same methodology employed in Sec. 4 of Ref. [18], which relies on the CAESAR formalism [35]. We introduce Σ(L), which is the probability an event shape has a value smaller than exp(−L). We have already seen in Sec. 2 that when we perform an expansion in the strong coupling, α s , at most 2 powers of L appear for each power of α s , i.e. where the LL terms are contained in g 1 (α s L), while the NLL terms are in g 2 (α s L).
The Herwig single emission probability can be written as dP Hw7 soft = 2α s C F π dq q d = 2α s C F π dp T p T dη =ᾱ dp T p T dη (B.3) whereᾱ = 2αsC F π , and p T = q and η = log( Q/p T ) are the Lund variables. The impact of the incorrect shower mapping can be written as where we have replaced the 1/2! multiplicity factor with the ordering condition i.e. either η 1 is in the left hemisphere and η 2 is in the right, or they are both in the same hemisphere and ordered with respect to each other. i = log(p T i /Q) and V is the shape observable expressed in terms of the Lund variables, the subscript "correct" means that V is calculated using the correct double-soft kinematics where p T 1 ≡ 1q1 is the transverse momentum of the first emitted gluon, while "PS" denotes the result obtained using the kinematics of the Herwig parton shower (in the double-soft limit). In Sec. 5 we have shown that the double-soft kinematics are correctly mapped if the transverse momenta or the dot products of the momenta of the emitted particles are preserved, so here we only need to consider the case of the q 2 -preserving scheme, which gives inaccurate kinematics when the two gluons are emitted from the same progenitor. We therefore only need to consider positive rapidities, provided we include a factor of 2 In the case of the q 2 -preserving scheme the contribution of the first gluon is modified: we label the new transverse momentum and rapidity as p T 1 and η 1 respectively, while we denote by p T 1 and η 1 the original values. Therefore from Eqn. (5.14) we can read that By observing that the recoil prescription does not change the light-cone momentum fraction of the first gluon, i.e.
we can write where in the first line we have removed the theta function coming from the angularordering condition, Θ(η 2 − η 1 ), and included a factor of 1/2 as the integrand is symmetric in the exchange 1 ↔ 2. In the second line we have defined x i = η i − i − L and reinserted an ordering x 2 > x 1 . Now, the only dependence on L is in the limits on the η integrals, which are trivial, and we can read off the leading power in L, This proves that this choice of the kinematic mapping introduces a NLL discrepancy at order α 2 s (while in the case of dipole showers, the first NLL discrepancy appears at order α 3 s [18]).