The rare semi-leptonic $B_c$ decays involving orbitally excited final mesons

The rare processes $B_c\to D_{(s)J} ^{(*)}\mu\bar{\mu}$, where $D_{(s)J}^{(*)}$ stands for the final meson $D_{s0}^*(2317)$, $D_{s1}(2460,2536)$,~$D_{s2}^*(2573)$, $D_0^*(2400)$, $D_{1}(2420,2430)$ or~$D_{2}^*(2460)$, are studied within the Standard Model. The hadronic matrix elements are evaluated in the Bethe-Salpeter approach and furthermore a discussion on the gauge-invariant condition of the annihilation hadronic currents is presented. Considering the penguin, box, annihilation, color-favored cascade and color-suppressed cascade contributions, the observables $\text{d}Br/\text{d}Q^2$, $A_{LPL}$, $A_{FB}$ and $P_L$ are calculated.

In the previous works [23,24], the process B c → D * s0 (2317)ll was calculated including only the b → sll effects, whose typical Feynman diagrams are Box and Penguin (BP) diagrams, as plotted in Figs. (1) (a, b). However, besides the BP effects, the Annihilation (Ann) diagrams, as shown in Figs. (1) (c), also make un-negligible contributions. On one hand, both BP and Ann diagrams are of order O(α em G f ) and the ratio of their CKM matrix elements is |V * cb V cs(d) |/|V * ts(d) V tb | ∼ 1. On the other hand, from Fig. (1) (c), we see that the color factors of Ann diagrams are 3 times larger than those of BP diagrams. Thus, when the decay B c → D * s0 (2317)ll is analyzed, it is necessary to include the Ann effects.
In addition to the BP and Ann effects, the process B c → D * s0 (2317)ll is also influenced by resonance cascade processes, such as B c → D *  should be experimentally removed. In Ref. [23], the regions [25], which are defined through comparing the BP and color-suppressed (CS) cascade contributions, are employed. However, in the B c → D * s0 (2317)ll process, both the color-favored (CF) and CS diagrams exist. Furthermore, the CF transition amplitudes are expected to be larger than the CS ones by a 3 times larger color factor approximately. Thus, it is necessary to redefine these regions with both CF and CS cascade influences.
So in this paper, we investigate B c → D * s0 (2317)ll transition including BP, Ann, CS and CF contributions. In addition, in order to give a more comprehensive discussion on the semi-leptonic rare decays of B c , the processes B c →D s1 (2460, 2536)ll, B c → D * s2 (2573)ll and B c → D where Q = P i −P f and P i(f ) stands for the momentum of the initial (finial) meson. V tb and V ts(d) denote the CKM matrix elements. C 10 is the Wilson coefficient. C ef f 7,9 are the combinations of the Wilson coefficients which are multiplied by the same hadronic matrix elements. The numerical value of C 10 and the explicit expressions of C ef f 7,9 can be found in Ref. [32]. The hadronic matrix elements W µ and W T µ are defined as where the definition σ µν = (i/2)[γ µ , γ ν ] is used.
Based on the effective theory [26] and the factorization hypothesis [33], the transition amplitude describing the Ann effects is [31] M Ann =V cb V * cs(d) where C 1,2 are the Wilson coefficients, whose values can be found in Ref. [32]. The annihilation hadronic current W µ ann is defined as W µ ann = W µ 1ann + W µ 2ann + W µ 3ann + W µ 4ann , where p q 1−4 and m q 1−4 are momenta and masses of the propagated quarks, respectively.
For the CS and CF cascade resonance effects, the transition amplitudes are [31] where M V and Γ V are the mass and full width of the resonance meson, respectively. Γ(V →ll) denotes the branching width of the transition V →ll. The resonance meson V stands for the particle J/ψ or ψ(2S). The CF hadronic current W µ CF is defined as Consequently, the total transition amplitude is 3

Hadronic Transition Matrix Elements in the BS Method
In Sec. 2, the transition amplitudes of the B c → D ( * ) (s)J ll processes are introduced and the hadronic matrix elements W (T ) , W ann and W CF are defined. In this section, within the BS method, we show how to calculate these hadronic matrix elements. In Sec. 3.1, we express the hadronic currents as the integrals of the wave functions. Sec. 3.2 is devoted to showing the wave functions of the mesons which are involved in this paper. Using these wave functions, we calculate the hadronic currents in Sec. 3

General Arguments on Hadronic Currents
In this part, we rewrite the hadronic currents as the integrals of the wave functions and present some general arguments.
According to the Mandelstam formalism [30], W (T ) can be expressed as the integrals of the 4-dimensional BS wave functions. In the spirit of the instantaneous approximation [34], the integrations with respect to q 0 i , where q i represents the relative momentum between the quark and anti-quark of the initial meson, can be performed first. And then we have [27,31] where the hadronic tensors Y µν V,A are defined as The term ϕ ++ i(f ) in Eqs. (8)(9) denotes the positive energy part of the initial (finial) wave function [34] and will be specified in the next subsection. In this paper we ignore the negative-energy parts since they give negligible contributions.
where q a is defined as where m b,c,s,d are masses of the constituent quarks. The parameters F i0,i± (i → f ) and Using Eqs. (10)(11)(12)(13)(14), we now discuss the gauge invariant condition of the Ann hadronic currents calculated in BS method. One may note that examining whether W ann satisfies the gauge invariant condition is equivalent to checking whether W ann · Q is zero. If we multiply Eqs. (10)(11)(12)(13) by Q µ , it is obvious that (W 1ann · Q)+(W 2ann · Q) cancels (W 4ann · Q)+(W 3ann · Q). Hence, we have W ann · Q = 0. This implies that the Ann hadronic currents in BS method indeed satisfy the gauge invariant condition. We stress that there is no need to specify the initial or final state in the process of obtaining W ann · Q = 0. Thus, our conclusion is quite general.
The decay constants of the scalar and axial-vector mesons can be found in Ref. [35].

Wave Functions in BS Method
In BS method, the meson is considered to be a bound state of two constituent quarks and can be described by the BS wave functions [28]. In the framework of instantaneous approximation [34], the time component of the BS wave functions' arguments can be integrated out and the BS equations are reduced to the Salpeter equations. By means of solving the Salpeter equations, we obtain the wave function [35][36][37][38] for each meson.
From Eq. (17), the wave functions of D s1 (2460, 2536) and D 1 (2420, 2430) can be constructed from the ones of 1 P 1 and 3 P 1 states. In the BS method, the positive energy wave functions of 1 P 1 and 3 P 1 states [39] are where ǫ A µ is the polarization vector of the axial-vector meson. The explicit expressions of b 1−4 and c 1−4 can be found in Ref. [39] and their numerical values can be obtained by solving the Salpeter equations [35]. In the processes of solving the Salpeter equations, the masses of 1 P 1 and 3 P 1 states, namely, M D (s) 1 P 1 and M D (s) 3 P 1 , are required. In analogy to the case of η 1 − η 8 mixing [44], we determine them from the following relationships [45,46], where M D 1 (2420,2430) and M D s1 (2460,2536) stand for the physical masses and we take them from Ref. [22].

(4)Wave Function of B c
The B c meson is considered as a 1 S 0 state, whose the positive energy wave function can be written as [36], where the parameters e 1−4 can be found in Ref. [36].

Calculations of Hadronic Matrix Elements
In this part, we calculate the hadronic currents through the formalism introduced above. Since W µ s have been investigated extensively in our previous papers [39,[47][48][49][50][51], here we do not introduce the W µ calculations but pay more attentions to W µ T,ann s. Please recall that W µ T s have been expressed in combinations of Y µν V,A s within Eq. (8), while in Eqs. (10)(11)(12)(13), W µ ann s are written in terms of F i,f 0(±) s. Hence, in order to obtain W µ T,ann , it is convenient to compute Y µν V,A s and F i,f 0(±) s first of all. From their definitions in Eq. (9) and Eq. (14), we see that the calculations of Y µν V,A s and F i,f 0(±) s are channel-dependent and the channels under our consideration include P → S, T, A transitions, where P, S, T, A are the abbreviations for pseudo-scalar, scalar, tensor, axial-vector mesons, respectively.

Hadronic Matrix Elements of P → S processes
First, we introduce the details of the Y µν V,A (P → S) estimations. We have expressed Y µν V,A s as the overlapping integrals of ϕ ++ i,f s in Eq. (9). In the P → S processes, the initial wave function ϕ ++ i corresponds to ϕ ++ 1 S 0 , while ϕ ++ f should be ϕ ++ 3 P 0 . The expressions of ϕ ++ 1 S 0 and ϕ ++ 3 P 0 are given in Eq. (20) and Eq. (15), respectively. Substituting Eqs. (15,20) into Eq. (9), the hadronic matrix elements Y µν V,A s can be obtained. In light of the forbidden parity, we have Y µν where the definition of q a has been given in Sec. 3.1, while q b is the relative momentum of the final meson. Due to the spectator approximation, the retarded relationship between q a and q b reads [27] Now we turn to the discussions of F i,f 0(±) (P → S)s. In Eq. (14), F i0(±) s are written in terms of ϕ ++ i s, while F f 0(±) s are shown in the integrals of ϕ ++ f s. Similar to the calculations of Y µν V,A (P → S)s, ϕ ++ i(f ) corresponds to ϕ ++ 1 S 0 ( 3 P 0 ) . So we have

Hadronic Matrix Elements of P → T processes
Here we deal with Y µν V,A in the P → T precesses. The calculations of Y µν V,A (P → T ) are similar to the ones of Y µν V,A (P → S), except replacing the final wave function ϕ ++ 3 P 0 by ϕ ++ 3 P 2 . The expression of ϕ ++ 3 P 2 can be found in Eq. (16). Hence, we have The expressions of F αµν V l and F αµν Ak , where l = 1, . . . , 7 and k = 1, 2, 3, are presented in Appendix. A.
Next, we pay attentions to F i0(±) (P → T )s. From Eq. (14), we see that F i0(±) (P → T )s are the same as F i0(±) (P → S)s, due to the identical initial meson B c in the decays P → S, T .
the situations are different. They should be calculated through Eq. (14), with the final wave After factoring the polarization tensor out, we have where

Hadronic Matrix Elements of P → A processes
Due to the mixing nature of the final mesons as formulated in Eq. (17), the calculations of And then, based on the mixing relationships in Eq. (17), we combine the results of P → A3 P 1 and P → A1 P 1 .
For Y µν V,A (P → A3 P 1 , 1 P 1 )s, we calculate them from Eq. (9), with the initial wave function For F i0(±) (P → A3 P 1 , 1 P 1 )s, we see that they are identical to F i0(±) (P → S)s. But as to F f 0(±) (P → A3 P 1 , 1 P 1 )s, we need to compute them by substituting ϕ ++ 3 P 1 , 1 P 1 into Eq. (14). The results read Finally, with the results above and the mixing relationship in Eq. (17), we can calculate the hadronic matrix elements of the physical processes from , During our calculations of Eq. (29), to avoid the kinematic confusion, we consider M f in Eqs. (27)(28) as the physical mass of the finial meson. (In this paper, the masses of 1 P 1 and 3 P 1 states introduced in Eq (19) are used only in solving the BS equations.) This approximation can also be found in the investigations of B → K 1 (1270, 1400)ll [52][53][54][55][56].

The Definitions of Form Factors
In the previous parts, we show how to calculate the hadronic currents. In order to show their results conveniently, here we parameterize the hadronic matrix elements in terms of the form factors. In this paper, we do not define the form factors of W CF s, because as introduced in Sec. 3.1, W µ CF (P → S, A) can be obtained from W µ CF (P → P, V ) by some trivial replacements, while W µ CF (P → T ) = 0. Hence, in the following paragraphs, we pay more attentions to the form factors of W (T ) and W ann s.
In the case of the P → Sll transitions, according to the Lorentz symmetry and the gauge invariant condition of the Ann currents discussed in Sec. 3.1, we have where Similarly, for P → T ll transitions, the definitions are shown as where T T zann and V T ann are the form factors.
As to P → All decays, the definitions take the following forms, where ann are the form factors.

Numerical Results of Form Factors
In this part, we present the numerical results of form factors and the according discussions.

Parameters in the Calculations
Here we specify the involved parameters. First, the masses and the lifetimes of B c and D ( * ) (s)J are required in our calculations and we take their values from Ref. [22]. Second, the BS-inputs are also needed, which include the Cornell-Potential-Parameters (CPPs) and the masses of the constituent quarks. The CPPs can be found in Ref. [57]. The masses of the constituent quarks are taken as m b = 4.96 GeV, m c = 1.62 GeV, m s = 0.5 GeV and m d = 0.311 GeV [47].

Results and Discussions on Form Factors
From the aforementioned parameters and the derivations in Sec. 3.3, the form factors can be evaluated. In the following paragraphs, we will show and discuss them. (2317)) are presented. These form factors are all positively related to Q 2 . This behavior can be understood from the facts that 1) as shown in Eqs. (8)(9), our hadronic currents W µ (T ) s are obtained from the integrals over the overlapping regions of the initial and final wave functions and 2) due to the retarded relationship in Eq. (22), the overlapping regions grow with increase in the variable Q 2 . (2317)) have also been calculated in the three-point QCD sum rules [23] and light-cone quark model [24]. The definitions of the W µ (T ) form factors in Refs. [23,24] are different from the ones in this paper. But if the same definitions are taken, the absolute values of our form factors are comparable with theirs. Fig. 3 (b) shows the form factors of W µ ann (B c → D * s0 (2317)). We see that B S z are complex. The reason is that in the calculations of the W ann , the quark propagators are involved, as shown in Eqs. (10)(11)(12)(13). In order to deal with these propagators, we separate them into two parts: the principal value terms and δ function ones. The real part of B S z comes from the principal value terms, while its imaginary part is caused by δ function terms. 2 (2573)) also increase monotonically as Q 2 grows. This similarity comes from the facts that both W µ (T ) (B c → D * s0 (2317)) and W µ (T ) (B c → D * s2 (2573)) are evaluated by Eqs. (8)(9).
In Figs. 4 (c, d), the Ann form factors of B c → D * s2 (2573)ll process are plotted. One may note that the absolute values of these form factors are quite smaller than the ones of W µ ann (B c → D * s0 (2317)). To see how this happens, one should recall that the Ann currents W ann are the sums of the terms W ann1,...,ann4 s. In the case of W µ ann (B c → D * s0 (2317)), the four terms all contribute. But as to W µ ann (B c → D * s2 (2573)), the vanishing decay constant of the final meson forbids the W ann1,ann2 contributions and leaves only W ann3,ann4 terms. Compared with the sums of W ann1 and W ann2 , the contributions of W ann3 and W ann4 are fairly suppressed. 3 Thus, we see the smaller W µ ann (B c → D * s2 (2573)) form factors in Figs. 4 (c, d).
First, we see that the form factors of W µ (T ) (B c → D s1 (2460, 2536)) are not of the same sign. To understand this feature, recall that in order to calculate W µ (T ) (B c → D s1 (2460, 2536)), the hadronic currents W (T ) (B c → D s 1 P 1 , 3 P 1 ) are first evaluated and then we mix the results according to the mixing relationship in Eq. (29). The form factors of W (T ) (B c → D s 1 P 1 , 3 P 1 ) are all of the same sign. But in the mixing step, we need to evaluate the sums and differences of the

The Calculations of Observables
In this part, we employ the helicity amplitude method [32] to calculate observables.
First of all, we need to split the total transition amplitudes as where M µ 1(2) can be determined by matching Eq. (7) to the equation above.
Finally, according to the derivations in Ref. [32], the differential branching fractions dBr/dQ 2 , the forward-backward asymmetries A F B , the longitudinal polarizations of the final mesons P L and the leptonic longitudinal polarization asymmetries A LP L can be expressed in terms of helicity amplitudes, which are where h denotes the helicity of l − , while the denotation λ = ( Plugging the helicity amplitudes H In this part, we lay stress on the introductions of numerical results of the observables. And in next section, the systematic discussions on the theoretical uncertainties will be shown. When the numerical values of observables are calculated in this paper, we have considered the BP, Ann, CS and CF diagrams. In order to show their influences clearly, for each channel, we plot 1) the observables where only BP contributions are considered, 2) the ones where BP and CS effects are contained, 3) the ones with BP and Ann influences and 4) the ones including the BP, Ann, CS and CF diagrams. In the following paragraphes, their comparisons and discussions will be presented.

The Observables of
In Figs. 11 (a, b), the differential branching fractions of B c → D * s0 (2317)µμ process are illustrated.
For dBr/dQ 2 which includes only BP contributions, as shown in the dash-dot line of Fig. 11 (a), we see that dBr/dQ 2 is biggest around Q 2 ∼ 10.5 GeV 2 and suppressed considerably at the end points. This is similar to the result in Ref. [24] but quite different from the one in Ref. [23].
If the Ann effects are added, as plotted in the dash-dot line of Fig. 11 (b), dBr/dQ 2 is enhanced un-negligibly around Q ∼ 12.5 GeV 2 .
For dBr/dQ 2 which contains BP and CS effects, as plotted in the solid line of Fig. 11 (a), because of the Breit-Wigner propagators in C CS 9 , the significant enlargements emerge around the resonance regions. If the Ann and CF diagrams are included, as displayed in solid line of Fig. 11 (b), dBr/dQ 2 around Q 2 ∼ M 2 J/ψ continues enlarging. But in light of the node structure of the ψ(2S) wave function, which leads to the cancelations in the W CF (B c → D * s0 (2317)ψ(2S) → D * s0 (2317)µμ) calculation, dBr/dQ 2 around Q 2 ∼ M 2 ψ(2S) changes imperceptibly. This feature can also be found in the processes B c → D (s) µμ [31].
In Figs. 11 (c, d), we illustrate A LP L s of the B c → D * s0 (2317)µμ process. If the Ann effects are added, as given in dash-dot line of Fig. 11. (d), A LP L deviates from −1 strongly over the low Q 2 area, while in the high Q 2 region, this kind of deviation becomes weaker.
To understand this feature, recall that the real part of Ann form factor ℜ[B S zann ] is positive within the low Q 2 domain but turns negative when Q 2 ≥ 12 GeV 2 , as shown in Fig. 3 (b). When   a, b), one may note that dBr/dQ 2 including the BP and Ann effects deviates imperceptibly from the one with only BP contribution. This is because that as plotted in Figs. 4 (c, d), the Ann form factors are quite small, which suppresses M ann considerably so that the Ann contributions are much less than the BP ones. Hence, as illustrated in Figs. 12 (a, b), dBr/dQ 2 s show the insensitivities to the Ann diagrams.

The Experimentally Excluded Regions and Integrated Branching Fractions
Using the results of dBr/dQ 2 s, as shown in Figs. 11-18 (a, b), now we define the experi- The second category contains B c → D * 2 (2460)(D * s2 (2573))µμ, B c → D s1 (2536)µμ and B c → D 1 (2420)µμ transitions, which are not sensitive to the CF contributions. So their experimentally excluded area is defined as Based on the experimentally excluded regions introduced above, the integrated branching fractions are calculated and shown in Table. 1. As seen in Table. 1, the branching fractions including BP and Ann effects are comparable with the ones containing both BP, Ann, CF and CS contributions. This implies that our experimentally excluded regions defined in Eqs. (37,38) are workable. Modes

Estimations of the Theoretical Uncertainties
In the previous section, the numerical results of the B c → D ( * ) (s)J µμ observables are discussed. In this part, we discuss their theoretical uncertainties.
In this paper, we estimate the theoretical uncertainties of the observables including two aspects. First, the theoretical errors from hadronic matrix elements are considered. Recall that our hadronic currents are calculated in the BS method and the obtained form factors are dependent on the numerical values of the BS inputs. In order to estimate the according systematic uncertainties, we calculate the observables with changing the BS inputs by ±5%.
Second, the systematic errors aroused by the factorization hypothesis are included. In the derivations of M Ann,CS,CF , the factorization hypothesis [33] is employed. In this method, in order to include the non-factorizable contributions, the number of colors N c in the expression (C 1 /N c + C 2 ) or (C 1 + C 2 /N c ) is treated as an adjustable parameter which should be determined by fitting the experimental data [58][59][60][61]. But since that the present experimental data on B c meson is still rare so that this parameter can not be obtained at the moment, we calculate the observables with N c = 3 but change the numerical values of N c within the region [2, ∞] for estimating systematic uncertainties brought by factorization hypothesis.
Actually, in recent years, several methods, dealing with the non-factorizable contributions more systematically, have been devoted to investigating the B c decays, such as perturbative QCD approach(PQCD) [62,63] and QCD factorization (QCDF) [64]. However, the channels in which the PQCD and QCDF are workable must have energetic final particles. Moreover as to B c → D ( * ) (s)J ll, the finial mesons have small recoil momenta in the high Q 2 domain. Hence, in this paper, we choose to employ the factorization method [33]. Similar situations can also be found in the calculations of B c → D ( * ) (s) ll [32,[65][66][67][68][69][70][71][72][73][74][75][76] in which the factorization method has to be used extensively to account for the non-factorizable effects.
Here we stress that using the factorization assumption to deal with the non-factorizable effects is a temporary way in the early stage of investigating the rare B c decays. A more systematical method is important and necessary. Hence, more work in the future is required.

Testing the Hadronic Matrix Elements
In the previous subsection, by changing the BS inputs within ±5%, we estimate the theoretical uncertainties from hadronic currents. Strictly speaking, this only measures parts of the uncertainties, because the systematic uncertainties from the approximations made within the BS method are not considered. Considering that this kind of uncertainties are rather difficult to be systematically estimated, in fact, we do not control the hadronic uncertainties confidently. 4 Hence, testing whether the hadronic currents are properly evaluated is important. First, we pay attentions to the decays B c → D Second, we turn to investigating B c → l AlA l BνB , whose typical diagrams are illustrated in Fig. 2. For Fig. 2 (a), the according hadronic matrix element is 0|cγ µ (1 − γ 5 )b|B c , which can be obtained from the future experimental data on pure leptonic decays B c → lν l . As to  Hence, we attempt to test them in an indirect way: we use the same framework and the same set of inputs as the ones, which are used to calculate W T , W 3ann , W 4ann and W (B c → D ( * ) sJ ), to investigate the processes B s → D * sJ µν, B → D * J µν and B c → χ cJ µν. The reasons for choosing these channels are that 1) these channels are induced by b → c(u)µν transitions, which are dominated by SM contributions from experiences of B (s) decays [22]; 2) unlike the non-leptonic decays, these semi-leptonic processes do not suffer from the theoretical uncertainties from the factorization problem. In our previous papers [50,77], the processes B s → D * sJ µν, B → D * J µν were calculated, while in Ref. [51], B c → χ cJ µν were analyzed.
In the paragraphs above, the channels B c → D In order to examine the BS inputs and the approximations of the BS method, we should pay attentions to the B c,s,u,d → D ( * ) s,d,u (η c , J/ψ)µν decays whose finial mesons are of S-wave states. In our previous papers [48,78], the observables of the processes B (s) → D ( * ) (s) µν are estimated and the results are in good agreements with the experimental observations [22]. In Ref. [79], the B c → J/ψ(η c )µν are analyzed and we expect that these channels can be tested by the future experimental data. If our results deviate from the future data, constraining our BS inputs or modifying BS method is required.
In this work, we take all the D

Conclusion
In this paper, including the BP, Ann, CS and CF contributions, we re-analyze the process  [24] but quite different from the ones in Ref. [23]. Once Ann, CS and CF Feynman diagrams are contained, the B c → D * s0 (2317)µμ observables change considerably, as shown in Figs. 11 (a-d).
2. As plotted in Figs. 14, 18 (a-h), the observables of the B c → D s1 (2536)(D 1 (2430))µμ processes behave quite sensitively to the Ann and CF influences. This makes that when these channels are analyzed, besides the BP and CS diagrams, it is necessary to include the Ann and CF ones. Here we present the explicit expressions of F α V 1−7 and F α A1−3 .