AliRoot Core  3dc7879 (3dc7879)
AliTPCDistortionFun.C
Go to the documentation of this file.
1 
46 #include "TCanvas.h"
47 #include "TF1.h"
48 #include "TLegend.h"
49 #include "AliTPCDistortions.h"
50 
52 
53 Double_t GetIFCDistortion(Double_t lx, Double_t ly, Double_t lz, Int_t icoord, Double_t shift=1.){
55 
56  static AliTPCDistortions transform;
57  static Bool_t doInit=kTRUE;
58  if (doInit){
59  transform.SetIFCShift(1);
60  transform.InitIFCShiftDistortion();
61  doInit=kFALSE;
62  }
63 
64  Double_t xyzIn[3]={lx, ly, lz};
65  Double_t xyzOut[3]={lx, ly, lz};
66  Int_t dummyROC=0;
67  if (lz<0) { dummyROC=36;}
68  transform.UndoIFCShiftDistortion(xyzIn,xyzOut,dummyROC);
69  Double_t result=0;
70  if (icoord<3) result=xyzOut[icoord]-xyzIn[icoord];
71  return result*shift;
72 }
73 
74 
75 
76 Double_t GetGGDistortion(Double_t lx, Double_t ly, Double_t lz, Int_t icoord, Double_t deltaVGGA=1., Double_t deltaVGGC=1.){
78 
79  static AliTPCDistortions transform;
80  static Bool_t doInit=kTRUE;
81  if (doInit){
82  transform.SetDeltaVGGA(1.);
83  transform.SetDeltaVGGC(1.);
84  transform.InitGGVoltErrorDistortion();
85  doInit=kFALSE;
86  }
87  Double_t xyzIn[3]={lx, ly, lz};
88  Double_t xyzOut[3]={lx, ly, lz};
89  Int_t dummyROC=0;
90  if (lz<0) { dummyROC=36;}
91  transform.UndoGGVoltErrorDistortion(xyzIn,xyzOut,dummyROC);
92  Double_t result=0;
93  if (icoord<3) result=xyzOut[icoord]-xyzIn[icoord];
94  if (lz<0) result*=deltaVGGA;
95  if (lz>0) result*=deltaVGGC;
96  return result;
97 }
98 
99 
102 
103  TCanvas *canvasIFC1D= new TCanvas("IFC radial shift", "IFC radial shift");
104  TLegend *legend = new TLegend(0.1,0.60,0.5,0.9, "Radial distortion due IFC shift 1 mm ");
105 
106  for (Int_t iradius=0; iradius<=10; iradius++){
107  Double_t radius=85+iradius*(245-85)/10.;
108  TF1 *f1= new TF1(Form("fIFC_ZX%f",radius),Form("0.1*GetIFCDistortion(%f,0,x,0)*sign(x)",radius),-250,250);
109  f1->SetMaximum( 0.6);
110  f1->SetMinimum(-0.6);
111  f1->SetNpx(200);
112  f1->SetLineColor(1+((20+iradius)%20));
113  f1->SetLineWidth(1);
114  if (iradius==0) f1->Draw("p");
115  f1->Draw("samep");
116  legend->AddEntry(f1,Form("R=%f",radius));
117  }
118  legend->Draw();
119 }
120 
123 
124  TCanvas *canvasIFC1D= new TCanvas("IFC radial shift - angle", "IFC radial shift angle");
125  TLegend *legend = new TLegend(0.1,0.60,0.9,0.9, "Radial distortion due IFC shift 1 mm ");
126 
127  for (Int_t iradius=0; iradius<=10; iradius++){
128  Double_t radius=85+iradius*(245-85)/10.;
129  TF1 *f1= new TF1(Form("fIFC_ZX%f",radius),Form("0.1*1000*(GetIFCDistortion(%f,0,x,0)-GetIFCDistortion(%f,0,x+1,0))",radius,radius),-250,250);
130  f1->SetMaximum( 10);
131  f1->SetMinimum(-10);
132  f1->SetNpx(200);
133  f1->SetLineColor(1+(iradius%10));
134  f1->SetLineWidth(1);
135  //f1->SetLineStyle(2+(iradius%6));
136  if (iradius==0) f1->Draw("");
137  f1->Draw("same");
138  legend->AddEntry(f1,Form("R=%f",radius));
139  }
140  legend->Draw();
141 }
142 
#define TObjArray
Double_t GetGGDistortion(Double_t lx, Double_t ly, Double_t lz, Int_t icoord, Double_t deltaVGGA=1., Double_t deltaVGGC=1.)
void MakePicIFCDXangle()
Double_t GetIFCDistortion(Double_t lx, Double_t ly, Double_t lz, Int_t icoord, Double_t shift=1.)
TObjArray * arrayPic
void MakePicIFCDX()