Removing obsolete includes of TGeant3.h
[u/mrichter/AliRoot.git] / VZERO / AliVZEROv4.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
16 /* $Id$ */
17
18 //////////////////////////////////////////////////////////////////////
19 //                                                                  //
20 //  (V-zero) detector  version 4  as designed by the Lyon group     //
21 //   All comments should be sent to Brigitte CHEYNIS :              //
22 //                                  b.cheynis@ipnl.in2p3.fr         // 
23 //   Geometry of the   24th of february 2004                        //
24 //  (now 4 rings instead of 3 rings as in previous version)         // 
25 //   V0R (now V0C) sits between Z values  -89.4 and  -84.9 cm       //
26 //   V0L (now V0A) sits between Z values +350.0 and +352.0 cm       //
27 //   New coordinate system has been implemented in october 2003     //
28 //                                                                  //
29 //////////////////////////////////////////////////////////////////////
30
31 // --- Standard libraries ---
32 #include <Riostream.h>
33 #include <stdlib.h>
34 #include <string.h>
35
36 // --- ROOT libraries ---
37 #include <TClonesArray.h>
38 #include <TGeometry.h>
39 #include <TLorentzVector.h>
40 #include <TMath.h>
41 #include <TNode.h>
42 #include <TObjectTable.h>
43 #include <TTUBE.h>
44 #include <TTUBS.h>
45 #include <TVirtualMC.h>
46 #include <TParticle.h>
47
48 // --- AliRoot header files ---
49 #include "AliRun.h"
50 #include "AliMC.h"
51 #include "AliConst.h"
52 #include "AliMagF.h"
53 #include "AliVZEROLoader.h"
54 #include "AliVZEROdigit.h"
55 #include "AliVZEROhit.h"
56 #include "AliVZEROv4.h"
57  
58 ClassImp(AliVZEROv4)
59
60 //_____________________________________________________________________________
61 AliVZEROv4:: AliVZEROv4():AliVZERO()
62 {
63 // Standard default constructor 
64 }
65
66 //_____________________________________________________________________________
67 AliVZEROv4::AliVZEROv4(const char *name, const char *title):
68  AliVZERO(name,title)
69 {
70
71 // Standard constructor for V-zero Detector  version 2
72
73   Int_t i;
74
75   printf("\n");
76   for(i=0;i<26;i++) printf("*");
77   printf(" Create VZERO object ");
78   for(i=0;i<26;i++) printf("*");
79   printf("\n");
80       
81   fLightYield              =  93.75; // Light yield in BC408 (93.75 eV per photon)
82   fLightAttenuation        =   0.05; // Light attenuation in fiber (0.05 per meter)
83   fnMeters                 =   15.0; // Number of meters of clear fibers to PM
84   fFibToPhot               =    0.3; // Attenuation at fiber-photocathode interface
85 }
86
87 //_____________________________________________________________________________
88 void AliVZEROv4::CreateGeometry()
89 {
90
91 // Creates the GEANT geometry of the V-zero Detector  version 3
92   
93   Int_t i;
94   
95   printf("\n");
96   for(i=0;i<26;i++) printf("*");
97   printf(" Create VZERO Geometry ");
98   for(i=0;i<26;i++) printf("*");
99   printf("\n");
100       
101   Int_t    *idtmed = fIdtmed->GetArray()-2999;
102
103   Int_t    ndetR = 1;
104   Int_t    ndetL = 1;
105  
106   Int_t    ncellsR = 1;
107   Int_t    ncellsL = 1;
108   
109   Int_t    idrotm[999];
110  
111   Float_t  height1Right, height2Right, height3Right, height4Right; 
112   Float_t  heightRight;
113   Float_t  theta;  
114   
115   Float_t  halfThickQua;
116   
117   Float_t  zdet;
118   Float_t  r0Right, r4Right;
119   Float_t  pi = TMath::Pi();
120     
121   height1Right    =     3.42;  // height of cell 1, in cm
122   height2Right    =     6.78;  // height of cell 2, in cm
123   height3Right    =     7.11;  // height of cell 3, in cm
124   height4Right    =    10.91;  // height of cell 4, in cm
125   
126   theta       = pi/6.0/2.0;       // half angular opening = 15 degrees
127     
128   halfThickQua= fThickness1/2.0;  // half thickness of elementary cell (inner ring)
129   
130 // distance 0.6 cm in zdet accounts for the fact V0R box back lid sits 0.6 away from 
131 // absorber nose sitting at 90 cm. Will use -zdet later...
132 // size of V0R box  (fThickness) is increased by 1 mm as compared to version v2 
133
134   fThickness  =    fThickness + 0.1;  
135   zdet        =    90.0 - 0.6 - fThickness/2.0;  // distance to vertex (along Z axis)   
136   r0Right     =    4.05;          // closest distance to center of the beam pipe
137   heightRight =    height1Right + height2Right + height3Right + height4Right;
138   r4Right     =    r0Right + heightRight;
139
140 // Creation of mother volume v0LE - left part - :
141 // Entrance face at  +350.0 cm  (new coordinate system) ...
142
143    Float_t   partube[3];
144    
145    partube[0] =  4.3;
146    partube[1] = 45.0;
147    partube[2] = fThickness1/2.0;   
148     
149    gMC->Gsvolu("V0LE","TUBE",idtmed[3005],partube,3);
150      
151 // Creation of five rings - left part - :
152 // Entrance face at +350.0 cm  (new coordinate system) ... 
153
154 // Mother volume v0L0 in which will be set 5 scintillator cells 
155
156   Float_t   partubs[5];  
157     
158   Float_t   r0Left      =   4.3;   
159   Float_t   height1Left =   3.6;  
160   Float_t   height2Left =   6.4;  
161   Float_t   height3Left =  14.9; 
162   Float_t   height4Left =  14.4;  
163   Float_t   heightLeft  =   height1Left + height2Left + height3Left + height4Left; 
164                                       
165   Float_t   r4Left      =   r0Left  + heightLeft; 
166   
167   partubs[0]     =  r0Left;
168   partubs[1]     =  r4Left;
169   partubs[2]     =  fThickness1/2.0;
170   partubs[3]     =  90.0-15.0;
171   partubs[4]     = 120.0-15.0;
172
173   gMC->Gsvolu("V0L0","TUBS",idtmed[3010],partubs,5);  // air volume
174   
175   Float_t  r1Left =  r0Left + height1Left;        
176      
177   partubs[0]     =  r0Left;
178   partubs[1]     =  r1Left;
179
180   gMC->Gsvolu("V0L1","TUBS",idtmed[3005],partubs,5);  // quartz volume
181   gMC->Gspos("V0L1",1,"V0L0", 0.0, 0.0 , 0.0, 0,"ONLY"); 
182
183   Float_t  r2Left  =  r1Left + height2Left;       
184   
185   partubs[0]     =  r1Left;
186   partubs[1]     =  r2Left;
187
188   gMC->Gsvolu("V0L2","TUBS",idtmed[3005],partubs,5);  // quartz volume
189   gMC->Gspos("V0L2",1,"V0L0", 0.0, 0.0 , 0.0, 0,"ONLY"); 
190   
191   Float_t  r3Left  =  r2Left + height3Left;       
192       
193   partubs[0]     =  r2Left;
194   partubs[1]     =  r3Left;
195
196   gMC->Gsvolu("V0L3","TUBS",idtmed[3005],partubs,5);  // quartz volume
197   gMC->Gspos("V0L3",1,"V0L0", 0.0, 0.0 , 0.0, 0,"ONLY");
198
199   partubs[0]     =  r3Left;
200   partubs[1]     =  r4Left;
201
202   gMC->Gsvolu("V0L4","TUBS",idtmed[3005],partubs,5);  // quartz volume
203   gMC->Gspos("V0L4",1,"V0L0", 0.0, 0.0 , 0.0, 0,"ONLY");  
204
205 // Creation of mother volume v0RI - right part - :
206   
207   partube[0] = r0Right - 0.2;
208   partube[1] = r4Right + 1.0;
209   partube[2] = fThickness/2.0; 
210       
211   gMC->Gsvolu("V0RI","TUBE",idtmed[3010],partube,3);
212   
213 // Creation of  carbon lids (3.5 mm thick) to keep v0RI box shut...
214
215   Float_t  lidThickness = 0.35;
216  
217   partube[0] =   r0Right;
218   partube[1] =   r4Right;
219   partube[2] =   +lidThickness/2.0;
220     
221   gMC->Gsvolu("V0CA","TUBE",idtmed[3001],partube,3); 
222   gMC->Gspos("V0CA",1,"V0RI",0.0,0.0, fThickness/2.0-partube[2],0,"ONLY");
223   gMC->Gspos("V0CA",2,"V0RI",0.0,0.0,-fThickness/2.0+partube[2],0,"ONLY");
224   
225 // Creation of aluminum rings to maintain the v0RI pieces ...
226
227   partube[0] =   r0Right - 0.2;
228   partube[1] =   r0Right;
229   partube[2] =   +fThickness/2.0;
230    
231   gMC->Gsvolu("V0IR","TUBE",idtmed[3003],partube,3);    
232   gMC->Gspos("V0IR",1,"V0RI",0.0,0.0,0.0,0,"ONLY");
233
234   partube[0] =   r4Right;
235   partube[1] =   r4Right + 1.0;
236   partube[2] =   +fThickness/2.0;
237  
238   gMC->Gsvolu("V0ER","TUBE",idtmed[3003],partube,3);    
239   gMC->Gspos("V0ER",1,"V0RI",0.0,0.0,0.0,0,"ONLY");
240   
241 // Mother volume v0R0 in which will be set 4  scintillator cells 
242   
243   partubs[0]      =  r0Right;
244   partubs[1]      =  r4Right;
245   partubs[2]      =  fThickness/2.0;
246   partubs[3]      =  90.0-15.0;
247   partubs[4]      = 120.0-15.0;
248
249   gMC->Gsvolu("V0R0","TUBS",idtmed[3010],partubs,5);  // air volume 
250
251 // Elementary cell of ring 1 :
252 // (cells of ring 1  will be shifted by 1.7 cm towards vertex to output fibers) 
253    
254   Float_t   offsetFibers  =  1.7;
255   Float_t   offset        =  fThickness/2.0 - lidThickness - fThickness1/2.0; 
256   Float_t   r1Right       =  r0Right + height1Right;
257       
258   partubs[0]     =  r0Right;
259   partubs[1]     =  r1Right;
260   partubs[2]     =  fThickness1/2.0;
261   
262   gMC->Gsvolu("V0R1","TUBS",idtmed[3005],partubs,5);  // scintillator volume
263   gMC->Gspos("V0R1",1,"V0R0", 0.0, 0.0 , -offset + offsetFibers, 0,"ONLY"); 
264
265 // Elementary cell of ring 2 :
266
267   Float_t   r2Right   =  r1Right + height2Right;       
268   
269   partubs[0]     =  r1Right;
270   partubs[1]     =  r2Right;
271
272   gMC->Gsvolu("V0R2","TUBS",idtmed[3005],partubs,5);  // scintillator volume
273   gMC->Gspos("V0R2",1,"V0R0", 0.0, 0.0 , -offset, 0,"ONLY"); 
274
275
276 // Elementary cell of ring 3 :
277
278   Float_t   r3Right   =  r2Right + height3Right;   
279      
280   partubs[0]     =  r2Right;
281   partubs[1]     =  r3Right;
282
283   gMC->Gsvolu("V0R3","TUBS",idtmed[3005],partubs,5);  // scintillator volume
284   gMC->Gspos("V0R3",1,"V0R0", 0.0, 0.0 , -offset, 0,"ONLY");
285
286 // Elementary cell of ring 4 :
287    
288   partubs[0]     =  r3Right;
289   partubs[1]     =  r4Right;
290
291   gMC->Gsvolu("V0R4","TUBS",idtmed[3005],partubs,5);  // scintillator volume
292   gMC->Gspos("V0R4",1,"V0R0", 0.0, 0.0 , -offset, 0,"ONLY");
293
294   Float_t  phiDeg = 180./6.; 
295
296 // Right part : 
297  
298   for(Float_t  phi = 15.0 ; phi < 360.0; phi = phi + phiDeg)
299       {        
300         AliMatrix(idrotm[902], 90.0, phi, 90.0, 90.0 +phi, 0.0 , 0.0);
301         gMC->Gspos("V0R0",ndetR,"V0RI",0.0,
302                           0.0,0.0,idrotm[902],"ONLY");
303         ndetR++;
304        }
305
306   gMC->Gspos("V0RI",1,"ALIC",0.0,0.0,-zdet,0,"ONLY");
307  
308   ncellsR = (ndetR - 1) * 4;  
309   printf("    Number of cells on Right side =   %d\n",  ncellsR);    
310
311 // Left part :
312   
313   for(Float_t  phi = 15.0 ; phi < 360.0; phi = phi + phiDeg)
314       {       
315         AliMatrix(idrotm[902], 90.0, phi, 90.0, 90.0 +phi, 0.0 , 0.0);
316         gMC->Gspos("V0L0",ndetL,"V0LE",0.0,
317                           0.0,0.0,idrotm[902],"ONLY");
318         ndetL++;
319        }
320
321   gMC->Gspos("V0LE",1,"ALIC",0.0,0.0,350.0+fThickness1/2.0,0,"ONLY");
322  
323   ncellsL = (ndetL - 1) * 4;
324   printf("    Number of cells on Left side  =   %d\n",  ncellsL);    
325   for(i=0;i<75;i++) printf("*");
326   printf("\n");
327            
328 }
329             
330 //_____________________________________________________________________________
331 void AliVZEROv4::BuildGeometry()
332 {
333   
334 // Builds simple ROOT TNode geometry for event display
335
336   Int_t i;
337
338   printf("\n");
339   for(i=0;i<30;i++) printf("*");
340   printf(" VZERO BuildGeometry ");
341   for(i=0;i<30;i++) printf("*");
342   printf("\n");
343
344   TNode *top; 
345
346   TNode *v0Rnode, *v0Rnode0, *v0Rnode6 , *v0Rnode7, *v0Rnode8, *v0Rnode9;
347   TNode *v0Rnode1, *v0Rnode2, *v0Rnode3, *v0Rnode4;
348   TNode *v0Lnode, *v0Lnode0;
349   TNode *v0Lnode1, *v0Lnode2, *v0Lnode3, *v0Lnode4;
350    
351   const int kColorVZERO  = kGreen;
352  
353   top = gAlice->GetGeometry()->GetNode("alice");
354
355   Float_t  height1Right, height2Right, height3Right, height4Right; 
356   Float_t  heightRight;
357   Float_t  theta;  
358   
359   Float_t  halfThickQua;
360   Float_t  zdet;
361   Float_t  r0Right, r4Right;
362   Float_t  pi = TMath::Pi();
363
364   height1Right      =     3.42;  // height of cell 1, in cm
365   height2Right      =     6.78;  // height of cell 2, in cm
366   height3Right      =     7.11;  // height of cell 3, in cm
367   height4Right      =    10.91;  // height of cell 4, in cm
368     
369   theta             =    pi/6.0/2.0;    
370    
371   halfThickQua      =    fThickness1/2.0; 
372   
373   zdet              =    90.0 - 0.6 - fThickness/2.0;   
374   r0Right           =    4.05;         
375   heightRight       =    height1Right + height2Right + height3Right + height4Right;
376   r4Right           =    r0Right + heightRight;
377   
378   Int_t     ndiv    =     1;  
379
380   Float_t   partube[3];
381
382   partube[0] =  r0Right - 0.2;
383   partube[1] =  r4Right + 1.0;
384   partube[2] = fThickness/2.0;   
385   
386   TTUBE *v0RI = new TTUBE("V0RI", "V0RI", "void", partube[0], partube[1], partube[2]);
387                 
388   top->cd();
389   
390   v0Rnode = new TNode("V0RI","V0RI",v0RI,0.0,0.0,-zdet,0);
391   
392   v0Rnode->SetLineColor(kYellow);
393   fNodes->Add(v0Rnode);  
394   v0Rnode->SetVisibility(2);     
395  
396 // Rondelles de carbone (epaisseur 3.5 mm) de maintien des cellules ...
397
398   Float_t  lidThickness = 0.35;
399   
400   partube[0] =   r0Right;
401   partube[1] =   r4Right;
402   partube[2] =   +lidThickness/2.0;
403   
404   TTUBE  *v0CA = new TTUBE("V0CA", "V0CA", "void",partube[0], partube[1], partube[2]);
405   
406   v0Rnode->cd();
407   v0Rnode6 = new TNode("V0CA", "V0CA",v0CA,0.0,0.0, fThickness/2.0-partube[2],0);        
408   v0Rnode6->SetLineColor(kYellow);
409   fNodes->Add(v0Rnode6); 
410   v0Rnode->cd();
411   v0Rnode7 = new TNode("V0CA", "V0CA",v0CA,0.0,0.0,-fThickness/2.0+partube[2],0);        
412   v0Rnode7->SetLineColor(kYellow);
413   fNodes->Add(v0Rnode7);
414   
415   partube[0] =   r0Right - 0.2;
416   partube[1] =   r0Right;
417   partube[2] =   +fThickness/2.0;
418   
419   TTUBE *v0IR = new TTUBE("V0IR","V0IR","void", partube[0], partube[1], partube[2]);
420  
421   v0Rnode->cd();
422   v0Rnode8 = new TNode("V0IR", "V0IR",v0IR,0.0,0.0,0.0,0);
423   v0Rnode8->SetLineColor(kYellow);
424   fNodes->Add(v0Rnode8);
425
426   partube[0] =   r4Right;
427   partube[1] =   r4Right + 1.0; 
428   partube[2] =   +fThickness/2.0;
429
430   TTUBE  *v0ER = new TTUBE("V0ER","V0ER","void", partube[0], partube[1], partube[2]);
431   
432   v0Rnode->cd();
433   v0Rnode9 = new TNode("V0ER", "V0ER",v0ER,0.0,0.0,0.0,0);
434   v0Rnode9->SetLineColor(kYellow);
435   fNodes->Add(v0Rnode9);
436   
437   Float_t   partubs[5];
438  
439   partubs[0]      =  r0Right;
440   partubs[1]      =  r4Right;
441   partubs[2]      =  fThickness/2.0;
442   partubs[3]      =  90.0-15.0;
443   partubs[4]      = 120.0-15.0;
444
445   TTUBS  *v0R0 = new TTUBS("V0R0", "V0R0", "void",partubs[0], partubs[1], partubs[2], 
446                                                   partubs[3], partubs[4]); 
447                                                   
448   v0R0->SetNumberOfDivisions(ndiv);                                               
449
450   Float_t   r1Right =  r0Right + height1Right;
451   Float_t   offset  = fThickness/2.0 - lidThickness - fThickness1/2.0; 
452   Float_t   offsetFibers = 1.7;
453     
454   partubs[0]     =  r0Right;
455   partubs[1]     =  r1Right;
456   partubs[2]     =  fThickness1/2.0;
457
458   TTUBS *v0R1 = new TTUBS("V0R1", "V0R1", "void", partubs[0], partubs[1], partubs[2], 
459                                                   partubs[3], partubs[4]);
460
461   v0R1->SetNumberOfDivisions(ndiv);
462   
463   Float_t   r2Right  =  r1Right + height2Right;       
464   
465   partubs[0]     =  r1Right;
466   partubs[1]     =  r2Right;
467
468   TTUBS *v0R2 = new TTUBS("V0R2", "V0R2", "void", partubs[0], partubs[1], partubs[2], 
469                                                   partubs[3], partubs[4]);
470
471   v0R2->SetNumberOfDivisions(ndiv);
472
473   Float_t   r3Right  =  r2Right + height3Right;     
474    
475   partubs[0]     =  r2Right;
476   partubs[1]     =  r3Right;
477   
478   TTUBS *v0R3 = new TTUBS("V0R3", "V0R3", "void", partubs[0], partubs[1], partubs[2], 
479                                                   partubs[3], partubs[4]);
480   v0R3->SetNumberOfDivisions(ndiv);
481
482   partubs[0]     =  r3Right;
483   partubs[1]     =  r4Right;
484   
485   TTUBS *v0R4 = new TTUBS("V0R4", "V0R4", "void", partubs[0], partubs[1], partubs[2], 
486                                                   partubs[3], partubs[4]);
487   v0R4->SetNumberOfDivisions(ndiv);
488         
489   Float_t  phi;
490   Float_t  phiDeg= 180./6.;
491     
492   Int_t    ndetR = 1; 
493
494   char     nameNode[12];  
495  
496   for (phi = 15.0; phi < 360.0; phi = phi + phiDeg)
497   
498   {     
499     TRotMatrix* mat920 = new TRotMatrix("rot920","rot920", 90.0, +phi, 90., 90.+phi, 0.0, 0.0 );        
500      
501     sprintf(nameNode,"SUBDER%d",ndetR);
502     
503     v0Rnode->cd();
504     v0Rnode0 = new TNode(nameNode,nameNode,v0R0,0.0,0.0, 0.0,mat920);    
505     v0Rnode0->SetLineColor(kYellow);
506     fNodes->Add(v0Rnode0);
507     ndetR++;
508     
509     sprintf(nameNode,"SUBDER%d",ndetR);
510     v0Rnode0->cd();    
511     v0Rnode1 = new TNode(nameNode,nameNode,v0R1,0.0,0.0, -offset+ offsetFibers ,0);      
512     v0Rnode1->SetLineColor(kColorVZERO);
513     fNodes->Add(v0Rnode1);
514     ndetR++;
515     
516     sprintf(nameNode,"SUBDER%d",ndetR);
517     v0Rnode0->cd();    
518     v0Rnode2 = new TNode(nameNode,nameNode,v0R2,0.0,0.0, -offset,0);     
519     v0Rnode2->SetLineColor(kColorVZERO);
520     fNodes->Add(v0Rnode2);
521     ndetR++;
522
523     sprintf(nameNode,"SUBDER%d",ndetR);
524     v0Rnode0->cd();    
525     v0Rnode3 = new TNode(nameNode,nameNode,v0R3,0.0,0.0, -offset,0);     
526     v0Rnode3->SetLineColor(kColorVZERO);
527     fNodes->Add(v0Rnode3);
528     ndetR++;
529
530     sprintf(nameNode,"SUBDER%d",ndetR);
531     v0Rnode0->cd();    
532     v0Rnode4 = new TNode(nameNode,nameNode,v0R4,0.0,0.0, -offset,0);     
533     v0Rnode4->SetLineColor(kColorVZERO);
534     fNodes->Add(v0Rnode4);
535     ndetR++;
536        
537     v0Rnode0->SetVisibility(2);    
538   }    
539
540 // Left side of VZERO :
541      
542   Float_t   r0Left      =    4.3;   
543   Float_t   height1Left =    3.6; 
544   Float_t   height2Left =    6.4;
545   Float_t   height3Left =   14.9; 
546   Float_t   height4Left =   14.4;     
547   Float_t   heightLeft  = height1Left + height2Left + height3Left + height4Left; 
548                                         
549   Float_t   r4Left      = r0Left  + heightLeft; 
550
551   partube[0] =  r0Left;
552   partube[1] =  r4Left;
553   partube[2] =  fThickness1/2.0; 
554   
555   TTUBE *v0LE = new TTUBE("V0LE", "V0LE", "void", partube[0], partube[1], partube[2]);
556                 
557   top->cd();
558   
559   v0Lnode = new TNode("V0LE","V0LE",v0LE,0.0,0.0,350.0+fThickness1/2.0,0);
560   
561   v0Lnode->SetLineColor(kBlue);
562   fNodes->Add(v0Lnode);
563   
564   v0Lnode->SetVisibility(2);
565
566   partubs[0]      =  r0Left;
567   partubs[1]      =  r4Left;
568   partubs[2]      =  fThickness1/2.0;
569   partubs[3]      =  90.0-15.0;
570   partubs[4]      = 120.0-15.0;
571   
572   TTUBS *v0L0 = new TTUBS("V0L0", "V0L0", "void", partubs[0], partubs[1], partubs[2], 
573                                                   partubs[3], partubs[4]);
574
575   v0L0->SetNumberOfDivisions(ndiv); 
576   v0L0->SetLineColor(7);
577   
578   Float_t   offsetLeft;
579   offsetLeft    = - fThickness1/2.0; 
580
581   Float_t   r1Left =  r0Left + height1Left;        
582       
583   partubs[0]     =  r0Left;
584   partubs[1]     =  r1Left;
585
586   TTUBS *v0L1 = new TTUBS("V0L1", "V0L1", "void", partubs[0], partubs[1], partubs[2], 
587                                                   partubs[3], partubs[4]);
588   v0L1->SetNumberOfDivisions(ndiv);
589   
590   Float_t   r2Left =  r1Left + height2Left;       
591   
592   partubs[0]     =  r1Left;
593   partubs[1]     =  r2Left;
594
595   TTUBS *v0L2 = new TTUBS("V0L2", "V0L2", "void", partubs[0], partubs[1], partubs[2], 
596                                                   partubs[3], partubs[4]);
597   v0L2->SetNumberOfDivisions(ndiv);
598
599   Float_t   r3Left =  r2Left + height3Left; 
600    
601   partubs[0]     =  r2Left;
602   partubs[1]     =  r3Left;
603   
604   TTUBS *v0L3 = new TTUBS("V0L3", "V0L3", "void", partubs[0], partubs[1], partubs[2], 
605                                                   partubs[3], partubs[4]);
606   v0L3->SetNumberOfDivisions(ndiv);
607
608   partubs[0]     =  r3Left;
609   partubs[1]     =  r4Left;
610   
611   TTUBS *v0L4 = new TTUBS("V0L4", "V0L4", "void", partubs[0], partubs[1], partubs[2], 
612                                                   partubs[3], partubs[4]);
613   v0L4->SetNumberOfDivisions(ndiv);
614  
615   Int_t    ndetL   = 1;
616  
617   for (phi = 15.0; phi < 360.0; phi = phi + phiDeg)
618   
619   {
620      
621     TRotMatrix* mat920 = new TRotMatrix("rot920","rot920", 90.0, +phi, 90., 90.+phi, 0.0, 0.0 );        
622      
623     sprintf(nameNode,"SUBDEL%d",ndetL);
624     
625     v0Lnode->cd();
626     v0Lnode0 = new TNode(nameNode,nameNode,v0L0,0.0,0.0, offsetLeft + halfThickQua,mat920);      
627     v0Lnode0->SetLineColor(kColorVZERO);
628     fNodes->Add(v0Lnode0);
629     ndetL++;
630     
631     sprintf(nameNode,"SUBDEL%d",ndetL);
632     v0Lnode0->cd();    
633     v0Lnode1 = new TNode(nameNode,nameNode,v0L1,0.0,0.0, 0.0,0);         
634     v0Lnode1->SetLineColor(kColorVZERO);
635     fNodes->Add(v0Lnode1);
636     ndetL++;
637     
638     sprintf(nameNode,"SUBDEL%d",ndetL);
639     v0Lnode0->cd();    
640     v0Lnode2 = new TNode(nameNode,nameNode,v0L2,0.0,0.0, 0.0,0);         
641     v0Lnode2->SetLineColor(kColorVZERO);
642     fNodes->Add(v0Lnode2);
643     ndetL++;
644
645     sprintf(nameNode,"SUBDEL%d",ndetL);
646     v0Lnode0->cd();    
647     v0Lnode3 = new TNode(nameNode,nameNode,v0L3,0.0,0.0, 0.0,0);         
648     v0Lnode3->SetLineColor(kColorVZERO);
649     fNodes->Add(v0Lnode3);
650     ndetL++;
651
652     sprintf(nameNode,"SUBDEL%d",ndetL);
653     v0Lnode0->cd();    
654     v0Lnode4 = new TNode(nameNode,nameNode,v0L4,0.0,0.0, 0.0,0);         
655     v0Lnode4->SetLineColor(kColorVZERO);
656     fNodes->Add(v0Lnode4);
657     ndetL++; 
658           
659     v0Lnode0->SetVisibility(2);    
660   }         
661 }  
662     
663 //_____________________________________________________________________________
664 void AliVZEROv4::CreateMaterials()
665 {
666
667 // Creates materials used for geometry 
668
669     Int_t i;
670
671     printf("\n");
672     for(i=0;i<25;i++) printf("*");
673     printf(" VZERO create materials ");
674     for(i=0;i<26;i++) printf("*");
675     printf("\n");
676     
677 /*
678     Float_t ppckov[14] = { 5.5e-9, 5.7e-9, 5.9e-9, 6.1e-9, 6.3e-9, 6.5e-9, 6.7e-9, 
679                            6.9e-9, 7.1e-9, 7.3e-9, 7.5e-9, 7.7e-9, 7.9e-9, 8.1e-9 };
680
681            
682     Float_t ppckov_alu[14] = { 5.5e-9, 5.7e-9, 5.9e-9, 6.1e-9, 6.3e-9, 6.5e-9, 6.7e-9, 
683                                6.9e-9, 7.1e-9, 7.3e-9, 7.5e-9, 7.7e-9, 7.9e-9, 8.1e-9 };
684                            
685     Float_t rindex_quarz[14] = { 1.52398,  1.53090, 1.53835, 1.54641, 1.55513, 1.56458, 
686                                  1.57488,  1.58611, 1.59842, 1.61197, 1.62696, 1.64362, 
687                                  1.662295, 1.68337 };
688                                  
689     Float_t absco_quarz[14] = { 105.8,  45.656, 35.665, 28.598, 25.007, 21.04, 17.525, 
690                                 14.177, 9.282, 4.0925, 1.149, 0.3627, 0.1497, 0.05 };   
691                                                                                         
692     Float_t effic_all[14]   = { 1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1. };
693     
694         
695     Float_t rindex_alu[14] = { 1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1. }; 
696     
697     
698     Float_t absco_alu[14]  = { 1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,1e-4,
699                                1e-4,1e-4,1e-4,1e-4 };
700     Float_t effic_alu[14]  = { 1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1.,1. };
701
702 */        
703
704     Int_t *idtmed = fIdtmed->GetArray()-2999;
705     
706     
707 //  Parameters related to Quarz (SiO2) :
708  
709     Float_t aqua[2], zqua[2], densqua, wmatqua[2];
710     Int_t nlmatqua;
711     
712     aqua[0]    = 28.09;
713     aqua[1]    = 16.;
714     zqua[0]    = 14.;
715     zqua[1]    = 8.;
716     densqua    = 2.64;
717     nlmatqua   = -2;
718     wmatqua[0] = 1.;
719     wmatqua[1] = 2.;
720
721 // Parameters  related to aluminum sheets :
722     
723     Float_t  aal   = 26.98;
724     Float_t  zal   = 13.00; 
725     Float_t  densal=   2.7; 
726     Float_t  radlal=   8.9;
727        
728 // Parameters  related to scintillator CH :
729     
730     Float_t ascin[2] = {1.00794,12.011};
731     Float_t zscin[2] = {1.,6.};
732     Float_t wscin[2] = {1.,1.};
733     Float_t denscin  = 1.032;
734     
735 // AIR
736
737     Float_t aAir[4]={12.,14.,16.,36.};
738     Float_t zAir[4]={6.,7.,8.,18.};
739     Float_t wAir[4]={0.000124,0.755267,0.231781,0.012827};
740     Float_t dAir = 1.20479E-3;
741     
742 //  Definition of materials :
743
744 // AIR
745
746     AliMixture( 1, "AIR A$", aAir,zAir,dAir,4,wAir);
747     AliMixture(11, "AIR I$", aAir,zAir,dAir,4,wAir);
748     AliMaterial( 2, "CARBON$"  , 12.01, 6.0, 2.265, 18.8, 49.9, 0, 0);
749     AliMixture(  3, "QUA", aqua, zqua, densqua, nlmatqua, wmatqua);
750     AliMaterial( 4, "ALUMINIUM1$", 26.98, 13., 2.7, 8.9, 37.2, 0, 0);
751     AliMaterial( 5, "ALUMINIUM2$", aal, zal, densal, radlal, 0, 0, 0);
752  
753     AliMixture( 6, "Scintillator$",ascin,zscin,denscin,-2,wscin);
754     
755      
756     Int_t   iSXFLD = gAlice->Field()->Integ();
757     Float_t sXMGMX = gAlice->Field()->Max();
758     
759     Float_t tmaxfd, stemax, deemax, epsil, stmin;
760         
761     tmaxfd = 10.;
762     stemax = 0.1;
763     deemax = 0.1;     
764     epsil  = 0.001;
765     stmin  = 0.001;
766   
767 //  Active Air :    
768     AliMedium(1, "ACTIVE AIR$", 1, 1, iSXFLD, sXMGMX,
769               10.0, 1.0, 0.1, 0.1, 10.0, 0, 0) ;
770
771 //  Inactive air : 
772   
773     AliMedium(11, "INACTIVE AIR$", 11, 0, iSXFLD, sXMGMX,
774               10.0, 1.0, 0.1, 0.1, 10.0, 0, 0) ;
775     
776     AliMedium(2, "CARBON$ ", 2,  1, iSXFLD, sXMGMX,
777               tmaxfd, stemax, deemax, epsil, stmin, 0, 0);   
778
779     AliMedium(3, "QUARZ$", 3, 1, iSXFLD, sXMGMX,
780               tmaxfd, fMaxStepQua, fMaxDestepQua, epsil, stmin, 0, 0);
781     
782     AliMedium(4,"ALUMINUM1$",4, 1, iSXFLD, sXMGMX,
783               tmaxfd, fMaxStepAlu, fMaxDestepAlu, epsil, stmin, 0, 0);
784               
785
786     AliMedium(5,"ALUMINUM2$",5, 1, iSXFLD, sXMGMX,
787               tmaxfd, fMaxStepAlu, fMaxDestepAlu, epsil, stmin, 0, 0);    
788
789     AliMedium(6,"SCINTILLATOR$",6, 1, iSXFLD, sXMGMX, 10.0, 0.1, 0.1, 0.003, 0.003, 0, 0);
790
791     gMC->Gstpar(idtmed[3000], "LOSS", 1.);  //  [3000] = air ACTIF  [3010] = air INACTIF
792     gMC->Gstpar(idtmed[3000], "HADR", 1.);
793     gMC->Gstpar(idtmed[3000], "DCAY", 1.);
794     gMC->Gstpar(idtmed[3000], "DRAY", 1.);
795     
796     gMC->Gstpar(idtmed[3001], "LOSS", 1.);  //  [3001] = carbon
797     gMC->Gstpar(idtmed[3001], "HADR", 1.);
798     gMC->Gstpar(idtmed[3001], "DCAY", 1.);
799     gMC->Gstpar(idtmed[3001], "DRAY", 1.);
800
801     gMC->Gstpar(idtmed[3002], "LOSS", 1.);  //  [3002] = quartz
802     gMC->Gstpar(idtmed[3002], "HADR", 1.);
803     gMC->Gstpar(idtmed[3002], "DCAY", 1.);
804     gMC->Gstpar(idtmed[3002], "DRAY", 1.);  
805     gMC->Gstpar(idtmed[3002], "CUTGAM",0.5E-4) ; 
806     gMC->Gstpar(idtmed[3002], "CUTELE",1.0E-4) ;
807     
808     gMC->Gstpar(idtmed[3003], "LOSS", 1.);  //  [3003] = normal aluminum
809     gMC->Gstpar(idtmed[3003], "HADR", 1.);
810     gMC->Gstpar(idtmed[3003], "DCAY", 1.);
811     gMC->Gstpar(idtmed[3003], "DRAY", 1.);
812     
813     gMC->Gstpar(idtmed[3004], "LOSS", 1.);  //  [3004] = reflecting aluminum
814     gMC->Gstpar(idtmed[3004], "HADR", 1.);
815     gMC->Gstpar(idtmed[3004], "DCAY", 1.);
816     gMC->Gstpar(idtmed[3004], "DRAY", 1.);
817     gMC->Gstpar(idtmed[3004], "CUTGAM",0.5E-4) ; 
818     gMC->Gstpar(idtmed[3004], "CUTELE",1.0E-4) ;
819     
820     gMC->Gstpar(idtmed[3005], "LOSS", 1.);  //  [3005] = scintillator
821     gMC->Gstpar(idtmed[3005], "HADR", 1.);
822     gMC->Gstpar(idtmed[3005], "DCAY", 1.);
823     gMC->Gstpar(idtmed[3005], "DRAY", 1.); 
824     gMC->Gstpar(idtmed[3005], "CUTGAM",0.5E-4) ; 
825     gMC->Gstpar(idtmed[3005], "CUTELE",1.0E-4) ;
826       
827     
828 //    geant3->Gsckov(idtmed[3002], 14, ppckov, absco_quarz, effic_all,rindex_quarz);    
829 //    geant3->Gsckov(idtmed[3004], 14, ppckov_alu, absco_alu, effic_alu, rindex_alu);
830
831 //    gMC->SetCerenkov(idtmed[3002], 14, ppckov, absco_quarz, effic_all,rindex_quarz);    
832 //    gMC->SetCerenkov(idtmed[3004], 14, ppckov_alu, absco_alu, effic_alu, rindex_alu);
833                                        
834 }
835
836 //_____________________________________________________________________________
837 void AliVZEROv4::DrawModule()
838 {
839
840 //  Drawing is done in DrawVZERO.C
841
842    Int_t i;
843
844    printf("\n");
845    for(i=0;i<30;i++) printf("*");
846    printf(" VZERO DrawModule ");
847    for(i=0;i<30;i++) printf("*");
848    printf("\n");
849 }
850
851 //_____________________________________________________________________________
852 void AliVZEROv4::Init()
853 {
854 // Initialises version 2 of the VZERO Detector
855 // Just prints an information message
856   
857    printf(" VZERO version %d initialized \n",IsVersion());
858    
859 //   gMC->SetMaxStep(fMaxStepAlu);
860 //   gMC->SetMaxStep(fMaxStepQua);
861    
862     AliVZERO::Init();  
863 }
864
865 //_____________________________________________________________________________
866 void AliVZEROv4::StepManager()
867 {
868  
869 // Step Manager, called at each step 
870  
871      Int_t     copy;
872      static    Int_t   vol[4];
873      static    Float_t hits[21];
874      static    Float_t eloss, tlength;
875      static    Int_t   nPhotonsInStep;
876      static    Int_t   nPhotons; 
877      static    Int_t   numStep;
878      Float_t   ringNumber;
879      Float_t   destep, step;
880      
881      numStep += 1; 
882           
883 //   We keep only charged tracks :
884      
885      if ( !gMC->TrackCharge() || !gMC->IsTrackAlive() ) return; 
886
887      vol[0]    = gMC->CurrentVolOffID(1, vol[1]);
888      vol[2]    = gMC->CurrentVolID(copy);
889      vol[3]    = copy;
890      
891      if      ( gMC->CurrentVolID(copy) == gMC->VolId("V0R1") ||
892                gMC->CurrentVolID(copy) == gMC->VolId("V0L1") )
893                ringNumber = 1.0;
894      else if ( gMC->CurrentVolID(copy) == gMC->VolId("V0R2") ||
895                gMC->CurrentVolID(copy) == gMC->VolId("V0L2") ) 
896                ringNumber = 2.0;  
897      else if ( gMC->CurrentVolID(copy) == gMC->VolId("V0R3") ||
898                gMC->CurrentVolID(copy) == gMC->VolId("V0L3") )
899                ringNumber = 3.0;
900      else if ( gMC->CurrentVolID(copy) == gMC->VolId("V0R4") ||
901                gMC->CurrentVolID(copy) == gMC->VolId("V0L4") )
902                ringNumber = 4.0;               
903      else
904                ringNumber = 0.0;
905
906      if  (  ringNumber > 0.5  ) { 
907      
908         destep    = gMC->Edep();
909         step      = gMC->TrackStep();
910         
911         nPhotonsInStep  = Int_t(destep / (fLightYield *1e-9) ); 
912         nPhotonsInStep  = gRandom->Poisson(nPhotonsInStep);
913         
914         eloss    += destep;
915         tlength  += step;        
916         
917         if  ( gMC->IsTrackEntering()  )  { 
918          
919             nPhotons  =  nPhotonsInStep;       
920             gMC->TrackPosition(fTrackPosition);
921             gMC->TrackMomentum(fTrackMomentum);
922             
923             Float_t pt  = TMath::Sqrt( fTrackMomentum.Px() * fTrackMomentum.Px() +
924                                        fTrackMomentum.Py() * fTrackMomentum.Py() );
925                
926             hits[0]  = fTrackPosition.X();
927             hits[1]  = fTrackPosition.Y();
928             hits[2]  = fTrackPosition.Z();               
929             hits[3]  = Float_t (gMC->TrackPid()); 
930
931             hits[4]  = gMC->TrackTime();
932             hits[5]  = gMC->TrackCharge();
933             hits[6]  = fTrackMomentum.Theta()*TMath::RadToDeg();
934             hits[7]  = fTrackMomentum.Phi()*TMath::RadToDeg();
935             hits[8]  = ringNumber;
936          
937             hits[9]  = pt;
938             hits[10] = fTrackMomentum.P();
939             hits[11] = fTrackMomentum.Px();
940             hits[12] = fTrackMomentum.Py();
941             hits[13] = fTrackMomentum.Pz();
942             
943             TParticle *par = gAlice->GetMCApp()->Particle(gAlice->GetMCApp()->GetCurrentTrackNumber());
944             hits[14] = par->Vx();
945             hits[15] = par->Vy();
946             hits[16] = par->Vz();
947             
948             tlength  = 0.0;
949             eloss    = 0.0;         
950          }
951          
952          nPhotons  = nPhotons + nPhotonsInStep;
953          
954          if( gMC->IsTrackExiting() || gMC->IsTrackStop() || gMC->IsTrackDisappeared()){
955          
956          nPhotons  = nPhotons - Int_t((Float_t(nPhotons) * fLightAttenuation * fnMeters));       
957          nPhotons  = nPhotons - Int_t( Float_t(nPhotons) * fFibToPhot);  
958          
959          hits[17] =   eloss;
960          hits[18] = tlength;
961          hits[19] = nPhotons;
962          hits[20] = GetCellId (vol, hits); 
963                  
964          AddHit(gAlice->GetMCApp()->GetCurrentTrackNumber(), vol, hits);
965                  
966          tlength         = 0.0;
967          eloss           = 0.0; 
968          nPhotons        =   0;
969          nPhotonsInStep  =   0;
970          
971          numStep         =   0;  
972          } 
973     }
974       
975 }
976
977 //_____________________________________________________________________________
978 void AliVZEROv4::AddHit(Int_t track, Int_t *vol, Float_t *hits)
979 {
980   
981 //  Adds a VZERO hit
982   
983   TClonesArray &lhits = *fHits;
984   new(lhits[fNhits++]) AliVZEROhit(fIshunt,track,vol,hits);
985 }
986
987 //_____________________________________________________________________________
988 void AliVZEROv4::AddDigits(Int_t *tracks, Int_t* digits) 
989 {
990
991 //  Adds a VZERO digit
992
993    TClonesArray  &ldigits = *fDigits;
994    new(ldigits[fNdigits++]) AliVZEROdigit(tracks, digits);
995 }
996
997 //_____________________________________________________________________________
998 void AliVZEROv4::MakeBranch(Option_t *option)
999 {
1000   
1001 // Creates new branches in the current Root Tree
1002     
1003   char branchname[10];
1004   sprintf(branchname,"%s",GetName());
1005   printf(" fBufferSize = %d \n",fBufferSize);
1006   
1007   const char *cH = strstr(option,"H");
1008   
1009   if (fHits   && TreeH() && cH) {
1010     TreeH()->Branch(branchname,&fHits, fBufferSize);
1011     printf("* AliDetector::MakeBranch * Making Branch %s for hits\n",branchname);
1012   }     
1013
1014   const char *cD = strstr(option,"D");
1015   
1016   if (fDigits   && fLoader->TreeD() && cD) {
1017     fLoader->TreeD()->Branch(branchname,&fDigits, fBufferSize);
1018     printf("* AliDetector::MakeBranch * Making Branch %s for digits\n",branchname);
1019   }  
1020    
1021 }
1022
1023 //_____________________________________________________________________________
1024 Int_t AliVZEROv4::GetCellId(Int_t *vol, Float_t *hits) 
1025 {
1026
1027   //   Returns Id of scintillator cell
1028   //   Right side from  0 to 47
1029   //   Left  side from 48 to 95
1030
1031    Int_t index = vol[1];
1032    fCellId     = 0;
1033    
1034    if (index < 10) index = index + 12; 
1035    
1036    if (hits[2] < 0.0) { 
1037       index = (index - 10) + ( ( Int_t(hits[8]) - 1 ) * 12);
1038       fCellId   = index;
1039    }
1040    else if (hits[2] > 0.0)
1041    {
1042       index = (index + 38) + ( ( Int_t(hits[8]) - 1 ) * 12);
1043       fCellId   = index;}
1044            
1045    return fCellId;
1046 }