1 // ProcessLevel.h is a part of the PYTHIA event generator.
2 // Copyright (C) 2012 Torbjorn Sjostrand.
3 // PYTHIA is licenced under the GNU GPL version 2, see COPYING for details.
4 // Please respect the MCnet Guidelines, see GUIDELINES for details.
6 // This file contains the main class for process-level event generation.
7 // ProcessLevel: administrates the selection of "hard" process.
9 #ifndef Pythia8_ProcessLevel_H
10 #define Pythia8_ProcessLevel_H
13 #include "BeamParticle.h"
16 #include "ParticleData.h"
17 #include "PartonDistributions.h"
18 #include "ProcessContainer.h"
19 #include "PythiaStdlib.h"
20 #include "ResonanceDecays.h"
22 #include "SigmaTotal.h"
23 #include "SusyCouplings.h"
24 #include "SusyLesHouches.h"
25 #include "StandardModel.h"
26 #include "UserHooks.h"
30 //==========================================================================
32 // The ProcessLevel class contains the top-level routines to generate
33 // the characteristic "hard" process of an event.
40 ProcessLevel() : iLHACont(-1) {}
42 // Destructor to delete processes in containers.
46 bool init( Info* infoPtrIn, Settings& settings,
47 ParticleData* particleDataPtrIn, Rndm* rndmPtrIn,
48 BeamParticle* beamAPtrIn, BeamParticle* beamBPtrIn,
49 Couplings* couplingsPtrIn, SigmaTotal* sigmaTotPtrIn, bool doLHAin,
50 SusyLesHouches* slhaPtrIn, UserHooks* userHooksPtrIn,
51 vector<SigmaProcess*>& sigmaPtrs, ostream& os = cout);
53 // Store or replace Les Houches pointer.
54 void setLHAPtr( LHAup* lhaUpPtrIn) {lhaUpPtr = lhaUpPtrIn;
55 if (iLHACont >= 0) containerPtrs[iLHACont]->setLHAPtr(lhaUpPtr);}
57 // Generate the next "hard" process.
58 bool next( Event& process);
60 // Special case: LHA input of resonance decay only.
61 bool nextLHAdec( Event& process);
63 // Accumulate and update statistics (after possible user veto).
66 // Print statistics on cross sections and number of events.
67 void statistics(bool reset = false, ostream& os = cout);
70 void resetStatistics();
72 // Add any junctions to the process event record list.
73 void findJunctions( Event& junEvent);
75 // Initialize and call resonance decays separately.
76 void initDecays( Info* infoPtrIn, ParticleData* particleDataPtrIn,
77 Rndm* rndmPtrIn, LHAup* lhaUpPtrIn) { infoPtr = infoPtrIn;
78 resonanceDecays.init( infoPtrIn, particleDataPtrIn, rndmPtrIn);
79 containerLHAdec.setLHAPtr(lhaUpPtrIn, particleDataPtrIn); }
80 bool nextDecays( Event& process) { return resonanceDecays.next( process);}
84 // Constants: could only be changed in the code itself.
85 static const int MAXLOOP;
87 // Generic info for process generation.
88 bool doSecondHard, doSameCuts, allHardSame, noneHardSame,
89 someHardSame, cutsAgree, cutsOverlap, doResDecays;
90 int nImpact, startColTag;
91 double mHatMin1, mHatMax1, pTHatMin1, pTHatMax1, mHatMin2, mHatMax2,
92 pTHatMin2, pTHatMax2, sigmaND, sumImpactFac, sum2ImpactFac;
94 // Vector of containers of internally-generated processes.
95 vector<ProcessContainer*> containerPtrs;
96 int iContainer, iLHACont;
99 // Ditto for optional choice of a second hard process.
100 vector<ProcessContainer*> container2Ptrs;
104 // Single half-dummy container for LHA input of resonance decay only.
105 ProcessContainer containerLHAdec;
107 // Pointer to various information on the generation.
110 // Pointer to the particle data table.
111 ParticleData* particleDataPtr;
113 // Pointer to the random number generator.
116 // Pointers to the two incoming beams.
117 BeamParticle* beamAPtr;
118 BeamParticle* beamBPtr;
120 // Pointer to Standard Model couplings, including alphaS and alphaEM.
121 Couplings* couplingsPtr;
123 // Pointer to SigmaTotal object needed to handle soft QCD processes.
124 SigmaTotal* sigmaTotPtr;
126 // Pointer to SusyLesHouches object for interface to SUSY spectra.
127 SusyLesHouches* slhaPtr;
129 // Pointer to userHooks object for user interaction with program.
130 UserHooks* userHooksPtr;
132 // Pointer to LHAup for generating external events.
135 // Initialization of some SLHA blocks.
138 // ResonanceDecay object does sequential resonance decays.
139 ResonanceDecays resonanceDecays;
141 // Generate the next event with one interaction.
142 bool nextOne( Event& process);
144 // Generate the next event with two hard interactions.
145 bool nextTwo( Event& process);
147 // Append the second to the first process list.
148 void combineProcessRecords( Event& process, Event& process2);
150 // Check that colours match up.
151 bool checkColours( Event& process);
153 // Print statistics when two hard processes allowed.
154 void statistics2(bool reset, ostream& os = cout);
156 // Statistics for Les Houches event classification.
157 vector<int> codeLHA, nEvtLHA;
161 //==========================================================================
163 } // end namespace Pythia8
165 #endif // Pythia8_ProcessLevel_H