Container class for VEvents deriving from VEvent.
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 12 Mar 2009 10:25:21 +0000 (10:25 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 12 Mar 2009 10:25:21 +0000 (10:25 +0000)
Can be used for mixing analysis.

STEER/AliMixedEvent.cxx [new file with mode: 0644]
STEER/AliMixedEvent.h [new file with mode: 0644]
STEER/STEERBaseLinkDef.h
STEER/libSTEERBase.pkg

diff --git a/STEER/AliMixedEvent.cxx b/STEER/AliMixedEvent.cxx
new file mode 100644 (file)
index 0000000..d9ec2af
--- /dev/null
@@ -0,0 +1,106 @@
+/**************************************************************************
+ * Copyright(c) 1998-2007, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/* $Id$ */
+
+
+//-------------------------------------------------------------------------
+//                          Class AliMixedEvent
+// VEvent which is the container of several VEvents 
+// Use Case: Event Mixing     
+// Origin: Andreas Morsch, CERN, Andreas.Morsch@cern.ch 
+//-------------------------------------------------------------------------
+
+
+#include "AliMixedEvent.h"
+#include <TMath.h>
+
+ClassImp(AliMixedEvent)
+
+
+AliMixedEvent::AliMixedEvent() :
+    AliVEvent(),
+    fNEvents(0),       
+    fNumberOfTracks(0),
+    fNTracksCumul(0)
+{
+    // Default constructor
+}
+
+
+AliMixedEvent::AliMixedEvent(const AliMixedEvent& Evnt) :
+    AliVEvent(Evnt) { } // Copy constructor
+
+AliMixedEvent& AliMixedEvent::operator=(const AliMixedEvent& vEvnt)
+{ if (this!=&vEvnt) { 
+    AliVEvent::operator=(vEvnt); 
+  }
+  
+  return *this; 
+}
+
+
+void AliMixedEvent::AddEvent(AliVEvent* evt)
+{
+    // Add a new event to the list
+    fEventList.Add(evt);
+}
+
+
+void AliMixedEvent::Init()
+{
+    // Initialize meta information
+    fNEvents = fEventList.GetEntries();
+    fNTracksCumul = new Int_t[fNEvents];
+    fNumberOfTracks = 0;
+    TIter next(&fEventList);
+    AliVEvent* event;
+    Int_t iev = 0;
+    
+    while((event = (AliVEvent*)next())) {
+       fNTracksCumul[iev++] = fNumberOfTracks;
+       fNumberOfTracks += (event->GetNumberOfTracks());
+    }
+}
+
+
+AliVParticle* AliMixedEvent::GetTrack(Int_t i) const
+{
+    // Return track # i
+    Int_t iEv  = TMath::BinarySearch(fNEvents, fNTracksCumul, i);
+
+    Int_t irel = i - fNTracksCumul[iEv];
+    AliVEvent* evt = (AliVEvent*) (fEventList.At(iEv));
+    return (evt->GetTrack(irel));
+}
+
+
+void AliMixedEvent::Reset()
+{
+    // Reset the event
+    fEventList.Clear();
+    fNEvents = 0;
+    fNumberOfTracks = 0;
+    if (fNTracksCumul) {
+       delete[]  fNTracksCumul;
+       fNTracksCumul = 0;
+    }
+}
+
+Int_t AliMixedEvent::EventIndex(Int_t itrack)
+{
+  // Return the event index for track #itrack
+  return  TMath::BinarySearch(fNEvents, fNTracksCumul, itrack);
+}
diff --git a/STEER/AliMixedEvent.h b/STEER/AliMixedEvent.h
new file mode 100644 (file)
index 0000000..9e933f4
--- /dev/null
@@ -0,0 +1,98 @@
+// -*- mode: C++ -*- 
+#ifndef ALIMIXEDEVENT_H
+#define ALIMIXEDEVENT_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+
+/* $Id$ */
+
+//-------------------------------------------------------------------------
+//                          Class AliMixedEvent
+// VEvent which is the container of several VEvents 
+// Use Case: Event Mixing     
+// Origin: Andreas Morsch, CERN, Andreas.Morsch@cern.ch 
+//-------------------------------------------------------------------------
+
+#include "AliVEvent.h"
+#include <TList.h>
+
+class AliMixedEvent : public AliVEvent {
+
+public:
+    
+    AliMixedEvent();
+    virtual ~AliMixedEvent() {Reset();} 
+    AliMixedEvent(const AliMixedEvent& vEvnt); 
+    AliMixedEvent& operator=(const AliMixedEvent& vEvnt);
+    
+    // Services from VEvent
+    virtual void AddObject(TObject* /*obj*/) {;}
+    virtual TObject* FindListObject(const char* /*name*/) {return 0;}
+    virtual TList* GetList() const {return 0;}
+    virtual void CreateStdContent() {;}
+    virtual void GetStdContent() {;}
+    virtual void ReadFromTree(TTree* /*tree*/, Option_t* /*opt*/) {;} 
+    virtual void WriteToTree(TTree* /*tree*/) const {;}
+    virtual void SetStdNames() {;} 
+    virtual void Print(Option_t * /*option*/) const  {;} 
+    // Specific Services
+    virtual void AddEvent(AliVEvent* evt);
+    virtual void Reset();
+    virtual void Init();
+    
+    // Header
+    virtual AliVHeader* GetHeader() const {return 0;}
+
+    // Delegated methods for fESDRun or AODHeader
+  
+    virtual void     SetRunNumber(Int_t /*n*/)         {;} 
+    virtual void     SetPeriodNumber(UInt_t /*n*/)     {;} 
+    virtual void     SetMagneticField(Double_t /*mf*/) {;} 
+    
+    virtual Int_t    GetRunNumber()     const  {return -999 ;} 
+    virtual UInt_t   GetPeriodNumber()  const  {return    0 ;} 
+    virtual Double_t GetMagneticField() const  {return -999.;} 
+    
+    virtual Double_t GetDiamondX() const {return -999.;}
+    virtual Double_t GetDiamondY() const {return -999.;}
+    virtual void     GetDiamondCovXY(Float_t cov[3]) const
+       {cov[0]=-999.; return;}
+
+  // Delegated methods for fHeader
+    virtual void      SetOrbitNumber(UInt_t /*n*/)        {;} 
+    virtual void      SetBunchCrossNumber(UShort_t /*n*/) {;}
+    virtual void      SetEventType(UInt_t /*eventType*/)  {;}
+    virtual void      SetTriggerMask(ULong64_t /*n*/)     {;}
+    virtual void      SetTriggerCluster(UChar_t /*n*/)    {;}
+
+    virtual UInt_t    GetOrbitNumber()      const  {return    0;}
+    virtual UShort_t  GetBunchCrossNumber() const  {return    0;}
+    virtual UInt_t    GetEventType()        const  {return    0;}
+    virtual ULong64_t GetTriggerMask()      const  {return    0;}
+    virtual UChar_t   GetTriggerCluster()   const  {return    0;}
+    
+    virtual Double_t  GetZDCN1Energy()            const {return -999.;}
+    virtual Double_t  GetZDCP1Energy()            const {return -999.;}
+    virtual Double_t  GetZDCN2Energy()            const {return -999.;}
+    virtual Double_t  GetZDCP2Energy()            const {return -999.;}
+    virtual Double_t  GetZDCEMEnergy(Int_t /*i*/) const {return -999.;}
+  // Tracks
+    virtual AliVParticle *GetTrack(Int_t i)  const;
+    virtual Int_t        GetNumberOfTracks() const {return fNumberOfTracks;}
+    virtual Int_t        GetNumberOfV0s()    const {return -999;}
+    virtual Int_t        EventIndex(Int_t itrack);
+
+  // Primary vertex
+    virtual const AliVVertex   *GetPrimaryVertex() const {return 0;}
+private:
+    TList   fEventList;         //! List of Events
+    Int_t   fNEvents;           //! Number of Events 
+    Int_t   fNumberOfTracks;    //! Total number of tracks
+    Int_t*  fNTracksCumul;      //! Cumulant
+    
+    ClassDef(AliMixedEvent, 0)  // Container for mixed events
+};
+#endif 
+
index 3c9b05f..69752cd 100644 (file)
@@ -18,6 +18,9 @@
 #pragma link C++ class AliVEventHandler+;
 #pragma link C++ class AliVEventPool+;
 
+
+#pragma link C++ class AliMixedEvent+;
+
 #pragma link C++ class AliPID+;
 #pragma link C++ class AliLog+;
 
index fce132d..ccaface 100644 (file)
@@ -4,6 +4,7 @@ SRCS = AliVParticle.cxx \
        AliVTrack.cxx \
        AliVVertex.cxx \
        AliVEvent.cxx \
+       AliMixedEvent.cxx \
        AliVHeader.cxx \
        AliVEventHandler.cxx \
        AliVEventPool.cxx \