updating apd and fee parameters
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALDigitizer.cxx
CommitLineData
61e0abb5 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
16/* $Id$ */
17
18//_________________________________________________________________________
ffa6d63b 19//
61e0abb5 20//////////////////////////////////////////////////////////////////////////////
ffa6d63b 21// Class performs digitization of Summable digits
22//
61e0abb5 23// In addition it performs mixing of summable digits from different events.
24//
25// For each event two branches are created in TreeD:
26// "EMCAL" - list of digits
27// "AliEMCALDigitizer" - AliEMCALDigitizer with all parameters used in digitization
28//
29// Note, that one cset title for new digits branch, and repeat digitization with
30// another set of parameters.
31//
32// Examples of use:
33// root[0] AliEMCALDigitizer * d = new AliEMCALDigitizer() ;
34// root[1] d->ExecuteTask()
35// Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
36// //Digitizes SDigitis in all events found in file galice.root
37//
38// root[2] AliEMCALDigitizer * d1 = new AliEMCALDigitizer("galice1.root") ;
39// // Will read sdigits from galice1.root
40// root[3] d1->MixWith("galice2.root")
41// Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
42// // Reads another portion of sdigits from galice2.root
43// root[3] d1->MixWith("galice3.root")
44// // Reads another portion of sdigits from galice3.root
45// root[4] d->ExecuteTask("deb timing")
46// // Reads SDigits from files galice1.root, galice2.root ....
47// // mixes them and stores produced Digits in file galice1.root
48// // deb - prints number of produced digits
49// // deb all - prints list of produced digits
50// // timing - prints time used for digitization
ffa6d63b 51////////////////////////////////////////////////////////////////////////////////////
61e0abb5 52//
ffa6d63b 53//*-- Author: Sahal Yacoob (LBL)
814ad4bf 54// based on : AliEMCALDigitizer
05a92d59 55// Modif:
56// August 2002 Yves Schutz: clone PHOS as closely as possible and intoduction
57// of new IO (à la PHOS)
1963b290 58// November 2003 Aleksei Pavlinov : adopted for Shish-Kebab geometry
59//_________________________________________________________________________________
61e0abb5 60
61// --- ROOT system ---
e939a978 62#include <TROOT.h>
63#include <TTree.h>
64#include <TSystem.h>
65#include <TBenchmark.h>
66#include <TList.h>
67#include <TH1.h>
68#include <TBrowser.h>
69#include <TObjectTable.h>
70#include <TRandom.h>
21e7df44 71#include <cassert>
61e0abb5 72
73// --- AliRoot header files ---
b42d4197 74#include "AliLog.h"
e4dcc484 75#include "AliRun.h"
ea5f3389 76#include "AliRunDigitizer.h"
5dee926e 77#include "AliRunLoader.h"
53f1c378 78#include "AliCDBManager.h"
79#include "AliCDBEntry.h"
61e0abb5 80#include "AliEMCALDigit.h"
05a92d59 81#include "AliEMCAL.h"
5dee926e 82#include "AliEMCALLoader.h"
61e0abb5 83#include "AliEMCALDigitizer.h"
84#include "AliEMCALSDigitizer.h"
85#include "AliEMCALGeometry.h"
05a92d59 86#include "AliEMCALTick.h"
315d1c64 87#include "AliEMCALHistoUtilities.h"
05a92d59 88
61e0abb5 89ClassImp(AliEMCALDigitizer)
90
91
92//____________________________________________________________________________
18a21c7c 93AliEMCALDigitizer::AliEMCALDigitizer()
94 : AliDigitizer("",""),
95 fDefaultInit(kTRUE),
96 fDigitsInRun(0),
97 fInit(0),
98 fInput(0),
99 fInputFileNames(0x0),
100 fEventNames(0x0),
101 fDigitThreshold(0),
102 fMeanPhotonElectron(0),
103 fPedestal(0),
104 fSlope(0),
105 fPinNoise(0),
106 fTimeResolution(0),
107 fTimeThreshold(0),
108 fTimeSignalLength(0),
109 fADCchannelEC(0),
110 fADCpedestalEC(0),
111 fNADCEC(0),
112 fEventFolderName(""),
113 fFirstEvent(0),
114 fLastEvent(0),
115 fControlHists(0),
53f1c378 116 fHists(0),fCalibData(0x0)
61e0abb5 117{
118 // ctor
88cb7938 119 InitParameters() ;
05a92d59 120 fManager = 0 ; // We work in the standalong mode
839828a6 121}
61e0abb5 122
839828a6 123//____________________________________________________________________________
18a21c7c 124AliEMCALDigitizer::AliEMCALDigitizer(TString alirunFileName, TString eventFolderName)
125 : AliDigitizer("EMCAL"+AliConfig::Instance()->GetDigitizerTaskName(), alirunFileName),
126 fDefaultInit(kFALSE),
127 fDigitsInRun(0),
128 fInit(0),
129 fInput(0),
130 fInputFileNames(0),
131 fEventNames(0),
132 fDigitThreshold(0),
133 fMeanPhotonElectron(0),
134 fPedestal(0),
135 fSlope(0),
136 fPinNoise(0),
137 fTimeResolution(0),
138 fTimeThreshold(0),
139 fTimeSignalLength(0),
140 fADCchannelEC(0),
141 fADCpedestalEC(0),
142 fNADCEC(0),
143 fEventFolderName(eventFolderName),
144 fFirstEvent(0),
145 fLastEvent(0),
146 fControlHists(0),
49f91537 147 fHists(0),fCalibData(0x0)
839828a6 148{
05a92d59 149 // ctor
a9485c37 150 InitParameters() ;
88cb7938 151 Init() ;
88cb7938 152 fManager = 0 ; // We work in the standalong mode
61e0abb5 153}
839828a6 154
61e0abb5 155//____________________________________________________________________________
18a21c7c 156AliEMCALDigitizer::AliEMCALDigitizer(const AliEMCALDigitizer & d)
157 : AliDigitizer(d.GetName(),d.GetTitle()),
158 fDefaultInit(d.fDefaultInit),
159 fDigitsInRun(d.fDigitsInRun),
160 fInit(d.fInit),
161 fInput(d.fInput),
162 fInputFileNames(d.fInputFileNames),
163 fEventNames(d.fEventNames),
164 fDigitThreshold(d.fDigitThreshold),
165 fMeanPhotonElectron(d.fMeanPhotonElectron),
166 fPedestal(d.fPedestal),
167 fSlope(d.fSlope),
168 fPinNoise(d.fPinNoise),
169 fTimeResolution(d.fTimeResolution),
170 fTimeThreshold(d.fTimeThreshold),
171 fTimeSignalLength(d.fTimeSignalLength),
172 fADCchannelEC(d.fADCchannelEC),
173 fADCpedestalEC(d.fADCpedestalEC),
174 fNADCEC(d.fNADCEC),
175 fEventFolderName(d.fEventFolderName),
176 fFirstEvent(d.fFirstEvent),
177 fLastEvent(d.fLastEvent),
178 fControlHists(d.fControlHists),
53f1c378 179 fHists(d.fHists),fCalibData(d.fCalibData)
88cb7938 180{
181 // copyy ctor
88cb7938 182 }
183
184//____________________________________________________________________________
18a21c7c 185AliEMCALDigitizer::AliEMCALDigitizer(AliRunDigitizer * rd)
186 : AliDigitizer(rd,"EMCAL"+AliConfig::Instance()->GetDigitizerTaskName()),
187 fDefaultInit(kFALSE),
188 fDigitsInRun(0),
189 fInit(0),
190 fInput(0),
191 fInputFileNames(0),
192 fEventNames(0),
193 fDigitThreshold(0.),
194 fMeanPhotonElectron(0),
195 fPedestal(0),
196 fSlope(0.),
197 fPinNoise(0),
198 fTimeResolution(0.),
199 fTimeThreshold(0),
200 fTimeSignalLength(0),
201 fADCchannelEC(0),
202 fADCpedestalEC(0),
203 fNADCEC(0),
204 fEventFolderName(0),
205 fFirstEvent(0),
206 fLastEvent(0),
207 fControlHists(0),
49f91537 208 fHists(0),fCalibData(0x0)
61e0abb5 209{
45fa49ca 210 // ctor Init() is called by RunDigitizer
88cb7938 211 fManager = rd ;
e0e18b6d 212 fEventFolderName = fManager->GetInputFolderName(0) ;
88cb7938 213 SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
05a92d59 214 InitParameters() ;
839828a6 215}
814ad4bf 216
839828a6 217//____________________________________________________________________________
218 AliEMCALDigitizer::~AliEMCALDigitizer()
219{
14ce0a6e 220 //dtor
5dee926e 221 if (AliRunLoader::GetRunLoader()) {
222 AliLoader *emcalLoader=0;
223 if ((emcalLoader = AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL")))
224 emcalLoader->CleanDigitizer();
225 }
226 else
227 AliDebug(1," no runloader present");
f5f384f4 228 delete [] fInputFileNames ;
229 delete [] fEventNames ;
64942713 230
231 if(fHists) delete fHists;
61e0abb5 232}
233
234//____________________________________________________________________________
09884213 235void AliEMCALDigitizer::Digitize(Int_t event)
05a92d59 236{
61e0abb5 237
238 // Makes the digitization of the collected summable digits
239 // for this it first creates the array of all EMCAL modules
240 // filled with noise (different for EMC, CPV and PPSD) and
241 // after that adds contributions from SDigits. This design
242 // helps to avoid scanning over the list of digits to add
243 // contribution of any new SDigit.
1963b290 244 static int nEMC=0; //max number of digits possible
61e0abb5 245
5dee926e 246 AliRunLoader *rl = AliRunLoader::GetRunLoader();
247 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(rl->GetDetectorLoader("EMCAL"));
14ce0a6e 248 Int_t readEvent = event ;
1963b290 249 // fManager is data member from AliDigitizer
45fa49ca 250 if (fManager)
14ce0a6e 251 readEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetCurrentEventNumber() ;
cf24bdc4 252 AliDebug(1,Form("Adding event %d from input stream 0 %s %s",
14ce0a6e 253 readEvent, GetTitle(), fEventFolderName.Data())) ;
254 rl->GetEvent(readEvent);
5dee926e 255
256 TClonesArray * digits = emcalLoader->Digits() ;
30bc5594 257 digits->Delete() ;
61e0abb5 258
e4dcc484 259 // Load Geometry
30bc5594 260 AliEMCALGeometry *geom = 0;
261 if (rl->GetAliRun()) {
262 AliEMCAL * emcal = (AliEMCAL*)rl->GetAliRun()->GetDetector("EMCAL");
263 geom = emcal->GetGeometry();
264 }
265 else
266 AliFatal("Could not get AliRun from runLoader");
e4dcc484 267
fe93d365 268 nEMC = geom->GetNCells();
269 AliDebug(1,Form("nEMC %i (number cells in EMCAL) | %s \n", nEMC, geom->GetName()));
270
88cb7938 271 Int_t absID ;
61e0abb5 272
88cb7938 273 digits->Expand(nEMC) ;
05a92d59 274
88cb7938 275 // get first the sdigitizer from the tasks list (must have same name as the digitizer)
5dee926e 276 if ( !emcalLoader->SDigitizer() )
277 emcalLoader->LoadSDigitizer();
278 AliEMCALSDigitizer * sDigitizer = dynamic_cast<AliEMCALSDigitizer*>(emcalLoader->SDigitizer());
88cb7938 279
280 if ( !sDigitizer )
281 Fatal("Digitize", "SDigitizer with name %s %s not found", fEventFolderName.Data(), GetTitle() ) ;
282
283 //take all the inputs to add together and load the SDigits
284 TObjArray * sdigArray = new TObjArray(fInput) ;
5dee926e 285 sdigArray->AddAt(emcalLoader->SDigits(), 0) ;
88cb7938 286 Int_t i ;
5dee926e 287
88cb7938 288 for(i = 1 ; i < fInput ; i++){
289 TString tempo(fEventNames[i]) ;
290 tempo += i ;
3880b0fb 291
292 AliRunLoader * rl2 = AliRunLoader::GetRunLoader(tempo) ;
293
294 if (rl2==0)
295 rl2 = AliRunLoader::Open(fInputFileNames[i], tempo) ;
296
45fa49ca 297 if (fManager)
14ce0a6e 298 readEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(i))->GetCurrentEventNumber() ;
299 Info("Digitize", "Adding event %d from input stream %d %s %s", readEvent, i, fInputFileNames[i].Data(), tempo.Data()) ;
3880b0fb 300 rl2->LoadSDigits();
14ce0a6e 301 rl2->GetEvent(readEvent);
3880b0fb 302 AliEMCALLoader *emcalLoader2 = dynamic_cast<AliEMCALLoader*>(rl2->GetDetectorLoader("EMCAL"));
303 sdigArray->AddAt(emcalLoader2->SDigits(), i) ;
0c219753 304 }
5dee926e 305
d968cee0 306 //Find the first tower with signal
a9485c37 307 Int_t nextSig = nEMC + 1 ;
88cb7938 308 TClonesArray * sdigits ;
309 for(i = 0 ; i < fInput ; i++){
310 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
814ad4bf 311 if ( !sdigits->GetEntriesFast() )
312 continue ;
88cb7938 313 Int_t curNext = dynamic_cast<AliEMCALDigit *>(sdigits->At(0))->GetId() ;
814ad4bf 314 if(curNext < nextSig)
315 nextSig = curNext ;
cf24bdc4 316 AliDebug(1,Form("input %i : #sdigits %i \n",
317 i, sdigits->GetEntriesFast()));
814ad4bf 318 }
cf24bdc4 319 AliDebug(1,Form("FIRST tower with signal %i \n", nextSig));
61e0abb5 320
88cb7938 321 TArrayI index(fInput) ;
814ad4bf 322 index.Reset() ; //Set all indexes to zero
61e0abb5 323
05a92d59 324 AliEMCALDigit * digit ;
325 AliEMCALDigit * curSDigit ;
ffa6d63b 326
814ad4bf 327 TClonesArray * ticks = new TClonesArray("AliEMCALTick",1000) ;
ffa6d63b 328
5dee926e 329 //Put Noise contribution
d25f2c54 330 for(absID = 0; absID < nEMC; absID++){ // Nov 30, 2006 by PAI; was from 1 to nEMC
88cb7938 331 Float_t amp = 0 ;
a9485c37 332 // amplitude set to zero, noise will be added later
d25f2c54 333 new((*digits)[absID]) AliEMCALDigit( -1, -1, absID, 0, TimeOfNoise() ); // absID-1->absID
814ad4bf 334 //look if we have to add signal?
d25f2c54 335 digit = dynamic_cast<AliEMCALDigit *>(digits->At(absID)); // absID-1->absID
88cb7938 336
814ad4bf 337 if(absID==nextSig){
54cbf5a6 338 //Add SDigits from all inputs
814ad4bf 339 ticks->Clear() ;
340 Int_t contrib = 0 ;
341 Float_t a = digit->GetAmp() ;
342 Float_t b = TMath::Abs( a /fTimeSignalLength) ;
88cb7938 343 //Mark the beginning of the signal
814ad4bf 344 new((*ticks)[contrib++]) AliEMCALTick(digit->GetTime(),0, b);
fdebddeb 345 //Mark the end of the signal
814ad4bf 346 new((*ticks)[contrib++]) AliEMCALTick(digit->GetTime()+fTimeSignalLength, -a, -b);
347
88cb7938 348 // loop over input
349 for(i = 0; i< fInput ; i++){ //loop over (possible) merge sources
350 if(dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
351 curSDigit = dynamic_cast<AliEMCALDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
814ad4bf 352 else
353 curSDigit = 0 ;
354 //May be several digits will contribute from the same input
fdebddeb 355 while(curSDigit && (curSDigit->GetId() == absID)){
814ad4bf 356 //Shift primary to separate primaries belonging different inputs
357 Int_t primaryoffset ;
358 if(fManager)
359 primaryoffset = fManager->GetMask(i) ;
360 else
361 primaryoffset = i ;
362 curSDigit->ShiftPrimary(primaryoffset) ;
363
364 a = curSDigit->GetAmp() ;
365 b = a /fTimeSignalLength ;
366 new((*ticks)[contrib++]) AliEMCALTick(curSDigit->GetTime(),0, b);
367 new((*ticks)[contrib++]) AliEMCALTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b);
368
369 *digit = *digit + *curSDigit ; //add energies
370
371 index[i]++ ;
88cb7938 372 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
1963b290 373 curSDigit = dynamic_cast<AliEMCALDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
814ad4bf 374 else
375 curSDigit = 0 ;
376 }
377 }
a9485c37 378 // add fluctuations for photo-electron creation
379 amp = sDigitizer->Calibrate(digit->GetAmp()) ; // GeV
a9485c37 380 amp *= static_cast<Float_t>(gRandom->Poisson(fMeanPhotonElectron)) / static_cast<Float_t>(fMeanPhotonElectron) ;
88cb7938 381
05a92d59 382 //calculate and set time
814ad4bf 383 Float_t time = FrontEdgeTime(ticks) ;
384 digit->SetTime(time) ;
385
386 //Find next signal module
a9485c37 387 nextSig = nEMC + 1 ;
88cb7938 388 for(i = 0 ; i < fInput ; i++){
389 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
814ad4bf 390 Int_t curNext = nextSig ;
391 if(sdigits->GetEntriesFast() > index[i] ){
88cb7938 392 curNext = dynamic_cast<AliEMCALDigit *>(sdigits->At(index[i]))->GetId() ;
814ad4bf 393 }
394 if(curNext < nextSig) nextSig = curNext ;
395 }
396 }
a9485c37 397 // add the noise now
fdebddeb 398 amp += TMath::Abs(gRandom->Gaus(0., fPinNoise)) ;
a9485c37 399 digit->SetAmp(sDigitizer->Digitize(amp)) ;
cf24bdc4 400 AliDebug(10,Form(" absID %5i amp %f nextSig %5i\n",
401 absID, amp, nextSig));
1963b290 402 } // for(absID = 1; absID <= nEMC; absID++)
814ad4bf 403
404 ticks->Delete() ;
405 delete ticks ;
61e0abb5 406
05a92d59 407 delete sdigArray ; //We should not delete its contents
61e0abb5 408
814ad4bf 409 //remove digits below thresholds
88cb7938 410 for(i = 0 ; i < nEMC ; i++){
411 digit = dynamic_cast<AliEMCALDigit*>( digits->At(i) ) ;
fdebddeb 412 Float_t threshold = fDigitThreshold ;
88cb7938 413 if(sDigitizer->Calibrate( digit->GetAmp() ) < threshold)
414 digits->RemoveAt(i) ;
415 else
814ad4bf 416 digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ;
a9485c37 417 }
61e0abb5 418
814ad4bf 419 digits->Compress() ;
61e0abb5 420
05a92d59 421 Int_t ndigits = digits->GetEntriesFast() ;
61e0abb5 422
1963b290 423 //Set indexes in list of digits and fill hists.
315d1c64 424 AliEMCALHistoUtilities::FillH1(fHists, 0, Double_t(ndigits));
1963b290 425 Float_t energy=0., esum=0.;
05a92d59 426 for (i = 0 ; i < ndigits ; i++) {
427 digit = dynamic_cast<AliEMCALDigit *>( digits->At(i) ) ;
814ad4bf 428 digit->SetIndexInList(i) ;
1963b290 429 energy = sDigitizer->Calibrate(digit->GetAmp()) ;
430 esum += energy;
e4dcc484 431 digit->SetAmp(DigitizeEnergy(energy, digit->GetId()) ) ; // for what ??
315d1c64 432 AliEMCALHistoUtilities::FillH1(fHists, 2, double(digit->GetAmp()));
433 AliEMCALHistoUtilities::FillH1(fHists, 3, double(energy));
434 AliEMCALHistoUtilities::FillH1(fHists, 4, double(digit->GetId()));
d25f2c54 435 // if(digit->GetId() == nEMC) {
436 // printf(" i %i \n", i );
437 // digit->Dump();
438 // assert(0);
439 //}
61e0abb5 440 }
315d1c64 441 AliEMCALHistoUtilities::FillH1(fHists, 1, esum);
814ad4bf 442}
61e0abb5 443
d62a891a 444// //_____________________________________________________________________
e4dcc484 445Int_t AliEMCALDigitizer::DigitizeEnergy(Float_t energy, Int_t AbsId)
814ad4bf 446{
e4dcc484 447 // Returns digitized value of the energy in a cell absId
53f1c378 448
449 // Load Geometry
450 const AliEMCALGeometry * geom = AliEMCALGeometry::GetInstance();
e4dcc484 451
452 if (geom==0)
d25f2c54 453 AliFatal("Did not get geometry from EMCALLoader");
e4dcc484 454
455 Int_t iSupMod = -1;
2bb3725c 456 Int_t nModule = -1;
e4dcc484 457 Int_t nIphi = -1;
458 Int_t nIeta = -1;
459 Int_t iphi = -1;
460 Int_t ieta = -1;
fdebddeb 461 Int_t channel = -999;
e4dcc484 462
2bb3725c 463 Bool_t bCell = geom->GetCellIndex(AbsId, iSupMod, nModule, nIphi, nIeta) ;
e4dcc484 464 if(!bCell)
d25f2c54 465 Error("DigitizeEnergy","Wrong cell id number : AbsId %i ", AbsId) ;
2bb3725c 466 geom->GetCellPhiEtaIndexInSModule(iSupMod,nModule,nIphi, nIeta,iphi,ieta);
d62a891a 467
53f1c378 468 if(fCalibData) {
469 fADCpedestalEC = fCalibData->GetADCpedestal(iSupMod,ieta,iphi);
470 fADCchannelEC = fCalibData->GetADCchannel(iSupMod,ieta,iphi);
e4dcc484 471 }
5eb74a51 472
473 channel = static_cast<Int_t> (TMath::Ceil( (energy + fADCpedestalEC)/fADCchannelEC )) ;
474
fdebddeb 475 if(channel > fNADCEC )
7d74eaca 476 channel = fNADCEC ;
814ad4bf 477 return channel ;
d62a891a 478
61e0abb5 479}
480
481//____________________________________________________________________________
05a92d59 482void AliEMCALDigitizer::Exec(Option_t *option)
483{
45fa49ca 484 // Steering method to process digitization for events
485 // in the range from fFirstEvent to fLastEvent.
486 // This range is optionally set by SetEventRange().
487 // if fLastEvent=-1, then process events until the end.
488 // by default fLastEvent = fFirstEvent (process only one event)
05a92d59 489
88cb7938 490 if (!fInit) { // to prevent overwrite existing file
491 Error( "Exec", "Give a version name different from %s", fEventFolderName.Data() ) ;
492 return ;
493 }
494
814ad4bf 495 if (strstr(option,"print")) {
b42d4197 496
88cb7938 497 Print();
814ad4bf 498 return ;
499 }
500
61e0abb5 501 if(strstr(option,"tim"))
502 gBenchmark->Start("EMCALDigitizer");
5dee926e 503
504 AliRunLoader *rl = AliRunLoader::GetRunLoader();
505 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(rl->GetDetectorLoader("EMCAL"));
506
5b4d2c50 507 // Post Digitizer to the white board
5dee926e 508 emcalLoader->PostDigitizer(this) ;
5b4d2c50 509
30bc5594 510 if (fLastEvent == -1) {
5dee926e 511 fLastEvent = rl->GetNumberOfEvents() - 1 ;
30bc5594 512 }
396a348e 513 else if (fManager)
1963b290 514 fLastEvent = fFirstEvent ; // what is this ??
396a348e 515
45fa49ca 516 Int_t nEvents = fLastEvent - fFirstEvent + 1;
cf24bdc4 517 Int_t ievent;
518
5dee926e 519 rl->LoadSDigits("EMCAL");
cf24bdc4 520 for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
64942713 521
5dee926e 522 rl->GetEvent(ievent);
05a92d59 523
814ad4bf 524 Digitize(ievent) ; //Add prepared SDigits to digits and add the noise
88cb7938 525
9e5d2067 526 WriteDigits() ;
88cb7938 527
61e0abb5 528 if(strstr(option,"deb"))
529 PrintDigits(option);
1963b290 530 if(strstr(option,"table")) gObjectTable->Print();
88cb7938 531
814ad4bf 532 //increment the total number of Digits per run
5dee926e 533 fDigitsInRun += emcalLoader->Digits()->GetEntriesFast() ;
61e0abb5 534 }
64942713 535
5dee926e 536 emcalLoader->CleanDigitizer() ;
5b4d2c50 537
61e0abb5 538 if(strstr(option,"tim")){
539 gBenchmark->Stop("EMCALDigitizer");
b42d4197 540 AliInfo(Form("Exec: took %f seconds for Digitizing %f seconds per event",
541 gBenchmark->GetCpuTime("EMCALDigitizer"), gBenchmark->GetCpuTime("EMCALDigitizer")/nEvents )) ;
88cb7938 542 }
61e0abb5 543}
544
05a92d59 545//____________________________________________________________________________
546Float_t AliEMCALDigitizer::FrontEdgeTime(TClonesArray * ticks)
88cb7938 547{
548 // Returns the shortest time among all time ticks
549
05a92d59 550 ticks->Sort() ; //Sort in accordance with times of ticks
551 TIter it(ticks) ;
552 AliEMCALTick * ctick = (AliEMCALTick *) it.Next() ;
553 Float_t time = ctick->CrossingTime(fTimeThreshold) ;
554
555 AliEMCALTick * t ;
556 while((t=(AliEMCALTick*) it.Next())){
557 if(t->GetTime() < time) //This tick starts before crossing
558 *ctick+=*t ;
559 else
560 return time ;
561
562 time = ctick->CrossingTime(fTimeThreshold) ;
563 }
564 return time ;
565}
566
567//____________________________________________________________________________
568Bool_t AliEMCALDigitizer::Init()
569{
570 // Makes all memory allocations
88cb7938 571 fInit = kTRUE ;
5dee926e 572 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
573
574 if ( emcalLoader == 0 ) {
575 Fatal("Init", "Could not obtain the AliEMCALLoader");
05a92d59 576 return kFALSE;
577 }
05a92d59 578
45fa49ca 579 fFirstEvent = 0 ;
580 fLastEvent = fFirstEvent ;
581
88cb7938 582 if(fManager)
583 fInput = fManager->GetNinputs() ;
584 else
585 fInput = 1 ;
586
587 fInputFileNames = new TString[fInput] ;
588 fEventNames = new TString[fInput] ;
589 fInputFileNames[0] = GetTitle() ;
590 fEventNames[0] = fEventFolderName.Data() ;
591 Int_t index ;
592 for (index = 1 ; index < fInput ; index++) {
593 fInputFileNames[index] = dynamic_cast<AliStream*>(fManager->GetInputStream(index))->GetFileName(0);
594 TString tempo = fManager->GetInputFolderName(index) ;
595 fEventNames[index] = tempo.Remove(tempo.Length()-1) ; // strip of the stream number added bt fManager
05a92d59 596 }
597
88cb7938 598 //to prevent cleaning of this object while GetEvent is called
5dee926e 599 emcalLoader->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
e52475ed 600
1c6f613c 601 //PH Print();
53f1c378 602 //Calibration instance
603 fCalibData = emcalLoader->CalibData();
88cb7938 604 return fInit ;
05a92d59 605}
606
607//____________________________________________________________________________
608void AliEMCALDigitizer::InitParameters()
14ce0a6e 609{
29b5d492 610 // Parameter initialization for digitizer
611 // Tune parameters - 24-nov-04; Apr 29, 2007
fe93d365 612 // New parameters JLK 14-Apr-2008
1963b290 613
fe93d365 614 fMeanPhotonElectron = 4400; // electrons per GeV
615 fPinNoise = 0.037; // pin noise in GEV from analysis test beam data
88cb7938 616 if (fPinNoise == 0. )
617 Warning("InitParameters", "No noise added\n") ;
1963b290 618 fDigitThreshold = fPinNoise * 3; // 3 * sigma
619 fTimeResolution = 0.3e-9 ; // 300 psc
05a92d59 620 fTimeSignalLength = 1.0e-9 ;
05a92d59 621
97075cd8 622 // These defaults are normally not used.
623 // Values are read from calibration database instead
624 fADCchannelEC = 0.0153; // Update 24 Apr 2007: 250./16/1024 - width of one ADC channel in GeV
625 fADCpedestalEC = 0.0 ; // GeV
1963b290 626 fNADCEC = (Int_t) TMath::Power(2,16) ; // number of channels in Tower ADC - 65536
a9485c37 627
1963b290 628 fTimeThreshold = 0.001*10000000 ; // Means 1 MeV in terms of SDigits amplitude ??
629 // hists. for control; no hists on default
630 fControlHists = 0;
631 fHists = 0;
05a92d59 632}
61e0abb5 633
61e0abb5 634//__________________________________________________________________
09884213 635void AliEMCALDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
05a92d59 636{
637 // Allows to produce digits by superimposing background and signal event.
61e0abb5 638 // It is assumed, that headers file with SIGNAL events is opened in
05a92d59 639 // the constructor.
640 // Sets the BACKGROUND event, with which the SIGNAL event is to be mixed
641 // Thus we avoid writing (changing) huge and expensive
61e0abb5 642 // backgound files: all output will be writen into SIGNAL, i.e.
643 // opened in constructor file.
644 //
05a92d59 645 // One can open as many files to mix with as one needs.
646 // However only Sdigits with the same name (i.e. constructed with the same SDigitizer)
647 // can be mixed.
61e0abb5 648
05a92d59 649 if( strcmp(GetName(), "") == 0 )
61e0abb5 650 Init() ;
814ad4bf 651
05a92d59 652 if(fManager){
9859bfc0 653 Error("MixWith", "Cannot use this method under AliRunDigitizer") ;
61e0abb5 654 return ;
05a92d59 655 }
88cb7938 656 // looking for file which contains AliRun
657 if (gSystem->AccessPathName(alirunFileName)) {// file does not exist
658 Error("MixWith", "File %s does not exist!", alirunFileName.Data()) ;
659 return ;
61e0abb5 660 }
88cb7938 661 // looking for the file which contains SDigits
5dee926e 662 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
663 TString fileName( emcalLoader->GetSDigitsFileName() ) ;
e191bb57 664 if ( eventFolderName != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
88cb7938 665 fileName = fileName.ReplaceAll(".root", "") + "_" + eventFolderName + ".root" ;
666 if ( (gSystem->AccessPathName(fileName)) ) {
667 Error("MixWith", "The file %s does not exist!", fileName.Data()) ;
668 return ;
814ad4bf 669 }
88cb7938 670 // need to increase the arrays
30bc5594 671 // MvL: This code only works when fInput == 1, I think.
88cb7938 672 TString tempo = fInputFileNames[fInput-1] ;
673 delete [] fInputFileNames ;
674 fInputFileNames = new TString[fInput+1] ;
675 fInputFileNames[fInput-1] = tempo ;
676
677 tempo = fEventNames[fInput-1] ;
678 delete [] fEventNames ;
679 fEventNames = new TString[fInput+1] ;
680 fEventNames[fInput-1] = tempo ;
681
682 fInputFileNames[fInput] = alirunFileName ;
683 fEventNames[fInput] = eventFolderName ;
684 fInput++ ;
685}
1963b290 686
fe93d365 687//__________________________________________________________________
1963b290 688void AliEMCALDigitizer::Print1(Option_t * option)
689{ // 19-nov-04 - just for convinience
690 Print();
691 PrintDigits(option);
692}
693
61e0abb5 694//__________________________________________________________________
eb0b1051 695void AliEMCALDigitizer::Print(Option_t*)const
88cb7938 696{
697 // Print Digitizer's parameters
fdebddeb 698 printf("Print: \n------------------- %s -------------", GetName() ) ;
88cb7938 699 if( strcmp(fEventFolderName.Data(), "") != 0 ){
700 printf(" Writing Digits to branch with title %s\n", fEventFolderName.Data()) ;
701
702 Int_t nStreams ;
703 if (fManager)
704 nStreams = GetNInputStreams() ;
705 else
706 nStreams = fInput ;
707
708 Int_t index = 0 ;
5dee926e 709
710 AliRunLoader *rl=0;
711
88cb7938 712 for (index = 0 ; index < nStreams ; index++) {
713 TString tempo(fEventNames[index]) ;
714 tempo += index ;
5dee926e 715 if ((rl = AliRunLoader::GetRunLoader(tempo)) == 0)
716 rl = AliRunLoader::Open(fInputFileNames[index], tempo) ;
717 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(rl->GetDetectorLoader("EMCAL"));
718 TString fileName( emcalLoader->GetSDigitsFileName() ) ;
e191bb57 719 if ( fEventNames[index] != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
88cb7938 720 fileName = fileName.ReplaceAll(".root", "") + "_" + fEventNames[index] + ".root" ;
721 printf ("Adding SDigits from %s %s\n", fInputFileNames[index].Data(), fileName.Data()) ;
722 }
5dee926e 723
724 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
725
726 printf("\nWriting digits to %s", emcalLoader->GetDigitsFileName().Data()) ;
88cb7938 727
728 printf("\nWith following parameters:\n") ;
729
730 printf(" Electronics noise in EMC (fPinNoise) = %f\n", fPinNoise) ;
731 printf(" Threshold in EMC (fDigitThreshold) = %f\n", fDigitThreshold) ;
732 printf("---------------------------------------------------\n") ;
ea5f3389 733 }
734 else
fdebddeb 735 printf("Print: AliEMCALDigitizer not initialized") ;
61e0abb5 736}
814ad4bf 737
61e0abb5 738//__________________________________________________________________
14ce0a6e 739void AliEMCALDigitizer::PrintDigits(Option_t * option)
740{
741 //utility method for printing digit information
742
5dee926e 743 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
744 TClonesArray * digits = emcalLoader->Digits() ;
745 TClonesArray * sdigits = emcalLoader->SDigits() ;
88cb7938 746
1963b290 747 printf("\n #Digits: %d : sdigits %d ", digits->GetEntriesFast(), sdigits->GetEntriesFast()) ;
5dee926e 748 printf("\n event %d", emcalLoader->GetRunLoader()->GetEventNumber());
88cb7938 749
750 if(strstr(option,"all")){
61e0abb5 751 //loop over digits
752 AliEMCALDigit * digit;
88cb7938 753 printf("\nEMC digits (with primaries):\n") ;
754 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
61e0abb5 755 Int_t index ;
88cb7938 756 for (index = 0 ; index < digits->GetEntries() ; index++) {
757 digit = dynamic_cast<AliEMCALDigit *>(digits->At(index)) ;
758 printf("\n%6d %8d %6.5e %4d %2d : ",
11f9c5ff 759 digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
61e0abb5 760 Int_t iprimary;
9859bfc0 761 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
88cb7938 762 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
9859bfc0 763 }
764 }
61e0abb5 765 }
1963b290 766 printf("\n");
61e0abb5 767}
05a92d59 768
769//__________________________________________________________________
770Float_t AliEMCALDigitizer::TimeOfNoise(void)
04f0bda3 771{
772 // Calculates the time signal generated by noise
26a2ef9d 773 //PH Info("TimeOfNoise", "Change me") ;
04f0bda3 774 return gRandom->Rndm() * 1.28E-5;
05a92d59 775}
776
88cb7938 777//__________________________________________________________________
778void AliEMCALDigitizer::Unload()
779{
fdebddeb 780 // Unloads the SDigits and Digits
5dee926e 781 AliRunLoader *rl=0;
782
88cb7938 783 Int_t i ;
784 for(i = 1 ; i < fInput ; i++){
785 TString tempo(fEventNames[i]) ;
5dee926e 786 tempo += i;
787 if ((rl = AliRunLoader::GetRunLoader(tempo)))
788 rl->GetDetectorLoader("EMCAL")->UnloadSDigits() ;
88cb7938 789 }
5dee926e 790 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
791 emcalLoader->UnloadDigits() ;
88cb7938 792}
793
814ad4bf 794//_________________________________________________________________________________________
9e5d2067 795void AliEMCALDigitizer::WriteDigits()
814ad4bf 796{
61e0abb5 797
814ad4bf 798 // Makes TreeD in the output file.
799 // Check if branch already exists:
800 // if yes, exit without writing: ROOT TTree does not support overwriting/updating of
801 // already existing branches.
802 // else creates branch with Digits, named "EMCAL", title "...",
803 // and branch "AliEMCALDigitizer", with the same title to keep all the parameters
804 // and names of files, from which digits are made.
61e0abb5 805
5dee926e 806 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
cf24bdc4 807
5dee926e 808 const TClonesArray * digits = emcalLoader->Digits() ;
809 TTree * treeD = emcalLoader->TreeD();
810 if ( !treeD ) {
811 emcalLoader->MakeDigitsContainer();
812 treeD = emcalLoader->TreeD();
813 }
88cb7938 814
815 // -- create Digits branch
816 Int_t bufferSize = 32000 ;
5dee926e 817 TBranch * digitsBranch = 0;
30bc5594 818 if ((digitsBranch = treeD->GetBranch("EMCAL"))) {
5dee926e 819 digitsBranch->SetAddress(&digits);
30bc5594 820 AliWarning("Digits Branch already exists. Not all digits will be visible");
821 }
5dee926e 822 else
823 treeD->Branch("EMCAL","TClonesArray",&digits,bufferSize);
824 //digitsBranch->SetTitle(fEventFolderName);
825 treeD->Fill() ;
88cb7938 826
5dee926e 827 emcalLoader->WriteDigits("OVERWRITE");
828 emcalLoader->WriteDigitizer("OVERWRITE");
88cb7938 829
830 Unload() ;
831
814ad4bf 832}
61e0abb5 833
fe93d365 834//__________________________________________________________________
1963b290 835void AliEMCALDigitizer::Browse(TBrowser* b)
836{
837 if(fHists) b->Add(fHists);
838 TTask::Browse(b);
839}
840
fe93d365 841//__________________________________________________________________
eb0b1051 842TList *AliEMCALDigitizer::BookControlHists(int var)
14ce0a6e 843{
844 // 22-nov-04
845 // histograms for monitoring digitizer performance
846
1963b290 847 Info("BookControlHists"," started ");
848 gROOT->cd();
5dee926e 849 const AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance();
1963b290 850 if(var>=1){
851 new TH1F("hDigiN", "#EMCAL digits with fAmp > fDigitThreshold",
852 fNADCEC+1, -0.5, Double_t(fNADCEC));
5dee926e 853 new TH1F("HDigiSumEnergy","Sum.EMCAL energy from digi", 1000, 0.0, 200.);
1963b290 854 new TH1F("hDigiAmp", "EMCAL digital amplitude", fNADCEC+1, -0.5, Double_t(fNADCEC));
855 new TH1F("hDigiEnergy","EMCAL cell energy", 2000, 0.0, 200.);
856 new TH1F("hDigiAbsId","EMCAL absId cells with fAmp > fDigitThreshold ",
857 geom->GetNCells(), 0.5, Double_t(geom->GetNCells())+0.5);
858 }
64942713 859
315d1c64 860 fHists = AliEMCALHistoUtilities::MoveHistsToList("EmcalDigiControlHists", kFALSE);
64942713 861 fHists = 0; //huh? JLK 03-Mar-2006
1963b290 862 return fHists;
863}
864
fe93d365 865//__________________________________________________________________
1963b290 866void AliEMCALDigitizer::SaveHists(const char* name, Bool_t kSingleKey, const char* opt)
867{
315d1c64 868 AliEMCALHistoUtilities::SaveListOfHists(fHists, name, kSingleKey, opt);
1963b290 869}