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