]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TPC/CalibMacros/MakeLookup.C
M AliTPCcalibBase.h - SetRun function
[u/mrichter/AliRoot.git] / TPC / CalibMacros / MakeLookup.C
CommitLineData
ddd07cdf 1//
2// Make lookup of distortion TPC distortion +(ITS and TRD)
3// Input: Residual histograms obtained in the AliTPCcalibTime
4
5// Residual histograms:
6// 1. TPC-ITS - entrance of the TPC
7// 2. TPC-ITS - at the vertex
8// 3. TPC-TRD - outer wall of the TPC
9
10// Histogram binning:
11// 1. Theta - fP3
12// 2. Phi - global phi at the entrance (case 1,2) and at the outer wall of TPC (case 3)
13// 3. snp(phi) - fP2 - local inclination angle at reference X
14
15// Output value:
16// Mean residuals, rms and number of entries in each bing
17//
18//
19
20
21/*
22 Example usage:
7f4cb119 23 gROOT->Macro("~/rootlogon.C")
ddd07cdf 24 gSystem->AddIncludePath("-I$ALICE_ROOT/STAT")
25 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC")
26 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros")
27 gSystem->Load("libANALYSIS");
28 gSystem->Load("libTPCcalib");
29 .L $ALICE_ROOT/TPC/CalibMacros/MakeLookup.C+
7f4cb119 30 Int_t run=115888
31 // .x $ALICE_ROOT/ANALYSIS/CalibMacros/Pass0/ConfigCalibTrain.C(run,"local:///lustre/alice/alien/alice/data/2010/OCDB")
32
ddd07cdf 33 MakeLookup(run,0);
34
35 //Local check of procedure
36 TTreeSRedirector * pcstream = new TTreeSRedirector("mean.root");
37 TFile f("CalibObjects.root");
38 AliTPCcalibTime *calibTime= f.Get("calibTime");
39 THnSparse * his = calibTime->GetResHistoTPCITS(0);
40 his->GetAxis(1)->SetRangeUser(-1,1);
41 his->GetAxis(2)->SetRange(0,1000000);
42 his->GetAxis(3)->SetRangeUser(-0.3,0.3);
43 //
44
45
46
47*/
48
49#if !defined(__CINT__) || defined(__MAKECINT__)
50#include "THnSparse.h"
51#include "TLatex.h"
52#include "TCanvas.h"
53#include "TLegend.h"
54#include "TSystem.h"
55#include "TFile.h"
56#include "TChain.h"
57#include "TCut.h"
58#include "TH3.h"
59#include "TProfile3D.h"
60#include "TMath.h"
61#include "TVectorD.h"
62#include "TMatrixD.h"
63#include "TStatToolkit.h"
64#include "TTreeStream.h"
65#include "AliExternalTrackParam.h"
66#include "AliESDfriend.h"
67#include "AliTPCcalibTime.h"
68#include "TROOT.h"
69#include "AliXRDPROOFtoolkit.h"
70#include "AliTPCCorrection.h"
71#include "AliTPCExBTwist.h"
72#include "AliTPCGGVoltError.h"
73#include "AliTPCComposedCorrection.h"
74#include "AliTPCExBConical.h"
7f4cb119 75#include "TPostScript.h"
76#include "TStyle.h"
ddd07cdf 77#include "AliTrackerBase.h"
7f4cb119 78#include "AliTPCCalibGlobalMisalignment.h"
79#include "AliTPCExBEffective.h"
80#include "AliTPCExBBShape.h"
ddd07cdf 81#endif
82
83TChain * chain=0;
7f4cb119 84void MakeLookup(Int_t run, Int_t mode);
ddd07cdf 85//void MakeLookupHisto(THnSparse * his, TTreeSRedirector *pcstream, const char* hname, Int_t run);
86void FitLookup(TChain *chainIn, const char *prefix, TVectorD &vecA, TVectorD &vecC, TVectorD& vecStatA, TVectorD &vecStatD, TCut cut, TObjArray *picArray);
87void MakeFits(Int_t run);
88void MakeFitTree();
89void DrawDistortionDy(TCut cutUser="", Double_t ymin=-0.6, Double_t ymax=0.6);
90void DrawDistortionMaps(const char *fname="mean.root");
7f4cb119 91TCanvas * DrawDistortionMaps(TTree * treedy, TTree * treedsnp, TTree* treed1Pt, const char * name, const char *title);
92AliTPCComposedCorrection * MakeComposedCorrection();
93void MakeLaserTree();
94void AddEffectiveCorrection(AliTPCComposedCorrection* comp);
ddd07cdf 95
96void MakeLookup(Int_t run, Int_t mode){
97 //
98 // make a lookup tree with mean values
7f4cb119 99 // 5. make laser tree
100 // 4.
ddd07cdf 101 gSystem->AddIncludePath("-I$ALICE_ROOT/STAT");
102 gSystem->Load("libANALYSIS");
103 gSystem->Load("libTPCcalib");
7f4cb119 104 if (mode==5) {MakeLaserTree();return;};
105 if (mode==4) {DrawDistortionMaps();return;}
ddd07cdf 106 if (mode==1){
107 DrawDistortionDy("abs(snp)<0.25");
108 DrawDistortionMaps();
109 return;
110 }
111 if (mode==2) return MakeFitTree();
7f4cb119 112 TFile f("TPCTimeObjects.root");
ddd07cdf 113 AliTPCcalibTime *calibTime= (AliTPCcalibTime*)f.Get("calibTime");
7f4cb119 114 if (!calibTime){
115 TFile* f2 = new TFile("CalibObjects.root");
116 calibTime= (AliTPCcalibTime*)f2->Get("calibTime");
117 if (!calibTime){
118 TFile* f3 = new TFile("../CalibObjects.root");
119 calibTime= (AliTPCcalibTime*)f3->Get("calibTime");
120 }
121 }
ddd07cdf 122 //
123 TTreeSRedirector * pcstream = new TTreeSRedirector("mean.root");
124 THnSparse * his = 0;
125 const char * hname[5]={"dy","dz","dsnp","dtheta","d1pt"};
126 if (calibTime){
127 for (Int_t ihis=0; ihis<5;ihis++){
128 his = calibTime->GetResHistoTPCITS(ihis);
129 his->GetAxis(1)->SetRangeUser(-1.1,1.1);
130 his->GetAxis(2)->SetRange(0,1000000);
131 his->GetAxis(3)->SetRangeUser(-0.5,0.5);
132 AliTPCCorrection::MakeDistortionMap(his,pcstream, Form("ITS%s",hname[ihis]),run);
133 //
134 his = calibTime->GetResHistoTPCvertex(ihis);
135 his->GetAxis(1)->SetRangeUser(-1.1,1.1);
136 his->GetAxis(2)->SetRange(0,1000000);
137 his->GetAxis(3)->SetRangeUser(-0.5,0.5);
138 AliTPCCorrection::MakeDistortionMap(his,pcstream, Form("Vertex%s",hname[ihis]),run);
139 //
140 his = calibTime->GetResHistoTPCTRD(ihis);
141 his->GetAxis(1)->SetRangeUser(-1.1,1.1);
142 his->GetAxis(2)->SetRange(0,1000000);
143 his->GetAxis(3)->SetRangeUser(-0.5,0.5);
144 AliTPCCorrection::MakeDistortionMap(his,pcstream, Form("TRD%s",hname[ihis]),run);
7f4cb119 145 his = calibTime->GetResHistoTPCTOF(ihis);
146 if (his){
147 his->GetAxis(1)->SetRangeUser(-1.1,1.1);
148 his->GetAxis(2)->SetRange(0,1000000);
149 his->GetAxis(3)->SetRangeUser(-0.5,0.5);
150 AliTPCCorrection::MakeDistortionMap(his,pcstream, Form("TOF%s",hname[ihis]),run);
151 }
ddd07cdf 152 }
153 }
154 delete pcstream;
155 gSystem->Exec("echo `pwd`/mean.root > distort.txt");
156 MakeFits(run);
157}
158
159
160void MakeFits(Int_t run){
161 //
162 // Make the fits of distortion
163 // store fit results and QA pictures in the file distortFit.root
164 //
165 TCut cut="entries>50&&rms>0";
166 TTreeSRedirector *pcstream = new TTreeSRedirector("distortFit.root");
167 AliXRDPROOFtoolkit tool;
168 TVectorD vecA[100],vecC[100], vecStatA[100],vecStatC[100];
169 TObjArray * picArray= new TObjArray();
170 //
171 TChain *chainITSdy = tool.MakeChain("distort.txt","ITSdy",0,100000);
172 TChain *chainITSdsnp = tool.MakeChain("distort.txt","ITSdsnp",0,100000);
173 TChain *chainTRDdy = tool.MakeChain("distort.txt","TRDdy",0,100000);
174 TChain *chainTRDdsnp = tool.MakeChain("distort.txt","TRDdsnp",0,100000);
175 TChain *chainVertexdy = tool.MakeChain("distort.txt","Vertexdy",0,100000);
176 TChain *chainVertexdsnp = tool.MakeChain("distort.txt","Vertexdsnp",0,100000);
177
178 (*pcstream)<<"fits"<<"run="<<run;
179
180 FitLookup(chainITSdy, "yits-tpc",vecA[0],vecC[0],vecStatA[0],vecStatC[0],cut, picArray);
181 (*pcstream)<<"fits"<<"itsdyA.="<<&vecA[0]<<"itsdyC.="<<&vecC[0];
182 (*pcstream)<<"fits"<<"itsdyAS.="<<&vecStatA[0]<<"itsdyCS.="<<&vecStatC[0];
183 FitLookup(chainITSdsnp, "snpits-tpc",vecA[1],vecC[1],vecStatA[1],vecStatC[1],cut, picArray);
184 (*pcstream)<<"fits"<<"itsdsnpA.="<<&vecA[1]<<"itsdsnpC.="<<&vecC[1];
185 (*pcstream)<<"fits"<<"itsdsnpAS.="<<&vecStatA[1]<<"itsdsnpCS.="<<&vecStatC[1];
186 //
187 FitLookup(chainTRDdy, "ytrd-tpc",vecA[5],vecC[5],vecStatA[5],vecStatC[5],cut, picArray);
188 (*pcstream)<<"fits"<<"trddyA.="<<&vecA[5]<<"trddyC="<<&vecC[5];
189 (*pcstream)<<"fits"<<"trddyAS.="<<&vecStatA[5]<<"trddyCS.="<<&vecStatC[5];
190 FitLookup(chainTRDdsnp, "snptrd-tpc",vecA[6],vecC[6],vecStatA[6],vecStatC[6],cut, picArray);
191 (*pcstream)<<"fits"<<"trddsnpA.="<<&vecA[6]<<"trddsnpC.="<<&vecC[6];
192 (*pcstream)<<"fits"<<"trddsnpAS.="<<&vecStatA[6]<<"trddsnpCS.="<<&vecStatC[6];
193 //
194 FitLookup(chainVertexdy, "yvertex-tpc",vecA[10],vecC[10],vecStatA[10],vecStatC[10],cut, picArray);
195 (*pcstream)<<"fits"<<"vertexdyA.="<<&vecA[10]<<"vertexdyC.="<<&vecC[10];
196 (*pcstream)<<"fits"<<"vertexdyAS.="<<&vecStatA[10]<<"vertexdyCS.="<<&vecStatC[10];
197 FitLookup(chainVertexdsnp, "snpvertex-tpc",vecA[11],vecC[11],vecStatA[11],vecStatC[11],cut, picArray);
198 (*pcstream)<<"fits"<<"vertexdsnpA.="<<&vecA[11]<<"vertexdsnpC.="<<&vecC[11];
199 (*pcstream)<<"fits"<<"vertexdsnpAS.="<<&vecStatA[11]<<"vertexdsnpCS.="<<&vecStatC[11];
200 //
201 (*pcstream)<<"fits"<<"\n";
202
203 pcstream->GetFile()->cd();
204 for (Int_t i=0;i<picArray->GetEntries(); i++){
205 TObject * obj = picArray->At(i);
206 if (obj) obj->Write(obj->GetName());
207 }
208 //
209 //
210 //
211 chainITSdy->AddFriend(chainITSdy,"ITSdy");
212 chainITSdy->AddFriend(chainITSdsnp,"ITSdsnp");
213 chainITSdy->AddFriend(chainTRDdy,"TRDdy");
214 chainITSdy->AddFriend(chainTRDdsnp,"TRDdsnp");
215 chainITSdy->AddFriend(chainVertexdy,"Vertexdy");
216 chainITSdy->AddFriend(chainVertexdsnp,"Vertexdsnp");
217 TTree * tree = chainITSdy->CloneTree();
218 tree->Write("distortionTree");
219 delete pcstream;
220
221 //
222 //
223}
224
225
226
227
228
229void FitLookup(TChain *chainIn, const char *prefix, TVectorD &vecA, TVectorD &vecC, TVectorD& vecStatA, TVectorD &vecStatC, TCut cut, TObjArray *picArray){
230 // TCut cut="entries>100&&rms>0";
231 vecStatA.ResizeTo(6);
232 vecStatC.ResizeTo(6);
233 vecA.ResizeTo(10);
234 vecC.ResizeTo(10);
235 Int_t npointsMax=30000000;
236 TStatToolkit toolkit;
237 Double_t chi2A=0;
238 Double_t chi2C=0;
239 Int_t npointsA=0;
240 Int_t npointsC=0;
241 TVectorD paramA;
242 TVectorD paramC;
243 TMatrixD covar;
244 TString *strFitYA;
245 TString *strFitYC;
246
247 TString fstring=""; // magnetic part
248 //fstring+="(1)++"; //0 - offset value
249 fstring+="(cos(phi))++"; //1 - cos part
250 fstring+="(sin(phi))++"; //2 - sin part
251 //
252 fstring+="(theta)++"; //3
253 fstring+="(theta*cos(phi))++"; //4
254 fstring+="(theta*sin(phi))++"; //5
255 //
256 fstring+="(snp)++"; //6 - delta X(radial) coeficients - offset
257 fstring+="(snp*cos(phi))++"; //7 - delta X(radial) coeficient
258 fstring+="(snp*sin(phi))++"; //8
259 //
260 //
261 strFitYA = TStatToolkit::FitPlane(chainIn,"mean", fstring.Data(),cut+"theta>0", chi2A,npointsA,paramA,covar,-1,0, npointsMax, kFALSE);
262 strFitYC = TStatToolkit::FitPlane(chainIn,"mean", fstring.Data(),cut+"theta<0", chi2C,npointsC,paramC,covar,-1,0, npointsMax,kFALSE);
263 strFitYA->Tokenize("++")->Print();
264 strFitYC->Tokenize("++")->Print();
265 chainIn->SetAlias("dyA",strFitYA->Data());
266 chainIn->SetAlias("dyC",strFitYC->Data());
267 //
268 TH1 * his=0;
269 vecA.ResizeTo(paramA.GetNrows());
270 vecC.ResizeTo(paramC.GetNrows());
271 vecA=paramA;
272 vecC=paramC;
273 // A side stat
274 vecStatA[0]=npointsA;
275 vecStatA[1]=TMath::Sqrt(chi2A/npointsA);
276 chainIn->Draw("mean",cut+"theta>0");
277 his=chainIn->GetHistogram();
278 his->SetName(Form("Orig #Delta%sA",prefix));
279 his->SetTitle(Form("Orig #Delta%sA",prefix));
280 picArray->AddLast(his->Clone());
281 vecStatA[2] = his->GetMean();
282 vecStatA[3] = his->GetRMS();
283 chainIn->Draw("mean-dyA",cut+"theta>0");
284 his=chainIn->GetHistogram();
285 his->SetName(Form("Corr. #Delta%sA",prefix));
286 his->SetTitle(Form("Corr. #Delta%sA",prefix));
287 picArray->AddLast(his->Clone());
288 vecStatA[4] = his->GetMean();
289 vecStatA[5] = his->GetRMS();
290 //
291 // C side stat
292 vecStatC[0]=npointsC;
293 vecStatC[1]=TMath::Sqrt(chi2C/npointsC);
294 chainIn->Draw("mean",cut+"theta<0");
295 his=chainIn->GetHistogram();
296 his->SetName(Form("Orig #Delta%sC",prefix));
297 his->SetTitle(Form("Orig #Delta%sC",prefix));
298 picArray->AddLast(his->Clone());
299 vecStatC[2] = his->GetMean();
300 vecStatC[3] = his->GetRMS();
301 chainIn->Draw("mean-dyC",cut+"theta<0");
302 his=chainIn->GetHistogram();
303 his->SetName(Form("Corr. #Delta%sC",prefix));
304 his->SetTitle(Form("Corr. #Delta%sC",prefix));
305 picArray->AddLast(his->Clone());
306 vecStatC[4] = his->GetMean();
307 vecStatC[5] = his->GetRMS();
308
309 vecStatA.Print();
310 vecStatC.Print();
311}
312
313
314
315void DrawDistortionDy(TCut cutUser, Double_t ymin, Double_t ymax){
316 //
317 //
318 //
319 TFile fplus("meanBplus.root");
320 TFile fminus("meanBminus.root");
321 TTree * titsDyPlus= (TTree*)fplus.Get("ITSdy");
322 TTree * titsDyMinus= (TTree*)fminus.Get("ITSdy");
323 TTree * ttrdDyPlus= (TTree*)fplus.Get("TRDdy");
324 TTree * ttrdDyMinus= (TTree*)fminus.Get("TRDdy");
325 //
326 TTree * titsDsnpPlus= (TTree*)fplus.Get("ITSdsnp");
327 TTree * titsDsnpMinus= (TTree*)fminus.Get("ITSdsnp");
328 TTree * ttrdDsnpPlus= (TTree*)fplus.Get("TRDdsnp");
329 TTree * ttrdDsnpMinus= (TTree*)fminus.Get("TRDdsnp");
330 TTree * titsDthetaPlus= (TTree*)fplus.Get("ITSdtheta");
331 TTree * titsDthetaMinus= (TTree*)fminus.Get("ITSdtheta");
332 TTree * ttrdDthetaPlus= (TTree*)fplus.Get("TRDdtheta");
333 TTree * ttrdDthetaMinus= (TTree*)fminus.Get("TRDdtheta");
334 TCut cut="rms>0&&entries>100&&abs(theta)<1&&abs(snp)<0.3";
335 cut+=cutUser;
336 //
337 titsDyPlus->AddFriend(titsDyMinus,"TM.");
338 titsDyPlus->AddFriend(titsDsnpMinus,"Tsnp.");
339 titsDyPlus->AddFriend(titsDthetaMinus,"Ttheta.");
340 ttrdDyPlus->AddFriend(ttrdDyMinus,"TM.");
341 titsDsnpPlus->AddFriend(titsDsnpMinus,"TM.");
342 ttrdDsnpPlus->AddFriend(ttrdDsnpMinus,"TM.");
343
344 //
345 //
346 //
347 Int_t entries=0;
348 TGraph * graphITSYC[3];
349 TGraph * graphITSYA[3];
350 TGraph * graphTRDYC[3];
351 TGraph * graphTRDYA[3];
352 //
353 //
354 entries=titsDyPlus->Draw("mean-TM..mean:phi",cut+"theta<-0.1","goff");
355 graphITSYC[0] = new TGraph(entries, titsDyPlus->GetV2(),titsDyPlus->GetV1());
356 titsDyMinus->Draw("mean:phi",cut+"theta<-0.1","goff");
357 graphITSYC[2] = new TGraph(entries, titsDyMinus->GetV2(),titsDyMinus->GetV1());
358 titsDyPlus->Draw("mean:phi",cut+"theta<-0.1","goff");
359 graphITSYC[1] = new TGraph(entries, titsDyPlus->GetV2(),titsDyPlus->GetV1());
360 //
361 entries=titsDyPlus->Draw("mean-TM..mean:phi",cut+"theta>0.1","goff");
362 graphITSYA[0] = new TGraph(entries, titsDyPlus->GetV2(),titsDyPlus->GetV1());
363 titsDyMinus->Draw("mean:phi",cut+"theta>0.1","goff");
364 graphITSYA[2] = new TGraph(entries, titsDyMinus->GetV2(),titsDyMinus->GetV1());
365 titsDyPlus->Draw("mean:phi",cut+"theta>0.1","goff");
366 graphITSYA[1] = new TGraph(entries, titsDyPlus->GetV2(),titsDyPlus->GetV1());
367 //
368 entries=ttrdDyPlus->Draw("mean-TM..mean:phi",cut+"theta<-0.1&&TM..entries>100","goff");
369 graphTRDYC[0] = new TGraph(entries, ttrdDyPlus->GetV2(),ttrdDyPlus->GetV1());
370 ttrdDyMinus->Draw("mean:phi",cut+"theta<-0.1","goff");
371 graphTRDYC[2] = new TGraph(entries, ttrdDyMinus->GetV2(),ttrdDyMinus->GetV1());
372 ttrdDyPlus->Draw("mean:phi",cut+"theta<-0.1","goff");
373 graphTRDYC[1] = new TGraph(entries, ttrdDyPlus->GetV2(),ttrdDyPlus->GetV1());
374 //
375 entries=ttrdDyPlus->Draw("mean-TM..mean:phi",cut+"theta>0.1&&TM..entries>100","goff");
376 graphTRDYA[0] = new TGraph(entries, ttrdDyPlus->GetV2(),ttrdDyPlus->GetV1());
377 ttrdDyMinus->Draw("mean:phi",cut+"theta>0.1","goff");
378 graphTRDYA[2] = new TGraph(entries, ttrdDyMinus->GetV2(),ttrdDyMinus->GetV1());
379 ttrdDyPlus->Draw("mean:phi",cut+"theta>0.1","goff");
380 graphTRDYA[1] = new TGraph(entries, ttrdDyPlus->GetV2(),ttrdDyPlus->GetV1());
381 //
382
383
384 //
385 {for (Int_t i=0; i<3; i++){
386 graphITSYC[i]->SetMaximum(ymax+0.2); graphITSYC[i]->SetMinimum(ymin);
387 graphITSYC[i]->GetXaxis()->SetTitle("#phi"); graphITSYC[i]->GetYaxis()->SetTitle("#Delta r#phi (cm)");
388 graphITSYC[i]->SetMarkerColor(i+1); graphITSYC[i]->SetMarkerStyle(20+i);
389 graphITSYA[i]->SetMaximum(ymax+0.2); graphITSYA[i]->SetMinimum(ymin);
390 graphITSYA[i]->GetXaxis()->SetTitle("#phi"); graphITSYA[i]->GetYaxis()->SetTitle("#Delta r#phi (cm)");
391 graphITSYA[i]->SetMarkerColor(i+1); graphITSYA[i]->SetMarkerStyle(20+i);
392 graphTRDYC[i]->SetMaximum(ymax+0.2); graphTRDYC[i]->SetMinimum(ymin);
393 graphTRDYC[i]->GetXaxis()->SetTitle("#phi"); graphTRDYC[i]->GetYaxis()->SetTitle("#Delta r#phi (cm)");
394 graphTRDYC[i]->SetMarkerColor(i+1); graphTRDYC[i]->SetMarkerStyle(20+i);
395 graphTRDYA[i]->SetMaximum(ymax+0.2); graphTRDYA[i]->SetMinimum(ymin);
396 graphTRDYA[i]->GetXaxis()->SetTitle("#phi"); graphTRDYA[i]->GetYaxis()->SetTitle("#Delta r#phi (cm)");
397 graphTRDYA[i]->SetMarkerColor(i+1); graphTRDYA[i]->SetMarkerStyle(20+i);
398 }
399 }
400 TLegend *legend=0;
401 TString cname="cdeltaY";
402 TString dirName=gSystem->GetFromPipe("dirname `pwd` | xargs basename ");
403 TString splus=gSystem->GetFromPipe("ls -al meanBplus.root | gawk '{print \"File: \"$8\" \"$6\" \"$7\" \";}'");
404 TString sminus=gSystem->GetFromPipe("ls -al meanBminus.root | gawk '{print \"File: \"$8\" \"$6\" \"$7\" \";}'");
405
406 TCanvas *cdeltaY = new TCanvas("cdeltaY","cdeltaY",1300,800);
407 cdeltaY->Divide(2,2);
408
409 cdeltaY->cd(1);
410 graphITSYC[0]->Draw("ap");
411 graphITSYC[1]->Draw("p");
412 graphITSYC[2]->Draw("p");
413 legend = new TLegend(0.6,0.6,1.0,1.0,"ITS-TPC C side #Delta r-#phi");
414 legend->AddEntry("",dirName.Data());
415 legend->AddEntry("",splus.Data());
416 legend->AddEntry("",sminus.Data());
417 legend->AddEntry(graphITSYC[0],"B_{plus}-B_{minus}");
418 legend->AddEntry(graphITSYC[1],"B_{plus}");
419 legend->AddEntry(graphITSYC[2],"B_{minus}");
420 legend->Draw();
421
422 cdeltaY->cd(2);
423 graphITSYA[0]->Draw("ap");
424 graphITSYA[1]->Draw("p");
425 graphITSYA[2]->Draw("p");
426 legend = new TLegend(0.6,0.6,1.0,1.0,"ITS-TPC A side #Delta r-#phi");
427 legend->AddEntry(graphITSYA[0],"B_{plus}-B_{minus}");
428 legend->AddEntry(graphITSYA[1],"B_{plus}");
429 legend->AddEntry(graphITSYA[2],"B_{minus}");
430 legend->Draw();
431
432 //
433 cdeltaY->cd(3);
434 graphTRDYC[0]->Draw("ap");
435 graphTRDYC[1]->Draw("p");
436 graphTRDYC[2]->Draw("p");
437 legend = new TLegend(0.6,0.6,1.0,1.0,"TRD-TPC C side #Delta r-#phi");
438 legend->AddEntry(graphTRDYC[0],"B_{plus}-B_{minus}");
439 legend->AddEntry(graphTRDYC[1],"B_{plus}");
440 legend->AddEntry(graphTRDYC[2],"B_{minus}");
441 legend->Draw();
442
443 cdeltaY->cd(4);
444 graphTRDYA[0]->Draw("ap");
445 graphTRDYA[1]->Draw("p");
446 graphTRDYA[2]->Draw("p");
447 legend = new TLegend(0.6,0.6,1.0,1.0,"TRD-TPC A side #Delta r-#phi");
448 legend->AddEntry(graphTRDYA[0],"B_{plus}-B_{minus}");
449 legend->AddEntry(graphTRDYA[1],"B_{plus}");
450 legend->AddEntry(graphTRDYA[2],"B_{minus}");
451 legend->Draw();
452 cdeltaY->SaveAs(dirName+"_deltaY.pdf");
453
454}
455
456
ddd07cdf 457
7f4cb119 458
ddd07cdf 459
460
461void MakeFitTree(){
462 //
463 // 1. Initialize ocdb e.g
464 // .x $ALICE_ROOT/ANALYSIS/CalibMacros/Pass0/ConfigCalibTrain.C(114972,)
465 // .x $ALICE_ROOT/ANALYSIS/CalibMacros/Pass0/ConfigCalibTrain.C(114972,"local:///lustre/alice/alien/alice/data/2010/OCDB")
466 //
467 //
ddd07cdf 468
7f4cb119 469 AliTPCComposedCorrection *cc= MakeComposedCorrection();
470 TObjArray * corr = (TObjArray*)(cc->GetCorrections());
471 // corr->AddLast(cc);
0b736a46 472 const Int_t kskip=23;
ddd07cdf 473 //
474 TFile f("mean.root");
475 TTree * tree= 0;
476 tree = (TTree*)f.Get("ITSdy");
7f4cb119 477 printf("ITSdy\n");
478 AliTPCCorrection::MakeTrackDistortionTree(tree,0,0,corr,kskip);
ddd07cdf 479 tree = (TTree*)f.Get("TRDdy");
7f4cb119 480 printf("TRDdy\n");
481 AliTPCCorrection::MakeTrackDistortionTree(tree,1,0,corr,kskip);
ddd07cdf 482 tree = (TTree*)f.Get("Vertexdy");
7f4cb119 483 printf("Vertexdy\n");
484 AliTPCCorrection::MakeTrackDistortionTree(tree,2,0,corr,kskip);
485 tree = (TTree*)f.Get("TOFdy");
486 printf("TOFdy\n");
487 if (tree) AliTPCCorrection::MakeTrackDistortionTree(tree,3,0,corr,kskip);
488
ddd07cdf 489 //
490 tree = (TTree*)f.Get("ITSdsnp");
7f4cb119 491 printf("ITSdsnp\n");
492 AliTPCCorrection::MakeTrackDistortionTree(tree,0,2,corr,kskip);
ddd07cdf 493 tree = (TTree*)f.Get("TRDdsnp");
7f4cb119 494 printf("TRDdsnp\n");
495 AliTPCCorrection::MakeTrackDistortionTree(tree,1,2,corr,kskip);
ddd07cdf 496 tree = (TTree*)f.Get("Vertexdsnp");
7f4cb119 497 printf("Vertexdsnp\n");
498 AliTPCCorrection::MakeTrackDistortionTree(tree,2,2,corr,kskip);
ddd07cdf 499 //
7f4cb119 500
501 tree = (TTree*)f.Get("ITSd1pt");
502 printf("ITSd1pt\n");
503 if (tree) AliTPCCorrection::MakeTrackDistortionTree(tree,0,4,corr,kskip);
504
505 tree = (TTree*)f.Get("TOFdz");
506 printf("TOFdz\n");
507 if (tree) AliTPCCorrection::MakeTrackDistortionTree(tree,3,0,corr,kskip);
508
509
ddd07cdf 510 tree = (TTree*)f.Get("ITSdz");
7f4cb119 511 printf("ITSdz\n");
512 AliTPCCorrection::MakeTrackDistortionTree(tree,0,1,corr,kskip);
ddd07cdf 513 tree = (TTree*)f.Get("Vertexdz");
7f4cb119 514 printf("Vertexdz\n");
515 AliTPCCorrection::MakeTrackDistortionTree(tree,2,1,corr,kskip);
ddd07cdf 516
517 //
518 tree = (TTree*)f.Get("ITSdtheta");
7f4cb119 519 printf("ITSdtheta\n");
520 AliTPCCorrection::MakeTrackDistortionTree(tree,0,3,corr,kskip);
ddd07cdf 521 //
522 tree = (TTree*)f.Get("Vertexdtheta");
7f4cb119 523 printf("Vertexdtheta\n");
524 AliTPCCorrection::MakeTrackDistortionTree(tree,2,3,corr,kskip);
ddd07cdf 525
526
527}
528
529
530
531void MakeGlobalFit(){
532 AliXRDPROOFtoolkit tool;
533 TChain *chain = tool.MakeChain("distortion.txt","fit",0,100000);
534
535 TCut cutS="rms>0&&entries>100";
536 TCut cut=cutS+"(ptype==0||ptype==2||ptype==3||(ptype==4&&dtype==0))";
537 Int_t npointsMax=30000000;
538 TStatToolkit toolkit;
539 Double_t chi2=0;
540 Int_t npoints=0;
541 TVectorD param;
542 TMatrixD covar;
543 chain->SetAlias("side","(1+(theta>0)*2)");
544 TString fstring=""; // magnetic part
545 fstring+="tX1++"; //1 - twist x in mrad
546 fstring+="tY1++"; //2 - twist y in mrad
547 fstring+="ExBCon++"; //3 - custom
548 fstring+="ExBCon*cos(phi)++"; //3 - custom
549 fstring+="ExBCon*sin(phi)++"; //3 - custom
550 fstring+="ggA++"; //4 - gating grid
551 fstring+="ggA*cos(phi)++"; //4 - gating grid
552 fstring+="ggA*sin(phi)++"; //4 - gating grid
553 fstring+="ggC++"; //5 - gating grid
554 fstring+="ggC*cos(phi)++"; //5 - gating grid
555 fstring+="ggC*sin(phi)++"; //5 - gating grid
556 fstring+="(ptype==2)++";
557 fstring+="(ptype==3)++";
558 fstring+="(ptype==4)++";
559 fstring+="(ptype==4)*bz++";
560
561 TString *strDelta = TStatToolkit::FitPlane(chain,"mean:rms", fstring.Data(),cut+"(Entry$%13==0)", chi2,npoints,param,covar,-1,0, npointsMax, kTRUE);
562 chain->SetAlias("delta",strDelta->Data());
563 strDelta->Tokenize("++")->Print();
564
565}
566
567
568void MakeGlobalFitRelative(Int_t highFrequency=0){
569 //
570 // Make a global fit of ExB
571 // To get rid of the misalignment errors -
572 // Use relative change of deltas for 2 different filed settings
573 //
574
575 AliXRDPROOFtoolkit tool;
7f4cb119 576 // TChain *chain = tool.MakeChain("distortion.txt","fit",0,100000);
ddd07cdf 577 TChain *chainPlus = tool.MakeChain("distortionPlus.txt","fit",0,100000);
578 TChain *chainMinus = tool.MakeChain("distortionMinus.txt","fit",0,100000);
7f4cb119 579 TChain *chain0 = tool.MakeChain("distortionField0.txt","fit",0,100000);
580
ddd07cdf 581 chainPlus->AddFriend(chainMinus,"M");
7f4cb119 582 chainPlus->AddFriend(chain0,"F0");
ddd07cdf 583
584 TCut cut="rms>0&&M.rms>0&&entries>100&&dtype==0&&(ptype==0||ptype==2||ptype==4)";
585 chainPlus->SetAlias("deltaM","mean-M.mean");
586 Int_t npointsMax=30000000;
587 TStatToolkit toolkit;
588 Double_t chi2=0;
589 Int_t npoints=0;
590 TVectorD param;
591 TMatrixD covar;
592 //
593 TString fstring=""; // fit string
594 fstring+="(ptype==0)*cos(phi)++"; // - primitive gx movement
595 fstring+="(ptype==0)*sin(phi)++"; // - primitive gy movement
596
597 fstring+="(tX1-M.tX1)++"; // 1 - twist x in mrad
598 fstring+="(tY1-M.tY1)++"; // 2 - twist y in mrad
599 //
600 fstring+="(ExBConA-M.ExBConA)++"; // 3 conical shape A side
601 fstring+="(ExBConC-M.ExBConC)++"; // 4 conical shape C side
602 //
603 fstring+="(ggA-M.ggA)++"; // 3 conical shape A side
604 fstring+="(ggC-M.ggC)++"; // 4 conical shape C side
605
606 if (highFrequency) {for (Int_t i=1; i<highFrequency; i++){
607 fstring+=Form("((ggA-M.ggA)*cos(%d*phi))++",i);
608 fstring+=Form("((ggA-M.ggA)*sin(%d*phi))++",i);
609 //
610 fstring+=Form("((ggC-M.ggC)*cos(%d*phi))++",i);
611 fstring+=Form("((ggC-M.ggC)*sin(%d*phi))++",i);
612 if (i>1){
613 fstring+=Form("((ExBConA-M.ExBConA)*cos(%d*phi))++",i);
614 fstring+=Form("((ExBConA-M.ExBConA)*sin(%d*phi))++",i);
615 fstring+=Form("((ExBConC-M.ExBConC)*cos(%d*phi))++",i);
616 fstring+=Form("((ExBConC-M.ExBConC)*sin(%d*phi))++",i);
617 }
618 }}
619
620 TString *strDelta = TStatToolkit::FitPlane(chainPlus,"deltaM:rms", fstring.Data(),cut+"(Entry$%7==0)", chi2,npoints,param,covar,-1,0, npointsMax, kTRUE);
621 chainPlus->SetAlias("delta",strDelta->Data());
622 TObjArray *fitArray = strDelta->Tokenize("++");
623 fitArray->Print();
624
625 //
626 // Draw results
627 //
628 TH1::AddDirectory(0);
629 TLatex *lfit=new TLatex;
630
631
632 TCanvas * canvasdY= new TCanvas(Form("deltaY%d",highFrequency),Form("deltaY%d",highFrequency),1200,800);
633 TCanvas * canvasdSnp= new TCanvas(Form("deltaSnp%d",highFrequency),Form("deltaSnp%d",highFrequency),1200,800);
634 TCanvas * canvasd1pt= new TCanvas(Form("delta1pt%d",highFrequency),Form("delta1pt%d",highFrequency),1200,800);
635
636 canvasdY->Divide(3,2);
637 chainPlus->SetMarkerStyle(25);
638 chainPlus->SetMarkerSize(0.5);
639 chainPlus->SetMarkerColor(1);
640 canvasdY->cd(1);
641 chainPlus->Draw("deltaM:delta",cut+"theta>0&&abs(snp)<0.2&&ptype==0","");
642 canvasdY->cd(2);
643 chainPlus->SetMarkerColor(1);
644 chainPlus->Draw("deltaM-delta>>deltaCorr(50,-0.4,0.4)",cut+"theta>0&&abs(snp)<0.2&&ptype==0","err");
645 chainPlus->SetMarkerColor(3);
646 chainPlus->Draw("deltaM",cut+"theta>0&&abs(snp)<0.2&&ptype==0","errsame");
647 canvasdY->cd(3);
648 chainPlus->SetMarkerColor(1);
649 chainPlus->Draw("deltaM:phi",cut+"theta>0&&abs(snp)<0.2&&ptype==0","");
650 chainPlus->SetMarkerColor(2);
651 chainPlus->Draw("deltaM-delta:phi",cut+"theta>0&&abs(snp)<0.2&&ptype==0","same");
652 canvasdY->cd(6);
653 chainPlus->SetMarkerStyle(26);
654 chainPlus->SetMarkerColor(1);
655 chainPlus->Draw("deltaM:phi",cut+"theta<0&&abs(snp)<0.2&&ptype==0","");
656 chainPlus->SetMarkerColor(4);
657 chainPlus->Draw("deltaM-delta:phi",cut+"theta<0&&abs(snp)<0.2&&ptype==0","same");
658
659
660
661 canvasdSnp->Divide(3,2);
662 chainPlus->SetMarkerStyle(25);
663 chainPlus->SetMarkerSize(0.5);
664 chainPlus->SetMarkerColor(1);
665 canvasdSnp->cd(1);
666 chainPlus->Draw("1000*deltaM:1000*delta",cut+"theta>0&&abs(snp)<0.2&&ptype==2","");
667 canvasdSnp->cd(2);
668 chainPlus->SetMarkerColor(1);
669 chainPlus->Draw("1000*(deltaM-delta)>>deltaCorr(50,-4.,4)",cut+"theta>0&&abs(snp)<0.2&&ptype==2","err");
670 chainPlus->SetMarkerColor(3);
671 chainPlus->Draw("1000*deltaM",cut+"theta>0&&abs(snp)<0.2&&ptype==2","errsame");
672 canvasdSnp->cd(3);
673 chainPlus->SetMarkerColor(1);
674 chainPlus->Draw("1000*(deltaM):phi",cut+"theta>0&&abs(snp)<0.2&&ptype==2","");
675 chainPlus->SetMarkerColor(2);
676 chainPlus->Draw("1000*(deltaM-delta):phi",cut+"theta>0&&abs(snp)<0.2&&ptype==2","same");
677 canvasdSnp->cd(6);
678 chainPlus->SetMarkerStyle(26);
679 chainPlus->SetMarkerColor(1);
680 chainPlus->Draw("1000*(deltaM):phi",cut+"theta<0&&abs(snp)<0.2&&ptype==2","");
681 chainPlus->SetMarkerColor(4);
682 chainPlus->Draw("1000*(deltaM-delta):phi",cut+"theta<0&&abs(snp)<0.2&&ptype==2","same");
683
684
685
686 canvasd1pt->Divide(3,2);
687 chainPlus->SetMarkerStyle(25);
688 chainPlus->SetMarkerSize(0.5);
689 chainPlus->SetMarkerColor(1);
690 canvasd1pt->cd(1);
691 chainPlus->Draw("deltaM:delta",cut+"theta>0&&abs(snp)<0.2&&ptype==4","");
692 canvasd1pt->cd(2);
693 chainPlus->SetMarkerColor(1);
694 chainPlus->Draw("deltaM-delta>>deltaCorr(50,-0.15,0.15)",cut+"theta>0&&abs(snp)<0.2&&ptype==4","err");
695 chainPlus->GetHistogram()->Fit("gaus");
696 chainPlus->SetMarkerColor(3);
697 chainPlus->Draw("deltaM",cut+"theta>0&&abs(snp)<0.2&&ptype==4","errsame");
698 canvasd1pt->cd(3);
699 chainPlus->SetMarkerColor(1);
700 chainPlus->Draw("deltaM:phi",cut+"theta>0&&abs(snp)<0.2&&ptype==4","");
701 chainPlus->SetMarkerColor(2);
702 chainPlus->Draw("deltaM-delta:phi",cut+"theta>0&&abs(snp)<0.2&&ptype==4","same");
703 canvasd1pt->cd(6);
704 chainPlus->SetMarkerStyle(26);
705 chainPlus->SetMarkerColor(1);
706 chainPlus->Draw("deltaM:phi",cut+"theta<0&&abs(snp)<0.2&&ptype==4","");
707 chainPlus->SetMarkerColor(4);
708 chainPlus->Draw("deltaM-delta:phi",cut+"theta<0&&abs(snp)<0.2&&ptype==4","same");
709
710 {for (Int_t ipad=0; ipad<3;ipad++){
711 if (ipad==0) canvasdY->cd(5);
712 if (ipad==1) canvasdSnp->cd(5);
713 if (ipad==2) canvasd1pt->cd(5);
714 lfit->SetTextAlign(12); lfit->SetTextSize(0.04);
715 {for (Int_t i=1; i<fitArray->GetEntries(); i++){
716 lfit->DrawLatex(0.1,1-i/9.,fitArray->At(i)->GetName());
717 }}
718 if (ipad==0) canvasdY->cd(4);
719 if (ipad==1) canvasdSnp->cd(4);
720 if (ipad==2) canvasd1pt->cd(4);
721 lfit->DrawLatex(0.1,0.9,"Global fit - TPC ITS matching in r-#phi");
722 lfit->DrawLatex(0.1,0.8,"Residual minimization:");
723 lfit->DrawLatex(0.1,0.7,"r#phi_{0.5T}-r#phi_{-0.5T} (cm)");
724 lfit->DrawLatex(0.1,0.6,"sin(r#phi)_{0.5T}-sin(r#phi)_{-0.5T} (mrad)");
725 lfit->DrawLatex(0.1,0.5,"1/pt_{0.5T}-1/pt_{-0.5T} (1/GeV)");
726 }}
727
728
729 TFile *fpic = new TFile(Form("fitPictures%d.root",highFrequency),"update");
730 canvasd1pt->Write();
731 canvasdY->Write();
732 canvasdSnp->Write();
733 //
734 canvasd1pt->SaveAs(Form("fitd1pt%d.pdf",highFrequency));
735 canvasdY->SaveAs(Form("fitdy%d.pdf",highFrequency));
736 canvasdSnp->SaveAs(Form("fitdsnp%d.pdf",highFrequency));
737 delete fpic;
738}
7f4cb119 739
740void DrawDistortionMaps(const char *fname){
741 TFile f(fname);
742 TTree *ITSdy=(TTree*)f.Get("ITSdy");
743 TTree *ITSdsnp=(TTree*)f.Get("ITSdsnp");
744 TTree *ITSd1pt=(TTree*)f.Get("ITSd1pt");
745 TTree *TRDdy=(TTree*)f.Get("TRDdy");
746 TTree *TRDdsnp=(TTree*)f.Get("TRDdsnp");
747 TTree *TRDd1pt=(TTree*)f.Get("TRDd1pt");
748 TTree *Vertexdy=(TTree*)f.Get("Vertexdy");
749 TTree *Vertexdsnp=(TTree*)f.Get("Vertexdsnp");
750 TTree *Vertexd1pt=(TTree*)f.Get("Vertexd1pt");
751 TCanvas *cits = DrawDistortionMaps(ITSdy,ITSdsnp,ITSd1pt,"ITS","ITS");
752 cits->SaveAs("matchingITS.pdf");
753 TCanvas *ctrd = DrawDistortionMaps(TRDdy,TRDdsnp,TRDd1pt,"TRD","TRD");
754 ctrd->SaveAs("matchingTRD.pdf");
755 TCanvas *cvertex = DrawDistortionMaps(Vertexdy,Vertexdsnp,Vertexd1pt,"Vertex","Vertex");
756 cvertex->SaveAs("matchingVertex.pdf");
757 //
758 //
759 TPostScript *ps = new TPostScript("matching.ps", 112);
760 gStyle->SetOptTitle(1);
761 gStyle->SetOptStat(0);
762 ps->NewPage();
763 cits->Draw();
764 ps->NewPage();
765 ctrd->Draw();
766 ps->NewPage();
767 cvertex->Draw();
768 ps->Close();
769 delete ps;
770}
771
772
773TCanvas * DrawDistortionMaps(TTree * treedy, TTree * treedsnp, TTree* treed1pt, const char * name, const char *title){
774 //
775 //
776 //
777 TH1::AddDirectory(0);
778 TCanvas *cdist = new TCanvas(name,name,1200,800);
779 cdist->Divide(3,2);
780 cdist->Draw("");
781 TH1 * his=0;
782
783 cdist->cd(1);
784 treedy->Draw("10*mean","rms>0&&abs(snp)<0.25&&entries>100","");
785 his = (TH1*)treedy->GetHistogram()->Clone();
786 his->SetName("dY"); his->SetXTitle("#Delta_{r#phi} (cm)");
787 his->Draw();
788 //
789 cdist->cd(2);
790 treedsnp->Draw("1000*mean","rms>0&&abs(snp)<0.25&&entries>200","");
791 his = (TH1*)treedsnp->GetHistogram()->Clone();
792 his->SetName("dsnp"); his->SetXTitle("#Delta_{sin(#phi)}");
793 his->Draw();
794 //
795 cdist->cd(3);
796 treed1pt->Draw("mean","rms>0&&abs(snp)<0.15&&entries>400","");
797 his = (TH1*)treed1pt->GetHistogram()->Clone();
798 his->SetName("d1pt"); his->SetXTitle("#Delta_{1/pt} (1/GeV)");
799 his->Draw();
800 //
801 treedy->SetMarkerSize(1);
802 treedsnp->SetMarkerSize(1);
803 treed1pt->SetMarkerSize(1);
804 {for (Int_t type=0; type<3; type++){
805 cdist->cd(4+type);
806 TTree * tree =treedy;
807 if (type==1) tree=treedsnp;
808 if (type==2) tree=treed1pt;
809 Int_t counter=0;
810 for (Double_t theta=-1; theta<=1; theta+=0.2){
811 if (TMath::Abs(theta)<0.11) continue;
812 TCut cut=Form("rms>0&&abs(snp)<0.25&&entries>20&&abs(theta-%f)<0.1",theta);
813 tree->SetMarkerStyle(20+counter);
814 Option_t *option=(counter==0)? "": "same";
815 if (theta>0) tree->SetMarkerColor(2);
816 if (theta<0) tree->SetMarkerColor(4);
817 if (type==0) tree->Draw("10*mean:phi>>his(45,-pi,pi)",cut,"profgoff");
818 if (type==1) tree->Draw("1000*mean:phi>>his(45,-pi,pi)",cut,"profgoff");
819 if (type==2) tree->Draw("mean:phi>>his(45,-pi,pi)",cut,"profgoff");
820 his = (TH1*)tree->GetHistogram()->Clone();
821 his->SetName(Form("%d_%d",counter,type));
822 his->SetXTitle("#phi");
823 his->SetMaximum(4);
824 his->SetMinimum(-4);
825 if (type==2){
826 his->SetMaximum(0.06);
827 his->SetMinimum(-0.06);
828 }
829 his->Draw(option);
830 counter++;
831 }
832 }
833 }
834 // cdist->SaveAs(Form("%s.pdf"),name);
835 return cdist;
836}
837
838
839
840
841
842AliTPCComposedCorrection * MakeComposedCorrection(){
843
844 //
845 Double_t bzField=AliTrackerBase::GetBz();
846 AliMagF* magF= (AliMagF*)(TGeoGlobalMagField::Instance()->GetField());
847 Double_t vdrift = 2.6; // [cm/us] // to be updated: per second (ideally)
848 Double_t ezField = 400; // [V/cm] // to be updated: never (hopefully)
849 Double_t wt = -10.0 * (bzField*10) * vdrift / ezField ;
0b736a46 850 Double_t T1 = 1.0;
851 Double_t T2 = 1.0;
7f4cb119 852 //
853 //
854
855
856 AliTPCExBTwist *twistX1= new AliTPCExBTwist;
857 twistX1->SetXTwist(0.001); // 1 mrad twist in x
858 twistX1->SetName("tX1");
859 //
860 AliTPCExBTwist *twistY1= new AliTPCExBTwist;
861 twistY1->SetYTwist(0.001); // 1 mrad twist in Y
862 twistY1->SetName("tY1");
863 //
864 AliTPCGGVoltError* ggErrorA = new AliTPCGGVoltError;
865 ggErrorA->SetDeltaVGGA(40.); // delta GG - 1 mm
866 ggErrorA->SetName("ggA");
867 //
868 AliTPCGGVoltError* ggErrorC = new AliTPCGGVoltError;
869 ggErrorC->SetDeltaVGGC(40.); // delta GG - 1 mm
870 ggErrorC->SetName("ggC");
871 // conical free param
872 //
873 AliTPCCalibGlobalMisalignment *shiftX=new AliTPCCalibGlobalMisalignment;
874 shiftX->SetXShift(0.1);
875 shiftX->SetName("shiftX");
876 AliTPCCalibGlobalMisalignment *shiftY=new AliTPCCalibGlobalMisalignment;
877 shiftY->SetYShift(0.1);
878 shiftY->SetName("shiftY");
879
880 AliTPCExBBShape *exb11 = new AliTPCExBBShape;
881 exb11->SetBField(magF);
882 exb11->SetName("exb11");
883 AliTPCExBBShape *exbT1 = new AliTPCExBBShape;
884 exbT1->SetBField(magF);
885 exbT1->SetName("exbT1");
886 AliTPCExBBShape *exbT2 = new AliTPCExBBShape;
887 exbT2->SetBField(magF);
888 exbT2->SetName("exbT2");
889
890 TObjArray * corr = new TObjArray;
891 corr->AddLast(twistX1);
892 corr->AddLast(twistY1);
893 corr->AddLast(ggErrorA);
894 corr->AddLast(ggErrorC);
895 corr->AddLast(shiftX);
896 corr->AddLast(shiftY);
897 corr->AddLast(exb11);
898 corr->AddLast(exbT1);
899 corr->AddLast(exbT2);
900 //
901 AliTPCComposedCorrection *cc= new AliTPCComposedCorrection ;
902 cc->SetCorrections((TObjArray*)(corr));
903 AddEffectiveCorrection(cc);
904 cc->SetOmegaTauT1T2(wt,T1,T2);
905 exbT1->SetOmegaTauT1T2(wt,T1+0.1,T2);
906 exbT2->SetOmegaTauT1T2(wt,T1,T2+0.1);
907 cc->Init();
908 // cc->SetMode(1);
909 cc->Print("DA"); // Print used correction classes
910 cc->SetName("Comp");
911 return cc;
912}
913
914
915
916void AddEffectiveCorrection(AliTPCComposedCorrection* comp){
917 //
918 //
919 //
920 TMatrixD polA(18,4);
921 TMatrixD polC(18,4);
922 TMatrixD valA(18,1);
923 TMatrixD valC(18,1);
924 Int_t counter=0;
925 { for (Int_t iside=0; iside<=1; iside++){
926 {for (Int_t ir=0; ir<3; ir++)
927 for (Int_t id=0; id<3; id++) {
928 polA(counter,0)=0;
929 polA(counter,1)=ir;
930 polA(counter,2)=id;
931 polC(counter,0)=0;
932 polC(counter,1)=ir;
933 polC(counter,2)=id;
934 valA(counter,0)=0;
935 valC(counter,0)=0;
936 counter++;
937 }
938 }
939 }
940 }
941 counter=0;
942 TObjArray * array = (TObjArray*) comp->GetCorrections();
943 { for (Int_t iside=0; iside<=1; iside++)
944 for (Int_t ir=0; ir<3; ir++)
945 for (Int_t id=0; id<3; id++) {
946 AliTPCExBEffective *eff=new AliTPCExBEffective;
947 valA*=0;
948 valC*=0;
949 eff->SetName(Form("eff%d_%d_%d",iside,ir,id));
950 if (iside==0) valA(counter,0)=0.1;
951 if (iside==1) valC(counter,0)=0.1;
952 eff->SetPolynoms(&polA,&polC);
953 eff->SetCoeficients(&valA,&valC);
954 valA.Print();
955 valC.Print();
956 array->AddLast(eff);
957 counter++;
958 }
959 }
960}
961
962
963
964void MakeLaserTree(){
965 //
966 //
967 //
968 AliTPCComposedCorrection *cc= MakeComposedCorrection();
969 TObjArray * corr = (TObjArray*)(cc->GetCorrections());
970 // corr->AddLast(cc);
971 TFile f("laserMean.root");
972 TTree* tree =(TTree*)f.Get("Mean");
973 AliTPCCorrection::MakeLaserDistortionTree(tree,corr,0);
974
975}