X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=TPC%2FAliTPC.cxx;h=b7aa8998cdd3b9277f01948c68a08500adc14092;hb=dd4bfd89c346ca0dd412e1327c6122a7e34d53de;hp=2a71077004e76be189d8348587c14a0b94c4c0f9;hpb=9ff4af8121de7f178a8d9a97e1c7c452f5a894ed;p=u%2Fmrichter%2FAliRoot.git

diff --git a/TPC/AliTPC.cxx b/TPC/AliTPC.cxx
index 2a71077004e..b7aa8998cdd 100644
--- a/TPC/AliTPC.cxx
+++ b/TPC/AliTPC.cxx
@@ -36,28 +36,26 @@
 #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"
@@ -68,46 +66,81 @@
 #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),
+		   fIsGEM(0)
+
 {
   //
   // Default constructor
   //
   fIshunt   = 0;
-  fHits     = 0;
-  fDigits   = 0;
-  fNsectors = 0;
-  fDigitsArray = 0;
-  fDefaults = 0;
-  fTrackHits = 0; 
+  for(Int_t i=0;i<4;i++) fCurrentIndex[i]=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),
+    fIsGEM(0)
+                  
 {
   //
   // Standard constructor
@@ -117,8 +150,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 +161,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;
+
+  for(Int_t i=0;i<4;i++) fCurrentIndex[i]=0;
 
   //
   fIshunt     =  0;
   //
   // Initialise color attributes
-  SetMarkerColor(kYellow);
+  //PH SetMarkerColor(kYellow);
 
   //
   //  Set TPC parameters
@@ -154,21 +182,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,14 +206,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
   
   fDigitsArray = 0x0;
   delete [] fNoiseTable;
   delete [] fActiveSectors;
-
+  if (fDebugStreamer) delete fDebugStreamer;
 }
 
 //_____________________________________________________________________________
@@ -202,91 +230,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()
 {
@@ -298,8 +241,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
@@ -327,7 +270,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);
@@ -351,30 +294,35 @@ void AliTPC::CreateMaterials()
   // drift gases 
   //----------------------------------------------------------------
 
-
+  //
   // Drift gases 1 - nonsensitive, 2 - sensitive
-  // Ne-CO2-N (85-10-5)
+  // Ne-CO2 (90-10) (volume) values at 20deg and 1 atm.
+  // rho(Ne) = 0.839 g/cm^3, rho(CO2) = 1.842 g/cm^3
+  
 
   amat[0]= 20.18;
   amat[1]=12.011;
   amat[2]=15.9994;
-  amat[3]=14.007;
+  // amat[3]=14.007;
 
   zmat[0]= 10.; 
   zmat[1]=6.;
   zmat[2]=8.;
-  zmat[3]=7.;
-
-  wmat[0]=0.7707;
-  wmat[1]=0.0539;
-  wmat[2]=0.1438;
-  wmat[3]=0.0316;
+  // zmat[3]=7.;
+
+  //wmat[0]=0.756992632;
+  wmat[0]=0.8038965;
+  //wmat[1]=0.056235789;
+  wmat[1]= 0.053519;
+  //wmat[2]=0.128469474;
+  wmat[2]= 0.1425743;
+  // wmat[3]=0.058395789;
  
-  density=0.0010252;
-
-  AliMixture(12,"Ne-CO2-N-1",amat,zmat,density,4,wmat);
-  AliMixture(13,"Ne-CO2-N-2",amat,zmat,density,4,wmat);
+  density=0.0009393;
 
+  AliMixture(12,"Ne-CO2-1",amat,zmat,density,3,wmat);
+  AliMixture(13,"Ne-CO2-2",amat,zmat,density,3,wmat);
+  AliMixture(35,"Ne-CO2-3",amat,zmat,density,3,wmat);
   //----------------------------------------------------------------------
   //               solid materials
   //----------------------------------------------------------------------
@@ -511,7 +459,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;
 
@@ -529,12 +477,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);
 
@@ -586,6 +534,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;
@@ -603,6 +569,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
 
@@ -639,16 +627,68 @@ 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
   //----------------------------------------------------------
 
   AliMedium(0, "Air", 11, 0, iSXFLD, sXMGMX, 10., 999., .1, .01, .1);
-  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(1, "Ne-CO2-1", 12, 0, iSXFLD, sXMGMX, 10., 999.,.1,.001, .001);
+  AliMedium(2, "Ne-CO2-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-3", 35, 1, iSXFLD, sXMGMX, 10., 999.,.1,.001, .001);
   //-----------------------------------------------------------  
   // tracking media for solids
   //-----------------------------------------------------------
@@ -670,7 +710,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)
@@ -715,8 +759,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;
@@ -730,41 +773,51 @@ 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)
+  //TBranch * branch=0;
+  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");
+  if (fHitType>1) fLoader->TreeH()->GetBranch("TPC2");
+  //else branch = fLoader->TreeH()->GetBranch("TPC");
+  else 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++){
@@ -791,6 +844,7 @@ void AliTPC::Digits2Raw()
     return;
   }
 
+  //
   AliSimDigits digarr;
   AliSimDigits* digrow = &digarr;
   digits->GetBranch("Segment")->SetAddress(&digrow);
@@ -871,11 +925,182 @@ void AliTPC::Digits2Raw()
 }
 
 
+//_____________________________________________________________________________
+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);
+	  ((AliSimDigits*)dig)->SetTrackIDFast( 3141593, iTimeBin,iPad,0); 
+	}
+      }
+      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
+  }
+  // get LHC clock phase from the digits tree
+
+  TParameter<float> *ph; 
+  Float_t phase;
+  TTree *digtree = fLoader->TreeD();
+  //
+  if(digtree){ // if TreeD exists
+    ph = (TParameter<float>*)digtree->GetUserInfo()->FindObject("lhcphase0");
+    phase = ph->GetVal();
+  }
+  else{ //TreeD does not exist
+    phase = 0.; 
+  }
+    //
+    // store lhc clock phase in S-digits tree
+    //
+    fLoader->TreeS()->GetUserInfo()->Add(new TParameter<float>("lhcphase0",phase));
+   //
+   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
+AliDigitizer* AliTPC::CreateDigitizer(AliDigitizationInput* digInput) const
 {
-  return new AliTPCDigitizer(manager);
+  return new AliTPCDigitizer(digInput);
 }
 //__
 void AliTPC::SDigits2Digits2(Int_t /*eventnumber*/)  
@@ -992,51 +1217,94 @@ 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");
-  }
+  //AliTPCParamSR *param=(AliTPCParamSR*)gDirectory->Get("75x40_100x60");
+  //gDirectory->Get("75x40_100x60");
+  AliTPCParamSR *param = (AliTPCParamSR*)AliTPCcalibDB::Instance()->GetParameters();
   if(!param){
     AliFatal("No TPC parameters found");
+    return;
   }
+  if (!param->IsGeoRead()){
+      //
+      // read transformation matrices for gGeoManager
+      //
+      param->ReadGeoMatrices();
+    }
+
 
 
   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());
+  
+  snprintf(strgamma4,1000,"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()) 
-    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();
+
+  if (fIsGEM==0){
+    printf ("TPC MWPC readout\n");
+    TFile *f=TFile::Open("$ALICE_ROOT/TPC/AliTPCprf2d.root");
+    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();
+  }
+
+  if (fIsGEM==1){
+    printf ("TPC GEM readout\n");
+    TFile *f=TFile::Open("$ALICE_ROOT/TPC/AliTPCprf2dGEM.root");
+    if (!f->IsOpen()) 
+      AliFatal("Can't open $ALICE_ROOT/TPC/AliTPCprf2dGEM.root !");
+    
+    TString s;
+    prfinner->Read("prf0");
+    //PH Set different names
+    s=prfinner->GetGRF()->GetName();
+    s+="in";
+    prfinner->GetGRF()->SetName(s.Data());
+    
+    prfouter1->Read("prf1");
+    s=prfouter1->GetGRF()->GetName();
+    s+="out1";
+    prfouter1->GetGRF()->SetName(s.Data());
+    
+    prfouter2->Read("prf2");  
+    s=prfouter2->GetGRF()->GetName();
+    s+="out2";
+    prfouter2->GetGRF()->SetName(s.Data());    
+    f->Close();
+  }
   savedir->cd();
 
   param->SetInnerPRF(prfinner);
@@ -1070,8 +1338,7 @@ void AliTPC::Hits2Digits()
   AliRunLoader* runLoader = fLoader->GetRunLoader(); 
 
   for (Int_t iEvent = 0; iEvent < runLoader->GetNumberOfEvents(); iEvent++) {
-    runLoader->GetEvent(iEvent);
-    SetActiveSectors();   
+    //PH    runLoader->GetEvent(iEvent);
     Hits2Digits(iEvent);
   }
 
@@ -1086,6 +1353,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.");
@@ -1116,14 +1384,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("Hits2Digits: Sector %d is active.",isec));
       Hits2DigitsSector(isec);
     }
     else {
-      AliDebug(1,Form("Hits2Digits","Sector %d is NOT active.",isec));
+      AliDebug(1,Form("Hits2Digits: Sector %d is NOT active.",isec));
     }
   
   fLoader->WriteDigits("OVERWRITE"); 
@@ -1183,6 +1469,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);
 
@@ -1226,9 +1531,13 @@ void AliTPC::Hits2SDigits()
     SetActiveSectors();
     Hits2SDigits2(iEvent);
   }
-
+  
   fLoader->UnloadHits();
   fLoader->UnloadSDigits();
+  if (fDebugStreamer) {
+    delete fDebugStreamer;
+    fDebugStreamer=0;
+  }    
 }
 //_____________________________________________________________________________
 
@@ -1251,7 +1560,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;
@@ -1259,17 +1568,10 @@ void AliTPC::Hits2DigitsSector(Int_t isec)
 
   Stat_t ntracks = tH->GetEntries();
 
-  if( ntracks > 0){
-
-  //------------------------------------------- 
-  //  Only if there are any tracks...
-  //-------------------------------------------
-
-    TObjArray **row;
-    
     Int_t nrows =fTPCParam->GetNRow(isec);
 
-    row= new TObjArray* [nrows+2]; // 2 extra rows for cross talk
+    TObjArray **row=new TObjArray* [nrows+2]; // 2 extra rows for cross talk
+    for(Int_t j=0;j<nrows+2;j++) row[j]=0;
     
     MakeSector(isec,nrows,tH,ntracks,row);
 
@@ -1312,7 +1614,6 @@ void AliTPC::Hits2DigitsSector(Int_t isec)
       
     delete [] row; // delete the array of pointers to TObjArray-s
         
-  } // ntracks >0
 
 } // end of Hits2DigitsSector
 
@@ -1331,6 +1632,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;
   
@@ -1384,8 +1690,12 @@ 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() - 1;  // saturation
@@ -1422,7 +1732,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
   //
@@ -1472,6 +1789,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++){
@@ -1479,7 +1805,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));
@@ -1607,11 +1934,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]; 
 
@@ -1672,7 +2034,7 @@ void AliTPC::MakeSector(Int_t isec,Int_t nrows,TTree *TH,
       }
 
       // Remove hits which arrive before the TPC opening gate signal
-      if(((fTPCParam->GetZLength()-TMath::Abs(tpcHit->Z()))
+      if(((fTPCParam->GetZLength(isec)-TMath::Abs(tpcHit->Z()))
 	  /fTPCParam->GetDriftV()+tpcHit->Time())<fTPCParam->GetGateDelay()) {
 	tpcHit = (AliTPChit*) NextHit();
 	continue;
@@ -1713,7 +2075,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()*
@@ -1727,28 +2089,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);
@@ -1760,7 +2175,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]++;	  
 	//----------------------------------
@@ -1818,32 +2233,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()
 {
@@ -1893,7 +2282,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];
@@ -1911,14 +2300,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
@@ -1935,7 +2342,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
@@ -1944,21 +2351,22 @@ void AliTPC::MakeBranch2(Option_t *option,const char */*file*/)
   AliDebug(1,"");
   if (fHitType<2) return;
   char branchname[10];
-  sprintf(branchname,"%s2",GetName());  
+  //sprintf(branchname,"%s2",GetName()); 
+  snprintf(branchname,10,"%s2",GetName()); 
   //
   // 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);
 //   }	
 }
 
@@ -1981,10 +2389,11 @@ void AliTPC::SetTreeAddress2()
   
   TBranch *branch;
   char branchname[20];
-  sprintf(branchname,"%s2",GetName());
+  //sprintf(branchname,"%s2",GetName());
+  snprintf(branchname,20,"%s2",GetName());
   //
   // Branch address for hit tree
-  TTree *treeH = TreeH();
+  TTree *treeH = fLoader->TreeH();
   if ((treeH)&&(fHitType&4)) {
     branch = treeH->GetBranch(branchname);
     if (branch) {
@@ -2004,14 +2413,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()
@@ -2024,7 +2425,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);
@@ -2083,8 +2485,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) {
@@ -2117,16 +2519,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)
 {
   //
@@ -2146,7 +2538,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()]);
     }
   }
 }
@@ -2159,7 +2551,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++){
@@ -2168,8 +2560,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();
@@ -2183,7 +2575,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;
@@ -2193,169 +2585,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)
 {
@@ -2373,13 +2602,22 @@ AliTPCParam* AliTPC::LoadTPCParam(TFile *file) {
 // and one has to decide where to store the TPC parameters
 // M.Kowalski
   char paramName[50];
-  sprintf(paramName,"75x40_100x60_150x60");
+  //sprintf(paramName,"75x40_100x60_150x60");
+  snprintf(paramName,50,"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;
+    //paramTPC = new AliTPCParamSR;
+    paramTPC = AliTPCcalibDB::Instance()->GetParameters();
+    if (!paramTPC->IsGeoRead()){
+      //
+      // read transformation matrices for gGeoManager
+      //
+      paramTPC->ReadGeoMatrices();
+    }
+  
   }
   return paramTPC;