22-dec-2003 NvE Event selection introduced in AliCollider.
authornick <nick@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 12 Jan 2004 08:23:22 +0000 (08:23 +0000)
committernick <nick@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 12 Jan 2004 08:23:22 +0000 (08:23 +0000)
24-dec-2003 NvE Some cosmetics in printout of AliVertex, AliEvent and AliCollider.
                Support for spectator tracks introduced in AliCollider.
28-dec-2003 NvE Facility introduced in AliCollider to set minimal momentum for spectator tracks
                to be stored.

RALICE/AliCollider.cxx
RALICE/AliCollider.h
RALICE/AliEvent.cxx
RALICE/AliVertex.cxx
RALICE/history.txt

index 4f1bcea2fc073411ba49f21eb767d961a9044177..7d570706edc176cd1ebc5345b261ef6706d6dd4a 100644 (file)
@@ -13,7 +13,7 @@
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
 
-// $Id: AliCollider.cxx,v 1.7 2003/10/26 14:53:44 nick Exp $
+// $Id: AliCollider.cxx,v 1.8 2003/12/18 09:28:06 nick Exp $
 
 ///////////////////////////////////////////////////////////////////////////
 // Class AliCollider
@@ -53,6 +53,8 @@
 //
 //  gen->Init("fixt",zp,ap,zt,at,158);
 //
+//  gen->SetTitle("SPS Pb-Pb collision at 158A GeV/c beam energy");
+//
 //  Int_t nevents=5;
 //
 //  AliRandom rndm;
@@ -88,6 +90,8 @@
 //
 //  gen->Init("fixt","nu_mu","p",1e11);
 //
+//  gen->SetTitle("Atmospheric nu_mu-p interaction at 1e20 eV");
+//
 //  Int_t nevents=10;
 //
 //  for (Int_t i=0; i<nevents; i++)
 //
 //
 //--- Author: Nick van Eijndhoven 22-nov-2002 Utrecht University
-//- Modified: NvE $Date: 2003/10/26 14:53:44 $ Utrecht University
+//- Modified: NvE $Date: 2003/12/18 09:28:06 $ Utrecht University
 ///////////////////////////////////////////////////////////////////////////
 
 #include "AliCollider.h"
@@ -124,6 +128,8 @@ AliCollider::AliCollider() : TPythia6()
 
  fEvent=0;
 
+ fSpecpmin=0;
+
  fFrame="none";
  fWin=0;
 
@@ -139,6 +145,13 @@ AliCollider::AliCollider() : TPythia6()
 
  fOutFile=0;
  fOutTree=0;
+
+ fSelections=0;
+ fSelect=0;
+
+ TString s=GetName();
+ s+=" (AliCollider)";
+ SetName(s.Data());
 }
 ///////////////////////////////////////////////////////////////////////////
 AliCollider::~AliCollider()
@@ -159,6 +172,11 @@ AliCollider::~AliCollider()
   delete fOutTree;
   fOutTree=0;
  }
+ if (fSelections)
+ {
+  delete fSelections;
+  fSelections=0;
+ }
 }
 ///////////////////////////////////////////////////////////////////////////
 void AliCollider::SetOutputFile(TString s)
@@ -373,6 +391,13 @@ 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.
+// Normally also the spectator tracks will be stored into the event structure.
+// The spectator tracks have a negative user Id to distinguish them from the
+// ordinary generated tracks.
+// In case the user has selected the creation of vertex structures, the spectator
+// tracks will be linked to the primary vertex.
+// However, specification of npt<0 will suppress the storage of spectator tracks.
+// In the latter case abs(npt) will be taken as the number of participants.  
 // In case of a standard Pythia run for 'elementary' particle interactions,
 // the value of npt is totally irrelevant.
 //
@@ -391,12 +416,27 @@ void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit)
 // In case of a nucleus-nucleus interaction, the proper A and Z values for 
 // the projectile and target particles are set in the event structure.
 // However, in this case both particle ID's are set to zero.
+//
+// Note : Only in case an event passed the selection criteria as specified
+//        via SelectEvent(), the event will appear on the output file.
 
  fEventnum++; 
 
+ Int_t specmode=1;
+ if (npt<0)
+ {
+  specmode=0;
+  npt=abs(npt);
+ }
+
  // Counters for the various (proj,targ) combinations : p+p, n+p, p+n and n+n
  Int_t ncols[4]={0,0,0,0};
 
+ Int_t zp=0;
+ Int_t ap=0;
+ Int_t zt=0;
+ Int_t at=0;
+
  Int_t ncol=1;
  if (fNucl)
  {
@@ -412,10 +452,10 @@ void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit)
   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;
+  zp=fZproj;
+  ap=fAproj;
+  zt=fZtarg;
+  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];
@@ -497,6 +537,8 @@ void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit)
  {
   fEvent=new AliEvent();
   fEvent->SetOwner();
+  fEvent->SetName(GetName());
+  fEvent->SetTitle(GetTitle());
  }
 
  fEvent->Reset();
@@ -508,6 +550,7 @@ void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit)
  AliPosition r,rx;
  Float_t v[3];
  AliVertex vert;
+ Ali3Vector pproj,ptarg;
 
  if (fVertexmode)
  {
@@ -530,6 +573,7 @@ void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit)
 
  Int_t kf=0;
  Float_t charge=0,mass=0;
+ char* name="";
 
  Int_t ntypes=4;
 
@@ -541,7 +585,9 @@ void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit)
  }
 
  // Generate all the various collisions
- Int_t first=1; // Flag to indicate the first collision process
+ fSelect=0;      // Flag to indicate whether the total event is selected or not 
+ Int_t select=0; // Flag to indicate whether the sub-event is selected or not 
+ Int_t first=1;  // Flag to indicate the first collision process
  Double_t pnucl;
  Int_t npart=0,ntk=0;
  Double_t dist=0;
@@ -558,6 +604,9 @@ void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit)
   {
    GenerateEvent();
 
+   select=IsSelected();
+   if (select) fSelect=1;
+
    if (first) // Store projectile and target information in the event structure
    {
     if (fNucl)
@@ -565,12 +614,14 @@ void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit)
      v[0]=GetP(1,1);
      v[1]=GetP(1,2);
      v[2]=GetP(1,3);
-     pnucl=sqrt(v[0]*v[0]+v[1]*v[1]+v[2]*v[2]);
+     pproj.SetVector(v,"car");
+     pnucl=pproj.GetNorm();
      fEvent->SetProjectile(fAproj,fZproj,pnucl);
      v[0]=GetP(2,1);
      v[1]=GetP(2,2);
      v[2]=GetP(2,3);
-     pnucl=sqrt(v[0]*v[0]+v[1]*v[1]+v[2]*v[2]);
+     ptarg.SetVector(v,"car");
+     pnucl=ptarg.GetNorm();
      fEvent->SetTarget(fAtarg,fZtarg,pnucl);
     }
     else
@@ -593,13 +644,14 @@ void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit)
 
    if (medit >= 0) Pyedit(medit); // Define which particles are to be kept
 
-   if (mlist >= 0) Pylist(mlist);
+   if (mlist>=0 && select) Pylist(mlist);
 
    npart=GetN();
    for (Int_t jpart=1; jpart<=npart; jpart++)
    {
     kf=GetK(jpart,2);
     charge=Pychge(kf)/3.;
+    Pyname(kf,name);
     mass=GetP(jpart,5);
 
     // 3-momentum in GeV/c
@@ -619,6 +671,7 @@ void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit)
     t.Reset();
     t.SetId(ntk);
     t.SetParticleCode(kf);
+    t.SetName(name);
     t.SetCharge(charge);
     t.SetMass(mass);
     t.Set3Momentum(p);
@@ -694,14 +747,152 @@ void AliCollider::MakeEvent(Int_t npt,Int_t mlist,Int_t medit)
   }
  }
 
+ // Include the spectator tracks in the event structure.
+ if (fNucl && specmode)
+ {
+  Float_t pmass=0.938272;
+  Float_t nmass=0.93956533;
+  v[0]=0;
+  v[1]=0;
+  v[2]=0;
+  r.SetPosition(v,"car");
+
+  zp=fZproj-(ncols[0]+ncols[2]);
+  if (zp<0) zp=0;
+  ap=fAproj-(ncols[0]+ncols[1]+ncols[2]+ncols[3]);
+  if (ap<0) ap=0;
+  zt=fZtarg-(ncols[0]+ncols[1]);
+  if (zt<0) zt=0;
+  at=fAtarg-(ncols[0]+ncols[1]+ncols[2]+ncols[3]);
+  if (at<0) at=0;
+
+  Int_t nspec=0;
+
+  if (pproj.GetNorm() > fSpecpmin)
+  {
+   kf=2212; // Projectile spectator protons
+   charge=Pychge(kf)/3.;
+   mass=pmass;
+   Pyname(kf,name);
+   for (Int_t iprojp=1; iprojp<=zp; iprojp++)
+   {
+    nspec++;
+    t.Reset();
+    t.SetId(-nspec);
+    t.SetParticleCode(kf);
+    t.SetName(name);
+    t.SetTitle("Projectile spectator proton");
+    t.SetCharge(charge);
+    t.SetMass(mass);
+    t.Set3Momentum(pproj);
+    t.SetBeginPoint(r);
+
+    fEvent->AddTrack(t);
+   }
+
+   kf=2112; // Projectile spectator neutrons
+   charge=Pychge(kf)/3.;
+   mass=nmass;
+   Pyname(kf,name);
+   for (Int_t iprojn=1; iprojn<=(ap-zp); iprojn++)
+   {
+    nspec++;
+    t.Reset();
+    t.SetId(-nspec);
+    t.SetParticleCode(kf);
+    t.SetName(name);
+    t.SetTitle("Projectile spectator neutron");
+    t.SetCharge(charge);
+    t.SetMass(mass);
+    t.Set3Momentum(pproj);
+    t.SetBeginPoint(r);
+
+    fEvent->AddTrack(t);
+   }
+  }
+
+  if (ptarg.GetNorm() > fSpecpmin)
+  {
+   kf=2212; // Target spectator protons
+   charge=Pychge(kf)/3.;
+   mass=pmass;
+   Pyname(kf,name);
+   for (Int_t itargp=1; itargp<=zt; itargp++)
+   {
+    nspec++;
+    t.Reset();
+    t.SetId(-nspec);
+    t.SetParticleCode(kf);
+    t.SetName(name);
+    t.SetTitle("Target spectator proton");
+    t.SetCharge(charge);
+    t.SetMass(mass);
+    t.Set3Momentum(ptarg);
+    t.SetBeginPoint(r);
+
+    fEvent->AddTrack(t);
+   }
+
+   kf=2112; // Target spectator neutrons
+   charge=Pychge(kf)/3.;
+   mass=nmass;
+   Pyname(kf,name);
+   for (Int_t itargn=1; itargn<=(at-zt); itargn++)
+   {
+    nspec++;
+    t.Reset();
+    t.SetId(-nspec);
+    t.SetParticleCode(kf);
+    t.SetName(name);
+    t.SetTitle("Target spectator neutron");
+    t.SetCharge(charge);
+    t.SetMass(mass);
+    t.Set3Momentum(ptarg);
+    t.SetBeginPoint(r);
+
+    fEvent->AddTrack(t);
+   }
+  }
+
+ // Link the spectator tracks to the primary vertex.
+ if (fVertexmode)
+ {
+  AliVertex* vp=fEvent->GetIdVertex(1);
+  if (vp)
+  {
+   for (Int_t ispec=1; ispec<=nspec; ispec++)
+   {
+    AliTrack* tx=fEvent->GetIdTrack(-ispec);
+    if (tx) vp->AddTrack(tx);
+   }
+  }
+ }
+}
+
  if (mlist && !(fEventnum%fPrintfreq)) cout << endl; // Create empty output line after the event
- if (fOutTree) fOutTree->Fill();
+
+ if (fOutTree && fSelect) fOutTree->Fill();
 }
 ///////////////////////////////////////////////////////////////////////////
-AliEvent* AliCollider::GetEvent()
+AliEvent* AliCollider::GetEvent(Int_t select)
 {
 // Provide pointer to the generated event structure.
- return fEvent;
+//
+// select = 0 : Always return the pointer to the generated event.
+//          1 : Only return the pointer to the generated event in case
+//              the event passed the selection criteria as specified via
+//              SelectEvent(). Otherwise the value 0 will be returned.
+//
+// By invoking GetEvent() the default of select=0 will be used.
+
+ if (!select || fSelect)
+ {
+  return fEvent;
+ }
+ else
+ {
+  return 0;
+ }
 }
 ///////////////////////////////////////////////////////////////////////////
 void AliCollider::EndRun()
@@ -759,3 +950,121 @@ void AliCollider::SetStable(Int_t id,Int_t mode)
  }
 }
 ///////////////////////////////////////////////////////////////////////////
+void AliCollider::SelectEvent(Int_t id)
+{
+// Add a particle to the event selection list.
+// The parameter "id" indicates the Pythia KF particle code, which
+// basically is the PDG particle identifier code.
+// In case the user has built a selection list via this procedure, only the
+// events in which one of the particles specified in the list was generated
+// will be kept. 
+// The investigation of the generated particles takes place when the complete
+// event is in memory, including all (shortlived) mother particles and resonances.
+// So, the settings of the various particle decay modes have no influence on
+// the event selection described here.
+//
+// If no list has been specified, all events will be accepted.  
+//
+// Note : id=0 will delete the selection list.
+//
+// Be aware of the fact that severe selection criteria (i.e. selecting only
+// rare events) may result in long runtimes before an event sample has been
+// obtained.
+//
+ if (!id)
+ {
+  if (fSelections)
+  {
+   delete fSelections;
+   fSelections=0;
+  }
+ }
+ else
+ {
+  Int_t kc=Pycomp(id);
+  if (!fSelections)
+  {
+   fSelections=new TArrayI(1);
+   fSelections->AddAt(kc,0);
+  }
+  else
+  {
+   Int_t exist=0;
+   Int_t size=fSelections->GetSize();
+   for (Int_t i=0; i<size; i++)
+   {
+    if (kc==fSelections->At(i))
+    {
+     exist=1;
+     break;
+    }
+   }
+  
+   if (!exist)
+   {
+    fSelections->Set(size+1);
+    fSelections->AddAt(kc,size);
+   }
+  }
+ }
+}
+///////////////////////////////////////////////////////////////////////////
+Int_t AliCollider::GetSelectionFlag()
+{
+// Return the value of the selection flag for the total event.
+// When the event passed the selection criteria as specified via
+// SelectEvent() the value 1 is returned, otherwise the value 0 is returned.
+ return fSelect;
+}
+///////////////////////////////////////////////////////////////////////////
+Int_t AliCollider::IsSelected()
+{
+// Check whether the generated (sub)event contains one of the particles
+// specified in the selection list via SelectEvent().
+// If this is the case or when no selection list is present, the value 1
+// will be returned, indicating the event is selected to be kept.
+// Otherwise the value 0 will be returned.
+
+ if (!fSelections) return 1; 
+
+ Int_t nsel=fSelections->GetSize();
+ Int_t npart=GetN();
+ Int_t kf,kc;
+
+ Int_t select=0;
+ for (Int_t jpart=1; jpart<=npart; jpart++)
+ {
+  kf=GetK(jpart,2);
+  kc=Pycomp(kf);
+  for (Int_t i=0; i<nsel; i++)
+  {
+   if (kc==fSelections->At(i))
+   {
+    select=1;
+    break;
+   }
+  }
+  if (select) break;
+ }
+ return select;
+}
+///////////////////////////////////////////////////////////////////////////
+void AliCollider::SetSpectatorPmin(Float_t pmin)
+{
+// Set minimal momentum in GeV/c for spectator tracks to be stored.
+// Spectator tracks with a momentum below this threshold will not be stored
+// in the (output) event structure.
+// This facility allows to minimise the output file size.
+// Note that when the user wants to boost the event into another reference
+// frame these spectator tracks might have got momenta above the threshold.
+// However, when the spectator tracks were not stored in the event structure
+// in the original frame, there is no way to retreive them anymore. 
+ fSpecpmin=pmin;
+}
+///////////////////////////////////////////////////////////////////////////
+Float_t AliCollider::GetSpectatorPmin()
+{
+// Provide the minimal spectator momentum in GeV/c.
+ return fSpecpmin;
+}
+///////////////////////////////////////////////////////////////////////////
index 968e4e1f9301d32123ddb21c5dfce35cecab2af9..53bb0016d52f730de95ffcfdc77a3eba8cfbca05 100644 (file)
@@ -3,12 +3,13 @@
 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  * See cxx source for full Copyright notice                               */
 
-// $Id: AliCollider.h,v 1.5 2003/10/26 14:53:44 nick Exp $
+// $Id: AliCollider.h,v 1.6 2003/12/18 09:28:06 nick Exp $
 
 #include "TPythia6.h"
 #include "TString.h"
 #include "TFile.h"
-#include "TTree.h" 
+#include "TTree.h"
+#include "TArrayI.h" 
 
 #include "AliEvent.h"
 #include "AliRandom.h"
@@ -30,9 +31,13 @@ class AliCollider : public TPythia6
   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 SetStable(Int_t id,Int_t mode=1);                // Declare a certain particle as stable or not
+  void SelectEvent(Int_t id);                           // Select only events containing specified particles
+  Int_t GetSelectionFlag();                             // Return the selection flag for this event
   void MakeEvent(Int_t npt=0,Int_t mlist=-1,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
+  AliEvent* GetEvent(Int_t select=0);                   // Provide pointer to the generated event structure
+  void SetSpectatorPmin(Float_t pmin);                  // Set minimal momentum for spectator track to be stored
+  Float_t GetSpectatorPmin();                           // Provide the minimal momentum for spectator tracks
 
  protected:
   Int_t fVertexmode;    // The vertex structure creation mode
@@ -53,12 +58,17 @@ class AliCollider : public TPythia6
   Float_t fFracnn;      // Fraction of n+n collisions
   AliRandom fRan;       // Random number generator
   AliEvent* fEvent;     // The produced event structure
+  Float_t fSpecpmin;    // The minimal momentum for spectator tracks to be stored
 
   TFile* fOutFile;      // The user defined output data file 
   TTree* fOutTree;      // The standard ROOT output tree
 
+  TArrayI* fSelections; // The particle KC codes for event selection
+  Int_t fSelect;        // Flag to indicate whether the total event is selected (1) or not (0)
+
+  Int_t IsSelected();   // Check whether (sub)event passed the selection criteria
   void GetFractions(Float_t zp,Float_t ap,Float_t zt,Float_t at); // Determine various N-N collision fractions
 
- ClassDef(AliCollider,5) // Pythia based universal physics event generator
+ ClassDef(AliCollider,6) // Pythia based universal physics event generator
 };
 #endif
index 375b1f8f991b54ce79528f4f4bf33dce70c8caa0..0914cf9f9f9ee8f9746bb1e2517f59636c8a4bec 100644 (file)
@@ -13,7 +13,7 @@
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
 
-// $Id: AliEvent.cxx,v 1.16 2003/12/03 14:30:26 nick Exp $
+// $Id: AliEvent.cxx,v 1.17 2003/12/18 09:28:06 nick Exp $
 
 ///////////////////////////////////////////////////////////////////////////
 // Class AliEvent
 // Note : All quantities are in GeV, GeV/c or GeV/c**2
 //
 //--- Author: Nick van Eijndhoven 27-may-2001 UU-SAP Utrecht
-//- Modified: NvE $Date: 2003/12/03 14:30:26 $ UU-SAP Utrecht
+//- Modified: NvE $Date: 2003/12/18 09:28:06 $ UU-SAP Utrecht
 ///////////////////////////////////////////////////////////////////////////
 
 #include "AliEvent.h"
@@ -365,7 +365,7 @@ void AliEvent::SetDayTime(TTimeStamp& stamp)
 {
 // Set the date and time stamp for this event.
 // An exact copy of the entered date/time stamp will be saved with an
-// accuracy of 1 second.
+// accuracy of 1 nanosecond.
  fDaytime=stamp;
 }
 ///////////////////////////////////////////////////////////////////////////
@@ -483,12 +483,18 @@ Int_t AliEvent::GetTargetId()
 void AliEvent::HeaderData()
 {
 // Provide event header information
- cout << " *" << ClassName() << "::Data* Name : " << GetName()
-      << " Title : " << GetTitle() << endl;
+ const char* name=GetName();
+ const char* title=GetTitle();
+ Int_t ndevs=GetNdevices();
+ cout << " *" << ClassName() << "::Data*";
+ if (strlen(name))  cout << " Name : " << GetName();
+ if (strlen(title)) cout << " Title : " << GetTitle();
+ cout << endl;
  cout << "  " << fDaytime.AsString() << endl;
- cout << "  Run : " << fRun << " Event : " << fEvent << endl;
+ cout << "  Run : " << fRun << " Event : " << fEvent
+      << " Number of devices : " << ndevs << endl;
 
- ShowDevices();
if (ndevs) ShowDevices();
 }
 ///////////////////////////////////////////////////////////////////////////
 void AliEvent::Data(TString f)
index e202dc9eb892f6206060b6b0314298a938ff6356..d5fc6cb7b24125cb436cacf33348eb22d4b716f3 100644 (file)
@@ -559,7 +559,12 @@ void AliVertex::AddVertex(AliVertex& v,Int_t connect)
 void AliVertex::Data(TString f)
 {
 // Provide vertex information within the coordinate frame f
- cout << " *AliVertex::Data* Name : " << GetName() << " Title : " << GetTitle() << endl;
+ const char* name=GetName();
+ const char* title=GetTitle();
+ cout << " *AliVertex::Data*";
+ if (strlen(name))  cout << " Name : " << GetName();
+ if (strlen(title)) cout << " Title : " << GetTitle();
+ cout << endl;
  cout << " Id : " << fUserId << " Invmass : " << GetInvmass()
       << " Charge : " << GetCharge() << " Momentum : " << GetMomentum()
       << " Ntracks : " << GetNtracks() << " Nvertices : " << fNvtx 
index 2d49a9bb8b53a4ebcb21b48f6e794e4227fcfa76..ad3bb5f4a87c597383a32092ffedfc307b5caefd 100644 (file)
                 in case the corresponding AliTrack object is deleted.
 19-dec-2003 NvE Slight modification in AliVertex.h and AliAttribObj.cxx to prevent gcc compiler
                 warning concerning an unused variable.
+22-dec-2003 NvE Event selection introduced in AliCollider.
+24-dec-2003 NvE Some cosmetics in printout of AliVertex, AliEvent and AliCollider.
+                Support for spectator tracks introduced in AliCollider.
+28-dec-2003 NvE Facility introduced in AliCollider to set minimal momentum for spectator tracks
+                to be stored.