Transverse momentum dependent distributions in dijet and heavy hadron pair production at EIC

JHEP03(2022)047 Published for SISSA by Springer Received:November 18, 2021 Accepted:February 21, 2022 Published:March 8, 2022 Transverse momentum dependent distributions in dijet and heavy hadron pair production at EIC Rafael F. del Castillo,aMiguel G. Echevarria,b,c Yiannis Makrisdand Ignazio Scimemia aDpto. de Física Teórica & IPARCOS, Universidad Complutense de Madrid, E-28040 Madrid, Spain bDepartment of Physics, University of the Basque Country UPV/EHU, Apartado 644, 48080 Bilbao, Spain cUniversity of Alcalá, Dep. of Physics and Mathematics, 28805 Alcalá de Henares (Madrid), Spain dINFN — Sezione di Pavia, via Bassi 6, I-27100 Pavia, Italy E-mail: [email protected],[email protected], [email protected],[email protected] Abstract: We discuss the measurement of gluon transverse momentum distribution (TMD) in dijet and heavy hadron pair (HHP) production in semi-inclusive deep inelastic scattering. The factorization of these processes in position space shows the appearance of a specific new soft factor matrix element on top of angular and complex valued anomalous dimensions. We show in detail how these features can be treated consistently and we discuss a scale prescription for the evolution kernel of the dijet soft function. As a result we obtain phenomenological predictions for unpolarized and angular modulated cross-sections for the electron-ion collider (EIC) using current available information on unpolarized TMD. Keywords: Jets, QCD Phenomenology ArXiv ePrint: 2111.03703 Open Access,c The Authors. Article funded by SCOAP3.https://doi.org/10.1007/JHEP03(2022)047 JHEP03(2022)047 Contents 1 Introduction 1 2 Factorization theorem, frame choice and modulations 3 2.1 Notation and kinematics 3 2.2 Factorization theorem for dijet and heavy hadron pair production 4 3 Cross-sections used in phenomenology 7 3.1 Extracting the Born-level cross-sections 7 3.2 Angle integrated and azimuthally modulated cross-section 8 4 Evolution kernels with angular dependent anomalous dimensions 9 4.1 Dijet soft function and angle dependent anomalous dimensions 9 4.2 Treatment of angular dependent anomalous dimensions and resummation 10 5 Evolution kernels and scale choices 16 5.1 ζ-prescription for dijet evolution kernel 17 6 Dijet and heavy hadron pair (HHP) production at EIC 20 6.1 Results 22 6.1.1 Results for dijet production 22 6.1.2 Results for heavy hadron production 23 7 Conclusions 25 A Hard prefactors 26 B Anomalous dimensions 27 1 Introduction The access to non-perturbative gluon distributions from experiments is notoriously challenging. This is also the case of gluon transverse momentum distributions (TMDs). Gluons enter directly in Higgs production in hadronic colliders [1–5] that has a relatively high mass and low production rates, and quarkonium production both at EIC and LHC [3,6–25] that is sensitive also to the heavy quark hadronization effects [20,21]. Recent studies (see for example [26,27]) suggest that the experimental observation of the dijet imbalance is possible at the future EIC. In a recent work [28] we have proposed the dijet and hadron pair production at electron-ion colliders (EIC) to probe gluon TMD. Previous studies on these – 1 – JHEP03(2022)047 dijet LO process: heavy meson pair at LO: (⇤g) <latexit sha1_base64="J8coPN5hkT51PA5EvkXIytYySMA=">AAAB83icbVBNSwMxEM3Wr1q/qh69BItSPZTdKuix4MVjBfsB3bXMptltaJJdkqxQSv+GFw+KePXPePPfmLZ70OqDgcd7M8zMC1POtHHdL6ewsrq2vlHcLG1t7+zulfcP2jrJFKEtkvBEdUPQlDNJW4YZTrupoiBCTjvh6Gbmdx6p0iyR92ac0kBALFnECBgr+VU/BiHg4RzHZ/1yxa25c+C/xMtJBeVo9suf/iAhmaDSEA5a9zw3NcEElGGE02nJzzRNgYwgpj1LJQiqg8n85ik+scoAR4myJQ2eqz8nJiC0HovQdgowQ73szcT/vF5moutgwmSaGSrJYlGUcWwSPAsAD5iixPCxJUAUs7diMgQFxNiYSjYEb/nlv6Rdr3kXtfrdZaVxmsdRREfoGFWRh65QA92iJmohglL0hF7Qq5M5z86b875oLTj5zCH6BefjG4NckJk=</latexit> (⇤f) <latexit sha1_base64="hxhjwwKSLwfAqWdYud89cC8OwLk=">AAAB83icbVBNSwMxEM3Wr1q/qh69BItSPZTdKuix4MVjBfsB3bXMptk2NMkuSVYoS/+GFw+KePXPePPfmLZ70OqDgcd7M8zMCxPOtHHdL6ewsrq2vlHcLG1t7+zulfcP2jpOFaEtEvNYdUPQlDNJW4YZTruJoiBCTjvh+Gbmdx6p0iyW92aS0EDAULKIETBW8qv+EISAh3McnfXLFbfmzoH/Ei8nFZSj2S9/+oOYpIJKQzho3fPcxAQZKMMIp9OSn2qaABnDkPYslSCoDrL5zVN8YpUBjmJlSxo8V39OZCC0nojQdgowI73szcT/vF5qousgYzJJDZVksShKOTYxngWAB0xRYvjEEiCK2VsxGYECYmxMJRuCt/zyX9Ku17yLWv3ustI4zeMooiN0jKrIQ1eogW5RE7UQQQl6Qi/o1UmdZ+fNeV+0Fpx85hD9gvPxDYHXkJg=</latexit> + + . . . + . . . pµ 1 <latexit sha1_base64="tGojqvU+NYcCLWarxz1q7RgcSio=">AAAB8HicdVBNSwMxEJ2tX7V+VT16CRbF07JbRT0WvHisYD+kXUs2zbahSXZJskJZ+iu8eFDEqz/Hm//GtF1Biz4YeLw3w8y8MOFMG8/7dApLyyura8X10sbm1vZOeXevqeNUEdogMY9VO8SaciZpwzDDaTtRFIuQ01Y4upr6rQeqNIvlrRknNBB4IFnECDZWukt6/n3WFemkV654rjcD8tzzBeLnVgVy1Hvlj24/Jqmg0hCOte74XmKCDCvDCKeTUjfVNMFkhAe0Y6nEguogmx08QUdW6aMoVrakQTP150SGhdZjEdpOgc1QL3pT8S+vk5roMsiYTFJDJZkvilKOTIym36M+U5QYPrYEE8XsrYgMscLE2IxKNoTvT9H/pFl1/VO3enNWqR3ncRThAA7hBHy4gBpcQx0aQEDAIzzDi6OcJ+fVeZu3Fpx8Zh9+wXn/AtvnkFw=</latexit> pµ 2 <latexit sha1_base64="ZTWNniBWanf5aeqn0rGMTiqUh0U=">AAAB8HicdVBNSwMxEJ2tX7V+VT16CRbF07JbRT0WvHisYD+kXUs2zbahSXZJskJZ+iu8eFDEqz/Hm//GtF1Biz4YeLw3w8y8MOFMG8/7dApLyyura8X10sbm1vZOeXevqeNUEdogMY9VO8SaciZpwzDDaTtRFIuQ01Y4upr6rQeqNIvlrRknNBB4IFnECDZWukt61fusK9JJr1zxXG8G5LnnC8TPrQrkqPfKH91+TFJBpSEca93xvcQEGVaGEU4npW6qaYLJCA9ox1KJBdVBNjt4go6s0kdRrGxJg2bqz4kMC63HIrSdApuhXvSm4l9eJzXRZZAxmaSGSjJfFKUcmRhNv0d9pigxfGwJJorZWxEZYoWJsRmVbAjfn6L/SbPq+qdu9easUjvO4yjCARzCCfhwATW4hjo0gICAR3iGF0c5T86r8zZvLTj5zD78gvP+Bd1xkF0=</latexit> qµ=Q p2(nµ¯nµ)=(0,0,0,Q) <latexit sha1_base64="7+faBURMonmc9B3AUsQjfCH/RUI=">AAACJHicdVDLSgMxFM34rPVVdekmWJQKWmZqtYIUBDcuLVhb6NSSSTNtaCYzTTJCGeZj3Pgrblz4wIUbv8XMOIKKnnDh5Jx7Se5xAkalMs03Y2p6ZnZuPreQX1xaXlktrK1fST8UmDSxz3zRdpAkjHLSVFQx0g4EQZ7DSMsZnSV+64YISX1+qSYB6XpowKlLMVJa6hVOxte2F8K67QqEo0Yc2XIsVFSJY1jiqbUPbQeJiMfJbRfWYcncS05jt1coWmUzBTTLh7Wjas3SJFO+rCLIcNErPNt9H4ce4QozJGXHMgPVjZBQFDMS5+1QkgDhERqQjqYceUR2o3TJGG5rpQ9dX+jiCqbq94kIeVJOPEd3ekgN5W8vEf/yOqFyj7sR5UGoCMefD7khg8qHSWKwTwXBik00QVhQ/VeIh0jHpXSu+e8h/E+uKmXroFxpVIunO1kcObAJtkAJWKAGTsE5uABNgMEtuAeP4Mm4Mx6MF+P1s3XKyGY2wA8Y7x9xE6Kh</latexit> kµ=⇠ p2xQ¯nµ=⇠ 2x(Q, 0,0,Q) <latexit sha1_base64="xn1t+/hv2JPZ1N8LvPoygrkyK8E=">AAACMHicbVBLSwMxGMzWV62vqkcvwaJU0LJbBb0IBQ967IJ9QLeWbJqtodnsmmSlZdmf5MWfohcFRbz6K0zbPdTWCYFhZj6Sb9yQUalM893ILCwuLa9kV3Nr6xubW/ntnboMIoFJDQcsEE0XScIoJzVFFSPNUBDku4w03P7VyG88EiFpwG/VMCRtH/U49ShGSkud/HX/zvEjeOl4AuHYGdAkduSDUHE5GSQ2hI6LRMyTSQhOpcpwkMCifWzqc2IfdfIFs2SOAeeJlZICSFHt5F+cboAjn3CFGZKyZZmhasdIKIoZSXJOJEmIcB/1SEtTjnwi2/F44QQeaKULvUDoyxUcq9MTMfKlHPquTvpI3ctZbyT+57Ui5V20Y8rDSBGOJw95EYMqgKP2YJcKghUbaoKwoPqvEN8jXYrSHed0CdbsyvOkXi5Zp6WyfVaoHKZ1ZMEe2AdFYIFzUAE3oApqAIMn8Ao+wKfxbLwZX8b3JJox0pld8AfGzy9tjKhe</latexit> Figure 1. Example LO diagrams for the two processes. The momenta qµand kµ(corresponding to the photon and incoming parton momenta respectively) are expressed in the Breit frame. processes were performed in refs. [29,30] for dijet and in refs. [31–33] for heavy-meson pair production in an electron-hadron collider. The relevant processes for our case are `+h→`0+J1+J2+X, and `+h→`0+H+¯ H+X , (1.1) where `and `0are the initial and final state leptons, his the colliding hadron, Jiand H/ ¯ Hare the jets and heavy mesons respectively and Xrepresent undetected particles. At leading order (LO), and ignoring the intrinsic momentum of partons inside the target hadron, the two hard-scattering processes are schematically shown in figure 1. We have shown that these processes are factorizable when considering the cross-section dσ dxdη1dη2dpTdrT ,(1.2) where xis the Bjorken variable, and ηi,rTand pTare respectively the rapidity, the sum of the transverse momenta (with respect to the beam axis) and the average scalar transverse momenta of the two final jets. The factorization condition is |rT|  pTin the Breit frame, as the virtual photon and target-hadron directions are back-to-back. These conditions are certainly fulfilled in the Breit frame for pT∈[5,40] GeV and in the central rapidity region. We also demand that there are no hierarchies among partonic Mandelstam variables, ˆs∼ |ˆ t| ∼ |ˆu|. On the experimental side recent studies [26] suggest that the measurement of dijet imbalance is possible at EIC. For the heavy-meson case monte-carlo generator studies suggest that charmed mesons can be reconstructed [34,35]. The charm production rates have been considered at the LO and NLO QCD for ep →c/¯c+Xin ref. [36]. For our case in principle we require the transverse momenta of the heavy mesons, pH/ ¯ H T, be parametrically larger than their mass, mH, i.e. pH/ ¯ H TmH. At leading order (LO) the jet processes can be initiated by either a gluon or a quark, while in the heavy-meson case only the gluon initial state is relevant. The factorized cross-section results in products/convolutions of several fundamental functions as TMDs, jet/heavy quark distributions, soft function and a new evolution kernel, detailed in [28]. The purpose of the present paper is to provide a phenomenological study of these processes, – 2 – JHEP03(2022)047 including theoretical errors. In order to achieve this, we have used the code Artemide [37, 38], introducing new moduli necessary to describe the present cases. The code already includes quark TMDPDF and the TMD evolution kernel extracted from Drell-Yan and semi-inclusive DIS experiments [39]. The factorization of the cross-section is produced in position space, described by the variable b, conjugate of the momentum rT. In order to Fourier transform the factorized cross-section from bspace to momentum space one has to perform an angular integration on φb(i.e. d2b=b db dφband vJ·b=vJbcos φb) that results non trivial because the anomalous dimensions of several functions depend on this angle and are complex valued. We show that this integration can be performed in resummed perturbation theory and that this addresses the problem of complex values in all anomalous dimensions. As a result the φb-angle integrated cross-section is then factorized into distributions that are b-rotation invariant. The φb-angle integrated dijet evolution kernel is derived integrating a system of coupled differential equations similar to the TMD case. We show and discuss here a specific scale choice prescription that is analogue to the ζ-prescription already discussed in [40], which was already implemented in Artemide. The factorization theorem and the detailed definition of the observables is provided in section 2. The cross-sections subtleties discovered when comparing to other groups are discussed in section 3. The resummation of logs involves angular integrations as explained in section 4, which leads to the evolution kernels detailed in section 5. The phenomenological results obtained with the codes that we have developed are summarized in section 6, after which the conclusions are drawn. 2 Factorization theorem, frame choice and modulations 2.1 Notation and kinematics In order to define angles and to deduce a factorized cross-section we need to establish some kinematics. The direction of the beam is fixed along the ˆzaxis. It is useful to define the four-vectors nµ=1 √2(1,0,0,1),¯nµ=1 √2(1,0,0,−1), so that n2= ¯n2= 0,¯n·n= 1. Then any other four vector can be decomposed into its light-cone components, pµ=p+¯nµ+p−nµ+pµ ⊥= (p+, p−, p⊥)n,(2.1) with p+=n·p, p−= ¯n·p, p2= 2p+p−+p2 ⊥= 2p+p−−p2.(2.2) The direction of the two jets are v1and v2, normalized as v2 J= ¯v2 J= 0, vJ·¯vJ= 1,with J= 1,2,(2.3) and ¯vJare defined by reversing the sign of the spacial components. We have then the standard Lorentz-invariants, Q2=−q2, x =Q2 2P·q,(2.4) – 3 – JHEP03(2022)047 where qµis the momentum of the virtual photon, Pµis the momentum of the target hadron. In the Breit frame we have qµ= (0,0,0, Q)and neglecting mass corrections Pµ=1 2x(Q, 0,0,−Q).The ratio of the longitudinal momenta of the incoming parton and the target hadron is ξ=k+ P+.with kµthe momentum of the parton entering hard process. We can then express the variables Qand ξin terms of the Born level kinematics using the pseudo-rapidities, η1and η2, and the transverse momentum, pT, of the two outgoing partons, Q= 2pTcosh(η−) exp(η+), ξ = 2xcosh(η+) exp(−η+),(2.5) where, neglecting corrections from the target hadron mass, η±=η1±η2 2.The partonic Mandelstam variables can be written using the same variables, ˆs= (q+k)2= 4p2 Tcosh2(η−), ˆ t= (q−p2)2=−4p2 Tcosh(η−) cosh(η+) exp(η1), ˆu= (q−p1)2=−4p2 Tcosh(η−) cosh(η+) exp(η2),(2.6) with pµ 1and pµ 2the momenta of the outgoing partons. At partonic level they satisfy ˆs+ˆ t+ ˆu=−Q2.(2.7) Finally, the transverse momentum imbalance of the two jets, rT, and the hard transverse momentum, pT, are defined through rT=p1T+p2T,pT=p1T−p2T 2,(2.8) where the sub-index 1, 2 refers to the final jets. At Born level p1T=−p2Tand thus rT= 0. It must be taken into account that the hadronization of the outgoing partons will form jet-like configurations along similar directions and wide angle radiation that can escape the jet clustering algorithm, affecting the imbalance. 2.2 Factorization theorem for dijet and heavy hadron pair production The factorization of dijet and heavy hadron pair production at leading power (LP) for semiinclusive deep inelastic experiments has already been provided in [28]. The cross-sections reported here do not take into account any leptonic fiducial cuts, which however could be implemented once the experimental conditions are established (especially at EIC). In this section we recall the main formulas that are used in our phenomenological description. We start with the dijet cross-section which can be written as a sum of terms depending on the parton that initiates the hard process (quark or gluon) dσ2J=dσ(γ∗g) + dσU(γ∗f),(2.9) dσ(γ∗g) = dσU(γ∗g) + dσL(γ∗g).(2.10) – 4 – JHEP03(2022)047 For an unpolarized hadronic process, the cross-section is made out of hard contributions from unpolarized initial quarks dσU(γ∗f), unpolarized initial gluons dσU(γ∗g), and linearly polarized gluons dσL(γ∗f). The quark contribution to the cross-section is dσU(γ∗f) dxdη1dη2dpTdrT =X f σfU 0HU γ∗f→gf (ˆs, ˆ t, ˆu, µ)Zd2b (2π)2exp(ib·rT)ff 1(ξ, b, µ, ζ1)(2.11) ×Sγf (b, ζ2, µ)Cg(b, R, µ)Jg(pT, R, µ)Cf(b, R, µ)Jf(pT, R, µ). In this formula ff 1is the unpolarized quark TMDPDF for flavor f,HUthe hard factor for the unpolarized quark case. The perturbative calculations of TMDPDF has been performed recently at NNLO [5,41–44] and N3LO [45,46]. The jets are described by the product of a collinear-soft function C(f,g)and a jet shape function J(f,g)specific for each partonic flavor. The calculation at NLO of these functions can be found in [47,48] for generic kT-type and cone jet algorithms. The factor Sγf is the dijet soft function for the fundamental representation of SU(3)C and calculated in [28] up to NLO. A corresponding soft factor, Sγg, for the adjoint representation of SU(3)Cis also necessary for the incoming gluon contribution. dσ(γ∗g) dxdη1dη2dpTdrT =X f Hµν γ∗g→f¯ f(ˆs, ˆ t, ˆu, µ)Zd2b (2π)2exp(ib·rT)Fg,µν(ξ, b, µ, ζ1)(2.12) ×Sγg(b, η1, η2, µ, ζ2)Cf(b, R, µ)Jf(pT, R, µ)C¯ f(b, R, µ)J¯ f(pT, R, µ), The hard factor Hµν(µ)accounts for contributions of unpolarized and linearly polarized gluons, Hµν γ∗g→f¯ f=σgU 0HU γ∗g→f¯ f gµν T d−2+σgL 0HL γ∗g→f¯ f−gµν T d−2+vµ 1Tvν 2T+vµ 2Tvν 1T 2v1T·v2T.(2.13) The TMD tensor Fg,µν can be also decomposed in terms of unpolarized and linearly polarized parts, Fµν g(ξ, b) = fg 1(ξ, b)gµν T d−2+h⊥ 1(ξ, b)gµν T d−2+bµbν b2,(2.14) with gµν T=gµν −nµ¯nν−¯nµnν. The hard factors are evaluated up to NNLO in the unpolarized case in [49,50] and at LO for the linearly polarized case [51]. In these equations fg 1and h⊥ 1represent the unpolarized and linearly polarized gluon TMD. Both of them are known perturbatively up to NNLO [5,42,44]. Combining eq. (2.10), (2.12), (2.13), (2.14) – 5 – JHEP03(2022)047 one obtains dσU(γ∗g) dxdη1dη2dpTdrT =σgU 0X f HU γ∗g→f¯ f(ˆs, ˆ t, ˆu, µ)Zd2b (2π)2exp(ib·rT)fg 1(ξ, b, µ, ζ1)(2.15) ×Sγg(b, ζ2, µ)Cf(b, R, µ)Jf(pT, R, µ)C¯ f(b, R, µ)J¯ f(pT, R, µ), dσL(γ∗g) dxdη1dη2dpTdrT =σgL 0X f HL γ∗g→f¯ f(ˆs, ˆ t, ˆu, µ)Zd2b (2π)2exp(ib·rT)h⊥ 1(ξ, b, µ, ζ1)(2.16) ×s2 b−c2 b 2Sγg(b, ζ2, µ)Cf(b, R, µ)Jf(pT, R, µ)C¯ f(b, R, µ)J¯ f(pT, R, µ). We use sb= sin φband cb= cos φbfor the sine and cosine of the angle φbbetween the vectors band v1T, respectively. Each of dσ has a hard factor that describes the initiating interaction. The coefficients σ(f,g),(U,L) 0are introduced such that the leading order hard functions are normalized to the unity, i.e. HU(L) LO = 1 + O(αs). The case of heavy hadron pair is very similar. The measured imbalance rTis rT=pH T+p¯ H T,(2.17) where the superscript Hindicates a generic heavy meson and ¯ Hthe corresponding antiparticle. The imbalance is measured in the Breit frame and assuming the TMD factorization scaling, i.e., |rT|  pH, ¯ H T. We also assume that the two heavy mesons are fragmented near the kinematic end-point and carry most of the energy of the heavy quark coming from the hard process. The cross-section reads dσ(γ∗g) dxdηHdη ¯ HdpTdrT =Hµν γ∗g→Q¯ Q(ˆs, ˆ t, ˆu, µ)Zdb (2π)2exp(ib·rT)Fg,µν(ξ, b, µ, ζ1) ×Sγg(b, µ, ζ2)JQ→H(b, pT, mQ, µ)J¯ Q→¯ H(b, pT, mQ, µ).(2.18) with ηHand η¯ Hthe pseudo-rapidities of the heavy mesons, JQ→Hthe heavy quark jetfunctions [52,53]. The hard, soft, and beam functions are the same as in the dijet case. In the hard function we do not consider corrections due to the quark mass and we define pT=|pH T|+|p¯ H T| 2,(2.19) The heavy quark jet functions, JQ→H, can be partially evaluated in perturbation theory as shown in [28]. We work in the limit pTmHΛQCD and the heavy quark jet function can be re-factorized using bHQET. We also have rTpTso that it is possible to find large logs of two parametrically different scales in the fragmentation process, µ+=mQ,and µJ=mQ rT pT ,(2.20) – 6 – JHEP03(2022)047 that need to be resummed to ensure the convergence of the expansion. Following [28] the jet function can be firstly factorized into a short distance matching coefficient and a bHQET matrix element, JQ→H(b, pT, mQ, µ) = H+(mQ, µ)JQ→Hb,mQ pT , µ,(2.21) where the coefficient H+is H+(mQ, µ) = |C+(mQ, µ)|2.(2.22) and the two-dimensional shape function is defined in momentum space as JQ→H(r) = 1 2p− HNCX Xh0|δ(2)r−iv(¯v·∂)W† vhvβ+|XHihXH|¯ hv,β+Wv/ ¯v|0i.(2.23) Notice that vis a Euclidean, two dimensional, transverse component of the light-like fourvector vµpointing along the direction of the boosted heavy meson. In position space JQ→H is obtained by Fourier transformation JQ→Hb,mQ pT , µ=Zdrexp(ib·r)JQ→H(r).(2.24) The one-loop expression for these quantities are calculated in [28]. 3 Cross-sections used in phenomenology The cross-sections presented in previous section are usually partially integrated in phenomenological observables. We discuss here these integrations, which also allow us to relate the normalization of our cross-section with the ones obtained in the literature. 3.1 Extracting the Born-level cross-sections The tree level cross-sections for the dijet and hadron pair production were considered at tree level in ref. [54]. We start considering the gluon case, from which one can easily deduce also the quark case. The gluon hard contribution to the cross-section is described by dσ(γ∗g) dxdη1dη2dpTdrT =N xshA0+A1cos 2φ0 p+···+B0cos 2φ0 r+···i,(3.1) and the azimuthal angles (φ0 r,φ0 p) of vectors pT,rTare measured with respect to the lepton plane. However, our preferred frame is the one where the φ`angle is measured in the plane defined by pTand qT, the sum of the lepton momenta, and φris the azimuthal angle between rTand pT. In this frame and integrating over the angle φ`we are left with: dσ(γ∗g) dxdη1dη2dpTdrT = 2πpTN xshA0+B2cos(2φr)i,(3.2) – 7 – JHEP03(2022)047 with the factor 2πcoming from φ`integration. The LO expressions are obtained by separating the unpolarized and linearly polarized gluon contributions and Fourier transforming. The unpolarized partonic part has a similar form also for quarks, so that we find dσU(γ∗g) dxdη1dη2dpTdrTLO =σgU 0Zd2b (2π)2exp(ib·rT)fg 1(ξ, b) = σgU 0fg 1(ξ, rT),(3.3) dσU(γ∗f) dxdη1dη2dpTdrTLO =σfU 0Zd2b (2π)2exp(ib·rT)ff 1(ξ, b) = σfU 0ff 1(ξ, rT),(3.4) The same for the linearly polarized gluons gives dσL(γ∗g) dxdη1dη2dpTdrTLO =σgL 0Zd2b (2π)2exp(irT·b)sin2φb−cos2φb 2h⊥ 1(ξ, b). =−σgL 0Zb db dφb 8π2exp irTbcos(φb−φr)cos(2φb)h⊥ 1(ξ, b) = cos(2φr)σgL 0Zb db 4πJ2(rTb)h⊥ 1(ξ, b) =−cos(2φr) 2σgL 0h⊥ 1(ξ, rT),(3.5) where notice that h⊥ 1(ξ, rT)is not the direct Fourier transform of h⊥ 1(ξ, b)and both functions can be related through eq. (2.20) in [3]. We obtain the σ(g,f)(U,L) 0prefactors from the structure functions given in eqs. (3.3, 3.5) in [54] and we list them in appendix A. 3.2 Angle integrated and azimuthally modulated cross-section The scalar cross-section that we finally consider in the phenomenological studies is obtained by integrating over the φrangle dσ dΠdrT =rTZ+π −π dφr dσ dΠdrT ,(3.6) where dΠ = dxdη1dη2dpT. Because the factorized cross-section is always expressed in position space one can write (here J0,2are Bessel functions) dσ dΠdrT =rTZ+π −π dφrZdb (2π)2exp hirTbcos(φb−φr)id˜σ(b) dΠdb =rTZ∞ 0 b db 2πJ0(rTb)Z+π −π dφb d˜σ(b) dΠdb =rTZ∞ 0 b db 2πJ0(rTb)Z+π −π dφbd˜σU(b) dΠdb−cos 2φb 2 d˜σL(b) dΠdb,(3.7) where dσU=dσU(γ∗f) + dσU(γ∗g),dσL=dσL(γ∗g)for the dijet case and dσU,L = dσU,L(γ∗g)for the heavy hadron pair case. In our phenomenological analysis we consider also the azimuthal angle average hcos 2φri ≡ "Z+π −π dφrcos 2φr dσ dΠdrT#, dσ dΠdrT .(3.8) – 8 – JHEP03(2022)047 The part of the cross-section relative to linearly polarized gluons can be treated similarly. In this case we need to incorporate an additional cos 2φbterm in the integrals but this is the only change since the soft and collinear-soft functions that appear in the two contributions are the same as for the dijet case. Using the trigonometric identity cos 2φb= 2 cos2φb−1the integrals for this case can be deduced from the discussion of the unpolarized cross-section by the replacement In(A)−→ −In(A+ 1) + 1 2In(A).(4.38) Equivalently for the case of angular modulation in eq. (3.8) the following transformations have to be performed, d˜σU(b) : In(A)−→ −In(A)+2In(A+ 1) d˜σL(b) : In(A)−→ −In(A)+2In(A+ 1) −2In(A+ 2) .(4.39) In all these cases one obtains a cancellation of the imaginary part of the cross-section. Treating perturbatively the angular integration as discussed in this section leads to write eq. (3.7) for the dijet case as dσ dΠdrT =dσU(γ∗g) dΠdrT +dσU(γ∗f) dΠdrT +dσL(γ∗g) dΠdrT ,(4.40) where dσU(γ∗g) dΠdrT =X f σgU 0HU γ∗g→f¯ f(ˆs, ˆ t, ˆu, µ =pT)Jf(pT, R, µJ)J¯ f(pT, R, µJ) ×Z+∞ 0 bdb J0(brT)fg 1(ξ, b)Rg({µk}, ζ1,0, ζ2,0)→(pT, p2 T,1)ˆσU g(b, R, {µi}), (4.41) dσU(γ∗f) dΠdrT =X f, ¯ f σfU 0HU γ∗f→gf (ˆs, ˆ t, ˆu, µ =pT)Jf(pT, R, µJ)Jg(pT, R, µJ) ×Z+∞ 0 bdb J0(brT)ff 1(ξ, b)Rq({µk}, ζ1,0, ζ2,0)→(pT, p2 T,1)ˆσU f(b, R, {µi}), (4.42) dσL(γ∗g) dΠdrT =X f σgL 0HL γ∗g→f¯ f(ˆs, ˆ t, ˆu, µ =pT)Jf(pT, R, µJ)J¯ f(pT, R, µJ) ×Z+∞ 0 bdb J0(brT)h⊥ 1(ξ, b)Rg({µk}, ζ1,0, ζ2,0)→(pT, p2 T,1)ˆσL g(b, R, {µi}), (4.43) where Rf,g are products of evolution kernels to be described in the next section, and ˆσU,L f,g are the result of φbangular integration and can be written as ˆσU g=IgU const.+as(µC)CU f(b, R, µC) + as(µC)CU ¯ f(b, R, µC) + as(µ0)SU γg(b, ζ2, µ0),(4.44) ˆσU f=IfU const.+as(µC)CU f(b, R, µC) + as(µC)CU g(b, R, µC) + as(µ0)SU γf (b, ζ2, µ0),(4.45) ˆσL g=IgL const.+as(µC)CL f(b, R, µC) + as(µC)CL ¯ f(b, R, µC) + as(µ0)SL γg(b, ζ2, µ0).(4.46) – 15 – JHEP03(2022)047 The functions Cand Sin eq. (4.44)–(4.46) are the result of the φbintegration in collinearsoft and dijet soft functions. For the heavy meson case we have just contributions from gluon scattering, dσ dΠdrT =dσU(γ∗g) dΠdrT +dσL(γ∗g) dΠdrT ,(4.47) and we have to change Jf, ¯ f→H+and Cf→ JQ→Hin eq. (4.41)–(4.43). In the case of angular modulation the cross-sections can also be written as in eq. (4.40)–(4.47), with the correct values of the functions Iconst.,Cand S. The non-perturbative effects are in all cases encoded in the evolution kernels, TMD and jet functions. In the next section we describe how the evolution kernels are defined. 5 Evolution kernels and scale choices The evolution kernels appearing in eq. (4.40)–(4.47) are Rg({µk}, ζ1,0, ζ2,0)→(pT, p2 T,1) =RJf(µJ→pT)2RCf(µC→pT)2 ×Rg F(µ0, ζ1,0)→(pT, p2 T)Rq S(µ0, ζ2,0)→(pT,1),(5.1) Rq({µk}, ζ1,0, ζ2,0)→(pT, p2 T,1) =RJf(µJ→pT)RJg(µJ→pT)RCf(µC→pT)RCg(µC→pT) ×Rq F(µ0, ζ1,0)→(pT, p2 T)Rg S(µ0, ζ2,0)→(pT,1),(5.2) where RJf,g is a jet function kernel, RCf,g is the one of collinear-soft functions, Rq,g Fthe one of TMD and finally Rqg Sis the one of the dijet soft function. In the heavy quark case the evolution kernels are parameterized like in eq. (5.1) with the usual changes Jf→H+ and Cf→ JQ→H. The kernels for single-scale evolution have a standard form and a review up to NLL is given in [58], Ri(µi→pT) = eKi(µi→pT)µi miωi(µi→pT) , i ={Cf,Cg, Jf, Jg,JQ→H, H+}(5.3) where ωi(µi→pT)NLL =−Γ0 i β0ln r+Γ1 Γ0−β1 β0αs(µi) 4π(r−1),(5.4) Ki(µi→pT)NLL =−γ0 i 2β0 ln r−2πΓ0 i (β0)2r−1−rln r αs(pT) +Γ1 Γ0−β1 β01−r+ ln r 4π+β1 8πβ0 ln2r,(5.5) – 16 – JHEP03(2022)047 with r=αs(pT)/αs(µi)and Γ0 Cf=−4CF,Γ0 Cg=−4CA, γ0 Cf/g = 0, mCf/g =Re−γE b, Γ0 Jf= 4CF,Γ0 Jg= 4CA, γ0 Jf= 6CF, γ0 Jg= 2β0, mJf/g =pTR, Γ0 J=−4CF, γ0 J= 4CF, mJ=mQ/pTe−γE b, Γ0 += 4CF, γ0 += 2CF, m+=mQ,(5.6) Initial scales µichoice is given in section 6. The TMD kernel is considered here in the ζ-prescription described in [40] and implemented in the code Artemide [37,38] that we use, Rq,g F({µ0, ζ0}→{µf, ζf}) = ζf ζµ(b, µf)−Dq,g(b,µf) .(5.7) In the next paragraph we define a ζ-prescription also for the dijet evolution kernel Rg S(µ0, ζ2,0)→(pT,1), which is the only missing part. 5.1 ζ-prescription for dijet evolution kernel The angular independent kernel of the dijet soft function is obtained as a solution of a coupled system of differential equations, reported in eq. (4.13), that are formally very similar to the TMD ones [59,60]. The anomalous dimensions are given by ¯γSγg (µ, ζ) = γcusp2CFln µ2 µ2 0−CAln ζ ζγg 2,0+δγγg S,(5.8) ¯γSγf (µ, ζ) = γcusp(CF+CA) ln µ2 µ2 0−CFln ζ ζγf 2,0+δγγf S,(5.9) where µ0=2 beγE, ζγg 2,0=4p2 T ˆs 2CF CA, ζγf 2,0=4p2 T ˆs CF+CA CFˆ t ˆu CF−CA CF,(5.10) and δγSare the non-cusp SF anomalous dimension, which is known up to three-loops for the gluon-channel and up to one-loop for the quark-channel and are reported in appendix. The anomalous dimension and the rapidity anomalous dimension (RAD) in eq. (4.6), (4.7) satisfy also −d dln ζ¯γSγi (µ, ζ) = d dln µDi(µ, b)=Γcusp(µ)(5.11) The evolution for the SF takes the general form Ri S({µi, ζi}→{µf, ζf}) = exp ZP¯γSγi(µ, ζ)dln µ−Di(µ, b)dln ζ(5.12) with i=q, g and {µi, ζi}and {µf, ζf}being the initial and final points of factorization and rapidity scales. The integration path Pis an arbitrary path in the {µ, ζ}-plane. Eq. (5.11) – 17 – JHEP03(2022)047 ensures that the evolution kernel is path only independent when one knows the complete perturbative expansion of the anomalous dimensions. Since this is not the case the path independence is broken. In order to partially restore the path independence we proceed as in [40] defining a ζ-prescription also for the dijet soft function evolution kernel. The ζ-prescription provides a way to choose the initial scale ζiof the evolution kernel as a function of µand bso that the SF does not depend on the initial scale µi. This is done by taking the integration path through a null-evolution line in the {µ, ζ}-plane and then taking a fixed-µevolution. To find the null-evolution line we interpret the pair of differential equations (5.11) as a two-dimensional gradient equation ∇F=EF, where E= (γS(µ, ζ),−DS(µ, b)). The null-evolution line is then an equipotential line of the field E. In particular, there is a special null-evolution line that passes through the saddle-point {µsaddle, ζsaddle}of the evolution field. We find that the saddle point is exactly µsaddle =µ0and ζγi saddle =ζγi 0. If we parameterize the null-evolution line as {µ, ζµ(b)}, the value of ζµis given by γSγi (µ, ζµ(b)) = 2DSγi (µ, b)dln ζµ(b) dln µ2,(5.13) which is solved perturbatively order by order in αs. The perturbative solution takes the form ζγg µ, pert(b) = µ µ0 2CF CAζγg 0evγg(µ,b),(5.14) ζγf µ, pert(b) = µ µ0 CF+CA CFζγf 0evγf (µ,b),(5.15) where vγi(µ, b) = ∞ X n=0 an s(µ)vγi n(Lµ),Lµ= ln(Bµ2e2γE), vγg 0(Lµ)=0,(5.16) vγg 1(Lµ) = 2CF CA"−β0 12L2 µ+ γ2 2CF−d(2,0) Γ0#,(5.17) vγg 2(Lµ) = 2CF CA"−β2 0 24L3 µ−β1 12 +β0Γ1 12Γ0L2 µ+ β0γ2 2CF 2Γ0−4β0d(2,0) 3Γ0!Lµ − γ2 2CFΓ1−d(2,0)Γ1 Γ2 0 + γ3 2CF−d(3,0) Γ0#,(5.18) vγf 0(Lµ)=0,(5.19) – 18 – JHEP03(2022)047 and we are using the following notation Di(µ, b) = Ci ∞ X n=1 an s(µ) n X k=0 Lk µd(n,k), δγS(µ) = ∞ X n=1 an s(µ)γn,(5.20) β(as) = −∞ X n=0 an+2 sβn,Γcusp(µ) = Ciγcusp(µ) = Ci ∞ X n=0 an+1 s(µ)Γn,(5.21) with Ci=CF, CAfor quark and gluon channel respectively. Notice that v0vanishes as it is proportional to the LO non-cusp AD, which is zero for the SF. The non-cusp AD is not known beyond LO for the quark-channel. The RAD is a function of band therefore has important non-perturbative corrections in the large-bregion. These corrections can be implemented as a model. The way to proceed is to solve (5.13) for a generic non-perturbative RAD. The equation is solvable but it is difficult to obtain the cancellation of perturbative logarithms in the small-bregion. Following [40] we use the perturbative solution for the small-bregion and move to the exact (generic RAD) solution for large-b: ζµ(b) = ζpert µ(b)e−b2/B2 NP +ζexact µ(b)1−e−b2/B2 NP ,(5.22) with BNP being the bvalue where non-perturbative (NP) effects become important (∼2.5 GeV−1). We have already discussed the perturbative solution to eq. (5.13). For the exact solution we find ζγg µ, exact(b) = µ2 µ2 0 2CF CAζγg 0e−gγg(as,DS)/DS,(5.23) ζγf µ, exact(b) = µ2 µ2 0 CF+CA CFζγf 0e−gγf (as,DS)/DS,(5.24) where gγi(as,DS) = 1 as Γ0 2β2 0 ∞ X n=0 an sgγi n(DS), gγg 0=2CF CAhe−p−1 + pi,(5.25) gγg 1=2CF CAβ1 β0e−p−1 + p−p2 2−Γ1 Γ0e−p−1 + p,(5.26) gγf 0=CF+CA CFhe−p−1 + pi,(5.27) gγf 1=CF+CA CFβ1 β0e−p−1 + p−p2 2−Γ1 Γ0e−p−1 + p,(5.28) and p= 2β0DS/Γ0. Finally, the evolution kernel that provides the evolution from the null-evolution line and that passes through the saddle-point to the final ζpoint is given by Rq,g S({µ0, ζ0}→{µf, ζf}) = ζf ζµ(b, µf)−Dq,g(b,µf) ,(5.29) – 19 – JHEP03(2022)047 and if we consider the evolution from an arbitrary initial scale we take Ri S({µi, ζi}→{µf, ζf}) = Ri S({µ0, ζ0}→{µf, ζf}) Ri S({µ0, ζ0}→{µi, ζi}).(5.30) with i=q, g. This discussion concludes the analysis of all terms that appear in the factorization theorem and the scale prescription. We are now ready for the implementation in the code Artemide [37,38]. 6 Dijet and heavy hadron pair (HHP) production at EIC In order to test the phenomenology developed in the previous sections we consider the case of the EIC. In [28] we already studied the coverage of the EIC and we concluded that the most favourable case is given for a value of mass energy for dijet production around √s= 140 GeV and central rapidity, η1=η2= 0. Typical values for jet radii and momenta at EIC are respectively R∼0.7and pT∼Q/2∼20 GeV. In order to simplify the discussion we show plots integrated over Bjorken variable x(the longitudinal fraction of momentum ξthat enters in the TMDPDFs is ξ∼2x) in the allowed kinematic intervals. For the case of central rapidity we have x∈(0.0859,0.5). The cross-sections that we plot are Zxmax xmin dx dσ dxdη1dη2dpTdrTη1, η2, pT (6.1) and its value is presented as a function of the small transverse momentum rT. The crosssections and the error bands are obtained by using and preparing specific moduli for the code Artemide [37,38]. In particular we use the TMD and the TMD evolution kernels already coded in Artemide, that come from the fit [39], while the new functions studied in this work are included in this code for the first time. The gluon TMD is not fitted yet, however in Artemide there is a parameterization for it. The code takes into account that the contribution of linearly polarized gluons is highly suppressed because in the small-b regime the matching of the linearly polarized gluon TMD onto the gluon PDF starts at order α1 sand not at order α0 slike other distributions. In ref. [5] the cross section obtained in this way agrees with Pythya 8 and current experimental results for the Higgs transverse momentum spectrum, which are however not very precise. The non-perturbative effects are expected to be important in the high-bregion and they should not alter the small-bbehavior of this distribution. Notice also that the non-perturbative effects play a role to control the behavior of the distribution around the Landau pole at large-b, which means a further suppression effect at large-b(as we also observe in the case of unpolarized distributions). Summing up, given the current perturbative and non-perturbative knowledge of TMDs, at this stage we prefer not to push for a hypothetical non-perturbative enhancement of the contribution of linearly polarized gluon TMD. The factorization that we propose in general needs information of the non-perturbative effects in several functions. For the dijet case we have C(b, R;pT) = RC(b, R;pT, µC)Cpert (b, R;µC)fNP C(b, R),(6.2) Sγi(b;pT,1) = RS({µ0, ζ0}→{pT,1})Spert γi (b;µ0, ζ0)fNP S(b),(6.3) – 20 – JHEP03(2022)047 C J S Bi NP (GeV−1) 2.5 2.5 2.5 C J bmax (GeV−1) 0.5 0.3 Table 1. Values of non-perturbative parameter BNP and bmax prescription chosen for collinear-soft function, heavy meson jet function and dijet soft function. Impact of the varition of BNP is shown in figure 2. where the functions with suffix pert refer to their perturbative part in the MS scheme which is currently known at one loop. Similarly for the HHP case we need J(b, mQ/pT;pT) = RJ(b, mQ/pT;pT, µJ)Jpert (b, mQ/pT;µJ)fNP J(b;mQ).(6.4) The non-perturbative effects are parameterized as fNP i(b) = exp −b2 (Bi NP)2, i =C,J, S. (6.5) The values of Bi NP define thee non-perturbative model and we have tested several combinations as shown in figure 2. Higher values of Bi NP are more sensitive to the perturbative series in the low transverse momentum spectrum, and in general provide higher values of the observables. In unpolarized TMD cases we have usually that typical values of Bi NP are around 1-3 GeV−1so we have found reasonable to fix their values as in table 1. The factorization scales µCfor dijet and µJfor HHP are chosen to minimize perturbative logarithms and to not hit the Landau pole of the strong coupling constant, µC= 2e−γE1 b+1 bmax ,(6.6) µJ=1 2e−γE1 b+1 bmax .(6.7) This scale choice deserves some comments. In the dijet case the scale choice does not include the dependence on the jet radius R. Similarly, the mass of the ratio mQ/pTdoes not enter the collinear-soft function and heavy meson jet. In all cases, this means that there is not a complete cancellation of the logarithms of these functions. The reason is that the φbintegration imposes some constraints on the choice of scales. In fact, the function A({µi})defined in eq. (4.21), which depends on the initial scale choice for the soft function, collinear-soft function and heavy meson jet, needs to be A>−1/2in order to have a well defined angular integration. Because of this constraint some scale choices which could be considered like for instance µC=Re−γE b, µJ=mQ/pTe−γE b,(6.8) can not be used. As a result in our approach we only partially resum the logs in the collinear-soft function and the heavy meson jet in order to maintain the structure of ζprescription and double scale evolution in the soft function that is described in section 5. This leads to the initial scales in eq. (6.6), (6.7). – 21 – JHEP03(2022)047 Finally for the dijet soft function we use the b∗-prescription in the same way as for the TMDPDF: µS=2e−γE b∗, b∗=b p1 + b2/b2 max .(6.9) Concerning the theoretical errors, the scale variations in collinear-soft and heavy meson jet function are the main source. This is due to the non-cancellation of logs in the functions by the choice of the initial scales. The choice of the values bmax for collinear-soft function and heavy meson jet account for the convergence of our perturbative result. A more consistent way to treat the resummation of these scales is left for a future work, involving the refactorization of these functions. For functions that do not depend on bthe initial scale choice does not require a prescription or NP-model and it is dictated by the cancellation of the logarithms. For the jet function and the H+matching coefficient we have µJ=pTR, µ+=mQ.(6.10) We use a NP-model for the rapidity anomalous dimension that enters the exact solution for the null-evolution ζµline as it is explained in section 5. In particular, we use the same model that has been used for TMDPDF in [39] DNP F,S =c0bb∗, c0= 2.5·10−2.(6.11) This model dictates how the rapidity anomalous dimension behaves in the large-bregion and is used for both dijet soft function and TMDPDF when performing double scale evolution. While a color re-scaling of the non-perturbative models for gluon TMDPDF and gluon channel soft function with respect to their quark analogues is possible, we observe that this change does not have a significant impact on the cross-section and, therefore, we choose to keep the same model for both quarks and gluons. 6.1 Results In this section we show our results for the differential cross-section for both dijet and heavy hadron pair production processes. Differential cross-sections are shown with error bands coming from scale variation of the different final and initial scales of the functions appearing in our factorization formulas. Scale variation bands are obtained by changing the considered scale by a factor of 2 up and down relative to its central value. 6.1.1 Results for dijet production In figure 3we show the impact of the change of jet radius, jet transverse momentum (hard scale) and jet pseudorapidity over total dijet cross-section. We show that for the variation of the jet radius we see a change of around 20% on the cross-section from the central value when taking the jet radius to be ±0.2from R= 0.7. For pTthere is a variation of an order of magnitude in the total cross-section when taking ±5GeV from 20 GeV. This corresponds to Q= 30,40,50 GeV respectively. Finally, for pseudorapidity variation we obtain an order of magnitude difference above and below when compared to the central – 22 – JHEP03(2022)047 � � � � � �� � � � � � �� �� �� ������� � � � � � Figure 2. Impact of BNP variation over dijets and heavy meson total cross-section. Legend correspond to (BS NP, BC NP)and (BS NP, BJ NP)for dijet and HHP production respectively � � � � � �� ��� ��� ��� ��� (a) Rvariation � � � � � �� ���� ���� � �� (b) pTvariation � � � � � �� ���� ���� ���� ���� � � �� (c) ηvariation (η1, η2) Figure 3. Impact of the variation of the jet radius (R), hard scale jet transverse momentum (pT) and jet pseudorapidity (ηi) for dijet production. For pseudorapidity variation legend is shown referring to (η1, η2)pair, dashed and dotted lines correspond to negative and positive rapidity respectively. rapidity case. Positive rapidities (ηi= 0.5) correspond to Q≃53 GeV while negative rapidities (ηi=−0.5) correspond to Q≃32 GeV, so both pTand ηplots are consistent. Notice that total dijet cross-section is not symmetrical for both jet rapidities as for quark channel we have both a quark and gluon jet in the final state. Every other plot is obtained taking R= 0.7,pT= 20 GeV and ηi= 0. In figure 4the result for the cross-section including quark and gluon channels is shown. We consider the contribution of linearly polarized gluons in a separate panel to show that their contribution is completely negligible, being a factor 103-104smaller. This leads to the conclusion that the contribution from the linearly polarized gluons can be neglected when considering the unpolarized cross-section. The angular modulation asymmetry is shown in figure 6, being around 5%. 6.1.2 Results for heavy hadron production The analysis for HHP has followed similar steps of the dijet case when possible. The differential cross-section including all channels is plotted in figure 5. A separate analysis of – 23 – JHEP03(2022)047 ������� ��� ��� ��� ��� ��� ��� � � � � � � � ��� ��� ��� ��� ��� ��� ������� ��� ��� ��� ��� ��� ��� ������� ��� ��� ��� ��� ��� ��� -� -� -� � � � � � -� -� -� � � � � � -� -� -� � � � � � -� -� -� � � � � � Figure 4. Cross-sections for dijet production at EIC with error-bands coming from scale dependence in collinear-soft factor (CSF), hard factor (Hard), jet distributions (Jet) and Wilson coefficients (OPE). Rows correspond to contributions from linearly polarized gluons (top) and total cross-section (bottom). √s= 140 GeV, R= 0.7,pT= 20 GeV, η1=η2= 0. � � � � � �� � � � � � �� � � � � � �� � � � � � �� ������� -� -� � � � � � � � � � � � � -� -� � � � � � ������� -� -� � � � � � ������� -� -� � � � � � Figure 5. Cross-sections for HHP production at EIC with error-bands coming from scale dependence in hard factor (Hard), heavy meson jet function (bHQET), heavy meson jet function matching coefficient (Hard+) and Wilson coefficients (OPE). The rows correspond to contributions from total cross-section (top) and linearly polarized gluons (bottom). √s= 140 GeV, pT= 20 GeV, η1=η2= 0. the contribution of linearly polarized gluons show also in this case that they are completely negligible being suppressed by a factor 102-103. The angular modulation asymmetry is shown in figure 7, being around 5%. – 24 – JHEP03(2022)047 [37] artemide web-page, https://teorica.fis.ucm.es/artemide/ [38] artemide repository, https://github.com/vladimirovalexey/artemide-public. [39] I. Scimemi and A. Vladimirov, Non-perturbative structure of semi-inclusive deep-inelastic and Drell-Yan scattering at small transverse momentum,JHEP 06 (2020) 137 [arXiv:1912.06532] [INSPIRE]. [40] I. Scimemi and A. Vladimirov, Systematic analysis of double-scale evolution,JHEP 08 (2018) 003 [arXiv:1803.11089] [INSPIRE]. [41] T. Gehrmann, T. Luebbert and L.L. Yang, Calculation of the transverse parton distribution functions at next-to-next-to-leading order,JHEP 06 (2014) 155 [arXiv:1403.6451] [INSPIRE]. [42] M.G. Echevarria, I. Scimemi and A. Vladimirov, Unpolarized transverse momentum dependent parton distribution and fragmentation functions at next-to-next-to-leading order, JHEP 09 (2016) 004 [arXiv:1604.07869] [INSPIRE]. [43] M.-X. Luo, X. Wang, X. Xu, L.L. Yang, T.-Z. Yang and H.X. Zhu, Transverse parton distribution and fragmentation functions at NNLO: the quark case,JHEP 10 (2019) 083 [arXiv:1908.03831] [INSPIRE]. [44] M.-X. Luo, T.-Z. Yang, H.X. Zhu and Y.J. Zhu, Transverse parton distribution and fragmentation functions at NNLO: the gluon case,JHEP 01 (2020) 040 [arXiv:1909.13820] [INSPIRE]. [45] M.-x. Luo, T.-Z. Yang, H.X. Zhu and Y.J. Zhu, Quark transverse parton distribution at the next-to-next-to-next-to-leading order,Phys. Rev. Lett. 124 (2020) 092001 [arXiv:1912.05778] [INSPIRE]. [46] M.-x. Luo, T.-Z. Yang, H.X. Zhu and Y.J. Zhu, Unpolarized quark and gluon TMD PDFs and FFs at N3LO,JHEP 06 (2021) 115 [arXiv:2012.03256] [INSPIRE]. [47] A. Hornig, Y. Makris and T. Mehen, Jet shapes in dijet events at the LHC in SCET,JHEP 04 (2016) 097 [arXiv:1601.01319] [INSPIRE]. [48] M.G.A. Buffing, Z.-B. Kang, K. Lee and X. Liu, A transverse momentum dependent framework for back-to-back photon+jet production,arXiv:1812.07549 [INSPIRE]. [49] T. Becher and M.D. Schwartz, Direct photon production with effective field theory,JHEP 02 (2010) 040 [arXiv:0911.0681] [INSPIRE]. [50] T. Becher, C. Lorentzen and M.D. Schwartz, Precision direct photon and W-boson spectra at high pTand comparison to LHC data,Phys. Rev. D 86 (2012) 054026 [arXiv:1206.6115] [INSPIRE]. [51] Y.-T. Chien et al., Recoil-free azimuthal angle for precision boson-jet correlation,Phys. Lett. B815 (2021) 136124 [arXiv:2005.12279] [INSPIRE]. [52] R.L. Jaffe and L. Randall, Heavy quark fragmentation into heavy mesons,Nucl. Phys. B 412 (1994) 79 [hep-ph/9306201] [INSPIRE]. [53] M. Fickinger, S. Fleming, C. Kim and E. Mereghetti, Effective field theory approach to heavy quark fragmentation,JHEP 11 (2016) 095 [arXiv:1606.07737] [INSPIRE]. [54] D. Boer, P.J. Mulders, C. Pisano and J. Zhou, Asymmetries in heavy quark pair and dijet production at an EIC,JHEP 08 (2016) 001 [arXiv:1605.07934] [INSPIRE]. – 31 – JHEP03(2022)047 [55] Z.-B. Kang, K. Lee, D.Y. Shao and J. Terry, The Sivers asymmetry in hadronic dijet production,JHEP 02 (2021) 066 [arXiv:2008.05470] [INSPIRE]. [56] P. Sun, C.P. Yuan and F. Yuan, Soft gluon resummations in dijet azimuthal angular correlations in hadronic collisions,Phys. Rev. Lett. 113 (2014) 232001 [arXiv:1405.1105] [INSPIRE]. [57] P. Sun, C.P. Yuan and F. Yuan, Transverse momentum resummation for dijet correlation in hadronic collisions,Phys. Rev. D 92 (2015) 094007 [arXiv:1506.06170] [INSPIRE]. [58] A. Hornig, Y. Makris and T. Mehen, Jet shapes in dijet events at the LHC in SCET,JHEP 04 (2016) 097 [arXiv:1601.01319] [INSPIRE]. [59] J. Collins, Foundations of perturbative QCD, Cambridge University Press, Cambridge U.K. (2013). [60] M.G. Echevarria, A. Idilbi and I. Scimemi, Factorization theorem for Drell-Yan at low qT and transverse momentum distributions on-the-light-cone,JHEP 07 (2012) 002 [arXiv:1111.4996] [INSPIRE]. [61] M.G. Echevarria, I. Scimemi and A. Vladimirov, Universal transverse momentum dependent soft function at NNLO,Phys. Rev. D 93 (2016) 054004 [arXiv:1511.05590] [INSPIRE]. – 32 –