]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TPC/CalibMacros/AnalyzeLaserCE.C
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Update master to aliroot
[u/mrichter/AliRoot.git] / TPC / CalibMacros / AnalyzeLaserCE.C
CommitLineData
35d81915 1/// \file AnalyzeLaserCE.C
2///
3/// Macro to perform fits of the Laser Central electrode data
4/// Several fit methods implemented
5///
6/// 0. RebuildCE("ce.root","pul.root"); - rebuild data from the scratch
7/// - the data will be storered in file inname
8///
9/// 1. RebuildData() - transform arbitrary layout of the Input data to the internal format
10/// StoreData(); - The data tree expected in file inname (see variable bellow)
11/// StoreTree(); - Modify inname and xxside and tcor in order to transform data
12///
13///
14/// 2. MakeFit(); - Make a fit of the data - already in internal format
15/// StoreData(); - Store
16/// StoreTree();
17///
18/// 3. MakeRes(); - Make the final calibration + combination of different components
19///
20/// 4. LoadViewer(); - Browse the fit parameters
21///
22/// ~~~
23/// .x ~/rootlogon.C
24/// gSystem->AddIncludePath("-I$ALICE_ROOT/TPC -I$ALICE_ROOT/STAT");
25/// gSystem->Load("libSTAT");
26/// .L $ALICE_ROOT/TPC/CalibMacros/AnalyzeLaserCE.C+
27///
28/// // setup aliases
29///
30/// qaside="CE_Q";
31/// taside="CE_T";
32/// raside="CE_RMS";
33/// qcside="CE_Q";
34/// tcside="CE_T";
35/// rcside="CE_RMS";
36/// ~~~
37///
38/// Calibration viewer variables:
39///
40/// Result - resulting correction
41/// out - outlyers not used for fit
42/// tcor - offset specified by user before fitting
43/// timeF1 - sector local fit - plane
44/// timeF2 - sector local fit - parabola
45/// timeIn - input times
46/// qIn - input charge
47/// out - outlyers not used for fit
48/// tcor - offset specified by user before fitting
49/// timeF1 - sector time local fit - plane
50/// timeF2 - sector time local fit - parabola
51/// qF1 - sector q local fit - plane
52/// qF2 - sector q local fit - parabola
53///
54/// fitted values
55/// ffit0 - base fit
56/// ffit1 - adding common shifts - alpha dependendent
57/// ffit2 - adding opposite shifts - alpha dependent
58///
59/// fGXY - global fit parameter - XY
60/// fInOut - global fit parameter - inner-outer sector matching
61/// fLX - global LX dependence
62///
63/// Gloabl fit o consist of
64/// -fGXY~-fLX~-fTL~-fOff~:ffit0~
65///
66///
67/// Control variable - check results
68///
69/// ffit2~-(timeIn~):lx~ - fit value minus input time
70///
71/// result cosntruction:
72/// (timeF2~-ffit2~+fTL~+fInOut~):Result~
73///
74/// timeF2~-Result~:ffit2~-fTL~-fInOut~
9e7da9da 75
9e7da9da 76
2c55475e 77#include <fstream>
9e7da9da 78#include "TString.h"
79#include "TSystem.h"
80#include "TTree.h"
2c55475e 81#include "TEntryList.h"
73e9c7cf 82#include "TCut.h"
9e7da9da 83#include "TStatToolkit.h"
84#include "AliTPCCalibViewer.h"
85#include "AliTPCCalibViewerGUI.h"
86#include "AliTPCPreprocessorOnline.h"
87#include "AliTPCCalibCE.h"
73e9c7cf 88#include "AliTPCCalibPulser.h"
0803cf01 89#include "TStopwatch.h"
9e7da9da 90//
91//Define interesting variables - file names
92//
93char * inname = "treeCE.root"; // input file with tree
94//
95// variable name definition in input tree - change it according the input
96//
73e9c7cf 97TString qaside("CE_Q");
98TString taside("CE_T");
99TString raside("CE_RMS");
100TString qcside("CE_Q");
101TString tcside("CE_T");
102TString rcside("CE_RMS");
103
9e7da9da 104//
105// correction variable - usually Pulser time
106//
73e9c7cf 107TString tcor("-pulCorr");
9e7da9da 108
109//
110char * fname = "treefitCE.root"; // output file with tree
111char * oname = "fitCE.root"; // output file with CalPads fit
112
113//
114//
115// Input CalPads
116//
117AliTPCCalPad *calPadIn = 0; // original time pad
118AliTPCCalPad *calPadF1 = 0; // original time pad - fit plane
119AliTPCCalPad *calPadF2 = 0; // original time pad - fit parabola
120AliTPCCalPad *calPadQIn = 0; // original Q pad
121AliTPCCalPad *calPadQF1 = 0; // original Q pad
122AliTPCCalPad *calPadQF2 = 0; // original Q pad
123
124AliTPCCalPad *calPadCor = 0; // base correction CalPad
125AliTPCCalPad *calPadOut = 0; // outlyer CalPad
126//
127// cuts
128//
0803cf01 129const Float_t tThr=0.3; // max diff in the sector
9e7da9da 130const Float_t qThr0=0.5; // max Q diff in the sector
131const Float_t qThr1=2; // max Q diff in the sector
132
133//
134//
135// fit Cal Pads
136AliTPCCalPad *calPad0 = 0; // global fit 0 - base
137AliTPCCalPad *calPad1 = 0; // global fit 1 - common behavior rotation -A-C
138AliTPCCalPad *calPad2 = 0; // gloabl fit 2 - CE missalign rotation A-C
139//
140AliTPCCalPad *calPadInOut = 0; // misaalign in-out
141AliTPCCalPad *calPadLX = 0; // local x missalign
142AliTPCCalPad *calPadTL = 0; // tan
143AliTPCCalPad *calPadQ = 0; // time (q) correction
144AliTPCCalPad *calPadGXY = 0; // global XY missalign (drift velocity grad)
145AliTPCCalPad *calPadOff = 0; // normalization offset fit
146AliTPCCalPad *calPadRes = 0; // final calibration
147
0803cf01 148TObjString strFit0=""; // string fit 0
149TObjString strFit1=""; // string fit 1
150TObjString strFit2=""; // string fit 2
151TVectorD vecFit0; // parameters fit 0
152TVectorD vecFit1; // parameters fit 1
153TVectorD vecFit2; // parameters fit 2
2c55475e 154TEntryList *elist=0; // event list -slection criterias
9e7da9da 155//
156// working variables
157//
158AliTPCCalibViewerGUI * viewer=0; //viewerGUI
159AliTPCCalibViewer *makePad=0; //viewer
160TTree * tree=0; // working tree
161
162void LoadViewer();
163void RebuildData(); // transform the input data to the fit format
164void MakeFit(); // make fits
73e9c7cf 165void MakeFitPulser(); // make fit for pulser correction data
9e7da9da 166//
167// internal functions
168//
169void MakeAliases0(); // Make Aliases 0 - for user tree
170void MakeAliases1(); // Make Aliases 1 - for default tree
171void LoadData(); // Load data
172void StoreData(); // store current data
173void StoreTree(); // store fit data in the output tree
174
175
0803cf01 176
9e7da9da 177void AnalyzeLaser(){
35d81915 178 ///
179
9e7da9da 180 LoadViewer();
181 MakeAliases1();
182}
183
73e9c7cf 184void MakeFitPulser(){
185 TStatToolkit stat;
186 Int_t npoints;
187 Double_t chi2;
188 TVectorD vec0,vec1,vec2;
189 TMatrixD mat;
190 TString fitvar="P_T.fElements";
191 TCut out("abs(P_T.fElements/P_T_Mean.fElements-1.001)<.002");
192 TCut fitadd("P_T.fElements>446");
193 TString fstring="";
194 fstring+="(sector>=0&&sector<9)++";
195 fstring+="(sector>=9&&sector<18)++";
196 fstring+="(sector>=18&&sector<27)++";
197 fstring+="(sector>=27&&sector<36)++";
198 fstring+="(sector>=36&&sector<45)++";
199 fstring+="(sector>=45&&sector<54)++";
200 fstring+="(sector>=54&&sector<63)++";
201 fstring+="(sector>=63&&sector<72)";
202
203 TString *pulOff =stat.FitPlane(tree,fitvar.Data(),fstring.Data(),(out+fitadd).GetTitle(),chi2,npoints,vec0,mat);
204
205 tree->SetAlias("pul0",pulOff->Data());
206 tree->SetAlias("pulCorr","P_T.fElements-pul0");
207 tree->SetAlias("pulOut",out.GetTitle());
208}
9e7da9da 209
210void MakeFit(){
35d81915 211 ///
212
9e7da9da 213 LoadData();
0803cf01 214 //LoadViewer();
215 //
216 makePad = new AliTPCCalibViewer(fname);
217 tree = makePad->GetTree();
218 MakeAliases1();
219 //
220 // MakeFitPulser();
9e7da9da 221 TStatToolkit stat;
222 Int_t npoints;
223 Double_t chi2;
9e7da9da 224 TMatrixD mat;
225 TString fstring="";
226 //
227 //Basic correction
228 //
229 fstring+="side++"; // offset on 2 different sides //1
230 //fstring+="(1/qp)++"; // Q -threshold effect correction //2
231 //fstring+="(qp)++"; // Q -threshold effect correction //3
232 fstring+="(inn)++"; // inner outer misalign - common //4
233 fstring+="(side*inn)++"; // - opposite //5
234 //
235 fstring+="(gyr)++"; // drift velocity gradient - common //6
236 fstring+="(side*gyr)++"; // - opposite //7
237 fstring+="(gxr)++"; // global x tilt - common //8
238 fstring+="(side*gxr)++"; // - opposite //9
239 //
240 fstring+="tl^2++"; // local phi correction //10
241 //
242 fstring+="(lxr)++"; // zr angle - common //11
243 fstring+="(side*lxr)++"; // - opposite //12
244 fstring+="(inn*lxr)++"; // inner outer angle - common //13
245 fstring+="(side*inn*lxr)++";// - opposite //14
246 fstring+="(lxr^2)++"; // zr second - common //15
247 fstring+="(side*lxr^2)++"; // - opposite //16
0803cf01 248 fstring+="(inn*lxr^2)++"; // zr second - common //15
249 fstring+="(inn*side*lxr^2)++"; // - opposite //16
9e7da9da 250 //
0803cf01 251 printf("Fit0\t start\n");
252 TString *fit0 =stat.FitPlane(tree,"dt",fstring.Data(),"cutF&&cutCE",chi2,npoints,vecFit0,mat,0.9);
9e7da9da 253 tree->SetAlias("f0",fit0->Data());
0803cf01 254 strFit0 = fit0->Data();
255 printf("Fit0\t end\n");
256 printf("Chi2/npoints\t=%f\n",TMath::Sqrt(chi2/npoints));
257 fit0->Tokenize("++")->Print();
9e7da9da 258 //
0803cf01 259 printf("Global tendencies extraction\n");
9e7da9da 260 //
261 // Extract variables
262 //
263 TString tmpstr = fstring;
264 TObjArray *arr = tmpstr.Tokenize("++");
265 TString fitQ("0"); // q correction
266 TString fitLX("0"); // lx correction
267 TString fitInOut("0"); // inner-outer - match
268 TString fitGXY("0"); // global xy fit
269 TString fitOff("0"); // side offsets
270 TString fitTL("0"); // side offsets
271 //
272 fitOff+="+";
0803cf01 273 fitOff+=vecFit0[0];
9e7da9da 274 fitOff+="+side*";
0803cf01 275 fitOff+=vecFit0[1];
9e7da9da 276 {
277 for(Int_t i=0;i<arr->GetEntriesFast();i++){
278 if (!arr->At(i)) continue;
279 TString *fitstr = new TString(arr->At(i)->GetName());
280 //
281 //Bool_t isQ = fitstr->Contains("qp)");
282 Bool_t isRot = fitstr->Contains("sin(")+fitstr->Contains("cos(");
283 Bool_t isLX = fitstr->Contains("lxr");
284 Bool_t isIn = fitstr->Contains("inn");
285 Bool_t isGXY = fitstr->Contains("gxr")+fitstr->Contains("gyr");
286 if (fitstr->Contains("tl^2")){
287 fitTL+="+";
288 fitTL+=(*fitstr)+"*";
0803cf01 289 fitTL+=vecFit0[i+1];
9e7da9da 290 }
291 if (isGXY){
292 fitGXY+="+";
293 fitGXY+=(*fitstr)+"*";
0803cf01 294 fitGXY+=vecFit0[i+1];
9e7da9da 295 }
9e7da9da 296 //
297 if (isLX&&!isRot&&!isIn){
298 fitLX+="+";
299 fitLX+=(*fitstr)+"*";
0803cf01 300 fitLX+=vecFit0[i+1];
9e7da9da 301 }
302 //
303 if (!isRot&&isIn){
304 fitInOut+="+";
305 fitInOut+=(*fitstr)+"*";
0803cf01 306 fitInOut+=vecFit0[i+1];
9e7da9da 307 }
308 }
309 }
310 //
311 tree->SetAlias("fInOut",fitInOut.Data());
312 tree->SetAlias("fLX",fitLX.Data());
313 tree->SetAlias("fGXY",fitGXY.Data());
314 tree->SetAlias("fOff",fitOff.Data());
315 //tree->SetAlias("fQ",fitQ.Data());
316 tree->SetAlias("fTL",fitTL.Data());
0803cf01 317
9e7da9da 318 //
319 //
0803cf01 320 // Common "deformation" tendencies
321 //
322 fstring+="(sin(atan2(gy.fElements,gx.fElements)))++";
323 fstring+="(cos(atan2(gy.fElements,gx.fElements)))++";
324 //
325 fstring+="(sin(atan2(gy.fElements,gx.fElements)*2))++";
326 fstring+="(cos(atan2(gy.fElements,gx.fElements)*2))++";
327 fstring+="(sin(atan2(gy.fElements,gx.fElements)*3))++";
328 fstring+="(cos(atan2(gy.fElements,gx.fElements)*3))++";
329 //
330 fstring+="(sin(atan2(gy.fElements,gx.fElements)*2))*lxr++";
331 fstring+="(cos(atan2(gy.fElements,gx.fElements)*2))*lxr++";
332 fstring+="(sin(atan2(gy.fElements,gx.fElements)*3))*lxr++";
333 fstring+="(cos(atan2(gy.fElements,gx.fElements)*3))*lxr++";
334 //
335
336 TString *fit1 =stat.FitPlane(tree,"dt",fstring.Data(),"cutF&&cutCE",chi2,npoints,vecFit1,mat,0.9);
337 tree->SetAlias("f1",fit1->Data());
338 strFit1 = fit1->Data();
339 printf("Fit1\t end\n");
340 printf("Chi2/npoints\t=%f\n",TMath::Sqrt(chi2/npoints));
341 fit1->Tokenize("++")->Print();
342
343 //
344 // Central electrode "deformation"
345 //
346 fstring+="(side*sin(atan2(gy.fElements,gx.fElements)))++";
347 fstring+="(side*cos(atan2(gy.fElements,gx.fElements)))++";
348 //
349 fstring+="(side*sin(atan2(gy.fElements,gx.fElements)*2))++";
350 fstring+="(side*cos(atan2(gy.fElements,gx.fElements)*2))++";
351 fstring+="(side*sin(atan2(gy.fElements,gx.fElements)*3))++";
352 fstring+="(side*cos(atan2(gy.fElements,gx.fElements)*3))++";
9e7da9da 353 //
0803cf01 354 fstring+="(side*sin(atan2(gy.fElements,gx.fElements)*2))*lxr++";
355 fstring+="(side*cos(atan2(gy.fElements,gx.fElements)*2))*lxr++";
356 fstring+="(side*sin(atan2(gy.fElements,gx.fElements)*3))*lxr++";
357 fstring+="(side*cos(atan2(gy.fElements,gx.fElements)*3))*lxr++";
358
359 TString *fit2 =stat.FitPlane(tree,"dt",fstring.Data(),"cutF&&abs(dt-f0)<0.7&&cutCE",chi2,npoints,vecFit2,mat,0.9);
360 strFit2 = fit2->Data();
361 printf("Fit2\t end\n");
362 printf("Chi2/npoints\t=%f\n",TMath::Sqrt(chi2/npoints));
363 fit2->Tokenize("++")->Print();
364 tree->SetAlias("f2",fit2->Data());
365 //
366 //
367 //
9e7da9da 368 calPad0 = makePad->GetCalPad("f0","1", "ffit0");
369 calPad1 = makePad->GetCalPad("f1","1", "ffit1");
370 calPad2 = makePad->GetCalPad("f2","1", "ffit2");
371 calPadInOut = makePad->GetCalPad("fInOut","1", "fInOut");
372 calPadLX = makePad->GetCalPad("fLX","1", "fLX");
373 calPadTL = makePad->GetCalPad("fTL","1", "fTL");
374 //calPadQ = makePad->GetCalPad("fQ","1", "fQ");
375 calPadGXY = makePad->GetCalPad("fGXY","1", "fGXY");
376 calPadOff = makePad->GetCalPad("fOff","1", "fOff");
377}
378
73e9c7cf 379
380
9e7da9da 381void LoadViewer(){
35d81915 382 /// Load calib Viewer
383
9e7da9da 384 TObjArray * array = AliTPCCalibViewerGUI::ShowGUI(fname);
385 viewer = (AliTPCCalibViewerGUI*)array->At(0);
386 makePad = viewer->GetViewer();
387 tree = viewer->GetViewer()->GetTree();
388 MakeAliases1();
389}
390
391
392
393
394
395
396void RebuildData(){
35d81915 397 /// transform the input data to the fit format
398
0803cf01 399 TStopwatch timer;
9e7da9da 400 makePad = new AliTPCCalibViewer(inname);
401 tree = makePad->GetTree();
0803cf01 402 //
403
404 timer.Print();
405 timer.Continue();
406 printf("-1. MaQkeFitPulser\n");
73e9c7cf 407 MakeFitPulser();
0803cf01 408 timer.Print();
409 timer.Continue();
9e7da9da 410 MakeAliases0(); //
411 //
412 printf("0. GetCalPads\n");
0803cf01 413 timer.Print();
9e7da9da 414 calPadCor = makePad->GetCalPad("tcor","1", "tcor");
0803cf01 415 calPadOut = makePad->GetCalPad("1","!((cutA||cutC)&&abs(ly.fElements/lx.fElements)<0.16)", "out");
416 calPadIn = makePad->GetCalPad("dt","1","timeIn");
417
418 calPadF1 = calPadIn->GlobalFit("timeF1", calPadOut,kFALSE,0,0,0.8); // use robust fit
9e7da9da 419 calPadQIn = makePad->GetCalPad("qa*(sector%36<18)+qc*(sector%36>17)","1","qIn");
420 //
421 //
422 printf("1. Create outlyer map\n");
0803cf01 423 timer.Print();
424 timer.Continue();
9e7da9da 425 //
426 // Update outlyer maps
427 //
428 for (Int_t isector=0;isector<72; isector++){
429 for (UInt_t ich=0;ich<calPadIn->GetCalROC(isector)->GetNchannels();ich++){
430 Float_t val0= calPadIn->GetCalROC(isector)->GetValue(ich);
431 Float_t val1= calPadF1->GetCalROC(isector)->GetValue(ich);
432 if (TMath::Abs(val0-val1)>tThr) calPadOut->GetCalROC(isector)->SetValue(ich,1);
433 }
434 }
0803cf01 435 printf("2. Fit sector parabolas\n");
436 timer.Print();
437 timer.Continue();
438 calPadF1 = calPadIn->GlobalFit("timeF1", calPadOut,kFALSE,0,0,0.9);
439 calPadF2 = calPadIn->GlobalFit("timeF2", calPadOut,kFALSE,1,0,0.9);
440 calPadQF1 = calPadQIn->GlobalFit("qF1", calPadOut,kFALSE,1);
9e7da9da 441 calPadQF2 = calPadQIn->GlobalFit("qF2", calPadOut,kFALSE,1);
442
443 //
444 // Update outlyer maps
445 //
0803cf01 446 printf("3. Update outlyer map\n");
9e7da9da 447 for (Int_t isector=0;isector<72; isector++){
448 for (UInt_t ich=0;ich<calPadIn->GetCalROC(isector)->GetNchannels();ich++){
449 Float_t val0= calPadQIn->GetCalROC(isector)->GetValue(ich);
450 Float_t val1= calPadQF2->GetCalROC(isector)->GetValue(ich);
0803cf01 451 if (val1<=0) {
9e7da9da 452 continue;
453 }
454 if (TMath::Abs(val0/val1)<qThr0) calPadOut->GetCalROC(isector)->SetValue(ich,1);
455 if (TMath::Abs(val0/val1)>qThr1) calPadOut->GetCalROC(isector)->SetValue(ich,1);
456 }
457 }
0803cf01 458 printf("4. Redo fit of the of parabola \n");
459 timer.Print();
460 timer.Continue();
73e9c7cf 461 //
462 //
463 AliTPCCalPad *calPadInCor = makePad->GetCalPad("dt","1","timeIn");
0803cf01 464 calPadF1 = calPadInCor->GlobalFit("timeF1", calPadOut,kFALSE,0,0.9);
465 calPadF2 = calPadInCor->GlobalFit("timeF2", calPadOut,kFALSE,1,0.9);
466 calPadQF1 = calPadQIn->GlobalFit("qF1", calPadOut,kFALSE,0,0,0.9);
467 calPadQF2 = calPadQIn->GlobalFit("qF2", calPadOut,kFALSE,1,0,0.9);
468 printf("5. End\n");
469 timer.Print();
470
9e7da9da 471}
472
473void LoadData(){
35d81915 474 /// Get Data
475
9e7da9da 476 TFile f(oname);
477 calPadIn = (AliTPCCalPad*)f.Get("timeIn"); // original time pad
478 calPadF1 = (AliTPCCalPad*)f.Get("timeF1"); // original time pad - fit plane
479 calPadF2 = (AliTPCCalPad*)f.Get("timeF2"); // original time pad - fit parabola
480 //
481 calPadQIn = (AliTPCCalPad*)f.Get("qIn"); // original time pad
482 calPadQF1 = (AliTPCCalPad*)f.Get("qF1"); // original time pad - fit plane
483 calPadQF2 = (AliTPCCalPad*)f.Get("qF2"); // original time pad - fit parabola
484 //
485 calPadCor = (AliTPCCalPad*)f.Get("tcor"); // base correction CalPad
486 calPadOut = (AliTPCCalPad*)f.Get("out"); // outlyer CalPad
487 //
488 calPad0 = (AliTPCCalPad*)f.Get("ffit0"); // global fit 0 - base
489 calPad1 = (AliTPCCalPad*)f.Get("ffit1"); // global fit 1 - common behavior rotation -A-C
490 calPad2 = (AliTPCCalPad*)f.Get("ffit2"); // gloabl fit 2 - CE missalign rotation A-C
491 calPadInOut = (AliTPCCalPad*)f.Get("fInOut");// misaalign in-out
492 calPadLX = (AliTPCCalPad*)f.Get("fLX"); // local x missalign
493 calPadTL = (AliTPCCalPad*)f.Get("fTL"); // local y/x missalign
494 calPadQ = (AliTPCCalPad*)f.Get("fQ"); // time (q) correction
495 calPadGXY = (AliTPCCalPad*)f.Get("fGXY"); // global XY missalign (drift velocity grad)
496 calPadOff = (AliTPCCalPad*)f.Get("fOff"); // normalization offset fit
497 calPadRes = (AliTPCCalPad*)f.Get("Result"); //resulting calibration
498}
499
500void StoreData(){
35d81915 501 /// Store data
502
9e7da9da 503 TFile * fstore = new TFile(oname,"recreate");
504 if (calPadIn) calPadIn->Write("timeIn"); // original time pad
505 if (calPadF1) calPadF1->Write("timeF1"); // original time pad - fit plane
506 if (calPadF2) calPadF2->Write("timeF2"); // original time pad - fit parabola
507 //
508 if (calPadQIn) calPadQIn->Write("qIn"); // original time pad
509 if (calPadQF1) calPadQF1->Write("qF1"); // original time pad - fit plane
510 if (calPadQF2) calPadQF2->Write("qF2"); // original time pad - fit parabola
511 //
512 if (calPadCor) calPadCor->Write("tcor"); // base correction CalPad
513 if (calPadOut) calPadOut->Write("out"); // outlyer CalPad
514 //
515 if (calPad0) calPad0->Write("ffit0"); // global fit 0 - base
516 if (calPad1) calPad1->Write("ffit1"); // global fit 1 - common behavior rotation -A-C
517 if (calPad2) calPad2->Write("ffit2"); // gloabl fit 2 - CE missalign rotation A-C
518 if (calPadInOut)calPadInOut->Write("fInOut"); // misaalign in-out
519 if (calPadLX) calPadLX->Write("fLX"); // local x missalign
520 if (calPadTL) calPadTL->Write("fTL"); // local y/x missalign
521 if (calPadQ) calPadQ->Write("fQ"); // time (q) correction
522 if (calPadGXY) calPadGXY->Write("fGXY"); // global XY missalign (drift velocity grad)
523 if (calPadOff) calPadOff->Write("fOff"); // normalization offset fit
524 if (calPadRes) calPadRes->Write("Result"); //resulting calibration
525 fstore->Close();
526 delete fstore;
527}
528
529void StoreTree(){
35d81915 530 ///
531
9e7da9da 532 AliTPCPreprocessorOnline * preprocesor = new AliTPCPreprocessorOnline;
533 //
534 if (calPadIn) preprocesor->AddComponent(calPadIn->Clone());
535 if (calPadF1) preprocesor->AddComponent(calPadF1->Clone());
536 if (calPadF2) preprocesor->AddComponent(calPadF2->Clone());
537 //
538 if (calPadQIn) preprocesor->AddComponent(calPadQIn->Clone());
539 if (calPadQF1) preprocesor->AddComponent(calPadQF1->Clone());
540 if (calPadQF2) preprocesor->AddComponent(calPadQF2->Clone());
541 //
542 if (calPadCor) preprocesor->AddComponent(calPadCor->Clone());
543 if (calPadOut) preprocesor->AddComponent(calPadOut->Clone());
544 //
545 if (calPad0) preprocesor->AddComponent(calPad0->Clone());
546 if (calPad1) preprocesor->AddComponent(calPad1->Clone());
547 if (calPad2) preprocesor->AddComponent(calPad2->Clone());
548 if (calPadInOut)preprocesor->AddComponent(calPadInOut->Clone());
549 if (calPadLX) preprocesor->AddComponent(calPadLX->Clone());
550 if (calPadTL) preprocesor->AddComponent(calPadTL->Clone());
551 if (calPadQ) preprocesor->AddComponent(calPadQ->Clone());
552 if (calPadGXY) preprocesor->AddComponent(calPadGXY->Clone());
553 if (calPadOff) preprocesor->AddComponent(calPadOff->Clone());
554 if (calPadRes) preprocesor->AddComponent(calPadRes->Clone());
555 preprocesor->DumpToFile(fname);
556 delete preprocesor;
557}
558
559
560void MakeAliases0(){
35d81915 561 /// Define variables and selection of outliers - for user defined tree
562
9e7da9da 563 tree->SetAlias("tcor",tcor.Data()); // correction variable
564 tree->SetAlias("ta",taside+".fElements");
565 tree->SetAlias("tc",tcside+".fElements");
566 tree->SetAlias("qa",qaside+".fElements");
567 tree->SetAlias("qc",qcside+".fElements");
568 tree->SetAlias("ra",raside+".fElements");
569 tree->SetAlias("rc",rcside+".fElements");
570 tree->SetAlias("side","1-(sector%36>17)*2");
571 tree->SetAlias("dt","(ta)*(sector%36<18)+(tc)*(sector%36>17)+tcor");
572 tree->SetAlias("cutA","qa>30&&qa<400&&abs(ta)<2&&ra>0.5&&ra<2");
573 tree->SetAlias("cutC","qc>30&&qc<400&&abs(tc)<2&&rc>0.5&&rc<2");
0803cf01 574 tree->SetAlias("cutF","((row.fElements+pad.fElements+sector)%19==0)"); // use just part of the statistic - 5 %
9e7da9da 575 tree->SetAlias("cutCE","V.out.fElements");
576 //
577 // fit param aliases
578 //
579 tree->SetAlias("inn","sector<36");
580 tree->SetAlias("gxr","(gx.fElements/250.)"); //
581 tree->SetAlias("gyr","(gy.fElements/250.)"); //
582 tree->SetAlias("lxr","(lx.fElements-133.41)/250.");
583 tree->SetAlias("qp","((sector%36<18)*sqrt(qa)/10.+(sector%36>17)*sqrt(qc)/10.)"); //
584 tree->SetAlias("tl","(ly.fElements/lx.fElements)/0.17");
585}
586
587
588void MakeAliases1(){
35d81915 589 /// Define variables and selection of outliers -for default usage
590
9e7da9da 591 tree->SetAlias("tcor","tcor.fElements"); // correction variable
592 tree->SetAlias("side","1-(sector%36>17)*2");
0803cf01 593 tree->SetAlias("dt","timeIn.fElements");
9e7da9da 594 //
595 tree->SetAlias("cutA","out.fElements==1");
596 tree->SetAlias("cutC","out.fElements==1");
0803cf01 597 tree->SetAlias("cutF","((row.fElements+pad.fElements+sector.fElements)%19==0)"); // use just part of the statistic - 5 %
598
9e7da9da 599 tree->SetAlias("cutCE","out.fElements<0.5");
600 //
601 // fit param aliases
602 //
603 tree->SetAlias("inn","sector<36");
604 tree->SetAlias("gxr","(gx.fElements/250.)"); //
605 tree->SetAlias("gyr","(gy.fElements/250.)"); //
606 tree->SetAlias("lxr","(lx.fElements-133.41)/250.");
607 tree->SetAlias("qp","(sqrt(qIn.fElements)/10.+(out.fElements>0.5))"); //
608 tree->SetAlias("tl","(ly.fElements/lx.fElements)/0.17");
609}
610
611
612void MakeRes()
613{
35d81915 614 /// make final calibration
615
9e7da9da 616 AliTPCCalPad * calPadRes0 =new AliTPCCalPad(*calPadIn);
9e7da9da 617 calPadRes0->Add(calPad2,-1); // remove global fit
618 calPadRes = calPadRes0->GlobalFit("Result", calPadOut,kTRUE,1,0.9);
619 //
620 //
621 {
622 Float_t tlmedian = calPadTL->GetMedian();
623 for (Int_t isector=0;isector<72; isector++){
624 for (UInt_t ich=0;ich<calPadIn->GetCalROC(isector)->GetNchannels();ich++){
625 //
626 //
627 Float_t val0 = calPadRes->GetCalROC(isector)->GetValue(ich);
628 if (TMath::Abs(val0)>0.5) calPadRes->GetCalROC(isector)->SetValue(ich,0);
629 Float_t tl = calPadTL->GetCalROC(isector)->GetValue(ich);
630 Float_t inOut = calPadInOut->GetCalROC(isector)->GetValue(ich);
631 calPadRes->GetCalROC(isector)->SetValue(ich,calPadRes->GetCalROC(isector)->GetValue(ich)+tl+inOut);
632 //
633 }
634 }
635 }
636 calPadRes->Add(calPadCor,-1); // remove back correction (e.g Pulser time 0)
9e7da9da 637}
638
639
640
641
73e9c7cf 642void RebuildCE(char *finname, char *pulname){
35d81915 643 /// Transformation from the CE to the visualization-analisys output
644 ///
645 /// finname = CE_Vscan_Run_61684-50_170.root;
646
9e7da9da 647 TFile f(finname);
648 AliTPCCalibCE * ce = (AliTPCCalibCE*)f.Get("AliTPCCalibCE");
649 //
650 AliTPCCalPad *padtime = new AliTPCCalPad((TObjArray*)ce->GetCalPadT0());
651 AliTPCCalPad *padRMS = new AliTPCCalPad((TObjArray*)ce->GetCalPadRMS());
652 AliTPCCalPad *padq = new AliTPCCalPad((TObjArray*)ce->GetCalPadQ());
653 padtime->SetName("CE_T");
654 padRMS->SetName("CE_RMS");
655 padq->SetName("CE_Q");
73e9c7cf 656
657 TFile f2(pulname);
658 AliTPCCalibPulser *pul = (AliTPCCalibPulser*)f2.Get("AliTPCCalibPulser");
659 AliTPCCalPad *pultime = new AliTPCCalPad((TObjArray*)pul->GetCalPadT0());
660 AliTPCCalPad *pulRMS = new AliTPCCalPad((TObjArray*)pul->GetCalPadRMS());
661 AliTPCCalPad *pulq = new AliTPCCalPad((TObjArray*)pul->GetCalPadQ());
662 pultime->SetName("P_T");
663 pulRMS->SetName("P_RMS");
664 pulq->SetName("P_Q");
665
666
9e7da9da 667 AliTPCPreprocessorOnline * preprocesor = new AliTPCPreprocessorOnline;
668 preprocesor->AddComponent(padtime);
669 preprocesor->AddComponent(padq);
670 preprocesor->AddComponent(padRMS);
73e9c7cf 671 preprocesor->AddComponent(pultime);
672 preprocesor->AddComponent(pulq);
673 preprocesor->AddComponent(pulRMS);
9e7da9da 674 preprocesor->DumpToFile(inname);
675}
0803cf01 676
677
678void AnalyzeLaserCE(){
679 RebuildCE("ce.root","pul.root");
680 StoreData();
681 StoreTree();
682 RebuildData();
683 StoreData();
684 StoreTree();
685 MakeFit();
686 StoreData();
687 StoreTree();
688 MakeRes();
689}
690
2c55475e 691
692
693void AddFiles(char *list){
35d81915 694 /// prepare viewer for data sets
695
2c55475e 696 fstream finput;
697 finput.open(list, ios_base::in);
698 //
699 TString currentFile;
700 Int_t counter=0;
701 if (!elist ) {
702 elist = new TEntryList("elist","elist");
703 tree->Draw(Form(">>elist%d",counter),"1","entrylist");
704 }
705 while(finput.good()) {
706 finput >> currentFile;
707 printf("Getting file%s\n",currentFile.Data());
708 TFile * fC = new TFile(currentFile.Data());
709 TTree * treeC = (TTree*) fC->Get("calPads");
710 makePad->AddFriend(treeC,Form("T%d",counter));
711 //
712 TString cutStrF1=Form("abs(T%d.timeIn.fElements-T%d.timeF1.fElements)<1.5",counter,counter);
713 TEntryList *celist = new TEntryList(Form("listF1%d",counter),Form("listF1%d",counter));
714 tree->Draw(Form(">>listF1%d",counter),cutStrF1,"entrylist");
715 if (celist->GetN()>0) elist->Add(celist);
716 counter++;
717 }
718 viewer->Reload();
719}
720
c33ec98a 721//void AnalyzeSectors(){
722// //
723// // Analyze sector by sector
724// //
725// TTreeSRedirector cstream("fits.root");
726// for (Int_t
727//}