1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 // Author: Andreas Morsch 27/10/2007
17 ////////////////////////////////////////////////////////////////////////////////
21 // TPythia is an interface class to C++ version of Pythia 8.1 //
22 // event generators, written by T.Sjostrand. //
24 // The user is assumed to be familiar with the Pythia package. //
25 // This class includes only a basic interface to Pythia8. Because Pythia8 is //
26 // also written in C++, its functions/classes can be called directly from a //
27 // compiled C++ script. //
28 // To call Pythia functions not available in this interface a dictionary must //
30 // see $ROOTSYS/tutorials/pythia/pythia8.C for an example of use from CINT. //
32 ////////////////////////////////////////////////////////////////////////////////
34 *------------------------------------------------------------------------------------*
36 | *------------------------------------------------------------------------------* |
39 | | PPP Y Y TTTTT H H III A Welcome to the Lund Monte Carlo! | |
40 | | P P Y Y T H H I A A This is PYTHIA version 8.100 | |
41 | | PPP Y T HHHHH I AAAAA Last date of change: 20 Oct 2007 | |
42 | | P Y T H H I A A | |
43 | | P Y T H H III A A Now is 27 Oct 2007 at 18:26:53 | |
45 | | Main author: Torbjorn Sjostrand; CERN/PH, CH-1211 Geneva, Switzerland, | |
46 | | and Department of Theoretical Physics, Lund University, Lund, Sweden; | |
47 | | phone: + 41 - 22 - 767 82 27; e-mail: torbjorn@thep.lu.se | |
48 | | Author: Stephen Mrenna; Computing Division, Simulations Group, | |
49 | | Fermi National Accelerator Laboratory, MS 234, Batavia, IL 60510, USA; | |
50 | | phone: + 1 - 630 - 840 - 2556; e-mail: mrenna@fnal.gov | |
51 | | Author: Peter Skands; CERN/PH, CH-1211 Geneva, Switzerland, | |
52 | | and Theoretical Physics Department, | |
53 | | Fermi National Accelerator Laboratory, MS 106, Batavia, IL 60510, USA; | |
54 | | phone: + 41 - 22 - 767 24 59; e-mail: skands@fnal.gov | |
56 | | The main program reference is the 'Brief Introduction to PYTHIA 8.1', | |
57 | | T. Sjostrand, S. Mrenna and P. Skands, arXiv:0710.3820 | |
59 | | The main physics reference is the 'PYTHIA 6.4 Physics and Manual', | |
60 | | T. Sjostrand, S. Mrenna and P. Skands, JHEP05 (2006) 026 [hep-ph/0603175]. | |
62 | | An archive of program versions and documentation is found on the web: | |
63 | | http://www.thep.lu.se/~torbjorn/Pythia.html | |
65 | | This program is released under the GNU General Public Licence version 2. | |
66 | | Please respect the MCnet Guidelines for Event Generator Authors and Users. | |
68 | | Disclaimer: this program comes without any guarantees. | |
69 | | Beware of errors and use common sense when interpreting results. | |
71 | | Copyright (C) 2007 Torbjorn Sjostrand | |
74 | *------------------------------------------------------------------------------* |
76 *------------------------------------------------------------------------------------*
79 #include "AliTPythia8.h"
81 #include "TClonesArray.h"
82 #include "TParticle.h"
83 #include "TDatabasePDG.h"
84 #include "TLorentzVector.h"
88 AliTPythia8* AliTPythia8::fgInstance = 0;
89 char* AliTPythia8::fgXmldocPath = 0;
91 //___________________________________________________________________________
92 AliTPythia8::AliTPythia8():
93 TGenerator("AliTPythia8", "AliTPythia8"),
99 Fatal("AliTPythia8", "There's already an instance of AliTPythia8");
101 delete fParticles; // was allocated as TObjArray in TGenerator
103 fParticles = new TClonesArray("TParticle",50);
104 if (fgXmldocPath != 0) {
105 fPythia = new Pythia8::Pythia(fgXmldocPath);
107 fPythia = new Pythia8::Pythia();
111 //___________________________________________________________________________
112 AliTPythia8::AliTPythia8(const char *xmlDir):
113 TGenerator("AliTPythia8", "AliTPythia8"),
115 fNumberOfParticles(0)
117 // Constructor with an xmlDir (eg "../xmldoc"
119 Fatal("AliTPythia8", "There's already an instance of AliTPythia8");
121 delete fParticles; // was allocated as TObjArray in TGenerator
123 fParticles = new TClonesArray("TParticle",50);
124 fPythia = new Pythia8::Pythia(xmlDir);
127 //___________________________________________________________________________
128 AliTPythia8::~AliTPythia8()
132 fParticles->Delete();
139 //___________________________________________________________________________
140 AliTPythia8* AliTPythia8::Instance()
142 // Return an instance of AliTPythia8
143 return fgInstance ? fgInstance : (fgInstance = new AliTPythia8()) ;
146 //___________________________________________________________________________
147 Bool_t AliTPythia8::Initialize(Int_t idAin, Int_t idBin, Double_t ecms)
150 AddParticlesToPdgDataBase();
151 return fPythia->init(idAin, idBin, ecms);
154 //___________________________________________________________________________
155 void AliTPythia8::GenerateEvent()
157 // Generate the next event
159 fNumberOfParticles = fPythia->event.size() - 1;
163 //___________________________________________________________________________
164 Int_t AliTPythia8::ImportParticles(TClonesArray *particles, Option_t *option)
166 // Import particles from Pythia stack
167 if (particles == 0) return 0;
168 TClonesArray &clonesParticles = *particles;
169 clonesParticles.Clear();
173 fNumberOfParticles = fPythia->event.size();
174 if (fPythia->event[0].id() == 90) {
178 if (!strcmp(option,"") || !strcmp(option,"Final")) {
179 for (i = 0; i < fNumberOfParticles; i++) {
180 if (fPythia->event[i].id() == 90) continue;
181 if (fPythia->event[i].isFinal()) {
182 new(clonesParticles[nparts]) TParticle(
183 fPythia->event[i].id(),
184 fPythia->event[i].isFinal(),
185 fPythia->event[i].mother1() + ioff,
186 fPythia->event[i].mother2() + ioff,
187 fPythia->event[i].daughter1() + ioff,
188 fPythia->event[i].daughter2() + ioff,
189 fPythia->event[i].px(), // [GeV/c]
190 fPythia->event[i].py(), // [GeV/c]
191 fPythia->event[i].pz(), // [GeV/c]
192 fPythia->event[i].e(), // [GeV]
193 fPythia->event[i].xProd(), // [mm]
194 fPythia->event[i].yProd(), // [mm]
195 fPythia->event[i].zProd(), // [mm]
196 fPythia->event[i].tProd()); // [mm/c]
198 } // final state partice
200 } else if (!strcmp(option,"All")) {
201 for (i = 0; i < fNumberOfParticles; i++) {
202 if (fPythia->event[i].id() == 90) continue;
203 new(clonesParticles[nparts]) TParticle(
204 fPythia->event[i].id(),
205 fPythia->event[i].isFinal(),
206 fPythia->event[i].mother1() + ioff,
207 fPythia->event[i].mother2() + ioff,
208 fPythia->event[i].daughter1() + ioff,
209 fPythia->event[i].daughter2() + ioff,
210 fPythia->event[i].px(), // [GeV/c]
211 fPythia->event[i].py(), // [GeV/c]
212 fPythia->event[i].pz(), // [GeV/c]
213 fPythia->event[i].e(), // [GeV]
214 fPythia->event[i].xProd(), // [mm]
215 fPythia->event[i].yProd(), // [mm]
216 fPythia->event[i].zProd(), // [mm]
217 fPythia->event[i].tProd()); // [mm/c]
221 if(ioff==-1) fNumberOfParticles--;
225 //___________________________________________________________________________
226 TObjArray* AliTPythia8::ImportParticles(Option_t* /* option */)
228 // Import particles from Pythia stack
232 fNumberOfParticles = fPythia->event.size();
233 if (fPythia->event[0].id() == 90) {
238 TClonesArray &a = *((TClonesArray*)fParticles);
239 for (Int_t i = 0; i < fNumberOfParticles; i++) {
240 if (fPythia->event[i].id() == 90) continue;
241 new(a[nparts]) TParticle(
242 fPythia->event[i].id(),
243 fPythia->event[i].isFinal(),
244 fPythia->event[i].mother1() + ioff,
245 fPythia->event[i].mother2() + ioff,
246 fPythia->event[i].daughter1() + ioff,
247 fPythia->event[i].daughter2() + ioff,
248 fPythia->event[i].px(), // [GeV/c]
249 fPythia->event[i].py(), // [GeV/c]
250 fPythia->event[i].pz(), // [GeV/c]
251 fPythia->event[i].e(), // [GeV]
252 fPythia->event[i].xProd(), // [mm]
253 fPythia->event[i].yProd(), // [mm]
254 fPythia->event[i].zProd(), // [mm]
255 fPythia->event[i].tProd()); // [mm/c]
258 if (ioff == -1) fNumberOfParticles--;
262 //___________________________________________________________________________
263 Int_t AliTPythia8::GetN() const
266 return (fPythia->event.size() - 1);
269 //___________________________________________________________________________
270 void AliTPythia8::ReadString(const char* string) const
273 fPythia->readString(string);
276 //___________________________________________________________________________
277 void AliTPythia8::ReadConfigFile(const char* string) const
280 fPythia->readFile(string);
283 //___________________________________________________________________________
284 void AliTPythia8::PrintStatistics() const
286 // Print end of run statistics
287 fPythia->statistics();
290 //___________________________________________________________________________
291 void AliTPythia8::EventListing() const
294 fPythia->event.list();
297 //___________________________________________________________________________
298 void AliTPythia8::AddParticlesToPdgDataBase() const
300 // Add some pythia specific particle code to the data base
302 TDatabasePDG *pdgDB = TDatabasePDG::Instance();
303 pdgDB->AddParticle("string","string", 0, kTRUE,
304 0, 0, "QCD string", 90);
305 pdgDB->AddParticle("rho_diff0", "rho_diff0", 0, kTRUE,
306 0, 0, "QCD diffr. state", 9900110);
307 pdgDB->AddParticle("pi_diffr+", "pi_diffr+", 0, kTRUE,
308 0, 1, "QCD diffr. state", 9900210);
309 pdgDB->AddParticle("omega_di", "omega_di", 0, kTRUE,
310 0, 0, "QCD diffr. state", 9900220);
311 pdgDB->AddParticle("phi_diff","phi_diff", 0, kTRUE,
312 0, 0, "QCD diffr. state", 9900330);
313 pdgDB->AddParticle("J/psi_di", "J/psi_di", 0, kTRUE,
314 0, 0, "QCD diffr. state", 9900440);
315 pdgDB->AddParticle("n_diffr0","n_diffr0",0,kTRUE,
316 0, 0, "QCD diffr. state", 9902110);
317 pdgDB->AddParticle("p_diffr+","p_diffr+", 0, kTRUE,
318 0, 1, "QCD diffr. state", 9902210);