Implementing the wrapper volumes (Mario)

[u/mrichter/AliRoot.git] / ITS / UPGRADE / AliITSUGeomTGeo.cxx

diff --git a/ITS/UPGRADE/AliITSUGeomTGeo.cxx b/ITS/UPGRADE/AliITSUGeomTGeo.cxx
index 9eddadbca14abf6e9dbce33b47553c0ea77a298a..f87ce3ebc5d07c0b1cc3fa794635dd9e2ac96e1a 100644 (file)
--- a/ITS/UPGRADE/AliITSUGeomTGeo.cxx
+++ b/ITS/UPGRADE/AliITSUGeomTGeo.cxx
@@ -32,10 +32,15 @@
 #include <TGeoManager.h>
 #include <TGeoPhysicalNode.h>
 #include <TDatime.h>
+#include <TMath.h>
+#include <TSystem.h>
 
 #include "AliITSUGeomTGeo.h"
 #include "AliLog.h"
 #include "AliAlignObj.h"
+#include "AliITSsegmentation.h"
+#include "AliITSUSegmentationPix.h"
+using namespace TMath;
 
 ClassImp(AliITSUGeomTGeo)
 
@@ -50,7 +55,7 @@ const char* AliITSUGeomTGeo::fgkITSDetTypeName[AliITSUGeomTGeo::kNDetTypes] = {"
 TString     AliITSUGeomTGeo::fgITSsegmFileName = "itsSegmentations.root";
 
 //______________________________________________________________________
-AliITSUGeomTGeo::AliITSUGeomTGeo(Bool_t build)
+AliITSUGeomTGeo::AliITSUGeomTGeo(Bool_t build, Bool_t loadSegm)
   :fVersion(kITSVNA)
   ,fNLayers(0)
   ,fNModules(0)
@@ -60,9 +65,10 @@ AliITSUGeomTGeo::AliITSUGeomTGeo(Bool_t build)
   ,fLastModIndex(0)
   ,fMatSens(0)
   ,fMatT2L(0)
+  ,fSegm(0)
 {
   // default c-tor
-  if (build) BuildITS();
+  if (build) BuildITS(loadSegm);
 }
 
 //______________________________________________________________________
@@ -77,6 +83,7 @@ AliITSUGeomTGeo::AliITSUGeomTGeo(const AliITSUGeomTGeo &src)
   ,fLastModIndex(0)
   ,fMatSens(0)
   ,fMatT2L(0)
+  ,fSegm(0)
 {
   // copy c-tor
   if (fNLayers) {
@@ -103,7 +110,17 @@ AliITSUGeomTGeo::AliITSUGeomTGeo(const AliITSUGeomTGeo &src)
       fMatT2L->SetOwner(kTRUE);
       for (int i=0;i<fNModules;i++) {
        const TGeoHMatrix* mat =(TGeoHMatrix*) src.fMatT2L->At(i);
-       fMatSens->AddAt(new TGeoHMatrix(*mat),i);
+       fMatT2L->AddAt(new TGeoHMatrix(*mat),i);
+      }
+    }
+    if (src.fSegm) {
+      int sz = src.fSegm->GetEntriesFast();
+      fSegm = new TObjArray(sz);
+      fSegm->SetOwner(kTRUE);
+      for (int i=0;i<sz;i++) {
+       AliITSsegmentation* sg = (AliITSsegmentation*)src.fSegm->UncheckedAt(i);
+       if (!sg) continue;
+       fSegm->AddAt(sg->Clone(),i);
       }
     }
   }
@@ -119,6 +136,7 @@ AliITSUGeomTGeo::~AliITSUGeomTGeo()
   delete[] fLastModIndex;
   delete fMatT2L;
   delete fMatSens;
+  delete fSegm;
 }
 
 
@@ -153,6 +171,16 @@ AliITSUGeomTGeo& AliITSUGeomTGeo::operator=(const AliITSUGeomTGeo &src)
        fMatT2L->AddAt(new TGeoHMatrix(*mat),i);
       }
     }
+    if (src.fSegm) {
+      int sz = src.fSegm->GetEntriesFast();
+      fSegm = new TObjArray(sz);
+      fSegm->SetOwner(kTRUE);
+      for (int i=0;i<sz;i++) {
+       AliITSsegmentation* sg = (AliITSsegmentation*)src.fSegm->UncheckedAt(i);
+       if (!sg) continue;
+       fSegm->AddAt(sg->Clone(),i);
+      }
+    }
     //
     if (fNLayers) {
       fNLadders   = new Int_t[fNLayers];
@@ -318,6 +346,7 @@ TGeoHMatrix* AliITSUGeomTGeo::GetMatrix(Int_t index)  const
 {
   // Get the transformation matrix for a given module 'index'
   // by quering the TGeoManager
+  static TGeoHMatrix matTmp;
   TGeoPNEntry *pne = GetPNEntry(index);
   if (!pne) return NULL;
 
@@ -331,9 +360,9 @@ TGeoHMatrix* AliITSUGeomTGeo::GetMatrix(Int_t index)  const
     AliError(Form("Volume path %s not valid!",path));
     return NULL;
   }
-  TGeoHMatrix *mat = gGeoManager->GetCurrentMatrix();
+  matTmp = *gGeoManager->GetCurrentMatrix();
   gGeoManager->PopPath();
-  return mat;
+  return &matTmp;
 }
 
 //______________________________________________________________________
@@ -446,6 +475,7 @@ TGeoHMatrix* AliITSUGeomTGeo::ExtractMatrixSens(Int_t index) const
 {
   // Get the transformation matrix of the SENSOR (not ncessary the same as the module) 
   // for a given module 'index' by quering the TGeoManager
+  static TGeoHMatrix matTmp;
   const TString kPathBase = Form("/ALIC_1/%s_2/",AliITSUGeomTGeo::GetITSVolPattern());
   const TString kNames = Form("%%s%s%%d_1/%s%%d_%%d/%s%%d_%%d/%s%%d_%%d"
                              ,AliITSUGeomTGeo::GetITSLayerPattern()
@@ -463,13 +493,13 @@ TGeoHMatrix* AliITSUGeomTGeo::ExtractMatrixSens(Int_t index) const
     AliError(Form("Error in cd-ing to %s",path.Data()));
     return 0;
   } // end if !gGeoManager
-  TGeoHMatrix* mat = gGeoManager->GetCurrentMatrix();
+  matTmp = *gGeoManager->GetCurrentMatrix(); // matrix may change after cd
   //RSS
   //  printf("%d/%d/%d %s\n",lay,ladd,detInLad,path.Data());
   //  mat->Print();
   // Retstore the modeler state.
   gGeoManager->PopPath();
-  return mat;
+  return &matTmp;
 }
 
 
@@ -498,7 +528,7 @@ TGeoPNEntry* AliITSUGeomTGeo::GetPNEntry(Int_t index) const
 }
 
 //______________________________________________________________________
-void AliITSUGeomTGeo::BuildITS()
+void AliITSUGeomTGeo::BuildITS(Bool_t loadSegm)
 {
   // exract upg ITS parameters from TGeo
   if (fVersion!=kITSVNA) {AliWarning("Already built"); return; // already initialized}
@@ -523,6 +553,11 @@ void AliITSUGeomTGeo::BuildITS()
   FetchMatrices();
   fVersion = kITSVUpg;
   //
+  if (loadSegm) {  // fetch segmentations
+    fSegm = new TObjArray();
+    AliITSUSegmentationPix::LoadSegmentations(fSegm,GetITSsegmentationFileName());
+  }
+  //
 }
 
 //______________________________________________________________________
@@ -634,10 +669,64 @@ void AliITSUGeomTGeo::FetchMatrices()
   if (!gGeoManager) AliFatal("Geometry is not loaded");
   fMatSens = new TObjArray(fNModules);
   fMatSens->SetOwner(kTRUE);
+  for (int i=0;i<fNModules;i++) fMatSens->AddAt(new TGeoHMatrix(*ExtractMatrixSens(i)),i);
+  CreateT2LMatrices();
+}
+
+//______________________________________________________________________
+void AliITSUGeomTGeo::CreateT2LMatrices()
+{
+  // create tracking to local (Sensor!) matrices
   fMatT2L  = new TObjArray(fNModules);  
   fMatT2L->SetOwner(kTRUE);
-  for (int i=0;i<fNModules;i++) {
-    fMatSens->AddAt(new TGeoHMatrix(*ExtractMatrixSens(i)),i);
-    fMatT2L->AddAt(new TGeoHMatrix(*ExtractMatrixT2L(i)),i);
+  TGeoHMatrix matLtoT;
+  double loc[3]={0,0,0},glo[3];
+  const double *rotm;
+  for (int isn=0;isn<fNModules;isn++) {
+    const TGeoHMatrix* matSens = GetMatrixSens(isn);
+    if (!matSens) {AliFatal(Form("Failed to get matrix for sensor %d",isn)); return;}
+    matSens->LocalToMaster(loc,glo);
+    rotm = matSens->GetRotationMatrix();
+    Double_t al = -ATan2(rotm[1],rotm[0]);
+    double sn=Sin(al), cs=Cos(al), r=glo[0]*sn-glo[1]*cs, x=r*sn, y=-r*cs; // sensor plane PCA to origin
+    TGeoHMatrix* t2l = new TGeoHMatrix();
+    t2l->RotateZ(ATan2(y,x)*RadToDeg()); // rotate in direction of normal to the sensor plane
+    t2l->SetDx(x);
+    t2l->SetDy(y);
+    t2l->MultiplyLeft(&matSens->Inverse());
+    fMatT2L->AddAt(t2l,isn);
+    /*
+    const double *gtrans = matSens->GetTranslation();
+    memcpy(&rotMatrix[0], matSens->GetRotationMatrix(), 9*sizeof(Double_t));
+    Double_t al = -ATan2(rotMatrix[1],rotMatrix[0]);
+    Double_t rSens = Sqrt(gtrans[0]*gtrans[0] + gtrans[1]*gtrans[1]);
+    Double_t tanAl = ATan2(gtrans[1],gtrans[0]) - Pi()/2; //angle of tangent
+    Double_t alTr = tanAl - al;
+    //
+    // The X axis of tracking frame must always look outward
+    loc[1] = rSens/2;
+    matSens->LocalToMaster(loc,glo);
+    double rPos = Sqrt(glo[0]*glo[0] + glo[1]*glo[1]);
+    Bool_t rotOutward = rPos>rSens ? kFALSE : kTRUE;
+    //
+    // Transformation matrix
+    matLtoT.Clear();
+    matLtoT.SetDx(-rSens*Sin(alTr)); // translation
+    matLtoT.SetDy(0.);
+    matLtoT.SetDz(gtrans[2]);
+    // Rotation matrix
+    rotMatrix[0]= 0;  rotMatrix[1]= 1;  rotMatrix[2]= 0; // + rotation
+    rotMatrix[3]=-1;  rotMatrix[4]= 0;  rotMatrix[5]= 0;
+    rotMatrix[6]= 0;  rotMatrix[7]= 0;  rotMatrix[8]= 1;
+    //
+    TGeoRotation rot;
+    rot.SetMatrix(rotMatrix);
+    matLtoT.MultiplyLeft(&rot);
+    if (rotOutward) matLtoT.RotateZ(180.);
+    // Inverse transformation Matrix
+    fMatT2L->AddAt(new TGeoHMatrix(matLtoT.Inverse()),isn);
+    */
   }
+  //
 }
+