K to ll' revision (fix to #138)

flavio/physics/kdecays/kll.py

@@ -9,15 +9,13 @@
 from math import pi, sqrt
 
 
-def amplitudes(par, wc, l1, l2):
-    r"""Amplitudes P and S entering the $K\to\ell_1^+\ell_2^-$ observables.
+def amplitudes_weak_eigst(par, wc, l1, l2):
+    r"""Amplitudes P and S for the decay $\bar K^0\to\ell_1^+\ell_2^-$.
 
     Parameters
     ----------
-
     - `par`: parameter dictionary
     - `wc`: Wilson coefficient dictionary
-    - `K`: should be `'KL'` or `'KS'`
     - `l1` and `l2`: should be `'e'` or `'mu'`
     """
     # masses
@@ -33,38 +31,73 @@ def amplitudes(par, wc, l1, l2):
     CSm = wc['CS_'+qqll] - wc['CSp_'+qqll]
     P = (ml2 + ml1)/mK * C10m + mK * CPm  # neglecting mu, md
     S = (ml2 - ml1)/mK * C9m  + mK * CSm  # neglecting mu, md
+    # Include complex part of the eff. operator prefactor. Phases matter.
     xi_t = ckm.xi('t', 'sd')(par)
     return xi_t * P, xi_t * S
 
 
+def amplitudes(par, wc, K, l1, l2):
+    r"""Amplitudes P and S entering the $K_{L,S}\to\ell_1^+\ell_2^-$ observables.
+
+    Parameters
+    ----------
+
+    - `par`: parameter dictionary
+    - `wc`: Wilson coefficient dictionary
+    - `K`: should be `'KL'` or `'KS'`
+    - `l1` and `l2`: should be `'e'` or `'mu'`
+    """
+    # KL, KS are linear combinations of K0, K0bar. So are the amplitudes.
+    S_K0bar, P_K0bar = amplitudes_weak_eigst(par, wc, l1, l2)
+    if l1 != l2:
+        S_aux, P_aux = amplitudes_weak_eigst(par, wc, l2, l1)
+        S_K0 = -S_aux.conjugate()
+        P_K0 = P_aux.conjugate()
+        if K == 'KL':
+            sig = +1
+        elif K == 'KS':
+            sig = -1
+        S = (S_K0 + sig * S_K0bar) / sqrt(2)
+        P = (P_K0 + sig * P_K0bar) / sqrt(2)
+    # Save computing time for special cases. See also arXiv:1711.11030.
+    elif l1 == l2:
+        if K == 'KL':
+            S = -1j * S_K0bar.imag * sqrt(2)
+            P = P_K0bar.real * sqrt(2)
+        elif K == 'KS':
+            S = -S_K0bar.real * sqrt(2)
+            P = -1j * P_K0bar.imag * sqrt(2)
+    return P, S
+
+
 def amplitudes_LD(par, K, l):
     r"""Long-distance amplitudes entering the $K\to\ell^+\ell^-$ observables."""
     ml = par['m_' + l]
-    mK = par['m_' + K]
+    mK = par['m_K0']
     s2w = par['s2w']
-    pre = 2 * ml / mK / s2w
+    pre = sqrt(2) * 2 * ml / mK / s2w
     # numbers extracted from arXiv:1711.11030
-    ASgaga = 2.49e-4 * (-2.821 + 1.216j)
-    ALgaga = 2.02e-4 * (par['chi_disp(KL->gammagamma)'] - 5.21j)
-    S = pre * ASgaga
-    P = -pre * ALgaga
-    return S, P
+    if K == 'KS':
+        ASgaga = 2.49e-4 * (-2.821 + 1.216j)
+        SLD = +pre * ASgaga
+        PLD = 0
+    elif K == 'KL':
+        ALgaga = 2.02e-4 * (par['chi_disp(KL->gammagamma)'] - 5.21j)
+        SLD = 0
+        PLD = -pre * ALgaga
+    return SLD, PLD
 
 
 def amplitudes_eff(par, wc, K, l1, l2, ld=True):
     r"""Effective amplitudes entering the $K\to\ell_1^+\ell_2^-$ observables."""
-    P, S = amplitudes(par, wc, l1, l2)
+    P, S = amplitudes(par, wc, K, l1, l2)
     if l1 != l2 or not ld:
         SLD = 0
         PLD = 0
     else:
         SLD, PLD = amplitudes_LD(par, K, l1)
-    if K == 'KS' and l1 == l2:
-        Peff = P.imag
-        Seff = S.real + SLD
-    if K == 'KL':
-        Peff = P.real + PLD
-        Seff = S.imag
+    Peff = P + PLD
+    Seff = S + SLD
     return Peff, Seff
 
 
@@ -95,12 +128,12 @@ def br_kll(par, wc_obj, K, l1, l2, ld=True):
     mK = par['m_K0']
     tauK = par['tau_'+K]
     fK = par['f_K0']
-    # appropriate CKM elements
+    # CKM part of the eff. operator prefactor included in Peff and Seff
     N = 4 * GF / sqrt(2) * alphaem / (4 * pi)
     beta = sqrt(lambda_K(mK**2, ml1**2, ml2**2)) / mK**2
     beta_p = sqrt(1 - (ml1 + ml2)**2 / mK**2)
     beta_m = sqrt(1 - (ml1 - ml2)**2 / mK**2)
-    prefactor = 2 * abs(N)**2 / 32. / pi * mK**3 * tauK * beta * fK**2
+    prefactor = abs(N)**2 / 32. / pi * mK**3 * tauK * beta * fK**2
     Peff, Seff = amplitudes_eff(par, wc, K, l1, l2, ld=ld)
     return prefactor * (beta_m**2 * abs(Peff)**2 + beta_p**2 * abs(Seff)**2)