+++ /dev/null
-#include "AliHBTReader.h"
-//_________________________________________________________________________
-///////////////////////////////////////////////////////////////////////////
-//
-// class AliHBTReader
-//
-// Reader Base class (reads particles and tracks and
-// puts it to the AliHBTRun objects
-//
-// Provides functionality for both buffering and non-buffering reading
-// This can be switched on/off via method SetEventBuffering(bool)
-// The main method that inheriting classes need to implement is ReadNext()
-// that read next event in queue.
-// The others are:
-// Bool_t ReadsTracks() const; specifies if reader is able to read simulated particles
-// Bool_t ReadsParticles() const; specifies if reader is able to read reconstructed tracks
-// void Rewind(); rewind reading to the beginning
-//
-// reading of next event is triggered via method Next()
-//
-// This class provides full functionality for reading from many sources
-// User can provide TObjArray of TObjString (SetDirs method or via parameter
-// in constructor) which desribes paths of directories to search data in.
-// If none specified current directory is searched.
-//
-// Piotr.Skowronski@cern.ch
-///////////////////////////////////////////////////////////////////////////
-
-#include <TString.h>
-#include <TObjString.h>
-#include <TObjArray.h>
-#include <TClass.h>
-#include <TRandom.h>
-#include <TH1.h>
-
-#include "AliHBTParticleCut.h"
-#include "AliHBTEvent.h"
-#include "AliHBTRun.h"
-
-ClassImp(AliHBTReader)
-//pure virtual
-
-/*************************************************************************************/
-
-AliHBTReader::AliHBTReader():
- fCuts(new TObjArray()),
- fDirs(0x0),
- fCurrentEvent(0),
- fCurrentDir(0),
- fNEventsRead(0),
- fTracksEvent(0x0),
- fParticlesEvent(0x0),
- fParticles(0x0),
- fTracks(0x0),
- fIsRead(kFALSE),
- fBufferEvents(kFALSE),
- fBlend(kFALSE),
- fFirst(0),
- fLast(0),
- fTrackCounter(0x0)
-{
-//constructor
-}
-/*************************************************************************************/
-
-AliHBTReader::AliHBTReader(TObjArray* dirs):
- fCuts(new TObjArray()),
- fDirs(dirs),
- fCurrentEvent(0),
- fCurrentDir(0),
- fNEventsRead(0),
- fTracksEvent(0x0),
- fParticlesEvent(0x0),
- fParticles(0x0),
- fTracks(0x0),
- fIsRead(kFALSE),
- fBufferEvents(kFALSE),
- fBlend(kFALSE),
- fFirst(0),
- fLast(0),
- fTrackCounter(0x0)
-{
-//ctor with array of directories to read as parameter
-}
-/*************************************************************************************/
-AliHBTReader::AliHBTReader(const AliHBTReader& in):
- TNamed(in),
- fCuts((in.fCuts)?(TObjArray*)in.fCuts->Clone():0x0),
- fDirs((in.fDirs)?(TObjArray*)in.fDirs->Clone():0x0),
- fCurrentEvent(0),
- fCurrentDir(0),
- fNEventsRead(0),
- fTracksEvent(0x0),
- fParticlesEvent(0x0),
- fParticles(0x0),
- fTracks(0x0),
- fIsRead(kFALSE),
- fBufferEvents(in.fBufferEvents),
- fBlend(in.fBlend),
- fFirst(in.fFirst),
- fLast(in.fLast),
- fTrackCounter(0x0)
-{
- //cpy constructor
-}
-
-AliHBTReader::~AliHBTReader()
-{
-//destructor
- if(fCuts)
- {
- fCuts->SetOwner();
- delete fCuts;
- }
- delete fParticlesEvent;
- delete fTracksEvent;
- delete fTrackCounter;
-}
-/*************************************************************************************/
-
-AliHBTReader& AliHBTReader::operator=(const AliHBTReader& in)
-{
- //Assigment operator
- if (this == &in) return *this;
- TNamed::operator=( (const TNamed&)in );
-
- fCuts = (in.fCuts)?(TObjArray*)in.fCuts->Clone():0x0;
- fDirs = (in.fDirs)?(TObjArray*)in.fDirs->Clone():0x0;
- fCurrentEvent = 0;
- fCurrentDir = 0;
- fNEventsRead = 0;
- fTracksEvent = 0x0;
- fParticlesEvent = 0x0;
- fParticles = 0x0;
- fTracks = 0x0;
- fIsRead = kFALSE;
- fBufferEvents = in.fBufferEvents;
- fBlend = in.fBlend;
- fFirst = in.fFirst;
- fLast = in.fLast;
- fTrackCounter = 0x0;
- return *this;
-}
-/*************************************************************************************/
-
-Int_t AliHBTReader::Next()
-{
-//moves to next event
-
- //if asked to read up to event nb. fLast, and it is overcome, report no more events
- if ((fNEventsRead > fLast) && (fLast > 0) ) return kTRUE;
-
- if (fTrackCounter == 0x0)//create Track Counter
- {
- fTrackCounter = new TH1I("trackcounter","Track Counter",20000,0,20000);
- fTrackCounter->SetDirectory(0x0);
- }
-
- do //if asked to read from event fFirst, rewind to it
- {
- if ( ReadNext() == kTRUE) //if no more evets, return it
- return kTRUE;
- }while (fNEventsRead < fFirst);
-
- //here we have event
-
- if (fBlend) Blend();//Mix particles order
-
- if (fBufferEvents)//store events if buffering is on
- {
- if ( ReadsTracks() && fTracksEvent)
- fTracks->SetEvent(fNEventsRead-1-fFirst,fTracksEvent);
- if ( ReadsParticles() && fParticlesEvent)
- fParticles->SetEvent(fNEventsRead-1-fFirst,fParticlesEvent);
- }
- return kFALSE;
-}
-/*************************************************************************************/
-
-void AliHBTReader::AddParticleCut(AliHBTParticleCut* cut)
-{
- //sets the new cut
-
- if (!cut) //if cut is NULL return with error
- {
- Error("AddParticleType","NULL pointers are not accepted any more.\nIf You want to accept all particles of this type, set an empty cut ");
- return;
- }
- AliHBTParticleCut *c = (AliHBTParticleCut*)cut->Clone();
- fCuts->Add(c);
-}
-/********************************************************************/
-
-AliHBTEvent* AliHBTReader::GetParticleEvent(Int_t n)
- {
- //returns Nth event with simulated particles
- if (ReadsParticles() == kFALSE)
- {
- Error("GetParticleEvent","This reader is not able to provide simulated particles.");
- return 0;
- }
-
- if (!fIsRead)
- {
- if (ReadsParticles() && (fParticles == 0x0)) fParticles = new AliHBTRun();
- if (ReadsTracks() && (fTracks == 0x0)) fTracks = new AliHBTRun();
-
- if (Read(fParticles,fTracks))
- {
- Error("GetParticleEvent","Error in reading");
- return 0x0;
- }
- else fIsRead = kTRUE;
- }
- return fParticles->GetEvent(n);
- }
-/********************************************************************/
-
-AliHBTEvent* AliHBTReader::GetTrackEvent(Int_t n)
- {
- //returns Nth event with reconstructed tracks
- if (ReadsTracks() == kFALSE)
- {
- Error("GetTrackEvent","This reader is not able to provide recosntructed tracks.");
- return 0;
- }
- if (!fIsRead)
- {
- if (ReadsParticles() && (fParticles == 0x0)) fParticles = new AliHBTRun();
- if (ReadsTracks() && (fTracks == 0x0)) fTracks = new AliHBTRun();
-
- if(Read(fParticles,fTracks))
- {
- Error("GetTrackEvent","Error in reading");
- return 0x0;
- }
- else fIsRead = kTRUE;
- }
- return fTracks->GetEvent(n);
- }
-/********************************************************************/
-
-Int_t AliHBTReader::GetNumberOfPartEvents()
- {
- //returns number of events of particles
- if (ReadsParticles() == kFALSE)
- {
- Error("GetNumberOfPartEvents","This reader is not able to provide simulated particles.");
- return 0;
- }
-
- if (!fIsRead)
- {
- if (ReadsParticles() && (fParticles == 0x0)) fParticles = new AliHBTRun();
- if (ReadsTracks() && (fTracks == 0x0)) fTracks = new AliHBTRun();
-
- if (Read(fParticles,fTracks))
- {
- Error("GetNumberOfPartEvents","Error in reading");
- return 0;
- }
- else fIsRead = kTRUE;
- }
- return fParticles->GetNumberOfEvents();
- }
-/********************************************************************/
-
-Int_t AliHBTReader::GetNumberOfTrackEvents()
- {
- //returns number of events of tracks
- if (ReadsTracks() == kFALSE)
- {
- Error("GetNumberOfTrackEvents","This reader is not able to provide recosntructed tracks.");
- return 0;
- }
- if (!fIsRead)
- {
- if (ReadsParticles() && (fParticles == 0x0)) fParticles = new AliHBTRun();
- if (ReadsTracks() && (fTracks == 0x0)) fTracks = new AliHBTRun();
-
- if(Read(fParticles,fTracks))
- {
- Error("GetNumberOfTrackEvents","Error in reading");
- return 0;
- }
- else fIsRead = kTRUE;
- }
- return fTracks->GetNumberOfEvents();
- }
-/********************************************************************/
-
-Int_t AliHBTReader::Read(AliHBTRun* particles, AliHBTRun *tracks)
-{
- //reads data and puts put to the particles and tracks objects
- //reurns 0 if everything is OK
- //
- Info("Read","");
-
- if ( ReadsParticles() && (particles == 0x0) ) //check if an object is instatiated
- {
- Error("Read"," particles object must be instatiated before passing it to the reader");
- return 1;
- }
- if ( ReadsTracks() && (tracks == 0x0) ) //check if an object is instatiated
- {
- Error("Read"," tracks object must be instatiated before passing it to the reader");
- return 1;
- }
-
- if (ReadsParticles()) particles->Reset();//clear runs == delete all old events
- if (ReadsTracks()) tracks->Reset();
-
- Rewind();
-
- Int_t i = 0;
- while(Next() == kFALSE)
- {
- if (ReadsTracks()) tracks->SetEvent(i,fTracksEvent);
- if (ReadsParticles()) particles->SetEvent(i,fParticlesEvent);
- i++;
- }
- return 0;
-}
-/*************************************************************************************/
-
-Bool_t AliHBTReader::Rejected(AliHBTParticle* p)
-{
- //Method examines whether particle meets all cut and particle type criteria
-
- if(p==0x0)//of corse we not pass NULL pointers
- {
- Warning("Rejected()","No Pasaran! We never accept NULL pointers");
- return kTRUE;
- }
- //if no particle is specified, we pass all particles
- //excluding NULL pointers, of course
- if ( fCuts->GetEntriesFast() == 0 ) return kFALSE; //if no cut specified accept all particles
- for(Int_t i=0; i<fCuts->GetEntriesFast(); i++)
- {
- AliHBTParticleCut &cut = *((AliHBTParticleCut*)fCuts->At(i));
- if(!cut.Rejected(p)) return kFALSE; //accepted
- }
-
- return kTRUE;//not accepted
-}
-/*************************************************************************************/
-
-Bool_t AliHBTReader::Rejected(Int_t pid)
-{
-//this method checks if any of existing cuts accepts this pid particles
-//or any cuts accepts all particles
-
- if(pid == 0)
- return kTRUE;
-
- if ( fCuts->GetEntriesFast() == 0 ) return kFALSE; //if no cut specified accept all particles
-
- for(Int_t i=0; i<fCuts->GetEntriesFast(); i++)
- {
- AliHBTParticleCut &cut = *((AliHBTParticleCut*)fCuts->At(i));
- //if some of cuts accepts all particles or some accepts particles of this type, accept
- if ( (cut.GetPID() == 0) || (cut.GetPID() == pid) ) return kFALSE;
- }
- return kTRUE;
-}
-/*************************************************************************************/
-
-TString& AliHBTReader::GetDirName(Int_t entry)
-{
-//returns directory name of next one to read
- TString* retval;//return value
- if (fDirs == 0x0)
- {
- retval = new TString(".");
- 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;
- }
-
- if (fDirs->GetEntries() == 0)
- {
- retval = new TString(".");
- return *retval;
- }
-
- TClass *objclass = fDirs->At(entry)->IsA();
- TClass *stringclass = TObjString::Class();
-
- TObjString *dir = (TObjString*)objclass->DynamicCast(stringclass,fDirs->At(entry));
-
- if(dir == 0x0)
- {
- Error("GetDirName","Object in TObjArray is not a TObjString or its descendant");
- retval = new TString();
- return *retval;
- }
- if (gDebug > 0) Info("GetDirName","Returned ok %s",dir->String().Data());
- return dir->String();
-}
-/*************************************************************************************/
-
-void AliHBTReader::Blend()
-{
- //randomly change positions of the particles after reading
- //is used to check if some distr depends on order of particles
- //(tracking gives particles Pt sorted)
-
- if (fParticlesEvent == 0x0) return;
-
- for (Int_t i = 2; i < fParticlesEvent->GetNumberOfParticles(); i++)
- {
- Int_t with = gRandom->Integer(i);
- fParticlesEvent->SwapParticles(i,with);
- if (fTracksEvent) fTracksEvent->SwapParticles(i,with);
- }
-}
-/*************************************************************************************/
-
-void AliHBTReader::WriteTrackCounter() const
-{
- //writes the counter histogram
-
- if (fTrackCounter) fTrackCounter->Write(0,TObject::kOverwrite);
- else
- {
- Warning("WriteTrackCounter","Counter is NULL");
- }
-}
+++ /dev/null
-#ifndef ALIHBTREADER_H
-#define ALIHBTREADER_H
-//_________________________________________________________________________
-///////////////////////////////////////////////////////////////////////////
-//
-// class AliHBTReader
-//
-// Reader Base class (reads particles and tracks and
-// puts it to the AliHBTRun objects
-//
-// Piotr.Skowronski@cern.ch
-///////////////////////////////////////////////////////////////////////////
-
-#include <TNamed.h>
-#include <TObjArray.h>
-
-class AliHBTRun;
-class AliHBTEvent;
-class AliHBTParticleCut;
-class AliHBTParticle;
-class TString;
-class TH1I;
-
-class AliHBTReader: public TNamed
-{
- public:
- AliHBTReader();
- AliHBTReader(TObjArray*);
- AliHBTReader(const AliHBTReader& in);
- virtual ~AliHBTReader();
-
- AliHBTReader& operator=(const AliHBTReader& in);
-
- virtual Int_t Next();
- virtual void Rewind() = 0; //
-
- virtual Bool_t ReadsTracks() const = 0; //specifies if reader is able to read simulated particles
- virtual Bool_t ReadsParticles() const = 0;//specifies if reader is able to read reconstructed tracks
-
- void AddParticleCut(AliHBTParticleCut* cut);
-
- virtual Int_t Read(AliHBTRun* particles, AliHBTRun *tracks);
-
- virtual AliHBTEvent* GetParticleEvent() {return fParticlesEvent;}//can not be const because position randomizer overloads it
- virtual AliHBTEvent* GetTrackEvent() {return fTracksEvent;}//
-
- virtual AliHBTEvent* GetParticleEvent(Int_t n);
- virtual AliHBTEvent* GetTrackEvent(Int_t n);
-
- virtual Int_t GetNumberOfPartEvents();
- virtual Int_t GetNumberOfTrackEvents();
-
- void SetDirs(TObjArray* dirs){fDirs = dirs;} //sets array directories names
- void SetEventBuffering(Bool_t flag){fBufferEvents = flag;}//switches on/off buffering - read data are kept in local buffer
- void SetBlend(Bool_t flag = kTRUE){fBlend=flag;} //set blending - randomizing particle order
- virtual Int_t GetNumberOfDirs() const {return (fDirs)?fDirs->GetEntries():0;}
- void ReadEventsFromTo(Int_t first,Int_t last){fFirst = first; fLast = last;}
- virtual TH1I* GetTrackCounter() const {return fTrackCounter;}
- virtual void WriteTrackCounter() const;
-
- protected:
-
- TObjArray* fCuts;//array with particle cuts
- TObjArray* fDirs;//arry with directories to read data from
-
- Int_t fCurrentEvent;//! number of current event in current directory
- Int_t fCurrentDir;//! number of current directory
-
- Int_t fNEventsRead;//!total
-
- AliHBTEvent* fTracksEvent; //! tracks read from current event
- AliHBTEvent* fParticlesEvent; //! particles read from current event
-
- AliHBTRun* fParticles; //!simulated particles
- AliHBTRun* fTracks; //!reconstructed tracks (particles)
-
- Bool_t fIsRead;//!flag indicating if the data are already read
- Bool_t fBufferEvents;//flag indicating if the data should be bufferred
-
- Bool_t fBlend;// flag indicating if randomly change positions of the particles after reading
-
- Int_t fFirst;//first event to return (all are before are skipped)
- Int_t fLast;//last
-
- TH1I* fTrackCounter; //histogram with number of tracks read
-
- virtual Int_t ReadNext() = 0; //this methods reads next event and put result in fTracksEvent and/or fParticlesEvent
- Bool_t Rejected(AliHBTParticle* p);
- Bool_t Rejected(Int_t pid);
- void Blend();
-
- TString& GetDirName(Int_t entry);
-
- private:
-
- ClassDef(AliHBTReader,4)//version 2 - TNamed as parental class
- //version 3 - Blending added
-};
-
-#endif
+++ /dev/null
-#include "AliHBTReaderInternal.h"
-//_________________________________________________________________________
-///////////////////////////////////////////////////////////////////////////
-// //
-// class AliHBTReaderInternal //
-// //
-// Multi file reader for Internal Data Format //
-// //
-// This reader reads data created by itself //
-// (method AliHBTReaderInternal::Write) //
-// Data are stored in form of tree of TClonesArray of AliHBTParticle's //
-// //
-// Piotr.Skowronski@cern.ch //
-// more info: http://aliweb.cern.ch/people/skowron/analyzer/index.html //
-// //
-///////////////////////////////////////////////////////////////////////////
-
-#include <TTree.h>
-#include <TFile.h>
-#include <TParticle.h>
-#include <TError.h>
-#include <AliRun.h>
-#include <AliMagF.h>
-
-#include "AliHBTRun.h"
-#include "AliHBTEvent.h"
-#include "AliHBTParticle.h"
-#include "AliHBTParticleCut.h"
-
-ClassImp(AliHBTReaderInternal)
-/********************************************************************/
-
-AliHBTReaderInternal::AliHBTReaderInternal():
- fFileName(),
- fPartBranch(0x0),
- fTrackBranch(0x0),
- fTree(0x0),
- fFile(0x0),
- fPartBuffer(0x0),
- fTrackBuffer(0x0)
-{
-//Defalut constructor
-}
-/********************************************************************/
-
-AliHBTReaderInternal::AliHBTReaderInternal(const char *filename):
- fFileName(filename),
- fPartBranch(0x0),
- fTrackBranch(0x0),
- fTree(0x0),
- fFile(0x0),
- fPartBuffer(0x0),
- fTrackBuffer(0x0)
-{
-//constructor
-//filename - name of file to open
-}
-/********************************************************************/
-
-AliHBTReaderInternal::AliHBTReaderInternal(TObjArray* dirs, const char *filename):
- AliHBTReader(dirs),
- fFileName(filename),
- fPartBranch(0x0),
- fTrackBranch(0x0),
- fTree(0x0),
- fFile(0x0),
- fPartBuffer(0x0),
- fTrackBuffer(0x0)
-{
-//ctor
-//dirs contains strings with directories to look data in
-//filename - name of file to open
-}
-/********************************************************************/
-
-AliHBTReaderInternal::AliHBTReaderInternal(const AliHBTReaderInternal& in):
- AliHBTReader(in),
- fFileName(in.fFileName),
- fPartBranch(0x0),
- fTrackBranch(0x0),
- fTree(0x0),
- fFile(0x0),
- fPartBuffer(0x0),
- fTrackBuffer(0x0)
-{
- //cpy constructor
-}
-/********************************************************************/
-
-AliHBTReaderInternal& AliHBTReaderInternal::operator=(const AliHBTReaderInternal& in)
-{
- //Assigment operator
- if (this == &in) return *this;
- Rewind();//close current session
- AliHBTReader::operator=((const AliHBTReader&)in);
- fFileName = in.fFileName;
- return *this;
-}
-/********************************************************************/
-AliHBTReaderInternal::~AliHBTReaderInternal()
- {
- //desctructor
- delete fTree;
- delete fFile;
- }
-/********************************************************************/
-
-void AliHBTReaderInternal::Rewind()
-{
- delete fTree;
- fTree = 0x0;
- delete fFile;
- fFile = 0x0;
- fCurrentDir = 0;
- fNEventsRead= 0;
-}
-/********************************************************************/
-
-Int_t AliHBTReaderInternal::ReadNext()
- {
- //reads data and puts put to the particles and tracks objects
- //reurns 0 if everything is OK
- //
- Info("ReadNext","");
- Int_t i; //iterator and some temprary values
- Int_t counter;
- AliHBTParticle* tpart = 0x0, *ttrack = 0x0;
-
- TDatabasePDG* pdgdb = TDatabasePDG::Instance();
- if (pdgdb == 0x0)
- {
- Error("ReadNext","Can not get PDG Particles Data Base");
- return 1;
- }
-
- if (fParticlesEvent == 0x0) fParticlesEvent = new AliHBTEvent();
- if (fTracksEvent == 0x0) fTracksEvent = new AliHBTEvent();
-
- fParticlesEvent->Reset();
- fTracksEvent->Reset();
-
- do //do{}while; is OK even if 0 dirs specified. In that case we try to read from "./"
- {
- if (fTree == 0x0)
- if( (i=OpenNextFile()) )
- {
- Error("ReadNext","Skipping directory due to problems with opening files. Errorcode %d",i);
- fCurrentDir++;
- continue;
- }
- if (fCurrentEvent == (Int_t)fTree->GetEntries())
- {
- delete fTree;
- fTree = 0x0;
- delete fFile;
- fFile = 0x0;
- fPartBranch = 0x0;
- fTrackBranch= 0x0;
- fCurrentDir++;
- continue;
- }
- /***************************/
- /***************************/
- /***************************/
-
-
- Info("ReadNext","Event %d",fCurrentEvent);
- fTree->GetEvent(fCurrentEvent);
-
- counter = 0;
- if (fPartBranch && fTrackBranch)
- {
- Info("ReadNext","Found %d tracks in total.",fTrackBuffer->GetEntries());
- Info("ReadNext","Found %d particles in total.",fPartBuffer->GetEntries());
- for(i = 0; i < fPartBuffer->GetEntries(); i++)
- {
- tpart = dynamic_cast<AliHBTParticle*>(fPartBuffer->At(i));
- ttrack = dynamic_cast<AliHBTParticle*>(fTrackBuffer->At(i));
-
- if( tpart == 0x0 ) continue; //if returned pointer is NULL
- if( ttrack == 0x0 ) continue; //if returned pointer is NULL
-
- if (AliHBTParticle::GetDebug() > 9)
- {
- Info("ReadNext","Particle:");
- tpart->Print();
- Info("ReadNext","Track:");
- ttrack->Print();
- }
- if (ttrack->GetUID() != tpart->GetUID())
- {
- Error("ReadNext","Sth. is wrong: Track and Particle has different UID.");
- Error("ReadNext","They probobly do not correspond to each other.");
- }
-
- for (Int_t s = 0; s < ttrack->GetNumberOfPids(); s++)
- {
- //check if we are intersted with particles of this type
- //if not take next partilce
- if( Rejected(ttrack->GetNthPid(s)) )
- {
- if (AliHBTParticle::GetDebug() > 9)
- Info("ReadNext","Track Incarnation %d did not pass PID cut.",ttrack->GetNthPid(s));
- continue;
- }
- TParticlePDG* pdgp = pdgdb->GetParticle(ttrack->GetNthPid(s));
- if (pdgp == 0x0)//PDG part corresponding to new incarnation
- {
- Error("ReadNext","Particle code unknown to PDG DB.");
- continue;
- }
-
- AliHBTParticle* track = new AliHBTParticle(*ttrack);
-
- //apart of setting PDG code of an incarnation
- //it is necessary tu recalculate energy on the basis of
- //new PDG code (mass) hypothesis
- Double_t mass = pdgp->Mass();//mass of new incarnation
- Double_t tEtot = TMath::Sqrt( ttrack->Px()*ttrack->Px() +
- ttrack->Py()*ttrack->Py() +
- ttrack->Pz()*ttrack->Pz() +
- mass*mass);//total energy of the new incarnation
- //update Energy and Calc Mass
- track->SetMomentum(ttrack->Px(),ttrack->Py(),ttrack->Pz(),tEtot);
- track->SetCalcMass(mass);
- track->SetPdgCode(ttrack->GetNthPid(s),ttrack->GetNthPidProb(s));
-
- if( Rejected(track) )
- {
- if (AliHBTParticle::GetDebug() > 9)
- Info("ReadNext","Track Incarnation %d did not pass cut.",ttrack->GetNthPid(s));
- delete track;
- continue;
- }
- AliHBTParticle* part = new AliHBTParticle(*tpart);//particle has only 1 incarnation (real)
-
- fParticlesEvent->AddParticle(part);
- fTracksEvent->AddParticle(track);
-
- counter++;
- }
- }
- Info("ReadNext"," Read: %d particles and tracks.",counter);
- }
- else
- {
- if (fPartBranch)
- {
- Info("ReadNext","Found %d particles in total.",fPartBuffer->GetEntries());
- for(i = 0; i < fPartBuffer->GetEntries(); i++)
- {
- tpart = dynamic_cast<AliHBTParticle*>(fPartBuffer->At(i));
- if(tpart == 0x0) continue; //if returned pointer is NULL
-
- for (Int_t s = 0; s < tpart->GetNumberOfPids(); s++)
- {
- if( pdgdb->GetParticle(tpart->GetNthPid(s)) == 0x0 ) continue; //if particle has crazy PDG code (not known to our database)
- if( Rejected(tpart->GetNthPid(s)) ) continue; //check if we are intersted with particles of this type
- AliHBTParticle* part = new AliHBTParticle(*tpart);
- part->SetPdgCode(tpart->GetNthPid(s),tpart->GetNthPidProb(s));
- if( Rejected(part) )
- {
- delete part;
- continue;
- }
- fParticlesEvent->AddParticle(part);
- counter++;
- }
- }
- Info("ReadNext"," Read: %d particles.",counter);
- }
- else Info("ReadNext"," Read: 0 particles.");
-
- if (fTrackBranch)
- {
- Info("ReadNext","Found %d tracks in total.",fTrackBuffer->GetEntries());
- for(i = 0; i < fTrackBuffer->GetEntries(); i++)
- {
- ttrack = dynamic_cast<AliHBTParticle*>(fTrackBuffer->At(i));
- if( ttrack == 0x0 ) continue; //if returned pointer is NULL
-
- for (Int_t s = 0; s < ttrack->GetNumberOfPids(); s++)
- {
- if( Rejected(ttrack->GetNthPid(s)) ) continue; //check if we are intersted with particles of this type
- //if not take next partilce
- TParticlePDG* pdgp = pdgdb->GetParticle(ttrack->GetNthPid(s));
- if (pdgp == 0x0)//PDG part corresponding to new incarnation
- {
- Error("ReadNext","Particle code unknown to PDG DB.");
- continue;
- }
- AliHBTParticle* track = new AliHBTParticle(*ttrack);
-
- //apart of setting PDG code of an incarnation
- //it is necessary tu recalculate energy on the basis of
- //new PDG code (mass) hypothesis
- Double_t mass = pdgp->Mass();//mass of new incarnation
- Double_t tEtot = TMath::Sqrt( ttrack->Px()*ttrack->Px() +
- ttrack->Py()*ttrack->Py() +
- ttrack->Pz()*ttrack->Pz() +
- mass*mass);//total energy of the new incarnation
- //update Energy and Calc Mass
- track->SetMomentum(ttrack->Px(),ttrack->Py(),ttrack->Pz(),tEtot);
- track->SetCalcMass(mass);
- track->SetPdgCode(ttrack->GetNthPid(s),ttrack->GetNthPidProb(s));
-
- if( Rejected(track) )
- {
- delete track;
- continue;
- }
- fTracksEvent->AddParticle(track);
-
- counter++;
- }
- }
- Info("ReadNext"," Read: %d tracks",counter);
- }
- else Info("ReadNext"," Read: 0 tracks.");
- }
-
- if (fPartBuffer) fPartBuffer->Delete();
- if (fTrackBuffer) fTrackBuffer->Delete();
-
- fCurrentEvent++;
- fNEventsRead++;
-
- return 0;
-
- }while(fCurrentDir < GetNumberOfDirs());
-
- return 1;//no more directories to read
-}
-/********************************************************************/
-
-Int_t AliHBTReaderInternal::OpenNextFile()
-{
- //open file in current directory
- const TString& dirname = GetDirName(fCurrentDir);
- if (dirname == "")
- {
- Error("OpenNextFile","Can not get directory name");
- return 4;
- }
-
- TString filename = dirname +"/"+ fFileName;
- fFile = TFile::Open(filename.Data());
- if ( fFile == 0x0 )
- {
- Error("OpenNextFile","Can't open file named %s",filename.Data());
- return 1;
- }
- if (fFile->IsOpen() == kFALSE)
- {
- Error("OpenNextFile","Can't open filenamed %s",filename.Data());
- return 1;
- }
-
- fTree = (TTree*)fFile->Get("data");
- if (fTree == 0x0)
- {
- Error("OpenNextFile","Can not get the tree.");
- return 1;
- }
-
-
- fTrackBranch = fTree->GetBranch("tracks");//get the branch with tracks
- if (fTrackBranch == 0x0) ////check if we got the branch
- {//if not return with error
- Info("OpenNextFile","Can't find a branch with tracks !\n");
- }
- else
- {
- if (fTrackBuffer == 0x0)
- fTrackBuffer = new TClonesArray("AliHBTParticle",15000);
- fTrackBranch->SetAddress(&fTrackBuffer);
- }
-
- fPartBranch = fTree->GetBranch("particles");//get the branch with particles
- if (fPartBranch == 0x0) ////check if we got the branch
- {//if not return with error
- Info("OpenNextFile","Can't find a branch with particles !\n");
- }
- else
- {
- if (fPartBuffer == 0x0)
- fPartBuffer = new TClonesArray("AliHBTParticle",15000);
- fPartBranch->SetAddress(&fPartBuffer);
- }
-
- Info("OpenNextFile","________________________________________________________");
- Info("OpenNextFile","Found %d event(s) in directory %s",
- (Int_t)fTree->GetEntries(),GetDirName(fCurrentDir).Data());
-
- fCurrentEvent = 0;
-
- return 0;
-}
-/********************************************************************/
-
-Int_t AliHBTReaderInternal::Write(AliHBTReader* reader,const char* outfile, Bool_t multcheck)
- {
- //reads tracks from reader and writes runs to file
- //reader - provides data for writing in internal format
- //name of output file
- //multcheck - switches of checking if particle was stored with other incarnation
- // usefull e.g. when using kine data, where all particles have 100% pid prob.and saves a lot of time
-
- Int_t i,j;
-
- ::Info("AliHBTReaderInternal::Write","________________________________________________________");
- ::Info("AliHBTReaderInternal::Write","________________________________________________________");
- ::Info("AliHBTReaderInternal::Write","________________________________________________________");
-
- TFile *histoOutput = TFile::Open(outfile,"recreate");
-
- if (!histoOutput->IsOpen())
- {
- ::Error("AliHBTReaderInternal::Write","File is not opened");
- return 1;
- }
-
- TTree *tracktree = new TTree("data","Tree with tracks");
-
- TClonesArray* pbuffer = new TClonesArray("AliHBTParticle",15000);
- TClonesArray* tbuffer = new TClonesArray("AliHBTParticle",15000);
-
- TClonesArray &particles = *pbuffer;
- TClonesArray &tracks = *tbuffer;
-
- TString name("Tracks");
-
- Int_t nt = reader->GetNumberOfTrackEvents();
- Int_t np = reader->GetNumberOfPartEvents();
-
- if (AliHBTParticle::GetDebug() > 0)
- ::Info("Write","Reader has %d track events and %d particles events.",nt,np);
-
- Bool_t trck = (nt > 0) ? kTRUE : kFALSE;
- Bool_t part = (np > 0) ? kTRUE : kFALSE;
-
- TBranch *trackbranch = 0x0, *partbranch = 0x0;
-
- if (trck) trackbranch = tracktree->Branch("tracks",&tbuffer,32000,0);
- if (part) partbranch = tracktree->Branch("particles",&pbuffer,32000,0);
-
- if ( (trck) && (part) && (np != nt))
- {
- ::Warning("AliHBTReaderInternal::Write","Number of track and particle events is different");
- }
-
- Int_t n;
- if (nt >= np ) n = nt; else n = np;
-
- if (AliHBTParticle::GetDebug() > 0)
- ::Info("Write","Will loop over %d events",n);
-
- for ( i =0;i< n; i++)
- {
- ::Info("AliHBTReaderInternal::Write","Event %d",i+1);
- Int_t counter = 0;
- if (trck && (i<=nt))
- {
- AliHBTEvent* trackev = reader->GetTrackEvent(i);
- for ( j = 0; j< trackev->GetNumberOfParticles();j++)
- {
- const AliHBTParticle& t = *(trackev->GetParticle(j));
- if (multcheck)
- {
- if (FindIndex(tbuffer,t.GetUID()))
- {
- if (AliHBTParticle::GetDebug()>4)
- {
- ::Info("Write","Track with Event UID %d already stored",t.GetUID());
- }
- continue; //not to write the same particles with other incarnations
- }
- }
- new (tracks[counter++]) AliHBTParticle(t);
- }
- ::Info("AliHBTReaderInternal::Write"," Tracks: %d",tracks.GetEntries());
- }else ::Info("AliHBTReaderInternal::Write","NO TRACKS");
-
- counter = 0;
- if (part && (i<=np))
- {
-// ::Warning("AliHBTReaderInternal::Write","Find index switched off!!!");
-
- AliHBTEvent* partev = reader->GetParticleEvent(i);
- for ( j = 0; j< partev->GetNumberOfParticles();j++)
- {
- const AliHBTParticle& part= *(partev->GetParticle(j));
- if (multcheck)
- {
- if (FindIndex(pbuffer,part.GetUID()))
- {
- if (AliHBTParticle::GetDebug()>4)
- {
- ::Info("Write","Particle with Event UID %d already stored",part.GetUID());
- }
- continue; //not to write the same particles with other incarnations
- }
- }
- new (particles[counter++]) AliHBTParticle(part);
- }
- ::Info("AliHBTReaderInternal::Write"," Particles: %d",particles.GetEntries());
- }else ::Info("AliHBTReaderInternal::Write","NO PARTICLES");
-
- histoOutput->cd();
- tracktree->Fill();
- tracktree->AutoSave();
- tbuffer->Delete();
- pbuffer->Delete();
- }
-
- histoOutput->cd();
- tracktree->Write(0,TObject::kOverwrite);
- delete tracktree;
-
- tbuffer->SetOwner();
- pbuffer->SetOwner();
- delete pbuffer;
- delete tbuffer;
-
- histoOutput->Close();
- return 0;
- }
-/********************************************************************/
-
-Bool_t AliHBTReaderInternal::FindIndex(TClonesArray* arr,Int_t idx)
-{
-//Checks if in the array exists already partilce with Unique ID idx
- if (arr == 0x0)
- {
- ::Error("FindIndex","Array is 0x0");
- return kTRUE;
- }
- TIter next(arr);
- AliHBTParticle* p;
- while (( p = (AliHBTParticle*)next()))
- {
- if (p->GetUID() == idx) return kTRUE;
- }
- return kFALSE;
-}
+++ /dev/null
-#ifndef ALIHBTREADERINTERNAL_H
-#define ALIHBTREADERINTERNAL_H
-//_________________________________________________________________________
-///////////////////////////////////////////////////////////////////////////
-// //
-// class AliHBTReaderInternal //
-// //
-// Multi file reader for Internal Data Format //
-// //
-// This reader reads data created by itself //
-// (method AliHBTReaderInternal::Write) //
-// //
-// Piotr.Skowronski@cern.ch //
-// more info: http://alisoft.cern.ch/people/skowron/analyzer/index.html //
-// //
-///////////////////////////////////////////////////////////////////////////
-
-#include "AliHBTReader.h"
-#include <TString.h>
-
-class TFile;
-class TTree;
-class TBranch;
-class TClonesArray;
-
-class AliHBTReaderInternal: public AliHBTReader
-{
- public:
- AliHBTReaderInternal();
- AliHBTReaderInternal(const char *filename);
- AliHBTReaderInternal(TObjArray* dirs, const char *filename);
- AliHBTReaderInternal(const AliHBTReaderInternal& in);
-
- virtual ~AliHBTReaderInternal();
-
- AliHBTReaderInternal& operator=(const AliHBTReaderInternal& in);
-
- void Rewind();
-
- Bool_t ReadsTracks() const {return kTRUE;}
- Bool_t ReadsParticles() const {return kTRUE;}
-
- static Int_t Write(AliHBTReader* reader,const char* outfile = "data.root", Bool_t multcheck = kFALSE);//reads tracks from runs and writes them to file
-
- protected:
-
- TString fFileName;//name of the file with tracks
- TBranch* fPartBranch;//!branch with particles
- TBranch* fTrackBranch;//!branch with tracks
- TTree* fTree;//!tree
- TFile* fFile;//!file
- TClonesArray* fPartBuffer;//!buffer array that tree is read to
- TClonesArray* fTrackBuffer;//!
-
- Int_t ReadNext();//reads next event
- Int_t OpenNextFile();//opens file to be read for given event
- static Bool_t FindIndex(TClonesArray* arr,Int_t idx);
-
- ClassDef(AliHBTReaderInternal,2)
-};
-#endif
+++ /dev/null
-#include "AliHBTReaderTPC.h"
-//______________________________________________
-//
-// class AliHBTReaderTPC
-//
-// reader for TPC tracks
-// needs galice.root
-// just to shut up coding conventions checker
-//
-// more info: http://aliweb.cern.ch/people/skowron/analyzer/index.html
-// Piotr.Skowronski@cern.ch
-//
-///////////////////////////////////////////////////////////////////////////
-#include <TTree.h>
-#include <TParticle.h>
-#include <TH1.h>
-
-
-#include <AliRun.h>
-#include <AliLoader.h>
-#include <AliStack.h>
-#include <AliMagF.h>
-#include <AliTPCtrack.h>
-#include <AliTPCLoader.h>
-
-#include "AliHBTEvent.h"
-#include "AliHBTParticle.h"
-#include "AliHBTTrackPoints.h"
-#include "AliHBTClusterMap.h"
-
-
-ClassImp(AliHBTReaderTPC)
-
-AliHBTReaderTPC::AliHBTReaderTPC():
- fFileName("galice.root"),
- fRunLoader(0x0),
- fTPCLoader(0x0),
- fMagneticField(0.0),
- fUseMagFFromRun(kTRUE),
- fNTrackPoints(0),
- fdR(0.0),
- fClusterMap(kFALSE),
- fNClustMin(0),
- fNClustMax(150),
- fNChi2PerClustMin(0.0),
- fNChi2PerClustMax(10e5),
- fC00Min(0.0),
- fC00Max(10e5),
- fC11Min(0.0),
- fC11Max(10e5),
- fC22Min(0.0),
- fC22Max(10e5),
- fC33Min(0.0),
- fC33Max(10e5),
- fC44Min(0.0),
- fC44Max(10e5)
-{
- //constructor
-
-}
-/********************************************************************/
-
-AliHBTReaderTPC::AliHBTReaderTPC(const Char_t* galicefilename):
- fFileName(galicefilename),
- fRunLoader(0x0),
- fTPCLoader(0x0),
- fMagneticField(0.0),
- fUseMagFFromRun(kTRUE),
- fNTrackPoints(0),
- fdR(0.0),
- fNClustMin(0),
- fNClustMax(150),
- fNChi2PerClustMin(0.0),
- fNChi2PerClustMax(10e5),
- fC00Min(0.0),
- fC00Max(10e5),
- fC11Min(0.0),
- fC11Max(10e5),
- fC22Min(0.0),
- fC22Max(10e5),
- fC33Min(0.0),
- fC33Max(10e5),
- fC44Min(0.0),
- fC44Max(10e5)
-{
- //constructor,
- //Defaults:
-}
-/********************************************************************/
-
-AliHBTReaderTPC::AliHBTReaderTPC(TObjArray* dirs, const Char_t* galicefilename):
- AliHBTReader(dirs),
- fFileName(galicefilename),
- fRunLoader(0x0),
- fTPCLoader(0x0),
- fMagneticField(0.0),
- fUseMagFFromRun(kTRUE),
- fNTrackPoints(0),
- fdR(0.0),
- fClusterMap(kFALSE),
- fNClustMin(0),
- fNClustMax(150),
- fNChi2PerClustMin(0.0),
- fNChi2PerClustMax(10e5),
- fC00Min(0.0),
- fC00Max(10e5),
- fC11Min(0.0),
- fC11Max(10e5),
- fC22Min(0.0),
- fC22Max(10e5),
- fC33Min(0.0),
- fC33Max(10e5),
- fC44Min(0.0),
- fC44Max(10e5)
-{
- //constructor,
- //Defaults:
- // galicefilename = "" - this means: Do not open gAlice file -
- // just leave the global pointer untached
-
-}
-/********************************************************************/
-AliHBTReaderTPC::AliHBTReaderTPC(const AliHBTReaderTPC& in):
- AliHBTReader(in),
- fFileName(in.fFileName),
- fRunLoader(0x0),
- fTPCLoader(0x0),
- fMagneticField(in.fMagneticField),
- fUseMagFFromRun(in.fUseMagFFromRun),
- fNTrackPoints(in.fNTrackPoints),
- fdR(in.fdR),
- fNClustMin(in.fNClustMin),
- fNClustMax(in.fNClustMax),
- fNChi2PerClustMin(in.fNChi2PerClustMin),
- fNChi2PerClustMax(in.fNChi2PerClustMax),
- fC00Min(in.fC00Min),
- fC00Max(in.fC00Max),
- fC11Min(in.fC11Min),
- fC11Max(in.fC11Max),
- fC22Min(in.fC22Min),
- fC22Max(in.fC22Max),
- fC33Min(in.fC33Min),
- fC33Max(in.fC33Max),
- fC44Min(in.fC44Min),
- fC44Max(in.fC44Max)
-{
- //cpy constructor,
-}
-/********************************************************************/
-
-AliHBTReaderTPC::~AliHBTReaderTPC()
-{
- //desctructor
- if (AliHBTParticle::GetDebug())
- {
- Info("~AliHBTReaderTPC","deleting Run Loader");
- AliLoader::SetDebug(AliHBTParticle::GetDebug());
- }
-
- delete fRunLoader;
-
- if (AliHBTParticle::GetDebug())
- {
- Info("~AliHBTReaderTPC","deleting Run Loader Done");
- }
-}
-/********************************************************************/
-
-AliHBTReaderTPC& AliHBTReaderTPC::operator=(const AliHBTReaderTPC& in)
-{
-//Assigment operator
-
- delete fRunLoader;
-
- fFileName = in.fFileName;
- fRunLoader = 0x0;
- fTPCLoader = 0x0;
- fMagneticField = in.fMagneticField;
- fUseMagFFromRun = in.fUseMagFFromRun;
- fNTrackPoints = in.fNTrackPoints;
- fdR = in.fdR;
- fNClustMin = in.fNClustMin;
- fNClustMax = in.fNClustMax;
- fNChi2PerClustMin = in.fNChi2PerClustMin;
- fNChi2PerClustMax = in.fNChi2PerClustMax;
- fC00Min = in.fC00Min;
- fC00Max = in.fC00Max;
- fC11Min = in.fC11Min;
- fC11Max = in.fC11Max;
- fC22Min = in.fC22Min;
- fC22Max = in.fC22Max;
- fC33Min = in.fC33Min;
- fC33Max = in.fC33Max;
- fC44Min = in.fC44Min;
- fC44Max = in.fC44Max;
- return *this;
-}
-/********************************************************************/
-
-void AliHBTReaderTPC::Rewind()
-{
-//Rewind reading to the beginning
- delete fRunLoader;
- fRunLoader = 0x0;
- fCurrentDir = 0;
- fNEventsRead= 0;
- if (fTrackCounter) fTrackCounter->Reset();
-}
-/********************************************************************/
-
-Int_t AliHBTReaderTPC::ReadNext()
- {
- //reads data and puts put to the particles and tracks objects
- //reurns 0 if everything is OK
- //
- Info("Read","");
-
- TObjArray *tarray = new TObjArray(5000); //cotainer for tpc tracks
- tarray->SetOwner(); //set the ownership of the objects it contains
- //when array is is deleted or cleared all objects
- //that it contains are deleted
-
- if (fParticlesEvent == 0x0) fParticlesEvent = new AliHBTEvent();
- if (fTracksEvent == 0x0) fTracksEvent = new AliHBTEvent();
-
- fParticlesEvent->Reset();
- fTracksEvent->Reset();
-
- do //do{}while; is OK even if 0 dirs specified. In that case we try to read from "./"
- {
-
- if (fRunLoader == 0x0)
- if (OpenNextSession()) continue;//directory counter is increased inside in case of error
-
- if (fCurrentEvent == fRunLoader->GetNumberOfEvents())
- {
- //read next directory
- delete fRunLoader;//close current session
- fRunLoader = 0x0;//assure pointer is null
- fCurrentDir++;//go to next dir
- continue;//directory counter is increased inside in case of error
- }
-
- Info("ReadNext","Reading Event %d",fCurrentEvent);
-
- fRunLoader->GetEvent(fCurrentEvent);
- TTree *tracktree = fTPCLoader->TreeT();//get the tree
- if (!tracktree) //check if we got the tree
- {//if not return with error
- Error("ReadNext","Can't get a tree with TPC tracks !\n");
- fCurrentEvent++;//go to next dir
- continue;
- }
-
- TBranch *trackbranch=tracktree->GetBranch("tracks");//get the branch with tracks
- if (!trackbranch) ////check if we got the branch
- {//if not return with error
- Error("ReadNext","Can't get a branch with TPC tracks !\n");
- fCurrentEvent++;//go to next dir
- continue;
- }
- Int_t ntpctracks=(Int_t)tracktree->GetEntries();//get number of TPC tracks
- Info("ReadNext","Found %d TPC tracks.",ntpctracks);
- //Copy tracks to array
-
- AliTPCtrack *iotrack=0;
-
- for (Int_t i=0; i<ntpctracks; i++) //loop over all tpc tracks
- {
- iotrack=new AliTPCtrack; //create new tracks
- trackbranch->SetAddress(&iotrack); //tell the branch ehere to put track data from tree(file)
- tracktree->GetEvent(i); //stream track i to the iotrack
- tarray->AddLast(iotrack); //put the track in the array
- }
-
- Double_t xk;
- Double_t par[5];
- Float_t phi, lam, pt;//angles and transverse momentum
- Int_t label; //label of the current track
-
- AliStack* stack = fRunLoader->Stack();
- if (stack == 0x0)
- {
- Error("ReadNext","Can not get stack for current event",fCurrentEvent);
- fCurrentEvent++;
- continue;
- }
- stack->Particles();
-
- for (Int_t i=0; i<ntpctracks; i++) //loop over all good tracks
- {
- iotrack = (AliTPCtrack*)tarray->At(i);
- label = iotrack->GetLabel();
-
- if (label < 0) continue;
-
- if (CheckTrack(iotrack)) continue;//Checks the cuts on track parameters cov. mtx etc
-
- TParticle *p = (TParticle*)stack->Particle(label);
-
- if(p == 0x0) continue; //if returned pointer is NULL
- if(p->GetPDG() == 0x0) continue; //if particle has crezy PDG code (not known to our database)
-
- if(Rejected(p->GetPdgCode())) continue; //check if we are intersted with particles of this type
- //if not take next partilce
-
- AliHBTParticle* part = new AliHBTParticle(*p,i);
- if(Rejected(part)) { delete part; continue;}//check if meets all criteria of any of our cuts
- //if it does not delete it and take next good track
-
-// iotrack->PropagateTo(3.,0.0028,65.19);
-// iotrack->PropagateToVertex(36.66,1.2e-3);//orig
- iotrack->PropagateToVertex(50.,0.0353);
-
- iotrack->GetExternalParameters(xk,par); //get properties of the track
-
- phi=TMath::ASin(par[2]) + iotrack->GetAlpha();
- if (phi<-TMath::Pi()) phi+=2*TMath::Pi();
- if (phi>=TMath::Pi()) phi-=2*TMath::Pi();
- lam=par[3];
- pt=1.0/TMath::Abs(par[4]);
-
- Double_t tpx = pt * TMath::Cos(phi); //track x coordinate of momentum
- Double_t tpy = pt * TMath::Sin(phi); //track y coordinate of momentum
- Double_t tpz = pt * lam; //track z coordinate of momentum
-
- Double_t mass = p->GetMass();
- Double_t tEtot = TMath::Sqrt( tpx*tpx + tpy*tpy + tpz*tpz + mass*mass);//total energy of the track
-
- AliHBTParticle* track = new AliHBTParticle(p->GetPdgCode(),i, tpx, tpy , tpz, tEtot, 0., 0., 0., 0.);
- if(Rejected(track))//check if meets all criteria of any of our cuts
- //if it does not delete it and take next good track
- {
- delete track;
- delete part;
- continue;
- }
-
- if (fNTrackPoints > 0)
- {
- AliHBTTrackPoints* tpts = new AliHBTTrackPoints(fNTrackPoints,iotrack,fdR);
- track->SetTrackPoints(tpts);
- }
- if ( fClusterMap )
- {
- AliHBTClusterMap* cmap = new AliHBTClusterMap(iotrack);
- track->SetClusterMap(cmap);
- }
-
- fParticlesEvent->AddParticle(part);
- fTracksEvent->AddParticle(track);
- }
-
- Info("ReadNext","Read %d tracks and %d particles from event %d (event %d in dir %d).",
- fParticlesEvent->GetNumberOfParticles(), fTracksEvent->GetNumberOfParticles(),
- fNEventsRead,fCurrentEvent,fCurrentDir);
- fTrackCounter->Fill(fTracksEvent->GetNumberOfParticles());
- fCurrentEvent++;
- fNEventsRead++;
- delete tarray;
- return 0;
- }while(fCurrentDir < GetNumberOfDirs());
-
- delete tarray;
- return 1;
- }
-/********************************************************************/
-
-Int_t AliHBTReaderTPC::OpenNextSession()
-{
-//Opens session thats from fCurrentDir
- TString filename = GetDirName(fCurrentDir);
- if (filename.IsNull())
- {
- DoOpenError("Can not get directory name");
- return 1;
- }
- filename = filename +"/"+ fFileName;
- fRunLoader = AliRunLoader::Open(filename,AliConfig::GetDefaultEventFolderName());
- if( fRunLoader == 0x0)
- {
- DoOpenError("Can not open session.");
- return 1;
- }
-
- fRunLoader->LoadHeader();
- if (fRunLoader->GetNumberOfEvents() <= 0)
- {
- DoOpenError("There is no events in this directory.");
- return 1;
- }
-
- if (fRunLoader->LoadKinematics())
- {
- DoOpenError("Error occured while loading kinematics.");
- return 1;
- }
- fTPCLoader = (AliTPCLoader*)fRunLoader->GetLoader("TPCLoader");
- if ( fTPCLoader == 0x0)
- {
- DoOpenError("Exiting due to problems with opening files.");
- return 1;
- }
- Info("OpenNextSession","________________________________________________________");
- Info("OpenNextSession","Found %d event(s) in directory %s",
- fRunLoader->GetNumberOfEvents(),GetDirName(fCurrentDir).Data());
- Float_t mf;
- if (fUseMagFFromRun)
- {
- fRunLoader->LoadgAlice();
- mf = fRunLoader->GetAliRun()->Field()->SolenoidField();
- Info("OpenNextSession","Setting Magnetic Field from run: B=%fT",mf/10.);
- fRunLoader->UnloadgAlice();
- }
- else
- {
- Info("OpenNextSession","Setting Own Magnetic Field: B=%fT",fMagneticField);
- if (fMagneticField == 0x0)
- {
- Fatal("OpenNextSession","Magnetic field can not be 0.");
- return 1;//pro forma
- }
- mf = fMagneticField*10.;
- }
- AliKalmanTrack::SetConvConst(1000/0.299792458/mf);
-
-
- if (fTPCLoader->LoadTracks())
- {
- DoOpenError("Error occured while loading TPC tracks.");
- return 1;
- }
-
- fCurrentEvent = 0;
- return 0;
-}
-/********************************************************************/
-
-void AliHBTReaderTPC::DoOpenError( const char *va_(fmt), ...)
-{
- // Does error display and clean-up in case error caught on Open Next Session
-
- va_list ap;
- va_start(ap,va_(fmt));
- Error("OpenNextSession", va_(fmt), ap);
- va_end(ap);
-
- delete fRunLoader;
- fRunLoader = 0x0;
- fTPCLoader = 0x0;
- fCurrentDir++;
-}
-
-/********************************************************************/
-Bool_t AliHBTReaderTPC::CheckTrack(AliTPCtrack* t) const
-{
- //Performs check of the track
- if ( (t->GetNumberOfClusters() > fNClustMax) || (t->GetNumberOfClusters() < fNClustMin) ) return kTRUE;
-
- Float_t chisqpercl = t->GetChi2()/((Double_t)t->GetNumberOfClusters());
- if ( (chisqpercl < fNChi2PerClustMin) || (chisqpercl > fNChi2PerClustMax) ) return kTRUE;
-
- Double_t cc[15];
- t->GetExternalCovariance(cc);
-
- if ( (cc[0] < fC00Min) || (cc[0] > fC00Max) ) return kTRUE;
- if ( (cc[2] < fC11Min) || (cc[2] > fC11Max) ) return kTRUE;
- if ( (cc[5] < fC22Min) || (cc[5] > fC22Max) ) return kTRUE;
- if ( (cc[9] < fC33Min) || (cc[9] > fC33Max) ) return kTRUE;
- if ( (cc[14] < fC44Min) || (cc[14] > fC44Max) ) return kTRUE;
-
-
- return kFALSE;
-
-}
-/********************************************************************/
-
-void AliHBTReaderTPC::SetNClustersRange(Int_t min,Int_t max)
-{
- //sets range of Number Of Clusters that tracks have to have
- fNClustMin = min;
- fNClustMax = max;
-}
-/********************************************************************/
-
-void AliHBTReaderTPC::SetChi2PerCluserRange(Float_t min, Float_t max)
-{
- //sets range of Chi2 per Cluster
- fNChi2PerClustMin = min;
- fNChi2PerClustMax = max;
-}
-/********************************************************************/
-
-void AliHBTReaderTPC::SetC00Range(Float_t min, Float_t max)
-{
- //Sets range of C00 parameter of covariance matrix of the track
- //it defines uncertainty of the momentum
- fC00Min = min;
- fC00Max = max;
-}
-/********************************************************************/
-
-void AliHBTReaderTPC::SetC11Range(Float_t min, Float_t max)
-{
- //Sets range of C11 parameter of covariance matrix of the track
- //it defines uncertainty of the momentum
- fC11Min = min;
- fC11Max = max;
-}
-/********************************************************************/
-
-void AliHBTReaderTPC::SetC22Range(Float_t min, Float_t max)
-{
- //Sets range of C22 parameter of covariance matrix of the track
- //it defines uncertainty of the momentum
- fC22Min = min;
- fC22Max = max;
-}
-/********************************************************************/
-
-void AliHBTReaderTPC::SetC33Range(Float_t min, Float_t max)
-{
- //Sets range of C33 parameter of covariance matrix of the track
- //it defines uncertainty of the momentum
- fC33Min = min;
- fC33Max = max;
-}
-/********************************************************************/
-
-void AliHBTReaderTPC::SetC44Range(Float_t min, Float_t max)
-{
- //Sets range of C44 parameter of covariance matrix of the track
- //it defines uncertainty of the momentum
- fC44Min = min;
- fC44Max = max;
-}
-
-/********************************************************************/
-/********************************************************************/
-/********************************************************************/
-
-/*
-void (AliTPCtrack* iotrack, Double_t curv)
-{
- Double_t x[5];
- iotrack->GetExternalPara
- //xk is a
- Double_t fP4=iotrack->GetC();
- Double_t fP2=iotrack->GetEta();
-
- Double_t x1=fX, x2=x1+(xk-x1), dx=x2-x1, y1=fP0, z1=fP1;
- Double_t c1=fP4*x1 - fP2, r1=sqrt(1.- c1*c1);
- Double_t c2=fP4*x2 - fP2, r2=sqrt(1.- c2*c2);
-
- fP0 += dx*(c1+c2)/(r1+r2);
- fP1 += dx*(c1+c2)/(c1*r2 + c2*r1)*fP3;
-
-}
-*/
-
+++ /dev/null
-#ifndef AliHBTReaderTPC_H
-#define AliHBTReaderTPC_H
-//______________________________________________
-//
-// class AliHBTReaderTPC
-//
-// reader for TPC tracks
-// needs galice.root
-//
-// more info: http://aliweb.cern.ch/people/skowron/analyzer/index.html
-// Piotr.Skowronski@cern.ch
-//
-///////////////////////////////////////////////////////////////////////////
-
-#include "AliHBTReader.h"
-#include <TString.h>
-
-class TFile;
-class TArrayF;
-class AliRunLoader;
-class AliTPCLoader;
-class AliTPCtrack;
-
-class AliHBTReaderTPC: public AliHBTReader
-{
- public:
- AliHBTReaderTPC();
- AliHBTReaderTPC(const Char_t* galicefilename);
- AliHBTReaderTPC(TObjArray* dirs, const Char_t* galicefilename = "galice.root");
- AliHBTReaderTPC(const AliHBTReaderTPC& in);
-
- virtual ~AliHBTReaderTPC();
-
- AliHBTReaderTPC& operator=(const AliHBTReaderTPC& in);
-
- void Rewind();
-
- Bool_t ReadsTracks() const {return kTRUE;}
- Bool_t ReadsParticles() const {return kTRUE;}
-
- void SetMagneticField(Float_t mf){fMagneticField=mf;}
- void UseMagneticFieldFromRun(Bool_t flag = kTRUE){fUseMagFFromRun=flag;}
-
- void SetNClustersRange(Int_t min,Int_t max);
- void SetChi2PerCluserRange(Float_t min, Float_t max);
- void SetC00Range(Float_t min, Float_t max);
- void SetC11Range(Float_t min, Float_t max);
- void SetC22Range(Float_t min, Float_t max);
- void SetC33Range(Float_t min, Float_t max);
- void SetC44Range(Float_t min, Float_t max);
- void SetNumberOfTrackPoints(Int_t n = 5,Float_t dr = 30.0) {fNTrackPoints = n; fdR = dr;}
- Int_t GetNumberOfTrackPoints() const {return fNTrackPoints;}
- void SetClusterMap(Bool_t flag = kTRUE){fClusterMap = flag;}
-
- protected:
- Int_t ReadNext();
- Int_t OpenNextSession();
- void DoOpenError(const char* msgfmt, ...);
-
- TString fFileName;//name of the file with galice.root
- AliRunLoader* fRunLoader;//!RL
- AliTPCLoader* fTPCLoader;//!TPCLoader
- Float_t fMagneticField;//magnetic field value that was enforced while reading
- Bool_t fUseMagFFromRun;//flag indicating if using field specified in gAlice (kTRUE)
- // or enforece other defined by fMagneticField
-
- Int_t fNTrackPoints;//number of track points; if==0 track points are not created
- Float_t fdR;//spacing between points (along radius) in cm
- //Track Points are needed for Anti-Merging Cut
-
- Bool_t fClusterMap;//Flag indicating if Claster Map should be created for each track
- //Claster map is needed for Anti-Splitting Cut
-
- //Cut Parameters specific to TPC tracks
-
- Int_t fNClustMin;//Number of clusters min value
- Int_t fNClustMax;//Number of clusters max value
-
- Float_t fNChi2PerClustMin;//Chi^2 per number of clusters min value
- Float_t fNChi2PerClustMax;//Chi^2 per number of clusters max value
-
- Float_t fC00Min;//C00 (0th diagonal element of covariance matrix) min value
- Float_t fC00Max;//C00 (0th diagonal element of covariance matrix) max value
-
- Float_t fC11Min;//C11 (1th diagonal element of covariance matrix) min value
- Float_t fC11Max;//C11 (1th diagonal element of covariance matrix) max value
-
- Float_t fC22Min;//C22 (2th diagonal element of covariance matrix) min value
- Float_t fC22Max;//C22 (2th diagonal element of covariance matrix) max value
-
- Float_t fC33Min;//C33 (3th diagonal element of covariance matrix) min value
- Float_t fC33Max;//C33 (3th diagonal element of covariance matrix) max value
-
- Float_t fC44Min;//C44 (4th diagonal element of covariance matrix) min value
- Float_t fC44Max;//C44 (4th diagonal element of covariance matrix) max value
-
- private:
-
- Bool_t CheckTrack(AliTPCtrack* t) const;
-
- ClassDef(AliHBTReaderTPC,3)
-};
-
-
-#endif
+++ /dev/null
-#if 0
- #include "/afs/cern.ch/user/s/skowron/aliroot/my/TPC/alles.h"
- #include "AliHBTReader.h"
- #include "AliHBTReaderKineTree.h"
- #include "AliHBTReaderITSv2.h"
- #include "AliHBTReaderITSv1.h"
- #include "AliHBTReaderTPC.h"
- #include "AliHBTParticleCut.h"
- #include "AliHBTEvent.h"
- #include "AliAODPairCut.h"
- #include "AliHBTQResolutionFctns.h"
- #include "AliHBTTwoTrackEffFctn.h"
- #include "AliHBTCorrelFctn.h"
- #include "TSystem.h"
- #include "TObjString.h"
- #include "TString.h"
- #include "AliPDG.h"
-#endif
-
-
-void AliHBTWriteInternFormat(Option_t* datatype, Option_t* processopt="TracksAndParticles",
- Int_t first = -1,Int_t last = -1,
- char *outfile = "data.root")
- {
-//HBT Anlysis Macro
-//Anlyzes TPC recontructed tracks and simulated particles that corresponds to them
-
-//datatype defines type of data to be read
-// Kine - analyzes Kine Tree: simulated particles
-// TPC - analyzes TPC tracking + particles corresponding to these tracks
-// ITSv1 - analyzes ITSv1 ----------------------//--------------------------
-// ITSv2 - analyzes ITSv2 ----------------------//--------------------------
-
-//processopt defines option passed to AliHBTAnlysis::Process method
-// default: TracksAndParticles - process both recontructed tracks and sim. particles corresponding to them
-// Tracks - process only recontructed tracks
-// Particles - process only simulated particles
-
-//Reads data from diroctories from first to last(including)
-// For examples if first=3 and last=5 it reads from
-// ./3/
-// ./4/
-// ./5/
-//if first or last is negative (or both), it reads from current directory
-//
-//these names I use when analysis is done directly from CASTOR files via RFIO
-
- const char* basedir=".";
- const char* serie="";
- const char* field = "";
- cout<<"AliHBTWriteInternFormat.C: datatype is "<<datatype<<" dir is basedir"<<endl;
- // Dynamically link some shared libs
-
- cout<<"AliHBTWriteInternFormat.C: Loading HBTAN .....\n";
- gSystem->Load("$(ALICE_ROOT)/lib/tgt_$(ALICE_TARGET)/libHBTAN");
- cout<<"AliHBTWriteInternFormat.C: ..... Loaded\n";
-
- Bool_t multcheck = kTRUE;
- /***********************************************************/
-
- AliHBTReader* reader;
- Int_t kine = strcmp(datatype,"Kine");
- Int_t ESD = strcmp(datatype,"ESD");
- Int_t TPC = strcmp(datatype,"TPC");
- Int_t ITSv1 = strcmp(datatype,"ITSv1");
- Int_t ITSv2 = strcmp(datatype,"ITSv2");
- Int_t intern = strcmp(datatype,"Intern");
-
- if(!kine)
- {
- reader = new AliHBTReaderKineTree();
- processopt="Particles"; //this reader by definition reads only simulated particles
- multcheck = kFALSE;
- }
- else if(!ESD)
- {
- AliHBTReaderESD* esdreader = new AliHBTReaderESD();
- esdreader->ReadParticles(kTRUE);
- reader = esdreader;
- multcheck = kTRUE;
- }
- else if(!TPC)
- {
- cout<<"AliHBTWriteInternFormat.C: Creating Reader TPC .....\n";
- reader = new AliHBTReaderTPC();
- multcheck = kFALSE;
- cout<<"AliHBTWriteInternFormat.C: ..... Created\n";
- }
- else if(!ITSv1)
- {
- reader = new AliHBTReaderITSv1();
- multcheck = kFALSE;
- }
- else if(!ITSv2)
- {
- cout<<"AliHBTWriteInternFormat.C: Creating Reader ITSv2 .....\n";
- reader = new AliHBTReaderITSv2();
- multcheck = kFALSE;
- cout<<"AliHBTWriteInternFormat.C: ..... Created\n";
- }
- else if(!intern)
- {
- reader = new AliHBTReaderInternal("data.root");
- multcheck = kTRUE;
- }
- else
- {
- cerr<<"Option "<<datatype<<" not recognized. Exiting"<<endl;
- return;
- }
-
- TObjArray* dirs=0;
- if ( (first >= 0) && (last>=0) && ( (last-first)>=0 ) )
- {//read from many dirs dirs
- cout<<"AliHBTWriteInternFormat.C: ..... Setting dirs first="<<first<<" last="<<last<<"\n";
- char buff[50];
- dirs = new TObjArray(last-first+1);
- dirs->SetOwner();
- for (Int_t i = first; i<=last; i++)
- {
- sprintf(buff,"%s/%s/%s/%d",basedir,field,serie,i);
- TObjString *odir= new TObjString(buff);
- dirs->Add(odir);
- }
- }
-
- reader->SetDirs(dirs);
-
- cout<<"AliHBTWriteInternFormat.C: P R O C S E S S I N G .....\n\n";
- AliHBTReaderInternal::Write(reader,outfile,multcheck);
- cout<<"\n\nAliHBTWriteInternFormat.C: F I N I S H E D\n";
-
- if (dirs) delete dirs;
- delete reader;
- }
-
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff