Optimizing Higgs factories by modifying the recoil mass

It is difficult to measure the $WW$-fusion Higgs production process ($e^+e^- \to \nu \bar{\nu} h$) at a lepton collider with a center of mass energy of 240-250 GeV due to its small rate and the large background from the Higgsstrahlung process with an invisible $Z$ ($e^+e^- \to hZ,\,Z\to \nu\bar{\nu}$). We construct a modified recoil mass variable, $m^p_{\rm recoil}$, defined using only the 3-momentum of the reconstructed Higgs particle, and show that it can better separate the $WW$-fusion and Higgsstrahlung events than the original recoil mass variable $m_{\rm recoil}$. Consequently, the $m^p_{\rm recoil}$ variable can be used to improve the overall precisions of the extracted Higgs couplings, in both the conventional framework and the effective-field-theory framework. We also explore the application of the $m^p_{\rm recoil}$ variable in the inclusive cross section measurements of the Higgsstrahlung process, while a quantitive analysis is left for future studies.


Introduction
A lepton collider running at a center of mass energy of around 240 to 250 GeV is ideal for studying the properties of the Higgs boson. It can collect a large amount of Higgsstrahlung (e + e − → hZ) events, which has a cross section maximized at around 250 GeV. At higher energies, the W W -fusion process of Higgs production (e + e − → ννh) can be better measured as its cross section increases with energy. It is important to have good measurements of the W W -fusion process, which provides information complementary to the one from the Higgsstrahlung process. In the conventional kappa framework, the W Wfusion process can constrain the hW W coupling and is also an important input for the determination of the Higgs total width. 1 In the effective-field-theory (EFT) framework, the W W -fusion and the Higgsstrahlung processes probe different combinations of EFT parameters. The inclusion of both processes, as well as the diboson one (e + e − → W W ), is crucial for discriminating different EFT parameters and obtaining robust constraints on all of them [2][3][4]. However, it is not guaranteed that the runs at energies higher than 240-250 GeV will be available. The proposed Circular Electron Positron Collider (CEPC) in China does not have plans for the 350 GeV run at the current moment [5]. For the Future Circular Collider (FCC)-ee at CERN [6] and the International Linear Collider (ILC) in Japan [7], a significant amount of time may also be spent at the 240 GeV/250 GeV run before moving on to higher energies. The measurements of the W W -fusion process at 240-250 GeV is therefore of great relevance to the study of Higgs physics. It is difficult to measure the W W -fusion process at 240-250 GeV for the following two reasons. First, it has a small rate at lower energy, with a cross section of 6.72 fb at 250 GeV assuming unpolarized beams, while the total cross section of the Higgsstrahlung process is 212 fb at the same energy [5]. Second, the Higgsstrahlung process with Z decaying invisibly (e + e − → hZ, Z → νν) is the dominate background of W W fusion. Both contributing to the channel e + e − → ννh as shown in Fig. 1, the cross section of the former is more than six times of the latter one. With longitudinal beam polarizations, the situation is slightly better. If the signs of the polarizations are sgn(P (e − ), P (e + )) = (−, +), the W W -fusion cross section is enhanced, and by a larger factor than the one of e + e − → hZ. The method of recoil mass can also be used to separate the W W -fusion and hZ events, as the reconstructed mass of the neutrino pair should center around the Z mass for hZ events. However, the discriminating power is limited by the detector resolutions (especially for hadronic Higgs decays) and other effects [5,8,9]. A consistent treatment of the hZZ and hW W couplings is also required since they are related by gauge invariance. In the EFT framework, the relation is complicated by the inclusion of Dimension-6 operators which generate anomalous couplings with Lorentz structures different from the standard model (SM) ones.
In this paper, we try to address the issues mentioned above and further optimize the measurements in e + e − → ννh at 240-250 GeV. We first perform a collider study in Section 2 with a comparison of the recoil mass variable and its variations. We try to validate our study by following closely Ref. [9]. We point out that the variable m p recoil , defined using only the 3-momentum of the reconstructed Higgs particle, could provide a discriminating power better than the original recoil mass variable m recoil does. We then implement the m recoil and m p recoil distributions in the EFT global analysis in Section 3, using the framework in Ref. [3]. We point out the importance of fitting the EFT parameters directly to the binned m recoil or m p recoil distribution instead of fitting them to the extracted cross sections of the W W -fusion process. We also apply m p recoil to the inclusive e + e − → hZ process in Section 4 and comment on its potential use in the inclusive hZ cross section measurements. Finally, we conclude in Section 5. We provide a short summary of the EFT framework used in our analysis in Appendix A and the numerical expressions of the EFT dependence of the (modified) recoil mass distributions in Appendix B.
The following collider scenarios are considered in our study: • CEPC with 5 ab −1 data collected at 240 GeV, with unpolarized beams [5]. This scenario can also be thought as the earlier stage of the FCC-ee, which also plans to collect 5 ab −1 at 240 GeV and eventually 1.5 ab −1 data at 350 GeV as well [6].
• ILC with 2 ab −1 data collected at 250 GeV and beam polarizations of P (e − , e + ) = (±0.8, ±0.3), which could be considered as the first stage of a full program with center of mass energies up to 500 GeV [10].
For the W W -fusion measurements, we focus on the channel with the Higgs decaying to a pair of bottom quarks (e + e − → ννh, h → bb), which has the largest branching ratio. The measurements of W W fusion at 240-250 GeV with other Higgs decay channels are not reported in the official documents due to the poor constraints (see e.g., Refs. [5,10,11]).

The modified recoil mass of e + e − → ννh, h → bb
At lepton colliders, the method of recoil mass can be used to reconstruct the mass of a particle without measuring its decay products. One of its most important applications is the measurement of the inclusive rate of the Higgsstrahlung process, e + e − → hZ. Assuming both the Higgs and Z are on mass shell, one could write the relation where the total center of mass energy √ s is fixed up to corrections from beam energy spread and initial state radiations. By measuring the energies and momenta of the Z decaying products one could reconstruct the mass of the Higgs particle. This can be used to select e + e − → hZ signal events without tagging the Higgs decay products, which makes it possible to measure the inclusive cross section of this channel. It also provides the best Higgs mass measurement. For example, a precision of 5.9 MeV can be achieved with the leptonic Z decay channels of the inclusive hZ measurements at the CEPC [5].
If the Higgs decay products are measured, the recoil mass can be turned around to reconstruct the Z mass, since the following relation also holds for an e + e − → hZ event, The recoil mass can then be defined as where E rec h and m rec h are the reconstructed Higgs energy and mass. For Higgs decaying to a pair of bottom quarks, they are the total energy and the invariant mass of the two b-jets. This offers a way to separate the Higgsstrahlung events with an invisible Z (e + e − → hZ, Z → νν) from the W W -fusion events, both contributing to the channel e + e − → ννh, as shown in Fig. 1. However, due to finite jet resolutions, beam energy spread and other effects, the recoil mass distribution of the hZ events has a rather large spread. This limits it discriminating power especially at the energy 240-250 GeV, for which the recoil mass distribution of the W W -fusion events spreads around the same region. We also find that for Higgs decaying to a pair of b-jets (h → bb) the uncertainty on the recoil mass is dominated by the energy and momentum resolutions of the b-jets. This is also obvious from the observation that the recoil mass distribution for the Higgs mass reconstruction in Eq. (2.1) is much narrower for the leptonic Z decay channel than for the hadronic one (see e.g. Ref. [5]).
While the recoil mass defined in Eq. (2.3) makes use of all the kinematic information, it does not make any assumption on the value of the Higgs mass. Both the Higgs width and the projected uncertainty of its mass are at the MeV level and can be neglected compared with the effects of jet resolution. Using the information of the Higgs mass, two modifications of the recoil mass can be constructed. The first, using only the reconstructed Higgs energy as the measurement input, is defined as where m h is fixed to be the Higgs mass, 125.09 GeV. The other, using only the reconstructed 3-momentum ( p rec h ) of the Higgs, is defined as where m h is again fixed to be 125.09 GeV. At the truth level, m recoil , m E recoil and m p recoil are all equivalent. However, the uncertainties in the energy and momentum measurements certainty have different impacts on the three variables. To illustrate this impact, we define, for a given event, a set of five parameters {δ m , δ E m , δ p m , δ E , δ p } which parameterize the differences between the reconstructed quantities and the true ones, with and where m true recoil is the true parton level recoil mass, E h and p h are the true energy and 3momentum of the Higgs. For hZ events, m true recoil = m Z (assuming it is on shell), and the three parameters δ m , δ E m and δ p m can be written in terms of δ E and δ p . At leading order, they are given (for hZ events) by  Note that δ E and δ p can be either positive or negative. The overall negative coefficients in Eq. (2.8) indicates that if the measured energy or 3-momentum of the Higgs is larger than its actual value, the recoil mass variables will be smaller than the Z mass, and vice versa.
where the small coefficients of δ p come from a suppression factor of ∼ | p h | 2 /E 2 h relative to the ones of δ E , shown in Eq. (2.8). The distributions of δ E and δ p for e + e − → hZ, Z → νν, h → bb at CEPC 240 GeV are shown on the left panel of Fig. 2, after applying the selection cuts which include a Higgs-mass-window cut of 105 GeV < m rec h < 135 GeV on the b-jet pair. While δ p has a slightly larger spread than δ E , its coefficients in Eq. (2.9) are much smaller. We therefore expect the distribution m p recoil to have the smallest spread, and the one of m E recoil to have the largest. This is verified on the right panel of Fig. 2 where the distributions of m recoil , m E recoil and m p recoil are shown. For the W W -fusion events, we expect a less significant difference among the distributions of the three variables (which are shown later in Fig. 3 & 4), since they do not have a Z in the event. The corresponding distributions for ILC 250 GeV are very similar to the ones in Fig. 2.
It should be noted that the distribution of δ E in Fig. 2 is asymmetric, suggesting that on average the measured energy of the b-jet pair is smaller than its actual value. This is due to the fact that in our simulation we do not apply any jet energy corrections that are widely used in the LHC experiments [12,13]. As a result, the central values of the m recoil and m E recoil distributions are also shifted to be larger than m Z . Assuming a jet energy correction mechanism will be implemented in the future lepton collider(s), one would expect a corrected central value and also some improvements on the energy measurement, and the m recoil (and m E recoil ) distribution will have a peak value around m Z and a slightly smaller spread. We do not expect the lack of jet energy correction to have a significant impact on our results since the m p recoil distribution still has a much smaller spread due to the parametric suppression of it uncertainty near the hZ threshold as discussed above. 2 Having found that the variable m p recoil could better reconstruct the Z mass than m recoil , we perform an analysis based on a fast simulation to explicitly exam their discriminating powers on the W W -fusion and hZ events at the CEPC 240 GeV (with unpolarized beams) and ILC 250 GeV (assuming P (e − , e + ) = (−0.8, +0.3)). We generate events for both processes using Madgraph5 [14], which are showered with Pythia [15] before passing to Delphes [16] with ILD cards (using the detector geometry and flavor tagging efficiencies given in Ref. [17]) for detector simulations. The interference term between the W W -fusion and hZ processes are ignored. It should be noted that the effects of ISR photons are not considered in the simulation with Madgraph5. However, we expect their effects to be much smaller compared with the ones of jet resolutions. We use the ILC analysis in Ref. [9] as a guide to validate our results from the simple simulation. While the Durham jet clustering algorithm is used in Ref. [9], it pointed out that the anti-k t jet algorithm with jet radius R = 1.5 has a similar performance in the Higgs invariant mass reconstruction, which is used in our simulation. We also follow closely the selection cuts in Ref. [9]. In particular, each event is required to have exactly two b-jets and a cut on the invariant mass of the b-jet pair, 105 GeV < m rec h < 135 GeV, is applied to reduce the backgrounds. The cuts related to variables in the Durham jet clustering algorithm are replaced by the simple requirement on jet number (= 2) in each event. After event selections, we scale the number of signal events of ILC 250 GeV to the ones in Ref. [9] (normalized to 2 ab −1 ). A similar scaling is also applied for CEPC, taking count of the differences in cross sections and selection efficiencies between CEPC and ILC.
The composition of the background in the e + e − → ννh channel is also listed in Ref. [9]. The major ones are ννbb and qq, which contributes to 42% and 34% of the total background after selection cuts. The qq background is difficult to simulate due to its huge cross section and tiny selection efficiency. For simplicity, we simulate only the ννbb background, apply the selection cuts and scale it up to match the total background number, given in Ref. [9] and normalized to our run scenarios. We expect this simple treatment to provide a reasonable estimation on the effects of the backgrounds.
After selection cuts, the m recoil and m p recoil distributions of hZ (Z → νν), W W -fusion and background events are shown in Fig. 3 for ILC 250 GeV (2 ab −1 data with P (e − , e + ) = (−0.8, +0.3)) and Fig. 4 for CEPC 240 GeV (5 ab −1 data, unpolarized beams). In Fig. 3 and top panel. To check the validity of our results, we compared our m recoil distributions to the ones in Refs. [9,5] and found a reasonable agreement in terms of spread ranges and distribution shapes. Comparing the m recoil and m p recoil distributions in Fig. 3 and 4, it is clear that m p recoil provides a better discrimination between hZ and W W -fusion events, with the m p recoil distribution of hZ having a much sharper peak around the Z mass than the m recoil one. On the other hand, the m p recoil distribution of the ννbb background has a more even spread than the m recoil one. This is because that for the ννbb background the bb pair does not come from the Higgs decay. Using the wrong mass assumption therefore makes the reconstruction of Z mass worse.
Following Ref. [9], we apply a χ 2 fit to the binned m recoil and m p recoil distributions to extract the precisions (one-sigma uncertainties) of the hZ and W W cross sections, denoted as σ hZ and σ W W →h . This is done by treating the overall rates, σ hZ and σ W W →h as free parameters while assuming a perfect knowledge on the the shapes of the distributions. Ref. [9] also treats the overall cross section of the background (σ bg ) as a free parameter. We consider two cases with σ bg treated as a free parameter as well as fixing it to the predicted value. The total χ 2 of the m recoil or m p recoil distribution is given by (2.10) where for each bin i, n i exp is the expected number of events from simulation and n i theory is a function of σ hZ and σ W W →h (and σ bg ). To ensure enough statistics in each bin, we choose a bin size of 5 GeV except for the first and last bin, which are chosen to include all the events below 65 GeV (60 GeV) and above 130 GeV (115 GeV) for the m recoil (m p recoil ) distribution at ILC, and all the events below 65 GeV (60 GeV) and above 120 GeV (110 GeV) for the m recoil (m p recoil ) distribution at CEPC. The results of the χ 2 fits are presented in Table 1 for ILC 250 GeV (2 ab −1 data with P (e − , e + ) = (−0.8, +0.3)) and Table 2 for CEPC 240 GeV (5 ab −1 data, unpolarized beams).
Overall, our results on the precision of the W W -fusion cross section is slightly worse than the ones in Ref. [9] (if normalized to the same luminosity) and the CEPC preCDR [5]. This is not surprising, since the results could depend on details of the simulation, for which we only performed a simplified study. In what follows, we shall focus on the relative difference between the results from m recoil and m p recoil . While the distributions of m p recoil clearly better separates the hZ and W W -fusion events as shown in Fig. 3 and Fig. 4, in the 3-parameter fit (shown on the left panels of Table 1 and   of σ W W →h from the m p recoil distribution is similar to (or even worse than) the one of m recoil , due to the fact that W W -fusion and background events have a larger overlap in the m p recoil distribution than in the m recoil one. m p recoil nevertheless significantly improves the precisions of σ hZ and σ bg in the 3-parameter fit. Assuming a good knowledge of the background, one may also fix the background cross section to the predicted value. In the 2-parameter fit with σ hZ and σ W W →h , the m p recoil distribution indeed provide a significantly better constraint on σ W W →h , with an improvement of about 30% at ILC and 20% at CEPC, compared with the constraint from the m recoil distribution.
In the kappa framework with the hZZ and hW W couplings treated as independent parameters, the W W -fusion measurement is an important input for constraining the hW W coupling (with the other being the Higgs decay, h → W W * ). With a 20-30% improvement on the precision of the W W -fusion cross section, a sizable improvement for the constraint on the hW W coupling is expected. The W W -fusion cross section is also an important input for the determination of the Higgs total width, following the relation [  to it. With a 20-30% improvement on the precision of σ(ννh, h → bb), as well as some possible improvement on the determination of BR(h → bb) from a better measurement of σ(hZ, Z → νν, h → bb), the precision of the Higgs total width obtained using Eq. (2.11) could be improved by at least 20-30% using the m p recoil variable. To conclude this section, we would like to emphasize that, while we try our best to validate our results, they do rely on simple simulations and should be explicitly tested by experimental groups with proper simulation tools. The fits performed in obtaining the results in Table 1 and 2 also assume a perfect knowledge of the distribution shapes for each process, which may not be a good assumption in an actual experiment. We also include only the ννbb background, while other backgrounds may have different kinematic features. Nevertheless, we expect our results to still hold qualitatively due to the simple reasoning that m p recoil has a smaller uncertainty and better reconstruct the Z mass.

Improving Higgs coupling constraints in the EFT framework
Having explored the capability of the m p recoil variable in improving the measurement of the W W -fusion process, we are now ready to exam its impact on the determination of the Higgs couplings. We choose to study the Higgs coupling constraints in a global effectivefield-theory (EFT) framework with dimension six (D6) operators. 4 Such a framework has several advantages. First, assuming the scale of new physics is high, the EFT with D6 operators gives a good parameterization of the effects of new physics and the results can be mapped to any specific model that satisfies the assumptions of the framework. Second, it takes count of the connections among different measurements. For instance, some operators contribute to both Higgs processes and the dibson process, and the triple gauge coupling (TGC) measurements from the diboson process can thus help the overall constraints on the Higgs couplings [24]. Gauge invariance is also imposed by construction in the EFT framework. We focus on the CEPC 240 GeV and follow Ref. [3] in terms of the basis choice and measurement inputs. In particular, focusing on the Higgs and diboson measurements at 240 GeV, and making reasonable assumptions, a total number of 11 parameters are sufficient to describe the contributions from the D6 operators. A short summary of the framework in Ref. [3] is provided in Appendix A. The methods we propose should nevertheless be applicable to other collider scenarios and frameworks.
A few important differences between the cross-section fit in Section 2 and the EFT analysis should be noted. While not specifically mentioned, the cross-section fit does make assumptions on the new physics, in particular that it only modifies the overall rates, not the differential distributions, of the hZ and W W -fusion processes. In the cross-section fit, the hW W and hZZ couplings are also assumed to be independent, regardless of their relation from gauge invariance. The EFT analysis, while imposing gauge invariance, contain anomalous couplings of the form hZ µν Z µν and hZ µ ∂ ν Z µν (and the same for W ) which have different momentum dependences from the SM couplings due to the extra derivatives. The potential new physics contribution to the hZγ vertex is also included in the EFT analysis, which could contribute to e + e − → hZ via an s-channel photon. It is thus beter to directly fit the EFT parameters to the m recoil or m p recoil distribution of the inclusive e + e − → ννh process instead of fitting them to the extracted precisions of cross sections from Section 2. By fitting to the inclusive e + e − → ννh process we also include the interference term of W W -fusion and hZ processes, which is usually ignored in crosssection fits. The expressions for the total cross section of e + e − → ννh and the (binned) differential ones of the m recoil and m p recoil distributions in terms of the EFT parameters are listed in Appendix B.
The measurement inputs of the Higgsstrahlung (e + e − → hZ) and diboson (e + e − → W W ) processes at CEPC 240 GeV are listed in Table 3. The estimations of hZ measurements are taken from Ref. [25], which updates the ones in the CEPC preCDR [5]. In addition, The angular observables of e + e − → hZ in Ref. [19] are included, for which we use only the channel e + e − → hZ, h → bb, Z → + − and assume a fixed 60% signal selection efficiency, following Refs. [21,3]. For the TGC measurements, we follow the treatment in Ref. [3] which adopts the one in Ref. [26] with the addition of a universal 1% systematical uncertainty in each bin of all differential distributions. We directly list the resultant one-sigma constraints of the anomalous TGC parameters and their correlations in Table 3. We construct the total χ 2 by summing over the χ 2 s of all measurements and perform global fits to obtain the precision reaches (one-sigma bounds) of the relevant EFT parameters.
We consider three scenarios in the global analysis. All three use the inputs on Higgsstrahlung and TGC measurements in Table 3, but different information on the measurement of e + e − → ννh. The first one uses only the total rate of e + e − → ννh. The second (third) uses the information in the m recoil (m p recoil ) distribution, with the EFT parameters CEPC 240 GeV, 5 ab −1 , unpolarized beams e + e − → hZ e + e − → W W σ(e + e − → hZ) 0.50% uncertainty correlation matrix  [3].
directly fitted to the binned distributions. We also compare the reach with the more conventional method of fitting the EFT parameters to the extracted precisions of hZ and W W -fusion cross sections in Section 2, ignoring the correlation between the two cross sections (which is often not reported). Since the background is assumed to be SM-like in the EFT analysis, 5 for the extracted precision of cross sections we also use the results of the two-parameter fit with fixed background, shown on the right panel of Table 2. The results of the 11-parameter fit are presented in Fig. 5.
Comparing the reaches of the three scenarios (shown by the orange, blue and cyan columns in Fig. 5), we indeed observe a non-negligible improvement on the overall reach if the information in the m recoil or m p recoil distribution is used. In particular, the reach on the parameter δc Z (corresponding to a shift in the SM hZZ and hW W couplings) is improved by more than 10%. Consequently, the reach onc eff gg , δy c , δy b and δy τ , which contribute to the Higgs decay to gg, cc, bb and ττ , have also been improved by a similar (or slightly less) factor. The reach with the m p recoil distribution is better than the one with The second (third) one uses the m recoil (m p recoil ) distribution, and contain two sub-scenarios. The one shown by the light shaded columns uses the extracted precisions of σ hZ and σ W W →h in Section 2 (the two-parameter fit on the right panel of Table 2, correlation ignored). The one shown by the solid columns is obtained from a direct fit to the binned m recoil (m p recoil ) distribution.
m recoil as we expected. However, the relative improvement from m recoil to m p recoil turns out to be very marginal. We also find that, if fitting the EFT parameters to the extracted cross sections σ hZ and σ W W →h without taking count of their correlation (the results are shown with light shades for the 2nd and 3rd columns), the reaches are worse than the ones from direct fitting the EFT parameters to the distributions, in particular for the m recoil distribution. This is because, as shown in Section 2, the uncertainties of σ hZ and σ W W →h have a large correlation between them due to the difficulty in separating the two, in particular for the m recoil distrubiton. This correlation is usually not reported in official documents, and the omission of it could lead to a considerable impact on the overall reach. It should be noted that, our results using only the total rate of e + e − → ννh is worse than the corresponding ones in Ref. [3]. This is because our estimations on the rate measurement of e + e − → ννh is more conservative than the one in Ref. [3], which is derived from the CEPC preCDR [5]. If the overall cross section measurement of e + e − → ννh can be improved (e.g. by optimizing the selection cuts), we expect the use of m recoil and m p recoil distributions would also bring a more significant improvement on the overall reach of the EFT fit. We have also chosen very conservative bin sizes to control the uncertainties in each bin from simulation. Further optimizations of the analysis may also provide substantial improvements on the reach with the m recoil and m p recoil distributions.
We also find that, if the TGCs can be measured with much better precisions such as the ones in Ref. [4] for ILC 250 GeV (which are one order of magnitude better than the ones in Table 3), or if multiple runs with different beam polarizations are available (also likely to be the case for ILC), the improvement from the m recoil and m p recoil distributions with respect to using only the total rate of e + e − → ννh becomes rather insignificant. This is not surprising, since very precise TGC measurements can effectively remove two degrees of freedom in the fit so that there is less need for additional handles to discriminate the parameters. The interference term of the e + e − → hZ diagram with an s-channel Z and the one with an s-channel photon is also sensitive to the beam polarization, which can help probe the operators that contribute to this interference [3,4]. In general, once a sufficient number of constraints are included in a global analysis, the overall precision reach is expected to be less sensitive to the impact of a single measurement, such as the one of e + e − → ννh. It is nevertheless important to optimize the measurements in order to maximize the sensitivity to new physics.

Applications on the inclusive hZ measurements
While we have focused on the W W -fusion measurements in the previous two sections, a question of great interest is whether the variable m p recoil could be applied to improve the inclusive measurement of e + e − → hZ, where the decay product of Z are tagged instead. While suffering from the jet resolution, the hadronic Z decay channel provides a slightly better measurement on σ(Zh) than the leptonic one, thanks to its large branching ratio. 6 An improved measurement of the hadronic Z channel could thus have a significant impact on the overall precision reach of σ(hZ). It is straightforward to write down the recoil mass and its two variations for the reconstruction of the Higgs mass in e + e − → hZ, which are where E rec Z , p rec Z and m rec Z are the reconstructed energy, 3-momentum and invariant mass of the Z from the two jets, while m Z is the true Z mass, fixed to be 91. 19 GeV. Similar to Eq. (2.9), we derive the deviations in the measured m recoil , m E recoil , m p recoil as a function of the deviations in the measured Z energy and 3-momentum to be where δ m , δ E m and δ p m are defined as with m true recoil = m h . For δ E and δ p , the definitions are where E Z and p Z are the true energy and 3-momentum of the Z. Similar to Eq. (2.9), in Eq. (4.2) the coefficients of δ p is also smaller than the ones of δ E , but with a suppression factor of ∼ | p Z | 2 /E 2 Z instead. The distributions of δ E and δ p for the reconstructed Z are shown in Fig. 6, with the details of simulation stated later in this section. Note that the cut on the Z-mass window has a strong impact on the distributions of δ E and δ p . For larger deviations of the measured energy and momentum from the true ones, the invariant mass also tends to be further away from its true value.
To compare the reconstruction power of m recoil and m p recoil on the Higgs mass, we perform a simple analysis using the simulation tools listed in Section 2. One important difference here is that for Higgs inclusive measurement with Z → qq, the final states could contain additional jets from Higgs decay, making it more difficult to reconstruct the Z. Due to the additional jets, we set the jet radius to R = 0.5 in order to reduce the contamination among the jets. For an event with more than two jets, we choose the pair of jets with an invariant mass that is closest to the value of Z mass. We then apply a Z-mass-window cut on the invariant mass of the jet pair, m qq , intended for removing backgrounds. The difference between m recoil and m p recoil is strongly correlated with the size of the Z window -in the limit that the invariant mass equals the actual Z mass, m recoil and m p recoil become equivalent. We therefore consider both a larger window, 70 GeV < m qq < 110 GeV, and a smaller one, 80 GeV < m qq < 95 GeV. The distributions of m recoil and m p recoil for e + e − → hZ, Z → qq after the selection cuts are shown in Fig. 7 for CEPC 240 GeV. To estimate the impact of the combinatorial problem in the reconstruction of Z, we first consider a case in which the Higgs is forced to decay invisibly in the   simulation. The only purpose of the invisible decay is to avoid having additional jets from the Higgs decay and ensure a clear identification of the Z jet-pair. The results are shown on the top panels of Fig. 7 for the two choices of Z-mass-window cuts. For this ideal case, it is clearly that m p recoil has a significantly narrower spread and provides a much better reconstruction of the Higgs mass than m recoil does. The improvement with m p recoil is more significant if a large Z-mass window cut is applied as we expected. For the realistic case with Higgs inclusive decays, the distributions are shown in the bottom panels of Fig. 7. The reconstruction of the Higgs mass is worse for both m recoil and m p recoil distributions due to the wrong jet-pairing. However, m p recoil still has a better performance than m recoil , so its usefulness is not washed out by the combinatorial problem. We also note that, due to the lack of jet energy correction mentioned in Section 2, our distributions of m recoil peak around 130 GeV rather than 125 GeV. While the central values of the distributions can be corrected, we expect m p recoil to still have a better performance than m recoil after the implementation of jet energy corrections due to the parametric suppression on the Our study shows that the m p recoil variable could better reconstruct the Higgs mass for the signal and does not have the same effect on backgrounds. As such, we expect it to provide a significant improvement on the inclusive cross section measurements of the Higgsstrahlung process compared with the conventional recoil mass variable m recoil . Needless to say, such an improvement is crucially relevant to the studies of the Higgs boson properties. We also find a similar behavior of the signal and background distributions at the ILC 250 GeV, the results of which are not specifically shown. Since we have only performed a simplified simulation analysis and have not considered some of the important backgrounds, we will restrain ourselves from doing any quantitative analysis on the inclusive σ(hZ) measurements and leave it for the experimental groups who have better tools for such an analysis.

Conclusions
In this paper, we explore the use of the recoil mass and its variations in the measurements of the W W -fusion process at a lepton collider with a center of mass energy of 240-250 GeV. We found the variable m p recoil , constructed using only the 3-momenta of the Higgs decay products, can better separate Higgsstrahlung events with an invisible Z from the W W -fusion events than the original recoil mass m recoil does, with an improvement up to 20-30% on the precision of the W W -fusion cross section. We study its impact in both the conventional framework and the effective-field-theory one. In the conventional framework, a better precision on the W W -fusion cross section leads to a significant improvement on the constraints of the hW W coupling and the total Higgs width. In a global analysis under the effective-field-theory framework, using the information in the m recoil or m p recoil distributions could improve the reach on some of the EFT parameters by more than 10% compared with just using the total rate of the e + e − → ννh channel. We find that fitting the EFT parameters directly to the binned distributions gives the best precision reach. On the other hand, if the EFT parameters are fitted to the precisions of the W Wfusion and hZ, Z → νν cross sections extracted from the m recoil distribution, the precision reach could suffer from the large correlation between the two cross sections if it is not taken count of. We also explore the use of m p recoil in the inclusive measurements of the Higgsstrahlung process (e + e − → hZ) with hadronic Zs and find that it can significantly improve the reconstruction of the Higgs at a center of mass energy of 240-250 GeV. The use of m p recoil could therefore potentially lead to an improvement on the overall precision of the inclusive hZ cross section measurements. The construction of m p recoil is also extremely simple and does not require any additional measurements. It should be straightforward to implement m p recoil in any studies that make use of the recoil mass distribution.
where the Higgs boson couplings to a pair of SM gauge bosons are given by The parameters in Eq. (A.2) are not all independent. Four constraints can be written down by imposing gauge invariances, which we choose to rewrite δc W , c W W , c W and c γ as where δm can only be induced by custodial symmetry breaking effects and is set to zero in our framework. For the Yukawa couplings, we focus on the ones of t, c, b, τ, µ, parameterized as The possible flavor violating Yukawa couplings from new physics are not considered. The anomalous triple gauge couplings (aTGCs) are parameterized as Gauge invariance further imposes the relations δκ Z = δg 1,Z − t 2 θ W δκ γ and λ Z = λ γ . This leaves three independent aTGC parameters, which we choose to be δg 1,Z , δκ γ and λ Z . Two of them, δg 1,Z and δκ γ , are related to the Higgs parameters and can be written as To summarize, in our framework the contribution from D6 operators to the Lagrangian in Eq. (A.1) can be parametrized by the following 12 parameters: δc Z , c ZZ , c Z , c γγ , c Zγ , c gg , δy t , δy c , δy b , δy τ , δy µ , λ Z . (A.7) Also following Ref. [3,28,24], we consider the EFT contribution to the hγγ and hZγ vertices at the tree level, in which case the only EFT parameter that contributes to the decay rate of h → γγ (h → Zγ) is c γγ (c Zγ ). For the decay h → gg, we include, in addition to c gg , the contributions of δy t and δy b , which enter the hgg vertex by modifying the Yukawa couplings in the fermion loops. It is also convenient to normalize c γγ , c Zγ and c gg with respect to the SM 1-loop contributions. We follow Ref. [3] and define the following parameters Furthermore, without measuring the tth process at high energies ( √ s 500 GeV) or at the LHC, the parameters c gg and δy t can not be independently constrained. Since we focus on the 240-250 GeV run at lepton colliders, we replace c gg and δy t byc eff gg in Eq. (A.9) which parametrize the total contribution to the hgg vertex. The number of parameters is thus reduced to 11, and the parameters are δc Z , c ZZ , c Z ,c γγ ,c Zγ ,c eff gg , δy c , δy b , δy τ , δy µ , λ Z , (A. 11) which are used in our EFT global analysis in Section 3.

B EFT expressions for the e + e − → ννh cross sections
We obtain the cross section of e + e − → ννh as a function of the EFT parameters by generating events using Madgraph5 [14] with the BSMC package [29,30]. The events are showered in Pythia [15] and passed to Delphes [16] with the ILD card for detector simulations, after which the selection cuts in Section 2 are applied. Here we do not impose the gauge invariance condition from Eq. (A.3) in order to show the different dependences on the Z and W parameters. For the binned differential distributions of m recoil (m p recoil ), the numerical coefficients in Eq. (B.1) are replaced by the ones in Table 4 (Table 5). We have also checked that the statistical uncertainties from sim-   Table 5: Same as Table 4 but for the m p recoil distribution.
ulation are under control. 7 We then impose the gauge invariance condition in Eq. (A.3) and construct χ 2 of the e + e − → ννh measurement, assuming the events follow a poisson distribution. For the total rate, we have where N sig and N bg are the number of signal and backgrounds after cuts, normalized to 5 ab −1 for CEPC. For the binned distributions, we use Eq. (B.2) to construct the χ 2 of each bin where N sig and N bg are the number of signal and backgrounds in the bin. We then sum over the χ 2 of all the bins, assuming no correlation among them. The χ 2 is then combined with the ones of other measurements for the global analysis in Section 3. We refer the readers to Ref. [3] for a complete set of expressions for the other relevant observables.