1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 /* History of cvs commits:
19 * Revision 1.111 2007/02/18 15:21:47 kharlov
20 * Corrections for memory leak in Digits2Raw due to AliAltroMapping
22 * Revision 1.110 2007/02/13 10:52:08 policheh
23 * Raw2SDigits() implemented
25 * Revision 1.109 2007/02/05 10:43:25 hristov
26 * Changes for correct initialization of Geant4 (Mihaela)
28 * Revision 1.108 2007/02/01 10:34:47 hristov
29 * Removing warnings on Solaris x86
31 * Revision 1.107 2007/01/29 16:29:37 kharlov
32 * Digits2Raw(): special workaround for digits with time out of range
34 * Revision 1.106 2007/01/17 17:28:56 kharlov
35 * Extract ALTRO sample generation to a separate class AliPHOSPulseGenerator
37 * Revision 1.105 2007/01/12 21:44:29 kharlov
38 * Simulate and reconstruct two gains simulaneouslsy
40 * Revision 1.104 2006/11/23 13:40:44 hristov
41 * Common class for raw data reading and ALTRO mappiing for PHOS and EMCAL (Gustavo, Cvetan)
43 * Revision 1.103 2006/11/14 17:11:15 hristov
44 * 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
46 * Revision 1.102 2006/10/27 17:14:27 kharlov
47 * Introduce AliDebug and AliLog (B.Polichtchouk)
49 * Revision 1.101 2006/10/13 06:47:29 kharlov
50 * Simulation of RAW data applies real mapping (B.Polichtchouk)
52 * Revision 1.100 2006/08/11 12:36:26 cvetan
53 * Update of the PHOS code needed in order to read and reconstruct the beam test raw data (i.e. without an existing galice.root)
55 * Revision 1.99 2006/06/28 11:36:09 cvetan
56 * 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.
58 * Revision 1.98 2006/05/11 11:30:48 cvetan
59 * 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
61 * Revision 1.97 2006/04/22 10:30:17 hristov
62 * Add fEnergy to AliPHOSDigit and operate with EMC amplitude in energy units (Yu.Kharlov)
64 * Revision 1.96 2006/04/07 08:41:59 hristov
65 * Follow AliAlignObj framework and remove AliPHOSAlignData (Yu.Kharlov)
67 * Revision 1.95 2006/03/14 19:40:41 kharlov
68 * Remove De-digitizing of raw data and digitizing the raw data fit
70 * Revision 1.94 2006/03/07 18:56:25 kharlov
71 * CDB is passed via environment variable
73 * Revision 1.93 2005/11/22 08:45:11 kharlov
74 * Calibration is read from CDB if any (Boris Polichtchouk)
76 * Revision 1.92 2005/11/03 13:09:19 hristov
77 * Removing meaningless const declarations (linuxicc)
79 * Revision 1.91 2005/07/27 15:08:53 kharlov
80 * Mixture ArCO2 is corrected
82 * Revision 1.90 2005/06/17 07:39:07 hristov
83 * Removing GetDebug and SetDebug from AliRun and AliModule. Using AliLog for the messages
85 * Revision 1.89 2005/05/28 12:10:07 schutz
86 * Copy constructor is corrected (by T.P.)
90 //_________________________________________________________________________
91 // Base Class for PHOS description:
92 // PHOS consists of a PbWO4 calorimeter (EMCA) and a gazeous charged
93 // particles detector (CPV or PPSD).
94 // The only provided method here is CreateMaterials,
95 // which defines the materials common to all PHOS versions.
97 //*-- Author: Laurent Aphecetche & Yves Schutz (SUBATECH)
98 //////////////////////////////////////////////////////////////////////////////
101 // --- ROOT system ---
105 #include <TVirtualMC.h>
110 // --- Standard library ---
112 // --- AliRoot header files ---
115 #include "AliPHOSGetter.h"
117 #include "AliPHOSDigitizer.h"
118 #include "AliPHOSSDigitizer.h"
119 #include "AliPHOSDigit.h"
120 #include "AliAltroBuffer.h"
121 #include "AliAltroMapping.h"
122 #include "AliCaloAltroMapping.h"
124 #include "AliCDBManager.h"
125 #include "AliCDBEntry.h"
126 #include "AliCDBStorage.h"
127 #include "AliPHOSCalibData.h"
128 #include "AliPHOSPulseGenerator.h"
130 #include "AliPHOSRawDecoder.h"
131 #include "AliPHOSRawDigiProducer.h"
135 //____________________________________________________________________________
136 AliPHOS:: AliPHOS() : AliDetector()
143 //____________________________________________________________________________
144 AliPHOS::AliPHOS(const char* name, const char* title): AliDetector(name, title)
146 // ctor : title is used to identify the layout
149 //____________________________________________________________________________
154 //____________________________________________________________________________
155 AliDigitizer* AliPHOS::CreateDigitizer(AliRunDigitizer* manager) const
157 return new AliPHOSDigitizer(manager);
160 //____________________________________________________________________________
161 void AliPHOS::CreateMaterials()
163 // Definitions of materials to build PHOS and associated tracking media.
164 // media number in idtmed are 699 to 798.
166 // --- The PbWO4 crystals ---
167 Float_t aX[3] = {207.19, 183.85, 16.0} ;
168 Float_t zX[3] = {82.0, 74.0, 8.0} ;
169 Float_t wX[3] = {1.0, 1.0, 4.0} ;
172 AliMixture(0, "PbWO4$", aX, zX, dX, -3, wX) ;
175 // --- The polysterene scintillator (CH) ---
176 Float_t aP[2] = {12.011, 1.00794} ;
177 Float_t zP[2] = {6.0, 1.0} ;
178 Float_t wP[2] = {1.0, 1.0} ;
181 AliMixture(1, "Polystyrene$", aP, zP, dP, -2, wP) ;
184 AliMaterial(2, "Al$", 26.98, 13., 2.7, 8.9, 999., 0, 0) ;
185 // --- Absorption length is ignored ^
187 // --- Tyvek (CnH2n) ---
188 Float_t aT[2] = {12.011, 1.00794} ;
189 Float_t zT[2] = {6.0, 1.0} ;
190 Float_t wT[2] = {1.0, 2.0} ;
193 AliMixture(3, "Tyvek$", aT, zT, dT, -2, wT) ;
195 // --- Polystyrene foam ---
196 Float_t aF[2] = {12.011, 1.00794} ;
197 Float_t zF[2] = {6.0, 1.0} ;
198 Float_t wF[2] = {1.0, 1.0} ;
201 AliMixture(4, "Foam$", aF, zF, dF, -2, wF) ;
204 Float_t aTIT[3] = {47.88, 26.98, 54.94} ;
205 Float_t zTIT[3] = {22.0, 13.0, 25.0} ;
206 Float_t wTIT[3] = {69.0, 6.0, 1.0} ;
209 AliMixture(5, "Titanium$", aTIT, zTIT, dTIT, -3, wTIT);
212 AliMaterial(6, "Si$", 28.0855, 14., 2.33, 9.36, 42.3, 0, 0) ;
216 // --- Foam thermo insulation ---
217 Float_t aTI[2] = {12.011, 1.00794} ;
218 Float_t zTI[2] = {6.0, 1.0} ;
219 Float_t wTI[2] = {1.0, 1.0} ;
222 AliMixture(7, "Thermo Insul.$", aTI, zTI, dTI, -2, wTI) ;
225 Float_t aTX[4] = {16.0, 28.09, 12.011, 1.00794} ;
226 Float_t zTX[4] = {8.0, 14.0, 6.0, 1.0} ;
227 Float_t wTX[4] = {292.0, 68.0, 462.0, 736.0} ;
230 AliMixture(8, "Textolit$", aTX, zTX, dTX, -4, wTX) ;
233 Float_t aFR[4] = {16.0, 28.09, 12.011, 1.00794} ;
234 Float_t zFR[4] = {8.0, 14.0, 6.0, 1.0} ;
235 Float_t wFR[4] = {292.0, 68.0, 462.0, 736.0} ;
238 AliMixture(9, "FR4$", aFR, zFR, dFR, -4, wFR) ;
240 // --- The Composite Material for micromegas (so far polyetylene) ---
241 Float_t aCM[2] = {12.01, 1.} ;
242 Float_t zCM[2] = {6., 1.} ;
243 Float_t wCM[2] = {1., 2.} ;
244 Float_t dCM = 0.935 ;
246 AliMixture(10, "Compo Mat$", aCM, zCM, dCM, -2, wCM) ;
249 AliMaterial(11, "Cu$", 63.546, 29, 8.96, 1.43, 14.8, 0, 0) ;
251 // --- G10 : Printed Circuit material ---
252 Float_t aG10[4] = { 12., 1., 16., 28.} ;
253 Float_t zG10[4] = { 6., 1., 8., 14.} ;
254 Float_t wG10[4] = { .259, .288, .248, .205} ;
257 AliMixture(12, "G10$", aG10, zG10, dG10, -4, wG10);
260 AliMaterial(13, "Pb$", 207.2, 82, 11.35, 0.56, 0., 0, 0) ;
262 // --- The gas mixture ---
264 Float_t aCO[2] = {12.0, 16.0} ;
265 Float_t zCO[2] = {6.0, 8.0} ;
266 Float_t wCO[2] = {1.0, 2.0} ;
267 Float_t dCO = 0.001977 ;
269 AliMixture(14, "CO2$", aCO, zCO, dCO, -2, wCO);
272 Float_t dAr = 0.001782 ;
273 AliMaterial(15, "Ar$", 39.948, 18.0, dAr, 14.0, 0., 0, 0) ;
275 // Ar+CO2 Mixture (80% / 20%)
276 Float_t arContent = 0.80 ; // Ar-content of the ArCO2-mixture
277 Float_t aArCO[3] = {39.948, 12.0, 16.0} ;
278 Float_t zArCO[3] = {18.0 , 6.0, 8.0} ;
280 wArCO[0] = arContent;
281 wArCO[1] = (1-arContent)*1;
282 wArCO[2] = (1-arContent)*2;
283 Float_t dArCO = arContent*dAr + (1-arContent)*dCO ;
284 AliMixture(16, "ArCO2$", aArCO, zArCO, dArCO, -3, wArCO) ;
286 // --- Stainless steel (let it be pure iron) ---
287 AliMaterial(17, "Steel$", 55.845, 26, 7.87, 1.76, 0., 0, 0) ;
290 // --- Fiberglass ---
291 Float_t aFG[4] = {16.0, 28.09, 12.011, 1.00794} ;
292 Float_t zFG[4] = {8.0, 14.0, 6.0, 1.0} ;
293 Float_t wFG[4] = {292.0, 68.0, 462.0, 736.0} ;
296 AliMixture(18, "Fibergla$", aFG, zFG, dFG, -4, wFG) ;
298 // --- Cables in Air box ---
301 Float_t aCA[4] = { 1.,12.,55.8,63.5 };
302 Float_t zCA[4] = { 1.,6.,26.,29. };
303 Float_t wCA[4] = { .014,.086,.42,.48 };
304 Float_t dCA = 0.8 ; //this density is raw estimation, if you know better - correct
306 AliMixture(19, "Cables $", aCA, zCA, dCA, -4, wCA) ;
310 Float_t aAir[4]={12.0107,14.0067,15.9994,39.948};
311 Float_t zAir[4]={6.,7.,8.,18.};
312 Float_t wAir[4]={0.000124,0.755267,0.231781,0.012827};
313 Float_t dAir = 1.20479E-3;
315 AliMixture(99, "Air$", aAir, zAir, dAir, 4, wAir) ;
317 // DEFINITION OF THE TRACKING MEDIA
319 // for PHOS: idtmed[699->798] equivalent to fIdtmed[0->100]
320 Int_t isxfld = gAlice->Field()->Integ() ;
321 Float_t sxmgmx = gAlice->Field()->Max() ;
323 // The scintillator of the calorimeter made of PBW04 -> idtmed[699]
324 AliMedium(0, "PHOS Xtal $", 0, 1,
325 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
327 // The scintillator of the CPV made of Polystyrene scintillator -> idtmed[700]
328 AliMedium(1, "CPV scint. $", 1, 1,
329 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
331 // Various Aluminium parts made of Al -> idtmed[701]
332 AliMedium(2, "Al parts $", 2, 0,
333 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.001, 0.001, 0, 0) ;
335 // The Tywek which wraps the calorimeter crystals -> idtmed[702]
336 AliMedium(3, "Tyvek wrapper$", 3, 0,
337 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.001, 0.001, 0, 0) ;
339 // The Polystyrene foam around the calorimeter module -> idtmed[703]
340 AliMedium(4, "Polyst. foam $", 4, 0,
341 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
343 // The Titanium around the calorimeter crystal -> idtmed[704]
344 AliMedium(5, "Titan. cover $", 5, 0,
345 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.0001, 0.0001, 0, 0) ;
347 // The Silicon of the pin diode to read out the calorimeter crystal -> idtmed[705]
348 AliMedium(6, "Si PIN $", 6, 0,
349 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.01, 0.01, 0, 0) ;
351 // The thermo insulating material of the box which contains the calorimeter module -> idtmed[706]
352 AliMedium(7, "Thermo Insul.$", 7, 0,
353 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
355 // The Textolit which makes up the box which contains the calorimeter module -> idtmed[707]
356 AliMedium(8, "Textolit $", 8, 0,
357 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
359 // FR4: The Plastic which makes up the frame of micromegas -> idtmed[708]
360 AliMedium(9, "FR4 $", 9, 0,
361 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.0001, 0, 0) ;
364 // The Composite Material for micromegas -> idtmed[709]
365 AliMedium(10, "CompoMat $", 10, 0,
366 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
368 // Copper -> idtmed[710]
369 AliMedium(11, "Copper $", 11, 0,
370 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.0001, 0, 0) ;
372 // G10: Printed Circuit material -> idtmed[711]
374 AliMedium(12, "G10 $", 12, 0,
375 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.01, 0, 0) ;
377 // The Lead -> idtmed[712]
379 AliMedium(13, "Lead $", 13, 0,
380 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
382 // The gas mixture: ArCo2 -> idtmed[715]
384 AliMedium(16, "ArCo2 $", 16, 1,
385 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.01, 0, 0) ;
387 // Stainless steel -> idtmed[716]
388 AliMedium(17, "Steel $", 17, 0,
389 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.0001, 0, 0) ;
391 // Fibergalss -> idtmed[717]
392 AliMedium(18, "Fiberglass$", 18, 0,
393 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
395 // Cables in air -> idtmed[718]
396 AliMedium(19, "Cables $", 19, 0,
397 isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
399 // Air -> idtmed[798]
400 AliMedium(99, "Air $", 99, 0,
401 isxfld, sxmgmx, 10.0, 1.0, 0.1, 0.1, 10.0, 0, 0) ;
404 //_____________________________________________________________________________
408 // Initialises cuts for PHOS
410 // --- Set decent energy thresholds for gamma and electron tracking
411 Int_t * idtmed = fIdtmed->GetArray() - 699 ;
413 // Tracking threshold for photons and electrons in the scintillator crystal
414 gMC->Gstpar(idtmed[699], "CUTGAM",0.5E-4) ;
415 gMC->Gstpar(idtmed[699], "CUTELE",1.0E-4) ;
417 // --- Generate explicitly delta rays in the titan cover ---
418 gMC->Gstpar(idtmed[704], "LOSS",3.) ;
419 gMC->Gstpar(idtmed[704], "DRAY",1.) ;
420 // --- and in aluminium parts ---
421 gMC->Gstpar(idtmed[701], "LOSS",3.) ;
422 gMC->Gstpar(idtmed[701], "DRAY",1.) ;
423 // --- and in PIN diode
424 gMC->Gstpar(idtmed[705], "LOSS",3) ;
425 gMC->Gstpar(idtmed[705], "DRAY",1) ;
426 // --- and in the passive convertor
427 gMC->Gstpar(idtmed[712], "LOSS",3) ;
428 gMC->Gstpar(idtmed[712], "DRAY",1) ;
429 // Tracking threshold for photons and electrons in the gas ArC02
430 gMC->Gstpar(idtmed[715], "CUTGAM",1.E-5) ;
431 gMC->Gstpar(idtmed[715], "CUTELE",1.E-5) ;
432 gMC->Gstpar(idtmed[715], "CUTNEU",1.E-5) ;
433 gMC->Gstpar(idtmed[715], "CUTHAD",1.E-5) ;
434 gMC->Gstpar(idtmed[715], "CUTMUO",1.E-5) ;
435 gMC->Gstpar(idtmed[715], "BCUTE",1.E-5) ;
436 gMC->Gstpar(idtmed[715], "BCUTM",1.E-5) ;
437 gMC->Gstpar(idtmed[715], "DCUTE",1.E-5) ;
438 gMC->Gstpar(idtmed[715], "DCUTM",1.E-5) ;
439 gMC->Gstpar(idtmed[715], "PPCUTM",1.E-5) ;
440 gMC->Gstpar(idtmed[715], "LOSS",2.) ;
441 gMC->Gstpar(idtmed[715], "DRAY",0.) ;
442 gMC->Gstpar(idtmed[715], "STRA",2.) ;
445 //____________________________________________________________________________
446 void AliPHOS::Digits2Raw()
448 // convert digits of the current event to raw data
450 AliPHOSLoader * loader = dynamic_cast<AliPHOSLoader*>(fLoader) ;
453 loader->LoadDigits();
454 TClonesArray* digits = loader->Digits() ;
457 AliError(Form("No digits found !"));
462 AliPHOSGeometry* geom = GetGeometry();
464 AliError(Form("No geometry found !"));
468 // some digitization constants
469 const Float_t kThreshold = 0.001; // skip digits below 1 MeV
470 const Int_t kAdcThreshold = 1; // Lower ADC threshold to write to raw data
474 // Create a shaper pulse object
475 AliPHOSPulseGenerator *pulse = new AliPHOSPulseGenerator();
477 Int_t *adcValuesLow = new Int_t[pulse->GetRawFormatTimeBins()];
478 Int_t *adcValuesHigh= new Int_t[pulse->GetRawFormatTimeBins()];
480 const Int_t maxDDL = 20;
481 AliAltroBuffer *buffer[maxDDL];
482 AliAltroMapping *mapping[maxDDL];
484 for(Int_t jDDL=0; jDDL<maxDDL; jDDL++) {
496 // loop over digits (assume ordered digits)
497 for (Int_t iDigit = 0; iDigit < digits->GetEntries(); iDigit++) {
498 AliPHOSDigit* digit = dynamic_cast<AliPHOSDigit *>(digits->At(iDigit)) ;
499 if (digit->GetEnergy() < kThreshold)
502 geom->AbsToRelNumbering(digit->GetId(), relId);
503 Int_t module = relId[0];
507 continue; // ignore digits from CPV
510 Int_t row = relId[2]-1;
511 Int_t col = relId[3]-1;
516 if(0<=row&&row<32 && 0<=col&&col<28) iRCU=0;
519 if(0<=row&&row<32 && 28<=col&&col<56) iRCU=1;
522 if(32<=row&&row<64 && 0<=col&&col<28) iRCU=2;
525 if(32<=row&&row<64 && 28<=col&&col<56) iRCU=3;
528 // PHOS EMCA has 4 DDL per module. Splitting is based on the (row,column) numbers.
529 // PHOS internal convention: 1<module<5.
530 Int_t iDDL = 4 * (module - 1) + iRCU;
533 if (iDDL != prevDDL) {
534 if (buffer[iDDL] == 0) {
535 // open new file and write dummy header
536 TString fileName = AliDAQ::DdlFileName("PHOS",iDDL);
538 TString path = gSystem->Getenv("ALICE_ROOT");
539 path += "/PHOS/mapping/RCU";
543 mapping[iDDL] = new AliCaloAltroMapping(path.Data());
544 buffer[iDDL] = new AliAltroBuffer(fileName.Data(),mapping[iDDL]);
545 buffer[iDDL]->WriteDataHeader(kTRUE, kFALSE); //Dummy;
550 AliDebug(2,Form("digit E=%.4f GeV, t=%g s, (mod,col,row)=(%d,%d,%d)\n",
551 digit->GetEnergy(),digit->GetTimeR(),
552 relId[0]-1,relId[3]-1,relId[2]-1));
553 // if a signal is out of time range, write only trailer
554 if (digit->GetTimeR() > pulse->GetRawFormatTimeMax()*0.5 ) {
555 AliInfo("Signal is out of time range.\n");
556 buffer[iDDL]->FillBuffer(0);
557 buffer[iDDL]->FillBuffer(pulse->GetRawFormatTimeBins() ); // time bin
558 buffer[iDDL]->FillBuffer(3); // bunch length
559 buffer[iDDL]->WriteTrailer(3, relId[3]-1, relId[2]-1, 0); // trailer
561 // calculate the time response function
563 Double_t energy = 0 ;
565 if (digit->GetId() <= geom->GetNModules() * geom->GetNCristalsInModule()) {
566 energy=digit->GetEnergy();
567 if(energy>eMax) {eMax=energy; modMax=module; colMax=col; rowMax=row;}
570 energy = 0; // CPV raw data format is now know yet
572 pulse->SetAmplitude(energy);
573 pulse->SetTZero(digit->GetTimeR());
574 pulse->MakeSamples();
575 pulse->GetSamples(adcValuesHigh, adcValuesLow) ;
576 buffer[iDDL]->WriteChannel(relId[3]-1, relId[2]-1, 0,
577 pulse->GetRawFormatTimeBins(), adcValuesLow , kAdcThreshold);
578 buffer[iDDL]->WriteChannel(relId[3]-1, relId[2]-1, 1,
579 pulse->GetRawFormatTimeBins(), adcValuesHigh, kAdcThreshold);
583 // write real header and close last file
584 for (Int_t iDDL=0; iDDL<maxDDL; iDDL++) {
586 buffer[iDDL]->Flush();
587 buffer[iDDL]->WriteDataHeader(kFALSE, kFALSE);
589 if (mapping[iDDL]) delete mapping[iDDL];
593 AliDebug(1,Form("Digit with max. energy: modMax %d colMax %d rowMax %d eMax %f\n",
594 modMax,colMax,rowMax,eMax));
597 loader->UnloadDigits();
600 //____________________________________________________________________________
601 void AliPHOS::Hits2SDigits()
603 // create summable digits
605 AliPHOSSDigitizer phosDigitizer(fLoader->GetRunLoader()->GetFileName().Data()) ;
606 phosDigitizer.SetEventRange(0, -1) ; // do all the events
607 phosDigitizer.ExecuteTask("all") ;
610 //____________________________________________________________________________
611 AliLoader* AliPHOS::MakeLoader(const char* topfoldername)
613 //different behaviour than standard (singleton getter)
614 // --> to be discussed and made eventually coherent
615 fLoader = new AliPHOSLoader(GetName(),topfoldername);
619 //____________________________________________________________________________
620 void AliPHOS::SetTreeAddress()
622 // Links Hits in the Tree to Hits array
625 sprintf(branchname,"%s",GetName());
626 // Branch address for hit tree
627 TTree *treeH = TreeH();
629 branch = treeH->GetBranch(branchname);
632 if (fHits == 0x0) fHits= new TClonesArray("AliPHOSHit",1000);
633 //AliInfo(Form("<%s> Setting Hits Address",GetName()));
634 branch->SetAddress(&fHits);
639 //____________________________________________________________________________
640 Bool_t AliPHOS::Raw2SDigits(AliRawReader* rawReader)
643 AliPHOSLoader * loader = dynamic_cast<AliPHOSLoader*>(fLoader) ;
646 tree = loader->TreeS() ;
648 loader->MakeTree("S");
649 tree = loader->TreeS() ;
652 TClonesArray * sdigits = loader->SDigits() ;
654 loader->MakeSDigitsArray();
655 sdigits = loader->SDigits();
659 AliPHOSRawDecoder dc(rawReader);
660 AliPHOSRawDigiProducer pr;
661 pr.MakeDigits(sdigits,&dc);
663 Int_t bufferSize = 32000 ;
664 TBranch * sdigitsBranch = tree->Branch("PHOS",&sdigits,bufferSize);
667 fLoader->WriteSDigits("OVERWRITE");