Adding a SetTreeAddress in AliPHOS and AliPHOSv0 classes for the re-analysis of the...
[u/mrichter/AliRoot.git] / PHOS / AliPHOSv0.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 // Implementation version v0 of PHOS Manager class 
20 // Layout EMC + PPSD has name GPS2  
21 //                  
22 //*-- Author: Yves Schutz (SUBATECH)
23
24
25 // --- ROOT system ---
26
27 #include "TBRIK.h"
28 #include "TNode.h"
29 #include "TRandom.h"
30
31
32 // --- Standard library ---
33
34 #include <stdio.h>
35 #include <string.h>
36 #include <stdlib.h>
37 #include <strstream.h>
38
39 // --- AliRoot header files ---
40
41 #include "AliPHOSv0.h"
42 #include "AliPHOSHit.h"
43 #include "AliPHOSDigit.h"
44 #include "AliPHOSReconstructioner.h"
45 #include "AliRun.h"
46 #include "AliConst.h"
47
48 ClassImp(AliPHOSv0)
49
50 //____________________________________________________________________________
51 AliPHOSv0::AliPHOSv0()
52 {
53   // ctor
54   fNTmpHits = 0 ; 
55   fTmpHits  = 0 ; 
56 }
57
58 //____________________________________________________________________________
59 AliPHOSv0::AliPHOSv0(const char *name, const char *title):
60   AliPHOS(name,title)
61 {
62   // ctor : title is used to identify the layout
63   //        GPS2 = 5 modules (EMC + PPSD)   
64   // We use 2 arrays of hits :
65   //
66   //   - fHits (the "normal" one), which retains the hits associated with
67   //     the current primary particle being tracked
68   //     (this array is reset after each primary has been tracked).
69   //
70   //   - fTmpHits, which retains all the hits of the current event. It 
71   //     is used for the digitization part.
72
73   fPinElectronicNoise = 0.010 ;
74   fDigitThreshold      = 1. ;   // 1 GeV 
75
76   // We do not want to save in TreeH the raw hits
77   // fHits   = new TClonesArray("AliPHOSHit",100) ;
78   // gAlice->AddHitList(fHits) ; 
79
80   // But save the cumulated hits instead (need to create the branch myself)
81   // It is put in the Digit Tree because the TreeH is filled after each primary
82  // and the TreeD at the end of the event (branch is set in FinishEvent() ).
83   
84   fTmpHits= new TClonesArray("AliPHOSHit",1000) ;
85
86   fNTmpHits = fNhits = 0 ;
87
88   fDigits = new TClonesArray("AliPHOSDigit",1000) ;
89
90
91   fIshunt     =  1 ; // All hits are associated with primary particles
92  
93   // gets an instance of the geometry parameters class  
94    
95   fGeom =  AliPHOSGeometry::GetInstance(title, "") ; 
96
97   if (fGeom->IsInitialized() ) 
98     cout << "AliPHOSv0 : PHOS geometry intialized for " << fGeom->GetName() << endl ;
99   else
100    cout << "AliPHOSv0 : PHOS geometry initialization failed !" << endl ;   
101
102   // For reconstruction
103   // fEmcRecPoints    = new  TObjArray(1000);
104   //fPpsdRecPoints   = new  TObjArray(1000);
105   //fRecParticles    = new  TClonesArray("AliPHOSRecParticle", 1000) ;
106   //fTrackSegments   = new  TClonesArray("AliPHOSTrackSegment", 1000) ;
107 }
108
109 //____________________________________________________________________________
110 AliPHOSv0::AliPHOSv0(AliPHOSReconstructioner * Reconstructioner, const char *name, const char *title):
111   AliPHOS(name,title)
112 {
113   // ctor : title is used to identify the layout
114   //        GPS2 = 5 modules (EMC + PPSD)   
115   // We use 2 arrays of hits :
116   //
117   //   - fHits (the "normal" one), which retains the hits associated with
118   //     the current primary particle being tracked
119   //     (this array is reset after each primary has been tracked).
120   //
121   //   - fTmpHits, which retains all the hits of the current event. It 
122   //     is used for the digitization part.
123
124   fPinElectronicNoise = 0.010 ;
125
126   // We do not want to save in TreeH the raw hits
127   //fHits   = new TClonesArray("AliPHOSHit",100) ;
128
129   fDigits = new TClonesArray("AliPHOSDigit",1000) ;
130   fTmpHits= new TClonesArray("AliPHOSHit",1000) ;
131
132   fNTmpHits = fNhits = 0 ;
133
134   fIshunt     =  1 ; // All hits are associated with primary particles
135  
136   // gets an instance of the geometry parameters class  
137   fGeom =  AliPHOSGeometry::GetInstance(title, "") ; 
138
139   if (fGeom->IsInitialized() ) 
140     cout << "AliPHOSv0 : PHOS geometry intialized for " << fGeom->GetName() << endl ;
141   else
142    cout << "AliPHOSv0 : PHOS geometry initialization failed !" << endl ;   
143
144   // Defining the PHOS Reconstructioner
145  
146  fReconstructioner = Reconstructioner ;
147
148
149 }
150
151 //____________________________________________________________________________
152 AliPHOSv0::~AliPHOSv0()
153 {
154   // dtor
155
156   if ( fTmpHits) {
157     fTmpHits->Delete() ; 
158     delete fTmpHits ;
159     fTmpHits = 0 ; 
160   }
161
162   if ( fEmcRecPoints ) {
163     fEmcRecPoints->Delete() ; 
164     delete fEmcRecPoints ; 
165     fEmcRecPoints = 0 ; 
166   }
167
168   if ( fPpsdRecPoints ) { 
169     fPpsdRecPoints->Delete() ;
170     delete fPpsdRecPoints ;
171     fPpsdRecPoints = 0 ; 
172   }
173   
174   if ( fTrackSegments ) {
175     fTrackSegments->Delete() ; 
176     delete fTrackSegments ;
177     fTrackSegments = 0 ; 
178   }
179
180 }
181
182 //____________________________________________________________________________
183 void AliPHOSv0::AddHit(Int_t primary, Int_t Id, Float_t * hits)
184 {
185   // Add a hit to the hit list.
186   // A PHOS hit is the sum of all hits in a single crystal
187   //   or in a single PPSD gas cell
188
189   Int_t hitCounter ;
190   TClonesArray &ltmphits = *fTmpHits ;
191   AliPHOSHit *newHit ;
192   AliPHOSHit *curHit ;
193   //  AliPHOSHit *curHit2 ;
194   Bool_t deja = kFALSE ;
195
196   // In any case, fills the fTmpHit TClonesArray (with "accumulated hits")
197
198   newHit = new AliPHOSHit(primary, Id, hits) ;
199
200   // We do not want to save in TreeH the raw hits 
201   //  TClonesArray &lhits = *fHits;
202
203   for ( hitCounter = 0 ; hitCounter < fNTmpHits && !deja ; hitCounter++ ) {
204     curHit = (AliPHOSHit*) ltmphits[hitCounter] ;
205   if( *curHit == *newHit ) {
206     *curHit = *curHit + *newHit ;
207     deja = kTRUE ;
208     }
209   }
210          
211   if ( !deja ) {
212     new(ltmphits[fNTmpHits]) AliPHOSHit(*newHit) ;
213     fNTmpHits++ ;
214   }
215
216   // We do not want to save in TreeH the raw hits 
217   //   new(lhits[fNhits]) AliPHOSHit(*newHit) ;    
218   //   fNhits++ ;
219
220   // Please note that the fTmpHits array must survive up to the
221   // end of the events, so it does not appear e.g. in ResetHits() (
222   // which is called at the end of each primary).  
223
224   delete newHit;
225
226 }
227
228
229 //____________________________________________________________________________
230 void AliPHOSv0::BuildGeometry()
231 {
232   // Build the PHOS geometry for the ROOT display
233   //BEGIN_HTML
234   /*
235     <H2>
236      PHOS in ALICE displayed by root
237     </H2>
238     <UL>
239     <LI> All Views
240     <P>
241     <CENTER>
242     <IMG Align=BOTTOM ALT="All Views" SRC="../images/AliPHOSv0AllViews.gif"> 
243     </CENTER></P></LI>
244     <LI> Front View
245     <P>
246     <CENTER>
247     <IMG Align=BOTTOM ALT="Front View" SRC="../images/AliPHOSv0FrontView.gif"> 
248     </CENTER></P></LI>
249      <LI> 3D View 1
250     <P>
251     <CENTER>
252     <IMG Align=BOTTOM ALT="3D View 1" SRC="../images/AliPHOSv03DView1.gif"> 
253     </CENTER></P></LI>
254     <LI> 3D View 2
255     <P>
256     <CENTER>
257     <IMG Align=BOTTOM ALT="3D View 2" SRC="../images/AliPHOSv03DView2.gif"> 
258     </CENTER></P></LI>
259     </UL>
260   */
261   //END_HTML  
262
263   this->BuildGeometryforPHOS() ; 
264   if ( ( strcmp(fGeom->GetName(), "GPS2" ) == 0 ) ) 
265     this->BuildGeometryforPPSD() ;
266   else
267     cout << "AliPHOSv0::BuildGeometry : no charged particle identification system installed" << endl; 
268
269 }
270
271 //____________________________________________________________________________
272 void AliPHOSv0:: BuildGeometryforPHOS(void)
273 {
274  // Build the PHOS-EMC geometry for the ROOT display
275
276   const Int_t kColorPHOS = kRed ;
277   const Int_t kColorXTAL = kBlue ;
278
279   Double_t const kRADDEG = 180.0 / kPI ;
280  
281   new TBRIK( "OuterBox", "PHOS box", "void", fGeom->GetOuterBoxSize(0)/2, 
282                                              fGeom->GetOuterBoxSize(1)/2, 
283                                              fGeom->GetOuterBoxSize(2)/2 );
284
285   // Textolit Wall box, position inside PHOS 
286   
287   new TBRIK( "TextolitBox", "PHOS Textolit box ", "void", fGeom->GetTextolitBoxSize(0)/2, 
288                                                           fGeom->GetTextolitBoxSize(1)/2, 
289                                                           fGeom->GetTextolitBoxSize(2)/2);
290
291   // Polystyrene Foam Plate
292
293   new TBRIK( "UpperFoamPlate", "PHOS Upper foam plate", "void", fGeom->GetTextolitBoxSize(0)/2, 
294                                                                 fGeom->GetSecondUpperPlateThickness()/2, 
295                                                                 fGeom->GetTextolitBoxSize(2)/2 ) ; 
296
297   // Air Filled Box
298  
299   new TBRIK( "AirFilledBox", "PHOS air filled box", "void", fGeom->GetAirFilledBoxSize(0)/2, 
300                                                             fGeom->GetAirFilledBoxSize(1)/2, 
301                                                             fGeom->GetAirFilledBoxSize(2)/2 );
302
303   // Crystals Box
304
305   Float_t xtlX = fGeom->GetCrystalSize(0) ; 
306   Float_t xtlY = fGeom->GetCrystalSize(1) ; 
307   Float_t xtlZ = fGeom->GetCrystalSize(2) ; 
308
309   Float_t xl =  fGeom->GetNPhi() * ( xtlX + 2 * fGeom->GetGapBetweenCrystals() ) / 2.0 + fGeom->GetModuleBoxThickness() ;
310   Float_t yl =  ( xtlY + fGeom->GetCrystalSupportHeight() + fGeom->GetCrystalWrapThickness() + fGeom->GetCrystalHolderThickness() ) / 2.0 
311              + fGeom->GetModuleBoxThickness() / 2.0 ;
312   Float_t zl =  fGeom->GetNZ() * ( xtlZ + 2 * fGeom->GetGapBetweenCrystals() ) / 2.0 +  fGeom->GetModuleBoxThickness() ;
313   
314   new TBRIK( "CrystalsBox", "PHOS crystals box", "void", xl, yl, zl ) ;
315
316 // position PHOS into ALICE
317
318   Float_t r = fGeom->GetIPtoOuterCoverDistance() + fGeom->GetOuterBoxSize(1) / 2.0 ;
319   Int_t number = 988 ; 
320   Float_t pphi =  TMath::ATan( fGeom->GetOuterBoxSize(0)  / ( 2.0 * fGeom->GetIPtoOuterCoverDistance() ) ) ;
321   pphi *= kRADDEG ;
322   TNode * top = gAlice->GetGeometry()->GetNode("alice") ;
323  
324   char * nodename = new char[20] ;  
325   char * rotname  = new char[20] ; 
326
327   for( Int_t i = 1; i <= fGeom->GetNModules(); i++ ) { 
328    Float_t angle = pphi * 2 * ( i - fGeom->GetNModules() / 2.0 - 0.5 ) ;
329    sprintf(rotname, "%s%d", "rot", number++) ;
330    new TRotMatrix(rotname, rotname, 90, angle, 90, 90 + angle, 0, 0);
331    top->cd();
332    sprintf(nodename,"%s%d", "Module", i) ;    
333    Float_t x =  r * TMath::Sin( angle / kRADDEG ) ;
334    Float_t y = -r * TMath::Cos( angle / kRADDEG ) ;
335    TNode * outerboxnode = new TNode(nodename, nodename, "OuterBox", x, y, 0, rotname ) ;
336    outerboxnode->SetLineColor(kColorPHOS) ;
337    fNodes->Add(outerboxnode) ;
338    outerboxnode->cd() ; 
339    // now inside the outer box the textolit box
340    y = ( fGeom->GetOuterBoxThickness(1) -  fGeom->GetUpperPlateThickness() ) / 2.  ;
341    sprintf(nodename,"%s%d", "TexBox", i) ;  
342    TNode * textolitboxnode = new TNode(nodename, nodename, "TextolitBox", 0, y, 0) ; 
343    textolitboxnode->SetLineColor(kColorPHOS) ;
344    fNodes->Add(textolitboxnode) ;
345    // upper foam plate inside outre box
346    outerboxnode->cd() ; 
347    sprintf(nodename, "%s%d", "UFPlate", i) ;
348    y =  ( fGeom->GetTextolitBoxSize(1) - fGeom->GetSecondUpperPlateThickness() ) / 2.0 ;
349    TNode * upperfoamplatenode = new TNode(nodename, nodename, "UpperFoamPlate", 0, y, 0) ; 
350    upperfoamplatenode->SetLineColor(kColorPHOS) ;
351    fNodes->Add(upperfoamplatenode) ;  
352    // air filled box inside textolit box (not drawn)
353    textolitboxnode->cd();
354    y = ( fGeom->GetTextolitBoxSize(1) - fGeom->GetAirFilledBoxSize(1) ) / 2.0 -  fGeom->GetSecondUpperPlateThickness() ;
355    sprintf(nodename, "%s%d", "AFBox", i) ;
356    TNode * airfilledboxnode = new TNode(nodename, nodename, "AirFilledBox", 0, y, 0) ; 
357    fNodes->Add(airfilledboxnode) ; 
358    // crystals box inside air filled box
359    airfilledboxnode->cd() ; 
360    y = fGeom->GetAirFilledBoxSize(1) / 2.0 - yl 
361        - ( fGeom->GetIPtoCrystalSurface() - fGeom->GetIPtoOuterCoverDistance() - fGeom->GetModuleBoxThickness() 
362        -  fGeom->GetUpperPlateThickness() -  fGeom->GetSecondUpperPlateThickness() ) ; 
363    sprintf(nodename, "%s%d", "XTBox", i) ; 
364    TNode * crystalsboxnode = new TNode(nodename, nodename, "CrystalsBox", 0, y, 0) ;    
365    crystalsboxnode->SetLineColor(kColorXTAL) ; 
366    fNodes->Add(crystalsboxnode) ; 
367   }
368 }
369
370 //____________________________________________________________________________
371 void AliPHOSv0:: BuildGeometryforPPSD(void)
372 {
373  //  Build the PHOS-PPSD geometry for the ROOT display
374  //BEGIN_HTML
375   /*
376     <H2>
377      PPSD displayed by root
378     </H2>
379     <UL>
380     <LI> Zoom on PPSD: Front View
381     <P>
382     <CENTER>
383     <IMG Align=BOTTOM ALT="PPSD Front View" SRC="../images/AliPHOSv0PPSDFrontView.gif"> 
384     </CENTER></P></LI>
385     <LI> Zoom on PPSD: Perspective View
386     <P>
387     <CENTER>
388     <IMG Align=BOTTOM ALT="PPSD Prespective View" SRC="../images/AliPHOSv0PPSDPerspectiveView.gif"> 
389     </CENTER></P></LI>
390     </UL>
391   */
392   //END_HTML  
393   Double_t const kRADDEG = 180.0 / kPI ;
394
395   const Int_t kColorPHOS = kRed ;
396   const Int_t kColorPPSD = kGreen ;
397   const Int_t kColorGas  = kBlue ;  
398   const Int_t kColorAir  = kYellow ; 
399
400   // Box for a full PHOS module
401
402   new TBRIK( "PPSDBox", "PPSD box", "void",  fGeom->GetPPSDBoxSize(0)/2, 
403                                              fGeom->GetPPSDBoxSize(1)/2, 
404                                              fGeom->GetPPSDBoxSize(2)/2 );
405
406   // Box containing one micromegas module 
407
408   new TBRIK( "PPSDModule", "PPSD module", "void",  fGeom->GetPPSDModuleSize(0)/2, 
409                                                    fGeom->GetPPSDModuleSize(1)/2, 
410                                                    fGeom->GetPPSDModuleSize(2)/2 );
411  // top lid
412
413   new TBRIK ( "TopLid", "Micromegas top lid", "void",  fGeom->GetPPSDModuleSize(0)/2,
414                                                        fGeom->GetLidThickness()/2,
415                                                        fGeom->GetPPSDModuleSize(2)/2 ) ; 
416  // composite panel (top and bottom)
417
418   new TBRIK ( "TopPanel", "Composite top panel", "void",  ( fGeom->GetPPSDModuleSize(0) - fGeom->GetMicromegasWallThickness() )/2,
419                                                             fGeom->GetCompositeThickness()/2,
420                                                           ( fGeom->GetPPSDModuleSize(2) - fGeom->GetMicromegasWallThickness() )/2 ) ;  
421   
422   new TBRIK ( "BottomPanel", "Composite bottom panel", "void",  ( fGeom->GetPPSDModuleSize(0) - fGeom->GetMicromegasWallThickness() )/2,
423                                                                   fGeom->GetCompositeThickness()/2,
424                                                                 ( fGeom->GetPPSDModuleSize(2) - fGeom->GetMicromegasWallThickness() )/2 ) ; 
425  // gas gap (conversion and avalanche)
426
427   new TBRIK ( "GasGap", "gas gap", "void",  ( fGeom->GetPPSDModuleSize(0) - fGeom->GetMicromegasWallThickness() )/2,
428                                             ( fGeom->GetConversionGap() +  fGeom->GetAvalancheGap() )/2,
429                                             ( fGeom->GetPPSDModuleSize(2) - fGeom->GetMicromegasWallThickness() )/2 ) ; 
430
431  // anode and cathode 
432
433   new TBRIK ( "Anode", "Anode", "void",  ( fGeom->GetPPSDModuleSize(0) - fGeom->GetMicromegasWallThickness() )/2,
434                                            fGeom->GetAnodeThickness()/2,
435                                          ( fGeom->GetPPSDModuleSize(2) - fGeom->GetMicromegasWallThickness() )/2 ) ; 
436
437   new TBRIK ( "Cathode", "Cathode", "void",  ( fGeom->GetPPSDModuleSize(0) - fGeom->GetMicromegasWallThickness() )/2,
438                                                fGeom->GetCathodeThickness()/2,
439                                              ( fGeom->GetPPSDModuleSize(2) - fGeom->GetMicromegasWallThickness() )/2 ) ; 
440  // PC  
441
442   new TBRIK ( "PCBoard", "Printed Circuit", "void",  ( fGeom->GetPPSDModuleSize(0) - fGeom->GetMicromegasWallThickness() )/2,
443                                                        fGeom->GetPCThickness()/2,
444                                                      ( fGeom->GetPPSDModuleSize(2) - fGeom->GetMicromegasWallThickness() )/2 ) ; 
445  // Gap between Lead and top micromegas
446
447   new TBRIK ( "LeadToM", "Air Gap top", "void", fGeom->GetPPSDBoxSize(0)/2,
448                                                 fGeom->GetMicro1ToLeadGap()/2,
449                                                 fGeom->GetPPSDBoxSize(2)/2  ) ;  
450  
451 // Gap between Lead and bottom micromegas
452
453   new TBRIK ( "MToLead", "Air Gap bottom", "void", fGeom->GetPPSDBoxSize(0)/2,
454                                                    fGeom->GetLeadToMicro2Gap()/2,
455                                                    fGeom->GetPPSDBoxSize(2)/2  ) ; 
456  // Lead converter
457    
458   new TBRIK ( "Lead", "Lead converter", "void", fGeom->GetPPSDBoxSize(0)/2,
459                                                 fGeom->GetLeadConverterThickness()/2,
460                                                 fGeom->GetPPSDBoxSize(2)/2  ) ; 
461
462      // position PPSD into ALICE
463
464   char * nodename = new char[20] ;  
465   char * rotname  = new char[20] ; 
466
467   Float_t r = fGeom->GetIPtoTopLidDistance() + fGeom->GetPPSDBoxSize(1) / 2.0 ;
468   Int_t number = 988 ; 
469   TNode * top = gAlice->GetGeometry()->GetNode("alice") ;
470  
471   for( Int_t i = 1; i <= fGeom->GetNModules(); i++ ) { // the number of PHOS modules
472     Float_t angle = fGeom->GetPHOSAngle(i) ;
473     sprintf(rotname, "%s%d", "rotg", number++) ;
474     new TRotMatrix(rotname, rotname, 90, angle, 90, 90 + angle, 0, 0);
475     top->cd();
476     sprintf(nodename, "%s%d", "Moduleg", i) ;    
477     Float_t x =  r * TMath::Sin( angle / kRADDEG ) ;
478     Float_t y = -r * TMath::Cos( angle / kRADDEG ) ;
479     TNode * ppsdboxnode = new TNode(nodename , nodename ,"PPSDBox", x, y, 0, rotname ) ;
480     ppsdboxnode->SetLineColor(kColorPPSD) ;
481     fNodes->Add(ppsdboxnode) ;
482     ppsdboxnode->cd() ;
483     // inside the PPSD box: 
484     //   1.   fNumberOfModulesPhi x fNumberOfModulesZ top micromegas
485     x = ( fGeom->GetPPSDBoxSize(0) - fGeom->GetPPSDModuleSize(0) ) / 2. ;  
486     {
487       for ( Int_t iphi = 1; iphi <= fGeom->GetNumberOfModulesPhi(); iphi++ ) { // the number of micromegas modules in phi per PHOS module
488         Float_t z = ( fGeom->GetPPSDBoxSize(2) - fGeom->GetPPSDModuleSize(2) ) / 2. ;
489         TNode * micro1node ; 
490         for ( Int_t iz = 1; iz <= fGeom->GetNumberOfModulesZ(); iz++ ) { // the number of micromegas modules in z per PHOS module
491           y = ( fGeom->GetPPSDBoxSize(1) - fGeom->GetMicromegas1Thickness() ) / 2. ; 
492           sprintf(nodename, "%s%d%d%d", "Mic1", i, iphi, iz) ;
493           micro1node  = new TNode(nodename, nodename, "PPSDModule", x, y, z) ;
494           micro1node->SetLineColor(kColorPPSD) ;  
495           fNodes->Add(micro1node) ; 
496           // inside top micromegas
497           micro1node->cd() ; 
498           //      a. top lid
499           y = ( fGeom->GetMicromegas1Thickness() - fGeom->GetLidThickness() ) / 2. ; 
500           sprintf(nodename, "%s%d%d%d", "Lid", i, iphi, iz) ;
501           TNode * toplidnode = new TNode(nodename, nodename, "TopLid", 0, y, 0) ;
502           toplidnode->SetLineColor(kColorPPSD) ;  
503           fNodes->Add(toplidnode) ; 
504           //      b. composite panel
505           y = y - fGeom->GetLidThickness() / 2. - fGeom->GetCompositeThickness() / 2. ; 
506           sprintf(nodename, "%s%d%d%d", "CompU", i, iphi, iz) ;
507           TNode * compupnode = new TNode(nodename, nodename, "TopPanel", 0, y, 0) ;
508           compupnode->SetLineColor(kColorPPSD) ;  
509           fNodes->Add(compupnode) ; 
510           //      c. anode
511           y = y - fGeom->GetCompositeThickness() / 2. - fGeom->GetAnodeThickness()  / 2. ; 
512           sprintf(nodename, "%s%d%d%d", "Ano", i, iphi, iz) ;
513           TNode * anodenode = new TNode(nodename, nodename, "Anode", 0, y, 0) ;
514           anodenode->SetLineColor(kColorPHOS) ;  
515           fNodes->Add(anodenode) ; 
516           //      d.  gas 
517           y = y - fGeom->GetAnodeThickness() / 2. - ( fGeom->GetConversionGap() +  fGeom->GetAvalancheGap() ) / 2. ; 
518           sprintf(nodename, "%s%d%d%d", "GGap", i, iphi, iz) ;
519           TNode * ggapnode = new TNode(nodename, nodename, "GasGap", 0, y, 0) ;
520           ggapnode->SetLineColor(kColorGas) ;  
521           fNodes->Add(ggapnode) ;          
522           //      f. cathode
523           y = y - ( fGeom->GetConversionGap() +  fGeom->GetAvalancheGap() ) / 2. - fGeom->GetCathodeThickness()  / 2. ; 
524           sprintf(nodename, "%s%d%d%d", "Cathode", i, iphi, iz) ;
525           TNode * cathodenode = new TNode(nodename, nodename, "Cathode", 0, y, 0) ;
526           cathodenode->SetLineColor(kColorPHOS) ;  
527           fNodes->Add(cathodenode) ;        
528           //      g. printed circuit
529           y = y - fGeom->GetCathodeThickness() / 2. - fGeom->GetPCThickness()  / 2. ; 
530           sprintf(nodename, "%s%d%d%d", "PC", i, iphi, iz) ;
531           TNode * pcnode = new TNode(nodename, nodename, "PCBoard", 0, y, 0) ;
532           pcnode->SetLineColor(kColorPPSD) ;  
533           fNodes->Add(pcnode) ;        
534           //      h. composite panel
535           y = y - fGeom->GetPCThickness() / 2. - fGeom->GetCompositeThickness()  / 2. ; 
536           sprintf(nodename, "%s%d%d%d", "CompDown", i, iphi, iz) ;
537           TNode * compdownnode = new TNode(nodename, nodename, "BottomPanel", 0, y, 0) ;
538           compdownnode->SetLineColor(kColorPPSD) ;  
539           fNodes->Add(compdownnode) ;   
540           z = z - fGeom->GetPPSDModuleSize(2) ;
541           ppsdboxnode->cd() ;
542         } // end of Z module loop     
543         x = x -  fGeom->GetPPSDModuleSize(0) ; 
544         ppsdboxnode->cd() ;
545       } // end of phi module loop
546     }
547     //   2. air gap      
548     ppsdboxnode->cd() ;
549     y = ( fGeom->GetPPSDBoxSize(1) - 2 * fGeom->GetMicromegas1Thickness() - fGeom->GetMicro1ToLeadGap() ) / 2. ; 
550     sprintf(nodename, "%s%d", "GapUp", i) ;
551     TNode * gapupnode = new TNode(nodename, nodename, "LeadToM", 0, y, 0) ;
552     gapupnode->SetLineColor(kColorAir) ;  
553     fNodes->Add(gapupnode) ;        
554     //   3. lead converter
555     y = y - fGeom->GetMicro1ToLeadGap() / 2. - fGeom->GetLeadConverterThickness() / 2. ; 
556     sprintf(nodename, "%s%d", "LeadC", i) ;
557     TNode * leadcnode = new TNode(nodename, nodename, "Lead", 0, y, 0) ;
558     leadcnode->SetLineColor(kColorPPSD) ;  
559     fNodes->Add(leadcnode) ;        
560     //   4. air gap
561     y = y - fGeom->GetLeadConverterThickness() / 2. - fGeom->GetLeadToMicro2Gap()  / 2. ; 
562     sprintf(nodename, "%s%d", "GapDown", i) ;
563     TNode * gapdownnode = new TNode(nodename, nodename, "MToLead", 0, y, 0) ;
564     gapdownnode->SetLineColor(kColorAir) ;  
565     fNodes->Add(gapdownnode) ;        
566     //    5.  fNumberOfModulesPhi x fNumberOfModulesZ bottom micromegas
567     x = ( fGeom->GetPPSDBoxSize(0) - fGeom->GetPPSDModuleSize(0) ) / 2. - fGeom->GetPhiDisplacement() ;  
568     {
569       for ( Int_t iphi = 1; iphi <= fGeom->GetNumberOfModulesPhi(); iphi++ ) { 
570         Float_t z = ( fGeom->GetPPSDBoxSize(2) - fGeom->GetPPSDModuleSize(2) ) / 2.  - fGeom->GetZDisplacement() ;;
571         TNode * micro2node ; 
572         for ( Int_t iz = 1; iz <= fGeom->GetNumberOfModulesZ(); iz++ ) { 
573           y = - ( fGeom->GetPPSDBoxSize(1) - fGeom->GetMicromegas2Thickness() ) / 2. ; 
574           sprintf(nodename, "%s%d%d%d", "Mic2", i, iphi, iz) ;
575           micro2node  = new TNode(nodename, nodename, "PPSDModule", x, y, z) ;
576           micro2node->SetLineColor(kColorPPSD) ;  
577           fNodes->Add(micro2node) ; 
578           // inside bottom micromegas
579           micro2node->cd() ; 
580           //      a. top lid
581           y = ( fGeom->GetMicromegas2Thickness() - fGeom->GetLidThickness() ) / 2. ; 
582           sprintf(nodename, "%s%d", "Lidb", i) ;
583           TNode * toplidbnode = new TNode(nodename, nodename, "TopLid", 0, y, 0) ;
584           toplidbnode->SetLineColor(kColorPPSD) ;  
585           fNodes->Add(toplidbnode) ; 
586           //      b. composite panel
587           y = y - fGeom->GetLidThickness() / 2. - fGeom->GetCompositeThickness() / 2. ; 
588           sprintf(nodename, "%s%d", "CompUb", i) ;
589           TNode * compupbnode = new TNode(nodename, nodename, "TopPanel", 0, y, 0) ;
590           compupbnode->SetLineColor(kColorPPSD) ;  
591           fNodes->Add(compupbnode) ; 
592           //      c. anode
593           y = y - fGeom->GetCompositeThickness() / 2. - fGeom->GetAnodeThickness()  / 2. ; 
594           sprintf(nodename, "%s%d", "Anob", i) ;
595           TNode * anodebnode = new TNode(nodename, nodename, "Anode", 0, y, 0) ;
596           anodebnode->SetLineColor(kColorPPSD) ;  
597           fNodes->Add(anodebnode) ; 
598           //      d. conversion gas
599           y = y - fGeom->GetAnodeThickness() / 2. - ( fGeom->GetConversionGap() +  fGeom->GetAvalancheGap() )  / 2. ; 
600           sprintf(nodename, "%s%d", "GGapb", i) ;
601           TNode * ggapbnode = new TNode(nodename, nodename, "GasGap", 0, y, 0) ;
602           ggapbnode->SetLineColor(kColorGas) ;  
603           fNodes->Add(ggapbnode) ;           
604           //      f. cathode
605           y = y - ( fGeom->GetConversionGap() + fGeom->GetAvalancheGap() ) / 2. - fGeom->GetCathodeThickness()  / 2. ; 
606           sprintf(nodename, "%s%d", "Cathodeb", i) ;
607           TNode * cathodebnode = new TNode(nodename, nodename, "Cathode", 0, y, 0) ;
608           cathodebnode->SetLineColor(kColorPPSD) ;  
609           fNodes->Add(cathodebnode) ;        
610           //      g. printed circuit
611           y = y - fGeom->GetCathodeThickness() / 2. - fGeom->GetPCThickness()  / 2. ; 
612           sprintf(nodename, "%s%d", "PCb", i) ;
613           TNode * pcbnode = new TNode(nodename, nodename, "PCBoard", 0, y, 0) ;
614           pcbnode->SetLineColor(kColorPPSD) ;  
615           fNodes->Add(pcbnode) ;        
616           //      h. composite pane
617           y = y - fGeom->GetPCThickness() / 2. - fGeom->GetCompositeThickness()  / 2. ; 
618           sprintf(nodename, "%s%d", "CompDownb", i) ;
619           TNode * compdownbnode = new TNode(nodename, nodename, "BottomPanel", 0, y, 0) ;
620           compdownbnode->SetLineColor(kColorPPSD) ;  
621           fNodes->Add(compdownbnode) ;        
622           z = z - fGeom->GetPPSDModuleSize(2) ;
623           ppsdboxnode->cd() ;
624         } // end of Z module loop     
625         x = x -  fGeom->GetPPSDModuleSize(0) ; 
626         ppsdboxnode->cd() ;
627       } // end of phi module loop
628     }
629   } // PHOS modules
630  
631   delete[] rotname ;  
632   delete[] nodename ; 
633
634 }
635
636 //____________________________________________________________________________
637 void AliPHOSv0::CreateGeometry()
638 {
639   // Create the PHOS geometry for Geant
640
641   AliPHOSv0 *phostmp = (AliPHOSv0*)gAlice->GetModule("PHOS") ;
642
643   if ( phostmp == NULL ) {
644     
645     fprintf(stderr, "PHOS detector not found!\n") ;
646     return;
647     
648   }
649   // Get pointer to the array containing media indeces
650   Int_t *idtmed = fIdtmed->GetArray() - 699 ;
651
652   Float_t bigbox[3] ; 
653   bigbox[0] =   fGeom->GetOuterBoxSize(0) / 2.0 ;
654   bigbox[1] = ( fGeom->GetOuterBoxSize(1) + fGeom->GetPPSDBoxSize(1) ) / 2.0 ;
655   bigbox[2] =   fGeom->GetOuterBoxSize(2) / 2.0 ;
656   
657   gMC->Gsvolu("PHOS", "BOX ", idtmed[798], bigbox, 3) ;
658   
659   this->CreateGeometryforPHOS() ; 
660   if ( strcmp( fGeom->GetName(), "GPS2") == 0  ) 
661     this->CreateGeometryforPPSD() ;
662   else
663     cout << "AliPHOSv0::CreateGeometry : no charged particle identification system installed" << endl; 
664   
665   // --- Position  PHOS mdules in ALICE setup ---
666   
667   Int_t idrotm[99] ;
668   Double_t const kRADDEG = 180.0 / kPI ;
669   
670   for( Int_t i = 1; i <= fGeom->GetNModules(); i++ ) {
671     
672     Float_t angle = fGeom->GetPHOSAngle(i) ;
673     AliMatrix(idrotm[i-1], 90.0, angle, 90.0, 90.0+angle, 0.0, 0.0) ;
674  
675     Float_t r = fGeom->GetIPtoOuterCoverDistance() + ( fGeom->GetOuterBoxSize(1) + fGeom->GetPPSDBoxSize(1) ) / 2.0 ;
676
677     Float_t xP1 = r * TMath::Sin( angle / kRADDEG ) ;
678     Float_t yP1 = -r * TMath::Cos( angle / kRADDEG ) ;
679
680     gMC->Gspos("PHOS", i, "ALIC", xP1, yP1, 0.0, idrotm[i-1], "ONLY") ;
681  
682   } // for GetNModules
683
684 }
685
686 //____________________________________________________________________________
687 void AliPHOSv0::CreateGeometryforPHOS()
688 {
689   // Create the PHOS-EMC geometry for GEANT
690     //BEGIN_HTML
691   /*
692     <H2>
693     Geant3 geometry tree of PHOS-EMC in ALICE
694     </H2>
695     <P><CENTER>
696     <IMG Align=BOTTOM ALT="EMC geant tree" SRC="../images/EMCinAlice.gif"> 
697     </CENTER><P>
698   */
699   //END_HTML  
700   
701   // Get pointer to the array containing media indexes
702   Int_t *idtmed = fIdtmed->GetArray() - 699 ;
703
704   // ---
705   // --- Define PHOS box volume, fPUFPill with thermo insulating foam ---
706   // --- Foam Thermo Insulating outer cover dimensions ---
707   // --- Put it in bigbox = PHOS
708
709   Float_t dphos[3] ; 
710   dphos[0] =  fGeom->GetOuterBoxSize(0) / 2.0 ;
711   dphos[1] =  fGeom->GetOuterBoxSize(1) / 2.0 ;
712   dphos[2] =  fGeom->GetOuterBoxSize(2) / 2.0 ;
713
714   gMC->Gsvolu("EMCA", "BOX ", idtmed[706], dphos, 3) ;
715
716   Float_t yO =  - fGeom->GetPPSDBoxSize(1)  / 2.0 ;
717
718   gMC->Gspos("EMCA", 1, "PHOS", 0.0, yO, 0.0, 0, "ONLY") ; 
719
720   // ---
721   // --- Define Textolit Wall box, position inside EMCA ---
722   // --- Textolit Wall box dimentions ---
723  
724  
725   Float_t dptxw[3];
726   dptxw[0] = fGeom->GetTextolitBoxSize(0) / 2.0 ;
727   dptxw[1] = fGeom->GetTextolitBoxSize(1) / 2.0 ;
728   dptxw[2] = fGeom->GetTextolitBoxSize(2) / 2.0 ;
729
730   gMC->Gsvolu("PTXW", "BOX ", idtmed[707], dptxw, 3);
731
732   yO =   (  fGeom->GetOuterBoxThickness(1) -   fGeom->GetUpperPlateThickness() ) / 2.  ;
733    
734   gMC->Gspos("PTXW", 1, "EMCA", 0.0, yO, 0.0, 0, "ONLY") ;
735
736   // --- 
737   // --- Define Upper Polystyrene Foam Plate, place inside PTXW ---
738   // --- immediately below Foam Thermo Insulation Upper plate ---
739
740   // --- Upper Polystyrene Foam plate thickness ---
741  
742   Float_t  dpufp[3] ;
743   dpufp[0] = fGeom->GetTextolitBoxSize(0) / 2.0 ; 
744   dpufp[1] = fGeom->GetSecondUpperPlateThickness() / 2. ;
745   dpufp[2] = fGeom->GetTextolitBoxSize(2) /2.0 ; 
746
747   gMC->Gsvolu("PUFP", "BOX ", idtmed[703], dpufp, 3) ;
748   
749   yO = ( fGeom->GetTextolitBoxSize(1) -  fGeom->GetSecondUpperPlateThickness() ) / 2.0 ;
750   
751   gMC->Gspos("PUFP", 1, "PTXW", 0.0, yO, 0.0, 0, "ONLY") ;
752   
753   // ---
754   // --- Define air-filled box, place inside PTXW ---
755   // --- Inner AIR volume dimensions ---
756  
757
758   Float_t  dpair[3] ;
759   dpair[0] = fGeom->GetAirFilledBoxSize(0) / 2.0 ;
760   dpair[1] = fGeom->GetAirFilledBoxSize(1) / 2.0 ;
761   dpair[2] = fGeom->GetAirFilledBoxSize(2) / 2.0 ;
762
763   gMC->Gsvolu("PAIR", "BOX ", idtmed[798], dpair, 3) ;
764   
765   yO = ( fGeom->GetTextolitBoxSize(1) -  fGeom->GetAirFilledBoxSize(1) ) / 2.0 -   fGeom->GetSecondUpperPlateThickness() ;
766   
767   gMC->Gspos("PAIR", 1, "PTXW", 0.0, yO, 0.0, 0, "ONLY") ;
768
769 // --- Dimensions of PbWO4 crystal ---
770
771   Float_t xtlX =  fGeom->GetCrystalSize(0) ; 
772   Float_t xtlY =  fGeom->GetCrystalSize(1) ; 
773   Float_t xtlZ =  fGeom->GetCrystalSize(2) ; 
774
775   Float_t dptcb[3] ;  
776   dptcb[0] =  fGeom->GetNPhi() * ( xtlX + 2 *  fGeom->GetGapBetweenCrystals() ) / 2.0 + fGeom->GetModuleBoxThickness() ;
777   dptcb[1] = ( xtlY +  fGeom->GetCrystalSupportHeight() +  fGeom->GetCrystalWrapThickness() + fGeom->GetCrystalHolderThickness() ) / 2.0 
778              + fGeom->GetModuleBoxThickness() / 2.0 ;
779   dptcb[2] = fGeom->GetNZ() * ( xtlZ + 2 * fGeom->GetGapBetweenCrystals() ) / 2.0 +  fGeom->GetModuleBoxThickness() ;
780   
781   gMC->Gsvolu("PTCB", "BOX ", idtmed[706], dptcb, 3) ;
782
783   yO =  fGeom->GetAirFilledBoxSize(1) / 2.0 - dptcb[1] 
784        - ( fGeom->GetIPtoCrystalSurface() - fGeom->GetIPtoOuterCoverDistance() - fGeom->GetModuleBoxThickness() 
785        -  fGeom->GetUpperPlateThickness() -  fGeom->GetSecondUpperPlateThickness() ) ;
786   
787   gMC->Gspos("PTCB", 1, "PAIR", 0.0, yO, 0.0, 0, "ONLY") ;
788
789   // ---
790   // --- Define Crystal BLock filled with air, position it inside PTCB ---
791   Float_t dpcbl[3] ; 
792   
793   dpcbl[0] = fGeom->GetNPhi() * ( xtlX + 2 * fGeom->GetGapBetweenCrystals() ) / 2.0 ;
794   dpcbl[1] = ( xtlY + fGeom->GetCrystalSupportHeight() + fGeom->GetCrystalWrapThickness() + fGeom->GetCrystalHolderThickness() ) / 2.0 ;
795   dpcbl[2] = fGeom->GetNZ() * ( xtlZ + 2 * fGeom->GetGapBetweenCrystals() ) / 2.0 ;
796   
797   gMC->Gsvolu("PCBL", "BOX ", idtmed[798], dpcbl, 3) ;
798   
799   // --- Divide PCBL in X (phi) and Z directions --
800   gMC->Gsdvn("PROW", "PCBL", Int_t (fGeom->GetNPhi()), 1) ;
801   gMC->Gsdvn("PCEL", "PROW", Int_t (fGeom->GetNZ()), 3) ;
802
803   yO = -fGeom->GetModuleBoxThickness() / 2.0 ;
804   
805   gMC->Gspos("PCBL", 1, "PTCB", 0.0, yO, 0.0, 0, "ONLY") ;
806
807   // ---
808   // --- Define STeel (actually, it's titanium) Cover volume, place inside PCEL
809   Float_t  dpstc[3] ; 
810   
811   dpstc[0] = ( xtlX + 2 * fGeom->GetCrystalWrapThickness() ) / 2.0 ;
812   dpstc[1] = ( xtlY + fGeom->GetCrystalSupportHeight() + fGeom->GetCrystalWrapThickness() + fGeom->GetCrystalHolderThickness() ) / 2.0 ;
813   dpstc[2] = ( xtlZ + 2 * fGeom->GetCrystalWrapThickness()  + 2 *  fGeom->GetCrystalHolderThickness() ) / 2.0 ;
814   
815   gMC->Gsvolu("PSTC", "BOX ", idtmed[704], dpstc, 3) ;
816
817   gMC->Gspos("PSTC", 1, "PCEL", 0.0, 0.0, 0.0, 0, "ONLY") ;
818
819   // ---
820   // --- Define Tyvek volume, place inside PSTC ---
821   Float_t  dppap[3] ;
822
823   dppap[0] = xtlX / 2.0 + fGeom->GetCrystalWrapThickness() ;
824   dppap[1] = ( xtlY + fGeom->GetCrystalSupportHeight() + fGeom->GetCrystalWrapThickness() ) / 2.0 ;
825   dppap[2] = xtlZ / 2.0 + fGeom->GetCrystalWrapThickness() ;
826   
827   gMC->Gsvolu("PPAP", "BOX ", idtmed[702], dppap, 3) ;
828   
829   yO = ( xtlY + fGeom->GetCrystalSupportHeight() + fGeom->GetCrystalWrapThickness() ) / 2.0 
830               - ( xtlY +  fGeom->GetCrystalSupportHeight() +  fGeom->GetCrystalWrapThickness() + fGeom->GetCrystalHolderThickness() ) / 2.0 ;
831    
832   gMC->Gspos("PPAP", 1, "PSTC", 0.0, yO, 0.0, 0, "ONLY") ;
833
834   // ---
835   // --- Define PbWO4 crystal volume, place inside PPAP ---
836   Float_t  dpxtl[3] ; 
837
838   dpxtl[0] = xtlX / 2.0 ;
839   dpxtl[1] = xtlY / 2.0 ;
840   dpxtl[2] = xtlZ / 2.0 ;
841   
842   gMC->Gsvolu("PXTL", "BOX ", idtmed[699], dpxtl, 3) ;
843
844   yO = ( xtlY + fGeom->GetCrystalSupportHeight() + fGeom->GetCrystalWrapThickness() ) / 2.0 - xtlY / 2.0 - fGeom->GetCrystalWrapThickness() ;
845   
846   gMC->Gspos("PXTL", 1, "PPAP", 0.0, yO, 0.0, 0, "ONLY") ;
847
848   // ---
849   // --- Define crystal support volume, place inside PPAP ---
850   Float_t dpsup[3] ; 
851
852   dpsup[0] = xtlX / 2.0 + fGeom->GetCrystalWrapThickness()  ;
853   dpsup[1] = fGeom->GetCrystalSupportHeight() / 2.0 ;
854   dpsup[2] = xtlZ / 2.0 +  fGeom->GetCrystalWrapThickness() ;
855
856   gMC->Gsvolu("PSUP", "BOX ", idtmed[798], dpsup, 3) ;
857
858   yO =  fGeom->GetCrystalSupportHeight() / 2.0 - ( xtlY +  fGeom->GetCrystalSupportHeight() + fGeom->GetCrystalWrapThickness() ) / 2.0 ;
859
860   gMC->Gspos("PSUP", 1, "PPAP", 0.0, yO, 0.0, 0, "ONLY") ;
861
862   // ---
863   // --- Define PIN-diode volume and position it inside crystal support ---
864   // --- right behind PbWO4 crystal
865
866   // --- PIN-diode dimensions ---
867
868  
869   Float_t dppin[3] ;
870   dppin[0] = fGeom->GetPinDiodeSize(0) / 2.0 ;
871   dppin[1] = fGeom->GetPinDiodeSize(1) / 2.0 ;
872   dppin[2] = fGeom->GetPinDiodeSize(2) / 2.0 ;
873  
874   gMC->Gsvolu("PPIN", "BOX ", idtmed[705], dppin, 3) ;
875  
876   yO = fGeom->GetCrystalSupportHeight() / 2.0 - fGeom->GetPinDiodeSize(1) / 2.0 ;
877  
878   gMC->Gspos("PPIN", 1, "PSUP", 0.0, yO, 0.0, 0, "ONLY") ;
879
880   // ---
881   // --- Define Upper Cooling Panel, place it on top of PTCB ---
882   Float_t dpucp[3] ;
883  // --- Upper Cooling Plate thickness ---
884  
885   dpucp[0] = dptcb[0] ;
886   dpucp[1] = fGeom->GetUpperCoolingPlateThickness() ;
887   dpucp[2] = dptcb[2] ;
888   
889   gMC->Gsvolu("PUCP", "BOX ", idtmed[701], dpucp,3) ;
890   
891   yO = (  fGeom->GetAirFilledBoxSize(1) -  fGeom->GetUpperCoolingPlateThickness() ) / 2. 
892        - ( fGeom->GetIPtoCrystalSurface() - fGeom->GetIPtoOuterCoverDistance() - fGeom->GetModuleBoxThickness()
893            - fGeom->GetUpperPlateThickness() - fGeom->GetSecondUpperPlateThickness() - fGeom->GetUpperCoolingPlateThickness() ) ; 
894   
895   gMC->Gspos("PUCP", 1, "PAIR", 0.0, yO, 0.0, 0, "ONLY") ;
896
897   // ---
898   // --- Define Al Support Plate, position it inside PAIR ---
899   // --- right beneath PTCB ---
900  // --- Al Support Plate thickness ---
901  
902   Float_t dpasp[3] ;
903   dpasp[0] =  fGeom->GetAirFilledBoxSize(0) / 2.0 ;
904   dpasp[1] = fGeom->GetSupportPlateThickness() / 2.0 ;
905   dpasp[2] =  fGeom->GetAirFilledBoxSize(2) / 2.0 ;
906   
907   gMC->Gsvolu("PASP", "BOX ", idtmed[701], dpasp, 3) ;
908   
909   yO = (  fGeom->GetAirFilledBoxSize(1) - fGeom->GetSupportPlateThickness() ) / 2. 
910        -  ( fGeom->GetIPtoCrystalSurface() - fGeom->GetIPtoOuterCoverDistance()
911            - fGeom->GetUpperPlateThickness() - fGeom->GetSecondUpperPlateThickness() + dpcbl[1] * 2 ) ;
912   
913   gMC->Gspos("PASP", 1, "PAIR", 0.0, yO, 0.0, 0, "ONLY") ;
914
915   // ---
916   // --- Define Thermo Insulating Plate, position it inside PAIR ---
917   // --- right beneath PASP ---
918   // --- Lower Thermo Insulating Plate thickness ---
919   
920   Float_t dptip[3] ;
921   dptip[0] = fGeom->GetAirFilledBoxSize(0) / 2.0 ;
922   dptip[1] = fGeom->GetLowerThermoPlateThickness() / 2.0 ;
923   dptip[2] = fGeom->GetAirFilledBoxSize(2) / 2.0 ;
924
925   gMC->Gsvolu("PTIP", "BOX ", idtmed[706], dptip, 3) ;
926
927   yO =  ( fGeom->GetAirFilledBoxSize(1) - fGeom->GetLowerThermoPlateThickness() ) / 2. 
928        -  ( fGeom->GetIPtoCrystalSurface() - fGeom->GetIPtoOuterCoverDistance() - fGeom->GetUpperPlateThickness() 
929             - fGeom->GetSecondUpperPlateThickness() + dpcbl[1] * 2 + fGeom->GetSupportPlateThickness() ) ;
930
931   gMC->Gspos("PTIP", 1, "PAIR", 0.0, yO, 0.0, 0, "ONLY") ;
932
933   // ---
934   // --- Define Textolit Plate, position it inside PAIR ---
935   // --- right beneath PTIP ---
936   // --- Lower Textolit Plate thickness ---
937  
938   Float_t dptxp[3] ;
939   dptxp[0] = fGeom->GetAirFilledBoxSize(0) / 2.0 ;
940   dptxp[1] = fGeom->GetLowerTextolitPlateThickness() / 2.0 ;
941   dptxp[2] = fGeom->GetAirFilledBoxSize(2) / 2.0 ;
942
943   gMC->Gsvolu("PTXP", "BOX ", idtmed[707], dptxp, 3) ;
944
945   yO =  ( fGeom->GetAirFilledBoxSize(1) - fGeom->GetLowerTextolitPlateThickness() ) / 2. 
946        -  ( fGeom->GetIPtoCrystalSurface() - fGeom->GetIPtoOuterCoverDistance() - fGeom->GetUpperPlateThickness() 
947             - fGeom->GetSecondUpperPlateThickness() + dpcbl[1] * 2 + fGeom->GetSupportPlateThickness() 
948             +  fGeom->GetLowerThermoPlateThickness() ) ;
949
950   gMC->Gspos("PTXP", 1, "PAIR", 0.0, yO, 0.0, 0, "ONLY") ;
951
952 }
953
954 //____________________________________________________________________________
955 void AliPHOSv0::CreateGeometryforPPSD()
956 {
957   // Create the PHOS-PPSD geometry for GEANT
958
959   //BEGIN_HTML
960   /*
961     <H2>
962     Geant3 geometry tree of PHOS-PPSD in ALICE
963     </H2>
964     <P><CENTER>
965     <IMG Align=BOTTOM ALT="PPSD geant tree" SRC="../images/PPSDinAlice.gif"> 
966     </CENTER><P>
967   */
968   //END_HTML  
969
970   // Get pointer to the array containing media indexes
971   Int_t *idtmed = fIdtmed->GetArray() - 699 ;
972   
973   // The box containing all ppsd's for one PHOS module filled with air 
974   Float_t ppsd[3] ; 
975   ppsd[0] = fGeom->GetPPSDBoxSize(0) / 2.0 ;  
976   ppsd[1] = fGeom->GetPPSDBoxSize(1) / 2.0 ; 
977   ppsd[2] = fGeom->GetPPSDBoxSize(2) / 2.0 ;
978
979   gMC->Gsvolu("PPSD", "BOX ", idtmed[798], ppsd, 3) ;
980
981   Float_t yO =  fGeom->GetOuterBoxSize(1) / 2.0 ;
982
983   gMC->Gspos("PPSD", 1, "PHOS", 0.0, yO, 0.0, 0, "ONLY") ; 
984
985   // Now we build a micromegas module
986   // The box containing the whole module filled with epoxy (FR4)
987
988   Float_t mppsd[3] ;  
989   mppsd[0] = fGeom->GetPPSDModuleSize(0) / 2.0 ;  
990   mppsd[1] = fGeom->GetPPSDModuleSize(1) / 2.0 ;  
991   mppsd[2] = fGeom->GetPPSDModuleSize(2) / 2.0 ;
992
993   gMC->Gsvolu("MPPS", "BOX ", idtmed[708], mppsd, 3) ;  
994  
995   // Inside mppsd :
996   // 1. The Top Lid made of epoxy (FR4) 
997
998   Float_t tlppsd[3] ; 
999   tlppsd[0] = fGeom->GetPPSDModuleSize(0) / 2.0 ; 
1000   tlppsd[1] = fGeom->GetLidThickness() / 2.0 ;
1001   tlppsd[2] = fGeom->GetPPSDModuleSize(2) / 2.0 ;
1002
1003   gMC->Gsvolu("TLPS", "BOX ", idtmed[708], tlppsd, 3) ; 
1004
1005   Float_t  y0 = ( fGeom->GetMicromegas1Thickness() - fGeom->GetLidThickness() ) / 2. ; 
1006
1007   gMC->Gspos("TLPS", 1, "MPPS", 0.0, y0, 0.0, 0, "ONLY") ; 
1008  
1009   // 2. the upper panel made of composite material
1010
1011   Float_t upppsd[3] ; 
1012   upppsd[0] = ( fGeom->GetPPSDModuleSize(0) - fGeom->GetMicromegasWallThickness() ) / 2.0 ;
1013   upppsd[1] = fGeom->GetCompositeThickness() / 2.0 ;
1014   upppsd[2] = ( fGeom->GetPPSDModuleSize(2) - fGeom->GetMicromegasWallThickness() ) / 2.0 ;
1015  
1016   gMC->Gsvolu("UPPS", "BOX ", idtmed[709], upppsd, 3) ; 
1017   
1018   y0 = y0 - fGeom->GetLidThickness() / 2. - fGeom->GetCompositeThickness() / 2. ; 
1019
1020   gMC->Gspos("UPPS", 1, "MPPS", 0.0, y0, 0.0, 0, "ONLY") ; 
1021
1022   // 3. the anode made of Copper
1023   
1024   Float_t anppsd[3] ; 
1025   anppsd[0] = ( fGeom->GetPPSDModuleSize(0) - fGeom->GetMicromegasWallThickness() ) / 2.0 ; 
1026   anppsd[1] = fGeom->GetAnodeThickness() / 2.0 ; 
1027   anppsd[2] = ( fGeom->GetPPSDModuleSize(2) - fGeom->GetMicromegasWallThickness() ) / 2.0  ; 
1028
1029   gMC->Gsvolu("ANPS", "BOX ", idtmed[710], anppsd, 3) ; 
1030   
1031   y0 = y0 - fGeom->GetCompositeThickness() / 2. - fGeom->GetAnodeThickness()  / 2. ; 
1032   
1033   gMC->Gspos("ANPS", 1, "MPPS", 0.0, y0, 0.0, 0, "ONLY") ; 
1034
1035   // 4. the conversion gap + avalanche gap filled with gas
1036
1037   Float_t ggppsd[3] ; 
1038   ggppsd[0] = ( fGeom->GetPPSDModuleSize(0) - fGeom->GetMicromegasWallThickness() ) / 2.0 ;
1039   ggppsd[1] = ( fGeom->GetConversionGap() +  fGeom->GetAvalancheGap() ) / 2.0 ; 
1040   ggppsd[2] = ( fGeom->GetPPSDModuleSize(2) - fGeom->GetMicromegasWallThickness() ) / 2.0 ;
1041
1042   gMC->Gsvolu("GGPS", "BOX ", idtmed[715], ggppsd, 3) ; 
1043   
1044   // --- Divide GGPP in X (phi) and Z directions --
1045   gMC->Gsdvn("GROW", "GGPS", fGeom->GetNumberOfPadsPhi(), 1) ;
1046   gMC->Gsdvn("GCEL", "GROW", fGeom->GetNumberOfPadsZ() , 3) ;
1047
1048   y0 = y0 - fGeom->GetAnodeThickness() / 2.  - ( fGeom->GetConversionGap() +  fGeom->GetAvalancheGap() ) / 2. ; 
1049
1050   gMC->Gspos("GGPS", 1, "MPPS", 0.0, y0, 0.0, 0, "ONLY") ; 
1051
1052
1053   // 6. the cathode made of Copper
1054
1055   Float_t cappsd[3] ;
1056   cappsd[0] = ( fGeom->GetPPSDModuleSize(0) - fGeom->GetMicromegasWallThickness() ) / 2.0 ;
1057   cappsd[1] = fGeom->GetCathodeThickness() / 2.0 ; 
1058   cappsd[2] = ( fGeom->GetPPSDModuleSize(2) - fGeom->GetMicromegasWallThickness() ) / 2.0  ;
1059
1060   gMC->Gsvolu("CAPS", "BOX ", idtmed[710], cappsd, 3) ; 
1061
1062   y0 = y0 - ( fGeom->GetAvalancheGap() +  fGeom->GetAvalancheGap() ) / 2. - fGeom->GetCathodeThickness()  / 2. ; 
1063
1064   gMC->Gspos("CAPS", 1, "MPPS", 0.0, y0, 0.0, 0, "ONLY") ; 
1065
1066   // 7. the printed circuit made of G10       
1067
1068   Float_t pcppsd[3] ; 
1069   pcppsd[0] = ( fGeom->GetPPSDModuleSize(0) - fGeom->GetMicromegasWallThickness() ) / 2,.0 ; 
1070   pcppsd[1] = fGeom->GetPCThickness() / 2.0 ; 
1071   pcppsd[2] = ( fGeom->GetPPSDModuleSize(2) - fGeom->GetMicromegasWallThickness() ) / 2.0 ;
1072
1073   gMC->Gsvolu("PCPS", "BOX ", idtmed[711], cappsd, 3) ; 
1074
1075   y0 = y0 - fGeom->GetCathodeThickness() / 2. - fGeom->GetPCThickness()  / 2. ; 
1076
1077   gMC->Gspos("PCPS", 1, "MPPS", 0.0, y0, 0.0, 0, "ONLY") ; 
1078
1079   // 8. the lower panel made of composite material
1080                                                     
1081   Float_t lpppsd[3] ; 
1082   lpppsd[0] = ( fGeom->GetPPSDModuleSize(0) - fGeom->GetMicromegasWallThickness() ) / 2.0 ; 
1083   lpppsd[1] = fGeom->GetCompositeThickness() / 2.0 ; 
1084   lpppsd[2] = ( fGeom->GetPPSDModuleSize(2) - fGeom->GetMicromegasWallThickness() ) / 2.0 ;
1085
1086   gMC->Gsvolu("LPPS", "BOX ", idtmed[709], lpppsd, 3) ; 
1087  
1088   y0 = y0 - fGeom->GetPCThickness() / 2. - fGeom->GetCompositeThickness()  / 2. ; 
1089
1090   gMC->Gspos("LPPS", 1, "MPPS", 0.0, y0, 0.0, 0, "ONLY") ; 
1091
1092   // Position the  fNumberOfModulesPhi x fNumberOfModulesZ modules (mppsd) inside PPSD to cover a PHOS module
1093   // the top and bottom one's (which are assumed identical) :
1094
1095    Float_t yt = ( fGeom->GetPPSDBoxSize(1) - fGeom->GetMicromegas1Thickness() ) / 2. ; 
1096    Float_t yb = - ( fGeom->GetPPSDBoxSize(1) - fGeom->GetMicromegas2Thickness() ) / 2. ; 
1097
1098    Int_t copyNumbertop = 0 ; 
1099    Int_t copyNumberbot = fGeom->GetNumberOfModulesPhi() *  fGeom->GetNumberOfModulesZ() ; 
1100
1101    Float_t x  = ( fGeom->GetPPSDBoxSize(0) - fGeom->GetPPSDModuleSize(0) ) / 2. ;  
1102
1103    for ( Int_t iphi = 1; iphi <= fGeom->GetNumberOfModulesPhi(); iphi++ ) { // the number of micromegas modules in phi per PHOS module
1104       Float_t z = ( fGeom->GetPPSDBoxSize(2) - fGeom->GetPPSDModuleSize(2) ) / 2. ;
1105
1106       for ( Int_t iz = 1; iz <= fGeom->GetNumberOfModulesZ(); iz++ ) { // the number of micromegas modules in z per PHOS module
1107         gMC->Gspos("MPPS", ++copyNumbertop, "PPSD", x, yt, z, 0, "ONLY") ;
1108         gMC->Gspos("MPPS", ++copyNumberbot, "PPSD", x, yb, z, 0, "ONLY") ; 
1109         z = z - fGeom->GetPPSDModuleSize(2) ;
1110       } // end of Z module loop   
1111       x = x -  fGeom->GetPPSDModuleSize(0) ; 
1112     } // end of phi module loop
1113
1114    // The Lead converter between two air gaps
1115    // 1. Upper air gap
1116
1117    Float_t uappsd[3] ;
1118    uappsd[0] = fGeom->GetPPSDBoxSize(0) / 2.0 ;
1119    uappsd[1] = fGeom->GetMicro1ToLeadGap() / 2.0 ; 
1120    uappsd[2] = fGeom->GetPPSDBoxSize(2) / 2.0 ;
1121
1122   gMC->Gsvolu("UAPPSD", "BOX ", idtmed[798], uappsd, 3) ; 
1123
1124   y0 = ( fGeom->GetPPSDBoxSize(1) - 2 * fGeom->GetMicromegas1Thickness() - fGeom->GetMicro1ToLeadGap() ) / 2. ; 
1125
1126   gMC->Gspos("UAPPSD", 1, "PPSD", 0.0, y0, 0.0, 0, "ONLY") ; 
1127
1128    // 2. Lead converter
1129  
1130   Float_t lcppsd[3] ; 
1131   lcppsd[0] = fGeom->GetPPSDBoxSize(0) / 2.0 ;
1132   lcppsd[1] = fGeom->GetLeadConverterThickness() / 2.0 ; 
1133   lcppsd[2] = fGeom->GetPPSDBoxSize(2) / 2.0 ;
1134  
1135   gMC->Gsvolu("LCPPSD", "BOX ", idtmed[712], lcppsd, 3) ; 
1136   
1137   y0 = y0 - fGeom->GetMicro1ToLeadGap() / 2. - fGeom->GetLeadConverterThickness() / 2. ; 
1138
1139   gMC->Gspos("LCPPSD", 1, "PPSD", 0.0, y0, 0.0, 0, "ONLY") ; 
1140
1141   // 3. Lower air gap
1142
1143   Float_t lappsd[3] ; 
1144   lappsd[0] = fGeom->GetPPSDBoxSize(0) / 2.0 ; 
1145   lappsd[1] = fGeom->GetLeadToMicro2Gap() / 2.0 ; 
1146   lappsd[2] = fGeom->GetPPSDBoxSize(2) / 2.0 ;
1147
1148   gMC->Gsvolu("LAPPSD", "BOX ", idtmed[798], lappsd, 3) ; 
1149     
1150   y0 = y0 - fGeom->GetLeadConverterThickness() / 2. - fGeom->GetLeadToMicro2Gap()  / 2. ; 
1151   
1152   gMC->Gspos("LAPPSD", 1, "PPSD", 0.0, y0, 0.0, 0, "ONLY") ; 
1153    
1154 }
1155
1156 //___________________________________________________________________________
1157 Int_t AliPHOSv0::Digitize(Float_t Energy)
1158 {
1159   // Applies the energy calibration
1160   
1161   Float_t fB = 100000000. ;
1162   Float_t fA = 0. ;
1163   Int_t chan = Int_t(fA + Energy*fB ) ;
1164   return chan ;
1165 }
1166
1167 //___________________________________________________________________________
1168 void AliPHOSv0::FinishEvent()
1169 {
1170   // Makes the digits from the sum of summed hit in a single crystal or PPSD gas cell
1171   // Adds to the energy the electronic noise
1172   // Keeps digits with energy above fDigitThreshold
1173
1174   // Save the cumulated hits instead of raw hits (need to create the branch myself)
1175   // It is put in the Digit Tree because the TreeH is filled after each primary
1176   // and the TreeD at the end of the event.
1177   
1178   
1179   Int_t i ;
1180   Int_t relid[4];
1181   Int_t j ; 
1182   TClonesArray &lDigits = *fDigits ;
1183   AliPHOSHit  * hit ;
1184   AliPHOSDigit * newdigit ;
1185   AliPHOSDigit * curdigit ;
1186   Bool_t deja = kFALSE ; 
1187   
1188   for ( i = 0 ; i < fNTmpHits ; i++ ) {
1189     hit = (AliPHOSHit*)fTmpHits->At(i) ;
1190     newdigit = new AliPHOSDigit( hit->GetPrimary(), hit->GetId(), Digitize( hit->GetEnergy() ) ) ;
1191     deja =kFALSE ;
1192     for ( j = 0 ; j < fNdigits ;  j++) { 
1193       curdigit = (AliPHOSDigit*) lDigits[j] ;
1194       if ( *curdigit == *newdigit) {
1195         *curdigit = *curdigit + *newdigit ; 
1196         deja = kTRUE ; 
1197       }
1198     }
1199     if ( !deja ) {
1200       new(lDigits[fNdigits]) AliPHOSDigit(* newdigit) ;
1201       fNdigits++ ;  
1202     }
1203  
1204     delete newdigit ;    
1205   } 
1206   
1207   // Noise induced by the PIN diode of the PbWO crystals
1208
1209   Float_t energyandnoise ;
1210   for ( i = 0 ; i < fNdigits ; i++ ) {
1211     newdigit =  (AliPHOSDigit * ) fDigits->At(i) ;
1212     fGeom->AbsToRelNumbering(newdigit->GetId(), relid) ;
1213
1214     if (relid[1]==0){   // Digits belong to EMC (PbW0_4 crystals)
1215       energyandnoise = newdigit->GetAmp() + Digitize(gRandom->Gaus(0., fPinElectronicNoise)) ;
1216
1217       if (energyandnoise < 0 ) 
1218         energyandnoise = 0 ;
1219
1220       if ( newdigit->GetAmp() < fDigitThreshold ) // if threshold not surpassed, remove digit from list
1221         fDigits->RemoveAt(i) ; 
1222     }
1223   }
1224   
1225   fDigits->Compress() ;  
1226
1227   fNdigits =  fDigits->GetEntries() ; 
1228   for (i = 0 ; i < fNdigits ; i++) { 
1229     newdigit = (AliPHOSDigit *) fDigits->At(i) ; 
1230     newdigit->SetIndexInList(i) ; 
1231   }
1232   
1233 }
1234
1235 //____________________________________________________________________________
1236 void AliPHOSv0::Init(void)
1237 {
1238   // Just prints an information message
1239   
1240   Int_t i;
1241
1242   printf("\n");
1243   for(i=0;i<35;i++) printf("*");
1244   printf(" PHOS_INIT ");
1245   for(i=0;i<35;i++) printf("*");
1246   printf("\n");
1247
1248   // Here the PHOS initialisation code (if any!)
1249
1250   for(i=0;i<80;i++) printf("*");
1251   printf("\n");
1252   
1253 }
1254
1255 //___________________________________________________________________________
1256 void AliPHOSv0::MakeBranch(Option_t* opt)
1257 {  
1258   // Create new branche in the current Root Tree in the digit Tree
1259
1260   AliDetector::MakeBranch(opt) ;
1261   
1262   char branchname[10];
1263   sprintf(branchname,"%s",GetName());
1264   char *cdD = strstr(opt,"D");
1265   if (fDigits && gAlice->TreeD() && cdD) {
1266     gAlice->TreeD()->Branch(branchname, &fDigits, fBufferSize);
1267   }
1268
1269   // Create new branche PHOSCH in the current Root Tree in the digit Tree for accumulated Hits
1270   if ( ! (gAlice->IsLegoRun()) ) { // only when not in lego plot mode 
1271     if ( fTmpHits && gAlice->TreeD()  && cdD) {
1272       char branchname[10] ;
1273       sprintf(branchname, "%sCH", GetName()) ;
1274       gAlice->TreeD()->Branch(branchname, &fTmpHits, fBufferSize) ;
1275     }   
1276   }
1277
1278 }
1279
1280
1281 //_____________________________________________________________________________
1282 void AliPHOSv0::Reconstruction(AliPHOSReconstructioner * Reconstructioner)
1283
1284   // 1. Reinitializes the existing RecPoint, TrackSegment, and RecParticles Lists and 
1285   // 2. Creates TreeR with a branch for each list
1286   // 3. Steers the reconstruction processes
1287   // 4. Saves the 3 lists in TreeR
1288   // 5. Write the Tree to File
1289   
1290   fReconstructioner = Reconstructioner ;
1291   
1292   char branchname[10] ;
1293   
1294   // 1.
1295
1296   //  gAlice->MakeTree("R") ; 
1297   Int_t splitlevel = 0 ; 
1298   
1299   if (fEmcRecPoints) { 
1300     fEmcRecPoints->Delete() ; 
1301     delete fEmcRecPoints ;
1302     fEmcRecPoints = 0 ; 
1303   }
1304
1305   //  fEmcRecPoints= new AliPHOSRecPoint::RecPointsList("AliPHOSEmcRecPoint", 1000) ; if TClonesArray
1306   fEmcRecPoints= new AliPHOSRecPoint::RecPointsList(1000) ; 
1307
1308   if ( fEmcRecPoints && gAlice->TreeR() ) {
1309     sprintf(branchname,"%sEmcRP",GetName()) ;
1310     
1311     // gAlice->TreeR()->Branch(branchname, &fEmcRecPoints, fBufferSize); if TClonesArray
1312     gAlice->TreeR()->Branch(branchname, "TObjArray", &fEmcRecPoints, fBufferSize, splitlevel) ; 
1313   }
1314
1315   if (fPpsdRecPoints) { 
1316     fPpsdRecPoints->Delete() ; 
1317     delete fPpsdRecPoints ; 
1318     fPpsdRecPoints = 0 ; 
1319   }
1320
1321   //  fPpsdRecPoints = new AliPHOSRecPoint::RecPointsList("AliPHOSPpsdRecPoint", 1000) ; if TClonesArray
1322   fPpsdRecPoints = new AliPHOSRecPoint::RecPointsList(1000) ;
1323
1324   if ( fPpsdRecPoints && gAlice->TreeR() ) {
1325     sprintf(branchname,"%sPpsdRP",GetName()) ;
1326      
1327      // gAlice->TreeR()->Branch(branchname, &fPpsdRecPoints, fBufferSize); if TClonesArray
1328     gAlice->TreeR()->Branch(branchname, "TObjArray", &fPpsdRecPoints, fBufferSize, splitlevel) ;
1329   }
1330
1331   if (fTrackSegments) { 
1332    fTrackSegments->Delete() ; 
1333     delete fTrackSegments ; 
1334     fTrackSegments = 0 ; 
1335   }
1336
1337   fTrackSegments = new AliPHOSTrackSegment::TrackSegmentsList("AliPHOSTrackSegment", 1000) ;
1338   if ( fTrackSegments && gAlice->TreeR() ) { 
1339     sprintf(branchname,"%sTS",GetName()) ;
1340     gAlice->TreeR()->Branch(branchname, &fTrackSegments, fBufferSize) ;
1341   }
1342
1343   if (fRecParticles) {  
1344     fRecParticles->Delete() ; 
1345     delete fRecParticles ; 
1346     fRecParticles = 0 ; 
1347   }
1348   fRecParticles = new AliPHOSRecParticle::RecParticlesList("AliPHOSRecParticle", 1000) ;
1349   if ( fRecParticles && gAlice->TreeR() ) { 
1350      sprintf(branchname,"%sRP",GetName()) ;
1351      gAlice->TreeR()->Branch(branchname, &fRecParticles, fBufferSize) ;
1352   }
1353   
1354   // 3.
1355
1356   fReconstructioner->Make(fDigits, fEmcRecPoints, fPpsdRecPoints, fTrackSegments, fRecParticles);
1357
1358   // 4. Expand or Shrink the arrays to the proper size
1359   
1360   Int_t size ;
1361   
1362   size = fEmcRecPoints->GetEntries() ;
1363   fEmcRecPoints->Expand(size) ;
1364  
1365   size = fPpsdRecPoints->GetEntries() ;
1366   fPpsdRecPoints->Expand(size) ;
1367
1368   size = fTrackSegments->GetEntries() ;
1369   fTrackSegments->Expand(size) ;
1370
1371   size = fRecParticles->GetEntries() ;
1372   fRecParticles->Expand(size) ;
1373
1374   gAlice->TreeR()->Fill() ;
1375   cout << "filled" << endl ;
1376   // 5.
1377
1378   gAlice->TreeR()->Write() ;
1379   cout << "writen" << endl ;
1380  
1381   // Deleting reconstructed objects
1382   ResetReconstruction();
1383
1384   
1385 }
1386
1387 //____________________________________________________________________________
1388 void AliPHOSv0::ResetDigits() 
1389
1390   // May sound strange, but cumulative hits are store in digits Tree
1391   AliDetector::ResetDigits();
1392   if(  fTmpHits ) {
1393     fTmpHits->Delete();
1394     fNTmpHits = 0 ;
1395   }
1396 }  
1397 //____________________________________________________________________________
1398 void AliPHOSv0::ResetReconstruction() 
1399
1400   // Deleting reconstructed objects
1401
1402   if ( fEmcRecPoints )   fEmcRecPoints->Delete();
1403   if ( fPpsdRecPoints )  fPpsdRecPoints->Delete();
1404   if ( fTrackSegments )  fTrackSegments->Delete();
1405   if ( fRecParticles )   fRecParticles->Delete();
1406   
1407 }
1408 //____________________________________________________________________________
1409
1410 //____________________________________________________________________________
1411 void AliPHOSv0::SetTreeAddress()
1412
1413   //  TBranch *branch;
1414   AliPHOS::SetTreeAddress();
1415
1416  //  //Branch address for TreeR: RecPpsdRecPoint
1417 //   TTree *treeR = gAlice->TreeR();
1418 //   if ( treeR && fPpsdRecPoints ) {
1419 //     branch = treeR->GetBranch("PHOSPpsdRP");
1420 //     if (branch) branch->SetAddress(&fPpsdRecPoints) ;
1421   //  }
1422 }
1423
1424 //____________________________________________________________________________
1425
1426 void AliPHOSv0::StepManager(void)
1427 {
1428   // Accumulates hits as long as the track stays in a single crystal or PPSD gas Cell
1429
1430   Int_t          relid[4] ;      // (box, layer, row, column) indices
1431   Float_t        xyze[4] ;       // position wrt MRS and energy deposited
1432   TLorentzVector pos ;
1433   Int_t copy ;
1434
1435   Int_t primary =  gAlice->GetPrimary( gAlice->CurrentTrack() ); 
1436   TString name = fGeom->GetName() ; 
1437   if ( name == "GPS2" ) { // the CPV is a PPSD
1438     if( gMC->CurrentVolID(copy) == gMC->VolId("GCEL") ) // We are inside a gas cell 
1439     {
1440       gMC->TrackPosition(pos) ;
1441       xyze[0] = pos[0] ;
1442       xyze[1] = pos[1] ;
1443       xyze[2] = pos[2] ;
1444       xyze[3] = gMC->Edep() ; 
1445
1446       if ( xyze[3] != 0 ) { // there is deposited energy 
1447         gMC->CurrentVolOffID(5, relid[0]) ;  // get the PHOS Module number
1448         gMC->CurrentVolOffID(3, relid[1]) ;  // get the Micromegas Module number 
1449       // 1-> Geom->GetNumberOfModulesPhi() *  fGeom->GetNumberOfModulesZ() upper                         
1450       //  >  fGeom->GetNumberOfModulesPhi()  *  fGeom->GetNumberOfModulesZ() lower
1451         gMC->CurrentVolOffID(1, relid[2]) ;  // get the row number of the cell
1452         gMC->CurrentVolID(relid[3]) ;        // get the column number 
1453
1454         // get the absolute Id number
1455
1456         Int_t absid ; 
1457         fGeom->RelToAbsNumbering(relid, absid) ; 
1458
1459         // add current hit to the hit list      
1460         AddHit(primary, absid, xyze);
1461
1462       } // there is deposited energy 
1463      } // We are inside the gas of the CPV  
1464    } // GPS2 configuration
1465   
1466    if(gMC->CurrentVolID(copy) == gMC->VolId("PXTL") )  //  We are inside a PBWO crystal
1467      {
1468        gMC->TrackPosition(pos) ;
1469        xyze[0] = pos[0] ;
1470        xyze[1] = pos[1] ;
1471        xyze[2] = pos[2] ;
1472        xyze[3] = gMC->Edep() ;
1473
1474        if ( xyze[3] != 0 ) {
1475           gMC->CurrentVolOffID(10, relid[0]) ; // get the PHOS module number ;
1476           relid[1] = 0   ;                    // means PBW04
1477           gMC->CurrentVolOffID(4, relid[2]) ; // get the row number inside the module
1478           gMC->CurrentVolOffID(3, relid[3]) ; // get the cell number inside the module
1479
1480       // get the absolute Id number
1481
1482           Int_t absid ; 
1483           fGeom->RelToAbsNumbering(relid, absid) ; 
1484  
1485       // add current hit to the hit list
1486
1487           AddHit(primary, absid, xyze);
1488     
1489        } // there is deposited energy
1490     } // we are inside a PHOS Xtal
1491 }
1492