Documentation of TauSpinner algorithms -- program for simulating spin effects in tau-lepton production at LHC

The tau-lepton plays an important role in the physics program at the Large Hadron Collider (LHC). It offers a probe in searches for New Physics and to measure parameters of the Standard Model. Spin of tau lepton is an interesting phenomenological quantity for separation of signal from background or in measuring properties of particles. A proper treatment of tau spin effects in the Monte Carlo simulations is important for understanding the detector acceptance as well as for the measurements/use of tau polarization and tau spin correlations. The TauSpinner is a tool to modify tau$spin effects in any sample. Also production matrix elements can be modified with its help. The general principle of TauSpinner algorithm rely on the kinematic of outgoing particles (or jets) only, and to average over all possible initial states accordingly to assumed parton distributions and cross-sections. Incoming on mass-shell parton four-momenta are reconstructed in approximation, but energy-momentum conservation is obeyed. Samples of events featuring tau leptons with their decay products, can be used as an input. The information on the polarization and spin correlations, and production/decay matrix elements is reconstructed from the kinematics of the tau lepton(s) (also nu_tau in case of W-mediated processes) and tau decay products (also jets). No other information is required. By calculating spin (or production/decay matrix element) weights, attributed to each event, it enables evaluation of modifications without regenerating events. Elements of program functionalities together with explanation of the use are given in Appendices. Present publication collect program documentation, previously scattered over several publiications.


Introduction 2 Theoretical basis
We start discussion from presenting exact formula on the cross section for pp → ττ X or pp → τ ν τ X processes, which requires knowledge of the whole process matrix elements, including decays of τ leptons. With about 20 τ decay channels we end-up with 400 possible configurations for each τ-pair production process. The exact formula can not be therfore used in practice. Nonetheless, it can serve as a starting point for the more practical to use ones, also those involving approximations. The formulae presented below for the H(Z/γ * ) → ττ can be applied also to the H ± (W ± ) → τ ν τ case. Then the second τ of the formulae should be replaced by ν τ . Consequently, part of formula which describes its decay is dropped out.

Exact formula
The basic formula for cross-section for the process pp → τ + τ − X; τ + → Y +ν ; τ − → Y ′− ν reads: Where M denotes complete matrix element. Obvious terms due to average over initial state spin degrees of freedom and flux factor are dropped out. Because of the τ lepton width being very small, the τ-lepton propagator module squared, to a perfect approximation, reduces to the Dirac δ function in lepton virtuality and so the infinitesimal element in the kinematical phase-space factorizes into dΩ prod dΩ τ + dΩ τ − . The formula (1) is not very convenient, because of over 20 possible τ decay channels which lead to 400 distinct processes and matrix elements. Fortunately, as for the phase space parametrization, because of a narrow τ width, matrix element can be also expressed as product of amplitudes for production and decays, summed over each τ spin projections (λ = 1, 2) on at this point unspecified axes The cross-section can be expressed then with the rewritten exact formula, which thanks to Fierz identity decouples summations of the decay and production matrix elements. It introduces the correlations back in a form of spin weight, This formula is a core of TauSpinner algorithm, but not only. It is used in Tauola τ decay library interface as foundation of the "after-burn" approach [6] and also in KKMC Monte Carlo [7]. For the details of geometrical formulation and naming conventions Ref. [8] may be useful. The spin weight wt spin is dimensionless, and contains information of all spin effects transmitted from the production to the decay of τ leptons, The indices i, j, of the adjoint (to 1/2 representations of fermions) rotation group representation run over t, x, y, z. The polarimetric vector h i τ ± is calculated from formula where σ i λ,λ stands for Pauli matrices. The h i τ ± are normalised further to set their time-like component to 1, namely . We will often use shortened form for the notation h i τ or h i ± . The superscript i will be dropped if we will have all components in mind.
The spin correlation matrix R i, j is calculated from formula and also normalized to set its time-like component to 1, namely R i, j = R i, j R t,t . Formula (4) is universal, it does not depend on the production process, and there are no approximations introduced. It features useful properties: 0 < wt spin < 4 and weight average < wt spin >= 1. Use of expressions (5) and (6) moves our description from language of spin amplitudes toward the language of cross-sections and probability distributions; the amplitudes squared contribute only.
All complexities are in calculation of R i, j spin correlation matrix for production and h i τ + h j τ − polarimetric vectors for the decay of τ leptons. The R i, j matrix describe the full spin correlation between the two τ's as well as the individual spin states of the τ's. The R i, j depends on kinematics of production process only, and h i τ + h j τ − respectively on kinematics of τ + and τ − . The explicit definition of the matrix R i, j and vectors h i τ + , h j τ − is rather lengthy and also well known, see [9,10] for detailed definitions. It is decay channel dependent.
So far we discussed case of τ pair in the final state, that is the case of H → ττ or Z/γ * → ττ processes. Case of only one τ in final state, that is the case of W ± → τ ν τ or H ± → τ ν τ is much simpler. The R i, j matrix is replaced by a vector of components R t , R x , R y , R z = 1, 0, 0, ∓1 respectively for W and charged Higgs decays. Sum over two indices is thus reduced to sum over one index only, axis z is the direction of τ ν from W ± /H ± decay as seen from τ-lepton rest-frame, spin weight read thus as wt spin = 1 ± h z τ . The frames choices (presented in the next Section) for that case could have been simpler than for Z/H decays, no details of transverse directions are then needed.
At this point we would like to decouple calculation of wt spin as much as possible from the PDF's dependence. Note, that this is needed for the Z/γ * mediated processes due to V − A structure of Z boson and γ * couplings to fermions. For the H or W mediated processes spin polarization state does not depend on flavours of incoming partons. For the Z/γ * the R i, j used in calculation of wt spin , is taken as weighted average (with PDFs and production matrix elements squared) over all flavour configurations, as in the following equation: No approximation is introduced in this way, denominator of Eq. (8) cancels explicitly the corresponding factor of Eq. (7) (that is also the integrand of the following Eq. (9), first line). To introduce corrections due to different spin effects, modified production process or decay model, in the generated sample (i.e. without re-generation of events) one can define the weight W T , representing a ratio of new to old cross-sections at each point in the phase-space.
The Eq. (7) for modified cross-section, takes then the form × wt spin × W T. 1 The double sum, runs in principle over all possible flavours for the two incoming partons: g,u,ū,d,d,c,c,s,s,b,b.
Variant of Eq. (9), described in the footnote 2 , is used in Tauola.
In general the weight W T , factorises into multiplicative components/weights: production (wt prod ), decay (wt τ ± decay ) and ratio of spin correlation/polarization weights, new to old one wt spin new /wt spin old : W T = wt prod wt τ + decay wt τ − decay wt spin new /wt spin old .
The first three terms of the weight, represent modification of the matrix elements for production and decays The combined weight as in Eq. (10), is nothing else than ratios of spin averaged amplitudes squared for the whole process; new to old. The eventual changes in the PDFs parametrizations should be taken into account in calculation of wt prod and also in ratio wt spin new /wt spin old in formula (10). In case production process or decay models are not modified, respectively wt prod and wt τ ± decay are equal to 1. Ratio of spin weights wt spin new /wt spin old allows for introduction of new spin effects. In case of originally unpolarized sample, the wt spin old = 1 and wt spin new alone allows to introduce desired spin effects.
Formula (9) can be rewritten to equivalent form, which underlines separation useful for the organization of program, with segments for production of τ leptons, their decays and the spin weights. Let us recall that the spin weight wt spin is the only part depending simultaneously on production and decay kinematics and for production of τ lepton pair has always a form of Eq. (4). Ratio of this weight calculated from new and old matrix elements has then to be taken. Reorganizing some terms we can finally write: dσ = dσ prod dΓ τ + decay dΓ τ − decay wt spin old wt spin new wt spin old , where All four component of the overall weight (10) are given separately to underline, wherever possible, their dependence on all or only sub-set of phase-space dimensions: production and/or decay. Note that for W ± /H ± → τ ν τ we have not prepared so far any option which allows to modify production processes.

WT calculations in TauSpinner
The formulae listed above are used in TauSpinner for calculation of distinct components of W T weight, see Eq. (10). Let us start with the main weight of the spin effects, wt spin defined by Eq. (4). With this weight, for samples where spin effects of τ production are absent, they can be inserted into decays distributions. Alternatively, its inverse can be used to remove spin effects from the sample where they are taken into account. A building block of this weight introducing the dependence on the production process, the matrix R i, j , is given by the formula (8). The other weights, w prod and w τ ± decay are defined by formulae (11,12). 2 It would seem natural to make a choice of the initial flavours, on the basis of probabilities obtained from production density 1/N normalization f (x 1 ,...) f (x 2 ,...) ∑ λ 1 ,λ 2 |M prod parton level | 2 over the flavours. This would also seem quite similar to what is used in Tauola interface [6]. Such solution would not be correct. The dependence on the decay configuration, already fixed and stored for the event, would be then ignored. In case of Tauola interface this information can be obtained from the event record.
It is useful to introduce notation P z τ = R t,z = R z,t , which represents longitudinal polarization of the single τ (if one integrates out possible configurations of the other τ). In usually sufficient approximation of helicity states, that is when transverse momenta of τ decay products are neglected, as is the case of ultrarelativistic τ leptons, the P z τ is the only non-trivial (dynamic dependent) element of matrix R. All others are equal ±1, or can be set to 0. We will return to the details later.
In calculation of wt spin and w prod physics uncertainty depend on accuratness of factorization assumptions for separation of patron level matrix elements and PDFs; also on choice of particular PDFs parametrization. Nothing of principle would change if instead of 2 → 2 body production matrix element one would use 2 → 2 + n where n denotes additional partons/jets or other particles (of summed spin-states). Except the choice of PDFs (eventually also model of underlying event interactions), for a parton level matrix elements one has to make a careful choice of how the hard scattering kinematics is reconstructed from information available in the event.
Calculation of wt τ ± decay involves modelling of τ decays only. It is of interest for studies of the resulting systematic errors. This part of the code is covered by the Tauola library.
We will return to these points later, but one has to keep it in mind already now, while reading the following.

Kinematical frames
Components necessary for the W T weight are calculated in different frames. This is a correct approach as long as details of boosts and rotations connecting frames are meticulously followed. Also spin states of the τ leptons may be defined in different frames than four momenta of the hard process. Bremsstrahlung photons and properties of the matrix element in their presence, require dedicated treatment. Basic formula (1) can in principle provide exact results. In practice, we have to introduce approximations to adapt to the fact that matrix elements are calculated from the parton level amplitudes, and to the way how these amplitudes are reconstructed by TauSpinner algorithm where only four momenta of outgoing τ leptons and partons can be used. Definition of quantization frames is an essential element of the arrangements, let us provide now the necessary details.
For calculating all components of event weights (to prepare kinematical configurations of the events necessary for that purpose) we define in total four frames (group of frames): A : Basic frame (starting point) for kinematic transformations and all other frames definitions. The rest frame of τ-pair is used (and not of τ-pair with the final state bremsstrahlung photons combined). Such a choice is possible thanks to properties of bremsstrahlung amplitudes. Emitted photons do not carry out the spin. This non-trivial observation is exploited also in Photos Monte Carlo phase space parametrization, see Refs. [11,12]. The definition of the rest frame of the τ lepton pairs, is completed with incoming partons set along the z-axis.
A ′ : The rest frame of the τ lepton pairs and final state bremsstrahlung photons included. Incoming partons are again set along the z-axis. This frame is used for calculating production weight wt prod and spin correlation matrix R i, j ; if no bremsstrahlung photons are present it is equivalent 3 to frame A. In every case we reconstruct x 1 , x 2 , the arguments of PDF function, from virtuality M and longitudinal to the beam direction component (p L ) of the intermediate state (sum of momenta of τ + τ − nγ) momentum in lab frame. The center of mass energy of pp scattering is used to calculate . The nγ is assumed to correspond to final state bremsstrahlung associated with τ-pair production and its momentum has to be taken into account, but only in calculation of virtuality for intermediate Z/γ * state in case of 2 → 2 amplitudes 4 . Note that A ′ is used for x 1 , x 2 calculation only. The 2 → 2 hard process scattering angle is calculated in frame A.
B : The rest frame of the τ lepton pairs, with τ-leptons along the z-axis. It is used for common orientation of production and decay coordinate systems, in this frame defined are polarimetric vectors h ± . Transformation A → B requires rotation.
C τ + ,C τ − : The rest frames of individual τ leptons, with direction of the boost to τ pair rest-frame along the z-axis.
Transformation B → C τ + or C τ − requires boost only.

Object
Frame of object calculation Frame of object use Comment Table 1: Frames used for calculation of components used in weights calculation, Eq. (10). Specified are also frames to which these components need to be transformed.
The rest frames of individual τ leptons, with ν τ along the z-axis. It is used for calculating decay weights. In fact in the frames D τ ± also polarimetric vectors h τ ± are calculated and then rotated back to frame C τ ± .
In Table 1, we summarize which frames are used for calculations of different variables and then to which frames they are boosted (rotated) before being used in formula (7). Let us stress, that it is not only important to define all frames of type A − D but also Lorentz transformations between them.
The frame of type A is used for calculation of hard process amplitudes. If only longitudinal components of spin density matrix are taken into account, we do not need to define in detail transverse directions (quantization frame versors) with respect of τ momenta, and control the corresponding details of boost methods to the rest frames of τ's is not needed. However, the code of TauSpinner controls such details of the boosts, because it is prepared to handle also cases when complete spin density matrix is used.
Frame of type B and frames of type C have common direction of the z-axis, also x-and y-axes coincide. Frame B is used for R i, j calculation if transverse degree of freedom are taken into account. Calculation of the spin weight, that is contraction of indices in wt spin = ∑ i j R i, j h i τ + h j τ − , is performed in frames C τ + , C τ − . It enable control on transverse spin effects. Frames of type D are used for calculation of decay matrix elements and h i τ + , h i τ − polarimetric vectors with methods from Tauola library, see Ref. [10] for details. Let us note, that we have used such frames in many applications already at a time of LEP analyzes, they were useful for comparisons and tests. These choices once established, can be easily modified, whenever necessary, e.g. to fulfil constraints of particular conventions of spin amplitudes calculations, see eg. in Ref. [8]. With the definition of frames completed, we can continue our discussion.

Variants for definition of τ lepton pair rest frame.
Let us comment on the possible future improvement for the definition of the frame A, potentially useful freedom of choice is available. At present, in TauSpinner definiton of lepton pair rest-frame A, beam direction, necessary to define Born-level scattering angle, the simplest possible ansatz, that in a given event there are no jets of substantial p T , is used. This can be improved without a need to re-do the design of the presented above tree of frames definitions. The necessary studies for matrix elements featuring one or two high p T jets are already documented in Refs. [13,14]. The corresponding modifications of the code were not yet introduced into TauSpinner, as there are no clear indications of their numerical importance, but it is an input for the possible forthcoming improvement.

Production weight wt prod and spin correlation matrix R i, j
The R i, j matrix is not a Lorentz invariant tensor. In fact, its indices run over coordinates in two frames, C τ + and C τ − at the same time as the input for matrix element calculation the four momenta of A frame are used. For details of relations between frames and elements used in weights calculation, see also Table 1.
For calculation of production weight wt prod and spin correlation matrix R i, j , centre-of-mass frame of the τ lepton pair with incoming partons being along z-axis (frame of type A) is not the optimal one. It is nonetheless used if there is no interest in transverse spin degrees of freedom. Frame of type B, with τ-leptons along z-axis, is more convenient and with x and y axes parallel to the ones of frames C τ ± completes the most convenient setup. In fact all these three frames B and C τ ± are used simultaneously for calculation of spin weight of Eq. (4).
The kinematical frames for τ-pair production and τ decays need to be explicitly related by the set of Lorentz transformations with the frames of matrix elements calculations. Let us recall details in a form of algorithmic steps, which need to be performed on all decay products of τ leptons before calculation of R i, j and h i/ j τ ± (for h i/ j τ ± further steps will be still needed): • Boost both τ leptons and their decay products from the laboratory frame to rest-frame of ττ pair.
• Rotate all to the frame, where direction of incoming partons is along the z-axis (frame A).
• Rotate so that τ leptons are set along the z-axis and incoming partons remain in z − y plane (frames C τ ± have to be correlated with B by a boost along z axis). These frames are used for R i, j spin states definition too.
-Presence of the final state bremstrahlung from τ pair brings complication, because the rest frame of the τ pair is not the rest frame of the resonance which decayed into τ pair. For the sake of kinematical transformations, frames A and B defined as above are used. However, for calculation of matrix elements we use frame A ′ , where efffectively bremsstrahlung photons are absorbed into τ lepton momenta. For calculation of x 1 , x 2 we use invariant mass, which can be easily calculated in frame A ′ .

Decay weight wt τ ± decay and polarimetric vectors h τ ±
The decay weight wt τ ± decay and polarimetric vectors h τ ± are calculated for each outgoing τ lepton in its rest frame, denoted as frame D τ ± that is a frame where neutrino from τ lepton decay, ν τ is along z-axis. The Lorentz transformations of the previous sub-section have to be supplemented with the one relating frames C τ ± with D τ ± .
• Boost all decay products along z-axis to the rest frame of the τ ± lepton (frames C τ ± ).
• Calculate in frames D τ ± polarimetic vectors h i τ ± , rotate them back from D τ ± to C τ ± . Only then, they can be contracted with R i, j for wt spin calculation.
Decay weights and polarimetic vectors for τ + and τ − are calculated in frame type D τ ± , using code from Tauola, which sets requirement of ν τ along z-axis.

Exact and approximate spin weight wt spin
Basic formula (3) of Section 2.1 is exact. Because of approximated calculations (introduction of parton level amplitudes) or as explained later in the Section, to obtain easier to interpret physical picture (valid in ultrarelativistic case only) we can introduce simplifications. No such simplification was introduced yet. Let us recall first that according to R i, j and h i τ ± definitions, R t,t = 1 and h t ± = 1 (we can use shorter notation h i ± = h i τ ± ).

Neglecting m 2 τ terms
With the choices of kinematical frames as defined in Section 3.1 and neglecting m 2 τ terms following the ultrarelativistic limit, one gets that R z,z = sign = ±1 (the sign depends whether decaying object is a scalar or a vector boson). Also, non-diagonal terms R t,x = 0, R t,y = 0, R x,t = 0, R y,t = 0, R z,x = 0, R z,y = 0, R x,z = 0, R y,z = 0, independently whether the decay of spin zero (scalar) or spin one (vector) state is considered. The formula (4) reduces to: or if notation of individual τ longitudinal polarization P z τ is used (for Higgs decays P z τ = 0), Let us stress, that the formulae (15,16) are quantization frame dependent. One has to be careful especially with signs, if the kinematical transformations, as explained in Section 3.1 are not followed. Distinct conventions for defining those frames are used in the literature.

Neglecting transverse spin correlations and spin state probabilities
Further approximation (not always used) is to include longitudinal spin effects only, means that terms R i, j are set to zero for i, j = x, y. In general, those terms can be large, but as they result in dependencies of transverse with respect to τ direction components of τ decay products momenta, often do not lead to measurable effects and can be dropped. With this approximation, formula (16) becomes Using formula (17), and resulting approximation, one will be able to introduce probabilities for spin states p z τ which will be defined later. For decay of scalar resonance (eg. Higgs) we define (equal) probabilities for τ pair configurations being left-right (+) or right-left (−). We rewrite formula (17) and decompose wt spin into sum of wt + spin and wt − spin : For decay of scalar into ττ pair, the spin state probability p z τ = 1 2 (the subscript z denotes spin projection on z axes of frames type C). For ultrarelativistic cases, the p z τ means also probability of τ helicity states. For the decay of vector resonance, the probability p z τ depends on intermediate boson virtuality and τ direction, we can rewrite formula (17) as Note that p z τ can be interpreted as a probability of creating first τ of the pair in a spin state (+) along z direction, whereas P z τ as its longitudinal polarization. Obviously, from comparison of above eqs. (17) and (19), following relation occur: which can be rewritten as

Longitudinal versus transverse spin correlations
In the formalism discussed above, there is seemingly rather minor difference between situation when complete spin effects are taken into account or only longitudinal ones. The case of longitudinal spin effects only, means that all R i, j components are set to zero except R t,z , R z,t and R t,t = 1, R z,z = ±1; the sign depends whether decaying object is a scalar or a vector and R t,z is equal to longitudinal τ polarisation P z τ . No difference in the re-weighting algorithm is needed in our solution. The transverse spin correlations are introduced with non-zero R x,x , R y,y and also off-diagonal R x,y , R y,x , without any modifications of the t and z components.

The 2 → 2 parton level processes
The complete event weight W T , is factorised into multiplicative components: production (wt prod ), decay (wt τ ± decay ) and spin correlation/polarization (wt spin ): If the sample is polarized the last factor should read wt spin new /wt spin old . Implemented methods calculate internally complete set of weights but as a default one, return the spin correlation/polarization weight, as most often used in the applications 5 . Other weights are returned by supplementary 5 Calculation is invoked with the call to calculateWeightFromParticlesH method. The code organization is a consequence of implemented tests, availability of tested methods and available benchmarks. Not all features were expected to be included at the beginning: At first, only longitudinal spin effects were expected to be taken into account and for 2 → 2 at parton level processes only. That is why, methods from KORALZ [15] could have been used. Later processes 2 → 4 with two jets in final state made original methods focused around calculation of τ polarization i.e. R t,z less intuitive/convenient. Introduction of options featuring complete spin effects made such organization even less suitable. Still we preserve backward compatibility, frame for functioning code and stability of results for accumulated over years tests and examples for use. methods 6 .

Production weight
The production weight wt prod which is a sum over all configurations of τ leptons spin and flavours of incoming partons of matrix element squared multiplied by PDF's is calculated in frame A ′ . The parton level angular kinematic configuration are usually taken in A frame. The R i, j are calculated in A or A ′ frame, except the case when transverse spin effects are taken into account and frame B has to be used.

Decay weights
Calculation of the decay weights wt τ + decay and wt τ − decay is a by-product of polarimetric vectors of h i τ + and h j τ − calculation. This calculation is performed in frames D τ ± . The methods of Tauola used for h j τ − calculation, at the same time calculate decay matrix elements squared summed over the τ spin which are used in wt τ ± decay .

Spin correlation matrix R i, j
In most cases the spin correlation matrix R i, j is calculated in frame A or A ′ . The frame B is used if transverse spin effects are taken into account. Calculation of production weight and of spin correlation matrix R i, j are related with each other. It is important to realize that in the collinear approximation where σ ± prod denotes the sum of matrix element squared multiplied by PDF's, the contribution, where first τ has spin ±. In this approximation and for Eq. (23) difference between frames A and A ′ can be overcomed.
Note, that the weight wt prod is obtaines as a by-product of R i, j calculation, in fact of its R t,z , R z,t components 7 . In case of using effective Born cross-section or directly formula using vector/axial couplings and Z/γ * propagators (implemented in TauSpinner) this connection is less transparent.

Case of spin = 0 resonance
The R z,z = −1 if τ mass is small in comparison to τ-pair virtuality. The R x,x = R y,y = 1 for the scalar Higgs boson and R x,x = R y,y = −1 for the pseudoscalar, while the mixed scalar-pseudoscalar represents only a slightly more complicated case; the formula includes scalar-pseudoscalar mixing angle and non-zero off-diagonal R x,y , R y,x terms. This extension is important for simulation of Higgs CP parity signatures.

Case of spin = 1 resonance (Drell-Yan)
The R z,z = 1 for neutral vector boson, if τ mass is small in comparison of τ-pair virtuality. The non-zero offdiagonal terms are At present implementation, if electroweak corrections are neglected, the transverse spin correlations for neutral vector boson/photon exchange are neglected as well and R x,x = R y,y = R x,y = R y,x = 0. We are simply missing the code for the Born level expressions taking into account transverse spin effects. Fortunately, if EW corrections calculation is activated, thanks to ported from Tauola tables of SANC library [16] results, not only first order EW corrections become available, but transverse spin effects as well.

Case of spin = 2 resonance X (non SM Higgs)
In this case the R z,z = 1 as in Drell-Yan case. The non-zero off-diagonal terms are also the difference is that p z τ feature different dependence on τ lepton directions defined in τ lepton pair rest frame B. It is described by spherical polynomial of the fourth order, and not of the second order as for Drell-Yan. This case was introduced with the paper [2] and corresponding example of TauSpinner application.

Calculating wt spin weight
For the wt spin calculation formulae (15,16) are available, however unless required by the user, the default version is to omit transverse spin correlations, hence put R x,x = R y,y = R x,y = R y,x = 0, effectively reducing calculations to formula (17). In the following, we will exploit prepared in the previous sub-section components of R i, j for the individual cases.

Case of spin = 2 resonance X (non SM Higgs)
where sign = +1. Note, that this form depends on Born level distribution provided at the execution time with the user function, which is used in calculating R i, j , at present in our examples only P z τ is non-zero.

Case of spin = 1 resonance (Drell-Yan)
(28) where sign = +1 and R EW 1,1 , R EW 1,2 , R EW 2,1 and R EW 2,2 denote averaged (with production cross-section) pre-tabulated transverse spin correlations terms, calculated using O(α) EW corrections from SANC program 8 . The R x,x = R y,y = R x,y = R y,x = 1 should be set, in case we want to use O(α) EW corrections from SANC. Or one can also set some of them to 1 and remaining ones to 0 in case one is interested in studying only some components of transverse-spin effects.
The R x,x , R y,y , R x,y , R y,x are used directly as transverse components for R i, j in case of non-standard model calculation. If one is not interested in the actual size of transverse effects, as predicted by O(α) EW corrections, one can also modify hard coded pre initialized to zero values, or use discussed previously non SM Higgs options. Generated Drell-Yan samples may feature only partial spin effects, e.g. just (1.0 + sign · h z + h z − ) part of the spin effects. Then wt spin has to take this into account, the denominator of such form has to be introduced, to remove this partial effect. The numerator, given in eq. (28), will introduce the required (full) spin effect. Let us list the formulae for the typical cases: • If in the generated sample spin correlations between two τ leptons are included but no other effects, the weight should take the form: Such choice for the sample introduces dominant spin effect (the correlation), which at the same time is free of systematic errors. It is a consequence of the spin 1 state decaying to τ lepton pair.
• If in the generated sample spin correlations and partial (no angular dependence and no incoming quark flavour dependence) polarization effects are included, the weight should take the form: where the denominator again represents spin correlations present in the sample. One can expect the P z,part τ = P z τ (s, cos(θ = 0)) calculated for incoming electrons (incoming electron, because it is the easiest way to obtain average correction) reproduces correctly average individual τ polarization for events of τ pair virtuality close to the Z boson peak. Such a choice for the generated sample could be motivated as follows. It features all spin effects except those which depend on incoming partons PDF's, thus is free of related systematic errors, making it convenient for studying systematic errors from PDFs. The corrections introduced with wt no angular spin are even smaller than those of wt no pol spin .

Use of spin correlation/polarisation weight in special cases
If generated sample feature spin effects then the wt spin of formula (3) can be used as weight W T = 1/wt spin to remove spin effects from the generated sample.
The cases when only part of spin effects is taken into account, like in [17], more specifically the spin correlation but no effects due to vector and axial couplings to the intermediate Z/γ * state, can be corrected with the help of the appropriate weights as in examples for Eqs. (29, 30). On the other hand, such partial spin effects can be also removed completely, e.g. for consistency checks. Also the case, when spin correlations and average over the τ lepton angles within rest-frame of lepton pair are present, but no angular dependence is taken into account, can be checked in an analogous way.

Born cross-section
Until now, we were expressing formulae in the language of spin amplitudes and spin weights (spin correlation matrix R i, j ). However, often cross-sections calculated at the parton level are available for the explicit helicity states of the outgoing τ leptons. Several such lowest order formulae for calculating cross-sections of an analytical form, are implemented in the code of the TauSpinner itself, or of its use examples. Some of them represent legacy code, ported from different projects, others were coded for TauSpinner application as benchmarks on the EW parameters setting and are used in the published tests. For completness of documentation we decribe them in Appendix B.9.

Polarimetric vectors, polarization and helicity states
Calculation of the polarimetric vectors h τ is performed using code from Tauola library. The conventions is that last component is time-like) following convention used in FORTRAN. The same FORTRAN code calculates decay matrix element squared, which are returned as well, as auxiliary information.
Depending on which decay channel is present, calculation is performed using respective functions from Tauola library. Not for all τ decay channels, see Table 2, h τ is nonetheless calculated for 97.5% of the decay width is covered, i.e. polarimetric vector h is calculated, otherwise its spatial components are set to 0. This is not the case for Tauola library, if used for generating τ decays. Then only for channels featuring more than 5 pions polarimetric vector is calculated with approximation. Still then, it is not set to zero.
To calculate decay matrix elements, τ decay products have to be boosted into τ rest-frame: we first boost all products to τ − τ pair rest frame, then rotate from A to B so that τ of Z/H decay is along z-axis, intermediate step is boost from the τ − τ pair rest-frame to C τ , then rotation to have neutrino from τ decay along z axis (D τ frame); calculated for that purpose angles φ 2 , θ 2 are stored for rotating back of the polarimetric vector to the τ rest-frame B or C τ .
As a by-product of h τ calculations, wt τ ± decay of Eq. (12) is obtained.

Longitudinal polarization
Calculation of components of R i, j spin density matrix of τ-pair production in ultrarelativistic limit simplifies, it is then equivalet to calculation of modules for matrix elements of τ-pair production in a given helicity configurations.
Depending on the approach, either p z τ or P z τ is calculated from the other one, using respectively formula (21) or (20). The case, when electroweak effects are included, technically differs only slightly 9 , despite the fact that then the transverse spin effects can be taken into account. In all cases the relation R t,z = P z τ = holds, no pretabulated electroweak results are used for this component of R i, j . This can be easily activated though, e.g. if results for large s are needed. Around the Z peak the presently implemented electroweak library may be less suitable. Formula (19) is used for attribution of helicity states.

Case of spin = 0 resonance (scalar (Higgs)
In H → ττ decays, the probability of the helicity state denoted p z τ is equal to 0.5 for P z τ = 1 and for P z τ = −1 as well, see Table 3.

Case of spin = 1 resonance (Drell-Yan)
In Z/γ * → ττ decays, the probability of the helicity state denoted p z τ is a function of the τ scattering angle, θ, and the center of mass squared of the hard process, s. It is true at the Born level, and in ultrarelativistic limit [17].
Note that the rest frame of the production process and the rest frame of the τ-pair might not be the same due to photon bremsstrahlung in the τ pair production vertex. We have discussed this already in Section 3.1, but let us complete some more technical details now.
The center-of-mass system of τ-pairs is used for calculating cos θ angle and the direction of incoming partons is set along z-direction in the laboratory frame. Following definitions are used: where b + , b − denote 3-vectors of incoming beams, in the laboratory frame defined along the z-axis. And the definition of effective cos θ is constructed following Ref. [19] cos θ = cos θ 1 · sin θ 2 + cosθ 2 · sin θ 1 sin θ 1 + sin θ 2 .
At the Born level, for a given flavour of incoming partons and in the ultrarelativistic limit, the p z τ is defined [17] Table 3: Probability for the helicities of the of τ's from different origins [17].
where dσ d cosθ (s, cos θ, P z τ ) = (1 + cos 2 θ)F 0 (s) and F i (s) are four form factors which depend on the initial and the final state fermion couplings to the Z boson and the propagator [17].
In case of TauSpinner algorithm flavour of the initial partons is assumed to be unknown, the probability p z τ is calculated as a weighted average over all possible initial state quark configurations. This requires averaging (summing) over all production matrix elements square convoluted with the structure functions where indices i,j denoting flavours of incoming partons are now explicitly written. For calculating the spin weight, the polarisation P z τ of the single τ in a mixed quantum state, is calculated neglecting transverse spin degrees of freedom (ultrarelativistic limit). The P z τ is then linear function of probability for the helicity state 10 , Eq.(21).

Case of spin = 2 resonance X (not a Higgs)
In this case user provided function is required to replace dσ d cos θ , the Non-Standard physics effective Born parametrizations are not included in the TauSpinner library, dummy function is provided only. It can be replaced at the execution time with the user one, see Ref. [2] for details and in particular TAUOLA/TauSpinner/examples/ tau-reweight-test.cxx directory. This 2 → 2 parton level process function nonSM_adoptH (or nonSM_adopt) has the following arguments (int ID, double S, double cost, int H1, int H2, int key). It can be activated in the user code at the execution time, see Appendix B.2.

Attributing τ (or τ pair) helicity states
For the individual event, the τ helicity (+/−) is attributed stochastically. The probability for given configurations of the longitudinal polarisation of τ's of different origins [17] is shown in Table 3. These probabilities can be used, when decay kinematical configurations of the decays are not yet constructed.
However, once the configurations of decay products are generated, event-by-event calculation of the probability for helicity states, has to take into account τ decay's polarimetric vectors h + and h − . The (+) is chosen if probabilities, calculated with the formulae below, are smaller than generated random number RRR, otherwise (−) is taken: 10 Note that in paper [1] formula P z τ = (2p z τ − 1) is used instead of Eq. (21). There, different (rotated by π angle with respect to axis perpendicular to reaction frame is used. Similar convention matching may be necessary in programs as well, see e.g. nonSM_adopt(); nonSM_adoptH() methods of tau-reweight-test.cxx.

The 2 → 4 process
At this point we have completed discussion for the case when TauSpinner is identifying and working with parton level processes of 2 → 2 type. Let us turn now to the case when 2 → 4 parton level processes are used. This extension has been introduced very recently, see [5] for detailed description. A baseline algorithm remains essentially the same. The differences can be grouped into few parts

• Kinematics
For calculation of factorisation scale, instead of virtuality of τ lepton pair, virtuality of lepton pair combined with jets is used. In fact few options to define the Q 2 factorisation scale have been implemented. The definitions of frames type A, B, C τ ± and D τ ± as explained for case of 2 → 2 processes are used. In principle an additional rotation by an angle π around axis x would be again needed for the τ − in transformation from B, C τ − frames. This is due to convention used for spin amplitudes 11 Definition of the previously discussed frame A ′ needs to be modified. Not only bremsstrahlung photons need to be integrated into τ leptons to obtain kinematic configurations for matrix element calculation, but outgoing partons (jets) four-momenta have to be taken into account as well.
• Matrix Element Calculation Matrix element squared are calculated for given helicity states of τ leptons and flavours of incoming/outgoing partons in frame A'. In particular formula (23) is used for calculation of R z,t and R t,z , which are then used as of frame B, without any modifications. The EW scheme used and α s (Q 2 ) can be configured using several options defined in Ref. [5].
• Extension due to flavours of accompanying partons/jets All sums ∑ f lavour as defined in Section 2.2 should read now as sums over flavours of incoming and outgoing partons. Of course PDFs distributions are used for incoming partons only. This modification is extending re-weighting algoritm dependence on final state components consisting of partons (hadronic jets). This technically simple change may be thus of a concern because of additional types of the systematic errors involved.
Except these changes, algorithms remain not modified. However only longitudinal spin effects can be simulated because we use helicity states, that is ultrarelativistic approximation in matrix element calculation. In future, if necessary, we can extend the case to complete spin effects treatment. Only minor extension of the program will be needed. Approximation is in the matrix element calculation only. It may even be possible to use for that purpose presently available MadGraph amplitudes discussed in [5]. Better control of the complex phases will be however needed.

Summary
In this paper we have presented theoretical basis and description of algorithms used in the TauSpinner program for simulating spin effects in the production and decay of τ leptons in proton-proton collisions at LHC. The primary source of τ leptons are Drell-Yan processes of single W and Z boson production, often accompanied by the additional jets. Since Higgs boson discovery by ATLAS and CMS Collaboration in 2012, the H → ττ decay became the most promising channel for studying CP properties of its couplings to fermions. With more data available from Run II of LHC also τ decays in case of multi-boson production will become part of interesting signatures.
Discussed here algorithms allow for separation of the τ lepton production and decay in phenomenological studies. We have discussed options used for simplifying descriptions and approximations used.
Although the framework of TauSpinner is prepared for implementation of the non Standard Model couplings, we have not given too much attention to this possible development path. In the paper appendices, instead, we have recalled all available in the distribution examples and applications and pointed to documenting them publications, wherever possible. The technical information given here is very brief, detailes are delegated to the README files available with the distributed code.

A Project software organization
The code of TauSpinner is distributed as a part of the Tauolapp • testsERW , see App. B.8, • applications, general methods, useful to deal with event record contents and README file is included only. • paper-application-studies, see App. B.1, and Ref. [3] (redirects to TAUOLA/TauSpinner/examples/applications/applications-plots-and-paper), • paper-CP see App. B.3 and Ref. [4].
There is no need for the user to read details of all this multitude of files, except the one of the interest. The code of the basic examples is stored in the following two files: TAUOLA/TauSpinner/examples/tau-reweight-test.cxx for 2 → 2 processes, TAUOLA/TauSpinner/examples/example-VBF/example-VBF.cxx for 2 → 4 processes.

B General information on user interface methods
Over the years several methods of initializing and operating TauSpinner were developped [1,2,3,4,5]. Technical descriptions and discussions on systematic errors were provided, but remained scattered over several publications. Collecting them or providing references from one place is now of convenience.
The rather nontrivial component of analysis code is a correct reading and interpretation of information in the event record, namely finding outgoing τ leptons and their decay products 13 . The standard of the HepMC format [21] of the event record is a usual input format for TauSpinner, but some helping classes are prepared in the examples, to allow also for reading LHEF format [20]. These interface classes are usually kept separately from the examples code and described in the appropriate README files of the particular example directory. In several cases also input files with small number of generated events are provided together with the examples.
The technical implementation of the methods used in examples, evolved during the last few years. Most of them use basic (sometimes also additional developped during the project) functionalites of the MCTester package [22,23] and its native analysis function MC_Analyze. This allows to standarize not only production of benchmark histograms, but also validation of class reading event record, readParticlesFromTAUOLA_HepMC method. This method is expected to be adapted or replaced by the user own, depending on the particular event record form.
Let us list some advantages of MCTester: 1. MCTester interface to HepMC is profoundly tested for many conventions on how to handle HepMC event record exceptions in the way it was filled in. MCTester offers thus debugging facility for communication too.

The only information exchanged between
MCTester and TauSpinner modules is the event weight which is passed to MC_Analyze function together with the event itself MC_Analyze(temp_event,WT); 3. Particular test specific plots, are defined in C scripts with names following MCTester analysis template script UserTreeAnalysis.C. These user-defined plots are then generated by MC_Analyze function.
The TauSpinner and MCTester codes are independent, use of the latter is not mandatory, see Appendix B.8.
We are now ready to return to documentation of the TauSpinner application methods. There are three categories of methods prepared for the user: for initialization, evaluation of particular event and finally for gathering auxiliary information. Let us start from general methods. Methods for the specific applications will be discussed later, case after case in sub-appendices.

Initialization
Before TauSpinner is activated, and individual events analyzed, general information about sample on which it is going to be used has to be provided. Initialization of TauSpinner has to be adjusted to the sample specific. The initialization is divided into two categories, the general one discussed now and the one specific for particular application which will be discussed in the sub-appendices. In every case the following initialization routine has to be invoked: to intialize τ decay library. Generation of decays will not be performed but constants necessary for calculation of decay matrix elements have to be set. The library of PDFs has to be initialized as well, string name="MSTW2008nnlo90cl.LHgrid"; LHAPDF::initPDFSetByName(name); for every application and finally the basic initialization method of TauSpinner can be invoked, initialize_spinner(Ipp, Ipol, nonSM2, nonSMN, CMSENE); The input parameters of this method are as follows: • bool Ipp = true; denote that program is to be used for pp collisions. At present it is the only available option.
• int Ipol; states of the polarization in the input sample -0 -events generated without spin effects -1 -events generated with all spin effects -2 -events generated with spin correlations only, no individual τ polarization.
-3 -events generated with spin correlations and polarization but missing angular dependence of individual τ polarization.
• int nonSM2 = 0 -Are we using nonSM calculations? y/n (1/0), • int nonSMN = 0 -If we are using nonSM calculations we may consider corrections to shapes only: y/n (1/0), this option was introduced for 2 → 2 algorithm only, • double CMSENE -Center of mass system energy, used in reconstruction of x 1 , x 2 , the arguments for PDF distributions.
There are two random number generators used in typical applications. If re-decay of τ leptons by Tauola is requested, then one may need to re-seed its native random generator, especially if parallel generation is performed. For this end, the sequence int ijklin=..., int ntotin=..., int ntot2n=...; Tauola::setSeed(ijklin,ntotin,ntot2n); or Tauola::setSeed(time(NULL), 0, 0); should be considered. The TauSpinner algorithm for attribution of τ helicities use no-sophisticate random number generator. It can be replaced with the command Tauola::setRandomGenerator( randomik ); Other initialization methods, specific only for some applications will be presented in sub-appendices.

Event loop
In a loop, for each consecutive event read from the datafile 14 , the four momenta and pdgID of the outgoing τ leptons and their 15 decay products are retrieved. The intermediate state boson X has to be known too (from the content of the event record or set manually).
For each event, information can be read from datafile with the help of int status = readParticlesFromTAUOLA_HepMC(input_file, X, tau, tau2, tau_daughters, tau_daughters2); method, present in file /TauSpinner/examples/read_particles_from_TAUOLA.cxx but we assume that user will replace ti with his own, matching the format of the datafile of the actual use.

Single event analysis
On the basis of the intermediate boson X pdgID algorithm decides which matrix element version is used for calculating R i, j (and production weight as well). The four momentum of X does not need to be equal to the sum of the outgoing τ leptons (or τ lepton plus ν τ in case of W ± or H ± mediated processes). The non-conservation is expected, and will be attributed to the presence of final state emitted photons (prior to τ decays). Then, the weight defined as in Eq. (4) is calculated wt spin = calculateWeightFromParticlesH(X, tau, tau2, tau_daughters,tau_daughters2); If there is only a single τ among X decay products then it is assumed that tau2=neutrino wt spin = calculateWeightFromParticlesWorHpn(X, tau, tau2, tau_daughters); The formula Eq. (4) is used in both cases, for the second one where instead of τ, the ν τ is present, a polarimetric vector (1, 0, 0, 0) is used.
The methods prepared for 2 → 4 processes have similar input parameters, only the four momenta corresponding to outgoing partons (jets) are added. The event can be read from datafile, e.g. with the help of int status = read_particles_for_VBF(input_file,p1,p2,X,p3,p4,tau1,tau2,tau1_daughters,tau2_daughters); and later the weight is calculated wt spin = calculateWeightFromParticlesVBF(p3, p4, X, tau1, tau2, tau1_daughters, tau2_daughters); Auxiliary information 14 In most cases event record format of HepMC is used as an input. 15  As a by-product of preparing input for calculation of spin weight, Eq. (4), auxiliary physics information is available and can be retrieved. This information is available through the following get functions that can be called after each calculateWeightFrom... function call. double wtnonSM = getWtNonSM() -returns production weight wt prod given by. Eq. (11), double WTp = getWtamplitP(); double WTm = getWtamplitM() return weights wt τ + decay , wt τ − decay of Eq. (12), corresponding to variation of τ ± decay matrix elements. Finally double getTauSpin() returns the τ lepton helicity, attributed stochastically from the p z τ probability given by Eq. (20), which is calculated from production and τ decay products kinematics. Note that such attribution is often subject of a further approximation than used for spin weight calculation: transverse spin effects are then not included, and 100% spin correlation is assumed. That is why, providing helicity only for the first τ is sufficient. The getEWwt(), getEWwt0() methods are described in Appendix B.4 Before we explain further details of examples, let us list in the Table 4, all weights and their methods, and where they are adressed in Appendix.

B.1 Longitudinal spin effects in 2 → 2 processes
From the first version of TauSpinner [1], the purpose of calculating spin weight wt spin given by formula (4), was to evaluate how spin correlations between produced τ + and τ − affect different observables. Applications for the cases of Z, γ * , mW ± , H, H ± mediated processes were prepared. Longitudinal spin effects were studied only. Other than t, z components of the R i, j matrix, Eq. (4), were set to zero. Source of the paper [1] is available in directory Tauola/TauSpinner/paper. Introduced at that time methods form the core of the TauSpinner code. Later, more refined matrix elements and accompanying methods were introduced. The general principles remained unchanged. The parameters of the Born amplitudes in case of the Higgs mediated process can be modified with: void setHiggsParameters(int jak, double mass, double width, double normalization).
At the time, these parameters were not yet known experimentally. For the Drell-Yan process the defaults of Tauola [6] are used.
Detailed discussion of examples is presented with multitude plots of Ref. [3]. Option depicted by the initialization parameter Ipol necessary to adopt to the level of spin effects embeded in, to be reweighted events, is discussed.

B.2 Beyond Standard Model Higgs and Drell-Yand 2 → 2 processes
With Ref. [2] re-weighting with the help of the production matrix element was introduced and wt prod of Eq. (11) became available. The methods setNonSMkey(1); double getWtNonSM(); were introduced. The facility was prepared for handling the 'new' matrix element introduced by the user. The methods set_nonSM_born( nonSM_adopt ); set_nonSM_bornH( nonSM_adoptH ); were provided. Numerical values of the coupling constants for the default matrix elements of TauSpinner can be changed in that way as well. We will exploit that point in context of 2 → 4 processes. The user provided matrix elements can be defined for the Higgs mediated processes or for the Dell-Yan ones. The algorithm will make a choice which one to use, on the basis of pdgID of the X. With wt prod one can change the initialization parameters of the Standard Model amplitudes as well.
If no user provided matrix elemens are attributed, program will print warning and exit(-1) command will be activated.

B.3 Examples for studies of Higgs parity observables.
Example program and reference numerical results for CP sensitive observables as of Ref. [4] are given in directory Tauola/TauSpinner/examples/CP-tests. The wt spin is used. Sub-directories with names starting with H-and Z-, are respectively for cases where the signal (Higgs boson) or background (Drell-Yan, Z-production) distributions are studied. The second part of the directory name, denotes decay channels of τ leptons, pi for τ → πν and rho for τ → ππ 0 ν.
Sub-directory generate-datafiles includes the program to read and generate τ decays for events stored in input Pythia HepMC files. The README file with instructions is prepared, but was not updated now in every detail. Source of the paper [4] is in directory: TAUOLA/TauSpinner/paper-CP.

B.3.1 CP Higgs
From the algorithmic point of view introduction of CP-sensititive transverse spin correlations into R i, j matrix was straightforward and completed in Ref. [4]. The task was simple, because due to the zero spin of the Higgs there is no dependence on the Higgs production process. To activate transverse spin correlations, one should initialize one of the below: • setHiggsParametersTR(-1.0, 1.0, 0.0, 0.0); for scalar Higgs, • setHiggsParametersTR( 1.0,-1.0, 0.0, 0.0); for pseudo-scalar Higgs, • double theta = 0.2; for mixed parity state choose theta for mixing angle and then setHiggsParametersTR(-cos(2*theta),cos(2*theta), -sin(2*theta),-sin(2*theta)); For the code, see TAUOLA/TauSpinner/examples/CP-tests and its sub-directories H-pi and H-rho and for an example of results obtained with the methods Refs. [24,25].

B.3.2 Transverse spin correlations in the background; 2 → 2 Drell-Yan process
With the help of the tables of electroweak corrections, the transverse spin effects of Z/γ * mediated processes can be introduced into the sample. As default those terms are switched off. To activate as explained in Ref. [4], initialization method setZgamMultipliersTR(1., 1., 1., 1. ); has to be used. User may be interested to introduce transverse spin correlations in part, then some of these coefficients can be set to zero or even one can reverse their sign, by setting multipliers to −1. If any component of transverse density matrix for Drell-Yan has been turned on using setZgamMultipliersTR, spin weight wt spin will also take into account transverse spin correlations R i, j , i, j = x, y (calculated from tables table1-1.txt table2-2.txt which must be stored in the same directory as the program) by calling 16  If tables are present in the directory of the executable, then they will be read in and used. Then, if double wt1=Tauola::getEWwt(); double wt0=Tauola::getEWwt0(); is executed, one can also improve the distribution and normalization of the cross section with the wt1/wt0 factor added in the wt prod event weight and not only take transverse spin effects discussed above.

B.5 Drell-Yan ττ pair production with two jets
Necessary methods were introduced/discussed in Ref. [5]: it turned out of interest to extend calculation of matrix elements to the case where four-momenta of two jets accompany outgoing τ lepton pair, are also used explicitely. This option is invoked, see TAUOLA/TauSpinner/examples/example-VBF by call to WT = calculateWeightFromParticlesVBF(p3, p4, X, tau1, tau2, tau1_daughters, tau2_daughters); instead of discussed previously WT = calculateWeightFromParticlesH(...).
The matrix elements replacement necessary for the calculation of the corresponing wt prod is set by TauSpinner::set_vbfdistrModif(vbfdistrModif); and then initialized with vbfinit_(&ref,&variant); For the non-standard calculation it is then available with setNonSMkey(1);. In [5] it was found interesting to use this interface to redefine initialization of the parameters in the default matrix elements. Distinct electroweak schemes can be initialized for standard and non-standard calculations. Variants for running α S can be also distinct for standard and non-standard calculation. Note that the default for non-standard calculation is to use the same set of matrix elements as for Standard Model calculation. The electroweak scheme can be nonetheless different (the one corresponding to variant). Also option for distinct scale in calculation of PDFs is prepared. It can be set with setPDFOpt(QCDdefault,QCDvariant); respecively for the default and non-standard calculation. Also, new parametrization of running strong coupling can be introduced with this interface. The default function can be replaced with the user own: TauSpinner::set_alphasModif(alphasModif);. These options were already explored in Ref. [5].

B.6 Spin2
So far we have discussed production weight wt prod as a tool to reweight to slightly modified initialization parameters of the Standard Model production model. As demonstrated in Ref. [26], matrix elements calculated for the non-standard type (eg. for non-standard resonance) may differ substantially from the default-standard mode ones. Invoking TauSpinner::set_vbfdistrModif(SPIN2::spin2distr); will activate user function, as of our example, the SPIN2::spin2distr(...,KEY). In this case, user matrix element prototype invokes a whole library which need initialization, see TAUOLA/TauSpinner/examples/example-VBF/SPIN2. Its organization and initialization is quite similar to the code of the TauSpinner library of ττ j j processes, but it is nonetheless completely independent. The spin2init_(&ref,&variant); has the same parameters as the vbfinit_, but not all of its functionality is at this moment active. The parameter ref is dummy and all constants and coupling are independently initialized, in particular hard coded own value for not-running α s is used. The parameter variant can be used to vary electroweak scheme though.
Ambiguity of choice for the hadronic currents used in simmulations may be a source of systematic errors, which is of the interest to be estimated. An attempt for such estimation is shown in an example TAUOLA/TauSpinner/examples/CLEO-to-RCHL/tau-reweight-CLEO-to-RCHL.cxx The calculation with two versions of hadronic currents is possible without any intervention into the code of libraries and interchanged thanks to the method Tauola::setNewCurrents( int key); where key=0 or key=1 is possible. Once hadronic current model is set, the usual weight calculation can be performed and double getWtamplitP(); double getWtamplitM(); can be used to obtain matrix elements squared respectively for τ + and τ − . Ratio of the results for the two versions of hadronic currents will constitute the wt τ + decay and wt τ − decay weights.

B.8 Not published tests programs
Set of selected test programs, coded without MCTester user analysis helper, are collected in directory TAUOLA/TauSpinner/examples/testsERW.

B.9 Function sigborn
Input parameters are pdgID of incoming partons, virtuality squared of the resonance Q 2 and cos θ of the outgoing lepton, the intermediate particle identifier is not passed as an argument. Appropriate choice of the function used is performed by the TauSpinner method calculating spin and production weight. In case of the SM Higgs boson production, assumed is gg → H production process only, the formula for the lowest-order cross-section summed over all τ helicity configurations using wrapper to disth function from Tauola code read as follows (disth(Q 2 , cos θ, 1, 1) + disth(Q 2 , cos θ, 1, −1) (39) +disth(Q 2 , cos θ, −1, 1) + disth(Q 2 , cos θ, −1, −1)).
In case of Drell-Yan process, respective formula reads, using wrapper to t_born function from Tauola code.
where Q 2 is the virtuality of the resonance, cos θ is the cosine of angle between τ + and the first beam and last two arguments are helicity states of τ + and τ − respectively. The 123231 = ∼ 2 · π/α 2 QED represent normalisation factor. Details of normalization were not studied so far, this was out of scope for the studied applications, as in the wt spin ratios of σ Born were used and for wt prod required precision was sufficent, in cases when required precision was higher the ratio of σ Born was used for wt prod as well.
The electroweak parameters are initialised with the setting as in Table 5, using wrappers to initwk_ function from Tauola library code. In Appendix B.4, it is explained how to activate electroweak one loop corrections, of physics content as described in Refs. [16,27] into effective Born angular distributions. Appendices C.7 and D of Ref. [6] document how the scheme of these electroweak effects calculation can be changed. The corresponding example code is provided in TAUOLA/TauSpinner/examples/testsERW/ex4