]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - EMCAL/AliEMCALGeometry.cxx
Coding conventions...
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALGeometry.cxx
index b2abc8443dad32a9045ca08059a6b4d535ce6879..b6cba1b5b5ea3f848ca40190c8ff696f54fd9c35 100644 (file)
@@ -36,7 +36,7 @@
 
 #include <assert.h>
 
-// --- AliRoot header files ---
+// --- Root header files ---
 #include <Riostream.h>
 #include <TBrowser.h>
 #include <TClonesArray.h>
@@ -47,6 +47,7 @@
 #include <TMatrixD.h>
 #include <TObjArray.h>
 #include <TObjString.h>
+#include <TVector2.h>
 #include <TVector3.h>
 
 // -- ALICE Headers.
@@ -237,7 +238,6 @@ void AliEMCALGeometry::Init(void){
   fArm1EtaMax     = +0.7;      // pseudorapidity, Ending EMCAL Eta position
   fIPDistance     = 454.0;      // cm, Radial distance to inner surface of EMCAL
   fPhiGapForSM    = 0.;         // cm, only for final TRD1 geometry
-  for(int i=0; i<12; i++) fMatrixOfSM[i] = 0;
 
   // geometry
   if(fGeoName.Contains("SHISH")){ // Only shahslyk now
@@ -334,7 +334,7 @@ void AliEMCALGeometry::Init(void){
       fECScintThick  = fECPbRadThickness = 0.5;
     }
     if(fGeoName.Contains("WSUC")){ // 18-may-05 - about common structure
-      fShellThickness = 30.; // should be change 
+      fShellThickness = 30.;       // should be change 
       fNPhi = fNZ = 4; 
     }
 
@@ -409,38 +409,9 @@ void AliEMCALGeometry::Init(void){
   fNCellsInTRUEta = 16 ;
   fNCellsInTRUPhi = 24 ;
 
-      // Define TGeoMatrix of SM - Jan 19, 2007 (just fro TRD1)
-  if(fGeoName.Contains("TRD1")) { // copy code from  AliEMCALv0::CreateSmod()
-    int nphism  = GetNumberOfSuperModules()/2;
-    double dphi = (GetArm1PhiMax() - GetArm1PhiMin())/nphism;
-    double rpos = (GetEnvelop(0) + GetEnvelop(1))/2.;
-    double phi, phiRad, xpos, ypos, zpos;
-    for(int i=0; i<nphism; i++){
-       phi    = GetArm1PhiMin() + dphi*(2*i+1)/2.; // phi= 90, 110, 130, 150, 170, 190
-       phiRad = phi*TMath::Pi()/180.;
-       xpos = rpos * TMath::Cos(phiRad);
-       ypos = rpos * TMath::Sin(phiRad);
-       zpos = fParSM[2];
-       if(i==5) {
-         xpos += (fParSM[1]/2. * TMath::Sin(phiRad)); 
-         ypos -= (fParSM[1]/2. * TMath::Cos(phiRad));
-       }
-       // pozitive z
-       int ind = 2*i;
-       TGeoRotation *geoRot0 = new TGeoRotation("geoRot0", 90.0, phi, 90.0, 90.0+phi, 0.0, 0.0);
-       fMatrixOfSM[ind] = new TGeoCombiTrans(Form("EmcalSM%2.2i",ind),
-                                                xpos,ypos, zpos, geoRot0);
-       // negaive z
-       ind++;
-       double phiy = 90. + phi + 180.;
-       if(phiy>=360.) phiy -= 360.;
-       TGeoRotation *geoRot1 = new TGeoRotation("geoRot1", 90.0, phi, 90.0, phiy, 180.0, 0.0);
-       fMatrixOfSM[ind] = new TGeoCombiTrans(Form("EmcalSM%2.2i",ind),
-                                          xpos,ypos,-zpos, geoRot1);
-    } // for
-  }
+  if(fGeoName.Contains("WSUC")) fNumberOfSuperModules = 1; // Jul 12, 2007
+
   fgInit = kTRUE; 
-  AliInfo(" is ended");  
 }
 
 void AliEMCALGeometry::PrintGeometry()
@@ -519,14 +490,6 @@ void AliEMCALGeometry::PrintGeometry()
         printf("\n");
 
       }
-      printf(" Matrix transformation\n");
-      for(Int_t i=0; i<12; i++) {
-        TGeoMatrix *m = fMatrixOfSM[i];
-        if(m==0) continue;
-        const double *xyz = m->GetTranslation();
-        printf(" %2.2i %s %s x %7.2f y %7.2f z %7.2f\n", 
-              i, m->GetName(), m->ClassName(), xyz[0],xyz[1],xyz[2]); 
-      }
 
       printf("\n Cells grid in phi directions : size %i\n", fCentersOfCellsPhiDir.GetSize());
       for(Int_t i=0; i<fCentersOfCellsPhiDir.GetSize(); i++) {
@@ -536,7 +499,6 @@ void AliEMCALGeometry::PrintGeometry()
       }
     }
   }
-  cout<<endl;
 }
 
 void AliEMCALGeometry::PrintCellIndexes(Int_t absId, int pri, char *tit)
@@ -710,6 +672,7 @@ Bool_t AliEMCALGeometry::IsInEMCAL(Double_t x, Double_t y, Double_t z) const {
        return 0;
  
      Double_t phi = TMath::ATan2(y,x) * 180./TMath::Pi();
+     if (phi < 0) phi += 360;  // phi should go from 0 to 360 in this case
      if (phi > fArm1PhiMin && phi < fArm1PhiMax)
        return 1;
   }
@@ -890,7 +853,7 @@ Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t &xr, Double_t
 
   // Shift index taking into account the difference between standard SM 
   // and SM of half size in phi direction
-  const Int_t phiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
+  const Int_t kphiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
   static Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta;
   if(!CheckAbsCellId(absId)) return kFALSE;
 
@@ -903,7 +866,7 @@ Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t &xr, Double_t
   if(nSupMod<10) {
     yr = fCentersOfCellsPhiDir.At(iphi);
   } else {
-    yr = fCentersOfCellsPhiDir.At(iphi + phiIndexShift);
+    yr = fCentersOfCellsPhiDir.At(iphi + kphiIndexShift);
   }
   AliDebug(1,Form("absId %i nSupMod %i iphi %i ieta %i xr %f yr %f zr %f ",absId,nSupMod,iphi,ieta,xr,yr,zr));
 
@@ -933,6 +896,110 @@ Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, TVector3 &vloc) const
   // Alice numbering scheme - Jun 03, 2006
 }
 
+Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t distEff, Double_t &xr, Double_t &yr, Double_t &zr) const
+{
+  // Jul 30, 2007 - taking into account position of shower max
+  // Look to see what the relative
+  // position inside a given cell is
+  // for a recpoint.
+  // In:
+  // absId   - cell is as in Geant,     0<= absId   < fNCells;
+  // e       - cluster energy
+  // OUT:
+  // xr,yr,zr - x,y,z coordinates of cell with absId inside SM 
+
+  // Shift index taking into account the difference between standard SM 
+  // and SM of half size in phi direction
+  const  Int_t kphiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
+  static Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta;
+  static Int_t iphim, ietam;
+  static AliEMCALShishKebabTrd1Module *mod = 0;
+  static TVector2 v;
+  if(!CheckAbsCellId(absId)) return kFALSE;
+
+  GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
+  GetModulePhiEtaIndexInSModule(nSupMod, nModule, iphim, ietam);
+  GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta); 
+  mod = GetShishKebabModule(ietam);
+  mod->GetPositionAtCenterCellLine(nIeta, distEff, v); 
+  xr = v.Y() - fParSM[0];
+  zr = v.X() - fParSM[2];
+
+  if(nSupMod<10) {
+    yr = fCentersOfCellsPhiDir.At(iphi);
+  } else {
+    yr = fCentersOfCellsPhiDir.At(iphi + kphiIndexShift);
+  }
+  AliDebug(1,Form("absId %i nSupMod %i iphi %i ieta %i xr %f yr %f zr %f ",absId,nSupMod,iphi,ieta,xr,yr,zr));
+
+  return kTRUE;
+}
+
+Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Int_t maxAbsId, Double_t distEff, Double_t &xr, Double_t &yr, Double_t &zr) const
+{
+  // Jul 31, 2007 - taking into account position of shower max and apply coor2.
+  // Look to see what the relative
+  // position inside a given cell is
+  // for a recpoint.
+  // In:
+  // absId     - cell is as in Geant,     0<= absId   < fNCells;
+  // maxAbsId  - abs id of cell with highest energy
+  // e         - cluster energy
+  // OUT:
+  // xr,yr,zr - x,y,z coordinates of cell with absId inside SM 
+
+  // Shift index taking into account the difference between standard SM 
+  // and SM of half size in phi direction
+  const  Int_t kphiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
+  static Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta;
+  static Int_t iphim, ietam;
+  static AliEMCALShishKebabTrd1Module *mod = 0;
+  static TVector2 v;
+
+  static Int_t nSupModM, nModuleM, nIphiM, nIetaM, iphiM, ietaM;
+  static Int_t iphimM, ietamM, maxAbsIdCopy=-1;
+  static AliEMCALShishKebabTrd1Module *modM = 0;
+  static Double_t distCorr;
+
+  if(!CheckAbsCellId(absId)) return kFALSE;
+
+  GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
+  GetModulePhiEtaIndexInSModule(nSupMod, nModule, iphim, ietam);
+  GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta); 
+  mod = GetShishKebabModule(ietam);
+
+  if(absId != maxAbsId) {
+    distCorr = 0.;
+    if(maxAbsIdCopy != maxAbsId) {
+      GetCellIndex(maxAbsId, nSupModM, nModuleM, nIphiM, nIetaM);
+      GetModulePhiEtaIndexInSModule(nSupModM, nModuleM, iphimM, ietamM);
+      GetCellPhiEtaIndexInSModule(nSupModM,nModuleM,nIphiM,nIetaM, iphiM, ietaM); 
+      modM = GetShishKebabModule(ietamM); // do I need this ?
+      maxAbsIdCopy = maxAbsId;
+    }
+    if(ietamM !=0) {
+      distCorr = GetEtaModuleSize()*(ietam-ietamM)/TMath::Tan(modM->GetTheta()); // Stay here
+      //printf(" distCorr %f | dist %f | ietam %i -> etamM %i\n", distCorr, dist, ietam, ietamM);  
+    }
+    // distEff += distCorr;
+  }
+  // Bad resolution in this case, strong bias vs phi
+  // distEff = 0.0; 
+  mod->GetPositionAtCenterCellLine(nIeta, distEff, v); // Stay here
+  xr = v.Y() - fParSM[0];
+  zr = v.X() - fParSM[2];
+
+  if(nSupMod<10) {
+    yr = fCentersOfCellsPhiDir.At(iphi);
+  } else {
+    yr = fCentersOfCellsPhiDir.At(iphi + kphiIndexShift);
+  }
+  AliDebug(1,Form("absId %i nSupMod %i iphi %i ieta %i xr %f yr %f zr %f ",absId,nSupMod,iphi,ieta,xr,yr,zr));
+
+  return kTRUE;
+}
+
 void AliEMCALGeometry::CreateListOfTrd1Modules()
 {
   // Generate the list of Trd1 modules
@@ -981,7 +1048,7 @@ void AliEMCALGeometry::CreateListOfTrd1Modules()
   fCentersOfCellsPhiDir.Set(fNPhi*fNPHIdiv);
   fPhiCentersOfCells.Set(fNPhi*fNPHIdiv);
 
-  Double_t R0 = GetIPDistance() + GetLongModuleSize()/2.;
+  Double_t r0 = GetIPDistance() + GetLongModuleSize()/2.;
   for(Int_t it=0; it<fNPhi; it++) { // cycle on modules
     ytCenterModule = -fParSM[1] + fPhiModuleSize*(2*it+1)/2;  // center of module
     for(Int_t ic=0; ic<fNPHIdiv; ic++) { // cycle on cells in module
@@ -996,7 +1063,7 @@ void AliEMCALGeometry::CreateListOfTrd1Modules()
       // Define grid on phi direction
       // Grid is not the same for different eta bin;
       // Effect is small but is still here
-      phi = TMath::ATan2(ytCenterCell, R0);
+      phi = TMath::ATan2(ytCenterCell, r0);
       fPhiCentersOfCells.AddAt(phi, ind);
 
       AliDebug(2,Form(" ind %2.2i : y %8.3f ", ind, fCentersOfCellsPhiDir.At(ind))); 
@@ -1045,57 +1112,32 @@ void AliEMCALGeometry::CreateListOfTrd1Modules()
 
 }
 
-void  AliEMCALGeometry::GetTransformationForSM()
-{
-  //Uses the geometry manager to
-  //load the transformation matrix
-  //for the supermodules
-  // Unused after 19 Jan, 2007 - keep for compatibility; 
-
-  return;
-  static Bool_t transInit=kFALSE;
-  if(transInit) return;
-
-  int i=0;
-  if(gGeoManager == 0) {
-    Info("CreateTransformationForSM() "," Load geometry : TGeoManager::Import()");
-    assert(0);
-  }
-  TGeoNode *tn = gGeoManager->GetTopNode();
-  TGeoNode *node=0, *xen1 = 0;
-  for(i=0; i<tn->GetNdaughters(); i++) {
-    node = tn->GetDaughter(i);
-    TString ns(node->GetName());
-    if(ns.Contains(GetNameOfEMCALEnvelope())) {
-      xen1 = node;
-      break;
-    }
-  }
-  if(!xen1) {
-    Info("CreateTransformationForSM() "," geometry has not EMCAL envelope with name %s", 
-    GetNameOfEMCALEnvelope());
-    assert(0);
-  }
-  printf(" i %i : EMCAL Envelope is %s : #SM %i \n", i, xen1->GetName(), xen1->GetNdaughters());
-  for(i=0; i<xen1->GetNdaughters(); i++) {
-    TGeoNodeMatrix *sm = (TGeoNodeMatrix*)xen1->GetDaughter(i);
-    fMatrixOfSM[i] = sm->GetMatrix();
-    //Compiler doesn't like this syntax...
-    //    printf(" %i : matrix %x \n", i, fMatrixOfSM[i]);
-  }
-  transInit = kTRUE;
-}
-
 void AliEMCALGeometry::GetGlobal(const Double_t *loc, Double_t *glob, int ind) const
 {
   // Figure out the global numbering
   // of a given supermodule from the
-  // local numbering
-  // Alice numbering - Jun 03,2006
-  //  if(fMatrixOfSM[0] == 0) GetTransformationForSM();
+  // local numbering and the transformation
+  // matrix stored by the geometry manager (allows for misaligned
+  // geometry)
 
   if(ind>=0 && ind < GetNumberOfSuperModules()) {
-    fMatrixOfSM[ind]->LocalToMaster(loc, glob);
+    TString volpath = "ALIC_1/XEN1_1/SMOD_";
+    volpath += ind+1;
+
+    if(GetKey110DEG() && ind>=10) {
+      volpath = "ALIC_1/XEN1_1/SM10_";
+      volpath += ind-10+1;
+    }
+
+    if(!gGeoManager->cd(volpath.Data()))
+      AliFatal(Form("AliEMCALGeometry::GeoManager cannot find path %s!",volpath.Data()));
+
+    TGeoHMatrix* m = gGeoManager->GetCurrentMatrix();
+    if(m) {
+      m->LocalToMaster(loc, glob);
+    } else {
+      AliFatal("Geo matrixes are not loaded \n") ;
+    }
   }
 }
 
@@ -1120,10 +1162,27 @@ void AliEMCALGeometry::GetGlobal(Int_t absId , double glob[3]) const
   glob[0]=glob[1]=glob[2]=0.0; // bad case
   if(RelPosCellInSModule(absId, loc)) {
     GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
-    fMatrixOfSM[nSupMod]->LocalToMaster(loc, glob);
+
+    TString volpath = "ALIC_1/XEN1_1/SMOD_";
+    volpath += (nSupMod+1);
+
+    if(GetKey110DEG() && nSupMod>=10) {
+      volpath = "ALIC_1/XEN1_1/SM10_";
+      volpath += (nSupMod-10+1);
+    }
+    if(!gGeoManager->cd(volpath.Data()))
+      AliFatal(Form("GeoManager cannot find path %s!",volpath.Data()));
+
+    TGeoHMatrix* m = gGeoManager->GetCurrentMatrix();
+    if(m) {
+      m->LocalToMaster(loc, glob);
+    } else {
+      AliFatal("Geo matrixes are not loaded \n") ;
+    }
   }
 }
 
+//___________________________________________________________________
 void AliEMCALGeometry::GetGlobal(Int_t absId , TVector3 &vglob) const
 { 
   // Alice numbering scheme - Jun 03, 2006
@@ -1134,7 +1193,14 @@ void AliEMCALGeometry::GetGlobal(Int_t absId , TVector3 &vglob) const
 
 }
 
-void AliEMCALGeometry::GetGlobal(const AliRecPoint *rp, TVector3 &vglob) const
+//____________________________________________________________________________
+void AliEMCALGeometry::GetGlobal(const AliRecPoint* /*rp*/, TVector3& /* vglob */) const
+{
+  AliFatal(Form("Please use GetGlobalEMCAL(recPoint,gpos) instead of GetGlobal!"));
+}
+
+//_________________________________________________________________________________
+void AliEMCALGeometry::GetGlobalEMCAL(const AliEMCALRecPoint *rp, TVector3 &vglob) const
 {
   // Figure out the global numbering
   // of a given supermodule from the
@@ -1143,9 +1209,8 @@ void AliEMCALGeometry::GetGlobal(const AliRecPoint *rp, TVector3 &vglob) const
   static TVector3 vloc;
   static Int_t nSupMod, nModule, nIphi, nIeta;
 
-  AliRecPoint *rpTmp = (AliRecPoint*)rp; // const_cast ??
-  if(!rpTmp) return;
-  AliEMCALRecPoint *rpEmc = (AliEMCALRecPoint*)rpTmp;
+  const AliEMCALRecPoint *rpTmp = rp;
+  const AliEMCALRecPoint *rpEmc = rpTmp;
 
   GetCellIndex(rpEmc->GetAbsId(0), nSupMod, nModule, nIphi, nIeta);
   rpTmp->GetLocalPosition(vloc);
@@ -1267,7 +1332,7 @@ Bool_t AliEMCALGeometry::GetAbsCellIdFromEtaPhi(Double_t eta, Double_t phi, Int_
   return kFALSE;
 }
 
-AliEMCALShishKebabTrd1Module* AliEMCALGeometry::GetShishKebabModule(Int_t neta)
+AliEMCALShishKebabTrd1Module* AliEMCALGeometry::GetShishKebabModule(Int_t neta) const
 {
   //This method was too long to be
   //included in the header file - the
@@ -1282,16 +1347,23 @@ AliEMCALShishKebabTrd1Module* AliEMCALGeometry::GetShishKebabModule(Int_t neta)
   return trd1;
 }
 
-void AliEMCALGeometry::Browse(TBrowser* b)
+void AliEMCALGeometry::Browse(TBrowser* b) const
 {
+  //Browse the modules
   if(fShishKebabTrd1Modules) b->Add(fShishKebabTrd1Modules);
-  for(int i=0; i<fNumberOfSuperModules; i++) {
-    if(fMatrixOfSM[i])  b->Add(fMatrixOfSM[i]);
-  }
 }
 
 Bool_t AliEMCALGeometry::IsFolder() const
 {
+  //Check if fShishKebabTrd1Modules is in folder
   if(fShishKebabTrd1Modules) return kTRUE;
   else                       return kFALSE;
 }
+
+Double_t AliEMCALGeometry::GetPhiCenterOfSM(Int_t nsupmod) const
+{
+  //returns center of supermodule in phi 
+  static int i = nsupmod/2;
+  return fPhiCentersOfSM[i];
+
+}