Merge c2e173c into 011236d

flavio/__init__.py

@@ -7,6 +7,7 @@
 from . import classes
 from .classes import Measurement, Parameter, ParameterConstraints, Observable, NamedInstanceClass
 from .config import config
+from .citations import Citations, default_citations, register_citation
 from flavio.physics.eft import WilsonCoefficients
 from flavio.parameters import default_parameters
 from flavio.functions import sm_prediction, sm_uncertainty, np_uncertainty, sm_error_budget, np_prediction, sm_covariance, combine_measurements

flavio/citations.py

@@ -0,0 +1,88 @@
+"""Citation class for handling theory calculation citations"""
+# citations.py is based on code from the PyBaMM project
+
+# Copyright (c) 2018-2020, the PyBaMM team.
+# All rights reserved.
+
+# Redistribution and use in source and binary forms, with or without
+# modification, are permitted provided that the following conditions are met:
+
+# * Redistributions of source code must retain the above copyright notice, this
+#   list of conditions and the following disclaimer.
+
+# * Redistributions in binary form must reproduce the above copyright notice,
+#   this list of conditions and the following disclaimer in the documentation
+#   and/or other materials provided with the distribution.
+
+# * Neither the name of the copyright holder nor the names of its
+#   contributors may be used to endorse or promote products derived from
+#   this software without specific prior written permission.
+
+# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
+# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
+# DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
+# FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+# DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+# SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
+# CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
+# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
+# OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+
+import flavio
+
+
+class Citations:
+
+    """Entry point to citations management.
+    This object may be used to register a citation is relevant for a particular
+    implementation. It can then also then print out a list of papers to be cited.
+
+    Examples
+    --------
+    >>> import flavio
+    >>> flavio.sm_prediction("DeltaM_s")
+    >>> print(flavio.default_citations)
+    """
+
+    def __init__(self):
+        self._reset()
+
+    def __iter__(self):
+        for citation in self._papers_to_cite:
+            yield citation
+
+    def __str__(self):
+        return ",".join(self._papers_to_cite)
+
+    def _reset(self):
+        "Reset citations to empty"
+        self._papers_to_cite = set()
+
+    @property
+    def tostring(self):
+        return str(self)
+
+    @property
+    def toset(self):
+        return self._papers_to_cite
+
+    def register(self, key):
+        """Register a paper to be cited. The intended use is that this method
+        should be called only when the referenced functionality is actually being used.
+
+        Parameters
+        ----------
+        key : str
+            The INSPIRE texkey for the paper to be cited
+        """
+        self._papers_to_cite.add(key)
+
+
+
+default_citations = Citations()
+# Register the flavio paper
+default_citations.register("Straub:2018kue")
+
+def register_citation(inspire_key):
+    flavio.default_citations.register(inspire_key)

flavio/classes.py

@@ -697,6 +697,20 @@ def get_measurements(self):
                             break
         return ms
 
+    def theory_citations(self, *args, **kwargs):
+        """Return a set of theory papers (in the form of INSPIRE texkeys) to
+        cite for the theory prediction for an observable.
+        
+        Arguments are passed to the observable and are necessary,
+        depending on the observable (e.g. $q^2$-dependent observables).
+        """
+        old_citations = flavio.default_citations
+        flavio.default_citations = flavio.Citations()
+        flavio.sm_prediction(self.name, *args, **kwargs)
+        SM_citations = flavio.default_citations.toset
+        flavio.default_citations = old_citations
+        return SM_citations
+
 
 
 class AuxiliaryQuantity(NamedInstanceClass):

flavio/physics/bdecays/bllgamma.py

@@ -52,7 +52,8 @@ def getfft(s, par, B, ff, ff0, lep, wc):
 
     #Add light meson resonances
     resonances = {'Bs': ['phi'], 'B0': ['rho0', 'omega']}
-    fVtofVemfactors = {'phi': -1/3, 'rho0': 1/sqrt(2), 'omega': 1/(3*sqrt(2))}
+    flavio.register_citation("Kozachuk:2017mdk")
+    fVtofVemfactors = {'phi': -1/3, 'rho0': 1/sqrt(2), 'omega': 1/(3*sqrt(2))} # Given in Sec. 8.A.3 of 1712.07926
     #We use the general B->rho and B->omega parameters, hence the isotopic factors
     resgV = {('Bs','phi'): -par['Bs->phi BSZ a0_T1'],
              ('B0','rho0'): 1/sqrt(2)*par['B->rho BSZ a0_T1'],
@@ -66,6 +67,7 @@ def getfft(s, par, B, ff, ff0, lep, wc):
     return (fftv, ffta)
 
 def getF1(s, par, B, ff, ff0, lep, wc):
+    flavio.register_citation("Melikhov:2004mk")
     scale = config['renormalization scale']['bllgamma']
     mb = running.get_mb(par, scale, nf_out=5)
     mbh = mb/par['m_'+B]
@@ -75,6 +77,7 @@ def getF1(s, par, B, ff, ff0, lep, wc):
     return (abs(wc['C9_'+label])**2 + abs(wc['C10_'+label])**2)*ff['v']**2 + 4*mbh**2/s**2*abs(wc['C7_'+bq]*fftv)**2 + 4*mbh/s*ff['v']*_Re(wc['C7_'+bq]*ffta*_Co(wc['C9_'+label]))
 
 def getF2(s, par, B, ff, ff0, lep, wc):
+    flavio.register_citation("Melikhov:2004mk")
     scale = config['renormalization scale']['bllgamma']
     mb = running.get_mb(par, scale, nf_out=5)
     mbh = mb/par['m_'+B]
@@ -85,6 +88,8 @@ def getF2(s, par, B, ff, ff0, lep, wc):
 
 
 def B10(s, par, B, ff, ff0, lep, wc):
+    flavio.register_citation("Melikhov:2004mk")
+    flavio.register_citation("Guadagnoli:2016erb")
     F1 = getF1(s, par, B, ff, ff0, lep, wc)
     F2 = getF2(s, par, B, ff, ff0, lep, wc)
     mlh = par['m_'+lep]/par['m_'+B]
@@ -94,11 +99,15 @@ def B10(s, par, B, ff, ff0, lep, wc):
 # B11 vanishes for the BR estimation
 
 def B12(s, par, B, ff, ff0, lep, wc):
+    flavio.register_citation("Melikhov:2004mk")
+    flavio.register_citation("Guadagnoli:2016erb")
     F1 = getF1(s, par, B, ff, ff0, lep, wc)
     F2 = getF2(s, par, B, ff, ff0, lep, wc)
     return s*(F1+F2)
 
 def B120(s, par, B, ff, ff0, lep, wc):
+    flavio.register_citation("Melikhov:2004mk")
+    flavio.register_citation("Guadagnoli:2016erb")
     scale = config['renormalization scale']['bllgamma']
     mb = running.get_mb(par, scale, nf_out=5)
     mbh = mb/par['m_'+B]
@@ -110,18 +119,21 @@ def B120(s, par, B, ff, ff0, lep, wc):
 
 
 def dG1dsMN(s, par, B, ff, ff0, lep, wc):
+    flavio.register_citation("Melikhov:2004mk")
     mlh = par['m_'+lep]/par['m_'+B]
     pref = prefactor(s, par, B, ff, lep, wc)*(1-s)**3*sqrt(1-4*mlh**2/s)
     return pref*(B10(s, par, B, ff, ff0, lep, wc) + (s - 4*mlh**2)/(3*s)*B12(s, par, B, ff, ff0, lep, wc))
 
 def dG2dsMN(s, par, B, ff, ff0, lep, wc):
+    flavio.register_citation("Melikhov:2004mk")
     mlh = par['m_'+lep]/par['m_'+B]
     bq = meson_quark[B]
     label = bq+lep+lep
     pref = 16*prefactor(s, par, B, ff, lep, wc)*(par['f_'+B]/par['m_'+B]*abs(wc['C10_'+label])*mlh)**2/(1-s)
     return pref*(-8*sqrt(s-4*mlh**2)*sqrt(s)+4*atanh(sqrt(1-4*mlh**2/s))*(1+s-mlh**2*(1-s)**2))
 
 def dG12dsMN(s, par, B, ff, ff0, lep, wc):
+    flavio.register_citation("Melikhov:2004mk")
     mlh = par['m_'+lep]/par['m_'+B]
     pref = 4*prefactor(s, par, B, ff, lep, wc)*par['f_'+B]/par['m_'+B]*(1-s)
     return pref*atanh(sqrt(1-4*mlh**2/s))*B120(s, par, B, ff, ff0, lep, wc)

flavio/physics/bdecays/bvll/observables_bs.py

@@ -34,6 +34,7 @@ def bsvll_obs(function, q2, wc_obj, par, B, V, lep):
 
 def S_theory_num_Bs(y, J, J_bar, J_h, i):
     # (42) of 1502.05509
+    flavio.register_citation("Descotes-Genon:2015hea")
     return 1/(1-y**2) * (J[i] + J_bar[i]) - y/(1-y**2) * J_h[i]
 
 def S_experiment_num_Bs(y, J, J_bar, J_h, i):
@@ -51,6 +52,7 @@ def S_experiment_Bs(y, J, J_bar, J_h, i):
 
 def dGdq2_ave_Bs(y, J, J_bar, J_h):
     # (48) of 1502.05509
+    flavio.register_citation("Descotes-Genon:2015hea")
     return (1/(1-y**2) * (observables.dGdq2(J) + observables.dGdq2(J_bar))
             - y/(1-y**2) * observables.dGdq2(J_h))/2.
 

flavio/physics/bdecays/bvll/qcdf.py

@@ -37,6 +37,8 @@ def get_input(par, B, V, scale):
 # see eqs. (18), (50), (79) of hep-ph/0106067v2
 # and eqs. (47), (48) of hep-ph/0412400
 def Cq34(q2, par, wc, B, V, scale):
+    flavio.register_citation("Beneke:2001at")
+    flavio.register_citation("Beneke:2004dp")
     # this is -C_q^12 (for q=u) or C_q^34 - eps_u * C_q^12 (for q=d,s) of hep-ph/0412400
     mB, mb, mc, alpha_s, q, eq, ed, eu, eps_u, qiqj = get_input(par, B, V, scale)
     T_t = wc['C3_'+qiqj] + 4/3.*(wc['C4_'+qiqj] + 12*wc['C5_'+qiqj] + 16*wc['C6_'+qiqj])
@@ -63,16 +65,18 @@ def T_para_minus_WA(q2, par, wc, B, V, scale):
     return -eq * 4*mB/mb * Cq34(q2, par, wc, B, V, scale)
 
 
-# B->V, weak annihaltion at O(1/mb)
+# B->V, weak annihilation at O(1/mb)
 
 # (51) of hep-ph/0412400
 
 def T_perp_WA_PowC_1(q2, par, wc, B, V, scale):
+    flavio.register_citation("Beneke:2004dp")
     mB, mb, mc, alpha_s, q, eq, ed, eu, eps_u, qiqj = get_input(par, B, V, scale)
     # NB, the remaining prefactors are added below in the function T_perp
     return -eq * 4/mb * (wc['C3_'+qiqj] + 4/3.*(wc['C4_'+qiqj] + 3*wc['C5_'+qiqj] + 4*wc['C6_'+qiqj]))
 
 def T_perp_WA_PowC_2(q2, par, wc, B, V, scale):
+    flavio.register_citation("Beneke:2004dp")
     mB, mb, mc, alpha_s, q, eq, ed, eu, eps_u, qiqj = get_input(par, B, V, scale)
     # NB, the remaining prefactors are added below in the function T_perp
     return eq * 2/mb * Cq34(q2, par, wc, B, V, scale)
@@ -84,11 +88,13 @@ def T_perp_WA_PowC_2(q2, par, wc, B, V, scale):
 
 # chromomagnetic dipole contribution
 def T_perp_plus_O8(q2, par, wc, B, V, u, scale):
+    flavio.register_citation("Beneke:2001at")
     mB, mb, mc, alpha_s, q, eq, ed, eu, eps_u, qiqj = get_input(par, B, V, scale)
     ubar = 1 - u
     return - (alpha_s/(3*pi)) * 4*ed*wc['C8eff_'+qiqj]/(u + ubar*q2/mB**2)
 
 def T_para_minus_O8(q2, par, wc, B, V, u, scale):
+    flavio.register_citation("Beneke:2001at")
     mB, mb, mc, alpha_s, q, eq, ed, eu, eps_u, qiqj = get_input(par, B, V, scale)
     ubar = 1 - u
     return (alpha_s/(3*pi)) * eq * 8 * wc['C8eff_'+qiqj]/(ubar + u*q2/mB**2)
@@ -148,6 +154,7 @@ def En_V(mB, mV, q2):
 
 # (27) of hep-ph/0106067v2
 def t_perp(q2, u, mq, par, B, V):
+    flavio.register_citation("Beneke:2001at")
     mB = par['m_'+B]
     mV = par['m_'+V]
     EV = En_V(mB, mV, q2)
@@ -162,6 +169,7 @@ def t_perp(q2, u, mq, par, B, V):
 
 # (28) of hep-ph/0106067v2
 def t_para(q2, u, mq, par, B, V):
+    flavio.register_citation("Beneke:2001at")
     mB = par['m_'+B]
     mV = par['m_'+V]
     EV = En_V(mB, mV, q2)
@@ -177,6 +185,7 @@ def B0diffBFS(q2, u, mq, mB):
 
 # (29) of hep-ph/0106067v2
 def B0(s, mq):
+    flavio.register_citation("Beneke:2001at")
     if s==0.:
         return -2.
     if 4*mq**2/s == 1.:
@@ -188,6 +197,7 @@ def B0(s, mq):
 
 # (30), (31) of hep-ph/0106067v2
 def i1_bfs(q2, u, mq, mB):
+    flavio.register_citation("Beneke:2001at")
     ubar = 1 - u
     iepsilon = 1e-8j
     mq2 = mq**2 - iepsilon
@@ -201,6 +211,7 @@ def i1_bfs(q2, u, mq, mB):
 
 # (32) of hep-ph/0106067v2
 def L1(x):
+    flavio.register_citation("Beneke:2001at")
     if x == 0.:
         return -(pi**2/6.)
     elif x == 1.:
@@ -219,12 +230,14 @@ def phiV(u, a1, a2):
 # moments of the B meson light-cone distribution amplitude as in
 # eq. (54) and surroundings of hep-ph/0106067v2
 def lB_minus(q2, par, B):
+    flavio.register_citation("Beneke:2001at")
     mB = par['m_'+B]
     mb = running.get_mb_pole(par)
     LambdaBar = mB - mb
     w0 = 2*LambdaBar/3
     return 1/(exp(-q2/mB/w0)/w0 * (-ei(q2/mB/w0) + 1j*pi))
 def lB_plus(par, B):
+    flavio.register_citation("Beneke:2001at")
     mB = par['m_'+B]
     mb = running.get_mb_pole(par)
     LambdaBar = mB - mb

flavio/physics/bdecays/bxgamma.py

@@ -40,6 +40,7 @@ def BRBXgamma(wc_obj, par, q, E0):
     cutoff $E_0$ in GeV
     (currently works only for `E0=1.6`).
     See arXiv:1503.01789 and references therein."""
+    flavio.register_citation("Misiak:2015xwa")
     scale = flavio.config['renormalization scale']['bxgamma']
     alphaem = 1/137.035999139 # this is alpha_e(0), a constant for our purposes
     bq = 'b' + q
@@ -55,6 +56,7 @@ def BRBXgamma(wc_obj, par, q, E0):
     # central value of non-perturbative correction (M. Misiak, private communication)
     P_nonpert = 0.00381745
     # eq. (2.1) of hep-ph/0609241
+    flavio.register_citation("Misiak:2006ab")
     return BRSL * abs(xi_t)**2/abs(Vcb)**2 * 6*alphaem/pi/C * (PE0 + P_nonpert + P_uncertainty)
 
 def PE0_BR_BXgamma(wc, par, q, E0):

flavio/physics/bdecays/bxll.py

@@ -287,6 +287,7 @@ def Phill_logqed(I, sh, mb, ml, alpha_s, alpha_e, wc, q, lep, mc):
     # Note the different powers of ast (\tilde\alpha_s) and k (\kappa) due to
     # the different normalisation of 1) C9 and C10 and 2) the overall
     # normalisation of phi_ll
+    flavio.register_citation("Huber:2015sra")
     e = np.zeros((11, 11), dtype=complex) # from 0 to 10
     if I == 'T' or I == 'L' or I == 'BR':
         e[7,7] = 8 * ast    * k    * sigmaij_I(I, 7, 7, sh)*wem(I, 7, 7, sh, mb, ml, scale, mc)
@@ -360,6 +361,7 @@ def f_u(alpha_e, alpha_s, alpha_s_mu0, mb, l1, l2):
     # contribution to the branching ratio, as C9 and C10 are formally of
     # O(\tilde \alpha_s \kappa). This is a numerically relevant choice (around
     # 15% change of the low-q^2 BR).
+    flavio.register_citation("Huber:2015sra")
     ast = alpha_s/4./pi
     k = alpha_e/alpha_s
     eta = alpha_s_mu0/alpha_s

flavio/physics/bdecays/bxll_qed.py

@@ -3,7 +3,7 @@
 
 See arXiv:1503.04849."""
 
-
+import flavio
 from flavio.math.functions import li2
 from math import log, pi, sqrt
 import numpy as np
@@ -368,6 +368,7 @@ def wem_79_low(sh, mb, ml, scale, mc):
 wem_dict_high['BR', 9, 9] = wem_99_high
 
 def wem(I, i, j, sh, mb, ml, scale, mc):
+    flavio.register_citation("Huber:2015sra")
     if sh < 0.5:
         return wem_dict_low[I, i, j](sh, mb, ml, scale, mc)
     else:

flavio/physics/bdecays/bxlnu.py

@@ -50,6 +50,7 @@ def _BR_BXclnu(par, wc_obj, lep, nu):
              + 12*gVS(xc, xl) * (wc['VR']*wc['T']).real
            )
     # eq. (26) of arXiv:1107.3100 + corrections (P. Gambino, private communication)
+    flavio.register_citation("Gambino:2011cq")
     r_BLO = ( 1
                  # NLO QCD
                  + alpha_s/pi * pc1(xc, mb)

flavio/physics/bdecays/formfactors/b_gamma/bgamma.py

@@ -1,3 +1,4 @@
+import flavio
 from flavio.classes import AuxiliaryQuantity, Implementation
 from flavio.config import config
 
@@ -7,6 +8,7 @@ def ff(q2, par, B):
 
     See hep-ph/0208256.pdf.
     """
+    flavio.register_citation("Kruger:2002gf")
     fB = par['f_'+B]
     mB = par['m_'+B]
     name = 'B->gamma KM '

flavio/physics/bdecays/formfactors/b_p/bcl.py

@@ -1,5 +1,6 @@
 from math import sqrt
 import numpy as np
+import flavio
 from flavio.physics.bdecays.formfactors.common import z
 from flavio.physics.bdecays.formfactors.b_p.isgurwise import isgur_wise
 
@@ -35,6 +36,7 @@ def ff(process, q2, par, n=3, t0=None):
 
     The standard convention defines the form factors $f_+$, $f_0$, and $f_T$.
     """
+    flavio.register_citation("Bourrely:2008za")
     pd = process_dict[process]
     mpl = par[process + ' BCL m+']
     m0 = par[process + ' BCL m0']
@@ -62,6 +64,7 @@ def ff_isgurwise(process, q2, par, scale, n=3, t0=None):
     and BCL parametrization (arXiv:0807.2722) for $f_0$ and $f_+$, but using
     an improved Isgur-Wise relation in the heavy quark limit for $f_T$.
     """
+    flavio.register_citation("Bourrely:2008za")
     pd = process_dict[process]
     mpl = par[process + ' BCL m+']
     m0 = par[process + ' BCL m0']