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) :
64 fPartDB = new ParticleDB();
66 TTherminator::~TTherminator()
69 if (sRPInstance) delete sRPInstance;
70 if (fEvent) delete fEvent;
71 if (fCalka) delete fCalka;
72 if (fPartDB) delete fPartDB;
75 void TTherminator::ReadParameters()
77 // Read in global model parameters
82 sNumEvents = atoi(sRPInstance->getPar("NumberOfEvents").Data());
83 tModel = sRPInstance->getPar("FreezeOutModel");
84 if (tModel.Contains("SingleFreezeOut"))
87 else if (tModel.Contains("Lhyquid3D"))
89 else if (tModel.Contains("Lhyquid2D"))
92 else if (tModel.Contains("BlastWaveVTDelay"))
94 else if (tModel.Contains("BlastWaveVT"))
96 else if (tModel.Contains("BlastWaveVLinearFormation"))
98 else if (tModel.Contains("BlastWaveVLinearDelay"))
100 else if (tModel.Contains("BlastWaveVLinear"))
102 else if (tModel.Contains("FiniteHydro"))
105 PRINT_MESSAGE("Unknown model type: " << tModel.Data());
106 PRINT_MESSAGE("Please provide the proper name");
109 if (atoi(sRPInstance->getPar("Randomize").Data()) == 1)
113 tTable = sRPInstance->getPar("TableType");
114 if (tTable.Contains("Mathematica"))
116 else if (tTable.Contains("SHARE"))
119 PRINT_MESSAGE("Unknown table type: " << tTable.Data());
122 sIntegrateSample = atoi(sRPInstance->getPar("NumberOfIntegrateSamples").Data());
123 fInputDir = sRPInstance->getPar("InputDirSHARE");
126 PRINT_DEBUG_1("RunBFPW::ReadParameters - Caught exception " << tError);
127 PRINT_MESSAGE("Did not find one of the neccessary parameters in the parameters file.");
132 void TTherminator::Initialize(){
133 // Initialize global variables
134 // and particle and decay tables
137 sRPFileName = "therminator.in";
138 sRPInstance = new ReadPar(sRPFileName.Data());
142 tParser = new Parser(fPartDB);
144 tParser->ReadShare();
146 fCalka = new Integrator(sIntegrateSample);
147 fCalka->ReadMultInteg(fPartDB);
150 // Read in particle table from share
151 // and add missing particle type to TDatabasePDG
156 // ParticleType *tPartBuf;
158 TDatabasePDG *tInstance = TDatabasePDG::Instance();
159 TParticlePDG *tParticleType;
161 // AliWarning(Form("Reading particle types from particles.data"));
163 ifstream in((fInputDir+"/"+"particles.data").Data());
164 // ifstream in("particles.data");
169 if ((in) && (in.is_open()))
179 if (/*(*str == '#')||*/(*str<65)||(*str>122))
181 in.getline(str1,200);
184 double mass, gamma, spin, tI3, tI, q, s, aq, as, c, ac, mc;
186 in>>mass>>gamma>>spin>>tI>>tI3>>q>>s>>aq>>as>>c>>ac>>mc;
188 tParticleType = tInstance->GetParticle((int) mc);
189 if (!tParticleType) {
191 if (strstr(str, "plu")) charge = 1;
192 if (strstr(str, "min")) charge = -1;
193 if (strstr(str, "plb")) charge = -1;
194 if (strstr(str, "mnb")) charge = 1;
195 if (strstr(str, "plp")) charge = 2;
196 if (strstr(str, "ppb")) charge = -2;
197 tInstance->AddParticle(str, str, mass, gamma == 0.0 ? 1:0, gamma, charge , "meson", (int) mc);
198 // printf("Added particle %s with PDG PID %d charge %d", str, (int) mc, charge);
199 //AliWarning(Form("Added particle %s with PDG PID %d charge %d", str, (int) mc, charge));
200 // 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));
208 void TTherminator::GenerateEvent()
210 // Run single event generation
211 if ((sRunType == 3) || (sRunType == 4))
216 if (!fEvent) fEvent = new Event(fPartDB, fCalka);
217 if (sRunType == 4) fEvent->Randomize();
219 fEvent->GenerateEvent();
220 fEvent->DecayParticles();
221 // fEvent->WriteEvent(0);
225 Int_t TTherminator::ImportParticles(TClonesArray *particles, Option_t */*option*/)
227 // Import particles from a generated event into an external array
228 const double kFmToGev = 0.197327;
230 if (particles == 0) return 0;
231 TClonesArray &particlesR = *particles;
234 if (!fEvent) return 0;
235 Int_t numpart = fEvent->GetParticleCount();
236 printf("\n TTherminator: Therminator stack contains %d particles.", numpart);
237 for (Int_t iPart=0; iPart<numpart; iPart++) {
238 Particle *tPart = fEvent->GetParticleOfCount(iPart);
240 tFather = tPart->GetFather();
243 TParticle *mother = (TParticle*) (particlesR.UncheckedAt(tFather+1));
244 mother->SetLastDaughter(iPart);
245 if (mother->GetFirstDaughter()==-1)
246 mother->SetFirstDaughter(iPart);
249 // printf("Putting %d %d %lf %d\n", tPart->GetParticleType()->GetPDGCode(), iPart, tPart->px, tPart);
250 new (particlesR[iPart]) TParticle(tPart->GetParticleType()->GetPDGCode(), tPart->HadDecayed(),
252 tPart->px, tPart->py, tPart->pz, tPart->GetEnergy() ,
253 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);
254 particlesR[iPart]->SetUniqueID(iPart);
260 TObjArray* TTherminator::ImportParticles(Option_t */*option*/)
262 // Import generated particles into an internal array
263 const double kFmToGev = 0.197327;
267 if (!fEvent) return 0;
268 Int_t numpart = fEvent->GetParticleCount();
269 printf("\n TTherminator: Therminator stack contains %d particles.", numpart);
270 fEvent->InitParticleScan();
271 for (Int_t iPart=0; iPart<numpart; iPart++) {
272 Particle *tPart = fEvent->GetNextParticle();
274 tFather = tPart->GetFather();
277 TParticle *mother = (TParticle*) (fParticles->UncheckedAt(tFather));
278 mother->SetLastDaughter(iPart);
279 if (mother->GetFirstDaughter()==-1)
280 mother->SetFirstDaughter(iPart);
282 TParticle* p = new TParticle(tPart->GetParticleType()->GetPDGCode(), tPart->HadDecayed(),
284 tPart->px, tPart->py, tPart->pz, tPart->GetEnergy() ,
285 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);
286 p->SetUniqueID(iPart);
289 nump = fEvent->GetParticleCount();