Array of track references added to AliMCParticle.
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 10 Oct 2007 14:12:23 +0000 (14:12 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 10 Oct 2007 14:12:23 +0000 (14:12 +0000)
The implementation uses a TRefArrays to objects cashed in AliMCEvent.

STEER/AliMCEvent.cxx
STEER/AliMCEvent.h
STEER/AliMCParticle.cxx
STEER/AliMCParticle.h

index eeec9656f00661033721514d82c317fb7884ce1b..d22440c421b969d6265385541782891c97b365fe 100644 (file)
@@ -26,7 +26,7 @@
 #include <TFile.h>
 #include <TParticle.h>
 #include <TClonesArray.h>
-#include <TObjArray.h>
+#include <TRefArray.h>
 
 #include "AliLog.h"
 #include "AliMCEvent.h"
@@ -42,7 +42,8 @@ AliMCEvent::AliMCEvent():
     fMCParticles(new TClonesArray("AliMCParticle",1000)),
     fMCParticleMap(0),
     fHeader(0),
-    fTrackReferences(0),
+    fTRBuffer(0),
+    fTrackReferences(new TClonesArray("AliTrackReference", 1000)),
     fTreeTR(0),
     fTmpTreeTR(0),
     fTmpFileTR(0),
@@ -58,6 +59,7 @@ AliMCEvent::AliMCEvent(const AliMCEvent& mcEvnt) :
     fMCParticles(0),
     fMCParticleMap(0),
     fHeader(0),
+    fTRBuffer(0),
     fTrackReferences(0),
     fTreeTR(0),
     fTmpTreeTR(0),
@@ -103,8 +105,8 @@ void AliMCEvent::ConnectTreeK (TTree* tree)
        fMCParticleMap->Clear();
        if (fNparticles>0) fMCParticleMap->Expand(fNparticles);}
     else
-       fMCParticleMap = new TObjArray(fNparticles);
-    fMCParticles->Clear();
+       fMCParticleMap = new TRefArray(fNparticles);
+
 }
 
 void AliMCEvent::ConnectTreeTR (TTree* tree)
@@ -121,7 +123,7 @@ void AliMCEvent::ConnectTreeTR (TTree* tree)
        ReorderAndExpandTreeTR();
     } else {
        // New format 
-       fTreeTR->SetBranchAddress("TrackReferences", &fTrackReferences);
+       fTreeTR->SetBranchAddress("TrackReferences", &fTRBuffer);
     }
 }
 
@@ -137,7 +139,7 @@ Int_t AliMCEvent::GetParticleAndTR(Int_t i, TParticle*& particle, TClonesArray*&
     particle = fStack->Particle(i);
     if (fTreeTR) {
        fTreeTR->GetEntry(fStack->TreeKEntry(i));
-       trefs    = fTrackReferences;
+       trefs    = fTRBuffer;
        return trefs->GetEntries();
     } else {
        trefs = 0;
@@ -158,10 +160,10 @@ void AliMCEvent::Clean()
     delete fStack;
 
     // Clear TR
-    if (fTrackReferences) {
-       fTrackReferences->Clear();
-       delete fTrackReferences;
-       fTrackReferences = 0;
+    if (fTRBuffer) {
+       fTRBuffer->Clear();
+       delete fTRBuffer;
+       fTRBuffer = 0;
     }
 }
 
@@ -169,9 +171,12 @@ void AliMCEvent::FinishEvent()
 {
     // Clean-up after event
      fStack->Reset(0);
+     fMCParticles->Clear();
+     fTrackReferences->Clear();
 }
 
 
+
 void AliMCEvent::DrawCheck(Int_t i, Int_t search)
 {
     //
@@ -187,14 +192,14 @@ void AliMCEvent::DrawCheck(Int_t i, Int_t search)
        AliWarning("AliMCEvent::GetEntry: Index out of range");
     }
     
-    Int_t nh = fTrackReferences->GetEntries();
+    Int_t nh = fTRBuffer->GetEntries();
     
     
     if (search) {
        while(nh <= search && i < fNparticles - 1) {
            i++;
            fTreeTR->GetEntry(fStack->TreeKEntry(i));
-           nh =  fTrackReferences->GetEntries();
+           nh =  fTRBuffer->GetEntries();
        }
        printf("Found Hits at %5d\n", i);
     }
@@ -214,7 +219,7 @@ void AliMCEvent::DrawCheck(Int_t i, Int_t search)
     a->Draw();
     
     for (Int_t ih = 0; ih < nh; ih++) {
-       AliTrackReference* ref = (AliTrackReference*) fTrackReferences->At(ih);
+       AliTrackReference* ref = (AliTrackReference*) fTRBuffer->At(ih);
        TMarker* m = new TMarker(ref->X(), ref->Y(), 20);
        m->Draw();
        m->SetMarkerSize(0.4);
@@ -233,8 +238,8 @@ void AliMCEvent::ReorderAndExpandTreeTR()
 
     fTmpFileTR = new TFile("TrackRefsTmp.root", "recreate");
     fTmpTreeTR = new TTree("TreeTR", "TrackReferences");
-    if (!fTrackReferences)  fTrackReferences = new TClonesArray("AliTrackReference", 100);
-    fTmpTreeTR->Branch("TrackReferences", "TClonesArray", &fTrackReferences, 32000, 0);
+    if (!fTRBuffer)  fTRBuffer = new TClonesArray("AliTrackReference", 100);
+    fTmpTreeTR->Branch("TrackReferences", "TClonesArray", &fTRBuffer, 32000, 0);
     
 
 //
@@ -358,14 +363,14 @@ void AliMCEvent::ReorderAndExpandTreeTR()
                        if (label == id) {
                            // secondary found
                            tr->SetDetectorId(ib-1);
-                           Int_t nref =  fTrackReferences->GetEntriesFast();
-                           TClonesArray &lref = *fTrackReferences;
+                           Int_t nref =  fTRBuffer->GetEntriesFast();
+                           TClonesArray &lref = *fTRBuffer;
                            new(lref[nref]) AliTrackReference(*tr);
                        }
                    } // hits
                } // branches
                fTmpTreeTR->Fill();
-               fTrackReferences->Clear();
+               fTRBuffer->Clear();
                ifills++;
            } // daughters
        } // has hits
@@ -393,8 +398,8 @@ void AliMCEvent::ReorderAndExpandTreeTR()
                    
                    if (label == ip) {
                        tr->SetDetectorId(ib-1);
-                       Int_t nref = fTrackReferences->GetEntriesFast();
-                       TClonesArray &lref = *fTrackReferences;
+                       Int_t nref = fTRBuffer->GetEntriesFast();
+                       TClonesArray &lref = *fTRBuffer;
                        new(lref[nref]) AliTrackReference(*tr);
                    }
                } // hits
@@ -402,7 +407,7 @@ void AliMCEvent::ReorderAndExpandTreeTR()
        } // entries
        it++;
        fTmpTreeTR->Fill();
-       fTrackReferences->Clear();
+       fTRBuffer->Clear();
        ifills++;
     } // tracks
     // Check
@@ -430,9 +435,12 @@ void AliMCEvent::ReorderAndExpandTreeTR()
 AliMCParticle* AliMCEvent::GetTrack(Int_t i) const
 {
     // Get MC Particle i
-    AliMCParticle *mcParticle;
-    TParticle     *particle;
-
+    AliMCParticle *mcParticle = 0;
+    TParticle     *particle   = 0;
+    TClonesArray  *trefs      = 0;
+    Int_t          ntref      = 0;
+    TRefArray     *rarray     = 0;
+    
     // Out of range check
     if (i < 0 || i >= fNparticles) {
        AliWarning(Form("AliMCEvent::GetEntry: Index out of range"));
@@ -441,12 +449,31 @@ AliMCParticle* AliMCEvent::GetTrack(Int_t i) const
     }
     
 
-    particle   = fStack->Particle(i);
+
     //
     // First check of the MC Particle has been already cached
     if(!fMCParticleMap->At(i)) {
+       // Get particle from the stack
+       particle   = fStack->Particle(i);
+       // Get track references from Tree TR
+       if (fTreeTR) {
+           fTreeTR->GetEntry(fStack->TreeKEntry(i));
+           trefs     = fTRBuffer;
+           ntref     = trefs->GetEntriesFast();
+           rarray    = new TRefArray(ntref);
+           Int_t nen = fTrackReferences->GetEntriesFast();
+           for (Int_t j = 0; j < ntref; j++) {
+               // Save the track references in a TClonesArray
+               AliTrackReference* ref = dynamic_cast<AliTrackReference*>((*fTRBuffer)[j]);
+               // Save the pointer in a TRefArray
+               new ((*fTrackReferences)[nen]) AliTrackReference(*ref);
+               rarray->AddAt((*fTrackReferences)[nen], j);
+               nen++;
+           } // loop over track references for entry i
+       } // if TreeTR available
+
        Int_t nentries = fMCParticles->GetEntriesFast();
-       new ((*fMCParticles)[nentries]) AliMCParticle(particle);
+       new ((*fMCParticles)[nentries]) AliMCParticle(particle, rarray);
        mcParticle = dynamic_cast<AliMCParticle*>((*fMCParticles)[nentries]);
        fMCParticleMap->AddAt(mcParticle, i);
     } else {
index 124db6a5b9b256451281f55beafaa99b744f41f7..b076edff2376960e6fb93728fc080e2383d48d29 100644 (file)
@@ -15,7 +15,7 @@
 
 
 #include <TTree.h>
-#include <TObjArray.h>
+#include <TRefArray.h>
 
 #include <AliVEvent.h>
 #include "AliVHeader.h"
@@ -107,16 +107,17 @@ private:
     virtual void      ReorderAndExpandTreeTR();
     
 private:
-    AliStack         *fStack;            //! Current pointer to stack
-    TClonesArray     *fMCParticles;      //! Pointer to list of particles
-    TObjArray        *fMCParticleMap;    //! Map of MC Particles
-    AliHeader        *fHeader;           //! Current pointer to header
-    TClonesArray     *fTrackReferences;  //! Current list of track references    
-    TTree            *fTreeTR;           //! Pointer to Track Reference Tree
-    TTree            *fTmpTreeTR;        //! Temporary tree TR to read old format
-    TFile            *fTmpFileTR;        //! Temporary file with TreeTR to read old format
-    Int_t             fNprimaries;       //! Number of primaries
-    Int_t             fNparticles;       //! Number of particles
+    AliStack         *fStack;           //! Current pointer to stack
+    TClonesArray     *fMCParticles;     //! Pointer to list of particles
+    TRefArray        *fMCParticleMap;   //! Map of MC Particles
+    AliHeader        *fHeader;          //! Current pointer to header
+    TClonesArray     *fTRBuffer;        //! Track reference buffer    
+    TClonesArray     *fTrackReferences; //! Array of track references
+    TTree            *fTreeTR;          //! Pointer to Track Reference Tree
+    TTree            *fTmpTreeTR;       //! Temporary tree TR to read old format
+    TFile            *fTmpFileTR;       //! Temporary file with TreeTR to read old format
+    Int_t             fNprimaries;      //! Number of primaries
+    Int_t             fNparticles;      //! Number of particles
     ClassDef(AliMCEvent, 0)  // AliVEvent realisation for MC data
 };
 
index 4bf7f7f0e5f62076b54542ad7b2ed86415c7c240..3f92b34488ca869e9ccd81c91a186c3e2a06b39c 100644 (file)
@@ -21,6 +21,8 @@
 //     Author: Andreas Morsch, CERN
 //-------------------------------------------------------------------------
 
+#include <TRefArray.h>
+
 #include "AliMCParticle.h"
 
 
@@ -28,16 +30,23 @@ ClassImp(AliMCParticle)
 
 AliMCParticle::AliMCParticle():
     AliVParticle(),
-    fParticle(0)
+    fParticle(0),
+    fTrackReferences(0),
+    fNTrackRef(0)
 {
     // Constructor
 }
 
-AliMCParticle::AliMCParticle(TParticle* part):
+AliMCParticle::AliMCParticle(TParticle* part, TRefArray* rarray):
     AliVParticle(),
-    fParticle(part)
+    fParticle(part),
+    fTrackReferences(rarray),
+    fNTrackRef(0)
 {
     // Constructor
+    if (rarray != 0) {
+       fNTrackRef = fTrackReferences->GetEntriesFast();
+    }
 }
     
     
index abe3b4f6ac8449c2f02adf9ee2c35896359b12ba..c47bbcd5396b83dee0f9e5d6acbe71579fc8a6f2 100644 (file)
@@ -7,24 +7,26 @@
 
 //-------------------------------------------------------------------------
 //     AliVParticle realisation for MC Particles
-//     Author: Markus Oldenburg, CERN
+//     Author: Andreas Morsch, CERN
 //-------------------------------------------------------------------------
 
 #include <Rtypes.h>
-#include <AliVParticle.h>
 #include <TParticle.h>
 #include <TParticlePDG.h>
+#include <TRefArray.h>
 
-class AliMCParticle: public AliVParticle {
+#include "AliTrackReference.h"
+#include "AliVParticle.h"
 
+class AliMCParticle: public AliVParticle {
 public:
     AliMCParticle();
-    AliMCParticle(TParticle* part);
+    AliMCParticle(TParticle* part, TRefArray* rarray = 0);
     virtual ~AliMCParticle() {}
     AliMCParticle(const AliMCParticle& mcPart); 
     AliMCParticle& operator=(const AliMCParticle& mcPart);
     
-    // kinematics
+    // Kinematics
     virtual Double_t Px()        const;
     virtual Double_t Py()        const;
     virtual Double_t Pz()        const;
@@ -46,9 +48,15 @@ public:
     
     // PID
     virtual const Double_t *PID() const {return 0;} // return PID object (to be defined, still)
-    
+
+    // Track References
+    Int_t              GetNumberOfTrackReferences() {return fNTrackRef;}
+    AliTrackReference* GetTrackReference(Int_t i)
+       {return dynamic_cast<AliTrackReference*>((*fTrackReferences)[i]);}
  private:
-    TParticle *fParticle; // The TParticle
+    TParticle *fParticle;             // The wrapped TParticle
+    TRefArray *fTrackReferences;      // Reference array to track references
+    Int_t      fNTrackRef;            // Number of track references
     
   ClassDef(AliMCParticle,0)  // AliVParticle realisation for MCParticles
 };