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