git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
fstream.h included instead of iostream.h
[u/mrichter/AliRoot.git] / THbtp / THBTprocessor.cxx
CommitLineData
6671656e 1
2#include <TGenerator.h>
3#include "THBTprocessor.h"
4#include <TParticle.h>
5#include <TMath.h>
6
7#include <iostream.h>
8#ifndef WIN32
9# define hbtprocessor hbtprocessor_
10# define ctest ctest_
11# define type_of_call
12#else
13# define hbtprocessor HBTPROCESSOR
14# define ctest CTEST
15# define type_of_call _stdcall
16#endif
17
18
19ClassImp(THBTprocessor)
20
21extern "C" void type_of_call hbtprocessor();
22extern "C" void type_of_call ctest();
23
24//Wrapper class to HBT processor fortran code.
25//For more information see AliGenHBTprocessor class
26//HBT processor is written by Lanny Ray
27//Wrapper class written by Piotr Krzysztof Skowronski (Piotr.Skowronski@cern.ch)
28//
29
30/*************************************************/
31THBTprocessor::THBTprocessor()// it is better not to intialize it:TGenerator("THBTprocessor","THBTprocessor")
32 //because it allocates memmory for TObjArray::fParticles which is not used in our case
33 //and we are using TClonesArray in import paerticles
34 {
35 //constructor
36 PARAMETERS.ALICE = 1; //flag that we are working in AliRoot (0==STAR stand alone mode)
37 PARAMETERS.pi = TMath::Pi();//3.141592654;
38 PARAMETERS.rad = 180.0/TMath::Pi();
39 PARAMETERS.hbc = 0.19732891;
40
41 fParticles = 0;
42 Initialize(); //Enforce initialization (cleaning all commons)
43 }
44
45void THBTprocessor::GenerateEvent()
46{
47//Starts processing
48
49 cout<<endl<<"Entering Fortran"<<endl;
50 ctest();
51 hbtprocessor();
52 cout<<endl<<"Exited Fortran"<<endl;
53 if(PARAMETERS.errorcode)
54 {
55 Fatal("THBTprocessor::GenerateEvent()","HBT Processor (fortran part) finished with errors\n \
56 See hbt_simulation.out file for more detailed information");
57 }
58 else
59 cout<<endl<<"GOOD ! HBT Processor finished without errors"<<endl;
60 cout<<"PARAMETERS.errorcode= "<<PARAMETERS.errorcode<<endl;
61
62}
63/*****************************************************************************************/
64void THBTprocessor::Initialize()
65{
66 //Call this method to set default values in PARAMETERS & MESH
67 //and zero other common block
68
69 cout<<"\nHBT PROCESSOR::Initialize() Setting Default valuses in all COMMON BLOCKS"<<endl;
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/*********************************************/
166
167 for (Int_t i =0; i<TRK_MAXLEN;i++)
168 {
169 TRACK1.trk_id[i] = 0;
170 TRACK1.trk_px_sec[i] = 0;
171 TRACK1.trk_py_sec[i] = 0;
172 TRACK1.trk_pz_sec[i] = 0;
173 TRACK1.trk_sector[i] = 0;
174 TRACK1.trk_flag[i] = 0;
175 TRACK1.trk_out_flag[i] = 0;
176 TRACK1.trk_merge_flag[i] = 0;
177 TRACK1.trk_ge_pid[i] = 0;
178 TRACK1.trk_start_vertex[i] = 0;
179 TRACK1.trk_stop_vertex[i] = 0;
180 TRACK1.trk_event_line[i] = 0;
181 TRACK1.trk_px[i] = 0;
182 TRACK1.trk_py[i] = 0;
183 TRACK1.trk_pz[i] = 0;
184 TRACK1.trk_E[i] = 0;
185 TRACK1.trk_pt[i] = 0;
186 TRACK1.trk_phi[i] = 0;
187 TRACK1.trk_eta[i] = 0;
188 }
189
190/*********************************************/
191
192 for (Int_t i =0; i<TRK2_MAXLEN;i++)
193 {
194 TRACK2.trk2_id[i] = 0;
195 TRACK2.trk2_px_sec[i] = 0;
196 TRACK2.trk2_py_sec[i] = 0;
197 TRACK2.trk2_pz_sec[i] = 0;
198 TRACK2.trk2_sector[i] = 0;
199 TRACK2.trk2_flag[i] = 0;
200 TRACK2.trk2_out_flag[i] = 0;
201 TRACK2.trk2_merge_flag[i] = 0;
202 TRACK2.trk2_ge_pid[i] = 0;
203 TRACK2.trk2_start_vertex[i] = 0;
204 TRACK2.trk2_stop_vertex[i] = 0;
205 TRACK2.trk2_event_line[i] = 0;
206 TRACK2.trk2_px[i] = 0;
207 TRACK2.trk2_py[i] = 0;
208 TRACK2.trk2_pz[i] = 0;
209 TRACK2.trk2_E[i] = 0;
210 TRACK2.trk2_pt[i] = 0;
211 TRACK2.trk2_phi[i] = 0;
212 TRACK2.trk2_eta[i] = 0;
213 }
214
215/*********************************************/
216
217 for (Int_t i =0; i<SEC_MAXLEN;i++)
218 {
219 SEC_TRK_MAP.stm_sec_id [i] = 0;
220 SEC_TRK_MAP.stm_n_trk_sec[i] = 0;
221 SEC_TRK_MAP.stm_flag[i] = 0;
222
223 for (Int_t j=0; j<MAX_TRK_SEC;j++)
224 SEC_TRK_MAP.stm_track_id[i][j] = 0;
225 }
226
227/*********************************************/
228
229 for (Int_t i =0; i<SEC_MAXLEN2;i++)
230 {
231 SEC_TRK_MAP2.stm_sec_id2[i] = 0;
232 SEC_TRK_MAP2.stm_n_trk_sec2[i] = 0;
233 SEC_TRK_MAP2.stm_flag2[i] = 0;
234
235 for (Int_t j=0; j<MAX_TRK_SEC;j++)
236 SEC_TRK_MAP2.stm_track_id2[i][j] = 0;
237 }
238
239/*********************************************/
240
241 for (Int_t i =0; i<PART_MAXLEN;i++)
242 {
243 PARTICLE.part_id[i] = 0;
244 PARTICLE.part_charge[i] = 0;
245 PARTICLE.part_mass[i] = 0;
246 PARTICLE.part_lifetime[i] = 0;
247 }
248
249
250/*********************************************/
251 for (Int_t i =0; i<MAX_C2_1D;i++)
252 {
253 CORRELATIONS.c2mod_like_1d[i] = 0;
254 CORRELATIONS.c2mod_unlike_1d[i] = 0;
255 CORRELATIONS.c2fit_like_1d[i] = 0;
256 CORRELATIONS.c2fit_unlike_1d[i] = 0;
257 CORRELATIONS.c2err_like_1d[i] = 0;
258 CORRELATIONS.c2err_unlike_1d[i] = 0;
259 }
260/*********************************************/
261 for (Int_t i =0; i<MAX_C2_3D;i++)
262 for (Int_t j =0; j<MAX_C2_3D;j++)
263 for (Int_t k =0; k<MAX_C2_3D;k++)
264 {
265 CORRELATIONS.c2mod_like_3d_fine[i][j][k] = 0.0;
266 CORRELATIONS.c2mod_unlike_3d_fine[i][j][k] = 0.0;
267 CORRELATIONS.c2mod_like_3d_coarse[i][j][k] = 0.0;
268 CORRELATIONS.c2mod_unlike_3d_coarse[i][j][k] = 0.0;
269 CORRELATIONS.c2fit_like_3d_fine[i][j][k] = 0.0;
270 CORRELATIONS.c2fit_unlike_3d_fine[i][j][k] = 0.0;
271 CORRELATIONS.c2fit_like_3d_coarse[i][j][k] = 0.0;
272 CORRELATIONS.c2fit_unlike_3d_coarse[i][j][k] = 0.0;
273 CORRELATIONS.c2err_like_3d_fine[i][j][k] = 0.0;
274 CORRELATIONS.c2err_unlike_3d_fine[i][j][k] = 0.0;
275 CORRELATIONS.c2err_like_3d_coarse[i][j][k] = 0.0;
276 CORRELATIONS.c2err_unlike_3d_coarse[i][j][k] = 0.0;
277 }
278/*********************************************/
279
280 EVENT_SUMMARY.niter_mean = 0.0;
281 EVENT_SUMMARY.niter_rms = 0.0;
282
283 EVENT_SUMMARY.npart1_mean = 0.0;
284 EVENT_SUMMARY.npart1_rms = 0.0;
285
286 EVENT_SUMMARY.npart2_mean = 0.0;
287 EVENT_SUMMARY.npart2_rms = 0.0;
288
289 EVENT_SUMMARY.npart_tot_mean = 0.0;
290 EVENT_SUMMARY.npart_tot_rms = 0.0;
291
292 EVENT_SUMMARY.nsec_flag_mean = 0.0;
293 EVENT_SUMMARY.nsec_flag_rms = 0.0;
294
295 EVENT_SUMMARY.frac_trks_out_mean = 0.0;
296 EVENT_SUMMARY.frac_trks_out_rms = 0.0;
297
298 EVENT_SUMMARY.frac_trks_flag_mean = 0.0;
299 EVENT_SUMMARY.frac_trks_flag_rms = 0.0;
300
301 EVENT_SUMMARY.chi_l1d_mean = 0.0;
302 EVENT_SUMMARY.chi_l1d_rms = 0.0;
303
304 EVENT_SUMMARY.chi_u1d_mean = 0.0;
305 EVENT_SUMMARY.chi_u1d_rms = 0.0;
306
307 for (Int_t i =0; i<MAX_EVENTS;i++)
308 {
309 EVENT_SUMMARY.n_part_used_1_store[i] = 0.0;
310 EVENT_SUMMARY.n_part_used_2_store[i] = 0.0;
311 EVENT_SUMMARY.n_part_tot_store[i] = 0.0;
312 EVENT_SUMMARY.num_sec_flagged_store[i] = 0.0;
313 EVENT_SUMMARY.frac_trks_out[i] = 0.0;
314 EVENT_SUMMARY.frac_trks_flag[i] = 0.0;
315 EVENT_SUMMARY.chisq_like_1d_store[i] = 0.0;
316 EVENT_SUMMARY.num_iter[i] = 0.0;
317 EVENT_SUMMARY.chisq_unlike_1d_store[i] = 0.0;
318 }
319/*********************************************/
320 for (Int_t i =0; i<MAX_C2_COUL;i++)
321 {
322 COULMB.c2_coul_like[i] = 0.0;
323 COULMB.c2_coul_unlike[i] = 0.0;
324 COULMB.q_coul[i] = 0.0;
325
326 }
327/*********************************************/
328 for (Int_t i =0; i<MAX_H_1D;i++)
329 {
330 HISTOGRAMS.hist1_pt_1[i] = 0;
331 HISTOGRAMS.hist1_phi_1[i] = 0;
332 HISTOGRAMS.hist1_eta_1[i] = 0;
333 HISTOGRAMS.hist1_pt_2[i] = 0;
334 HISTOGRAMS.hist1_phi_2[i] = 0;
335 HISTOGRAMS.hist1_eta_2[i] = 0;
336 HISTOGRAMS.htmp1_pt_1[i] = 0;
337 HISTOGRAMS.htmp1_phi_1[i] = 0;
338 HISTOGRAMS.htmp1_eta_1[i] = 0;
339 HISTOGRAMS.htmp1_pt_2[i] = 0;
340 HISTOGRAMS.htmp1_phi_2[i] = 0;
341 HISTOGRAMS.htmp1_eta_2[i] = 0;
342 HISTOGRAMS.href1_pt_1[i] = 0;
343 HISTOGRAMS.href1_phi_1[i] = 0;
344 HISTOGRAMS.href1_eta_1[i] = 0;
345 HISTOGRAMS.href1_pt_2[i] = 0;
346 HISTOGRAMS.href1_phi_2[i] = 0;
347 HISTOGRAMS.href1_eta_2[i] = 0;
348 HISTOGRAMS.hinc1_pt_1[i] = 0;
349 HISTOGRAMS.hinc1_phi_1[i] = 0;
350 HISTOGRAMS.hinc1_eta_1[i] = 0;
351 HISTOGRAMS.hinc1_pt_2[i] = 0;
352 HISTOGRAMS.hinc1_phi_2[i] = 0;
353 HISTOGRAMS.hinc1_eta_2[i] = 0;
354 HISTOGRAMS.hist_like_1d[i] = 0;
355 HISTOGRAMS.hist_unlike_1d[i] = 0;
356 HISTOGRAMS.htmp_like_1d[i] = 0;
357 HISTOGRAMS.htmp_unlike_1d[i] = 0;
358 HISTOGRAMS.href_like_1d[i] = 0;
359 HISTOGRAMS.href_unlike_1d[i] = 0;
360 HISTOGRAMS.hinc_like_1d[i] = 0;
361 HISTOGRAMS.hinc_unlike_1d[i] = 0;
362 }
363
364 for (Int_t i =0; i<MAX_H_3D;i++)
365 for (Int_t j =0; j<MAX_H_3D;j++)
366 for (Int_t k =0; k<MAX_H_3D;k++)
367 {
368 HISTOGRAMS.hist_like_3d_fine[i][j][k] = 0;
369 HISTOGRAMS.hist_unlike_3d_fine[i][j][k] = 0;
370 HISTOGRAMS.hist_like_3d_coarse[i][j][k] = 0;
371 HISTOGRAMS.hist_unlike_3d_coarse[i][j][k] = 0;
372 HISTOGRAMS.htmp_like_3d_fine[i][j][k] = 0;
373 HISTOGRAMS.htmp_unlike_3d_fine[i][j][k] = 0;
374 HISTOGRAMS.htmp_like_3d_coarse[i][j][k] = 0;
375 HISTOGRAMS.htmp_unlike_3d_coarse[i][j][k] = 0;
376 HISTOGRAMS.href_like_3d_fine[i][j][k] = 0;
377 HISTOGRAMS.href_unlike_3d_fine[i][j][k] = 0;
378 HISTOGRAMS.href_like_3d_coarse[i][j][k] = 0;
379 HISTOGRAMS.href_unlike_3d_coarse[i][j][k] = 0;
380 HISTOGRAMS.hinc_like_3d_fine[i][j][k] = 0;
381 HISTOGRAMS.hinc_unlike_3d_fine[i][j][k] = 0;
382 HISTOGRAMS.hinc_like_3d_coarse[i][j][k] = 0;
383 HISTOGRAMS.hinc_unlike_3d_coarse[i][j][k] = 0;
384 }
385/*********************************************/
386
387
388/*********************************************/
389
390// cout<<" FINISHED"<<endl;
391
392}
393
394
395/*****************************************************************************************/
396
397
398Int_t THBTprocessor::ImportParticles(TClonesArray *particles, Option_t *option)
399 {
400 //Copy particle data into TClonesArray
401 if (particles == 0) return 0;
402 TClonesArray &Particles = *particles;
403 Particles.Clear();
404
405 Int_t nrpart = 0;
406 for (Int_t i=0; i < TRK_MAXLEN; i++)
407 {
408
409
410 if (TRACK1.trk_E[i] == 0.) continue;
411
412 Float_t px = TRACK1.trk_px[i];
413 Float_t py = TRACK1.trk_py[i];
414 Float_t pz = TRACK1.trk_pz[i];
415// Float_t pE = TRACK.trk_E[i];
416 Float_t mass = PARTICLE.part_mass[TRACK1.trk_ge_pid[i]];
417
418 new(Particles[nrpart++]) TParticle(0,0,0,0,0,0,px,py,pz,
419 TMath::Sqrt(mass*mass+px*px+py*py+pz*pz),
420 0,0,0,0);
421 }
422 return nrpart;
423 }
424
425/*****************************************************************************************/
426
427void THBTprocessor::PrintEvent()
428 {
429 //Prints all particles (common block data)
430 cout<<"Print Event"<<endl;
431 for (Int_t i=0; i<TRK_MAXLEN;i++)
432 {
433 if(TRACK1.trk_E[i]==0.) continue;
434
435 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;
436 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;
437 }
438 }
439
440
441/*****************************************************************************************/