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