whizard is hosted by Hepforge, IPPP Durham
Previous Up Next

Chapter 4  Running WHIZARD

4.1  Setting up the process list

WHIZARD is able to cope with arbitrary scattering processes allowed by the selected physical model. The program is not a library of processes, but the user can compile a list of processes he is interested in, write a file whizard.prc and put it in the subdirectory conf. If no file is provided by the user, the file whizard.prc.default (Fig. 4.1) will be used. Note that as per default CompHEP and MadGraph are disabled all processes by these generators will be skipped. Only O’Mega and trivial test matrix elements will be produced. The user may take this file as a template, insert or delete processes as desired, and save the new version as whizard.prc.1


                 # WHIZARD configuration file
                 
                 # The selected model
                 model   SM
                 
                 # Processes
                 #  Methods: chep=CompHEP, mad=MadGraph, omega=O'Mega, test=trivial)
                 #  Options: s      selected diagrams (CompHEP/MadGraph)
                 #           r      restricted intermediate state (O'Mega)
                 #           c      apply exact color algebra (O'Mega)
                 #           n:XXX  coupling order (MadGraph)
                 #           w:XXX  width scheme (O'Mega)
                 #           p      transfer polarization (test)
                 #           u      unit matrix element (test)
                 #
                 # Tag           In      Out             Method  Option
                 #=====================================================
                 ee_c            e1,E1   e1,E1           chep
                 ee_m            e1,E1   e1,E1           mad
                 ee_o            e1,E1   e1,E1           omega
                 ww_c            e1,E1   W+,W-           chep
                 ww_m            e1,E1   W+,W-           mad
                 ww_o            e1,E1   W+,W-           omega
                 zh_c            e1,E1   Z,H             chep
                 zh_m            e1,E1   Z,H             mad
                 zh_o            e1,E1   Z,H             omega
                 nnh_c           e1,E1   n1,N1,H         chep
                 nnh_m           e1,E1   n1,N1,H         mad
                 nnh_o           e1,E1   n1,N1,H         omega
                 nnbb_m          e1,E1   n1,N1,b,B       mad
                 nnbb_o          e1,E1   n1,N1,b,B       omega
                 
Figure 4.1: Default process configuration file. Only O’Mega processes will be generated per default, the others will be skipped.

4.1.1  Model selection

WHIZARD is not restricted to the Standard Model. The physics model is selected in the header of whizard.prc, respectively.

For each model, all particles, parameters, and the vertices which go into the phase space setup are defined in a single file with extension .mdl. The vertex definitions include the complete specification of Feynman rules for MadGraph and CompHEP, while the O’Mega vertices are currently hard-coded in the corresponding executable. The model files can be found in the conf/models subdirectory. Along with the master files SM.mdl, MSSM.mdl, etc., for each model XXX there is also a file parameters.XXX.omega.f90 which contains the translation to the parameter definitions used internally by O’Mega.

The distribution contains the following models (see Tab. 4.2):

  • QED: QED with three lepton generations.
  • QCD: QCD with three quark generations.
  • SM: The Standard Model with e, sinθW and MZ as independent parameters in the electroweak sector. Some remarks are in order2:
    • The u and d masses are zero, all other masses are non-zero by default and user-definable.
    • All fermions, except for the top quark, have zero width.
    • By default, the CKM matrix is the unit matrix, but for each model where it matters, a variant with nontrivial CKM matrix is available. Keep in mind that keeping the full CKM matrix may result in less efficient calculations.
    • The Higgs boson is coupled to all massive fermions, with one exception: O’Mega has no Higgs couplings to the strange quark.
  • SM_ac: The Standard Model with anomalous couplings. This is usually considered in the context of dynamical electroweak symmetry breaking, hence, by default, the Higgs mass is set to a large value. O’Mega supports anomalous trilinear and quartic gauge couplings3.
  • SM_km: This model has been included in version 1.90. It is intended for LHC or CLIC physics of Higgsless models or models without unitarizing UV completion. A K-matrix unitarization guarantees a physically sane description for LHC and CLIC energies. More details can be found in [13].
  • MSSM: The Minimal Supersymmetric Standard Model (MSSM)4. The input data may be provided in the SUSY Les Houches Accord format (see Sec. 4.4.7). Typically, such a file is generated by automatic SUSY spectrum and decay packages.

    For convenience, there are two MSSM input files provided in the conf directory: sps1a.in implements the mass spectrum of the standard SUSY point SPS1a [14] with zero particle widths; this file has been used for generating the comparison results of Ref. [12]. The file sps1ap.in corresponds to the reference point SPS1a’ as it has been defined for the SPA convention [15]; this file includes width values for all particles that are computed, in some cases, with higher-order corrections and running couplings incorporated. Be aware that, for a consistent tree-level calculation, it may be necessary to compute widths in a different convention. This caveat is discussed further in Ref. [12].

  • MSSM_Grav: This variant of the MSSM has been included in version 1.93. It contains the MSSM with an additional gravitino and is intended for the studies of GMSB models.
  • NMSSM: This is the Next-to-Minimal Supersymmetric SM, which has been contributed to WHIZARD by Felix Braam [16].
  • PSSSM: A non-minimal MSSM version motivated by an E6 → Pati-Salam GUT group. This is a very extended model which has been implemented by Daniel Wiesler.

Model typewith CKM matrixtrivial CKM
QED with e,µ,τ,γQED
QCD with d,u,s,c,b,t,gQCD
Standard ModelSM_CKMSM
SM with anomalous couplingsSM_ac
SM with K matrixSM_km
MSSMMSSM_CKMMSSM
MSSM with gravitinoMSSM_Grav
NMSSMNMSSM_CKMNMSSM
E6 SSMPSSSM
Littlest Higgs modelLittlest
Littlest Higgs model (ungauged U(1))Littlest_Eta
Simplest Little Higgs modelSimplest
Simplest Little Higgs (universal coupl.)Simplest_univ
UEDUED
SM with ZZprime
SM with graviton resonanceXdim
SUSY toy model with gravitinoGravTest
SM as templateTemplate
[SM, all unitary gauge, CompHEP modelSM_ug ]
Figure 4.2: Summary of models currently supported by WHIZARD.

There are more BSM models, whose inclusion into WHIZARD started from version 1.51 on:

  • Littlest: The Littlest Higgs model, which contains in addition to the SM a B′ (or γ′), a W′/Z′ triplet, a complex isospin-2 scalar multiplet ψ, and a t′.
  • Littlest_Eta: The variant of the Littlest Higgs model which contains a pseudoscalar η instead of the B′ vector boson.
  • Simplest: A different Little Higgs model which has recurrences of all SM fermions and new isodoublet gauge bosons.
  • Simplest_univ: A variant of the previous model with universal, but not anomaly-free, fermion assignments.
  • Xdim: The SM augmented by a spin-2 graviton as a toy model
  • GravTest: The SM augmented by a photino and a gravitino as a toy model.
  • UED: Universal Extra Dimensions.
  • Zprime: The SM with a completely generic Z′ resonance.
  • Template: This is just a copy of the SM. It might serve as a template if the user wants to implement his/her own model.

The user might wonder whether this list can be extended. Actually, it is trivial to set up, e.g., a SM version with a different input parameter scheme or different particle naming conventions, etc.: Copy SM.mdl to, e.g., SM_new.mdl, make the changes in that file, and provide copies of the corresponding files for use of the matrix element generators5. There is a model Template as described above, which is just the SM and might serve as a template to implement the user’s own model. However, for this task one needs to modify the module Template in omega-src/bundle/src/models5.ml (which requires programming in O’Caml) and then the corresponding files Template.mdl and parameters.Template.omega.f90 in conf/models. This is a rather cumbersome task. There is fortunately an easier way: using the interface of WHIZARD to the program FeynRules [17]. Here, one can write down a Mathematica file containing a Lagrangian in a very intuitive language, and then FeynRules derives the mass eigenstates, mixing and Feynman rules and imports the model into WHIZARD. Although this has been developped already for WHIZARD version 1 and is fully supported and functional for that version, we refer to the manual of version 2 for more information, or alternatively the WHIZARD Wiki page http://projects.hepforge.org/whizard/trac/wiki.

4.1.2  Process list

WHIZARD is capable of computing 2→ n scattering processes and 1→ n decay processes. In the former case, the overall result (the integral) displayed by the program will be the total cross section in fb. In the latter case, it will be the partical decay width in GeV.

For scattering processes, n=1 (single-particle production) is allowed67. Obviously, this requires at least one beam to be structured, such that a energy spectrum is available for the initial state.

Each process is defined by a single line in whizard.prc, for instance:

h_production    e1,E1   n1,N1,H         omega
z_decay         Z       e1,E1           omega
z_production    e1,E1   Z               omega

(see Fig. 4.1). The meaning of the entries within a line is as follows:

The first entry is a unique alphanumeric tag which will be used to identify the process. It can be arbitrarily chosen, but should be a valid Fortran 95 token, lower case only. In particular, one may use lower-case characters and digits as well as the underscore, where the first character is not a digit.

The second entry defines the (partonic) initial state. It consists of two particle names, separated by a comma, with no space in between. The particle names are those given in the model file mentioned above, where typically for each particle several possible names are listed. As an example, in Fig. 4.3, the Standard Model particles are displayed as they appear in SM.mdl. This file also contains alternative names, quantum number assignments, and the physical parameters of the model.


ParticleNameAntiparticleMassWidth
d-quarkdD  
u-quarkuU  
s-quarksSms 
c-quarkcCmc 
b-quarkbBmb 
t-quarktTmtopwtop
electrone1E1me 
muone2E2mmu 
tau-leptone3E3mtau 
e-neutrinon1N1  
µ-neutrinon2N2  
τ-neutrinon3N3  
gluonG   
photonA   
Z-bosonZ mZwZ
W-bosonW+W-mWwW
HiggsH mHwH
Figure 4.3: Particle names (CompHEP naming scheme) and parameters for the Standard Model.

The third entry defines the final state. It consists of one, two or more particle names, separated by commas, with no space in between.

Flavor sums8 are defined by particle names separated by colons. For instance, the line

  eeqq  e1,E1   u:d:s,U:D:S     omega

will create a sum of all processes ee+qq, with q being any of the light quarks. Similarly,

  qqmm  u:d:s,U:D:S   e2,E2     omega

creates a sum of all processes qq→ µµ+. This is especially suitable for hadronic collisions. One should note that the mirror process qq→ µµ+ is not automatically included. To get the correct results for a hadron collider, one can either double cross section values and symmetrize final-state distributions9 (for a symmetric initial state, e.g. pp) or define the process as10

  qqmm u:d:s:U:D:S,u:d:s:U:D:S e2,E2 omega

Note that in the diagnostics messages at runtime, only a single flavor assignment will be shown, which is used as a template for generating phase space parameterizations. For the matrix element evaluation and event generation, however, all valid flavor assignments are used.

To make this more readable, there is the possibility to specify aliases: The above definition is equivalent to

  alias  q  u:d:s:U:D:S
  qqmm   q,q   e2,E2    omega

Alias definitions can be arbitrarily mixed with process definitions. They are applied consecutively to all in- and out-state specifications below their definition.

The particles to be summed over must have identical spin, mass and color representation; otherwise the summation will not work. This feature is therefore useful mainly for summing over fermions with identical quantum numbers (quarks or leptons).

The order of particles in the final state can matter if PYTHIA/JETSET is used for fragmenting events. If no information on color connections is available in simulated events (for O’Mega matrix elements if the ’c’ option is not set 11), the (dominant) color configuration has to be inferred from the particle ordering, if at all possible 12. To transfer the color connections to JETSET, particles that make up a color-singlet string should be put together, beginning with the quark, possibly intermediate gluons, and ending with the antiquark. This suffices for well-defined flavor states. If flavor sums are used and more than one quark pair is present, all possible orders of particles will show up in the final state.

The fourth entry selects the method for the generation of the matrix element, while the fifth entry is either absent or consists of a string of code letters separated by commas (blanks in between are allowed). Some code letters are followed by a colon ’:’ that indicates a value. As a default there are two options for matrix elements, either the method omega using the intrinsic generator O’Mega or the method test which produces dummy matrix elements for certain consistency tests described below.

For cross checks and debugging there are also (when enabled during the configuration of WHIZARD) the methods chep and mad for using matrix elements generated by CompHEP or MadGraph, respectively. These features have become obsolete for the user at latest with version 1.50 or 1.90, when O’Mega become the main or default generator WHIZARD. For completeness reasons the details and options for these two generators will be described at the end of this subsection.

O’Mega can generate matrix elements with an arbitrary number of final-state particles, limited only by the host resources (WHIZARD has extensively been tested only for processes up to 2→ 6, some very special cases of 2→ 7 or 2→ 8 have been studied). For multiparticle processes, code generated by O’Mega is expected to be the fastest, since more redundancies are eliminated from the matrix element evaluation.

There are several options (entries in the fifth column) for the O’Mega matrix element generator:

  • r: The letter r indicates restrictions on intermediate states for the amplitude. Effectively, this also selects classes of Feynman diagrams. Example:
        nnh        e1,E1   n1,N1,H   omega
        nnh_zh     e1,E1   n1,N1,H   omega  r: 3+4~Z
      
    For more details, see Sec. 5.2.
  • w: The letter w modifies the width treatment for unstable particles. The default behavior is the expected one: the width is inserted only for s-channel propagators. This may lead to gauge-invariance issues. The option w:fudged (or abbreviated w:f) applies the ’fudge-factor’ prescription for unstable particles where all terms are collected over the respective propagator denominators (in contrast to the usual fixed-width prescription). This prescription is useful if the process includes the exchange of photons and W bosons in the t-channel, for instance:
        enw        e1,E1  n1,E1,W-    omega    w:fudged
      
    A fixed-width option can be enforced by specifying w:constant (or w:c). Finally, all widths can be forced to zero both in the s- and t-channel by the option w:zero (or w:z).
  • O’Mega comes in two versions: By default, only the leading-1/Nc color factor is generated, and no color-flow information is available. If the option flag c is given, WHIZARD uses the tensor-product feature for process definitions of O’Mega originally intended for flavor summation to treat color exactly and to generate color-flow information. Example:
        uudd_leading   u,U  d,D      omega
        uudd_full      u,U  d,D      omega   c
      
    The modified version also eliminates some redundancies in summing over flavor.

    In the current implementation, the default (leading-color) version can be more efficient (unless many flavor states need to be summed over), so it is still useful if color is not an issue. This is the reason why we decided to keep the leading-color version as the default. Note that the O’Mega generator within WHIZARD 2 takes care of the color-flow information completely by itself and delivers matrix elements that contain the full color correlations as a default.


The test method generates a constant matrix element, independent of the particles supplied for the process. There are two choices for the normalization:

  1. If the option ’u’ is given, the matrix element is unity. With zero masses and no cuts, the integral should be exactly
    σn(s) = 
    (2π)4
    s
    Φn(s) = 
    1
    16π s
     
    Φn(s)
    Φ2(s)
    ,     (1)
    where the volume of the massless n-particle phase space is given by
    Φn(s) = 
    1
    4(2π)5
     


    s
    16π2



    n−2



     
    1
    (n−1)!(n−2)!
    .     (2)
    For n≠2 the phase space volume is dimensionful, so the units of the integral are fb×GeV2(n−2).
  2. Without option (default), the result of the integration is not a cross section, but the phase space volume, normalized by the massless phase space volume (2). Thus, if all final-state particles are massless and no cuts are specified, the integral is unity regardless of the number of particles.

Cuts and nonvanishing masses will affect the result accordingly. If events are generated, the final-state particles are distributed uniformly in phase space.

Note that the phase-space integration for the test matrix element is organized in the same way as it would be for the real 2→ n process. Since such a phase space parameterization is not optimized for the constant matrix element that is supplied instead, good convergence is not guaranteed. (Setting stratified = F may be helpful here.)

The possibility to call a dummy matrix element in this way allows to histogram spectra or structure functions: Choose a trivial process such as uudd, select the test method, switch on structure functions for one (or both) beams, and generate events. The distribution of the final-state mass squared reflects the x dependence of the selected structure function.

For the test method there are two more options:

  • p:: The letter p indicates that beam polarization should be transferred to the final state unchanged. (Otherwise, it would be ignored.) This is only possible if initial and final state coincide, e.g, e+ee+e. The option is useful for testing polarized spectra and structure functions.
        ee_test_pol    e1,E1  e1,E1   test   p
      
  • u:: The letter u indicates that the matrix element is not just constant, but unity. Such code may serve as the initial setup for inserting some nontrivial matrix element by hand.




For historical and completeness reasons we mention here the technical details and available options when using CompHEP or MadGraph matrix elements.

The CompHEP version included in WHIZARD can only generate matrix elements for up to four particles in the final state. No polarization or flavor summation is possible (in contrast to the official CompHEP version), and color-flow information is not explicitly generated in the version used by WHIZARD. However, arbitrary models can be defined and the predefined models can be modified to include anomalous couplings etc. CompHEP matrix elements, evaluated by the trace method, are most efficient for simple 2→ 2 or 2→ 3 processes, while for 2→ 4 processes the other two packages typically generate faster code.

The MadGraph version inside WHIZARD has the Standard Model implemented with the exception of the trilinear or quartic Higgs couplings. The user may choose the maximal order of the strong (QCD), electromagnetic (QED), or weak (QFD) coupling for the matrix element evaluation. The implementation included in WHIZARD can generate matrix elements with up to six final-state particles. No flavor summation is possible.

There are several options in the fifth entry of the input file for the two generators: The letter s indicates that only a selection of particular Feynman diagrams is taken. Example:

    nnh        e1,E1   n1,N1,H   chep
    nnh_zh     e1,E1   n1,N1,H   chep   s
  

The setup of the necessary diagram selection file is described in Sec. 5.1.

For MadGraphonly it is possible to define the maximal order in each of the coupling constants: ‘QCD’ counts powers of gs (gluon couplings), ‘QED’ counts powers of e (photon couplings), and ‘QFD’ counts powers of the remaining electroweak couplings. As an example, the following specifications generate the QCD and the electroweak contribution to uū→ dd separately:

    uudd_qcd   u,U   d,D    mad    n:QED<=0, n:QFD<=0
    uudd_ew    u,U   d,D    mad    n:QCD<=0

There is only an upper bound on the number of couplings of a given species supported, so <= is the only available relation.

MadGraph always uses a fixed-width prescription, both in the s- and the t-channel, so there is no option. The “fudged” width option mentioned above can be activated for CompHEP by setting the flag gwidth. However, this is a runtime flag that enters the input file whizard.in in the block parameter_input.

4.2  Compilation, installation, and bundles

After the processes have been selected, the matrix elements and executable code are compiled and linked by the simple command

make prg

This will generate an executable named whizard in the bin subdirectory.

An alternative is

make bundle-src
make bundle

which will create a restricted bundle containing the Fortran 95 source code of the selected matrix elements and the libraries which have not been explicitly disabled.13 This bundle will be put as a gzipped tar file in the main directory under the name whizard-bundle-yymmdd-hhmm.tgz, where yymmdd-hhmm specifies the current date and time. The bundle can be unpacked on a different platform. There, it can be configured, compiled and run without the need for a working O’Mega (i.e. O’Caml) installation. Important: In the bundle, you have to use the commands make prg_bundle install_bundle instead of make prg install, otherwise it won’t work. (The necessary O’Mega library will be contained in the bundle, of course.) The configuration and input files (whizard.prc, whizard.in, etc.) will be copied from the conf directory of the distribution (not from the results directory). The site-specific configuration config.site will not be copied (precisely because it is site-specific). The restricted bundle has all capabilities of the WHIZARD system except for the possibility to re-generate the matrix elements.

If only the source code for the process list should be generated, say

make proc-src

The files can be found in the subdirectory processes-src. If desired, the source code can be modified before a further make call compiles and links them together.

A run of the Monte Carlo generator requires, apart from the executable, the following input and configuration files. In their initial form they are found in the conf subdirectory:

  • whizard.in containing all input data and parameters.
  • whizard.mdl containing vertex definitions used in phase space setup. If an appropriate phase space configuration file whizard.phs exists, the model file is not necessary.
  • whizard.cut1 containing a-priori cut definitions. This file may be empty.
  • whizard.cut5 containing cut and histogram definitions for the data analysis step (the 5th pass of the program). This file may be empty.

When the command

make install

is issued, all necessary configuration files, together with the whizard executable, are copied into the results subdirectory. If no user-defined files are found, default ones will be used. The contents of the results directory may be copied anywhere without affecting the functionality of the whizard executable within, leaving room for the simultaneous use of multiple WHIZARD copies with different process content. (’make install’ implies ’make prg’, so the latter is in fact redundant.)

Alternatively, one can omit the install step and copy the WHIZARD executable to any place in the system, or just leave it in the bin subdirectory. For a successful run of WHIZARD, one needs, of course, copies of the input files in the working directory.

Dealing with compilation problems

Although the authors have tried to make the program compile and run in various environments, it may happen that something goes wrong. The following problems are, at least, known:

  • The O’Caml compiler is missing on the host computer. This is not a fatal problem since the compiler is part of most Linux distributions and available via MAC OS X port. Furthermore, it is easily available via http://pauillac.inria.fr/ocaml. Note that the O'Caml programming language can be installed on virtually any UNIX and MAC OS X system, and that it is perfectly possible to do this on a user account without superuser access.
  • One of the PERL scripts which organize matrix element generation fails. This may be due to inconsistent shell return codes on some systems (IRIX, for instance). Since the OS world seems to converge towards Linux, this problem rarely occurs in reality. (Please notify the author if you really need a fix for a different OS. However, for WHIZARD 1 we basically consider Linux, MAC OS X and ALPHA as supported architectures.)
  • A matrix element generator (O’Mega, CompHEP, MadGraph) fails. The reason is probably an attempt to generate a process that is not supported or which is physically impossible. 14 If this happens, there could be some trace in the temporary directory that has been created (in /tmp if TMPDIR has not been set). Note that you should remove this directory by hand if the program terminates abnormally during process creation.
  • The Fortran 95 compiler fails. You may have found an actual WHIZARD bug, or you may have discovered a bug in the Fortran 95 compiler (this has happened more than once during the development of WHIZARD!). Make sure that you have the most recent compiler version and complain to the compiler vendor, if nothing else helps. Note however, that WHIZARD 1 explores no Fortran 03 features and most Fortran 95 are very mature nowadays).
  • The linker fails due to incompatible object file formats. The reason could be that the WHIZARD configuration has chosen a Fortran 77 compiler different from the one used in generating the precompiled Fortran 77 libraries (PYTHIA, for instance). Set the correct compiler using the environment variable F77, reconfigure and recompile.
  • The program crashes with a runtime error. This can be due to an incorrect input syntax, e.g., in a NAMELIST block. A frequent error is the specification of parameters that are misspelled or unknown to the current model. Another source of runtime errors can be incompatible external library formats (PYTHIA etc.). Otherwise, in particular if there is no reasonable error message, it may be a WHIZARD bug. Please report this to one of the email addresses below.
  • Floating underflow messages after successful program execution are usually harmless.
  • CompHEP cannot be compiled, so the executable x_comphep is missing when WHIZARD tries to generate matrix elements. This is likely due to missing C libraries needed by CompHEP. In particular, although WHIZARD will call CompHEP in command-line mode, the CompHEP executable needs the X libraries and include files installed on the host system. Note that CompHEP is just for self-testing purposes and is considered to be no longer maintained by us.

In any case, whether you could not solve a configuration or compilation problem or you found a workaround, please contact the author

kilian@physik.uni-siegen.de, ohl@physik.uni-wuerzburg.de, juergen.reuter@desy.de, christian.speckner@physik.uni-freiburg.de

so that the issue can be dealt with in a future release. (Please make sure that you are running the most recent version of WHIZARD, or at least look in the current CHANGES file whether a bug has been fixed in the meantime.)

4.3  Process selection

For running WHIZARD, the user can freely choose among the processes he has defined in the configuration file and compiled in the executable. The parameter process_id in the main input file (see below) is a string which contains the list of process tags, e.g.

  process_id = "zh, zz, zw"

If there is more than one process, integration of the cross sections will be performed in sequence. If a subsequent simulation step is requested, the events will be mixed according to their relative total cross sections. If fragmentation is enabled, one can even mix WHIZARD processes with processes generated by PYTHIA (see Sec. 4.8).

Since the integration grids are written to separate files, they can be reused later if events are to be simulated for the same or a different collection of processes, provided the input parameters and cuts are identical. This may save a lot of adaptation and integration time.

If the parameter process_id is left empty, WHIZARD will just print the list of processes, write a logfile containing the complete parameter list (including user settings), and exit.

4.4  Input data

When the whizard executable is started, it will use the parameters from the copy of the file whizard.in in the working directory. In this file, the physical input parameters and runtime switches are set. There are very many parameters, and all of them are optional. WHIZARD will try to insert sensible default values if the user has not specified a value. Of course, if no process tag is given, WHIZARD can do nothing, so it will just display the list of processes and exit.

The file whizard.in is written using the Fortran 95 NAMELIST conventions. It consists of several blocks marked by a keyword beginning with the & character and closed by a slash /. Each block contains a list of assignments of the form variable = value. Variables left out are assigned their default value. The assignments are separated by commas or newlines. Whitespace can be inserted anywhere. Comments can be inserted in this file; they begin with ! and extend up to the end of the line.

An example for the input file is shown in Fig. 4.4. In many cases, it suffices to fill the process_input block and leave all others empty. If polarization, beamstrahlung etc. are intended, the beam_input blocks must also be considered.

As an alternative to the NAMELIST format, WHIZARD supports the SUSY Les Houches Accord format for its input file. See Sec. 4.4.7 for details.


                      &process_input
                      process_id = "nnh"
                      sqrts = 500
                      luminosity = 100
                      /
                     
                      &integration_input
                      calls = 
                        1  10000
                        5  10000
                        2  20000
                      /
                     
                      &simulation_input
                      write_events = T
                      /
                     
                      &diagnostics_input /
                     
                      &parameter_input
                      Me = 0
                      Ms = 0
                      Mc = 0
                      MH = 115
                      wH = 0.3228E-02
                      /
                     
                      &beam_input 
                      /
                     
                      &beam_input 
                      /
Figure 4.4: Sample input file.

4.4.1  Files

The process_input block allows to specify a generic filename. This can be used to organize WHIZARD run data by giving specific filenames. For instance, to indicate the process energy in the filename for all output files, one may write

&process_input
 process_id = "nnh"
 sqrts = 800
 filename = "whizard_800"
/

Then, the output files will be whizard_800.out and whizard.800.nnh.out, etc. You do not need to write a separate file whizard_800.in, however: if a specific file is not found, the program will always look for the default one, in this case whizard.in.

Similarly, the working directory can be specified by the directory parameter. Of course, the first input file to be read must reside in the initial directory, otherwise it cannot be found (unless you define it on the command line, see Sec. 5.6).

The filename setting is overridden by specific filenames such as input_file, or e.g. in the block integration_input the variables read_phase_space_file or write_phase_space_file (cf. below). Such filenames will be interpreted relative to the chosen directory (if set), unless they begin with a slash, e.g. for read_model_file within the integration_input block 15,

  read_model_file = "/home/whizard/standard_model/whizard"

or with “./”, e.g.,

  read_model_file = "./whizard"

The latter form refers to the initial working directory. Note that in all cases the filename extension (here, .mdl) will be appended to the filename you specify.

It is possible to specify an additional input_file name in the process_input block. This input file will be read after the current one. New settings of parameters override old ones. In principle, any number of input files can be read consecutively and merged in this way. Reading input is finished when an input file does not exist, or neither directory nor input_file is changed in the last reading. Obviously, one should avoid loops in file reading.

All external file names used by WHIZARD can be similarly redefined. As a special rule, the filename for reading something (cuts, grids, events etc.) will automatically be set equal to the corresponding filename for writing, unless it is specified explicitly. The reason is that those files will be rewritten if the corresponding input file is not found, which could lead to accidentally overwriting files which is not intended.

4.4.2  The process_input block

These parameters should always be inspected.

ParameterValueDefaultDescription
process_idstringemptyProcess tag(s) as defined in whizard.prc. It should contain the list of processes to activate, separated by commas or blanks, enclosed in quotes.
cm_frameT/FTIf true, the c.m. frame is the lab frame, the beams are in ± z directions, and the total c.m. energy is given by sqrts. If false, the beam energies and directions must be specified below in the blocks beam_input.
sqrtsnumber0If this number is greater than the sum of the incoming particle masses, it specifies the c.m. energy of the initial state in GeV. Applies only if cm_frame is true, and is ignored for decay processes.
luminositynumber0Integrated luminosity in fb−1. A nonzero value will activate event generation. (Alternatively, the number of events can be specified below in the block simulation_input.) The luminosity value is ignored for decay processes.
polarized_beamsT/FFIf true, the helicity content of the beams must be specified below in the blocks beam_input.
structured_beamsT/FFIf true, the nature of the incoming beams must be specified below in the blocks beam_input.
beam_recoilT/FFIf true, and if structure functions (e.g., ISR) are selected, the recoil of the partons against the beam remnant (e.g., the emitted photons) is taken into account. The pT distribution is computed within the approximation valid for the emission, hence it will be accurate at low pT.
recoil_conserve_momentumT/FFApplies only if beam_recoil is set: if true, keep momentum balance between parton and recoil momenta at the expense of energy balance. If false, keep energy balance at the expense of momentum balance.
filenamestringemptyBase filename (w/o extension) to be used for all input/output files instead of the string "whizard".
directorystringemptyWorking directory for all further reading/writing of files.
input_filestringemptyIf nonempty, read the specified input file after the current one. The extension .in will be appended to the filename.
input_slha_formatT/FFIf true, assume that the next input file is in SUSY Les Houches Accord format (see Sec. 4.4.7). If false, determine the format from the first line.

4.4.3  The integration_input block

These parameters will not affect the cross section value (with the possible exception of the default cut parameters) but can improve or disprove the running time, the stability and accuracy of the result, and the efficiency of event generation.

ParameterValueDefaultDescription
calls6 numbers(yes)Array describing the number of iterations and number of calls per integration pass. Default values depend on the selected process. See below in subsection 4.6 for details.
seedintegerundefinedRandom number generator seed (integer). When omitted, the time counter will be used, resulting in a different value each run.
reset_seed_each_processT/FFReset the random number generator seed to seed not just once, but each time a process is integrated. This is useful for comparing matrix elements which should be identical, but, e.g., have been generated by different programs.
accuracy_goalnumber0Goal for the accuracy estimate (6th column in the output). When this goal is reached and the efficiency goal is either also reached or unset, further grid adaptation iterations will be skipped (Sec. 4.9).
efficiency_goalpercentage100Goal for the reweighting efficiency estimate (7th column in the output). When this goal is reached and the accuracy goal is either also reached or unset, further grid adaptation iterations will be skipped (Sec. 4.9).
time_limit_adaptationinteger0If nonzero, grid adaptation for the current process will be stopped after the specified number of minutes, and the final integration pass started (Sec. 4.9).
stratifiedT/FTUse stratified (T) / importance (F) sampling.
use_efficiencyT/FFUse efficiency (T) / accuracy (F) as the criterion for adapting the channel weights.
weights_powernumber0.25Power used for adapting the channel weights. Lower value means slower adaptation (to suppress fluctuations).
min_binsinteger3Minimal number of bins per integration dimension.
max_binsinteger20Maximal number of bins per integration dimension. This number will be used as long as there are enough sampling points, otherwise the number of bins will be decreased.
min_calls_per_bininteger10Minimal number of points per bin, integration dimension, and integration channel. If this limit cannot be satisfied, the total number of points will be increased.
min_calls_per_channelinteger0All integration channel will get at least (approximately) this number of points. Prevents channels from being dropped during adaptation.
write_gridsT/FTWrite grids to files whizard.grb (best grid) and whizard.grc (current grid), to be reused later.
write_grids_rawT/FFUse binary format for writing grids. Saves memory at the expense of portability.
write_grids_filestringemptyIf nonempty, use the specified filename for writing grids instead of the default. The file extensions are appended to string.
write_all_gridsT/FFWrite, in addition, after each iteration the current grid to file whizard.grXXX, where XXX is the iteration number.
write_all_grids_filestringemptyIf nonempty, use the specified filename for writing the extra grids instead of the default. The file extensions are appended to string.
read_gridsT/FFRead existing grids whizard.grb and whizard.grc if they have been written by a previous run. This avoids the time-consuming adaptation step. Makes sense only if no physical parameters have been changed.
read_grids_rawT/FFIf true, search first for binary grid files, then for ASCII grids. If false, do search first for ASCII..
read_grids_forceT/FFSet this to T if you want to read the grids from file even if some parameters have changed. Use with care! This may result in a program crash if the grid structures are incompatible.
read_grids_filestringemptyIf nonempty, use the specified filename for reading grids instead of the default. The file extensions are appended to string.
generate_phase_spaceT/FTGenerate a phase space configuration appropriate for the current process and write it to whizard.phx.
read_model_filestringemptyIf nonempty, read vertex definitions for phase space setup from string.mdl instead of the default whizard.mdl.
write_phase_space_filestringemptyWrite phase space configuration to string.phx instead of the default.
read_phase_spaceT/FTRead phase space configuration from whizard.phs or a previously generated file whizard.phx if possible.
read_phase_space_filestringemptyRead phase space configuration from string.phs or string.phx instead of the default.
phase_space_onlyT/FFStop the program after phase space generation.
use_equivalencesT/FTIf true, use permutation symmetry when updating grids to improve the quality of the results.
azimuthal_dependenceT/FFIf false, it is assumed that the scattering does not depend on the overall azimuthal angle. This will be automatically T if general beam polarization is switched on, therefore the user need only access this parameter in the case of azimuthal-dependent cuts.
write_phase_space_channels_filestringemptyShow phase space channels in string.ps instead of the default file whizard-channels.ps. Note that you need to call CHANNELS= string make -e channels in order to generate string.ps

The following parameters affect the details of the algorithm by which WHIZARD generates the phase space setup. Under normal circumstances, there should be no need to change them.

ParameterValueDefaultDescription
off_shell_linesinteger1Maximum number of off-shell-lines allowed for Feynman graphs which are initially taken into account for the phase space configuration. Log-enhanced (massless) propagators are not counted as off-shell.
extra_off_shell_linesinteger1Use configurations with more off-shell lines, if they happen to be maximally resonant.
splitting_depthinteger1Up to this number of branchings, a (massless) propagator will be considered as log-enhanced and mapped like a photon propagator.
exchange_linesinteger3Up to this number of t-channel propagators, a multiperipheral graph will be taken into account.
show_deleted_channelsT/FFWith extra_off_shell_lines, extra channels will be generated which are deleted if they do not contain enough resonances. With this flag, they are just commented out, so they could be manually activated.
single_off_shell_decaysT/FTWhether single-off-shell decays are relevant for the phase space configuration.
double_off_shell_decaysT/FFWhether double-off-shell decays are relevant for the phase space configuration.
single_off_shell_branchingsT/FTWhether single-off-shell branchings are relevant for the phase space configuration.
double_off_shell_branchingsT/FTWhether double-off-shell branchings are relevant for the phase space configuration.
massive_fsrT/FTWhether the radiation of a massive particle in the final state is relevant for the phase space configuration.
threshold_massnumber50A particle with a mass up to this value will be considered as massless for the purpose of phase-space setup. (But the true mass is taken into account when the particle appears as a resonant intermediate state.)
threshold_mass_tnumber200A particle with a mass up to this value will be considered as massless for the purpose of phase-space setup, when it appears as a t-channel propagator.
initial_decays_fatalT/FTAs the phase space maps cannot describe on-shell decays of beam particles properly, WHIZARD normally gives a fatal error when such configurations are encountered. This option changes this to a warning at the price of a potentially screwed phase space setup.

The following parameters affect the setup of kinematical cuts. Note that the default cuts are taken into account only if there is no appropriate entry in the user cut file whizard.cut1. However, even if such an entry exists, the cut parameters below affect the phase-space mappings, so adjusting them to order-of-magnitude may improve the convergence of the result.

ParameterValueDefaultDescription
default_jet_cutnumber10The default invariant mass cut in GeV applied to pairs of massless colored particles.
default_mass_cutnumber10The default invariant mass cut in GeV applied to pair production of massless colorless charged particles and to photon emission.
default_energy_cutnumber10The default energy cut in GeV applied to photon and gluon emission.
default_q_cutnumber10The default Q cut in GeV applied to photon and gluon exchange.
write_default_cuts_filestringfileIf nonempty, write the list of default cuts to this file (augmented by the file extension) instead of the default. Note that the settings in this file are overwritten by a user-defined cut configuration, if present.
read_cuts_filestringemptyLook for user-defined cut configurations in string.cut1 instead of whizard.cut1.
user_cut_modeinteger0Set this nonzero to activate a user-defined cut function (Sec. 5.4).
user_weight_modeinteger0Set this nonzero to activate a user-defined weight function (Sec. 5.5).

4.4.4  The simulation_input block

These parameters control event generation, hadronization and output. If the luminosity (above) is left undefined or zero, one may still generate a fixed number of events by setting n_events nonzero. By default, the events are written to file in raw format and can be reused in that way. Writing events in another format (which uses up more disk space) must be requested explicitly. However, if this has not been done, one may use WHIZARD afterwards as a translator to any file format by reading pre-generated events and immediately writing them in another format:

&integration_input  read_grids = T /
&simulation_input   read_events = T  write_events = T /

WHIZARD also contains an experimental version of a kT-ordered parton shower, written by Sebastian Schmidt. In the version it is implemented here, it is highly experimental and does work only for LHEF event formats (cf. below, this is the only one which allows for a varying number of particles in the event). A much more advanced version can be found in WHIZARD 2. For more details about the shower confer the WHIZARD 2 manual. The shower can be switched on with the shower flag.

ParameterValueDefaultDescription
n_eventsinteger0Number of (unweighted) events to generate at least, irrespective of the luminosity setting.
n_callsinteger0Number of matrix-element calls (weighted events) to execute at least, irrespective of the luminosity setting.
n_events_warmupinteger0Number of extra warmup events (see below, Sec. 4.7).
unweightedT/FTReweight events to generate an unweighted event sample.
normalize_weightT/FTIf true, normalize the event weight to unity. If false, normalize to the total cross section.
write_weightsT/FFIf unweighted=F, write weight distribution data to whizard.wgt.
write_weights_filestringemptyWrite weight distribution to string.wgt instead.
safety_factornumber1Multiply the estimate for the highest weight by this factor before starting event generation.
write_eventsT/FFWrite generated events to file whizard.evt to be used by an external analysis package.
write_events_formatinteger1The format to be used for writing events, where the file extension depends on the format (.evt for format = 1, see Sec. 4.7).
write_events_filestringemptyIf nonempty, use string as filename for writing events, where the file extension will be appended.
events_per_fileinteger0If positive, begin a new event file once the number of entries exceeds this number. The event file counter is appended to each event file name, separated with a dot (before the file extension). This feature applies only to non-binary event formats.
bytes_per_fileinteger0If positive, begin a new event file once the number of bytes in the file exceeds this number. The event file counter is appended to each event file name, separated with a dot (before the file extension). This feature applies only to non-binary event formats. See Sec. 4.7.1.
min_file_countinteger1If event files are split, use this index for the first event file. Increase the counter by one for each successive event file.
max_file_countinteger999Limit for the event file counter; if this limit is exceeded, event generation is terminated. (For weighted events only, this is an error condition since the event sample must be complete for being usable.)
write_events_rawT/FTWrite events to whizard.evx in condensed binary format, so they can be internally reused in another run.
write_events_raw_filestringemptyWrite raw events to string.evx instead.
read_eventsT/FFRead events from file whizard.evx (raw format) instead of generating them. This is equivalent to read_events_raw.
read_events_forceT/FTThis was intended to force WHIZARD to read in events from file even if some parameters have changed. However, the MD5 checksum implemented to check for parameter changes has some deficiencies. Hence, we always enforce to read in events, because this feature could otherwise not be used at all. Use with great care! (Note that this problem has been solved in WHIZARD 2.)
read_events_raw_filestringemptyRead raw events from string.evx instead.
keep_beam_remnantsT/FFKeep the beam remnants in the event record when applying structure functions. See Sec. 4.4.8.
keep_initialsT/FFKeep the beam particles and the partons which initiate the hard scattering in the event record. See Sec. 4.4.8.
guess_color_flowT/FFInfer the color flow for hadronization from the particle ordering, if it is nontrivial and not available directly.
recalculateT/FFRecalculate the matrix element value for each event of a previously generated sample. Setting this flag automatically turns on reading grids and events from file.
fragmentT/FFFragment the events depending on the value of fragmentation_method (see Sec. 4.8).
fragmentation_methodinteger1Method used for fragmentation if fragment is true: 0=no fragmentation; 1=JETSET; 2=PYTHIA; 3=user.
user_fragmentation_modeinteger0When user-defined fragmentation routines are called, this parameter may select different modes.
pythia_parametersstringemptyString to be given to PYTHIA’s pygive call before starting event generation. This allows to modify PYTHIA/JETSET properties, set particle masses, etc. The string is also available within user-defined fragmentation routines and can there be abused for different purposes.
pythia_processesstringemptyPYTHIA background processes to be simulated in addition to the WHIZARD processes: A list of integers separated by blanks, enclosed in quotation marks. Refer to the PYTHIA manual for the list of processes.
showerT/FFSwitches the internal shower on or off. As a default it is off. Note that the shower only works with LHEF event format (foramt type 6). This shower is highly experimental. Mainly intended for testing purposes.
shower_nfinteger5Number of light flavors in the shower.
shower_running_alpha_sT/FFWhether to use a running strong coupling αs or not in the shower. As a default it is fixed.
shower_alpha_snumber0.2The value of the strong coupling constant αs as used by the internal shower. The default is 0.2.
shower_lambdanumber0.29The value of the QCD scale ΛQCD as used by the internal shower. The default is 0.29 GeV.
shower_t_minnumber1.0The infrared cutoff for the evolution parameter tmin as used by the internal shower. The default is 1 GeV.
shower_mdnumber0.330The constituent d quark mass md as used by the internal shower. The default is 0.330 GeV.
shower_munumber0.330The constituent u quark mass mu as used by the internal shower. The default is 0.330 GeV.
shower_msnumber0.500The constituent s quark mass ms as used by the internal shower. The default is 0.5 GeV.
shower_mcnumber1.5The constituent c quark mass mc as used by the internal shower. The default is 1.5 GeV.
shower_mbnumber4.8The constituent b quark mass mb as used by the internal shower. The default is 4.8 GeV.

4.4.5  The diagnostics_input block

These parameters do not affect the result, but the information displayed on screen and stored in files.

ParameterValueDefaultDescription
chattinessinteger4How much information to show on screen: (0) only fatal errors, (1) and non-fatal errors, (2) and warnings, (3) and messages, (4) and results, (5) and debugging messages (if any).
catch_signalsT/FTIf the compiler supports it, try to catch external signals such as SIGINT and SIGXCPU and exit gracefully, closing files first.
time_limitinteger0If nonzero, exit gracefully after the given number of minutes has passed. This is useful to prevent an external kill within a batch environment.
warn_empty_channelT/FFIssue a warning whenever the integral within a phase space channel is zero.
screen_eventsT/FFWhether to show generated events on screen.
screen_histogramsT/FFWhether to show histograms on screen.
screen_diagnosticsT/FFWhether to repeat the input parameters on screen.
show_pythia_bannerT/FTWhether to display the PYTHIA banner page if fragmentation is enabled.
show_pythia_initializationT/FTWhether to display the PYTHIA initialization messages if fragmentation is enabled.
show_pythia_statisticsT/FTWhether to display the PYTHIA statistics summary after event generation is completed.
write_logfileT/FTWhether to write the (process-specific) output file(s) whizard.XXX.out.
write_logfile_filestringemptyUse this as the filename for the logfile.
show_inputT/FTWhether to repeat the input parameters in the logfile.
show_resultsT/FTWhether to show the integration results in namelist format in the logfile.
show_phase_spaceT/FFWhether to show the phase space configuration in the logfile.
show_cutsT/FTWhether to show the cut configuration in the logfile.
show_historiesT/FFWhether to show the individual VAMP channel histories in the logfile.
show_historyT/FTWhether to show the overall VAMP history in the logfile.
show_weightsT/FTWhether to show the weight adaptation in the logfile.
show_eventT/FFWhether to show the last event in the logfile.
show_histogramsT/FFWhether to show histograms in the logfile.
show_overflowT/FFWhether to show events beyond the first or last bin in histogram listings.
show_excessT/FTWhether to show a summary of events with weight exceeding one.
read_analysis_filestringemptyUse this (string.cut5) as the filename for the analysis setup instead of whizard.cut5
plot_widthnumber130The width in mm of the plots if online analysis is enabled.
plot_heightnumber90The height in mm of the plots if online analysis is enabled.
plot_excessT/FTIn the plots, display excess events in red.
plot_historyT/FTIf this is enabled, write a graphics driver file for displaying the integration history, i.e., the integral with error bars for each iteration. Use make history to generate the graphics file whizard-history.ps.
plot_grids_channelsstringemptyThe string is a list of phase-space channels (integers) for which the bin distribution will be histogrammed. Use make grids to generate the graphics file whizard-grids.ps.
plot_grids_logscalenumber10Use logarithmic scale for the grid plots if the bin width varies over more than this ratio.
slha_rewrite_inputT/FTIf SUSY Les Houches Accord data have been used, whether to repeat this input in the process-specific logfiles (including comments) or to rewrite it there using the data which have actually been used.
slha_ignore_errorsT/FFIf this is false, an error signaled in the SLHA input file (in the SPINFO or DCINFO block) will cause WHIZARD to stop before calculating anything. If true, such errors will be displayed, but the run continues.

4.4.6  The parameter_input block

These are the physical constants used in evaluating the matrix elements. If this block is left empty, default values (see below) will be inserted.

Which constants are actually present, and which ones are dependent or derived, depends on the physical model. Consult the model files in conf/models for the corresponding parameters. For example, the constants relevant for the Standard Model SM_CKM are shown below together with their default values:

ParameterValueDefaultDescription
GFnumber1.16639× 10−5Fermi constant
mZnumber91.1882Z-boson mass
mWnumber80.419W-boson mass
mHnumber200Higgs mass
alphasnumber0.1178Strong coupling constant αs(MZ)
menumber0.000511electron mass
mmunumber0.1057muon mass
mtaunumber1.777τ-lepton mass
msnumber0.12s-quark mass
mcnumber1.25c-quark mass
mbnumber4.2b-quark mass
mtopnumber174t-quark mass
wtopnumber1.523t-quark width
wZnumber2.443Z-boson width
wWnumber2.049W-boson width
wHnumber1.419Higgs width
vckm11number0.97383Vud
vckm12number0.2272Vus
vckm13number0.00396Vub
vckm21number-0.2271Vcd
vckm22number0.97296Vcs
vckm23number0.04221Vcb
vckm31number0.00814Vtd
vckm32number-0.04161Vts
vckm33number0.99910Vtb
khgaznumber1.000anomaly Higgs coupling K factors
khgaganumber1.000anomaly Higgs coupling K factors
khggnumber1.000anomaly Higgs coupling K factors

The dependent parameters (such as e, sinθW, Vcb etc.) will be shown in the output file whizard.out, but one should not try to reset them in the input file. The dependencies take into account the tree approximation only, such that gauge invariance is automatically respected. This may lead to unwanted side-effects in particular cases, e.g., the Higgs coupling to b quarks is proportional to the b mass, so setting mb=0 removes all Hbb subprocesses.

Setting fermion masses to zero will considerably speed up the matrix element evaluation since certain helicity combinations vanish identically (with O’Mega and MadGraph). In that case, cuts may be needed for a finite answer. However, if necessary, finite masses can be kept for all particles except the light quarks u and d.

Concerning quark masses in particular, they depend in fact on the chosen scale. The default values for the parameters include the strong coupling constant defined at the scale µ=MZ, but the quark masses evaluated at the low scale µ=2 GeV. These default values are taken from the PDG 2000 compilation 16. For a consistent scheme, running quark masses should be inserted instead wherever this is of importance. For instance, the b mass is only about 2.9 GeV at µ=MZ, thus affecting the Higgs coupling to b quarks.

Although they can be calculated from the other parameters in principle, the particle widths are retained as independent input parameters. Thus, the branching ratio of a particular resonance decay can be tuned by modifying the resonance width. When the Higgs mass is changed, the Higgs width has to be changed accordingly if Higgs decays are considered. The default values for the particle widths are given by the sum of the tree-level 1→ 2 decay channels, using the default masses and couplings as input. This choice is questionable as well. Since O’Mega(MadGraph and CompHEP, too) can only calculate cross sections at tree level, this is at least consistent; however, for each particular problem one should re-investigate the underlying assumption that 2-particle decays are dominant, and recalculate the widths accordingly. This cannot be done automatically by WHIZARD 1.

Finally, even if a parameter is defined and appears in the output, this does not necessarily mean that its value is used by all three programs CompHEP, MadGraph and O’Mega. For instance, the CKM matrix is not taken into account by the MadGraph version of the Standard Model.

For using CompHEP, there are two switches within the parameter_input block:

ParameterValueDefaultDescription
gwidthT/FFUse the gauge-invariant width prescription for unstable particles in CompHEP matrix elements. This is the so-called overall factor scheme. For O’Mega, the corresponding ’f’ option has to be specified in the process configuration file already.
rwidthT/FFUse a running-width prescription for unstable particles in CompHEP matrix elements.

4.4.7  Input in SUSY Les Houches Accord format

The SUSY Les Houches Accord (SLHA) [10, 11] describes a format to pass MSSM/NMSSM physics data between computer programs. WHIZARD accepts its input file in SLHA format. This requires that either

  • The first line of the input file begins with the string
    # SUSY Les Houches Accord
    
  • or the parameter input_slha_format in the process_input block (see above) is .true. This flag can be set on the command line, or in an input file which is read before the SLHA data.

It is possible to

  • read SLHA input after the usual input file: In the process_input block of the standard input file, set input_file to the name of the SLHA file (without the extension .in which each input file must have).
  • embed standard WHIZARD input within the SLHA file: Enclose it in a block WHIZARD_INPUT. The text of this block is read by WHIZARD only and copied verbatim into an internal file which is read after the SLHA data. Note that, according to the SLHA standard, the first character of each line within this block has to be blank.

In this way, several input files in both formats may be read consecutively, which could be useful to set defaults and override some parameters computed by a SUSY spectrum calculator. However, the recommended procedure is to supply a standard whizard.in file where the process selection, collider information, and WHIZARD switches and runtime parameters are set, while containing all physical data in a separate file written in SLHA format which is read afterwards. The latter may be produced by an automatic spectrum and decay calculator.

In the output, if SLHA data have been used, the data will be repeated. WHIZARD reorders the blocks, where the blocks that actually have been used (which may include MODSEL, SMINPUTS, MINPAR, MASS, and mixing data) come first. The parameter slha_rewrite_input in the diagnostics_input block controls whether the input data are just copied into the output (preserving all comments, but not necessarily the ordering of blocks) or rewritten by WHIZARD. In the latter case, the relevant parameters in the output are made consistent between the standard WHIZARD format and the SLHA format, where applicable. Furthermore, for each process a block PROCESS will be appended which contains the integration results in a format which follows the SLHA conventions. The block CSINFO (cross section calculator information) is added, analogous to the SPINFO and DCINFO blocks described in the standard (see Fig. 4.5).


                 ! Begin SLHA data
                # SUSY Les Houches Accord Version 1.0
                Block MODSEL  # Select model
                    1    1    # mSugra
                Block SMINPUTS   # Standard Model inputs
                    1   1.279340000e+02  # alpha_em^(-1)(MZ)
                    2   1.166370000e-05  # G_Fermi
                    3   1.172000000e-01  # alpha_s(MZ)MSbar
                    4   9.118760000e+01  # MZ(pole)
                    5   4.250000000e+00  # Mb(mb)
                    6   1.743000000e+02  # Mtop(pole)
                    7   1.777000000e+00  # Mtau(pole)
                
                # WHIZARD: The following blocks have not been used:
                Block EXTPAR          # non-universal SUSY breaking parameters
                
                # WHIZARD: The following blocks have been added:
                BLOCK PROCESS   # Result of cross section calculation
                    1          2         # Number of incoming particles
                    2         11         #   Incoming particle: e
                    2        -11         #   Incoming particle: a-e
                   11          2         # Number of outgoing particles
                   12    1000024         #   Outgoing particle: ch1+
                   12   -1000024         #   Outgoing particle: a-ch1+
                   21    1.62815594E+02  # Cross section in fb
                   22    8.00542154E-02  # Cross section error in fb
                   23    1.56582525E+00  # Chi-squared value
                #
                BLOCK CSINFO   # Cross section calculator information
                    1    WHIZARD    # Calculator name
                    2    1.97       # Calculator version
                #
                 ! End SLHA data
Figure 4.5: Sample output in SLHA format

WHIZARD will use the MASS and DECAY information in the SLHA file to set the physical masses and widths of the particles it recognizes. Branching ratio information is not used. Note that the masses of some SM particles (Z, t, b, and τ) are defined in SMINPUTS. For completeness, they may also be set in the MASS block (where the W mass is located). If the values disagree, the MASS information takes precedence for WHIZARD’s needs. The widths of Z, W, and t may be set using the DECAY format. In all cases, default values (which may have been set in a previous input file) will be taken without warning where data are missing.

While the SLHA is intended for the MSSM and its extensions, it is possible to use this format to pass physical parameters for a SM calculation, too. If the SM is chosen when configuring WHIZARD, only the blocks MODSEL, SMINPUTS, MASS and DECAY information will be used. By convention, we define the model parameter (in the MODSEL block) to be negative for a non-SUSY model.

4.4.8  The beam_input blocks

The input file is finished by two blocks which describe the properties of the first and second incoming particle beam, respectively. It depends on the master flags in the process_input block above whether particular information within these blocks is actually used:

  • If cm_frame is false, the beam_input blocks specify information on the beam energies and directions. If cm_frame is true (default), this information is discarded, the beam energies are given by the sqrts value, and the beam directions are assumed along the positive and negative z axis, respectively.
  • If polarized_beams is true, the beam_input blocks specify information on the beam polarization. If polarized_beams is false (default), this information is discarded and the beams are assumed to be unpolarized.
  • If structured_beams is true, the beam_input blocks specify information on the beam particle types, beam energy spectra and/or structure functions. If structured_beams is false (default), the beam particles are assumed to be the incoming particles of the hard process, and the partonic energies are set equal to the beam energies.

If all of the master flags are kept at their default values, the content of the beam_input blocks is completely discarded, so they may be left empty.

For decay processes, there are some obvious restrictions: There is only one “beam”, so the second block is ignored. While cm_frame and polarized_beams retain their meaning (i.e., one can simulate the decay of a moving and/or polarized particle), structured_beams is ignored, and no structure functions can be selected.

Energy and direction

If the master flag cm_frame is false, the beam energies and directions are set by

ParameterValueDefaultDescription
energynumber0If greater than the beam particle mass, this specifies the beam energy in the lab frame. Otherwise, the beam energy is set equal to the particle mass (fixed target).
anglenumber0If direction is not set, this specifies a rotation of the beam axis in the lab frame around the positive y axis. (By default, the beam directions are along the positive/negative z axis, so a rotation by the angle π/2 turns them into the positive/negative x axis.) If direction is set, this parameter is ignored.
directionthree numbers0 0 0If any component is nonzero, this vector explicitly specifies the direction of the given beam in the lab frame.

Polarization

If the master flag polarized_beams is true, the polarization is set for each beam by

ParameterValueDefaultDescription
vector_polarizationT/FFIf false (default), use the standard helicity basis (left-/right-handed). Set this flag if you need another basis, in particular transversal polarization.
polarizationtwo (three) numbers0 0 0Fraction of left/right polarization (fermions, photons, gluons), resp. left/longitudinal/right polarization (massive vector bosons). If the vector polarization model is selected, the three numbers denote the polarization vector.

For fermions, the default (helicity basis) polarization model describes longitudinally polarized beams with the specified fraction of left- and right-handed particles. In the case of massless vector bosons, this corresponds to left and right circular polarization along the beam axis. Massive vector bosons have three entries, where the middle entry specifies the fraction of longitudinal polarization. If the numbers sum up to less than 1, the remainder is filled with unpolarized particles. For instance, .5 0 and .75 .25 both stand for 50 % left-handed polarized electrons.

The quantization axis of a polarized particle is the momentum vector in the lab frame. For particles at rest, the quantization axis is the positive z axis.

The vector polarization model currently applies to fermions only. The three numbers denote the components of a three-vector P with length less or equal to 1. The helicity density matrix for fermions is given by

ρ = 1/2(1 + 
P
 
·
σ
 
),     (3)

where the polarization vector has components (Px,Py,Pz). For antifermions this is replaced by

ρ = 1/2(1 + 
P
 
·
σ
 
*),     (4)

so the y component changes sign.

While the z component of P describes longitudinal polarization, the x and y components correspond to transversal polarization. Longitudinally polarized states can therefore be specified in both polarization models, while transversal polarization requires the vector model.

Spectra and structure functions

If the master flag structured_beams is true, WHIZARD will check if any of the structure functions defined below is applicable for this beam. For each beam, the beam particle is defined by either one of

ParameterValueDefaultDescription
particle_codeinteger0PDG code of the incoming beam particle.
particle_namestring"UNKNOWN"Name of the incoming beam particle.

If both entries are specified, the given PDG code takes precedence.

The further parameters specify which structure functions are used for each beam and fix structure function parameters where necessary. Structure functions are applied consecutively (Tab. 4.1), choosing a sensible chain of beam particle splittings: For instance, if the beam particle is an electron and the process is initiated by a quark, switching on both EPA and PDF will take the quark content of a photon which is emitted (in the Weizsäcker-Williams approximation) by the incoming electron. WHIZARD will not check whether the selected combination makes physically sense.

If it is desired to have the beam remnants in the event record, one should set the parameter keep_beam_remnants in the simulation_input block. The beam remnants will be photons (beamstrahlung, ISR), charged particles (EPA) or undefined hadron remnants (LHAPDF,PDF). A beam remnant is always massless, regardless which mass the particle has been assigned somewhere else, and it will have no transverse momentum. There will be exactly one beam remnant per structure function even though the number of photons radiated in beamstrahlung and ISR is not well-defined in reality. Beam remnants cannot be used for the built-in analysis since they are assigned no binary codes. Their main virtue is to make possible a reconstruction of individual splittings in the initial state in an external analysis.

Furthermore, one can set keep_initials. In that case, the beam particles and the partons which initiate the hard scattering are kept in the beam record. Together with the beam remnants, this allows for full mother-daughter relationships to be displayed. For fixed-energy beams, there are no beam remnants, and the initial partons are identical with the beam particles.


Structure functionin-particleout-particleuser beam pol.
File eventsany= inkeep
User spectrum (Sec. 5.3)any anykeep/ignore
CIRCE (beamstrahlung)e±e±keep
 e±γignore
CIRCE2 (beamstr.,γ spectra)e±e±ignore
ISR (photon radiation)charged particleunchanged keep
EPA (effective photon approximation)charged particleγ ignore
EWA (effective W approximation)weak doubletweak doublet ignore
LHAPDF (parton distribution functions)p,p or γq or gignore
PDF (parton distribution functions)p,p or γq or gignore
User structure function (Sec. 5.3)any anykeep/ignore
Table 4.1: Available structure functions. Each transforms an in-particle into a collinear out particle with a certain momentum fraction. They may be concatenated in the order given in the table. In the present implementation, the user-specified beam polarization is either kept unchanged for all events, or is overwritten.

Events from file

WHIZARD is able to read a sample of events from file and use it for integration and event generation. If both beams are read from file in this way and the filename is identical for both beams, the program expects a pair of real numbers in each line of the file. (Otherwise, it expects separate files with a single real number per line.) Each number corresponds to the energy carried by the beam particle. Empty lines and comment lines (marked by # or !) are skipped.

When the end of the beam event file is reached, a warning is issued and the file is reread from the beginning. Note that reusing the event sample in this way leads to systematic errors in the calculation. For instance, the fluctuations of the computed cross section between subsequent iterations can exceed the error estimate by orders of magnitude. In the final integration step, this could be visible as a χ2 value which is anomalously large. To avoid this problem, it is recommended to use an event sample which is statistically distributed and contains a sufficient number of events so that no rewinding is necessary.

ParameterValueDefaultDescription
FILE_events_onT/FFWhether to read the beam events from file.
FILE_events_filestringemptyFile name. The directory prefix (as specified by directory above, if any) will be prepended and the extension .bev will be appended to this string. If this file does not exist, the default file whizard.bev is searched instead.
FILE_events_energiesT/FTWhether the numbers in the file should be interpreted as absolute energies (T) or as energy fractions with respect to the nominal beam energy, i.e., √s/2 in the c.m. frame (F).

Parton distribution functions

There are two possibilities two access parton distribution functions inside hadrons: either via the LHAPDF library or via the PDFLIB inside the CERNLIB. Meanwhile we consider LHAPDF to be the default option, while we mention the PDFLIB setup just for backwards compatiblity and completeness. Distribution functions exist for all quarks and for the gluon. For LHAPDF, the parameters LHAPDF_on, PDF_scale (if PDF_running_scale is false) and LHAPDF_file for the PDF set are mandatory. For an overview over available sets and parameters confer the LHAPDF manual.

ParameterValueDefaultDescription
LHAPDF_onT/FFWhether to use PDFs from LHAPDF.
LHAPDF_filestringemptyFile name for PDF set, e.g. cteq6ll.LHpdf
LHAPDF_setinteger0Number of the structure function set within the given file. Usually, 0 is the central value, if values other than zero are available they correspond to error sets.
PDF_running_scaleT/FFWhether to use a running scale (the total partonic invariant mass) for the structure function evaluation.
PDF_scalenumber0If PDF_running_scale is false, this parameter must be set to a fixed energy scale for the structure function evaluation.

For PDFLIB inside the CERNLIB setup, the parameters are listed below (for a more detailed description, refer to the PDFLIB manual). For the numerical parameters, a zero value will cause the PDFLIB default value to be taken, so only PDF_on and PDF_scale (if PDF_running_scale is false) are mandatory.

ParameterValueDefaultDescription
PDF_onT/FFWhether to use PDFs.
PDF_ngroupinteger0Number of author group. If left 0, the PDFLIB default value (MRS) will be used.
PDF_nsetinteger0Number of the structure function set by the given group. If left 0, the PDFLIB default value will be used.
PDF_nflinteger0Number of active flavors.
PDF_lointeger0Order of αs calculation.
PDF_running_scaleT/FFWhether to use a running scale (the total partonic invariant mass) for the structure function evaluation.
PDF_scalenumber0If PDF_running_scale is false, this parameter must be set to a fixed energy scale for the structure function evaluation.
PDF_QCDL4number0QCD scale in GeV for four active flavors.

Beamstrahlung (CIRCE1)

For beamstrahlung, the accelerator type should always be defined, and it is not recommended to use an energy different from the values for which valid CIRCE [6] parameterizations exist. There is the option to either use the generator mode of CIRCE, or to accept the distributions calculated by CIRCE as analytical spectra. For efficiency reasons, the generator mode is typically preferred. An example input file whizard.in.circe_isr for the default process ee_o can be found the conf directory.

ParameterValueDefaultDescription
CIRCE_onT/FFWhether to apply beamstrahlung (electron, positron or photon beam).
CIRCE_generateT/FTWhether to use the CIRCE event generator or the analytical beamstrahlung spectra.
CIRCE_mapT/FTWhether to apply a mapping to improve convergence (relevant for analytical spectra only).
CIRCE_sqrtsnumbersThe reference energy. Note that CIRCE spectra are defined only for the discrete energy values (cf. the CIRCE write-up for more info) .
CIRCE_verinteger0The CIRCE version number required.
CIRCE_revinteger0The CIRCE revision number required.
CIRCE_accinteger0The accelerator type as needed by CIRCE.
CIRCE_chatinteger0The level of CIRCE diagnostics.

In general, more information can be found in the CIRCE manual on the HepForge page:

http://projects.hepforge.org/whizard/circe1.pdf

Photon spectra (CIRCE2)

CIRCE2 is intended as a replacement for CIRCE1. While the previous implementation provided smooth spectra, the new version, which has to be used for photon collider spectra where no CIRCE1 version is available, provides histogrammed spectra which are read from a data file. These data files are partially contained in the WHIZARD distribution, in the subdirectory circe2-src/data. For more information on the CIRCE2 generator and possible updates on the spectra confer the manual/documented source code available from the HepForge page:

http://projects.hepforge.org/whizard/circe2.pdf

There is the option to either use the generator mode of CIRCE2, or to accept the distributions calculated by CIRCE2 as (histogrammed) spectra. For efficiency reasons, the generator mode is typically preferred.

The following parameters can be defined (identical for both beams):

ParameterValueDefaultDescription
CIRCE2_onT/FFWhether to use CIRCE2 for simulating spectra.
CIRCE2_generateT/FTWhether to use the CIRCE2 generator or the (histogrammed) spectra.
CIRCE2_verboseT/FTWhether to print CIRCE2 messages on screen.
CIRCE2_filestringemptyThe data file to read by CIRCE2.
CIRCE2_designstring"∗"The accelerator design for CIRCE2.
CIRCE2_polarizedT/FTUse a polarized (T) or polarization-averaged (F) spectrum in the datafile.
CIRCE2_map-1/0/1-1Apply a power mapping for a singularity at x=0 or x=1, or no mapping (-1). Relevant only for histogrammed spectra, irrelevant for the CIRCE2 generator mode.
CIRCE2_powernumber2Parameter for the mapping, if any.

If CIRCE2 is switched on, the CIRCE2_file parameter must always be given. The filename will be interpreted relative to the current working directory. The accelerator design should be specified if the data file contains more than one dataset. (The default matches the last dataset within the data file.)

It is important to note that polarization treatment is different from the standard WHIZARD conventions when CIRCE2 is used. The CIRCE2 data files introduce fixed polarization distributions of the beams which originate from realistic simulations. Therefore, the polarized_beams switch and the polarization settings should not be used together with CIRCE2 under normal circumstances. Instead, the appropriate data file with polarization should be taken and CIRCE2_polarized switched on. This will generate the appropriate polarized spectra, properly normalized.

On the other hand, if a polarization-averaged dataset is specified, CIRCE2_polarized should also be switched off.

Despite the fact that the polarization distributions are a fixed property of any given dataset, one can reweight the polarization contributions by hand. To this effect one can set the parameters polarized_beams and polarization. For instance, by switching on only one helicity combination, the properties of the corresponding cross section component can be studied. This feature should be used with care, since the results may be unphysical.

In the spectra mode (CIRCE2_generate=false), if no mapping is specified by the user, the program inserts a mapping at x=0 (power = 3) for photons, and a mapping at x=1 (power = 12) for electrons. This assumes that all particles originate from an incident e beam, where for electrons there is a strong singularity (in fact, a smeared-out delta function) at x=1. For a different setup, e.g., a positron beam, the mapping parameters have to be adjusted by hand. Note that without the x=1 mapping, the delta-function singularity for electrons will almost certainly be missed by the adaptive integration, and the result may stabilize for a wrong cross section.

A sample input file, whizard.in.photoncollider, for the assumed hard process γγ → e e+ can be found in the conf subdirectory.

ISR spectrum

The ISR spectrum for charged particles accounts for the leading-logarithmic effect of multiple photon emission. In the leading-logarithmic approximation, universal terms of order єn with є=α/πlns/m2 can be absorbed in a structure function. The leading singularity at x=1 can be calculated exactly, while the other terms of the structure function are known only to finite order.

The implementation in WHIZARD takes into account the LLA structure function ([7]), including terms up to order єn, where n=0,1,2,3, depending on the value of ISR_LLA_order. The default value is n=3. The scale which enters the ISR spectrum is given by ISR_Q_max; the default is to take the collider c.m. energy.

Since αQED is usually calculated from GF as far as the hard scattering is concerned, it makes sense to reset it to 1/137 here. The incoming particle mass (e.g, me) must be defined here only if it has been set to zero in the parameter input block (e.g., to speed up the calculation).

ParameterValueDefaultDescription
ISR_onT/FFWhether to apply ISR (electron or positron beam).
ISR_alphanumbere2/4πThe value of αQED to be used for the spectrum.
ISR_m_innumbermThe mass of the incoming particle.
ISR_Q_maxnumbersThe hard scale which cuts off photon radiation.
ISR_LLA_order0/1/2/33The order of the leading-logarithmic approximation.
ISR_mapT/FTWhether to use a mapping of the singularity at x=1 when evaluating the structure function (recommended; note that switching this off might even lead to an uncaught arithmetic exception).

You can find an example, whizard.in.circe_isr for the sample process ee_o in the conf subdirectory.

EPA spectrum

The EPA spectrum is the Weizsäcker-Williams approximation for a real photon radiated off the incoming beam. Here, Qmax and mX must be supplied, the other parameters are optional. (If the produced system is massive, mX should be set equal or less to the mass sum of all produced particles.)

ParameterValueDefaultDescription
EPA_onT/FFWhether to use the EPA spectrum (photon beam).
EPA_mapT/FTWhether to apply a mapping to improve convergence.
EPA_alphanumbere2/4πThe value of αQED to be used for the spectrum.
EPA_m_innumbermThe mass of the incoming beam particle.
EPA_mXnumber0The lower cutoff for the produced invariant mass.
EPA_Q_maxnumber0The upper cutoff on the virtuality of the photon (Qmax>0).
EPA_x0number0The lower cutoff on the energy fraction of the incoming photon
EPA_x1number1The upper cutoff on the energy fraction of the incoming photon

Again, an example, whizard.in.epa for the sample process γγ → e e+ named aaee can be found in the conf subdirectory.

EWA spectrum

The EWA spectrum is equivalent of the Weizsäcker-Williams approximation for massive W and Z bosons radiated off the incoming beam. For details about the EWA and its implementation into WHIZARD confer [13].

ParameterValueDefaultDescription
EWA_onT/FFWhether to use the EWA spectrum (W/Z beam).
EWA_mapT/FTWhether to apply a mapping to improve convergence.
EWA_pT_maxnumber0The upper cutoff on the pT of the W/Z (pT,max>0).
EWA_x0number0The lower cutoff on the energy fraction of the incoming W/Z
EWA_x1number1The upper cutoff on the energy fraction of the incoming W/Z

Again, an example, whizard.in.ewa for the sample process W+WZZ named wwzz in the leptonic setup of a 3 TeV CLIC collider can be found in the conf subdirectory.

4.5  Cuts

WHIZARD allows for cuts in some standard kinematical variables. Cuts may be necessary in case the matrix element is infinite without them (which happens if massless particles are either exchanged or produced). In this case, WHIZARD will apply default cuts on the invariant mass of colored or charged particle pairs, on the energy of emitted photons or gluons, and on the momentum transfer to exchanged photons or gluons. The values of those cuts are controlled by parameters in the input file (see above).


! e-  e+  ->  e-  e+  gamma
! 16   8       1   2    4
process eeg
  cut Q of 10 within -99999 -1
  cut Q of 17 within -99999 -1
  cut M of  3 within 10 99999
  cut E of  4 within  5 99999
  cut PT of 4 within 19 99999
  cut THETA(DEG) of  4 1 within 5 180
  cut THETA(DEG) of  4 2 within 5 180
Figure 4.6: Cut configuration file

Alternatively, the user may specify his own set of cuts by supplying an entry in the file whizard.cut1. If such an entry exists, the program will ignore the default cuts, and the user is responsible for setting up his cuts in such a way that the cross section is still finite. The format of such an entry for the process ee+e e+ γ is shown in Fig. 4.6. The contents of whizard.cut0, where the default cuts are logged by WHIZARD (if any), may be used as a template.

The entry is begun by the keyword process followed by the process tag (eeg in this case). If the same set of cuts should be applied to several processes, one may specify more than one tag per entry, e.g.

process uudd uuss uucc
    cut ...

The cut configuration file may contain any number of process entries. Comments are marked by the letter ! or #.


CodeAlternative code(s)# ArgsDescription
- 0−2No cut
MQ1(Signed) invariant mass M=sgn(p2)√|p2|
LMLQ1log10 |M|
MSQQSQ S T U1Squared invariant mass M2 = p2
E 1Energy in the lab frame
LE 1log10E
PT 1Transverse momentum p
LPT 1log10 p
PL 1Longitudinal momentum pL
P 1Absolute value of momentum |p|
YRAP RAPIDITY1Rapidity y
ETA 1Pseudorapidity η
DETADELTA-ETA2Pseudorapidity distance Δη
PHPHI1Azimuthal angle φ (lab frame) in radians
PHDPHID PHI(DEG)1Azimuthal angle φ (lab frame) in degrees
DPHDPHI DELTA-PHI2Azimuthal distance Δφ (lab frame) in radians
DPHDDPHID DELTA-PHI(DEG)2Azimuthal distance Δφ (lab frame) in degrees
AAANGLE-ABS TH-ABS THETA-ABS1Absolute polar angle θabs (lab frame) in radians. Reference axis is the z-axis.
AADANGLE-ABS(DEG)  
 TH-ABS(DEG) THETA-ABS(DEG)1Absolute polar angle θabs (lab frame) in degrees
CTACOS(TH-ABS) COS(THETA-ABS)1cosθabs
AANGLE TH THETA2Relative polar angle θ (lab frame) in radians
ADANGLE(DEG)  
 TH(DEG) THETA(DEG)2Relative polar angle θ (lab frame) in degrees
CTCOS(TH) COS(THETA)2cosθ
A*ANGLE* TH* THETA*2Relative polar angle θ* (rest frame of part.#2) in radians
AD*ANGLE*(DEG)  
 TH*(DEG) THETA*(DEG)2Relative polar angle θ* (rest frame of part.#2) in degrees
CT*COS(TH*) COS(THETA*)2cosθ*
DRDELTA-R CONE2Distance in η-φ space, i.e. √Δη2 + Δφ2
LDRLOG-DELTA-R LOG-CONE2log10Δη2 + Δφ2
VBFVBF_ETA HEMI HEMISPHERE2η1 × η2
Table 4.2: Cut and histogram code letters. Masses, energies and momenta are measured in GeV.

For a given process, each cut is indicated by the keyword cut, followed by a code letter, the keyword of, one or two binary codes, and cut window specifications. The code letters are listed in Table 4.2. The binary codes indicate the combination of particle momenta to which the cut should be applied. As the comment lines in Fig. 4.6 suggest, the outgoing particles are given codes 1,2,4,…. Two additional codes are assigned to the incoming particles, where the first incoming particle has the highest number. Momentum sums (particle combinations) are defined by adding binary codes, such that the code 5, for instance, stands for the sum of momenta of particles 1 and 4.

The initial momenta are always counted negative, such that (in Fig. 4.6) the number 17 stands for the difference of the momenta of particle 1 (outgoing) and 16 (incoming). The same applies to angles: The polar angle between 1 and 16 is not A of 1 16 but A of 1 8 since the initial particle 16 goes into the direction of −p(8).

By default, the inserted momenta are those of the elementary partonic process. Thus, one can assume momentum conservation. For instance, in a 2→ 2 scattering process, the invariant mass of 1+2=3 (the outgoing particles) is the same as the invariant mass of 4+8=12 (the incoming particles). In some circumstances, one may need the incoming beam momenta instead. To insert those in place of the incoming parton momenta, use a negative number. In the previous example, the invariant mass of −12 would be the invariant mass of the incoming beams which is equal to the nominal c.m. energy. Of course, there is no difference unless the process is a 2→ n scattering process with structure functions enabled.

Cut window specifications consist of the keyword within and two real numbers, optionally followed by additional specifications like

  cut PT of  4 within 0 150
  cut M  of  3 within 80 100 or 180 200 or 500 99999

The absence of an upper bound should be marked by a number outside the kinematical range.

4.6  Integration

4.6.1  Running WHIZARD

After the input parameters and cuts have been specified, the integration can be started by

make run

or, if we are in the working directory (results),

./whizard

After each iteration, the result will be displayed on screen and logged in the file whizard.XXX.out, where XXX is the process tag. An example is shown in Fig. 4.7. As it is expected, the integration converges towards a stable result, and at the same time the estimated reweighting efficiency increases (allowing for some fluctuations).


         ! WHIZARD 1.97 (May 31 2011)
         ! Reading process data from file whizard.in
         ! Wrote whizard.out
         !
         ! Process nnbb:
         !    e a-e -> nu_e a-nu_e   b a-b
         !   32  16 ->    1      2   4   8
        ....

         ! WHIZARD run for process nnbb:
         !=============================================================================
         ! It      Calls  Integral[fb]  Error[fb]   Err[%]    Acc  Eff[%]   Chi2 N[It]
         !-----------------------------------------------------------------------------
         ! Reading cut configuration data from file whizard.cut1
         ! No cut data found for process nnbb
         ! Preparing (fixed weights):   2 samples of      10000 calls ...
            1      20000  5.6761132E+01  1.07E+00    1.88    2.65*  1.45    1.52   2
         !-----------------------------------------------------------------------------
         ! Adapting (variable wgts.):  10 samples of      20000 calls ...
            2      20000  5.9620359E+01  8.15E-01    1.37    1.93*  6.58
            3      20000  5.8344071E+01  4.62E-01    0.79    1.12*  8.99
            4      20000  5.8023577E+01  4.17E-01    0.72    1.02*  9.04
            5      20000  5.7938879E+01  3.83E-01    0.66    0.94* 11.95
            6      20000  5.8152958E+01  3.73E-01    0.64    0.91* 11.85
            7      20000  5.8132996E+01  3.66E-01    0.63    0.89* 12.02
            8      20000  5.8252891E+01  3.50E-01    0.60    0.85* 13.61
            9      20000  5.8296069E+01  3.42E-01    0.59    0.83* 12.74
           10      20000  5.8391405E+01  3.32E-01    0.57    0.81* 15.93
           11      20000  5.8168510E+01  3.19E-01    0.55    0.78* 16.66
         !-----------------------------------------------------------------------------
         ! Integrating (fixed wgts.):   3 samples of      50000 calls ...
           12     150000  5.8232347E+01  1.17E-01    0.20    0.78  10.11    0.49   3
         !-----------------------------------------------------------------------------
Figure 4.7: Sample output of running WHIZARD.

The number of iterations and the number of calls per iterations are as given in the input file whizard.in (see Fig. 4.4). Adaptation proceeds in three passes which are separated by a horizontal line in the output:

In the first integration pass, the grids are adapted after each iteration, but the relative weight of the integration channels is kept fixed. In the example of Fig. 4.7, this pass consists of two iterations with 10,000 calls each. This pass is useful in particular for many-particle processes, where it takes some time before VAMP has found the physical region in each channel.

In the second pass, the grids within each phase space channel and the relative weights of the phase space channels are adapted after each iteration. The example of Fig. 4.7 shows how the integral converges, the error estimate improves, and the estimated reweighting efficiency increases, until some saturation is reached. Here, the sixth column is of main interest, which is given by the estimator for the relative error multiplied by the square root of the number of calls. One expects a number of the order of unity (or smaller) here for a well-adapted grid, independent of the number of calls. Each time this number improves, it is marked by an asterisk. The best grid will be used later for the final integration and event generation steps.

In the third pass, the relative weights are fixed again, and the remaining iterations are averaged for the final cross-section result.

4.6.2  The phase space configuration

The critical point for the integration is the phase space configuration. This is set up automatically by WHIZARD before starting the integration. When integrating a given process for the first time, a phase space configuration is generated and written to the file whizard.phx. In subsequent runs for the same process, this file is reused. If you would like to edit and modify this file or if you want to set up your own configuration, the new configuration should be stored in the file whizard.phs which takes precedence over whizard.phx and will not get overwritten. The automatically generated file whizard.phx will be deleted when the process list changes or, in any case, by a make clean command.

For each process, WHIZARD phase space consists of a set of channels which formally correspond to Feynman diagrams. In the automatically generated phase space configuration file whizard.phx, these channels are listed in order, using binary codes for the internal lines. The channels that have been generated can be drawn and viewed on screen or printed:

make channels

This generates a LATEX file whizard-channels.tex which uses the FeynMF package for drawing Feynman graphs, and transforms it into a PostScript file whizard-channels.ps. These files are located in the results directory. Note that this is currently available only for the automatically generated phase space, not for a phase space setup supplied by the user. (Note further that the Feynman diagrams drawn are only sample diagrams for the corresponding phase space channels. Especially if the user edits phase space by removing channels from the file whizard.phx by hand, still the full matrix element is taken into account!)

Similarly, the bin distributions of the grids as they have been adapted during the WHIZARD run can be converted into histograms by

make grids

and viewed as the file whizard-grids.ps (note that this only works if the user has set the variable plot_grids_channels in the diagnostics_input block).

Again, it should be kept in mind that the choice of phase space channels does not influence the integration results (all Feynman graphs are included in the cross section calculation and event generation, but for the phase space setup only some of them are used — and some are even used twice). However, the accuracy and the speed of convergence depend on this choice. In many cases it will not be necessary to modify anything, but if you are not satisfied with the results of the default configuration, the following remarks may help.

  • The channels are grouped in groves: Sets of channels that are equivalent in some sense and therefore should get equal weight. All channels within a grove are sampled with the same number of points per iteration which may increase or decrease during the adaptation procedure. To prevent channels from being dropped altogether, one may set min_calls_per_channel to a positive value, e.g. 100.
  • When a Feynman graph is interpreted as a phase-space parameterization channel, the internal lines are considered either as resonant, log-enhanced, or off-shell. A resonant line either corresponds to a resonance which is kinematically allowed to go on-shell. A log-enhanced line is an (approximately) massless propagator, either in the s-channel (e.g., a photon which splits into massless quarks), or in the t channel (e.g., a Coulomb photon which is exchanged between the initial particles).
  • The first criterion for selecting a channel is the number of off-shell lines. The default value is off_shell_lines=1. Since in dominant graphs most internal lines are either resonant or log-enhanced, this already covers many cases.
  • By this criterion, some important graphs are missed, e.g., the Higgs radiation off a top quark in e+e annihiliation, which has two off-shell lines. The rule is that additional channels are selected with extra_off_shell_lines (default 1). For those, the log-enhanced lines count as off-shell. It turns out that this rule is sufficient for essentially all physically relevant cases. However, it may fail in more obscure cases (e.g., if non-resonant regions of phase space are selected by cuts), so increasing these parameters manually, which adds more Feynman diagram classes to the phase-space setup, is a possible option for the user.
  • At high energies, multiperipheral graphs with particle exchange become important even though they may be subleading according to the former rules. This fact is accounted for by setting a threshold_mass_t value equal to 200 GeV, so W and Z bosons in the t-channel are considered light, and W/Z exchange graphs counted log-enhanced. For the s-channel, the threshold_mass value is 50 GeV, so W and Z are considered heavy. Of course, these values may be adjusted if the kinematics is completely different.
  • If no phase space channels are found which satisfy all criteria, WHIZARD will give it another try with off_shell_lines increased.
  • Five additional flags give finer control: single_off_shell_branchings, double_off _shell_branchings, single_off_shell_decays, double_off_shell_decays, massive _fsr. With the exception of the second flag, all are .true. by default, so switching them off will remove phase space channels.
  • In particularly complicated situations you might wish to have full control and write your own configuration in whizard.phs, taking the pre-generated configuration whizard.phx as the starting point. To make this easier, for nonzero extra_off_shell_lines you will find extra channels in the output which are commented out. Instead of finding the binary coding of Feynman graphs yourself, you may just remove the comment characters for a inactive set of phase space channels.

The number of channels in the default configuration tends to be somewhat too large. Many channels turn out to be irrelevant and will be dropped during the self-adaptation procedure. If you can guess which ones those are, you may use the available parameters and options for manual intervention to optimize the channel selection. As an example, the effective Zee coupling is almost zero because of sin2θw∼ 1/4, so graphs with Z emission from electron lines are suppressed. You could remove such graphs, either manually, or by removing the Zee vertex from the whizard.mdl file before starting the program. (You have to be careful not to drop resonant Zee decays, then.) However, in most cases, while this may save some running time for the adaptation, it is not worth the effort.

On the other hand, it may happen that important channels are missed by the default configuration. In this case, the restrictions should be weakened.

4.6.3  What if it does not converge?

Although WHIZARD is able to produce stable results with the default settings for a large variety of scattering processes, it is not possible to guarantee fast convergence of the cross section estimate in all circumstances. There is no simple rule for the achievable accuracy and event generation efficiency for a given process; the examples in Sec. 6 should provide some guideline. In any case the screen output of the adaptation procedure should be inspected, and questionable results should not be used for simulation.

Experience shows that in well-behaved 2→ 4 processes an accuracy (Err/Exp in the output) below 1 and an estimated reweighting efficiency above 10 % is reasonable. Due to the increased dimensionality of phase space, for 2→ 6 processes the numbers are typically worse by a factor 2 to 4.

If convergence is apparently not reached, i.e., there are sizable fluctuations in the result, this may be due to one of the following causes:

  • The integral is actually infinite. One should check the cuts if they are sufficient for a finite result. Note that a singularity with a small coefficient is likely to produce an instability rather than an obvious divergence, if the number of calls per iteration is finite.

    A notorious example are processes where the final state contains µ+eνeνµ or µ+ūdνµ: For zero muon mass, there is no problem, but if the muon mass is given the physical value, this final state can result from µ+µ production and muon decay. The muon width is zero, which is reasonable since the muon decay occurs outside the interaction point, but with cuts on visible particles only, real muon decay makes the integral divergent.

  • The number of iterations is insufficient. For simple 2→ 3 processes, one reaches a stable result after 2−5 iterations, but in complicated cases (many particles with photon emission with weak cuts, etc.) 20−30 iterations may be needed. WHIZARD will make a guess, but this may be not sufficient in some cases.
  • The number of calls per iteration is insufficient. With increasing dimension and complexity of the phase space, this number should be increased as well.
  • The phase space parameterization is inadequate (see above). One may try to include or remove channels by modifying parameters such as the number of off-shell lines allowed in generating parameterizations.
  • The parameterization is adequate in principle, but some subdominant (but still important) channels got dropped during the adaptation procedure. This may be prevented by setting min_calls_per_channel. The channel weight history can be inspected in the process-specific output file.
  • The matrix element is numerically unstable. This is probably not an issue for O’Mega or MadGraph matrix elements with realistic kinematics, but cannot be excluded. CompHEP matrix elements, which result from squared amplitudes, are susceptible for numerical instabilities which occur for extreme kinematics. Actually, if the compiler supports it, one may switch to quadruple precision (configure with --enable-quadruple) in order to ameliorate this problem.
  • One should be careful to avoid the effects of QED gauge invariance violation, which can be large if unstable particles and photons are around (e.g., single W production in e+e collisions). With CompHEP matrix elements, the gwidth switch should be set in such cases, while O’Mega matrix elements should be generated with the ‘f’ option. This enables the overall factor scheme, where the resonant contributions to the matrix element are evaluated for zero width of charged particles, and the Breit-Wigner correction factor is applied to the complete matrix element at the end. Another gauge-invariant prescription (the complex-mass scheme) is foreseen to be implemented in a future version of O’Mega within the context of WHIZARD 2.

4.6.4  Output files

All results are logged in a generic file whizard.out and in process-specific files whizard.XXX.out, where XXX is the process tag. The output files are rewritten after each iteration. Depending on the settings of the diagnostics parameters, they will contain additional information: The contents of the input file are repeated, including also the default values which have been inserted by the program. Furthermore, it shows the detailed history of the individual integration channels and of their relative weights.

After each iteration, the current VAMP grids are written to file as well: whizard.XXX.grc for the current grid and whizard.XXX.grb for the best grid. These may be reused in another run.

4.7  Event generation

4.7.1  Unweighted events

Event generation is activated in the input file, either by specifying a nonzero luminosity

&process_input
...
luminosity = 100
/

or by a nonzero number of events

&simulation_input
n_events = 10000
/

In the former case, the actual number of generated events depends on the calculated cross section. If both numbers are given, the one leading to the larger event sample will be used. It is also possible to specify n_calls, which is the number of matrix element calls. However, this is more useful for weighted event samples (see below).

The simulation makes sense only if WHIZARD could find a stable result for the cross section. This should be judged from inspecting the results displayed on screen and in the log file whizard.out.

With event generation activated, WHIZARD will not stop after the integration step, but start event generation. By default, WHIZARD will generate unweighted events, which result from reweighting the events produced by the VAMP sampling algorithm. In the event generation step, importance sampling is used, such that the events are random and have the proper distribution as required for a realistic simulation.

WHIZARD is able to write events to file in various formats (note that the manual of WHIZARD 2, http://projects.hepforge.org/whizard/manual.pdf), does contain examples for all the different event formats presented here, though the debug format does look different there):

  • The file whizard.evx contains the events in a binary machine-dependent format. By setting read_events=T they can be re-read, bypassing the time-consuming event generation. Writing this file is suppressed if write_events_raw=F.
  • The file whizard.evt (or, alternatively, whizard.pythia or whizard.stdhep, depending on the format) is written if write_events=T. The event file format is controlled by the parameter write_events_format. The supported formats are:
    • 0 Verbose format (ASCII), useful mainly for debugging.
    • 1 HEPEVT format (ASCII). The format corresponds to the HEPEVT standard, such as if the contents of a HEPEVT common block were written in their natural order (Fig. 4.8, 4.9). As an extension, the helicity (−1, 0, or 1) is written for each particle in the same line after JDAHEP. For unfragmented events, the first line will contain further information: the number of beam remnants and outgoing partons listed separately, and the event weight.
    • 2 SHORT format (ASCII). This is similar to HEPEVT, but unnecessary entries are suppressed (Fig. 4.10).
    • 3 LONG format (ASCII). The same as SHORT, but there is an additional line per record which contains the value of the integrand and the integrand ratio (Fig. 4.10). In a normal run (recalculate=F), the ratio is unity (see Sec.5.5).
    • 4 ATHENA format (ASCII). This is the input format used by the ATHENA software (ATLAS). A slight variation of the HEPEVT format above. (cf. Fig. 4.12)
    • 5 Les Houches Accord format LHA (ASCII). This format implements the Les Houches Accord as in the original MadGraph implementation for transferring parton-level events to hadronization packages. To use it, you must set fragment=T, and you may set fragmentation_method=0 since any built-in fragmentation is irrelevant for the output. (cf. Fig. 4.13)
    • 6 Les Houches Accord format LHEF (ASCII). The new version containing XML headers [18], cf. Fig. 4.14. For more details cf. the WHIZARD 2 manual,
      http://projects.hepforge.org/whizard/manual.pdf.
    • 11/12/13 PYLIST formats 1,2,3 (ASCII). This is useful in particular for fully fragmented events, but can be used as well if fragmentation is turned off (see the PYTHIA manual [1] for details).17
    • 20 STDHEP format for writing the HEPEVT common block (binary, machine-independent)
    • 21 STDHEP_UP format for writing the HEPRUP/HEPEUP common blocks (binary, machine-independent)

Note that the formats that require external libraries (PYTHIA, STDHEP) do not work if the corresponding libraries have not been linked to the main program.

Within each event, the particles will appear in the order determined by the process definition /footnotefor O’Mega and MadGraph; CompHEP may reorder them.. Hadronization packages such as JETSET require a particular order of colored particles to properly assign showers, be may necessary for the user to rearrange them before feeding them into the hadronization step. In general, if color-flow information is absent 18, one may set guess_color_flow to tell the program that the particles are to be combined as singlets in the order in which they appear in the event record.

Reweighting Monte-Carlo events requires the knowledge of the highest weight. This can never be known with certainty in advance, so there may be some fraction of events with weight greater than one. By default, WHIZARD will sum up the excess weight of those events and show the induced error in the total cross section, which should be small compared to the statistical error of event generation. If desired, the excess events can be further analyzed by WHIZARD’s own analysis capabilities (see below), to make sure that they have no significant impact on final results. By setting the safety_factor in the simulation_input block to a number greater than one, the amount and weight of excess events (if any) will be reduced. However, the time needed for event generation will be increased by the same factor.

WHIZARD will use the adaptation step to determine the highest weight beforehand. This requires sufficient statistics in the adaptation procedure. If many events have to be generated, the default choice for the number of calls per iteration may be too low, and too many excess events are produced. As a rule of thumb, the number of calls per adaptation iteration should be of the order of the number of unweighted events divided by the estimated efficiency.

Alternatively, a “test” run can be made by specifying a nonzero value of n_events_warmup in the input file. For instance, one may generate 200,000 (weighted) events to determine the highest weight before generating the actual sample to be used for the reweighting procedure:

&simulation_input
n_events = 10000
n_events_warmup = 200000
/

Fortunately, this is not necessary in most cases since the pre-determined highest weight is usually quite accurate.


  integer, parameter :: nmxhep = 4000
  integer :: nevhep, nhep
  integer, dimension(nmxhep) :: isthep, idhep
  integer, dimension(2, nmxhep) :: jmohep, jdahep
  real(kind=double), dimension(5, nmxhep) :: phep
  real(kind=double), dimension(4, nmxhep) :: vhep
  common /HEPEVT/ nevhep, nhep, isthep, idhep, &
       & jmohep, jdahep, phep, vhep
Figure 4.8: Definition of the HEPEVT common block in Fortran 95 format


 write(u, 11) nhep
 do i=1, nhep
    write (u,13) isthep(I), idhep(I), &
         & (jmohep(j,i), j=1,2), (jdahep(j,i), j=1,2)
    write (u,12) (phep(j,i), j=1,5)
    write (u,12) (vhep(j,i), j=1,min(5,ubound(vhep,1)))
 end do

 11 format(10I5)
 12 format(10F17.10)
 13 format(I9,I9,5I5)
Figure 4.9: Output routine which would produce the event file whizard.evt from the contents of a HEPEVT common block. u is an output unit.


write(u,15) evt%n_out + evt%n_beam_remnants, &
     &      evt%n_out, evt%n_beam_remnants, weight
do i=1, evt%n_out
   p = array(particle_momentum(evt%prt(i)))
   write (u,13) particle_code(evt%prt(i)), evt%hel_out(i)
   write (u,12) p(1:3), p(0), particle_mass(evt%prt(i))
end do
12 format(10F17.10)
13 format(I9,I9,5I5)
15 format(3I5,1X,F17.10)
Figure 4.10: Output routine for the SHORT event file format. Note that the number of beam remnants and the event weight are not part of the HEPEVT standard. (We here and in the sequel do not include the code for the beam remnants and/or initials.)


write(u,15) evt%n_out + evt%n_beam_remnants, &
     &      evt%n_out, evt%n_beam_remnants, weight
do i=1, evt%n_out
   p = array(particle_momentum(evt%prt(i)))
   write(u,13) particle_code(evt%prt(i)), evt%hel_out(i)
   write(u,12) p(1:3), p(0), particle_mass(evt%prt(i))
end do
write(u, 14) evt%mc_function_value, evt%mc_function_ratio
12 format(10F17.10)
13 format(I9,I9,5I5)
14 format(ES17.10,F17.10)
15 format(3I5,1X,F17.10)
Figure 4.11: Output routine for the LONG event file format. Note that the number of beam remnants and the event weight are not part of the HEPEVT standard.


write(u,11) evt%count, event_n_entries (evt)
k = 0
do i=1, evt%n_out
   k = k + 1
   p = array(particle_momentum(evt%prt(i)))
   write (u,13) k, 1, particle_code(evt%prt(i)), &
        & evt%n_in+1, 2*evt%n_in, (0, j=1,2)
   write (u,12) p(1:3), p(0), particle_mass(evt%prt(i))
   write (u,12) (0._default, j=1,4)
end do
11 format(10I5)
12 format(10F17.10)
13 format(I9,I9,5I5)
Figure 4.12: Output routine for the ATHENA event file format. (If the initial particles and/or the beam remnants are to be kept, these are written before the outgoing particles.)


write (u,16) nup, idprup, xwgtup, scalup, aqedup, aqcdup
write (u,17) idup(:nup)
write (u,17) mothup(1,:nup)
write (u,17) mothup(2,:nup)
write (u,17) icolup(1,:nup)
write (u,17) icolup(2,:nup)
write (u,17) istup(:nup)
write (u,17) ispinup(:nup)
do i=1, nup
   write (u,18) i, pup((/4,1,2,3/),i)
end do
16 format(2(1x,I5),1x,F17.10,3(1x,F13.6))
17 format(500(1x,I5))
18 format(1x,I5,4(1x,F17.10))
Figure 4.13: Output routine for the (old) Les Houches Accord (LHA) event file format.


  subroutine lhef_write_init (unit, process_id)
    integer, intent(in) :: unit
    integer :: i
    character(*), dimension(:), intent(in) :: process_id
    write (unit, "(A)")  '<LesHouchesEvents version="1.0">'
    write (unit, "(A)")  '<header>'
    write (unit, "(2x,A)")  '<generator>'
    write (unit, "(4x,A)")  &
         '<name version="<<Version>>" date="<<Date>>">' // &
         'WHIZARD</name>'
    do i = 1, size (process_id)
       write (unit, "(4x,A)") &
            '<process_id>' // trim (adjustl (process_id(i))) &
            // '</process_id>'
    end do
    write (unit, "(2x,A)")  '</generator>'
    write (unit, "(A)")  '</header>'
  end subroutine lhef_write_init

  subroutine lhef_write_heprup (unit, process_id)
    character(*), dimension(:), intent(in) :: process_id
    integer, intent(in) :: unit
    integer :: i
    write (unit, "(A)") '<init>'
    write (unit, "(2(2x,I10),2(2x,G22.15E3),8x,6(2x,I10))") &
         idbmup(1:2), ebmup(1:2), &
         pdfgup(1:2), pdfsup(1:2), &
         idwtup, nprup
    do i = 1, size (process_id)
       write (unit, "(3(2x,G22.15E3),8x,1(2x,I10))") &
            xsecup(i), xerrup(i), xmaxup(i), lprup(i)
    end do
    write (unit, "(A)") '</init>'
  end subroutine lhef_write_heprup

  subroutine lhef_write_hepeup (unit, pro_id)
    integer, intent(in) :: unit, pro_id
    integer :: i
    idprup = pro_id
    write (unit, "(A)") '<event>'
    write (unit, "(2(2x,I10),2(2x,G22.15E3),8x,2(2x,G22.15E3))") &
         nup, idprup, xwgtup, scalup, aqedup, aqcdup
    do i = 1, nup
       write (unit, "(6(2x,I10),8x,3(2x,G22.15E3),8x,2(2x,G22.15E3),32x,2(2x,G22.15E3))") &
            idup(i), istup(i), &
            mothup(1:2,i), icolup(1:2,i), &
            pup(1:5,i), vtimup(i), spinup(i)
    end do
    write (unit, "(A)") '</event>'
  end subroutine lhef_write_hepeup

  subroutine lhef_write_final (unit)
    integer, intent(in) :: unit
    write (unit, "(A)")  '</LesHouchesEvents>'
  end subroutine lhef_write_final
Figure 4.14: Output routine for the (new) Les Houches Accord (LHEF) event file format.

4.7.2  Weighted events

The default setup will generate unweighted events, which is the proper way of simulating actual physics runs. The cross section error of the simulation run is equal to the square root of the number of events, as it should be.

This requires an unweighting step where, depending on the ratio of the current event’s weight and the maximum weight, a large fraction of the generated event sample is thrown away. Clearly, information is lost in this step. To keep all available information, one can use weighted events instead:

&simulation_input
unweighted = F
...
/

The number of weighted events is calculated such that (approximately) the same unweighted event sample would result from unweighting them afterwards, as it would from directly generating unweighting events. If luminosity or n_events do not require a higher number, the total number of weighted events generated is approximately equal to n_calls.

If there is more than one process, the luminosity is distributed among them proportional to the individual cross sections. In contrast to unweighted generation, the events are ordered, such that all events belonging to the first process are generated first, then the next one, and so on. The exact number of weighted events for a given process depends on both the cross section and the efficiency. For a process with low efficiency, a larger number of weighted events is needed for the same cross section than for a high-efficiency process. (In unweighted event samples, the individual event numbers depend only on the cross sections.)

These events can also be written to file and analyzed by an external program. For each event, the event weight is a number between zero and one (if no excess events are present), and it is printed as the last number in the first line of the output. Hence, in this definition the weight is the probability with which the event would have been kept in the sample if an unweighting step had taken place. As a result, differential (binned) distributions should be identical for weighted and unweighted event samples, apart from fluctuations.

Unfortunately, the HEPEVT standard does not allow for event weights. Therefore, the PYLIST and STDHEP output formats which read out the HEPEVT common block cannot be used for weighted events. By contrast, event weights can be transferred using the HEPEUP common block, so the STDHEP_UP format does work for weighted events. Keep in mind that weighted events are ordered, so files containing those must not be truncated or extended.

4.7.3  Built-in analysis

WHIZARD has capabilities to analyze the events without reference to external packages. To this end, a configuration file whizard.cut5 should be set up which has a similar format as whizard.cut1. For instance, in Fig. 4.15 a configuration file is shown which tells WHIZARD to analyze the generated event sample in three different ways: to calculate the cross section within a window 80<M(e+e)<100 GeV, then within 180<M(e+e)<200 GeV with an additional cut on p(γ), and finally to generate a histogram of p(e) with a cut on the photon energy.


! e-  e+  ->  e-  e+  gamma
! 16   8       1   2    4
process eeg
  cut M of 3 within  80 100
and
  cut M  of 3 within 180 200
  cut PT of 4 within 100 99999
and
  cut E of 4 within 0 100
  histogram PT of 1 within 0 500
Figure 4.15: Analysis configuration file.

The histogram specification may take either one of the form

histogram PT of 1 within 0 500
histogram PT of 1 within 0 500 nbin 50
histogram PT of 1 within 0 500 step 10

Instead of the number of bins (20 by default), the bin width can be specified by the keyword step.

The cross section results are displayed on screen. In addition, they are written together with the histogram data to the file whizard-plots.dat. If show_excess=T has been set in the input file, the summary and the histograms will contain also information about any events with weight greater than one, such that their properties can be examined.

Clearly, for a quick analysis it is not necessary to simulate proper unweighted events. Therefore, one may set unweighted=F in order to analyze a sample of weighted events, which saves a lot of execution time.

The contents of whizard-plots.dat can be visualized directly by means of the gamelan graphical analysis package which is contained in the distribution. This requires the existence of the MetaPost program, which is part of standard LATEX distributions. Plots can then be made by typing

make plots

This will result in a Postscript file whizard-plots.ps which can be viewed and printed. Note that this works only if you have generated events, otherwise WHIZARD might fail with the message

No rule to make target plots.

The event sample can further be analyzed in terms of color flow, helicity and flavor content. Example

! e-  e+  ->  q  qbar  b  bbar
! 16   8      1   2    4   8
process qqbb
  cut M of 3 within 85 95
  color 0 0 4 0 6 0
  helicity  1 -1 -1 1
  helicity -1  1 -1 1
  flavor    2 -2 5 -5
  histogram M of 12 within 0 500

This would plot a histogram of the invariant mass of the bb pair with the additional requirement that the other quark pair must be uū, that the b (anti)quarks must be left- resp. right-handed, and that the uū and bb pairs are each in a color-singlet state.

Color flow is defined in a natural way: The particles are numbered consecutively beginning from the initial ones (do not confuse this with the binary codes used for specifying cuts and histograms). For an incoming quark or gluon, the index of the particle its color flows to is given, while for an outgoing quark or gluon it is the index of the particle its color comes from. (This differs from WHIZARD’s internal convention where incoming and outgoing particles are equally treated, but it is easier to use.) For colorless particles and for antiquarks, the entry is zero. So, in the above example, the q color comes from the q, and the b color from the b.

Helicities and flavors are indicated as usual: Fermions and massless vector bosons have helicity ± 1, massive vector bosons can come with helicity 0 in addition, and for scalars the helicity is always zero. Flavor is denoted by standard PDG codes.

The order should be as in the above example: first cut specifications (the event must pass all cuts), then color, helicity and flavor specifications (the event must match any of the given helicity combinations and any of the given flavor combinations), finally histogram specifications. Any of these may be absent; an empty entry will just display the total cross section again. More analysis configurations can be added for the same process, they must be separated by the keyword and.

4.8  Fragmentation

The Les Houches standard [9] for transferring events together with weight and color information allows for a consistent interface of partonic event generators and programs that account for showering, fragmentation, and hadronization. In order to use such an external program, WHIZARD has to write events to file in a format that supports this information. Currently, there is the LHA format (code: 5) that writes, for each event, the contents of the HEPEUP block in ASCII format, and the STDHEP_UP format (code: 21) that writes this information into a portable binary STDHEP file. The first record of this file will contain the HEPRUP common block; then, for each event, the HEPEUP common block is written.

Alternatively, WHIZARD provides an internal interface to routines for fragmentation, decays and hadronization. For each event, WHIZARD calls the corresponding external program that fills the HEPEVT common block and writes this to file. To activate this, the parameter fragment in the simulation_input block must be set to true. Of course, the corresponding library has to be properly linked. Only unweighted events can be fragmented. (This is because the HEPEVT standard does not allow for the event weight to be transferred.)

Once fragment is set, various options are available, depending on the value of the parameter fragmentation_method:

  1. No fragmentation within WHIZARD (for completeness).
  2. JETSET fragmentation. The final state is fragmented and hadronized by a call to the PYEXEC fragmentation routine inside PYTHIA, i.e., by the former JETSET program. If available, color flow information is used. Otherwise, i.e. for CompHEP matrix elements and O’Mega without the ’c’ option, the color connection properties of the final state are guessed by JETSET from the particle ordering, which must therefore carefully be chosen by the user when defining the process. Note that even for processes like e+ejj where O’Mega without the option ’c’ gives the correct color factors, it only gives the correct color flows for (JETSET) fragmentation when using the ’c’ option.
  3. PYTHIA fragmentation. The final state is fragmented and hadronized by making the process generated by WHIZARD an external process within PYTHIA. If possible and depending on the PYTHIA parameter settings, beam remnants are constructed, multiple interactions and pile-up events are added, and more.
  4. This parameter value is intended to transfer control to fragmentation routines written by the user and included in the user module. Of course, such code can make use of the HEPRUP and HEPEUP common block contents. For additional flexibility, there is also a parameter user_fragmentation_mode (integer) that allows for selecting different versions of user fragmentation.

Full color information is available when using the ’c’ option for O’Mega matrix elements (or for MadGraph matrix elements). However, in the other cases, for PYTHIA fragmentation WHIZARD will provide missing color information if can be uniquely determined. Otherwise, it tries to guess color connections from the particle ordering in the event record (just as JETSET does). If color connection information cannot be determined, colored states cannot be fragmented in this way. If fragmentation is nevertheless desired (even though the results may not be accurate), one has to revert to JETSET fragmentation.

No effort is done to streamline WHIZARD and PYTHIA parameters, so changing a value in the WHIZARD input file will have no effect on the fragmentation step. However, PYTHIA parameters can be explicitly set by an interface to PYTHIA’s pygive routine. For instance, the line

  pythia_parameters = "MSTP(41)=0; MSTP(81)=0"

in the simulation_input block disables resonance decays and multiple interactions within PYTHIA. Similarly, the Higgs mass may be set by a call like

  pythia_parameters = "PMAS(25,1)=115"

The detailed list of parameters available is described in the PYTHIA manual [1].

The internal analysis routines described in Sec. 4.7.3 always refer to the event before fragmentation. To analyze the fragmented event, the user must use an external package. Note that WHIZARD 2 meanwhile allows for analyzing fragmented events with the internal analysis tools, however.

4.8.1  Parton showers and matching

In WHIZARD there is a version of a kT-ordered parton shower included written by Sebastian Schmidt. This is, however, a highly experimental version mainly intended for testing purposes. It is quite similar to the PYTHIA parton shower [1], but does contain only very basic tunes to LEP Run 2 jet observables. Furthermore, it is only a shower for final state radiation (FSR); there is no implementation of an initial state shower in WHIZARD 1. The internal shower has been improved a lot and augmented by a different (so-called analytical) parton shower within WHIZARD 2. For more details about the shower algorithm confer the WHIZARD 2 manual http://projects.hepforge.org/whizard/manual.pdf. The internal parton shower can be switched on for event generation in the simulation_input block. Note that the only event format which supports showered events from the WHIZARD-internal parton shower is the (new) Les Houches event format (LHEF, event format 6) as this is the only format that allows for a variable number of particles.

There is, of course, also the possibility to use showering (and hadronization) using the PYTHIA interface. In the default settings, PYTHIA events are constructed with parton showers developing along colored and charged partons. In WHIZARD, fixed-order matrix elements are evaluated. The matching of both is nontrivial. For this reason, PYTHIA parton showers are disabled in the current WHIZARD version. In general, the availability of exact higher-order matrix elements (albeit tree-level) within WHIZARD makes matching without double-counting a delicate task. This problem of potential double counting has been (partially) solved in the framework of so called matching algorithms by either vetoing shower events entering hard matrix elements regimes or reweighting with the corresponding Sudakov factors according to the shower history. None of these algorithms have been implemented within WHIZARD 1, but they are (at least partially) available inside WHIZARD 2.

There is also an experimental Herwig interface, but we do not discuss the details here.

4.8.2  Mixing WHIZARD and PYTHIA events

Since WHIZARD provides an external process to PYTHIA when PYTHIA fragmentation is enabled, it is perfectly possible to mix PYTHIA and WHIZARD events. To this end, one should specify PYTHIA processes in the simulation_input block. For instance,

  pythia_processes = "22, 24"

will simulate e+eZH and ZZ in addition to the WHIZARD-generated processes. There are some caveats, however:

  • PYTHIA processes are generated only if PYTHIA fragmentation is enabled. Consequently, the PYTHIA events are mixed with the fully fragmented WHIZARD events and optionally written to file. There is no way to use the built-in analysis for those events.
  • Since parton showers (ISR and FSR) are disabled globally, this applies to the PYTHIA processes as well. However, one could switch the showering on by modifying the inputs.
  • When the number of events to generate is derived from the luminosity input parameter, it is assumed that no PYTHIA events are present. To obtain the actual luminosity corresponding to an event sample in the presence of PYTHIA events, this number has to be multiplied by the fraction of WHIZARD events in this sample. The parameter n_events always refers to the total number of events. When show_pythia_statistics is true, a summary of the generated channels and cross sections is printed.
  • One should be careful not to double-count processes, except if this is actually desired (e.g., for comparing WHIZARD and PYTHIA event samples for related processes).

4.8.3  PYTHIA event display

When fragmentation is enabled, the conventions for writing events to file are changed. The verbose format is then not available since some information (used internally by WHIZARD) is lost after fragmentation. The format codes 1 (HEPEVT), 2 (SHORT), 3 (LONG), and 20/21 (STDHEP) are available. In addition, there are three new format codes (11, 12, 13) which allow for displaying events using PYTHIA’s PYLIST call. The respective argument given to PYLIST is 1, 2, or 3, corresponding to three different verbosity levels. (For automatic processing, the more condensed HEPEVT or STDHEP formats are better suited.) Note again that these output formats might not work properly if you use a CERN library compiled with a different compiler than the WHIZARD core program, because in this case the output units of the two compiler libraries are independent of each other. The PYTHIA events then just appear in a file fort.<unit>, where <unit> is the current output unit.

4.9  Controlling WHIZARD runs

For a new process, it is not easy to predict how long it will take to reach a sufficiently precise result or to generate a sufficiently large event sample. The achievable efficiency and accuracy and the overall running time strongly depend on the complexity of the process, the singularity structure of the matrix element and how well WHIZARD is able to map it, the presence or absence of structure functions, and, clearly, the evaluation time of the matrix element for a single event.

Ideally, computer time and space is abundant, one calcaulates cross sections to the highest possible accuracy, and infinitely large event samples are generated. In practice, one does not want to wait too long for useable results, and one is limited by CPU time and space restrictions. Priorities have to be defined. To meet these needs, various switches and levers are available that control WHIZARD runs. All of them are set in the main input file and are listed in the tables of Sec. 4.4.

4.9.1  Grid adaptation

Parameters that control the grid adaptation and integration passes of WHIZARD are collected in the integration_input block in the input file, listed in Sec. 4.4.3

Obviously, the array calls is the principal place for improving the accuracy and efficiency of the adaptation and integration runs. In no case more iterations or calls than specified in calls will be executed, where WHIZARD uses default values for entries left empty by the user.

Increasing the number of adaptation iterations normally improves (decreases) the constant α in the error estimate for the total cross section,

Δσ
σ
 = 
α
Nevents
,     (5)

until a limit is reached which depends on the matrix element singularity structure and on the choice of phase space integration channels. The estimate for α is printed as the accuracy in the screen output table. Increasing the number of points per iteration tends to decrease the number of iterations needed to reach this limit. In particular, this improves the precision in estimating α and reduces inevitable fluctuations around the ideal grid configuration. Simultaneously the reweighting efficiency є increases, which is the ratio of the average to the maximum event weight. Note that fluctuations in the є estimate are significantly stronger than for α.

If the program has been killed for any reason, it can be restarted with read_grids set, so adaptation is resumed after the last complete iteration. This mechanism can also be used to increase the number of iterations afterwards for a job that has already been completed.

However, one wants to avoid a unnecessary large number of iterations and calls. Adaptation iterations are skipped if any one of the following parameters is set:

accuracy_goal:
If this number is greater than zero, adaptation will stop when the accuracy estimate of the last iteration is smaller.
efficiency_goal:
If this number is greater than zero, adaptation will stop when the efficiency estimate of the last iteration is higher. Note that the efficiency is measured in percent.
time_limit_adaptation:
If this number is greater than zero, adaptation will stop when the real time after the last iteration is later than the specified time limit, measured in minutes from the program start.

In any case, the program will not stop but the third pass (integration) will be executed after the adaptation pass, where the number of integration steps and calls per iteration are taken from the corresponding calls entry.

4.9.2  Simulation

The reweighting efficiency in the final event generation pass is determined by the grid quality, measured by the efficiency є. Since the precision in estimating є is not very high, it may be appropriate to insert a safety_factor f>1 to avoid excess events with weight >1. The time needed for generating a fixed number of unweighted events is then proportional to 1/(fє). As a rule, investing more time in grid adaptation results in a better є value, and thus less time for event generation, up to a certain limit.

The number of unweighted events is specified by either setting the parameter n_events in the simulation_input block, or by a luminosity value in the process_input block. The number that leads to a larger event sample will take priority.

If the program is killed during event generation, the event sample is incomplete. The switch catch_signals can be used to catch terminate signals if the OS and compiler supports it, so it is guaranteed the the event files are properly closed before the program terminates. Then, WHIZARD can be restarted with read_events set (which implies read_grids). The events already generated will be reread from the raw event file whizard.evx until the end is reached, and further events will be appended to this file. Event output files in other formats are re-generated on the fly. Thisy method can also be used to increase an event sample afterwards, spending CPU time essentially on newly generated events only.

4.9.3  Hard limits

There are situations where one has to control the amount of memory, CPU time, and disk space which the program uses. In particular, batch installations put limits on these resources and kill processes which exceed them.

Some of these issues can be handled by WHIZARD:

time_limit:
This parameter (in minutes) sets an internal limit on the wallclock time the run may use. When it is exceeded, the program terminates. Depending whether the termination occured during integration or simulation, the program can be restarted such that it will resume after the last grid or event written, respectively.
events_per_file:
This parameter takes care of event files becoming too large: once it is exceeded, a new file is begun. An external analysis would then process the files consecutively. Note that this does not apply to the internal event file whizard.evx which cannot be split; if this file becomes too large, one has to switch it off by write_events_raw=F.
bytes_per_file:
The analogous parameter which counts the number of bytes allowed for an event file. This relies on certain assumptions on the number of bytes an event takes in a particular format.

1
A better method is to give the user-created file a different name and use a symbolic link to whizard.prc instead. The same applies to whizard.in.
2
The MadGraph version implemented in WHIZARD does not treat multi-scalar vertices correctly. Don’t use it if you are interested in double (or triple) Higgs production. CompHEP uses Feynman gauge (i.e., includes ghosts as external states) for a more efficient evaluation of electroweak processes, while the other generators use unitarity gauge. This should not be visible to the user. Note also that CompHEP provides polarization-averaged matrix elements only. There is a CompHEP model using unitarity gauge, SM_ug, which turns out to be less efficient in many cases).
3
The CompHEP version currently has only the anomalous quartic couplings. The MadGraph version is just the Standard Model.
4
The MSSM is supported only by the O’Mega matrix element generator, while the included MadGraph and CompHEP versions will fall back to the Standard Model.
5
Furthermore, by writing the Feynman rules into a Lagrangian file, one may easily modify or set up a new model for CompHEP, which then can be used by WHIZARD. Modifying the model for O’Mega requires changing program sources, however. The same holds for MadGraph, if you are interested in Lorentz or Dirac structures that are not present in the SM.
6
It does not work for CompHEP matrix elements, however.
7
The current implementation is not well suited for single-particle processes if the initial-state energy distributions are strongly peaked. A similar caveat holds for processes with an s-channel resonance, for instance, e+e→µ+µ with a continuous beamstrahlung/ISR spectrum for the incoming electrons, integrated over the whole energy range. Nevertheless, with large statistics reasonable results can be obtained. Some improvement has been achieved by using a specific mapping for s-channel like topologies.
8
Not available for CompHEP or MadGraph matrix elements)
9
Naively summing over identical particles in the final state may result in double-counting, and becomes very inefficient if many combinations are allowed. The ‘c’ option takes care of both issues (see below).
10
If the list of flavor sums or the number of particles in the process becomes very long, it may happen that the compiler complains for the number of continuation lines becoming too large. This problem is easily fixed by increasing the number of continuation lines accepted by the compiler. You should consult your Fortran compiler documentation for the appropriate flag which may then be set in the config.site file. Don’t forget to clean up and configure before restarting the compilation.
11
And also for CompHEP.
12
CompHEP, unfortunately, reorders particles according its own conventions.
13
Those that are part of the WHIZARD system. Precompiled libraries such as the LHAPDF have to be present on the target system and can be enabled there if they are needed.
14
For MadGraph in particular, if there is an ambiguity you may have to set the QCD order as an option.
15
Note that if WHIZARD does not find the corresponding file it takes its default file.
16
Note that we do not keep this updated. The parameter file of version 2, however, contain the most recent PDG values
17
This works only if PYTHIA has been linked. Also it does not work (properly) if WHIZARD and the CERN library containing PYTHIA have been compiled by different compilers, because PYTHIA or JETSET cannot connect the proper output unit in that case. So then a file like fort.<unit> will appear containing the PYTHIA/JETSET events.
18
For CompHEP matrix elements, and for O’Mega without the ’c’ option.

Previous Up Next