1.Guizhou Key Laboratory in Physics and Related Areas, Guizhou University of Finance and Economics, Guiyang 550025, China 2.Department of Physics, Guizhou University, Guiyang 550025, China 3.Key Laboratory of Quark and Lepton Physics (MOE) and Institute of Particle Physics, Central China Normal University, Wuhan 430079, China Received Date:2020-03-03 Available Online:2020-07-01 Abstract:The exclusive photoproduction of vector mesons ($J/\psi$ and $\phi$) is investigated by considering the next-to-leading order corrections in the framework of the color glass condensate. We compare the next-to-leading order modified dipole amplitude with the HERA data, finding a good agreement. Our studies show that the $\chi^2/d.o.f$ from the leading order, running coupling, and collinearly improved next-to-leading order dipole amplitudes are 2.159, 1.097, and 0.932 for the elastic cross-section, and 2.056, 1.449, and 1.357 for the differential cross-section, respectively. The results indicate that the higher-order corrections contribute significantly to the vector meson productions, and the description of the experimental data is dramatically improved once the higher order corrections are included. We extend the next-to-leading order exclusive vector meson production model to LHC energies using the same parameters obtained from HERA. We find that our model provides a rather good description of the $J/\psi$ and $\phi$ data in proton-proton collisions at 7 TeV and 13 TeV in LHCb experiments.
HTML
--> --> --> -->
2.1.Exclusive photoproduction of vector mesons in dipole model at non-zero momentum transfer
In terms of the dipole model, the vector meson production in an exclusive diffractive $ \gamma^* p\rightarrow Vp $ scattering can be viewed as three separated subprocesses [48], as shown in Fig. 1. The first subprocess is the formation of a dipole (a quark-antiquark pair) derived from a virtual photon fluctuation. The second subprocess is the interaction between the dipole and the proton via exchanging gluons. The last subprocess is the recombination of the outgoing quark-antiquark pair to generate a final vector meson. Therefore, the scattering amplitude of the diffractive process can be factorized into three ingredients: the photon wave function, the dipole-proton scattering amplitude, and the vector meson wave function. Placing all ingredients together, one can write the imaginary part of the scattering amplitude for a vector meson production as Figure1. Schematic diagram of a vector meson production in $\gamma^* p\rightarrow Vp$ within the dipole model. Three separated subprocesses were denoted by I, II, and III, respectively.
where z is the longitudinal momentum fraction of the incoming photon carried by a quark, x is the Bjorken variable, and $ Q^2 $ is the photon virtuality. The variable q denotes the momentum transfer, whose relationship with the squared momentum transfer is $ t = -{{q}}^2 $. The remaining two dimensional vectors r and b are the transverse size of the quark-antiquark dipole and the impact parameter, respectively. $ \Psi $ is the wave function of the incoming photon, which can be accurately calculated by QED, and it is well known in the literature [49, 50]. $ \Psi_{V}^{*} $ denotes the final vector meson wave function, unlike the photon wave function, it has various prescriptions as we shall discuss at the end of this section. $ (\Psi_{V}^{*}\Psi)_{T,L} $ represent the transverse and longitudinal overlap function between the photon and vector meson, respectively. We note that Eq. (1) is a scattering amplitude containing only the forward component. To obtain the nonforward scattering amplitude, one can multiply the forward wave functions by a phase factor $ \exp[\pm{\rm{i}}(1-z){{r}}\cdot{{q}}/2] $, as was done in Ref. [51]. Using this approach and assuming that the S-matrix is purely real (or the amplitude is purely imaginary), the scattering amplitude can be written as
where $ T(x,{{r}},{{b}}) = 1-S(x,{{r}},{{b}}) $ describes the scattering amplitude between the dipole and proton, which contains all basic information regarding the strong interactions between the dipole and proton. By taking into account the corrections from the real part of the scattering amplitude and the skewness effect, the differential cross-section of an exclusive vector meson photoproduction can be written as [18]:
where $ \beta $ is the ratio of the real to imaginary part of the scattering amplitude, and the factor $ (1+\beta^{2}) $ is to include the correction from the missing real part of the scattering amplitude due to the amplitude, $ {\cal{A}}^{\gamma^* p\rightarrow Vp}_{T,L} $, in Eq. (2), only considering the contribution from the imaginary part. The skewness effect factor $ R_{g} $ is derived from the fact that the momentum fraction of the exchanging gluons between the proton and dipole legs can be different. The parameters associated with these two corrections can be expressed by the imaginary part as follows:
The dipole-proton scattering amplitude originates from the solution to the evolution equations, such as the IIM model [52] inspired by the LO BK equation. In most cases, the impact parameter dependence is disregarded in the BK equation, as the dipole amplitude develops a power-like b behaviour, called Coulomb tails, which yield unphysical results, i.e., the total cross-section violation of the Froissart unitarity bound. To avoid the above-mentioned difficulty, a general strategy is to build an impact parameter independent dipole amplitude inspired by the BK equation, then a model is employed to include the impact parameter dependence, such as two typical models IP-Sat [31] and b-CGC [35]. In this study, we use almost the same scheme as described, albeit with an impact-parameter-independent dipole amplitude resulting from a numerical solution to the LO, rc, ci BK evolution equations. We introduce the impact parameter via multiplying the numerical dipole amplitude with a Gaussian b dependence. In view of the advantage of the method②, which was proposed in Ref. [33] by Marquet, Peschanski, and Soyez (MPS), in study the t-distribution of differential cross sections of photoproduction of vector mesons, we shall follow the MPS strategy in this study. Following Ref. [33], the dipole-proton scattering amplitude can be rewritten in terms of the momentum transfer q instead of the impact parameter b by using the Fourier transform
where the factor $ {\rm{e}}^{-B{{q}}^{2}} $ originates from the nonperturbative effects, R can be interpreted as the radius of proton and $ {\cal{N}}(r,x) $ is an impact parameter independent dipole amplitude. We would like to point out that B and R are free parameters in our fit, which shall be determined by fitting to HERA data. 22.2.Dipole evolution equations -->
2.2.Dipole evolution equations
A key ingredient to calculate the differential cross-section is the dipole-proton scattering amplitude. It is known that almost all past studies on the differential cross-section of vector meson production in the framework of CGC were hovered on the LO level in the literature [53, 54]. Although the LO dipole amplitude can describe the diffractive vector meson production experimental data at HERA at certain uncertainties [33, 55], the precision of the model has to be improved to distinguish the dynamic mechanism of the CGC evolution from the DGLAP evolution, as the DGLAP formulism also provides a good description of the data [34]. Indeed numerous efforts have been made to improve the accuracy of the CGC theory by including other higher order corrections, such as quark loops [10, 11], gluon loops [12], and pomeron loops [56]. The running coupling effects dramatically slow down the evolution of the gluon system, which yield a good description of the latest data from HERA on reduced cross-sections [45]. Similarly, the direct numerical solution of the full NLO BK equation also shows that it slows down the evolution [40]. Based on the significance of the NLO corrections, we extend the LO vector meson production formalism to the NLO of this study. In the next section, we demonstrate that the descriptions of the experimental data are dramatically improved once the NLO corrections are included. The LO BK equation describes the evolution of a quark-antiquark (with a quark at $ x_{\bot} $ and an antiquark at $ y_{\bot} $) dipole with the rapidity Y by the emission of a soft gluon. In the large $ N_c $ limit, this can be written as
where $ \bar{\alpha}_s = \alpha_sN_c/\pi $. Here, $ z_{\bot} $ denotes the transverse coordinate of emitted gluon in the evolution. In Eq. (9), we used the notation $ {{r}} = x_{\bot}-y_{\bot} $, $ {{r}}_1 = x_{\bot}-z_{\bot} $, and $ {{r}}_2 = z_{\bot}-y_{\bot} $ to denote the transverse size of the parent and new daughter dipoles, respectively. The BK equation is obtained at a leading logarithmic approximation, and it has been found that it is insufficient when compared with experimental data [44, 45, 57]. Therefore, significant efforts have been made to improve the understanding of the dipole's evolution at NLO accuracy. The first improvement to the LO BK equation was performed by including quark loops. After resumming $ \alpha_s N_f $ to all orders, one can obtain an evolution equation with running coupling corrections [10, 11], which is called an rcBK equation. The rcBK equation is given by
The numerical solution of the rcBK equation was obtained by Albacete et al. [44, 45]. They found that the proton structure function can be efficiently described under this evolution equation. However, the quark loops corrections are not the only source of the higher order corrections, the complete NLO corrections should also include the contributions from gluon loops and the tree gluon diagrams with quadratic and cubic nonlinearities [12]. Considering all these contributions, we obtain the full NLO BK evolution equation
In Eq. (13), we employed the notation $ {{r}}_1' = x_{\bot}-z_{\bot}' $, $ {{r}}_2' = y_{\bot}-z_{\bot}' $, and $ {{r}}_3 = z_{\bot}-z_{\bot}' $ to denote the transverse size of dipoles. From Eq. (14), we can see there is a double logarithmic term $ \ln\dfrac{r_1^2}{r^2}\ln\dfrac{r_2^2}{r^2} $ in the evolution kernel, which renders the full NLO BK equation unstable [40]. The solution can turn to a negative value for some region due to the double logarithmic term. Thus, one needs to make a resummation of these double logarithms under the double logarithmic approximation (DLA), as it has done by Iancu et al. in Ref. [41]. When this resummation is applied to the full NLO BK equation, the double logarithmic term is removed from kernel $ K_1 $, and the resummation will modify kernel $ K_1 $ by multiplying it with kernel
with $\rho = \sqrt{\ln\dfrac{r_1^2}{r^2}\ln\dfrac{r_2^2}{r^2}}$. In addition to the double logarithmic term, the single transverse logarithms (STL) will also generate large logarithmic corrections to the evolution equation, as shown in Ref. [36]. The effect of the single transverse logarithm resummation will also modify kernel $ K_1 $ by multiplying it with kernel
with anomalous dimension $ A_1 = \dfrac{11}{12} $. By resumming the large single and double transverse logarithms as in Ref. [36], the collinearly-improved version of BK evolution equation reads
Notably, Eqs. (9), (11), and (19) shall be numerically solved, and their solutions shall be used as dipole amplitudes to calculate the elastic and differential cross-sections in the following section. 22.3.Wavefunctions for vector meson -->
2.3.Wavefunctions for vector meson
Another ingredient to compute the differential cross-section for vector meson production is the overlap function $ (\Psi_{V}^{*}\Psi)_{T,L} $, which depends on the quark momentum fraction z, the dipole transverse size $ {{r}} $, and the photon virtuality $ Q^{2} $. The overlap function has various prescriptions, such as the boosted Gaussian, Gauss-LC, and DGKP [18]. Ref. [33] showed that for an identified meson, not all overlap functions provide an equally good description of the experimental data, and a meson has its own favorite wavefunction. We focus on studying the higher order effects for vector meson production in this study. Thus, we shall use an unified formalism of the wavefunction for different mesons to gain a better insight into the higher order effects. The overlap function between the photon and the vector meson has transverse and longitudinal components and can be written as [18]
In Eqs. (21) and (22), $ \phi(r,z) $ is the scalar function. In our study, the boosted Gaussian scalar functions are employed, since they works well for both light and heavy mesons [58]. In the boosted Gaussian formalism, the scalar functions are given by
The variable $ \epsilon $ in the Bessel functions in Eqs. (21) and (22) is $ \epsilon^{2} = z(1-z)Q^{2}+m_f^2 $. The values of the parameters $ M_V $, $ m_f $, $ N_{T,L} $, and $ R_{T,L} $ in the above equations are given in Table 1. It is worth noting that the longitudinal component is ignored in most studies due to its small contribution [54, 59], which is safe in a very small photon virtuality regime, like quasi-real photoproduction. However, the longitudinal component can provide a significant contribution in a large photon virtuality region, as we discussed in Ref. [60]. Thus, the longitudinal component is included in this study, as we compare with the data at various photon virtualities.
meson
$M_V/{\rm{GeV}}$
$m_f/{\rm{GeV}}$
$N_{T}$
$N_{L}$
$R_{T}/{\rm{GeV}}^{-2}$
$R_{L}/{\rm{GeV}}^{-2}$
$J/\psi$
3.097
1.4
0.578
0.575
2.3
2.3
$\phi$
1.019
0.14
0.919
0.825
11.2
11.2
Table1.Parameters of boosted Gaussian formalism for $J/\psi$ and $\phi$ [18].
-->
3.1.Numerical setup and data selection
The LO, rc and ci BK evolution equations are integral-differential equations. To obtain their numerical solutions, we can solve them on a lattice. In this study, we discretize the variable r into 256 points ($ r_{\rm min} = 2.06 \times 10^{-9} $ and $ r_{\rm max} = 54.6 $). Throughout this numerical study, the unit of dipole size r is $ {\rm{GeV}}^{-1} $. For the rapidity, the number of points are set to 100 with the step size $ \triangle Y = 0.2 $. This setup can ensure that the grid is sufficiently small for this purpose. To perform the numerical simulations, we employ the GNU scientific library (GSL). The main GSL subroutines we have used are Runge-Kutta for solving ordinary differential equations, the adaptive integration for numerical integrals, and the cubic spline interpolation for interpolating data points. To solve these integro-differential equations, initial conditions are required. There are several kinds of initial conditions in the literature, such as GBW [30] and MV [61] models. Refs. [44] and [45] use both of these models as initial conditions for the rcBK equation in the fit of the reduced cross-sections, they showed that the MV initial condition is significantly more favorable by the experimental data than the GBW model. Thus, we adopt the MV model as initial condition in this study [61],
with $ \gamma = 1.13 $, $ Q^2_{s0} = 0.15\, {\rm{GeV^2}} $, and $ \Lambda_{\rm QCD} = 0.241\, {\rm{GeV}} $. In our analysis, we use Eq. (3) to fit the differential cross-section and the elastic cross-section for $ J/\psi $ and $ \phi $ productions. The experimental data are taken from the ZEUS Collaboration ($ J/\psi $ [24], $ \phi $ [23]) and H1 Collaboration ($ J/\psi $ [26], $ \phi $ [25]). Notably, our studies are within the framework of the CGC, which is valid in the range $ x \leq x_{0} $ with $ x_0 = 10^{-2} $. Therefore, the data points with x larger than $ x_{0} $ are automatically excluded in the data set. Further, we exclude the data with large error bars at large photon virtuality. After selection, the total number amounts to 177 data points, which are used for the fit. For the details, the elastic cross-section data of $ J/\psi $ and $ \phi $ productions are 58 and 61 points, and the differential cross-section data of $ J/\psi $ and $ \phi $ productions are 24 and 34 points, respectively. 23.2.Fitting results with HERA data -->
3.2.Fitting results with HERA data
To demonstrate the significance of high order corrections in the description of the HERA data, one needs to compute the vector meson productions with the LO and NLO dipole amplitudes and compare these calculations. We concentrate on the study of the higher order effects for the vector meson production in this study. Therefore, there are only two free parameters B and R, as indicated by Eq. (8). The other parameters, such as $ \gamma $, $ Q_{s0}^2 $ in the initial condition in Eq. (25), are directly taken from Ref. [44], since it has been shown in Ref. [57] that the initial condition effects are eventually washed-out as the evolution developing with rapidity, so the solutions of the LO and NLO BK equations are insensitive to the choice of the initial condition at high rapidities. Tables 2 and 3 show these two parameters, and $ \chi^{2}/d.o.f $ results from our fit. From the values of the $ \chi^{2}/d.o.f $ in the last columns of Tables 2 and 3, the NLO descriptions of the vector meson productions are better than the LO case, which indicate that the NLO corrections play an important role in the diffractive process. In particular, the $ \chi^{2}/d.o.f $ resulting from the fit to the elastic cross-section, $ \sigma $, there is a large improvement as compared to the LO description once the NLO corrections are included. By global analysis, the values of $ \chi^{2}/d.o.f $ calculated from the rc and ci dipole amplitudes are closer to unity than those calculated from the LO amplitude.
${\cal{N}}(r,x)$
$B/{\rm{GeV}}^{-2}$
$R/{\rm{GeV}}^{-1}$
$\chi^{2}/d.o.f$
LO
2.500
3.800
2.159
rc
1.954
3.791
1.097
ci
2.060
3.737
0.932
Table2.Parameters and $\chi^{2}/d.o.f$ results for elastic cross-section with different dipole amplitudes.
${\cal{N}}(r,x)$
$B/{\rm{GeV}}^{-2}$
$R/{\rm{GeV}}^{-1}$
$\chi^{2}/d.o.f$
LO
2.253
3.223
2.056
rc
2.200
3.480
1.449
ci
2.175
3.349
1.357
Table3.Parameters and $\chi^{2}/d.o.f$ results for differential cross-section with different dipole amplitudes.
Figure 2 shows the elastic cross-sections $ \sigma $ for $ J/\psi $ and $ \phi $ productions as a function of the photon virtuality $ Q^{2} $. The dotted blue, dashed red, and solid black lines represent the results calculatied using the LO, rc, and ci dipole amplitudes, respectively (similarly hereinafter). For each meson, we consider the experimental data both from H1 and ZEUS collaborations. For $ J/\psi $, one can see that the higher order dipole amplitudes are in good agreement with the experimental measurement. For $ \phi $ production, it seems that all dipole amplitudes provide a similarly good description of the data in moderate $ Q^{2} $, however only the rc and ci amplitudes can provide a precise description of the data in low $ Q^{2} $. From Fig. 2, it is almost clear that the NLO amplitudes are more favored by the experimental data. Figure2. (color online) Elastic cross-section $\sigma$ for $J/\psi$ and $\phi$ as a function of $Q^{2}$.
The elastic cross-section $ \sigma $ for $ J/\psi $ and $ \phi $ productions as a function of photon-hadron center of mass energies $ W_{\gamma p} $ at different photon virtuality $ Q^{2} $ are shown in Fig. 3. The left panels of Fig. 3 show that the theoretical calculations from the NLO amplitudes are more consistent with the $ J/\psi $ data. From the right panels of Fig. 3, one can see that the LO calculations have a rather poor description of the experimental data, while the NLO computations provide a relatively good description of the data although the quality is not as good as the $ J/\psi $ case, since the experimental data for the $ \phi $ meson have large uncertainties. From Fig. 3 shows that the NLO calculations have a better agreement with experimental data than the LO BK equation for both $ J/\psi $ and $ \phi $. Figure3. (color online) Elastic cross-section $\sigma$ for $J/\psi$ and $\phi$ as a function of $W_\gamma p$ at different $Q^{2}$.
The differential cross-section $ {\rm d}\sigma/{\rm d}t $ for $ J/\psi $ and $ \phi $ as a function of the squared momentum transfer t at different photon virtuality $ Q^{2} $ are shown in Fig. 4. From Fig. 4, it seems that the LO and NLO calculations provide a similar quality description of the experimental data. This is due to the small dataset with large error bars. However, one can clearly see from the last column in Table 3 that the $ \chi^{2}/d.o.f $ computed from the NLO dipole amplitudes are significantly smaller than the ones from LO cases, which indicate that the NLO corrections take an effective role in the diffractive vector meson productions. Figure4. (color online) Differential cross-section ${\rm d}\sigma/{\rm d}t$ for $J/\psi$ and $\phi$ as a function of $t$ at different $Q^{2}$.
Furthermore, it should be noted that there is no significant difference between the description of the experimental data from the rcBK and ci BK calculations. To better interpret the underlying reasons, we have plotted the dipole amplitude, $ {\cal{N}}(r,x) $, as a function of the dipole size, r, for three different rapidities in Fig. 5. From Fig. 5, we find that the difference between the LO and the NLO (rcBK and ci BK) dipole amplitudes is evident. However, the difference between the amplitudes from rcBK and ci BK equations is miniscule up to large rapidities, i.e., $ Y = 5 $. The largest available rapidity at HERA is approximately $ Y = 5 $; therefore, it is almost impossible to discriminate the NLO running coupling effect from collinear resummations with current HERA data. This is the reason why we cannot see a remarkable difference between $ \chi^{2}/d.o.f $ resulting from running coupling and resummation improved dipole amplitudes. Figure5. (color online) Dipole amplitude for LO, rc, and ci BK evolution equations at three different rapidities.
23.3.Predictions for LHC -->
3.3.Predictions for LHC
The experimental data from LHC offer a peculiar way to test the hadronic structure, as the higher energy collision will touch even the small-x region. Fig. 5 shows the difference between the dipole amplitudes from rc and ci BK equations at larger rapidities, $ Y>5 $ (smaller-x region). Thus, the predictions for the LHC energies are meaningful as higher precision and rapidity data will be released by the LHCb collaboration. In the high energy proton-proton collisions, there are events involving interactions at large impact parameters, where the electromagnetic interaction is dominant. In these photon-induced processes, the two protons are kept intact after the interaction. For the total cross-section, this can be written in terms of a convolution of the equivalent photon flux and the photon-proton production cross-section. Therefore, the rapidity distribution for the exclusive vector meson production is given by
where y is the rapidity of the produced vector meson, $ \sigma_{\gamma p \rightarrow V p} $ is the total photon-proton cross section, and $ \omega $ is the photon energy ($ \omega_{\rm L} = \dfrac{M_{V}}{2} \exp(-y) $ and $ \omega_{\rm R} = \dfrac{M_{V}}{2} \exp(y) $). Note that there are two terms on the right hand side of the rapidity distribution equation. This is because the photon can be emitted either from the left or the right proton. In Eq. (26), $ \displaystyle\frac{{\rm d}N}{{\rm d}\omega} $ is the equivalent photon spectrum of the relativistic proton. In the Weisz$ \ddot{a} $cker-Williams approximation, this can be written as
where $ \xi = 1 + [\,(0.71 \,{\rm{GeV}}^2)/Q_{{\rm{min}}}^2\,] $ with $ Q_{{\rm{min}}}^2 \approx (\omega/\gamma_L)^2 $ at the high energy limit, $ \sqrt{s} $ is the proton-proton center of mass energy, and $ \gamma_L $ is the Lorentz factor. Moreover, the total photon-proton cross-section $ \sigma_{\gamma p \rightarrow V p} $ can be integrated from the differential cross-section in Eq. (3). The integral over t can be rewritten as follows
$ \sigma_{\gamma p \rightarrow V p} = \int_{-\infty}^{0}\frac{{\rm d}\sigma^{\gamma p\rightarrow Vp}}{{\rm d}t}{\rm d}t. $
(28)
Using the above formalism and the parameters obtained from fitting the HERA data, we can predict the rapidity distributions for diffractive $ J/\psi $ and $ \phi $ productions in proton-proton collisions at LHC energies. Figures 6 and 7 show our predictions for the rapidity distributions of exclusive $ J/\psi $ and $ \phi $ in proton-proton collisions at 7 TeV and 13 TeV, respectively. For the LHC kinematics region, there are possible rapidities whose corresponding Bjorken-x can be larger than $ x_0 $ for one of the protons, but still smaller than $ x_0 $ for the other proton. To obtain a smooth curve, we make a linear extrapolation for the dipole amplitudes when $ x > 0.01 $. We consider three kinds of dipole amplitudes (LO, rc, ci amplitudes) to calculate the rapidity distributions for exclusive vector meson productions and compare with the released data from LHCb [27, 28]. The numerical results in Figs. 6 and 7 show that the NLO dipole amplitudes provide better agrement with experimental data points, as expected. For completeness, we present our predictions of the total cross-section with different kinds of dipole amplitudes in Tables 4 and 5. From the tables, one can see that the production rates of the vector mesons ($ J/\psi $ and $ \phi $) are suppressed by the NLO effect, which satisfy theoretical expectations.
LO/nb
rc/nb
ci/nb
7 TeV
37.181
23.204
24.229
13 TeV
61.217
31.016
37.608
Table4.J/ψ total cross-section with different dipole evolution equations in pp collisions.
LO/nb
rc/nb
ci/nb
7 TeV
331.567
301.806
317.551
13 TeV
419.120
374.290
398.908
Table5.? total cross section with different dipole evolution equations in pp collisions.
Figure6. (color online) Predictions for rapidity distributions of $J/\psi$ and $\phi$ mesons in pp collisions at 7 TeV as a function of y.
Figure7. (color online) Predictions for rapidity distributions of $J/\psi$ and $\phi$ mesons in pp collisions at 13 TeV as a function of y.
In summary, we investigated the exclusive vector meson photoproduction for $ J/\psi $ and $ \phi $ at HERA in the framework of color glass condensate. By comparing the results from the rcBK and ci BK equations with those from the LO BK equation, we find that the results from NLO equations are more consistent with experimental data than the LO BK equation. We also present our predictions for the rapidity distributions in pp collisions by using parameters obtained from fitting the HERA data. These results indicate that the NLO effects are significant in the calculation of the vector meson production at LHC energies. Furthermore, the higher order corrections considered in this work are part of the NLO corrections to the BK evolution equation. As we have studied in Ref. [42], the rare fluctuations also require large corrections to the evolution equation once the gluon loop contributions are included into the rcBK equation. Therefore, the exclusive vector meson production with a rare fluctuation corrections is worth exploring in a future study.