2 Unit test for some functions classes used in the $ALICE_ROOT/TPC/Base directory:
4 gSystem->SetIncludePath("-I$ROOTSYS/include -I$ALICE_ROOT/ -I$ALICE_ROOT/install/include -I$ALICE_ROOT/STEER -I$ALICE_ROOT/ITS -I$ALICE_ROOT/TRD -I$ALICE_ROOT/TOF -I$ALICE_ROOT/RAW -I$ALICE_ROOT/STAT -I$ALICE_ROOT/TPC/TPCbase -I$ALICE_ROOT/TPCcalib");
7 .L $ALICE_ROOT/TPC/Base/test/UnitTest.C+
8 UnitTestAliTPCCalPadTree();
9 TestCorrection_AliTPCCorrection_AddCorrectionCompact();
15 #include "TLinearFitter.h"
17 #include "AliSysInfo.h"
22 #include "AliTPCCalPad.h"
23 #include "AliTPCCalibViewer.h"
24 #include "AliTPCcalibDButil.h"
25 #include "AliTPCCorrection.h"
26 #include "AliTPCComposedCorrection.h"
27 #include "AliTPCExBTwist.h"
28 #include "AliTPCFCVoltError3D.h"
29 #include "AliTPCROCVoltError3D.h"
30 #include "AliTPCBoundaryVoltError.h"
31 #include "AliTPCCalibGlobalMisalignment.h"
32 #include "AliCDBEntry.h"
33 #include "TStopwatch.h"
34 #include "TGeoMatrix.h"
35 #include "TGeoGlobalMagField.h"
38 // PARAMETERS to set from outside:
40 TString baseDir="/hera/alice/wiechula/calib/guiTrees"; // TO FIX specification of inout data
44 Bool_t TestCorrection_AliTPCCorrection_AddCorrectionCompact();
45 Bool_t TestCorrection_AliTPCExBTwistAddCorrectionCompact();
46 Bool_t TestCorrection_AliTPCFCVoltError3DAddCorrectionCompact();
47 Bool_t TestCorrection_AliTPCRocVoltError3DAddCorrectionCompact();
48 Bool_t TestCorrection_AliTPCBoundaryVoltErrorAddCorrectionCompact();
49 Bool_t TestCorrection_AliTPCCalibGlobalMisalignmentAddCorrectionCompact();
50 Bool_t TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact();
51 Bool_t TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact_TPCCalibCorrection(Bool_t fast=kFALSE);
53 void UnitTestAliTPCCalPadTree(){
55 // Make a UnitTest of the AliTPCCalPad
56 // a.) TTree functionaility
57 // b.) MedianFilterFunctionality
58 // c.) LTMFilterFunctionality
60 TObjArray *fArray = new TObjArray(100);
61 TTree * treePad=AliTPCcalibDButil::ConnectGainTrees(baseDir);
62 for (Int_t i=0; i<5; i+=2){
63 AliTPCCalPad * padLx = AliTPCCalPad::MakePadFromTree(treePad,"lx.fElements","Lx",kTRUE);
64 AliTPCCalPad * padLy = AliTPCCalPad::MakePadFromTree(treePad,"ly.fElements","Ly",kTRUE);
65 AliTPCCalPad * padLLx = AliTPCCalPad::MakePadFromTree(treePad,"lx.fElements","LLx",kTRUE);
66 AliTPCCalPad * padLLy = AliTPCCalPad::MakePadFromTree(treePad,"ly.fElements","LLy",kTRUE);
67 AliTPCCalPad * padMax = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MaxCharge.fElements","QMax",kTRUE);
68 AliTPCCalPad * padMean = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MeanCharge.fElements","QTot",kTRUE);
69 AliTPCCalPad * padMaxL = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MaxCharge.fElements","QMax",kTRUE);
70 AliTPCCalPad * padMeanL = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MeanCharge.fElements","QTot",kTRUE);
72 padLx->MedianFilter(i,2*i);
73 padLy->MedianFilter(i,2*i);
74 padLLx->LTMFilter(i,2*i,1.00, 0);
75 padLLy->LTMFilter(i,2*i,1.00, 0);
76 padMax->MedianFilter(i,2*i);
77 padMean->MedianFilter(i,2*i);
78 padMaxL->LTMFilter(i,2*i,0.8,0);
79 padMeanL->LTMFilter(i,2*i,0.8,0);
81 padLx->SetName(TString::Format("Lx%d",i).Data());
82 padLy->SetName(TString::Format("Ly%d",i).Data());
83 padLLx->SetName(TString::Format("LLx%d",i).Data());
84 padLLy->SetName(TString::Format("LLy%d",i).Data());
85 padMax->SetName(TString::Format("QMax%d",i).Data());
86 padMean->SetName(TString::Format("QTot%d",i).Data());
87 padMaxL->SetName(TString::Format("QMaxL%d",i).Data());
88 padMeanL->SetName(TString::Format("QTotL%d",i).Data());
89 fArray->AddLast(padLx);
90 fArray->AddLast(padLy);
91 fArray->AddLast(padLLx);
92 fArray->AddLast(padLLy);
93 fArray->AddLast(padMax);
94 fArray->AddLast(padMean);
95 fArray->AddLast(padMaxL);
96 fArray->AddLast(padMeanL);
98 AliTPCCalibViewer::MakeTree("QAtest.root", fArray,0);
100 // 2.) Check invariants
102 TFile*fout= TFile::Open("QAtest.root");
103 TTree * tree = (TTree*)fout->Get("calPads");
104 Int_t isOutM0 = tree->Draw("(Ly2.fElements-Ly0.fElements)>>his0(100,-10,10)","abs((Ly2.fElements-Ly0.fElements))>2","goff");
105 Int_t isOutM1=tree->Draw("(Lx2.fElements-Lx0.fElements)/0.75>>his1(100,-10,10)","abs((Lx2.fElements-Lx0.fElements))>0","goff");
106 printf("IsOut=%d\t%d\n",isOutM0,isOutM1);
107 if ((isOutM0+isOutM1)==0) ::Info("UnitTestAliTPCCalPadTree","MedianTest OK");
108 if (isOutM0||isOutM1) ::Fatal("UnitTestAliTPCCalPadTree","MedianTest FAILED");
110 Int_t isOutL0 = tree->Draw("(LLy2.fElements-Ly0.fElements)>>his0(100,-10,10)","abs((LLy2.fElements-LLy0.fElements))>0","goff");
111 Int_t isOutL1=tree->Draw("(LLx2.fElements-Lx0.fElements)/0.75>>his1(100,-10,10)","abs((LLx2.fElements-LLx0.fElements))>0","goff");
112 printf("IsOut=%d\t%d\n",isOutL0,isOutL1);
113 if ((isOutL0+isOutL1)==0) ::Info("UnitTestAliTPCCalPadTree","LTMTest OK");
114 if (isOutL0||isOutL1) ::Fatal("UnitTestAliTPCCalPadTree","LTMTest FAILED");
118 Bool_t TestCorrection_AliTPCExBTwistAddCorrectionCompact(){
121 // 1.) Test ExB twist AddCorrectionCompact
123 Bool_t isOK[10]={kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
124 AliTPCComposedCorrection *compCorrTwist = new AliTPCComposedCorrection();
125 AliTPCExBTwist *twistX = new AliTPCExBTwist;
126 AliTPCExBTwist *twistY = new AliTPCExBTwist;
127 twistX->SetXTwist(0.001); // 1 mrad twist in x
128 twistY->SetYTwist(0.001); // 1 mrad twist in x
129 isOK[0]&=compCorrTwist->AddCorrectionCompact(twistX,0.5);
130 isOK[0]&=compCorrTwist->AddCorrectionCompact(twistY,0.5);
131 isOK[0]&=compCorrTwist->AddCorrectionCompact(twistX,0.5);
132 isOK[0]&=compCorrTwist->AddCorrectionCompact(twistY,0.5);
133 isOK[0]&=compCorrTwist->AddCorrectionCompact(twistY,-1);
134 isOK[0]&=compCorrTwist->AddCorrectionCompact(twistX,-1);
135 isOK[1]=compCorrTwist->GetCorrections()->GetEntries()==1;
137 AliTPCExBTwist *twistRes=0;
138 if (isOK[1]==kFALSE){
143 twistRes= dynamic_cast<AliTPCExBTwist *>(compCorrTwist->GetSubCorrection(0));
149 isOK[3] &= (twistRes->GetXTwist()==0);
150 isOK[4] &= (twistRes->GetYTwist()==0);
154 for (Int_t i=0; i<5; i++) res&=isOK[i];
156 if (isOK[0]==kFALSE){
157 ::Error("TestCorrection_AddCorrectionCompact","AliTPCExBTwist -ADD FAILED");
159 ::Info("TestCorrection_AddCorrectionCompact","AliTPCExBTwist -ADD OK");
161 if (isOK[1]==kFALSE){
162 ::Error("TestCorrection_AddCorrectionCompact","AliTPCExBTwist - wrong entries FAILED");
164 ::Info("TestCorrection_AddCorrectionCompact","AliTPCExBTwist - entries OK");
166 if (isOK[2]==kFALSE || isOK[3]==kFALSE ||isOK[4]==kFALSE ){
167 ::Error("TestCorrection_AddCorrectionCompact","AliTPCExBTwist - inconsitent entries FAILED");
169 ::Info("TestCorrection_AddCorrectionCompact","AliTPCExBTwist - consistent entries OK");
176 Bool_t TestCorrection_AliTPCFCVoltError3DAddCorrectionCompact(){
178 // TestCorrection_AliTPCFCVoltError3DAddCorrectionCompact
180 const Float_t kEpsilon=0.000001;
181 Bool_t isOK[10]={kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
182 AliTPCComposedCorrection *compCorrComp = new AliTPCComposedCorrection();
183 AliTPCFCVoltError3D *corr0 = new AliTPCFCVoltError3D;
184 AliTPCFCVoltError3D *corr1 = new AliTPCFCVoltError3D;
185 for (Int_t isec=0; isec<36; isec++){
186 corr0->SetRodVoltShiftA(isec,TMath::Cos(TMath::Pi()*isec/36),kFALSE);
187 corr0->SetRodVoltShiftC(isec,TMath::Cos(TMath::Pi()*isec/36),kFALSE);
188 corr1->SetRodVoltShiftA(isec,TMath::Sin(TMath::Pi()*isec/36),kFALSE);
189 corr1->SetRodVoltShiftC(isec,TMath::Sin(TMath::Pi()*isec/36),kFALSE);
190 corr1->SetCopperRodShiftA(isec,TMath::Sin(TMath::Pi()*isec/36),kFALSE);
191 corr1->SetCopperRodShiftC(isec,TMath::Sin(TMath::Pi()*isec/36),kFALSE);
194 isOK[0]&=compCorrComp->AddCorrectionCompact(corr0,0.5);
195 isOK[0]&=compCorrComp->AddCorrectionCompact(corr1,0.5);
196 isOK[0]&=compCorrComp->AddCorrectionCompact(corr0,0.5);
197 isOK[0]&=compCorrComp->AddCorrectionCompact(corr1,0.5);
198 isOK[0]&=compCorrComp->AddCorrectionCompact(corr1,-1);
199 isOK[0]&=compCorrComp->AddCorrectionCompact(corr0,-1);
200 isOK[1]=compCorrComp->GetCorrections()->GetEntries()==1;
201 AliTPCFCVoltError3D *corrRes=0;
202 if (isOK[1]==kFALSE){
207 corrRes= dynamic_cast<AliTPCFCVoltError3D *>(compCorrComp->GetSubCorrection(0));
213 for (Int_t isec=0; isec<36; isec++){
214 isOK[3] &=( TMath::Abs(corrRes->GetRodVoltShiftA(isec))<kEpsilon);
215 isOK[4] &=( TMath::Abs(corrRes->GetRodVoltShiftC(isec))<kEpsilon);
216 isOK[5] &=( TMath::Abs(corrRes->GetCopperRodShiftA(isec))<kEpsilon);
217 isOK[6] &=( TMath::Abs(corrRes->GetCopperRodShiftC(isec))<kEpsilon);
222 for (Int_t i=0; i<5; i++) res&=isOK[i];
224 if (isOK[0]==kFALSE){
225 ::Error("TestCorrection_AddCorrectionCompact","AliTPCFCVoltError3D -ADD FAILED");
227 ::Info("TestCorrection_AddCorrectionCompact","AliTPCFCVoltError3D -ADD OK");
229 if (isOK[1]==kFALSE){
230 ::Error("TestCorrection_AddCorrectionCompact","AliTPCFCVoltError3D - wrong entries FAILED");
232 ::Info("TestCorrection_AddCorrectionCompact","AliTPCFCVoltError3D - entries OK");
234 if (isOK[2]==kFALSE || isOK[3]==kFALSE ||isOK[4]==kFALSE ){
235 ::Error("TestCorrection_AddCorrectionCompact","AliTPCFCVoltError3D - inconsitent entries FAILED");
237 ::Info("TestCorrection_AddCorrectionCompact","AliTPCFCVoltError3D - consistent entries OK");
245 Bool_t TestCorrection_AliTPCRocVoltError3DAddCorrectionCompact(){
247 // AliTPCRocVoltError3DAddCorrectionCompact
249 const Float_t kEpsilon=0.00000001;
250 Bool_t isOK[10]={kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
251 AliTPCComposedCorrection *compCorrROCVoltError3D = new AliTPCComposedCorrection();
252 AliTPCROCVoltError3D *corr0 = new AliTPCROCVoltError3D;
253 AliTPCROCVoltError3D *corr1 = new AliTPCROCVoltError3D;
254 TMatrixD matrixDz(72,3);
255 for (Int_t isec=0; isec<72; isec++){
256 matrixDz(isec,0)=gRandom->Rndm()*0.1;
257 matrixDz(isec,1)=gRandom->Rndm()*0.001;
258 matrixDz(isec,2)=gRandom->Rndm()*0.001;
260 corr0->SetROCData(&matrixDz);
262 corr1->SetROCData(&matrixDz);
264 isOK[0]&=compCorrROCVoltError3D->AddCorrectionCompact(corr0,0.5);
265 isOK[0]&=compCorrROCVoltError3D->AddCorrectionCompact(corr1,0.5);
266 isOK[0]&=compCorrROCVoltError3D->AddCorrectionCompact(corr0,0.5);
267 isOK[0]&=compCorrROCVoltError3D->AddCorrectionCompact(corr1,0.5);
268 isOK[0]&=compCorrROCVoltError3D->AddCorrectionCompact(corr1,-1);
269 isOK[0]&=compCorrROCVoltError3D->AddCorrectionCompact(corr0,-1);
270 isOK[1]=compCorrROCVoltError3D->GetCorrections()->GetEntries()==1;
271 AliTPCROCVoltError3D *corrRes=0;
272 if (isOK[1]==kFALSE){
277 corrRes= dynamic_cast<AliTPCROCVoltError3D *>(compCorrROCVoltError3D->GetSubCorrection(0));
283 isOK[3]=TMath::Abs(corrRes->GetMatrix()->Sum())<kEpsilon;
287 for (Int_t i=0; i<5; i++) res&=isOK[i];
289 if (isOK[0]==kFALSE){
290 ::Error("TestCorrection_AddCorrectionCompact","AliTPCROCVoltError3D -ADD FAILED");
292 ::Info("TestCorrection_AddCorrectionCompact","AliTPCROCVoltError3D -ADD OK");
294 if (isOK[1]==kFALSE){
295 ::Error("TestCorrection_AddCorrectionCompact","AliTPCROCVoltError3D - wrong entries FAILED");
297 ::Info("TestCorrection_AddCorrectionCompact","AliTPCROCVoltError3D - entries OK");
299 if (isOK[2]==kFALSE || isOK[3]==kFALSE ||isOK[4]==kFALSE ){
300 ::Error("TestCorrection_AddCorrectionCompact","AliTPCROCVoltError3D - inconsitent entries FAILED");
302 ::Info("TestCorrection_AddCorrectionCompact","AliTPCROCVoltError3D - consistent entries OK");
312 Bool_t TestCorrection_AliTPCBoundaryVoltErrorAddCorrectionCompact(){
314 // AliTPCBoundaryVoltErrorAddCorrectionCompact
316 const Float_t kEpsilon=0.00000001;
317 Bool_t isOK[10]={kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
318 AliTPCComposedCorrection *compCorrBoundaryVoltError = new AliTPCComposedCorrection();
319 AliTPCBoundaryVoltError *corr0 = new AliTPCBoundaryVoltError;
320 AliTPCBoundaryVoltError *corr1 = new AliTPCBoundaryVoltError;
321 Float_t boundaries[8];
322 for (Int_t ibound=0; ibound<8; ibound++){
323 boundaries[ibound]=gRandom->Rndm()-0.5;
325 corr0->SetBoundariesA(boundaries);
326 corr1->SetBoundariesA(boundaries);
328 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr0,0.5);
329 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr1,0.5);
330 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr0,0.5);
331 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr1,0.5);
332 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr1,-1);
333 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr0,-1);
334 isOK[1]=compCorrBoundaryVoltError->GetCorrections()->GetEntries()==1;
335 AliTPCBoundaryVoltError *corrRes=0;
336 if (isOK[1]==kFALSE){
341 corrRes= dynamic_cast<AliTPCBoundaryVoltError *>(compCorrBoundaryVoltError->GetSubCorrection(0));
347 for (Int_t ibound=0; ibound<8; ibound++){
348 isOK[3]&=TMath::Abs(corrRes->GetBoundariesA(ibound))<kEpsilon;
349 isOK[3]&=TMath::Abs(corrRes->GetBoundariesC(ibound))<kEpsilon;
354 for (Int_t i=0; i<5; i++) res&=isOK[i];
356 if (isOK[0]==kFALSE){
357 ::Error("TestCorrection_AddCorrectionCompact","AliTPCBoundaryVoltError -ADD FAILED");
359 ::Info("TestCorrection_AddCorrectionCompact","AliTPCBoundaryVoltError -ADD OK");
361 if (isOK[1]==kFALSE){
362 ::Error("TestCorrection_AddCorrectionCompact","AliTPCBoundaryVoltError - wrong entries FAILED");
364 ::Info("TestCorrection_AddCorrectionCompact","AliTPCBoundaryVoltError - entries OK");
366 if (isOK[2]==kFALSE || isOK[3]==kFALSE ||isOK[4]==kFALSE ){
367 ::Error("TestCorrection_AddCorrectionCompact","AliTPCBoundaryVoltError - inconsitent entries FAILED");
369 ::Info("TestCorrection_AddCorrectionCompact","AliTPCBoundaryVoltError - consistent entries OK");
379 Bool_t TestCorrection_AliTPCCalibGlobalMisalignmentAddCorrectionCompact(){
381 // AliTPCCalibGlobalMisalignmentAddCorrectionCompact
382 // Invariant used in test is not exact it is only approximate - as matrix multiplication is not comulative
383 // !!!! BUG FOUND ????
384 // hmatrix1->GetTranslation()[idelta]=xxx; // does not work as expected Translation is set, visible in Print but not used later
385 const Float_t kEpsilon=0.0001;
386 Bool_t isOK[10]={kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
387 Double_t delta[3]={0.01,0.02,0.03};
389 AliTPCComposedCorrection *compCorrBoundaryVoltError = new AliTPCComposedCorrection();
390 AliTPCCalibGlobalMisalignment *corr0 = new AliTPCCalibGlobalMisalignment;
391 AliTPCCalibGlobalMisalignment *corr1 = new AliTPCCalibGlobalMisalignment;
392 AliTPCCalibGlobalMisalignment *corr2 = new AliTPCCalibGlobalMisalignment;
393 TObjArray sectorAlign(72);
395 TGeoHMatrix *hmatrix0 = new TGeoHMatrix;
396 TGeoHMatrix *hmatrix1 = new TGeoHMatrix;
397 hmatrix0->RotateX(TMath::RadToDeg()*0.0001);
398 hmatrix0->RotateY(TMath::RadToDeg()*0.0002);
399 hmatrix0->RotateZ(TMath::RadToDeg()*0.0003);
400 hmatrix1->SetTranslation(delta);
401 for (Int_t isec=0; isec<72; isec++){
402 if ((isec%2)==0) sectorAlign.AddAt(hmatrix0,isec);
403 if ((isec%2)==1) sectorAlign.AddAt(hmatrix1,isec);
405 corr0->SetAlignGlobal(hmatrix0);
406 corr1->SetAlignGlobal(hmatrix1);
407 corr0->SetAlignGlobalDelta(hmatrix1);
408 corr1->SetAlignGlobalDelta(hmatrix0);
409 corr2->SetAlignSectors(§orAlign);
411 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr0,0.5);
412 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr1,0.5);
413 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr0,0.5);
414 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr1,0.5);
415 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr0,-1);
416 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr1,-1);
418 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr2,1);
419 isOK[0]&=compCorrBoundaryVoltError->AddCorrectionCompact(corr2,-1);
422 isOK[1]=compCorrBoundaryVoltError->GetCorrections()->GetEntries()==1;
423 AliTPCCalibGlobalMisalignment *corrRes=0;
424 if (isOK[1]==kFALSE){
429 corrRes= dynamic_cast<AliTPCCalibGlobalMisalignment *>(compCorrBoundaryVoltError->GetSubCorrection(0));
435 for (Int_t itrans=0; itrans<3; itrans++){
436 isOK[2+itrans]&=TMath::Abs(corrRes->GetAlignGlobal()->GetTranslation()[itrans])<kEpsilon;
437 for (Int_t isec=0; isec<72; isec++){
438 isOK[2+itrans]&=TMath::Abs(((TGeoHMatrix*)(corrRes->GetAlignSectors()->At(isec)))->GetTranslation()[itrans])<kEpsilon;
441 corrRes->GetAlignGlobal()->Print();
442 corrRes->GetAlignSectors()->At(0)->Print();
443 corrRes->GetAlignSectors()->At(1)->Print();
447 for (Int_t i=0; i<5; i++) res&=isOK[i];
449 if (isOK[0]==kFALSE){
450 ::Error("TestCorrection_AddCorrectionCompact","AliTPCCalibGlobalMisalignment -ADD FAILED");
452 ::Info("TestCorrection_AddCorrectionCompact","AliTPCCalibGlobalMisalignment -ADD OK");
454 if (isOK[1]==kFALSE){
455 ::Error("TestCorrection_AddCorrectionCompact","AliTPCCalibGlobalMisalignment - wrong entries FAILED");
457 ::Info("TestCorrection_AddCorrectionCompact","AliTPCCalibGlobalMisalignment - entries OK");
459 if (isOK[2]==kFALSE || isOK[3]==kFALSE ||isOK[4]==kFALSE ){
460 ::Error("TestCorrection_AddCorrectionCompact","AliTPCCalibGlobalMisalignment - inconsitent entries FAILED");
462 ::Info("TestCorrection_AddCorrectionCompact","AliTPCCalibGlobalMisalignment - consistent entries OK");
468 Bool_t TestCorrection_AliTPCCorrection_AddCorrectionCompact(){
472 TestCorrection_AliTPCExBTwistAddCorrectionCompact();
473 TestCorrection_AliTPCFCVoltError3DAddCorrectionCompact();
474 TestCorrection_AliTPCRocVoltError3DAddCorrectionCompact();
475 TestCorrection_AliTPCBoundaryVoltErrorAddCorrectionCompact();
476 TestCorrection_AliTPCCalibGlobalMisalignmentAddCorrectionCompact();
479 Bool_t TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact_TPCCalibCorrection(Bool_t fast){
483 const Int_t npointsTest=10000;
484 const Float_t kEpsilon=0.001; //10 microns
485 TGeoGlobalMagField::Instance()->SetField(new AliMagF("Maps","Maps", -1., -1., AliMagF::k5kG));
487 // 0.) Read an input OCDB entry
489 TFile * f = TFile::Open("$ALICE_OCDB/alice/data/2010/OCDB/TPC/Calib/Correction/Run0_999999999_v8_s0.root");
490 AliCDBEntry * entry=(AliCDBEntry*)f->Get("AliCDBEntry");
491 TObjArray * corrArray = (TObjArray *)entry->GetObject();
492 AliTPCComposedCorrection *compInput = (AliTPCComposedCorrection *)corrArray->At(0);
493 AliTPCComposedCorrection *compInputFast = new AliTPCComposedCorrection;
494 Int_t ncorrs = compInput->GetCorrections()->GetEntries();
495 TObjArray arrayInputFast(ncorrs);
497 // 1.) Test each individual correction
499 for (Int_t icorr=0; icorr<ncorrs; icorr++){
500 TString clName=compInput->GetSubCorrection(icorr)->IsA()->GetName();
501 if (fast){ // skip slow correction
502 if ( clName.Contains("AliTPCFCVoltError3D"))continue;
503 if ( clName.Contains("AliTPCROCVoltError3D"))continue;
505 // if ( clName.Contains("AliTPCExBBShape"))continue;
506 AliTPCCorrection *corrInput=compInput->GetSubCorrection(icorr);
508 ::Info("TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact_TPCCalibCorrection",TString::Format("%s\t%s",corrInput->IsA()->GetName(),corrInput->GetName()).Data());
509 AliTPCComposedCorrection *compTest0= new AliTPCComposedCorrection;
510 AliTPCComposedCorrection *compTest1= new AliTPCComposedCorrection;
511 compTest0->AddCorrectionCompact(corrInput,0.5);
512 compTest0->AddCorrectionCompact(corrInput,0.5);
513 compTest1->AddCorrectionCompact(corrInput,1);
514 compTest1->AddCorrectionCompact(corrInput,-1);
515 corrInput->AddVisualCorrection(corrInput,10);
516 compTest0->AddVisualCorrection(compTest0,11);
517 compTest1->AddVisualCorrection(compTest1,12);
518 compTest0->SetOmegaTauT1T2(0.35,1,1);
519 compTest1->SetOmegaTauT1T2(0.35,1,1);
520 corrInput->SetOmegaTauT1T2(0.35,1,1);
521 for (Int_t icoord=0; icoord<3; icoord++){
522 TVectorD dvecTest0(npointsTest);
523 TVectorD dvecTest1(npointsTest);
524 for (Int_t ipoint=0; ipoint<npointsTest; ipoint++){
525 Double_t r= 85.+gRandom->Rndm()*150;
526 Double_t phi= gRandom->Rndm()*TMath::TwoPi();
527 Double_t z=500*(gRandom->Rndm()-0.5);
528 dvecTest0[ipoint]=AliTPCCorrection::GetCorrXYZ(r*TMath::Cos(phi),r*TMath::Sin(phi),z, icoord, 11)-AliTPCCorrection::GetCorrXYZ(r*TMath::Cos(phi),r*TMath::Sin(phi),z, icoord, 10);
529 dvecTest1[ipoint]=AliTPCCorrection::GetCorrXYZ(r*TMath::Cos(phi),r*TMath::Sin(phi),z, icoord, 12);
531 Double_t mean0 = TMath::Mean(npointsTest, dvecTest0.GetMatrixArray());
532 Double_t rms0 = TMath::RMS(npointsTest, dvecTest0.GetMatrixArray());
533 Double_t mean1 = TMath::Mean(npointsTest, dvecTest1.GetMatrixArray());
534 Double_t rms1 = TMath::RMS(npointsTest, dvecTest1.GetMatrixArray());
535 if (TMath::Abs(rms0)>kEpsilon){
536 ::Error("TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact_TPCCalibCorrection",TString::Format("Test0:\t%s\t%3.5f\t%3.5f FAILED",clName.Data(),mean0,rms0).Data());
538 ::Info("TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact_TPCCalibCorrection",TString::Format("Test0:\t%s\t%3.5f\t%3.5f OK",clName.Data(),mean0,rms0).Data());
540 if (TMath::Abs(rms1)>kEpsilon){
541 ::Error("TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact_TPCCalibCorrection",TString::Format("Test1:\t%s\t%3.5f\t%3.5f FAILED",clName.Data(),mean1,rms1).Data());
543 ::Info("TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact_TPCCalibCorrection",TString::Format("Test1:\t%s\t%3.5f\t%3.5f OK",clName.Data(),mean1,rms1).Data());
550 Bool_t TestCorrection_AliTPCComposedCorrectionAddCorrectionCompact(){
552 // Tests of AliTPCComposedCorrection
553 // 1.) Make linear combination correction example using weights.
554 // Test correction checking invariant inverse x orig (there are simpler way to do inversion using AliTPCInverseCorrection)
556 // 2.) Make compact for of the Composed correction. Test correction checking invariant inverse x orig
558 const Int_t npointsTest=10000;
559 const Float_t kEpsilon=0.0001; // using Floating point precission
561 // 0.) Read an input OCDB entry
563 TFile * f = TFile::Open("$ALICE_OCDB/alice/data/2010/OCDB/TPC/Calib/Correction/Run0_999999999_v8_s0.root");
564 AliCDBEntry * entry=(AliCDBEntry*)f->Get("AliCDBEntry");
565 TObjArray * corrArray = (TObjArray *)entry->GetObject();
566 AliTPCComposedCorrection *compInput = (AliTPCComposedCorrection *)corrArray->At(0);
567 AliTPCComposedCorrection *compInputFast = new AliTPCComposedCorrection;
568 Int_t ncorrs = compInput->GetCorrections()->GetEntries();
569 TObjArray arrayInputFast(ncorrs);
570 for (Int_t icorr=0; icorr<ncorrs; icorr++){
571 TString clName=compInput->GetSubCorrection(icorr)->IsA()->GetName();
572 if ( clName.Contains("AliTPCFCVoltError3D"))continue;
573 if ( clName.Contains("AliTPCROCVoltError3D"))continue;
574 if ( clName.Contains("AliTPCExBBShape"))continue;
575 arrayInputFast.AddLast(compInput->GetSubCorrection(icorr));
577 compInputFast->SetCorrections(&arrayInputFast);
582 // 1.) Make linear combination correction example using weights.
583 // Test correction checking invariant inverse x orig (there are simpler way to do inversion using AliTPCInverseCorrection)
584 AliTPCComposedCorrection *compInverse = new AliTPCComposedCorrection();
585 Bool_t isOK1[10]={kTRUE};
586 TObjArray * collection= dynamic_cast<TObjArray*>(compInput->GetCorrections());
587 Int_t entries = collection->GetEntries();
588 TVectorD weights(entries+1);
589 for (Int_t i=0; i<entries+1; i++) weights[i]=-1.0;
591 TObjArray * arrayInvariant = new TObjArray(entries+1);
592 arrayInvariant->AddLast(compInput);
593 for (Int_t i=0; i<entries; i++) arrayInvariant->AddLast( collection->At(i));
594 compInverse->SetCorrections( arrayInvariant);
595 compInverse->SetWeights(&weights);
596 compInverse->AddVisualCorrection(compInverse,1);
597 compInput->AddVisualCorrection(compInput,2);
598 TF1 finv1("finv1","AliTPCCorrection::GetCorrXYZ(x,x,100,0,1)",85,245);
600 TVectorD vecCompInverse(npointsTest);
601 for (Int_t icoord=0; icoord<3; icoord++){
602 for (Int_t ipoint=0; ipoint<npointsTest; ipoint++){
603 Double_t r= 85.+gRandom->Rndm()*150;
604 Double_t phi= gRandom->Rndm()*TMath::TwoPi();
605 Double_t z=500*(gRandom->Rndm()-0.5);
606 vecCompInverse[ipoint]=AliTPCCorrection::GetCorrXYZ(r*TMath::Cos(phi),r*TMath::Sin(phi),z, icoord, 1);
608 Double_t rms=TMath::RMS(npointsTest,vecCompInverse.GetMatrixArray());
609 Double_t mean=TMath::Mean(npointsTest,vecCompInverse.GetMatrixArray());
610 isOK1[icoord]=TMath::Abs(rms)<kEpsilon;
611 isOK1[icoord]&=TMath::Abs(mean)<kEpsilon;
613 if (isOK1[0]==kFALSE || isOK1[1]==kFALSE ||isOK1[2]==kFALSE ){
614 ::Error("TestCorrection_AddCorrectionCompact",TString::Format("AliTPCComposedCorrection - Test1 (%d,%d,%d) FAILED",isOK1[0], isOK1[1],isOK1[2]).Data());
616 ::Info("TestCorrection_AddCorrectionCompact","AliTPCComposedCorrection - Test1 OK");
619 // 2.) Make compact for of the Composed correction. Test correction checking invariant inverse x orig
620 // This take time - dostortion has to be recalculated
621 AliTPCComposedCorrection *compOutInverseCompact = new AliTPCComposedCorrection();
622 compOutInverseCompact->AddCorrectionCompact(compInputFast,1);
623 compOutInverseCompact->AddCorrectionCompact(compInputFast,-1);
624 compOutInverseCompact->SetOmegaTauT1T2(0,1,1);
625 compInputFast->SetOmegaTauT1T2(0,1,1);
626 compOutInverseCompact->AddVisualCorrection(compOutInverseCompact,10);
627 compInputFast->AddVisualCorrection(compInput,3);
630 TF1 fcomp("fcomp","AliTPCCorrection::GetCorrXYZ(x,x,100,0,10)",85,245);
631 TF1 forig("forig","-AliTPCCorrection::GetCorrXYZ(x,x,100,0,3)",85,245);
632 TF1 fdiff("fdiff","AliTPCCorrection::GetCorrXYZ(x,x,100,0,10)+AliTPCCorrection::GetCorrXYZ(x,x,100,0,2)",85,245);