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