RALICE/AliCollider.cxx

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