-
#include "AliHBTReaderITSv2.h"
-#include <iostream.h>
-#include <fstream.h>
+#include <Varargs.h>
+
#include <TString.h>
#include <TObjString.h>
#include <TTree.h>
-#include <TFile.h>
#include <TParticle.h>
#include <AliRun.h>
+#include <AliRunLoader.h>
+#include <AliLoader.h>
+#include <AliStack.h>
#include <AliMagF.h>
#include <AliITStrackV2.h>
-//#include <AliITSParam.h>
-#include <AliITStrackerV2.h>
-#include <AliITSgeom.h>
#include "AliHBTRun.h"
#include "AliHBTEvent.h"
ClassImp(AliHBTReaderITSv2)
-AliHBTReaderITSv2::
- AliHBTReaderITSv2(const Char_t* trackfilename, const Char_t* clusterfilename,
- const Char_t* galicefilename)
- :fTrackFileName(trackfilename),fClusterFileName(clusterfilename),
- fGAliceFileName(galicefilename)
+AliHBTReaderITSv2::AliHBTReaderITSv2():
+ fFileName("galice.root"),
+ fRunLoader(0x0),
+ fITSLoader(0x0),
+ fMagneticField(0.0),
+ fUseMagFFromRun(kTRUE)
{
//constructor,
//Defaults:
- // trackfilename = "AliITStracksV2.root"
- // clusterfilename = "AliITSclustersV2.root"
// galicefilename = "galice.root"
- fParticles = new AliHBTRun();
- fTracks = new AliHBTRun();
- fIsRead = kFALSE;
}
/********************************************************************/
-AliHBTReaderITSv2::
- AliHBTReaderITSv2(TObjArray* dirs, const Char_t* trackfilename,
- const Char_t* clusterfilename, const Char_t* galicefilename)
- : AliHBTReader(dirs),
- fTrackFileName(trackfilename),fClusterFileName(clusterfilename),
- fGAliceFileName(galicefilename)
+AliHBTReaderITSv2::AliHBTReaderITSv2(const Char_t* galicefilename):
+ fFileName(galicefilename),
+ fRunLoader(0x0),
+ fITSLoader(0x0),
+ fMagneticField(0.0),
+ fUseMagFFromRun(kTRUE)
{
//constructor,
//Defaults:
- // trackfilename = "AliITStracksV2.root"
- // clusterfilename = "AliITSclustersV2.root"
// galicefilename = "galice.root"
-
- fParticles = new AliHBTRun();
- fTracks = new AliHBTRun();
- fIsRead = kFALSE;
}
/********************************************************************/
-AliHBTReaderITSv2::~AliHBTReaderITSv2()
- {
- if (fParticles) delete fParticles;
- if (fTracks) delete fTracks;
- }
-/********************************************************************/
-/********************************************************************/
-
-AliHBTEvent* AliHBTReaderITSv2::GetParticleEvent(Int_t n)
- {
- //returns Nth event with simulated particles
- if (!fIsRead)
- if(Read(fParticles,fTracks))
- {
- Error("GetParticleEvent","Error in reading");
- return 0x0;
- }
-
- return fParticles->GetEvent(n);
- }
-/********************************************************************/
-
-AliHBTEvent* AliHBTReaderITSv2::GetTrackEvent(Int_t n)
- {
- //returns Nth event with reconstructed tracks
- if (!fIsRead)
- if(Read(fParticles,fTracks))
- {
- Error("GetTrackEvent","Error in reading");
- return 0x0;
- }
- return fTracks->GetEvent(n);
- }
+AliHBTReaderITSv2::AliHBTReaderITSv2(TObjArray* dirs, const Char_t* galicefilename):
+ AliHBTReader(dirs),
+ fFileName(galicefilename),
+ fRunLoader(0x0),
+ fITSLoader(0x0),
+ fMagneticField(0.0),
+ fUseMagFFromRun(kTRUE)
+{
+ //constructor,
+ //Defaults:
+ // galicefilename = "galice.root"
+
+}
/********************************************************************/
-
-Int_t AliHBTReaderITSv2::GetNumberOfPartEvents()
- {
- //returns number of events of particles
- if (!fIsRead)
- if(Read(fParticles,fTracks))
- {
- Error("GetNumberOfPartEvents","Error in reading");
- return 0;
- }
- return fParticles->GetNumberOfEvents();
- }
-
+
+void AliHBTReaderITSv2::Rewind()
+{
+ //rewinds reading
+ delete fRunLoader;
+ fRunLoader = 0x0;
+ fCurrentDir = 0;
+ fNEventsRead= 0;
+}
/********************************************************************/
-Int_t AliHBTReaderITSv2::GetNumberOfTrackEvents()
- {
- //returns number of events of tracks
- if (!fIsRead)
- if(Read(fParticles,fTracks))
- {
- Error("GetNumberOfTrackEvents","Error in reading");
- return 0;
- }
- return fTracks->GetNumberOfEvents();
- }
-
-/********************************************************************/
+AliHBTReaderITSv2::~AliHBTReaderITSv2()
+{
+ //dtor
+ delete fRunLoader;
+}
/********************************************************************/
-Int_t AliHBTReaderITSv2::Read(AliHBTRun* particles, AliHBTRun *tracks)
+
+Int_t AliHBTReaderITSv2::ReadNext()
{
- Int_t Nevents = 0; //number of events found in given directory
- Int_t Ndirs; //number of the directories to be read
- Int_t Ntracks; //number of tracks in current event
- Int_t currentdir = 0; //number of events in the current directory
- Int_t totalNevents = 0; //total number of events read from all directories up to now
- register Int_t i; //iterator
+//reads data from next event
+
+ register Int_t i = 0; //iterator
- TFile *aTracksFile;//file with tracks
- TFile *aClustersFile;//file with clusters
- TFile *aGAliceFile;//file name with galice
-
// AliITStrackerV2 *tracker; // ITS tracker - used for cooking labels
- TTree *tracktree; // tree for tracks
+ TTree *tracktree = 0x0; // tree for tracks
+ AliITStrackV2 *iotrack = 0x0;
Double_t xk;
Double_t par[5]; //Kalman track parameters
Float_t phi, lam, pt;//angles and transverse momentum
Int_t label; //label of the current track
-
- char tname[100]; //buffer for tree name
- AliITStrackV2 *iotrack= new AliITStrackV2(); //buffer track for reading data from tree
-
- if (!particles) //check if an object is instatiated
- {
- Error("Read"," particles object must instatiated before passing it to the reader");
- }
- if (!tracks) //check if an object is instatiated
- {
- Error("Read"," tracks object must instatiated before passing it to the reader");
- }
- particles->Reset();//clear runs == delete all old events
- tracks->Reset();
-
- if (fDirs) //if array with directories is supplied by user
- {
- Ndirs = fDirs->GetEntries(); //get the number if directories
- }
- else
- {
- Ndirs = 0; //if the array is not supplied read only from current directory
- }
-// cout<<"Found "<<Ndirs<<" directory entries"<<endl;
+ if (fParticlesEvent == 0x0) fParticlesEvent = new AliHBTEvent();
+ if (fTracksEvent == 0x0) fTracksEvent = new AliHBTEvent();
+
+ fParticlesEvent->Reset();
+ fTracksEvent->Reset();
do //do while is good even if Ndirs==0 (than read from current directory)
{
- if( (i=OpenFiles(aTracksFile,aClustersFile,aGAliceFile,currentdir)) )
+ if (fRunLoader == 0x0)
+ if (OpenNextFile()) continue;//directory counter is increased inside in case of error
+
+ if (fCurrentEvent == fRunLoader->GetNumberOfEvents())
{
- Error("Read","Exiting due to problems with opening files. Errorcode %d",i);
- delete iotrack;
- return i;
+ //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);
- if (gAlice->TreeE())//check if tree E exists
+ tracktree=fITSLoader->TreeT();
+ if (!tracktree)
{
- Nevents = (Int_t)gAlice->TreeE()->GetEntries();//if yes get number of events in gAlice
- cout<<"________________________________________________________\n";
- cout<<"Found "<<Nevents<<" event(s) in directory "<<GetDirName(currentdir)<<endl;
- cout<<"Setting Magnetic Field. Factor is "<<gAlice->Field()->Factor()<<endl;
- AliKalmanTrack::SetConvConst(100/0.299792458/0.2/gAlice->Field()->Factor());
+ Error("ReadNext","Can't get a tree with ITS tracks");
+ fCurrentEvent++;
+ continue;
}
- else
- {//if not return an error
- Error("Read","Can not find Header tree (TreeE) in gAlice");
- delete iotrack;
- return 1;
- }
-
- AliITSgeom *geom=(AliITSgeom*)aClustersFile->Get("AliITSgeom");
- if (!geom)
- {
- Error("Read","Can't get the ITS geometry!");
- delete iotrack;
- return 2;
+
+ TBranch *tbranch=tracktree->GetBranch("tracks");
+ if (!tbranch)
+ {
+ Error("ReadNext","Can't get a branch with ITS tracks");
+ fCurrentEvent++;
+ continue;
}
- for(Int_t currentEvent =0; currentEvent<Nevents;currentEvent++)//loop over all events
+ AliStack* stack = fRunLoader->Stack();
+ if (stack == 0x0)
{
- cout<<"Reading Event "<<currentEvent<<endl;
-
- aGAliceFile->cd();
- gAlice->GetEvent(currentEvent);
-
-// TParticle * part = gAlice->Particle(0);
-
- aClustersFile->cd();
-// Double_t orz=part->Vz();
-// tracker = new AliITStrackerV2(geom,currentEvent,orz); //<---- this is for Massimo version
-// tracker = new AliITStrackerV2(geom,currentEvent);
- sprintf(tname,"TreeT_ITS_%d",currentEvent);
-
- tracktree=(TTree*)aTracksFile->Get(tname);
- if (!tracktree)
- {
- Error("Read","Can't get a tree with ITS tracks");
- delete iotrack;
- // delete tracker;
- return 4;
- }
- TBranch *tbranch=tracktree->GetBranch("tracks");
-
- Ntracks=(Int_t)tracktree->GetEntries();
-
- Int_t accepted = 0;
- Int_t tpcfault = 0;
- Int_t itsfault = 0;
- for (i=0; i<Ntracks; i++) //loop over all tpc tracks
+ Error("ReadNext","Can not get stack for current event",fCurrentEvent);
+ fCurrentEvent++;
+ continue;
+ }
+
+ //must be here because on the beginning conv. const. is not set yet
+ if (iotrack == 0x0) iotrack = new AliITStrackV2(); //buffer track for reading data from tree
+
+ Int_t ntr = (Int_t)tracktree->GetEntries();
+
+ for (i=0; i < ntr; i++) //loop over all tpc tracks
+ {
+ tbranch->SetAddress(&iotrack);
+ tracktree->GetEvent(i);
+
+ label=iotrack->GetLabel();
+ if (label < 0)
+ {
+ continue;
+ }
+
+ TParticle *p = 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(Pass(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(Pass(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();
+
+ 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(Pass(track))//check if meets all criteria of any of our cuts
+ //if it does not delete it and take next good track
{
- if(i%100 == 0)cout<<"all: "<<i<<" accepted: "<<accepted<<" tpc faults: "<<tpcfault<<"\r";
-
- tbranch->SetAddress(&iotrack);
- tracktree->GetEvent(i);
-
- label=iotrack->GetLabel();
- if (label < 0)
- {
- tpcfault++;
- continue;
- }
-
- TParticle *p = (TParticle*)gAlice->Particle(label);
- if(Pass(p->GetPdgCode())) continue; //check if we are intersted with particles of this type
- //if not take next partilce
-
- AliHBTParticle* part = new AliHBTParticle(*p);
- if(Pass(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->Propagate(iotrack->GetAlpha(),3.,0.1/65.19*1.848,0.1*1.848);
- iotrack->PropagateToVertex();
-
- 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(), tpx, tpy , tpz, tEtot, 0., 0., 0., 0.);
- if(Pass(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;
- }
- particles->AddParticle(totalNevents,part);//put track and particle on the run
- tracks->AddParticle(totalNevents,track);
- accepted++;
- }//end of loop over tracks in the event
+ delete track;
+ delete part;
+ continue;
+ }
- aTracksFile->Delete(tname);
- aTracksFile->Delete("tracks");
-// delete tracker;
+ fParticlesEvent->AddParticle(part);
+ fTracksEvent->AddParticle(track);
+ }//end of loop over tracks in the event
- totalNevents++;
- cout<<"all: "<<i<<" accepted: "<<accepted<<" tpc faults: "<<tpcfault<<" its faults: "<<itsfault<<endl;
+ Info("ReadNext","Read %d tracks and %d particles from event %d (event %d in dir %d).",
+ fParticlesEvent->GetNumberOfParticles(), fTracksEvent->GetNumberOfParticles(),
+ fNEventsRead,fCurrentEvent,fCurrentDir);
- }//end of loop over events in current directory
- CloseFiles(aTracksFile,aClustersFile,aGAliceFile);
- currentdir++;
- }while(currentdir < Ndirs);//end of loop over directories specified in fDirs Obj Array
+ fCurrentEvent++;
+ fNEventsRead++;
+ delete iotrack;
+ return 0;
+ }while(fCurrentDir < GetNumberOfDirs());//end of loop over directories specified in fDirs Obj Array
delete iotrack;
- fIsRead = kTRUE;
- return 0;
+ return 1;
}
/********************************************************************/
-Int_t AliHBTReaderITSv2::OpenFiles
-(TFile*& aTracksFile, TFile*& aClustersFile, TFile*& agAliceFile,Int_t event)
+Int_t AliHBTReaderITSv2::OpenNextFile()
{
- //opens all the files
-
-
- const TString& dirname = GetDirName(event);
- if (dirname == "")
- {
- Error("OpenFiles","Can not get directory name");
- return 4;
- }
-
- TString filename = dirname +"/"+ fTrackFileName;
- aTracksFile = TFile::Open(filename.Data());
- if ( aTracksFile == 0x0 )
- {
- Error("OpenFiles","Can't open file with tacks named %s",filename.Data());
- return 1;
- }
- if (!aTracksFile->IsOpen())
- {
- Error("OpenFiles","Can't open file with tacks named %s",filename.Data());
- return 1;
- }
-
- filename = dirname +"/"+ fClusterFileName;
- aClustersFile = TFile::Open(filename.Data());
- if ( aClustersFile == 0x0 )
- {
- Error("OpenFiles","Can't open file with TPC clusters named %s",filename.Data());
- return 2;
- }
- if (!aClustersFile->IsOpen())
- {
- Error("OpenFiles","Can't open file with TPC clusters named %s",filename.Data());
- return 2;
- }
+ //opens next file
+ TString filename = GetDirName(fCurrentDir);
+ if (filename.IsNull())
+ {
+ DoOpenError("Can not get directory name");
+ return 1;
+ }
+ filename = filename +"/"+ fFileName;
+ fRunLoader = AliRunLoader::Open(filename,AliConfig::fgkDefaultEventFolderName);
+ if( fRunLoader == 0x0)
+ {
+ DoOpenError("Can not open session.");
+ return 1;
+ }
+
+ if (fRunLoader->GetNumberOfEvents() <= 0)
+ {
+ DoOpenError("There is no events in this directory.");
+ return 1;
+ }
- filename = dirname +"/"+ fGAliceFileName;
- agAliceFile = TFile::Open(filename.Data());
- if ( agAliceFile== 0x0)
- {
- Error("OpenFiles","Can't open file with TPC clusters named %s",filename.Data());
- return 3;
- }
- if (!agAliceFile->IsOpen())
- {
- Error("OpenFiles","Can't open file with TPC clusters named %s",filename.Data());
- return 3;
- }
+ if (fRunLoader->LoadKinematics())
+ {
+ DoOpenError("Error occured while loading kinematics.");
+ return 1;
+ }
+ fITSLoader = fRunLoader->GetLoader("ITSLoader");
+ if ( fITSLoader == 0x0)
+ {
+ DoOpenError("Exiting due to problems with opening files.");
+ return 1;
+ }
- if (!(gAlice=(AliRun*)agAliceFile->Get("gAlice")))
- {
- Error("OpenFiles","gAlice have not been found on %s !\n",filename.Data());
- return 5;
- }
+ Info("OpenNextSession","________________________________________________________");
+ Info("OpenNextSession","Found %d event(s) in directory %s",
+ fRunLoader->GetNumberOfEvents(),GetDirName(fCurrentDir).Data());
+ Float_t mf;
+ if (fUseMagFFromRun)
+ {
+ if (fRunLoader->LoadgAlice())
+ {
+ DoOpenError("Error occured while loading AliRun.");
+ return 1;
+ }
+ 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);
- return 0;
+ if (fITSLoader->LoadTracks())
+ {
+ DoOpenError("Error occured while loading TPC tracks.");
+ return 1;
+ }
+
+ fCurrentEvent = 0;
+ return 0;
}
/********************************************************************/
-/********************************************************************/
-
-void AliHBTReaderITSv2::CloseFiles(TFile*& tracksFile, TFile*& clustersFile, TFile*& gAliceFile)
+void AliHBTReaderITSv2::DoOpenError( const char *va_(fmt), ...)
{
- //closes the files
- tracksFile->Close();
- delete tracksFile;
- tracksFile = 0x0;
-
- clustersFile->Close();
- delete clustersFile;
- clustersFile = 0x0;
-
- delete gAlice;
- gAlice = 0;
+ // Does error display and clean-up in case error caught on Open Next Session
- if (gAliceFile)
- {
- gAliceFile->Close();
- delete gAliceFile;
- gAliceFile = 0x0;
- }
+ va_list ap;
+ va_start(ap,va_(fmt));
+ Error("OpenNextFile", va_(fmt), ap);
+ va_end(ap);
+
+ delete fRunLoader;
+ fRunLoader = 0x0;
+ fITSLoader = 0x0;
+ fCurrentDir++;
}
/********************************************************************/
-
-
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff