]>
Commit | Line | Data |
---|---|---|
2e967919 | 1 | /////////////////////////////////////////////////////////////////////////////// |
2 | // // | |
3 | // TTherminator: global interface class to the Therminator model. // | |
4 | // Initialized global variables and runs the event generation // | |
5 | // // | |
6 | /////////////////////////////////////////////////////////////////////////////// | |
7 | /****************************************************************************** | |
8 | * T H E R M I N A T O R * | |
9 | * THERMal heavy-IoN generATOR * | |
10 | * version 1.0 * | |
11 | * * | |
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 * | |
18 | * * | |
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 * | |
22 | * * | |
23 | * Homepage: http://hirg.if.pw.edu.pl/en/therminator/ * | |
24 | * * | |
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 * | |
30 | * available. * | |
31 | * * | |
32 | *****************************************************************************/ | |
33 | #include <TTherminator.h> | |
34 | #include <fstream> | |
35 | #include <TDatabasePDG.h> | |
36 | #include <TParticle.h> | |
37 | #include <TClonesArray.h> | |
38 | ||
39 | ReadPar *sRPInstance; | |
40 | STR sRPFileName; | |
41 | int sNumEvents; | |
42 | int sRunType; | |
43 | int sTables; | |
44 | int sModel; | |
45 | int sIntegrateSample; | |
46 | ||
47 | ClassImp(TTherminator) | |
48 | ||
49 | TTherminator::TTherminator(): | |
50 | fCalka(0), | |
51 | fEvent(0), | |
52 | fPartDB(0) | |
53 | { | |
54 | // Default constructor | |
55 | fPartDB = new ParticleDB(); | |
56 | } | |
57 | TTherminator::TTherminator(const TTherminator & therm) : | |
58 | fCalka(0), | |
59 | fEvent(0), | |
60 | fPartDB(0) | |
61 | { | |
62 | // Copy constructor | |
63 | fPartDB = new ParticleDB(); | |
64 | } | |
65 | TTherminator::~TTherminator() | |
66 | { | |
67 | // Destructor | |
68 | if (sRPInstance) delete sRPInstance; | |
69 | if (fEvent) delete fEvent; | |
70 | if (fCalka) delete fCalka; | |
71 | if (fPartDB) delete fPartDB; | |
72 | } | |
73 | ||
74 | void TTherminator::ReadParameters() | |
75 | { | |
76 | // Read in global model parameters | |
77 | STR tModel; | |
78 | STR tTable; | |
79 | ||
80 | try { | |
81 | sNumEvents = atoi(sRPInstance->getPar("NumberOfEvents").Data()); | |
82 | tModel = sRPInstance->getPar("FreezeOutModel"); | |
83 | if (tModel.Contains("SingleFreezeOut")) | |
84 | sModel = 0; | |
85 | /*MCH begin*/ | |
86 | else if (tModel.Contains("Lhyquid3D")) | |
87 | sModel = 10; | |
88 | else if (tModel.Contains("Lhyquid2D")) | |
89 | sModel = 11; | |
90 | /*MCH end*/ | |
91 | else if (tModel.Contains("BlastWaveVTDelay")) | |
92 | sModel = 6; | |
93 | else if (tModel.Contains("BlastWaveVT")) | |
94 | sModel = 2; | |
95 | else if (tModel.Contains("BlastWaveVLinearFormation")) | |
96 | sModel = 7; | |
97 | else if (tModel.Contains("BlastWaveVLinearDelay")) | |
98 | sModel = 8; | |
99 | else if (tModel.Contains("BlastWaveVLinear")) | |
100 | sModel = 4; | |
101 | else if (tModel.Contains("FiniteHydro")) | |
102 | sModel = 5; | |
103 | else { | |
104 | PRINT_MESSAGE("Unknown model type: " << tModel.Data()); | |
105 | PRINT_MESSAGE("Please provide the proper name"); | |
106 | exit(0); | |
107 | } | |
108 | if (atoi(sRPInstance->getPar("Randomize").Data()) == 1) | |
109 | sRunType = 4; | |
110 | else | |
111 | sRunType = 3; | |
112 | tTable = sRPInstance->getPar("TableType"); | |
113 | if (tTable.Contains("Mathematica")) | |
114 | sTables = 0; | |
115 | else if (tTable.Contains("SHARE")) | |
116 | sTables = 1; | |
117 | else { | |
118 | PRINT_MESSAGE("Unknown table type: " << tTable.Data()); | |
119 | exit(0); | |
120 | } | |
121 | sIntegrateSample = atoi(sRPInstance->getPar("NumberOfIntegrateSamples").Data()); | |
122 | } | |
123 | catch (STR tError) { | |
124 | PRINT_DEBUG_1("RunBFPW::ReadParameters - Caught exception " << tError); | |
125 | PRINT_MESSAGE("Did not find one of the neccessary parameters in the parameters file."); | |
126 | exit(0); | |
127 | } | |
128 | } | |
129 | ||
130 | void TTherminator::Initialize(){ | |
131 | // Initialize global variables | |
132 | // and particle and decay tables | |
133 | Parser *tParser; | |
134 | ||
135 | sRPFileName = "therminator.in"; | |
136 | sRPInstance = new ReadPar(sRPFileName.Data()); | |
137 | ||
138 | ReadParameters(); | |
139 | ||
140 | tParser = new Parser(fPartDB); | |
141 | ||
142 | tParser->ReadShare(); | |
143 | ||
144 | fCalka = new Integrator(sIntegrateSample); | |
145 | fCalka->ReadMultInteg(fPartDB); | |
146 | ||
147 | ||
148 | // Read in particle table from share | |
149 | // and add missing particle type to TDatabasePDG | |
150 | ||
151 | char str[50]; | |
152 | char str1[200]; | |
153 | ||
154 | // ParticleType *tPartBuf; | |
155 | ||
156 | TDatabasePDG *tInstance = TDatabasePDG::Instance(); | |
157 | TParticlePDG *tParticleType; | |
158 | ||
159 | // AliWarning(Form("Reading particle types from particles.data")); | |
160 | ifstream in("particles.data"); | |
161 | ||
162 | int charge; | |
163 | ||
164 | int number=0; | |
165 | if ((in) && (in.is_open())) | |
166 | { | |
167 | //START OF HEAD-LINE | |
168 | in.ignore(200,'\n'); | |
169 | in.ignore(200,'\n'); | |
170 | in.ignore(200,'\n'); | |
171 | //END OF HEAD-LINE | |
172 | ||
173 | while (in>>str) | |
174 | { | |
175 | if (/*(*str == '#')||*/(*str<65)||(*str>122)) | |
176 | { | |
177 | in.getline(str1,200); | |
178 | continue; | |
179 | } | |
180 | double mass, gamma, spin, tI3, tI, q, s, aq, as, c, ac, mc; | |
181 | ||
182 | in>>mass>>gamma>>spin>>tI>>tI3>>q>>s>>aq>>as>>c>>ac>>mc; | |
183 | number++; | |
184 | tParticleType = tInstance->GetParticle((int) mc); | |
185 | if (!tParticleType) { | |
186 | charge = 0; | |
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)); | |
197 | ||
198 | } | |
199 | } | |
200 | in.close(); | |
201 | } | |
202 | } | |
203 | ||
204 | void TTherminator::GenerateEvent() | |
205 | { | |
206 | // Run single event generation | |
207 | if ((sRunType == 3) || (sRunType == 4)) | |
208 | { | |
209 | if (sRunType==4) | |
210 | fCalka->Randomize(); | |
211 | ||
212 | if (!fEvent) fEvent = new Event(fPartDB, fCalka); | |
213 | if (sRunType == 4) fEvent->Randomize(); | |
214 | ||
215 | fEvent->GenerateEvent(); | |
216 | fEvent->DecayParticles(); | |
217 | // fEvent->WriteEvent(0); | |
218 | } | |
219 | } | |
220 | ||
221 | Int_t TTherminator::ImportParticles(TClonesArray *particles, Option_t *option) | |
222 | { | |
223 | // Import particles from a generated event into an external array | |
224 | if (particles == 0) return 0; | |
225 | TClonesArray &particlesR = *particles; | |
226 | particlesR.Clear(); | |
227 | Int_t nump = 0; | |
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); | |
233 | Int_t tFather; | |
234 | tFather = tPart->GetFather(); | |
235 | ||
236 | if (tFather> -1) { | |
237 | TParticle *mother = (TParticle*) (particlesR.UncheckedAt(tFather+1)); | |
238 | mother->SetLastDaughter(iPart); | |
239 | if (mother->GetFirstDaughter()==-1) | |
240 | mother->SetFirstDaughter(iPart); | |
241 | } | |
242 | nump++; | |
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(), | |
245 | tFather, -1, -1, -1, | |
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 | |
248 | ); | |
249 | particlesR[iPart]->SetUniqueID(iPart); | |
250 | } | |
251 | ||
252 | return nump; | |
253 | } | |
254 | ||
255 | TObjArray* TTherminator::ImportParticles(Option_t *option) | |
256 | { | |
257 | // Import generated particles into an internal array | |
258 | Int_t nump = 0; | |
259 | fParticles->Clear(); | |
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(); | |
266 | Int_t tFather; | |
267 | tFather = tPart->GetFather(); | |
268 | ||
269 | if (tFather> -1) { | |
270 | TParticle *mother = (TParticle*) (fParticles->UncheckedAt(tFather)); | |
271 | mother->SetLastDaughter(iPart); | |
272 | if (mother->GetFirstDaughter()==-1) | |
273 | mother->SetFirstDaughter(iPart); | |
274 | } | |
275 | TParticle* p = new TParticle(tPart->GetParticleType()->GetPDGCode(), tPart->HadDecayed(), | |
276 | tFather, -1, -1, -1, | |
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 | |
279 | ); | |
280 | p->SetUniqueID(iPart); | |
281 | fParticles->Add(p); | |
282 | } | |
283 | nump = fEvent->GetParticleCount(); | |
284 | ||
285 | return fParticles; | |
286 | ||
287 | } |