qqbarthreshold is hosted by Hepforge, IPPP Durham
QQbar_threshold
QQbar_threshold Namespace Reference

Classes

struct  contribution_factors
 Multiplicative coefficients for the potentials and current matching coefficients. More...
 
struct  Et_Gammat_range
 Grid range in rescaled energy and width coordinates. More...
 
class  integrator
 Class for integrating a function of type Integrand. More...
 
struct  options
 Options for \(q\bar{q}\) functions. More...
 
struct  scheme
 Mass scheme. More...
 
struct  threshold
 Center-of-mass energy of a quark threshold. More...
 

Enumerations

enum  scheme_name {
  PS, PSshift, OneS, OneSshift,
  Pole, MSshift
}
 Names of mass schemes. More...
 
enum  SM { none = 0, QED = 1 << 0, Higgs = 1 << 1 , all = (1 << 3) - 1 }
 Standard Model effects beyond QCD. More...
 
enum  production_channel { photon_only = 0, S_wave_only, all }
 Production channels. More...
 

Functions

double alpha_s_top (double mu, double alpha_s_mZ=QQbar_threshold::alpha_s_mZ)
 The strong coupling used for top-related functions. More...
 
double alpha_s_top (double mu, alpha_s_options alpha_s_mu0)
 The strong coupling used for top-related functions. More...
 
double alpha_s_bottom (double mu, double alpha_s_mZ=QQbar_threshold::alpha_s_mZ)
 The strong coupling used for bottom-related functions. More...
 
double alpha_s_bottom (double mu, alpha_s_options alpha_s_mu0)
 The strong coupling used for bottom-related functions. More...
 
double axial_Z_coupling (double T3, double)
 Axial vector coupling of fermion to Z boson. More...
 
double vector_Z_coupling (double T3, double e)
 Vector coupling of fermion to Z boson. More...
 
double ttbar_energy_level (int n, double mu, double m, order o, options opt=top_options())
 Computes the binding energy \(E_n\) of an S-wave \(t\bar{t}\) bound state. More...
 
double ttbar_energy_level (int n, double mu, top_properties M, order o, options opt=top_options())
 Computes the binding energy \(E_n\) of an S-wave \(t\bar{t}\) bound state. More...
 
double bbbar_energy_level (int n, double mu, double m, order o, options opt=bottom_options())
 Computes the binding energy \(E_n\) of an S-wave \(b\bar{b}\) bound state. More...
 
template<class Integrand >
double integrate (Integrand I, double x_min, double x_max)
 Compute the integral of a function. More...
 
template<class Integrand >
integrator< Integrand > make_integrator (Integrand I)
 helper function to construct an integrator object More...
 
void load_grid (char const *filename)
 Loads a precomputed grid. More...
 
void load_grid (std::string const &filename)
 
void load_grid (std::ifstream &&in)
 
void load_nonresonant_grid (char const *filename)
 Loads a precomputed grid for the non-resonant contribution. More...
 
void load_nonresonant_grid (std::string const &filename)
 
void load_nonresonant_grid (std::ifstream &&in)
 
std::string grid_directory ()
 Get directory with default grids. More...
 
Et_Gammat_range grid_range ()
 Get grid range. More...
 
bool is_shifted_mass (scheme_name s)
 Determines whether the given scheme is a shifted scheme. More...
 
scheme_name to_shifted_mass (scheme_name s)
 Returns the corresponding shifted scheme. More...
 
bool is_shifted_mass (scheme const &s)
 Determines whether the given scheme is a shifted scheme. More...
 
double bottom_pole_mass (double mu, double m, order o, options opt=bottom_options())
 Computes the pole mass from a bottom quark mass in a different scheme. More...
 
double top_pole_mass (double mu, double m, order o, options opt=top_options())
 Computes the pole mass from a top quark masses in a different scheme. More...
 
double top_pole_mass (double mu, top_properties M, order o, options opt=top_options())
 Computes the pole mass from a top quark masses in a different scheme. More...
 
double top_ISR_log (double sqrt_s)
 Computes the logarithm generated through initial state radiation. More...
 
double bottom_ISR_log (double sqrt_s)
 Computes the logarithm generated through initial state radiation. More...
 
double ISR_log (double sqrt_s, double alpha)
 Computes the logarithm generated through initial state radiation. More...
 
double structure_function (double x, double beta)
 Computes the electron (or positron) structure function. More...
 
double modified_structure_function (double x, double beta)
 Modified electron (or positron) structure function. More...
 
double luminosity_function (double x, double beta)
 Computes the electron (or positron) luminosity function. More...
 
double modified_luminosity_function (double t, double beta)
 Modified electron (or positron) luminosity function. More...
 
std::string bottom_internal_settings (double mu, double m, order o, options opt)
 Describe internal settings for bottom-related functions. More...
 
std::string bottom_internal_settings (double sqrt_s, double mu, double m, order o, options opt)
 Describe internal settings for bottom-related functions. More...
 
std::string top_internal_settings (double mu, double m, order o, options opt)
 Describe internal settings for top-related functions. More...
 
std::string top_internal_settings (double sqrt_s, scales mu, top_properties M, order o, options opt)
 Describe internal settings for top-related functions. More...
 
options bottom_options ()
 Default option settings for bottom-related functions. More...
 
options top_options ()
 Default option settings for top-related functions. More...
 
double ttbar_residue (int n, double mu, double m, order o, options opt=top_options())
 Computes the residue \(Z_n\) of an S-wave \(t\bar{t}\) bound state. More...
 
complex ttbar_residue (int n, double mu, top_properties M, order o, options opt=top_options())
 Computes the residue \(Z_n\) of an S-wave \(t\bar{t}\) bound state. More...
 
double bbbar_residue (int n, double mu, double m, order o, options opt=bottom_options())
 Computes the residue \(Z_n\) of an S-wave \(b\bar{b}\) bound state. More...
 
double top_width (double mu, double m, order o, options opt=top_options())
 Computes the top quark decay width. More...
 
double ttbar_width (int n, double mu, top_properties M, order o, options opt=top_options())
 Computes the decay width of a top-antitop bound state. More...
 
double ttbar_xsection (double sqrt_s, scales mu, top_properties M, order o, options opt=top_options())
 Computes the cross section \(e^+e^- \to W^+W^-b\bar{b}\). More...
 
double ttbar_xsection (double sqrt_s, double mu, top_properties M, order o, options opt=top_options())
 Computes the cross section \(e^+e^- \to W^+W^-b\bar{b}\). More...
 
double bbbar_xsection (double sqrt_s, double mu, double m, order o, options opt=bottom_options())
 Computes the cross section \(e^+e^- \to b\bar{b}\). More...
 
double ttbar_threshold (double mu, double m, order o, options opt=top_options())
 Computes the naive threshold for \(t\bar{t}\) production. More...
 
double bbbar_threshold (double mu, double m, order o, options opt=bottom_options())
 Computes the naive threshold for \(b\bar{b}\) production. More...
 
double bbbar_R_ratio (double sqrt_s, double mu, double m, order o, options opt=bottom_options())
 Bottom contribution to the hadronic R ratio. More...
 
double ttbar_R_ratio (double sqrt_s, scales mu, top_properties M, order o, options opt=top_options())
 Top contribution to the hadronic R ratio. More...
 

Variables

constexpr double alpha_s_mZ = 0.1184
 Default value for strong coupling at the scale mZ.
 
constexpr double alpha_mZ = 1/128.944
 QED coupling at the scale mZ.
 
constexpr double alpha_Y = 1/132.274
 QED coupling at the scale mu_alpha_Y.
 
constexpr double mu_alpha_Y = 10.2
 Typical scale for Y resonances.
 
constexpr double mZ = 91.1876
 Mass of the Z boson.
 
constexpr double mW = 80.385
 Mass of the W boson.
 
constexpr double me = 0.0005109989461
 Electron mass.
 
constexpr double G_F = 0.000011663787
 Fermi constant - this is exclusively used for calculating the top width.
 
constexpr double m_Higgs = 125.
 Default value for mass of the Higgs boson.
 
constexpr double alpha_QED = 1/137.035999074
 Fine structure constant.
 
constexpr double e_u = 2./3.
 Electric charge of the top quark in units of the positron charge.
 
constexpr double e_d = -1./3.
 Electric charge of the bottom quark in units of the positron charge.
 
constexpr double e = -1
 Electric charge of the electron in units of the positron charge.
 
constexpr double cw2 = mW*mW/(mZ*mZ)
 cosine of the weak mixing angle squared
 
constexpr double sw2 = 1 - cw2
 sine of the weak mixing angle squared
 
constexpr double T3_nu = 1./2.
 Weak isospin of neutrino.
 
constexpr double T3_e = -1./2.
 Weak isospin of electron.
 
constexpr double T3_u = 1./2.
 Weak isospin of top.
 
constexpr double T3_d = -1./2.
 Weak isospin of bottom.
 
constexpr double mb_SI = 4.203
 Reference scale-invariant mass for bottom quarks.
 
constexpr double mu_thr = 2*mb_SI
 Decoupling threshold for bottom quarks.
 
constexpr int nl_bottom = 4
 Number of light flavours for bottom-related functions.
 
constexpr int nl_top = 5
 Number of light flavours for top-related functions.
 
constexpr double mu_f_bottom = 2.
 Default PS scale for bottom.
 
constexpr double mu_f_top = 20.
 Default PS scale for top.
 
constexpr double invGeV2_to_pb = 389379300.
 Conversion factor from 1/GeV^2 to picobarn.
 

Detailed Description

Enumeration Type Documentation

◆ scheme_name

Names of mass schemes.

If a shifted scheme is used in the evaluation of an expression the mass is first converted to the pole scheme at the given order and then inserted into the expression for the pole mass.

For insertion schemes the leading-order conversion to the pole scheme is used. Higher-order corrections in the relation to the pole mass are treated as insertions of counterterms [Beneke:2014].

See also
is_shifted_mass
Enumerator
PS 

The potential-subtracted insertion scheme.

PSshift 

The potential-subtracted shift scheme.

OneS 

The 1S insertion scheme.

OneSshift 

The 1S shift scheme.

Pole 

The pole scheme.

MSshift 

The \(\overline{\rm MS}\) scheme.

◆ SM

enum SM
strong

Standard Model effects beyond QCD.

Effects can be combined with the | operator e.g.

options opt;
opt.beyond_QCD = SM::QED | SM::Higgs;

activates both QED and Higgs effects.

Enumerator
none 

No corrections beyond QCD.

QED 

Corrections due to photons.

Higgs 

Corrections due to the physical Higgs boson.

all 

All Standard Model corrections.

◆ production_channel

enum production_channel
strong

Production channels.

Enumerator
photon_only 

Production only via a virtual photon.

S_wave_only 

Production only via S-wave photon or Z.

all 

All possible production channels.

Function Documentation

◆ alpha_s_top() [1/2]

double QQbar_threshold::alpha_s_top ( double  mu,
double  alpha_s_mZ = QQbar_threshold::alpha_s_mZ 
)

The strong coupling used for top-related functions.

Parameters
murenormalisation scale
alpha_s_mZstrong coupling at the scale mZ
Returns
The value of the strong coupling theory at the given scale

Both input and output value for the strong coupling are in the theory with five active flavours.

◆ alpha_s_top() [2/2]

double QQbar_threshold::alpha_s_top ( double  mu,
alpha_s_options  alpha_s_mu0 
)

The strong coupling used for top-related functions.

Parameters
murenormalisation scale
alpha_s_mu0a pair {mu0, alpha_s} of a reference scale mu0 and the strong coupling at the scale mu0
Returns
The value of the strong coupling at the given renormalisation scale

Both input and output value for the strong coupling are in the theory with five active flavours.

◆ alpha_s_bottom() [1/2]

double QQbar_threshold::alpha_s_bottom ( double  mu,
double  alpha_s_mZ = QQbar_threshold::alpha_s_mZ 
)

The strong coupling used for bottom-related functions.

Parameters
murenormalisation scale
alpha_s_mZstrong coupling at the scale mZ
Returns
The value of the strong coupling theory at the given scale

For the input value five active flavours are assumed. The output value is given in the four-flavour theory. Decoupling is performed at the scale mu_thr assuming a scale-invariant bottom mass mb_SI.

◆ alpha_s_bottom() [2/2]

double QQbar_threshold::alpha_s_bottom ( double  mu,
alpha_s_options  alpha_s_mu0 
)

The strong coupling used for bottom-related functions.

Parameters
murenormalisation scale
alpha_s_mu0a pair {mu0, alpha_s} of a reference scale mu0 and the strong coupling at the scale mu0
Returns
The value of the strong coupling at the given renormalisation scale

For the input value, five active flavours are assumed if mu0 is above the bottom decoupling threshold mu_thr, and four active flavours if mu0 is below mu_thr. The output value is always given in the four-flavour theory. If necessary, decoupling is performed at the scale mu_thr assuming a scale-invariant bottom mass mb_SI.

◆ axial_Z_coupling()

double QQbar_threshold::axial_Z_coupling ( double  T3,
double   
)
inline

Axial vector coupling of fermion to Z boson.

Parameters
T3Weak isospin
Returns
a_f, the axial vector coupling to the Z boson

◆ vector_Z_coupling()

double QQbar_threshold::vector_Z_coupling ( double  T3,
double  e 
)
inline

Vector coupling of fermion to Z boson.

Parameters
T3Weak isospin
eElectric charge in units of positron charge
Returns
v_f, the vector coupling to the Z boson

◆ ttbar_energy_level() [1/2]

double QQbar_threshold::ttbar_energy_level ( int  n,
double  mu,
double  m,
order  o,
options  opt = top_options() 
)

Computes the binding energy \(E_n\) of an S-wave \(t\bar{t}\) bound state.

Parameters
nprincipal quantum number
murenormalisation scale
mtop quark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The binding energy \(E_n\).

The bound state mass is given by \(2m_t+E_n\), where \(m_t\) is the mass of the top quark in the given scheme.

The supported options and their default values are:

See also
options, bbbar_energy_level, ttbar_residue

◆ ttbar_energy_level() [2/2]

double QQbar_threshold::ttbar_energy_level ( int  n,
double  mu,
top_properties  M,
order  o,
options  opt = top_options() 
)

Computes the binding energy \(E_n\) of an S-wave \(t\bar{t}\) bound state.

Parameters
nprincipal quantum number
murenormalisation scale
Mtop quark mass and width
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The binding energy \(E_n\).

This function includes corrections due to a non-zero top-quark width.

The bound state mass is given by \(2m_t+E_n\), where \(m_t\) is the mass of the top quark in the given scheme.

The supported options and their default values are:

See also
options, bbbar_energy_level, ttbar_residue

◆ bbbar_energy_level()

double QQbar_threshold::bbbar_energy_level ( int  n,
double  mu,
double  m,
order  o,
options  opt = bottom_options() 
)

Computes the binding energy \(E_n\) of an S-wave \(b\bar{b}\) bound state.

Parameters
nprincipal quantum number
murenormalisation scale
mbottom quark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The binding energy \(E_n\).

The bound state mass is given by \(M_{\Upsilon(nS)} = 2m_b+E_n\), where \(m_b\) is the mass of the bottom quark in the given scheme.

The supported options and their default values are:

See also
options, ttbar_energy_level, bbbar_residue

◆ integrate()

double integrate ( Integrand  I,
double  x_min,
double  x_max 
)
inline

Compute the integral of a function.

Parameters
integrandIntegrand function
x_minLower bound of the integration.
x_maxUpper bound of the integration.
Returns
The value of \(\int_{\texttt{x\_min}}^{\texttt{x\_max}} dx {\texttt{integrand(x)}}\)

◆ make_integrator()

integrator<Integrand> QQbar_threshold::make_integrator ( Integrand  I)
inline

helper function to construct an integrator object

Parameters
IThe integrand
Returns
An integrator object for I

◆ load_grid() [1/3]

void QQbar_threshold::load_grid ( char const *  filename)

Loads a precomputed grid.

Parameters
filenameName of the grid to be loaded

For order > 1, the cross section functions require a precomputed grid. Some grids are already packaged with this library. New grids can be generated with the QQbarGridCalc.m Mathematica package

See also
bbbar_xsection, ttbar_xsection, QQbarGridCalc.m, QQbarThreshold

◆ load_grid() [2/3]

void QQbar_threshold::load_grid ( std::string const &  filename)

This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.

◆ load_grid() [3/3]

void QQbar_threshold::load_grid ( std::ifstream &&  in)

This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.

Parameters
instream to the grid to be loaded

◆ load_nonresonant_grid() [1/3]

void QQbar_threshold::load_nonresonant_grid ( char const *  filename)

Loads a precomputed grid for the non-resonant contribution.

Parameters
filenameName of the grid to be loaded

Explicitly loading a grid is only required for exotic values of the top mass or the Z mass. For these cases grids can be generated with the QQbarGridCalc.m Mathematica package

◆ load_nonresonant_grid() [2/3]

void QQbar_threshold::load_nonresonant_grid ( std::string const &  filename)

This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.

◆ load_nonresonant_grid() [3/3]

void QQbar_threshold::load_nonresonant_grid ( std::ifstream &&  in)

This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.

Parameters
instream to the grid to be loaded

◆ grid_directory()

std::string QQbar_threshold::grid_directory ( )

Get directory with default grids.

Returns
The directory containing the default grids

◆ grid_range()

Et_Gammat_range QQbar_threshold::grid_range ( )

Get grid range.

Returns
The coordinate range of the currently loaded grid.
See also
Et_Gammat_range

◆ is_shifted_mass() [1/2]

bool QQbar_threshold::is_shifted_mass ( scheme_name  s)
inline

Determines whether the given scheme is a shifted scheme.

Parameters
sThe name of the mass scheme
Returns
true for shifted mass schemes, false otherwise

◆ to_shifted_mass()

scheme_name QQbar_threshold::to_shifted_mass ( scheme_name  s)
inline

Returns the corresponding shifted scheme.

Parameters
sThe name of the original mass scheme
Returns
The corresponding shifted scheme

For example, this returns PSshift if the input scheme is PS. If the input scheme is a shifted scheme or the pole scheme the input scheme is returned.

◆ is_shifted_mass() [2/2]

bool QQbar_threshold::is_shifted_mass ( scheme const &  s)
inline

Determines whether the given scheme is a shifted scheme.

Parameters
sThe mass scheme
Returns
true for shifted mass schemes, false otherwise

◆ bottom_pole_mass()

double QQbar_threshold::bottom_pole_mass ( double  mu,
double  m,
order  o,
options  opt = bottom_options() 
)

Computes the pole mass from a bottom quark mass in a different scheme.

Parameters
murenormalisation scale
mbottom quark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
the bottom pole mass

The supported options and their default values are:

Example

This computes the pole mass for a PS bottom quark mass of 4.5 GeV at a scale of 5 GeV at NNNLO:

options opt = bottom_options();
opt.mass_scheme = {PSshift, 2};
double mb_PS = 4.5;
double mb_pole = pole_bottom_mass(5, mb_PS, N3LO, opt);
See also
mass_schemes.hpp

◆ top_pole_mass() [1/2]

double QQbar_threshold::top_pole_mass ( double  mu,
double  m,
order  o,
options  opt = top_options() 
)

Computes the pole mass from a top quark masses in a different scheme.

Parameters
murenormalisation scale
mtop quark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
the top pole mass

The supported options and their default values are:

Example

This computes the pole mass for a \(\overline{\rm MS}\) top quark mass \(m_t(m_t) = 162\) GeV at N3LO:

options opt = top_options();
double mt_mt = 162;
opt.mass_scheme = {MSbar, mt_mt};
double mt_pole = pole_top_mass(mt_mt, mt_mt, N3LO, opt);
See also
mass_schemes.hpp

◆ top_pole_mass() [2/2]

double QQbar_threshold::top_pole_mass ( double  mu,
top_properties  M,
order  o,
options  opt = top_options() 
)

Computes the pole mass from a top quark masses in a different scheme.

Parameters
murenormalisation scale
Mtop quark mass and width
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
the top pole mass

This function includes corrections from a non-zero width in the mass conversion. Currently, this is only relevant for the conversion from the 1S scheme.

The supported options and their default values are:

See also
mass_schemes.hpp

◆ top_ISR_log()

double QQbar_threshold::top_ISR_log ( double  sqrt_s)

Computes the logarithm generated through initial state radiation.

Parameters
sqrt_scenter of mass energy
Returns
\(\beta = -2*\alpha(m_Z)/\pi [\log(m_e^2/s) + 1]\)
See also
alpha_mZ, me, structure_function, luminosity_function, modified_structure_function, modified_luminosity_function

◆ bottom_ISR_log()

double QQbar_threshold::bottom_ISR_log ( double  sqrt_s)

Computes the logarithm generated through initial state radiation.

Parameters
sqrt_scenter of mass energy
Returns
\(\beta = -2*\alpha(m_\Upsilon)/\pi [\log(m_e^2/s) + 1]\)
See also
alpha_Y, me, structure_function, luminosity_function, modified_structure_function, modified_luminosity_function

◆ ISR_log()

double QQbar_threshold::ISR_log ( double  sqrt_s,
double  alpha 
)

Computes the logarithm generated through initial state radiation.

Parameters
sqrt_scenter of mass energy
alphaelectromagnetic coupling constant
Returns
\(\beta = -2*\alpha/\pi [\log(m_e^2/s) + 1]\)
See also
me, structure_function, luminosity_function, modified_structure_function, modified_luminosity_function

◆ structure_function()

double QQbar_threshold::structure_function ( double  x,
double  beta 
)

Computes the electron (or positron) structure function.

Parameters
xfraction of center-of-mass energy squared between 0 and 1
betavalue of the logarithm to be resummed, see ISR_log
Returns
the value of the structure function
See also
me, structure_function, luminosity_function, modified_structure_function, modified_luminosity_function

◆ modified_structure_function()

double QQbar_threshold::modified_structure_function ( double  x,
double  beta 
)

Modified electron (or positron) structure function.

Parameters
tconvolution parameter between 0 and 1
betavalue of the logarithm to be resummed, see ISR_log
Returns
the value of the modified structure function

The modified structure function \(\tilde{\Gamma}_{ee}(t, \beta)\) is given by \(\tilde{\Gamma}_{ee}(t, \beta) = 2/\beta t^{2/\beta-1} \Gamma_{ee}(1-t^{2/\beta}, \beta)\), where \(\Gamma_{ee}(1-t^{2/\beta}, \beta)\) is the electron structure function.

See also
structure_function, luminosity_function, modified_structure_function, modified_luminosity_function

◆ luminosity_function()

double QQbar_threshold::luminosity_function ( double  x,
double  beta 
)

Computes the electron (or positron) luminosity function.

Parameters
xfraction of center-of-mass energy squared between 0 and 1
betavalue of the logarithm to be resummed, see ISR_log
Returns
the value of the luminosity function

The luminosity function \(L\) is defined as \(L(x) = \int_x^1 \frac{dy}{y} f(y)f(x/y)\), where \(f\) is the electron structure function

See also
structure_function, luminosity_function, modified_structure_function, modified_luminosity_function

◆ modified_luminosity_function()

double QQbar_threshold::modified_luminosity_function ( double  t,
double  beta 
)

Modified electron (or positron) luminosity function.

Parameters
tconvolution parameter between 0 and 1
betavalue of the logarithm to be resummed, see ISR_log
Returns
the value of the modified luminosity function

The modified luminosity function \(\tilde{L}(t, \beta)\) is defined as \(\tilde{L}(t, \beta) = t^{1/\beta-1}/\beta L(1-t^{1/\beta}, \beta)\), where \(L\) is the luminosity function

See also
structure_function, luminosity_function, modified_structure_function, modified_luminosity_function

◆ bottom_internal_settings() [1/2]

std::string QQbar_threshold::bottom_internal_settings ( double  mu,
double  m,
order  o,
options  opt 
)

Describe internal settings for bottom-related functions.

Parameters
murenormalisation scale
mquark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
A string containing the names and values of internal variables

This function is intended for debugging purposes.

◆ bottom_internal_settings() [2/2]

std::string QQbar_threshold::bottom_internal_settings ( double  sqrt_s,
double  mu,
double  m,
order  o,
options  opt 
)

Describe internal settings for bottom-related functions.

Parameters
sqrt_scenter-of-mass energy
murenormalisation scale
mbottom quark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
A string containing the names and values of internal variables

This function is intended for debugging purposes.

◆ top_internal_settings() [1/2]

std::string QQbar_threshold::top_internal_settings ( double  mu,
double  m,
order  o,
options  opt 
)

Describe internal settings for top-related functions.

Parameters
murenormalisation scale
mquark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
A string containing the names and values of internal variables

This function is intended for debugging purposes.

◆ top_internal_settings() [2/2]

std::string QQbar_threshold::top_internal_settings ( double  sqrt_s,
scales  mu,
top_properties  M,
order  o,
options  opt 
)

Describe internal settings for top-related functions.

Parameters
sqrt_scenter-of-mass energy
murenormalisation scales
Mtop quark mass and width
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
A string containing the names and values of internal variables

This function is intended for debugging purposes.

◆ bottom_options()

options QQbar_threshold::bottom_options ( )

◆ top_options()

options QQbar_threshold::top_options ( )

◆ ttbar_residue() [1/2]

double QQbar_threshold::ttbar_residue ( int  n,
double  mu,
double  m,
order  o,
options  opt = top_options() 
)

Computes the residue \(Z_n\) of an S-wave \(t\bar{t}\) bound state.

Parameters
nprincipal quantum number
murenormalisation scale
mtop quark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The residue \(Z_n\).

The residue in the pole scheme is defined as \(Z_N = \frac{4m^2}{(2m+E_N)^2}c_v \bigg[c_v - \frac{E_N}{m}\frac{d_v}{3}\bigg]|\psi_N(0)|^2\,.\) \(E_N\) is the binding energy, \(\psi_N(0)\) is the wave function at the origin, and \(c_v\) and \(d_v\) are the hard current matching coefficients. Deviating from the general conversion rules, the residue \(Z_N^{\rm RS}\) in an insertion scheme is given by \(Z_N^{\rm RS} = m^{\rm RS} \times [Z_N/m^2]^{\rm RS}\), where the usual prescription for scheme conversion is applied to the term in square brackets.

The supported options and their default values are:

See also
options, ttbar_energy_level, bbbar_residue

◆ ttbar_residue() [2/2]

complex QQbar_threshold::ttbar_residue ( int  n,
double  mu,
top_properties  M,
order  o,
options  opt = top_options() 
)

Computes the residue \(Z_n\) of an S-wave \(t\bar{t}\) bound state.

Parameters
nprincipal quantum number
murenormalisation scale
Mtop quark mass and width
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The residue \(Z_n\).

This function includes corrections due to a non-zero top quark width.

The residue in the pole scheme is defined as \(Z_N = \frac{4m^2}{(2m+E_N)^2}c_v \bigg[c_v - \frac{E_N}{m}\frac{d_v}{3}\bigg]|\psi_N(0)|^2\,.\) \(E_N\) is the binding energy, \(\psi_N(0)\) is the wave function at the origin, and \(c_v\) and \(d_v\) are the hard current matching coefficients. Deviating from the general conversion rules, the residue \(Z_N^{\rm RS}\) in an insertion scheme is given by \(Z_N^{\rm RS} = m^{\rm RS} \times [Z_N/m^2]^{\rm RS}\), where the usual prescription for scheme conversion is applied to the term in square brackets.

The supported options and their default values are:

See also
options, ttbar_energy_level, bbbar_residue

◆ bbbar_residue()

double QQbar_threshold::bbbar_residue ( int  n,
double  mu,
double  m,
order  o,
options  opt = bottom_options() 
)

Computes the residue \(Z_n\) of an S-wave \(b\bar{b}\) bound state.

Parameters
nprincipal quantum number
murenormalisation scale
mbottom quark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The binding energy \(E_n\).

The residue in the pole scheme is defined as \(Z_N = \frac{4m^2}{(2m+E_N)^2}c_v \bigg[c_v - \frac{E_N}{m}\frac{d_v}{3}\bigg]|\psi_N(0)|^2\,.\) \(E_N\) is the binding energy, \(\psi_N(0)\) is the wave function at the origin, and \(c_v\) and \(d_v\) are the hard current matching coefficients. Deviating from the general conversion rules, the residue \(Z_N^{\rm RS}\) in an insertion scheme is given by \(Z_N^{\rm RS} = m^{\rm RS} \times [Z_N/m^2]^{\rm RS}\), where the usual prescription for scheme conversion is applied to the term in square brackets.

The supported options and their default values are:

See also
options, bbbar_energy_level, ttbar_residue

◆ top_width()

double QQbar_threshold::top_width ( double  mu,
double  m,
order  o,
options  opt = top_options() 
)

Computes the top quark decay width.

Parameters
murenormalisation scale
mtop quark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The decay width

The supported options and their default values are:

See also
options

◆ ttbar_width()

double QQbar_threshold::ttbar_width ( int  n,
double  mu,
top_properties  M,
order  o,
options  opt = top_options() 
)

Computes the decay width of a top-antitop bound state.

Parameters
nprincipal quantum number
murenormalisation scale
Mtop quark mass and width
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The decay width

The supported options and their default values are:

See also
options

◆ ttbar_xsection() [1/2]

double QQbar_threshold::ttbar_xsection ( double  sqrt_s,
scales  mu,
top_properties  M,
order  o,
options  opt = top_options() 
)

Computes the cross section \(e^+e^- \to W^+W^-b\bar{b}\).

Parameters
sqrt_scenter-of-mass energy
murenormalisation scales {soft, weak}
Mtop quark mass and width
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The cross section in picobarn

The supported options and their default values are:

See also
options, ttbar_threshold, threshold

◆ ttbar_xsection() [2/2]

double QQbar_threshold::ttbar_xsection ( double  sqrt_s,
double  mu,
top_properties  M,
order  o,
options  opt = top_options() 
)

Computes the cross section \(e^+e^- \to W^+W^-b\bar{b}\).

Parameters
sqrt_scenter-of-mass energy
muoverall value for all renormalisation scales
Mtop quark mass and width
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The cross section in picobarn

The supported options and their default values are:

See also
options, ttbar_threshold, threshold

◆ bbbar_xsection()

double QQbar_threshold::bbbar_xsection ( double  sqrt_s,
double  mu,
double  m,
order  o,
options  opt = bottom_options() 
)

Computes the cross section \(e^+e^- \to b\bar{b}\).

Parameters
sqrt_scenter-of-mass energy
musoft renormalisation scale
mbottom quark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The cross section in picobarn

The supported options and their default values are:

See also
options, bbbar_threshold, threshold

◆ ttbar_threshold()

double QQbar_threshold::ttbar_threshold ( double  mu,
double  m,
order  o,
options  opt = top_options() 
)

Computes the naive threshold for \(t\bar{t}\) production.

Parameters
murenormalisation scale
mtop quark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The center-of-mass energy of the \(t\bar{t}\) production threshold

The naive production threshold is twice the pole quark mass. If the mass is given in a different mass scheme the corresponding pole mass is used.

The supported options and their default values are:

◆ bbbar_threshold()

double QQbar_threshold::bbbar_threshold ( double  mu,
double  m,
order  o,
options  opt = bottom_options() 
)

Computes the naive threshold for \(b\bar{b}\) production.

Parameters
murenormalisation scale
mbottom quark mass
operturbative order (LO, NLO, N2LO, or N3LO)
optfurther options
Returns
The center-of-mass energy of the \(b\bar{b}\) production threshold

The naive production threshold is twice the pole quark mass. If the mass is given in a different mass scheme the corresponding pole mass is used.

The supported options and their default values are:

◆ bbbar_R_ratio()

double QQbar_threshold::bbbar_R_ratio ( double  sqrt_s,
double  mu,
double  m,
order  o,
options  opt = bottom_options() 
)

Bottom contribution to the hadronic R ratio.

The bottom contribution to the hadronic R ratio is defined by bbbar_R_ratio = \(\sigma_b/\sigma_{pt} \) with \(\sigma_{pt} = \frac{4 \pi \alpha(\mu_\alpha)^2}{3s}\), where \(\sigma_b\) is the bottom production cross section bbbar_xsection converted to natural units. \(\alpha(\mu_\alpha)\) is the running on-shell QED coupling at the scale \(\mu_\alpha\) set by the options.

The supported options and their default values are:

◆ ttbar_R_ratio()

double QQbar_threshold::ttbar_R_ratio ( double  sqrt_s,
scales  mu,
top_properties  M,
order  o,
options  opt = top_options() 
)

Top contribution to the hadronic R ratio.

The bottom contribution to the hadronic R ratio is defined by ttbar_R_ratio = \(\sigma_t/\sigma_{pt} \) with \(\sigma_{pt} = \frac{4 \pi \alpha(\mu_\alpha)^2}{3s}\), where \(\sigma_t\) is the top production cross section bbbar_xsection converted to natural units. \(\alpha(\mu_\alpha)\) is the running on-shell QED coupling at the scale \(\mu_\alpha\) set by the options.

The supported options and their default values are: