1 #ifndef ALI_TPC_CORRECTION_H
2 #define ALI_TPC_CORRECTION_H
4 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
5 * See cxx source for full Copyright notice */
7 // _________________________________________________________________
10 // <h2> AliTPCCorrection class </h2>
12 // This class provides a general framework to deal with space point
13 // distortions. An correction class which inherits from here is for example
14 // AliTPCExBBShape or AliTPCExBTwist
16 // General functions are (for example): <br>
17 // CorrectPoint(x,roc) where x is the vector of inital positions in
18 // cartesian coordinates and roc represents the Read Out chamber number
19 // according to the offline naming convention. The vector x is overwritten
20 // with the corrected coordinates.
22 // An alternative usage would be CorrectPoint(x,roc,dx), which leaves the
23 // vector x untouched, put returns the distortions via the vector dx
25 // Several plot functionalities (see example below), general solvers as well as simplified interpolation techniques are implemented.
27 // The class allows "effective Omega Tau" corrections to be shifted to the
28 // single distortion classes.
30 // Note: This class is normally used via the class AliTPCComposedCorrection
33 // Begin_Macro(source)
35 // gROOT->SetStyle("Plain"); gStyle->SetPalette(1);
36 // TCanvas *c2 = new TCanvas("c2","c2",800,1200); c2->Divide(2,3);
37 // AliTPCROCVoltError3D roc; // EXAMPLE PLOTS - SEE BELOW
38 // roc.SetOmegaTauT1T2(0,1,1); // B=0
39 // Float_t z0 = 1; // at +1 cm -> A side
40 // c2->cd(1); roc.CreateHistoDRinXY(1.,300,300)->Draw("cont4z");
41 // c2->cd(3);roc.CreateHistoDRPhiinXY(1.,300,300)->Draw("cont4z");
42 // c2->cd(5);roc.CreateHistoDZinXY(1.,300,300)->Draw("cont4z");
44 // c2->cd(2);roc.CreateHistoDRinZR(phi0)->Draw("surf2");
45 // c2->cd(4);roc.CreateHistoDRPhiinZR(phi0)->Draw("surf2");
46 // c2->cd(6);roc.CreateHistoDZinZR(phi0)->Draw("surf2");
53 // Date: 27/04/2010 <br>
54 // Authors: Magnus Mager, Stefan Rossegger, Jim Thomas
56 // _________________________________________________________________
58 ////////////////////////////////////////////////////////////////////////////////
59 // AliTPCCorrection class //
60 ////////////////////////////////////////////////////////////////////////////////
68 class TTreeSRedirector;
69 class AliExternalTrackParam;
75 class AliTPCCorrection : public TNamed {
77 enum CompositionType {kParallel,kQueue};
80 AliTPCCorrection(const char *name,const char *title);
81 virtual ~AliTPCCorrection();
84 // functions to correct a space point
85 void CorrectPoint ( Float_t x[],const Short_t roc);
86 void CorrectPoint (const Float_t x[],const Short_t roc,Float_t xp[]);
87 virtual void GetCorrection(const Float_t x[],const Short_t roc,Float_t dx[]);
89 // functions to distort a space point
90 void DistortPoint ( Float_t x[],const Short_t roc);
91 void DistortPoint (const Float_t x[],const Short_t roc,Float_t xp[]);
92 virtual void GetDistortion(const Float_t x[],const Short_t roc,Float_t dx[]);
94 // initialization and update functions
96 virtual void Update(const TTimeStamp &timeStamp);
98 // convenience functions
99 virtual void Print(Option_t* option="") const;
101 TH2F* CreateHistoDRinXY (Float_t z=10.,Int_t nx=100,Int_t ny=100);
102 TH2F* CreateHistoDRPhiinXY(Float_t z=10.,Int_t nx=100,Int_t nphi=100);
103 TH2F* CreateHistoDZinXY (Float_t z=10.,Int_t nx=100,Int_t ny=100);
105 TH2F* CreateHistoDRinZR (Float_t phi=0.,Int_t nZ=100,Int_t nR=100);
106 TH2F* CreateHistoDRPhiinZR(Float_t phi=0.,Int_t nZ=100,Int_t nR=100);
107 TH2F* CreateHistoDZinZR (Float_t phi=0.,Int_t nZ=100,Int_t nR=100);
110 TTree* CreateDistortionTree(Double_t step=5);
111 static void MakeDistortionMap(THnSparse * his0, TTreeSRedirector *pcstream, const char* hname, Int_t run);
112 // normally called directly in the correction classes which inherit from this class
113 virtual void SetOmegaTauT1T2(Float_t omegaTau,Float_t t1,Float_t t2);
114 AliExternalTrackParam * FitDistortedTrack(AliExternalTrackParam & trackIn, Double_t refX, Int_t dir,TTreeSRedirector *pcstream);
115 void StoreInOCDB(Int_t startRun, Int_t endRun, const char *comment=0);
116 static void MakeTrackDistortionTree(TTree *tinput, Int_t dtype, Int_t ptype, const TObjArray * corrArray, Int_t step=1, Bool_t debug=0);
117 static void MakeLaserDistortionTree(TTree* tree, TObjArray *corrArray, Int_t itype);
119 void FastSimDistortedVertex(Double_t orgVertex[3], Int_t nTracks, AliESDVertex &aV, AliESDVertex &avOrg, AliESDVertex &cV, AliESDVertex &cvOrg, TTreeSRedirector * const pcstream, Double_t etaCuts);
121 static void AddVisualCorrection(AliTPCCorrection* corr, Int_t position);
122 static Double_t GetCorrSector(Double_t sector, Double_t r, Double_t kZ, Int_t axisType, Int_t corrType=0);
123 static Double_t GetCorrXYZ(Double_t gx, Double_t gy, Double_t gz, Int_t axisType, Int_t corrType=0);
126 TH2F* CreateTH2F(const char *name,const char *title,
127 const char *xlabel,const char *ylabel,const char *zlabel,
128 Int_t nbinsx,Double_t xlow,Double_t xup,
129 Int_t nbinsy,Double_t ylow,Double_t yup);
131 static const Double_t fgkTPCZ0; // nominal gating grid position
132 static const Double_t fgkIFCRadius; // Mean Radius of the Inner Field Cage ( 82.43 min, 83.70 max) (cm)
133 static const Double_t fgkOFCRadius; // Mean Radius of the Outer Field Cage (252.55 min, 256.45 max) (cm)
134 static const Double_t fgkZOffSet; // Offset from CE: calculate all distortions closer to CE as if at this point
135 static const Double_t fgkCathodeV; // Cathode Voltage (volts)
136 static const Double_t fgkGG; // Gating Grid voltage (volts)
137 static const Double_t fgkdvdE; // [cm/V] drift velocity dependency on the E field (from Magboltz for NeCO2N2 at standard environment)
138 static const Double_t fgkEM; // charge/mass in [C/kg]
139 static const Double_t fgke0; // vacuum permittivity [A·s/(V·m)]
142 enum {kNR= 72 }; // Number of R points in the table for interpolating distortion data
143 enum {kNPhi= 18*10+1}; // Number of Phi points in the table for interpolating distortion data ( plus one extra for 360 == 0 )
144 enum {kNZ= 166}; // Number of Z points in the table for interpolating distortion data
146 Double_t fgkRList[kNR]; // points in the radial direction (for the lookup table)
147 Double_t fgkPhiList[kNPhi]; // points in the phi direction (for the lookup table)
148 Double_t fgkZList[kNZ]; // points in the z direction (for the lookup table)
150 // Simple Interpolation functions: e.g. with tricubic interpolation (not yet in TH3)
151 Int_t fILow, fJLow, fKLow; // variable to help in the interpolation
152 void Interpolate2DEdistortion( const Int_t order, const Double_t r, const Double_t z,
153 const Double_t er[kNZ][kNR], Double_t &erValue );
154 void Interpolate3DEdistortion( const Int_t order, const Double_t r, const Float_t phi, const Double_t z,
155 const Double_t er[kNZ][kNPhi][kNR], const Double_t ephi[kNZ][kNPhi][kNR],
156 const Double_t ez[kNZ][kNPhi][kNR],
157 Double_t &erValue, Double_t &ephiValue, Double_t &ezValue);
158 Double_t Interpolate2DTable( const Int_t order, const Double_t x, const Double_t y,
159 const Int_t nx, const Int_t ny, const Double_t xv[], const Double_t yv[],
160 const TMatrixD &array );
161 Double_t Interpolate3DTable( const Int_t order, const Double_t x, const Double_t y, const Double_t z,
162 const Int_t nx, const Int_t ny, const Int_t nz,
163 const Double_t xv[], const Double_t yv[], const Double_t zv[],
164 TMatrixD **arrayofArrays );
165 Double_t Interpolate( const Double_t xArray[], const Double_t yArray[],
166 const Int_t order, const Double_t x );
167 void Search( const Int_t n, const Double_t xArray[], const Double_t x, Int_t &low );
168 virtual Int_t IsPowerOfTwo ( Int_t i ) const ;
172 // Algorithms to solve the laplace or poisson equation
173 void PoissonRelaxation2D(TMatrixD &arrayV, TMatrixD &chargeDensity,
174 TMatrixD &arrayErOverEz, TMatrixD &arrayDeltaEz,
175 const Int_t rows, const Int_t columns, const Int_t iterations,
176 const Bool_t rocDisplacement = kTRUE);
178 void PoissonRelaxation3D( TMatrixD **arrayofArrayV, TMatrixD **arrayofChargeDensities,
179 TMatrixD **arrayofEroverEz, TMatrixD **arrayofEPhioverEz, TMatrixD **arrayofEz,
180 const Int_t rows, const Int_t columns, const Int_t phislices,
181 const Float_t deltaphi, const Int_t iterations, const Int_t summetry,
182 const Bool_t rocDisplacement = kTRUE);
185 Double_t fT1; // tensor term of wt - T1
186 Double_t fT2; // tensor term of wt - T2
187 static TObjArray *fgVisualCorrection; // array of orrection for visualization
189 AliTPCCorrection(const AliTPCCorrection &); // not implemented
190 AliTPCCorrection &operator=(const AliTPCCorrection &); // not implemented
192 void InitLookUpfulcrums(); // to initialize the grid of the look up table
194 ClassDef(AliTPCCorrection,4);