AliRoot Core  v5-06-30 (35d6c57)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliTPCComposedCorrection Class Reference

AliTPCComposedCorrection class. More...

#include <AliTPCComposedCorrection.h>

Inheritance diagram for AliTPCComposedCorrection:

Public Types

enum  CompositionType { kParallel, kQueue, kQueueResidual }
 
- Public Types inherited from AliTPCCorrection
enum  CompositionType { kParallel, kQueue }
 

Public Member Functions

 AliTPCComposedCorrection ()
 
 AliTPCComposedCorrection (TCollection *corrections, CompositionType mode)
 
virtual ~AliTPCComposedCorrection ()
 
virtual Bool_t AddCorrectionCompact (AliTPCCorrection *corr, Double_t weight)
 
void SetOmegaTauT1T2 (Float_t omegaTau, Float_t t1, Float_t t2)
 
TCollection * GetCorrections () const
 
void SetCorrections (const TCollection *corrections)
 
CompositionType GetMode () const
 
void SetMode (CompositionType mode)
 
virtual void GetCorrection (const Float_t x[], const Short_t roc, Float_t dx[])
 
virtual void GetDistortion (const Float_t x[], const Short_t roc, Float_t dx[])
 
virtual AliTPCCorrectionGetSubCorrection (Int_t ipos)
 
virtual AliTPCCorrectionGetSubCorrection (const char *cname)
 
virtual void Print (Option_t *option="") const
 
virtual void Init ()
 
virtual void Update (const TTimeStamp &timeStamp)
 
void SetWeights (TVectorD *weights)
 
const TVectorD * GetWeights () const
 
- Public Member Functions inherited from AliTPCCorrection
 AliTPCCorrection ()
 
 AliTPCCorrection (const char *name, const char *title)
 
virtual ~AliTPCCorrection ()
 
void CorrectPoint (Float_t x[], Short_t roc)
 
void CorrectPointLocal (Float_t x[], Short_t roc)
 
void CorrectPoint (const Float_t x[], Short_t roc, Float_t xp[])
 
virtual void GetCorrectionDz (const Float_t x[], Short_t roc, Float_t dx[], Float_t delta)
 
virtual void GetCorrectionIntegralDz (const Float_t x[], Short_t roc, Float_t dx[], Float_t delta)
 
void DistortPoint (Float_t x[], Short_t roc)
 
void DistortPointLocal (Float_t x[], Short_t roc)
 
void DistortPoint (const Float_t x[], Short_t roc, Float_t xp[])
 
virtual void GetDistortionDz (const Float_t x[], Short_t roc, Float_t dx[], Float_t delta)
 
virtual void GetDistortionIntegralDz (const Float_t x[], Short_t roc, Float_t dx[], Float_t delta)
 
virtual void SetCorrScaleFactor (Float_t)
 
virtual Float_t GetCorrScaleFactor () const
 
TH2F * CreateHistoDRinXY (Float_t z=10., Int_t nx=100, Int_t ny=100)
 
TH2F * CreateHistoDRPhiinXY (Float_t z=10., Int_t nx=100, Int_t nphi=100)
 
TH2F * CreateHistoDZinXY (Float_t z=10., Int_t nx=100, Int_t ny=100)
 
TH2F * CreateHistoDRinZR (Float_t phi=0., Int_t nZ=100, Int_t nR=100)
 
TH2F * CreateHistoDRPhiinZR (Float_t phi=0., Int_t nZ=100, Int_t nR=100)
 
TH2F * CreateHistoDZinZR (Float_t phi=0., Int_t nZ=100, Int_t nR=100)
 
TTree * CreateDistortionTree (Double_t step=5)
 
AliExternalTrackParam * FitDistortedTrack (AliExternalTrackParam &trackIn, Double_t refX, Int_t dir, TTreeSRedirector *pcstream)
 
void StoreInOCDB (Int_t startRun, Int_t endRun, const char *comment=0)
 
void FastSimDistortedVertex (Double_t orgVertex[3], Int_t nTracks, AliESDVertex &aV, AliESDVertex &avOrg, AliESDVertex &cV, AliESDVertex &cvOrg, TTreeSRedirector *const pcstream, Double_t etaCuts)
 

Private Member Functions

AliTPCComposedCorrectionoperator= (const AliTPCComposedCorrection &)
 
 AliTPCComposedCorrection (const AliTPCComposedCorrection &)
 

Private Attributes

TCollection * fCorrections
 The corrections this one is composed of. More...
 
CompositionType fMode
 The way to apply the corrections (see general class documentation) More...
 
TVectorD * fWeights
 optional vector with weights - used for fit benchmarking More...
 

Additional Inherited Members

- Static Public Member Functions inherited from AliTPCCorrection
static void MakeDistortionMap (THnSparse *his0, TTreeSRedirector *pcstream, const char *hname, Int_t run, Float_t refX, Int_t type, Int_t integ=1)
 
static void MakeDistortionMapCosmic (THnSparse *his0, TTreeSRedirector *pcstream, const char *hname, Int_t run, Float_t refX, Int_t type)
 
static void MakeDistortionMapSector (THnSparse *his0, TTreeSRedirector *pcstream, const char *hname, Int_t run, Int_t type)
 
static void MakeTrackDistortionTree (TTree *tinput, Int_t dtype, Int_t ptype, const TObjArray *corrArray, Int_t step=1, Int_t offset=0, Bool_t debug=0)
 
static void MakeSectorDistortionTree (TTree *tinput, Int_t dtype, Int_t ptype, const TObjArray *corrArray, Int_t step=1, Int_t offset=0, Bool_t debug=0)
 
static void MakeLaserDistortionTreeOld (TTree *tree, TObjArray *corrArray, Int_t itype)
 
static void MakeLaserDistortionTree (TTree *tree, TObjArray *corrArray, Int_t itype)
 
static void AddVisualCorrection (AliTPCCorrection *corr, Int_t position)
 
static AliTPCCorrectionGetVisualCorrection (Int_t position)
 
static AliTPCCorrectionGetVisualCorrection (const char *corName)
 
static TObjArrayGetVisualCorrections ()
 
static Double_t GetCorrSector (Double_t sector, Double_t r, Double_t kZ, Int_t axisType, Int_t corrType=0)
 
static Double_t GetCorrXYZ (Double_t gx, Double_t gy, Double_t gz, Int_t axisType, Int_t corrType=0)
 
static Double_t GetCorrXYZDz (Double_t gx, Double_t gy, Double_t gz, Int_t axisType, Int_t corrType=0, Double_t delta=5)
 
static Double_t GetCorrXYZIntegrateZ (Double_t gx, Double_t gy, Double_t gz, Int_t axisType, Int_t corrType=0, Double_t delta=5)
 
static Double_t GetDistXYZ (Double_t gx, Double_t gy, Double_t gz, Int_t axisType, Int_t corrType=0)
 
static Double_t GetDistXYZDz (Double_t gx, Double_t gy, Double_t gz, Int_t axisType, Int_t corrType=0, Double_t delta=5)
 
static Double_t GetDistXYZIntegrateZ (Double_t gx, Double_t gy, Double_t gz, Int_t axisType, Int_t corrType=0, Double_t delta=5)
 
- Protected Types inherited from AliTPCCorrection
enum  { kNR = 72 }
 
enum  { kNPhi = 18*10+1 }
 
enum  { kNZ = 166 }
 
- Protected Member Functions inherited from AliTPCCorrection
TH2F * CreateTH2F (const char *name, const char *title, const char *xlabel, const char *ylabel, const char *zlabel, Int_t nbinsx, Double_t xlow, Double_t xup, Int_t nbinsy, Double_t ylow, Double_t yup)
 
void Interpolate2DEdistortion (Int_t order, Double_t r, Double_t z, const Double_t er[kNZ][kNR], Double_t &erValue)
 
void Interpolate3DEdistortion (Int_t order, Double_t r, Float_t phi, Double_t z, const Double_t er[kNZ][kNPhi][kNR], const Double_t ephi[kNZ][kNPhi][kNR], const Double_t ez[kNZ][kNPhi][kNR], Double_t &erValue, Double_t &ephiValue, Double_t &ezValue)
 
Double_t Interpolate2DTable (Int_t order, Double_t x, Double_t y, Int_t nx, Int_t ny, const Double_t xv[], const Double_t yv[], const TMatrixD &array)
 
Double_t Interpolate3DTable (Int_t order, Double_t x, Double_t y, Double_t z, Int_t nx, Int_t ny, Int_t nz, const Double_t xv[], const Double_t yv[], const Double_t zv[], TMatrixD **arrayofArrays)
 
Double_t Interpolate (const Double_t xArray[], const Double_t yArray[], Int_t order, Double_t x)
 
void Search (Int_t n, const Double_t xArray[], Double_t x, Int_t &low)
 
Float_t Interpolate2DTable (Int_t order, Double_t x, Double_t y, Int_t nx, Int_t ny, const Double_t xv[], const Double_t yv[], const TMatrixF &array)
 
Float_t Interpolate3DTable (Int_t order, Double_t x, Double_t y, Double_t z, Int_t nx, Int_t ny, Int_t nz, const Double_t xv[], const Double_t yv[], const Double_t zv[], TMatrixF **arrayofArrays)
 
Float_t Interpolate (const Double_t xArray[], const Float_t yArray[], Int_t order, Double_t x)
 
virtual Int_t IsPowerOfTwo (Int_t i) const
 
void PoissonRelaxation2D (TMatrixD &arrayV, TMatrixD &chargeDensity, TMatrixD &arrayErOverEz, TMatrixD &arrayDeltaEz, Int_t rows, Int_t columns, Int_t iterations, Bool_t rocDisplacement=kTRUE)
 
void PoissonRelaxation3D (TMatrixD **arrayofArrayV, TMatrixD **arrayofChargeDensities, TMatrixD **arrayofEroverEz, TMatrixD **arrayofEPhioverEz, TMatrixD **arrayofEz, Int_t rows, Int_t columns, Int_t phislices, Float_t deltaphi, Int_t iterations, Int_t summetry, Bool_t rocDisplacement=kTRUE)
 
void SetIsLocal (Bool_t isLocal)
 
Bool_t IsLocal () const
 
- Protected Attributes inherited from AliTPCCorrection
Double_t fgkRList [kNR]
 points in the radial direction (for the lookup table) More...
 
Double_t fgkPhiList [kNPhi]
 points in the phi direction (for the lookup table) More...
 
Double_t fgkZList [kNZ]
 points in the z direction (for the lookup table) More...
 
Int_t fILow
 
Int_t fJLow
 
Int_t fKLow
 variable to help in the interpolation More...
 
Double_t fT1
 tensor term of wt - T1 More...
 
Double_t fT2
 tensor term of wt - T2 More...
 
Bool_t fIsLocal
 switch to indicate that the distortion is a local vector drphi/dz, dr/dz More...
 
- Static Protected Attributes inherited from AliTPCCorrection
static const Double_t fgkTPCZ0 = 249.7
 nominal gating grid position More...
 
static const Double_t fgkIFCRadius = 83.5
 Mean Radius of the Inner Field Cage ( 82.43 min, 83.70 max) (cm) More...
 
static const Double_t fgkOFCRadius = 254.5
 Mean Radius of the Outer Field Cage (252.55 min, 256.45 max) (cm) More...
 
static const Double_t fgkZOffSet = 0.2
 Offset from CE: calculate all distortions closer to CE as if at this point. More...
 
static const Double_t fgkCathodeV = -100000.0
 Cathode Voltage (volts) More...
 
static const Double_t fgkGG = -70.0
 Gating Grid voltage (volts) More...
 
static const Double_t fgkdvdE = 0.0024
 [cm/V] drift velocity dependency on the E field (from Magboltz for NeCO2N2 at standard environment) More...
 
static const Double_t fgkEM = -1.602176487e-19/9.10938215e-31
 charge/mass in [C/kg] More...
 
static const Double_t fgke0 = 8.854187817e-12
 vacuum permittivity [A·s/(V·m)] More...
 
static TObjArrayfgVisualCorrection =0
 array of orrection for visualization More...
 

Detailed Description

AliTPCComposedCorrection class.

This class is creating a correction that is composed out of smaller corrections. There are two ways the sub-corrections can be combined into this one:

  1. kParallel: All corrections are applied at the given position x and the dx terms are summed up (this commutes).
  2. kQueue: The corrections are called in order. The first one at the given position x resulting in dx1, the second one is called at the corrected position (x+dx1) resulting in dx2, the third one is then called at position (x+dx1+dx2) and so forth. dx=dx1+dx2+... is returned.
  3. kQueueResidual: like kQueue with the exception that in case of a positive weight the 'Distortion' is called and in case of a negative weight the 'Correction' is called, where the absolute of the weight will be applied to the correction For the inverse of the correction this is taken into account by reversing the order the corrections are applied in the kQueue case (no issue for kParallel).

Example usage:

AliMagF mag("mag","mag");
AliTPCExBBShape exb; // B field shape distortions
exb.SetBField(&mag);
AliTPCExBTwist twist; // ExB Twist distortions
twist.SetXTwist(0.001);
TObjArray cs; cs.Add(&exb); cs.Add(&twist);
cc.SetOmegaTauT1T2(wt,T1,T2);
cc.Print("DA");
cc.CreateHistoDRPhiinZR(0,100,100)->Draw("surf2");
Author
Magnus Mager, Stefan Rossegger, Jim Thomas
Date
27/04/2010

This class is creating a correction that is composed out of smaller corrections. There are two ways the sub-corrections can be combined into this one:

  1. kParallel: All corrections are applied at the given position x and the dx terms are summed up (this commutes).
  2. kQueue: The corrections are called in order. The first one at the given position x resulting in dx1, the second one is called at the corrected position (x+dx1) resulting in dx2, the third one is then called at position (x+dx1+dx2) and so forth. dx=dx1+dx2+... is returned. For the inverse of the correction this is taken into account by reversing the order the corrections are applied in the kQueue case (no issue for kParallel).
Author
Magnus Mager, Stefan Rossegger, Jim Thomas
Date
27/04/2010

Definition at line 33 of file AliTPCComposedCorrection.h.

Member Enumeration Documentation

Enumerator
kParallel 
kQueue 
kQueueResidual 

Definition at line 35 of file AliTPCComposedCorrection.h.

Constructor & Destructor Documentation

AliTPCComposedCorrection::AliTPCComposedCorrection ( )

default constructor

Definition at line 70 of file AliTPCComposedCorrection.cxx.

AliTPCComposedCorrection::AliTPCComposedCorrection ( TCollection *  corrections,
AliTPCComposedCorrection::CompositionType  mode 
)

Constructor that defines the set of corrections, this one is composed of.

Definition at line 81 of file AliTPCComposedCorrection.cxx.

AliTPCComposedCorrection::~AliTPCComposedCorrection ( )
virtual

destructor

Definition at line 93 of file AliTPCComposedCorrection.cxx.

AliTPCComposedCorrection::AliTPCComposedCorrection ( const AliTPCComposedCorrection )
private

Member Function Documentation

Bool_t AliTPCComposedCorrection::AddCorrectionCompact ( AliTPCCorrection corr,
Double_t  weight 
)
virtual

Add correction - better name needed (left/right) - for the moment I assumme they commute Why not to just use array of corrections - CPU consideration Assumptions:

  • origin of distortion/correction are additive
  • corrections/distortion are small and they commute
  • only correction ot the same type supported

Reimplemented from AliTPCCorrection.

Definition at line 110 of file AliTPCComposedCorrection.cxx.

Referenced by TestCorrection_AliTPCBoundaryVoltErrorAddCorrectionCompact(), TestCorrection_AliTPCCalibGlobalMisalignmentAddCorrectionCompact(), TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact(), TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact_TPCCalibCorrection(), TestCorrection_AliTPCExBTwistAddCorrectionCompact(), TestCorrection_AliTPCFCVoltError3DAddCorrectionCompact(), and TestCorrection_AliTPCRocVoltError3DAddCorrectionCompact().

void AliTPCComposedCorrection::GetCorrection ( const Float_t  x[],
const Short_t  roc,
Float_t  dx[] 
)
virtual

This applies all correction and the specified manner (see general class description for details).

Reimplemented from AliTPCCorrection.

Definition at line 168 of file AliTPCComposedCorrection.cxx.

void AliTPCComposedCorrection::GetDistortion ( const Float_t  x[],
const Short_t  roc,
Float_t  dx[] 
)
virtual

This applies all distortions and the specified manner (see general class descxiption for details).

Reimplemented from AliTPCCorrection.

Definition at line 211 of file AliTPCComposedCorrection.cxx.

Referenced by GetCorrection().

CompositionType AliTPCComposedCorrection::GetMode ( ) const
inline

Definition at line 45 of file AliTPCComposedCorrection.h.

AliTPCCorrection * AliTPCComposedCorrection::GetSubCorrection ( const char *  cname)
virtual

Definition at line 159 of file AliTPCComposedCorrection.cxx.

const TVectorD* AliTPCComposedCorrection::GetWeights ( ) const
inline

Definition at line 59 of file AliTPCComposedCorrection.h.

void AliTPCComposedCorrection::Init ( )
virtual

Initialization funtion (not used at the moment)

Reimplemented from AliTPCCorrection.

Definition at line 295 of file AliTPCComposedCorrection.cxx.

Referenced by MakeComposedCorrection(), RegisterAliTPCROCVoltError3D(), RegisterAliTPCROCVoltError3DSector(), and AliTPCcalibDB::Update().

AliTPCComposedCorrection& AliTPCComposedCorrection::operator= ( const AliTPCComposedCorrection )
private
void AliTPCComposedCorrection::Print ( Option_t *  option = "") const
virtual

Print function to check which correction classes are used option=="d" prints details regarding the setted magnitude option=="a" prints the C0 and C1 coefficents for calibration purposes

Reimplemented from AliTPCCorrection.

Definition at line 266 of file AliTPCComposedCorrection.cxx.

Referenced by AliTPCCorrectionDemo(), LoadModels(), MakeComposedCorrection(), MakeComposedCorrectionExB(), RegisterAliTPCBoundaryVoltError(), RegisterAliTPCFCVoltError3D(), RegisterAliTPCFCVoltError3DRodFCSideRadiusType(), RegisterAliTPCROCVoltError3D(), and RegisterAliTPCROCVoltError3DSector().

void AliTPCComposedCorrection::SetMode ( CompositionType  mode)
inline

Definition at line 46 of file AliTPCComposedCorrection.h.

void AliTPCComposedCorrection::SetOmegaTauT1T2 ( Float_t  omegaTau,
Float_t  t1,
Float_t  t2 
)
virtual

Gives the possibility to set the OmegaTau plus Tensor corrections T1 and T2 (effective omega Tau) to each subcorrection (since they might become event specific due to changing drift velocity)

The omegaTau comes idealy from the Database, since it is a function of drift velocity, B and E field e.g. omegaTau = -10.0 * Bz * vdrift / Ez ; // with Bz in kG and Ez in V/cm omegaTau = -0.325 for Bz=5kG, Ez=400V/cm and vdrift = 2.6cm/muSec The T1 and T2 tensors were measured in a dedicated calibration run

Note: overwrites previously set values!

Reimplemented from AliTPCCorrection.

Definition at line 328 of file AliTPCComposedCorrection.cxx.

Referenced by AliTPCCorrectionDemo(), LoadModels(), MakeComposedCorrection(), MakeDistortionMap(), RegisterAliTPCBoundaryVoltError(), RegisterAliTPCFCVoltError3D(), RegisterAliTPCFCVoltError3DRodFCSideRadiusType(), RegisterAliTPCROCVoltError3D(), RegisterAliTPCROCVoltError3DSector(), TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact(), and TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact_TPCCalibCorrection().

void AliTPCComposedCorrection::SetWeights ( TVectorD *  weights)
inline
void AliTPCComposedCorrection::Update ( const TTimeStamp &  timeStamp)
virtual

Update function

Reimplemented from AliTPCCorrection.

Definition at line 310 of file AliTPCComposedCorrection.cxx.

Member Data Documentation

TCollection* AliTPCComposedCorrection::fCorrections
private
CompositionType AliTPCComposedCorrection::fMode
private

The way to apply the corrections (see general class documentation)

Definition at line 63 of file AliTPCComposedCorrection.h.

Referenced by GetCorrection(), GetDistortion(), GetMode(), and SetMode().

TVectorD* AliTPCComposedCorrection::fWeights
private

optional vector with weights - used for fit benchmarking

Definition at line 64 of file AliTPCComposedCorrection.h.

Referenced by AddCorrectionCompact(), GetCorrection(), GetDistortion(), GetWeights(), SetWeights(), and ~AliTPCComposedCorrection().


The documentation for this class was generated from the following files: