Azimuthal angle correlations at large rapidities: revisiting density variation mechanism

In the paper we discuss the angular correlation present in hadron-hadron collisions at large rapidity difference ($\bas\,y_{12}\,\gg\,1$). We find that in the CGC/saturation approach the largest contribution stems from the density variation mechanism. Our principal results are that the odd Fourier harmonics($v_{2n+1}$), decrease substantially as function of $y_{12}$, while the even harmonics ($v_{2n}$ ), increase considerably with a growth of $y_{12}$.


Introduction
In this paper we address the problem of the azimuthal angle correlations of two hadrons with transverse momenta p T 1 and p T 2 and rapidities y 1 and y 2 , at large values of y 12 ≡ |y 1 − y 1 | 1/ᾱ S . Our main theoretical assumption is that these correlations stem from interactions in the initial state. We are aware that, unlike rapidity correlations which at large rapidities are originated from the initial state interactions due to causality reasons [1], a substantial part of these correlations could be due to the interactions in the final states [2]. On the other hand, it has been demonstrated that at small rapidity differenceᾱ S y 12 < 1 the interactions in the initial state [3][4][5][6][7][8][9][10][11] gives the value of the correlations, which describe the major part of the experimentally observed correlations [12][13][14][15][16][17][18][19][20][21][22].
In this paper we concentrate our efforts, on calculating the long range rapidity part of angular correlations with large value of the rapidity difference y 12 . All previous calculations, assumed thatᾱ S y 12 < 1 [3][4][5][6][7][8][9][10][11]. It turns out, that in this kinematic region, the main source of the azimuthal angle correlations, is the Bose-Einstein correlations of identical gluons, which corresponds to the interference diagram in the production of two partonic showers. Intuitively, we expect that the correlations in the process, where two different gluons are produced from two different partonic showers, should not depend on the difference of rapidities (y 12 ), as well as on the values of y 1 and y 2 . Using the AGK cutting rules [23] * one can prove that the two gluon correlations can be calculated using the Mueller diagrams [24] of Fig. 1 [24] for two partonic showers production. Fig. 1-a describes the square of the production amplitudes, while Fig. 1-b corresponds to the interference diagram which leads to the Bose-Einstein correlations. The wavy lines show the BFKL Pomeron [25,26], while the helical lines denote gluons. Fig. 1-c shows the example of a more complicated structure of the partonic cascades, than the exchange of the BFKL Pomeron. The colour of the lines indicates the parton shower.
The diagrams of Fig. 1 lead to correlations which do not depend on y 1 and y 2 , but only forᾱ S y 12 1. For large y 12 the contributions of Fig. 1 decrease. The main goal of this paper to find the contributions which survive at large y 12 (ᾱ S y 12 1).
At large y 12 , we have to take into account the emission gluons, with rapidities y 2 < y i < y 1 which transform the Mueller diagram of Fig. 1-b to more general diagrams of Fig. 2. The general features of Fig. 1b is that the lower Pomerons carry momenta Q T + p 12 and − Q T − p 12 with p 12 = p T 1 − p T 2 . Q T denotes the momentum along the BFKL Pomeron. After integration over Q T , we obtain p 12 ∼ 1/R h , where R h is the size of the target(projectile), which has a non-perturbative origin. Roughly speaking, the correlation function turns out to be proportional to G (p 12 ), where G denotes the non-perturbative form factor of the target or projectile [9]. This conclusion stems from the value of the typical Q T for the BFKL Pomeron, which is determined by the size of the largest dipoles in the Pomeron. Fig. 2 does not have these features. We will show that the azimuthal angle correlations originate from the integration over Q T (see Fig. 2), due to the structure of the vertices of emission of the gluons with p T 1 and p T 2 , which have contributions * In the framework of perturbative QCD for the inclusive cross sections, the AGK cutting rules were discussed and proven in Refs. [27][28][29][30][31][32][33][34]. However, in Ref. [35] it was shown that the AGK cutting rules are violated for double inclusive production. This violation is intimately related to the enhanced diagrams [31,34,35] and to the production of gluon from triple Pomeron vertex. It reflects the fact that different cuts of the triple BFKL Pomeron vertex with produced gluon, lead to different contributions. Recall, that we do not consider such diagrams.
proportional to ( p T 1 · Q T ) n ( p T 2 · Q T ) n . Recall, that these kind of vertices, are the only possibilities to obtain angular correlations in the classical Regge analysis [36]. This mechanism for azimuthal angular correlations was suggested in Ref. [37] (see also Refs. [7,[38][39][40]), and in the review of Ref. [40], was called, the density variation mechanism.
The paper is organized as follows. In the next section we discuss the contribution of the diagram of Fig. 2 in the momentum representation. In the remainder of the paper, we will use the mixed representation: the dipole sizes and momentum transferred (Q T ), which will be introduced in section 3. Section 4 is devoted to the discussion of the single inclusive production in the Colour Glass Condensate (CGC)/saturation approach. The double inclusive production is considered in section 5, in which the rapidity dependence of the master diagram of Fig. 2 will be calculated. In section 6, we estimate the angular correlation function and Fourier harmonics v n , and we present our prediction for dependence of v n on the difference of rapidities (y 12 ). In section 7 we draw our conclusions and outline problems for future investigation.   Figure 2. The generalization of Fig. 1-b forᾱ S y 12 1. The wavy lines show the BFKL Pomeron [25,26], while the helical lines denote gluons. The colour (blue and black) of the lines indicates the parton shower.

Correlations in the momentum representation
The double inclusive cross section of Fig. 2 takes the following form as well as all other functions φ of this type, are the correlation functions which at Q T = 0, give the probability to find a gluon with transverse momentum k T in the hadron(nucleus) of the projectile (target). φ k T , − k + Q T ; k T , − k T + Q T describes the interaction of two gluons with momenta k T and k T , which scatter at momentum transferred Q" T . N (Q T ) is a pure phenomenological form factor that describes the probability to find two Pomerons in the projectile or target, with transferred moment Q T and − Q T . C F = N 2 c − 1 /2N c where N c is the number of colours. The Lipatov vertex Γ µ (k T , p T 1 ) has the following form At high energies the parton densities φ(. . . ; Y ) in Eq. (2.1) and in Eq. (2.5), are proportional to exp (∆ BFKL Y ) for the BFKL Pomeron, where ∆ BF KL = 2.8ᾱ S is the intercept of the BFKL Pomeron. Bearing this in mind, one can see, that the interference diagram for the double inclusive cross section does not depend on y 1 , y 2 or on y 12 .
The main diagram of Fig. 1-a, also does not depend on rapidities, and its expression has the following form:

BFKL Pomeron in the mixed representation
For a more convenient presentation, it turns out that the most economical way of calculating the diagram of Fig. 2, is to use the mixed representation of the BFKL Pomeron Green's function, G r, R, Q T , Y , where r and R are the sizes of two interacting dipoles , Q T denotes the momentum transferred by the Pomeron, and Y the rapidity between the two dipoles. This Green's function is well known [26] and we discuss it here for the completeness of presentation, referring to Refs. [26,41] for all details. It has the following form: where ψ(z) is the Euler ψ-function (see Ref. [42] formulae 8.36) and ∆ BFKL =ᾱ S 4 ln 2, D =ᾱ S 14ζ(3), ξ = ln r 2 1 /r 2 2 . Each term in Eq. (3.1) has a very simple structure, being the typical contribution of the Regge pole exchange: the product of two vertices and Regge-pole propagator. From Eq. (3.2) one can see that at large Y the main contribution comes from the term with n = 0, and in what follows we will concentrate on this particular term. The vertices with n = 0 have been determined in Refs. [26,41], and they have an elegant form in the complex number representation for the point on the two dimensional plane: viz.
For r(x, y) : ρ = x + iy; ρ * = x − iy; For Q T (Q x , Q y ) : q = Q x + iQ y ; q * = Q x − iQ y ; (3.3) Using this notation the vertices have the following structure: At Q T → 0 this vertex takes the form: For small values of ν (which are related to the region of largeᾱ S Y 1), Eq. (3.5) can be simplified and reduced to the form: Using that The contribution of the first term in Eq. (3.1) can be reduced to the following form for the scattering amplitude of two dipoles with sizes r 1 and r 2 : where ∆ BFKL and D are defined in Eq. (3.2).
In the derivation of Eq. (3.9) we neglected the contributions that are proportional to since this contribution will be the same as in Eq. (3.9), but with, ξ = ln To integrate over ν, we use the method of steepest descent, and the expansion of ω (ν, 0) at small ν (diffusion approximation, see the second equation in Eq. (3.2)).
N (r 1 , r 2 ; Y ) denotes the imaginary part of the dipole-dipole sacttering amplitude at Q T = 0, which is related to the cross section. One can check that Eq. (3.9) has the correct dimension.
4 Single inclusive production in a one parton shower 4.1 BFKL Pomeron: the simplest approach for a one parton shower The single inclusive cross section resulting from the one BFKL Pomeron is known, and it is equal to The relation between the parton densities φ and the Green's function of the BFKL Pomeron has been given in Ref. [29]: where N BFKL (r, r 1 ; Y ) is given by Eq. (3.1) or by Eq. (3.9), in the high energy limit. Eq. (4.2) can be re-written as follow We have Plugging in Eq. (4.2) and Eq. (4.4) into Eq. (4.1) we obtain [29] where N pr and N tr denote the probability to find a dipole in the projectile and target, respectively. r 1 and r 2 are the typical dipoles sizes in the projectile and target.
As can be seen from Eq. (2.1) we need to generalize Eq. (4.5) for the case Q T = 0. Eq. (4.1) has to be replaced by we re-write Eq. (4.5) in the form

General estimates
It should be stressed that the single inclusive production has the form of Eq. (4.5) and Eq. (4.8) as was shown in Ref. [29] for the general structure of the single parton shower. For example, for the process shown in Fig. 1-c. We need only to substitute N G tr (r, N tr (r, r 2 ; y, Q T ) is a solution to the non-linear evolution equation. For the case of inclusive production, we can considerably simplify estimates noting that where Q s (y) denotes the saturation momentum.
In other words, the main contribution to inclusive production comes from the vicinity of the saturation scale, where r 2 Q 2 s ≈ 1. Fortunately, the behavior of N in this kinematic region is determined by the linear BFKL evolution equation [43][44][45] and has the following form [46] N tr (r, We have seen in Eq. (3.8) that for Q T = 0 , the scattering amplitude decreases at Q 4 T r 2 r 2 2 1. Therefore, we need to consider rather small values of Q T : Q 4 T r 2 r 2 2 ≤ 1. The product of vertices that determines the amplitude has two terms ( see Eq. (3.5) ) which are proportional to r 2 /r 2 2 iν and to Iν . Therefore, the maximum of ∇ 2 r N can be reached if r 2 /r 2 2 e κ y ∼ 1 and Q 4 T r 2 r 2 2 e κ y ∼ 1 and the amplitude then has the following form The first term does not depend on Q T and, therefore, the upper limit of the integral over Q T , goes up to (Q max T ) 2 ≈ 1/(r r 2 ). The second term both for Q 2 T r r 2 < e − 1 2 κ y and for Q 2 T r r 2 > e − 1 2 κ y turns out to be small. Indeed, in the first region the amplitude is small, while in the second region we are deep in the saturation domain where ∇ 2 r N → 0. Hence, we expect that in the integral over Q T , the first term gives a larger contribution than the second term and we will keep only this contribution in our estimates. 5 Double inclusive cross section for two parton shower production

The simplest diagram
In this section we calculate the simplest diagram of Fig. 2. We need to integrate the product of two BFKL Pomerons over Q T (see Eq. (2.5)): From Eq. (2.5) in the momentum representation, we see that r 1 = r 1 (r 2 = r 2 ) but they are close to each other, being determined by the same momentum k T . We assume that p T 1 < k T , since k T ∼ Q s (Y − y 1 ) µ soft . Considering r 1 ≈ r 1 r 2 ≈ r 2 we will show that in the integral over Q T , the typical Q T ∼ 1/r 2 . In other words, the dependence of Q T is determined by the largest of interacting dipoles.
Assuming that both ν 1 and ν 2 are small, we see that all four terms are equal to each other, and the integral can be written as follows.
The appearance of the pole ν 1 = −ν 2 indicates that the contribution from this kinematic region is large.
Closing the contour of integration on ν 2 over the pole, we obtain I = 2 6 π ν 2 1 1 r 2 2 (5.4) Actually, the double inclusive cross section depends on ∇ 2 N as we argued in the previous section. Repeating the procedure for we obtain for small ν 1 and ν 2 : Taking the integral over ν 1 , using the method of steepest descent, we obtain the following contribution: where ∆ BFKL and D are defined in Eq. (3.2).
Rewriting Eq. (2.5) in the coordinate representation we obtain: In Eq. (5.8) we neglected the terms which are proportional to Q 2 T (see Eq. (2.5)) as the typical Q T are small as we have argued, and because these terms do not lead to additional correlations in the azimuthal angles.
We have discussed the integral over Q T , and it has the form of Eq. (5.7). The extra e i r 1 · Q T give an additional numerical factor, replacing 2 5 by 2 7 in Eq. (5.7). To integrate over k T and k T we replace Now we can take the integrals over r i bearing in mind Eq. (5.7) and The integrals overr 1 andr 1 have the following form (see Ref. [42] formulae 6.511 (6)) Using Eq. (5.10) we obtain Collecting Eq. (5.10),Eq. (5.1) and Eq. (5.11) we see that the main contribution stems from the region kr 2 1 and the integral over k T has the form The integral over k T has the same structure while the integration in Eq. (5.1) goes to infinity. As the result we can reduce the integral to the form; Finally, (5.14)

The CGC/saturation approach
The integral over k T in Eq. (5.13) has an infra-red singularity with a cutoff at p T 2 since we assume that p T 2 is the smallest momentum. This reflects the principle feature of the BFKL Pomeron parton cascade, which has diffusion, both in the region of small and large transverse momenta. On the other hand we know that the CGC/saturation approach suppressed the diffusion in the small momenta [27], providing the natural cutoff for the infrared divergency. We expect that such a cutoff will be the value of the smallest saturation momenta: Q s (Y − y 1 ) or Q s (y 2 ), which will replace one of p 2 T 2 in the dominator of Eq. (5.14). Therefore, we anticipate that for the realistic structure of the one parton shower cascade, (see Fig. 1-c for example), the contribution for the double inclusive cross section will be different.
We need to specify the behavior of the scattering amplitude in the vicinity of the saturation scale. We have discussed the basic formulae [46] of Eq. (4.11), but for integration over the dipole sizes we need to know the size of this region. The scattering amplitude can be written in the form: where ω(γ, 0) is given by Eq. (3.2), replacing 1 2 +iν ≡ γ and ξ = ln r 2 1 /r 2 2 . The saturation scale is deterring by the line on which the amplitude is a constant (C) of the order one. This leads to the following equation for the saturation scale [43,46]: which results in the value of γ cr given by the equation: which gives γ cr = 0.37 and the equation for the saturation momentum: Expanding the phase ω (γ, 0) Y − (1 − γ)ξ in the vicinity ∆ξ = ξ − ξ s and ∆γ = γ − γ cr we obtain At first sight Eq. (5.19) shows that the amplitude has a maximum at τ = r 2 1 Q 2 s = 1. However, this is not correct. Eq. (5.19) gives the correct behavior for τ < 1 while for τ > 1 we need to take into account the interaction of the BFKL Pomerons and the non-linear evolution, generated by these interactions. The general result of this evolution is the fact that the amplitude depends on one variable [47] τ , i.e. N (τ ) (as it has geometric scaling behavior). The peak at τ = 1 appears in From Eq. (5.20) we can conclude that the width of the distribution in r 2 1 is of the order of Q 2 s , but depends crucially on the model for the Pomeron interaction. In Fig. 3-a we plot this value for the behavior of the scattering amplitude deep in the saturation domain (see Ref. [48]).  This approach is not correct for τ → 1 and −∇ 2 N = 1.58 at τ = 1, but it starts to be small at τ > 2, which could be large enough to trust the formulae of Ref. [48]. At least such a conclusion can be justified considering the fit of the DIS data in the saturation model of Refs. [49,50], which based on the idea of Ref. [51], and has the correct behavior of the scattering amplitude, both deep in the saturation domain, and near τ = 1. Hence, we expect that ∇ 2 N decreases faster than we can see from Eq. (5.19). Bearing these conclusions in mind, we will calculated the contribution of Fig. 2, keeping all N in Eq. (5.8) in the vicinities of the saturation scales, by replacing We will show in the folowing, that we cannot take all six Pomerons in the vicinities of the saturation scale. We have to take two of them, either deep in the saturation domain, or in the perturbative QCD region. We choose to take two Pomerons between rapidities y 1 and y 2 , (see Fig. 3-b) i.e. in the perturbative QCD region. Unfortunately, we cannot use the AGK cutting rules [23], which state that these Pomerons will be not affected by the Pomeron interaction, and the contributions of these interactions (see red Pomeron in Fig. 3-b) are canceled. Indeed, it has been proven that for the double inclusive production [52], that they do not work in perturbative QCD. On the other hand, these Pomerons carry transverse momentum Q T , unlike the others one in the diagram, which is larger than the saturation scale Q s (y 2 ) and therefore, their contributions are suppressed in comparison with the other Pomerons in Fig. 2. In addition our choice provides the natural matching with the regionᾱ S y 12 < 1.
The integration over Q T will produce the same result as Eq. (5.7), as in the previous section. We re-write the integration over φ i (see Eq. (5.9)) in the following way: We see that the integrals over r 1 and r 1 leads to r 1 ∼ 1/Q s (Y − y 1 ) and r 1 ∼ 1/Q s (Y − y 1 ). The same holds for the integrals over r 2 and r 2 , leading to r 2 ∼ 1/Q s (y 2 ) and r 2 ∼ 1/Q s (y 2 ) . Assuming that Q s (Y − y 1 ) > Q s (y 2 ) we conclude that r i and r i are much smaller than r 2 and r 2 . Replacing whereγ = 1 − γ cr , we obtain from Eq. (5.21) that integration over r takes the form Recall that we consider Q s = Q s (Y − Y 1 ) in Eq. (5.23). For p T 1 Q s (Y − y 1 ) we can replace e −i p T 1 · r 1 = 1 in Eq. (5.21). In this case the integral has the form where τ = k 2 /Q 2 s . The integral over r in the lower part of the diagram takes the form: Using Eq. (5.25) for p T 2 Q s (y 2 ) the integral over k T can be reduced to Finally, collecting all numerical coefficients, we obtain where constant C is the value of the amplitude at τ = 1.

This contribution is proportional to
. We need to estimate the diagram of Fig. 1-a (see Eq. (2.7)). This diagram can be re-written as Examining Eq. (4.5), one can see that in general case when Y − y 1 = y 1 and Y − y 2 = y 2 all four Pomerons cannot be in the vicinity of the saturation scale. Actually we have two kinematic regions which give the maximal contributions (assuming Q s (Y − y 1 ) > Q 2 s (y 1 )): 1.
In the region 1 the upper Pomeron is in the vicinity of the saturation scale, while the lower Pomeron is in the perturbative QCD region. In region 2 the lower Pomeron is in the vicinity of the saturation scale, and the upper Pomeron is deep inside the saturation domain. As we have discussed (see Fig. 3-a) ∇ 2 N decreases in the saturation region much faster than in the perturbation QCD region and, therefore, we bf assume that the kinematic region 1 gives the largest contribution. Hence, for p T 1 Q 2 s (y 1 ) we obtain In Eq. (5.29) we used backward evolution, from the saturation boundary where N = C.
The ratio of two contributions takes the following form: One can see that Eq. (5.30) demonstrates the additional suppression in comparison with the calculation of the simplest diagram, due to infrared cutoff at Q s (y 2 ) instead of p T 2 . The factor exp (2∆ BFKL y 12 ) reflects the fact that two BFKL Pomerons between rapidities y 1 and y 2 are taken in the perturbative QCD region. It should be stressed that we can trust our estimates only for values of y 12 at which the exchange of the BFKL Pomeron with rapidity y 12 give the contribution smaller than C. This condition means that 1 (2 D y 12 ) e 2∆ BFKL y 12 < C (5.31) Taking ∆ BFKL = 0.25 and Q 2 s (Y ) ∝ exp (λY ) with λ = 0.25 (these values correspond to the BFKL phenomenology) we obtain that the l.h.s. of Eq. (5.31) is smaller than 0.15 for y 12 ≤ 7. Therefore, we can trust our estimates shown in Fig. 4 for C > 0.15. We are used to take C = 0.3 which lead to the contribution of the shadowing corrections of the order of 30%.   Fig. 4-a is the estimates for the y 12 dependence of the Bose-Einstein contribution at small y 12 [8,53]. Fig. 4-a and Fig. 4-b give the estimates in the leading order of perturbative QCD withᾱ S = 0.25. In Fig. 4-c we take ∆ BFKL = 0.25 and Q 2 s (Y ) ∝ exp (λY ) with λ = 0.25. These numbers correspond to the BFKL phenomenology.
In Fig. 4 we plotted ratio R as function of y 12 for y 12 ≤ 7 (see Eq. (5.31). One can see that the ratio increases for large y 12 .

Azimuthal angle correlations
The azimuthal angle correlations stem from terms Q T · r i n in the vertices (see Eq. (3.5) and Eq. (3.6)). Indeed, after integrating over r i these terms transform to expressions of the following type [37]: , which lead to term of ( p T 1 · p T 2 ) m . We have illustrated in Eq. (3.5) and Eq. (3.6) how these originate from the general form of BFKL the Pomeron vertices in the coordinate representation. From Eq. (3.5) and Eq. (3.6) only terms proportional to Q T · r i n with even n appear in the expansion. Therefore, the azimuthal angle (φ) correlation function contains only terms cos 2n (φ), and it is invariant with respect to φ → π − φ. In other words, v n with odd n, turn out to be zero. Hence, we have the first prediction: the value v n with odd n should decreases with y 12 , and their dependence should follow the dotted lines in Fig. 4-a.
We return to Eq. (5.1) and integrate over Q T , collecting terms that depend on the angles between Q T and r i , which we have neglected in the previous section. As we have learned the typical values of Q T ∝ 1/r 2 ∼ 1/r 2 where r 2 and r 2 are larger than r 1 and r 1 . In other words , we showed that the main contributions stem from the kinematic regions: we conclude that r 1 (r 1 ) r 2 (r 2 ). The typical Q T is determined by the largest dipoles and, therefore, we expect Q T ≈ 1/r 2 (1/r 2 ), as has been demonstrated above. Bearing these estimate in mind, we can replace vertices V ν 1 r 1 , Q T and V ν 2 r 1 , Q T in Eq. (5.1) by Eq. (3.6) in which we put Q T = 1/r 2 and Q T = 1/r 2 , respectively. Taking into account that r 1 /r 2 1(r 1 /r 2 1) we obtain At first sight Eq. (6.1) should enter two angles between Q T and r 1 and r 1 , respectively. However, in the integrand for integration over r i (see Eq. (5.9)) depends only on one vector p T 1 . Therefore, after integration over all angles, we find that the angle φ in Eq. (6.1) is the angle between Q T and p T 1 .
For vertices V * ν 1 r 2 , Q T and V * ν 2 r 2 , Q T in Eq. (5.1) we use Eq. (3.8). Finally, we need to evaluate the integral Eq. (6.1) Each term in Eq. (6.6) can be factorized as a product of two functions which depend on r i 1 and on r i 2 . Bearing this feature in mind we calculate each term going to the momentum representation using Eq. (5.21). We obtain a product of functions of k T . Each of these function has the following general form: As we have seen the dependence on r i stem from the integration over Q T or, in other words, from I Q In I Q dependence on r 1 and r 1 can be extracted explicitly, leading to F (r) ∝ 1/r. Hence the momentum image for Eq. (6.7) has a simple form: For j = 2 and j = 4 which we need to calculate Eq. (6.6) we have Note that for integration over r 1 , Eq. (6.8) takes the form The term r 2 1 ( n 1 · n 2 ) 2 + r 2 1 ( n 1 · n 2 ) 2 can be re-written as r 1,µ r 1,ν + r 1,µ r 1,ν r 2,µ r 2,ν and in the momentum representation it looks as The expressions for A and B can be written in a general form. Assuming that both p T 1 and p T 2 are smaller than Q s (y 2 ), we can expand the answer, only taking into account terms that are proportional to v n,n = 2π 0 dφ C p T , φ cos (n φ) In Eq. (5.30) we have calculated the part of C (p T , φ) which does not depend on φ, which coincides with C (p T , φ = 0) = R of Eq. (5.30) for Q s (Y − y 1 ) Q s (y 2 ). To calculate the contribution to C, which depends on φ, we need to take the separate integrals over ν 1 and ν 2 since the terms, which are proportional to cos 2 (φ) and cos 4 (φ), do not have a pole at ν 1 = −ν 2 (see Eq. (6.5)). These integrations lead to the following extra factor in C (p T , φ) − C (p T , φ = 0) . We took factors proportional to p T from the expression for A (k T , p T 1 ) and A (4) (k T , p T 1 ) putting p T 1 = p T 2 = p T . To find the final correlation function and v 2,2 and v 4,4 , we need to collect all numerical factors that come from A (k T , p T 1 ) , A (4) (k T , p T 1 ) and Eq. (6.6), and to integrate over φ, as given in Eq. (6.17).
Note, that in the symmetric kinematics, where Y − y 1 = y 2 = 1 2 (Y − y 12 ), ξ = 0 and Eq. (6.19) vanishes. In this case, we have to use Eq. (3.5) instead of Eq. (3.6), keeping track of the corrections, which are proportional to ν i . As the result, we can consider ξ = 0 in Eq. (6.19), but we need to replace factor ξ 2 by 1.
Eq. (5.30) and Eq. (6.19) suffer la numerical uncertainties, which stem both from the values of soft parametersμ soft and µ soft as well as the values of the saturation scale at low energies, and from the integration in Eq. (5.23) and Eq. (5.25), which were taken neglecting contribution from the region τ < 1. On the other hand, the contribution to the double inclusive cross sections of the diagram of Fig. 2 at α S y 12 1 coincide with the contribution of Fig. 1 Therefore, to obtain the realistic estimate we use the following procedure of matching where v 2 (p T = 0.5 GeV ) | F ig. 1−b and v 4 (p T = 0.5 GeV ) | F ig. 1−b are taken from Ref. [11] where the estimates were performed based on the model for soft interaction which describes all features of soft interaction at high energy and provides the interface with the hard processes.   5 shows the p T and y dependence of the v 2 and v 4 using Eq. (6.21) for normalization. In addition we take ∆ BFKL = 0.25 and Q 2 s (y) ∝ exp (λ y) with λ = 0.25. These values correspond to the BFKL Pomeron phenomenology. We believe that this figure illustrates the scale of rapidity dependence and will be instructive for future experimental observations.

Conclusions
In this paper we generalize the interference diagram, that described the Bose-Einstein correlation for small rapidity differenceᾱ S y 12 1, to include the emission of the gluons with rapidities (y i ) between y 1 and y 2 (y 1 , y i < y 2 ). We calculate the resulting diagram in CGC/saturation approach and make two observations which we consider as the main result of this paper. The first one is a substantial decrease of the odd Fourier harmonics v 2n+1 as a function of the rapidity difference y 12 ( see Fig. 4-c). The second result is, that even Fourier harmonic v 2n has a rather strong dependence on y 12 , showing a considerable increase in the region of large y 12 ( see Fig. 5). We believe that our calculations, that have been performed both for the simplest diagrams and for the CGC/saturation approach, will be instructive for further development of the approach especially in the part that is related to the integration of the momenta transferred by the BFKL Pomerons.
We demonstrated in this paper the general origin of the density variation mechanism, whose nature does not depend on the technique that has been used. This mechanism has to be taken into account, since it leads to the values of Fouriers harmonics that are large and of the order of v n that have been observed experimentally.
We hope that the paper will be useful in the clarification of the origin of the angular correlation, especially for hadron-hadron scattering at high energy. We firmly believe that the experimental observation of both phenomena: the sharp decrease of v n with odd n and the substantial increase of v n with even n as a function of y 12 , will be a strong argument for CGC/saturation nature of the angular correlations.