]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - HBTAN/AliHBTReaderTPC.cxx
Comilation Warning removal
[u/mrichter/AliRoot.git] / HBTAN / AliHBTReaderTPC.cxx
index e0094115432f350e9e749479b17ccc380d06ab73..1cab431c8481ac6ad97d9551fcc79973e7fe73e8 100644 (file)
@@ -1,13 +1,30 @@
 #include "AliHBTReaderTPC.h"
+//______________________________________________
+//
+// class AliHBTReaderTPC
+//
+// reader for TPC tracking
+// needs galice.root, AliTPCtracks.root, AliTPCclusters.root
+//
+// more info: http://aliweb.cern.ch/people/skowron/analyzer/index.html
+// Piotr.Skowronski@cern.ch
+//
+///////////////////////////////////////////////////////////////////////////
 #include <TTree.h>
 #include <TFile.h>
 #include <TParticle.h>
+#include <TH1.h>
+
+#include <Varargs.h>
 
 #include <AliRun.h>
+#include <AliLoader.h>
+#include <AliStack.h>
 #include <AliMagF.h>
 #include <AliTPCtrack.h>
 #include <AliTPCParam.h>
+#include <AliTPCtracker.h>
+#include <AliTPCLoader.h>
 
 #include "AliHBTRun.h"
 #include "AliHBTEvent.h"
 #include "AliHBTParticleCut.h"
 
 
+
 ClassImp(AliHBTReaderTPC)
-//______________________________________________
-//
-// class AliHBTReaderTPC
-//
-//reader for TPC tracking
-//needs galice.root, AliTPCtracks.root
-//
-//more info: http://alisoft.cern.ch/people/skowron/analyzer/index.html
-//Piotr.Skowronski@cern.ch
-//
 
+AliHBTReaderTPC::AliHBTReaderTPC():
+ fFileName("galice.root"),
+ fRunLoader(0x0),
+ fTPCLoader(0x0),
+ fMagneticField(0.0),
+ fUseMagFFromRun(kTRUE),
+ fNClustMin(0),
+ fNClustMax(190),
+ fNChi2PerClustMin(0.0),
+ fNChi2PerClustMax(10e5),
+ fC44Min(0.0),
+ fC44Max(10e5)
+{
+  //constructor, 
+  //Defaults:
+  //  galicefilename = ""  - this means: Do not open gAlice file - 
+  //                         just leave the global pointer untouched
+  
+}
 
-AliHBTReaderTPC:: AliHBTReaderTPC(const Char_t* trackfilename,
-                                  const Char_t* clusterfilename,
-                                  const Char_t* galicefilename):
- fParticles(new AliHBTRun()),
- fTracks(new AliHBTRun()),
- fTrackFileName(trackfilename),
- fClusterFileName(clusterfilename),
- fGAliceFileName(galicefilename),
- fIsRead(kFALSE),
- fMagneticField(0.4),
- fUseMagFFromRun(kTRUE)
+AliHBTReaderTPC::AliHBTReaderTPC(const Char_t* galicefilename):
+ fFileName(galicefilename),
+ fRunLoader(0x0),
+ fTPCLoader(0x0),
+ fMagneticField(0.0),
+ fUseMagFFromRun(kTRUE),
+ fNClustMin(0),
+ fNClustMax(190),
+ fNChi2PerClustMin(0.0),
+ fNChi2PerClustMax(10e5),
+ fC44Min(0.0),
+ fC44Max(10e5)
 {
   //constructor, 
   //Defaults:
-  //  trackfilename = "AliTPCtracks.root"
-  //  clusterfilename = "AliTPCclusters.root"
   //  galicefilename = ""  - this means: Do not open gAlice file - 
-  //                         just leave the global pointer untached
-  //this class is not supposed to be written, that is why I let myself and user
-  //for comfort of having default constructor which allocates dynamic memmory
-  //in case it is going to be written, it should be changed
+  //                         just leave the global pointer untouched
   
 }
 /********************************************************************/
-AliHBTReaderTPC::AliHBTReaderTPC(TObjArray* dirs,
-                  const Char_t* trackfilename, const Char_t* clusterfilename,
-                  const Char_t* galicefilename):
+AliHBTReaderTPC::AliHBTReaderTPC(TObjArray* dirs, const Char_t* galicefilename):
  AliHBTReader(dirs), 
- fParticles(new AliHBTRun()),
- fTracks(new AliHBTRun()),
- fTrackFileName(trackfilename),
- fClusterFileName(clusterfilename),
- fGAliceFileName(galicefilename),
- fIsRead(kFALSE),
- fMagneticField(0.4),
- fUseMagFFromRun(kTRUE)
+ fFileName(galicefilename),
+ fRunLoader(0x0),
+ fTPCLoader(0x0),
+ fMagneticField(0.0),
+ fUseMagFFromRun(kTRUE),
+ fNClustMin(0),
+ fNClustMax(190),
+ fNChi2PerClustMin(0.0),
+ fNChi2PerClustMax(10e5),
+ fC44Min(0.0),
+ fC44Max(10e5)
 {
   //constructor, 
   //Defaults:
-  //  trackfilename = "AliTPCtracks.root"
-  //  clusterfilename = "AliTPCclusters.root"
   //  galicefilename = ""  - this means: Do not open gAlice file - 
   //                         just leave the global pointer untached
+  
 }
 /********************************************************************/
 
 AliHBTReaderTPC::~AliHBTReaderTPC()
- {
+{
  //desctructor
-   delete fParticles;
-   delete fTracks;
- }
-/********************************************************************/
-
-AliHBTEvent* AliHBTReaderTPC::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* AliHBTReaderTPC::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);
- }
-/********************************************************************/
-
-Int_t AliHBTReaderTPC::GetNumberOfPartEvents()
- {
- //returns number of events of particles
-   if (!fIsRead) 
-    if ( Read(fParticles,fTracks))
-     {
-       Error("GetNumberOfPartEvents","Error in reading");
-       return 0;
-     }
-   return fParticles->GetNumberOfEvents();
- }
-
+   if (AliHBTParticle::GetDebug()) 
+    {
+      Info("~AliHBTReaderTPC","deleting Run Loader");
+      AliLoader::SetDebug(AliHBTParticle::GetDebug());
+    }
+   
+   delete fRunLoader;
+   
+   if (AliHBTParticle::GetDebug()) 
+    {
+      Info("~AliHBTReaderTPC","deleting Run Loader Done");
+    }
+}
 /********************************************************************/
-Int_t AliHBTReaderTPC::GetNumberOfTrackEvents()
- {
- //returns number of events of tracks
-  if (!fIsRead)
-    if(Read(fParticles,fTracks))
-     {
-       Error("GetNumberOfTrackEvents","Error in reading");
-       return 0;
-     }
-  return fTracks->GetNumberOfEvents();
- }
+void AliHBTReaderTPC::Rewind()
+{
+  delete fRunLoader;
+  fRunLoader = 0x0;
+  fCurrentDir = 0;
+  fNEventsRead= 0;
+  if (fTrackCounter) fTrackCounter->Reset();
+}
 /********************************************************************/
 
-
-Int_t AliHBTReaderTPC::Read(AliHBTRun* particles, AliHBTRun *tracks)
+Int_t AliHBTReaderTPC::ReadNext()
  {
  //reads data and puts put to the particles and tracks objects
  //reurns 0 if everything is OK
  //
   Info("Read","");
-  Int_t i; //iterator and some temprary values
-  Int_t Nevents = 0;
-  Int_t totalNevents = 0;
-  TFile *aTracksFile;//file with tracks
-  TFile *aClustersFile;//file with clusters
-  TFile *aGAliceFile;//!ile name with galice
-
-  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();
+  
   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
-  Int_t currentdir = 0;
 
-  Int_t Ndirs;
-  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
-   }
+  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( (i=OpenFiles(aTracksFile,aClustersFile,aGAliceFile,currentdir)) )
+      
+    if (fRunLoader == 0x0) 
+      if (OpenNextSession()) 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);
-       currentdir++;
-       continue;
+       //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);
     
-    if (gAlice->TreeE())//check if tree E exists
+    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
      {
-      Nevents = (Int_t)gAlice->TreeE()->GetEntries();//if yes get number of events in gAlice
-      Info("Read","________________________________________________________");
-      Info("Read","Found %d event(s) in directory %s",Nevents,GetDirName(currentdir).Data());
-      Float_t mf;
-      if (fUseMagFFromRun)
-       {
-         mf = gAlice->Field()->SolenoidField();
-         Info("Read","Setting Magnetic Field from run: B=%fT",mf/10.);
-       }
-      else
-       {
-         Info("Read","Setting Own Magnetic Field: B=%fT",fMagneticField);
-         mf = fMagneticField*10.;
-       }
-      AliKalmanTrack::SetConvConst(1000/0.299792458/mf);
+       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
      }
-    else
-     {//if not return an error
-       Error("Read","Can not find Header tree (TreeE) in gAlice");
-       currentdir++;
+
+    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;
      }
-  
-    aGAliceFile->cd();//set cluster file active 
-    AliTPCParam *TPCParam= (AliTPCParam*)aGAliceFile->Get("75x40_100x60_150x60");
-    if (!TPCParam) 
-      { 
-       TPCParam= (AliTPCParam*)aGAliceFile->Get("75x40_100x60");
-       if (!TPCParam) 
+    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;
+       
+       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(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(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(Pass(track))//check if meets all criteria of any of our cuts
+                    //if it does not delete it and take next good track
         { 
-          Error("Read","TPC parameters have not been found !\n");
-          currentdir++;
+          delete track;
+          delete part;
           continue;
         }
+       fParticlesEvent->AddParticle(part);
+       fTracksEvent->AddParticle(track);
       }
-
-  
-    for(Int_t currentEvent =0; currentEvent<Nevents;currentEvent++)//loop over all events
-     {
-       Info("Read","Reading Event %d",currentEvent);
-       /**************************************/
-        /**************************************/
-         /**************************************/ 
-         
-         aTracksFile->cd();//set track file active
-          
-         Char_t  treename[100];
-         sprintf(treename,"TreeT_TPC_%d",currentEvent);//prepare name of the tree
-   
-         TTree *tracktree=0;
-         
-         tracktree=(TTree*)aTracksFile->Get(treename);//get the tree 
-         if (!tracktree) //check if we got the tree
-          {//if not return with error
-            Error("Read","Can't get a tree with TPC tracks !\n"); 
-            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("Read","Can't get a branch with TPC tracks !\n"); 
-            continue;
-          }
-         Int_t NTPCtracks=(Int_t)tracktree->GetEntries();//get number of TPC tracks 
-         Info("Read","Found %d TPC tracks.",NTPCtracks);
-         //Copy tracks to array
-         
-         AliTPCtrack *iotrack=0;
-         
-//         aClustersFile->cd();//set cluster file active 
-   
-         for (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
-          }
-         
-         aTracksFile->Delete(treename);//delete tree from memmory (and leave untached on disk)- we do not need it any more
-         aTracksFile->Delete("tracks");//delete branch from memmory
-         
-         trackbranch = 0x0;
-         tracktree = 0x0;
-   
-         Double_t xk;
-         Double_t par[5];
-         Float_t phi, lam, pt;//angles and transverse momentum
-         Int_t label; //label of the current track
-         
-         aGAliceFile->cd();
-         gAlice->GetEvent(currentEvent); 
-
-         gAlice->Particles();
-         
-         for (i=0; i<NTPCtracks; i++) //loop over all good tracks
-          { 
-            iotrack = (AliTPCtrack*)tarray->At(i);
-            label = iotrack->GetLabel();
-            
-            if (label < 0) continue;
-            
-            
-            TParticle *p = (TParticle*)gAlice->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);
-            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->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);
-          }
-         tarray->Clear(); //clear the array
-         
-        /**************************************/
-       /**************************************/
-      /**************************************/  
-     totalNevents++;
-    }
-    //save environment (resouces) --
-    //clean your place after the work
-    CloseFiles(aTracksFile,aClustersFile,aGAliceFile); 
-    currentdir++;
-   }while(currentdir < Ndirs);
-
+      
+    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;
-  fIsRead = kTRUE;
-  return 0;
+  return 1;
  }
-
 /********************************************************************/
-Int_t AliHBTReaderTPC::OpenFiles
-(TFile*& aTracksFile, TFile*& aClustersFile, TFile*& agAliceFile,Int_t event)
+
+Int_t AliHBTReaderTPC::OpenNextSession()
 {
- //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());
+  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;
+   }
+
+  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);
+
+
+  fRunLoader->CdGAFile();
+  AliTPCParam *TPCParam= (AliTPCParam*)gDirectory->Get("75x40_100x60");
+  if (!TPCParam) 
+   {
+    TPCParam=(AliTPCParam *)gDirectory->Get("75x40_100x60_150x60");
+    if (!TPCParam) 
+     { 
+       DoOpenError("TPC parameters have not been found !\n");
        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;
-    }
-*/
-   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 (!(gAlice=(AliRun*)agAliceFile->Get("gAlice"))) 
-    {
-      Error("OpenFiles","gAlice have not been found on %s !\n",filename.Data());
-      return 5;
-    }
-
-   return 0; 
+  if (fTPCLoader->LoadTracks())
+   {
+     DoOpenError("Error occured while loading TPC tracks.");
+     return 1;
+   }
+  
+  fCurrentEvent = 0;
+  return 0;
 }
 /********************************************************************/
 
-/********************************************************************/
-  
-void AliHBTReaderTPC::CloseFiles(TFile*& tracksFile, TFile*& clustersFile, TFile*& gAliceFile)
+void AliHBTReaderTPC::DoOpenError( const char *va_(fmt), ...)
 {
-  //closes the files
-  tracksFile->Close();
-  delete tracksFile;
-  tracksFile = 0x0;
+  // Does error display and clean-up in case error caught on Open Next Session
 
-//  clustersFile->Close();
-//  delete clustersFile;
-//  clustersFile = 0x0;
+   va_list ap;
+   va_start(ap,va_(fmt));
+   Error("OpenNextSession", va_(fmt), ap);
+   va_end(ap);
+   
+   delete fRunLoader;
+   fRunLoader = 0x0;
+   fTPCLoader = 0x0;
+   fCurrentDir++;
+}
 
-  delete gAlice;
-  gAlice = 0;
+/********************************************************************/
+Bool_t AliHBTReaderTPC::CheckTrack(AliTPCtrack* t)
+{
+  //Performs check of the track
+  if ( (t->GetNumberOfClusters() > fNClustMax) || (t->GetNumberOfClusters() < fNClustMin) ) return kTRUE;
 
-  if (gAliceFile) 
-   {
-     gAliceFile->Close();
-     delete gAliceFile;
-     gAliceFile = 0x0;
-   }
+  Double_t cc[15];
+  t->GetCovariance(cc);
+  if ( (cc[9] < fC44Min) || (cc[9] > fC44Max) ) return kTRUE;
+  
+  Float_t chisqpercl = t->GetChi2()/((Double_t)t->GetNumberOfClusters());
+  if ( (chisqpercl < fNChi2PerClustMin) || (chisqpercl > fNChi2PerClustMax) ) 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::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;
+
+}
+*/
+