]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPC.cxx
Merged with pythia6_address.c
[u/mrichter/AliRoot.git] / TPC / AliTPC.cxx
index a71d5c712ba4957a5caf85a77245d217ca44fbe9..b5fc5e1f94dcfc2a028ab3816ba0211a6d8bbf43 100644 (file)
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
 
-/*
-$Log$
-Revision 1.57  2002/05/07 17:23:11  kowal2
-Linear gain inefficiency instead of the step one at the wire edges.
-
-Revision 1.56  2002/04/04 16:26:33  kowal2
-Digits (Sdigits) go to separate files now.
-
-Revision 1.55  2002/03/29 06:57:45  kowal2
-Restored backward compatibility to use the hits from Dec. 2000 production.
-
-Revision 1.54  2002/03/18 17:59:13  kowal2
-Chnges in the pad geometry - 3 pad lengths introduced.
-
-Revision 1.53  2002/02/25 11:02:56  kowal2
-Changes towards speeding up the code. Thanks to Marian Ivanov.
-
-Revision 1.52  2002/02/18 09:26:09  kowal2
-Removed compiler warning
-
-Revision 1.51  2002/01/21 17:13:21  kowal2
-New track hits using root containers. Setting active sectors added.
-
-Revision 1.50  2001/12/06 14:16:19  kowal2
-meaningfull printouts
-
-Revision 1.49  2001/11/30 11:55:37  hristov
-Noise table created in Hits2SDigits (M.Ivanov)
-
-Revision 1.48  2001/11/24 16:10:21  kowal2
-Faster algorithms.
-
-Revision 1.47  2001/11/19 10:25:34  kowal2
-Nearest integer instead of integer when converting to ADC counts
-
-Revision 1.46  2001/11/07 06:47:12  kowal2
-Removed printouts
-
-Revision 1.45  2001/11/03 13:33:48  kowal2
-Updated algorithms in Hits2SDigits, SDigits2Digits,
-Hits2ExactClusters.
-Added method Merge
-
-Revision 1.44  2001/08/30 09:28:48  hristov
-TTree names are explicitly set via SetName(name) and then Write() is called
-
-Revision 1.43  2001/07/28 12:02:54  hristov
-Branch split level set to 99
-
-Revision 1.42  2001/07/28 11:38:52  hristov
-Loop variable declared once
-
-Revision 1.41  2001/07/28 10:53:50  hristov
-Digitisation done according to the general scheme (M.Ivanov)
-
-Revision 1.40  2001/07/27 13:03:14  hristov
-Default Branch split level set to 99
-
-Revision 1.39  2001/07/26 09:09:34  kowal2
-Hits2Reco method added
-
-Revision 1.38  2001/07/20 14:32:43  kowal2
-Processing of many events possible now
-
-Revision 1.37  2001/06/12 07:17:18  kowal2
-Hits2SDigits method implemented (summable digits)
-
-Revision 1.36  2001/05/16 14:57:25  alibrary
-New files for folders and Stack
-
-Revision 1.35  2001/05/08 16:02:22  kowal2
-Updated material specifications
-
-Revision 1.34  2001/05/08 15:00:15  hristov
-Corrections for tracking in arbitrary magnenetic field. Changes towards a concept of global Alice track. Back propagation of reconstructed tracks (Yu.Belikov)
-
-Revision 1.33  2001/04/03 12:40:43  kowal2
-Removed printouts
-
-Revision 1.32  2001/03/12 17:47:36  hristov
-Changes needed on Sun with CC 5.0
-
-Revision 1.31  2001/03/12 08:21:50  kowal2
-Corrected C++ bug in the material definitions
-
-Revision 1.30  2001/03/01 17:34:47  kowal2
-Correction due to the accuracy problem
-
-Revision 1.29  2001/02/28 16:34:40  kowal2
-Protection against nonphysical values of the avalanche size,
-10**6 is the maximum
-
-Revision 1.28  2001/01/26 19:57:19  hristov
-Major upgrade of AliRoot code
-
-Revision 1.27  2001/01/13 17:29:33  kowal2
-Sun compiler correction
-
-Revision 1.26  2001/01/10 07:59:43  kowal2
-Corrections to load points from the noncompressed hits.
-
-Revision 1.25  2000/11/02 07:25:31  kowal2
-Changes due to the new hit structure.
-Memory leak removed.
-
-Revision 1.24  2000/10/05 16:06:09  kowal2
-Forward declarations. Changes due to a new class AliComplexCluster.
-
-Revision 1.23  2000/10/02 21:28:18  fca
-Removal of useless dependecies via forward declarations
-
-Revision 1.22  2000/07/10 20:57:39  hristov
-Update of TPC code and macros by M.Kowalski
-
-Revision 1.19.2.4  2000/06/26 07:39:42  kowal2
-Changes to obey the coding rules
-
-Revision 1.19.2.3  2000/06/25 08:38:41  kowal2
-Splitted from AliTPCtracking
-
-Revision 1.19.2.2  2000/06/16 12:59:28  kowal2
-Changed parameter settings
-
-Revision 1.19.2.1  2000/06/09 07:15:07  kowal2
-
-Defaults loaded automatically (hard-wired)
-Optional parameters can be set via macro called in the constructor
-
-Revision 1.19  2000/04/18 19:00:59  fca
-Small bug fixes to TPC files
-
-Revision 1.18  2000/04/17 09:37:33  kowal2
-removed obsolete AliTPCDigitsDisplay.C
-
-Revision 1.17.2.2  2000/04/10 08:15:12  kowal2
-
-New, experimental data structure from M. Ivanov
-New tracking algorithm
-Different pad geometry for different sectors
-Digitization rewritten
-
-Revision 1.17.2.1  2000/04/10 07:56:53  kowal2
-Not used anymore - removed
-
-Revision 1.17  2000/01/19 17:17:30  fca
-Introducing a list of lists of hits -- more hits allowed for detector now
-
-Revision 1.16  1999/11/05 09:29:23  fca
-Accept only signals > 0
-
-Revision 1.15  1999/10/08 06:26:53  fca
-Removed ClustersIndex - not used anymore
-
-Revision 1.14  1999/09/29 09:24:33  fca
-Introduction of the Copyright and cvs Log
-
-*/
+/* $Id$ */
 
 ///////////////////////////////////////////////////////////////////////////////
 //                                                                           //
@@ -189,81 +33,52 @@ Introduction of the Copyright and cvs Log
 
 //
 
+#include <Riostream.h>
+#include <stdlib.h>
+
+#include <TFile.h>  
+#include <TGeometry.h>
+#include <TInterpreter.h>
 #include <TMath.h>
-#include <TRandom.h>
-#include <TVector.h>
 #include <TMatrix.h>
-#include <TGeometry.h>
 #include <TNode.h>
-#include <TTUBS.h>
 #include <TObjectTable.h>
-#include "TParticle.h"
-#include "AliTPC.h"
-#include <TFile.h>  
+#include <TParticle.h>
 #include <TROOT.h>
+#include <TRandom.h>
 #include <TSystem.h>     
-#include "AliRun.h"
-#include <iostream.h>
-#include <stdlib.h>
-#include <fstream.h>
-#include "AliMC.h"
-#include "AliMagF.h"
-#include "AliTrackReference.h"
-
+#include <TTUBS.h>
+#include <TTree.h>
+#include <TVirtualMC.h>
+#include <TString.h>
+#include <TF2.h>
+#include <TStopwatch.h>
 
-#include "AliTPCParamSR.h"
-#include "AliTPCPRF2D.h"
-#include "AliTPCRF1D.h"
+#include "AliArrayBranch.h"
 #include "AliDigits.h"
+#include "AliMagF.h"
+#include "AliPoints.h"
+#include "AliRun.h"
+#include "AliRunLoader.h"
 #include "AliSimDigits.h"
+#include "AliTPC.h"
+#include "AliTPC.h"
+#include "AliTPCDigitsArray.h"
+#include "AliTPCLoader.h"
+#include "AliTPCPRF2D.h"
+#include "AliTPCParamSR.h"
+#include "AliTPCRF1D.h"
 #include "AliTPCTrackHits.h"
 #include "AliTPCTrackHitsV2.h"
-#include "AliPoints.h"
-#include "AliArrayBranch.h"
-
-
-#include "AliTPCDigitsArray.h"
-#include "AliComplexCluster.h"
-#include "AliClusters.h"
-#include "AliTPCClustersRow.h"
-#include "AliTPCClustersArray.h"
-
-#include "AliTPCcluster.h"
-#include "AliTPCclusterer.h"
-#include "AliTPCtracker.h"
-
-#include <TInterpreter.h>
-#include <TTree.h>
-
+#include "AliTrackReference.h"
+#include "AliMC.h"
+#include "AliTPCDigitizer.h"
+#include "AliTPCBuffer.h"
+#include "AliTPCDDLRawData.h"
+#include "AliLog.h"
 
 
 ClassImp(AliTPC) 
-
-//_____________________________________________________________________________
-// helper class for fast matrix and vector manipulation - no range checking
-// origin - Marian Ivanov
-
-class AliTPCFastMatrix : public TMatrix {
-public :
-  AliTPCFastMatrix(Int_t row_lwb, Int_t row_upb, Int_t col_lwb, Int_t col_upb);
-  inline Float_t & UncheckedAt(Int_t rown, Int_t coln) const  {return  (fIndex[coln])[rown];} //fast acces   
-  inline Float_t   UncheckedAtFast(Int_t rown, Int_t coln) const  {return  (fIndex[coln])[rown];} //fast acces   
-};
-
-AliTPCFastMatrix::AliTPCFastMatrix(Int_t row_lwb, Int_t row_upb, Int_t col_lwb, Int_t col_upb):
-  TMatrix(row_lwb, row_upb,col_lwb,col_upb)
-   {
-   };
-//_____________________________________________________________________________
-class AliTPCFastVector : public TVector {
-public :
-  AliTPCFastVector(Int_t size);
-  inline Float_t & UncheckedAt(Int_t index) const  {return  fElements[index];} //fast acces  
-};
-
-AliTPCFastVector::AliTPCFastVector(Int_t size):TVector(size){
-};
-
 //_____________________________________________________________________________
 AliTPC::AliTPC()
 {
@@ -275,11 +90,14 @@ AliTPC::AliTPC()
   fDigits   = 0;
   fNsectors = 0;
   fDigitsArray = 0;
-  fClustersArray = 0;
   fDefaults = 0;
   fTrackHits = 0; 
   fTrackHitsOld = 0;   
-  fHitType = 2; //default CONTAINERS - based on ROOT structure 
+#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
+  fHitType = 4; // ROOT containers
+#else
+  fHitType = 2; //default CONTAINERS - based on ROOT structure
+#endif 
   fTPCParam = 0;    
   fNoiseTable = 0;
   fActiveSectors =0;
@@ -297,9 +115,8 @@ AliTPC::AliTPC(const char *name, const char *title)
   //
   // Initialise arrays of hits and digits 
   fHits     = new TClonesArray("AliTPChit",  176);
-  gAlice->AddHitList(fHits); 
+  gAlice->GetMCApp()->AddHitList(fHits); 
   fDigitsArray = 0;
-  fClustersArray= 0;
   fDefaults = 0;
   //
   fTrackHits = new AliTPCTrackHitsV2;  
@@ -314,7 +131,11 @@ AliTPC::AliTPC(const char *name, const char *title)
 
   fNoiseTable =0;
 
+#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
+  fHitType = 4; // ROOT containers
+#else
   fHitType = 2;
+#endif
   fActiveSectors = 0;
   //
   // Initialise counters
@@ -334,13 +155,18 @@ AliTPC::AliTPC(const char *name, const char *title)
   if (!strcmp(title,"Default")) {       
     fTPCParam = new AliTPCParamSR;
   } else {
-    cerr<<"AliTPC warning: in Config.C you must set non-default parameters\n";
+    AliWarning("In Config.C you must set non-default parameters.");
     fTPCParam=0;
   }
 
 }
 
 //_____________________________________________________________________________
+AliTPC::AliTPC(const AliTPC& t):AliDetector(t){
+  //
+  // dummy copy constructor
+  //
+}
 AliTPC::~AliTPC()
 {
   //
@@ -363,31 +189,14 @@ void AliTPC::AddHit(Int_t track, Int_t *vol, Float_t *hits)
   //
   // Add a hit to the list
   //
-  //  TClonesArray &lhits = *fHits;
-  //  new(lhits[fNhits++]) AliTPChit(fIshunt,track,vol,hits);
   if (fHitType&1){
     TClonesArray &lhits = *fHits;
     new(lhits[fNhits++]) AliTPChit(fIshunt,track,vol,hits);
   }
   if (fHitType>1)
-   AddHit2(track,vol,hits);
+    AddHit2(track,vol,hits);
 }
 
-void  AliTPC::AddTrackReference(Int_t lab, TLorentzVector p, TLorentzVector x){
-  //
-  // add a trackrefernce to the list
-  if (!fTrackReferences) {
-    cerr<<"Container trackrefernce not active\n";
-    return;
-  }
-  Int_t nref = fTrackReferences->GetEntriesFast();
-  TClonesArray &lref = *fTrackReferences;
-  AliTrackReference * ref =  new(lref[nref]) AliTrackReference;
-  ref->SetMomentum(p[0],p[1],p[2]);
-  ref->SetPosition(x[0],x[1],x[2]);
-  ref->SetTrack(lab);
-}
 //_____________________________________________________________________________
 void AliTPC::BuildGeometry()
 {
@@ -474,7 +283,7 @@ void AliTPC::BuildGeometry()
 }    
 
 //_____________________________________________________________________________
-Int_t AliTPC::DistancetoPrimitive(Int_t , Int_t )
+Int_t AliTPC::DistancetoPrimitive(Int_t , Int_t ) const
 {
   //
   // Calculate distance from TPC to mouse on the display
@@ -483,13 +292,6 @@ Int_t AliTPC::DistancetoPrimitive(Int_t , Int_t )
   return 9999;
 }
 
-void AliTPC::Clusters2Tracks(TFile *of) {
-  //-----------------------------------------------------------------
-  // This is a track finder.
-  //-----------------------------------------------------------------
-  AliTPCtracker tracker(fTPCParam);
-  tracker.Clusters2Tracks(gFile,of);
-}
 
 //_____________________________________________________________________________
 void AliTPC::CreateMaterials()
@@ -603,42 +405,40 @@ void AliTPC::CreateMaterials()
   // gases - mixtures, ID >= 20 pure gases, <= 10 ID < 20 -compounds
   //----------------------------------------------------------------
 
-  char namate[21]; 
+  char namate[21]=""
   density = 0.;
   Float_t am=0;
   Int_t nc;
-  Float_t rho,absl,X0,buf[1];
+  Float_t rho,absl,x0,buf[1];
   Int_t nbuf;
   Float_t a,z;
 
-  for(nc = 0;nc<fNoComp;nc++)
-    {
+  for(nc = 0;nc<fNoComp;nc++) {
     
-      // retrive material constants
-      
-      gMC->Gfmate((*fIdmate)[fMixtComp[nc]],namate,a,z,rho,X0,absl,buf,nbuf);
-
-      amat[nc] = a;
-      zmat[nc] = z;
-
-      Int_t nnc = (fMixtComp[nc]>=20) ? fMixtComp[nc]%20 : fMixtComp[nc]%10;
-      am += fMixtProp[nc]*((fMixtComp[nc]>=20) ? apure[nnc] : amol[nnc]); 
-      density += fMixtProp[nc]*rho;  // density of the mixture
+    // retrive material constants
       
-    }
+    gMC->Gfmate((*fIdmate)[fMixtComp[nc]],namate,a,z,rho,x0,absl,buf,nbuf);
 
+    amat[nc] = a;
+    zmat[nc] = z;
+    
+    Int_t nnc = (fMixtComp[nc]>=20) ? fMixtComp[nc]%20 : fMixtComp[nc]%10;
+    
+    am += fMixtProp[nc]*((fMixtComp[nc]>=20) ? apure[nnc] : amol[nnc]); 
+    density += fMixtProp[nc]*rho;  // density of the mixture
+    
+  }
+  
   // mixture proportions by weight!
+  
+  for(nc = 0;nc<fNoComp;nc++) {
 
-  for(nc = 0;nc<fNoComp;nc++)
-    {
-
-      Int_t nnc = (fMixtComp[nc]>=20) ? fMixtComp[nc]%20 : fMixtComp[nc]%10;
+    Int_t nnc = (fMixtComp[nc]>=20) ? fMixtComp[nc]%20 : fMixtComp[nc]%10;
 
-      wmat[nc] = fMixtProp[nc]*((fMixtComp[nc]>=20) ? 
-                 apure[nnc] : amol[nnc])/am;
+    wmat[nc] = fMixtProp[nc]*((fMixtComp[nc]>=20) ? 
+                             apure[nnc] : amol[nnc])/am;
 
-    
+  } 
 
   // Drift gases 1 - nonsensitive, 2 - sensitive
 
@@ -846,7 +646,7 @@ void AliTPC::CreateMaterials()
 
   // get epoxy
 
-  gMC->Gfmate((*fIdmate)[45],namate,amat[1],zmat[1],rho,X0,absl,buf,nbuf);
+  gMC->Gfmate((*fIdmate)[45],namate,amat[1],zmat[1],rho,x0,absl,buf,nbuf);
 
   // Carbon fiber
 
@@ -859,7 +659,7 @@ void AliTPC::CreateMaterials()
 
   // get SiO2
 
-  gMC->Gfmate((*fIdmate)[38],namate,amat[0],zmat[0],rho,X0,absl,buf,nbuf); 
+  gMC->Gfmate((*fIdmate)[38],namate,amat[0],zmat[0],rho,x0,absl,buf,nbuf); 
 
   wmat[0]=0.725; // by weight!
   wmat[1]=0.275;
@@ -929,7 +729,7 @@ Float_t AliTPC::GetNoise()
 }
 
 
-Bool_t  AliTPC::IsSectorActive(Int_t sec)
+Bool_t  AliTPC::IsSectorActive(Int_t sec) const
 {
   //
   // check if the sector is active
@@ -940,6 +740,7 @@ Bool_t  AliTPC::IsSectorActive(Int_t sec)
 void    AliTPC::SetActiveSectors(Int_t * sectors, Int_t n)
 {
   // activate interesting sectors
+  SetTreeAddress();//just for security
   if (fActiveSectors) delete [] fActiveSectors;
   fActiveSectors = new Bool_t[fTPCParam->GetNSector()];
   for (Int_t i=0;i<fTPCParam->GetNSector();i++) fActiveSectors[i]=kFALSE;
@@ -953,6 +754,7 @@ void    AliTPC::SetActiveSectors(Int_t flag)
   //
   // activate sectors which were hitted by tracks 
   //loop over tracks
+  SetTreeAddress();//just for security
   if (fHitType==0) return;  // if Clones hit - not short volume ID information
   if (fActiveSectors) delete [] fActiveSectors;
   fActiveSectors = new Bool_t[fTPCParam->GetNSector()];
@@ -962,16 +764,23 @@ void    AliTPC::SetActiveSectors(Int_t flag)
   }
   for (Int_t i=0;i<fTPCParam->GetNSector();i++) fActiveSectors[i]=kFALSE;
   TBranch * branch=0;
-  if (fHitType>1) branch = gAlice->TreeH()->GetBranch("TPC2");
-  else branch = gAlice->TreeH()->GetBranch("TPC");
-  Stat_t ntracks = gAlice->TreeH()->GetEntries();
+  if (TreeH() == 0x0)
+   {
+     AliFatal("Can not find TreeH in folder");
+     return;
+   }
+  if (fHitType>1) branch = TreeH()->GetBranch("TPC2");
+  else branch = TreeH()->GetBranch("TPC");
+  Stat_t ntracks = TreeH()->GetEntries();
   // loop over all hits
-  for(Int_t track=0;track<ntracks;track++){
+  AliDebug(1,Form("Got %d tracks",ntracks));
+  
+  for(Int_t track=0;track<ntracks;track++) {
     ResetHits();
     //
     if (fTrackHits && fHitType&4) {
-      TBranch * br1 = gAlice->TreeH()->GetBranch("fVolumes");
-      TBranch * br2 = gAlice->TreeH()->GetBranch("fNVolumes");    
+      TBranch * br1 = TreeH()->GetBranch("fVolumes");
+      TBranch * br2 = TreeH()->GetBranch("fNVolumes");
       br1->GetEvent(track);
       br2->GetEvent(track);
       Int_t *volumes = fTrackHits->GetVolumes();
@@ -981,7 +790,7 @@ void    AliTPC::SetActiveSectors(Int_t flag)
     
     //
     if (fTrackHitsOld && fHitType&2) {
-      TBranch * br = gAlice->TreeH()->GetBranch("fTrackHitsInfo");
+      TBranch * br = TreeH()->GetBranch("fTrackHitsInfo");
       br->GetEvent(track);
       AliObjectArray * ar = fTrackHitsOld->fTrackHitsInfo;
       for (UInt_t j=0;j<ar->GetSize();j++){
@@ -989,415 +798,146 @@ void    AliTPC::SetActiveSectors(Int_t flag)
       } 
     }    
   }
-  
 }  
 
 
 
 
-void AliTPC::Digits2Clusters(TFile *of, Int_t eventnumber)
-{
-  //-----------------------------------------------------------------
-  // This is a simple cluster finder.
-  //-----------------------------------------------------------------
-  AliTPCclusterer::Digits2Clusters(fTPCParam,of,eventnumber);
-}
-
-extern Double_t SigmaY2(Double_t, Double_t, Double_t);
-extern Double_t SigmaZ2(Double_t, Double_t);
 //_____________________________________________________________________________
-void AliTPC::Hits2Clusters(TFile *of, Int_t eventn)
+void AliTPC::Digits2Raw()
 {
-  //--------------------------------------------------------
-  // TPC simple cluster generator from hits
-  // obtained from the TPC Fast Simulator
-  // The point errors are taken from the parametrization
-  //--------------------------------------------------------
-
-  //-----------------------------------------------------------------
-  // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
-  //-----------------------------------------------------------------
-  // Adopted to Marian's cluster data structure by I.Belikov, CERN,
-  // Jouri.Belikov@cern.ch
-  //----------------------------------------------------------------
-  
-  /////////////////////////////////////////////////////////////////////////////
-  //
-  //---------------------------------------------------------------------
-  //   ALICE TPC Cluster Parameters
-  //--------------------------------------------------------------------
-       
-  
+// convert digits of the current event to raw data
 
-  // Cluster width in rphi
-  const Float_t kACrphi=0.18322;
-  const Float_t kBCrphi=0.59551e-3;
-  const Float_t kCCrphi=0.60952e-1;
-  // Cluster width in z
-  const Float_t kACz=0.19081;
-  const Float_t kBCz=0.55938e-3;
-  const Float_t kCCz=0.30428;
+  static const Int_t kThreshold = 0;
+  static const Bool_t kCompress = kTRUE;
 
-  TDirectory *savedir=gDirectory; 
-
-  if (!of->IsOpen()) {
-     cerr<<"AliTPC::Hits2Clusters(): output file not open !\n";
-     return;
+  fLoader->LoadDigits();
+  TTree* digits = fLoader->TreeD();
+  if (!digits) {
+    AliError("No digits tree");
+    return;
   }
 
-  //if(fDefaults == 0) SetDefaults();
-
-  Float_t sigmaRphi,sigmaZ,clRphi,clZ;
-  //
-  TParticle *particle; // pointer to a given particle
-  AliTPChit *tpcHit; // pointer to a sigle TPC hit
-  Int_t sector;
-  Int_t ipart;
-  Float_t xyz[5];
-  Float_t pl,pt,tanth,rpad,ratio;
-  Float_t cph,sph;
-  
-  //---------------------------------------------------------------
-  //  Get the access to the tracks 
-  //---------------------------------------------------------------
-  
-  TTree *tH = gAlice->TreeH();
-  Stat_t ntracks = tH->GetEntries();
-
-  //Switch to the output file
-  of->cd();
-
-  char   cname[100];
-
-  sprintf(cname,"TreeC_TPC_%d",eventn);
-
-  fTPCParam->Write(fTPCParam->GetTitle());
-  AliTPCClustersArray carray;
-  carray.Setup(fTPCParam);
-  carray.SetClusterType("AliTPCcluster");
-  carray.MakeTree();
-
-  Int_t nclusters=0; //cluster counter
-  
-  //------------------------------------------------------------
-  // Loop over all sectors (72 sectors for 20 deg
-  // segmentation for both lower and upper sectors)
-  // Sectors 0-35 are lower sectors, 0-17 z>0, 17-35 z<0
-  // Sectors 36-71 are upper sectors, 36-53 z>0, 54-71 z<0
-  //
-  // First cluster for sector 0 starts at "0"
-  //------------------------------------------------------------
-   
-  for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++){
-    //MI change
-    fTPCParam->AdjustCosSin(isec,cph,sph);
-    
-    //------------------------------------------------------------
-    // Loop over tracks
-    //------------------------------------------------------------
-    
-    for(Int_t track=0;track<ntracks;track++){
-      ResetHits();
-      tH->GetEvent(track);
-      //
-      //  Get number of the TPC hits
-      //     
-       tpcHit = (AliTPChit*)FirstHit(-1);
-
-      // Loop over hits
-      //
-       while(tpcHit){
-        if (tpcHit->fQ == 0.) {
-           tpcHit = (AliTPChit*) NextHit();
-           continue; //information about track (I.Belikov)
-        }
-       sector=tpcHit->fSector; // sector number
-
-       if(sector != isec){
-        tpcHit = (AliTPChit*) NextHit();
-        continue; 
-       }
-       ipart=tpcHit->Track();
-       particle=gAlice->Particle(ipart);
-       pl=particle->Pz();
-       pt=particle->Pt();
-       if(pt < 1.e-9) pt=1.e-9;
-       tanth=pl/pt;
-       tanth = TMath::Abs(tanth);
-       rpad=TMath::Sqrt(tpcHit->X()*tpcHit->X() + tpcHit->Y()*tpcHit->Y());
-       ratio=0.001*rpad/pt; // pt must be in MeV/c - historical reason
-
-       //   space-point resolutions
-       
-       sigmaRphi=SigmaY2(rpad,tanth,pt);
-       sigmaZ   =SigmaZ2(rpad,tanth   );
-       
-       //   cluster widths
-       
-       clRphi=kACrphi-kBCrphi*rpad*tanth+kCCrphi*ratio*ratio;
-       clZ=kACz-kBCz*rpad*tanth+kCCz*tanth*tanth;
-       
-       // temporary protection
-       
-       if(sigmaRphi < 0.) sigmaRphi=0.4e-3;
-       if(sigmaZ < 0.) sigmaZ=0.4e-3;
-       if(clRphi < 0.) clRphi=2.5e-3;
-       if(clZ < 0.) clZ=2.5e-5;
-       
-       //
-       
-       //
-       // smearing --> rotate to the 1 (13) or to the 25 (49) sector,
-       // then the inaccuracy in a X-Y plane is only along Y (pad row)!
-       //
-        Float_t xprim= tpcHit->X()*cph + tpcHit->Y()*sph;
-       Float_t yprim=-tpcHit->X()*sph + tpcHit->Y()*cph;
-       xyz[0]=gRandom->Gaus(yprim,TMath::Sqrt(sigmaRphi));   // y
-          Float_t alpha=(isec < fTPCParam->GetNInnerSector()) ?
-         fTPCParam->GetInnerAngle() : fTPCParam->GetOuterAngle();
-          Float_t ymax=xprim*TMath::Tan(0.5*alpha);
-          if (TMath::Abs(xyz[0])>ymax) xyz[0]=yprim; 
-       xyz[1]=gRandom->Gaus(tpcHit->Z(),TMath::Sqrt(sigmaZ)); // z
-          if (TMath::Abs(xyz[1])>fTPCParam->GetZLength()) xyz[1]=tpcHit->Z(); 
-       xyz[2]=sigmaRphi;                                     // fSigmaY2
-       xyz[3]=sigmaZ;                                        // fSigmaZ2
-       xyz[4]=tpcHit->fQ;                                    // q
-
-        AliTPCClustersRow *clrow=carray.GetRow(sector,tpcHit->fPadRow);
-        if (!clrow) clrow=carray.CreateRow(sector,tpcHit->fPadRow);    
-
-        Int_t tracks[3]={tpcHit->Track(), -1, -1};
-       AliTPCcluster cluster(tracks,xyz);
-
-        clrow->InsertCluster(&cluster); nclusters++;
-
-        tpcHit = (AliTPChit*)NextHit();
-        
-
-      } // end of loop over hits
-
-    }   // end of loop over tracks
-
-    Int_t nrows=fTPCParam->GetNRow(isec);
-    for (Int_t irow=0; irow<nrows; irow++) {
-        AliTPCClustersRow *clrow=carray.GetRow(isec,irow);
-        if (!clrow) continue;
-        carray.StoreRow(isec,irow);
-        carray.ClearRow(isec,irow);
+  AliSimDigits digarr;
+  AliSimDigits* digrow = &digarr;
+  digits->GetBranch("Segment")->SetAddress(&digrow);
+
+  const char* fileName = "AliTPCDDL.dat";
+  AliTPCBuffer* buffer  = new AliTPCBuffer(fileName);
+  //Verbose level
+  // 0: Silent
+  // 1: cout messages
+  // 2: txt files with digits 
+  //BE CAREFUL, verbose level 2 MUST be used only for debugging and
+  //it is highly suggested to use this mode only for debugging digits files
+  //reasonably small, because otherwise the size of the txt files can reach
+  //quickly several MB wasting time and disk space.
+  buffer->SetVerbose(0);
+
+  Int_t nEntries = Int_t(digits->GetEntries());
+  Int_t previousSector = -1;
+  Int_t subSector = 0;
+  for (Int_t i = 0; i < nEntries; i++) {
+    digits->GetEntry(i);
+    Int_t sector, row;
+    fTPCParam->AdjustSectorRow(digarr.GetID(), sector, row);
+    if(previousSector != sector) {
+      subSector = 0;
+      previousSector = sector;
     }
 
-  } // end of loop over sectors  
-
-  cerr<<"Number of made clusters : "<<nclusters<<"                        \n";
-  carray.GetTree()->SetName(cname);
-  carray.GetTree()->Write();
-
-  savedir->cd(); //switch back to the input file
-  
-} // end of function
-
-//_________________________________________________________________
-void AliTPC::Hits2ExactClustersSector(Int_t isec)
-{
-  //--------------------------------------------------------
-  //calculate exact cross point of track and given pad row
-  //resulting values are expressed in "digit" coordinata
-  //--------------------------------------------------------
-
-  //-----------------------------------------------------------------
-  // Origin: Marian Ivanov  GSI Darmstadt, m.ivanov@gsi.de
-  //-----------------------------------------------------------------
-  //
-  if (fClustersArray==0){    
-    return;
-  }
-  //
-  TParticle *particle; // pointer to a given particle
-  AliTPChit *tpcHit; // pointer to a sigle TPC hit
-  //  Int_t sector,nhits;
-  Int_t ipart;
-  const Int_t kcmaxhits=30000;
-  AliTPCFastVector * xxxx = new AliTPCFastVector(kcmaxhits*4);
-  AliTPCFastVector & xxx = *xxxx;
-  Int_t maxhits = kcmaxhits;
-  //construct array for each padrow
-  for (Int_t i=0; i<fTPCParam->GetNRow(isec);i++) 
-    fClustersArray->CreateRow(isec,i);
-  
-  //---------------------------------------------------------------
-  //  Get the access to the tracks 
-  //---------------------------------------------------------------
-  
-  TTree *tH = gAlice->TreeH();
-  Stat_t ntracks = tH->GetEntries();
-  Int_t npart = gAlice->GetNtrack();
-  //MI change
-  TBranch * branch=0;
-  if (fHitType>1) branch = tH->GetBranch("TPC2");
-  else branch = tH->GetBranch("TPC");
-
-  //------------------------------------------------------------
-  // Loop over tracks
-  //------------------------------------------------------------
+    if (sector < 36) { //inner sector [0;35]
+      if (row != 30) {
+       //the whole row is written into the output file
+       buffer->WriteRowBinary(kThreshold, digrow, 0, 0, 0, 
+                              sector, subSector, row);
+      } else {
+       //only the pads in the range [37;48] are written into the output file
+       buffer->WriteRowBinary(kThreshold, digrow, 37, 48, 1, 
+                              sector, subSector, row);
+       subSector = 1;
+       //only the pads outside the range [37;48] are written into the output file
+       buffer->WriteRowBinary(kThreshold, digrow, 37, 48, 2, 
+                              sector, subSector, row);
+      }//end else
+
+    } else { //outer sector [36;71]
+      if (row == 54) subSector = 2;
+      if ((row != 27) && (row != 76)) {
+       buffer->WriteRowBinary(kThreshold, digrow, 0, 0, 0,
+                              sector, subSector, row);
+      } else if (row == 27) {
+       //only the pads outside the range [43;46] are written into the output file
+       buffer->WriteRowBinary(kThreshold, digrow, 43, 46, 2,
+                                sector, subSector, row);
+       subSector = 1;
+       //only the pads in the range [43;46] are written into the output file
+       buffer->WriteRowBinary(kThreshold, digrow, 43, 46, 1,
+                                sector, subSector, row);
+      } else if (row == 76) {
+       //only the pads outside the range [33;88] are written into the output file
+       buffer->WriteRowBinary(kThreshold, digrow, 33, 88, 2,
+                              sector, subSector, row);
+       subSector = 3;
+       //only the pads in the range [33;88] are written into the output file
+       buffer->WriteRowBinary(kThreshold, digrow, 33, 88, 1,
+                                sector, subSector, row);
+      }
+    }//end else
+  }//end for
 
-  for(Int_t track=0;track<ntracks;track++){ 
-    Bool_t isInSector=kTRUE;
-    ResetHits();
-     isInSector = TrackInVolume(isec,track);
-    if (!isInSector) continue;
-    //MI change
-    branch->GetEntry(track); // get next track
-    //
-    //  Get number of the TPC hits and a pointer
-    //  to the particles
-    // Loop over hits
-    //
-    Int_t currentIndex=0;
-    Int_t lastrow=-1;  //last writen row
+  delete buffer;
+  fLoader->UnloadDigits();
 
-    //M.I. changes
+  AliTPCDDLRawData rawWriter;
+  rawWriter.SetVerbose(0);
 
-    tpcHit = (AliTPChit*)FirstHit(-1);
-    while(tpcHit){
-      
-      Int_t sector=tpcHit->fSector; // sector number
-      if(sector != isec){
-       tpcHit = (AliTPChit*) NextHit();
-       continue; 
-      }
+  rawWriter.RawData(fileName);
+  gSystem->Unlink(fileName);
 
-      ipart=tpcHit->Track();
-      if (ipart<npart) particle=gAlice->Particle(ipart);
-      
-      //find row number
-
-      Float_t  x[3]={tpcHit->X(),tpcHit->Y(),tpcHit->Z()};
-      Int_t    index[3]={1,isec,0};
-      Int_t    currentrow = fTPCParam->GetPadRow(x,index) ;    
-      if (currentrow<0) {tpcHit = (AliTPChit*)NextHit(); continue;}
-      if (lastrow<0) lastrow=currentrow;
-      if (currentrow==lastrow){
-       if ( currentIndex>=maxhits){
-         maxhits+=kcmaxhits;
-         xxx.ResizeTo(4*maxhits);
-       }     
-       xxx(currentIndex*4)=x[0];
-       xxx(currentIndex*4+1)=x[1];
-       xxx(currentIndex*4+2)=x[2];     
-       xxx(currentIndex*4+3)=tpcHit->fQ;
-       currentIndex++; 
-      }
-      else 
-       if (currentIndex>2){
-         Float_t sumx=0;
-         Float_t sumx2=0;
-         Float_t sumx3=0;
-         Float_t sumx4=0;
-         Float_t sumy=0;
-         Float_t sumxy=0;
-         Float_t sumx2y=0;
-         Float_t sumz=0;
-         Float_t sumxz=0;
-         Float_t sumx2z=0;
-         Float_t sumq=0;
-         for (Int_t index=0;index<currentIndex;index++){
-           Float_t x,x2,x3,x4;
-           x=x2=x3=x4=xxx(index*4);
-           x2*=x;
-           x3*=x2;
-           x4*=x3;
-           sumx+=x;
-           sumx2+=x2;
-           sumx3+=x3;
-           sumx4+=x4;
-           sumy+=xxx(index*4+1);
-           sumxy+=xxx(index*4+1)*x;
-           sumx2y+=xxx(index*4+1)*x2;
-           sumz+=xxx(index*4+2);
-           sumxz+=xxx(index*4+2)*x;
-           sumx2z+=xxx(index*4+2)*x2;   
-           sumq+=xxx(index*4+3);
-         }
-         Float_t centralPad = (fTPCParam->GetNPads(isec,lastrow)-1)/2;
-         Float_t det=currentIndex*(sumx2*sumx4-sumx3*sumx3)-sumx*(sumx*sumx4-sumx2*sumx3)+
-           sumx2*(sumx*sumx3-sumx2*sumx2);
-         
-         Float_t detay=sumy*(sumx2*sumx4-sumx3*sumx3)-sumx*(sumxy*sumx4-sumx2y*sumx3)+
-           sumx2*(sumxy*sumx3-sumx2y*sumx2);
-         Float_t detaz=sumz*(sumx2*sumx4-sumx3*sumx3)-sumx*(sumxz*sumx4-sumx2z*sumx3)+
-           sumx2*(sumxz*sumx3-sumx2z*sumx2);
-         
-         Float_t detby=currentIndex*(sumxy*sumx4-sumx2y*sumx3)-sumy*(sumx*sumx4-sumx2*sumx3)+
-           sumx2*(sumx*sumx2y-sumx2*sumxy);
-         Float_t detbz=currentIndex*(sumxz*sumx4-sumx2z*sumx3)-sumz*(sumx*sumx4-sumx2*sumx3)+
-           sumx2*(sumx*sumx2z-sumx2*sumxz);
-         
-         if (TMath::Abs(det)<0.00001){
-            tpcHit = (AliTPChit*)NextHit();
-           continue;
-         }
-       
-         Float_t y=detay/det+centralPad;
-         Float_t z=detaz/det;  
-         Float_t by=detby/det; //y angle
-         Float_t bz=detbz/det; //z angle
-         sumy/=Float_t(currentIndex);
-         sumz/=Float_t(currentIndex);
-
-         AliTPCClustersRow * row = (fClustersArray->GetRow(isec,lastrow));
-         if (row!=0) {
-           AliComplexCluster* cl = new((AliComplexCluster*)row->Append()) AliComplexCluster ;
-           //    AliComplexCluster cl;
-           cl->fX=z;
-           cl->fY=y;
-           cl->fQ=sumq;
-           cl->fSigmaX2=bz;
-           cl->fSigmaY2=by;
-           cl->fTracks[0]=ipart;
-         }
-         currentIndex=0;
-         lastrow=currentrow;
-       } //end of calculating cluster for given row
-       
-       
-      tpcHit = (AliTPChit*)NextHit();
-    } // end of loop over hits
-  }   // end of loop over tracks 
-  //write padrows to tree 
-  for (Int_t ii=0; ii<fTPCParam->GetNRow(isec);ii++) {
-    fClustersArray->StoreRow(isec,ii);    
-    fClustersArray->ClearRow(isec,ii);        
+  if (kCompress) {
+    AliInfo("Compressing raw data");
+    rawWriter.RawDataCompDecompress(kTRUE);
+    gSystem->Unlink("Statistics");
   }
-  xxxx->Delete();
 }
 
 
 
+//______________________________________________________________________
+AliDigitizer* AliTPC::CreateDigitizer(AliRunDigitizer* manager) const
+{
+  return new AliTPCDigitizer(manager);
+}
 //__
-void AliTPC::SDigits2Digits2(Int_t eventnumber)  
+void AliTPC::SDigits2Digits2(Int_t /*eventnumber*/)  
 {
   //create digits from summable digits
   GenerNoise(500000); //create teble with noise
-  char  sname[100];
-  char  dname[100];
-  sprintf(sname,"TreeS_%s_%d",fTPCParam->GetTitle(),eventnumber);
-  sprintf(dname,"TreeD_%s_%d",fTPCParam->GetTitle(),eventnumber);
 
   //conect tree with sSDigits
-  TTree *t;
-  if (gAlice->GetTreeDFile()) {
-    t = (TTree *)gAlice->GetTreeDFile()->Get(sname); 
-  } else {
-    t = (TTree *)gDirectory->Get(sname); 
-  }
-  if (!t) {
-    cerr<<"TPC tree with sdigits not found"<<endl;
-    return;
+  TTree *t = fLoader->TreeS();
+
+  if (t == 0x0) {
+    fLoader->LoadSDigits("READ");
+    t = fLoader->TreeS();
+    if (t == 0x0) {
+      AliError("Can not get input TreeS");
+      return;
+    }
   }
+  
+  if (fLoader->TreeD() == 0x0) fLoader->MakeTree("D");
+  
   AliSimDigits digarr, *dummy=&digarr;
-  t->GetBranch("Segment")->SetAddress(&dummy);
+  TBranch* sdb = t->GetBranch("Segment");
+  if (sdb == 0x0) {
+    AliError("Can not find branch with segments in TreeS.");
+    return;
+  }  
+
+  sdb->SetAddress(&dummy);
+      
   Stat_t nentries = t->GetEntries();
 
   // set zero suppression
@@ -1413,14 +953,9 @@ void AliTPC::SDigits2Digits2(Int_t eventnumber)
   AliTPCDigitsArray *arr = new AliTPCDigitsArray; 
   arr->SetClass("AliSimDigits");
   arr->Setup(fTPCParam);
-// Note that methods arr->MakeTree have different signatures
-  if (gAlice->GetTreeDFile()) {
-    arr->MakeTree(gAlice->GetTreeDFile());
-  } else {
-    arr->MakeTree(fDigitsFile);
-  }
+  arr->MakeTree(fLoader->TreeD());
   
-  AliTPCParam * par =fTPCParam;
+  AliTPCParam * par = fTPCParam;
 
   //Loop over segments of the TPC
 
@@ -1428,10 +963,11 @@ void AliTPC::SDigits2Digits2(Int_t eventnumber)
     t->GetEvent(n);
     Int_t sec, row;
     if (!par->AdjustSectorRow(digarr.GetID(),sec,row)) {
-      cerr<<"AliTPC warning: invalid segment ID ! "<<digarr.GetID()<<endl;
+      AliWarning(Form("Invalid segment ID ! %d",digarr.GetID()));
       continue;
     }
     if (!IsSectorActive(sec)) continue;
+    
     AliSimDigits * digrow =(AliSimDigits*) arr->CreateRow(sec,row);
     Int_t nrows = digrow->GetNRows();
     Int_t ncols = digrow->GetNCols();
@@ -1452,18 +988,11 @@ void AliTPC::SDigits2Digits2(Int_t eventnumber)
     //if you cahnge implementation
     //of the Alidigits - it must be rewriten -
     for (Int_t i= 0; i<nelems; i++){
-      //      Float_t q = *pamp0;
-      //q/=16.;  //conversion faktor
-      //Float_t noise= GetNoise(); 
-      //q+=noise;      
-      //q= TMath::Nint(q);
       Float_t q = TMath::Nint(Float_t(*pamp0)/16.+GetNoise());
       if (q>zerosup){
        if (q>saturation) q=saturation;      
        *pamp1=(Short_t)q;
-       //if (ptracks0[0]==0)
-       //  ptracks1[0]=1;
-       //else
+
        ptracks1[0]=ptracks0[0];        
        ptracks1[nelems]=ptracks0[nelems];
        ptracks1[2*nelems]=ptracks0[2*nelems];
@@ -1476,157 +1005,28 @@ void AliTPC::SDigits2Digits2(Int_t eventnumber)
 
     arr->StoreRow(sec,row);
     arr->ClearRow(sec,row);   
-    // cerr<<sec<<"\t"<<row<<"\n";   
   }  
 
     
   //write results
-
-  
-  arr->GetTree()->SetName(dname);  
-  arr->GetTree()->AutoSave();  
-  delete arr;
-}
-//_________________________________________
-void AliTPC::Merge(TTree * intree, Int_t *mask, Int_t nin, Int_t outid)
-{
-  
-  //intree - pointer to array of trees with s digits
-  //mask   - mask for each 
-  //nin    - number of inputs
-  //outid  - event  number of the output event
-
-  //create digits from summable digits 
-  //conect tree with sSDigits
-
-    
-  AliSimDigits ** digarr = new AliSimDigits*[nin];
-  for (Int_t i1=0;i1<nin; i1++){
-    digarr[i1]=0;
-    intree[i1].GetBranch("Segment")->SetAddress(&digarr[i1]);
-  }
-  Stat_t nentries = intree[0].GetEntries();
-  
-  //make tree with digits   
-  char  dname[100];
-  sprintf(dname,"TreeD_%s_%d",fTPCParam->GetTitle(),outid);
-  AliTPCDigitsArray *arr = new AliTPCDigitsArray; 
-  arr->SetClass("AliSimDigits");
-  arr->Setup(fTPCParam);
-  arr->MakeTree(fDigitsFile);  
-
-  // set zero suppression
-
-  fTPCParam->SetZeroSup(2);
-
-  // get zero suppression
-
-  Int_t zerosup = fTPCParam->GetZeroSup();
-
-
-  AliTPCParam * par =fTPCParam;
-
-  //Loop over segments of the TPC
-  for (Int_t n=0; n<nentries; n++) {
-    
-    for (Int_t i=0;i<nin; i++){ 
-      //connect proper digits
-      intree[i].GetEvent(n);      
-      digarr[i]->ExpandBuffer();
-      digarr[i]->ExpandTrackBuffer();
-    }      
-    Int_t sec, row;
-    if (!par->AdjustSectorRow(digarr[0]->GetID(),sec,row)) {
-      cerr<<"AliTPC warning: invalid segment ID ! "<<digarr[0]->GetID()<<endl;
-      continue;
-    }
-
-    AliSimDigits * digrow =(AliSimDigits*) arr->CreateRow(sec,row);
-    Int_t nrows = digrow->GetNRows();
-    Int_t ncols = digrow->GetNCols();
-
-    digrow->ExpandBuffer();
-    digrow->ExpandTrackBuffer();
-   
-    for (Int_t rows=0;rows<nrows; rows++){
-      for (Int_t col=0;col<ncols; col++){
-        Float_t q=0;
-        Int_t label[1000]; // i hope no more than 300 events merged
-        Int_t labptr = 0;
-        // looop over digits
-        for (Int_t i=0;i<nin; i++){ 
-          q  += digarr[i]->GetDigitFast(rows,col);
-          for (Int_t tr=0;tr<3;tr++) {
-            Int_t lab = digarr[i]->GetTrackIDFast(rows,col,tr);
-            if ( lab > 1) {
-              label[labptr]=lab+mask[i]-2;
-              labptr++;
-            }
-          }
-        }
-        //add noise
-        q = gRandom->Gaus(q,fTPCParam->GetNoise()*fTPCParam->GetNoiseNormFac()*16); 
-        
-        q/=16;  //conversion faktor
-        q=(Short_t)q;
-
-        if (q> zerosup){
-
-          if(q > fTPCParam->GetADCSat()) q = fTPCParam->GetADCSat();
-          digrow->SetDigitFast((Short_t)q,rows,col);  
-          for (Int_t tr=0;tr<3;tr++){
-            if (tr<labptr)
-              ((AliSimDigits*)digrow)->SetTrackIDFast(label[tr],rows,col,tr);
-            else
-              ((AliSimDigits*)digrow)->SetTrackIDFast(-1,rows,col,tr);
-          }
-        }
-      } 
-     }         
-      arr->StoreRow(sec,row);
-
-      arr->ClearRow(sec,row);   
-  }  
-  
-  delete digarr;
-  arr->GetTree()->SetName(dname); 
-  arr->GetTree()->Write(); 
-  delete arr;  
-  
-}
-
-/*_________________________________________
-void AliTPC::SDigits2Digits(Int_t eventnumber)
-{
-
-
-  cerr<<"Digitizing TPC...\n";
-
-  Hits2Digits(eventnumber);
+  fLoader->WriteDigits("OVERWRITE");
    
-    
-  //write results
-
-  //  char treeName[100];
-
-  //  sprintf(treeName,"TreeD_%s_%d",fTPCParam->GetTitle(),eventnumber);
-  
-  //  GetDigitsArray()->GetTree()->Write(treeName);  
+  delete arr;
 }
-*/
 //__________________________________________________________________
 void AliTPC::SetDefaults(){
-
+  //
+  // setting the defaults
+  //
    
-   cerr<<"Setting default parameters...\n";
-
   // Set response functions
 
+  //
+  AliRunLoader* rl = (AliRunLoader*)fLoader->GetEventFolder()->FindObject(AliRunLoader::GetRunLoaderName());
+  rl->CdGAFile();
   AliTPCParamSR *param=(AliTPCParamSR*)gDirectory->Get("75x40_100x60");
   if(param){
-    printf("You are using 2 pad-length geom hits with 3 pad-lenght geom digits...\n");
+    AliInfo("You are using 2 pad-length geom hits with 3 pad-lenght geom digits...");
     delete param;
     param = new AliTPCParamSR();
   }
@@ -1634,8 +1034,7 @@ void AliTPC::SetDefaults(){
     param=(AliTPCParamSR*)gDirectory->Get("75x40_100x60_150x60");
   }
   if(!param){
-    printf("No TPC parameters found\n");
-    exit(4);
+    AliFatal("No TPC parameters found");
   }
 
 
@@ -1649,13 +1048,26 @@ void AliTPC::SetDefaults(){
   
   TDirectory *savedir=gDirectory;
   TFile *f=TFile::Open("$ALICE_ROOT/TPC/AliTPCprf2d.root");
-  if (!f->IsOpen()) 
-    cerr<<"Can't open $ALICE_ROOT/TPC/AliTPCprf2d.root !\n" ;
-     exit(3);
-  }
+  if (!f->IsOpen()) 
+    AliFatal("Can't open $ALICE_ROOT/TPC/AliTPCprf2d.root !");
+
+  TString s;
   prfinner->Read("prf_07504_Gati_056068_d02");
+  //PH Set different names
+  s=prfinner->GetGRF()->GetName();
+  s+="in";
+  prfinner->GetGRF()->SetName(s.Data());
+
   prfouter1->Read("prf_10006_Gati_047051_d03");
+  s=prfouter1->GetGRF()->GetName();
+  s+="out1";
+  prfouter1->GetGRF()->SetName(s.Data());
+
   prfouter2->Read("prf_15006_Gati_047051_d03");  
+  s=prfouter2->GetGRF()->GetName();
+  s+="out2";
+  prfouter2->GetGRF()->SetName(s.Data());
+
   f->Close();
   savedir->cd();
 
@@ -1673,12 +1085,47 @@ void AliTPC::SetDefaults(){
 
 }
 //__________________________________________________________________  
+void AliTPC::Hits2Digits()  
+{
+  //
+  // creates digits from hits
+  //
+
+  fLoader->LoadHits("read");
+  fLoader->LoadDigits("recreate");
+  AliRunLoader* runLoader = fLoader->GetRunLoader(); 
+
+  for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
+    runLoader->GetEvent(iEvent);
+    SetActiveSectors();   
+    Hits2Digits(iEvent);
+  }
+
+  fLoader->UnloadHits();
+  fLoader->UnloadDigits();
+} 
+//__________________________________________________________________  
 void AliTPC::Hits2Digits(Int_t eventnumber)  
 { 
  //----------------------------------------------------
  // Loop over all sectors for a single event
  //----------------------------------------------------
-
+  AliRunLoader* rl = (AliRunLoader*)fLoader->GetEventFolder()->FindObject(AliRunLoader::GetRunLoaderName());
+  rl->GetEvent(eventnumber);
+  if (fLoader->TreeH() == 0x0) {
+    if(fLoader->LoadHits()) {
+      AliError("Can not load hits.");
+    }
+  }
+  SetTreeAddress();
+  
+  if (fLoader->TreeD() == 0x0 ) {
+    fLoader->MakeTree("D");
+    if (fLoader->TreeD() == 0x0 ) {
+      AliError("Can not get TreeD");
+      return;
+    }
+  }
 
   if(fDefaults == 0) SetDefaults();  // check if the parameters are set
   GenerNoise(500000); //create teble with noise
@@ -1686,38 +1133,36 @@ void AliTPC::Hits2Digits(Int_t eventnumber)
   //setup TPCDigitsArray 
 
   if(GetDigitsArray()) delete GetDigitsArray();
-
+  
   AliTPCDigitsArray *arr = new AliTPCDigitsArray; 
   arr->SetClass("AliSimDigits");
   arr->Setup(fTPCParam);
-// Note that methods arr->MakeTree have different signatures
-  if (gAlice->GetTreeDFile()) {
-    arr->MakeTree(gAlice->GetTreeDFile());
-  } else {
-    arr->MakeTree(fDigitsFile);
-  }
+
+  arr->MakeTree(fLoader->TreeD());
   SetDigitsArray(arr);
 
   fDigitsSwitch=0; // standard digits
 
-  cerr<<"Digitizing TPC -- normal digits...\n";
-
-  for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++) if (IsSectorActive(isec)) Hits2DigitsSector(isec); 
-   
-  // write results
-
-  char treeName[100];
-
-  sprintf(treeName,"TreeD_%s_%d",fTPCParam->GetTitle(),eventnumber);
+  for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++) 
+    if (IsSectorActive(isec)) {
+      AliDebug(1,Form("Hits2Digits","Sector %d is active.",isec));
+      Hits2DigitsSector(isec);
+    }
+    else {
+      AliDebug(1,Form("Hits2Digits","Sector %d is NOT active.",isec));
+    }
+  
+  fLoader->WriteDigits("OVERWRITE"); 
+  
+//this line prevents the crash in the similar one
+//on the beginning of this method
+//destructor attempts to reset the tree, which is deleted by the loader
+//need to be redesign
+  if(GetDigitsArray()) delete GetDigitsArray();
+  SetDigitsArray(0x0);
   
-  GetDigitsArray()->GetTree()->SetName(treeName);  
-  GetDigitsArray()->GetTree()->AutoSave();  
-
-
 }
 
-
-
 //__________________________________________________________________
 void AliTPC::Hits2SDigits2(Int_t eventnumber)  
 { 
@@ -1726,29 +1171,40 @@ void AliTPC::Hits2SDigits2(Int_t eventnumber)
   //   summable digits - 16 bit "ADC", no noise, no saturation
   //-----------------------------------------------------------
 
- //----------------------------------------------------
- // Loop over all sectors for a single event
- //----------------------------------------------------
+  //----------------------------------------------------
+  // Loop over all sectors for a single event
+  //----------------------------------------------------
+
+  AliRunLoader* rl = fLoader->GetRunLoader();
 
+  rl->GetEvent(eventnumber);
+  if (fLoader->TreeH() == 0x0) {
+    if(fLoader->LoadHits()) {
+      AliError("Can not load hits.");
+      return;
+    }
+  }
+  SetTreeAddress();
 
+
+  if (fLoader->TreeS() == 0x0 ) {
+    fLoader->MakeTree("S");
+  }
+  
   if(fDefaults == 0) SetDefaults();
+  
   GenerNoise(500000); //create table with noise
   //setup TPCDigitsArray 
 
   if(GetDigitsArray()) delete GetDigitsArray();
 
+  
   AliTPCDigitsArray *arr = new AliTPCDigitsArray; 
   arr->SetClass("AliSimDigits");
   arr->Setup(fTPCParam);
-// Note that methods arr->MakeTree have different signatures
-  if (gAlice->GetTreeSFile()) {
-    arr->MakeTree(gAlice->GetTreeSFile());
-  } else {
-    arr->MakeTree(fDigitsFile);
-  }
-  SetDigitsArray(arr);
+  arr->MakeTree(fLoader->TreeS());
 
-  cerr<<"Digitizing TPC -- summable digits...\n"; 
+  SetDigitsArray(arr);
 
   fDigitsSwitch=1; // summable digits
   
@@ -1756,23 +1212,22 @@ void AliTPC::Hits2SDigits2(Int_t eventnumber)
 
   fTPCParam->SetZeroSup(0);
 
- for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++) if (IsSectorActive(isec)) Hits2DigitsSector(isec);
-
-
-  // write results
-
-  char treeName[100];
+  for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++) 
+    if (IsSectorActive(isec)) {
+      Hits2DigitsSector(isec);
+    }
 
-  sprintf(treeName,"TreeS_%s_%d",fTPCParam->GetTitle(),eventnumber);
-  
-  GetDigitsArray()->GetTree()->SetName(treeName); 
-  GetDigitsArray()->GetTree()->AutoSave(); 
+  fLoader->WriteSDigits("OVERWRITE");
 
+//this line prevents the crash in the similar one
+//on the beginning of this method
+//destructor attempts to reset the tree, which is deleted by the loader
+//need to be redesign
+  if(GetDigitsArray()) delete GetDigitsArray();
+  SetDigitsArray(0x0);
 }
-
-
-
 //__________________________________________________________________
+
 void AliTPC::Hits2SDigits()  
 { 
 
@@ -1780,49 +1235,22 @@ void AliTPC::Hits2SDigits()
   //   summable digits - 16 bit "ADC", no noise, no saturation
   //-----------------------------------------------------------
 
- //----------------------------------------------------
- // Loop over all sectors for a single event
- //----------------------------------------------------
-  //MI change - for pp run
-  Int_t eventnumber = gAlice->GetEvNumber();
-
-  if(fDefaults == 0) SetDefaults();
-  GenerNoise(500000); //create table with noise
-
-  //setup TPCDigitsArray 
-
-  if(GetDigitsArray()) delete GetDigitsArray();
+  fLoader->LoadHits("read");
+  fLoader->LoadSDigits("recreate");
+  AliRunLoader* runLoader = fLoader->GetRunLoader(); 
 
-  AliTPCDigitsArray *arr = new AliTPCDigitsArray; 
-  arr->SetClass("AliSimDigits");
-  arr->Setup(fTPCParam);
-// Note that methods arr->MakeTree have different signatures
-  if (gAlice->GetTreeSFile()) {
-    arr->MakeTree(gAlice->GetTreeSFile());
-  } else {
-    arr->MakeTree(fDigitsFile);
+  for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
+    runLoader->GetEvent(iEvent);
+    SetTreeAddress();
+    SetActiveSectors();
+    Hits2SDigits2(iEvent);
   }
-  SetDigitsArray(arr);
-
-  cerr<<"Digitizing TPC -- summable digits...\n"; 
-
-  //  fDigitsSwitch=1; // summable digits  -for the moment direct
-
-  for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++) if (IsSectorActive(isec)) Hits2DigitsSector(isec);
-
-
-  // write results
-  char treeName[100];
-
-  sprintf(treeName,"TreeD_%s_%d",fTPCParam->GetTitle(),eventnumber);
-  
-  GetDigitsArray()->GetTree()->SetName(treeName); 
-  GetDigitsArray()->GetTree()->AutoSave(); 
 
+  fLoader->UnloadHits();
+  fLoader->UnloadSDigits();
 }
-
-
 //_____________________________________________________________________________
+
 void AliTPC::Hits2DigitsSector(Int_t isec)
 {
   //-------------------------------------------------------------------
@@ -1842,7 +1270,12 @@ void AliTPC::Hits2DigitsSector(Int_t isec)
 
   if(fDefaults == 0) SetDefaults();
 
-  TTree *tH = gAlice->TreeH(); // pointer to the hits tree
+  TTree *tH = TreeH(); // pointer to the hits tree
+  if (tH == 0x0) {
+    AliFatal("Can not find TreeH in folder");
+    return;
+  }
+
   Stat_t ntracks = tH->GetEntries();
 
   if( ntracks > 0){
@@ -1853,47 +1286,50 @@ void AliTPC::Hits2DigitsSector(Int_t isec)
 
     TObjArray **row;
     
-    //printf("*** Processing sector number %d ***\n",isec);
-
-      Int_t nrows =fTPCParam->GetNRow(isec);
+    Int_t nrows =fTPCParam->GetNRow(isec);
 
-      row= new TObjArray* [nrows];
+    row= new TObjArray* [nrows+2]; // 2 extra rows for cross talk
     
-      MakeSector(isec,nrows,tH,ntracks,row);
+    MakeSector(isec,nrows,tH,ntracks,row);
 
-      //--------------------------------------------------------
-      //   Digitize this sector, row by row
-      //   row[i] is the pointer to the TObjArray of AliTPCFastVectors,
-      //   each one containing electrons accepted on this
-      //   row, assigned into tracks
-      //--------------------------------------------------------
+    //--------------------------------------------------------
+    //   Digitize this sector, row by row
+    //   row[i] is the pointer to the TObjArray of TVectors,
+    //   each one containing electrons accepted on this
+    //   row, assigned into tracks
+    //--------------------------------------------------------
 
-      Int_t i;
+    Int_t i;
 
-      if (fDigitsArray->GetTree()==0) fDigitsArray->MakeTree(fDigitsFile);
-
-      for (i=0;i<nrows;i++){
+    if (fDigitsArray->GetTree()==0) {
+      AliFatal("Tree not set in fDigitsArray");
+    }
 
-       AliDigits * dig = fDigitsArray->CreateRow(isec,i); 
+    for (i=0;i<nrows;i++){
+      
+      AliDigits * dig = fDigitsArray->CreateRow(isec,i); 
 
-       DigitizeRow(i,isec,row);
+      DigitizeRow(i,isec,row);
 
-       fDigitsArray->StoreRow(isec,i);
+      fDigitsArray->StoreRow(isec,i);
 
-       Int_t ndig = dig->GetDigitSize(); 
-       if (gDebug > 10) printf("*** Sector, row, compressed digits %d %d %d ***\n",isec,i,ndig);        
+      Int_t ndig = dig->GetDigitSize(); 
+       
+      AliDebug(10,
+              Form("*** Sector, row, compressed digits %d %d %d ***\n",
+                   isec,i,ndig));        
        
-        fDigitsArray->ClearRow(isec,i);  
+      fDigitsArray->ClearRow(isec,i);  
 
    
-       } // end of the sector digitization
+    } // end of the sector digitization
 
-      for(i=0;i<nrows;i++){
-        row[i]->Delete();  
-        delete row[i];   
-      }
+    for(i=0;i<nrows+2;i++){
+      row[i]->Delete();  
+      delete row[i];   
+    }
       
-       delete [] row; // delete the array of pointers to TObjArray-s
+    delete [] row; // delete the array of pointers to TObjArray-s
         
   } // ntracks >0
 
@@ -1913,9 +1349,8 @@ void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
   // Modified: Marian Ivanov GSI Darmstadt, m.ivanov@gsi.de
   //-----------------------------------------------------------------
  
-
   Float_t zerosup = fTPCParam->GetZeroSup();
-  Int_t nrows =fTPCParam->GetNRow(isec);
+
   fCurrentIndex[1]= isec;
   
 
@@ -1927,10 +1362,10 @@ void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
   //  and a single track signal
   //    
 
-  AliTPCFastMatrix *m1 = new AliTPCFastMatrix(0,nofPads,0,nofTbins); // integrated
-  AliTPCFastMatrix *m2 = new AliTPCFastMatrix(0,nofPads,0,nofTbins); // single
+  TMatrix *m1 = new TMatrix(0,nofPads,0,nofTbins); // integrated
+  TMatrix *m2 = new TMatrix(0,nofPads,0,nofTbins); // single
   //
-  AliTPCFastMatrix &total  = *m1;
+  TMatrix &total  = *m1;
 
   //  Array of pointers to the label-signal list
 
@@ -1943,14 +1378,14 @@ void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
   //
   //calculate signal 
   //
-  Int_t row1 = TMath::Max(irow-fTPCParam->GetNCrossRows(),0);
-  Int_t row2 = TMath::Min(irow+fTPCParam->GetNCrossRows(),nrows-1);
+  Int_t row1=irow;
+  Int_t row2=irow+2; 
   for (Int_t row= row1;row<=row2;row++){
     Int_t nTracks= rows[row]->GetEntries();
     for (i1=0;i1<nTracks;i1++){
       fCurrentIndex[2]= row;
-      fCurrentIndex[3]=irow;
-      if (row==irow){
+      fCurrentIndex[3]=irow+1;
+      if (row==irow+1){
        m2->Zero();  // clear single track signal matrix
        Float_t trackLabel = GetSignal(rows[row],i1,m2,m1,indexRange); 
        GetList(trackLabel,nofPads,m2,indexRange,pList);
@@ -1965,11 +1400,9 @@ void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
   Int_t gi=-1;
   Float_t fzerosup = zerosup+0.5;
   for(Int_t it=0;it<nofTbins;it++){
-    Float_t *pq = &(total.UncheckedAt(0,it));
     for(Int_t ip=0;ip<nofPads;ip++){
       gi++;
-      Float_t q=*pq;      
-      pq++;
+      Float_t q=total(ip,it);      
       if(fDigitsSwitch == 0){
        q+=GetNoise();
         if(q <=fzerosup) continue; // do not fill zeros
@@ -1979,10 +1412,10 @@ void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
       }
 
       else {
-       if(q <= 0.) continue; // do not fill zeros
-       if(q>2000.) q=2000.;
-       q *= 16.;
-       q = TMath::Nint(q);
+       if(q <= 0.) continue; // do not fill zeros
+       if(q>2000.) q=2000.;
+       q *= 16.;
+       q = TMath::Nint(q);
       }
 
       //
@@ -2000,13 +1433,11 @@ void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
 */
 //End_Html
       dig->SetDigitFast((Short_t)q,it,ip);
-      if (fDigitsArray->IsSimulated())
-       {
-        ((AliSimDigits*)dig)->SetTrackIDFast(tracks[0],it,ip,0);
-        ((AliSimDigits*)dig)->SetTrackIDFast(tracks[1],it,ip,1);
-        ((AliSimDigits*)dig)->SetTrackIDFast(tracks[2],it,ip,2);
-       }
-     
+      if (fDigitsArray->IsSimulated()) {
+       ((AliSimDigits*)dig)->SetTrackIDFast(tracks[0],it,ip,0);
+       ((AliSimDigits*)dig)->SetTrackIDFast(tracks[1],it,ip,1);
+       ((AliSimDigits*)dig)->SetTrackIDFast(tracks[2],it,ip,2);
+      }
     
     } // end of loop over time buckets
   }  // end of lop over pads 
@@ -2023,14 +1454,13 @@ void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
 
   delete m1;
   delete m2;
-  //  delete m3;
 
 } // end of DigitizeRow
 
 //_____________________________________________________________________________
 
 Float_t AliTPC::GetSignal(TObjArray *p1, Int_t ntr, 
-             AliTPCFastMatrix *m1, AliTPCFastMatrix *m2,Int_t *indexRange)
+             TMatrix *m1, TMatrix *m2,Int_t *indexRange)
 {
 
   //---------------------------------------------------------------
@@ -2044,13 +1474,13 @@ Float_t AliTPC::GetSignal(TObjArray *p1, Int_t ntr,
   // Modified: Marian Ivanov 
   //-----------------------------------------------------------------
 
-  AliTPCFastVector *tv;
+  TVector *tv;
 
-  tv = (AliTPCFastVector*)p1->At(ntr); // pointer to a track
-  AliTPCFastVector &v = *tv;
+  tv = (TVector*)p1->At(ntr); // pointer to a track
+  TVector &v = *tv;
   
   Float_t label = v(0);
-  Int_t centralPad = (fTPCParam->GetNPads(fCurrentIndex[1],fCurrentIndex[3])-1)/2;
+  Int_t centralPad = (fTPCParam->GetNPads(fCurrentIndex[1],fCurrentIndex[3]-1)-1)/2;
 
   Int_t nElectrons = (tv->GetNrows()-1)/4;
   indexRange[0]=9999; // min pad
@@ -2058,8 +1488,8 @@ Float_t AliTPC::GetSignal(TObjArray *p1, Int_t ntr,
   indexRange[2]=9999; //min time
   indexRange[3]=-1; // max time
 
-  AliTPCFastMatrix &signal = *m1;
-  AliTPCFastMatrix &total = *m2;
+  TMatrix &signal = *m1;
+  TMatrix &total = *m2;
   //
   //  Loop over all electrons
   //
@@ -2074,23 +1504,23 @@ Float_t AliTPC::GetSignal(TObjArray *p1, Int_t ntr,
     Float_t *weight = & (fTPCParam->GetResWeight(0));
 
     if (n>0) for (Int_t i =0; i<n; i++){       
-       Int_t pad=index[1]+centralPad;  //in digit coordinates central pad has coordinate 0
-
-         if (pad>=0){
-        Int_t time=index[2];    
-         Float_t qweight = *(weight)*eltoadcfac;
-        
-        if (m1!=0) signal.UncheckedAt(pad,time)+=qweight;
-         total.UncheckedAt(pad,time)+=qweight;
-        if (indexRange[0]>pad) indexRange[0]=pad;
-        if (indexRange[1]<pad) indexRange[1]=pad;
-        if (indexRange[2]>time) indexRange[2]=time;
-        if (indexRange[3]<time) indexRange[3]=time;
-
-        index+=3;
-        weight++;      
-
-       }        
+      Int_t pad=index[1]+centralPad;  //in digit coordinates central pad has coordinate 0
+
+      if (pad>=0){
+       Int_t time=index[2];     
+       Float_t qweight = *(weight)*eltoadcfac;
+       
+       if (m1!=0) signal(pad,time)+=qweight;
+       total(pad,time)+=qweight;
+       if (indexRange[0]>pad) indexRange[0]=pad;
+       if (indexRange[1]<pad) indexRange[1]=pad;
+       if (indexRange[2]>time) indexRange[2]=time;
+       if (indexRange[3]<time) indexRange[3]=time;
+       
+       index+=3;
+       weight++;       
+
+      }         
     }
   } // end of loop over electrons
   
@@ -2098,7 +1528,7 @@ Float_t AliTPC::GetSignal(TObjArray *p1, Int_t ntr,
 }
 
 //_____________________________________________________________________________
-void AliTPC::GetList(Float_t label,Int_t np,AliTPCFastMatrix *m,
+void AliTPC::GetList(Float_t label,Int_t np,TMatrix *m,
                      Int_t *indexRange, Float_t **pList)
 {
   //----------------------------------------------------------------------
@@ -2109,7 +1539,7 @@ void AliTPC::GetList(Float_t label,Int_t np,AliTPCFastMatrix *m,
   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
   //-----------------------------------------------------------------
 
-  AliTPCFastMatrix &signal = *m;
+  TMatrix &signal = *m;
 
   // lop over nonzero digits
 
@@ -2117,77 +1547,73 @@ void AliTPC::GetList(Float_t label,Int_t np,AliTPCFastMatrix *m,
     for(Int_t ip=indexRange[0];ip<indexRange[1]+1;ip++){
 
 
-        // accept only the contribution larger than 500 electrons (1/2 s_noise)
+      // accept only the contribution larger than 500 electrons (1/2 s_noise)
 
-        if(signal(ip,it)<0.5) continue; 
+      if(signal(ip,it)<0.5) continue; 
 
-
-        Int_t globalIndex = it*np+ip; // globalIndex starts from 0!
+      Int_t globalIndex = it*np+ip; // globalIndex starts from 0!
         
-        if(!pList[globalIndex]){
+      if(!pList[globalIndex]){
         
-          // 
-         // Create new list (6 elements - 3 signals and 3 labels),
-         //
-
-          pList[globalIndex] = new Float_t [6];
-
-         // set list to -1 
-
-          *pList[globalIndex] = -1.;
-          *(pList[globalIndex]+1) = -1.;
-          *(pList[globalIndex]+2) = -1.;
-          *(pList[globalIndex]+3) = -1.;
-          *(pList[globalIndex]+4) = -1.;
-          *(pList[globalIndex]+5) = -1.;
-
-
-          *pList[globalIndex] = label;
-          *(pList[globalIndex]+3) = signal(ip,it);
-        }
-        else{
-
-         // check the signal magnitude
-
-          Float_t highest = *(pList[globalIndex]+3);
-          Float_t middle = *(pList[globalIndex]+4);
-          Float_t lowest = *(pList[globalIndex]+5);
-
-         //
-         //  compare the new signal with already existing list
-         //
+       // 
+       // Create new list (6 elements - 3 signals and 3 labels),
+       //
 
-          if(signal(ip,it)<lowest) continue; // neglect this track
+       pList[globalIndex] = new Float_t [6];
 
-         //
+       // set list to -1 
+       
+       *pList[globalIndex] = -1.;
+       *(pList[globalIndex]+1) = -1.;
+       *(pList[globalIndex]+2) = -1.;
+       *(pList[globalIndex]+3) = -1.;
+       *(pList[globalIndex]+4) = -1.;
+       *(pList[globalIndex]+5) = -1.;
+
+       *pList[globalIndex] = label;
+       *(pList[globalIndex]+3) = signal(ip,it);
+      }
+      else {
 
-          if (signal(ip,it)>highest){
-            *(pList[globalIndex]+5) = middle;
-            *(pList[globalIndex]+4) = highest;
-            *(pList[globalIndex]+3) = signal(ip,it);
+       // check the signal magnitude
 
-            *(pList[globalIndex]+2) = *(pList[globalIndex]+1);
-            *(pList[globalIndex]+1) = *pList[globalIndex];
-            *pList[globalIndex] = label;
-         }
-          else if (signal(ip,it)>middle){
-            *(pList[globalIndex]+5) = middle;
-            *(pList[globalIndex]+4) = signal(ip,it);
+       Float_t highest = *(pList[globalIndex]+3);
+       Float_t middle = *(pList[globalIndex]+4);
+       Float_t lowest = *(pList[globalIndex]+5);
+       
+       //
+       //  compare the new signal with already existing list
+       //
+       
+       if(signal(ip,it)<lowest) continue; // neglect this track
 
-            *(pList[globalIndex]+2) = *(pList[globalIndex]+1);
-            *(pList[globalIndex]+1) = label;
-         }
-          else{
-            *(pList[globalIndex]+5) = signal(ip,it);
-            *(pList[globalIndex]+2) = label;
-         }
-        }
+       //
 
+       if (signal(ip,it)>highest){
+         *(pList[globalIndex]+5) = middle;
+         *(pList[globalIndex]+4) = highest;
+         *(pList[globalIndex]+3) = signal(ip,it);
+         
+         *(pList[globalIndex]+2) = *(pList[globalIndex]+1);
+         *(pList[globalIndex]+1) = *pList[globalIndex];
+         *pList[globalIndex] = label;
+       }
+       else if (signal(ip,it)>middle){
+         *(pList[globalIndex]+5) = middle;
+         *(pList[globalIndex]+4) = signal(ip,it);
+         
+         *(pList[globalIndex]+2) = *(pList[globalIndex]+1);
+         *(pList[globalIndex]+1) = label;
+       }
+       else{
+         *(pList[globalIndex]+5) = signal(ip,it);
+         *(pList[globalIndex]+2) = label;
+       }
+      }
+      
     } // end of loop over pads
   } // end of loop over time bins
 
-
-
 }//end of GetList
 //___________________________________________________________________
 void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
@@ -2217,14 +1643,14 @@ void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
  
   //----------------------------------------------
   // Create TObjArray-s, one for each row,
-  // each TObjArray will store the AliTPCFastVectors
-  // of electrons, one AliTPCFastVectors per each track.
+  // each TObjArray will store the TVectors
+  // of electrons, one TVectors per each track.
   //---------------------------------------------- 
     
-  Int_t *nofElectrons = new Int_t [nrows]; // electron counter for each row
-  AliTPCFastVector **tracks = new AliTPCFastVector* [nrows]; //pointers to the track vectors
+  Int_t *nofElectrons = new Int_t [nrows+2]; // electron counter for each row
+  TVector **tracks = new TVector* [nrows+2]; //pointers to the track vectors
 
-  for(i=0; i<nrows; i++){
+  for(i=0; i<nrows+2; i++){
     row[i] = new TObjArray;
     nofElectrons[i]=0;
     tracks[i]=0;
@@ -2235,7 +1661,7 @@ void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
   //--------------------------------------------------------------------
   //  Loop over tracks, the "track" contains the full history
   //--------------------------------------------------------------------
-
+  
   Int_t previousTrack,currentTrack;
   previousTrack = -1; // nothing to store so far!
 
@@ -2246,7 +1672,7 @@ void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
     if (!isInSector) continue;
     //MI change
     branch->GetEntry(track); // get next track
-
+    
     //M.I. changes
 
     tpcHit = (AliTPChit*)FirstHit(-1);
@@ -2264,133 +1690,131 @@ void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
        continue; 
       }
 
-       currentTrack = tpcHit->Track(); // track number
-
+      currentTrack = tpcHit->Track(); // track number
 
-        if(currentTrack != previousTrack){
+      if(currentTrack != previousTrack){
                           
-           // store already filled fTrack
+       // store already filled fTrack
               
-          for(i=0;i<nrows;i++){
-             if(previousTrack != -1){
-              if(nofElectrons[i]>0){
-                 AliTPCFastVector &v = *tracks[i];
-                v(0) = previousTrack;
-                 tracks[i]->ResizeTo(4*nofElectrons[i]+1); // shrink if necessary
-                row[i]->Add(tracks[i]);                     
-              }
-               else{
-                 delete tracks[i]; // delete empty AliTPCFastVector
-                 tracks[i]=0;
-              }
-            }
-
-             nofElectrons[i]=0;
-             tracks[i] = new AliTPCFastVector(481); // AliTPCFastVectors for the next fTrack
-
-          } // end of loop over rows
+       for(i=0;i<nrows+2;i++){
+         if(previousTrack != -1){
+           if(nofElectrons[i]>0){
+             TVector &v = *tracks[i];
+             v(0) = previousTrack;
+             tracks[i]->ResizeTo(4*nofElectrons[i]+1); // shrink if necessary
+             row[i]->Add(tracks[i]);                     
+           }
+           else {
+             delete tracks[i]; // delete empty TVector
+             tracks[i]=0;
+           }
+         }
+
+         nofElectrons[i]=0;
+         tracks[i] = new TVector(481); // TVectors for the next fTrack
+
+       } // end of loop over rows
               
-           previousTrack=currentTrack; // update track label 
-       }
+       previousTrack=currentTrack; // update track label 
+      }
           
-       Int_t qI = (Int_t) (tpcHit->fQ); // energy loss (number of electrons)
+      Int_t qI = (Int_t) (tpcHit->fQ); // energy loss (number of electrons)
 
-       //---------------------------------------------------
-       //  Calculate the electron attachment probability
-       //---------------------------------------------------
+      //---------------------------------------------------
+      //  Calculate the electron attachment probability
+      //---------------------------------------------------
 
 
-        Float_t time = 1.e6*(fTPCParam->GetZLength()-TMath::Abs(tpcHit->Z()))
-                                                        /fTPCParam->GetDriftV(); 
-       // in microseconds!     
-       Float_t attProb = fTPCParam->GetAttCoef()*
-         fTPCParam->GetOxyCont()*time; //  fraction! 
+      Float_t time = 1.e6*(fTPCParam->GetZLength()-TMath::Abs(tpcHit->Z()))
+       /fTPCParam->GetDriftV(); 
+      // in microseconds!      
+      Float_t attProb = fTPCParam->GetAttCoef()*
+       fTPCParam->GetOxyCont()*time; //  fraction! 
    
-       //-----------------------------------------------
-       //  Loop over electrons
-       //-----------------------------------------------
-       Int_t index[3];
-       index[1]=isec;
-        for(Int_t nel=0;nel<qI;nel++){
-          // skip if electron lost due to the attachment
-          if((gRandom->Rndm(0)) < attProb) continue; // electron lost!
-         xyz[0]=tpcHit->X();
-         xyz[1]=tpcHit->Y();
-         xyz[2]=tpcHit->Z();   
-         //
-         // protection for the nonphysical avalanche size (10**6 maximum)
-         //  
-          Double_t rn=TMath::Max(gRandom->Rndm(0),1.93e-22);
-         xyz[3]= (Float_t) (-gasgain*TMath::Log(rn)); 
-         index[0]=1;
-         
-         TransportElectron(xyz,index); //MI change -august       
-         Int_t rowNumber;
-         fTPCParam->GetPadRow(xyz,index); //MI change august
-         rowNumber = index[2];
-         //transform position to local digit coordinates
-         //relative to nearest pad row 
-         if ((rowNumber<0)||rowNumber>=fTPCParam->GetNRow(isec)) continue;
-  Float_t x1,y1;
-         if (isec <fTPCParam->GetNInnerSector()) {
-           x1 = xyz[1]*fTPCParam->GetInnerPadPitchWidth();
-           y1 = fTPCParam->GetYInner(rowNumber);
-         }
-         else{
-           x1=xyz[1]*fTPCParam->GetOuterPadPitchWidth();
-           y1 = fTPCParam->GetYOuter(rowNumber);
-         }
-
-         // gain inefficiency at the wires edges - linear
-
-         x1=TMath::Abs(x1);
-         y1-=1.;
-          if(x1>y1) xyz[3]*=TMath::Max(1.e-6,(y1-x1+1.));      
-       
-         nofElectrons[rowNumber]++;      
-         //----------------------------------
-         // Expand vector if necessary
-         //----------------------------------
-         if(nofElectrons[rowNumber]>120){
-           Int_t range = tracks[rowNumber]->GetNrows();
-           if((nofElectrons[rowNumber])>(range-1)/4){
-        
-             tracks[rowNumber]->ResizeTo(range+400); // Add 100 electrons
-           }
+      //-----------------------------------------------
+      //  Loop over electrons
+      //-----------------------------------------------
+      Int_t index[3];
+      index[1]=isec;
+      for(Int_t nel=0;nel<qI;nel++){
+       // skip if electron lost due to the attachment
+       if((gRandom->Rndm(0)) < attProb) continue; // electron lost!
+       xyz[0]=tpcHit->X();
+       xyz[1]=tpcHit->Y();
+       xyz[2]=tpcHit->Z();     
+       //
+       // protection for the nonphysical avalanche size (10**6 maximum)
+       //  
+       Double_t rn=TMath::Max(gRandom->Rndm(0),1.93e-22);
+       xyz[3]= (Float_t) (-gasgain*TMath::Log(rn)); 
+       index[0]=1;
+       
+       TransportElectron(xyz,index);    
+       Int_t rowNumber;
+       fTPCParam->GetPadRow(xyz,index); 
+       // row 0 - cross talk from the innermost row
+       // row fNRow+1 cross talk from the outermost row
+       rowNumber = index[2]+1; 
+       //transform position to local digit coordinates
+       //relative to nearest pad row 
+       if ((rowNumber<0)||rowNumber>fTPCParam->GetNRow(isec)+1) continue;
+       Float_t x1,y1;
+       if (isec <fTPCParam->GetNInnerSector()) {
+         x1 = xyz[1]*fTPCParam->GetInnerPadPitchWidth();
+         y1 = fTPCParam->GetYInner(rowNumber);
+       }
+       else{
+         x1=xyz[1]*fTPCParam->GetOuterPadPitchWidth();
+         y1 = fTPCParam->GetYOuter(rowNumber);
+       }
+       // gain inefficiency at the wires edges - linear
+       x1=TMath::Abs(x1);
+       y1-=1.;
+       if(x1>y1) xyz[3]*=TMath::Max(1.e-6,(y1-x1+1.)); 
+       
+       nofElectrons[rowNumber]++;        
+       //----------------------------------
+       // Expand vector if necessary
+       //----------------------------------
+       if(nofElectrons[rowNumber]>120){
+         Int_t range = tracks[rowNumber]->GetNrows();
+         if((nofElectrons[rowNumber])>(range-1)/4){
+           
+           tracks[rowNumber]->ResizeTo(range+400); // Add 100 electrons
          }
-         
-          AliTPCFastVector &v = *tracks[rowNumber];
-         Int_t idx = 4*nofElectrons[rowNumber]-3;
-         Real_t * position = &(((AliTPCFastVector&)v).UncheckedAt(idx)); //make code faster
-         memcpy(position,xyz,4*sizeof(Float_t));
-       } // end of loop over electrons
+       }
+       
+       TVector &v = *tracks[rowNumber];
+       Int_t idx = 4*nofElectrons[rowNumber]-3;
+       Real_t * position = &(((TVector&)v)(idx)); //make code faster
+       memcpy(position,xyz,4*sizeof(Float_t));
+       
+      } // end of loop over electrons
 
-        tpcHit = (AliTPChit*)NextHit();
-        
-      } // end of loop over hits
-    } // end of loop over tracks
+      tpcHit = (AliTPChit*)NextHit();
+      
+    } // end of loop over hits
+  } // end of loop over tracks
 
     //
     //   store remaining track (the last one) if not empty
     //
-
-     for(i=0;i<nrows;i++){
-       if(nofElectrons[i]>0){
-          AliTPCFastVector &v = *tracks[i];
-         v(0) = previousTrack;
-          tracks[i]->ResizeTo(4*nofElectrons[i]+1); // shrink if necessary
-         row[i]->Add(tracks[i]);  
-       }
-       else{
-          delete tracks[i];
-          tracks[i]=0;
-       }  
-      }  
-
-          delete [] tracks;
-          delete [] nofElectrons;
+  
+  for(i=0;i<nrows+2;i++){
+    if(nofElectrons[i]>0){
+      TVector &v = *tracks[i];
+      v(0) = previousTrack;
+      tracks[i]->ResizeTo(4*nofElectrons[i]+1); // shrink if necessary
+      row[i]->Add(tracks[i]);  
+    }
+    else{
+      delete tracks[i];
+      tracks[i]=0;
+    }  
+  }  
+  
+  delete [] tracks;
+  delete [] nofElectrons;
 
 } // end of MakeSector
 
@@ -2401,40 +1825,33 @@ void AliTPC::Init()
   //
   // Initialise TPC detector after definition of geometry
   //
-  Int_t i;
-  //
-  if(fDebug) {
-    printf("\n%s: ",ClassName());
-    for(i=0;i<35;i++) printf("*");
-    printf(" TPC_INIT ");
-    for(i=0;i<35;i++) printf("*");
-    printf("\n%s: ",ClassName());
-    //
-    for(i=0;i<80;i++) printf("*");
-    printf("\n");
-  }
+  AliDebug(1,"*********************************************");
 }
 
 //_____________________________________________________________________________
-void AliTPC::MakeBranch(Option_t* option, const char *file)
+void AliTPC::MakeBranch(Option_t* option)
 {
   //
   // Create Tree branches for the TPC.
   //
+  AliDebug(1,"");
   Int_t buffersize = 4000;
   char branchname[10];
   sprintf(branchname,"%s",GetName());
-
-  AliDetector::MakeBranch(option,file);
-
+  
+  const char *h = strstr(option,"H");
+  
+  if ( h && (fHitType<=1) && (fHits == 0x0)) fHits = new TClonesArray("AliTPChit", 176);//skowron 20.06.03
+  
+  AliDetector::MakeBranch(option);
+  
   const char *d = strstr(option,"D");
-
-  if (fDigits   && gAlice->TreeD() && d) {
-      MakeBranchInTree(gAlice->TreeD(), 
-                       branchname, &fDigits, buffersize, file);
+  
+  if (fDigits   && fLoader->TreeD() && d) {
+    MakeBranchInTree(gAlice->TreeD(), branchname, &fDigits, buffersize, 0);
   }    
 
-  if (fHitType>1) MakeBranch2(option,file); // MI change 14.09.2000
+  if (fHitType>1) MakeBranch2(option,0); // MI change 14.09.2000
 }
  
 //_____________________________________________________________________________
@@ -2457,7 +1874,6 @@ void AliTPC::SetSecAL(Int_t sec)
   //-----------------------------------------------------------------
   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
   //-----------------------------------------------------------------
-
   fSecAL = sec;
 }
 
@@ -2471,7 +1887,6 @@ void AliTPC::SetSecAU(Int_t sec)
   //-----------------------------------------------------------------
   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
   //-----------------------------------------------------------------
-
   fSecAU = sec;
 }
 
@@ -2552,17 +1967,15 @@ void AliTPC::SetGasMixt(Int_t nc,Int_t c1,Int_t c2,Int_t c3,Float_t p1,
 
   // gax mixture definition
 
- fNoComp = nc;
 fNoComp = nc;
  
- fMixtComp[0]=c1;
- fMixtComp[1]=c2;
- fMixtComp[2]=c3;
 fMixtComp[0]=c1;
 fMixtComp[1]=c2;
 fMixtComp[2]=c3;
 
- fMixtProp[0]=p1;
- fMixtProp[1]=p2;
- fMixtProp[2]=p3; 
+  fMixtProp[0]=p1;
+  fMixtProp[1]=p2;
+  fMixtProp[2]=p3; 
 }
 //_____________________________________________________________________________
 
@@ -2592,32 +2005,12 @@ void AliTPC::TransportElectron(Float_t *xyz, Int_t *index)
   // ExB
   
   if (fTPCParam->GetMWPCReadout()==kTRUE){
-    Float_t x1=xyz[0];
-    fTPCParam->Transform2to2NearestWire(xyz,index);
-    Float_t dx=xyz[0]-x1;
+    Float_t dx = fTPCParam->Transform2to2NearestWire(xyz,index);
     xyz[1]+=dx*(fTPCParam->GetOmegaTau());
   }
-  //add nonisochronity (not implemented yet)
-  
+  //add nonisochronity (not implemented yet)  
 }
   
-ClassImp(AliTPCdigit)
-//_____________________________________________________________________________
-AliTPCdigit::AliTPCdigit(Int_t *tracks, Int_t *digits):
-  AliDigit(tracks)
-{
-  //
-  // Creates a TPC digit object
-  //
-  fSector     = digits[0];
-  fPadRow     = digits[1];
-  fPad        = digits[2];
-  fTime       = digits[3];
-  fSignal     = digits[4];
-}
-
 ClassImp(AliTPChit)
  
 //_____________________________________________________________________________
@@ -2635,16 +2028,16 @@ AliHit(shunt,track)
   fQ          = hits[3];
 }
  
-
 //________________________________________________________________________
 // Additional code because of the AliTPCTrackHitsV2
 
-void AliTPC::MakeBranch2(Option_t *option,const char *file)
+void AliTPC::MakeBranch2(Option_t *option,const char */*file*/)
 {
   //
   // Create a new branch in the current Root Tree
   // The branch of fHits is automatically split
   // MI change 14.09.2000
+  AliDebug(1,"");
   if (fHitType<2) return;
   char branchname[10];
   sprintf(branchname,"%s2",GetName());  
@@ -2652,65 +2045,26 @@ void AliTPC::MakeBranch2(Option_t *option,const char *file)
   // Get the pointer to the header
   const char *cH = strstr(option,"H");
   //
-  if (fTrackHits   && gAlice->TreeH() && cH && fHitType&4) {    
-    //    AliObjectBranch * branch = new AliObjectBranch(branchname,"AliTPCTrackHitsV2",&fTrackHits, 
-    //                                            gAlice->TreeH(),fBufferSize,99);
-    //TBranch * branch = 
-    gAlice->TreeH()->Branch(branchname,"AliTPCTrackHitsV2",&fTrackHits, 
-                                                  fBufferSize,99);
-
-    // gAlice->TreeH()->GetListOfBranches()->Add(branch);
-    if (GetDebug()>1) 
-      printf("* AliDetector::MakeBranch * Making Branch %s for trackhits\n",branchname);
-    const char folder [] = "RunMC/Event/Data";
-    if (GetDebug())
-      printf("%15s: Publishing %s to %s\n",ClassName(),branchname,folder);
-    Publish(folder,&fTrackHits,branchname);
-    if (file) {
-        TBranch *b = gAlice->TreeH()->GetBranch(branchname);
-        TDirectory *wd = gDirectory;
-        b->SetFile(file);
-        TIter next( b->GetListOfBranches());
-        while ((b=(TBranch*)next())) {
-         b->SetFile(file);
-        }
-        wd->cd(); 
-        if (GetDebug()>1) 
-             cout << "Diverting branch " << branchname << " to file " << file << endl;  
-    }
-  }    
+  if (fTrackHits   && TreeH() && cH && fHitType&4) {
+    AliDebug(1,"Making branch for Type 4 Hits");
+    TreeH()->Branch(branchname,"AliTPCTrackHitsV2",&fTrackHits,fBufferSize,99);
+  }
 
-  if (fTrackHitsOld   && gAlice->TreeH() && cH && fHitType&2) {    
+  if (fTrackHitsOld   && TreeH() && cH && fHitType&2) {    
+    AliDebug(1,"Making branch for Type 2 Hits");
     AliObjectBranch * branch = new AliObjectBranch(branchname,"AliTPCTrackHits",&fTrackHitsOld, 
-                                                  gAlice->TreeH(),fBufferSize,99);
-    //TBranch * branch = gAlice->TreeH()->Branch(branchname,"AliTPCTrackHitsV2",&fTrackHits, 
-    //                                                    fBufferSize,99);
-
-    gAlice->TreeH()->GetListOfBranches()->Add(branch);
-    if (GetDebug()>1) 
-      printf("* AliDetector::MakeBranch * Making Branch %s for trackhits\n",branchname);
-    const char folder [] = "RunMC/Event/Data";
-    if (GetDebug())
-      printf("%15s: Publishing %s to %s\n",ClassName(),branchname,folder);
-    Publish(folder,&fTrackHitsOld,branchname);
-    if (file) {
-        TBranch *b = gAlice->TreeH()->GetBranch(branchname);
-        TDirectory *wd = gDirectory;
-        b->SetFile(file);
-        TIter next( b->GetListOfBranches());
-        while ((b=(TBranch*)next())) {
-         b->SetFile(file);
-        }
-        wd->cd(); 
-        if (GetDebug()>1) 
-             cout << "Diverting branch " << branchname << " to file " << file << endl;  
-    }
+                                                   TreeH(),fBufferSize,99);
+    TreeH()->GetListOfBranches()->Add(branch);
   }    
 }
 
 void AliTPC::SetTreeAddress()
 {
-  if (fHitType<=1) AliDetector::SetTreeAddress();
+  //Sets tree address for hits  
+  if (fHitType<=1) {
+    if (fHits == 0x0 ) fHits = new TClonesArray("AliTPChit", 176);//skowron 20.06.03
+    AliDetector::SetTreeAddress();
+  }
   if (fHitType>1) SetTreeAddress2();
 }
 
@@ -2719,19 +2073,39 @@ void AliTPC::SetTreeAddress2()
   //
   // Set branch address for the TrackHits Tree
   // 
+  AliDebug(1,"");
+  
   TBranch *branch;
   char branchname[20];
   sprintf(branchname,"%s2",GetName());
   //
   // Branch address for hit tree
-  TTree *treeH = gAlice->TreeH();
+  TTree *treeH = TreeH();
   if ((treeH)&&(fHitType&4)) {
     branch = treeH->GetBranch(branchname);
-    if (branch) branch->SetAddress(&fTrackHits);
+    if (branch) {
+      branch->SetAddress(&fTrackHits);
+      AliDebug(1,"fHitType&4 Setting");
+    }
+    else 
+      AliDebug(1,"fHitType&4 Failed (can not find branch)");
+    
   }
   if ((treeH)&&(fHitType&2)) {
     branch = treeH->GetBranch(branchname);
-    if (branch) branch->SetAddress(&fTrackHitsOld);
+    if (branch) {
+      branch->SetAddress(&fTrackHitsOld);
+      AliDebug(1,"fHitType&2 Setting");
+    }
+    else
+      AliDebug(1,"fHitType&2 Failed (can not find branch)");
+  }
+  //set address to TREETR
+  
+  TTree *treeTR = TreeTR();
+  if (treeTR && fTrackReferences) {
+    branch = treeTR->GetBranch(GetName());
+    if (branch) branch->SetAddress(&fTrackReferences);
   }
 
 }
@@ -2749,27 +2123,24 @@ void AliTPC::AddHit2(Int_t track, Int_t *vol, Float_t *hits)
   // add hit to the list  
   Int_t rtrack;
   if (fIshunt) {
-    int primary = gAlice->GetPrimary(track);
-    gAlice->Particle(primary)->SetBit(kKeepBit);
+    int primary = gAlice->GetMCApp()->GetPrimary(track);
+    gAlice->GetMCApp()->Particle(primary)->SetBit(kKeepBit);
     rtrack=primary;
   } else {
     rtrack=track;
-    gAlice->FlagTrack(track);
+    gAlice->GetMCApp()->FlagTrack(track);
   }  
-  //AliTPChit *hit = (AliTPChit*)fHits->UncheckedAt(fNhits-1);
-  //if (hit->fTrack!=rtrack)
-  //  cout<<"bad track number\n";
   if (fTrackHits && fHitType&4) 
     fTrackHits->AddHitKartez(vol[0],rtrack, hits[0],
-                            hits[1],hits[2],(Int_t)hits[3]);
+                             hits[1],hits[2],(Int_t)hits[3]);
   if (fTrackHitsOld &&fHitType&2 ) 
     fTrackHitsOld->AddHitKartez(vol[0],rtrack, hits[0],
-                            hits[1],hits[2],(Int_t)hits[3]);
+                                hits[1],hits[2],(Int_t)hits[3]);
   
 }
 
 void AliTPC::ResetHits()
-{
+{ 
   if (fHitType&1) AliDetector::ResetHits();
   if (fHitType>1) ResetHits2();
 }
@@ -2790,6 +2161,9 @@ AliHit* AliTPC::FirstHit(Int_t track)
 }
 AliHit* AliTPC::NextHit()
 {
+  //
+  // gets next hit
+  //
   if (fHitType>1) return NextHit2();
   
   return AliDetector::NextHit();
@@ -2806,7 +2180,7 @@ AliHit* AliTPC::FirstHit2(Int_t track)
   // 
   if(track>=0) {
     gAlice->ResetHits();
-    gAlice->TreeH()->GetEvent(track);
+    TreeH()->GetEvent(track);
   }
   //
   if (fTrackHits && fHitType&4) {
@@ -2843,19 +2217,17 @@ void AliTPC::LoadPoints(Int_t)
 {
   //
   Int_t a = 0;
-  /*  if(fHitType==1) return AliDetector::LoadPoints(a);
-  LoadPoints2(a);
-  */
+
   if(fHitType==1) AliDetector::LoadPoints(a);
   else LoadPoints2(a);
-   
-  // LoadPoints3(a);
-
 }
 
 
 void AliTPC::RemapTrackHitIDs(Int_t *map)
 {
+  //
+  // remapping
+  //
   if (!fTrackHits) return;
   
   if (fTrackHitsOld && fHitType&2){
@@ -2882,7 +2254,7 @@ Bool_t   AliTPC::TrackInVolume(Int_t id,Int_t track)
   //
   //  return kTRUE;
   if (fTrackHitsOld && fHitType&2) {
-    TBranch * br = gAlice->TreeH()->GetBranch("fTrackHitsInfo");
+    TBranch * br = TreeH()->GetBranch("fTrackHitsInfo");
     br->GetEvent(track);
     AliObjectArray * ar = fTrackHitsOld->fTrackHitsInfo;
     for (UInt_t j=0;j<ar->GetSize();j++){
@@ -2891,14 +2263,14 @@ Bool_t   AliTPC::TrackInVolume(Int_t id,Int_t track)
   }
 
   if (fTrackHits && fHitType&4) {
-    TBranch * br1 = gAlice->TreeH()->GetBranch("fVolumes");
-    TBranch * br2 = gAlice->TreeH()->GetBranch("fNVolumes");    
+    TBranch * br1 = TreeH()->GetBranch("fVolumes");
+    TBranch * br2 = TreeH()->GetBranch("fNVolumes");    
     br2->GetEvent(track);
     br1->GetEvent(track);    
     Int_t *volumes = fTrackHits->GetVolumes();
     Int_t nvolumes = fTrackHits->GetNVolumes();
     if (!volumes && nvolumes>0) {
-      printf("Problematic track\t%d\t%d",track,nvolumes);
+      AliWarning(Form("Problematic track\t%d\t%d",track,nvolumes));
       return kFALSE;
     }
     for (Int_t j=0;j<nvolumes; j++)
@@ -2906,7 +2278,7 @@ Bool_t   AliTPC::TrackInVolume(Int_t id,Int_t track)
   }
 
   if (fHitType&1) {
-    TBranch * br = gAlice->TreeH()->GetBranch("fSector");
+    TBranch * br = TreeH()->GetBranch("fSector");
     br->GetEvent(track);
     for (Int_t j=0;j<fHits->GetEntriesFast();j++){
       if (  ((AliTPChit*)fHits->At(j))->fSector==id) return kTRUE;
@@ -2929,7 +2301,7 @@ void AliTPC::LoadPoints2(Int_t)
   if (fHitType&2) nhits = fTrackHitsOld->GetEntriesFast();
   
   if (nhits == 0) return;
-  Int_t tracks = gAlice->GetNtrack();
+  Int_t tracks = gAlice->GetMCApp()->GetNtrack();
   if (fPoints == 0) fPoints = new TObjArray(tracks);
   AliHit *ahit;
   //
@@ -3004,7 +2376,7 @@ void AliTPC::LoadPoints3(Int_t)
   //
   Int_t nhits = fTrackHits->GetEntriesFast();
   if (nhits == 0) return;
-  Int_t tracks = gAlice->GetNtrack();
+  Int_t tracks = gAlice->GetMCApp()->GetNtrack();
   if (fPoints == 0) fPoints = new TObjArray(2*tracks);
   fPoints->Expand(2*tracks);
   AliHit *ahit;
@@ -3026,11 +2398,9 @@ void AliTPC::LoadPoints3(Int_t)
   // Loop over all the hits and store their position
   //
   ahit = FirstHit2(-1);
-  //for (Int_t hit=0;hit<nhits;hit++) {
 
   Int_t lastrow = -1;
   while (ahit){
-    //    ahit = (AliHit*)fHits->UncheckedAt(hit);
     trk=ahit->GetTrack(); 
     Float_t  x[3]={ahit->X(),ahit->Y(),ahit->Z()};
     Int_t    index[3]={1,((AliTPChit*)ahit)->fSector,0};
@@ -3068,7 +2438,6 @@ void AliTPC::LoadPoints3(Int_t)
       points->SetMarkerSize(0.2);
       points->SetDetector(this);
       points->SetParticle(trk);
-      //      points->SetPolyMarker(ntrk[trk],coor[trk],GetMarkerStyle()20);
       points->SetPolyMarker(ntrk[trk],coor[trk],30);
       fPoints->AddAt(points,tracks+trk);
       delete [] coor[trk];
@@ -3082,173 +2451,56 @@ void AliTPC::LoadPoints3(Int_t)
 
 
 
-void AliTPC::FindTrackHitsIntersection(TClonesArray * arr)
+AliLoader* AliTPC::MakeLoader(const char* topfoldername)
 {
-
-  //
-  //fill clones array with intersection of current point with the
-  //middle of the row
-  Int_t sector;
-  Int_t ipart;
-  
-  const Int_t kcmaxhits=30000;
-  AliTPCFastVector * xxxx = new AliTPCFastVector(kcmaxhits*4);
-  AliTPCFastVector & xxx = *xxxx;
-  Int_t maxhits = kcmaxhits;
-      
-  //
-  AliTPChit * tpcHit=0;
-  tpcHit = (AliTPChit*)FirstHit2(-1);
-  Int_t currentIndex=0;
-  Int_t lastrow=-1;  //last writen row
-
-  while (tpcHit){
-    if (tpcHit==0) continue;
-    sector=tpcHit->fSector; // sector number
-    ipart=tpcHit->Track();
-    
-    //find row number
-    
-    Float_t  x[3]={tpcHit->X(),tpcHit->Y(),tpcHit->Z()};
-    Int_t    index[3]={1,sector,0};
-    Int_t    currentrow = fTPCParam->GetPadRow(x,index) ;      
-    if (currentrow<0) continue;
-    if (lastrow<0) lastrow=currentrow;
-    if (currentrow==lastrow){
-      if ( currentIndex>=maxhits){
-       maxhits+=kcmaxhits;
-       xxx.ResizeTo(4*maxhits);
-      }     
-      xxx(currentIndex*4)=x[0];
-      xxx(currentIndex*4+1)=x[1];
-      xxx(currentIndex*4+2)=x[2];      
-      xxx(currentIndex*4+3)=tpcHit->fQ;
-      currentIndex++;  
-    }
-    else 
-      if (currentIndex>2){
-       Float_t sumx=0;
-       Float_t sumx2=0;
-       Float_t sumx3=0;
-       Float_t sumx4=0;
-       Float_t sumy=0;
-       Float_t sumxy=0;
-       Float_t sumx2y=0;
-       Float_t sumz=0;
-       Float_t sumxz=0;
-       Float_t sumx2z=0;
-       Float_t sumq=0;
-       for (Int_t index=0;index<currentIndex;index++){
-         Float_t x,x2,x3,x4;
-         x=x2=x3=x4=xxx(index*4);
-         x2*=x;
-         x3*=x2;
-         x4*=x3;
-         sumx+=x;
-         sumx2+=x2;
-         sumx3+=x3;
-         sumx4+=x4;
-         sumy+=xxx(index*4+1);
-         sumxy+=xxx(index*4+1)*x;
-         sumx2y+=xxx(index*4+1)*x2;
-         sumz+=xxx(index*4+2);
-         sumxz+=xxx(index*4+2)*x;
-         sumx2z+=xxx(index*4+2)*x2;     
-         sumq+=xxx(index*4+3);
-       }
-       Float_t centralPad = (fTPCParam->GetNPads(sector,lastrow)-1)/2;
-       Float_t det=currentIndex*(sumx2*sumx4-sumx3*sumx3)-sumx*(sumx*sumx4-sumx2*sumx3)+
-         sumx2*(sumx*sumx3-sumx2*sumx2);
-       
-       Float_t detay=sumy*(sumx2*sumx4-sumx3*sumx3)-sumx*(sumxy*sumx4-sumx2y*sumx3)+
-         sumx2*(sumxy*sumx3-sumx2y*sumx2);
-       Float_t detaz=sumz*(sumx2*sumx4-sumx3*sumx3)-sumx*(sumxz*sumx4-sumx2z*sumx3)+
-         sumx2*(sumxz*sumx3-sumx2z*sumx2);
-       
-       Float_t detby=currentIndex*(sumxy*sumx4-sumx2y*sumx3)-sumy*(sumx*sumx4-sumx2*sumx3)+
-         sumx2*(sumx*sumx2y-sumx2*sumxy);
-       Float_t detbz=currentIndex*(sumxz*sumx4-sumx2z*sumx3)-sumz*(sumx*sumx4-sumx2*sumx3)+
-         sumx2*(sumx*sumx2z-sumx2*sumxz);
-       
-       Float_t y=detay/det+centralPad;
-       Float_t z=detaz/det;    
-       Float_t by=detby/det; //y angle
-       Float_t bz=detbz/det; //z angle
-       sumy/=Float_t(currentIndex);
-       sumz/=Float_t(currentIndex);
-       
-       AliComplexCluster cl;
-       cl.fX=z;
-       cl.fY=y;
-       cl.fQ=sumq;
-       cl.fSigmaX2=bz;
-       cl.fSigmaY2=by;
-       cl.fTracks[0]=ipart;
-       
-       AliTPCClustersRow * row = (fClustersArray->GetRow(sector,lastrow));
-       if (row!=0) row->InsertCluster(&cl);
-       currentIndex=0;
-       lastrow=currentrow;
-      } //end of calculating cluster for given row
-               
-  } // end of loop over hits
-  xxxx->Delete();
-
+  //Makes TPC loader
+  fLoader = new AliTPCLoader(GetName(),topfoldername);
+  return fLoader;
 }
-//_______________________________________________________________________________
-void AliTPC::Digits2Reco(Int_t firstevent,Int_t lastevent)
-{
-  // produces rec points from digits and writes them on the root file
-  // AliTPCclusters.root
-
-  TDirectory *cwd = gDirectory;
-
 
-  AliTPCParamSR *dig=(AliTPCParamSR *)gDirectory->Get("75x40_100x60");
-  if(dig){
-    printf("You are running 2 pad-length geom hits with 3 pad-length geom digits\n");
-    delete dig;
-    dig = new AliTPCParamSR();
-  }
-  else
-  {
-   dig=(AliTPCParamSR *)gDirectory->Get("75x40_100x60_150x60"); 
-  }
-  if(!dig){
-   printf("No TPC parameters found\n");
-   exit(3);
-  }
-   
-  SetParam(dig);
-  cout<<"AliTPC::Digits2Reco: TPC parameteres have been set"<<endl; 
-  TFile *out;
-  if(!gSystem->Getenv("CONFIG_FILE")){
-    out=TFile::Open("AliTPCclusters.root","recreate");
-  }
-  else {
-    const char *tmp1;
-    const char *tmp2;
-    char tmp3[80];
-    tmp1=gSystem->Getenv("CONFIG_FILE_PREFIX");
-    tmp2=gSystem->Getenv("CONFIG_OUTDIR");
-    sprintf(tmp3,"%s%s/AliTPCclusters.root",tmp1,tmp2);
-    out=TFile::Open(tmp3,"recreate");
+////////////////////////////////////////////////////////////////////////
+AliTPCParam* AliTPC::LoadTPCParam(TFile *file) {
+//
+// load TPC paarmeters from a given file or create new if the object
+// is not found there
+// 12/05/2003 This method should be moved to the AliTPCLoader
+// and one has to decide where to store the TPC parameters
+// M.Kowalski
+  char paramName[50];
+  sprintf(paramName,"75x40_100x60_150x60");
+  AliTPCParam *paramTPC=(AliTPCParam*)file->Get(paramName);
+  if (paramTPC) {
+    AliDebugClass(1,Form("TPC parameters %s found.",paramName));
+  } else {
+    AliWarningClass("TPC parameters not found. Create new (they may be incorrect)");
+    paramTPC = new AliTPCParamSR;
   }
+  return paramTPC;
 
-  TStopwatch timer;
-  cout<<"AliTPC::Digits2Reco - determination of rec points begins"<<endl;
-  timer.Start();
-  cwd->cd();
-  for(Int_t iev=firstevent;iev<lastevent+1;iev++){
+// the older version of parameters can be accessed with this code.
+// In some cases, we have old parameters saved in the file but 
+// digits were created with new parameters, it can be distinguish 
+// by the name of TPC TreeD. The code here is just for the case 
+// we would need to compare with old data, uncomment it if needed.
+//
+//  char paramName[50];
+//  sprintf(paramName,"75x40_100x60");
+//  AliTPCParam *paramTPC=(AliTPCParam*)in->Get(paramName);
+//  if (paramTPC) {
+//    cout<<"TPC parameters "<<paramName<<" found."<<endl;
+//  } else {
+//    sprintf(paramName,"75x40_100x60_150x60");
+//    paramTPC=(AliTPCParam*)in->Get(paramName);
+//    if (paramTPC) {
+//     cout<<"TPC parameters "<<paramName<<" found."<<endl;
+//    } else {
+//     cerr<<"TPC parameters not found. Create new (they may be incorrect)."
+//         <<endl;    
+//     paramTPC = new AliTPCParamSR;
+//    }
+//  }
+//  return paramTPC;
 
-    printf("Processing event %d\n",iev);
-    Digits2Clusters(out,iev);
-  }
-  cout<<"AliTPC::Digits2Reco - determination of rec points ended"<<endl;
-  timer.Stop();
-  timer.Print();
-  out->Close();
-  cwd->cd(); 
+}
 
 
-}