EMCAL/AliEMCALGeometry.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id$*/
  17
  18 //_________________________________________________________________________
  19 // Geometry class  for EMCAL : singleton
  20 // EMCAL consists of layers of scintillator and lead
  21 // Places the the Barrel Geometry of The EMCAL at Midrapidity
  22 // between 80 and 180(or 190) degrees of Phi and
  23 // -0.7 to 0.7 in eta
  24 // Number of Modules and Layers may be controlled by
  25 // the name of the instance defined
  26 //     EMCAL geometry tree:
  27 //     EMCAL -> superModule -> module -> tower(cell)
  28 //     Indexes
  29 //     absId -> nSupMod     -> nModule -> (nIphi,nIeta)
  30 //
  31 //*-- Author: Sahal Yacoob (LBL / UCT)
  32 //     and  : Yves Schutz (SUBATECH)
  33 //     and  : Jennifer Klay (LBL)
  34 //     SHASHLYK : Aleksei Pavlinov (WSU)
  35 //
  36
  37 #include <assert.h>
  38
  39 // --- Root header files ---
  40 #include <Riostream.h>
  41 #include <TBrowser.h>
  42 #include <TClonesArray.h>
  43 #include <TGeoManager.h>
  44 #include <TGeoMatrix.h>
  45 #include <TGeoNode.h>
  46 #include <TList.h>
  47 #include <TMatrixD.h>
  48 #include <TObjArray.h>
  49 #include <TObjString.h>
  50 #include <TVector2.h>
  51 #include <TVector3.h>
  52
  53 // -- ALICE Headers.
  54 #include "AliLog.h"
  55
  56 // --- EMCAL headers
  57 #include "AliEMCALGeometry.h"
  58 #include "AliEMCALShishKebabTrd1Module.h"
  59 #include "AliEMCALRecPoint.h"
  60 #include "AliEMCALDigit.h"
  61 #include "AliEMCALHistoUtilities.h"
  62
  63 ClassImp(AliEMCALGeometry)
  64
  65 // these initialisations are needed for a singleton
  66 AliEMCALGeometry  *AliEMCALGeometry::fgGeom      = 0;
  67 Bool_t             AliEMCALGeometry::fgInit      = kFALSE;
  68 Char_t*            AliEMCALGeometry::fgDefaultGeometryName = "SHISH_77_TRD1_2X2_FINAL_110DEG scTh=0.176 pbTh=0.144";
  69 //
  70 // Usage:
  71 //        You can create the AliEMCALGeometry object independently from anything.
  72 //        You have to use just the correct name of geometry. If name is empty string the
  73 //        default name of geometry will be used.
  74 //
  75 //  AliEMCALGeometry* g = AliEMCALGeometry::GetInstance(name,title); // first time
  76 //  ..
  77 //  g = AliEMCALGeometry::GetInstance();                             // after first time
  78 //
  79 //  MC:   If you work with MC data you have to get geometry the next way:
  80 //  ==                                      =============================
  81 //  AliRunLoader    *rl   = AliRunLoader::GetRunLoader();
  82 // AliEMCALGeometry *geom = dynamic_cast<AliEMCAL*>(rl->GetAliRun()->GetDetector("EMCAL"))->GetGeometry();
  83 // TGeoManager::Import("geometry.root");
  84
  85 AliEMCALGeometry::AliEMCALGeometry()
  86   : AliGeometry(),
  87     fGeoName(0),fArrayOpts(0),fAlFrontThick(0.),fECPbRadThickness(0.),fECScintThick(0.),
  88     fNECLayers(0),fArm1PhiMin(0.),fArm1PhiMax(0.),fArm1EtaMin(0.),fArm1EtaMax(0.),fIPDistance(0.),
  89     fShellThickness(0.),fZLength(0.),fGap2Active(0.),fNZ(0),fNPhi(0),fSampling(0.),fNumberOfSuperModules(0),
  90     fSteelFrontThick(0.),fFrontSteelStrip(0.),fLateralSteelStrip(0.),fPassiveScintThick(0.),fPhiModuleSize(0.),
  91     fEtaModuleSize(0.),fPhiTileSize(0.),fEtaTileSize(0.),fLongModuleSize(0.),fNPhiSuperModule(0),fNPHIdiv(0),fNETAdiv(0),
  92     fNCells(0),fNCellsInSupMod(0),fNCellsInModule(0),fNTRUEta(0),fNTRUPhi(0),
  93     fNCellsInTRUEta(0), fNCellsInTRUPhi(0), fTrd1Angle(0.),f2Trd1Dx2(0.),
  94     fPhiGapForSM(0.),fKey110DEG(0),fPhiBoundariesOfSM(0), fPhiCentersOfSM(0),fEtaMaxOfTRD1(0),
  95     fTrd2AngleY(0.),f2Trd2Dy2(0.),fEmptySpace(0.),fTubsR(0.),fTubsTurnAngle(0.),fCentersOfCellsEtaDir(0),
  96     fCentersOfCellsXDir(0),fCentersOfCellsPhiDir(0),fEtaCentersOfCells(0),fPhiCentersOfCells(0),
  97     fShishKebabTrd1Modules(0), fNAdditionalOpts(0),
  98     fILOSS(-1), fIHADR(-1)
  99 {
 100   // default ctor only for internal usage (singleton)
 101   // must be kept public for root persistency purposes, but should never be called by the outside world
 102   //  CreateListOfTrd1Modules();
 103   AliDebug(2, "AliEMCALGeometry : default ctor ");
 104 }
 105 //______________________________________________________________________
 106 AliEMCALGeometry::AliEMCALGeometry(const Text_t* name, const Text_t* title)
 107   : AliGeometry(name, title),
 108     fGeoName(0),fArrayOpts(0),fAlFrontThick(0.),fECPbRadThickness(0.),fECScintThick(0.),
 109     fNECLayers(0),fArm1PhiMin(0.),fArm1PhiMax(0.),fArm1EtaMin(0.),fArm1EtaMax(0.),fIPDistance(0.),
 110     fShellThickness(0.),fZLength(0.),fGap2Active(0.),fNZ(0),fNPhi(0),fSampling(0.),fNumberOfSuperModules(0),
 111     fSteelFrontThick(0.),fFrontSteelStrip(0.),fLateralSteelStrip(0.),fPassiveScintThick(0.),fPhiModuleSize(0.),
 112     fEtaModuleSize(0.),fPhiTileSize(0.),fEtaTileSize(0.),fLongModuleSize(0.),fNPhiSuperModule(0),fNPHIdiv(0),fNETAdiv(0),
 113     fNCells(0),fNCellsInSupMod(0),fNCellsInModule(0),fNTRUEta(0),fNTRUPhi(0),
 114     fNCellsInTRUEta(0), fNCellsInTRUPhi(0), fTrd1Angle(0.),f2Trd1Dx2(0.),
 115     fPhiGapForSM(0.),fKey110DEG(0),fPhiBoundariesOfSM(0), fPhiCentersOfSM(0), fEtaMaxOfTRD1(0),
 116     fTrd2AngleY(0.),f2Trd2Dy2(0.),fEmptySpace(0.),fTubsR(0.),fTubsTurnAngle(0.),fCentersOfCellsEtaDir(0),
 117     fCentersOfCellsXDir(0),fCentersOfCellsPhiDir(0),fEtaCentersOfCells(0),fPhiCentersOfCells(0),
 118     fShishKebabTrd1Modules(0),fNAdditionalOpts(0),
 119     fILOSS(-1), fIHADR(-1)
 120 {
 121   // ctor only for internal usage (singleton)
 122   AliDebug(2, Form("AliEMCALGeometry(%s,%s) ", name,title));
 123
 124   Init();
 125
 126   CreateListOfTrd1Modules();
 127
 128   if (AliDebugLevel()>=2) {
 129     PrintGeometry();
 130   }
 131
 132 }
 133 //______________________________________________________________________
 134 AliEMCALGeometry::AliEMCALGeometry(const AliEMCALGeometry& geom)
 135   : AliGeometry(geom),
 136     fGeoName(geom.fGeoName),
 137     fArrayOpts(geom.fArrayOpts),
 138     fAlFrontThick(geom.fAlFrontThick),
 139     fECPbRadThickness(geom.fECPbRadThickness),
 140     fECScintThick(geom.fECScintThick),
 141     fNECLayers(geom.fNECLayers),
 142     fArm1PhiMin(geom.fArm1PhiMin),
 143     fArm1PhiMax(geom.fArm1PhiMax),
 144     fArm1EtaMin(geom.fArm1EtaMin),
 145     fArm1EtaMax(geom.fArm1EtaMax),
 146     fIPDistance(geom.fIPDistance),
 147     fShellThickness(geom.fShellThickness),
 148     fZLength(geom.fZLength),
 149     fGap2Active(geom.fGap2Active),
 150     fNZ(geom.fNZ),
 151     fNPhi(geom.fNPhi),
 152     fSampling(geom.fSampling),
 153     fNumberOfSuperModules(geom.fNumberOfSuperModules),
 154     fSteelFrontThick(geom.fSteelFrontThick),
 155     fFrontSteelStrip(geom.fFrontSteelStrip),
 156     fLateralSteelStrip(geom.fLateralSteelStrip),
 157     fPassiveScintThick(geom.fPassiveScintThick),
 158     fPhiModuleSize(geom.fPhiModuleSize),
 159     fEtaModuleSize(geom.fEtaModuleSize),
 160     fPhiTileSize(geom.fPhiTileSize),
 161     fEtaTileSize(geom.fEtaTileSize),
 162     fLongModuleSize(geom.fLongModuleSize),
 163     fNPhiSuperModule(geom.fNPhiSuperModule),
 164     fNPHIdiv(geom.fNPHIdiv),
 165     fNETAdiv(geom.fNETAdiv),
 166     fNCells(geom.fNCells),
 167     fNCellsInSupMod(geom.fNCellsInSupMod),
 168     fNCellsInModule(geom.fNCellsInModule),
 169     fNTRUEta(geom.fNTRUEta),
 170     fNTRUPhi(geom.fNTRUPhi),
 171     fNCellsInTRUEta(geom.fNCellsInTRUEta),
 172     fNCellsInTRUPhi(geom.fNCellsInTRUPhi),
 173     fTrd1Angle(geom.fTrd1Angle),
 174     f2Trd1Dx2(geom.f2Trd1Dx2),
 175     fPhiGapForSM(geom.fPhiGapForSM),
 176     fKey110DEG(geom.fKey110DEG),
 177     fPhiBoundariesOfSM(geom.fPhiBoundariesOfSM),
 178     fPhiCentersOfSM(geom.fPhiCentersOfSM),
 179     fEtaMaxOfTRD1(geom.fEtaMaxOfTRD1),
 180     fTrd2AngleY(geom.fTrd2AngleY),
 181     f2Trd2Dy2(geom.f2Trd2Dy2),
 182     fEmptySpace(geom.fEmptySpace),
 183     fTubsR(geom.fTubsR),
 184     fTubsTurnAngle(geom.fTubsTurnAngle),
 185     fCentersOfCellsEtaDir(geom.fCentersOfCellsEtaDir),
 186     fCentersOfCellsXDir(geom.fCentersOfCellsXDir),
 187     fCentersOfCellsPhiDir(geom.fCentersOfCellsPhiDir),
 188     fEtaCentersOfCells(geom.fEtaCentersOfCells),
 189     fPhiCentersOfCells(geom.fPhiCentersOfCells),
 190     fShishKebabTrd1Modules(geom.fShishKebabTrd1Modules),
 191     fNAdditionalOpts(geom.fNAdditionalOpts),
 192     fILOSS(geom.fILOSS), fIHADR(geom.fIHADR)
 193 {
 194   //copy ctor
 195 }
 196
 197 //______________________________________________________________________
 198 AliEMCALGeometry::~AliEMCALGeometry(void){
 199     // dtor
 200 }
 201 //______________________________________________________________________
 202 void AliEMCALGeometry::Init(void){
 203   // Initializes the EMCAL parameters
 204   // naming convention : GUV_WX_N_ gives the composition of a tower
 205   // WX inform about the composition of the EM calorimeter section:
 206   //   thickness in mm of Pb radiator (W) and of scintillator (X), and number of scintillator layers (N)
 207   // New geometry: EMCAL_55_25
 208   // 24-aug-04 for shish-kebab
 209   // SHISH_25 or SHISH_62
 210   // 11-oct-05   - correction for pre final design
 211   // Feb 06,2006 - decrease the weight of EMCAL
 212   //
 213   // Oct 30,2006 - SHISH_TRD1_CURRENT_1X1, SHISH_TRD1_CURRENT_2X2 or SHISH_TRD1_CURRENT_3X3;
 214   //
 215
 216   fAdditionalOpts[0] = "nl=";       // number of sampling layers (fNECLayers)
 217   fAdditionalOpts[1] = "pbTh=";     // cm, Thickness of the Pb   (fECPbRadThick)
 218   fAdditionalOpts[2] = "scTh=";     // cm, Thickness of the Sc    (fECScintThick)
 219   fAdditionalOpts[3] = "latSS=";    // cm, Thickness of lateral steel strip (fLateralSteelStrip)
 220   fAdditionalOpts[4] = "allILOSS="; // = 0,1,2,3,4 (4 - energy loss without fluctuation)
 221   fAdditionalOpts[5] = "allIHADR="; // = 0,1,2 (0 - no hadronic interaction)
 222
 223   fNAdditionalOpts = sizeof(fAdditionalOpts) / sizeof(char*);
 224
 225   fgInit = kFALSE; // Assume failed until proven otherwise.
 226   fGeoName   = GetName();
 227   fGeoName.ToUpper();
 228   fKey110DEG = 0;
 229   if(fGeoName.Contains("110DEG") || fGeoName.Contains("CURRENT")) fKey110DEG = 1; // for GetAbsCellId
 230   fShishKebabTrd1Modules = 0;
 231   fTrd2AngleY = f2Trd2Dy2 = fEmptySpace = fTubsR = fTubsTurnAngle = 0;
 232
 233   fNZ             = 114;        // granularity along Z (eta)
 234   fNPhi           = 168;        // granularity in phi (azimuth)
 235   fArm1PhiMin     = 80.0;       // degrees, Starting EMCAL Phi position
 236   fArm1PhiMax     = 190.0;      // degrees, Ending EMCAL Phi position
 237   fArm1EtaMin     = -0.7;       // pseudorapidity, Starting EMCAL Eta position
 238   fArm1EtaMax     = +0.7;       // pseudorapidity, Ending EMCAL Eta position
 239   fIPDistance     = 454.0;      // cm, Radial distance to inner surface of EMCAL
 240   fPhiGapForSM    = 0.;         // cm, only for final TRD1 geometry
 241
 242   // geometry
 243   if(fGeoName.Contains("SHISH")){ // Only shahslyk now
 244     // 7-sep-05; integration issue
 245     fArm1PhiMin     = 80.0;     // 60  -> 80
 246     fArm1PhiMax     = 180.0;    // 180 -> 190
 247
 248     fNumberOfSuperModules = 10; // 12 = 6 * 2 (6 in phi, 2 in Z);
 249     fSteelFrontThick = 2.54;    //  9-sep-04
 250     fIPDistance      = 460.0;
 251     fFrontSteelStrip = fPassiveScintThick = 0.0; // 13-may-05
 252     fLateralSteelStrip = 0.025; // before MAY 2005
 253     fPhiModuleSize   = fEtaModuleSize   = 11.4;
 254     fPhiTileSize = fEtaTileSize      = 5.52; // (11.4-5.52*2)/2. = 0.18 cm (wall thickness)
 255     fNPhi            = 14;
 256     fNZ              = 30;
 257     fAlFrontThick    = fGap2Active = 0;
 258     fNPHIdiv = fNETAdiv = 2;
 259
 260     fNECLayers       = 62;
 261     fECScintThick    = fECPbRadThickness = 0.2;
 262     fSampling        = 1.;  // 30-aug-04 - should be calculated
 263     if(fGeoName.Contains("TWIST")) { // all about EMCAL module
 264       fNZ             = 27;  // 16-sep-04
 265     } else if(fGeoName.Contains("TRD")) {
 266       fIPDistance      = 428.0;  //  11-may-05
 267       fSteelFrontThick = 0.0;    // 3.17 -> 0.0; 28-mar-05 : no stell plate
 268       fNPhi            = 12;
 269       fSampling       = 12.327;
 270       fPhiModuleSize = fEtaModuleSize = 12.26;
 271       fNZ            = 26;     // 11-oct-04
 272       fTrd1Angle     = 1.3;    // in degree
 273 // 18-nov-04; 1./0.08112=12.327
 274 // http://pdsfweb01.nersc.gov/~pavlinov/ALICE/SHISHKEBAB/RES/linearityAndResolutionForTRD1.html
 275       if(fGeoName.Contains("TRD1")) {       // 30-jan-05
 276         // for final design
 277         fPhiGapForSM    = 2.;         // cm, only for final TRD1 geometry
 278         if(fGeoName.Contains("MAY05") || fGeoName.Contains("WSUC") || fGeoName.Contains("FINAL") || fGeoName.Contains("CURRENT")){
 279           fNumberOfSuperModules = 12; // 20-may-05
 280           if(fGeoName.Contains("WSUC")) fNumberOfSuperModules = 1; // 27-may-05
 281           fNECLayers     = 77;       // (13-may-05 from V.Petrov)
 282           fPhiModuleSize = 12.5;     // 20-may-05 - rectangular shape
 283           fEtaModuleSize = 11.9;
 284           fECScintThick  = fECPbRadThickness = 0.16;// (13-may-05 from V.Petrov)
 285           fFrontSteelStrip   = 0.025;// 0.025cm = 0.25mm  (13-may-05 from V.Petrov)
 286           fLateralSteelStrip = 0.01; // 0.01cm  = 0.1mm   (13-may-05 from V.Petrov) - was 0.025
 287           fPassiveScintThick = 0.8;  // 0.8cm   = 8mm     (13-may-05 from V.Petrov)
 288           fNZ                = 24;
 289           fTrd1Angle         = 1.5;  // 1.3 or 1.5
 290
 291           if(fGeoName.Contains("FINAL") || fGeoName.Contains("CURRENT")) { // 9-sep-05
 292             fNumberOfSuperModules = 10;
 293             if(GetKey110DEG()) {
 294               fNumberOfSuperModules = 12;// last two modules have size 10 degree in phi (180<phi<190)
 295               fArm1PhiMax = 200.0;       // for XEN1 and turn angle of super modules
 296             }
 297             if(fGeoName.Contains("FINAL")) {
 298               fPhiModuleSize = 12.26 - fPhiGapForSM / Float_t(fNPhi); // first assumption
 299             } else if(fGeoName.Contains("CURRENT")) {
 300               fECScintThick      = 0.176; // 10% of weight reduction
 301               fECPbRadThickness  = 0.144; //
 302               fLateralSteelStrip = 0.015; // 0.015cm  = 0.15mm (Oct 30, from Fred)
 303               fPhiModuleSize     = 12.00;
 304               fPhiGapForSM       = (12.26 - fPhiModuleSize)*fNPhi; // have to check
 305             }
 306             fEtaModuleSize = fPhiModuleSize;
 307             if(fGeoName.Contains("HUGE")) fNECLayers *= 3; // 28-oct-05 for analysing leakage
 308           }
 309         }
 310       } else if(fGeoName.Contains("TRD2")) {       // 30-jan-05
 311         fSteelFrontThick = 0.0;         // 11-mar-05
 312         fIPDistance+= fSteelFrontThick; // 1-feb-05 - compensate absence of steel plate
 313         fTrd1Angle  = 1.64;             // 1.3->1.64
 314         fTrd2AngleY = fTrd1Angle;       //  symmetric case now
 315         fEmptySpace    = 0.2; // 2 mm
 316         fTubsR         = fIPDistance; // 31-jan-05 - as for Fred case
 317
 318         fPhiModuleSize  = fTubsR*2.*TMath::Tan(fTrd2AngleY*TMath::DegToRad()/2.);
 319         fPhiModuleSize -= fEmptySpace/2.; // 11-mar-05
 320         fEtaModuleSize  = fPhiModuleSize; // 20-may-05
 321         fTubsTurnAngle  = 3.;
 322       }
 323       fNPHIdiv = fNETAdiv  = 2;   // 13-oct-04 - division again
 324       if(fGeoName.Contains("3X3")) {   // 23-nov-04
 325         fNPHIdiv = fNETAdiv  = 3;
 326       } else if(fGeoName.Contains("4X4")) {
 327         fNPHIdiv = fNETAdiv  = 4;
 328       } else if(fGeoName.Contains("1X1")) {
 329         fNPHIdiv = fNETAdiv  = 1;
 330       }
 331     }
 332     if(fGeoName.Contains("25")){
 333       fNECLayers     = 25;
 334       fECScintThick  = fECPbRadThickness = 0.5;
 335     }
 336     if(fGeoName.Contains("WSUC")){ // 18-may-05 - about common structure
 337       fShellThickness = 30.;       // should be change
 338       fNPhi = fNZ = 4;
 339     }
 340
 341     CheckAdditionalOptions();
 342     DefineSamplingFraction();
 343
 344     fPhiTileSize = fPhiModuleSize/double(fNPHIdiv) - fLateralSteelStrip; // 13-may-05
 345     fEtaTileSize = fEtaModuleSize/double(fNETAdiv) - fLateralSteelStrip; // 13-may-05
 346
 347     // constant for transition absid <--> indexes
 348     fNCellsInModule  = fNPHIdiv*fNETAdiv;
 349     fNCellsInSupMod = fNCellsInModule*fNPhi*fNZ;
 350     fNCells         = fNCellsInSupMod*fNumberOfSuperModules;
 351     if(GetKey110DEG()) fNCells -= fNCellsInSupMod;
 352
 353     fLongModuleSize = fNECLayers*(fECScintThick + fECPbRadThickness);
 354     if(fGeoName.Contains("MAY05")) fLongModuleSize += (fFrontSteelStrip + fPassiveScintThick);
 355
 356     // 30-sep-04
 357     if(fGeoName.Contains("TRD")) {
 358       f2Trd1Dx2 = fEtaModuleSize + 2.*fLongModuleSize*TMath::Tan(fTrd1Angle*TMath::DegToRad()/2.);
 359       if(fGeoName.Contains("TRD2")) {  // 27-jan-05
 360         f2Trd2Dy2 = fPhiModuleSize + 2.*fLongModuleSize*TMath::Tan(fTrd2AngleY*TMath::DegToRad()/2.);
 361       }
 362     }
 363   } else Fatal("Init", "%s is an undefined geometry!", fGeoName.Data()) ;
 364
 365   fNPhiSuperModule = fNumberOfSuperModules/2;
 366   if(fNPhiSuperModule<1) fNPhiSuperModule = 1;
 367
 368   fShellThickness = fAlFrontThick + fGap2Active + fNECLayers*GetECScintThick()+(fNECLayers-1)*GetECPbRadThick();
 369   if(fGeoName.Contains("SHISH")) {
 370     fShellThickness = fSteelFrontThick + fLongModuleSize;
 371     if(fGeoName.Contains("TWIST")) { // 13-sep-04
 372       fShellThickness  = TMath::Sqrt(fLongModuleSize*fLongModuleSize + fPhiModuleSize*fEtaModuleSize);
 373       fShellThickness += fSteelFrontThick;
 374     } else if(fGeoName.Contains("TRD")) { // 1-oct-04
 375       fShellThickness  = TMath::Sqrt(fLongModuleSize*fLongModuleSize + f2Trd1Dx2*f2Trd1Dx2);
 376       fShellThickness += fSteelFrontThick;
 377       // Local coordinates
 378       fParSM[0] = GetShellThickness()/2.;
 379       fParSM[1] = GetPhiModuleSize() * GetNPhi()/2.;
 380       fParSM[2] = 350./2.;
 381     }
 382   }
 383
 384   fZLength        = 2.*ZFromEtaR(fIPDistance+fShellThickness,fArm1EtaMax); // Z coverage
 385   fEnvelop[0]     = fIPDistance; // mother volume inner radius
 386   fEnvelop[1]     = fIPDistance + fShellThickness; // mother volume outer r.
 387   fEnvelop[2]     = 1.00001*fZLength; // add some padding for mother volume.
 388
 389   fNumberOfSuperModules = 12;
 390
 391   // SM phi boundaries - (0,1),(2,3) .. (10,11) - has the same boundaries; Nov 7, 2006
 392   fPhiBoundariesOfSM.Set(fNumberOfSuperModules);
 393   fPhiCentersOfSM.Set(fNumberOfSuperModules/2);
 394   fPhiBoundariesOfSM[0] = TMath::PiOver2() - TMath::ATan2(fParSM[1] , fIPDistance); // 1th and 2th modules)
 395   fPhiBoundariesOfSM[1] = TMath::PiOver2() + TMath::ATan2(fParSM[1] , fIPDistance);
 396   fPhiCentersOfSM[0]     = TMath::PiOver2();
 397   for(int i=1; i<=4; i++) { // from 2th ro 9th
 398     fPhiBoundariesOfSM[2*i]   = fPhiBoundariesOfSM[0] + 20.*TMath::DegToRad()*i;
 399     fPhiBoundariesOfSM[2*i+1] = fPhiBoundariesOfSM[1] + 20.*TMath::DegToRad()*i;
 400     fPhiCentersOfSM[i]         = fPhiCentersOfSM[0]     + 20.*TMath::DegToRad()*i;
 401   }
 402   fPhiBoundariesOfSM[11] = 190.*TMath::DegToRad();
 403   fPhiBoundariesOfSM[10] = fPhiBoundariesOfSM[11] - TMath::ATan2((fParSM[1]) , fIPDistance);
 404   fPhiCentersOfSM[5]      = (fPhiBoundariesOfSM[10]+fPhiBoundariesOfSM[11])/2.;
 405
 406   //TRU parameters. These parameters values are not the final ones.
 407   fNTRUEta = 3 ;
 408   fNTRUPhi = 1 ;
 409   fNCellsInTRUEta = 16 ;
 410   fNCellsInTRUPhi = 24 ;
 411
 412   if(fGeoName.Contains("WSUC")) fNumberOfSuperModules = 1; // Jul 12, 2007
 413
 414   fgInit = kTRUE;
 415 }
 416
 417 void AliEMCALGeometry::PrintGeometry()
 418 {
 419   // Separate routine is callable from broswer; Nov 7,2006
 420   printf("\nInit: geometry of EMCAL named %s :\n", fGeoName.Data());
 421   if(fArrayOpts) {
 422     for(Int_t i=0; i<fArrayOpts->GetEntries(); i++){
 423       TObjString *o = (TObjString*)fArrayOpts->At(i);
 424       printf(" %i : %s \n", i, o->String().Data());
 425     }
 426   }
 427   printf("Granularity: %d in eta and %d in phi\n", GetNZ(), GetNPhi()) ;
 428   printf("Layout: phi = (%7.1f, %7.1f), eta = (%5.2f, %5.2f), IP = %7.2f -> for EMCAL envelope only\n",
 429            GetArm1PhiMin(), GetArm1PhiMax(),GetArm1EtaMin(), GetArm1EtaMax(), GetIPDistance() );
 430
 431   printf( "               ECAL      : %d x (%f cm Pb, %f cm Sc) \n",
 432   GetNECLayers(), GetECPbRadThick(), GetECScintThick() ) ;
 433   printf("                fSampling %5.2f \n",  fSampling );
 434   if(fGeoName.Contains("SHISH")){
 435     printf(" fIPDistance       %6.3f cm \n", fIPDistance);
 436     if(fSteelFrontThick>0.)
 437     printf(" fSteelFrontThick  %6.3f cm \n", fSteelFrontThick);
 438     printf(" fNPhi %i   |  fNZ %i \n", fNPhi, fNZ);
 439     printf(" fNCellsInModule %i : fNCellsInSupMod %i : fNCells %i\n",fNCellsInModule, fNCellsInSupMod, fNCells);
 440     if(fGeoName.Contains("MAY05")){
 441       printf(" fFrontSteelStrip         %6.4f cm (thickness of front steel strip)\n",
 442       fFrontSteelStrip);
 443       printf(" fLateralSteelStrip       %6.4f cm (thickness of lateral steel strip)\n",
 444       fLateralSteelStrip);
 445       printf(" fPassiveScintThick  %6.4f cm (thickness of front passive Sc tile)\n",
 446       fPassiveScintThick);
 447     }
 448     printf(" X:Y module size     %6.3f , %6.3f cm \n", fPhiModuleSize, fEtaModuleSize);
 449     printf(" X:Y   tile size     %6.3f , %6.3f cm \n", fPhiTileSize, fEtaTileSize);
 450     printf(" #of sampling layers %i(fNECLayers) \n", fNECLayers);
 451     printf(" fLongModuleSize     %6.3f cm \n", fLongModuleSize);
 452     printf(" #supermodule in phi direction %i \n", fNPhiSuperModule );
 453   }
 454   printf(" fILOSS %i : fIHADR %i \n", fILOSS, fIHADR);
 455   if(fGeoName.Contains("TRD")) {
 456     printf(" fTrd1Angle %7.4f\n", fTrd1Angle);
 457     printf(" f2Trd1Dx2  %7.4f\n",  f2Trd1Dx2);
 458     if(fGeoName.Contains("TRD2")) {
 459       printf(" fTrd2AngleY     %7.4f\n", fTrd2AngleY);
 460       printf(" f2Trd2Dy2       %7.4f\n", f2Trd2Dy2);
 461       printf(" fTubsR          %7.2f cm\n", fTubsR);
 462       printf(" fTubsTurnAngle  %7.4f\n", fTubsTurnAngle);
 463       printf(" fEmptySpace     %7.4f cm\n", fEmptySpace);
 464     } else if(fGeoName.Contains("TRD1")){
 465       printf("SM dimensions(TRD1) : dx %7.2f dy %7.2f dz %7.2f (SMOD, BOX)\n",
 466       fParSM[0],fParSM[1],fParSM[2]);
 467       printf(" fPhiGapForSM  %7.4f cm (%7.4f <- phi size in degree)\n",
 468       fPhiGapForSM, TMath::ATan2(fPhiGapForSM,fIPDistance)*TMath::RadToDeg());
 469       if(GetKey110DEG()) printf(" Last two modules have size 10 degree in  phi (180<phi<190)\n");
 470       printf(" phi SM boundaries \n");
 471       for(int i=0; i<fPhiBoundariesOfSM.GetSize()/2.; i++) {
 472         printf(" %i : %7.5f(%7.2f) -> %7.5f(%7.2f) : center %7.5f(%7.2f) \n", i,
 473         fPhiBoundariesOfSM[2*i], fPhiBoundariesOfSM[2*i]*TMath::RadToDeg(),
 474                fPhiBoundariesOfSM[2*i+1], fPhiBoundariesOfSM[2*i+1]*TMath::RadToDeg(),
 475                fPhiCentersOfSM[i], fPhiCentersOfSM[i]*TMath::RadToDeg());
 476       }
 477       printf(" fShishKebabTrd1Modules has %i modules : max eta %5.4f \n",
 478                fShishKebabTrd1Modules->GetSize(),fEtaMaxOfTRD1);
 479
 480       printf("\n Cells grid in eta directions : size %i\n", fCentersOfCellsEtaDir.GetSize());
 481       for(Int_t i=0; i<fCentersOfCellsEtaDir.GetSize(); i++) {
 482         printf(" ind %2.2i : z %8.3f : x %8.3f \n", i,
 483                fCentersOfCellsEtaDir.At(i),fCentersOfCellsXDir.At(i));
 484         int ind=0; // Nov 21,2006
 485         for(Int_t iphi=0; iphi<fCentersOfCellsPhiDir.GetSize(); iphi++) {
 486           ind = iphi*fCentersOfCellsEtaDir.GetSize() + i;
 487           printf("%6.4f ", fEtaCentersOfCells[ind]);
 488           if((iphi+1)%12 == 0) printf("\n");
 489         }
 490         printf("\n");
 491
 492       }
 493
 494       printf("\n Cells grid in phi directions : size %i\n", fCentersOfCellsPhiDir.GetSize());
 495       for(Int_t i=0; i<fCentersOfCellsPhiDir.GetSize(); i++) {
 496         double phi=fPhiCentersOfCells.At(i);
 497         printf(" ind %2.2i : y %8.3f : phi %7.5f(%6.2f) \n", i, fCentersOfCellsPhiDir.At(i),
 498                phi, phi*TMath::RadToDeg());
 499       }
 500     }
 501   }
 502 }
 503
 504 void AliEMCALGeometry::PrintCellIndexes(Int_t absId, int pri, char *tit)
 505 {
 506   // Service methods
 507   Int_t nSupMod, nModule, nIphi, nIeta;
 508   Int_t iphi, ieta;
 509   TVector3 vg;
 510
 511   GetCellIndex(absId,  nSupMod, nModule, nIphi, nIeta);
 512   printf(" %s | absId : %i -> nSupMod %i nModule %i nIphi %i nIeta %i \n", tit, absId,  nSupMod, nModule, nIphi, nIeta);
 513   if(pri>0) {
 514     GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi,ieta);
 515     printf(" local SM index : iphi %i : ieta %i \n", iphi,ieta);
 516     GetGlobal(absId, vg);
 517     printf(" vglob : mag %7.2f : perp %7.2f : z %7.2f : eta %6.4f : phi %6.4f(%6.2f) \n",
 518            vg.Mag(), vg.Perp(), vg.Z(), vg.Eta(), vg.Phi(), vg.Phi()*TMath::RadToDeg());
 519   }
 520 }
 521
 522 //______________________________________________________________________
 523 void AliEMCALGeometry::CheckAdditionalOptions()
 524 {
 525   // Feb 06,2006
 526   // Additional options that
 527   // can be used to select
 528   // the specific geometry of
 529   // EMCAL to run
 530   // Dec 27,2006
 531   // adeed allILOSS= and allIHADR= for MIP investigation
 532   fArrayOpts = new TObjArray;
 533   Int_t nopt = AliEMCALHistoUtilities::ParseString(fGeoName, *fArrayOpts);
 534   if(nopt==1) { // no aditional option(s)
 535     fArrayOpts->Delete();
 536     delete fArrayOpts;
 537     fArrayOpts = 0;
 538     return;
 539   }
 540   for(Int_t i=1; i<nopt; i++){
 541     TObjString *o = (TObjString*)fArrayOpts->At(i);
 542
 543     TString addOpt = o->String();
 544     Int_t indj=-1;
 545     for(Int_t j=0; j<fNAdditionalOpts; j++) {
 546       TString opt = fAdditionalOpts[j];
 547       if(addOpt.Contains(opt,TString::kIgnoreCase)) {
 548           indj = j;
 549         break;
 550       }
 551     }
 552     if(indj<0) {
 553       AliDebug(2,Form("<E> option |%s| unavailable : ** look to the file AliEMCALGeometry.h **\n",
 554                       addOpt.Data()));
 555       assert(0);
 556     } else {
 557       AliDebug(2,Form("<I> option |%s| is valid : number %i : |%s|\n",
 558                       addOpt.Data(), indj, fAdditionalOpts[indj]));
 559       if       (addOpt.Contains("NL=",TString::kIgnoreCase))   {// number of sampling layers
 560         sscanf(addOpt.Data(),"NL=%i", &fNECLayers);
 561         AliDebug(2,Form(" fNECLayers %i (new) \n", fNECLayers));
 562       } else if(addOpt.Contains("PBTH=",TString::kIgnoreCase)) {//Thickness of the Pb(fECPbRadThicknes)
 563         sscanf(addOpt.Data(),"PBTH=%f", &fECPbRadThickness);
 564       } else if(addOpt.Contains("SCTH=",TString::kIgnoreCase)) {//Thickness of the Sc(fECScintThick)
 565         sscanf(addOpt.Data(),"SCTH=%f", &fECScintThick);
 566       } else if(addOpt.Contains("LATSS=",TString::kIgnoreCase)) {// Thickness of lateral steel strip (fLateralSteelStrip)
 567         sscanf(addOpt.Data(),"LATSS=%f", &fLateralSteelStrip);
 568         AliDebug(2,Form(" fLateralSteelStrip %f (new) \n", fLateralSteelStrip));
 569       } else if(addOpt.Contains("ILOSS=",TString::kIgnoreCase)) {// As in Geant
 570         sscanf(addOpt.Data(),"ALLILOSS=%i", &fILOSS);
 571         AliDebug(2,Form(" fILOSS %i \n", fILOSS));
 572       } else if(addOpt.Contains("IHADR=",TString::kIgnoreCase)) {// As in Geant
 573         sscanf(addOpt.Data(),"ALLIHADR=%i", &fIHADR);
 574         AliDebug(2,Form(" fIHADR %i \n", fIHADR));
 575       }
 576     }
 577   }
 578 }
 579
 580 void AliEMCALGeometry::DefineSamplingFraction()
 581 {
 582   // Jun 05,2006
 583   // Look http://rhic.physics.wayne.edu/~pavlinov/ALICE/SHISHKEBAB/RES/linearityAndResolutionForTRD1.html
 584   // Keep for compatibilty
 585   //
 586   if(fNECLayers == 69) {        // 10% layer reduction
 587     fSampling = 12.55;
 588   } else if(fNECLayers == 61) { // 20% layer reduction
 589     fSampling = 12.80;
 590   } else if(fNECLayers == 77) {
 591     if       (fECScintThick>0.175 && fECScintThick<0.177) { // 10% Pb thicknes reduction
 592       fSampling = 10.5; // fECScintThick = 0.176, fECPbRadThickness=0.144;
 593     } else if(fECScintThick>0.191 && fECScintThick<0.193) { // 20% Pb thicknes reduction
 594       fSampling = 8.93; // fECScintThick = 0.192, fECPbRadThickness=0.128;
 595     }
 596   }
 597 }
 598
 599 //______________________________________________________________________
 600 void AliEMCALGeometry::GetCellPhiEtaIndexInSModuleFromTRUIndex(const Int_t itru, const Int_t iphitru, const Int_t ietatru, Int_t &iphiSM, Int_t &ietaSM) const
 601 {
 602
 603   // This method transforms the (eta,phi) index of cells in a
 604   // TRU matrix into Super Module (eta,phi) index.
 605
 606   // Calculate in which row and column where the TRU are
 607   // ordered in the SM
 608
 609   Int_t col = itru/ fNTRUPhi ;
 610   Int_t row = itru - col*fNTRUPhi ;
 611
 612   iphiSM = fNCellsInTRUPhi*row + iphitru  ;
 613   ietaSM = fNCellsInTRUEta*col + ietatru  ;
 614 }
 615
 616 //______________________________________________________________________
 617 AliEMCALGeometry *  AliEMCALGeometry::GetInstance(){
 618   // Returns the pointer of the unique instance
 619
 620   AliEMCALGeometry * rv = static_cast<AliEMCALGeometry *>( fgGeom );
 621   return rv;
 622 }
 623
 624 //______________________________________________________________________
 625 AliEMCALGeometry* AliEMCALGeometry::GetInstance(const Text_t* name,
 626                                                 const Text_t* title){
 627     // Returns the pointer of the unique instance
 628
 629     AliEMCALGeometry * rv = 0;
 630     if ( fgGeom == 0 ) {
 631       if ( strcmp(name,"") == 0 ) { // get default geometry
 632          fgGeom = new AliEMCALGeometry(fgDefaultGeometryName, title);
 633       } else {
 634          fgGeom = new AliEMCALGeometry(name, title);
 635       }  // end if strcmp(name,"")
 636       if ( fgInit ) rv = (AliEMCALGeometry * ) fgGeom;
 637       else {
 638          rv = 0;
 639          delete fgGeom;
 640          fgGeom = 0;
 641       } // end if fgInit
 642     }else{
 643         if ( strcmp(fgGeom->GetName(), name) != 0) {
 644           printf("\ncurrent geometry is %s : ", fgGeom->GetName());
 645           printf(" you cannot call %s ",name);
 646         }else{
 647           rv = (AliEMCALGeometry *) fgGeom;
 648         } // end
 649     }  // end if fgGeom
 650     return rv;
 651 }
 652
 653 Bool_t AliEMCALGeometry::IsInEMCAL(Double_t x, Double_t y, Double_t z) const {
 654   // Checks whether point is inside the EMCal volume, used in AliEMCALv*.cxx
 655   //
 656   // Code uses cylindrical approximation made of inner radius (for speed)
 657   //
 658   // Points behind EMCAl, i.e. R > outer radius, but eta, phi in acceptance
 659   // are considered to inside
 660
 661   Double_t r=sqrt(x*x+y*y);
 662
 663   if ( r > fEnvelop[0] ) {
 664      Double_t theta;
 665      theta  =    TMath::ATan2(r,z);
 666      Double_t eta;
 667      if(theta == 0)
 668        eta = 9999;
 669      else
 670        eta    =   -TMath::Log(TMath::Tan(theta/2.));
 671      if (eta < fArm1EtaMin || eta > fArm1EtaMax)
 672        return 0;
 673
 674      Double_t phi = TMath::ATan2(y,x) * 180./TMath::Pi();
 675      if (phi < 0) phi += 360;  // phi should go from 0 to 360 in this case
 676      if (phi > fArm1PhiMin && phi < fArm1PhiMax)
 677        return 1;
 678   }
 679   return 0;
 680 }
 681 // ==
 682
 683 //
 684 // == Shish-kebab cases ==
 685 //
 686 Int_t AliEMCALGeometry::GetAbsCellId(Int_t nSupMod, Int_t nModule, Int_t nIphi, Int_t nIeta) const
 687 {
 688   // 27-aug-04;
 689   // corr. 21-sep-04;
 690   //       13-oct-05; 110 degree case
 691   // May 31, 2006; ALICE numbering scheme:
 692   // 0 <= nSupMod < fNumberOfSuperModules
 693   // 0 <= nModule  < fNPHI * fNZ ( fNPHI * fNZ/2 for fKey110DEG=1)
 694   // 0 <= nIphi   < fNPHIdiv
 695   // 0 <= nIeta   < fNETAdiv
 696   // 0 <= absid   < fNCells
 697   static Int_t id=0; // have to change from 0 to fNCells-1
 698   if(fKey110DEG == 1 && nSupMod >= 10) { // 110 degree case; last two supermodules
 699     id  = fNCellsInSupMod*10 + (fNCellsInSupMod/2)*(nSupMod-10);
 700   } else {
 701     id  = fNCellsInSupMod*nSupMod;
 702   }
 703   id += fNCellsInModule *nModule;
 704   id += fNPHIdiv *nIphi;
 705   id += nIeta;
 706   if(id<0 || id >= fNCells) {
 707 //     printf(" wrong numerations !!\n");
 708 //     printf("    id      %6i(will be force to -1)\n", id);
 709 //     printf("    fNCells %6i\n", fNCells);
 710 //     printf("    nSupMod %6i\n", nSupMod);
 711 //     printf("    nModule  %6i\n", nModule);
 712 //     printf("    nIphi   %6i\n", nIphi);
 713 //     printf("    nIeta   %6i\n", nIeta);
 714     id = -TMath::Abs(id); // if negative something wrong
 715   }
 716   return id;
 717 }
 718
 719 Bool_t  AliEMCALGeometry::CheckAbsCellId(Int_t absId) const
 720 {
 721   // May 31, 2006; only trd1 now
 722   if(absId<0 || absId >= fNCells) return kFALSE;
 723   else                            return kTRUE;
 724 }
 725
 726 Bool_t AliEMCALGeometry::GetCellIndex(Int_t absId,Int_t &nSupMod,Int_t &nModule,Int_t &nIphi,Int_t &nIeta) const
 727 {
 728   // 21-sep-04; 19-oct-05;
 729   // May 31, 2006; ALICE numbering scheme:
 730   //
 731   // In:
 732   // absId   - cell is as in Geant,     0<= absId   < fNCells;
 733   // Out:
 734   // nSupMod - super module(SM) number, 0<= nSupMod < fNumberOfSuperModules;
 735   // nModule  - module number in SM,     0<= nModule  < fNCellsInSupMod/fNCellsInSupMod or(/2) for tow last SM (10th and 11th);
 736   // nIphi   - cell number in phi driection inside module; 0<= nIphi < fNPHIdiv;
 737   // nIeta   - cell number in eta driection inside module; 0<= nIeta < fNETAdiv;
 738   //
 739   static Int_t tmp=0, sm10=0;
 740   if(!CheckAbsCellId(absId)) return kFALSE;
 741
 742   sm10 = fNCellsInSupMod*10;
 743   if(fKey110DEG == 1 && absId >= sm10) { // 110 degree case; last two supermodules
 744     nSupMod = (absId-sm10) / (fNCellsInSupMod/2) + 10;
 745     tmp     = (absId-sm10) % (fNCellsInSupMod/2);
 746   } else {
 747     nSupMod = absId / fNCellsInSupMod;
 748     tmp     = absId % fNCellsInSupMod;
 749   }
 750
 751   nModule  = tmp / fNCellsInModule;
 752   tmp     = tmp % fNCellsInModule;
 753   nIphi   = tmp / fNPHIdiv;
 754   nIeta   = tmp % fNPHIdiv;
 755
 756   return kTRUE;
 757 }
 758
 759 void AliEMCALGeometry::GetModulePhiEtaIndexInSModule(Int_t nSupMod, Int_t nModule,  int &iphim, int &ietam) const
 760 {
 761   // added nSupMod; - 19-oct-05 !
 762   // Alice numbering scheme        - Jun 01,2006
 763   // ietam, iphi - indexes of module in two dimensional grid of SM
 764   // ietam - have to change from 0 to fNZ-1
 765   // iphim - have to change from 0 to nphi-1 (fNPhi-1 or fNPhi/2-1)
 766   static Int_t nphi;
 767
 768   if(fKey110DEG == 1 && nSupMod>=10) nphi = fNPhi/2;
 769   else                               nphi = fNPhi;
 770
 771   ietam = nModule/nphi;
 772   iphim = nModule%nphi;
 773 }
 774
 775 void AliEMCALGeometry::GetCellPhiEtaIndexInSModule(Int_t nSupMod, Int_t nModule, Int_t nIphi, Int_t nIeta,
 776 int &iphi, int &ieta) const
 777 {
 778   //
 779   // Added nSupMod; Nov 25, 05
 780   // Alice numbering scheme  - Jun 01,2006
 781   // IN:
 782   // nSupMod - super module(SM) number, 0<= nSupMod < fNumberOfSuperModules;
 783   // nModule  - module number in SM,     0<= nModule  < fNCellsInSupMod/fNCellsInSupMod or(/2) for tow last SM (10th and 11th);
 784   // nIphi   - cell number in phi driection inside module; 0<= nIphi < fNPHIdiv;
 785   // nIeta   - cell number in eta driection inside module; 0<= nIeta < fNETAdiv;
 786   //
 787  // OUT:
 788   // ieta, iphi - indexes of cell(tower) in two dimensional grid of SM
 789   // ieta - have to change from 0 to (fNZ*fNETAdiv-1)
 790   // iphi - have to change from 0 to (fNPhi*fNPHIdiv-1 or fNPhi*fNPHIdiv/2-1)
 791   //
 792   static Int_t iphim, ietam;
 793
 794   GetModulePhiEtaIndexInSModule(nSupMod,nModule, iphim, ietam);
 795   //  ieta  = ietam*fNETAdiv + (1-nIeta); // x(module) = -z(SM)
 796   ieta  = ietam*fNETAdiv + (fNETAdiv - 1 - nIeta); // x(module) = -z(SM)
 797   iphi  = iphim*fNPHIdiv + nIphi;     // y(module) =  y(SM)
 798
 799   if(iphi<0 || ieta<0)
 800   AliDebug(1,Form(" nSupMod %i nModule %i nIphi %i nIeta %i => ieta %i iphi %i\n",
 801   nSupMod, nModule, nIphi, nIeta, ieta, iphi));
 802 }
 803
 804 Int_t  AliEMCALGeometry::GetSuperModuleNumber(Int_t absId)  const
 805 {
 806   // Return the number of the  supermodule given the absolute
 807   // ALICE numbering id
 808
 809   static Int_t nSupMod, nModule, nIphi, nIeta;
 810   GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
 811   return nSupMod;
 812 }
 813
 814 void  AliEMCALGeometry::GetModuleIndexesFromCellIndexesInSModule(Int_t nSupMod, Int_t iphi, Int_t ieta,
 815                         Int_t &iphim, Int_t &ietam, Int_t &nModule) const
 816 {
 817   // Transition from cell indexes (ieta,iphi) to module indexes (ietam,iphim, nModule)
 818   static Int_t nphi;
 819   nphi  = GetNumberOfModuleInPhiDirection(nSupMod);
 820
 821   ietam  = ieta/fNETAdiv;
 822   iphim  = iphi/fNPHIdiv;
 823   nModule = ietam * nphi + iphim;
 824 }
 825
 826 Int_t  AliEMCALGeometry::GetAbsCellIdFromCellIndexes(Int_t nSupMod, Int_t iphi, Int_t ieta) const
 827 {
 828   // Transition from super module number(nSupMod) and cell indexes (ieta,iphi) to absId
 829   static Int_t ietam, iphim, nModule;
 830   static Int_t nIeta, nIphi; // cell indexes in module
 831
 832   GetModuleIndexesFromCellIndexesInSModule(nSupMod, iphi, ieta, ietam, iphim, nModule);
 833
 834   nIeta = ieta%fNETAdiv;
 835   nIeta = fNETAdiv - 1 - nIeta;
 836   nIphi = iphi%fNPHIdiv;
 837
 838   return GetAbsCellId(nSupMod, nModule, nIphi, nIeta);
 839 }
 840
 841
 842 // Methods for AliEMCALRecPoint - Feb 19, 2006
 843 Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t &xr, Double_t &yr, Double_t &zr) const
 844 {
 845   // Look to see what the relative
 846   // position inside a given cell is
 847   // for a recpoint.
 848   // Alice numbering scheme - Jun 08, 2006
 849   // In:
 850   // absId   - cell is as in Geant,     0<= absId   < fNCells;
 851   // OUT:
 852   // xr,yr,zr - x,y,z coordinates of cell with absId inside SM
 853
 854   // Shift index taking into account the difference between standard SM
 855   // and SM of half size in phi direction
 856   const Int_t kphiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
 857   static Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta;
 858   if(!CheckAbsCellId(absId)) return kFALSE;
 859
 860   GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
 861   GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta);
 862
 863   xr = fCentersOfCellsXDir.At(ieta);
 864   zr = fCentersOfCellsEtaDir.At(ieta);
 865
 866   if(nSupMod<10) {
 867     yr = fCentersOfCellsPhiDir.At(iphi);
 868   } else {
 869     yr = fCentersOfCellsPhiDir.At(iphi + kphiIndexShift);
 870   }
 871   AliDebug(1,Form("absId %i nSupMod %i iphi %i ieta %i xr %f yr %f zr %f ",absId,nSupMod,iphi,ieta,xr,yr,zr));
 872
 873   return kTRUE;
 874 }
 875
 876 Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t loc[3]) const
 877 {
 878   // Alice numbering scheme - Jun 03, 2006
 879   loc[0] = loc[1] = loc[2]=0.0;
 880   if(RelPosCellInSModule(absId, loc[0],loc[1],loc[2])) {
 881     return kTRUE;
 882   }
 883   return kFALSE;
 884 }
 885
 886 Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, TVector3 &vloc) const
 887 {
 888   static Double_t loc[3];
 889   if(RelPosCellInSModule(absId,loc)) {
 890     vloc.SetXYZ(loc[0], loc[1], loc[2]);
 891     return kTRUE;
 892   } else {
 893     vloc.SetXYZ(0,0,0);
 894     return kFALSE;
 895   }
 896   // Alice numbering scheme - Jun 03, 2006
 897 }
 898
 899 Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t distEff, Double_t &xr, Double_t &yr, Double_t &zr) const
 900 {
 901   // Jul 30, 2007 - taking into account position of shower max
 902   // Look to see what the relative
 903   // position inside a given cell is
 904   // for a recpoint.
 905   // In:
 906   // absId   - cell is as in Geant,     0<= absId   < fNCells;
 907   // e       - cluster energy
 908   // OUT:
 909   // xr,yr,zr - x,y,z coordinates of cell with absId inside SM
 910
 911   // Shift index taking into account the difference between standard SM
 912   // and SM of half size in phi direction
 913   const  Int_t kphiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
 914   static Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta;
 915   static Int_t iphim, ietam;
 916   static AliEMCALShishKebabTrd1Module *mod = 0;
 917   static TVector2 v;
 918   if(!CheckAbsCellId(absId)) return kFALSE;
 919
 920   GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
 921   GetModulePhiEtaIndexInSModule(nSupMod, nModule, iphim, ietam);
 922   GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta);
 923
 924   mod = GetShishKebabModule(ietam);
 925   mod->GetPositionAtCenterCellLine(nIeta, distEff, v);
 926   xr = v.Y() - fParSM[0];
 927   zr = v.X() - fParSM[2];
 928
 929   if(nSupMod<10) {
 930     yr = fCentersOfCellsPhiDir.At(iphi);
 931   } else {
 932     yr = fCentersOfCellsPhiDir.At(iphi + kphiIndexShift);
 933   }
 934   AliDebug(1,Form("absId %i nSupMod %i iphi %i ieta %i xr %f yr %f zr %f ",absId,nSupMod,iphi,ieta,xr,yr,zr));
 935
 936   return kTRUE;
 937 }
 938
 939 Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Int_t maxAbsId, Double_t distEff, Double_t &xr, Double_t &yr, Double_t &zr) const
 940 {
 941   // Jul 31, 2007 - taking into account position of shower max and apply coor2.
 942   // Look to see what the relative
 943   // position inside a given cell is
 944   // for a recpoint.
 945   // In:
 946   // absId     - cell is as in Geant,     0<= absId   < fNCells;
 947   // maxAbsId  - abs id of cell with highest energy
 948   // e         - cluster energy
 949   // OUT:
 950   // xr,yr,zr - x,y,z coordinates of cell with absId inside SM
 951
 952   // Shift index taking into account the difference between standard SM
 953   // and SM of half size in phi direction
 954   const  Int_t kphiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
 955   static Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta;
 956   static Int_t iphim, ietam;
 957   static AliEMCALShishKebabTrd1Module *mod = 0;
 958   static TVector2 v;
 959
 960   static Int_t nSupModM, nModuleM, nIphiM, nIetaM, iphiM, ietaM;
 961   static Int_t iphimM, ietamM, maxAbsIdCopy=-1;
 962   static AliEMCALShishKebabTrd1Module *modM = 0;
 963   static Double_t distCorr;
 964
 965   if(!CheckAbsCellId(absId)) return kFALSE;
 966
 967   GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
 968   GetModulePhiEtaIndexInSModule(nSupMod, nModule, iphim, ietam);
 969   GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta);
 970   mod = GetShishKebabModule(ietam);
 971
 972   if(absId != maxAbsId) {
 973     distCorr = 0.;
 974     if(maxAbsIdCopy != maxAbsId) {
 975       GetCellIndex(maxAbsId, nSupModM, nModuleM, nIphiM, nIetaM);
 976       GetModulePhiEtaIndexInSModule(nSupModM, nModuleM, iphimM, ietamM);
 977       GetCellPhiEtaIndexInSModule(nSupModM,nModuleM,nIphiM,nIetaM, iphiM, ietaM);
 978       modM = GetShishKebabModule(ietamM); // do I need this ?
 979       maxAbsIdCopy = maxAbsId;
 980     }
 981     if(ietamM !=0) {
 982       distCorr = GetEtaModuleSize()*(ietam-ietamM)/TMath::Tan(modM->GetTheta()); // Stay here
 983       //printf(" distCorr %f | dist %f | ietam %i -> etamM %i\n", distCorr, dist, ietam, ietamM);
 984     }
 985     // distEff += distCorr;
 986   }
 987   // Bad resolution in this case, strong bias vs phi
 988   // distEff = 0.0;
 989   mod->GetPositionAtCenterCellLine(nIeta, distEff, v); // Stay here
 990   xr = v.Y() - fParSM[0];
 991   zr = v.X() - fParSM[2];
 992
 993   if(nSupMod<10) {
 994     yr = fCentersOfCellsPhiDir.At(iphi);
 995   } else {
 996     yr = fCentersOfCellsPhiDir.At(iphi + kphiIndexShift);
 997   }
 998   AliDebug(1,Form("absId %i nSupMod %i iphi %i ieta %i xr %f yr %f zr %f ",absId,nSupMod,iphi,ieta,xr,yr,zr));
 999
1000   return kTRUE;
1001 }
1002
1003 void AliEMCALGeometry::CreateListOfTrd1Modules()
1004 {
1005   // Generate the list of Trd1 modules
1006   // which will make up the EMCAL
1007   // geometry
1008
1009   AliDebug(2,Form(" AliEMCALGeometry::CreateListOfTrd1Modules() started "));
1010
1011   AliEMCALShishKebabTrd1Module *mod=0, *mTmp=0; // current module
1012   if(fShishKebabTrd1Modules == 0) {
1013     fShishKebabTrd1Modules = new TList;
1014     fShishKebabTrd1Modules->SetName("ListOfTRD1");
1015     for(int iz=0; iz< GetNZ(); iz++) {
1016       if(iz==0) {
1017         mod  = new AliEMCALShishKebabTrd1Module(TMath::Pi()/2.,this);
1018       } else {
1019         mTmp  = new AliEMCALShishKebabTrd1Module(*mod);
1020         mod   = mTmp;
1021       }
1022       fShishKebabTrd1Modules->Add(mod);
1023     }
1024   } else {
1025     AliDebug(2,Form(" Already exits : "));
1026   }
1027   mod = (AliEMCALShishKebabTrd1Module*)fShishKebabTrd1Modules->At(fShishKebabTrd1Modules->GetSize()-1);
1028   fEtaMaxOfTRD1 = mod->GetMaxEtaOfModule(0);
1029
1030   AliDebug(2,Form(" fShishKebabTrd1Modules has %i modules : max eta %5.4f \n",
1031                   fShishKebabTrd1Modules->GetSize(),fEtaMaxOfTRD1));
1032   // Feb 20,2006;
1033   // Jun 01, 2006 - ALICE numbering scheme
1034   // define grid for cells in eta(z) and x directions in local coordinates system of SM
1035   // Works just for 2x2 case only -- ?? start here
1036   //
1037   //
1038   // Define grid for cells in phi(y) direction in local coordinates system of SM
1039   // as for 2X2 as for 3X3 - Nov 8,2006
1040   //
1041   AliDebug(2,Form(" Cells grid in phi directions : size %i\n", fCentersOfCellsPhiDir.GetSize()));
1042   Int_t ind=0; // this is phi index
1043   Int_t ieta=0, nModule=0, iphiTemp;
1044   Double_t xr, zr, theta, phi, eta, r, x,y;
1045   TVector3 vglob;
1046   Double_t ytCenterModule=0.0, ytCenterCell=0.0;
1047
1048   fCentersOfCellsPhiDir.Set(fNPhi*fNPHIdiv);
1049   fPhiCentersOfCells.Set(fNPhi*fNPHIdiv);
1050
1051   Double_t r0 = GetIPDistance() + GetLongModuleSize()/2.;
1052   for(Int_t it=0; it<fNPhi; it++) { // cycle on modules
1053     ytCenterModule = -fParSM[1] + fPhiModuleSize*(2*it+1)/2;  // center of module
1054     for(Int_t ic=0; ic<fNPHIdiv; ic++) { // cycle on cells in module
1055       if(fNPHIdiv==2) {
1056         ytCenterCell = ytCenterModule + fPhiTileSize *(2*ic-1)/2.;
1057       } else if(fNPHIdiv==3){
1058         ytCenterCell = ytCenterModule + fPhiTileSize *(ic-1);
1059       } else if(fNPHIdiv==1){
1060         ytCenterCell = ytCenterModule;
1061       }
1062       fCentersOfCellsPhiDir.AddAt(ytCenterCell,ind);
1063       // Define grid on phi direction
1064       // Grid is not the same for different eta bin;
1065       // Effect is small but is still here
1066       phi = TMath::ATan2(ytCenterCell, r0);
1067       fPhiCentersOfCells.AddAt(phi, ind);
1068
1069       AliDebug(2,Form(" ind %2.2i : y %8.3f ", ind, fCentersOfCellsPhiDir.At(ind)));
1070       ind++;
1071     }
1072   }
1073
1074   fCentersOfCellsEtaDir.Set(fNZ *fNETAdiv);
1075   fCentersOfCellsXDir.Set(fNZ *fNETAdiv);
1076   fEtaCentersOfCells.Set(fNZ *fNETAdiv * fNPhi*fNPHIdiv);
1077   AliDebug(2,Form(" Cells grid in eta directions : size %i\n", fCentersOfCellsEtaDir.GetSize()));
1078   for(Int_t it=0; it<fNZ; it++) {
1079     AliEMCALShishKebabTrd1Module *trd1 = GetShishKebabModule(it);
1080     nModule = fNPhi*it;
1081     for(Int_t ic=0; ic<fNETAdiv; ic++) {
1082       if(fNPHIdiv==2) {
1083         trd1->GetCenterOfCellInLocalCoordinateofSM(ic, xr, zr);      // case of 2X2
1084         GetCellPhiEtaIndexInSModule(0, nModule, 0, ic, iphiTemp, ieta);
1085       } if(fNPHIdiv==3) {
1086         trd1->GetCenterOfCellInLocalCoordinateofSM_3X3(ic, xr, zr);  // case of 3X3
1087         GetCellPhiEtaIndexInSModule(0, nModule, 0, ic, iphiTemp, ieta);
1088       } if(fNPHIdiv==1) {
1089         trd1->GetCenterOfCellInLocalCoordinateofSM_1X1(xr, zr);      // case of 1X1
1090         GetCellPhiEtaIndexInSModule(0, nModule, 0, ic, iphiTemp, ieta);
1091       }
1092       fCentersOfCellsXDir.AddAt(float(xr) - fParSM[0],ieta);
1093       fCentersOfCellsEtaDir.AddAt(float(zr) - fParSM[2],ieta);
1094       // Define grid on eta direction for each bin in phi
1095       for(int iphi=0; iphi<fCentersOfCellsPhiDir.GetSize(); iphi++) {
1096         x = xr + trd1->GetRadius();
1097         y = fCentersOfCellsPhiDir[iphi];
1098         r = TMath::Sqrt(x*x + y*y + zr*zr);
1099         theta = TMath::ACos(zr/r);
1100         eta   = AliEMCALShishKebabTrd1Module::ThetaToEta(theta);
1101         //        ind   = ieta*fCentersOfCellsPhiDir.GetSize() + iphi;
1102         ind   = iphi*fCentersOfCellsEtaDir.GetSize() + ieta;
1103         fEtaCentersOfCells.AddAt(eta, ind);
1104       }
1105       //printf(" ieta %i : xr + trd1->GetRadius() %f : zr %f : eta %f \n", ieta, xr + trd1->GetRadius(), zr, eta);
1106     }
1107   }
1108   for(Int_t i=0; i<fCentersOfCellsEtaDir.GetSize(); i++) {
1109     AliDebug(2,Form(" ind %2.2i : z %8.3f : x %8.3f", i+1,
1110                     fCentersOfCellsEtaDir.At(i),fCentersOfCellsXDir.At(i)));
1111   }
1112
1113 }
1114
1115 void AliEMCALGeometry::GetGlobal(const Double_t *loc, Double_t *glob, int ind) const
1116 {
1117   // Figure out the global numbering
1118   // of a given supermodule from the
1119   // local numbering and the transformation
1120   // matrix stored by the geometry manager (allows for misaligned
1121   // geometry)
1122
1123   if(ind>=0 && ind < GetNumberOfSuperModules()) {
1124     TString volpath = "ALIC_1/XEN1_1/SMOD_";
1125     volpath += ind+1;
1126
1127     if(GetKey110DEG() && ind>=10) {
1128       volpath = "ALIC_1/XEN1_1/SM10_";
1129       volpath += ind-10+1;
1130     }
1131
1132     if(!gGeoManager->cd(volpath.Data()))
1133       AliFatal(Form("AliEMCALGeometry::GeoManager cannot find path %s!",volpath.Data()));
1134
1135     TGeoHMatrix* m = gGeoManager->GetCurrentMatrix();
1136     if(m) {
1137       m->LocalToMaster(loc, glob);
1138     } else {
1139       AliFatal("Geo matrixes are not loaded \n") ;
1140     }
1141   }
1142 }
1143
1144 void AliEMCALGeometry::GetGlobal(const TVector3 &vloc, TVector3 &vglob, int ind) const
1145 {
1146   //Figure out the global numbering
1147   //of a given supermodule from the
1148   //local numbering given a 3-vector location
1149
1150   static Double_t tglob[3], tloc[3];
1151   vloc.GetXYZ(tloc);
1152   GetGlobal(tloc, tglob, ind);
1153   vglob.SetXYZ(tglob[0], tglob[1], tglob[2]);
1154 }
1155
1156 void AliEMCALGeometry::GetGlobal(Int_t absId , double glob[3]) const
1157 {
1158   // Alice numbering scheme - Jun 03, 2006
1159   static Int_t nSupMod, nModule, nIphi, nIeta;
1160   static double loc[3];
1161
1162   if (!gGeoManager || !gGeoManager->IsClosed()) {
1163     AliError("Can't get the global coordinates! gGeoManager doesn't exist or it is still open!");
1164     return;
1165   }
1166
1167   glob[0]=glob[1]=glob[2]=0.0; // bad case
1168   if(RelPosCellInSModule(absId, loc)) {
1169     GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
1170
1171     TString volpath = "ALIC_1/XEN1_1/SMOD_";
1172     volpath += (nSupMod+1);
1173
1174     if(GetKey110DEG() && nSupMod>=10) {
1175       volpath = "ALIC_1/XEN1_1/SM10_";
1176       volpath += (nSupMod-10+1);
1177     }
1178     if(!gGeoManager->cd(volpath.Data()))
1179       AliFatal(Form("GeoManager cannot find path %s!",volpath.Data()));
1180
1181     TGeoHMatrix* m = gGeoManager->GetCurrentMatrix();
1182     if(m) {
1183       m->LocalToMaster(loc, glob);
1184     } else {
1185       AliFatal("Geo matrixes are not loaded \n") ;
1186     }
1187   }
1188 }
1189
1190 //___________________________________________________________________
1191 void AliEMCALGeometry::GetGlobal(Int_t absId , TVector3 &vglob) const
1192 {
1193   // Alice numbering scheme - Jun 03, 2006
1194   static Double_t glob[3];
1195
1196   GetGlobal(absId, glob);
1197   vglob.SetXYZ(glob[0], glob[1], glob[2]);
1198
1199 }
1200
1201 //____________________________________________________________________________
1202 void AliEMCALGeometry::GetGlobal(const AliRecPoint* /*rp*/, TVector3& /* vglob */) const
1203 {
1204   AliFatal(Form("Please use GetGlobalEMCAL(recPoint,gpos) instead of GetGlobal!"));
1205 }
1206
1207 //_________________________________________________________________________________
1208 void AliEMCALGeometry::GetGlobalEMCAL(const AliEMCALRecPoint *rp, TVector3 &vglob) const
1209 {
1210   // Figure out the global numbering
1211   // of a given supermodule from the
1212   // local numbering for RecPoints
1213
1214   static TVector3 vloc;
1215   static Int_t nSupMod, nModule, nIphi, nIeta;
1216
1217   const AliEMCALRecPoint *rpTmp = rp;
1218   const AliEMCALRecPoint *rpEmc = rpTmp;
1219
1220   GetCellIndex(rpEmc->GetAbsId(0), nSupMod, nModule, nIphi, nIeta);
1221   rpTmp->GetLocalPosition(vloc);
1222   GetGlobal(vloc, vglob, nSupMod);
1223 }
1224
1225 void AliEMCALGeometry::EtaPhiFromIndex(Int_t absId,Double_t &eta,Double_t &phi) const
1226 {
1227   // Nov 16, 2006- float to double
1228   // version for TRD1 only
1229   static TVector3 vglob;
1230   GetGlobal(absId, vglob);
1231   eta = vglob.Eta();
1232   phi = vglob.Phi();
1233 }
1234
1235 void AliEMCALGeometry::EtaPhiFromIndex(Int_t absId,Float_t &eta,Float_t &phi) const
1236 {
1237   // Nov 16,2006 - should be discard in future
1238   static TVector3 vglob;
1239   GetGlobal(absId, vglob);
1240   eta = float(vglob.Eta());
1241   phi = float(vglob.Phi());
1242 }
1243
1244 Bool_t AliEMCALGeometry::GetPhiBoundariesOfSM(Int_t nSupMod, Double_t &phiMin, Double_t &phiMax) const
1245 {
1246   // 0<= nSupMod <=11; phi in rad
1247   static int i;
1248   if(nSupMod<0 || nSupMod >11) return kFALSE;
1249   i = nSupMod/2;
1250   phiMin = fPhiBoundariesOfSM[2*i];
1251   phiMax = fPhiBoundariesOfSM[2*i+1];
1252   return kTRUE;
1253 }
1254
1255 Bool_t AliEMCALGeometry::GetPhiBoundariesOfSMGap(Int_t nPhiSec, Double_t &phiMin, Double_t &phiMax) const
1256 {
1257   // 0<= nPhiSec <=4; phi in rad
1258   // 0;  gap boundaries between  0th&2th  | 1th&3th SM
1259   // 1;  gap boundaries between  2th&4th  | 3th&5th SM
1260   // 2;  gap boundaries between  4th&6th  | 5th&7th SM
1261   // 3;  gap boundaries between  6th&8th  | 7th&9th SM
1262   // 4;  gap boundaries between  8th&10th | 9th&11th SM
1263   if(nPhiSec<0 || nPhiSec >4) return kFALSE;
1264   phiMin = fPhiBoundariesOfSM[2*nPhiSec+1];
1265   phiMax = fPhiBoundariesOfSM[2*nPhiSec+2];
1266   return kTRUE;
1267 }
1268
1269 Bool_t AliEMCALGeometry::SuperModuleNumberFromEtaPhi(Double_t eta, Double_t phi, Int_t &nSupMod) const
1270 {
1271   // Return false if phi belongs a phi cracks between SM
1272
1273   static Int_t i;
1274
1275   if(TMath::Abs(eta) > fEtaMaxOfTRD1) return kFALSE;
1276
1277   phi = TVector2::Phi_0_2pi(phi); // move phi to (0,2pi) boundaries
1278   for(i=0; i<6; i++) {
1279     if(phi>=fPhiBoundariesOfSM[2*i] && phi<=fPhiBoundariesOfSM[2*i+1]) {
1280       nSupMod = 2*i;
1281       if(eta < 0.0) nSupMod++;
1282       AliDebug(1,Form("eta %f phi %f(%5.2f) : nSupMod %i : #bound %i", eta,phi,phi*TMath::RadToDeg(), nSupMod,i));
1283       return kTRUE;
1284     }
1285   }
1286   return kFALSE;
1287 }
1288
1289 Bool_t AliEMCALGeometry::GetAbsCellIdFromEtaPhi(Double_t eta, Double_t phi, Int_t &absId) const
1290 {
1291   // Nov 17,2006
1292   // stay here - phi problem as usual
1293   static Int_t nSupMod, i, ieta, iphi, etaShift, nphi;
1294   static Double_t absEta=0.0, d=0.0, dmin=0.0, phiLoc;
1295   absId = nSupMod = - 1;
1296   if(SuperModuleNumberFromEtaPhi(eta, phi, nSupMod)) {
1297     // phi index first
1298     phi    = TVector2::Phi_0_2pi(phi);
1299     phiLoc = phi - fPhiCentersOfSM[nSupMod/2];
1300     nphi   = fPhiCentersOfCells.GetSize();
1301     if(nSupMod>=10) {
1302       phiLoc = phi - 190.*TMath::DegToRad();
1303       nphi  /= 2;
1304     }
1305
1306     dmin   = TMath::Abs(fPhiCentersOfCells[0]-phiLoc);
1307     iphi   = 0;
1308     for(i=1; i<nphi; i++) {
1309       d = TMath::Abs(fPhiCentersOfCells[i] - phiLoc);
1310       if(d < dmin) {
1311         dmin = d;
1312         iphi = i;
1313       }
1314       //      printf(" i %i : d %f : dmin %f : fPhiCentersOfCells[i] %f \n", i, d, dmin, fPhiCentersOfCells[i]);
1315     }
1316     // odd SM are turned with respect of even SM - reverse indexes
1317     AliDebug(2,Form(" iphi %i : dmin %f (phi %f, phiLoc %f ) ", iphi, dmin, phi, phiLoc));
1318     // eta index
1319     absEta   = TMath::Abs(eta);
1320     etaShift = iphi*fCentersOfCellsEtaDir.GetSize();
1321     dmin     = TMath::Abs(fEtaCentersOfCells[etaShift]-absEta);
1322     ieta     = 0;
1323     for(i=1; i<fCentersOfCellsEtaDir.GetSize(); i++) {
1324       d = TMath::Abs(fEtaCentersOfCells[i+etaShift] - absEta);
1325       if(d < dmin) {
1326         dmin = d;
1327         ieta = i;
1328       }
1329     }
1330     AliDebug(2,Form(" ieta %i : dmin %f (eta=%f) : nSupMod %i ", ieta, dmin, eta, nSupMod));
1331
1332     if(eta<0) iphi = (nphi-1) - iphi;
1333     absId = GetAbsCellIdFromCellIndexes(nSupMod, iphi, ieta);
1334
1335     return kTRUE;
1336   }
1337   return kFALSE;
1338 }
1339
1340 AliEMCALShishKebabTrd1Module* AliEMCALGeometry::GetShishKebabModule(Int_t neta) const
1341 {
1342   //This method was too long to be
1343   //included in the header file - the
1344   //rule checker complained about it's
1345   //length, so we move it here.  It returns the
1346   //shishkebabmodule at a given eta index point.
1347
1348   static AliEMCALShishKebabTrd1Module* trd1=0;
1349   if(fShishKebabTrd1Modules && neta>=0 && neta<fShishKebabTrd1Modules->GetSize()) {
1350     trd1 = (AliEMCALShishKebabTrd1Module*)fShishKebabTrd1Modules->At(neta);
1351   } else trd1 = 0;
1352   return trd1;
1353 }
1354
1355 void AliEMCALGeometry::Browse(TBrowser* b) const
1356 {
1357   //Browse the modules
1358   if(fShishKebabTrd1Modules) b->Add(fShishKebabTrd1Modules);
1359 }
1360
1361 Bool_t AliEMCALGeometry::IsFolder() const
1362 {
1363   //Check if fShishKebabTrd1Modules is in folder
1364   if(fShishKebabTrd1Modules) return kTRUE;
1365   else                       return kFALSE;
1366 }
1367
1368 Double_t AliEMCALGeometry::GetPhiCenterOfSM(Int_t nsupmod) const
1369 {
1370   //returns center of supermodule in phi
1371   int i = nsupmod/2;
1372   return fPhiCentersOfSM[i];
1373
1374 }