]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TPC/CalibMacros/FitRodShift.C
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
doxy: TPC/CalibMacros
[u/mrichter/AliRoot.git] / TPC / CalibMacros / FitRodShift.C
CommitLineData
35d81915 1/// \file FitRodShift.C
2///
3/// \author marian.ivanov@cern.ch, Stefan.Rossegger@cern.ch
4///
5/// Macro to fit alignment/E field distortion maps
6/// As a input output cluster distortion maps (produced by AliTPCcalibAling class)
7/// are used. Cluster distrotion maps granularity (180 *phi x44 *theta x 53 R)
8///
9/// In total 440 parameters to fit using global fit:
10///
11/// 1. Rotation and translation for each sector (72x2)
12/// 2. Rotation and translation for each quadrant of OROC (36x4x2)
13/// 3. Rod/strip shifts in IFC and OFC (18 sectors x 2 sides x 2 )
14/// 4. Rotated clips in IFC and OFC
15///
16///
17/// Input file mean.root with distortion tree expected to be in directory:
18/// ../mergeField0/clusterDY.root
19/// The ouput file fitRodShift.root contains:
20/// 1. Resulting (residual) AliTPCCalibMisalignment and AliTPCFCVoltError3D classes
21/// 2. All important temporary results are stored in the workspace (TTreeSRedirector) associated
22/// with the file - fitrodShift
23/// 2.a Fit parameters with errors
24/// 2.b QA default graphs
25/// 2.c QA defualt histograms
26///
27///
28///
29/// Functions:
30/// 1. LoadModels() - load models to fit - Rod shift (10,20)
31/// 2. RegisterAlignFunction() - register align functions (10-16)
32/// 3. LoadTrees - load trees and make aliases
33/// 4. PrintFit - helper function to print substring of the fit
34/// 5. DeltaLookup - function to calulate the distortion for given distortion cunction
35/// -
36/// 6. MakeAliases - Make tree aliases shortcuts for the fitting function
37///
38/// 7. MakeAlignCorrection - Crete alignment entry to be stored in the OCDB
39/// - Dump fit parameters and errors
40/// 8. MakeQuadrantCorrection -
41/// - Dump fit parameters and errors
42///
43/// 9. FitRodShifts - main minimization function
44///
45/// ~~~
46/// .x ~/rootlogon.C
47/// .x ~/NimStyle.C
48/// .L $ALICE_ROOT/TPC/CalibMacros/FitRodShift.C+
49/// FitRodShift(kFALSE);
50/// ~~~
32b5322b 51
52#if !defined(__CINT__) || defined(__MAKECINT__)
53#include "TFile.h"
54#include "TTreeStream.h"
55#include "TMath.h"
56#include "TGraph.h"
57#include "TRandom.h"
58#include "TTree.h"
59#include "TF1.h"
60#include "TH2F.h"
61#include "TH3F.h"
62#include "TAxis.h"
63#include "TPad.h"
64#include "TCanvas.h"
65#include "TStyle.h"
66#include "AliTPCParamSR.h"
67#include "TDatabasePDG.h"
68#include "AliTPCFCVoltError3D.h"
69#include "AliTPCROCVoltError3D.h"
70#include "AliTPCComposedCorrection.h"
71#include "AliTPCBoundaryVoltError.h"
72#include "AliCDBEntry.h"
73#include "AliTPCROC.h"
74#include <TStatToolkit.h>
75#include "TCut.h"
76#include "TGraphErrors.h"
77#include "AliTPCCalibGlobalMisalignment.h"
78#endif
79
80TTreeSRedirector *pcWorkspace= 0; // Workspace to store the fit parameters and QA histograms
81TTree * treeDY= 0; //distortion tree - Dy
82TTree * treeDZ= 0; //distortion tree - Dz
83// fit models ...
84AliTPCFCVoltError3D *rotOFC =0; // fit models
85AliTPCFCVoltError3D *rodOFC1 =0;
86AliTPCFCVoltError3D *rodOFC2 =0;
87AliTPCFCVoltError3D *rotIFC =0;
88AliTPCFCVoltError3D *rodIFC1 =0;
89AliTPCFCVoltError3D *rodIFC2 =0;
90AliTPCFCVoltError3D *rodAll =0; //all lookup initialized
91//
92AliTPCComposedCorrection *corFit0 = 0; // composed correction
93void FitFunctionQA();
94void MakeQA();
95void DrawAlignParam();
96
97void PrintFit(TString fitString, TString filter){
35d81915 98 /// helper function to print substring of the fit
99 /// TString fitString =*strFitGA;
100 /// TString filter="rot";
101
32b5322b 102 TObjArray *arr = fitString.Tokenize("++");
103 Int_t entries=arr->GetEntries();
104 for (Int_t i=0; i<entries; i++){
105 TString aaa = arr->At(i)->GetName();
106 if (aaa.Contains(filter)==1){
107 printf("%d\t%s\n",i,aaa.Data());
108 }
109 }
110}
111
112
113void LoadModels(){
35d81915 114 /// Load the models from the file
115 /// Or create it
116
32b5322b 117 Int_t volt = 1;
118 TFile f("OCDB-RodShifts.root");
119 AliCDBEntry *entry = (AliCDBEntry*) f.Get("AliCDBEntry");
120 if (entry) { // load file
121 AliTPCComposedCorrection *cor = (AliTPCComposedCorrection*) entry->GetObject();
122 cor->Print();
123 TCollection *iter = cor->GetCorrections();
124 rotOFC = (AliTPCFCVoltError3D*)iter->FindObject("rotOFC");
125 rodOFC1 = (AliTPCFCVoltError3D*)iter->FindObject("rodOFC1");
126 rodOFC2 = (AliTPCFCVoltError3D*)iter->FindObject("rodOFC2");
127 rotIFC = (AliTPCFCVoltError3D*)iter->FindObject("rotIFC");
128 rodIFC1 = (AliTPCFCVoltError3D*)iter->FindObject("rodIFC1");
129 rodIFC2 = (AliTPCFCVoltError3D*)iter->FindObject("rodIFC2");
130 rodAll = (AliTPCFCVoltError3D*)iter->FindObject("rodAll");
131 // rotOFC->CreateHistoDZinXY(1,500,500)->Draw("surf2");
132 } else {
133 // OFC
134 rotOFC = new AliTPCFCVoltError3D();
135 rotOFC->SetOmegaTauT1T2(0,1,1);
136 rotOFC->SetRotatedClipVoltA(1,volt);
137 rotOFC->SetRotatedClipVoltC(1,volt);
138 //
139 rodOFC1 = new AliTPCFCVoltError3D();
140 rodOFC1->SetOmegaTauT1T2(0,1,1);
141 rodOFC1->SetRodVoltShiftA(18,volt);
142 rodOFC1->SetRodVoltShiftC(18,volt);
143 //
144 rodOFC2 = new AliTPCFCVoltError3D();
145 rodOFC2->SetOmegaTauT1T2(0,1,1);
146 rodOFC2->SetCopperRodShiftA(18,volt);
147 rodOFC2->SetCopperRodShiftC(18,volt);
148 // IFC
149 rotIFC = new AliTPCFCVoltError3D();
150 rotIFC->SetOmegaTauT1T2(0,1,1);
151 rotIFC->SetRotatedClipVoltA(0,volt);
152 rotIFC->SetRotatedClipVoltC(0,volt);
153 //
154 rodIFC1 = new AliTPCFCVoltError3D();
155 rodIFC1->SetOmegaTauT1T2(0,1,1);
156 rodIFC1->SetRodVoltShiftA(0,volt);
157 rodIFC1->SetRodVoltShiftC(0,volt);
158 //
159 rodIFC2 = new AliTPCFCVoltError3D();
160 rodIFC2->SetOmegaTauT1T2(0,1,1);
161 rodIFC2->SetCopperRodShiftA(0,volt);
162 rodIFC2->SetCopperRodShiftC(0,volt);
163 // dummy object with all inits
164 rodAll = new AliTPCFCVoltError3D();
165 Double_t volt0=0.0000000001;
166 rodAll->SetRotatedClipVoltA(1,volt0);
167 rodAll->SetRotatedClipVoltC(1,volt0);
168 rodAll->SetRodVoltShiftA(18,volt0);
169 rodAll->SetRodVoltShiftC(18,volt0);
170 rodAll->SetCopperRodShiftA(18,volt0);
171 rodAll->SetCopperRodShiftC(18,volt0);
172 rodAll->SetRotatedClipVoltA(0,volt0);
173 rodAll->SetRotatedClipVoltC(0,volt0);
174 rodAll->SetRodVoltShiftA(0,volt0);
175 rodAll->SetRodVoltShiftC(0,volt0);
176 rodAll->SetCopperRodShiftA(0,volt0);
177 rodAll->SetCopperRodShiftC(0,volt0);
178 rodAll->SetOmegaTauT1T2(0,1,1);
179 //
180 //
181 // Initialization of the lookup tables
182 //
183 printf(" ------- OFC rotated clip:\n"); rotOFC->InitFCVoltError3D();
184 printf(" ------- OFC rod & strip:\n"); rodOFC1->InitFCVoltError3D();
185 printf(" ------- OFC copper rod:\n"); rodOFC2->InitFCVoltError3D();
186 printf(" ------- IFC rotated clip:\n"); rotIFC->InitFCVoltError3D();
187 printf(" ------- IFC rod & strip:\n"); rodIFC1->InitFCVoltError3D();
188 printf(" ------- IFC copper rod:\n"); rodIFC2->InitFCVoltError3D();
189 printf(" ------- Dummy all:\n"); rodAll->InitFCVoltError3D();
190 // give names
191 rotOFC->SetName("rotOFC");rotOFC->SetTitle("rotOFC");
192 rodOFC1->SetName("rodOFC1");rodOFC1->SetTitle("rodOFC1");
193 rodOFC2->SetName("rodOFC2");rodOFC2->SetTitle("rodOFC2");
194 rotIFC->SetName("rotIFC");rotIFC->SetTitle("rotIFC");
195 rodIFC1->SetName("rodIFC1");rodIFC1->SetTitle("rodIFC1");
196 rodIFC2->SetName("rodIFC2");rodIFC2->SetTitle("rodIFC2");
197 rodAll->SetName("rodAll");rodAll->SetTitle("rodAll");
198 //
199 // save in file
200 AliTPCComposedCorrection *corFit0 = new AliTPCComposedCorrection();
201 TObjArray *cs = new TObjArray();
202 cs->Add(rotIFC); cs->Add(rotOFC);
203 cs->Add(rodIFC1); cs->Add(rodOFC1);
204 cs->Add(rodIFC2); cs->Add(rodOFC2);
205 cs->Add(rodAll);
206 corFit0->SetCorrections(cs);
207 corFit0->SetOmegaTauT1T2(0,1.,1.);
208 corFit0->Print();
209 corFit0->StoreInOCDB(0,1000000,"testForFit");
210 }
211 //
212 AliTPCCorrection::AddVisualCorrection(rotOFC,0);
213 AliTPCCorrection::AddVisualCorrection(rodOFC1,1);
214 AliTPCCorrection::AddVisualCorrection(rodOFC2,2);
215 AliTPCCorrection::AddVisualCorrection(rotIFC,3);
216 AliTPCCorrection::AddVisualCorrection(rodIFC1,4);
217 AliTPCCorrection::AddVisualCorrection(rodIFC2,5);
218 AliTPCCorrection::AddVisualCorrection(rodAll,6);
219}
220
221
222void RegisterAlignFunction(){
35d81915 223 /// Register primitive alignment components.
224 /// Linear conbination of primitev forulas used for fit
225 /// The nominal delta 1 mm in shift and 1 mrad in rotation
226 /// Primitive formulas registeren in AliTPCCoreection::AddvisualCorrection
227 /// 10 - deltaX
228 /// 11 - deltaY
229 /// 12 - deltaZ
230 /// 13 - rot0 (phi)
231 /// 14 - rot1 (theta)
232 /// 15 - rot2
233
32b5322b 234 TGeoHMatrix matrixX;
235 TGeoHMatrix matrixY;
236 TGeoHMatrix matrixZ;
237 TGeoRotation rot0;
238 TGeoRotation rot1;
239 TGeoRotation rot2; //transformation matrices
240 TGeoRotation rot90; //transformation matrices
241 matrixX.SetDx(0.1); matrixY.SetDy(0.1); matrixZ.SetDz(0.1); //1 mm translation
242 rot0.SetAngles(0.001*TMath::RadToDeg(),0,0);
243 rot1.SetAngles(0,0.001*TMath::RadToDeg(),0);
244 rot2.SetAngles(0,0,0.001*TMath::RadToDeg());
245 //how to get rot02 ?
246 rot90.SetAngles(0,90,0);
247 rot2.MultiplyBy(&rot90,kTRUE);
248 rot90.SetAngles(0,-90,0);
249 rot2.MultiplyBy(&rot90,kFALSE);
250 AliTPCCalibGlobalMisalignment *alignRot0 =new AliTPCCalibGlobalMisalignment;
251 alignRot0->SetAlignGlobal(&rot0);
252 AliTPCCalibGlobalMisalignment *alignRot1=new AliTPCCalibGlobalMisalignment;
253 alignRot1->SetAlignGlobal(&rot1);
254 AliTPCCalibGlobalMisalignment *alignRot2=new AliTPCCalibGlobalMisalignment;
255 alignRot2->SetAlignGlobal(&rot2);
256 //
257 AliTPCCalibGlobalMisalignment *alignTrans0 =new AliTPCCalibGlobalMisalignment;
258 alignTrans0->SetAlignGlobal(&matrixX);
259 AliTPCCalibGlobalMisalignment *alignTrans1=new AliTPCCalibGlobalMisalignment;
260 alignTrans1->SetAlignGlobal(&matrixY);
261 AliTPCCalibGlobalMisalignment *alignTrans2=new AliTPCCalibGlobalMisalignment;
262 alignTrans2->SetAlignGlobal(&matrixZ);
263 AliTPCCorrection::AddVisualCorrection(alignTrans0 ,10);
264 AliTPCCorrection::AddVisualCorrection(alignTrans1 ,11);
265 AliTPCCorrection::AddVisualCorrection(alignTrans2 ,12);
266 AliTPCCorrection::AddVisualCorrection(alignRot0 ,13);
267 AliTPCCorrection::AddVisualCorrection(alignRot1 ,14);
268 AliTPCCorrection::AddVisualCorrection(alignRot2 ,15);
269 TObjArray arrAlign(6);
270 arrAlign.AddAt(alignTrans0->Clone(),0);
271 arrAlign.AddAt(alignTrans1->Clone(),1);
272 arrAlign.AddAt(alignTrans2->Clone(),2);
273 arrAlign.AddAt(alignRot0->Clone(),3);
274 arrAlign.AddAt(alignRot1->Clone(),4);
275 arrAlign.AddAt(alignRot2->Clone(),5);
276 //combAlign.SetCorrections((TObjArray*)arrAlign.Clone());
277}
278
279
280void LoadTrees(){
35d81915 281 /// 1. Load trees
282 /// 2. Make standard cuts - filter the tree
283 /// 3. makeStandard aliases
284
32b5322b 285 TFile *fy = new TFile("clusterDY.root");
286 TFile *fz = new TFile("clusterDZ.root");
287 treeDY= (TTree*)fy->Get("delta");
288 treeDZ= (TTree*)fz->Get("delta");
289 TCut cutAll = "entries>3000&&abs(kZ)<1&&localX>80&&localX<246&&abs(sector-int(sector)-0.5)<0.41&&abs(localX-134)>2&&rmsG<0.5&&abs(localX-189)>3";
290
291 treeDY->Draw(">>outListY",cutAll,"entryList");
292 TEntryList *elistOutY = (TEntryList*)gDirectory->Get("outListY");
293 treeDY->SetEntryList(elistOutY);
294 treeDZ->Draw(">>outListZ",cutAll,"entryList");
295 TEntryList *elistOutZ = (TEntryList*)gDirectory->Get("outListZ");
296 treeDZ->SetEntryList(elistOutZ);
297 //
298 // Make aliases
299 //
300 treeDY->SetAlias("dsec","(sector-int(sector)-0.5)");
301 treeDY->SetAlias("dsec0","(sector-int(sector))");
302 treeDY->SetAlias("signy","(-1.+2*(sector-int(sector)>0.5))");
303 treeDY->SetAlias("dx","(localX-165.5)"); // distance X to ref. plane
304 treeDY->SetAlias("dxm","0.001*(localX-165.5)");
305 treeDY->SetAlias("rx","(localX/166.5)"); //ratio distrnace to reference plane
306 treeDY->SetAlias("dq1","(((q1==0)*(-rx)+q1*(1-rx))*signy)");
307 //
308 {for (Int_t isec=0; isec<18; isec++){ // sectors
309 treeDY->SetAlias(Form("sec%d",isec), Form("(abs(sector-%3.1lf)<0.5)",isec+0.5));
310 treeDZ->SetAlias(Form("sec%d",isec), Form("(abs(sector-%3.1lf)<0.5)",isec+0.5));
311 }}
312 treeDY->SetAlias("gy","localX*sin(pi*sector/9)");
313 treeDY->SetAlias("gx","localX*cos(pi*sector/9)");
314 treeDY->SetAlias("side","(-1.+2.*(kZ>0))");
315 treeDY->SetAlias("drphi","mean");
316 treeDY->SetMarkerStyle(25);
317}
318
319Double_t DeltaLookup(Double_t sector, Double_t localX, Double_t kZ, Double_t xref, Int_t value, Int_t corr){
35d81915 320 /// Distortion maps are calculated relative to the reference X plane
321 /// The same procedure applied for given correction
322 /// fd(sector, localX, kZ) ==> fd(sector, localX, kZ)-fd(sector, xref, kZ)*localX/xref
323
32b5322b 324 Double_t distortion = AliTPCCorrection::GetCorrSector(sector,localX,kZ,value,corr);
325 Double_t distortionRef = AliTPCCorrection::GetCorrSector(sector,xref,kZ,value,corr)*localX/xref;
326 return distortion-distortionRef;
327}
328
329void MakeAliases(){
35d81915 330 /// make alias names
331
32b5322b 332 AliTPCROC * roc = AliTPCROC::Instance();
333 Double_t xref = ( roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5;
334 treeDY->SetAlias("iroc","(localX<134)"); // IROC
335 treeDY->SetAlias("oroc","(localX>134)"); // OROC
336 treeDY->SetAlias("q1","abs(localX-161)<28"); // OROC first haf
337 treeDY->SetAlias("q2","(localX>189)"); // OROC second half
338 //
339 treeDY->SetAlias("dIFC","abs(localX-83)"); // distance to IFC
340 treeDY->SetAlias("dOFC","abs(localX-250)"); // distance to OFC
341 treeDY->SetAlias("errY","(1.+1/(1+(dIFC/2.))+1/(1+dOFC/2.))");
342
343 //
344 // rotated clip - IFC --------------------
345 treeDY->SetAlias("rotClipI","DeltaLookup(sector,localX,kZ,165.6,1,3+0)");
346 // rotated clip - OFC --------------------
347 treeDY->SetAlias("rotClipO","DeltaLookup(sector,localX,kZ,165.6,1,0+0)");
348 for (Int_t isec=0;isec<18; isec++){
349 //
350 // alignment IROC - shift cm - rotation mrad
351 treeDY->SetAlias(Form("dyI_%d",isec),Form("(iroc*sec%d)",isec));
352 treeDY->SetAlias(Form("drotI_%d",isec),Form("(iroc*sec%d*(localX-%3.1lf)/1000.)",isec,xref));
353 // alignment OROC
354 treeDY->SetAlias(Form("dyO_%d",isec),Form("(sec%d*oroc-sec%d*localX/165.6)",isec,isec)); // fix local Y shift - do not use it
355 treeDY->SetAlias(Form("drotO_%d",isec),Form("(sec%d*oroc*(localX-%3.1lf)/1000.)",isec,xref));
356 //
357 // quadrant shift OROC
358 treeDY->SetAlias(Form("dqO0_%d",isec),Form("(sec%d*sign(dsec)*(q1-localX/165.6))",isec));
359 // quadrant delta rotation OFC inner
360 treeDY->SetAlias(Form("dqOR0_%d",isec),Form("(sec%d*sign(dsec)*q1*(localX-165.6)/1000.)",isec));
361 //
362 treeDY->SetAlias(Form("dqO1_%d",isec),Form("(q2*sec%d*sign(dsec))",isec)); // delta
363 treeDY->SetAlias(Form("dqOR1_%d",isec),Form("(q2*sec%d*sign(dsec)*(localX-165.6)/1000.)",isec));
364 treeDY->SetAlias(Form("dqO2_%d",isec),Form("(q2*sec%d)",isec)); //common
365 treeDY->SetAlias(Form("dqOR2_%d",isec),Form("(q2*sec%d*(localX-165.6)/1000.)",isec));
366 //
367 //
368 // shifted rod - OFC
369 treeDY->SetAlias(Form("rodStripO_%d",isec),Form("DeltaLookup(sector-%d,localX,kZ,165.6,1,1+0)",isec));
370 // Shifted cooper rod OFC
371 treeDY->SetAlias(Form("coppRodO_%d",isec),Form("DeltaLookup(sector-%d,localX,kZ,165.6,1,2+0)",isec));
372 // shifted rod - IFC
373 treeDY->SetAlias(Form("rodStripI_%d",isec),Form("DeltaLookup(sector-%d,localX,kZ,165.6,1,4+0)",isec));
374 // Shifted cooper rod IFC
375 treeDY->SetAlias(Form("coppRodI_%d",isec),Form("DeltaLookup(sector-%d,localX,kZ,165.6,1,5+0)",isec));
376 }
377 //
378 // fitted correction - Using the saved coposed corrections
379 //
380 treeDY->SetAlias("dAlign","DeltaLookup(sector,localX,kZ,165.6,1,1000+0)");
381 treeDY->SetAlias("dQuadrant","DeltaLookup(sector,localX,kZ,165.6,1,1100+0)");
382 treeDY->SetAlias("dField","DeltaLookup(sector,localX,kZ,165.6,1,100+0)");
383 treeDY->SetAlias("dAll","(dAlign+dQuadrant+dField)");
384}
385
386AliTPCCalibGlobalMisalignment *MakeAlignCorrection(TVectorD paramA, TVectorD paramC, TMatrixD covar, Double_t chi2){
35d81915 387 /// Make a global alignmnet
388 /// Take a fit parameters and make a combined correction:
389 /// Only delta local Y and delta phi are fitted - not sensitivity for other parameters
390 /// GX and GY shift extracted per side.
391 ///
392 /// Algorithm:
393 /// 1. Loop over sectors
394 /// 2. Make combined AliTPCCalibGlobalMisalignment
395
32b5322b 396 AliTPCCalibGlobalMisalignment *alignLocal =new AliTPCCalibGlobalMisalignment;
397 Int_t offset=3;
398 AliTPCROC * roc = AliTPCROC::Instance();
399 Double_t xref = ( roc->GetPadRowRadii(0,0)+roc->GetPadRowRadii(36,roc->GetNRows(36)-1))*0.5;
400 //
401 pcWorkspace->GetFile()->cd();
402 TObjArray * array = new TObjArray(72);
403
404 for (Int_t side=-1; side<2; side+=2){
405 TVectorD &param = (side==1) ? paramA : paramC;
406 TGeoHMatrix matrixGX;
407 TGeoHMatrix matrixGY;
408 matrixGX.SetDx(0.1*param[1]);
409 matrixGY.SetDy(0.1*param[2]);
410 //
411 for (Int_t isec=0; isec<18; isec++){
412 TGeoHMatrix matrixOROC;
413 TGeoHMatrix matrixIROC;
414 TGeoRotation rotIROC;
415 TGeoRotation rotOROC;
416 Double_t phi= (Double_t(isec)+0.5)*TMath::Pi()/9.;
417 //
418 Double_t drotIROC = -param[offset+isec+18]*0.001;
419 Double_t drotOROC = -param[offset+isec+36]*0.001;
420 Double_t dlyIROC = param[offset+isec];
421 Double_t dlyIROC0 = param[offset+isec]+drotIROC*xref; // shift at x=0
422 Double_t dlyOROC = 0;
423 Double_t dlyOROC0 = 0+drotOROC*xref; // shift at x=0
424 //
425 Double_t dgxIROC = TMath::Cos(phi)*0+TMath::Sin(phi)*dlyIROC0;
426 Double_t dgyIROC = TMath::Sin(phi)*0-TMath::Cos(phi)*dlyIROC0;
427 Double_t dgxOROC = TMath::Cos(phi)*0+TMath::Sin(phi)*dlyOROC0;
428 Double_t dgyOROC = TMath::Sin(phi)*0-TMath::Cos(phi)*dlyOROC0;
429 Double_t errYIROC=TMath::Sqrt(covar(offset+isec, offset+isec)*chi2);
430 Double_t errPhiIROC=TMath::Sqrt(covar(offset+isec+18, offset+isec+18)*chi2);
431 Double_t errPhiOROC=TMath::Sqrt(covar(offset+isec+36, offset+isec+36)*chi2);
432 matrixIROC.SetDx(dgxIROC);
433 matrixIROC.SetDy(dgyIROC);
434 matrixOROC.SetDx(dgxOROC);
435 matrixOROC.SetDy(dgyOROC);
436 rotIROC.SetAngles(drotIROC*TMath::RadToDeg(),0,0);
437 matrixIROC.Multiply(&matrixGX);
438 matrixIROC.Multiply(&matrixGY);
439 matrixIROC.Multiply(&rotIROC);
440 rotOROC.SetAngles(drotOROC*TMath::RadToDeg(),0,0);
441 matrixOROC.Multiply(&matrixGX);
442 matrixOROC.Multiply(&matrixGY);
443 matrixOROC.Multiply(&rotOROC);
444 if (side>0){
445 array->AddAt(matrixIROC.Clone(),isec);
446 array->AddAt(matrixOROC.Clone(),isec+36);
447 }
448 if (side<0){
449 array->AddAt(matrixIROC.Clone(),isec+18);
450 array->AddAt(matrixOROC.Clone(),isec+36+18);
451 }
452 Double_t rms=TMath::Sqrt(chi2); //chi2 normalized 1/npoints before
453 (*pcWorkspace)<<"align"<<
454 "isec="<<isec<< // sector
455 "side="<<side<< // side
456 "phi="<<phi<< // phi
457 "rms="<<rms<< // rms in cm
458 "cov.="<<&covar<<
459 // errors
460 "errYIROC="<<errYIROC<< //error in local y
461 "errPhiIROC="<<errPhiIROC<< //error in phi IROC
462 "errPhiOROC="<<errPhiOROC<< //error in phi OROC
463 //
464 "dlyIROC="<<dlyIROC<< //delta Local Y at refX=165 -IROC
465 "dlyIROC0="<<dlyIROC0<< //delta Local Y at refX=0 -IROC
466 "dgxIROC="<<dgxIROC<< //
467 "dgyIROC="<<dgyIROC<<
468 "drotIROC="<<drotIROC<<
469 //
470 "dlyOROC="<<dlyOROC<< // assuming 0 at
471 "dlyOROC0="<<dlyOROC0<<
472 "dgxOROC="<<dgxOROC<<
473 "dgyOROC="<<dgyOROC<<
474 "drotOROC="<<drotOROC<<
475 "\n";
476 }
477 }
478 alignLocal->SetAlignSectors(array);
479 AliTPCCorrection::AddVisualCorrection(alignLocal,1000);
480 /*
481 Check:
482 alignLocal->CreateHistoDRPhiinXY(10,100,100)->Draw("colz"); // OK
483 treeDY->Draw("DeltaLookup(sector,localX,kZ,165.6,1,1001):localX:int(sector)","kZ>0","colz");
484 //
485 (( (*pcWorkspace)<<"align"))->GetTree()->Draw("dlyOROC-dlyOROC:isec","","*")
486 treeDY->Draw("tfitA:DeltaLookup(sector,localX,kZ,165.6,1,1001):sector","kZ>0","colz");
487
488 */
489 return alignLocal;
490}
491
492AliTPCCalibGlobalMisalignment *MakeQuadrantCorrection(TVectorD paramA, TVectorD paramC, TMatrixD covar, Double_t chi2){
35d81915 493 /// Make a global alignmnet
494 /// side= 1 - A side
495 /// side=-1 - C side
496 /// Take a fit parameters and make a combined correction:
497 /// Only delta local Y and delta phi are fitted - not sensitivity for other parameters
498 /// GX and GY shift extracted per side.
499 ///
500 /// Algorithm:
501 /// 1. Loop over sectors
502 /// 2. Make combined AliTPCCalibGlobalMisalignment
503
32b5322b 504 AliTPCCalibGlobalMisalignment *alignLocalQuadrant =new AliTPCCalibGlobalMisalignment;
505 //
506 Int_t offset=3+3*18;
507 TVectorD quadrantQ0(36); //dly+=sign*(*fQuadrantQ0)[isec]; // shift in cm
508 TVectorD quadrantRQ0(36); //dly+=sign*(*fQuadrantRQ0)[isec]*(pos[0]-xref);
509 TVectorD quadrantQ1(36); //dly+=sign*(*fQuadrantQ1)[isec]; // shift in cm
510 TVectorD quadrantRQ1(36); //dly+=sign*(*fQuadrantRQ1)[isec]*(pos[0]-xref);
511 TVectorD quadrantQ2(36); //dly+=(*fQuadrantQ2)[isec]; // shift in cm
512 TVectorD quadrantRQ2(36); //dly+=(*fQuadrantRQ2)[isec]*(pos[0]-xref);
513
514 for (Int_t side=-1; side<2; side+=2){
515 Int_t offsetSec= (side==1) ? 0:18;
516 TVectorD &param = (side==1) ? paramA : paramC;
517 for (Int_t isec=0; isec<18; isec++){
518 Double_t q0=-param[offset+isec+0*18];
519 Double_t rq0=-param[offset+isec+1*18]*0.001;
520 Double_t q1=-param[offset+isec+2*18];
521 Double_t rq1=-param[offset+isec+3*18]*0.001;
522 Double_t q2=-param[offset+isec+4*18];
523 Double_t rq2=-param[offset+isec+5*18]*0.001;
524 //
525 Double_t sq0=TMath::Sqrt(covar(offset+isec+0*18,offset+isec+0*18)*chi2);
526 Double_t srq0=TMath::Sqrt(covar(offset+isec+1*18,offset+isec+1*18)*chi2)*0.001;
527 Double_t sq1=TMath::Sqrt(covar(offset+isec+2*18,offset+isec+2*18)*chi2);
528 Double_t srq1=TMath::Sqrt(covar(offset+isec+3*18,offset+isec+3*18)*chi2)*0.001;
529 Double_t sq2=TMath::Sqrt(covar(offset+isec+4*18,offset+isec+4*18)*chi2);
530 Double_t srq2=TMath::Sqrt(covar(offset+isec+5*18,offset+isec+5*18)*chi2)*0.001;
531 //
532 quadrantQ0[offsetSec+isec]=q0;
533 quadrantRQ0[offsetSec+isec]=rq0;
534 quadrantQ1[offsetSec+isec]=q1;
535 quadrantRQ1[offsetSec+isec]=rq1;
536 quadrantQ2[offsetSec+isec]=q2;
537 quadrantRQ2[offsetSec+isec]=rq2;
538 Double_t rms=TMath::Sqrt(chi2); //chi2 normalized 1/npoints before
539 (*pcWorkspace)<<"quadrant"<<
540 "isec="<<isec<< // sector
541 "side="<<side<< // side
542 "rms="<<rms<< // rms in cm
543 "cov.="<<&covar<<
544 //
545 "q0="<<q0<< // quadrant alignment
546 "rq0="<<rq0<<
547 "q1="<<q1<<
548 "rq1="<<rq1<<
549 "q2="<<q2<<
550 "rq2="<<rq2<<
551 "sq0="<<sq0<< //rms of quadrant parameters
552 "srq0="<<srq0<<
553 "sq1="<<sq1<<
554 "srq1="<<srq1<<
555 "sq2="<<sq2<<
556 "srq2="<<srq2<<
557 "\n";
558 }
559 }
560 alignLocalQuadrant->SetQuadranAlign(&quadrantQ0, &quadrantRQ0, &quadrantQ1, &quadrantRQ1, &quadrantQ2, &quadrantRQ2);
561 AliTPCCorrection::AddVisualCorrection(alignLocalQuadrant,1100);
562 /*
563 (( (*pcWorkspace)<<"quadrant"))->GetTree()->Draw("q0:isec","","*")
564
565 alignLocalQuadrant->CreateHistoDRPhiinXY(10,100,100)->Draw("colz"); //OK
566
567 treeDY->Draw("tfitA-DeltaLookup(sector,localX,kZ,165.6,1,1000):DeltaLookup(sector,localX,kZ,165.6,1,1100):int(sector)","kZ>0&&oroc","colz");
568 */
569 return alignLocalQuadrant;
570}
571
572AliTPCFCVoltError3D* MakeEfieldCorrection(TVectorD paramA, TVectorD paramC, TMatrixD covar, Double_t chi2){
35d81915 573 /// Make a global AliTPCFCVoltError3D object
574 ///
575 /// Take a fit parameters and make a combined correction:
576 /// Only delta local Y and delta phi are fitted - not sensitivity for other parameters
577 /// GX and GY shift extracted per side.
578 ///
579 /// Algorithm:
580 /// 1. Loop over sectors
581 /// 2. Make combined AliTPCCalibGlobalMisalignment
582
32b5322b 583 Int_t offset=3+9*18;
584 //
585 AliTPCFCVoltError3D* corrField = new AliTPCFCVoltError3D;
586 //
587 Double_t rotIROCA=paramA[offset+0];
588 Double_t rotIROCC=paramC[offset+0];
589 Double_t rotOROCA=paramA[offset+1];
590 Double_t rotOROCC=paramC[offset+1];
591 //
592 corrField->SetRotatedClipVoltA(0,paramA[offset+0]);
593 corrField->SetRotatedClipVoltC(0,paramC[offset+0]);
594 corrField->SetRotatedClipVoltA(1,paramA[offset+1]);
595 corrField->SetRotatedClipVoltC(1,paramC[offset+1]);
596 {for (Int_t isec=0; isec<18; isec++){
597 corrField->SetRodVoltShiftA(isec,paramA[offset+2+36+isec]); // rod shift IFC
598 corrField->SetRodVoltShiftA(18+isec,paramA[offset+2+isec]); // rod shift OFC
599 corrField->SetRodVoltShiftC(isec,paramC[offset+2+36+isec]); // rod shift IFC
600 corrField->SetRodVoltShiftC(18+isec,paramC[offset+2+isec]); // rod shift OFC
601 Double_t rodIROCA=paramA[offset+2+36+isec];
602 Double_t rodIROCC=paramC[offset+2+36+isec];
603 Double_t rodOROCA=paramA[offset+2+isec];
604 Double_t rodOROCC=paramC[offset+2+isec];
605 //
606 Double_t srodIROCA=TMath::Sqrt(covar(offset+2+36+isec,offset+2+36+isec)*chi2);
607 Double_t srodOROCA=TMath::Sqrt(covar(offset+2+isec,offset+2+isec)*chi2);
608 Double_t phi= (Double_t(isec)+0.5)*TMath::Pi()/9.;
609 Double_t rms=TMath::Sqrt(chi2); //chi2 normalized 1/npoints before
610 (*pcWorkspace)<<"field"<<
611 "isec="<<isec<< // sector
612 //"side="<<side<< // side
613 "phi="<<phi<< // phi
614 "rms="<<rms<< // rms in cm
615 "cov.="<<&covar<<
616 //
617 "rotIROCA="<<rotIROCA<<
618 "rotIROCC="<<rotIROCC<<
619 "rotOROCA="<<rotOROCA<<
620 "rotOROCC="<<rotOROCC<<
621 //
622 "rodIROCA="<<rodIROCA<<
623 "rodIROCC="<<rodIROCC<<
624 "rodOROCA="<<rodOROCA<<
625 "rodOROCC="<<rodOROCC<<
626 //
627 "srodIROCA="<<srodIROCA<<
628 "srodOROCA="<<srodOROCA<<
629 "srodIROCC="<<srodIROCA<<
630 "srodOROCC="<<srodOROCA<<
631 "\n";
632 }
633 }
634 corrField->SetOmegaTauT1T2(0,1.,1.);
635 corrField->InitFCVoltError3D();
636 AliTPCCorrection::AddVisualCorrection(corrField,100);
637 return corrField;
638}
639
640
641void FitRodShift(Bool_t flagIFCcopper = kTRUE) {
35d81915 642 /// Main fit function
643 /// In total 440 parameters to fit using global fit:
644 /// 1. Rotation and translation for each sector (72x2)
645 /// 2. Rotation and translation for each quadrant of OROC (36x4x2)
646 /// 3. Rod/strip shifts in IFC and OFC (18 sectors x 2 sides x 2 )
647 /// 4. Rotated clips in IFC and OFC
648
32b5322b 649 LoadTrees();
650 LoadModels();
651 RegisterAlignFunction();
652 MakeAliases();
653
654 pcWorkspace= new TTreeSRedirector("fitAlignLookup.root");
655 TFormula::SetMaxima(10000);
656 //
657 // 1. fit Global
658 //
659 Double_t chi2G=0; Int_t npointsG=0; TVectorD paramG; TMatrixD covarG;
660 TString fstringGlobal="";
661 fstringGlobal+="DeltaLookup(sector,localX,kZ,165.6,1,10-0)++"; // deltaGX
662 fstringGlobal+="DeltaLookup(sector,localX,kZ,165.6,1,11-0)++"; // deltaGY
663 fstringGlobal+="DeltaLookup(sector,localX,kZ,165.6,1,10-0)*side++"; // deltaGX - sides
664 fstringGlobal+="DeltaLookup(sector,localX,kZ,165.6,1,11-0)*side++"; // deltaGY -sides
665 //
666 TString *strFitGlobal = TStatToolkit::FitPlane(treeDY,"drphi", fstringGlobal.Data(),"1", chi2G,npointsG,paramG,covarG,-1,0, 10000000, kTRUE);
667 treeDY->SetAlias("fitYGlobal",strFitGlobal->Data());
668 strFitGlobal->Tokenize("++")->Print();
669 printf("chi2=%f\n",TMath::Sqrt(chi2G/npointsG));
670 // testplots - if necessary - e.g.
671 // treeDY->Draw("AliTPCCorrection::GetCorrSector(sector,localX,kZ,1,2):sector:localX","Cut","colz")
672
673 // FIT PREPARATION +++++++++++++++++++++++++++++++++
674 // define the fit function
675
676 TString fstring="";
677
678 //
679 // Alignment
680 //
681 {
682 fstring+="DeltaLookup(sector,localX,kZ,165.6,1,10-0)++"; // deltaGX
683 fstring+="DeltaLookup(sector,localX,kZ,165.6,1,11-0)++"; // deltaGY
684 // alignment IROC
685 for (Int_t i=0;i<18;i++){
686 fstring+=Form("dyI_%d++",i); // alignment - linear shift (in cm)
687 }
688 for (Int_t i=0;i<18;i++){
689 fstring+=Form("drotI_%d++",i); // alignment - rotation (in mrad)
690 }
691 // alignment OROC
692 for (Int_t i=0;i<18;i++){ // shift of OROC is not allowed
693 //fstring+=Form("dyO_%d++",i); // alignment - linear shift (in cm)
694 }
695 for (Int_t i=0;i<18;i++){
696 fstring+=Form("drotO_%d++",i); // alignment - rotation (in mrad)
697 }
698 //
699 for (Int_t i=0;i<18;i++){
700 fstring+=Form("dqO0_%d++",i); // alignment - quadrant shift OROC in
701 }
702 for (Int_t i=0;i<18;i++){
703 fstring+=Form("dqOR0_%d++",i); // alignment - quadrant rotation OROC in
704 }
705 for (Int_t i=0;i<18;i++){
706 fstring+=Form("dqO1_%d++",i); // alignment - quadrant shift OROC out
707 }
708 for (Int_t i=0;i<18;i++){
709 fstring+=Form("dqOR1_%d++",i); // alignment - quadrant rotation OROC out
710 }
711 for (Int_t i=0;i<18;i++){
712 fstring+=Form("dqO2_%d++",i); // alignment - quadrant shift OROC out
713 }
714 for (Int_t i=0;i<18;i++){
715 fstring+=Form("dqOR2_%d++",i); // alignment - quadrant rotation OROC out
716 }
717 }
718 //
719 // Field distortion
720 //
721 {
722 // rotated clip - IFC --------------------
723 treeDY->SetAlias("rotClipI","DeltaLookup(sector,localX,kZ,165.6,1,3+0)");
724 fstring+="rotClipI++";
725 // rotated clip - OFC --------------------
726 treeDY->SetAlias("rotClipO","DeltaLookup(sector,localX,kZ,165.6,1,0+0)");
727 fstring+="rotClipO++";
728
729 // shifted rod - OFC
730 for (Int_t i=0;i<18;i++){
731 fstring+=Form("rodStripO_%d++",i);
732 }
733 // Shifted cooper rod OFC
734 for (Int_t i=0;i<18;i++){
735 fstring+=Form("coppRodO_%d++",i);
736 }
737 // shifted rod - IFC
738 for (Int_t i=0;i<18;i++){
739 fstring+=Form("rodStripI_%d++",i);
740 }
741 // Shifted cooper rod IFC
742 if (flagIFCcopper) for (Int_t i=0;i<18;i++){
743 fstring+=Form("coppRodI_%d++",i);
744 }
745 }
746 //fstring+=fstringGlobal;
747
748 // FIT ++++++++++++++++++++++++++++++++++++++++++++++++++++++++
749
750 Double_t chi2A=0; Int_t npointsA=0; TVectorD paramA; TMatrixD covarA;
751 Double_t chi2C=0; Int_t npointsC=0; TVectorD paramC; TMatrixD covarC;
752
753
754 printf("Fitting A side\n");
755 TCut cutAllA = "Cut&&kZ>0";
756 TString *strFitGA = TStatToolkit::FitPlane(treeDY,"drphi:errY", fstring.Data(),"kZ>0", chi2A,npointsA,paramA,covarA,-1,0, 10000000, kTRUE);
757 treeDY->SetAlias("tfitA",strFitGA->Data());
758 // strFitGA->Tokenize("++")->Print();
759 printf("chi2=%f\n",TMath::Sqrt(chi2A/npointsA));
760 printf("Sigma Y:\t%f (mm)\n",10.*TMath::Sqrt(covarA(3,3)*chi2A/npointsA));
761 printf("IROC Sigma Angle:\t%f (mrad)\n",TMath::Sqrt(covarA(3+18,3+18)*chi2A/npointsA));
762 printf("OROC Sigma Angle:\t%f (mrad)\n",TMath::Sqrt(covarA(3+36,3+36)*chi2A/npointsA));
763
764 printf("Fitting C side\n");
765 TCut cutAllC = "Cut&&kZ<0";
766 TString *strFitGC = TStatToolkit::FitPlane(treeDY,"drphi:errY", fstring.Data(),"kZ<0", chi2C,npointsC,paramC,covarC,-1,0, 10000000, kTRUE);
767 treeDY->SetAlias("tfitC",strFitGC->Data());
768 // strFitGC->Tokenize("++")->Print();
769 printf("chi2=%f\n",TMath::Sqrt(chi2C/npointsC));
770 printf("Sigma Y:\t%f (mm)\n",10.*TMath::Sqrt(covarC(3,3)*chi2C/npointsC));
771 printf("IROC Sigma Angle:\t%f (mrad)\n",TMath::Sqrt(covarC(3+18,3+18)*chi2C/npointsC));
772 printf("OROC Sigma Angle:\t%f (mrad)\n",TMath::Sqrt(covarC(3+36,3+36)*chi2C/npointsC));
773
774 //
775 // make combined correction
776 //
777 TH2 *hdist =0;
778 //
779 AliTPCCalibGlobalMisalignment * alignLocal = MakeAlignCorrection(paramA, paramC, covarA, chi2A/npointsA);
780 alignLocal->SetName("alignLocal");
781 alignLocal->Write("alignLocal");
782 hdist = alignLocal->CreateHistoDRPhiinXY(10,250,250);
783 hdist->SetName("AlignLocalAside"); hdist->SetTitle("Alignment map (A side)");
784 hdist->Write("AlignLocalAside");
785 hdist = alignLocal->CreateHistoDRPhiinXY(-10,250,250);
786 hdist->SetName("AlignLocalCside"); hdist->SetTitle("Alignment map (C side)");
787 hdist->Write("AlignLocalCside");
788 //
789 AliTPCCalibGlobalMisalignment *alignQuadrant = MakeQuadrantCorrection(paramA, paramC, covarA, chi2A/npointsA);
790 alignQuadrant->SetName("alignQuadrant");
791 alignQuadrant->Write("alignQuadrant");
792 hdist = alignQuadrant->CreateHistoDRPhiinXY(10,250,250);
793 hdist->SetName("AlignQuadrantAside"); hdist->SetTitle("Quadrant Alignment map (A side)");
794 hdist->Write("AlignQuadrantAside");
795 hdist = alignQuadrant->CreateHistoDRPhiinXY(-10,250,250);
796 hdist->SetName("AlignQuadrantCside"); hdist->SetTitle("Quadrant Alignment map (C side)");
797 hdist->Write("AlignQuadrantCside");
798 //
799 AliTPCFCVoltError3D* corrField = MakeEfieldCorrection(paramA, paramC, covarA, chi2A/npointsA);
800 corrField->SetName("corrField");
801 corrField->Write("corrField");
802
803 hdist = corrField->CreateHistoDRPhiinXY(10,250,250);
804 hdist->SetName("AlignEfieldAside"); hdist->SetTitle("Efield Alignment map (A side)");
805 hdist->Write("AlignEfieldAside");
806 hdist = corrField->CreateHistoDRPhiinXY(-10,250,250);
807 hdist->SetName("AlignEfieldCside"); hdist->SetTitle("Efield Alignment map (C side)");
808 hdist->Write("AlignEfieldCside");
809
810
811 //
812 delete pcWorkspace;
813 MakeQA();
814 return;
815}
816
817
818void MakeQA(){
819 LoadTrees();
820 LoadModels();
821 RegisterAlignFunction();
822 MakeAliases();
823 TFile f("fitAlignLookup.root","update");
824 AliTPCCalibGlobalMisalignment* alignLocal = (AliTPCCalibGlobalMisalignment*)f.Get("alignLocal");
825 AliTPCCalibGlobalMisalignment* alignQuadrant = (AliTPCCalibGlobalMisalignment*)f.Get("alignQuadrant");
826 AliTPCFCVoltError3D *corrField= (AliTPCFCVoltError3D*)f.Get("corrField");
827 AliTPCCorrection::AddVisualCorrection(alignLocal,1000);
828 AliTPCCorrection::AddVisualCorrection(alignQuadrant,1100);
829 AliTPCCorrection::AddVisualCorrection(corrField,100);
830 FitFunctionQA();
831 DrawAlignParam();
832 f.Close();
833}
834
835void FitFunctionQA(){
35d81915 836 ///
837
32b5322b 838 TH1 *his=0;
839 TCanvas *canvasDist= new TCanvas("FitQA","fitQA",1200,800);
840 canvasDist->Divide(2,2);
841 //
842 canvasDist->cd(1)->SetLogy(kFALSE); canvasDist->cd(1)->SetRightMargin(0.15);
843 treeDY->Draw("10*mean:sector:localX","kZ<0&&localX<134","colz");
844 his= (TH1*)treeDY->GetHistogram()->Clone(); his->SetName("DeltaRPhi1");
845 his->GetXaxis()->SetTitle("Sector");
846 his->GetZaxis()->SetTitle("R (cm)");
847 his->GetYaxis()->SetTitle("#Delta_{r#phi} (mm)");
848 his->Draw("cont2");
849 his->Draw("colz");
850 //
851 canvasDist->cd(2)->SetLogy(kFALSE); canvasDist->cd(2)->SetRightMargin(0.15);
852 treeDY->Draw("10*mean:sector:localX","kZ<0&&localX>160","colz");
853 his= (TH1*)treeDY->GetHistogram()->Clone(); his->SetName("DeltaRPhi2");
854 his->GetXaxis()->SetTitle("Sector");
855 his->GetZaxis()->SetTitle("R (cm)");
856 his->GetYaxis()->SetTitle("#Delta_{r#phi} (mm)");
857 his->Draw("cont2");
858 his->Draw("colz");
859 //
860 canvasDist->cd(3)->SetLogy(kFALSE); canvasDist->cd(3)->SetRightMargin(0.15);
861 treeDY->Draw("10*(mean-dAll):sector:localX","kZ<0&&localX<134","colz");
862 his= (TH1*)treeDY->GetHistogram()->Clone(); his->SetName("DeltaRPhi3"); his->SetTitle("Delta #RPhi");
863 his->GetXaxis()->SetTitle("Sector");
864 his->GetZaxis()->SetTitle("R (cm)");
865 his->GetYaxis()->SetTitle("#Delta_{r#phifit}-#Delta_{r#phi} (mm)");
866 his->Draw("cont2");
867 his->Draw("colz");
868 //
869 canvasDist->cd(4)->SetLogy(kFALSE); canvasDist->cd(4)->SetRightMargin(0.15);
870 treeDY->SetMarkerColor(1);
871 treeDY->Draw("10*mean:10*dAll:localX","","colz");
872 his= (TH1*)treeDY->GetHistogram()->Clone(); his->SetName("DeltaRPhi4");
873 his->GetXaxis()->SetTitle("Fit value #Delta_{r#phi} (mm)");
874 his->GetYaxis()->SetTitle("#Delta_{r#phi} (mm)");
875 his->GetZaxis()->SetTitle("R (cm)");
876 his->Draw("cont2");
877 his->Draw("colz");
878 //
879 canvasDist->Write("FitQA");
880}
881
882void DrawAlignParam(){
35d81915 883 ///
884
32b5322b 885 TFile f("fitAlignLookup.root");
886 TTree * treeAlign=(TTree*)f.Get("align");
887 TTree * treeQuadrant=(TTree*)f.Get("quadrant");
888 TCanvas *canvasAlign=new TCanvas("align","align",1000,800);
889 TH1 *his=0;
890 canvasAlign->Divide(2,2);
891 treeAlign->SetMarkerStyle(25);
892 treeQuadrant->SetMarkerStyle(25);
893 gStyle->SetOptStat(kTRUE);
894 //
895 canvasAlign->cd(1); canvasAlign->cd(1)->SetRightMargin(0.15);
896 treeAlign->Draw("10*dlyIROC*side:isec:1+side","","colz");
897 his= (TH1*)treeAlign->GetHistogram()->Clone(); his->SetName("dlyIROC");his->SetTitle("IROC Alignment #Delta_{r#phi}");
898 his->GetXaxis()->SetTitle("sector");
899 his->GetYaxis()->SetTitle("#Delta_{r#phi} (mm)");
900 his->GetZaxis()->SetTitle("side");
901 his->Draw("cont2");
902 his->Draw("colz");
903 //
904 canvasAlign->cd(2); canvasAlign->cd(2)->SetRightMargin(0.15);
905 treeAlign->Draw("1000*drotIROC*side:isec:1+side","","colz");
906 his= (TH1*)treeAlign->GetHistogram()->Clone(); his->SetName("drotIROC");his->SetTitle("IROC Angular Alignment #Delta_{r#phi}");
907 his->GetXaxis()->SetTitle("sector");
908 his->GetYaxis()->SetTitle("#Delta_{r#phi} (mrad)");
909 his->GetZaxis()->SetTitle("side");
910 his->Draw("cont2");
911 his->Draw("colz");
912 //
913 canvasAlign->cd(4);canvasAlign->cd(4)->SetRightMargin(0.15);
914 treeAlign->Draw("1000*drotOROC*side:isec:1+side","","colz");
915 his= (TH1*)treeAlign->GetHistogram()->Clone(); his->SetName("drotOROC");his->SetTitle("OROC Angular Alignment #Delta_{r#phi}");
916 his->GetXaxis()->SetTitle("sector");
917 his->GetYaxis()->SetTitle("#Delta_{r#phi} (mrad)");
918 his->GetZaxis()->SetTitle("side");
919 his->Draw("cont2");
920 his->Draw("colz");
921 //
922 canvasAlign->cd(3);canvasAlign->cd(3)->SetRightMargin(0.15);
923 treeQuadrant->Draw("10*q2*side:isec:1+side","","colz");
924 his= (TH1*)treeQuadrant->GetHistogram()->Clone(); his->SetName("drphiOROC");his->SetTitle("OROC Alignment Outer Quadrant #Delta_{r#phi}");
925 his->GetXaxis()->SetTitle("sector");
926 his->GetYaxis()->SetTitle("#Delta_{r#phi} (mm)");
927 his->GetZaxis()->SetTitle("side");
928 his->Draw("cont2");
929 his->Draw("colz");
930
931
932}