Phenomenological studies of initial state effects and jet quenching in High-Energy Nuclear Collisions at LHC

Departamento de F´ısica de Part´ıculas e Instituto Galego de F´ısica de Altas Enerx´ıas Phenomenological studies of initial state effects and jet quenching in High-Energy Nuclear Collisions at LHC Carlota Andr´es Casas Santiago de Compostela, setembro de 2017. UNIVERSIDADE DE SANTIAGO DE COMPOSTELA Departamento de F´ısica de Part´ıculas e Instituto Galego de F´ısica de Altas Enerx´ıas Phenomenological studies of initial state effects and jet quenching in High-Energy Nuclear Collisions at LHC Carlota Andr´es Casas Santiago de Compostela, setembro de 2017. UNIVERSIDADE DE SANTIAGO DE COMPOSTELA Departamento de F´ısica de Part´ıculas e Instituto Galego de F´ısica de Altas Enerx´ıas Phenomenological studies of initial state effects and jet quenching in High-Energy Nuclear Collisions at LHC Tese presentada para optar ao grao de Doutora en F´ısica por: Carlota Andr´es Casas Setembro, 2017 UNIVERSIDADE DE SANTIAGO DE COMPOSTELA Departamento de F´ısica de Part´ıculas e Instituto Galego de F´ısica de Altas Enerx´ıas Carlos Alberto Salgado L´opez, Profesor Titular de F´ısica Te´orica da Universidade de Santiago de Compostela e Carlos Pajares Vales, Catedr´atico Em´erito de F´ısica Te´orica da Universidade de Santiago de Compostela, CERTIFICAN: que a memoria titulada Phenomenological studies of initial state effects and jet quenching in High-Energy Nuclear Collisions at LHC foi realizada, baixo a nosa direcci´on, por Carlota Andr´es Casas, no departamento de F´ısica de Part´ıculas e Instituto Galego de F´ısica de Altas Enerx´ıas desta Universidade e constit´ue o traballo de Tese que presenta para optar ao grao de Doutora en F´ısica. Asinado: Asinado: Carlos Alberto Salgado L´opez. Carlos Pajares Vales. Santiago de Compostela, setembro de 2017. Santiago de Compostela, setembro de 2017. A mi familia ”Si en mi tarjeta pusiera Emilio, en lugar de Emilia, qu´e distinta habr´ıa sido mi vida”. Emilia Pardo Baz´an Abstract Heavy ion collisions (HICs) are the appropriate tools to study Quantum Chromodynamics (QCD) under extreme conditions of energy and density, which are very different from those inside the atomic nucleus. In high-energy nuclear collisions critical temperatures and densities that allow the formation of the so-called quark-gluon-plasma (QGP) are reached. This thesis is focused on the analysis of two types of effects arising in HICs: Initial state effects (IS) and Final state effects (FS). Among the former, there are the nuclear parton distribution functions (nPDFs), whose precise determination is crucial for the correct interpretation of any observable used in HICs. In this thesis a global analysis of nPDFs at next-to-next-to-leading order in perturbative QCD (pQCD) is performed. Another very interesting issue regarding the initial stage of HICs are the collective phenomena that give rise to the QGP. These are addressed using an approach denominated percolation of strings (SPM). The results obtained in this framework for different observables are compared to available experimental data. Amongst the FS it is worth stressing hard probes, which are observables characterized by a high energy or mass. An analysis of single-inclusive suppression of hard particles at different center of mass energies and centralities in presented. The main result of this work is the extraction of the so-called jet quenching coefficient, ˆq. KEYWORDS Heavy ion collisions (HICs), QGP, nPDFs, percolation of strings, jet quenching List of Papers List of the papers that were pusblished, or in process of being published, with the work present in this thesis: 1. Universal geometrical scaling for hadronic interactions, C. Andr´es, A. Moscoso and C. Pajares, Nucl. Phys. A 901 14 (2013). 2. Onset of the ridge structure in AA, pA, and pp collisions, C. Andr´es, A. Moscoso and C. Pajares, Phys. Rev. C 90, 054902 (2014). 3. Universal geometrical scaling of the elliptic flow, C. Andr´es, A. Moscoso, J. Dias de Deus, C. Pajares and Carlos A. Salgado, Phys. Rev. C 90, 034901 (2015). 4. Energy versus centrality dependence of the jet quenching parameter ˆqat RHIC and LHC: a new puzzle?, Carlota Andr´es, N´estor Armesto, Matthew Luzum, Carlos A. Salgado and P´ıa Zurita, Eur. Phys. J. C 76, 475 (2016). 5. Energy loss as the origin of a universal scaling law of the elliptic flow, C. Andr´es, Mikhail Braun, C. Pajares, Eur. Phys. J. A 53 no. 3, 41 (2017). List of proceedings: 1. High-density QCD and the new LHC data, A. Moscoso, C. Andr´es and C. Pajares, Theor. Math. Phys. 176 937 (2013). 2. Universal geometrical scaling in pp, pA and AA collisions and saturation of gluons, C. Andr´es, A. Moscoso and C. Pajares, AIP Conf. Proc. 1606 283 (2014). 3. Universal geometrical scaling of the elliptic flow, C. Andr´es, A. Moscoso, J. Dias de Deus, C. Pajares and Carlos A. Salgado, EPJ Web Conf. 90 08003 (2015). 4. Modelling jet quenching with Quenching Weights, C. Andr´es, N. Armesto, Carlos A. Salgado and Yan Zhu, J. Phys. Conf. Ser. 612 012001 (2015). 5. The onset of the ridge structure in AA, pA and pp collisions, C. Andr´es, A. Moscoso and C. Pajares, Nuclear and Particle Physics Proceedings, 273–275, 1513 (2016). 6. Extracting ˆqfrom single inclusive data at RHIC and at the LHC for different centralities: a new puzzle?, C. Andr´es, N. Amesto, M. Luzum, C. A. Salgado and P. Zurita, arXiv:1612.06781 [nucl-th], accepted in Nuclear and Particle Physics Proceedings by Elsevier. 7. Extracting ˆqin event-by-event hydrodynamics and the centrality/energy puzzle, C. Andr´es, N. Armesto, H. Niemi, R. Paatelainen, C. A. Salgado and P. Zurita, arXiv:1705.01493 [nucl-th], accepted in Nuclear Physics A. I acknowledge the financial support given by the Ministerio Espa˜nol de Educaci´on, Cultura y Deportes under the grant FPU2013-03558. Contents Abstract 3 1 Introduction 5 1.1 Quantum Chromodynamics . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 1.2 Factorization and parton distributions . . . . . . . . . . . . . . . . . . . . . 8 1.3 Hard Probes in Heavy Ion Collisions . . . . . . . . . . . . . . . . . . . . . . 9 1.4 Percolationofstrings............................... 11 2 Nuclear parton distribution functions at NNLO 13 2.1 DGLAP evolution equations . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 2.2 Deep inelastic scattering . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.2.1 Deep inelastic scattering at LO . . . . . . . . . . . . . . . . . . . . . 18 2.2.2 Deep inelastic scattering at NLO . . . . . . . . . . . . . . . . . . . . 20 2.2.3 Deep inelastic scattering at NNLO . . . . . . . . . . . . . . . . . . . 24 2.3 Drell-Yan ..................................... 27 2.4 Treatment of heavy-flavor . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 2.4.1 The FFN and ZM-VFN schemes . . . . . . . . . . . . . . . . . . . . . 30 2.4.2 GM-VFNschemes ............................ 31 2.5 Mellinevolution.................................. 33 2.6 Nuclear parton distribution functions . . . . . . . . . . . . . . . . . . . . . . 42 2.7 A-ZGlobalanalysis................................ 44 2.7.1 Parametrization of nPDFs . . . . . . . . . . . . . . . . . . . . . . . . 45 2.7.2 Experimentaldata ............................ 47 2.7.3 Analysismethod ............................. 51 2.7.4 Results................................... 54 3 Percolation of strings 59 3.1 Geometric scaling for hadronic interactions . . . . . . . . . . . . . . . . . . . 63 3.1.1 The saturation momentum . . . . . . . . . . . . . . . . . . . . . . . . 64 3.1.2 Comparison with experimental data . . . . . . . . . . . . . . . . . . . 65 3.2 The near-side ridge structure . . . . . . . . . . . . . . . . . . . . . . . . . . . 70 3.2.1 The near-side ridge in the SPM . . . . . . . . . . . . . . . . . . . . . 71 3.2.2 Results................................... 74 3.3 Geometric scaling of elliptic flow . . . . . . . . . . . . . . . . . . . . . . . . . 78 3.3.1 Universal scaling law . . . . . . . . . . . . . . . . . . . . . . . . . . . 79 1 2CONTENTS 3.3.2 Discussion................................. 80 3.4 Energy loss as the origin of the scaling law of v2................ 84 3.4.1 EnergyLoss................................ 84 3.4.2 Discussion................................. 87 4 Suppression of high-pTparticles in HICs 91 4.1 Jetquenching................................... 91 4.1.1 Collisional energy loss . . . . . . . . . . . . . . . . . . . . . . . . . . 91 4.1.2 Radiative energy loss . . . . . . . . . . . . . . . . . . . . . . . . . . . 92 4.1.3 Jet quenching models . . . . . . . . . . . . . . . . . . . . . . . . . . . 94 4.2 Particleobservables................................ 97 4.2.1 Nuclear modification factor . . . . . . . . . . . . . . . . . . . . . . . 97 4.2.2 High-pTflowharmonics ......................... 98 4.3 Energy loss in the ASW framework . . . . . . . . . . . . . . . . . . . . . . . 100 4.4 Medium induced gluon radiation . . . . . . . . . . . . . . . . . . . . . . . . . 102 4.4.1 Opacity expansion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 104 4.4.2 Multiple soft scattering . . . . . . . . . . . . . . . . . . . . . . . . . . 104 4.5 Quenchingweights ................................ 105 4.5.1 Expandingmedium............................ 106 4.6 Energy and centrality dependence of ˆq..................... 106 4.6.1 Energy loss approach . . . . . . . . . . . . . . . . . . . . . . . . . . . 107 4.6.2 From hydrodynamics to the transport coefficient . . . . . . . . . . . . 108 4.6.3 Hydrodynamic modeling of the medium . . . . . . . . . . . . . . . . 110 4.6.4 Results................................... 111 4.6.5 Limitations and conclusions . . . . . . . . . . . . . . . . . . . . . . . 116 4.6.6 One step forward: EbyE hydrodynamics . . . . . . . . . . . . . . . . 118 Conclusions 123 Resumen 127 Appendix A About the Mellin technique in global analysis 133 List of Figures 138 List of Tables 139 Bibliography 141 Abstract Heavy ion collisions (HICs) are the fundamental tool to study Quantum Chromodynamics (QCD) under extreme conditions of energy and density, which are very different from those inside the atomic nucleus. In high-energy nuclear collisions critical temperatures and densities that allow the formation of the so-called quark-gluon-plasma (QGP) are reached. The QGP is composed by deconfined quarks and gluons – the degrees of freedom of the QCD Lagrangian density. This plasma is an almost perfect liquid where its constituents interact strongly. Consequently, its study provides a way to better constrain the non-perturbative regime of QCD, which currently in not very well understood. Moreover, the QGP is described by present cosmological models as the state of matter during the first microseconds after the Big-Bang. Therefore, the study of the QGP in heavy ion collisions allows us to analyze the evolution from the deconfined state of the matter existing in the early stages of the Universe to the normal confined matter. In other words, they may be helpful to disclose the origin of mass and confinement. High-energy nuclear collisions are, definitely, a big window to observe, explore, and understand the origin and evolution of our Universe. Nonetheless, the mean life of the hot nuclear matter formed in HICs is very small – of around 100 ys = 10−22 s – and, hence, it cannot be directly detected. Its properties must be indirectly studied in the final hadronic state of the collision. Usually, distributions of soft particles are used to obtain signs of the collective behavior of the QGP and to try to describe it hydrodynamically; while hard probes are employed to study the effect of the medium on processes that can be computed perturbatively. This thesis is focused on the analysis of two types of effects: •Initial stage effects (IS). These are previous to thermalization. A basic ingredient towards the understanding of heavy ion collisions and, hence, the description of the QGP, is the nuclear parton distribution functions (nPDFs). They contain the information about the partonic structure (quarks and gluons) of protons and neutrons bound in the colliding nuclei. They are, consequently, long-distance contributions which are not part of the perturbative domain of QCD. Their precise determination is crucial for the correct interpretation of any observable used in HICs. Thanks to their universality and to the fact that their evolution with respect to a 3 4CONTENTS concrete initial scale can be addressed by perturbative QCD, a technique called global analysis has been developed to extract them. The obtaining of nPDFs at next-to-next- to-leading order in pQCD using global analysis is the subject of chapter 2. Another very interesting issue regarding HICs are the collective phenomena that give rise to the QGP. These are usually examined within the Color Glass Condensate (CGC) framework. A simplified model which encodes some of the properties of the QGP is percolation of strings (SPM). This approach and some of its results compared to available experimental data are presented in chapter 3. •Final state effects (FS). Amongst them it is worth stressing hard probes, which are observables characterized by a high energy or mass. Therefore, they are part of the perturbative sector of QCD. Their behavior in vacuum, i.e, in collisions where a plasma is not created, is well known. Thus, their modifications – with respect to the vacuum case – due to the presence of the nuclear medium are analyzed in order to extract the properties of the QGP. Among these hard probes there is jet quenching. The suppression of jets, and, particularly, the single-inclusive suppression of particles with high transverse momentum in HICs is addressed in chapter 4. Chapter 1 Introduction The existence of a new form of nuclear matter that would be created at very high energies and densities was predicted for the first time in the 1970’s [1–5]. Under these extreme conditions, short range interactions dominate over the long range ones which begin to be screened by the color sources nearby. Short-range interactions are characterized by a small coupling constant – due to asymptotic freedom. Therefore, the QCD matter at high energy and density, usually called quark-gluon plasma (QGP), is composed by deconfined quarks and gluons. This allows us to study QCD under exceptional conditions which cannot be reached in elementary particle interactions. In addition, the QGP was formed in the very early stages of the Universe, which means that the phase transition from confined to deconfined nuclear matter may provide information about the dynamics of the early Universe. This led to the development of the Heavy-Ion Collisions (HICs) programs: the Alternating Gradient Synchrotron (AGS) at the Brookhaven National Laboratory (BNL), the Super Proton Synchrotron (SPS) at Conseil Europ´een pour la Recherche Nucl´eaire (CERN), the Relativistic Heavy Ion Collider (RHIC) at BNL, and the Large Hadron Collider (LHC) at CERN. According to Lattice QCD calculations, this phase transition is a fast and continuous cross-over at temperatures of Tc≃154 MeV (critical temperature), which corresponds to a critical energy density of c∼1 GeV/fm3[6]. Experimental data have shown since long time ago that the energy density achieved in ultra-relativistic Heavy-Ion collisions (HICs) is higher than the critical one [7–9]. Therefore, Heavy-Ion Collisions give a window of opportunity to study the properties of the QGP. The systems produced in this kind of collisions expand in time scales O(10 fm) which are much larger than the typical time scales of individual processes in QCD. During such a long period of time, collective phenomena occur and medium effects are experimentally accessible. The study of processes that are sensitive to the degree of collectiveness of the system is the main purpose of HICs. As it was just mentioned, the mean life of the QGP is really small, of the order of the transverse size of the atomic nucleus: τmean ∼A1/3; that is, τmean ∼6 fm for lead, the nucleus used in HICs at the LHC. Therefore, a direct observation of it is impossible. Thus, only indirect signals can be used to characterize the properties of the QGP. Usually, these signatures are classified in: hard probes, involving scales in the perturbative region of QCD, 5 12 CHAPTER 1. INTRODUCTION percolation threshold. The value of the critical density ηc= 1.13 has been computed by several numerical studies [48–50] using an homogeneous distribution for the colliding nuclei. However, in HICs the distribution of strings in highly non-uniform since there are more nucleons – and, hence, more strings – in the center than in the edge of the nuclei. The results on ηcare higher in this case [51]. This threshold was already reached in Pb-Pb central collisions at the SPS and in semi-peripheral collisions at RHIC energies. The connection between percolation and QCD has been formally studied. These studies are based on the relation between SU(N) gauge field systems and spin [52], and on Polyakov loops [53]. The addition of some dynamics to this geometric picture is required in order to compare with experimental data. Specifically, a model for the decays of clusters into hadrons is needed. This will be addressed in Chapter 3, where the expressions of the average momentum and average multiplicity of the clusters from the composition of the color fields of the strings will be derived. Chapter 2 Nuclear parton distribution functions at NNLO 2.1 DGLAP evolution equations It was already mentioned in the Introduction that parton distribution functions are extracted from experimental data. However, PDFs need to be evaluated at the relevant scale, Q2, of the experiment. The Q2-dependence of parton distribution functions can be determined perturbatively as long as Q2is large enough. The transverse resolution of these processes is set by their virtuality, Q2. Therefore, the Q2-dependence of parton distribution functions is strictly related to the collinear divergence, a singularity due to the emission of gluons at low transverse momentum. It was already remarked that at leading order (nuclear) parton distributions, fi(x, Q2) represent the probability of finding a parton iwith fraction of momentum xof the hadron (free or bound in nucleus) at the energy scale Q2. The Q2-evolution of the distribution is calculated by considering the probability of emitting a parton with transverse momentum close to Q2. The collinear divergences appearing when computing these parton emissions give rise to logarithmic divergences, ln(Q2/ΛQCD). As αsis small (Q2is large), parton emissions are suppressed by powers of the strong coupling constant. Nonetheless, this logarithmic terms can compensate this suppression. These large logarithms are in general obtained from the phase space region where multiple emissions are ordered by transverse momenta, with subsequent emissions having less momenta 1. This brings on a resummation of large logarithms, as the typical ones in the renormalization group equations approach. The Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) evolution equations are [22–24] d dln Q2fi(x, Q) = αs(Q) 2πZ1 x dz zPij (z, αs(Q)) fjx z, Q,(2.1.1) 1This is the case of PDFs, which are space-like (DIS, q2<0). On the contrary, FFs are time-like (q2>0 for e+e−annihilation). 13 14 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO where iand jrun over all the possible parton flavors (2nf+ 1 equations, nfbeing the number of quark flavors). The splitting functions are denoted by Pij. They can be computed in perturbation theory and they describe the probability of the evolution of a parton jinto a parton i. In consequence, in general, they are not symmetric, Pij 6=Pji. The splitting functions themselves can be expressed as a perturbative expansion in powers of αs: Pij(x, αs(Q)) = P( 0) ij (x) + αs(Q) 2πP(1) ij (x) + αs(Q) 2π2 P(2) ij (x) + ... (2.1.2) Due to charge conjugation invariance and flavor symmetry their number can be reduced significantly: Pqiqj=P¯qi¯qj=δikPv qq +Ps qq , Pqi¯qj=P¯qiqj=δikPv q¯q +Ps q¯q , 2nfPqig= 2nfP¯qig=Pqg , Pgqi=Pg¯qi=Pgq , (2.1.3) where v stands for the flavor-diagonal valence quantity, and sfor the flavor-independent sea contribution. The LO splitting functions are [23] Pqq =CF 1 + x2 (1 −x)+ + 2δ(1 −x), Pqg =1 2x2+ (1 −x)2, Pgq =CF"1 + (1 −x)2 x#, Pgg = 2Ncx (1 −x)+ +1−x x+x(1 −x)+11Nc−2nf 6δ(1 −x), (2.1.4) where nfstands for the number of flavors, Ncfor the number of colors and CF= (N2 c− 1)/2Nc. The “+” prescription is defined by Z1 0 dx f(x) (1 −x)+ =Z1 0 dx f(x)−f(1) 1−x.(2.1.5) 2.1. DGLAP EVOLUTION EQUATIONS 15 Higher order contributions to the splitting functions have already been calculated. The splitting functions at NLO can be found in [54,55] . For the NNLO splitting functions see [56,57]. Often, the physical basis of the PDFs is modified to a more appropriate one, which decouples the non-singlet distributions from the gluon. This basis is formed by the aforementioned singlet, non-singlet, valence, and gluon, the three former being, suppressing the functional dependencies, respectively [58] Σ = nf X i (qi+ ¯qi),(2.1.6) q± ij = (qi±¯qi)−(¯qj±qj),(2.1.7) qv ns = nf X i (qi−¯qi).(2.1.8) It can be demonstrated that in this basis the splitting functions are for the non-singlet and the valence, respectively [58]: P± ns =Pv qq ±Pv q¯q ,(2.1.9) Pv ns =Pv qq −Pv q¯q +nf(Ps qq −Ps q¯q)≡P− ns +Ps ns .(2.1.10) Regarding the singlet and the gluon, the expression for the quark-quark splitting function is needed [57] Pqq =P+ ns +nf(Ps qq +Ps ¯qq)≡P+ ns +Pps .(2.1.11) The DGLAP equations in this basis, suppressing the functional dependencies, are d dln Q2qi ns =X n=0 αs 4πn+1 P(n)i ns ⊗qi ns (i=±,v) ,(2.1.12) d dln Q2Σ g=αs 2πPqq Pqg Pgq Pgg ⊗Σ g,(2.1.13) where ⊗stands for the Mellin convolution, [a⊗b](x)≡Z1 x dy ya(y)bx y.(2.1.14) Using these equations, the Q2-evolution with respect to an initial scale Q0of the parton distribution functions is predicted. However, PDFs are unknown a priori at any initial scale, hence, they have to be extracted from a global analysis of (DIS, Drell-Yan...) data. The procedure is as follows, first of all, among all the experimental data available, those which are in the kinematic region where pQCD is valid are selected. Then, all the selected observables are computed at parton level at the order desired in pQCD, in this case, at NNLO. Next, the parton distribution functions are parametrized at the initial scale, Q0. For 16 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO each data point, the PDFs are evolved using the DGLAP equations (including the splitting functions up to NNLO) from the initial scale to the scale of the process, Q > Q0, and, then, the soft part is convoluted with the hard one to obtain the theoretical prediction for the observable. Then a χ2-distribution is built and the value of the parameters is found minimizing iteratively the χ2. It is important to highlight that if the observables are computed including NNLO corrections in QCD, the splitting function up to order α2 shave to be included in the evolution giving rise to terms of order α3 sin the DGLAP equations. In principle, the number of parameters given to the functional form of the PDFs at the initial scale is arbitrary. The value of Q0is also arbitrary and each nPDFs collaboration chooses it according to different criteria. At the initial scale, some of the parameters can be fixed by the proton sum rules, related with the baryon number conservation, Z1 0 dx uvx, Q2 0= 2 Z1 0 dx dvx, Q2 0= 2 ,Z1 0 dx svx, Q2 0= 0 ,(2.1.15) and the momentum sum rule, Z1 0 dx x (Σ (x, Q2 0) + g(x, Q2 0)) = 1 .(2.1.16) 2.2 Deep inelastic scattering As it has been already outlined in the previous sections both PDFs and nPDFs are obtained from several experimental data and then evolved using the DGLAP evolution equations – global analysis. In this section we will discuss the case of lepton-proton(nucleus) collisions involving a large momentum exchange: Deep Inelastic Scattering (DIS) experiments. These experiments usually provide the strongest constraints to PDFs in global PDF fits and one of the most rigorous tests of perturbative QCD, in general. In this process, illustrated in Fig. 2.1, the projectile is a high-energy lepton which scatters off a hadron target: lN →l0X. The four-momentum of the incoming lepton is denoted by kµ= (E,~ k) and the momentum of the outgoing lepton by k0µ= (E0,~ k0) . The momentum of the hadron is pµ= (EN, ~p), which can be written as pµ= (M,~ 0) in the hadron rest frame. The momentum transfer is given by qµ=kµ−k0µ. In the simplest case both the incoming and the outgoing leptons are charged leptons (an electron or a muon) and the interaction is dominated by the exchange of a virtual photon, i.e, the interaction is mainly electromagnetic. However, corrections by the exchange of a Z0boson, i.e, neutral current (NC) interaction, might be relevant, see left panel of Fig. 2.1. (Anti)neutrinos have no electric charge so do not interact electromagnetically, but they do interact with quarks through week interactions. Neutrino DIS refers to the scattering process that involves the conversion of a(n) (anti)neutrino to its associated charged lepton, exchanging a virtual massive vector boson, the W±. In this case, the interaction is called charged current (CC) and it is represented, at leading order (LO), on the right panel of Fig. 2.1 . 2.2. DEEP INELASTIC SCATTERING 17 Figure 2.1: Feynman diagram of the electron-proton NC (ep →eX, left panel) and CC (νep→eX, right panel) at leading order. The standard DIS invariant variables are defined by: Q2≡ −q2, x ≡Q2 2p·q=Q2 2Mν , y ≡p·q p·k=ν E0,(2.2.17) where the latter equalities refer to the target rest frame. Mdenotes the rest mass of the nucleon, ν=E−E0is the energy transfer to the target nucleon. The invariant Q2is called virtuality and it is related to the scale of distances involved in the interaction. The Bj¨orken scaling variable xis the fraction of momentum from the proton (free or bound) carried by the interacting quark. The invariant yrepresents the fraction of the energy of the ingoing lepton transferred to the virtual photon. The unpolarized cross section for both charged-lepton and neutrino DIS can be written in terms of the structure functions. These are adimensional functions which parametrize the structure of the target as ‘seen’ by the virtual photon (or, in general, gauge boson). These functions are F2,FLand F32. In the electromagnetic case – interaction mediated by a virtual photon – only F2and FLappear. Fem 2dominates the cross section and, for this reason, it is the best known of the structure functions, as Fem Lappears at higher orders in perturbation theory. In the case of charged current (or neutrino) DIS the structure functions involved are F2and F3, which describe the interaction with (anti)neutrinos. Due to the weak nature of the neutrino interaction the use of heavy nuclear targets is unavoidable. 2Sometimes instead of FL,F1is used which is just a combination of F2and FL 18 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO 2.2.1 Deep inelastic scattering at LO The charged-lepton electromagnetic DIS cross section, neglecting the contribution of the exchange of a Z0boson, can be written as: d2σ dxdQ2=4πα2 em xQ41−y+y2 2Fem 2(x, Q2)−y2 2Fem L(x, Q2),(2.2.18) where αem ≡e2/4πstands for the fine structure constant and FLis given by Fem L(x, Q2) = Fem 2(x, Q2)1 + 4x2M2 Q2−2xF em 1(x, Q2).(2.2.19) At very high energies, in particular for Q2∼M2 Z, the contribution to the chargedlepton DIS from the boson Zexchange cannot be neglected. Eq. (2.2.18) can be generalized including both γand Zexchanges to: d2σ dxdQ2=4πα2 em xQ4hxy2F1(x, Q2) + (1 −y)F2(x, Q2) + y(1 −y 2)F3(x, Q2)i,(2.2.20) where now the parity-violating structure function F3appears. The unpolarized neutrino (antineutrino) DIS cross section for νp →lX can be written as: d2σν(¯ν) dxdy =G2 FME π1−y−M 2ExyFν(¯ν) 2+y2xF ν(¯ν) 1+ (−)y1−y 2xF ν(¯ν) 3, (2.2.21) where GFis the Fermi constant defined as GF √2=g2 8m2 W .(2.2.22) Let us discuss now the electromagnetic DIS interaction in the naive parton model. The parton model approximation states that for hard enough interactions, the virtual photon interacts only with a single point-like parton inside the target proton and these partons can be treated approximately as free particles. In this model, considering point-like (massless) and spin 1/2 partons (quarks) moving parallel to the hadron (with pT= 0) the EM longitudinal structure functions are exactly zero: Fem L= 0. Hence, in the limit M2/Q2→0, the Callan- Gross relation is satisfied [59]: 2xF1(x) = F2(x).(2.2.23) This relation illustrates a fundamental property of spin-1/2 particles, that they cannot absorb a longitudinally polarized photon [60]. The electromagnetic F1and F2in the parton model are, respectively, Fem 1(x, Q2) = X i e2 iqi(x, Q2) Fem 2(x, Q2) = X i 2xe2 iqi(x, Q2),(i= 1,2, L) (2.2.24) 2.2. DEEP INELASTIC SCATTERING 19 Figure 2.2: Bj¨orken scaling: F2versus Q2for different values of x. Courtesy of the Particle Data Group [65]. where irepresents the sum over all partons and eiare their fractional electric charges. Some comments are in order. First of all, in this approximation the structure functions are independent of the virtuality of the process, Q2. They only depend on one dimensionless variable x. This phenomenon is known as Bj¨orken-scaling [61] and was a remarkable success of the original parton model. It was observed for the first time in the SLAC experiments [62] giving direct evidence of the constituents of the nucleon. The data showing this approximate behavior span two decades of experiments, from the early SLAC measurements to the more recent experiments at HERA [63,64], see Fig. 2.2. The lack of any scale dependence in the structure functions is a consequence of the parton model’s assumptions: partons are treated as point-like particles, and consequently having no characteristic length scale. The simple parton model provided a good phenomenological description of early DIS measurements. The asymptotic freedom of QCD allows for a consistent description of Bj¨orkenscaling, where the constituents of the hadron can be treated as independent, non-interacting point-like particles at high virtuality Q2. The partons in this model were therefore rapidly 20 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO associated with the quarks and gluons of QCD. However, as it can be seen in Fig. 2.2, the Q2-independence of the structure functions is a good approximation, but it is not completely exact. Some deviations from the scaling law appear. These are called scaling violations and they can be perfectly understood by the DGLAP equations. These higher-order corrections to the point vertex cross section introduce logarithms of Q2which break the Bj¨orken scaling of the structure functions. In fact, the observation of such scaling violations was considered as one of the most robust experimental validations of QCD. For the neutrino scattering, the interaction is mediated by the electroweak chargedcurrent, CC, so the quark distributions (PDFs) are weighted by their corresponding weak charges. It has to be taken into account whether the interacting lepton is a neutrino or an antineutrino and the mixture of flavors in the quark sector given by the CKM matrix. When the interacting lepton is a neutrino, νp →lX, the structure function Fν 2at LO is [66]: Fν 2(x, Q2)=2xd0(x, Q2) + s0(x, Q2) + ¯u(x, Q2) + ¯c(x, Q2),(2.2.25) with d0(x, Q2) = |Vud|2d(x, Q2) + |Vus|2s(x, Q2) s0(x, Q2) = |Vcd|2d(x, Q2) + |Vcs|2s(x, Q2),(2.2.26) where Vij is the element ij of the CKM matrix 3. The structure function F3at LO is given by [66]: Fν 3(x, Q2)=2d0(x, Q2) + s0(x, Q2)−¯u(x, Q2)−¯c(x, Q2).(2.2.27) 2.2.2 Deep inelastic scattering at NLO In this section the diagrams that contribute at next-to-leading order (in QCD) to the one studied in the previous section will be summarized. At O(αs) in QCD the structure functions are modified by two kind of processes: initial (or final) state radiation and virtual corrections. The initial (final) state radiation consists of having extra QCD particles in the initial (final) state; i.e, in the radiation of an extra real gluon, as it can be seen in Fig. 2.3. The virtual corrections correspond to the emission and then absorption of a virtual gluon. This can be due to the self-energy corrections of the external legs (see right panel of Fig. 2.4) or to the loop correction to the quark-photon vertex (see left panel of Fig. 2.4). 3When the interacting lepton is an antineutrino, ¯νp →lX,F¯ν 2can be obtained from Eq. (2.2.25) changing the PDF of each particle by the PDF or its antiparticle. 2.2. DEEP INELASTIC SCATTERING 21 Figure 2.3: Feynman diagrams corresponding to the real-gluon emission process γ∗q→qg. Figure 2.4: Virtual corrections to DIS. Collinear,infrared (IR) and ultraviolet (UV), divergences show up when trying to compute the QCD corrections to DIS (at any order in perturbation theory). For instance, in the case of initial state radiation, see Fig. 2.3, collinear divergences appear when the emission is parallel to the direction of the incoming or outgoing quark and infrared divergences happen if the momentum of the emitted particle goes to zero. These divergences can be regulated in a gauge-invariant way using the dimensional regularization methods [67]. All these calculations can be performed in N= 4 −space-time dimensions. Then, all the divergences appear as 1/n-poles. When regularization is appropriately performed, all the divergences, except of the collinear ones, turn out to cancel as a consequence of the IR safety of QCD yielding a finite contribution to the cross section. The latter singularities need to be removed by absorbing them into the parton densities. Including all the contributions both at LO and NLO, the full Fem 2structure function can be written as: Fem 2(x, Q2) = X i xe2 ifi(x) +αs 2πZ1 0 dξ ξX jPij x ξlog Q2 κ2+Wij(x)fj(ξ) + O(α2 s)#.(2.2.28) 28 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO where x1≡√τeyR, x2≡√τe−yR, τ ≡M2 s,(2.3.48) √sbeing the center of mass energy of the two hadrons system. As usual, q(i)(xi) is the parton distribution function. NLO QCD corrections are added to the Born-level Drell-Yan cross section. There are two type of corrections: having a quark and antiquark in the initial state or having a quark (antiquark) and a gluon. In the first case, O(αs) QCD corrections come from gluon radiation, the emission of an additional unobserved gluon off the quark or the antiquark before the annihilation, see Fig 2.6, or from virtual corrections, one-loop diagrams equivalent to the ones in DIS, see Fig 2.7. On the other hand, we have initial state gluons: instead of having a quark and a antiquark in the initial state, we have a quark (antiquark) and a gluon. This gluon decays into a quark-antiquark pair, and the quark (antiquark) is annihilated with the antiquark (quark) from the other hadron into a virtual photon which decays into the dilepton pair, see Fig. 2.8. All the subprocesses corresponding to the NLO QCD corrections to Drell-Yan dilepton production are summarized in Table 2.3. Figure 2.6: NLO gluon radiation in Drell-Yan dilepton production. Process q+ ¯q→µ++ µ−+g. Figure 2.7: NLO Virtual corrections to the Born term: q+ (¯q)→γ∗. The LO cross section is finite, however NLO QCD corrections give rise to infrared, ultraviolet, and collinear singularities. To handle them a regularization procedure is needed. The detailed computation of the NLO Drell-Yan cross section can be found in [81]. The full expression for the Drell-Yan cross section at NLO is given in [82]. 2.3. DRELL-YAN 29 Figure 2.8: Gluon-quark subprocess at NLO in pQCD, g+q(¯q)→µ++µ+q(¯q), in Drell-Yan dilepton production. Drell-Yan NLO subprocesses αs q+ (¯q)→γ∗one-loop (virtual) correction q+ ¯q→γ∗+g g+q(¯q)→γ∗+q(¯q) Table 2.3: List of subprocesses in Drell-Yan dilepton production at NLO. Drell-Yan dilepton production is an essential process to test the Standard Model (SM) in an accurate way at hadron colliders, since it has a large cross section, a clean experimental signature, and it is very sensitive to the properties of the gauge bosons. It is also a fundamental tool for the extraction of the (n)PDFs. Thus, the QCD corrections at NNLO to Drell-Yan are essential to have the precision required for present and future colliders. The processes which contribute to the Drell-Yan cross section at NNLO in pQCD are listed in Table 2.4, where Vstands for the vector boson exchanged, γor Z0boson. As in the case of the NLO corrections, infrared, ultraviolet, and collinear poles show up. A regularization procedure is used and the divergent part is absorbed into the parton densities. For the results for Drell-Yan dilepton production at NNLO see [83,84]. Drell-Yan NNLO subprocesses α2 s q+ (¯q)→Vtwo-loop (virtual) correction q+ ¯q→V+gone-loop correction q+ ¯q→V+g+g q+ ¯q→V+q+ ¯q g+q(¯q)→V+q(¯q) one-loop correction g+q(¯q)→V+q(¯q) + g q(¯q) + q(¯q)→V+q(¯q) + q(¯q) g+g+→V+q+ ¯q Table 2.4: List of subprocesses in Drell-Yan dilepton production at NNLO. 30 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO 2.4 Treatment of heavy-flavor So far in this thesis, the assumption that all the quarks are massless has been made. Nonetheless, the appropriate treatment of heavy quarks in global analysis is essential for the precision demanded by current hadron colliders [95]. The zero-mass limit for all the partons is considered a good approximation at virtualities far above all quark mass thresholds, denoted by mi:Qmi. This obviously does not hold when Q≤mi. A consistent treatment of heavy quark mass effects in pQCD over the whole energy range, from Q≤mito Qmi, is mandatory. Dealing with heavy quarks in pQCD is a delicate task. Traditionally, there have been two kinematical regimes where this treatment could be simplified: the regime where Q2≤mQ and the one where Q2≥mQ, being mQthe mass of the heavy quark. These two regions correspond to the two bounded regimes between which all the current heavy quark schemes try to interpolate. These are, respectively the Fixed Flavor Number scheme (FFNS) and the Zero-Mass Variable Flavor Number scheme (ZM-VFNS). For simplicity, a theory with nflight flavors and a single heavy quark with mass mQwill be considered. 2.4.1 The FFN and ZM-VFN schemes The simplest scheme adopted to include heavy quark effects in global analysis is the Fixed Flavor Number scheme (FFNS). The assumption made is that all the partons in the theory are the nflight quarks, considered as massless, and the gluons. Indeed, the initial state nucleon is assumed to not to have heavy quark component, which is treated as a final state particle. Unfortunately, this model becomes unreliable when the scale of the hard process, Q, becomes larger than the mass of the heavy quark, mQ. In this scheme, when the factorization and the renormalization scales are µr=µf=µ, the DIS structure function in the MS factorization scheme can be written as, Fλnf, m2 H, Q2= nf X i Cλ i nf,Q2 m2 Q ,µ2 m2 Q ,Q2 µ2!,(2.4.49) where the sum is only over light quark flavors. This structure function can be divided in two components: FL, the one where only light flavors are present, and FH, which includes the contribution of the final state heavy quark. F(nf, Q2, m2 Q) = FL(nf, Q2) + FH(nf, Q2, m2 Q),(2.4.50) where FL(nf, Q2) = nf X i Linl,Q2 µ2⊗fi(nf, µ2),(2.4.51) FH(nf, Q2, m2 Q) = nf X i Hi nf,Q2 m2 Q ,µ2 m2 Q ,Q2 µ2!⊗fi(nf, µ2).(2.4.52) 2.4. TREATMENT OF HEAVY-FLAVOR 31 Here the Wilson coefficients including heavy quark lines are denoted by H, and Lare the diagrams that do not. The expression (2.4.49) is reliable in the heavy quark mass threshold region and below. In this scheme, the Wilson coefficients contain unresummed logarithms of the ratio Q2/m2 Q which can become large at scales much larger than the heavy quark mass. In the Zero-Mass Variable Flavor Number (ZM-VFN) prescription these problems are solved by treating the heavy quark as a massless parton above its mass threshold. That is that the parton distribution of the heavy quark remains zero if Q2≤m2 Q, but it follows DGLAP evolution when Q2> m2 Q. The renormalization of the heavy quark PDF resums all the logarithms via the DGLAP equations. Thus, the number of flavors is increased by one when crossing the heavy quark mass threshold in this approach. In the ZM-VFN the structure function is F(nf+ 1, x, Q2) = nf+1 X i Cinf+ 1,Q2 µ2⊗fi(nf+ 1, µ2).(2.4.53) This method solves the problems that arises at large scales in FFNS. However, the treatment of heavy quarks is made only in terms of massless partons: it completely ignores the massive contributions to the Wilson coefficients. Hence, it is not accurate when the powers m2 Q/Q2are significant. Despite the simplicity of the procedure, it has been used in several global analysis [90,96]. 2.4.2 GM-VFN schemes The General Mass Variable Flavor Number schemes (GM-VFNS) have been developed in order to avoid the shortcomings of both the FFN and ZM-VFN schemes. In such approaches the treatment of heavy quarks is usually reduced to the FFN scheme at low scales and the ZM-VFNS procedure at high scales, with some interpolation between these two regions. The GM-VFN schemes require that the ZM-VFN and FFN calculations coincide at large scales, where the heavy quark mass dependence of the FFN Wilson coefficients can be neglected, FL(nf, Q2) + lim Q2m2 QFH(nf, Q2, m2 Q)=F(nf+ 1, x, Q2).(2.4.54) This constraint implies that the PDFs in the two schemes need to be related by a perturbatively computable transformation: fi(nf+ 1, µ2) = nf X j Aij nf,µ2 m2 Q!⊗fj(nf, µ2),(2.4.55) where the A’s are not square matrices: iruns over the nf+ 1 partons in the ZM-VFNS, and jruns over the nfpartons in the FFNS. These matrices have been determined up to NNLO 32 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO in αSin Refs. [97,98]. In practice, the GM-VFN schemes are chains of FFN-type schemes with increasing nf as the scale increases over each quark mass threshold, requiring the physical observables to be continuous across these thresholds. The matching condition at the heavy quark mass threshold for the structure function is FGM(m2 Q) = nl X j CFF jnf, m2 Q⊗fj(nf) = nf+1 X i CVF inf+ 1, m2 Q⊗fi(nf+ 1) .(2.4.56) Introducing Eq. (2.4.55) in Eq. (2.4.56) one gets: CFF jnf, m2 Q= nf+1 X i CVF inf+ 1, m2 Q⊗Aij nf, m2 Q,(2.4.57) which is the minimal GM-VFN scheme [99]. This last relation has to be satisfied order by order in αs. To illustrate this, let us consider Eq. (2.4.57) at NLO for the gluon: CNLO g(nf+ 1, mQ) = CNLO g(nf, mQ)−CLO H(nf+ 1, mQ)⊗ALO Hg(nf, m2 Q).(2.4.58) This is the expression used to define the original ACOT (M. A. G Aivazis, J. C. Collins, F. I. Olness, and W. K. Tung) scheme [100]. The rightmost term on the r.h.s. of Eq. (2.4.58) when moved to the l.h.s. is known as the subtraction term and allows to define the coefficient CNLO g(nf, mQ) in the ACOT scheme. The subtraction term allows the elimination of the IR-unsafe logarithms present in the FFN approach. However, there is ambiguity to define different GM-VFNS, which has resulted in the existence of several prescriptions. This freedom of definition comes from the existence of terms proportional to powers of mQ/Q that can be interchanged between the Wilson coefficients in Eq. (2.4.58) without changing the final value of the physical observable, i.e., of the structure function. The already-mentioned ACOT scheme does not try to exploit this degeneracy. The Simplified ACOT scheme, or S-ACOT [101, 102], chooses a much simpler prescription of the GM-VFNS. In the S-ACOT calculation the heavy-quark-initiated subprocesses are chosen to be equal to the ZM formulae. Therefore, only the full mQ-dependent contributions from the light-parton-initiated contributions have to be computed. Another choice is the TR type schemes by R. S. Thorne and R. G. Roberts [103, 104], which also require the scale derivatives of heavy flavor structure functions, dF2/d ln Q2, to be continuous at the transition point. More recent refinements to this prescription can be found in Refs. [105–107]. Let us summarize the TR prescription in [107] since it is the one implemented in our global analysis. It has been shown that the choice of coefficient functions removes the ambiguity in defining a GM-VFNS. However, there are additional ambiguities: the ordering of FH 2(x, Q2) is different for the nfand nf+ 1 regions. This is illustrated in the following table where the 2.5. MELLIN EVOLUTION 33 expressions order by order both below and above the heavy quark mass threshold are shown. nf−flavor nf+ 1 −flavor LO αs 4πCFF,1 2,Hg ⊗gnfCV F,0 2,HH ⊗(h+¯ h) NLO αs 4π2(CFF,2 2,Hg ⊗gnf+CF F,2 2,Hq ⊗Σnf)αs 4π(CV F,1 2,HH ⊗(h+¯ h) + CV F,1 2,Hg ⊗gnf+1 ) NNLO αs 4π3PiCFF,3 2,Hi ⊗fnf iαs 4π2PjCV F,2 2,Hj ⊗fnf+1 j Being hthe parton distribution function of the heavy quark. This means that switching directly from a fixed order with nfactive quarks to a fixed order with nf+ 1 active quarks leads to a discontinuity in FH 2(x, Q2). ACOT scheme uses the same order of αsabove and below the transition point, for instance, at NLO αs 4πCFF,1 2,Hg ⊗gnf→CV F,0 2,HH ⊗(h+¯ h) + αs 4π(CV F,1 2,HH ⊗(h+¯ h) + CFF,1 2,Hg ⊗gnf+1).(2.4.59) Then, the structure function is automatically continuous. However, as CF F,1 2,Hg contains only information on P0 qg not on P1 qg, there is effectively LO evolution below the heavy quark mass threshold and NLO evolution above it. Hence, the slope d FH 2(x, Q2)/d ln Q2is discontinuous. In the TR scheme [103, 104] – and all its revisions – has the same order above and below the transition point is used, but it adds a term independent of Q2above the transition point to preserve the continuity of the structure function αs(Q2) 4πCFF,1 2,Hg (Q2/m2 Q)⊗gnf(Q2)→αs(M2) 4πCFF,1 2,Hg (1) ⊗gnf(M2) (2.4.60) +CV F,0 2,HH (Q2/m2 Q)⊗(h+¯ h)(Q2). The O(αs) term is frozen through Q2=m2 Q. Consequently, the definition of the ordering is consistent within each region, except for the addition of a constant term above Q2=m2 Q which does not affect the evolution. This implies that in order to implement a GM-VFNS at NNLO ,O(α3 s) heavy-flavor coefficient functions are needed for Q2≤m2 Qand that their contribution will be frozen for Q2≥m2 Q. The heavy-flavor Wilson coefficients for DIS have been computed up to two-loop order in [108–110]. Recently, O(α3 s) contributions to these coefficients have been computed in the asymptotic region of large momentum transfer, i.e. for Q2m2 Q[111, 112]. The full expressions up to two-loop corrections and the asymptotic ones at three loops have been implemented in our code. The massive Wilson coefficients of F2,FLand F3for chargedcurrent DIS in the same limit can be found in [113,114]. Most of the expressions regarding the heavy-quark treatment cannot be analytically transformed into Mellin space, hence, we have pre-evaluated their momenta and saved them in grids to be read as needed. 2.5 Mellin evolution In Section 2.1, the DGLAP evolution equations, Eqs. (2.1.13) and (2.1.12), were presented. These integro-differential equations cannot be solved analytically in the x-space. Because of 34 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO this, several techniques to accomplish a numerical solution have been developed. Some of the most common methods employed to solve numerically the DGLAP equations in the x-space will be introduced in the following, playing special attention to the Mellinevolution, that is the approach used in this thesis. As it was already indicated, solving the DGLAP equations in the x-space analytically is impossible. One way to deal with this in the force-brute method [85, 86]. This consists in discretizing the equation by a finite difference method. The code is very simple, however, the computation time is large due to the large number of steps (in xand in t) required in order to obtain an accurate evolution. This method was very popular in the 90’s. The Laguerre technique [87] is one of the fastest methods. Splitting functions and PDFs are expanded in terms of Laguerre polynomials. The DGLAP equations turn out to be a simple summation of Laguerre-expansion coefficients, so, it is a very efficient numerical method. This method has some disadvantages. First, obtaining an accurate evolution at large xis not easy. And, second, since the Laguerre polynomials Ln(−ln x) diverge when x goes to zero, it could have convergence issues in the small xregion, which is an important region in high-energy physics. One of the most established approaches nowadays in the x-space is given by the QCD evolution package HOPPET [88]. The splitting functions are represented by their convolution with a set of piecewise polynomial basis functions. Then, Runge-Kutta techniques are used for the evolution in Q2. Good speed and accuracy are obtained by this method. Another popular method is a semianalytical one involving the construction of an evolution operator expressed in form of a rapidly convergent series of matrices, depending only on the splitting functions [89]. This technique is used, for example, in EPS09 NLO nPDFs [90]. The Mellin transformation method is one of the most popular evolution techniques. As it was said in the previous section, the DGLAP equations, Eqs. (2.1.13) and (2.1.12), contain a Mellin convolution. When a Mellin transform is applied to these equations, their right-hand side becomes a simple multiplication of two (Mellin) moments: the moment of the splitting functions, and the moment of the distribution functions (see Appendix A). Therefore, using the Mellin technique, the DGLAP integro-differential equations in x-space are converted into ordinary differential equations that can be solved analytically to all orders. The main inconvenient of this method is that the functional form of the parton distribution functions at the initial scale, Q0, is restricted to have an easily computed Mellin moment. In general lines the application of this method to global analysis is the following: •The (n)PDFs at the initial scale are parametrized in x-space and then the Mellin transform is performed. Flexible polynomial forms are used for the parton distributions at the initial scale in order to be able to calculate their Mellin moments easily. 2.5. MELLIN EVOLUTION 35 •The evolution of the (n)PDFs from the initial scale, Q0, to the scale of the experiment, Q>Q0, is also done in the moment-space. The moment of the splitting functions – anomalous dimensions – are well known up to NNLO [56,57]. Then, obtaining an analytic solution for the DGLAP evolution equations in the Mellin-space is straightforward. •Once the data to be used have been chosen, the corresponding observables (crosssections or structure functions) are computed in Mellin space. In the particular case of contributions that cannot be analytically transformed into momentum space, we preevaluated the momenta and saved them in grids to be read as needed. The principal advantage of this is that all the time-consuming integrations are already dealt with in the determination of the perturbative cross sections. Then, thanks to factorization (see Section 1.2) the theoretical observable can be calculated (in moment-space). Usually, DIS experimental data are in terms of the structure functions, so, instead of the cross sections, the structure functions are used. As it was shown in previous sections, the structure functions are written in terms of the Wilson coefficients and the parton distribution functions. Because of this, the expressions at NNLO for the Wilson coefficients in Mellin-space have to be computed. •A Mellin inversion of the total cross section is performed to go back to x-space and obtain the theoretical prediction of the observable. •Finally, a χ2-distribution is built and the value of the parameters is found iteratively minimizing the χ2. In Appendix A, further details about the Wilson coefficients in moment-space and the anomalous dimensions can be found. The Mellin inversion is also explained in this appendix. As it was mentioned, the DGLAP equations can be solved analytically in Mellin-space. The non-singlet DGLAP equation has an analytical exact solution at NNLO, which is derived in Section 8 of [91]. Regarding the singlet (and gluon) equations, there are no exact solutions in moment-space beyond those known at LO. However, a solution can be constructed using a logarithmic ansatz. This logarithmic series can be improved in order to capture higher order contributions in the truncated solution, a feature that can be very appealing for phenomenological purposes [91]. In our case, the higher order logarithmic approximation of the NNLO singlet solution is used, see Section 10 of [91]. This solution has corrections of O(α3 s) in the evolution, but the anomalous dimensions used go up only to O(α2 s). Using O(α3 s)-terms in the evolution improves the result comparing to the one truncated at O(α2 s). To see the expression of this solution see Section 10 of [91]. In summary, we have implemented from scratch in our code the exact solution of the DGLAP equations at NNLO in pQCD for the non-singlet and the truncated evolution of the singlet including these O(α3 s)-terms. This implies that all the expressions of the anomalous dimensions at NNLO in pQCD have also 36 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO been implemented. Our evolution, both for the non-singlet and the singlet sectors, has been successfully crosschecked with the evolution code QCD-PEGASUS [92]. QCD-PEGASUS results have been carefully cross-checked up to Q2= 10000 GeV2with [93]. In Figs. 2.9, 2.10, 2.11, and 2.12, free-proton PDFs are plotted versus xfor Q2= 1 GeV2,Q2= 10 GeV2,Q2= 100 GeV2, and Q2= 10000 GeV2, respectively. Black points correspond to MMHT2014 NNLO free-proton PDFs [94], whose evolution is given by QCD-PEGASUS [92]. The different curves represent our evolution for A= 1, that is, for the proton, whose PDFs have been parametrized according to MMHT2014 at the initial scale Q2 0= 1 GeV2. Dashed lines correspond to a Mellin contour given by 144 moments and φ= 3π/4. The different colors stand for the different values of the parameter c: black c= 1, red c= 0.9, green c= 0.8, and blue c= 0.7 (see Appendix A). For continuous lines, the same criteria apply, but, in this case, the contour is given by 136 Mellin-moments. As it can be seen in these figures, there is a good agreement between our evolution and the one given by QCD-PEGASUS [92] in the wide range in virtuality considered 4. In order to choose the most appropriate Mellin contour, we plot in Fig. 2.13 the ratio of the curves in Fig. 2.9. The same color criterion has been employed. As it can be seen in Fig. 2.13 the agreement between our evolution and that of QCD-PEGASUS [92] is excellent. Only at very large and very low x, some discrepancies arise. On the one hand, in the (very) large-xregion, parton densities are numerically very small, typically falling as (1 −x)βwith β≥3. Moreover, in this region there are not any experimental data which could constrain the PDFs. On the other hand, at low-x, deviations show up for x≤10−4. However, there are no data constrains at all in nPDFs fits for x≤10−3. Hence, all the available sets in this region are model-dependent, since the nPDFs at low-xare mostly extrapolations which depend on the chosen parametrization. Taking into account the current level of precision of nPDF analysis, the differences that arise in the evolution in this region are negligible. The Mellin contour that was finally implemented in our code is that corresponding to the dashed blue line, i.e., the values φ= 3π/4 and c= 0.8 and 144 Mellin-moments are taken. 4In Figs. 2.9, 2.10, 2.11, and 2.12, our evolution of charm and bottom PDFs are plotted compared to that obtained by MMHT2014 global analysis [94], where the GM-VFNS in the TR scheme was also employed. 2.5. MELLIN EVOLUTION 37 0 0.2 0.4 0.6 0.8 10 -6 10 -4 10 -2 1 -0.1 0 0.1 0.2 0.3 0.4 10 -6 10 -4 10 -2 1 -0.005 0 0.005 0.01 0.015 10 -6 10 -4 10 -2 1 0 0.2 0.4 0.6 0.8 10 -6 10 -4 10 -2 1 0 0.2 0.4 0.6 0.8 10 -6 10 -4 10 -2 1 0 0.1 0.2 0.3 0.4 10 -6 10 -4 10 -2 1 -0.1 -0.05 0 0.05 0.1 10 -6 10 -4 10 -2 1 -0.1 -0.05 0 0.05 0.1 10 -6 10 -4 10 -2 1 Nmom=136, C=1 Nmom=136, C=0.9 Nmom=136, C=0.8 Nmom=136, C=0.7 Nmom=144, C=1 Nmom=144, C=0.9 Nmom=144, C=0.8 Nmom=144, C=0.7 uv x x x dvsv us x x x proton PDFs at Q2=1 GeV2 dsss c x x x b g -30 -20 -10 0 10 -6 10 -4 10 -2 1 Figure 2.9: NNLO-evolution at Q2= 1 GeV2for the different PDFs versus x. Dashed lines correspond to 144 Mellin-moments and continuous lines to 136. Different colors correspond to the different values of the parameter c: black c= 1, red c= 0.9, green c= 0.8, and blue c= 0.7 (see Appendix A). Black points correspond to MMHT2014 NNLO free-proton PDFs [94]. 44 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO •EMC effect. For 0.3.x.0.8, FA 2/Fd 2decreases reaching a minimum at x∼0.6 and then increases again. The EMC effect suggests that in the region where 0.3.x.0.8 – dominated by valence quark distributions for DIS – the valence quark distribution of free nucleons is larger than the one of bound nucleons. Many models successfully explain this phenomenon: nuclear binding or pion cloud effects [126,127]; combinations of various different models [128]; and two-nucleon short-range correlations [129]. •Fermi-motion. For x&0.8, FA 2/Fd 2>1. Nucleons in nuclei move with an average momentum of kF. For nuclear experimental data, the xvariable is determined in the approximation that the nucleon is at rest. In consequence, the experimental measured structure function F2in case of nuclear DIS is the convolution of the structure function, F2, of a free nucleon and the momentum distribution of the nucleon in the nucleus. However, since there are no data in the region of high x, global analysis impose to the nuclear parton distributions, fA i:fA i→0 and fA i/fi→L, where L1 when x→1. 2.7 A-Z Global analysis In this section my work in collaboration with Dr. Pia Zurita (Brookhaven National Laboratory, Upton, NY) is presented. We have performed a global analysis in QCD of collinearly factorized nuclear parton distribution functions and their uncertainties. In this analysis, available data of charged-lepton deep inelastic scattering as well as neutrino deep inelastic scattering are included. The global fit is done at NNLO in pQCD, which means that all the observables are calculated at NNLO, as explained in Section 2.2, and that also the evolution is performed at NNLO in QCD, see Sections 2.1 and 2.5. LHC data are not included in this study. However, several articles have indicated the constraining power of these recent measurements at the LHC [130, 131] in the context of reweighting [132–135], and a new global analysis of nPDFs at NLO, EPPS16, has recently incorporated them [136]. In the present work, Drell-Yan experimental results on the structure functions are being implemented and will be presented in a forthcoming publication. The treatment of heavy quarks follows the GM-VFNS approach explained in Section 2.4. This is the first global analysis of nPDFs performed at NNLO within the GM-VFN scheme 5. The analysis is performed using the Hessian method, that will be explained, together with the estimation of the uncertainties, in Subsection 2.7.3. The extraction of nPDFs from experimental data is fundamental for the understanding of nuclear modifications outlined in the previous section. On the other hand, nuclear parton distributions are essential for the analysis of a wide variety of nuclear experiments, such as heavy ion collisions at the Large Hadron Collider (LHC) and at the Relativistic Heavy-Ion 5There is a previous set of NNLO nPDFs [96], but the simplistic ZM-VFNS is used. 2.7. A-Z GLOBAL ANALYSIS 45 Collider (RHIC). Moreover, the precise knowledge of nPDFs will be mandatory for future colliders, as the EIC and the LHeC. In the last few years, there has been a lot of improvement on the determination of nPDFs, both from the experimental and the theoretical sides. Furthermore, global fits also permit to study the range of applicability of the factorization theorems [19,20] and the universality of nPDFs. Finally, free-proton global fits make use of neutrino-nucleus DIS data in order to study flavor decomposition, hence, nPDFs are needed also in free-proton PDFs analysis. However, nPDFs are much less known than free proton PDFs, especially due to the wider diversity of data – that cover a larger kinematical range – available for the proton community, as it can be seen in Fig. 2.16. Nowadays, the available nuclear (charged-lepton) DIS data go from x≈0.01 to x≈1 and they are often a ratio of the structure function of the studied nucleus, A, with respect to the structure function of deuterium – or sometimes from another nucleus B. These data are basically sensitive to the valence quarks, though for x.0.01 some sensitivity to the sea quarks arises. Gluons are almost not constrained by DIS – and DY – experiments. Many analysis such as HKN07 [137], EPS09 [90], DSSZ12 [138], nCTEQ15 [139], and the recent EPPS16 [136], have included among their data inclusive pion production in d-Au collisions at RHIC [140, 141], since this observable has a potential to better constrain the gluon distribution. However, inclusive pion data are different to the previous mentioned sets – DY, and neutral and charged current DIS – in the sense that they have an additional dependence on the fragmentation functions. As the parton-to-pion fragmentation functions may experience a nuclear modification [21], the interpretation of the inclusive pion production observable is nowadays ambiguous. That is the reason why this observable has not been included in our analysis so far. With respect to neutrino DIS experiments, they may be helpful to constrain light quark flavor, since they provide an electroweak observable. These data have been controversial in the past, as it appeared to show some tension with (charged-lepton NC) DIS data [143]. Nonetheless, this seemed to be mostly a normalization issue and neutrino data have been employed in DSSZ12 [138] and EPPS16 [136] with no difficulties. This section is organized as follows: the parametrization employed for the nuclear modifications of PDFs is described in the next section. Then, the experimental data used in our analysis are presented. In Subsection 2.7.3 the analysis procedure is summarized and, finally, in subsection 2.7.4 our results are presented and discussed. 2.7.1 Parametrization of nPDFs Following the typical framework used in most of the nPDFs global fits, the bound proton PDF for a mass number Aand parton species i,fA i(x, Q2), is defined relative to the corresponding free proton PDF, fp i(x, Q2), as fA i(x, Q2 0) = RA i(x, Q2 0)fp i(x, Q2 0),(2.7.61) 46 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO Figure 2.16: (Left panel) Typical kinematic range of data used in nPDFs global analysis. (Right panel) Typical kinematic range employed in free-proton global analysis. Figures taken from [142]. where RA i(x, Q2 0) denotes the nuclear modification. Our free proton baseline is MMHT2014 NNLO [94], which is defined in a GM-VFNS to deal with heavy quark mass effects. Consistently, we adopt the same prescription for the treatment of heavy quark effects. The nPDFs are then obtained by Eq. (2.7.61) at an initial scale of Q0by determining the nuclear modification factors, RA i(x, Q2 0), from the experimental data. In order to parametrize the RA i(x, Q2 0) in Eq. (2.7.61), we assume isospin invariance for bound protons and neutrons, that is, up=dn, and dp=un. Therefore, the uquark density in a nucleus Awith Zprotons – and A−Znucleons – at a scale Q2is given by uA(x, Q2) = Z AfA u(x, Q2) + Z−A AfA d(x, Q2),(2.7.62) and similarly for dA, ¯uA, and ¯ dA. In contrast to all the nPDFs sets available, except for HKN07 [137], nuclear effects in deuterium are considered here. Due to the lack of diversity of the experimental data available, nuclear modifications for each parton flavor cannot be independently determined. At this stage, as only charged lepton and neutrino DIS data are employed in our fit, we define only three nuclear corrections at the parametrization scale Q0:RA Vfor both valence quark distributions, RA Sfor all the sea quarks, and RA gfor gluons. This is a standard way to proceed that has been used in many nuclear global fits, such as EPS09 [90], and DSSZ12 [138]. However, the most recent analyses, nCTEQ15 [139] and EPPS16 [136], allow some flavor separation for the valence quark distributions in the case of nCTEQ15 [139], and complete separation for the valence and some separation for the sea quarks in the case of EPPS16 [136]. Considering flavor 2.7. A-Z GLOBAL ANALYSIS 47 separation in our fit would imply a larger number of parameters and, at this stage, due to the lack of sensitivity of the implemented data, they would not be well-determined making the uncertainties huge. Our parametrization of the nuclear corrections at the initial scale Q0is given by RA V=NV(A) + (1 −AaV) (1 −x)βVb(A)xαV+c(A)x2αV+d(A)x3αV, RA S=NS(A) + (1 −AaS) (1 −x)βSb(A)xαS+c(A)x2αS+d(A)x3αS, RA g=Ng(A) + (1 −Aag) (1 −x)βgbgx+cgx2,(2.7.63) where the dependencies on Aof the parameters are indicated and βgis fixed to βg= 0.1. In the case of NV(A), NS(A), and Ng(A), the Adependencies are Ni(A) = Aγi(i=V, S, g).(2.7.64) The coefficients b(A), c(A) and d(A) depend on A. The dependence on Aof d(A) is given by d(A) = d1Ad2.(2.7.65) The baryon number and momentum sum rules given by Eq. (2.1.15) and Eq. (2.1.16), respectively, allows us to determine b,c,NV, and Ngfor each nucleus separately. Therefore, the latter quantities are not parametrized. Thus, we have in total 12 free parameters. Above the initial scale, Q > Q0, nuclear PDFs are obtained solving the DGLAP evolution equations with 3-loop splitting functions in the Mellin space, as explained in Sections 2.1 and 2.5. Our initial scale is Q2 0= 1 GeV2, consistently with that of MMHT2014 NNLO [94]. The charm and bottom quark masses are respectively: mb= 4.75 GeV, and mc= 1.4 GeV. The factorization, µF, and renormalization, µR, scales are taken to be µ2 F=µ2 R=Q2. The value of the strong coupling constant is set at the renormalization scale as αs(Q0) = 0.46797. 2.7.2 Experimental data In this first study we restrict ourselves to fixed-target neutral and charged current leptonnuclei deep-inelastic scattering structure functions (or cross sections). Regarding charged-lepton DIS, measurements from a wide variety of experiments, such as NMC [119, 120, 144, 145], SLAC139 [117], and EMC [146] are incorporated in this analysis. Guided by the free-proton baseline fit MMHT2014 [94], the kinematical cuts applied on these data are: Q2≥2 GeV2and W2≥15 GeV2. These data sets – except of those referring to deuteron targets – are listed in Table 2.5, as well as the number of data points that survive our kinematical cuts and the corresponding publication reference. As it was said in the previous section, deuterium is usually considered as approximately free. In fact, only HKN07 nPDFs global analysis [137] have considered any nuclear correction to deuterium. 48 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO However, some nuclear effects exist [116] and, deuterium data have been incorporated in several free-proton PDF global fits, such as MMHT2014 [94], and NNPDF3.0 and subsequent analysis [147]. These analysis have shown the power of using deuterium data in order to separate the uand ddistributions at moderate values of x[94]. Therefore, in this study, DIS on deuteron targets measurements have been incorporated. These data are shown in Table 2.6. On the other hand, results on neutrino DIS using either lead or iron targets from CDHS [148] and CHORUS [149] collaborations, respectively, are also included. The kinematical cuts set on these data are: Q2≥5 GeV2and W2≥25 GeV2, according to MMHT2014 [94]. These data are summarized in Table 2.5. In Tables 2.5 and 2.6 the data sets employed are listed. In most of the cases the observable used is a ratio of the structure function F2of the nucleus A– or cross-section – with respect to that of deuterium or that of a lighter nucleus B. However, for neutrino DIS and, occasionally, for charged-lepton DIS, only measurements of the absolute structure function F2– or cross section – are available. This is the reason why the cuts imposed here are more restrictive – especially for charged-lepton DIS – than those used in other nPDFs global analysis which only include observables in form of ratios where higher-twist effects may cancel. The notations RAand R0Aused in these tables designate the following observables RA=FA L+FA 24M2x2 Q2 FA 2−FA L ,(2.7.66) R0A=FA L FA 2(1 + 4M2x2 Q2)−FA L ,(2.7.67) where Mis the mass of the nucleon. Table 2.5: Data sets included in the present analysis (except deuterium data), listed in order of growing nuclear mass number. In the second column the observable used in our analysis is indicated, when needed the nuclear mass number is as well specified. In the third and fourth columns the number of data points that survive our kinematical cuts and their contribution to the χ2are – respectively – shown. In the last column the corresponding publication reference is indicated. Experiment Observable Ndat χ2Reference DESY HERMES σHe(3)/σd43 48.6 [150] SLAC E-139 σHe(4)/σd2 0.7 [117] CERN NMC 95, re. FHe(4) 2/Fd 215 13.9 [119] CERN NMC 95 FLi(6) 2/Fd 214 15.8 [144] CERN NMC 95, Q2dependence FLi(6) 2/Fd 2132 151.0 [144] SLAC E-139 σBe(9)/σd2 0.3 [117] CERN NMC 96 FBe(9) 2/FC(12) 215 4.5 [120] CERN BCDMS FC(12) 2162 422.0 [151] CERN EMC FC(12) 2/Fd 29 8.5 [152] 2.7. A-Z GLOBAL ANALYSIS 49 Experiment Observable Ndat χ2Reference CERN EMC σC(12)/σd8 20.5 [153] CERN EMC FC(12) 231 37.3 [154] SLAC E-139 σC(12)/σd1 0.1 [117] Fermilab E665 σC(12)/σd4 8.6 [118] CERN NMC 95, re. FC(12) 2/Fd 215 23.0 [119] CERN NMC 95 FC(12) 2/Fd 215 20.6 [144] CERN NMC 95, Q2dependence FC(12) 2/Fd 2145 149.8 [144] CERN NMC 95, re. FC(12) 2/FLi(6) 218 18.7 [119] SLAC E-143 R0C(12) 5 2.8 [155] CERN BCDMS FN(14) 2/Fd 29 26.8 [156] DESY HERMES σN(14)/σd42 36.0 [150] SLAC E-049 σAl(27)/σd3 2.0 [157] SLAC E-139 σAl(27)/σd2 0.1 [117] CERN NMC 96 FAl(27) 2/FC(12) 215 7.1 [120] CERN EMC σCa(40)/σd8 8.5 [153] CERN EMC FCa(40) 232 24.2 [154] SLAC E-139 σCa(40)/σd1 0.001 ∼0 [117] Fermilab E665 σCa(40)/σd4 7.7 [118] CERN NMC 92 RCa(40) −RC(12) 4 1.4 [158] CERN NMC 95, re. FCa(40) 2/Fd 214 17.6 [119] CERN NMC 95, re. FCa(40) 2/FLi 218 11.0 [119] CERN NMC 95, re. FCa(40) 2/FC(12) 218 20.5 [119] CERN NMC 96 FCa(40) 2/FC(12) 215 9.2 [120] CERN BCDMS FFe(56) 2/Fd 26 14.7 [156] CERN BCDMS FFe(56) 2/Fd 210 25.0 [159] CERN EMC FFe(56) 2195 352.5 [161] SLAC E-087 σFe(56)/σd1 0.3 [160] SLAC E-139 σFe(56)/σd5 2.3 [117] CERN NMC 96 FFe(56) 2/FC(12) 215 14.3 [120] CERN EMC FCu(64) 2/Fd 29 4.7 [152] CERN EMC FCu(64) 2/Fd 219 13.8 [146] DESY HERMES σKr(84)/σd34 35.8 [150] SLAC E-139 σAg(108)/σd1 0.1 [117] CERN EMC FSn(117) 2/Fd 28 17.6 [152] CERN NMC 96 FSn(117) 2/FC(12) 215 6.2 [120] CERN NMC 96, Q2dependence FSn(117) 2/FC(12) 2138 108.0 [145] CERN NMC 96 R0Sn(117) −R0C(12) 13 4.7 [145] Fermilab E665 σXe(132)/σd3 1.4 [162] SLAC E-139 σAu(197)/σd2 0.5 [117] 50 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO Experiment Observable Ndat χ2Reference CERN CDHS σνFe(56) 300 256.1 [148] CERN CDHS σ¯νF e(56) 302 256.3 [148] CERN NMC 96 FPb(208) 2/FC(12) 215 6.1 [120] Fermilab E665 σPb(208)/σd4 15.9 [118] CERN CHORUS σνP b(208) 292 196.7 [149] CERN CHORUS σ¯νP b(208) 292 287.5 [149] Table 2.6: DIS on deuteron targets data sets included in the present analysis. In the second column the observable used in our fit is indicated. In the third and fourth columns the number of data points that survive our kinematical cuts and their contribution to the χ2are – respectively – shown. In the last column the corresponding publication reference is indicated. Experiment Observable Ndat χ2Reference CERN EMC 2Fd 2/Fp 266 19.2 [163] CERN EMC 2Fd 235 34.5 [154] CERN BCDMS Rd9 6.8 [164] CERN BCDMS Fd 2246 114.5 [164] CERN BCDMS 2Fd 2/Fp 2−1 11 21.3 [165] Fermilab E665 σd/σp4 3.9 [166] Fermilab E665 2Fd 2/Fp 2−1 7 3.4 [167] Fermilab E665 Fd 2/Fp 253 44.3 [168] SLAC E-140 R0d1 2.8 [169] SLAC E-140 Fd 21 0.005 ∼0 [169] CERN NMC Fd 2/Fp 2148 175.6 [170] CERN NMC Rd/Rp13 13.3 [170] DESY HERMES σd/σp21 9.5 [171] Total (of both tables) 3115 3188.4 Charged-lepton DIS The structure function FA 2of a nucleus with atomic number Zand mass number Ais given by FA 2=Z AFp,A 2+A−Z AFn,A 2,(2.7.68) where Fp,A 2and Fn,A 2are the structure functions of the bound protons and neutrons. However, this is not the observable reported in the original publications of charged-lepton DIS. Instead, the isoscalar corrected structure function was defined as that containing equal number of protons and neutrons: ˆ FA 2=1 2Fp,A 2+1 2Fn,A 2.(2.7.69) 2.7. A-Z GLOBAL ANALYSIS 51 The leading cause of using this definition was to eliminate the effects arising for the different number of protons and neutrons of the considered nucleus Awhen comparing to deuteron structure functions in such a way that ˆ FA 2/Fd 2directly illustrated nuclear effects on the PDFs. The structure function F2can be expressed in terms of the isoscalar corrected structure function as FA 2=ˆ FA 2δ , (2.7.70) where δ=2 A Z+ (A−Z)Fn,A 2 Fp,A 2 1 + Fn,A 2 Fp,A 2 .(2.7.71) The ratio Fn,A 2/Fp,A 2is considered to be free of nuclear effects and parametrized by the different collaborations according to the DIS measurements on proton and deuterium. Using Eq. (2.7.70) we compute from the isoscalar corrected structure function, ˆ FA 2, published by the different collaborations, the structure function FA 2, which is the one employed in this work. Finally, in the case of neutrino DIS, both target-mass [172] and radiative [173] corrections to the cross-section have been included. 2.7.3 Analysis method The typical procedure to extract the optimal values of the parameters that fit the experimental data is to minimize the global χ2-function. We use for this purpose the routine MINUIT [174] from CERNLIB library. The χ2is defined as χ2({aj})≡X i [Di−Ti({aj})]2 σ2 i ,(2.7.72) where Diare the measured experimental values, Tiare the corresponding theoretical predictions, and σ2 iare the experimental errors. In fact, σ2 i. are the systematic and statistical uncertainties added in quadrature, since in most of the cases the correlation matrices are not available. The parameters {aj}are the set of parameters – 12 in our case – that define the nuclear modification at the initial scale. We define our central fit as that corresponding to the minimum value of the global χ2 obtained using our set of free parameters {aj}, χ2({a0 j})≡min χ2({aj}).(2.7.73) 52 CHAPTER 2. NUCLEAR PARTON DISTRIBUTION FUNCTIONS AT NNLO Estimation of the uncertainties The set of parameters {a0 j}represent our best estimate of nPDFs. However, due to the existence of experimental uncertainties, we can move along the neighborhood of {a0 j}having still a good agreement with data. The main goal of this subsection is to quantify the uncertainties of the nPDFs and also their propagation to any hard process quantity X. The Hessian approach [175], which will be briefly described in the following, is used for this purpose. The basic assumption of the Hessian method is that the dominant behavior of the χ2- function near the fitted minimum can be approximated by a quadratic form of the fitting parameters, {aj}, χ2({aj})≈χ2 0+X i,j Hij yiyj,(2.7.74) where yi=ai−a0 iare the parameters shifts from their best-fit values, χ2 0=χ2({a0 j}) is the value of the χ2-function in the minimum, and Hij is the Hessian matrix, defined as Hij ≡1 2∂2χ2 ∂yi∂yjai=a0 i .(2.7.75) Being a real and symmetric matrix, Hij, has a complete set of n– 12 in our case – orthonormal eigenvectors v(k) iwith eigenvalues ksuch that X j Hijv(k) j=kv(k) i. X i v(j) iv(k) i=δjk (2.7.76) We can define a new set of parameters {zj}replacing our parameters {yj}by zi≡√iX j v(i) jyj.(2.7.77) In these new coordinates, Eq. (2.7.74) reduces to ∆χ2≡χ2−χ2 0≈X i z2 i.(2.7.78) In other words, the surfaces of constant χ2are – in the quadratic approximation – hyperspheres in z-space, with distance to the minimum given by Eq. (2.7.78). Now let us consider any physical quantity Xwhich depends on the PDFs, i.e, which is a function of the parameters {aj}. Assuming that the linear term of Taylor expansion of X around its central value X0≡X(a0 i)≡X(zi= 0) gives an adequate approximation, one has ∆X≡X−X0≈X j ∂X ∂zi zi,(2.7.79) 2.7. A-Z GLOBAL ANALYSIS 53 where the z-gradient of Xis evaluated at the global minimum, that is, at the origin in z-space. Since χ2increases uniformly in all z-space directions, the vector in this space that maximizes ∆Xfor a given ∆χ2is that in the direction of the gradient of Xwith length of p∆χ2. For the square deviation we obtain (∆X)2≈∆χ2X i∂X ∂zi2 .(2.7.80) In order to compute the partial derivatives in Eq. (2.7.80), a set of 2n– 24 in our case – auxiliary PDFs, S± i, in the z-space are defined as S± 1=±p∆χ2(1,0, ...0) S± 2=±p∆χ2(0,1, ...0) . . . S± N=±p∆χ2(0,0, ...N),(2.7.81) where nis the original number of parameters aj. These S± iPDFs, usually known as error sets, are those in which the fit parameters are changed by a fixed amount in the z-space direction separately. The central set, S0= (0,0, ..., 0), is the set giving the minimum χ2. Using these PDF error sets, the derivatives in Eq. (2.7.80) can be approximated by ∂X ∂zi≈X(S+ i)−X(S− i) 2p∆χ2(2.7.82) Inserting this expression in Eq. (2.7.80) we get ∆X≈1 2sX iX(S+ i)−X(S− i)2.(2.7.83) How to choose ∆χ2? Ideally, one would expect the errors to be given by ∆χ2= 1 for one standard deviation (90 % C.L. limit). This is adequate when fitting consisting data with ideal Gaussian uncertainties to a well-defined theory. However, global analysis combine data from a wide variety of independent experiments, where there are unknown experimental and theoretical uncertainties. Therefore, the situation is far from being ideal. So, how to determine ∆χ2in global analysis? The 90 % C.L. limit can be computed for each data set kwith Nkdata points. Then, ∆χ2can be chosen to ensure that each data set is described within its 90 % C.L. limit. The 90 % C.L. limit for each data set kis defined by [176] Zξk 0 dχ2 2 Γ(Nk/2)χ2 2Nk/2−1 e−χ2/2= 0.90 ,(2.7.84) 60 CHAPTER 3. PERCOLATION OF STRINGS Figure 3.1: Overlapping discs up to the percolation phase transition. Figure taken from [52]. discs in percolation theory. A cluster of strings has an homogeneous color field which is the result of the vectorial sum of those of the original strings. At a given critical density, ηc, a macroscopic cluster appears across the collision surface. This marks the percolation phase transition. Hence, the nature of this transition is geometrical. Consider a flat large 2-dimensional surface, SA=πR2, which is the total nuclear overlapping area. On this surface, ndiscs of area S1=πr2 0are randomly distributed, allowing overlap between them. If the number of discs grows, clusters of overlapping discs start to form. The density of discs is given by ρ=n/SA=n/πR2. When this density increases, the average cluster size increases and at a certain critical density, ρc, the cluster occupies the whole surface. This is known as percolation and it is shown in Fig. 3.1. The percolation threshold, ηc, that is, the onset of continuum percolation, is related to the critical density, ρc, by ηc=ρcS1=ρcπr2 0.(3.0.2) ηchas been computed using numerical simulations for different systems [48–50]. The results are in the range ηc= 1.12 −1.5 depending on the profile function used for the colliding nuclei – homogeneous or a 3-parameter Fermi distribution. When the percolation point has been reached, i.e., when the largest cluster expands over the entire surface, there is still a considerable fraction of the surface which is empty. Indeed, at the threshold, only 1−exp(ηc)∼2/3 of the surface is covered by discs. The expressions of the average transverse momentum, hpTi, and average multiplicity, hµi, of the particles produced in a cluster can be derived from the composition of the color fields of its strings. Let us consider a cluster made of nstrings, with area Snand color charge ~ Qn. ~ Qnis the vectorial sum of the color charges of each individual string, ~ Qi, ~ Q2 n= n X 1 ~ Qi!2 .(3.0.3) Since strings colors are arbitrarily oriented, in the limit of large n, the average of ~ Qi·~ Qj 61 (with i6=j) is zero. Hence, ~ Q2 n=n~ Q2 1. Due to the Gauss theorem, Qn=rnSn S1 Q1.(3.0.4) This approximate equation will be used as a smooth interpolation between the non-overlapping and the total overlapping extreme cases. From the color charge in Eq. (3.0.4) and applying the Schwinger formula Eq. (3.0.1) the expressions for the average multiplicity and the average transverse momentum of a cluster can be obtained: hµni=rnSn S1hµ1i,hp2 Tni=rnS1 Snhp2 T1i,(3.0.5) where hµ1iand hpT1istand, respectively, for the average multiplicity and transverse momentum of the particles produced in a single string. Note that a kind of conservation law of the total transverse momentum produced, holds 1 nhµnihp2 Tni=hµ1ihp2 T1i.(3.0.6) In low energy peripheral HICs the number of nucleon-nucleon interactions is small and so is the number of strings formed. Thus, in this case, strings act as independent sources and, so, Sn=nS1. Then, from Eq. (3.0.5), the total average multiplicity and transverse momentum are: hµni=nhµ1iand hpTni=hpT1i. On the contrary, in the total overlapping limit case, Sn=S1and, then, Qn=√n Q1,hµni=√nhµ1i, and, hp2 Tni=√nhp2 T1i. To obtain the mean pTand the mean multiplicity of a collision at a given centrality, one needs to sum over all formed clusters and to average over all events, hµi=PNevents i=1 Pjhµnji Nevents ,hpTi=PNevents i=1 PjhµnjihpTnji PNevents i=1 Pjhµnji,(3.0.7) where the sum over jgoes over all individual clusters (j), each one constituted by njstrings and occupying an area Snj. These expressions show a good agreement with experimental data. For instance, as it can be seen in Fig 3.2, when clustering is included, there is a prefect agreement of the multiplicity of negatively charged particles in Pb-Pb collisions with the data from NA49 Collaboration [189]. Let us consider now the “thermodynamic limit”, i.e., the limit where the number of strings goes to infinity, n→ ∞, keeping the percolation parameter,η=ρS1, fixed. In this limit, the distribution of the overlaps is Poissonian, Pn=ηne−η/n!, and [191] hµni=nF(η)hµ1i,hp2 Tni=hp2 T1i F(η),(3.0.8) 62 CHAPTER 3. PERCOLATION OF STRINGS Figure 3.2: Mean multiplicity of negatively charged particles in Pb-Pb collisions at Plab = 158AGeV/c not including clustering (dashed line) and including it (continuous line) compared to NA49 experimental data [189]. Figure taken from [190]. where F(η) = s1−e−η η(3.0.9) is a continuous geometric saturation function whose value is close to unity at low ηand to zero at very high ones. In Eq. (3.0.9), the numerator, 1 −e−η, represents the fraction of the area covered by the strings. A more realistic implementation implies a modification of this area in Eq. (3.0.9). In HICs experimental data on the multiplicity are known to be well described by a negative binomial distribution P(µ, k) = γk Γ(k)µ! Γ(µ+k) (1 + γ)µ+k,(3.0.10) where kstands for 1 k≡hµ2i−hµi2−hµi hµi2,(3.0.11) and γ=k/hµi. P(µ, k) can be written as a convolution of the probability of having a cluster composed by nstrings, W(n), and the probability for that cluster to fragment into µparticles, G(n, µ), [192] P(µ, k) = Z∞ 0 dn W(n)G(n, µ).(3.0.12) 3.1. GEOMETRIC SCALING FOR HADRONIC INTERACTIONS 63 G(n, µ) is taken Poissonian, G(n, µ) = e−nnµ µ!,(3.0.13) and W(n) is taken as a Gamma distribution, W(n) = γ Γ(k)(γn)k−1e−γn ,(3.0.14) with kgiven by 1 k=hn2i−hni2 hni2.(3.0.15) There are several reasons for this choice. First, the gamma distribution reproduces to a good approximation the cluster size distributions at different centralities. Indeed, in peripheral collisions the density of strings is small and there are only very few overlapping strings. In this case, the cluster size is peaked at low values of the number of strings of the cluster. As the centrality increases, so does the density of strings, and there are more and more overlapping strings. The cluster size distribution becomes strongly modified. Second, there is a more technical reason related to the renormalization group [193]. R ≡ 1/k can be interpreted as the fluctuations on the number of strings in the clusters. Note that, on the one hand, in the large density limit, ηlarge, hn2i−hni2∼ hni, then k→ hni → ∞ [192,194]. On the other hand, in the low density limit, η1, the multiplicity, µ, is Poisson-like, hence, its variance coincides whith its mean value, and k→ ∞, see Eq. (3.0.11). At intermediate energy densities, kmust have a minimum close to the critical density. In the following sections, some of the results of this model and their comparison with the available experimental data will be discussed. 3.1 Geometric scaling for hadronic interactions It has been recently shown that the pT-spectra of charged particles in p-p collisions exhibit geometric scaling [195, 196]. Indeed, the pT-spectra in p-p collisions in the broad range of energies from 0.9 to 7 TeV, scale in a single variable τ≡p2 T/Qp s2, where the proton saturation momentum Qp sis given by (Qp s)2≡Q2 0W pTλ ,(3.1.16) where W=√sand λ= 0.27. An extension of this geometric scaling to A-A collisions, for both RHIC and LHC energies, at different centralities, and different nuclei is considered here. It is shown that this scaling is not only valid for each collision at fixed centrality separately, but also for any 64 CHAPTER 3. PERCOLATION OF STRINGS centrality for τ < 1. Furthermore, the hard spectrum, i.e., the spectrum for τ > 1 is also analyzed. The hard multiplicity turns out to decrease with the size of the participant nuclei – as it would be expected from jet quenching. 3.1.1 The saturation momentum The saturation momentum, QA s, in A-A collisions is defined as, QA s2= (Qp s)2Nβ(s)/2 AA1/6A NA1/3 ,(3.1.17) where β(s) = 1 3 1−1 1 + ln ps/s0+ 1 ,(3.1.18) NAbeing the number of participant nucleons divided by two and Athe mass number. At high energy, β(s) = 1/3. With this value, for NA=Athe well-known behavior of QA s for central collisions QA s2= (Qp s)2A1/3is recovered. The parametrization of Eq. (3.1.17) is based on the description of the experimental data on dN/dy for p-p and A-A collisions at all centralities, rapidities, and energies using the framework of percolation of strings [197, 198]. The dependence of the charged particle multiplicity on the center of mass energy is the same in p-p and A-A collisions, as it is shown in reference [197]. The observed differences are due to the energy conservation effect that can be incorporated in β(s). In fact, dnAA dy y=0 ∼NANβ(s) A−1dnpp dy y=0 ,(3.1.19) where β(s) is given by Eq. (3.1.18). As it was already explained, the string percolation can be regarded as a simpler implementation of the CGC [199]. Therefore, the number of color flux tubes of the glasma, (QA s)2R2 A, corresponds to the number of clusters of strings (effective number of sources), η1/2R2 A, in SPM. In this way, the dependence of the string density, η, on s,Aand NAis translated into (QA s)2, resulting in Eq. (3.1.17). The values of the parameters λ= 0.27 and Q0= 1 GeV/c are taken from the reference [195]. A slightly worse scaling has been obtained with λ= 0.30. The value of √s0= 245 GeV was obtained in references [197, 198] and indicates the energy scale of the energy conservation effect. 3.1. GEOMETRIC SCALING FOR HADRONIC INTERACTIONS 65 (QA s)2in Eq. (3.1.17) does not depend on the rapidity. Nevertheless, only the central pseudorapidity region will be used. In most of the string models, strings stretching among sea quarks and antiquarks are expanded only in the central pseudorapidity range. This fact changes the dependence of strings on NA. For central pseudorapidity, Ns∼N4/3 Aand outside this range Ns∼NA. This change gives rise to a smoother dependence of (QA s)2on NA outside the central pseudorapidity region. 3.1.2 Comparison with experimental data In Fig. 3.3, the multiplicities for Cu-Cu 0 −6 % at 62.4 and 200 GeV, and Au-Au 0 −6 % at 62.4 and 200 GeV [203,204], together with Pb-Pb 40 −50 %, 5 −10 %, and 0 −5 % at 2.76 TeV [205] as a function of τ≡p2 T/QA s2, in the pseudorapidity range 0.2< η < 1.4 are plotted. The number of participants used is the mean value corresponding to the given centrality. Experimental data taken from ALICE [205] correspond to a pseudorapidity range – including the η= 0 region –, in which dn/dη is smaller than in the pseudorapidity range considered here. Because of this a 15 % correction was applied to the normalization [198]. It can be seen that Pb-Pb 0−5 %, and 5−10 % data points lie in the same line for τ < 1. Also the Cu-Cu, and the Au-Au data at both energies are approximately in the same curve for τ < 1. Only Pb-Pb 40 −50 % present some departure around 1. For τ > 1 a suppression for the heavy nuclei is shown. In order to see the differences between the different sizes of projectile and target, Fig. 3.4 in presented. In this figure, data for p-p collisions at 0.9, 2.36, and 7 TeV [206, 207] are plotted together with the recent p-Pb data at 5.02 TeV [208]. Au-Au central data at 62.4 GeV [203], and Pb-Pb central data at 2.76 TeV [205], already shown in Fig. 3.3, are also plotted here. In the case of p-Pb, hNparti= 7.9 is used, according to reference [209]. It is observed that p-p data at different energies are all in the same line, satisfying geometric scaling, as it was shown in reference [195]. Moreover, p-Pb, Pb-Pb and Au-Au central data at low τare very close to these p-p data. When τbecomes larger the difference between these sets of data increases. For τ > 1 the suppression is larger for Pb-Pb 0 −5 % central than for p-Pb, and the latter are more suppressed than the p-p sets. The scaling for Q2< Q2 swas predicted by CGC and by phenomenological saturation models [210–213]. The HERA data on DIS at low xshow scaling even for very high Q2, Q2<400 GeV2/c2[214]. Indeed, the solution of the BFKL evolution equation shows that the scaling can be extended to intermediate Q2, 1 .ln(Q2/Q2 s)ln(Q2 s/Λ2 QCD) [215]. Our comparison with data shows that p-p collisions present geometric scaling even for τ > 1. However, this is not true for A-A collisions. Notice that in the p-p case, jet quenching is not expected, since a high density medium is not formed. On the contrary, in A-A collisions, at high LHC energies and high multiplicity events, jet quenching has been predicted [216,217]; but the weight of these events compared to minimum bias is negligible. 66 CHAPTER 3. PERCOLATION OF STRINGS 10 − 310 − 210 − 1 1.0 10 102 τ 10 − 8 10 − 6 10 − 4 10 − 2 1.0 102 N − 1 A (2 πpT ) − 1 d 2 Nch/dηdpT ( GeV/c ) − 2 0 . 2 <η < 1 . 4 Au-Au 62.4 GeV 0-6% PHOBOS Au-Au 200 GeV 0-6% PHOBOS Cu-Cu 62.4 GeV 0-6% PHOBOS Cu-Cu 200 GeV 0-6% PHOBOS Pb-Pb 2.76 TeV 0-5% ALICE Pb-Pb 2.76 TeV 5-10% ALICE Pb-Pb 2.76 TeV 40-50% ALICE Figure 3.3: Charged particle multiplicity per participant at pseudorapidity 0.2< η < 1.4 for Au-Au and Cu-Cu central collisions at two RHIC energies 62.4 and 200 GeV [203,204], and for Pb-Pb collisions at 2.76 TeV [205] versus τ. 10 − 210 − 1 1.0 10 τ 10 − 7 10 − 5 10 − 3 10 − 1 10 N − 1 A (2 πpT ) − 1 d 2 Nch/dηdpT ( GeV/c ) − 2 0 . 2 <η < 1 . 4 p-p 0.9 TeV CMS p-p 2.36 TeV CMS p-p 7 TeV CMS Au-Au 62.4 GeV 0-6% PHOBOS Pb-Pb 2.76 TeV 0-5% ALICE p-Pb 5.02 TeV ALICE Figure 3.4: Charged particle multiplicity per at pseudorapidity 0.2< η < 1.4 for p-p collisions [206, 207], Au-Au 0 −6 % central collisions at 62.4 GeV [203], Pb-Pb 0 −5 % collisions at 2.76 TeV [205], and p-Pb data at 5.02 TeV [208] versus τ. 3.1. GEOMETRIC SCALING FOR HADRONIC INTERACTIONS 67 10 − 310 − 210 − 1 1.0 10 102103 τ 10 − 11 10 − 8 10 − 5 102 10 N − 1 A (2 πpT ) − 1 d 2 Nch/dηdpT ( GeV/c ) − 2 0 . 2 <η < 1 . 4 p-p 0.9 TeV CMS p-p 2.36 TeV CMS p-p 7 TeV CMS Cu-Cu 62.4 GeV 0-6% PHOBOS Cu-Cu 200 GeV 0-6% PHOBOS Au-Au 62.4 GeV 0-6% PHOBOS Au-Au 200 GeV 0-6% PHOBOS Pb-Pb 2.76 TeV 0-5% ALICE Pb-Pb 2.76 TeV 5-10% ALICE Pb-Pb 2.76 TeV 40-50% ALICE p-Pb 5.02 TeV ALICE Figure 3.5: Charged particle multiplicity per participant in the pseudorapidity range 0.2< η < 1.4 for all the heavy ion collisions considered versus τ. In Fig. 3.5, all previous data were plotted together to show that for τ < 1 there are not many differences among them. This does not occur for τ > 1 region, where suppression is larger for larger sizes of nuclei. It is worth emphasizing that the parametrization in Eq. (3.1.17) allows us to deal with all type of collisions at any energy. For that purpose, the introduction of the β(s) function is required, which implies a different power for different number of participants. In reference [196] the geometric scaling was shown for p-p collisions at LHC and A-A collisions at RHIC using an energy independent exponent λ, but this exponent for p-p was different that the one from A-A. If an energy independent exponent is used for both, p-p and A-A collisions, the scaling is spoiled as it can be seen in Fig. 3.6, where p-p and several A-A transverse momentum distributions using β(s) = 1/3 are plotted. In order to see the quality of this extended scaling, the ratio of Pb-Pb 0 −5 % at 2.76 TeV, Au-Au 0−6 % at 200 GeV, Cu-Cu 0−6 % at 62.4 GeV, p-Pb at 5.02 TeV over Cu-Cu 0−6 % at 200 GeV as a function of τis shown in Fig. 3.7. Even for such different projectiles and targets and also in the broad range of energies considered, an approximate scaling at low τis observed. Notice that for τ < 1 the data extend over three orders of magnitude. Excluding Cu-Cu data, ratios vary between 0.7 and 1.3 in the range 0.2< τ < 11. The violation of the scaling is clear for τ > 1, showing a higher suppression for heavier nuclei. The hierarchy of the scaling violations for τ > 1 agree with the expected suppression due to 1Note that for pT<ΛQCD or, equivalently, τ < 0.1−0.2, there is no reason to expect geometric scaling. 68 CHAPTER 3. PERCOLATION OF STRINGS 10 − 210 − 1 1.0 10 τ 10 − 7 10 − 5 10 − 3 10 − 1 10 N − 1 A (2 πpT ) − 1 d 2 Nch/dηdpT ( GeV/c ) − 2 0 . 2 <η < 1 . 4 p-p 0.9 TeV CMS p-p 2.36 TeV CMS p-p 7 TeV CMS Au-Au 62.4 GeV 0-6% PHOBOS Cu-Cu 200 GeV 0-6% PHOBOS Pb-Pb 2.76 TeV 40-50% ALICE Figure 3.6: Charged particle multiplicity per participant in the pseudorapidity range 0.2< η < 1.4 for p-p collisions [206,207], Au-Au 0 −6 % at 62.4 GeV [203], Cu-Cu 0 −6 % at 200 GeV [204], and Pb-Pb 40 −50 % at 2.76 TeV [205] versus τusing β(s) = 1/3. jet quenching. Based on this approximate scaling, the multiplicity of soft and hard particles per participant, defining soft particles as those with pT< Qsand hard as those with pT> Qs, can be computed. In fact, 1 NA dNsoft ch dη =1 NAZQ2 s 0 dp2 T dN2 ch dηdp2 T =1 NAZQ2 s 0 dp2 T 1 Q2 0 F(τ),(3.1.20) and with the following change of variables dp2 T Q2 0 =2 2 + λW Q02λ 2+λ τ−λ 2+λNβ(s)/2 AA1/6A NA1/3 dτ , (3.1.21) the fraction of soft and hard multiplicities over the total multiplicity results Rs≡dNsoft ch /dη dNtot ch /dη =R1 0dτ τ−λ 2+λF(τ) R∞ 0dτ τ−λ 2+λF(τ), Rh≡dNhard ch /dη dNtot ch /dη =R∞ 1dτ τ−λ 2+λF(τ) R∞ 0dτ τ−λ 2+λF(τ). (3.1.22) In Table 3.1, the fractions of soft and hard multiplicities for the centralities, energies, and collisions considered are presented. In order to do the integration in Eq. (3.1.22), first a fit to each pT-distribution separately is done, and, afterwards, this is integrated. 3.1. GEOMETRIC SCALING FOR HADRONIC INTERACTIONS 69 Figure 3.7: Ratio of Pb-Pb 0 −5%at2.76 TeV [205], Au-Au 0 −6 % at 200 GeV [203], Cu-Cu 0 −6 % at 62.4 GeV [204], and p-Pb at 5.02 TeV [208] w.r.t. Cu-Cu 0 −6 % at 200 GeV [204] versus τ. The hard fraction decreases slowly with center of mass energy and with the size of the participant nuclei, as it can be seen in Table 3.1. In the case of p-p collisions the dependence of the energy is very weak, varying from 9 % to 8 % in the broad range 0.9−7 TeV. However, larger differences between Au-Au central collisions at 62.4 GeV and Pb-Pb central collisions at 2.76 TeV are observed: their hard fractions are, respectively, 7 % and 2 %. The hard fraction of Pb-Pb is the same for peripheral (40−50 %) and central (0−5 %) collisions. In summary, the weight of the hard collisions – and, hence, the weight of the hard multiplicities – is believed to increase when the energy or the size of the participant nuclei grow; however, the opposite is observed here, within our definition of hard and soft collisions. This result was also pointed out by ALICE Collaboration data on the sphericity as a function of the energy and the charge multiplicity in p-p collisions at 0.9, 2.36, and 7 TeV [218]. The sphericity measures the jet activity, in such a way that one event with dijet back-to-back implies sphericity zero. On the contrary, an event where all the produced particles are distributed isotropically in phase space implies sphericity one. Most of the Monte Carlo codes (Pythia 8, Perugia-0, Phojet, Atlas-CSC) predict a decreasing of the sphericity with energy and charged multiplicity, whereas the data show the opposite trend. In the CGC [15,16,219] the more the saturation momentum grows with the energy and the number of participants, the smaller is the room for hard collisions. This also happens in the framework of percolation 76 CHAPTER 3. PERCOLATION OF STRINGS 0 50 100 150 200 250 300 350 400 450 Nch 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 A 1 Au-Au 200 GeV STAR Collaboration Au-Au 62 GeV STAR Collaboration Figure 3.9: Comparison between our results on the strength of the near-side ridge for Au-Au collisions at two RHIC energies √s= 200 GeV (red line) and √s= 62 GeV (blue line) with experimental data [246] versus Nch. the nucleus radius should be similar to 1/R, as it is obtained. The pseudorapidity width, obtained from Eq. (3.2.26), compared with the experimental data on Au-Au at √s= 200 GeV and √s= 62 GeV [246] is shown in Fig. 3.11. The value of c1is 0.23. NAand NCare taken from the quoted experimental analysis [256]. It is observed that our result for √s= 200 GeV is slightly larger than the corresponding one at √s= 62 GeV. Experimental data are very close at both energies. In Fig. 3.12 experimental data on the azimuthal width for Au-Au collisions at √s= 200 GeV and √s= 62 GeV [246] are compared to our model. The values of c2are: c2= 0.866 for Au-Au at √s= 200 GeV, and c2= 0.890 for Au-Au at √s= 62 GeV. Increasing energy and centrality, the azimuthal width decreases, in agreement with the trend of the experimental data. For both widths a qualitative agreement is obtained. The dependence of A1,σ∆η, and σ∆φon energy and centrality, resulting from equations (3.2.29) and (3.2.30), is very similar to the one obtained in the glasma picture. In this approach, RdN/dy and σ∆ηare proportional to 1/αs(Qs) and σ∆φis proportional to 1/Qs. Hence, both RdN/dy and σ∆ηgrow with ln sand ln NA. In the high density limit, the same dependence on the center of mass energy and on the number of participants for both observables is obtained in percolation . In the case of σ∆φ, as 1/Qs∼N1/6 As∆/2and r0F(η)1/2∼r0η1/4∼r0N1/6 As∆/2. 3.2. THE NEAR-SIDE RIDGE STRUCTURE 77 0 50 100 150 200 250 300 350 Nch 0.00 0.05 0.10 0.15 0.20 A 1 p-Pb 5.02 TeV CMS Collaboration pp 7 TeV CMS Collaboration Figure 3.10: Comparison between our results on the strength of the near-side ridge for p-p collisions at √s= 7 TeV (blue line) and p-Pb collisions √s= 5.02 TeV (red line) with experimental data [229] versus Nch. 1.0 1.5 2.0 2.5 η 0.5 1.0 1.5 2.0 2.5 σ ∆ η Au-Au 62 GeV STAR Collaboration Au-Au 200 GeV STAR Collaboration Figure 3.11: Pseudorapidity width of the near-side ridge for Au-Au collisions at two RHIC energies √s= 200 GeV and √s= 62 GeV [246] versus η. Curves obtained from Eq. (3.2.26) for Au-Au collisions at √s= 200 GeV (red) and √s= 62 GeV (blue). 78 CHAPTER 3. PERCOLATION OF STRINGS 1.0 1.5 2.0 2.5 3.0 η 0.6 0.7 0.8 0.9 1.0 1.1 σ ∆ φ Au-Au 62GeV STAR Collaboration Au-Au 200 GeV STAR Collaboration Figure 3.12: Azimuthal width of the near-side ridge for Au-Au collisions at √s= 200 GeV and √s= 62 GeV [246] versus η. Curves obtained from Eq. (3.2.32) for Au-Au collisions at √s= 200 GeV (red) and √s= 62 GeV (blue). To conclude, it has been shown that percolation of strings may naturally explain the anomalous dependence of the near-side ridge structure correlation on the multiplicity observed in Au-Au collisions at √s= 200 GeV and √s= 62 GeV. The onset of the ridge structure in high multiplicity p-p events at √s= 7 TeV and in high multiplicity p-Pb at √s= 5.02 TeV is also explained. Furthermore, our model qualitatively describes the dependence of the azimuthal and pseudorapidity widths on multiplicity. Our framework can be regarded as a complementary picture to the glasma in the description of the initial state, able to explore the transition from low to high density. Most of the ingredients used here can be considered initial state effects, although the quenching of the partons produced in the cluster of strings is a final state effect. 3.3 Geometric scaling of elliptic flow The discovery of a sizable elliptic flow in A-A collisions, first observed at RHIC [7, 8] and later at the LHC [257], turned up as an experimental major breakthrough. The observed anisotropic flow can be exclusively understood if the measured particles in the final state depend not only on the physical conditions realized locally at their production point, but also on the global geometry of the event. This non-local information can solely emerge as a collective effect, requiring strong interactions among the relevant degrees of freedom, i.e. quarks 3.3. GEOMETRIC SCALING OF ELLIPTIC FLOW 79 and gluons. The study of higher harmonics has also shown very interesting features, including the ridge structure seen in A-A collisions [221, 222, 225, 226], p-Pb collisions [229, 230] and also in high multiplicity p-p collisions [227], which was analyzed in Section 3.2. Along these lines, it is pointed out that some scaling laws satisfied by the elliptic flow could be very useful to determine some properties of the initial stage of the collision which should be preserved by the hydrodynamic evolution [258]. The experimental data on the elliptic flow of charged particles showed up a universal scaling law related to the gluon saturation momentum. This scaling law is also satisfied by the photon data, suggesting that the elliptic flow of charged particles and photons should have a common origin. 3.3.1 Universal scaling law The experimental data for v2at RHIC and at the LHC normalized to the saturation momentum, eccentricity, and radius of the collision area satisfy geometric scaling: v2(pT) 1QA sL=f(τ),(3.3.33) where 1=2 πZπ/2 0 dϕ cos 2ϕR2−R2 ϕ R2, Rϕ=RAsin(ϕ−α) sin ϕ,(3.3.34) α= arcsin( b 2RA sin ϕ), R2=hR2 ϕi=2 πZπ/2 0 dϕ R2 ϕ,(3.3.35) and τ=p2 T (QA s)2,(3.3.36) being QA sthe saturation momentum, RAthe radius of the nucleus, and Lthe length associated to the size of the collision area at a given impact parameter and energy. Indeed, the product QA sLis the inverse of the Knudsen number, i.e., the mean free path normalized to the length measured as the number of scattering centers. The scaling law (3.3.33) is tested in the range 0 < τ < 1. 1is a measurement of the eccentricity of the collision. It does not depend on the distribution of scattering centers – partons or nucleons in the transverse plane. It is only determined by the almond shape of the collision at a given impact parameter. The scaling variable τis known from the geometric scaling verified in DIS, p-p, p-A, and A-A collisions [195,196,259–261]. That is, τ=p2 T/QA s2, where, now, QA s2= (Qp s)2Aβ(s)/2N1/6 A,(3.3.37) being NAthe number of wounded nucleons. β(s) and the proton saturation momentum, Qp s, are given, respectively, by equations (3.1.18) and (3.1.16)4. 4Note that the definition of QA shas been slightly modified from the one given in 3.1.17. 80 CHAPTER 3. PERCOLATION OF STRINGS 3.3.2 Discussion In Fig. 3.13 (a), v2(pT) for Au-Au collisions for different centralities at RHIC [262], and for Pb-Pb collisions at LHC [257] divided by the product 1QA sLcomputed for each centrality and energy is plotted. The usual values of band NAfor each centrality are taken, in order to compute 1and QA susing Eqs. (3.3.34) and (3.3.35), and Eq. (3.3.37), respectively. The length, L, measures the number of longitudinal scatterings, which in the Glauber model is proportional to N1/3 A. Nevertheless, (1 + N1/3 A)/2 is employed, as in most of the strings models, DPM [247,248], quark gluon string model [249], Venus [185], and EPOS [250]. The values of b,NAand 1for each centrality and energy are shown in Table 3.3. √s200 GeV (PHENIX) 2.76 TeV (ALICE) Centrality 10-20 % 20-30 % 30-40 % 40-50 % 10-20 % 20-30 % 30-40 % 40-50 % b(fm) 5.7 7.4 8.7 9.9 5.6 7.4 8.9 10.1 NA117.3 83.3 57.1 37.2 130.05 92.9 64.25 42.35 10.208 0.286 0.356 0.436 0.172 0.238 0.300 0.357 Table 3.3: Values of the impact parameter, NA=Npart/2 and 1for PHENIX [256] and ALICE [263] at different centralities. The solid black line in this figures corresponds to a fit to these data, given by v2 1QA sL=aτb,(3.3.38) where a= 0.1264±0.0076 and b= 0.404±0.025. Fig 3.13 shows that this scaling is satisfied. In order to see the quality of this scaling, the ratio of Pb-Pb 10 −20 % at 2.76 TeV, Pb-Pb 40−50 % at 2.76 TeV, Au-Au 20−30 % at 200 GeV and Au-Au 30−40 % at 200 GeV over Pb-Pb 30−40 % at 2.76 TeV as a function of τis shown in Fig. 3.13 (b) . All the ratios lie in the range 0.8−1.15 for the whole τconsidered, showing that the scaling is quite good 5. The experimental data used in Fig. 3.13 correspond to event plane [262] and 4-particle cumulant measurements and, therefore, include some amount of fluctuations. As in the scaling law of Eq. (3.3.33) the quantities ε1,QA s, and Lhave nothing to do with fluctuations, the latter could give rise to the residual differences between the experimental data and the function f(τ) of Eq. (3.3.33). Changing the eccentricity, 1, by the usual eccentricity, =hy2−x2i/hy2+x2i, or by the participant eccentricity, the scaling is not satisfied for both Monte-Carlo Glauber and Color Glass distributions. This fact does not mean that the initial state should give the 5Most of the experimental error data are of the order of 10 %. 3.3. GEOMETRIC SCALING OF ELLIPTIC FLOW 81 0.0 0.2 0.4 0.6 0.8 1.0 τ 0.00 0.02 0.04 0.06 0.08 0.10 0.12 0.14 2 v 2 / ( ² 1 Qs (1 + N 1 / 3 A )) (a) fit to data Pb-Pb 2.76 TeV ALICE 10-20% Pb-Pb 2.76 TeV ALICE 20-30% Pb-Pb 2.76 TeV ALICE 30-40% Pb-Pb 2.76 TeV ALICE 40-50% Au-Au 200 GeV PHENIX 10-20% Au-Au 200 GeV PHENIX 20-30% Au-Au 200 GeV PHENIX 30-40% Au-Au 200 GeV PHENIX 40-50% 0.0 0.2 0.4 0.6 0.8 τ 0.7 0.8 0.9 1.0 1.1 1.2 R (b) Au-Au 200 GeV PHENIX 20-30% Au-Au 200 GeV PHENIX 30-40% Ratio Pb-Pb 2.76 TeV ALICE 10-20% Ratio Pb-Pb 2.76 TeV ALICE 40-50% Figure 3.13: (a) v2divided by the product 1QA sLfor 10 −20 %, 20 −30 %, 30 −40 % and 40 −50 % Au-Au collisions at 200 GeV [262], and for 10 −20 %, 20 −30 %, 30 −40 % and 40 −50 % Pb-Pb collisions at 2.76 TeV [257] versus τ. The solid black line is a fit to data according to Eq. (3.3.38). (b) Ratio of Pb-Pb 10−20 %, Pb-Pb 40−50 % at 2.76 TeV [257] , Au-Au 20 −30 %, and Au-Au 30 −40 % at 200 GeV [262] over Pb-Pb 30 −40 % at 2.76 TeV [257] versus τ. corresponding eccentricity of a hard profile, such as it is defined in Eq. (3.3.34). Probably, the scaling law could be preserved using other eccentricities, but in this case some changes in the dependence of Land NAare necessary. 82 CHAPTER 3. PERCOLATION OF STRINGS Figure 3.14: (a) v2of π,kand pdivided by the product 1QA sLfor 10 −20 %, 20 −30 %, 30 −40 % and 40 −50% Au-Au collisions at 200 GeV [269], and for 10 −20 %, 20 −30 %, 30 −40 % and 40 −50 % Pb-Pb collisions at 2.76 TeV [271] versus τ, for τ < 1. The solid line is a fit to Eq. (3.3.39). (b) Ratio of the experimental points over the fitting function Eq. (3.3.39) versus τ, for 0.1< τ < 1. The scaling law is also satisfied for specified particles, such as pions, kaons, and protons, as it can be seen in Fig. 3.14. In the proton case, an effective transverse momentum, Q0A s, instead of QA s,Q0A s2=N0.045 sQA s2, where Nsis the number of strings, has been used. It is known that for central collisions the ratio baryon/meson increases with pTup to a moderate value of transverse momentum. In central collisions, due to the strong color field formed, the color flux tubes in the glasma picture or the cluster of strings in the SPM have a larger string tension, producing high mass particles more efficiently. In addition to that, inside a cluster of many strings the flavor of each single string recombines with the flavor 3.3. GEOMETRIC SCALING OF ELLIPTIC FLOW 83 0.0 0.2 0.4 0.6 0.8 1.0 τ 0.05 0.00 0.05 0.10 0.15 2 v 2 / ( ² 1 Qs (1 + N 1 / 3 A )) (b) fit to data Direct photons for Pb-Pb 2.76 TeV ALICE 0-40% Direct photons for Au-Au 200 GeV PHENIX 0-20%. RXN (up) and BBC (down) Direct photons for Au-Au 200 GeV PHENIX 20-40%. RXN (up) and BBC (down) Figure 3.15: v2divided by the product 1QA sLfor direct photons at 0 −20 % and 20 −40 % Au-Au collisions at 200 GeV [269] and direct photons at 0 −40 % Pb-Pb collisions at 2.76 TeV [270] versus τ. The solid black line is the scaling curve of Fig. (3.13). All the data with pT< Qshave been included. of other strings, producing also baryons more efficiently. Concerning the pT-distributions, these two effects can be taken into account in an effective way, defining a Q0A sfor baryon production, related to QA sby a factor which was obtained in [264]. This was done by fitting the dependence of the ratio of the pT-integrated nucleon distribution over the pT-integrated pion distribution on the number of collisions, and, so, on the number of strings, Ns. For this study data from PHENIX are used [265], resulting in a N0.09 sdependence which, in terms of Q2 s, is N0.045 s[199]. It is observed in Fig. 3.14 that all the data lie, approximately, in the same curve, parametrized by v2 1QA sL=τ a+bτ +c√τ,(3.3.39) with a= 0.573 ±0.011, b= 4.76 ±0.23 and c= 1.52 ±0.34. Although in this case the scaling is not as accurate as the obtained for charged particles, a good agreement for τ > 0.1 is obtained. The discrepancies are not higher than 20 % for most of the experimental points. For τ < 0.1 a great departure occurs, probably motivated by the precision of our fit in this low τregion and the proximity with ΛQCD. Points in this region are not shown in Fig. 3.14 (b) in order to keep a good visibility of the remaining points to evaluate the scaling law. Data on v2for p-Pb collisions have not been included in these analysis due to the uncer- 84 CHAPTER 3. PERCOLATION OF STRINGS tainties in the values of NAat a given impact parameter. In addition, since direct photon production satisfies geometric scaling [261], its elliptic flow may be of the same size and pTshape of the rest of particles. In order to check this point, in Fig. 3.15, ALICE data [270] and the PHENIX data [269] at different centralities are plotted. PHENIX Collaboration quote two different points at the same pTand centrality obtained by different analysis methods (BBC and RXN detectors). In any case, data are close to the scaling curve. The scaling law of elliptic flow cannot be derived from the geometric scaling of the multiplicities in a simple way [272]. Thus, it would be interesting to know the origin of the v2-scaling law and the role played by saturation on it. Moreover, it would be also worth looking at the possibility of a scaling law similar to Eq. (3.3.33) for higher harmonics. These will be the subjects of the next section. In summary, it has been shown that experimental data on the elliptic flow of charged particles for Au-Au and Pb-Pb collisions for different centralities at RHIC and LHC energies satisfy a scaling law. The elliptic flow for identified particles, π,k, and p, lies in the same curve. The photon data, despite their large uncertainties, also verify this scaling. Other than the eccentricity, this scaling law involves the number of scatterings and a function which depends only on p2 T/QA s2. The number of scatterings is the only involved quantity related to final state effects, since the rest of variables have to do with the geometry and the gluon saturation. 3.4 Energy loss as the origin of the scaling law of v2 In the previous section, Section 3.3, a scaling law for the elliptic flow was obtained, see Eq. (3.3.33). Nonetheless, this scaling could not be derived from the geometric scaling of the transverse momentum distributions described in Section 3.1. The origin of this scaling is the subject of this section. It will be shown here that the interaction of the partons produced in the collision with the color field gives rise to this scaling. Furthermore, a detailed functional form of the scaling which shows a very good agreement with data is obtained. 3.4.1 Energy Loss To describe the dynamics of particle production at high energies in the soft domain the model of color strings with fusion and percolation [17] is used. As it was already explained, string decays are assumed to follow the Schwinger mechanism of pair production in the strong external field. The momentum distribution of these initial partons is azimuthally isotropic, P(p0) = Ce−p2 0/σ ,(3.4.40) 3.4. ENERGY LOSS AS THE ORIGIN OF THE SCALING LAW OF V285 where p0is the initial transverse momentum, σis the string tension, and Cthe normalization factor. It is important to stress again that p0is different from the observed particle transverse momentum, pT, since the parton has to traverse the cluster area emitting gluons on its way out. Hence, in fact, the momentum distribution of the observed particles has the following form P(p, φ) = C e−p2 0(p,l(φ))/σ ,(3.4.41) where φis the azimuthal angle and l(φ) is the path length inside the nuclear overlap through which the observed particle has passed before being detected. Note that due to string tension fluctuations, the distribution in Eq. (3.4.40) is transformed into the thermal one P(p0) = Ce−p0/T ,(3.4.42) where the temperature is T=pσ/2 [273, 274]. In this calculation, the latter thermal distribution is employed. Radiative energy loss has been extensively studied in the framework of perturbative QCD for a parton traversing the QGP as a result of multiple collisions with the medium scattering centers [275–279] and it will be explained in detail and used in a phenomenological analysis of single-inclusive production in the following chapter. In this case, the physical picture is different: the created parton, with a relatively small transverse momentum, moves in the external gluon field of the string or cluster of strings, which, approximately, can be taken as constant and orthogonal to the direction of the parton. This precludes from application of the perturbative QCD. In the same vein as the mechanism of pair creation, one may assume that the reaction force due to radiation is similar to the QED case, where a charged particle is moving in an external electromagnetic field E. For an ultra-relativistic particle in a very strong field, this force causes an energy loss given by [280] dp(l) dl =−012 e2(eEp(l))2/3,(3.4.43) By integrating this equation, the dependence p(l) on the path length, l, is found. This leads to our quenching formula, p0(p, l(φ)) = p1 + κp−1/3T4/3l(φ)3,(3.4.44) where the string tension, eE/π =σ, has been identified by comparison with Eq. (3.4.40) and the Schwinger equation for pair creation in a constant electric field. It has also been introduced the dimensionless quenching coefficient κ, taking the role of 0.12e2. The validity of using the QED equation for the QCD case may raise certain doubts. However, it was found using the AdS/CFT correspondence that for the N= 4 SUSY Yang-Mills theory the energy loss of a colored charge moving in the external chromodynamic field is essentially given by the same equation as in QED case (see reference [281]). 92 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS Assuming that the energy of the propagating parton is much larger that the energy of the plasma constituents, E, E0k, Bj¨orken obtained dE coll dx =−2π α2 s2 3±11 + nf 6ln 2hkiE M2∝ −α2 sT2ln E . (4.1.3) Here ±stands for an effective color charge of quarks an gluons, respectively. The only difference between the energy loss of quarks and gluons is the prefactor. Note that this expression is obtained in the fixed coupling approximation. This was later improved, including running coupling, finite energy kinematics, and quarkmass effects by various authors [291–293]. Implementing the running coupling, for E M2/T, where Mis the mass of the heavy quark, the collisional energy loss is given by [292]: •Light quark, gluon: dE coll dx q,g =−CR 4αs(ET)m2 Dln ET m2 D,(4.1.4) •Heavy quark: dE coll dx Q =dE coll dx q +2π CR 9T2αs(M2)αs(ET) ln ET M2 .(4.1.5) Here mDis the Debye mass, the inverse of the screening length, m2 D= 4π αsT21 + nf 6 [294] and CR= 4/3 (3) is the quark (gluon) color factor. As a numerical example, taking E= 20 GeV, M= 1.3 GeV (charm quark) and a medium with T= 0.4 GeV and mD= 1 GeV , the elastic energy losses per unit-length are dE coll/dx|q=−2.3 GeV/fm and dE coll/dx|Q= −2.6 GeV/fm. 4.1.2 Radiative energy loss The dominant mechanism of energy loss of a fast parton in a QCD environment is radiative energy loss. The hard parton traversing the QGP suffers multiple scatterings with it, inducing extra gluon radiation with respect to vacuum – medium-induced gluon radiation. Let us start by thin media, Lλ, being Lthe length of the media and λthe mean free path. In this case, the propagating particle suffers at most one single scattering and the QCD radiation spectrum is just given by the Bethe-Heitler (BH) bremsstrahlung expression, obtained in [295] 1. The independent Bethe-Heitler gluon spectrum is ωdIrad/dω ∝ω−1L2, where ωis the energy of the radiated gluon. When the medium is thick, Lλ,coherence is a relevant effect that needs to be taken into account. When a high-energy parton traverses a medium, coherence effects between 1The BH formula was first obtained for QED bremsstrahlung in [296] 4.1. JET QUENCHING 93 emitter and emitted quanta due to successive scattering centers with the medium show up. This leads to a destructive interference of gluon radiation with respect to the incoherent sum of scattering centers [297]. This effect was first derived in QED and it is known as LPM effect (Landau-Pomeranchuk-Migdal) [298,299]. Coherence becomes important when the formation time of the radiated gluon, τf, is much larger than the mean free path, τfλ. The formation time is the time needed by the radiated gluon to decohere from the projectile. This time can be computed using the uncertainty principle, τf=ω k2 ⊥ ,(4.1.6) where ωand k⊥are, respectively, the energy and the transverse momentum of the gluon. The decoherence between the projectile and the gluon is achieved when the gluon becomes on-shell (k2= 0). In other words, when the phase, ϕ, – acquired through multiple scatterings by the gluon – is close to the unity [300]: ϕ=k2 ⊥ 2ω∆l,(4.1.7) where ∆lis the distance travelled by the gluon during its formation time. The transport coefficient, which describes the medium-induced transverse momentum squared, k2 ⊥, that the medium transfers to the parton per mean free path λ, can be written (for a static medium) as: ˆq=k2 ⊥ λ.(4.1.8) The limiting energy of the gluons that decohere from its parent parton is called the characteristic gluon frequency,ωc. For a hard parton propagating through a finite path length Lin the medium, ωccan be easily obtained from Eq. (4.1.7) ωc=1 2ˆqL2.(4.1.9) In medium, a coherence time,τdec, can be also defined as the kicks that the projectile receives from the medium until its decoherence from the parent parton, k2 ⊥= ˆqτdec .(4.1.10) Using Eq. (4.1.6) and Eq. (4.1.10) and assuming that the coherence and the formation times are approximately the same, the coherence time can be obtained: τdec =rω ˆq.(4.1.11) The gluon energy spectrum per unit length can be estimated qualitatively by taking just the coherent sum of all the scatterings, ωd2Irad dzdω ∝αsrˆq ω,(4.1.12) 94 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS for ω < ωc. Note, that due to the destructive interference, the LPM spectrum – proportional to ω−1/2– is suppressed in the infrared, i.e. for small ω’s, compared to the independent Bethe-Heitler gluon spectrum, which was proportional to ω−1. The total energy loss can be obtained by integrating Eq. (4.1.12), ∆Erad =ZL 0 dz Zωc 0 dω ω d2I dzdω ∝αsωc∝L2.(4.1.13) The QCD energy loss shows a characteristic L2-dependence which is also present in the case of Abelian (QED) plasmas, see, for instance, [301]. This is a general feature of the medium-induced energy loss of any in-medium newborn particle. The total collisional energy loss can be obtained, in first approximation, by integrating Eq. (4.1.3), which gives ∆Ecoll ∝Lln(L). Therefore, for an extended medium, the radiative energy loss is dominant with respect to the collisional (elastic) energy loss. A numerical illustration of this, to be compared with the elastic losses of O(2 GeV/fm), would be to consider a gluon with E= 20 GeV in a medium with ˆq= 2 GeV2/fm and L= 6 fm, then dE rad/dx is O(10 GeV/fm). Gluon emission off a heavy quark differs from that off a light quark, already in vacuum. Due to kinematics constraints, the radiation is suppressed at angles smaller than M/E, being Mthe mass of the heavy quark. This is known as the dead cone effect [302] and it results into a suppression of the total gluon radiation emitted off heavy quarks. In medium, this reduction is non-trivial [301] and the corrections are O(M/E) [303]. In summary, radiative energy loss is the dominant process of energy loss in the QGP. However, it has its limitations. Due to the dead cone effect, it predicts different amounts of energy loss for light and heavy quarks [302,303]. 4.1.3 Jet quenching models The energy-loss expressions presented in the previous subsections refer to a static and uniform QGP traversed by an infinite-energy parton. However, the real situation is much more complicated. For instance, the temperature of the QGP, and therefore its Debye mass and transport coefficient, is position-dependent; the QGP evolves dynamically and, consequently, the medium properties (ˆqand mD) are also time dependent; and energy loss fluctuations show up due to the finite size of the medium. All these effects may produce significant deviations from the above-mentioned analytical energy loss formulas. Thus, more general approaches of jet quenching which compute parton energy loss within pQCD (regardless whether the properties of the medium itself can be treated perturbatively) have been developed. These main frameworks are: BDMPS-Z [275–277,307, 308], ASW [278,300,303,309], DGLV [279,304,310,311], AMY [312–314] and Higher-Twist (HT) [315–318]. All these models are based on factorization. They all compute the radiative branching of a hard parton propagating through a colored dense medium. In order to do so, they all assume that hadronization takes place in vacuum and that only a global energy loss affects 4.1. JET QUENCHING 95 the production of the fragmenting particles. Moreover, some other approximations are made in these models: •The eikonal approximation: The energy of the radiated quanta, ω, and the energy of the emitter, E, are much larger than the transverse momentum exchanged with the medium, q⊥, i.e., ω, E q⊥. In addition, in general, only soft gluons are taken into account: ωE. •The transverse momentum of the radiated parton, k⊥, is much smaller than its energy, ω:k⊥ω. That is, the gluon radiation is collinear with respect to the parent parton. •The scattering in the medium occurs locally. This means that the mean free path of the propagating particle, λ, is much larger than the Debye screening length, m−1 D: λm−1 D. These analysis differ in their assumptions about the relationships between the relevant scales (emitter energy E, virtuality Q2, typical momentum µ≈mD, and extent Lof the medium), as well as by how they approximate the space-time profile of the medium. A brief description of each one is given here. •BDMPS-Z/ASW. The BDMPS-Z model was the first to be developed. It uses the path integral formalism to calculate gluon radiation off an emitter. Energy loss in a colored medium is usually computed in a multiple soft-scattering approximation. A hard parton traversing the medium interacts with various scattering centers and splits into an outgoing parton as well as a radiated gluon. The ASW framework follows the same formalism, but it includes the interference effects between vacuum and mediuminduced radiation. Both approaches define a density distribution of scattering centers, ρ, which depends on the time (to account for the characteristic dilution of a medium in expansion). Collisional energy loss is not considered. There are two possible approximations of these models: multiple soft scattering approximation (also called dipole or saddle point approximation), used by both formalisms, and the opacity expansion, only employed by ASW [309]. In the first one, the medium is described by only one parameter, the transport coefficient, ˆq. On the other hand, the medium is characterized in the second approximation by two parameters: the density of scattering centers, n, (or the mean free path λ) and the Debye screening mass, mD. An expansion in terms of the number of scattering centers is made. Usually, this series is truncated at N= 1, which corresponds to a single hard scattering. This second approximation results in exactly the same gluon spectra as the DGLV approach (see next). Multiple gluon emission is computed (in both approximations) by an incoherent sum of single emission giving rise to a Poisson distribution. This distribution is usually called Quenching Weights (QW) [300] and will be used in the present thesis to perform a study of the experimental data. 96 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS •DGLV. In the DGLV approach the medium consists of, as in the BDMPS-Z/ASW approach, almost static (heavy) scattering centers producing a screened Coulomb (Yukawa) potential. DGLV considers the single-hard radiation spectrum (opacity expansion). Despite most of the calculations only include the leading term (N= 1), the behavior of gluon radiation at larger opacities were explored [279, 319]. Independent gluon emission is also considered to compute multiple gluon emission in a fashion similar to the QW. •AMY describes parton energy loss through field theory considering a weakly-coupled QGP in equilibrium. The medium properties are encoded in the temperature, T, and the chemical potential, µB, and a hierarchy TgT g2Tis assumed. The hard parton scatters off other partons in the medium, leading to momentum transfers of O(gT) and inducing collinear radiation. Rate equations between quarks and gluons are used to evolve the branching of the leading parton. This procedure allows to evolve the emission probability distribution as the leading parton loses its energy. It does not take into account neither vacuum radiation nor vacuum-medium interference. •Higher-Twist describes the multiple scattering of a parton as power corrections to the leading-twist cross section. The medium properties are given by the higher-twist matrix element, that factorizes from nPDFs. The evolution of the parton branching is directed by the DGLAP equations [22–24] considering in this way the interference between vacuum and medium-induced radiation. These corrections are enhanced by the medium length Land suppressed by the power of the hard scale Q2. Therefore, it is more appropriate for thin than thick media. Originally, this approach computed the medium corrections to nuclear DIS. For further details and for a direct comparison among these formalisms, see reference [320]. Monte Carlo codes to simulate modification of the jet shower and, consequently, jet quenching, are also available. Some of these generators are: PYQUEN [321], which includes both collisional and radiative energy loss, but it does not modify the whole branching process. JEWEL [322,323], which implements elastic scattering in DGLAP evolution plus radiative energy loss through a multiplicative constant in the infrared part of the splitting functions. Q-PYTHIA [324] that includes radiative energy loss of the ASW type. Elastic energy loss is neglected. And MARTINI [325], where the evolution is based on AMY rates. These all use PYTHIA [326] (or its more recent versions) as high energy generator of p-p collisions. During the last five years, several extensions of the standard computations of energy loss have been developed. The q¯q-antenna set-up [327–332] has been derived to account for interference between the different parton emitters. A new formalism based on soft collinear effective theory (SCET) that computes radiative energy loss has also started to be developed [333,334]. 4.2. PARTICLE OBSERVABLES 97 4.2 Particle observables 4.2.1 Nuclear modification factor Hard Probes are a successful tool to study the properties of the QGP. The simplest consequence of jet quenching in HICs is the suppression of the single inclusive high-pThadron spectrum relative to that in proton-proton collisions. The observable used to quantify this is the nuclear modification factor, RAA, which characterizes how the number of hadrons, h, produced in the collision of two nuclei A-A varies with respect to the equivalent number of proton-proton, p-p, collisions. It is given by the ratio of the hadron spectrum in A-A, over, the same quantity in p-p, but scaled by the number of binary nucleon collisions, hNcolli: RAA =dNAA/d2pTdy hNcollidNpp/dp2 Tdy .(4.2.14) This observable has information both from initial and final state, that is why nPDFs should be precisely determined. In the absence of nuclear modifications – of initial and final state – this observable would be equal to unity. By controlling nPDFs from pA collisions, where final state effects should be absent, we can disentangle both effects. Figure 4.1: RAA versus pTmeasured in central Au-Au collisions at 200 GeV for π0[335], η mesons [336], charged hadrons [35], and photons [337,338]. The nuclear modification factor depends, in general, on the identity, transverse momentum, pT, and pseudorapity, η, of the particle, as well as on the energy and centrality of the collision. As it can be seen in Fig. 4.1 for RHIC top energies, above pT∼5 GeV, π0[335], η[336], and charged hadrons [34,35] show all a common RAA ≈0.2. This suppression is much larger 98 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS than the one expected from nPDFs, which indicates that jet quenching effects are present in the QGP. Direct photons, which do not interact with the medium, have an RAA compatible with one [337]. The fact that RAA ≈0.2 irrespective of the nature of the produced hadron is consistent with the scenario where fragmentation into hadrons takes place in vacuum – but it is energy-rescaled. The nuclear modification factor at RHIC was fitted using ASW, HT and AMY energy loss models using a common three-dimensional relativistic fluid dynamics [339]. However, the value of the jet quenching coefficient, ˆq, obtained in the three approaches by fitting the data differs significantly. Let us move now from RHIC to the LHC, i.e., to O(10) times higher center of mass energy. As expected, the nuclear modification factor is now a bit smaller. For instance, for particles with pT∼6−7 GeV, RAA ≈0.13 in the most central Pb-Pb collisions at the LHC . A more interesting feature of RAA is its centrality- and pT-dependence. At RHIC top energies, the nuclear modification factor remains relatively constant from 5 GeV up to the highest transverse momenta measured so far, pT∼20 GeV (see Fig. 4.1). The much larger kinematical range opened at the LHC, allows us to observe an increase on RAA for higher pT’s up to RAA ≈0.4 for pT∼50 GeV in the most central collisions, see Fig. 4.2. As it is shown in the same figure, by moving from the most central to the most peripheral collisions the suppression is reduced. In a more peripheral collision, the overlap region where the QGP is formed is smaller and has lower density. This leads to a smaller path length, which is translated into a smaller amount of energy loss. Hence, the centrality dependence of the nuclear modification factor supports the picture of energy loss. Single inclusive measurements and in particular the nuclear modification factor are very useful. However, they have their limitations. As it was explained before, a factorization between medium effects and fragmentation is needed. A more general approach, where this separation is not a requirement, such as fully reconstructed jets, is desirable. Jet observables are very sensitive to jet quenching effects, but they will not be presented in this thesis, as they are not part of our analysis. Just to mention that the LHC has given us access to many jet measurements. Among them, there are the Dijet asymmetry for back-to- back jet pairs, AJ, measured in ATLAS [29] and CMS [30,42]; the Dijet azimuthal distribution [30], ∆φ; the average missing transverse momentum,h/ pk Ti[30]; the jet fragmentation function [340]; and jet shapes [341]. 4.2.2 High-pTflow harmonics The almost perfect fluid QGP may be understood through two key experimental signatures: collective flow and jet quenching. Collective flow observables have been successfully described by means of event-by-event (EbyE) relativistic viscous hydrodynamics [342–346]. More recently, the effects of event-by-event fluctuations have been analyzed in the context of hard probes. Measuring the azimuthal asymmetry of hard particles was proposed for the first time 4.2. PARTICLE OBSERVABLES 99 Figure 4.2: RAA versus pTfor charged hadrons measured for different centralities in Pb-Pb collisions at 2.76 TeV from ALICE [205]. in [347, 348]. It was suggested in there that the path length dependence due to an initial anisotropy would lead to a non-zero elliptic flow, v2, at high-pT. The first measurement of high-pTv2was published in 2006 by PHENIX Collaboration [349]. However, a simultaneous description of RAA and high-pTv2has been a puzzle during the last ten years. In fact, RAA was reasonably described by all the energy loss models, but the computed high-pTelliptic flow under-predicted the data [350]. Once EbyE fluctuations and more realistic initial conditions were implemented, a new theoretical high-pTv2definition was proposed that might solve this puzzle [351]. First of all, EbyE relativistic hydrodynamics is mandatory to study the flow harmonics. For instance, all hydrodynamic simulations with averaged smooth initial conditions predict a triangular flow, v3, equal to zero. Only in the case that event-by-event fluctuations are considered a non-zero v3can be obtained across all pT’s. However, until 2016 all the model calculations – energy loss approach plus relativistic hydrodynamic medium modeling – were incapable of computing both the nuclear modification factor and the high-pTflow harmonics. It was suggested in [351] that this may be due to the fact that the theoretically computed 100 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS quantity, vhard 2, was not the appropriate quantity to compare with the experimental data. In fact, the experimental pT>10 GeV flow coefficients vexp n(pT) are measured using the scalar product [352] vexp n(pT) = hvsoft nvhard n(pT) cos nψsoft n−ψhard n(pT)i svsoft n2,(4.2.15) where vsoft nand ψsoft nare the integrated soft flow harmonic and the corresponding event plane angle, and vhard n(pT) is the second Fourier coefficient of RAA(pT, φ), vhard n(pT) = 1 2πR2π 0dφ cos nφ −nψ hard nRAA(pT, φ) RAA(pT)(4.2.16) where ψhard n(pT) is given by ψhard n(pT) = 1 narctan R2π 0dφ sin(nφ)RAA(pT, φ) R2π 0dφ cos(nφ)RAA(pT, φ)!.(4.2.17) Hence, the calculation of high-pTflow harmonics vn(pT), requires modeling of both the soft and hard sectors of HICs. Note that in the case of smooth-averaged hydrodynamics there is only one “event” and, then, Eq. (4.2.15) is reduced to vhard n(pT), Eq. (4.2.16). In summary, EbyE fluctuations combined with the proper definition of vexp 2at high-pTcould provide a solution to the RAA ⊗v2puzzle. However, the energy loss model was implemented ad hoc in reference [351]. Therefore, a reliable energy loss model – combined with EbyE relativistic hydrodynamics – is mandatory to study the solution proposed in [351]. 4.3 Energy loss in the ASW framework In this section a brief summary of the path-integral formalism used to describe the propagation of hard particles through the QGP is given. Then the single-inclusive particle spectrum in the ASW approach (already introduced in Subsection 4.1.3) will be presented, as it is the formalism used in our analysis. Jet quenching is usually described by considering an elementary hard collision, with a cross section computed by pQCD, which produces a high-energy parton with large transverse momentum, pT, relative to the beam direction. In the high energy limit, as it was already explained, the dominant energy loss process is radiative (inelastic) energy loss. In this limit, the propagation time of the parton through the QGP is much smaller than the time scale of modifications of the medium. Consequently, the medium can be treated as a background field that interacts with the hard probe by means of very soft gluons [275, 277, 353], see Fig. 4.3. This limit can be summarized as Eω |k|,|p⊥|  T, ΛQCD ,(4.3.18) 4.3. ENERGY LOSS IN THE ASW FRAMEWORK 101 where Eis the (high) energy of the parton, ωand kare respectively the energy and momentum of the radiated gluons, and p⊥is the transverse momentum (with respect to the initial parton) accumulated by the projectile due to the radiative interactions with the medium. T and ΛQCD are the characteristic energy scales of the dense medium. Figure 4.3: Representation of multiple scattering of a high energetic quark with static medium components, represented as a small dark blob. The main effect of these inelastic interactions with the medium is the color rotation of the parton wave function, usually called eikonal phase [354]. As the energy of the parton is very high, in the reference frame where the parton is at rest, it traverses the target in such a small time that its transverse position does not vary during its propagation. This is known as eikonal approximation. Consequently, at ultrarelativistic energies, a color rotation on the parton wave function, due to the color field of the target, is induced. This effect is described by Wilson lines W(x0+, L+;x⊥) = Pexp (ig ZL+ x0+ dx+A−(x+,x⊥)),(4.3.19) where x⊥is the transverse position of the projectile, x0+and L+are the light-cone medium boundaries 2and A−≡TaAa −designs the medium color field components 3in a given light-cone ordering P. Wilson lines, therefore, are only valid to describe the propagation of partons that follow a straight line. A simple derivation of the Wilson line can be obtained in terms of multiple scatterings [355], as shown in Fig. 4.3. When the energy of the projectile is not as large, some subdominating terms with respect to the eikonal approximation can be incorporated, giving rise to a two dimensional pathintegral. This is the Green function, which replaces the Wilson line, and that allows some Brownian perturbations in the transverse plane of the propagating parton, 2The light-cone variables are given by x±= (x0±x3)/√2 and x⊥= (x1, x2). 3Tais given in the fundamental or adjoint representation whether the radiated quanta is, respectively, a quark or a gluon. 108 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS We use the quenching weights, P(), tabulated in [360], to model the amount of energy loss of highly energetic partons which will eventually fragment into a given hadron hin vacuum. The corresponding cross section reads dσAA→h dydpT =ZdqTdz dσAA→k dydqT P()Dk→h(z, µ2 F)δ(pT−z(1 −)qT),(4.6.28) where the cross section for producing parton kis dσAA→k dydqT =Zdx1dx2x1fA i(x1, µ2 F)x2fA j(x2, µ2 F)dˆσij→k dˆ t.(4.6.29) Note that any difference between parton and hadron rapidities is neglected here, and that all renormalization, factorization and fragmentation scales are taken to be equal, µF=pT. The partonic cross section is computed at NLO using the code in [372]. The free proton PDFs set CTEQ6.6M [373] with EPS09 [90] nuclear corrections is employed. Vacuum fragmentation functions DSS [374,375] are used. The quenching weights, P(), are defined as in Eq. (4.5.25) and they are computed in the multiple soft scattering approximation, see Subsection 4.4.2. It is worth emphasizing that this formalism is based on two main assumptions: i) the subsequent medium-induced gluon emissions are independent and ii) the fragmentation functions are not modified, i.e., fragmentation takes place in vacuum. Both of them find solid theoretical support in recent analyses of coherence effects in the medium. Beginning with the simplified setup of a QCD antenna [327,330,361,376], a pair of color-correlated partons with opening angle Θ emitting a soft gluon, a simple picture of jet quenching arises [377]: a medium of length Land jet quenching parameter ˆqhas a typical transverse momentum scale for color correlations Λ⊥∼1/√ˆqL; when the typical transverse size of the jet, r⊥∼ΘL, is smaller than this scale, r⊥<Λ⊥, the medium cannot resolve the inner structure of the jet, which remains unmodified and radiates medium-induced gluons as a whole with the total charge of the jet. This is the totally coherent case. Clearly, this picture implies that the fragmentation function remains basically unchanged if color coherence is maintained. However, it still depends on the fraction of momentum, z, and only a global energy loss affects the production of the fragmenting particles. This picture of jet quenching dictated by color coherence is in qualitative agreement with the experimental findings at the LHC [30,43,341,378–380] — see [381] for a quantitative analysis of some data. 4.6.2 From hydrodynamics to the transport coefficient The quenching weights are tabulated in [360] for the case of a static medium of finite length Land transport coefficient ˆq, where ωc=1 2ˆqL2, R =ωcL . (4.6.30) 4.6. ENERGY AND CENTRALITY DEPENDENCE OF ˆ Q109 For a medium in dynamic evolution, when the jet quenching parameter can be written as ˆq(τ)∼1/τα, a dynamical scaling law was found [382] which relates the resulting spectra with an equivalent static scenario. Based on this scaling law, effective ωeff cand Reff for an hydrodynamical medium profile are computed as ωeff c(x0, y0, τprod, φ) = Zdξ ξ ˆq(ξ),(4.6.31) Reff (x0, y0, τprod, φ) = 3 2Zdξ ξ2ˆq(ξ),(4.6.32) that reproduces Eqs. (4.6.30) for the static case. Similar implementations of the hydrodynamic model have been used before [339,357,366,367], so it is a rather standard procedure4. The production point of the parton at time τprod is distributed according to an Ncollscaling in the transverse plane and the azimuthal angle φis taken as a random number in [0,2π]. Each parton traverses the medium in a straight-line parameterized by the proper time ξat each point in the transverse plane. Then, it is only needed to specify a relation between the local value of the hydrodynamical variables at (x⊥(ξ), y⊥(ξ)) and the local value of the transport coefficient ˆq(ξ). Following [357], ˆq(ξ) = K·23/4(ξ),(4.6.33) where K≃1 corresponds to the ideal QGP, see the estimate in [368]. Other relations between the transport coefficient and the local thermodynamical quantities have been explored for instance in Ref. [339]. The local energy density (ξ) is taken from a hydrodynamic simulation of the medium. Several different options will be considered in the next section. In an expanding medium like this one, there is an ambiguity on the value of the transport coefficient, defined by Eq. (4.6.33), for values smaller than the proper time, τ0, when relativistic hydrodynamics is started. To quantify this uncertainty, three different extrapolations for the time from the production time to the proper time are considered: •Case (i): ˆq(ξ) = 0 for ξ < τ0. •Case (ii): ˆq(ξ) = ˆq(τ0) for ξ < τ0. •Case (iii): ˆq(ξ) = ˆq(τ0)/ξ3/4for ξ < τ0. The ˆq(ξ) = 0-extrapolation is an extreme case. It considers no energy loss at all before the thermalization time. This is a strong assumption as neither thermalization nor isotropization is necessary in the quenching weights approach. The second case assumes a continuous interaction from the production time (taken to be τprod ≃0.04 fm/c) and the third one considers a free-streaming medium with energy density decreasing as (ξ)∼1/ξ. 4Note that the prescription of Reff is slightly modified here. Now, Reff is the second moment of ˆq(ξ). The results are similar with the old prescription (see Eqs. (4.2)-(4.4) in Ref. [357]) but with improved stability for functional dependencies of ˆq(ξ) that are divergent in 1/ξ. 110 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS The production point of the hard scattering is characterized by a production weight w(x0, y0) calculated as w(x0, y0) = TA(x0, y0)TA(~ b−(x0, y0)) ,(4.6.34) where TAare the profile functions computed from a 3-parameter Fermi distribution at a given impact parameter ~ btaken from [383]. The average fragmentation functions for a parton k which has propagated through the medium and hadronizes in the vacuum to a hadron hcan be computed, in terms of these weights, as Dmed k→h(z, µ2 F) = 1 NZdφ dx0dy0w(x0, y0)Zdζ 1−ζPk(x0, y0, φ, ζ)Dvac k→hz 1−ζ, µ2 F, (4.6.35) where Pk(x0, y0, φ, ζ) is the quenching weight for parton kand the normalization is N= 2πRdx0dy0w(x0, y0). With this definition, the cross section, Eq. (4.6.28), is simply dσAA→h dydpT =ZdqTdz dσAA→k dydqT Dmed k→h(z, µ2 F)δ(pT−zqT).(4.6.36) 4.6.3 Hydrodynamic modeling of the medium The space-time distribution of the local energy density is obtained by solving the relativistic hydrodynamic equations. Such simulations require the specification of initial values for the energy momentum tensor, as well as parameters that describe medium properties, neither of which are accurately known. In order to test the robustness of our results and conclusions with respect to these uncertainties, we repeat all calculations using space-time profiles from various different smooth-averaged hydrodynamic descriptions. The first, which we refer to as “Hirano”, corresponds to the calculation described in [384–386], to which we refer the reader for details. In short, this simulation uses an optical Glauber model where the initial entropy density at initial proper time τ0= 0.6 fm is given by a linear combination of the number density of participant nucleons, ρpart, and binary collisions, ρbin: s∝(1 −x)ρpart +xρcoll ,(4.6.37) with binary collision fraction x= 0.15. A bag model equation of state is used, with chemical freeze out enforced at Tch = 170 MeV, and kinetic freeze out at Tf= 100 MeV, below which temperature the medium has frozen out and no energy loss occurs. This is an ideal hydrodynamic calculation, with vanishing viscosity. The other two hydrodynamical models correspond exactly to the calculations in [387] (for 200 GeV Au-Au collisions at RHIC) and [388] (for 2.76 TeV Pb-Pb collisions at the LHC), to which we again refer the reader for all relevant details. One calculation, which we refer to as “Glauber”, uses as initial condition an energy density proportional to the density of binary collisions, ρbin, while the ratio of shear viscosity to 4.6. ENERGY AND CENTRALITY DEPENDENCE OF ˆ Q111 entropy density is fixed to a constant value of η/s = 0.08. The final calculation is referred to as “fKLN”. This simulation takes its initial condition from a factorized Kharzeev-Levin-Nardi model [389], with the shear viscosity set to η/s = 0.16. Both of the latter simulations begin at an initial proper time of τ0= 1 fm and use an equation of state inspired by lattice QCD calculations. Each system is assumed to be in chemical equilibrium until it reaches a freeze out temperature of Tf= 140 MeV. All of these calculations have been successfully tested against several experimental data, but use different choices for initial conditions, thermalization time, viscosity, equation of state, etc. Thus, we expect that the variation in our results from using these different models should give a reasonable indication of the uncertainty coming from the hydrodynamic background. It will be shown that such uncertainty is negligible with respect to our main conclusions. 4.6.4 Results This analysis is restricted to single-inclusive suppression at RHIC [390] and the LHC [205]. Other observables are not considered here, as they may involve other effects related with fragmentation, mass effects on the energy loss mechanism, etc. Something new is, the centrality dependence of both RHIC and LHC data. Most of previous studies have analyzed the most central collisions or the centrality dependence only for one energy [339, 366,367]. Therefore, this is the first study of both centrality and (center of mass) energy dependence of the nuclear modification factor. In Fig. 4.6 our results for different values of Ktogether with the experimental data from the PHENIX Collaboration [390] on suppression of inclusive neutral pions on Au-Au collisions at √sNN = 200 GeV (72 data points) are plotted. In Fig. 4.7 we plot the corresponding results for LHC Pb-Pb collisions at √sNN = 2.76 TeV from ALICE [205], 156 data points 5. We restrict to pT>5 GeV/c to stay in a region where pQCD can be applied and no large contribution from other effects like flow is expected. In both figures “Hirano” hydrodynamic simulation was used. The results with the other two above-mentioned hydrodynamic profiles are very similar to these ones. We have performed a χ2fit to the best value of Kfor each energy and centrality, and for each assumption of hydrodynamical profile or behavior of ˆqat values of proper time smaller than the thermalization time, τ0, assumed in each hydrodynamical simulation. For the case of ALICE data [205] we add the systematic and statistical errors in quadrature, as no particular instructions of how to include them in a fit are provided. For the case of RHIC, the latest analysis includes the contribution from several different error sources. The two methods lead 5Note that the plotted values of Kare different in both figures. 112 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS 0 2 46 8 10 12 14 16 18 20 pT (GeV) 0.2 0.4 0.6 0.8 1.0 RAA Au-Au 200 GeV 0-5% PHENIX 0 2 46 8 10 12 14 16 18 20 pT (GeV) 0.2 0.4 0.6 0.8 1.0 RAA Au-Au 200 GeV 0-10% PHENIX 0 2 46 8 10 12 14 16 18 20 pT (GeV) 0.2 0.4 0.6 0.8 1.0 RAA Au-Au 200 GeV 10-20% PHENIX 0 2 46 8 10 12 14 16 18 20 pT (GeV) 0.2 0.4 0.6 0.8 1.0 RAA Au-Au 200 GeV 20-30% PHENIX 0 2 46 8 10 12 14 16 18 20 pT (GeV) 0.2 0.4 0.6 0.8 1.0 RAA Au-Au 200 GeV 30-40% PHENIX Figure 4.6: Suppression of inclusive π0in Au-Au collisions at √sNN = 200 GeV for different values of the parameter K(see Eq. (4.6.33)) compared with PHENIX data at different centralities [390]. Curves from top to bottom correspond to K=K0/1.46, with K0= 2,2.25,2.5,...,6, using the “Hirano” hydrodynamical model and the energy density prior to the start of hydrodynamical evolution taken as constant, see the previous subsection. 4.6. ENERGY AND CENTRALITY DEPENDENCE OF ˆ Q113 0 10 20 30 40 50 pT (GeV) 10 − 1 1.0 RAA PbPb 2.76 TeV 0-5% ALICE 0 10 20 30 40 50 pT (GeV) 10 − 1 1.0 RAA PbPb 2.76 TeV 5-10% ALICE 0 10 20 30 40 50 pT (GeV) 10 − 1 1.0 RAA PbPb 2.76 TeV 10-20% ALICE 0 10 20 30 40 50 pT (GeV) 10 − 1 1.0 RAA PbPb 2.76 TeV 20-30% ALICE 0 10 20 30 40 50 pT (GeV) 10 − 1 1.0 RAA PbPb 2.76 TeV 30-40% ALICE 0 10 20 30 40 50 pT (GeV) 10 − 1 1.0 RAA PbPb 2.76 TeV 40-50% ALICE Figure 4.7: Suppression of inclusive charged particles in Pb-Pb collisions at √sNN = 2.76 TeV for different values of the parameter K(see Eq. (4.6.33)) compared to ALICE data at different centralities [205]. Curves from top to bottom correspond to K=K0/1.46, with K0= 0.5,0.7,0.9,...,3.1, using the “Hirano” hydrodynamical model and the energy density prior to the start of hydrodynamical evolution taken as constant, see the previous subsection. 114 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS 246 8 10 12 b (fm) 2 4 6 8 10 12 14 16 18 K =ˆ q / 2 ² 3 / 4 ˆ q( τ ) =0 , τ <τ 0 RHIC 200 GeV Hirano fKLN Glauber 246 8 10 12 b (fm) 2 4 6 8 10 12 14 16 18 K =ˆ q / 2 ² 3 / 4 ˆ q( τ ) =0 , τ <τ 0 LHC 2.76 TeV Hirano fKLN Glauber Figure 4.8: K-factors obtained from fits to PHENIX RAA data [390] (left panel) and to ALICE RAA data [205] (right panel) using different hydrodynamical profiles as a function of the average impact parameter for each centrality class and for ˆq(ξ) = 0 before thermalization, see the previous subsections. to comparable values of K(differences ∼5%) except for the most peripheral bins, for which the Kvalues in the case of errors added in quadrature are ∼30% smaller. The uncertainty band is determined by ∆χ2= 1. In order to make the comparison between RHIC and the LHC, these issues need to be taken into account, although the conclusions do not change at the qualitative level. In the left panels of Fig. 4.8, Fig. 4.9 and Fig. 4.10 we plot the different values of the K-parameter fitted to the PHENIX data [390] for different combinations of hydrodynamic profiles and behavior before the proper time. The corresponding values for the LHC [205] are plotted in the right panels of the same figures. Let us comment now on these results. First, the extracted values of Kare compatible for the cases of either frozen energy density or free streaming before τ0, and the results for the three different hydrodynamic models are similar. This is not the case when no quenching is assumed before τ0; for this assumption, the two viscous hydrodynamic implementations which use a common (larger) τ0require a larger Kthan the ideal hydrodynamic model that considers a smaller τ0, with actual values which become unrealistically large. Therefore, we do not consider the results obtained for this assumption for the discussion of the values of K, but the qualitative behavior that we find is in agreement with the two other assumptions. In any case they clearly illustrate the importance of the treatment of early times in jet quenching computations. Second, for the most peripheral collisions at the LHC, model “Glauber” demands a much larger Kthan the others, while model “Hirano” returns a rather flat value of Kfor all centralities. Third, the trend of the results at RHIC is a slight decrease with 4.6. ENERGY AND CENTRALITY DEPENDENCE OF ˆ Q115 246 8 10 12 b (fm) 1.0 1.5 2.0 2.5 3.0 K =ˆ q / 2 ² 3 / 4 ˆ q( τ ) =ˆ q( τ 0) , τ <τ 0 RHIC 200 GeV Hirano fKLN Glauber 246 8 10 12 b (fm) 1.0 1.5 2.0 2.5 3.0 K =ˆ q / 2 ² 3 / 4 ˆ q( τ ) =ˆ q( τ 0) , τ <τ 0 LHC 2.76 TeV Hirano fKLN Glauber Figure 4.9: K-factors obtained from fits to PHENIX RAA data [390] (left panel) and to ALICE RAA data [205] (right panel) using different hydrodynamical profiles as a function of the average impact parameter for each centrality class and the energy density prior to the start of hydrodynamical evolution taken as constant, see the previous subsections. decreasing centrality, although compatible with constant, while at the LHC the behavior is constant except for the smaller centralities, where the behavior, as it was mentioned above, depends very much on the hydrodynamic profile employed. Finally, we would like to understand the systematics and relation of LHC and RHIC results for the K-factor that we obtain. First, we notice that, in principle, Eq. (4.6.33) determines how far or close the perturbative estimate ˆq≃23/4is from our value fitted to experimental data. In this sense, we note that there is a clear departure from unity of this value for the case of RHIC. This fact was found several times [339,357]6. We also find that the corresponding value of Kis smaller at the LHC, a fact which has been already found before by other groups [366] but with a smaller decrease (a factor ∼25 % compared to our factor 2–3)7. The study of the centrality dependence is, nonetheless, more surprising. The extracted value of Kseems to depend mainly on the energy of the collision and much less (if any) on the centrality. This is not the behavior one would expect from a naive interpreta- 6Note that the difference of the present extraction K∼2–3 and K∼4 from [357] comes mainly from the new definition of Rin (4.6.32), indicating, again, the important role of the geometry in the extraction of ˆq. 7Nevertheless this comparison needs to be taken with caution as the values of ˆq/T3quoted in [366] are performed at a given temperature and no systematics with temperature is presented. Moreover, ˆqis not the natural fitting parameter in the models studied in that reference but a derived quantity once the parameters of the different models are extracted from the data. In our case, ˆqis the natural parameter, given by (4.4.24), and the K-factor has a well defined meaning. 116 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS 246 8 10 12 b (fm) 1.0 1.5 2.0 2.5 K =ˆ q / 2 ² 3 / 4 free streaming, τ <τ 0 RHIC 200 GeV Hirano fKLN Glauber 246 8 10 12 b (fm) 1.0 1.5 2.0 2.5 K =ˆ q / 2 ² 3 / 4 free streaming, τ <τ 0 LHC 2.76 TeV Hirano fKLN Glauber Figure 4.10: K-factors obtained from fits to PHENIX RAA data [390] (left panel) and to ALICE RAA data [205] (right panel) using different hydrodynamical profiles as a function of the average impact parameter for each centrality class and for the free-streaming extrapolation, see the previous subsections. tion in which the K-factor only indicates the departure from the leading order perturbative estimate determined by temperature. In this naive interpretation, a medium with a smaller temperature (RHIC) would need higher orders of the perturbative series to be included, while a medium at higher temperature would be closer to the ideal limit. This simple interpretation does not correspond, notwithstanding, to the present findings as there is an overlap on typical energy densities between central Au-Au at RHIC and semi-peripheral Pb-Pb at the LHC, so their values of Kshould coincide in this naive interpretation. In order to provide an estimate of this overlap, we plot, in Figure 4.11, the K-factors obtained for different centralities and energies versus an energy density times formation time τ0extracted from the experimental data using Bj¨orken estimates [391,392] – we have checked that the overlap is similar if we plot as a function of the maximum energy density of the hydrodynamical profiles that we have used to perform the fits. 4.6.5 Limitations and conclusions One-particle inclusive suppression of particles produced at high transverse momenta at RHIC and the LHC as a function of centrality has been studied. By defining a constant K-factor with respect to the perturbative estimate ˆq≃23/4we fit the corresponding experimental data at RHIC and LHC for different centralities. The fitted value at RHIC confirms previous estimates [339,357] of large corrections to the ideal case, although the actual numerical 4.6. ENERGY AND CENTRALITY DEPENDENCE OF ˆ Q117 246 8 10 12 ²τ 0 (GeV / fm2 / c) 1.0 1.5 2.0 2.5 3.0 K =ˆ q / 2 ² 3 / 4 ˆ q( τ ) =ˆ q( τ 0) , τ <τ 0 Hirano RHIC fKLN RHIC Glauber RHIC Hirano LHC fKLN LHC Glauber LHC Figure 4.11: K-factor obtained from fits to RAA data at RHIC and LHC energies for different centrality classes plotted as a function of an estimate of the energy density times formation time τ0of the QCD medium formed in each case. The τ0estimates are taken from [391,392]. value is a bit smaller, due to a new, more stable, definition of the effective values of the static scenario equivalent to the evolving medium, Eq. (4.6.32). For the case of the LHC, instead, the extracted value of Kis close to unity. One would be tempted to make the naive interpretation that the medium created at the LHC, having a larger temperature, is closer to the ideal case than the one at RHIC, for which larger corrections or even a strongly coupling treatment, could be needed. This naive interpretation finds difficulties to be accommodated, however, with the fact that the centrality dependencies at RHIC and the LHC separately are rather flat, that is, the change in the value of Kis not simply due to the different temperature (or energy density), as there is a large region of overlap between RHIC and the LHC for different centralities. At this moment we do not have an interpretation for this finding which, in any case, should be checked by other model implementations of jet quenching. It is also worth noticing that the extraction of the value of Kin the case of RHIC depends on a single set of experimental data, namely inclusive π0suppression measured by PHENIX. The corresponding results from STAR on π++π−suppression [393] show a smaller suppression but the smaller range of transverse momentum studied makes our analysis to be not very reliable. For this reason we have chosen not to include this set of data in the fit. For the LHC, on the other hand, CMS [36] and ATLAS [394] have measured the suppression of inclusive charged particles with results almost identical to the ones from the ALICE collaboration8. 8On the other hand, ALICE data are restricted to mid-rapidities where the boost invariant picture of 124 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS Percolation of strings The string percolation model (SPM) is a non-perturbative approach which captures some of the main features of the Color Glass Condensate (CGC). In the SPM, the multiplicity and the mean pT-distributions result from the formation of clusters of strings which introduce correlations among the produced particles. Several observables have been analyzed in chapter 3 within this framework. It has been shown that the geometric scaling observed both in HICs and in p-p collisions at different energies and centralities can be derived in the SPM. The anomalous dependence of the near-side ridge structure in Au-Au collisions at RHIC can be also naturally explained in this model. Moreover, the onset of this structure in high multiplicity events in proton-proton collisions is also described. Finally, experimental data of elliptic flow are shown to exhibit geometric scaling in Pb-Pb collisions at the LHC and Au-Au collisions at RHIC. This scaling law is also extended to the elliptic flow of identified particles: pions, kaons, and protons. It is suggested that this scaling may be due to the energy lost by the interaction of the parton with the color field of the strings. •Final state effects Regarding final state effects, hard probes and, more concretely, jet quenching has been studied in chapter 4. We have in there performed one of the most ample analysis of experimental data and its consequences for the QGP. Single-particle inclusive suppression of particles produced at high-pTboth at RHIC and at the LHC as a function of centrality has been presented. By defining a constant K-factor with respect to the perturbative estimate, ˆq=K23/4–being the energy density – the corresponding nuclear modification factor data at RHIC and LHC for different centralities have been fitted. This is done by combining our energy loss model, the ASW Quenching Weights, with different smooth-averaged hydrodynamic simulations, to check that the outcome is independent of the hydrodynamic model employed. The K-factor obtained at the LHC is approximately 1. The value at RHIC is ∼2−3 times larger, confirming large corrections to the ideal case. However, its dependence on the centrality at RHIC and at the LHC separately is rather flat. Whether the K- factor was determined, by, for instance the temperature, the most central collisions at RHIC should present a value similar to semi-peripheral LHC data. Consequently, the K-factor would not depend on the local properties of the medium as energy density or temperature, but on global collision quantities such as the center of mass energy. This is a very unexpected result for which we cannot yet provide a clear interpretation. Finally, an equivalent work to the one just described, but using now an event-by-event hydrodynamic profile, EKRT EbyE hydrodynamics, is shown. This analysis leads to the same conclusions that our previous one, where smooth-averaged hydrodynamic 4.6. ENERGY AND CENTRALITY DEPENDENCE OF ˆ Q125 simulations were used. It is worth stressing that this new study can be regarded as a first step towards the simultaneous description of the nuclear modification factor and the flow harmonics at high-pT. Nowadays, there is not any realistic and reliable model of jet quenching which enable to compute both observables. However, recent works show the possibility of combining an event-by-event hydrodynamics and a complete formalism of jet quenching, as those implemented in our analysis. Therefore, this last study is a good starting point in this direction. 126 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS Resumen La fuerza fuerte es una de las cuatro interacciones fundamentales de la naturaleza. Dicha fuerza es la responsable de que los nucleones – protones y neutrones – permanezcan unidos formando los n´ucleos at´omicos. La teor´ıa que describe dicha interacci´on data de los a˜nos 1970 y se denomina Cromodin´amica Cu´antica (QCD, por sus siglas en ingl´es). La QCD describe las interacciones entre los quarks y gluones, part´ıculas fundamentales que constituyen los hadrones. Actualmente la QCD es considerada uno de los pilares fundamentales del denominado Modelo Est´andar de las part´ıculas elementales y sus interacciones. Entre las propiedades de la QCD cabe destacar dos de gran relevancia: la libertad asint´otica y el confinamiento. La libertad asint´otica es la reducci´on de la fuerza de las interacciones entre los quarks y gluones a medida que la escala de la energ´ıa de dichas interacciones aumenta, o lo que es lo mismo, la de distancia decrece. Por tanto, los quarks se mueven en el interior de los hadrones como part´ıculas libres, lo que permite el uso de lo que se conoce como t´ecnicas perturbativas. Sin embargo, a distancias largas, la fuerza de esta interacci´on aumenta confinando los quarks y gluones en el interior de los hadrones. Esta es la raz´on por la que quarks y gluones no existen de forma aislada en la naturaleza, sino que est´an formando parte de los hadrones. Las colisiones de iones pesados (HICs, por sus siglas en ingl´es) son la herramienta fundamental para estudio de la Cromodin´amica Cu´antica bajo condiciones extremas de temperatura y densidad, muy distintas a las que existen en el interior del n´ucleo at´omico. En estas colisiones nucleares de alta energ´ıa se alcanzan temperaturas y densidades cr´ıticas que permiten la formaci´on del denominado quark-gluon plasma (QGP). La existencia de esta nueva forma de materia nuclear fue predicha por primera vez en los a˜nos 1970. Bajo estas condiciones, las interacciones de distancia corta comienzan a dominar sobre las de distancia larga que empiezan a verse apantalladas por las fuentes de color a su alrededor. Las interacciones de distancia corta se caracterizan por una constante de acoplamiento peque˜na debido a la libertad asint´otica caracter´ıstca de la QCD. Por tanto, el QGP, est´a constituido por quarks y gluones – grados de libertad del lagrangiano de QCD – deconfinados. Este plasma es un l´ıquido casi perfecto, cuyos componentes est´an estrechamente unidos. Por ello, su estudio nos puede ayudar a entender mejor la parte no perturbativa de la QCD, de la cual el conocimiento hoy en d´ıa es limitado. Adem´as, el QGP es descrito por los modelos cosmol´ogicos actuales como el estado de la 127 128 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS materia durante los primeros microsegundos tras el Big-Bang. Por consiguiente, el estudio del QGP en colisiones de iones pesados, nos permite analizar la evoluci´on de este estado de la materia (deconfinada) que exist´ıa en los inicios del Universo hacia la materia normal confinada. O lo que es lo mismo, nos pueden revelar informaci´on sobre el origen de la masa y del confinamiento. Las colisiones nucleares de alta energ´ıa son, en definitiva, una gran oportunidad de observar y entender el origen y la evoluci´on de nuestro Universo. Esto ha llevado al desarrollo de diversos programas de f´ısica de iones pesados: en el Alternating Gradient Synchrotron (AGS) en Brookhaven National Laboratory (BNL), en el Super Proton Synchrotron (SPS) en el Conseil Europ´een pour la Recherche Nucl´eaire (CERN), en el Relativistic Heavy Ion Collider (RHIC) en BNL y en el Large Hadron Collider (LHC) en CERN. Seg´un c´alculos en Lattice QCD la temperatura cr´ıtica para la formaci´on del QGP es Tc≃154 MeV, que se corresponde a una densidad cr´ıtica de energ´ıa de c∼1 GeV/fm3. Los datos experimentales muestran desde hace mucho tiempo que la temperatura alcanzada en las colisiones relativistas de iones pesados es mayor que dicha temperatura cr´ıtica. En consecuencia, este tipo de colisiones son las propicias para estudiar las propiedades del QGP. Los sistemas producidos en estas colisiones se caracterizan por la presencia de fen´omenos colectivos y efectos del medio (del QGP) que son accesibles experimentalmente. El estudio de procesos que son sensibles al grado de colectividad del sistema es el objetivo fundamental de las colisiones de iones pesados. Sin embargo, la vida media del plasma de quarks y gluones producido en HICs es muy peque˜na – del orden de 100 ys = 10−24 s – y, por tanto, no puede ser detectado directamente. Las propiedades de este estado de la materia han de ser estudiadas indirectamente en el estado hadr´onico final de la colisi´on. Generalmente, se utilizan distribuciones de part´ıculas soft para ver indicios de comportamiento colectivo del QGP y tratar de dar una posible descripci´on hidrodin´amica del mismo; al mismo tiempo que se usan hard probes para ver el efecto del medio sobre procesos que puedan ser calculados perturbativamente. Esta memoria se centra en el an´alisis de dos tipos de efectos: •Efectos de estado inicial (IS). Son efectos previos a la termalizaci´on. –Funciones de distribuci´on part´onicas nucleares (nPDFs). Un ingrediente b´asico para entender las colisiones de iones pesados y poder describir, por tanto, el QGP formado en ellas, son las funciones de distribuci´on part´onicas nucleares (nPDFs, por sus siglas en ingl´es). Estas contienen la informaci´on de la estructura part´onica (quarks y gluones) de los protones y neutrones que forman los n´ucleos que colisionan. Son, por tanto, contribuciones de distancia larga que no pertenecen al dominio perturbativo de la QCD. Su determinaci´on precisa es crucial para la correcta interpretaci´on de todos los observables utilizados en HICs. Gracias a su universalidad y al hecho de que su evoluci´on respecto a una determinada escala inicial s´ı es perturbativa, las nPDFs son obtenidas mediante una t´ecnica denominada an´alisis global. La extracci´on de estas distribuciones 4.6. ENERGY AND CENTRALITY DEPENDENCE OF ˆ Q129 mediante un an´alisis global a next-to-next-to leading order (NNLO) en QCD es el objeto del cap´ıtulo 2 de esta tesis. La t´ecnica del an´alisis global – o fit global – consiste en extraer las nPDFs a partir de diversos experimentos: dispersi´on inel´astica profunda – DIS, por sus siglas en ingl´es –, Drell-Yan, etc. El procedimiento es el que sigue. Primero, se seleccionan los observables experimentales a usar y se computan a nivel part´onico al orden deseado en QCD perturbativa, en nuestro caso a NNLO. A continuaci´on, las nPDFs son parametrizadas a una escala inicial, Q0. Para cada set de datos, las nPDFs son evolucionadas usando las conocidas como ecuaciones de Dokshitzer-Gribov- Lipatov-Altarelli-Parisi (DGLAP) (tambi´en a NNLO) desde la escala inicial a la escala del proceso, Q > Q0. En nuestro caso esta evoluci´on es realizada en el espacio Mellin. Despu´es, se convoluciona la parte de peque˜na distancia (hard) con la de larga distancia (soft) para obtener el observable en espacio Mellin. Posteriormente se realiza una inversi´on de Mellin para obtener la predicci´on te´orica del observable en espacio x. Finalmente, se construye una distribuci´on χ2y se extrae el valor de los par´ametros minimizando iterativamente el χ2. En el cap´ıtulo 2 de esta tesis se presenta dicho fit global de nPDFs a NNLO. Este an´alisis incorpora un tratamiento ´ıntegro de los efectos de quarks pesados siguiendo el denominado Esquema General con N´umero Variable de Sabores (GMVFNS, por sus siglas en ingl´es). Este es el primer fit global hecho a NNLO en QCD con un tratamiento exhaustivo de los efectos de masa de los quarks. En este fit se utilizan datos de DIS con leptones cargados y de DIS con neutrinos. Cabe destacar que entre el primer tipo de datos, se incluyen datos donde el blanco es deuterio, tratando as´ı los efectos nucleares del deuterio normalmente obviados en los an´alisis globales de funciones de distribuci´on part´onicas nucleares. Por otra parte, se han cuantificado las incertidumbres experimentales dando lugar a 24 sets de errores, adem´as del set central de nPDFs. Adem´as de ser imprescendibles para entender las ya mencionadas colisiones de iones pesados que tiene lugar en el Large Hadron Collider (LHC) en CERN (Suiza) y en Relativistic Heavy-Ion Collider (RHIC) en Brookhaven (USA), el conocimiento preciso de las nPDFs ser´a crucial para los futuros aceleradores, como el Electron Ion Collider (EIC) y el Large Hadron-Electron Collider (LHeC). Por ello su determinaci´on cada vez m´as precisa es una l´ınea de investigaci´on imperativa no solo para la comunidad de iones pesados, sino para la f´ısica de altas energ´ıas actual y futura. –Percolaci´on de cuedas. Otro tema de gran inter´es son los fen´omenos colectivos que dan lugar a la producci´on del QGP. Estos son habitualmente analizados bajo el marco del Clolor 130 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS Glass Condensate (CGC) en HICs. Un modelo simplificado que captura varias de las propiedades del CGC es la percolaci´on de cuerdas. Este enfoque y algunos de sus resultados comparados con los datos experimentales se presentan en el cap´ıtulo 3. Las colisiones de iones pesados son descritas en este modelo mediante la formaci´on de cuerdas que unen los partones del proyectil y del blanco de la colisi´on. Estas cuerdas de color pueden ser vistas como cilindros expandidos en la direcci´on longitudinal cuyas intersecciones con el plano transverso son peque˜nas ´areas llenas del campo de color de los partones colisionantes. Al crecer el n´umero at´omico o la energ´ıa de la colisi´on, el n´umero de cuerdas crece y estas empiezan a solaparse formando agregados, usualmente denominados clusters. Esto da lugar a tres importantes consecuencias: -La carga de color de un cluster de ncuerdas es √n– y no n– veces la carga elemental de las cuerdas. Como consecuencia, la multiplicidad se reduce respecto al caso de cuerdas independientes. -Las part´ıculas son creadas v´ıa producci´on de pares quark-antiquark en el campo de color mediante un mecanismo tipo Schwinger. -La percolaci´on es un fen´omeno cr´ıtico: alcanzada una cierta densidad de cuerdas, denominada densidad cr´ıtica o de percolaci´on, se forma un ´unico cluster macrosc´opico que ocupa toda la superficie de solapamiento nuclear. Diversos observables en HICs ser´an analizados en el cap´ıtulo 3 mediante dicho formalismo. Entre estos, se encuentra el denominado escaleo geom´etrico observado no solo en colisiones de iones persados, sino tambi´en en colisiones prot´on-prot´on. En dicho cap´ıtulo, se mostrar´a que este escaleo puede ser explicado mediante la percolaci´on de cuerdas para distintas colisiones, energ´ıas y centralidades. Otro fen´omeno analizado es la llamada near-side ridge structure. La dependencia an´omala de esta estructura encontrada en colisiones de n´ucleos de oro ser´a explicada de forma natural en el contexto de la percolaci´on de cuerdas. Adem´as, la aparici´on de dicha estructura en eventos de alta multiplicidad en colisiones prot´on-prot´on tambi´en ser´a descrita. Por ´ultimo, se mostrar´a que los datos experimetales de flujo el´ıptico en colisiones de plomo-plomo en el LHC y de oro-oro en RHIC tambi´en satisfacen una ley de escaleo. Dicha ley se extiende tambi´en al flujo el´ıptico de piones, kaones y protones. Se mostrar´a que dicho escaleo puede tener como origen la p´erdida de energ´ıa debida a la interacci´on del part´on emitido con el campo de color de las cuerdas. •Efectos de estado final (FS). Entre ellos, cabe destacar las hard probes, esto es, part´ıculas caracterizadas por un alta 4.6. ENERGY AND CENTRALITY DEPENDENCE OF ˆ Q131 energ´ıa o masa cuyo estudio pertenece al r´egimen perturbativo de la QCD. El estudio de procesos tipo hard (hard probes) se ha convertido en uno de los campos m´as activos de investigaci´on en colisiones de iones pesados debido al amplio rango de escalas de energ´ıa que engloba. –Supresi´on de jets Su comportamiento en vac´ıo, esto es, en colisiones en las que no se forma un plasma, es conocido. En consequencia, se analizan las modificaciones debidas a la presencia de un medio nuclear, es decir, debidas a la fomaci´on del QGP, permiti´endonos extraer as´ı las propiedades de ´este. Entre estas hard probes cabe destacar la supresi´on de jets (jet quenching). La supresi´on de jets, y m´as concretamente, la supresi´on inclusiva de part´ıculas de alto momento en HICs ser´a analizada en el cap´ıtulo 4. El fen´omeno supresi´on de jets se refiere m´as concretamente al conjunto de efectos que sufren las part´ıculas tipo hard al propagarse por el QGP. Fue propuesto por Bjorken en 1982 y m´as tarde observado en RHIC y LHC. Hoy en d´ıa se trata de un fen´omeno bien establecido del que se tienen muchos datos experimentales tanto de RHIC como del LHC. No obstante, desde el punto de vista te´orico no ha sido todav´ıa completamente entendido. En una colisi´on de iones pesados los partones hard se forman con el QGP, el cual tienen que atravesar antes de dar lugar a las part´ıculas finales detectadas en los experimentos. Durante dicha propagaci´on sufren jet quenching. El objetivo fundamental de los anal´ısis de dicho fen´omeneno es extraer las propiedades del QGP a trav´es del estudio de la modificaci´on de la propagaci´on de los partones debido a la presencia del plasma. Uno de los estudios m´as extensos de datos experimentales y sus consecuencias para las propiedades del medio ser´a presentado en cap´ıtulo 4. En dicho estudio, se analizan datos tanto de RHIC como de LHC del denomiado factor de modificaci´on nuclear para producci´on single-inclusive de part´ıculas de alto momento transverso. La p´erdida de energ´ıa en el medio es analizada mediante los denominados Quenching Weights (QW) y la evoluci´on din´amica del medio es tratada mediante distintas simulaciones hidrodin´amicas. El resultado fundamental de este an´alisis es la extracci´on del denominado coeficiente de tranporte de jets, ˆq. En nuestro estudio el coeficiente de tranporte de jets, ˆq, es determinadado por la densisdad de energ´ıa local, , extraida de diversos modelos hidrodin´amicos. El coeficiente de transporte es definido como ˆq=K23/4y los valores de Kson fiteados a los datos experimentales del factor de modificaci´on nuclear a distintas centralidades y energ´ıas. El resultado obtenido es que este factor Kes ∼2−3 mayor para RHIC que para el LHC y que, sorprendentemente, este factor Kparece no 132 CHAPTER 4. SUPPRESSION OF HIGH-PTPARTICLES IN HICS depender de las propiedades locales del medio, como la temperatura. De hecho, Kdepende principalmente de la energ´ıa del centro de masa de la colisi´on y es b´asicamente independiente de la centralidad de la misma. El formalismo de p´erdida de energ´ıa utilizado en este trabajo, son los ya mencionados Quenching Weights. Estos son unas probabilidades de p´erdida de energ´ıa tipo Poisson basadas en la emisi´on independiente de gluones inducida por el medio y calculados en el formalismo de Armesto-Salgado-Wiedemann (ASW). Para modelar la parte soft se utilizan primero tres perfiles hidrodin´amicos smooth-averaged que usan diferentes ecuaciones de estado, condiciones iniciales, tiempo de termalizaci´on, etc., demostrando as´ı que nuestras conclusiones son independientes del modelo hidrdin´amico empleado. En el cap´ıtulo 4 este an´alisis, sus resultados, limitaciones y conclusiones son presentados en detalle. Por ´ultimo, se presenta un trabajo equivalente al anterior, pero usando una simulaci´on hidrodin´amica event-by-event. Este estudio lleva a las mismas conclusiones que el previo, en el que los modelos hidrodin´amicos empleados eran smoothaveraged. Cabe destacar, adem´as, que este nuevo an´alisis sirve como un paso previo hacia la descripci´on simult´anea del factor de modificaci´on nuclear y los arm´onicos a alto momento transverso. A d´ıa de hoy, ning´un modelo realista de jet quenching ha sido capaz de determinar correctamente ambos observables. Sin embargo, estudios recientes muestran la posibilidad de calcular ambos mediante la combinaci´on de un modelo hidrodin´amico event-by-event y un formalismo completo para el jet quenching, como los que tenemos implementados en nuestro an´alisis. Por ello, este ´ultimo trabajo sirve de punto de partida en esta direcci´on. En resumen, en esta tesis se consideran diversos aspectos relacionados con los estados inicial y final en colisiones de iones pesados. En la introducci´on se presentan algunos aspectos generales de la Cromodin´amica Cu´antica motivando su estudio en colisiones de iones pesados. El el cap´ıtulo dedicado a las funciones de distribuci´on part´onicas nucleares, cap´ıtulo 2, se presenta un fit global a next-to-next-to leading order en QCD perturbativa. En el cap´ıtulo 3, se describe el formalismo de percolaci´on de cuerdas y diversos observables en HICs son analizados mediante el uso de dicho modelo. Por ´ultimo, en el cap´ıtulo 4 se presenta un an´alisis de la supresi´on inclusiva de part´ıculas de alto momento a diversas energ´ıas y centralidades, cuyo resultado principal es la determinaci´on del denominado coeficiente de tranporte de jets, ˆq. Appendix A About the Mellin technique in global analysis The Mellin transform of a function f(x) is given by ˆ f(N) = Z1 0 dx xN−1f(x),(A.1) where Nis an integer number. An analytic continuation to complex values of the argument can be done [398]. The DGLAP equations, Eqs. (2.1.12) and (2.1.13), where written in terms of a convolution of the PDFs (soft part) and the splitting functions (hard part). Here, it will be proven than this kind of convolution in the x-space, when a Mellin transform is applied, can be written as a product. The Mellin convolution is f⊗g(x)≡Z1 x dz zf(z)gx z,(A.2) a Mellin transform is applied, [ f⊗g(N) = Z1 0 dx xN−1Z1 x dz zf(z)gx z.(A.3) Changing the order of the integrals, [ f⊗g(N) = Z1 0 dz zf(z)Zz 0 dx xN−1gx z(A.4) Denoting by y=x/z we have [ f⊗g(N) = Z1 0 dz zN−1f(z)Z1 0 dy yN−1g(y) = ˆ f(N) ˆg(N) (A.5) This is what is done for the DGLAP equations. It is numerically useful, as instead of working with convolutions, it is possible to work with products in the Mellin-space. Therefore, the products of the Mellin moments of the PDFs and the splitting functions need to be 133