]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Add support for the physics event selection
authorakisiel <akisiel@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 14 Jan 2010 15:51:15 +0000 (15:51 +0000)
committerakisiel <akisiel@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 14 Jan 2010 15:51:15 +0000 (15:51 +0000)
PWG2/FEMTOSCOPY/AliFemto/AliFemtoBasicEventCut.cxx
PWG2/FEMTOSCOPY/AliFemto/AliFemtoBasicEventCut.h
PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorEventMult.cxx
PWG2/FEMTOSCOPY/AliFemto/AliFemtoCutMonitorEventMult.h
PWG2/FEMTOSCOPY/AliFemto/AliFemtoEvent.cxx
PWG2/FEMTOSCOPY/AliFemto/AliFemtoEvent.h
PWG2/FEMTOSCOPY/AliFemto/AliFemtoEventReaderESDChain.cxx
PWG2/FEMTOSCOPY/AliFemto/AliFemtoEventReaderESDChain.h

index 3a0bff2fae34f13a1472ebd67bc1afad3412bbd4..362aa4f369a4172d55bf7f95fc0043632c08ad25 100644 (file)
@@ -17,6 +17,7 @@ AliFemtoBasicEventCut::AliFemtoBasicEventCut() :
   fEventMult(),
   fVertZPos(),
   fAcceptBadVertex(false), 
+  fAcceptOnlyPhysics(true), 
   fNEventsPassed(0), 
   fNEventsFailed(0)
 {
@@ -44,7 +45,8 @@ bool AliFemtoBasicEventCut::Pass(const AliFemtoEvent* event){
      (mult <= fEventMult[1]) && 
      (vertexZPos > fVertZPos[0]) &&
      (vertexZPos < fVertZPos[1]) &&
-     (fAcceptBadVertex || (event->PrimVertCov()[4] > -1000.0)));
+     (fAcceptBadVertex || (event->PrimVertCov()[4] > -1000.0)) &&
+     ((!fAcceptOnlyPhysics) || (event->IsCollisionCandidate())));
   goodEvent ? fNEventsPassed++ : fNEventsFailed++ ;
 //   cout << "AliFemtoBasicEventCut:: return : " << goodEvent << endl;
   return (goodEvent);
@@ -71,3 +73,11 @@ bool AliFemtoBasicEventCut::GetAcceptBadVertex()
 {
   return fAcceptBadVertex;
 }
+void AliFemtoBasicEventCut::SetAcceptOnlyPhysics(bool b)
+{
+  fAcceptOnlyPhysics = b;
+}
+bool AliFemtoBasicEventCut::GetAcceptOnlyPhysics()
+{
+  return fAcceptOnlyPhysics;
+}
index 4790b3770cf41cf9961e6b29a420e6b18bffd10b..65d2807bc493692ba1b377ae725f47a3243b53e7 100644 (file)
@@ -26,9 +26,11 @@ public:
   void SetEventMult(const int& lo,const int& hi);
   void SetVertZPos(const float& lo, const float& hi);
   void SetAcceptBadVertex(bool b);
+  void SetAcceptOnlyPhysics(bool b);
   int NEventsPassed() const;
   int NEventsFailed() const;
   bool GetAcceptBadVertex();
+  bool GetAcceptOnlyPhysics();
 
   virtual AliFemtoString Report();
   virtual bool Pass(const AliFemtoEvent* event);
@@ -42,6 +44,7 @@ private:   // here are the quantities I want to cut on...
   bool fAcceptBadVertex;  // Set to true to accept events with bad vertex
   long fNEventsPassed;    // Number of events checked by this cut that passed
   long fNEventsFailed;    // Number of events checked by this cut that failed
+  bool fAcceptOnlyPhysics;// Accept only physics events
 
 #ifdef __ROOT__
   ClassDef(AliFemtoBasicEventCut, 1)
index df3b752bf4626e7393d89d37ea988448d23cf6f0..d29423aec03d7e81fb28de19582819e58f6310a0 100644 (file)
@@ -12,7 +12,8 @@
 #include <TList.h>
 
 AliFemtoCutMonitorEventMult::AliFemtoCutMonitorEventMult():
-  fEvMult(0)
+  fEvMult(0),
+  fNormEvMult(0)
 {
   // Default constructor
   fEvMult = new TH1D("EvMult", "Event Multiplicity", 5001, -0.5, 5000.5);
@@ -20,27 +21,36 @@ AliFemtoCutMonitorEventMult::AliFemtoCutMonitorEventMult():
 
 AliFemtoCutMonitorEventMult::AliFemtoCutMonitorEventMult(const char *aName):
   AliFemtoCutMonitor(),
-  fEvMult(0)
+  fEvMult(0),
+  fNormEvMult(0)
 {
   // Normal constructor
   char name[200];
   snprintf(name, 200, "EvMult%s", aName);
   fEvMult = new TH1D(name, "Event Multiplicity", 5001, -0.5, 5000.5);
+
+  snprintf(name, 200, "NormEvMult%s", aName);
+  fNormEvMult = new TH1D(name, "Normalized Event Multiplicity", 5001, -0.5, 5000.5);
 }
 
 AliFemtoCutMonitorEventMult::AliFemtoCutMonitorEventMult(const AliFemtoCutMonitorEventMult &aCut):
   AliFemtoCutMonitor(),
-  fEvMult(0)
+  fEvMult(0),
+  fNormEvMult(0)
 {
   // copy constructor
   if (fEvMult) delete fEvMult;
   fEvMult = new TH1D(*aCut.fEvMult);
+
+  if (fNormEvMult) delete fNormEvMult;
+  fNormEvMult = new TH1D(*aCut.fNormEvMult);
 }
 
 AliFemtoCutMonitorEventMult::~AliFemtoCutMonitorEventMult()
 {
   // Destructor
   delete fEvMult;
+  delete fNormEvMult;
 }
 
 AliFemtoCutMonitorEventMult& AliFemtoCutMonitorEventMult::operator=(const AliFemtoCutMonitorEventMult& aCut)
@@ -52,6 +62,9 @@ AliFemtoCutMonitorEventMult& AliFemtoCutMonitorEventMult::operator=(const AliFem
   if (fEvMult) delete fEvMult;
   fEvMult = new TH1D(*aCut.fEvMult);
   
+  if (fNormEvMult) delete fNormEvMult;
+  fNormEvMult = new TH1D(*aCut.fNormEvMult);
+  
   return *this;
 }
 
@@ -66,18 +79,21 @@ void AliFemtoCutMonitorEventMult::Fill(const AliFemtoEvent* aEvent)
 {
   // Fill in the monitor histograms with the values from the current track
   fEvMult->Fill(aEvent->NumberOfTracks());
+  fNormEvMult->Fill(aEvent->UncorrectedNumberOfPrimaries());
 }
 
 void AliFemtoCutMonitorEventMult::Write()
 {
   // Write out the relevant histograms
   fEvMult->Write();
+  fNormEvMult->Write();
 }
 
 TList *AliFemtoCutMonitorEventMult::GetOutputList()
 {
   TList *tOutputList = new TList();
   tOutputList->Add(fEvMult);
+  tOutputList->Add(fNormEvMult);
 
   return tOutputList;
 }
index b961ebe1b7a501415ffb237ffee62f23badd68a3..1289386026c6de2501a69a0461529e7cb8e2194e 100644 (file)
@@ -21,7 +21,7 @@ class TList;
 
 class AliFemtoCutMonitorEventMult : public AliFemtoCutMonitor{
   
-public:
+ public:
   AliFemtoCutMonitorEventMult();
   AliFemtoCutMonitorEventMult(const char *aName);
   AliFemtoCutMonitorEventMult(const AliFemtoCutMonitorEventMult &aCut);
@@ -43,8 +43,9 @@ public:
 
   virtual TList *GetOutputList();
 
-private:
-  TH1D *fEvMult;    // Multiplicity distribution
+ private:
+  TH1D *fEvMult;     // Multiplicity distribution
+  TH1D *fNormEvMult; // Normalized event multiplicity distribution
 };
 
 #endif
index 1a4ac50fece72fb99fe7393d41b2125eefc24abe..004edf006a29822feb63dd576abda098081798d2 100644 (file)
@@ -29,6 +29,7 @@ AliFemtoEvent::AliFemtoEvent():
   fRunNumber(0),
   fNumberOfTracks(0),
   fMagneticField(0),
+  fIsCollisionCandidate(kTRUE),
   fPrimVertPos(0,0,0),
   fPrimVertCov(),
   fTrackCollection(0),
@@ -67,6 +68,7 @@ AliFemtoEvent::AliFemtoEvent(const AliFemtoEvent& ev, AliFemtoTrackCut* tCut, Al
   fRunNumber(0),
   fNumberOfTracks(0),
   fMagneticField(0),
+  fIsCollisionCandidate(kTRUE),
   fPrimVertPos(0,0,0),
   fPrimVertCov(),
   fTrackCollection(0),
@@ -97,7 +99,8 @@ AliFemtoEvent::AliFemtoEvent(const AliFemtoEvent& ev, AliFemtoTrackCut* tCut, Al
   fZDCParticipants=ev.fZDCParticipants;
   fNumberOfTracks = ev.fNumberOfTracks;
   fMagneticField= ev.fMagneticField;
-  
+  fIsCollisionCandidate = ev.fIsCollisionCandidate;
+
   fTriggerMask=ev.fTriggerMask;     // Trigger Type (mask)
   fTriggerCluster=ev.fTriggerCluster;
   fReactionPlaneAngle=ev.fReactionPlaneAngle;
@@ -143,6 +146,7 @@ AliFemtoEvent::AliFemtoEvent(const AliFemtoEvent& ev):
   fRunNumber(0),
   fNumberOfTracks(0),
   fMagneticField(0),
+  fIsCollisionCandidate(kTRUE),
   fPrimVertPos(0,0,0),
   fPrimVertCov(),
   fTrackCollection(0),
@@ -173,7 +177,7 @@ AliFemtoEvent::AliFemtoEvent(const AliFemtoEvent& ev):
   fZDCParticipants=ev.fZDCParticipants;
   fNumberOfTracks = ev.fNumberOfTracks;
   fMagneticField= ev.fMagneticField;
-  
+  fIsCollisionCandidate = ev.fIsCollisionCandidate;
   fTriggerMask=ev.fTriggerMask;     // Trigger Type (mask)
   fTriggerCluster=ev.fTriggerCluster;
   fReactionPlaneAngle=ev.fReactionPlaneAngle;
@@ -222,7 +226,8 @@ AliFemtoEvent& AliFemtoEvent::operator=(const AliFemtoEvent& aEvent)
   fZDCParticipants=aEvent.fZDCParticipants;
   fNumberOfTracks = aEvent.fNumberOfTracks;
   fMagneticField= aEvent.fMagneticField;
-  
+  fIsCollisionCandidate = aEvent.fIsCollisionCandidate;
+
   fTriggerMask=aEvent.fTriggerMask;     // Trigger Type (mask)
   fTriggerCluster=aEvent.fTriggerCluster;
   fReactionPlaneAngle=aEvent.fReactionPlaneAngle;
@@ -315,6 +320,7 @@ void AliFemtoEvent::SetPrimVertCov(const double* v){
   fPrimVertCov[5] = v[5];
 }
 void AliFemtoEvent::SetMagneticField(const double& magF){fMagneticField = magF;}
+void AliFemtoEvent::SetIsCollisionCandidate(const bool& is){fIsCollisionCandidate = is;}
 
 void AliFemtoEvent::SetTriggerMask(const unsigned long int& aTriggerMask) {fTriggerMask=aTriggerMask;}
 void AliFemtoEvent::SetTriggerCluster(const unsigned char& aTriggerCluster) {fTriggerCluster=aTriggerCluster;}
@@ -336,6 +342,7 @@ const double* AliFemtoEvent::PrimVertCov() const {return fPrimVertCov;}
 double AliFemtoEvent::MagneticField() const {return fMagneticField;}
 unsigned long int AliFemtoEvent::TriggerMask() const {return fTriggerMask;}
 unsigned char AliFemtoEvent::TriggerCluster() const {return fTriggerCluster;}
+bool AliFemtoEvent::IsCollisionCandidate() const {return fIsCollisionCandidate;}
 
 
 float AliFemtoEvent::ZDCN1Energy() const {return fZDCN1Energy;}       
@@ -357,5 +364,19 @@ double AliFemtoEvent::UncorrectedNumberOfNegativePrimaries() const
 
 double AliFemtoEvent::UncorrectedNumberOfPrimaries() const
 {
-  return NumberOfTracks();
+  // Count number of normalized charged tracks 
+  Int_t tNormTrackCount = 0;
+  for (AliFemtoTrackIterator iter=fTrackCollection->begin();iter!=fTrackCollection->end();iter++){
+    if (!((*iter)->Flags()&(AliFemtoTrack::kTPCrefit))) continue;
+    if ((*iter)->TPCncls() < 50) continue;
+    if ((*iter)->TPCchi2()/(*iter)->TPCncls() > 60.0) continue;
+    if ((*iter)->ImpactD() > 6.0) continue;
+    if ((*iter)->ImpactZ() > 6.0) continue;
+    if (fabs((*iter)->P().pseudoRapidity()) > 0.9) continue;
+
+    tNormTrackCount++;
+  }
+
+  return tNormTrackCount;
+  //  return NumberOfTracks();
 }
index f4f48792f75358b68f71fb68703f01060c8bb18d..1d63559d978f6e17b42fa5a978d51cd17cd5247e 100644 (file)
@@ -55,6 +55,7 @@ public:
   AliFemtoKinkCollection* KinkCollection() const;
   AliFemtoTrackCollection* TrackCollection() const;
   double MagneticField() const;
+  bool IsCollisionCandidate() const;
 
   //functions for alice variables
   float ZDCN1Energy() const;      
@@ -75,7 +76,8 @@ public:
   void SetPrimVertPos(const AliFemtoThreeVector& v);
   void SetPrimVertCov(const double* v);
   void SetMagneticField(const double& x);
-  
+  void SetIsCollisionCandidate(const bool& is);
+
    //functions for alice variables
   void SetZDCN1Energy(const float& x);      
   void SetZDCP1Energy(const float& x);      
@@ -97,7 +99,8 @@ private:
   unsigned short fRunNumber;             // run number the event belong to
   unsigned short fNumberOfTracks;        // total number of TPC tracks
   double fMagneticField;                 // magnetic field in Z direction
-
+  bool fIsCollisionCandidate;            // is collision candidate
+  
   AliFemtoThreeVector fPrimVertPos;      // primary vertex position
   double fPrimVertCov[6];                // primary vertex covariances
   AliFemtoTrackCollection* fTrackCollection; // collection of tracks
index 0e9f29eea760d5ff8053476ec93f392ad14bf0b0..ca84992482cdd3fec4fcdc5558b00f807deee956 100644 (file)
@@ -38,7 +38,10 @@ AliFemtoEventReaderESDChain::AliFemtoEventReaderESDChain():
   fNumberofEvent(0),
   fCurEvent(0),
   fCurFile(0),
-  fEvent(0x0)
+  fEvent(0x0),
+  fUsePhysicsSel(kFALSE),
+  fSelect(0x0),
+  fTrackType(kGlobal)
 {
   //constructor with 0 parameters , look at default settings 
 //   fClusterPerPadrow = (list<Int_t> **) malloc(sizeof(list<Int_t> *) * AliESDfriendTrack::kMaxTPCcluster);
@@ -61,7 +64,10 @@ AliFemtoEventReaderESDChain::AliFemtoEventReaderESDChain(const AliFemtoEventRead
   fNumberofEvent(0),
   fCurEvent(0),
   fCurFile(0),
-  fEvent(0x0)
+  fEvent(0x0),
+  fUsePhysicsSel(kFALSE),
+  fSelect(0x0),
+  fTrackType(kGlobal)
 {
   // Copy constructor
   fConstrained = aReader.fConstrained;
@@ -72,6 +78,10 @@ AliFemtoEventReaderESDChain::AliFemtoEventReaderESDChain(const AliFemtoEventRead
   fCurFile = aReader.fCurFile;
   //  fEvent = new AliESD(*aReader.fEvent);
   fEvent = new AliESDEvent();
+  fUsePhysicsSel = aReader.fUsePhysicsSel;
+  if (aReader.fUsePhysicsSel)
+    fSelect = new AliPhysicsSelection();
+  fTrackType = aReader.fTrackType;
 //   fEventFriend = aReader.fEventFriend;
 //   fClusterPerPadrow = (list<Int_t> **) malloc(sizeof(list<Int_t> *) * AliESDfriendTrack::kMaxTPCcluster);
 //   for (int tPad=0; tPad<AliESDfriendTrack::kMaxTPCcluster; tPad++) {
@@ -106,6 +116,7 @@ AliFemtoEventReaderESDChain::~AliFemtoEventReaderESDChain()
 //     delete fSharedList[tPad];
 //   }
 //   delete [] fSharedList;
+  if (fSelect) delete fSelect;
 }
 
 //__________________
@@ -123,7 +134,11 @@ AliFemtoEventReaderESDChain& AliFemtoEventReaderESDChain::operator=(const AliFem
   fCurFile = aReader.fCurFile;
   if (fEvent) delete fEvent;
   fEvent = new AliESDEvent();
+  fTrackType = aReader.fTrackType;
 
+  fUsePhysicsSel = aReader.fUsePhysicsSel;
+  if (aReader.fUsePhysicsSel)
+    fSelect = new AliPhysicsSelection();
   //  fEventFriend = aReader.fEventFriend;
   
 //   if (fClusterPerPadrow) {
@@ -203,6 +218,12 @@ bool AliFemtoEventReaderESDChain::GetUseTPCOnly() const
   return fUseTPCOnly;
 }
 
+void AliFemtoEventReaderESDChain::SetUsePhysicsSelection(const bool usephysics)
+{
+  fUsePhysicsSel = usephysics;
+  if (!fSelect) fSelect = new AliPhysicsSelection();
+}
+
 AliFemtoEvent* AliFemtoEventReaderESDChain::ReturnHbtEvent()
 {
   // Get the event, read all the relevant information
@@ -215,8 +236,17 @@ AliFemtoEvent* AliFemtoEventReaderESDChain::ReturnHbtEvent()
   cout<<"starting to read event "<<fCurEvent<<endl;
   //  fEvent->SetESDfriend(fEventFriend);
   if(fEvent->GetAliESDOld())fEvent->CopyFromOldESD();
-       
+  
   hbtEvent = new AliFemtoEvent;
+
+  if (fUsePhysicsSel) {
+    hbtEvent->SetIsCollisionCandidate(fSelect->IsCollisionCandidate(fEvent));
+    if (!(fSelect->IsCollisionCandidate(fEvent)))
+      printf("Event not a collision candidate\n");
+  }
+  else
+    hbtEvent->SetIsCollisionCandidate(kTRUE);
+
   //setting basic things
   //  hbtEvent->SetEventNumber(fEvent->GetEventNumber());
   hbtEvent->SetRunNumber(fEvent->GetRunNumber());
@@ -290,10 +320,24 @@ AliFemtoEvent* AliFemtoEventReaderESDChain::ReturnHbtEvent()
     {
       bool  tGoodMomentum=true; //flaga to chcek if we can read momentum of this track
                
-      AliFemtoTrack* trackCopy = new AliFemtoTrack();  
       const AliESDtrack *esdtrack=fEvent->GetTrack(i);//getting next track
       //      const AliESDfriendTrack *tESDfriendTrack = esdtrack->GetFriendTrack();
 
+      // If reading ITS-only tracks, reject all with TPC
+      if (fTrackType == kITSOnly) {
+       if (esdtrack->GetStatus() & AliESDtrack::kTPCrefit) continue;
+       if (!(esdtrack->GetStatus() & AliESDtrack::kITSrefit)) continue;
+       if (esdtrack->GetStatus() & AliESDtrack::kTPCin) continue;
+       UChar_t iclm = esdtrack->GetITSClusterMap();
+       Int_t incls = 0;
+       for (int iter=0; iter<6; iter++) if (iclm&(1<<iter)) incls++;
+       if (incls<=3) {
+         cout << "Rejecting track with " << incls << " clusters" << endl;
+         continue;
+       }
+      }
+
+      AliFemtoTrack* trackCopy = new AliFemtoTrack();  
       trackCopy->SetCharge((short)esdtrack->GetSign());
 
       //in aliroot we have AliPID 
@@ -526,6 +570,10 @@ Float_t AliFemtoEventReaderESDChain::GetSigmaToVertex(double *impact, double *co
   return d;
 }
 
+void AliFemtoEventReaderESDChain::SetReadTrackType(ReadTrackType aType)
+{
+  fTrackType = aType;
+}
 
 
 
index 80d1f442b0158d1a5c6ab00c5f7a7f1939f9074f..517d052bcfb23c591a64e0319b6d4b7833205614 100644 (file)
@@ -17,6 +17,7 @@
 #include "TTree.h"
 #include "AliESDEvent.h"
 #include "AliESDfriend.h"
+#include "AliPhysicsSelection.h"
 #include <list>
 
 class AliFemtoEvent;
@@ -24,6 +25,9 @@ class AliFemtoEvent;
 class AliFemtoEventReaderESDChain : public AliFemtoEventReader 
 {
  public:
+  enum TrackType {kGlobal=0, kTPCOnly=1, kITSOnly=2, kSPDTracklet=3};
+  typedef enum TrackType ReadTrackType;
+
   AliFemtoEventReaderESDChain();
   AliFemtoEventReaderESDChain(const AliFemtoEventReaderESDChain& aReader);
   ~AliFemtoEventReaderESDChain();
@@ -35,10 +39,13 @@ class AliFemtoEventReaderESDChain : public AliFemtoEventReader
   void SetConstrained(const bool constrained);
   void SetReadTPCInner(const bool readinner);
   void SetUseTPCOnly(const bool usetpconly);
+  void SetUsePhysicsSelection(const bool usephysics);
 
   bool GetConstrained() const;
   bool GetReadTPCInner() const;
   bool GetUseTPCOnly() const;
+  
+  void SetReadTrackType(ReadTrackType aType);
 
   void SetESDSource(AliESDEvent *aESD);
   //  void SetESDfriendSource(AliESDfriend *aFriend);
@@ -57,6 +64,9 @@ class AliFemtoEventReaderESDChain : public AliFemtoEventReader
   unsigned int   fCurFile;       //number of current file
   AliESDEvent*   fEvent;         //ESD event
   //  AliESDfriend*  fEventFriend;
+  bool           fUsePhysicsSel; //if true the physics selection class will be used
+  AliPhysicsSelection *fSelect;  //Class to select only physics events
+  ReadTrackType  fTrackType;     // Type of track read
 
 /*   list<Int_t>  **fSharedList;       //! Table (one list per padrow) of clusters which are shared */
 /*   list<Int_t>  **fClusterPerPadrow; //! Table (one list per padrow) of clusters in each padrow */