Module flavio.physics.bdecays.bpll_lfv

Functions for exclusive $B\to P\ell^+\ell^-$ decays.

Functions

def bpll_dbrdq2(q2, wc, par, B, P, l1, l2)

Expand source code

def bpll_dbrdq2(q2, wc, par, B, P, l1, l2):
    tauB = par['tau_'+B]
    dBR = tauB * bpll_obs(dGdq2, q2, wc, par, B, P, l1, l2)
    if P == 'pi0':
        # factor of 1/2 for neutral pi due to pi = (uubar-ddbar)/sqrt(2)
        return dBR / 2.
    else:
        return dBR

def bpll_dbrdq2_int(q2min, q2max, wc_obj, par, B, P, l1, l2, epsrel=0.005)

Expand source code

def bpll_dbrdq2_int(q2min, q2max, wc_obj, par, B, P, l1, l2, epsrel=0.005):
    scale = config['renormalization scale']['bpll']
    label = meson_quark[(B,P)] + l1 + l2 # e.g. bsmumu, bdtautau
    wc = wc_obj.get_wc(label, scale, par)
    if all([abs(v) < 1e-12 for v in wc.values()]):
        # if all WCs are essentially zero, return BR=0
        return 0
    def obs(q2):
        return bpll_dbrdq2(q2, wc, par, B, P, l1, l2)
    return flavio.math.integrate.nintegrate(obs, q2min, q2max, epsrel=epsrel)/(q2max-q2min)

def bpll_dbrdq2_tot_func(B, P, l1, l2)

Expand source code

def bpll_dbrdq2_tot_func(B, P, l1, l2):
    def fct(wc_obj, par):
        mB = par['m_'+B]
        mP = par['m_'+P]
        ml1 = par['m_'+l1]
        ml2 = par['m_'+l2]
        q2max = (mB-mP)**2
        q2min = (ml1+ml2)**2
        return bpll_dbrdq2_int(q2min, q2max, wc_obj, par, B, P, l1, l2)*(q2max-q2min)
    return fct

def bpll_dbrdq2_tot_lfv_comb_func(B, P, l1, l2)

Expand source code

def bpll_dbrdq2_tot_lfv_comb_func(B, P, l1, l2):
    def fct(wc_obj, par):
        mB = par['m_'+B]
        mP = par['m_'+P]
        ml1 = par['m_'+l1]
        ml2 = par['m_'+l2]
        q2max = (mB-mP)**2
        q2min = (ml1+ml2)**2
        return (
            + bpll_dbrdq2_int(q2min, q2max, wc_obj, par, B, P, l1, l2)
            + bpll_dbrdq2_int(q2min, q2max, wc_obj, par, B, P, l2, l1)
        )*(q2max-q2min)
    return fct

def bpll_obs(function, q2, wc, par, B, P, l1, l2)

Expand source code

def bpll_obs(function, q2, wc, par, B, P, l1, l2):
    ml1 = par['m_'+l1]
    ml2 = par['m_'+l2]
    mB = par['m_'+B]
    mP = par['m_'+P]
    if q2 <= (ml1+ml2)**2 or q2 > (mB-mP)**2:
        return 0
    scale = config['renormalization scale']['bpll']
    mb = running.get_mb(par, scale)
    h     = helicity_amps(q2, wc, par, B, P, l1, l2)
    J     = angular.angularcoeffs_general_p(h, q2, mB, mP, mb, 0, ml1, ml2)
    return function(J)

def dGdq2(J)

Expand source code

def dGdq2(J):
    return 2 * (J['a'] + J['c']/3.)

def helicity_amps(q2, wc, par_dict, B, P, l1, l2)

Expand source code

def helicity_amps(q2, wc, par_dict, B, P, l1, l2):
    par = par_dict.copy()
    scale = config['renormalization scale']['bpll']
    wc_eff = get_wceff_lfv(q2, wc, par, B, P, l1, l2, scale)
    ml1 = par['m_'+l1]
    ml2 = par['m_'+l2]
    mB = par['m_'+B]
    mP = par['m_'+P]
    mb = running.get_mb(par, scale)
    N = prefactor(q2, par, B, P)
    ff = get_ff(q2, par, B, P)
    h = angular.helicity_amps_p(q2, mB, mP, mb, 0, ml1, ml2, ff, wc_eff, N)
    return h