1 Introduction

Angular correlations between two particles have been established as a powerful tool to study the properties of the system created in high energy collisions of hadrons and nuclei [116]. These measurements are usually performed in a two dimensional space as a function of \(\Delta \eta \) and \(\Delta \varphi \). Here \(\Delta \eta \) and \(\Delta \varphi \) are the differences in pseudorapidity \(\eta = -{\mathrm {ln}}[\tan (\theta /2)]\) (where \(\theta \) is the polar angle of a particle relative to the beam axis) and in azimuthal angle \(\varphi \) of the two particles.

In heavy-ion collisions at both the Relativistic Heavy Ion Collider (RHIC) [311] and at the Large Hadron Collider (LHC) [1216], these correlations exhibit characteristic structures: (a) a peak at (\(\Delta \eta \),\(\Delta \varphi \)\(=\) (0, 0), usually referred to as the near-side jet peak, resulting from intra-jet correlations as well as correlation due to decay of resonances and quantum statistics correlations, (b) an elongated structure over \(\Delta \eta \) at \(\Delta \varphi \) \(=\) \(\pi \) originating partially from correlations between particles from back-to-back jets and from collective effects such as anisotropic flow, and (c) a similar component at \(\Delta \varphi \) \(=\) 0 extending to large values of \(\Delta \eta \), usually called the near-side ridge, whose origin was subject of a theoretical debate [1731]. Although initially the near-side ridge was also attributed to jet–medium interactions [1720], it is now believed to be associated to the development of collective motion [2431] and to initial state density fluctuations, including the initial state effects within the framework of the Color Glass Condensate (CGC) [2123].

Similar structures have recently been reported in two-particle correlation analyses in smaller systems. In particular, the CMS Collaboration, by studying angular correlations between two particles in \(\Delta \eta \) and \(\Delta \varphi \), reported the development of an enhancement of correlations on the near-side (i.e. \(\Delta \varphi \) \(=\) 0) in high- compared to low-multiplicity pp collisions at \(\sqrt{s} = 7\) TeV that persists over large values of \(\Delta \eta \) [32]. In the subsequent data taking periods at the LHC, similar ridge structures were observed on both the near- and the away-side in high-multiplicity p–Pb collisions at \(\sqrt{s_{\mathrm {NN}}}\) \(=\) 5.02 TeV [3338]. The origin of these effects, appearing in small systems, is still debated theoretically. In particular, it was suggested in [3941] that in high-multiplicity collisions the small system develops collective motion during a short hydrodynamic expansion phase. On the other hand, in [4244] the authors suggested that the ridge structure can be understood within the CGC framework.

The ALICE Collaboration also reported a particle mass ordering in the extracted \(v_{2}\) (i.e. the second coefficient of the Fourier expansion of the azimuthal distribution of particles relative to the symmetry plane) values for \(\pi ^{\pm }\), \({\mathrm {K}}^{\pm }\), and p(\(\overline{{\mathrm {p}}}\)) in high-multiplicity p–Pb collisions [45]. This mass ordering becomes evident once the correlations observed in the lowest multiplicity class are subtracted from the ones recorded in the highest multiplicity class. The ordering is less pronounced, yet still present, if this subtraction procedure is not applied. Similar mass ordering in Pb–Pb collisions [46] is usually attributed to the interplay between radial and elliptic flow induced by the collective motion of the system. These observations in p–Pb collisions were reproduced by models incorporating a hydrodynamic expansion of the system [47, 48]. Recently, it was suggested in [49] that the signatures of collective effects observed in experiments could be partially described by models that couple the hot QCD matter created in these small systems, described as an ensemble of non-interacting particles, to a late stage hadronic cascade model. More recently, the CMS Collaboration demonstrated that the effects responsible for the observed correlations in high-multiplicity p–Pb events are of multiparticle nature [50]. This strengthens the picture of the development of collective effects even in these small systems.

The charge-dependent part of two-particle correlations is traditionally studied with the balance function (BF) [51], described in detail in Sect. 4. Such studies have emerged as a powerful tool to probe the properties of the system created in high energy collisions. Particle production is governed by conservation laws, such as local charge conservation. The latter ensures that each charged particle is balanced by an oppositely-charged partner, created at the same location in space and time. The BF reflects the distribution of balancing charges in momentum space. It is argued to be a sensitive probe of both the time when charges are created [51, 52] and of the collective motion of the system [26, 53]. In particular, the width of the balance function is expected to be small in the case of a system consisting of particles that are created close to the end of its evolution and are affected by radial flow [26, 5153]. On the other hand, a wide balance function distribution might signal the creation of balancing charges at the first stages of the system’s evolution [26, 5153] and the reduced contribution or absence of radial flow.

In this article, we extend the previous measurements [54] by reporting results on the balance function in pp, p–Pb, and Pb–Pb collisions at \(\sqrt{s_{\mathrm {NN}}} = 7\), 5.02, and 2.76 TeV, respectively. The data were recorded with the ALICE detector [5557]. The results are presented as a function of multiplicity and transverse momentum (\(p_{{\mathrm {T}}}\)) to investigate potential scaling properties and similarities or differences between the three systems. The article is organized as follows: Sect. 2 briefly describes the experimental setup, while details about the data sample and the selection criteria are introduced in Sect. 3. In Sect. 4, the analysis technique and the applied corrections are illustrated. In Sect. 5, the specifics about the estimation of the systematic uncertainties are described. Section 6 discusses the results followed by a detailed comparison with models to investigate the influence of different mechanisms (e.g. unrelated to hydrodynamic effects) on the balance functions. In the same section, the comparison of the results among the three systems is presented.

2 Experimental setup

ALICE [57] is one of the four major detectors at the LHC. It is designed to efficiently reconstruct and identify particles in the high-particle density environment of central Pb–Pb collisions [58, 59]. The experiment consists of a number of central barrel detectors positioned inside a solenoidal magnet providing a 0.5 T field parallel to the beam direction, and a set of forward detectors. The central detector systems of ALICE provide full azimuthal coverage for track reconstruction within a pseudorapidity window of \(|\eta | < 0.9\). The experimental setup is also optimized to provide good momentum resolution (about \(1~\%\) at \(p_{{\mathrm {T}}}~< 1\) GeV/c) and particle identification (PID) over a broad momentum range [60].

For this analysis, charged particles were reconstructed using the Time Projection Chamber (TPC) [61] and the Inner Tracking System (ITS) [57]. The TPC is the main tracking detector of the central barrel [61], consisting of 159 pad rows grouped into 18 sectors that cover the full azimuth within \(|\eta | < 0.9\). The inner and outer radii of the detector are 85 and 247 cm, respectively. The ITS consists of six layers of silicon detectors employing three different technologies. The two innermost layers, positioned at \(r = 3.9\) and 7.6 cm, are Silicon Pixel Detectors (SPD), followed by two layers of Silicon Drift Detectors (SDD) at \(r = 15\) and 23.9 cm. Finally, the two outermost layers are double-sided Silicon Strip Detectors (SSD) at \(r = 38\) and 43 cm.

A set of forward detectors, the V0 scintillator arrays [62], were used in the trigger logic and the multiplicity determination. The V0 consists of two systems, the V0A and the V0C, positioned on both sides of the interaction point along the beam. They cover the pseudorapidity ranges \(2.8 < \eta < 5.1\) and \(-3.7 < \eta < -1.7\) for the V0A and the V0C, respectively.

For more details on the ALICE detector setup and its performance in the LHC run 1, see [57, 60].

3 Analysis details

This analysis is based on data from pp, p–Pb, and Pb–Pb collisions. The data were recorded for pp collisions during the 2010 run at \(\sqrt{s} = 7\) TeV, for p–Pb collisions during the 2013 run at \(\sqrt{s_{\mathrm {NN}}} = 5.02\) TeV, and for Pb–Pb collisions during the 2010 and 2011 runs at \(\sqrt{s_{\mathrm {NN}}} = 2.76\) TeV. In p–Pb collisions, the nucleon–nucleon centre-of-mass system was shifted with respect to the ALICE laboratory system by a rapidity of \(-\)0.465 in the direction of the proton beam. For simplicity, the pseudorapidity in the laboratory frame is denoted, throughout this article, with \(\eta \) for all systems (note that for pp and Pb–Pb collisions the laboratory and the centre-of-mass systems coincide).

Minimum-bias p–Pb and Pb–Pb events were triggered by the coincidence between signals from the two sides of the V0 detector. For the pp run, the minimum-bias trigger definition was modified to require at least one hit in the SPD  or either of the V0 detectors. In addition, for Pb–Pb, an online selection based on the V0 detectors was used to increase the number of events with high multiplicity. An offline event selection exploiting the signal arrival time in V0A and V0C, with a 1 ns resolution, was used to discriminate background (e.g. beam-gas) from collision events. This led to a reduction of background events in the analyzed samples to a negligible fraction (\({<}0.1~\%\)) for all systems [60]. All events retained in the analysis had a reconstructed primary vertex position along the beam axis (\(z_{vtx}\)) within 10 cm from the nominal interaction point. Finally, events with multiple reconstructed vertices were rejected, leading to a negligible amount of pile-up events for all systems [60].

After all the selection criteria, approximately \(240\times 10^6\), \(100\times 10^6\), and \(35\times 10^6\) events were analyzed for pp, p–Pb, and Pb–Pb, respectively.

Tracks are reconstructed from a collection of space points (clusters) inside the TPC. The tracking algorithm, based on the Kalman filter, provides the quality of the fit by calculating its \(\chi ^2\) value. Each space-point is reconstructed at one of the TPC padrows, where the deposited ionization energy is also measured. The specific ionization energy loss (\(\langle {\mathrm {d}}E/{\mathrm {d}}x \rangle \)) is estimated by averaging this ionization over all clusters associated to the track. The procedure has an uncertainty, which we later refer to as \(\sigma _{{\mathrm {d}}E/{\mathrm {d}}x}\).

To select primary tracks with high efficiency and to minimize the contribution from background tracks (i.e. secondary particles originating either from weak decays or from the interaction of particles with the detector material), all selected tracks were required to have at least 70 reconstructed space points out of the maximum of 159 possible in the TPC. In addition, the \(\chi ^2\) per degree of freedom per TPC space point of the momentum fit was required to be below 2. To further reduce the contamination from background tracks, only tracks with a distance of closest approach (DCA) to the primary vertex in both the xy-plane (\({\mathrm {DCA}}_{\mathrm {xy}}\)) and the z coordinate (\({\mathrm {DCA}}_{\mathrm {z}}\)) below a threshold value (i.e. \({\mathrm {DCA}}_{\mathrm {xy}} < 2.4\) cm and \({\mathrm {DCA}}_{\mathrm {z}} < 3.0\) cm) were analyzed. These requirements lead to a reconstruction efficiency of about \(80~\%\) for primary particles and a contamination from secondaries of about \(5~\%\) at \(p_{{\mathrm {T}}} = 1\) GeV/c [63] in pp collisions. The efficiency is similar in p–Pb collisions and it is lower by about 3–5 \(\%\) in central Pb–Pb collisions, according to detailed Monte Carlo simulations. In addition, electrons originating from \(\gamma \)-conversion and \(\pi ^0\)-Dalitz decays were removed based on the energy loss \(({\mathrm {d}}E/{\mathrm {d}}x)\) measured by the TPC. Tracks for which the measured \({\mathrm {d}}E/{\mathrm {d}}x\) lied within \(3\sigma _{{\mathrm {d}}E/{\mathrm {d}}x}\) of the Bethe–Bloch parametrization of \(\langle {\mathrm {d}}E/{\mathrm {d}}x \rangle \) for electrons and at least \(3\sigma _{{\mathrm {d}}E/{\mathrm {d}}x}\) away from the relevant parametrizations for pions, kaons, and protons, were removed.

All particles were reconstructed within \(|\eta | < 0.8\). This selection excludes possible biases from the tracking efficiency that becomes lower for \(|\eta | > 0.8\) as compared to \(|\eta | < 0.8\). The particles selected in this analysis have a transverse momentum in the range \(0.2 < p_{{\mathrm {T}}} < 15.0\) GeV/c.

In order to reduce the contribution from track splitting (i.e. incorrect reconstruction of a signal produced by one track as two tracks) and merging (i.e. two nearby tracks being reconstructed as one track) in the active volume of the TPC, a selection based on the closest distance of two tracks in the TPC volume was applied when forming particle pairs. This was done by excluding pairs with a minimum pseudorapidity difference of \(|\Delta \eta |<0.02\) and angular distance \(|\Delta \varphi ^{*}|<0.02\) rad. Here \(\Delta \varphi ^*\) is the angular distance between two tracks, accounting also for their curvature due to their charge, according to:

$$\begin{aligned} \Delta \varphi ^{*}=\varphi _{1}-\varphi _{2}-\alpha _1+\alpha _2, \end{aligned}$$
(1)

where \(\varphi _{1}\) and \(\varphi _{2}\) are the azimuthal angles of the two tracks at the vertex, and \(\alpha _{i}\) (with \(i = 1,2\)) is given by

$$\begin{aligned} \alpha _{i} = q_i\left| \arcsin \left( \frac{0.0075B_{z}(\mathrm {T}) r(\mathrm {cm})}{p_{\mathrm {T}i}(\mathrm {GeV/}c)}\right) \right| \end{aligned}$$
(2)

In Eq. 2, \(q_{1}\) and \(q_{2}\) stand for the charge of each track, \(B_{z}\) is the magnetic field in the z direction, r corresponds to the radius of the smallest distance of the tracks in the detector used (\(0.8 < r < 2.5\) m with a step of \(\Delta r = 0.2\) cm, for the TPC) and \(p_{\mathrm {T}1}\) and \(p_{\mathrm {T}2}\) are the transverse momentum values of the two particles forming the pair.

3.1 Multiplicity classes in pp, p–Pb, and Pb–Pb collisions

The analyzed events were divided into multiplicity classes using the V0A detector. Since this detector does not provide any tracking information, the amplitude of the signal from each cell, which is proportional to the number of particles that hit a cell, was used as a proxy for multiplicity [64]. The choice of the V0A as the default multiplicity estimator was driven by the fact that in p–Pb collisionsFootnote 1 this detector is located in the direction of the Pb–ion and thus is sensitive to its fragmentation [64]. In addition, this choice allowed for reducing autocorrelation biases introduced when the multiplicity class was estimated in the same \(\eta \) range as the one used to measure correlations. For consistency, the same multiplicity estimator was used for the other two systems. For the V0 detectors, a calibration procedure [60, 62] (i.e. gain equalization) was performed to account for fluctuations induced by the hardware performance, and for the different conditions of the LHC machine for each running period.

For each multiplicity class, the raw transverse momentum spectrum for charged particles with \(p_{{\mathrm {T}}} > 0.2\) GeV/c reconstructed in \(|\eta | < 0.8\) was extracted. These raw spectra were corrected for detector acceptance and efficiency using Monte Carlo simulations with PYTHIA [65], DPMJET [66], and HIJING [67] event generators for pp, p–Pb, and Pb–Pb, respectively. The ALICE detector response for these events was determined using a GEANT3 [68] simulation. In addition to the reconstruction efficiency, a correction related to the contamination from secondaries originating from weak decays and from the interaction of particles with the material of the detector was applied. This correction was estimated with both the aforementioned simulations and also using a data-driven method, based on fitting the DCA distributions with templates extracted from Monte Carlo for primary particles and secondaries originating either from weak decays or from the interaction of other particles with the detector material, as described in [69]. The resulting corrected charged-particle multiplicity was calculated by integrating the corrected transverse momentum spectrum over the region with \(p_{{\mathrm {T}}} > 0.2\) GeV/c.

Table 1 presents the multiplicity classes in terms of percentage of the multiplicity distribution, and the corresponding number of charged particles with \(p_{{\mathrm {T}}}\) \(~> 0.2\) GeV/c reconstructed at \(|\eta | < 0.8\) for all three systems. The resulting values for \(\mathrm {N}_{\mathrm {charged}}\) are subject to an overall tracking efficiency uncertainty of \(4~\%\) [70].

Table 1 Corrected mean charged particle multiplicities (for \(p_{{\mathrm {T}}}\) \(> 0.2\) GeV/c, and \(|\eta | < 0.8\)) for event classes defined by the percentage of the V0A multiplicity distribution for pp, p–Pb, and Pb–Pb collisions at \(\sqrt{s_{\mathrm {NN}}} = 7\), 5.02, and 2.76 TeV, respectively

4 Balance function

The charge-dependent correlations are studied using the balance function [51] for pairs of charged particles with angular differences \(\Delta \eta \) and \(\Delta \varphi \). For each pair, the first (“trigger”) particle has a transverse momentum \(p_{\mathrm {T,trig}}\), while the second (“associated”) charged particle has a transverse momentum \(p_{\mathrm {T,assoc}}\).

The associated yield per trigger particle is then calculated for different charge combinations. For one charge combination \((+,-),\) it is defined as

$$\begin{aligned} c_{(+,-)}=\frac{1}{N_{{\mathrm {trig}},+}}\frac{{\mathrm {d}}^2 N_{\mathrm {assoc},-}}{{\mathrm {d}}\Delta \eta {\mathrm {d}}\Delta \varphi } = S_{(+,-)}/f_{(+,-)} \end{aligned}$$
(3)

and similarly for the other charge combinations. The signal \(S_{(+,-)}= 1/N_{\mathrm {trig},+}{\mathrm {d}}^2N_{\mathrm {same},(+,-)}/{\mathrm {d}}\Delta \eta {\mathrm {d}}\Delta \varphi \) is constructed from the number of positive trigger particles \(N_{\mathrm {trig},+}\) and the particle pair distribution \({\mathrm {d}}^2N_{\mathrm {same},(+,-)}/{\mathrm {d}}\Delta \eta {\mathrm {d}}\Delta \varphi \), formed in \(\Delta \eta \)-\(\Delta \varphi \) with positive and negative particles from the same event. Both terms are corrected for detector inefficiencies and contamination from secondary particles on a track-by-track basis, using the corrections described in Sect. 3.1 as an inverse weight. \(S_{(+,-)}\) is computed after summing separately over all events the two components \(N_{\mathrm {trig},+}\) and \({\mathrm {d}}^2N_{\mathrm {same},(+,-)}/{\mathrm {d}}\Delta \eta {\mathrm {d}}\Delta \varphi \).

The background distribution \(f_{(+,-)}= \alpha {\mathrm {d}}^2N_{mixed,+-}/{\mathrm {d}}\Delta \eta {\mathrm {d}}\Delta \varphi \) corrects for particle pair-acceptance. It is constructed by combining a trigger particle from one event with associated particles from other events. This procedure is known as the event mixing technique. These mixed pairs are formed from events having the same multiplicity classes and \(z_{vtx}\) within \(\pm 2\) cm of each other. Each trigger particle is mixed with associated particles from at least 5 events. The coefficient \(\alpha \) in Eq. 3 is used to normalize the mixed-event distribution to unity in the \(\Delta \eta \) region of maximal pair acceptance. Finally, the associated yield per trigger particle is computed by calculating the weighted-average of the corresponding yields for several intervals of \(V_z\). This is done to account for the different pair acceptance and efficiency as a function of \(V_z\).

The balance function is then defined as the difference of the associated yields per trigger particle for unlike and like-sign combinations [51], according to

$$\begin{aligned} B(\Delta \eta ,\Delta \varphi ) = \frac{1}{2}[c_{(+,-)} + c_{(-,+)} - c_{(+,+)} - c_{(-,-)}] \end{aligned}$$
(4)

The resulting two-dimensional distributions are projected separately onto \(\Delta \eta \) and \(\Delta \varphi \) and the widths, \(\sigma _{\Delta \eta }\) and \(\sigma _{\Delta \varphi }\), are calculated as the standard deviation of the distributions. In this analysis, the projection in \(\Delta \eta \) is done on the near- (\(-\pi /2 < \Delta \varphi < \pi /2\)) and on the away-side (\(\pi /2<\Delta \varphi <3\pi /2\)), separately.

Three transverse momentum intervals are used in the analysis: the low (\(0.2 < p_{\mathrm {T,assoc}} < p_{\mathrm {T,trig}} < 2.0\) GeV/\(c\)), intermediate (\(2.0 < p_{\mathrm {T,assoc}} < 3.0 < p_{\mathrm {T,trig}} < 4.0\) GeV/\(c\)), and high (\(3.0 < p_{\mathrm {T,assoc}} < 8.0 < p_{\mathrm {T,trig}} < 15.0\) GeV/\(c\)) \(p_{{\mathrm {T}}}\) regions. Note that the integral of the balance function reported in this article does not reach unity but rather 0.5 due to the requirement imposed on the \(p_{{\mathrm {T}}}\) of the “trigger” and the “associated” particles.

For \(0.2 < p_{\mathrm {T,assoc}} < p_{\mathrm {T,trig}} < 2.0\) GeV/\(c\), the width in \(\Delta \eta \) and \(\Delta \varphi \) is calculated in \(|\Delta \eta |<1.6\) and \(-\pi /2<\Delta \varphi <\pi /2\). For higher values of transverse momentum, the balance function distributions are fitted with a sum of a Gaussian and a constant. The width is then calculated within \(3\sigma _{{\mathrm {Gauss}}}\), with \(\sigma _{{\mathrm {Gauss}}}\) extracted from the Gaussian of the aforementioned fit. The statistical error of the width is calculated using the subsample method [71, 72]. The values of \(\sigma _{\Delta \eta }\) and \(\sigma _{\Delta \varphi }\) are calculated for each subsample (maximum 10 subsamples were used) and the statistical uncertainty is estimated from the spread of these independent results.

5 Systematic uncertainty

In all figures except Fig. 1, the data points are plotted with their statistical and systematic uncertainties indicated by error bars and open boxes around each point, respectively. The systematic uncertainty was obtained by varying the event, track, and pair selection criteria, as will be explained in the following paragraphs. The contribution of each source was calculated as the spread of the values of each data point, extracted from variations of the selection criteria. If statistically significant, each contribution was added in quadrature to obtain the final systematic uncertainty. Following this procedure, the resulting maximum values of the systematic uncertainty over all multiplicity classes and systems for the balance function projections in \(\Delta \eta \) and \(\Delta \varphi \) were less than 5 \(\%\). In what follows, we report the maximum systematic uncertainties over all multiplicity classes for each system for \(\sigma _{\Delta \eta }\) and \(\sigma _{\Delta \varphi }\).

The Pb–Pb data samples were analyzed separately for two magnetic field configurations. The difference of \(1.5~\%\) in the results was taken as a systematic uncertainty. For all systems, different LHC periods, reflecting different machine conditions and detector configurations (e.g. non-working channels), were analyzed separately. The corresponding maximum systematic uncertainties over all multiplicity classes was \(1.1~\%\). Furthermore, the influence on the results of different tracking strategies was studied by repeating the analysis using tracks reconstructed by the combination of signals from the TPC and the ITS. The relevant maximum systematic uncertainties from this source were 1.2, 0.2, and \(1.2~\%\) for pp, p–Pb, and Pb–Pb, respectively. Finally, the contribution coming from the V0 gain equalization in pp collisions was investigated by equalizing the signal per V0 ring, per channel, and per detector. The study did not reveal any systematic differences in the obtained results.

Fig. 1
figure 1

The balance function \(B(\Delta \eta ,\Delta \varphi \)) for charged particles with \(0.2 < p_{\mathrm {T,assoc}} < p_{\mathrm {T,trig}} < 2.0\) GeV/\(c\) in Pb–Pb, p–Pb, and pp collisions at \(\sqrt{s_{\mathrm {NN}}} = 2.76\), 5.02, and 7 TeV, respectively. From top to bottom the \(0{-}5~\%\) for Pb–Pb and \(0{-}10~\%\) for p–Pb and pp collisions, \(30{-}40~\%\), and the \(70{-}80~\%\) multiplicity classes are shown

In addition, several of the track quality criteria defined by the tracking algorithm described in Sect. 3 were varied. The uncertainty related to the electron rejection criterium was studied by varying the requirement on the expected Bethe–Bloch parameterization of the momentum dependence of \({\mathrm {d}}E/{\mathrm {d}}x\) for electrons from \(3\sigma \) to \(5\sigma \). This contribution was negligible in the pp system, while it was 0.1 and \(0.2~\%\) for p–Pb and Pb–Pb, respectively. The requirement on the closest distance of two tracks of a pair in the TPC was varied from \(\Delta \eta = 0.01\) to \(\Delta \eta = 0.03\) and from \(\Delta \varphi ^{*} = 0.01\) rad to \(\Delta \varphi ^{*} = 0.03\) rad. This source was found to yield negligible systematic uncertainty for the pp system, while the maximum contribution for p–Pb and Pb–Pb systems were 0.2 and \(0.7~\%\), respectively. The systematic uncertainty of the track-by-track correction for efficiency and contamination was estimated from Monte Carlo simulations. For this, the results of the analysis of a sample at the event generator level (i.e. without invoking either the detector geometry or the reconstruction algorithm) were compared with the results of the analysis over the output of the full reconstruction chain, using the corrections for detector inefficiencies and acceptance discussed in Sect. 3. This source resulted into a partially correlated uncertainty of around \(0.4~\%\) for the case of pp and p–Pb, and \(1.1~\%\) for the Pb–Pb system.

The resulting values for the systematics are summarized in Table 2, for all systems. The table provides the maximum value for every source over all multiplicity classes and transverse momentum ranges.

Table 2 The maximum value of the systematic uncertainties on the width of the balance function over all multiplicity classes for each of the sources studied for the pp, p–Pb and Pb–Pb systems

Finally, different multiplicity estimators were used to study the variations coming from the multiplicity class definition. There was no systematic uncertainty assigned for this contribution. The results obtained with the two forward detectors (e.g. V0A and V0C) show no significant difference. On the other hand, a slightly weaker narrowing of the balance function with increasing multiplicity is observed when the central barrel detector is used for both measuring the correlations and the multiplicity class definition, in the pp and p–Pb systems. These differences are coming from physics processes (e.g. back-to-back jets), whose contribution is reduced if one defines multiplicity classes using a detector located further away from mid-rapidity. This also justifies the reason why the V0A detector was chosen as the multiplicity estimator in this analysis.

6 Results

6.1 Balance function in the low transverse momentum region

Figure 1 presents the balance function for charged particles in \(\Delta \eta \) and \(\Delta \varphi \) for three multiplicity classes of Pb–Pb, p–Pb, and pp collisions at \(\sqrt{s_{\mathrm {NN}}} = 2.76\), 5.02, and 7 TeV, respectively. From top to bottom the results for the highest (i.e. 0–5 \(\%\) for Pb–Pb collisions and 0–10 \(\%\) for p–Pb and pp collisions), intermediate (i.e. 30–40 \(\%\)), and lowest (i.e. 70–80 \(\%\)) multiplicity classes are shown. The trigger and associated particles are selected from the low transverse momentum region \(0.2 < p_{\mathrm {T,assoc}} < p_{\mathrm {T,trig}} < 2.0\) GeV/\(c\). The bulk of the charge-dependent correlation yield is located on the near-side (\(-\pi /2<\Delta \varphi <\pi /2\)). In this region, the balance function becomes narrower with increasing multiplicity for all three collision systems. The peak values of the balance function also change with multiplicity, with higher values corresponding to collisions with higher multiplicity. On the away-side (\(\pi /2 < \Delta \varphi < 3\pi /2\)), the balance function has a larger magnitude for lower multiplicity events. In addition, a depletion in the correlation pattern around \((\Delta \eta ,\Delta \varphi )=(0,0)\) starts to emerge in mid-central (e.g. 30–40 \(\%\) multiplicity class) events in Pb–Pb collisions and becomes more pronounced in p–Pb and pp collisions with decreasing multiplicity. The origin of this structure will be discussed later.

Fig. 2
figure 2

The balance function for charged particles with \(0.2 < p_{\mathrm {T,assoc}} < p_{\mathrm {T,trig}} < 2.0\) GeV/\(c\) as a function of \(\Delta \eta \) on the near-side (upper row) and away-side (middle row) and \(\Delta \varphi \) (lower row) in different multiplicity classes of Pb–Pb in panels a, d and g, p–Pb in panels b, e and h, and pp collisions in panels c, f and i at \(\sqrt{s_{\mathrm {NN}}} = 2.76\), 5.02, and 7 TeV, respectively

The integral of the balance function over the acceptance is related to measures of charge fluctuations as argued in [52], and is between 0.25 and 0.35 (i.e. 0.5 and 0.7 in case the \(p_{{\mathrm {T}}}\) requirement between the “trigger” and the “associated” particles is not imposed) for all systems and multiplicity classes. For each system it reveals a mild multiplicity class dependence which, for Pb–Pb, could explain the increase of multiplicity fluctuations for central compared to peripheral events reported in [73].

6.1.1 Balance function projections

Figure 2 presents for Pb–Pb, p–Pb, and pp collisions the projections of the two-dimensional balance function in \(\Delta \eta \)  on the near-side (panels (a), (b), (c) ) and away-side (in panels (d), (e), (f)), and \(\Delta \varphi \) in panels (g), (h), (i), respectively. The statistical uncertainty, usually smaller than the marker size, is represented by the error bar while the systematic uncertainty, calculated as the quadratic sum of the correlated and the uncorrelated part, by the box around each data point. The balance function as a function of the relative pseudorapidity difference \(\Delta \eta \) on the near-side exhibits a strong multiplicity dependence for all collision systems. In particular, the distribution narrows and the peak value becomes larger for high- compared to low-multiplicity events. As a function of the relative azimuthal angle \(\Delta \varphi \) on the near-side, the balance function exhibits the same qualitative features as for \(\Delta \eta \), i.e. narrower distributions with larger magnitude for increasing event multiplicity in all three systems. However, the magnitude of the balance function on the away-side exhibits a different trend, with larger values of B(\(\Delta \eta \)) and B(\(\Delta \varphi \)) measured for low- compared to high-multiplicity events.

As already discussed in Sect. 3, in p–Pb collisions, the nucleon–nucleon centre-of-mass system shifts by a rapidity of \(-\)0.465 with respect to the ALICE laboratory system in the direction of the proton beam. The influence of this shift was studied with simulations and, although the balance function is not translationally-invariant, the shift does not lead to any significant difference in either the projections of the balance function or the extracted widths.

As indicated previously, starting from mid-central events in Pb–Pb collisions a distinct depletion is observed in the two-dimensional distribution around \((\Delta \eta ,\Delta \varphi )=(0,0)\) that becomes more pronounced in events with low multiplicities, and in particular in p–Pb and pp collisions. The fact that the aforementioned depletion does not seem to be restricted to a very narrow window in either \(\Delta \eta \) or \(\Delta \varphi \) (the structures extend to \(-0.4 < \Delta \eta < 0.4\) and \(-\pi /6 < \Delta \varphi < \pi /6\)) indicates that the origin is not due to detector effects, as was confirmed by independent studies involving modification of cuts controlling track splitting and merging. One possible mechanism that could create such a structure is the charge-dependent short-range correlations such as Coulomb attraction and repulsion, or quantum statistics correlations [7476]. To test this hypothesis, a criterium on the minimum transverse momentum difference \(\Delta p_{{\mathrm {T}}}\) between two particles of a pair was applied. The value was varied from \(\Delta p_{{\mathrm {T}}}~=~0\) GeV/c to \({\Delta }p_{{\mathrm {T}}}~=~0.2\) GeV/c. The choice for the selected values is driven by the fact that the bulk of short-range correlations are expected to have \({\Delta }p_{{\mathrm {T}}}<0.1\) GeV/c [77]. The depletion is less pronounced with increasing value of \({\Delta }p_{{\mathrm {T}}}\) and vanishes for \({\Delta }p_{{\mathrm {T}}}~=~0.2\) GeV/c. The disappearance of the depletion was also achieved by increasing the lower transverse momentum threshold for both the trigger and the associated particle to \(p_{{\mathrm {T}}} > 0.5\) GeV/c. Both these observations are inline with the hypothesis that the depletion originates from (mainly) quantum statistics correlations and Coulomb effects. The physics conclusion, i.e. narrower distributions with increasing event multiplicity, does not change applying one of these criteria.

6.1.2 Comparison with models

In Fig. 3a, d, g the balance function in \(\Delta \eta \) on the near- (a) and away-side (d), and in \(\Delta \varphi \) (g) are compared with Monte Carlo calculations using the HIJING [67] and AMPT [78, 79] event generators. The figures show the 0–5 \(\%\) multiplicity class of Pb–Pb collisions. In AMPT simulations, the string melting option was used, with parameters tuned to describe the experimental data on anisotropic flow at LHC energies [80, 81]. The centrality classes were defined based on the module of the impact parameter. It is seen that neither AMPT nor HIJING are able to describe the balance function projections in \(\Delta \eta \) on the near-side (see Fig. 3a), since they expect not only much broader distributions but they also underestimate the magnitude of the balance function. On the other hand, the projection of the balance function in \(\Delta \eta \) on the away-side (Fig. 3d) indicates that AMPT is in qualitative agreement with the data points, contrary to HIJING that predicts a significantly larger magnitude of the balance function. Finally, the \(\Delta \varphi \) projection of the balance function in Fig. 3g shows that while HIJING is still not able to describe the data points, AMPT predicts narrower distributions on the near-side but with a much smaller magnitude than the one experimentally measured.

Fig. 3
figure 3

The balance function for charged particles with \(0.2 < p_{\mathrm {T,assoc}} < p_{\mathrm {T,trig}} < 2.0\) GeV/\(c\) as a function of \(\Delta \eta \) on the near-side (upper row) and away-side (middle row) and as a function of \(\Delta \varphi \) (lower row) for Pb–Pb (panels a, d and g), p–Pb (panels b, e and h) and pp collisions (panels c, f and i) compared with results from various event generators. Only the highest multiplicity class is shown, i.e. 0–5 % for Pb–Pb and 0–10 % for p–Pb and pp collisions

The comparison of the experimental results for the 0–10 \(\%\) multiplicity class in p–Pb collisions with model predictions is presented in Fig. 3b, e, h. For this comparison, results from Monte Carlo calculations using the DPMJET [66] and AMPT [78, 79] event generators were used. DPMJET is a model based on independent pp collisions, describing hard processes, hadron–hadron interactions, and hadronic interactions involving photons, without any collective effects. This model fails to describe the experimental data points in either of the two projections in \(\Delta \eta \), i.e. on the near- and the away-side in Fig. 3b, e, respectively, expecting much broader distributions with smaller (larger) magnitude on the near-(away-)side. In addition, for the balance function projection in \(\Delta \varphi \) presented in Fig. 3h, DPMJET predicts broader distributions with a smaller magnitude compared to the measured data points on the near-side, but also exhibits a correlation peak on the away-side contrary to what is observed experimentally. On the other hand, AMPT, as in the case of the Pb–Pb collisions, seems to describe better the balance function projections in both \(\Delta \eta \) and \(\Delta \varphi \).

For pp collisions, the experimental results are compared with two variants of calculations using PYTHIA8 tune 4C [82] in Fig. 3c, f, i. This tune contains modified multi-parton interaction (MPI) parameters that allow it to describe the multiplicity dependence of \(\langle p_{\mathrm {T}} \rangle \) [63]. The default calculation includes the color reconnection mechanism, which is switched off in the second configuration. The version of PYTHIA8 without the inclusion of color reconnection expects a broader balance function near-side projection in \(\Delta \eta \) with a smaller magnitude than the one measured as observed in Fig. 3c. On the other hand, the same tune predicts larger magnitude than the one measured for the balance function away-side projection in \(\Delta \eta \) (see Fig. 3f). Finally, for the projection in \(\Delta \varphi \), this tune expects significantly broader distributions on the near-side than the measured ones, with an extra correlation peak developing on the away-side which is not observed experimentally. On the other hand, the tune of PYTHIA8 with the inclusion of color reconnection describes the experimental measurement fairly well in both \(\Delta \eta \) and \(\Delta \varphi \) projections.

As discussed in the previous paragraphs, there are models that exhibit a correlation peak on the away-side contrary to what is supported by the data. For this reason, the width of the balance function distribution in \(\Delta \eta \) and \(\Delta \varphi \) will be extracted and compared with models on the near-side only.

Fig. 4
figure 4

The multiplicity-class dependence of \(\sigma _{\Delta \eta }\) in Pb–Pb, p–Pb, and pp collisions at \(\sqrt{s_{\mathrm {NN}}} = 2.76\), 5.02, and 7 TeV compared with results from various event generators in panels a, c, and e. Panels b, d, and f show the relative decrease of \(\sigma _{\Delta \eta }\) calculated with respect to \(\sigma _{\Delta \eta }^{70{-}80\,\%}\), as a function of the multiplicity class. The transverse momentum values for both the trigger and the associated particles satisfy the condition \(0.2 < p_{\mathrm {T,assoc}} < p_{\mathrm {T,trig}} < 2.0\) GeV/\(c\)

6.1.3 Balance function width

To quantify the narrowing of the balance function width as a function of multiplicity, the standard deviation \(\sigma \) is calculated as described in Sect. 4. The panels (a), (c), and (e) of Fig. 4 present the evolution of \(\sigma _{\Delta \eta }\) on the near-side with multiplicity class, expressed by the multiplicity percentile for Pb–Pb, p–Pb, and pp collisions, respectively. Note that the multiplicity decreases from left to right along the horizontal axis. The statistical uncertainties of the data points are represented by the error bars and are usually smaller than the marker size. For all collision systems, a significant narrowing of the balance function in \(\Delta \eta \) with increasing multiplicity is observed.

The panels (b), (d), and (f) of Fig. 4 show the relative decrease of \(\sigma _{\Delta \eta }\), expressed by the ratio of \(\sigma _{\Delta \eta }\) for each multiplicity class over the value in the lowest multiplicity class, i.e. 70–80 % for all collision systems. The narrowing of the balance function with increasing multiplicity is most prominent in Pb–Pb collisions where the relative decrease between the largest and lowest multiplicity class is \(21.2 \pm 2.4 ({\mathrm {stat.}}) \pm 2.4 (\mathrm {syst.})~\%\). A significant relative decrease is also observed for the other two systems with values of \(6.7 \pm 0.2 ({\mathrm {stat.}}) \pm 0.4 ({\mathrm {syst.}})~\%\) and \(7.0 \pm 0.3 ({\mathrm {stat.}}) \pm 1.4 ({\mathrm {syst.}})~\%\) in p–Pb and pp collisions, respectively. Note though that the multiplicities in these three systems are significantly different (see e.g. Table 1)

In Fig. 4a the width in \(\Delta \eta \) for Pb–Pb collisions is compared with the results from HIJING and AMPT. Neither model describes the experimentally observed narrowing of the balance function with increasing multiplicity. This is also reflected in Fig. 4b where the relative decrease for both models is around 4 \(\%\).

Figure 4c shows the comparison of \(\sigma _{\Delta \eta }\)  in p–Pb collisions with model calculations. It is seen that DPMJET results in broader balance function distributions compared to AMPT. In addition, both models expect narrower balance function distributions compared to experimental measurements for low multiplicity classes (starting from 60 \(\%\) for DPMJET and 40 \(\%\) for AMPT). However, with increasing multiplicity (i.e. below 60 \(\%\) for DPMJET and 30 \(\%\) for AMPT) the balance function distributions are significantly narrower in the experiment compared to either of the models. Similar to the Pb–Pb case, neither of the models is able to reproduce the significant decrease of the width with increasing multiplicity observed in data. This is also reflected in Fig. 4d, where the relative decrease of the width between the highest and lowest multiplicity class for DPMJET and AMPT is marginal and not larger than 2 \(\%\).

The experimental results for pp collisions are compared with model predictions in Fig. 4e. PYTHIA8 without color reconnection, represented by the solid line, fails to describe the significant narrowing of the balance function with increasing multiplicity. The values of \(\sigma _{\Delta \eta }\) for this calculation are comparable within uncertainties to the ones obtained for the lowest multiplicity class in data. On the other hand, the inclusion of color reconnection, see the dashed line in Fig. 4e, results in a qualitatively similar narrowing as the one observed in the measurements. The absolute value of \(\sigma _{\Delta \eta }\) is lower than the experimental results for almost all multiplicity classes. Quantum statistics correlations are not included in the simulation, which might be the reason for this difference. Figure 4f that presents the relative decrease of \(\sigma _{\Delta \eta }\)  quantifies the previous observations. It is seen that PYTHIA8 without color reconnection shows a rather weak (i.e. around 2 \(\%\)) narrowing of the balance function with increasing multiplicity. This narrowing may result from the increased resonance yield for high- compared to low-multiplicity pp events [54]. The version of PYTHIA8 with the inclusion of color reconnection expects a relative reduction of around 7 \(\%\), in quantitative agreement with the measurement.

Fig. 5
figure 5

The multiplicity-class dependence of \(\sigma _{\Delta \varphi }\) in Pb–Pb, p–Pb, and pp collisions at \(\sqrt{s_{\mathrm {NN}}} = 2.76\), 5.02, and 7 TeV compared with results from various event generators in panels a, c, and e. Panels b, d, and f show the relative decrease of \(\sigma _{\Delta \varphi }\) calculated with respect to \(\sigma _{\Delta \varphi }^{70{-}80\,\%}\) as a function of the multiplicity class. The transverse momentum values for both the trigger and the associated particles satisfy the condition \(0.2 < p_{\mathrm {T,assoc}} < p_{\mathrm {T,trig}} < 2.0\) GeV/\(c\)

Figure 5 presents the multiplicity dependence of \(\sigma _{\Delta \varphi }\) in Pb–Pb, p–Pb, and pp collisions in panels (a), (c), and (e), respectively. All three systems exhibit a significant multiplicity-dependent narrowing of the balance function in \(\Delta \varphi \). Panels (b), (d), and (f) quantify this narrowing by presenting the decrease of the width in \(\Delta \varphi \) for each multiplicity class relative to the lowest multiplicity class. The data exhibit a narrowing of \(26.5 \pm 1.0 ({\mathrm {stat.}}) \pm 1.4 ({\mathrm {syst.}})~\%\), \(10.2 \pm 0.3 ({\mathrm {stat.}}) \pm 0.2 ({\mathrm {syst.}})~\%\), and \(10.8 \pm 0.4 ({\mathrm {stat.}}) \pm 1.4 ({\mathrm {syst.}})~\%\) in Pb–Pb, p–Pb, and pp collisions.

The multiplicity dependence of the width in \(\Delta \varphi \) in Pb–Pb collisions is compared with expectations from HIJING and AMPT in Fig. 5a. HIJING fails to describe the experimental measurements while AMPT expects a significant decrease of \(\sigma _{\Delta \varphi }\) with increasing multiplicity. The relative decrease in AMPT is about \(18~\%\), see Fig. 5b, and can be attributed to a rather strong multiplicity-dependent radial flow in the model that acts over the balancing pairs, retaining their initial correlations in \(\Delta \varphi \).

Fig. 6
figure 6

The balance function for charged particles with \(2.0 < p_{\mathrm {T,assoc}} < 3.0 < p_{\mathrm {T,trig}} < 4.0\) GeV/\(c\) as a function of \(\Delta \eta \) (upper row) and \(\Delta \varphi \) (lower row) in different multiplicity classes of Pb–Pb, in panels a and d, p–Pb, in panels b and e, and pp collisions, in panels c and f, at \(\sqrt{s_{\mathrm {NN}}} = 2.76\), 5.02, and 7 TeV, respectively

The measurements in p–Pb collisions are compared with the results from DPMJET and AMPT in Fig. 5c. Neither DPMJET, which does not exhibit a significant dependence on the event multiplicity, nor AMPT, which exhibits a relative decrease of around \(4~\%\), can quantitatively describe the experimental findings, as demonstrated in Fig. 5d.

Finally, the values of \(\sigma _{\Delta \varphi }\) in pp collisions are compared in Fig. 5e with the two variants of PYTHIA8 calculations described before. Similarly to the picture that emerged from the comparison of \(\sigma _{\Delta \eta }\), the variant of PYTHIA8 calculation without the inclusion of color reconnection does not describe the strong multiplicity dependence reported in pp collision data. However, the calculation with color reconnection exhibits a qualitatively similar decrease of \(\sigma _{\Delta \varphi }\) with increasing multiplicity. The relative decrease for this model is around 10 \(\%\), in quantitative agreement with the experimental results, as indicated in Fig. 5f.

The comparison between the data and the corresponding expectations from models like PYTHIA, illustrates the potentially significant role of color reconnection on charge-dependent correlations for small systems such as pp collisions. The effect of color reconnection in PYTHIA8 is strongly connected to MPIs, whose number increases with increasing multiplicity. In high-multiplicity pp events, MPIs lead to many color strings that will overlap in physical space. Within PYTHIA8 approach, these strings are given a probability to be reconnected and hence hadronize not independently, but rather in a process that resembles collective final-state effects. This process results in a transverse boost of the fragments that leads to the development of final-state correlations between charged particles in a similar way as a collective radial boost does.

6.2 Balance function at high transverse momentum

In order to study if the narrowing of the balance function is restricted to the bulk particle production at low \(p_{{\mathrm {T}}}\) or is also connected to hard processes, the balance function was also measured in all collision systems for higher values of transverse momentum for both trigger and associated particles. Figure 6 presents the projections of the two- dimensional balance functions in \(\Delta \eta \) in panels (a), (b), (c), and \(\Delta \varphi \) in panels (d), (e), (f) for \(2.0 < p_{\mathrm {T,assoc}} < 3.0 < p_{\mathrm {T,trig}} < 4.0\) GeV/\(c\) in Pb–Pb, p–Pb, and pp collisions, respectively. The analysis of Pb–Pb and p–Pb collisions was also extended to higher transverse momenta, \(3.0 < p_{\mathrm {T,assoc}} < 8.0 < p_{\mathrm {T,trig}} < 15.0\) GeV/\(c\) , shown in Fig. 7.

Fig. 7
figure 7

The balance function for charged particles with \(3.0 < p_{\mathrm {T,assoc}} < 8.0 < p_{\mathrm {T,trig}} < 15.0\) GeV/\(c\) as a function of \(\Delta \eta \) (upper row) and \(\Delta \varphi \) (lower row) in different multiplicity classes of Pb–Pb in panels a and c and p–Pb collisions in panels b and d at \(\sqrt{s_{\mathrm {NN}}} = 2.76\) and 5.02 TeV, respectively

The charge-dependent correlations exhibit little if any multiplicity dependence, in contrast to the findings from the lower transverse momentum region. In addition, the distributions in the intermediate and high-\(p_{{\mathrm {T}}}\) range are significantly narrower than the corresponding distributions at lower values of \(p_{{\mathrm {T}}}\) for each multiplicity class.

Fig. 8
figure 8

The multiplicity-class dependence of \(\sigma _{\Delta \eta }\) (a) and \(\sigma _{\Delta \varphi }\) (b) in Pb–Pb collisions at \(\sqrt{s_{\mathrm {NN}}} = 2.76\) TeV. The different markers represent the low (i.e. \(0.2 < p_{\mathrm {T,assoc}} < p_{\mathrm {T,trig}} < 2.0\) GeV/\(c\) with red circles), intermediate (i.e. \(2.0 < p_{{\mathrm {T,assoc}}} < 3.0 < p_{{\mathrm {T,trig}}} < 4.0\) GeV/c with blue squares), and high (i.e. \(3.0 < p_{\mathrm {T,assoc}} < 8.0 < p_{\mathrm {T,trig}} < 15.0\) GeV/\(c\) with green triangles) transverse momentum regions used in this analysis

Fig. 9
figure 9

The multiplicity-class dependence of \(\sigma _{\Delta \eta }\) (a) and \(\sigma _{\Delta \varphi }\) (b) in p–Pb collisions at \(\sqrt{s_{\mathrm {NN}}} = 5.02\) TeV. The different markers represent the low (i.e. \(0.2 < p_{\mathrm {T,assoc}} < p_{\mathrm {T,trig}} < 2.0\) GeV/\(c\) with red circles), intermediate (i.e. \(2.0 < p_{{\mathrm {T,assoc}}} < 3.0 < p_{{\mathrm {T,trig}}} < 4.0\) GeV/c with blue squares), and high (i.e. \(3.0 < p_{\mathrm {T,assoc}} < 8.0 < p_{\mathrm {T,trig}} < 15.0\) GeV/\(c\) with green triangles) transverse momentum regions used in this analysis

Fig. 10
figure 10

The multiplicity-class dependence of the width of the balance function in \(\Delta \eta \) (a) and in \(\Delta \varphi \) (b) in pp collisions at \(\sqrt{s} = 7\) TeV. The results correspond to the intermediate transverse momentum region (i.e. \(2.0 < p_{\mathrm {T,assoc}} < 3.0 < p_{\mathrm {T,trig}} < 4.0\) GeV/c). The data points are compared with two versions of PYTHIA8 calculations

The widths of the balance function, \(\sigma _{\Delta \eta }\) and \(\sigma _{\Delta \varphi }\) for the different transverse momentum regions, are presented in Figs. 8 and 9 as a function of the multiplicity class, for Pb–Pb and p–Pb collisions, respectively. The observed narrowing of the balance function with increasing multiplicity is restricted to the lower transverse momentum region, i.e. where the bulk of particles are produced. For higher transverse momenta, the multiplicity class dependence is significantly reduced, or even vanishes. In addition, the values of \(\sigma _{\Delta \eta }\) and \(\sigma _{\Delta \varphi }\) decrease with increasing \(p_{{\mathrm {T}}}\) for a given multiplicity class. This decrease can be attributed to the transition to a region where initial hard-scattering processes and parton fragmentation become the dominant particle production mechanism. The emerging hadrons are thus correlated within a cone whose angular size decreases with increasing \(p_{{\mathrm {T}}}\).

For pp collisions, the widths of the balance function \(\sigma _{\Delta \eta }\) and \(\sigma _{\Delta \varphi }\) are compared with results from PYTHIA in Fig. 10. The tune of PYTHIA8 without the inclusion of color reconnection is found to describe the data at a quantitative level, for both \(\sigma _{\Delta \eta }\) and \(\sigma _{\Delta \varphi }\). On the other hand, PYTHIA8 with the inclusion of color reconnection shows a broadening of the distributions with increasing multiplicity in both \(\Delta \eta \) and \(\Delta \varphi \), which is not supported by the data.

6.3 Comparison between the three systems

Fig. 11
figure 11

The width of the balance function in \(\Delta \eta \) (a) and in \(\Delta \varphi \) (b) for the three systems analyzed (pp, p–Pb, and Pb–Pb), as a function of the charged-particle multiplicity, estimated with the V0A for \(|\eta | < 0.8\) and \(p_\mathrm{{T}} > 0.2\) GeV/c. The low-, intermediate-, and high-\(p_{{\mathrm {T}}}\) intervals correspond to \(0.2 < p_{\mathrm {T,assoc}} < p_{\mathrm {T,trig}} < 2.0\) GeV/\(c\), \(2.0 < p_{\mathrm {T,assoc}} < 3.0 < p_{\mathrm {T,trig}} < 4.0\) GeV/\(c\), and \(3.0 < p_{\mathrm {T,assoc}} < 8.0 < p_{\mathrm {T,trig}} < 15.0\) GeV/\(c\), respectively

A comparison of the widths of the balance function in pp, p–Pb, and Pb–Pb as a function of particle multiplicity can provide direct information about differences and similarities between these systems in e.g. particle production mechanisms. It is important to note though, that this is performed for different center-of-mass energies which could complicate the comparison.

Figure 11 presents the charged-particle multiplicity dependence of the width of the balance function in \(\Delta \eta \) (a) and \(\Delta \varphi \) (b) for all three systems. The results for the low-, intermediate- and high-\(p_{{\mathrm {T}}}\) intervals are shown in the same plot. Multiplicity is defined as the number of charged particles reconstructed in \(|\eta | < 0.8\) and \(p_{{\mathrm {T}}}\) \(> 0.2\) GeV/c, as described in Sect. 3. It is seen that between the pp and the p–Pb systems, and for overlapping multiplicities in the low-\(p_{{\mathrm {T}}}\) region, the width in both \(\Delta \eta \) and \(\Delta \varphi \) has similar values. This could indicate that the charge-dependent correlations have similar origin in these two systems. On the other hand, the comparison of the results between p–Pb and Pb–Pb at the overlapping multiplicities indicate differences for both \(\sigma _{\Delta \eta }\) and (to a smaller extent in) \(\sigma _{\Delta \varphi }\). The origin of the charge-dependent correlations probed with the balance function in Pb–Pb collisions is believed to be related to radial flow and/or to a delayed hadronization scenario. The differences observed in the results of the Pb–Pb system compared with the ones in pp and p–Pb collisions at similar multiplicities could be explained by a different mechanism that drives the charge-dependent correlations in smaller systems.

With increasing values of transverse momentum, the balance functions become narrower and exhibit no significant multiplicity dependence for all systems, as discussed previously. The origin of these correlations at these transverse momentum ranges could be connected to initial hard parton scattering and subsequent fragmentation. The agreement of the values of both \(\sigma _{\Delta \eta }\) and \(\sigma _{\Delta \varphi }\) for all multiplicities over all three systems clearly indicates that the dynamics responsible for the high-\(p_{{\mathrm {T}}}\) charge-dependent correlations do not change significantly between pp, p–Pb, and Pb–Pb.

Fig. 12
figure 12

The multiplicity-class dependence of the width of the balance function in \(\Delta \eta \) (a) and in \(\Delta \varphi \) (b) for the three systems analyzed (pp, p–Pb, and Pb–Pb) relative to the 70–80 \(\%\) multiplicity class

The narrowing of the balance function in both \(\Delta \eta \) and \(\Delta \varphi \) is a distinct characteristic of the low transverse momentum region. Figure 12 visually illustrates the relative decrease of the different systems in this region. It is interesting that the relative decrease in this representation is similar between the two small systems (around 7 and 10.5 \(\%\) in \(\Delta \eta \) and \(\Delta \varphi \), respectively) while different for Pb–Pb. The results from the analysis of Pb–Pb collisions illustrate a significantly larger relative decrease of 21.2 \(\%\) in \(\Delta \eta \) (26.5 \(\%\) for \(\Delta \varphi \)). A direct comparison of the width at the same multiplicity class can not be done because, for the same class, the physics conditions are quite different for pp, p–Pb, and Pb–Pb collisions. However, the comparison of the relative decrease and the agreement of the results in both \(\Delta \eta \) and \(\Delta \varphi \) between the two small systems could indicate that they share a similar mechanism which is responsible for the decrease of the width with increasing multiplicity.

7 Summary

This article reports the first measurements of the balance function for charged particles in pp, p–Pb, and Pb–Pb collisions at the LHC measured with the ALICE detector. For all three systems, the balance function in both relative pseudorapidity (\(\Delta \eta \)) and relative azimuthal angle (\(\Delta \varphi \)) was studied for up to 9 multiplicity classes, and different trigger and associated particle transverse momentum. The widths of the balance functions in \(\Delta \eta \) and \(\Delta \varphi \) were found to decrease with increasing multiplicity for all systems only in the low-\(p_{{\mathrm {T}}}\) region (for \(p_{{\mathrm {T}}}\) \(< 2.0\) GeV/c). For higher values of \(p_{{\mathrm {T}}}\), the multiplicity-class dependence is significantly reduced, if not vanished, and the correlations of balancing partners are stronger with respect to the low-\(p_{{\mathrm {T}}}\) region. Models incorporating collective effects, such as AMPT, reproduce the narrowing of the experimental points qualitatively in \(\Delta \varphi \), but fail to reproduce the dependence in \(\Delta \eta \). On the other hand, models based on independent pp collisions such as DPMJET and HIJING do not show any narrowing in p–Pb and Pb–Pb. The comparison of the results in pp collisions with different PYTHIA8 tunes indicates the importance of MPIs and of the color reconnection mechanism, whose inclusion within this model allows for a qualitative description of the experimentally measured narrowing with increasing multiplicity at low values of transverse momentum.