]>
Commit | Line | Data |
---|---|---|
57f936af | 1 | /* |
43e55a68 | 2 | Simple compiled macro for declaration of static distortion function |
3 | on top of the AliTPCDistortion class. | |
4 | Why: | |
5 | 1. Use static function in the fitting procedure | |
6 | 2. Usage in TFormual, TF1, Tf2 ... for visualization. | |
7 | 3. Usage in tree->Draw() for visualization | |
8 | 4. Simple visualization of fit residuals in multidemension - using tree Draw functionality | |
9 | ||
10 | ||
11 | ||
12 | Usage: | |
57f936af | 13 | gSystem->AddIncludePath("-I$ALICE_ROOT/TPC"); |
43e55a68 | 14 | .L $ALICE_ROOT/TPC/CalibMacros/AliTPCDistortions.cxx+ |
57f936af | 15 | .L $ALICE_ROOT/TPC/CalibMacros/AliTPCDistortionFun.C+ |
43e55a68 | 16 | |
17 | Example: | |
57f936af | 18 | // |
43e55a68 | 19 | // Draw integrated distortion in local x |
57f936af | 20 | // |
7f550ed6 | 21 | TF2 fdistIFCXvZX("fdistIFCXvZX","GetIFCDistortion(y,0,x,0)*sign(x)",-250,250,80,250); |
22 | fdistIFCXvZX.SetNpx(200); | |
23 | fdistIFCXvZX.SetNpy(200); | |
57f936af | 24 | fdistIFCXvZX.GetXaxis()->SetTitle("Z (cm)"); |
25 | fdistIFCXvZX.GetYaxis()->SetTitle("local X (cm)"); | |
26 | fdistIFCXvZX->Draw("colz"); | |
27 | // | |
43e55a68 | 28 | // Draw local distortion angle dx/dz in mrad |
57f936af | 29 | // |
7f550ed6 | 30 | TF2 fangleIFCXvZX("fangleIFCXvZX","1000*(GetIFCDistortion(y,0,x,0)-GetIFCDistortion(y,0,x-1,0))",-250,250,85,245); |
43e55a68 | 31 | fangleIFCXvZX.SetNpx(200); |
32 | fangleIFCXvZX.SetNpy(200); | |
57f936af | 33 | fangleIFCXvZX.GetXaxis()->SetTitle("Z (cm)"); |
34 | fangleIFCXvZX.GetYaxis()->SetTitle("local X (cm)"); | |
35 | fangleIFCXvZX->Draw("colz"); | |
36 | ||
7f550ed6 | 37 | // |
38 | // | |
39 | // | |
40 | TF2 fangleGGXvZX("fangleGGXvZX","1000*(GetGGDistortion(y,0,x,0,1*400,1*400)-GetGGDistortion(y,0,x-1,0,1*400,1*400))*sign(x)",-250,250,85,245); | |
41 | fangleGGXvZX.SetNpx(200); | |
42 | fangleGGXvZX.SetNpy(200); | |
43 | fangleGGXvZX.GetXaxis()->SetTitle("Z (cm)"); | |
44 | fangleGGXvZX.GetYaxis()->SetTitle("local X (cm)"); | |
45 | fangleGGXvZX->Draw("colz"); | |
46 | ||
47 | ||
48 | ||
57f936af | 49 | */ |
50 | ||
51 | ||
7f550ed6 | 52 | #include "TCanvas.h" |
53 | #include "TF1.h" | |
54 | #include "TLegend.h" | |
57f936af | 55 | #include "AliTPCDistortions.h" |
56 | ||
7f550ed6 | 57 | TObjArray *arrayPic=new TObjArray; |
57f936af | 58 | |
59 | Double_t GetIFCDistortion(Double_t lx, Double_t ly, Double_t lz, Int_t icoord, Double_t shift=1.){ | |
60 | // | |
61 | static AliTPCDistortions transform; | |
7f550ed6 | 62 | static Bool_t doInit=kTRUE; |
63 | if (doInit){ | |
64 | transform.SetIFCShift(1); | |
65 | transform.InitIFCShiftDistortion(); | |
66 | doInit=kFALSE; | |
67 | } | |
68 | ||
57f936af | 69 | Double_t xyzIn[3]={lx, ly, lz}; |
70 | Double_t xyzOut[3]={lx, ly, lz}; | |
57f936af | 71 | Int_t dummyROC=0; |
72 | if (lz<0) { dummyROC=36;} | |
73 | transform.UndoIFCShiftDistortion(xyzIn,xyzOut,dummyROC); | |
74 | Double_t result=0; | |
75 | if (icoord<3) result=xyzOut[icoord]-xyzIn[icoord]; | |
76 | return result*shift; | |
77 | } | |
7f550ed6 | 78 | |
79 | ||
80 | ||
81 | Double_t GetGGDistortion(Double_t lx, Double_t ly, Double_t lz, Int_t icoord, Double_t deltaVGGA=1., Double_t deltaVGGC=1.){ | |
82 | // | |
83 | // GG distortion induced distortions | |
84 | // | |
85 | static AliTPCDistortions transform; | |
86 | static Bool_t doInit=kTRUE; | |
87 | if (doInit){ | |
88 | transform.SetDeltaVGGA(1.); | |
89 | transform.SetDeltaVGGC(1.); | |
90 | transform.InitGGVoltErrorDistortion(); | |
91 | doInit=kFALSE; | |
92 | } | |
93 | Double_t xyzIn[3]={lx, ly, lz}; | |
94 | Double_t xyzOut[3]={lx, ly, lz}; | |
95 | Int_t dummyROC=0; | |
96 | if (lz<0) { dummyROC=36;} | |
97 | transform.UndoGGVoltErrorDistortion(xyzIn,xyzOut,dummyROC); | |
98 | Double_t result=0; | |
99 | if (icoord<3) result=xyzOut[icoord]-xyzIn[icoord]; | |
100 | if (lz<0) result*=deltaVGGA; | |
101 | if (lz>0) result*=deltaVGGC; | |
102 | return result; | |
103 | } | |
104 | ||
105 | ||
106 | void MakePicIFCDX(){ | |
107 | // | |
108 | // | |
109 | // | |
110 | TCanvas *canvasIFC1D= new TCanvas("IFC radial shift", "IFC radial shift"); | |
111 | TLegend *legend = new TLegend(0.1,0.60,0.5,0.9, "Radial distortion due IFC shift 1 mm "); | |
112 | ||
113 | for (Int_t iradius=0; iradius<=10; iradius++){ | |
114 | Double_t radius=85+iradius*(245-85)/10.; | |
115 | TF1 *f1= new TF1(Form("fIFC_ZX%f",radius),Form("0.1*GetIFCDistortion(%f,0,x,0)*sign(x)",radius),-250,250); | |
116 | f1->SetMaximum( 0.6); | |
117 | f1->SetMinimum(-0.6); | |
118 | f1->SetNpx(200); | |
119 | f1->SetLineColor(1+((20+iradius)%20)); | |
120 | f1->SetLineWidth(1); | |
121 | if (iradius==0) f1->Draw("p"); | |
122 | f1->Draw("samep"); | |
123 | legend->AddEntry(f1,Form("R=%f",radius)); | |
124 | } | |
125 | legend->Draw(); | |
126 | } | |
127 | ||
128 | void MakePicIFCDXangle(){ | |
129 | // | |
130 | // | |
131 | // | |
132 | TCanvas *canvasIFC1D= new TCanvas("IFC radial shift - angle", "IFC radial shift angle"); | |
133 | TLegend *legend = new TLegend(0.1,0.60,0.9,0.9, "Radial distortion due IFC shift 1 mm "); | |
134 | ||
135 | for (Int_t iradius=0; iradius<=10; iradius++){ | |
136 | Double_t radius=85+iradius*(245-85)/10.; | |
137 | 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); | |
138 | f1->SetMaximum( 10); | |
139 | f1->SetMinimum(-10); | |
140 | f1->SetNpx(200); | |
141 | f1->SetLineColor(1+(iradius%10)); | |
142 | f1->SetLineWidth(1); | |
143 | //f1->SetLineStyle(2+(iradius%6)); | |
144 | if (iradius==0) f1->Draw(""); | |
145 | f1->Draw("same"); | |
146 | legend->AddEntry(f1,Form("R=%f",radius)); | |
147 | } | |
148 | legend->Draw(); | |
149 | } | |
150 |