Higher Order QCD predictions for Associated Higgs production with anomalous couplings to gauge bosons

We present predictions for the associated production of a Higgs boson at NLO+PS accuracy, including the effect of anomalous interactions between the Higgs and gauge bosons. We present our results in different frameworks, one in which the interaction vertex between the Higgs boson and Standard Model $W$ and $Z$ bosons is parameterized in terms of general Lorentz structures, and one in which Electroweak symmetry breaking is manifestly linear and the resulting operators arise through a six-dimensional effective field theory framework. We present analytic calculations of the Standard Model and Beyond the Standard Model contributions, and discuss the phenomenological impact of the higher order pieces. Our results are implemented in the NLO Monte Carlo program MCFM, and interfaced to shower Monte Carlos through the {\sc Powheg} box framework.


Introduction
The LHC's discovery of a particle consistent with the predicted Standard Model (SM) Higgs boson has opened the door to a full understanding of electroweak symmetry breaking in nature. One of the key aims of Run II of the LHC is to study the properties and interactions of the Higgs in as much detail as possible, with the ultimate goal of confirming, or seriously constraining, the possibility of new physics and/or anomalous interactions.
One of the most interesting electroweak processes to study at the LHC is the interaction of the Higgs boson with massive vector bosons (W, Z). The primary role of the Higgs is to generate masses for these particles and ensure perturbative unitarity in vector boson scattering and any deviation from the SM Higgs-vector boson vertex could be indicative of new physics contributions. At the LHC the dominant Higgs production mechanism occurs through the fusion of gluons via a top quark loop. Therefore the total inclusive Higgs cross section at the LHC is more sensitive to the top Yukawa coupling than potential anomalous interactions of the Higgs with vector bosons. An obvious place to study the interaction between the Higgs and vector bosons are the decays H → V V * . However since the Higgs is considerably lighter than the 2m V threshold the decay phase space is restricted, forcing one of the final state vector bosons off-shell. Consequently, anomalous interactions that modify the high energy behaviour of the vertex, are suppressed due to the kinematic requirements. Accordingly, the best places to constrain anomalous interactions of the Higgs and vector bosons are those sensitive to said vertex in production, namely Vector Boson Fusion (VBF), Higgs in association with a hard jet, and associated production (V H). Of these, associated production -which occurs through an s-channel production mechanism -is particularly appealing, since one can directly probe the high energy behaviour of the interaction through, for instance, the invariant mass of the Higgs Vector system, m V H .
A simple way to encode effects of new physics in the Higgs sector is to study Higgs anomalous couplings (HAC) [1]. This parametrization does not rely on assumptions about whether EWSB is linearly or non-linearly realized, as it only relies on the Higgs as a scalar degree of freedom and the preservation of U (1) EM , i.e. by saturating all possible Lorentz structures in the vertex with the lowest number of derivatives. This parametrization was successfully used at LEP in the study of anomalous trilinear gauge couplings [2] and adopted in the study of BSM effects in the Higgs couplings.
An alternative way to describe indirect effects of new physics is to use an Effective Field Theory (EFT) approach. Within this approach, one could assume a linear realization of EWSB with the Higgs as a doublet of SU (2), and write down all the relevant operators which satisfy SU (2) L × U (1) Y [3]. This effective Lagrangian can be written in several equivalent ways which account for the choice of a basis. In this paper we will be using the proposal in Refs. [4,5]. A translation into other choices of basis can be done using, e.g. the tool Rosetta [6]. Also, one could write an EFT for a non-linear realization of EWSB as in Refs. [7]. In either case, there is a correspondence between the HAC and EFT approaches, see e.g. [5].
In this paper we focus on searching for BSM effects in Higgs production in association with a massive vector boson. The Higgs associated production process is defined through the following reaction, Where V represents an electroweak vector boson. In the SM V is constrained to be either a W or a Z, whilst including the higher dimensional operators allow for potential exchange of an off-shell virtual photon. The massive bosons are unstable and their decay products are measured at collider experiments. Leptonic decays of the vector boson are the cleanest experimentally, whilst the decay H → bb corresponds to the maximal Higgs branching ratio. Therefore, unless otherwise stated, the process we study in this paper corresponds to q(p 1 ) + q(p 2 ) → V (→ 1 (p 3 ) + 2 (p 4 )) + H(→ b(p 5 ) + b(p 6 )) (1.2) 1 and 2 correspond to either two charged leptons (V = Z) or a charged lepton and neutrino (V = W ). Associated production of a Higgs boson and a Z includes the following production process g(p 1 ) + g(p 2 ) → Z(→ − (p 3 ) + + (p 4 )) + H(→ b(p 5 ) + b(p 6 )) (1.3) This process is formally O(α 2 S ) and therefore a consistent treatment in a perturbative expansion would first include this piece at NNLO. However, the large gluon flux at the LHC, coupled with the boosted topology typical of experimental searches results in a significant contribution from this initial state. The SM contribution corresponds to two types of topologies, one in which the top loop radiates the Z, and one in which the Z is produced at the HZZ vertex. In the SM these two terms destructively interfere.
Given its phenomenological relevance significant theoretical effort has been invested in the V H processes. The NNLO corrections for the inclusive (i.e. on-shell V H) cross section were presented in ref. [8] based on the previous calculation of corrections to the Drell-Yan process [9]. The results of ref. [8] included the on-shell contributions from the gg pieces described in the proceeding paragraphs (a study of the gg pieces at NLO in the heavy top EFT was presented in ref [10]). A second type of heavy quark initiated contribution arises at NNLO and contains a qq pair. The leading contributions were calculated terms of an asymptotic expansion in m 2 t in ref. [11]. A fully differential calculation for the Drell-Yan type W H process was presented in [12], and was extended to include NLO H → bb decays in [13]. A similar calculation for ZH, including the gg diagrams was presented in ref [14] while a complete, fully differential calculation of the NNLO production (Drell-Yan and heavy quark) with decays for bb at NLO can be found in ref. [15]. The calculation of the NLO EW corrections for W H were presented in [16]. Fully differential predictions to the H → bb decay were computed at NNLO in [17] (the total inclusive H → bb rate is known to O(α 4 s ) [18]). In addition to the parton level calculations discussed above there has been significant progress matching parton level predictions to parton showers, allowing for full simulation of the LHC collisions. An NLO matched prediction for V H using the Powheg [19] framework was first presented in ref. [20] and using the MC@NLO [21] setup [22]. Studies using merged NLO samples of V H + 0 and 1 jet were in Powheg [23] and SHERPA [24]. A study for W H including anomalous couplings was presented at parton level in VBFNLO in ref . [25]. EW corrections have been implemented in the HAWK Monte Carlo code [26], including a study of anomalous HV V interactions at NLO in QCD.
The aim of this paper is to provide a framework to combine the precision predictions described above with the anomalous coupling prediction in a general HAC or EFT framework. We will do this be by calculating the V H processes at NLO including a general parameterization of the HV V vertex. We then interface our results to the POHWEG-BOX [19,20,27] allowing for full event simulation. This paper proceeds as follows, in section 2 we discuss the implementation of the anomalous couplings through the EFT Lagrangian. Section 3 presents our calculation in detail and provides the amplitudes for V H production at NLO including the anomalous couplings. Section 4 includes predictions at fixed order and NLO+PS accuracy and we present phenomenological results for LHC Run II. We draw our conclusions in section 5.

The Effective Standard Model
In this paper we focus on the effects of heavy New Physics in production of a Higgs boson in association with a vector boson. We are interested, hence, in three-point functions involving the Higgs and two vector bosons [1,5] as well as, possibly, couplings of the Higgs to a vector boson and two fermions. These Higgs anomalous couplings (HAC) are a model-independent parametrization which respects the fundamental symmetries of the SM at energies below electroweak symmetry breaking (EWSB), namely Lorentz and U (1) EM invariance, assuming that the Higgs is a neutral, scalar particle. The HAC can be related to an Effective Field Theory approach (EFT), where new resonances participating in EWSB are integrated out. The relation between HAC and EFT depends on assumptions of how EWSB occurs, e.g. whether the symmetry is linearly or non-linearly realized. In this paper we will match results in terms of HAC with a linearly realized EFT in which the Higgs h is part of a doublet of SU (2) L . We follow the conventions for EFT operators in [4,5], which are based on the work in Ref. [28]. The relevant part of the Lagrangian is as follows, where Φ is the Higgs doublet, and the dual field strength tensors are defined by Table 1: Translation between the HAC and EFT coefficients.
In Table 1 we show the relation between HAC and coefficients of the EFT. The basis we have chosen in this paper is not a unique choice and one can use, for instance, Rosetta [6] as a tool to translate among different basis.
The effect of HAC/EFT on the Higgs coupling to vector bosons is to introduce a nontrivial momentum dependence in the vertex, as one can see by inspecting the Feynman rules of the Higgs with W W , ZZ and Zγ vector bosons, which are presented in Figure 1. The tree-level SM piece is independent of the momentum, whereas New Physics potentially introduces new momentum-dependent Lorentz structures. New Physics can thus affect both rates of production and decay, as well as differential distributions. Exploiting differences in shape due to these new effects is one of the main avenues to look for indirect signals of New Physics. Here understanding SM higher-order corrections is especially important. This is primarily due to changes in shape arising from higher order matrix elements which are particularly relevant in the high-p T region.
We observe that in possible models which may generate these anomalous couplings, i.e. UV completions, not all Lorentz structures may be simultaneously generated. Indeed, in a large class of models, HAC of the type g (2) hvv do not occur, e.g. in 2HDMs [29], radion/dilaton exchange [29] or supersymmetric loops involving sfermions or gauginos [30]. This makes the study of phenomenology in which g (2) hvv = 0 and g (2) hvv = 0 particularly interesting. Finally, we comment on the current bounds for these operators from a global analysis including LEP and LHC Run 1 performed in Ref. [31], see Refs. [32] for other examples of fits in this context. This analysis took into account all the CP-even operators including pure gauge and operators involving fermions, but not the CP-odd couplingsg. In Table 2 Z ⌫ (q) haz p µ p ⌫ # Figure 1: Feynman rules for the Higgs-Z-Z, Higgs-W + -W and Higgs--Z vertices in the Higgs Anomalous Coupling description. All momenta are flowing into the vertex.
we show the result of a global fit. When g (2) hvv = 0, i.e.c W = c HW , the limits from a global fit of Run1 data are more stringent [29] c HW = c W = (0.0004, 0.02) . (2.5)

Calculation
In this section, we present the calculation of the NLO QCD predictions for a Higgs produced in association with a vector boson in the presence of anomalous couplings. We first present the calculation of the fixed order results needed for the NLO implementation in MCFM,

Calculation
In this section, we present the calculation of the NLO QCD predictions for a Higgs produced in association with a vector boson in the presence of anomalous couplings. We first present the calculation of the fixed order results needed for the NLO implementation in MCFM, and then proceed to discuss the implementation in the Powheg-Box. Our results are computed in terms of the modified Feynman rules presented in Figure 1 such that the anomalous couplings are a function of g hV V . At NLO accuracy the production and decay for pp → V H →leptons +bb completely factorize due to SU (3) color structure. This is because gluon radiation linking the initial state to the final state has no contribution at NLO since its interference with the LO amplitude results in a contribution proportional to Tr(T a ), where T a is an SU (3) generator. We therefore present amplitudes for pp → V H →leptons H, and allow the subsequent MC code to decay the Higgs boson (PYTHIA, or MCFM). In this paper the MCFM prediction corresponds to a LO decay, whilst the PYTHIA prediction includes effects from the parton shower.

Amplitudes for W H production
The LO amplitude for the associated production of a W and Higgs boson has the following form, where we have defined the full amplitude A (0) 5 in terms of a color stripped primitive amplitude A (0) 5 ; at LO the color factor is simply δ i 1 j 2 . In addition to the extraction of the overall color factor in eq. 3.1 we have also extracted the electroweak pre-factors from the fermionic W vertices. We note however that we have not extracted the g from the HW W vertex since this will be modified during our calculation. Although we will consider decays of the Higgs to bb, for simplicity we suppress the decay of the Higgs in this section.
We will also require the amplitudes needed to construct the NLO corrections. This consists of two new amplitudes, the one-loop virtual amplitude A (1) 5 and the real emission amplitude containing an additional parton A (0) 6 (in this case an additional gluon). The virtual primitive amplitude is defined as follows, The real emission amplitude including the parton 7 g 1 is thus The color stripping for the real emission amplitude is slightly more complicated than the LO and depends on the color matrix T a 7 , however there is still only one unique color ordering which simplifies the calculation significantly.
Since we are interested in associated production we are able to factorize the QCD corrections which affect the initial state, from the modifications to the HV V vertex. The factorization proceeds as follows for all of the primitive amplitudes we have considered, In the above equation J µ SM (p 1 , p 2 , V 1 , X) represents the production of a (chiral) current in the SM, if j = 5 then X = 0, whilst for the real emission amplitude j = 6 and X corresponds to the emission of an additional gluon. The second current V ν Λ (V 2 , H) corresponds to the splitting of the initial vector boson V 1 into V 2 and H, with the subsequent decays of V 2 to leptons included. Finally the two currents are connected by vector boson propagator, which in the unitary gauge is defined as follows, Since the W bosons decay we work in the complex mass scheme. However, in this calculation the longitudinal k µ k ν pieces do not contribute since J µ SM (p 1 , p 2 , V 1 , X)P 12X µ = 0 for massless initial and final states. We also frequently use the following function in our calculations In order to complete our discussion of the amplitudes for W H production we present results for the currents defined above. Since we are discussing W bosons, it is natural to relate the chiral currents to helicity amplitudes. Therefore we will use the language and notation of the spinor helicity formalism in the following definitions. We refer readers unfamiliar with the spinor-helicity formalism to one of the many reviews. For instance a detailed introduction can be found in ref. [33].
The first current we define is the modified decay current including the effect of the dimension-6 operators W µ where P 2 is the inflowing momenta, H 234 is the outgoing Higgs boson, and p 3 and p 4 represent the momenta of the final state leptons. The explicit form for this current is as follows, Next we consider the chiral currents in the SM. The current needed for the construction of the LO amplitude is rather simple, The virtual current corresponding to a vertex correction, is also very simple, Finally the current corresponding to the emission of an additional gluon, necessary in the real calculation has two possible helicity configurations, Contracting these various currents together results in the amplitudes for the production of W H at NLO in QCD including the effects of the dimension-6 operators. For example the LO contraction is,

Amplitudes for ZH production
Next we consider the production of a Z boson in association with a Higgs boson. The situation is slightly more complicated than the W H example considered in the previous section as the internal boson can be either a Z or a virtual photon γ * . In the SM the latter case does not occur, but the full general anomalous coupling parametrization allows for this contribution. We therefore parametrize the LO amplitude as follows, Unlike the case for W H there are now four allowed helicity configurations, corresponding to the selection of L and R helicities for particles 1 and 3. The left and right handed couplings (in the SM) are defined as follows, The sign in the v q L distinguishes between up (+) and down (-) type quarks. The amplitudes involving Z exchange can be obtained from the results presented in the previous section. The results for the analogous case in which W → Z can be obtained by simply swapping g → g cos θ 2 W and g hZZ in eq. 3.7. The results then correspond to the LL configuration, other configurations can be obtained from fermion line reversal symmetries. The current for a virtual photon exchange is given by, The SM currents are related to those described in the previous section.

Implementation in Monte Carlo codes
The amplitudes calculated above have been implemented into the parton level Monte Carlo code MCFM. Infrared divergences are regulated using the Catani Seymour Dipole subtraction method [34]. We use the default MCFM electroweak parameters, which correspond to the following choices, The remaining EW parameters are defined in terms of the above input parameters. Since we make the choice of defining the input of our program in terms of the dimension-6 Wilson coefficients of eq. (2.2), some additional effects are taken into account to fully map the physical effects of the EFT description into our HAC Lagrangian. Of the operators that we consider in our implementation -O W , O B and O γ -lead to non-canonical SU (2) × U (1) gauge boson kinetic terms after electroweak symmetry breaking. The field and gauge coupling redefinitions necessary to canonically normalise the theory lead to O(Λ −2 ) modifications of both the EW parameters in terms of the inputs as well as the couplings of gauge bosons to fermions as compared to the usual SM expectations. There is some freedom in how these redefinitions are performed and therefore the places in which these corrections appear although physical observables are naturally independent of such choices at this order in the EFT expansion. Appendix A describes the choices we make and therefore the origin of the relations and corrections that follow.
For the EW parameters, we work in the m W , m Z , G F scheme, and define the SM values for the Weinberg angle, electric charge and Higgs v.e.v as These are corrected due to a relative shift in the Z-boson mass which redefines the mixing angle and α(m Z ) as follows, While the definition of the Higgs v.e.v in terms of G F is unchanged. All other EW parameters are derived from the modified valuesĉ 2 W and e using SM relations. The coupling between a photon and fermion is corrected by the term, while the left and right handed couplings of the Z to a fermion, f , with weak isospin T f 3 and electric charge Q f are shifted as follows, Phenomenologically speaking, the effect of these re-definitions is minor, and typically results in changes which are sub-precent compared to predictions which do not alter the EW parameters by the EFT operators. For the choice of parameters simulated in this paper, only c W affects the EW parameters and neutral gauge boson couplings. The Higgs width is also modified as a result of the anomalous interactions, and we use the eHDECAY implementation described in ref. [4] to define the modifications to the Higgs width.
The Powheg-Box [19,27] provides a mechanism to match fixed order results to parton shower level Monte Carlo codes. In our case the implementation is rather straightforward, in particular since associated production in the SM has already been considered in the literature [20,35] . Therefore to include our results in the Powheg-Box we have simply updated the existing matrix elements with our own calculations described earlier in this section. The only technical task is to ensure that all of the variables required in the MCFM matrix element routines, are appropriately initialized by to the values assigned in Powheg. This is achieved through an interface to the MCFM routines which is called once at runtime.

Fixed order results
In this section we present results obtained at fixed order in perturbation theory. Specifically we study the dependence on the total rate at NLO as a function of the anomalous couplings. In order to simplify our results we do not vary all parameters continuously. Instead we focus on parameters which are representative of the phenomenology at the LHC. Ignoring for now the CP-violating operators we note that the in our basis the variables c W and c HW are sufficient to probe the Lorentz structures of the Feynman rule associated with g (1) hww and g (2) hww . In particular c HW = 0 and c W = 0 probes the regime in which g (2) hww modifies the vertex, and if c W = −c HW then the regime in which g (1) hww modifies the vertex is selected. In order to study the impact of the anomalous couplings we calculate cross sections for the LHC at the 13 TeV in which the final state particles have to satisfy the following phase space selection criteria, W H : We refer to this selection as our basic-cuts. Since the effects of the EFT operators are more apparent at higher energies we also define a high-p V T selection cut in which we impose an additional cut on p V T > 200 GeV. We note that p V T is a well defined experimental observable for both W and Z processes.
Our results for the basic and high-p V T cuts are presented in Figures 2 and 3 for W H and ZH processes respectively. The results have been obtained using MCFM, with the parameters described in the previous section. We use the CT10 PDF set [36] for NLO predictions and CTEQ6L1 for LO predictions. The renormalization and factorization scale have been set to µ = m V H . In each of Figures 2 and 3 the plot on the left side represents the cross section obtained using the basic cuts, whilst those on the right hand side correspond to the cross sections obtained with the additional high-p V T selection requirement applied. In both figures we plot the cross section as a function of c W and c HW . We present contours which correspond to the values of c W and c HW needed to obtain a 10, 20 or 30 % deviation from the SM prediction. For reference using our cuts described above the SM predictions are : 9.7 fb and 1.8 fb for W + H with the basic and high-p V T cuts respectively, and 5.1 and 0.54 for ZH. Our results contain terms up to order c 2 X as can be clearly seen from the figures, since the results for constant cross section form ellipses. Including these terms is somewhat dangerous, since in general they correspond to regions in which the EFT may be breaking down. This is because 8 dimensional operators also contribute first at this level, and therefore should not be ignored in the calculation. Therefore in an attempt to quantify the range of validity of our results we present a contour which corresponds to the regime in which the linear part of the cross section corresponds to 95% (solid) and 90% (dashed) of the total, or in other words the higher order pieces in the EFT should be small (both from 8 dimensional operators and loop corrections to the 6 dimensional operators which go like quadratic pieces). Experimental results can then be used to set reliable EFT bounds inside of this contour. We stress that values which lie outside of this contour (i.e. larger absolute values of c X ) cannot be reliably excluded given our theoretical accuracy, and given the form of our results, it is clear that there are regions outside of the EFT validity which conspire to produce small changes in the total cross section.
The hallmark of EFT operators is a lack of suppression at high energies due to poor high-energy behaviour. Therefore, a natural place to search for the impact of the higher dimensional operators is the region which is sensitive to larger values ofŝ. Since m V H cannot be directly cut upon in the experiment for leptonic W H final states, we use p V T as a proxy forŝ. The plots on the right hand side of Figures 2 and 3 present these results. As expected we see a significant increase in sensitivity in the (c HW , c W ) plane compared to the more inclusive analysis. By looking at high p V T one improves the sensitivity from around |c W + c HW | ∼ < 0.005 to around |c W + c HW | ∼ < 0.002. To demonstrate the flexibility of our code, we present results in the HAC basis, rather than the EFT approach, in Figure 4. In this setup the anomalous couplings are parameterized in terms of general Lorentz invariant operators in the Lagrangian. In this approach corrections from higher dimensional operators are not a concern so we present the full ellipses, and do not present EFT validity contours for these plots. We note that, since the HAC basis does not necessitate deviations in the EW parameters due to kinetic mixing of operators, we use the SM EW parameter choices (corresponding to the MCFM defaults) for these plots. In Figure 5 we investigate the impact of the NLO corrections to the anomalous couplings. We define the following ratio, Figure 4: Contours representing deviations as a function of the Higgs anomalous couplings g (1) hvv v and g (2) hvv v, where v is the Higgs vev, from the NLO SM prediction for W H (left) and ZH (right) production at the 13 TeV LHC.
Here R N LO (c W , c HW ) is defined as the full NLO result, divided by the NLO SM piece plus the LO anomalous coupling pieces. The results for W + H and ZH are presented in the figure. As might be expected from the inclusive K-factor the deviations are around ±20% depending on the position in the (c W , c HW ) plane. Around (−0.015, 0.01) the NLO corrections suppress the result one would obtain if a LO prediction were used, and previous limits in this region of phase space (using the LO prediction) may be too optimistic. On the other hand, away from this region the corrections tend to be positive and will improve existing limits. We note that the region which corresponds to that in which our EFT calculation is valid intersects the region in which the impact of the NLO corrections is most rapidly changing. This suggests that using a flat K factor to re-weight the anomalous coupling part of the calculation is not advisable.
Finally in Figure 6 we present contours of constant cross section in the (c HW , c HW ) plane, i.e. we study the impact of including CP odd operators. The CP odd operators do not interfere with the SM amplitude, so they first enter the cross section at O(c 2 HW ). This can be seen in the figures via thec HW ↔ −c HW symmetry in the figure. In order to have a relatively small deviation from the SM therefore requires a negative value of c HW which can compensate for the positive definite correction arising from the CP odd operator. In Figure 6 we present results for the basic cuts only, however improved results could be obtained by optimizing the analysis to look at high p V T observables. In addition angular distributions between final state particles are particularly sensitive to the CP structure of operators and represent a promising avenue of study.

NLO + Parton shower results
We now turn to the presentation of the fixed order plus parton shower results making use of the MCFM/Powheg-Box interface. The two processes considered are the production of the Higgs in association with a Z or W boson, where the Higgs decays to bb and the weak boson decays leptonically to e or µ. Events were generated for some characteristic values of the EFT Wilson coefficients and showered/hadronised with Pythia8. The decay of the Higgs was also performed by Pythia8 [37], with the total rate normalised to the NLO production cross section times the branching fraction as calculated by eHDECAY [38] Process Cuts 2 b-jets, p T > 25 GeV, η b < 2.5 1 lepton, ± (e or µ) 2 leptons, + , − (e or µ) p T < 25 GeV, |η | < 2.5 Table 3: Table summarizing the selection cuts performed on events in the two vector boson associated production modes.
to NLO accuracy in both α S and α EW . Event reconstruction and the implementation of basic selection cuts, summarised in Table 3, was performed using MadAnalysis5 [39] which makes use of Fastjet [40]. Benchmark EFT scenarios are selected to be within the allowed regions of recently performed global fits. In each case, an estimate of scale uncertainty is evaluated by varying the renormalisation and factorisation scales between half and twice the central scale, which is the invariant mass of the V H system. This is combined in quadrature with the usual Monte Carlo statistics uncertainty arising from the finite number of events generated. NLO and LO samples were generated with the CTEQ10 and CTEQ6L1 PDF sets respectively and PDF uncertainties were not estimated.

Gluon initiated contribution to HZ
In order to highlight the importance of the gg-initiated contribution to ZH production, a sample was generated separately and compared to the pure qq-initiated piece at NLO. In general, due to the 2m t thresholds, the kinematics of the box configuration will prefer a significantly harder region of phase space than the Drell-Yan like topology and, if it is not adequately taken into account, could show up as a fictitious EFT-like signal. Figure 7 overlays the two contributions in several differential distributions, showing the relative size of the would-be 'signal' that one may observe if the gg piece were not factored into the SM prediction. The contribution of the sub-process to the inclusive cross section is minor, of order 3-4%, but as it populates a high p T region, where the SM cross section is also quite small, the gg contribution can show up as an O(10 − 15%) effect in the tails of differential distributions, mimicking a potential EFT-like signature. The effect of this contribution on the N j distribution is even more striking, given the increased emission probabilities of the initial state, reaching around 40%. In Sections 4.2.2, this contribution is taken into account as part of the SM prediction for ZH associated production. The final panel in Figure 7 (and similar panels in future figures) encapsulates the emission of radiation in the form of the number of jets in the process. These jets can arise from either the matrix element or the subsequent parton shower. Since the matrix elements can provide at most one additional jet at NLO, the parton shower provides the additional jet multiplies beyond those of the NLO matching, which for these processes corresponds to N j > 3 (qq), and N j > 2 (gg) (and hence these rates are subject to larger theoretical uncertainties.)

LO +PS vs NLO +PS
To asses the impact of taking higher order effects into account, we now compare the new implementation at NLO to a LO one in MCFM, post showering and hadronisation. Figure 8 depicts a selection of differential distributions in the SM and for one of our benchmark points of c W = 0.004 (see discussion in Section 4.2.2) for HZ and HW production respectively.
Here the predominant effect is that of a relatively flat K-factor that is not sensitive to the presence of the new EFT interactions which, in any case, are colour neutral. The rightmost distributions in the upper half of Figure 8 show the p T of the ZH system and are therefore sensitive to the 'kick' that it receives from additional radiation. We see a mild rise in the tail between the NLO case that we would expect given that it captures the full phase space of one additional emission compared to the LO case which remains within the soft and collinear approximation of the parton shower.

EFT effects
We now turn to examining the effect of switching on one or more of the previously defined Wilson coefficients that affect associated production in both the ZH and W H channels. We limit ourselves to the c W and c HW coefficients as they are sufficient to capture the additional non SM-like momentum dependence brought about by dimension-6 operators. Figure 7 displays a number of characteristic differential distributions evaluated using the values of c W = −0.02 and c HW = 0.015, which saturate the bounds set by the most recent global fits [31]. In general, very large effects are expected for such sizable values of the coefficients and considering the discussion in Section 4.1, the validity of such an EFT description is called into question in the phase space regions where the BSM effects are important. Considering Figure 3, it is clear that the values of the coefficients lie well outside of the regions in which the quadratic piece of the EFT contribution makes up less than 10% of the overall contribution. The p b T distribution in Figure 9, for example, highlights very clearly the onset of a breakdown of the EFT in the high p T tail, where the relative contribution of the c W = −0.02 point changes sign, suggesting the dominance of the (c W /Λ) 2 term. Therefore, although these values of Wilson coefficients are technically 'allowed', the evidence in this section as well as in Section 4.1, suggests that we are not yet at a point in which the sensitivity of experiments can provide meaningful information about the coefficients affecting this Higgs production process.
Since we have not included effects from dimension 8 operators, which may be as large as the aforementioned squared EFT contributions, we prefer to present results using more conservative values of the coefficients, where the EFT interpretation is better motivated. These values are chosen from the criteria delineated in Section 4.1, i.e., the requirement that the squared terms do not make up more than 10% of the overall contribution. This leads us to choose |c W |, |c HW | = 0.004. Our two benchmark points derived from this arec W = 0.004 andc W = −c HW = −0.004. The relationship imposed in the latter choice is motivated by the results of previous works that calculated the low energy EFT coefficients predicted by a number of UV scenarios [29]. In the Two-Higgs Doublet Model, for example, this relationship is always satisfied at the matching scale. From a phenomenological perspective, this relation is also special because it corresponds to the elimination of one of the two momentum structures, g (2) hvv , present in the extended Higgs-gauge vertices (see Table 1 and discussion in Section 2).
For ZH production, Figure 10 shows differential distributions with respect to the Higgs-Z invariant mass, Z-boson p T and the number of jets (N j ) normalised to the 0-jet bin comparing the SM to the two EFT benchmarks. For the N j distribution, an additional cut on the Higgs p T of 200 GeV is applied in order to isolate the region where the EFT contributions are most important. In the case of W H production, the leptonic decay of the W + includes a neutrino which contributes to real missing energy on the event, preventing the construction of some of the kinematic variables available to the Z-boson associated production process, namely m V H and p V T , the total invariant mass and the vector boson transverse momentum. We trade these two observables for the total transverse mass of the system and the transverse momentum of the Higgs boson. Here, the total transverse mass of the HW system is defined including the two b-jets, the lepton and the missing transverse energy, The observable is the analogue of M V H in the ZH case and is an approximation of the momentum flowing through the W H vertex. These variables are shown in Figure 10 along with the normalised N j distributions, as in the ZH case. We see that these more conservative choices for the Wilson coefficients still permit O(20-50%) deviations in the tails of the various distributions for thec W = 0.004 benchmark with a clear preference for large momentum flow through the vertex. We also observe that the size of the deviation in the M V H distribution for ZH correspond roughly to the size of the deviation in p Z T at the corresponding energy scale, i.e, M V H ∼ p Z T /2, demonstrating the expected correlation between the two observables. The second benchmark ofc W = −c HW = −0.004 does not exhibit such large deviations, instead contributing a relatively flatter depletion of the differential rate. This can be traced to the different Lorentz structure governing the effective vertex. The difference between the two benchmarks shows that g (1) hvv leads to much more striking 'EFT-like' deviations than g (1) hvv . Looking more closely at the Feynman rules of Figure 1, we see that g (2) hvv goes as the square of the individual momenta of the Z bosons, while g (1) hvv goes as the product of the two Z boson momenta. As a consequence, in high centre of mass energy limit of the ZH production matrix element, g (2) hvv leads to a richer energy dependence, containing terms proportional to higher powers of Mandelstam variables ∝ st/M 2 Z , t 2 /M 2 Z that are not present when only considering g (1) hvv contributions. The N j distributions -although suffering from somewhat low MC statistics due to the p V T > 200 GeV requirement -appear to follow a similar trend.

Conclusions
Physics Beyond the Standard Model is likely to be connected to the Higgs sector, generically leading to deviations in the Higgs behaviour with respect to SM predictions. These indirect probes of new physics require a precise understanding of the SM contributions as well as the interplay between the SM and New Physics in observables. Among the different LHC Higgs observables, the production in association with a vector boson is specially sensitive to effects of new heavy particles in kinematic distributions and ratios of cross sections [31,41].
In this paper we have presented predictions for the associated production of a Higgs boson in association with a W or Z vector boson, including anomalous couplings between the Higgs and vector boson, not present in the Standard Model. Our predictions include effects in QCD beyond the Leading Order in perturbation theory. We presented predictions at fixed order (NLO) and matched to parton showers using the Powheg formalism (NLO+PS).
Anomalous couplings in the HV V vertex (HAC) can arise in many extensions of the SM. A general model independent parameterization can be obtained by saturating the Lorentz structures of the four-dimensional HV V interactions in the Lagrangian. Particular models then correspond to some (or all) of the new couplings acquiring non-zero values. An interesting class of models arise when the scale of new physics is large and can be integrated out of the Standard Model. In these scenarios the SM is treated as an effective field theory (EFT). We matched our results from the EFT to a linearly realized breaking of the EW symmetry, in which the Higgs is a doublet of SU (2) L . Transitioning between the two calculations setups is straightforward, and we presented results in both the HAC and EFT frameworks.
In order to maximize the physics potential of the LHC it is essential that precise theoretical predictions are used to compare theoretical predictions and experimental data. Matched parton showers, which combine the normalization and matrix elements of a Next-to-Leading Order calculation, and a leading logarithmic resummation of soft collinear logarithms provide a good framework for comparing theoretical predictions to data. The Powheg-Box provides a public format to match results obtained at fixed order to parton showers, allowing for full event simulation. We calculated the NLO corrections to V H production including the effects of anomalous couplings using analytic amplitudes and the spinor helicity formalism. We implemented this calculation into MCFM and modified the existing V H processes in Powheg to incorporate our new matrix elements.
We used our results to study the phenomenological impact of our calculation at the Run II of the LHC operating at √ s = 13 TeV. We demonstrated the capabilities of our code both at fixed order and NLO+PS accuracy generating events and showering them with the PYTHIA parton shower. We focused on parameter selections which are consistent with limits obtained during Run I. In this region NLO effects change the differential distributions by around O(20%). Our results will be made publicly available in the released versions of the MCFM and Powheg codes.    Q Z = ê c WŝW T 3 1 + δT Z 3 − Qŝ 2 W 1 + δQ Z , (A.10) are given in eqs. (3.24) and (3.25). We note that the difference between using the SM and EFT definitions of parameters multiplying a dimension-6 O(Λ −2 ) contribution is higher order in the EFT expansion. The extraction of the EW parameters from the {m Z , m W , G F } set of inputs discussed in sec. 3.3 follows from these definitions. It is important to stress that these definitions are valid for the subset of operators that are implemented in our code, namely O W , O HW , O B , O HB and O γ . In general, the presence of a complete dimension-6 basis of operators will lead to more modifications, such as with the O H and O T operators modifying the canonical normalisation of the Higgs field and therefore its couplings. A consequence of this can be seen in Table 1, where these two Wilson coefficients appear in the SM-like g (3) hzz structure.