]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - EMCAL/AliEMCALGeoUtils.cxx
Added some methods for testing purpose
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALGeoUtils.cxx
index da5b2d8e6bbb37180db4ed7a8c49da5deac71797..0ea50c146c4b5c565912a299c4a4aed998db9b98 100644 (file)
@@ -26,7 +26,7 @@
 //        You have to use just the correct name of geometry. If name is empty string the
 //        default name of geometry will be used.
 //         
-//  AliEMCALGeoUtils* geom = new AliEMCALGeoUtils("EMCAL_COMPLETE","EMCAL");
+//  AliEMCALGeoUtils* geom = new AliEMCALGeoUtils("EMCAL_COMPLETEV1","EMCAL");
 //  TGeoManager::Import("geometry.root");
 //
 //  MC:   If you work with MC data you have to get geometry the next way: 
@@ -42,6 +42,7 @@
 #include <TParticle.h>
 #include <TGeoManager.h>
 #include <TGeoMatrix.h>
+#include <TGeoBBox.h>
 #include <TList.h>
 #include <TBrowser.h>
 
@@ -173,7 +174,7 @@ AliEMCALGeoUtils::AliEMCALGeoUtils(const Text_t* name, const Text_t* title)
        
   if (AliDebugLevel()>=2) {
     fEMCGeometry->Print();
-    PrintGeometry();
+    PrintGeometryGeoUtils();
   }
 
   for (Int_t ix = 0; ix < 48; ix++)
@@ -184,7 +185,7 @@ AliEMCALGeoUtils::AliEMCALGeoUtils(const Text_t* name, const Text_t* title)
 
 //____________________________________________________________________________
 AliEMCALGeoUtils & AliEMCALGeoUtils::operator = (const AliEMCALGeoUtils  & /*rvalue*/) { 
-
+  //assing operator
   Fatal("assignment operator", "not implemented") ; 
     return *this ;
 }
@@ -193,12 +194,12 @@ AliEMCALGeoUtils & AliEMCALGeoUtils::operator = (const AliEMCALGeoUtils  & /*rva
 AliEMCALGeoUtils::~AliEMCALGeoUtils(void)
 {
   // dtor
-  for(Int_t smod = 0 ; smod < fEMCGeometry->GetNumberOfSuperModules(); smod++){
-    if(fkSModuleMatrix[smod])
-       delete fkSModuleMatrix[smod] ;
-      fkSModuleMatrix[smod]=0 ;
-  }
-  if(fEMCGeometry){
+  if (fEMCGeometry){ 
+    for(Int_t smod = 0 ; smod < fEMCGeometry->GetNumberOfSuperModules(); smod++){
+      if(fkSModuleMatrix[smod])
+        delete fkSModuleMatrix[smod] ;
+        fkSModuleMatrix[smod]=0 ;
+    }
     delete fEMCGeometry; fEMCGeometry = 0 ;
   }
 }
@@ -300,6 +301,16 @@ void AliEMCALGeoUtils::PrintCellIndexes(Int_t absId, int pri, const char *tit) c
   }
 }
 
+void AliEMCALGeoUtils::PrintLocalTrd1(Int_t pri) const
+{
+  // For comparing with numbers from drawing
+  for(Int_t i=0; i<GetShishKebabTrd1Modules()->GetSize(); i++){
+    printf(" %s | ", GetShishKebabModule(i)->GetName());
+    if(i==0 && pri<1) GetShishKebabModule(i)->PrintShish(1);
+    else     GetShishKebabModule(i)->PrintShish(pri);
+  }
+}
+
 //________________________________________________________________________________________________
 void AliEMCALGeoUtils::EtaPhiFromIndex(Int_t absId,Double_t &eta,Double_t &phi) const
 {
@@ -659,6 +670,7 @@ void AliEMCALGeoUtils::CreateListOfTrd1Modules()
   // Generate the list of Trd1 modules
   // which will make up the EMCAL
   // geometry
+  // key: look to the AliEMCALShishKebabTrd1Module::
 
   AliDebug(2,Form(" AliEMCALGeometry::CreateListOfTrd1Modules() started "));
 
@@ -785,7 +797,7 @@ AliEMCALShishKebabTrd1Module* AliEMCALGeoUtils::GetShishKebabModule(Int_t neta)
 }
 
 //___________________________________________________________________
-void AliEMCALGeoUtils::PrintGeometry()
+void AliEMCALGeoUtils::PrintGeometryGeoUtils()
 {
   //Print information from geometry
   fEMCGeometry->PrintGeometry();
@@ -968,7 +980,7 @@ Bool_t AliEMCALGeoUtils::GetAbsFastORIndexFromTRU(const Int_t iTRU, const Int_t
 {
        //Trigger mapping method, get  FastOr Index from TRU
 
-    if (iTRU > 31 || iTRU < 0 || iADC > 95 || iADC < 0) 
+  if (iTRU > 31 || iTRU < 0 || iADC > 95 || iADC < 0) 
        {
                AliError("TRU out of range!");
                return kFALSE;
@@ -1053,6 +1065,8 @@ Bool_t AliEMCALGeoUtils::GetPositionInSMFromAbsFastORIndex(const Int_t id, Int_t
 //________________________________________________________________________________________________
 Bool_t AliEMCALGeoUtils::GetPositionInEMCALFromAbsFastORIndex(const Int_t id, Int_t& iEta, Int_t& iPhi) const
 {
+  //Trigger mapping method, get position in EMCAL from FastOR index
+
        Int_t iSM=-1;
        
        if (GetPositionInSMFromAbsFastORIndex(id, iSM, iEta, iPhi))
@@ -1086,7 +1100,8 @@ Bool_t AliEMCALGeoUtils::GetAbsFastORIndexFromPositionInTRU(const Int_t iTRU, co
 //________________________________________________________________________________________________
 Bool_t AliEMCALGeoUtils::GetAbsFastORIndexFromPositionInSM(const Int_t  iSM, const Int_t iEta, const Int_t iPhi, Int_t& id) const
 {
-       //
+  //Trigger mapping method, from position in SM Index get FastOR index 
+
        if (iSM < 0 || iSM > 11 || iEta < 0 || iEta > 23 || iPhi < 0 || iPhi > 11) 
        {
                AliError("Out of range!");
@@ -1110,7 +1125,8 @@ Bool_t AliEMCALGeoUtils::GetAbsFastORIndexFromPositionInSM(const Int_t  iSM, con
 //________________________________________________________________________________________________
 Bool_t AliEMCALGeoUtils::GetAbsFastORIndexFromPositionInEMCAL(const Int_t iEta, const Int_t iPhi, Int_t& id) const
 {
-       //
+  //Trigger mapping method, from position in EMCAL Index get FastOR index 
+
        if (iEta < 0 || iEta > 47 || iPhi < 0 || iPhi > 63 ) 
        {
                AliError("Out of range!");
@@ -1131,6 +1147,8 @@ Bool_t AliEMCALGeoUtils::GetAbsFastORIndexFromPositionInEMCAL(const Int_t iEta,
 //________________________________________________________________________________________________
 Bool_t AliEMCALGeoUtils::GetFastORIndexFromCellIndex(const Int_t id, Int_t& idx) const
 {
+  //Trigger mapping method, from cell index get FastOR index 
+
        Int_t iSupMod, nModule, nIphi, nIeta, iphim, ietam;
        
        Bool_t isOK = GetCellIndex( id, iSupMod, nModule, nIphi, nIeta );
@@ -1148,7 +1166,9 @@ Bool_t AliEMCALGeoUtils::GetFastORIndexFromCellIndex(const Int_t id, Int_t& idx)
 //________________________________________________________________________________________________
 Bool_t AliEMCALGeoUtils::GetCellIndexFromFastORIndex(const Int_t id, Int_t idx[4]) const
 {
-       Int_t iSM=-1, iEta=-1, iPhi=-1;
+  //Trigger mapping method, from FASTOR index get cell index 
+
+  Int_t iSM=-1, iEta=-1, iPhi=-1;
        if (GetPositionInSMFromAbsFastORIndex(id, iSM, iEta, iPhi))
        {
                Int_t ix = 2 * iEta;
@@ -1171,6 +1191,8 @@ Bool_t AliEMCALGeoUtils::GetCellIndexFromFastORIndex(const Int_t id, Int_t idx[4
 //________________________________________________________________________________________________
 Bool_t AliEMCALGeoUtils::GetTRUIndexFromSTUIndex(const Int_t id, Int_t& idx) const
 {
+  //Trigger mapping method, from STU index get TRU index 
+
        if (id > 31 || id < 0) 
        {
                AliError(Form("TRU index out of range: %d",id));
@@ -1185,6 +1207,8 @@ Bool_t AliEMCALGeoUtils::GetTRUIndexFromSTUIndex(const Int_t id, Int_t& idx) con
 //________________________________________________________________________________________________
 Int_t AliEMCALGeoUtils::GetTRUIndexFromSTUIndex(const Int_t id) const
 {
+  //Trigger mapping method, from STU index get TRU index 
+
        if (id > 31 || id < 0) 
        {
                AliError(Form("TRU index out of range: %d",id));
@@ -1222,6 +1246,7 @@ void AliEMCALGeoUtils::BuildFastOR2DMap()
 //________________________________________________________________________________________________
 Bool_t AliEMCALGeoUtils::GetFastORIndexFromL0Index(const Int_t iTRU, const Int_t id, Int_t idx[], const Int_t size) const
 {
+  //Trigger mapping method, from L0 index get FastOR index 
        if (size <= 0 ||size > 4)
        {
                AliError("Size not supported!");
@@ -1233,7 +1258,7 @@ Bool_t AliEMCALGeoUtils::GetFastORIndexFromL0Index(const Int_t iTRU, const Int_t
        switch (size)
        {
                case 1: // Cosmic trigger
-                       if (!GetAbsFastORIndexFromTRU(iTRU, id, idx[0])) return kFALSE;
+                       if (!GetAbsFastORIndexFromTRU(iTRU, id, idx[1])) return kFALSE;
                        break;
                case 4: // 4 x 4
                        for (Int_t k = 0; k < 4; k++)
@@ -1270,13 +1295,14 @@ const TGeoHMatrix * AliEMCALGeoUtils::GetMatrixForSuperModule(Int_t smod) const
        //    TGeoHMatrix* m = gGeoManager->GetCurrentMatrix();
        
        if(gGeoManager){
-               char path[255] ;
-               sprintf(path,"/ALIC_1/XEN1_1/SMOD_%d",smod+1) ;
+    const Int_t buffersize = 255;
+               char path[buffersize] ;
+               snprintf(path,buffersize,"/ALIC_1/XEN1_1/SMOD_%d",smod+1) ;
                //TString volpath = "ALIC_1/XEN1_1/SMOD_";
            //volpath += smod+1;
 
                if(fKey110DEG && smod >= 10){
-                         sprintf(path,"/ALIC_1/XEN1_1/SM10_%d",smod-10+1) ;
+                         snprintf(path,buffersize,"/ALIC_1/XEN1_1/SM10_%d",smod-10+1) ;
                        //volpath = "ALIC_1/XEN1_1/SM10_";
                        //volpath += smod-10+1;
                }
@@ -1318,3 +1344,140 @@ void AliEMCALGeoUtils::GetModulePhiEtaIndexInSModuleFromTRUIndex(Int_t itru, Int
   //printf(" GetModulePhiEtaIndexInSModuleFromTRUIndex : itru %2i iphitru %2i ietatru %2i iphiSM %2i ietaSM %2i \n", 
   // itru, iphitru, ietatru, iphiSM, ietaSM);
 }
+
+//__________________________________________________________________________________________________________________
+void AliEMCALGeoUtils::RecalculateTowerPosition(Float_t drow, Float_t dcol, const Int_t sm, const Float_t depth,
+                                                const Float_t misaligTransShifts[15], const Float_t misaligRotShifts[15], Float_t global[3]) const
+{ //Transform clusters cell position into global with alternative method, taking into account the depth calculation.
+  //Input are: the tower indeces, 
+  //           supermodule, 
+  //           particle type (photon 0, electron 1, hadron 2 )
+  //           misalignment shifts to global position in case of need.
+  // Federico.Ronchetti@cern.ch
+  
+  // To use in a print later
+  //Int_t iphi = drow;
+  //Int_t ieta = dcol;
+  
+  if(gGeoManager){
+    //Recover some stuff
+    
+    gGeoManager->cd("ALIC_1/XEN1_1");
+    TGeoNode        *geoXEn1 = gGeoManager->GetCurrentNode();
+    TGeoNodeMatrix  *geoSM[4];        
+    TGeoVolume      *geoSMVol[4];     
+    TGeoShape       *geoSMShape[4];    
+    TGeoBBox        *geoBox[4];        
+    TGeoMatrix      *geoSMMatrix[4];       
+    
+    for(int iSM = 0; iSM < 4; iSM++) {  
+      geoSM[iSM]       = dynamic_cast<TGeoNodeMatrix *>(geoXEn1->GetDaughter(iSM));
+      geoSMVol[iSM]    = geoSM[iSM]->GetVolume(); 
+      geoSMShape[iSM]  = geoSMVol[iSM]->GetShape();
+      geoBox[iSM]      = dynamic_cast<TGeoBBox *>(geoSMShape[iSM]);
+      geoSMMatrix[iSM] = geoSM[iSM]->GetMatrix();
+    }
+    
+    if(sm % 2 == 0) {
+      dcol = 47. - dcol;
+      drow = 23. - drow;
+    }
+    
+    Int_t i      = 0; // one always needs "i"
+    Int_t istrip = 0;
+    Float_t z0   = 0;
+    Float_t zb   = 0;
+    Float_t z_is = 0;
+    
+    Float_t x,y,z; // return variables in terry's RF
+    
+    //***********************************************************
+    //Do not like this: too many hardcoded values, is it not already stored somewhere else?
+    //                : need more comments in the code 
+    //***********************************************************
+    
+    Float_t dz = 6.0;   // base cell width in eta
+    Float_t dx = 6.004; // base cell width in phi
+    
+    
+    //Float_t L = 26.04; // active tower length for hadron (lead+scint+paper)
+    // we use the geant numbers 13.87*2=27.74
+    Float_t teta1 = 0.;
+    
+    i = sm;
+    
+    if (dcol >= 47.5 || dcol<-0.5) {
+      exit(0);
+    }
+    
+    if (drow >= 23.5 || drow<-0.5) {
+      exit(0);
+    }
+    if (sm > 13 || sm <0) {
+      exit(0);
+    }
+        
+    istrip = int ((dcol+0.5)/2);
+    
+    // tapering angle
+    teta1 = TMath::DegToRad() * istrip * 1.5;
+    
+    // calculation of module corner along z 
+    // as a function of strip
+    
+    for (int is=0; is<= istrip; is++) {
+      
+      teta1 = TMath::DegToRad() * (is*1.5 + 0.75);
+      if(is==0)
+        z_is = z_is + 2*dz*TMath::Cos(teta1);
+      else
+        z_is = z_is + 2*dz*TMath::Cos(teta1) + 2*dz*TMath::Sin(teta1)*TMath::Tan(teta1-0.75*TMath::DegToRad());
+      
+    }
+    
+    z0 = dz*(dcol-2*istrip+0.5);
+    zb = (2*dz-z0-depth*TMath::Tan(teta1));
+    
+    z = z_is - zb*TMath::Cos(teta1);
+    y = depth/TMath::Cos(teta1) + zb*TMath::Sin(teta1);
+    
+    x = (drow + 0.5)*dx;
+    
+    // moving the origin from terry's RF
+    // to the GEANT one
+    
+    double xx =  y - geoBox[i]->GetDX();
+    double yy = -x + geoBox[i]->GetDY(); 
+    double zz =  z - geoBox[i]->GetDZ(); 
+    const double localIn[3] = {xx, yy, zz};
+    double dglobal[3];
+    //geoSMMatrix[i]->Print();
+    //printf("TFF Local    (row = %d, col = %d, x = %3.2f,  y = %3.2f, z = %3.2f)\n", iphi, ieta, localIn[0], localIn[1], localIn[2]);
+    geoSMMatrix[i]->LocalToMaster(localIn, dglobal);
+    //printf("TFF Global   (row = %2.0f, col = %2.0f, x = %3.2f,  y = %3.2f, z = %3.2f)\n", drow, dcol, dglobal[0], dglobal[1], dglobal[2]);
+    
+    //apply global shifts
+    if(sm == 2 || sm == 3) {//sector 1
+      global[0] = dglobal[0] + misaligTransShifts[3] + misaligRotShifts[3]*TMath::Sin(TMath::DegToRad()*20) ; 
+      global[1] = dglobal[1] + misaligTransShifts[4] + misaligRotShifts[4]*TMath::Cos(TMath::DegToRad()*20) ; 
+      global[2] = dglobal[2] + misaligTransShifts[5];
+    }
+    else if(sm == 0 || sm == 1){//sector 0
+      global[0] = dglobal[0] + misaligTransShifts[0]; 
+      global[1] = dglobal[1] + misaligTransShifts[1]; 
+      global[2] = dglobal[2] + misaligTransShifts[2];
+    }
+    else {
+      AliInfo("Careful, correction not implemented yet!");
+      global[0] = dglobal[0] ;
+      global[1] = dglobal[1] ;
+      global[2] = dglobal[2] ;
+    }
+    
+    
+  }
+  else{
+    AliFatal("Geometry boxes information, check that geometry.root is loaded\n");
+  }
+  
+}