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