Update timestamps for new AMANDA simulation (17/02/2015)
[u/mrichter/AliRoot.git] / TPC / CalibMacros / MakeGlobalFit.C
CommitLineData
7f4cb119 1 /*
2 gROOT->Macro("~/rootlogon.C")
3 gSystem->AddIncludePath("-I$ALICE_ROOT/STAT")
4 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC")
5 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros")
6 gSystem->Load("libANALYSIS");
7 gSystem->Load("libTPCcalib");
8
9 .L $ALICE_ROOT/TPC/CalibMacros/MakeGlobalFit.C+
10 MakeGlobalFit();
11
12 rm matching.ps
13 aliroot -b -q $ALICE_ROOT/TPC/CalibMacros/MakeGlobalFit.C | tee globalFit.log
14
15*/
16
17#if !defined(__CINT__) || defined(__MAKECINT__)
18#include "THnSparse.h"
19#include "TLatex.h"
20#include "TCanvas.h"
21#include "TLegend.h"
22#include "TSystem.h"
23#include "TFile.h"
24#include "TChain.h"
25#include "TCut.h"
26#include "TH3.h"
27#include "TH2F.h"
28#include "TProfile3D.h"
29#include "TMath.h"
30#include "TVectorD.h"
31#include "TMatrixD.h"
32#include "TStatToolkit.h"
33#include "TTreeStream.h"
34#include "AliExternalTrackParam.h"
35#include "AliESDfriend.h"
36#include "AliTPCcalibTime.h"
37#include "TROOT.h"
38#include "AliXRDPROOFtoolkit.h"
39#include "AliTPCCorrection.h"
40#include "AliTPCExBTwist.h"
41#include "AliTPCGGVoltError.h"
42#include "AliTPCComposedCorrection.h"
43#include "AliTPCExBConical.h"
44#include "TPostScript.h"
45#include "TStyle.h"
46#include "AliTrackerBase.h"
47#include "AliTPCCalibGlobalMisalignment.h"
48#include "AliTPCExBEffective.h"
49#include "TEntryList.h"
0b736a46 50#include "AliCDBMetaData.h"
51#include "AliCDBId.h"
52#include "AliCDBManager.h"
53#include "AliCDBStorage.h"
54#include "AliTPCcalibDB.h"
7f4cb119 55#endif
56
57const char* chDetName[5]={"ITS","TRD", "Vertex", "TOF","Laser"};
58const char* chParName[5]={"rphi","z", "snp", "tan","1/pt"};
59const char* chSideName[2]={"A side","C side"};
0b736a46 60Bool_t enableDet[5]={1,1,1,0,1}; // detector for fit ITS=0, TRD=1, Vertex=2, Laser=4
7f4cb119 61Bool_t enableParam[5]={1,0,1,0,1}; //
62Bool_t enableSign=kFALSE;
63Bool_t useEff0=kFALSE;
64Bool_t useEffD=kFALSE;
65Bool_t useEffR=kFALSE;
35d81915 66/// \file MakeGlobalFit.C
67
7f4cb119 68TChain *chain = 0;
69TChain *chainRef = 0;
70Bool_t printMatrix=kFALSE;
71TString *fitString=0;
72//
73// cuts
74//
75TCut cut="1";
76TCut cutD="1";
77
78void MakeChain();
79void MakeAliases();
80void MakeCuts();
81void MakeFit(TCut cutCustom);
0b736a46 82void MakeOCDBEntry(Int_t refRun);
7f4cb119 83TCanvas* DrawFitITS(const char *name);
84TCanvas* DrawFitVertex(const char *name);
85TCanvas* DrawFitLaser(const char *cname);
86TCanvas* DrawCorrdY();
87TCanvas* DrawCorrdSnp();
88TCanvas * DrawFitdY(const char *name);
89TCanvas * DrawFitdSnp(const char *name);
90void PrintMatch();
91TCanvas * MakeComposedCorrection(const char * name);
92
93void MakeAliases(){
94 chain->SetAlias("isITS","(dtype==0||dtype==2)");
95 chain->SetAlias("isTRD","(dtype==1)");
96 chain->SetAlias("isLaser","(dtype==4)");
97 chain->SetAlias("r","sqrt(gx*gx+gy*gy)");
98 chain->SetAlias("r250","(sqrt(gx*gx+gy*gy)/250.)");
99 chain->SetAlias("weight","((dtype==4)*rms*10+rms)"); // downscale laser
100 chain->SetAlias("side","(1+(theta>0)*2)");
101 chain->SetAlias("mdelta","(mean-R.mean-isLaser*((dRrec-R.dRrec)*tan(asin(snp))))");
102 chain->SetAlias("drift","(1-abs(gz/250))");
103 chain->SetAlias("r0","(r-160)/80");
104 chain->SetAlias("delta","(0*gx)");
105}
106
107
108void MakeGlobalFit(){
35d81915 109 ///
110
7f4cb119 111 gROOT->Macro("~/rootlogon.C");
0b736a46 112 //gROOT->Macro("NimStyle.C");
7f4cb119 113 gSystem->AddIncludePath("-I$ALICE_ROOT/STAT");
114 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC");
115 gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
116 gSystem->Load("libANALYSIS");
117 gSystem->Load("libTPCcalib");
118
119 MakeChain();
7f4cb119 120 gStyle->SetOptTitle(1);
121 gStyle->SetOptStat(0);
0b736a46 122 TPostScript *ps = new TPostScript("exbFit.eps", 112);
123 ps->NewPage();
124 TCanvas *c=0;
125 useEff0=kFALSE; useEffD=kFALSE; useEffR=kFALSE;
126 ps->NewPage();
7f4cb119 127 MakeFit("1");
0b736a46 128 //
7f4cb119 129 ps->NewPage();
0b736a46 130 c=DrawFitdY("dY-Physical");
131 c->Update();
7f4cb119 132 ps->NewPage();
0b736a46 133 c->Update();
134
7f4cb119 135 ps->NewPage();
0b736a46 136 c=DrawFitdSnp("dSnp-Physical");
137 c->Update();
138
7f4cb119 139 ps->NewPage();
0b736a46 140 c=DrawFitITS("ITS Physical");
141 c->Update();
142
7f4cb119 143 ps->NewPage();
0b736a46 144 c=DrawFitVertex("Vertex Physical");
145 c->Update();
146
147 ps->NewPage();
148 c=DrawFitLaser("Laser Physical");
149 c->Update();
150
7f4cb119 151 ps->NewPage();
0b736a46 152 c=MakeComposedCorrection("Correction physical");
153 c->Update();
154
7f4cb119 155 //
156 //
157 //
158 //
159 useEff0=kTRUE; useEffD=kTRUE; useEffR=kTRUE; enableSign=kFALSE;
160 ps->NewPage();
161 MakeFit("1");
162 ps->NewPage();
0b736a46 163 c=DrawFitdY("dY-Physical+Effective ");
164 c->Update();
165
7f4cb119 166 ps->NewPage();
0b736a46 167 c=DrawFitdSnp("dSnp-Physical+Effective ");
168 c->Update();
169
7f4cb119 170 ps->NewPage();
0b736a46 171 c=DrawFitITS("ITS Physical+Effective");
172 c->Update();
173
7f4cb119 174 ps->NewPage();
0b736a46 175 c=DrawFitVertex("Vertex Physical+Effective");
176 c->Update();
177
7f4cb119 178 ps->NewPage();
0b736a46 179 c=DrawFitLaser("Laser Physical +Effective ");
180 c->Update();
181
7f4cb119 182 ps->NewPage();
0b736a46 183 c=MakeComposedCorrection("Correction physical+Effective");
184 c->Update();
7f4cb119 185 //
0b736a46 186 useEff0=kTRUE; useEffD=kTRUE; useEffR=kTRUE; enableSign=kTRUE;
7f4cb119 187 ps->NewPage();
188 MakeFit("1");
0b736a46 189 //
7f4cb119 190 ps->NewPage();
0b736a46 191 c=DrawFitdY("dY-Physical+Effective Sign");
192 c->Update();
193
7f4cb119 194 ps->NewPage();
0b736a46 195 c=DrawFitdSnp("dSnp-Physical+Effective Sign");
196 c->Update();
197
7f4cb119 198 ps->NewPage();
0b736a46 199 c=DrawFitITS("ITS Physical+Effective Sign");
200 c->Update();
201
7f4cb119 202 ps->NewPage();
0b736a46 203 c=DrawFitVertex("Vertex Physical+Effective Sign");
204 c->Update();
205
7f4cb119 206 ps->NewPage();
0b736a46 207 c=DrawFitLaser("Laser Physical +Effective Sign");
208 c->Update();
209
7f4cb119 210 ps->NewPage();
0b736a46 211 c=MakeComposedCorrection("Correction physical+Effective Sign");
212 c->Update();
213
7f4cb119 214 ps->Close();
215 delete ps;
216
0b736a46 217 //
7f4cb119 218}
219
220void MakeChain(){
35d81915 221 ///
222
7f4cb119 223 TH1::AddDirectory(0);
0b736a46 224 TFile * f0 =0; // file 0 field
225 TFile * fp =0; // file plus
226 TFile * fm =0; // file minus
227 TTree * tree0=0;
228 TTree * treeP=0;
229 TTree * treeM=0;
7f4cb119 230 //
231 chain = new TChain("fit","fit");
232 chainRef = new TChain("fit","fit");
0b736a46 233 //
234 //
235 //
7f4cb119 236 for (Int_t idet=0; idet<5; idet++){
237 for (Int_t ipar=0; ipar<5; ipar++){
0b736a46 238 f0= TFile::Open(Form("../mergeField0/distortion%d_%d.root",idet,ipar));
239 fp= TFile::Open(Form("../mergePlus/distortion%d_%d.root",idet,ipar));
240 fm= TFile::Open(Form("../mergeMinus/distortion%d_%d.root",idet,ipar));
241 tree0 = (f0) ? (TTree*)f0->Get("fit"):0;
242 treeP = (fp) ? (TTree*)fp->Get("fit"):0;
243 treeM = (fm) ? (TTree*)fm->Get("fit"):0;
244 //
245 if ( ipar==0 || ipar==2){
246 if (tree0 && treeP){
247 chain->Add(Form("../mergePlus/distortion%d_%d.root",idet,ipar));
248 chainRef->Add(Form("../mergeField0/distortion%d_%d.root",idet,ipar));
7f4cb119 249 }
0b736a46 250 if (tree0 && treeM){
251 chain->Add(Form("../mergeMinus/distortion%d_%d.root",idet,ipar));
252 chainRef->Add(Form("../mergeField0/distortion%d_%d.root",idet,ipar));
7f4cb119 253 }
0b736a46 254 }
255 if ( ipar==1 || ipar==3 || ipar==4){
256 if (treeP && treeM){
257 chain->Add(Form("../mergePlus/distortion%d_%d.root",idet,ipar));
258 chainRef->Add(Form("../mergeMinus/distortion%d_%d.root",idet,ipar));
7f4cb119 259 }
260 }
261 }
262 }
263 chain->AddFriend(chainRef,"R");
264 MakeAliases();
265 MakeCuts();
266}
267
268
269void MakeCuts(){
35d81915 270 ///
271
7f4cb119 272 TCut cutS="((rms>0&&R.rms>0&&entries>0&&R.entries>0))"; // statistic cuts
273 TCut cutType="((dtype==R.dtype)&&(ptype==R.ptype))"; // corresponding types
274 TCut cutOut="(ptype==0)*abs(mdelta)<(0.3+rms)||(ptype==0&&abs(mdelta*85)<(0.3+rms*85))"; // corresponding types
275 //
276
277 chain->Draw(">>outList",cutS+cutType+cutOut+"abs(snp)<0.5","entryList");
278 TEntryList *elistOut = (TEntryList*)gDirectory->Get("outList");
279 chain->SetEntryList(elistOut);
280
281
282}
283
284
285TMatrixD * MakeCorrelation(TMatrixD &matrix){
35d81915 286 ///
287
7f4cb119 288 Int_t nrows = matrix.GetNrows();
289 TMatrixD * mat = new TMatrixD(nrows,nrows);
290 for (Int_t irow=0; irow<nrows; irow++)
291 for (Int_t icol=0; icol<nrows; icol++){
292 (*mat)(irow,icol)= matrix(irow,icol)/TMath::Sqrt(matrix(irow,irow)*matrix(icol,icol));
293 }
294 return mat;
295}
296
297void MakeDetCut(){
298 cutD=Form("((dtype==%d)||(dtype==%d)||(dtype==%d)||(dtype==%d)||(dtype==%d))",enableDet[0]?0:0,enableDet[1]?1:0,enableDet[2]?2:0,enableDet[3]?3:0,enableDet[4]?4:0);
299 cutD+=Form("((ptype==%d)||(ptype==%d)||(ptype==%d)||(ptype==%d)||(ptype==%d))",enableParam[0]?0:0,enableParam[1]?1:0,enableParam[2]?2:0,enableParam[3]?3:0,enableParam[4]?4:0);
300}
301
302void MakeFit(TCut cutCustom){
303 MakeDetCut();
304 Int_t npointsMax=30000000;
305 TStatToolkit toolkit;
306 Double_t chi2=0;
307 Int_t npoints=0;
308 TVectorD param;
309 TMatrixD covar;
310 //
311 TString fstring=""; // magnetic part
312 fstring+="(tX1-R.tX1)++"; // twist X
313 fstring+="(tY1-R.tY1)++"; // twist Y
314 fstring+="(sign(bz)*(tX1-R.tX1))++"; // twist X
315 fstring+="(sign(bz)*(tY1-R.tY1))++"; // twist Y
7f4cb119 316
317 {if (enableDet[0]){
318 fstring+="(isITS*shiftX)++"; // shift X - ITS
319 fstring+="(isITS*shiftY)++"; // shift Y - ITS
320 fstring+="(isITS*sign(bz)*shiftX)++"; // shift X - ITS
321 fstring+="(isITS*sign(bz)*shiftY)++"; // shift Y - ITS
322 }}
323 {if (enableDet[1]){
324 fstring+="(isTRD*shiftX)++"; // shift X - TRD
325 fstring+="(isTRD*shiftY)++"; // shift Y - TRD
326 fstring+="(isTRD*sign(bz)*shiftX)++"; // shift X - TRD
327 fstring+="(isTRD*sign(bz)*shiftY)++"; // shift Y - TRD
328 }}
329 //
330 if (enableDet[4]){
331 }
332 TString fstringCustom="";
333 if (useEff0){
0b736a46 334 fstring+="(exbT1-exb11-(R.exbT1-R.exb11))++"; // T1 adjustment
335 fstring+="(exbT2-exb11-(R.exbT2-R.exb11))++"; // T2 adjustment
7f4cb119 336 fstringCustom+="(eff0_0_0-R.eff0_0_0)++"; // effective correction constant part
337 fstringCustom+="(eff1_0_0-R.eff1_0_0)++"; //
338 }
339 if (useEffD){
340 fstringCustom+="(eff0_0_1-R.eff0_0_1)++"; // effective correction drift part
341 fstringCustom+="(eff1_0_1-R.eff1_0_1)++"; //
342 fstringCustom+="(eff0_0_2-R.eff0_0_2)++"; // effective correction drift part
343 fstringCustom+="(eff1_0_2-R.eff1_0_2)++"; //
344 }
345 if (useEffR){
346 fstringCustom+="(eff0_1_0-R.eff0_1_0)++"; // effective correction radial part
347 fstringCustom+="(eff1_1_0-R.eff1_1_0)++"; //
348 fstringCustom+="(eff0_1_1-R.eff0_1_1)++"; // effective correction radial part
349 fstringCustom+="(eff1_1_1-R.eff1_1_1)++"; //
350 fstringCustom+="(eff0_1_2-R.eff0_1_2)++"; // effective correction radial part
351 fstringCustom+="(eff1_1_2-R.eff1_1_2)++"; //
352 fstringCustom+="(eff0_2_0-R.eff0_2_0)++"; // effective correction radial part
353 fstringCustom+="(eff1_2_0-R.eff1_2_0)++"; //
354 fstringCustom+="(eff0_2_1-R.eff0_2_1)++"; // effective correction radial part
355 fstringCustom+="(eff1_2_1-R.eff1_2_1)++"; //
356 fstringCustom+="(eff0_2_2-R.eff0_2_2)++"; // effective correction radial part
357 fstringCustom+="(eff1_2_2-R.eff1_2_2)++"; //
358 }
359 if (enableSign){
360 if (useEff0){
361 fstringCustom+="sign(bz)*(eff0_0_0-R.eff0_0_0)++"; // effective correction constant part
362 fstringCustom+="sign(bz)*(eff1_0_0-R.eff1_0_0)++"; //
363 }
364 if (useEffD){
365 fstringCustom+="sign(bz)*(eff0_0_1-R.eff0_0_1)++"; // effective correction drift part
366 fstringCustom+="sign(bz)*(eff1_0_1-R.eff1_0_1)++"; //
367 fstringCustom+="sign(bz)*(eff0_0_2-R.eff0_0_2)++"; // effective correction drift part
368 fstringCustom+="sign(bz)*(eff1_0_2-R.eff1_0_2)++"; //
369 }
370 if (useEffR){
371 fstringCustom+="sign(bz)*(eff0_1_0-R.eff0_1_0)++"; // effective correction radial part
372 fstringCustom+="sign(bz)*(eff1_1_0-R.eff1_1_0)++"; //
373 fstringCustom+="sign(bz)*(eff0_1_1-R.eff0_1_1)++"; // effective correction radial part
374 fstringCustom+="sign(bz)*(eff1_1_1-R.eff1_1_1)++"; //
375 fstringCustom+="sign(bz)*(eff0_1_2-R.eff0_1_2)++"; // effective correction radial part
376 fstringCustom+="sign(bz)*(eff1_1_2-R.eff1_1_2)++"; //
377 fstringCustom+="sign(bz)*(eff0_2_0-R.eff0_2_0)++"; // effective correction radial part
378 fstringCustom+="sign(bz)*(eff1_2_0-R.eff1_2_0)++"; //
379 fstringCustom+="sign(bz)*(eff0_2_1-R.eff0_2_1)++"; // effective correction radial part
380 fstringCustom+="sign(bz)*(eff1_2_1-R.eff1_2_1)++"; //
381 fstringCustom+="sign(bz)*(eff0_2_2-R.eff0_2_2)++"; // effective correction radial part
382 fstringCustom+="sign(bz)*(eff1_2_2-R.eff1_2_2)++"; //
383 }
384 }
385 //
386 fstring+=fstringCustom;
387 //
388 TString *strDelta = TStatToolkit::FitPlane(chain,"mdelta:weight", fstring.Data(),cut+cutD+cutCustom, chi2,npoints,param,covar,-1,0, npointsMax, kTRUE);
389 if (printMatrix) MakeCorrelation(covar)->Print();
390 chain->SetAlias("delta",strDelta->Data());
391 strDelta->Tokenize("++")->Print();
392 fitString = strDelta;
393 PrintMatch();
394}
395
396
397
398
399
400void PrintMatch(){
35d81915 401 /// Print detector matching info
402
7f4cb119 403 for (Int_t ipar=0; ipar<5; ipar++){
404 for (Int_t idet=0; idet<5; idet++){
405 Double_t mean0,rms0,mean1,rms1;
0b736a46 406 Int_t entries = chain->Draw("delta-mdelta>>rhis",cut+Form("dtype==%d&&ptype==%d",idet,ipar),"goff");
7f4cb119 407 if (entries==0) continue;
408 TH1 * his = (TH1*)(chain->GetHistogram()->Clone());
409 mean1=his->GetMean();
410 rms1 =his->GetRMS();
411 delete his;
412 //
0b736a46 413 entries = chain->Draw("mdelta>>rhis",cut+Form("dtype==%d&&ptype==%d",idet,ipar),"goff");
7f4cb119 414 if (entries==0) continue;
415 his = (TH1*)(chain->GetHistogram()->Clone());
416 //
417 mean0=his->GetMean();
418 rms0 =his->GetRMS();
419
420 printf("ptype==%s\tdtype==%s\tMean=%f -> %f\tRMS=%f -> %f\n",chParName[ipar],chDetName[idet], mean0,mean1,rms0,rms1);
421 delete his;
422 }
423 }
424}
425
426
427
428
429TCanvas* DrawFitITS(const char *name){
35d81915 430 ///
431
7f4cb119 432 TLegend *legend=0;
433 TCanvas *canvas = new TCanvas(name,name,800,800);
434 canvas->Divide(1,2);
435 Int_t entries=0;
436 const char * grname[10]={"A side (B=0.5T)","A side (B=-0.5T)", "C side (B=0.5T)", "C side (B=-0.5T)",0};
437
438 TGraph *graphsdyITS[10];
439 TGraph *graphsdyITSc[10];
440 entries=chain->Draw("mdelta:phi",cut+"ptype==0&&dtype==0&&bz>0&&theta>0","goff");
441 graphsdyITS[0]=new TGraph(entries,chain->GetV2(), chain->GetV1());
442 entries=chain->Draw("mdelta:phi",cut+cut+"ptype==0&&dtype==0&&bz<0&&theta>0","goff");
443 graphsdyITS[1]=new TGraph(entries,chain->GetV2(), chain->GetV1());
444 entries=chain->Draw("mdelta:phi",cut+cut+"ptype==0&&dtype==0&&bz>0&&theta<0","goff");
445 graphsdyITS[2]=new TGraph(entries,chain->GetV2(), chain->GetV1());
446 entries=chain->Draw("mdelta:phi",cut+cut+"ptype==0&&dtype==0&&bz<0&&theta<0","goff");
447 graphsdyITS[3]=new TGraph(entries,chain->GetV2(), chain->GetV1());
448 //
449 entries=chain->Draw("mdelta-delta:phi",cut+cut+"ptype==0&&dtype==0&&bz>0&&theta>0","goff");
450 graphsdyITSc[0]=new TGraph(entries,chain->GetV2(), chain->GetV1());
451 entries=chain->Draw("mdelta-delta:phi",cut+cut+"ptype==0&&dtype==0&&bz<0&&theta>0","goff");
452 graphsdyITSc[1]=new TGraph(entries,chain->GetV2(), chain->GetV1());
453 entries=chain->Draw("mdelta-delta:phi",cut+cut+"ptype==0&&dtype==0&&bz>0&&theta<0","goff");
454 graphsdyITSc[2]=new TGraph(entries,chain->GetV2(), chain->GetV1());
455 entries=chain->Draw("mdelta-delta:phi",cut+cut+"ptype==0&&dtype==0&&bz<0&&theta<0","goff");
456 graphsdyITSc[3]=new TGraph(entries,chain->GetV2(), chain->GetV1());
457
458 canvas->cd(1);
459 legend = new TLegend(0.6,0.6,1.0,1.0,"#Delta R#phi TPC-ITS");
460 for (Int_t i=0; i<4; i++){
461 graphsdyITS[i]->SetMaximum(0.5);
462 graphsdyITS[i]->SetMinimum(-0.5);
463 graphsdyITS[i]->SetMarkerColor(i+1);
464 graphsdyITS[i]->SetMarkerStyle(25+i);
465 // graphsdyITS[i]->SetName(grname[i]);
466 graphsdyITS[i]->GetXaxis()->SetTitle("#phi");
467 graphsdyITS[i]->GetYaxis()->SetTitle("#Delta R#Phi (cm)");
468 if (i==0) graphsdyITS[i]->Draw("ap");
469 graphsdyITS[i]->Draw("p");
470 legend->AddEntry(graphsdyITS[i],grname[i]);
471 }
472 legend->Draw();
473 canvas->cd(2);
474 legend = new TLegend(0.6,0.6,1.0,1.0,"#Delta R#phi TPC-ITS corrected");
475 for (Int_t i=0; i<4; i++){
476 graphsdyITSc[i]->SetMaximum(0.5);
477 graphsdyITSc[i]->SetMinimum(-0.5);
478 graphsdyITSc[i]->SetMarkerColor(i+1);
479 graphsdyITSc[i]->SetMarkerStyle(25+i);
480 // graphsdyITS[i]->SetName(grname[i]);
481 graphsdyITSc[i]->GetXaxis()->SetTitle("#phi");
482 graphsdyITSc[i]->GetYaxis()->SetTitle("#Delta R#Phi (cm)");
483 if (i==0) graphsdyITS[i]->Draw("ap");
484 graphsdyITSc[i]->Draw("p");
485 legend->AddEntry(graphsdyITSc[i],grname[i]);
486 }
487 legend->Draw();
488 return canvas;
489}
490
491
492TCanvas* DrawFitLaser(const char *cname){
35d81915 493 ///
494
7f4cb119 495 TH1::AddDirectory(0);
496 TCut cutLaser=cut+"isLaser&&bz<0";
497 TCanvas *canvas= new TCanvas(cname, cname,800,800);
498 canvas->Divide(2,2);
499 canvas->Draw();
500 TLegend*legend=0;
501 //
502 TH1 *his[16]={0};
503 {for (Int_t icorr=0; icorr<2; icorr++){
504 for (Int_t iside=0; iside<2; iside++){
505 canvas->cd(iside+1);
506 for (Int_t ib=0; ib<4; ib++){
507 TCut cutB="1";
508 if (iside==0) cutB=Form("(id==%d&&gz>0)",ib);
509 if (iside==1) cutB=Form("(id==%d&&gz<0)",ib);
510 //cutB.Print();
511 if (icorr==0) chain->Draw("10*mdelta:r",cutLaser+cutB,"prof");
512 if (icorr==1) chain->Draw("10*(mdelta-delta):r",cutLaser+cutB,"prof");
513 his[icorr*8+iside*4+ib] = (TH1*)(chain->GetHistogram()->Clone());
514 his[icorr*8+iside*4+ib]->SetName(Form("B%d%d%d",icorr,iside,ib));
515 his[icorr*8+iside*4+ib]->SetTitle(Form("Bundle %d",ib));
516 his[icorr*8+iside*4+ib]->SetMarkerColor(ib+1);
517 his[icorr*8+iside*4+ib]->SetMarkerStyle(ib+25);
518 his[icorr*8+iside*4+ib]->SetMarkerSize(0.4);
519 his[icorr*8+iside*4+ib]->SetMaximum(3);
520 his[icorr*8+iside*4+ib]->SetMinimum(-3);
521 his[icorr*8+iside*4+ib]->GetXaxis()->SetTitle("r (cm)");
522 his[icorr*8+iside*4+ib]->GetYaxis()->SetTitle("#Delta r#phi (mm)");
523 }
524 }
525 }}
526 //
527 for (Int_t icorr=0; icorr<2; icorr++){
528 for (Int_t iside=0; iside<2; iside++){
529 canvas->cd(icorr*2+iside+1);
530 legend = new TLegend(0.6,0.6,1.0,1.0,Form("#Delta R#phi Laser-%s",chSideName[iside]));
531 for (Int_t ib=0; ib<4; ib++){
532 if (ib==0) his[icorr*2+iside*4+ib]->Draw();
533 his[icorr*8+iside*4+ib]->Draw("same");
534 legend->AddEntry(his[icorr*8+iside*4+ib]);
535 }
536 legend->Draw();
537 }
538 }
539 return canvas;
540}
541
542
543
544
545
546
547
548
549
550
551
552TCanvas* DrawFitVertex(const char *name){
35d81915 553 ///
554
7f4cb119 555 TLegend *legend=0;
556 TCanvas *canvas = new TCanvas(name,name,800,800);
557 canvas->Divide(1,2);
558 Int_t entries=0;
559 const char * grname[10]={"A side (B=0.5T)","A side (B=-0.5T)", "C side (B=0.5T)", "C side (B=-0.5T)",0};
560
561 TGraph *graphsdyVertex[10];
562 TGraph *graphsdyVertexc[10];
563 entries=chain->Draw("mdelta:phi",cut+cut+"ptype==0&&dtype==2&&bz>0&&theta>0","goff");
564 graphsdyVertex[0]=new TGraph(entries,chain->GetV2(), chain->GetV1());
565 entries=chain->Draw("mdelta:phi",cut+cut+"ptype==0&&dtype==2&&bz<0&&theta>0","goff");
566 graphsdyVertex[1]=new TGraph(entries,chain->GetV2(), chain->GetV1());
567 entries=chain->Draw("mdelta:phi",cut+cut+"ptype==0&&dtype==2&&bz>0&&theta<0","goff");
568 graphsdyVertex[2]=new TGraph(entries,chain->GetV2(), chain->GetV1());
569 entries=chain->Draw("mdelta:phi",cut+cut+"ptype==0&&dtype==2&&bz<0&&theta<0","goff");
570 graphsdyVertex[3]=new TGraph(entries,chain->GetV2(), chain->GetV1());
571 //
572 entries=chain->Draw("mdelta-delta:phi",cut+cut+"ptype==0&&dtype==2&&bz>0&&theta>0","goff");
573 graphsdyVertexc[0]=new TGraph(entries,chain->GetV2(), chain->GetV1());
574 entries=chain->Draw("mdelta-delta:phi",cut+cut+"ptype==0&&dtype==2&&bz<0&&theta>0","goff");
575 graphsdyVertexc[1]=new TGraph(entries,chain->GetV2(), chain->GetV1());
576 entries=chain->Draw("mdelta-delta:phi",cut+cut+"ptype==0&&dtype==2&&bz>0&&theta<0","goff");
577 graphsdyVertexc[2]=new TGraph(entries,chain->GetV2(), chain->GetV1());
578 entries=chain->Draw("mdelta-delta:phi",cut+cut+"ptype==0&&dtype==2&&bz<0&&theta<0","goff");
579 graphsdyVertexc[3]=new TGraph(entries,chain->GetV2(), chain->GetV1());
580
581 canvas->cd(1);
582 legend = new TLegend(0.6,0.6,1.0,1.0,"#Delta R#phi TPC-Vertex");
583 for (Int_t i=0; i<4; i++){
584 graphsdyVertex[i]->SetMaximum(1);
585 graphsdyVertex[i]->SetMinimum(-1);
586 graphsdyVertex[i]->SetMarkerColor(i+1);
587 graphsdyVertex[i]->SetMarkerStyle(25+i);
588 // graphsdyVertex[i]->SetName(grname[i]);
589 graphsdyVertex[i]->GetXaxis()->SetTitle("#phi");
590 graphsdyVertex[i]->GetYaxis()->SetTitle("#Delta R#Phi (cm)");
591 if (i==0) graphsdyVertex[i]->Draw("ap");
592 graphsdyVertex[i]->Draw("p");
593 legend->AddEntry(graphsdyVertex[i],grname[i]);
594 }
595 legend->Draw();
596 canvas->cd(2);
597 legend = new TLegend(0.6,0.6,1.0,1.0,"#Delta R#phi TPC-Vertex corrected");
598 for (Int_t i=0; i<4; i++){
599 graphsdyVertexc[i]->SetMaximum(1);
600 graphsdyVertexc[i]->SetMinimum(-1);
601 graphsdyVertexc[i]->SetMarkerColor(i+1);
602 graphsdyVertexc[i]->SetMarkerStyle(25+i);
603 // graphsdyVertex[i]->SetName(grname[i]);
604 graphsdyVertexc[i]->GetXaxis()->SetTitle("#phi");
605 graphsdyVertexc[i]->GetYaxis()->SetTitle("#Delta R#Phi (cm)");
606 if (i==0) graphsdyVertex[i]->Draw("ap");
607 graphsdyVertexc[i]->Draw("p");
608 legend->AddEntry(graphsdyVertexc[i],grname[i]);
609 }
610 legend->Draw();
611 return canvas;
612}
613
614
615
616
617
618
619TCanvas* DrawCorrdY(){
620
621 TLegend *legend=0;
622 TCanvas *canvas = new TCanvas("Corrections","Corrections",800,800);
623 canvas->Divide(1,3);
624
625 TGraph *graphsITS[10];
626 TGraph *graphsTRD[10];
627 TGraph *graphsVertex[10];
628 Int_t entries=0;
629 const char * grname[10]={"X twist (1 mrad)","Y twist (1 mrad)", "X shift (1 mm)", "Y shift (1 mm)", "drPhi (1 mm)"};
630 //
631 entries=chain->Draw("tX1:phi",cut+cut+"ptype==0&&dtype==0&&bz>0","goff");
632 graphsITS[0]=new TGraph(entries,chain->GetV2(), chain->GetV1());
633 entries=chain->Draw("tY1:phi",cut+cut+"ptype==0&&dtype==0&&bz>0","goff");
634 graphsITS[1]=new TGraph(entries,chain->GetV2(), chain->GetV1());
635 entries=chain->Draw("shiftX:phi",cut+cut+"ptype==0&&dtype==0&&bz>0","goff");
636 graphsITS[2]=new TGraph(entries,chain->GetV2(), chain->GetV1());
637 entries=chain->Draw("shiftY:phi",cut+cut+"ptype==0&&dtype==0&&bz>0","goff");
638 graphsITS[3]=new TGraph(entries,chain->GetV2(), chain->GetV1());
639 entries=chain->Draw("drPhiA+drPhiC:phi",cut+cut+"ptype==0&&dtype==0&&bz>0","goff");
640 graphsITS[4]=new TGraph(entries,chain->GetV2(), chain->GetV1());
641 //
642 entries=chain->Draw("tX1:phi",cut+cut+"ptype==0&&dtype==1&&bz>0","goff");
643 graphsTRD[0]=new TGraph(entries,chain->GetV2(), chain->GetV1());
644 entries=chain->Draw("tY1:phi",cut+cut+"ptype==0&&dtype==1&&bz>0","goff");
645 graphsTRD[1]=new TGraph(entries,chain->GetV2(), chain->GetV1());
646 entries=chain->Draw("shiftX:phi",cut+cut+"ptype==0&&dtype==1&&bz>0","goff");
647 graphsTRD[2]=new TGraph(entries,chain->GetV2(), chain->GetV1());
648 entries=chain->Draw("shiftY:phi",cut+cut+"ptype==0&&dtype==1&&bz>0","goff");
649 graphsTRD[3]=new TGraph(entries,chain->GetV2(), chain->GetV1());
650 entries=chain->Draw("drPhiA+drPhiC:phi",cut+cut+"ptype==0&&dtype==1&&bz>0","goff");
651 graphsTRD[4]=new TGraph(entries,chain->GetV2(), chain->GetV1());
652 //
653 entries=chain->Draw("tX1:phi",cut+cut+"ptype==0&&dtype==2&&bz>0","goff");
654 graphsVertex[0]=new TGraph(entries,chain->GetV2(), chain->GetV1());
655 entries=chain->Draw("tY1:phi",cut+cut+"ptype==0&&dtype==2&&bz>0","goff");
656 graphsVertex[1]=new TGraph(entries,chain->GetV2(), chain->GetV1());
657 entries=chain->Draw("shiftX:phi",cut+cut+"ptype==0&&dtype==2&&bz>0","goff");
658 graphsVertex[2]=new TGraph(entries,chain->GetV2(), chain->GetV1());
659 entries=chain->Draw("shiftY:phi",cut+cut+"ptype==0&&dtype==2&&bz>0","goff");
660 graphsVertex[3]=new TGraph(entries,chain->GetV2(), chain->GetV1());
661 entries=chain->Draw("drPhiA+drPhiC:phi",cut+cut+"ptype==0&&dtype==2&&bz>0","goff");
662 graphsVertex[4]=new TGraph(entries,chain->GetV2(), chain->GetV1());
663
664 canvas->cd(1);
665 //
666 legend = new TLegend(0.6,0.6,1.0,1.0,"#Delta R#Phi TPC-Vertex");
667 for (Int_t i=0; i<5; i++){
668 graphsVertex[i]->SetMaximum(0.2);
669 graphsVertex[i]->SetMinimum(-0.2);
670 graphsVertex[i]->SetMarkerColor(i+1);
671 graphsVertex[i]->SetMarkerStyle(25+i);
672 graphsVertex[i]->SetName(grname[i]);
673 graphsVertex[i]->GetXaxis()->SetTitle("#phi");
674 graphsVertex[i]->GetYaxis()->SetTitle("#DeltaR#Phi (cm)");
675 if (i==0) graphsVertex[i]->Draw("ap");
676 graphsVertex[i]->Draw("p");
677 legend->AddEntry(graphsITS[i],grname[i]);
678 }
679 legend->Draw();
680
681 canvas->cd(2);
682 legend = new TLegend(0.6,0.6,1.0,1.0,"#Delta R#Phi ITS-TPC");
683 for (Int_t i=0; i<5; i++){
684 graphsITS[i]->SetMaximum(0.2);
685 graphsITS[i]->SetMinimum(-0.2);
686 graphsITS[i]->SetMarkerColor(i+1);
687 graphsITS[i]->SetMarkerStyle(25+i);
688 graphsITS[i]->SetName(grname[i]);
689 graphsITS[i]->GetXaxis()->SetTitle("#phi");
690 graphsITS[i]->GetYaxis()->SetTitle("#Delta R#Phi (cm)");
691 if (i==0) graphsITS[i]->Draw("ap");
692 graphsITS[i]->Draw("p");
693 legend->AddEntry(graphsITS[i],grname[i]);
694 }
695 legend->Draw();
696 //
697 canvas->cd(3);
698 legend = new TLegend(0.6,0.6,1.0,1.0,"#Delta R#Phi TPC-TRD");
699 for (Int_t i=0; i<5; i++){
700 graphsTRD[i]->SetMaximum(0.2);
701 graphsTRD[i]->SetMinimum(-0.2);
702 graphsTRD[i]->SetMarkerColor(i+1);
703 graphsTRD[i]->SetMarkerStyle(25+i);
704 graphsTRD[i]->SetName(grname[i]);
705 graphsTRD[i]->GetXaxis()->SetTitle("#phi");
706 graphsTRD[i]->GetYaxis()->SetTitle("#DeltaR#Phi (cm)");
707 if (i==0) graphsTRD[i]->Draw("ap");
708 graphsTRD[i]->Draw("p");
709 legend->AddEntry(graphsITS[i],grname[i]);
710 }
711 legend->Draw();
712 return canvas;
713}
714
715
716TCanvas * DrawCorrdSnp(){
717
718 TLegend *legend=0;
719 TCanvas *canvas = new TCanvas("Corrections dSnp","Corrections dSnp",800,800);
720 canvas->Divide(1,3);
721
722 TGraph *graphsITS[10];
723 TGraph *graphsTRD[10];
724 TGraph *graphsVertex[10];
725 Int_t entries=0;
726 const char * grname[10]={"X twist (1 mrad)","Y twist (1 mrad)", "X shift (1 mm)", "Y shift (1 mm)",0};
727 //
728 entries=chain->Draw("1000*tX1:phi",cut+cut+"ptype==2&&dtype==0&&bz>0","goff");
729 graphsITS[0]=new TGraph(entries,chain->GetV2(), chain->GetV1());
730 entries=chain->Draw("1000*tY1:phi",cut+cut+"ptype==2&&dtype==0&&bz>0","goff");
731 graphsITS[1]=new TGraph(entries,chain->GetV2(), chain->GetV1());
732 entries=chain->Draw("1000*shiftX:phi",cut+cut+"ptype==2&&dtype==0&&bz>0","goff");
733 graphsITS[2]=new TGraph(entries,chain->GetV2(), chain->GetV1());
734 entries=chain->Draw("1000*shiftY:phi",cut+cut+"ptype==2&&dtype==0&&bz>0","goff");
735 graphsITS[3]=new TGraph(entries,chain->GetV2(), chain->GetV1());
736 //
737 entries=chain->Draw("1000*tX1:phi",cut+cut+"ptype==2&&dtype==1&&bz>0","goff");
738 graphsTRD[0]=new TGraph(entries,chain->GetV2(), chain->GetV1());
739 entries=chain->Draw("1000*tY1:phi",cut+cut+"ptype==2&&dtype==1&&bz>0","goff");
740 graphsTRD[1]=new TGraph(entries,chain->GetV2(), chain->GetV1());
741 entries=chain->Draw("1000*shiftX:phi",cut+cut+"ptype==2&&dtype==1&&bz>0","goff");
742 graphsTRD[2]=new TGraph(entries,chain->GetV2(), chain->GetV1());
743 entries=chain->Draw("1000*shiftY:phi",cut+cut+"ptype==2&&dtype==1&&bz>0","goff");
744 graphsTRD[3]=new TGraph(entries,chain->GetV2(), chain->GetV1());
745 //
746 entries=chain->Draw("1000*tX1:phi",cut+cut+"ptype==2&&dtype==2&&bz>0","goff");
747 graphsVertex[0]=new TGraph(entries,chain->GetV2(), chain->GetV1());
748 entries=chain->Draw("1000*tY1:phi",cut+cut+"ptype==2&&dtype==2&&bz>0","goff");
749 graphsVertex[1]=new TGraph(entries,chain->GetV2(), chain->GetV1());
750 entries=chain->Draw("1000*shiftX:phi",cut+cut+"ptype==2&&dtype==2&&bz>0","goff");
751 graphsVertex[2]=new TGraph(entries,chain->GetV2(), chain->GetV1());
752 entries=chain->Draw("1000*shiftY:phi",cut+cut+"ptype==2&&dtype==2&&bz>0","goff");
753 graphsVertex[3]=new TGraph(entries,chain->GetV2(), chain->GetV1());
754
755 canvas->cd(1);
756 //
757 legend = new TLegend(0.6,0.6,1.0,1.0,"#Delta sin(#phi) TPC-Vertex");
758 for (Int_t i=0; i<4; i++){
759 graphsVertex[i]->SetMaximum(2);
760 graphsVertex[i]->SetMinimum(2);
761 graphsVertex[i]->SetMarkerColor(i+1);
762 graphsVertex[i]->SetMarkerStyle(25+i);
763 graphsVertex[i]->SetName(grname[i]);
764 graphsVertex[i]->GetXaxis()->SetTitle("#phi");
765 graphsVertex[i]->GetYaxis()->SetTitle("#Deltasin(#phi) (mrad)");
766 if (i==0) graphsVertex[i]->Draw("ap");
767 graphsVertex[i]->Draw("p");
768 legend->AddEntry(graphsITS[i],grname[i]);
769 }
770 legend->Draw();
771
772 canvas->cd(2);
773 legend = new TLegend(0.6,0.6,1.0,1.0,"#Delta sin(#phi) ITS-TPC");
774 for (Int_t i=0; i<4; i++){
775 graphsITS[i]->SetMaximum(1);
776 graphsITS[i]->SetMinimum(-1);
777 graphsITS[i]->SetMarkerColor(i+1);
778 graphsITS[i]->SetMarkerStyle(25+i);
779 graphsITS[i]->SetName(grname[i]);
780 graphsITS[i]->GetXaxis()->SetTitle("#phi");
781 graphsITS[i]->GetYaxis()->SetTitle("#Delta sin(#phi) (mrad)");
782 if (i==0) graphsITS[i]->Draw("ap");
783 graphsITS[i]->Draw("p");
784 legend->AddEntry(graphsITS[i],grname[i]);
785 }
786 legend->Draw();
787 //
788 canvas->cd(3);
789 legend = new TLegend(0.6,0.6,1.0,1.0,"#Delta sin(#phi) TPC-TRD");
790 for (Int_t i=0; i<4; i++){
791 graphsTRD[i]->SetMaximum(1);
792 graphsTRD[i]->SetMinimum(-1);
793 graphsTRD[i]->SetMarkerColor(i+1);
794 graphsTRD[i]->SetMarkerStyle(25+i);
795 graphsTRD[i]->SetName(grname[i]);
796 graphsTRD[i]->GetXaxis()->SetTitle("#phi");
797 graphsTRD[i]->GetYaxis()->SetTitle("#Deltasin(#phi) (mrad)");
798 if (i==0) graphsTRD[i]->Draw("ap");
799 graphsTRD[i]->Draw("p");
800 legend->AddEntry(graphsITS[i],grname[i]);
801 }
802 legend->Draw();
803 return canvas;
804 }
805
806
807TCanvas * DrawFitdY(const char *name){
35d81915 808 ///
809
7f4cb119 810 TH1::AddDirectory(0);
811 TCanvas *canvas = new TCanvas(name,name,800,800);
812 canvas->Divide(3,5);
813 for (Int_t idet=0; idet<5; idet++){
814 chain->SetMarkerColor(4);
815 chain->SetMarkerStyle(25);
816 chain->SetMarkerSize(0.3);
817 chain->SetLineColor(2);
818 //
819 canvas->cd(idet*3+1);
820 chain->Draw("mdelta:delta",cut+Form("ptype==0&&dtype==%d",idet),"");
821 //
822 canvas->cd(idet*3+2);
823 chain->SetLineColor(2);
824 chain->Draw("mdelta-delta",cut+Form("ptype==0&&dtype==%d",idet),"");
825 chain->SetLineColor(1);
826 chain->Draw("mdelta",cut+Form("ptype==0&&dtype==%d",idet),"same");
827 //
828 canvas->cd(idet*3+3);
829 chain->SetMarkerColor(1);
830 chain->Draw("mdelta:phi",cut+Form("ptype==0&&dtype==%d&&theta>0&&bz>0",idet),"");
831 chain->SetMarkerColor(2);
832 chain->Draw("mdelta-delta:phi",cut+Form("ptype==0&&dtype==%d&&theta>0&&bz>0",idet),"same");
833
834 }
835 return canvas;
836}
837
838TCanvas * DrawFitdSnp(const char *name){
35d81915 839 ///
840
7f4cb119 841 TH1::AddDirectory(0);
842 TCanvas *canvas = new TCanvas(name,name,800,800);
843 canvas->Divide(3,5);
844 {for (Int_t idet=0; idet<5; idet++){
845 chain->SetMarkerColor(4);
846 chain->SetMarkerStyle(25);
847 chain->SetMarkerSize(0.3);
848 chain->SetLineColor(2);
849 //
850 canvas->cd(idet*3+1);
851 chain->Draw("1000*mdelta:1000*delta",cut+Form("ptype==2&&dtype==%d",idet),"");
852 //
853 canvas->cd(idet*3+2);
854 chain->SetLineColor(2);
855 chain->Draw("1000*(mdelta-delta)",cut+Form("ptype==2&&dtype==%d",idet),"");
856 chain->SetLineColor(1);
857 chain->Draw("1000*mdelta",cut+Form("ptype==2&&dtype==%d",idet),"same");
858 //
859 canvas->cd(idet*3+3);
860 chain->SetMarkerColor(1);
861 chain->Draw("1000*(mdelta):phi",cut+Form("ptype==2&&dtype==%d&&theta>0&&bz>0",idet),"");
862 chain->SetMarkerColor(2);
863 chain->Draw("1000*(mdelta-delta):phi",cut+Form("ptype==2&&dtype==%d&&theta>0&&bz>0",idet),"same");
864 }}
865 return canvas;
866}
867
0b736a46 868
869
870
871
872
7f4cb119 873TCanvas * MakeComposedCorrection(const char *name){
874
875 TString fit = chain->GetAlias("delta");
876 TObjArray * array = fit.Tokenize("++");
877 Int_t nfun=array->GetEntries();
878 Double_t wt = 0.3 ;
0b736a46 879 Double_t T1 = 1.0;
880 Double_t T2 = 1.0;
7f4cb119 881 //
882 // sign independent correction
883 //
884 AliTPCExBEffective *eff = new AliTPCExBEffective;
885 AliTPCCalibGlobalMisalignment *shiftITS = new AliTPCCalibGlobalMisalignment;
886 AliTPCExBTwist *twist = new AliTPCExBTwist;
887 //
888 // sign dependent correction
889 //
890 AliTPCExBEffective *effS = new AliTPCExBEffective;
891 AliTPCCalibGlobalMisalignment *shiftITSS = new AliTPCCalibGlobalMisalignment;
892 AliTPCExBTwist *twistS = new AliTPCExBTwist;
893
894 TMatrixD polA(100,4);
895 TMatrixD polC(100,4);
896 TMatrixD valA(100,1);
897 TMatrixD valC(100,1);
0b736a46 898 TMatrixD valAS(100,1);
899 TMatrixD valCS(100,1);
7f4cb119 900 Int_t counterA=0;
901 Int_t counterC=0;
0b736a46 902 Int_t counterAS=0;
903 Int_t counterCS=0;
7f4cb119 904
905 {
906 for (Int_t i=1; i<nfun;i++){
907 TObjString fitName=array->At(i)->GetName();
908 TObjString fitVal= &(fitName.String()[1+fitName.String().Last('(')]);
909 Double_t value= fitVal.String().Atof();
910 if (fitName.String().Contains("sign")) continue;
911 if (fitName.String().Contains("tX1")){
912 twist->SetXTwist(value*0.001);
913 }
914 if (fitName.String().Contains("tY1")){
915 twist->SetYTwist(value*0.001);
916 }
917
918 if (fitName.String().Contains("shiftX")&&fitName.String().Contains("isITS")){
919 shiftITS->SetXShift(value*0.1);
920 }
921 if (fitName.String().Contains("shiftY")&&fitName.String().Contains("isITS")){
922 shiftITS->SetYShift(value*0.1);
923 }
924
925 if (fitName.String().Contains("eff")){
926 Int_t index=fitName.String().First("_")-1;
927 Int_t side=atoi(&(fitName.String()[index]));
928 Int_t px =atoi(&(fitName.String()[index+2]));
929 Int_t pd =atoi(&(fitName.String()[index+4]));
930 Int_t pp =0;
931 //printf("%s\t%d\t%d\t%d\t%d\t%f\n",fitName.GetName(),side,px,pd,pp, value);
932 if (side==0){
933 polA(counterA,0)=0;
934 polA(counterA,1)=px;
935 polA(counterA,2)=pd;
936 polA(counterA,3)=pp;
937 valA(counterA,0)=value*0.1;
938 counterA++;
939 }
940 if (side==1){
941 polC(counterC,0)=0;
942 polC(counterC,1)=px;
943 polC(counterC,2)=pd;
944 polC(counterC,3)=pp;
945 valC(counterC,0)=value*0.1;
946 counterC++;
947 }
948 }
949 }
950 }
951 polA.ResizeTo(counterA,4);
952 polC.ResizeTo(counterC,4);
953 valA.ResizeTo(counterA,1);
954 valC.ResizeTo(counterC,1);
955 eff->SetPolynoms(&polA,&polC);
956 eff->SetCoeficients(&valA,&valC);
957 eff->SetOmegaTauT1T2(wt,T1,T2);
958 shiftITS->SetOmegaTauT1T2(wt,T1,T2);
959 twist->SetOmegaTauT1T2(wt,T1,T2);
960 //
961 //
962 //
963 {
0b736a46 964 counterAS=0;
965 counterCS=0;
7f4cb119 966 for (Int_t i=1; i<nfun;i++){
967 TObjString fitName=array->At(i)->GetName();
968 TObjString fitVal= &(fitName.String()[1+fitName.String().Last('(')]);
969 if (!fitName.String().Contains("sign")) continue;
970 Double_t value= fitVal.String().Atof();
971 if (fitName.String().Contains("tX1")){
972 twistS->SetXTwist(value*0.001);
973 }
974 if (fitName.String().Contains("tY1")){
975 twistS->SetYTwist(value*0.001);
976 }
977
978 if (fitName.String().Contains("shiftX")&&fitName.String().Contains("isITS")){
979 shiftITSS->SetXShift(value*0.1);
980 }
981 if (fitName.String().Contains("shiftY")&&fitName.String().Contains("isITS")){
982 shiftITSS->SetYShift(value*0.1);
983 }
984
985 if (fitName.String().Contains("eff")){
986 Int_t index=fitName.String().First("_")-1;
987 Int_t side=atoi(&(fitName.String()[index]));
988 Int_t px =atoi(&(fitName.String()[index+2]));
989 Int_t pd =atoi(&(fitName.String()[index+4]));
990 Int_t pp =0;
991 //printf("%s\t%d\t%d\t%d\t%d\t%f\n",fitName.GetName(),side,px,pd,pp, value);
992 if (side==0){
0b736a46 993 polA(counterAS,0)=0;
994 polA(counterAS,1)=px;
995 polA(counterAS,2)=pd;
996 polA(counterAS,3)=pp;
997 valAS(counterAS,0)=value*0.1;
998 counterAS++;
7f4cb119 999 }
1000 if (side==1){
0b736a46 1001 polC(counterCS,0)=0;
1002 polC(counterCS,1)=px;
1003 polC(counterCS,2)=pd;
1004 polC(counterCS,3)=pp;
1005 valCS(counterCS,0)=value*0.1;
1006 counterCS++;
7f4cb119 1007 }
1008 }
1009 }
1010 }
1011 polA.ResizeTo(counterA,4);
1012 polC.ResizeTo(counterC,4);
1013 valA.ResizeTo(counterA,1);
1014 valC.ResizeTo(counterC,1);
1015 effS->SetPolynoms(&polA,&polC);
0b736a46 1016 effS->SetCoeficients(&valAS,&valCS);
7f4cb119 1017 effS->SetOmegaTauT1T2(wt,T1,T2);
1018 shiftITSS->SetOmegaTauT1T2(wt,T1,T2);
1019 twistS->SetOmegaTauT1T2(wt,T1,T2);
0b736a46 1020 //
1021 // Make combined correction
1022 //
7f4cb119 1023
0b736a46 1024 TObjArray * corr0 = new TObjArray;
1025 TObjArray * corrP = new TObjArray;
1026 TObjArray * corrM = new TObjArray;
1027 AliTPCExBEffective *eff0 = new AliTPCExBEffective;
1028 AliTPCCalibGlobalMisalignment *shiftITSP = new AliTPCCalibGlobalMisalignment;
1029 AliTPCExBTwist *twistP = new AliTPCExBTwist;
1030 AliTPCExBEffective *effP = new AliTPCExBEffective;
1031 AliTPCCalibGlobalMisalignment *shiftITSM = new AliTPCCalibGlobalMisalignment;
1032 AliTPCExBTwist *twistM = new AliTPCExBTwist;
1033 AliTPCExBEffective *effM = new AliTPCExBEffective;
1034 //
1035 shiftITSP->SetXShift(shiftITS->GetXShift()+shiftITSS->GetXShift()); // shift due to the B field
1036 shiftITSP->SetYShift(shiftITS->GetXShift()+shiftITSS->GetYShift());
1037 shiftITSM->SetXShift(shiftITS->GetXShift()-shiftITSS->GetXShift());
1038 shiftITSM->SetYShift(shiftITS->GetXShift()-shiftITSS->GetYShift());
1039 //
1040 twistP->SetXTwist(twist->GetXTwist()+twistS->GetXTwist()); // twist between field - both used
1041 twistP->SetYTwist(twist->GetYTwist()+twistS->GetYTwist()); //
1042 twistM->SetXTwist(twist->GetXTwist()-twistS->GetXTwist());
1043 twistM->SetYTwist(twist->GetYTwist()-twistS->GetYTwist());
1044 //
1045 effP->SetPolynoms(&polA,&polC); // effective correction
1046 effP->SetCoeficients(&valA,&valC);
1047 effM->SetPolynoms(&polA,&polC);
1048 effM->SetCoeficients(&valA,&valC);
1049 //
1050 eff0->SetPolynoms((TMatrixD*)&polA,(TMatrixD*)&polC);
1051 eff0->SetCoeficients((TMatrixD*)&valA,(TMatrixD*)&valC);
1052 //
1053 //
1054 //
1055
1056 //
1057 corrP->AddLast(shiftITSP);
1058 corrP->AddLast(twistP);
1059 corrP->AddLast(effP);
1060 corrM->AddLast(shiftITSM);
1061 corrM->AddLast(twistM);
1062 corrM->AddLast(effM);
1063 corr0->AddLast(eff0);
1064 //
7f4cb119 1065
0b736a46 1066 AliTPCComposedCorrection *comp0= new AliTPCComposedCorrection ;
1067 AliTPCComposedCorrection *compP= new AliTPCComposedCorrection ;
1068 AliTPCComposedCorrection *compM= new AliTPCComposedCorrection ;
1069 comp0->SetCorrections((TObjArray*)(corr0));
1070 compP->SetCorrections((TObjArray*)(corrP));
1071 compM->SetCorrections((TObjArray*)(corrM));
1072 //
1073 comp0->SetOmegaTauT1T2(0*wt,T1,T2);
1074 compP->SetOmegaTauT1T2(wt,T1,T2);
1075 compM->SetOmegaTauT1T2(-wt,T1,T2);
1076
1077 TFile f("correctionsExB.root","update");
1078 f.mkdir(name);
1079 f.cd(name);
1080 printf("\nDump correction B=+0.5T:\n");
1081 compP->Print("da");
1082 printf("\nDump correction B=-0.5T:\n");
1083 compM->Print("da");
1084 //
1085 comp0->Write("ExB-Field0");
1086 compP->Write("ExB-Bplus");
1087 compM->Write("ExB-Bminus");
1088 //
7f4cb119 1089 //
1090 TCanvas *c = new TCanvas(name,name,800,800);
0b736a46 1091 c->Divide(4,2);
7f4cb119 1092 c->cd(1);
1093 shiftITS->CreateHistoDRPhiinXY()->Draw("surf2");
1094 c->cd(2);
1095 twist->CreateHistoDRPhiinXY()->Draw("surf2");
1096 c->cd(3);
1097 eff->CreateHistoDRPhiinZR()->Draw("surf2");
7f4cb119 1098 c->cd(5);
0b736a46 1099 shiftITSS->CreateHistoDRPhiinXY()->Draw("surf2");
7f4cb119 1100 c->cd(6);
0b736a46 1101 twistS->CreateHistoDRPhiinXY()->Draw("surf2");
1102 c->cd(7);
7f4cb119 1103 effS->CreateHistoDRPhiinZR()->Draw("surf2");
0b736a46 1104 //
1105 c->cd(4);
1106 compP->CreateHistoDRPhiinZR()->Draw("surf2");
1107 c->cd(8);
1108 compM->CreateHistoDRPhiinZR()->Draw("surf2");
1109
7f4cb119 1110
1111 return c;
1112}
0b736a46 1113
1114
1115void MakeOCDBEntry(Int_t refRun){
35d81915 1116 /// make a Correction OCDB entry
1117 /// take the fit values writen in config file
1118 ///
1119 /// 1. Read previous value used in calibration
1120 /// OCDB has to be initialized before
1121
0b736a46 1122 gROOT->Macro(Form("ConfigCalibTrain.C(%d)",refRun)); // configuring calib db
1123 gROOT->LoadMacro("AddTaskTPCCalib.C");
1124 gROOT->ProcessLine(Form("ConfigOCDB(%d);",refRun));
1125 AliTPCCorrection *corr = AliTPCcalibDB::Instance()->GetTPCComposedCorrection();
1126 //
1127 TFile f("correctionsExB.root");
1128 AliTPCComposedCorrection * corrP0=(AliTPCComposedCorrection *)f.Get("/Correction physical/ExB-Field0");
1129 AliTPCComposedCorrection * corrPP=(AliTPCComposedCorrection *)f.Get("/Correction physical/ExB-Bplus");
1130 AliTPCComposedCorrection * corrPM=(AliTPCComposedCorrection *)f.Get("/Correction physical/ExB-Bminus");
1131 //
1132 AliTPCComposedCorrection * corrPE0=(AliTPCComposedCorrection *)f.Get("/Correction physical+Effective/ExB-Field0");
1133 AliTPCComposedCorrection * corrPEP=(AliTPCComposedCorrection *)f.Get("/Correction physical+Effective/ExB-Bplus");
1134 AliTPCComposedCorrection * corrPEM=(AliTPCComposedCorrection *)f.Get("/Correction physical+Effective/ExB-Bminus");
1135 //
1136 corrP0->SetName("Field0 correction");
1137 corrPP->SetName("FieldP correction");
1138 corrPM->SetName("FieldM correction");
1139 corrPE0->SetName("Field0 correction");
1140 corrPEP->SetName("FieldP correction");
1141 corrPEM->SetName("FieldM correction");
1142
1143 TObjArray corrPhysical;
1144 TObjArray corrPhysicalEffective;
1145 //
1146 // add the base correction
1147 corrP0->GetCorrections()->Add(corr);
1148 corrPP->GetCorrections()->Add(corr);
1149 corrPM->GetCorrections()->Add(corr);
1150 corrPE0->GetCorrections()->Add(corr);
1151 corrPEP->GetCorrections()->Add(corr);
1152 corrPEM->GetCorrections()->Add(corr);
1153 //
1154 corrPhysical.AddLast(corrP0); corrPhysical.AddLast(corrPP); corrPhysical.AddLast(corrPM);
1155 corrPhysicalEffective.AddLast(corrPE0); corrPhysicalEffective.AddLast(corrPEP); corrPhysicalEffective.AddLast(corrPEM);
1156 //
1157 // make OCDB entries
1158 TString userName=gSystem->GetFromPipe("echo $USER");
1159 TString date=gSystem->GetFromPipe("date");
1160 TString ocdbStorage="local:////";
1161 ocdbStorage+=gSystem->GetFromPipe("pwd")+"/OCDB";
1162 //
1163 AliCDBMetaData *metaData= new AliCDBMetaData();
1164 metaData->SetObjectClassName("TObjArray");
1165 metaData->SetResponsible("Marian Ivanov");
1166 metaData->SetBeamPeriod(1);
1167 metaData->SetAliRootVersion("05-25-01"); //root version
1168 metaData->SetComment(Form("Correction calibration. User: %s\n Data%s",userName.Data(),date.Data()));
1169 AliCDBId* id1=NULL;
1170 id1=new AliCDBId("TPC/Calib/Correction", 0, AliCDBRunRange::Infinity());
1171 AliCDBStorage* gStorage = 0;
1172
1173 gStorage=AliCDBManager::Instance()->GetStorage((ocdbStorage+"Physical").Data());
1174 gStorage->Put(&corrPhysical, (*id1), metaData);
1175 gStorage = AliCDBManager::Instance()->GetStorage((ocdbStorage+"PhysicalEffective").Data());
1176 gStorage->Put(&corrPhysicalEffective, (*id1), metaData);
1177}
1178