]>
Commit | Line | Data |
---|---|---|
45575004 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | ||
1c01b4f8 | 16 | // $Id: AliCollider.cxx,v 1.6 2003/08/29 09:05:11 nick Exp $ |
fdbea0ce | 17 | |
18 | /////////////////////////////////////////////////////////////////////////// | |
19 | // Class AliCollider | |
20 | // Pythia based universal physics event generator. | |
4b570fab | 21 | // This class is derived from TPythia6 and has some extensions to |
fdbea0ce | 22 | // support also generation of nucleus-nucleus interactions and to allow |
23 | // investigation of the effect of detector resolving power. | |
24 | // Furthermore, the produced event information is provided in a format | |
25 | // using the AliEvent structure. | |
26 | // For the produced AliTrack objects, the particle ID code is set to the | |
27 | // Pythia KF value, which is compatible with the PDG identifier. | |
28 | // This will allow a direct analysis of the produced data using the | |
29 | // Ralice physics analysis tools. | |
30 | // | |
31 | // For further details concerning the produced output structure, | |
32 | // see the docs of the memberfunctions SetVertexMode and SetResolution. | |
33 | // | |
34 | // Example job of minimum biased Pb+Pb interactions : | |
35 | // -------------------------------------------------- | |
36 | // { | |
37 | // gSystem->Load("libEG"); | |
38 | // gSystem->Load("libEGPythia6"); | |
39 | // gSystem->Load("ralice"); | |
40 | // | |
41 | // AliCollider* gen=new AliCollider(); | |
42 | // | |
43 | // gen->SetOutputFile("test.root"); | |
44 | // gen->SetVertexMode(3); | |
45 | // gen->SetResolution(1e-4); // 1 micron vertex resolution | |
46 | // | |
47 | // gen->SetRunNumber(1); | |
48 | // | |
49 | // Int_t zp=82; | |
50 | // Int_t ap=208; | |
51 | // Int_t zt=82; | |
52 | // Int_t at=208; | |
53 | // | |
54 | // gen->Init("fixt",zp,ap,zt,at,158); | |
55 | // | |
56 | // Int_t nevents=5; | |
57 | // | |
58 | // AliRandom rndm; | |
59 | // Float_t* rans=new Float_t[nevents]; | |
60 | // rndm.Uniform(rans,nevents,2,ap+at); | |
61 | // Int_t npart; | |
62 | // for (Int_t i=0; i<nevents; i++) | |
63 | // { | |
64 | // npart=rans[i]; | |
65 | // gen->MakeEvent(npart); | |
66 | // | |
67 | // AliEvent* evt=gen->GetEvent(); | |
68 | // | |
69 | // evt->List(); | |
70 | // } | |
71 | // | |
72 | // gen->EndRun(); | |
73 | // } | |
74 | // | |
75 | // | |
76 | // Example job of a cosmic nu+p atmospheric interaction. | |
77 | // ----------------------------------------------------- | |
78 | // { | |
79 | // gSystem->Load("libEG"); | |
80 | // gSystem->Load("libEGPythia6"); | |
81 | // gSystem->Load("ralice"); | |
82 | // | |
83 | // AliCollider* gen=new AliCollider(); | |
84 | // | |
85 | // gen->SetOutputFile("test.root"); | |
86 | // | |
87 | // gen->SetRunNumber(1); | |
88 | // | |
89 | // gen->Init("fixt","nu_mu","p",1e11); | |
90 | // | |
91 | // Int_t nevents=10; | |
92 | // | |
93 | // for (Int_t i=0; i<nevents; i++) | |
94 | // { | |
95 | // gen->MakeEvent(0,1); | |
96 | // | |
97 | // AliEvent* evt=gen->GetEvent(); | |
98 | // | |
84bb7c66 | 99 | // evt->Data(); |
fdbea0ce | 100 | // } |
101 | // | |
102 | // gen->EndRun(); | |
103 | // } | |
104 | // | |
105 | // | |
106 | //--- Author: Nick van Eijndhoven 22-nov-2002 Utrecht University | |
1c01b4f8 | 107 | //- Modified: NvE $Date: 2003/08/29 09:05:11 $ Utrecht University |
fdbea0ce | 108 | /////////////////////////////////////////////////////////////////////////// |
109 | ||
110 | #include "AliCollider.h" | |
c72198f1 | 111 | #include "Riostream.h" |
fdbea0ce | 112 | |
113 | ClassImp(AliCollider) // Class implementation to enable ROOT I/O | |
114 | ||
c72198f1 | 115 | AliCollider::AliCollider() : TPythia6() |
fdbea0ce | 116 | { |
117 | // Default constructor. | |
118 | // All variables initialised to default values. | |
119 | fVertexmode=0; // No vertex structure creation | |
120 | fResolution=1e-5; // Standard resolution is 0.1 micron | |
121 | fRunnum=0; | |
122 | fEventnum=0; | |
123 | fPrintfreq=1; | |
124 | ||
125 | fEvent=0; | |
126 | ||
127 | fFrame="none"; | |
128 | fWin=0; | |
129 | ||
130 | fNucl=0; | |
131 | fZproj=0; | |
132 | fAproj=0; | |
133 | fZtarg=0; | |
134 | fAtarg=0; | |
135 | fFracpp=0; | |
136 | fFracnp=0; | |
137 | fFracpn=0; | |
138 | fFracnn=0; | |
139 | ||
140 | fOutFile=0; | |
141 | fOutTree=0; | |
142 | } | |
143 | /////////////////////////////////////////////////////////////////////////// | |
144 | AliCollider::~AliCollider() | |
145 | { | |
146 | // Default destructor | |
147 | if (fEvent) | |
148 | { | |
149 | delete fEvent; | |
150 | fEvent=0; | |
151 | } | |
152 | if (fOutFile) | |
153 | { | |
154 | delete fOutFile; | |
155 | fOutFile=0; | |
156 | } | |
157 | if (fOutTree) | |
158 | { | |
159 | delete fOutTree; | |
160 | fOutTree=0; | |
161 | } | |
162 | } | |
163 | /////////////////////////////////////////////////////////////////////////// | |
164 | void AliCollider::SetOutputFile(TString s) | |
165 | { | |
166 | // Create the output file containing all the data in ROOT output format. | |
167 | if (fOutFile) | |
168 | { | |
169 | delete fOutFile; | |
170 | fOutFile=0; | |
171 | } | |
172 | fOutFile=new TFile(s.Data(),"RECREATE","AliCollider data"); | |
173 | ||
174 | if (fOutTree) | |
175 | { | |
176 | delete fOutTree; | |
177 | fOutTree=0; | |
178 | } | |
179 | fOutTree=new TTree("T","AliCollider event data"); | |
180 | ||
181 | Int_t bsize=32000; | |
182 | Int_t split=0; | |
183 | fOutTree->Branch("Events","AliEvent",&fEvent,bsize,split); | |
184 | } | |
185 | /////////////////////////////////////////////////////////////////////////// | |
186 | void AliCollider::SetVertexMode(Int_t mode) | |
187 | { | |
188 | // Set the mode of the vertex structure creation. | |
189 | // | |
190 | // By default all generated tracks will only appear in the AliEvent | |
191 | // structure without any primary (and secondary) vertex structure. | |
192 | // The user can build the vertex structure if he/she wants by means | |
193 | // of the beginpoint location of each AliTrack. | |
194 | // | |
195 | // However, one can also let AliCollider automatically create | |
196 | // the primary (and secondary) vertex structure(s). | |
197 | // In this case the primary vertex is given Id=1 and all sec. vertices | |
198 | // are given Id's 2,3,4,.... | |
199 | // All vertices are created as standalone entities in the AliEvent structure | |
200 | // without any linking between the various vertices. | |
201 | // For this automated process, the user-selected resolution | |
202 | // (see SetResolution) is used to decide whether or not certain vertex | |
203 | // locations can be resolved. | |
204 | // In case no vertex creation is selected (i.e. the default mode=0), | |
205 | // the value of the resolution is totally irrelevant. | |
206 | // | |
207 | // The user can also let AliCollider automatically connect the sec. vertices | |
208 | // to the primary vertex (i.e. mode=3). This process will also automatically | |
209 | // generate the tracks connecting the vertices. | |
210 | // Note that the result of the mode=3 operation may be very sensitive to | |
211 | // the resolution parameter. Therefore, no attempt is made to distinguish | |
212 | // between secondary, tertiary etc... vertices. All sec. vertices are | |
213 | // linked to the primary one. | |
214 | // | |
215 | // Irrespective of the selected mode, all generated tracks can be obtained | |
216 | // directly from the AliEvent structure. | |
217 | // In case (sec.) vertex creation is selected, all generated vertices can | |
218 | // also be obtained directly from the AliEvent structure. | |
219 | // These (sec.) vertices contain only the corresponding pointers to the various | |
220 | // tracks which are stored in the AliEvent structure. | |
221 | // | |
222 | // Overview of vertex creation modes : | |
223 | // ----------------------------------- | |
224 | // mode = 0 ==> No vertex structure will be created | |
225 | // 1 ==> Only primary vertex structure will be created | |
226 | // 2 ==> Unconnected primary and secondary vertices will be created | |
227 | // 3 ==> Primary and secondary vertices will be created where all the | |
228 | // sec. vertices will be connected to the primary vertex. | |
229 | // Also the vertex connecting tracks will be automatically | |
230 | // generated. | |
231 | // | |
232 | if (mode<0 || mode >3) | |
233 | { | |
234 | cout << " *AliCollider::SetVertexMode* Invalid argument mode : " << mode << endl; | |
235 | fVertexmode=0; | |
236 | } | |
237 | else | |
238 | { | |
239 | fVertexmode=mode; | |
240 | } | |
241 | } | |
242 | /////////////////////////////////////////////////////////////////////////// | |
243 | Int_t AliCollider::GetVertexMode() | |
244 | { | |
245 | // Provide the current mode for vertex structure creation. | |
246 | return fVertexmode; | |
247 | } | |
248 | /////////////////////////////////////////////////////////////////////////// | |
249 | void AliCollider::SetResolution(Double_t res) | |
250 | { | |
251 | // Set the resolution (in cm) for resolving (sec.) vertices. | |
252 | // By default this resolution is set to 0.1 micron. | |
253 | // Note : In case no vertex creation has been selected, the value of | |
254 | // the resolution is totally irrelevant. | |
255 | fResolution=fabs(res); | |
256 | } | |
257 | /////////////////////////////////////////////////////////////////////////// | |
258 | Double_t AliCollider::GetResolution() | |
259 | { | |
260 | // Provide the current resolution (in cm) for resolving (sec.) vertices. | |
261 | return fResolution; | |
262 | } | |
263 | /////////////////////////////////////////////////////////////////////////// | |
264 | void AliCollider::SetRunNumber(Int_t run) | |
265 | { | |
266 | // Set the user defined run number. | |
267 | // By default the run number is set to 0. | |
268 | fRunnum=run; | |
269 | } | |
270 | /////////////////////////////////////////////////////////////////////////// | |
271 | Int_t AliCollider::GetRunNumber() | |
272 | { | |
273 | // Provide the user defined run number. | |
274 | return fRunnum; | |
275 | } | |
276 | /////////////////////////////////////////////////////////////////////////// | |
277 | void AliCollider::SetPrintFreq(Int_t n) | |
278 | { | |
279 | // Set the print frequency for every 'n' events. | |
280 | // By default the printfrequency is set to 1 (i.e. every event). | |
281 | fPrintfreq=n; | |
282 | } | |
283 | /////////////////////////////////////////////////////////////////////////// | |
284 | Int_t AliCollider::GetPrintFreq() | |
285 | { | |
286 | // Provide the user selected print frequency. | |
287 | return fPrintfreq; | |
288 | } | |
289 | /////////////////////////////////////////////////////////////////////////// | |
290 | void AliCollider::Init(char* frame,char* beam,char* target,Float_t win) | |
291 | { | |
292 | // Initialisation of the underlying Pythia generator package. | |
293 | // This routine just invokes TPythia6::Initialize(...) and the arguments | |
294 | // have the corresponding meaning. | |
295 | // The event number is reset to 0. | |
296 | fEventnum=0; | |
297 | fNucl=0; | |
298 | fFrame=frame; | |
299 | fWin=win; | |
300 | Initialize(frame,beam,target,win); | |
c72198f1 | 301 | |
302 | cout << " *AliCollider::Init* Standard Pythia initialisation." << endl; | |
303 | cout << " Beam particle : " << beam << " Target particle : " << target | |
304 | << " Frame = " << frame << " Energy = " << win | |
305 | << endl; | |
fdbea0ce | 306 | } |
307 | /////////////////////////////////////////////////////////////////////////// | |
308 | void AliCollider::Init(char* frame,Int_t zp,Int_t ap,Int_t zt,Int_t at,Float_t win) | |
309 | { | |
310 | // Initialisation of the underlying Pythia generator package for the generation | |
311 | // of nucleus-nucleus interactions. | |
312 | // In addition to the Pythia standard arguments 'frame' and 'win', the user | |
da17f667 | 313 | // can specify here (Z,A) values of the projectile and target nuclei. |
314 | // | |
315 | // Note : The 'win' value denotes either the cms energy per nucleon-nucleon collision | |
316 | // (i.e. frame="cms") or the momentum per nucleon in all other cases. | |
317 | // | |
fdbea0ce | 318 | // The event number is reset to 0. |
319 | fEventnum=0; | |
320 | fNucl=1; | |
321 | fFrame=frame; | |
322 | fWin=win; | |
323 | fZproj=0; | |
324 | fAproj=0; | |
325 | fZtarg=0; | |
326 | fAtarg=0; | |
327 | fFracpp=0; | |
328 | fFracnp=0; | |
329 | fFracpn=0; | |
330 | fFracnn=0; | |
331 | ||
332 | if (ap<1 || at<1 || zp>ap || zt>at) | |
333 | { | |
334 | cout << " *AliCollider::Init* Invalid input value(s). Zproj = " << zp | |
335 | << " Aproj = " << ap << " Ztarg = " << zt << " Atarg = " << at << endl; | |
336 | return; | |
337 | } | |
338 | ||
339 | fZproj=zp; | |
340 | fAproj=ap; | |
341 | fZtarg=zt; | |
342 | fAtarg=at; | |
343 | ||
344 | cout << " *AliCollider::Init* Nucleus-Nucleus generator initialisation." << endl; | |
345 | cout << " Zproj = " << zp << " Aproj = " << ap << " Ztarg = " << zt << " Atarg = " << at | |
346 | << " Frame = " << frame << " Energy = " << win | |
347 | << endl; | |
348 | } | |
349 | /////////////////////////////////////////////////////////////////////////// | |
350 | void AliCollider::GetFractions(Float_t zp,Float_t ap,Float_t zt,Float_t at) | |
351 | { | |
352 | // Determine the fractions for the various N-N collision processes. | |
353 | // The various processes are : p+p, n+p, p+n and n+n. | |
354 | if (zp<0) zp=0; | |
355 | if (zt<0) zt=0; | |
356 | ||
357 | fFracpp=0; | |
358 | fFracnp=0; | |
359 | fFracpn=0; | |
360 | fFracnn=0; | |
361 | ||
362 | if (ap>0 && at>0) | |
363 | { | |
364 | fFracpp=(zp/ap)*(zt/at); | |
365 | fFracnp=(1.-zp/ap)*(zt/at); | |
366 | fFracpn=(zp/ap)*(1.-zt/at); | |
367 | fFracnn=(1.-zp/ap)*(1.-zt/at); | |
368 | } | |
369 | } | |
370 | /////////////////////////////////////////////////////////////////////////// | |
371 | void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit) | |
372 | { | |
373 | // Generate one event. | |
374 | // In case of a nucleus-nucleus interaction, the argument 'npt' denotes | |
375 | // the number of participant nucleons. | |
376 | // In case of a standard Pythia run for 'elementary' particle interactions, | |
377 | // the value of npt is totally irrelevant. | |
da17f667 | 378 | // |
fdbea0ce | 379 | // The argument 'mlist' denotes the list mode used for Pylist(). |
da17f667 | 380 | // Note : mlist<0 suppresses the invokation of Pylist(). |
381 | // By default, no listing is produced (i.e. mlist=-1). | |
382 | // | |
c72198f1 | 383 | // The argument 'medit' denotes the edit mode used for Pyedit(). |
384 | // Note : medit<0 suppresses the invokation of Pyedit(). | |
385 | // By default, only 'stable' final particles are kept (i.e. medit=1). | |
386 | // | |
da17f667 | 387 | // In the case of a standard Pythia run concerning 'elementary' particle |
388 | // interactions, the projectile and target particle ID's for the created | |
389 | // event structure are set to the corresponding Pythia KF codes. | |
390 | // All the A and Z values are in that case set to zero. | |
391 | // In case of a nucleus-nucleus interaction, the proper A and Z values for | |
392 | // the projectile and target particles are set in the event structure. | |
393 | // However, in this case both particle ID's are set to zero. | |
fdbea0ce | 394 | |
395 | fEventnum++; | |
396 | ||
397 | // Counters for the various (proj,targ) combinations : p+p, n+p, p+n and n+n | |
398 | Int_t ncols[4]={0,0,0,0}; | |
399 | ||
c72198f1 | 400 | Int_t ncol=1; |
fdbea0ce | 401 | if (fNucl) |
402 | { | |
403 | if (npt<1 || npt>(fAproj+fAtarg)) | |
404 | { | |
405 | cout << " *AliCollider::MakeEvent* Invalid input value. npt = " << npt | |
406 | << " Aproj = " << fAproj << " Atarg = " << fAtarg << endl; | |
407 | return; | |
408 | } | |
409 | ||
410 | // Determine the number of nucleon-nucleon collisions | |
c72198f1 | 411 | ncol=npt/2; |
fdbea0ce | 412 | if (npt%2 && fRan.Uniform()>0.5) ncol+=1; |
413 | ||
414 | // Determine the number of the various types of N+N interactions | |
415 | Int_t zp=fZproj; | |
416 | Int_t ap=fAproj; | |
417 | Int_t zt=fZtarg; | |
418 | Int_t at=fAtarg; | |
419 | Int_t maxa=2; // Indicator whether proj (1) or target (2) has maximal A left | |
420 | if (ap>at) maxa=1; | |
421 | Float_t* rans=new Float_t[ncol]; | |
422 | fRan.Uniform(rans,ncol); | |
423 | Float_t rndm=0; | |
424 | for (Int_t i=0; i<ncol; i++) | |
425 | { | |
426 | GetFractions(zp,ap,zt,at); | |
427 | rndm=rans[i]; | |
428 | if (rndm<=fFracpp) // p+p interaction | |
429 | { | |
430 | ncols[0]++; | |
4b570fab | 431 | if (maxa==2) |
fdbea0ce | 432 | { |
433 | at--; | |
434 | zt--; | |
435 | } | |
436 | else | |
437 | { | |
438 | ap--; | |
439 | zp--; | |
440 | } | |
441 | } | |
442 | if (rndm>fFracpp && rndm<=(fFracpp+fFracnp)) // n+p interaction | |
443 | { | |
444 | ncols[1]++; | |
4b570fab | 445 | if (maxa==2) |
fdbea0ce | 446 | { |
447 | at--; | |
448 | zt--; | |
449 | } | |
450 | else | |
451 | { | |
452 | ap--; | |
453 | } | |
454 | } | |
455 | if (rndm>(fFracpp+fFracnp) && rndm<=(fFracpp+fFracnp+fFracpn)) // p+n interaction | |
456 | { | |
457 | ncols[2]++; | |
4b570fab | 458 | if (maxa==2) |
fdbea0ce | 459 | { |
460 | at--; | |
461 | } | |
462 | else | |
463 | { | |
464 | ap--; | |
465 | zp--; | |
466 | } | |
467 | } | |
468 | if (rndm>(fFracpp+fFracnp+fFracpn)) // n+n interaction | |
469 | { | |
470 | ncols[3]++; | |
4b570fab | 471 | if (maxa==2) |
fdbea0ce | 472 | { |
473 | at--; | |
474 | } | |
475 | else | |
476 | { | |
477 | ap--; | |
478 | } | |
479 | } | |
480 | } | |
481 | delete [] rans; | |
c72198f1 | 482 | } |
fdbea0ce | 483 | |
1c01b4f8 | 484 | if (!(fEventnum%fPrintfreq)) |
485 | { | |
486 | cout << " *AliCollider::MakeEvent* Run : " << fRunnum << " Event : " << fEventnum | |
487 | << endl; | |
488 | if (fNucl) | |
fdbea0ce | 489 | { |
1c01b4f8 | 490 | cout << " npart = " << npt << " ncol = " << ncol |
491 | << " ncolpp = " << ncols[0] << " ncolnp = " << ncols[1] | |
492 | << " ncolpn = " << ncols[2] << " ncolnn = " << ncols[3] << endl; | |
fdbea0ce | 493 | } |
1c01b4f8 | 494 | } |
fdbea0ce | 495 | |
fdbea0ce | 496 | if (!fEvent) |
497 | { | |
498 | fEvent=new AliEvent(); | |
499 | fEvent->SetOwner(); | |
500 | } | |
501 | ||
502 | fEvent->Reset(); | |
503 | fEvent->SetRunNumber(fRunnum); | |
504 | fEvent->SetEventNumber(fEventnum); | |
505 | ||
da17f667 | 506 | AliTrack t; |
507 | Ali3Vector p; | |
508 | AliPosition r,rx; | |
509 | Float_t v[3]; | |
fdbea0ce | 510 | AliVertex vert; |
da17f667 | 511 | |
fdbea0ce | 512 | if (fVertexmode) |
513 | { | |
514 | // Make sure the primary vertex gets correct location and Id=1 | |
da17f667 | 515 | v[0]=0; |
516 | v[1]=0; | |
517 | v[2]=0; | |
518 | r.SetPosition(v,"car"); | |
519 | v[0]=fResolution; | |
520 | v[1]=fResolution; | |
521 | v[2]=fResolution; | |
522 | r.SetPositionErrors(v,"car"); | |
523 | ||
fdbea0ce | 524 | vert.SetId(1); |
525 | vert.SetTrackCopy(0); | |
526 | vert.SetVertexCopy(0); | |
da17f667 | 527 | vert.SetPosition(r); |
fdbea0ce | 528 | fEvent->AddVertex(vert,0); |
529 | } | |
530 | ||
c72198f1 | 531 | Int_t kf=0; |
fdbea0ce | 532 | Float_t charge=0,mass=0; |
533 | ||
fdbea0ce | 534 | Int_t ntypes=4; |
535 | ||
536 | // Singular settings for a normal Pythia elementary particle interation | |
537 | if (!fNucl) | |
538 | { | |
539 | ntypes=1; | |
540 | ncols[0]=1; | |
541 | } | |
542 | ||
543 | // Generate all the various collisions | |
da17f667 | 544 | Int_t first=1; // Flag to indicate the first collision process |
545 | Double_t pnucl; | |
fdbea0ce | 546 | Int_t npart=0,ntk=0; |
547 | Double_t dist=0; | |
548 | for (Int_t itype=0; itype<ntypes; itype++) | |
549 | { | |
550 | if (fNucl) | |
551 | { | |
552 | if (itype==0 && ncols[itype]) Initialize(fFrame,"p","p",fWin); | |
553 | if (itype==1 && ncols[itype]) Initialize(fFrame,"n","p",fWin); | |
554 | if (itype==2 && ncols[itype]) Initialize(fFrame,"p","n",fWin); | |
555 | if (itype==3 && ncols[itype]) Initialize(fFrame,"n","n",fWin); | |
556 | } | |
557 | for (Int_t jcol=0; jcol<ncols[itype]; jcol++) | |
558 | { | |
559 | GenerateEvent(); | |
560 | ||
da17f667 | 561 | if (first) // Store projectile and target information in the event structure |
562 | { | |
563 | if (fNucl) | |
564 | { | |
565 | v[0]=GetP(1,1); | |
566 | v[1]=GetP(1,2); | |
567 | v[2]=GetP(1,3); | |
568 | pnucl=sqrt(v[0]*v[0]+v[1]*v[1]+v[2]*v[2]); | |
569 | fEvent->SetProjectile(fAproj,fZproj,pnucl); | |
570 | v[0]=GetP(2,1); | |
571 | v[1]=GetP(2,2); | |
572 | v[2]=GetP(2,3); | |
573 | pnucl=sqrt(v[0]*v[0]+v[1]*v[1]+v[2]*v[2]); | |
574 | fEvent->SetTarget(fAtarg,fZtarg,pnucl); | |
575 | } | |
576 | else | |
577 | { | |
578 | v[0]=GetP(1,1); | |
579 | v[1]=GetP(1,2); | |
580 | v[2]=GetP(1,3); | |
581 | pnucl=sqrt(v[0]*v[0]+v[1]*v[1]+v[2]*v[2]); | |
582 | kf=GetK(1,2); | |
583 | fEvent->SetProjectile(0,0,pnucl,kf); | |
584 | v[0]=GetP(2,1); | |
585 | v[1]=GetP(2,2); | |
586 | v[2]=GetP(2,3); | |
587 | pnucl=sqrt(v[0]*v[0]+v[1]*v[1]+v[2]*v[2]); | |
588 | kf=GetK(2,2); | |
589 | fEvent->SetTarget(0,0,pnucl,kf); | |
590 | } | |
591 | first=0; | |
592 | } | |
593 | ||
594 | if (medit >= 0) Pyedit(medit); // Define which particles are to be kept | |
fdbea0ce | 595 | |
da17f667 | 596 | if (mlist >= 0) Pylist(mlist); |
fdbea0ce | 597 | |
c72198f1 | 598 | npart=GetN(); |
599 | for (Int_t jpart=1; jpart<=npart; jpart++) | |
fdbea0ce | 600 | { |
c72198f1 | 601 | kf=GetK(jpart,2); |
602 | charge=Pychge(kf)/3.; | |
603 | mass=GetP(jpart,5); | |
fdbea0ce | 604 | |
605 | // 3-momentum in GeV/c | |
c72198f1 | 606 | v[0]=GetP(jpart,1); |
607 | v[1]=GetP(jpart,2); | |
608 | v[2]=GetP(jpart,3); | |
fdbea0ce | 609 | p.SetVector(v,"car"); |
610 | ||
611 | // Production location in cm. | |
c72198f1 | 612 | v[0]=GetV(jpart,1)/10; |
613 | v[1]=GetV(jpart,2)/10; | |
614 | v[2]=GetV(jpart,3)/10; | |
da17f667 | 615 | r.SetPosition(v,"car"); |
fdbea0ce | 616 | |
617 | ntk++; | |
618 | ||
619 | t.Reset(); | |
620 | t.SetId(ntk); | |
621 | t.SetParticleCode(kf); | |
622 | t.SetCharge(charge); | |
623 | t.SetMass(mass); | |
624 | t.Set3Momentum(p); | |
625 | t.SetBeginPoint(r); | |
626 | ||
627 | fEvent->AddTrack(t); | |
628 | ||
629 | // Build the vertex structures if requested | |
630 | if (fVertexmode) | |
631 | { | |
632 | // Check if track belongs within the resolution to an existing vertex | |
633 | Int_t add=0; | |
634 | for (Int_t jv=1; jv<=fEvent->GetNvertices(); jv++) | |
635 | { | |
636 | AliVertex* vx=fEvent->GetVertex(jv); | |
637 | if (vx) | |
638 | { | |
639 | rx=vx->GetPosition(); | |
640 | dist=rx.GetDistance(r); | |
641 | if (dist < fResolution) | |
642 | { | |
643 | AliTrack* tx=fEvent->GetIdTrack(ntk); | |
644 | if (tx) | |
645 | { | |
646 | vx->AddTrack(tx); | |
647 | add=1; | |
648 | } | |
649 | } | |
650 | } | |
651 | if (add) break; // No need to look further for vertex candidates | |
652 | } | |
653 | ||
654 | // If track was not close enough to an existing vertex | |
655 | // a new secondary vertex is created | |
656 | if (!add && fVertexmode>1) | |
657 | { | |
658 | AliTrack* tx=fEvent->GetIdTrack(ntk); | |
659 | if (tx) | |
660 | { | |
da17f667 | 661 | v[0]=fResolution; |
662 | v[1]=fResolution; | |
663 | v[2]=fResolution; | |
664 | r.SetPositionErrors(v,"car"); | |
fdbea0ce | 665 | vert.Reset(); |
666 | vert.SetTrackCopy(0); | |
667 | vert.SetVertexCopy(0); | |
668 | vert.SetId((fEvent->GetNvertices())+1); | |
669 | vert.SetPosition(r); | |
670 | vert.AddTrack(tx); | |
671 | fEvent->AddVertex(vert,0); | |
672 | } | |
673 | } | |
674 | } | |
675 | } // End of loop over the produced particles for each collision | |
676 | } // End of loop over number of collisions for each type | |
677 | } // End of loop over collision types | |
678 | ||
679 | // Link sec. vertices to primary if requested | |
680 | // Note that also the connecting tracks are automatically created | |
681 | if (fVertexmode>2) | |
682 | { | |
683 | AliVertex* vp=fEvent->GetIdVertex(1); // Primary vertex | |
684 | if (vp) | |
685 | { | |
686 | for (Int_t i=2; i<=fEvent->GetNvertices(); i++) | |
687 | { | |
688 | AliVertex* vx=fEvent->GetVertex(i); | |
689 | if (vx) | |
690 | { | |
691 | if (vx->GetId() != 1) vp->AddVertex(vx); | |
692 | } | |
693 | } | |
694 | } | |
695 | } | |
696 | ||
c72198f1 | 697 | if (mlist && !(fEventnum%fPrintfreq)) cout << endl; // Create empty output line after the event |
fdbea0ce | 698 | if (fOutTree) fOutTree->Fill(); |
699 | } | |
700 | /////////////////////////////////////////////////////////////////////////// | |
701 | AliEvent* AliCollider::GetEvent() | |
702 | { | |
703 | // Provide pointer to the generated event structure. | |
704 | return fEvent; | |
705 | } | |
706 | /////////////////////////////////////////////////////////////////////////// | |
707 | void AliCollider::EndRun() | |
708 | { | |
709 | // Properly close the output file (if needed). | |
710 | if (fOutFile) | |
711 | { | |
712 | fOutFile->Write(); | |
713 | fOutFile->Close(); | |
714 | cout << " *AliCollider::EndRun* Output file correctly closed." << endl; | |
715 | } | |
716 | } | |
717 | /////////////////////////////////////////////////////////////////////////// |