Search for CP violation effects in the $h\to \tau\tau$ decay with $e^+e^-$ colliders

A new method is proposed to reconstruct the neutrinos in the $e^+e^-\to Zh$ process followed by the $h\to\tau\tau$ decay. With the help of a refined Higgs momentum reconstruction from the recoiling system and the impact parameters, high precision in the determination of the momentum of neutrinos can be achieved. The prospect of measuring the Higgs CP mixing angle with the $h\to\tau\tau$ decay at future $e^+e^-$ colliders is studied with the new method. The analysis is based on a detailed detector simulation of the signal and backgrounds. The fully reconstructed neutrinos and also other visible products from the tau decay are used to build matrix element (ME) based CP observables. With 5 $ab^{-1}$ of data at $E_{\text{CM}}=250$ GeV, a precision of $2.9^\circ$ can be achieved for the CP mixing angle with three main one-prong decay modes of the taus. The precision is found to be about 35\% better than the other methods.


Introduction
To explain the matter-antimatter asymmetry in the current observed Universe, C and CP symmetry violation is one of the three necessary conditions listed by Sakharov [1]. In the electroweak baryogenesis models, the CP violation phase from the CKM matrix is too small to accommodate the observed baryon-antibaryon imbalance. Hence the Standard Model has been extended to introduce extra new CP violating sources, such as in the 2HDM, Left-Right Symmetric and SUSY models.
On the other hand, the recently discovered Higgs [2][3][4] has opened a new window toward the new physics. The precision measurement of the Higgs properties will be one of the most important targets of the next generation collider including high-luminosity upgrade of LHC [5] and e + e − colliders [6][7][8][9]. Besides the determination of the signal strengths of the Higgs for different channels [10][11][12], the large data that will be accumulated in future colliders would allow us to measure the CP structure of some certain interaction with extraordinary precision.
Using 8 TeV LHC data, ATLAS and CMS have already measured the CP property of the Higgs coupling with the Z boson using the four-lepton decay channel (h → ZZ → ( )( )). A pure CP odd assumption is excluded at 99% confidence level [13][14][15]. On the other hand, a CP violating Higgs sector can also be searched in the fermionic decays of Higgs, using an effective Lagrangian of the following a e-mail: xin.chen@cern.ch b e-mail: ycwu@physics.carleton.ca form: where φ is the mixing angle between the CP even and odd terms. The advantage of this decay is that, unlike the bosonic decays of Higgs where the pseudoscalar interaction may arise from a dimension-6 operator such as hZ µνZ µν , the pseudoscalar term in Eq. 1 is dimension-4 and the coefficient can be large, so that a large CP violating effect is possible when the Higgs is a CP mixture. Although difficult to reconstruct, the tau decays are complex enough to yield information of the tau spin correlation [16,17]. Hence h → τ + τ − is a golden channel for the measurement of the CP nature of the Higgs which has led to an enormous amount of studies on the prospects of the measurement of the CP mixing angle(φ) at LHC using either the ggF/ZH or the VBF production mode [18][19][20][21][22][23][24][25][26][27] as well as at e + e − colliders [28][29][30].
The tau decay will always produce at least one neutrino which can hardly be reconstructed at hadron colliders. 1 Hence, the proposed methods at LHC utilize either only the visible products (such as the four pions in the ρρ mode) or the impact parameter of the tau [23], all of these can be reconstructed without knowing the full information of the neutrino. However, the resolution of the hadron collider will significantly affect the reconstruction of the CP sensitive observable [26]. Furthermore, besides the clean background, the lepton collider can also provide the possibility to fully reconstruct the neutrinos. Thus the lepton collider will be a better facility for the study of the CP properties of the Higgs.
In this paper, we propose a new method that can be used to reconstruct the momentum of the neutrinos from tau decay at future e + e − machine. In this method we utilize all the information including mass peaks of all involving particles as well as the impact parameters of each tau decay which will be measured with higher precision at lepton collider than at hadron collider. We find that using the new proposed method we can precisely determine the momentum of the missing neutrinos. With these fully reconstructed neutrinos as well as the visible decay products we construct a CP sensitive observable based on event by event matrix element to test the CP structure of the coupling between Higgs and τ leptons. After performing a detailed Monte Carlo (MC) simulation of the signal and corresponding backgrounds at e + e − collider, we have checked capability of the matrix element-based observable and also compared with the usually used CP observable defined as the angle between some certain planes.
The rest of the paper is organized as follows. In Sect. 2 we present in detail the method that we used to reconstruct the neutrino. We first refine the Higgs momentum from recoiling Z boson, and then determine the momentum of neutrinos by a global fitting utilizing impact parameter and on-shell conditions. A detailed MC simulation study of this method is given in Sect. 3. We use a matrix element-based CP observable to give the prospect of the measurement of the CP mixing angle. A comparison with previously used method is also given in Sect. 3. We summarize our discussion in Sect. 4.
2 The reconstruction of neutrinos at e + e − collider The leading Higgs production channel at e + e − collider is In order to study the CP property of the hτ τ coupling, after the Higgs decays into τ pair (including the CP violation effect as stated in Eq. (1)), and the Z boson decays either leptonically (Z → , = e/µ) or hadronically (Z → jj), the following one-prong decay channels for τleptons will be used: For clarity, before introducing the reconstructed method, we will first introduce details of the simulation and basic reconstruction for signal and background processes. The charged particle in each process leaves a track inside the detector which is used to construct the impact parameters. The τ mass will be used in the reconstruction of the neutrino, as well as the ρ mass in the ρ decay channel. Further, the mass of the Higgs will also be used to provide more information. Note that one can construct the Higgs peak from either the momentum of τ pair or the recoil mass constructed by the momentum of the Z boson, and how well the Higgs momentum can be reconstructed will significantly affect the further determination of the neutrino momentum. Thus the Higgs momentum reconstruction from the recoil system will also be investigated.

Signal and background simulation and reconstruction
The signal and background events are generated with Mad-Graph5 [31] and passed to Pythia8 [32] for the resonance decay (CP mixing Higgs, tau and rho mesons) and parton shower. The spin correlation between two taus is retained during the whole procedure. The signal process as demonstrated above is e + e − → Zh (σ = 212 fb), at a center of mass energy E CM = 250 GeV, with h → τ τ and Z → or Z → jj(2 jets). The dominant backgrounds for this signal are following: The events are afterwards passed through the DELPHES [33] simulation using the ILD card based on the TDR of ILC [7], in which the tracking efficiency for charged tracks is 99% with fiducial pseudorapidity range up to |η| = 2.4; the identification efficiency for electrons and muons is 95%, the calorimeter towers are simulated with a particle flow algorithm with a coverage up to |η| = 3.0. The track momentum is smeared with a resolution of 0.01 2 + (10 −4 p T ) 2 (p T in GeV) in a magnetic field of 3.5 T, and a resolution of 0.001 in η and φ (azimuthal angle) for its direction. The calorimeter energy is smeared with a resolution of √ A 2 E 2 + B 2 E, where the coefficients for the constant and stochastic terms are A = 1.0% (1.5%) and B = 15% (50%), respectively, for the EM (hadronic) calorimeter. To reject non-prompt leptons from jet fragmentation and heavy flavor meson decay, a lepton isolation cut of p PF T /p T < 0.7 is applied, where p PF T is the sum over particle flow objects (tracks and calorimeter clusters) around the lepton with p T > 0.5 GeV and ∆R = ∆η 2 + ∆φ 2 < 0.4 with respect to the lepton. This efficiency of this cut ranges from 91% in the low p T (5 GeV< p T < 10 GeV) to 99% in the high p T (> 30 GeV) region, and mostly affects the soft leptons from tau decays.
To have the best impact parameter resolutions which will be used to reconstruct the neutrinos, a minimum p T of 5 GeV is applied on the leptons and charged pions from taus. The resolution of the impact parameter in the transverse plane (in the z-axis) is set to σ d = 5 µm (σ z = 10 µm) as in Eq. 7, which can be achieved in the next generation e + e − colliders [7,8].
The hadronic taus (τ h ) are clustered by the Anti-k t jet algorithm [34] with a cone parameter of 0.4. A simple flat τ -tagging efficiency of 60% (0.5%) is applied on the real (fake) taus, taking into account the current τ -tagging performance at the LHC experiments [35,36]. With this fake tau rate, the W + W − → 4j background with two jets faking two taus can be neglected. When the Z boson decays hadronically, the particle flow objects associated with the selected electrons, muons and hadronic taus are firstly removed from the event, and then a coneless exclusive (k t ) algorithm [37] is applied on the remaining objects to cluster up to two jets, which are assumed to originate from Z → jj. The four-momentum of Z can be best reconstructed in this way, which further helps to calculate the Higgs mass recoiling against Z with the given CM energy.
The efficiency of finding two jets with p T > 5 GeV and |η| < 3 is about 95%. The net efficiencies after the above object selections in different di-tau decay modes and the requirement that leptons from the h/Z decay should have opposite charge are listed in Tab. 1. To measure the CP mixing angle, it is essential to identify the different tau decay modes efficiently. The development of tau substructure algorithms in the ATLAS and CMS experiments has recently made this possible based on the particle flow method [38,39]. The neutral pion energy can also be measured with a precision of 15% with the substructure. In this work, we assume that different tau decay channels can be efficiently identified in the future e + e − colliders, with no crosstalk for simplicity, and a precision of 10% on the neutral pion energy measurement can be achieved.

Refined Higgs four-momentum reconstruction
In general, the reconstructed Higgs peak will always have a long tail as we ignored the initial state radiation (ISR) photons. To best reconstruct the Higgs four-momentum from the recoiling Z boson, taking into account the known Higgs mass at 125 GeV, and assuming that the ISR photons are mostly collinear with the beam, a quantity x can be solved which represents the fraction of momentum carried away by the ISR photon from the positron beam traveling in the positive z-direction, where E, m and p z are the energy, mass and z-component of the momentum of the recoiling Z boson. The value of x can also be negative, which means that the ISR photon carries momentum from the electron beam along the negative z-axis. 2 With the x solved, the four-momentum of the Zh system can be expressed as p tot = (E, p x , p y , p z ) = (250-125|x|, 0, 0, -125x) GeV, from which the Higgs recoil four-momentum can be calculated as p RC h = p totp Z . Because the tau mass is much smaller than the Higgs mass, the tau from Higgs decay is highly boosted and the neutrinos from tau decay are almost collinear with its visible decay products. With this assumption, the Higgs four-momentum can also be expressed as p h = p vis1 /x 1 + p vis2 /x 2 , where p vis1,2 are the four-momenta of visible decay products of the two taus, and x 1,2 are the fractions of momentum carried by the visible products. With these quantities defined, a χ 2 can be minimized per event to get the best momentum resolution: where m Z is the Z boson mass constructed from the decay products of the Z boson (two leptons or two jets), f j1,2 are the jet energy correction factors on the two jets from Z → jj decay, and variables with energy dimensionality are in GeV. The jet energy resolution used here is about 6% and this is a conservative choice and can be further improved in future electron colliders, such as ILC [6,7], CEPC [8] and CLIC [9] (4%). By minimizing the χ 2 in Eq. 5 with x 1,2 and f j1,2 freely floating, not only the twofold ambiguity of x in Eq. 4 is resolved (the solution with a smaller χ 2 is chosen), but also the jet resolution, and hence the resolutions of p RC h are improved. The improvement is shown in Fig. 1. The denominator values in Eq. 5 are chosen such that the best resolutions as in Fig. 1 are achieved. In the Z → channels, the last two terms in Eq. 5 are not present since the leptons' momenta are precisely measured when compared to the jets. Eq. 5 fulfills two purposes in this case. One is to resolve the x ambiguity, and the other is to resolve the two-fold ambiguity of the e + e − e + τ h and µ + µ − µ + τ h final states. There are two ways to assign the two same-sign leptons to the Z and h decays. With current setup, the fraction of wrong assignments is negligible.

Impact parameter
Ideally, except for the lepton channels, using the missing momentum, the Higgs four-momentum and the tau mass, one already has enough constraints to reconstruct the two neutrinos. However, in practice, the reconstructed momentum of the neutrino still has large uncertainty. Thus the usage of the further information from the impact parameters [40, 41] will be absolutely helpful in the reconstruction of the neutrinos especially in lepton channel. In order to make use of the impact parameters, a χ 2 IP is reconstructed for each tau which will be minimized in further global fittings; this will be explained in the next section.
The method to calculate χ 2 IP for one tau is as follows: given the transverse impact parameter d and the direction of tau momentum in the transverse plane, the intersection point between the tau flight direction and the track trajectory in the transverse plan can be found. The transverse impact parameter d is defined in the transverse plane as the minimum distance from the interaction point to the charged track. As demonstrated in Fig. 2, the tau is produced at the collision point O, flight length from O to the point of decay P in the transverse plane is L, the point of closest approach to O with a backward extrapolation of the trajectory is D, and the extrapolated arc length in transverse plane is S. There is just one intersection point when O is inside the circle of the trajectory (Fig. 2(a)), but when O is outside, the number of intersection points can be 2, 1 or 0 ( Fig. 2(b,c)). The distance between O and D is the impact parameter d 0 . The impact parameter in the z-axis can be calculated from P : The impact parameters will get updated values (d fit 0 and z fit 0 ) when χ 2 IP is minimized. When just one intersection point exists, the relevant contribution reads When two intersection points exist, the χ 2 for each possibility is calculated according to Eq. 7 and the smaller one is taken. When no intersection point exists, as demonstrated by the dashed arrow in Fig. 2(c), it is assumed to be due to the uncertainty in d 0 measurement, and the best-fit d fit 0 will be around the distance between O and D (denoted d C 0 ). The O is determined by translating the dashed vector to become tangential to the trajectory (the blue vector), and z fit 0 is calculated from the tangential point P . The relevant contribution then reads When the dashed vector in Fig. 2(c) points away from the track curvature, no tangential point can be found, and O coincides with D. Although the tracks in Fig. 2 all travel anti-clockwise, but depending on the particle charge, they may also travel clockwise, i.e., from P to D. Situations with tracks traveling clockwise can be obtained by taking a mirror image of the plots in Fig. 2, with essentially no change to Eq. 7 and 8. The cases where the angle between the tau flight direction and the track trajectory is an obtuse angle are also considered in our method. However, since the tau is boosted and its decay products are highly collinear, this rarely happens.

Reconstruction of the neutrino momenta
In order to best estimate the neutrino momenta, with the benefit of accurate particle momentum reconstruction and clean background from a e + e − collider, the direction (in η, φ) and the magnitude of each of the two neutrinos from tau decays are scanned globally for each event.
The available information such as the tau mass, the Higgs four-momentum from Z recoil, and the impact parameters of the charged tracks from tau decays, are all used to achieved this goal. With the four-momenta of all final state particles reconstructed, a matrix element-based method is used to fully extract the information of CP, which can result in an improved sensitivity with respect to the usually used method [23].
The global χ 2 of the likelihood that is minimized for every event is expressed as where p h = p vis1 + p mis1 + p vis2 + p mis2 is calculated from the four-momenta of visible and invisible decay products of the two taus (without any collinear assumptions), p RC h is obtained by minimizing Eq. 5 in the previous step, m τ 1,2 are the masses of the taus, and χ 2 IP is the term accounting for the contribution from the impact parameters of the tracks that can help the fit to find the correct neutrino directions as stated in Sect. 2.3. For Z → (jj) channels, σ RC = 0.5 (4.0) GeV and σ τ = 0.1 (0.2) GeV. The resolution parameters are set so as to achieve the best reconstructed neutrino momenta (magnitude and direction close to the true values). For the states with an intermediate ρ meson, extra terms of (m ρ − 0.775) 2 /σ 2 ρ + (f ρ − 1) 2 /0.10 2 are added to Eq. 9, where f ρ is the factor multiplied to the ρ meson's energy for a better resolution. Correspondingly, σ τ is doubled for this case, and σ ρ = 0.15 (0.30) for Z → (jj).
In the per-event χ 2 minimization, the (η, φ) of one neutrino is firstly scanned, from which the magnitudes of the neutrinos' momenta and the direction of the other neutrino can be obtained from the tau mass and recoil constraints in Eq. 9. Conversely, the scan is repeated from the (η, φ) of the other neutrino. Finally, a fit using MI-NUIT [42] is performed around the minimal point found by the scans for a better estimation. The ∆R and momentum ratios between the reconstructed and true neutrinos in the π + ρ channel are shown in Fig. 3.
In the channels with a lepton (e/µ) from tau, the momenta of the two neutrinos from the same tau cannot be fully reconstructed. Instead, the combined four-momenta of them is reconstructed, with the di-neutrino mass as an extra degree of freedom (the relative angles between the two neutrinos are "integrated" out). Although there are eight constraints in Eq. 9, it is still difficult to achieve. In this case, another extra term, −2 ln P(∆R, m mis ), is added. It is the joint probability distribution of ∆R between the lepton and di-neutrino momenta, and m mis the di-neutrino invariant mass. This function is obtained depending on the tau momentum from the Monte Carlo simulation, which has been used in the di-tau MMC mass by ATLAS [43]. The quality of the neutrino momentum reconstruction in the +π/ρ and Z → channels are shown in Fig. 4. Due to m mis and worse Higgs recoil momentum resolutions, the +π/ρ channels with Z → jj have limited h,fit > 122 GeV 80 GeV< m fit Z <100 GeV 120 GeV< m h <130 GeV 1.5 GeV< mτ <2.0 GeV mρ >0.3 GeV (for channels with ρ)  and Z → jj decays before the cuts are applied. With 5 ab −1 of e + e − collision data at E CM = 250 GeV, the expected numbers of events entering the CP sensitivity test are listed in Tab. 3. Note that the whole per-event fitting procedure not only helps to determine the momentum of neutrino and also provides better resolutions for these mass observables and thus helps to further suppress the backgrounds.

Matrix element-based measurement of CP
With the fully reconstructed momentum of final states from the Higgs decay, a method based on the matrix ele-  Table 3. The expected numbers of signal and background events with 5 ab −1 of data in each channel after the cuts in Tab. 2 with Z → , jj combined. Note that the + π/ρ and Z → jj channels are excluded.
ments could be used to probe the CP mixing angle. From Eq. 1, the matrix element squared with φ being the CP mixing angle can be expressed as: (10) where I 1 = A + B, I 2 = 2C and I 3 = A − B. For example, in the π +π channel, the coefficients have relatively simple expressions: where p τ − and p τ + are the momenta of the two taus, which include the information of the reconstructed neutrinos and p1p2p3p4 is short for µνρσ p µ 1 p ν 2 p ρ 3 p σ 4 . The full expression of the coefficients for the other decay mode combinations can be found in Appendix A where we also list the effective Lagrangian used for tau decay. For channels involving leptonic decays of tau, a phase space integral is performed on the internal degrees of freedom between the two neutrinos, at the cost of a bit loss of sensitivity to CP. In the ρ + ρ and π + ρ channels, the neutral pion(s) is assigned to the corresponding decay chain without any ambiguity, since the neutral pion from different decay chain is highly collimated with the original tau. An optimal obsevable [44,45], defined as OO = I 2 /I 1 for each event, is used to distinguish signals with different CP mixing angles. For signals with a positive (negative) φ, the mean of the OO distribution will be shifted to negative (positive) values, as shown in Fig. 6(a).
Template probability density functions (PDF) for different CP angles are first obtained from the simulations. A binned likelihood function (L) is then calculated from the  Fig. 4. The ∆R between the reconstructed and true neutrino (or di-neutrino) momenta for the τ → νν (a) and τ → πν + ρν (b) decays, the ratio of the reconstructed to true di-neutrino momentum for the τ → νν decay (c), and the 2-D distribution of the true versus reconstructed di-neutrino mass (d), in the + π/ρ and Z → channels.
PDF as a function of the CP mixing angle φ. The bestfit CP angle is found with the least negative logarithm of the likelihood (NLL), and the 1σ confidence interval is obtained by finding the angles whose NLLs are 0.5 higher than the minimum NLL. Fig. 6(b) shows the expected ∆NLL as a function of the CP angle φ, with all five channels in Tab. 3 combined. With 5 ab −1 of e + e − collision data (can be achieved by the CEPC [8]), a precision of 0.05 radians can be reached. On the other hand, if the integrated luminosity is 2 ab −1 (can be achieved by the ILC [46]), a precision of about 0.09 can be reached.
Actually, one can also built an angle ∆φ M E (−π ≤ ∆φ M E ≤ π) to probe the CP mixing angle which satisfies: With this definition, the matrix element (which also determines the distribution of ∆φ M E ) becomes The distribution of ∆φ M E is presented in Fig. 7(a) with two different choice of CP mixing angle φ. ∆φ M E makes full use of the information inside the matrix element and performs better when φ is large. Fig. 7(b) shows the NLL comparison between OO and ∆φ M E , we find that when φ is large, ∆φ M E is better than OO (as in the OO construction, the term I 3 is omitted). But in the small φ region, the two methods give quite similar results.
Further, the ρ+ρ channel is used to compare the sensitivities of the ∆φ IP method (using the impact parameters to approximately reconstruct the tau decay planes), the ∆φ CP method (using the visible ρ meson decay products to reconstruct the CP sensitive angle ∆φ CP ) and the OO method, the details of the definitions for ∆φ IP and ∆φ CP are given in Appendix B. Figure 8(a, b) shows the distributions of ∆φ IP , ∆φ CP in this channel. It is seen that the IP method have much worse separation power between CP even and CP mixed states than the ∆φ CP method. Figure 8(c) shows the ∆NLL variations with the ∆φ CP and the OO, from which one can read off that the latter is about 35% better than the former. (a1, b1, c1) and Z → jj (a2, b2, c2) decays with the three tau decay modes considered and assuming 5 ab −1 of data. Note that the + π/ρ and Z → jj channels are excluded.

Conclusion
Understanding the CP property of the Higgs is one of the primary goals in future e + e − colliders or Higgs factory. If the Higgs is found to be a CP mixture, a door to matterantimatter imbalance and New Physics would be opened.
The H → τ τ decay has a unique position in the Higgs CP search, as the CP odd contribution can enter the Lagrangian at the tree level.
In this paper, we proposed a new method utilizing the impact parameters and the on-shell conditions to fully reconstruct the momentum of the neutrinos from the tau decay (in the leptonic decay channel, the sum of the two neutrinos). With fully reconstructed the momenta of final state particles from the tau decays using a global likeli- hood scan and fit, the Higgs CP information can be best estimated based on a matrix element method. In the reconstruction, the resolution degradation due to ISR photons can be reduced by a pre-fit of the Higgs' momentum from the Z recoil, which takes into account the Higgs mass and the special H → τ τ decay topology.
We have performed a complete MC simulation taking into account the performance of current LHC detectors for the signal process e + e − → Zh → ( /jj)(τ + τ − ) followed by three main one-prong decay channels of tau and also the corresponding backgrounds. The refined Higgs momentum can further help to reduce the background leaving enough signal rates. By virtue of fully reconstructed neutrinos, the matrix element is calculated event by event for all di-tau decay modes, and matrix element based observables are built (OO and ∆φ ME ) to probe the CP mixing angle φ of the hτ τ coupling. At larger value of φ, ∆φ ME is definitely better than OO, but for small values of φ, the performances of them are found to be similar.
Template PDFs are obtained from the simulation, and the NLL analysis based on the PDF shows that with 5 ab −1 (2 ab −1 ) of data at E CM = 250 GeV center of mass energy from e + e − collisions, a CP mixing angle can be measured to a precision of 0.05 (0.09) radians, or 2.9 • (5.2 • ). The comparison with other CP observables (∆φ IP and ∆φ CP ) has also been presented, and we find that the CP mixing angle sensitivity reach based on our method is about at least 35% better than the other methods. The method can also be extended to other collider analyses that are sensitive to the neutrino momentum reconstruction.

A Matrix element expressions for the channels considered in this work
The effective Lagrangian we used for the calculation of the matrix elements is +C ρ (τ γ µ P L ν τ (π 0 ∂ µ π − − π − ∂ µ π 0 ) + h.c.) For each Channel, the matrix element square has following form: Note that in the + π/ρ channels, when constructing the Optimal Observable, the internal degrees of freedom between the two neutrinos are integrated out in the corresponding coefficients, leaving only the combined fourmomentum of the di-neutrino in the expressions. In the ρ + ρ channel, the two neutral pions are taken to be distinguishable.