Signature of 2 HDM at Higgs Factories

The full one-loop corrections, both the weak and QED corrections, to the process e+e− → Zφ (φ = h0,H0) at the Higgs factories are presented. Up to the O(αEM ) level, the virtual corrections are evaluated in a Feynman-diagrammatic approach. The real emission corrections are computed using Feynman Diagram Calculation (FDC) system and the collinear divergences are regularized by the structure functions of electron. Using the FDC system, we study the corrections both in the SM and the two Higgs doublet model (2HDM), respectively. After taking into account experimental constraints from the current LHC data, we propose four interesting benchmark scenarios of 2HDM for future colliders. By using these benchmark scenarios, we evaluate the deviation of ∆σ(e+e− → Zφ) from their Standard Model (SM) values. We also examine φ → bb̄ and φ → τ+τ−, which may receive large EW contribution from triple Higgs couplings which are absent in the SM. It is found that for these benchmark scenarios, both EW and real emission corrections are sizable and could be measured at a future e+e− collider such as the ILC, CLIC and CEPC. xiewh@ihep.ac.cn r.benbrik@uca.ac.ma ahabjia@gmail.com twain@ihep.ac.cn yanqishu@ucas.ac.cn


Introduction
A Higgs-like boson (h) has been discovered in the first Run of the Large Hadron Collider (LHC) in 2012 [1,2]. Based on the data of Run-1, ATLAS and CMS collaborations has established Higgs boson mass with m h = 125.09 ± 0.21 (stat.) ± 0.11 (syst.) GeV [3]. LHC also performed several Higgs coupling measurements to a precision around 10 − 20%, such as Higgs couplings to dibosons V V with V = W ± , Z, γ. Very recently, LHC have also measured Higgs couplings to the fermions of the third generation over 5σ, via the process pp → V h(h → bb) [4,5], via the VBF process pp → jjh(h → τ + τ − ) [6,7], and pp → tth with combined Higgs boson decay final states [8,9]. These measurements demonstrated that the SM works quite well in the current Higgs data.
One of the tasks of the LHC Run-2 at 13 TeV (and 14 TeV) with High Luminosity (HL) would be to improve all the aforementioned measurements and to perform new ones such as accessing h → γZ as well as the triple self-coupling of the Higgs boson. It is expected that the new LHC Run-2 will pin down the uncertainty in h → bb and h → τ + τ − to 10-13% and 6-8% for bottom quarks and tau leptons, respectively. These measurements will be further ameliorated by the High Luminosity option for the LHC (HL-LHC) down to uncertainties of 4-7% for b quarks and 2-5% for τ leptons [12]. Moreover, in the clean environment of the e + e − Linear Collider (LC), which can act as a Higgs factory, the uncertainties on h → bb and h → τ + τ − would be much smaller reaching 0.6% for the couplings in h 0 → bb and 1.3% for those in h 0 → τ + τ − [13,14].
The above accuracies on fermionic Higgs decay measurements, if reached, are of the size comparable to the effects of radiative corrections to some Higgs decays. Therefore, one can use these radiative correction effects to distinguish between the Standard Model (SM) and various beyond-standard models (BSM). In this respect, precise calculations of Higgs-boson production and decay rates have been performed already quite some time ago with great achievements (see e.g. [15,16]). QCD corrections to Higgs decays into quarks are very well known up to O(α 3 s ) as well as additional corrections at O(α 2 s ) that involve logarithms of the light-quark masses and also heavy top contributions [15]. Electroweak radiative corrections to fermionic decays (bb and τ + τ − ) of the Higgs boson in the SM are also well established [17][18][19] in the on-shell scheme. In the framework of the Two-Higgs-Doublet Model (2HDM), several studies have been carried out to evaluate the electroweak corrections to fermionic Higgs decays [20,21]. The calculation of Ref. [20] is done in the on-shell scheme except for the Higgs field renormalization where the MS subtraction has been used, while the one of Ref. [21] is performed using the on-shell renormalization scheme of [22].
At the LHC, due to the large theoretical uncertainties like the PDF and large experimental background, the precision measurements of the Higgs boson are rather challenging. In contrast, e + e − colliders can offer us precision measurements on the production and decay properties of Higgs boson. At an 240 GeV e + e − collider, the Higgs-strahlung process e + e − → Zh is the dominant production channel for the Higgs boson, for which properties of Higgs boson can be reconstructed by using the recoiled Z boson via its leptonic decay. For a center-of-mass energy of 240-250 GeV and an integrated luminosity of 250 fb −1 , O(10 5 ) Higgs bosons per year will be produced which can lead to a measurement of the Higgs couplings at percent level [23][24][25]. At the international linear collider (ILC) expriments [26], the luminosity is expected to be high and a quite small experimental error is expected. Correspondingly, precise theoretical predictions to the physical observables related to Higgs boson are required.
In this work, we will propose a few benchmark scenarios of 2HDM after taking into account theoretical constraints as well as experimental restrictions from recent LHC data. By using these benchmark scenarios, we will study the effects of electroweak radiative corrections to the production process e + e − → Zφ and the decay processes φ → bb and φ → τ + τ − . We will also report our full next leading order calculation on the cross section of e + e − → Zφ with φ = h 0 , H 0 in 2HDM by including both virtual correction up to one-loop and real emission of a photon. Our results are consistent with those given in the previous work [27] under the same setting of 2HDM parameters. The complete SM 1-loop correction of σ(e + e − → Zh SM ) have already been calculated using the GRACE [28] and we reproduced the results for comparison. We noticed that the EW corrections of new physics in our benchmark scenarios can be of order −10% ∼ −20% even after taking into account the recent constraints LHC data. While the contribution of real emission has positive sign.
To examine the deviation caused by the new physics, we will also evaluate the ratio of branching fractions of Higgs decays in the 2HDM [29,30], This ratio of Higgs boson decay widths is independent of the production process, but it is sensitive to the electroweak corrections of new physics that may not affect the production rate of Higgs bosons. This ratio is less sensitive to the systematic errors (which drop out in the ratio) so it could be useful to discriminate the SM against other models such as 2HDM or supersymmetric models.
The paper is organized as follows. In section 2, we review the Yukawa textures, scalar potential and Higgs self-couplings of the 2HDM model, as well as the theoretical and experimental constraints on the model. In section 3, we outline the framework of our calculation and specifies the renormalization scheme we will use. In section 4, we present the numerical results. We conclude this work in section 5.

The 2HDM model 2.1 Yukawa textures
In the 2HDM, fermion and gauge boson masses are generated from two Higgs doublets Φ 1,2 where both of them acquire vacuum expectation values v 1,2 . If both Higgs fields couple to all fermions, Flavor Changing Neutral Currents (FCNC) are generated which can invalidate some low energy observables in B, D and K physics. In order to avoid such FCNC, the Paschos-Glashow-Weinberg theorem [31] proposes a Z 2 symmetry that forbids FCNC couplings at the tree level. Depending on the Z 2 assignment, we have four type of models [16,32].
where Φ d,l (d, l = 1, 2) represents Φ 1 or Φ 2 and Y f (f = u, d or ℓ) stands for Yukawa matrices. The two complex scalar SU(2) doublets can be decomposed according to where v 1,2 are the vacuum expectation values of Φ 1,2 . The mass eigenstates for the Higgs bosons are obtained by orthogonal transformations,  with the generic orthogonal matrix From the eight fields initially present in the two scalar doublets, three of them, namely the Goldstone bosons G ± and G 0 , are eaten by the longitudinal components of W ± and Z, respectively. The remaining five are physical Higgs fields, two CP-even H 0 and h 0 , a CP-odd A 0 , and a pair of charged scalars H ± .
Writing the Yukawa interactions eq. (2) in terms of mass eigenstates of the neutral and charged Higgs bosons yields with the W boson mass m W and the SU(2) gauge coupling constant g. We have used the notation s x and c x as short-hand notations for sin(x) and cos(x), respectively. The mixing angle β is defined by via tan β = v 2 /v 1 .
It has been shown that the 2HDM has a decoupling limit which is reached for cos(β−α) = 0 and m H 0 ,A 0 ,H ± ≫ m Z [33]. In this limit, the coupling of the CP-even h 0 to SM particles completely mimic the SM Higgs couplings including the triple coupling h 0 h 0 h 0 . Moreover, the model possesses also an alignment limit [34], in which one of the CP-even Higgs bosons h 0 or H 0 looks like SM Higgs particle if sin(β − α) → 1 or cos(β − α) → 1. In the decoupling or alignment limit with α = β − π/2, the above triple Higgs couplings reduce to the simplified form

Theoretical and experimental constraints
The 2HDM has several theoretical constraints which we briefly address here. In order to ensure vacuum stability of the 2HDM, the scalar potential must satisfy conditions that guarantee that its bounded from below, i.e. that the requirement V 2HDM ≥ 0 is satisfied for all directions of Φ 1 and Φ 2 components. This requirement imposes the following conditions on the coefficients λ i [35,36]: In addition to the constraints from positivity of the scalar potential, there is another set of constraints by requiring perturbative tree-level unitarity for scattering of Higgs bosons and longitudinally polarized gauge bosons. These constraints are taken from [37,38]. Moreover, we also force the potential to be perturbative by imposing that all quartic coefficients of the scalar potential satisfy |λ i | ≤ 8π (i = 1, ..., 5). Besides these theoretical bounds, we have indirect experimental constraints from B physics observables on 2HDM parameters such as tan β and the charged Higgs boson mass. It is well known that in the framework of 2HDM-II and IV, for example, the measurement of the b → sγ branching ratio requires the charged Higgs boson mass to be heavier than 580 GeV [39,40] for any value of tan β ≥ 1. Such a limit is much lower for the other 2HDM types [41]. In 2HDM-I and III, as long as tan β ≥ 2, it is possible to have charged Higgs bosons as light as 100 GeV [41,42] while being consistent with all B physics constraints as well as with LEP and LHC limits [43][44][45][46][47][48].
We stress in passing that after the Higgs-like particle discovery, several theoretical studies have performed global-fit analysis for the 2HDM to pin down the allowed regions of parameter space both for a SM-like Higgs h 0 [49] as well as for a SM-like Higgs boson H 0 [50].
The parameter space of our benchmark scenarios in this paper is already constrained by the limits obtained from various searches for additional Higgs bosons at the LHC, and by the requirement that one of the neutral scalars should match the properties of the observed Higgs-like boson. We evaluate the former constraints with the public code HiggsBounds [51], and the latter with the code HiggsSignals [52].
Before presenting our results, we would like to mention that we have performed cross checks with the results [27] for the subclass of Yukawa corrections considered there and found perfect agreement. In order to illustrate the effect of the new physics, we have chosen four bencmhark poins scenarios as showen in Table.2 based on the best-fit from the latest results of Higgs data using HiggsBounds and HiggsSignals public codes.  It is noticed that in the scenarios of both BP1-h and BP1-H, a light charged Higgs boson less than 200 GeV can be consistent with the current LHC Higgs data.  In Table. (3), triple Higgs couplings are tabulated. It is noticed that the triple coupling λ h 0 H ± H ∓ for the BP1-h and BP2-h case could be quite large. While for the case BP1-H and BP2-H, the self coupling of λ h 0 h 0 h 0 could be sizable.
3 One-loop calculation and renormalization scheme Calculations of higher order corrections in perturbation theory in general lead to ultra-violet (UV) divergences. The standard procedure to eliminate these UV divergences consists in renormalization of the bare Lagrangian by redefinition of couplings and fields. In the SM, the on-shell renormalization scheme is well elaborated [53][54][55]. For the 2HDM, several extensions of the SM renormalization scheme exist in the literature [20][21][22]56]. Recently, gauge independent renormalization schemes have been proposed [57,58], with e.g. MS renormalization for the mixing angles and the soft Z 2 breaking term in the Higgs sector [58]. In the present study, we adopt the on-shell renormalization scheme used also in [20], which is an extension of the on-shell scheme of the SM: the gauge sector is renormalized in analogy to [53,54] concerning vector-boson masses and field renormalization; also fermion mass and field renormalization is treated in an analogous way (see also [55]).
and for the vector-boson fields By the on-shell definition of the electroweak mixing angle (s W = sin θ W , c W = cos θ W ) according to s 2 W = 1 − m 2 W /m 2 Z , the corresponding countreterm is determined through the W and Z mass countreterms, As shown in Ref. [55], the renormalization constant for charge δZ e is obtained as with This corresponds to the running coupling constant α(0) obtained at Thomson limit. As the first term in δZ e is sensitive to the hadronic contribution, usually a non-perturbative parameter ∆α (5) hadron (m Z ) is used to to absorb the hadronic contribution, namely δZ e is rewritten as Here we call it α(0). Another two schemes, α(m z ) and α( √ s), are more preferred, in which the large logarithm from leptons are also absorbed into the redefinition of running coupling constant [55,60]. Corresponding renormalization constant is And the running coupling constant is replaced with Results under these two schemes will be independent of log(m e ). In the following, we use α(m Z ) scheme as our default choice.
For renormalization of the Higgs sector we take over the approach used in [20], which means on-shell renormalization • for the h 0 , H 0 tadpoles, yielding zero for the renormalized tadpoles and thus v 1,2 at the minimum of the potential also at one-loop order, • for all physical masses from the Higgs potential, defining the masses whereas Higgs field renormalization is done in the MS scheme. We assign renormalization constants Z Φ i for the two Higgs doublets in (3) and counter-terms for the v i within the doublets, according to and expand the Z factors Z Φ i = 1 + δZ Φ i to one-loop order. The MS condition yields the field renormalization constants as follows for all types of models listed in Table 1: with ∆ = 2/(4 − D) − γ + log4π from dimensional regularization, the color factor N C = 3 for quarks and N C = 1 for leptons, and the gauge couplings g and g ′ . The factors ξ A 0 u,d,l can be found in Table 1. Eq. (21) is a generalization of the work of [20] with respect to the various Yukawa structures of the 2HDM.
The renormalized self-energy of the Higgs field φ is the following finite combination of the unrenormalized self-energy and counter-terms, with the on-shell mass counter-term Owing to the MS field renormalization, a finite wave function renormalization has to be assigned to each external φ in a physical amplitude. This quantity is determined by the derivative of the renormalized self-energy Σ ′ h 0 on the mass shell, given by Application to the one-loop calculation for the fermionic Higgs boson decay φ → ff yields the decay amplitude which can be written as follows, where ∆M 1 is the sum of the one-loop vertex diagrams V φff 1 and the vertex counter-term δ(φff ), Σ φφ ′ is the φ-φ ′ mixing (φ, φ ′ ) = (h 0 , H 0 ), δα represents the counter-term for the mixing angle α, andẐ φ is the finite wave function renormalization of the external φ fixed by the derivative of the renormalized self-energy specified above in Eq. (23). Given the fact that the mixing angle α is an independent parameter, it can be renormalized in a way independent of all the other renormalization conditions. A simple renormalization condition for α is to require that δα absorbs the transition φ -φ ′ in the non-diagonal part ∆M 12 of the fermionic Higgs decay amplitude. Therefore, the angle α is hence the CP-even Higgs-boson mixing angle also at the one-loop level, and the decay amplitude M 1 simplifies to the ∆M 1 term only.
The amplitude given in Eq. (24) together with its ingredients is a generalization of the work in [20], extended to all charged fermions and for the various 2HDM types. ∆M 1 contains besides the genuine vertex corrections the counter-term δ(φff ) for Higgs-fermionfermion vertex, which reads as follows, where PSfrag replacements l k l m l k l m Figure 1: Generic one-loop 2HDM Feynman diagrams contributing to Γ 1 (h 0 → bb) and can be expressed in terms of the scalar functions of the fermion self-energy, and the universal part This universality is a consequence of the renormalization condition (see the discussion in [20]), which is also used in the Minimal Supersymmetric SM (MSSM), see e.g. [62][63][64]. It is important that the singular part of the difference in the lhs. of (32) vanishes. The singular part of δv is (δv/v) MS = − 1 64π 2 (3g 2 + g ′2 )∆, which is equal to the expression found in the MSSM and constitutes a check of our calculation. We end this section by showing in Fig. (1) the one-loop Feynman diagrams in 2HDM for simplicity we draw h 0 → bb and h 0 → τ + τ − , where S stands for (H ± , A 0 , H 0 , G ± ) for both decays while F represents (b,t) for h 0 → bb and (τ, ν τ ) for h 0 → τ + τ − . In the SM limit [33], diagrams (4, 5, 10, 11) and (2,8) with (S,S)=(H ± , G ± ) vanish. Consequently, the important effects come from diagrams (1, 2) and (7,8) respectively for h 0 → bb and h 0 → τ + τ − .
In the present work, computation of all the one-loop amplitudes and counter-terms is done with the help of FeynArts and FormCalc [65] packages. Numerical evaluations of the scalar integrals are done with LoopTools [66]. We have also tested the cancellation of UV divergences both analytically and numerically. Before illustrating our findings, we first present the one-loop quantities that we are interested in. At one-loop order the decay width of the Higgs-boson into bb and τ + τ − is given by the following expressions, where β 2 = 1 − 4m 2 f /m 2 φ . We will parameterize the tree level width by the Fermi constant G F , i.e. we use the relation where ∆r incorporates higher-order corrections. According to the above relation, the oneloop decay width Eq. (33) becomes To parameterize the quantum corrections, we define the following one-loop ratios: where we also take the SM decay width Γ SM 1 (h → ff ) with the one-loop electroweak corrections. The two ratios defined above will take the following form: 3.2 Full correction of e + e − → Zh SM in SM In our calculation, a small photon mass λ is introduced to regularize the soft divergence. Meanwhile, two cutoffs, ∆E and ∆θ, are introduced to deal with the IR singularities in real correction processes. Thus NLO corrections can be expressed as where the soft part is given by and the HC + CT part is given by Both the soft and virtual parts are obtained with FormCalc, while the other parts are obtained with the help of FDC [67]. The λ independence have been checked when we combine the soft and virtual parts. The dependence on ∆E and ∆θ are shown in Tables 4 and 5. In Table 4, the dependence is seen in a wide range and we choose δ s = 10 −3 as our default choice. In Table 5, the result becomes cut dependent when ∆θ is smaller than 10 −4 . It is because the approximation in Eq. (40) demands ∆θ ≫ m e / √ s ∼ 2 × 10 −6 . Thus we choose ∆θ = 10 −3 as our choice.    One-loop radiation correction includes collinear singularities when m e goes to zero. In our calculation electron has nonzero mass, but the singularities will not disappear. They become terms proportional to log(m e ). Some of them are cancelled when summing up virtual and real corrections, some of them are absorbed into the redefinition of running coupling constant as mentioned above, but some are remained. To deal with this, we used following fixed order electron structure function which can be derived from Eq. (11) of Ref. [68]: with is regularized Altarelli-Parisi splitting function. One-loop structure function is used here instead of most-commonly-used resummed one to ensure the cancellation of collinear singularities. The cancellation is shown in Table 6. We vary the mass of electron with a factor of k from 2 −4 to 2 8 and find the result unchanged. Also, we can see that singular terms only appear in σ V +S and σ SC parts.   In this section we calculate the one-loop radiative corrections to the Higgs-strahlung process e + e − → Zφ in the SM and 2HDM. It tourns out to deal with corrections of the following corrections φZZ which are part of the correction to double Higgs strahlung process. Those couplings are given at the tree-level by: where s W = sin θ W and c W = cos θ W for the electroweak mixing angle θ W . It is clear that the processes e + e − → Zh 0 , ZH 0 are complementary to the Higgs pair production ones. At the tree-level no difference in the 2HDM as compared to the MSSM. Not like in Ref [27] where they studied e + e − → Zh, ZH in the alignement limit with sin(β − α) = 1 and small tan β, we instead based on the latest beft-fit of α and β from LHC Higgs data. We have calculated the radiative corrections to the tree level e + e − → Zh 0 , ZH 0 processes in 2HDM in the Feynman gauge including all the particles of the model in the loops.The Feynman diagrams from 2HDM contributing to e + e − → Zh 0 process (for simplicity we take φ = h 0 ) are shown in Figs. 2, 3, 4 and 5 for illustrations. Since we are dealing with a processes at the one-loop level, a systematic treatment of the UV divergences have to be considered. We will use the on-shell renormalization scheme in which the input parameters coincide with the physical masses and couplings. As performed before, the counter-terms are calculated by specific renormalization conditions which allow us to cancel the UV divergences of the diagrams with loops. Inserting these redefinitions into the Lagrangian, we find the following counter term for h 0 ZZ and H 0 ZZ: The Feynman diagrams for the one-loop virtual contributions are displayed in Figures In the limit of vanishing electron mass, The differential LO cross section for unpolarized particles is given by ) and the amplitudes A µ are given in [59]. The part of the amplitude following from the one-loop diagrams, M 1 , can be projected on the two invariants P L,R as The omitted terms are of the typev(1 ± γ 5 )u; they vanish in the interference with M 0 and are hence not required at one-loop order. For the counterterm part of the amplitude, δM 1 , the projection on these two invariants is exact, The NLO differential cross section can be written as follows, (spin summation to be understood) with the Born cross section from Eq. (47). At one-loop order, only the interference term contributes, on which our numerical analysis is based. Nevertheless, we will comment also on the purelyquadratic term, which may be useful to give some partial information on the size of higher-order corrections of O(α 2 ).
The integrated cross section at the one-loop level, σ 1 , derived from Eq. (50) with the interference term only, can be written in the following way, The relative correction δ in Eq.(53) can be decomposed into the following parts, indicating their origin, The electroweak terms δ vertex and δ boxes do not contain any virtual-photon diagram. According to the discussion at the end of section.3.2, all virtual-photon diagrams for vertex, box, and external wave-function contributions have been properly separated and combined with real bremsstrahlung to form the subclass of QED corrections, described by δ QED . The analytical expression for the soft photon Bremsstrahlung can be found in [55]

Numerical results and discussions
In this section, we describe in detail the numerical results that are obtained from the analysis of the processes e + e − → Z + h/H at one loop level in the SM and 2HDM. Basically, we concern with the following two quantities: the ratio of the weak corrections and the ratio of the full one loop corrections (including real emissions) to the leading order results are shown. The ratio of the weak corrections is defined as: where σ weak denotes the 1-loop weak corrections and σ (0) is the tree-level contribution of the SM and the 2HDM, respectively. The ratio of the full next-to-leading order results is defined as: where σ (1) denotes the full 1-loop corrections including real emissions The total cross section for e + e − → Zh SM at NLO including the real emissions can be found in Table 7. The result under α(0) scheme still depends on log(m e ) while the other two are independent. We can see that scale dependence has been greatly improved at NLO. And our results for weak part are consistent with those in Ref. [60].
In Fig. (6), the dependence on the renormalization scale µ and the collision energy of the cross section are shown. In Fig. (6(a)), it is noticed that the leading order results can  Figure 4: Generic self-energies one-loop 2HDM Feynman diagrams contributing to e + e − → Zh 0 (we take φ = h 0 ).
change drastically when the unphysical renormalization scale µ varies from 0 to √ s, i.e. the difference between the results of α(0) and α( √ s) can reach up to 15% as given in Table (7).
In contrast, the next-leading order results are improved a lot. For example, the difference between the results of α(0) and α( √ s) can only reach 1% or so. In Fig. (6(b)), the lineshapes of the cross section at the leading order varying with the collision energy √ s are shown. The largest difference between different renormalization scales occur near the threshold. The larger the collision energy √ s goes, the less the difference betweem the results of different renormalization scales. In Fig. (7), we compare the results of leading order and the full next leading order. In Fig. (7(a)), the distribution of the transverse momentum of Higgs boson in the process e + e − → Zh is shown. It is noticed that the lineshapes of the leading and full next leading order results are similar except a global scaling factor. In Fig. (7(b)), the cross section with the vary of collision energy is shown. It is noticed that when √ s < 700 GeV, the cross section of the full next leading order is typically smaller than that of the leading order, but when √ s > 700 GeV, the cross section of the full next leading order becomes larger, though the difference is not dramatic.  Due to the constraints from the conditions of vacuum stability and the unitarity conditions of perturbation to the theoretical parameters, it is noticed that the parameter λ 5 is allowed to vary in the range from -1 to 1 or so in type-I, and from 1 to 9 or so in type-II at the process e + e − → Zh, respectively. For the same reasons, the parameter λ 5 is allowed to vary in the range from -0.5 to 0.3 or so in type-I and from -0.6 to 0.3 or so in type-II at the process e + e − → ZH. The parameter of λ 5 is related to the triple Higgs boson couplings, which is given as [33] Once λ 5 is fixed, then triple Higgs boson couplings are fixed, as shown in Eq. (7). To    demonstrate the effects of triple Higgs boson self-couplings, we take a few typical values of λ 5 = −0.5, λ 5 = 0 and λ 5 = 0.3 in Fig. (8(a)) and λ 5 = 1, λ 5 = 2 and λ 5 = 5 Fig. (8(b)), respectively. In Fig.(8), results of the ratios in Eq. (55) and Eq. (56) are shown, where the result of the SM at the leading order is computed with the renormalization scheme α(M z ) and the h = h SM . The NLO results of SM include both the weak correction and the real emission. It is noticed that signs of the real emission and the weak interaction are different in the SM. The real emission tends to increase the total cross section by a factor +0.5%. From the Fig.  (8(b)), it is noticed that all the lines have a bump near the collision energy 350-400 GeV, which can be attributed to the contribution of top quark pair in the loop functions.
The contribution of new physics tends to further decrease the contribution of weak correction of the SM, as shown in Fig. (8) with three typical values of the parameter λ 5 . It it noticed from Fig. (8) that the larger the value of λ 5 , the smaller is the correction, which is due to the cancellation in the triple Higgs couplings, i.e. a larger the value of λ 5 corresponds a smaller triple Higgs boson couplings given in Eq. (7), as demonstrated in table (3). In Fig. (9), two benchmark points for an alternative interpretation of the SM-like Higgs boson are considered, in this scenario the heavy CP even Higgs boson has a mass near 125 GeV. The detailed information on the mass spectra and parameters of these two bench mark points are presented in Table (2) and labelled as BP1-H and BP2-H. We take three typical values of λ 5 = −0.5, λ 5 = 0, and λ 5 = 0.3, respectively. It it noticed that the contributions of the first benchmark for different typical values of λ 5 increase from −10% to −20%, while those of the second benchmark changed from −20% to −30% or so due to the large contribution of Higgs boson's loop.
Furthermore, it is noticed that the weak corrections of 2HDM for these two benchmark points and the typical values of λ 5 have the same sign to those of SM, which is hold for both the type-I and type-II as well as for the process e + e − → Zh and e + e − → ZH.
The last but not the least, it is noticed that all the lines have a bump near the collision energy 350-400 GeV, which can be attributed to the contribution of top quark pair in the loop functions.
In Fig. (10) and (11) we present the contribution of new physics ∆ EW Zh and ∆ W EAK ZH , respectively for different values λ 5 . The quantity ∆ EW is defined as given below which is to describe the contribution of new physics compared with the results of the SM but the real emission is not taken into account. Obviously, the pure contribution of the SM is subtracted in this quantity, but the interferences between the new physics and the SM are counted. In Fig. (10), the benchmark points BP1-h and BP2-h are considered. In Fig. (10(a)), there is a bump near the 380 GeV or so, which is due to the contribution of charged Higgs bosons. In Fig. (10(b)), there is a little bump near the 350 GeV, which is due to the contribution of interference of the SM and the new physics. For BP2-h, the contribution of Higgs sector is almost constant and does not change with collision energy.
In Fig. (11), the bench mark points BP1-H and BP2-H are considered. In Fig. (11(a)), there is a bump near the 270 GeV or so, which is due to the contribution of charged Higgs bosons. In Fig. (11(b)), there are a little bump near the 350 GeV, which is due to the contribution of interference of the SM and the new physics. The change of λ 5 can change ∆ EW a few percent. Furthermore, there is a little dip near 730 GeV or so, which is due to the contribution of H − A via the diagram s 12 and s 19 in Fig. (4). For BP2-H, the contribution of Higgs sector is almost constant and does not change with collision energy.
In Figure.(12), we discuss numerically the correlations of ∆ bb (h 0 ) and ∆ τ τ (h 0 ) in four types of 2HDMs for λ 5 and tan β in the ranges [−2, 2], [1,5] respectively, the results are given for cos(β − α) ≤ 0 (blue point) and ≥ 0 (black points) are shown. Our results are in good agreements with [10,11]. Red and yellow ellipse indicate the expected prediction at ±1σ at HL-LHC and CEPC, respectively. In each plot, the point (0,0) corresponds to the SM. Predictions of all the types get close to the expected precision except for 2HDMtype-II. It turns out that if cos(β − α) change the sign, ∆ bb (h 0 ) and ∆ τ τ (h 0 ) for each type lead to deviate in different directions. Therefore we can detemine the sign of cos(β − α) by using a precise measurements of bb and τ τ . However it difficult to discriminate the types of 2HDMs by evaluating on h 0 → bb and h 0 → ττ because behaviors of ∆ bb (h 0 ) and ∆ τ τ (h 0 ) depend on the sign of cos(β − α). At HL-LHC, h 0 τ τ and h 0 bb couplings are expected to be measured with about 8% and 11%, respectively dependung on tan β. In that case, we can discriminate the types of Yukawa interactions by using HL-LHC. At the CEPC, however, the Higgs couplins measurements have typically O(1)% level resolution: e.g h 0 couplings to b and τ can determine in type-III and IV, respectively. In Figure.(13), we focus on the heavy Higgs H is SM-like and whether we can distinguish it from h 0 SM-like case using the Higgs couplings measurements. Areas of the 1σ accuracy from the fiducial points are shown. In this Figure, it is clear that by using precision measurements of the Higgs boson couplingss, we can distinguish H-SM-like from h-SM-like in all types of 2HDM. We note that the fiducial points in all types have a large radiative corrections which may be excluded by the sensitivity for Hbb and Hτ τ because the ratios of decay rate for the fermion are differents, therefore we may be able to distinguish the H-SM-like case by the very precise measurement of the Higgs boson couplings if the integratred luminosity is large enough.

Discussions and Conclusion
We proposed 4 benchmark scenarios of the 2HDM after taking into account the current Higgs data from the LHC. And we have evaluated the radiative corrections to the process e + e − → Zφ in the SM and in these 4 benchmark scenarios up to one-loop level. It is noticed that the calculation given by Whizard the initial state radiation [69] can cause the cross section of e + e − → Zh in the SM decreasing by a factor 10%, where the contribution of high order logarithms are resumed. While here we noticed that the real emission from the initial state can increase the cross section by a factor +0.5%. The difference can be attributed to the fact that our results are the fixed order ones. Moreover, the lineshape of P t (h) given in Fig. (7(a)) is different from the Fig. (2b) given in [70], where the distribution is a normalized one.
We have used an on-shell renormalization scheme for all parameters except for wave PSfrag replacements l k l m l k l m Figure 12: Correlation between relative precisions ∆ τ τ (h 0 ) and ∆ bb (h 0 ) in 2HDM. The ellipse show the 68% confidence regions for these couplings expected from the HL-LHC and CEPC [12].
function renormalization of the Higgs doublet which has been done in the MS scheme. Using it, we also compute the decays φ → bb and φ → τ + τ − with φ = h 0 and H 0 in the four types of 2HDM given in Table (1) by including the EW corrections. We have shown that in type II and IV the electroweak radiative corrections in these two decay processes are rather small due to the heavy states A 0 , H 0 and H ± have a mass of order of 600 GeV while they could be sizable for 2HDM type I and III, as shown in Fig. (12) and Fig. (13). Our results demonstrate that e + e − colliders (epecially the Higgs factories with √ s = 250 GeV) can offer us the potential capability to distinguish between various 2HDM models by looking at these quantum effects in Higgs observables. As demonstrated in the Table (2), except performing precision measurements, linear colliders also have the potential to discover new physics. For example, it is possible to discover the light charged Higgs boson in BP1-h and BP1-H via e + e − → H + H − processes. Obviously, the energy scan of e + e − colliders can also help to detect the mass spectra of Higgs bosons and even triple couplings of Higgs bosons, as shown in Fig. (10) and Fig. (11). It is hopeful to reconstruct the sector of Higgs potential after combining all the future data of production and decay. Our benchmark points show that it also possible to explore new physics beyond the SM, even if some new particles could not be directly produced in the future experiments. PSfrag replacements l k l m l k l m Figure 13: Correlation between relative precisions ∆ τ τ (H 0 ) and ∆ bb (H 0 ) in 2HDM. The coding color is the same as in Fig.12.