Chapter 10 Implemented physics
10.1 The hard interaction modelsIn this section, we give a brief overview over the different incarnations of models for the description of the realm of subatomic particles and their interactions inside WHIZARD. In Sec. 10.1.1, the Standard Model (SM) itself and straightforward extensions and modifications thereof in the gauge, fermionic and Higgs sector are described. Then, Sec. 10.1.2 gives a list and short description of all genuine beyond the SM models (BSM) that are currently implemented in WHIZARD and its matrix element generator O’Mega. Additional models beyond that can be integrated and handled via the interfaces to external tools like SARAH and FeynRules, or the universal model format UFO, cf. Chap. 17. 10.1.1 The Standard Model and friends10.1.2 Beyond the Standard Model
Strongly Interacting Models and Composite ModelsHiggsless models have been studied extensively before the Higgs boson discovery at the LHC Run I in 2012 in order to detect possible loopholes in the electroweak Higgs sector discovery potential of this collider. The Threesite Higgsless Model is one of the simplest incarnations of these models, and was one of the first BSM models beyond SUSY and Little Higgs models that have been implemented in WHIZARD [39]. It is also called the Minimal Higgsless Model (MHM) [40] is a minimal deconstructed Higgsless model which contains only the first resonance in the tower of KaluzaKlein modes of a Higgsless extradimensional model. It is a nonrenormalizable, effective theory whose gauge group is an extension of the SM with an extra SU(2) gauge group. The breaking of the extended electroweak gauge symmetry is accomplished by a set of nonlinear sigma fields which represent the effects of physics at a higher scale and make the theory nonrenormalizable. The physical vector boson spectrum contains the usual photon, W^{±} and Z bosons as well as a W′^{±} and Z′ boson. Additionally, a new set of heavy fermions are introduced to accompany the new gauge group “site” which mix to form the physical eigenstates. This mixing is controlled by the small mixing parameter є_{L} which is adjusted to satisfy constraints from precision observables, such as the S parameter [41]. Here, additional weak gauge boson production at the LHC was one of the focus of the studies with WHIZARD [42]. Supersymmetric ModelsWHIZARD/O’Mega was the first multileg matrixelement/event generator to include the full Minimal Supersymmetric Standard Model (MSSM), and also the NMSSM. The SUSY implementations in WHIZARD have been extensively tested [43, 44], and have been used for many theoretical and experimental studies (some prime examples being [45, 46, 56]. Little Higgs ModelsInofficial modelsThere have been several models that have been included within the WHIZARD/O’Mega framework but never found their way into the official release series. One famous example is the noncommutative extension of the SM, the NCSM. There have been several studies, e.g. simulations on the schannel production of a Z boson at the photon collider option of the ILC [49]. Also, the production of electroweak gauge bosons at the LHC in the framework of the NCSM have been studied [50]. 10.2 The SUSY Les Houches Accord (SLHA) interfaceTo be filled in ... [52, 53, 54]. The neutralino sector deserves special attention. After diagonalization of the mass matrix expresssed in terms of the gaugino and higgsino eigenstates, the resulting mass eigenvalues may be either negative or positive. In this case, two procedures can be followed. Either the masses are rendered positive and the associated mixing matrix gets purely imaginary entries or the masses are kept signed, the mixing matrix in this case being real. According to the SLHA agreement, the second option is adopted. For a specific eigenvalue, the phase is absorbed into the definition of the relevant eigenvector, rendering the mass negative. However, WHIZARD has not yet officially tested for negative masses. For external SUSY models (cf. Chap. 17) this means, that one must be careful using a SLHA file with explicit factors of the complex unity in the mixing matrix, and on the other hand, real and positive masses for the neutralinos. For the hardcoded SUSY models, this is completely handled internally. Especially Ref. [56] discusses the details of the neutralino (and chargino) mixing matrix. 10.3 Lepton Collider Beam SpectraFor the simulation of lepton collider beam spectra there are two dedicated tools, CIRCE1 and CIRCE2 that have been written as in principle independent tools. Both attempt to describe the details of electron (and positron) beams in a realistic lepton collider environment. Due to the quest for achieving high peak luminosities at e^{+}e^{−} machines, the goal is to make the spatial extension of the beam as small as possible but keeping the area of the beam roughly constant. This is achieved by forcing the beams in the final focus into the shape of a quasi2D bunch. Due to the high charge density in that bunch, the bunch electron distribution is modified by classical electromagnetic radiation, so called beamstrahlung. The two CIRCE packages are intended to perform a simulation of this beamstrahlung and its consequences on the electron beam spectrum as realistic as possible. More details about the two packages can be found in their standalone documentations. We will discuss the basic features of leptoncollider beam simulations in the next two sections, including the technicalities of passing simulations of the machine beam setup to WHIZARD. This will be followed by a section on the simulation of photon collider spectra, included for historical reasons. 10.3.1 CIRCE1While the bunches in a linear collider cross only once, due to their small size they experience a strong beambeam effect. There is a code to simulate the impact of this effect on luminosity and background, called GuineaPig++ [10, 11, 12]. This takes into account the details of the accelerator, the final focus etc. on the structure of the beam and the main features of the resulting energy spectrum of the electrons and positrons. It offers the stateoftheart simulation of leptoncollider beam spectra as close as possible to reality. However, for many highluminosity simulations, event files produced with GuineaPig++ are usually too small, in the sense that not enough independent events are available for physics simulations. Lepton collider beam spectra do peak at the nominal beam energy (√s/2) of the collider, and feature very steeply falling tails. Such steeply falling distributions are very poorly mapped by histogrammed distributions with fixed bin widths. The main working assumption to handle such spectra are being followed within CIRCE1:
The two powers α and β are the main coefficients that can be tuned in order to describe the spectrum with CIRCE1 as close as possible as the original GuineaPig++ spectrum. More details about how CIRCE1 works and what it does can be found in its own writeup in circe1/share/doc. 10.3.2 CIRCE2The two conditions listed in 10.3.1 are too restrictive and hence insufficient to describe more complicated leptoncollider beam spectra, as they e.g. occur in the CLIC drivebeam design. Here, the two beams are highly correlated and also a powerlaw description does not give good enough precision for the tails. To deal with these problems, CIRCE2 starts with a twodimensional histogram featuring factorized, but variable bin widths in order to simulate the steep parts of the distributions. The limited statistics from too small GuineaPig++ event output files leads to correlated fluctuations that would leave strange artifacts in the distributions. To abandon them, Gaussian filters are applied to smooth out the correlated fluctuations. Here care has to be taken when going from the continuum in x momentum fraction space to the corresponding boundaries: separate smoothing procedures are being applied to the bins in the continuum region and those in the boundary in order to avoid artificial unphysical beam energy spreads. Fig. 10.1 shows the smoothing of the distribution for the bin at the x_{e+} = 1 boundary. The blue dots show the direct GuineaPig++ output comprising the fluctuations due to the low statistics. Gaussian filters with widths of 3 and 10 bins, respectively, have been applied (orange and green dots, resp.). While there is still considerable fluctuation for 3 bin width Gaussian filtering, the distribution is perfectly smooth for 10 bin width. Hence, five bin widths seem a reasonable compromise for histograms with a total of 100 bins. Note that the bins are not equidistant, but shrink with a power law towards the x_{e−} = 1 boundary on the right hand side of Fig. 10.1. WHIZARD ships (inside its subpackage CIRCE2) with prepared beam spectra ready to be used within CIRCE2 for the ILC beam spectra used in the ILC TDR [13, 14, 15, 16, 17]. These comprise the designed staging energies of 200 GeV, 230 GeV, 250 GeV, 350 GeV, and 500 GeV. Note that all of these spectra up to now do not take polarization of the original beams on the beamstrahlung into account, but are polarizationaveraged. For backwards compatibility, also the 500 GeV spectra for the TESLA design [27, 28], here both for polarized and polarizationaveraged cases, are included. Correlated spectra for CLIC staging energies like 350 GeV, 1400 GeV and 3000 GeV are not yet (as of version 2.2.4) included in the WHIZARD distribution. In the following we describe how to obtain such files with the tools included in WHIZARD(resp. CIRCE2). The procedure is equivalent to the socalled lumilinker construction used by Timothy Barklow (SLAC) together with the legacy version WHIZARD 1.95. The workflow to produce such files is to run GuineaPig++ with the following input parameters:
This demands from GuineaPig++ the generation of distributions for the e^{−}e^{+}, e^{∓}γ, and γγ components of the beamstrahlung’s spectrum, respectively. These are the files lumi.ee.out, lumi.eg.out, lumi.ge.out, and lumi.gg.out, respectively. These contain pairs (E_{1}, E_{2}) of beam energies, not fractions of the original beam energy. Huge event numbers are out in here, as GuineaPig++ will produce only a small fraction due to a very low generation efficiency. The next step is to transfer these output files from GuineaPig++ into input files used with CIRCE2. This is done by means of the tool circe_tool.opt that is installed together with the WHIZARD main binary and libraries. The user should run this executable with the following input file:
The first line defines the output file, that later can be read in into the beamstrahlung’s description of WHIZARD (cf. below). Then, in the second line the design of the collider (here: ILC for 500 GeV centerofmass energy, with the number of bins) is specified. The next line tells the tool to take the unpolarized case, then the GuineaPig++ parameters (event file and luminosity) are set. In the last three lines, details concerning the adaptation of the simulation as well as the smoothing procedure are being specified: the number of iterations in the adaptation procedure, and for the smoothing with the Gaussian filter first in the continuum and then at the two edges of the spectrum. For more details confer the documentation in the CIRCE2 subpackage. This produces the corresponding input files that can be used within WHIZARD to describe beamstrahlung for lepton colliders, using a SINDARIN input file like:
10.3.3 Photon Collider SpectraFor details confer the complete writeup of the CIRCE2 subpackage. 10.4 Transverse momentum for ISR photonsThe structure functions that describe the splitting of a beam particle into a particle pair, of which one enters the hard interaction and the other one is radiated, are defined and evaluated in the strict collinear approximation. In particular, this holds for the ISR structure function which describes the radiation of photons off a charged particle in the initial state. The ISR structure function that is used by WHIZARD is understood to be inclusive, i.e., it implicitly contains an integration over transverse momentum. This approach is to be used for computing a total cross section via integrate. In WHIZARD, it is possible to unfold this integration, as a transformation that is applied by simulate step, event by event. The resulting modified events will show a proper logarithmic momentumtransfer (Q^{2}) distribution for the radiated photons. The recoil is applied to the hardinteraction system, such that fourmomentum and √ŝ are conserved. The distribution is cut off by Q_{max}^{2} (cf. isr_q_max) for large momentum transfer, and smoothly by the parton mass (cf. isr_mass) for small momentum transfer. To activate this modification, set
before, or as an option to, the simulate command. Limitations: the current implementation of the p_{T} modification works only for the symmetric doubleISR case, i.e., both beams have to be charged particles with identical mass (e.g., e^{+}e^{−}). The mode recoil generates exactly one photon per beam, i.e., it modifies the momentum of the single collinear photon that the ISR structure function implementation produces, for each beam. (It is foreseen that further modes or options will allow to generate multiple photons. Alternatively, the PYTHIA shower can be used to simulate multiple photons radiated from the initial state.) 10.5 Transverse momentum for the EPA approximationFor the equivalentphoton approximation (EPA), which is also defined in the collinear limit, recoil momentum can be inserted into generated events in an entirely analogous way. The appropriate settings are
Limitations: as for ISR, the current implementation of the p_{T}
modification works only for the symmetric doubleEPA case. Both
incoming particles of the hard process must be photons, while both
beams must be charged particles with identical mass (e.g., e^{+}e^{−}).
Furthermore, the current implementation does not respect the
kinematical limit parameter It is possible to combine the ISR and EPA handlers, for processes where ISR is active for one of the beams, EPA for the other beam. For this scenario to work, both handler switches must be on, and both mode strings must coincide. The parameters are set separately for ISR and EPA, as described above. 10.6 Resonances and continuum10.6.1 Complete matrix elementsMany elementary physical processes are composed of contributions that can be qualified as (multiply) resonant or continuum. For instance, the amplitude for the process e^{+}e^{−}→ qq qq, evaluated at tree level in perturbation theory, contains Feynman diagrams with zero, one, or two W and Z bosons as virtual lines. If the kinematical constraints allow this, two vector bosons can become simultaneously onshell in part of phase space. To a first approximation, this situation is understood as W^{+}W^{−} or ZZ production with subsequent decay. The kinematical distributions show distinct resonances in the quarkpair spectra. Other graphs contain only one schannel W/Z boson, or none at all, such as graphs with qq production and subsequent gluon radiation, splitting into another qq pair. A WHIZARD declaration of the form
produces the full set of graphs for the selected final state, which after squaring and integrating yields the exact treelevel result for the process. The result contains all doubly and singly resonant parts, with correct resonance shapes, as well as the continuum contribution and all interference. This is, to given order in perturbation theory, the best possible approximation to the true result. 10.6.2 Processes restricted to resonancesFor an intuitive separation of a twoboson “signal” contribution, it is possible to restrict the set of graphs to a certain intermediate state. For instance, the declaration
generates an amplitude that contains only those Feynman graphs where the specified quarks are connected to a Z virtual line. The result may be understood as ZZ production with subsequent decay, where the Z resonances exhibit a BreitWigner shape. Combining this with the analogous W^{+}W^{−} restricted process, the user can generate “signal” processes. Adding both “signal” cross sections WW and ZZ will result in a reasonable approximation to the exact treelevel cross section. The amplitude misses the singleresonant and continuum contributions, and the squared amplitude misses the interference terms, however. More importantly, the restricted processes as such are not gaugeinvariant (with respect to the electroweak gauge group), and they are no longer dominant away from resonant kinematics. We therefore strongly recommend that such restricted processes are always accompanied by a cut setup that restricts the kinematics to an approximately onshell pattern for both resonances. For instance:
In this region, the gaugedependent and continuum contributions are strictly subdominant. Away from the resonance(s), the results for a restricted process are meaningless, and the full process has to be computed instead. 10.6.3 Factorized processesAnother method for obtaining the signal contribution is a proper factorization into resonance production and decay. We would have to generate a production process and two decay processes:
All three processes must be integrated. The integration results are partial decay widths and the ZZ production cross section, respectively. (Note that cut expressions in SINDARIN apply to all integrations, so make sure that no productionprocess cuts are active when integrating the decay processes.) During a later eventgeneration step, the Z decays can then be activated by declaring the Z as unstable,
and then simulating the production process
The generated events will consist of fourfermion final states, including all combinations of both decay modes. It is important to note that in this setup, the invariant uū and dd masses will be always exactly equal to the Z mass. There is no BreitWigner shape involved. However, in this approximation the results are gaugeinvariant, as there is no offshell contribution involved. For further details on factorized processes and spin correlations, cf. Sec. 5.8.2. 10.6.4 Resonance insertion in the event recordFrom the above discussion, we may conclude that it is always preferable to compute the complete process for a given final state, as long as this is computationally feasible. However, in the simulation step this approach also has a drawback. Namely, if a partonshower module (see below) is switched on, the partonshower algorithm relies on event details in order to determine the radiation pattern of gluons and further splitting. In the generated event records, the fullprocess events carry the signature of nonresonant continuum production with no intermediate resonances. The parton shower will thus start the evolution at the process energy scale, the total available energy. By contrast, for an electroweak production and decay process, the evolution should start only at the vector boson mass, m_{Z}. In effect, even though the resonant contribution of WW and ZZ constitutes the bulk of the cross section, the radiation pattern follows the dynamics of fourquark continuum production. In general, the number of radiated hadrons will be too high. To overcome this problem, there is a refinement of the process description available in WHIZARD. By modifying the process declaration to
we advise the program to produce not just the complete matrix element, but also all possible restricted matrix elements containing resonant intermediate states. This has no effect at all on the integration step, and thus on the total cross section. However, when subsequently events are generated with this setting, the program
checks, for each event, the kinematics and determines the set of potentially
resonant contributions. The criterion is whether the offshellness of a
particular wouldbe resonance is less than the resonance width multiplied by
the value of It has to be mentioned that generating the matrixelement code for all
possible resonance histories takes additional computing resources. In the
current default setup, this feature is switched off. It has to be explicitly
activated via the Also, the feature can be activated or deactivated individually for each process, such as in
If the flag is There are two additional parameters that can finetune the conditions for
resonance insertion in the event record. Firstly, the parameter
will reduce the probability for the event to be qualified as resonant by
e^{−1}= 37 % if the kinematics is offshell by four units of the width,
and by e^{−4}=2 % at eight units of the width. Beyond this point, the
setting of the Note that if an event, by the above mechanism, is identified as following a certain resonance history, the assigned color flow will be chosen to match the resonance history, not the complete matrix element. This may result in a reassignment of color flow with respect to the original partonic event. Finally, we mention the order of execution: any additional
matrix element code is compiled and linked when 10.7 Parton showers and HadronizationIn order to produce sensible events, final state QCD (and also QED) radiation has to be considered as well as the binding of strongly interacting partons into mesons and baryons. Furthermore, final state hadronic resonances undergo subsequent decays into those particles showing up in (or traversing) the detector. The latter are mostly pions, kaons, photons, electrons and muons. The physics associated with these topics can be divided into the perturbative part which is the regime of the parton shower, and the nonperturbative part which is the regime for the hadronization. WHIZARD comes with its own two different parton shower implementations, an analytic and a socalled k_{T}ordered parton shower that will be detailed in the next section. Note that in general it is not advisable to use different shower and hadronization methods, or in other words, when using shower and hadronization methods from different programs these would have to be tuned together again with the corresponding data. Parton showers are approximations to full matrix elements taking only the leading color flow into account, and neglecting all interferences between different amplitudes leading to the same exclusive final state. They rely on the QCD (and QED) splitting functions to describe the emissions of partons off other partons. This is encoded in the socalled Sudakov form factor [29]:
This gives the probability for a parton to evolve from scale t_{2} to t_{1} without any further emissions of partons. t is the evolution parameter of the shower, which can be a parton energy, an emission angle, a virtuality, a transverse momentum etc. The variable z relates the two partons after the branching, with the most common choice being the ratio of energies of the parton after and before the branching. For finalstate radiation brachings occur after the hard interaction, the evolution of the shower starts at the scale of the hard interaction, t ∼ ŝ, down to a cutoff scale t = t_{cut} that marks the transition to the nonperturbative regime of hadronization. In the spacelike evolution for the initialstate shower, the evolution is from a cutoff representing the factorization scale for the parton distribution functions (PDFs) to the inverse of the hard process scale, −ŝ. Technically, this evolution is then backwards in (shower) time [30], leading to the necessity to include the PDFs in the Sudakov factors. The main switches for the shower and hadronization which are realized as transformations on the partonic events within WHIZARD are ?allow_shower and ?allow_hadronization, which are true by default and only there for technical reasons. Next, different shower and hadronization methods can be chosen within WHIZARD:
The snippet above shows the default choices in WHIZARD namely WHIZARD’s intrinsic parton shower, but PYTHIA6 as hadronization tool. (Note that WHIZARD does not have its own hadronization module yet.) The usage of PYTHIA6 for showering and hadronization will be explained in Sec. 10.7.3, while the two different implementations of the WHIZARD homebrew parton showers are discussed in Sec. 10.7.1 and 10.7.2, respectively. 10.7.1 The k_{T}ordered parton shower10.7.2 The analytic parton shower10.7.3 Parton shower and hadronization from PYTHIA6Development of the PYTHIA6 generator for parton shower and hadronization (the Fortran version) has been discontinued by the authors several years ago. Hence, the final release of that program is frozen. This allowed to ship this final version, v6.427, with the WHIZARD distribution without the need of updating it all the time. One of the main reasons for that inclusion – besides having the standard tool for showering and hadronization for decays at hand – is to allow for backwards validation within WHIZARD particularly for the event samples generated for the development of linear collider physics: first for TESLA, JLC and NLC, and later on for the Conceptual and Technical Design Report for ILC, for the Conceptual Design Report for CLIC as well as for the Letters of Intent for the LC detectors, ILD and SiD. Usually, an external parton shower and hadronization program (PS) is steered via the transfer of event files that are given to the PS via LHE events, while the PS program then produces hadron level events, usually in HepMC format. These can then be directed towards a full or fast detector simulation program. As PYTHIA6 has been completely integrated inside the WHIZARD framework, the showered or more general hadron level events can be returned to and kept inside WHIZARD’s internal event record, and hence be used in WHIZARD’s internal event analysis. In that way, the events can be also written out in event formats that are not supported by PYTHIA6, e.g. LCIO via the output capabilities of WHIZARD. There are several switches to directly steer PYTHIA6 (the values in brackets correspond to the PYTHIA6 variables):
The values given above are the default values. The first value corresponds to the PYTHIA6 parameter PARJ(82), its squared being the minimal virtuality that is allowed for the parton shower, i.e. the crossover to the hadronization. The same parameter is used also for the WHIZARD showers. ps_fsr_lambda is the equivalent of PARP(72) and is the Λ_{QCD} for the final state shower. The corresponding variable for the initial state shower is called PARP(61) in PYTHIA6. By the next variable (MSTJ(45)), the maximal number of flavors produced in splittings in the shower is given, together with the number of active flavors in the running of α_{s}. ?ps_isr_alphas_running which corresponds to MSTP(64) in PYTHIA6 determines whether or net a running α_{s} is taken in the spacelike initial state showers. The same variable for the final state shower is MSTJ(44). For fixed α_{s}, the default value is given by ps_fixed_alpha, corresponding to PARU(111). MSTP(62) determines whether the ISR shower is angular order, i.e. whether angles are increasing towards the hard interaction. This is per default true, and set in the variable ?ps_isr_angular_ordered. The width of the distribution for the primordial (intrinsic) k_{T} distribution (which is a nonperturbative quantity) is the PYTHIA6 variable PARP(91), while in WHIZARD it is given by pythia_isr_primordial_kt_width. The next variable (PARP(93) gives the upper cutoff for that distribution, which is 5 GeV per default. For splitting in spacelike showers, there is a cutoff on the z variable named ps_isr_z_cutoff in WHIZARD. This corresponds to one minus the value of the PYTHIA6 parameter PARP(66). PARP(65), on the other hand, gives the minimal (effective) energy for a timelike or onshell emitted parton on a spacelike QCD shower, given by the SINDARIN parameter ps_isr_minenergy. Whether or not partons emitted from spacelike showers are allowed to be only onshell is given by ?ps_isr_only_onshell_emitted_partons, MSTP(63) in PYTHIA6 language. For more details confer the PYTHIA6 manual [31]. Any other nonstandard PYTHIA6 parameter can be fed into the parton shower via the string variable
Variables set here get preference over the ones set explicitly by dedicated SINDARIN commands. For example, the OPAL tune for hadronic final states can be set via:
A very common error that appears quite often when using PYTHIA6 for SUSY or any other model having a stable particle that serves as a possible Dark Matter candidate, is the following warning/error message:
In that case, PYTHIA6 gets a stable particle (here the lightest neutralino with the PDG code 1000022) handed over and does not know what to do with it. Particularly, it wants to treat it as a heavy resonance which should be decayed, but does not know how do that. After a certain number of tries (in the example abobe 10k), WHIZARD ends with a fatal error telling the user that the event transformation for the parton shower in the simulation has failed without producing a valid event. The solution to work around that problem is to let PYTHIA6 know that the neutralino (or any other DM candidate) is stable by means of
Here, 1000022 has to be replaced by the stable dark matter candidate or longlived particle in the user’s favorite model. Also note that with other options being passed to PYTHIA6 the MDCY option above has to be added to an existing $ps_PYTHIA_PYGIVE command separated by a semicolon. 10.7.4 Parton shower and hadronization from PYTHIA810.7.5 Other tools for parton shower and hadronization10.7.6 Loopinduced processesIn order to steer loopinduced processes the usage of the OLP OpenLoopsis required. Information on the interface and setting up this program can be found in Sec. 9.5 and Sec. 5.11.1. Furthermore the following settings should be observed
