Z⁰ boson associated b-jet production in high-energy nuclear collisions

The formation of quark-gluon plasma (QGP), which was produced in the early stages of the high-energy nucleus-nucleus collisions at the Relativistic Heavy Ion Collider (RHIC) and the Large Hadron Collider (LHC), offers a new possibility to test quantum chromodynamics (QCD) under an extremely hot and dense deconfined state of nuclear matter. The high-$ p_T $ partons (quarks and gluons) produced in the initial hard scattering strongly interact with the QGP and dissipate their energy to the medium, which is referred to as the jet quenching effect [1–5]. Consequently, the "quenched jet" observables are used to quantify the properties [6] of the hot and dense QCD matter by investigating their medium modifications in heavy-ion collisions (HIC) relative to their p+p baselines.

Recently, the associated production of a vector boson (photon γ or electroweak boson such as $ Z^0 $ and $ W^{+/-} $) and jets (V + jet) has been extensively studied both theoretically [7–11] and experimentally [12–20] to test the fundamental properties of QCD and improve the constraints on the parton distribution function (PDF) of a proton. More importantly, because the vector boson would not involve the strong interaction with the medium and gauge the initial energy of the tagged jets, V+jet is recognized as the ideal probe of the properties of QGP [21–38].

In particular, new measurements of the associated production of the $ Z^0 $ boson and b-jet (denoted as $ Z^0$ + b-jet) in p+p collisions at the LHC have been performed by ATLAS and CMS [39–45], since the final state b-jet associated with the $ Z^0 $ boson is the dominant background of the associated production of Higgs and $ Z^0 $ bosons ($ Z^0+H\rightarrow Z+b\bar{b} $) within the standard model (SM) [46] and can test many physics scenarios beyond the SM that predict new generation mechanisms of b quarks and $ Z^0 $ bosons [41]. It is noted that in HIC, the $ Z^0 $ boson tagged b-jet (as the initial energy of the b quark is well gauged by the vector boson and thus its energy loss can be directly obtained) is particularly suitable for exploring the quenching of the heavy flavor jet [47]. The "dead-cone" effect [48] of heavy quarks in QGP may lead to a smaller energy loss compared to light flavors, which is known as the mass effect of jet quenching and has attracted intensive investigations [49–54]. At the particle level, the latest measurements indicate that the yield of the B meson appears to be less suppressed than that of the D meson in nucleus-nucleus collisions [55–57]. However, no clear evidence was found at the full-jet level in the previous experimental measurements, such as the comparison of the $ R_{AA} $ between the inclusive jet and b-jet [58, 59] and the comparison of the $ p_T $ imbalance ($ x_J $) between inclusive dijets and $ b\bar{b} $ dijets [60], except for some preliminary indication in the recent measurement implemented by the ATLAS collaboration [61]. The possible reasons for this can be manifold; e.g., the large contribution of gluon-initiated b-parton processes may thwart the attempts at solving the problem. Previous studies [21, 47] have indicated that the dominant contribution of the $ Z^0 $ tagged jet is a quark-initiated jet, and the study of $ Z^0$+ b-jet in HIC, especially its different medium modifications compared with that of $ Z^0$ + jet, will provide a very useful tool to directly address the mass effect between the light-quark jet and massive bottom jet. Nevertheless, thus far, studies on the associated production of the b-jet and $ Z^0 $ boson in nucleus-nucleus collisions are lacking.

With this in mind, in this work, we present a Monte Carlo transport simulation including elastic (collisional) [62] and inelastic (radiative) [63–66] interactions of the energetic parton in the hot/dense QCD medium, while taking the next-to-leading order (NLO) plus parton shower (PS) generated initial hard parton spectrum [67] as input, to study the in-medium modification of the vector boson $ Z^0 $ tagged b-jets. This framework was employed to describe the heavy-flavor jet production of high-energy nuclear collisions in our previous studies [68–73]. We first present our numerical results for $ Z^0$ + jet and compare them with the available experimental data to test the applicability of our model. Then, we calculate the angular correlations of $ Z^0$ + b-jet in A+A collisions and demonstrate that the modifications of these correlations are sensitive to the initial b-jet $ p_T $ distribution instead of the azimuthal angle. In contrast to the case of $ Z^0$ + jet, the requirement of b-tagging excludes the contribution from multiple jets, so that the azimuthal angular correlations of $ Z^0$ + b-jet show distinct pattern modifications. With the high purity of the light-quark jet in $ Z^0$ + jet events [21, 47], we expect to address the mass dependence of the jet quenching effects between $ Z^0$ + jet and $ Z^0$ + b-jet.

The remainder of this paper is organized as follows. In Sec. II, we present the productions of $ Z^0$ + b-jet in p+p collisions calculated via the Monte Carlo event generator and comparisons with experimental data. In Sec. III, we discuss our treatments of the jet in-medium evolution in A+A collisions. In Sec. IV, we present the simulated results and discussions of the azimuthal angular correlation, transverse momentum imbalance, and nuclear modification factor of $ Z^0$ + b-jet in HIC. Finally, Sec. V summarizes the study.

Before we move into the study on $ Z^0$ + b-jet production in HIC, we should address its production in p+p collisions. Figure 1 shows a few processes [39] contributing to the associated production of a $ Z^0 $ boson and b-jet. In Figs. 1(a) and (b), an initial bottom quark from the parton distribution function (PDF) derived from the gluon distribution of one beam particle suffers the hard scattering and then turns into a b-jet and a $ Z^0 $ boson in the final state. In the bottom two diagrams, the $ b\bar{b} $ pairs originate from the hard scattering and then turn into two b-jets associated with an emitted $ Z^0 $ boson in the final state.

Figure 1.  (color online) Feynman diagrams contributing to the associated production of the b-jet with the $ Z^0 $ boson.

In this work, we use the MC@NLO event generator SHERPA [67] to obtain the initial $ Z^0$ + b-jet production in a p+p collision. The tree-level matrix elements are calculated using the internal modules Amegic [74] and Comix [75], and the one-loop virtual correction is calculated using the external program BlackHat [76]. The parton shower based on the Catani-Seymour [77] subtraction method is matched with NLO QCD matrix elements via the MC@NLO method [78]. The NLO PDF from NNPDF3.0 [79] with a 5-flavor scheme has been chosen in the calculations. FastJet [80] with an anti-$ k_{T} $ algorithm is used in the final-state jet reconstruction.

To compare our calculation results based on SHERPA with experimental data for p+p collisions, the same configurations implemented by the CMS collaboration [45] were used in our simulations. The $ Z^0 $ boson is reconstructed according to its decay channels $ Z^0\rightarrow e^{+}e^{-} $ and $ Z^0\rightarrow \mu^{+}\mu^{-} $. The transverse momentum of the electron and muon candidates is required to be larger than 20 GeV. To exclude the barrel-endcap transition region, the electrons are selected within the pseudorapidity region $ |\eta|<1.44 $ or $ 1.57<|\eta|<2.4 $, while muons are selected within $ |\eta|<2.4 $. According to the requirement of the experiment, the events are considered only when the invariant mass of the electron or muon pairs lies in the region $ 70<M_{ll}<111 $ GeV. The jets associated with the $ Z^0 $ boson are reconstructed by FastJet using the anti-$ k_T $ algorithm with a cone size of $ R=0.5 $. To reduce the contribution from the underlying event, the reconstructed jets must be in the pseudorapidity region $ |\eta^{{\rm{jet}}}|<2.4 $ and have $ p_{T,{\rm{jet}}}>30 $ GeV. The contribution from the underlying event is less than 5%, as estimated by CMS [45], because the production of softer jets is significantly suppressed by this requirement.

For events with at least one b-jet, we show the differential cross sections calculated via SHERPA as a function of the leading b-jet $ p_T $ and $ Z^0 $ boson $ p_T $ in Fig. 2. In the upper panel of Fig. 3, the azimuthal angular correlation between the leading b-jet and $ Z^0 $ boson ($ \Delta\phi_{bZ}=|\phi_{b}-\phi_{Z}| $) is compared between the calculation and the CMS data. In the lower panel of Fig. 3, for events with two b-jets (hereinafter $ Z^0\,+\,2$ b-jets), we plot the differential cross sections as a function of the azimuthal angle between the two b-jets ($ \Delta\phi_{bb}=|\phi_{b1}-\phi_{b2}| $). To quantify the deviation of our calculations from the experimental data points, we estimate the $ \chi^2/{\rm{d.o.f}} $ of these observables and present their values in the figures, where $ \chi^2=\sum_{i}\dfrac{[D_i-T_i]^2}{\delta_i^2} $; $ D_i $ and $ \delta_i $ are the center value and uncertainties, respectively, of the i-th experimental data point; $ T_i $ represents the theoretical value; and $ {\rm{d.o.f}} $ represents the number of compared data points. Our calculations based on SHERPA agree well with the experimental measurements.

Figure 2.  (color online) Differential cross section of $ Z^0$ + b-jet simulated by SHERPA (blue line) in the p+p collision at $ \sqrt{s}=8 $ TeV as a function of the transverse momentum of the highest-$ p_T $ b-jet (upper panel) and transverse momentum of the $ Z^0 $ boson (bottom panel) compared with the CMS data [45]; the $ \chi^2/{\rm{d.o.f}} $ of our fit to the CMS data is also presented in the plots.

Figure 3.  (color online) Differential cross section of $ Z^0$ + b-jet simulated by SHERPA (blue line) in the p+p collision at $ \sqrt{s}=8 $ TeV as a function of the azimuthal angular difference $ \Delta\phi_{bZ}=|\phi_{b}-\phi_{Z}| $ of the $ Z^0 $ boson and b-jet (upper panel) and azimuthal angular difference $ \Delta\phi_{bb}=|\phi_{b1}-\phi_{b2}| $ of the two b-jets (bottom panel) compared with the CMS data [45]; the $ \chi^2/{\rm{d.o.f}} $ of our fit to the CMS data is also presented in the plots.

In high-energy nuclear collisions, a droplet of an exotic state of nuclear matter, i.e., QGP, is expected to be formed. The high $ p_T $ partons produced in the hard scattering propagating in the QGP suffer both collisional and radiative energy loss as a result of the in-medium interaction. Numerous theoretical approaches [81–88] and Monte Claro models [89–96] have been developed in the last two decades to describe the heavy flavor meson production in HIC at the RHIC and LHC. Among them, the Langevin transport equations have been employed effectively to describe the heavy quark evolution in the expanding QCD medium [62, 68, 69, 97–100]. Because complete treatment of the heavy quark jets propagating in the QGP usually needs a simultaneous description of space-time evolution for both light and heavy partons [52, 87, 101, 102], we use the modified Langevin equations [62, 68, 69] to describe the propagation of heavy quarks in QGP and take into account the in-medium energy loss of light partons [63–66, 103, 104]. Additionally, owing to the lack of a unified theoretical approach that covers the full phase space of in-medium jet evolution in HIC, there are usually two methods to combine the vacuum parton shower with the medium-induced radiation [105, 106]. The first one is to introduce medium modifications on the vacuum parton splitting at a higher virtuality scale $ t\ge\hat{q}\tau_f $ (where $ \hat{q} $ is the jet transport parameter and $ \tau_f $ is the formation length of the radiated gluon), as in the Q-PYTHIA [107], MATTER [108], and JEWEL [109] models. In this work, we employ the alternative treatment implemented in LBT [110] and MARTINI [111] to simulate the jet energy loss at low virtuality $ t\le\hat{q}\tau_f $ and high energy $ E\gg q_{\perp} $ (where $ q_{\perp} $ represents the momentum exchange between the hard parton and the hot medium). According to this strategy, we take the p+p events produced by SHERPA with a full vacuum parton shower as the input, sample their initial spatial positions using the MC-Glauber model [112], and then simulate the subsequent in-medium evolution.

The movement of a heavy quark with large mass ($ M\gg T $) propagating in the hot/dense nuclear matter and suffering a large number of random kicks from the medium can be modeled as a Brownian motion [97]. Hence, a discrete Langevin equation can be utilized to describe the propagation of heavy quarks in the QCD medium [62, 68, 69, 90]:

where $ \Delta t $ represents the timestep in the Monte Carlo simulation, and $ \Gamma(p) $ is the drag coefficient representing the dissipation effect and controlling the strength of quasi-elastic scattering. $ \vec{\xi}(t) $ is the stochastic term that obeys a Gaussian probability distribution, i.e.,

and leads to

The diffusion coefficient κ is related to the drag coefficient Γ by the fluctuation-dissipation relation [113]:

where $ D_s $ is the spacial diffusion coefficient. The last term $ -\vec{p}_g $ is the recoil momentum due to the medium-induced gluon radiation, which is discussed in the following section. At each timestep, we boost partons to the local rest frame of the expanding medium to update the four-momentum and then boost them back to the laboratory frame to update the spatial position. Note that the procedure is done for $ T>T_c $, where $ T_c=165 $ MeV is the QCD transition temperature [114, 115]. The space-time evolution profile of the bulk medium in Pb+Pb collision is provided by the smoothed VISHNU [116] code. Even though the event-by-event fluctuation effects on the jet energy loss are small, the initial geometry fluctuation may be non-negligible for other observables such as particle collective flow ($ v_n $) [117–119].

Meanwhile, the calculation for leading logarithmic accuracy at Hard-Thermal-Loop approximation [103, 104] is employed in our framework to take into account the collisional energy loss of light quarks and gluons:

where L represents the transport path of the partons along the propagating direction, $ \alpha_s $ is the strong coupling constant, $ C_s $ is the quadratic Casimir in color representation, and $ \mu_{D} $ represents the Debye screening mass in the QCD medium. Note that Eq. (7) is only employed for light partons to consider their collisional energy loss because they cannot be treated as massive particles to evolve in the medium with the Langevin equations. During each timestep, the amount of collisional energy loss of a light parton can be calculated by integrating Eq. (7). Because medium-induced gluon radiation is the dominant energy loss mechanism for light partons, this treatment can be regarded as an effective approximation.

The inelastic scattering also plays an important role in the in-medium energy loss of energetic partons [120, 121]. In our work, the Higher-Twist (HT) radiated gluon spectra [63–66] is implemented to simulate the medium-induced gluon radiation when a parton propagates in the dense and hot QCD matter:

where x and $ k_\perp $ represent the energy fraction and the transverse momentum of the radiated gluon, respectively. $ \alpha_s $ is the strong coupling constant, which is fixed at $ \alpha_s=0.3 $ in our calculations, $ C_s $ is the quadratic Casimir in color representation, and $ P(x) $ is the splitting function [122] for the splitting processes $ q\rightarrow q+g $ and $ g\rightarrow g+g $ ($ g\rightarrow q+\bar{q} $ process is negligible owing to its low probability [100]).

$ \tau_f $ is the radiated gluon formation time defined as $ \tau_f=2Ex(1-x)/(k^2_\perp+x^2m^2) $, and $ t-t_i $ is the time interval between two instances of inelastic scattering. In the rigorous Higher-Twist calculations, the splitting function $ P(x) $ of heavy quarks should be mass dependent [65]. Nevertheless, as a simplified treatment in many transport models, such as LBT [93, 110], QLBT [123], and BAMPS [124, 125], it's convenient to uniformly write $ P(x) $ as a massless form for heavy and light quarks, while the dominant mass effect of the gluon radiation spectra of heavy quarks can be approximately presented as an overall factor $ (k_{\perp}^2/(k_{\perp}^2+x^2m^2))^4 $. Additionally, $ \hat{q} $ is the jet transport coefficient [126]:

where $ \hat{q}_{0} $ denotes the value of $ \hat{q} $ at the center of the bulk medium at the initial time $ \tau_0=0.6\; {\rm{fm/c}} $, and $\rho^{\rm QGP}(\tau,r)$ is the parton number density where the parton is probed. To take into account the radial flow effect [127], the four-momentum of parton $ p^{\mu} $ and the four flow velocity of the medium in the collision frame $ u^{\mu} $ act as a modification for $ \hat{q} $ in an expanding nuclear medium. We still hold that $ \hat{q} $/$ T^3={\rm{const}}. $ in the current framework as an effective approximation. More studies focusing on the puzzle of the temperature dependence of $ \hat{q}/T^3 $ can be found in the Refs. [6, 92, 128–130]. The last term of Eq. (8) represents the dead-cone effect [51, 65] which suppresses the gluon radiation of heavy quarks at a small angle ($ \theta\sim M/E $) owing to their large mass. We also noticed the up-to-data development of the HT approach in Ref. [131] achieved by including both transverse and longitudinal momentum exchanges between hard partons and the QCD medium and will consider their results in our framework in the future.

An imposed cut-off of the radiated gluon energy fraction $ x_{{\rm{min}}}=\dfrac{\mu_{D}}{E} $ is taken to avoid the divergence near $ x\rightarrow 0 $, where $ \mu_D^2=4\pi\alpha_s(1+\dfrac{n_f}{6}) T^2 $ is the Debye screening mass induced by the QGP medium. Therefore, following the method introduced in [93], one can estimate the mean number of radiated gluons$ \langle N(t,\Delta t)\rangle $ during a timestep $ \Delta t $ by integrating the phase space of x, $ k_{\perp} $, and t in Eq. (8):

By assuming that the multiple gluon radiation is a Poisson process, we obtain the probability distribution of the radiation number $ P(n,t,\Delta t) $ during a timestep, as well as the total inelastic scattering probability $P_{\rm rad}(t,\Delta t)$:

In our Monte Carlo simulation, during every timestep, the $P_{\rm rad}(t,\Delta t)$ is first evaluated to determine whether the radiation occurs. If accepted, the Possion distribution function $ P(n,t,\Delta t) $ is used for the sampling of the radiated gluon number. Finally, the four-momentum of the radiated gluon can be sampled according to the spectrum ${\rm d}N/{\rm d}x{\rm d}k_{\perp}^{2}$ expressed in Eq. (8). In Fig. 4, for a consistent comparison between our Monte Carlo simulation and the analytical calculation, we estimate the collisional and radiative energy losses of the gluons, light quarks, charm, and bottom in a static medium ($ T=400 $ MeV). Here we fix the parton energy (50 GeV) at each evolution timestep and also restore the initial time $ t_i $ in Eq. (8) to be 0 (same as the treatment in Ref. [91]) because their variations during the Monte Carlo simulation are not automatically included in the analytical calculation. We find that the MC results agree well with the analytical calculations. For the collisional energy loss, the $ \Delta E_{\rm coll} $ of the gluon is $ 9/4 $ times that of the light quark owing to the large color factor, and in our framework, the $\Delta E_{\rm coll}$ of the heavy quark is comparable to that of the light quark. For the radiative energy loss, a clear mass hierarchy for different parton species can be found: $\Delta E^g_{\rm{rad}} > \Delta E^q_{\rm{rad}} > \Delta E^c_{\rm{rad}} > $ $ \Delta E^b_{\rm{rad}}$. For a long propagation time $ t=4 $ fm, we find that the radiative energy loss dominates the total parton energy loss because of its quadratic dependence on the path length.

Figure 4.  (color online) Collisional (a) and radiative (b) energy losses of gluons, light quarks, charm, and bottom with initial energy $ E_0=50$ GeV in a static medium with temperature $ T=400$ MeV. The Monte Carlo simulations are compared with the semi-analytical calculations.

In general, there are two parameters in our framework that need to be determined: the jet transport coefficient $ \hat{q} $ and the diffusion coefficient $ D_s $. We treat $ \hat{q} $ and κ as two independent parameters to be constrained by experimental data. First, the value of $ \hat{q} $ is determined via global extraction of the single hadron production in Pb+Pb collisions [132], in which $ q_0=1.2 $ GeV$ ^2 $/fm is obtained at the LHC energy. After $ \hat{q} $ is fixed, we extract the best value $ D_s(2\pi T)\sim4 $ via $ \chi^2 $ fitting to the D meson $ R_{AA} $ data [133, 134], which is consistent with the results of $ D_s(2\pi T)=3.7\sim7 $ reported by the Lattice QCD [135].

In this section, to estimate the medium modification of jet observables in nucleus-nucleus collisions, we use the p+p events provided by SHERPA as the input of our simulation within the hydrodynamic background to study the in-medium jet evolution. Before proceeding to the $ Z^0 $ tagged b-jet, we calculate the azimuthal angular correlation ($ \Delta\phi_{jZ}=|\phi_{{\rm{jet}}}-\phi_Z| $) and transverse momentum imbalance ($ x_{jZ}=p_T^{{\rm{jet}}}/p_T^Z $) of $ Z^0$ + jet, as well as the nuclear modification factor $ R_{AA} $ of the inclusive b-jet, and compare our theoretical results with the available experimental data. Then, we calculate the $ Z^0$ + b-jet observables, including the azimuthal angular correlation between the $ Z^0 $ boson and b-jet ($\Delta\phi_{bZ}=|\phi_{b{\rm{-jet}}}-\phi_Z|$), angle separation between the two Z-tagged b-jets ($ \Delta\phi_{bb}=|\phi_{b1}-\phi_{b2}| $), transverse momentum ($x_{bZ}=p_T^{b{\rm{-jet}}}/p_T^Z$), and nuclear modification factor $ I_{AA} $ [21, 35] defined as

Here, $ \langle N_{{\rm{bin}}} \rangle $ denotes the average number of binary nucleon-nucleon collisions among A+A collisions calculated via the Glauber model [112].

In Fig. 5, we show our calculated $ \Delta\phi_{jZ} $ distributions for both p+p and 0%–30% Pb+Pb collisions at $ \sqrt{s_{NN}}$ = 5.02 TeV compared with the CMS experimental data [31]. Additionally, the $ \chi^2/{\rm{d.o.f}} $ of our model fitting to the CMS data for both p+p and Pb+Pb are presented. The same configurations in the jet reconstruction used by the CMS are employed. All final-state jets are reconstructed by FastJet using the anti-$ k_T $ algorithm with $ R=0.3 $ and require $ p_T^{{\rm{jet}}}>30 $ GeV. The selected $ Z^0 $ bosons are reconstructed by the electron or muon pairs based on their decay channels ($ Z^0\rightarrow e^{+}e^{-} $ and $ Z^0\rightarrow \mu^{+}\mu^{-} $) and require $ p_T^{Z}>60 $ GeV. Note that these distributions are normalized by the number of $ Z^0 $events, and the transverse momentum imbalances are subjected to Gaussian smearing [31] to take into account the detector resolution effects. This reveals that the distribution of azimuthal angular correlation in Pb+Pb collisions suffers a suppression in a small $ \Delta\phi_{jZ} $ region relative to the p+p baseline, which is consistent with the CMS measurement. However, in a large angle region ($ \Delta\phi_{jZ}\sim \pi $, where the $ Z^0 $ boson and jet are almost back-to-back), this suppression is not very apparent. The reason for this behavior has been discussed in detail [26, 34]; i.e., the small $ \Delta\phi_{jZ} $ region is dominated by the multiple-jet processes and the large $ \Delta\phi_{jZ} $ region is dominated by soft/collinear radiation. Usually, the jet energy of multiple-jet processes is relatively low and easier to be shifted below the jet selection threshold (30 GeV) because of parton energy loss [27]. We also notice two significant differences between the p+p and Pb+Pb CMS data near $ \Delta\phi_{jZ}=\pi/2 $ and $ \Delta\phi_{jZ}=7\pi/8 $, and our results cannot fit these points well. However, similar abnormal behaviors of the CMS data are not found in the measurements of the γ-jet [136]; hence, we guess that they may be caused by the statistical fluctuations in the experiment.

Figure 5.  (color online) Distributions of the azimuthal angle difference $ \Delta\phi_{jZ} $ between the $ Z^0 $ boson and the jet both in p+p and 0%–30% Pb+Pb collisions at 5.02 TeV compared with CMS data [31]. The distributions are scaled by the number of $ Z^0 $ events N in p+p collisions. The values of $ \chi^2/{\rm{d.o.f}} $ of our fit to the CMS data for both p+p and Pb+Pb are also presented in the plots.

In Fig. 6, we compute $ x_{jZ} $ distribution for $ Z^0$ + jet in both p+p and 0%–30% Pb+Pb collisions at $ \sqrt{s_{NN}}$ = 5.02 TeV compared with the CMS data. The $ \chi^2/{\rm{d.o.f}} $ are also presented in the figures, indicating that our calculations are consistent with the experimental data, but the p+p baseline needs to be improved. Note here that selected $ Z^0$ + jet pairs are required to be almost back-to-back ($ \Delta\phi_{jZ}>7\pi/8 $). Relative to the p+p baseline, in Pb+Pb collisions, we find that the $ x_{jZ} $ distribution is shifted toward smaller values, exhibiting an enhancement at $ 0<x_{jZ}<0.7 $ and suppression at $ 0.7<x_{jZ}<2 $. As $ x_{jZ} $ represents the transverse momentum imbalance of $ Z^0 $ and the jet, it can be understood that for each $ Z^0$ + jet pair, the values of $ x_{jZ} $ decrease owing to the jet energy loss and thus are shifted to a smaller $ x_{jZ} $ observed in the final state.

Figure 6.  (color online) Distributions of the transverse momentum balance $ x_{jZ} $ of $ Z^0$ + jet in both p+p and 0%–30% Pb+Pb collisions at 5.02 TeV compared with CMS data [31]. The distribution is normalized by the number of $ Z^0 $ events and $ Z^0$ + jet pairs are required with $ \Delta\phi_{jZ}>7\pi/8 $. The values of $ \chi^2/{\rm{d.o.f}} $ of our fit to the CMS data for both p+p and Pb+Pb are also presented in the plots.

In Fig. 7, we investigate the nuclear modification factor $ R_{AA} $ of the inclusive b-jet in Pb+Pb collisions at $ \sqrt{s_{NN}}=2.76 $ TeV in comparison with the CMS measurements [58] to test our model calculations. The values of $ \chi^2/{\rm{d.o.f}} $ are shown in each panel of Fig. 7 to quantify the deviation of our calculations from the CMS data for different centrality bins. We find that our calculations are essentially consistent with the experimental data, but we may need to pay attention to some discrepancies. First, our theoretical results show weak $ p_T $ dependence, but the experimental $ R_{AA} $ appears to increase visibly with the jet $ p_T $ in the 10%–30% and 30%–50% centrality bins. Second, some data points ($ p_T=100 $ GeV at 0%−10% and 10%−30%, $ p_T=210 $ GeV at 30%–50% and 50%–100%) cannot be well described by our calculations.

Figure 7.  (color online) Nuclear modification factor $ R_{AA} $ of the b-jet. The p+p baseline is provided by SHERPA, and the theoretical calculations are compared with CMS data [58] at centralities of 0%−10%, 10%−30%, 30%−50%, and 50%−100%. The values of $ \chi^2/{\rm{d.o.f}} $ of our fit to the CMS data are also presented in the plots.

The good agreement between our model calculations and the data of the $ Z^0$ + jet and inclusive b-jet makes it possible to study the medium modification of the $ Z^0$ + b-jet in nuclear-nuclear collisions. In Fig. 8, we calculate the azimuthal angular correlation of the $ Z^0 $ boson and b-jet in p+p and 0%–10% Pb+Pb collisions at $ \sqrt{s_{NN}}$ = 5.02 TeV. The b-jets associated with the $ Z^0 $ boson are reconstructed by FastJet using the anti-$ k_T $ algorithm with a cone size of $ \Delta R=0.5 $, $ |\eta^{{\rm{jet}}}|<2.4 $, and $ p_T^{{\rm{jet}}}>30 $ GeV for both p+p and Pb+Pb collisions. These distributions are normalized by the initial $ Z^0$ + b-jet event number (in p+p collision) to address the medium modification. We observe an overall suppression in Pb+Pb collisions relative to the p+p baseline. We show their ratio PbPb/pp in the middle panel of Fig. 8 and find that the suppression for the $ Z^0$ + b-jet has a far weaker dependence on $ \Delta\phi_{bZ} $ than that for the $ Z^0$ + jet, which exhibits stronger suppression in a small $ \Delta\phi_{jZ} $ region where multiple jets dominate. All the selected jets must first be b quark tagged, and this requirement significantly reduces the contribution from multiple-jet processes when we consider the azimuthal angular ($ \Delta\phi_{bZ} $) distribution.

Figure 8.  (color online) Upper panel: distributions of the azimuthal angular correlation of the $ Z^0 $ boson and b-jet in both p+p and 0%−10% Pb+Pb collisions at 5.02 TeV; the distributions are scaled by the number of $ Z^0$ + b-jet events N in p+p collisions. Middle panel: ratio of the azimuthal angle correlations in Pb+Pb to p+p. Bottom panel: averaged b-jet transverse momentum $ \langle p_{T} \rangle $ as a function of $ \Delta\phi_{bZ} $

To address the key factor that leads to the flat suppression on $ \Delta\phi_{bZ} $ distribution, we estimate the averaged b-jet transverse momentum $ \langle p_{T} \rangle $ as a function of $ \Delta\phi_{bZ} $, which can be calculated by:

The decrease of the selected event number in A+A collisions results from the in-medium energy loss, which shifts lower $ p_T $ jet below the kinematic selection cut. The initial $ \langle p_{T} \rangle $ distribution actually reflects the $ \Delta\phi_{bZ} $ dependence of this shift. It turns out that the distribution of $ \langle p_{T} \rangle $ versus $ \Delta\phi_{bZ} $ is nearly a constant value at 55 GeV as shown in the bottom panel of Fig. 8, which leads to the rather flat suppression on $ \Delta\phi_{bZ} $ distribution. Here the band of $ \langle p_{T} \rangle $ distribution represents the statistical standard errors in the simulations.

As mentioned in Sec. II, the azimuthal angular separation $ \Delta\phi_{bb} $ of the two b-jets tagged by $ Z^0 $ boson is also a useful observable to distinguish the contribution from subprocesses where $ Z^0 $ boson is emitted from one of the final state b quark or gluon splitting ($ g\rightarrow b\bar{b} $) [41], as shown in diagrams (c) and (d) of Fig. 1. Note that these two categories of contributions corresponding to the cases that the two b-jets are almost back-to-back or collinear. What interests us is how the $ \Delta\phi_{bb} $ distribution of these two categories of $ Z^0\,+\,2$ b-jets would be modified in the QGP. As shown in the top panel of Fig. 9, we plot the $ \Delta\phi_{bb} $ distributions both in p+p and 0−10% Pb+Pb collisions at $ \sqrt{s_{NN}} $=5.02 TeV and also plot the ratio PbPb/pp in the middle panel. In the upper panel, we can find a kink at $ \Delta\phi_{bb}\sim\pi/5 $ and a peak at $ \Delta\phi_{bb}\sim\pi $ in Fig. 9, which present the two contributions of $ Z^0\,+\,2$ b-jet production: the two b-jets are almost back-to-back or collinear. The two-peaks distribution is similar to what we have observed in the angular correlations of $ b\bar{b} $ dijets [68], and the peaks at smaller and larger $ \Delta\phi_{bb} $ regions corresponding to the gluon splitting (GSP) processes and flavor creation (FCR) processes. In the middle panel, we observe an upward trend of the ratio from 0.5 to 0.7 as $ \Delta\phi_{bb} $ increases. To figure out the $ \Delta\phi_{bb} $ dependence of the ratio PbPb/pp, we estimate the initial jet $ \langle p_{T} \rangle $ of the leading and sub-leading b-jet in p+p collisions, as shown in the bottom panel of Fig. 9. We find that, for both the leading and subleading one, $ \langle p_{T} \rangle $ is increasing with $ \Delta\phi_{bb} $, which shows a similar trend with that of the ratio PbPb/pp versus $ \Delta\phi_{bb} $. It may indicate that the medium modification of $ \Delta\phi_{bb} $ in Pb+Pb has a close connection with the initial b-jet $ \langle p_T \rangle $ distribution versus $ \Delta\phi_{bb} $ in p+p collisions.

Figure 9.  (color online) Distributions of the azimuthal angular separation $ \Delta\phi_{bb} $ of the two b-jets tagged by $ Z^0 $ boson both in p+p and 0−10% Pb+Pb collisions at 5.02 TeV, the distributions are scaled by the $ Z^0\,+\,2$ b-jets event number N in p+p collisions. Middle panel: ratio of the azimuthal angular separation in Pb+Pb to that in p+p. Bottom panel: averaged transverse momentum of the leading and the sub-leading b-jet as a function of $ \Delta\phi_{bb} $.

The associated production of $ Z^0$ + jet may shed new light on the mass dependence of the jet quenching effect in nuclear matter, owing to the high purity of light-quark-initiated jets. In Ref. [47], the contributions from light-quark-jets and gluon-jets in the $ Z^0$ + jet production are approximately 70% and 30%, respectively. To verify this point, in Fig. 10, we estimate the gluon-jet fraction in four categories of jets in p+p collisions at 5.02 TeV: inclusive jet, inclusive b-jet, $ Z^0 $ tagged jet, and $ Z^0 $ tagged b-jet. The gluon-jet in the events can be identified by requiring that the gluon is the leading parton in the jets, which can be easily implemented in the FastJet program. We find that at $ p_T^{{\rm{jet}}}\sim 50 $ GeV, the gluon-jet fraction is approximately 50% for the inclusive jet and 30% for the $ Z^0 $ tagged jet. The $ Z^0 $-tagging requirement considerably decreases the gluon jet contribution (by 40%), especially at a low $ p_T $. More importantly, for the inclusive b-jet and $ Z^0 $ tagged b-jet, the contributions from the gluon-jet are significantly suppressed owing to the requirement of b-quark tagging and show almost equal values.

Figure 10.  (color online) Gluon initiated jet fraction as a function of the transverse momentum of the inclusive jet (red dash), inclusive b-jet (blue dash), $ Z^0$ + jet (red solid), and $ Z^0$ + b-jet (blue solid).

We plot the scaled $ x_{jZ} $ distributions of $ Z^0$ + jet in both p+p and 0%–10% Pb+Pb collisions at $ \sqrt{s_{NN}} $=5.02 TeV in Fig. 11(a) and, with the same configuration, the $ x_{bZ} $ distributions of $ Z^0$ + b-jet in Fig. 11(b). All the selected jets (b-jets) must satisfy $ \Delta\phi_{jZ}>7\pi/8 $ ($ \Delta\phi_{bZ}>7\pi/8 $) to guarantee that they are almost back-to-back with the $ Z^0 $ boson. We observe that the distributions of $ x_J $ shift toward smaller values for both $ Z^0$ + jet and $ Z^0$ + b-jet in Pb+Pb collisions relative to their p+p baselines, owing to the energy loss of the tagged jets. To perform a more intuitive comparison between $ Z^0$ + jet and $ Z^0$ + b-jet, we show the difference of the $ x_J $ distribution in Pb+Pb from that in p+p (see the bottom panels of Figs. 11(a) and (b)), which have positive values at $ 0.2<x_{J}<0.8 $ and negative values at $ 0.8<x_{J}<1.6 $. Then, we find that the absolute value of the difference (PbPb-pp) of $ Z^0$ + jet is larger than that of $ Z^0$ + b-jet. Furthermore, we estimate the shifting of the mean value of momentum imbalance $ \Delta {\langle x_{jZ} \rangle } ={\langle x_{jZ} \rangle }_{pp}-{\langle x_{jZ} \rangle }_{PbPb} $ for $ Z^0$ + jet, and $\Delta {\langle x_{bZ} \rangle }= {\langle x_{bZ} \rangle }_{pp}-{\langle x_{bZ} \rangle }_{PbPb}$ for $ Z^0$ + b-jet, as shown in Table 1.

Figure 11.  (color online) (a) $ x_{jZ} $ distributions of $ Z^0$ + jet in both p+p and 0%−10% Pb+Pb collisions at 5.02 TeV, as well as their difference in Pb+Pb to p+p. The distributions are normalized by the $ Z^0$ + jet number in both p+p and Pb+Pb collisions. (b): $ x_{bZ} $ distributions of $ Z^0$ + b-jet in both p+p and 0%−10% Pb+Pb collisions at 5.02 TeV, as well as their difference in Pb+Pb to p+p. The distributions are normalized by the $ Z^0$ + b-jet number in both p+p and Pb+Pb collisions.

	$ Z^0 $+ jet	$ Z^0$ + b-jet
$ {\left\langle x_{J} \right\rangle }_{pp} $	0.987±0.0047	0.941±0.0056
${\left\langle x_{J} \right\rangle }_{\rm PbPb}$	0.851±0.0061	0.849±0.0064
$ \Delta \left\langle x_{J} \right\rangle $	0.136±0.0108	0.092±0.012

Table 1.  Mean value of momentum imbalance $ x_J $ of $ Z^0 $ + jet and $ Z^0 $ + b-jet in both p+p and 0%–10% Pb+Pb collisions at $ \sqrt{s_{NN}}=5.02 $ TeV, as well as the shifting of the mean value of momentum imbalance $\Delta x_{J}={\left\langle x_{J} \right\rangle }_{pp}-{\left\langle x_{J} \right\rangle }_{\rm PbPb}$. The standard errors of $ x_J $ in the simulations are also presented.

Here, J denotes different processes. Note that the standard errors of $ x_J $ in our simulations are also presented in Table 1. We find that within the statistical uncertainties, $ \Delta\langle x_{jZ} \rangle $ ($ \sim0.136\pm0.0108 $) for $ Z^0$ + jet is visibly larger than $ \Delta\langle x_{bZ} \rangle $ ($ \sim0.092\pm0.012 $) for $ Z^0$ + b-jet, indicating stronger modifications to the light-quark jet compared with the b-jet. Note that $ x_{jZ} $ ($ x_{bZ} $) is defined by the ratio of jet (b-jet) $ p_T $ to $ Z^0 $ boson $ p_T $. The shifting of $ x_{jZ} $ ($ x_{bZ} $) toward smaller values is directly related to the amount of jet (b-jet) energy loss; hence, $ \Delta\langle x_{jZ} \rangle>\Delta\langle x_{bZ} \rangle $ indicates that light quark jets lose more energy than b-jets.

The nuclear modification factor $ I_{AA} $ is practically another good observable to address the mass hierarchy and flavor dependence of the jet quenching effect. Comparisons of $ I_{AA} $ between $ Z^0$ + jet and $ Z^0$ + b-jet would provide more reliable evidence of the mass effect of jet quenching. For this purpose, we present the calculations of $ I_{AA} $ of $ Z^0$ + jet and $ Z^0$ + b-jet in 0%−10% Pb+Pb collisions at $ \sqrt{s_{NN}} $=5.02 TeV as a function of jet $ p_T $ within three $ p^Z_T $ windows in Fig. 12. First, we find different shapes of $ I_{AA} $ in these three panels. In the left panel ($ 40<p_T^Z<60$ GeV), the two curves of $ I_{AA} $ for $ Z^0$ + jet and $ Z^0$ + b-jet are flat, but in the right panel ($ 80<p_T^Z<120$ GeV), the curves show enhancement at $ p_T^{{\rm{jet}}}<80 $ GeV. This is because if we constrain $ 80<p_T^Z< 120 $ GeV in the event selection, the cross section at $ p_T^{{\rm{jet}}}<80 $ GeV is steeper, falling with jet $ p_T $. Then, the jets shifting from a higher $ p_T $ to a lower $ p_T $ because of the in-medium energy loss naturally lead to the relatively large $ I_{AA} $ values in the lower $ p_T^{\rm jet} $ region, even larger than one. Additionally, we observe smaller values of $ I_{AA} $ for Z + jet relative to those for Z + b-jet in the three panels of Fig. 12, which indicates that the yield of $ Z^0 $ tagged light-quark jets suffers stronger suppression after traversing the QCD matter than that of $ Z^0 $ tagged b-jets. Note that $ I_{AA} $ is directly related to the $ p_T $ shifting of the tagged jets due to in-medium energy loss. The smaller value of $ I_{AA} $ for $ Z^0 $ + jet indicates a larger jet $ p_T $ shift compared with that of $ Z^0$ + b-jet, suggesting that the energy loss of $ Z^0 $ tagged b-jets is smaller than that of $ Z^0 $ tagged light-quark jets. Our conclusions are consistent with the recent ATLAS measurements [61], which indicate that the $ R_{AA} $ of the b-jet is higher than that of the inclusive jet. We hope that our complementary predictions can be tested via future measurements at the LHC, which may be helpful for solving the puzzle of the mass hierarchy of jet quenching.

Figure 12.  (color online) Nuclear modification factor as a function of the transverse momentum of $Z^0 $ + jet (blue line) and $Z^0$ + b-jet (red line) within three $p_T^Z$ ranges, i.e., 40−60 GeV, 60−80 GeV, and 80−120 GeV, in 0%−10% centrality Pb+Pb collisions at $\sqrt{s_{NN}}=5.02$ TeV.

Vector boson-tagged heavy quark jets are promising new tools for studying the jet quenching effect. In this work, we present a Monte Carlo transport simulation, which takes into account the elastic and inelastic jet interactions within a hydrodynamic background, to study the in-medium modification of $ Z^0 $ tagged b-jets. The NLO+PS event generator SHERPA was used to provide the p+p baseline of $ Z^0$ + b-jet production, which agreed well with the CMS measurements. This framework has been proven to give good descriptions of medium modifications of the $ \Delta\phi_{jZ} $ and $ x_{jZ} $ of $ Z^0$ + jet, as well as the $ R_{AA} $ of the inclusive b-jet, measured in Pb+Pb collisions at the LHC.

The angular correlation between the vector boson and heavy quark-tagged jets may be a new promising observable for studying the in-medium jet interaction. We present the first calculation of the azimuthal angular correlation $ \Delta\phi_{bZ} $ of $ Z^0$ + b-jet in both p+p and 0%−10% Pb+Pb collisions at $ \sqrt{s_{NN}}$ = 5.02 TeV. We observe a flat suppression factor versus $ \Delta\phi_{bZ} $, in contrast to the case of $ Z^0$ + jet, because the requirement of b-tagging excludes the contribution from multiple-jet processes. Additionally, we calculate the medium modification of the azimuthal angular correlation $ \Delta\phi_{bb} $ in central Pb+Pb collisions at $ \sqrt{s_{NN}}=5.02 $ TeV and observe stronger suppression in a smaller $ \Delta\phi_{bb} $ region of the distribution relative to that at $ \Delta\phi_{bb}\sim\pi $. By analyzing the $ \langle p_T \rangle $ distribution of the tagged b-jets, we find that the medium modification pattern on $ \Delta\phi_{bZ} $ ($ \Delta\phi_{bb} $) in Pb+Pb has a close connection with the initial $ \langle p_T \rangle $ distribution versus $ \Delta\phi_{bZ} $ ($ \Delta\phi_{bb} $) in p+p collisions. These investigations may help us to understand the experimental measurements of jet angular correlations at the LHC in recent years, e.g., for the γ + jet and $ Z^0$ + jet.

Finally, we predict that the mass effect of jet quenching can be addressed by comparing the medium modifications of $ Z^0$ + jet and $ Z^0$ + b-jet. With the high purity of the quark jet in $ Z^0$ + (b-) jet events, we estimate the medium modification of the transverse momentum imbalance $ x_{jZ} $ ($ x_{bZ} $) and the nuclear modification factor $ I_{AA} $ for both $ Z^0$ + jet and $ Z^0$ + b-jet in Pb+Pb collisions. We find a larger shift of $ x_{jZ} $ and stronger suppression of $ I_{AA} $ for $ Z^0$ + jet than for $ Z^0$ + b-jet, indicating that b-jets lose less energy than light quark jets. These predictions can be tested via future measurements at the LHC and may provide a key to solving the puzzle of the mass hierarchy of jet quenching.