ATO-98 AddCorrectionCompact - to make corrections compact in memory + UnitTests...
authormivanov <marian.ivanov@cern.ch>
Sat, 13 Dec 2014 23:21:19 +0000 (00:21 +0100)
committermivanov <marian.ivanov@cern.ch>
Sat, 13 Dec 2014 23:24:59 +0000 (00:24 +0100)
TPC/Base/AliTPCComposedCorrection.cxx
TPC/Base/AliTPCComposedCorrection.h
TPC/Base/AliTPCCorrection.cxx
TPC/Base/AliTPCCorrection.h
TPC/Base/AliTPCExBTwist.cxx
TPC/Base/AliTPCExBTwist.h
TPC/Base/AliTPCFCVoltError3D.cxx
TPC/Base/AliTPCFCVoltError3D.h
TPC/Base/test/UnitTest.C

index 00804e6..d7d396c 100644 (file)
@@ -109,6 +109,48 @@ AliTPCComposedCorrection::~AliTPCComposedCorrection() {
   if (fWeights) delete fWeights;
 }
 
+
+Bool_t AliTPCComposedCorrection::AddCorrectionCompact(AliTPCCorrection* corr, Double_t weight){
+  //
+  // Add correction  - better name needed (left/right) - for the moment I assumme they commute
+  // Why not to just use array of corrections - CPU consideration
+  // Assumptions:
+  //  - origin of distortion/correction are additive
+  //  - corrections/distortion are small   and they commute
+  //  - only correction ot the same type supported 
+  const Int_t knCorr=100;
+  if (corr==NULL) {
+    AliError("Zerro pointer - correction");
+    return kFALSE;
+  }
+  if (!fCorrections) fCorrections= new TObjArray(knCorr);
+  // 1.) Case of  Composed correction
+  AliTPCComposedCorrection * corrC = dynamic_cast<AliTPCComposedCorrection *>(corr);
+  if (corrC){
+    Int_t ncorrs= corrC->fCorrections->GetEntries();
+    Bool_t isOK=kTRUE;
+    for (Int_t icorr=0; icorr<ncorrs; icorr++){
+      isOK&=AddCorrectionCompact(corrC->GetSubCorrection(icorr),weight*(*((corrC->fWeights)))[icorr]);
+    }
+    return isOK;
+  }
+  // 2.) Find subcorrection of the same type
+  AliTPCCorrection * toAdd=0;
+  Int_t ncorr=fCorrections->GetEntries();
+  for (Int_t icorr=0; icorr<ncorr; icorr++){
+    if (GetSubCorrection(icorr)==NULL)  continue;
+    if (GetSubCorrection(icorr)->IsA()==corr->IsA())  toAdd=GetSubCorrection(icorr);
+  }
+  // 3.) create of givent type  if does not exist 
+  if (toAdd==NULL){
+    toAdd= (AliTPCCorrection*)((corr->IsA())->New());
+    fCorrections->Add(toAdd);
+  }
+  // 4.) add to object of given type 
+  return toAdd->AddCorrectionCompact(corr, weight);
+}
+
+
 AliTPCCorrection * AliTPCComposedCorrection::GetSubCorrection(Int_t ipos){
   //
   //
index 4e6c0a8..f1ba877 100644 (file)
@@ -39,7 +39,7 @@ public:
   AliTPCComposedCorrection();
   AliTPCComposedCorrection(TCollection *corrections,CompositionType mode);
   virtual ~AliTPCComposedCorrection();
-
+  virtual Bool_t AddCorrectionCompact(AliTPCCorrection* corr, Double_t weight);
   void SetOmegaTauT1T2(Float_t omegaTau,Float_t t1,Float_t t2);
 
   TCollection* GetCorrections() const {return fCorrections;}
index 20059bf..fc7d120 100644 (file)
@@ -176,6 +176,21 @@ AliTPCCorrection::~AliTPCCorrection() {
   //
 }
 
+Bool_t AliTPCCorrection::AddCorrectionCompact(AliTPCCorrection* corr, Double_t weight){
+  //
+  // Add correction  and make them compact
+  // Assumptions:
+  //  - origin of distortion/correction are additive
+  //  - only correction ot the same type supported ()
+  if (corr==NULL) {
+    AliError("Zerro pointer - correction");
+    return kFALSE;
+  }  
+  AliError(TString::Format("Correction %s not implementend",IsA()->GetName()).Data());
+  return kFALSE;
+}
+
+
 void AliTPCCorrection::CorrectPoint(Float_t x[], Short_t roc) {
   //
   // Corrects the initial coordinates x (cartesian coordinates)
index db1a245..0067556 100644 (file)
@@ -30,7 +30,7 @@ public:
   AliTPCCorrection();
   AliTPCCorrection(const char *name,const char *title);
   virtual ~AliTPCCorrection();
-  
+  virtual Bool_t AddCorrectionCompact(AliTPCCorrection* corr, Double_t weight);
 
   // functions to correct a space point
           void CorrectPoint (      Float_t x[], Short_t roc);
index 4c120f9..c60d880 100644 (file)
@@ -42,6 +42,23 @@ AliTPCExBTwist::~AliTPCExBTwist() {
   //
 }
 
+Bool_t AliTPCExBTwist::AddCorrectionCompact(AliTPCCorrection* corr, Double_t weight){
+  //
+  // Add correction  and make them compact
+  // Assumptions:
+  //  - origin of distortion/correction are additive
+  //  - only correction ot the same type supported ()
+  if (corr==NULL) {
+    AliError("Zerro pointer - correction");
+    return kFALSE;
+  }  
+  AliTPCExBTwist * corrC = dynamic_cast< AliTPCExBTwist*>(corr);
+  if (corrC == NULL) return kFALSE;
+  fXTwist+=weight*corrC->fXTwist;        // Twist of E to B field in X-Z [rad]
+  fYTwist+=weight*corrC->fYTwist;        // Twist of E to B field in Y-Z [rad]
+  return kTRUE;
+}
+
 
 
 void AliTPCExBTwist::Init() {
index a203cb2..0ec6e46 100644 (file)
@@ -35,6 +35,7 @@
 //   <p>
 // Date: 27/04/2010  <br>
 // Authors: Jim Thomas, Magnus Mager, Stefan Rossegger 
+// Support since 2010: mrain.ivanov@cern.ch
 // End_Html 
 // _________________________________________________________________
 
@@ -44,7 +45,7 @@ class AliTPCExBTwist : public AliTPCCorrection {
 public:
   AliTPCExBTwist();
   virtual ~AliTPCExBTwist();
-
+  virtual Bool_t AddCorrectionCompact(AliTPCCorrection* corr, Double_t weight);
   // initialization and update functions
   virtual void Init();
   virtual void Update(const TTimeStamp &timeStamp);
index 033c80d..9aa6c8a 100644 (file)
@@ -188,6 +188,36 @@ AliTPCFCVoltError3D::~AliTPCFCVoltError3D() {
   }
 }
 
+
+Bool_t AliTPCFCVoltError3D::AddCorrectionCompact(AliTPCCorrection* corr, Double_t weight){
+  //
+  // Add correction  and make them compact
+  // Assumptions:
+  //  - origin of distortion/correction are additive
+  //  - only correction ot the same type supported ()
+  if (corr==NULL) {
+    AliError("Zerro pointer - correction");
+    return kFALSE;
+  }  
+  AliTPCFCVoltError3D * corrC = dynamic_cast<AliTPCFCVoltError3D *>(corr);
+  if (corrC == NULL) return kFALSE;
+  //
+  for (Int_t isec=0; isec<36; isec++){
+    fRodVoltShiftA[isec]+= weight*corrC->fRodVoltShiftA[isec];      // Rod (&strips) shift in Volt (40V~1mm) 
+    fRodVoltShiftC[isec]+=weight*corrC->fRodVoltShiftC[isec];      // Rod (&strips) shift in Volt (40V~1mm) 
+    fCopperRodShiftA[isec]+=weight*corrC->fCopperRodShiftA[isec];    // only Rod shift 
+    fCopperRodShiftC[isec]+=weight*corrC->fCopperRodShiftC[isec];    // only Rod shift 
+  }
+  for (Int_t isec=0; isec<2; isec++){  
+    fRotatedClipVoltA[isec]+= weight*corrC->fRotatedClipVoltA[isec];    // rotated clips at HV rod
+    fRotatedClipVoltC[isec]+= weight*corrC-> fRotatedClipVoltC[isec];    // rotated clips at HV rod
+  }
+
+  return kTRUE;
+}
+
+
+
 void AliTPCFCVoltError3D::Init() {
   //
   // Initialization funtion
index 8575bd6..1ca4a0c 100644 (file)
@@ -16,7 +16,7 @@ class AliTPCFCVoltError3D : public AliTPCCorrection {
 public:
   AliTPCFCVoltError3D();
   virtual ~AliTPCFCVoltError3D();
-
+  virtual Bool_t AddCorrectionCompact(AliTPCCorrection* corr, Double_t weight);
   // initialization and update functions
   virtual void Init();
   virtual void Update(const TTimeStamp &timeStamp);
index 780f372..cba8bdd 100644 (file)
@@ -5,6 +5,8 @@
    
   .L $ALICE_ROOT/TPC/Base/test/UnitTest.C+ 
   UnitTestAliTPCCalPadTree();
+  TestCorrection_AliTPCExBTwistAddCorrectionCompact();
+  TestCorrection_AliTPCFCVoltError3DAddCorrectionCompact()
 
 */
 
 #include "AliTPCCalPad.h"
 #include "AliTPCCalibViewer.h"
 #include "AliTPCcalibDButil.h"
+#include "AliTPCCorrection.h"
+#include "AliTPCComposedCorrection.h"
+#include "AliTPCExBTwist.h"
+#include "AliTPCFCVoltError3D.h"
+//#include "AliTPCROCVoltError3D.h"
+//#include "AliTPCBoundaryVoltError.h"
 
 //
 // PARAMETERS to set from outside:
 //
-TString baseDir="/hera/alice/wiechula/calib/guiTrees";
+TString baseDir="/hera/alice/wiechula/calib/guiTrees";  // TO  FIX specification of inout data
 //
 //
 
@@ -39,14 +47,14 @@ void  UnitTestAliTPCCalPadTree(){
   TObjArray *fArray = new TObjArray(100);
   TTree * treePad=AliTPCcalibDButil::ConnectGainTrees(baseDir);
   for (Int_t i=0; i<5; i+=2){
-    AliTPCCalPad * padLx = AliTPCCalPad::MakePadFromTree(treePad,"lx.fElements","Lx");
-    AliTPCCalPad * padLy = AliTPCCalPad::MakePadFromTree(treePad,"ly.fElements","Ly");
-    AliTPCCalPad * padLLx = AliTPCCalPad::MakePadFromTree(treePad,"lx.fElements","LLx");
-    AliTPCCalPad * padLLy = AliTPCCalPad::MakePadFromTree(treePad,"ly.fElements","LLy");
-    AliTPCCalPad * padMax = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MaxCharge.fElements","QMax");
-    AliTPCCalPad * padMean = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MeanCharge.fElements","QTot");
-    AliTPCCalPad * padMaxL = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MaxCharge.fElements","QMax");
-    AliTPCCalPad * padMeanL = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MeanCharge.fElements","QTot");
+    AliTPCCalPad * padLx = AliTPCCalPad::MakePadFromTree(treePad,"lx.fElements","Lx",kTRUE);
+    AliTPCCalPad * padLy = AliTPCCalPad::MakePadFromTree(treePad,"ly.fElements","Ly",kTRUE);
+    AliTPCCalPad * padLLx = AliTPCCalPad::MakePadFromTree(treePad,"lx.fElements","LLx",kTRUE);
+    AliTPCCalPad * padLLy = AliTPCCalPad::MakePadFromTree(treePad,"ly.fElements","LLy",kTRUE);
+    AliTPCCalPad * padMax = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MaxCharge.fElements","QMax",kTRUE);
+    AliTPCCalPad * padMean = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MeanCharge.fElements","QTot",kTRUE);
+    AliTPCCalPad * padMaxL = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MaxCharge.fElements","QMax",kTRUE);
+    AliTPCCalPad * padMeanL = AliTPCCalPad::MakePadFromTree(treePad,"QA.2010.LHC10d.MeanCharge.fElements","QTot",kTRUE);
     if (i>0) {
       padLx->MedianFilter(i,2*i);
       padLy->MedianFilter(i,2*i);
@@ -92,3 +100,129 @@ void  UnitTestAliTPCCalPadTree(){
   if ((isOutL0+isOutL1)==0) ::Info("UnitTestAliTPCCalPadTree","LTMTest OK");
   if (isOutL0||isOutL1) ::Fatal("UnitTestAliTPCCalPadTree","LTMTest FAILED");
 }
+
+
+Bool_t  TestCorrection_AliTPCExBTwistAddCorrectionCompact(){
+  //
+  // 
+  //#include "AliTPCFCVoltError3D.h"
+  //#include "AliTPCROCVoltError3D.h"
+  //#include "AliTPCBoundaryVoltError.h"
+  // 1.) Test ExB twist AddCorrectionCompact
+  //
+  // 1.) Test ExB twist AddCorrectionCompact
+  //
+  Bool_t isOKTwist[10]={kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
+  AliTPCComposedCorrection *compCorrTwist = new AliTPCComposedCorrection();
+  AliTPCExBTwist  *twistX    = new  AliTPCExBTwist;
+  AliTPCExBTwist  *twistY    = new  AliTPCExBTwist;
+  twistX->SetXTwist(0.001);  // 1 mrad twist in x
+  twistY->SetYTwist(0.001);  // 1 mrad twist in x
+  isOKTwist[0]&=compCorrTwist->AddCorrectionCompact(twistX,1);
+  isOKTwist[0]&=compCorrTwist->AddCorrectionCompact(twistY,1);
+  isOKTwist[0]&=compCorrTwist->AddCorrectionCompact(twistY,-1);
+  isOKTwist[0]&=compCorrTwist->AddCorrectionCompact(twistX,-1);
+  isOKTwist[1]=compCorrTwist->GetCorrections()->GetEntries()==1;
+  AliTPCExBTwist  *twistRes=0;
+  if (isOKTwist[1]==kFALSE){
+    isOKTwist[2]=kFALSE;
+    isOKTwist[3]=kFALSE;
+    isOKTwist[4]=kFALSE;
+  }else{
+    twistRes=  dynamic_cast<AliTPCExBTwist *>(compCorrTwist->GetSubCorrection(0));
+    if (twistRes==NULL){
+      isOKTwist[2]=kFALSE;
+      isOKTwist[3]=kFALSE;
+      isOKTwist[4]=kFALSE;
+    }else{
+      isOKTwist[3] &= (twistRes->GetXTwist()==0);
+      isOKTwist[4] &= (twistRes->GetYTwist()==0);
+    }
+  }
+  Bool_t res=kTRUE;
+  for (Int_t i=0; i<5; i++) res&=isOKTwist[i];
+  {
+    if (isOKTwist[0]==kFALSE){
+      ::Error("TestCorrection_AddCorrectionCompact","AliTPCExBTwist -ADD FAILED");
+    }else{
+      ::Info("TestCorrection_AddCorrectionCompact","AliTPCExBTwist -ADD OK");
+    }
+    if (isOKTwist[1]==kFALSE){
+      ::Error("TestCorrection_AddCorrectionCompact","AliTPCExBTwist - wrong entries  FAILED");
+    }else{
+      ::Info("TestCorrection_AddCorrectionCompact","AliTPCExBTwist - entries  OK");
+    }
+    if (isOKTwist[2]==kFALSE || isOKTwist[3]==kFALSE ||isOKTwist[4]==kFALSE ){
+      ::Error("TestCorrection_AddCorrectionCompact","AliTPCExBTwist - inconsitent entries  FAILED");    
+    }else{
+      ::Info("TestCorrection_AddCorrectionCompact","AliTPCExBTwist - consistent entries  OK");    
+    }
+  }    
+  return res;
+} 
+
+
+Bool_t  TestCorrection_AliTPCFCVoltError3DAddCorrectionCompact(){
+  const Float_t kEpsilon=0.000001;
+  //
+  // 1.) Test ExB  AddCorrectionCompact
+  //
+  Bool_t isOK[10]={kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
+  AliTPCComposedCorrection *compCorrComp = new AliTPCComposedCorrection();
+  AliTPCFCVoltError3D  *corr0    = new  AliTPCFCVoltError3D;
+  AliTPCFCVoltError3D  *corr1    = new  AliTPCFCVoltError3D;
+  for (Int_t isec=0; isec<36; isec++){
+    corr0->SetRodVoltShiftA(isec,TMath::Cos(TMath::Pi()*isec/36),kFALSE);
+    corr0->SetRodVoltShiftC(isec,TMath::Cos(TMath::Pi()*isec/36),kFALSE);
+    corr1->SetRodVoltShiftA(isec,TMath::Sin(TMath::Pi()*isec/36),kFALSE);
+    corr1->SetRodVoltShiftC(isec,TMath::Sin(TMath::Pi()*isec/36),kFALSE);
+    corr1->SetCopperRodShiftA(isec,TMath::Sin(TMath::Pi()*isec/36),kFALSE);
+    corr1->SetCopperRodShiftC(isec,TMath::Sin(TMath::Pi()*isec/36),kFALSE);
+  }
+  //
+  isOK[0]&=compCorrComp->AddCorrectionCompact(corr0,1);
+  isOK[0]&=compCorrComp->AddCorrectionCompact(corr1,1);
+  isOK[0]&=compCorrComp->AddCorrectionCompact(corr1,-1);
+  isOK[0]&=compCorrComp->AddCorrectionCompact(corr0,-1);
+  isOK[1]=compCorrComp->GetCorrections()->GetEntries()==1;
+  AliTPCFCVoltError3D  *corrRes=0;
+  if (isOK[1]==kFALSE){
+    isOK[2]=kFALSE;
+    isOK[3]=kFALSE;
+    isOK[4]=kFALSE;
+  }else{
+    corrRes=  dynamic_cast<AliTPCFCVoltError3D *>(compCorrComp->GetSubCorrection(0));
+    if (corrRes==NULL){
+      isOK[2]=kFALSE;
+      isOK[3]=kFALSE;
+      isOK[4]=kFALSE;
+    }else{
+      for (Int_t isec=0; isec<36; isec++){
+       isOK[3] &=( TMath::Abs(corrRes->GetRodVoltShiftA(isec))<kEpsilon);
+       isOK[4] &=( TMath::Abs(corrRes->GetRodVoltShiftC(isec))<kEpsilon);
+       isOK[5] &=( TMath::Abs(corrRes->GetCopperRodShiftA(isec))<kEpsilon);
+       isOK[6] &=( TMath::Abs(corrRes->GetCopperRodShiftC(isec))<kEpsilon);
+      }
+    }
+  }
+  Bool_t res=kTRUE;
+  for (Int_t i=0; i<5; i++) res&=isOK[i];
+  {
+    if (isOK[0]==kFALSE){
+      ::Error("TestCorrection_AddCorrectionCompact","AliTPCFCVoltError3D -ADD FAILED");
+    }else{
+      ::Info("TestCorrection_AddCorrectionCompact","AliTPCFCVoltError3D -ADD OK");
+    }
+    if (isOK[1]==kFALSE){
+      ::Error("TestCorrection_AddCorrectionCompact","AliTPCFCVoltError3D - wrong entries  FAILED");
+    }else{
+      ::Info("TestCorrection_AddCorrectionCompact","AliTPCFCVoltError3D - entries  OK");
+    }
+    if (isOK[2]==kFALSE || isOK[3]==kFALSE ||isOK[4]==kFALSE ){
+      ::Error("TestCorrection_AddCorrectionCompact","AliTPCFCVoltError3D - inconsitent entries  FAILED");    
+    }else{
+      ::Info("TestCorrection_AddCorrectionCompact","AliTPCFCVoltError3D - consistent entries  OK");    
+    }
+  }    
+  return res;
+}