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>
48 ClassImp(TTherminator)
50 TTherminator::TTherminator():
55 // Default constructor
56 fPartDB = new ParticleDB();
58 TTherminator::TTherminator(const TTherminator & therm) :
65 fPartDB = new ParticleDB();
67 TTherminator::~TTherminator()
70 if (sRPInstance) delete sRPInstance;
71 if (fEvent) delete fEvent;
72 if (fCalka) delete fCalka;
73 if (fPartDB) delete fPartDB;
76 void TTherminator::ReadParameters()
78 // Read in global model parameters
83 sNumEvents = atoi(sRPInstance->getPar("NumberOfEvents").Data());
84 tModel = sRPInstance->getPar("FreezeOutModel");
85 if (tModel.Contains("SingleFreezeOut"))
88 else if (tModel.Contains("Lhyquid3D"))
90 else if (tModel.Contains("Lhyquid2D"))
93 else if (tModel.Contains("BlastWaveVTDelay"))
95 else if (tModel.Contains("BlastWaveVT"))
97 else if (tModel.Contains("BlastWaveVLinearFormation"))
99 else if (tModel.Contains("BlastWaveVLinearDelay"))
101 else if (tModel.Contains("BlastWaveVLinear"))
103 else if (tModel.Contains("FiniteHydro"))
106 PRINT_MESSAGE("Unknown model type: " << tModel.Data());
107 PRINT_MESSAGE("Please provide the proper name");
110 if (atoi(sRPInstance->getPar("Randomize").Data()) == 1)
114 tTable = sRPInstance->getPar("TableType");
115 if (tTable.Contains("Mathematica"))
117 else if (tTable.Contains("SHARE"))
120 PRINT_MESSAGE("Unknown table type: " << tTable.Data());
123 sIntegrateSample = atoi(sRPInstance->getPar("NumberOfIntegrateSamples").Data());
124 fInputDir = sRPInstance->getPar("InputDirSHARE");
127 PRINT_DEBUG_1("RunBFPW::ReadParameters - Caught exception " << tError);
128 PRINT_MESSAGE("Did not find one of the neccessary parameters in the parameters file.");
133 void TTherminator::Initialize(){
134 // Initialize global variables
135 // and particle and decay tables
138 sRPFileName = "therminator.in";
139 sRPInstance = new ReadPar(sRPFileName.Data());
143 tParser = new Parser(fPartDB);
145 tParser->ReadShare();
147 fCalka = new Integrator(sIntegrateSample);
148 fCalka->ReadMultInteg(fPartDB);
151 // Read in particle table from share
152 // and add missing particle type to TDatabasePDG
157 // ParticleType *tPartBuf;
159 TDatabasePDG *tInstance = TDatabasePDG::Instance();
160 TParticlePDG *tParticleType;
162 // AliWarning(Form("Reading particle types from particles.data"));
164 ifstream in((fInputDir+"/"+"particles.data").Data());
165 // ifstream in("particles.data");
170 if ((in) && (in.is_open()))
180 if (/*(*str == '#')||*/(*str<65)||(*str>122))
182 in.getline(str1,200);
185 double mass, gamma, spin, tI3, tI, q, s, aq, as, c, ac, mc;
187 in>>mass>>gamma>>spin>>tI>>tI3>>q>>s>>aq>>as>>c>>ac>>mc;
189 tParticleType = tInstance->GetParticle((int) mc);
190 if (!tParticleType) {
192 if (strstr(str, "plu")) charge = 1;
193 if (strstr(str, "min")) charge = -1;
194 if (strstr(str, "plb")) charge = -1;
195 if (strstr(str, "mnb")) charge = 1;
196 if (strstr(str, "plp")) charge = 2;
197 if (strstr(str, "ppb")) charge = -2;
198 tInstance->AddParticle(str, str, mass, gamma == 0.0 ? 1:0, gamma, charge , "meson", (int) mc);
199 // printf("Added particle %s with PDG PID %d charge %d", str, (int) mc, charge);
200 //AliWarning(Form("Added particle %s with PDG PID %d charge %d", str, (int) mc, charge));
201 // 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));
211 void TTherminator::GenerateEvent()
213 // Run single event generation
214 if ((sRunType == 3) || (sRunType == 4))
219 if (!fEvent) fEvent = new Event(fPartDB, fCalka);
220 if (sRunType == 4) fEvent->Randomize();
222 fEvent->GenerateEvent();
223 fEvent->DecayParticles();
224 // fEvent->WriteEvent(0);
228 Int_t TTherminator::ImportParticles(TClonesArray *particles, Option_t */*option*/)
230 // Import particles from a generated event into an external array
231 const double kFmToGev = 0.197327;
233 if (particles == 0) return 0;
234 TClonesArray &particlesR = *particles;
237 if (!fEvent) return 0;
238 Int_t numpart = fEvent->GetParticleCount();
239 printf("\n TTherminator: Therminator stack contains %d particles.\n", numpart);
240 for (Int_t iPart=0; iPart<numpart; iPart++) {
241 Particle *tPart = fEvent->GetParticleOfCount(iPart);
243 tFather = tPart->GetFather();
246 TParticle *mother = (TParticle*) (particlesR.UncheckedAt(tFather+1));
247 mother->SetLastDaughter(iPart);
248 if (mother->GetFirstDaughter()==-1)
249 mother->SetFirstDaughter(iPart);
252 // printf("Putting %d %d %lf %d\n", tPart->GetParticleType()->GetPDGCode(), iPart, tPart->px, tPart);
253 new (particlesR[iPart]) TParticle(tPart->GetParticleType()->GetPDGCode(), tPart->HadDecayed(),
255 tPart->px, tPart->py, tPart->pz, tPart->GetEnergy() ,
256 tPart->rx*1.e-13*kFmToGev, tPart->ry*1.e-13*kFmToGev, tPart->rz*1.e-13*kFmToGev, tPart->rt*1.e-13*kFmToGev/3e10);
257 particlesR[iPart]->SetUniqueID(iPart);
263 TObjArray* TTherminator::ImportParticles(Option_t */*option*/)
265 // Import generated particles into an internal array
266 const double kFmToGev = 0.197327;
270 if (!fEvent) return 0;
271 Int_t numpart = fEvent->GetParticleCount();
272 printf("\n TTherminator: Therminator stack contains %d particles.", numpart);
273 fEvent->InitParticleScan();
274 for (Int_t iPart=0; iPart<numpart; iPart++) {
275 Particle *tPart = fEvent->GetNextParticle();
277 tFather = tPart->GetFather();
280 TParticle *mother = (TParticle*) (fParticles->UncheckedAt(tFather));
281 mother->SetLastDaughter(iPart);
282 if (mother->GetFirstDaughter()==-1)
283 mother->SetFirstDaughter(iPart);
285 TParticle* p = new TParticle(tPart->GetParticleType()->GetPDGCode(), tPart->HadDecayed(),
287 tPart->px, tPart->py, tPart->pz, tPart->GetEnergy() ,
288 tPart->rx*1.e-13*kFmToGev, tPart->ry*1.e-13*kFmToGev, tPart->rz*1.e-13*kFmToGev, tPart->rt*1.e-13*kFmToGev/3e10);
289 p->SetUniqueID(iPart);
292 nump = fEvent->GetParticleCount();