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());
124 PRINT_DEBUG_1("RunBFPW::ReadParameters - Caught exception " << tError);
125 PRINT_MESSAGE("Did not find one of the neccessary parameters in the parameters file.");
130 void TTherminator::Initialize(){
131 // Initialize global variables
132 // and particle and decay tables
135 sRPFileName = "therminator.in";
136 sRPInstance = new ReadPar(sRPFileName.Data());
140 tParser = new Parser(fPartDB);
142 tParser->ReadShare();
144 fCalka = new Integrator(sIntegrateSample);
145 fCalka->ReadMultInteg(fPartDB);
148 // Read in particle table from share
149 // and add missing particle type to TDatabasePDG
154 // ParticleType *tPartBuf;
156 TDatabasePDG *tInstance = TDatabasePDG::Instance();
157 TParticlePDG *tParticleType;
159 // AliWarning(Form("Reading particle types from particles.data"));
160 ifstream in("particles.data");
165 if ((in) && (in.is_open()))
175 if (/*(*str == '#')||*/(*str<65)||(*str>122))
177 in.getline(str1,200);
180 double mass, gamma, spin, tI3, tI, q, s, aq, as, c, ac, mc;
182 in>>mass>>gamma>>spin>>tI>>tI3>>q>>s>>aq>>as>>c>>ac>>mc;
184 tParticleType = tInstance->GetParticle((int) mc);
185 if (!tParticleType) {
187 if (strstr(str, "plu")) charge = 1;
188 if (strstr(str, "min")) charge = -1;
189 if (strstr(str, "plb")) charge = -1;
190 if (strstr(str, "mnb")) charge = 1;
191 if (strstr(str, "plp")) charge = 2;
192 if (strstr(str, "ppb")) charge = -2;
193 tInstance->AddParticle(str, str, mass, gamma == 0.0 ? 1:0, gamma, charge , "meson", (int) mc);
194 // printf("Added particle %s with PDG PID %d charge %d", str, (int) mc, charge);
195 //AliWarning(Form("Added particle %s with PDG PID %d charge %d", str, (int) mc, charge));
196 // 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));
204 void TTherminator::GenerateEvent()
206 // Run single event generation
207 if ((sRunType == 3) || (sRunType == 4))
212 if (!fEvent) fEvent = new Event(fPartDB, fCalka);
213 if (sRunType == 4) fEvent->Randomize();
215 fEvent->GenerateEvent();
216 fEvent->DecayParticles();
217 // fEvent->WriteEvent(0);
221 Int_t TTherminator::ImportParticles(TClonesArray *particles, Option_t *option)
223 // Import particles from a generated event into an external array
224 if (particles == 0) return 0;
225 TClonesArray &particlesR = *particles;
228 if (!fEvent) return 0;
229 Int_t numpart = fEvent->GetParticleCount();
230 printf("\n TTherminator: Therminator stack contains %d particles.", numpart);
231 for (Int_t iPart=0; iPart<numpart; iPart++) {
232 Particle *tPart = fEvent->GetParticleOfCount(iPart);
234 tFather = tPart->GetFather();
237 TParticle *mother = (TParticle*) (particlesR.UncheckedAt(tFather+1));
238 mother->SetLastDaughter(iPart);
239 if (mother->GetFirstDaughter()==-1)
240 mother->SetFirstDaughter(iPart);
243 // printf("Putting %d %d %lf %d\n", tPart->GetParticleType()->GetPDGCode(), iPart, tPart->px, tPart);
244 new (particlesR[iPart]) TParticle(tPart->GetParticleType()->GetPDGCode(), tPart->HadDecayed(),
246 tPart->px, tPart->py, tPart->pz, tPart->GetEnergy() ,
247 tPart->rx*1.e-13, tPart->ry*1.e-13, tPart->rz*1.e-13, tPart->rt*1.e-13/300000000
249 particlesR[iPart]->SetUniqueID(iPart);
255 TObjArray* TTherminator::ImportParticles(Option_t *option)
257 // Import generated particles into an internal array
260 if (!fEvent) return 0;
261 Int_t numpart = fEvent->GetParticleCount();
262 printf("\n TTherminator: Therminator stack contains %d particles.", numpart);
263 fEvent->InitParticleScan();
264 for (Int_t iPart=0; iPart<numpart; iPart++) {
265 Particle *tPart = fEvent->GetNextParticle();
267 tFather = tPart->GetFather();
270 TParticle *mother = (TParticle*) (fParticles->UncheckedAt(tFather));
271 mother->SetLastDaughter(iPart);
272 if (mother->GetFirstDaughter()==-1)
273 mother->SetFirstDaughter(iPart);
275 TParticle* p = new TParticle(tPart->GetParticleType()->GetPDGCode(), tPart->HadDecayed(),
277 tPart->px, tPart->py, tPart->pz, tPart->GetEnergy() ,
278 tPart->rx*1.e-13, tPart->ry*1.e-13, tPart->rz*1.e-13, tPart->rt*1.e-13/300000000
280 p->SetUniqueID(iPart);
283 nump = fEvent->GetParticleCount();