| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // 26 //
27 /// \file eventgenerator/pythia/decayer6/src/P 27 /// \file eventgenerator/pythia/decayer6/src/Pythia6.cc
28 /// \brief Implementation of the Pythia6 class 28 /// \brief Implementation of the Pythia6 class
29 29
30 // ------------------------------------------- 30 // ----------------------------------------------------------------------------
31 // According to TPythia6 class from Root: 31 // According to TPythia6 class from Root:
32 // (The TPythia6 class is an interface class t 32 // (The TPythia6 class is an interface class to F77 routines in Pythia6 //
33 // CERNLIB event generators, written by T.Sjos << 33 // CERNLIB event generators, written by T.Sjostrand.)
34 // http://root.cern.ch/ 34 // http://root.cern.ch/
35 // see http://root.cern.ch/root/License.html 35 // see http://root.cern.ch/root/License.html
36 // 36 //
37 // The complete Pythia6 documentation can be f 37 // The complete Pythia6 documentation can be found at:
38 // http://home.thep.lu.se/~torbjorn/pythiaaux/ 38 // http://home.thep.lu.se/~torbjorn/pythiaaux/recent.html
39 // ------------------------------------------- 39 // ----------------------------------------------------------------------------
40 40
41 // ******************************************* 41 // ******************************************************************************
42 // ******************************************* 42 // ******************************************************************************
43 // ** 43 // ** **
44 // ** 44 // ** **
45 // ** *......* W 45 // ** *......* Welcome to the Lund Monte Carlo! **
46 // ** *:::!!:::::::::::* 46 // ** *:::!!:::::::::::* **
47 // ** *::::::!!::::::::::::::* P 47 // ** *::::::!!::::::::::::::* PPP Y Y TTTTT H H III A **
48 // ** *::::::::!!::::::::::::::::* P 48 // ** *::::::::!!::::::::::::::::* P P Y Y T H H I A A **
49 // ** *:::::::::!!:::::::::::::::::* P 49 // ** *:::::::::!!:::::::::::::::::* PPP Y T HHHHH I AAAAA **
50 // ** *:::::::::!!:::::::::::::::::* P 50 // ** *:::::::::!!:::::::::::::::::* P Y T H H I A A **
51 // ** *::::::::!!::::::::::::::::*! P 51 // ** *::::::::!!::::::::::::::::*! P Y T H H III A A **
52 // ** *::::::!!::::::::::::::* !! 52 // ** *::::::!!::::::::::::::* !! **
53 // ** !! *:::!!:::::::::::* !! T 53 // ** !! *:::!!:::::::::::* !! This is PYTHIA version 6.418 **
54 // ** !! !* -><- * !! L 54 // ** !! !* -><- * !! Last date of change: 9 Jun 2008 **
55 // ** !! !! !! 55 // ** !! !! !! **
56 // ** !! !! !! N 56 // ** !! !! !! Now is 0 Jan 2000 at 0:00:00 **
57 // ** !! !! 57 // ** !! !! **
58 // ** !! lh !! D 58 // ** !! lh !! Disclaimer: this program comes **
59 // ** !! !! w 59 // ** !! !! without any guarantees. Beware **
60 // ** !! hh !! o 60 // ** !! hh !! of errors and use common sense **
61 // ** !! ll !! w 61 // ** !! ll !! when interpreting results. **
62 // ** !! !! 62 // ** !! !! **
63 // ** !! C 63 // ** !! Copyright T. Sjostrand (2008) **
64 // ** 64 // ** **
65 // ** An archive of program versions and docum 65 // ** An archive of program versions and documentation is found on the web: **
66 // ** http://www.thep.lu.se/~torbjorn/Pythia.h 66 // ** http://www.thep.lu.se/~torbjorn/Pythia.html **
67 // ** 67 // ** **
68 // ** When you cite this program, the official 68 // ** When you cite this program, the official reference is to the 6.4 manual: **
69 // ** T. Sjostrand, S. Mrenna and P. Skands, J 69 // ** T. Sjostrand, S. Mrenna and P. Skands, JHEP05 (2006) 026 **
70 // ** (LU TP 06-13, FERMILAB-PUB-06-052-CD-T) 70 // ** (LU TP 06-13, FERMILAB-PUB-06-052-CD-T) [hep-ph/0603175]. **
71 // ** 71 // ** **
72 // ** Also remember that the program, to a lar 72 // ** Also remember that the program, to a large extent, represents original **
73 // ** physics research. Other publications of 73 // ** physics research. Other publications of special relevance to your **
74 // ** studies may therefore deserve separate m 74 // ** studies may therefore deserve separate mention. **
75 // ** 75 // ** **
76 // ** Main author: Torbjorn Sjostrand; Departm 76 // ** Main author: Torbjorn Sjostrand; Department of Theoretical Physics, **
77 // ** Lund University, Solvegatan 14A, S-223 77 // ** Lund University, Solvegatan 14A, S-223 62 Lund, Sweden; **
78 // ** phone: + 46 - 46 - 222 48 16; e-mail: 78 // ** phone: + 46 - 46 - 222 48 16; e-mail: torbjorn@thep.lu.se **
79 // ** Author: Stephen Mrenna; Computing Divisi 79 // ** Author: Stephen Mrenna; Computing Division, GDS Group, **
80 // ** Fermi National Accelerator Laboratory, 80 // ** Fermi National Accelerator Laboratory, MS 234, Batavia, IL 60510, USA; **
81 // ** phone: + 1 - 630 - 840 - 2556; e-mail: 81 // ** phone: + 1 - 630 - 840 - 2556; e-mail: mrenna@fnal.gov **
82 // ** Author: Peter Skands; Theoretical Physic 82 // ** Author: Peter Skands; Theoretical Physics Department, **
83 // ** Fermi National Accelerator Laboratory, 83 // ** Fermi National Accelerator Laboratory, MS 106, Batavia, IL 60510, USA; **
84 // ** and CERN/PH, CH-1211 Geneva, Switzerla 84 // ** and CERN/PH, CH-1211 Geneva, Switzerland; **
85 // ** phone: + 41 - 22 - 767 24 59; e-mail: 85 // ** phone: + 41 - 22 - 767 24 59; e-mail: skands@fnal.gov **
86 // ** 86 // ** **
87 // ** 87 // ** **
88 // ******************************************* 88 // ******************************************************************************
89 89
90 #include "Pythia6.hh" 90 #include "Pythia6.hh"
91 91
>> 92 #include <iostream>
92 #include <cstdlib> 93 #include <cstdlib>
93 #include <cstring> 94 #include <cstring>
94 #include <iostream> <<
95 95
96 #ifndef WIN32 96 #ifndef WIN32
97 # define pycomp pycomp_ << 97 # define pycomp pycomp_
98 # define py1ent py1ent_ << 98 # define py1ent py1ent_
99 # define type_of_call << 99 # define type_of_call
100 #else 100 #else
101 # define pycomp PYCOMP << 101 # define pycomp PYCOMP
102 # define py1ent PY1ENT << 102 # define py1ent PY1ENT
103 # define type_of_call _stdcall << 103 # define type_of_call _stdcall
104 #endif 104 #endif
105 105
106 // pythia6 functions 106 // pythia6 functions
107 extern "C" << 107 extern "C" {
108 { << 108 int type_of_call pycomp(int *kf);
109 int type_of_call pycomp(int* kf); <<
110 void type_of_call py1ent(int&, int&, double& 109 void type_of_call py1ent(int&, int&, double&, double&, double&);
111 void* pythia6_common_address(const char*); << 110 void* pythia6_common_address(const char*);
112 } 111 }
113 112
114 // Direct declaration of pythia6 common blocks 113 // Direct declaration of pythia6 common blocks
115 // extern "C" { 114 // extern "C" {
116 // extern Pyjets_t pyjets_; 115 // extern Pyjets_t pyjets_;
117 // extern Pydat1_t pydat1_; 116 // extern Pydat1_t pydat1_;
118 // extern Pydat3_t pydat3_; 117 // extern Pydat3_t pydat3_;
119 // } 118 // }
120 119
121 Pythia6* Pythia6::fgInstance = 0; << 120 Pythia6* Pythia6::fgInstance = 0;
122 121
123 //....oooOO0OOooo........oooOO0OOooo........oo 122 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
124 123
125 Pythia6* Pythia6::Instance() << 124 Pythia6* Pythia6::Instance()
126 { 125 {
127 /// Static access method << 126 /// Static access method
128 127
129 if (!fgInstance) fgInstance = new Pythia6(); << 128 if ( ! fgInstance ) fgInstance = new Pythia6();
130 129
131 return fgInstance; << 130 return fgInstance;
132 } 131 }
133 132
134 //....oooOO0OOooo........oooOO0OOooo........oo 133 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
135 134
136 Pythia6::Pythia6() : fParticles(0), fPyjets(0) << 135 Pythia6::Pythia6()
>> 136 : fParticles(0),
>> 137 fPyjets(0),
>> 138 fPydat1(0),
>> 139 fPydat3(0)
137 { 140 {
138 /// Pythia6 constructor: creates a vector of << 141 /// Pythia6 constructor: creates a vector of Pythia6Particle in which it will
139 /// store all particles. Note that there may << 142 /// store all particles. Note that there may be only one functional Pythia6
140 /// object at a time, so it's not use to cre << 143 /// object at a time, so it's not use to create more than one instance of it.
141 << 144
142 // Protect against multiple objects. All a << 145 // Protect against multiple objects. All access should be via the
143 // Instance member function. << 146 // Instance member function.
144 if (fgInstance) { << 147 if ( fgInstance ) {
145 std::cerr << "There's already an instance << 148 std::cerr << "There's already an instance of Pythia6" << std::endl;
146 exit(1); << 149 exit (1);
147 } << 150 }
148 << 151
149 fParticles = new ParticleVector(); << 152 fParticles = new ParticleVector();
150 << 153
151 // Initialize common-blocks << 154 // Initialize common-blocks
152 fPyjets = (Pyjets_t*)pythia6_common_address( << 155 fPyjets = (Pyjets_t*) pythia6_common_address("PYJETS");
153 fPydat1 = (Pydat1_t*)pythia6_common_address( << 156 fPydat1 = (Pydat1_t*) pythia6_common_address("PYDAT1");
154 fPydat3 = (Pydat3_t*)pythia6_common_address( << 157 fPydat3 = (Pydat3_t*) pythia6_common_address("PYDAT3");
155 << 158
156 // Alternative way to initialize common-bloc << 159 // Alternative way to initialize common-blocks
157 // usind direct declaration of pythia6 commo << 160 // usind direct declaration of pythia6 common blocks
158 // fPyjets = &pyjets_; << 161 // fPyjets = &pyjets_;
159 // fPydat1 = &pydat1_; << 162 // fPydat1 = &pydat1_;
160 // fPydat3 = &pydat3_; << 163 // fPydat3 = &pydat3_;
161 } 164 }
162 165
163 //....oooOO0OOooo........oooOO0OOooo........oo 166 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
164 167
165 Pythia6::~Pythia6() 168 Pythia6::~Pythia6()
166 { 169 {
167 /// Destroy the object, delete and dispose a << 170 /// Destroy the object, delete and dispose all Pythia6Particles currently on
168 /// list. << 171 /// list.
169 172
170 if (fParticles) { << 173 if ( fParticles ) {
171 ParticleVector::const_iterator it; << 174 ParticleVector::const_iterator it;
172 for (it = fParticles->begin(); it != fPart << 175 for ( it = fParticles->begin(); it != fParticles->end(); it++ )
173 delete *it; << 176 delete *it;
174 delete fParticles; << 177 delete fParticles;
175 } << 178 }
176 } 179 }
177 180
178 //....oooOO0OOooo........oooOO0OOooo........oo 181 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
179 182
180 int Pythia6::Pycomp(int kf) << 183 int Pythia6::Pycomp(int kf)
181 { 184 {
182 /// Interface with fortran routine pycomp << 185 /// Interface with fortran routine pycomp
183 186
184 return pycomp(&kf); << 187 return pycomp(&kf);
185 } 188 }
186 189
187 //....oooOO0OOooo........oooOO0OOooo........oo 190 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
188 191
189 void Pythia6::Py1ent(int ip, int kf, double pe 192 void Pythia6::Py1ent(int ip, int kf, double pe, double theta, double phi)
190 { 193 {
191 /// Add one entry to the event record, i.e. << 194 /// Add one entry to the event record, i.e. either a parton or a
192 /// particle. << 195 /// particle.
193 /// << 196 ///
194 /// IP: normally line number for the part << 197 /// IP: normally line number for the parton/particle. There are two
195 /// exceptions: << 198 /// exceptions:
196 /// << 199 ///
197 /// If IP = 0: line number 1 is used << 200 /// If IP = 0: line number 1 is used and PYEXEC is called.
198 /// If IP < 0: line -IP is used, with << 201 /// If IP < 0: line -IP is used, with status code K(-IP,2)=2
199 /// rather than 1; thus a << 202 /// rather than 1; thus a parton system may be built
200 /// up by filling all but << 203 /// up by filling all but the last parton of the
201 /// system with IP < 0. << 204 /// system with IP < 0.
202 /// KF: parton/particle flavour code (PDG << 205 /// KF: parton/particle flavour code (PDG code)
203 /// PE: parton/particle energy. If PE is << 206 /// PE: parton/particle energy. If PE is smaller than the mass,
204 /// the parton/particle is taken to b << 207 /// the parton/particle is taken to be at rest.
205 /// THETA: << 208 /// THETA:
206 /// PHI: polar and azimuthal angle for the << 209 /// PHI: polar and azimuthal angle for the momentum vector of the
207 /// parton/particle. << 210 /// parton/particle.
208 211
209 py1ent(ip, kf, pe, theta, phi); << 212 py1ent(ip, kf, pe, theta, phi);
210 } 213 }
211 214
212 //....oooOO0OOooo........oooOO0OOooo........oo 215 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
213 216
214 int Pythia6::ImportParticles(ParticleVector* p 217 int Pythia6::ImportParticles(ParticleVector* particles, const char* option)
215 { 218 {
216 /// Default primary creation method. It rea << 219 /// Default primary creation method. It reads the /HEPEVT/ common block which
217 /// has been filled by the GenerateEvent me << 220 /// has been filled by the GenerateEvent method. If the event generator does
218 /// not use the HEPEVT common block, This r << 221 /// not use the HEPEVT common block, This routine has to be overloaded by
219 /// the subclasses. << 222 /// the subclasses.
220 /// The function loops on the generated par << 223 /// The function loops on the generated particles and store them in
221 /// the TClonesArray pointed by the argumen << 224 /// the TClonesArray pointed by the argument particles.
222 /// The default action is to store only the << 225 /// The default action is to store only the stable particles (ISTHEP = 1)
223 /// This can be demanded explicitly by sett << 226 /// This can be demanded explicitly by setting the option = "Final"
224 /// If the option = "All", all the particle << 227 /// If the option = "All", all the particles are stored.
225 << 228
226 if (particles == 0) return 0; << 229 if ( particles == 0 ) return 0;
227 << 230
228 ParticleVector::const_iterator it; << 231 ParticleVector::const_iterator it;
229 for (it = particles->begin(); it != particle << 232 for ( it = particles->begin(); it != particles->end(); it++ )
230 delete *it; << 233 delete *it;
231 particles->clear(); << 234 particles->clear();
232 << 235
233 int numpart = fPyjets->N; << 236 int numpart = fPyjets->N;
234 int nparts = 0; << 237 int nparts=0;
235 if (!strcmp(option, "") || !strcmp(option, " << 238 if (!strcmp(option,"") || !strcmp(option,"Final")) {
236 for (int i = 0; i < numpart; i++) { << 239 for (int i = 0; i<numpart; i++) {
237 if (fPyjets->K[0][i] == 1) { << 240
238 // << 241 if (fPyjets->K[0][i] == 1) {
239 // Use the common block values for th << 242 //
240 // << 243 // Use the common block values for the TParticle constructor
241 particles->push_back(new Pythia6Partic << 244 //
242 fPyjets->K[0][i], fPyjets->K[1][i], << 245 particles->push_back(
243 fPyjets->P[0][i], fPyjets->P[1][i], << 246 new Pythia6Particle(
244 fPyjets->V[0][i], fPyjets->V[1][i], << 247 fPyjets->K[0][i] ,
245 fPyjets->V[4][i])); << 248 fPyjets->K[1][i] ,
246 << 249 fPyjets->K[2][i] ,
247 // if(gDebug) printf("%d %d %d! ", << 250 fPyjets->K[3][i] ,
248 nparts++; << 251 fPyjets->K[4][i] ,
>> 252 fPyjets->P[0][i] ,
>> 253 fPyjets->P[1][i] ,
>> 254 fPyjets->P[2][i] ,
>> 255 fPyjets->P[3][i] ,
>> 256 fPyjets->P[4][i] ,
>> 257 fPyjets->V[0][i] ,
>> 258 fPyjets->V[1][i] ,
>> 259 fPyjets->V[2][i] ,
>> 260 fPyjets->V[3][i] ,
>> 261 fPyjets->V[4][i]));
>> 262
>> 263 // if(gDebug) printf("%d %d %d! ",i,fPyjets->K[1][i],numpart);
>> 264 nparts++;
>> 265 }
>> 266 }
>> 267 }
>> 268 else if (!strcmp(option,"All")) {
>> 269 for (int i = 0; i<numpart; i++) {
>> 270 particles->push_back(
>> 271 new Pythia6Particle(
>> 272 fPyjets->K[0][i] ,
>> 273 fPyjets->K[1][i] ,
>> 274 fPyjets->K[2][i] ,
>> 275 fPyjets->K[3][i] ,
>> 276 fPyjets->K[4][i] ,
>> 277 fPyjets->P[0][i] ,
>> 278 fPyjets->P[1][i] ,
>> 279 fPyjets->P[2][i] ,
>> 280 fPyjets->P[3][i] ,
>> 281 fPyjets->P[4][i] ,
>> 282 fPyjets->V[0][i] ,
>> 283 fPyjets->V[1][i] ,
>> 284 fPyjets->V[2][i] ,
>> 285 fPyjets->V[3][i] ,
>> 286 fPyjets->V[4][i]));
249 } 287 }
250 } << 288 nparts=numpart;
251 } << 289 }
252 else if (!strcmp(option, "All")) { <<
253 for (int i = 0; i < numpart; i++) { <<
254 particles->push_back(new Pythia6Particle <<
255 fPyjets->K[0][i], fPyjets->K[1][i], fP <<
256 fPyjets->P[0][i], fPyjets->P[1][i], fP <<
257 fPyjets->V[0][i], fPyjets->V[1][i], fP <<
258 } <<
259 nparts = numpart; <<
260 } <<
261 290
262 return nparts; << 291 return nparts;
263 } 292 }
264 293
265 //....oooOO0OOooo........oooOO0OOooo........oo 294 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
266 295