1 Introduction

The search for Higgs bosons [16] is, and has been, a major cornerstone of the physics programmes of past, present and future high energy colliders. This has become even more important in view of the recent discovery of a Higgs signal by ATLAS [7, 8] and CMS [911]. Determining the properties of this newly observed state and comparing the measurements to explicit theories beyond the Standard Model (SM) is one of the present challenges. These theories often contain enlarged Higgs sectors with multiple Higgs bosons. Even in the presence of a signal, it is therefore important that the LEP, Tevatron and LHC experiments present exclusion limits from the non-observation of Higgs bosons in various channels. These are very useful for constraining the available parameter space of those models which are able to fit correctly the observed Higgs signal. Such constraints will need to be taken into account also in the future interpretation of the Higgs results in the context of models of new physics.

In this paper we describe new developments in version 4 of the publicly available Fortran code HiggsBounds [1214], which has been designed for exactly this purpose. For the complementary approach, to test whether a model is compatible with the observed LHC Higgs signal (and possible future signals of additional Higgs bosons), we have also developed the sister program HiggsSignals, which has been described in Ref. [15]. It is highly recommended to use these two programs in parallel to obtain the most complete test for extensions of the SM Higgs sector.

The experimental analyses implemented in HiggsBounds usually take one of two forms. Dedicated analyses have been carried out in order to constrain some of the most popular models, such as the SM [711, 16] and various benchmark scenarios in the Minimal Supersymmetric Standard Model (MSSM) [1720]. In addition, model-independent limits on the cross sections of individual signal topologies (such as \(e^+e^-\rightarrow h_iZ\rightarrow b\bar{b}Z\)) have been published. In the former type of analyses several search channels (or signal topologies) are typically combined in order to maximize the discovery/exclusion reach. However, the re-interpretation of these results in the context of different models than those already investigated by the search analysis requires detailed knowledge of the individual efficiencies (or signal contaminations) of the investigated search channels. In contrast, the latter type of analysis can be used easily to test a wider class of models.

HiggsBounds-4 has been designed to facilitate the task of comparing Higgs sector predictions with existing exclusion limits, thus allowing the user to quickly and conveniently check a wide variety of models against the state-of-the-art results from Higgs searches. Version 4 differs significantly from previous versions of the code (described in [12, 13, 21]) in several respects. The code now fully supports testing models against exclusion limits from the LHC, which are implemented for analyses performed at center-of-mass energies of both \(\sqrt{s}=7 \) and 8 TeV. The main algorithm of HiggsBounds has been extended to ensure a reliable application of exclusion limits in the presence of a signal (as is now observed in the LHC data). The model-likeness test, which tests whether a given model fulfills the assumptions of a particular Higgs search to a sufficient degree, has been fully rewritten to enable in particular the limits from SM Higgs searches at the LHC to be applied in a wider context. We introduce an option to take into account theory uncertainties on the Higgs mass predictions, which are relevant, for instance, for the lightest Higgs boson mass in the MSSM. An alternative statistical treatment for the LEP constraints (in the form of a \(\chi ^2\) output) is provided. Finally, we describe an improved input/output for supersymmetric (SUSY) models that can now be given in the SLHA format [22, 23]. The main focus of this updated documentation is to provide a detailed description of these new developments, to show relevant physics examples of where improvements can be expected, and to introduce the user to how the improved HiggsBounds code can be used in practice.

This paper is structured as follows. In Sect. 2 we give a general introduction to the statistical approach employed in the HiggsBounds code, and describe, in particular, the way in which this approach has been extended in HiggsBounds-4. Section 3 gives a thorough description of the different methods of providing theory input for HiggsBounds, and their extension to LHC7/8 predictions. This is followed by Sect. 4, which contains a discussion of the major new developments, including numerical examples. Finally, Sect. 5 contains the technical details on how these new features can be used in practice, extending the original HiggsBounds manual [12, 13] with a description of the new subroutines, data files and example programs that are provided. In an Appendix we list and provide references to all the experimental analyses that provide results implemented in the code, including the analyses added in the latest public version (HiggsBounds-4.1).

2 General approach of HiggsBounds-4

The general concept of HiggsBounds, including details on the treatment of limits from LEP and the Tevatron, has already been published in [1214] (see also Ref. [21]). From a conceptual point of view, the extension of HiggsBounds to include LHC limits is straightforward. The technicalities of this implementation, and how it modifies the user input, is discussed in Sect. 3. Our aim here is to give a brief introduction to the purpose of the code and the methods it uses. We also introduce one conceptual change with respect to previous versions, which has been prompted by the application of HiggsBounds to models which feature a Higgs boson with a mass close to the observed LHC Higgs signal.

The basic input for HiggsBounds (which the user has to provide) are the relevant physical quantities predicted for the Higgs sector of the model under consideration. The necessary predictions for each Higgs boson \(H_i \; (i=1,\ldots , n_{H^0}+ n_{H^{\pm }})\) in the model are, schematically,

$$\begin{aligned}&M_{H_i},\quad \Gamma _{{\text {tot}}}(H_i),\quad \text {BR}_{\text {model}}(H_i\rightarrow ...),\quad \frac{\sigma _{\text {model}}(P(H_i))}{\sigma _{{\text {ref}}}(P(H))}, \end{aligned}$$

i.e. the Higgs boson mass, its total decay width (it is assumed that the narrow width approximation holds), its decay branching ratios, and the production cross sections, normalized to a particular reference value. Here, \(P\) denotes a specific Higgs production process. If \(P\) exists in the SM, its cross section, \(\sigma ^\text {SM}(P(H))\), evaluated at the same mass value, \(M_H = M_{H_i}\), is typically used as the reference cross section, \(\sigma _{{\text {ref}}}\). In some cases it can also be necessary to supply additional predictions, such as the \(\text {BR}(t\rightarrow bH^+)\), or the \({\mathcal {CP}}\) properties of the neutral Higgs bosons. Variations on the input format are offered, which allow the user to specify a simpler set of input quantities when certain basic approximations are valid. A complete list of the options for giving model input is given in Sect. 3.

In addition to the model predictions, the other important ingredient of HiggsBounds is the experimental data. Exclusion limits from (negative) Higgs searches are collected from the experimental publications, with the aim of keeping the code as up-to-date as possible with the latest developments. Currently the code includes results from LEP, the Tevatron and the LHC experiments. More information on which analyses are available in HiggsBounds is provided in Appendix 7.1. The data for these analyses is contained in tables of expected exclusion limits at 95 % C.L. in the absence of a signal (based on Monte Carlo simulations), and the corresponding observed limits, as a function of the Higgs boson mass. The list consists both of analyses for which model-independent limits were published, and of dedicated analyses carried out specifically under the assumption of the SM (like most LHC searches to date), or for Higgs bosons with certain \({\mathcal {CP}}\) properties. These limits can be applied to models with Higgs bosons which show these characteristics to a sufficient degree Footnote 1 at the parameter point being considered.

We call the application of the limit from a particular Higgs search to one of the Higgs bosons of the model under study (or to two of the Higgs bosons, for searches involving two Higgs bosons) an “analysis application”, which we denote by \(X\) in the following.Footnote 2 Each analysis application has a corresponding signal cross section prediction \(\sigma (X)\), which HiggsBounds uses to calculate the relevant quantity \(Q_{\text {model}}(X)\) for which the experimental limit is given; typically this is a conveniently normalized cross section times a branching ratio. The corresponding experimental quantities are denoted \(Q_{\text {expec}} (X)\) and \(Q_\text {obs}(X)\) for the expected and observed limits, respectively. If two Higgs bosons have a narrow mass separation, \(\delta M=M_{h_i}-M_{h_j}\), then their predicted cross sections are added for certain analyses where the mass resolution is limited and interference effects are expected to be negligible. The settings for the maximal \(\delta M_h\) can be varied by the user separately for LEP, Tevatron, and LHC analyses (the default values are \(0\) GeV for LEP and \(10\) GeV for Tevatron/LHC).

HiggsBounds operates by considering, for each analysis application, the ratio of the model predictions, \(Q_{\text {model}}(X)\), to the experimental limits. To ensure that the result can be interpreted as an exclusion at \(95\,\%\) C.L. (which is the same confidence level as adopted by the individual analyses), it is crucial that the model prediction is only compared to the experimentally observed limit for one particular analysis application. In a first step, HiggsBounds therefore uses the expected experimental limits to determine the analysis application \(X_0\) with the highest statistical sensitivity to exclude the model point under consideration,

$$\begin{aligned} X_0=X:\max \frac{Q_{\text {model}}(X)}{Q_\text {expec}(X)}. \end{aligned}$$
(1)

In the second step, HiggsBounds then performs the exclusion test for the Higgs boson and analysis combination represented by \(X_0\), by computing the ratio to the observed limit

$$\begin{aligned} k_0=\frac{Q_{\text {model}}(X_0)}{Q_\text {obs}(X_0)}. \end{aligned}$$
(2)

If \(k_0>1\), HiggsBounds concludes that this parameter point of the tested model is excluded at \(95 \,\%\) C.L.Footnote 3

The statistical method as described here (in the following referred to as the classic method) has been the only mode of operation available in previous HiggsBounds versions. For HiggsBounds-4, we have extended this method to perform better in situations where a Higgs boson signal is present (as is now the case in the LHC data). The problem of the classic method arises for models with multiple Higgs bosons. If one of these has a mass close to that of the observed signal (which is likely, since any reasonable model should also explain this signal), its analysis applications will test the model predictions against limits (for various channels) in the signal region. In this region, the expected limits (based on the background-only hypothesis) will continue to improve with more experimental data and optimized analysis methods, whereas the observed limits can never be expected to reach exclusion at the SM level (provided a true signal of near-SM strength is what is observed). For model points where the most sensitive analysis application \(X_0\) is a test of the signal-like Higgs boson, the classic HiggsBounds method would therefore never yield exclusion. Moreover, constraints on the remaining Higgs spectrum (with less expected sensitivity) are not applied. Even if the exclusion remains formally valid at \(95\,\%\) C.L., it could be anticipated that this problem would eventually become serious enough to limit the usability of the code.

Among the several possible ways that the HiggsBounds algorithm could be extended to address this problem, all involving different compromises, we have opted for a solution which involves a slight violation of the strict testing of only one experimental limit. We call this the full HiggsBounds method. In summary, this method performs the original HiggsBounds test separately for each Higgs boson in the model. In the full HiggsBounds method, the first step is to evaluate the most sensitive analysis application \(X_i\) for each Higgs boson \(H_i\) according to

$$\begin{aligned} X_i=X(H_i):\max \frac{Q_{\text {model}}\left( X(H_i)\right) }{Q_\text {expec}\left( X(H_i)\right) }. \end{aligned}$$
(3)

This is followed by a straightforward exclusion test on the individually most sensitive analysis applications

$$\begin{aligned} k_i=\frac{Q_{\text {model}}(X_i)}{Q_\text {obs}(X_i)}. \end{aligned}$$
(4)

The result of these tests contains more information than the single test of HiggsBounds classic (such as exclusion/non-exclusion by individual Higgs bosons), which is now made available to the user (see Sect. 5 for details). A combined HiggsBounds exclusion is also calculated, where the result is interpreted as model exclusion if \(k_i>1\) for any of the \(k_i\). The combined (single-number) output is then calculated as

$$\begin{aligned}&k_0=\max _i k_i,\end{aligned}$$
(5)
$$\begin{aligned}&X_0=X_i : \max _i k_i . \end{aligned}$$
(6)

By the construction of the full method, it follows directly that the two methods are equivalent for models with a single Higgs boson. It is also clear that the full method can only give stronger exclusion than the classic method. This is consistent with the fact that the exclusion of the full method will correspond to a limit at somewhat lower statistical confidence level than \(95\,\%\). Still, the deviation from the strict \(95\,\%\) C.L. should be minor in this approach compared to the alternative (naive) testing of all Higgs bosons versus all observed limits, since the number of Higgs bosons in a model in general is much smaller than the number of implemented experimental analyses. Furthermore, a non-negligible dilution of the \(95\,\%\) C.L. interpretation of the combined result is only expected in the case where more than one test \(X_i\) leads to a ratio \(k_i \approx 1\).

To illustrate the difference between the classic and full methods of HiggsBounds, we show in Fig. 1 three versions of the excluded regions in an MSSM benchmark scenario, the so-called \(M_h^\mathrm{mod+}\) scenario [24]. The MSSM has three neutral Higgs bosons (\(h,H,A\)), where in this scenario the \(h\) boson can have a mass close to the LHC signal around \(M_h\sim 125\,\, \mathrm {GeV}\) (this region, considering a \(2\) \((3)\) GeV total uncertainty on \(M_h\) is indicated by dark (light) green colour in the figure). The exclusion bounds, as evaluated by HiggsBounds, are shown separately for LEP exclusion (blue) and the LHC (red). When evaluating the limits in this figure, a theory uncertainty of \(3\) GeV is taken into account in the evaluation of the lightest Higgs mass, see Sect. 4.2 for details on how this is done. As can be seen from this figure, the full method gives the strongest exclusion, corresponding to the most accurate application of the existing limits in this scenario (as also used in [24]). The difference to the classic method can be seen in particular for high \(M_A\) and high \(\tan \beta \) (the decoupling regime). Here the applicability of the classic method is limited, since the globally most sensitive channel is a search for the lightest (SM-like) Higgs boson, which cannot be excluded when is mass its in the signal region, \(M_h\simeq 125\,\, \mathrm {GeV}\). This is in contrast to the results in the full method, which can be further illustrated by looking at the contribution of individual Higgs bosons as shown in Fig. 2 for the same MSSM example. The first panel shows the exclusion contributed by \(h\). The narrow unexcluded region around \(M_A=135\) GeV results from a particular channel (\(pp\rightarrow VH\), \(H\rightarrow b\bar{b}\)) being the most sensitive. For this channel, the observed limit is not strong enough here to exclude the lightest Higgs. The second panel shows the exclusion for \(H/A\). They are treated together, since their masses are close to degenerate over most of the parameter space. The dominant exclusion therefore comes from the same search channels and their signal rates are added. Finally, the last plot shows the exclusion from \(H^\pm \). The exclusion region presented for the full method in Fig. 1 consists of the union of the three different exclusion regions shown here. In the HiggsBounds distribution we provide an updated example program, HBwithFH, which can be used to test MSSM parameter points for exclusion using either the full or the classic HiggsBounds methods.

Fig. 1
figure 1

Exclusion regions in the MSSM parameter space for the \(M_h^\mathrm{mod+}\) benchmark scenario [24]. Results from HiggsBounds full (left) are compared to the results from HiggsBounds classic (right). The colours show exclusion by the LHC (red), LEP (blue), and the favored region where \(M_h=125.7\pm 2\) GeV (dark green), \(M_h=125.7\pm 3\) GeV (light green)

Fig. 2
figure 2

Contribution to the full HiggsBounds exclusion in the MSSM parameter space for the \(M_h^\mathrm{mod+}\) benchmark scenario [24] from exclusion of the individual Higgs bosons: \(h\) (left), \(H/A\) (center), and \(H^\pm \) (right). The colour coding is the same as in Fig. 1

3 Theoretical predictions

The theoretical model predictions, which are compared to the experimental data in the analysis applications, are computed from the user input. A detailed explanation of this input has been given for previous version of HiggsBounds [12, 13]. For completeness, we include the full input specification here, describing both the original input and the updates in HiggsBounds-4. For the input of the required theory predictions to HiggsBounds, the user can choose between four different formats. In the command line version of the program, these are labelled by the variable whichinput, which can take the values hadr (hadronic), part (partonic), effC (effective couplings), or SLHA (SUSY Les Houches Accord). The type of input required for each of these formats is summarized in Table 1, and the following subsections contain more detailed information about each type of input. In particular, we describe the internal treatment of LHC cross sections (at \(7\) and \(8\) TeV). Technical details on the subroutines which should be used for the different input options can be found in Sect. 5.

Table 1 The possible values of the variable whichinput, which indicates the format of the theoretical predictions provided by the user for the neutral Higgs sector

3.1 Hadronic input

The hadronic input option whichinput=hadr requires model input in the most general form and therefore contains the lowest degree of approximation. The complete input in this format involves the specification of:

  1. (1)

    The masses for the neutral Higgs bosons, \(h_k\) (\(k=1,\ldots , n_{H^0})\), and the (singly) charged Higgs bosons, \(H_k^{\pm }\) \((k=1,\ldots , n_{H^{\pm }})\), in the model.

    $$\begin{aligned} M_{h_k},\, M_{H_k^{\pm }}, \end{aligned}$$
  2. (2)

    their total decay widths,

    $$\begin{aligned}&\Gamma _{\text {tot}}(h_k),\,\Gamma _{\text {tot}}(H_k^{\pm }), \end{aligned}$$
  3. (3)

    whether the neutral Higgs bosons are \({\mathcal {CP}}\)-even, \({\mathcal {CP}}\)-odd or states of mixed \({\mathcal {CP}}\),

  4. (4)

    the neutral Higgs branching ratios which have SM equivalents

    $$\begin{aligned}&\text {BR}_{\text {model}}(h_k\rightarrow \text {SM}),\\&\text {with SM } = \{s\bar{s}, c\bar{c}, b\bar{b}, \mu ^+ \mu ^- , \tau ^+\tau ^- , W^+W^- , ZZ,\\&\qquad \qquad \qquad \quad Z\gamma , \gamma \gamma , gg\} , \end{aligned}$$
  5. (5)

    the neutral Higgs branching ratios without SM equivalents

    $$\begin{aligned}&\text {BR}_{\text {model}}(h_k\rightarrow h_i h_j), \\&\text {BR}_{\text {model}}(h_k\rightarrow \mathrm{invisible}), \end{aligned}$$
  6. (6)

    the charged Higgs branching ratios to SM particles

    $$\begin{aligned}&\text {BR}(H_j^{+} \rightarrow {\text {SM}}),\,\mathrm{where}\, \text {SM}= \{c\bar{s}, c\bar{b}, \tau ^+ \nu _{\tau }\}, \end{aligned}$$
  7. (7)

    the top quark branching ratios

    $$\begin{aligned}&\text {BR}(t \rightarrow W^{+}b),\\&\text {BR}(t \rightarrow H_j^{+}b), \end{aligned}$$
  8. (8)

    normalized cross sections ratios \(R_{\sigma }(P)\) (see below for the definition), for the LEP Higgs production processes

    $$\begin{aligned}&e^+e^- \rightarrow h_j Z ,\\&e^+e^- \rightarrow b \bar{b} h_j ,\\&e^+e^- \rightarrow \tau ^+ \tau ^- h_j ,\\&e^+e^- \rightarrow h_i h_j ,\\&e^+e^- \rightarrow H^+_j H^-_j , \end{aligned}$$
  9. (9)

    normalized ratios \(R_{\sigma }(P)\) (see below) of hadronic Higgs production cross sections at the Tevatron for the processes

    $$\begin{aligned}&p\bar{p} \rightarrow h_j ,\\&p\bar{p} \rightarrow b h_j ,\\&p\bar{p} \rightarrow h_j W ,\\&p\bar{p} \rightarrow h_j Z ,\\&p\bar{p} \rightarrow h_j \mathrm{\, via \, VBF} ,\\&p\bar{p} \rightarrow t \bar{t} h_j , \end{aligned}$$
  10. (10)

    and, finally, normalized ratios \(R_{\sigma }(P)\) of hadronic Higgs production cross sections at the LHC (both for \(\sqrt{s}=7\) and \(8\,\, \mathrm {TeV}\), given as separate inputs) for the processes

    $$\begin{aligned}&pp \rightarrow h_j ,\\&pp \rightarrow b h_j ,\\&pp \rightarrow h_j W ,\\&pp \rightarrow h_j Z ,\\&pp \rightarrow h_j \mathrm{\, via \, VBF} ,\\&pp \rightarrow t \bar{t} h_j . \end{aligned}$$

It is important to note that this corresponds to an exhaustive list of the possible input. In certain (most) cases only a subset of these inputs may be required. For example, if the user is only interested in limits from the LHC, no LEP or Tevatron cross sections need to be given as input (the corresponding values can be set to zero). Likewise, for models with only neutral Higgs bosons, no input involving either the charged Higgs sector or top decays will be used. The meaning of most of these quantities should be pretty clear already from the notation; for those that require further clarifications we provide some more details below.

For input items (8), (9) and (10), the normalized cross section of a Higgs production process \(P\) is simply defined by

$$\begin{aligned} R_{\sigma }(P)=\frac{\sigma _\mathrm{model} (P)}{\sigma _\mathrm{ref} (P)}. \end{aligned}$$
(7)

Where a SM equivalent exists, the reference cross section \(\sigma _\mathrm{ref} (P)\) for Higgs boson \(h_k\) is \(\sigma _\mathrm{ref} (P)=\sigma _\mathrm{SM} (P)\), evaluated for \(M^\mathrm{SM}_H=M_{h_k}\). The only neutral Higgs production process considered in HiggsBounds up to now which does not have a SM equivalent is \(e^+e^- \rightarrow h_j h_i\). In this case, we choose as the reference cross section a fictitious production process for two scalar particles (\(H'\), \(H\)) with masses \(M_{H'}=M_{h_j}\) and \(M_H = M_{h_i}\) that proceeds via a virtual \(Z\) exchange with a standardized squared coupling constant \(g_{H'HZ}^2 =e^2/(4 s_\mathrm {w}^2 c_\mathrm {w}^2)\), where \(e\) denotes the electromagnetic coupling constant, and \(s_\mathrm {w}\) (\(c_\mathrm {w}\)) the sine (cosine) of the weak mixing angle, respectively.Footnote 4 Similarly, for the process \(P=e^+e^- \rightarrow H^+_j H^-_j\), the reference cross section is the cross section of the process \(e^+e^- \rightarrow H^+ H^-\) in a two-Higgs-doublet model (e.g. the MSSM) at tree-level (i.e. with \(s\)-channel \(\gamma \) and \(Z\) exchange). This reference cross section depends solely on the mass of the charged Higgs boson and SM quantities. As a consequence, the leading-order prediction in the MSSM is \(R_{\sigma }(e^+e^- \rightarrow H^+_j H^-_j)=1\).

The branching ratio to “invisible”, \(\text {BR}_{\text {model}}(h_k\rightarrow \mathrm invisible)\), is defined as the branching ratio of a neutral Higgs boson into particles which are only infered in the detector by their contribution to the missing transverse energy. Examples of this includes the lightest neutralino in the MSSM [25, 26], inert scalars [2729], or majorons in supersymmetric models with spontaneous breaking of \(\mathcal {R}\)-parity [3032].

The hadronic input is the most versatile, since it allows the user to provide the predictions in the form of the most precise calculations available. It is also the only input format which allows for studying e.g. effects of parton distribution functions or hadronic uncertainties on the Higgs exclusion bounds. On the other hand, this input format is also the most demanding, and in order to make it more convenient for the user to correctly normalize the rate predictions, the HiggsBounds library provides a series of Fortran functions which allow the user to access the predictions of certain SM quantities, including the hadronic SM Higgs production cross sections, total decay width, and branching ratios as a function of the Higgs mass.

3.2 Parton-level input

The second possibility for specifying the HiggsBounds input (with whichinput=part) is using ratios of partonic cross sections as far as is possible. This input format is in many cases more convenient for the user to calculate than the hadronic option. It requires (at most) the following model predictions

(1)–(8) the same as for whichinput=hadr,

(9) normalized ratios \(R_{\sigma }(P)\) (as defined by Eq. 7) of hadronic Higgs production cross sections at the Tevatron (\(\sqrt{s}=1.96 \,\, \mathrm {TeV}\)) and the LHC (at \(\sqrt{s}=7/8\,\, \mathrm {TeV}\)) for the neutral Higgs production processes

$$\begin{aligned}&p\bar{p} \rightarrow h_j \mathrm{\, via \, VBF}\\&p\bar{p} \rightarrow t \bar{t} h_j \\&p p \rightarrow h_j \mathrm{\, via \, VBF}\\&p p \rightarrow t \bar{t} h_j \end{aligned}$$

(10) normalized ratios \(R^{h_j+y}_{nm}\) (defined below) of parton-level cross sections for neutral Higgs production, which are assumed to be valid both at the Tevatron and the LHC, for the following processes

$$\begin{aligned}&gg,b\bar{b} \rightarrow h_j,\\&u\bar{d}, c\bar{s} \rightarrow h_j W^+, \\&\bar{u}d,\bar{c}s \rightarrow h_j W^-,\\&gg,q\bar{q} \rightarrow h_j Z\; (q=d,u,s,c,b), \\&bg \rightarrow h_j b. \end{aligned}$$

The normalized cross section ratio \(R^{h_j+y}_{nm}\) for a partonic neutral Higgs production process, \(nm\rightarrow h_j + y\), is defined by

$$\begin{aligned} R^{h_j+y}_{nm}=\frac{\hat{\sigma }^\mathrm{model}_{nm\rightarrow h_{j}+y} }{\hat{\sigma }^\mathrm{SM}_{nm\rightarrow H+y} }. \end{aligned}$$
(8)

It should be calculated for a parton-system center-of-mass energy squared \(\hat{s}=\hat{s}_0\), where \(\hat{s}_0\) denotes the partonic production threshold \(\hat{s}_0=(M_{h_j}+M_y)^2\). For this approximation to be valid the dependence of the partonic cross section on \(\hat{s}\) is required to be mild. For the case of single Higgs boson production, \(M_y=0\).

The partonic cross section ratios \(R^{h_j+y}_{nm}\) can be a lot easier to calculate than hadronic cross section ratios for a wide range of models. In addition, it is often (at least to a good approximation) the case that

$$\begin{aligned} R^{h_j +W^+}_{nm}=R^{h_j +W^-}_{nm}=R^{h_j +Z}_{q\bar{q}} \end{aligned}$$

for all \(nm\). This reduces the number of partonic cross section ratios which need to be provided by the user from \(13\) to \(5\). An example of a model of this type is given by the MSSM with real parameters, where the common ratio can be calculated approximately from the normalized (squared) effective coupling of the Higgs boson to a pair of \(Z\) bosons, \((g^{\text {model}}_{h_jZZ}/g^\text {SM}_{HZZ})^2\). For a more complete discussion of the use of the efffective coupling approximation as input to HiggsBounds, see Sect. 3.3.

Internally, HiggsBounds uses the approximate relations

$$\begin{aligned} R^{\mathrm {TEV}}_\sigma (P)&\displaystyle \simeq \sum \limits _{\{n,m\}} R_{nm}^{h_j+y} \frac{ \sigma _\text {SM}(p\bar{p} \rightarrow n m \rightarrow H+y) }{\sigma _\text {SM}(p\bar{p} \rightarrow H+y)}\,, \end{aligned}$$
(9)
$$\begin{aligned} R^{\mathrm {LHC}}_\sigma (P)&\displaystyle \simeq \sum \limits _{\{n,m\}} R_{nm}^{h_j+y} \frac{ \sigma _\text {SM}(p p \rightarrow n m \rightarrow H+y) }{\sigma _\text {SM}(p p \rightarrow H+y)}\,, \end{aligned}$$
(10)

to calculate the missing hadronic cross section ratios from the partonic cross section ratios. Here, \(\sigma _\text {SM}(p p \rightarrow n m \rightarrow H+y)\) denotes the contribution from the partonic initial state \(nm\) to the total hadronic cross section for the process \(p p \rightarrow H+y\) in the SM. The hadronic LHC ratios are evaluated separately for \(7\) and \(8\,\, \mathrm {TeV}\) but, as already mentioned, using the same values for \(R_{nm}^{h_j+y}\).

The hadronic cross section ratios for the SM appearing in Eqs. (9), (10) are provided within HiggsBounds. Further discussion of the applicability of this approximation, and details of how the ratios \(\sigma _\text {SM}(p\bar{p} \rightarrow n m \rightarrow H+y)/\sigma _\text {SM}(p\bar{p} \rightarrow H+y)\) are calculated for the Tevatron have been presented in [12]. In HiggsBounds-4, these cross section ratios are provided for the LHC with \(\sqrt{s}=7\,\, \mathrm {TeV}\) and \(\sqrt{s}=8\,\, \mathrm {TeV}\) based on the prediction of the LHC Higgs Cross Section Working Group [3335] for the gluon fusion cross section and the \(HZ\) and \(HW\) cross sections, and the program bbh@nnlo 1.3 [36]. Results on the relative contributions from different parton configurations (as indicated in the figure) to the total hadronic cross section for single Higgs production, \(HZ\) production, \(HW^\pm \) production (at \(\sqrt{s}=8\,\, \mathrm {TeV}\)) are shown in Fig. 3.

Fig. 3
figure 3

Relative SM contributions from different partonic subprocesses to the total hadronic cross sections for \(pp\rightarrow H\) (left), \(pp\rightarrow HZ\) (center), and \(pp\rightarrow HW^\pm \) (right) at \(\sqrt{s}=8\) TeV

3.3 The effective coupling approximation

With the effective coupling option (whichinput=effC), the user input is greatly simplified and reduced to a smaller number of quantities. From this input, physical predictions corresponding to input with the partonic option whichinput=part are calculated. The effective couplings involve specifying (at most)

  1. (1)–(2)

    the same as for whichinput=hadr,

  2. (3)

    normalized scalar and pseudoscalar (squared) effective Higgs couplings to fermions

    $$\begin{aligned}&\left( \frac{g^{\text {model}}_{s,h_k(\text {OP})}}{g^\text {SM}_{H(\text {OP})}}\right) ^2,\, \left( \frac{g^{\text {model}}_{p,h_k(\text {OP})}}{g^\text {SM}_{H(\text {OP})}}\right) ^2,\,\\&\text { with OP } = \{s\bar{s}, c\bar{c}, b\bar{b}, t\bar{t}, \mu ^+ \mu ^-, \tau ^+\tau ^-\}, \end{aligned}$$
  3. (4)

    normalized (squared) effective Higgs couplings to bosons

    $$\begin{aligned}&\left( \frac{g^{\text {model}}_{h_ih_jZ}}{g^{{\text {ref}}}_{H'HZ}}\right) ^2,\, \left( \frac{g^{\text {model}}_{h_k(\text {OP})}}{g^\text {SM}_{H(\text {OP})}}\right) ^2,\,\\&\text { with OP } = \{W^+W^-, ZZ, Z\gamma , \gamma \gamma , gg, ggZ\}, \end{aligned}$$
  4. (5)

    neutral Higgs branching ratios without SM equivalents, charged Higgs branching ratios to SM particles, and top quark branching ratios according to (5)–(7) of whichinput=hadr.

The scalar and pseudoscalar components of a Higgs coupling to a fermion pair are defined in the usual way, via the Feynman rule for the coupling of a generic neutral Higgs boson \(h\) to fermions:

$$\begin{aligned} g_{h f \bar{f}} = \mathrm {i}(g_s \mathbf{1} + g_p \gamma _5), \end{aligned}$$
(11)

where \(g_s\) and \(g_p\) are real-valued scalar and pseudoscalar coupling constants respectively, and \(\mathbf{1}\) and \(\gamma _5\) are the usual matrices in Dirac space. A CP-even scalar, like the SM Higgs boson, has \(g_p=0\) and a CP-odd scalar has \(g_s=0\).

Where these exist, the reference couplings are taken as the SM tree-level equivalents:

$$\begin{aligned}&\left( g^\text {SM}_{HZZ}\right) ^2 = \left( \frac{e M_Z^2}{s_\mathrm {w}M_W}\right) ^2, \end{aligned}$$
(12)
$$\begin{aligned}&\left( g^\text {SM}_{HWW}\right) ^2 = \left( \frac{eM_W}{s_\mathrm {w}}\right) ^2, \end{aligned}$$
(13)
$$\begin{aligned}&\left( g^\text {SM}_{Hf\bar{f}}\right) ^2 = \left( \frac{e m_f}{2s_\mathrm {w}M_W}\right) ^2, \end{aligned}$$
(14)

where \(M_Z\) is the \(Z\) boson mass, \(M_W\) the \(W\) boson mass, and \(m_f\) the mass of fermion \(f\). The reference coupling \((g^{{\text {ref}}}_{H'HZ})^2\), that does not have a SM equivalent, is defined as

$$\begin{aligned} (g^{{\text {ref}}}_{H'HZ})^2= \frac{e^2}{4s_\mathrm {w}^2 c_\mathrm {w}^2}. \end{aligned}$$
(15)

The effective couplings \((g^{\text {model}}_{h_k \gamma \gamma }/g^\text {SM}_{H \gamma \gamma })^2\) (and similarly for \(\gamma Z\)) are loop-induced. They can most conveniently be defined via

$$\begin{aligned} \left( \frac{g^{\text {model}}_{h_k \gamma \gamma }}{g^\text {SM}_{H \gamma \gamma }}\right) ^2 = \frac{\Gamma ^{\text {model}}_{h_k\rightarrow \gamma \gamma }(M_{h_k})}{\Gamma ^\text {SM}_{H\rightarrow \gamma \gamma }(M_H=M_{h_k})}. \end{aligned}$$
(16)

For the Higgs-gluon-gluon effective coupling, \((g^{\text {model}}_{h_kgg}/g^\text {SM}_{Hgg})^2\), there is a choice of definition. It can be defined either via decay widths as

$$\begin{aligned} \left( \frac{g^{\text {model}}_{h_k gg}}{g^\text {SM}_{Hgg}}\right) ^2 =\frac{\Gamma ^{\text {model}}_{h_k\rightarrow gg}(M_{h_k})}{\Gamma ^\text {SM}_{H\rightarrow gg}(M_H=M_{h_k})}, \end{aligned}$$
(17)

or via partonic cross sections:

$$\begin{aligned} \left( \frac{g^{\text {model}}_{h_k gg}}{g^\text {SM}_{Hgg}}\right) ^2=R_{g g}^{h_k}. \end{aligned}$$
(18)

It has to be understood that both these definitions will involve approximations. It is therefore only recommended to use the input of effective couplings when both definitions result in similar values for \((g^{\text {model}}_{h_kgg}/g^\text {SM}_{Hgg})^2\). However, under certain circumstances, this condition can be relaxed. For example, in a model in which the LEP searches for Higgs bosons decaying into hadrons are not relevant, the effective \(h_kgg\) coupling can be defined solely by Eq. (18). Conversely, if for some reason the gluon fusion Higgs production mechanism is not relevant, the effective coupling can be defined solely by Eq. (17).

The calculation of LEP and Tevatron cross section ratios from the effective couplings has been described in [22, 23]. Here we shall focus on the extension of this procedure to LHC cross sections. Partonic cross section ratios for the LHC are calculated as

$$\begin{aligned}&R_{g g}^{h_k}= \left( \frac{g^{\text {model}}_{h_kgg}}{g^\text {SM}_{Hgg}}\right) ^2 \,,\end{aligned}$$
(19)
$$\begin{aligned}&R_{g g}^{h_kZ}=\left( \frac{g^{\text {model}}_{h_kggZ}}{g^\text {SM}_{HggZ}}\right) ^2\,,\end{aligned}$$
(20)
$$\begin{aligned}&R_{b \bar{b} }^{ h_k} = R_{b g,\bar{b} g}^{h_kb,h_k \bar{b}} =\left( \frac{g^{\text {model}}_{s,h_k b\bar{b}}}{g^\text {SM}_{H b\bar{b}}}\right) ^2+ \left( \frac{g^{\text {model}}_{p,h_k b\bar{b}}}{g^\text {SM}_{H b\bar{b}}}\right) ^2 \,, \end{aligned}$$
(21)
$$\begin{aligned}&R_{q \bar{q}' }^{ h_k W^+} = R_{q' \bar{q} }^{ h_k W^-} = \left( \frac{g^{\text {model}}_{h_kWW}}{g^\text {SM}_{HWW}}\right) ^2 \,,\end{aligned}$$
(22)
$$\begin{aligned}&R_{q'' \bar{q}'' }^{ h_k Z} = \left( \frac{g^{\text {model}}_{h_kZZ}}{g^\text {SM}_{HZZ}}\right) ^2 \,, \end{aligned}$$
(23)

where \((q,q')\in \{(u,d),(c,s)\}\) and \(q''\in \{u,d,c,s,b\}\). The normalized hadronic cross section ratio for \(t \bar{t}\) production together with a \({\mathcal {CP}}\)-even Higgs boson for the LHC is obtained using

$$\begin{aligned} R_{\sigma }(p p \rightarrow t \bar{t} h^{{\mathcal {CP}}\mathrm {-even}}_k ) = \left( \frac{g^{\text {model}}_{s,h_k t \bar{t}}}{g^\text {SM}_{H t \bar{t}}}\right) ^2. \end{aligned}$$
(24)

For \(h_k\) production via VBF, the normalized hadronic cross section ratio is calculated using the approximate relation

$$\begin{aligned}&R^{\text {LHC}}_{\sigma }(p p \rightarrow h_k \text { via VBF})\nonumber \\&\quad = R_{\text {VBF,LHC}}^{WW}\, \left( \frac{g^{\text {model}}_{h_kWW}}{g^\text {SM}_{HWW}}\right) ^2+R_{\text {VBF,LHC}}^{ZZ}\, \left( \frac{g^{\text {model}}_{h_kZZ}}{g^\text {SM}_{HZZ}}\right) ^2. \end{aligned}$$
(25)

The \(M_H\) dependence of the relative fractions of VBF events induced from \(WW\) and \(ZZ\) fusion, denoted as \(R_{\text {VBF,LHC}}^{WW}\) and \(R_{\text {VBF,LHC}}^{ZZ}\), respectively, is mild and, for the Tevatron case, can be approximated by constant values. For the LHC case, we obtain these functions by fitting to SM results produced with HAWK 1.1 [37, 38] for \(p p\) collisions at \(\sqrt{s}=7\,\, \mathrm {TeV}\). Note that, for models where \((g^{\text {model}}_{h_kWW}/g^\text {SM}_{HWW})^2=(g^{\text {model}}_{h_kZZ}/g^\text {SM}_{HZZ})^2=(g^{\text {model}}_{h_kVV}/g^\text {SM}_{HVV})^2\) (which is often the case), Eq. (25) reduces simply to

$$\begin{aligned} R^{\text {LHC}}_{\sigma }(p p \rightarrow h_k \text { via VBF}) = \left( \frac{g^{\text {model}}_{h_kVV}}{g^\text {SM}_{HVV}}\right) ^2, \end{aligned}$$
(26)

independent of \(R_{\text {VBF,LHC}}^{WW}\) and \(R_{\text {VBF,LHC}}^{ZZ}\).

The input scheme for decay widths and branching ratios is not affected by the extension to include LHC results in HiggsBounds-4, and the calculation from the effective couplings follows what has been published earlier [12, 13].Footnote 5 The only difference compared to previous versions is that we have updated the SM reference values to agree with those recommended by the LHC Higgs cross section working group [3335] for \(M_H\) between \(80\,\, \mathrm {GeV}\) and \(1\,\, \mathrm {TeV}\).

In order to constrain Higgs bosons with masses below \(\sim 10 \,\, \mathrm {GeV}\), the effective coupling input option is usually not appropriate because the implemented SM reference branching ratios are rather inaccurate for such low masses. It is therefore strongly recommended to use one of the other input formats and enter the branching ratios for such a light Higgs boson directly. It can also be relevant in this mass region to consider constraints from other colliders than LEP, which are not included in HiggsBounds.

3.4 Input using the SUSY LesHouches Accord

For the convenience of users interested in supersymmetric models, an input option using the SUSY Les Houches Accord (SLHA) [22, 23] is now offered. This option (available by setting whichinput=SLHA) uses the calculated decay information from the SLHA file and cross sections approximated through the effective couplings approach. It requires the following input to be read from an SLHA file:

  1. (i)

    The masses for the neutral Higgs bosons \(h_k \; (k=1\ldots n_{H^0})\) and singly charged Higgs bosons \(H_j^{\pm } \; (j=1\ldots n_{H^{\pm }})\),

    $$\begin{aligned} M_{h_k},\, M_{H_j^{\pm }}, \end{aligned}$$
  2. (ii)

    the Higgs total decay widths,

    $$\begin{aligned} \Gamma _{\text {tot}}(h_k),\,\Gamma _{\text {tot}}(H_k^{\pm }), \end{aligned}$$
  3. (iii)

    neutral Higgs branching ratios with SM equivalents

    $$\begin{aligned}&\text {BR}_{\text {model}}(h_k\rightarrow \text {SM}),\\&\text {with SM } = \{s\bar{s} , c\bar{c} , b\bar{b} , \mu ^+ \mu ^- , \tau ^+\tau ^- , W^+W^- , ZZ,\\&\qquad \qquad \qquad \quad Z\gamma , \gamma \gamma , gg\} , \end{aligned}$$

    Note that the decay modes have to be specified as two-body decays, irrespectively of whether they are on- or off-shell.

  4. (iv)

    the neutral Higgs branching ratios without SM equivalents

    $$\begin{aligned}&\text {BR}_{\text {model}}(h_k\rightarrow h_i h_j), \\&\text {BR}_{\text {model}}(h_k\rightarrow \mathrm{invisible}), \end{aligned}$$
  5. (v)

    the charged Higgs branching ratios to SM particles

    $$\begin{aligned} \text {BR}(H_j^{+} \rightarrow {\text {SM}}),\,\mathrm{where}\, \text {SM } = \{c\bar{s}, c\bar{b}, \tau ^+ \nu _{\tau }\}, \end{aligned}$$
  6. (vi)

    the top quark branching ratios

    $$\begin{aligned}&\text {BR}(t \rightarrow W^{+}b),\\&\text {BR}(t \rightarrow H_j^{+}b), \end{aligned}$$
  7. (vii)

    normalized scalar and pseudoscalar (squared) effective Higgs couplings to fermions

    $$\begin{aligned}&\left( \frac{g^{\text {model}}_{s,h_k(\text {OP})}}{g^\text {SM}_{H(\text {OP})}}\right) ^2,\, \left( \frac{g^{\text {model}}_{p,h_k(\text {OP})}}{g^\text {SM}_{H(\text {OP})}}\right) ^2,\\&\quad \text { with OP } = \{b\bar{b}, t\bar{t}, \tau ^+\tau ^-\}, \end{aligned}$$
  8. (viii)

    normalized (squared) effective Higgs couplings to bosons

    $$\begin{aligned}&\left( \frac{g^{\text {model}}_{h_ih_jZ}}{g^{{\text {ref}}}_{H'HZ}}\right) ^2,\, \left( \frac{g^{\text {model}}_{h_k(\text {OP})}}{g^\text {SM}_{H(\text {OP})}}\right) ^2,\,\\&\quad \text { with OP } = \{W^+W^-, ZZ, gg, ggZ\}, \end{aligned}$$

In the SLHA input the effective couplings are only used to calculate the Higgs production cross section ratios (unlike the effective coupling option, where they are also used to calculate the branching ratios). The Higgs decay branching ratios are taken directly from the corresponding decay blocks in the SLHA file. In the case of incomplete input, Higgs masses which are not specified in the SLHA file will be set equal to minus one (such that these Higgs bosons are not tested against any limits), whereas any other input that is not specified will be set equal to zero.

Table 2 lists the PDG codes of the particles that can be considered as neutral Higgs bosons by HiggsBounds. The setting of nHzero determines how many of these are used, starting from the top of Table 2. For example, with nHzero=3, the properties of particles with the PDG numbers 25, 35, and 36 are read from the SLHA file. Note that no \({\mathcal {CP}}\) properties are inferred from the PDG numbers of the neutral Higgs bosons. The invisible Higgs branching ratios are obtained from the branching ratios of Higgs bosons into a weakly-interacting lightest supersymmetric particle (LSP). HiggsBounds finds the weakly-interacting candidate with the lowest mass (considering neutralinos and sneutrinos as candidates) and confirms that this particle is indeed the LSP by comparing its mass against the masses of the charged leptons, the lightest chargino, and the gluino. If the LSP is not neutral, the invisible Higgs branching ratio is set to zero.

Table 2 PDG particle codes for particles that represent neutral Higgs bosons. A number nHzero of these will be considered by HiggsBounds, starting from the top row

To specify the required effective couplings, as described by points (vii) and (viii) above, HiggsBounds requires two extra input blocks which are not part of the normal SLHA. An example of these blocks is shown in Table 3. There are some cases when HiggsBounds is unable to use an SLHA file, including any of the following:

Table 3 Examples of the two new SLHA blocks which are required by HiggsBounds when using the SLHA input option
  • The Block MODSEL indicates that there is \(\mathcal {R}\)-parity violation,

  • either Block SPINFO or DCINFO contains an entry with the label ‘4’ (which is used to indicate an unphysical parameter point),

  • the number of neutral Higgs bosons, nHzero, is greater than 5,

  • the number of charged Higgs bosons, nHplus, is greater than 1.

The settings for nHzero and nHplus are given as input, either as arguments to the subroutine initialize_HiggsBounds or on the command line, and are not read from the SLHA file. If HiggsBounds is unable to use an SLHA file (i.e. if one of the situations listed above applies), it might still be possible to run HiggsBounds with one of the other input options. It is nevertheless recommended that the user investigates and understands the reason behind the SLHA failure before proceeding.

The supersymmetric spectrum calculator SPheno [39, 40] can directly write the HiggsBounds SLHA input blocks to its output SLHA file.Footnote 6 For FeynHiggs [4145] we provide a stand-alone program, HBSLHAinputblocks fromFH, which creates the necessary SLHA blocks from the FeynHiggs output.

4 New developments in HiggsBounds-4

4.1 Applying exclusion limits to arbitrary Higgs models

The aim of HiggsBounds is to apply limits derived in Higgs collider searches to models which have not been directly investigated by the experimental analyses. These models can be arbitrary in the sense that they may contain any number of neutral or (singly) charged Higgs bosons,Footnote 7 or particles which behave like (elementary) Higgs bosons in Higgs collider searches. Examples of the latter include theories with composite Higgs bosons [46] or dilatons [47]. More specifically, the requirements on the theory in order for the results of HiggsBounds to be reliable are the following:

  1. (i)

    The narrow width approximation must be applicable, such that the predictions for each process can be factorized into Higgs production and decay.

  2. (ii)

    The investigated model should not change the signature of the background processes considerably. Usually, new physics models which show strong deviations from the SM in the background processes of Higgs searches are not considered in the literature, since this would often put them in conflict with SM electroweak precision data [48, 49]. Hence, they would most likely not be interesting for HiggsBounds anymore. The presence of such backgrounds would rather correspond to an opportunity for the discovery of physics beyond the SM in other areas.

  3. (iii)

    The investigated model should not significantly change the kinematical distributions of the signal topology \(X\) (e.g. the \(\eta \) and \(p_T\) distributions of the final state particles) from that assumed in the corresponding analysis. For a more detailed discussion of this requirement, see Ref. [12, 13].

The above requirements are typically sufficient to ensure the applicability of model-independent exclusion limits, i.e. limits on a cross section of a certain Higgs signal topology, composed of one production and one decay process. If further model assumptions have been made in the experimental analysis, for instance on the \({\mathcal {CP}}\)-properties of the Higgs boson or on the top quark branching ratios, HiggsBounds checks whether the investigated model fulfills them before applying the analysis.

The application of exclusion limits to arbitrary Higgs models becomes less trivial if the experimental analysis combines several Higgs signal topologies under the assumption of a specific model. This is the case for most of the Tevatron and LHC Higgs searches, where a SM Higgs boson is assumed. The exclusion limit is then set on a common signal scale factor for all considered SM Higgs topologies (also called signal strength modifier), usually denoted by \(\mu \) (sometimes also \(\sigma /\sigma _\text {SM}\) is used). For an analysis considering \(i=1,\ldots , N\) signal topologies (each consisting of a production mode \(P\) and a final state \(F\)), the prediction for this quantity can be computed for Higgs boson \(h_k\) of the investigated model as

$$\begin{aligned} \mu = \frac{\sum _{i=1}^{N} \epsilon _i [\sigma _\mathrm {model} (P(h_k)) \times \text {BR}_\mathrm {model}(h_k\rightarrow F)]_i}{\sum _{i=1}^{N} \epsilon _j [\sigma _\mathrm {SM} (P(H)) \times \text {BR}_\mathrm {SM}(H\rightarrow F)]_j}. \end{aligned}$$
(27)

The channel efficiencies, \(\epsilon _i\), are assumed to be the same for the model and the SM (see requirements (ii) and (iii) above). If these channel efficiencies were published together with the exclusion limit posed by an experimental analysis, the signal strength modifier \(\mu \) could be computed for a given model without further assumptions. However, these efficiencies have so far been made publicly available only in a very few cases.Footnote 8 In HiggsBounds we therefore neglect the channel efficiencies in Eq. (27), leading to an unavoidable model-dependence of the resulting limit, since the calculation of \(\mu \) via Eq. (27) with all \(\epsilon _i \equiv 1\), is strictly speaking only valid if the model predictions for all signal topologies of the analysis contribute to the total signal rate in (approximately) the same proportions as in the SM.

In order to ensure that an analysis is applied only when this last requirement is fulfilled by the model, HiggsBounds performs a SM-likeness test for every Higgs analysis performed under SM assumptions. A test of this kind has been present in all versions of HiggsBounds [12, 13]. However, this test was significantly improved in HiggsBounds-3.8.0 onwards, as described in Ref. [21], and it is this improved version which we describe here. Neglecting the channel efficiencies, the predicted signal strength modifier \(\mu \) given in Eq. (27) can be obtained as \(\mu \approx \sum _{i=1}^N \omega _i c_i\), where

$$\begin{aligned} c_i = \frac{[\sigma _\mathrm {model}(P(h_k)) \text {BR}_\mathrm {model}(h_k\rightarrow F)]_i}{[\sigma _\mathrm {SM}(P(H)) \text {BR}_\mathrm {SM}(H\rightarrow F)]_i} \end{aligned}$$
(28)

and

$$\begin{aligned} \omega _i = \frac{[\sigma _\mathrm {SM}(P(H)) \text {BR}_\mathrm {SM}(H\rightarrow F)]_i}{\sum _{j=1}^N [\sigma _\mathrm {SM}(P(H)) \text {BR}_\mathrm {SM}(H\rightarrow F)]_j} \end{aligned}$$
(29)

are the (SM normalized) channel signal strengths and the SM channel weights, respectively. The requirement that the signal topologies contribute in similar proportions to the total signal rate as in the SM is fulfilled if all channel signal strengths \(c_i\) are similar to the total signal strength modifier \(\mu \) (and thus similar to each other). A possible SM-likeness criterion would therefore be to require

$$\begin{aligned} \Delta \equiv \max _{i}~\left| \frac{\delta c_i}{\mu } \right| < \zeta \end{aligned}$$
(30)

with \(\delta c_i = c_i - \mu \) and \(\zeta \sim \mathcal {O}\,(\mathrm {few}~\%)\), i.e. that the maximal relative deviation of the channel signal strength modifiers from the total signal strength modifier is less than a few percent. In fact, this criterion is very similar to what was used in earlier versions of HiggsBounds. However, this choice was found to be too restrictive in some cases, since it may reject an analysis application which is actually justifiable, leading to overly conservative results. In particular, it is reasonable that channels contributing only a few percent to the total signal rate should be allowed to deviate more from their SM expectations, since their influence on \(\mu \) is subdominant. We therefore introduce the SM channel weights \(\omega _i\) in an improved SM-likeness test criterion,

$$\begin{aligned} \Delta \equiv \max _{i}~\omega _i \left| \frac{\delta c_i}{\mu } \right| < \zeta . \end{aligned}$$
(31)

The default setting in HiggsBounds is \(\zeta = 2\,\%\). This is a conservative choice, considering that the uncertainties on the rate predictions for individual channels (even in the SM) are generally larger. With the improved SM-likeness test, the maximal weighted deviation of an individual signal strength modifier from the total signal strength modifier is required to be less than \(2\,\%\). Models fulfilling this SM-likeness test for a SM analysis can be safely tested against its exclusion limit.

To illustrate the inclusion of the SM weights \(\omega _i\) in the SM-likeness test criterion, we consider as an example the ATLAS \(H\rightarrow \gamma \gamma \) search [50] and test a model with a single Higgs boson with mass \(m=125\,\, \mathrm {GeV}\). We depart from the SM by modifying either the squared effective Higgs coupling ratio to gluons (normalized to the SM), \(g_{Hgg}^2\), or the coupling to vector bosons, \(g_{HVV}^2\) (\(V=W,Z\)). All other effective Higgs couplings, in particular the \(H\gamma \gamma \) coupling, are set to their SM values. At \(m=125\,\, \mathrm {GeV}\), the SM weights for the LHC at \(\sqrt{s}=7\,\, \mathrm {TeV}\) are

$$\begin{aligned} \omega \approx (87.7\,\%,~6.8\,\%,~3.2\,\%,~1.8\,\%,~0.5\,\%) \end{aligned}$$
(32)

for the production processes \((gg\rightarrow H,~\mathrm {VBF},~HW,~HZ, Ht\bar{t})\). In Fig. 4 we show how the total signal strength modifier \(\mu \) and the \(c_i\) for the signal topologies are influenced by the modified effective Higgs couplings. Varying \(g^2_{Hgg}\) influences only the \(gg\rightarrow H\) (ggf) cross section. However, due to its large SM weight, \(\omega _\mathrm {ggf} \approx 87.7\,\%\), the total signal strength modifier \(\mu \) follows \(c (\mathrm {ggf})\) closely. The failure of the SM-likeness test at \(g^2_{Hgg} = 0.835\) and \(1.225\) is therefore eventually caused by the ggf signal topology, although the deviation \(\delta c_i\) for the remaining signal topologies VBF, \(HW\), \(HZ\) and \(Ht\bar{t}\) is much larger here. However, the SM weights of these channels are much smaller. The same effects can be seen when varying \(g_{HVV}^2\) (\(V=W,Z\)). Now, the \(c_i\) of the VBF, \(HW\), \(HZ\) signal topologies are affected by the modified effective coupling, but the total signal strength modifier \(\mu \) is only slightly influenced due the small weight of these channels. Again, the deviation between \(\mu \) and \(c(\mathrm {ggf})\) eventually causes the SM-likeness test to fail. Due to the inclusion of the SM weights in Eq. (31), subdominant signal topologies are allowed to deviate further from \(\mu \).

Fig. 4
figure 4

Performance of the SM-likeness test. Total signal strength modifier \(\mu \) and the relevant individual signal strength modifiers \(c_i\) for the ATLAS \(H\rightarrow \gamma \gamma \) search [50] with modified effective Higgs couplings (relative to the SM) \(g^2_{Hgg}\) (left) and \(g_{HVV}^2\) (\(V=W,Z\)) (right) for a Higgs boson with mass \(m=125\,\, \mathrm {GeV}\). The gray regions indicate the parameters for which the SM-likeness test fails

In comparison with the old SM-likeness test (which was used in HiggsBounds up to version 3.7.0), the new criterion leads to a wider applicability of SM Higgs search results to other Higgs sectors, and thus to a significant improvement of the performance of HiggsBounds. This is shown in Fig. 5 for the \(M_h^\mathrm{mod+}\) benchmark scenario of the MSSM [24]. Without SM weights (left panel), the LHC exclusion approximately follows the results from the dedicated MSSM search for \(H/A\rightarrow \tau \tau \) [51], and no additional exclusion can be set. In particular there is no LHC exclusion from the SM-like Higgs boson, \(h\). With the full weighted criterion active (the default setting in HiggsBounds-4), the lightest MSSM Higgs boson can become sufficiently SM-like at large \(M_A\) and small \(\tan \beta \) for the combined SM Higgs searches of ATLAS and CMS to be applied, giving additional areas of exclusion (right panel). Exclusion by individual Higgs bosons for the same scenario can be seen in Fig. 2, which has also been produced using the weighted SM-likeness criterion.

Fig. 5
figure 5

Combined HiggsBounds exclusion in the the \(M_h^\mathrm{mod+}\) benchmark scenario of the MSSM using a SM-likeness test without weights, Eq. (30) (left), and using the new SM-likeness test with weights included, Eq. (31) (right). The colour coding is the same as in Fig. 1

4.2 Including Higgs mass uncertainties

In several theories with extended Higgs sectors, the Higgs boson masses are not free parameters but can be predicted as a function of the other model parameters up to a certain theoretical accuracy. This is the case, for example, in the MSSM where out of the five physical Higgs states typically only one mass, \(M_A\) or \(M_{H^\pm }\), is used as an (on-shell) input parameter. The remaining Higgs masses then become predictions of the model, with a theoretical uncertainty that varies within the parameter space and with the sophistication of the theoretical prediction.

We have extended HiggsBounds to be able to take this type of theoretical uncertainty into account when evaluating the Higgs exclusion limits. For theories that have no Higgs mass uncertainties, or where they are negligibly small, this new feature can be left deactivated. In HiggsBounds-4, the Higgs mass uncertainties are taken into account approximately by varying each mass within a user-defined interval.Footnote 9 If the tested Higgs boson is unexcluded by the probed limit (in the normal HiggsBounds sense) for any mass in this interval, the tested parameter point of the model is regarded as being unexcluded. This leads to an overall conservative (weaker) result for the exclusion limit when uncertainties are included.

Technically, the number of mass points \(N\) considered in the variation can be specified by a variable. The default setting is \(N=3\) (this corresponds to testing the central mass value, \(M_H\), and the two endpoints, \(M_H\pm \Delta M_H\), of the specified uncertainty interval). When a sensitive limit varies rapidly with \(M_H\), it is advisable to increase \(N\) for best results. The mass variation is performed for each Higgs boson independently. In the classic HiggsBounds method this variation is also simultaneous, which leads to a multi-dimensional computation grid with a complexity growing as \(\mathcal {O}(N^{n_H})\), where \(n_H\) is the number of Higgs bosons with a non-zero mass uncertainty.Footnote 10 For the full method, since the limit from each Higgs boson is already considered independently of the others, the complexity remains managable, i.e. \(\mathcal {O}(n_H N )\). Nevertheless, the user is strongly encouraged to specify uncertainties only for those Higgs bosons where this is numerically relevant.

The effects of a theoretical mass uncertainty on the resulting HiggsBounds limits are demonstrated in Fig. 6, which shows the combined exclusion for a SM-like Higgs boson with \(\Delta M_H=0\) GeV (solid lines), and similarly for a Higgs boson with SM-like couplings but a theory mass uncertainty of \(\Delta M_H=2\) GeV (dashed lines). In this figure, the mass range excluded at \(95\,\%\) C.L. corresponds to where the limit on \(\sigma _{\mathrm {model}}/\sigma _{\mathrm {SM}}<1\). Including the mass uncertainties can be seen here as a broadening of the allowed range for the Higgs mass prediction in the model by \(\pm 2\) GeV around the signal region. It can also be seen that for a given mass point the resulting upper limit on the signal rate is always weaker or equal to the upper limit obtained without theoretical mass uncertainty. Including a theory mass uncertainty in HiggsBounds therefore produces overall more conservative limits, which is as expected.

Fig. 6
figure 6

Upper limits from HiggsBounds (at \(95\,\%\) C.L.) on the relative signal strength versus the Higgs boson mass in the SM, which has zero theoretical uncertainty (dashed lines), and in a model with a SM-like Higgs boson with a theoretical mass uncertainty of \(2\) GeV (solid lines). The two colours indicate mass ranges where the most sensitive limit is from either LEP (black) or the LHC (red)

This point is further illustrated in Fig. 7, which shows the resulting limits from the light (SM-like) MSSM Higgs boson, \(h\), when running HiggsBounds full for the \(M_h^\mathrm{max}\) benchmark scenario [24]. Similar to above, the green band shows the region of parameter space (in this scenario) where \(M_h=125.7\pm 2(3)\,\, \mathrm {GeV}\). For large values of \(M_A\) and \(\tan \beta \,\!\!,\) the \(M_h^\mathrm{max}\) scenario gives rise to values of \(M_h\) that are too high compared to the measured LHC signal. The predicted value for \(M_h\) increases with \(\tan \beta \). \(M_h\gtrsim 128\,\, \mathrm {GeV}\) is excluded when no theory uncertainty is applied (cf. Fig. 6). The three panels in Fig. 7 show, from left to right, the results when using a mass uncertainty (resulting from the calculation of \(M_h\) in the model [45]) of \(\Delta M_h=0\) GeV, \(\Delta M_h=1\) GeV, and \(\Delta M_h=2\) GeV. It can be seen that the exclusion at high \(M_A\) from the limit on the lightest Higgs boson goes down to lower \(\tan \beta \,\) values when \(\Delta M_h\) is small. This illustrates the importance of taking Higgs mass uncertainties into account when interpreting exclusion limits (and compatibility with observed signals, see [15]) in the MSSM and other scenarios for physics beyond the SM.

Fig. 7
figure 7

Contribution from the lightest MSSM Higgs boson, \(h\), to the full HiggsBounds exclusion in the parameter space for the \(M_h^\mathrm{max}\) benchmark scenario [24]. The results are shown for a theory mass uncertainty of \(\Delta M_h=0\) GeV (left), \(\Delta M_h=1\) GeV (center), and \(\Delta M_h=2\) GeV (right). The colour coding is the same as used in Fig. 1

4.3 LEP \(\chi ^2\) extension

An unfortunate limitation of both the model-independent limits implemented in HiggsBounds, as well as the model-dependent search limits presented by the experiments, is that they are available only at one fixed confidence interval, which is \(95\,\%\) C.L. for all searches implemented here. The result of an exclusion provided by HiggsBounds based on these searches therefore has a confidence limit of at least \(95\,\%\) C.L., and in many cases higher. However, the exact level of confidence at which a signal with the properties given to HiggsBounds is either excluded or allowed, is generally unknown.

This has unfortunate consequences for the use of these limits in applications like global fits (see e.g. Refs. [5254] for examples of such fits in the MSSM). There, a model point where the predicted Higgs signal is excluded for example at \(96\,\%\) C.L., i.e. with a significance of slightly more than \(2\,\sigma \), might still be a very good fit if the other properties of the model point in the global fit match the data well. However, the conventional HiggsBounds output only contains information about whether the parameter point is experimentally excluded at at least \(95\,\%\) C.L. and thus can only be treated as a “hard cut” on the validity of a parameter point.

In order to circumvent this problem, at least for the LEP Higgs searches, the full information on \(CL_{s+b}\) and \(CL_s\) for all Higgs mass combinations in the model-independent LEP searches from [17] have been re-calculated for varying cross sections [55]. These can be written as \(\sigma _i=f_i\, \sigma _{i, \mathrm {ref}}\), where \(\sigma _{i, \mathrm {ref}}\) is the reference cross section times branching fraction for search \(i\), motivated by the SM Higgs boson or the corresponding cross section for non-SM-Higgs bosons (see [56] for details), and \(f_i\) is an arbitrary scaling parameter. A logarithmic grid in the scaling parameters \(f_i\) with \(100\) points between \(10^{-3}\) and \(1\) is used. Using an interpolation, the actual \(CL\) can be calculated for every Higgs production mode at LEP for every physically allowed cross section.

This \(CL\) can then be transferred into a quantity whose properties closely follow that of a \(\chi ^2\) function. This is achieved by assuming that the distribution of \(-2\ln Q\) [56] is Gaussian in the asymptotic limit. Transferring the one-sided \(CL\) into the two-sided calculation of a \(\chi ^2\), the following formula can then be used

$$\begin{aligned} \chi ^2_{H} = 2\,\mathrm {InvErf}^{2}(1-2CL_{s+b}). \end{aligned}$$

The resulting \(\chi ^2_{H}\) can be used as a continuous expression of the agreement between the result of the LEP Higgs boson searches and the model predictions. Note that, in the case of a strong excess in one of the searches, \(\chi ^2_{H}\) is not only large for models whose predicted cross section times branching fraction is above the observed limit, but also for predictions much smaller than the observed rate in data.

In cases when the predicted cross section is lower than the minimal (rescaled) value available in the table, the corresponding \(\chi ^2\) value is set to zero. When the predicted cross section exceeds the tabulated values, no reliable \(\chi ^2\) value can be calculated, and the value \(\chi ^2=-999\) is returned to indicate that a problem has occured. This default behavior can be changed (by setting a flag in usefulbits.f90) to use instead the \(\chi ^2\) value for the maximal (rescaled) cross section available for that combination of Higgs masses.

An example of the relation between the LEP \(CL_{s+b}\) and \(\chi ^2_{H}\), also for different values of \(f\), is given in Fig. 8. It can be seen that for \(CL_{s+b}\approx 0.5\), indicating very good agreement of the signal plus background prediction with the data, fluctuations of \(\chi ^2_H\) around 0 are unavoidable, but numerically irrelevant. In addition, the possibility exists to follow a prescription from [57] to include a mass uncertainty into \(\chi ^2_H\) by folding the full \(\chi ^2\) distribution with a gaussian \(G_{\Delta M_H}\) with a mass uncertainty \(\Delta M_H\) given by the user,Footnote 11 instead of evaluating \(\chi ^2_H\) just at the given \(M_H\):

$$\begin{aligned}&\chi ^2_{H,\mathrm {bare}}(M_H,\Delta M_H)\\&\quad =-2\ln \left( \,\,\int _{-\infty }^{+\infty }e^{-\frac{1}{2}\chi ^2_H(M')}\, G_{\Delta M_H}(M_H-M')\mathrm {d}M'\right) . \end{aligned}$$

Since the folding introduces small, but non-zero values \(\chi ^2_{H,\mathrm {bare}}(M_H,\Delta M_H)\) for \(M_H>116.4\) GeV, where no sensitivity is expected for the SM-like Higgs search channels, the final \(\chi ^2_H(M_H,\Delta M_H)\) is obtained by subtracting \(\chi ^2_H(116.4\,\mathrm {GeV},\Delta M_H)\) from \(\chi ^2_{H,\mathrm {bare}}(M_H,\Delta M_H)\) for \(M_H\le 116.4\) GeV, and by setting \(\chi ^2_H(M_H,\Delta M_H)=0\) above. A similar procedure is adapted for non-SM like searches, where the point of vanishing sensitivity is determined for each search prior to the folding.

Fig. 8
figure 8

Examples for transferring the LEP \(CL\) into a value for \(\chi ^2_{H}\) using three different values of the scale factor \(f=(0.25, 0.5, 1.0)\). The upper row shows a the LEP \(CL_{s+b}\) result for \(e^+e^-\rightarrow hZ\) in the SM, and b the corresponding \(\chi ^2_H\) (note the logarithmic scale). In the lower row, similar results are shown for \(e^+e^-\rightarrow hA\rightarrow b\bar{b}b\bar{b}\), with c again being the \(CL_{s+b}\) values and d the result for \(\chi ^2_{H}\)

This implementation has already been used in global fits of constrained SUSY models [5254]. A non-trivial example of how the LEP \(\chi ^2\) information can be applied is given in Fig. 9. This figure shows the MSSM low-\(M_H\) benchmark scenario [24], where the heavier \({\mathcal {CP}}\)-even Higgs boson is interpreted as the LHC signal around \(M_H\sim 126\) GeV. In that case, the lightest Higgs boson, \(h\), is usually below the SM LEP limit and has suppressed couplings to gauge bosons. This is reflected in the figure, where a sizeable \(\chi ^2\) penalty can be seen to result in parts of the parameter space, corresponding to regions of low \(M_h\) (an uncertainty of \(2\,\, \mathrm {GeV}\) was used here) where the couplings to gauge bosons is such that the LEP Higgs searches are sensitive to the production of such a state. The sharp edge in the \(\chi ^2\) distribution in Fig. 9 is obtained at the boundary between two regions of parameter space where the \(\chi ^2\) contribution comes from the channels \(e^+e^-\rightarrow hZ\), \(h\rightarrow b\bar{b}\) and \(e^+e^-\rightarrow hA \rightarrow 4b\), respectively. Using the LEP \(\chi ^2\) information together with the HiggsBounds exclusion at \(95\,\%\) C.L. from Tevatron/LHC, Fig. 9 (right) gives the most complete information available from direct Higgs search limits.

Fig. 9
figure 9

HiggsBounds results for the LEP \(\chi ^2\) (colours) in the low-\(M_H\) scenario of the MSSM [24]. The LEP \(\chi ^2\) information is shown both on its own (left), and with with the combined LHC exclusion bounds in gray (right). The latest limit from ATLAS charged Higgs searches [58] is not applied here. These results, which are included from HiggsBounds-4.1, lead to exclusion over the whole parameter space of this benchmark scenario

Even after the discovery of a SM-like Higgs boson, Higgs boson exclusions still plays, and will continue to play, a vital role in fitting models of physics beyond the SM with an extended Higgs sector. It would therefore be to great advantage if the Tevatron and LHC collaborations could follow the example of the LEP Higgs WG and provide exclusion limits for varying values of \(f\,\sigma _{\mathrm {ref}}\) in addition to the results that are presented at \(95\,\%\) C.L..

5 User operating instructions

In this section we describe in detail the two main methods to use HiggsBounds-4: The command line version and the library of subroutines. There is also an online version that provides quick access to all the functionality of HiggsBounds, without the need to install the code.

5.1 Installation

The HiggsBounds source code, the online version, and the documentation can all be obtained at the URL http://higgsbounds.hepforge.org The HiggsBounds-4 code is mostly written in Fortran 90 but includes also a few Fortran 2003 features. It has been tested with a variety of Fortran compilers, including the free GNU compilerFootnote 12 (gfortran) which accompanies most Linux distributions.

Before compiling the HiggsBounds code, the user should first make changes to the configure script to appropriately reflect the compiler and path settings on the user’s system. The code can then be compiled by running

$$\begin{aligned}&\mathtt{{./configure}}&\\&\mathtt{{make}}&\end{aligned}$$

which creates the main HiggsBounds executable and the library of subroutines, libHB.a. Any program for which the HiggsBounds subroutines should be used can be compiled and linked to the library by adding \(\mathtt{{-L<HBpath> -lHB}}\) to the command line, for example,

$$\begin{aligned}&\mathtt{{gfortran}}\ \mathtt{{myprog.f90 -o}}\ \mathtt{{myprog}}\\&\quad \mathtt{{-L<HBpath> -lHB}} \end{aligned}$$

where \(\mathtt{{<\!\!HBpath\!\!>}}\) is the location of the HiggsBounds library. The HiggsBounds subroutines make use of the Fortran file handles 10, 11, 44, 45 and 87, which means that users should avoid these file handles in programs calling the subroutines.

The default behavior of HiggsBounds-4 is to use the full (new) method to generate combined exclusion. To set the classic method as the default, the user can modify the flag run_HB_classic in the file usefulbits.f90 before compiling HiggsBounds. When running the subroutine version of HiggsBounds, it is also possible to access the results from both methods without changing the default behavior, see below.

The library of subroutines and the command-line version share a common set of features, which we will describe first. We will then give the proper operating instructions to use each of these HiggsBounds formats individually.

5.2 Common features: Input

Regardless of the operation mode, HiggsBounds requires five basic types of user input:

  1. (i)

    the number of neutral Higgs bosons in the model under study (nHzero)

  2. (ii)

    the number of (singly, positively) charged Higgs bosons in the model (nHplus)

  3. (iii)

    the set of experimental analyses which should be considered (whichanalyses)

  4. (iv)

    the desired input format for the theoretical predictions (whichinput)

  5. (v)

    the theoretical predictions for the Higgs sector of the model (given as arrays)

The variables nHzero and nHplus are currently both limited to the range \(0\)\(9\), but if necessary this could easily be extended in the future. The possible values for the choice of experimental analyses (whichanalyses) are described in Table 4.

Table 4 Possible settings for the variable whichanalyses, which determines the subset of experimental analyses to be considered by HiggsBounds

HiggsBounds expects the theoretical input to be provided in one of three formats labelled by the variable whichinput. These formats are described in detail in Sect. 3, and their required input is briefly summarized in Table 1. In Appendix 8 (Tables 11, 12, 13, 14) we assign names and list the full contents of all the possible input arrays for the theory predictions. These names will be used below to describe the input requirements of each version of HiggsBounds individually.

5.3 Common features: Output

HiggsBounds provides the user with four types of basic output:

  1. (i)

    whether or not the model under study is excluded by Higgs searches at the 95 % C.L. (HBresult)

  2. (ii)

    the reference number of the analysis application (\(X_0\)) with the highest statistical sensitivity (chan)

  3. (iii)

    the number of Higgs bosons that contributed to the theoretical rate for the corresponding process (ncombined)

  4. (iv)

    the ratio \(k_0=Q_\mathrm{model}/Q_\mathrm{obs}\) for the process with highest statistical sensitivity (obsratio).

As discussed in Sect. 2, the extended HiggsBounds algorithm now offers similar quantities to be calculated individually for each Higgs boson in the model. When making use of the full method, the corresponding output quantities are promoted to arrays of length \(n+1\), where \(n\) is the total number of (neutral and charged) Higgs bosons in the model. The combined result (contained in element \(0\) of these arrays) from this extended test can also be used in analogy to the result of HiggsBounds classic. When several Higgs bosons exclude the same point through different searches, the values for chan, obsratio, and ncombined in the combined result refers to the channel giving the strongest exclusion.

Table 5 shows how to interpret the possible values of HBresult and obsratio (one entry in the case of arrays), which are complementary. When using either the library of subroutines or the command-line version, the keys associating the reference numbers (as given by chan) with the analysis applications is written in human-readable format in the file Key.dat. In the online version, this information appears directly on the screen. When the SLHA option is used for input, the HiggsBounds results can be added to SLHA files in the form of a new block, called HiggsBoundsResults. An example of this block is shown in Table 6.

Table 5 The possible values of the output variables HBresult and obsratio, which indicate whether a parameter point has been excluded at the 95 % C.L. by the experimental results under consideration
Table 6 Example of the output SLHA Block HiggsBoundsResults. Strings appearing in the output are wrapped with ‘\(||\)

5.4 Library of subroutines

In this section we list all the user subroutines available through the HiggsBounds library.

figure a

In each run, this subroutine must be called before any other subroutine of the HiggsBounds package, and it must be called only once. It performs preparatory operations such as initialization of arrays and reading in the tables of experimental data. If the neutral Higgs sector should not be tested with HiggsBounds, the user should set nHzero=0. Similarly, if the user does not wish to test the charged Higgs sector, set nHplus=0. The possible settings for whichanalyses are shown in Table 4.

figure b

This is an alternative version of initialize_HiggsBounds, which takes an integer for the third argument instead of a string constant. This code specifies the set of experimental data that is used by HiggsBounds according to the first column of Table 4.

figure c

This subroutine sets the model input for the neutral Higgs sector using the effective couplings (whichinput=effC), as defined in Sect. 3.3. Using this method also excludes the use of either parton-level or hadron-level input. The meaning of the input arrays (all of length \(n=\mathtt{nHzero}>0\)) is summarized in Appendix 8, Table 11. If any of the effective couplings are deemed to be irrelevant, the corresponding array may be filled with zeros. However, this needs to be exercised with some caution for quantities relevant in a SM-like Higgs search (as most of the limits reported from the LHC are). It is possible that setting certain couplings artificially to zero could lead to the model failing the SM-likeness test, cf. Sect. 4.1.

figure d

This routine is used to set the input for the neutral Higgs sector using parton-level cross sections (whichinput=part), as defined in Sect. 3.2. Using this method excludes the simultaneous use of effective couplings or hadron-level input. The meaning of the input arrays (of length \(n=\mathtt{nHzero}>0\)) are summarized in more detail in Appendix 8, Table 12 (cross sections) and Table 13 (branching ratios). As for the effective coupling case, quantities which are not required by any channel that has a competitive sensitivity can be set to zero to simplify the input (and the same caveats about searches for SM-like Higgs bosons apply).

figure e

This subroutine sets the input for the neutral Higgs sector using hadron-level cross sections (whichinput=hadr), as defined in Sect. 3.1. Using this method excludes the use of effective couplings or parton-level input. The names for the input arrays (of length \(n=\mathtt{nHzero}>0\)) are described in Appendix 8, Table 12 (cross sections) and Table 13 (branching ratios). Similarly to the effective coupling case, quantities which are not required by any channel that has a competitive sensitivity can be set to zero to simplify the input (and the same caveats about searches for SM-like Higgs bosons apply).

figure f

The subroutine HiggsBounds_charged_input gives the charged Higgs sector input to HiggsBounds. The use of this subroutine is only required if \(k=\mathtt{nHplus}\) is non-zero (recall that nHplus is set in subroutine initialize_HiggsBounds). Currently, only results from searches for light charged Higgs bosons (\(M_{H^\pm }< m_t\)) are available. Once results from heavy charged Higgs searches are presented, this interface will be extended with input of the necessary cross sections. The names used for the input arrays are described in Appendix 8, Table 13.

figure g

This subroutine can be used for supersymmetric theories as an alternative to the other routines which provide model input to HiggsBounds. It is called with a string-type variable, SLHAfilename, which gives the name of an SLHA file (full path should be included if not in the current working directory). The model predictions are then read in from this file, which should contain the two HiggsBounds-specific blocks as described in Sect. 3.4. Furthermore, it will set the mass uncertainties of the neutral and charged Higgs bosons according to the values given in the SLHA block DMASS (if available).

figure h

This subroutine allows the user to specify theory mass uncertainties for the neutral and charged Higgs bosons of the model. The implementation and use of these uncertainties when the exclusion limits are evaluated is discussed in detail in Sect. 4.2. The default is for all the uncertainties to be zero. The treatment of mass uncertainties in the limit setting is invoked automatically by setting any of them to a non-zero value. The routine takes two arrays as arguments: dMhneut(n) (of length \(n=1\ldots \) nHzero), which specifies the (absolute) uncertainties for the neutral Higgs boson masses in GeV, and dMhch(k) (\(k=1\ldots \) nHplus) which does the same for the charged Higgs bosons. If either \(\mathtt{nHzero}=0\) or \(\mathtt{nHplus}=0\), the corresponding uncertainty array will not be used (and can therefore be set to arbitrary values).

figure i

After initializating and setting the model input using one of the methods discussed above, this subroutine is called to perform the main part of the HiggsBounds calculations. The results from the run is given as output. The combined result, HBresult, is reported according to the description in Table 5. The channel with the highest exclusion power is identified by its code, chan (the channel codes are written to the file Key.dat), and the corresponding ratio of the model prediction to the observed limit in this channel is given by obsratio. Finally, the number of Higgs bosons combined in this prediction is ncombined. The default behavior of this subroutine (which can be controlled by setting a flag in usefulbits.f90) is to use the full exclusion method of HiggsBounds, rather than the classic method employed in previous versions.

figure j

This subroutine can be used to run HiggsBounds directly in the classic mode, without changing any flag. As discussed in Sect. 2, the classic method tests for exclusion using only the globally most sensitive analysis (considering all the Higgs boson). This corresponds to the behavior of HiggsBounds prior to version 4. The output variables have the same definitions as for run_HiggsBounds.

figure k

This subroutine runs HiggsBounds in the full mode. This is similar to the default behavior of run_HiggsBounds, but with the important difference that when running this subroutine the results from each individual Higgs boson can be accessed in the output. Each of the output variables is therefore an array (with elements \(n=0\ldots \) nHzero+nHplus), where element \(0\) contains the combined result (the same as obtained from run_HiggsBounds) and the remaining entries hold the individual results: first the entries for the neutral Higgs bosons, followed by the results for the charged Higgs bosons.

figure l

This subroutine produces the results for a single Higgs boson, which should be indexed by h. The indexing is such that the neutral Higgs bosons correspond to \(\mathtt{h}=1\ldots \mathtt{nHzero}\), followed by the charged Higgs bosons of the model, \(\mathtt{h}=\mathtt{nHzero}+1\ldots \mathtt{nHzero+nHplus}\). To get the results for more than one individual Higgs boson, it is recommended to instead use the subroutine run_HiggsBounds_full for better performance.

figure m

When using the SLHA input, the subroutine HiggsBounds_SLHA_output can be called after using (any of the different) run_HiggsBounds routines in order to write the block HiggsBoundsResults to the SLHA file. The results are written in terms of the combined exclusion, see Table 6 for an example.

figure n

The subroutine finish_HiggsBounds should be called once at the end of the program, after all other HiggsBounds subroutines. This deallocates the allocatable arrays used within HiggsBounds.

5.5 Command-line version

When using HiggsBounds from the command line, the run options are specified in the program call, which should be of the form

$$\begin{aligned}&./\,\, \mathtt {HiggsBounds <\!whichanalyses\!> <\!whichinput\!>}\\&\quad \mathtt {<\!nHzero\!> <\!nHplus\!> <\!prefix\!>} \end{aligned}$$

The variable \(<\!\!\!\mathtt{prefix}\!\!\!>\) is a string which is added to the front of input and output file names. It may include directory names or other information identifying the run files. If \(\mathtt{whichinput=SLHA}\), \(\mathtt{<\!prefix\!>}\) should contain the full name of the SLHA file to use, including the path if it is not in the current working directory. When running HiggsBounds from the command line, the program behaviour (full/classic) is determined by a flag specified in the file usefulbits.f90 (the same as for the subroutine run_HiggsBounds). The default setting is that the full method is used.

5.5.1 Input file format

The arrays containing the theoretical model predictions are read in from text files, with each value given in a separate column (separated by whitespace). The contents of each input file is described in Tables 7 and 8. Note that all these files will not be necessary at the same time. This will be specified below. Each row in the input files starts with a line number, \(\mathtt{k}\), which identifies predictions belonging to the same parameter point in different files. The input files must not contain any comments or blank lines. Care should be taken with the order of the array elements in the files. The elements of a one-dimensional array, e.g. Mh for \(\mathtt{nHzero}=3\), is given in the order

Table 7 File names and data format for the contents of HiggsBounds input files (part I). The right column shows the order of the input data arrays within one row of the input file (k is the line number). For the order of elements within the arrays, see the text for details

Mh(1), Mh(2), Mh(3).

The correspondence between the array elements and the physical input quantities is clarified in Appendix B. Not all of the elements of the two dimensional arrays are required. From the array g2hjhiZ only the lower left triangle (including the diagonal) is required (and similarly for lepCS_hjhi_ratio below), since they are symmetric matrices. From the general matrix

$$\begin{aligned} \begin{pmatrix} \mathtt{g2hjhiZ(1,1)} &{}\quad \mathtt{g2hjhiZ(1,2)} &{}\quad \mathtt{g2hjhiZ(1,3)}\\ \mathtt{g2hjhiZ(2,1)} &{}\quad \mathtt{g2hjhiZ(2,2)} &{}\quad \mathtt{g2hjhiZ(2,3)}\\ \mathtt{g2hjhiZ(3,1)} &{}\quad \mathtt{g2hjhiZ(3,2)} &{}\quad \mathtt{g2hjhiZ(3,3)}\\ \end{pmatrix}\!\!, \end{aligned}$$

the required elements should be written in the input file using the order

g2hjhiZ(1,1), g2hjhiZ(2,1), g2hjhiZ(2,2),

g2hjhiZ(3,1), g2hjhiZ(3,2), g2hjhiZ(3,3).

The branching ratios for the Higgs decays to lighter Higgs bosons, \(h_j \rightarrow h_i h_i\), are given via the matrix BR_hjhihi(j,i):

$$\begin{aligned} \begin{pmatrix} \mathtt{BR\_hjhihi(1,1)} &{}\;\; \mathtt{BR\_hjhihi(1,2)}&{}\;\; \mathtt{BR\_hjhihi(1,3)}\\ \mathtt{BR\_hjhihi(2,1)}&{}\;\; \mathtt{BR\_hjhihi(2,2)} &{}\;\; \mathtt{BR\_hjhihi(2,3)}\\ \mathtt{BR\_hjhihi(3,1)}&{}\;\; \mathtt{BR\_hjhihi(3,2)}&{}\;\; \mathtt{BR\_hjhihi(3,3)}\\ \end{pmatrix} \end{aligned}$$

Here, only the off-diagonal components are required since the diagonal elements are not physical quantities. The required elements should be given in the order

BR_hjhihi(1,2), BR_hjhihi(1,3),

BR_hjhihi(2,1), BR_hjhihi(2,3),

BR_hjhihi(3,1), BR_hjhihi(3,2).

Table 8 File names and data format for the contents of HiggsBounds input files (part II). The right column shows the order of the input data arrays within one row of the input file (k is the line number). For the order of elements within the arrays, see the text for details. Note that several arrays appear in two different input files. These files are never required simultaneously in one run of HiggsBounds

The file MHall_uncertainties.dat is optional. If it is provided, HiggsBounds will automatically include the theoretical Higgs mass uncertainties given in the file. This has been described in more detail in Sect. 4.2. The file additional.dat is also listed as optional. If this file is included, it can have any number of columns greater than \(1\) (as for the previous files, the first entry on each line should still be the line number). This input is particularly useful to keep track of variables which are not required by HiggsBounds, but which are helpful when plotting the results from a parameter scan. For example, in the case of the MSSM, additional.dat could be used to store the values of \(\tan \beta \).

As in the subroutine version, the command line version of HiggsBounds expects a subset of the total list of possible input arrays, which depends on the chosen setting of whichinput. The maximal list of files used for each value of whichinput is given in Table 9. Furthermore, some of the arrays will not be relevant for some of the choices for whichanalyses. The command line version of HiggsBounds will consider the list of input files appropriate to the settings of whichinput and whichanalyses, and then attempt to read only those input files where at least one of the arrays contained in the file will be used. Table 10 contains a list of which input files are actually relevant for each setting of whichanalyses. Finally, the model predictions for the neutral and charged Higgs sectors are independently specified in different input files (except for the common optional input file MHall_uncertainties.dat). Therefore, the files Mhplus_GammaTot.dat, LEP_HpHm_CS_ratios.dat, BR_t.dat and BR_Hplus.dat are only required if the user wants to test the charged Higgs sector (nHplus > 0). On the other hand, if the user is only interested in the constraints from charged Higgs boson searches, it is sufficient to give only these files while setting nHzero = 0. As for the subroutine version, if the user does not require processes involving a particular branching ratio or cross section ratio to be checked by HiggsBounds, that particular array can be filled with zeros.

Table 9 List of possible input files for each setting of whichinput. Optional files are marked with (*). The files required can also depend on the setting of whichanalyses, see Table 10
Table 10 List of input files relevant to each setting of whichanalyses (marked by ‘y’). The required files also depend on the settings of whichinput, nHzero and nHplus; see Table 9 and the text for details

For supersymmetric models, one possible way of generating HiggsBounds input files is to use the model building tool SARAH [5961] in conjunction with the spectrum generator SPheno [39, 40]. These codes can directly write out the HiggsBounds input files required for the effective coupling approximation.

5.5.2 Output file format

When the command line version of HiggsBounds is used with whichinput=hadr, part or effC, the output is written to the file \(<\!\!\!\mathtt{prefix}\!\!\!>\) HiggsBounds_results.dat. A sample of this output is shown in Fig. 10. The key to the process numbering is written to the file \(<\!\mathtt{prefix}\!>\) Key.dat.

Fig. 10
figure 10

Sample output file of the type \(<\mathtt{prefix}>\) HiggsBounds_Results.dat obtained by running HiggsBounds from the command line. The results for the first four parameter points are shown

When the command-line version of HiggsBounds is used with whichinput=SLHA, the HiggsBounds results are added to the SLHA file in the form of the block HiggsBoundsResults, see Table 6. It should be noted that it is not efficient to use the command-line version of HiggsBounds with SLHA input for large parameter scans, since the experimental data tables must be read in again for each SLHA file. If this is a concern, a better option is to use the HiggsBounds subroutines to create a program which can be called from the command line. An example program, HBwithSLHA, demonstrating this is included.

5.5.3 Example input files

The HiggsBounds package includes a full set of sample input files for the case \(n_H=3\), \(n_{H^+}=1\), contained in the folder example_data. Each filename is prefixed with HB_randomtest50points_. To run the command-line version of HiggsBounds with these files as input, use, for example,

$$\begin{aligned}&\mathtt {./HiggsBounds}\,\,\mathtt {LandH}\,\,\mathtt {effC}\,\,\mathtt {3}\,\,\mathtt {1\,\, 'example\_data/}\\&\quad \mathtt {HB\_randomtest50points\_'} \end{aligned}$$

where the values of whichanalyses and whichinput can be varied as desired. The setting \(\mathtt{nHplus}=0\) can also be used if the user does not wish to test the charged Higgs sector,

$$\begin{aligned}&\mathtt {./HiggsBounds}\,\,\mathtt {LandH}\,\,\mathtt {effC}\,\,\mathtt {3}\,\,\mathtt {0\,\, 'example\_data/}\\&\quad \mathtt {HB\_randomtest50points\_'} \end{aligned}$$

5.6 Example programs

We provide a number of example programs which demonstrate the different features of the HiggsBounds subroutines. These are available in the subfolder /example_programs/ of the main installation directory. After the HiggsBounds library (libHB.a) has been compiled (using ./configure; make as described previously), each of the examples can be compiled and with the command

$$\begin{aligned} \mathtt{{make <\!\!program}\,\,\mathtt {name\!\!>}} \end{aligned}$$

More generally, any program linking with the HiggsBounds libraries can be compiled (assuming gfortran is used) as follows

$$\begin{aligned}&\mathtt{{gfortran}}\,\,\mathtt {example\_program.F} \\&\quad \mathtt { -o}\,\, \mathtt {example\_program}\,\mathtt{{-L<\!HBpath\!> -lHB}} \end{aligned}$$

The following example programs are provided with HiggsBounds-4:

  • SM_vs_4thGen This first example compares the Higgs exclusion limits in the SM to those in a model with a heavy fourth fermion generation. The program demonstrates the use of the HiggsBounds functions for the SM branching ratios and total decay widths to calculate the Higgs decay widths and the effective normalized squared couplings. This information is then used as input, and for each mass point HiggsBounds is called once with the SM input and once for the 4th generation model. The results are written to two separate output files.

  • HBwithFH This program demonstrate the use of the subroutine version of HiggsBounds to test exclusion of MSSM model points. The theory predictions are provided by linking to the publicly available MSSM Higgs spectrum calculator FeynHiggs [4145]. The model parameters should be specified in the source file (see the code for details). The results are written directly to the screen.

  • HBwithCPsuperH This example is similar to the HBwithFH example above, but uses the spectrum generator CPSuperH [6264] for the theory predictions instead of FeynHiggs. As above, the model parameters should be specified directly in the source file and the results are written directly to the screen.

  • HBwithFH_dm This is an updated version of the HBwithFH example, which demonstrates the use of several new features in HiggsBounds-4. It makes it use of FeynHiggs-calculated Higgs mass uncertainties when evaluating the exclusion limits, and the output is provided both as a combined result and in terms of exclusion information from the individual MSSM Higgs bosons. Also in this case the model parameters should be specified directly in the source file, and the output is written to the screen.

  • HBwithSLHA When using this example, the user can provide input in the SLHA format with one or more input files. The program run settings are fixed as \(\mathtt{|<\!\!whichanalyses\!\!>=}\mathtt{LandH|}\), \(\mathtt{|<\!\!nHzero\!\!>=3|}\) and \(\mathtt{|<\!\!nHplus\!\!>=1|}\), but this can be changed in the code. The set of SLHA files to be used as input should be named \(\mathtt{|<\!\!stem\!\!>.i|}\) where \(\mathtt{|i|}\) \(=1\ldots \) \(\mathtt{|npoints|}\). The program can be called from the command line as:

    $$\begin{aligned}&\mathtt{{./example\_programs/HBwithSLHA}} \\&\quad \mathtt{{<\!\!npoints\!\!> <\!\!stem\!\!>}} \end{aligned}$$

    The block \(\mathtt{|HiggsBoundsResults|}\) will be added to each SLHA file. In addition, the HiggsBounds results (\(\mathtt{|i|, |HBresult|, |chan| ,|obsratio|,} \mathtt{|ncombined|}\)) are collected in a summary output file called \(\mathtt{|<\!stem\!>-fromHB|}\). As in the command-line version, the HiggsBounds results are obtained using either the full (default) or classic method, following the setting of the corresponding flag in usefulbits.f90.

In addition to the programs listed here, there are two more example codes specifically for the LEP \(\chi ^2\) extension. These are discussed in the next section.

5.7 Installing and using the LEP \(\varvec{\chi }^2\) extension

The usage of the LEP \(\chi ^2\) information is restricted to the subroutine version of HiggsBounds. In order to enable this feature, the user first needs to download a separate package containing the binary files with the relevant experimental information from the URL

http://higgsbounds.hepforge.org.

These files are contained in the tarball \(\mathtt{csboutput\_}\mathtt{trans\_binary.tar.gz}\), which should be extracted to a user-defined directory \(\mathtt{<\!clsbtablesdir\!>}\) (not exceeding 80 characters), such that the following file structure is obtained:

$$\begin{aligned}&\mathtt{{<\!clsbtablesdir\!>/csboutput\_trans\_binary/}}\\&\mathtt{{*.binary}} \end{aligned}$$

A convenient choice for \(\mathtt{<\!\!clsbtablesdir\!\!>}\) might be the HiggsBounds main directory. In the next step, \(\mathtt{<\!\!clsbtablesdir\!\!>}\) has to be specified in the script configure-with-chisq, in addition to the usual compiler settings, cf. Sect. 5.1. Then, the HiggsBounds library can be built with:

$$\begin{aligned}&\mathtt{{./configure-with-chisq}}\\&\mathtt{{make}} \end{aligned}$$

After a successful compilation, new subroutines for the LEP \(\chi ^2\) extension are available. These are described in the following.

figure o

This subroutine initializes the new arrays and tables needed for the LEP \(\chi ^2\) extension. It reads in all the relevant experimental information from the binary files installed in \(\mathtt{<\!\!clsbtablesdir\!\!>}\).

figure p

This routine is run to calculate the LEP \(\chi ^2\) value. The user can specify a theoretical mass uncertainty (in GeV), theory_uncertainty. Note that this value is only used here, and not in “standard” mass uncertainties for the limits (which can be different). The resulting \(\chi ^2\) value is reported both including (chisq_withtheory) and without (chisq_withouttheory) this Higgs mass uncertainty. The channel code for the experimental analysis from which the \(\chi ^2\) value is derived is also given. This subroutine requires a preceeding call to the subroutine run_HiggsBounds_classic, in order to determine the most sensitive analyses for the model. Therefore, the usage of the LEP \(\chi ^2\) extension always requires a simultaneous run of the standard HiggsBounds program.

$$\begin{aligned}&\mathtt {finish\_HiggsBounds\_chisqtables}(\,) \end{aligned}$$

This deallocates the new arrays and tables, and should be called at the end of a run.

5.7.1 Usage

The typical sequence of subroutine calls when using the LEP \(\chi ^2\) extension is the following:

figure q

Note that the LEP \(\chi ^2\) functionality requires a classic HiggsBounds run to determine the most sensitive channel. For a consistent combination of LEP \(\chi ^2\) extension with Tevatron and LHC limits we recommend to initialize the LEP \(\chi ^2\) functionality with the option whichanalyses=“onlyL” and performing a separate HiggsBounds run to consider the hadronic collider limits, i.e. whichanalyses=“onlyH”. Two example programs are provided for the LEP \(\chi ^2\) extension. They are called HBchisq and HBchisqwithSLHA, respectively, and are both contained in the directory example_programs of the main HiggsBounds directory. After setting up HiggsBounds to use the LEP \(\chi ^2\) extensions, these examples can be compiled with

$$\begin{aligned} \mathtt{{make~HBchisq}} \end{aligned}$$

The first example, HBchisq, simply scans over the SM Higgs boson mass within the range \(M_H \in {[100,\ 120]}\ \text {GeV}\) and evaluates the LEP exclusion \(\chi ^2\) value. The second program, HBchisqwithSLHA, runs HiggsBounds on a set of \(n\) SLHA files, named \(\mathtt {<}SLHA-filename{>}. \mathtt i \) with \(\mathtt i =1\ldots n\). It is called as

$$\begin{aligned}&\mathtt{{/HBchisqwithSLHA}}~\mathtt{{<\!number~of~files\!>}}\\&\quad \mathtt{{<\!SLHA-filename\!>}} \end{aligned}$$

The results, including the LEP \(\chi ^2\) values, are printed for all parameter points to the new file \(\mathtt{<SLHA-filename>-} \mathtt{fromHB}\).

6 Conclusions

We have presented HiggsBounds-4, an extension of the HiggsBounds program which can be used to study exclusion bounds on arbitrary Higgs sectors using experimental results from LEP, the Tevatron and the LHC. It includes the latest LHC results presented in 2013, many of which are based on the full 8 TeV dataset.

We briefly reviewed the options for user input, including a new option (in the case of a supersymmetric Higgs sector) which allows an SLHA file to be used as input. Several improvements and updates of the code have been presented. This includes in particular an improved SM-likeness test that now takes into account the relative weight of the contributing channels, and in this way substantially enlarges the parameter space in which SM analyses can be applied. We have included the option of a theoretical Higgs mass uncertainty, which can be relevant, e.g., in the MSSM. Taking the theory uncertainties into account conservatively broadens the range of non-excluded Higgs mass values. Concerning the LEP limits, we include an option to obtain the full \(\chi ^2\) information, i.e. not “only” a hard 95 % C.L. cut. This is particularly useful for fits in the Higgs sector.

In view of the discovery of a Higgs signal at the LHC at \(\sim 125.5 \,\, \mathrm {GeV}\) we have included the option to test every Higgs boson in the model under consideration individually. In this way we slightly deviate from the pure 95 % C.L. exclusion limit, but we ensure that models do not falsely pass the HiggsBounds test because the spectrum contains one (SM-like) Higgs boson at a mass of \(\sim 125.5 \,\, \mathrm {GeV}\).

HiggsBounds can now readily be used together with its new sister code, HiggsSignals [15]. HiggsSignals performs a \(\chi ^2\) evaluation of the compatibility between the predictions of arbitrary Higgs sectors to measured signal rates. This includes in particular the possibility to test the model predictions against the observed signal at \(\sim 125.5 \,\, \mathrm {GeV}\), but also future, hypothetical, signals of extended Higgs sectors. A combined analysis using both codes exploits all the public information on the Higgs signal and the Higgs exclusion bounds obtained at LEP, the Tevatron and the LHC.