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