Diffractive and non-diffractive wounded nucleons and final states in pA collisions

We review the state-of-the-art of Glauber-inspired models for estimating the distribution of the number of participating nucleons in pA and AA collisions. We argue that there is room for improvement in these model when it comes to the treatment of diffractive excitation processes, and present a new simple Glauber-like model where these processes are better taken into account. We also suggest a new way of using the number of participating, or wounded, nucleons to extrapolate event characteristics from pp collisions, and hence get an estimate of basic hadronic final-state properties in pA collisions, which may be used to extract possible nuclear effects. The new method is inspired by the Fritiof model, but based on the full, semi-hard multiparton interaction model of Pythia8.


Introduction
An important topic in the studies of the strong interaction is the understanding of the features of hot and dense nuclear matter. To correctly interpret signals for collective behaviour in high energy nucleus-nucleus collisions, it is necessary to have a realistic extrapolation of the dynamics in pp collisions. Here experiments on pA collisions have been regarded as an important intermediate step.
As an example refs. [1,2] have discussed the possibility to discriminate between the dynamics of the wounded nucleon model and that of the Color Glass Condensate formalism in pPb collisions at the LHC. An extrapolation of results from pp to pA and AA collisions is generally performed using the Glauber formalism [3,4]. This model is based on the eikonal approximation, where the interaction is driven by absorption into inelastic channels. Elastic scattering is then the shadow of absorption, and determined by the optical theorem. The projectile nucleon(s) are assumed to travel along straight lines and undergo multiple sub-collisions with nucleons in the target. The Glauber model has been commonly used in experiments at RHIC and LHC, e.g. to estimate the number of participant nucleons, N part , and the number of binary nucleon-nucleon collisions, N coll , as a function of centrality. A basic assumption JHEP10(2016)139 is then that one can compare a pA or an AA collision, at a certain centrality with, e.g., N part /2 or N coll times the corresponding result in pp collisions (for which N part = 2). A comparison with a fit to pp collision data, folded by the distribution in N part /2 or N coll , can then be used to investigate nuclear effects on various observables.
There are several problems related to such analyses, and in this paper we will concentrate on two of them: • Since the actual impact parameter is not a physical observable, the experiments typically select an observable, which is expected to be strongly correlated with the impact parameter (such as a forward energy or particle flow). This implies that the definition of centrality becomes detector dependent, which, among other problems, also implies difficulties when comparing experimental results with each other and with theoretical calculations.
• When the interaction is driven by absorption, shadow scattering (meaning diffraction) can contain elastic as well as single and double diffractive excitation. This is important since experiments at high energy colliders show, that diffractive excitation is a significant fraction of the total cross section, and not limited to low masses (see e.g. [5][6][7]). Thus the driving force in Glauber's formalism should be the absorptive, meaning the non-diffractive inelastic cross section, and not the total inelastic cross section.
In the following we will argue that the approximations normally used in this procedure are much too crude, and we will present a number of suggestions for how they can be improved, both in the way N part and N coll are calculated and the way pp event characteristics are extrapolated to get reference distributions. In both cases we will show that diffractive processes play an important role. In Glauber's original analysis only elastic scattering was taken into account, but it was early pointed out by Gribov [8], that diffractive excitation of the intermediate nucleons gives a significant contribution. However, problems encountered when taking diffractive excitation into account have implied, that this has frequently been neglected, also in recent applications (see e.g. the review by Miller et al. [4]). Thus the "black disk" approximation, and other simplifying treatments, are still frequently used in analyses of experimental results. 1 A way to include diffractive excitation in a Glauber analysis, using the Good-Walker formalism, was formulated by Heiselberg et al. [10]. It was further developed in several papers (see refs. [11][12][13][14] and further references in there) and is often called the "Glauber-Gribov" (GG) model. In the Good-Walker formalism [15], diffractive excitation is described as the result of fluctuations in the nucleon's partonic substructure. When used in impact parameter space, it has the advantage that saturation effects can easily be taken into account, which makes it particularly suited for applications in collisions with nuclei.
The "Glauber-Gribov" model has been applied both to data from RHIC and in recent analyses of data from the LHC, e.g. in refs. [16,17] However, although this formalism implies JHEP10(2016)139 a significant improvement of the data analyses, also in this formulation the treatment of diffractive excitation is simplified, as the full structure of single excitation of either the projectile or the target, and of double diffraction, is not taken into account. As we will show in this paper, this simplification causes important problems, and we will here present a very simple model which separates the fluctuations in the projectile and the target nucleons.
To guide us in our investigation of conventional Glauber models we use the dipsy Monte Carlo program [18][19][20], which is based on Mueller's dipole approach to BFKL evolution [21,22], but also includes important non-leading effects, saturation and confinement. It reproduces fairly well both total, elastic, and diffractive pp cross sections, and has also recently been applied to pA collisions [9]. The dipsy model gives a very detailed picture of correlations and fluctuations in the initial state of a nucleon, and by combining it with a simple geometrical picture of the distribution of nucleons in a nucleus in its ground state, we can build up an equally detailed picture of the initial states in pA and AA collisions. This allows us to gain new insights into the pros and cons of the approximations made in conventional Glauber Models.
The dipsy program is also able to produce fully exclusive hadronic final states in pp collisions, giving a reasonable description of minimum bias data from e.g. the LHC [20]. It could, in principle also be used to directly model final states in pA and AA, but due to some shortcomings, we will in this paper instead only use general features of these final states to motivate a revival of the old Fritiof model [23,24] with great similarities with the original "wounded nucleon" model [25]. (For a more recent update of the wounded nucleon model see ref. [26].) For energies up to (and including) those at fixed target experiments at CERN, the particle density at mid-rapidity in pp collisions is almost energy independent. For higher energies the density increases, and the p ⊥ distribution gets a tail to larger values. However, for minimum bias events with lower p ⊥ , the wounded nucleon model still works with the multiplicity scaling with the number of participating (wounded) nucleons, both at RHIC [27,28] and LHC [29]. For higher p ⊥ the distributions scale, however, better with the number of binary N N collisions, indicating the effect of hard parton-parton sub-collisions [17].
We will here argue that, due to the relatively flat distribution in rapidity of high-mass diffractive processes, absorbed and diffractively excited nucleons will contribute to the pA (and in principle also AA) final states in very similar ways, as wounded nucleons. We will also present preliminary results where we use our modified GG model to calculate the number distribution of wounded nucleons in pA, and from that construct hadronic final states by stacking diffractive excitation events, on top of a primary non-diffractive scattering, using Pythia8 with its semi-hard multi-parton interaction picture of hadronic collisions.
Although this remarkably simple picture gives very promising results, we find that there is a need for differentiating between diffractively and non-diffractively wounded nucleons. We will here be helped by the simple model mentioned above, in which fluctuations in the projectile and the target nucleon are treated separately. The model involves treating both the projectile and target as semi-transparent disks, separately fluctuating between two sizes according to a given probability. The radii, the transparency and the fluctuation probability is then adjusted to fit the non-diffractive nucleon-nucleon cross section, as well JHEP10(2016)139 as the elastic, single diffractive and double diffractive cross sections. Even though this is a rather crude model, it will allow us to investigate effects of the difference between diffractively and non-diffractively wounded nucleons.
We will begin this article by establishing in section 2 the framework we will use to describe high energy nucleon-nucleon scattering, with special emphasis on the Good-Walker formalism for diffractive excitation. In section 3 we will then use this framework to analyse the Glauber formalism in general and define the concept of a wounded target cross section. In section 4 we dissect the conventional Glauber models and the Glauber-Gribov model together with the dipsy model and present some comparisons of the resulting number distributions of wounded nucleons in pA. In section 5 we then go on to present our proposed model for constructing fully exclusive hadronic final states, and compare the procedure to recent results on particle distributions in pA collisions from the LHC, before we present conclusions and an outlook in section 6.
2 Dynamics of high energy pp scattering

Multiple sub-collisions and perturbative parton-parton interaction
As mentioned in the introduction, at energies up to those at fixed target experiments and the ISR at CERN, the pp cross sections and particle density, dn/dy are relatively independent of energy. For collisions with nuclei the wounded nucleon model works quite well [25], which formed the basis for the development of the Fritiof model [23]. This model worked very well within that energy range, but at higher energies it could not in a satisfactory way reproduce the development of a high p ⊥ tail caused by hard parton-parton interactions. Nevertheless the wounded nucleon model works well for minimum bias events even at LHC energies, if the rising rapidity plateau in pp collisions is taken into account, although the production of high p ⊥ particles appear to scale better with the number of N N collisions. These features may be interpreted as signals for dominance of soft interactions, and were the basis for the development of the Fritiof model [23]. This model worked very well within that energy range, but at higher energies, available atpp colliders at CERN and Fermilab, the effects of (multiple) hard parton-parton sub-collisions became increasingly important, and not so easily incorporated in the Fritiof model. Today high energy collisions (above √ s ∼ 100 GeV) are more often described as the result of multiple partonic sub-collisions, described by perturbative QCD. This picture was early proposed by Sjöstrand and van Zijl [30], and is implemented in the Pythia8 event generator [31]. This picture has also been applied in other generators such as herwig [32], Sherpa [33], dipsy [18,20], and others. The dominance of perturbative effects can here be understood from the suppression of low-p ⊥ partons due to saturation, as expressed e.g. in the Color Glass Condensate formalism [34].

The eikonal approximation
The large cross sections in hadronic collisions imply that unitarity constraints are important, and the elastic amplitude has to satisfy the optical theorem, which with convenient JHEP10(2016)139 normalisation reads Here the sum runs over all inelastic channels j. In high energy pp collisions the real part of the elastic amplitude is small, which indicates that the interaction is dominated by absorption into inelastic channels, with elastic scattering formed as the diffractive shadow of this absorption. This diffractive scattering is dominated by small p ⊥ , and the scattered proton continues essentially along its initial direction. At high energies and small transverse momenta, multiple scattering corresponds to a convolution in transverse momentum space, which is represented by a product in transverse coordinate space. This implies that diffraction and rescattering is more easily described in impact parameter space. In a situation where all inelastic channels correspond to absorption (meaning no diffractive excitation), the optical theorem in eq. (2.1) implies that the elastic amplitude in impact parameter space is given by Here P abs (b) = j |A j (b)| 2 represents the probability for absorption into inelastic channels.
If the absorption probability in the Born approximation is given by 2F (b), then unitarity is restored by rescattering effects, which exponentiates in b-space and give the eikonal approximation: 3) To simplify the notation we introduce the nearly real amplitude T = −iA el = 1 − S. The relation in eq. (2.2) then gives S(b) = e −F (b) and T (b) = 1 − e −F (b) . The optical theorem then gives We note that the possibility of diffractive excitation is not included here. Therefore the absorptive cross section in eq. (2.3) is the same as the inelastic cross section. How to include diffractive excitation and its relation to fluctuations will be discussed below in section 2.3. We then also note that diffractive excitation is very sensitive to saturation effects, as the fluctuations go to zero when saturation drives the interaction towards the black limit.
That rescattering exponentiates in transverse coordinate space also makes this formulation suitable for generalisations to collisions with nuclei.

Dipole models in transverse coordinate space
In this paper we will use our implementation of Mueller's dipole model, called dipsy, in order to have a model which gives a realistic picture of correlations and fluctuations in the colliding nucleons. In this way we can evaluate to what extent Glauber-like models  Figure 1. A colour dipole cascade in transverse coordinate space. A dipole can radiate a gluon. The gluon carries away colour, which implies that the dipole is split in two dipoles, which in the large N c limit radiate further gluons independently.  Figure 2. In a collision between two dipole cascades, two dipoles can interact via gluon exchange. As the exchanged gluon carries colour, the two dipole chains become recoupled. are able to take such effects into account. The dipsy model has been described in a series of papers [18][19][20] and we will here only give a very brief description. Mueller's dipole model [21,22] is a formulation of LL BFKL evolution in impact parameter space. A colour charge is always screened by an accompanying anti-charge. A charge-anti-charge pair can emit bremsstrahlung gluons in the same way as an electric dipole, with a probability per unit rapidity for a dipole (r 0 , r 1 ) to emit a gluon in the point r 2 , given by (cf. (2.5) The important difference from electro-magnetism is that the emitted gluon carries away colour, which implies that the dipole splits in two dipoles. These dipoles can then emit further gluons in a cascade, producing a chain of dipoles as illustrated in figure 1. When two such chains, accelerated in opposite directions, meet, they can interact via gluon exchange. This implies exchange of colour, and thus a reconnection of the chains as shown in figure 2.
The elastic scattering amplitude for gluon exchange is in the Born approximation given by BFKL evolution is a stochastic process, and many sub-collisions may occur independently. Summing over all possible pairs gives the total Born amplitude The unitarised amplitude then becomes

The Lund dipole model dipsy
The dipsy model [18][19][20] is a generalisation of Mueller's cascade, which includes a set of corrections: • Important non-leading effects in BFKL evolution. Most essential are those related to energy conservation and running α s .
• Saturation from Pomeron loops in the evolution. Dipoles with identical colours form colour quadrupoles, which give Pomeron loops in the evolution. These are not included in Mueller's model or in the BK equation.
• Confinement via a gluon mass satisfies t-channel unitarity.
• It can be applied to collisions between electrons, protons, and nuclei.
Some results for pp total and elastic cross sections are shown in refs. [35,36]. We note that there is no input structure functions in the model; the gluon distributions are generated within the model. We also note that the elastic cross section goes to zero in the dip of the t-distribution, as the real part of the amplitude is neglected.

Diffractive excitation and the Good-Walker formalism
In his analysis of the Glauber formalism, Gribov considered low mass excitation in the resonance region, but experiments at high energy colliders have shown, that diffractive excitation is not limited to low masses, and that high mass diffraction is a significant fraction of the pp cross section also at high energies (see e.g. [5][6][7]). Diffractive excitation is often described within the Mueller-Regge formalism [37], where high-mass diffraction is given by a triple-Pomeron diagram. Saturation effects imply, however, that complicated diagrams with Pomeron loops have to be included, which leads to complicated resummation schemes, see e.g. refs. [38][39][40]. These effects make the application in Glauber calculations quite difficult.
High mass diffraction can also be described, within the Good-Walker formalism [15], as the result of fluctuations in the nucleon's partonic substructure. Diffractive excitation is here obtained when the projectile is a linear combination of states with different absorption probabilities. This formalism was first applied to pp collisions by Miettinen and Pumplin [41], and later within the formalism for QCD cascades by Hatta et al. [42] and by Avsar and coworkers [36,43]. When used in impact parameter space, this formulation has the advantage that saturation effects can easily be taken into account, and this feature makes it particularly suited in applications for collisions with nuclei. (For a BFKL Pomeron, the Good-Walker and the Mueller-Regge formalisms describe the same physics, seen from different sides [44].) As an illustration of the Good-Walker mechanism, we can study a photon in an optically active medium. For a photon beam passing a black absorber, the waves around the JHEP10(2016)139 absorber are scattered elastically, within a narrow forward cone. In the optically active medium, right-handed and left-handed photons move with different velocities, meaning that they propagate as particles with different mass. Study a beam of right-handed photons hitting a polarised target, which absorbs photons polarised in the x-direction. The diffractively scattered beam is then a mixture of right-and left-handed photons. If the right-handed photons have lower mass, this means that the diffractive beam contains also photons excited to a state with higher mass.

A projectile with substructure colliding with a structureless target
For a projectile with a substructure, the mass eigenstates can differ from the eigenstates of diffraction. Call the diffractive eigenstates Φ k , with elastic scattering amplitudes T k . The mass eigenstates Ψ i are linear combinations of the states Φ k : (2.10) The elastic scattering amplitude is given by 11) and the elastic cross section The amplitude for diffractive transition to the mass eigenstate Ψ k is given by which gives a total diffractive cross section (including elastic scattering) (2.14) Consequently the cross section for diffractive excitation is given by the fluctuations: We note in particular that in this case the absorptive cross section equals the inelastic non-diffractive cross section. Averaging over different eigenstates eq. (2.3) gives l be the diffractive eigenstates for the projectile and the target respectively, and T kl the corresponding eigenvalue. (We here make the assumption that the set of eigenstates for the projectile are the same, for all possible target states. This assumption is also made in the dipsy model discussed above.) The total diffractive cross section, including elastic scattering, is then obtained by taking the average of T 2 kl over all possible states for the projectile and the target. Subtracting the elastic scattering then gives the total cross section for diffractive excitation: Here the subscripts p and t denote averages over the projectile and target respectively. Taking the average over target states before squaring gives the probability for an elastic interaction for the target. Subtracting single diffraction of the projectile and the target from the total in eq. (2.17) will finally give the double diffraction. Thus we get the following relations: where σ Dp and σ Dt is single diffractive excitation of the projectile and target respectively and σ DD is double diffractive excitation. Also here the absorptive cross section, which will be important in the following discussion of the Glauber model, corresponds to the non-diffractive inelastic cross section:

Diffractive eigenstates at high energies
In the early work by Miettinen and Pumplin [41], the authors suggested that the diffractive eigenstates correspond to different geometrical configurations of the valence quarks, as a result of their relative motion within a hadron. At higher energies the proton's partonic structure is dominated by gluons. The BFKL evolution is a stochastic process, and it is then natural to interpret the perturbative parton cascades as the diffractive eigenstates (which may also depend on the positions of the emitting valence partons). This was the assumption in the work by Hatta et al. [42] and in the dipsy model. Within the dipsy model, based on BFKL dynamics, it was possible to obtain a fair description of both the experimental cross section [36,43] and final state properties [45] for diffractive JHEP10(2016)139 excitation. In the GG model two sources to fluctuations are considered; first fluctuations in the geometric distribution of valence quarks, and secondly fluctuations in the emitted gluon cascades, called colour fluctuations or flickering. In ref. [13] it was concluded that the latter is expected to dominate at high energies. We here also note that at very high energies, when saturation drives the interaction towards the black limit, the fluctuations go to zero. This implies that diffractive excitation is largest in peripheral collisions, where saturation is less effective. This is true both for pp collisions and collisions with nuclei. (Although diffractive excitation of the projectile is almost zero in central pA collisions, this is not the case for nucleons in the target.) 3 Glauber formalism for collisions with nuclei

General formalism
High energy nuclear collisions are usually analysed within the Glauber formalism [3] (for a more recent overview see [4]). In this formalism, target nucleons are treated as independent, and any interaction between them is neglected. 2 The projectile nucleon(s) travel along straight lines, and undergo multiple diffractive sub-collisions with small transverse momenta. As mentioned in the introduction, multiple scattering, which in transverse momentum space corresponds to a convolution of the scattering S-matrices, corresponds to a product in transverse coordinate space. Thus the matrices S (pNν ) , for the encounters of the proton with the different nucleons in the target nucleus, factorise: (3.1) We denote the impact parameters for the projectile and for the different nucleons in the target nucleus by b and b ν respectively, and defineb ν ≡ b − b ν . Using the notation in eq. (2.4), we then get the following elastic scattering amplitude for a proton hitting a nucleus with A nucleons: 2 In the dipsy model gluons with the same colour can interfere, also when they come from different nucleons. This so-called inter-nucleon swing mechanism was shown [9] to have noticeable effects in photonnucleus collisions, but in pA, especially for heavy nuclei, the effects were less that 5%. We have therefore chosen to ignore such effects in this paper, but may return to the issue in a future publication.

JHEP10(2016)139
The inelastic cross section (now equal to the absorptive) would be equal to the difference between these two, in accordance with eq. (2.3).
Fluctuations in the pp interaction are discussed in the following subsection. Fluctuations and correlations in the nucleon distribution within the nucleus are difficult to treat analytically, and therefore most easily studied by means of a Monte Carlo, as discussed further in sections 3.4, 4 and 5 below. Valuable physical insight can, however, be gained in an approximation where all correlations between target nucleons are neglected. Such an approximation, called the optical limit, is discussed in section 3.5.

Gribov corrections. Fluctuations in the pp interaction
Gribov pointed out that the original Glauber model gets significant corrections due to possible diffractive excitation. In the literature it is, however, common to take only diffractive excitation of the projectile into account, disregarding possible excitation of the target nucleons. In this section we will develop the formalism to account for excitations of nucleons in both projectile and target. We will then see that in many cases fluctuations in the target nucleons will average out, while in other cases they may give important effects. (Fluctuations in both projectile and target will, however, be even more essential in nucleus-nucleus collisions, which we plan to discuss in a future publication.)

Total and elastic cross sections
When the nucleons can be in different diffractive eigenstates, the amplitudes T (pNν ) in eq. (3.2) are matrices T (pNν ) k,lν , depending on the states k for the projectile and l ν for the target nucleon ν. The elastic pA amplitude, T (pA) (b) , can then still be calculated from eq. (3.2), by averaging over all values for k and l ν , with ν running from 1 to A. Thus When evaluating the averages in these equations, it is essential that the projectile proton stays in the same diffractive eigenstate, Φ k , throughout the whole passage through the target nucleus, while the states, Φ lν , for the nucleons in the target nucleus are uncorrelated from each other. This implies that for a fixed projectile state k, the average of the S-matrix over different states, l ν , for the target nucleons factorise in eq. (3.1) or (3.5). Thus we have Here · · · k ( · · · lν ) denotes average over projectile (target nucleon) substructures k (l ν ), while · · · bν denotes average over the target nucleon positions b ν , an as beforeb ν ≡ b−b ν . We introduce the following notation for the average of the pp amplitude over target states:

JHEP10(2016)139
The pA amplitude can then be written in the form where the average is taken over the target nucleon positions b ν and the projectile states, k. The total and elastic cross sections in eqs. (3.5) and (3.6) are finally obtained from eq. (2.4). We want here to emphasise that these expressions only contain the first moment with respect to the fluctuations in the target states, l ν , but also all higher moments of the fluctuations in the projectile states, k.
To evaluate the b-integrated cross sections, we must know both the distribution of the (correlated) nucleon positions, b ν , and the b-dependence of the pp amplitude T The distribution of nucleon positions is normally handled by a Monte Carlo, as will be discussed in section 3.4. When fluctuations and diffractive excitation was neglected in section 3.1, the b-dependence of T (pp) (b) could be well approximated by a Gaussian distribution C exp (−b 2 /2B), corresponding to an exponential elastic cross section dσ/dt ∝ exp (Bt). With fluctuations it is necessary to take the unitarity constraint T ≤ 1 into account, which implies that a large cross section must be associated with a wider distribution. One should then check that after averaging the differential elastic cross section reproduces the observed slope. 3

Specification of "wounded" nucleons
The notion of "wounded" nucleons was introduced by Bia las, Bleszyński, and Czyż in 1976 [25], based on the idea that inelastic pA or AA collisions can be described as a sum of independent contributions from the different participating nucleons. 4 In ref. [25] diffractive excitation was neglected, and thus "wounded nucleons" was identical to inelastically interacting nucleons. 5 Although the importance of diffractive excitation was pointed out by Gribov already in 1968 [8], it has, as far as we know, never been discussed whether or not diffractively excited nucleons should be regarded as wounded. These nucleons contribute to the inelastic, but not to the absorptive cross section, as defined in eq. (2.19).
Diffractive excitation is usually fitted to a distribution proportional to dM 2 is the intercept of the Pomeron trajectory, estimated to around 1.2 from the HERA structure functions at small x. More complicated diagrams tend, however, to reduce ǫ. (In ref. [40] it is shown JHEP10(2016)139 that the largest correction is a four-Pomeron diagram, which gives a contribution with ǫ = 0.) Fits to LHC data [6,7] give ǫ ≈ 0.1, but with rather large uncertainties.
If ǫ is small, diffractively excited target nucleons can contribute to particle production both in the forward and in the central region. If ǫ instead is large, diffraction would contribute mainly close to the nucleus fragmentation region. For ǫ ≈ 0.1, the experimentally favoured value, the contribution in the central region would be suppressed by a factor exp(−0.1 · ∆η) ∼ 1/2 for pPb collisions at LHC. We conclude that the definition of wounded nucleons should depend critically upon both the experimental observable studied in a certain analyses, and upon the still uncertain M X -dependence of diffractive excitation at LHC energies. (In section 5.1 we will show that a simple model, assuming similar contributions from absorbed and diffractively excited nucleons actually quite successfully describes the final state in pPb collisions at LHC.) Below we present first results for the absorbed, non-diffractive, nucleons, followed by results when diffractively excited nucleons are included.

Wounded nucleon cross sections
Absorptive cross section. We first assume that wounded nucleons correspond to nucleons absorbed via gluon exchange, which for large values of ǫ would be relevant for observables in the central region, away from the nucleus fragmentation region. Due to the relation T = 1 − S, the absorptive cross section in eq. (2.19) can also be written dσ abs /d 2 b = 1 − S 2 . We here note that, as the S-matrix factorises in the elastic amplitude in eqs. (3.1) and (3.5), this is also the case for S 2 . This implies that (3.10) In analogy with eq. (3.8) for σ tot , also here, when taking the average over the target states l ν , the factors in the product depend only on the projectile state k and the positionsb ν . We here introduce the notation This quantity represents the probability that nucleon ν is absorbed by a projectile in state k.
Averaging over all values for k and b ν , it gives the total pA absorptive, meaning inelastic non-diffractive, cross section This expression equals the probability that at least one target nucleon is absorbed.
Cross section including diffractively excited target nucleons. We now discuss the situation when also diffractively excited target nucleons should be counted as wounded.
(The case with an excited projectile proton is discussed below.) The probability for a JHEP10(2016)139 nucleon, ν, in the nucleus to be diffractively excited is obtained from eq. (2.18) by adding single and double diffraction: Adding the absorptive cross section in eq. (2.19) we obtain the total probability that a target nucleon, ν, is excited or broken up by either diffraction or absorption, and we will call such nucleons inclusively wounded (w inc ), as opposed to absorptively wounded (w abs ).
In analogy with eq. (3.11) we define W which gives the probability that the target nucleon ν is either absorbed or diffractively excited, by a projectile in state k. Thus, if these target nucleons are counted as wounded, the cross section is also given by eq. (3.12), when W . We note that the expression for the wounded nucleon cross section resembles the total one in eqs. (3.9) and (3.6), with T is known including its b-dependence. This is not the case for W Elastically scattered projectile protons. We should note that the probabilities given above include events, where the projectile is elastically scattered, and thus not regarded as a wounded nucleon. The probability for this to happen in an event with diffractively excited target nucleons, is given by the relation ( · · · p and · · · t denote averages over projectile and target states respectively) (3. 16) In case these events do not contribute to the observable under study, this contribution should thus be removed. For a large target nucleus, this is generally a small contribution.

Wounded nucleon multiplicity
In the following we let W k denote either W , depending upon whether or not diffractively excited target nucleons should be counted as wounded.

JHEP10(2016)139
Average number of wounded nucleons. As W k (b ν ) denotes the probability that target nucleon ν is wounded, the average number of wounded nucleons in the target is then (for fixed b) given by ν W k (b ν ) k,bν , obtained by summing over target nucleons ν, and averaging also over projectile states k and all target positions b ν . Averaging over impact parameters, b, is only meaningful, when we calculate the average number of wounded target nucleons per event with at least one wounded nucleon, which we denote N t w . This is obtained by dividing by the probability in eq. (3.12). Integrating over b, weighting by the same absorptive probability, and normalising by the total absorptive cross section (also integrated over b) we get Note that the total number of wounded nucleons is given by N w = N t w + 1, as the projectile proton should be added, provided the projectile proton is not elastically scattered (in which case all wounded target nucleons have to be diffractively excited).
Multiplicity distribution for wounded nucleons. It is also possible to calculate the probability distribution in the number of wounded target nucleons N t w . For fixed projectile states k and target nucleon positionsb ν , the probability for target nucleon ν to be wounded, or not wounded, is W k (b ν ) and 1 − W k (b ν ) respectively. For fixed k the probability distribution in the number of absorbed target nucleons is then given by Here the sum goes over all subsets C N t w of t w wounded target nucleons, and C N t w is the set of the remaining A − N t w target nucleons, which thus are not wounded. The states of the target nucleons can be assumed to be uncorrelated, and the averages could therefore be taken separately, as in eq. (3.11). The state k and positions b ν or b µ give, however, correlations between the different factors, and these averages must be taken after the multiplication, which gives the result The distribution in eq. (3.19) includes the possibility for N t w = 0. As for the average number of wounded nucleons above, to get the normalised multiplicity distribution for events, with N t w ≥ 1, we should divide by the probability in eq. (3.12). The final distribution is then obtained by integrating over b, with a weight given by the same absorption probability. This gives the result where dP/d N t w (b) and dσ pA w (b)/d 2 b are the expressions in eqs. (3.19)) and (3.12).

JHEP10(2016)139
We want here to emphasise that the quantity W (w abs ) k contains the average of the square of the amplitude T (pp) , and is therefore not simply determined from the average l , which appears in the expression for the total and elastic cross sections in eqs. (3.9) and (3.6). This contrasts to the situation for inclusively wounded nucleons, where W

Nucleus geometry and quasi-elastic scattering
In a real nucleus the nucleons are subject to forces with a hard repulsive core, and their different points r ν are therefore not uncorrelated. In Glauber's original papers this correlation was neglected, and this approximation is discussed in the subsequent section.
In addition to the suppression of nucleons at small separations, the geometrical structure will fluctuate from event to event. These fluctuations are not only a computational problem, but have also physical consequences. Just as fluctuations in the nucleon substructure can induce diffractive excitation of the nucleon, fluctuations in the nucleus substructure induces diffractive excitation of the nucleus. If the projectile is elastically scattered these events are called quasi-elastic. The fluctuations in the target nucleon positions are also directly reproduced by the Monte Carlo programs mentioned above, and within the Good-Walker formalism the quasi-elastic cross section, σ el * , is given by (cf. eq. (2.18)): (3.21) The average over the target states here includes averaging over all geometric distributions of nucleons in the nucleus, and all partonic states of these nucleons. Note that this expression includes the elastic proton-nucleus scattering (given by T 2 p,t ). Some results for quasielastic pPb collisions are presented in ref. [9,46].

Optical limit -Uncorrelated nucleons and large nucleus approximations
Even though the averages in eqs. (3.5) and (3.12) factorise, they are still complicated by the fact that all factors S (pNν ) are different, due to the different values for the impact parameters. It is interesting to study simplifying approximations, assuming uncorrelated nucleon positions and large nuclei. This is generally called the optical limit. It was used by Glauber in his initial study [3], and is also described in the review by Miller et al. [4], for a situation when diffractive excitation is neglected. We here discuss the modifications necessary when diffractive excitation is included, also separating single excitation of projectile and target, and double diffraction.

Uncorrelated nucleons
Neglecting the correlations between the nucleon positions in the target nucleus, the individual nucleons can be described by a smooth density A · ρ(b ν ) (normalised so that which enter the total pA cross section in eq. (3.5), are uncorrelated and give the same result, depending only on projectile state and impact parameter k and b: In the same way all factors W k (b ν ), entering the wounded nucleon cross sections in eqs. (3.11) and (3.15), give equal contributions: (3.23)

Large nucleus
If, in addition to the approximations in eqs. (3.22) and (3.23), the width of the nucleus (specified by ρ) is much larger than the extension of the pp interaction (specified by T (pp) ), further simplifications are possible. For the amplitude in eq. (3.22) we can integrate over b ν , and get the approximation We have here introduced the notation σ pp tot,k for the total cross section for a projectile proton in state k, averaged over all states for a target proton.
In the same way we get , and σ pp w,k is the corresponding pp cross section for a projectile in state k.

Total cross section
Inserting eq. (3.24) into eqs. (3.5) -(3.6) gives the total cross section for a projectile in state k hitting a nucleus: The total pA cross section is then finally obtained by averaging over projectile states, k, and integrating over impact parameters, b: We note here in particular, that in this approximation the b-dependence of T (pp) k (b) is unimportant, and the result depends only on its integral σ pp tot,k /2. We also note that to calculate the elastic pA cross section ∼ d 2 b (T (b)) 2 , which has a steeper b-dependence, a knowledge about this dependence is also needed.
Proton-deuteron cross section

JHEP10(2016)139
Neglecting fluctuations, eqs. (3.26) and (3.27) would give the simpler result For the special case with a deuteron target we then get the result 6 and with the estimate For a non-fluctuating amplitude, the optical theorem gives a direct connection between the total and elastic cross sections. As the integral over d 2 b gives the Fourier transform at q = 0, we have . (3.30) HereT (q) denotes the Fourier transform of the amplitude T (b). For a Gaussian interaction profile we get (σ where the slope B is a measure of the width of the interaction. As σ el is determined by the squared amplitude, the ratio σ el /σ tot will be larger for a strong interaction with a short range, than for a weaker interaction with a wider range. For the general case with fluctuating amplitudes, we can using the results in eq. (2.18), in an analogous way rewrite (σ Here σ pp Dp denotes the cross section for single diffractive excitation of the projectile proton (i.e. on one side only). For a fluctuating amplitude we then get instead of eq. (3.29) . (3.33) The negative term in eq. (3.33) represents a shadowing effect, which for a deuteron target has one contribution from the elastic proton-nucleon cross section, and another from diffractive excitation. Note in particular, that it is only single diffraction which enters, with an excited projectile but an elastically scattered target nucleon. (This would be particularly important in case of a photon or a pion projectile.) Larger target nuclei For a larger target higher moments, T (pp) n be determined from the total cross section and the cross section for diffractive excitation. They can be calculated if we know the full probability distribution, dP/d T (pp) t , for the pp amplitude averaged over target states, but for varying projectile states. 7 In addition also higher moments of the nucleus density, d 2 b ρ n (b), are needed.
We also note here that the factorisation feature in eq. (3.2) is not realised in AA collisions. This implies that also in the optical limit, the AA-results cannot be directly expressed in terms of the moments T (pp) n t p .

Wounded nucleon cross sections
Also for cross sections corresponding to wounded (absorptively or inclusively) nucleons, approximations analogous to eqs.
The average in eq. (3.34) includes averages of all possible powers (σ pp w,k ) n k . For n = 1 this is just equal to the pp cross section σ pp w for (with w denoting either absorptively or inclusively wounded), but for higher moments a knowledge of the full probability distribution for σ pp w,k is needed, in analogy with eq. (3.27) for the total pA cross section. Note, however, that a similar relation is not satisfied for the elastic or total inelastic cross sections, σ el and σ in = σ tot − σ el , which as seen in eq. (2.18) contain the average over projectile states k before squaring.

Average number of wounded nucleons
In eq. (3.34) ρ(b)σ pp w,k represents the probability that a specific target nucleon is wounded, in a collision with a projectile in state k at an impact parameter b. In the optical limit this probability is the same for all A target nucleons. Averaging over projectile states k then gives the average number of wounded target nucleons for an encounter at this b-value. Dividing by the probability for a "wounded" event, we get the average number of wounded target nucleons per wounded event for this b: Normalising by the probability for absorption in eq. (3.34), and integrating over b with a weight given by the same probability, then gives with dσ pA w /d 2 b given by eq. (3.34). As noted above, this needs knowledge of the full probability distribution for σ pp w,k .

JHEP10(2016)139
3.5.6 Multiplicity distribution for wounded nucleons As in section 3.3, when calculating the full distribution in N t w (b), it is important to take the average over projectile states k after multiplication of the different nucleon absorption probabilities, which gives Similar to the general result in eq. (3.19), this expression includes the probability for zero target participants. Normalising by the probability for absorption in eq. (3.34), and integrating over b with a weight given by the same probability, gives finally .

Models for pp scattering used in Glauber calculations
As mentioned in section 3.4, most analyses today use a Monte Carlo simulation to generate a realistic distribution of nucleons within the nucleus, including fluctuations which cause quasi-elastic scattering of the nucleus [9,46] as well as initial state anisotropies (e.g. [47]).
In contrast most Glauber Monte Carlos use a rather simple model for the pp interaction. In this section we discuss some models which have been used in analyses of experimental data. We will also comment on the pros and cons, when these models are applied to pA collisions.
In the optical approximation, where the extension of the nucleus is much larger than the range of the pp interaction, the results for pA collisions can be expressed in terms of integrated pp amplitudes, without knowledge of their respective impact parameter dependence (see eq. (3.24)). It is therefore most essential to use a model, where the integrated pp cross sections are well reproduced. Note, however, that although the total pA cross section is most sensitive to the integrated total pp cross section, the b-dependence is very important for the ratio between the elastic and total cross sections (see eq. (3.31)). This feature naturally also affects the ratio between the inelastic and the total cross sections.
As mentioned in the introduction, the problems encountered when taking fluctuations and diffractive excitation of the nucleons properly into account in the Glauber model, have implied that these effects are neglected or severely approximated in many applications, see e.g. ref. [4]. However also in models which do include fluctuations, as far as we know no published analysis uses a model which can separate single excitation of the projectile from that of the target, and from double excitation. This is a problem as the various pA cross sections in section 3 contain powers of pp amplitudes averaged in different ways over projectile and target fluctuation.
We first discuss some simple models determined by just a few parameters, and then the more ambitious approach by Strikman and coworkers, using a continuous distribution for the fluctuations. (i) Black disk model. The simplest approximation is the "black disk model" with a fixed radius. Here diffractive excitation is completely neglected, and the target in a nucleon-nucleon collision acts as a black absorber. The projectile nucleon travels along a straight line, and interacts inelastically if the transverse distance to a nucleon in the target is smaller than a distance R, which gives This results in the following cross sections: Here σ D denotes the cross section for diffractive excitation. (See eq. (2.4) with F = ∞.) This is in clear contrast to the experimental result σ el ≈ σ tot /4 and the total diffractive excitation of the same order of magnitude as σ el . This again illustrates how a short range amplitude gives a large σ el /σ tot ratio. The radius can therefore be adjusted to reproduce the experimental value for one of these three cross sections, at the cost of not reproducing the other two.
As discussed in ref. [9], choosing to reproduce σ pp tot , the simple black-disk result for pPb collisions agrees rather well with the dipsy model for σ pPb tot , but not so well for σ pPb el or σ pPb in . Similarly adjusting R to reproduce σ pp in or σ pp abs gives results which agree with dipsy for the corresponding pPb cross section, but not for the other.
The black disk model is implemented in many Monte Carlos, e.g. in the PHOBOS Monte Carlo [48,49]. It is also used in refs. [46,47] where the authors study fluctuations in the distribution of nucleons within the nucleus, but do not address the fluctuations in the pp interaction.
(ii) Grey disk and Gaussian profile. Also other shapes for a non-fluctuating pp interaction have been used in the literature [4,24,50]. The simplest example is a fixed semi-transparent "grey disk", with opacity given by the parameter α: which gives σ el : σ tot : σ in = α : 2 : 2 − α.
These models contain two parameters (with α ≤ 1 to satisfy the unitarity constraint T ≤ 1), and it is therefore possible to fit e.g. the total and the elastic cross sections, with the inelastic (non-diffractive) cross section given by the difference between these two. The lower ratio σ el /σ tot is a consequence of the wider interaction range. We note, however,

JHEP10(2016)139
that even if typical events are well reproduced it is often interesting to study rare events in the tail of a distribution. As an example the tail of the pp amplitude out to large b-values may be important for the probability to produce rare events with many pN sub-collisions at large impact separation. The Gaussian profile may e.g. thus give a larger tail than the gray disk, also when they give very similar averages.

Models including fluctuations
To account for diffractive excitation, we must allow the pp amplitude to fluctuate. Models used in the literature do, however, not separate fluctuations in the projectile and the target. From eqs. (3.9) and (3.6) we see that if the amplitude is adjusted to reproduce the amplitude averaged over target states, T (pp) (b) t , then the correct result for the pA total cross section will be obtained. The fluctuations included in the model should then only describe fluctuations in the projectile state, and should thus reproduce the cross section for single excitation of the projectile. Such a model will, however, not reproduce cross sections for absorptively wounded nucleons properly, as will be discussed further below.
(iii) Fluctuating grey disk. The simplest model accounting for diffractive excitation is the fluctuating "grey disk model". Here it is assumed that within a radius R the projectile is absorbed with probability a, with 0 < a < 1. This implies that T (b) 2 = T (b) , and the resulting pp cross sections are here The two parameters R and a can now be adjusted to reproduce e.g. the total and the elastic pp cross sections. At LHC this would give a ≈ 1/2. The cross section for diffractive excitation should here be interpreted as representing only the single excitation of the projectile, while target excitation is part of the absorptive cross section. With a = 1/2 this is quite an overestimate. It corresponds rather to the total diffractive excitation, which implies that the results for σ abs /σ tot is close to the experimental value. The relation between the absorptive and diffractive cross section, which together make up the inelastic cross section, is also fixed in this model. The agreement of the fluctuating gray disk with dipsy results for pPb collisions are not superior to those of the black disk model [9].
(iv) Fluctuating Gaussian profile. Here the profile in eq. (4.4) gives the probability for absorption. Thus T = 1 with probability α exp(−b 2 /2B) while T = 0 with probability 1 − α exp(−b 2 /2B). As for the fluctuating gray disk this implies that T (b) 2 = T (b) . This does not change the total and elastic cross sections, but it splits the inelastic one JHEP10(2016)139 into relative fractions a non-diffractive (absorptive) and a diffractive part, with the result σ tot : σ el : σ D : σ abs = 4 : α : 2 − α : 2. For α ≈ 1 this gives the same result as the fluctuating gray disk. (Although this implies that the results will be very similar in the optical limit, it does not mean that the results are identical for a more realistic nucleus. This will be particularly true for the tail at very large numbers of wounded nucleons.) This model is also an option in the PHENIX Monte Carlo.
(v) Fluctuating black disk model. It is possible to let the radius of the black disk fluctuate. As for the two previous model, the fact that T (b) is either 1 or 0 implies that T (b) 2 = T (b) , which gives σ abs = σ tot /2. Such a fluctuating black disk model has sometimes been used in connection with the GG model described in section 4.2, and will be further described below.
(vi) A new simple model allowing for separate projectile and target excitations.
The main reason neither of the above models are able to properly take into account the diffractive aspects of nucleon collisions, is that the fluctuations in the cross sections are not treated in terms of fluctuations in the projectile and target separately. Interpreting the amplitude T as the average over target states, which as mentioned above can give a correct total cross section, excitation of target nucleons will not be separated from the absorptive cross section.
To redeem this we have constructed a new model, which in some sense is the minimal possible extension needed to reproduce all relevant semi-inclusive cross sections. The basis of the model is having fluctuating sizes of the colliding nucleons. With some probability, c, a nucleon with radius r 1 , can fluctuate into a larger radius r 2 . This will then give us the elastic amplitude for a projectile with radius R p colliding with a target with radius, R t , Here α is again an opacity parameter between zero and one, which together with c, r 1 and r 2 gives us four parameters which can be adjusted to reproduce the relevant nucleonnucleon cross sections σ abs , σ el , σ Dp = σ Dt and σ DD . Below we will refer to this model as 2×2-disk.

The approach by Strikman and coworkers
An ambitious approach to describe fluctuations in pp scattering, for use in the Glauber model, was presented in refs. [10,11]. This model has been further extended in several papers by Alvioli, Strikman and coworkers; for a general overview see ref. [12] and further references in there. Recent studies, with applications to the LHC, discuss effects of colour fluctuations (or flickering) [13], and evidence for x-dependent proton colour fluctuations [14]. The model does not take into account the possibility of separate excitations of the projectile and the target, and the fluctuations in the target are not considered. In this section we will also discuss how the model can be modified to take the target fluctuations into account. (Note that our amplitude T is in ref. [12] denoted Γ.)

pp total cross section
The basic feature of the model is a description of the fluctuations in the NN total cross section, as a smooth function, which has the form Here σ is regarded as the total pp cross section in a single event, with the probability distribution P tot (σ), while the observed total cross section σ tot is given by the average in eq. (4.8). 8 For the functional form in eq. (4.7), the average and the width of the distribution are related to (but not identical to) the parameters σ 0 and Ω, while ρ is a normalisation constant. In eqs. (3.6) and (3.9) we see that the total pA cross section contains all possible moments with respect to the fluctuations in the projectile state, but only the average (the first moment) with respect to the fluctuations in the target state. Thus, although target fluctuations are not considered explicitely, we conclude that if σ is interpreted as the average over target states and only the average over projectile states is described by the distribution in eq. (4.7), then the total pA cross sections will (in the large nucleus approximation) be determined in terms of all possible moments σ N , obtained from the distribution P tot (σ), and the average in eq. (4.8) will correctly give the total pp cross section. With this interpretation the width of the distribution can also be determined from eq. (3.32), which gives the second moment . (4.10) Here the first term in the parenthesis would give σ 2 corresponding to Glauber's result, while the second, determined by single excitation of the projectile, is the result of fluctuations in the projectile state. Eq. (4.10) has been used by Blaettel et al. [11] together with eq. (3.33) to estimate the width from shadowing in pd collisions at fixed target energies. They also estimated the width from diffractive excitation data at the CERNpp collider. With data from TOTEM [51][52][53] and ALICE [54] for elastic and single diffractive cross sections and elastic forward slope, supplemented by the assumption that the diffractive slope is approximately half the elastic (as is the case at 560 GeV [55]), we get for 7 TeV the estimated width σ 2 − σ 2 ≈ 0.4 σ . As mentioned above, the amplitude for larger nuclei the amplitude in eq. (3.26) contains also higher moments of the pp amplitude. Blaettel et al. estimated also the third moment, σ 3 from data for diffractive excitation in pd scattering, and JHEP10(2016)139 they also studied other analytic forms. Most recent applications use, however, the form in eq. (4.7), in which the higher moments are fixed by a determination of the width.

Elastic cross section
We use the notation to describe the b-dependence of the fluctuating cross section σ in eq. (4.9). This gives (4.12) As pointed out earlier, the relation between σ el and σ tot depends on the width of the interaction. Thus, although the elastic and total cross sections for fixed b are given by the same average over target fluctuations, the elastic cross section is not determined by the σ-distribution in eq. (4.7), unless it is supplemented by a knowledge of the b-dependence (for all values of σ). We here note that the distribution in eq. (4.7) has a tail out to large cross sections. The unitarity constraint T (b) < 1, or dσ tot /d 2 b < 2, therefore implies that a large value for σ must be associated with a wider b-distribution. The effect of different assumptions about the b-dependence will be discussed in section 4.2.4. This feature implies of course that also the inelastic cross section cannot be directly determined from eq. (4.7).

Wounded nucleon cross section
As discussed in section 3.3.2, the definition of a wounded nucleon may depend upon the specific observables under consideration. As pointed out earlier, in cases where the absorptive cross section is the most relevant, this is given by which cannot be determined without knowing how the separate fluctuations in the projectile and target result in single and double diffractive excitation. We see that in contrast to the expressions entering the total and elastic cross sections in eq. (4.12), this expression contains also the second moment with respect to the target fluctuations. In ref. [12] Alvioli and Strikman identify the differential wounded nucleon cross section with the total inelastic pp cross section (which includes diffractive excitation). In the hypothetical situation where the target did not fluctuate, after averaging over projectile fluctuations this also gives the absorptive (inelastic non-diffractive) cross section.

JHEP10(2016)139
However, if T is identified with the amplitude averaged over target states, as in eq. (4.11) (which gives the correct result for the total cross section), then we get instead (4.14) From eq. (2.18) we see that this corresponds exactly to the inclusively wounded nucleon cross section dσ w inc /d 2 b, where σ w inc now includes diffractively excited target nucleons: We here note that, although dσ w /d 2 b in eq. (4.14) contains the same average of T (pp) k,l (b) over target states, to integrate this expression over b we also need to know the b-dependence of T (pp) (b, σ) for all σ. We also note that the integral The distribution P w (σ w ) is consequently not easily related to the distribution in the total cross section P tot (σ tot ). Lacking a detailed description, Strikman et al. use an approximation assuming the proportional distribution which for the absorption probability would mean P abs (σ) ∝ P tot (σ/λ abs ), (4.16) where λ abs = σ abs /σ tot . This approximation may be less accurate, since for non-peripheral while for peripheral collisions with small T we have dP abs (b)/dT (b) = 2. Also, even if the analytic form in eq. (4.7) may give a satisfying result, there is no obvious reason why the same value of the width parameter Ω, should be applicable as the one determined from shadowing or diffractive excitation.
Monte Carlo implementations. The GG model has been implemented in Monte Carlo simulations in many applications to pA collisions, e.g. in ref. [12][13][14]. In experimental analyses it has been combined with earlier Monte Carlos, where the parameters in one of the simple models described in section 4.1 are allowed to vary according to eq. (4.7) (or using a scaled version as in eq. (4.16), but typically using the total inelastic cross section rather than the absorptive), in a way reproducing the total (or the inelastic) cross section respectively. The PHOBOS Monte Carlo [48] with a black disk with a variable radius, which is also used by e.g. ATLAS [16]. Fitting to the inelastic pp cross section here overestimates the number of absorbed nucleons. In ref. [13] it is argued that this is a small effect, as the cross section for diffractive excitation of the projectile proton is small in pA collisions. However, the cross section for target nucleon excitation is not small, and although the cross section for diffractive projectile excitations is small, it may have a significant effect on the tail of the wounded nucleon distribution at high multiplicities. These problems will be further discussed in section 5.

Impact-parameter profile
To investigate further, we need to make assumptions about the impact-parameter dependence, T (pp) (b, σ) in eq. where B is proportional to σ. The proportionality factor could then be fit together with the σ 0 and Ω parameters of eq. (4.7) to the total and elastic pp cross sections from eq. (4.12) and the inclusively wounded cross section in eqs. (4.14) and (4.15): In this case we find that the wounded cross section distribution can indeed be written as a simple scaling of the total, P w inc (σ) = P tot (σ/λ w inc ), however, the same would still not be true for P w abs . We also find that for the Gaussian profile, the unitarity constraint, T (pp) (b, σ) ≤ 1, gives a hard limit on σ tot − σ w inc = σ el + σ Dp < σ tot /4, which is not found experimentally. To proceed we therefore decided to choose a different form of the b distribution. What is used by ATLAS in e.g. [16] is a black disk approximation: We will instead use a semi-transparent disk with (cf. eq. (4.4)) where the unitarity constraint gives us σ el + σ Dp < σ tot /2, which can easily accommodate experimental data.

Conclusion on the GG formalism
We conclude that for the pA total cross sections, it is straight forward to use the GG formalism by Strikman et al., interpreting σ as the total cross section averaged over target states. The distribution P tot (σ) then describes the fluctuations in the projectile states. The average and the variance of P tot are given by eqs. (4.8) and (4.10). However, to get the elastic or inclusively wounded nucleon cross section (including target excitation), we also need to know the b-dependence If wounded nucleons are interpreted as absorbed nucleons, we also need to know (T (pp) (b)) 2 t . To estimate these quantities in a way consistent with eq. (2.18), we believe it is better to use a formalism which include individual excitation of both projectile and target.

Consequences of fluctuating pp cross section
Adding fluctuations to the pp cross section dramatically changes the distributions of the number of wounded nucleons. But since pp data only offers inclusive and semi-inclusive cross sections to compare models to, one is given little guidance to why one model works better than another. Although the dipsy model is less than perfect in reproducing experimental data, it includes those fluctuations in the nucleon wave function, which we argue are important when considering the number of participating nucleons in pA and AA collisions. Thus although it only works at high energies due to lack of quarks in the proton, the description of high-p ⊥ particles is poor, and generation of exclusive diffractive final states is difficult, we believe these deficiencies are less important when describing the fluctuations.

Comparison with dipsy
When comparing the GG results in eq. (4.18) with results from dipsy, we first look at the total cross section. As descussed above, we interpret the GG fluctuating total cross section σ in eq. (4.18) as describing fluctuations in the projectile, averaged over target states: The parameters in the GG distribution P tot (σ) can then be tuned to reproduce the corresponding distribution in dipsy, which is obtained by generating a large ensemble of targets for each projectile, and for each target calculate T at a large number of impact parameters. To get the corresponding results for the elastic and the "wounded" cross sections we have to make an assumption about the b-distribution of the amplitude T (b) t = T (pp) (b, σ), appearing in eq. (4.18). We here make the simple approximation in eq. (4.19), and calculate the cross section σ w inc from eq. (4.18).
Tuning the parameters in P tot (σ) to the dipsy results is now done fitting the cross sections σ tot , σ el , and σw inc using a χ 2 fit. The values obtained with dipsy are here assigned weights corresponding to the relative error one would expect from experiment (taken from the analysis in ref. [56]). The result of the fit is shown in the first line of table 1.
The result of the fit is compared with the dipsy results for the distributions P tot (σ) and P w inc (σ) in figure 3. We note here that the b-dependence assumed in eq. (4.19) implies that the distribution P (σ w inc ) is given by a scaled P tot -distribution P w inc (σ) ∝ P tot (σ/λ), where λ = σ w inc /σ tot .
It is clearly seen, that the high-σ tails of the dipsy distributions are not reproduced by the functional form for P tot in eq. (4.7). Since dipsy provides a picture of the fluctuations built upon a full dynamical model, it is reasonable to believe that the shape of the dipsy distributions are closer to reality than eq. (4.7). We therefore try a new parametrisation which makes it easier to obtain a large high-σ tail, namely a log-normal distribution:  The fit to the dipsy cross sections with the log-normal distribution is also shown in table 1. The corresponding distributions are shown in figure 3 for two different width parameters, labeled Ω = 0.25 and Ω = 0.33. We see that the larger value matches the dipsy distribution perfectly, while the lower value is close to the GG curve below the maximum but has a higher tail for larger σ.
We note, however, that for technical reasons the diffractive cross section in dipsy is calculated demanding a central rapidity gap, restricting the masses to M 2 X ≤ √ s · (1 GeV).
This implies that the fluctuations are somewhat underestimated. We therefore believe that the functional form is quite realistic, while the width is underestimated. Results obtained when tuning instead to the experimental cross sections are presented in the following subsection.

Comparison to data
We now repeat the same procedure as in the previous section, but with experimental results for the relevant cross sections. There is no experimental access to the distributions in cross section, but the integrated inclusive and semi-inclusive cross sections are measured, and we here use values from ref. [56], extrapolated to √ s NN = 5 TeV: σ tot = 93.2 ± 2.3 mb, σ el = 23.2 ± 1.2 mb and σ w inc = 63.0 ± 1.8 mb. (4.22) Note that the diffractive cross sections here have been extrapolated into unmeasured M X regions to the full 0 < M X < √ S interval. As mentioned above this was not done in for the DIPSY diffractive cross sections in section 4.3.1, where by construction a rapidity gap is required at mid rapidity. Hence we expect that the fluctuations for DIPSY are underestimated as compared to data. The parameter values obtained by minimising the χ 2 are listed in table 2.
In figure 4 we compare the fits of the two parametrisations of P w inc (σ) shown above. We see that the new parametrisation provides larger fluctuations in the high-σ tail, as expected. It should be noted that both fits reproduce the experimental cross sections  Table 2. GG parameter values obtained by fit to inclusive and semi-inclusive cross-section from data.  well within the experimental errors. For comparison we also show log-normal distribution P w inc (σ), when it's width and mean are fitted to dipsy by eye (denoted Ω = 0.33), which is significantly more narrow. We suspect that while the log-normal parametrisation probably gives a more realistic description of the high-σ fluctuations in the GG formalism, it is far from the whole story. The GG results presented here are obtained assuming that all fluctuations are ascribed to fluctuations in projectile size, as described in section 4.2.4. In dipsy, however, the cross section fluctuations arise from a combination of fluctuations in size and fluctuations in gluon density. We believe that updating the profile functions from simple disks or Gaussian distributions to more realistic ones, could provide a better handle on the parametrisations of the cross section from pp data, this will be investigated in a future publication. So far we have described a prescription which seems to both catch the necessary physics to calculate the inclusive wounded cross section, with all parameters being obtainable from pp data. We will now apply this to pA collisions.

Distributions of wounded nucleons
Using the considerations about fluctuations in the wounded cross section, we will now turn to generation of distributions of wounded nucleons. Normally, in inelastic, non-diffractive pA collisions, the number of wounded nucleons is always one plus the number of inelastic, non-diffractive NN interactions. In the following we will make the distinction between JHEP10(2016)139 r 1 r 2 α c 0.15 fm 1.07 fm 0.97 0.42 Table 3. Table of parameters of the 2×2-disk model fitted to pp data. diffractively and absorptively wounded nucleons. In order to avoid situations where the projectile should sometimes be counted twice as a wounded nucleon, we will solely talk about the number of wounded nucleons in the target, which we denote N t w . We note also that since the number of sub-collisions and the number of wounded nucleons are trivially connected, the question whether a specific observable scales better with wounded nucleons or with NN sub-collisions, is much more relevant for nucleus-nucleus collisions.

Inclusively wounded nucleons
We will describe the nucleus' transverse structure using a Woods-Saxon distribution in the GLISSANDO parameterisation [57,58], where the density is given by: where R is the nuclear radius, a is "skin width", and ρ 0 is the central density. The parameter w describes a possible non-constant density, but is zero for lead. The nucleons are generated with a hard core, which thus introduces short range correlations among the nucleons. As shown by Rybczynski and Broniowski [59], the correct two-particle correlation can be obtained if the nucleons are generated with a minimum distance equal to 2r core . Using r core = 0.45 fm and a skin width of a = 0.459, the radius of the Lead nucleus becomes R Pb = 6.406 according to the parameterisation in [58].
For each nucleus state we generate a random impact parameter wrt. the projectile proton and proceed to determine which nucleons will be wounded, following the previously outlined models.
In figure 5, we show the distribution in the number of inclusively wounded nucleons (using σ pp w inc = 63.0 mb) for: a black disk model without any pp cross section fluctuations; GG with parameters fitted to data in section 4.3.2; GG with P tot (σ) given by eq. (4.21), also fitted to data; and the new simplified model outlined in section 4.1 (here called 2×2disk) Fitting the latter to the cross sections in eq. (4.22) as well as to the double diffractive cross section σ DD = 3.2 mb, we obtain the parameters listed in table 3.
Looking at the individual distributions in figure 5, we see that all three inclusions of additional fluctuations in the cross section, significantly increases the tail of the distribution compared to the black disk. The 2×2-disk model has fewer fluctuations to very large N t w inc numbers, and the dip in the distribution around N t w inc = 10, also indicates that the fluctuations are too crude. The difference between GG with the original parametrisation and the log-normal distribution is visible in the tail above N t w inc ≈ 35, as expected. One would therefore expect only an effect in the central events.

Distinguishing between absorptively and diffractively wounded nucleons
In our interpretation of the GG model in section 4.2.3, it can be used to calculate the sum of absorptively and diffractively wounded nucleons. In the Monte Carlo one would, however, like to have an impact parameter dependent recipe for each sub-collision to decide whether or not a target nucleon is diffractively or absorptively wounded, when hit by a projectile in a definite state p. This amounts to calculating the ratio of the absorptive to the inclusively wounded cross sections for a given sub-collision, and compare it to a random number

(4.24)
For the 2×2-disk model this is done easily, as the above ratio reduces to: (4.25) The GG model on the other hand, implies averaging over target nucleon states, and provides thus no distinction. Instead we follow the 2×2-disk model to calculate the the conditional probability to be diffractively wounded, if a nucleon is already inclusively wounded in the GG model. This is: where the first term is a requirement that the two nucleons are separated by an amount such that a fluctuation in size is necessary to be wounded. In figure 6 we show distributions of N t w abs for the 2×2-disk model and for the corrected GG model, using both parametrisations of P (σ).

Modelling final states in pA collisions
In this section we will take the knowledge about distributions of wounded nucleons and investigate the consequences for final states in pA collisions. We will discuss a few views on modelling particle production in such collisions, all assuming that a full final state of a pA collision can be adequately modelled by stacking pp events on top of each other, here modelled using Pythia8. Following the introduction of the models, we will compare to data, both multiplicity as function of centrality, and inclusive p ⊥ spectra. Finally we will give an estimate of the theoretical uncertainties present at this early stage of the model.

Generating final states with Pythia8
The general methodology for generating final states, which will be pursued here, will have the following ingredients: • For each collision, a Glauber calculation is performed as outlined in section 4.4, setting up the nuclear geometry.
• The total number of inclusively wounded target nucleons is calculated, as well as the number of absorptively wounded targets, if the two differ in the considered approach.
• Sub-collisions are generated as pp collisions, according to two separate approaches, which will be outlined in the following.
• Each sub-collision is treated separately in terms of colour reconnection and hadronisation. Efforts to include cross talk between sub-collisions will be the subject of a future publication.

JHEP10(2016)139
Cross talk between sub-collisions is, however, included in one respect by accounting for energy-momentum conservation in all approaches. As before, we will concentrate on pPb collisions at √ s NN = 5 TeV. The methodology is, however, not limited to this, and generalisation to AA collisions will be the subject of a future publication.

Wounded nucleons and multi-parton interactions
In ref. [25] Bia las et al. noticed that the central particle density in pA collisions scales approximately with the number of "wounded" or "participating" nucleons, dN (pA) /dη ≈ (N wounded /2) dN (pp) /dη. The projectile proton was here included as one of the wounded nucleons, and the distribution in rapidity could be described if each wounded target nucleon gives a contribution proportional to (η − η t )/(η p − η t ) where η t,p are the rapidities of the target and projectile respectively. The wounded projectile proton gives a similar contribution with p exchanged for t. The wounded nucleon model worked well for minimum bias events and low p ⊥ particles, while high p ⊥ particles scale better with the number of N N collisions, which can be understood if the high-p ⊥ particles originate from independent partonic subcollisions. (See e.g. ref. [60].) A model with this feature, called G-Pythia, has been used in analyses by ALICE [17].
These results can be given a heuristic interpretation in terms of the Landau-Pomeranchuk formation time. The formation time for a hadron is, in a frame where p L = 0: This implies that a produced pion will resolve the nucleus at a length scale given roughly by 1/p ⊥ . For p ⊥ < 1 GeV the resolution scale is larger than that of the individual nuclei, while for p ⊥ larger than ∼ 1 GeV, constituents of individual nucleons can be resolved. Below we will compare two models, generated with the help of Pythia 8. The model denoted "Absorptive" is similar to G-Pythia. Here each N N subcollision is treated as a pp collision. 9 The second model, explained in the next subsection, is called FritiofP8 and is more similar to the wounded nucleon model. We note that the models should not be compared on equal footing. From the above arguments, the "Absorptive" model is expected to describe the high-p ⊥ part of the spectrum better, while FritiofP8, being similar to the wounded nucleon model, is expected to describe the low-p ⊥ part, and thus also the total multiplicity, best.
Technically, the subcollisions are in both models generated with Pythia8. This means that for each sub-collision, multiple partonic interactions are created in decreasing order 9 Note that the G-Pythia approach really uses a black disk Glauber calculation with σNN = σin = σ abs + σDp + σDt + σDD, and lets the collisions be a mixture of the four corresponding processes, while we we consider only absorptive collisions, as we believe this is more in line with the original model by Glauber. There is, however, no notable difference for observables in the near central rapidity range, taken with a minimum bias trigger.

JHEP10(2016)139
of p ⊥ with the probability: starting from a maximum scale related to the impact parameter of the sub-collision. The cross section is obtained by treating everything as perturbative QCD 2 → 2 scatterings, but since the cross section diverges at low p ⊥ , it is regulated at low p ⊥ using: Here p ⊥0 is a tunable parameter. Aside from momentum conservation, Pythia8 also rescales the PDF every time a quark has been used in an MPI. When using this MPI model for generating pA collisions we maintain momentum conservation, but do not maintain the rescaling of the PDF between separate NN collisions.

The revived Fritiof model
A very different approach was used in the Fritiof model [23]. Where the Pythia8 MPI model assumes everything can be described by perturbative scatterings, Fritiof imposed a soft model for everything, specifically limiting it's range of validity to low-p ⊥ processes. 10 In the Fritiof model it is assumed that a soft min-bias interaction causes a momentum exchange, which in light-cone variables has the form This produces two excited states assumed to decay like strings stretching the rapidity range between the initial beam rapidities and a point distributed evenly within the kinematically allowed region The result is approximately reproducing the original wounded nucleon model [25], but it is in the Fritiof model also assumed that a secondary encounter with another nucleon will increase the excitation, thus leading to a logarithmic scale breaking. In ref. [61] it was suggested to extend the Fritiof model to include the possibility for a hard scattering and associated bremsstrahlung when the energy is high enough. At LHC collision energies, the necessity for including the possibility for at least one such interaction is apparent.
In figure 7 we show a schematic picture of a projectile proton wounding a number of nucleons. The picture is strongly oversimplified, showing only the main gluon propagators, i.e. no initial/final-state radiation (ISR/FSR) or multi-parton interactions (MPI) in the individual sub-collision. Nucleons ν 1 , . . . ν 4 are wounded absorptively with ν 4 being the hardest or "primary" wounded nucleon, which contributes to the hadronic multiplicity the full rapidity region. The other absorptively wounded nucleons, ν 1 . . . ν 3 , contributes only to the parts of the rapidity range in the nucleus direction. As indicated by the exchanged Pomeron, IP , ν 5 is only diffractively excited, and will also only contribute in the nucleus direction. Thinking in terms of cut Pomeron diagramsà la AGK [62] we show in figure 8 the similarity between the diagram describing diffractive excitation in proton-proton scattering and a fully absorptive proton-deuteron scattering. It is not far fetched to assume that the triple-Pomeron vertex in both cases are distributed in approximately the same way in rapidity, i.e., that the gap size in the single diffractive excitation in pp would be distributed in the same way as the the size of the region of rapidity populated by hadrons from both wounded nucleons in a pd collision.

JHEP10(2016)139
As discussed in section 3.3 the distribution in diffractive masses indicates a fairly flat distribution in rapidity of the triple-Pomeron vertex as ǫ is close to zero. 11 We will therefore assume as a first approximation, that the secondary absorptive collisions in a pA collision can be approximated by single diffractive collisions.
Treating secondary absorptive collisions as single diffractive excitation has an additional added benefit. In the Pythia8 implementation, one can model high mass soft excitation using a perturbative approach where the exited proton can undergo multiple partonic interactions, as in eq. (5.2), and ISR is included. It is thus possible to treat absorptively wounded nucleons differently, depending on whether the mass of the excited system is larger or smaller than a pre-set threshold mass scale.
Finally, in figure 7, we also have ν 5 , which is an standard diffractively excited nucleon, and will be modelled as such.

Comparison to data
We now compare the two methods for particle production, which were introduced above. Stacking absorptive events on top of each other is labelled "Absorptive", we use a black disk Glauber model with σ abs = 67.9 mb to calculate the number of absorptive sub-collisions event by event. The model including both diffractive excitation and the Fritiof-inspired absorptive sub-collisions is labelled "FritiofP8". To calculate the amount of wounded nucleons we use the modified GG model with cross section fluctuations described by the log-normal distribution in eq. (4.21) and including the modifications introduced in section 4.2.4, as well as distinguishing between absorptive and diffractive events using the 2×2-disk modification, introduced in section 4.4.2. All parameters are fitted to pp data.

Centrality estimation and multiplicity
The primary observable we wish to discuss, is the charged particle pseudo-rapidity distribution at different centralities, as measured by ATLAS [16]. When comparing Monte Carlo predictions to data in pp, the work flow has matured greatly over the past years, with the advent of automated frameworks for performing such tasks, such as Rivet [63]. In this framework, equal treatment of theory and unfolded data is ensured by publishing measurements along with an implementation of the analysis. This is not yet tradition in the heavy ion community, and the data comparisons shown here, is the result of our own Rivet implementation on the analysis, based on the paper, with data obtained from HepData [64].
In the experimental analysis by ATLAS, event centrality is calculated by taking fractiles of the distribution in E ⊥ of charged particles in the interval 12 3.1 < η < 4.9. For this particular observable, unfolded data has not been published, but we will still compare theoretical curves for the two previously outlined particle production models. In figure 9 we show the model stacking N t w abs absorptive events spanning the whole rapidity region (denoted "Absorptive") reaches a much higher E ⊥ than the Fritiof inspired model (denoted "FritiofP8") With one absorptive event spanning the whole rapidity region, N t w abs −1 absorptive events modelled as diffractive excitation, and N t w inc − N t w abs events from diffractive excitation. We note that the "FritiofP8" results agree almost perfectly with the data from Atlas, while the Absorptive model reaches significantly higher E ⊥ values. In figure 10 we show pseudorapidity distributions for different centralities, where we have used the same cuts as ATLAS, but reconstructed fractiles from our own distribution. 13 We see that while both models describe the peripheral events reasonably well (which is expected), the new FritiofP8 model based on diffractive excitation does a much better job describing both average multiplicity and the forward-backward asymmetry, as expected.

Inclusive transverse momentum
We now compare to centrality-inclusive charged particle p ⊥ spectra in different ranges of η as measured by CMS [65]. In figure 11(a) we show the transverse momentum distribution for −1.0 < η < 1.0. We see that the FritiofP8 model performs well at low p ⊥ , while Absorptive performs well at high p ⊥ , as expected.
The same picture is seen when going to large negative η ( figure 11(c)), but performance of the FritiofP8 model improves slightly when going to large positive η shown in figure 11(b) for 1.3 < η < 1.8 (the proton side). The "absorptive" model performs as before, but it is rather surprising that the FritiofP8 model performs poorly at high p ⊥ here. One explanation could be that the parton distribution function used for for secondary absorptive events is a Pomeron PDF and not a proton PDF, due to the fact that secondary absorptive events are modelled as single diffractive events. This is also a possible explanation for the poor performance at high p ⊥ in figure 11(a) and figure 11(c).

Uncertainties
The method presented here for generating final states in pA is interesting, as it gives qualitatively correct description of the multiplicity in both the proton and the nucleus direction. It is still mostly a proof-of-principle since, as we will demonstrate here, using the Pythia8 default settings introduces several hidden assumptions. We will discuss these assumptions by giving a rough estimate of the uncertainty associated with each of them. That uncertainty will decrease, or vanish entirely, when the assumptions are dealt with more carefully, one by one, which will be done in one or more future publications.

PDFs and MPI activity
In the previous section we described how the secondary absorptive sub-collisions are approximated as single diffractive excitation events. The perturbative handling of single diffractive events at high masses in Pythia8 relies on a factorised Pomeron approach, where the diffractive state is modelled by a Pomeron-proton collision, including MPI, and we will study two important uncertainties here.
• The Pomeron PDFs used in the MPI machinery are not really appropriate in our model of the secondary absorptive sub-collisions, since it is still really the parton density in the proton which should drive the MPI. To see possible effects, we have tried to make the Pomeron more proton-like, by modifying the PDF used in Pythia8 to have much more small-x gluons. 14 This will increase MPI activity. • Another way of modifying the MPI activity is to change the Pomeron-proton cross section used in Pythia8. This is not a physical cross section, but rather a free parameter in the program which only affects MPI activity, and is adjusted to fit data. The default value of this parameter is 10 mb. We here increase it to it's maximal allowed value, 40 mb, to better reflect a pp absorptive cross section.
In figure 12 we show the pseudo-rapidity distribution of charged particles for the highest centrality bin. The two variations above are labelled "PDF" and "MPI" respectively while FritiofP8 is the same as before.
In total, the envelope of the three lines gives what we believe to be a reasonable, albeit conservative, estimate of the uncertainty so far associated with the approximations regarding the parton densities and the amount of MPIs.

GG uncertainty
In section 4.3.1 we described how the Glauber-Gribov cross section fluctuations could be parameterised with a log-normal distribution. This new parameterisation was used in all the previous data comparison plots, here we show how the traditional parameterisation (also fitted to pp data) compared to the new one. Since distributions of wounded nucleons for the two models differs most in the tail, we are most sensitive in central collisions, In figure 13 we show the uncertainty in central particle production arising from changing the parametrisation of cross section fluctuations. We see that the uncertainty covers data well, but is smaller than the uncertainties in the handling of secondary absorptive events above.  Figure 13. Pseudo-rapidity distribution of charged particle multiplicity for centrality 0-1% compared to two parameterisations for the GG model.

Conclusions and outlook
Collisions of heavy nuclei with protons or each other, are often understood in terms of their deviation from minimum bias pp collisions, which after many years of work on models for multiple partonic interactions, are fairly well understood. The extrapolation from pp to pA or AA, which is usually based on Glauber models, is therefore a crucial step toward understanding the deviations, and thus crucial for understanding how and if a deconfined plasma of quarks and gluons are created in heavy ion collisions, and what kind of impact it has on final state observables. In this article we have discussed Glauber models, their extensions including fluctuations and the impact these extensions on the number of wounded nucleons, followed by a possible way of using that knowledge to generate full, exclusive final states as a sum of pp collisions. We will conclude on these two parts separately, followed by an outlook primarily focused on extending the model for full final states, to also include microscopic models for collective effects.

Fluctuations in the Glauber formalism
Following Good and Walker, we have discussed diffractive excitation as a manifestation of fluctuations in the substructure of the nucleons. As a result of this discussion we identified the inclusively wounded cross section to be the relevant pp cross section for calculations of the number of wounded nucleons. The wounded cross section has contributions from absorptive processes plus diffractively excited target nucleons, in brief form we have for the semi-inclusive pp cross sections: σ w inc = σ abs + σ DD + σ Dt = σ tot − σ el − σ Dp . (6.1)

JHEP10(2016)139
We have discussed the developments in calculating the number distribution of inclusively wounded nucleons in the so-called Glauber-Gribov approach by Strikman and coworkers. By comparing distributions of the cross section to ones calculated with the dipsy model and measurements, we find that the parametrisation of the cross section suggested in the Glauber-Gribov approach does not fully include all fluctuations necessary to describe the ones in dipsy.
The simplest Glauber calculations using a fixed black disk, can clearly not describe any fluctuations, and we want to emphasize that just setting σ NN = σ w inc in such a calculation, is not enough if one wants to calculate the contribution from the wounded nucleons to a centrality defining observable.
We have shown that a Glauber-Gribov calculation with a black disk fluctuating in size, can be used to predict the distribution of inclusively wounded nucleons, if the fluctuating black disk is attributed to fluctuations in the projectile, while averaging over fluctuations in the target nucleon. However, if one wants to separate absorptively wounded nucleons from diffractively wounded ones, it is necessary to also consider fluctuations in the individual target nucleons.
We have suggested a new functional form for the fluctuations in the pp interaction, with a higher tail to a larger cross section. Instead of the parametrization introduced by Strikman et al., and used in several experimental analyses, we suggest a log-normal distribution, which is believed to give a somewhat better description of the inclusive distribution, though not necessarily a more realistic picture of the cross section fluctuations. We have also included fluctuations in the target nucleons by introducing a crude Glauber-like model, where the radii of the projectile and target are allowed to fluctuate independently between two values. The model includes four parameters, and can be fitted to describe four independent semi-inclusive pp cross sections, including the inclusive wounded one. By using this model to include projectile fluctuations in the Glauber-Gribov approach, we separate the inclusively wounded nucleon into absorptively and diffractively wounded ones.
The parametrisations and toy-models studied are built on the assumption that distributions of wounded nucleons can be described solely by fluctuations in the proton size; projectile size for the inclusive distributions, and target nucleon size for distinguishing between absorptively and diffractively wounded nucleons. Thinking about the dynamics of more involved calculations, like the dipsy model used in this paper, it is clear that size fluctuations cannot account for all the relevant physics. In the dipsy model, fluctuations in cross section will also arise when e.g. a small projectile is very dense, and therefore gives rise to a larger cross section. Similarly, a large projectile can be dilute, giving rise to a smaller cross section. In the language of this article, such effects would need to be accounted for by a fluctuating, b-dependent opacity, resulting in a profile function which would go beyond the simple Gaussian, grey -or black disks, but still not be as calculationally involved as the full, dynamical models. This will be the subject of a future publication.

Full final states
We have given a proof-of-principle for an approach to model exclusive final states in pA collisions as a sum of several pp collisions. The approach uses Pythia8 to calculate JHEP10(2016)139 the hardest absorptive sub-collision as a normal non-diffractive pp collision, while the subsequent absorptive collisions are modelled as single diffractive events. This was inspired by the old Fritiof model, which is valid at lower energies, but adds another dimension, as high-mass single diffractive exchanges can now be treated perturbatively and allows for multiple parton-parton scatterings.
We have shown that this approach, in a quite crude implementation, is able to give a reasonable description of some recently published final-state measurements of pPb collisions at the LHC, but the uncertainties in our approach are quite large. In future studies we will try to eliminate these uncertainties, in the hope to get a theoretically well motivated and more accurate description of data.
To do this it is helpful to have data published in a usable form for comparison with event generators. The LHC pp community has come a long way in this respect by publishing many of their results in the form of Rivet routines. In this way the measurements are presented in a clearly reproducible form, including all relevant kinematical cuts and unfolding of detector effects but free from dependence of theoretical models. In order to allow for the development of better event generators for heavy ion collisions, it is imperative that the experimental heavy ion community adapts a similar way of presenting their results.

Outlook: modelling collective effects
In both the theoretical and experimental communities around heavy ion collisions, much attention is given to observables thought to convey information about a possible plasma state created in the collision. A direct application of the work presented in this article, is to use the final state extrapolation as a baseline for calculations of collective effects using microscopic, QCD based models. In the pp community, much attention have been given to models of final state interactions including colour reconnections, rope hadronisation and junction formation [67][68][69]. These models have, in pp been shown to reproduce: the enhancement of strange hadrons to non-strange hadrons in dense environments [70]; flowlike effects in hadron ratios as function of p ⊥ ; and preliminary studies have shown that interactions between strings can produce a ridge [71].
To implement such models in a framework like the one presented here, a necessary component is a good understanding of the sub-collisions in transverse space. Such a picture is not included in e.g. the Pythia8 MPI model, and must therefore be added a posteriori. Guidance can then be had from the dipsy model, which includes detailed information about the transverse space structure, but does not produce final states describing data as well as Pythia8.
Finally, the final state model introduced here should be developed further to model fully exclusive hadronic final states also in AA collisions. Since every projectile here also becomes a target, we suspect that the model cannot be transferred one-to-one, but that some modifications may be needed. Here modelling collective effects using microscopic models is also highly desirable.