6671656e |
1 | #include <TGenerator.h> |
2 | #include "THBTprocessor.h" |
3 | #include <TParticle.h> |
4 | #include <TMath.h> |
5 | |
6 | #include <iostream.h> |
7 | #ifndef WIN32 |
8 | # define hbtprocessor hbtprocessor_ |
9 | # define ctest ctest_ |
10 | # define type_of_call |
11 | #else |
12 | # define hbtprocessor HBTPROCESSOR |
13 | # define ctest CTEST |
14 | # define type_of_call _stdcall |
15 | #endif |
16 | |
17 | |
18 | ClassImp(THBTprocessor) |
19 | |
20 | extern "C" void type_of_call hbtprocessor(); |
21 | extern "C" void type_of_call ctest(); |
22 | |
23 | //Wrapper class to HBT processor fortran code. |
24 | //For more information see AliGenHBTprocessor class |
25 | //HBT processor is written by Lanny Ray |
26 | //Wrapper class written by Piotr Krzysztof Skowronski (Piotr.Skowronski@cern.ch) |
27 | // |
28 | |
29 | /*************************************************/ |
30 | THBTprocessor::THBTprocessor()// it is better not to intialize it:TGenerator("THBTprocessor","THBTprocessor") |
31 | //because it allocates memmory for TObjArray::fParticles which is not used in our case |
32 | //and we are using TClonesArray in import paerticles |
33 | { |
34 | //constructor |
35 | PARAMETERS.ALICE = 1; //flag that we are working in AliRoot (0==STAR stand alone mode) |
36 | PARAMETERS.pi = TMath::Pi();//3.141592654; |
37 | PARAMETERS.rad = 180.0/TMath::Pi(); |
38 | PARAMETERS.hbc = 0.19732891; |
39 | |
40 | fParticles = 0; |
41 | Initialize(); //Enforce initialization (cleaning all commons) |
42 | } |
43 | |
44 | void THBTprocessor::GenerateEvent() |
45 | { |
46 | //Starts processing |
47 | |
48 | cout<<endl<<"Entering Fortran"<<endl; |
49 | ctest(); |
50 | hbtprocessor(); |
51 | cout<<endl<<"Exited Fortran"<<endl; |
52 | if(PARAMETERS.errorcode) |
53 | { |
54 | Fatal("THBTprocessor::GenerateEvent()","HBT Processor (fortran part) finished with errors\n \ |
55 | See hbt_simulation.out file for more detailed information"); |
56 | } |
57 | else |
58 | cout<<endl<<"GOOD ! HBT Processor finished without errors"<<endl; |
59 | cout<<"PARAMETERS.errorcode= "<<PARAMETERS.errorcode<<endl; |
60 | |
61 | } |
62 | /*****************************************************************************************/ |
63 | void THBTprocessor::Initialize() |
4c90f2a0 |
64 | { |
65 | //IT RESETS ALL THE PREVIOUS SETTINGS |
6671656e |
66 | //Call this method to set default values in PARAMETERS & MESH |
67 | //and zero other common block |
68 | |
4c90f2a0 |
69 | if(gDebug) cout<<"\nHBT PROCESSOR::Initialize() Setting Default valuses in all COMMON BLOCKS"<<endl; |
6671656e |
70 | |
71 | PARAMETERS.ref_control = 2; |
72 | PARAMETERS.switch_1d = 3; |
73 | PARAMETERS.switch_3d = 1; |
74 | PARAMETERS.switch_type = 3; |
75 | PARAMETERS.switch_coherence = 0; |
76 | PARAMETERS.switch_coulomb = 3; |
77 | PARAMETERS.switch_fermi_bose= 1; |
78 | PARAMETERS.trk_accep = 1.0; |
79 | PARAMETERS.print_full = 1; |
80 | PARAMETERS.print_sector_data= 1; |
81 | |
82 | PARAMETERS.n_pid_types = 2; |
83 | PARAMETERS.pid[0] = 8; |
84 | PARAMETERS.pid[1] = 9; |
85 | PARAMETERS.deltap = 0.1; |
86 | PARAMETERS.maxit = 50; |
87 | PARAMETERS.delchi = 1.0; |
88 | PARAMETERS.irand = 76564; |
89 | PARAMETERS.lambda = 0.6; |
90 | PARAMETERS.R_1d = 7.0; |
91 | PARAMETERS.Rside = 6.0; |
92 | |
93 | PARAMETERS.Rout = 7.0; |
94 | PARAMETERS.Rlong = 4.0; |
95 | PARAMETERS.Rperp = 6.0; |
96 | PARAMETERS.Rparallel = 4.0; |
97 | PARAMETERS.R0 = 4.0; |
98 | PARAMETERS.Q0 = 9.0; |
99 | |
100 | PARAMETERS.n_part_1_trk = 0; |
101 | PARAMETERS.n_part_2_trk = 0; |
102 | PARAMETERS.n_part_tot_trk = 0; |
103 | PARAMETERS.n_part_used_1_trk = 0; |
104 | |
105 | PARAMETERS.n_part_used_2_trk = 0; |
106 | PARAMETERS.n_part_1_trk2 = 0; |
107 | PARAMETERS.n_part_2_trk2 = 0; |
108 | PARAMETERS.n_part_tot_trk2 = 0; |
109 | PARAMETERS.n_part_used_1_trk2 = 0; |
110 | PARAMETERS.n_part_used_2_trk2 = 0; |
111 | PARAMETERS.n_part_used_1_ref = 0; |
112 | PARAMETERS.n_part_used_2_ref = 0; |
113 | PARAMETERS.n_part_used_1_inc = 0; |
114 | PARAMETERS.n_part_used_2_inc = 0; |
115 | |
116 | PARAMETERS.num_pairs_like = 0; |
117 | PARAMETERS.num_pairs_unlike = 0; |
118 | PARAMETERS.num_pairs_like_ref = 0; |
119 | PARAMETERS.num_pairs_like_inc = 0; |
120 | PARAMETERS.num_pairs_unlike_inc = 0; |
121 | PARAMETERS.event_line_counter = 0; |
122 | PARAMETERS.file10_line_counter = 0; |
123 | |
124 | PARAMETERS.chisq_wt_like_1d = 1.0; |
125 | PARAMETERS.chisq_wt_unlike_1d = 1.0; |
126 | PARAMETERS.chisq_wt_like_3d_fine = 1.0; |
127 | |
128 | PARAMETERS.chisq_wt_unlike_3d_fine = 1.0; |
129 | PARAMETERS.chisq_wt_like_3d_coarse = 1.0; |
130 | PARAMETERS.chisq_wt_unlike_3d_coarse= 1.0; |
131 | PARAMETERS.chisq_wt_hist1_1 = 1.0; |
132 | PARAMETERS.chisq_wt_hist1_2 = 1.0; // /////internal comment 25 fields |
133 | |
134 | /*********************************************/ |
135 | |
136 | |
137 | MESH.n_pt_bins = 50; //OK |
138 | MESH.pt_min = 0.1; //Pt in GeV/c //OK |
139 | MESH.pt_max = 0.98; //Pt in GeV/c |
140 | MESH.n_phi_bins = 50; //OK |
141 | MESH.phi_min = 0.0; //OK |
142 | MESH.phi_max = 360.0; //OK |
143 | MESH.n_eta_bins = 50; //OK |
144 | MESH.eta_min =-1.5; //OK |
145 | MESH.eta_max = 1.5; //OK |
146 | MESH.n_1d_fine = 10; //OK |
147 | MESH.binsize_1d_fine = 0.01; //ok |
148 | MESH.n_1d_coarse = 2; //O |
149 | MESH.binsize_1d_coarse= 0.05; //ok |
150 | MESH.n_3d_fine = 8; //OK |
151 | MESH.binsize_3d_fine = 0.01; //ok |
152 | MESH.n_3d_coarse = 2; //OK |
153 | MESH.binsize_3d_coarse= 0.08; //ok |
154 | MESH.n_3d_fine_project= 3; //OK |
155 | MESH.n_px_bins = 20; //OK |
156 | MESH.px_min =-1.0; //OK |
157 | MESH.px_max = 1.0; //OK |
158 | MESH.n_py_bins = 20; //OK |
159 | MESH.py_min =-1.0; //OK |
160 | MESH.py_max = 1.0; //OK |
161 | MESH.n_pz_bins = 70; //OK |
162 | MESH.pz_min =-3.6; //OK |
163 | MESH.pz_max = 3.6; //OK |
164 | |
165 | /*********************************************/ |
450f427d |
166 | |
167 | Int_t i; //loop variable |
6671656e |
168 | |
450f427d |
169 | for (i =0; i<TRK_MAXLEN;i++) |
6671656e |
170 | { |
171 | TRACK1.trk_id[i] = 0; |
172 | TRACK1.trk_px_sec[i] = 0; |
173 | TRACK1.trk_py_sec[i] = 0; |
174 | TRACK1.trk_pz_sec[i] = 0; |
175 | TRACK1.trk_sector[i] = 0; |
176 | TRACK1.trk_flag[i] = 0; |
177 | TRACK1.trk_out_flag[i] = 0; |
178 | TRACK1.trk_merge_flag[i] = 0; |
179 | TRACK1.trk_ge_pid[i] = 0; |
180 | TRACK1.trk_start_vertex[i] = 0; |
181 | TRACK1.trk_stop_vertex[i] = 0; |
182 | TRACK1.trk_event_line[i] = 0; |
183 | TRACK1.trk_px[i] = 0; |
184 | TRACK1.trk_py[i] = 0; |
185 | TRACK1.trk_pz[i] = 0; |
186 | TRACK1.trk_E[i] = 0; |
187 | TRACK1.trk_pt[i] = 0; |
188 | TRACK1.trk_phi[i] = 0; |
189 | TRACK1.trk_eta[i] = 0; |
190 | } |
191 | |
192 | /*********************************************/ |
193 | |
450f427d |
194 | for (i =0; i<TRK2_MAXLEN;i++) |
6671656e |
195 | { |
196 | TRACK2.trk2_id[i] = 0; |
197 | TRACK2.trk2_px_sec[i] = 0; |
198 | TRACK2.trk2_py_sec[i] = 0; |
199 | TRACK2.trk2_pz_sec[i] = 0; |
200 | TRACK2.trk2_sector[i] = 0; |
201 | TRACK2.trk2_flag[i] = 0; |
202 | TRACK2.trk2_out_flag[i] = 0; |
203 | TRACK2.trk2_merge_flag[i] = 0; |
204 | TRACK2.trk2_ge_pid[i] = 0; |
205 | TRACK2.trk2_start_vertex[i] = 0; |
206 | TRACK2.trk2_stop_vertex[i] = 0; |
207 | TRACK2.trk2_event_line[i] = 0; |
208 | TRACK2.trk2_px[i] = 0; |
209 | TRACK2.trk2_py[i] = 0; |
210 | TRACK2.trk2_pz[i] = 0; |
211 | TRACK2.trk2_E[i] = 0; |
212 | TRACK2.trk2_pt[i] = 0; |
213 | TRACK2.trk2_phi[i] = 0; |
214 | TRACK2.trk2_eta[i] = 0; |
215 | } |
216 | |
217 | /*********************************************/ |
218 | |
450f427d |
219 | for (i =0; i<SEC_MAXLEN;i++) |
6671656e |
220 | { |
221 | SEC_TRK_MAP.stm_sec_id [i] = 0; |
222 | SEC_TRK_MAP.stm_n_trk_sec[i] = 0; |
223 | SEC_TRK_MAP.stm_flag[i] = 0; |
224 | |
225 | for (Int_t j=0; j<MAX_TRK_SEC;j++) |
226 | SEC_TRK_MAP.stm_track_id[i][j] = 0; |
227 | } |
228 | |
229 | /*********************************************/ |
230 | |
450f427d |
231 | for (i =0; i<SEC_MAXLEN2;i++) |
6671656e |
232 | { |
233 | SEC_TRK_MAP2.stm_sec_id2[i] = 0; |
234 | SEC_TRK_MAP2.stm_n_trk_sec2[i] = 0; |
235 | SEC_TRK_MAP2.stm_flag2[i] = 0; |
236 | |
237 | for (Int_t j=0; j<MAX_TRK_SEC;j++) |
238 | SEC_TRK_MAP2.stm_track_id2[i][j] = 0; |
239 | } |
240 | |
241 | /*********************************************/ |
242 | |
450f427d |
243 | for (i =0; i<PART_MAXLEN;i++) |
6671656e |
244 | { |
245 | PARTICLE.part_id[i] = 0; |
246 | PARTICLE.part_charge[i] = 0; |
247 | PARTICLE.part_mass[i] = 0; |
248 | PARTICLE.part_lifetime[i] = 0; |
249 | } |
250 | |
251 | |
252 | /*********************************************/ |
450f427d |
253 | for (i =0; i<MAX_C2_1D;i++) |
6671656e |
254 | { |
255 | CORRELATIONS.c2mod_like_1d[i] = 0; |
256 | CORRELATIONS.c2mod_unlike_1d[i] = 0; |
257 | CORRELATIONS.c2fit_like_1d[i] = 0; |
258 | CORRELATIONS.c2fit_unlike_1d[i] = 0; |
259 | CORRELATIONS.c2err_like_1d[i] = 0; |
260 | CORRELATIONS.c2err_unlike_1d[i] = 0; |
261 | } |
262 | /*********************************************/ |
450f427d |
263 | for (i =0; i<MAX_C2_3D;i++) |
6671656e |
264 | for (Int_t j =0; j<MAX_C2_3D;j++) |
265 | for (Int_t k =0; k<MAX_C2_3D;k++) |
266 | { |
267 | CORRELATIONS.c2mod_like_3d_fine[i][j][k] = 0.0; |
268 | CORRELATIONS.c2mod_unlike_3d_fine[i][j][k] = 0.0; |
269 | CORRELATIONS.c2mod_like_3d_coarse[i][j][k] = 0.0; |
270 | CORRELATIONS.c2mod_unlike_3d_coarse[i][j][k] = 0.0; |
271 | CORRELATIONS.c2fit_like_3d_fine[i][j][k] = 0.0; |
272 | CORRELATIONS.c2fit_unlike_3d_fine[i][j][k] = 0.0; |
273 | CORRELATIONS.c2fit_like_3d_coarse[i][j][k] = 0.0; |
274 | CORRELATIONS.c2fit_unlike_3d_coarse[i][j][k] = 0.0; |
275 | CORRELATIONS.c2err_like_3d_fine[i][j][k] = 0.0; |
276 | CORRELATIONS.c2err_unlike_3d_fine[i][j][k] = 0.0; |
277 | CORRELATIONS.c2err_like_3d_coarse[i][j][k] = 0.0; |
278 | CORRELATIONS.c2err_unlike_3d_coarse[i][j][k] = 0.0; |
279 | } |
280 | /*********************************************/ |
281 | |
282 | EVENT_SUMMARY.niter_mean = 0.0; |
283 | EVENT_SUMMARY.niter_rms = 0.0; |
284 | |
285 | EVENT_SUMMARY.npart1_mean = 0.0; |
286 | EVENT_SUMMARY.npart1_rms = 0.0; |
287 | |
288 | EVENT_SUMMARY.npart2_mean = 0.0; |
289 | EVENT_SUMMARY.npart2_rms = 0.0; |
290 | |
291 | EVENT_SUMMARY.npart_tot_mean = 0.0; |
292 | EVENT_SUMMARY.npart_tot_rms = 0.0; |
293 | |
294 | EVENT_SUMMARY.nsec_flag_mean = 0.0; |
295 | EVENT_SUMMARY.nsec_flag_rms = 0.0; |
296 | |
297 | EVENT_SUMMARY.frac_trks_out_mean = 0.0; |
298 | EVENT_SUMMARY.frac_trks_out_rms = 0.0; |
299 | |
300 | EVENT_SUMMARY.frac_trks_flag_mean = 0.0; |
301 | EVENT_SUMMARY.frac_trks_flag_rms = 0.0; |
302 | |
303 | EVENT_SUMMARY.chi_l1d_mean = 0.0; |
304 | EVENT_SUMMARY.chi_l1d_rms = 0.0; |
305 | |
306 | EVENT_SUMMARY.chi_u1d_mean = 0.0; |
307 | EVENT_SUMMARY.chi_u1d_rms = 0.0; |
308 | |
450f427d |
309 | for (i =0; i<MAX_EVENTS;i++) |
6671656e |
310 | { |
311 | EVENT_SUMMARY.n_part_used_1_store[i] = 0.0; |
312 | EVENT_SUMMARY.n_part_used_2_store[i] = 0.0; |
313 | EVENT_SUMMARY.n_part_tot_store[i] = 0.0; |
314 | EVENT_SUMMARY.num_sec_flagged_store[i] = 0.0; |
315 | EVENT_SUMMARY.frac_trks_out[i] = 0.0; |
316 | EVENT_SUMMARY.frac_trks_flag[i] = 0.0; |
317 | EVENT_SUMMARY.chisq_like_1d_store[i] = 0.0; |
318 | EVENT_SUMMARY.num_iter[i] = 0.0; |
319 | EVENT_SUMMARY.chisq_unlike_1d_store[i] = 0.0; |
320 | } |
321 | /*********************************************/ |
450f427d |
322 | for (i =0; i<MAX_C2_COUL;i++) |
6671656e |
323 | { |
324 | COULMB.c2_coul_like[i] = 0.0; |
325 | COULMB.c2_coul_unlike[i] = 0.0; |
326 | COULMB.q_coul[i] = 0.0; |
327 | |
328 | } |
329 | /*********************************************/ |
450f427d |
330 | for (i =0; i<MAX_H_1D;i++) |
6671656e |
331 | { |
332 | HISTOGRAMS.hist1_pt_1[i] = 0; |
333 | HISTOGRAMS.hist1_phi_1[i] = 0; |
334 | HISTOGRAMS.hist1_eta_1[i] = 0; |
335 | HISTOGRAMS.hist1_pt_2[i] = 0; |
336 | HISTOGRAMS.hist1_phi_2[i] = 0; |
337 | HISTOGRAMS.hist1_eta_2[i] = 0; |
338 | HISTOGRAMS.htmp1_pt_1[i] = 0; |
339 | HISTOGRAMS.htmp1_phi_1[i] = 0; |
340 | HISTOGRAMS.htmp1_eta_1[i] = 0; |
341 | HISTOGRAMS.htmp1_pt_2[i] = 0; |
342 | HISTOGRAMS.htmp1_phi_2[i] = 0; |
343 | HISTOGRAMS.htmp1_eta_2[i] = 0; |
344 | HISTOGRAMS.href1_pt_1[i] = 0; |
345 | HISTOGRAMS.href1_phi_1[i] = 0; |
346 | HISTOGRAMS.href1_eta_1[i] = 0; |
347 | HISTOGRAMS.href1_pt_2[i] = 0; |
348 | HISTOGRAMS.href1_phi_2[i] = 0; |
349 | HISTOGRAMS.href1_eta_2[i] = 0; |
350 | HISTOGRAMS.hinc1_pt_1[i] = 0; |
351 | HISTOGRAMS.hinc1_phi_1[i] = 0; |
352 | HISTOGRAMS.hinc1_eta_1[i] = 0; |
353 | HISTOGRAMS.hinc1_pt_2[i] = 0; |
354 | HISTOGRAMS.hinc1_phi_2[i] = 0; |
355 | HISTOGRAMS.hinc1_eta_2[i] = 0; |
356 | HISTOGRAMS.hist_like_1d[i] = 0; |
357 | HISTOGRAMS.hist_unlike_1d[i] = 0; |
358 | HISTOGRAMS.htmp_like_1d[i] = 0; |
359 | HISTOGRAMS.htmp_unlike_1d[i] = 0; |
360 | HISTOGRAMS.href_like_1d[i] = 0; |
361 | HISTOGRAMS.href_unlike_1d[i] = 0; |
362 | HISTOGRAMS.hinc_like_1d[i] = 0; |
363 | HISTOGRAMS.hinc_unlike_1d[i] = 0; |
364 | } |
365 | |
450f427d |
366 | for (i =0; i<MAX_H_3D;i++) |
6671656e |
367 | for (Int_t j =0; j<MAX_H_3D;j++) |
368 | for (Int_t k =0; k<MAX_H_3D;k++) |
369 | { |
370 | HISTOGRAMS.hist_like_3d_fine[i][j][k] = 0; |
371 | HISTOGRAMS.hist_unlike_3d_fine[i][j][k] = 0; |
372 | HISTOGRAMS.hist_like_3d_coarse[i][j][k] = 0; |
373 | HISTOGRAMS.hist_unlike_3d_coarse[i][j][k] = 0; |
374 | HISTOGRAMS.htmp_like_3d_fine[i][j][k] = 0; |
375 | HISTOGRAMS.htmp_unlike_3d_fine[i][j][k] = 0; |
376 | HISTOGRAMS.htmp_like_3d_coarse[i][j][k] = 0; |
377 | HISTOGRAMS.htmp_unlike_3d_coarse[i][j][k] = 0; |
378 | HISTOGRAMS.href_like_3d_fine[i][j][k] = 0; |
379 | HISTOGRAMS.href_unlike_3d_fine[i][j][k] = 0; |
380 | HISTOGRAMS.href_like_3d_coarse[i][j][k] = 0; |
381 | HISTOGRAMS.href_unlike_3d_coarse[i][j][k] = 0; |
382 | HISTOGRAMS.hinc_like_3d_fine[i][j][k] = 0; |
383 | HISTOGRAMS.hinc_unlike_3d_fine[i][j][k] = 0; |
384 | HISTOGRAMS.hinc_like_3d_coarse[i][j][k] = 0; |
385 | HISTOGRAMS.hinc_unlike_3d_coarse[i][j][k] = 0; |
386 | } |
387 | /*********************************************/ |
388 | |
389 | |
390 | /*********************************************/ |
391 | |
392 | // cout<<" FINISHED"<<endl; |
393 | |
394 | } |
395 | |
396 | |
397 | /*****************************************************************************************/ |
398 | |
399 | |
400 | Int_t THBTprocessor::ImportParticles(TClonesArray *particles, Option_t *option) |
401 | { |
402 | //Copy particle data into TClonesArray |
403 | if (particles == 0) return 0; |
404 | TClonesArray &Particles = *particles; |
405 | Particles.Clear(); |
406 | |
407 | Int_t nrpart = 0; |
408 | for (Int_t i=0; i < TRK_MAXLEN; i++) |
409 | { |
410 | |
411 | |
412 | if (TRACK1.trk_E[i] == 0.) continue; |
413 | |
414 | Float_t px = TRACK1.trk_px[i]; |
415 | Float_t py = TRACK1.trk_py[i]; |
416 | Float_t pz = TRACK1.trk_pz[i]; |
417 | // Float_t pE = TRACK.trk_E[i]; |
418 | Float_t mass = PARTICLE.part_mass[TRACK1.trk_ge_pid[i]]; |
419 | |
420 | new(Particles[nrpart++]) TParticle(0,0,0,0,0,0,px,py,pz, |
421 | TMath::Sqrt(mass*mass+px*px+py*py+pz*pz), |
422 | 0,0,0,0); |
423 | } |
424 | return nrpart; |
425 | } |
426 | |
427 | /*****************************************************************************************/ |
428 | |
429 | void THBTprocessor::PrintEvent() |
430 | { |
431 | //Prints all particles (common block data) |
432 | cout<<"Print Event"<<endl; |
433 | for (Int_t i=0; i<TRK_MAXLEN;i++) |
434 | { |
435 | if(TRACK1.trk_E[i]==0.) continue; |
436 | |
437 | cout<<"trk_id: "<<TRACK1.trk_id[i]<<" trk_px :"<<TRACK1.trk_px[i]<<" trk_py :"<<TRACK1.trk_py[i]<<" trk_pz :"<<TRACK1.trk_pz[i]<<endl; |
438 | cout<<" trk_E: "<<TRACK1.trk_E[i]<<" trk_pt: "<<TRACK1.trk_pt[i]<<" trk_phi: "<<TRACK1.trk_phi[i]<<" trk_eta: "<<TRACK1.trk_eta[i]<<endl; |
439 | } |
440 | } |
441 | |
442 | |
443 | /*****************************************************************************************/ |