Weak Radiative Decays of the B Meson and Bounds on $M_{H^\pm}$ in the Two-Higgs-Doublet Model

TL;DR

This paper updates the bounds on the charged Higgs boson mass in the Two-Higgs-Doublet Model based on recent Belle collaboration data on B-meson decays, highlighting the impact of experimental improvements and methodological choices.

Contribution

It provides revised bounds on $M_{H^}$ in the Two-Higgs-Doublet Model using the latest experimental data and clarifies ambiguities related to photon energy cutoff choices.

Findings

Lower bound on $M_{H^}$ is now 570-800 GeV at 95% C.L.

Updated bounds are sensitive to the method of analysis and photon energy cutoff.

The experimental branching ratio is below the Standard Model prediction but remains consistent.

Abstract

In a recent publication [6], the Belle collaboration updated their analysis of the inclusive weak radiative B-meson decay, including the full dataset of (772 +_ 11) x 10^6 BBbar pairs. Their result for the branching ratio is now below the Standard Model prediction [7,8], though it remains consistent with it. However, bounds on the charged Higgs boson mass in the Two-Higgs-Doublet Model get affected in a significant manner. In the so-called Model-II, the 95% C.L. lower bound on $M_{H^\pm}$ is now in the 570-800 GeV range, depending quite sensitively on the method applied for its determination. Our present note is devoted to presenting and discussing the updated bounds, as well as to clarifying several ambiguities that one might encounter in evaluating them. One of such ambiguities stems from the photon energy cutoff choice, which deserves re-consideration in view of the improved experimental accuracy.

Tables

Table 1. Table 1: Quantities Δ q subscript Δ 𝑞 \Delta_{q} from Eq. ( 7 ) evaluated using various approaches.

$E_0$[GeV]	${\Delta }_s^{\text{BF}}$	${\Delta }_s^{\text{Belle}}$	${\Delta }_s^{\text{fix}}$	${\Delta }_{s+d}^{\text{fix}}$	${\Delta }_d^{\text{fix}}$
$1.7$	$\left(1.5\pm 0.4\right)\%$	?	$1.3\%$	$1.5\%$	$5.3\%$
$1.8$	$\left(3.4\pm 0.6\right)\%$	$\left(3.69\pm 1.39\right)$%	$3.0\%$	$3.4\%$	$10.5\%$
$1.9$	$\left(6.8\pm 1.1\right)\%$	?	$5.5\%$	$6.0\%$	$15.7\%$
$2.0$	$\left(11.9\pm 2.0\right)\%$	?	$10.0\%$	$10.5\%$	$22.5\%$

Table 2. Table 2: Experimental results and their averages for ℬ s ​ γ × 10 6 subscript ℬ 𝑠 𝛾 superscript 10 6 {\mathcal{B}}_{s\gamma}\times 10^{6} (upper rows) and ℬ ( s + d ) ​ γ × 10 6 subscript ℬ 𝑠 𝑑 𝛾 superscript 10 6 {\mathcal{B}}_{(s+d)\gamma}\times 10^{6} (lower rows) at each value of E 0 subscript 𝐸 0 E_{0} . Each world average (w.a.) is first calculated at E 0 subscript 𝐸 0 E_{0} (10th column), and then extrapolated to 1.6 1.6 1.6\, GeV (11th column) using Δ s BF superscript subscript Δ 𝑠 BF \Delta_{s}^{\!\mbox{\tiny BF}} or Δ s + d BF superscript subscript Δ 𝑠 𝑑 BF \Delta_{s+d}^{\!\mbox{\tiny BF}} from Section  2 . In the last two columns, the ratios R γ ( × 10 5 ) R_{\gamma}(\times 10^{5}) are calculated from the corresponding averages for ℬ ( s + d ) ​ γ subscript ℬ 𝑠 𝑑 𝛾 {\mathcal{B}}_{(s+d)\gamma} using ℬ c ​ ℓ ​ ν = 0.1067 ​ ( 16 ) subscript ℬ 𝑐 ℓ 𝜈 0.1067 16 {\mathcal{B}}_{c\ell\nu}=0.1067(16) .

	Babar				Belle			CLEO	w.a.	w.a.	$R_{\gamma }$	$R_{\gamma }$
$E_0$	[9]	[44]	[45]	aver	[6]	[46]	aver	[47]	($E_0$)	(1.6)	($E_0$)	(1.6)
1.7					$306\left(28\right)$		$306\left(28\right)$		$306\left(28\right)$	$311\left(28\right)$
					$\mathrm{𝟑𝟐𝟎}\left(\mathrm{𝟐𝟗}\right)$		$320\left(29\right)$		$320\left(29\right)$	$326\left(30\right)$	$300\left(28\right)$	$305\left(28\right)$
1.8	$321\left(34\right)$			$321\left(34\right)$	$301\left(22\right)$		$301\left(22\right)$		$307\left(19\right)$	$318\left(19\right)$
	$\mathrm{𝟑𝟑𝟓}\left(\mathrm{𝟑𝟓}\right)$			$335\left(35\right)$	$\mathrm{𝟑𝟏𝟓}\left(\mathrm{𝟐𝟑}\right)$		$315\left(23\right)$		$321\left(19\right)$	$333\left(20\right)$	$301\left(19\right)$	$312\left(19\right)$
1.9	$300\left(24\right)$	$\mathrm{𝟑𝟐𝟗}\left(\mathrm{𝟓𝟐}\right)$	$366\left(104\right)$	$308\left(22\right)$	$294\left(18\right)$	$\mathrm{𝟑𝟓𝟏}\left(\mathrm{𝟑𝟕}\right)$	$305\left(16\right)$		$306\left(13\right)$	$327\left(14\right)$
	$\mathrm{𝟑𝟏𝟑}\left(\mathrm{𝟐𝟓}\right)$	$344\left(54\right)$	$\mathrm{𝟑𝟖𝟏}\left(\mathrm{𝟏𝟎𝟖}\right)$	$321\left(23\right)$	$\mathrm{𝟑𝟎𝟕}\left(\mathrm{𝟏𝟗}\right)$	$367\left(39\right)$	$319\left(17\right)$		$320\left(14\right)$	$343\left(15\right)$	$300\left(14\right)$	$322\left(15\right)$
2.0	$280\left(19\right)$		$339\left(79\right)$	$283\left(18\right)$	$279\left(15\right)$		$279\left(15\right)$	$293\left(46\right)$	$281\left(11\right)$	$315\left(14\right)$
	$\mathrm{𝟐𝟗𝟐}\left(\mathrm{𝟐𝟎}\right)$		$\mathrm{𝟑𝟓𝟑}\left(\mathrm{𝟖𝟑}\right)$	$296\left(19\right)$	$\mathrm{𝟐𝟗𝟐}\left(\mathrm{𝟏𝟓}\right)$		$292\left(15\right)$	$\mathrm{𝟑𝟎𝟔}\left(\mathrm{𝟒𝟗}\right)$	$294\left(11\right)$	$331\left(14\right)$	$276\left(11\right)$	$310\left(14\right)$

Table 3. Table 3: Bounds on M H ± subscript 𝑀 superscript 𝐻 plus-or-minus M_{H^{\pm}} obtained using different methods.

Model	$R_{\gamma }^{\exp }\times {10}^3$	$95\%$C.L. bounds			$99\%$C.L. bounds
		1-sided	2-sided	FC	1-sided	2-sided	FC
	$3.05\pm 0.28$	307	268	268	230	208	208
I	$3.12\pm 0.19$	401	356	356	313	288	288
($\tan \beta =1$)	$3.22\pm 0.15$	504	445	445	391	361	361
	$3.05\pm 0.28$	740	591	569	477	420	411
II	$3.12\pm 0.19$	795	645	628	528	468	461
(absolute)	$3.22\pm 0.15$	692	583	580	490	440	439

Equations

$$R_{\gamma }=\frac{{\mathcal{B}}_{s\gamma }+{\mathcal{B}}_{d\gamma }}{{\mathcal{B}}_{c\ell \nu }}\equiv \frac{{\mathcal{B}}_{(s+d)\gamma }}{{\mathcal{B}}_{c\ell \nu }}.$$

$${\mathcal{L}}_{\mathrm{weak}}\sim \sum _i{C}_i{Q}_i,$$

$$Q_1=({\bar{s}}_L{\gamma }_{\mu }{T}^a{c}_L)({\bar{c}}_L{\gamma }^{\mu }{T}^a{b}_L),Q_2=({\bar{s}}_L{\gamma }_{\mu }{c}_L)({\bar{c}}_L{\gamma }^{\mu }{b}_L),$$

$$Q_7=\frac{e}{16{\pi }^2}{m}_b({\bar{s}}_L{\sigma }^{\mu \nu }{b}_R)F_{\mu \nu },Q_8=\frac{g}{16{\pi }^2}{m}_b({\bar{s}}_L{\sigma }^{\mu \nu }{T}^a{b}_R)G_{\mu \nu }^a.$$

$$\Gamma (\bar{B}\to X_s\gamma )=\Gamma (b\to X_s^p\gamma )+\delta {\Gamma }_{\mathrm{nonp}},$$

$${\left[\frac{\delta {\Gamma }_{\mathrm{nonp}}}{\Gamma (b\to X_s^p\gamma )}\right]}_{\mathrm{only}{C}_7}=-\frac{{\mu }_{\pi }^2+3{\mu }_G^2}{2m_b^2}+\mathcal{O}\left(\frac{{\alpha }_{\mathrm{s}}{\Lambda }^2}{(m_b-2E_0{)}^2},\frac{{\Lambda }^3}{m_b^3}\right),$$

$${\Delta }_q\equiv \frac{{\mathcal{B}}_{q\gamma }(1.6)}{{\mathcal{B}}_{q\gamma }(E_0)}-1,$$

$$R_{\gamma }^{\mathrm{SM}}=(3.31\pm 0.22)\times 10^{-3},$$

$$R_{\gamma }^{\exp }=(3.22\pm 0.15)\times 10^{-3}$$

Full text

TTP17-004

IFT-1/2017

**Weak Radiative Decays of the $B$ Meson and

Bounds on $M_{H^{\pm}}$ in the Two-Higgs-Doublet Model

**

Mikołaj Misiak1,2 and  Matthias Steinhauser3

1 *Institute of Theoretical Physics, Faculty of Physics, University of Warsaw,

02-093 Warsaw, Poland.

2 *Theoretical Physics Department, CERN, CH-1211 Geneva 23, Switzerland.

3 *Institut für Theoretische Teilchenphysik, Karlsruhe Institute of Technology (KIT),

76128 Karlsruhe, Germany.

**Abstract

**

In a recent publication [6], the Belle collaboration updated their analysis of the inclusive weak radiative $B$-meson decay, including the full dataset of $(772\pm 11)\times 10^{6}~{}B\bar{B}$ pairs. Their result for the branching ratio is now below the Standard Model prediction [7, 8], though it remains consistent with it. However, bounds on the charged Higgs boson mass in the Two-Higgs-Doublet Model get affected in a significant manner. In the so-called Model-II, the $95\%\,$C.L. lower bound on $M_{H^{\pm}}$ is now in the $570$–$800\,$GeV range, depending quite sensitively on the method applied for its determination. Our present note is devoted to presenting and discussing the updated bounds, as well as to clarifying several ambiguities that one might encounter in evaluating them. One of such ambiguities stems from the photon energy cutoff choice, which deserves re-consideration in view of the improved experimental accuracy.

1 Introduction

In the absence of any new strongly-interacting particles discovered at the LHC, one observes growing interest in models where kinematically accessible exotic particles take part in the electroweak (EW) interactions only. The simplest of such models are constructed by extending the Higgs sector of the Standard Model (SM) via introduction of another $SU(2)_{\rm weak}$ doublet. There are several versions of the Two-Higgs-Doublet Model (2HDM) that differ in the Higgs boson couplings to fermions. They are usually arranged in such a way that no tree-level Flavour-Changing Neutral Current interactions arise [1, 2]. In the so-called Model-I, fermions receive their masses in the SM-like manner, from Yukawa couplings to only one of the two Higgs doublets. In Model-II, the Yukawa couplings are as in the Minimal Supersymmetric Standard Model, i.e. one of the doublets called $H_{u}$ gives masses to the up-type quarks, while the other doublet, $H_{d}$, gives masses to both the down-type quarks and the leptons.111 Couplings to leptons are irrelevant for our considerations throughout the paper. Thus, whatever we write about Models I and II, is also true for the models called “X” and “Y”, respectively – see Tabs. 1 and 6 of Ref. [3].

Within the 2HDM, the physical spin-zero boson spectrum consists of one charged scalar $H^{\pm}$, one neutral pseudoscalar $A^{0}$, and a pair of scalars, $H^{0}$ and $h^{0}$, the latter of which is identified with the recently discovered SM-like Higgs boson. If the beyond-SM (BSM) scalars become very heavy ($M_{H^{\pm}},M_{A^{0}},M_{H^{0}}\sim M\gg m_{h^{0}}$), they undergo decoupling, and the model reduces to the SM at scales much smaller than $M$. Thus, any claims [4] concerning exclusion of the 2HDM in its full parameter space imply claiming exclusion of the SM, too.

As it is well known (see, e.g., Ref. [5]), strong constraints on $M_{H^{\pm}}$ follow from measurements of the inclusive weak radiative $B$-meson decay branching ratio. The most precise results come from the Belle collaboration, especially from their recent analysis based on the full $(772\pm 11)\times 10^{6}~{}B\bar{B}$ pair dataset [6]. Their updated result is now below the Standard Model prediction [7, 8], though it remains consistent with it. On the other hand, the 2HDM effects in Model-II can only enhance the decay rate. In consequence, the lower bound on $M_{H^{\pm}}$ in this model becomes very strong, reaching the range of $570$–$800\,$GeV. At the same time, the bound becomes very sensitive to the method applied for its determination. Given the relevance of the considered bound for many popular BSM models with extended Higgs sectors, a detailed discussion of this issue is necessary. This is the main purpose of our present paper.

Following Eqs. (1.1) and (1.2) of Ref. [8], as well as Eq. (9) of Ref. [7], we shall use the CP- and isospin-averaged branching ratios ${\mathcal{B}}_{s\gamma}$ and ${\mathcal{B}}_{d\gamma}$ of the weak radiative decays, normalizing them to the analogously averaged branching ratio ${\mathcal{B}}_{c\ell\nu}$ of the semileptonic decay. The main observable for our considerations will be the ratio

[TABLE]

We prefer to use $R_{\gamma}$ rather than ${\mathcal{B}}_{s\gamma}$ for two reasons. First, the currently most precise experimental results come from the fully inclusive analyses of Belle [6] and Babar [9] where the actually measured quantity is ${\mathcal{B}}_{(s+d)\gamma}$. Second, a normalization to the semileptonic rate removes the main contribution to the parametric uncertainty of around 1.5% on the theory side, while it cannot introduce a larger uncertainty on the experimental side. Our treatment of the experimental results will be described in Section 4.

Constraints on $M_{H^{\pm}}$ from observables other than $R_{\gamma}$ (or ${\mathcal{B}}_{s\gamma}$) have been reviewed in several recent articles – see, e.g., Refs. [3, 10, 11]. Their common property in Model-II is that they become relevant either for small or for very large ratio of the vacuum expectation values $v_{u}/v_{d}\equiv\tan\beta$. On the other hand, $R_{\gamma}$ provides a bound that cannot be avoided for any $\tan\beta$, and turns out to be the strongest one in quite a wide range of $\tan\beta$.

Our paper is organized as follows. In the next section, the basic framework for analyzing the considered decays is outlined, with extended explanations concerning the photon energy cutoff issue. Section 3 is devoted to recalling the SM prediction for $R_{\gamma}$, and discussing the size of possible extra contributions in the 2HDM. In Section 4, we collect all the available experimental results for ${\mathcal{B}}_{s\gamma}$ and/or ${\mathcal{B}}_{(s+d)\gamma}$, calculate their weighted averages in several ways, and convert them to $R_{\gamma}$. The resulting bounds on $M_{H^{\pm}}$ are derived in Section 5. We conclude in Section 6.

2 Photon energy cutoff in ${\mathcal{B}}_{s\gamma}$ and ${\mathcal{B}}_{d\gamma}$

The radiative decays we are interested in proceed dominantly via quark-level transitions $b\to s\gamma$,  $b\to d\gamma$, and their CP-conjugates. A suppression by small Cabibbo-Kobayashi-Maskawa (CKM) angles makes ${\mathcal{B}}_{d\gamma}$ about 20 times smaller than ${\mathcal{B}}_{s\gamma}$. For definiteness, we shall discuss ${\mathcal{B}}_{s\gamma}$ in what follows, making separate comments on ${\mathcal{B}}_{d\gamma}$ only at points where anything beyond a trivial replacement of the quark flavours matters.

Theoretical analyses of rare $B$-meson decays are most conveniently performed in the framework of an effective theory that arises after decoupling the $W$-boson, the heavier SM particles, and all the (relevant) BSM particles. We assume here that the BSM particles being decoupled are much heavier than the $b$-quark, but their masses are not much above a TeV. In such a case, the decoupling can be performed in a single step, at a common renormalization scale $\mu_{0}\sim m_{t}$. In the effective theory below $\mu_{0}$, the weak interaction Lagrangian that matters for $b\to s\gamma$ takes the form

[TABLE]

where $Q_{i}$ are dimension-five and -six operators of either four-quark type, e.g.,

[TABLE]

or dipole type

[TABLE]

A complete list of $Q_{i}$ that matter in the SM or 2HDM at the Leading Order222 For the Next-to-Leading (NLO) EW corrections, several extra four-quark operators need to be included – see Eq. (2) of Ref. [12]. In that paper, such corrections were calculated within the SM. The 2HDM case is still pending. (LO) in $\alpha_{\rm em}$ can be found in Eq. (1.6) of Ref. [8]. Their Wilson coefficients $C_{i}(\mu_{0})$ are evaluated perturbatively in $\alpha_{\mathrm{s}}$ by matching several effective-theory Green’s functions with those of the SM or 2HDM. Such calculations have now reached the Next-to-Next-to-Leading Order (NNLO) accuracy in QCD, i.e. $C_{i}(\mu_{0})$ are known up to ${\mathcal{O}}(\alpha_{\mathrm{s}}^{2})$. In the dipole operator case, performing a three-loop matching [13, 14] was necessary at the NNLO.

In the next step, the Wilson coefficients are evolved according to their renormalization group equations down to the scale $\mu_{b}\sim m_{b}$, in order to resum large logarithms of the form $\left(\alpha_{\mathrm{s}}\ln(\mu_{0}^{2}/\mu_{b}^{2})\right)^{n}\sim\left(\alpha_{\mathrm{s}}\ln(m_{t}^{2}/m_{b}^{2})\right)^{n}$. To achieve this at the NNLO level, anomalous dimension matrices up to four loops [15] had to be determined. At present, all the Wilson coefficients $C_{i}(\mu_{b})$ are known with a precision that is sufficient for evaluating $R_{\gamma}$ at the NNLO in QCD.

While the calculations of $C_{i}(\mu_{b})$ are purely perturbative, one needs to take nonperturbative effects into account when determining the physical decay rates. For $\bar{B}\to X_{s}\gamma$ (with $\bar{B}$ denoting either $\bar{B}^{0}$ or $B^{-}$), the decay rate is a sum of the dominant perturbative contribution and a subdominant nonperturbative one $\delta\Gamma_{\rm nonp}$, i.e.

[TABLE]

where a photon energy cutoff $E_{\gamma}>E_{0}$ in the decaying particle rest frame is imposed on both sides.333 The rates would be ill-defined without such a cutoff. The partonic final state $X_{s}^{p}$ is assumed to consist of charmless quarks and gluons, while the hadronic state $X_{s}$ is assumed to contain no charmed or $c\bar{c}$ hadrons. The latter requirement is in principle stronger than demanding that $X_{s}$ as a whole is charmless. However, all the measurements to date have been performed with $E_{0}\geq 1.7\,$GeV, in which case the $c\bar{c}$ hadrons and/or pairs of charmed hadrons are kinematically forbidden in $X_{s}$ anyway. There is no experimental restriction on extra photons or lepton pairs in $X_{s}$, but their contribution corresponds to very small NLO QED corrections that are only partly included on the theory side.

The nonperturbative contribution $\delta\Gamma_{\rm nonp}$ in Eq. (5) is strongly dependent on $E_{0}$. For $E_{0}=1.6\,$GeV, it shifts the SM prediction for ${\mathcal{B}}_{s\gamma}$ by almost $+3\%$ [16],444 Such a central value of the shift corresponds actually to the effect of $N(E_{0})$ in Eq. (D.4) of Ref. [8] where a normalization to the semileptonic rate was used, and some of the nonperturbative effects were relegated to the semileptonic phase-space factor $C$ in Eqs. (D.2)–(D.3) there. while the corresponding uncertainty is estimated at the $\pm 5\%$ level [17]. For higher values of $E_{0}$, theoretical uncertainties grow (see below), while the experimental ones decrease thanks to lower background subtraction errors. To resolve this issue, it has become standard to perform a data-driven extrapolation of the experimental results down to $E_{0}=1.6\,$GeV, and compare with theory at that point.

A few comments about such an extrapolation need to be made. First, it is instructive to have a look at Fig. 1 which presents the background-subtracted photon energy ($E^{*}_{\gamma}$) spectrum in the $\Upsilon(4S)$ frame, as determined by Belle in their full-dataset measurement [6]. Photon energies $E_{\gamma}$ in the $B$-meson rest frame differ from $E^{*}_{\gamma}$ by boost factors that do not exceed 1.07. One can see that energies below $2\,$GeV are well in the tail of the spectrum. On the other hand, a large set of measurements that gives quite a precise weighted average for ${\mathcal{B}}_{(s+d)\gamma}$ is available already at $E_{0}=1.9\,$GeV (see Section 4). Thus, the extrapolation we need is really a short one, and only in the tail of the spectrum.

To understand the growth of theoretical uncertainties with $E_{0}$, one begins with considering the case when $C_{7}$ is assumed to be the only nonvanishing Wilson coefficient at the scale $\mu_{b}$. In such a case, the fixed-order Heavy Quark Effective Theory (HQET) formalism can be used to show that [18, 19, 20]

[TABLE]

provided  $m_{b}-2E_{0}\gg\Lambda$  with  $\Lambda\sim\Lambda_{\rm QCD}$.  The quantities $\mu_{\pi}^{2}$ and $\mu_{G}^{2}$ are of order $\Lambda^{2}$, and are currently quite well known from fits to the measured semileptonic decay spectra [21]. With growing $E_{0}$, at some point one enters into the region where  $m_{b}-2E_{0}\sim\Lambda$,  and the fixed-order HQET calculation is no longer applicable. Instead, the leading nonperturbative effect is parameterized in terms of a universal shape function [22, 23]. We need to rely on models for this function, which is the main reason why the theory uncertainties grow with $E_{0}$.

A number of shape-function models have been invented in the past, with their parameters constrained by measurements of the semileptonic and radiative $B$-meson decay spectra – see, e.g., Refs. [24, 25]. In Fig. 13 of Ref. [25], one can see that the $\bar{B}\to X_{s}\gamma$ photon energy spectrum becomes quite unique already at $E_{\gamma}=1.9\,$GeV, at least for the considered class of models. Such a uniqueness is indeed expected below the point where the shape-function description starts to overlap with the fixed-order HQET description. Future studies with precise Belle-II data should shed more light on the actual location of this point. Our present approach relies on the assumption that $E_{0}=1.6\,$GeV is definitely below this point.

Choosing $1.6\,$GeV as the default $E_{0}$ to compare the fixed-order HQET predictions with the (extrapolated) experimental results for ${\mathcal{B}}_{s\gamma}$ was first suggested in Ref. [26], at the time when no precise data on the spectrum were available, and one had to rely on a limited class of shape-function models. At present, one might wonder whether this default $E_{0}$ might be shifted upwards. However, since such a decision would need to be made on the basis of the experimental data, shifting the default $E_{0}$ could hardly improve anything with respect to the current extrapolation approach. Another question that one might ask is whether the extrapolation method (say, from $1.9$ to $1.6$) is indeed superior with respect to direct measurements at lower values of $E_{0}$ (but in the range $[1,6,1.9]$). The answer depends on the balance of uncertainties: the extrapolation ones and the background subtraction ones. We shall come back to this issue in Section 4.

Effects of extrapolations from $E_{0}$ to $1.6\,$GeV can be parameterized by

[TABLE]

with $q=s,d$ or $s+d$. Numerical values of this quantity obtained with the help of various methods are presented in Tab. 1. Those denoted by $\Delta_{s}^{\!\mbox{\tiny BF}}$ were evaluated in Ref. [27] where the measured semileptonic and radiative $B$-meson decay spectra (as available in 2005) were used to determine the $b$-quark mass $m_{b}$ and the parameter $\mu_{\pi}^{2}$ in three different renormalization schemes. Next, these parameters were inserted into the Kagan-Neubert shape function model (Eq. (24) of Ref. [24]). The shape function was then convoluted with the perturbatively calculated photon energy spectrum in the $b$-quark decay, which led to a prediction for the physical photon energy spectrum in the $B$-meson decay.

In the next column of Tab. 1, the quantities $\Delta_{s}^{\!\mbox{\tiny Belle}}$ were obtained in Ref. [6] using essentially the same method but with the radiative spectrum only, as measured in the very analysis of Ref. [6]. In that case, only the result for $E_{0}=1.8\,$GeV is publicly available at present. The shape function model was used in the experimental analysis not only for the extrapolation in $E_{0}$, but also for efficiency estimates and boosting between $E_{\gamma}^{*}$ and $E_{\gamma}$. The best fit for $m_{b}$ and $\mu_{\pi}^{2}$ in Ref. [6] leads to a good description of the measured spectrum (solid histogram in Fig. 1), and at the same time is consistent with the semileptonic fits [21].

The last three columns of Tab. 1 have been obtained using the approach of Refs. [7, 8] (perturbative & fixed-order HQET), in which case the photon energy spectrum is determined mainly by the perturbative gluon bremsstrahlung. In these cases, no uncertainties are quoted, as we do not know at which $E_{0}$ the fixed-order HQET description breaks down. The subleading ${\mathcal{O}}(\alpha_{\mathrm{s}}\Lambda^{2})$ nonperturbative corrections [28] begin to rapidly increase at $E_{0}$ around $1.8\,$GeV due to $(m_{b}-2E_{0})^{2}$ in their denominators, but their overall suppression factor is small, and they remain under control in the whole region of interest (up to $2\,$GeV).

The quantities $\Delta_{q}^{\!\mbox{\tiny fix}}$ involve effects of the photon bremsstrahlung in decays of the $b$ quark to three light (anti)quarks, as calculated in Refs. [29, 30]. Such effects are small in ${\mathcal{B}}_{s\gamma}$ (unless one goes well below $E_{0}=1.6\,$GeV) but become much more relevant in ${\mathcal{B}}_{d\gamma}$ where the tree-level $b\to du\bar{u}\gamma$ transitions are not CKM-suppressed with respect to the leading $b\to d\gamma$ one. In effect, $\Delta_{d}^{\!\mbox{\tiny fix}}$ are visibly different from $\Delta_{s}^{\!\mbox{\tiny fix}}$. However, $\Delta_{s+d}^{\!\mbox{\tiny fix}}$ is not much different from $\Delta_{s}^{\!\mbox{\tiny fix}}$ due to the dominance of ${\mathcal{B}}_{s\gamma}$ over ${\mathcal{B}}_{d\gamma}$. Such photon bremsstrahlung effects involve collinear singularities in the limit of vanishing quark masses, which signals the presence of important nonperturbative effects that need to be described in terms of fragmentation functions [31], and are poorly known. Fortunately, their overall suppression factors in ${\mathcal{B}}_{s\gamma}$ and ${\mathcal{B}}_{(s+d)\gamma}$ are strong enough, and the corresponding uncertainties are far below the dominant nonperturbative ones.

It is interesting to observe in Tab. 1 that $\Delta_{s}^{\!\mbox{\tiny BF}}$ and $\Delta_{s}^{\!\mbox{\tiny Belle}}$ are quite close to $\Delta_{s}^{\!\mbox{\tiny fix}}$ and $\Delta_{s+d}^{\!\mbox{\tiny fix}}$. It gives us a hope that the breakdown of the fixed-order HQET description, even if present, is not dramatic in the considered region of $E_{0}$. In effect, our sensitivity to ambiguities in modeling the shape functions is likely to be quite limited, at least for the purpose of the $1.9\to 1.6$ extrapolations. However, a devoted analysis with the most recent data and a wide class of shape function models is necessary to estimate the corresponding uncertainty in a reliable manner.555 As follows from Ref. [17], operators other than $Q_{7}$ give rise to relevant nonperturbative effects, which may increase the extrapolation uncertainties. Since such an analysis is still awaited, we shall proceed with using $\Delta_{s}^{\!\mbox{\tiny BF}}$ in what follows for the extrapolation of ${\mathcal{B}}_{s\gamma}$. As far as the extrapolation of ${\mathcal{B}}_{(s+d)\gamma}$ is concerned, we are going to rescale $\Delta_{s}^{\!\mbox{\tiny BF}}$ according to the fixed-order results, namely use $\Delta_{s+d}^{\!\mbox{\tiny BF}}\equiv\Delta_{s}^{\!\mbox{\tiny BF}}\times\Delta_{s+d}^{\!\mbox{\tiny fix}}/\Delta_{s}^{\!\mbox{\tiny fix}}$.

3 The ratio $R_{\gamma}$ in the SM and 2HDM

Although the perturbative decay rate $\Gamma(b\to X_{s}^{p}\gamma)$ in Eq. (5) may seem straightforward to evaluate, its determination to better than $\pm 5\%$ accuracy requires including the NNLO QCD corrections, which is a highly nontrivial task. While the Wilson coefficients are already known to sufficient accuracy both in the SM and 2HDM (as already mentioned in the previous section), our knowledge of the NNLO corrections is yet incomplete in the case of matrix elements, namely interferences among on-shell decay amplitudes generated at the scale $\mu_{b}$ by the operators $Q_{i}$. The matrix elements are the same in the SM and in the 2HDM.

At the NNLO level, we can restrict our attention to the operators listed in Eqs. (3)–(4), as the remaining ones can be neglected due to their small Wilson coefficients. The $Q_{7}$-$Q_{7}$ and $Q_{7}$-$Q_{8}$ interference terms are already known at ${\mathcal{O}}(\alpha_{\mathrm{s}}^{2})$ in a complete manner [32, 33, 34, 35, 36]. The NNLO interference terms not involving $Q_{7}$ can be separated into two-body final state contributions (trivially derived from the NLO results) or relatively small $(n\geq 3)$-body final state contributions that have been calculated so far [37, 38, 39] only in the Brodsky-Lepage-Mackenzie (BLM) [40] approximation. The main perturbative uncertainty comes from the $Q_{1,2}$-$Q_{7}$ interferences at ${\mathcal{O}}(\alpha_{\mathrm{s}}^{2})$. Their BLM parts, as well as effects of nonvanishing quark masses on the gluon lines were evaluated in Refs. [37, 41, 42] for arbitrary values of the charm quark mass $m_{c}$. The remaining parts were found only in the limits $m_{c}\gg m_{b}/2$ [43] or $m_{c}=0$ [8], and then an interpolation between these two limits was performed [8].

With all the NNLO QCD, NLO EW and nonperturbative corrections evaluated to date, the SM prediction for $R_{\gamma}$ at $E_{0}=1.6\,$GeV reads [7]

[TABLE]

where the overall uncertainty has been obtained by combining in quadrature the nonperturbative one $(\pm 5\%)$, the parametric one $(\pm 1.5\%)$, the one stemming from neglected higher-order effects $(\pm 3\%)$, and the one due to the above-mentioned interpolation in $m_{c}$ $(\pm 3\%)$.

In the 2HDM, additional contributions to the Wilson coefficient matching arise from diagrams with the physical charged scalar exchanges. The relevant couplings and sample diagrams can be found, e.g., in Sec. 2.3 of Ref. [14]. We evaluate $R_{\gamma}$ in Model-I and Model-II with the same accuracy as in the SM, up to the missing NLO EW corrections to the charged Higgs contributions. Apart from the SM parameters, the results depend only on $M_{H^{\pm}}$ and $\tan\beta$. They are plotted in Fig. 2 as functions of $M_{H^{\pm}}$ in two cases of particular interest: Model-I with $\tan\beta=1$ and Model-II with $\tan\beta=50$. The solid and dashed curves in these plots correspond to the 2HDM and SM cases, respectively. Dotted lines indicate the experimental average to be discussed in the next section.

In Model-I, the charged Higgs contribution to the decay amplitude is proportional to $\cot^{2}\beta$, and it interferes with the SM one in a destructive manner. In Model-II, the interference is always constructive, and the charged Higgs amplitude has the form666 The term proportional to $\tan^{2}\beta$ is suppressed by the strange quark mass, and we neglect it here. $A+B\cot^{2}\beta$. The quantities $A$ and $B$ depend on $M_{H^{\pm}}$ only, and they have the same sign. In consequence, an absolute bound on $M_{H^{\pm}}$ can be derived from $R_{\gamma}$ in Model-II by setting the $\cot^{2}\beta$ term to zero. In practice, $R_{\gamma}$ begins to be practically independent of $\tan\beta$ already around $\tan\beta\simeq 2$. The $\tan\beta=50$ case in Fig. 2 indicates that the absolute bound on $M_{H^{\pm}}$ is going to be in the few-hundred GeV region. In Model-I, sizeable deviations of $R_{\gamma}$ from its SM value occur only for moderate or small values of $\tan\beta$. The $\tan\beta=1$ case displayed in Fig. 2 shows that our sensitivity to $M_{H^{\pm}}$ in this case is almost as strong as in Model-II.

4 Determining the current experimental average for $R_{\gamma}$

All the available measurements of ${\mathcal{B}}_{(s+d)\gamma}$ and ${\mathcal{B}}_{s\gamma}$, as well as our averages of them are collected in Tab. 2. The results of Babar have been obtained using three methods: fully inclusive [9], semi-inclusive [44], and the hadronic-tag one [45]. Belle has used the fully inclusive [6] and semi-inclusive [46] approaches, while their hadronic-tag analysis is still awaited. In the measurement of CLEO [47], the fully inclusive method was used.

The most precise results come from the fully inclusive analyses where the actually measured quantity is ${\mathcal{B}}_{(s+d)\gamma}$. The same refers to the hadronic-tag result of Babar, which is actually also fully inclusive. In the semi-inclusive cases, a single kaon in the final state was required, so the measurements accounted directly for ${\mathcal{B}}_{s\gamma}$. We indicate this in Tab. 2 by typesetting the corresponding numbers in bold.

Belle and CLEO provided their ${\mathcal{B}}_{(s+d)\gamma}$ results explicitly, while Babar rescaled them to ${\mathcal{B}}_{s\gamma}$, quoting in each case the necessary CKM factor together with its uncertainty. In Tab. 2, we “undo” the rescaling using precisely the same factors. On the other hand, in the two semi-inclusive cases, we derive ${\mathcal{B}}_{(s+d)\gamma}$ from ${\mathcal{B}}_{s\gamma}$ using a rescaling factor $(1.047\pm 0.003)$ that we calculate at $E_{0}=1.9\,$GeV as in Refs [7, 8]. Our factor differs only slightly (by $0.2\%$) from $1+|V_{td}/V_{ts}|^{2}$, due to the $b\to du\bar{u}\gamma$ effects. Rescaling the semi-inclusive results is a minor issue anyway, as they come with considerably larger experimental errors.

The reader is referred to the original experimental papers [6, 9, 44, 45, 46, 47] for the decomposition of errors into the statistical, systematic and occasionally the spectrum-modeling ones. Here we have added them in quadrature for the purpose of determining our naive averages, in which no correlations have been taken into account. In several cases, we can compare our averages with the very recent ones of HFAG [48] where, we believe, the necessary correlations have been included. For instance, the two Belle results for ${\mathcal{B}}_{s\gamma}$ at $E_{0}=1.9\,$GeV lead to the naive average of  ${\tt av}[294(18),351(37)]=305(16)$, which perfectly agrees with Ref. [48]. In the same row of Tab. 2, the two less precise results of Babar give  ${\tt av}[329(52),366(104)]=336(46)$, which again overlaps with Ref. [48]. In this case, the most precise result of Babar has not been included in the HFAG average. We have been informed that this point is going to be corrected soon [49].

As far as the world average for ${\mathcal{B}}_{s\gamma}$ extrapolated to $E_{0}=1.6\,$GeV is concerned, Ref. [48] gives $(3.32\pm 0.15)\times 10^{-4}$, which is quite close to our $(3.27\pm 0.14)\times 10^{-4}$ in the row containing the semi-inclusive measurements. We do not know which inputs have been used in this average of HFAG. Concerning the extrapolation, they have indicated using the method of Ref. [27].

Comparing the uncertainties in the four alternative averages for $R_{\gamma}$ at $E_{0}=1.6\,$GeV in the last column of Tab. 2, one can see that the first two of them are less accurate. Thus, at the moment, the balance of the background subtraction and extrapolation uncertainties points towards using the results extrapolated from 1.9 or 2.0, at least when one takes the errors from Ref. [27] for granted. Since there is not much difference in the uncertainties of these two averages, we suggest discarding the 2.0 one, as it requires a longer extrapolation. Thus, we recommend adopting

[TABLE]

as the current experimental average for $R_{\gamma}$ at $E_{0}=1.6\,$GeV.

5 Bounds on $M_{H^{\pm}}$

In this section, we shall use $R_{\gamma}$ to derive bounds on $M_{H^{\pm}}$ in the 2HDM. We are going to treat all the uncertainties as stemming from Gaussian probability distributions, which is obviously an ad-hoc assumption, although consistent with combining various partial uncertainties in quadrature on the theory side, and in the experimental averages. In any case, the quoted confidence levels of our bounds should be taken with a grain of salt.

The left plot in Fig. 3 shows a Gaussian probability distribution for our average in Eq. (9). In Model-II, only enhancements of $R_{\gamma}$ with respect to the SM prediction (8) are possible. Thus, if there were no theoretical uncertainties, only $R_{\gamma}>3.31\times 10^{-3}$ would be accessible in Model-II. This is marked by the shaded regions (both light and dark) in the considered plot. The integrated probability over the dark-shaded region amounts to 5%. The border between the light- and dark-shaded regions corresponds to the central value for $R_{\gamma}$ obtained for $M_{H^{\pm}}\simeq 1276\,$GeV in the limit $\cot\beta\to 0$. Thus, one might expect that the $95\%\,$C.L. lower bound for $M_{H^{\pm}}$ should amount to $1276\,$GeV in the absence of theoretical uncertainties. We are not assuming here that Model-II is valid for sure. Instead, we are allowing for a possibility that it gets excluded (together with the SM) if $R_{\gamma}^{\exp}$ is sufficiently far below the SM prediction.

To include the theory uncertainties, one follows the standard confidence belt construction (see, e.g., Sec. 39.4.2.1 of Ref. [50]). For each $M_{H^{\pm}}$, one considers a Gaussian probability distribution around the theoretical central value, with its variance obtained by combining the experimental and theoretical uncertainties in quadrature. Next, a confidence interval corresponding to (say) $95\%$ integrated probability is determined. It can be placed either centrally (for a derivation of 2-sided bounds), or maximally shifted in either way (for 1-sided bounds), or in an intermediate way, like in the Feldman-Cousins (FC) approach [51]. This is illustrated in the right plot of Fig. 3, for Model-II with $\cot\beta\to 0$. The red, black, green, and blue intervals correspond to the 2-sided, upper 1-sided, lower 1-sided and FC cases, respectively. The experimental central value from Eq. (9) is marked by the vertical dashed line. On the vertical axis, we use $1\,{\rm TeV}/M_{H^{\pm}}$ that is restricted to be positive, which makes our case very similar to the one in Sec. IV-B of Ref. [51].

It is the freedom of the confidence interval placement that makes the resulting bounds on $M_{H^{\pm}}$ somewhat ambiguous. If we choose the FC (blue) intervals, low values of $R_{\gamma}^{\exp}$ can never lead to exclusion of Model-II in its whole parameter space. If we choose the upper 1-sided (black) intervals, our method is actually equivalent to using the experimental upper bound on $R_{\gamma}^{\exp}$ rather than the actual measurement. In this case, the previously discussed example from the left plot of Fig. 3 is recovered in the limit of no theory uncertainties. In the literature, bounds on $M_{H^{\pm}}$ have been derived using either the 1-sided (e.g., Refs. [7, 14, 26, 52]) or 2-sided (e.g., Refs. [6, 53]) approaches, and the method choice was not always explicitly spelled out.

In Tab. 3, we present the bounds we obtain following three different methods, and using three out of four777 We omit the one requiring the longest ($2.0\to 1.6$) extrapolation in $E_{0}$. averages for $R_{\gamma}^{\exp}$ from Tab. 2. The rows corresponding to our preferred choice (Eq. (9)) are displayed in bold. For Model-I we set $\tan\beta=1$, while the absolute bounds ($\cot\beta\to 0$) are shown for Model-II. In the Model-I case, the lower rather than the upper 1-sided intervals are employed.

It is interesting to observe that stronger bounds on $M_{H^{\pm}}$ in Model-II are found from the two less precise averages, just because their central values turn out to be lower. These averages are less sensitive to the $E_{0}$-extrapolation issues, which might be helpful in accepting the ones derived from Eq. (9) as conservative. The situation in Model-I is reverse – the most precise average gives the strongest bounds, as naively expected.

By coincidence, our 2-sided $95\%\,$C.L. bound of $583\,$GeV in Model-II practically overlaps with the $580\,$GeV one that has been obtained in Ref. [6] from their single measurement alone (giving ${\mathcal{B}}_{(s+d)\gamma}$ with a lower central value but larger uncertainty than the one corresponding to our Eq. (9)). Since this bound is also the most conservative one, we suggest choosing it for updated combinations with constraints from other observables.

In Fig. 4, the $95\%\,$C.L. bounds on $M_{H^{\pm}}$ are shown as functions of $\tan\beta$. Above $\tan\beta\simeq 2$, the Model-I bound becomes weaker than the LEP one ($\simeq 80\,$GeV [50]), while the Model-II one gets saturated by its $\tan\beta\to\infty$ limit ($\simeq 580\,$GeV). Our plot terminates on the left side at $\tan\beta=0.4$. For lower values of $\tan\beta$, the bound from $R_{b}$ becomes more important in Model-II (see Figs. 13 and 14 of Ref. [53]), while $R_{\gamma}$ alone in Model-I becomes insufficient due to possible changes of sign in the coefficient $C_{7}$. In the latter case, including the $b\to s\ell^{+}\ell^{-}$ observables becomes necessary – see Ref. [54] and references therein.

6 Conclusions

We derived updated constraints on $M_{H^{\pm}}$ in the 2HDM that get imposed by measurements of the inclusive weak radiative $B$-meson decay branching ratio. Although in principle straightforward, such a derivation faces several ambiguities stemming mainly from the photon energy cutoff choice. We presented an extended discussion of this issue, and updated the experimental averages. In Model-I, relevant constraints are obtained only for \tan\beta\;\hbox to0.0pt{\lower 3.5pt\hbox{\mathchar 0\relax\sim}\hss}\raise 1.0pt\hbox{<}\;2. In Model-II, the absolute ($\tan\beta$-independent) $95\%\,$C.L. bounds are in the $570$–$800\,$GeV range. We recommend one of the most conservative choices, namely $580\,$GeV, to be used for combinations with constraints from other observables. This value overlaps with the bound derived from the most recent single measurement alone [6].

Acknowledgments

We are grateful to Akimasa Ishikawa, Phillip Urquijo, Sławomir Tkaczyk and Aleksander Filip Żarnecki for helpful discussions, as well as to Paolo Gambino for reading the manuscript and useful comments. The research of M.S. has been supported by the BMBF grant 05H15VKCCA. M.M. acknowledges partial support from the National Science Centre (Poland) research project, decision no. DEC-2014/13/B/ST2/03969, as well as by the Munich Institute for Astro- and Particle Physics (MIAPP) of the DFG cluster of excellence “Origin and Structure of the Universe”.