Higgs decay into a lepton pair and a photon revisited

We present new calculations of the differential decay rates for $H\to \ell^+\ell^- \gamma$ with $\ell=e$ or $\mu$ in the Standard Model. The branching fractions and forward-backward asymmetries, defined in terms of the flight direction of the photon relative to the lepton momenta, depend on the cuts on energies and invariant masses of the final state particles. For typical choices of these cuts we find the branching ratios $B(H\to e \bar e \gamma)=6.1\cdot 10^{-5}$ and $B(H\to \mu \bar \mu \gamma)=6.7\cdot 10^{-5}$ and the forward-backward asymmetries $\mathcal{A}^{(e)}_{\text{FB}}=0.366$ and $\mathcal{A}^{(\mu)}_{\text{FB}}=0.280$. We provide compact analytic expressions for the differential decay rates for the use in experimental analyses.


I. INTRODUCTION
Since the discovery of a Higgs boson with a mass of 125 GeV [1,2] in 2012, the LHC experiments CMS and ATLAS put a major effort into the precise determination of its couplings. The Standard Model (SM) accomodates a minimal Higgs sector, with just one Higgs doublet, and it is natural to ask whether Nature foresees a richer Higgs sector than the SM. A possible imprint of an extended Higgs sector are deviations of the measured couplings from their SM predictions [3]. To date CMS and ATLAS have studied the couplings of the dicovered Higgs boson to W [4][5][6][7] and Z [7][8][9][10] bosons, τ leptons [11][12][13], b [14,15] and t [16,17] quarks, as well as photons [18,19]. The latter coupling is loop-suppressed and has been probed though the rare decay H → γγ. Rare Higgs decays are especially sensitive to physics beyond the SM and even probe scenarios with only one Higgs doublet as in the SM. For instance, H → γγ data were instrumental to rule out a fourth sequential fermion generation from a global analysis of Higgs signal strengths [20].
Phenomenologiocal analyses of two-Higgs-doublet models (2HDM) usually assume simple versions of the Yukawa sector (called type I,II,X, or Y), in which different observables become correlated and the largest imprints are on the heavy fermions of the third generation [21][22][23][24][25][26][27][28][29][30][31]. In such models the dynamics of light fermions of the first and second generation follow the pattern of the third generation and measurements of Higgs decay rates involving light fermions will only provide redundant information. However, as outlined in the following paragraph there are well-motivated phenomenological reasons to consider the possibility that physics beyond the SM shows imprints on the decays of the 125 GeV Higgs bosons into final states containing light fermions. It is therefore mandatory to measure the corresponding decay rates accurately and to compare the data with precise SM predictions.
In this paper we study the SM predictions for the rare decays H → + − γ with focus on = e and = µ. While the amplitude of H → + − is suppressed by one power of the Yukawa coupling y = m /v, where m is the lepton mass and v = 174 GeV is the vacuum expectation value (vev) of the Higgs field, the radiative analogue does not suffer from this suppression: Electroweak loop contributions instead involve the Higgs coupling to heavy gauge bosons or to the top quark and permit a non-zero decay amplitude even for y = 0, producing the lepton pair in a state with angular momentum j = 1. The decay rate Γ(H → e + e − γ) exceeds Γ(H → e + e − ) by far, while Γ(H → µ + µ − γ) and Γ(H → µ + µ − ) are comparable in size. In the case of Γ(H → τ + τ − γ) the electroweak loop contribution is much smaller than the tree-level contribution proportional to y 2 τ , which simply amounts to the bremsstrahlung contribution to H → τ + τ − . It is important to note that H → + − and H → + − γ probe different sectors of BSM models (chirality-flipping vs. chirality-conserving couplings to lepton fields) and are therefore complementary. A further motivation to study H → µ + µ − γ is the 3.7σ discrepancy between the measured anomalous magnetic moment of the muon, a µ , and its SM prediction [32]. a µ involves the magnetic operatorL µ Φσ αβ µ R F αβ , with the lepton doublet L µ = (ν µ , µ), the Higgs doublet Φ, and the electromagnetic field strength tensor F αβ . BSM models with loop contributions to the coefficient of the magnetic operator may as well affect the H-μ-µ-γ couplings. Another related topic are the hints of violation of lepton flavour universality encoded in the ratios R K ( * ) ≡ B(B → K ( * ) µ + µ − )/B(B → K ( * ) e + e − ), [33,34] which support BSM physics coupling to left-chiral leptons [35,36]. Also here the underlying BSM dynamics can eventually be tested with H → + − γ. None of the H → + − γ decays has been observed yet, cf. Ref. [37] for LHC limits. Analytic expressions for differential H → + − γ decay rates have been derived in Refs. [38] and [39].
After the discovery of the 125 GeV Higgs boson Ref. [38] was updated [40] and two new detailed analyses based on novel calculations of the decay rate have been presented in Refs. [42] and [43]. Comparing the decay rates dΓ(H → + − γ)/d √ s, where √ s is the invariant mass of the lepton pair, presented in these papers, we find significant discrepancies, which motivates the new calculation of this decay rate presented in this paper. We use a linear R ξ gauge, so that the we can use the vanishing of the W and Z gauge parameters as a check of our calculation. This check is especially valuable in the context of the Z width, which must be taken into account when the invariant mass of the lepton pair is close to the Z boson mass and special care is needed to ensure a gauge-independent result [44,45]. In this kinematic region our decay is indistinguishable from H → Zγ (for recent LHC search limits cf. [46,47]). In Ref. [42] it is argued that H → Zγ is not a properly defined physical process and one should instead discuss the full decay chain, including the Z decay, such as H → + − γ. Since we keep the gauge parameters arbitrary, we can track how the unphysical gauge-dependent pieces of the H → Z[→ + − ]γ sub-processe cancel with those of other diagrams. We will quote a compact formula for the differential decay rate d 2 Γ(H → + − γ)/(ds dt) with respect to the Mandelstam variables s and t, where √ t is the invariant mass of the ( − , γ) pair, and discuss the forward-backward asymmetry of the photon.
Our paper is organized as follows: In the following section we present our calculation and discuss our results, including a comparison with the literature. Sec. III contains our conclusions, followed by two appendices guiding through our analytic results.

A. Amplitudes
The amplitude for the tree-level photon emission process (see figure 1) is with our conventions for the kinematical variables as following: We denote four-momenta of photon, lepton and antilepton by k, p 1 , p 2 , respectively. Squared invariant masses are denoted by the Mandelstam variables where m H is the Higgs boson mass. e, m W , and θ W are the electromagnetic coupling constant, mass of the W boson, and weak mixing angle, respectively. u and v are the lepton and antilepton spinors and ε is the polarization vector of the photon.
In the one-loop contribution we can neglect y . The Feynman diagrams may be grouped into several classes as depicted in Fig. 2. Sample diagrams can be found in Fig. 3. The one-loop amplitude can be parametrized as (2) P L,R = (1 ∓ γ 5 )/2 are chiral projectors and the coefficients a 1,2 and b 1,2 are functions of s, t, u and the particle masses. We present the analytic results for a 1 and b 1 in Appendix A in Eqs. (A1) and (A2); a 2 , b 2 are obtained by interchanging t and u: a 2 (t, u) = a 1 (u, t) and b 2 (t, u) = b 1 (u, t). The compact results in Eqs. (A1) and (A2) involve the coefficient functions appearing in the Passarino-Veltman decomposition [48,49] of the tensor integrals. A result fully reduced to scalar one-loop functions [50,51]  We now comment on several differences with respect to existing results in literature. Refs. [39,52] contain additional terms of the form µνρσ k ρ (p 1 + p 2 ) σ ε ν * involving the Levi-Civita tensor in the final result for the loop amplitude. This contribution is absent in our result in equation (2) 1 . We have found that this term indeed appears in the top-quark triangle diagram (a) of Fig. 3, but cancels with the corresponding diagram with opposite fermion number flow. We have used 't Hooft-Veltman scheme for the treatment of Using a nonlinear gauge [53] the authors of Ref. [38] have identified classes of diagrams which separately satisfy the electromagnetic Ward identity. Using instead the usual linear R ξ -gauge we find a straightforward cancellation of the Z boson gauge parameter ξ Z , while the cancellation of the dependence on the Ref. [42]. We instead start with strictly real gauge boson masses, verify the ξ W independence of the result, and subsequently add the finite Z width Γ Z to the final, gauge-independent result.
In most phase space regions our loop functions are real; exceptions are kinematical situations such as t > M 2 W or u > M 2 W permitting on-shell cuts of the loops. Switching to the complex mass scheme makes the real loop expressions develop imaginary parts proportional to the gauge boson width; in the corresponding phase space regions the decay rates found in the two approaches differ by terms quadratic in Γ W,Z . We checked that the difference in dΓ(H → + − γ)/d √ s between the two approaches is numerically negligible. Therefore the treatment of Γ Z cannot be the reason for the numerical differences between our result and the various results in the literature.
We next shortly describe the tools used in our calculation. We have generated the Feynman diagrams with the FeynArts package [54]. For the evaluation of the loop integrals and the reduction to scalar basis  The tree-level contribution for the case of electrons is negligible. The only cut which we impose for these plots is E γ, min = 5 GeV which merely amounts to a lowering of the maximum value of m .
The tree contribution to the decay rate shown in Fig. 1 reads with The contribution from the one-loop diagrams is While we set m to zero in Eq. (5), we retain a nonvanishing value of m in the kinematical limits of the phase-space intergration. The limits for the variables s and t can be expressed as The tree-level contribution exhibits an infrared pole as s approaches its maximum value s max = m 2 H that corresponds to a vanishing photon energy E γ in the Higgs boson rest-frame, E γ = (m 2 H − s)/2m H . For the evaluation of the total decay rate we impose the cut E γ, min which lowers the maximum value of s to s cut = m 2 H − 2m H E γ, min . The resulting differential decay rate over m = √ s, for = e, µ, is shown in Fig. 4. One notes the enhancement from the Z-pole as well as the tail of photon pole starting at m , min .
The contribution of the interference between tree-level and one-loop contributions is negligible, as well as the effect of the tree-level diagrams in the case of electrons. For the evaluation of the full decay rates we employ the kinematical cuts of Ref. [40,42], namely: The results do not change if E 1 and E 2 are interchanged.
We use the following input for for the physical parameters: The difference between Γ (e) and Γ (µ) stems from the tree-level contribution. Fig. 5 shows the differential decay rate with respect to the invariant mass of the lepton-photon pair. With a total Higgs width of 4.1 MeV the rates in Eq. (9) correspond to the branching ratios B(H → eēγ) = 6.1 · 10 −5 , B(H → µμγ) = 6.7 · 10 −5 .
These branching ratios are roughly three times smaller than B(H → µμ).

C. Forward-backward asymmetry
Here we present the differential decay distribution with respect to cos θ ( ) , where θ ( ) is the angle between lepton and the photon in the rest frame of the Higgs boson, t = E γ (E 1 − | p 1 | cos θ ( ) ). The resulting distribution for the case of = µ is shown in Fig. 6. For this evaluation we apply the cuts m µµ > 0.1m H and E γ > 5 GeV (and no cuts on E 1,2 ). We define the forward-backward asymmetry with respect to θ ( ) as With the cuts m > 0.1m H and E γ,min = 5 GeV, applied to both the numerator and the denominator in Eq. (11), we obtain the numerical values: We note that the contribution of the tree-level diagrams to the numerators of the asymmetries in Eq. (11) can be neglected relatively to the dominant one-loop contribution and the numerators for the electron and muon case are essentially identical. Thus the difference between A FB is numerically more pronounced than the one between Γ (e) and Γ (µ) in Eq. (9), because different cuts are used in Eqs. (9) and (12).

D. Comparison with previous results
The main goal of our paper is the resolution of the discrepancies between the different results in the literature. Only Abbasabadi at al. [38] and Chen at al. [39] provided analytic result. The latter paper is the only one containing terms with the Levi-Civita tensor and the origin and cancellation of such terms is  [42] and [40] are denoted by red short-dashed and blue long-dashed lines, respectively. discussed above in Sec. II A. In the case of the former paper we have numerically evaluated the presented formula (taking into account the typo reported in Ref. [40]) and only find quantitative agreement in some regions, while we significantly disagree in others. In the next step we have digitalized the plots for dΓ(H → + − γ)/dm ee of Refs. [42] and [40], which have used the same cuts on the kinematic variables. We compare the two results and ours in Fig. 7. From this plot we can see that the difference between the previous works in Refs. [40,42] is up to 30%. Our result is close to the one of Ref. [42] for m ee 40 GeV, but significantly deviates for smaller values of m ee . One may speculate that the choice of the QED fine structure constant α, which fixes the e 4 term in the overall normalization constant N in Eq. only alleviates the tension in Fig. 7 for m ee 40 GeV, but does not fully resolve it. Furthermore, the agreement with the total decay rate becomes worse. Furthermore the shape of the distributions is different, and eventually the numerical integration over t is the reason for this discrepancy. We have cross-checked our result by using different Monte-Carlo generators, i.e. Vegas and the one implemented in Mathematica. Our result for the integrated rate Γ (e) = 0.249 keV is in reasonable agreement with the result Γ (e) = 0.233 keV given in Ref. [42]. We further remark that we disagree with Ref. [42] in the tree-level contribution to the integrated deacy rate of H → µμγ amplitude by a factor of 2. Reference [43] further presents results for dΓ(H → + − γ)/dm ee for a different choice of cuts, namely: where ∆R γf = (∆η 2 + ∆φ 2 ) 1/2 denotes the rapidity-azimuthal angle separation. Digitalizing the plot in this paper as well, we compare the presented results with ours in Fig. 8. We have found that the effect of the cuts in equation (13) does not alter the the loop-induced distribution by more than 2% in the region where we observe deviations from Ref. [43]. Therefore, we add in the same figure the digitalized result of Ref.
[40] that does not employ any cut. We observe good agreement between Refs. [40,43] in a region below the Z peak, where we agree with these results well for m µµ 70 GeV, while deviating otherwise. For m µµ > M Z we agree well with the result of Ref. [43].

III. CONCLUSIONS
The results in the literature for the differential decay rates dΓ(H → + − γ)/dm with = e, µ, differ substantially. We have performed a new calculation of the differential decay rate d 2 Γ(H → + − γ)/(ds dt), where s is the squared invariant mass of the lepton-antilepton pair and t is the corresponding quantity for the lepton-photon pair. We have performed our calculation in an R ξ gauge and have verified the gauge independence of the result. After presenting various differential decay distributions we have studied the forward-backward asymmetry defined in terms of the flight direction of the photon with respect to the lepton. These asymmetries, quoted in Eq. (12), are sizable.
In experimental studies one defines cuts in the laboratory reference frame rather than the Higgs rest frame and the comparison between data and SM prediction requires the use of the fully differential decay rate. For this purpose we present an analytic expression in a compact form in Appendix A and provide ancillary files for the use by experimental collaborations.