Gauge-independent Renormalization of the 2-Higgs-Doublet Model

The 2-Higgs-Doublet Model (2HDM) belongs to the simplest extensions of the Standard Model (SM) Higgs sector that are in accordance with theoretical and experimental constraints. In order to be able to properly investigate the experimental Higgs data and, in the long term to distinguish between possible models beyond the SM, precise predictions for the Higgs boson observables have to be made available on the theory side. This requires the inclusion of the higher order corrections. In this work, we investigate in detail the renormalization of the 2HDM, a pre-requisite for the computation of higher order corrections. We pay particular attention to the renormalization of the mixing angles $\alpha$ and $\beta$, which diagonalize the Higgs mass matrices and which enter all Higgs observables. The implications of various renormalization schemes in next-to-leading order corrections to the sample processes $H^\pm \to W^\pm h/H$ and $H \to ZZ$ are investigated. Based on our findings, we will present a renormalization scheme that is at the same time process independent, gauge independent and numerically stable.


Introduction
The discovery of a new scalar particle by the LHC experiments ATLAS [1] and CMS [2] in 2012 and its subsequent confirmation as being the Higgs boson [3][4][5][6] marked a milestone for particle physics. At the same time, it triggered a change of paradigm. The Higgs particle, which formerly was the object of experimental searches, has itself become a tool in the search for New Physics (NP). Although the Standard Model (SM) of particle physics has been tested in previous and present experiments at the highest accuracy, there remain many open questions that cannot be answered within this model. The SM is therefore regarded as the low-energy description of some more fundamental theory that becomes effective at higher energy scales. A plethora of NP models have been discussed, among them e.g. supersymmetry (SUSY) as one of the most popular and most intensely studied Beyond the SM (BSM) extensions. Supersymmetry requires the introduction of at least two complex Higgs doublets. The Higgs sector of the Minimal Supersymmetric extension of the SM (MSSM) [7][8][9][10] is a special case of the 2-Higgs-Doublet Model (2HDM) [7,11,12] type II. While the parameters of the SUSY Higgs potential are restricted due to SUSY relations, general 2HDMs allow for much more freedom in the choice of the parameters. They are therefore an ideal framework to study the implications of an extended Higgs sector for Higgs phenomenology at the LHC. This is reflected in the experimental analyses that interpret the results in various benchmark models, among them the 2HDM. The precise investigation of the Higgs sector aims at getting insights into the nature of electroweak symmetry breaking (EWSB) and at clarifying the question whether it is based on weakly or strongly interacting dynamics. Deviations in the properties of the discovered SM-like Higgs boson are hints towards NP. In particular, the higher precision in the Higgs couplings measurements at the LHC run 2 and in the high-luminosity option allows to search indirectly for BSM effects. This becomes increasingly important in view of the null results of direct searches for NP so far. 1 The precise measurements on the experimental side, however, call for precise predictions of parameters and observables from theory. Accurate theory predictions are indispensable not only for the proper interpretation of the experimental data, but also for the correct determination of the parameter space that is still allowed in the various models, and, finally, for the distinction between different BSM extensions.
With this paper we contribute to the effort of providing precise predictions for parameters and observables relevant for the phenomenology at the LHC and future e + e − colliders. We investigate higher order corrections in the framework of the 2HDM. While 2HDMs are interesting because they contain the MSSM Higgs sector as a special case, they also belong to the simplest SM extensions respecting basic experimental and theoretical constraints that are testable at the LHC. After EWSB they feature five physical Higgs bosons, two neutral CP-even, one neutral CPodd and two charged Higgs bosons. They represent an ideal benchmark framework to investigate the various possible NP effects to be expected at the LHC in multi-Higgs boson sectors. Finally, specific 2HDM versions also allow for a Dark Matter (DM) candidate [15][16][17][18][19][20][21]. In the past, numerous papers have provided higher order corrections to the 2HDM parameters, production cross sections and decay widths. Several papers have dealt with the renormalization of the 2HDM (see e.g. [22][23][24]). In particular, the renormalization of the mixing angles α and β is of interest. While α is introduced to diagonalize the mass matrix of the neutral CP-even Higgs sector, the angle β appears in the diagonalization of the CP-odd and the charged Higgs sector, respectively. These angles define the Higgs couplings to the SM particles and thus enter all Higgs observables like e.g. production cross sections and decay widths. For the MSSM it was stated in [25] that a renormalization scheme for the only mixing angle taken as an independent parameter from the scalar sector, β, cannot be simultaneously gauge independent, process independent and numerically stable. In the 2HDM also α needs to be renormalized, which has important consequences for the choice of the renormalization scheme. If the tadpoles are treated in the usual way, which we call the standard approach (cf. 3.1.1), a process-independent definition of the angular counterterms is prone to lead to gauge-dependent amplitudes and consequently to gauge-dependent physical observables. This is the case e.g. in the scheme presented in [23]. There are essentially two possibilities to circumvent the emergence of this gauge dependence. Either one gives up the requirement of process independence and fixes α and β in terms of a physical observable or one changes the treatment of the tadpoles. As we will see, this will decouple the issue of gauge dependence from the definition of δα and δβ and allow for process-and gaugeindependent angular counterterms leading to manifestly gauge-independent amplitudes.
In this paper we study in detail the renormalization of the 2HDM Higgs sector with the main focus on the investigation of the gauge dependence of the renormalization of the mixings angles α and β. We propose several schemes and compare them both to the ones in the literature and amongst each other. In sample decay processes we investigate the numerical differences and in particular the numerical stability of the various renormalization prescriptions. Our results presented here will serve as basis for the further computation of the one-loop electroweak corrections to all 2HDM Higgs boson decays.
The organization of the paper is as follows. In sec. 2 we introduce the model and set up our notation. The following sec. 3 is devoted to the detailed presentation and discussion of the various renormalization prescriptions that will be applied. Section 4 deals with the computation of the electroweak (EW) one-loop corrections to various decay processes and the discussion of the gauge dependence of the angular counterterms. In sec. 5 we present our numerical results. We finish with the conclusions in sec. 6. The paper is accompanied by an extensive Appendix to serve as starting point for further investigations of the 2HDM renormalization.

Description of the Model
We work in the framework of a general 2HDM with a global discrete Z 2 symmetry that is softly broken. The kinetic term of the two SU (2) Higgs doublets Φ 1 and Φ 2 is given by in terms of the covariant derivative where τ a denote the Pauli matrices, W a µ and B µ the SU (2) L and U (1) Y gauge bosons, respectively, and g and g the corresponding gauge couplings. The scalar potential that can be built from the two SU (2) Higgs doublets can be written as 3) The discrete Z 2 symmetry (Φ 1 → −Φ 1 , Φ 2 → Φ 2 ) ensures the absence of tree-level Flavour Changing Neutral Currents. Assuming CP conservation, the 2HDM potential depends on eight real parameters, three mass parameters, m 11 , m 22 and m 12 , and five coupling parameters λ 1 -λ 5 . Through the term proportional to m 2 12 the discrete Z 2 symmetry is softly broken. After EWSB the neutral components of the Higgs doublets develop vacuum expectation values (VEVs), which are real in the CP-conserving case. Expanding about the VEVs v 1 and v 2 and expressing each doublet Φ i (i = 1, 2) in terms of the charged complex field φ + i and the real neutral CP-even and CP-odd fields ρ i and η i , respectively, leads to the mass matrices, which are obtained from the terms bilinear in the Higgs fields in the potential. Due to charge and CP conservation they decompose into 2 × 2 matrices M S , M P and M C for the neutral CP-even, neutral CP-odd and charged Higgs sector. They are diagonalized by the following orthogonal transformations leading to the physical Higgs states, a neutral light CP-even, h, a neutral heavy CP-even, H, a neutral CP-odd, A, and two charged Higgs bosons, H ± . The massless pseudo-Nambu-Goldstone bosons G ± and G 0 are absorbed by the longitudinal components of the massive gauge bosons, the charged W ± and the Z boson, respectively. The rotation matrices in terms of the mixing angles ϑ = α and β, respectively, read R(ϑ) = cos ϑ − sin ϑ sin ϑ cos ϑ . (2.8) The mixing angle β is related to the two VEVs as with v 2 1 + v 2 2 = v 2 ≈ (246 GeV) 2 , while the mixing angle α is expressed through where (M S ) ij (i, j = 1, 2) denote the matrix elements of the neutral CP-even scalar mass matrix M S . With we have [23] tan 2α = where we have introduced the abbreviation and used short-hand notation s x ≡ sin x etc.
The minimization conditions of the Higgs potential require the terms linear in the Higgs fields to vanish in the vacuum, i.e.
where the brackets denote the vacuum. The corresponding coefficients, the tadpole parameters T 1 and T 2 , therefore have to be zero. The tadpole conditions at lowest order are given by There are various possibilities to choose the set of independent parameters that parametrizes the Higgs potential V . Thus, Eqs. (2.15) and (2.16) can be used to replace m 2 11 and m 2 22 by the tadpole parameters T 1 and T 2 . The VEV v can furthermore be expressed in terms of the physical gauge boson masses M W and M Z and the electric charge e. In the following, we will choose the set of independent parameters such that the parameters can be related to as many physical quantities as possible. Our set is given by the Higgs boson masses, the tadpole parameters, the two mixing angles, the soft breaking parameter, the massive gauge boson masses and the electric charge. Additionally, we will need the fermion masses m f for the Higgs decays into fermions which will be used for a process-dependent definition of the angular counterterms.

Renormalization
In this section we will present the various renormalization schemes that we will apply in the renormalization of the 2HDM and that will be investigated with respect to their gauge parameter dependence and their numerical stability. We will use these schemes in sample processes given by the EW one-loop corrected decays of the charged Higgs boson into a W ± and a CP-even Higgs boson, H ± → W ± h/H, and of the heavy H into a Z boson pair, H → ZZ. The computation of the EW one-loop corrections leads to ultraviolet (UV) divergences. In the charged Higgs decay we will furthermore encounter infrared (IR) divergences because of massless photons running in the loops. The UV divergences in the virtual corrections are canceled by the renormalization of the parameters involved in the EW corrections of the process, while the IR ones are subtracted by taking into account the real corrections. The renormalization of the above decay processes requires the renormalization of the electroweak sector and of the Higgs sector. We will also compute the EW one-loop corrections to the decays of H and A into τ leptons, H/A → τ τ . These processes will be exploited for a process-dependent definition of the angular counterterms, which will be presented as a possible renormalization scheme among others. The corrections to the decays into τ leptons also require the renormalization of the fermion sector. Note, that the renormalization of the CKM matrix, which we will assume to be real, will not play a role in our renormalization procedure. We start by replacing the relevant parameters by the renormalized ones and their corresponding counterterms: Gauge sector: The massive gauge boson masses and the electric charge are replaced by Equally, the VEV v, which will be expressed in terms of these parameters, is replaced by The gauge boson fields are renormalized by their field renormalization constants δZ, Fermion sector: The counterterms to the fermion masses m f are defined through The bare left-and right-handed fermion fields are replaced by their corresponding renormalized fields according to Higgs sector: The renormalization is performed in the mass basis and the mass counterterms are defined through and the mixing angles by While the tadpoles vanish at leading order, the terms linear in the Higgs fields get loop contributions at higher orders. Therefore, also the tadpole parameters T 1 and T 2 have to be renormalized in order to fulfill the tadpole conditions Eqs. (2.14). The tadpoles are hence replaced as

Renormalization conditions
The finite parts of the counterterms are fixed by the renormalization conditions. Throughout we will fix the renormalization constants for the masses and fields through on-shell (OS) conditions. The renormalization schemes differ, however, in the treatment of the tadpoles and of the mixing angles. We will describe two different approaches for the treatment of the tadpoles. Both of them apply the same renormalization conditions for the tadpoles. They differ, however, in the way the minimum conditions are applied when the mass counterterms are generated. As a consequence, the tadpole counterterms can either explicitly show up in the mass counterterms or not. The latter case, that we will call 'alternative tadpole' or in short 'tadpole' scheme, has the virtue that the mass counterterms are manifestly gauge independent, while in the former one, named 'standard tadpole' or simply 'standard' scheme, this is not the case. The authors of Ref. [23] have combined the standard tadpole scheme with the definition of the counterterms through off-diagonal wave function renormalization constants. This 'KOSY' scheme, denoted by the initials of the authors, leads to manifestly gauge-dependent decay amplitudes, as we will show. In the alternative tadpole scheme not only this problem does not occur, but in addition, the angular counterterms are explicitly gauge independent. If the angular counterterms are defined in a 'process-dependent' scheme via a physical process, the decay amplitude is gauge independent irrespective of the treatment of the tadpoles. The only difference lies in the gauge independence of the angular counterterms in case the alternative tadpole scheme is adopted. In the following, the renormalization conditions of the various schemes will be introduced.

Standard Tadpole Scheme
We start by presenting the usual, i.e. 'standard', approach in the renormalization of the 2HDM as also applied in [23,24]. The gauge bosons are renormalized through OS conditions, which implies the following counterterms for the masses, where the superscript T denotes the transverse part of the respective self-energy Σ. In order to guarantee the correct OS properties the wave function renormalization constants have to be introduced as The electric charge is defined to be the full electron-positron photon coupling for OS external particles in the Thomson limit, implying that all corrections to this vertex vanish OS and for zero momentum transfer. The counterterm for the electric charge in terms of the transverse photon-photon and photon-Z self-energies reads [26] δZ α(0) where s W /c W ≡ sin θ W / cos θ W and θ W denotes the Weinberg angle. Note that the sign in the second term of Eq. (3.23) differs from the one in [26] due to our sign conventions in the covariant derivative of Eq. (2.2). In our computation, however, we will use the fine structure constant at the Z boson mass α(M 2 Z ) as input, so that the results are independent of large logarithms due to light fermions f = t. The counterterm δZ e is therefore modified as [26] δZ where the transverse part of the photon self-energy Σ T γγ in Eq. (3.25) includes only the light fermion contributions. For the computation of the EW one-loop corrected Higgs decay widths we also need to renormalize the coupling g, which can be related to e and the gauge boson masses as so that its counterterm can be expressed in terms of the gauge boson mass counterterms through Defining the following structure for the fermion self-energies the fermion mass counterterms applying OS conditions are given by The fermion wave function renormalization constants are determined from The OS conditions for the physical Higgs bosons yield the following Higgs mass counterterms The appearance of the tadpole counterterms in Eqs. (3.31) and (3.32) can be understood by recalling that the parameters m 2 11 and m 2 22 , which enter the mass matrices, can be replaced by the tadpole coefficients T 1 and T 2 . Applying the shifts Eq. (3.19) and rotating into the mass eigenbasis yield the above conditions in the OS scheme. The relations between the tadpole counterterms in the mass basis and δT 1,2 are given by The tadpoles are renormalized such that the correct vacuum is reproduced at one-loop order, leading to the renormalization conditions The T 1,2 stand for the contributions coming from the corresponding genuine Higgs boson tadpole graphs in the gauge basis. For the wave function renormalization constants the OS renormalization implies the following conditions (3.39)

The KOSY Scheme
We now turn to the renormalization conditions for the mixing angles. The renormalization scheme chosen in [23], the 'KOSY' scheme, uses the standard tadpole scheme. For the renormalization of the mixing angles it is based on the idea of making the counterterms δα and δβ appear in the inverse propagator matrix and hence in the wave function renormalization constants, in a way that is consistent with the internal relations of the 2HDM. This can be achieved by renormalizing in the mass basis (f 1 , f 2 ) T , but temporarily switching to the gauge basis (γ 1 , γ 2 ) T , and back again, The fields f i and γ i (i = 1, 2) and the mixing angle ϑ stand here for any of the field pairs in the mass and gauge basis, respectively, defined in Eqs. (2.5)-(2.7), together with their corresponding mixing angle, i.e. (f i ; γ i ; ϑ) = (H, h; ρ i ; α), (G 0 , A; η i ; β) and (G ± , H ± ; φ ± i ; β). With the field renormalization matrix Z γ in the gauge basis being a real symmetric matrix the following parametrization of the field renormalization matrices in the mass basis can be chosen [23,24] The off-diagonal elements are identified with the off-diagonal wave function renormalization constants in the mass basis. For the CP-even scalar sector we obtain The superscript 'OS' indicates the OS renormalization scheme for the wave function constants. The counterterm δC h will not be used again. While the mixing angle β diagonalizes both the charged and the CP-odd mass matrices and we have altogether four off-diagonal wave function constants in the charged and CP-odd Higgs sector, Eq. (3.41) implies only three free parameters to be fixed, namely δβ, δC A and δC H ± . Consequently, one has to choose three out of four possible conditions and not all scalar fields can be OS at the same time. If we choose e.g. the OS renormalized δZ OS G 0 A , δZ OS G ± H ± and δZ OS H ± G ± to fix the counterterms, we ensure H ± to be OS. This scheme can hence be used in the process H ± → W ± h/H, where we have an external charged Higgs boson. 2 This yields the following possible first set of counterterms, Choosing on the other hand the set δZ OS G 0 A , δZ OS AG 0 and δZ OS H ± G ± we get a second possible set There are two more sets that can be chosen. However, we are not going to use them and hence they will not be repeated here. Replacing the OS conditions given in Eqs. (3.37), (3.38) and (3.39) in Eqs. (3.44), (3.46) and (3.49), respectively, yields the following counterterms for the mixing angles α and β (3.54) As already mentioned before and as we will demonstrate later in detail for the example of the charged Higgs boson decay, the application of this renormalization scheme not only makes a gauge-independent definition of the counterterms impossible, but more seriously, leads to unphysical gauge-dependent decay amplitudes. The computation of the loop-corrected amplitude in the general R ξ gauge shows that after including all counterterms but the ones for the angles, there remains a residual gauge dependence that is UV-divergent. The angular counterterms must therefore reveal exactly the same UV-divergent gauge dependence but with opposite sign. The counterterm δα is found to have exactly this UV-divergent ξ-dependent counterpart, needed to render the amplitude gauge independent. However, in addition, δα and δβ contain ξ-dependent finite terms, which reintroduce a gauge dependence into the amplitude. To get rid of these finite gauge-dependent terms in δβ, the authors of Ref. [24] suggest to drop the assumption that Z f is symmetric, thereby yielding additional renormalization conditions. These are then exploited to move the gauge dependence of δβ into δC f 3 . While this scheme would in principle allow to eliminate the gauge dependence of δβ, it cannot be applied in processes that involve the renormalization of α. The UV-divergent ξ-dependent counterterm δα is needed to cancel the UV-divergent ξ-dependent counterpart in the loop-corrected amplitude, that is encountered in the standard renormalization scheme. In practice, however, this procedure cannot be applied, as it lacks an unambiguous prescription on how to extract the truly gauge-independent parts from the loop-corrected amplitude and from the counterterms. The extraction of the gaugeindependent part is not straightforward as the loop functions A 0 and B 0 [27,28] which appear in the angular counterterms, can be rewritten in terms of higher n-point scalar integrals that contain the gauge parameter ξ besides additional gauge-independent components.

Alternative Tadpole Scheme
We now present a renormalization scheme that fulfills the requirements for a possible gaugeindependent definition of the angular counterterms. It relies on the application of the renormalization scheme worked out in Ref. [29]. In Appendix A we show in detail how this scheme works and in particular we present its extension from the SM case [29] to the 2HDM. The generic diagrams contributing to the self-energies defined in this 'alternative tadpole' scheme, called Σ tad in the following, are shown in Fig. 1. Besides the generic one-particle irreducible (1PI) diagrams depicted by the first two topologies in Fig. 1, they also contain the tadpole diagrams connected to the self-energies through the CP-even Higgs bosons h and H that are represented by the third topology. The application of the tadpole scheme alters the structure of the mass counterterms and of the off-diagonal wave function renormalization constants 4 such that now the loop-corrected amplitude including all counterterms but those for the angles does not encounter a UV-divergent ξ dependence any more. Hence, also the angular counterterms can and even have to be defined in a gauge-independent way by applying appropriate renormalization conditions. h/H Besides the angular counterterms, also the mass counterterms, defined via OS conditions become gauge independent in the tadpole scheme. This has been shown for the electroweak sector in [30]. All counterterms of the electroweak sector have exactly the same structure as in the standard scheme, but the self-energies Σ appearing in Eqs. (3.20)-(3.23) have to be replaced by the self-energies Σ tad containing the tadpole contributions. Note however, that there are no tadpole contributions for Σ T γZ so that Furthermore, due to the fact that the tadpoles are independent of the external momentum the derivatives of the self-energies do not change, The Higgs mass counterterms become And for the Higgs wave function renormalization constants we obtain keeping in mind that Eq. (3.56) holds. Applying the same procedure for the definition of the angular counterterms as in the standard scheme, but with the different treatment of the tadpoles, the angular counterterms in the tadpole scheme read (3.64) Compared to the standard scheme, the self-energies are replaced by the Σ tad and no tadpole counterterms appear any more.
The application of the tadpole scheme not only allows for a gauge-independent definition of the angular counterterms but also requires it in order to ensure a gauge-independent physical decay amplitude. Note that the counterterms (3.62)-(3.64) still contain a ξ dependence and hence, a ξ-independent definition has yet to be found. In the MSSM, several schemes for the renormalization of tan β have been proposed and used, see e.g. [25,[31][32][33][34][35][36][37][38]. The renormalization prescriptions have been discussed in detail in [25] with respect to their gauge dependence, process independence and numerical stability (see also [39]). Renormalization prescriptions making use of physical quantities like Higgs boson masses or physical processes clearly lead to a gaugeindependent prescription. However, they were found to be numerically unstable in the former case, while the latter case may be viewed as unsatisfactory, as the definition via a specific process makes tan β a non-universal, flavour-dependent quantity [25]. Finally, DR prescriptions lead in the R ξ gauge to gauge-independence of δ tan β in the MSSM at one-loop level, but not at two-loop level [25,40]. We now present a renormalization scheme that leads to ξ-independent δα and δβ and also addresses the problem of extracting the gauge-independent part unambiguously.

On-shell tadpole-pinched scheme
The scheme we propose here combines the virtues of the tadpole scheme with the unambiguous extraction of the truly gauge-independent parts of the angular counterterms. It is based on the renormalization schemes presented in [38] and in [41,42]. The former defines the angular counterterms in a physical way as residues of poles appearing in one-loop corrections, while in [41,42] the pinch technique 5 [52][53][54][55][56] is used to extract the truly gauge-independent parts of the angular counterterms. Both methods lead to the same gauge-independent definitions of the counterterms. With the help of the PT it is possible to define the pinched self-energies Σ. The self-energies are related to the tadpole self-energies evaluated in the Feynman gauge as where ξ stands for the gauge fixing parameters ξ Z , ξ W and ξ γ of the R ξ gauge. Note, that in order to apply the PT the tadpole scheme has to be used. 6 For better readability we omitted the superscript 'tad' in Σ. The self-energy Σ add in Eq. (3.65) is an additional contribution that is explicitly independent of the gauge fixing parameter ξ. Applying [42] we arrive at the following counterterms (3.68) These angular counterterms are different from the ones obtained in the KOSY scheme, so that the classification as an independent renormalization scheme is justified. The additional contribution Σ add Hh has been given in [42] for the MSSM. We have derived the remaining two contributions Σ add G 0 A and Σ add G ± H ± . Altogether we have where B 0 is the scalar two-point function [27,28].
p tadpole-pinched scheme As indicated by the name, this scheme differs from the OS tadpole-pinched scheme solely in the scale at which the self-energies, appearing in the definition of the angular counterterms, are evaluated. The self-energies are evaluated at the average of the particle momenta squared [42], , respectively, and we will henceforth refer to this scheme as the p -scheme. When the self-energies are evaluated at p 2 the additional self-energies Σ add vanish, as can easily be seen from Eqs. (3.69)-(3.71), and the pinched self-energies are given by the tadpole self-energies Σ tad computed in the Feynman gauge, i.e.
The angular counterterms then read (3.76)

Process-dependent Scheme
We will also investigate the renormalization of the mixing angles through a physical process. Provided the alternative tadpole scheme is applied, this leads to a manifestly gauge-independent definition of the mixing angle counterterms. In order to fix the respective angular counterterm we will require the next-to-leading order (NLO) Higgs decay width, in which the angle appears, to be equal to the leading order (LO) one, i.e.
where Γ virt denotes the contribution of all virtual one-loop corrections to the decay width and Γ c.t. the counterterm contributions. This implies (for NLO processes that do not encounter real corrections, see below) and allows to fix the angular counterterm via the decay process. This scheme has some drawbacks, however, cf. [25]. Conceptually, it is not satisfying as the definition of the mixing angles becomes non-universal and flavour-dependent. From a calculational point of view, it is involved as it requires the computation of loop-corrected three-particle vertices. Another problem is related to the choice of the process that defines the counterterm. The definition through a process receiving QED corrections that cannot be separated from the rest of the EW corrections would entail real radiative corrections in the counterterm. This is precluded, however, as this counterterm would inevitably depend on some detector sensitivity ∆E via the photon phase space cut and thereby introduce a dependence on the experimental setting. This forbids e.g. the definition of the angular counterterms appearing in the loop corrected decay H ± → W ± h through the process H ± → W ± H. Finally, care has to be taken to choose a process that is phenomenologically accessible. This eliminates e.g. the choice of H → ZZ. With the 125 GeV Higgs boson being very SM-like and hence coupling with full SM strength to the Z bosons, sum rules lead to a tiny coupling of the heavy Higgs boson to massive gauge bosons and hence a very small H → ZZ decay width. In this paper we choose, as proposed in [25], the decays H → τ τ and A → τ τ in order to define δβ via the latter and δα via the former. In both decays the QED corrections form a UV-finite subset of the full EW one-loop corrections.

One-Loop EW Corrected Decay Widths
In this section we present the EW one-loop corrections to the processes 7 The charged Higgs decays (4.1) will serve us to discuss in detail the renormalization of the mixing angles α and β in view of a gauge-independent definition. In this context, the fermionic decays (4.3) will be used for a process-dependent definition of the angular counterterms. Note that we could have equally well chosen h → τ τ instead of H → τ τ . The numerical implications of the different renormalization schemes shall be investigated in the subsequent section. This will be done not only for the charged Higgs decays, but also for another sample process, the heavy Higgs decay into a Z boson pair (4.2).

Electroweak One-Loop Corrections to H ± → W ± h/H
The decays of the charged Higgs boson into the charged W ± boson and a CP-even Higgs boson φ = h or H, depend through the couplings on the mixing angle combinations 5) and the LO decay width is given by The NLO decay width can be written as The one-loop correction Γ (1) consists of the virtual corrections, the counterterm contributions and the real corrections. The counterterms cancel the UV divergences and the real corrections the IR divergences encountered in the virtual corrections. The diagrams contributing to the latter are depicted in Fig. 2 and show the pure vertex corrections (a) and the corrections (b)-(e) to the external legs. The counterterm diagram is shown in (f). The vertex corrections comprise the 1PI diagrams given by the triangle diagrams with scalars, fermions and gauge bosons in the loops, as shown in the first two rows of Fig. 3, and the diagrams involving four-particle vertices (last four diagrams of Fig. 3). The corrections to the external legs in Fig. 2 (b) and (c) vanish due to the OS renormalization of the scalars, while the vanishing of the mixing contribution (d) is ensured by a Slavnov-Taylor identity [59] 8 and the one of (e) by the Ward identity for an OS W ± boson. The vertex contributions with a photon in the loop lead to IR divergences that need to be canceled by the real corrections. These are computed from the diagrams displayed in Fig. 4. They consist of the proper bremsstrahlung contributions (a)-(c), where a photon is radiated from the charged initial and final state particles, and the diagram (d) involving a fourparticle vertex with a photon. Note, that this last diagram leads to an IR-finite contribution.
(a)   Figure 3: Generic diagrams contributing to the vertex corrections in H ± → W ± h/H.
The NLO contributions factorize from the LO amplitude, so that the one-loop corrected decay width can be cast into the form type I II lepton-specific flipped The counterterm contribution ∆ ct is given in terms of the wave function renormalization constants, the coupling and angle counterterms. For φ ≡ h it reads and for φ ≡ H, As the expressions for the counterterm ∆ ct and the virtual and real contributions ∆ virt and ∆ real in terms of scalar one-, two-and three-point functions are rather lengthy, we do not display them explicitly here.

Electroweak One-Loop Corrections to H → τ τ and A → τ τ
The LO decay width for the process H → τ τ reads with the coupling modification factor g Hτ τ in the 2HDM, which depends on the 2HDM type. We give in Table 1 the coupling factors for all neutral Higgs bosons to τ leptons in the different realizations of the 2HDM. For the decay A → τ τ the LO decay width is with g Aτ τ given in Table 1. These two processes can hence be used to define the counterterms for α and β.
The EW NLO corrections to H → τ τ consist of the virtual corrections, the counterterms and the real corrections. The generic contributions to the virtual corrections are depicted in Fig. 5. The 1PI diagrams of the vertex corrections are shown in Fig. 6 and consist of the triangle diagrams with scalars, fermions, massive gauge bosons and photons in the loop. The corrections to the external legs in Fig. 5 (b), (d) and (e) vanish because of the OS renormalized H and τ , respectively. Diagram (c) is zero because of CP conservation. Diagram (f) finally vanishes because of a Slavnov-Taylor identity. The real corrections consist of the diagrams where a photon is radiated off either of the final state τ leptons. We explicitly checked that all NLO corrections factorize from the LO width so that the NLO decay width can be cast into the form (4.14)  Figure 6: Generic diagrams contributing to the vertex corrections in H → τ τ .
For ∆ ct we have ∆ ct = δZ HH + g hτ τ g Hτ τ δZ hH + δZ L τ τ + δZ R τ τ + 2 Note, that the pure QED contributions in ∆ virt and ∆ ct can be separated from the weak contributions in a gauge-invariant way and form a UV-finite subset by themselves. This is important as it allows to define the angular counterterm via this process through the purely weak NLO contributions, see also the discussion in section 3.1.4. Requiring the following renormalization condition for the process-dependent definition of δα, 16) and imposing this condition only on the weak part of the decay width we arrive at the processdependent counterterm definition δα H→τ τ = − g Hτ τ 2g hτ τ δZ HH + g hτ τ g Hτ τ δZ hH + δZ L,weak τ τ + δZ R,weak τ τ H→τ τ . (4.17) The superscript 'weak' indicates that in the respective counterterms and in the virtual correction only the purely weak contributions are taken into account. For example for ∆ virt,weak H→τ τ this means that corrections stemming from diagrams in Fig. 6 that involve photons are dropped.
The counterterm δ tan β or δβ, respectively, which is necessary in (4.17), can be defined in a process-dependent scheme via the NLO decay A → τ τ as outlined in the following. Again the NLO contributions consist of virtual, counterterm and real diagrams. The generic ones for the former two are shown in Fig. 7 and the 1PI diagrams of the vertex corrections are summarized in Fig. 8. The loops contain scalars, fermions, massive gauge bosons and photons. The loops with photons induce IR divergences that are canceled by the real corrections. The corrections to the external legs in Fig. 7 (b), (d) and (e) vanish due to OS renormalization conditions, those in (c) because of CP invariance and those in (f) because of a Slavnov-Taylor identity. Also in this process the pure QED corrections can be separated from the remainder in a gauge-invariant way  and form a UV-finite subset so that the NLO decay width can be used for the process-dependent definition of the counterterm δβ through the requirement With the factorization Γ NLO (A → τ τ ) = Γ LO 1 + ∆ virt + ∆ ct + ∆ real (4.19) and the counterterm we arrive by imposing the condition (4.18) at Again the superscript 'weak' denotes the purely weak contributions to the respective counterterms and to the virtual corrections. Thus, ∆ virt,weak A→τ τ is given by the purely weak virtual corrections to A → τ τ at NLO which are computed from the diagrams in Fig. 8 discarding those with photons in the loop.

The gauge (in)dependence of the angular counterterms
The question of gauge dependence in the standard scheme: In order to investigate the question whether the angular counterterms can be defined in a gauge-independent way, we have calculated the one-loop corrected decay width for the charged Higgs decays in the general R ξ gauge. When we apply the standard scheme, the computation of the NLO amplitude M H ± →W ± h including all counterterms but the one for the anglesi.e. δc β−α is set to zero -yields an amplitude that depends on the gauge parameters as follows, where we have introduced the abbreviation (V ≡ W, Z) in terms of the scalar one-point function A 0 [27,28]. With p 1 we denote the incoming fourmomentum of H ± and with * (p 3 ) the polarization vector of the outgoing W ± boson with four-momentum p 3 and Note that α V is UV-divergent. This result shows explicitly what we have already stated before: In the standard renormalization scheme, the NLO decay amplitude without the angular counterterms has a residual UV-divergent gauge dependence. This can only be canceled by the angular counterterms. Therefore, the counterterms cannot be defined in a gauge-independent way. This gauge dependence is independent of the renormalization scheme chosen for the angular counterterms. It is purely due to the treatment of the tadpoles.
in terms of the scalar two-point function B 0 , we find the following gauge-dependent results for the angular counterterms, and (4.27) Here the symbol ξ=1 represents the counterterm result obtained for ξ = ξ W = ξ Z = 1. The result for δβ (2) looks similar with the appropriate mass replacements and ξ W → ξ Z . The second line in Eq. (4.26) has the appropriate structure to cancel the remaining UV-divergent gauge dependence in the amplitude (4.22). However, the additional finite terms in (4.26) and (4.27) proportional to the β-integrals defined above, reintroduce a gauge dependence into the amplitude. In [24] it was argued that the gauge dependence of δβ can be moved into the unphysical counterterm δC f , see Eq. (3.41). Yet, lacking a method to define uniquely the gauge-dependent parts in the standard scheme, where the PT cannot be applied, it remains unclear, how this could be accomplished. The situation is even worse for δα, where we necessarily have to retain the gauge-dependent part proportional to the UV-divergent A 0 functions, but must move the rest into δC f . To summarize, this result shows that not only is it impossible to arrive at a gauge-independent definition of δα in the standard scheme, but it also explicitly demonstrates that the KOSY scheme leads to an unphysical gauge dependence of the decay amplitude, which one cannot be disposed of in a straightforward way. This is not only true for the charged Higgs bosons decays we are discussing. In fact, the investigation of the origin of this gauge dependence shows, that the standard tadpole scheme inevitably leads to gauge-dependent decay widths in case the KOSY scheme is applied for the mixing angles.
If we define the angular counterterms via a physical process, however, namely through the decay widths H → τ τ and A → τ τ , compute the contribution of the counterterm δc β−α , and extract the ξ-dependent parts we obtain the following, It is exactly the same as Eq. (4.22) but with opposite sign, so that altogether the EW oneloop corrected decay width is gauge independent and UV-finite as required. The standard treatment of the tadpoles combined with a process-dependent definition hence leads to a gaugeindependent physical result, as it should be. The counterterms, however, necessarily contain a gauge dependence.
Gauge-independent angular counterterms: For the angular counterterms to be gauge-independent the loop-corrected amplitude including all counterterms but the angular ones must be independent of ξ. This can be achieved by treating the tadpoles according to Ref. [29], cf. the discussion in section 3.1.3. It means that in the counterterms Eq. (4.10) and Eq. (4.11), respectively, the self-energies Σ and the tadpole counterterms δT , contained in the wave function constants, the scalar mass counterterms and the angular counterterms, have to be replaced by Σ tad and δT = 0. Note, that the change to this tadpole scheme in principle implies new vertices arising from constant tadpole contributions to the respective original vertices, cf. App. A. In the 2HDM, however, there is no quartic vertex between two scalars, a charged Higgs and a charged gauge boson, h/H − h/H − H ± − W ∓ , where one of the external h/H legs would carry the additional tadpole contribution. Therefore, the process H ± → W ± h/H does not receive additional tadpole diagrams. The counterterms δα, δβ, δZ hH , δZ Hh δZ G 0 A and δZ G ± H ± change however. With these modifications the gauge-dependent part of the amplitude with the angular counterterms set to zero, becomes The amplitude without the mixing angle counterterm is itself gauge-independent, so that it is possible to provide a gauge-independent renormalization of the angular counterterms.
a) Gauge-independent tadpole-pinched scheme: The pinch technique allows to extract from the Green functions the truly gauge-independent part. Combined with the tadpole scheme this leads to manifestly gauge-independent angular counterterms. Choosing the OS scale, they are given by Eqs. (3.66)-(3.71). In the p scheme the formulae simplify to (3.74)-(3.76). In the numerical analysis we will apply both choices.
b) Gauge-independent process-dependent definition of the angular counterterms: Another possibility to arrive at a truly gauge-independent definition of the angular counterterms is the definition via the physical processes H/A → τ τ , provided of course that the framework of the tadpole scheme is applied.
In the processes H/A → τ τ no new diagrams are introduced when switching to the tadpole scheme, while the counterterms do change. In the tadpole scheme the process-dependent definition of δα and δβ through the requirement Eq. (4.16) and Eq. (4.18), respectively, then indeed leads to gauge independence of both counterterms and hence also of δc β−α , i.e.
(δc β−α ) tad, proc-dep ξ = 0 . (4.30) We have seen in Eq. (4.28) that the treatment of the tadpoles in the standard scheme cannot lead to gauge-independent angular counterterms, although they are defined through a physical process. In detail, this gauge parameter dependence stems from δα, whereas δβ is gaugeindependent in the process-dependent definition also without applying the tadpole scheme. Thus we have This result shows two important things: First, the process-dependent definition of the angular counterterms leads to gauge-independent counterterms only if the tadpole scheme is applied. Second, Eqs. (4.31)-(4.33) demonstrate, that in a process-dependent definition of the counterterms the difference between the application of the tadpole and the standard scheme is a gauge-dependent expression that solely depends on A 0 functions, which are UV-divergent. As the 2HDM is renormalizable this implies that also in the amplitude the difference in the application of the two schemes must be UV-divergent and must have the same structure, since the divergences have to cancel. In conclusion, this means: The definition of the angular counterterms via any physical process leads for any NLO decay process to a gauge-independent result, independently of the treatment of the tadpoles.
In the following numerical analysis in section 5 we will apply all three types of renormalization schemes, the standard, the tadpole-pinched and the process-dependent scheme, and compare them to each other. We will do this for the sample processes H ± → W ± h/H and H → ZZ. In order to describe also for this latter process the implications of the tadpole scheme, required for a gauge-independent definition of the angular counterterms, we briefly repeat the ingredients of the EW one-loop corrections to H → ZZ.

Electroweak One-Loop Corrections to H → ZZ
The LO decay width for the process is given by and depends on the mixing angles through the coupling factor g HZZ = c β−α .  diagrams for the virtual corrections and the counterterm are depicted in Fig. 9. The 1PI diagrams contributing to the vertex corrections are given by the triangle diagrams with scalars, fermions, massive gauge bosons and ghost particles in the loops, as shown in the first three rows of Fig. 10, and by the diagrams involving four-particle vertices (last four diagrams of Fig. 10). The corrections to the external leg in Fig. 9 (b) vanish due to the OS renormalization of the H. The mixing contributions (c) and (d) vanish because of the Ward identity for the OS Z boson. The counterterm amplitude is given by where the µ * denote the polarization vectors of the outgoing Z bosons with four-momentum p 3 and p 4 , respectively. If the tadpole scheme is applied, the HZZ vertex is modified by additional tadpole contributions, which lead to further diagrams, that have to be taken into account in the computation of the decay width. They are shown in Fig. 11. As the formula for the vertex corrections and counterterms in terms of the scalar one-, two-and three-point functions are quite lengthy, we do not display them explicitly here.

Numerical Analysis
For the computation of the NLO EW corrections to the Higgs decay widths described in the previous section we have performed two independent calculations. Both of them employed the Mathematica package FeynArts 3.9 [62,63] to generate the amplitudes at LO and NLO in the general R ξ gauge. To this end, the model file for a CP-conserving 2HDM was used, which is already implemented in the package. Additionally, all tadpole and self-energy amplitudes, needed for the definition of the counterterms and wave function renormalization constants, have been generated in the general R ξ gauge. The contraction of the Dirac matrices and formulation of the results in terms of Passarino-Veltman functions has been done with FeynCalc 8.2.0 [64] in one calculation and with FormCalc [65] in the other. The dimensionally regularized [66,67] integrals have been evaluated numerically with the help of the C++ library LoopTools 2.12 [65].
For one of the two calculations the Python progam 2HDMCalc was developed that links FeynArts, generates the needed counterterms dynamically from the 2HDM Lagrangian by calling a Mathematica script and combines the LO, NLO and counterterms calculated by FeynCalc into the full partial decay widths. These are then evaluated numerically by linking LoopTools. Finally, the LO and NLO partial decay widths are written out for all renormalization schemes of the mixing angles introduced above. The outcome of this program was compared to the results of the second independent computation. All results agree within numerical errors.
In the following we specify the input parameters that we used for the numerical evaluation. As explained in section 3 we use the fine structure constant α at the Z boson mass scale, given by [68] α(M 2 Z ) = have only a small influence on our results. In order to be consistent with the ATLAS and CMS analyses, we follow the recommendation of the LHC Higgs Cross Section Working Group (HXSWG) [69,71] and use the following OS value for the top quark mass m t = 172.5 GeV . as recommended by [69]. Omitting CP violation we consider the CKM matrix to be real, with the CKM matrix elements given by [68]  The SM-like Higgs mass value, denoted by m H SM , has been set to [72] m H SM = 125.09 GeV . The IR divergences in the computation of the NLO corrections to the process H ± → W ± H/h require the inclusion of the real corrections to regularize the decay width. This introduces a dependence on the detector sensitivity ∆E for the resolution of the soft photons from the real corrections. We showed that this dependence is small [73]. For our analysis we fixed the value to ∆E = 10 GeV . (5.9) In the subsequently presented plots we only used 2HDM parameter sets that are not yet excluded by experiment and that fulfill certain theoretical constraints. These data sets have been generated with the tool ScannerS [74]. 9 The applied theoretical constraints require that the chosen CP-conserving vacuum is the global minimum [75], that the 2HDM potential is bounded from below [76] and that tree-level unitarity holds [77,78]. For consistency with experimental data the following conditions have been imposed. The electroweak precision constraints [79][80][81][82][83][84][85] have to be satisfied, i.e. the S, T, U variables [79] predicted by the model are within the 95% ellipsoid centered on the best fit point to the EW data. Indirect experimental constraints are due to loop processes involving charged Higgs bosons, that depend on tan β via the charged Higgs coupling to the fermions. They are mainly due to B physics observables [86][87][88] and the measurement of R b [89][90][91][92]. We have included the most recent bound of m H ± > ∼ 480 GeV for the type II and flipped 2HDM [93]. The results from LEP [94] and the recent ones from the LHC [95,96] 10 constrain the charged Higgs mass to be above O(100 GeV) depending on the model type. In order to check the compatibility with the LHC Higgs data ScannerS is interfaced with SusHi [98] which computes the Higgs production cross sections through gluon fusion and b-quark fusion at NNLO QCD. All other production cross sections are taken at NLO QCD [70]. The 2HDM decays were obtained from HDECAY [99,100]. Note that in the computation of these processes all EW corrections were consistently neglected, as they are not available for the 2HDM. The exclusion limits were checked by using HiggsBounds [101][102][103] and the compatibility with the observed signal for the 125 GeV Higgs boson was tested with HiggsSignals [104]. For further details we refer to [105].
In our numerical analysis we investigate the applicability of the various proposed renormalization schemes. The goal is to find a renormalization scheme for the 2HDM, that is process independent, gauge independent and numerically stable. All results that we show are for the 2HDM type II.

Gauge dependence of the KOSY scheme
We start by analyzing the gauge dependence of the partial decay width, introduced through the renormalization of the mixing angles α and β in the KOSY scheme. As an example we choose the charged Higgs boson decay into the W boson and the light CP-even scalar h corresponding to H SM , H ± → W ± h. For the renormalization of β we use the charged sector and call the renormalization scheme accordingly KOSY c . The corresponding angular counterterm δβ (1) is defined in Eqs. (3.53), while δα is given by Eq. (3.52). The size of the gauge dependence will be quantified by It parametrizes the deviation of the NLO partial decay width for an arbitrarily chosen gauge parameter ξ in the R ξ gauge from the reference decay width chosen to be the NLO width in the Feynman gauge, normalized to the reference value. For simplicity we only vary the gauge parameter ξ W and set ξ Z = 1. The 2HDM scenario Scen1 that we investigate is defined by the input parameters (5.11) Figure 12 shows the ξ W dependence of our process, , as a function of ξ W . The kinks in the figure are due to threshold effects in the B 0 functions entering the counterterms. In detail, the kinks are given by the following parameter configurations and counterterms Kink ξ W Kinematic point Origin With a relative variation of the NLO width of up to 20% due to the change of the gauge parameter, the figure clearly demonstrates the gauge dependence of the NLO decay width in the KOSY scheme. The explicit calculation shows that for large values of ξ W the partial decay width drops as −(m 2 H − m 2 h ) ln(ξ W ). This explicit gauge dependence makes a practical use of the KOSY scheme impossible as it leads to non-physical gauge dependences in the decay widths.

The processes Γ(H ± → W ± h/H) at NLO
We move on to the investigation of the size of the NLO corrections to the processes H ± → W ± h/H and their dependence on the renormalization scheme. In our scenarios h corresponds to the SM-like Higgs bosons. We define the quantity  (5.14) In Fig. 13  : on-shell tadpole-pinched, δβ (1) or δβ (2) KOSY c,o : gauge-dependent scheme, δβ (1) or δβ (2) .
The process-dependent renormalization refers to the renormalization of α via the process H → τ τ and of β via A → τ τ . The process-dependent renormalization can be performed by applying either the standard or the alternative tadpole scheme. The investigation of the decay widths shows, however, that all decays discussed in this analysis, i.e. H ± → W ± h/H and H → ZZ, are invariant with respect to a change of the tadpole scheme. 11 In the process-independent schemes we can choose to renormalize β either through the charged sector, with the counterterm given by δβ (1) , or through the CP-odd sector, with the counterterm given by δβ (2) . For the shown m H ± range the LO decay width varies from Γ LO = 0.0750 GeV at m H ± = 654 GeV to Γ LO = 0.1474 GeV at m H ± = 804 GeV.
In Fig. 13 (left) we show results for the process-dependent renormalization and for some representatives of the process-independent schemes, the pOS o , the p c and for comparison also the KOSY c scheme. As can be inferred from the left plot, the process-dependent renormalization leads to much larger NLO corrections than the other schemes. The NLO corrections can increase the LO width by more than a factor of three. For the process-independent renormalization schemes on the other hand, the NLO corrections are much milder and vary between about −11 to 20% depending on the renormalization scheme and the charged Higgs mass value (and discarding the unphysical KOSY scheme). This can be inferred from Fig. 13 (right) which displays the results for the process-independent schemes, where the β renormalization is performed both through the charged and through the CP-odd sector. 12 Provided that the same choice for the β renormalization is made, the OS tadpole-pinched scheme, pOS, leads to results closer to the KOSY scheme than the p tadpole-pinched scheme. This is due to the fact that the KOSY and the pOS scheme use the scale of the OS masses for the evaluation of the self-energies. Also note that the schemes which rely on the CP-odd sector for the renormalization of β, show a slightly weaker dependence on the mass of the charged Higgs boson, as the latter enters the counterterm δβ (2) only through a few diagrams (namely the tadpole contributions). An important conclusion, which can be drawn from the plots, is that the process-dependent renormalization scheme is less advisable due to the induced unnaturally large NLO corrections compared to the results in the other renormalization schemes.
Discarding the numerically unstable process-dependent scheme and the unphysical KOSY scheme, we can use the comparison of the results for p c and p o and the comparison of those for pOS c and pOS o to estimate the remaining theoretical uncertainty due to missing higher order corrections, based on a change of the renormalization scheme for β. In the same way we can estimate the uncertainty based on a variation of the renormalization scale by comparing the results for pOS o and p o or the results for pOS c and p c . In the investigated m H ± range from the lower to the upper end, the remaining uncertainty varies between 1% and 11%, when estimated from the scale change, and from close to 0 to 18%, when estimated from the change of the β renormalization scheme. Note also that the results in the tadpole-pinched scheme, when evaluated at the OS scale, are less affected by a change of the renormalization scheme for δβ than in the p scheme. The renormalization of β through the charged sector is less sensitive to the scale choice than δβ (2) , which uses the CP-odd sector, as can be inferred by comparing p c with pOS c on the one hand, and p o and pOS o on the other hand. Taking these as indicators for theoretical uncertainties, one might draw the conclusion that the pOS c scheme would be the best choice here. Finally, we note that the kinks, which are independent of the renormalization scheme, are due to the thresholds in the following counterterms and parameter configurations

Kink
Kinematic point Origin Fig. 14 we show the relative NLO corrections for H ± → W ± h as a function of the LO width for all generated scenarios compatible with the applied theoretical and experimental constraints. The colours indicate the results for the process-dependent scheme, the p tadpolepinched schemes, the OS tadpole-pinched schemes and the KOSY c scheme. The plots demonstrate that the process-dependent renormalization in general leads to relative NLO corrections that are one to two orders of magnitude above those obtained in the other schemes, which yield corrections of typically 13 a few percent up to 40%, as can be inferred from the right plot.
In Fig. 15 we show the relative NLO corrections for the process H ± → W ± H with the parameters given by Scen3, Eq. (5.14). In the plotted m A range the LO decay width, which does not depend on m A , is given by Γ LO = 4.0568 GeV. In the left plot we have included the results for the process-dependent renormalization, for pOS o , p c and KOSY c . The right plot includes all renormalization schemes but the process-dependent one. The relative corrections lie between about −7.70 to −7.97% in the investigated mass range. 14 Altogether the results for all schemes lie very close to each other, with the process-dependent scheme deviating the most from the remaining schemes, although the difference in ∆Γ is of maximally 0.16% only. This behaviour can be understood by looking at the counterterm for the NLO process, Eq. (4.11). The contributions from the angular counterterms δα and δβ come with the factor 1/t β−α , which is numerically very small in the SM-like limit h ≡ H SM . Therefore any difference in the renormalization schemes for the angles will barely manifest itself in the total NLO corrections. The zoomed in region in Fig. 15 (right) again shows that the KOSY scheme is closer to pOS than to the other schemes and that the usage of the OS scale in δβ is less sensitive to a change of the renormalization scheme, while the renormalization of β via the charged sector is less sensitive to a scale change than the one through the CP-odd sector.
14 The small mA mass range is due to the fact that all other parameter points for this scenario are excluded. Figure 15: Relative NLO corrections to H ± → W ± H for various renormalization schemes, with the 2HDM parameters given by Scen3, Eq. (5.14); left: with, right: without the process-dependent renormalization. In the right plot the lines for KOSY c and pOS c lie on top of each other.

The process Γ(H → ZZ) at NLO
We now turn to the discussion of the NLO corrections to the heavy Higgs boson decay into a pair of Z bosons, H → ZZ. The scenario we have chosen is given by In the left plot the process-dependent renormalization is included. Additionally we show representatives for process-independent schemes, the pOS o , the p c and the KOSY c scheme. Again the counterterm definition via tauonic heavy Higgs decays leads to much larger corrections than the other schemes. In the investigated mass range it can increase the LO decay width by more than a factor of two. The observed coincidence of the results for the process-independent and process-dependent renormalization schemes at m H = 690 GeV is accidental. The relative corrections in the process-dependent renormalization start to increase quickly again for different m H values. The NLO increase in the process-independent schemes, on the other hand, ranges from about -3 to 17% in the investigated parameter range. The right plot shows the same behaviour we have seen previously. The results in the KOSY and in the pOS scheme are closer to each other than to the p scheme. Furthermore, the change of the β renormalization scheme affects the pOS scheme less than the p scheme and the β renormalization through the charged sector is less sensitive to a change in the renormalization scale than the one through the CP-odd sector. Overall, in the investigated mass range, the theoretical uncertainty due to missing higher order corrections can be estimated to be of less than a percent to around 6% based on a scale  change, and it ranges from the permille level to about 4% when estimated from the change of the β renormalization scheme, discarding the numerically unstable process-dependent scheme. Figure 17 shows the relative NLO corrections ∆Γ H→ZZ for H → ZZ as a function of the LO width for all generated scenarios compatible with the applied theoretical and experimental constraints. The colours indicate the results for the various renormalization schemes. The plots clearly demonstrate the numerical instability of the process-dependent renormalization, which exceeds the relative corrections in the other schemes by one to two orders of magnitudes. For the process-independent schemes the relative corrections are typically of the order of a few percent to 40%, discarding the region with small LO widths.
Altogether we conclude, that the choice of the KOSY scheme for the renormalization of the angular counterterms is precluded due to its manifest gauge dependence. The choice of the process-dependent scheme is not advisable, as it leads to very large relative NLO corrections 15 . The process-independent tadpole-pinched schemes lead to results that are manifestly gaugeindependent and numerically stable. Among these schemes the OS tadpole-pinched scheme turns out to be more stable when changing the β renormalization scheme than the p scheme for our investigated scenarios.

Conclusions and Outlook
We have investigated the renormalization of the 2HDM with special focus on the mixing angles α and β which diagonalize the Higgs mass matrices. These angles are highly relevant for the phenomenology of the Higgs bosons as they enter the Higgs boson couplings and therefore all Higgs observables. We have shown that if the tadpoles are treated in the more usual approach, which we called 'standard tadpole', a process-independent definition of the angular counterterms leads to gauge-dependent decay amplitudes and thus to gauge-dependent physical observables. Therefore, the counterterms δα and δβ either have to be defined through a physical process, or the treatment of the tadpoles has to be changed. Following the 'alternative tadpole' scheme as proposed in [29] allows for a manifestly gauge-independent definition of the masses and in particular of the mixing angles.
In this work we presented several distinct renormalization schemes and investigated their implications by applying them to the NLO EW corrections in the decays H ± → W ± h, H ± → W ± H and H → ZZ. It was explicitly shown that the scheme presented in [23] leads to gauge-dependent decay widths. This scheme applies the standard tadpole scheme and relates the angular counterterms to the off-diagonal wave function renormalization constants. By using the alternative tadpole scheme together with the modified Higgs self-energies obtained from the application of the pinch technique we introduced the 'tadpole-pinched' scheme as a manifestly gauge-independent scheme for the angular counterterms. We furthermore investigated the process-dependent definition of δα and δβ through the decays H → τ τ and A → τ τ , respectively. In this scheme the angular counterterms are gauge dependent when the standard tadpole scheme is applied, they are gauge independent in case the alternative tadpole scheme is used. For the investigated decay processes and scenarios, the process-dependent scheme turned out to lead to unnaturally large relative NLO corrections. Based on the investigated parameter sets and decay widths this leads us to the conclusion to propose the tadpole-pinched scheme as the renormalization scheme for the mixing angles that is at the same time process independent, gauge independent and numerically stable.
In order to complete the renormalization of the 2HDM, also the renormalization of the softbreaking parameter m 2 12 has to be investigated. This parameter appears in the couplings of the Higgs self-interactions and hence impacts the Higgs-to-Higgs decay widths. The renormalization of m 2 12 and the phenomenological investigation of the implications of the higher order corrections for Higgs phenomenology will be the subject of a follow-up paper.

A.1 Derivation of the Tadpole Scheme
We start by setting the notation and by presenting the standard scheme before we move on to the derivation of the tadpole scheme in the 2HDM.
A.1.1 Setting of the notation and tadpole renormalization The expansion of the two Higgs doublets Φ 1 and Φ 2 about the VEVs, cf. Eq. (2.4), leads to the mass matrices that are obtained from the terms bilinear in the Higgs fields in the 2HDM potential. Due to CP-and charge conservation they decompose into 2 × 2 matrices for the neutral CP-even, neutral CP-odd and charged Higgs sector, respectively. As we have seen in Sec. 2 the minimum conditions of the potential require the tree-level tadpole parameters T 1 and T 2 to vanish. At lowest order they are given by Eqs. (2.15) and (2.16). These tadpole conditions can be exploited to eliminate m 11 and m 22 . Higher order corrections, however, lead to non-vanishing tadpole contributions that have to be taken into account. Applying Eqs. (2.15) and (2.16) we arrive at the following mass matrices Here we have explicitly kept the tadpole parameters although they vanish at tree level. This helps us to keep track of their non-vanishing contributions at higher orders when performing the renormalization program. The mass matrices are diagonalized by the rotation matrices R rotating the scalar fields from the gauge basis into the mass basis, cf. Eqs. (2.5)-(2.7), The scalar mass eigenstates with same quantum numbers, grouped into the doublets (H, h), A) and (G ± , H ± ), mix at higher orders. The wave function renormalization constants, introduced in Eqs. (3.14)-(3.16) for the three doublets, also develop non-vanishing mixing contributions and form 2 × 2 matrices with off-diagonal elements. In the following we will use a generic notation and denote with φ 1 and φ 2 the two scalars of the same doublet. With this notation we then have for Eqs. (3.14)-(3.16) For the diagonal mass matrices, denoted from now on generically by D 2 φ , we introduce the counterterm matrix δD 2 φ , which is a symmetric 2 × 2 matrix whose specific form will be determined below. With these definitions the renormalized self-energyΣ φ becomeŝ The self-energy Σ φ is a symmetric 2 × 2 matrix containing the 1PI self-energies of the scalar doublet (φ 1 , φ 2 ). We require OS renormalization conditions for the scalar Higgs fields yielding the following conditions for the counterterm δD 2 φ and the wave function renormalization constants δZ φ , (i = 1, 2) So far we have not specified δD 2 φ . Its exact form depends on the treatment of the tadpoles in the renormalization procedure and will be elaborated below. In order to guarantee the correct minimization conditions for the Higgs potential also at one-loop order, the tadpoles are renormalized asT where T 1 and T 2 are the sum of all one-loop tadpole contributions to the fields ρ 1 and ρ 2 , respectively, in the gauge basis. Applying the renormalization conditions we have for the tadpole counterterms the conditions In the mass basis we have and The renormalization conditions for the tadpoles are shown pictorially in Fig. 18.

A.1.2 Mass counterterms and wave function renormalization constants in the standard scheme
Regarding the renormalization of the masses, the bare mass of each particle in the 2HDM is split into a physical mass and a counterterm as specified in section 3. The VEVs v 1 and v 2 , respectively v, are fixed at one-loop level such that their values in the tree-level mass relations (a) for the scalars, derived by calculating explicitly Eqs. (A.4)-(A.6), lead to the OS physical masses at one-loop level. The shift from the bare parameter to the physical one-loop value is hence fully contained in the mass counterterms. In generic notation the diagonalized bare mass matrices read where the subscript 0 denotes the bare quantities and ϕ = α for the CP-even, ϕ = β for the CP-odd and charged doublets, respectively. We have explicitly kept the bare tadpole parameters to keep track of their renormalization. Taking into account the renormalization of the tadpole parameters given in Eq. (A.14) we arrive at the NLO counterterm for the mass matrix Insertion of Eq. (A.18) in the renormalization conditions (A.10)-(A.12) we get the field strength renormalization constants and mass counterterms in the standard scheme These formulae can easily be generalized to the fermion and gauge boson sector. There, however, no tadpole counterterms will be involved, as they are not part of the tree-level mass relations. The counterterms introduced in Eqs. (A.28)-(A.30) are in general gauge dependent, which is not a problem, as long as all gauge dependencies cancel in physical observables. Since the renormalized masses must be gauge independent, the bare masses must be gauge dependent as well.

A.1.3 Mass counterterms and wave function renormalization constants in the tadpole scheme
We have seen that in the standard tadpole scheme the correct vacuum is reproduced by renormalizing the VEVs at higher orders accordingly. Derived from the gauge-dependent loopcorrected potential, the VEVs themselves are gauge dependent. The counterterms and the bare masses, that are given in terms of the VEVs therefore become gauge dependent, as the physical OS masses are gauge independent. In the tadpole scheme [29] the same renormalization conditions as given in Eq. (A.14), respectively in Eq. (A. 16), are used. The crucial point, however, is the inclusion of the minimization conditions of the potential such that the mass and coupling counterterms can be defined in a gauge-independent way. This is achieved in the following way: In the alternative tadpole scheme the bare masses are expressed in terms of the tree-level VEVs. As the tree-level VEVs are gauge independent, the bare masses do not depend on the gauge choice either. In order to still reproduce the correct minimum at higher orders, the VEVs acquire a shift. This shift now affects the counterterms and not the bare masses, as the latter are expressed in terms of the tree-level VEVs. The gauge dependences related to the VEV shifts cancels those of the counterterms, so that the counterterms become gauge independent themselves. Together with the gauge-independent bare masses the OS renormalized masses are gauge independent as they should be. The VEVs are hence shifted when going from LO to NLO as We emphasize that v 1,2 represent the tree-level values of the VEVs. The shifts δv 1,2 are fixed by the minimization, that is, by the tadpole conditions. The tadpole parameters are given in terms of the VEVs, cf. Eqs. (2.15) and (2.16), so that a shift in the VEVs corresponds to a shift in the tadpole parameters. Note that we apply the term 'shift' here in order to describe the changes of the parameters due to the VEV shifts, and to differ these from the counterterms for the chosen set of independent parameters.
The shifts in the VEVs are propagated into all parameters that depend on the VEVs. These shifts are determined as follows: (i) Express the parameters in terms of v 1 and v 2 . (ii) Perform the shifts Eq. (A.31) of the VEVs. (iii) Apply the tree-level relations between the VEVs and the various parameters to remove the redundant parameters m 2 11 , m 2 22 and/or to simplify the expressions as convenient.
Thus by shifting and subsequently applying the tadpole conditions Eqs. (2.15) and (2.16) we obtain Since the VEVs are determined order by order by applying the VEV shifts such that the tadpole conditions (2.15) and (2.16) hold we identify on the right-hand side of both equations the shift of the tadpole parameters induced by the shift of the VEVs with the counterterms δT 1 and δT 2 . By comparing the coefficients of δv 1,2 in Eqs. (A.32) and (A.33) with the elements of the CP-even mass matrix given in Eq. (A.1) the following relation between the VEV shifts and the tadpole counterterms, that determine δv 1,2 can be derived Rotation to the mass basis yields By applying the renormalization condition depicted diagrammatically in Fig. 18, the shift can be interpreted as a connected tadpole diagram, containing the Higgs tadpole and its propagator at zero momentum transfer, where h i ∈ {H, h} stands for the physical Higgs particles. For the consistent application of the tadpole scheme the VEV shifts have to be applied wherever the VEVs appear explicitly. As the calculation of the tadpole diagrams is usually performed in the mass basis, but the VEV shifts are introduced most conveniently in the gauge basis, we give the relation between the two bases, For the illustration of the implications of the tadpole scheme we consider a specific example, namely the NLO effects of the VEV shifts on the CP-odd mass matrix given in Eq. (A.2). The application of the shifts requires the replacement of the tadpoles by T i + δT i , with the δT i given in Eqs. (A.32) and (A. 33), and the replacement of all occurring VEVs by v i + δv i so that we have Having applied the shifts, we can now use the tree-level relations again to eliminate the last matrix in Eq. (A.38), as the tadpole parameters vanish at tree-level. The rotation to the mass basis is performed by applying the rotation matrix R(β) which is defined as the matrix diagonalizing the tree-level mass matrix M 2 η . We get where we applied the definition of Λ 5 Eq. (4.24) and the tree-level relation for the mass of the pseudoscalar [12,23] We furthermore applied the definition of the tadpole matrix in the mass basis, Eq. (A.18). In the last line we defined the terms ∆ G 0 G 0 , ∆ G 0 A and ∆ AA that contain all effects of the VEV shifts on the physical mass matrix D η . These shifts can be further evaluated. In order to do so, we introduce the coupling constants for the trilinear Higgs couplings [23] g By using the explicit form of the tadpole counterterm δT G 0 G 0 given in Eq. (A.22) the vanishing the remaining shifts for the whole scalar sector. The total shift of the mass matrices, which is given by the shifts ∆D induced by the NLO shifts of the VEVs and by the mass counterterms, then reads with the explicit form of the additional mass shifts (i = 1, 2) we obtain by inserting Eq. (A.50) in Eq. (A.9) the following form of the renormalized self-energy, And finally the counterterms and wave function renormalization constants in the tadpole scheme read These results can be generalized to the gauge boson and fermion sectors. The application of the tadpole scheme hence requires a redefinition of the self-energies as depicted diagrammatically in Fig. 19. In the gauge and fermion sectors this implies that the tadpole diagrams of the scalar Higgs bosons that couple to the gauge boson and fermion, respectively, have to be included in their self-energy. Furthermore, in the scalar sector the tadpole counterterms drop out of the definition of the wave function renormalization constants and mass counterterms. 16 The VEV shifts introduced in Eq. (A.31) also have implications for the coupling constants of the vertices. Let us consider the example of the Higgs H coupling to a pair of Z µ Z ν bosons. 16 In the gauge and fermion sectors they do not appear anyway as the mass matrices do not depend on m 2 11 and m 2 22 that are traded for the tadpoles.
iΣ tad (p 2 ) := + + Figure 19: Modified self-energy iΣ tad (p 2 ) in the tadpole scheme, consisting of all 1PI self-energy diagrams together with the one-loop tadpole diagrams, indicated by a gray blob.
Defining the needed coupling constants through the Feynman rules HZ µ Z ν : ig HZZ g µν (A.57) HHZ µ Z ν : ig HHZZ g µν , (A. 58) we have (A.60) The shifts Eq. (A.31) introduce a shift in the coupling constants. In order to perform this shift consistently, the coupling constants must be expressed in terms of the VEVs v 1 and v 2 . When doing so, care has to be taken, to differentiate between the angles α and β in the sense of mixing angles and β in the sense of the ratio of the VEVs, cf. Eq. (2.9), and α in the sense of the ratio of the 2HDM parameters 17 given in Eq. (2.12). The VEV shifts only affect the latter two.
The quartic coupling obviously does not receive any shift, while g HZZ contains β as ratio of the VEVs so that it receives a shift. The angle α is a mixing angle here. At NLO we therefore have to make the replacement ig HZZ → ig HZZ + ig The above result can be generalized to the whole 2HDM. In the tadpole scheme additional virtual vertex corrections have to be taken into account that manifest themselves in form of tadpdole vertex diagrams. The rule to be applied here is, that every 2HDM trilinear vertex receives corrections if the resulting quartic coupling constant exists, that connects the original trilinear vertex to the CP-even Higgs, that is attached to the tadpole. In the case above the vertex g HHZZ exists, so that the vertex acquires a tadpole contribution with H, but not with h, as the vertex g hHZZ does not exist.
As last example we look at the coupling between W ± µ , H ± (p ) and h(p), where p (p) denotes the outgoing (incoming) momentum of H ± (h). The Feynman rule for the coupling is given by with the coupling constant Both angles in this coupling are true mixing angles, so that no VEV shift has to be applied. Therefore the vertex does not change in the tadpole scheme. This is in accordance with the rule given above: There exists no vertex g W ± H ± hh nor g W ± H ± hH that could connect a tadpole with h or H to the trilinear vertex.

A.2 Rules for the Tadpole Scheme in the 2HDM
In this Appendix we summarize all rules of the tadpole scheme for the 2HDM at NLO. The general rules are: Self-energies: The self-energies in the wave function renormalization constants and counterterms change such that they contain additional tadpole contributions: Σ(p 2 ) → Σ tad (p 2 ).
Vertex corrections: In the virtual vertex corrections additional tadpole contributions have to be taken into account if the resulting coupling exists in the 2HDM.
Explicitly, this means that the following counterterms are the same in the standard and the alternative tadpole scheme: Counterterms independent of the choice of the tadpole scheme: for all possible combinations of fermions F , gauge bosons V , ghosts U , scalars S and φ i,j ≡ H, h, G 0 , A, G ± , H ± within the 2HDM.
The following counterterms and wave function renormalization constants depend on the choice of the tadpole scheme. We give the relations between the standard tadpole scheme, denoted by the superscript 'stand', and the alternative tadpole scheme denoted by the superscript 'tad'. The subscript 'trunc' means, that all spinors, all Lorentz structure of the vector bosons and the Lorentz structure of the coupling has been suppressed where applicable.
Tadpole-scheme-dependent counterterms: Gauge sector: Fermion sector: Scalar sector: (δZ φ i φ j ) tad = (δZ φ i φ j ) stand (A.70) We encounter additional contributions to the vertices when changing from the standard to the tadpole scheme. Here below, the g denote the coupling constants, i.e. we have suppressed the Lorentz structure of the vertex where applicable.
Triple scalar vertices: 71) for all scalars φ i,j,k ≡ H, h, G 0 , A, G ± , H ± , wherever the resulting quartic couplings λ φ i φ j φ k h and λ φ i φ j φ k H exist in the 2HDM. Scalar-vector-vector vertices: for all scalars φ i,j,k ≡ H, h, G 0 , A, G ± , H ± , and gauge bosons V j,k ≡ γ, Z, W ± , wherever the resulting quartic couplings λ φ i V j V k h and λ φ i V j V k H exist in the 2HDM.