25-nov-2002 NvE User defined particle id code introduced in AliTrack.
authornick <nick@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 27 Nov 2002 21:25:52 +0000 (21:25 +0000)
committernick <nick@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 27 Nov 2002 21:25:52 +0000 (21:25 +0000)
27-nov-2002 NvE New class AliCollider introduced and RALICEHeaders.h and RALICELinkdef.h
                updated accordingly.

RALICE/AliCollider.cxx [new file with mode: 0644]
RALICE/AliCollider.h [new file with mode: 0644]
RALICE/AliTrack.cxx
RALICE/AliTrack.h
RALICE/RALICEHeaders.h
RALICE/RALICELinkDef.h
RALICE/history.txt

diff --git a/RALICE/AliCollider.cxx b/RALICE/AliCollider.cxx
new file mode 100644 (file)
index 0000000..daf5c15
--- /dev/null
@@ -0,0 +1,640 @@
+// $Id$
+
+///////////////////////////////////////////////////////////////////////////
+// Class AliCollider
+// Pythia based universal physics event generator.
+// This event class is derived from TPythia6 and has some extensions to
+// support also generation of nucleus-nucleus interactions and to allow
+// investigation of the effect of detector resolving power.
+// Furthermore, the produced event information is provided in a format
+// using the AliEvent structure.
+// For the produced AliTrack objects, the particle ID code is set to the
+// Pythia KF value, which is compatible with the PDG identifier.
+// This will allow a direct analysis of the produced data using the
+// Ralice physics analysis tools.
+//
+// For further details concerning the produced output structure,
+// see the docs of the memberfunctions SetVertexMode and SetResolution.
+//
+// Example job of minimum biased Pb+Pb interactions :
+// --------------------------------------------------
+// {
+//  gSystem->Load("libEG");
+//  gSystem->Load("libEGPythia6");
+//  gSystem->Load("ralice");
+//
+//  AliCollider* gen=new AliCollider();
+//
+//  gen->SetOutputFile("test.root");
+//  gen->SetVertexMode(3);    
+//  gen->SetResolution(1e-4); // 1 micron vertex resolution
+//
+//  gen->SetRunNumber(1);
+//
+//  Int_t zp=82;
+//  Int_t ap=208;
+//  Int_t zt=82;
+//  Int_t at=208;
+//
+//  gen->Init("fixt",zp,ap,zt,at,158);
+//
+//  Int_t nevents=5;
+//
+//  AliRandom rndm;
+//  Float_t* rans=new Float_t[nevents];
+//  rndm.Uniform(rans,nevents,2,ap+at);
+//  Int_t npart;
+//  for (Int_t i=0; i<nevents; i++)
+//  {
+//   npart=rans[i];
+//   gen->MakeEvent(npart);
+//
+//   AliEvent* evt=gen->GetEvent();
+//  
+//   evt->List();
+//  }
+//
+//  gen->EndRun();
+// }
+//
+//
+// Example job of a cosmic nu+p atmospheric interaction.
+// -----------------------------------------------------
+// {
+//  gSystem->Load("libEG");
+//  gSystem->Load("libEGPythia6");
+//  gSystem->Load("ralice");
+//
+//  AliCollider* gen=new AliCollider();
+//
+//  gen->SetOutputFile("test.root");
+//
+//  gen->SetRunNumber(1);
+//
+//  gen->Init("fixt","nu_mu","p",1e11);
+//
+//  Int_t nevents=10;
+//
+//  for (Int_t i=0; i<nevents; i++)
+//  {
+//   gen->MakeEvent(0,1);
+//
+//   AliEvent* evt=gen->GetEvent();
+//  
+//   evt->Info();
+//  }
+//
+//  gen->EndRun();
+// }
+//
+//
+//--- Author: Nick van Eijndhoven 22-nov-2002 Utrecht University
+//- Modified: NvE $Date$ Utrecht University
+///////////////////////////////////////////////////////////////////////////
+
+#include "AliCollider.h"
+ClassImp(AliCollider) // Class implementation to enable ROOT I/O
+AliCollider::AliCollider()
+{
+// Default constructor.
+// All variables initialised to default values.
+ fVertexmode=0;    // No vertex structure creation
+ fResolution=1e-5; // Standard resolution is 0.1 micron
+ fRunnum=0;
+ fEventnum=0;
+ fPrintfreq=1;
+
+ fEvent=0;
+
+ fFrame="none";
+ fWin=0;
+
+ fNucl=0;
+ fZproj=0;
+ fAproj=0;
+ fZtarg=0;
+ fAtarg=0;
+ fFracpp=0;
+ fFracnp=0;
+ fFracpn=0;
+ fFracnn=0;
+
+ fOutFile=0;
+ fOutTree=0;
+}
+///////////////////////////////////////////////////////////////////////////
+AliCollider::~AliCollider()
+{
+// Default destructor
+ if (fEvent)
+ {
+  delete fEvent;
+  fEvent=0;
+ }
+ if (fOutFile)
+ {
+  delete fOutFile;
+  fOutFile=0;
+ }
+ if (fOutTree)
+ {
+  delete fOutTree;
+  fOutTree=0;
+ }
+}
+///////////////////////////////////////////////////////////////////////////
+void AliCollider::SetOutputFile(TString s)
+{
+// Create the output file containing all the data in ROOT output format.
+ if (fOutFile)
+ {
+  delete fOutFile;
+  fOutFile=0;
+ }
+ fOutFile=new TFile(s.Data(),"RECREATE","AliCollider data");
+
+ if (fOutTree)
+ {
+  delete fOutTree;
+  fOutTree=0;
+ }
+ fOutTree=new TTree("T","AliCollider event data");
+
+ Int_t bsize=32000;
+ Int_t split=0;
+ fOutTree->Branch("Events","AliEvent",&fEvent,bsize,split);
+}
+///////////////////////////////////////////////////////////////////////////
+void AliCollider::SetVertexMode(Int_t mode)
+{
+// Set the mode of the vertex structure creation.
+//
+// By default all generated tracks will only appear in the AliEvent
+// structure without any primary (and secondary) vertex structure.
+// The user can build the vertex structure if he/she wants by means
+// of the beginpoint location of each AliTrack.
+//
+// However, one can also let AliCollider automatically create
+// the primary (and secondary) vertex structure(s).
+// In this case the primary vertex is given Id=1 and all sec. vertices
+// are given Id's 2,3,4,....
+// All vertices are created as standalone entities in the AliEvent structure
+// without any linking between the various vertices.
+// For this automated process, the user-selected resolution
+// (see SetResolution) is used to decide whether or not certain vertex
+// locations can be resolved.
+// In case no vertex creation is selected (i.e. the default mode=0),
+// the value of the resolution is totally irrelevant.
+//
+// The user can also let AliCollider automatically connect the sec. vertices
+// to the primary vertex (i.e. mode=3). This process will also automatically
+// generate the tracks connecting the vertices.
+// Note that the result of the mode=3 operation may be very sensitive to
+// the resolution parameter. Therefore, no attempt is made to distinguish
+// between secondary, tertiary etc... vertices. All sec. vertices are
+// linked to the primary one.
+//  
+// Irrespective of the selected mode, all generated tracks can be obtained
+// directly from the AliEvent structure.
+// In case (sec.) vertex creation is selected, all generated vertices can
+// also be obtained directly from the AliEvent structure. 
+// These (sec.) vertices contain only the corresponding pointers to the various
+// tracks which are stored in the AliEvent structure.
+//
+// Overview of vertex creation modes :
+// -----------------------------------
+// mode = 0 ==> No vertex structure will be created
+//        1 ==> Only primary vertex structure will be created
+//        2 ==> Unconnected primary and secondary vertices will be created
+//        3 ==> Primary and secondary vertices will be created where all the
+//              sec. vertices will be connected to the primary vertex.
+//              Also the vertex connecting tracks will be automatically
+//              generated. 
+//
+ if (mode<0 || mode >3)
+ {
+  cout << " *AliCollider::SetVertexMode* Invalid argument mode : " << mode << endl;
+  fVertexmode=0;
+ }
+ else
+ {
+  fVertexmode=mode;
+ }
+}
+///////////////////////////////////////////////////////////////////////////
+Int_t AliCollider::GetVertexMode()
+{
+// Provide the current mode for vertex structure creation.
+ return fVertexmode;
+}
+///////////////////////////////////////////////////////////////////////////
+void AliCollider::SetResolution(Double_t res)
+{
+// Set the resolution (in cm) for resolving (sec.) vertices.
+// By default this resolution is set to 0.1 micron.
+// Note : In case no vertex creation has been selected, the value of
+//        the resolution is totally irrelevant.
+ fResolution=fabs(res);
+}
+///////////////////////////////////////////////////////////////////////////
+Double_t AliCollider::GetResolution()
+{
+// Provide the current resolution (in cm) for resolving (sec.) vertices.
+ return fResolution;
+}
+///////////////////////////////////////////////////////////////////////////
+void AliCollider::SetRunNumber(Int_t run)
+{
+// Set the user defined run number.
+// By default the run number is set to 0.
+ fRunnum=run;
+}
+///////////////////////////////////////////////////////////////////////////
+Int_t AliCollider::GetRunNumber()
+{
+// Provide the user defined run number.
+ return fRunnum;
+}
+///////////////////////////////////////////////////////////////////////////
+void AliCollider::SetPrintFreq(Int_t n)
+{
+// Set the print frequency for every 'n' events.
+// By default the printfrequency is set to 1 (i.e. every event).
+ fPrintfreq=n;
+}
+///////////////////////////////////////////////////////////////////////////
+Int_t AliCollider::GetPrintFreq()
+{
+// Provide the user selected print frequency.
+ return fPrintfreq;
+}
+///////////////////////////////////////////////////////////////////////////
+void AliCollider::Init(char* frame,char* beam,char* target,Float_t win)
+{
+// Initialisation of the underlying Pythia generator package.
+// This routine just invokes TPythia6::Initialize(...) and the arguments
+// have the corresponding meaning.
+// The event number is reset to 0.
+ fEventnum=0;
+ fNucl=0;
+ fFrame=frame;
+ fWin=win;
+ Initialize(frame,beam,target,win);
+}
+///////////////////////////////////////////////////////////////////////////
+void AliCollider::Init(char* frame,Int_t zp,Int_t ap,Int_t zt,Int_t at,Float_t win)
+{
+// Initialisation of the underlying Pythia generator package for the generation
+// of nucleus-nucleus interactions.
+// In addition to the Pythia standard arguments 'frame' and 'win', the user
+// can specify here (Z,A) values of the projectile and target nuclei and the number
+// 'npart' of the participant nucleons for this collision.
+// The event number is reset to 0.
+ fEventnum=0;
+ fNucl=1;
+ fFrame=frame;
+ fWin=win;
+ fZproj=0;
+ fAproj=0;
+ fZtarg=0;
+ fAtarg=0;
+ fFracpp=0;
+ fFracnp=0;
+ fFracpn=0;
+ fFracnn=0;
+
+ if (ap<1 || at<1 || zp>ap || zt>at)
+ {
+  cout << " *AliCollider::Init* Invalid input value(s). Zproj = " << zp
+       << " Aproj = " << ap << " Ztarg = " << zt << " Atarg = " << at << endl;
+  return;
+ }
+
+ fZproj=zp;
+ fAproj=ap;
+ fZtarg=zt;
+ fAtarg=at;
+
+ cout << " *AliCollider::Init* Nucleus-Nucleus generator initialisation." << endl;
+ cout << " Zproj = " << zp << " Aproj = " << ap << " Ztarg = " << zt << " Atarg = " << at
+      << " Frame = " << frame << " Energy = " << win
+      << endl;
+}
+///////////////////////////////////////////////////////////////////////////
+void AliCollider::GetFractions(Float_t zp,Float_t ap,Float_t zt,Float_t at)
+{
+// Determine the fractions for the various N-N collision processes.
+// The various processes are : p+p, n+p, p+n and n+n.
+ if (zp<0) zp=0;
+ if (zt<0) zt=0;
+
+ fFracpp=0;
+ fFracnp=0;
+ fFracpn=0;
+ fFracnn=0;
+
+ if (ap>0 && at>0)
+ {
+  fFracpp=(zp/ap)*(zt/at);
+  fFracnp=(1.-zp/ap)*(zt/at);
+  fFracpn=(zp/ap)*(1.-zt/at);
+  fFracnn=(1.-zp/ap)*(1.-zt/at);
+ }
+}
+///////////////////////////////////////////////////////////////////////////
+void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit)
+{
+// Generate one event.
+// In case of a nucleus-nucleus interaction, the argument 'npt' denotes
+// the number of participant nucleons.
+// In case of a standard Pythia run for 'elementary' particle interactions,
+// the value of npt is totally irrelevant.
+// The argument 'medit' denotes the edit mode used for Pyedit().
+// By default, only 'stable' final particles are kept (i.e. medit=1). 
+// The argument 'mlist' denotes the list mode used for Pylist().
+// By default, no listing is produced (i.e. mlist=0).
+
+ fEventnum++; 
+
+ // Counters for the various (proj,targ) combinations : p+p, n+p, p+n and n+n
+ Int_t ncols[4]={0,0,0,0};
+
+ if (fNucl)
+ {
+  if (npt<1 || npt>(fAproj+fAtarg))
+  {
+   cout << " *AliCollider::MakeEvent* Invalid input value. npt = " << npt
+        << " Aproj = " << fAproj << " Atarg = " << fAtarg << endl;
+   return;
+  }
+
+  // Determine the number of nucleon-nucleon collisions
+  Int_t ncol=npt/2.;
+  if (npt%2 && fRan.Uniform()>0.5) ncol+=1;
+
+  // Determine the number of the various types of N+N interactions
+  Int_t zp=fZproj;
+  Int_t ap=fAproj;
+  Int_t zt=fZtarg;
+  Int_t at=fAtarg;
+  Int_t maxa=2; // Indicator whether proj (1) or target (2) has maximal A left
+  if (ap>at) maxa=1;
+  Float_t* rans=new Float_t[ncol];
+  fRan.Uniform(rans,ncol);
+  Float_t rndm=0;
+  for (Int_t i=0; i<ncol; i++)
+  {
+   GetFractions(zp,ap,zt,at);
+   rndm=rans[i];
+   if (rndm<=fFracpp) // p+p interaction
+   {
+    ncols[0]++;
+    if (maxa=2)
+    {
+     at--;
+     zt--;
+    } 
+    else
+    {
+     ap--;
+     zp--;
+    }
+   }
+   if (rndm>fFracpp && rndm<=(fFracpp+fFracnp)) // n+p interaction
+   {
+    ncols[1]++;
+    if (maxa=2)
+    {
+     at--;
+     zt--;
+    } 
+    else
+    {
+     ap--;
+    }
+   }
+   if (rndm>(fFracpp+fFracnp) && rndm<=(fFracpp+fFracnp+fFracpn)) // p+n interaction
+   {
+    ncols[2]++;
+    if (maxa=2)
+    {
+     at--;
+    } 
+    else
+    {
+     ap--;
+     zp--;
+    }
+   }
+   if (rndm>(fFracpp+fFracnp+fFracpn)) // n+n interaction
+   {
+    ncols[3]++; 
+    if (maxa=2)
+    {
+     at--;
+    } 
+    else
+    {
+     ap--;
+    }
+   }
+  }
+  delete [] rans;
+
+  if (!(fEventnum%fPrintfreq))
+  {
+   cout << " *AliCollider::MakeEvent* Run : " << fRunnum << " Event : " << fEventnum
+        << endl;
+   cout << " npart = " << npt << " ncol = " << ncol 
+        << " ncolpp = " << ncols[0] << " ncolnp = " << ncols[1]
+        << " ncolpn = " << ncols[2] << " ncolnn = " << ncols[3] << endl;
+  }
+
+ }
+
+ if (!fEvent)
+ {
+  fEvent=new AliEvent();
+  fEvent->SetOwner();
+ }
+
+ fEvent->Reset();
+ fEvent->SetRunNumber(fRunnum);
+ fEvent->SetEventNumber(fEventnum);
+
+ AliVertex vert;
+ if (fVertexmode)
+ {
+  // Make sure the primary vertex gets correct location and Id=1
+  vert.SetId(1);
+  vert.SetTrackCopy(0);
+  vert.SetVertexCopy(0);
+  fEvent->AddVertex(vert,0);
+ }
+
+ AliTrack t;
+ Ali3Vector p;
+ AliPosition r,rx;
+ Float_t v[3];
+
+ Int_t kf=0,kc=0;
+ Float_t charge=0,mass=0;
+
+ TMCParticle* part=0;
+
+ Int_t ntypes=4;
+
+ // Singular settings for a normal Pythia elementary particle interation 
+ if (!fNucl)
+ {
+  ntypes=1;
+  ncols[0]=1;
+ }
+
+ // Generate all the various collisions
+ Int_t npart=0,ntk=0;
+ Double_t dist=0;
+ for (Int_t itype=0; itype<ntypes; itype++)
+ {
+  if (fNucl)
+  {
+   if (itype==0 && ncols[itype]) Initialize(fFrame,"p","p",fWin);
+   if (itype==1 && ncols[itype]) Initialize(fFrame,"n","p",fWin);
+   if (itype==2 && ncols[itype]) Initialize(fFrame,"p","n",fWin);
+   if (itype==3 && ncols[itype]) Initialize(fFrame,"n","n",fWin);
+  }
+  for (Int_t jcol=0; jcol<ncols[itype]; jcol++)
+  {
+   GenerateEvent();
+
+   Pyedit(medit); // Define which particles are to be kept
+
+   if (mlist) Pylist(mlist);
+
+   ImportParticles();
+   npart=0;
+   if (fParticles) npart=fParticles->GetEntries();
+
+   for (Int_t jpart=0; jpart<npart; jpart++)
+   {
+    part=(TMCParticle*)fParticles->At(jpart);
+    if (!part) continue;
+
+    kf=part->GetKF();
+    kc=Pycomp(kf);
+
+    charge=GetKCHG(kc,1)/3.;
+    if (kf<0) charge*=-1;
+    mass=GetPMAS(kc,1);
+
+    // 3-momentum in GeV/c
+    v[0]=part->GetPx();
+    v[1]=part->GetPy();
+    v[2]=part->GetPz();
+    p.SetVector(v,"car");
+
+    // Production location in cm.
+    v[0]=(part->GetVx())/10;
+    v[1]=(part->GetVy())/10;
+    v[2]=(part->GetVz())/10;
+    r.SetVector(v,"car");
+
+    ntk++;
+
+    t.Reset();
+    t.SetId(ntk);
+    t.SetParticleCode(kf);
+    t.SetCharge(charge);
+    t.SetMass(mass);
+    t.Set3Momentum(p);
+    t.SetBeginPoint(r);
+
+    fEvent->AddTrack(t);
+
+    // Build the vertex structures if requested
+    if (fVertexmode)
+    {
+     // Check if track belongs within the resolution to an existing vertex
+     Int_t add=0;  
+     for (Int_t jv=1; jv<=fEvent->GetNvertices(); jv++)
+     {
+      AliVertex* vx=fEvent->GetVertex(jv);
+      if (vx)
+      {
+       rx=vx->GetPosition();
+       dist=rx.GetDistance(r);
+       if (dist < fResolution)
+       {
+        AliTrack* tx=fEvent->GetIdTrack(ntk);
+        if (tx)
+        { 
+         vx->AddTrack(tx);
+         add=1;
+        }
+       }
+      }
+      if (add) break; // No need to look further for vertex candidates
+     }
+
+     // If track was not close enough to an existing vertex
+     // a new secondary vertex is created      
+     if (!add && fVertexmode>1)
+     {
+      AliTrack* tx=fEvent->GetIdTrack(ntk);
+      if (tx)
+      {
+       vert.Reset();
+       vert.SetTrackCopy(0);
+       vert.SetVertexCopy(0);
+       vert.SetId((fEvent->GetNvertices())+1);
+       vert.SetPosition(r);
+       vert.AddTrack(tx);
+       fEvent->AddVertex(vert,0);
+      } 
+     }
+    }
+   } // End of loop over the produced particles for each collision
+  } // End of loop over number of collisions for each type
+ } // End of loop over collision types
+
+ // Link sec. vertices to primary if requested
+ // Note that also the connecting tracks are automatically created
+ if (fVertexmode>2)
+ {
+  AliVertex* vp=fEvent->GetIdVertex(1); // Primary vertex
+  if (vp)
+  {
+   for (Int_t i=2; i<=fEvent->GetNvertices(); i++)
+   {
+    AliVertex* vx=fEvent->GetVertex(i);
+    if (vx)
+    {
+     if (vx->GetId() != 1) vp->AddVertex(vx);
+    }
+   }
+  }
+ }
+
+ if (mlist) cout << endl; // Create empty output line after the event
+ if (fOutTree) fOutTree->Fill();
+}
+///////////////////////////////////////////////////////////////////////////
+AliEvent* AliCollider::GetEvent()
+{
+// Provide pointer to the generated event structure.
+ return fEvent;
+}
+///////////////////////////////////////////////////////////////////////////
+void AliCollider::EndRun()
+{
+// Properly close the output file (if needed).
+ if (fOutFile)
+ {
+  fOutFile->Write();
+  fOutFile->Close();
+  cout << " *AliCollider::EndRun* Output file correctly closed." << endl;
+ }
+}
+///////////////////////////////////////////////////////////////////////////
diff --git a/RALICE/AliCollider.h b/RALICE/AliCollider.h
new file mode 100644 (file)
index 0000000..112ff57
--- /dev/null
@@ -0,0 +1,63 @@
+#ifndef ALICOLLIDER_H
+#define ALICOLLIDER_H
+
+// $Id$
+
+#include "Riostream.h"
+#include "TPythia6.h"
+#include "TMCParticle.h"
+#include "TString.h"
+#include "TFile.h"
+#include "TTree.h" 
+
+#include "AliEvent.h"
+#include "AliRandom.h"
+class AliCollider : public TPythia6
+{
+ public:
+  AliCollider();                                        // Default constructor
+  ~AliCollider();                                       // Default destructor
+  void SetOutputFile(TString name);                     // Initialise the ROOT output data file
+  void SetVertexMode(Int_t mode);                       // Select mode for (sec.) vertex structure creation
+  Int_t GetVertexMode();                                // Provide vertex structure creation mode
+  void SetResolution(Double_t res);                     // Set resolution (in cm) for resolving (sec.) vertices
+  Double_t GetResolution();                             // Provide (sec.) vertex resolving resolution (in cm)
+  void SetRunNumber(Int_t run);                         // Set user defined run number
+  Int_t GetRunNumber();                                 // Provide the user defined run number
+  void SetPrintFreq(Int_t n);                           // Set print frequency for every 'n' events
+  Int_t GetPrintFreq();                                 // Provide the print frequency
+  void Init(char* frame,char* beam,char* target,Float_t win); // Standard Pythia initialisation
+  void Init(char* frame,Int_t zp,Int_t ap,Int_t zt,Int_t at,Float_t win); // Nucl-Nucl initialisation
+  void MakeEvent(Int_t npt,Int_t mlist=0,Int_t medit=1);// Generate a single event with npt participant nucleons
+  void EndRun();                                        // Properly close all buffers and output file
+  AliEvent* GetEvent();                                 // Provide pointer to the generated event structure
+
+ protected:
+  Int_t fVertexmode;    // The vertex structure creation mode
+  Double_t fResolution; // The resolution (in cm) for resolving (sec.) vertices 
+  Int_t fRunnum;        // The user defined run number
+  Int_t fEventnum;      // The automatically updated event number
+  Int_t fPrintfreq;     // The user selected print frequency
+  char* fFrame;         // The Pythia frame indicator
+  Float_t fWin;         // The Pythia energy indicator
+  Int_t fNucl;          // Flag to denote nucleus-nucleus (1) or standard Pythia (0) running
+  Int_t fZproj;         // Z of the projectile particle
+  Int_t fAproj;         // A of the projectile particle
+  Int_t fZtarg;         // Z of the target particle
+  Int_t fAtarg;         // A of the target particle
+  Float_t fFracpp;      // Fraction of p+p collisions
+  Float_t fFracnp;      // Fraction of n+p collisions
+  Float_t fFracpn;      // Fraction of p+n collisions
+  Float_t fFracnn;      // Fraction of n+n collisions
+  AliRandom fRan;       // Random number generator
+  AliEvent* fEvent;     // The produced event structure
+
+  TFile* fOutFile;      // The user defined output data file 
+  TTree* fOutTree;      // The standard ROOT output tree
+
+  void GetFractions(Float_t zp,Float_t ap,Float_t zt,Float_t at); // Determine various N-N collision fractions
+
+ ClassDef(AliCollider,1) // Pythia based universal physics event generator
+};
+#endif
index 4a64ec553fda46e23c0ca3b08589991f7a263eeb..2739a922ac9fa5e46fa48a546afe92ef22f3e27d 100644 (file)
@@ -136,6 +136,7 @@ AliTrack::AliTrack(AliTrack& t)
  fChi2=t.GetChi2();
  fNdf=t.GetNdf();
  fUserId=t.GetId();
+ fCode=t.GetParticleCode();
  fNdec=t.GetNdecay();
  fNsig=t.GetNsignals();
  fNmasses=t.GetNMassHypotheses();
@@ -195,6 +196,7 @@ void AliTrack::Reset()
  fChi2=0;
  fNdf=0;
  fUserId=0;
+ fCode=0;
  fNdec=0;
  fNsig=0;
  fNmasses=0;
@@ -271,8 +273,8 @@ void AliTrack::Info(TString f)
 // Provide track information within the coordinate frame f
  Double_t m=GetMass();
  Double_t dm=GetResultError();
- cout << " *AliTrack::Info* Id : " << fUserId << " Mass : " << m
-      << " error : " << dm << " Charge : " << fQ
+ cout << " *AliTrack::Info* Id : " << fUserId << " Code : " << fCode
+      << " Mass : " << m << " error : " << dm << " Charge : " << fQ
       << " Momentum : " << GetMomentum() << " Nmass hyp. : " << fNmasses
       << " Ntracks : " << fNdec << " Nsignals : " << fNsig << endl;
  for (Int_t i=0; i<fNmasses; i++)
@@ -932,13 +934,13 @@ AliPosition AliTrack::GetImpactPoint(TString q)
 ///////////////////////////////////////////////////////////////////////////
 void AliTrack::SetId(Int_t id)
 {
-// Set a user defined identifier for this track.
+// Set a user defined unique identifier for this track.
  fUserId=id;
 }
 ///////////////////////////////////////////////////////////////////////////
 Int_t AliTrack::GetId()
 {
-// Provide the user defined identifier of this track.
+// Provide the user defined unique identifier of this track.
  return fUserId;
 }
 ///////////////////////////////////////////////////////////////////////////
@@ -992,3 +994,15 @@ Int_t AliTrack::GetNdf()
  return fNdf;
 }
 ///////////////////////////////////////////////////////////////////////////
+void AliTrack::SetParticleCode(Int_t code)
+{
+// Set the user defined particle id code (e.g. the PDF convention).
+ fCode=code;
+}
+///////////////////////////////////////////////////////////////////////////
+Int_t AliTrack::GetParticleCode()
+{
+// Provide the user defined particle id code.
+ return fCode;
+}
+///////////////////////////////////////////////////////////////////////////
index 09c8445b31fee09a7adbd01d44847ae5fd3428c0..16d18e73b5c515d1adbe2c3f2914a17888832cb4 100644 (file)
@@ -60,14 +60,16 @@ class AliTrack : public TObject,public Ali4Vector
   Double_t GetRapidity();           // Provide rapidity value w.r.t. z-axis
   void SetImpactPoint(AliPosition p,TString q); // Set the impact-point in plane "q=0"
   AliPosition GetImpactPoint(TString q);        // Provide the impact-point in plane "q=0"
-  void SetId(Int_t id);             // Set the user defined identifier
-  Int_t GetId();                    // Provide the user defined identifier
+  void SetId(Int_t id);             // Set the user defined unique track identifier
+  Int_t GetId();                    // Provide the user defined unique track identifier
   void SetClosestPoint(AliPosition p); // Set position p as point of closest approach w.r.t. some reference
   AliPosition GetClosestPoint();       // Provide point of closest approach w.r.t. some reference
   void SetChi2(Float_t chi2);       // Set the chi-squared value of the track fit
   void SetNdf(Int_t ndf);           // Set the number of degrees of freedom for the track fit
   Float_t GetChi2();                // Provide the chi-squared value of the track fit
   Int_t GetNdf();                   // Provide the number of degrees of freedom for the track fit
+  void SetParticleCode(Int_t code); // Set the user defined particle id code (e.g. the PDF convention)
+  Int_t GetParticleCode();          // Provide the user defined particle id code
 
  
  protected:
@@ -90,10 +92,11 @@ class AliTrack : public TObject,public Ali4Vector
   AliPosition fClosest;  // The (extrapolated) point of closest approach w.r.t some reference
   Float_t fChi2;         // The Chi-squared of the track fit
   Int_t fNdf;            // The number of degrees of freedom of the track fit
+  Int_t fCode;           // The user defined particle id code
 
  private:
   void Dump(AliTrack* t,Int_t n,TString f); // Recursively print all decay levels
  
- ClassDef(AliTrack,1) // Handling of the attributes of a reconstructed particle track.
+ ClassDef(AliTrack,2) // Handling of the attributes of a reconstructed particle track.
 };
 #endif
index 5384aa675f31c9b9343ee8c4367bbffbb1ce84a6..f7a9218910653e5c0732cbd82cdbc3c865d4a762 100644 (file)
@@ -28,3 +28,4 @@
 #include "AliVertex.h"
 #include "AliInvmass.h"
 #include "AliEvent.h"
+#include "AliCollider.h"
index f00c8ee124a0dbfcf17c1b10cfc11e730ed46abf..5f223c0bf08ec5373d4ec66b965de6d17fd200c2 100644 (file)
@@ -33,5 +33,6 @@
  #pragma link C++ class AliVertex+;
  #pragma link C++ class AliInvmass+;
  #pragma link C++ class AliEvent+;
+ #pragma link C++ class AliCollider+;
 #endif
  
index 08f7277351a9ded10f77e48c20f2582aee9815d4..6d94b9bcd7ca8066fd37cd887d71c2c6c5885a5f 100644 (file)
 28-oct-2002 NvE "Riostream.h" introduced to replace the standard C++ includes and 
                 static_cast<int> introduced in AliVertex::Draw to prevent a warning
                 message for the new gcc compiler version.
+25-nov-2002 NvE User defined particle id code introduced in AliTrack.
+27-nov-2002 NvE New class AliCollider introduced and RALICEHeaders.h and RALICELinkdef.h
+                updated accordingly.