Next-to-leading order QCD corrections to paired $B_c$ production in $e^+e^-$ annihilation

We present theoretical analysis of paired $B_c$ mesons production in $e^+e^-$ annihilation at different energy scales taking into account full next-to-leading order QCD corrections. Both possible electroweak channels are considered: production via virtual photon and via virtual $Z$-boson. We study in detail the role of radiative QCD corrections, which were found to be significant especially at low energies. It is shown that the contribution from $Z$-boson is significant at high energies ($\sqrt{s}>M_Z/2$) especially in the case of paired production of pseudo-scalar and vector ($B_c + B_c^*$) mesons. Azimuthal asymmetry induced by a $P$-violating weak interaction with $Z$-boson is also analyzed.


Introduction
In the recent years bound systems of heavy quarks became a perfect laboratory for rigorous theoretical and experimental studies of Standard Model phenomena. Due to the large mass of quarks, production of quark-antiquark pairs can be described by the perturbation theory, while hadronization into observed mesons can be accounted for with the use of effective non-relativistic theories such as NRQCD [1] and pNRQCD [2][3][4]. Among existing heavy quarkonia systems, B cmesons (bc) are still the most poorly studied experimentally compared to for example charmed or beauty mesons. This is due to the fact that the production cross section is small compared to (cc) or (bb) systems since two pairs of heavy quark pairs must be produced simultaniously. The first observation of B c (0 − ) was done in hadronic production by CDF collaboration [5] and then it was extensively studied at Tevatron and LHC. This pseudoscalar 0 − state is the only one experimentally confirmed (bc) state (recently ATLAS has found evidence for its radial excitation B c (2S ) [6]). The most significant efforts in the experimental study of B c mesons nowadays are undertaking by the LHCb collaboration, which has studied about a dozen of decay channels of B c meson (see [7] and references therein) as well as its hadronic production [8] which they found in a good agreement with theory [9][10][11][12][13].
In the present work we study paired production of S -wave (bc)-mesons through virtual γ * and Z * -boson exchanges in e + e − -annihilation taking into account full next-to-leading order (NLO) QCD radiative corrections. We show that these one-loop corrections give sensitive contribution to the whole production cross sections. Additionally, theoretical errors induced by the choice of the renormalization scale in strong coupling α S (µ) become several times smaller at NLO compared to leading order (LO). In the case of the production of a pair of pseudoscalar mesons there is an interesting feature related to the cross section becoming zero (both within LO and NLO) at some value of center of mass energy near threshold, which can be potentially interesting for experimental search and study of B * c mesons, since all (bc)-mesons produced exclusively at this point are 1 − . We show that the account for virtual Z * -boson exchange changes significantly the behavior of B c + B * c production cross section, in particular it changes the asymptotic behavior at large energies and completely changes the cross section dependence on the renormalization scale µ. Such effects are absent in the case of B c + B c and B * c + B * c production. Parity violation induced by weak interaction is shown to be significant at the energies near m Z /2 for the production cross section of B c + B * c and B * c + B * c , while in the case of B c + B c there is no parity violation in all energy range since final state with such quantum numbers can not be produced from the axial part of Z-current and thus there is no V-A interference.

The method
To calculate paired B ( * ) c production cross sections we employed a systematic framework offered by NRQCD factorization [1]. This way cross sections are being expressed as the sum over bc and cb channels of products of hard perturbative bc, cb cross sections and soft nonperturbative NRQCD matrix elements.

NRQCD factorization and projectors
Production and bound state dynamics of quarkonia involves several different energy scales. First, we have hard scales given by heavy quark masses m b and m c . Next, there are so called soft and ultrasoft scales. The soft scale is given by the relative momentum of heavy quark-antiquark pair in the quarkonium rest frame |p| ∼ mv (∼ 1/r). Here v is the relative velocity of quarkantiquark motion, r is the radius of quarkonium and m is the reduced mass, which in the case of B c meson is given by m = m c m b m c +m b . The ultrasoft scale is given by kinetic energy of quarkantiquark motion inside quarkonium E ∼ mv 2 and of course we have nonperturbative confinement scale Λ QCD . The nonrelativistic QCD (NRQCD) is the effective theory obtained by integrating out hard degrees of freedom from QCQ. The general hierarchy of scales for quarkonia production and bound state dynamics is m b , m c mv, mv 2 , Λ QCD , as v 1. First, it sets a firm ground for constructing a self consistent effective theory for bound state dynamics and, second, it allows us to write down production cross sections for quarkonia in a factorized form: hard perturbative production of heavy quark pairs and subsequent soft nonperturbative formation of quarkonia bound states out of these heavy quark-antiquark pairs. To calculate matrix elements for e + e − → B +( * ) c (P 1 ) + B −( * ) c (P 2 ) production processes, we start from the matrix element for e + e − → c(p c )b(pb) + b(p b )c(pc) with heavy quarks and antiquarks on their mass shells: For each of the B c mesons the momenta of heavy quarks are written as

LO and NLO cross section calculation
So, to calculate production cross sections of B ( * ) c pairs we need to know expressions for hard production cross sections of two heavy quark-antiquark pairs projected on states with given quantum numbers. This problem requires automation already at tree level, where we have 8 diagrams, see Fig. 1. To generate Feynman diagrams and obtain analytical expressions for them we used FeynArts package [47]. The subsequent processing, such as color and Dirac algebra were performed using several different strategies based on usage of different computer algebra systems such as Mathematica with FeynCalc [48,49] package , FORM [50] or Redberry [51]. This way, for example, for A V P µ µ 1 matrix element accounting only for virtual photon exchange we get where k 1 , k 2 are incoming momenta of electron-positron pair, s = (k 1 + k 2 ) 2 , r = m c m c +m b and C F = To get an expression for cross section this matrix element should be squared, summed over B * c meson polarizations, averaged over polarizations of initial electron-positron pair and integrated over final phase space. In Appendix A we have gathered the results for leading order differential B ( * ) c + B ( * ) c cross sections including contributions both from virtual γ * and Z * -boson exchanges. It should be noted, that the results for the case of γ * contribution have already appeared in literature [22], while account for Z * -boson exchange is a new result.
The NLO calculation is more involved. At this order we have 172 nonzero diagrams, see Fig. 1 for some sample diagrams. To regularize both UV and IR singularities we used dimensional regularization (D = 4 − 2 ). In order to calculate hard perturbative production cross section we used so called method of regions [52]. This way we may directly deduce NRQCD short-distance coefficients (hard production cross sections in our case) by setting relative momenta of heavy quark-antiquark motion inside quarkonia to zero before carrying out loop integrations rather than resorting to much more expensive matching calculations. The diagram calculation itself was again performed using several computational strategies. To reduce diagrams to master integrals we used Laporta algorithm [53] and its implementation in FIRE [54,55]. Before this step however one should resolve a possible linear dependence of propagators and perform partial fractioning when needed. For that purpose we used $Apart [56] package. The tensor integrals were treated using Passarino-Veltman reduction [57]. Another strategy employed our own private Mathematica code, which performed topological analysis of the generated diagrams and produced FORM rules for diagram reduction to canonical prototypes involving numerator simplification, partial fractioning and Passarino-Veltman reduction together with subsequent reduction of canonical prototypes to master integrals using reduction database generated with FIRE. For numerical evaluation of master integrals we used PackageX [58]. In the case of two-point A 0 and B 0 master integrals it was required to hold terms up to 1 order, so we used known analytical results which are given for reader's convenience in Appendix B. For renormalization procedure we employed so called on-shell scheme 1 : , m is the heavy quark pole mass, n f -number of fermions taken into account in gluon self-energies and γ E -Euler's gamma constant. The resulting expressions for renormalized one-loop cross sections in terms of master integrals could be found in supplementary material. For γ 5 treatment we used two different prescriptions, namely West [59] and Larin [60] schemes . It turns out, that for the production of a pair of pseudoscalar (PP) mesons and production of vector and pseudoscalar meson (VP) both schemes give identical results. At the same time in the case of the production of a pair of vector mesons (VV) both schemes gave slightly (less then a percent) different results. In this case we regard this discrepancy as an error related to γ 5 scheme dependence. See [61] for discussion of problems present when dealing with γ 5 in dimensional regularization. It is also worth to mention that if we consider only γ * contribution both schemes give identical results for all final states.

Results
To present numerical estimates we used the following values of parameters: The running of strong coupling constant was taken at one-loop accuracy and is given by , 1 The gauge coupling renormalization was performed in the modified minimal subtraction MS scheme.
with α S (m Z ) = 0.1184. The renormalization scale µ R and strong coupling scale Q were chosen as 2 µ R = Q = µ = √ s. Fig. 2 shows tree level cross sections near production threshold (left) and near Z-boson peak (right). First of all, from the left plot we see that the production cross section of a pair of pseudoscalar (PP) mesons e + + e − → B + c + B − c is zero at some center of mass energy near threshold. This point is where r = m c /M and M = m c + m b . If we neglect Z-boson contribution then the cross section at this point is exactly zero. Accounting for Z-boson it is given by: which is negligibly small (as it is generally the case for Z-boson contribution at small energies). At the same time, production cross sections of GeV (moreover σ VV /σ PV ≈ 4 at this energy). So, one can use this energy point as a some kind of "B * c -factory" since all (bc)-mesons produced exclusively at this point are 1 − . While there is no yet any experimental evidence for B * c , such setup can be potentially interesting for experimental search and study of B * c mesons. At high energies the cross sections for all considered processes tend to zero as 1/s 3 . However, without the account of Z-boson contribution the production cross section of pseudo-scalar plus vector (PV) mesons goes to zero as 1/s 4 . This can be explained by the fact that PV production is forbidden in massless limit due to chirality conservation (each quark-antiquark pair is produced from either photon or gluon which both are 1 − and there is no way to combine them into 1 − + 0 − state combination). In the case of massive quarks it is additionally suppressed by M/ √ s factor giving extra power of s to the asymptotic behavior of cross section. Accounting for Z-boson contribution containing axial current this restriction is removed and σ PV → 1/s 3 at high energies. In the case of a pair of pseudoscalar (PP) and a pair of vector (VV) mesons an account of Z-boson contribution raises cross section at high energies by a factor: important at low energies and may even raise the cross sections by a factor of two as it is the case for PV mesons production. At high energies σ NLO /σ LO tends to ≈ 1.2. Next, from Fig. 3 it is seen that in the PP case loop contribution is negative at low energies (square of one-loop contribution is discarded at O(α 3 S )). This leads to unphysical result near s PP point defined before in (6), where tree cross section goes to zero. So, to correctly describe production cross section at this point at least two-loop calculation is required as well as account for relativistic O(v 2 ) corrections (such corrections to paired B c -production were recently studied in [23]).
To investigate the dependence of full NLO cross section on the choice of renormalization scale and to deduce a range of admissible values we employed the following strategy. We defined the measure of importance of O(α 3 S ) corrections as: In order to have perturbation series "convergent" 3 R should be small in all energy range. So, we first adopted the constraint R = z for a set of fixed z values in the range [0.1 ÷ 0.5] and solved the corresponding equation with respect to µ. Fig. 5 (left) shows obtained renormalization scales µ as a function of √ s in the case of production of a pair of pseudoscalar mesons (results for other states are the same and we omitted them in the figure). It is seen, that if one takes µ proportional to √ s, R will be constant in all energy range. In particular, taking µ = √ s we have R = 0.3, while varying µ in the range √ s/2 ≤ µ ≤ 2 √ s, R varies in the range 0.2 ≤ R ≤ 0.5. This result is quite predictable, since it is common to take renormalization scale equal to some characteristic energy of the process and in the case of the considered processes √ s appears in all propagators. It is interesting to note, that in the case of production of pseudoscalar plus vector (PV) mesons without account for virtual Z * boson the situation is completely different. Fig. 5 (right) shows, that in this case the renormalization scale, which ensures smallness of R value at all energies, is inversely proportional to √ s. This is probably related to the fact that in this case there is an additional M/ √ s suppression in the amplitude as discussed previously. Varying µ value in the range √ s/2 ≤ µ ≤ 2 √ s the NLO contribution for the production of different states varies as  Fig. 6 (right) shows this dependence for LO and NLO cross sections. It is seen, that the variation of µ leads to 35% uncertainty in the LO cross section and only to 10% uncertainty for the NLO cross section. So, as expected an account for NLO corrections stabilizes the dependence of cross section on µ and provides higher precision of theoretical predictions. Finally, Tab. 1 presents cross sections obtained at different energies and Fig. 7 shows results for NLO cross sections with the uncertainties induced by the variation of µ.
Weak interactions are known to be P-violating and thus we should observe this phenomenon in the considered processes. To study P-asymmetry we choose the following quantity: which is of the order of the difference (σ(θ = π) − σ(θ = −π)). In the case of PP production this asymmetry is exactly zero in all energy range. This is easily explained by the fact that axial part of Z-current which is 1 + can not produce two 0 − states and thus axial current gives no contribution in this case and there is no V-A interference and parity violation. In the case of PV and VV mesons production the dependence of A on the center of mass energy is shown in Fig. 8. It is seen, that the maximal asymmetry should be observed at √ s ≈ m Z /2 and it is almost zero at Z-boson peak. Specifically, the asymmetry becomes zero at the following points in the case of PV and VV mesons  Figure 6: Scale dependence of the NLO to LO ratio R defined in (7) (left) and the total cross section (right). Both plots are evaluated for VV final state at √ s = 20 GeV (results for PP and PV are almost the same and thus they were omitted). Filled areas correspond to interval √ s/2 ≤ µ ≤ 2 √ s. Two curves on the right plot correspond to leading order (dashed curve) and full NLO cross section (solid curve).    production:

Conclusions
In the present paper we have considered paired production of different S -wave B c states in e + e − -annihilation via virtual γ * and Z * bosons taking into account full NLO QCD corrections. We have found that one-loop corrections are quite sizable and may raise the LO cross section up to 2 times at low energies, while at high energies σ NLO /σ LO ≈ 1.2. The dependence on the renormalization scale µ stabilizes with the account of NLO corrections, so that errors induced by variation of the scale are reduced from 30% to only 10% at NLO compared to LO.
An interesting feature of the production cross section for PP state is that it becomes zero at some energy point (6)  exchange changes asymptotic behavior of the cross section at high energies from 1/s 4 to 1/s 3additional 1/s suppression in the case of vector current only can be explained by approximate chirality conservation. Another interesting feature related with this is that account for Z-boson completely changes PV cross section dependence on the renormalization scale µ.
Finally, parity asymmetry induced by weak current appears to be most significant near √ s = m Z /2 point, while it is almost zero at √ s = m Z . It is interesting to note, that in the case of PP mesons production there is no parity asymmetry since a pair of 0 − states can not be produced by axial current thus there is no V-A interference in the cross section.