Changes for correct initialization of Geant4 (Mihaela)
[u/mrichter/AliRoot.git] / PHOS / AliPHOS.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15 /* $Id$ */
16 /* History of cvs commits:
17  *
18  * $Log$
19  * Revision 1.108  2007/02/01 10:34:47  hristov
20  * Removing warnings on Solaris x86
21  *
22  * Revision 1.107  2007/01/29 16:29:37  kharlov
23  * Digits2Raw(): special workaround for digits with time out of range
24  *
25  * Revision 1.106  2007/01/17 17:28:56  kharlov
26  * Extract ALTRO sample generation to a separate class AliPHOSPulseGenerator
27  *
28  * Revision 1.105  2007/01/12 21:44:29  kharlov
29  * Simulate and reconstruct two gains simulaneouslsy
30  *
31  * Revision 1.104  2006/11/23 13:40:44  hristov
32  * Common class for raw data reading and ALTRO mappiing for PHOS and EMCAL (Gustavo, Cvetan)
33  *
34  * Revision 1.103  2006/11/14 17:11:15  hristov
35  * Removing inheritances from TAttLine, TAttMarker and AliRndm in AliModule. The copy constructor and assignment operators are moved to the private part of the class and not implemented. The corresponding changes are propagated to the detectors
36  *
37  * Revision 1.102  2006/10/27 17:14:27  kharlov
38  * Introduce AliDebug and AliLog (B.Polichtchouk)
39  *
40  * Revision 1.101  2006/10/13 06:47:29  kharlov
41  * Simulation of RAW data applies real mapping (B.Polichtchouk)
42  *
43  * Revision 1.100  2006/08/11 12:36:26  cvetan
44  * Update of the PHOS code needed in order to read and reconstruct the beam test raw data (i.e. without an existing galice.root)
45  *
46  * Revision 1.99  2006/06/28 11:36:09  cvetan
47  * New detector numbering scheme (common for DAQ/HLT/Offline). All the subdetectors shall use the AliDAQ class for the sim and rec of the raw data. The AliDAQ and raw reader classes now provide all the necessary interfaces to write and select the detector specific raw-data payload. Look into the AliDAQ.h and AliRawReader.h for more details.
48  *
49  * Revision 1.98  2006/05/11 11:30:48  cvetan
50  * Major changes in AliAltroBuffer. Now it can be used only for writing of raw data. All the corresponding read method are removed. It is based now on AliFstream in order to avoid endianess problems. The altro raw data is written always with little endian
51  *
52  * Revision 1.97  2006/04/22 10:30:17  hristov
53  * Add fEnergy to AliPHOSDigit and operate with EMC amplitude in energy units (Yu.Kharlov)
54  *
55  * Revision 1.96  2006/04/07 08:41:59  hristov
56  * Follow AliAlignObj framework and remove AliPHOSAlignData (Yu.Kharlov)
57  *
58  * Revision 1.95  2006/03/14 19:40:41  kharlov
59  * Remove De-digitizing of raw data and digitizing the raw data fit
60  *
61  * Revision 1.94  2006/03/07 18:56:25  kharlov
62  * CDB is passed via environment variable
63  *
64  * Revision 1.93  2005/11/22 08:45:11  kharlov
65  * Calibration is read from CDB if any (Boris Polichtchouk)
66  *
67  * Revision 1.92  2005/11/03 13:09:19  hristov
68  * Removing meaningless const declarations (linuxicc)
69  *
70  * Revision 1.91  2005/07/27 15:08:53  kharlov
71  * Mixture ArCO2 is corrected
72  *
73  * Revision 1.90  2005/06/17 07:39:07  hristov
74  * Removing GetDebug and SetDebug from AliRun and AliModule. Using AliLog for the messages
75  *
76  * Revision 1.89  2005/05/28 12:10:07  schutz
77  * Copy constructor is corrected (by T.P.)
78  *
79  */
80
81 //_________________________________________________________________________
82 // Base Class for PHOS description:
83 //   PHOS consists of a PbWO4 calorimeter (EMCA) and a gazeous charged 
84 //    particles detector (CPV or PPSD).
85 //   The only provided method here is CreateMaterials, 
86 //    which defines the materials common to all PHOS versions.   
87 // 
88 //*-- Author: Laurent Aphecetche & Yves Schutz (SUBATECH) 
89 //////////////////////////////////////////////////////////////////////////////
90
91
92 // --- ROOT system ---
93 class TFile;
94 #include <TFolder.h> 
95 #include <TTree.h>
96 #include <TVirtualMC.h> 
97 #include <TH1F.h> 
98 #include <TF1.h> 
99 #include <TRandom.h> 
100
101 // --- Standard library ---
102
103 // --- AliRoot header files ---
104 #include "AliMagF.h"
105 #include "AliPHOS.h"
106 #include "AliPHOSGetter.h"
107 #include "AliRun.h"
108 #include "AliPHOSDigitizer.h"
109 #include "AliPHOSSDigitizer.h"
110 #include "AliPHOSDigit.h"
111 #include "AliAltroBuffer.h"
112 #include "AliAltroMapping.h"
113 #include "AliCaloAltroMapping.h"
114 #include "AliLog.h"
115 #include "AliCDBManager.h"
116 #include "AliCDBEntry.h"
117 #include "AliCDBStorage.h"
118 #include "AliPHOSCalibData.h"
119 #include "AliPHOSPulseGenerator.h"
120 #include "AliDAQ.h"
121
122 ClassImp(AliPHOS)
123
124 //____________________________________________________________________________
125   AliPHOS:: AliPHOS() : AliDetector()
126 {
127   // Default ctor
128   fName   = "PHOS" ;
129
130 }
131
132 //____________________________________________________________________________
133 AliPHOS::AliPHOS(const char* name, const char* title): AliDetector(name, title)
134 {
135   //   ctor : title is used to identify the layout
136 }
137
138 //____________________________________________________________________________
139 AliPHOS::~AliPHOS() 
140 {  
141 }
142
143 //____________________________________________________________________________
144 AliDigitizer* AliPHOS::CreateDigitizer(AliRunDigitizer* manager) const
145 {
146   return new AliPHOSDigitizer(manager);
147 }
148
149 //____________________________________________________________________________
150 void AliPHOS::CreateMaterials()
151 {
152   // Definitions of materials to build PHOS and associated tracking media.
153   // media number in idtmed are 699 to 798.
154
155   // --- The PbWO4 crystals ---
156   Float_t aX[3] = {207.19, 183.85, 16.0} ;
157   Float_t zX[3] = {82.0, 74.0, 8.0} ;
158   Float_t wX[3] = {1.0, 1.0, 4.0} ;
159   Float_t dX = 8.28 ;
160
161   AliMixture(0, "PbWO4$", aX, zX, dX, -3, wX) ;
162
163
164   // --- The polysterene scintillator (CH) ---
165   Float_t aP[2] = {12.011, 1.00794} ;
166   Float_t zP[2] = {6.0, 1.0} ;
167   Float_t wP[2] = {1.0, 1.0} ;
168   Float_t dP = 1.032 ;
169
170   AliMixture(1, "Polystyrene$", aP, zP, dP, -2, wP) ;
171
172   // --- Aluminium ---
173   AliMaterial(2, "Al$", 26.98, 13., 2.7, 8.9, 999., 0, 0) ;
174   // ---         Absorption length is ignored ^
175
176  // --- Tyvek (CnH2n) ---
177   Float_t aT[2] = {12.011, 1.00794} ;
178   Float_t zT[2] = {6.0, 1.0} ;
179   Float_t wT[2] = {1.0, 2.0} ;
180   Float_t dT = 0.331 ;
181
182   AliMixture(3, "Tyvek$", aT, zT, dT, -2, wT) ;
183
184   // --- Polystyrene foam ---
185   Float_t aF[2] = {12.011, 1.00794} ;
186   Float_t zF[2] = {6.0, 1.0} ;
187   Float_t wF[2] = {1.0, 1.0} ;
188   Float_t dF = 0.12 ;
189
190   AliMixture(4, "Foam$", aF, zF, dF, -2, wF) ;
191
192  // --- Titanium ---
193   Float_t aTIT[3] = {47.88, 26.98, 54.94} ;
194   Float_t zTIT[3] = {22.0, 13.0, 25.0} ;
195   Float_t wTIT[3] = {69.0, 6.0, 1.0} ;
196   Float_t dTIT = 4.5 ;
197
198   AliMixture(5, "Titanium$", aTIT, zTIT, dTIT, -3, wTIT);
199
200  // --- Silicon ---
201   AliMaterial(6, "Si$", 28.0855, 14., 2.33, 9.36, 42.3, 0, 0) ;
202
203
204
205   // --- Foam thermo insulation ---
206   Float_t aTI[2] = {12.011, 1.00794} ;
207   Float_t zTI[2] = {6.0, 1.0} ;
208   Float_t wTI[2] = {1.0, 1.0} ;
209   Float_t dTI = 0.04 ;
210
211   AliMixture(7, "Thermo Insul.$", aTI, zTI, dTI, -2, wTI) ;
212
213   // --- Textolith ---
214   Float_t aTX[4] = {16.0, 28.09, 12.011, 1.00794} ;
215   Float_t zTX[4] = {8.0, 14.0, 6.0, 1.0} ;
216   Float_t wTX[4] = {292.0, 68.0, 462.0, 736.0} ;
217   Float_t dTX    = 1.75 ;
218
219   AliMixture(8, "Textolit$", aTX, zTX, dTX, -4, wTX) ;
220
221   //--- FR4  ---
222   Float_t aFR[4] = {16.0, 28.09, 12.011, 1.00794} ;
223   Float_t zFR[4] = {8.0, 14.0, 6.0, 1.0} ;
224   Float_t wFR[4] = {292.0, 68.0, 462.0, 736.0} ;
225   Float_t dFR = 1.8 ; 
226
227   AliMixture(9, "FR4$", aFR, zFR, dFR, -4, wFR) ;
228
229   // --- The Composite Material for  micromegas (so far polyetylene) ---                                       
230   Float_t aCM[2] = {12.01, 1.} ; 
231   Float_t zCM[2] = {6., 1.} ; 
232   Float_t wCM[2] = {1., 2.} ; 
233   Float_t dCM = 0.935 ; 
234
235   AliMixture(10, "Compo Mat$", aCM, zCM, dCM, -2, wCM) ;
236
237   // --- Copper ---                                                                    
238   AliMaterial(11, "Cu$", 63.546, 29, 8.96, 1.43, 14.8, 0, 0) ;
239  
240   // --- G10 : Printed Circuit material ---                                                  
241   Float_t aG10[4] = { 12., 1., 16., 28.} ;
242   Float_t zG10[4] = { 6., 1., 8., 14.} ;
243   Float_t wG10[4] = { .259, .288, .248, .205} ;
244   Float_t dG10  = 1.7 ;
245   
246   AliMixture(12, "G10$", aG10, zG10, dG10, -4, wG10);
247
248   // --- Lead ---                                                                     
249   AliMaterial(13, "Pb$", 207.2, 82, 11.35, 0.56, 0., 0, 0) ;
250
251  // --- The gas mixture ---                                                                
252  // Co2
253   Float_t aCO[2] = {12.0, 16.0} ; 
254   Float_t zCO[2] = {6.0, 8.0} ; 
255   Float_t wCO[2] = {1.0, 2.0} ; 
256   Float_t dCO = 0.001977 ; 
257
258   AliMixture(14, "CO2$", aCO, zCO, dCO, -2, wCO);
259
260  // Ar
261   Float_t dAr = 0.001782 ; 
262   AliMaterial(15, "Ar$", 39.948, 18.0, dAr, 14.0, 0., 0, 0) ;   
263  
264   // Ar+CO2 Mixture (80% / 20%)
265   Float_t arContent = 0.80 ;  // Ar-content of the ArCO2-mixture
266   Float_t aArCO[3]  = {39.948, 12.0, 16.0} ;
267   Float_t zArCO[3]  = {18.0  ,  6.0,  8.0} ;
268   Float_t wArCO[3];
269   wArCO[0] = arContent;
270   wArCO[1] = (1-arContent)*1;
271   wArCO[2] = (1-arContent)*2;
272   Float_t dArCO = arContent*dAr + (1-arContent)*dCO ;
273   AliMixture(16, "ArCO2$", aArCO, zArCO, dArCO,  -3, wArCO) ;
274
275   // --- Stainless steel (let it be pure iron) ---
276   AliMaterial(17, "Steel$", 55.845, 26, 7.87, 1.76, 0., 0, 0) ;
277
278
279   // --- Fiberglass ---
280   Float_t aFG[4] = {16.0, 28.09, 12.011, 1.00794} ;
281   Float_t zFG[4] = {8.0, 14.0, 6.0, 1.0} ;
282   Float_t wFG[4] = {292.0, 68.0, 462.0, 736.0} ;
283   Float_t dFG    = 1.9 ;
284
285   AliMixture(18, "Fibergla$", aFG, zFG, dFG, -4, wFG) ;
286
287   // --- Cables in Air box  ---
288   // SERVICES
289
290   Float_t aCA[4] = { 1.,12.,55.8,63.5 };
291   Float_t zCA[4] = { 1.,6.,26.,29. }; 
292   Float_t wCA[4] = { .014,.086,.42,.48 };
293   Float_t dCA    = 0.8 ;  //this density is raw estimation, if you know better - correct
294
295   AliMixture(19, "Cables  $", aCA, zCA, dCA, -4, wCA) ;
296
297
298   // --- Air ---
299   Float_t aAir[4]={12.0107,14.0067,15.9994,39.948};
300   Float_t zAir[4]={6.,7.,8.,18.};
301   Float_t wAir[4]={0.000124,0.755267,0.231781,0.012827};
302   Float_t dAir = 1.20479E-3;
303  
304   AliMixture(99, "Air$", aAir, zAir, dAir, 4, wAir) ;
305
306   // DEFINITION OF THE TRACKING MEDIA
307
308   // for PHOS: idtmed[699->798] equivalent to fIdtmed[0->100]
309   Int_t   isxfld = gAlice->Field()->Integ() ;
310   Float_t sxmgmx = gAlice->Field()->Max() ;
311
312   // The scintillator of the calorimeter made of PBW04                              -> idtmed[699]
313   AliMedium(0, "PHOS Xtal    $", 0, 1,
314             isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
315
316   // The scintillator of the CPV made of Polystyrene scintillator                   -> idtmed[700]
317   AliMedium(1, "CPV scint.   $", 1, 1,
318             isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
319
320   // Various Aluminium parts made of Al                                             -> idtmed[701]
321   AliMedium(2, "Al parts     $", 2, 0,
322              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.001, 0.001, 0, 0) ;
323
324   // The Tywek which wraps the calorimeter crystals                                 -> idtmed[702]
325   AliMedium(3, "Tyvek wrapper$", 3, 0,
326              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.001, 0.001, 0, 0) ;
327
328   // The Polystyrene foam around the calorimeter module                             -> idtmed[703]
329   AliMedium(4, "Polyst. foam $", 4, 0,
330              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
331
332   // The Titanium around the calorimeter crystal                                    -> idtmed[704]
333   AliMedium(5, "Titan. cover $", 5, 0,
334              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.0001, 0.0001, 0, 0) ;
335
336   // The Silicon of the pin diode to read out the calorimeter crystal               -> idtmed[705] 
337  AliMedium(6, "Si PIN       $", 6, 0,
338              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.01, 0.01, 0, 0) ;
339
340  // The thermo insulating material of the box which contains the calorimeter module -> idtmed[706]
341   AliMedium(7, "Thermo Insul.$", 7, 0,
342              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
343
344   // The Textolit which makes up the box which contains the calorimeter module      -> idtmed[707]
345   AliMedium(8, "Textolit     $", 8, 0,
346              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
347
348   // FR4: The Plastic which makes up the frame of micromegas                        -> idtmed[708]
349   AliMedium(9, "FR4 $", 9, 0,
350              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.0001, 0, 0) ; 
351
352
353   // The Composite Material for  micromegas                                         -> idtmed[709]
354   AliMedium(10, "CompoMat   $", 10, 0,
355              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
356
357   // Copper                                                                         -> idtmed[710]
358   AliMedium(11, "Copper     $", 11, 0,
359              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.0001, 0, 0) ;
360
361   // G10: Printed Circuit material                                                  -> idtmed[711]
362  
363   AliMedium(12, "G10        $", 12, 0,
364              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.01, 0, 0) ;
365
366   // The Lead                                                                       -> idtmed[712]
367  
368   AliMedium(13, "Lead      $", 13, 0,
369              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
370
371   // The gas mixture: ArCo2                                                         -> idtmed[715]
372  
373   AliMedium(16, "ArCo2      $", 16, 1,
374              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.01, 0, 0) ;
375  
376   // Stainless steel                                                                -> idtmed[716]
377   AliMedium(17, "Steel     $", 17, 0,
378              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.0001, 0, 0) ;
379
380   // Fibergalss                                                                     -> idtmed[717]
381   AliMedium(18, "Fiberglass$", 18, 0,
382              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
383
384   // Cables in air                                                                  -> idtmed[718]
385   AliMedium(19, "Cables    $", 19, 0,
386              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
387
388   // Air                                                                            -> idtmed[798] 
389   AliMedium(99, "Air          $", 99, 0,
390              isxfld, sxmgmx, 10.0, 1.0, 0.1, 0.1, 10.0, 0, 0) ;
391 }
392
393 //_____________________________________________________________________________
394 void AliPHOS::Init()
395 {
396   //
397   // Initialises cuts for PHOS
398   //
399   // --- Set decent energy thresholds for gamma and electron tracking
400   Int_t * idtmed = fIdtmed->GetArray() - 699 ; 
401
402   // Tracking threshold for photons and electrons in the scintillator crystal 
403   gMC->Gstpar(idtmed[699], "CUTGAM",0.5E-4) ; 
404   gMC->Gstpar(idtmed[699], "CUTELE",1.0E-4) ;
405  
406   // --- Generate explicitly delta rays in the titan cover ---
407   gMC->Gstpar(idtmed[704], "LOSS",3.) ;
408   gMC->Gstpar(idtmed[704], "DRAY",1.) ;
409   // --- and in aluminium parts ---
410   gMC->Gstpar(idtmed[701], "LOSS",3.) ;
411   gMC->Gstpar(idtmed[701], "DRAY",1.) ;
412   // --- and in PIN diode
413   gMC->Gstpar(idtmed[705], "LOSS",3) ;
414   gMC->Gstpar(idtmed[705], "DRAY",1) ;
415   // --- and in the passive convertor
416   gMC->Gstpar(idtmed[712], "LOSS",3) ;
417   gMC->Gstpar(idtmed[712], "DRAY",1) ;
418   // Tracking threshold for photons and electrons in the gas ArC02 
419   gMC->Gstpar(idtmed[715], "CUTGAM",1.E-5) ; 
420   gMC->Gstpar(idtmed[715], "CUTELE",1.E-5) ;
421   gMC->Gstpar(idtmed[715], "CUTNEU",1.E-5) ;
422   gMC->Gstpar(idtmed[715], "CUTHAD",1.E-5) ;
423   gMC->Gstpar(idtmed[715], "CUTMUO",1.E-5) ;
424   gMC->Gstpar(idtmed[715], "BCUTE",1.E-5) ;
425   gMC->Gstpar(idtmed[715], "BCUTM",1.E-5) ;
426   gMC->Gstpar(idtmed[715], "DCUTE",1.E-5) ;
427   gMC->Gstpar(idtmed[715], "DCUTM",1.E-5) ;
428   gMC->Gstpar(idtmed[715], "PPCUTM",1.E-5) ;
429   gMC->Gstpar(idtmed[715], "LOSS",2.) ;
430   gMC->Gstpar(idtmed[715], "DRAY",0.) ;
431   gMC->Gstpar(idtmed[715], "STRA",2.) ;
432 }
433
434 //____________________________________________________________________________
435 void AliPHOS::Digits2Raw()
436 {
437 // convert digits of the current event to raw data
438   
439   AliPHOSLoader * loader = dynamic_cast<AliPHOSLoader*>(fLoader) ; 
440
441   // get the digits
442   loader->LoadDigits();
443   TClonesArray* digits = loader->Digits() ;
444
445   if (!digits) {
446     AliError(Form("No digits found !"));
447     return;
448   }
449
450   // get the geometry
451   AliPHOSGeometry* geom = GetGeometry();
452   if (!geom) {
453     AliError(Form("No geometry found !"));
454     return;
455   }
456
457   // some digitization constants
458   const Float_t    kThreshold = 0.001; // skip digits below 1 MeV
459   const Int_t      kAdcThreshold = 1;  // Lower ADC threshold to write to raw data
460
461   AliAltroBuffer* buffer = NULL;
462   Int_t prevDDL = -1;
463
464   // Create a shaper pulse object
465   AliPHOSPulseGenerator *pulse = new AliPHOSPulseGenerator();
466
467   Int_t *adcValuesLow = new Int_t[pulse->GetRawFormatTimeBins()];
468   Int_t *adcValuesHigh= new Int_t[pulse->GetRawFormatTimeBins()];
469
470   //!!!!for debug!!!
471   Int_t modMax=-111;
472   Int_t colMax=-111;
473   Int_t rowMax=-111;
474   Float_t eMax=-333;
475   //!!!for debug!!!
476
477   // loop over digits (assume ordered digits)
478   for (Int_t iDigit = 0; iDigit < digits->GetEntries(); iDigit++) {
479     AliPHOSDigit* digit = dynamic_cast<AliPHOSDigit *>(digits->At(iDigit)) ;
480     if (digit->GetEnergy() < kThreshold) 
481       continue;
482     Int_t relId[4];
483     geom->AbsToRelNumbering(digit->GetId(), relId);
484     Int_t module = relId[0];
485  
486     // Begin FIXME 
487     if (relId[1] != 0) 
488       continue;    // ignore digits from CPV
489     // End FIXME 
490
491     Int_t row = relId[2]-1;
492     Int_t col = relId[3]-1;
493
494     Int_t iRCU = -111;
495
496     //RCU0
497     if(0<=row&&row<32 && 0<=col&&col<28) iRCU=0;
498
499     //RCU1
500     if(0<=row&&row<32 && 28<=col&&col<56) iRCU=1;
501
502     //RCU2
503     if(32<=row&&row<64 && 0<=col&&col<28) iRCU=2;
504
505     //RCU3
506     if(32<=row&&row<64 && 28<=col&&col<56) iRCU=3;
507
508
509     // PHOS EMCA has 4 DDL per module. Splitting is based on the (row,column) numbers.
510     // PHOS internal convention: 1<module<5.
511     Int_t iDDL = 4 * (module - 1) + iRCU;
512
513     // new DDL
514     if (iDDL != prevDDL) {
515       // write real header and close previous file
516       if (buffer) {
517         buffer->Flush();
518         buffer->WriteDataHeader(kFALSE, kFALSE);
519         delete buffer;
520       }
521
522       // open new file and write dummy header
523       TString fileName = AliDAQ::DdlFileName("PHOS",iDDL);
524
525       TString path = gSystem->Getenv("ALICE_ROOT");
526       path += "/PHOS/mapping/RCU";
527       path += iRCU;
528       path += ".data";
529
530       AliAltroMapping* mapping = new AliCaloAltroMapping(path.Data());
531       buffer = new AliAltroBuffer(fileName.Data(),mapping);
532       buffer->WriteDataHeader(kTRUE, kFALSE);  //Dummy;
533
534       prevDDL = iDDL;
535     }
536
537     AliDebug(2,Form("digit E=%.4f GeV, t=%g s, (mod,col,row)=(%d,%d,%d)\n",
538                     digit->GetEnergy(),digit->GetTimeR(),
539                     relId[0]-1,relId[3]-1,relId[2]-1));
540     // if a signal is out of time range, write only trailer
541     if (digit->GetTimeR() > pulse->GetRawFormatTimeMax()*0.5 ) {
542       AliInfo("Signal is out of time range.\n");
543       buffer->FillBuffer(0);
544       buffer->FillBuffer(pulse->GetRawFormatTimeBins() );  // time bin
545       buffer->FillBuffer(3);                               // bunch length
546       buffer->WriteTrailer(3, relId[3]-1, relId[2]-1, 0);  // trailer
547       
548     // calculate the time response function
549     } else {
550       Double_t energy = 0 ;
551       module = relId[0];
552       if (digit->GetId() <= geom->GetNModules() * geom->GetNCristalsInModule()) {
553         energy=digit->GetEnergy();
554         if(energy>eMax) {eMax=energy; modMax=module; colMax=col; rowMax=row;}
555       }
556       else {
557         energy = 0; // CPV raw data format is now know yet
558       }
559       pulse->SetAmplitude(energy);
560       pulse->SetTZero(digit->GetTimeR());
561       pulse->MakeSamples();
562       pulse->GetSamples(adcValuesHigh, adcValuesLow) ; 
563       buffer->WriteChannel(relId[3]-1, relId[2]-1, 0, 
564                            pulse->GetRawFormatTimeBins(), adcValuesLow , kAdcThreshold);
565       buffer->WriteChannel(relId[3]-1, relId[2]-1, 1, 
566                            pulse->GetRawFormatTimeBins(), adcValuesHigh, kAdcThreshold);
567     }
568   }
569   
570   // write real header and close last file
571   if (buffer) {
572     buffer->Flush();
573     buffer->WriteDataHeader(kFALSE, kFALSE);
574     delete buffer;
575   }
576   
577   AliDebug(1,Form("Digit with max. energy:  modMax %d colMax %d rowMax %d  eMax %f\n",
578          modMax,colMax,rowMax,eMax));
579
580   delete pulse;
581   loader->UnloadDigits();
582 }
583
584 //____________________________________________________________________________
585 void AliPHOS::Hits2SDigits()  
586
587 // create summable digits
588
589   AliPHOSSDigitizer phosDigitizer(fLoader->GetRunLoader()->GetFileName().Data()) ;
590   phosDigitizer.SetEventRange(0, -1) ; // do all the events
591   phosDigitizer.ExecuteTask("all") ; 
592 }
593
594 //____________________________________________________________________________
595 AliLoader* AliPHOS::MakeLoader(const char* topfoldername)
596 {
597 //different behaviour than standard (singleton getter)
598 // --> to be discussed and made eventually coherent
599  fLoader = new AliPHOSLoader(GetName(),topfoldername);
600  return fLoader;
601 }
602
603 //____________________________________________________________________________
604 void AliPHOS::SetTreeAddress()
605
606   // Links Hits in the Tree to Hits array
607   TBranch *branch;
608   char branchname[20];
609   sprintf(branchname,"%s",GetName());
610   // Branch address for hit tree
611     TTree *treeH = TreeH();
612   if (treeH) {
613     branch = treeH->GetBranch(branchname);
614     if (branch) 
615      { 
616        if (fHits == 0x0) fHits= new TClonesArray("AliPHOSHit",1000);
617        //AliInfo(Form("<%s> Setting Hits Address",GetName()));
618        branch->SetAddress(&fHits);
619      }
620   }
621 }
622