2nd phase of AOD implementation - reading w/o checking the cuts
authorskowron <skowron@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 17 Aug 2004 18:26:21 +0000 (18:26 +0000)
committerskowron <skowron@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 17 Aug 2004 18:26:21 +0000 (18:26 +0000)
ANALYSIS/AliAOD.cxx
ANALYSIS/AliAOD.h
ANALYSIS/AliReader.cxx
ANALYSIS/AliReader.h
ANALYSIS/AliReaderAOD.cxx
ANALYSIS/AliReaderAOD.h
ANALYSIS/WriteAOD.C
ANALYSIS/analysis.C

index 4c620e7..53925eb 100644 (file)
@@ -46,6 +46,37 @@ AliAOD::AliAOD():
 }
 /**************************************************************************/
 
+AliAOD::AliAOD(const AliAOD& in):
+ TObject(in),
+ fParticles((TClonesArray*)in.fParticles->Clone()),
+  fIsRandomized(in.fIsRandomized),
+  fPrimaryVertexX(fPrimaryVertexX),
+  fPrimaryVertexY(in.fPrimaryVertexY),
+  fPrimaryVertexZ(in.fPrimaryVertexZ),
+  fParticleClass(in.fParticleClass)
+{
+//copy constructor
+}
+/**************************************************************************/
+
+AliAOD& AliAOD::operator=(const AliAOD& in)
+{
+//assigment operator  
+
+  if (this == &in ) return *this;
+  
+  delete fParticles;
+  fParticles = (TClonesArray*)in.fParticles->Clone();
+  fIsRandomized = in.fIsRandomized ;
+  fPrimaryVertexX = in.fPrimaryVertexX ;
+  fPrimaryVertexY = in.fPrimaryVertexY ;
+  fPrimaryVertexZ = in.fPrimaryVertexZ ;
+  fParticleClass = in.fParticleClass ; //althought it is pointer, this points to object in class list of gROOT
+  return *this;
+}
+
+/**************************************************************************/
+
 AliAOD::~AliAOD()
 {
   //Destructor
@@ -55,6 +86,85 @@ AliAOD::~AliAOD()
 }
 /**************************************************************************/
 
+void AliAOD::CopyData(AliAOD* aod)
+{
+ //Copys all data from aod, but leaves local type of particles
+ if (aod == 0x0) return;
+ if (aod == this) return;
+ AliAOD& in = *this;
+ fIsRandomized = in.fIsRandomized ;
+ fPrimaryVertexX = in.fPrimaryVertexX ;
+ fPrimaryVertexY = in.fPrimaryVertexY ;
+ fPrimaryVertexZ = in.fPrimaryVertexZ ;
+ fParticleClass = in.fParticleClass ; //althought it is pointer, this points to object in class list of gROOT
+
+ if (in.fParticles == 0x0)
+  {//if in obj has null fParticles we delete ours
+    delete fParticles;
+    fParticles = 0x0;
+  }
+ else
+  { 
+    if (fParticles)
+     { //if ours particles were already created
+       if (fParticles->GetClass() != in.fParticles->GetClass())
+        {//if in obj has 
+          delete fParticles;
+          fParticles = (TClonesArray*)in.fParticles->Clone();
+        }
+       else
+        {
+         //it should be faster than cloning
+          Int_t inentr = in.fParticles->GetEntriesFast();
+          Int_t curentr = fParticles->GetEntriesFast();
+
+          TClonesArray& arr = *fParticles;
+
+          //we have to take care about different sizes of arrays
+          if ( curentr < inentr )
+           {
+             for (Int_t i = 0; i < curentr; i++)
+              {
+                TObject& inobj = *(in.fParticles->At(i));
+                TObject& obj = *(fParticles->At(i));
+                obj = inobj;
+              }
+
+             for (Int_t i = curentr; i < inentr; i++)
+              {
+                TObject& inobj = *(in.fParticles->At(i));
+                TObject& obj =  *((TObject*)(fParticleClass->New(arr[i])));
+                obj = inobj;
+              }
+           }
+          else 
+           {
+             for (Int_t i = 0; i < inentr; i++)
+              {
+                TObject& inobj = *(in.fParticles->At(i));
+                TObject& obj = *(fParticles->At(i));
+                obj = inobj;
+              }
+             
+             for (Int_t i = curentr ; i >= inentr ; i--)
+              {
+                fParticles->RemoveAt(i);
+              }
+           }
+        } 
+     }
+    else
+     {
+       fParticles = (TClonesArray*)in.fParticles->Clone();
+     } 
+  } 
+}
+/**************************************************************************/
+
 void AliAOD::SetParticleClassName(const char* classname)
 {
 //Sets type of particle that is going to be stored 
index 2688f3b..c4ae46c 100644 (file)
@@ -21,8 +21,12 @@ class AliAOD: public TObject {
 public:
   AliAOD();
   virtual ~AliAOD();
-
-  virtual TClonesArray*    GetParticles() {return fParticles;};
+  
+  AliAOD(const AliAOD& in);
+  virtual AliAOD& operator=(const AliAOD& in);
+  virtual void             CopyData(AliAOD* aod);//Copys all data from aod, but leaves local type of particles
+  
+  virtual TClonesArray*    GetParticles() const {return fParticles;} 
   virtual void             SetParticleClassName(const char* classname);
   virtual void             SetParticleClass(TClass* pclass);
   
@@ -34,6 +38,7 @@ public:
                                        Double_t vx, Double_t vy, Double_t vz, Double_t time);
   
   virtual void             Reset();
+  
   void                     SwapParticles(Int_t i, Int_t j);//swaps particles positions; used by AliReader::Blend
   Bool_t                   IsRandomized() const {return fIsRandomized;}
   void                     SetRandomized(Bool_t flag = kTRUE){fIsRandomized = flag;}
@@ -45,6 +50,7 @@ public:
   void                     Move(Double_t x, Double_t y, Double_t z);//moves all spacial coordinates about this vector
   virtual void             SetOwner(Bool_t owner);
   virtual void             Print(Option_t* /*option*/ = 0);
+  const TClass*            GetParticleClass() const {return fParticleClass;}
 private:
   TClonesArray            *fParticles;   // array of AOD particles, AliAOD is owner of particles
   Bool_t                   fIsRandomized;//flag indicating if positions of particles were randomized - used by HBTAN
index 41bbf11..230bad1 100644 (file)
@@ -388,27 +388,33 @@ Bool_t  AliReader::Rejected(Int_t pid)
 }
 /*************************************************************************************/
 
-TString& AliReader::GetDirName(Int_t entry)
+TString AliReader::GetDirName(Int_t entry)
 {
 //returns directory name of next one to read
-  TString* retval;//return value
+  TString  retval;//return value
   if (fDirs ==  0x0)
-   {
-     retval = new TString(".");
-     return *retval;
+   { 
+     if (entry == 0)
+      {
+       retval = ".";
+       return retval;
+      }
+     else
+      {
+        return retval;
+      }  
    }
 
   if ( (entry>fDirs->GetEntries()) || (entry<0))//if out of bounds return empty string
    {                                            //note that entry==0 is accepted even if array is empty (size=0)
      Error("GetDirName","Name out of bounds");
-     retval = new TString();
-     return *retval;
+     return retval;
    }
 
   if (fDirs->GetEntries() == 0)
    { 
-     retval = new TString(".");
-     return *retval;
+     retval = ".";
+     return retval;
    }
 
   TClass *objclass = fDirs->At(entry)->IsA();
@@ -419,11 +425,11 @@ TString& AliReader::GetDirName(Int_t entry)
   if(dir == 0x0)
    {
      Error("GetDirName","Object in TObjArray is not a TObjString or its descendant");
-     retval = new TString();
-     return *retval;
+     return retval;
    }
   if (gDebug > 0) Info("GetDirName","Returned ok %s",dir->String().Data());
-  return dir->String();
+  retval = dir->String();
+  return retval;
 }
 /*************************************************************************************/
 
index 2411aa8..fd04fcd 100644 (file)
@@ -129,7 +129,7 @@ class AliReader: public TNamed
     Bool_t               Rejected(Int_t pid);//Checks if a given pid passes cuts
     void                 Blend();//Mixes current events in a symmetric way so after mixing thy are consistent
     
-    TString&             GetDirName(Int_t entry);
+    TString              GetDirName(Int_t entry);
     
   private:
   
index 6575d6c..3e4b38f 100644 (file)
@@ -7,6 +7,151 @@ ClassImp(AliReaderAOD)
 #include <TTree.h>
 #include "AliAOD.h"
 
+
+const TString AliReaderAOD::fgkTreeName("TAOD");
+const TString AliReaderAOD::fgkRecosntructedDataBranchName("reconstructed.");
+const TString AliReaderAOD::fgkSimulatedDataBranchName("simulated.");
+
+AliReaderAOD::AliReaderAOD(const Char_t* aodfilename):
+ fFileName(aodfilename),
+ fReadSim(kFALSE),
+ fTree(0x0),
+ fFile(0x0),
+ fSimBuffer(0x0),
+ fRecBuffer(0x0)
+{
+  //ctor
+}
+/********************************************************************/
+
+AliReaderAOD::~AliReaderAOD()
+{
+//dtor
+  if (fEventSim == fSimBuffer )
+   {
+    fEventSim = 0x0;
+    fEventRec = 0x0;
+   }
+  delete fSimBuffer;
+  delete fRecBuffer;
+  
+  delete fTree;
+  delete fFile;
+}
+/********************************************************************/
+
+void AliReaderAOD::Rewind()
+{
+//Rewinds reading
+  delete fTree;
+  fTree = 0x0;
+  delete fFile;
+  fFile = 0x0;
+  fCurrentDir = 0;
+  fNEventsRead= 0;
+}
+/********************************************************************/
+Int_t AliReaderAOD::ReadNext()
+{
+//Reads next event
+  
+  do  //do{}while; is OK even if 0 dirs specified. In that case we try to read from "./"
+    {
+      if (fFile == 0x0)
+       {
+         Int_t opened = OpenFile(fCurrentDir);//rl is opened here
+         if (opened)
+           {
+             //Error("ReadNext","Error Occured while opening directory number %d",fCurrentDir);
+             fCurrentDir++;
+             continue;
+           }
+         fCurrentEvent = 0;
+       }
+      
+      fTree->GetEvent(fCurrentEvent);
+      
+      //Temporary testing sollution
+      fEventSim = fSimBuffer;
+      fEventRec = fRecBuffer;
+      
+      fCurrentEvent++;
+      fNEventsRead++;
+
+      if (fTree)
+       {
+         if ( fCurrentEvent >= fTree->GetEntries() )
+          {
+            delete fTree;
+            fTree = 0x0;
+            delete fFile;
+            fFile = 0x0;
+            fCurrentDir++;
+          } 
+       }
+
+
+      return 0;//success -> read one event
+      
+    }while(fCurrentDir < GetNumberOfDirs());//end of loop over directories specified in fDirs Obj Array  
+   
+  return 1; //no more directories to read
+  
+  
+}
+/********************************************************************/
+
+Int_t AliReaderAOD::OpenFile(Int_t n)
+{
+//opens fFile with  tree
+
+ const TString& dirname = GetDirName(n);
+ if (dirname == "")
+  {
+    if (AliVAODParticle::GetDebug() > 2 )
+      {
+        Info("OpenFile","Got empty string as a directory name."); 
+      }
+   return 1;
+  }
+ TString filename = dirname +"/"+ fFileName;
+ fFile = TFile::Open(filename.Data()); 
+ if ( fFile == 0x0)
+  {
+    Error("OpenFile","Can't open fFile %s",filename.Data());
+    return 2;
+  }
+ if (!fFile->IsOpen())
+  {
+    Error("OpenFile","Can't open fFile  %s",filename.Data());
+    delete fFile;
+    fFile = 0x0;
+    return 3;
+  }
+ fTree = dynamic_cast<TTree*>(fFile->Get(fgkTreeName));
+ if (fTree == 0x0)
+  {
+    if (AliVAODParticle::GetDebug() > 2 )
+      {
+        Info("ReadNext","Can not find TTree object named %s",fgkTreeName.Data());
+      }
+    fCurrentDir++;
+    delete fFile;//we have to assume there is no more ESD objects in the fFile
+    fFile = 0x0;
+    return 4;
+  }
+
+  fTree->SetBranchAddress(fgkSimulatedDataBranchName,&fSimBuffer);
+  fTree->SetBranchAddress(fgkRecosntructedDataBranchName,&fRecBuffer);
+  
+  return 0;
+}
+
+/********************************************************************/
+
 Int_t AliReaderAOD::WriteAOD(AliReader* reader, const char* outfilename, const char* pclassname,  Bool_t /*multcheck*/)
 {
 //reads tracks from runs and writes them to file
@@ -26,37 +171,52 @@ Int_t AliReaderAOD::WriteAOD(AliReader* reader, const char* outfilename, const c
      return -1;
    }
 
-  TTree *tree = new TTree("TAOD","Tree with tracks");
+  TTree *tree = new TTree(fgkTreeName,"Tree with tracks");
   
   TBranch *recbranch = 0x0, *simbranch = 0x0;
   
-  
   AliAOD* eventrec = new AliAOD();//must be created before Branch is called. Otherwise clones array is not splitted
   AliAOD* eventsim = new AliAOD();//AOD together with fParticles clones array knowing exact type of particles
   
   eventrec->SetParticleClassName(pclassname);
   eventsim->SetParticleClassName(pclassname);
+  AliAOD* recbuffer = eventrec;
+  AliAOD* simbuffer = eventsim;
   
-  if (reader->ReadsRec()) recbranch = tree->Branch("reconstructed","AliAOD",&eventrec,32000,99);
-  if (reader->ReadsSim()) simbranch = tree->Branch("simulated","AliAOD",&eventsim,32000,99);
+  if (reader->ReadsRec()) recbranch = tree->Branch(fgkRecosntructedDataBranchName,"AliAOD",&recbuffer,32000,99);
+  if (reader->ReadsSim()) simbranch = tree->Branch(fgkSimulatedDataBranchName,"AliAOD",&simbuffer,32000,99);
 
-  delete eventsim;
-  delete eventrec;
-  
   reader->Rewind();
   while (reader->Next() == kFALSE)
    {
      
      if (reader->ReadsRec())
-      {
-        eventrec = reader->GetEventRec();
-        recbranch->SetAddress(&eventrec);
+      {//here we can get AOD that has different particle type
+        AliAOD* event = reader->GetEventRec();
+        if ( eventrec->GetParticleClass() != event->GetParticleClass() )
+         {//if class type is not what what we whant we copy particles
+           eventrec->CopyData(event);
+           recbuffer = eventrec;
+         }
+        else
+         {//else just pointer to event from input reader is passed
+           recbuffer = event;
+         } 
       }
 
      if (reader->ReadsSim())
       {
-        eventsim = reader->GetEventSim();
-        simbranch->SetAddress(&eventsim);
+        AliAOD* event = reader->GetEventRec();
+        if ( eventsim->GetParticleClass() != event->GetParticleClass() )
+         {//if class type is not what what we whant we copy particles
+           eventsim->CopyData(event);
+           simbuffer = eventrec;
+         }
+        else
+         {//else just pointer to event from input reader is passed
+           simbuffer = event;
+         } 
       }
      eventrec->GetParticle(0)->Print();
      eventsim->GetParticle(0)->Print();
@@ -66,6 +226,10 @@ Int_t AliReaderAOD::WriteAOD(AliReader* reader, const char* outfilename, const c
   ::Info("AliReaderAOD::Write","Written %d events",tree->GetEntries());
   outfile->cd();
   tree->Write();
+
+  delete eventsim;
+  delete eventrec;
+  
   delete tree;
   delete outfile;
   return 0; 
index 52f8b34..b7e4592 100644 (file)
@@ -3,25 +3,45 @@
 
 #include "AliReader.h"
 
+class TTree;
+class TFile;
+
 class AliReaderAOD: public AliReader
 {
   public:
-    AliReaderAOD(const Char_t* aodfilename = "AliAOD.root"){}
-    virtual ~AliReaderAOD(){}
+    AliReaderAOD(const Char_t* aodfilename = "AOD.root");
+    virtual ~AliReaderAOD();
 
     void          ReadSimulatedData(Bool_t flag){fReadSim = flag;}//switches reading MC data
     Bool_t        ReadsRec() const {return kTRUE;}
     Bool_t        ReadsSim() const {return fReadSim;}
 
+    void          Rewind();
+
 
     static Int_t WriteAOD(AliReader* reader, const char* outfilename = "AliAOD.root", //reads tracks from runs and writes them to file
                           const char* pclassname = "AliAODParticle", Bool_t multcheck = kFALSE);
     
   protected:
+    virtual Int_t         ReadNext();
+    virtual Int_t        OpenFile(Int_t evno);//opens files to be read for given event
+
+    static const TString  fgkTreeName;//name of branch holding simulated data
+    static const TString  fgkRecosntructedDataBranchName;//name of branch holding reconstructed data
+    static const TString  fgkSimulatedDataBranchName;//name of branch holding simulated data
+    
+  
   private:
     TString fFileName;//File name
     
     Bool_t  fReadSim;//indicates if to read simulated data
+
+    TTree*        fTree;//!tree
+    TFile*        fFile;//!file
+    AliAOD*       fSimBuffer;//!buffer array that tree is read to
+    AliAOD*       fRecBuffer;//!
+    
+    
     ClassDef(AliReaderAOD,1)
 };
 
index bb85d39..be01e60 100644 (file)
@@ -12,8 +12,8 @@
 #endif
 
 
-void WriteAOD(Option_t* datatype, Option_t* processopt="TracksAndParticles",
-                Int_t first = -1,Int_t last = -1,
+void WriteAOD(Option_t* datatype, Int_t first = -1,Int_t last = -1,
+                Option_t* processopt="TracksAndParticles",
                 char *outfile = "AOD.root")
  {
 //datatype defines type of data to be read
index 0a62ca9..fb34393 100644 (file)
@@ -22,9 +22,11 @@ void analysis(Int_t first = -1, Int_t last = -1, const char* directory=".")
      }
    }
 
-  AliReaderESD* reader = new AliReaderESD(dirs);
-  reader->ReadSimulatedData(kTRUE);
-  reader->SetReadMostProbableOnly(kTRUE);
+  AliReaderAOD* reader = new AliReaderAOD("AOD.root");
+  reader->SetDirs(dirs);
+//  AliReaderESD* reader = new AliReaderESD(dirs);
+//  reader->ReadSimulatedData(kTRUE);
+//  reader->SetReadMostProbableOnly(kTRUE);
 
 /*  
    //example PID cuts