Trilinear Higgs coupling determination via single-Higgs differential measurements at the LHC

We study one-loop effects induced by an anomalous Higgs trilinear coupling on total and differential rates for the H→4ℓ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$H\rightarrow 4\ell $$\end{document} decay and some of the main single-Higgs production channels at the LHC, namely, VBF, VH, tt¯H\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$t{\bar{t}}H$$\end{document} and tHj. Our results are based on a public code that calculates these effects by simply reweighting samples of Standard-Model-like events for a given production channel. For VH and tt¯H\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$t{\bar{t}}H$$\end{document} production, where differential effects are particularly relevant, we include Standard Model electroweak corrections, which have similar sizes but different kinematic dependences. Finally, we study the sensitivity of future LHC runs to determine the trilinear coupling via inclusive and differential measurements, considering also the case where the Higgs couplings to vector bosons and the top quark is affected by new physics. We find that the constraints on the couplings and the relevance of differential distributions critically depend on the expected experimental and theoretical uncertainties.


Introduction
Since its discovery in 2012 [1,2], evidence has been steadily accumulating that the properties of the scalar particle at 125 GeV correspond to those of the Higgs boson predicted by the Standard Model (SM) of elementary particles and interactions. ATLAS and CMS experiments have analyzed data from several inverse femtobarns of integrated luminosity at different energies, providing already rather precise measurements of the Higgs couplings to the vector bosons and to the fermions of the third generation [3][4][5][6][7]. Prospects of the next LHC runs for further improving the precision on the value of a e-mail: fabio.maltoni@uclouvain.be b e-mail: davide.pagani@tum.de c e-mail: ambresh.shivaji@uclouvain.be d e-mail: xiaoran.zhao@uclouvain.be these couplings and for measuring also the couplings to the second-generation fermions are very good. On the contrary, the situation and prospects for determining the properties of the scalar potential, i.e., of the Higgs self-couplings at the LHC are less clear and therefore the subject of an intense theoretical and experimental activity.
At low energy, the potential for a scalar particle of mass m H can be parametrized as a polynomial, where v = ( √ 2G F ) −1/2 ∼ 246 GeV is the vacuum expectation value after electroweak-symmetry-breaking (EWSB). In the SM, renormalizability and gauge invariance dictate that the Higgs potential depends only on two parameters, μ and λ, where is the Higgs doublet. After EWSB the potential in Eq. (2) gives rise to the mass of the physical Higgs boson m H , which, together with the vacuum expectation value v, are related to μ and λ via μ 2 = m 2 H /2 and λ = m 2 H /(2v 2 ). As a result, by fixing μ and λ, in the SM the Higgs selfcouplings are completely determined, leading in Eq. (1) to On the other hand, new physics could modify the Higgs potential at low energy, by altering the value of λ 3 (or in general λ i for the i-point Higgs self-couplings) without affecting the value of m H and v. This can be realized either directly (e.g. by extending the scalar sector) or indirectly (due to the exchange of new virtual states). In addition, modifications in the self-interactions would be induced if the Higgs boson is a composite state.
Since double Higgs production directly depends on the Higgs trilinear coupling at LO, it is the standard process for studying λ 3 at the LHC. However, the cross section of its main production channel, the gluon fusion, is only about 35 fb at 13 TeV [8][9][10], so it is much smaller than single Higgs production cross section, which is about 50 pb [11]. Several phenomenological studies have been performed on the determination of λ 3 via the relevant experimental signatures emerging from this process: bbγ γ [12][13][14][15][16][17], bbτ τ [13,18], bbW + W − [19] and bbbb [20][21][22]. Also tt H H [23,24] and V H H [25] production processes have beeen considered. Nevertheless, given the complexity of a realistic experimental set-up, the final precision that could be achieved on the determination of λ 3 is still unclear. On the contrary, it is established that at the LHC perspectives of inferring information on λ 4 from the triple Higgs production are quite bleak [26,27], due to the smallness of the corresponding cross section [8,28], together with a rather weak dependence on this parameter. Even at a future 100 TeV proton-proton collider a considerable amount of integrated luminosity will be necessary in order to obtain rather loose bounds [29][30][31].
At the moment the strongest experimental bounds on nonresonant double-Higgs production have been obtained in the CMS analysis of the bbγ γ signature [32], where cross sections larger than about 19 times the predicted SM value have been excluded. However, exclusion limits on λ 3 are found to be strongly dependent on the value of the top-Higgs coupling and they are of the order λ 3 < − 9 λ SM 3 and λ 3 > 15 λ SM 3 for the SM-like case. These new limits, together with the slightly weaker ones from the ATLAS bbbb [33] and CMS bbτ τ [34] analyses at 13 TeV, improve the results from 8 TeV data (λ 3 < − 17.5 λ SM 3 and λ 3 > 22.5 λ SM 3 [35]), however, also with a high integrated luminosity (HL) of 3000 fb −1 , a further improvement may not be so tremendous. The most optimistic experimental studies for HL-LHC suggest that it could be possible to exclude values in the range λ 3 < − 1.3 λ SM 3 and λ 3 > 8.7 λ SM 3 via the bbγ γ signatures [36]. Additional and complementary strategies for the determination of λ 3 are thus desirable at the moment.
In Ref. [37] an indirect method of measuring λ 3 via EW radiative corrections in e + e − → Z H process was proposed. Recently, the same idea has also been extended to the LHC, by (globally) studying λ 3 -dependent EW corrections in single Higgs production and decay processes [38][39][40][41]. Moreover, limits on λ 3 can also be derived by two-loop effects in EW precision observables [42,43] such as the measurements of m W and of the S, T oblique parameters. These studies have confirmed that indirect bounds on λ 3 can be competitive with the direct ones inferred from the di-Higgs production channel. For example, a simple one-parameter fit to the signal-strengths' measurements at 8 TeV [7] gives − 9.4 λ SM 3 < λ 3 < 17.0 λ SM 3 [39], comparable to the current best constraints from the bbγ γ CMS measurement mentioned above. In both analyses no other deviations for the Higgs couplings are considered. The usefulness of a joint analysis of λ 3 indirect effects on single-Higgs production and direct effects on di-Higgs production has been discussed and quantified in Ref. [41]. As already suggested in Refs. [39,40], the role of differential distributions and their non-flat dependence on λ 3 has been proved to be crucial in Ref. [41], especially when not only anomalous λ 3 effects but also modifications of the couplings to the other particles are considered, as expected in a general effective-field-theory (EFT) approach.
The purpose of this work is threefold. First, we present an automated public code for generating events including λ 3 effects at one loop, thus allowing the study of differential effects in VBF, V H and tt H production 1 and all the relevant Higgs decays. The code is based on two independent and procedurally different implementations in the Mad-Graph5_aMC@NLO framework [44]. 2 Second, with the help of this code, we extend at the differential level the results of Ref. [39], where all the relevant single Higgs production (ggF, VBF, V H, tt H) and decay channels (γ γ , V V * , f f , gg) have been analyzed and included in a global fit only at the inclusive level. Indeed, in Ref. [39] the usefulness of differential distributions has already been explored but only for the case of V H and tt H production, providing the results on which the analysis in Ref. [41] rely. Differential information for VBF and V H production has been presented also in Ref. [40]. Here we scrutinize all the relevant distributions that are potentially affected by anomalous λ 3 effects, presenting for the first time detailed results at the differential level for tt H production, for the H → 4 decay and also for the t H j process, for which also inclusive results are new. Moreover, for the case of V H and tt H production, where loop-induced λ 3 effects are not flat, we repeat the calculation of NLO EW corrections [45][46][47][48][49][50] in the SM, which are also not flat, in order to check the robustness of our strategy. As expected, we find that NLO EW corrections are essential for a precise determination of anomalous λ 3 effects, but also that they do not jeopardize the sensitivity of indirect λ 3 determination. We use for the calculation the EW extension of the automated MadGraph5_aMC@NLO framework that has already been used and validated in Refs. [48,49,[51][52][53][54]. 3 Finally we perform a fit for λ 3 based on the future projections of ATLAS-HL for single-Higgs production and decay at 14 TeV [55,56]. We consider the effects induced on the fit by additional degrees of freedom, namely anomalous Higgs couplings with the vector bosons and/or the top quark. We investigate the impact on the fit of three different factors: the differential information, the experimental and theoretical uncertainties, and the inclusion of the two aforementioned additional degrees of freedom. We find that in a global fit, including all the possible production and decay channels, two additional degrees of freedom such as those considered here do not preclude the possibility of setting sensible λ 3 bounds, especially, they have a tiny impact on the upper bound for positive λ 3 values. On the contrary the role of differential information may be relevant, depending on the assumptions on experimental and theoretical uncertainties.
The structure of the paper is the following. In Sect. 2 we briefly repeat and comment on the main formulas that are relevant in the study of one-loop induced λ 3 effects. Differential results for the various processes that have been mentioned before are given in Sect. 3. In Sect. 4 we present the extension of the analysis of λ 3 effects including NLO EW corrections and we study the impact on differential distributions and inclusive rates. In Sect. 5 we introduce the framework for the fit and discuss the results obtained. Details as regards the statistical treatment of uncertainties in the fit are collected in Appendix A.

Technical set-up
2.1 Self-couplings effects in single Higgs production and decays at one loop Higgs self-couplings can be studied in a model-dependent approach, e.g., choosing a particular UV-complete scenario, or in a model-independent approach, as done in this work. However, there are different ways in which the modifications of trilinear and quartic couplings can be parametrized and they rely on different theoretical assumptions. If new physics is at scales sufficiently higher than the energies where measurements are performed, the SMEFT offers a consistent and model-independent way of organizing generic deformations to the Higgs interactions. Moreover, radiative corrections can be consistently performed within this framework. However, it has been shown [39] that in the case of the trilinear coupling and at the order we are considering, one-loop corrections for single Higgs processes, adding higher-dimensional operators that only affect the Higgs self-couplings (the ( † ) n operators with n > 2) or directly introducing an anomalous coupling, are two fully equivalent approaches for deforming the SM Higgs potential. In other words, regarding the Higgs selfcouplings, differences between an EFT and an anomalous coupling parametrization will arise only for final states fea-turing more Higgs bosons and/or at higher loop level, i.e., with the appearance of higher-point interactions (starting from quadrilinear). It is essential to note that the single Higgs production and the decay channels are not sensitive to λ 4 at one loop. For this reason, although results are written in terms of λ 3 , they can easily be translated in terms of the Wilson coefficient in front of the dimension-6 operator ( † ) 3 , or those for higher-dimension ( † ) n operators. In the following we briefly repeat and comment the main formulas that have been introduced and discussed at length in Ref. [39]; very minor modifications are present in the notation and in the definition of the corresponding quantities.
The λ 3 -dependent part of the NLO EW corrections to single Higgs processes are gauge invariant and finite. These contributions can be organized in two categories: a universal part proportional to (λ 3 ) 2 , which arises from the wavefunction renormalization of external Higgs boson and thus does not depend on the kinematics, and a process-dependent part linear in λ 3 , which is also sensitive to the kinematics. In the presence of modified trilinear coupling, the master formula for the λ 3 -dependence of a generic observable (total/differential cross section or decay width) can be written as where C 1 is the process-and kinematic-dependent component and LO stands for the LO prediction including any factorizable higher-order correction. In particular, we assume that QCD corrections do in general factorize, which, for the V H and VBF case, has been shown to be a correct approach up to NNLO in Ref. [40]. 4 Representative diagrams contributing to the C 1 for the different processes are depicted in Fig. 1 The This is particularly convenient for the discussion in Sect. 4, where we will analyze NLO EW corrections in the SM in conjunction with λ 3 -induced effects. In conclusion, the relative corrections due to the trilinear coupling can be expressed as which manifestly goes to zero in the κ 3 → 1 limit.
Numerical values of C 1 at the inclusive level for the processes considered in this work are reported in Table 1. The calculation of C 1 for single-top-Higgs production, which appears for the first time here, is non-trivial and discussed in Sect. 3.4. The range of validity of Eq. (9) has been identified in Ref. [39] as |κ 3 | < 20, given the values of δ Z H and C 1 in Table 1. As we will see, at the differential level this limit may be too loose since C 1 can receive large enhancements (see Sect. 3.3). On the other hand, we believe that the constraint |κ 3 | 6 identified in Ref. [57] is appropriate for inclusive double Higgs production, but it is too strong for the case of single-Higgs production. Indeed the violation of perturbativity for the H H H vertex is kinematic dependent and the condition |κ 3 | 6 arises from the configuration with two H bosons on-shell and the third one with virtuality slightly larger than 2m H . This is the kinematic configuration present above the threshold in double Higgs production, where the bulk of its cross section comes from, but is never present in single Higgs production, since only one Higgs boson can be on-shell in the H H H vertex appearing at one loop. While δ Z H is a universal quantity, C 1 is process and kinematics dependent. We have employed two independent methods for the computation of C 1 for differential cross sections and decay rates. They correspond to the two different codes publicly available, which agree within their numerical accuracy.
In the first we parametrize the finite one-loop corrections due to the trilinear Higgs self-coupling as form factors that are functions of the external momenta. One-loop integrals are computed using the LoopTools package [58] and the form factors are implemented as effective new vertices in a dedicated UFO model file [59], which is then used in Mad-Graph5_aMC@NLO [44]. As a result, parton-level events can be generated including O(λ 3 ) effects and any interesting observable analyzed. Our current implementation of form factors allows the computation of differential C 1 for VBF, V H and H → 4l processes. At the order of accuracy of our calculation, ggF production and all the other 1 → 2 decays do not have a kinematic dependence for C 1 ; results at the inclusive level are sufficient for any kinematic configuration and thus taken from Ref. [39].
On the other hand, the implementation of form factors for tt H and t H j processes would be quite cumbersome as there are many one-loop integrals that contribute. A different strategy, based on reweighting, has therefore been devised. With this second method, one starts by generating a sample of (unweighted) parton-level events at leading order. These events are then used as input for a code that computes the momentum-dependent weight where, following the notation of Ref. [39], M 0 refers to the tree-level amplitude and M 1 λ SM 3 to the SM virtual corrections depending on λ 3 . Then LO events are reweighted by multiplying the weight of each event i by the corresponding w i . In this way, also with this method it is possible to calculate C 1 for any differential distribution. The required one-loop matrix elements are computed using the capabilities of MadGraph5_aMC@NLO for evaluating loop diagrams [48,53,60,61]. For each process, we use diagram filters in order to select the relevant one-loop matrix elements featuring the Higgs self-coupling. We find this method much faster and efficient than the one based on form factors also for V H and VBF processes, thus we actually employ it for deriving all the results presented in this work and we suggest the usage of this version of the code. Note, however, that the other method offers at least in principle the possibility of explicitly including NLO QCD corrections on top of λ 3 -induced effects for VBF and V H production.

Results for differential distributions
Results are presented in this section and have been obtained with the following input parameters: which are taken from Ref. [62]. We use as PDF set the PDF4LHC2015 distributions with the factorization scale at where m(i) are the masses of the particles i in the final state. 5 In the following subsections we provide differential results for various relevant observables in VBF, V H, tt H and t H j production channels and in the H → 4l decay channel. Each plot has the layout that is described in the following. The upper panel displays the LO distribution (red) and O(λ 3 ) corrections alone (blue), both normalized by their value for the total cross section. In other words, we compare the shape of LO distributions with the shape of the contributions induced by C 1 in Eq. (5), which is thus independent on the value of κ 3 . The lower panel display C 1 both at differential level (green) and for the total cross section/decay (blue). The latter values are also summarized in Table 1 and will be used in Sect. 5 for the representative fit results.

VBF
Vector boson fusion is generated by requiring EW production of Higgs plus two jets, which includes also V H configurations with the vector boson V decaying into two jets. We effectively eliminate V H contributions by applying the following kinematic cuts [62] on the two final-state jets: In Fig. 2, we present C 1 for representative distributions, namely, p T (H ), p T ( j 1 ), m( j j) and m(H j j). In fact, we have checked that similar effects characterize other observables, which, however, we do not show. As already noticed in Refs. [39,40] the value of C 1 is not particularly large and rather flat for all the distributions shown here; C 1 = 0.63% for the total cross section and never exceeds 0.70% at the differential level. At variance with the case of V H and tt H considered in the following, loop corrections featuring trilinear Higgs self-couplings involve Higgs propagators connecting the final-state Higgs and internal V propagators. Thus, no Sommerfeld enhancement is present at threshold. In this respect, the interest of VBF for what concerns the indirect determination of λ 3 is mostly limited to the shift in the total rate, which, even though modest, is anyway relevant. Indeed, VBF is the channel with the second largest cross section and the smallest of the theoretical uncertainties [62], as can also be found in Table 4 in Appendix A.

VH
In Figs. 3 and 4 we show the differential C 1 for Z H and W H(W = W + , W − ), respectively. As discussed in Refs. [39,40]   The latter is particularly interesting because C 1 is enhanced w.r.t. the inclusive case in the region corresponding to the largest cross section.
We also looked at possible effects due to the Z polarization or, in other words, measurable via the angular distributions of Z (and H ) decay products. We did not see any enhancement or shape dependence for these distributions. Furthermore, one should bear in mind that also the loop-induced gg → H Z process gives a non-negligible amount of the NNLO cross section, order ∼ 1/6 at 13 TeV. This process also has a dependence on λ 3 , but only at two-loop level, and should exhibit a shape dependence. However, this calculation is not technically feasible yet.

ttH
Together with gluon-fusion production, the tt H channel plays a major role in providing information of the top-quark couplings to the Higgs. Its importance can be gauged by sim-ply considering its weight in a global κ-framework fit [63] or in the SMEFT framework [62]. The same importance should be ascribed to this process also from the point of view of the sensitivity to λ 3 : C 1 for tt H is the largest among all production channels and with the most significant kinematic dependence [39]. In Fig. 5, we show the most important kinematic distributions in this channel. C 1 for total cross section is 3.52% and can increase up to ∼ 5% in p T distributions. Similarly, with the binning chosen in Fig. 5, C 1 for the invariant mass distributions can be as large as 10% close to threshold, even though, once again, in the same region the cross section is suppressed by phase space. The origin of the large phase-space dependence of C 1 is again due to Sommerfeld enhancements in the threshold regions that are induced by interactions among the top (anti)quark and the Higgs boson.

tHj
Although it is characterized by a rather small cross section at the LHC, single associate production of a Higgs with a  single top is a particularly rich and interesting process, especially in searching for observables sensitive to relative phases among the Higgs couplings to fermions and bosons [64][65][66][67]. Naively, one would expect this process to have a sensitivity to the trilinear one between that of VBF and tt H; the t H j process features a top quark in the final state as well as W boson(s) in the propagators. The contribution of one-loop diagrams featuring the Higgs self-coupling to this process has not been considered in Ref. [39] for two major reasons. The first one was of phenomenological nature: in the SM this process is barely observable at the Run II of the LHC. The second one is of a technical nature: the calculation needs a careful check of EW gauge invariance and UV finiteness, since a few subtleties, which are not present for the other processes discussed in this work, arise. We describe them in the following. Similar to the case of the H → γ γ decay [38,39], Goldstone bosons appear in the Feynman diagrams contributing to the LO. Thus, H GG and H H GG interactions are present in one-loop EW corrections. While the former is not modi-fied by ( † ) n effective operators, the latter is indeed modified [38,39]. The calculation can be consistently performed in two different ways: either directly eliminating Goldstone bosons by employing the unitary gauge, as also done for other quantities in Refs. [39,42], or keeping track of H H GG effects in the intermediate calculation steps, as we explain in the following and as we actually will do in our calculation.
In a generic gauge, the on-shell renormalization of the EW sector [68] involves the counterterm for the Goldstone self-energy, which depends on the Higgs tadpole counter term δt, which in turn depends on the trilinear coupling λ 3 . Therefore, if we only modify the value of λ 3 , the Goldstone self-energy counterterm receives a UV-divergent contribution proportional to (κ 3 − 1), which is not cancelled by any divergence from loop diagrams. Instead, if we consistently take into account the modification of the H H GG vertex, loop diagrams featuring a seagull in the G propagator are also present; they exactly cancel the UV-divergent contribution proportional to (κ 3 − 1) in the Goldstone self-energy counter term, leading to the same result one would obtain in In our results we include both t H j andt H j channels and we do not apply cuts on the jet, since the result is infrared finite. We find the C 1 for the total cross section is about 0.91%. In Fig. 6, we show C 1 for kinematic distributions such as p T (H ), p T (t), m(t H) and m(t H j). We note that unlike the other variables p T (t) does not decrease monotonically as we move from low to high p T values. Near threshold m(t H) displays a quite impressive difference in shape.

H → 4
The Higgs decay into four fermions is the only Higgs decay channel with non-trivial final-state kinematics. Moreover, it is the only one where a priori also C 1 can have a shape dependence. Indeed, all the other decays correspond to a 1 → 2 process, and since the H boson is a scalar, there is not a preferred direction in its reference frame. In the previous study [39] the C 1 for H → Z Z * decay was calculated to be 0.83%. Although the full off-shell configuration was taken into account, possible angles between the decay products were not analyzed. Using the form-factor code mentioned above we calculate C 1 for H → e + e − μ + μ − channel. We analyzed C 1 for many observables involving the four leptons, but we found that it has in general almost no kinematic dependence. As an example, in Fig. 7, we display C 1 for leading and subleading lepton pair invariant masses. Since the Higgs boson interactions with the final-state fermions are negligible, this result can be extended to all the other decays into four leptons and in general into four fermions.

Anomalous trilinear effects and the NLO electroweak corrections
The set of one-loop corrections to single Higgs production and decays involving the trilinear Higgs self-coupling is gauge invariant and finite. However, in the SM, performing a perturbative expansion in power of α s and α, other contributions are present at the same order of accuracy. In other words, λ 3 -induced effects at one-loop level should be considered as a gauge-invariant and finite subset of the complete NLO EW corrections, which also include effects form virtual W, Z and photons as well as real emission contributions. 6 As shown in the previous section, the possibility of measuring anomalous λ 3 effects via precise predictions in single Higgs production relies both on the precision of experimental measurements and SM theory predictions. In particular, regarding the theoretical accuracy, while it is reasonable to assume that QCD corrections in general factorize λ 3 effects, as explained in Sect. 2, this is in general not true for NLO EW corrections. The purpose of this section is to provide a consistent extension of the master formula in Eq. (5) that includes also NLO EW corrections and to investigate their impact in the determination of anomalous λ 3 effects. All the calculations of the NLO EW corrections in the SM, with the exception of ggF taken from [62], are performed in a completely automated approach via an extension of the MadGraph5_aMC@NLO framework that has already been used and validated in Refs. [48,49,[51][52][53][54]. Concerning the renormalization, we use the G μ -scheme, consistently with the input parameters listed in Sect. 3.
At the differential level we limit ourselves to the study of the V H and tt H processes, where the C 1 dependence on the kinematics is large. We have also computed the differential EW corrections to the t H j production channel, but we do not report plots here, as its phenomenological relevance will be marginal at 13 TeV LHC for our purposes. The differential case is of particular interest since EW corrections in the SM, due to the Sudakov logarithms, are large in the boosted regime, i.e., exactly in the opposite phase-space region where λ 3 -induced effects are sizable, the production threshold, as already discussed before.
The master formula in Eq. (5) can be improved including NLO EW correction in the following way: where δ EW λ 3 =0 represents the part of the NLO EW K -factor in the SM 6 In the EW sector of the SM all the interactions are determined by the mass of the fermions, m H and three additional parameters, which are typically m W , m Z and α or G μ . In general, it is not possible to alter at NLO EW accuracy a derived quantity, such as λ = m 2 H /(2v 2 ), without spoiling the renormalizability of the theory. The special case of λ 3 in single Higgs production at one loop has been discussed in detail in Refs. [39,42].
that does not depend on λ 3 , namely, 7 In Eq. (14), SM NLO stands for the observable at LO + NLO EW accuracy. Thus, in the limit λ 3 → 1, BSM NLO → SM NLO . As can be noted, the Z BSM H term factorizes the NLO EW contributions in the SM, while C 1 does not. Indeed, in general, EW loop corrections on top of λ 3 -induced effects need a dedicated two-loop calculation and a full-fledged EFT approach in order to obtain UV-finite results; only the Z BSM H contribution is completely model-independent and factorizes the NLO EW corrections in the SM. However, it is worth to note that, assuming factorization also for the C 1 contributions, terms of the order κ 3 C 1 × δ EW λ 3 =0 would be anyway negligible, since either δ EW λ 3 =0 (Sudakov logarithms in the boosted regime) or C 1 (Sommerfeld enhancement in the threshold region) is sizable, but never both of them at the same time. This will be clear in the differential plots we display in the following.
The EW K -factor at the inclusive level can be found for all processes in Table 2, while relevant differential results for Z H, W H and tt H are displayed in Figs. 8, 9 and 10, respectively. In each figure, plots on the left show the p T (H ) distributions, while plots on the right those for the invariant mass of the final state. In the upper plots we display the ratio (σ BSM NLO − σ LO )/σ LO for different values of λ 3 , (−10, 0, 1, 2, 10). In practice, the case λ 3 = 1, directly denoted as SM in the plots, corresponds to the differential (K EW −1) in the SM. The lower plots display the ratio σ BSM NLO /σ NLO (solid lines) and the term 1 + δ κ 3 (dashed lines) for different values of λ 3 , (−10, 0, 1, 2, 10). In practice, the former is our prediction at NLO EW accuracy for the signal strengths μ i that will enter in the fits of the next section, 8 the latter is the definition at LO used in previous work.
First, we comment on the shape of the NLO EW corrections in the upper plots. The general trend in the SM is characterized by large negative Sudakov logarithms in the tails, 7 Here, in order to keep the notation simple, with the symbol δ Z H we still refer to only the λ 3 contributions to the Higgs-wave-function counterterm. Thus, δ EW λ3=0 contains further contributions to the Higgswave-function counterterm that do not depend on λ 3 . 8 The signal strengths μ i is better defined afterwards in Eq. (17). In this section κ i = 1. By looking at the upper plots in Figs. 8, 9 and 10 is evident that EW corrections have to be included in order to correctly identify anomalous λ 3 effects. On the other hand, NLO EW corrections do not largely affect the value of the signal strengths, i.e. the ratio of the BSM and SM prediction. This fact can be seen in the lower plots, where we display σ BSM NLO /σ NLO (solid lines) and 1 + δ κ 3 (dashed lines), which indeed corresponds to the aforementioned ratio with or without NLO EW corrections both in the numerator and denominator. 9 Solid and dashed lines are in general very close, especially for small values of κ 3 . It is interesting to note that for very small values of m(tt H) a value κ λ = −10 would lead to corrections that are negative and larger in absolute value than the LO. This is due to the very large C 1 (see Fig. 5) and points to the necessity of including higher-order λ 3 -induced effects for large values of κ 3 and for this specific phase-space region.
For the case of total cross sections, we plot, in Fig. 11, the ratio σ BSM NLO /σ SM NLO as a function of κ 3 in the range [−10, 10] for all the single Higgs production processes, including also t H j. Differences with the corresponding 1 + δ κ 3 ratios, 9 Actually, the ratio without NLO EW corrections should be σ BSM λ3 /σ SM λ3 , but its difference with 1 + δ κ3 , which is used in previous work and more useful for a direct comparison, is negligible. which have been presented in Ref. [39], are hardly visible and thus we do not show them here.
In conclusion, constraints on λ 3 from a global fit based on the value of the signal strengths μ i at inclusive [39] or differential level [40] will not be affected by NLO EW corrections. On the other hand, in the experimental analyses EW corrections have to be taken into account, especially at the differential level, for the determination of the value of the signal strengths, which is in general important for any BSM study and not peculiar for our case.

Constraining κ 3 through a global fit
In this section we discuss the role of differential distributions in the determination of λ 3 via single-Higgs production and decays measurements. The aim of this section is threefold. First, we show that the interplay of theoretical and experi-mental uncertainties have a large impact in the determination of the constraints on λ 3 , especially when differential information is exploited. Second, we discuss how bounds on λ 3 are affected by the presence of additional anomalous couplings in the fit. In particular, we progressively lift the assumption that the Higgs couplings to the top quark and to vector bosons are SM-like. Third, we include the EW effects discussed in the previous section in the fit analysis, providing consistent formulas for repeating the fit in conjunction with additional Higgs anomalous couplings. The reader who is only interested in the results may wish to go directly to Sect. 5.2 and skip Sect. 5.1, where formulas are given and a few technical details are discussed.

Combined parameterization of κ 3 , κ t and κ V effects
In order to parametrize the Higgs anomalous couplings to the top quark and to the vector bosons we use the coupling modi- fiers κ t and κ V , respectively (see Ref. [7] for definitions). We are interested in how additional BSM effects entering at LO may alter the determination of κ 3 and the relevance of differential measurements. The choice of the (κ 3 , κ t , κ V ) kappa framework is driven by simplicity; our main purpose is to add new degrees of freedom in the fit and identify in which configurations the differential information may be particularly relevant. On the other hand, while the cases of new physics entering only via κ 3 and a very general EFT parametrization (ten independent parameters) have already been explored [41], simplified intermediate parameterizations have not been considered yet. These parameterizations, such as the one used here, may be useful to identify relevant scenarios for which the determination of κ 3 is feasible. 10 10 In fact, the analysis carried here is a particular choice of two (linear combinations) of the 10 Wilson coefficients identified in Ref. [41]: κ t is related to δy t and κ V to c W = c Z . First of all we extended the framework and the notation introduced in Ref. [39] in order to take into account in the fit differential information, EW corrections in the production, and κ t and κ V dependence. The experimental inputs entering the fit are the signal strengths, which are defined for any particular combination i → H → f of production and decay channel as . (16) In Eq. (16), the quantities μ i and μ f are the production cross sections σ (i) (i = ggF, VBF, W H, Z H, tt H) and the BR( f ) ( f = γ γ, V V * , f f ) divided by their SM values, respectively. Assuming on-shell production, the product μ i × μ f is the measured rate for the i → H → f process divided by the corresponding SM prediction. This is valid also for differential distributions involving the reconstructed momentum of the Higgs boson. For simplicity, in the following we will refer with the symbol σ to both total cross sections or (bins in) differential distributions. The signal-strength productions μ i are given by where κ ggF = κ tt H = κ t and κ V H = κ VBF = κ V and the effect of κ 3 is parameterized via the quantity δμ i (κ 3 ). In presence of NLO EW corrections δμ i (κ 3 ) is given by where we have explicitly shown which quantities depend on the specific production process i. For differential distributions, differential K EW have to be used. As can be noted, we did not include κ t and κ V effects entering at one loop. As we will see in the results of the fit, we are going to probe deviations at the percent level in κ t and κ V . Thus, κ 2 t κ 3 and κ 2 V κ 3 effects are negligible for our purposes. On the other hand, terms of order κ 2 t κ 2 3 and κ 2 V κ 2 3 may be more important and in fact those from the Higgs wave-function can be consistently resummed; they are included in Eq. (17) via the Z BSM H κ 2 i term. 11 However, unless differently specified, we verified that also the inclusion of the κ 2 t κ 2 3 and κ 2 V κ 2 3 contributions has a negligible impact in the results presented in the following. It is also important to note that any further κ t or κ V dependence that may be introduced by NLO EW corrections, on top of those already present at LO, is negligible. Indeed, NLO EW corrections are per se at the percent level and their anomalous κ t component would be of the order of few percents of the corrections themselves. Thus, these kinds of effects, which similarly to those of order κ 2 t κ 2 3 and κ 2 V κ 2 3 can actually be calculated only in an EFT framework, are expected to be of the order ∼ α(κ i − 1) and therefore at the permille level or even smaller in our analysis. For this reason, we can safely ignore them. Similarly, the signal strength μ f for the Higgs decays H → f is given by NLO EW corrections in Higgs decays are small at inclusive level, therefore we can safely ignore them. The partial decay width in a given channel is given by where SM LO is the total width at LO in the SM. The SM widths SM ( f ) in Eq. (19) can be obtained by setting κ 3 = κ f = 1 in BSM ( f ). In order to ensure that the contribution of the Higgs-wave-function renormalization does not affect the branching ratios, in this case we resummed also the SM part (Z H = 1/(1−κ 2 3 δ Z H )) and factorize it to the κ 2 f dependence, as done in Ref. [39] for the LO analysis. For the γ γ decay channel κ γ γ depends on κ t and κ V , 12 κ V V * = κ V and κ f f = 1. Using Eq. (20) the signal strength for the decay becomes where in the last step we have assumed that C 1 is small, which is indeed true for the decay channels.

Results: comparison of differential and inclusive information in different scenarios
A first global fit on single Higgs channels has been performed in Ref. [39] using the 8 TeV LHC data, and a similar analysis has been applied to a future LHC scenario (CMS-HL-II) with 3000 fb −1 . Only total cross section information was used and especially the fit included only λ 3 as a variable. In Ref. [41] a first attempt to use differential rate information provided in Ref. [39] was made by extrapolating the projections on total cross section from ATLAS-HL [55,56,70] with 3000 fb −1 . 12 The relevant expression is where the functions A 1 (τ W ) and A 1/2 (τ t ) are defined in Ref. [69]. Since no differential information is available in the measured data at the moment, we focus on the same future scenario at 14 TeV (ATLAS-HL) considered in Ref. [41]. However, our results cannot be directly compared with those in Ref. [41], since there are a few differences in the treatment of the inputs from experimental projections. Details are reported in Appendix A, where we also carefully describe the procedure of the fit we performed and the assumptions on the uncertainties. In short, bounds on κ 3 , κ t and κ V are obtained by maximizing a log-likelihood function.
We perform the fit considering two very different scenarios for the uncertainties. In the first scenario (S1), only the statistical uncertainty is included. This crude assumption corresponds to the ideal (and rather unrealistic) situation where theoretical and experimental systematic uncertainties are negligible. On the other hand, we exploit it for a direct comparison with the second scenario (S2), where both theoretical and experimental systematic uncertainties are taken into account. At the differential level we performed the combination of the uncertainties via two different approaches that are described in detail in Appendix A. For this reason differential results for this second scenario always appear as bands rather than lines, accounting the uncertainty related to the different assumptions on the systematic and theoretical errors.
Before performing the global fit, we separately consider the different experimental inputs corresponding to ggF, VBF, V H and tt H production 13 and we restrict to the configura-tion with κ 3 only (κ t = κ V = 1). We remind the reader that different decay channels are entering for each of the production processes. Results are shown in Fig. 12, where the plot on the left refers to scenario S1 and the plot on the right to scenario S2. For the case of V H and tt H production dashed lines correspond to the fit of differential information; details of the binning are reported in Appendix A.
The different shapes of the curves for values smaller and larger than κ 3 = 1 can be understood from the behavior of κ 3 and κ 2 3 terms in Eqs. (6) and (13). While for κ 3 < 1 both the κ 3 and the κ 2 3 terms induce negative contributions in the production signal strengths, for κ 3 > 1 there are large cancellations that suppress the effect of κ 3 . If we only include the statistical uncertainty (S1) the ggF-like channel provides the best constraints for κ 3 both for the regions κ 3 > 1 and κ 3 < 1, where also tt H is giving strong constraints, which are not improved by the inclusion of differential information. A similar effect is visible also for V H; differential information does not lead to any significant improvement. On the other hand, in the region κ 3 > 1 we see a clear improvement due to differential information for tt H, although bounds from this single production process are not sufficient to set a constraint in the region for κ 3 > 1.
The plot on the right (S2) shows that including theoretical and experimental systematic uncertainties makes a difference. The tt H process is giving the strongest constraints in the region κ 3 < 1 and receive improvements from the differential information, with a tiny dependence on the assumption made for the combination of the uncertainties. This difference is induced by the change of the ggF result moving from scenario S1 to scenario S2 rather than by an improvement for tt H. Note, however, that the impact of the differential information for ggF production is not known and, while the exact calculation of the (two-)loop-induced effects from λ 3   Fig. 13 1σ and 2σ bounds on κ 3 including all production processes, based on future projections for ATLAS-HL at 14 TeV. Left: only statistical uncertainty (S1). Right: experimental systematic uncertainty and theoretical uncertainty included (S2) in pp → H j would be useful, it is currently out of reach. Although constraints from ggF becomes much weaker in scenario S2, in the region κ 3 > 1 they are still the strongest. At variance with ggF, tt H is in general very slightly affected by theoretical and systematic uncertainties since the dominant error is of statistical origin. Regarding the bounds on κ 3 from VBF-like and V H-like channels, they are always worse than those from ggF and tt H, even when the differential information is used for V H.
Next, we perform the global fit including all the experimental data as input and taking into account the anomalous couplings κ t and κ V . In Fig. 13 we present bounds after combining all the production channels, under different assumptions: i) only κ 3 is anomalous, ii) κ 3 , κ t or κ 3 , κ V are anomalous, iii) all three parameters κ 3 , κ V , κ t are anomalous. In the presence of anomalous couplings other than κ 3 , we marginalize over them. The plot on the left refers to scenario S1, only statistical uncertainties, and the one on the right to scenario S2, systematic and theoretical uncertainties included. As we expect, in scenario S1 the differential information (dashed line) does not noticeably improve any of the constraints, while in the scenario S2 in the region κ 3 < 1 and especially in the region κ 3 > 1 differential information from V H and tt H leads to a clear improvement of the constraints. What, instead, is not obvious, especially given the findings of Ref. [41], is the effect induced by anomalous κ t and/or κ V terms to the fit. While constraints in the region κ 3 < 1 are relaxed, although not washed out completely, by the inclusion of one or two new degrees of freedom, in the region κ 3 > 1 they are almost unaltered. In other words, in scenario S2, bounds in the region κ 3 > 1 are more affected by differential information than by the addition of the κ t or κ V parameters. Moreover, especially in the region κ 3 < 1, these two parameters alter the κ 3 constraints more in the unrealistic scenario S1 than S2. We describe rather in detail the observed features exploiting the information contained in Fig. 12.
In scenario S1 for κ 3 < 1 the constraints are strongly affected by the inclusion of κ t and/or κ V since the global fit with only κ 3 is completely dominated by ggF in that region. For this process only the total cross-section information is used in the fit, so that a flat direction appears, i.e., the fit is dominated by one input, 14 which is sufficient for setting constraints on only κ 3 but not at the same time on κ 3 and κ t , κ V . To resolve this degeneracy, more constraining information must be added to the fit. Indeed, the constraints with two parameters (κ 3 , κ t or κ 3 , κ V ) or three (κ 3 , κ t , κ V ) are in the region of the constraints from VBF and tt H in Fig. 12.
The previous argument cannot be applied to the region κ 3 > 1 for scenario S1, where the bounds in the global fit with only κ 3 are not completely dominated by ggF. Indeed the tt H (and in a smaller way the VBF) contribution is not negligible in that region, as can be seen from the left plot of Fig. 12. Moreover, at variance with ggF production, there is not a large background in tt H production for the experimental signatures involving the Higgs to μ + μ − decay, whose branching ratio has a different κ V and κ t dependence w.r.t. γ γ and V V * decays, and for values κ 3 ∼ 8 the impact of decays is more relevant. For this reason tt H and ggF are sufficient for constraining one, two or three parameters, with negligible difference when parameters other than κ 3 are marginalized. We explicitly verified this feature.
Moving to scenario S2, the plot on the right where all uncertainties are included, for κ 3 < 1 the bounds are dom- inated by tt H channel. For this reason there is a smaller dependence on the number of parameters considered in the fit and a larger sensitivity to the differential information, which is present for the same reason also in the region κ 3 > 1.
It is clear that the role of the ggF is essential when the impact of differential information is investigated in the global fit. When ggF is dominant, since there is no differential dependence, it masks the relevance of differential distributions. On the other hand, when tt H is dominant, the differential information becomes relevant. Above all, one should bear in mind that the impact of κ 3 on ggF distributions has not been calculated because of technical reasons; the exact two-loop calculation is beyond the current technology, but could be relevant too. To this purpose, in the following we look at constraints in the (κ 3 , κ t ) and (κ 3 , κ V ) plane with and without the contributions from VBF and ggF, which hides the impact of the differential information. We consider only scenario S2, which is more realistic.
In Fig. 14, we provide 1σ and 2σ contours in (κ 3 , κ t ) plane without (left) and with (right) anomalous κ V , which is anyway marginalized. Upper plots includes all the production channels, whereas in the lower ones only V H and tt H enter. Analogous plots are provided in Fig. 15 for the (κ 3 , κ V ), without and with anomalous κ t . First of all, one can note that due to the κ t dependence of the gluon-fusion channel and tt H channel, in the upper plots the constraints on κ 3 in presence of κ t (Fig. 14) are stronger than those in the presence of κ V (Fig. 15). 15 Also, in the upper plots, having two independent parameters (left) or marginalizing on an additional third one (right) does not lead to qualitatively significant differences. As also discussed before, the impact of the differential information is more important.
If we consider the lower plots the situation is very different. First, constraints with two or three parameters are qualitatively different. Second, the impact of the distributions is much more relevant. In the lower-left plot of Fig. 15 a flat direction is clearly resolved by differential information. The bottom-line is that by changing the number of free parameters and the number of inputs entering in the fit, the relevance of differential distributions and the sensitivity of the κ 3 -limits on additional parameters can be considerably altered. The range of the lower plots is much larger than in the upper plots; for this reason the exclusion of κ 2 3 κ 2 t and/or κ 2 3 κ 2 V terms from Eq. (17) would lead to visible effects to the 2σ contours, anyway without altering the qualitative information.

Conclusion
We have studied one-loop λ 3 effects for all the relevant single Higgs production modes at the LHC (ggF, VBF, V H, tt H, t H j) and decays (γ γ , V V * , f f , gg), extending and completing the results presented in Ref. [39]. In particular, we have calculated differential results for VBF, V H, tt H and t H j production and H → 4 decay. We have developed an automated code, which has been made public, for generating events including one-loop λ 3 effects. All the distributions that may be potentially affected by anomalous values of λ 3 have been scrutinized: differential level results for tt H production, H → 4 decay and also for the t H j process have been presented here for the first time.
We find that the production modes with a large kinematic dependence on λ 3 are V H, tt H and t H j. In particular, V H and tt H can provide additional sensitivity on λ 3 at differ-ential level. For these two channels we have consistently combined complete SM NLO EW corrections with anomalous λ 3 -induced effects at differential level. The same combination has been performed, at inclusive level, also for all the other production processes. We have verified the robustness of our strategy: NLO EW corrections are essential for a precise determination of anomalous λ 3 effects, but they do not jeopardize the efficiency of indirect λ 3 determination. We note that NLO EW corrections to t H j in the SM were unknown and have been calculated for the first time too.
Finally, we have performed a fit for κ 3 based on the future projections of ATLAS-HL for single-Higgs production and decay at 14 TeV [55,56]. We have considered the effects induced on the fit by additional degrees of freedom, in particular, anomalous Higgs couplings with the vector bosons and/or the top quark. We have found that, in a global fit, including all the possible production and decay channels, two additional degrees of freedom such as those considered here do not preclude the possibility of setting sensible λ 3 bounds, especially, they have a tiny impact on the upper bound for positive λ 3 values. On the contrary, the role of differential information may be relevant, critically depending on the assumptions on the future experimental and theoretical uncertainties. We have also shown that the relevance of differential distributions and the sensitivity on κ 3 can be considerably altered by varying the relation among the number of free parameters and the number of inputs entering in the fit.
Our results clearly illustrate the complementarity of precise single-Higgs measurements and double Higgs searches at the LHC for constraining λ 3 with the current and future accumulated luminosity. We therefore encourage experimental collaborations to use the MC tool provided here for performing λ 3 determination via single Higgs measurements, taking into account all the possible correlations among theoretical and experimental uncertainties of the different production and decay channels. (absolute) uncertainty for each channel obtained by summing in quadrature statistical (σ stat l, f ), theoretical (σ th l, f ) and experimental systematic (σ sys l, f ) uncertainties. The statistical and theoretical uncertainties are calculated as where th i is the relative theoretical uncertainty for each production mechanism i and we treat it as uncorrelated with the other different production mechanism. We list the theoretical uncertainty in Table 4, they are taken from the YR4 [62]. Concerning experimental systematic uncertainties, we list them directly in Table 3, based on an estimation from [55,56] and expressed directly as a number of events and not as relative numbers. In scenario S1, discussed in Sect. 5.2, we set σ th l, f = σ sys l, f = 0, while in scenario S2 we keep these uncertainties.
So far we discussed the case of the total cross section; numbers listed in Table 3 are for inclusive Higgs production. In the case of the differential distributions for V H and tt H, Eq. (24) can be generalized by independently considering each bin for these two processes. In practice, we split N l, f into several p T (H ) bins and for each bin j the number of events is given by where r i j (i = tt H, Z H, W H) is the ratio of the cross section of the bin j with the total cross section for process i. In other words, for each production-like mode we use r i j only form the dominant production process and decay, i.e., N X −like,i, f, j → r X j . The same assumption is made for the background. NLO EW K -factors at differential level are considered and used for computing the μ f i, j , which is simply the signalstrength prediction for each bin j. The chosen binning and the corresponding r i j , K EW and C 1 for each bin can be found in Table 5. For each bin, the statistical uncertainty is determined via its number of events and the relative theoretical uncertainty is assumed to be the same at the inclusive level. We may overestimate or underestimate the theoretical uncer-  tainty, since correlations are certainly present in the different bins but also in the different processes.
Concerning the experimental systematic uncertainty, we consider two cases. Either we scale it as r i j , so that in each bin the relative uncertainty that is present at the level of the total cross section is preserved, or as r i j , so that the sum in quadrature of the uncertainties for each bin is giving the value of the uncertainty for the total cross section, as in the case of the statistical uncertainty. The difference between the two approaches is small, as can be seen in Sect. 5; for the plots concerning scenario S2, the two approaches correspond to the border of the bands.