7562371111111a9251e5f206273ead453ec1f19f
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALv0.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15 /* $Id$ */
16
17 //_________________________________________________________________________
18 // Implementation version v0 of EMCAL Manager class 
19 // An object of this class does not produce hits nor digits
20 // It is the one to use if you do not want to produce outputs in TREEH or TREED
21 // This class places a Geometry of the EMCAL in the ALICE Detector as defined in AliEMCALGeometry.cxx                 
22 //*-- Author: Yves Schutz (SUBATECH)
23 //*-- and   : Sahal Yacoob (LBL / UCT)
24 //          : Aleksei Pavlinov (WSU)     SHASHLYK
25
26 // --- ROOT system ---
27 #include <cassert>
28
29 #include <TGeometry.h>
30 #include <TGeoPhysicalNode.h>
31 #include <TGeoManager.h>
32 #include <TGeoMatrix.h>
33 #include <TVirtualMC.h>
34 #include <TArrayI.h>
35 #include <TROOT.h>
36 #include <TList.h>
37 #include <TVector2.h>
38 #include <cassert>
39
40 //--- EMCAL system---
41 #include "AliEMCALShishKebabTrd1Module.h"
42
43 // --- Standard library ---
44
45 //#include <stdio.h>
46
47 // --- AliRoot header files ---
48
49 #include "AliEMCALv0.h"
50 #include "AliEMCALGeometry.h"
51 #include "AliRun.h"
52 #include "AliLog.h"
53 #include "AliGeomManager.h"
54 #include "AliEMCALSpaceFrame.h"
55
56 ClassImp(AliEMCALv0)
57
58 // EMCAL material
59 enum
60  {
61   kIdAIR   = 1599, 
62   kIdPB    = 1600, 
63   kIdSC    = 1601, 
64   kIdSTEEL = 1603
65  };
66
67
68 //______________________________________________________________________
69 AliEMCALv0::AliEMCALv0()
70   : AliEMCAL(),
71     fShishKebabModules(),fEnvelop1(0),fIdRotm(0),fIdTmedArr(0),
72     fSampleWidth(0),fSmodPar0(0),fSmodPar1(0),fSmodPar2(0),fCalFrame(0)
73 {
74   //default ctor
75   for(Int_t i = 0; i < 5 ; i++) fParEMOD[i]=0.0;
76 }
77
78 //______________________________________________________________________
79 AliEMCALv0::AliEMCALv0(const char *name, const char *title)
80   : AliEMCAL(name,title),
81     fShishKebabModules(),fEnvelop1(0),fIdRotm(0),fIdTmedArr(0),
82     fSampleWidth(0),fSmodPar0(0),fSmodPar1(0),fSmodPar2(0),fCalFrame(0)
83 {
84   // ctor : title is used to identify the layout
85   // Apr 25, 2006
86   // Nov 22, 2006 - case of 1X1  
87   
88   for(Int_t i = 0; i < 5 ; i++) fParEMOD[i]=0.0;
89   TString ntmp(GetTitle());
90   ntmp.ToUpper();
91
92   AliEMCALGeometry *geom = GetGeometry() ; 
93   fShishKebabModules = geom->GetShishKebabTrd1Modules(); 
94 }
95
96 //______________________________________________________________________
97 void AliEMCALv0::CreateGeometry()
98 {
99   // Create the EMCAL geometry for Geant
100   // Geometry of a tower
101   
102   AliEMCALGeometry * geom = GetGeometry() ; 
103   TString gn(geom->GetName());
104   gn.ToUpper(); 
105   
106   if(!(geom->IsInitialized())){
107     Error("CreateGeometry","EMCAL Geometry class has not been set up.");
108   } // end if
109
110   // Get pointer to the array containing media indices
111   fIdTmedArr = fIdtmed->GetArray() - 1599 ;
112   
113   fIdRotm = 1;
114   //  gMC->Matrix(nmat, theta1, phi1, theta2, phi2, theta3, phi3) - see AliModule
115   AliMatrix(fIdRotm, 90.0, 0., 90.0, 90.0, 0.0, 0.0) ; 
116   
117   // Create the EMCAL Mother Volume (a polygone) within which to place the Detector and named XEN1 
118   
119   Float_t envelopA[10];
120   if(gn.Contains("WSUC") ) { // TRD1 for WSUC facility
121     // 17-may-05 - just BOX
122     envelopA[0] = 26;
123     envelopA[1] = 15;
124     envelopA[2] = 30;
125     gMC->Gsvolu("XEN1", "BOX", fIdTmedArr[kIdSC], envelopA, 3) ;
126     fEnvelop1.Set(3);
127     for(int i=0; i<3; i++) fEnvelop1[i] = envelopA[i]; // 23-may-05  
128     // Position the EMCAL Mother Volume (XEN1) in WSUC  
129     gMC->Gspos("XEN1", 1, "WSUC", 0.0, 0.0, 0.0, fIdRotm, "ONLY") ;
130   } else { 
131     envelopA[0] = geom->GetArm1PhiMin();                         // minimum phi angle
132     envelopA[1] = geom->GetArm1PhiMax() - geom->GetArm1PhiMin(); // angular range in phi
133     envelopA[2] = geom->GetNPhiSuperModule();                    // number of sections in phi
134     envelopA[3] = 2;                                             // 2 z coordinates
135     envelopA[4] = -geom->GetEnvelop(2)/2.;                       // zmin - includes padding
136     envelopA[5] = geom->GetEnvelop(0) ;                          // rmin at z1 - includes padding
137     envelopA[6] = geom->GetEnvelop(1) ;                          // rmax at z1 - includes padding
138     envelopA[7] = geom->GetEnvelop(2)/2.;                        // zmax includes padding
139     
140     envelopA[8] = envelopA[5] ;                                  // radii are the same.
141     envelopA[9] = envelopA[6] ;                                  // radii are the same.
142
143     gMC->Gsvolu("XEN1", "PGON", fIdTmedArr[kIdAIR], envelopA, 10) ;   // Polygone filled with air 
144     fEnvelop1.Set(10, envelopA);
145     if (gDebug==2) {
146       printf("CreateGeometry: XEN1 = %f, %f\n", envelopA[5], envelopA[6]); 
147       printf("CreateGeometry: XU0 = %f, %f\n", envelopA[5], envelopA[6]); 
148     }
149     // Position the EMCAL Mother Volume (XEN1) in Alice (ALIC)  
150     gMC->Gspos(geom->GetNameOfEMCALEnvelope(), 1, "ALIC", 0.0, 0.0, 0.0, fIdRotm, "ONLY") ;
151   }
152
153   // COMPACT, TRD1
154   AliDebug(2,Form("Shish-Kebab geometry : %s", GetTitle())); 
155   CreateShishKebabGeometry();
156
157   //Space Frame
158   AliDebug(2,"Creating EMCAL Space Frame");
159   fCalFrame = new AliEMCALSpaceFrame();
160   fCalFrame->CreateGeometry();
161
162 }
163
164 //______________________________________________________________________
165 void AliEMCALv0::Init(void)
166 {
167     // Just prints an information message
168   AliEMCAL::Init();
169   if(AliLog::GetGlobalDebugLevel()>0) { 
170     TString message("\n") ; 
171     message += "*****************************************\n" ;
172     
173     // Here the EMCAL initialisation code (if any!)
174     
175     AliEMCALGeometry * geom = GetGeometry() ; 
176     
177     if (geom!=0) {   
178       message += "AliEMCAL " ; 
179       message += Version() ; 
180       message += "EMCAL geometry initialized for " ; 
181       message += geom->GetName()  ;
182     }
183     else {
184       message += "AliEMCAL " ; 
185       message += Version() ;  
186       message += "EMCAL geometry initialization failed !" ; 
187     }
188     message += "\n*****************************************" ;
189     printf("%s",message.Data() ) ; 
190   }
191 }
192
193 // 24-aug-04 by PAI
194 //______________________________________________________________________
195 void AliEMCALv0::CreateShishKebabGeometry()
196 {  
197   // TRD1
198   AliEMCALGeometry * g = GetGeometry(); 
199   TString gn(g->GetName()); gn.ToUpper(); 
200   // see AliModule::fFIdTmedArr
201   //  fIdTmedArr = fIdtmed->GetArray() - 1599 ; // see AliEMCAL::::CreateMaterials()
202   //  int kIdAIR=1599, kIdPB = 1600, kIdSC = 1601, kIdSTEEL = 1603;
203   //  idAL = 1602;
204   Double_t par[10], xpos=0., ypos=0., zpos=0.;
205
206   CreateSmod(g->GetNameOfEMCALEnvelope());
207
208   CreateEmod("SMOD","EMOD"); // 18-may-05
209
210   if(g->GetKey110DEG()) CreateEmod("SM10","EMOD"); // Nov 1,2006 
211
212   // Sensitive SC  (2x2 tiles)
213   double parSCM0[5]={0,0,0,0}, *dummy = 0, parTRAP[11];
214
215   //      double wallThickness = g->GetPhiModuleSize()/2. -  g->GetPhiTileSize(); // Need check
216   double wallThickness = g->GetPhiModuleSize()/g->GetNPHIdiv() -  g->GetPhiTileSize();
217   for(int i=0; i<3; i++) parSCM0[i] = fParEMOD[i] - wallThickness;
218   parSCM0[3] = fParEMOD[3];
219   gMC->Gsvolu("SCM0", "TRD1", fIdTmedArr[kIdAIR], parSCM0, 4);
220   gMC->Gspos("SCM0", 1, "EMOD", 0., 0., 0., 0, "ONLY") ;
221
222   if(g->GetNPHIdiv()==2 && g->GetNETAdiv()==2) {
223     // Division to tile size - 1-oct-04
224     AliDebug(2,Form(" Divide SCM0 on y-axis %i\n", g->GetNETAdiv()));
225     gMC->Gsdvn("SCMY","SCM0", g->GetNETAdiv(), 2); // y-axis
226     // Trapesoid 2x2
227     parTRAP[0] = parSCM0[3];    // dz
228     parTRAP[1] = TMath::ATan2((parSCM0[1]-parSCM0[0])/2.,2.*parSCM0[3])*180./TMath::Pi(); // theta
229     parTRAP[2] = 0.;           // phi
230     // bottom
231     parTRAP[3] = parSCM0[2]/2.; // H1
232     parTRAP[4] = parSCM0[0]/2.; // BL1
233     parTRAP[5] = parTRAP[4];    // TL1
234     parTRAP[6] = 0.0;           // ALP1
235     // top
236     parTRAP[7] = parSCM0[2]/2.; // H2
237     parTRAP[8] = parSCM0[1]/2.; // BL2
238     parTRAP[9] = parTRAP[8];    // TL2
239     parTRAP[10]= 0.0;           // ALP2
240     AliDebug(2,Form(" ** TRAP ** \n"));
241     for(int i=0; i<11; i++) AliDebug(3, Form(" par[%2.2i] %9.4f\n", i, parTRAP[i]));
242     
243     gMC->Gsvolu("SCMX", "TRAP", fIdTmedArr[kIdSC], parTRAP, 11);
244     xpos = +(parSCM0[1]+parSCM0[0])/4.;
245     gMC->Gspos("SCMX", 1, "SCMY", xpos, 0.0, 0.0, 0, "ONLY") ;
246
247     // Using rotation because SCMX should be the same due to Pb tiles
248     xpos = -xpos; 
249     AliMatrix(fIdRotm, 90.0,180., 90.0, 270.0, 0.0,0.0) ;
250     gMC->Gspos("SCMX", 2, "SCMY", xpos, 0.0, 0.0, fIdRotm, "ONLY");
251     // put LED to the SCM0 
252     AliEMCALShishKebabTrd1Module *mod = (AliEMCALShishKebabTrd1Module*)fShishKebabModules->At(0);
253     gMC->Gsvolu("PBTI", "BOX", fIdTmedArr[kIdPB], dummy, 0);
254     
255     par[1] = parSCM0[2]/2;            // y 
256     par[2] = g->GetECPbRadThick()/2.; // z
257
258     int nr=0;
259     ypos = 0.0; 
260     zpos = -fSampleWidth*g->GetNECLayers()/2. + g->GetECPbRadThick()/2.;
261     double xCenterSCMX =  (parTRAP[4] +  parTRAP[8])/2.;
262     AliDebug(2,Form(" Pb tiles \n"));
263     for(int iz=0; iz<g->GetNECLayers(); iz++){
264       par[0] = (parSCM0[0] + mod->GetTanBetta()*fSampleWidth*iz)/2.;
265       xpos   = par[0] - xCenterSCMX;
266       gMC->Gsposp("PBTI", ++nr, "SCMX", xpos, ypos, zpos, 0, "ONLY", par, 3) ;
267       AliDebug(3,Form(" %i xpos %f zpos %f par[0] %f \n", iz+1, xpos, zpos, par[0]));
268       zpos += fSampleWidth;
269     } 
270     AliDebug(2,Form(" Number of Pb tiles in SCMX %i \n", nr));
271     
272   } else if(g->GetNPHIdiv()==3 && g->GetNETAdiv()==3) {
273     printf(" before AliEMCALv0::Trd1Tower3X3() : parSCM0");
274     for(int i=0; i<4; i++) printf(" %7.4f ", parSCM0[i]);
275     printf("\n"); 
276     Trd1Tower3X3(parSCM0);
277   } else if(g->GetNPHIdiv()==1 && g->GetNETAdiv()==1) {
278     // no division in SCM0
279     Trd1Tower1X1(parSCM0);
280   } else if(g->GetNPHIdiv()==4 && g->GetNETAdiv()==4) {
281     Trd1Tower4X4();
282   }
283   
284 }
285
286 //______________________________________________________________________
287 void AliEMCALv0::CreateSmod(const char* mother)
288
289   // 18-may-05; mother="XEN1"; 
290   // child="SMOD" from first to 10th, "SM10" (11th and 12th) (TRD1 case)
291   AliEMCALGeometry * g = GetGeometry(); 
292   TString gn(g->GetName()); gn.ToUpper();
293
294   Double_t par[3], xpos=0., ypos=0., zpos=0., rpos=0., dphi=0., phi=0.0, phiRad=0.;
295   Double_t par1C = 0.;
296   //  ===== define Super Module from air - 14x30 module ==== ;
297   fSampleWidth = double(g->GetECPbRadThick()+g->GetECScintThick());
298   AliDebug(2,Form("\n ## Super Module | fSampleWidth %5.3f ## %s \n", fSampleWidth, gn.Data()));
299   par[0] = g->GetShellThickness()/2.;
300   par[1] = g->GetPhiModuleSize()*g->GetNPhi()/2.; 
301   par[2] = g->GetEtaModuleSize()*15.; 
302   fIdRotm=0;
303   int nphism = g->GetNumberOfSuperModules()/2; // 20-may-05
304   if(nphism>0) {
305     dphi = (g->GetArm1PhiMax() - g->GetArm1PhiMin())/nphism;
306     //    if(g->GetKey110DEG()) dphi = (g->GetArm1PhiMax() - g->GetArm1PhiMin())/(nphism-1);
307     rpos = (g->GetEnvelop(0) + g->GetEnvelop(1))/2.;
308     AliDebug(2,Form(" rpos %8.2f : dphi %6.1f degree \n", rpos, dphi));
309   }
310
311   if(gn.Contains("WSUC")) {
312     par[0] = g->GetPhiModuleSize()*g->GetNPhi()/2.; 
313     par[1] = g->GetShellThickness()/2.;
314     par[2] = g->GetEtaModuleSize()*g->GetNZ()/2. + 5; 
315
316     gMC->Gsvolu("SMOD", "BOX", fIdTmedArr[kIdAIR], par, 3);
317
318     AliDebug(2,Form("SMOD in WSUC : tmed %i | dx %7.2f dy %7.2f dz %7.2f (SMOD, BOX)\n", 
319                     fIdTmedArr[kIdAIR], par[0],par[1],par[2]));
320     fSmodPar0 = par[0]; 
321     fSmodPar1 = par[1];
322     fSmodPar2 = par[2];
323     nphism   =  g->GetNumberOfSuperModules();
324   } else {
325     par[2]  = 350./2.; // 11-oct-04 - for 26 division
326     AliDebug(2,Form(" par[0] %7.2f (old) \n",  par[0]));
327     Float_t *parSM = g->GetSuperModulesPars(); 
328     for(int i=0; i<3; i++) par[i] = parSM[i];
329   }
330   gMC->Gsvolu("SMOD", "BOX", fIdTmedArr[kIdAIR], par, 3);
331   AliDebug(2,Form("tmed %i | dx %7.2f dy %7.2f dz %7.2f (SMOD, BOX)\n", 
332                   fIdTmedArr[kIdAIR], par[0],par[1],par[2]));
333   fSmodPar0 = par[0]; 
334   fSmodPar2 = par[2];
335   if(g->GetKey110DEG()) { // 12-oct-05
336     par1C = par[1];
337     par[1] /= 2.;
338     gMC->Gsvolu("SM10", "BOX", fIdTmedArr[kIdAIR], par, 3);
339     AliDebug(2,Form(" Super module with name \"SM10\" was created too par[1] = %f\n", par[1]));
340     par[1] = par1C;
341   }
342   // Steel plate
343   if(g->GetSteelFrontThickness() > 0.0) { // 28-mar-05
344     par[0] = g->GetSteelFrontThickness()/2.;
345     gMC->Gsvolu("STPL", "BOX", fIdTmedArr[kIdSTEEL], par, 3);
346     printf("tmed %i | dx %7.2f dy %7.2f dz %7.2f (STPL) \n", fIdTmedArr[kIdSTEEL], par[0],par[1],par[2]);
347     xpos = -(g->GetShellThickness() - g->GetSteelFrontThickness())/2.;
348     gMC->Gspos("STPL", 1, "SMOD", xpos, 0.0, 0.0, 0, "ONLY") ;
349   }
350
351   int nr=0, nrsmod=0, i0=0;
352
353   // Turn whole super module
354   for(int i=i0; i<nphism; i++) {
355     if(gn.Contains("WSUC")) {
356       xpos = ypos = zpos = 0.0;
357       fIdRotm = 0;
358       gMC->Gspos("SMOD", 1, mother, xpos, ypos, zpos, fIdRotm, "ONLY") ;
359       printf(" fIdRotm %3i phi %6.1f(%5.3f) xpos %7.2f ypos %7.2f zpos %7.2f \n", 
360       fIdRotm, phi, phiRad, xpos, ypos, zpos);
361       nr++;
362     } else { // TRD1 
363       TString smName("SMOD"); // 12-oct-05
364       if(i==5 && g->GetKey110DEG()) {
365         smName = "SM10";
366         nrsmod = nr;
367         nr     = 0;
368       }
369       phi    = g->GetArm1PhiMin() + dphi*(2*i+1)/2.; // phi= 70, 90, 110, 130, 150, 170
370       phiRad = phi*TMath::Pi()/180.;
371
372       AliMatrix(fIdRotm, 90.0, phi, 90.0, 90.0+phi, 0.0, 0.0);
373
374       xpos = rpos * TMath::Cos(phiRad);
375       ypos = rpos * TMath::Sin(phiRad);
376       zpos = fSmodPar2; // 21-sep-04
377       if(i==5 && g->GetKey110DEG()) {
378         xpos += (par1C/2. * TMath::Sin(phiRad)); 
379         ypos -= (par1C/2. * TMath::Cos(phiRad)); 
380       }
381       
382       // 1th module in z-direction;
383       gMC->Gspos(smName.Data(), ++nr, mother, xpos, ypos, zpos, fIdRotm, "ONLY") ;
384       AliDebug(3, Form(" %s : %2i fIdRotm %3i phi %6.1f(%5.3f) xpos %7.2f ypos %7.2f zpos %7.2f : i %i \n", 
385                        smName.Data(), nr, fIdRotm, phi, phiRad, xpos, ypos, zpos, i));
386       // 2th module in z-direction;
387       // turn arround X axis; 0<phi<360
388       double phiy = 90. + phi + 180.;
389       if(phiy>=360.) phiy -= 360.;
390       
391       AliMatrix(fIdRotm, 90.0, phi, 90.0, phiy, 180.0, 0.0);
392       gMC->Gspos(smName.Data(), ++nr, mother, xpos, ypos, -zpos, fIdRotm, "ONLY");
393       AliDebug(3, Form(" %s : %2i fIdRotm %3i phiy %6.1f  xpos %7.2f ypos %7.2f zpos %7.2f \n", 
394                        smName.Data(), nr, fIdRotm, phiy, xpos, ypos, -zpos));
395     }
396   }
397   AliDebug(2,Form(" Number of Super Modules %i \n", nr+nrsmod));
398 }
399
400 //______________________________________________________________________
401 void AliEMCALv0::CreateEmod(const char* mother, const char* child)
402
403   // 17-may-05; mother="SMOD"; child="EMOD"
404   AliEMCALGeometry * g = GetGeometry(); 
405   TString gn(g->GetName()); gn.ToUpper(); 
406   // Module definition
407   Double_t xpos=0., ypos=0., zpos=0.;
408   Double_t parSCPA[5], zposSCPA=0.; // passive SC - 13-MAY-05, TRD1 case
409   Double_t trd1Angle = g->GetTrd1Angle()*TMath::DegToRad(), tanTrd1 = TMath::Tan(trd1Angle/2.);
410   int nr=0;
411   fIdRotm=0;
412   if(strcmp(mother,"SMOD")==0) {
413     fParEMOD[0] = g->GetEtaModuleSize()/2.;   // dx1
414     fParEMOD[1] = g->Get2Trd1Dx2()/2.;        // dx2
415     fParEMOD[2] = g->GetPhiModuleSize()/2.;;  // dy
416     fParEMOD[3] = g->GetLongModuleSize()/2.;  // dz
417     gMC->Gsvolu(child, "TRD1", fIdTmedArr[kIdSTEEL], fParEMOD, 4);
418     //      if(gn.Contains("WSUC") || gn.Contains("MAY05")){
419     if(0){ // Jul 12 - should be checked
420       parSCPA[0] = g->GetEtaModuleSize()/2. + tanTrd1*g->GetFrontSteelStrip();   // dx1
421       parSCPA[1] = parSCPA[0]               + tanTrd1*g->GetPassiveScintThick(); // dx2
422       parSCPA[2] = g->GetPhiModuleSize()/2.;     // dy
423       parSCPA[3] = g->GetPassiveScintThick()/2.; // dz
424       gMC->Gsvolu("SCPA", "TRD1", fIdTmedArr[kIdSC], parSCPA, 4);
425       zposSCPA   = -fParEMOD[3] + g->GetFrontSteelStrip() + g->GetPassiveScintThick()/2.;
426       gMC->Gspos ("SCPA", ++nr, child, 0.0, 0.0, zposSCPA, 0, "ONLY");
427     }
428   }
429   
430   nr   = 0;
431   // X->Z(0, 0); Y->Y(90, 90); Z->X(90, 0)
432   AliEMCALShishKebabTrd1Module *mod=0; // current module
433
434   for(int iz=0; iz<g->GetNZ(); iz++) {
435     Double_t  angle=90., phiOK=0;
436     mod = (AliEMCALShishKebabTrd1Module*)fShishKebabModules->At(iz);
437     angle = mod->GetThetaInDegree();
438     if(!gn.Contains("WSUC")) { // ALICE 
439       if(iz==0) AliMatrix(fIdRotm, 0.,0., 90.,90., 90.,0.); // z'(x); y'(y); x'(z)
440       else      AliMatrix(fIdRotm, 90.-angle,180., 90.0,90.0, angle, 0.);
441       phiOK = mod->GetCenterOfModule().Phi()*180./TMath::Pi(); 
442       //          printf(" %2i | angle | %6.3f - %6.3f = %6.3f(eta %5.3f)\n", 
443       //iz+1, angle, phiOK, angle-phiOK, mod->GetEtaOfCenterOfModule());
444       xpos = mod->GetPosXfromR() + g->GetSteelFrontThickness() - fSmodPar0;
445       zpos = mod->GetPosZ() - fSmodPar2;
446       
447       int iyMax = g->GetNPhi();
448       if(strcmp(mother,"SMOD") && g->GetKey110DEG()) {
449         iyMax /= 2;
450       }
451       for(int iy=0; iy<iyMax; iy++) { // flat in phi
452         ypos = g->GetPhiModuleSize()*(2*iy+1 - iyMax)/2.;
453         gMC->Gspos(child, ++nr, mother, xpos, ypos, zpos, fIdRotm, "ONLY") ;
454         //printf(" %2i xpos %7.2f ypos %7.2f zpos %7.2f fIdRotm %i\n", nr, xpos, ypos, zpos, fIdRotm);
455         AliDebug(3,Form("%3.3i(%2.2i,%2.2i) ", nr,iy+1,iz+1));
456       }
457       //PH          printf("\n");
458     } else { //WSUC
459       if(iz==0) AliMatrix(fIdRotm, 0.,0., 90.,0., 90.,90.); // (x')z; y'(x); z'(y)
460       else      AliMatrix(fIdRotm, 90-angle,270., 90.0,0.0, angle,90.);
461       phiOK = mod->GetCenterOfModule().Phi()*180./TMath::Pi(); 
462       printf(" %2i | angle -phiOK | %6.3f - %6.3f = %6.3f(eta %5.3f)\n", 
463              iz+1, angle, phiOK, angle-phiOK, mod->GetEtaOfCenterOfModule());
464       zpos = mod->GetPosZ()      - fSmodPar2;
465       ypos = mod->GetPosXfromR() - fSmodPar1;
466       printf(" zpos %7.2f ypos %7.2f fIdRotm %i\n xpos ", zpos, xpos, fIdRotm);
467       for(int ix=0; ix<g->GetNPhi(); ix++) { // flat in phi
468         xpos = g->GetPhiModuleSize()*(2*ix+1 - g->GetNPhi())/2.;
469         gMC->Gspos(child, ++nr, mother, xpos, ypos, zpos, fIdRotm, "ONLY") ;
470         printf(" %7.2f ", xpos);
471       }
472       printf("\n");
473     }
474   }
475   AliDebug(2,Form(" Number of modules in Super Module(%s) %i \n", mother, nr));
476 }
477
478 //______________________________________________________________________
479 void AliEMCALv0::Trd1Tower3X3(const double *parSCM0)
480 {
481   // Started Dec 8,2004 by PAI
482   // Fixed Nov 13,2006
483   printf(" AliEMCALv0::Trd1Tower3X3() : parSCM0");
484   for(int i=0; i<4; i++) printf(" %7.4f ", parSCM0[i]);
485   printf("\n"); 
486   // Nov 10, 2006 - different name of SCMX
487   double parTRAP[11], *dummy=0;
488   AliEMCALGeometry * g = GetGeometry(); 
489   TString gn(g->GetName()), scmx; 
490   gn.ToUpper(); 
491  // Division to tile size 
492   AliDebug(2,Form("Trd1Tower3X3() : Divide SCM0 on y-axis %i", g->GetNETAdiv()));
493   gMC->Gsdvn("SCMY","SCM0", g->GetNETAdiv(), 2); // y-axis
494   double dx1=parSCM0[0], dx2=parSCM0[1], dy=parSCM0[2], dz=parSCM0[3];
495   double ndiv=3., xpos=0.0;
496   // should be defined once
497   gMC->Gsvolu("PBTI", "BOX", fIdTmedArr[kIdPB], dummy, 0);
498   for(int ix=1; ix<=3; ix++) { // 3X3
499     scmx = "SCX"; // Nov 10,2006 
500     // ix=1
501     parTRAP[0] = dz;
502     double xCentBot = 2.*dx1/3.;
503     double xCentTop = 2.*(dx2/4. + dx1/12.);
504     parTRAP[1] = TMath::ATan2((xCentTop-xCentBot),2.*dz)*TMath::RadToDeg(); // theta
505     parTRAP[2] = 0.;           // phi
506     // bottom
507     parTRAP[3] = dy/ndiv;      // H1
508     parTRAP[4] = dx1/ndiv;     // BL1
509     parTRAP[5] = parTRAP[4];   // TL1
510     parTRAP[6] = 0.0;          // ALP1
511     // top
512     parTRAP[7] = dy/ndiv;      // H2
513     parTRAP[8] = dx2/2 - dx1/6.;// BL2
514     parTRAP[9] = parTRAP[8];   // TL2
515     parTRAP[10]= 0.0;          // ALP2
516     xpos = (xCentBot+xCentTop)/2.;
517
518     if      (ix==3) {
519       parTRAP[1] = -parTRAP[1];
520       xpos = -xpos;
521     } else if(ix==2) { // central part is box but we treat as trapesoid due to numbering
522       parTRAP[1] = 0.;
523       parTRAP[8] = dx1/ndiv;     // BL2
524       parTRAP[9] = parTRAP[8];   // TL2
525       xpos = 0.0;
526     }
527     AliDebug(2,Form(" ** TRAP ** xpos %9.3f\n", xpos));
528     for(int i=0; i<11; i++) AliDebug(2,Form(" par[%2.2i] %9.4f\n", i, parTRAP[i]));
529
530     scmx += ix;
531     gMC->Gsvolu(scmx.Data(), "TRAP", fIdTmedArr[kIdSC], parTRAP, 11);
532     gMC->Gspos(scmx.Data(), 1, "SCMY", xpos, 0.0, 0.0, 0, "ONLY") ;
533
534     PbInTrap(parTRAP, scmx);
535   }
536   AliDebug(2,"Trd1Tower3X3 - Ver. 1.0 : was tested.");
537 }
538
539 // 8-dec-04 by PAI
540 //______________________________________________________________________
541 void AliEMCALv0::PbInTrap(const double parTRAP[11], TString n)
542 {
543  // see for example CreateShishKebabGeometry(); just for case TRD1
544   static int nr=0;
545   AliDebug(2,Form(" Pb tiles : nrstart %i\n", nr));
546   AliEMCALGeometry * g = GetGeometry(); 
547
548   double par[3];
549   //  double fSampleWidth = double(g->GetECPbRadThick()+g->GetECScintThick());
550   double xpos = 0.0, ypos = 0.0;
551   double zpos = -fSampleWidth*g->GetNECLayers()/2. + g->GetECPbRadThick()/2.;
552
553   double coef = (parTRAP[8] -  parTRAP[4]) / (2.*parTRAP[0]);
554   double xCenterSCMX =  (parTRAP[4] +  parTRAP[8])/2.; // ??
555   //  double tan = TMath::Tan(parTRAP[1]*TMath::DegToRad());
556
557   par[1] = parTRAP[3];              // y 
558   par[2] = g->GetECPbRadThick()/2.; // z
559   for(int iz=0; iz<g->GetNECLayers(); iz++){
560     par[0] = parTRAP[4] + coef*fSampleWidth*iz;
561     xpos   = par[0] - xCenterSCMX;
562     if(parTRAP[1] < 0.) xpos = -xpos;
563     gMC->Gsposp("PBTI", ++nr, n.Data(), xpos, ypos, zpos, 0, "ONLY", par, 3) ;
564     AliDebug(2,Form(" %i xpos %9.3f zpos %9.3f par[0] %9.3f |", iz+1, xpos, zpos, par[0]));
565     zpos += fSampleWidth;
566     if(iz%2>0) printf("\n");
567   } 
568   AliDebug(2,Form(" Number of Pb tiles in SCMX %i coef %9.7f \n", nr, coef));
569   AliDebug(2,Form(" par[1] %9.3f  par[2] %9.3f ypos %9.3f \n", par[1], par[2], ypos)); 
570   AliDebug(2,Form(" PbInTrap Ver. 1.0 : was tested."));
571 }
572
573 // 8-dec-04 by PAI
574 //______________________________________________________________________
575 void AliEMCALv0::Trd1Tower4X4() const
576 {
577  // Not ready yet
578 }
579
580 //______________________________________________________________________
581 void AliEMCALv0::Trd1Tower1X1(double *parSCM0)
582 {
583   // Started Nov 22,2006 by PAI
584   AliDebug(1," AliEMCALv0::Trd1Tower1X1() : parSCM0");
585   for(int i=0; i<4; i++) printf(" %7.4f ", parSCM0[i]);
586   printf("\n"); 
587
588   // No division - keeping the same volume logic 
589   // and as consequence the same abs is scheme
590   AliDebug(2,"Trd1Tower1X1() : Create SCMX(SCMY) as SCM0");
591
592   gMC->Gsvolu("SCMY", "TRD1", fIdTmedArr[kIdAIR], parSCM0, 4);
593   gMC->Gspos("SCMY", 1, "SCM0", 0.0, 0.0, 0.0, 0, "ONLY");
594   gMC->Gsvolu("SCMX", "TRD1", fIdTmedArr[kIdSC], parSCM0, 4);
595   gMC->Gspos("SCMX", 1, "SCMY", 0.0, 0.0, 0.0, 0, "ONLY");
596
597   // should be defined once
598   double *dummy=0;
599   gMC->Gsvolu("PBTI", "BOX", fIdTmedArr[kIdPB], dummy, 0);
600
601   PbInTrd1(parSCM0, "SCMX");
602
603   AliDebug(1,"Trd1Tower1X1() : Ver. 0.1 : was tested.");
604 }
605
606 //______________________________________________________________________
607 void AliEMCALv0::PbInTrd1(const double *parTrd1, TString n)
608 {
609  // see PbInTrap(const double parTrd1[11], TString n)
610   static int nr=0, ndeb=2;
611   AliDebug(ndeb,Form(" Pb tiles : nrstart %i\n", nr));
612   AliEMCALGeometry * g = GetGeometry(); 
613
614   double par[3];
615   //  double fSampleWidth = double(g->GetECPbRadThick()+g->GetECScintThick());
616   double xpos = 0.0, ypos = 0.0;
617   double zpos = -fSampleWidth*g->GetNECLayers()/2. + g->GetECPbRadThick()/2.;
618   double coef = (parTrd1[1] -  parTrd1[0]) / (2.*parTrd1[3]);
619
620   par[1] = parTrd1[2];              // y 
621   par[2] = g->GetECPbRadThick()/2.; // z
622
623   for(int iz=0; iz<g->GetNECLayers(); iz++){
624     par[0] = parTrd1[0] + coef*fSampleWidth*iz;
625     gMC->Gsposp("PBTI", ++nr, n.Data(), xpos, ypos, zpos, 0, "ONLY", par, 3) ;
626     AliDebug(2,Form(" %i xpos %9.3f zpos %9.3f par[0] %9.3f |", iz+1, xpos, zpos, par[0]));
627     zpos += fSampleWidth;
628     if(iz%2>0) printf("\n");
629   } 
630   AliDebug(ndeb,Form(" Number of Pb tiles in SCMX %i coef %9.7f ", nr, coef));
631   AliDebug(ndeb,Form(" PbInTrd1 Ver. 0.1 : was tested."));
632 }
633
634 // 3-feb-05
635 //______________________________________________________________________
636 void AliEMCALv0::Scm0InTrd2(const AliEMCALGeometry * g, const Double_t emodPar[5], Double_t parSCM0[5])
637 {
638   // Passive material inside the detector
639   double wallThickness = g->GetPhiModuleSize()/2. -  g->GetPhiTileSize(); //Need check
640   AliDebug(2,Form(" wall thickness %7.5f \n", wallThickness));
641   for(int i=0; i<4; i++) { // on pictures sometimes I can not see 0 -> be carefull!!
642     parSCM0[i] = emodPar[i] - wallThickness;
643     AliDebug(2,Form(" %i parSCMO %7.3f emodPar %7.3f : dif %7.3f \n", 
644                     i, parSCM0[i],emodPar[i], parSCM0[i]-emodPar[i]));
645   }
646   parSCM0[4] = emodPar[4];
647   gMC->Gsvolu("SCM0", "TRD2", fIdTmedArr[kIdSC], parSCM0, 5); // kIdAIR -> kIdSC
648   gMC->Gspos("SCM0", 1, "EMOD", 0., 0., 0., 0, "ONLY") ;
649   // Division 
650   if(g->GetNPHIdiv()==2 && g->GetNETAdiv()==2) {
651     Division2X2InScm0(g, parSCM0);
652   } else {
653     Info("Scm0InTrd2"," no division SCM0 in this geometry |%s|\n", g->GetName());
654     assert(0);
655   }
656 }
657
658 //______________________________________________________________________
659 void AliEMCALv0::Division2X2InScm0(const AliEMCALGeometry * g, const Double_t parSCM0[5])
660 {
661   // Division 2X2
662   Double_t parTRAP[11], xpos=0.,ypos=0., dx1=0.0,dx2=0.,dy1=0.0,dy2=0.,dz=0.
663   ,dr1=0.0,dr2=0.;
664   fIdRotm=0;
665
666   Info("Division2X2InScm0","Divide SCM0 on y-axis %i\n", g->GetNETAdiv());
667   TString n("SCMX"), overLapFlagSCMY("ONLY"), overLapFlagSCMX("ONLY");
668   n = "SCM0"; // for testing - 14-mar-05
669   if(n=="SCM0"){
670     PbInTrapForTrd2(parSCM0, n);
671     // overLapFlagSCMY=overLapFlagSCMX="MANY"; // do not work
672     return;
673   }
674
675   dy1 = parSCM0[2] , dy2 = parSCM0[3], dz = parSCM0[4];
676
677   parTRAP[0] = parSCM0[4];    // dz
678   parTRAP[1] = TMath::ATan2((dy2-dy1)/2.,2.*dz)*TMath::RadToDeg();
679   parTRAP[2] = 90.;           // phi
680   // bottom
681   parTRAP[3] = parSCM0[2]/2.; // H1
682   parTRAP[4] = parSCM0[0];    // BL1
683   parTRAP[5] = parTRAP[4];    // TL1
684   parTRAP[6] = 0.0;           // ALP1
685   // top
686   parTRAP[7] = parSCM0[3]/2.; // H2
687   parTRAP[8] = parSCM0[1];    // BL2
688   parTRAP[9] = parTRAP[8];    // TL2
689   parTRAP[10]= 0.0;           // ALP2
690   AliDebug(2,Form(" ** SCMY ** \n"));
691            for(int i=0; i<11; i++) AliDebug(2,Form(" par[%2.2i] %9.4f\n", i, parTRAP[i]));
692
693   fIdRotm=0;
694   gMC->Gsvolu("SCMY", "TRAP", fIdTmedArr[kIdSC], parTRAP, 11); // kIdAIR -> kIdSC
695   ypos = +(parTRAP[3]+parTRAP[7])/2.; //
696   AliDebug(2,Form(" Y shift SCMY inside SCM0 : %7.3f : opt %s\n", ypos, overLapFlagSCMY.Data())); 
697   gMC->Gspos("SCMY", 1, "SCM0", 0.0, ypos, 0.0, fIdRotm,  overLapFlagSCMY.Data()) ;
698   // Rotation SCMY around z-axis on 180 degree; x'=-x; y'=-y and z=z
699   AliMatrix(fIdRotm, 90.0,180., 90.0, 270.0, 0.0,0.0) ;
700   // We may have problem with numeration due to rotation - 4-feb-05
701   gMC->Gspos("SCMY", 2, "SCM0", 0.0, -ypos, 0.0, fIdRotm,  overLapFlagSCMY.Data()); 
702
703   Info("Division2X2InScm0","Divide SCMY on x-axis %i\n", g->GetNPHIdiv());
704   dx1 = parSCM0[0]; 
705   dx2 = parSCM0[1]; 
706   dr1=TMath::Sqrt(dx1*dx1+dy1*dy1);
707   dr2=TMath::Sqrt(dx2*dx2+dy2*dy2);
708
709   parTRAP[0] = parSCM0[4];    // dz
710   parTRAP[1] = TMath::ATan2((dr2-dr1)/2.,2.*dz)*TMath::RadToDeg(); // 
711   parTRAP[2] = 45.;           // phi
712   // bottom
713   parTRAP[3] = parSCM0[2]/2.; // H1
714   parTRAP[4] = parSCM0[0]/2.; // BL1
715   parTRAP[5] = parTRAP[4];    // TL1
716   parTRAP[6] = 0.0;           // ALP1
717   // top
718   parTRAP[7] = parSCM0[3]/2.; // H2
719   parTRAP[8] = parSCM0[1]/2;  // BL2
720   parTRAP[9] = parTRAP[8];    // TL2
721   parTRAP[10]= 0.0;           // ALP2
722   AliDebug(2,Form(" ** SCMX ** \n"));
723   for(int i=0; i<11; i++) AliDebug(2,Form(" par[%2.2i] %9.4f\n", i, parTRAP[i]));
724
725   fIdRotm=0;
726   gMC->Gsvolu("SCMX", "TRAP", fIdTmedArr[kIdSC], parTRAP, 11);
727   xpos = (parTRAP[4]+parTRAP[8])/2.;
728   AliDebug(2,Form(" X shift SCMX inside SCMX : %7.3f : opt %s\n", xpos, overLapFlagSCMX.Data())); 
729   gMC->Gspos("SCMX", 1, "SCMY", xpos, 0.0, 0.0, fIdRotm,  overLapFlagSCMX.Data()) ;
730   //  AliMatrix(fIdRotm, 90.0,270., 90.0, 0.0, 0.0,0.0); // x'=-y; y'=x; z'=z
731   AliMatrix(fIdRotm, 90.0,90., 90.0, -180.0, 0.0,0.0);     // x'=y;  y'=-x; z'=z
732   gMC->Gspos("SCMX", 2, "SCMY", -xpos, 0.0, 0.0, fIdRotm,  overLapFlagSCMX.Data()) ;
733   // PB:
734   if(n=="SCMX" && overLapFlagSCMY == "ONLY") {
735     PbInTrapForTrd2(parTRAP, n);
736   }
737
738
739 // 4-feb-05 by PAI
740 //______________________________________________________________________
741 void AliEMCALv0::PbInTrapForTrd2(const double *parTRAP, TString name)
742 {
743  // TRD2 cases
744   Double_t *dummy=0;
745   TString pbShape("BOX"), pbtiChonly("ONLY");
746   if(name=="SCM0") {
747     pbShape    = "TRD2";
748     //    pbtiChonly = "MANY";
749   }
750   gMC->Gsvolu("PBTI", pbShape.Data(), fIdTmedArr[kIdPB], dummy, 0);
751
752   int nr=0;
753   Info("PbInTrapForTrd2"," Pb tiles inside %s: shape %s :pbtiChonly %s\n nrstart %i\n", 
754   name.Data(), pbShape.Data(), pbtiChonly.Data(), nr);
755   AliEMCALGeometry * g = GetGeometry(); 
756
757   double par[5], parPB[5], parSC[5];
758   //double fSampleWidth = double(g->GetECPbRadThick()+g->GetECScintThick());
759   double xpos = 0.0, ypos = 0.0;
760   double zpos = -fSampleWidth*g->GetNECLayers()/2. + g->GetECPbRadThick()/2.;
761   if(name == "SCMX") { // common trapezoid - 11 parameters
762     double coef = (parTRAP[8] -  parTRAP[4]) / (2.*parTRAP[0]);
763     double xCenterSCMX =  (parTRAP[4] +  parTRAP[8])/2.; // the same for y
764     AliDebug(2,Form(" xCenterSCMX %8.5f : coef %8.7f \n", xCenterSCMX, coef));
765
766     par[2] = g->GetECPbRadThick()/2.; // z
767     for(int iz=0; iz<g->GetNECLayers(); iz++){
768       par[0] = parTRAP[4] + coef*fSampleWidth*iz;
769       par[1] = par[0];
770       xpos   = ypos = par[0] - xCenterSCMX;
771     //if(parTRAP[1] < 0.) xpos = -xpos;
772       gMC->Gsposp("PBTI", ++nr, name.Data(), xpos, ypos, zpos, 0, "ONLY", par, 3) ;
773       AliDebug(2,Form(" %2.2i xpos %8.5f zpos %6.3f par[0,1] %6.3f |", iz+1, xpos, zpos, par[0]));
774       if(iz%2>0) AliDebug(2,Form("\n"));
775       zpos += fSampleWidth;
776     } 
777     AliDebug(2,Form(" Number of Pb tiles in SCMX %i coef %9.7f \n", nr, coef));
778     AliDebug(2,Form(" par[1] %9.5f  par[2] %9.5f ypos %9.5f \n", par[1], par[2], ypos)); 
779   } else if(name == "SCM0") { // 1-mar-05 ; TRD2 - 5 parameters
780     AliDebug(2,Form(" SCM0 par = "));
781     for(int i=0; i<5; i++) AliDebug(2,Form(" %9.5f ", parTRAP[i]));
782     AliDebug(2,Form("\n zpos %f \n",zpos));
783
784     double tanx = (parTRAP[1] -  parTRAP[0]) / (2.*parTRAP[4]); //  tanx =  tany now
785     double tany = (parTRAP[3] -  parTRAP[2]) / (2.*parTRAP[4]), ztmp=0.;
786     parPB[4] = g->GetECPbRadThick()/2.;
787     parSC[2] = g->GetECScintThick()/2.;
788     for(int iz=0; iz<g->GetNECLayers(); iz++){
789       ztmp     = fSampleWidth*double(iz);
790       parPB[0] = parTRAP[0] + tanx*ztmp;
791       parPB[1] = parPB[0]   + tanx*g->GetECPbRadThick();
792       parPB[2] = parTRAP[2] + tany*ztmp;
793       parPB[3] = parPB[2]   + tany*g->GetECPbRadThick();
794       gMC->Gsposp("PBTI", ++nr, name.Data(), xpos, ypos, zpos, 0, pbtiChonly.Data(), parPB, 5) ;
795       AliDebug(2,Form("\n PBTI %2i | zpos %6.3f | par = ", nr, zpos));
796       /*
797       for(int i=0; i<5; i++) printf(" %9.5f ", parPB[i]);
798       // individual SC tile
799       parSC[0] = parPB[0];
800       parSC[1] = parPB[1];
801       gMC->Gsposp("SCTI", nr, name.Data(), xpos, ypos, zpos+g->GetECScintThick(), 
802       0, pbtiChonly.Data(), parSC, 3) ;
803       printf("\n SCTI     zpos %6.3f | par = ", zpos+g->GetECScintThick());
804       for(int i=0; i<3; i++) printf(" %9.5f ", parPB[i]);
805       */
806       zpos  += fSampleWidth;
807     }
808     AliDebug(2,Form("\n"));
809   }
810   Info("PbInTrapForTrd2", "Ver. 0.03 : was tested.");
811 }
812
813 // 15-mar-05
814 //______________________________________________________________________
815 void AliEMCALv0::PbmoInTrd2(const AliEMCALGeometry * g, const Double_t emodPar[5], Double_t parPBMO[5])
816 {
817   // Pb inside Trd2
818   Info("PbmoInTrd2"," started : geometry %s ", g->GetName());
819   double wallThickness = g->GetPhiModuleSize()/2. -  g->GetPhiTileSize();
820   AliDebug(2,Form(" wall thickness %7.5f \n", wallThickness));
821   for(int i=0; i<4; i++) {
822     parPBMO[i] = emodPar[i] - wallThickness;
823     AliDebug(2,Form(" %i parPBMO %7.3f emodPar %7.3f : dif %7.3f \n", 
824                     i, parPBMO[i],emodPar[i], parPBMO[i]-emodPar[i]));
825   }
826   parPBMO[4] = emodPar[4];
827   gMC->Gsvolu("PBMO", "TRD2", fIdTmedArr[kIdPB], parPBMO, 5);
828   gMC->Gspos("PBMO", 1, "EMOD", 0., 0., 0., 0, "ONLY") ;
829   // Division 
830   if(g->GetNPHIdiv()==2 && g->GetNETAdiv()==2) {
831     Division2X2InPbmo(g, parPBMO);
832     AliDebug(2,Form(" PBMO, division 2X2 | geometry |%s|\n", g->GetName()));
833   } else {
834     AliDebug(2,Form(" no division PBMO in this geometry |%s|\n", g->GetName()));
835     assert(0);
836   }
837 }
838
839 //______________________________________________________________________
840 void AliEMCALv0::Division2X2InPbmo(const AliEMCALGeometry * g, const Double_t parPBMO[5]) 
841 {
842   // Division 2X2
843   Info("Division2X2InPbmo"," started : geometry %s ", g->GetName());
844   //Double_t *dummy=0;
845   //  gMC->Gsvolu("SCTI", "BOX", fIdTmedArr[kIdSC], dummy, 0);
846
847   double parSC[3];
848   //double fSampleWidth = double(g->GetECPbRadThick()+g->GetECScintThick());
849   double xpos = 0.0, ypos = 0.0, zpos = 0.0, ztmp=0;;
850   double tanx = (parPBMO[1] -  parPBMO[0]) / (2.*parPBMO[4]); //  tanx =  tany now
851   double tany = (parPBMO[3] -  parPBMO[2]) / (2.*parPBMO[4]);
852   const Int_t buffersize = 10;
853   char name[buffersize], named[buffersize], named2[buffersize];
854
855   AliDebug(2,Form(" PBMO par = "));
856   for(int i=0; i<5; i++) AliDebug(2,Form(" %9.5f ", parPBMO[i]));
857   AliDebug(2,Form("\n"));
858
859   parSC[2] = g->GetECScintThick()/2.;
860   zpos = -fSampleWidth*g->GetNECLayers()/2. + g->GetECPbRadThick() + g->GetECScintThick()/2.;
861   AliDebug(2,Form(" parSC[2] %9.5f \n", parSC[2]));
862   for(int iz=0; iz<g->GetNECLayers(); iz++){
863     ztmp     = g->GetECPbRadThick() + fSampleWidth*double(iz); // Z for previous PB
864     parSC[0] =  parPBMO[0] + tanx*ztmp;
865     parSC[1] =  parPBMO[2] + tany*ztmp;
866
867     snprintf(name,buffersize,"SC%2.2i", iz+1);
868     gMC->Gsvolu(name, "BOX", fIdTmedArr[kIdSC], parSC, 3);
869     gMC->Gspos(name, 1, "PBMO", xpos, ypos, zpos, 0, "ONLY") ;
870     AliDebug(2,Form("%s | zpos %6.3f | parSC[0,1]=(%7.5f,%7.5f) -> ", 
871                     name, zpos, parSC[0], parSC[1]));
872     
873     snprintf(named,buffersize,"SY%2.2i", iz+1);
874     printf(" %s -> ", named);
875     gMC->Gsdvn(named,name, 2, 2);
876
877     snprintf(named2,buffersize,"SX%2.2i", iz+1);
878     printf(" %s \n", named2);
879     gMC->Gsdvn(named2,named, 2, 1);
880
881     zpos    += fSampleWidth;
882   }
883 }
884
885 //______________________________________________________________________
886 AliEMCALShishKebabTrd1Module* AliEMCALv0::GetShishKebabModule(Int_t neta)
887
888   // 28-oct-05
889   AliEMCALShishKebabTrd1Module* trd1=0;
890   if(fShishKebabModules && neta>=0 && neta<fShishKebabModules->GetSize()) {
891     trd1 = (AliEMCALShishKebabTrd1Module*)fShishKebabModules->At(neta);
892   }
893   return trd1;
894 }
895
896 //_____________________________________________________________________________
897 void AliEMCALv0::AddAlignableVolumes() const
898 {
899   //Add volumes which are alignable (?)
900   TString ntmp(GetTitle()); // name of EMCAL geometry
901
902   if(ntmp.Contains("WSUC")) {
903     AddAlignableVolumesInWSUC();  // WSUC case
904   } else {
905     AddAlignableVolumesInALICE(); // ALICE case
906   }
907 }
908
909 //______________________________________________________________________
910 void AliEMCALv0::AddAlignableVolumesInALICE() const
911 {
912   //
913   // Create entries for alignable volumes associating the symbolic volume
914   // name with the corresponding volume path. Needs to be synchronized with
915   // eventual changes in the geometry.
916   //
917
918   Float_t * pars = GetGeometry()->GetSuperModulesPars();
919   double rpos = (GetGeometry()->GetEnvelop(0) + GetGeometry()->GetEnvelop(1))/2.;
920   double phi, phiRad, xpos, ypos, zpos;
921
922   AliGeomManager::ELayerID idEMCAL = AliGeomManager::kEMCAL;
923   Int_t modUID, modnum = 0;
924   TString volpath, symname;
925
926   Int_t nSMod = GetGeometry()->GetNumberOfSuperModules(); 
927   for (Int_t smodnum=0; smodnum < nSMod; smodnum++) {
928     modUID = AliGeomManager::LayerToVolUID(idEMCAL,modnum++);
929     volpath = "ALIC_1/XEN1_1/SMOD_";
930     volpath += (smodnum+1);
931     symname = "EMCAL/FullSupermodule";
932     symname += (smodnum+1);
933
934     if(GetGeometry()->GetKey110DEG() && smodnum>=10) {
935       volpath = "ALIC_1/XEN1_1/SM10_";
936       volpath += (smodnum-10+1);
937       symname = "EMCAL/HalfSupermodule";    
938       symname += (smodnum-10+1);
939     }
940
941     if(!gGeoManager->SetAlignableEntry(symname.Data(),volpath.Data(),modUID))
942       AliFatal("AliEMCALv0::Unable to set alignable entry!!");
943
944     // Creates the Tracking to Local transformation matrix for EMCAL
945     // modules                         
946     TGeoPNEntry *alignableEntry = gGeoManager->GetAlignableEntryByUID(modUID) ;
947
948     phiRad = GetGeometry()->GetPhiCenterOfSM(smodnum);  //comes in radians, not degrees
949     phi = phiRad*180./TMath::Pi();             //need degrees for rot. matrix
950     xpos = rpos * TMath::Cos(phiRad);
951     ypos = rpos * TMath::Sin(phiRad);
952     zpos = pars[2];
953     if(GetGeometry()->GetKey110DEG() && smodnum >= 10) {
954       xpos += (pars[1]/2. * TMath::Sin(phiRad));
955       ypos -= (pars[1]/2. * TMath::Cos(phiRad));
956     }
957
958     TGeoHMatrix *matTtoL;
959     TGeoHMatrix *globMatrix = alignableEntry->GetGlobalOrig();
960
961     if(smodnum%2 == 0) {
962       // pozitive z                                                        
963       TGeoTranslation geoTran0(xpos,ypos, zpos);                        
964       TGeoRotation geoRot0("geoRot0", 90.0, phi, 90.0, 90.0+phi, 0.0, 0.0);
965       TGeoCombiTrans mat0(geoTran0, geoRot0);
966       matTtoL = new TGeoHMatrix(mat0);
967
968       matTtoL->MultiplyLeft(&(globMatrix->Inverse()));
969       alignableEntry->SetMatrix(matTtoL);
970
971     } else {
972       // negative z                                                                                 
973       double phiy = 90. + phi + 180.;
974       if(phiy>=360.) phiy -= 360.;
975       TGeoTranslation geoTran1(xpos,ypos,-zpos);
976       TGeoRotation geoRot1("geoRot1", 90.0, phi, 90.0, phiy, 180.0, 0.0);
977       TGeoCombiTrans mat1(geoTran1, geoRot1);
978       matTtoL = new TGeoHMatrix(mat1);
979
980       matTtoL->MultiplyLeft(&(globMatrix->Inverse()));
981       alignableEntry->SetMatrix(matTtoL);
982
983     }
984
985   }
986
987 }
988
989 //______________________________________________________________________
990 void AliEMCALv0::AddAlignableVolumesInWSUC() const
991 {
992   //
993   // Create entries for alignable volumes associating the symbolic volume
994   // name with the corresponding volume path. Needs to be synchronized with
995   // eventual changes in the geometry.
996   //
997
998   TString vpstr1 = "WSUC_1/XEN1_1/SMOD_";
999   TString snstr1 = "EMCAL/CosmicTestSupermodule";
1000   TString volpath, symname;
1001
1002   // #SM is just one 
1003   for (Int_t smodnum=0; smodnum < 1; smodnum++) {
1004     symname = snstr1;
1005     symname += (smodnum+1);
1006     volpath = vpstr1;
1007     volpath += (smodnum+1);
1008     if(!gGeoManager->SetAlignableEntry(symname.Data(),volpath.Data()))
1009       AliFatal("Unable to set alignable entry!!");
1010   }
1011 }