]> git.uio.no Git - u/mrichter/AliRoot.git/blame - EMCAL/AliEMCALGeometry.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Additional protection in the destructor: when you have a chain of calls returning...
[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;
d297ef6e 82Char_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// == =============================
95// AliRunLoader *rl = AliRunLoader::GetRunLoader();
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
273 if(!(fGeoName.Contains("EMCAL_PDC06") || fGeoName.Contains("EMCAL_COMPLETE") || fGeoName.Contains("EMCAL_WSUC"))) {
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
d297ef6e 334 // constant for transition absid <--> indexes
335 fNCellsInModule = fNPHIdiv*fNETAdiv;
336 fNCellsInSupMod = fNCellsInModule*fNPhi*fNZ;
337 fNCells = fNCellsInSupMod*fNumberOfSuperModules;
338 if(GetKey110DEG()) fNCells -= fNCellsInSupMod;
937d0661 339
1963b290 340 fNPhiSuperModule = fNumberOfSuperModules/2;
341 if(fNPhiSuperModule<1) fNPhiSuperModule = 1;
d297ef6e 342
343 fPhiTileSize = fPhiModuleSize/double(fNPHIdiv) - fLateralSteelStrip; // 13-may-05
344 fEtaTileSize = fEtaModuleSize/double(fNETAdiv) - fLateralSteelStrip; // 13-may-05
1d46d1f6 345
d297ef6e 346 fLongModuleSize = fNECLayers*(fECScintThick + fECPbRadThickness);
347 f2Trd1Dx2 = fEtaModuleSize + 2.*fLongModuleSize*TMath::Tan(fTrd1Angle*TMath::DegToRad()/2.);
348 if(!fGeoName.Contains("WSUC")) fShellThickness = TMath::Sqrt(fLongModuleSize*fLongModuleSize + f2Trd1Dx2*f2Trd1Dx2);
facb5967 349
350 //These parameters are used to create the mother volume to hold the supermodules
351 //2cm padding added to allow for misalignments - JLK 30-May-2008
352 fEnvelop[0] = fIPDistance - 1.; // mother volume inner radius
353 fEnvelop[1] = fIPDistance + fShellThickness + 1.; // mother volume outer r.
354 fEnvelop[2] = fZLength + 2.; //mother volume length
9cff4509 355
d297ef6e 356 // Local coordinates
357 fParSM[0] = GetShellThickness()/2.;
358 fParSM[1] = GetPhiModuleSize() * GetNPhi()/2.;
facb5967 359 fParSM[2] = fZLength/4.; //divide by 4 to get half-length of SM
1d46d1f6 360
361 // SM phi boundaries - (0,1),(2,3) .. (10,11) - has the same boundaries; Nov 7, 2006
362 fPhiBoundariesOfSM.Set(fNumberOfSuperModules);
363 fPhiCentersOfSM.Set(fNumberOfSuperModules/2);
364 fPhiBoundariesOfSM[0] = TMath::PiOver2() - TMath::ATan2(fParSM[1] , fIPDistance); // 1th and 2th modules)
1d46d1f6 365 fPhiCentersOfSM[0] = TMath::PiOver2();
d297ef6e 366 if(fNumberOfSuperModules > 1)
367 fPhiBoundariesOfSM[1] = TMath::PiOver2() + TMath::ATan2(fParSM[1] , fIPDistance);
368 if(fNumberOfSuperModules > 2) {
369 for(int i=1; i<=4; i++) { // from 2th ro 9th
370 fPhiBoundariesOfSM[2*i] = fPhiBoundariesOfSM[0] + 20.*TMath::DegToRad()*i;
371 fPhiBoundariesOfSM[2*i+1] = fPhiBoundariesOfSM[1] + 20.*TMath::DegToRad()*i;
372 fPhiCentersOfSM[i] = fPhiCentersOfSM[0] + 20.*TMath::DegToRad()*i;
373 }
374 }
375 if(fNumberOfSuperModules > 10) {
376 fPhiBoundariesOfSM[11] = 190.*TMath::DegToRad();
377 fPhiBoundariesOfSM[10] = fPhiBoundariesOfSM[11] - TMath::ATan2((fParSM[1]) , fIPDistance);
378 fPhiCentersOfSM[5] = (fPhiBoundariesOfSM[10]+fPhiBoundariesOfSM[11])/2.;
1d46d1f6 379 }
d297ef6e 380
381 //called after setting of scintillator and lead layer parameters
382 DefineSamplingFraction();
1d46d1f6 383
85c25c2e 384 // TRU parameters - Apr 29,08 by PAI.
385 // These parameters values was updated at Nov 05, 2007
386 // As is on Olivier BOURRION (LPSC) ppt preasentation
387 // at ALICE trigger meeting at 13th-14th March
388 fNTRUEta = 1; // was 3
389 fNTRUPhi = 3; // was 1
390 fNModulesInTRUEta = 24; // was 8
391 fNModulesInTRUPhi = 4; // was 12
392 // Jet trigger
393 // 3*6*10 + 2*6*2 = 204 -> matrix (nphi(17), neta(12))
394 fNEtaSubOfTRU = 6;
1d46d1f6 395
89557f6d 396 fgInit = kTRUE;
2012850d 397}
173558f2 398
1d46d1f6 399void AliEMCALGeometry::PrintGeometry()
400{
401 // Separate routine is callable from broswer; Nov 7,2006
89557f6d 402 printf("\nInit: geometry of EMCAL named %s :\n", fGeoName.Data());
403 if(fArrayOpts) {
404 for(Int_t i=0; i<fArrayOpts->GetEntries(); i++){
405 TObjString *o = (TObjString*)fArrayOpts->At(i);
406 printf(" %i : %s \n", i, o->String().Data());
407 }
408 }
1d46d1f6 409 printf("Granularity: %d in eta and %d in phi\n", GetNZ(), GetNPhi()) ;
410 printf("Layout: phi = (%7.1f, %7.1f), eta = (%5.2f, %5.2f), IP = %7.2f -> for EMCAL envelope only\n",
411 GetArm1PhiMin(), GetArm1PhiMax(),GetArm1EtaMin(), GetArm1EtaMax(), GetIPDistance() );
c63c3c5d 412
1d46d1f6 413 printf( " ECAL : %d x (%f cm Pb, %f cm Sc) \n",
414 GetNECLayers(), GetECPbRadThick(), GetECScintThick() ) ;
415 printf(" fSampling %5.2f \n", fSampling );
d297ef6e 416 printf(" fIPDistance %6.3f cm \n", fIPDistance);
417 printf(" fNPhi %i | fNZ %i \n", fNPhi, fNZ);
418 printf(" fNCellsInModule %i : fNCellsInSupMod %i : fNCells %i\n",fNCellsInModule, fNCellsInSupMod, fNCells);
419 printf(" X:Y module size %6.3f , %6.3f cm \n", fPhiModuleSize, fEtaModuleSize);
420 printf(" X:Y tile size %6.3f , %6.3f cm \n", fPhiTileSize, fEtaTileSize);
421 printf(" #of sampling layers %i(fNECLayers) \n", fNECLayers);
422 printf(" fLongModuleSize %6.3f cm \n", fLongModuleSize);
423 printf(" #supermodule in phi direction %i \n", fNPhiSuperModule );
89557f6d 424 printf(" fILOSS %i : fIHADR %i \n", fILOSS, fIHADR);
d297ef6e 425 printf(" fTrd1Angle %7.4f\n", fTrd1Angle);
426 printf(" f2Trd1Dx2 %7.4f\n", f2Trd1Dx2);
427 printf("SM dimensions(TRD1) : dx %7.2f dy %7.2f dz %7.2f (SMOD, BOX)\n",
428 fParSM[0],fParSM[1],fParSM[2]);
429 printf(" fPhiGapForSM %7.4f cm (%7.4f <- phi size in degree)\n",
430 fPhiGapForSM, TMath::ATan2(fPhiGapForSM,fIPDistance)*TMath::RadToDeg());
431 if(GetKey110DEG()) printf(" Last two modules have size 10 degree in phi (180<phi<190)\n");
432 printf(" phi SM boundaries \n");
433 for(int i=0; i<fPhiBoundariesOfSM.GetSize()/2.; i++) {
434 printf(" %i : %7.5f(%7.2f) -> %7.5f(%7.2f) : center %7.5f(%7.2f) \n", i,
435 fPhiBoundariesOfSM[2*i], fPhiBoundariesOfSM[2*i]*TMath::RadToDeg(),
436 fPhiBoundariesOfSM[2*i+1], fPhiBoundariesOfSM[2*i+1]*TMath::RadToDeg(),
437 fPhiCentersOfSM[i], fPhiCentersOfSM[i]*TMath::RadToDeg());
438 }
439 printf(" fShishKebabTrd1Modules has %i modules : max eta %5.4f \n",
440 fShishKebabTrd1Modules->GetSize(),fEtaMaxOfTRD1);
441
442 printf("\n Cells grid in eta directions : size %i\n", fCentersOfCellsEtaDir.GetSize());
443 for(Int_t i=0; i<fCentersOfCellsEtaDir.GetSize(); i++) {
444 printf(" ind %2.2i : z %8.3f : x %8.3f \n", i,
445 fCentersOfCellsEtaDir.At(i),fCentersOfCellsXDir.At(i));
446 int ind=0; // Nov 21,2006
447 for(Int_t iphi=0; iphi<fCentersOfCellsPhiDir.GetSize(); iphi++) {
448 ind = iphi*fCentersOfCellsEtaDir.GetSize() + i;
449 printf("%6.4f ", fEtaCentersOfCells[ind]);
450 if((iphi+1)%12 == 0) printf("\n");
1d46d1f6 451 }
d297ef6e 452 printf("\n");
453
1d46d1f6 454 }
d297ef6e 455
456 printf("\n Cells grid in phi directions : size %i\n", fCentersOfCellsPhiDir.GetSize());
457 for(Int_t i=0; i<fCentersOfCellsPhiDir.GetSize(); i++) {
458 double phi=fPhiCentersOfCells.At(i);
459 printf(" ind %2.2i : y %8.3f : phi %7.5f(%6.2f) \n", i, fCentersOfCellsPhiDir.At(i),
460 phi, phi*TMath::RadToDeg());
461 }
462
1d46d1f6 463}
464
d297ef6e 465//______________________________________________________________________
1d46d1f6 466void AliEMCALGeometry::PrintCellIndexes(Int_t absId, int pri, char *tit)
467{
468 // Service methods
2bb3725c 469 Int_t nSupMod, nModule, nIphi, nIeta;
1d46d1f6 470 Int_t iphi, ieta;
471 TVector3 vg;
472
2bb3725c 473 GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
474 printf(" %s | absId : %i -> nSupMod %i nModule %i nIphi %i nIeta %i \n", tit, absId, nSupMod, nModule, nIphi, nIeta);
1d46d1f6 475 if(pri>0) {
2bb3725c 476 GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi,ieta);
1d46d1f6 477 printf(" local SM index : iphi %i : ieta %i \n", iphi,ieta);
478 GetGlobal(absId, vg);
479 printf(" vglob : mag %7.2f : perp %7.2f : z %7.2f : eta %6.4f : phi %6.4f(%6.2f) \n",
480 vg.Mag(), vg.Perp(), vg.Z(), vg.Eta(), vg.Phi(), vg.Phi()*TMath::RadToDeg());
481 }
482}
483
484//______________________________________________________________________
fc575e27 485void AliEMCALGeometry::CheckAdditionalOptions()
486{
487 // Feb 06,2006
89557f6d 488 // Additional options that
489 // can be used to select
490 // the specific geometry of
491 // EMCAL to run
492 // Dec 27,2006
493 // adeed allILOSS= and allIHADR= for MIP investigation
c63c3c5d 494 fArrayOpts = new TObjArray;
fc575e27 495 Int_t nopt = AliEMCALHistoUtilities::ParseString(fGeoName, *fArrayOpts);
c63c3c5d 496 if(nopt==1) { // no aditional option(s)
497 fArrayOpts->Delete();
498 delete fArrayOpts;
499 fArrayOpts = 0;
500 return;
501 }
502 for(Int_t i=1; i<nopt; i++){
503 TObjString *o = (TObjString*)fArrayOpts->At(i);
504
505 TString addOpt = o->String();
506 Int_t indj=-1;
fc575e27 507 for(Int_t j=0; j<fNAdditionalOpts; j++) {
508 TString opt = fAdditionalOpts[j];
c63c3c5d 509 if(addOpt.Contains(opt,TString::kIgnoreCase)) {
510 indj = j;
511 break;
512 }
513 }
514 if(indj<0) {
e5a93224 515 AliDebug(2,Form("<E> option |%s| unavailable : ** look to the file AliEMCALGeometry.h **\n",
516 addOpt.Data()));
c63c3c5d 517 assert(0);
518 } else {
e5a93224 519 AliDebug(2,Form("<I> option |%s| is valid : number %i : |%s|\n",
520 addOpt.Data(), indj, fAdditionalOpts[indj]));
c63c3c5d 521 if (addOpt.Contains("NL=",TString::kIgnoreCase)) {// number of sampling layers
522 sscanf(addOpt.Data(),"NL=%i", &fNECLayers);
e5a93224 523 AliDebug(2,Form(" fNECLayers %i (new) \n", fNECLayers));
b44d5aa4 524 } else if(addOpt.Contains("PBTH=",TString::kIgnoreCase)) {//Thickness of the Pb(fECPbRadThicknes)
c63c3c5d 525 sscanf(addOpt.Data(),"PBTH=%f", &fECPbRadThickness);
b44d5aa4 526 } else if(addOpt.Contains("SCTH=",TString::kIgnoreCase)) {//Thickness of the Sc(fECScintThick)
c63c3c5d 527 sscanf(addOpt.Data(),"SCTH=%f", &fECScintThick);
b44d5aa4 528 } else if(addOpt.Contains("LATSS=",TString::kIgnoreCase)) {// Thickness of lateral steel strip (fLateralSteelStrip)
529 sscanf(addOpt.Data(),"LATSS=%f", &fLateralSteelStrip);
530 AliDebug(2,Form(" fLateralSteelStrip %f (new) \n", fLateralSteelStrip));
89557f6d 531 } else if(addOpt.Contains("ILOSS=",TString::kIgnoreCase)) {// As in Geant
532 sscanf(addOpt.Data(),"ALLILOSS=%i", &fILOSS);
533 AliDebug(2,Form(" fILOSS %i \n", fILOSS));
534 } else if(addOpt.Contains("IHADR=",TString::kIgnoreCase)) {// As in Geant
535 sscanf(addOpt.Data(),"ALLIHADR=%i", &fIHADR);
536 AliDebug(2,Form(" fIHADR %i \n", fIHADR));
c63c3c5d 537 }
538 }
539 }
540}
541
25b033cf 542void AliEMCALGeometry::DefineSamplingFraction()
543{
544 // Jun 05,2006
545 // Look http://rhic.physics.wayne.edu/~pavlinov/ALICE/SHISHKEBAB/RES/linearityAndResolutionForTRD1.html
546 // Keep for compatibilty
547 //
548 if(fNECLayers == 69) { // 10% layer reduction
549 fSampling = 12.55;
550 } else if(fNECLayers == 61) { // 20% layer reduction
551 fSampling = 12.80;
552 } else if(fNECLayers == 77) {
d297ef6e 553 if (fECScintThick>0.159 && fECScintThick<0.161) { // original sampling fraction, equal layers
554 fSampling = 12.327; // fECScintThick = fECPbRadThickness = 0.160;
555 } else if (fECScintThick>0.175 && fECScintThick<0.177) { // 10% Pb thicknes reduction
25b033cf 556 fSampling = 10.5; // fECScintThick = 0.176, fECPbRadThickness=0.144;
557 } else if(fECScintThick>0.191 && fECScintThick<0.193) { // 20% Pb thicknes reduction
558 fSampling = 8.93; // fECScintThick = 0.192, fECPbRadThickness=0.128;
559 }
d297ef6e 560
25b033cf 561 }
562}
563
356fd0a9 564//______________________________________________________________________
85c25c2e 565void AliEMCALGeometry::GetModulePhiEtaIndexInSModuleFromTRUIndex(const Int_t itru, const Int_t iphitru, const Int_t ietatru, Int_t &iphiSM, Int_t &ietaSM) const
356fd0a9 566{
567
85c25c2e 568 // This method transforms the (eta,phi) index of module in a
356fd0a9 569 // TRU matrix into Super Module (eta,phi) index.
570
33d0b833 571 // Calculate in which row and column where the TRU are
356fd0a9 572 // ordered in the SM
573
85c25c2e 574 Int_t col = itru/ fNTRUPhi ; // indexes of TRU in SM
33d0b833 575 Int_t row = itru - col*fNTRUPhi ;
356fd0a9 576
85c25c2e 577 iphiSM = fNModulesInTRUPhi*row + iphitru ;
578 ietaSM = fNModulesInTRUEta*col + ietatru ;
579 //printf(" GetModulePhiEtaIndexInSModuleFromTRUIndex : itru %2i iphitru %2i ietatru %2i iphiSM %2i ietaSM %2i \n",
580 // itru, iphitru, ietatru, iphiSM, ietaSM);
356fd0a9 581}
f0377b23 582
b13bbe81 583//______________________________________________________________________
584AliEMCALGeometry * AliEMCALGeometry::GetInstance(){
05a92d59 585 // Returns the pointer of the unique instance
586
e52475ed 587 AliEMCALGeometry * rv = static_cast<AliEMCALGeometry *>( fgGeom );
588 return rv;
2012850d 589}
173558f2 590
b13bbe81 591//______________________________________________________________________
592AliEMCALGeometry* AliEMCALGeometry::GetInstance(const Text_t* name,
593 const Text_t* title){
594 // Returns the pointer of the unique instance
595
596 AliEMCALGeometry * rv = 0;
597 if ( fgGeom == 0 ) {
89557f6d 598 if ( strcmp(name,"") == 0 ) { // get default geometry
599 fgGeom = new AliEMCALGeometry(fgDefaultGeometryName, title);
600 } else {
601 fgGeom = new AliEMCALGeometry(name, title);
602 } // end if strcmp(name,"")
603 if ( fgInit ) rv = (AliEMCALGeometry * ) fgGeom;
604 else {
605 rv = 0;
606 delete fgGeom;
607 fgGeom = 0;
608 } // end if fgInit
b13bbe81 609 }else{
e5a93224 610 if ( strcmp(fgGeom->GetName(), name) != 0) {
87926378 611 printf("\ncurrent geometry is %s : ", fgGeom->GetName());
612 printf(" you cannot call %s ",name);
b13bbe81 613 }else{
9859bfc0 614 rv = (AliEMCALGeometry *) fgGeom;
e52475ed 615 } // end
b13bbe81 616 } // end if fgGeom
617 return rv;
2012850d 618}
173558f2 619
d297ef6e 620//_____________________________________________________________________________
ab37d09c 621Bool_t AliEMCALGeometry::IsInEMCAL(Double_t x, Double_t y, Double_t z) const {
dc7da436 622 // Checks whether point is inside the EMCal volume, used in AliEMCALv*.cxx
ab37d09c 623 //
624 // Code uses cylindrical approximation made of inner radius (for speed)
625 //
626 // Points behind EMCAl, i.e. R > outer radius, but eta, phi in acceptance
627 // are considered to inside
628
629 Double_t r=sqrt(x*x+y*y);
630
631 if ( r > fEnvelop[0] ) {
632 Double_t theta;
633 theta = TMath::ATan2(r,z);
634 Double_t eta;
635 if(theta == 0)
636 eta = 9999;
637 else
638 eta = -TMath::Log(TMath::Tan(theta/2.));
639 if (eta < fArm1EtaMin || eta > fArm1EtaMax)
640 return 0;
641
642 Double_t phi = TMath::ATan2(y,x) * 180./TMath::Pi();
2038074b 643 if (phi < 0) phi += 360; // phi should go from 0 to 360 in this case
ab37d09c 644 if (phi > fArm1PhiMin && phi < fArm1PhiMax)
645 return 1;
646 }
647 return 0;
648}
1963b290 649
650//
651// == Shish-kebab cases ==
652//
d297ef6e 653//________________________________________________________________________________________________
2bb3725c 654Int_t AliEMCALGeometry::GetAbsCellId(Int_t nSupMod, Int_t nModule, Int_t nIphi, Int_t nIeta) const
dc7da436 655{
656 // 27-aug-04;
d87bd045 657 // corr. 21-sep-04;
658 // 13-oct-05; 110 degree case
dc7da436 659 // May 31, 2006; ALICE numbering scheme:
660 // 0 <= nSupMod < fNumberOfSuperModules
2bb3725c 661 // 0 <= nModule < fNPHI * fNZ ( fNPHI * fNZ/2 for fKey110DEG=1)
dc7da436 662 // 0 <= nIphi < fNPHIdiv
663 // 0 <= nIeta < fNETAdiv
664 // 0 <= absid < fNCells
665 static Int_t id=0; // have to change from 0 to fNCells-1
666 if(fKey110DEG == 1 && nSupMod >= 10) { // 110 degree case; last two supermodules
667 id = fNCellsInSupMod*10 + (fNCellsInSupMod/2)*(nSupMod-10);
d87bd045 668 } else {
dc7da436 669 id = fNCellsInSupMod*nSupMod;
d87bd045 670 }
2bb3725c 671 id += fNCellsInModule *nModule;
dc7da436 672 id += fNPHIdiv *nIphi;
1963b290 673 id += nIeta;
dc7da436 674 if(id<0 || id >= fNCells) {
500aeccc 675// printf(" wrong numerations !!\n");
676// printf(" id %6i(will be force to -1)\n", id);
677// printf(" fNCells %6i\n", fNCells);
678// printf(" nSupMod %6i\n", nSupMod);
2bb3725c 679// printf(" nModule %6i\n", nModule);
500aeccc 680// printf(" nIphi %6i\n", nIphi);
681// printf(" nIeta %6i\n", nIeta);
dc7da436 682 id = -TMath::Abs(id); // if negative something wrong
1963b290 683 }
684 return id;
685}
686
d297ef6e 687//________________________________________________________________________________________________
dc7da436 688Bool_t AliEMCALGeometry::CheckAbsCellId(Int_t absId) const
fc575e27 689{
dc7da436 690 // May 31, 2006; only trd1 now
691 if(absId<0 || absId >= fNCells) return kFALSE;
692 else return kTRUE;
1963b290 693}
694
d297ef6e 695//________________________________________________________________________________________________
2bb3725c 696Bool_t AliEMCALGeometry::GetCellIndex(Int_t absId,Int_t &nSupMod,Int_t &nModule,Int_t &nIphi,Int_t &nIeta) const
fc575e27 697{
dc7da436 698 // 21-sep-04; 19-oct-05;
699 // May 31, 2006; ALICE numbering scheme:
4bba84bd 700 //
701 // In:
702 // absId - cell is as in Geant, 0<= absId < fNCells;
703 // Out:
704 // nSupMod - super module(SM) number, 0<= nSupMod < fNumberOfSuperModules;
2bb3725c 705 // nModule - module number in SM, 0<= nModule < fNCellsInSupMod/fNCellsInSupMod or(/2) for tow last SM (10th and 11th);
4bba84bd 706 // nIphi - cell number in phi driection inside module; 0<= nIphi < fNPHIdiv;
707 // nIeta - cell number in eta driection inside module; 0<= nIeta < fNETAdiv;
708 //
d87bd045 709 static Int_t tmp=0, sm10=0;
dc7da436 710 if(!CheckAbsCellId(absId)) return kFALSE;
711
d87bd045 712 sm10 = fNCellsInSupMod*10;
dc7da436 713 if(fKey110DEG == 1 && absId >= sm10) { // 110 degree case; last two supermodules
714 nSupMod = (absId-sm10) / (fNCellsInSupMod/2) + 10;
715 tmp = (absId-sm10) % (fNCellsInSupMod/2);
d87bd045 716 } else {
dc7da436 717 nSupMod = absId / fNCellsInSupMod;
718 tmp = absId % fNCellsInSupMod;
d87bd045 719 }
1963b290 720
2bb3725c 721 nModule = tmp / fNCellsInModule;
722 tmp = tmp % fNCellsInModule;
dc7da436 723 nIphi = tmp / fNPHIdiv;
724 nIeta = tmp % fNPHIdiv;
1963b290 725
726 return kTRUE;
727}
728
d297ef6e 729//________________________________________________________________________________________________
2bb3725c 730void AliEMCALGeometry::GetModulePhiEtaIndexInSModule(Int_t nSupMod, Int_t nModule, int &iphim, int &ietam) const
fc575e27 731{
1d46d1f6 732 // added nSupMod; - 19-oct-05 !
dc7da436 733 // Alice numbering scheme - Jun 01,2006
1d46d1f6 734 // ietam, iphi - indexes of module in two dimensional grid of SM
735 // ietam - have to change from 0 to fNZ-1
736 // iphim - have to change from 0 to nphi-1 (fNPhi-1 or fNPhi/2-1)
d87bd045 737 static Int_t nphi;
738
dc7da436 739 if(fKey110DEG == 1 && nSupMod>=10) nphi = fNPhi/2;
d87bd045 740 else nphi = fNPhi;
741
2bb3725c 742 ietam = nModule/nphi;
743 iphim = nModule%nphi;
d87bd045 744}
745
d297ef6e 746//________________________________________________________________________________________________
2bb3725c 747void AliEMCALGeometry::GetCellPhiEtaIndexInSModule(Int_t nSupMod, Int_t nModule, Int_t nIphi, Int_t nIeta,
e52475ed 748int &iphi, int &ieta) const
fc575e27 749{
1d46d1f6 750 //
751 // Added nSupMod; Nov 25, 05
752 // Alice numbering scheme - Jun 01,2006
4bba84bd 753 // IN:
754 // nSupMod - super module(SM) number, 0<= nSupMod < fNumberOfSuperModules;
2bb3725c 755 // nModule - module number in SM, 0<= nModule < fNCellsInSupMod/fNCellsInSupMod or(/2) for tow last SM (10th and 11th);
4bba84bd 756 // nIphi - cell number in phi driection inside module; 0<= nIphi < fNPHIdiv;
757 // nIeta - cell number in eta driection inside module; 0<= nIeta < fNETAdiv;
758 //
759 // OUT:
1d46d1f6 760 // ieta, iphi - indexes of cell(tower) in two dimensional grid of SM
761 // ieta - have to change from 0 to (fNZ*fNETAdiv-1)
762 // iphi - have to change from 0 to (fNPhi*fNPHIdiv-1 or fNPhi*fNPHIdiv/2-1)
763 //
dc7da436 764 static Int_t iphim, ietam;
765
2bb3725c 766 GetModulePhiEtaIndexInSModule(nSupMod,nModule, iphim, ietam);
1d46d1f6 767 // ieta = ietam*fNETAdiv + (1-nIeta); // x(module) = -z(SM)
768 ieta = ietam*fNETAdiv + (fNETAdiv - 1 - nIeta); // x(module) = -z(SM)
dc7da436 769 iphi = iphim*fNPHIdiv + nIphi; // y(module) = y(SM)
1d46d1f6 770
771 if(iphi<0 || ieta<0)
2bb3725c 772 AliDebug(1,Form(" nSupMod %i nModule %i nIphi %i nIeta %i => ieta %i iphi %i\n",
773 nSupMod, nModule, nIphi, nIeta, ieta, iphi));
1963b290 774}
e52475ed 775
d297ef6e 776//________________________________________________________________________________________________
e52475ed 777Int_t AliEMCALGeometry::GetSuperModuleNumber(Int_t absId) const
778{
1d46d1f6 779 // Return the number of the supermodule given the absolute
780 // ALICE numbering id
fc575e27 781
2bb3725c 782 static Int_t nSupMod, nModule, nIphi, nIeta;
783 GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
e52475ed 784 return nSupMod;
785}
786
d297ef6e 787//________________________________________________________________________________________________
1d46d1f6 788void AliEMCALGeometry::GetModuleIndexesFromCellIndexesInSModule(Int_t nSupMod, Int_t iphi, Int_t ieta,
2bb3725c 789 Int_t &iphim, Int_t &ietam, Int_t &nModule) const
1d46d1f6 790{
2bb3725c 791 // Transition from cell indexes (ieta,iphi) to module indexes (ietam,iphim, nModule)
1d46d1f6 792 static Int_t nphi;
793 nphi = GetNumberOfModuleInPhiDirection(nSupMod);
794
795 ietam = ieta/fNETAdiv;
796 iphim = iphi/fNPHIdiv;
2bb3725c 797 nModule = ietam * nphi + iphim;
1d46d1f6 798}
799
d297ef6e 800//________________________________________________________________________________________________
1d46d1f6 801Int_t AliEMCALGeometry::GetAbsCellIdFromCellIndexes(Int_t nSupMod, Int_t iphi, Int_t ieta) const
802{
803 // Transition from super module number(nSupMod) and cell indexes (ieta,iphi) to absId
2bb3725c 804 static Int_t ietam, iphim, nModule;
1d46d1f6 805 static Int_t nIeta, nIphi; // cell indexes in module
806
2bb3725c 807 GetModuleIndexesFromCellIndexesInSModule(nSupMod, iphi, ieta, ietam, iphim, nModule);
1d46d1f6 808
809 nIeta = ieta%fNETAdiv;
810 nIeta = fNETAdiv - 1 - nIeta;
811 nIphi = iphi%fNPHIdiv;
812
2bb3725c 813 return GetAbsCellId(nSupMod, nModule, nIphi, nIeta);
1d46d1f6 814}
815
816
e52475ed 817// Methods for AliEMCALRecPoint - Feb 19, 2006
d297ef6e 818//________________________________________________________________________________________________
14e75ea7 819Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t &xr, Double_t &yr, Double_t &zr) const
e52475ed 820{
1933eff2 821 // Look to see what the relative
822 // position inside a given cell is
823 // for a recpoint.
824 // Alice numbering scheme - Jun 08, 2006
4bba84bd 825 // In:
826 // absId - cell is as in Geant, 0<= absId < fNCells;
827 // OUT:
828 // xr,yr,zr - x,y,z coordinates of cell with absId inside SM
fc575e27 829
d25f2c54 830 // Shift index taking into account the difference between standard SM
831 // and SM of half size in phi direction
37890aaf 832 const Int_t kphiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
2bb3725c 833 static Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta;
e52475ed 834 if(!CheckAbsCellId(absId)) return kFALSE;
835
2bb3725c 836 GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
837 GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta);
e52475ed 838
1d46d1f6 839 xr = fCentersOfCellsXDir.At(ieta);
840 zr = fCentersOfCellsEtaDir.At(ieta);
e52475ed 841
1933eff2 842 if(nSupMod<10) {
1d46d1f6 843 yr = fCentersOfCellsPhiDir.At(iphi);
1933eff2 844 } else {
37890aaf 845 yr = fCentersOfCellsPhiDir.At(iphi + kphiIndexShift);
1933eff2 846 }
d25f2c54 847 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 848
e52475ed 849 return kTRUE;
850}
851
d297ef6e 852//________________________________________________________________________________________________
14e75ea7 853Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t loc[3]) const
854{
855 // Alice numbering scheme - Jun 03, 2006
856 loc[0] = loc[1] = loc[2]=0.0;
857 if(RelPosCellInSModule(absId, loc[0],loc[1],loc[2])) {
858 return kTRUE;
859 }
860 return kFALSE;
861}
862
d297ef6e 863//________________________________________________________________________________________________
14e75ea7 864Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, TVector3 &vloc) const
865{
866 static Double_t loc[3];
867 if(RelPosCellInSModule(absId,loc)) {
868 vloc.SetXYZ(loc[0], loc[1], loc[2]);
869 return kTRUE;
870 } else {
871 vloc.SetXYZ(0,0,0);
872 return kFALSE;
873 }
874 // Alice numbering scheme - Jun 03, 2006
875}
876
d297ef6e 877//________________________________________________________________________________________________
1ae500a2 878Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t distEff, Double_t &xr, Double_t &yr, Double_t &zr) const
879{
880 // Jul 30, 2007 - taking into account position of shower max
881 // Look to see what the relative
882 // position inside a given cell is
883 // for a recpoint.
884 // In:
885 // absId - cell is as in Geant, 0<= absId < fNCells;
886 // e - cluster energy
887 // OUT:
888 // xr,yr,zr - x,y,z coordinates of cell with absId inside SM
889
890 // Shift index taking into account the difference between standard SM
891 // and SM of half size in phi direction
37890aaf 892 const Int_t kphiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
1ae500a2 893 static Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta;
894 static Int_t iphim, ietam;
895 static AliEMCALShishKebabTrd1Module *mod = 0;
896 static TVector2 v;
897 if(!CheckAbsCellId(absId)) return kFALSE;
898
899 GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
900 GetModulePhiEtaIndexInSModule(nSupMod, nModule, iphim, ietam);
901 GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta);
902
903 mod = GetShishKebabModule(ietam);
904 mod->GetPositionAtCenterCellLine(nIeta, distEff, v);
905 xr = v.Y() - fParSM[0];
906 zr = v.X() - fParSM[2];
907
908 if(nSupMod<10) {
909 yr = fCentersOfCellsPhiDir.At(iphi);
910 } else {
37890aaf 911 yr = fCentersOfCellsPhiDir.At(iphi + kphiIndexShift);
1ae500a2 912 }
913 AliDebug(1,Form("absId %i nSupMod %i iphi %i ieta %i xr %f yr %f zr %f ",absId,nSupMod,iphi,ieta,xr,yr,zr));
914
915 return kTRUE;
916}
917
d297ef6e 918//________________________________________________________________________________________________
1ae500a2 919Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Int_t maxAbsId, Double_t distEff, Double_t &xr, Double_t &yr, Double_t &zr) const
920{
921 // Jul 31, 2007 - taking into account position of shower max and apply coor2.
922 // Look to see what the relative
923 // position inside a given cell is
924 // for a recpoint.
925 // In:
926 // absId - cell is as in Geant, 0<= absId < fNCells;
927 // maxAbsId - abs id of cell with highest energy
928 // e - cluster energy
929 // OUT:
930 // xr,yr,zr - x,y,z coordinates of cell with absId inside SM
931
932 // Shift index taking into account the difference between standard SM
933 // and SM of half size in phi direction
37890aaf 934 const Int_t kphiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
1ae500a2 935 static Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta;
936 static Int_t iphim, ietam;
937 static AliEMCALShishKebabTrd1Module *mod = 0;
938 static TVector2 v;
939
940 static Int_t nSupModM, nModuleM, nIphiM, nIetaM, iphiM, ietaM;
941 static Int_t iphimM, ietamM, maxAbsIdCopy=-1;
942 static AliEMCALShishKebabTrd1Module *modM = 0;
943 static Double_t distCorr;
944
945 if(!CheckAbsCellId(absId)) return kFALSE;
946
947 GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
948 GetModulePhiEtaIndexInSModule(nSupMod, nModule, iphim, ietam);
949 GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta);
950 mod = GetShishKebabModule(ietam);
951
952 if(absId != maxAbsId) {
953 distCorr = 0.;
954 if(maxAbsIdCopy != maxAbsId) {
955 GetCellIndex(maxAbsId, nSupModM, nModuleM, nIphiM, nIetaM);
956 GetModulePhiEtaIndexInSModule(nSupModM, nModuleM, iphimM, ietamM);
957 GetCellPhiEtaIndexInSModule(nSupModM,nModuleM,nIphiM,nIetaM, iphiM, ietaM);
958 modM = GetShishKebabModule(ietamM); // do I need this ?
959 maxAbsIdCopy = maxAbsId;
960 }
961 if(ietamM !=0) {
962 distCorr = GetEtaModuleSize()*(ietam-ietamM)/TMath::Tan(modM->GetTheta()); // Stay here
963 //printf(" distCorr %f | dist %f | ietam %i -> etamM %i\n", distCorr, dist, ietam, ietamM);
964 }
965 // distEff += distCorr;
966 }
967 // Bad resolution in this case, strong bias vs phi
968 // distEff = 0.0;
969 mod->GetPositionAtCenterCellLine(nIeta, distEff, v); // Stay here
970 xr = v.Y() - fParSM[0];
971 zr = v.X() - fParSM[2];
972
973 if(nSupMod<10) {
974 yr = fCentersOfCellsPhiDir.At(iphi);
975 } else {
37890aaf 976 yr = fCentersOfCellsPhiDir.At(iphi + kphiIndexShift);
1ae500a2 977 }
978 AliDebug(1,Form("absId %i nSupMod %i iphi %i ieta %i xr %f yr %f zr %f ",absId,nSupMod,iphi,ieta,xr,yr,zr));
979
980 return kTRUE;
981}
982
d297ef6e 983//________________________________________________________________________________________________
e52475ed 984void AliEMCALGeometry::CreateListOfTrd1Modules()
985{
1d46d1f6 986 // Generate the list of Trd1 modules
987 // which will make up the EMCAL
988 // geometry
fc575e27 989
e5a93224 990 AliDebug(2,Form(" AliEMCALGeometry::CreateListOfTrd1Modules() started "));
991
e52475ed 992 AliEMCALShishKebabTrd1Module *mod=0, *mTmp=0; // current module
993 if(fShishKebabTrd1Modules == 0) {
994 fShishKebabTrd1Modules = new TList;
1d46d1f6 995 fShishKebabTrd1Modules->SetName("ListOfTRD1");
e52475ed 996 for(int iz=0; iz< GetNZ(); iz++) {
997 if(iz==0) {
998 mod = new AliEMCALShishKebabTrd1Module(TMath::Pi()/2.,this);
999 } else {
1000 mTmp = new AliEMCALShishKebabTrd1Module(*mod);
1001 mod = mTmp;
1002 }
1003 fShishKebabTrd1Modules->Add(mod);
1004 }
1005 } else {
e5a93224 1006 AliDebug(2,Form(" Already exits : "));
e52475ed 1007 }
1d46d1f6 1008 mod = (AliEMCALShishKebabTrd1Module*)fShishKebabTrd1Modules->At(fShishKebabTrd1Modules->GetSize()-1);
1009 fEtaMaxOfTRD1 = mod->GetMaxEtaOfModule(0);
1010
1011 AliDebug(2,Form(" fShishKebabTrd1Modules has %i modules : max eta %5.4f \n",
1012 fShishKebabTrd1Modules->GetSize(),fEtaMaxOfTRD1));
e52475ed 1013 // Feb 20,2006;
dc7da436 1014 // Jun 01, 2006 - ALICE numbering scheme
e52475ed 1015 // define grid for cells in eta(z) and x directions in local coordinates system of SM
1d46d1f6 1016 // Works just for 2x2 case only -- ?? start here
1017 //
1018 //
1019 // Define grid for cells in phi(y) direction in local coordinates system of SM
1020 // as for 2X2 as for 3X3 - Nov 8,2006
1021 //
1022 AliDebug(2,Form(" Cells grid in phi directions : size %i\n", fCentersOfCellsPhiDir.GetSize()));
1023 Int_t ind=0; // this is phi index
85327f24 1024 Int_t ieta=0, nModule=0, iphiTemp;
1d46d1f6 1025 Double_t xr, zr, theta, phi, eta, r, x,y;
1026 TVector3 vglob;
85327f24 1027 Double_t ytCenterModule=0.0, ytCenterCell=0.0;
1d46d1f6 1028
1029 fCentersOfCellsPhiDir.Set(fNPhi*fNPHIdiv);
1030 fPhiCentersOfCells.Set(fNPhi*fNPHIdiv);
1031
37890aaf 1032 Double_t r0 = GetIPDistance() + GetLongModuleSize()/2.;
1d46d1f6 1033 for(Int_t it=0; it<fNPhi; it++) { // cycle on modules
1034 ytCenterModule = -fParSM[1] + fPhiModuleSize*(2*it+1)/2; // center of module
1035 for(Int_t ic=0; ic<fNPHIdiv; ic++) { // cycle on cells in module
1036 if(fNPHIdiv==2) {
1037 ytCenterCell = ytCenterModule + fPhiTileSize *(2*ic-1)/2.;
1038 } else if(fNPHIdiv==3){
1039 ytCenterCell = ytCenterModule + fPhiTileSize *(ic-1);
d25f2c54 1040 } else if(fNPHIdiv==1){
1041 ytCenterCell = ytCenterModule;
1d46d1f6 1042 }
1043 fCentersOfCellsPhiDir.AddAt(ytCenterCell,ind);
1044 // Define grid on phi direction
1045 // Grid is not the same for different eta bin;
1046 // Effect is small but is still here
37890aaf 1047 phi = TMath::ATan2(ytCenterCell, r0);
1d46d1f6 1048 fPhiCentersOfCells.AddAt(phi, ind);
1049
1050 AliDebug(2,Form(" ind %2.2i : y %8.3f ", ind, fCentersOfCellsPhiDir.At(ind)));
1051 ind++;
1052 }
1053 }
1054
1055 fCentersOfCellsEtaDir.Set(fNZ *fNETAdiv);
1056 fCentersOfCellsXDir.Set(fNZ *fNETAdiv);
1057 fEtaCentersOfCells.Set(fNZ *fNETAdiv * fNPhi*fNPHIdiv);
1058 AliDebug(2,Form(" Cells grid in eta directions : size %i\n", fCentersOfCellsEtaDir.GetSize()));
1059 for(Int_t it=0; it<fNZ; it++) {
e52475ed 1060 AliEMCALShishKebabTrd1Module *trd1 = GetShishKebabModule(it);
2bb3725c 1061 nModule = fNPhi*it;
1d46d1f6 1062 for(Int_t ic=0; ic<fNETAdiv; ic++) {
1063 if(fNPHIdiv==2) {
d25f2c54 1064 trd1->GetCenterOfCellInLocalCoordinateofSM(ic, xr, zr); // case of 2X2
2bb3725c 1065 GetCellPhiEtaIndexInSModule(0, nModule, 0, ic, iphiTemp, ieta);
1d46d1f6 1066 } if(fNPHIdiv==3) {
1067 trd1->GetCenterOfCellInLocalCoordinateofSM_3X3(ic, xr, zr); // case of 3X3
2bb3725c 1068 GetCellPhiEtaIndexInSModule(0, nModule, 0, ic, iphiTemp, ieta);
d25f2c54 1069 } if(fNPHIdiv==1) {
1070 trd1->GetCenterOfCellInLocalCoordinateofSM_1X1(xr, zr); // case of 1X1
2bb3725c 1071 GetCellPhiEtaIndexInSModule(0, nModule, 0, ic, iphiTemp, ieta);
1d46d1f6 1072 }
d25f2c54 1073 fCentersOfCellsXDir.AddAt(float(xr) - fParSM[0],ieta);
1074 fCentersOfCellsEtaDir.AddAt(float(zr) - fParSM[2],ieta);
1d46d1f6 1075 // Define grid on eta direction for each bin in phi
1076 for(int iphi=0; iphi<fCentersOfCellsPhiDir.GetSize(); iphi++) {
1077 x = xr + trd1->GetRadius();
1078 y = fCentersOfCellsPhiDir[iphi];
1079 r = TMath::Sqrt(x*x + y*y + zr*zr);
1080 theta = TMath::ACos(zr/r);
1081 eta = AliEMCALShishKebabTrd1Module::ThetaToEta(theta);
1082 // ind = ieta*fCentersOfCellsPhiDir.GetSize() + iphi;
1083 ind = iphi*fCentersOfCellsEtaDir.GetSize() + ieta;
1084 fEtaCentersOfCells.AddAt(eta, ind);
1085 }
1086 //printf(" ieta %i : xr + trd1->GetRadius() %f : zr %f : eta %f \n", ieta, xr + trd1->GetRadius(), zr, eta);
e52475ed 1087 }
1088 }
1d46d1f6 1089 for(Int_t i=0; i<fCentersOfCellsEtaDir.GetSize(); i++) {
e5a93224 1090 AliDebug(2,Form(" ind %2.2i : z %8.3f : x %8.3f", i+1,
1d46d1f6 1091 fCentersOfCellsEtaDir.At(i),fCentersOfCellsXDir.At(i)));
e52475ed 1092 }
e5a93224 1093
e52475ed 1094}
1095
d297ef6e 1096//________________________________________________________________________________________________
14e75ea7 1097void AliEMCALGeometry::GetGlobal(const Double_t *loc, Double_t *glob, int ind) const
e52475ed 1098{
14e75ea7 1099 // Figure out the global numbering
1100 // of a given supermodule from the
b8b0f8c2 1101 // local numbering and the transformation
1102 // matrix stored by the geometry manager (allows for misaligned
1103 // geometry)
14e75ea7 1104
e52475ed 1105 if(ind>=0 && ind < GetNumberOfSuperModules()) {
b8b0f8c2 1106 TString volpath = "ALIC_1/XEN1_1/SMOD_";
9aa6a5f6 1107 volpath += ind+1;
b8b0f8c2 1108
1109 if(GetKey110DEG() && ind>=10) {
1110 volpath = "ALIC_1/XEN1_1/SM10_";
9aa6a5f6 1111 volpath += ind-10+1;
b8b0f8c2 1112 }
1113
1114 if(!gGeoManager->cd(volpath.Data()))
9aa6a5f6 1115 AliFatal(Form("AliEMCALGeometry::GeoManager cannot find path %s!",volpath.Data()));
b8b0f8c2 1116
1117 TGeoHMatrix* m = gGeoManager->GetCurrentMatrix();
1118 if(m) {
1119 m->LocalToMaster(loc, glob);
1120 } else {
1121 AliFatal("Geo matrixes are not loaded \n") ;
1122 }
e52475ed 1123 }
1124}
1125
d297ef6e 1126//________________________________________________________________________________________________
25b033cf 1127void AliEMCALGeometry::GetGlobal(const TVector3 &vloc, TVector3 &vglob, int ind) const
1128{
1129 //Figure out the global numbering
1130 //of a given supermodule from the
1131 //local numbering given a 3-vector location
1132
1133 static Double_t tglob[3], tloc[3];
1134 vloc.GetXYZ(tloc);
1135 GetGlobal(tloc, tglob, ind);
1136 vglob.SetXYZ(tglob[0], tglob[1], tglob[2]);
1137}
1138
d297ef6e 1139//________________________________________________________________________________________________
14e75ea7 1140void AliEMCALGeometry::GetGlobal(Int_t absId , double glob[3]) const
1141{
1142 // Alice numbering scheme - Jun 03, 2006
1143 static Int_t nSupMod, nModule, nIphi, nIeta;
1144 static double loc[3];
1145
aaa3cb7c 1146 if (!gGeoManager || !gGeoManager->IsClosed()) {
1147 AliError("Can't get the global coordinates! gGeoManager doesn't exist or it is still open!");
1148 return;
1149 }
1150
14e75ea7 1151 glob[0]=glob[1]=glob[2]=0.0; // bad case
1152 if(RelPosCellInSModule(absId, loc)) {
1153 GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
b8b0f8c2 1154
1155 TString volpath = "ALIC_1/XEN1_1/SMOD_";
1156 volpath += (nSupMod+1);
1157
1158 if(GetKey110DEG() && nSupMod>=10) {
1159 volpath = "ALIC_1/XEN1_1/SM10_";
1160 volpath += (nSupMod-10+1);
1161 }
1162 if(!gGeoManager->cd(volpath.Data()))
1163 AliFatal(Form("GeoManager cannot find path %s!",volpath.Data()));
1164
1165 TGeoHMatrix* m = gGeoManager->GetCurrentMatrix();
1166 if(m) {
1167 m->LocalToMaster(loc, glob);
1168 } else {
1169 AliFatal("Geo matrixes are not loaded \n") ;
1170 }
14e75ea7 1171 }
e52475ed 1172}
1173
9aa6a5f6 1174//___________________________________________________________________
14e75ea7 1175void AliEMCALGeometry::GetGlobal(Int_t absId , TVector3 &vglob) const
1176{
1177 // Alice numbering scheme - Jun 03, 2006
1178 static Double_t glob[3];
1179
1180 GetGlobal(absId, glob);
1181 vglob.SetXYZ(glob[0], glob[1], glob[2]);
1182
1183}
1184
9aa6a5f6 1185//____________________________________________________________________________
1186void AliEMCALGeometry::GetGlobal(const AliRecPoint* /*rp*/, TVector3& /* vglob */) const
1187{
1188 AliFatal(Form("Please use GetGlobalEMCAL(recPoint,gpos) instead of GetGlobal!"));
1189}
1190
1191//_________________________________________________________________________________
1192void AliEMCALGeometry::GetGlobalEMCAL(const AliEMCALRecPoint *rp, TVector3 &vglob) const
e52475ed 1193{
664bfd66 1194 // Figure out the global numbering
1195 // of a given supermodule from the
1196 // local numbering for RecPoints
fc575e27 1197
e52475ed 1198 static TVector3 vloc;
14e75ea7 1199 static Int_t nSupMod, nModule, nIphi, nIeta;
e52475ed 1200
9aa6a5f6 1201 const AliEMCALRecPoint *rpTmp = rp;
1202 const AliEMCALRecPoint *rpEmc = rpTmp;
e52475ed 1203
14e75ea7 1204 GetCellIndex(rpEmc->GetAbsId(0), nSupMod, nModule, nIphi, nIeta);
e52475ed 1205 rpTmp->GetLocalPosition(vloc);
1206 GetGlobal(vloc, vglob, nSupMod);
1207}
1208
d297ef6e 1209//________________________________________________________________________________________________
1d46d1f6 1210void AliEMCALGeometry::EtaPhiFromIndex(Int_t absId,Double_t &eta,Double_t &phi) const
664bfd66 1211{
1d46d1f6 1212 // Nov 16, 2006- float to double
1213 // version for TRD1 only
664bfd66 1214 static TVector3 vglob;
1215 GetGlobal(absId, vglob);
1216 eta = vglob.Eta();
1217 phi = vglob.Phi();
1218}
1219
d297ef6e 1220//________________________________________________________________________________________________
1d46d1f6 1221void AliEMCALGeometry::EtaPhiFromIndex(Int_t absId,Float_t &eta,Float_t &phi) const
1222{
1223 // Nov 16,2006 - should be discard in future
1224 static TVector3 vglob;
1225 GetGlobal(absId, vglob);
1226 eta = float(vglob.Eta());
1227 phi = float(vglob.Phi());
1228}
1229
d297ef6e 1230//________________________________________________________________________________________________
1d46d1f6 1231Bool_t AliEMCALGeometry::GetPhiBoundariesOfSM(Int_t nSupMod, Double_t &phiMin, Double_t &phiMax) const
1232{
1233 // 0<= nSupMod <=11; phi in rad
1234 static int i;
1235 if(nSupMod<0 || nSupMod >11) return kFALSE;
1236 i = nSupMod/2;
1237 phiMin = fPhiBoundariesOfSM[2*i];
1238 phiMax = fPhiBoundariesOfSM[2*i+1];
1239 return kTRUE;
1240}
1241
d297ef6e 1242//________________________________________________________________________________________________
1d46d1f6 1243Bool_t AliEMCALGeometry::GetPhiBoundariesOfSMGap(Int_t nPhiSec, Double_t &phiMin, Double_t &phiMax) const
1244{
1245 // 0<= nPhiSec <=4; phi in rad
1246 // 0; gap boundaries between 0th&2th | 1th&3th SM
1247 // 1; gap boundaries between 2th&4th | 3th&5th SM
1248 // 2; gap boundaries between 4th&6th | 5th&7th SM
1249 // 3; gap boundaries between 6th&8th | 7th&9th SM
1250 // 4; gap boundaries between 8th&10th | 9th&11th SM
1251 if(nPhiSec<0 || nPhiSec >4) return kFALSE;
1252 phiMin = fPhiBoundariesOfSM[2*nPhiSec+1];
1253 phiMax = fPhiBoundariesOfSM[2*nPhiSec+2];
1254 return kTRUE;
1255}
1256
d297ef6e 1257//________________________________________________________________________________________________
1d46d1f6 1258Bool_t AliEMCALGeometry::SuperModuleNumberFromEtaPhi(Double_t eta, Double_t phi, Int_t &nSupMod) const
1259{
1260 // Return false if phi belongs a phi cracks between SM
1261
1262 static Int_t i;
1263
1264 if(TMath::Abs(eta) > fEtaMaxOfTRD1) return kFALSE;
1265
1266 phi = TVector2::Phi_0_2pi(phi); // move phi to (0,2pi) boundaries
1267 for(i=0; i<6; i++) {
1268 if(phi>=fPhiBoundariesOfSM[2*i] && phi<=fPhiBoundariesOfSM[2*i+1]) {
1269 nSupMod = 2*i;
1270 if(eta < 0.0) nSupMod++;
d25f2c54 1271 AliDebug(1,Form("eta %f phi %f(%5.2f) : nSupMod %i : #bound %i", eta,phi,phi*TMath::RadToDeg(), nSupMod,i));
1d46d1f6 1272 return kTRUE;
1273 }
1274 }
1d46d1f6 1275 return kFALSE;
1276}
1277
d297ef6e 1278//________________________________________________________________________________________________
1d46d1f6 1279Bool_t AliEMCALGeometry::GetAbsCellIdFromEtaPhi(Double_t eta, Double_t phi, Int_t &absId) const
1280{
1281 // Nov 17,2006
1282 // stay here - phi problem as usual
1283 static Int_t nSupMod, i, ieta, iphi, etaShift, nphi;
1284 static Double_t absEta=0.0, d=0.0, dmin=0.0, phiLoc;
1285 absId = nSupMod = - 1;
1286 if(SuperModuleNumberFromEtaPhi(eta, phi, nSupMod)) {
1287 // phi index first
1288 phi = TVector2::Phi_0_2pi(phi);
1289 phiLoc = phi - fPhiCentersOfSM[nSupMod/2];
1290 nphi = fPhiCentersOfCells.GetSize();
1291 if(nSupMod>=10) {
1292 phiLoc = phi - 190.*TMath::DegToRad();
1293 nphi /= 2;
1294 }
1295
1296 dmin = TMath::Abs(fPhiCentersOfCells[0]-phiLoc);
1297 iphi = 0;
1298 for(i=1; i<nphi; i++) {
1299 d = TMath::Abs(fPhiCentersOfCells[i] - phiLoc);
1300 if(d < dmin) {
1301 dmin = d;
1302 iphi = i;
1303 }
1304 // printf(" i %i : d %f : dmin %f : fPhiCentersOfCells[i] %f \n", i, d, dmin, fPhiCentersOfCells[i]);
1305 }
1306 // odd SM are turned with respect of even SM - reverse indexes
1307 AliDebug(2,Form(" iphi %i : dmin %f (phi %f, phiLoc %f ) ", iphi, dmin, phi, phiLoc));
1308 // eta index
1309 absEta = TMath::Abs(eta);
1310 etaShift = iphi*fCentersOfCellsEtaDir.GetSize();
1311 dmin = TMath::Abs(fEtaCentersOfCells[etaShift]-absEta);
1312 ieta = 0;
1313 for(i=1; i<fCentersOfCellsEtaDir.GetSize(); i++) {
1314 d = TMath::Abs(fEtaCentersOfCells[i+etaShift] - absEta);
1315 if(d < dmin) {
1316 dmin = d;
1317 ieta = i;
1318 }
1319 }
1320 AliDebug(2,Form(" ieta %i : dmin %f (eta=%f) : nSupMod %i ", ieta, dmin, eta, nSupMod));
1321
1322 if(eta<0) iphi = (nphi-1) - iphi;
1323 absId = GetAbsCellIdFromCellIndexes(nSupMod, iphi, ieta);
1324
1325 return kTRUE;
1326 }
1327 return kFALSE;
1328}
1329
d297ef6e 1330//________________________________________________________________________________________________
1ae500a2 1331AliEMCALShishKebabTrd1Module* AliEMCALGeometry::GetShishKebabModule(Int_t neta) const
fc575e27 1332{
1333 //This method was too long to be
1334 //included in the header file - the
1335 //rule checker complained about it's
1336 //length, so we move it here. It returns the
1337 //shishkebabmodule at a given eta index point.
1338
1339 static AliEMCALShishKebabTrd1Module* trd1=0;
1340 if(fShishKebabTrd1Modules && neta>=0 && neta<fShishKebabTrd1Modules->GetSize()) {
1341 trd1 = (AliEMCALShishKebabTrd1Module*)fShishKebabTrd1Modules->At(neta);
1342 } else trd1 = 0;
1343 return trd1;
1344}
1d46d1f6 1345
85c25c2e 1346//________________________________________________________________________________________________
1347Int_t AliEMCALGeometry::GetAbsTRUNumberFromNumberInSm(const Int_t row, const Int_t col, const Int_t sm)
1348{ // Nov 6, 2007
1349 Int_t itru = row + col*GetNModulesInTRUPhi() + sm*GetNTRU();
1350 // printf(" GetAbsTRUNumberFromNumberInSm : row %2i col %2i sm %2i -> itru %2i\n", row, col, sm, itru);
1351 return itru;
1352}
1353
d297ef6e 1354//________________________________________________________________________________________________
6f377f0c 1355void AliEMCALGeometry::Browse(TBrowser* b)
1d46d1f6 1356{
37890aaf 1357 //Browse the modules
1d46d1f6 1358 if(fShishKebabTrd1Modules) b->Add(fShishKebabTrd1Modules);
1359}
1360
d297ef6e 1361//________________________________________________________________________________________________
1d46d1f6 1362Bool_t AliEMCALGeometry::IsFolder() const
1363{
37890aaf 1364 //Check if fShishKebabTrd1Modules is in folder
1d46d1f6 1365 if(fShishKebabTrd1Modules) return kTRUE;
1366 else return kFALSE;
1367}
9aa6a5f6 1368
d297ef6e 1369//________________________________________________________________________________________________
9aa6a5f6 1370Double_t AliEMCALGeometry::GetPhiCenterOfSM(Int_t nsupmod) const
1371{
37890aaf 1372 //returns center of supermodule in phi
4238c618 1373 int i = nsupmod/2;
9aa6a5f6 1374 return fPhiCentersOfSM[i];
1375
1376}
225cd96d 1377//____________________________________________________________________________
1378Bool_t AliEMCALGeometry::Impact(const TParticle * particle) const
1379{
1380 // Tells if a particle enters EMCAL
1381 Bool_t in=kFALSE;
1382 Int_t AbsID=0;
1383 TVector3 vtx(particle->Vx(),particle->Vy(),particle->Vz());
1384 TVector3 vimpact(0,0,0);
1385 ImpactOnEmcal(vtx,particle->Theta(),particle->Phi(),AbsID,vimpact);
1386 if(AbsID!=0)
1387 in=kTRUE;
1388 return in;
1389}
1390//____________________________________________________________________________
1391void AliEMCALGeometry::ImpactOnEmcal(TVector3 vtx, Double_t theta, Double_t phi,
1392 Int_t & absId, TVector3 & vimpact) const
1393{
1394 // calculates the impact coordinates on EMCAL (centre of a tower/not on EMCAL surface)
1395 // of a neutral particle
1396 // emitted in the vertex vtx[3] with direction theta and phi in the ALICE global coordinate system
1397
1398 TVector3 p(TMath::Sin(theta)*TMath::Cos(phi),TMath::Sin(theta)*TMath::Sin(phi),TMath::Cos(theta)) ;
1399
1400 vimpact.SetXYZ(0,0,0);
1401 absId=-1;
1402 if(phi==0 || theta==0) return;
1403
1404 TVector3 direction;
1405 Double_t factor = (GetIPDistance()-vtx[1])/p[1];
1406 direction = vtx + factor*p;
1407
1408 if (!gGeoManager){
1409 AliFatal("Geo manager not initialized\n");
1410 }
1411 //from particle direction -> tower hitted
1412 GetAbsCellIdFromEtaPhi(direction.Eta(),direction.Phi(),absId);
1413
1414 //tower absID hitted -> tower/module plane (evaluated at the center of the tower)
1415 Int_t nSupMod, nModule, nIphi, nIeta;
1416 Double_t loc[3],loc2[3],loc3[3];
1417 Double_t glob[3]={},glob2[3]={},glob3[3]={};
1418
1419 if(!RelPosCellInSModule(absId,loc)) return;
1420
1421 //loc is cell center of tower
1422 GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
1423
1424 //look at 2 neighbours-s cell using nIphi={0,1} and nIeta={0,1}
1425 Int_t nIphi2,nIeta2,absId2,absId3;
1426 if(nIeta==0) nIeta2=1;
1427 else nIeta2=0;
1428 absId2=GetAbsCellId(nSupMod,nModule,nIphi,nIeta2);
1429 if(nIphi==0) nIphi2=1;
1430 else nIphi2=0;
1431 absId3=GetAbsCellId(nSupMod,nModule,nIphi2,nIeta);
1432
1433 //2nd point on emcal cell plane
1434 if(!RelPosCellInSModule(absId2,loc2)) return;
1435
1436 //3rd point on emcal cell plane
1437 if(!RelPosCellInSModule(absId3,loc3)) return;
1438
1439 TString volpath = "ALIC_1/XEN1_1/SMOD_";
1440 volpath += (nSupMod+1);
1441
1442 if(GetKey110DEG() && nSupMod>=10) {
1443 volpath = "ALIC_1/XEN1_1/SM10_";
1444 volpath += (nSupMod-10+1);
1445 }
1446 if(!gGeoManager->cd(volpath.Data())){
1447 AliFatal(Form("GeoManager cannot find path %s!",volpath.Data()))
1448 return;
1449 }
1450 TGeoHMatrix* m = gGeoManager->GetCurrentMatrix();
1451 if(m) {
1452 m->LocalToMaster(loc, glob);
1453 m->LocalToMaster(loc2, glob2);
1454 m->LocalToMaster(loc3, glob3);
1455 } else {
1456 AliFatal("Geo matrixes are not loaded \n") ;
1457 }
1458
1459 //Equation of Plane from glob,glob2,glob3 (Ax+By+Cz+D=0)
1460 Double_t A = glob[1]*(glob2[2]-glob3[2]) + glob2[1]*(glob3[2]-glob[2]) + glob3[1]*(glob[2]-glob2[2]);
1461 Double_t B = glob[2]*(glob2[0]-glob3[0]) + glob2[2]*(glob3[0]-glob[0]) + glob3[2]*(glob[0]-glob2[0]);
1462 Double_t C = glob[0]*(glob2[1]-glob3[1]) + glob2[0]*(glob3[1]-glob[1]) + glob3[0]*(glob[1]-glob2[1]);
1463 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]);
1464 D=-D;
1465
1466 //shift equation of plane from tower/module center to surface along vector (A,B,C) normal to tower/module plane
1467 Double_t dist = GetLongModuleSize()/2.;
1468 Double_t norm = TMath::Sqrt(A*A+B*B+C*C);
1469 Double_t glob4[3]={};
1470 TVector3 dir(A,B,C);
1471 TVector3 point(glob[0],glob[1],glob[2]);
1472 if(point.Dot(dir)<0) dist*=-1;
1473 glob4[0]=glob[0]-dist*A/norm;
1474 glob4[1]=glob[1]-dist*B/norm;
1475 glob4[2]=glob[2]-dist*C/norm;
1476 D = glob4[0]*A + glob4[1]*B + glob4[2]*C ;
1477 D = -D;
1478
1479 //Line determination (2 points for equation of line : vtx and direction)
1480 //impact between line (particle) and plane (module/tower plane)
1481 Double_t den = A*(vtx(0)-direction(0)) + B*(vtx(1)-direction(1)) + C*(vtx(2)-direction(2));
1482 if(den==0){
1483 printf("ImpactOnEmcal() No solution :\n");
1484 return;
1485 }
1486
1487 Double_t length = A*vtx(0)+B*vtx(1)+C*vtx(2)+D;
1488 length /=den;
1489
1490 vimpact.SetXYZ(vtx(0)+length*(direction(0)-vtx(0)),vtx(1)+length*(direction(1)-vtx(1)),vtx(2)+length*(direction(2)-vtx(2)));
1491
1492 //shift vimpact from tower/module surface to center along vector (A,B,C) normal to tower/module plane
1493 vimpact.SetXYZ(vimpact(0)+dist*A/norm,vimpact(1)+dist*B/norm,vimpact(2)+dist*C/norm);
1494
1495 return;
1496}