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