]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPC.cxx
Corrected raw->sdigits method. Using default altro raw-data reader class. Avoid alloc...
[u/mrichter/AliRoot.git] / TPC / AliTPC.cxx
index 0a419aeb6af2f42e68f9ca3e3ea9310b9a8183d7..af8bc32128b5d6871e480467dcd9767601026f5c 100644 (file)
 #include <Riostream.h>
 #include <stdlib.h>
 
+#include <TF2.h>
 #include <TFile.h>  
-#include <TGeometry.h>
+#include <TGeoGlobalMagField.h>
 #include <TInterpreter.h>
 #include <TMath.h>
 #include <TMatrixF.h>
-#include <TVector.h>
-#include <TNode.h>
 #include <TObjectTable.h>
 #include <TParticle.h>
 #include <TROOT.h>
 #include <TRandom.h>
+#include <TStopwatch.h>
+#include <TString.h>
 #include <TSystem.h>     
-#include <TTUBS.h>
 #include <TTree.h>
+#include <TVector.h>
 #include <TVirtualMC.h>
-#include <TString.h>
-#include <TF2.h>
-#include <TStopwatch.h>
+#include <TParameter.h>
 
 #include "AliDigits.h"
 #include "AliMagF.h"
-#include "AliPoints.h"
 #include "AliRun.h"
 #include "AliRunLoader.h"
 #include "AliSimDigits.h"
 #include "AliTPCPRF2D.h"
 #include "AliTPCParamSR.h"
 #include "AliTPCRF1D.h"
-//#include "AliTPCTrackHits.h"
 #include "AliTPCTrackHitsV2.h"
 #include "AliTrackReference.h"
 #include "AliMC.h"
+#include "AliStack.h"
 #include "AliTPCDigitizer.h"
 #include "AliTPCBuffer.h"
 #include "AliTPCDDLRawData.h"
 #include "AliLog.h"
-
+#include "AliTPCcalibDB.h"
+#include "AliTPCCalPad.h"
+#include "AliTPCCalROC.h"
+#include "AliTPCExB.h"
+#include "AliRawReader.h"
+#include "AliTPCRawStreamV3.h"
+#include "TTreeStream.h"
 
 ClassImp(AliTPC) 
 //_____________________________________________________________________________
-AliTPC::AliTPC()
+  AliTPC::AliTPC():AliDetector(),
+                  fDefaults(0),
+                  fSens(0),
+                  fNsectors(0),
+                  fDigitsArray(0),
+                  fTPCParam(0),
+                  fTrackHits(0),
+                  fHitType(0),
+                  fDigitsSwitch(0),
+                  fSide(0),
+                   fPrimaryIonisation(0),
+                  fNoiseDepth(0),
+                  fNoiseTable(0),
+                  fCurrentNoise(0),
+                  fActiveSectors(0),
+                   fGainFactor(1.),
+                   fDebugStreamer(0),
+                   fLHCclockPhaseSw(0)
+
 {
   //
   // Default constructor
   //
   fIshunt   = 0;
-  fHits     = 0;
-  fDigits   = 0;
-  fNsectors = 0;
-  fDigitsArray = 0;
-  fDefaults = 0;
-  fTrackHits = 0; 
   //  fTrackHitsOld = 0;   
 #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;
-  fSens = 0;
-
 }
  
 //_____________________________________________________________________________
 AliTPC::AliTPC(const char *name, const char *title)
-      : AliDetector(name,title)
+  : AliDetector(name,title),
+                   fDefaults(0),
+                  fSens(0),
+                  fNsectors(0),
+                  fDigitsArray(0),
+                  fTPCParam(0),
+                  fTrackHits(0),
+                  fHitType(0),
+                  fDigitsSwitch(0),
+                  fSide(0),
+                   fPrimaryIonisation(0),
+                  fNoiseDepth(0),
+                  fNoiseTable(0),
+                  fCurrentNoise(0),
+                   fActiveSectors(0),
+                   fGainFactor(1.),
+                   fDebugStreamer(0),
+                   fLHCclockPhaseSw(0)
+                  
 {
   //
   // Standard constructor
@@ -117,8 +147,6 @@ AliTPC::AliTPC(const char *name, const char *title)
   // Initialise arrays of hits and digits 
   fHits     = new TClonesArray("AliTPChit",  176);
   gAlice->GetMCApp()->AddHitList(fHits); 
-  fDigitsArray = 0;
-  fDefaults = 0;
   //
   fTrackHits = new AliTPCTrackHitsV2;  
   fTrackHits->SetHitPrecision(0.002);
@@ -130,23 +158,20 @@ AliTPC::AliTPC(const char *name, const char *title)
   //fTrackHitsOld->SetStepPrecision(0.003);  
   //fTrackHitsOld->SetMaxDistance(100); 
 
-  fNoiseTable =0;
 
 #if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
   fHitType = 4; // ROOT containers
 #else
   fHitType = 2;
 #endif
-  fActiveSectors = 0;
-  //
-  // Initialise counters
-  fNsectors = 0;
+
+
 
   //
   fIshunt     =  0;
   //
   // Initialise color attributes
-  SetMarkerColor(kYellow);
+  //PH SetMarkerColor(kYellow);
 
   //
   //  Set TPC parameters
@@ -154,21 +179,21 @@ AliTPC::AliTPC(const char *name, const char *title)
 
 
   if (!strcmp(title,"Default")) {       
-    fTPCParam = new AliTPCParamSR;
+    //fTPCParam = new AliTPCParamSR;
+    fTPCParam = AliTPCcalibDB::Instance()->GetParameters();
   } else {
     AliWarning("In Config.C you must set non-default parameters.");
     fTPCParam=0;
   }
-  fSens = 0;
 
 }
-
-//_____________________________________________________________________________
-AliTPC::AliTPC(const AliTPC& t):AliDetector(t){
-  //
-  // dummy copy constructor
+void AliTPC::CreateDebugStremer(){
   //
+  // Create Debug streamer to check simulation
+  // 
+  fDebugStreamer = new TTreeSRedirector("TPCSimdebug.root");
 }
+//_____________________________________________________________________________
 AliTPC::~AliTPC()
 {
   //
@@ -178,11 +203,14 @@ AliTPC::~AliTPC()
   fIshunt   = 0;
   delete fHits;
   delete fDigits;
-  delete fTPCParam;
+  //delete fTPCParam;
   delete fTrackHits; //MI 15.09.2000
   //  delete fTrackHitsOld; //MI 10.12.2001
-  if (fNoiseTable) delete [] fNoiseTable;
-
+  
+  fDigitsArray = 0x0;
+  delete [] fNoiseTable;
+  delete [] fActiveSectors;
+  if (fDebugStreamer) delete fDebugStreamer;
 }
 
 //_____________________________________________________________________________
@@ -199,91 +227,6 @@ void AliTPC::AddHit(Int_t track, Int_t *vol, Float_t *hits)
     AddHit2(track,vol,hits);
 }
 
-//_____________________________________________________________________________
-void AliTPC::BuildGeometry()
-{
-
-  //
-  // Build TPC ROOT TNode geometry for the event display
-  //
-  TNode *nNode, *nTop;
-  TTUBS *tubs;
-  Int_t i;
-  const int kColorTPC=19;
-  char name[5], title[25];
-  const Double_t kDegrad=TMath::Pi()/180;
-  const Double_t kRaddeg=180./TMath::Pi();
-
-
-  Float_t innerOpenAngle = fTPCParam->GetInnerAngle();
-  Float_t outerOpenAngle = fTPCParam->GetOuterAngle();
-
-  Float_t innerAngleShift = fTPCParam->GetInnerAngleShift();
-  Float_t outerAngleShift = fTPCParam->GetOuterAngleShift();
-
-  Int_t nLo = fTPCParam->GetNInnerSector()/2;
-  Int_t nHi = fTPCParam->GetNOuterSector()/2;  
-
-  const Double_t kloAng = (Double_t)TMath::Nint(innerOpenAngle*kRaddeg);
-  const Double_t khiAng = (Double_t)TMath::Nint(outerOpenAngle*kRaddeg);
-  const Double_t kloAngSh = (Double_t)TMath::Nint(innerAngleShift*kRaddeg);
-  const Double_t khiAngSh = (Double_t)TMath::Nint(outerAngleShift*kRaddeg);  
-
-
-  const Double_t kloCorr = 1/TMath::Cos(0.5*kloAng*kDegrad);
-  const Double_t khiCorr = 1/TMath::Cos(0.5*khiAng*kDegrad);
-
-  Double_t rl,ru;
-  
-
-  //
-  // Get ALICE top node
-  //
-
-  nTop=gAlice->GetGeometry()->GetNode("alice");
-
-  //  inner sectors
-
-  rl = fTPCParam->GetInnerRadiusLow();
-  ru = fTPCParam->GetInnerRadiusUp();
-
-  for(i=0;i<nLo;i++) {
-    sprintf(name,"LS%2.2d",i);
-    name[4]='\0';
-    sprintf(title,"TPC low sector %3d",i);
-    title[24]='\0';
-    
-    tubs = new TTUBS(name,title,"void",rl*kloCorr,ru*kloCorr,250.,
-                     kloAng*(i-0.5)+kloAngSh,kloAng*(i+0.5)+kloAngSh);
-    tubs->SetNumberOfDivisions(1);
-    nTop->cd();
-    nNode = new TNode(name,title,name,0,0,0,"");
-    nNode->SetLineColor(kColorTPC);
-    fNodes->Add(nNode);
-  }
-
-  // Outer sectors
-
-  rl = fTPCParam->GetOuterRadiusLow();
-  ru = fTPCParam->GetOuterRadiusUp();
-
-  for(i=0;i<nHi;i++) {
-    sprintf(name,"US%2.2d",i);
-    name[4]='\0';
-    sprintf(title,"TPC upper sector %d",i);
-    title[24]='\0';
-    tubs = new TTUBS(name,title,"void",rl*khiCorr,ru*khiCorr,250,
-                     khiAng*(i-0.5)+khiAngSh,khiAng*(i+0.5)+khiAngSh);
-    tubs->SetNumberOfDivisions(1);
-    nTop->cd();
-    nNode = new TNode(name,title,name,0,0,0,"");
-    nNode->SetLineColor(kColorTPC);
-    fNodes->Add(nNode);
-  }
-
-}    
-
 //_____________________________________________________________________________
 void AliTPC::CreateMaterials()
 {
@@ -295,8 +238,8 @@ void AliTPC::CreateMaterials()
   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
   //-----------------------------------------------------------------
 
-   Int_t iSXFLD=gAlice->Field()->Integ();
-  Float_t sXMGMX=gAlice->Field()->Max();
+   Int_t iSXFLD=((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->Integ();
+  Float_t sXMGMX=((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->Max();
 
   Float_t amat[5]; // atomic numbers
   Float_t zmat[5]; // z
@@ -324,7 +267,7 @@ void AliTPC::CreateMaterials()
   wmat[0]=0.2729;
   wmat[1]=0.7271;
 
-  density=0.001977;
+  density=0.001754609;
 
 
   AliMixture(10,"CO2",amat,zmat,density,2,wmat);
@@ -348,9 +291,11 @@ void AliTPC::CreateMaterials()
   // drift gases 
   //----------------------------------------------------------------
 
-
+  //
   // Drift gases 1 - nonsensitive, 2 - sensitive
-  // Ne-CO2-N (85-10-5)
+  // Ne-CO2-N2 (90-10-5) (volume) values at 20deg and 1 atm.
+  // rho(Ne) = 0.839 g/cm^3, rho(CO2) = 1.842 g/cm^3, rho(N2) = 1.165 g/cm^3
+  // for the calculation - everything is normalized to 1
 
   amat[0]= 20.18;
   amat[1]=12.011;
@@ -362,16 +307,16 @@ void AliTPC::CreateMaterials()
   zmat[2]=8.;
   zmat[3]=7.;
 
-  wmat[0]=0.7707;
-  wmat[1]=0.0539;
-  wmat[2]=0.1438;
-  wmat[3]=0.0316;
+  wmat[0]=0.756992632;
+  wmat[1]=0.056235789;
+  wmat[2]=0.128469474;
+  wmat[3]=0.058395789;
  
-  density=0.0010252;
+  density=0.000904929;
 
   AliMixture(12,"Ne-CO2-N-1",amat,zmat,density,4,wmat);
   AliMixture(13,"Ne-CO2-N-2",amat,zmat,density,4,wmat);
-
+  AliMixture(30,"Ne-CO2-N-3",amat,zmat,density,4,wmat);
   //----------------------------------------------------------------------
   //               solid materials
   //----------------------------------------------------------------------
@@ -508,7 +453,7 @@ void AliTPC::CreateMaterials()
   wmat[0]=0.194;
   wmat[1]=0.023;
   wmat[2]=0.443;
-  wmat[3]=0.340;
+  wmat[3]=0.34;
 
   density=1.82;
 
@@ -526,12 +471,12 @@ void AliTPC::CreateMaterials()
   zmat[2]=8.;
   zmat[3]=14.;
 
-  wmat[0]=0.225;
+  wmat[0]=0.257;
   wmat[1]=0.03;
-  wmat[2]=0.443;
+  wmat[2]=0.412;
   wmat[3]=0.3;
 
-  density=1.163;
+  density=1.725;
 
   AliMixture(21, "Prepreg3",amat,zmat,density,4,wmat);
 
@@ -583,6 +528,24 @@ void AliTPC::CreateMaterials()
 
   AliMaterial(25,"Cu",amat[0],zmat[0],density,999.,999.);
 
+  // brass
+
+  amat[0] = 63.546;
+  zmat[0] = 29.;
+  //
+  amat[1]= 65.409;
+  zmat[1]= 30.;
+  //
+  wmat[0]= 0.6;
+  wmat[1]= 0.4;
+
+  //
+  density = 8.23;
+  
+  //
+  AliMixture(33,"Brass",amat,zmat,density,2,wmat);
+  
   // Epoxy - C14 H20 O3
  
   amat[0]=12.011;
@@ -600,6 +563,28 @@ void AliTPC::CreateMaterials()
   density=1.25;
 
   AliMixture(26,"Epoxy",amat,zmat,density,-3,wmat);
+  //
+  // epoxy film - 90% epoxy, 10% glass fiber 
+  //
+  amat[0]=12.01;
+  amat[1]=1.;
+  amat[2]=15.994;
+  amat[3]=28.086;
+
+  zmat[0]=6.;
+  zmat[1]=1.;
+  zmat[2]=8.;
+  zmat[3]=14.;
+
+  wmat[0]=0.596;
+  wmat[1]=0.071;
+  wmat[2]=0.257;
+  wmat[3]=0.076;
+
+
+  density=1.345;
+
+  AliMixture(34, "Epoxy-film",amat,zmat,density,4,wmat);
 
   // Plexiglas  C5H8O2
 
@@ -636,6 +621,58 @@ void AliTPC::CreateMaterials()
   density=7.87;
 
   AliMaterial(29,"Fe",amat[0],zmat[0],density,999.,999.);
+  //
+  // Peek - (C6H4-O-OC6H4-O-C6H4-CO)n
+  amat[0]=12.011;
+  amat[1]=1.;
+  amat[2]=15.9994;
+
+  zmat[0]=6.;
+  zmat[1]=1.;
+  zmat[2]=8.;
+
+  wmat[0]=19.;
+  wmat[1]=12.;
+  wmat[2]=3.;
+  //
+  density=1.3;
+  //
+  AliMixture(30,"Peek",amat,zmat,density,-3,wmat);  
+  //
+  //  Ceramics - Al2O3
+  //
+  amat[0] = 26.98;
+  amat[1]= 15.9994;
+
+  zmat[0] = 13.;
+  zmat[1]=8.;
+  wmat[0]=2.;
+  wmat[1]=3.;
+  density = 3.97;
+
+  AliMixture(31,"Alumina",amat,zmat,density,-2,wmat);   
+
+  //
+  // liquids
+  //
+
+  // water
+
+  amat[0]=1.;
+  amat[1]=15.9994;
+
+  zmat[0]=1.;
+  zmat[1]=8.;
+
+  wmat[0]=2.;
+  wmat[1]=1.;
+
+  density=1.;
+
+  AliMixture(32,"Water",amat,zmat,density,-2,wmat);  
+
  
   //----------------------------------------------------------
   // tracking media for gases
@@ -645,7 +682,7 @@ void AliTPC::CreateMaterials()
   AliMedium(1, "Ne-CO2-N-1", 12, 0, iSXFLD, sXMGMX, 10., 999.,.1,.001, .001);
   AliMedium(2, "Ne-CO2-N-2", 13, 1, iSXFLD, sXMGMX, 10., 999.,.1,.001, .001);
   AliMedium(3,"CO2",10,0, iSXFLD, sXMGMX, 10., 999.,.1, .001, .001); 
-
+  AliMedium(20, "Ne-CO2-N-3", 30, 1, iSXFLD, sXMGMX, 10., 999.,.1,.001, .001);
   //-----------------------------------------------------------  
   // tracking media for solids
   //-----------------------------------------------------------
@@ -667,7 +704,11 @@ void AliTPC::CreateMaterials()
   AliMedium(16,"G10",22,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
   AliMedium(17,"Plexiglas",27,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
   AliMedium(18,"Steel",29,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001); 
-    
+  AliMedium(19,"Peek",30,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
+  AliMedium(21,"Alumina",31,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);    
+  AliMedium(22,"Water",32,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
+  AliMedium(23,"Brass",33,0, iSXFLD, sXMGMX, 10., 999., .1, .001, .001);
+  AliMedium(24,"Epoxyfm",34,0, iSXFLD, sXMGMX, 10., 999., .1, .0005, .001);  
 }
 
 void AliTPC::GenerNoise(Int_t tablesize)
@@ -712,8 +753,7 @@ 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()];
+  if (!fActiveSectors) fActiveSectors = new Bool_t[fTPCParam->GetNSector()];
   for (Int_t i=0;i<fTPCParam->GetNSector();i++) fActiveSectors[i]=kFALSE;
   for (Int_t i=0;i<n;i++) 
     if ((sectors[i]>=0) && sectors[i]<fTPCParam->GetNSector())  fActiveSectors[sectors[i]]=kTRUE;
@@ -727,41 +767,49 @@ void    AliTPC::SetActiveSectors(Int_t flag)
   //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()];
+  if (!fActiveSectors) fActiveSectors = new Bool_t[fTPCParam->GetNSector()];
   if (flag) {
     for (Int_t i=0;i<fTPCParam->GetNSector();i++) fActiveSectors[i]=kTRUE;
     return;
   }
   for (Int_t i=0;i<fTPCParam->GetNSector();i++) fActiveSectors[i]=kFALSE;
   TBranch * branch=0;
-  if (TreeH() == 0x0)
+  if (fLoader->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();
+  if (fHitType>1) branch = fLoader->TreeH()->GetBranch("TPC2");
+  else branch = fLoader->TreeH()->GetBranch("TPC");
+  Stat_t ntracks = fLoader->TreeH()->GetEntries();
   // loop over all hits
-  AliDebug(1,Form("Got %d tracks",ntracks));
+  AliDebug(1,Form("Got %d tracks", (Int_t) ntracks));
   
   for(Int_t track=0;track<ntracks;track++) {
     ResetHits();
     //
     if (fTrackHits && fHitType&4) {
-      TBranch * br1 = TreeH()->GetBranch("fVolumes");
-      TBranch * br2 = TreeH()->GetBranch("fNVolumes");
+      TBranch * br1 = fLoader->TreeH()->GetBranch("fVolumes");
+      TBranch * br2 = fLoader->TreeH()->GetBranch("fNVolumes");
       br1->GetEvent(track);
       br2->GetEvent(track);
       Int_t *volumes = fTrackHits->GetVolumes();
-      for (Int_t j=0;j<fTrackHits->GetNVolumes(); j++)
-       fActiveSectors[volumes[j]]=kTRUE;
+      for (Int_t j=0;j<fTrackHits->GetNVolumes(); j++) {
+       if (volumes[j]>-1 && volumes[j]<fTPCParam->GetNSector()) {
+         fActiveSectors[volumes[j]]=kTRUE;
+       }
+       else {
+           AliError(Form("Volume %d -> sector number %d is outside (0..%d)",
+                         j,
+                         volumes[j],
+                         fTPCParam->GetNSector()));
+       }
+      }
     }
     
     //
 //     if (fTrackHitsOld && fHitType&2) {
-//       TBranch * br = TreeH()->GetBranch("fTrackHitsInfo");
+//       TBranch * br = fLoader->TreeH()->GetBranch("fTrackHitsInfo");
 //       br->GetEvent(track);
 //       AliObjectArray * ar = fTrackHitsOld->fTrackHitsInfo;
 //       for (UInt_t j=0;j<ar->GetSize();j++){
@@ -780,7 +828,6 @@ void AliTPC::Digits2Raw()
 // convert digits of the current event to raw data
 
   static const Int_t kThreshold = 0;
-  static const Bool_t kCompress = kTRUE;
 
   fLoader->LoadDigits();
   TTree* digits = fLoader->TreeD();
@@ -866,14 +913,162 @@ void AliTPC::Digits2Raw()
   rawWriter.RawData(fileName);
   gSystem->Unlink(fileName);
 
-  if (kCompress) {
-    AliInfo("Compressing raw data");
-    rawWriter.RawDataCompDecompress(kTRUE);
-    gSystem->Unlink("Statistics");
-  }
 }
 
 
+//_____________________________________________________________________________
+Bool_t AliTPC::Raw2SDigits(AliRawReader* rawReader){
+  // Converts the TPC raw data into summable digits
+  // The method is used for merging simulated and
+  // real data events
+  if (fLoader->TreeS() == 0x0 ) {
+    fLoader->MakeTree("S");
+  }
+
+  if(fDefaults == 0) SetDefaults();  // check if the parameters are set
+
+  //setup TPCDigitsArray 
+  if(GetDigitsArray()) delete GetDigitsArray();
+
+  AliTPCDigitsArray *arr = new AliTPCDigitsArray; 
+  arr->SetClass("AliSimDigits");
+  arr->Setup(fTPCParam);
+  arr->MakeTree(fLoader->TreeS());
+
+  SetDigitsArray(arr);
+
+  // set zero suppression to "0"
+  fTPCParam->SetZeroSup(0);
+
+  // Loop over sectors
+  const Int_t kmaxTime = fTPCParam->GetMaxTBin();
+  const Int_t kNIS = fTPCParam->GetNInnerSector();
+  const Int_t kNOS = fTPCParam->GetNOuterSector();
+  const Int_t kNS = kNIS + kNOS;
+
+  // Setup storage
+  AliTPCROC * roc = AliTPCROC::Instance();
+  Int_t nRowsMax = roc->GetNRows(roc->GetNSector()-1);
+  Int_t nPadsMax = roc->GetNPads(roc->GetNSector()-1,nRowsMax-1);
+  Short_t** allBins = new Short_t*[nRowsMax];
+  for (Int_t iRow = 0; iRow < nRowsMax; iRow++) {
+    Int_t maxBin = kmaxTime*nPadsMax;
+    allBins[iRow] = new Short_t[maxBin];
+    memset(allBins[iRow],0,sizeof(Short_t)*maxBin);
+  }
+
+  for(Int_t iSector = 0; iSector < kNS; iSector++) {
+    
+    Int_t nRows = fTPCParam->GetNRow(iSector);
+    Int_t nDDLs = 0, indexDDL = 0;
+    if (iSector < kNIS) {
+      nDDLs = 2;
+      indexDDL = iSector * 2;
+    }
+    else {
+      nDDLs = 4;
+      indexDDL = (iSector-kNIS) * 4 + kNIS * 2;
+    }
+
+    // Load the raw data for corresponding DDLs
+    rawReader->Reset();
+
+    AliTPCAltroMapping** mapping =AliTPCcalibDB::Instance()->GetMapping();
+    AliTPCRawStreamV3 input(rawReader,(AliAltroMapping**)mapping);
+    rawReader->Select("TPC",indexDDL,indexDDL+nDDLs-1);
+
+    // Clean storage
+    for (Int_t iRow = 0; iRow < nRowsMax; iRow++) {
+      Int_t maxBin = kmaxTime*nPadsMax;
+      memset(allBins[iRow],0,sizeof(Short_t)*maxBin);
+    }
+
+    // Begin loop over altro data
+    while (input.NextDDL()) {
+
+      if (input.GetSector() != iSector)
+       AliFatal(Form("Sector index mismatch ! Expected (%d), but got (%d) !",iSector,input.GetSector()));
+
+      //loop over pads
+      while ( input.NextChannel() ) {
+
+        Int_t iRow = input.GetRow();
+        if (iRow < 0 || iRow >= nRows)
+          AliFatal(Form("Pad-row index (%d) outside the range (%d -> %d) !",
+                        iRow, 0, nRows -1));
+        Int_t iPad = input.GetPad();
+
+        Int_t maxPad = fTPCParam->GetNPads(iSector,iRow);
+
+        if (iPad < 0 || iPad >= maxPad)
+          AliFatal(Form("Pad index (%d) outside the range (%d -> %d) !",
+                        iPad, 0, maxPad -1));
+
+        //loop over bunches
+        while ( input.NextBunch() ){
+          Int_t  startTbin    = (Int_t)input.GetStartTimeBin();
+          Int_t  bunchlength  = (Int_t)input.GetBunchLength();
+          const UShort_t *sig = input.GetSignals();
+          for (Int_t iTime = 0; iTime<bunchlength; iTime++){
+            Int_t iTimeBin=startTbin-iTime;
+            if ( iTimeBin < 0 || iTimeBin >= kmaxTime) {
+              continue;
+              //AliFatal(Form("Timebin index (%d) outside the range (%d -> %d) !",
+              //               iTimeBin, 0, kmaxTime -1));
+            }
+
+            Int_t maxBin = kmaxTime*maxPad;
+            if (((iPad*kmaxTime+iTimeBin) >= maxBin) ||
+                ((iPad*kmaxTime+iTimeBin) < 0))
+              AliFatal(Form("Index outside the allowed range"
+                            " Sector=%d Row=%d Pad=%d Timebin=%d"
+                            " (Max.index=%d)",iSector,iRow,iPad,iTimeBin,maxBin));
+            allBins[iRow][iPad*kmaxTime+iTimeBin] = sig[iTime];
+          }
+        }
+      } // End loop over altro data
+    }
+
+    // Now fill the digits array
+    if (fDigitsArray->GetTree()==0) {
+      AliFatal("Tree not set in fDigitsArray");
+    }
+
+    for (Int_t iRow = 0; iRow < nRows; iRow++) {
+      AliDigits * dig = fDigitsArray->CreateRow(iSector,iRow);
+      Int_t maxPad = fTPCParam->GetNPads(iSector,iRow);
+      for(Int_t iPad = 0; iPad < maxPad; iPad++) {
+       for(Int_t iTimeBin = 0; iTimeBin < kmaxTime; iTimeBin++) {
+         Short_t q = allBins[iRow][iPad*kmaxTime + iTimeBin];
+         if (q <= 0) continue;
+         q *= 16;
+         dig->SetDigitFast((Short_t)q,iTimeBin,iPad);
+       }
+      }
+      fDigitsArray->StoreRow(iSector,iRow);
+      Int_t ndig = dig->GetDigitSize(); 
+       
+      AliDebug(10,
+              Form("*** Sector, row, compressed digits %d %d %d ***\n",
+                   iSector,iRow,ndig));        
+       
+      fDigitsArray->ClearRow(iSector,iRow);  
+
+    } // end of the sector digitization
+  }
+
+  fLoader->WriteSDigits("OVERWRITE");
+
+  if(GetDigitsArray()) delete GetDigitsArray();
+  SetDigitsArray(0x0);
+
+  // cleanup storage
+  for (Int_t iRow = 0; iRow < nRowsMax; iRow++)
+    delete [] allBins[iRow];
+  delete [] allBins;
+
+  return kTRUE;
+}
 
 //______________________________________________________________________
 AliDigitizer* AliTPC::CreateDigitizer(AliRunDigitizer* manager) const
@@ -954,7 +1149,7 @@ void AliTPC::SDigits2Digits2(Int_t /*eventnumber*/)
     Short_t * pamp1 = digrow->GetDigits();
     Int_t   * ptracks1 = digrow->GetTracks();
     Int_t  nelems =nrows*ncols;
-    Int_t saturation = fTPCParam->GetADCSat();
+    Int_t saturation = fTPCParam->GetADCSat() - 1;
     //use internal structure of the AliDigits - for speed reason
     //if you cahnge implementation
     //of the Alidigits - it must be rewriten -
@@ -996,14 +1191,21 @@ void AliTPC::SetDefaults(){
   AliRunLoader* rl = (AliRunLoader*)fLoader->GetEventFolder()->FindObject(AliRunLoader::GetRunLoaderName());
   rl->CdGAFile();
   AliTPCParamSR *param=(AliTPCParamSR*)gDirectory->Get("75x40_100x60");
-  if(param){
-    AliInfo("You are using 2 pad-length geom hits with 3 pad-lenght geom digits...");
-    delete param;
-    param = new AliTPCParamSR();
-  }
-  else {
-    param=(AliTPCParamSR*)gDirectory->Get("75x40_100x60_150x60");
-  }
+  // if(param){
+//     AliInfo("You are using 2 pad-length geom hits with 3 pad-lenght geom digits...");
+//     delete param;
+//     param = new AliTPCParamSR();
+//   }
+//   else {
+//     param=(AliTPCParamSR*)gDirectory->Get("75x40_100x60_150x60");
+//   }
+  param = (AliTPCParamSR*)AliTPCcalibDB::Instance()->GetParameters();
+  if (!param->IsGeoRead()){
+      //
+      // read transformation matrices for gGeoManager
+      //
+      param->ReadGeoMatrices();
+    }
   if(!param){
     AliFatal("No TPC parameters found");
   }
@@ -1012,11 +1214,24 @@ void AliTPC::SetDefaults(){
   AliTPCPRF2D    * prfinner   = new AliTPCPRF2D;
   AliTPCPRF2D    * prfouter1   = new AliTPCPRF2D;
   AliTPCPRF2D    * prfouter2   = new AliTPCPRF2D;  
-  AliTPCRF1D     * rf    = new AliTPCRF1D(kTRUE);
-  rf->SetGauss(param->GetZSigma(),param->GetZWidth(),1.);
-  rf->SetOffset(3*param->GetZSigma());
-  rf->Update();
+
   
+  //AliTPCRF1D     * rf    = new AliTPCRF1D(kTRUE);
+  //rf->SetGauss(param->GetZSigma(),param->GetZWidth(),1.);
+  //rf->SetOffset(3*param->GetZSigma());
+  //rf->Update();
+  //
+  // Use gamma 4
+  //
+  char  strgamma4[1000];
+  sprintf(strgamma4,"AliTPCRF1D::Gamma4((x-0.135+%f)*%f,55,160)",3*param->GetZSigma(), 1000000000*param->GetTSample()/param->GetZWidth());
+  
+  TF1 * fgamma4 = new TF1("fgamma4",strgamma4, -1,1);
+  AliTPCRF1D     * rf    = new AliTPCRF1D(kTRUE,1000);
+  rf->SetParam(fgamma4,param->GetZWidth(), 1,0.2);
+  rf->SetOffset(3*param->GetZSigma()); 
+  rf->Update();
+
   TDirectory *savedir=gDirectory;
   TFile *f=TFile::Open("$ALICE_ROOT/TPC/AliTPCprf2d.root");
   if (!f->IsOpen()) 
@@ -1061,14 +1276,19 @@ void AliTPC::Hits2Digits()
   //
   // creates digits from hits
   //
+  if (!fTPCParam->IsGeoRead()){
+    //
+    // read transformation matrices for gGeoManager
+    //
+    fTPCParam->ReadGeoMatrices();
+  }
 
   fLoader->LoadHits("read");
   fLoader->LoadDigits("recreate");
   AliRunLoader* runLoader = fLoader->GetRunLoader(); 
 
   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
-    runLoader->GetEvent(iEvent);
-    SetActiveSectors();   
+    //PH    runLoader->GetEvent(iEvent);
     Hits2Digits(iEvent);
   }
 
@@ -1083,6 +1303,7 @@ void AliTPC::Hits2Digits(Int_t eventnumber)
  //----------------------------------------------------
   AliRunLoader* rl = (AliRunLoader*)fLoader->GetEventFolder()->FindObject(AliRunLoader::GetRunLoaderName());
   rl->GetEvent(eventnumber);
+  SetActiveSectors();   
   if (fLoader->TreeH() == 0x0) {
     if(fLoader->LoadHits()) {
       AliError("Can not load hits.");
@@ -1113,14 +1334,32 @@ void AliTPC::Hits2Digits(Int_t eventnumber)
   SetDigitsArray(arr);
 
   fDigitsSwitch=0; // standard digits
-
+  // here LHC clock phase
+  Float_t lhcph = 0.;
+  switch (fLHCclockPhaseSw){
+  case 0: 
+    // no phase
+    lhcph=0.;
+    break;
+  case 1:
+    // random phase
+    lhcph = (Int_t)(gRandom->Rndm()/0.25);    
+    break;
+  case 2:
+    lhcph=0.;
+    // not implemented yet
+    break;
+  }
+  // adding phase to the TreeD user info 
+  fLoader->TreeD()->GetUserInfo()->Add(new TParameter<float>("lhcphase0",lhcph));
+  //
   for(Int_t isec=0;isec<fTPCParam->GetNSector();isec++) 
     if (IsSectorActive(isec)) {
-      AliDebug(1,Form("Hits2Digits","Sector %d is active.",isec));
+      AliDebug(1,Form("Hits2DigitsSector %d is active.",isec));
       Hits2DigitsSector(isec);
     }
     else {
-      AliDebug(1,Form("Hits2Digits","Sector %d is NOT active.",isec));
+      AliDebug(1,Form("Hits2DigitsSector %d is NOT active.",isec));
     }
   
   fLoader->WriteDigits("OVERWRITE"); 
@@ -1180,6 +1419,25 @@ void AliTPC::Hits2SDigits2(Int_t eventnumber)
   fDigitsSwitch=1; // summable digits
   
     // set zero suppression to "0"
+  // here LHC clock phase
+  Float_t lhcph = 0.;
+  switch (fLHCclockPhaseSw){
+  case 0: 
+    // no phase
+    lhcph=0.;
+    break;
+  case 1:
+    // random phase
+    lhcph = (Int_t)(gRandom->Rndm()/0.25);    
+    break;
+  case 2:
+    lhcph=0.;
+    // not implemented yet
+    break;
+  }
+  // adding phase to the TreeS user info 
+  
+  fLoader->TreeS()->GetUserInfo()->Add(new TParameter<float>("lhcphase0",lhcph));
 
   fTPCParam->SetZeroSup(0);
 
@@ -1206,6 +1464,13 @@ void AliTPC::Hits2SDigits()
   //   summable digits - 16 bit "ADC", no noise, no saturation
   //-----------------------------------------------------------
 
+  if (!fTPCParam->IsGeoRead()){
+    //
+    // read transformation matrices for gGeoManager
+    //
+    fTPCParam->ReadGeoMatrices();
+  }
+  
   fLoader->LoadHits("read");
   fLoader->LoadSDigits("recreate");
   AliRunLoader* runLoader = fLoader->GetRunLoader(); 
@@ -1216,9 +1481,13 @@ void AliTPC::Hits2SDigits()
     SetActiveSectors();
     Hits2SDigits2(iEvent);
   }
-
+  
   fLoader->UnloadHits();
   fLoader->UnloadSDigits();
+  if (fDebugStreamer) {
+    delete fDebugStreamer;
+    fDebugStreamer=0;
+  }    
 }
 //_____________________________________________________________________________
 
@@ -1241,7 +1510,7 @@ void AliTPC::Hits2DigitsSector(Int_t isec)
 
   if(fDefaults == 0) SetDefaults();
 
-  TTree *tH = TreeH(); // pointer to the hits tree
+  TTree *tH = fLoader->TreeH(); // pointer to the hits tree
   if (tH == 0x0) {
     AliFatal("Can not find TreeH in folder");
     return;
@@ -1249,11 +1518,7 @@ void AliTPC::Hits2DigitsSector(Int_t isec)
 
   Stat_t ntracks = tH->GetEntries();
 
-  if( ntracks > 0){
 
-  //------------------------------------------- 
-  //  Only if there are any tracks...
-  //-------------------------------------------
 
     TObjArray **row;
     
@@ -1302,7 +1567,6 @@ void AliTPC::Hits2DigitsSector(Int_t isec)
       
     delete [] row; // delete the array of pointers to TObjArray-s
         
-  } // ntracks >0
 
 } // end of Hits2DigitsSector
 
@@ -1321,6 +1585,11 @@ void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
   //-----------------------------------------------------------------
  
   Float_t zerosup = fTPCParam->GetZeroSup();
+  AliTPCCalPad * gainTPC = AliTPCcalibDB::Instance()->GetDedxGainFactor(); 
+  AliTPCCalPad * noiseTPC = AliTPCcalibDB::Instance()->GetPadNoise(); 
+  AliTPCCalROC * gainROC = gainTPC->GetCalROC(isec);  // pad gains per given sector
+  AliTPCCalROC * noiseROC = noiseTPC->GetCalROC(isec);  // noise per given sector
+
 
   fCurrentIndex[1]= isec;
   
@@ -1374,11 +1643,15 @@ void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
     for(Int_t ip=0;ip<nofPads;ip++){
       gi++;
       Float_t q=total(ip,it);      
-      if(fDigitsSwitch == 0){
-       q+=GetNoise();
+      if(fDigitsSwitch == 0){  
+       Float_t gain = gainROC->GetValue(irow,ip);  // get gain for given - pad-row pad 
+       Float_t noisePad = noiseROC->GetValue(irow,ip); 
+       //
+       q*=gain;
+       q+=GetNoise()*noisePad;
         if(q <=fzerosup) continue; // do not fill zeros
         q = TMath::Nint(q);
-        if(q > fTPCParam->GetADCSat()) q = fTPCParam->GetADCSat();  // saturation
+        if(q >= fTPCParam->GetADCSat()) q = fTPCParam->GetADCSat() - 1;  // saturation
 
       }
 
@@ -1412,7 +1685,14 @@ void AliTPC::DigitizeRow(Int_t irow,Int_t isec,TObjArray **rows)
     
     } // end of loop over time buckets
   }  // end of lop over pads 
+  //
+  // test
+  //
+  //
 
+  // glitch filters if normal simulated digits
+  //
+  if(!fDigitsSwitch) ((AliSimDigits*)dig)->GlitchFilter();
   //
   //  This row has been digitized, delete nonused stuff
   //
@@ -1451,7 +1731,7 @@ Float_t AliTPC::GetSignal(TObjArray *p1, Int_t ntr,
   TVector &v = *tv;
   
   Float_t label = v(0);
-  Int_t centralPad = (fTPCParam->GetNPads(fCurrentIndex[1],fCurrentIndex[3]-1)-1)/2;
+  Int_t centralPad = (fTPCParam->GetNPads(fCurrentIndex[1],fCurrentIndex[3]-1))/2;
 
   Int_t nElectrons = (tv->GetNrows()-1)/5;
   indexRange[0]=9999; // min pad
@@ -1462,6 +1742,15 @@ Float_t AliTPC::GetSignal(TObjArray *p1, Int_t ntr,
   TMatrixF &signal = *m1;
   TMatrixF &total = *m2;
   //
+  // Get LHC clock phase
+  //
+  TParameter<float> *ph;
+  if(fDigitsSwitch){// s-digits
+    ph = (TParameter<float>*)fLoader->TreeS()->GetUserInfo()->FindObject("lhcphase0");  
+  }
+  else{ // normal digits
+    ph = (TParameter<float>*)fLoader->TreeD()->GetUserInfo()->FindObject("lhcphase0");
+  } 
   //  Loop over all electrons
   //
   for(Int_t nel=0; nel<nElectrons; nel++){
@@ -1469,7 +1758,8 @@ Float_t AliTPC::GetSignal(TObjArray *p1, Int_t ntr,
     Float_t aval =  v(idx+4);
     Float_t eltoadcfac=aval*fTPCParam->GetTotalNormFac(); 
     Float_t xyz[4]={v(idx+1),v(idx+2),v(idx+3),v(idx+5)};
-    Int_t n = ((AliTPCParamSR*)fTPCParam)->CalcResponseFast(xyz,fCurrentIndex,fCurrentIndex[3]);
+    Int_t n = ((AliTPCParamSR*)fTPCParam)->CalcResponseFast(xyz,fCurrentIndex,
+                                                           fCurrentIndex[3],ph->GetVal());
 
     Int_t *index = fTPCParam->GetResBin(0);  
     Float_t *weight = & (fTPCParam->GetResWeight(0));
@@ -1597,11 +1887,46 @@ void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
   // contains the track label and the position of electrons.
   //-----------------------------------------------------------------
 
+  // 
+  // The trasport of the electrons through TPC drift volume
+  //    Drift (drift velocity + velocity map(not yet implemented)))
+  //    Application of the random processes (diffusion, gas gain)
+  //    Systematic effects (ExB effect in drift volume + ROCs)  
+  //
+  // Algorithm:
+  // Loop over primary electrons:
+  //    Creation of the secondary electrons
+  //    Loop over electrons (primary+ secondaries)
+  //        Global coordinate frame:
+  //          1. Skip electrons if attached  
+  //          2. ExB effect in drift volume
+  //             a.) Simulation   calib->GetExB()->CorrectInverse(dxyz0,dxyz1);
+  //             b.) Reconstruction -  calib->GetExB()->CorrectInverse(dxyz0,dxyz1);
+  //          3. Generation of gas gain (Random - Exponential distribution) 
+  //          4. TransportElectron function (diffusion)
+  //
+  //        5. Conversion to the local coordinate frame  pad-row, pad, timebin
+  //        6. Apply Time0 shift - AliTPCCalPad class 
+  //            a.) Plus sign in simulation
+  //            b.) Minus sign in reconstruction 
+  // end of loop          
+  //
   //-----------------------------------------------------------------
   // Origin: Marek Kowalski  IFJ, Krakow, Marek.Kowalski@ifj.edu.pl
+  // Origin: Marian Ivanov,  marian.ivanov@cern.ch
   //-----------------------------------------------------------------
+  AliTPCcalibDB* const calib=AliTPCcalibDB::Instance();
+  if (gAlice){ // Set correctly the magnetic field in the ExB calculation
+    if (!calib->GetExB()){
+      AliMagF * field = ((AliMagF*)TGeoGlobalMagField::Instance()->GetField()); 
+      if (field) {
+       calib->SetExBField(field);
+      }
+    }
+  }
 
   Float_t gasgain = fTPCParam->GetGasGain();
+  gasgain = gasgain/fGainFactor;
   Int_t i;
   Float_t xyz[5]; 
 
@@ -1661,6 +1986,13 @@ void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
        continue; 
       }
 
+      // Remove hits which arrive before the TPC opening gate signal
+      if(((fTPCParam->GetZLength(isec)-TMath::Abs(tpcHit->Z()))
+         /fTPCParam->GetDriftV()+tpcHit->Time())<fTPCParam->GetGateDelay()) {
+       tpcHit = (AliTPChit*) NextHit();
+       continue;
+      }
+
       currentTrack = tpcHit->Track(); // track number
 
       if(currentTrack != previousTrack){
@@ -1696,7 +2028,7 @@ void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
       //---------------------------------------------------
 
 
-      Float_t time = 1.e6*(fTPCParam->GetZLength()-TMath::Abs(tpcHit->Z()))
+      Float_t time = 1.e6*(fTPCParam->GetZLength(isec)-TMath::Abs(tpcHit->Z()))
        /fTPCParam->GetDriftV(); 
       // in microseconds!      
       Float_t attProb = fTPCParam->GetAttCoef()*
@@ -1710,28 +2042,81 @@ void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
       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();     
+       
+       //
+       // ExB effect
+       //
+       Double_t dxyz0[3],dxyz1[3];
+       dxyz0[0]=tpcHit->X();
+       dxyz0[1]=tpcHit->Y();
+       dxyz0[2]=tpcHit->Z();   
+       if (calib->GetExB()){
+         calib->GetExB()->CorrectInverse(dxyz0,dxyz1);
+       }else{
+         AliError("Not valid ExB calibration");
+         dxyz1[0]=tpcHit->X();
+         dxyz1[1]=tpcHit->Y();
+         dxyz1[2]=tpcHit->Z();         
+       }
+       xyz[0]=dxyz1[0];
+       xyz[1]=dxyz1[1];
+       xyz[2]=dxyz1[2];        
+       //
+       //
        //
        // 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); 
+       Int_t padrow = fTPCParam->GetPadRow(xyz,index); 
+       //
+       // Add Time0 correction due unisochronity
+       // xyz[0] - pad row coordinate 
+       // xyz[1] - pad coordinate
+       // xyz[2] - is in now time bin coordinate system
+       Float_t correction =0;
+       if (calib->GetPadTime0()){
+         if (!calib->GetPadTime0()->GetCalROC(isec)) continue;   
+         Int_t npads = fTPCParam->GetNPads(isec,padrow);
+         //      Int_t pad  = TMath::Nint(xyz[1]+fTPCParam->GetNPads(isec,TMath::Nint(xyz[0]))*0.5);
+         // pad numbering from -npads/2 .. npads/2-1
+         Int_t pad  = TMath::Nint(xyz[1]+npads/2);
+         if (pad<0) pad=0;
+         if (pad>=npads) pad=npads-1;
+         correction = calib->GetPadTime0()->GetCalROC(isec)->GetValue(padrow,pad);
+         //      printf("%d\t%d\t%d\t%f\n",isec,padrow,pad,correction);
+         if (fDebugStreamer){
+           (*fDebugStreamer)<<"Time0"<<
+             "isec="<<isec<<
+             "padrow="<<padrow<<
+             "pad="<<pad<<
+             "x0="<<xyz[0]<<
+             "x1="<<xyz[1]<<
+             "x2="<<xyz[2]<<
+             "hit.="<<tpcHit<<
+             "cor="<<correction<<
+             "\n";
+         }
+       }
+       xyz[2]+=correction;
+       xyz[2]+=fTPCParam->GetNTBinsL1();    // adding Level 1 time bin offset
+       //
        // Electron track time (for pileup simulation)
-       xyz[4] = tpcHit->Time()/fTPCParam->GetTSample();
+       xyz[2]+=tpcHit->Time()/fTPCParam->GetTSample(); // adding time of flight
+       xyz[4] =0;
+
+       //
        // 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;
+       /*      Float_t x1,y1;
        if (isec <fTPCParam->GetNInnerSector()) {
          x1 = xyz[1]*fTPCParam->GetInnerPadPitchWidth();
          y1 = fTPCParam->GetYInner(rowNumber);
@@ -1743,7 +2128,7 @@ void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
        // 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.)); 
+       if(x1>y1) xyz[3]*=TMath::Max(1.e-6,(y1-x1+1.)); */
        
        nofElectrons[rowNumber]++;        
        //----------------------------------
@@ -1801,32 +2186,6 @@ void AliTPC::Init()
   AliDebug(1,"*********************************************");
 }
 
-//_____________________________________________________________________________
-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());
-  
-  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   && fLoader->TreeD() && d) {
-    MakeBranchInTree(gAlice->TreeD(), branchname, &fDigits, buffersize, 0);
-  }    
-
-  if (fHitType>1) MakeBranch2(option,0); // MI change 14.09.2000
-}
 //_____________________________________________________________________________
 void AliTPC::ResetDigits()
 {
@@ -1876,7 +2235,7 @@ void AliTPC::TransportElectron(Float_t *xyz, Int_t *index)
   // xyz and index must be already transformed to system 1
   //
 
-  fTPCParam->Transform1to2(xyz,index);
+  fTPCParam->Transform1to2(xyz,index);  // mis-alignment applied in this step
   
   //add diffusion
   Float_t driftl=xyz[2];
@@ -1894,14 +2253,32 @@ void AliTPC::TransportElectron(Float_t *xyz, Int_t *index)
     Float_t dx = fTPCParam->Transform2to2NearestWire(xyz,index);
     xyz[1]+=dx*(fTPCParam->GetOmegaTau());
   }
-  //add nonisochronity (not implemented yet)  
+  //add nonisochronity (not implemented yet) 
+  
 }
   
 ClassImp(AliTPChit)
+  //______________________________________________________________________
+  AliTPChit::AliTPChit()
+            :AliHit(),
+            fSector(0),
+            fPadRow(0),
+            fQ(0),
+            fTime(0)
+{
+  //
+  // default
+  //
+
+}
 //_____________________________________________________________________________
-AliTPChit::AliTPChit(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits):
-AliHit(shunt,track)
+AliTPChit::AliTPChit(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits)
+          :AliHit(shunt,track),
+            fSector(0),
+            fPadRow(0),
+            fQ(0),
+            fTime(0)
 {
   //
   // Creates a TPC hit object
@@ -1918,7 +2295,7 @@ AliHit(shunt,track)
 //________________________________________________________________________
 // Additional code because of the AliTPCTrackHitsV2
 
-void AliTPC::MakeBranch2(Option_t *option,const char */*file*/)
+void AliTPC::MakeBranch(Option_t *option)
 {
   //
   // Create a new branch in the current Root Tree
@@ -1932,16 +2309,16 @@ void AliTPC::MakeBranch2(Option_t *option,const char */*file*/)
   // Get the pointer to the header
   const char *cH = strstr(option,"H");
   //
-  if (fTrackHits   && TreeH() && cH && fHitType&4) {
+  if (fTrackHits   && fLoader->TreeH() && cH && fHitType&4) {
     AliDebug(1,"Making branch for Type 4 Hits");
-    TreeH()->Branch(branchname,"AliTPCTrackHitsV2",&fTrackHits,fBufferSize,99);
+    fLoader->TreeH()->Branch(branchname,"AliTPCTrackHitsV2",&fTrackHits,fBufferSize,99);
   }
 
-//   if (fTrackHitsOld   && TreeH() && cH && fHitType&2) {    
+//   if (fTrackHitsOld   && fLoader->TreeH() && cH && fHitType&2) {    
 //     AliDebug(1,"Making branch for Type 2 Hits");
 //     AliObjectBranch * branch = new AliObjectBranch(branchname,"AliTPCTrackHits",&fTrackHitsOld, 
-//                                                    TreeH(),fBufferSize,99);
-//     TreeH()->GetListOfBranches()->Add(branch);
+//                                                    fLoader->TreeH(),fBufferSize,99);
+//     fLoader->TreeH()->GetListOfBranches()->Add(branch);
 //   } 
 }
 
@@ -1967,7 +2344,7 @@ void AliTPC::SetTreeAddress2()
   sprintf(branchname,"%s2",GetName());
   //
   // Branch address for hit tree
-  TTree *treeH = TreeH();
+  TTree *treeH = fLoader->TreeH();
   if ((treeH)&&(fHitType&4)) {
     branch = treeH->GetBranch(branchname);
     if (branch) {
@@ -1987,14 +2364,6 @@ void AliTPC::SetTreeAddress2()
 //     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);
-  }
-
 }
 
 void AliTPC::FinishPrimary()
@@ -2007,7 +2376,8 @@ void AliTPC::FinishPrimary()
 void AliTPC::AddHit2(Int_t track, Int_t *vol, Float_t *hits)
 { 
   //
-  // add hit to the list  
+  // add hit to the list
+
   Int_t rtrack;
   if (fIshunt) {
     int primary = gAlice->GetMCApp()->GetPrimary(track);
@@ -2066,8 +2436,8 @@ AliHit* AliTPC::FirstHit2(Int_t track)
   // track is returned
   // 
   if(track>=0) {
-    gAlice->ResetHits();
-    TreeH()->GetEvent(track);
+    gAlice->GetMCApp()->ResetHits();
+    fLoader->TreeH()->GetEvent(track);
   }
   //
   if (fTrackHits && fHitType&4) {
@@ -2100,16 +2470,6 @@ AliHit* AliTPC::NextHit2()
     return 0;
 }
 
-void AliTPC::LoadPoints(Int_t)
-{
-  //
-  Int_t a = 0;
-
-  if(fHitType==1) AliDetector::LoadPoints(a);
-  else LoadPoints2(a);
-}
-
-
 void AliTPC::RemapTrackHitIDs(Int_t *map)
 {
   //
@@ -2129,7 +2489,7 @@ void AliTPC::RemapTrackHitIDs(Int_t *map)
     TClonesArray * arr = fTrackHits->GetArray();;
     for (Int_t i=0;i<arr->GetEntriesFast();i++){
       AliTrackHitsParamV2 * info = (AliTrackHitsParamV2 *)(arr->At(i));
-      info->fTrackID = map[info->fTrackID];
+      info->SetTrackID(map[info->GetTrackID()]);
     }
   }
 }
@@ -2142,7 +2502,7 @@ Bool_t   AliTPC::TrackInVolume(Int_t id,Int_t track)
   //
   //  return kTRUE;
  //  if (fTrackHitsOld && fHitType&2) {
-//     TBranch * br = TreeH()->GetBranch("fTrackHitsInfo");
+//     TBranch * br = fLoader->TreeH()->GetBranch("fTrackHitsInfo");
 //     br->GetEvent(track);
 //     AliObjectArray * ar = fTrackHitsOld->fTrackHitsInfo;
 //     for (UInt_t j=0;j<ar->GetSize();j++){
@@ -2151,8 +2511,8 @@ Bool_t   AliTPC::TrackInVolume(Int_t id,Int_t track)
 //   }
 
   if (fTrackHits && fHitType&4) {
-    TBranch * br1 = TreeH()->GetBranch("fVolumes");
-    TBranch * br2 = TreeH()->GetBranch("fNVolumes");    
+    TBranch * br1 = fLoader->TreeH()->GetBranch("fVolumes");
+    TBranch * br2 = fLoader->TreeH()->GetBranch("fNVolumes");    
     br2->GetEvent(track);
     br1->GetEvent(track);    
     Int_t *volumes = fTrackHits->GetVolumes();
@@ -2166,7 +2526,7 @@ Bool_t   AliTPC::TrackInVolume(Int_t id,Int_t track)
   }
 
   if (fHitType&1) {
-    TBranch * br = TreeH()->GetBranch("fSector");
+    TBranch * br = fLoader->TreeH()->GetBranch("fSector");
     br->GetEvent(track);
     for (Int_t j=0;j<fHits->GetEntriesFast();j++){
       if (  ((AliTPChit*)fHits->At(j))->fSector==id) return kTRUE;
@@ -2176,169 +2536,6 @@ Bool_t   AliTPC::TrackInVolume(Int_t id,Int_t track)
 
 }
 
-//_____________________________________________________________________________
-void AliTPC::LoadPoints2(Int_t)
-{
-  //
-  // Store x, y, z of all hits in memory
-  //
-  //  if (fTrackHits == 0 && fTrackHitsOld==0) return;
-  if (fTrackHits == 0 ) return;
-  //
-  Int_t nhits =0;
-  if (fHitType&4) nhits = fTrackHits->GetEntriesFast();
-  //  if (fHitType&2) nhits = fTrackHitsOld->GetEntriesFast();
-  
-  if (nhits == 0) return;
-  Int_t tracks = gAlice->GetMCApp()->GetNtrack();
-  if (fPoints == 0) fPoints = new TObjArray(tracks);
-  AliHit *ahit;
-  //
-  Int_t *ntrk=new Int_t[tracks];
-  Int_t *limi=new Int_t[tracks];
-  Float_t **coor=new Float_t*[tracks];
-  for(Int_t i=0;i<tracks;i++) {
-    ntrk[i]=0;
-    coor[i]=0;
-    limi[i]=0;
-  }
-  //
-  AliPoints *points = 0;
-  Float_t *fp=0;
-  Int_t trk;
-  Int_t chunk=nhits/4+1;
-  //
-  // Loop over all the hits and store their position
-  //
-  ahit = FirstHit2(-1);
-  while (ahit){
-    trk=ahit->GetTrack();
-    if(ntrk[trk]==limi[trk]) {
-      //
-      // Initialise a new track
-      fp=new Float_t[3*(limi[trk]+chunk)];
-      if(coor[trk]) {
-       memcpy(fp,coor[trk],sizeof(Float_t)*3*limi[trk]);
-       delete [] coor[trk];
-      }
-      limi[trk]+=chunk;
-      coor[trk] = fp;
-    } else {
-      fp = coor[trk];
-    }
-    fp[3*ntrk[trk]  ] = ahit->X();
-    fp[3*ntrk[trk]+1] = ahit->Y();
-    fp[3*ntrk[trk]+2] = ahit->Z();
-    ntrk[trk]++;
-    ahit = NextHit2();
-  }
-
-
-
-  //
-  for(trk=0; trk<tracks; ++trk) {
-    if(ntrk[trk]) {
-      points = new AliPoints();
-      points->SetMarkerColor(GetMarkerColor());
-      points->SetMarkerSize(GetMarkerSize());
-      points->SetDetector(this);
-      points->SetParticle(trk);
-      points->SetPolyMarker(ntrk[trk],coor[trk],GetMarkerStyle());
-      fPoints->AddAt(points,trk);
-      delete [] coor[trk];
-      coor[trk]=0;
-    }
-  }
-  delete [] coor;
-  delete [] ntrk;
-  delete [] limi;
-}
-
-
-//_____________________________________________________________________________
-void AliTPC::LoadPoints3(Int_t)
-{
-  //
-  // Store x, y, z of all hits in memory
-  // - only intersection point with pad row
-  if (fTrackHits == 0) return;
-  //
-  Int_t nhits = fTrackHits->GetEntriesFast();
-  if (nhits == 0) return;
-  Int_t tracks = gAlice->GetMCApp()->GetNtrack();
-  if (fPoints == 0) fPoints = new TObjArray(2*tracks);
-  fPoints->Expand(2*tracks);
-  AliHit *ahit;
-  //
-  Int_t *ntrk=new Int_t[tracks];
-  Int_t *limi=new Int_t[tracks];
-  Float_t **coor=new Float_t*[tracks];
-  for(Int_t i=0;i<tracks;i++) {
-    ntrk[i]=0;
-    coor[i]=0;
-    limi[i]=0;
-  }
-  //
-  AliPoints *points = 0;
-  Float_t *fp=0;
-  Int_t trk;
-  Int_t chunk=nhits/4+1;
-  //
-  // Loop over all the hits and store their position
-  //
-  ahit = FirstHit2(-1);
-
-  Int_t lastrow = -1;
-  while (ahit){
-    trk=ahit->GetTrack(); 
-    Float_t  x[3]={ahit->X(),ahit->Y(),ahit->Z()};
-    Int_t    index[3]={1,((AliTPChit*)ahit)->fSector,0};
-    Int_t    currentrow = fTPCParam->GetPadRow(x,index) ;
-    if (currentrow!=lastrow){
-      lastrow = currentrow;
-      //later calculate intersection point           
-      if(ntrk[trk]==limi[trk]) {
-       //
-       // Initialise a new track
-       fp=new Float_t[3*(limi[trk]+chunk)];
-       if(coor[trk]) {
-         memcpy(fp,coor[trk],sizeof(Float_t)*3*limi[trk]);
-         delete [] coor[trk];
-       }
-       limi[trk]+=chunk;
-       coor[trk] = fp;
-      } else {
-       fp = coor[trk];
-      }
-      fp[3*ntrk[trk]  ] = ahit->X();
-      fp[3*ntrk[trk]+1] = ahit->Y();
-      fp[3*ntrk[trk]+2] = ahit->Z();
-      ntrk[trk]++;
-    }
-    ahit = NextHit2();
-  }
-  
-  //
-  for(trk=0; trk<tracks; ++trk) {
-    if(ntrk[trk]) {
-      points = new AliPoints();
-      points->SetMarkerColor(GetMarkerColor()+1);
-      points->SetMarkerStyle(5);
-      points->SetMarkerSize(0.2);
-      points->SetDetector(this);
-      points->SetParticle(trk);
-      points->SetPolyMarker(ntrk[trk],coor[trk],30);
-      fPoints->AddAt(points,tracks+trk);
-      delete [] coor[trk];
-      coor[trk]=0;
-    }
-  }
-  delete [] coor;
-  delete [] ntrk;
-  delete [] limi;
-}
-
-
 
 AliLoader* AliTPC::MakeLoader(const char* topfoldername)
 {
@@ -2362,7 +2559,15 @@ AliTPCParam* AliTPC::LoadTPCParam(TFile *file) {
     AliDebugClass(1,Form("TPC parameters %s found.",paramName));
   } else {
     AliWarningClass("TPC parameters not found. Create new (they may be incorrect)");
-    paramTPC = new AliTPCParamSR;
+    //paramTPC = new AliTPCParamSR;
+    paramTPC = AliTPCcalibDB::Instance()->GetParameters();
+    if (!paramTPC->IsGeoRead()){
+      //
+      // read transformation matrices for gGeoManager
+      //
+      paramTPC->ReadGeoMatrices();
+    }
+  
   }
   return paramTPC;