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