Coverity corrections
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALEMCGeometry.cxx
CommitLineData
0c5b726e 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id: AliEMCALEMCGeometry.cxx 29514 2008-10-26 10:24:38Z hristov $*/
17
18//_________________________________________________________________________
19// Geometry class for EMCAL : singleton
20// EMCAL consists of layers of scintillator and lead
21// with scintillator fiber arranged as "shish-kebab" skewers
22// Places the the Barrel Geometry of The EMCAL at Midrapidity
23// between 80 and 180(or 190) degrees of Phi and
24// -0.7 to 0.7 in eta
25//
26// EMCAL geometry tree:
27// EMCAL -> superModule -> module -> tower(cell)
28// Indexes
29// absId -> nSupMod -> nModule -> (nIphi,nIeta)
30//
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//
45//*-- Author: Sahal Yacoob (LBL / UCT)
46// and : Yves Schutz (SUBATECH)
47// and : Jennifer Klay (LBL)
48// and : Aleksei Pavlinov (WSU)
49// and : Magali Estienne (SUBATECH)
50
51// --- Root header files ---
52#include <TObjArray.h>
53#include <TObjString.h>
54#include <TRegexp.h>
55
56// -- ALICE Headers.
57#include "AliLog.h"
58
59// --- EMCAL headers
60#include "AliEMCALEMCGeometry.h"
61#include <cassert>
62
63ClassImp(AliEMCALEMCGeometry)
64
65// these initialisations are needed for a singleton
66Bool_t AliEMCALEMCGeometry::fgInit = kFALSE;
67const Char_t* AliEMCALEMCGeometry::fgkDefaultGeometryName = "EMCAL_COMPLETE";
68
69
70AliEMCALEMCGeometry::AliEMCALEMCGeometry()
71 : TNamed(),
72 fGeoName(0),fArrayOpts(0),fNAdditionalOpts(0),fECPbRadThickness(0.),fECScintThick(0.),
73 fNECLayers(0),fArm1PhiMin(0.),fArm1PhiMax(0.),fArm1EtaMin(0.),fArm1EtaMax(0.),fIPDistance(0.),
74 fShellThickness(0.),fZLength(0.),fNZ(0),fNPhi(0),fSampling(0.),fNumberOfSuperModules(0),
75 fFrontSteelStrip(0.),fLateralSteelStrip(0.),fPassiveScintThick(0.),fPhiModuleSize(0.),
76 fEtaModuleSize(0.),fPhiTileSize(0.),fEtaTileSize(0.),fLongModuleSize(0.),fNPhiSuperModule(0),
77 fNPHIdiv(0),fNETAdiv(0), fNCells(0),fNCellsInSupMod(0),fNCellsInModule(0),
78 // Trigger staff
79 fNTRUEta(0), fNTRUPhi(0), fNModulesInTRUEta(0), fNModulesInTRUPhi(0), fNEtaSubOfTRU(0),
80 //
81 fTrd1Angle(0.),f2Trd1Dx2(0.),
82 fPhiGapForSM(0.),fKey110DEG(0),fPhiBoundariesOfSM(0), fPhiCentersOfSM(0),fEtaMaxOfTRD1(0),
83 fCentersOfCellsEtaDir(0), fCentersOfCellsXDir(0),fCentersOfCellsPhiDir(0),
84 fEtaCentersOfCells(0),fPhiCentersOfCells(0),fShishKebabTrd1Modules(0),
85 fILOSS(-1), fIHADR(-1),
86 //obsolete member data
87 fAlFrontThick(0.), fGap2Active(0.), fSteelFrontThick(0.), fTrd2AngleY(0.),
88 f2Trd2Dy2(0.), fEmptySpace(0.), fTubsR(0.), fTubsTurnAngle(0.)
89{
90 // default ctor only for internal usage (singleton)
91 // must be kept public for root persistency purposes,
92 // but should never be called by the outside world
7e1d9a9b 93 fParSM[0]=0; fParSM[1]=0; fParSM[2]=0;
a5f0a2ba 94 fEnvelop[0] = 0; fEnvelop[1] = 0; fEnvelop[2] = 0;
a51e676d 95 for(Int_t i = 0; i < 6; i++) fkAdditionalOpts[i] = "";
a5f0a2ba 96
0c5b726e 97 AliDebug(2, "AliEMCALEMCGeometry : default ctor ");
98}
99//______________________________________________________________________
100AliEMCALEMCGeometry::AliEMCALEMCGeometry(const Text_t* name, const Text_t* title) :
101 TNamed(name,title),
102 fGeoName(0),fArrayOpts(0),fNAdditionalOpts(0),fECPbRadThickness(0.),fECScintThick(0.),
103 fNECLayers(0),fArm1PhiMin(0.),fArm1PhiMax(0.),fArm1EtaMin(0.),fArm1EtaMax(0.),fIPDistance(0.),
104 fShellThickness(0.),fZLength(0.),fNZ(0),fNPhi(0),fSampling(0.),fNumberOfSuperModules(0),
105 fFrontSteelStrip(0.),fLateralSteelStrip(0.),fPassiveScintThick(0.),fPhiModuleSize(0.),
106 fEtaModuleSize(0.),fPhiTileSize(0.),fEtaTileSize(0.),fLongModuleSize(0.),fNPhiSuperModule(0),
107 fNPHIdiv(0),fNETAdiv(0), fNCells(0),fNCellsInSupMod(0),fNCellsInModule(0),
108 // Trigger staff
109 fNTRUEta(0), fNTRUPhi(0), fNModulesInTRUEta(0), fNModulesInTRUPhi(0), fNEtaSubOfTRU(0),
110 //
111 fTrd1Angle(0.),f2Trd1Dx2(0.),
112 fPhiGapForSM(0.),fKey110DEG(0),fPhiBoundariesOfSM(0), fPhiCentersOfSM(0), fEtaMaxOfTRD1(0),
113 fCentersOfCellsEtaDir(0),fCentersOfCellsXDir(0),fCentersOfCellsPhiDir(0),
114 fEtaCentersOfCells(0),fPhiCentersOfCells(0),fShishKebabTrd1Modules(0),
115 fILOSS(-1), fIHADR(-1),
116 //obsolete member data
117 fAlFrontThick(0.), fGap2Active(0.), fSteelFrontThick(0.), fTrd2AngleY(0.),
118 f2Trd2Dy2(0.), fEmptySpace(0.), fTubsR(0.), fTubsTurnAngle(0.)
119{
120 // ctor only for internal usage (singleton)
121 AliDebug(2, Form("AliEMCALEMCGeometry(%s,%s) ", name,title));
122
123 Init();
124
125 // CreateListOfTrd1Modules();
126
127 if (AliDebugLevel()>=2) {
128 PrintGeometry();
129 }
130
131}
132//______________________________________________________________________
133AliEMCALEMCGeometry::AliEMCALEMCGeometry(const AliEMCALEMCGeometry& geom)
134 : TNamed(geom),
135 fGeoName(geom.fGeoName),
136 fArrayOpts(geom.fArrayOpts),
137 fNAdditionalOpts(geom.fNAdditionalOpts),
138 fECPbRadThickness(geom.fECPbRadThickness),
139 fECScintThick(geom.fECScintThick),
140 fNECLayers(geom.fNECLayers),
141 fArm1PhiMin(geom.fArm1PhiMin),
142 fArm1PhiMax(geom.fArm1PhiMax),
143 fArm1EtaMin(geom.fArm1EtaMin),
144 fArm1EtaMax(geom.fArm1EtaMax),
145 fIPDistance(geom.fIPDistance),
146 fShellThickness(geom.fShellThickness),
147 fZLength(geom.fZLength),
148 fNZ(geom.fNZ),
149 fNPhi(geom.fNPhi),
150 fSampling(geom.fSampling),
151 fNumberOfSuperModules(geom.fNumberOfSuperModules),
152 fFrontSteelStrip(geom.fFrontSteelStrip),
153 fLateralSteelStrip(geom.fLateralSteelStrip),
154 fPassiveScintThick(geom.fPassiveScintThick),
155 fPhiModuleSize(geom.fPhiModuleSize),
156 fEtaModuleSize(geom.fEtaModuleSize),
157 fPhiTileSize(geom.fPhiTileSize),
158 fEtaTileSize(geom.fEtaTileSize),
159 fLongModuleSize(geom.fLongModuleSize),
160 fNPhiSuperModule(geom.fNPhiSuperModule),
161 fNPHIdiv(geom.fNPHIdiv),
162 fNETAdiv(geom.fNETAdiv),
163 fNCells(geom.fNCells),
164 fNCellsInSupMod(geom.fNCellsInSupMod),
165 fNCellsInModule(geom.fNCellsInModule),
166 // Trigger staff
167 fNTRUEta(geom.fNTRUEta),
168 fNTRUPhi(geom.fNTRUPhi),
169 fNModulesInTRUEta(geom.fNModulesInTRUEta),
170 fNModulesInTRUPhi(geom.fNModulesInTRUPhi),
171 fNEtaSubOfTRU(geom.fNEtaSubOfTRU),
172 //
173 fTrd1Angle(geom.fTrd1Angle),
174 f2Trd1Dx2(geom.f2Trd1Dx2),
175 fPhiGapForSM(geom.fPhiGapForSM),
176 fKey110DEG(geom.fKey110DEG),
177 fPhiBoundariesOfSM(geom.fPhiBoundariesOfSM),
178 fPhiCentersOfSM(geom.fPhiCentersOfSM),
179 fEtaMaxOfTRD1(geom.fEtaMaxOfTRD1),
180 fCentersOfCellsEtaDir(geom.fCentersOfCellsEtaDir),
181 fCentersOfCellsXDir(geom.fCentersOfCellsXDir),
182 fCentersOfCellsPhiDir(geom.fCentersOfCellsPhiDir),
183 fEtaCentersOfCells(geom.fEtaCentersOfCells),
184 fPhiCentersOfCells(geom.fPhiCentersOfCells),
185 fShishKebabTrd1Modules(geom.fShishKebabTrd1Modules),
186 fILOSS(geom.fILOSS), fIHADR(geom.fIHADR),
187 //obsolete member data
188 fAlFrontThick(geom.fAlFrontThick),
189 fGap2Active(geom.fGap2Active),
190 fSteelFrontThick(geom.fSteelFrontThick),
191 fTrd2AngleY(geom.fTrd2AngleY),
192 f2Trd2Dy2(geom.f2Trd2Dy2),
193 fEmptySpace(geom.fEmptySpace),
194 fTubsR(geom.fTubsR),
195 fTubsTurnAngle(geom.fTubsTurnAngle)
196{
197 //copy ctor
7e1d9a9b 198 fParSM[0]=geom.fParSM[0];
199 fParSM[1]=geom.fParSM[1];
200 fParSM[2]=geom.fParSM[2];
a5f0a2ba 201 fEnvelop[0] = geom.fEnvelop[0];
202 fEnvelop[1] = geom.fEnvelop[1];
203 fEnvelop[2] = geom.fEnvelop[2];
a51e676d 204 for(Int_t i = 0; i < 6; i++) fkAdditionalOpts[i] = geom.fkAdditionalOpts[i];
a5f0a2ba 205
0c5b726e 206}
207
208//______________________________________________________________________
209AliEMCALEMCGeometry::~AliEMCALEMCGeometry(void){
210 // dtor
211}
212
213//______________________________________________________________________
214void AliEMCALEMCGeometry::Init(void){
215 //
216 // Initializes the EMCAL parameters based on the name
217 // Only Shashlyk geometry is available, but various combinations of
218 // layers and number of supermodules can be selected with additional
219 // options or geometry name
220 //
221
222 fkAdditionalOpts[0] = "nl="; // number of sampling layers (fNECLayers)
223 fkAdditionalOpts[1] = "pbTh="; // cm, Thickness of the Pb (fECPbRadThick)
224 fkAdditionalOpts[2] = "scTh="; // cm, Thickness of the Sc (fECScintThick)
225 fkAdditionalOpts[3] = "latSS="; // cm, Thickness of lateral steel strip (fLateralSteelStrip)
226 fkAdditionalOpts[4] = "allILOSS="; // = 0,1,2,3,4 (4 - energy loss without fluctuation)
227 fkAdditionalOpts[5] = "allIHADR="; // = 0,1,2 (0 - no hadronic interaction)
228
229 fNAdditionalOpts = sizeof(fkAdditionalOpts) / sizeof(char*);
230
231 // geometry
232 fgInit = kFALSE; // Assume failed until proven otherwise.
233 fGeoName = GetName();
234 fGeoName.ToUpper();
235
236 //Convert old geometry names to new ones
237 if(fGeoName.Contains("SHISH_77_TRD1_2X2_FINAL_110DEG")) {
238 if(fGeoName.Contains("PBTH=0.144") && fGeoName.Contains("SCTH=0.176")) {
239 fGeoName = "EMCAL_COMPLETE";
240 } else {
241 fGeoName = "EMCAL_PDC06";
242 }
243 }
244 if(fGeoName.Contains("WSUC")) fGeoName = "EMCAL_WSUC";
245
246 //check that we have a valid geometry name
9f467289 247 if(!(fGeoName.Contains("EMCAL_PDC06") || fGeoName.Contains("EMCAL_COMPLETE") || fGeoName.Contains("EMCAL_WSUC") || fGeoName.Contains("EMCAL_FIRSTYEAR"))) {
0c5b726e 248 Fatal("Init", "%s is an undefined geometry!", fGeoName.Data()) ;
249 }
250
251 // Option to know whether we have the "half" supermodule(s) or not
252 fKey110DEG = 0;
253 if(fGeoName.Contains("COMPLETE") || fGeoName.Contains("PDC06")) fKey110DEG = 1; // for GetAbsCellId
254 fShishKebabTrd1Modules = 0;
255
256 // JLK 13-Apr-2008
257 //default parameters are those of EMCAL_COMPLETE geometry
258 //all others render variations from these at the end of
259 //geometry-name specific options
260
261 fNumberOfSuperModules = 12; // 12 = 6 * 2 (6 in phi, 2 in Z)
262 fNPhi = 12; // module granularity in phi within smod (azimuth)
263 fNZ = 24; // module granularity along Z within smod (eta)
264 fNPHIdiv = fNETAdiv = 2; // tower granularity within module
265 fArm1PhiMin = 80.0; // degrees, Starting EMCAL Phi position
266 fArm1PhiMax = 200.0; // degrees, Ending EMCAL Phi position
267 fArm1EtaMin = -0.7; // pseudorapidity, Starting EMCAL Eta position
268 fArm1EtaMax = +0.7; // pseudorapidity, Ending EMCAL Eta position
269 fIPDistance = 428.0; // cm, radial distance to front face from nominal vertex point
270 fPhiGapForSM = 2.; // cm, only for final TRD1 geometry
271 fFrontSteelStrip = 0.025; // 0.025cm = 0.25mm (13-may-05 from V.Petrov)
272 fPassiveScintThick = 0.8; // 0.8cm = 8mm (13-may-05 from V.Petrov)
273 fLateralSteelStrip = 0.01; // 0.01cm = 0.1mm (13-may-05 from V.Petrov) - was 0.025
274 fTrd1Angle = 1.5; // in degrees
275
276 fSampling = 1.; // should be calculated with call to DefineSamplingFraction()
277 fNECLayers = 77; // (13-may-05 from V.Petrov) - can be changed with additional options
278 fECScintThick = 0.176; // scintillator layer thickness
279 fECPbRadThickness = 0.144; // lead layer thickness
280
281 fPhiModuleSize = 12.26 - fPhiGapForSM / Float_t(fNPhi); // first assumption
282 fEtaModuleSize = fPhiModuleSize;
283
284 fZLength = 700.; // Z coverage (cm)
285
286
287 //needs to be called for each geometry and before setting geometry
288 //parameters which can depend on the outcome
289 CheckAdditionalOptions();
290
291 //modifications to the above for PDC06 geometry
292 if(fGeoName.Contains("PDC06")){ // 18-may-05 - about common structure
293 fECScintThick = fECPbRadThickness = 0.16;// (13-may-05 from V.Petrov)
294 CheckAdditionalOptions();
295 }
296
297 //modifications to the above for WSUC geometry
298 if(fGeoName.Contains("WSUC")){ // 18-may-05 - about common structure
299 fPhiModuleSize = 12.5; // 20-may-05 - rectangular shape
300 fEtaModuleSize = 11.9;
301 fECScintThick = fECPbRadThickness = 0.16;// (13-may-05 from V.Petrov)
302 fNumberOfSuperModules = 1; // 27-may-05
303 fShellThickness = 30.; // should be change
304 fNPhi = fNZ = 4;
305 CheckAdditionalOptions();
306 }
307
9f467289 308 //In 2009-2010 data taking runs only 4 SM, in the upper position.
309 if(fGeoName.Contains("FIRSTYEAR")){
310 fNumberOfSuperModules = 4;
311 fArm1PhiMax = 120.0;
312 CheckAdditionalOptions();
313 }
314
0c5b726e 315 // constant for transition absid <--> indexes
9f467289 316 fNCellsInModule = fNPHIdiv*fNETAdiv;
0c5b726e 317 fNCellsInSupMod = fNCellsInModule*fNPhi*fNZ;
318 fNCells = fNCellsInSupMod*fNumberOfSuperModules;
319 if(GetKey110DEG()) fNCells -= fNCellsInSupMod;
320
321 fNPhiSuperModule = fNumberOfSuperModules/2;
322 if(fNPhiSuperModule < 1) fNPhiSuperModule = 1;
323
324 fPhiTileSize = fPhiModuleSize/double(fNPHIdiv) - fLateralSteelStrip; // 13-may-05
325 fEtaTileSize = fEtaModuleSize/double(fNETAdiv) - fLateralSteelStrip; // 13-may-05
326
327 fLongModuleSize = fNECLayers*(fECScintThick + fECPbRadThickness);
328 f2Trd1Dx2 = fEtaModuleSize + 2.*fLongModuleSize*TMath::Tan(fTrd1Angle*TMath::DegToRad()/2.);
329 if(!fGeoName.Contains("WSUC")) fShellThickness = TMath::Sqrt(fLongModuleSize*fLongModuleSize + f2Trd1Dx2*f2Trd1Dx2);
330
331 //These parameters are used to create the mother volume to hold the supermodules
332 //2cm padding added to allow for misalignments - JLK 30-May-2008
333 fEnvelop[0] = fIPDistance - 1.; // mother volume inner radius
334 fEnvelop[1] = fIPDistance + fShellThickness + 1.; // mother volume outer r.
335 fEnvelop[2] = fZLength + 2.; //mother volume length
336
337 // Local coordinates
338 fParSM[0] = GetShellThickness()/2.;
339 fParSM[1] = GetPhiModuleSize() * GetNPhi()/2.;
340 fParSM[2] = fZLength/4.; //divide by 4 to get half-length of SM
341
342 // SM phi boundaries - (0,1),(2,3) .. (10,11) - has the same boundaries; Nov 7, 2006
343 fPhiBoundariesOfSM.Set(fNumberOfSuperModules);
344 fPhiCentersOfSM.Set(fNumberOfSuperModules/2);
345 fPhiBoundariesOfSM[0] = TMath::PiOver2() - TMath::ATan2(fParSM[1] , fIPDistance); // 1th and 2th modules)
346 fPhiCentersOfSM[0] = TMath::PiOver2();
347 if(fNumberOfSuperModules > 1)
348 fPhiBoundariesOfSM[1] = TMath::PiOver2() + TMath::ATan2(fParSM[1] , fIPDistance);
349 if(fNumberOfSuperModules > 2) {
9f467289 350 Int_t maxPhiBlock =fNumberOfSuperModules/2-1;
351 if(fNumberOfSuperModules > 10) maxPhiBlock = 4;
352 for(int i=1; i<=maxPhiBlock; i++) { // from 2th ro 9th
0c5b726e 353 fPhiBoundariesOfSM[2*i] = fPhiBoundariesOfSM[0] + 20.*TMath::DegToRad()*i;
354 fPhiBoundariesOfSM[2*i+1] = fPhiBoundariesOfSM[1] + 20.*TMath::DegToRad()*i;
9f467289 355 fPhiCentersOfSM[i] = fPhiCentersOfSM[0] + 20.*TMath::DegToRad()*i;
0c5b726e 356 }
357 }
358 if(fNumberOfSuperModules > 10) {
359 fPhiBoundariesOfSM[11] = 190.*TMath::DegToRad();
360 fPhiBoundariesOfSM[10] = fPhiBoundariesOfSM[11] - TMath::ATan2((fParSM[1]) , fIPDistance);
9f467289 361 fPhiCentersOfSM[5] = (fPhiBoundariesOfSM[10]+fPhiBoundariesOfSM[11])/2.;
0c5b726e 362 }
363
364 //called after setting of scintillator and lead layer parameters
365 DefineSamplingFraction();
366
367
368 // TRU parameters - Apr 29,08 by PAI.
369 // These parameters values was updated at Nov 05, 2007
370 // As is on Olivier BOURRION (LPSC) ppt preasentation
371 // at ALICE trigger meeting at 13th-14th March
372 fNTRUEta = 1; // was 3
373 fNTRUPhi = 3; // was 1
374 fNModulesInTRUEta = 24; // was 8
375 fNModulesInTRUPhi = 4; // was 12
376 // Jet trigger
377 // 3*6*10 + 2*6*2 = 204 -> matrix (nphi(17), neta(12))
378 fNEtaSubOfTRU = 6;
0c5b726e 379
380 fgInit = kTRUE;
381}
382
383//___________________________________________________________________
384void AliEMCALEMCGeometry::PrintGeometry()
385{
386 // Separate routine is callable from broswer; Nov 7,2006
387 printf("\nInit: geometry of EMCAL named %s :\n", fGeoName.Data());
388 if(fArrayOpts) {
389 for(Int_t i=0; i<fArrayOpts->GetEntries(); i++){
390 TObjString *o = (TObjString*)fArrayOpts->At(i);
391 printf(" %i : %s \n", i, o->String().Data());
392 }
393 }
394 printf("Granularity: %d in eta and %d in phi\n", GetNZ(), GetNPhi()) ;
395 printf("Layout: phi = (%7.1f, %7.1f), eta = (%5.2f, %5.2f), IP = %7.2f -> for EMCAL envelope only\n",
396 GetArm1PhiMin(), GetArm1PhiMax(),GetArm1EtaMin(), GetArm1EtaMax(), GetIPDistance() );
397
398 printf( " ECAL : %d x (%f cm Pb, %f cm Sc) \n",
399 GetNECLayers(), GetECPbRadThick(), GetECScintThick() ) ;
400 printf(" fSampling %5.2f \n", fSampling );
401 printf(" fIPDistance %6.3f cm \n", fIPDistance);
402 printf(" fNPhi %i | fNZ %i \n", fNPhi, fNZ);
403 printf(" fNCellsInModule %i : fNCellsInSupMod %i : fNCells %i\n",fNCellsInModule, fNCellsInSupMod, fNCells);
404 printf(" X:Y module size %6.3f , %6.3f cm \n", fPhiModuleSize, fEtaModuleSize);
405 printf(" X:Y tile size %6.3f , %6.3f cm \n", fPhiTileSize, fEtaTileSize);
406 printf(" #of sampling layers %i(fNECLayers) \n", fNECLayers);
407 printf(" fLongModuleSize %6.3f cm \n", fLongModuleSize);
408 printf(" #supermodule in phi direction %i \n", fNPhiSuperModule );
409 printf(" fILOSS %i : fIHADR %i \n", fILOSS, fIHADR);
410 printf(" fTrd1Angle %7.4f\n", fTrd1Angle);
411 printf(" f2Trd1Dx2 %7.4f\n", f2Trd1Dx2);
412 printf("SM dimensions(TRD1) : dx %7.2f dy %7.2f dz %7.2f (SMOD, BOX)\n",
413 fParSM[0],fParSM[1],fParSM[2]);
414 printf(" fPhiGapForSM %7.4f cm (%7.4f <- phi size in degree)\n",
415 fPhiGapForSM, TMath::ATan2(fPhiGapForSM,fIPDistance)*TMath::RadToDeg());
416 if(GetKey110DEG()) printf(" Last two modules have size 10 degree in phi (180<phi<190)\n");
417 printf(" phi SM boundaries \n");
418 for(int i=0; i<fPhiBoundariesOfSM.GetSize()/2.; i++) {
419 printf(" %i : %7.5f(%7.2f) -> %7.5f(%7.2f) : center %7.5f(%7.2f) \n", i,
420 fPhiBoundariesOfSM[2*i], fPhiBoundariesOfSM[2*i]*TMath::RadToDeg(),
421 fPhiBoundariesOfSM[2*i+1], fPhiBoundariesOfSM[2*i+1]*TMath::RadToDeg(),
422 fPhiCentersOfSM[i], fPhiCentersOfSM[i]*TMath::RadToDeg());
423 }
424
425}
426
427//______________________________________________________________________
428void AliEMCALEMCGeometry::CheckAdditionalOptions()
429{
430 // Feb 06,2006
431 // Additional options that
432 // can be used to select
433 // the specific geometry of
434 // EMCAL to run
435 // Dec 27,2006
436 // adeed allILOSS= and allIHADR= for MIP investigation
437 fArrayOpts = new TObjArray;
438 Int_t nopt = ParseString(fGeoName, *fArrayOpts);
439 if(nopt==1) { // no aditional option(s)
440 fArrayOpts->Delete();
441 delete fArrayOpts;
442 fArrayOpts = 0;
443 return;
444 }
445 for(Int_t i=1; i<nopt; i++){
446 TObjString *o = (TObjString*)fArrayOpts->At(i);
447
448 TString addOpt = o->String();
449 Int_t indj=-1;
450 for(Int_t j=0; j<fNAdditionalOpts; j++) {
451 TString opt = fkAdditionalOpts[j];
452 if(addOpt.Contains(opt,TString::kIgnoreCase)) {
453 indj = j;
454 break;
455 }
456 }
457 if(indj<0) {
458 AliDebug(2,Form("<E> option |%s| unavailable : ** look to the file AliEMCALGeometry.h **\n",
459 addOpt.Data()));
460 assert(0);
461 } else {
462 AliDebug(2,Form("<I> option |%s| is valid : number %i : |%s|\n",
463 addOpt.Data(), indj, fkAdditionalOpts[indj]));
464 if (addOpt.Contains("NL=",TString::kIgnoreCase)) {// number of sampling layers
465 sscanf(addOpt.Data(),"NL=%i", &fNECLayers);
466 AliDebug(2,Form(" fNECLayers %i (new) \n", fNECLayers));
467 } else if(addOpt.Contains("PBTH=",TString::kIgnoreCase)) {//Thickness of the Pb(fECPbRadThicknes)
468 sscanf(addOpt.Data(),"PBTH=%f", &fECPbRadThickness);
469 } else if(addOpt.Contains("SCTH=",TString::kIgnoreCase)) {//Thickness of the Sc(fECScintThick)
470 sscanf(addOpt.Data(),"SCTH=%f", &fECScintThick);
471 } else if(addOpt.Contains("LATSS=",TString::kIgnoreCase)) {// Thickness of lateral steel strip (fLateralSteelStrip)
472 sscanf(addOpt.Data(),"LATSS=%f", &fLateralSteelStrip);
473 AliDebug(2,Form(" fLateralSteelStrip %f (new) \n", fLateralSteelStrip));
474 } else if(addOpt.Contains("ILOSS=",TString::kIgnoreCase)) {// As in Geant
475 sscanf(addOpt.Data(),"ALLILOSS=%i", &fILOSS);
476 AliDebug(2,Form(" fILOSS %i \n", fILOSS));
477 } else if(addOpt.Contains("IHADR=",TString::kIgnoreCase)) {// As in Geant
478 sscanf(addOpt.Data(),"ALLIHADR=%i", &fIHADR);
479 AliDebug(2,Form(" fIHADR %i \n", fIHADR));
480 }
481 }
482 }
483}
484
485//__________________________________________________________________
486void AliEMCALEMCGeometry::DefineSamplingFraction()
487{
488 // Jun 05,2006
489 // Look http://rhic.physics.wayne.edu/~pavlinov/ALICE/SHISHKEBAB/RES/linearityAndResolutionForTRD1.html
490 // Keep for compatibilty
491 //
492 if(fNECLayers == 69) { // 10% layer reduction
493 fSampling = 12.55;
494 } else if(fNECLayers == 61) { // 20% layer reduction
495 fSampling = 12.80;
496 } else if(fNECLayers == 77) {
497 if (fECScintThick>0.159 && fECScintThick<0.161) { // original sampling fraction, equal layers
498 fSampling = 12.327; // fECScintThick = fECPbRadThickness = 0.160;
499 } else if (fECScintThick>0.175 && fECScintThick<0.177) { // 10% Pb thicknes reduction
500 fSampling = 10.5; // fECScintThick = 0.176, fECPbRadThickness=0.144;
501 } else if(fECScintThick>0.191 && fECScintThick<0.193) { // 20% Pb thicknes reduction
502 fSampling = 8.93; // fECScintThick = 0.192, fECPbRadThickness=0.128;
503 }
504
505 }
506}
507
508//________________________________________________________________________________________________
509Double_t AliEMCALEMCGeometry::GetPhiCenterOfSM(Int_t nsupmod) const
510{
511 //returns center of supermodule in phi
512 int i = nsupmod/2;
513 return fPhiCentersOfSM[i];
514
515}
516
517//________________________________________________________________________________________________
518Bool_t AliEMCALEMCGeometry::GetPhiBoundariesOfSM(Int_t nSupMod, Double_t &phiMin, Double_t &phiMax) const
519{
520 // 0<= nSupMod <=11; phi in rad
521 static int i;
522 if(nSupMod<0 || nSupMod >11) return kFALSE;
523 i = nSupMod/2;
524 phiMin = (Double_t)fPhiBoundariesOfSM[2*i];
525 phiMax = (Double_t)fPhiBoundariesOfSM[2*i+1];
526 return kTRUE;
527}
528
529//________________________________________________________________________________________________
530Bool_t AliEMCALEMCGeometry::GetPhiBoundariesOfSMGap(Int_t nPhiSec, Double_t &phiMin, Double_t &phiMax) const
531{
532 // 0<= nPhiSec <=4; phi in rad
533 // 0; gap boundaries between 0th&2th | 1th&3th SM
534 // 1; gap boundaries between 2th&4th | 3th&5th SM
535 // 2; gap boundaries between 4th&6th | 5th&7th SM
536 // 3; gap boundaries between 6th&8th | 7th&9th SM
537 // 4; gap boundaries between 8th&10th | 9th&11th SM
538 if(nPhiSec<0 || nPhiSec >4) return kFALSE;
539 phiMin = fPhiBoundariesOfSM[2*nPhiSec+1];
540 phiMax = fPhiBoundariesOfSM[2*nPhiSec+2];
541 return kTRUE;
542}
543
544//________________________________________________________________________________________________
545int AliEMCALEMCGeometry::ParseString(const TString &topt, TObjArray &Opt)
546{
547 //Parse string, does what? GCB 08/09
548 Ssiz_t begin, index, end, end2;
549 begin = index = end = end2 = 0;
550 TRegexp separator("[^ ;,\\t\\s/]+");
551 while ( (begin < topt.Length()) && (index != kNPOS) ) {
552 // loop over given options
553 index = topt.Index(separator,&end,begin);
554 if (index >= 0 && end >= 1) {
555 TString substring(topt(index,end));
556 Opt.Add(new TObjString(substring.Data()));
557 }
558 begin += end+1;
559 }
560 return Opt.GetEntries();
561}
562