ROOT logo
/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

// _________________________________________________________________
//
// Begin_Html
//   <h2>  AliTPCSpaceCharge3D class   </h2>    
//   The class calculates the space point distortions due to an arbitrary space
//   charge distribution in 3D. 
//   <p>
//   The method of calculation is based on the analytical solution for the Poisson 
//   problem in 3D (cylindrical coordinates). The solution is used in form of 
//   look up tables, where the pre calculated solutions for different voxel 
//   positions are stored. These voxel solutions can be summed up according 
//   to the weight of the position of the applied space charge distribution.
//   Further details can be found in \cite[chap.5]{PhD-thesis_S.Rossegger}.
//   <p>
//   The class also allows a simple scaling of the resulting distortions
//   via the function SetCorrectionFactor. This speeds up the calculation 
//   time under the assumption, that the distortions scales linearly with the 
//   magnitude of the space charge distribution $\rho(r,z)$ and the shape stays 
//   the same at higher luminosities.
//   <p>
//   In contrast to the implementation in 2D (see the class AliTPCSpaceChargeabove), 
//   the input charge distribution can be of arbitrary character. An example on how 
//   to produce a corresponding charge distribution can be found in the function 
//   WriteChargeDistributionToFile. In there, a $\rho(r,z) = (A-B\,z)/r^2$, 
//   with slightly different magnitude on the A and C side (due to the muon absorber),
//   is superpositioned with a few leaking wires at arbitrary positions. 
//
//
//   Marian Ivanov change: 26.06.2013
//   Usage of the realy 3D space charge map as an optional input
//   SetInputSpaceCharge map.
//   In case given map is used 2 2D maps are ignored and  scaling functions  $\rho(r,z) = (A-B\,z)/r^2$, 
//   will not work
//


// End_Html
//
// Begin_Macro(source)
//   {
//   gROOT->SetStyle("Plain"); gStyle->SetPalette(1);
//   TCanvas *c2 = new TCanvas("cAliTPCSpaceCharge3D","cAliTPCSpaceCharge3D",500,400); 
//   AliTPCSpaceCharge3D sc;
//   sc.WriteChargeDistributionToFile("SC_zr2_GGleaks.root");
//   sc.SetSCDataFileName("SC_zr2_GGleaks.root");
//   sc.SetOmegaTauT1T2(0,1,1); // B=0
//   sc.InitSpaceCharge3DDistortion();
//   sc.CreateHistoDRinXY(15,300,300)->Draw("colz");
//   return c2;
//   } 
// End_Macro
//
// Begin_Html
//   <p>
//   Date: 19/06/2010  <br>                                                       
//   Authors: Stefan Rossegger                                                
// End_Html 
// _________________________________________________________________


#include "AliMagF.h"
#include "TGeoGlobalMagField.h"
#include "AliTPCcalibDB.h"
#include "AliTPCParam.h"
#include "AliLog.h"
#include "TH2F.h"
#include "TH3F.h"
#include "TFile.h"
#include "TVector.h"
#include "TMatrix.h"
#include "TMatrixD.h"

#include "TMath.h"
#include "AliTPCROC.h"
#include "AliTPCSpaceCharge3D.h"
#include "AliSysInfo.h"

ClassImp(AliTPCSpaceCharge3D)

AliTPCSpaceCharge3D::AliTPCSpaceCharge3D()
  : AliTPCCorrection("SpaceCharge3D","Space Charge - 3D"),
    fC0(0.),fC1(0.),
    fCorrectionFactor(1.),
    fInitLookUp(kFALSE),
    fSCDataFileName(""),
    fSCLookUpPOCsFileName3D(""),
    fSCLookUpPOCsFileNameRZ(""),
    fSCLookUpPOCsFileNameRPhi(""),
    fSCdensityInRZ(0),
    fSCdensityInRPhiA(0), 
    fSCdensityInRPhiC(0),
    fSpaceChargeHistogram3D(0),
    fSpaceChargeHistogramRPhi(0),
    fSpaceChargeHistogramRZ(0)
{
  //
  // default constructor
  //

  // Array which will contain the solution according to the setted charge density distribution
  // see InitSpaceCharge3DDistortion() function
  for ( Int_t k = 0 ; k < kNPhi ; k++ ) {
    fLookUpErOverEz[k]   =  new TMatrixF(kNR,kNZ);  
    fLookUpEphiOverEz[k] =  new TMatrixF(kNR,kNZ);
    fLookUpDeltaEz[k]    =  new TMatrixF(kNR,kNZ); 
    fSCdensityDistribution[k] = new TMatrixF(kNR,kNZ);
  }
  fSCdensityInRZ   = new TMatrixD(kNR,kNZ);
  fSCdensityInRPhiA = new TMatrixD(kNR,kNPhi);
  fSCdensityInRPhiC = new TMatrixD(kNR,kNPhi);

  // location of the precalculated look up tables

  fSCLookUpPOCsFileName3D="$(ALICE_ROOT)/TPC/Calib/maps/sc_3D_raw_18-18-26_17p-18p-25p-MN30.root"; // rough estimate
  fSCLookUpPOCsFileNameRZ="$(ALICE_ROOT)/TPC/Calib/maps/sc_radSym_35-01-51_34p-01p-50p_MN60.root";
  fSCLookUpPOCsFileNameRPhi="$(ALICE_ROOT)/TPC/Calib/maps/sc_cChInZ_35-144-26_34p-18p-01p-MN30.root";
  //  fSCLookUpPOCsFileNameRPhi="$(ALICE_ROOT)/TPC/Calib/maps/sc_cChInZ_35-36-26_34p-18p-01p-MN40.root";
 


  // standard location of the space charge distibution ... can be changes
  fSCDataFileName="$(ALICE_ROOT)/TPC/Calib/maps/sc_3D_distribution_Sim.root";

  //  SetSCDataFileName(fSCDataFileName.Data()); // should be done by the user


}

AliTPCSpaceCharge3D::~AliTPCSpaceCharge3D() {
  //
  // default destructor
  //
 
  for ( Int_t k = 0 ; k < kNPhi ; k++ ) {
    delete fLookUpErOverEz[k];
    delete fLookUpEphiOverEz[k];
    delete fLookUpDeltaEz[k];
    delete fSCdensityDistribution[k];
  }
  delete fSCdensityInRZ;
  delete fSCdensityInRPhiA;
  delete fSCdensityInRPhiC;
  delete fSpaceChargeHistogram3D;
  delete fSpaceChargeHistogramRPhi;
  delete fSpaceChargeHistogramRZ;
}


void AliTPCSpaceCharge3D::Init() {
  //
  // Initialization funtion
  //
  
  AliMagF* magF= (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
  if (!magF) AliError("Magneticd field - not initialized");
  Double_t bzField = magF->SolenoidField()/10.; //field in T
  AliTPCParam *param= AliTPCcalibDB::Instance()->GetParameters();
  if (!param) AliError("Parameters - not initialized");
  Double_t vdrift = param->GetDriftV()/1000000.; // [cm/us]   // From dataBase: to be updated: per second (ideally)
  Double_t ezField = 400; // [V/cm]   // to be updated: never (hopefully)
  Double_t wt = -10.0 * (bzField*10) * vdrift / ezField ; 
  // Correction Terms for effective omegaTau; obtained by a laser calibration run
  SetOmegaTauT1T2(wt,fT1,fT2);

  InitSpaceCharge3DDistortion(); // fill the look up table
}

void AliTPCSpaceCharge3D::Update(const TTimeStamp &/*timeStamp*/) {
  //
  // Update function 
  //
  AliMagF* magF= (AliMagF*)TGeoGlobalMagField::Instance()->GetField();
  if (!magF) AliError("Magneticd field - not initialized");
  Double_t bzField = magF->SolenoidField()/10.; //field in T
  AliTPCParam *param= AliTPCcalibDB::Instance()->GetParameters();
  if (!param) AliError("Parameters - not initialized");
  Double_t vdrift = param->GetDriftV()/1000000.; // [cm/us]  // From dataBase: to be updated: per second (ideally)
  Double_t ezField = 400; // [V/cm]   // to be updated: never (hopefully)
  Double_t wt = -10.0 * (bzField*10) * vdrift / ezField ; 
  // Correction Terms for effective omegaTau; obtained by a laser calibration run
  SetOmegaTauT1T2(wt,fT1,fT2);

  //  SetCorrectionFactor(1.); // should come from some database

}


void AliTPCSpaceCharge3D::GetCorrection(const Float_t x[],const Short_t roc,Float_t dx[]) {
  //
  // Calculates the correction due the Space Charge effect within the TPC drift volume
  //   

  if (!fInitLookUp) {
    AliInfo("Lookup table was not initialized! Performing the inizialisation now ...");
    InitSpaceCharge3DDistortion();
  }

  Int_t   order     = 1 ;    // FIXME: hardcoded? Linear interpolation = 1, Quadratic = 2         
                        
  Float_t intEr, intEphi, intdEz ;
  Double_t r, phi, z ;
  Int_t    sign;

  r      =  TMath::Sqrt( x[0]*x[0] + x[1]*x[1] ) ;
  phi    =  TMath::ATan2(x[1],x[0]) ;
  if ( phi < 0 ) phi += TMath::TwoPi() ;                   // Table uses phi from 0 to 2*Pi
  z      =  x[2] ;                                         // Create temporary copy of x[2]

  if ( (roc%36) < 18 ) {
    sign =  1;       // (TPC A side)
  } else {
    sign = -1;       // (TPC C side)
  }
  
  if ( sign==1  && z <  fgkZOffSet ) z =  fgkZOffSet;    // Protect against discontinuity at CE
  if ( sign==-1 && z > -fgkZOffSet ) z = -fgkZOffSet;    // Protect against discontinuity at CE
  

  if ( (sign==1 && z<0) || (sign==-1 && z>0) ) // just a consistency check
    AliError("ROC number does not correspond to z coordinate! Calculation of distortions is most likely wrong!");

  // Get the Er and Ephi field integrals plus the integral over DeltaEz 
  intEr      = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, 
				  fgkRList, fgkZList, fgkPhiList, fLookUpErOverEz  );
  intEphi    = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, 
				  fgkRList, fgkZList, fgkPhiList, fLookUpEphiOverEz);
  intdEz = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, 
				  fgkRList, fgkZList, fgkPhiList, fLookUpDeltaEz   );

  // Calculate distorted position
  if ( r > 0.0 ) {
    phi =  phi + fCorrectionFactor *( fC0*intEphi - fC1*intEr ) / r;      
    r   =  r   + fCorrectionFactor *( fC0*intEr   + fC1*intEphi );  
  }
  Double_t dz = intdEz * fCorrectionFactor * fgkdvdE;
 
  // Calculate correction in cartesian coordinates
  dx[0] = - (r * TMath::Cos(phi) - x[0]);
  dx[1] = - (r * TMath::Sin(phi) - x[1]); 
  dx[2] = - dz;  // z distortion - (scaled with driftvelocity dependency on the Ez field and the overall scaling factor)

}

void AliTPCSpaceCharge3D::InitSpaceCharge3DDistortion() {
 //
  // Initialization of the Lookup table which contains the solutions of the 
  // "space charge" (poisson) problem - Faster and more accureate
  //
  // Method: Weighted sum-up of the different fields within the look up table 
  // but using two lookup tables with higher granularity in the (r,z) and the (rphi)- plane to emulate
  // more realistic space charges. (r,z) from primary ionisation. (rphi) for possible Gating leaks

  if (fInitLookUp) {
    AliInfo("Lookup table was already initialized!  Doing it again anyway ...");
    return;
  }
  
  // ------------------------------------------------------------------------------------------------------
  // step 1: lookup table in rz, fine grid, radial symetric, to emulate primary ionization

  AliInfo("Step 1: Preparation of the weighted look-up tables.");
   
  // lookup table in rz, fine grid

  TFile *fZR = new TFile(fSCLookUpPOCsFileNameRZ.Data(),"READ");
  if ( !fZR ) {
    AliError("Precalculated POC-looup-table in ZR could not be found");
    return;
  } 

  // units are in [m]
  TVector *gridf1  = (TVector*) fZR->Get("constants"); 
  TVector &grid1 = *gridf1;
  TMatrix *coordf1  = (TMatrix*) fZR->Get("coordinates");
  TMatrix &coord1 = *coordf1;
  TMatrix *coordPOCf1  = (TMatrix*) fZR->Get("POCcoord");
  TMatrix &coordPOC1 = *coordPOCf1;
  
  Int_t rows      = (Int_t)grid1(0);   // number of points in r direction  - from RZ or RPhi table
  Int_t phiSlices = (Int_t)grid1(1);   // number of points in phi          - from RPhi table
  Int_t columns   = (Int_t)grid1(2);   // number of points in z direction  - from RZ table

  Float_t gridSizeR   =  (fgkOFCRadius-fgkIFCRadius)/(rows-1);   // unit in [cm]
  Float_t gridSizeZ   =  fgkTPCZ0/(columns-1);                  // unit in [cm]

  // temporary matrices needed for the calculation // for rotational symmetric RZ table, phislices is 1
  
  TMatrixD *arrayofErA[kNPhiSlices], *arrayofdEzA[kNPhiSlices]; 
  TMatrixD *arrayofErC[kNPhiSlices], *arrayofdEzC[kNPhiSlices]; 

  TMatrixD *arrayofEroverEzA[kNPhiSlices], *arrayofDeltaEzA[kNPhiSlices];
  TMatrixD *arrayofEroverEzC[kNPhiSlices], *arrayofDeltaEzC[kNPhiSlices];

 
  for ( Int_t k = 0 ; k < phiSlices ; k++ ) {
   
    arrayofErA[k]   =   new TMatrixD(rows,columns) ;
    arrayofdEzA[k]  =   new TMatrixD(rows,columns) ;
    arrayofErC[k]   =   new TMatrixD(rows,columns) ;
    arrayofdEzC[k]  =   new TMatrixD(rows,columns) ;

    arrayofEroverEzA[k]   =   new TMatrixD(rows,columns) ;
    arrayofDeltaEzA[k]    =   new TMatrixD(rows,columns) ;
    arrayofEroverEzC[k]   =   new TMatrixD(rows,columns) ;
    arrayofDeltaEzC[k]    =   new TMatrixD(rows,columns) ;

    // zero initialization not necessary, it is done in the constructor of TMatrix 

  }
 
  // list of points as used during sum up
  Double_t  rlist1[kNRows], zedlist1[kNColumns];// , philist1[phiSlices];
  for ( Int_t i = 0 ; i < rows ; i++ )    {
    rlist1[i] = fgkIFCRadius + i*gridSizeR ;
    for ( Int_t j = 0 ; j < columns ; j++ ) { 
      zedlist1[j]  = j * gridSizeZ ;
    }
  }
  
  TTree *treePOC = (TTree*)fZR->Get("POCall");

  TVector *bEr  = 0;   //TVector *bEphi= 0;
  TVector *bEz  = 0;
  
  treePOC->SetBranchAddress("Er",&bEr);
  treePOC->SetBranchAddress("Ez",&bEz);


  // Read the complete tree and do a weighted sum-up over the POC configurations
  // +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  
  Int_t treeNumPOC = (Int_t)treePOC->GetEntries(); // Number of POC conf. in the look-up table
  Int_t ipC = 0; // POC Conf. counter (note: different to the POC number in the tree!)

  for (Int_t itreepC=0; itreepC<treeNumPOC; itreepC++) { // ------------- loop over POC configurations in tree
  
    treePOC->GetEntry(itreepC);

    // center of the POC voxel in [meter]
    Double_t r0 = coordPOC1(ipC,0);
    Double_t phi0 = coordPOC1(ipC,1);
    Double_t z0 = coordPOC1(ipC,2);

    ipC++; // POC configuration counter

    // weights (charge density) at POC position on the A and C side (in C/m^3/e0)
    // note: coordinates are in [cm]
    Double_t weightA = GetSpaceChargeDensity(r0*100,phi0, z0*100, 1);  // partial load in r,z
    Double_t weightC = GetSpaceChargeDensity(r0*100,phi0,-z0*100, 1);  // partial load in r,z
  
    // Summing up the vector components according to their weight

    Int_t ip = 0;
    for ( Int_t j = 0 ; j < columns ; j++ ) { 
      for ( Int_t i = 0 ; i < rows ; i++ )    {
	for ( Int_t k = 0 ; k < phiSlices ; k++ ) {
		  
	  // check wether the coordinates were screwed
	  if (TMath::Abs((coord1(0,ip)*100-rlist1[i]))>1 || 
	      TMath::Abs((coord1(2,ip)*100-zedlist1[j])>1)) { 
	    AliError("internal error: coordinate system was screwed during the sum-up");
	    printf("sum-up: (r,z)=(%f,%f)\n",rlist1[i],zedlist1[j]);
	    printf("lookup: (r,z)=(%f,%f)\n",coord1(0,ip)*100,coord1(2,ip)*100);
	    AliError("Don't trust the results of the space charge calculation!");
	  }
	  
	  // unfortunately, the lookup tables were produced to be faster for phi symmetric charges
	  // This will be the most frequent usage (hopefully)
	  // That's why we have to do this here ...

	  TMatrixD &erA   =  *arrayofErA[k]  ;
	  TMatrixD &dEzA  =  *arrayofdEzA[k]   ;
   
	  TMatrixD &erC   =  *arrayofErC[k]  ;
	  TMatrixD &dEzC  =  *arrayofdEzC[k]   ;
   
	  // Sum up - Efield values in [V/m] -> transition to [V/cm]
	  erA(i,j) += ((*bEr)(ip)) * weightA /100;
	  erC(i,j) += ((*bEr)(ip)) * weightC /100;
	  dEzA(i,j)  += ((*bEz)(ip)) * weightA /100;
	  dEzC(i,j)  += ((*bEz)(ip)) * weightC /100;

	  // increase the counter
	  ip++;
	}
      }
    } // end coordinate loop
  } // end POC loop


  // -------------------------------------------------------------------------------
  // Division by the Ez (drift) field and integration along z

  //  AliInfo("Step 1: Division and integration");

  Double_t ezField = (fgkCathodeV-fgkGG)/fgkTPCZ0; // = Electric Field (V/cm) Magnitude ~ -400 V/cm;

  for ( Int_t k = 0 ; k < phiSlices ; k++ ) { // phi loop

    // matrices holding the solution - summation of POC charges // see above
    TMatrixD &erA   =  *arrayofErA[k]  ;
    TMatrixD &dezA  =  *arrayofdEzA[k]   ;
    TMatrixD &erC   =  *arrayofErC[k]  ;
    TMatrixD &dezC  =  *arrayofdEzC[k]   ;

    // matrices which will contain the integrated fields (divided by the drift field)
    TMatrixD &erOverEzA   =  *arrayofEroverEzA[k]  ;
    TMatrixD &deltaEzA    =  *arrayofDeltaEzA[k];
    TMatrixD &erOverEzC   =  *arrayofEroverEzC[k]  ;
    TMatrixD &deltaEzC    =  *arrayofDeltaEzC[k];    
    
    for ( Int_t i = 0 ; i < rows ; i++ )    { // r loop
      for ( Int_t j = columns-1 ; j >= 0 ; j-- ) {// z loop 
	// Count backwards to facilitate integration over Z

	Int_t index = 1 ; // Simpsons rule if N=odd.If N!=odd then add extra point 
	                  // by trapezoidal rule.  

	erOverEzA(i,j) = 0; //ephiOverEzA(i,j) = 0;
	deltaEzA(i,j) = 0;
	erOverEzC(i,j) = 0; //ephiOverEzC(i,j) = 0; 
	deltaEzC(i,j) = 0;

	for ( Int_t m = j ; m < columns ; m++ ) { // integration

	  erOverEzA(i,j)   += index*(gridSizeZ/3.0)*erA(i,m)/(-1*ezField) ;
	  erOverEzC(i,j)   += index*(gridSizeZ/3.0)*erC(i,m)/(-1*ezField)  ;
	  deltaEzA(i,j)    += index*(gridSizeZ/3.0)*dezA(i,m)/(-1) ;
	  deltaEzC(i,j)    += index*(gridSizeZ/3.0)*dezC(i,m)/(-1) ;

	  if ( index != 4 )  index = 4; else index = 2 ;

	}

	if ( index == 4 ) {
	  erOverEzA(i,j)   -= (gridSizeZ/3.0)*erA(i,columns-1)/(-1*ezField) ;
	  erOverEzC(i,j)   -= (gridSizeZ/3.0)*erC(i,columns-1)/(-1*ezField) ;
	  deltaEzA(i,j)    -= (gridSizeZ/3.0)*dezA(i,columns-1)/(-1) ;
	  deltaEzC(i,j)    -= (gridSizeZ/3.0)*dezC(i,columns-1)/(-1) ;
	}
	if ( index == 2 ) {
	  erOverEzA(i,j)   += (gridSizeZ/3.0)*(0.5*erA(i,columns-2)-2.5*erA(i,columns-1))/(-1*ezField) ;
	  erOverEzC(i,j)   += (gridSizeZ/3.0)*(0.5*erC(i,columns-2)-2.5*erC(i,columns-1))/(-1*ezField) ;
	  deltaEzA(i,j)    += (gridSizeZ/3.0)*(0.5*dezA(i,columns-2)-2.5*dezA(i,columns-1))/(-1) ;
	  deltaEzC(i,j)    += (gridSizeZ/3.0)*(0.5*dezC(i,columns-2)-2.5*dezC(i,columns-1))/(-1) ;
	}
	if ( j == columns-2 ) {
	  erOverEzA(i,j)   = (gridSizeZ/3.0)*(1.5*erA(i,columns-2)+1.5*erA(i,columns-1))/(-1*ezField) ;
	  erOverEzC(i,j)   = (gridSizeZ/3.0)*(1.5*erC(i,columns-2)+1.5*erC(i,columns-1))/(-1*ezField) ;
	  deltaEzA(i,j)    = (gridSizeZ/3.0)*(1.5*dezA(i,columns-2)+1.5*dezA(i,columns-1))/(-1) ;
	  deltaEzC(i,j)    = (gridSizeZ/3.0)*(1.5*dezC(i,columns-2)+1.5*dezC(i,columns-1))/(-1) ;
	}
	if ( j == columns-1 ) {
	  erOverEzA(i,j)   = 0;   
	  erOverEzC(i,j)   = 0;
	  deltaEzA(i,j)    = 0;  
	  deltaEzC(i,j)    = 0;
	}
      }
    }

  }
  
  //  AliInfo("Step 1: Interpolation to Standard grid");

  // -------------------------------------------------------------------------------
  // Interpolate results onto the standard grid which is used for all AliTPCCorrections classes

  const Int_t order  = 1  ;  // Linear interpolation = 1, Quadratic = 2  

  Double_t  r, z;//phi, z ;
  for ( Int_t k = 0 ; k < kNPhi ; k++ ) {
    //    phi = fgkPhiList[k] ;
	
    // final lookup table
    TMatrixF &erOverEzFinal   =  *fLookUpErOverEz[k]  ;
    TMatrixF &deltaEzFinal    =  *fLookUpDeltaEz[k]   ;
	
    // calculated and integrated tables - just one phi slice
    TMatrixD &erOverEzA   =  *arrayofEroverEzA[0]  ;
    TMatrixD &deltaEzA    =  *arrayofDeltaEzA[0];
    TMatrixD &erOverEzC   =  *arrayofEroverEzC[0]  ;
    TMatrixD &deltaEzC    =  *arrayofDeltaEzC[0];    
    
    
    for ( Int_t j = 0 ; j < kNZ ; j++ ) {

      z = TMath::Abs(fgkZList[j]) ;  // z position is symmetric
    
      for ( Int_t i = 0 ; i < kNR ; i++ ) { 
	r = fgkRList[i] ;

	// Interpolate Lookup tables onto standard grid
	if (fgkZList[j]>0) {
	  erOverEzFinal(i,j)   = Interpolate2DTable(order, r, z, rows, columns, rlist1, zedlist1, erOverEzA  );
	  deltaEzFinal(i,j)    = Interpolate2DTable(order, r, z, rows, columns, rlist1, zedlist1, deltaEzA   );
	} else {
	  erOverEzFinal(i,j)   = Interpolate2DTable(order, r, z, rows, columns, rlist1, zedlist1, erOverEzC  );
	  deltaEzFinal(i,j)    = - Interpolate2DTable(order, r, z, rows, columns, rlist1, zedlist1, deltaEzC   );
	  // negative coordinate system on C side
	}

      } // end r loop
    } // end z loop
  } // end phi loop

  
  // clear the temporary arrays lists
  for ( Int_t k = 0 ; k < phiSlices ; k++ )  {

    delete arrayofErA[k];  
    delete arrayofdEzA[k];
    delete arrayofErC[k];  
    delete arrayofdEzC[k];

    delete arrayofEroverEzA[k];  
    delete arrayofDeltaEzA[k];
    delete arrayofEroverEzC[k];  
    delete arrayofDeltaEzC[k];

  }

  fZR->Close();

  // ------------------------------------------------------------------------------------------------------
  // Step 2: Load and sum up lookup table in rphi, fine grid, to emulate for example a GG leak
 
  //  AliInfo("Step 2: Preparation of the weighted look-up table");
 
  TFile *fRPhi = new TFile(fSCLookUpPOCsFileNameRPhi.Data(),"READ");
  if ( !fRPhi ) {
    AliError("Precalculated POC-looup-table in RPhi could not be found");
    return;
  } 

  // units are in [m]
  TVector *gridf2  = (TVector*) fRPhi->Get("constants"); 
  TVector &grid2 = *gridf2;
  TMatrix *coordf2  = (TMatrix*) fRPhi->Get("coordinates");
  TMatrix &coord2 = *coordf2;
  TMatrix *coordPOCf2  = (TMatrix*) fRPhi->Get("POCcoord");
  TMatrix &coordPOC2 = *coordPOCf2;

  rows      = (Int_t)grid2(0);   // number of points in r direction   
  phiSlices = (Int_t)grid2(1);   // number of points in phi           
  columns   = (Int_t)grid2(2);   // number of points in z direction   

  gridSizeR   =  (fgkOFCRadius-fgkIFCRadius)/(rows-1);   // unit in [cm]
  Float_t gridSizePhi =  TMath::TwoPi()/phiSlices;         // unit in [rad]
  gridSizeZ   =  fgkTPCZ0/(columns-1);                  // unit in [cm]
 
  // list of points as used during sum up
  Double_t  rlist2[kNRows], philist2[kNPhiSlices], zedlist2[kNColumns]; 
  for ( Int_t k = 0 ; k < phiSlices ; k++ ) {
    philist2[k] =  gridSizePhi * k;
    for ( Int_t i = 0 ; i < rows ; i++ )    {
      rlist2[i] = fgkIFCRadius + i*gridSizeR ;
      for ( Int_t j = 0 ; j < columns ; j++ ) { 
	zedlist2[j]  = j * gridSizeZ ;
      }
    }
  } // only done once
 
  // temporary matrices needed for the calculation 

  TMatrixD *arrayofErA2[kNPhiSlices], *arrayofEphiA2[kNPhiSlices], *arrayofdEzA2[kNPhiSlices];
  TMatrixD *arrayofErC2[kNPhiSlices], *arrayofEphiC2[kNPhiSlices], *arrayofdEzC2[kNPhiSlices]; 

  TMatrixD *arrayofEroverEzA2[kNPhiSlices], *arrayofEphioverEzA2[kNPhiSlices], *arrayofDeltaEzA2[kNPhiSlices]; 
  TMatrixD *arrayofEroverEzC2[kNPhiSlices], *arrayofEphioverEzC2[kNPhiSlices], *arrayofDeltaEzC2[kNPhiSlices]; 

 
  for ( Int_t k = 0 ; k < phiSlices ; k++ ) {
   
    arrayofErA2[k]   =   new TMatrixD(rows,columns) ;
    arrayofEphiA2[k] =   new TMatrixD(rows,columns) ;
    arrayofdEzA2[k]  =   new TMatrixD(rows,columns) ;
    arrayofErC2[k]   =   new TMatrixD(rows,columns) ;
    arrayofEphiC2[k] =   new TMatrixD(rows,columns) ;
    arrayofdEzC2[k]  =   new TMatrixD(rows,columns) ;

    arrayofEroverEzA2[k]   =   new TMatrixD(rows,columns) ;
    arrayofEphioverEzA2[k] =   new TMatrixD(rows,columns) ; 
    arrayofDeltaEzA2[k]    =   new TMatrixD(rows,columns) ;
    arrayofEroverEzC2[k]   =   new TMatrixD(rows,columns) ;
    arrayofEphioverEzC2[k] =   new TMatrixD(rows,columns) ; 
    arrayofDeltaEzC2[k]    =   new TMatrixD(rows,columns) ;

    // zero initialization not necessary, it is done in the constructor of TMatrix 

  }
 
    
  treePOC = (TTree*)fRPhi->Get("POCall");

  //  TVector *bEr  = 0;   // done above
  TVector *bEphi= 0;
  //  TVector *bEz  = 0;   // done above

  treePOC->SetBranchAddress("Er",&bEr);
  treePOC->SetBranchAddress("Ephi",&bEphi);
  treePOC->SetBranchAddress("Ez",&bEz);

  // Read the complete tree and do a weighted sum-up over the POC configurations
  // +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  
  treeNumPOC = (Int_t)treePOC->GetEntries(); // Number of POC conf. in the look-up table
  ipC = 0; // POC Conf. counter (note: different to the POC number in the tree!)

  for (Int_t itreepC=0; itreepC<treeNumPOC; itreepC++) { // ------------- loop over POC configurations in tree
  
    treePOC->GetEntry(itreepC);

    // center of the POC voxel in [meter]
    Double_t r0 = coordPOC2(ipC,0);
    Double_t phi0 = coordPOC2(ipC,1);
    //    Double_t z0 = coordPOC2(ipC,2);

     // weights (charge density) at POC position on the A and C side (in C/m^3/e0)
    // note: coordinates are in [cm]
    Double_t weightA = GetSpaceChargeDensity(r0*100,phi0, 0.499, 2);  // partial load in r,phi
    Double_t weightC = GetSpaceChargeDensity(r0*100,phi0,-0.499, 2);  // partial load in r,phi
  
    //    printf("-----\n%f %f : %e %e\n",r0,phi0,weightA,weightC);

    // Summing up the vector components according to their weight

    Int_t ip = 0;
    for ( Int_t j = 0 ; j < columns ; j++ ) { 
      for ( Int_t i = 0 ; i < rows ; i++ )    {
	for ( Int_t k = 0 ; k < phiSlices ; k++ ) {
		  
	  // check wether the coordinates were screwed
	  if (TMath::Abs((coord2(0,ip)*100-rlist2[i]))>1 || 
	      TMath::Abs((coord2(1,ip)-philist2[k]))>1 || 
	      TMath::Abs((coord2(2,ip)*100-zedlist2[j]))>1) { 
	    AliError("internal error: coordinate system was screwed during the sum-up");
	    printf("lookup: (r,phi,z)=(%f,%f,%f)\n",coord2(0,ip)*100,coord2(1,ip),coord2(2,ip)*100);
	    printf("sum-up: (r,phi,z)=(%f,%f,%f)\n",rlist2[i],philist2[k],zedlist2[j]);
	    AliError("Don't trust the results of the space charge calculation!");
	  }
	  
	  // unfortunately, the lookup tables were produced to be faster for phi symmetric charges
	  // This will be the most frequent usage (hopefully)
	  // That's why we have to do this here ...

	  TMatrixD &erA   =  *arrayofErA2[k]  ;
	  TMatrixD &ephiA =  *arrayofEphiA2[k];
	  TMatrixD &dEzA  =  *arrayofdEzA2[k] ;
   
	  TMatrixD &erC   =  *arrayofErC2[k]  ;
	  TMatrixD &ephiC =  *arrayofEphiC2[k];
	  TMatrixD &dEzC  =  *arrayofdEzC2[k]   ;
   
	  // Sum up - Efield values in [V/m] -> transition to [V/cm]
	  erA(i,j) += ((*bEr)(ip)) * weightA /100;
	  erC(i,j) += ((*bEr)(ip)) * weightC /100;
	  ephiA(i,j) += ((*bEphi)(ip)) * weightA/100; // [V/rad]
	  ephiC(i,j) += ((*bEphi)(ip)) * weightC/100; // [V/rad]
	  dEzA(i,j)  += ((*bEz)(ip)) * weightA /100;
	  dEzC(i,j)  += ((*bEz)(ip)) * weightC /100;

	  // increase the counter
	  ip++;
	}
      }
    } // end coordinate loop
    
    
    // Rotation and summation in the rest of the dPhiSteps
    // which were not stored in the this tree due to storage & symmetry reasons

    
    Int_t phiPoints = (Int_t) grid2(1);
    Int_t phiPOC    = (Int_t) grid2(4);
    
    //   printf("%d %d\n",phiPOC,flagRadSym);
    
    for (Int_t phiiC = 1; phiiC<phiPOC; phiiC++) { // just used for non-radial symetric table 
      
      Double_t phi0R = phiiC*phi0*2 + phi0; // rotate further

      // weights (charge density) at POC position on the A and C side (in C/m^3/e0)
      // note: coordinates are in [cm] // ecxept z
      weightA = GetSpaceChargeDensity(r0*100,phi0R, 0.499, 2);  // partial load in r,phi
      weightC = GetSpaceChargeDensity(r0*100,phi0R,-0.499, 2);  // partial load in r,phi
    
      // printf("%f %f : %e %e\n",r0,phi0R,weightA,weightC);
         
      // Summing up the vector components according to their weight
      ip = 0;
      for ( Int_t j = 0 ; j < columns ; j++ ) { 
	for ( Int_t i = 0 ; i < rows ; i++ )    {
	  for ( Int_t k = 0 ; k < phiSlices ; k++ ) {
	    
	    // Note: rotating the coordinated during the sum up
	    
	    Int_t rotVal = (phiPoints/phiPOC)*phiiC;
	    Int_t ipR = -1;
	    
	    if ((ip%phiPoints)>=rotVal) {
	      ipR = ip-rotVal;
	    } else {
	      ipR = ip+(phiPoints-rotVal);
	    }
	    
	    // unfortunately, the lookup tables were produced to be faster for phi symmetric charges
	    // This will be the most frequent usage 
	    // That's why we have to do this here and not outside the loop ...
	    
	    TMatrixD &erA   =  *arrayofErA2[k]  ;
	    TMatrixD &ephiA =  *arrayofEphiA2[k];
	    TMatrixD &dEzA  =  *arrayofdEzA2[k]   ;
	    
	    TMatrixD &erC   =  *arrayofErC2[k]  ;
	    TMatrixD &ephiC =  *arrayofEphiC2[k];
	    TMatrixD &dEzC  =  *arrayofdEzC2[k]   ;
       
	    // Sum up - Efield values in [V/m] -> transition to [V/cm]
	    erA(i,j) += ((*bEr)(ipR)) * weightA /100;
	    erC(i,j) += ((*bEr)(ipR)) * weightC /100;
	    ephiA(i,j) += ((*bEphi)(ipR)) * weightA/100; // [V/rad]
	    ephiC(i,j) += ((*bEphi)(ipR)) * weightC/100; // [V/rad]
	    dEzA(i,j)  += ((*bEz)(ipR)) * weightA /100;
	    dEzC(i,j)  += ((*bEz)(ipR)) * weightC /100;

	    // increase the counter
	    ip++;
	  }
	}
      } // end coordinate loop

    } // end phi-POC summation (phiiC)

    ipC++; // POC configuration counter

    //   printf("POC: (r,phi,z) = (%f %f %f) | weight(A,C): %03.1lf %03.1lf\n",r0,phi0,z0, weightA, weightC);
    
  }




  // -------------------------------------------------------------------------------
  // Division by the Ez (drift) field and integration along z

  //  AliInfo("Step 2: Division and integration");


  for ( Int_t k = 0 ; k < phiSlices ; k++ ) { // phi loop

    // matrices holding the solution - summation of POC charges // see above
    TMatrixD &erA   =  *arrayofErA2[k]  ;
    TMatrixD &ephiA =  *arrayofEphiA2[k];
    TMatrixD &dezA  =  *arrayofdEzA2[k]   ;
    TMatrixD &erC   =  *arrayofErC2[k]  ;
    TMatrixD &ephiC =  *arrayofEphiC2[k];
    TMatrixD &dezC  =  *arrayofdEzC2[k]   ;

    // matrices which will contain the integrated fields (divided by the drift field)
    TMatrixD &erOverEzA   =  *arrayofEroverEzA2[k]  ;
    TMatrixD &ephiOverEzA =  *arrayofEphioverEzA2[k];
    TMatrixD &deltaEzA    =  *arrayofDeltaEzA2[k];
    TMatrixD &erOverEzC   =  *arrayofEroverEzC2[k]  ;
    TMatrixD &ephiOverEzC =  *arrayofEphioverEzC2[k];
    TMatrixD &deltaEzC    =  *arrayofDeltaEzC2[k];    
    
    for ( Int_t i = 0 ; i < rows ; i++ )    { // r loop
      for ( Int_t j = columns-1 ; j >= 0 ; j-- ) {// z loop 
	// Count backwards to facilitate integration over Z

	Int_t index = 1 ; // Simpsons rule if N=odd.If N!=odd then add extra point by trapezoidal rule.  

	erOverEzA(i,j) = 0; 
	ephiOverEzA(i,j) = 0;
	deltaEzA(i,j) = 0;
	erOverEzC(i,j) = 0; 
	ephiOverEzC(i,j) = 0; 
	deltaEzC(i,j) = 0;

	for ( Int_t m = j ; m < columns ; m++ ) { // integration

	  erOverEzA(i,j)   += index*(gridSizeZ/3.0)*erA(i,m)/(-1*ezField) ;
	  erOverEzC(i,j)   += index*(gridSizeZ/3.0)*erC(i,m)/(-1*ezField)  ;
	  ephiOverEzA(i,j) += index*(gridSizeZ/3.0)*ephiA(i,m)/(-1*ezField)  ;
	  ephiOverEzC(i,j) += index*(gridSizeZ/3.0)*ephiC(i,m)/(-1*ezField)  ;
	  deltaEzA(i,j)    += index*(gridSizeZ/3.0)*dezA(i,m)/(-1) ;
	  deltaEzC(i,j)    += index*(gridSizeZ/3.0)*dezC(i,m)/(-1) ;

	  if ( index != 4 )  index = 4; else index = 2 ;

	}

	if ( index == 4 ) {
	  erOverEzA(i,j)   -= (gridSizeZ/3.0)*erA(i,columns-1)/(-1*ezField) ;
	  erOverEzC(i,j)   -= (gridSizeZ/3.0)*erC(i,columns-1)/(-1*ezField) ;
	  ephiOverEzA(i,j) -= (gridSizeZ/3.0)*ephiA(i,columns-1)/(-1*ezField) ;
	  ephiOverEzC(i,j) -= (gridSizeZ/3.0)*ephiC(i,columns-1)/(-1*ezField) ;
	  deltaEzA(i,j)    -= (gridSizeZ/3.0)*dezA(i,columns-1)/(-1) ;
	  deltaEzC(i,j)    -= (gridSizeZ/3.0)*dezC(i,columns-1)/(-1) ;
	}
	if ( index == 2 ) {
	  erOverEzA(i,j)   += (gridSizeZ/3.0)*(0.5*erA(i,columns-2)-2.5*erA(i,columns-1))/(-1*ezField) ;
	  erOverEzC(i,j)   += (gridSizeZ/3.0)*(0.5*erC(i,columns-2)-2.5*erC(i,columns-1))/(-1*ezField) ;
	  ephiOverEzA(i,j) += (gridSizeZ/3.0)*(0.5*ephiA(i,columns-2)-2.5*ephiA(i,columns-1))/(-1*ezField) ;
	  ephiOverEzC(i,j) += (gridSizeZ/3.0)*(0.5*ephiC(i,columns-2)-2.5*ephiC(i,columns-1))/(-1*ezField) ;
	  deltaEzA(i,j)    += (gridSizeZ/3.0)*(0.5*dezA(i,columns-2)-2.5*dezA(i,columns-1))/(-1) ;
	  deltaEzC(i,j)    += (gridSizeZ/3.0)*(0.5*dezC(i,columns-2)-2.5*dezC(i,columns-1))/(-1) ;
	}
	if ( j == columns-2 ) {
	  erOverEzA(i,j)   = (gridSizeZ/3.0)*(1.5*erA(i,columns-2)+1.5*erA(i,columns-1))/(-1*ezField) ;
	  erOverEzC(i,j)   = (gridSizeZ/3.0)*(1.5*erC(i,columns-2)+1.5*erC(i,columns-1))/(-1*ezField) ;
	  ephiOverEzA(i,j) = (gridSizeZ/3.0)*(1.5*ephiA(i,columns-2)+1.5*ephiA(i,columns-1))/(-1*ezField) ;
	  ephiOverEzC(i,j) = (gridSizeZ/3.0)*(1.5*ephiC(i,columns-2)+1.5*ephiC(i,columns-1))/(-1*ezField) ;
	  deltaEzA(i,j)    = (gridSizeZ/3.0)*(1.5*dezA(i,columns-2)+1.5*dezA(i,columns-1))/(-1) ;
	  deltaEzC(i,j)    = (gridSizeZ/3.0)*(1.5*dezC(i,columns-2)+1.5*dezC(i,columns-1))/(-1) ;
	}
	if ( j == columns-1 ) {
	  erOverEzA(i,j)   = 0;   
	  erOverEzC(i,j)   = 0;
	  ephiOverEzA(i,j) = 0; 
	  ephiOverEzC(i,j) = 0;
	  deltaEzA(i,j)    = 0;  
	  deltaEzC(i,j)    = 0;
	}
      }
    }

  }
  
  AliInfo("Step 2: Interpolation to Standard grid");

  // -------------------------------------------------------------------------------
  // Interpolate results onto the standard grid which is used for all AliTPCCorrections classes


  for ( Int_t k = 0 ; k < kNPhi ; k++ ) {
    Double_t phi = fgkPhiList[k] ;
	
    // final lookup table
    TMatrixF &erOverEzFinal   =  *fLookUpErOverEz[k]  ;
    TMatrixF &ephiOverEzFinal =  *fLookUpEphiOverEz[k];
    TMatrixF &deltaEzFinal    =  *fLookUpDeltaEz[k]   ;
	
    for ( Int_t j = 0 ; j < kNZ ; j++ ) {

      z = TMath::Abs(fgkZList[j]) ;  // z position is symmetric
    
      for ( Int_t i = 0 ; i < kNR ; i++ ) { 
	r = fgkRList[i] ;

	// Interpolate Lookup tables onto standard grid
	if (fgkZList[j]>0) {
	  erOverEzFinal(i,j)   += Interpolate3DTable(order, r, z, phi, rows, columns, phiSlices, 
					       rlist2, zedlist2, philist2, arrayofEroverEzA2  );
	  ephiOverEzFinal(i,j) += Interpolate3DTable(order, r, z, phi, rows, columns, phiSlices,
					       rlist2, zedlist2, philist2, arrayofEphioverEzA2);
	  deltaEzFinal(i,j)    += Interpolate3DTable(order, r, z, phi, rows, columns, phiSlices,
					       rlist2, zedlist2, philist2, arrayofDeltaEzA2   );
	} else {
	  erOverEzFinal(i,j)   += Interpolate3DTable(order, r, z, phi, rows, columns, phiSlices, 
					       rlist2, zedlist2, philist2, arrayofEroverEzC2  );
	  ephiOverEzFinal(i,j) += Interpolate3DTable(order, r, z, phi, rows, columns, phiSlices,
					       rlist2, zedlist2, philist2, arrayofEphioverEzC2);
	  deltaEzFinal(i,j)  -=  Interpolate3DTable(order, r, z, phi, rows, columns, phiSlices,
					       rlist2, zedlist2, philist2, arrayofDeltaEzC2   );
	}

      } // end r loop
    } // end z loop
  } // end phi loop
  
  
  // clear the temporary arrays lists
  for ( Int_t k = 0 ; k < phiSlices ; k++ )  {

    delete arrayofErA2[k];  
    delete arrayofEphiA2[k];
    delete arrayofdEzA2[k];
    delete arrayofErC2[k];  
    delete arrayofEphiC2[k];
    delete arrayofdEzC2[k];

    delete arrayofEroverEzA2[k];  
    delete arrayofEphioverEzA2[k];
    delete arrayofDeltaEzA2[k];
    delete arrayofEroverEzC2[k];  
    delete arrayofEphioverEzC2[k];
    delete arrayofDeltaEzC2[k];

  }
 
  fRPhi->Close();
  
  // FINISHED

  fInitLookUp = kTRUE;

}

void AliTPCSpaceCharge3D::InitSpaceCharge3DDistortionCourse() {
  //
  // Initialization of the Lookup table which contains the solutions of the 
  // "space charge" (poisson) problem 
  //
  // The sum-up uses a look-up table which contains different discretized Space charge fields 
  // in order to calculate the corresponding field deviations due to a given (discretized)
  // space charge distribution ....
  //
  // Method of calculation: Weighted sum-up of the different fields within the look up table
  // Note: Full 3d version: Course and slow ...

  if (fInitLookUp) {
    AliInfo("Lookup table was already initialized!");
    //    return;
  }

  AliInfo("Preparation of the weighted look-up table");
   
  TFile *f = new TFile(fSCLookUpPOCsFileName3D.Data(),"READ");
  if ( !f ) {
    AliError("Precalculated POC-looup-table could not be found");
    return;
  }

  // units are in [m]
  TVector *gridf  = (TVector*) f->Get("constants"); 
  TVector &grid = *gridf;
  TMatrix *coordf  = (TMatrix*) f->Get("coordinates");
  TMatrix &coord = *coordf;
  TMatrix *coordPOCf  = (TMatrix*) f->Get("POCcoord");
  TMatrix &coordPOC = *coordPOCf;
  
  Bool_t flagRadSym = 0;
  if (grid(1)==1 && grid(4)==1) {
    //   AliInfo("LOOK UP TABLE IS RADIAL SYMETTRIC - Field in Phi is ZERO");
    flagRadSym=1;
  }

  Int_t rows      = (Int_t)grid(0);   // number of points in r direction 
  Int_t phiSlices = (Int_t)grid(1);   // number of points in phi         
  Int_t columns   = (Int_t)grid(2);   // number of points in z direction 

  const Float_t gridSizeR   =  (fgkOFCRadius-fgkIFCRadius)/(rows-1);   // unit in [cm]
  const Float_t gridSizePhi =  TMath::TwoPi()/phiSlices;         // unit in [rad]
  const Float_t gridSizeZ   =  fgkTPCZ0/(columns-1);                  // unit in [cm]
 
  // temporary matrices needed for the calculation
  TMatrixD *arrayofErA[kNPhiSlices], *arrayofEphiA[kNPhiSlices], *arrayofdEzA[kNPhiSlices];
  TMatrixD *arrayofErC[kNPhiSlices], *arrayofEphiC[kNPhiSlices], *arrayofdEzC[kNPhiSlices];

  TMatrixD *arrayofEroverEzA[kNPhiSlices], *arrayofEphioverEzA[kNPhiSlices], *arrayofDeltaEzA[kNPhiSlices];
  TMatrixD *arrayofEroverEzC[kNPhiSlices], *arrayofEphioverEzC[kNPhiSlices], *arrayofDeltaEzC[kNPhiSlices];

 
  for ( Int_t k = 0 ; k < phiSlices ; k++ ) {
   
    arrayofErA[k]   =   new TMatrixD(rows,columns) ;
    arrayofEphiA[k] =   new TMatrixD(rows,columns) ; // zero if radial symmetric
    arrayofdEzA[k]  =   new TMatrixD(rows,columns) ;
    arrayofErC[k]   =   new TMatrixD(rows,columns) ;
    arrayofEphiC[k] =   new TMatrixD(rows,columns) ; // zero if radial symmetric
    arrayofdEzC[k]  =   new TMatrixD(rows,columns) ;

    arrayofEroverEzA[k]   =   new TMatrixD(rows,columns) ;
    arrayofEphioverEzA[k] =   new TMatrixD(rows,columns) ; // zero if radial symmetric
    arrayofDeltaEzA[k]    =   new TMatrixD(rows,columns) ;
    arrayofEroverEzC[k]   =   new TMatrixD(rows,columns) ;
    arrayofEphioverEzC[k] =   new TMatrixD(rows,columns) ; // zero if radial symmetric
    arrayofDeltaEzC[k]    =   new TMatrixD(rows,columns) ;

    // Set the values to zero the lookup tables
    // not necessary, it is done in the constructor of TMatrix - code deleted

  }
 
  // list of points as used in the interpolation (during sum up)
  Double_t  rlist[kNRows], zedlist[kNColumns] , philist[kNPhiSlices];
  for ( Int_t k = 0 ; k < phiSlices ; k++ ) {
    philist[k] =  gridSizePhi * k;
    for ( Int_t i = 0 ; i < rows ; i++ )    {
      rlist[i] = fgkIFCRadius + i*gridSizeR ;
      for ( Int_t j = 0 ; j < columns ; j++ ) { 
	zedlist[j]  = j * gridSizeZ ;
      }
    }
  } // only done once
  
  
  TTree *treePOC = (TTree*)f->Get("POCall");

  TVector *bEr  = 0;   TVector *bEphi= 0;   TVector *bEz  = 0;
  
  treePOC->SetBranchAddress("Er",&bEr);
  if (!flagRadSym) treePOC->SetBranchAddress("Ephi",&bEphi);
  treePOC->SetBranchAddress("Ez",&bEz);


  // Read the complete tree and do a weighted sum-up over the POC configurations
  // +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  
  Int_t treeNumPOC = (Int_t)treePOC->GetEntries(); // Number of POC conf. in the look-up table
  Int_t ipC = 0; // POC Conf. counter (note: different to the POC number in the tree!)

  for (Int_t itreepC=0; itreepC<treeNumPOC; itreepC++) { // ------------- loop over POC configurations in tree
  
    treePOC->GetEntry(itreepC);

    // center of the POC voxel in [meter]
    Double_t r0 = coordPOC(ipC,0);
    Double_t phi0 = coordPOC(ipC,1);
    Double_t z0 = coordPOC(ipC,2);

    ipC++; // POC configuration counter

    // weights (charge density) at POC position on the A and C side (in C/m^3/e0)
    // note: coordinates are in [cm]
    Double_t weightA = GetSpaceChargeDensity(r0*100,phi0, z0*100,0);
    Double_t weightC = GetSpaceChargeDensity(r0*100,phi0,-z0*100,0);
  
    // Summing up the vector components according to their weight

    Int_t ip = 0;
    for ( Int_t j = 0 ; j < columns ; j++ ) { 
      for ( Int_t i = 0 ; i < rows ; i++ )    {
	for ( Int_t k = 0 ; k < phiSlices ; k++ ) {
		  
	  // check wether the coordinates were screwed
	  if (TMath::Abs((coord(0,ip)*100-rlist[i]))>1 || 
	      TMath::Abs((coord(1,ip)-philist[k]))>1 || 
	      TMath::Abs((coord(2,ip)*100-zedlist[j]))>1) { 
	    AliError("internal error: coordinate system was screwed during the sum-up");
	    printf("lookup: (r,phi,z)=(%f,%f,%f)\n",coord(0,ip)*100,coord(1,ip),coord(2,ip)*100);
	    printf("sum-up: (r,phi,z)=(%f,%f,%f)\n",rlist[i],philist[k],zedlist[j]);
	    AliError("Don't trust the results of the space charge calculation!");
	  }
	  
	  // unfortunately, the lookup tables were produced to be faster for phi symmetric charges
	  // This will be the most frequent usage (hopefully)
	  // That's why we have to do this here ...

	  TMatrixD &erA   =  *arrayofErA[k]  ;
	  TMatrixD &ephiA =  *arrayofEphiA[k];
	  TMatrixD &dEzA  =  *arrayofdEzA[k]   ;
   
	  TMatrixD &erC   =  *arrayofErC[k]  ;
	  TMatrixD &ephiC =  *arrayofEphiC[k];
	  TMatrixD &dEzC  =  *arrayofdEzC[k]   ;
   
	  // Sum up - Efield values in [V/m] -> transition to [V/cm]
	  erA(i,j) += ((*bEr)(ip)) * weightA /100;
	  erC(i,j) += ((*bEr)(ip)) * weightC /100;
	  if (!flagRadSym) {
	    ephiA(i,j) += ((*bEphi)(ip)) * weightA/100; // [V/rad]
	    ephiC(i,j) += ((*bEphi)(ip)) * weightC/100; // [V/rad]
	  }
	  dEzA(i,j)  += ((*bEz)(ip)) * weightA /100;
	  dEzC(i,j)  += ((*bEz)(ip)) * weightC /100;

	  // increase the counter
	  ip++;
	}
      }
    } // end coordinate loop
    
    
    // Rotation and summation in the rest of the dPhiSteps
    // which were not stored in the this tree due to storage & symmetry reasons

    Int_t phiPoints = (Int_t) grid(1);
    Int_t phiPOC    = (Int_t) grid(4);
    
    //   printf("%d %d\n",phiPOC,flagRadSym);
    
    for (Int_t phiiC = 1; phiiC<phiPOC; phiiC++) { // just used for non-radial symetric table 
      
      r0 = coordPOC(ipC,0);
      phi0 = coordPOC(ipC,1);
      z0 = coordPOC(ipC,2);
      
      ipC++; // POC conf. counter
      
      // weights (charge density) at POC position on the A and C side (in C/m^3/e0)
      // note: coordinates are in [cm]
      weightA = GetSpaceChargeDensity(r0*100,phi0, z0*100,0); 
      weightC = GetSpaceChargeDensity(r0*100,phi0,-z0*100,0);
      
      //     printf("%f %f %f: %e %e\n",r0,phi0,z0,weightA,weightC);
       
      // Summing up the vector components according to their weight
      ip = 0;
      for ( Int_t j = 0 ; j < columns ; j++ ) { 
	for ( Int_t i = 0 ; i < rows ; i++ )    {
	  for ( Int_t k = 0 ; k < phiSlices ; k++ ) {
	    
	    // Note: rotating the coordinated during the sum up
	    
	    Int_t rotVal = (phiPoints/phiPOC)*phiiC;
	    Int_t ipR = -1;
	    
	    if ((ip%phiPoints)>=rotVal) {
	      ipR = ip-rotVal;
	    } else {
	      ipR = ip+(phiPoints-rotVal);
	    }
	    
	    // unfortunately, the lookup tables were produced to be faster for phi symmetric charges
	    // This will be the most frequent usage 
	    // That's why we have to do this here and not outside the loop ...
	    
	    TMatrixD &erA   =  *arrayofErA[k]  ;
	    TMatrixD &ephiA =  *arrayofEphiA[k];
	    TMatrixD &dEzA  =  *arrayofdEzA[k]   ;
	    
	    TMatrixD &erC   =  *arrayofErC[k]  ;
	    TMatrixD &ephiC =  *arrayofEphiC[k];
	    TMatrixD &dEzC  =  *arrayofdEzC[k]   ;
       
	    // Sum up - Efield values in [V/m] -> transition to [V/cm]
	    erA(i,j) += ((*bEr)(ipR)) * weightA /100;
	    erC(i,j) += ((*bEr)(ipR)) * weightC /100;
	    if (!flagRadSym) {
	      ephiA(i,j) += ((*bEphi)(ipR)) * weightA/100; // [V/rad]
	      ephiC(i,j) += ((*bEphi)(ipR)) * weightC/100; // [V/rad]
	    }
	    dEzA(i,j)  += ((*bEz)(ipR)) * weightA /100;
	    dEzC(i,j)  += ((*bEz)(ipR)) * weightC /100;

	    // increase the counter
	    ip++;
	  }
	}
      } // end coordinate loop

    } // end phi-POC summation (phiiC)
   

    // printf("POC: (r,phi,z) = (%f %f %f) | weight(A,C): %03.1lf %03.1lf\n",r0,phi0,z0, weightA, weightC);
    
  }



  // -------------------------------------------------------------------------------
  // Division by the Ez (drift) field and integration along z

  AliInfo("Division and integration");

  Double_t ezField = (fgkCathodeV-fgkGG)/fgkTPCZ0; // = Electric Field (V/cm) Magnitude ~ -400 V/cm;

  for ( Int_t k = 0 ; k < phiSlices ; k++ ) { // phi loop

    // matrices holding the solution - summation of POC charges // see above
    TMatrixD &erA   =  *arrayofErA[k]  ;
    TMatrixD &ephiA =  *arrayofEphiA[k];
    TMatrixD &dezA  =  *arrayofdEzA[k]   ;
    TMatrixD &erC   =  *arrayofErC[k]  ;
    TMatrixD &ephiC =  *arrayofEphiC[k];
    TMatrixD &dezC  =  *arrayofdEzC[k]   ;

    // matrices which will contain the integrated fields (divided by the drift field)
    TMatrixD &erOverEzA   =  *arrayofEroverEzA[k]  ;
    TMatrixD &ephiOverEzA =  *arrayofEphioverEzA[k];
    TMatrixD &deltaEzA    =  *arrayofDeltaEzA[k];
    TMatrixD &erOverEzC   =  *arrayofEroverEzC[k]  ;
    TMatrixD &ephiOverEzC =  *arrayofEphioverEzC[k];
    TMatrixD &deltaEzC    =  *arrayofDeltaEzC[k];    
    
    for ( Int_t i = 0 ; i < rows ; i++ )    { // r loop
      for ( Int_t j = columns-1 ; j >= 0 ; j-- ) {// z loop 
	// Count backwards to facilitate integration over Z

	Int_t index = 1 ; // Simpsons rule if N=odd.If N!=odd then add extra point by trapezoidal rule.  

	erOverEzA(i,j) = 0; ephiOverEzA(i,j) = 0;  deltaEzA(i,j) = 0;
	erOverEzC(i,j) = 0; ephiOverEzC(i,j) = 0;  deltaEzC(i,j) = 0;

	for ( Int_t m = j ; m < columns ; m++ ) { // integration

	  erOverEzA(i,j)   += index*(gridSizeZ/3.0)*erA(i,m)/(-1*ezField) ;
	  erOverEzC(i,j)   += index*(gridSizeZ/3.0)*erC(i,m)/(-1*ezField)  ;
	  if (!flagRadSym) {
	    ephiOverEzA(i,j) += index*(gridSizeZ/3.0)*ephiA(i,m)/(-1*ezField)  ;
	    ephiOverEzC(i,j) += index*(gridSizeZ/3.0)*ephiC(i,m)/(-1*ezField)  ;
	  }
	  deltaEzA(i,j)    += index*(gridSizeZ/3.0)*dezA(i,m)/(-1) ;
	  deltaEzC(i,j)    += index*(gridSizeZ/3.0)*dezC(i,m)/(-1) ;

	  if ( index != 4 )  index = 4; else index = 2 ;

	}

	if ( index == 4 ) {
	  erOverEzA(i,j)   -= (gridSizeZ/3.0)*erA(i,columns-1)/(-1*ezField) ;
	  erOverEzC(i,j)   -= (gridSizeZ/3.0)*erC(i,columns-1)/(-1*ezField) ;
	  if (!flagRadSym) {
	    ephiOverEzA(i,j) -= (gridSizeZ/3.0)*ephiA(i,columns-1)/(-1*ezField) ;
	    ephiOverEzC(i,j) -= (gridSizeZ/3.0)*ephiC(i,columns-1)/(-1*ezField) ;
	  }
	  deltaEzA(i,j)    -= (gridSizeZ/3.0)*dezA(i,columns-1)/(-1) ;
	  deltaEzC(i,j)    -= (gridSizeZ/3.0)*dezC(i,columns-1)/(-1) ;
	}
	if ( index == 2 ) {
	  erOverEzA(i,j)   += (gridSizeZ/3.0)*(0.5*erA(i,columns-2)-2.5*erA(i,columns-1))/(-1*ezField) ;
	  erOverEzC(i,j)   += (gridSizeZ/3.0)*(0.5*erC(i,columns-2)-2.5*erC(i,columns-1))/(-1*ezField) ;
	  if (!flagRadSym) {
	    ephiOverEzA(i,j) += (gridSizeZ/3.0)*(0.5*ephiA(i,columns-2)-2.5*ephiA(i,columns-1))/(-1*ezField) ;
	    ephiOverEzC(i,j) += (gridSizeZ/3.0)*(0.5*ephiC(i,columns-2)-2.5*ephiC(i,columns-1))/(-1*ezField) ;
	  }
	  deltaEzA(i,j)    += (gridSizeZ/3.0)*(0.5*dezA(i,columns-2)-2.5*dezA(i,columns-1))/(-1) ;
	  deltaEzC(i,j)    += (gridSizeZ/3.0)*(0.5*dezC(i,columns-2)-2.5*dezC(i,columns-1))/(-1) ;
	}
	if ( j == columns-2 ) {
	  erOverEzA(i,j)   = (gridSizeZ/3.0)*(1.5*erA(i,columns-2)+1.5*erA(i,columns-1))/(-1*ezField) ;
	  erOverEzC(i,j)   = (gridSizeZ/3.0)*(1.5*erC(i,columns-2)+1.5*erC(i,columns-1))/(-1*ezField) ;
	  if (!flagRadSym) {
	    ephiOverEzA(i,j) = (gridSizeZ/3.0)*(1.5*ephiA(i,columns-2)+1.5*ephiA(i,columns-1))/(-1*ezField) ;
	    ephiOverEzC(i,j) = (gridSizeZ/3.0)*(1.5*ephiC(i,columns-2)+1.5*ephiC(i,columns-1))/(-1*ezField) ;
	  }
	  deltaEzA(i,j)    = (gridSizeZ/3.0)*(1.5*dezA(i,columns-2)+1.5*dezA(i,columns-1))/(-1) ;
	  deltaEzC(i,j)    = (gridSizeZ/3.0)*(1.5*dezC(i,columns-2)+1.5*dezC(i,columns-1))/(-1) ;
	}
	if ( j == columns-1 ) {
	  erOverEzA(i,j)   = 0;   
	  erOverEzC(i,j)   = 0;
	  if (!flagRadSym) {
	    ephiOverEzA(i,j) = 0; 
	    ephiOverEzC(i,j) = 0;
	  }
	  deltaEzA(i,j)    = 0;  
	  deltaEzC(i,j)    = 0;
	}
      }
    }

  }
  

 
  AliInfo("Interpolation to Standard grid");

  // -------------------------------------------------------------------------------
  // Interpolate results onto the standard grid which is used for all AliTPCCorrections classes

  const Int_t order  = 1  ;  // Linear interpolation = 1, Quadratic = 2  

  Double_t  r, phi, z ;
  for ( Int_t k = 0 ; k < kNPhi ; k++ ) {
    phi = fgkPhiList[k] ;
	
    TMatrixF &erOverEz   =  *fLookUpErOverEz[k]  ;
    TMatrixF &ephiOverEz =  *fLookUpEphiOverEz[k];
    TMatrixF &deltaEz    =  *fLookUpDeltaEz[k]   ;
	
    for ( Int_t j = 0 ; j < kNZ ; j++ ) {

      z = TMath::Abs(fgkZList[j]) ;  // z position is symmetric
    
      for ( Int_t i = 0 ; i < kNR ; i++ ) { 
	r = fgkRList[i] ;

	// Interpolate Lookup tables onto standard grid
	if (fgkZList[j]>0) {
	  erOverEz(i,j)   = Interpolate3DTable(order, r, z, phi, rows, columns, phiSlices, 
					       rlist, zedlist, philist, arrayofEroverEzA  );
	  ephiOverEz(i,j) = Interpolate3DTable(order, r, z, phi, rows, columns, phiSlices,
					       rlist, zedlist, philist, arrayofEphioverEzA);
	  deltaEz(i,j)    = Interpolate3DTable(order, r, z, phi, rows, columns, phiSlices,
					       rlist, zedlist, philist, arrayofDeltaEzA   );
	} else {
	  erOverEz(i,j)   = Interpolate3DTable(order, r, z, phi, rows, columns, phiSlices, 
					       rlist, zedlist, philist, arrayofEroverEzC  );
	  ephiOverEz(i,j) = Interpolate3DTable(order, r, z, phi, rows, columns, phiSlices,
					       rlist, zedlist, philist, arrayofEphioverEzC);
	  deltaEz(i,j)  = - Interpolate3DTable(order, r, z, phi, rows, columns, phiSlices,
					       rlist, zedlist, philist, arrayofDeltaEzC   );
	  // negative coordinate system on C side
	}

      } // end r loop
    } // end z loop
  } // end phi loop

 
  // clear the temporary arrays lists
  for ( Int_t k = 0 ; k < phiSlices ; k++ )  {

    delete arrayofErA[k];  
    delete arrayofEphiA[k];
    delete arrayofdEzA[k];
    delete arrayofErC[k];  
    delete arrayofEphiC[k];
    delete arrayofdEzC[k];

    delete arrayofEroverEzA[k];  
    delete arrayofEphioverEzA[k];
    delete arrayofDeltaEzA[k];
    delete arrayofEroverEzC[k];  
    delete arrayofEphioverEzC[k];
    delete arrayofDeltaEzC[k];

  }

  fInitLookUp = kTRUE;

}


void AliTPCSpaceCharge3D::SetSCDataFileName(TString fname) {
  //
  // Set & load the Space charge density distribution from a file 
  // (linear interpolation onto a standard grid)
  //

 
  fSCDataFileName = fname;

  TFile *f = new TFile(fSCDataFileName.Data(),"READ");
  if (!f) { 
    AliError(Form("File %s, which should contain the space charge distribution, could not be found",
		  fSCDataFileName.Data()));
    return;
  }
 
  TH2F *densityRZ = (TH2F*) f->Get("SpaceChargeInRZ");
  if (!densityRZ) { 
    AliError(Form("The indicated file (%s) does not contain a histogram called %s",
		  fSCDataFileName.Data(),"SpaceChargeInRZ"));
    return;
  }

  TH3F *densityRPhi = (TH3F*) f->Get("SpaceChargeInRPhi");
  if (!densityRPhi) { 
    AliError(Form("The indicated file (%s) does not contain a histogram called %s",
		  fSCDataFileName.Data(),"SpaceChargeInRPhi"));
    return;
  }
 

  Double_t  r, phi, z ;

  TMatrixD &scDensityInRZ   =  *fSCdensityInRZ;
  TMatrixD &scDensityInRPhiA   =  *fSCdensityInRPhiA;
  TMatrixD &scDensityInRPhiC   =  *fSCdensityInRPhiC;
  for ( Int_t k = 0 ; k < kNPhi ; k++ ) {
    phi = fgkPhiList[k] ;
    TMatrixF &scDensity   =  *fSCdensityDistribution[k]  ;
    for ( Int_t j = 0 ; j < kNZ ; j++ ) {
      z = fgkZList[j] ; 
      for ( Int_t i = 0 ; i < kNR ; i++ ) { 
	r = fgkRList[i] ;

	// partial load in (r,z)
	if (k==0) // do just once
	  scDensityInRZ(i,j) =  densityRZ->Interpolate(r,z);

	// partial load in (r,phi)
	if ( j==0 || j == kNZ/2 ) {
	  if (z>0) 
	    scDensityInRPhiA(i,k) =  densityRPhi->Interpolate(r,phi,0.499);  // A side
	  else 
	    scDensityInRPhiC(i,k) =  densityRPhi->Interpolate(r,phi,-0.499); // C side
	}

	// Full 3D configuration ...
	if (z>0) 
	   scDensity(i,j) = scDensityInRZ(i,j) + scDensityInRPhiA(i,k); 
	else
	   scDensity(i,j) = scDensityInRZ(i,j) + scDensityInRPhiC(i,k); 
      }
    }
  }

  f->Close();

  fInitLookUp = kFALSE;

  
}

void     AliTPCSpaceCharge3D::SetInputSpaceCharge(TH3 * hisSpaceCharge3D,  TH2 * hisRPhi, TH2* hisRZ, Double_t norm){
  //
  // Use 3D space charge map as an optional input
  // The layout of the input histogram is assumed to be: (phi,r,z)
  // Density histogram is expreseed is expected to bin in  C/m^3
  //
  // Standard histogram interpolation is used in order to use the density at center of voxel
  // 
  //
  fSpaceChargeHistogram3D = hisSpaceCharge3D;
  fSpaceChargeHistogramRPhi = hisRPhi;
  fSpaceChargeHistogramRZ = hisRZ;

  Double_t  r, phi, z ;
  TMatrixD &scDensityInRZ   =  *fSCdensityInRZ;
  TMatrixD &scDensityInRPhiA   =  *fSCdensityInRPhiA;
  TMatrixD &scDensityInRPhiC   =  *fSCdensityInRPhiC;
  //
  Double_t rmin=hisSpaceCharge3D->GetYaxis()->GetBinCenter(0);
  Double_t rmax=hisSpaceCharge3D->GetYaxis()->GetBinUpEdge(hisSpaceCharge3D->GetYaxis()->GetNbins());
  Double_t zmin=hisSpaceCharge3D->GetZaxis()->GetBinCenter(0);
  Double_t zmax=hisSpaceCharge3D->GetZaxis()->GetBinCenter(hisSpaceCharge3D->GetZaxis()->GetNbins());

  for ( Int_t k = 0 ; k < kNPhi ; k++ ) {
    phi = fgkPhiList[k] ;
    TMatrixF &scDensity   =  *fSCdensityDistribution[k]  ;
    for ( Int_t j = 0 ; j < kNZ ; j++ ) {
      z = fgkZList[j] ; 
      for ( Int_t i = 0 ; i < kNR ; i++ ) { 
	// Full 3D configuration ...
	r = fgkRList[i] ;
	if (r>rmin && r<rmax && z>zmin && z< zmax){	  	  
	  // partial load in (r,z)
	  if (k==0) {
	    if (fSpaceChargeHistogramRZ) scDensityInRZ(i,j) = norm*fSpaceChargeHistogramRZ->Interpolate(r,z) ;
	  }
	  // partial load in (r,phi)
	  if ( (j==0 || j == kNZ/2) && fSpaceChargeHistogramRPhi) {
	    if (z>0) 
	      scDensityInRPhiA(i,k) = norm*fSpaceChargeHistogramRPhi->Interpolate(phi,r);  // A side
	    else 
	      scDensityInRPhiC(i,k) = norm*fSpaceChargeHistogramRPhi->Interpolate(phi+TMath::TwoPi(),r); // C side
	  }
	  
	  if (z>0) 
	    scDensity(i,j) = norm*fSpaceChargeHistogram3D->Interpolate(phi,r,z); 
	  else
	    scDensity(i,j) = norm*fSpaceChargeHistogram3D->Interpolate(phi,r,z);
	}
      }
    }
  }
  
  fInitLookUp = kFALSE;

}


Float_t  AliTPCSpaceCharge3D::GetSpaceChargeDensity(Float_t r, Float_t phi, Float_t z, Int_t mode) {
  //
  // returns the (input) space charge density at a given point according 
  // Note: input in [cm], output in [C/m^3/e0] !!
  //

  while (phi<0) phi += TMath::TwoPi();
  while (phi>TMath::TwoPi()) phi -= TMath::TwoPi();


  // Float_t sc =fSCdensityDistribution->Interpolate(r0,phi0,z0);
  Int_t order = 1; //
  Float_t sc = 0;

  if (mode == 0) { // return full load
    sc = Interpolate3DTable(order, r, z, phi, kNR, kNZ, kNPhi, 
			    fgkRList, fgkZList, fgkPhiList, fSCdensityDistribution );
 
  } else if (mode == 1) { // return partial load in (r,z)
    TMatrixD &scDensityInRZ   =  *fSCdensityInRZ;
    sc = Interpolate2DTable(order, r, z, kNR, kNZ, fgkRList, fgkZList, scDensityInRZ );
    
  } else if (mode == 2) { // return partial load in (r,phi)

    if (z>0) {
      TMatrixD &scDensityInRPhi   =  *fSCdensityInRPhiA;
      sc = Interpolate2DTable(order, r, phi, kNR, kNPhi, fgkRList, fgkPhiList, scDensityInRPhi );
    } else {
      TMatrixD &scDensityInRPhi   =  *fSCdensityInRPhiC;
      sc = Interpolate2DTable(order, r, phi, kNR, kNPhi, fgkRList, fgkPhiList, scDensityInRPhi );
    }

  } else {
    // should i give a warning?
    sc = 0;
  }
  
  //  printf("%f %f %f: %f\n",r,phi,z,sc);
  
  return sc;
}


TH2F * AliTPCSpaceCharge3D::CreateHistoSCinXY(Float_t z, Int_t nx, Int_t ny, Int_t mode) {
  //
  // return a simple histogramm containing the space charge distribution (input for the calculation)
  //

  TH2F *h=CreateTH2F("spaceCharge",GetTitle(),"x [cm]","y [cm]","#rho_{sc} [C/m^{3}/e_{0}]",
		     nx,-250.,250.,ny,-250.,250.);

  for (Int_t iy=1;iy<=ny;++iy) {
    Double_t yp = h->GetYaxis()->GetBinCenter(iy);
    for (Int_t ix=1;ix<=nx;++ix) {
      Double_t xp = h->GetXaxis()->GetBinCenter(ix);
    
      Float_t r = TMath::Sqrt(xp*xp+yp*yp);
      Float_t phi = TMath::ATan2(yp,xp); 
   
      if (85.<=r && r<=250.) {
	Float_t sc = GetSpaceChargeDensity(r,phi,z,mode)/fgke0; // in [C/m^3/e0]
	h->SetBinContent(ix,iy,sc); 
      } else {
	h->SetBinContent(ix,iy,0.);
      }
    }
  }
  
  return h;
} 

TH2F * AliTPCSpaceCharge3D::CreateHistoSCinZR(Float_t phi, Int_t nz, Int_t nr,Int_t mode ) {
  //
  // return a simple histogramm containing the space charge distribution (input for the calculation)
  //

  TH2F *h=CreateTH2F("spaceCharge",GetTitle(),"z [cm]","r [cm]","#rho_{sc} [C/m^{3}/e_{0}]",
		     nz,-250.,250.,nr,85.,250.);

  for (Int_t ir=1;ir<=nr;++ir) {
    Float_t r = h->GetYaxis()->GetBinCenter(ir);
    for (Int_t iz=1;iz<=nz;++iz) {
      Float_t z = h->GetXaxis()->GetBinCenter(iz);
      Float_t sc = GetSpaceChargeDensity(r,phi,z,mode)/fgke0; // in [C/m^3/e0]
      h->SetBinContent(iz,ir,sc);
    }
  }

  return h;
} 

void AliTPCSpaceCharge3D::WriteChargeDistributionToFile(const char* fname) {
  //
  // Example on how to write a Space charge distribution into a File
  //  (see below: estimate from scaling STAR measurements to Alice)
  // Charge distribution is splitted into two (RZ and RPHI) in order to speed up
  // the needed calculation time
  //

  TFile *f = new TFile(fname,"RECREATE");
  
  // some grid, not too course
  Int_t nr = 350;
  Int_t nphi = 180;
  Int_t nz = 500;

  Double_t dr = (fgkOFCRadius-fgkIFCRadius)/(nr+1);
  Double_t dphi = TMath::TwoPi()/(nphi+1);
  Double_t dz = 500./(nz+1);
  Double_t safty = 0.; // due to a root bug which does not interpolate the boundary (first and last bin) correctly


  // Charge distribution in ZR (rotational symmetric) ------------------

  TH2F *histoZR = new TH2F("chargeZR","chargeZR",
			   nr,fgkIFCRadius-dr-safty,fgkOFCRadius+dr+safty,
			   nz,-250-dz-safty,250+dz+safty);
 
  for (Int_t ir=1;ir<=nr;++ir) {
    Double_t rp = histoZR->GetXaxis()->GetBinCenter(ir);
    for (Int_t iz=1;iz<=nz;++iz) {
      Double_t zp = histoZR->GetYaxis()->GetBinCenter(iz);
      
      // recalculation to meter
      Double_t lZ = 2.5; // approx. TPC drift length
      Double_t rpM = rp/100.; // in [m]
      Double_t zpM = TMath::Abs(zp/100.); // in [m]
      
      // setting of mb multiplicity and Interaction rate
      Double_t multiplicity = 950;
      Double_t intRate = 7800;

      // calculation of "scaled" parameters
      Double_t a = multiplicity*intRate/79175;
      Double_t b = a/lZ;
	
      Double_t charge = (a - b*zpM)/(rpM*rpM); // charge in [C/m^3/e0]
	
      charge = charge*fgke0; // [C/m^3]
	
      if (zp<0) charge *= 0.9; // e.g. slightly less on C side due to front absorber

      //  charge = 0; // for tests
      histoZR->SetBinContent(ir,iz,charge); 
    }
  }
    
  histoZR->Write("SpaceChargeInRZ");
  

  // Charge distribution in RPhi (e.g. Floating GG wire) ------------
  
  TH3F *histoRPhi = new TH3F("chargeRPhi","chargeRPhi",
			     nr,fgkIFCRadius-dr-safty,fgkOFCRadius+dr+safty,
			     nphi,0-dphi-safty,TMath::TwoPi()+dphi+safty,
			     2,-1,1); // z part - to allow A and C side differences
  
  // some 'arbitrary' GG leaks
  Int_t   nGGleaks = 5;
  Double_t secPosA[5]    = {3,6,6,11,13};         // sector
  Double_t radialPosA[5] = {125,100,160,200,230}; // radius in cm
  Double_t secPosC[5]    = {1,8,12,15,15};        // sector
  Double_t radialPosC[5] = {245,120,140,120,190}; // radius in cm

  for (Int_t ir=1;ir<=nr;++ir) {
    Double_t rp = histoRPhi->GetXaxis()->GetBinCenter(ir);
    for (Int_t iphi=1;iphi<=nphi;++iphi) {
      Double_t phip = histoRPhi->GetYaxis()->GetBinCenter(iphi);
      for (Int_t iz=1;iz<=2;++iz) {
	Double_t zp = histoRPhi->GetZaxis()->GetBinCenter(iz);
	
	Double_t charge = 0;
	
	for (Int_t igg = 0; igg<nGGleaks; igg++) { // loop over GG leaks
	  
	  // A side
	  Double_t secPos = secPosA[igg]; 
	  Double_t radialPos = radialPosA[igg];

	  if (zp<0) { // C side
	    secPos = secPosC[igg]; 
	    radialPos = radialPosC[igg];
	  }	    

	  // some 'arbitrary' GG leaks
	  if (  (phip<(TMath::Pi()/9*(secPos+1)) && phip>(TMath::Pi()/9*secPos) ) ) { // sector slice
	    if ( rp>(radialPos-2.5) && rp<(radialPos+2.5))  // 5 cm slice
	      charge = 300; 
	  }
	  
	}		
	
	charge = charge*fgke0; // [C/m^3]

	histoRPhi->SetBinContent(ir,iphi,iz,charge); 
      }
    }
  }

  histoRPhi->Write("SpaceChargeInRPhi");

  f->Close();
  
}


void AliTPCSpaceCharge3D::Print(const Option_t* option) const {
  //
  // Print function to check the settings of the boundary vectors
  // option=="a" prints the C0 and C1 coefficents for calibration purposes
  //

  TString opt = option; opt.ToLower();
  printf("%s\n",GetTitle());
  printf(" - Space Charge effect with arbitrary 3D charge density (as input).\n");
  printf("   SC correction factor: %f \n",fCorrectionFactor);

  if (opt.Contains("a")) { // Print all details
    printf(" - T1: %1.4f, T2: %1.4f \n",fT1,fT2);
    printf(" - C1: %1.4f, C0: %1.4f \n",fC1,fC0);
  } 
   
  if (!fInitLookUp) AliError("Lookup table was not initialized! You should do InitSpaceCharge3DDistortion() ...");

} 



void AliTPCSpaceCharge3D::InitSpaceCharge3DPoisson(Int_t kRows, Int_t kColumns, Int_t kPhiSlices, Int_t kIterations){
  //
  // MI extension  - calculate E field
  //               - inspired by  AliTPCROCVoltError3D::InitROCVoltError3D()
  // Initialization of the Lookup table which contains the solutions of the 
  // Dirichlet boundary problem
  // Calculation of the single 3D-Poisson solver is done just if needed
  // (see basic lookup tables in header file)
  //  
  Int_t kPhiSlicesPerSector = kPhiSlices/18; 
  //
  const Int_t   order       =    1  ;  // Linear interpolation = 1, Quadratic = 2  
  const Float_t gridSizeR   =  (fgkOFCRadius-fgkIFCRadius) / (kRows-1) ;
  const Float_t gridSizeZ   =  fgkTPCZ0 / (kColumns-1) ;
  const Float_t gridSizePhi =  TMath::TwoPi() / ( 18.0 * kPhiSlicesPerSector);

  // temporary arrays to create the boundary conditions
  TMatrixD *arrayofArrayV[kPhiSlices], *arrayofCharge[kPhiSlices] ; 
  TMatrixD *arrayofEroverEz[kPhiSlices], *arrayofEphioverEz[kPhiSlices], *arrayofDeltaEz[kPhiSlices] ; 

  for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) {
    arrayofArrayV[k]     =   new TMatrixD(kRows,kColumns) ;
    arrayofCharge[k]     =   new TMatrixD(kRows,kColumns) ;
    arrayofEroverEz[k]   =   new TMatrixD(kRows,kColumns) ;
    arrayofEphioverEz[k] =   new TMatrixD(kRows,kColumns) ;
    arrayofDeltaEz[k]    =   new TMatrixD(kRows,kColumns) ;
  }
  
  // list of point as used in the poisson relation and the interpolation (during sum up)
  Double_t  rlist[kRows], zedlist[kColumns] , philist[kPhiSlices];
  for ( Int_t k = 0 ; k < kPhiSlices ; k++ ) {
    philist[k] =  gridSizePhi * k;
    for ( Int_t i = 0 ; i < kRows ; i++ )    {
      rlist[i] = fgkIFCRadius + i*gridSizeR ;
      for ( Int_t j = 0 ; j < kColumns ; j++ ) { // Fill Vmatrix with Boundary Conditions
	zedlist[j]  = j * gridSizeZ ;
      }
    }
  }

  // ==========================================================================
  // Solve Poisson's equation in 3D cylindrical coordinates by relaxation technique
  // Allow for different size grid spacing in R and Z directions
  
  const Int_t   symmetry = 0;
 
  // Set bondaries and solve Poisson's equation --------------------------
  
  if ( !fInitLookUp ) {
    
    AliInfo(Form("Solving the poisson equation (~ %d sec)",2*10*(int)(kPhiSlices/10)));
    
    for ( Int_t side = 0 ; side < 2 ; side++ ) {  // Solve Poisson3D twice; once for +Z and once for -Z
      AliSysInfo::AddStamp("RunSide", 1,side,0);
      for ( Int_t k = 0 ; k < kPhiSlices ; k++ )  {
	TMatrixD &arrayV    =  *arrayofArrayV[k] ;
	TMatrixD &charge    =  *arrayofCharge[k] ;
	
	//Fill arrays with initial conditions.  V on the boundary and Charge in the volume.
// 	for ( Int_t i = 0 ; i < kRows ; i++ ) {
// 	  for ( Int_t j = 0 ; j < kColumns ; j++ ) {  // Fill Vmatrix with Boundary Conditions
// 	    arrayV(i,j) = 0.0 ; 
// 	    charge(i,j) = 0.0 ;

// // 	    Float_t radius0 = rlist[i] ;
// // 	    Float_t phi0    = gridSizePhi * k ;
	    
// 	    // To avoid problems at sector boundaries, use an average of +- 1 degree from actual phi location
// // 	    if ( j == (kColumns-1) ) {
// // 	      arrayV(i,j) = 0.5*  ( GetROCVoltOffset( side, radius0, phi0+0.02 ) + GetROCVoltOffset( side, radius0, phi0-0.02 ) ) ;

// // 	      if (side==1) // C side
// // 		arrayV(i,j) = -arrayV(i,j); // minus sign on the C side to allow a consistent usage of global z when setting the boundaries
// // 	    }
// 	  }
// 	}      
	
	for ( Int_t i = 1 ; i < kRows-1 ; i++ ) { 
	  for ( Int_t j = 1 ; j < kColumns-1 ; j++ ) {	    
	    Float_t radius0 = rlist[i] ;
	    Float_t phi0    = gridSizePhi * k ;
	    Double_t z0 = zedlist[j];
	    if (side==1) z0= -TMath::Abs(zedlist[j]);
	    arrayV(i,j) = 0.0 ; 
	    charge(i,j)  =  fSpaceChargeHistogram3D->Interpolate(phi0,radius0,z0);
	  }
	}
      }      
      AliSysInfo::AddStamp("RunPoisson", 2,side,0);
      
      // Solve Poisson's equation in 3D cylindrical coordinates by relaxation technique
      // Allow for different size grid spacing in R and Z directions
      
      //      PoissonRelaxation3D( arrayofArrayV, arrayofCharge, 
      // 			   arrayofEroverEz, arrayofEphioverEz, arrayofDeltaEz,
      // 			   kRows, kColumns, kPhiSlices, gridSizePhi, kIterations, 
      // 			   symmetry , fROCdisplacement) ;
      PoissonRelaxation3D( arrayofArrayV, arrayofCharge, 
			   arrayofEroverEz, arrayofEphioverEz, arrayofDeltaEz,
			   kRows, kColumns, kPhiSlices, gridSizePhi, kIterations, 
			   symmetry ) ;
      
      //Interpolate results onto a custom grid which is used just for these calculations.
      Double_t  r, phi, z ;
      for ( Int_t k = 0 ; k < kNPhi ; k++ ) {
	phi = fgkPhiList[k] ;
	
	TMatrixF &erOverEz   =  *fLookUpErOverEz[k]  ;
	TMatrixF &ephiOverEz =  *fLookUpEphiOverEz[k];
	TMatrixF &deltaEz    =  *fLookUpDeltaEz[k]   ;
	
	for ( Int_t j = 0 ; j < kNZ ; j++ ) {

	  z = TMath::Abs(fgkZList[j]) ;  // Symmetric solution in Z that depends only on ABS(Z)
  
	  if ( side == 0 &&  fgkZList[j] < 0 ) continue; // Skip rest of this loop if on the wrong side
	  if ( side == 1 &&  fgkZList[j] > 0 ) continue; // Skip rest of this loop if on the wrong side
	  
	  for ( Int_t i = 0 ; i < kNR ; i++ ) { 
	    r = fgkRList[i] ;

	    // Interpolate basicLookup tables; once for each rod, then sum the results
	    erOverEz(i,j)   = Interpolate3DTable(order, r, z, phi, kRows, kColumns, kPhiSlices, 
						 rlist, zedlist, philist, arrayofEroverEz  );
	    ephiOverEz(i,j) = Interpolate3DTable(order, r, z, phi, kRows, kColumns, kPhiSlices,
						 rlist, zedlist, philist, arrayofEphioverEz);
	    deltaEz(i,j)    = Interpolate3DTable(order, r, z, phi, kRows, kColumns, kPhiSlices,
						 rlist, zedlist, philist, arrayofDeltaEz  );

	    if (side == 1)  deltaEz(i,j) = -  deltaEz(i,j); // negative coordinate system on C side

	  } // end r loop
	}// end z loop
      }// end phi loop
      AliSysInfo::AddStamp("Interpolate Poisson", 3,side,0);      
      if ( side == 0 ) AliInfo(" A side done");
      if ( side == 1 ) AliInfo(" C side done");
    } // end side loop
  }
  
  // clear the temporary arrays lists
  for ( Int_t k = 0 ; k < kPhiSlices ; k++ )  {
    delete arrayofArrayV[k];
    delete arrayofCharge[k];
    delete arrayofEroverEz[k];  
    delete arrayofEphioverEz[k];
    delete arrayofDeltaEz[k];
  }
 

  fInitLookUp = kTRUE;

}



AliTPCSpaceCharge3D * AliTPCSpaceCharge3D::MakeCorrection22D(const char* fileName, Double_t multiplicity, Double_t intRate, Double_t epsIROC, Double_t epsOROC,
							     Double_t gasfactor, 
							     Double_t radialScaling){
  //
  // Origin: Christian Lippmann, CERN, Christian.Lippmann@cern.ch based on the internal note (xxx ...)
  // adopted by Marian Ivanov (different epsilon in IROC and OROC)
  //
  // Charge distribution is splitted into two (RZ and RPHI) in order to speed up
  // the needed calculation time.  
  //
  // Explanation of variables:
  // 1) multiplicity: charghed particle dn/deta for top 80% centrality (660 for 2011,
  //    expect 950 for full energy)
  // 2) intRate: Total interaction rate (e.g. 50kHz for the upgrade)
  // 3) eps: Number of backdrifting ions per primary electron (0 for MWPC, e.g.10 for GEM)
  // 4) gasfactor: Use different gas. E.g. Ar/CO2 has twice the primary ionization, ion drift
  //    velocity factor 2.5 slower, so  gasfactor = 5.
  //


  TFile *f = new TFile(fileName, "RECREATE");  
  // some grid, not too coarse
  const Int_t nr   = 350;
  const Int_t nphi = 180;
  const Int_t nz   = 500;
  const Double_t kROROC=134;  // hardwired OROC radius 

  Double_t dr = (fgkOFCRadius-fgkIFCRadius)/(nr+1);
  Double_t dphi = TMath::TwoPi()/(nphi+1);
  Double_t dz = 500./(nz+1);
  Double_t safty = 0.; // due to a root bug which does not interpolate the boundary ..
  // .. (first and last bin) correctly

  // Charge distribution in ZR (rotational symmetric) ------------------

  TH2F *histoZR = new TH2F("chargeZR", "chargeZR",
                           nr, fgkIFCRadius-dr-safty, fgkOFCRadius+dr+safty,
                           nz, -250-dz-safty, 250+dz+safty);

  // For the normalization to same integral as radial exponent = 2
  Double_t radialExponent             = -2.; // reference = 2
  Double_t radiusInner                = histoZR->GetXaxis()->GetBinCenter(1) / 100.;//in [m]
  Double_t radiusOuter                = histoZR->GetXaxis()->GetBinCenter(nr) / 100.;//in [m]
  Double_t integralRadialExponent2    = TMath::Power(radiusOuter,radialExponent+1) * 1./(radialExponent+1) 
    - TMath::Power(radiusInner,radialExponent+1) * 1./(radialExponent+1);
  
  radialExponent                      = -radialScaling; // user set   
  Double_t integralRadialExponentUser = 0.;
  if(radialScaling > 1 + 0.000001 || radialScaling < 1 - 0.000001 ) // to avoid n = -1
    integralRadialExponentUser = TMath::Power(radiusOuter,radialExponent+1) * 1./(radialExponent+1) 
      - TMath::Power(radiusInner,radialExponent+1) * 1./(radialExponent+1);
  else
    integralRadialExponentUser = TMath::Log(radiusOuter) - TMath::Log(radiusInner);
    
  Double_t normRadialExponent         = integralRadialExponent2 / integralRadialExponentUser;
 
  for (Int_t ir=1;ir<=nr;++ir) {
    Double_t rp = histoZR->GetXaxis()->GetBinCenter(ir);
    for (Int_t iz=1;iz<=nz;++iz) {
      Double_t zp = histoZR->GetYaxis()->GetBinCenter(iz);
      Double_t eps = (rp <kROROC) ? epsIROC:epsOROC; 
      // recalculation to meter
      Double_t lZ = 2.5; // approx. TPC drift length
      Double_t rpM = rp/100.; // in [m]
      Double_t zpM = TMath::Abs(zp/100.); // in [m]
 
      // calculation of "scaled" parameters
      Double_t a = multiplicity*intRate/76628;
      //Double_t charge = gasfactor * ( a / (rpM*rpM) * (1 - zpM/lZ) ); // charge in [C/m^3/e0], no IBF
      Double_t charge = normRadialExponent * gasfactor * ( a / (TMath::Power(rpM,radialScaling)) * (1 - zpM/lZ + eps) ); // charge in [C/m^3/e0], with IBF
      
      charge = charge*fgke0;          // [C/m^3]
      if (zp<0) charge *= 0.9; // Slightly less on C side due to front absorber

      histoZR->SetBinContent(ir, iz, charge); 
    }
  }
    
  histoZR->Write("SpaceChargeInRZ");
  //
  // Charge distribution in RPhi (e.g. Floating GG wire) ------------
  //
  TH3F *histoRPhi = new TH3F("chargeRPhi", "chargeRPhi",
                             nr, fgkIFCRadius-dr-safty, fgkOFCRadius+dr+safty,
                             nphi, 0-dphi-safty, TMath::TwoPi()+dphi+safty,
                             2, -1, 1); // z part - to allow A and C side differences  
  histoRPhi->Write("SpaceChargeInRPhi");
  f->Close();
  //
  //
  //
  AliTPCSpaceCharge3D *spaceCharge = new AliTPCSpaceCharge3D;
  spaceCharge->SetSCDataFileName(fileName);
  spaceCharge->SetOmegaTauT1T2(0.325,1,1); // Ne CO2
  spaceCharge->InitSpaceCharge3DDistortion();
  spaceCharge->AddVisualCorrection(spaceCharge,1);
  //
  f = new TFile(fileName, "update");
  spaceCharge->Write("spaceCharge");
  f->Close();
  return spaceCharge;
}

 AliTPCSpaceCharge3D.cxx:1
 AliTPCSpaceCharge3D.cxx:2
 AliTPCSpaceCharge3D.cxx:3
 AliTPCSpaceCharge3D.cxx:4
 AliTPCSpaceCharge3D.cxx:5
 AliTPCSpaceCharge3D.cxx:6
 AliTPCSpaceCharge3D.cxx:7
 AliTPCSpaceCharge3D.cxx:8
 AliTPCSpaceCharge3D.cxx:9
 AliTPCSpaceCharge3D.cxx:10
 AliTPCSpaceCharge3D.cxx:11
 AliTPCSpaceCharge3D.cxx:12
 AliTPCSpaceCharge3D.cxx:13
 AliTPCSpaceCharge3D.cxx:14
 AliTPCSpaceCharge3D.cxx:15
 AliTPCSpaceCharge3D.cxx:16
 AliTPCSpaceCharge3D.cxx:17
 AliTPCSpaceCharge3D.cxx:18
 AliTPCSpaceCharge3D.cxx:19
 AliTPCSpaceCharge3D.cxx:20
 AliTPCSpaceCharge3D.cxx:21
 AliTPCSpaceCharge3D.cxx:22
 AliTPCSpaceCharge3D.cxx:23
 AliTPCSpaceCharge3D.cxx:24
 AliTPCSpaceCharge3D.cxx:25
 AliTPCSpaceCharge3D.cxx:26
 AliTPCSpaceCharge3D.cxx:27
 AliTPCSpaceCharge3D.cxx:28
 AliTPCSpaceCharge3D.cxx:29
 AliTPCSpaceCharge3D.cxx:30
 AliTPCSpaceCharge3D.cxx:31
 AliTPCSpaceCharge3D.cxx:32
 AliTPCSpaceCharge3D.cxx:33
 AliTPCSpaceCharge3D.cxx:34
 AliTPCSpaceCharge3D.cxx:35
 AliTPCSpaceCharge3D.cxx:36
 AliTPCSpaceCharge3D.cxx:37
 AliTPCSpaceCharge3D.cxx:38
 AliTPCSpaceCharge3D.cxx:39
 AliTPCSpaceCharge3D.cxx:40
 AliTPCSpaceCharge3D.cxx:41
 AliTPCSpaceCharge3D.cxx:42
 AliTPCSpaceCharge3D.cxx:43
 AliTPCSpaceCharge3D.cxx:44
 AliTPCSpaceCharge3D.cxx:45
 AliTPCSpaceCharge3D.cxx:46
 AliTPCSpaceCharge3D.cxx:47
 AliTPCSpaceCharge3D.cxx:48
 AliTPCSpaceCharge3D.cxx:49
 AliTPCSpaceCharge3D.cxx:50
 AliTPCSpaceCharge3D.cxx:51
 AliTPCSpaceCharge3D.cxx:52
 AliTPCSpaceCharge3D.cxx:53
 AliTPCSpaceCharge3D.cxx:54
 AliTPCSpaceCharge3D.cxx:55
 AliTPCSpaceCharge3D.cxx:56
 AliTPCSpaceCharge3D.cxx:57
 AliTPCSpaceCharge3D.cxx:58
 AliTPCSpaceCharge3D.cxx:59
 AliTPCSpaceCharge3D.cxx:60
 AliTPCSpaceCharge3D.cxx:61
 AliTPCSpaceCharge3D.cxx:62
 AliTPCSpaceCharge3D.cxx:63
 AliTPCSpaceCharge3D.cxx:64
 AliTPCSpaceCharge3D.cxx:65
 AliTPCSpaceCharge3D.cxx:66
 AliTPCSpaceCharge3D.cxx:67
 AliTPCSpaceCharge3D.cxx:68
 AliTPCSpaceCharge3D.cxx:69
 AliTPCSpaceCharge3D.cxx:70
 AliTPCSpaceCharge3D.cxx:71
 AliTPCSpaceCharge3D.cxx:72
 AliTPCSpaceCharge3D.cxx:73
 AliTPCSpaceCharge3D.cxx:74
 AliTPCSpaceCharge3D.cxx:75
 AliTPCSpaceCharge3D.cxx:76
 AliTPCSpaceCharge3D.cxx:77
 AliTPCSpaceCharge3D.cxx:78
 AliTPCSpaceCharge3D.cxx:79
 AliTPCSpaceCharge3D.cxx:80
 AliTPCSpaceCharge3D.cxx:81
 AliTPCSpaceCharge3D.cxx:82
 AliTPCSpaceCharge3D.cxx:83
 AliTPCSpaceCharge3D.cxx:84
 AliTPCSpaceCharge3D.cxx:85
 AliTPCSpaceCharge3D.cxx:86
 AliTPCSpaceCharge3D.cxx:87
 AliTPCSpaceCharge3D.cxx:88
 AliTPCSpaceCharge3D.cxx:89
 AliTPCSpaceCharge3D.cxx:90
 AliTPCSpaceCharge3D.cxx:91
 AliTPCSpaceCharge3D.cxx:92
 AliTPCSpaceCharge3D.cxx:93
 AliTPCSpaceCharge3D.cxx:94
 AliTPCSpaceCharge3D.cxx:95
 AliTPCSpaceCharge3D.cxx:96
 AliTPCSpaceCharge3D.cxx:97
 AliTPCSpaceCharge3D.cxx:98
 AliTPCSpaceCharge3D.cxx:99
 AliTPCSpaceCharge3D.cxx:100
 AliTPCSpaceCharge3D.cxx:101
 AliTPCSpaceCharge3D.cxx:102
 AliTPCSpaceCharge3D.cxx:103
 AliTPCSpaceCharge3D.cxx:104
 AliTPCSpaceCharge3D.cxx:105
 AliTPCSpaceCharge3D.cxx:106
 AliTPCSpaceCharge3D.cxx:107
 AliTPCSpaceCharge3D.cxx:108
 AliTPCSpaceCharge3D.cxx:109
 AliTPCSpaceCharge3D.cxx:110
 AliTPCSpaceCharge3D.cxx:111
 AliTPCSpaceCharge3D.cxx:112
 AliTPCSpaceCharge3D.cxx:113
 AliTPCSpaceCharge3D.cxx:114
 AliTPCSpaceCharge3D.cxx:115
 AliTPCSpaceCharge3D.cxx:116
 AliTPCSpaceCharge3D.cxx:117
 AliTPCSpaceCharge3D.cxx:118
 AliTPCSpaceCharge3D.cxx:119
 AliTPCSpaceCharge3D.cxx:120
 AliTPCSpaceCharge3D.cxx:121
 AliTPCSpaceCharge3D.cxx:122
 AliTPCSpaceCharge3D.cxx:123
 AliTPCSpaceCharge3D.cxx:124
 AliTPCSpaceCharge3D.cxx:125
 AliTPCSpaceCharge3D.cxx:126
 AliTPCSpaceCharge3D.cxx:127
 AliTPCSpaceCharge3D.cxx:128
 AliTPCSpaceCharge3D.cxx:129
 AliTPCSpaceCharge3D.cxx:130
 AliTPCSpaceCharge3D.cxx:131
 AliTPCSpaceCharge3D.cxx:132
 AliTPCSpaceCharge3D.cxx:133
 AliTPCSpaceCharge3D.cxx:134
 AliTPCSpaceCharge3D.cxx:135
 AliTPCSpaceCharge3D.cxx:136
 AliTPCSpaceCharge3D.cxx:137
 AliTPCSpaceCharge3D.cxx:138
 AliTPCSpaceCharge3D.cxx:139
 AliTPCSpaceCharge3D.cxx:140
 AliTPCSpaceCharge3D.cxx:141
 AliTPCSpaceCharge3D.cxx:142
 AliTPCSpaceCharge3D.cxx:143
 AliTPCSpaceCharge3D.cxx:144
 AliTPCSpaceCharge3D.cxx:145
 AliTPCSpaceCharge3D.cxx:146
 AliTPCSpaceCharge3D.cxx:147
 AliTPCSpaceCharge3D.cxx:148
 AliTPCSpaceCharge3D.cxx:149
 AliTPCSpaceCharge3D.cxx:150
 AliTPCSpaceCharge3D.cxx:151
 AliTPCSpaceCharge3D.cxx:152
 AliTPCSpaceCharge3D.cxx:153
 AliTPCSpaceCharge3D.cxx:154
 AliTPCSpaceCharge3D.cxx:155
 AliTPCSpaceCharge3D.cxx:156
 AliTPCSpaceCharge3D.cxx:157
 AliTPCSpaceCharge3D.cxx:158
 AliTPCSpaceCharge3D.cxx:159
 AliTPCSpaceCharge3D.cxx:160
 AliTPCSpaceCharge3D.cxx:161
 AliTPCSpaceCharge3D.cxx:162
 AliTPCSpaceCharge3D.cxx:163
 AliTPCSpaceCharge3D.cxx:164
 AliTPCSpaceCharge3D.cxx:165
 AliTPCSpaceCharge3D.cxx:166
 AliTPCSpaceCharge3D.cxx:167
 AliTPCSpaceCharge3D.cxx:168
 AliTPCSpaceCharge3D.cxx:169
 AliTPCSpaceCharge3D.cxx:170
 AliTPCSpaceCharge3D.cxx:171
 AliTPCSpaceCharge3D.cxx:172
 AliTPCSpaceCharge3D.cxx:173
 AliTPCSpaceCharge3D.cxx:174
 AliTPCSpaceCharge3D.cxx:175
 AliTPCSpaceCharge3D.cxx:176
 AliTPCSpaceCharge3D.cxx:177
 AliTPCSpaceCharge3D.cxx:178
 AliTPCSpaceCharge3D.cxx:179
 AliTPCSpaceCharge3D.cxx:180
 AliTPCSpaceCharge3D.cxx:181
 AliTPCSpaceCharge3D.cxx:182
 AliTPCSpaceCharge3D.cxx:183
 AliTPCSpaceCharge3D.cxx:184
 AliTPCSpaceCharge3D.cxx:185
 AliTPCSpaceCharge3D.cxx:186
 AliTPCSpaceCharge3D.cxx:187
 AliTPCSpaceCharge3D.cxx:188
 AliTPCSpaceCharge3D.cxx:189
 AliTPCSpaceCharge3D.cxx:190
 AliTPCSpaceCharge3D.cxx:191
 AliTPCSpaceCharge3D.cxx:192
 AliTPCSpaceCharge3D.cxx:193
 AliTPCSpaceCharge3D.cxx:194
 AliTPCSpaceCharge3D.cxx:195
 AliTPCSpaceCharge3D.cxx:196
 AliTPCSpaceCharge3D.cxx:197
 AliTPCSpaceCharge3D.cxx:198
 AliTPCSpaceCharge3D.cxx:199
 AliTPCSpaceCharge3D.cxx:200
 AliTPCSpaceCharge3D.cxx:201
 AliTPCSpaceCharge3D.cxx:202
 AliTPCSpaceCharge3D.cxx:203
 AliTPCSpaceCharge3D.cxx:204
 AliTPCSpaceCharge3D.cxx:205
 AliTPCSpaceCharge3D.cxx:206
 AliTPCSpaceCharge3D.cxx:207
 AliTPCSpaceCharge3D.cxx:208
 AliTPCSpaceCharge3D.cxx:209
 AliTPCSpaceCharge3D.cxx:210
 AliTPCSpaceCharge3D.cxx:211
 AliTPCSpaceCharge3D.cxx:212
 AliTPCSpaceCharge3D.cxx:213
 AliTPCSpaceCharge3D.cxx:214
 AliTPCSpaceCharge3D.cxx:215
 AliTPCSpaceCharge3D.cxx:216
 AliTPCSpaceCharge3D.cxx:217
 AliTPCSpaceCharge3D.cxx:218
 AliTPCSpaceCharge3D.cxx:219
 AliTPCSpaceCharge3D.cxx:220
 AliTPCSpaceCharge3D.cxx:221
 AliTPCSpaceCharge3D.cxx:222
 AliTPCSpaceCharge3D.cxx:223
 AliTPCSpaceCharge3D.cxx:224
 AliTPCSpaceCharge3D.cxx:225
 AliTPCSpaceCharge3D.cxx:226
 AliTPCSpaceCharge3D.cxx:227
 AliTPCSpaceCharge3D.cxx:228
 AliTPCSpaceCharge3D.cxx:229
 AliTPCSpaceCharge3D.cxx:230
 AliTPCSpaceCharge3D.cxx:231
 AliTPCSpaceCharge3D.cxx:232
 AliTPCSpaceCharge3D.cxx:233
 AliTPCSpaceCharge3D.cxx:234
 AliTPCSpaceCharge3D.cxx:235
 AliTPCSpaceCharge3D.cxx:236
 AliTPCSpaceCharge3D.cxx:237
 AliTPCSpaceCharge3D.cxx:238
 AliTPCSpaceCharge3D.cxx:239
 AliTPCSpaceCharge3D.cxx:240
 AliTPCSpaceCharge3D.cxx:241
 AliTPCSpaceCharge3D.cxx:242
 AliTPCSpaceCharge3D.cxx:243
 AliTPCSpaceCharge3D.cxx:244
 AliTPCSpaceCharge3D.cxx:245
 AliTPCSpaceCharge3D.cxx:246
 AliTPCSpaceCharge3D.cxx:247
 AliTPCSpaceCharge3D.cxx:248
 AliTPCSpaceCharge3D.cxx:249
 AliTPCSpaceCharge3D.cxx:250
 AliTPCSpaceCharge3D.cxx:251
 AliTPCSpaceCharge3D.cxx:252
 AliTPCSpaceCharge3D.cxx:253
 AliTPCSpaceCharge3D.cxx:254
 AliTPCSpaceCharge3D.cxx:255
 AliTPCSpaceCharge3D.cxx:256
 AliTPCSpaceCharge3D.cxx:257
 AliTPCSpaceCharge3D.cxx:258
 AliTPCSpaceCharge3D.cxx:259
 AliTPCSpaceCharge3D.cxx:260
 AliTPCSpaceCharge3D.cxx:261
 AliTPCSpaceCharge3D.cxx:262
 AliTPCSpaceCharge3D.cxx:263
 AliTPCSpaceCharge3D.cxx:264
 AliTPCSpaceCharge3D.cxx:265
 AliTPCSpaceCharge3D.cxx:266
 AliTPCSpaceCharge3D.cxx:267
 AliTPCSpaceCharge3D.cxx:268
 AliTPCSpaceCharge3D.cxx:269
 AliTPCSpaceCharge3D.cxx:270
 AliTPCSpaceCharge3D.cxx:271
 AliTPCSpaceCharge3D.cxx:272
 AliTPCSpaceCharge3D.cxx:273
 AliTPCSpaceCharge3D.cxx:274
 AliTPCSpaceCharge3D.cxx:275
 AliTPCSpaceCharge3D.cxx:276
 AliTPCSpaceCharge3D.cxx:277
 AliTPCSpaceCharge3D.cxx:278
 AliTPCSpaceCharge3D.cxx:279
 AliTPCSpaceCharge3D.cxx:280
 AliTPCSpaceCharge3D.cxx:281
 AliTPCSpaceCharge3D.cxx:282
 AliTPCSpaceCharge3D.cxx:283
 AliTPCSpaceCharge3D.cxx:284
 AliTPCSpaceCharge3D.cxx:285
 AliTPCSpaceCharge3D.cxx:286
 AliTPCSpaceCharge3D.cxx:287
 AliTPCSpaceCharge3D.cxx:288
 AliTPCSpaceCharge3D.cxx:289
 AliTPCSpaceCharge3D.cxx:290
 AliTPCSpaceCharge3D.cxx:291
 AliTPCSpaceCharge3D.cxx:292
 AliTPCSpaceCharge3D.cxx:293
 AliTPCSpaceCharge3D.cxx:294
 AliTPCSpaceCharge3D.cxx:295
 AliTPCSpaceCharge3D.cxx:296
 AliTPCSpaceCharge3D.cxx:297
 AliTPCSpaceCharge3D.cxx:298
 AliTPCSpaceCharge3D.cxx:299
 AliTPCSpaceCharge3D.cxx:300
 AliTPCSpaceCharge3D.cxx:301
 AliTPCSpaceCharge3D.cxx:302
 AliTPCSpaceCharge3D.cxx:303
 AliTPCSpaceCharge3D.cxx:304
 AliTPCSpaceCharge3D.cxx:305
 AliTPCSpaceCharge3D.cxx:306
 AliTPCSpaceCharge3D.cxx:307
 AliTPCSpaceCharge3D.cxx:308
 AliTPCSpaceCharge3D.cxx:309
 AliTPCSpaceCharge3D.cxx:310
 AliTPCSpaceCharge3D.cxx:311
 AliTPCSpaceCharge3D.cxx:312
 AliTPCSpaceCharge3D.cxx:313
 AliTPCSpaceCharge3D.cxx:314
 AliTPCSpaceCharge3D.cxx:315
 AliTPCSpaceCharge3D.cxx:316
 AliTPCSpaceCharge3D.cxx:317
 AliTPCSpaceCharge3D.cxx:318
 AliTPCSpaceCharge3D.cxx:319
 AliTPCSpaceCharge3D.cxx:320
 AliTPCSpaceCharge3D.cxx:321
 AliTPCSpaceCharge3D.cxx:322
 AliTPCSpaceCharge3D.cxx:323
 AliTPCSpaceCharge3D.cxx:324
 AliTPCSpaceCharge3D.cxx:325
 AliTPCSpaceCharge3D.cxx:326
 AliTPCSpaceCharge3D.cxx:327
 AliTPCSpaceCharge3D.cxx:328
 AliTPCSpaceCharge3D.cxx:329
 AliTPCSpaceCharge3D.cxx:330
 AliTPCSpaceCharge3D.cxx:331
 AliTPCSpaceCharge3D.cxx:332
 AliTPCSpaceCharge3D.cxx:333
 AliTPCSpaceCharge3D.cxx:334
 AliTPCSpaceCharge3D.cxx:335
 AliTPCSpaceCharge3D.cxx:336
 AliTPCSpaceCharge3D.cxx:337
 AliTPCSpaceCharge3D.cxx:338
 AliTPCSpaceCharge3D.cxx:339
 AliTPCSpaceCharge3D.cxx:340
 AliTPCSpaceCharge3D.cxx:341
 AliTPCSpaceCharge3D.cxx:342
 AliTPCSpaceCharge3D.cxx:343
 AliTPCSpaceCharge3D.cxx:344
 AliTPCSpaceCharge3D.cxx:345
 AliTPCSpaceCharge3D.cxx:346
 AliTPCSpaceCharge3D.cxx:347
 AliTPCSpaceCharge3D.cxx:348
 AliTPCSpaceCharge3D.cxx:349
 AliTPCSpaceCharge3D.cxx:350
 AliTPCSpaceCharge3D.cxx:351
 AliTPCSpaceCharge3D.cxx:352
 AliTPCSpaceCharge3D.cxx:353
 AliTPCSpaceCharge3D.cxx:354
 AliTPCSpaceCharge3D.cxx:355
 AliTPCSpaceCharge3D.cxx:356
 AliTPCSpaceCharge3D.cxx:357
 AliTPCSpaceCharge3D.cxx:358
 AliTPCSpaceCharge3D.cxx:359
 AliTPCSpaceCharge3D.cxx:360
 AliTPCSpaceCharge3D.cxx:361
 AliTPCSpaceCharge3D.cxx:362
 AliTPCSpaceCharge3D.cxx:363
 AliTPCSpaceCharge3D.cxx:364
 AliTPCSpaceCharge3D.cxx:365
 AliTPCSpaceCharge3D.cxx:366
 AliTPCSpaceCharge3D.cxx:367
 AliTPCSpaceCharge3D.cxx:368
 AliTPCSpaceCharge3D.cxx:369
 AliTPCSpaceCharge3D.cxx:370
 AliTPCSpaceCharge3D.cxx:371
 AliTPCSpaceCharge3D.cxx:372
 AliTPCSpaceCharge3D.cxx:373
 AliTPCSpaceCharge3D.cxx:374
 AliTPCSpaceCharge3D.cxx:375
 AliTPCSpaceCharge3D.cxx:376
 AliTPCSpaceCharge3D.cxx:377
 AliTPCSpaceCharge3D.cxx:378
 AliTPCSpaceCharge3D.cxx:379
 AliTPCSpaceCharge3D.cxx:380
 AliTPCSpaceCharge3D.cxx:381
 AliTPCSpaceCharge3D.cxx:382
 AliTPCSpaceCharge3D.cxx:383
 AliTPCSpaceCharge3D.cxx:384
 AliTPCSpaceCharge3D.cxx:385
 AliTPCSpaceCharge3D.cxx:386
 AliTPCSpaceCharge3D.cxx:387
 AliTPCSpaceCharge3D.cxx:388
 AliTPCSpaceCharge3D.cxx:389
 AliTPCSpaceCharge3D.cxx:390
 AliTPCSpaceCharge3D.cxx:391
 AliTPCSpaceCharge3D.cxx:392
 AliTPCSpaceCharge3D.cxx:393
 AliTPCSpaceCharge3D.cxx:394
 AliTPCSpaceCharge3D.cxx:395
 AliTPCSpaceCharge3D.cxx:396
 AliTPCSpaceCharge3D.cxx:397
 AliTPCSpaceCharge3D.cxx:398
 AliTPCSpaceCharge3D.cxx:399
 AliTPCSpaceCharge3D.cxx:400
 AliTPCSpaceCharge3D.cxx:401
 AliTPCSpaceCharge3D.cxx:402
 AliTPCSpaceCharge3D.cxx:403
 AliTPCSpaceCharge3D.cxx:404
 AliTPCSpaceCharge3D.cxx:405
 AliTPCSpaceCharge3D.cxx:406
 AliTPCSpaceCharge3D.cxx:407
 AliTPCSpaceCharge3D.cxx:408
 AliTPCSpaceCharge3D.cxx:409
 AliTPCSpaceCharge3D.cxx:410
 AliTPCSpaceCharge3D.cxx:411
 AliTPCSpaceCharge3D.cxx:412
 AliTPCSpaceCharge3D.cxx:413
 AliTPCSpaceCharge3D.cxx:414
 AliTPCSpaceCharge3D.cxx:415
 AliTPCSpaceCharge3D.cxx:416
 AliTPCSpaceCharge3D.cxx:417
 AliTPCSpaceCharge3D.cxx:418
 AliTPCSpaceCharge3D.cxx:419
 AliTPCSpaceCharge3D.cxx:420
 AliTPCSpaceCharge3D.cxx:421
 AliTPCSpaceCharge3D.cxx:422
 AliTPCSpaceCharge3D.cxx:423
 AliTPCSpaceCharge3D.cxx:424
 AliTPCSpaceCharge3D.cxx:425
 AliTPCSpaceCharge3D.cxx:426
 AliTPCSpaceCharge3D.cxx:427
 AliTPCSpaceCharge3D.cxx:428
 AliTPCSpaceCharge3D.cxx:429
 AliTPCSpaceCharge3D.cxx:430
 AliTPCSpaceCharge3D.cxx:431
 AliTPCSpaceCharge3D.cxx:432
 AliTPCSpaceCharge3D.cxx:433
 AliTPCSpaceCharge3D.cxx:434
 AliTPCSpaceCharge3D.cxx:435
 AliTPCSpaceCharge3D.cxx:436
 AliTPCSpaceCharge3D.cxx:437
 AliTPCSpaceCharge3D.cxx:438
 AliTPCSpaceCharge3D.cxx:439
 AliTPCSpaceCharge3D.cxx:440
 AliTPCSpaceCharge3D.cxx:441
 AliTPCSpaceCharge3D.cxx:442
 AliTPCSpaceCharge3D.cxx:443
 AliTPCSpaceCharge3D.cxx:444
 AliTPCSpaceCharge3D.cxx:445
 AliTPCSpaceCharge3D.cxx:446
 AliTPCSpaceCharge3D.cxx:447
 AliTPCSpaceCharge3D.cxx:448
 AliTPCSpaceCharge3D.cxx:449
 AliTPCSpaceCharge3D.cxx:450
 AliTPCSpaceCharge3D.cxx:451
 AliTPCSpaceCharge3D.cxx:452
 AliTPCSpaceCharge3D.cxx:453
 AliTPCSpaceCharge3D.cxx:454
 AliTPCSpaceCharge3D.cxx:455
 AliTPCSpaceCharge3D.cxx:456
 AliTPCSpaceCharge3D.cxx:457
 AliTPCSpaceCharge3D.cxx:458
 AliTPCSpaceCharge3D.cxx:459
 AliTPCSpaceCharge3D.cxx:460
 AliTPCSpaceCharge3D.cxx:461
 AliTPCSpaceCharge3D.cxx:462
 AliTPCSpaceCharge3D.cxx:463
 AliTPCSpaceCharge3D.cxx:464
 AliTPCSpaceCharge3D.cxx:465
 AliTPCSpaceCharge3D.cxx:466
 AliTPCSpaceCharge3D.cxx:467
 AliTPCSpaceCharge3D.cxx:468
 AliTPCSpaceCharge3D.cxx:469
 AliTPCSpaceCharge3D.cxx:470
 AliTPCSpaceCharge3D.cxx:471
 AliTPCSpaceCharge3D.cxx:472
 AliTPCSpaceCharge3D.cxx:473
 AliTPCSpaceCharge3D.cxx:474
 AliTPCSpaceCharge3D.cxx:475
 AliTPCSpaceCharge3D.cxx:476
 AliTPCSpaceCharge3D.cxx:477
 AliTPCSpaceCharge3D.cxx:478
 AliTPCSpaceCharge3D.cxx:479
 AliTPCSpaceCharge3D.cxx:480
 AliTPCSpaceCharge3D.cxx:481
 AliTPCSpaceCharge3D.cxx:482
 AliTPCSpaceCharge3D.cxx:483
 AliTPCSpaceCharge3D.cxx:484
 AliTPCSpaceCharge3D.cxx:485
 AliTPCSpaceCharge3D.cxx:486
 AliTPCSpaceCharge3D.cxx:487
 AliTPCSpaceCharge3D.cxx:488
 AliTPCSpaceCharge3D.cxx:489
 AliTPCSpaceCharge3D.cxx:490
 AliTPCSpaceCharge3D.cxx:491
 AliTPCSpaceCharge3D.cxx:492
 AliTPCSpaceCharge3D.cxx:493
 AliTPCSpaceCharge3D.cxx:494
 AliTPCSpaceCharge3D.cxx:495
 AliTPCSpaceCharge3D.cxx:496
 AliTPCSpaceCharge3D.cxx:497
 AliTPCSpaceCharge3D.cxx:498
 AliTPCSpaceCharge3D.cxx:499
 AliTPCSpaceCharge3D.cxx:500
 AliTPCSpaceCharge3D.cxx:501
 AliTPCSpaceCharge3D.cxx:502
 AliTPCSpaceCharge3D.cxx:503
 AliTPCSpaceCharge3D.cxx:504
 AliTPCSpaceCharge3D.cxx:505
 AliTPCSpaceCharge3D.cxx:506
 AliTPCSpaceCharge3D.cxx:507
 AliTPCSpaceCharge3D.cxx:508
 AliTPCSpaceCharge3D.cxx:509
 AliTPCSpaceCharge3D.cxx:510
 AliTPCSpaceCharge3D.cxx:511
 AliTPCSpaceCharge3D.cxx:512
 AliTPCSpaceCharge3D.cxx:513
 AliTPCSpaceCharge3D.cxx:514
 AliTPCSpaceCharge3D.cxx:515
 AliTPCSpaceCharge3D.cxx:516
 AliTPCSpaceCharge3D.cxx:517
 AliTPCSpaceCharge3D.cxx:518
 AliTPCSpaceCharge3D.cxx:519
 AliTPCSpaceCharge3D.cxx:520
 AliTPCSpaceCharge3D.cxx:521
 AliTPCSpaceCharge3D.cxx:522
 AliTPCSpaceCharge3D.cxx:523
 AliTPCSpaceCharge3D.cxx:524
 AliTPCSpaceCharge3D.cxx:525
 AliTPCSpaceCharge3D.cxx:526
 AliTPCSpaceCharge3D.cxx:527
 AliTPCSpaceCharge3D.cxx:528
 AliTPCSpaceCharge3D.cxx:529
 AliTPCSpaceCharge3D.cxx:530
 AliTPCSpaceCharge3D.cxx:531
 AliTPCSpaceCharge3D.cxx:532
 AliTPCSpaceCharge3D.cxx:533
 AliTPCSpaceCharge3D.cxx:534
 AliTPCSpaceCharge3D.cxx:535
 AliTPCSpaceCharge3D.cxx:536
 AliTPCSpaceCharge3D.cxx:537
 AliTPCSpaceCharge3D.cxx:538
 AliTPCSpaceCharge3D.cxx:539
 AliTPCSpaceCharge3D.cxx:540
 AliTPCSpaceCharge3D.cxx:541
 AliTPCSpaceCharge3D.cxx:542
 AliTPCSpaceCharge3D.cxx:543
 AliTPCSpaceCharge3D.cxx:544
 AliTPCSpaceCharge3D.cxx:545
 AliTPCSpaceCharge3D.cxx:546
 AliTPCSpaceCharge3D.cxx:547
 AliTPCSpaceCharge3D.cxx:548
 AliTPCSpaceCharge3D.cxx:549
 AliTPCSpaceCharge3D.cxx:550
 AliTPCSpaceCharge3D.cxx:551
 AliTPCSpaceCharge3D.cxx:552
 AliTPCSpaceCharge3D.cxx:553
 AliTPCSpaceCharge3D.cxx:554
 AliTPCSpaceCharge3D.cxx:555
 AliTPCSpaceCharge3D.cxx:556
 AliTPCSpaceCharge3D.cxx:557
 AliTPCSpaceCharge3D.cxx:558
 AliTPCSpaceCharge3D.cxx:559
 AliTPCSpaceCharge3D.cxx:560
 AliTPCSpaceCharge3D.cxx:561
 AliTPCSpaceCharge3D.cxx:562
 AliTPCSpaceCharge3D.cxx:563
 AliTPCSpaceCharge3D.cxx:564
 AliTPCSpaceCharge3D.cxx:565
 AliTPCSpaceCharge3D.cxx:566
 AliTPCSpaceCharge3D.cxx:567
 AliTPCSpaceCharge3D.cxx:568
 AliTPCSpaceCharge3D.cxx:569
 AliTPCSpaceCharge3D.cxx:570
 AliTPCSpaceCharge3D.cxx:571
 AliTPCSpaceCharge3D.cxx:572
 AliTPCSpaceCharge3D.cxx:573
 AliTPCSpaceCharge3D.cxx:574
 AliTPCSpaceCharge3D.cxx:575
 AliTPCSpaceCharge3D.cxx:576
 AliTPCSpaceCharge3D.cxx:577
 AliTPCSpaceCharge3D.cxx:578
 AliTPCSpaceCharge3D.cxx:579
 AliTPCSpaceCharge3D.cxx:580
 AliTPCSpaceCharge3D.cxx:581
 AliTPCSpaceCharge3D.cxx:582
 AliTPCSpaceCharge3D.cxx:583
 AliTPCSpaceCharge3D.cxx:584
 AliTPCSpaceCharge3D.cxx:585
 AliTPCSpaceCharge3D.cxx:586
 AliTPCSpaceCharge3D.cxx:587
 AliTPCSpaceCharge3D.cxx:588
 AliTPCSpaceCharge3D.cxx:589
 AliTPCSpaceCharge3D.cxx:590
 AliTPCSpaceCharge3D.cxx:591
 AliTPCSpaceCharge3D.cxx:592
 AliTPCSpaceCharge3D.cxx:593
 AliTPCSpaceCharge3D.cxx:594
 AliTPCSpaceCharge3D.cxx:595
 AliTPCSpaceCharge3D.cxx:596
 AliTPCSpaceCharge3D.cxx:597
 AliTPCSpaceCharge3D.cxx:598
 AliTPCSpaceCharge3D.cxx:599
 AliTPCSpaceCharge3D.cxx:600
 AliTPCSpaceCharge3D.cxx:601
 AliTPCSpaceCharge3D.cxx:602
 AliTPCSpaceCharge3D.cxx:603
 AliTPCSpaceCharge3D.cxx:604
 AliTPCSpaceCharge3D.cxx:605
 AliTPCSpaceCharge3D.cxx:606
 AliTPCSpaceCharge3D.cxx:607
 AliTPCSpaceCharge3D.cxx:608
 AliTPCSpaceCharge3D.cxx:609
 AliTPCSpaceCharge3D.cxx:610
 AliTPCSpaceCharge3D.cxx:611
 AliTPCSpaceCharge3D.cxx:612
 AliTPCSpaceCharge3D.cxx:613
 AliTPCSpaceCharge3D.cxx:614
 AliTPCSpaceCharge3D.cxx:615
 AliTPCSpaceCharge3D.cxx:616
 AliTPCSpaceCharge3D.cxx:617
 AliTPCSpaceCharge3D.cxx:618
 AliTPCSpaceCharge3D.cxx:619
 AliTPCSpaceCharge3D.cxx:620
 AliTPCSpaceCharge3D.cxx:621
 AliTPCSpaceCharge3D.cxx:622
 AliTPCSpaceCharge3D.cxx:623
 AliTPCSpaceCharge3D.cxx:624
 AliTPCSpaceCharge3D.cxx:625
 AliTPCSpaceCharge3D.cxx:626
 AliTPCSpaceCharge3D.cxx:627
 AliTPCSpaceCharge3D.cxx:628
 AliTPCSpaceCharge3D.cxx:629
 AliTPCSpaceCharge3D.cxx:630
 AliTPCSpaceCharge3D.cxx:631
 AliTPCSpaceCharge3D.cxx:632
 AliTPCSpaceCharge3D.cxx:633
 AliTPCSpaceCharge3D.cxx:634
 AliTPCSpaceCharge3D.cxx:635
 AliTPCSpaceCharge3D.cxx:636
 AliTPCSpaceCharge3D.cxx:637
 AliTPCSpaceCharge3D.cxx:638
 AliTPCSpaceCharge3D.cxx:639
 AliTPCSpaceCharge3D.cxx:640
 AliTPCSpaceCharge3D.cxx:641
 AliTPCSpaceCharge3D.cxx:642
 AliTPCSpaceCharge3D.cxx:643
 AliTPCSpaceCharge3D.cxx:644
 AliTPCSpaceCharge3D.cxx:645
 AliTPCSpaceCharge3D.cxx:646
 AliTPCSpaceCharge3D.cxx:647
 AliTPCSpaceCharge3D.cxx:648
 AliTPCSpaceCharge3D.cxx:649
 AliTPCSpaceCharge3D.cxx:650
 AliTPCSpaceCharge3D.cxx:651
 AliTPCSpaceCharge3D.cxx:652
 AliTPCSpaceCharge3D.cxx:653
 AliTPCSpaceCharge3D.cxx:654
 AliTPCSpaceCharge3D.cxx:655
 AliTPCSpaceCharge3D.cxx:656
 AliTPCSpaceCharge3D.cxx:657
 AliTPCSpaceCharge3D.cxx:658
 AliTPCSpaceCharge3D.cxx:659
 AliTPCSpaceCharge3D.cxx:660
 AliTPCSpaceCharge3D.cxx:661
 AliTPCSpaceCharge3D.cxx:662
 AliTPCSpaceCharge3D.cxx:663
 AliTPCSpaceCharge3D.cxx:664
 AliTPCSpaceCharge3D.cxx:665
 AliTPCSpaceCharge3D.cxx:666
 AliTPCSpaceCharge3D.cxx:667
 AliTPCSpaceCharge3D.cxx:668
 AliTPCSpaceCharge3D.cxx:669
 AliTPCSpaceCharge3D.cxx:670
 AliTPCSpaceCharge3D.cxx:671
 AliTPCSpaceCharge3D.cxx:672
 AliTPCSpaceCharge3D.cxx:673
 AliTPCSpaceCharge3D.cxx:674
 AliTPCSpaceCharge3D.cxx:675
 AliTPCSpaceCharge3D.cxx:676
 AliTPCSpaceCharge3D.cxx:677
 AliTPCSpaceCharge3D.cxx:678
 AliTPCSpaceCharge3D.cxx:679
 AliTPCSpaceCharge3D.cxx:680
 AliTPCSpaceCharge3D.cxx:681
 AliTPCSpaceCharge3D.cxx:682
 AliTPCSpaceCharge3D.cxx:683
 AliTPCSpaceCharge3D.cxx:684
 AliTPCSpaceCharge3D.cxx:685
 AliTPCSpaceCharge3D.cxx:686
 AliTPCSpaceCharge3D.cxx:687
 AliTPCSpaceCharge3D.cxx:688
 AliTPCSpaceCharge3D.cxx:689
 AliTPCSpaceCharge3D.cxx:690
 AliTPCSpaceCharge3D.cxx:691
 AliTPCSpaceCharge3D.cxx:692
 AliTPCSpaceCharge3D.cxx:693
 AliTPCSpaceCharge3D.cxx:694
 AliTPCSpaceCharge3D.cxx:695
 AliTPCSpaceCharge3D.cxx:696
 AliTPCSpaceCharge3D.cxx:697
 AliTPCSpaceCharge3D.cxx:698
 AliTPCSpaceCharge3D.cxx:699
 AliTPCSpaceCharge3D.cxx:700
 AliTPCSpaceCharge3D.cxx:701
 AliTPCSpaceCharge3D.cxx:702
 AliTPCSpaceCharge3D.cxx:703
 AliTPCSpaceCharge3D.cxx:704
 AliTPCSpaceCharge3D.cxx:705
 AliTPCSpaceCharge3D.cxx:706
 AliTPCSpaceCharge3D.cxx:707
 AliTPCSpaceCharge3D.cxx:708
 AliTPCSpaceCharge3D.cxx:709
 AliTPCSpaceCharge3D.cxx:710
 AliTPCSpaceCharge3D.cxx:711
 AliTPCSpaceCharge3D.cxx:712
 AliTPCSpaceCharge3D.cxx:713
 AliTPCSpaceCharge3D.cxx:714
 AliTPCSpaceCharge3D.cxx:715
 AliTPCSpaceCharge3D.cxx:716
 AliTPCSpaceCharge3D.cxx:717
 AliTPCSpaceCharge3D.cxx:718
 AliTPCSpaceCharge3D.cxx:719
 AliTPCSpaceCharge3D.cxx:720
 AliTPCSpaceCharge3D.cxx:721
 AliTPCSpaceCharge3D.cxx:722
 AliTPCSpaceCharge3D.cxx:723
 AliTPCSpaceCharge3D.cxx:724
 AliTPCSpaceCharge3D.cxx:725
 AliTPCSpaceCharge3D.cxx:726
 AliTPCSpaceCharge3D.cxx:727
 AliTPCSpaceCharge3D.cxx:728
 AliTPCSpaceCharge3D.cxx:729
 AliTPCSpaceCharge3D.cxx:730
 AliTPCSpaceCharge3D.cxx:731
 AliTPCSpaceCharge3D.cxx:732
 AliTPCSpaceCharge3D.cxx:733
 AliTPCSpaceCharge3D.cxx:734
 AliTPCSpaceCharge3D.cxx:735
 AliTPCSpaceCharge3D.cxx:736
 AliTPCSpaceCharge3D.cxx:737
 AliTPCSpaceCharge3D.cxx:738
 AliTPCSpaceCharge3D.cxx:739
 AliTPCSpaceCharge3D.cxx:740
 AliTPCSpaceCharge3D.cxx:741
 AliTPCSpaceCharge3D.cxx:742
 AliTPCSpaceCharge3D.cxx:743
 AliTPCSpaceCharge3D.cxx:744
 AliTPCSpaceCharge3D.cxx:745
 AliTPCSpaceCharge3D.cxx:746
 AliTPCSpaceCharge3D.cxx:747
 AliTPCSpaceCharge3D.cxx:748
 AliTPCSpaceCharge3D.cxx:749
 AliTPCSpaceCharge3D.cxx:750
 AliTPCSpaceCharge3D.cxx:751
 AliTPCSpaceCharge3D.cxx:752
 AliTPCSpaceCharge3D.cxx:753
 AliTPCSpaceCharge3D.cxx:754
 AliTPCSpaceCharge3D.cxx:755
 AliTPCSpaceCharge3D.cxx:756
 AliTPCSpaceCharge3D.cxx:757
 AliTPCSpaceCharge3D.cxx:758
 AliTPCSpaceCharge3D.cxx:759
 AliTPCSpaceCharge3D.cxx:760
 AliTPCSpaceCharge3D.cxx:761
 AliTPCSpaceCharge3D.cxx:762
 AliTPCSpaceCharge3D.cxx:763
 AliTPCSpaceCharge3D.cxx:764
 AliTPCSpaceCharge3D.cxx:765
 AliTPCSpaceCharge3D.cxx:766
 AliTPCSpaceCharge3D.cxx:767
 AliTPCSpaceCharge3D.cxx:768
 AliTPCSpaceCharge3D.cxx:769
 AliTPCSpaceCharge3D.cxx:770
 AliTPCSpaceCharge3D.cxx:771
 AliTPCSpaceCharge3D.cxx:772
 AliTPCSpaceCharge3D.cxx:773
 AliTPCSpaceCharge3D.cxx:774
 AliTPCSpaceCharge3D.cxx:775
 AliTPCSpaceCharge3D.cxx:776
 AliTPCSpaceCharge3D.cxx:777
 AliTPCSpaceCharge3D.cxx:778
 AliTPCSpaceCharge3D.cxx:779
 AliTPCSpaceCharge3D.cxx:780
 AliTPCSpaceCharge3D.cxx:781
 AliTPCSpaceCharge3D.cxx:782
 AliTPCSpaceCharge3D.cxx:783
 AliTPCSpaceCharge3D.cxx:784
 AliTPCSpaceCharge3D.cxx:785
 AliTPCSpaceCharge3D.cxx:786
 AliTPCSpaceCharge3D.cxx:787
 AliTPCSpaceCharge3D.cxx:788
 AliTPCSpaceCharge3D.cxx:789
 AliTPCSpaceCharge3D.cxx:790
 AliTPCSpaceCharge3D.cxx:791
 AliTPCSpaceCharge3D.cxx:792
 AliTPCSpaceCharge3D.cxx:793
 AliTPCSpaceCharge3D.cxx:794
 AliTPCSpaceCharge3D.cxx:795
 AliTPCSpaceCharge3D.cxx:796
 AliTPCSpaceCharge3D.cxx:797
 AliTPCSpaceCharge3D.cxx:798
 AliTPCSpaceCharge3D.cxx:799
 AliTPCSpaceCharge3D.cxx:800
 AliTPCSpaceCharge3D.cxx:801
 AliTPCSpaceCharge3D.cxx:802
 AliTPCSpaceCharge3D.cxx:803
 AliTPCSpaceCharge3D.cxx:804
 AliTPCSpaceCharge3D.cxx:805
 AliTPCSpaceCharge3D.cxx:806
 AliTPCSpaceCharge3D.cxx:807
 AliTPCSpaceCharge3D.cxx:808
 AliTPCSpaceCharge3D.cxx:809
 AliTPCSpaceCharge3D.cxx:810
 AliTPCSpaceCharge3D.cxx:811
 AliTPCSpaceCharge3D.cxx:812
 AliTPCSpaceCharge3D.cxx:813
 AliTPCSpaceCharge3D.cxx:814
 AliTPCSpaceCharge3D.cxx:815
 AliTPCSpaceCharge3D.cxx:816
 AliTPCSpaceCharge3D.cxx:817
 AliTPCSpaceCharge3D.cxx:818
 AliTPCSpaceCharge3D.cxx:819
 AliTPCSpaceCharge3D.cxx:820
 AliTPCSpaceCharge3D.cxx:821
 AliTPCSpaceCharge3D.cxx:822
 AliTPCSpaceCharge3D.cxx:823
 AliTPCSpaceCharge3D.cxx:824
 AliTPCSpaceCharge3D.cxx:825
 AliTPCSpaceCharge3D.cxx:826
 AliTPCSpaceCharge3D.cxx:827
 AliTPCSpaceCharge3D.cxx:828
 AliTPCSpaceCharge3D.cxx:829
 AliTPCSpaceCharge3D.cxx:830
 AliTPCSpaceCharge3D.cxx:831
 AliTPCSpaceCharge3D.cxx:832
 AliTPCSpaceCharge3D.cxx:833
 AliTPCSpaceCharge3D.cxx:834
 AliTPCSpaceCharge3D.cxx:835
 AliTPCSpaceCharge3D.cxx:836
 AliTPCSpaceCharge3D.cxx:837
 AliTPCSpaceCharge3D.cxx:838
 AliTPCSpaceCharge3D.cxx:839
 AliTPCSpaceCharge3D.cxx:840
 AliTPCSpaceCharge3D.cxx:841
 AliTPCSpaceCharge3D.cxx:842
 AliTPCSpaceCharge3D.cxx:843
 AliTPCSpaceCharge3D.cxx:844
 AliTPCSpaceCharge3D.cxx:845
 AliTPCSpaceCharge3D.cxx:846
 AliTPCSpaceCharge3D.cxx:847
 AliTPCSpaceCharge3D.cxx:848
 AliTPCSpaceCharge3D.cxx:849
 AliTPCSpaceCharge3D.cxx:850
 AliTPCSpaceCharge3D.cxx:851
 AliTPCSpaceCharge3D.cxx:852
 AliTPCSpaceCharge3D.cxx:853
 AliTPCSpaceCharge3D.cxx:854
 AliTPCSpaceCharge3D.cxx:855
 AliTPCSpaceCharge3D.cxx:856
 AliTPCSpaceCharge3D.cxx:857
 AliTPCSpaceCharge3D.cxx:858
 AliTPCSpaceCharge3D.cxx:859
 AliTPCSpaceCharge3D.cxx:860
 AliTPCSpaceCharge3D.cxx:861
 AliTPCSpaceCharge3D.cxx:862
 AliTPCSpaceCharge3D.cxx:863
 AliTPCSpaceCharge3D.cxx:864
 AliTPCSpaceCharge3D.cxx:865
 AliTPCSpaceCharge3D.cxx:866
 AliTPCSpaceCharge3D.cxx:867
 AliTPCSpaceCharge3D.cxx:868
 AliTPCSpaceCharge3D.cxx:869
 AliTPCSpaceCharge3D.cxx:870
 AliTPCSpaceCharge3D.cxx:871
 AliTPCSpaceCharge3D.cxx:872
 AliTPCSpaceCharge3D.cxx:873
 AliTPCSpaceCharge3D.cxx:874
 AliTPCSpaceCharge3D.cxx:875
 AliTPCSpaceCharge3D.cxx:876
 AliTPCSpaceCharge3D.cxx:877
 AliTPCSpaceCharge3D.cxx:878
 AliTPCSpaceCharge3D.cxx:879
 AliTPCSpaceCharge3D.cxx:880
 AliTPCSpaceCharge3D.cxx:881
 AliTPCSpaceCharge3D.cxx:882
 AliTPCSpaceCharge3D.cxx:883
 AliTPCSpaceCharge3D.cxx:884
 AliTPCSpaceCharge3D.cxx:885
 AliTPCSpaceCharge3D.cxx:886
 AliTPCSpaceCharge3D.cxx:887
 AliTPCSpaceCharge3D.cxx:888
 AliTPCSpaceCharge3D.cxx:889
 AliTPCSpaceCharge3D.cxx:890
 AliTPCSpaceCharge3D.cxx:891
 AliTPCSpaceCharge3D.cxx:892
 AliTPCSpaceCharge3D.cxx:893
 AliTPCSpaceCharge3D.cxx:894
 AliTPCSpaceCharge3D.cxx:895
 AliTPCSpaceCharge3D.cxx:896
 AliTPCSpaceCharge3D.cxx:897
 AliTPCSpaceCharge3D.cxx:898
 AliTPCSpaceCharge3D.cxx:899
 AliTPCSpaceCharge3D.cxx:900
 AliTPCSpaceCharge3D.cxx:901
 AliTPCSpaceCharge3D.cxx:902
 AliTPCSpaceCharge3D.cxx:903
 AliTPCSpaceCharge3D.cxx:904
 AliTPCSpaceCharge3D.cxx:905
 AliTPCSpaceCharge3D.cxx:906
 AliTPCSpaceCharge3D.cxx:907
 AliTPCSpaceCharge3D.cxx:908
 AliTPCSpaceCharge3D.cxx:909
 AliTPCSpaceCharge3D.cxx:910
 AliTPCSpaceCharge3D.cxx:911
 AliTPCSpaceCharge3D.cxx:912
 AliTPCSpaceCharge3D.cxx:913
 AliTPCSpaceCharge3D.cxx:914
 AliTPCSpaceCharge3D.cxx:915
 AliTPCSpaceCharge3D.cxx:916
 AliTPCSpaceCharge3D.cxx:917
 AliTPCSpaceCharge3D.cxx:918
 AliTPCSpaceCharge3D.cxx:919
 AliTPCSpaceCharge3D.cxx:920
 AliTPCSpaceCharge3D.cxx:921
 AliTPCSpaceCharge3D.cxx:922
 AliTPCSpaceCharge3D.cxx:923
 AliTPCSpaceCharge3D.cxx:924
 AliTPCSpaceCharge3D.cxx:925
 AliTPCSpaceCharge3D.cxx:926
 AliTPCSpaceCharge3D.cxx:927
 AliTPCSpaceCharge3D.cxx:928
 AliTPCSpaceCharge3D.cxx:929
 AliTPCSpaceCharge3D.cxx:930
 AliTPCSpaceCharge3D.cxx:931
 AliTPCSpaceCharge3D.cxx:932
 AliTPCSpaceCharge3D.cxx:933
 AliTPCSpaceCharge3D.cxx:934
 AliTPCSpaceCharge3D.cxx:935
 AliTPCSpaceCharge3D.cxx:936
 AliTPCSpaceCharge3D.cxx:937
 AliTPCSpaceCharge3D.cxx:938
 AliTPCSpaceCharge3D.cxx:939
 AliTPCSpaceCharge3D.cxx:940
 AliTPCSpaceCharge3D.cxx:941
 AliTPCSpaceCharge3D.cxx:942
 AliTPCSpaceCharge3D.cxx:943
 AliTPCSpaceCharge3D.cxx:944
 AliTPCSpaceCharge3D.cxx:945
 AliTPCSpaceCharge3D.cxx:946
 AliTPCSpaceCharge3D.cxx:947
 AliTPCSpaceCharge3D.cxx:948
 AliTPCSpaceCharge3D.cxx:949
 AliTPCSpaceCharge3D.cxx:950
 AliTPCSpaceCharge3D.cxx:951
 AliTPCSpaceCharge3D.cxx:952
 AliTPCSpaceCharge3D.cxx:953
 AliTPCSpaceCharge3D.cxx:954
 AliTPCSpaceCharge3D.cxx:955
 AliTPCSpaceCharge3D.cxx:956
 AliTPCSpaceCharge3D.cxx:957
 AliTPCSpaceCharge3D.cxx:958
 AliTPCSpaceCharge3D.cxx:959
 AliTPCSpaceCharge3D.cxx:960
 AliTPCSpaceCharge3D.cxx:961
 AliTPCSpaceCharge3D.cxx:962
 AliTPCSpaceCharge3D.cxx:963
 AliTPCSpaceCharge3D.cxx:964
 AliTPCSpaceCharge3D.cxx:965
 AliTPCSpaceCharge3D.cxx:966
 AliTPCSpaceCharge3D.cxx:967
 AliTPCSpaceCharge3D.cxx:968
 AliTPCSpaceCharge3D.cxx:969
 AliTPCSpaceCharge3D.cxx:970
 AliTPCSpaceCharge3D.cxx:971
 AliTPCSpaceCharge3D.cxx:972
 AliTPCSpaceCharge3D.cxx:973
 AliTPCSpaceCharge3D.cxx:974
 AliTPCSpaceCharge3D.cxx:975
 AliTPCSpaceCharge3D.cxx:976
 AliTPCSpaceCharge3D.cxx:977
 AliTPCSpaceCharge3D.cxx:978
 AliTPCSpaceCharge3D.cxx:979
 AliTPCSpaceCharge3D.cxx:980
 AliTPCSpaceCharge3D.cxx:981
 AliTPCSpaceCharge3D.cxx:982
 AliTPCSpaceCharge3D.cxx:983
 AliTPCSpaceCharge3D.cxx:984
 AliTPCSpaceCharge3D.cxx:985
 AliTPCSpaceCharge3D.cxx:986
 AliTPCSpaceCharge3D.cxx:987
 AliTPCSpaceCharge3D.cxx:988
 AliTPCSpaceCharge3D.cxx:989
 AliTPCSpaceCharge3D.cxx:990
 AliTPCSpaceCharge3D.cxx:991
 AliTPCSpaceCharge3D.cxx:992
 AliTPCSpaceCharge3D.cxx:993
 AliTPCSpaceCharge3D.cxx:994
 AliTPCSpaceCharge3D.cxx:995
 AliTPCSpaceCharge3D.cxx:996
 AliTPCSpaceCharge3D.cxx:997
 AliTPCSpaceCharge3D.cxx:998
 AliTPCSpaceCharge3D.cxx:999
 AliTPCSpaceCharge3D.cxx:1000
 AliTPCSpaceCharge3D.cxx:1001
 AliTPCSpaceCharge3D.cxx:1002
 AliTPCSpaceCharge3D.cxx:1003
 AliTPCSpaceCharge3D.cxx:1004
 AliTPCSpaceCharge3D.cxx:1005
 AliTPCSpaceCharge3D.cxx:1006
 AliTPCSpaceCharge3D.cxx:1007
 AliTPCSpaceCharge3D.cxx:1008
 AliTPCSpaceCharge3D.cxx:1009
 AliTPCSpaceCharge3D.cxx:1010
 AliTPCSpaceCharge3D.cxx:1011
 AliTPCSpaceCharge3D.cxx:1012
 AliTPCSpaceCharge3D.cxx:1013
 AliTPCSpaceCharge3D.cxx:1014
 AliTPCSpaceCharge3D.cxx:1015
 AliTPCSpaceCharge3D.cxx:1016
 AliTPCSpaceCharge3D.cxx:1017
 AliTPCSpaceCharge3D.cxx:1018
 AliTPCSpaceCharge3D.cxx:1019
 AliTPCSpaceCharge3D.cxx:1020
 AliTPCSpaceCharge3D.cxx:1021
 AliTPCSpaceCharge3D.cxx:1022
 AliTPCSpaceCharge3D.cxx:1023
 AliTPCSpaceCharge3D.cxx:1024
 AliTPCSpaceCharge3D.cxx:1025
 AliTPCSpaceCharge3D.cxx:1026
 AliTPCSpaceCharge3D.cxx:1027
 AliTPCSpaceCharge3D.cxx:1028
 AliTPCSpaceCharge3D.cxx:1029
 AliTPCSpaceCharge3D.cxx:1030
 AliTPCSpaceCharge3D.cxx:1031
 AliTPCSpaceCharge3D.cxx:1032
 AliTPCSpaceCharge3D.cxx:1033
 AliTPCSpaceCharge3D.cxx:1034
 AliTPCSpaceCharge3D.cxx:1035
 AliTPCSpaceCharge3D.cxx:1036
 AliTPCSpaceCharge3D.cxx:1037
 AliTPCSpaceCharge3D.cxx:1038
 AliTPCSpaceCharge3D.cxx:1039
 AliTPCSpaceCharge3D.cxx:1040
 AliTPCSpaceCharge3D.cxx:1041
 AliTPCSpaceCharge3D.cxx:1042
 AliTPCSpaceCharge3D.cxx:1043
 AliTPCSpaceCharge3D.cxx:1044
 AliTPCSpaceCharge3D.cxx:1045
 AliTPCSpaceCharge3D.cxx:1046
 AliTPCSpaceCharge3D.cxx:1047
 AliTPCSpaceCharge3D.cxx:1048
 AliTPCSpaceCharge3D.cxx:1049
 AliTPCSpaceCharge3D.cxx:1050
 AliTPCSpaceCharge3D.cxx:1051
 AliTPCSpaceCharge3D.cxx:1052
 AliTPCSpaceCharge3D.cxx:1053
 AliTPCSpaceCharge3D.cxx:1054
 AliTPCSpaceCharge3D.cxx:1055
 AliTPCSpaceCharge3D.cxx:1056
 AliTPCSpaceCharge3D.cxx:1057
 AliTPCSpaceCharge3D.cxx:1058
 AliTPCSpaceCharge3D.cxx:1059
 AliTPCSpaceCharge3D.cxx:1060
 AliTPCSpaceCharge3D.cxx:1061
 AliTPCSpaceCharge3D.cxx:1062
 AliTPCSpaceCharge3D.cxx:1063
 AliTPCSpaceCharge3D.cxx:1064
 AliTPCSpaceCharge3D.cxx:1065
 AliTPCSpaceCharge3D.cxx:1066
 AliTPCSpaceCharge3D.cxx:1067
 AliTPCSpaceCharge3D.cxx:1068
 AliTPCSpaceCharge3D.cxx:1069
 AliTPCSpaceCharge3D.cxx:1070
 AliTPCSpaceCharge3D.cxx:1071
 AliTPCSpaceCharge3D.cxx:1072
 AliTPCSpaceCharge3D.cxx:1073
 AliTPCSpaceCharge3D.cxx:1074
 AliTPCSpaceCharge3D.cxx:1075
 AliTPCSpaceCharge3D.cxx:1076
 AliTPCSpaceCharge3D.cxx:1077
 AliTPCSpaceCharge3D.cxx:1078
 AliTPCSpaceCharge3D.cxx:1079
 AliTPCSpaceCharge3D.cxx:1080
 AliTPCSpaceCharge3D.cxx:1081
 AliTPCSpaceCharge3D.cxx:1082
 AliTPCSpaceCharge3D.cxx:1083
 AliTPCSpaceCharge3D.cxx:1084
 AliTPCSpaceCharge3D.cxx:1085
 AliTPCSpaceCharge3D.cxx:1086
 AliTPCSpaceCharge3D.cxx:1087
 AliTPCSpaceCharge3D.cxx:1088
 AliTPCSpaceCharge3D.cxx:1089
 AliTPCSpaceCharge3D.cxx:1090
 AliTPCSpaceCharge3D.cxx:1091
 AliTPCSpaceCharge3D.cxx:1092
 AliTPCSpaceCharge3D.cxx:1093
 AliTPCSpaceCharge3D.cxx:1094
 AliTPCSpaceCharge3D.cxx:1095
 AliTPCSpaceCharge3D.cxx:1096
 AliTPCSpaceCharge3D.cxx:1097
 AliTPCSpaceCharge3D.cxx:1098
 AliTPCSpaceCharge3D.cxx:1099
 AliTPCSpaceCharge3D.cxx:1100
 AliTPCSpaceCharge3D.cxx:1101
 AliTPCSpaceCharge3D.cxx:1102
 AliTPCSpaceCharge3D.cxx:1103
 AliTPCSpaceCharge3D.cxx:1104
 AliTPCSpaceCharge3D.cxx:1105
 AliTPCSpaceCharge3D.cxx:1106
 AliTPCSpaceCharge3D.cxx:1107
 AliTPCSpaceCharge3D.cxx:1108
 AliTPCSpaceCharge3D.cxx:1109
 AliTPCSpaceCharge3D.cxx:1110
 AliTPCSpaceCharge3D.cxx:1111
 AliTPCSpaceCharge3D.cxx:1112
 AliTPCSpaceCharge3D.cxx:1113
 AliTPCSpaceCharge3D.cxx:1114
 AliTPCSpaceCharge3D.cxx:1115
 AliTPCSpaceCharge3D.cxx:1116
 AliTPCSpaceCharge3D.cxx:1117
 AliTPCSpaceCharge3D.cxx:1118
 AliTPCSpaceCharge3D.cxx:1119
 AliTPCSpaceCharge3D.cxx:1120
 AliTPCSpaceCharge3D.cxx:1121
 AliTPCSpaceCharge3D.cxx:1122
 AliTPCSpaceCharge3D.cxx:1123
 AliTPCSpaceCharge3D.cxx:1124
 AliTPCSpaceCharge3D.cxx:1125
 AliTPCSpaceCharge3D.cxx:1126
 AliTPCSpaceCharge3D.cxx:1127
 AliTPCSpaceCharge3D.cxx:1128
 AliTPCSpaceCharge3D.cxx:1129
 AliTPCSpaceCharge3D.cxx:1130
 AliTPCSpaceCharge3D.cxx:1131
 AliTPCSpaceCharge3D.cxx:1132
 AliTPCSpaceCharge3D.cxx:1133
 AliTPCSpaceCharge3D.cxx:1134
 AliTPCSpaceCharge3D.cxx:1135
 AliTPCSpaceCharge3D.cxx:1136
 AliTPCSpaceCharge3D.cxx:1137
 AliTPCSpaceCharge3D.cxx:1138
 AliTPCSpaceCharge3D.cxx:1139
 AliTPCSpaceCharge3D.cxx:1140
 AliTPCSpaceCharge3D.cxx:1141
 AliTPCSpaceCharge3D.cxx:1142
 AliTPCSpaceCharge3D.cxx:1143
 AliTPCSpaceCharge3D.cxx:1144
 AliTPCSpaceCharge3D.cxx:1145
 AliTPCSpaceCharge3D.cxx:1146
 AliTPCSpaceCharge3D.cxx:1147
 AliTPCSpaceCharge3D.cxx:1148
 AliTPCSpaceCharge3D.cxx:1149
 AliTPCSpaceCharge3D.cxx:1150
 AliTPCSpaceCharge3D.cxx:1151
 AliTPCSpaceCharge3D.cxx:1152
 AliTPCSpaceCharge3D.cxx:1153
 AliTPCSpaceCharge3D.cxx:1154
 AliTPCSpaceCharge3D.cxx:1155
 AliTPCSpaceCharge3D.cxx:1156
 AliTPCSpaceCharge3D.cxx:1157
 AliTPCSpaceCharge3D.cxx:1158
 AliTPCSpaceCharge3D.cxx:1159
 AliTPCSpaceCharge3D.cxx:1160
 AliTPCSpaceCharge3D.cxx:1161
 AliTPCSpaceCharge3D.cxx:1162
 AliTPCSpaceCharge3D.cxx:1163
 AliTPCSpaceCharge3D.cxx:1164
 AliTPCSpaceCharge3D.cxx:1165
 AliTPCSpaceCharge3D.cxx:1166
 AliTPCSpaceCharge3D.cxx:1167
 AliTPCSpaceCharge3D.cxx:1168
 AliTPCSpaceCharge3D.cxx:1169
 AliTPCSpaceCharge3D.cxx:1170
 AliTPCSpaceCharge3D.cxx:1171
 AliTPCSpaceCharge3D.cxx:1172
 AliTPCSpaceCharge3D.cxx:1173
 AliTPCSpaceCharge3D.cxx:1174
 AliTPCSpaceCharge3D.cxx:1175
 AliTPCSpaceCharge3D.cxx:1176
 AliTPCSpaceCharge3D.cxx:1177
 AliTPCSpaceCharge3D.cxx:1178
 AliTPCSpaceCharge3D.cxx:1179
 AliTPCSpaceCharge3D.cxx:1180
 AliTPCSpaceCharge3D.cxx:1181
 AliTPCSpaceCharge3D.cxx:1182
 AliTPCSpaceCharge3D.cxx:1183
 AliTPCSpaceCharge3D.cxx:1184
 AliTPCSpaceCharge3D.cxx:1185
 AliTPCSpaceCharge3D.cxx:1186
 AliTPCSpaceCharge3D.cxx:1187
 AliTPCSpaceCharge3D.cxx:1188
 AliTPCSpaceCharge3D.cxx:1189
 AliTPCSpaceCharge3D.cxx:1190
 AliTPCSpaceCharge3D.cxx:1191
 AliTPCSpaceCharge3D.cxx:1192
 AliTPCSpaceCharge3D.cxx:1193
 AliTPCSpaceCharge3D.cxx:1194
 AliTPCSpaceCharge3D.cxx:1195
 AliTPCSpaceCharge3D.cxx:1196
 AliTPCSpaceCharge3D.cxx:1197
 AliTPCSpaceCharge3D.cxx:1198
 AliTPCSpaceCharge3D.cxx:1199
 AliTPCSpaceCharge3D.cxx:1200
 AliTPCSpaceCharge3D.cxx:1201
 AliTPCSpaceCharge3D.cxx:1202
 AliTPCSpaceCharge3D.cxx:1203
 AliTPCSpaceCharge3D.cxx:1204
 AliTPCSpaceCharge3D.cxx:1205
 AliTPCSpaceCharge3D.cxx:1206
 AliTPCSpaceCharge3D.cxx:1207
 AliTPCSpaceCharge3D.cxx:1208
 AliTPCSpaceCharge3D.cxx:1209
 AliTPCSpaceCharge3D.cxx:1210
 AliTPCSpaceCharge3D.cxx:1211
 AliTPCSpaceCharge3D.cxx:1212
 AliTPCSpaceCharge3D.cxx:1213
 AliTPCSpaceCharge3D.cxx:1214
 AliTPCSpaceCharge3D.cxx:1215
 AliTPCSpaceCharge3D.cxx:1216
 AliTPCSpaceCharge3D.cxx:1217
 AliTPCSpaceCharge3D.cxx:1218
 AliTPCSpaceCharge3D.cxx:1219
 AliTPCSpaceCharge3D.cxx:1220
 AliTPCSpaceCharge3D.cxx:1221
 AliTPCSpaceCharge3D.cxx:1222
 AliTPCSpaceCharge3D.cxx:1223
 AliTPCSpaceCharge3D.cxx:1224
 AliTPCSpaceCharge3D.cxx:1225
 AliTPCSpaceCharge3D.cxx:1226
 AliTPCSpaceCharge3D.cxx:1227
 AliTPCSpaceCharge3D.cxx:1228
 AliTPCSpaceCharge3D.cxx:1229
 AliTPCSpaceCharge3D.cxx:1230
 AliTPCSpaceCharge3D.cxx:1231
 AliTPCSpaceCharge3D.cxx:1232
 AliTPCSpaceCharge3D.cxx:1233
 AliTPCSpaceCharge3D.cxx:1234
 AliTPCSpaceCharge3D.cxx:1235
 AliTPCSpaceCharge3D.cxx:1236
 AliTPCSpaceCharge3D.cxx:1237
 AliTPCSpaceCharge3D.cxx:1238
 AliTPCSpaceCharge3D.cxx:1239
 AliTPCSpaceCharge3D.cxx:1240
 AliTPCSpaceCharge3D.cxx:1241
 AliTPCSpaceCharge3D.cxx:1242
 AliTPCSpaceCharge3D.cxx:1243
 AliTPCSpaceCharge3D.cxx:1244
 AliTPCSpaceCharge3D.cxx:1245
 AliTPCSpaceCharge3D.cxx:1246
 AliTPCSpaceCharge3D.cxx:1247
 AliTPCSpaceCharge3D.cxx:1248
 AliTPCSpaceCharge3D.cxx:1249
 AliTPCSpaceCharge3D.cxx:1250
 AliTPCSpaceCharge3D.cxx:1251
 AliTPCSpaceCharge3D.cxx:1252
 AliTPCSpaceCharge3D.cxx:1253
 AliTPCSpaceCharge3D.cxx:1254
 AliTPCSpaceCharge3D.cxx:1255
 AliTPCSpaceCharge3D.cxx:1256
 AliTPCSpaceCharge3D.cxx:1257
 AliTPCSpaceCharge3D.cxx:1258
 AliTPCSpaceCharge3D.cxx:1259
 AliTPCSpaceCharge3D.cxx:1260
 AliTPCSpaceCharge3D.cxx:1261
 AliTPCSpaceCharge3D.cxx:1262
 AliTPCSpaceCharge3D.cxx:1263
 AliTPCSpaceCharge3D.cxx:1264
 AliTPCSpaceCharge3D.cxx:1265
 AliTPCSpaceCharge3D.cxx:1266
 AliTPCSpaceCharge3D.cxx:1267
 AliTPCSpaceCharge3D.cxx:1268
 AliTPCSpaceCharge3D.cxx:1269
 AliTPCSpaceCharge3D.cxx:1270
 AliTPCSpaceCharge3D.cxx:1271
 AliTPCSpaceCharge3D.cxx:1272
 AliTPCSpaceCharge3D.cxx:1273
 AliTPCSpaceCharge3D.cxx:1274
 AliTPCSpaceCharge3D.cxx:1275
 AliTPCSpaceCharge3D.cxx:1276
 AliTPCSpaceCharge3D.cxx:1277
 AliTPCSpaceCharge3D.cxx:1278
 AliTPCSpaceCharge3D.cxx:1279
 AliTPCSpaceCharge3D.cxx:1280
 AliTPCSpaceCharge3D.cxx:1281
 AliTPCSpaceCharge3D.cxx:1282
 AliTPCSpaceCharge3D.cxx:1283
 AliTPCSpaceCharge3D.cxx:1284
 AliTPCSpaceCharge3D.cxx:1285
 AliTPCSpaceCharge3D.cxx:1286
 AliTPCSpaceCharge3D.cxx:1287
 AliTPCSpaceCharge3D.cxx:1288
 AliTPCSpaceCharge3D.cxx:1289
 AliTPCSpaceCharge3D.cxx:1290
 AliTPCSpaceCharge3D.cxx:1291
 AliTPCSpaceCharge3D.cxx:1292
 AliTPCSpaceCharge3D.cxx:1293
 AliTPCSpaceCharge3D.cxx:1294
 AliTPCSpaceCharge3D.cxx:1295
 AliTPCSpaceCharge3D.cxx:1296
 AliTPCSpaceCharge3D.cxx:1297
 AliTPCSpaceCharge3D.cxx:1298
 AliTPCSpaceCharge3D.cxx:1299
 AliTPCSpaceCharge3D.cxx:1300
 AliTPCSpaceCharge3D.cxx:1301
 AliTPCSpaceCharge3D.cxx:1302
 AliTPCSpaceCharge3D.cxx:1303
 AliTPCSpaceCharge3D.cxx:1304
 AliTPCSpaceCharge3D.cxx:1305
 AliTPCSpaceCharge3D.cxx:1306
 AliTPCSpaceCharge3D.cxx:1307
 AliTPCSpaceCharge3D.cxx:1308
 AliTPCSpaceCharge3D.cxx:1309
 AliTPCSpaceCharge3D.cxx:1310
 AliTPCSpaceCharge3D.cxx:1311
 AliTPCSpaceCharge3D.cxx:1312
 AliTPCSpaceCharge3D.cxx:1313
 AliTPCSpaceCharge3D.cxx:1314
 AliTPCSpaceCharge3D.cxx:1315
 AliTPCSpaceCharge3D.cxx:1316
 AliTPCSpaceCharge3D.cxx:1317
 AliTPCSpaceCharge3D.cxx:1318
 AliTPCSpaceCharge3D.cxx:1319
 AliTPCSpaceCharge3D.cxx:1320
 AliTPCSpaceCharge3D.cxx:1321
 AliTPCSpaceCharge3D.cxx:1322
 AliTPCSpaceCharge3D.cxx:1323
 AliTPCSpaceCharge3D.cxx:1324
 AliTPCSpaceCharge3D.cxx:1325
 AliTPCSpaceCharge3D.cxx:1326
 AliTPCSpaceCharge3D.cxx:1327
 AliTPCSpaceCharge3D.cxx:1328
 AliTPCSpaceCharge3D.cxx:1329
 AliTPCSpaceCharge3D.cxx:1330
 AliTPCSpaceCharge3D.cxx:1331
 AliTPCSpaceCharge3D.cxx:1332
 AliTPCSpaceCharge3D.cxx:1333
 AliTPCSpaceCharge3D.cxx:1334
 AliTPCSpaceCharge3D.cxx:1335
 AliTPCSpaceCharge3D.cxx:1336
 AliTPCSpaceCharge3D.cxx:1337
 AliTPCSpaceCharge3D.cxx:1338
 AliTPCSpaceCharge3D.cxx:1339
 AliTPCSpaceCharge3D.cxx:1340
 AliTPCSpaceCharge3D.cxx:1341
 AliTPCSpaceCharge3D.cxx:1342
 AliTPCSpaceCharge3D.cxx:1343
 AliTPCSpaceCharge3D.cxx:1344
 AliTPCSpaceCharge3D.cxx:1345
 AliTPCSpaceCharge3D.cxx:1346
 AliTPCSpaceCharge3D.cxx:1347
 AliTPCSpaceCharge3D.cxx:1348
 AliTPCSpaceCharge3D.cxx:1349
 AliTPCSpaceCharge3D.cxx:1350
 AliTPCSpaceCharge3D.cxx:1351
 AliTPCSpaceCharge3D.cxx:1352
 AliTPCSpaceCharge3D.cxx:1353
 AliTPCSpaceCharge3D.cxx:1354
 AliTPCSpaceCharge3D.cxx:1355
 AliTPCSpaceCharge3D.cxx:1356
 AliTPCSpaceCharge3D.cxx:1357
 AliTPCSpaceCharge3D.cxx:1358
 AliTPCSpaceCharge3D.cxx:1359
 AliTPCSpaceCharge3D.cxx:1360
 AliTPCSpaceCharge3D.cxx:1361
 AliTPCSpaceCharge3D.cxx:1362
 AliTPCSpaceCharge3D.cxx:1363
 AliTPCSpaceCharge3D.cxx:1364
 AliTPCSpaceCharge3D.cxx:1365
 AliTPCSpaceCharge3D.cxx:1366
 AliTPCSpaceCharge3D.cxx:1367
 AliTPCSpaceCharge3D.cxx:1368
 AliTPCSpaceCharge3D.cxx:1369
 AliTPCSpaceCharge3D.cxx:1370
 AliTPCSpaceCharge3D.cxx:1371
 AliTPCSpaceCharge3D.cxx:1372
 AliTPCSpaceCharge3D.cxx:1373
 AliTPCSpaceCharge3D.cxx:1374
 AliTPCSpaceCharge3D.cxx:1375
 AliTPCSpaceCharge3D.cxx:1376
 AliTPCSpaceCharge3D.cxx:1377
 AliTPCSpaceCharge3D.cxx:1378
 AliTPCSpaceCharge3D.cxx:1379
 AliTPCSpaceCharge3D.cxx:1380
 AliTPCSpaceCharge3D.cxx:1381
 AliTPCSpaceCharge3D.cxx:1382
 AliTPCSpaceCharge3D.cxx:1383
 AliTPCSpaceCharge3D.cxx:1384
 AliTPCSpaceCharge3D.cxx:1385
 AliTPCSpaceCharge3D.cxx:1386
 AliTPCSpaceCharge3D.cxx:1387
 AliTPCSpaceCharge3D.cxx:1388
 AliTPCSpaceCharge3D.cxx:1389
 AliTPCSpaceCharge3D.cxx:1390
 AliTPCSpaceCharge3D.cxx:1391
 AliTPCSpaceCharge3D.cxx:1392
 AliTPCSpaceCharge3D.cxx:1393
 AliTPCSpaceCharge3D.cxx:1394
 AliTPCSpaceCharge3D.cxx:1395
 AliTPCSpaceCharge3D.cxx:1396
 AliTPCSpaceCharge3D.cxx:1397
 AliTPCSpaceCharge3D.cxx:1398
 AliTPCSpaceCharge3D.cxx:1399
 AliTPCSpaceCharge3D.cxx:1400
 AliTPCSpaceCharge3D.cxx:1401
 AliTPCSpaceCharge3D.cxx:1402
 AliTPCSpaceCharge3D.cxx:1403
 AliTPCSpaceCharge3D.cxx:1404
 AliTPCSpaceCharge3D.cxx:1405
 AliTPCSpaceCharge3D.cxx:1406
 AliTPCSpaceCharge3D.cxx:1407
 AliTPCSpaceCharge3D.cxx:1408
 AliTPCSpaceCharge3D.cxx:1409
 AliTPCSpaceCharge3D.cxx:1410
 AliTPCSpaceCharge3D.cxx:1411
 AliTPCSpaceCharge3D.cxx:1412
 AliTPCSpaceCharge3D.cxx:1413
 AliTPCSpaceCharge3D.cxx:1414
 AliTPCSpaceCharge3D.cxx:1415
 AliTPCSpaceCharge3D.cxx:1416
 AliTPCSpaceCharge3D.cxx:1417
 AliTPCSpaceCharge3D.cxx:1418
 AliTPCSpaceCharge3D.cxx:1419
 AliTPCSpaceCharge3D.cxx:1420
 AliTPCSpaceCharge3D.cxx:1421
 AliTPCSpaceCharge3D.cxx:1422
 AliTPCSpaceCharge3D.cxx:1423
 AliTPCSpaceCharge3D.cxx:1424
 AliTPCSpaceCharge3D.cxx:1425
 AliTPCSpaceCharge3D.cxx:1426
 AliTPCSpaceCharge3D.cxx:1427
 AliTPCSpaceCharge3D.cxx:1428
 AliTPCSpaceCharge3D.cxx:1429
 AliTPCSpaceCharge3D.cxx:1430
 AliTPCSpaceCharge3D.cxx:1431
 AliTPCSpaceCharge3D.cxx:1432
 AliTPCSpaceCharge3D.cxx:1433
 AliTPCSpaceCharge3D.cxx:1434
 AliTPCSpaceCharge3D.cxx:1435
 AliTPCSpaceCharge3D.cxx:1436
 AliTPCSpaceCharge3D.cxx:1437
 AliTPCSpaceCharge3D.cxx:1438
 AliTPCSpaceCharge3D.cxx:1439
 AliTPCSpaceCharge3D.cxx:1440
 AliTPCSpaceCharge3D.cxx:1441
 AliTPCSpaceCharge3D.cxx:1442
 AliTPCSpaceCharge3D.cxx:1443
 AliTPCSpaceCharge3D.cxx:1444
 AliTPCSpaceCharge3D.cxx:1445
 AliTPCSpaceCharge3D.cxx:1446
 AliTPCSpaceCharge3D.cxx:1447
 AliTPCSpaceCharge3D.cxx:1448
 AliTPCSpaceCharge3D.cxx:1449
 AliTPCSpaceCharge3D.cxx:1450
 AliTPCSpaceCharge3D.cxx:1451
 AliTPCSpaceCharge3D.cxx:1452
 AliTPCSpaceCharge3D.cxx:1453
 AliTPCSpaceCharge3D.cxx:1454
 AliTPCSpaceCharge3D.cxx:1455
 AliTPCSpaceCharge3D.cxx:1456
 AliTPCSpaceCharge3D.cxx:1457
 AliTPCSpaceCharge3D.cxx:1458
 AliTPCSpaceCharge3D.cxx:1459
 AliTPCSpaceCharge3D.cxx:1460
 AliTPCSpaceCharge3D.cxx:1461
 AliTPCSpaceCharge3D.cxx:1462
 AliTPCSpaceCharge3D.cxx:1463
 AliTPCSpaceCharge3D.cxx:1464
 AliTPCSpaceCharge3D.cxx:1465
 AliTPCSpaceCharge3D.cxx:1466
 AliTPCSpaceCharge3D.cxx:1467
 AliTPCSpaceCharge3D.cxx:1468
 AliTPCSpaceCharge3D.cxx:1469
 AliTPCSpaceCharge3D.cxx:1470
 AliTPCSpaceCharge3D.cxx:1471
 AliTPCSpaceCharge3D.cxx:1472
 AliTPCSpaceCharge3D.cxx:1473
 AliTPCSpaceCharge3D.cxx:1474
 AliTPCSpaceCharge3D.cxx:1475
 AliTPCSpaceCharge3D.cxx:1476
 AliTPCSpaceCharge3D.cxx:1477
 AliTPCSpaceCharge3D.cxx:1478
 AliTPCSpaceCharge3D.cxx:1479
 AliTPCSpaceCharge3D.cxx:1480
 AliTPCSpaceCharge3D.cxx:1481
 AliTPCSpaceCharge3D.cxx:1482
 AliTPCSpaceCharge3D.cxx:1483
 AliTPCSpaceCharge3D.cxx:1484
 AliTPCSpaceCharge3D.cxx:1485
 AliTPCSpaceCharge3D.cxx:1486
 AliTPCSpaceCharge3D.cxx:1487
 AliTPCSpaceCharge3D.cxx:1488
 AliTPCSpaceCharge3D.cxx:1489
 AliTPCSpaceCharge3D.cxx:1490
 AliTPCSpaceCharge3D.cxx:1491
 AliTPCSpaceCharge3D.cxx:1492
 AliTPCSpaceCharge3D.cxx:1493
 AliTPCSpaceCharge3D.cxx:1494
 AliTPCSpaceCharge3D.cxx:1495
 AliTPCSpaceCharge3D.cxx:1496
 AliTPCSpaceCharge3D.cxx:1497
 AliTPCSpaceCharge3D.cxx:1498
 AliTPCSpaceCharge3D.cxx:1499
 AliTPCSpaceCharge3D.cxx:1500
 AliTPCSpaceCharge3D.cxx:1501
 AliTPCSpaceCharge3D.cxx:1502
 AliTPCSpaceCharge3D.cxx:1503
 AliTPCSpaceCharge3D.cxx:1504
 AliTPCSpaceCharge3D.cxx:1505
 AliTPCSpaceCharge3D.cxx:1506
 AliTPCSpaceCharge3D.cxx:1507
 AliTPCSpaceCharge3D.cxx:1508
 AliTPCSpaceCharge3D.cxx:1509
 AliTPCSpaceCharge3D.cxx:1510
 AliTPCSpaceCharge3D.cxx:1511
 AliTPCSpaceCharge3D.cxx:1512
 AliTPCSpaceCharge3D.cxx:1513
 AliTPCSpaceCharge3D.cxx:1514
 AliTPCSpaceCharge3D.cxx:1515
 AliTPCSpaceCharge3D.cxx:1516
 AliTPCSpaceCharge3D.cxx:1517
 AliTPCSpaceCharge3D.cxx:1518
 AliTPCSpaceCharge3D.cxx:1519
 AliTPCSpaceCharge3D.cxx:1520
 AliTPCSpaceCharge3D.cxx:1521
 AliTPCSpaceCharge3D.cxx:1522
 AliTPCSpaceCharge3D.cxx:1523
 AliTPCSpaceCharge3D.cxx:1524
 AliTPCSpaceCharge3D.cxx:1525
 AliTPCSpaceCharge3D.cxx:1526
 AliTPCSpaceCharge3D.cxx:1527
 AliTPCSpaceCharge3D.cxx:1528
 AliTPCSpaceCharge3D.cxx:1529
 AliTPCSpaceCharge3D.cxx:1530
 AliTPCSpaceCharge3D.cxx:1531
 AliTPCSpaceCharge3D.cxx:1532
 AliTPCSpaceCharge3D.cxx:1533
 AliTPCSpaceCharge3D.cxx:1534
 AliTPCSpaceCharge3D.cxx:1535
 AliTPCSpaceCharge3D.cxx:1536
 AliTPCSpaceCharge3D.cxx:1537
 AliTPCSpaceCharge3D.cxx:1538
 AliTPCSpaceCharge3D.cxx:1539
 AliTPCSpaceCharge3D.cxx:1540
 AliTPCSpaceCharge3D.cxx:1541
 AliTPCSpaceCharge3D.cxx:1542
 AliTPCSpaceCharge3D.cxx:1543
 AliTPCSpaceCharge3D.cxx:1544
 AliTPCSpaceCharge3D.cxx:1545
 AliTPCSpaceCharge3D.cxx:1546
 AliTPCSpaceCharge3D.cxx:1547
 AliTPCSpaceCharge3D.cxx:1548
 AliTPCSpaceCharge3D.cxx:1549
 AliTPCSpaceCharge3D.cxx:1550
 AliTPCSpaceCharge3D.cxx:1551
 AliTPCSpaceCharge3D.cxx:1552
 AliTPCSpaceCharge3D.cxx:1553
 AliTPCSpaceCharge3D.cxx:1554
 AliTPCSpaceCharge3D.cxx:1555
 AliTPCSpaceCharge3D.cxx:1556
 AliTPCSpaceCharge3D.cxx:1557
 AliTPCSpaceCharge3D.cxx:1558
 AliTPCSpaceCharge3D.cxx:1559
 AliTPCSpaceCharge3D.cxx:1560
 AliTPCSpaceCharge3D.cxx:1561
 AliTPCSpaceCharge3D.cxx:1562
 AliTPCSpaceCharge3D.cxx:1563
 AliTPCSpaceCharge3D.cxx:1564
 AliTPCSpaceCharge3D.cxx:1565
 AliTPCSpaceCharge3D.cxx:1566
 AliTPCSpaceCharge3D.cxx:1567
 AliTPCSpaceCharge3D.cxx:1568
 AliTPCSpaceCharge3D.cxx:1569
 AliTPCSpaceCharge3D.cxx:1570
 AliTPCSpaceCharge3D.cxx:1571
 AliTPCSpaceCharge3D.cxx:1572
 AliTPCSpaceCharge3D.cxx:1573
 AliTPCSpaceCharge3D.cxx:1574
 AliTPCSpaceCharge3D.cxx:1575
 AliTPCSpaceCharge3D.cxx:1576
 AliTPCSpaceCharge3D.cxx:1577
 AliTPCSpaceCharge3D.cxx:1578
 AliTPCSpaceCharge3D.cxx:1579
 AliTPCSpaceCharge3D.cxx:1580
 AliTPCSpaceCharge3D.cxx:1581
 AliTPCSpaceCharge3D.cxx:1582
 AliTPCSpaceCharge3D.cxx:1583
 AliTPCSpaceCharge3D.cxx:1584
 AliTPCSpaceCharge3D.cxx:1585
 AliTPCSpaceCharge3D.cxx:1586
 AliTPCSpaceCharge3D.cxx:1587
 AliTPCSpaceCharge3D.cxx:1588
 AliTPCSpaceCharge3D.cxx:1589
 AliTPCSpaceCharge3D.cxx:1590
 AliTPCSpaceCharge3D.cxx:1591
 AliTPCSpaceCharge3D.cxx:1592
 AliTPCSpaceCharge3D.cxx:1593
 AliTPCSpaceCharge3D.cxx:1594
 AliTPCSpaceCharge3D.cxx:1595
 AliTPCSpaceCharge3D.cxx:1596
 AliTPCSpaceCharge3D.cxx:1597
 AliTPCSpaceCharge3D.cxx:1598
 AliTPCSpaceCharge3D.cxx:1599
 AliTPCSpaceCharge3D.cxx:1600
 AliTPCSpaceCharge3D.cxx:1601
 AliTPCSpaceCharge3D.cxx:1602
 AliTPCSpaceCharge3D.cxx:1603
 AliTPCSpaceCharge3D.cxx:1604
 AliTPCSpaceCharge3D.cxx:1605
 AliTPCSpaceCharge3D.cxx:1606
 AliTPCSpaceCharge3D.cxx:1607
 AliTPCSpaceCharge3D.cxx:1608
 AliTPCSpaceCharge3D.cxx:1609
 AliTPCSpaceCharge3D.cxx:1610
 AliTPCSpaceCharge3D.cxx:1611
 AliTPCSpaceCharge3D.cxx:1612
 AliTPCSpaceCharge3D.cxx:1613
 AliTPCSpaceCharge3D.cxx:1614
 AliTPCSpaceCharge3D.cxx:1615
 AliTPCSpaceCharge3D.cxx:1616
 AliTPCSpaceCharge3D.cxx:1617
 AliTPCSpaceCharge3D.cxx:1618
 AliTPCSpaceCharge3D.cxx:1619
 AliTPCSpaceCharge3D.cxx:1620
 AliTPCSpaceCharge3D.cxx:1621
 AliTPCSpaceCharge3D.cxx:1622
 AliTPCSpaceCharge3D.cxx:1623
 AliTPCSpaceCharge3D.cxx:1624
 AliTPCSpaceCharge3D.cxx:1625
 AliTPCSpaceCharge3D.cxx:1626
 AliTPCSpaceCharge3D.cxx:1627
 AliTPCSpaceCharge3D.cxx:1628
 AliTPCSpaceCharge3D.cxx:1629
 AliTPCSpaceCharge3D.cxx:1630
 AliTPCSpaceCharge3D.cxx:1631
 AliTPCSpaceCharge3D.cxx:1632
 AliTPCSpaceCharge3D.cxx:1633
 AliTPCSpaceCharge3D.cxx:1634
 AliTPCSpaceCharge3D.cxx:1635
 AliTPCSpaceCharge3D.cxx:1636
 AliTPCSpaceCharge3D.cxx:1637
 AliTPCSpaceCharge3D.cxx:1638
 AliTPCSpaceCharge3D.cxx:1639
 AliTPCSpaceCharge3D.cxx:1640
 AliTPCSpaceCharge3D.cxx:1641
 AliTPCSpaceCharge3D.cxx:1642
 AliTPCSpaceCharge3D.cxx:1643
 AliTPCSpaceCharge3D.cxx:1644
 AliTPCSpaceCharge3D.cxx:1645
 AliTPCSpaceCharge3D.cxx:1646
 AliTPCSpaceCharge3D.cxx:1647
 AliTPCSpaceCharge3D.cxx:1648
 AliTPCSpaceCharge3D.cxx:1649
 AliTPCSpaceCharge3D.cxx:1650
 AliTPCSpaceCharge3D.cxx:1651
 AliTPCSpaceCharge3D.cxx:1652
 AliTPCSpaceCharge3D.cxx:1653
 AliTPCSpaceCharge3D.cxx:1654
 AliTPCSpaceCharge3D.cxx:1655
 AliTPCSpaceCharge3D.cxx:1656
 AliTPCSpaceCharge3D.cxx:1657
 AliTPCSpaceCharge3D.cxx:1658
 AliTPCSpaceCharge3D.cxx:1659
 AliTPCSpaceCharge3D.cxx:1660
 AliTPCSpaceCharge3D.cxx:1661
 AliTPCSpaceCharge3D.cxx:1662
 AliTPCSpaceCharge3D.cxx:1663
 AliTPCSpaceCharge3D.cxx:1664
 AliTPCSpaceCharge3D.cxx:1665
 AliTPCSpaceCharge3D.cxx:1666
 AliTPCSpaceCharge3D.cxx:1667
 AliTPCSpaceCharge3D.cxx:1668
 AliTPCSpaceCharge3D.cxx:1669
 AliTPCSpaceCharge3D.cxx:1670
 AliTPCSpaceCharge3D.cxx:1671
 AliTPCSpaceCharge3D.cxx:1672
 AliTPCSpaceCharge3D.cxx:1673
 AliTPCSpaceCharge3D.cxx:1674
 AliTPCSpaceCharge3D.cxx:1675
 AliTPCSpaceCharge3D.cxx:1676
 AliTPCSpaceCharge3D.cxx:1677
 AliTPCSpaceCharge3D.cxx:1678
 AliTPCSpaceCharge3D.cxx:1679
 AliTPCSpaceCharge3D.cxx:1680
 AliTPCSpaceCharge3D.cxx:1681
 AliTPCSpaceCharge3D.cxx:1682
 AliTPCSpaceCharge3D.cxx:1683
 AliTPCSpaceCharge3D.cxx:1684
 AliTPCSpaceCharge3D.cxx:1685
 AliTPCSpaceCharge3D.cxx:1686
 AliTPCSpaceCharge3D.cxx:1687
 AliTPCSpaceCharge3D.cxx:1688
 AliTPCSpaceCharge3D.cxx:1689
 AliTPCSpaceCharge3D.cxx:1690
 AliTPCSpaceCharge3D.cxx:1691
 AliTPCSpaceCharge3D.cxx:1692
 AliTPCSpaceCharge3D.cxx:1693
 AliTPCSpaceCharge3D.cxx:1694
 AliTPCSpaceCharge3D.cxx:1695
 AliTPCSpaceCharge3D.cxx:1696
 AliTPCSpaceCharge3D.cxx:1697
 AliTPCSpaceCharge3D.cxx:1698
 AliTPCSpaceCharge3D.cxx:1699
 AliTPCSpaceCharge3D.cxx:1700
 AliTPCSpaceCharge3D.cxx:1701
 AliTPCSpaceCharge3D.cxx:1702
 AliTPCSpaceCharge3D.cxx:1703
 AliTPCSpaceCharge3D.cxx:1704
 AliTPCSpaceCharge3D.cxx:1705
 AliTPCSpaceCharge3D.cxx:1706
 AliTPCSpaceCharge3D.cxx:1707
 AliTPCSpaceCharge3D.cxx:1708
 AliTPCSpaceCharge3D.cxx:1709
 AliTPCSpaceCharge3D.cxx:1710
 AliTPCSpaceCharge3D.cxx:1711
 AliTPCSpaceCharge3D.cxx:1712
 AliTPCSpaceCharge3D.cxx:1713
 AliTPCSpaceCharge3D.cxx:1714
 AliTPCSpaceCharge3D.cxx:1715
 AliTPCSpaceCharge3D.cxx:1716
 AliTPCSpaceCharge3D.cxx:1717
 AliTPCSpaceCharge3D.cxx:1718
 AliTPCSpaceCharge3D.cxx:1719
 AliTPCSpaceCharge3D.cxx:1720
 AliTPCSpaceCharge3D.cxx:1721
 AliTPCSpaceCharge3D.cxx:1722
 AliTPCSpaceCharge3D.cxx:1723
 AliTPCSpaceCharge3D.cxx:1724
 AliTPCSpaceCharge3D.cxx:1725
 AliTPCSpaceCharge3D.cxx:1726
 AliTPCSpaceCharge3D.cxx:1727
 AliTPCSpaceCharge3D.cxx:1728
 AliTPCSpaceCharge3D.cxx:1729
 AliTPCSpaceCharge3D.cxx:1730
 AliTPCSpaceCharge3D.cxx:1731
 AliTPCSpaceCharge3D.cxx:1732
 AliTPCSpaceCharge3D.cxx:1733
 AliTPCSpaceCharge3D.cxx:1734
 AliTPCSpaceCharge3D.cxx:1735
 AliTPCSpaceCharge3D.cxx:1736
 AliTPCSpaceCharge3D.cxx:1737
 AliTPCSpaceCharge3D.cxx:1738
 AliTPCSpaceCharge3D.cxx:1739
 AliTPCSpaceCharge3D.cxx:1740
 AliTPCSpaceCharge3D.cxx:1741
 AliTPCSpaceCharge3D.cxx:1742
 AliTPCSpaceCharge3D.cxx:1743
 AliTPCSpaceCharge3D.cxx:1744
 AliTPCSpaceCharge3D.cxx:1745
 AliTPCSpaceCharge3D.cxx:1746
 AliTPCSpaceCharge3D.cxx:1747
 AliTPCSpaceCharge3D.cxx:1748
 AliTPCSpaceCharge3D.cxx:1749
 AliTPCSpaceCharge3D.cxx:1750
 AliTPCSpaceCharge3D.cxx:1751
 AliTPCSpaceCharge3D.cxx:1752
 AliTPCSpaceCharge3D.cxx:1753
 AliTPCSpaceCharge3D.cxx:1754
 AliTPCSpaceCharge3D.cxx:1755
 AliTPCSpaceCharge3D.cxx:1756
 AliTPCSpaceCharge3D.cxx:1757
 AliTPCSpaceCharge3D.cxx:1758
 AliTPCSpaceCharge3D.cxx:1759
 AliTPCSpaceCharge3D.cxx:1760
 AliTPCSpaceCharge3D.cxx:1761
 AliTPCSpaceCharge3D.cxx:1762
 AliTPCSpaceCharge3D.cxx:1763
 AliTPCSpaceCharge3D.cxx:1764
 AliTPCSpaceCharge3D.cxx:1765
 AliTPCSpaceCharge3D.cxx:1766
 AliTPCSpaceCharge3D.cxx:1767
 AliTPCSpaceCharge3D.cxx:1768
 AliTPCSpaceCharge3D.cxx:1769
 AliTPCSpaceCharge3D.cxx:1770
 AliTPCSpaceCharge3D.cxx:1771
 AliTPCSpaceCharge3D.cxx:1772
 AliTPCSpaceCharge3D.cxx:1773
 AliTPCSpaceCharge3D.cxx:1774
 AliTPCSpaceCharge3D.cxx:1775
 AliTPCSpaceCharge3D.cxx:1776
 AliTPCSpaceCharge3D.cxx:1777
 AliTPCSpaceCharge3D.cxx:1778
 AliTPCSpaceCharge3D.cxx:1779
 AliTPCSpaceCharge3D.cxx:1780
 AliTPCSpaceCharge3D.cxx:1781
 AliTPCSpaceCharge3D.cxx:1782
 AliTPCSpaceCharge3D.cxx:1783
 AliTPCSpaceCharge3D.cxx:1784
 AliTPCSpaceCharge3D.cxx:1785
 AliTPCSpaceCharge3D.cxx:1786
 AliTPCSpaceCharge3D.cxx:1787
 AliTPCSpaceCharge3D.cxx:1788
 AliTPCSpaceCharge3D.cxx:1789
 AliTPCSpaceCharge3D.cxx:1790
 AliTPCSpaceCharge3D.cxx:1791
 AliTPCSpaceCharge3D.cxx:1792
 AliTPCSpaceCharge3D.cxx:1793
 AliTPCSpaceCharge3D.cxx:1794
 AliTPCSpaceCharge3D.cxx:1795
 AliTPCSpaceCharge3D.cxx:1796
 AliTPCSpaceCharge3D.cxx:1797
 AliTPCSpaceCharge3D.cxx:1798
 AliTPCSpaceCharge3D.cxx:1799
 AliTPCSpaceCharge3D.cxx:1800
 AliTPCSpaceCharge3D.cxx:1801
 AliTPCSpaceCharge3D.cxx:1802
 AliTPCSpaceCharge3D.cxx:1803
 AliTPCSpaceCharge3D.cxx:1804
 AliTPCSpaceCharge3D.cxx:1805
 AliTPCSpaceCharge3D.cxx:1806
 AliTPCSpaceCharge3D.cxx:1807
 AliTPCSpaceCharge3D.cxx:1808
 AliTPCSpaceCharge3D.cxx:1809
 AliTPCSpaceCharge3D.cxx:1810
 AliTPCSpaceCharge3D.cxx:1811
 AliTPCSpaceCharge3D.cxx:1812
 AliTPCSpaceCharge3D.cxx:1813
 AliTPCSpaceCharge3D.cxx:1814
 AliTPCSpaceCharge3D.cxx:1815
 AliTPCSpaceCharge3D.cxx:1816
 AliTPCSpaceCharge3D.cxx:1817
 AliTPCSpaceCharge3D.cxx:1818
 AliTPCSpaceCharge3D.cxx:1819
 AliTPCSpaceCharge3D.cxx:1820
 AliTPCSpaceCharge3D.cxx:1821
 AliTPCSpaceCharge3D.cxx:1822
 AliTPCSpaceCharge3D.cxx:1823
 AliTPCSpaceCharge3D.cxx:1824
 AliTPCSpaceCharge3D.cxx:1825
 AliTPCSpaceCharge3D.cxx:1826
 AliTPCSpaceCharge3D.cxx:1827
 AliTPCSpaceCharge3D.cxx:1828
 AliTPCSpaceCharge3D.cxx:1829
 AliTPCSpaceCharge3D.cxx:1830
 AliTPCSpaceCharge3D.cxx:1831
 AliTPCSpaceCharge3D.cxx:1832
 AliTPCSpaceCharge3D.cxx:1833
 AliTPCSpaceCharge3D.cxx:1834
 AliTPCSpaceCharge3D.cxx:1835
 AliTPCSpaceCharge3D.cxx:1836
 AliTPCSpaceCharge3D.cxx:1837
 AliTPCSpaceCharge3D.cxx:1838
 AliTPCSpaceCharge3D.cxx:1839
 AliTPCSpaceCharge3D.cxx:1840
 AliTPCSpaceCharge3D.cxx:1841
 AliTPCSpaceCharge3D.cxx:1842
 AliTPCSpaceCharge3D.cxx:1843
 AliTPCSpaceCharge3D.cxx:1844
 AliTPCSpaceCharge3D.cxx:1845
 AliTPCSpaceCharge3D.cxx:1846
 AliTPCSpaceCharge3D.cxx:1847
 AliTPCSpaceCharge3D.cxx:1848
 AliTPCSpaceCharge3D.cxx:1849
 AliTPCSpaceCharge3D.cxx:1850
 AliTPCSpaceCharge3D.cxx:1851
 AliTPCSpaceCharge3D.cxx:1852
 AliTPCSpaceCharge3D.cxx:1853
 AliTPCSpaceCharge3D.cxx:1854
 AliTPCSpaceCharge3D.cxx:1855
 AliTPCSpaceCharge3D.cxx:1856
 AliTPCSpaceCharge3D.cxx:1857
 AliTPCSpaceCharge3D.cxx:1858
 AliTPCSpaceCharge3D.cxx:1859
 AliTPCSpaceCharge3D.cxx:1860
 AliTPCSpaceCharge3D.cxx:1861
 AliTPCSpaceCharge3D.cxx:1862
 AliTPCSpaceCharge3D.cxx:1863
 AliTPCSpaceCharge3D.cxx:1864
 AliTPCSpaceCharge3D.cxx:1865
 AliTPCSpaceCharge3D.cxx:1866
 AliTPCSpaceCharge3D.cxx:1867
 AliTPCSpaceCharge3D.cxx:1868
 AliTPCSpaceCharge3D.cxx:1869
 AliTPCSpaceCharge3D.cxx:1870
 AliTPCSpaceCharge3D.cxx:1871
 AliTPCSpaceCharge3D.cxx:1872
 AliTPCSpaceCharge3D.cxx:1873
 AliTPCSpaceCharge3D.cxx:1874
 AliTPCSpaceCharge3D.cxx:1875
 AliTPCSpaceCharge3D.cxx:1876
 AliTPCSpaceCharge3D.cxx:1877
 AliTPCSpaceCharge3D.cxx:1878
 AliTPCSpaceCharge3D.cxx:1879
 AliTPCSpaceCharge3D.cxx:1880
 AliTPCSpaceCharge3D.cxx:1881
 AliTPCSpaceCharge3D.cxx:1882
 AliTPCSpaceCharge3D.cxx:1883
 AliTPCSpaceCharge3D.cxx:1884
 AliTPCSpaceCharge3D.cxx:1885
 AliTPCSpaceCharge3D.cxx:1886
 AliTPCSpaceCharge3D.cxx:1887
 AliTPCSpaceCharge3D.cxx:1888
 AliTPCSpaceCharge3D.cxx:1889
 AliTPCSpaceCharge3D.cxx:1890
 AliTPCSpaceCharge3D.cxx:1891
 AliTPCSpaceCharge3D.cxx:1892
 AliTPCSpaceCharge3D.cxx:1893
 AliTPCSpaceCharge3D.cxx:1894
 AliTPCSpaceCharge3D.cxx:1895
 AliTPCSpaceCharge3D.cxx:1896
 AliTPCSpaceCharge3D.cxx:1897
 AliTPCSpaceCharge3D.cxx:1898
 AliTPCSpaceCharge3D.cxx:1899
 AliTPCSpaceCharge3D.cxx:1900
 AliTPCSpaceCharge3D.cxx:1901
 AliTPCSpaceCharge3D.cxx:1902
 AliTPCSpaceCharge3D.cxx:1903
 AliTPCSpaceCharge3D.cxx:1904
 AliTPCSpaceCharge3D.cxx:1905
 AliTPCSpaceCharge3D.cxx:1906
 AliTPCSpaceCharge3D.cxx:1907
 AliTPCSpaceCharge3D.cxx:1908
 AliTPCSpaceCharge3D.cxx:1909
 AliTPCSpaceCharge3D.cxx:1910
 AliTPCSpaceCharge3D.cxx:1911
 AliTPCSpaceCharge3D.cxx:1912
 AliTPCSpaceCharge3D.cxx:1913
 AliTPCSpaceCharge3D.cxx:1914
 AliTPCSpaceCharge3D.cxx:1915
 AliTPCSpaceCharge3D.cxx:1916
 AliTPCSpaceCharge3D.cxx:1917
 AliTPCSpaceCharge3D.cxx:1918
 AliTPCSpaceCharge3D.cxx:1919
 AliTPCSpaceCharge3D.cxx:1920
 AliTPCSpaceCharge3D.cxx:1921
 AliTPCSpaceCharge3D.cxx:1922
 AliTPCSpaceCharge3D.cxx:1923
 AliTPCSpaceCharge3D.cxx:1924
 AliTPCSpaceCharge3D.cxx:1925
 AliTPCSpaceCharge3D.cxx:1926
 AliTPCSpaceCharge3D.cxx:1927
 AliTPCSpaceCharge3D.cxx:1928
 AliTPCSpaceCharge3D.cxx:1929
 AliTPCSpaceCharge3D.cxx:1930
 AliTPCSpaceCharge3D.cxx:1931
 AliTPCSpaceCharge3D.cxx:1932
 AliTPCSpaceCharge3D.cxx:1933
 AliTPCSpaceCharge3D.cxx:1934
 AliTPCSpaceCharge3D.cxx:1935
 AliTPCSpaceCharge3D.cxx:1936
 AliTPCSpaceCharge3D.cxx:1937
 AliTPCSpaceCharge3D.cxx:1938
 AliTPCSpaceCharge3D.cxx:1939
 AliTPCSpaceCharge3D.cxx:1940