Measurement of the $B_{s}^{0} \rightarrow D_{s}^{(*)+}D_{s}^{(*)-}$ branching fractions

The branching fraction of the decay $B_{s}^{0} \rightarrow D_{s}^{(*)+}D_{s}^{(*)-}$ is measured using $pp$ collision data corresponding to an integrated luminosity of $1.0fb^{-1}$, collected using the LHCb detector at a centre-of-mass energy of $7$TeV. It is found to be \begin{align*} {\mathcal{B}}(B_{s}^{0}\rightarrow~D_{s}^{(*)+}D_{s}^{(*)-}) = (3.05 \pm 0.10 \pm 0.20 \pm 0.34)\%, \end{align*} where the uncertainties are statistical, systematic, and due to the normalisation channel, respectively. The branching fractions of the individual decays corresponding to the presence of one or two $D^{*\pm}_{s}$ are also measured. The individual branching fractions are found to be \begin{align*} {\mathcal{B}}(B_{s}^{0}\rightarrow~D_{s}^{*\pm}D_{s}^{\mp}) = (1.35 \pm 0.06 \pm 0.09 \pm 0.15)\%, \newline{\mathcal{B}}(B_{s}^{0}\rightarrow~D_{s}^{*+}D_{s}^{*-}) = (1.27 \pm 0.08 \pm 0.10 \pm 0.14)\%. \end{align*} All three results are the most precise determinations to date.


Introduction
Because of B 0 s − B 0 s oscillations, the mass and flavour eigenstates of the B 0 s system do not coincide. The B 0 s meson mass eigenstates have a relative decay width difference ∆Γ s /Γ s , where ∆Γ s (Γ s ) is the difference (average) of the decay widths between the heavy and light states. The relative decay width difference is one of the key parameters of the B The data used in the analysis presented in this paper correspond to an integrated luminosity of 1.0 fb −1 , collected by the LHCb experiment during the 2011 run period. The branching fraction of the full B 0 s → D ( * )+ s D ( * )− s decay is determined relative to the B 0 → D + s D − decay, which has a similar final state and a precisely measured branching fraction. The charm daughters are reconstructed using the D + s → K + K − π + and D − → K + π − π − final states. Throughout the paper, unless stated otherwise, charge-conjugate modes are implied and summed over. The branching fraction ratio is determined as where f d (f s ) is the fraction of B 0 (B 0 s ) mesons produced in the fragmentation of a b quark, B 0 / B 0 s is the relative efficiency of the B 0 to the B 0 s selections, B(D − → K + π − π − ) and B(D + s → K − K + π + ) are the branching fractions of the charm daughter decays, and N B 0 s /N B 0 is the relative yield of B 0 s and B 0 candidates. The branching fraction of the exclusive B 0 s → D * + s D * − s decay is determined in the same way, along with the branching fraction of B 0 s → D * ± s D ∓ s . The branching fraction of the B 0 s → D + s D − s decay has been previously measured by LHCb using the same data as this analysis [8], and is therefore not determined in this study. However, the selection efficiency and yield in the B 0 s → D + s D − s channel are determined in this analysis, as both are needed for the calculation of the total B 0 s → D ( * )+ s D ( * )− s selection efficiency.

Detector and simulation
The LHCb detector [9,10] is a single-arm forward spectrometer covering the pseudorapidity range 2 < η < 5, designed for the study of particles containing b or c quarks. The detector includes a high-precision tracking system consisting of a silicon-strip vertex detector surrounding the pp interaction region, a large-area silicon-strip detector located upstream of a dipole magnet with a bending power of about 4 Tm, and three stations of silicon-strip detectors and straw drift tubes placed downstream of the magnet. The tracking system provides a measurement of momentum, p, of charged particles with a relative uncertainty that varies from 0.5% at low momentum to 1.0% at 200 GeV/c. The minimum distance of a track to a primary vertex, the impact parameter, is measured with a resolution of (15 + 29/p T ) µm, where p T is the component of the momentum transverse to the beam, in GeV/c. Different types of charged hadrons are distinguished using information from two ring-imaging Cherenkov detectors. Photons, electrons and hadrons are identified by a calorimeter system consisting of scintillating-pad and preshower detectors, an electromagnetic calorimeter and a hadronic calorimeter. Muons are identified by a system composed of alternating layers of iron and multiwire proportional chambers. The event selection is performed in two stages, with an initial online selection followed by a tighter offline selection. The online event selection is performed by a trigger [11], which consists of a hardware stage, based on information from the calorimeter and muon systems, followed by a software stage, which performs a full event reconstruction.
In the simulation, pp collisions are generated using Pythia 6 [12] with a specific LHCb configuration [13]. Decays of hadronic particles are described by EvtGen [14], in which final-state radiation is generated using Photos [15]. The interaction of the generated particles with the detector, and its response, are implemented using the Geant4 toolkit [16] as described in Ref. [17].

Signal selection
The D * + s meson decays to a D + s meson and either a photon or a neutral pion (93.5 ± 0.7)% and (5.8 ± 0.7)% of the time, respectively, nearly saturating the total branching fraction. The remainder of the decays are ignored in this analysis. Neither of the neutral particles is reconstructed in the decay chain, and the individual B 0 s → D * ± s D ∓ s and B 0 s → D * + s D * − s decays are identified through the reconstructed invariant mass of the D + s D − s system. The individual peaks from B 0 s → D * ± s (→ D ± s γ)D ∓ s and B 0 s → D * ± s (→ D ± s π 0 )D ∓ s are not resolved. Therefore the reconstructed D + s D − s mass distribution has three separate peaks, corresponding to decays containing zero, one, or two D * ± s particles.
At the hardware trigger stage, events are required to have a muon with high p T or a hadron, photon or electron with high transverse energy in the calorimeters. For hadrons, the transverse energy threshold is 3.5 GeV. Candidate B 0 s and B 0 mesons are used in the analysis if at least one of the associated tracks is selected by the hardware trigger, or if the event is triggered independently of the particles in the signal decay. The software trigger considers all charged particles with p T > 500 MeV/c and constructs two-, three-, or four-track secondary vertices which require a significant displacement from the primary pp interaction vertices. At least one charged particle must have a transverse momentum p T > 1.7 GeV/c and be inconsistent with originating from a primary vertex. A multivariate algorithm [18] is used for the identification of secondary vertices consistent with the decay of a b hadron. The selection to this point is hereafter referred to as the initial selection.
Signal B 0 s and normalisation B 0 → D + s D − candidates are required to satisfy a number of additional conditions in order to be included in the final samples. Kaons and pions are required to be identified by the particle identification (PID) system. All D + s and D − candidates must have an invariant mass within ±30 MeV/c 2 of their known values [19]. Signal B 0 s candidates are required to have a reconstructed mass in the range 4750 − 5800 MeV/c 2 , whereas B 0 candidates must have a mass in the range 5050 − 5500 MeV/c 2 . After these requirements are applied there are still contributions from other b-hadron decays into final states with two charm particles. The decays As these backgrounds accumulate in reconstructed mass close to the signal peaks, candidates consistent with any one of these background decay hypotheses are rejected in the selection by applying a veto based on the invariant mass of the candidate under the alternative particle type hypotheses. Candidate D + s mesons are vetoed if they have a reconstructed mass in the range 2271 − 2301 MeV/c 2 when the K + candidate is assumed to be a proton, or a mass in the range 1835−1905 MeV/c 2 when the K + candidate is assigned the π + mass. Candidate D − mesons are vetoed if they have a reconstructed mass in the range 1950 − 1990 MeV/c 2 when a π − candidate is assigned the kaon mass. In a simulated sample of where h is either a kaon or pion, are examined as other potential background sources, but are all disregarded because of either a small selection efficiency or small branching fraction relative to the signal channels.
In order to further improve the purity of the signal and normalisation samples, a boosted decision tree (BDT) classifier is used to distinguish real B 0 (s) decays from combinatorial background [20]. The BDT is trained using the AdaBoost algorithm [21] to distinguish simulated B 0 s signal decays from background candidates obtained from mass sidebands in the data. Background candidates must contain a B 0 s candidate with a mass greater than 5600 MeV/c 2 and two D ± s candidates with masses less than 1930 MeV/c 2 or greater than , and for B 0 → D + s D − . Each efficiency is presented relative to the previous cut and measured using simulated events, except for the PID efficiency which is obtained from data. The D + s veto was only applied to the normalisation mode,  and N b is the total number of combinatorial background events as taken from the fit. The same BDT classifier and selection criteria are also applied to the normalisation sample.
The efficiencies of the selection criteria in both the signal and normalisation channels are listed in Table 1. The efficiencies of the background vetoes, trigger, reconstruction, and BDT selection are determined using simulated signal samples. The efficiencies of identifying K + and π + mesons are determined using a calibration data sample of D * + → D 0 (→ K − π + )π + decays, with kinematic quantities reweighted to match those of the signal candidates. decay is determined by calculating a weighted average of the individual signal channels, with weights given by the relative yields in data. The relative efficiencies of the B 0 decay to the three individual channels and the full decay are given in Table 2.

Signal and background shapes
The B 0 s and B 0 yields in the signal channels and the normalisation mode are extracted by performing a three-dimensional extended unbinned maximum likelihood fit to the mass distributions of the B 0 (s) meson and the two charm daughters. In order to determine the yields for the individual signal peaks, the B 0 s candidate mass distribution in each channel is modelled using simulated signal events. The B 0 s → D + s D − s peak is parameterised as the sum of a Crystal Ball function [22] and a Gaussian function.
The tail parameters of the Crystal Ball function, the ratio of the width of its Gaussian core to the width of the Gaussian function, and the relative weight of each function in the full distribution, are taken from simulation. The mean and width of the Gaussian core are allowed to float. The two D ± s distributions are also parameterised using this model, with all shape parameters fixed to the values found in simulation.
Because of the kinematic differences between the D * + s → D + s γ and D * + s → D + s π 0 decays, the peak of the B 0 s → D * ± s D ∓ s mass distribution is parameterised by a superposition of two Gaussian functions. The individual mean values, the ratio of the widths, and the fraction of each Gaussian function in the full distribution are fixed to values taken from simulation. The peak corresponding to B 0 s → D * + s D * − s decays is modelled using a single Gaussian function, with the mean fixed to the value found from simulated events.
There is also a component in the fit to describe the presence of background decays of the form B 0 s → D + sJ D − s , where the D + sJ can be either a D s1 (2460) + or a D s0 (2317) + meson that decays to a D + s along with some combination of photons and neutral or charged pions. As some decay products are missed, this background is present only at the low mass region of the signal distribution. The shape of the distribution is determined by fitting to B 0 s → D s1 (2460) + D − s simulated events, as the contribution from D s1 (2460) + is currently the best understood among the D + sJ decays. It is found to be well modelled by an Argus function [23], all shape parameters for which are fixed to the values found in simulation.
The combinatorial background shape in the B 0 s candidate mass distribution is parameterised by a second-order polynomial, and the model is validated with candidates passing a wrong-sign version of the selection. The wrong-sign selection is identical to the signal selection but instead looks for events containing two D + s mesons. The parameters of the combinatorial background distribution are allowed to float in the full fit to data, and are found to be compatible with the values obtained from the fit to the wrong-sign sample. The

Decay Mode Yield
combinatorial background shape in the D ± s distribution is determined using events taken from the high-mass sideband region of the B 0 s distribution, and is found to be consistent with a first-order polynomial. The impact of adding a small Gaussian contribution to account for the presence of real D ± s mesons in the combinatorial background was found to be minimal, with the observed deviations from the nominal signal yields being smaller than the statistical uncertainty in each case.
The B 0 distribution is modelled using the same parameterisation as for the full B 0 s distribution, with one exception. The peak where either the D + s or D − comes from the decay of an excited state is modelled by a superposition of three Gaussian functions, rather than the two-Gaussian model used in the B 0 s case, to account for the difference in distributions from D * + s and D * − decays, as the D * − decay contains a π 0 in the final state more frequently than D * + s decays. There is also a small contribution from the decay B 0 s → D − s D + , which is modelled with the same distribution as for the signal B 0 s candidates.

Fit results
The fit to the signal data samples is shown in Fig. 1, where the triple peaked structure of the full decay is clearly visible. The yields for the individual signal channels and the two backgrounds are given in Table 3. The total B 0 s → D ( * )+ s D ( * )− s yield is the sum of the individual signal channel yields, with the uncertainty calculated using the correlation coefficients between the individual yields, and is found to be 2230 ± 63. The full fit to the data sample for the normalisation mode is shown in Fig. 2, and the yields are given in Table 4. Almost all B 0 → D * ± s D * ∓ decays are reconstructed with a mass lower than the 5050 MeV/c 2 mass cut imposed on the B 0 candidates. There is thus a relatively small yield from this channel. Only the main B 0 → D + s D − peak is used for normalisation purposes.

Systematic uncertainties
A number of systematic uncertainties affect the measurements of the ratios of branching fractions; the sources and magnitudes of these uncertainties are summarised in Table 5.

Decay Mode Yield
Combinatorial background 1542 ± 56 b fragmentation fraction ratio, f s /f d = 0.259 ± 0.015 [24]. Part of the uncertainty on this ratio is due to the ratio of the charm branching fractions B(D + s → K + K − π + )/B(D − → K + π − π − ) = 0.594 ± 0.020 [24], the inverse of which is used in the measurements presented in this paper, as shown in Eq. 1. With the two values from Ref. [24], the part of the uncertainty on f s /f d due to the charm branching fractions cancels, leading to a total uncertainty for the product 7%. The fit model used for the yield extraction is validated using pseudoexperiments and is found to be unbiased. The uncertainty due to the imperfect knowledge of the shape of the full mass distribution is investigated by measuring the yields using alternative models for each of the peaks. The D + s D − s peak is modelled with an Apollonios function [25] or a Cruijff function [26], the D * ± s D ∓ s peak is modelled using a single Gaussian function, and the D * + s D * − s peak is modelled using a combination of two Gaussian functions. The  branching fraction ratio. The uncertainty on the combinatorial background yield is determined by considering the differences when instead fitting this background with an exponential function, and is of the order of 1.5% for all of the branching fraction ratios.
The dominant uncertainty for the B 0 s → D * + s D * − s decay channel results from the lack of knowledge of the B 0 s → D + sJ D − s background decays. The shape of this background overlaps mostly with the B 0 s → D * + s D * − s signal decays, and therefore the systematic uncertainty due to this background shape is much larger for this channel (5.0%) than for the other two exclusive branching fractions (0.2% − 0.4%). The uncertainty is measured by repeating the fit with the cut-off point of the Argus function varied from 5050 MeV/c 2 to 5200 MeV/c 2 , where the upper limit is chosen in order to account for the presence of decays containing D s0 (2317) + mesons. The changes to the yields from the values found in the nominal fit are calculated in each case. The systematic uncertainty in each channel is then assigned as the RMS of the full set of deviations. The uncertainty on the overall branching fraction ratio is also determined in this way, and is found to be 1.9%.
The uncertainties on the overall efficiencies due to the limited size of the simulated   There is a systematic uncertainty arising from the calculation of the efficiencies of the PID cuts. The calibration of the data samples is performed in bins of momentum and pseudorapidity, which results in an uncertainty on the calculated efficiencies owing to the finite size of the D * + → D 0 π + calibration samples and the binning scheme used. The uncertainty resulting from the calibration sample size and binning scheme is determined by redoing the calibration using different binning schemes. Another systematic uncertainty is due to the presence of a small combinatorial background component in the samples that are used to determine the PID efficiencies. The systematic uncertainty due to this contamination is estimated by comparing the efficiencies found in data to those found when calibrating simulated signal events. The total uncertainties due to the PID efficiency calculation for the three branching fraction ratios presented in this paper are shown in Table 5. The value for B 0 s → D ( * )+ s D ( * )− s is again the weighted average of the contributing channels, with the uncertainty for the D + s D − s contribution being 1.1%. The uncertainty of 1.5% from the trigger response is assessed by considering variations in the response between data and simulation. The individual uncertainties are combined in quadrature to give the total relative systematic uncertainties for each measurement given in Table 5.
Using the current world average measurement of the B 0 → D + s D − branching fraction of (7.2 ± 0.8) × 10 −3 [19], gives where the uncertainties are statistical, systematic, and due to the branching fraction of the normalisation channel, respectively. Figure 3 shows the LHCb measurement of the total B 0 s → D ( * )+ s D ( * )− s branching fraction, along with the previous measurements by Belle [5], CDF [6], and D0 [7], the average of these previous measurements as calculated by HFAG [27], and the theoretical value from Ref. [3]. The theoretical prediction is for a decay time t = 0, while the measurements integrate over the B 0 s meson lifetime; the correction factor for mixing is known [28], but has not been applied. The LHCb result is consistent with all previous measurements and calculations, and is the most precise determination to date. In addition, the LHCb measurements of the individual B 0 s → D * ± s D ∓ s and B 0 s → D * + s D * − s branching fractions are consistent with, and more precise than, all previous measurements.
Using this measurement of the branching fraction of B 0 s → D [2] A. Lenz, Theoretical update of B-mixing and lifetimes, arXiv:1205.1444.