1 ///////////////////////////////////////////////////////////////////////////////
3 // TTherminator: global interface class to the Therminator model. //
4 // Initialized global variables and runs the event generation //
6 ///////////////////////////////////////////////////////////////////////////////
7 /******************************************************************************
8 * T H E R M I N A T O R *
9 * THERMal heavy-IoN generATOR *
12 * Authors of the model: Wojciech Broniowski, Wojciech.Broniowski@ifj.edu.pl, *
13 * Wojciech Florkowski, Wojciech.Florkowski@ifj.edu.pl *
14 * Authors of the code: Adam Kisiel, kisiel@if.pw.edu.pl *
15 * Tomasz Taluc, ttaluc@if.pw.edu.pl *
16 * Code designers: Adam Kisiel, Tomasz Taluc, Wojciech Broniowski, *
17 * Wojciech Florkowski *
19 * For the detailed description of the program and furhter references *
20 * to the description of the model plesase refer to: nucl-th/0504047, *
21 * accessibile at: http://www.arxiv.org/nucl-th/0504047 *
23 * Homepage: http://hirg.if.pw.edu.pl/en/therminator/ *
25 * This code can be freely used and redistributed. However if you decide to *
26 * make modifications to the code, please contact the authors, especially *
27 * if you plan to publish the results obtained with such modified code. *
28 * Any publication of results obtained using this code must include the *
29 * reference to nucl-th/0504047 and the published version of it, when *
32 *****************************************************************************/
33 #include <TTherminator.h>
35 #include <TDatabasePDG.h>
36 #include <TParticle.h>
37 #include <TClonesArray.h>
47 ClassImp(TTherminator)
49 TTherminator::TTherminator():
54 // Default constructor
55 fPartDB = new ParticleDB();
57 TTherminator::TTherminator(const TTherminator & therm) :
63 fPartDB = new ParticleDB();
65 TTherminator::~TTherminator()
68 if (sRPInstance) delete sRPInstance;
69 if (fEvent) delete fEvent;
70 if (fCalka) delete fCalka;
71 if (fPartDB) delete fPartDB;
74 void TTherminator::ReadParameters()
76 // Read in global model parameters
81 sNumEvents = atoi(sRPInstance->getPar("NumberOfEvents").Data());
82 tModel = sRPInstance->getPar("FreezeOutModel");
83 if (tModel.Contains("SingleFreezeOut"))
86 else if (tModel.Contains("Lhyquid3D"))
88 else if (tModel.Contains("Lhyquid2D"))
91 else if (tModel.Contains("BlastWaveVTDelay"))
93 else if (tModel.Contains("BlastWaveVT"))
95 else if (tModel.Contains("BlastWaveVLinearFormation"))
97 else if (tModel.Contains("BlastWaveVLinearDelay"))
99 else if (tModel.Contains("BlastWaveVLinear"))
101 else if (tModel.Contains("FiniteHydro"))
104 PRINT_MESSAGE("Unknown model type: " << tModel.Data());
105 PRINT_MESSAGE("Please provide the proper name");
108 if (atoi(sRPInstance->getPar("Randomize").Data()) == 1)
112 tTable = sRPInstance->getPar("TableType");
113 if (tTable.Contains("Mathematica"))
115 else if (tTable.Contains("SHARE"))
118 PRINT_MESSAGE("Unknown table type: " << tTable.Data());
121 sIntegrateSample = atoi(sRPInstance->getPar("NumberOfIntegrateSamples").Data());
122 fInputDir = sRPInstance->getPar("InputDirSHARE");
125 PRINT_DEBUG_1("RunBFPW::ReadParameters - Caught exception " << tError);
126 PRINT_MESSAGE("Did not find one of the neccessary parameters in the parameters file.");
131 void TTherminator::Initialize(){
132 // Initialize global variables
133 // and particle and decay tables
136 sRPFileName = "therminator.in";
137 sRPInstance = new ReadPar(sRPFileName.Data());
141 tParser = new Parser(fPartDB);
143 tParser->ReadShare();
145 fCalka = new Integrator(sIntegrateSample);
146 fCalka->ReadMultInteg(fPartDB);
149 // Read in particle table from share
150 // and add missing particle type to TDatabasePDG
155 // ParticleType *tPartBuf;
157 TDatabasePDG *tInstance = TDatabasePDG::Instance();
158 TParticlePDG *tParticleType;
160 // AliWarning(Form("Reading particle types from particles.data"));
162 ifstream in((fInputDir+"/"+"particles.data").Data());
163 // ifstream in("particles.data");
168 if ((in) && (in.is_open()))
178 if (/*(*str == '#')||*/(*str<65)||(*str>122))
180 in.getline(str1,200);
183 double mass, gamma, spin, tI3, tI, q, s, aq, as, c, ac, mc;
185 in>>mass>>gamma>>spin>>tI>>tI3>>q>>s>>aq>>as>>c>>ac>>mc;
187 tParticleType = tInstance->GetParticle((int) mc);
188 if (!tParticleType) {
190 if (strstr(str, "plu")) charge = 1;
191 if (strstr(str, "min")) charge = -1;
192 if (strstr(str, "plb")) charge = -1;
193 if (strstr(str, "mnb")) charge = 1;
194 if (strstr(str, "plp")) charge = 2;
195 if (strstr(str, "ppb")) charge = -2;
196 tInstance->AddParticle(str, str, mass, gamma == 0.0 ? 1:0, gamma, charge , "meson", (int) mc);
197 // printf("Added particle %s with PDG PID %d charge %d", str, (int) mc, charge);
198 //AliWarning(Form("Added particle %s with PDG PID %d charge %d", str, (int) mc, charge));
199 // AliWarning(Form("Quantum numbers q s c aq as ac tI3 %lf %lf %lf %lf %lf %lf %lf", q, s, c, aq, as, ac, tI3));
207 void TTherminator::GenerateEvent()
209 // Run single event generation
210 if ((sRunType == 3) || (sRunType == 4))
215 if (!fEvent) fEvent = new Event(fPartDB, fCalka);
216 if (sRunType == 4) fEvent->Randomize();
218 fEvent->GenerateEvent();
219 fEvent->DecayParticles();
220 // fEvent->WriteEvent(0);
224 Int_t TTherminator::ImportParticles(TClonesArray *particles, Option_t *option)
226 // Import particles from a generated event into an external array
227 if (particles == 0) return 0;
228 TClonesArray &particlesR = *particles;
231 if (!fEvent) return 0;
232 Int_t numpart = fEvent->GetParticleCount();
233 printf("\n TTherminator: Therminator stack contains %d particles.", numpart);
234 for (Int_t iPart=0; iPart<numpart; iPart++) {
235 Particle *tPart = fEvent->GetParticleOfCount(iPart);
237 tFather = tPart->GetFather();
240 TParticle *mother = (TParticle*) (particlesR.UncheckedAt(tFather+1));
241 mother->SetLastDaughter(iPart);
242 if (mother->GetFirstDaughter()==-1)
243 mother->SetFirstDaughter(iPart);
246 // printf("Putting %d %d %lf %d\n", tPart->GetParticleType()->GetPDGCode(), iPart, tPart->px, tPart);
247 new (particlesR[iPart]) TParticle(tPart->GetParticleType()->GetPDGCode(), tPart->HadDecayed(),
249 tPart->px, tPart->py, tPart->pz, tPart->GetEnergy() ,
250 tPart->rx*1.e-13, tPart->ry*1.e-13, tPart->rz*1.e-13, tPart->rt*1.e-13/300000000
252 particlesR[iPart]->SetUniqueID(iPart);
258 TObjArray* TTherminator::ImportParticles(Option_t *option)
260 // Import generated particles into an internal array
263 if (!fEvent) return 0;
264 Int_t numpart = fEvent->GetParticleCount();
265 printf("\n TTherminator: Therminator stack contains %d particles.", numpart);
266 fEvent->InitParticleScan();
267 for (Int_t iPart=0; iPart<numpart; iPart++) {
268 Particle *tPart = fEvent->GetNextParticle();
270 tFather = tPart->GetFather();
273 TParticle *mother = (TParticle*) (fParticles->UncheckedAt(tFather));
274 mother->SetLastDaughter(iPart);
275 if (mother->GetFirstDaughter()==-1)
276 mother->SetFirstDaughter(iPart);
278 TParticle* p = new TParticle(tPart->GetParticleType()->GetPDGCode(), tPart->HadDecayed(),
280 tPart->px, tPart->py, tPart->pz, tPart->GetEnergy() ,
281 tPart->rx*1.e-13, tPart->ry*1.e-13, tPart->rz*1.e-13, tPart->rt*1.e-13/300000000
283 p->SetUniqueID(iPart);
286 nump = fEvent->GetParticleCount();