]>
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()); | |
ae89e57a | 122 | fInputDir = sRPInstance->getPar("InputDirSHARE"); |
2e967919 | 123 | } |
124 | catch (STR tError) { | |
125 | PRINT_DEBUG_1("RunBFPW::ReadParameters - Caught exception " << tError); | |
126 | PRINT_MESSAGE("Did not find one of the neccessary parameters in the parameters file."); | |
127 | exit(0); | |
128 | } | |
129 | } | |
130 | ||
131 | void TTherminator::Initialize(){ | |
132 | // Initialize global variables | |
133 | // and particle and decay tables | |
134 | Parser *tParser; | |
135 | ||
136 | sRPFileName = "therminator.in"; | |
137 | sRPInstance = new ReadPar(sRPFileName.Data()); | |
138 | ||
139 | ReadParameters(); | |
140 | ||
141 | tParser = new Parser(fPartDB); | |
142 | ||
143 | tParser->ReadShare(); | |
144 | ||
145 | fCalka = new Integrator(sIntegrateSample); | |
146 | fCalka->ReadMultInteg(fPartDB); | |
147 | ||
148 | ||
149 | // Read in particle table from share | |
150 | // and add missing particle type to TDatabasePDG | |
151 | ||
152 | char str[50]; | |
153 | char str1[200]; | |
154 | ||
155 | // ParticleType *tPartBuf; | |
156 | ||
157 | TDatabasePDG *tInstance = TDatabasePDG::Instance(); | |
158 | TParticlePDG *tParticleType; | |
159 | ||
160 | // AliWarning(Form("Reading particle types from particles.data")); | |
ae89e57a | 161 | |
162 | ifstream in((fInputDir+"/"+"particles.data").Data()); | |
163 | // ifstream in("particles.data"); | |
2e967919 | 164 | |
165 | int charge; | |
166 | ||
167 | int number=0; | |
168 | if ((in) && (in.is_open())) | |
169 | { | |
170 | //START OF HEAD-LINE | |
171 | in.ignore(200,'\n'); | |
172 | in.ignore(200,'\n'); | |
173 | in.ignore(200,'\n'); | |
174 | //END OF HEAD-LINE | |
175 | ||
176 | while (in>>str) | |
177 | { | |
178 | if (/*(*str == '#')||*/(*str<65)||(*str>122)) | |
179 | { | |
180 | in.getline(str1,200); | |
181 | continue; | |
182 | } | |
183 | double mass, gamma, spin, tI3, tI, q, s, aq, as, c, ac, mc; | |
184 | ||
185 | in>>mass>>gamma>>spin>>tI>>tI3>>q>>s>>aq>>as>>c>>ac>>mc; | |
186 | number++; | |
187 | tParticleType = tInstance->GetParticle((int) mc); | |
188 | if (!tParticleType) { | |
189 | charge = 0; | |
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)); | |
200 | ||
201 | } | |
202 | } | |
203 | in.close(); | |
204 | } | |
205 | } | |
206 | ||
207 | void TTherminator::GenerateEvent() | |
208 | { | |
209 | // Run single event generation | |
210 | if ((sRunType == 3) || (sRunType == 4)) | |
211 | { | |
212 | if (sRunType==4) | |
213 | fCalka->Randomize(); | |
214 | ||
215 | if (!fEvent) fEvent = new Event(fPartDB, fCalka); | |
216 | if (sRunType == 4) fEvent->Randomize(); | |
217 | ||
218 | fEvent->GenerateEvent(); | |
219 | fEvent->DecayParticles(); | |
220 | // fEvent->WriteEvent(0); | |
221 | } | |
222 | } | |
223 | ||
224 | Int_t TTherminator::ImportParticles(TClonesArray *particles, Option_t *option) | |
225 | { | |
226 | // Import particles from a generated event into an external array | |
227 | if (particles == 0) return 0; | |
228 | TClonesArray &particlesR = *particles; | |
229 | particlesR.Clear(); | |
230 | Int_t nump = 0; | |
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); | |
236 | Int_t tFather; | |
237 | tFather = tPart->GetFather(); | |
238 | ||
239 | if (tFather> -1) { | |
240 | TParticle *mother = (TParticle*) (particlesR.UncheckedAt(tFather+1)); | |
241 | mother->SetLastDaughter(iPart); | |
242 | if (mother->GetFirstDaughter()==-1) | |
243 | mother->SetFirstDaughter(iPart); | |
244 | } | |
245 | nump++; | |
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(), | |
248 | tFather, -1, -1, -1, | |
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 | |
251 | ); | |
252 | particlesR[iPart]->SetUniqueID(iPart); | |
253 | } | |
254 | ||
255 | return nump; | |
256 | } | |
257 | ||
258 | TObjArray* TTherminator::ImportParticles(Option_t *option) | |
259 | { | |
260 | // Import generated particles into an internal array | |
261 | Int_t nump = 0; | |
262 | fParticles->Clear(); | |
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(); | |
269 | Int_t tFather; | |
270 | tFather = tPart->GetFather(); | |
271 | ||
272 | if (tFather> -1) { | |
273 | TParticle *mother = (TParticle*) (fParticles->UncheckedAt(tFather)); | |
274 | mother->SetLastDaughter(iPart); | |
275 | if (mother->GetFirstDaughter()==-1) | |
276 | mother->SetFirstDaughter(iPart); | |
277 | } | |
278 | TParticle* p = new TParticle(tPart->GetParticleType()->GetPDGCode(), tPart->HadDecayed(), | |
279 | tFather, -1, -1, -1, | |
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 | |
282 | ); | |
283 | p->SetUniqueID(iPart); | |
284 | fParticles->Add(p); | |
285 | } | |
286 | nump = fEvent->GetParticleCount(); | |
287 | ||
288 | return fParticles; | |
289 | ||
290 | } |