Precision luminosity measurements at LHCb

Measuring cross-sections at the LHC requires the luminosity to be determined accurately at each centre-of-mass energy $\sqrt{s}$. In this paper results are reported from the luminosity calibrations carried out at the LHC interaction point 8 with the LHCb detector for $\sqrt{s}$ = 2.76, 7 and 8 TeV (proton-proton collisions) and for $\sqrt{s_{NN}}$ = 5 TeV (proton-lead collisions). Both the"van der Meer scan"and"beam-gas imaging"luminosity calibration methods were employed. It is observed that the beam density profile cannot always be described by a function that is factorizable in the two transverse coordinates. The introduction of a two-dimensional description of the beams improves significantly the consistency of the results. For proton-proton interactions at $\sqrt{s}$ = 8 TeV a relative precision of the luminosity calibration of 1.47% is obtained using van der Meer scans and 1.43% using beam-gas imaging, resulting in a combined precision of 1.12%. Applying the calibration to the full data set determines the luminosity with a precision of 1.16%. This represents the most precise luminosity measurement achieved so far at a bunched-beam hadron collider.


Introduction
The determination of the cross-section of a given subatomic process at high energy collidingbeam experiments is generally performed by the measurement of an interaction rate. To determine such a cross-section on an absolute scale, a measurement of the colliding-beam luminosity must be performed. The requirement for the accuracy on the value of the cross-section is usually driven by the precision of theoretical predictions for the process. At the LHCb experiment [1] the cross-section measurements for the production of vector bosons (Z and W ) [2,3] and the exclusive two-photon production of muon pairs [4] motivate an accuracy of order 1-2% for the luminosity calibration.
The instantaneous luminosity L is defined by the relation between the reaction rate R and the process cross-section σ R = L σ .
The instantaneous luminosity for a colliding bunch pair can be written as [5][6][7] where N 1 and N 2 are the populations of the colliding bunches of beam 1 and beam 2, ν rev is the revolution frequency and the beam overlap integral Ω embodies the passage of the two bunches with spatial particle density distributions ρ 1 (x, y, z, t) and ρ 2 (x, y, z, t) accross each other. In the limit of ultra-relativistic particles (velocity close to the speed of light, v ≈ c), crossing at small angle, the beam overlap integral is given by Ω = 2c ρ 1 (x, y, z, t) ρ 2 (x, y, z, t) dx dy dz dt .
Methods for absolute luminosity determination can be classified as being either direct or indirect. Indirect methods are e.g. the use of the optical theorem to make a simultaneous measurement of the elastic and total cross-sections [8,9], or the comparison with a process for which the absolute cross-section is known, either from theory or by a previous direct measurement. Direct methods derive the luminosity from the measurements of the colliding beam parameters. The analysis described in this paper relies on two direct methods to determine the absolute luminosity calibration: the "van der Meer scan" method (VDM) [10][11][12][13][14] and the "beam-gas imaging" method (BGI) [10,15], the latter making use of unique capabilities of the LHCb experiment. The VDM method exploits the ability to move the beams in both transverse coordinates with high precision and thus to scan the overlap integral of the colliding beams at different relative beam positions while measuring a relative rate. This method, which was first applied at the CERN ISR [11], is also being used by the other LHC experiments [16][17][18]. The BGI method is based on reconstructing vertices of interactions between beam particles and gas nuclei in the beam vacuum to measure the angles, positions and shapes of the individual beams without displacing them. The shapes obtained with these data are constrained by the distribution of vertices measured with beam-beam interactions. In both methods, data taken with the LHCb detector located at interaction point (IP) 8 are used in conjunction with data from the LHC beam instrumentation.
At the LHC, from 2009 to 2013, several luminosity calibration measurements were performed with a gradually improving precision. Different nucleon-nucleon centre-of-mass energies √ s and different beam species were used: protons on protons (pp), lead on lead (Pb-Pb) and protons on lead (pPb or Pbp, where the first/second beam species applies to beam 1/beam 2 in the standard LHC definition [19], see Fig. 1). First LHC luminosity calibrations were obtained by LHCb using pp collision data collected at the end of 2009 at √ s = 900 GeV [20] and in 2010 at √ s = 7 TeV [10,21,22] with an accuracy that was limited by the systematic uncertainties associated with the normalization of the colliding bunch populations [23,24]. Recent detailed studies of the LHC beam current transformers (BCTs) significantly reduced these uncertainties [25][26][27], thus facilitating an improvement of the final precision of the luminosity calibration. In this paper results are reported from luminosity calibration experiments carried out at the LHC IP8 with the LHCb detector from 2011 to 2013, for √ s = 2.76, 7 and 8 TeV in pp collisions and for √ s NN = 5 TeV in pPb and Pbp collisions. In addition to performing luminosity calibration measurements, LHCb provided related beam-gas interaction measurements as a service to the other LHC experiments. This included the measurement of the total charge outside the nominally filled slots ("ghost charge", see Sec. 3) and of the single beam size as a function of time during the VDM scans of these other experiments.
The precision of the luminosity calibration in the LHCb experiment is now limited by the systematic uncertainties of the beam overlap determination. These systematic uncertainties are different, to a large extent, for the VDM and BGI methods. Therefore, the comparison provides an important cross check of the results. The calibration measurements obtained with the VDM and BGI methods are found to be consistent and are averaged for the final result.
Since the absolute calibration can only be performed during specific running periods, a relative normalization method is needed to transport the results of the absolute luminosity calibration to the complete data-taking period. To this end, several observables are used, each one corresponding to an effective visible cross-section σ vis . The corresponding crosssection is calibrated for each variable using the measurements of the absolute luminosity during specific data-taking periods. The integrated luminosity for an arbitrary period of data taking is then obtained from the accumulated counts of a calibrated visible cross-section.
In the present paper we first describe briefly the LHCb experimental setup and datataking conditions in Sec. 2, emphasizing the aspects relevant to the analysis presented here. Section 3 is devoted to the normalization of the bunch population, while the methods used for the relative normalization technique are given in Sec. 4. In Sec. 5 we introduce the luminosity formalism for colliding beams. The determination of the luminosity with the BGI method is detailed in Sec. 6 and with the VDM scan method in Sec. 7. The combination of the results and conclusions are given in Sec. 8. (muon stations 1, 2, 3, 4, and 5) (drawing from Ref. [28]).

Experimental setup and data-taking conditions
The LHCb detector (Fig. 1) is a single-arm forward spectrometer with a polar angular coverage of approximately 15 to 300 mrad in the horizontal (bending) plane, and 15 to 250 mrad in the vertical plane. It is designed for the study of particles containing b or c quarks and is described in detail elsewhere [1].
The apparatus contains tracking detectors, ring-imaging Cherenkov detectors, calorimeters, and a muon identification system. The tracking system comprises the vertex locator (VELO) surrounding the beam interaction region, a tracking station upstream of the dipole magnet and three tracking stations located downstream of the magnet. Particles traversing the spectrometer experience a bending-field integral of around 4 Tm.
The VELO plays an essential role in the application of the VDM and BGI luminosity calibration methods at LHCb. It consists of two horizontally retractable halves, each action point, cross at least three stations. This requirement ensures that the track will be properly reconstructed. The resulting arrangement of the stations which respects the requirements, while being close to the beam for precision, and introducing a minimum amount of material to traversing particles, is shown in figure 3.7. An additional two VELO stations, located more upstream, are called the pile-up system. This identifies bunch crossings with multiple interactions and through the first-level hardware trigger vetoes such events, as detailed in subsection 3.5.1. The VELO uses semi-circular silicon sensors in a 10 4 mbar vacuum, separated from the machine vacuum by a corrugated 300 µm thick Aluminium foil. A corrugated design minimises the interaction length encountered by particles, allows the sensors to overlap and o↵ers greater mechanical strength compared to a flat foil. The foil protects the machine vacuum from the lower quality vacuum inside the VELO and shields the sensors from the RF currents induced by the beams. On the sensor side, the foil is coated to electrically insulate it from the sensors. Both the sensors and foil can be moved to and from the beam line within a range from 5 mm Figure 2: Sketch of the VELO sensor positions. The luminous region is schematically depicted with a filled ellipse. Its longitudinal extent, RMS σ = 53 mm, is indicative. Sensors measuring the R (φ) coordinates are shown as blue (red) lines. The LHC beam of ring 1 (2) enters from the left (right) on this sketch. The coordinate system is defined in Sec. 5 (drawing from Ref. [28]).
having 21 modules of radial and azimuthal silicon-strip sensors in a half-circle shape (Fig. 2). Two additional stations (Pile-Up System, PU) upstream of the VELO tracking stations are mainly used in the hardware trigger. The VELO has a large acceptance for beam-beam interactions owing to its many layers of silicon sensors and their close proximity to the beam line. During nominal operation, the distance between the closest sensor strip and the beams is only 8.2 mm. During injection and beam adjustments, the two VELO halves are kept apart in a retracted position 30 mm away from the beams. They are brought to their nominal position close to the beams during stable beam periods only. More details about the VELO can be found in Ref. [29]. The LHCb trigger system [30] consists of two separate levels: a hardware trigger, which is implemented in custom electronics, and a software trigger, executed on a farm of commercial processors. The hardware trigger is designed to have an accept rate of 1 MHz and uses information from the PU sensors of the VELO, the calorimeters and the muon system. These detectors send information to the hardware decision unit, where selection algorithms are run synchronously with the 40 MHz LHC bunch crossing. For every nominal bunch-crossing slot (i.e. each 25 ns) the hardware decision unit sends its information to the LHCb readout supervisor, which distributes the synchronous hardware trigger decision to all front-end electronics. For every positive hardware decision the full event information of all sub-detectors is sent to the processor farm and is made available to the software trigger algorithms.
For luminosity calibration and monitoring, a trigger strategy is adopted to select beam-beam inelastic interactions and interactions of the beams with the residual gas in the vacuum chamber. Events are collected for the four bunch-crossing types: two colliding bunches (bb), one beam 1 bunch with no beam 2 bunch (be), one beam 2 bunch with no beam 1 bunch (eb) and nominally empty bunch slots (ee). Here "b" stands for "bunch" and "e" stands for "empty". The first two categories of crossings produce particles in the forward direction and are triggered using calorimeter information. An additional PU veto is applied for be crossings. Crossings of the type eb produce particles in the backward direction, are triggered by demanding a minimal hit multiplicity in the PU, and are vetoed by calorimeter activity. The trigger for ee crossings is defined as the logical OR of the conditions used for the be and eb crossings in order to be sensitive to background from both beams. In addition to these specific triggers, a decision based on a hardware trigger sensitive to any activity in the PU and calorimeter is available. The latter hardware trigger configuration is used for most measurements described in this paper. Events are then further selected by the software trigger based on track and vertex reconstruction using VELO hits. During VDM scans specialized trigger configurations are defined that optimize the data taking for these measurements (see Sec. 7).
The reconstruction of interaction vertices (also called "primary vertices", PVs) is performed using standard LHCb algorithms [31]. The initial estimate of the PV position is based on an iterative clustering of tracks. For each track the distance of closest approach (DOCA) with respect to all other tracks is calculated and tracks are clustered into a PV candidate if their DOCA is less than 1 mm. An initial position of the PV is obtained from the weighted average of the points of closest approach between all track pairs, after removing outliers. The final PV coordinates are determined by iteratively improving the position determination with an adaptive, weighted, least-squares fit. Participating tracks are assigned weights depending on their impact parameter with respect to the PV. The procedure is repeated for all possible track clusters, excluding tracks from previously reconstructed PVs, retaining only those with at least five tracks. For the analysis described here only PVs with a larger number of tracks are used since they provide better position resolution. For the study of beam-gas interactions only PVs with at least ten tracks are used and at least 25 tracks are required for the study of beam-beam interactions. For specific studies different criteria are applied as described below.
The full list of luminosity calibrations discussed in this paper is summarized in Table 1. The table is divided into five sections following the different nucleon-nucleon centre-of-mass energies and beam species involved. A first measurement with intentionally enlarged beta functions at the IP (β * = 10 m) was performed in October 2011 with pp collisions at √ s = 7 TeV. Several fills in 2012 were dedicated to luminosity calibration for pp collisions at √ s = 8 TeV, although only the measurements in July and November were performed with large β * . The April measurements were performed in non-optimal conditions, with focused beams (β * = 3 m) and with a tilted crossing plane (a non-zero vertical half crossing angle φ y ), and are therefore primarily used for the VDM calibration method and to cross-check the effects on the BGI method of the finite vertex resolution. Calibrations for pPb and Pbp were conducted in January 2013 at √ s NN = 5 TeV with VDM scans only.
Further pp calibrations were performed at √ s = 2.76 TeV in February 2013, exclusively using the BGI method. The number of bunches per beam is also given in the table. No active gas injection was used to enhance the beam-gas rates and the end of 2011, though a first rate increase was obtained in October 2011 by degrading the beam vacuum by switching off the VELO ion pumps. Thus, three configurations of the VELO vacuum In most calibration measurements the number of bunches per beam was the same for beam 1 and beam 2. For the pPb and Pbp fills where this was not the case two numbers are given, the first for the number of beam 1 bunches, the second for the number of beam 2 bunches. The number of colliding bunches at LHCb is indicated in parentheses (fifth column). Half crossing angles φ x and φ y , and β * are given as nominal values. The VELO vacuum state during BGI measurements is indicated in the column "Gas injection". A state "off" means that gas injection was turned off and the VELO ion pumps were turned off, which resulted in a residual vacuum pressure about a factor four higher than nominal. A state "on" indicates that neon gas was being injected into the beam vacuum. During VDM measurements the state was always "off".

Period
Fill where the VELO ion pumps were switched off, and one where, in addition to running with pumps off, neon gas was injected into the VELO vacuum chamber (see Sec. 6). All pp BGI calibration measurements of 2012 and 2013 took advantage of gas injection. During VDM calibration scans, gas injection was always off. In all pp calibration runs discussed here the initial bunch populations ranged between 0.6 and 1.1 × 10 11 particles. For the pPb and Pbp runs they varied between 1 and 2 × 10 10 elementary charges (for both beam species). Calibration experiments with the VDM method included a variety of beam displacement sequences. The details of these individual experiments are given in the section devoted to the VDM analysis (Sec. 7). In fills 3503, 3537 and 3540, no luminosity calibration was performed at IP8, though the LHCb experiment provided ghost charge and beam size measurements for the benefit of the luminosity calibrations conducted in other LHC experiments.

Bunch current normalization
Various detector systems are used to determine with high precision the population of particles in each colliding and non-colliding bunch in the LHC. The longitudinal structure of the LHC beams is shaped by the 400 MHz radio frequency (RF) system. Both LHC rings are filled with bunches at locations ("RF buckets") defined by the RF system and are organized in "slots", which contain each ten consecutive buckets. Ideally, only one of these buckets is filled with a bunch, called the "main bunch", the other nine are nominally empty. Only a subset of the slots are filled in a given filling scheme. In each filled slot, the main bunch occupies the same bucket number. In practice, a small fraction (typically < 10 −3 ) of the charge in a slot occupies nominally empty buckets and are called "satellite" bunches. Additionally, also the nominally empty slots may contain charges. The total charge outside the filled slots is called "ghost charge".

Bunch population measurement
To measure the population in the main bunches, specific instruments are used to determine the overall circulating charge, the relative charge in the filled bunches, the fraction of the charge in the satellite buckets of the filled slots and the fraction of ghost charge.
Four independent direct-current current-transformers (DCCTs), two per ring, are used to measure the total beam current circulating in each LHC ring. The DCCT is designed to be insensitive to the time structure of the beam [32]. Two fast bunch current transformers (FBCTs), one per ring, provide a relative measure of the individual charges on a slot-by-slot basis [33]. The FBCT is designed to produce a signal proportional to the charge in each 25 ns LHC bunch slot. The captured particles of an LHC bunch are contained within an RF bucket of 1-1.5 ns length at ±2 standard deviations [34]. Since 2012, one longitudinal density monitor (LDM) [35,36] per LHC ring is available for detecting synchrotron radiation photons emitted by particles deflected in a magnetic field. The LDMs are used to obtain the longitudinal beam charge distribution with a time resolution of about 90 ps to resolve the charge distribution in individual RF buckets. Finally, the ghost charge fraction is obtained by counting beam-gas interactions with the LHCb detector in nominally empty (ee) compared to the rates in nominally filled (bb, be and eb) bunch crossings. Previous LHC luminosity calibration experiments showed that one of the dominant uncertainties arises from the normalization of the bunch population product N 1 N 2 . As a consequence, a detailed study of the normalization was carried out using data from the LHC beam current transformers (BCTs) and from the LHC experiments. A dedicated analysis procedure was defined and bunch population uncertainties were quantified for the 2010 LHC luminosity calibration measurements [23,24]. The precision was limited by the understanding of the BCT data at that stage. Since then, a number of additional tests were carried out that significantly improved the understanding of the bunch current measurements. Careful calibration measurements and systematic studies of the DCCTs improved the dominant uncertainty by an order of magnitude [25,37]. Uncertainties on the beam current product for the 2011-2013 measurements are well below 1% and are given in more detail below.
The accuracy of the relative bunch populations determined with the FBCT is crosschecked against results from other measurements, such as those obtained from the ATLAS BPTX button pick-up [38] and those derived from the LHCb beam-gas interaction rates [26]. The sum of the FBCT signals of all nominally filled bunch slots is normalized to the total number of particles measured by the DCCTs after subtraction of the ghost charge and satellite charges, defining N j,i as the bunch population of the nominally filled RF bucket of bunch slot i of beam j, and I DCCT,j as the current measured by the DCCTs and Z j e the charge of a beam particle (82 e for Pb beams). The sum runs over all nominally filled slots and the S FBCT,j,i are the signals measured by the FBCT of ring j. The ghost charge fraction is denoted f ghost,j and the fraction of the charge in satellite bunches f sat,j,i for beam j and slot i. Ghost charge fractions for the 2011-2013 LHC luminosity calibration fills range up to about 2.5%. As mentioned above, these measurements are performed with the LHCb detector. The results and methods are described in detail in Sec. 3.2.
Satellite charges have been observed in various ways with the LHC detectors by detecting longitudinally displaced collisions (see for example Ref. [23]). The total satellite population fraction (f sat,j,i ) in a bunch slot is usually less than a percent compared to the associated main bunch population. Nevertheless, it needs to be quantified to obtain a precise measurement of the bunch population that actually contributes to the luminosity.

Ghost charge
The determination of the ghost charge from the beam-gas interaction rate measurements was pioneered in a previous LHCb luminosity calibration [39]. The results presented here benefit from the larger number of beam-gas events obtained with neon gas injection in the beam vacuum chamber, which allows the uncertainty to be reduced and provides a more

Beam
Efficiency average j including double-counting excluding double-counting 1 1.05 ± 0.03 0.93 ± 0.02 2 0.90 ± 0.01 0.86 ± 0.01 detailed determination of the charge distribution over the LHC ring in a shorter time. Systematic uncertainties are further reduced by a better trigger efficiency calibration. The ghost charge measurement is based on the same data sample as used for the BGI analysis. The trigger requirements are described in Sec. 6. To ensure that each vertex is a result of a beam-gas interaction and is assigned to the correct beam, several selection criteria are applied [37] that are based on the track directions (all forward for beam 1, all backward for beam 2), on the transverse position (to exclude interactions with material in the vicinity of the beams), on the longitudinal position and on the vertex track multiplicity.
The LHCb data acquisition is synchronised with the LHC RF system with a granularity of 25 ns. The sampling phase of the detectors relative to the LHC clock is optimized to provide the highest efficiency for nominally filled RF buckets, but the trigger efficiency may vary across the 25 ns bunch slot. Since the ghost charge is distributed over all RF buckets inside the 25 ns slots, the trigger efficiency must be known for all possible phases. A first efficiency measurement was performed in 2010 [10,39], resulting in a ghost charge uncertainty of about 20% per beam. A new dedicated measurement was performed in 2012 with the aim of reducing this uncertainty by acquiring data for more clock phases and by using neon gas injection to increase the statistical accuracy. The efficiency is determined by measuring dead-time corrected beam-gas interaction rates from non-colliding bunches at different clock phases and comparing them with the standard phase (zero clock shift). The absolute rate is measured as function of clock shift in 2.5 ns steps. The beam intensity decay observed during the measurement is taken into account.
If a beam-gas interaction occurs near the bunch slot edges, that is, the originating charge is near the previous or next clock cycle, the resulting VELO sensor signals may be sufficiently long that they are also seen in the neighbouring clock cycle. Therefore, depending on where the charges are located within the 25 ns bunch slot, some vertices are counted twice and thus bias the ghost charge or trigger efficiency measurement. To take this double-counting effect into account, the efficiency is measured including all beam-gas events or, alternatively, excluding double-counted vertices. In addition, the efficiency is measured for different vertex track multiplicity thresholds (from 8 to 12 tracks) to account for the slightly different trigger conditions used for this measurement as compared to later BGI measurements. The results of the trigger efficiency calibration are shown in Fig. 3 and the values averaged over the 25 ns clock cycle are summarized in Table 2.
The increased rate of beam-gas interaction data acquired with neon gas injection   Figure 3: Relative beam-gas trigger efficiency as function of LHCb detector clock shift with respect to the LHC reference timing, (left) including or (right) excluding double-counted beam-gas interaction vertices. The efficiency is shown relative to the value at the nominal clock setting (i.e. zero shift). The shaded areas indicate the variation between the results for thresholds corresponding to 8 and 12 tracks. The data points, appearing in groups of three, indicate measurements applying the 8, 10 and 12 track thresholds. enables a measurement of the charge distribution over the ring circumference. In Fig. 4 the ghost charge per 25 ns slot is shown as function of slot number (BCID) using data from fill 2520 as an example. Ghost charges are observed around the nominally filled bunches and are mostly absent further than about 20 slots away from filled bunches.
Ghost charge fractions during LHC luminosity calibration fills are measured in fourminute time bins. For each time bin the ghost charge fraction is evaluated with both counting methods: including and excluding double-counted vertices and applying the corresponding average trigger efficiency of Table 2. If all charges are evenly spread within their bunch slot, each evaluation would provide a different result before efficiency correction, but the same result after efficiency correction. After efficiency correction the differences between the two evaluations are small. This observation is in agreement with the LDM measurements [40], which show that the ghost charge tends to be spread evenly over all RF buckets of a bunch slot. The LDM information on the charge distribution within the nominally empty bunch slots is not used in the results except for fill 3542 during which the trigger was not configured to perform this measurement. The average of the two efficiency-corrected evaluations is taken as final value for the ghost fraction, while their difference is taken as systematic uncertainty. The trigger efficiency uncertainty taken from Table 2 is added in quadrature with the systematic uncertainty. A summary of all ghost charge measurements performed for the special luminosity fills in 2011, 2012 and 2013 is provided in Table 3.
With good accuracy using an average value for a fill. In this case the RMS over the fill, given in Table 3, should be taken into account in the uncertainty. On the contrary, for the intermediate-energy fills, an increase in the ghost charge fraction over time warrants a time dependent correction to the total beam intensity. As an example, the difference in ghost charge evolution seen between high-and intermediate-energy fills is shown in Fig. 5

Total uncertainty
A summary of the bunch population product uncertainties is given in Table 4 for each luminosity calibration fill. The systematic uncertainties for the ghost charge corrections of the two beams described in the previous section are assumed to be fully correlated with each other, i.e. the final ghost charge uncertainty on the bunch population product is the linear sum of the ghost charge systematic uncertainty of each beam. The satellite fractions provided by the LDM [40] are measured at the beginning and at the end of the fill. Here, the average of these two measurements is used. The average satellite fractions for all colliding bunches and fills with β * = 10 m at √ s = 8 TeV are 0.25% and 0.18% for beam 1 and beam 2, respectively. The uncertainty on the satellite fraction correction is taken as the full difference between the fractions measured at the beginning and end of fill. Assuming the uncertainties are fully correlated between the two beams, the uncertainty on the population product due to the satellite fraction correction is taken as the linear sum of the average uncertainties per beam, and is given as the average per fill in Table 4. The beam population product normalization uncertainty is dominated by the DCCT measurement. All fills listed in Table 4 are subject to the same procedure to evaluate the beam population product uncertainty. For fills with β * = 10 m and √ s = 8 TeV, the average uncertainty on the bunch population product weighted with the number of measurements amounts to 0.22% at 68% confidence level.

Relative luminosity calibration
Absolute luminosity calibrations are performed during short periods of data-taking. To be able to determine the integrated luminosity for any data sample obtained during long periods, the interaction rate of standard processes is measured continuously. The effective cross-section corresponding to these standard processes is determined by counting the visible interaction rates during the specific periods when the absolute luminosity is

Interaction rate determination
The luminosity is proportional to the average number of visible proton-proton interactions per beam-beam crossing, µ vis . The subscript "vis" is used to indicate that this particular definition of interaction rate does not need to have a simple physics interpretation. Any interaction rate that can be measured under stable conditions can be used as such a relative luminosity monitor. The interaction rates are acquired and stored together with the physics data as "luminosity data". During further processing of the data the relevant luminosity information is kept in the same storage entity. Thus, it remains possible to select only part of the full data set for analysis and still keep the capability to determine the corresponding integrated luminosity. Triggers, which initiate the full readout of the LHCb detector, are created for a random choice of beam crossings at a fixed average frequency. These are called "luminosity Table 5: Definition of luminosity observables used in the analysis. The fiducial volume used here is a cylinder of radius < 4 mm around the z axis and bound by |z| < 300 mm. It is used to cut either on the point of closest approach of a track relative to the z axis or on the position of a vertex.

Name
Description Origin triggers". During normal physics data-taking, the overall rate is chosen to be 1000 Hz. Of this rate, 70% is assigned to slots where two bunches cross (bb), 15% to slots with only a beam-1 bunch (be), 10% to those with only a beam-2 bunch (eb) and the remaining 5% to slots that are empty (ee). The events taken for crossing types other than bb are used for background subtraction and beam monitoring. Interaction rates are measured by processing the random luminosity triggers and these rates are stored in a small number of "luminosity observables". The set of luminosity observables comprises the number of vertices and tracks reconstructed in the VELO, the number of muons reconstructed in the muon system, the number of hits in the PU and in the SPD in front of the calorimeters, and the transverse energy deposition in the calorimeters. The number of vertices in the VELO that fall within a limited region around the nominal interaction point and VELO tracks crossing this region are counted separately. Some of these observables are directly obtained from the hardware trigger decision unit, others are the result of partial event reconstruction in the software trigger or in the off-line software. Observables used in this analysis are summarized in Table 5.
The luminosity for a given data set can be determined by integrating the values of observables that are proportional to the instantaneous luminosity and by applying the corresponding absolute calibration constant. However, this procedure sets stringent requirements on the stability of the observable and on its linearity in the presence of multiple interactions. Alternatively, one may determine the relative luminosity from the fraction of "empty" or invisible events in bb crossings which we denote by P 0 . An invisible event is defined by applying an observable-specific threshold, below which it is considered that no pp interaction is seen in the corresponding bunch crossing. For a colliding bunch pair, the number of interactions per bunch crossing follows a Poisson distribution with mean value proportional to the luminosity, hence the luminosity is proportional to − ln P 0 . In the absence of backgrounds, the average number of visible pp interactions per crossing can be obtained from the fraction of empty bb crossings by µ vis = − ln P bb 0 . This "zero-count" method is both robust and easy to implement [41]. The choice of a low visibility threshold ensures a better behaviour under gain or efficiency variations of the observable than the straightforward linear summing method. In addition, any non-linearity encountered with multiple events does not play a role when counting empty slots.
Assuming equal particle populations in bb, be, and eb bunches and no particles in ee slots, backgrounds are subtracted using where P i 0 (i = bb, ee, be, eb) are the probabilities to find an empty event in a bunch-crossing slot for the four different bunch-crossing types. In Eq. (5) it is implicitly assumed that all bunches of the same type have the same properties. The consequences of this approximation will be discussed in Sec. 4.2. The P ee 0 contribution is added because it is also contained in the P be 0 and P eb 0 terms. The purpose of the background subtraction, Eq. (5), is to correct the count-rate in the bb crossings for the detector response, which is due to beam-gas interactions and detector noise. In principle, the noise background is measured during ee crossings. In the presence of parasitic beam protons in ee bunch positions (ghost charge), it is not correct to evaluate the noise from P ee 0 . In addition, the detector signals are not fully confined within one 25 ns bunch-crossing slot for some of the observables. The empty (ee) bunch-crossing slots immediately following a bb, be or eb crossing slot contain detector signals from interactions occurring in the preceding slot ("spill-over"). However, because the filling schemes used for the data-taking described here did not contain adjacent filled slots, the spill-over background is negligible in the bb, be and eb crossings. Since the detector noise for the selected observables is small (see Sec. 4.2) the term ln P ee 0 in Eq. (5) is neglected.
The results of the zero-count method based on the number of tracks and vertices reconstructed in the VELO are found to be the most stable. An empty event is defined to have < 2 tracks in the VELO. A VELO track is defined by at least three R clusters and three φ clusters on a straight line in the VELO detector. The number of tracks reconstructed in the VELO restricted to a fiducial region is chosen as the reference observable.

Systematic uncertainties
The zero-count method is valid if an event is considered empty when the value of the observable is exactly zero. However, if the observable is affected by noise such that its value is never zero, the threshold discriminating empty events has to be increased. This is the case for the ECalo and Calo observables, used as a cross-check, for which a positive threshold must be chosen. The introduced bias depends on the noise distribution, the  one-interaction spectrum and the average number of interactions per crossing. While the latter was kept approximately constant during the 2012 data-taking period, the Calo noise distribution was changing due to ageing of the hadron calorimeter. In the second half of 2012, the HCAL gain was adjusted more frequently, thus keeping the noise distribution more stable.
The noise distribution is measured in ee crossings. Histograms of the mean value of the noise are shown for the ECalo and the Track observable in Fig. 6. For the ECalo observable, two peaks are observed in the pedestal distribution. This is attributed to a change of operating conditions, which is not easily corrected for. Therefore, for cross-checks using the ECalo and Calo observables we discard the runs for which the ECalo pedestal mean is lower than 4.75. The remaining larger fraction of runs spans the full year and is subsequently used for assessment of systematic uncertainties. The Track observable has typically less than 2.5 tracks per 100 ee crossings, which induces a negligible bias. The systematic uncertainty due to noise is negligible. Equation 5 assumes that the proton populations in the be and eb crossings are the same as in the bb crossings. With a population spread of typically 10% and a beam-gas background fraction for the reference observable < 0.3% compared to the pp interactions (see Fig. 7) the effect of the spread is small and therefore neglected.
The measured µ vis values can be contaminated by other backgrounds than beam-gas interactions, e.g. collisions between satellite and main bunches, and interactions with material in the VELO. We reduce such effects by applying a fiducial volume cut to the Track observable. To assess the magnitude of potentially unaccounted background, a comparison is made between the µ vis values measured with (Track ) and without (TrackNR) the fiducial volume cut, see Fig. 8. The observed discrepancy is used as an estimate of the systematic uncertainty due to beam-beam background.   The longitudinal position of the luminous region depends on the transverse beam separation in the crossing plane. In 2012 the instantaneous luminosity at LHCb was kept approximately constant by varying the separation of the two beams at the crossing point [42] (so called "luminosity levelling"). Due to imperfections in the procedure, the beams were displaced in a direction not exactly orthogonal to the crossing plane. As a result, there was a significant variation in the longitudinal position of the luminous region ξ lz , as shown in Fig. 9. The reconstruction efficiency for tracks and vertices is not uniform along z at the scale of the observed variations. Therefore, a correction needs to be applied to the observed µ vis values that are measured using VELO observables. The Calo observable is not affected.
From simulation we determine I(0|z), the probability to obtain an empty event while having one interaction at z, I(0|z) ≡ P (empty event | one interaction at z) ,  distribution, the probabilityĪ(0|f ) to have an empty event while one interaction occurred isĪ (0|f ) = I(0|z) f (z) dz .
The absolute normalization of I(0|z) andĪ(0|f ) depends on the underlying interaction generator. However, the normalization does not affect the luminosity measurement if used consistently. To avoid scaling µ vis with factors largely different from unity, the correction is made with respect to a reference valueĪ(0|ref) that corresponds to a Gaussian probability density g(z) of the longitudinal vertex distribution centred at zero and having an RMS of σ lz = 50 mm,Ī (0|ref) = I(0|z) g(z) dz .
The observed values, µ raw vis , are proportional to 1 −Ī(0|f ). Thus, the corrected values are given by The z distribution of the vertices is well approximated with a Gaussian function. Examples of the correction factors for a Gaussian vertex distribution with σ lz = 50 mm are shown in  An efficiency correction is implemented as a factor, according to Eq. (9), evaluated as an average over about one-month running periods. To take into account a possible inaccuracy of the efficiency obtained by simulation, a comparison is made with an unaffected observable (Calo). The systematic uncertainty due to the residual dependence of the ratio µ Track /µ Calo on the longitudinal position of the luminous region ξ lz is estimated as follows. First, the data are divided in 5 mm bins in ξ lz and the median of the ratio in each bin is calculated. The difference of the median with respect to that at ξ lz = 0 is assumed to be due to imperfect correction of the Track observable. The relative difference is then averaged over the full data set taking into account the luminosity content. Finally, the resulting difference of 0.19% is taken as the systematic uncertainty due to the imperfect correction for the efficiency of the observable.
The numbers of protons, beam sizes and transverse offsets at the interaction point vary across bunches. Thus, the µ vis value varies across bb crossings. An estimate of the spread of µ vis values is the RMS divided by the mean across bunch crossings, as shown in Fig. 11. Due to the non-linearity of the logarithmic function, ideally one first needs to compute µ vis values for different bunch crossings and then take the average. However, for short time intervals the number of samples is insufficient to make an unbiased measurement per bunch crossing using the zero-count method, while µ vis may not be constant when the intervals are too long due to e.g. loss of bunch population and emittance growth.
During physics data-taking, the bandwidth reserved for luminosity triggers is limited. Therefore, a statistically significant measurement of the luminosity cannot be obtained for each bunch crossing individually by integrating over periods shorter than about 30 minutes. The bias and systematic uncertainty introduced by this limitation is evaluated with a simulation. To reflect the luminosity integration for physics data, the number of visible interactions is counted in short time intervals ignoring the spread of µ vis values across bunch crossings. A set of 30 consecutive short intervals is used to accumulate a sufficient number of events per BCID. Then, a correction for the spread is calculated and applied as described below. For the following discussion backgrounds are not considered since their effect on the correction is negligible. Let n ti and k ti denote the number of random triggers and the number of empty events, respectively, in bb crossing slot i and short time interval t. Where the index t is omitted, an implicit sum is assumed over the set T of consecutive short intervals. Similarly, in case the index i is omitted, a summation over all bb bunch crossings is assumed. A correction factor is calculated for every set T and is applied to each short period t ∈ T where the average in the numerator is taken over i. The corrected estimate of the number of visible interactions is where f trig is the probability that a bb crossing is randomly triggered. A simulation study is performed to compare the bias of the estimated number of interactions before and after the correction procedure. The rate of triggers and the number of bunch crossings is chosen to reflect the typical running conditions. The µ values across bunch crossings are sampled from a normal distribution. A luminosity half-life of two hours is assumed. The bias is calculated as function of mean µ value and relative RMS, and is shown in Fig. 12. To estimate the residual bias of the correction technique on the data, we perform a simulation for each long period T . First, the µ vis value is estimated for each short period t and each bunch crossing i with which has the desirable property that it coincides with the projection estimates (µ t and µ i ) when the true µ value does not change over time or across colliding bunches. Then, for each (t, i) pair, k ti is sampled from a binomial distribution with success probability e −µ i and number of trials equal to n ti . As for the actual data, Eqs. (10) and (11) are used to estimate the number of visible interactions. Finally, the bias is obtained from the difference of the estimated and the true number of visible interactions, averaged over 25 independent repetitions of the simulation. Histograms of the average values of κ T and the residual relative bias for each run are shown in Fig. 13. The relative integrated bias over the full data set is assigned as a systematic uncertainty (0.14%).
In addition, a cross-check is made using the Muon observable, which is less sensitive to the spread owing to its low µ vis values ranging from 0.07 to 0.15. 2 The ratio µ Track /µ Muon as a function of the relative RMS of µ vis across bunch crossings is fitted with an even quadratic polynomial. The 0.5% of runs with extreme values of the ratio are excluded from the fit. The maximum relative difference between the predicted value at the mean spread and at zero spread gives an estimate of the residual bias. Since it has an opposite sign with respect to the residual bias obtained from the simulation, the result is taken as an additional systematic uncertainty (0.09%). The stability of the reference observable is demonstrated in Fig. 14, which shows the ratio of the relative luminosities determined with the zero-count method using the Track and the Calo observables. These two observables use different sub-detectors and have different systematic uncertainties. The variation of the ratio unexplained by statistical fluctuations is assigned as a systematic uncertainty to the relative luminosity measured using the Track observable. A similar cross-check with the ratio of the relative luminosities using the Track and the Vertex observables shows negligible discrepancy.
The systematic uncertainties of the relative luminosity measurement are summarized in Table 6. By summing the effect of the different sources in quadrature, we conclude that the relative luminosity measurement introduces a systematic uncertainty of 0.31% for the pp run at 8 TeV. The quoted uncertainty applies when the full dataset is used; for specific choices of partial datasets a different value may apply. In the case of the 2013 running conditions (proton-lead and pp at 2.76 TeV), the corrections due to the movement of the luminous region and the bunch spread are small. Since the available data sample size is insufficient to reliably perform the corresponding cross-checks, the full amount of each correction is assigned as an uncertainty. The beam-beam background uncertainty is estimated to be up to 1% for the proton-lead data taking, owing to the very low µ values (0.01-0.02) of these runs. A higher uncertainty of about 0.5% due to the bunch spread is estimated for the 2011 data taking. This is explained by worse conditions in the beginning of the year, when the spread of µ across bunches reached 30%, which leads to a correction of up to 7% for some fills. The 2011 data taking at 7 TeV was affected by parasitic collisions due to a vanishing net crossing angle for one of the magnet polarity settings. This background ranges between 0.2% and 0.7% and a correction is applied averaging over time intervals of a few weeks each, during which data were taken under similar conditions. The average correction amounts to about 0.4% and since only about half of the 2011 running period is affected, an uncertainty due to parasitic collisions of 0.2% is assigned on the full period. In addition, the estimated uncertainty due to beam-beam background from 2012 is added in quadrature to obtain 0.24% uncertainty for 2011. The stability of the effective process is estimated using only data that is not affected by parasitic collisions.

Formalism for the luminosity of colliding beams
In a cyclical collider, such as the LHC, the average instantaneous luminosity of one pair of colliding bunches can be expressed as [5] where we have introduced the velocities v 1 and v 2 of the particles (in the approximation of zero emittance the velocities are the same within one bunch). The particle densities ρ j (r, t) (j = 1, 2) at position r = (x, y, z) and time t are normalized such that their individual integrals over all space are unity at all times. For highly relativistic beams colliding with a small half crossing-angle φ, the Møller factor 2c cos 2 φ 2c and one recovers Eqs. (2) and (3). The LHCb system of coordinates, which is used here, is chosen as a right-handed cartesian coordinate system with its origin at the nominal interaction point IP8. The z axis points towards the LHCb dipole magnet along the nominal average beam-line, the x axis lies in the horizontal plane, with x > 0 pointing approximately toward the centre of the LHC ring, and the y axis completes the right-handed system. This system almost coincides with the LHC coordinate system. Small angles due to the known LHC plane inclination and other magnetic lattice imperfections have negligible influence on the measurement of the overlap integral as only the crossing angles are relevant, not the individual beam directions. Up to a normalization factor, ρ bb (x, y, z, t) = ρ 1 (x, y, z, t) ρ 2 (x, y, z, t) is the distribution of interactions from the luminous region in the laboratory frame. If both ρ 1 and ρ 2 factorize as a product of a longitudinal and a transverse density (relative to the direction of motion of the bunch), the spatial distribution integrated over time 3 can be expressed as where n(z) is a shape factor which depends on z only. This relation between the distributions of beam-beam and beam-gas interactions is used in the BGI analysis. Determining the luminosity or the reference cross-section requires measuring the bunch population products N 1 N 2 , as discussed in Sec. 3, and evaluating the overlap integral Ω. We briefly describe the principles of the two methods that are used in this paper to determine the latter.

Beam overlap measurement methods
The first method was introduced by van der Meer to measure the luminosity of the coasting beams at the Intersecting Storage Rings (ISR) [11]. The method was further extended to measure the luminosity of a collider with bunched beams [12] and is the main method used to determine the luminosity at the other LHC experiments. The key principle of the VDM scan method is to express the overlap integral in terms of rates that are experimental observables as opposed to measuring the bunch density functions. Experimentally, the method consists in moving the beams across each other in two orthogonal directions. The overlap integral can be inferred from the rates measured at different beam separations, provided the beam displacements are calibrated as absolute distances.
A reaction rate R per bunch crossing is measured that is proportional to the luminosity and depends on the two orthogonal transverse separations of the two beams ∆x and ∆y. Measuring this rate relative to the revolution frequency ν rev (approximately 11245 Hz at the LHC) defines the parameter µ, which is the average number of reactions per bunch crossing. In the case where the spatial distributions of the beams can be factorized in the two coordinates x and y, it is sufficient to measure µ (and thus R) as a function of ∆x (at a fixed ∆y 0 ) and as a function of ∆y (at a fixed ∆x 0 ). One can show that the interaction cross-section is then given by The pair of separation values (∆x 0 ,∆y 0 ) is called the working point and is typically chosen to be as close as possible to the point where the luminosity is at its maximum. However, Eq. (15) is valid for any values of ∆x 0 and ∆y 0 . It can be shown that it is also valid in the presence of non-zero crossing angles [13].
The VDM method has the advantage of using a measured rate as its only observable, which is experimentally simple. The experimental difficulties of the VDM method arise mostly from the fact that the beams must be moved to perform the measurement. The exact displacements ∆x and ∆y in Eq. (15) steered by the LHC magnets are calibrated at each interaction point in a so-called length scale calibration (LSC). While the resulting corrections are typically of the order of 1%, some non-reproducibilities have been observed between two consecutive scans without being able to identify the cause. Another difficulty originates from beam-beam effects. When the beams are displaced, a change in β * (dynamic beta effect) and a beam deflection may be produced, which both influence the observed rate. The resulting corrections to the visible cross-section depend on the LHC optics, the beam parameters and filling scheme, and must be evaluated at each interaction point (see Sec. 7.6).
In addition, when performed with one vertical and one horizontal scan, the VDM method is valid only under the assumption that the distributions along the transverse variables x and y are independent, i.e. that the x (y) shape measured at a working point ∆y 0 (∆x 0 ) does not depend on the working point position. As will be shown in the analysis described here, this assumption is not valid at the required precision.
An alternative to the VDM scan method for measuring the luminosity is provided by the BGI method [15], which was first applied at the LHCb experiment in 2009 [20] and 2010 [10]. The principle of this method is to evaluate the overlap integral by measuring all required observables in Eq. (3) using the spatial distribution of beam-gas and beam-beam interaction vertices. The details of the measurement are discussed in Sec. 6. Measuring the shapes of stationary beams avoids changes due to beam-beam effects and other, non reproducible, effects due to beam steering. Furthermore, at the LHC the BGI measurements at a given IP (here at LHCb in IP8) can be made in parallalel to the VDM scans of other LHC experiments and can therefore be made more frequently.
On the other hand, while the β * and crossing angles used at the LHC do not impact the VDM method to first order, the BGI measurement relies on the vertex measurement to determine the bunch shape. Therefore, an increased β * is preferable to avoid limitations introduced by the detector resolution. At LHCb, in 2012, pp physics data were acquired at β * = 3 m, while the most precise BGI luminosity calibrations fills were carried out with β * = 10 m. The knowledge of the crossing angle is also important since the luminosity reduction due to the crossing angle has been as large as 20%. A non-vanishing crossing angle is necessary to avoid interactions between the main bunch and out-of-time charges captured in the next RF bucket, which occur near z = ±37.5 cm. Such displaced collisions, if present, must be disentangled from beam-gas interactions. They can be completely avoided by introducing a sufficiently large crossing angle. The VDM measurement can exclude interactions occurring away from the interaction point and is therefore less affected by these satellite collisions.
The VDM and BGI methods are complementary, in the sense that their systematic uncertainties on the overlap integral are highly uncorrelated, and a luminosity calibration performed with both methods in the same fill permits their systematic uncertainties to be constrained further. At present this can only be done at the LHCb experiment.
The analyses of the VDM and BGI luminosity calibration measurements presented here indicate that the observed luminosity profiles and vertex distributions are not consistent with Gaussian bunch distributions. It is found that a sum of two Gaussian functions ("double Gaussian" shape model) is sufficient to describe the x and y shapes of each bunch as well as the resulting luminous region. However, the joint two-dimensional transverse distribution of the bunches is found to be non-factorizable in the transverse coordinates. Therefore, as explained in Sec. 5.3, the transverse shape of the bunches is modelled with a sum of four two-dimensional Gaussian functions, which is in general non-factorizable.
In order to explain the full analysis of the present work, which involves a detailed fit model with a sum of Gaussian terms, it is useful to consider first the formalism for the ideal case of pure Gaussian beams and then describe the two-dimensional (non-factorizable) Gaussian model used in this work.

Luminosity in the case of purely Gaussian beams
The overlap integral in Eq. (13) can be calculated analytically when the single beam distributions ρ j (j = 1, 2) are the product of three Gaussian functions, each one depending on a single spatial coordinate m =x j ,ŷ j , orẑ j . The beam reference framesx j ,ŷ j ,ẑ j are right-handed systems and the longitudinal axisẑ j is assumed to be parallel to the velocity vector of the bunch v j . It is also assumed that theŷ j axes of the two colliding bunches are parallel to the y axis of the laboratory frame. The beam crossing plane, defined by the velocity vectors v 1 and v 2 , is here assumed to coincide with the xz plane. This condition was not respected only for the April 2012 fills. The relevant modifications of the formulae below are discussed in Sec. 7.2. We assume the bunches are centred at r j = (ξ xj , ξ yj , ξ zj ) at time t = 0, with a particle density function described by a normalized Gaussian function for beam j = 1, 2 and coordinate m =x,ŷ,ẑ, (16) where σ mj denotes the RMS of the corresponding Gaussian function.
Assuming that ρ j (r, t) = ρ j (r − v j t, 0), one can show that the overlap integral becomes where the following quantities have been introduced and ∆m = ξ m1 − ξ m2 (with m = x, y) are the transverse beam separations evaluated at the moment t = 0 when the colliding bunches are at the same z position. In the LHCb experiment, this z position (called z rf ) is defined by the LHC RF timing and needs not coincide with the location z = 0 of the LHCb laboratory frame nor with the geometrical crossing point of the two beam trajectories. The longitudinal position ξ lz of the luminous region is related to the beam separation ∆x and longitudinal bunch crossing point z rf with The index l indicates here a property of the luminous region, as opposed to a single beam property.
One can also show that the longitudinal size σ lz of the luminous region is related to the convolved bunch length σ z by Therefore, if one has a measurement of the transverse bunch size of the individual beams, of the crossing angle and of the longitudinal size of the luminous region, one can evaluate the longitudinal bunch convolution.

Double Gaussian shape model
A factorizable transverse beam distribution with double Gaussian projections has the density where g(m; µ, σ) indicates a normalized Gaussian function of the variable m with parameters µ and σ. By convention, the narrow (n) and wide (w) components in each projection have widths σ mn and σ mw , and weights w m and 1 − w m . The weights w ixiy in the sum representation are defined as The wide and narrow components are assumed to have the same mean, as supported by the data. Moreover, it is assumed that the 3-dimensional bunch distribution factorizes in a transverse (ρ(x,ŷ)) and a longitudinal component, where the latter is modelled with a Gaussian function. Non-factorizability can be introduced into the model in Eq. (21) by modifying the weights w ixiy from Eq. (22). For instance, in an extreme case, one can have w nw = w wn = 0 and w nn + w ww = 1, which corresponds to a sum of two 2-dimensional Gaussian functions. To allow for a gradual transition between this extreme case and the case of factorizable beams, it is useful to define the weights as a linear combination where the coefficient f parametrizes the factorizability. In the fully non-factorizable case (f = 0) there is no distinction between the x and y weights, thus the parameters w x and w y are (arbitrarily) combined in a single weight. As a result of the single beam model from Eq. (21), the shape of the luminous region and the overlap integral are described by a weighted sum of 16 components. Explicitly, the beam overlap integral is given by where I denotes the set of indices i x , i y , j x , j y , while w ixiy,1 and w jxjy,2 are the weights from Eq. (23) for beam 1 and beam 2. Each partial overlap integral Ω I is evaluated with Eq. (17).

Beam-gas imaging method
In this section, the BGI methodology and calibration results are presented in detail. A description of the data taking conditions (trigger settings, vacuum conditions) and of the event selection are given. Studies of the vertex position resolution and the unfolding method are presented. The resolution is determined from data, separately for beambeam and beam-gas interactions. An analysis of the resolution-corrected vertex position distribution is then performed, which uses both beam-gas vertices and beam-beam vertices to perform a global fit of the beam parameters (angles, luminous region length, longitudinal crossing position, transverse beam shapes). A double Gaussian model is used which also allows for a non-factorizability of the x and y distributions. Simulation is used to verify the soundness of the fit procedure. Several checks are made, based on data, to quantify systematic uncertainties. The list of dedicated luminosity calibration fills discussed in this paper for 2011, 2012 and 2013 can be found in Table 1. We focus on the 8 TeV pp data set taken in 2012, because it gives the most precise results. For the other data sets the analysis is similar and only the differences with the former are discussed in Sec. 6.2.

Data-taking conditions and event selection
For dedicated luminosity calibrations the LHC is filled with a low number of bunches, of the order of 50 per beam or less, and a large gap between bunches (∼ 1 µs) is maintained. Under these conditions the vacuum pressure at the interaction point is ∼ 10 −9 mbar, producing a beam-gas trigger rate of about 0.5 Hz per 10 11 protons. 4 Performing a BGI measurement with such low rates necessitates integration of a measurement over a period of up to 8 hours. Significant limitations in the precision are caused by the low event rate, beam drift and by emittance growth over the integration time. In order to mitigate this, in 2011 the VELO vacuum pumps located close to the interaction point were switched off, thus increasing the beam-gas rate by about a factor of four. To increase the rate further and to take full advantage of the BGI capabilities, the use of a gas injection system was proposed [15], developed and commissioned in the LHCb experiment [37]. A first gas injection test with circulating beams was performed in November 2011. When activating the system, neon gas is injected in the VELO, thus raising the pressure from about 10 −9 mbar to slightly above 10 −7 mbar. Once the injection is stopped, the nominal pressure of ∼ 10 −9 mbar is recovered within 20 minutes. The pumps are switched off during the gas injection. The effect of gas injection on the pressure and beam-gas interaction rate is shown in Fig. 15. For pp collisions at √ s = 8 TeV, the recorded beam-gas event rates per 10 11 protons were 98 Hz for beam 1 and 82 Hz for beam 2. The corresponding rates of the hardware trigger were about 2.1 kHz and 1.3 kHz for beam 1 and beam 2, respectively.
All luminosity calibrations acquired in 2012 and 2013 took advantage of a simplified activity trigger in all bunch crossing types. For these fills the hardware trigger used information from the SPD and PU sensors and from the calorimeters. The software trigger dedicated to the BGI measurement accepted events with a multiplicity larger than ten tracks. The vertex position z was required to lie within a range of −2000 < z < 400 mm for beam 1-gas events and 0 < z < 2000 mm beam 2-gas events. Of the triggered events in bb crossings with a vertex in the range −300 < z < 300 mm only a fixed fraction was accepted ("prescaled") to keep their total rate below 15 kHz, close to the maximum rate that can be recorded. Interactions with a transverse vertex position more than 4 mm from the nominal beam line were not accepted in order to reject interactions with the material of the RF foil of the VELO [29]. value when neon is injected is to be multiplied by approximately a factor four to account for the lesser gauge sensitivity to neon gas. The absolute pressure increase cannot be accurately determined as the gas composition in unknown before the start of injection. The beam-gas rate for triggered events increases from about 2 to 72 Hz/bunch for beam 1 and 1.3 to 51 Hz/bunch for beam 2.
All events acquired are reconstructed offline to determine their vertex position. A vertex has a number of tracks N Tr associated with it, each track having either a forward or backward direction with respect to beam 1. Their multiplicities are defined as N fwd and N bwd , respectively. Forward tracks correspond to particles moving towards the LHCb spectrometer and benefit from additional information such as energy and transverse momentum, which can improve the vertex resolution. The VELO acceptance for forward tracks vanishes when they originate from z 600 mm. Backward tracks are only recorded by the VELO and can be detected when they are produced from a vertex with longitudinal position z −95 mm. The longitudinal vertex distribution for all bunch-crossing types before applying further selection criteria is shown in Fig. 16. The acceptance limits for beam-gas events from beam 1 and beam 2 are visible as a sharp drop in rate. It can be seen that the acceptance drops in two stages for beam 1 at about 200 mm and 400 mm; this effect is related to the positions of VELO sensors in the forward region. The distribution of beam-beam interactions in the luminous region is reduced due to the prescale factor, which also affects beam-gas interactions in bb crossings located in |z| < 300 mm. The xz and yz vertex distributions for non-colliding bunches are shown in Fig. 17 for the first 1000 vertices per beam for two fills.
The BGI method relies on differentiation of beam-gas and beam-beam vertices. The reconstructed vertex is required to contain at least 10 tracks, N Tr ≥ 10. To exclude Crossing types ee (green), be (blue) and eb (red) contain only beam-gas events while the bb crossing type (black) contains beam-beam vertices in the central region and beam-gas vertices over the whole range. The effect of the prescale factor on the beam-beam interaction rate is visible. The exclusion region of ±300 mm for beam-gas vertices in bb crossings is visible as the step given by the prescale factor applied to these vertices inside this range.
interactions with material in the VELO RF foil, the vertex position must be within a radial distance of 2 mm from the beam line. The longitudinal position of a beamgas vertex must satisfy −1000 mm ≤ z ≤ 500 mm for beam 1-gas vertices and 0 mm ≤ z ≤ 1000 mm for beam 2-gas vertices. In bb crossings, these ranges are reduced to −1000 mm ≤ z ≤ −250 mm and 250 mm ≤ z ≤ 1000 mm, respectively. For beam 1 (beam 2)-gas vertices all tracks are required to be in the forward (backward) direction. Beam-beam interactions are selected requiring at least two tracks in both directions in addition to a minimum requirement of 25 on the vertex track multiplicity.
The measurement of the beam angles combines beam-gas interaction vertices from colliding and non-colliding bunches. For the measurements of the overlap integral, beamgas and beam-beam events from colliding bunch pairs are used. Beam-gas interactions in non-colliding and empty bunch crossings determine the ghost charge fractions and are used in the beam-gas vertex resolution determination.

Vertex position resolution
The knowledge of the vertex position resolution is a central ingredient for the measurement of the absolute luminosity as the observed vertex distribution is a convolution of the physical beam or luminous region with the detector resolution. To reduce the impact of the resolution on the measurement of the overlap integral, the beam optics for dedicated luminosity calibration fills had a β * value of 10 m (compared to 3 m used for physics production runs) and an increased transverse emittance. This resulted in a beam width about twice as large as the transverse vertex resolution.
The vertex resolution used for the BGI analysis is understood as the standard deviation of the distribution of the residual distance in one coordinate between the true vertex position and its measured position. The resolution depends on the number of tracks associated with a vertex, the longitudinal position and whether the vertex originates from a beam-gas collision with only forward or only backward tracks, or from a beam-beam collision with both forward and backward tracks. Although the value of the resolution in z is about ten times worse than that in the transverse directions, its effect can be neglected owing to the much larger luminous region length (∼ 60 mm). Therefore, only the resolution in the transverse x and y directions is considered here. The resolution measurement method (described below) has been verified with simulated events. Our studies show that a better resolution is predicted by simulation compared to the measurement from data. The difference could be explained by the imperfect alignment of the VELO sensors.

Resolution for beam-beam interaction vertices
The longitudinal distribution of selected beam-beam interaction vertices and the corresponding distribution of the number of tracks per vertex are shown in Fig. 18. Without external knowledge of the true position of primary vertices, the residual distance to the true position and therefore the resolution cannot be measured directly. Instead, one can measure the residual distance between two reconstructed vertices originating from the same collision by dividing the tracks forming one vertex into two randomly chosen samples ("split vertex method"). Defining the absolute vertex positions v 1 and v 2 (v = x or y) resulting from the primary vertex splitting and ∆v = v 1 − v 2 the distance between the two measurements, the Gaussian width σ ∆v of all ∆v distances is a convolution of both vertex resolutions and depends on the vertex track multiplicity N Tr,1 and N Tr,2 of each vertex where σ res,v (N Tr ) is the resolution for the vertex track multiplicity N Tr . Indices 1, 2 denote here the two vertices resulting from the splitting. It is observed that the position of the original primary vertex is not identical to the average position of the two sub-vertices. To minimize any systematic bias due to the resolution measurement procedure, the analysis of the beam parameters is performed using the average position of the two sub-vertices rather than the position of the original primary vertex. The dependence on N Tr is determined using the full z range of the luminous region. In a second step, the variation of the resolution as function of z is addressed by introducing distribution of all ∆v i measurements are fitted with a Gaussian function to determine the width of the distribution σ ∆v i . An example of such a distribution is shown in Fig. 19 (left) together with the result of the fit. In this particular case, the number of tracks in both sets being equal, the measured distribution width σ ∆v i is directly related to the vertex resolution σ res,v (N Tr ) using Eq. (25). All measured distribution widths σ ∆x i as function of (N Tr,1 , N Tr,2 ) combinations are shown in Fig. 19 (right). The resolution for a given vertex track multiplicity σ res,v (N Tr ) is obtained by fitting all distribution widths σ ∆v i with a least squares minimization constrained by Eq. (25). Results for the resolution as a function of vertex track multiplicity are shown in Fig. 20. As can be seen in Fig. 18 (right), the number of vertices with more than 120 tracks vanishes, limiting the resolution measurement up to about 60 tracks per vertex. The resolution for vertices with a larger number of tracks is obtained by extrapolation of a parametrization, The factor A, the power B and constant term C are measured by fitting all σ ∆v i measurements. Typical values are given for fill 2520 in the caption of The variation of the vertex resolution for beam-beam interactions as a function of z is determined by comparing the residual distribution with the parametrization. A correction factor F z is introduced for each z bin such that the residuals are minimized. The resulting values of F z range between 0.98 and 1.04 [37].

Resolution for beam-gas interaction vertices
In previous measurements [20], the beam-gas vertex resolution had to be extrapolated from beam-beam resolution measurements due to an insufficient number of beam-gas events. With neon gas injection, the increased rate allows a direct measurement of the beam-gas vertex resolution to be made for both beams on a fill-by-fill basis. The measurement principle is similar to that used for the beam-beam vertex resolution. The main differences reside in the vertex selection and in the z dependence of the resolution measurement.
The detector acceptance and extrapolation distance vary considerably within the ±1 m z range used for the BGI analysis, leading to different distributions of vertex track multiplicity and vertex resolution. Therefore, beam-gas interaction vertex resolutions were measured separately in five z ranges for beam 1 ( [800, 1200] mm). Results are shown in Fig. 21. The distributions of the residuals between the direct resolution measurements and the resolution parametrization are shown in Fig. 22.
There is no significant bias in the parametrization and the statistical spread of about 0.2 µm can be neglected. The z dependence is evaluated with a procedure similar to that used for beam-beam interactions. The beam-gas vertex resolution can be tested by measuring the single beam width at different z positions. While the physical beam width is unknown, its relative change as a function of z can be predicted from the machine optics and is expected to behave as assuming the waist position is at z = 0. This is called the "hourglass" effect. In fill 2520, measurements were performed with β * = 3 m optics, providing smaller beam sizes and a stronger hourglass effect than in other luminosity calibration fills. Therefore, the beam width measurements for this fill are more sensitive to the resolution. A deviation from the expected z dependence of the beam width was observed and additional resolution factors f z are introduced based on this deviation. These corrections factors are measured for each fill separately. With a β * value of 10 m the hourglass effect is negligible (less than 1% at 1 m). Similarly the beam width to resolution ratio is also larger and the correction factors have a smaller impact on the measured width. A final verification of single beam width measurements as function of z is shown in Fig. 23 after applying the correction factors. The resolution has been measured independently for all dedicated fills.

Resolution function for a sample of vertices
The average position of the two sub-vertices from a split primary vertex is used to measure the beam shapes. The observed vertex distribution is a convolution of the density function such as Eq. (21) The resolution is unbiased, i.e. the functions are centred at the origin. By dividing the full range of the distribution of the resolution estimates of all vertices in the sample into N g equal-sized bins, the widths and weights are determined by taking the centre and population of the bins, respectively. The weights c kx and c ky are the relative resolution weights of the effective resolution functions for x and y. Choosing N g = 3 gives a good description of the sample resolution-function; a larger number of Gaussian functions does not change the results.
Since the transverse distributions of single beams and of the luminous region are expressed as a superposition of Gaussian functions, describing the resolution function also as a sum of Gaussian functions results in an analytical expression for their convolution. In this approximation, the observed transverse distribution for one of the beam components is written as a superposition of convolved Gaussian density functions. Each intervening Gaussian width (σ mn , σ mw ) is replaced by a resolution-convolved width in the coordinate m = x, y: σ * kmn = σ 2 res,km + σ 2 mn and σ * kmw = σ 2 res,km + σ 2 mw .
A similar treatment can be applied to the transverse distribution of the luminous region, using in this case the products of single beam density functions as given by Eq. (14) [37].

Measurement of the overlap integral
The knowledge of the three-dimensional bunch shapes ρ j (x, y, z) is required to evaluate the overlap integral defined in Eq. (3). To determine the value of the parameters described in Eq. (23), a fitting procedure is performed that proceeds in several steps. First, the directions of the single beams are obtained from the positions of beam-gas interaction vertices using all bunches in the filling scheme simultaneously. These beam directions determine the crossing angles and are used to project vertex positions in relatively large z ranges onto a reference plane to ease the fitting procedure. The rest of the procedure is applied only to colliding bunches, taking each bunch-pair individually. To avoid problems with beam drifts and emittance growth, the data are grouped into data-taking periods of about 20 minutes. In a second step, transverse properties of single beams are obtained using the assumption that the shapes in the two transverse coordinates are factorizable. The parameter values obtained are used as initial values for the following step. The third step consists of a fit to beam-gas vertices of bunch pairs of both beams simultaneously, together with beam-beam interaction vertices in their luminous region. This fit is performed in two transverse coordinates separately, still assuming factorizability. In the following, fourth, step the two beams and their luminous region are fitted in both transverse coordinates simultaneously with the full two-dimensional model. The initial parameters of the latter fit are provided by the fit performed in the previous step. Finally, the z positions of the beam-beam vertices are used to determine the longitudinal properties of the bunches. At this point, all parameters needed to evaluate the overlap integral are determined. The individual steps of the procedure are described in more detail below.
The beam angles α m,j for beam j = 1, 2 and axis m = x, y in the laboratory reference frame are measured using beam-gas vertex positions. While the luminosity measurement is based on vertices in the colliding bunch pairs, vertices originating from non-colliding bunch crossings are valuable to measure the crossing angles as they cover the full z range owing to the absence of beam-beam background. In a first pass, a straight-line least-square fit is performed to all beam-gas interaction vertices weighting the positions according to their resolution and to an initial estimate of the beam width. An example of a crossing angle measurement using vertex positions directly is shown in Fig. 17. In a second pass, events are binned in 50 mm z intervals with centre z c . Their transverse vertex position v j,m is projected to z c using as initial estimates values of α m,j obtained in the first pass. A weighted straight-line fit is then performed through the transverse positions obtained by Gaussian fits to the distributions in the z intervals. The statistical uncertainty in the angles is about 10 −2 µrad. The half crossing angle, which is the angle of interest to measure the overlap integral, is defined as φ m ≡ (α m,1 − α m,2 )/2.
In the second step of the fitting procedure, the transverse shapes of each individual bunch are analysed using beam-gas interaction data from both beams. The data are divided into different z ranges to combine only data of similar resolution. Three slices for the beam-gas samples in the ranges −1000 mm < z < −250 mm (beam 1) and 250 mm < z < 1000 mm (beam 2) are chosen. The vertices within a slice are projected along the beam direction onto a plane perpendicular to the z axis. This coordinate translation neglects the hourglass effect by assuming a constant beam shape along the beam axis. This is justified since, with β * = 10 m, the beam is broadened by 0.5% over 1 m. The transverse beam-shape model is fitted to the data in the three slices simultaneously using the double-Gaussian density model of Eq. (21) and convolution with the resolution model.
For the third and fourth steps, the beam-beam data are divided into z slices and a new fit is performed using beam-gas and beam-beam data simultaneously. The transverse distribution of the luminous region changes as function of z. Since the number of events is sufficient, the analysis is simplified by taking many slices so that in each of these slices the fit model can be approximated by the function value at the centre of the z range. Thus, for the beam-beam vertices, the range −100 mm < z < 110 mm is divided into 18 slices. The properties of the luminous region follow from the single beam properties, thus no new shape parameters are introduced, while the values obtained from the single beam fits are used as initial values. Free amplitude parameters are introduced for each single beam and luminous region z slice to take into account the combined effect of trigger and reconstruction efficiency and absolute rates. The third step consists of a 1-d fit assuming a factorizable description in the x and y coordinate. As is shown below, a twodimensional model is needed for the transverse shapes of the beams that can accommodate non-factorizable distributions in the x and y coordinates. Such a fit is performed in the fourth step, where a 2-d fit is performed to both coordinates simultaneously.
The model for the transverse distributions in the two fit passes is similar. The singlebeam density function is defined in terms of two Gaussian functions for each coordinate.
A factorizability parameter f j (for beam j) is used as defined in Eq. (23). For the 1-d fits the values of f j (j = 1, 2) are fixed to 1, decoupling the two transverse coordinates. Following this model, described in Sec. 5.3, the observed transverse vertex distribution per beam for a given z range is fitted with the resolution-convolved width defined in Eq. (29). The fit parameters are, per beam j, the Gaussian parameters w m,j , ξ m,j , σ mn,j and σ mw,j for both axes m = x, y, the factorizability parameter f j , and a free amplitude A j,k per z slice k. Figure 24 shows 1-d global fit results for the first measurement of the first bunch pair of fill 2852 as example.
The next step introduces values of the factorizability parameters f j different from unity, and therefore a combined fit coupled in the two transverse coordinates is mandatory. The large number of parameters (18 beam parameters + 24 amplitudes) for the 2-d global fit requires good starting values. The results of the 1-d global fits are used as starting values for the final global fits except for the starting values for the factorizability parameters f j , which are set to 0.5. An example of a global fit result displaying only one z slice per beam and one luminous region slice is shown in Fig. 25 for the first bunch and first measurement performed with gas injection and a β * = 10 m lattice (fill 2853). Evidence for a significant non-factorizability of the beam shape is discussed further in Sec. 6.5. The χ 2 /ndf, with ndf the number of degrees of freedom of the fit, is typically between 1 and 1.1 for the 2-d fit. The non-factorizability of the beams in the transverse coordinates can affect the overlap integral by up to 3%.
Finally, to be able to calculate the overlap integral from Eq. must be known. The transverse offsets ∆m of Eq. (17) have to be evaluated at z rf , which is defined by the LHC RF phase. The values of the parameters σ z and z rf can be obtained from an analysis of the longitudinal vertex distribution of beam-beam interactions using the relations given in Eqs. (19) and (20). A fit is performed to the luminous region distribution of beam-beam interaction vertices for the same data sets as the transverse fits. Following Eqs. (14) and (21), the luminous region distribution ρ bb (z) is represented by the sum of sixteen Gaussian contributions. Complete factorization of the z dependence of the bunch distribution is assumed. The uncertainty introduced by this assumption is discussed further in Sec. 6.6. Under this assumption each contribution to the luminous region has a length, amplitude and longitudinal offset, which depend only on the transverse single beam parameters and the quantities z rf , the combination σ 2 z = σ 2 z1 + σ 2 z2 in Eq. (18) and an arbitrary overall amplitude parameter A l . The quantities z rf and σ 2 z1 + σ 2 z2 are common to all Gaussian contributions. The relative weights of the sixteen contributions follow from Eq. (24) according to the fraction of luminosity they carry and do not introduce new parameters. Thus, the longitudinal distribution of the luminous region is fitted with only three free parameters and makes use of transverse fit parameters fixed by the global transverse fit. Because the reconstruction efficiency in the VELO is not constant over the full luminous region, the z distribution is corrected for the efficiency obtained from simulation.
At this point, all parameters are measured and the overlap integral can be calculated following Eqs. (17) and (18). The statistical uncertainty is evaluated by sampling the multivariate normal distributions of the parameters using the fit results as mean values and the covariance matrix provided by the last fitting step. The resulting statistical uncertainty on the overlap integral evaluated per bunch pair and in 20 minute periods is typically less than 0.5%.

Generic simulation
The BGI method relies on an accurate beam shape description. To test whether the fitting procedure described in the previous section gives unbiased results, simulated data sets are created with a Monte Carlo method. The development of the more complex 2-d fit model described in Sec. 5.3 was motivated by the possible non-factorizability of the beams in the x and y directions and the observation that the 2-d properties of the beams could, in principle, be measured with beam-gas interactions. This capability was first tested with simulated data. Results showing evidence for beam non-factorizability are presented in Sec. 6.5.
Datasets of simulated vertices are generated for single beams and the luminous region as follows. Single beam vertex positions v m (in the axis m = x, y) are generated by sampling Eq. (21) at a fixed z = 0 position for both beams. A random z position v z,j is then assigned to each vertex in the range −1000 < z < −250 mm for beam 1 and 250 < z < 1000 mm for beam 2. The transverse coordinates of the vertices are translated to v z,j according to the beam direction. The z dependence of the reconstruction efficiency is implemented with a linear reduction of vertices as function of z. Per simulated dataset, about 5 × 10 4 vertices are generated per beam, similar to the number of events acquired per bunch pair during 20 minutes of data taking.
The shape of the luminous region and the overlap integral are sampled over space and time using the single beam transverse distributions ρ j (j = 1, 2) defined in Sec. 5 and the longitudinal bunch shape ρ zj , which is assumed to be a single Gaussian function. A sufficiently large sampling volume S vol = ∆x ∆y ∆z ∆t is chosen (±0.8 mm in x and y, ±250 mm in z and ±1.2 ns in t). A number N s of random samples of (x, y, z, t) ∈ S vol are generated uniformly over the volume. Vertices are retained according to the probability density ρ s bb (x, y, z, t) = ρ 1 (x, y, z, t)ρ z1 (z − ct) ρ 2 (x, y, z, t)ρ z2 (z + ct) , where the z and t dependence in ρ 1,2 just expresses a translation along the beam direction. Following the usual rejection sampling method, a randomly sampled vertex is retained if a uniform random number u assigned to it in a range [0, A bb ] satisfies u ≤ ρ s bb (x, y, z, t), with the arbitrary constant A bb ≥ max(ρ bb ). The numerical value of a generated overlap integral is calculated using the fraction N bb of vertices retained compared to the total number of samples N s generated in the volume S vol with Each simulated primary vertex generated with the above method is assigned a track multiplicity according to the distributions found in data. Using the resolution parametriza- tion measured with data as described in Sec. 6.2, each vertex is assigned a measurement deviation in x and y by sampling a normal Gaussian distribution. The generated datasets are stored with the same format and are processed with the same algorithms as used for the data. The fitting algorithms are tested with different beam parameters and are validated with simulated datasets before being applied to the data. An example generated with beam parameters similar to those found in the data is shown in Fig. 26. Detailed studies validate that the simulation input parameters can consistently be recovered with the fitting procedure.

Evidence of non-factorizable beam shapes
Discrepancies of the order of 3% are observed in visible cross-section measurements performed with the BGI method in the four July 2012 fills when fitted with a model factorizable in the x and y coordinates. Since the beam-gas interaction vertices provide a complete transverse view of the beams, the factorizability hypothesis can be tested.
A set of simulated data samples with non-factorizable beams (f 1,2 = 0) is fitted with the 2-d global fit model. A factorizable version of the model is obtained by fixing the factorizability parameters at unity. A fit with these parameters left free can describe also non-factorizable distributions. The effect of the pulls showing the difference between a factorizable and non-factorizable model is shown in Fig. 27. The fit converges towards the correct value of f j = 0 when left free. In addition, the (x, y) distribution of pulls reveals a  clearly visible cross-like structure when fitted with a factorizable model, which cannot fully describe the beam shape. The measurement performed on one data set acquired in July 2012 is also shown in Fig. 27 as an example. The fit correctly describes the beam shapes and the fit pulls are more uniformly distributed. The 2-d fit model converges towards non-factorizable beams and the pulls of the fit assuming a factorizable beam display the same structure as in the distributions simulated with non-factorizable beams. Data sets have been generated with the parameter f j = 0 and another set of samples with f j = 1. Results for the measured factorizability parameter f j are shown in Fig. 28 both for simulated samples and for data. In general, the fit algorithm can reliably measure the value of f 1,2 except when the beams are close to a single Gaussian shape, where the parameters f 1,2 have little or no meaning. The distribution of values of f j in the fits to the data shows a clear dominance of non-factorizable bunch shapes. The factorizability parameter is only meaningful if the beam has a double Gaussian shape in both transverse coordinates. If the beam shape is single Gaussian in one plane only, the beam is by definition factorizable in the model used and f j can not be measured. The ability to measure f j depends thus on the "strength" of the double Gaussian shape of the beams defined here as S j,m = 1 − σ main,m,j σ rms,m,j , with σ 2 rms,m,j = w m,j σ 2 mn,j + (1 − w m,j )σ 2 mw,j for beam j = 1, 2 and plane m = x, y. The indices n and w denote the narrow and wide width Gaussian component, respectively, while the width σ main is the width of the Gaussian (n or w) which carries the largest weight. A single Gaussian shape has therefore a vanishing strength parameter.
With few exceptions, all bunch pairs measured in April and July 2012 have a double Gaussian shape and S j,m is significantly larger than zero for both beams and planes. The November fills 3311 and 3316 are clearly different from the earlier fills, as all bunches have S j,m values smaller than 0.04.

Results and systematic uncertainties
Instantaneous luminosity values for each colliding bunch pair are evaluated using Eq. (1) and (2) with the overlap integral Ω measured with data integrated over intervals of about 20 minutes. The luminosity measurements per colliding bunch pair are used to evaluate the visible cross-section for specific observables (Sec. 4) with where µ vis is the visible average interaction rate for the reference observable. In the BGI measurement, the Vertex observable is chosen as reference owing to its time stability and low background. Cross-section results for the Vertex observable for all dedicated luminosity calibration fills in 2012 at √ s = 8 TeV with β * = 10 m and nominally head-on beams are shown in Fig. 29, together with a comparison between 1-d and 2-d fits of the transverse bunch properties. One observes that the statistical uncertainties are small and that the results using the 2-d fits are consistent between fills. Comparing results from both fit methods shows the importance of measuring the shapes in two dimensions to take the non-factorizability of the description between the two transverse coordinates into account; the 1-d fits do not display the same consistency. One notes also that the difference between the 1-d and 2-d fit results ranges between 1 and 3% in the July fills, for which the bunches clearly had double Gaussian shapes, while it is < 1% in the November fills, for which the bunches were almost single Gaussian (see Sec. 6.5).
The bunch shapes change over the course of a fill due to emittance growth and other factors, such as beam-beam effects or beam position drifts. Any change in the beam shape influences the overlap integral. Furthermore, the bunch population product decays over time, reducing the luminosity. In contrast, the cross-section is a physical observable and must be stable over time. The consistency of the measured cross-section, together with the corresponding overlap integral, rate and bunch population product is shown in Fig. 30 for one colliding bunch pair for two different fills. While the intensity product decay is typically smooth, the rate fluctuations follow the variations of the overlap integral. The figure also shows that the variation of the overlap integral in adjacent 20 minute intervals is very small (much less than a percent), indicating that the effect of e.g. emittance growth and beam drifts during the short intervals is negligible. In the following, sources of systematic uncertainties and their effect on the measurement precision will be described.
The impact of the beam-beam interaction vertex resolution on the cross-section depends on the transverse size of the luminous region. Comparing the results obtained with different is representative of the importance of the beam-beam resolution in the cross-section measurement. Here, the σ lm represent the measured values of σ lm and σ res,lm the beambeam resolutions (m = x, y). A value of R = 0 means that the resolution is negligible compared to the beam size. Cross-section measurements with different R values are obtained by using different cuts on the vertex track multiplicity and by using data acquired with a β * value of 3 and 10 m. This procedure gives four sets of results with different values of R.
The results for the six fills with β * = 10 m and two fills with β * = 3 m are combined in Fig. 31. Measurements performed with the larger β * and the best resolution (high track multiplicity cut) provide the smallest R value, while a worse resolution combined with the smaller β * results in a larger R value. The clear correlation between the cross-section Fill 3311 lasted about 4h30 with a beam intensity product decrease of less than 2%; in this fill the luminosity reduction is mostly caused by emittance growth. and the R value visible in Fig. 31 is an indication that the resolution is not perfectly understood. The cross-section obtained for similar R values and different β * are similar. This shows that the difference between the measurements at β * = 3 m and β * = 10 m can be attributed to the resolution description. The leftmost group of data points are the results presented in Fig. 29. Those results are obtained with the high track multiplicity cut and with a β * = 10 m beam optics. This combination provides the best measurement conditions and those measurements, called here σ c are used as central value for the final results. The cross-section σ e obtained by extrapolating the cross-section to R = 0 based on the β * = 10 m fills with the low and high track multiplicity cuts (blue and magenta leftmost measurement groups) is used to evaluate the uncertainty due to the beam-beam resolution. The difference of ∆σ = σ c − σ e = 0.93% between the central value σ c and the extrapolated value σ e is taken as systematic uncertainty. This is the largest single contribution to the uncertainty of the cross-section result.
As discussed in Sec. 6.2, a set of correction factors to the beam-gas interaction vertex resolution have been determined to reach consistent beam width measurements at different z positions. The necessity to include correction factors is an indication of additional systematic uncertainties. The overlap integrals have also been evaluated without the resolution correction factors. The full difference between the results with and without correction factors amounts to 0.55%, which is used as the corresponding systematic uncertainty.
A misalignment of the VELO detector can correlate the positions of the sub-vertices and degrade the vertex resolution. This effect is included in the systematic uncertainties assigned to the resolution as described above. However, detector misalignment can also affect the crossing-angle measurement. Various versions of the detector alignment are produced, all of which provide acceptable and comparable results in the overall alignment quality. The different alignment versions arise from introducing different sets of constraints to satisfy internal consistency checks. The same data set is reconstructed with all alignment versions and measurements of the half crossing angle are performed. In the x direction, differences of the order of 10 µrad are found, an order of magnitude larger than the statistical uncertainty. The crossing angle uncertainty in the y axis is about 3 µrad and has a negligible impact on the luminosity. The luminosity uncertainty from the crossing angle correction depends on the bunch width and length and is different for each bunch pair. The full range of 0.45% in the various measured cross-sections for different alignment versions is taken as systematic uncertainty related to this source. The different alignment versions have negligible impact on the bunch shape measurements.
The bunch shape model defined in Sec. 5.3 can describe all observed bunch shapes with a χ 2 /ndf close to one. Nevertheless, the accuracy of the fit procedure as well as the capacity to describe different shapes is verified. The simulation (Sec. 6.4) is used to generate datasets of vertex distributions with different input parameters in order to test the capability to reproduce the input values by the fit procedure. In addition, the systematic uncertainty arising from fitting a double Gaussian model to different simulated shapes is estimated. Single, double or triple Gaussian shapes, or a "Super Gaussian shape" 5 [43], are tried.
A given set of input parameters is used to generate a sample of statistically independent datasets, each of which is analysed using the same algorithms as for real data. All shapes give a χ 2 /ndf close to unity and can be considered to be well described by the fit model. The results of the measurements of the simulated data are compared to the input parameters on the basis of the value of the overlap integral. The difference in the results indicates a 0.5% systematic uncertainty due to the fit model and accuracy of the fitting procedure. The double Gaussian fit model used to describe the transverse bunch shapes does not allow a description of a possible third component. For example, a fraction of the protons measured in a bunch could be present in a wider Gaussian shape that would not be measured with a double Gaussian function. In this case, the tails of the measured distributions would have a larger population fraction than expected. The fraction of vertices in the tails beyond the double-sided 99 percentile predicted by the fit, is checked for all measurements with β * = 10 m and for the simulated datasets. The tail population is about 2% for the single beams and about 1.5% for the luminous region (while 1% is expected). This tail population, however, is also observed in the simulated datasets, indicating that the higher tail population is a result of the fitting procedure. Therefore, the corresponding bias is already included in the fit model uncertainty given above, and no additional systematic uncertainty is assigned.
Measurements of cross-sections for all colliding bunch pairs in the fills of 2012 with β * = 10 m are shown in Fig. 32. The measurements use the Vertex observable as reference and have an RMS spread of 0.54%. Cross-section results are also shown on a fill-by-fill basis in Fig. 32; the indicated error bar is the RMS of all measurements of the corresponding fill. Since the statistical uncertainty per measurement is significantly less than 0.5%, the observed spread of 0.54% is due to the combination of statistical fluctuations and additional systematic effects. The full RMS is taken as systematic uncertainty on the cross-section. This uncertainty covers uncorrelated bunch-by-bunch uncertainties such as the shape description, which is influenced by the fit model and detector resolutions, or uncertainties in the bunch population measurements.
Some structure remains visible in the measurement shown in Fig. 29, pointing towards additional systematic uncertainties or a less than perfect bunch shape description. Correlations between the cross-sections and the major variables entering into the cross-section measurement (interaction rate, bunch population product, overlap integral, luminous region z position, crossing angle corrections) have been checked and are found to be negligible.
As described in Sec. 6.3, the convolved bunch length σ 2 z1 + σ 2 z2 and bunch crossing position z rf are measured with a fit to the longitudinal distribution of the luminous region. For each colliding bunch pair measurement, all vertices in the range |z| < 250 mm are selected, regardless of the track directions. This selection reduces the effect of the z dependence in the VELO region and limits the distortion of the luminous region shape. On the other hand, some beam-gas interaction vertices can be included in the distribution. An additional measurement of σ 2 z1 + σ 2 z2 and z rf is performed by requiring that vertices contain at least two forward and two backward tracks to exclude beam-gas interactions. This requirement distorts its longitudinal shape as backward tracks with an origin z −95 mm are not detected by the VELO. The reconstruction efficiency vtx obtained from simulation and corresponding to this requirement is shown in Fig. 33 and is used to correct the raw distribution. An example of a fit to the longitudinal distribution of the luminous region with this requirement is also shown in the figure for one bunch-pair measurement. All cross-section measurements have been analysed with both track requirements. The difference between the results obtained with the two methods is 0.04% and is taken as systematic uncertainty related to the reconstruction efficiency.
Two dedicated BGI measurements were performed with beams displaced vertically with repect to each other by 180 µm in fills 2852 and 2853 with a β * = 10 m optics and were preceded and followed by periods with no displacement. The beam offset reduced the overlap integral by a factor of about four. The convolved bunch length σ 2 z1 + σ 2 z2 measured during those fills is shown in Fig. 34. Periods where the beams were offset show the same σ 2 z1 + σ 2 z2 values as when the beams were head-on. With a relative displacement of 180 µm, the luminosity is dominated by the interaction of the wide beam components of both beams. On the other hand, the luminosity is dominated by the interaction of the narrow beam components when the beams collide head-on. The equality of the convolved bunch length values for head-on and offset collisions is an indication that the wide and narrow bunch components share the same length, justifying the assumption that the z coordinate is factorizable. The observed difference of 0.05% in the cross-section for the periods with offset beams compared to the head-on beams is taken as systematic uncertainty related to the convolved bunch length σ 2 z1 + σ 2 z2 measurement. The beam-gas interaction rate is proportional to the residual gas pressure at the interaction point. The BGI method measures the beam shapes with beam-gas interaction vertices, and therefore assumes that the pressure is uniform in the transverse directions. The relative error induced by a pressure gradient is estimated by evaluating the effect of a distortion on the overlap integral considering a constant pressure along the y axis and a pressure gradient along the x axis at the experimentally determined limit. A measurement was performed in 2010 during fill 1422 to verify the homogeneity of the pressure in the x direction by displacing the beams by 0.3 mm and has been used in the luminosity measurement for 2010 [10]. The relative uncertainty on Ω introduced by the pressure gradient is at most 0.03%.
The VELO transverse dimensions, which fix the absolute scale of the vertex transverse distribution measurements, were checked in the laboratory at different temperatures on individual silicon sensors. The relative uncertainty on the cross-section determination due to the transverse scale uncertainty is estimated to be at most 0.05%. The FBCT response is approximatively linear with respect to bunch intensities. However, small deviations from linearity lead to a non zero offset when extrapolating to zero bunch intensity. This offset can be inferred from the combination of all cross-sections measured during a fill, taking advantage of the fact that the cross-section is independent of the choice of the bunch pair. An analysis of the FBCT non-linearity is performed as outlined in Ref. [26]. A fit is performed with the offsets and the improved cross-section value as free parameters. As expected, the χ 2 /ndf values are improved by this procedure. All cross-section results obtained after applying the FBCT offset corrections for fills with β * = 10 m and √ s = 8 TeV give a distribution with an RMS of 0.48% instead of 0.54% without the correction. The central value is changed by 0.04%. This small deviation is not applied to the final result, but is taken as systematic uncertainty related to a potential FBCT offset. The systematic uncertainties introduced by the bunch population determination are described in Sec. 3. A difference of 0.2% is observed between the background-subtracted interaction rate measurement with the restricted Vertex observable and the corresponding unrestricted observable. This difference is attributed to the background subtraction and assigned as systematic uncertainty. The reference cross-section used for physics data taking is based on the Track observable, which is affected by beam-gas interaction background when the neon gas injection system is used for the BGI calibration data taking. To transport the visible cross-section based on the Vertex observable to that based on the Track observable, their ratio is measured in periods without gas injection. A variation of 0.2% is observed, which is taken as systematic uncertainty. Using the relation between the Vertex and Track observable of µ Track /µ Vertex = 1.106, the final calibration result is σ Track = 60.62 ± 0.87 mb. A summary of all uncertainties is provided in Table 7. The values shown for the other energies will be discussed below.  the neon gas injection system was active. The luminosity measurement and evaluation of systematic uncertainties follows the same procedure as with √ s = 8 TeV data. Crosssection results for √ s = 2.76 TeV are shown in Fig. 35. The RMS of the measurements is 1.3%. A half crossing angle of 885 µrad was chosen to avoid collisions between satellite and main bunches. Differences observed in the data as compared to 2012 data are discussed below. The beams have a double Gaussian shape and are significantly non-factorizable (f j is close to zero). However, the factorizability is more difficult to measure than with √ s = 8 TeV data and the uncertainties on the f j parameter are larger. The beam-beam resolution has a small impact on the cross-section compared to 2012 (uncertainty of 0.40%), owing to the comparatively large transverse bunch size. On the other hand, the uncertainty related to the beam-gas vertex position resolution has significantly larger impact than in 2012, due to the lower number of vertices in the luminous region. The cross-section difference measured with and without all beam-gas resolution correction factors amounts to 1.31%. The uncertainty related to the detector alignment of 0.9% is estimated at twice the value obtained in 2012 because the crossing angle correction is about twice as large. The convolved bunch length σ 2 z1 + σ 2 z2 measurement plays a more important role compared to the √ s = 8 TeV data due to the larger crossing angle correction. An uncertainty of 0.1% is assigned to the determination of σ 2 z1 + σ 2 z2 and z rf and 0.17% to the reconstruction efficiency. The FBCT offset fit changes the cross-section by 0.05% and reduces the overall RMS to 1.1%. The systematic uncertainty related to the ghost charge amounts to 0.07%. Uncertainties for the observable background subtraction and fit model are taken from the 2012 measurements. The reference cross-section for the Vertex observable is 46.4 ± 1.0 mb. A summary of all uncertainties is provided in Table 7. As for the √ s = 8 TeV data, the reference cross-section used for physics data taking is based on the Track observable. Using the relation between the Vertex and Track observable of µ Track /µ Vertex = 1.135, the final calibration result is σ Track = 52.7 ± 1.2 mb. In the 2011 pp calibration, the main differences with the situation in 2012 are the absence of the neon gas-injection system and a different trigger configuration. Therefore, the beam-gas interaction rate is a factor 20 lower than in 2012. To partly compensate for the lower rate, measurements are performed in one-hour periods, potentially introducing effects of emittance growth and beam drifts. In addition, data are available in one dedicated calibration fill only. Cross-sections results for √ s = 7 TeV are shown in Fig. 36. The RMS of the measurements is 2.0%.
Collisions between satellite and main bunches were observed at z = ±375 mm due to the low half crossing angle of 270 µrad. Therefore, beam-gas vertices in the region −550 mm < z < 550 mm are discarded to exclude beam-beam vertices in the single beam selection. This requirement, however, reduces the number of vertices available for the single beam measurement. Furthermore, the remaining vertices are measured with a worse resolution leading to a bigger impact on the beam size measurement. The limited amount of beam-gas vertices reduces the accuracy of the resolution determination and the correction factors described in Sec. 6.2 cannot be determined. An uncertainty of 10% is assumed for the beam-gas vertex resolution leading to an uncertainty of 2.8% on the cross-section calibration.
The ratio of the resolution width to the beam size is smaller in 2011, potentially reducing the uncertainty resulting from the beam-beam resolution. However, the same uncertainty as for √ s = 8 TeV is conservatively assumed for the impact of the beambeam vertex resolution. Given the smaller crossing angle, the uncertainty related to the detector alignment effects is estimated at 0.3%. The beam factorizablility is more difficult to measure in 2011 due to the lower number of vertices. The beams are mostly non-factorizable (f j is typically close to zero). Owing to the low beam-gas induced background, the reference cross-section can be based on the Track observable directly without using the ratio Vertex to Track. A global satellite fraction correction of 0.78% is applied to the result shown in Fig. 36. The correction cannot be applied on a per bunch basis as the LDM instrument was not operational yet for this fill. A summary of all uncertainties is provided in This value is consistent with the result obtained with the BGI method described in a previous LHCb publication [10]. The improvement in the present result is due to the better bunch population measurement, while the uncertainties in the overlap integral are similar.  Table 1 and the scan parameters are listed in Table 8. Four x-y scan pairs were performed in April and six in July. In all scans the two beams were moved symmetrically, typically covering a ±6σ n range of beam separations, where σ n is the nominal beam width assuming nominal values of β * and transverse normalized emittance n = 3.75 µm rad. The last scan pair in each fill had a nominal working point at a relatively large offset ∆x 0 , ∆y 0 ≈ +2σ n . These offset scan pairs are only used for cross-checks and are not considered in the cross-section determination because of their high sensitivity to systematic effects (e.g. beam orbit drift, factorizability, linear correlations). The first and fourth scans in April were performed along x and y axes that are rotated with respect to the principal axes of the LHC optics. This particular scan pair is analysed, but not used in the final result.
Beam movements recorded with LHC beam position monitors (BPMs) upstream and downstream of LHCb are shown in Fig. 37. The BPM measurements are not used in the analysis as their time stability is insufficiently accurate. During a scan the beam separation values are calculated on the basis of a detailed model of the LHC magnets and are set by the accelerator control system. Dedicated length scale calibration scans are performed in each VDM calibration fill to experimentally verify and calibrate the beam displacements using the precisely known geometry of the VELO (see Sec. 7.5).
During VDM scan sessions, a large fraction of the available trigger bandwidth was allocated to randomly triggered bunch crossings; 20 kHz were devoted to the crossings with collisions, 2 kHz to the crossing slots where only one of two beams was present, and 0.5 kHz to the empty crossing slots. In addition, starting from the July 2012 session,  beam-gas events were recorded simultaneously. Due to the small beam-gas interaction rate of a few hundred Hz in total, the collected events are only used for a cross-check of the beam positions. The average decay time of the bunch population product N 1 N 2 was 36 (70) hours in the April (July) session. The value of N 1 N 2 changed by about 1% during a scan pair. Therefore, the rates are normalized by N 1 N 2 of each colliding bunch pair at every scan point by defining a specific average number of interactions per crossing To reduce the noise associated with the N 1,2 measurements, the data for each beam are approximated by a smoothing spline [44]. In addition to the bunch population changes, the luminosity stability may be limited by changes in the bunch profiles, e.g. by emittance growth. The luminosity stability is checked several times during the scans when the beams were brought back to their nominal position. The evolution of µ sp , averaged over bunch pairs, is shown in Fig. 38. The average luminosity decay time is measured to be 29 (58) hours in the April (July) session, which is largely due to the decay in the bunch population product. The average luminosity drop caused by emittance growth (drop of µ sp ) amounts to 1.2% (0.5%) during the entire calibration session in April (July). The scan points have been taken from lower to higher ∆x, ∆y values, therefore, the luminosity drop due to emittance growth effectively enhances the rate at negative values and reduces the rate at positive values, so that the net effect on the integrals in Eq. (15) cancels to first order since the curves are symmetric. Therefore, the systematic error due to the emittance growth is considered to be negligible.
As will be shown in Sec. 7.8, the measurement using the fifth scan pair in July is more sensitive to beam orbit drifts, because it was performed at a slightly offset working point and took longer than usual (42 instead of 24 minutes). Since the orbit drifts cannot be precisely corrected for, this scan pair is not used for the central value of the result but is accounted for in the systematic uncertainty. In total, for the cross-section measurement, two reference scan pairs are used from the April session and four from the July session.

Overlap integral model
The rate measurements from the orthogonal VDM scans performed provide no information on the factorizability of the single beams. However, the cross-section measurement is sensitive to the latter. In order to impose constraints on the factorizability, the BGI measurements performed in the same fills on the same bunch pairs can be used. To facilitate such approach, it is advantageous to use an identical model of the single beam densities for both BGI and VDM analyses. In this section, a model for the overlap integral as function of beam separation is developed based on the double Gaussian beam shape description from Sec. 5.3. It should be noted that the application of the classic VDM method in the case of factorizable beam distributions (see Eq. (15)) only requires a good empirical description of the VDM scan data.
The formalism discussed in Sec. 5.2 is only valid in the case of a horizontal beam crossing plane. However, in the April calibration the crossing plane was rotated along z, which necessitates the following extended treatment. In this more general case, the exponent in the overlap integral for Gaussian beams from Eq. (17) gains a dependence on the product of the transverse beam separations ∆x∆y and the normalization factor is modified. It is useful to extend the definition of the effective convolved widths of the luminous region Σ m (m = x, y) from Eq. (18) (taking cos φ m ≈ 1) as so that Σ y also has a term that depends on the corresponding half crossing angle in the yz plane φ y . The general expression for the overlap integral (in the case of pure Gaussian beams) is then given by where g is the normalized Gaussian probability density function and λ is the correlation coefficient, which is related to the bunch lengths and the half crossing angles by The variables λ 1,2 are the individual bunch correlation coefficients, which take into account a possible rotation of the principal axes of the bunches (around the direction of motion). It is seen from Eq. (39) that C xy is only non-zero when both φ x and φ y are non-zero, which occurs when the beam crossing plane is tilted (i.e. neither strictly horizontal nor vertical). The overlap integral in the double Gaussian model, which is discussed in Sec. 5.3, is given by Eq. (24) with each partial overlap integral expressed as Ω ix,iy,jx,jy = Ω(∆x, ∆y | Σ x,ixjx , Σ y,iyjy , λ(Σ x,ixjx , Σ y,iyjy )) , where i x , i y , j x , j y take the values n and w. Here, it is implicitly assumed that all Gaussian components have the same centres, same bunch lengths and same correlation coefficients (i.e. C xy andλ). The weights in Eq. (24) depend on the factorizability parameter f j and two weight parameters (w x,j and w y,j ) for each beam j = 1, 2.
The rate that is measured in VDM scans only provides direct information on the beam overlap integral and not on the single beam densities. Therefore, it is convenient to parametrize the model of the overlap integral using the effective convolved widths rather than the underlying beam width parameters. Only three of the four effective convolved widths per plane are linearly independent, thus two parameters are sufficient to quantify their relative magnitudes, with the choice for m = x, y. In terms of beam width parameters we have thus A m describes the asymmetry between the differences in the widths of the wide and narrow components of the two beams. The special case A m = 1 corresponds to σ m1,n = σ m1,w , i.e. the distribution of beam 1 is Gaussian in the m axis. In the parametrization of Eqs. (41) and (42) where the ratios r m,ij are defined as r m,ij ≡ Σ m,ij /Σ m,nn . For the actual fitting procedure, it is advantageous to use scale (or width) parameters that have a model-independent meaning. For example, the RMS of the luminosity profiles at ∆x = 0 and ∆y = 0 are such parameters, which can be easily estimated from the scan data to obtain starting values for the fit. Using the parameters defined above, the RMS S x of the luminosity profile at ∆y = 0 is given by x,nn I w I r y,iy jy r 2 x,ixjx (1 − λ 2 (Σ x,nn r x,ixjx , Σ y,nn r y,iyjy )) I w I r y,iy jy (45) and similarly for the other coordinate. The values Σ x,nn and Σ y,nn are obtained by solving the system of equations defined by the above equation. Finally, in the double Gaussian model, the shape of the overlap integral as function of beam separation is parametrized with the following 14 parameters S x , S y , R x , R y , A x , A y , f 1 , f 2 , w x,1 , w y,1 , w x,2 , w y,2 , C xy ,λ .
Only some of those parameters remain free for the cross-section determination as explained in the following section.

Cross-section determination
All colliding bunch pairs are analysed individually. The data are fitted simultaneously for pairs of x and y scans. The value of µ vis at each step k is described with where N 1,2 are the bunch intensities and µ sp be(eb) is the specific µ value of the beam-gas background for beam 1 (2). Two position parameters, ∆x 0 and ∆y 0 , were introduced to account for the fact that the luminosity may reach a maximum at a non-zero nominal separation (∆x, ∆y) due to an imperfect alignment of the beams.
The last two terms in Eq. (47) are due to beam-gas interactions and are proportional to the beam intensity and the residual gas pressure. The value of µ sp be(eb) is determined using events from data taken simultaneously for non-colliding bunches for each scan independently and is typically 2 × 10 −14 (1 × 10 −14 ). In most cases the pressure was very stable during VDM scans. During the first pair of scans in the July session a drop of about 10% was observed due to the neon gas injection system being used beforehand. In order to take this into account, an exponential dependence on time is assumed for the first pair of scans in July, while in all other cases constant specific background is assumed.
As already discussed above, due to the pattern of two orthogonal movements used in the scans, the data provide no information on the factorizability parameters and the linear correlation parameters. Therefore, in order to obtain the required additional information, we use the BGI measurements performed in the same fills as the VDM scans. The BGI analysis gives f 1,2 values compatible with zero. In the VDM analysis the f 1,2 parameters are fixed to zero. In this fully non-factorizable case, the two weight parameters w x and w y for each beam are perfectly anticorrelated as seen from Eq. (23). Therefore, only one weight parameter per beam is used for the VDM analysis. The value of C xy is computed using the BGI measurements of the convolved bunch length and crossing angles. No significant linear correlation in the beam distributions is observed, thus the value ofλ is set to zero. The available number of counts per scan point and the relatively coarse scanning grid do not allow measuring both R m and A m . The parameter R m describes the main features of Σ m,ij , Eq. (44), and is determined in the fit, while A m is fixed to zero for the determination of the central value of σ vis . The systematic uncertainties arising from the assumptions on the fit parameters are discussed in Sec. 7.7. In total nine parameters remain free, including the visible cross-section.
Initially, the fit is performed using the uncertainty estimates obtained from the data δµ vis = 1/N 0 − 1/N , where N 0 is the number of empty events in a total of N events. For small values of µ vis , these uncertainty estimates are correlated with the µ vis values themselves, thus biasing the data weights and the fit result. To mitigate this problem, the data are fitted a second time using uncertainties δµ vis = exp(μ vis )/N − 1/N , which are based on the predictionsμ vis of the first fit.
For presentation purposes we define a corrected specific interaction rate per crossing Then, the function σ vis Ω(∆x, ∆y) represents the fit to the corrected data points µ sp vis,i . An example of a fit for a single bunch pair is shown in Fig. 39 for one scan pair. The measured cross-section values and the values of χ 2 per degree of freedom, χ 2 /ndf, obtained for all bunch and scan pairs are shown in Fig. 40. Higher values of σ vis are obtained for the fifth scan pair in July. The associated uncertainty is discussed in detail in Sec. 7.8. For the April scan session, the values of χ 2 /ndf are on average higher than one. This can be explained by the smaller beam size and the fact that the uncertainty of the beam separation is not taken into account. Finally, the cross-section is obtained for each calibration session by calculating the weighted average of all measurements from reference scan pairs.
The sources of systematic effects that influence the result of the fit described with Eq. (47) can be grouped into two major categories. First, there can be systematic effects related to the inputs of the fit, namely bunch intensities (described in Sec. 3), rates (of signal and background) and beam positions (Secs. 7.4-7.6). Second, the influence of the VDM profile model, the technical aspects of the fitting procedure and discrepancies between repeated measurements are described in Secs. 7.7 and 7.8. The variation in the obtained cross-section among bunch pairs is used to estimate the uncertainty due to the relative bunch intensity measurement. The uncertainty due to a potential non-linearity of the FBCT device is estimated to be about 0.05% using the same method as described in Sec. 6.6. Moreover, by using the independent measurements from the ATLAS BPTX system [38], a discrepancy in the final result of about 0.1% is found, which is assigned as an additional systematic uncertainty.

Rate measurement
For the absolute calibration, the beam-gas related backgrounds are subtracted taking into account the difference in bunch populations. There is a statistical uncertainty associated with the measured specific background per proton used for subtraction. The value of the specific background is shared among cross-section fits of bunch pairs and introduces a correlation, which is taken into account when combining the measurements. This uncertainty is estimated to be less than 0.1%. The beam-beam related background is estimated to be 0.1-0.2% by taking the difference between the visible cross-sections of the restricted and the non-restricted Track observables. For the central value of the visible cross-section we use the less affected Track observable, while the full difference is taken as systematic uncertainty.
In the presence of a non-zero beam crossing angle, the luminous region position varies with transverse beam separation. The track reconstruction efficiency is not uniform as function of the primary vertex longitudinal position. Therefore, a correction to the µ vis values is applied using the same principle as described in Sec. 4.
The longitudinal position and size of the luminous region are measured for each step using a single Gaussian fit to the selected beam-beam vertices, taking into account the vertex reconstruction efficiency as a function of z. The quantities are linearly extrapolated to large separations, where the luminosity vanishes and not enough vertices are available (< 1000). The measured longitudinal parameters of the luminous region are shown in Fig. 41 and the calculated correction factors to the observed µ vis for each scan and step are shown in Fig. 42. It can be noted that the correction factors are approximately linear as function of separation. Therefore, no large effect on the VDM profile width is expected. The effect of the correction amounts to +0.32% (−0.03%) in the April (July) calibration. An estimate of the associated systematic uncertainty is made using a comparison with an unaffected observable (Calo). In the approximation of pure Gaussian beams, the rate at zero beam separation is proportional to σ vis /(Σ x Σ y ), see Eqs. (1), (2) and (17). Equivalently, the measured cross-section is proportional to the rate at zero separation and the product of the widths of the VDM profiles. It is useful to separate the correction factor dependence on ξ lz (Fig. 10) in a linear and a higher order part. The former only affects the rate at zero separation, while the latter mainly affects the widths. The residual (after correction) slope of the ratio µ Track /µ Calo as function of ξ lz is measured and normalized to the ratio at ξ lz = 0. The ξ lz value at the working point is multiplied by that slope to obtain the uncertainty on the cross-section due to the linear part of the correction. The uncertainty due to the higher order part of the correction is estimated by comparing the product Σ x Σ y measured using the Track and the Calo observable. The two uncertainties are summed linearly to obtain the systematic uncertainty of 0.13% (0.07%) for the April (July) calibration.

Length scale
The nominal beam separation values ∆x and ∆y are calculated from the LHC magnet currents at every scan step. An absolute calibration of the beam separation is made using the more precise VELO scale. In principle, each beam can have different calibration constants. However, since for all scans the beams were moved symmetrically (as opposed to one beam at a time), only the average (or common) length scale matters to first order for calibrating the separation. During the length scale calibration scans the beams were moved in several steps, with each measurement lasting from one to four minutes. The time intervals during which the beams were stationary are determined using the recorded states of the deflection magnets in the LHC logging data. Two methods are used in order to obtain the absolute length scale calibration. Beam-beam vertices are selected using randomly triggered events in bb crossings. In the 2013 calibrations the average number of interactions per crossing was very low. Therefore, the standard trigger requirement for beam-beam events in bb crossings was used in order to collect a sufficient amount of data. To remove background from material interactions, only vertices that lie within 2 mm from the median in each axis are retained. Potential beam-gas interaction background is reduced by applying a loose cut on the longitudinal position of the vertices |z| < 350 mm and by requiring at least one track in each direction. The vertex coordinates are binned in x and y. The bin width is chosen automatically based on the data according to the Freedman-Diaconis rule. 6 Data that deviates more than 5σ from the median is discarded.
An empirical model is fitted to the histogram for each coordinate and the mean of the model distribution is used as an estimator for the position of the luminous region. A sum of two Gaussian functions, where the two means are not required to be equal, is found to fit the data well in all cases. The distribution of the luminous region is not expected to change during the scans if the beam separation is kept constant. Using this property as a constraint, an additional global fit per scan is performed, which has a single set of shape parameters with the exception of the mean, which is independent for each step.
A simultaneous drift of both beams in the same direction affects the length scale measurement. The drift is estimated using the steps before and after the scans, when the beams are nominally centred. Differences of the luminous region position before and after a scan are attributed entirely to simultaneous beam drift. A correction is applied for each step during the scan by using a linear interpolation. Since the beam separation is kept constant, beam-beam effects do not influence the measurement.
The measured luminous region position for each step is fitted against the average of the nominal beam positions using a weighted linear fit. The fit is performed simultaneously for all colliding bunch pairs using independent intercepts and a common slope. The results of this fit are shown in Fig. 43. The slope parameter of the fit is an estimate of the length scale calibration factor.
Beam-gas imaging method: Individual beam positions are measured using beamgas interactions. The neon gas injection system was used to enhance the beam-gas interaction rate. The beams were moved similarly to VDM scans, but in a few large steps. Beam-gas vertices for each beam are selected using the corresponding standard trigger requirements for beam-empty and beam-beam crossings. In order to remove background from material interactions, only vertices that lie within 2 mm from the median in each axis are retained. Potential beam-beam interaction background is reduced by requiring no tracks in the opposite direction of the corresponding beam. Consistent comparison  Beam-gas imaging 0.29 ± 0.06 0.69 ± 0.05 between measurements from colliding and non-colliding bunches is ensured by imposing a cut on the longitudinal position of the vertices |z| > 300 mm. To obtain the single-beam offsets, the x and y coordinates of the beam-gas vertices are fitted with a straight line as function of z. The weight assigned to the individual vertices takes into account the beam width and the vertex resolution. Individual beam positions are measured before and after the scans, when the beams are nominally centred. No significant beam drifts are observed at the 0.5 µm level, thus no correction is made. The beam-beam deflection effect influences beam positions in bb bunch crossings. We take this into account as described in Sec. 7.6, where Σ values are taken from the VDM scan measurements. Taking this effect into account improves the consistency between measurements from colliding and non-colliding bunches.
The measured beam position for each step is fitted against the nominal setting using a weighted linear fit. The fit is performed simultaneously for all bunches in each beam using one slope and one intercept parameter. The results of this fit are shown in Fig. 43. To obtain an effective common length scale, the average between both beams is taken.
The length scale calibration constants are summarized in Table 9. A statistically significant discrepancy of about 0.5% is observed between the calibration factors for the x axis obtained using the two methods in July. The origin of this difference is not understood. Therefore, we assign 0.5% uncertainty to all length scale calibrations. It is assumed that this uncertainty is correlated between calibrations. The statistical uncertainties are small compared to the discrepancy and are neglected. The calibration constants used for July are those obtained by the constant separation method. Thus, we maintain consistency of treatment with other calibration sessions, for which we have not performed a beam-gas imaging calibration.

Beam-beam effects
The electromagnetic interaction between charged particles of two colliding bunches is called beam-beam effect. There are two aspects of this interaction that affect the VDM scan measurements. The first, called dynamic β effect, is the result of the mutual (de)focusing of the two colliding bunches. It leads to a change in β * (and thus beam size), which depends on the transverse beam separation. Therefore, the transverse size of the bunches is not constant during the VDM scans. Secondly, the closed orbits of the bunches are distorted by the angular kick induced by their electromagnetic repulsion. This beam-beam deflection effect has a different magnitude depending on transverse separation, thus distorting the scan profiles.
The so-called beam-beam parameter that quantifies the strength of the dynamic β effect in the axis m, for a normally distributed bunch in beam j is given by [46] ξ m,j = α c 2π where j denotes the bunch in the opposite beam, N j Z j is the bunch charge in units of elementary charge, β * m,j is the value of the β function at the IP, E j is the ring energy setting (particle energy divided by particle charge) and σ m,j are the bunch sizes. The function ξ rel m was modelled [47] using the MAD-X optics software [48] and is shown in Fig. 44. It is independent of the machine optics parameters and the bunch properties. For small separations ∆x and ∆y the value of the ξ rel m is close to unity and approximately constant. Therefore, the beam-beam force is approximately linear, resembling the force of a quadrupole field. However, the force becomes non-linear for large separations. The above equation is also valid in the case of ion beams.
The ratio between the perturbed β * and the unperturbed β * 0 as function of separation is given by where Q m,j is the machine tune of beam j in the axis m. While there is a recursive dependence between ξ and β * , this is only a second order effect and can be neglected by using the nominal β * 0 in Eq. (49). Possible collisions at other IPs also contribute to the dynamic β independently from the local beam separation. Such contributions, provided beams are not moving at other IPs, are constant and effectively modify the nominal value of β * . Since these modifications are much smaller than the uncertainty in the value of the nominal β * itself, they are neglected.
To calculate the effect on the luminosity, the relative change of β * is computed with respect to a reference value. Applying a correction to the rate (or the overlap integral) effectively removes the dependence of bunch shapes on separation. The choice of the reference β * value is arbitrary and determines the values of the effective separationindependent bunch parameters. It can be shown that this choice has a negligible impact on the net correction. It is beneficial to use β * (0, 0) as the reference, in order to enable consistent comparisons of VDM profile parameters from different scans. The relative β * change (see Fig. 45) is calculated for each colliding bunch pair using the parameters of the nominal optics, the measured bunch intensities and the convolved widths. In this case, the shape of the bunches is approximated with a single Gaussian distribution and the size of the colliding bunches in each pair is assumed equal. The relative change of the luminosity is obtained using the perturbed bunch widths. Finally, the correction factor to the rate is obtained, which is shown in Fig. 46.
The beam-beam angular kick causes a deflection of beam j, which is calculated numerically using the formalism of Ref. [49] where N j Z j is the bunch charge in the opposite beam in units of elementary charge, i is the imaginary unit, F 0 is a complex function and Gaussian bunch profiles are assumed. The shift of the closed orbit in the m coordinate for beam j is given by while the effect on the beam separation is the sum of the individual beam shifts x scan y scan The corrections that are added to the nominal beam separation are shown in Fig. 47. The procedure outlined above was verified by the direct observation of beam-beam deflections with LHC orbit data [50], which were found to agree well with the expectations. The systematic uncertainty of the correction for beam-beam effects is estimated by numerically propagating the uncertainties of the input parameters. In addition to β * and machine tune, uncertainties are assigned to the ratio of bunch widths and the convolved bunch width to account for the assumption of identical bunch shapes and Gaussian bunches,  respectively. The assumed uncertainties on input parameters and their correlations are listed in Table 10. The uncertainty is similar for the two 8 TeV pp calibrations and amounts to about 0.3%.

Fit model uncertainty
The choice of a particular model to describe the single beam shapes gives rise to a systematic uncertainty. The bias in the case of the BGI method is studied by simulating various deviations of the single-beam distributions from the assumed model as described in Sec. 6.6. We assume that the same uncertainty of 0.5% applies to the VDM calibration and that it is fully correlated between the methods. In addition, the fitting procedure applied to VDM scan data might introduce a bias. The potential fit bias is estimated to be less than 0.2% by applying the analysis to simulated VDM scans with experimental conditions and beam properties similar to those of data. The linear correlations in VDM scan profiles can be non-zero due to a tilt of the beam crossing plane or due to non-zero linear correlations with respect to the scanning axes of the transverse bunch distributions in Eq. (39). The latter are found to be small and are neglected. The bunch crossing plane is nominally horizontal in all calibration fills, except in April, when the tilt was approximately −21.5°. Taking into account the correlation, the April calibration result changes by about +0.3%. The effect of the correlation is most pronounced if the scan pair crossing point is not centred at zero, as seen from the exponent factor in Eq. (38). For the nominal scans, the correction from the exponent factor is negligible. However, the correction to the fourth (offset) scan pair in April amounts to −8%, improving significantly the consistency with the reference scans. Rather than propagating the uncertainty of the fixed parameter C xy to the fit result, half of the effect of the correction is used as an estimate of the systematic uncertainty.
The major uncertainty of the VDM calibration is due to the assumptions made on the parameter values as described in Sec. 7.3. In order to estimate the corresponding bias, the fits are redone without fixing parameters but using constraints obtained from BGI measurements. The information that is used from the BGI is limited to the factorizability of the beams, thus the VDM calibration is as independent as possible. In particular, constraints are not applied directly using dimensional parameters (widths) from BGI measurements, which are subject to resolution and alignment uncertainties. Moreover, for each scan pair, all bunch pair fits are subject to the same constraints, with their width taking into account the spread of the parameter values among bunch pairs.
First of all, Gaussian constraints on the factorizability parameters f 1,2 are added to the fit. However, these constraints are not sufficient to ensure that the beams are non-factorizable. This is easily seen by considering sets of degenerate parameters for which the beam distribution is Gaussian in at least one of the axes. Therefore, a generic measure of the transverse factorizability is introduced as where ρ j is the transverse beam distribution of beam j = 1, 2. The measure v j is zero if ρ j is factorizable. The values of f 1,2 and v 1,2 measured with BGI, as well as the predicted values at the time of the VDM scans are shown in Fig. 48. The values of v 1,2 , being functions of the individual beam distributions, cannot be computed using only the parameters describing the overlap integral from Eq. (46). Therefore, two additional parameters for each coordinate m = x, y are added to the overlap fit model, which are sufficient to compute v 1,2 : Gaussian constraints of 10% around the measured value from BGI are applied to C x,y . Additionally, to avoid unphysical sets of beam parameters, beams are required to have similar transverse sizes in each coordinate. This requirement is ensured by imposing a weak Gaussian constraint of 0 ± 0.1 on the asymmetry of the RMS of the two beam distributions. Such constraints are justified by the equal design β * and emittance of the two beams. The asymmetries estimated from the BGI measurements for all bunch pairs lie within ±1σ, indicating that the width of the constraints is sufficiently large. The VDM data for each bunch and scan pair are fitted using the additional parameters and the constraints described above. Moreover, profiles of χ 2 /ndf are obtained as function of σ vis . The minima of the χ 2 /ndf profiles are very close to the values of σ vis obtained with the fit where σ vis is a free parameter as seen in Fig. 49. The new values of σ vis are averaged for each scan pair and the obtained value is compared to the baseline cross-section for that pair. A difference of 0.6% and 0.9% is found for the April and the July calibration, respectively. The difference is attributed to the fact that f 1,2 and A x,y are fixed for the baseline fits and is assigned as a systematic uncertainty.
Another systematic uncertainty arises from the uncertainty in the central values of the constraints obtained from BGI analysis. To estimate this uncertainty, the constraints on f 1,2 and v 1,2 are varied by one standard deviation. The cross-section obtained with the modified constraints is compared to that obtained with the nominal constraints. The difference of 0.3% (for both April and July calibrations) is taken as an additional systematic uncertainty, which is considered fully correlated with the BGI result. The average χ 2 /ndf profiles for the reference scan pairs are shown in Fig. 50 together with the relative contributions of the VDM data and the constraints.

Reproducibility
An important source of non-reproducibility of VDM scan measurements is the drift of the beam orbits. In addition, in the presence of a non-zero crossing angle, a drift in the longitudinal bunch crossing point, z rf , can influence the beam separation. While the effect of the latter is found to be negligible, the former has potentially a sizeable effect, which cannot be reliably corrected. For the following discussion data are averaged over bunch pairs since these effects are common for all bunches. The relation between the position of the luminous region and the beam separation from Eq. (19) can be used to estimate orbit drifts and z rf values in the case of non-zero crossing angle. In the general case of a crossing angle in both xz and yz planes, the luminous region centre for Gaussian beams is where k contains the beam widths and angles and ∆x is the separation in the crossing plane, which is given by ∆x = ∆x cos ψ + ∆y sin ψ , tan ψ = tan φ y tan φ x .
Here ψ is the tilt of the crossing plane, which is approximately −21.5°in April 2012 and negligible for all other calibrations.
In accordance with Eq. (19), the value of z rf is measured by interpolating the longitudinal position of the luminous region ξ lz at the separation that corresponds to the maximum  luminosity. It is implicitly assumed that z rf does not change and that there is no orbit drift during a scan pair. When the scanning axis is orthogonal to the beam crossing plane (e.g. y scans in July) this measurement cannot be made, as ξ lz is not expected to change during the corresponding scan. The measured z rf values, which are shown in Fig. 51, are found to be consistent between scans and with those obtained with the BGI analysis in the same fills. Therefore, it is assumed that z rf values do not vary during a scan and the associated uncertainty is negligible.
The overall drift of the beam separation between scan pairs is automatically taken into account in the VDM analysis, since it is effectively only a shift of the nominal separation where the luminosity is maximal. However, the drift during a scan pair can introduce a bias to the measurement. It is useful to split the drift into two components: a "slow" one that corresponds to the time scale of a scan pair and a "fast" one that corresponds to the time scale of a step.
The slow beam separation drift can be estimated from the fitted position of the VDM profile maximum. One estimate per coordinate is obtained for each scan pair, as shown in Fig. 52.
The fast component of the beam separation drift is more difficult to estimate. The presence of a crossing angle enables an estimation of the drift in the beam separation in the crossing plane ∆x . This is possible because of the correlation between the ∆x and the z position of the luminous region as seen in Eq. (19). While the latter is strictly true It is assumed that there is no significant orbit and z rf drift during the reference scans. The data for each scan direction are fitted with a smooth function. In July the crossing plane is orthogonal to the y axis, thus ξ lz is not expected to change during the y scans and the data are omitted. The non-linearity of the curves is due to second order effects, which are not expected for pure Gaussian beams. only for Gaussian beams, the exact form of the function ξ lz (∆x ) can be estimated by fitting the measured ξ lz as function of ∆x with a smooth function. The data and the estimated dependence are shown in Fig. 53. The deviations of the data points ξ lz from the curve give an estimate of the beam drift. This approach is only reliable when the reference curve is obtained by averaging enough measurements. Therefore, the drift is estimated only for the three nominally head-on steps for each scan, as shown in Fig. 54.
To estimate the effect of the orbit drifts on the cross-section measurement a simulation approach is used. The fast component of the drift in each coordinate is modelled with a probabilistic process according to a Brownian motion. The diffusion coefficient, which is the parameter of the model, is estimated from the measured fast component of the drift in ∆x . The value is found to be 0.005 µm 2 s −1 and it is assumed to be equal for both coordinates and constant during a fill. The separation drift model is constrained to the   Fig. 55. The systematic uncertainty is evaluated by analysing 400 statistically independent simulations for each scan session, where the simulated drift is added to the nominal beam separation. The average bias on the cross-section and its RMS are summarized on Fig. 56. It is seen that the average bias, which is mainly due to the slow component of the drift, is small. On the other hand, the uncertainty in the bias, which is represented by the RMS of all simulations and is driven by the fast drift component, is sizeable.
Since the fast drifts are modelled to be independent for each scan pair, the uncertainty on the average cross-section is reduced. The uncertainty due to the slow component is assumed to be correlated, thus the estimates for individual scan pairs are averaged. The two uncertainties are added in quadrature to obtain a total uncertainty of 0.54% and 0.21% in April and July, respectively. The non-reproducibility can be estimated directly by the deviations from the average of the cross-section from individual scan pairs. The maximum deviations for reference scan pairs observed in April and July are 0.03% and 0.30%, respectively. Assuming that the deviations are mainly due to the drift, and in order to avoid double counting, only the larger value of the drift estimate and the maximum deviation from the average is taken as a systematic uncertainty. It is seen in Fig. 56 that the scan pairs with an offset working point are more sensitive to random beam drifts, which can be explained by the large derivative of the VDM profile at non-zero beam separation. This is the main reason for excluding offset scans from the cross-section measurement. Moreover, since scan pair number five in July had a working point offset by about 30 µm in x, it is also more sensitive to beam drifts. Therefore, it is not included in the determination of the central value of the cross-section. Considering the reference scan pairs (one to four) in July, the uncertainties of the drift bias are too large compared to the observed fluctuations. This suggests that the value used for the simulation parameter, on which these uncertainties depend directly, is overestimated or the employed model does not describe well the beam drifts. Therefore, the higher cross-section measured from scan pair number five (see Fig. 40) cannot be explained by random beam drifts and the deviation from the average cross-section (using pairs one to four) of 0.8% is taken as a systematic uncertainty to account for a potential unknown source of non-reproducibility. The latter systematic uncertainty is assigned to all pp calibrations and is considered to be correlated among calibrations.
The effects of the applied corrections are summarized in Table 11. The values for the other energies will be discussed below (Sec. 7.10).

Results
The reference cross-section for pp collisions at √ s = 8 TeV for the Track observable is determined by computing a weighted average of the results from the calibrations in the April and July fills. It is assumed that most of the systematic uncertainties are fully correlated in order to avoid underestimating the uncertainty on the combined result. The individual calibration results and the average reference cross-section are shown in Table 12.
A list of all uncertainties for the July and April calibrations is provided in Table 13 along with estimates of their correlation. The values for the other energies are discussed below.

Summary of other van der Meer scan calibrations
The VDM analysis presented in this paper focuses on the calibration of the reference cross-section at  Tables 1 and 8. Three symmetric x-y scan pairs were performed. The average decay time of the bunch population product amounts to 70 hours and the luminosity drop caused by emittance growth is negligible.
The VDM analysis at 7 TeV is performed in exactly the same manner as for the 8 TeV calibrations, which are discussed in detail. Given the similarity of the conditions, the systematic uncertainties on the fit model, the non-reproducibility and the length scale calibration are directly translated from the 2012 calibration. The corrections and the systematic uncertainty due to beam-beam effects are found to be similar, owing to the similar beam energy and bunch intensities. The dominating systematic uncertainties are related to the constraints on the beam factorizability as obtained from the BGI analysis (Sec. 6.7). By using these constraints instead of making parameter assumptions (see Sec. 7.3), the measurement is shifted by 0.9%. The systematic uncertainty arising from the uncertainty in the central values of the constraints is found to be 0.8% and it is fully Individual calibrations of the pPb and Pbp reference cross-sections were performed in separate fills (3505, 3542). Due to the limited time of the ion runs, normal fills for physics data taking were used. Therefore, the filling schemes (see Table 1) were not optimized as for dedicated luminosity calibration fills. The number of colliding bunches was 38, instead of usually 16 in dedicated fills. Moreover, the peak µ vis value was very low, < 0.02. In order to decrease statistical uncertainties, the random trigger rate was doubled to 45 kHz during the pPb calibration. The statistical uncertainties on the cross-section measurement are significantly larger compared to pp calibrations and amount to 0.25%.
The statistical uncertainty on the cross-section measurement per colliding bunch pair and per scan pair is about 2%. Therefore, a systematic structure in the fit residuals, which leads to a measurement bias of that order, may not be revealed by the χ 2 /ndf values of individual fits. The fit quality of a model is estimated from the sum of the residuals of all bunch pairs as function of beam separation. The aggregated residuals show a statistically significant structure if a pure Gaussian model is employed. It is found that two factorizable empirical models (a double Gaussian function with a negative weight for the narrow component and a Gaussian function multiplied by an even fourth order polynomial) fit the VDM scan profiles well and give similar cross-section values. In addition, two non-factorizable empirical models for the two-dimensional luminosity profiles have similar fit quality. The first non-factorizable model is a sum of two twodimensional Gaussian functions with a negative weight for the narrow term. The second non-factorizable model is a xy rotationally symmetric Gaussian function multiplied by an even fourth order polynomial and scaled by the x and y profile widths. The two nonfactorizable models give similar cross-section values, which are about 2% lower than the values from the factorizable models. The VDM scans performed provide no information on the factorizability. Therefore, the full range of obtained cross-section values is considered by taking the central value at the middle and assigning half of the span as a fit model uncertainty.
An independent length scale calibration was performed for both fills. The LHC optics setup was almost identical in the two proton-lead beam configurations, thus the calibration constants are expected to be identical. However, a difference of about 1% is observed and is assigned as a systematic uncertainty to both calibrations.
The uncertainty due to the beam orbit drift is estimated to be 0.67% and 0.88% for the pPb and the Pbp calibration, respectively. The deviation of the visible cross-section measurement from the average in repeated VDM scans amounts to 0.23% and 1.31% for the pPb and the Pbp calibration, respectively. Following the procedure from Sec. 7.8, only the larger value of the uncertainty due to orbit drift and the deviation from the average is taken as a systematic uncertainty. Beam-beam effects have a small impact on the pPb and Pbp calibrations owing to the low bunch intensities of (1.3-1.7)×10 10 elementary charges. The low µ values effectively increase the fraction of beam-gas and beam-beam induced backgrounds. The corresponding uncertainties are evaluated as described in Sec. 7.4 and amount to 0.3% and 0.7% for beam-gas and beam-beam background respectively. It is assumed that there is no correlation between the uncertainties of the pPb and Pbp calibrations when it is partial. This is done to avoid underestimating the uncertainty on the ratio of the luminosities of the two data samples, which enters into ratios of cross-section measurements. A summary of all uncertainties is provided in Table 13. The reference cross-section for the Track observable is 2.126 ± 0.049 b and 2.120 ± 0.053 b for the pPb and Pbp mode, respectively. The compatibility of the two results indicates that the Track observable has similar efficiency for the two beam modes.

Summary and conclusion
Since some of the luminosity calibrations have been performed with both the VDM and BGI method (pp cross-section at √ s = 8 TeV and √ s = 7 TeV), the best result is obtained by computing an average of the two methods, taking into account the correlation between the systematic error sources. A summary of the final reference cross-section results is presented in Table 15.
The BGI and VDM calibrations of the visible pp cross-section at √ s = 8 TeV achieve very similar precision. Therefore, for simplicity the two results are combined with equal weights; the correlated components of the uncertainties are averaged linearly and the others are averaged in quadrature. The different sources of uncertainty are compared in Table 14, indicating also whether they are correlated or not. The result of the combination is given in Table 15. If the uncertainty on the propagation to physics data (0.31%, see Sec. 4) is included, the total uncertainty on the luminosity is 1.16%. The latter uncertainty is valid if the complete 2012 data set or a major part of it is used. In some cases, for small partial sets, the uncertainty may be different. Table 15: Results of the luminosity calibration measurements. The total uncertainty on the luminosity calibration (last column) is the sum in quadrature of the absolute calibration uncertainty (fourth column) and the relative calibration uncertainty (fifth column). The weights used to obtain the average absolute calibration at 8 and 7 TeV (pp) are given in the third column. The part of the uncertainty that is correlated between VDM and BGI calibrations is shown in parentheses (fourth column).  Figure 57: Measurements of the visible cross-section for the Vertex observable as function of centre-of-mass energy. The measurement at 2.76 TeV [51] performed using the VDM method (solid triangle) is corrected for the reduced efficiency due to the VELO being not fully closed. The visible cross-section for proton-lead collisions at 5 TeV is scaled by A −2/3 . A comparison is made with the luminosity-independent measurements of the pp inelastic cross-section by the TOTEM collaboration [8,52]  The measurements from other experiments are scaled with the LHCb efficiency for inelastic events η Vertex,LHCb , obtained from simulation. The uncertainties of the direct measurements from ALICE and ATLAS are dominated by the extrapolation of the visible cross-section to the total inelastic cross-section and are not to be compared with the uncertainties of the LHCb measurements. The tick marks represent the uncertainty due to the luminosity calibration only. Data points at the same centre-of-mass energy are displaced horizontally for clarity.

Method
A weighted average of the two pp cross-section measurements at √ s = 7 TeV is computed, taking into account that the VDM calibration is more precise. The correlated part of the uncertainty due to the bunch population measurements and the factorizability is estimated to be 1%. A weight of 0.87 is given to the VDM measurement based on the uncorrelated part of the uncertainties. The final result given in Table 15 has a precision of 1.63%. If the uncertainty on the propagation to physics data (0.53%, see Sec. 4) is also included, the total uncertainty on the luminosity is 1.71%. As for the 2012 data set, the total uncertainty on the luminosity is valid if a major part of the 2011 data set is used. Measurements of the visible cross-section for the Vertex observable as function of centre-of-mass energy are shown in Fig. 57. The visible cross-section for proton-lead collisions at 5 TeV is scaled by A −2/3 (A = 208 for lead). A comparison is made with the luminosity-independent measurements of the pp inelastic cross-section by the TOTEM collaboration [8,52] at 7 and 8 TeV, and direct measurements by the ALICE [53] and the ATLAS [54] experiments. The measurements from other experiments are scaled with the LHCb efficiency for inelastic events η Vertex,LHCb , which is obtained from simulation with a negligible statistical uncertainty. The values of η Vertex,LHCb are 0.729, 0.768 and 0.758 for 2.76, 7 and 8 TeV, respectively. A measurement from CMS [55] is not shown as it does not include an extrapolation to the total inelastic cross-section. No systematic uncertainties are included to account for the η Vertex,LHCb scaling, nor for the scaling with A −2/3 . Figure 57 also shows the result of a first luminosity calibration of the pp cross-section at 2.76 TeV which was performed using the VDM method in 2011 [51]. During this data taking period, the VELO was positioned with its sensitive area at a minimum distance of 13 mm from the beam instead of the nominal 8 mm. The corresponding drop in efficiency is estimated to be 5% from simulation. Taking this difference into account, as well as the unaccounted potential non-factorizability of the beams, a good agreement is found when comparing to the more precise BGI measurement from 2013.
An earlier measurement of the pp visible cross-section based on the VDM and BGI methods for the Track observable at √ s = 7 TeV was reported in Ref. [10]. The data available in the older publication were not sufficient to measure the factorizability of the beams and complete factorizability was assumed. It was shown here that neglecting effects of non-factorizability may cause an underestimate of the visible cross-section at the few percent level. Nevertheless, the value (58.8 ± 2.0 mb) is consistent with the significantly more precise result reported here. The results reported here supersede those of Ref. [10].
In conclusion, several luminosity calibration measurements were performed at the LHC using the LHCb detector and two experimental methods, the VDM scan method and the BGI method. Ghost charge fractions were also measured using beam-gas interactions and the results are used in several luminosity calibrations by other LHC experiments. The LHCb luminosity calibrations were made for proton-proton collisions at three different centre-of-mass energies (2.76, 7 and 8 TeV) and for proton-lead collisions at the equivalent nucleon-nucleon centre-of-mass energy of 5 TeV. The analysis strategies and the results of the calibrations were presented in detail. Compared to the calibration performed in 2010, an improvement by an order of magnitude was achieved in the bunch population normalization, chiefly obtained by means of a thorough study of the LHC beam current measuring devices. This achievement opened the way to a global reduction of the systematic and statistical uncertainties for both the VDM scan method and the BGI method. A controlled gas injection into the LHC vacuum was employed to increase the beam-gas interaction rate by almost two orders of magnitude and detailed systematic studies (including offset scans and reproducibility checks) were conducted with the VDM scan method. Modelling the non-factorizability in the transverse distribution of the bunch particles is required and, if neglected, would have changed the calibration results by up to 3%. In the case of proton-proton collisions at 8 TeV, a precision of 1.47% is obtained with the VDM scan method and 1.43% with the BGI method. When combining the results, the precision obtained on the reference visible cross-section is 1.12%, which constitutes to date the most precise luminosity calibration at a bunched-beam hadron collider. The precision of the calibrations for the other beams and beam energies is close to 2%. The luminosity calibration results are used to determine a reference cross-section, which is employed in the LHCb physics data analysis to measure absolute cross-sections of various processes.