updating apd and fee parameters
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALDigitizer.cxx
... / ...
CommitLineData
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//_________________________________________________________________________
19//
20//////////////////////////////////////////////////////////////////////////////
21// Class performs digitization of Summable digits
22//
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
51////////////////////////////////////////////////////////////////////////////////////
52//
53//*-- Author: Sahal Yacoob (LBL)
54// based on : AliEMCALDigitizer
55// Modif:
56// August 2002 Yves Schutz: clone PHOS as closely as possible and intoduction
57// of new IO (à la PHOS)
58// November 2003 Aleksei Pavlinov : adopted for Shish-Kebab geometry
59//_________________________________________________________________________________
60
61// --- ROOT system ---
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>
71#include <cassert>
72
73// --- AliRoot header files ---
74#include "AliLog.h"
75#include "AliRun.h"
76#include "AliRunDigitizer.h"
77#include "AliRunLoader.h"
78#include "AliCDBManager.h"
79#include "AliCDBEntry.h"
80#include "AliEMCALDigit.h"
81#include "AliEMCAL.h"
82#include "AliEMCALLoader.h"
83#include "AliEMCALDigitizer.h"
84#include "AliEMCALSDigitizer.h"
85#include "AliEMCALGeometry.h"
86#include "AliEMCALTick.h"
87#include "AliEMCALHistoUtilities.h"
88
89ClassImp(AliEMCALDigitizer)
90
91
92//____________________________________________________________________________
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),
116 fHists(0),fCalibData(0x0)
117{
118 // ctor
119 InitParameters() ;
120 fManager = 0 ; // We work in the standalong mode
121}
122
123//____________________________________________________________________________
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),
147 fHists(0),fCalibData(0x0)
148{
149 // ctor
150 InitParameters() ;
151 Init() ;
152 fManager = 0 ; // We work in the standalong mode
153}
154
155//____________________________________________________________________________
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),
179 fHists(d.fHists),fCalibData(d.fCalibData)
180{
181 // copyy ctor
182 }
183
184//____________________________________________________________________________
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),
208 fHists(0),fCalibData(0x0)
209{
210 // ctor Init() is called by RunDigitizer
211 fManager = rd ;
212 fEventFolderName = fManager->GetInputFolderName(0) ;
213 SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
214 InitParameters() ;
215}
216
217//____________________________________________________________________________
218 AliEMCALDigitizer::~AliEMCALDigitizer()
219{
220 //dtor
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");
228 delete [] fInputFileNames ;
229 delete [] fEventNames ;
230
231 if(fHists) delete fHists;
232}
233
234//____________________________________________________________________________
235void AliEMCALDigitizer::Digitize(Int_t event)
236{
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.
244 static int nEMC=0; //max number of digits possible
245
246 AliRunLoader *rl = AliRunLoader::GetRunLoader();
247 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(rl->GetDetectorLoader("EMCAL"));
248 Int_t readEvent = event ;
249 // fManager is data member from AliDigitizer
250 if (fManager)
251 readEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetCurrentEventNumber() ;
252 AliDebug(1,Form("Adding event %d from input stream 0 %s %s",
253 readEvent, GetTitle(), fEventFolderName.Data())) ;
254 rl->GetEvent(readEvent);
255
256 TClonesArray * digits = emcalLoader->Digits() ;
257 digits->Delete() ;
258
259 // Load Geometry
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");
267
268 nEMC = geom->GetNCells();
269 AliDebug(1,Form("nEMC %i (number cells in EMCAL) | %s \n", nEMC, geom->GetName()));
270
271 Int_t absID ;
272
273 digits->Expand(nEMC) ;
274
275 // get first the sdigitizer from the tasks list (must have same name as the digitizer)
276 if ( !emcalLoader->SDigitizer() )
277 emcalLoader->LoadSDigitizer();
278 AliEMCALSDigitizer * sDigitizer = dynamic_cast<AliEMCALSDigitizer*>(emcalLoader->SDigitizer());
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) ;
285 sdigArray->AddAt(emcalLoader->SDigits(), 0) ;
286 Int_t i ;
287
288 for(i = 1 ; i < fInput ; i++){
289 TString tempo(fEventNames[i]) ;
290 tempo += i ;
291
292 AliRunLoader * rl2 = AliRunLoader::GetRunLoader(tempo) ;
293
294 if (rl2==0)
295 rl2 = AliRunLoader::Open(fInputFileNames[i], tempo) ;
296
297 if (fManager)
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()) ;
300 rl2->LoadSDigits();
301 rl2->GetEvent(readEvent);
302 AliEMCALLoader *emcalLoader2 = dynamic_cast<AliEMCALLoader*>(rl2->GetDetectorLoader("EMCAL"));
303 sdigArray->AddAt(emcalLoader2->SDigits(), i) ;
304 }
305
306 //Find the first tower with signal
307 Int_t nextSig = nEMC + 1 ;
308 TClonesArray * sdigits ;
309 for(i = 0 ; i < fInput ; i++){
310 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
311 if ( !sdigits->GetEntriesFast() )
312 continue ;
313 Int_t curNext = dynamic_cast<AliEMCALDigit *>(sdigits->At(0))->GetId() ;
314 if(curNext < nextSig)
315 nextSig = curNext ;
316 AliDebug(1,Form("input %i : #sdigits %i \n",
317 i, sdigits->GetEntriesFast()));
318 }
319 AliDebug(1,Form("FIRST tower with signal %i \n", nextSig));
320
321 TArrayI index(fInput) ;
322 index.Reset() ; //Set all indexes to zero
323
324 AliEMCALDigit * digit ;
325 AliEMCALDigit * curSDigit ;
326
327 TClonesArray * ticks = new TClonesArray("AliEMCALTick",1000) ;
328
329 //Put Noise contribution
330 for(absID = 0; absID < nEMC; absID++){ // Nov 30, 2006 by PAI; was from 1 to nEMC
331 Float_t amp = 0 ;
332 // amplitude set to zero, noise will be added later
333 new((*digits)[absID]) AliEMCALDigit( -1, -1, absID, 0, TimeOfNoise() ); // absID-1->absID
334 //look if we have to add signal?
335 digit = dynamic_cast<AliEMCALDigit *>(digits->At(absID)); // absID-1->absID
336
337 if(absID==nextSig){
338 //Add SDigits from all inputs
339 ticks->Clear() ;
340 Int_t contrib = 0 ;
341 Float_t a = digit->GetAmp() ;
342 Float_t b = TMath::Abs( a /fTimeSignalLength) ;
343 //Mark the beginning of the signal
344 new((*ticks)[contrib++]) AliEMCALTick(digit->GetTime(),0, b);
345 //Mark the end of the signal
346 new((*ticks)[contrib++]) AliEMCALTick(digit->GetTime()+fTimeSignalLength, -a, -b);
347
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])) ;
352 else
353 curSDigit = 0 ;
354 //May be several digits will contribute from the same input
355 while(curSDigit && (curSDigit->GetId() == absID)){
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]++ ;
372 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
373 curSDigit = dynamic_cast<AliEMCALDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
374 else
375 curSDigit = 0 ;
376 }
377 }
378 // add fluctuations for photo-electron creation
379 amp = sDigitizer->Calibrate(digit->GetAmp()) ; // GeV
380 amp *= static_cast<Float_t>(gRandom->Poisson(fMeanPhotonElectron)) / static_cast<Float_t>(fMeanPhotonElectron) ;
381
382 //calculate and set time
383 Float_t time = FrontEdgeTime(ticks) ;
384 digit->SetTime(time) ;
385
386 //Find next signal module
387 nextSig = nEMC + 1 ;
388 for(i = 0 ; i < fInput ; i++){
389 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
390 Int_t curNext = nextSig ;
391 if(sdigits->GetEntriesFast() > index[i] ){
392 curNext = dynamic_cast<AliEMCALDigit *>(sdigits->At(index[i]))->GetId() ;
393 }
394 if(curNext < nextSig) nextSig = curNext ;
395 }
396 }
397 // add the noise now
398 amp += TMath::Abs(gRandom->Gaus(0., fPinNoise)) ;
399 digit->SetAmp(sDigitizer->Digitize(amp)) ;
400 AliDebug(10,Form(" absID %5i amp %f nextSig %5i\n",
401 absID, amp, nextSig));
402 } // for(absID = 1; absID <= nEMC; absID++)
403
404 ticks->Delete() ;
405 delete ticks ;
406
407 delete sdigArray ; //We should not delete its contents
408
409 //remove digits below thresholds
410 for(i = 0 ; i < nEMC ; i++){
411 digit = dynamic_cast<AliEMCALDigit*>( digits->At(i) ) ;
412 Float_t threshold = fDigitThreshold ;
413 if(sDigitizer->Calibrate( digit->GetAmp() ) < threshold)
414 digits->RemoveAt(i) ;
415 else
416 digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ;
417 }
418
419 digits->Compress() ;
420
421 Int_t ndigits = digits->GetEntriesFast() ;
422
423 //Set indexes in list of digits and fill hists.
424 AliEMCALHistoUtilities::FillH1(fHists, 0, Double_t(ndigits));
425 Float_t energy=0., esum=0.;
426 for (i = 0 ; i < ndigits ; i++) {
427 digit = dynamic_cast<AliEMCALDigit *>( digits->At(i) ) ;
428 digit->SetIndexInList(i) ;
429 energy = sDigitizer->Calibrate(digit->GetAmp()) ;
430 esum += energy;
431 digit->SetAmp(DigitizeEnergy(energy, digit->GetId()) ) ; // for what ??
432 AliEMCALHistoUtilities::FillH1(fHists, 2, double(digit->GetAmp()));
433 AliEMCALHistoUtilities::FillH1(fHists, 3, double(energy));
434 AliEMCALHistoUtilities::FillH1(fHists, 4, double(digit->GetId()));
435 // if(digit->GetId() == nEMC) {
436 // printf(" i %i \n", i );
437 // digit->Dump();
438 // assert(0);
439 //}
440 }
441 AliEMCALHistoUtilities::FillH1(fHists, 1, esum);
442}
443
444// //_____________________________________________________________________
445Int_t AliEMCALDigitizer::DigitizeEnergy(Float_t energy, Int_t AbsId)
446{
447 // Returns digitized value of the energy in a cell absId
448
449 // Load Geometry
450 const AliEMCALGeometry * geom = AliEMCALGeometry::GetInstance();
451
452 if (geom==0)
453 AliFatal("Did not get geometry from EMCALLoader");
454
455 Int_t iSupMod = -1;
456 Int_t nModule = -1;
457 Int_t nIphi = -1;
458 Int_t nIeta = -1;
459 Int_t iphi = -1;
460 Int_t ieta = -1;
461 Int_t channel = -999;
462
463 Bool_t bCell = geom->GetCellIndex(AbsId, iSupMod, nModule, nIphi, nIeta) ;
464 if(!bCell)
465 Error("DigitizeEnergy","Wrong cell id number : AbsId %i ", AbsId) ;
466 geom->GetCellPhiEtaIndexInSModule(iSupMod,nModule,nIphi, nIeta,iphi,ieta);
467
468 if(fCalibData) {
469 fADCpedestalEC = fCalibData->GetADCpedestal(iSupMod,ieta,iphi);
470 fADCchannelEC = fCalibData->GetADCchannel(iSupMod,ieta,iphi);
471 }
472
473 channel = static_cast<Int_t> (TMath::Ceil( (energy + fADCpedestalEC)/fADCchannelEC )) ;
474
475 if(channel > fNADCEC )
476 channel = fNADCEC ;
477 return channel ;
478
479}
480
481//____________________________________________________________________________
482void AliEMCALDigitizer::Exec(Option_t *option)
483{
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)
489
490 if (!fInit) { // to prevent overwrite existing file
491 Error( "Exec", "Give a version name different from %s", fEventFolderName.Data() ) ;
492 return ;
493 }
494
495 if (strstr(option,"print")) {
496
497 Print();
498 return ;
499 }
500
501 if(strstr(option,"tim"))
502 gBenchmark->Start("EMCALDigitizer");
503
504 AliRunLoader *rl = AliRunLoader::GetRunLoader();
505 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(rl->GetDetectorLoader("EMCAL"));
506
507 // Post Digitizer to the white board
508 emcalLoader->PostDigitizer(this) ;
509
510 if (fLastEvent == -1) {
511 fLastEvent = rl->GetNumberOfEvents() - 1 ;
512 }
513 else if (fManager)
514 fLastEvent = fFirstEvent ; // what is this ??
515
516 Int_t nEvents = fLastEvent - fFirstEvent + 1;
517 Int_t ievent;
518
519 rl->LoadSDigits("EMCAL");
520 for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
521
522 rl->GetEvent(ievent);
523
524 Digitize(ievent) ; //Add prepared SDigits to digits and add the noise
525
526 WriteDigits() ;
527
528 if(strstr(option,"deb"))
529 PrintDigits(option);
530 if(strstr(option,"table")) gObjectTable->Print();
531
532 //increment the total number of Digits per run
533 fDigitsInRun += emcalLoader->Digits()->GetEntriesFast() ;
534 }
535
536 emcalLoader->CleanDigitizer() ;
537
538 if(strstr(option,"tim")){
539 gBenchmark->Stop("EMCALDigitizer");
540 AliInfo(Form("Exec: took %f seconds for Digitizing %f seconds per event",
541 gBenchmark->GetCpuTime("EMCALDigitizer"), gBenchmark->GetCpuTime("EMCALDigitizer")/nEvents )) ;
542 }
543}
544
545//____________________________________________________________________________
546Float_t AliEMCALDigitizer::FrontEdgeTime(TClonesArray * ticks)
547{
548 // Returns the shortest time among all time ticks
549
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
571 fInit = kTRUE ;
572 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
573
574 if ( emcalLoader == 0 ) {
575 Fatal("Init", "Could not obtain the AliEMCALLoader");
576 return kFALSE;
577 }
578
579 fFirstEvent = 0 ;
580 fLastEvent = fFirstEvent ;
581
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
596 }
597
598 //to prevent cleaning of this object while GetEvent is called
599 emcalLoader->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
600
601 //PH Print();
602 //Calibration instance
603 fCalibData = emcalLoader->CalibData();
604 return fInit ;
605}
606
607//____________________________________________________________________________
608void AliEMCALDigitizer::InitParameters()
609{
610 // Parameter initialization for digitizer
611 // Tune parameters - 24-nov-04; Apr 29, 2007
612 // New parameters JLK 14-Apr-2008
613
614 fMeanPhotonElectron = 4400; // electrons per GeV
615 fPinNoise = 0.037; // pin noise in GEV from analysis test beam data
616 if (fPinNoise == 0. )
617 Warning("InitParameters", "No noise added\n") ;
618 fDigitThreshold = fPinNoise * 3; // 3 * sigma
619 fTimeResolution = 0.3e-9 ; // 300 psc
620 fTimeSignalLength = 1.0e-9 ;
621
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
626 fNADCEC = (Int_t) TMath::Power(2,16) ; // number of channels in Tower ADC - 65536
627
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;
632}
633
634//__________________________________________________________________
635void AliEMCALDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
636{
637 // Allows to produce digits by superimposing background and signal event.
638 // It is assumed, that headers file with SIGNAL events is opened in
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
642 // backgound files: all output will be writen into SIGNAL, i.e.
643 // opened in constructor file.
644 //
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.
648
649 if( strcmp(GetName(), "") == 0 )
650 Init() ;
651
652 if(fManager){
653 Error("MixWith", "Cannot use this method under AliRunDigitizer") ;
654 return ;
655 }
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 ;
660 }
661 // looking for the file which contains SDigits
662 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
663 TString fileName( emcalLoader->GetSDigitsFileName() ) ;
664 if ( eventFolderName != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
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 ;
669 }
670 // need to increase the arrays
671 // MvL: This code only works when fInput == 1, I think.
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}
686
687//__________________________________________________________________
688void AliEMCALDigitizer::Print1(Option_t * option)
689{ // 19-nov-04 - just for convinience
690 Print();
691 PrintDigits(option);
692}
693
694//__________________________________________________________________
695void AliEMCALDigitizer::Print(Option_t*)const
696{
697 // Print Digitizer's parameters
698 printf("Print: \n------------------- %s -------------", GetName() ) ;
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 ;
709
710 AliRunLoader *rl=0;
711
712 for (index = 0 ; index < nStreams ; index++) {
713 TString tempo(fEventNames[index]) ;
714 tempo += index ;
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() ) ;
719 if ( fEventNames[index] != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
720 fileName = fileName.ReplaceAll(".root", "") + "_" + fEventNames[index] + ".root" ;
721 printf ("Adding SDigits from %s %s\n", fInputFileNames[index].Data(), fileName.Data()) ;
722 }
723
724 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
725
726 printf("\nWriting digits to %s", emcalLoader->GetDigitsFileName().Data()) ;
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") ;
733 }
734 else
735 printf("Print: AliEMCALDigitizer not initialized") ;
736}
737
738//__________________________________________________________________
739void AliEMCALDigitizer::PrintDigits(Option_t * option)
740{
741 //utility method for printing digit information
742
743 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
744 TClonesArray * digits = emcalLoader->Digits() ;
745 TClonesArray * sdigits = emcalLoader->SDigits() ;
746
747 printf("\n #Digits: %d : sdigits %d ", digits->GetEntriesFast(), sdigits->GetEntriesFast()) ;
748 printf("\n event %d", emcalLoader->GetRunLoader()->GetEventNumber());
749
750 if(strstr(option,"all")){
751 //loop over digits
752 AliEMCALDigit * digit;
753 printf("\nEMC digits (with primaries):\n") ;
754 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
755 Int_t index ;
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 : ",
759 digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
760 Int_t iprimary;
761 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
762 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
763 }
764 }
765 }
766 printf("\n");
767}
768
769//__________________________________________________________________
770Float_t AliEMCALDigitizer::TimeOfNoise(void)
771{
772 // Calculates the time signal generated by noise
773 //PH Info("TimeOfNoise", "Change me") ;
774 return gRandom->Rndm() * 1.28E-5;
775}
776
777//__________________________________________________________________
778void AliEMCALDigitizer::Unload()
779{
780 // Unloads the SDigits and Digits
781 AliRunLoader *rl=0;
782
783 Int_t i ;
784 for(i = 1 ; i < fInput ; i++){
785 TString tempo(fEventNames[i]) ;
786 tempo += i;
787 if ((rl = AliRunLoader::GetRunLoader(tempo)))
788 rl->GetDetectorLoader("EMCAL")->UnloadSDigits() ;
789 }
790 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
791 emcalLoader->UnloadDigits() ;
792}
793
794//_________________________________________________________________________________________
795void AliEMCALDigitizer::WriteDigits()
796{
797
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.
805
806 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
807
808 const TClonesArray * digits = emcalLoader->Digits() ;
809 TTree * treeD = emcalLoader->TreeD();
810 if ( !treeD ) {
811 emcalLoader->MakeDigitsContainer();
812 treeD = emcalLoader->TreeD();
813 }
814
815 // -- create Digits branch
816 Int_t bufferSize = 32000 ;
817 TBranch * digitsBranch = 0;
818 if ((digitsBranch = treeD->GetBranch("EMCAL"))) {
819 digitsBranch->SetAddress(&digits);
820 AliWarning("Digits Branch already exists. Not all digits will be visible");
821 }
822 else
823 treeD->Branch("EMCAL","TClonesArray",&digits,bufferSize);
824 //digitsBranch->SetTitle(fEventFolderName);
825 treeD->Fill() ;
826
827 emcalLoader->WriteDigits("OVERWRITE");
828 emcalLoader->WriteDigitizer("OVERWRITE");
829
830 Unload() ;
831
832}
833
834//__________________________________________________________________
835void AliEMCALDigitizer::Browse(TBrowser* b)
836{
837 if(fHists) b->Add(fHists);
838 TTask::Browse(b);
839}
840
841//__________________________________________________________________
842TList *AliEMCALDigitizer::BookControlHists(int var)
843{
844 // 22-nov-04
845 // histograms for monitoring digitizer performance
846
847 Info("BookControlHists"," started ");
848 gROOT->cd();
849 const AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance();
850 if(var>=1){
851 new TH1F("hDigiN", "#EMCAL digits with fAmp > fDigitThreshold",
852 fNADCEC+1, -0.5, Double_t(fNADCEC));
853 new TH1F("HDigiSumEnergy","Sum.EMCAL energy from digi", 1000, 0.0, 200.);
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 }
859
860 fHists = AliEMCALHistoUtilities::MoveHistsToList("EmcalDigiControlHists", kFALSE);
861 fHists = 0; //huh? JLK 03-Mar-2006
862 return fHists;
863}
864
865//__________________________________________________________________
866void AliEMCALDigitizer::SaveHists(const char* name, Bool_t kSingleKey, const char* opt)
867{
868 AliEMCALHistoUtilities::SaveListOfHists(fHists, name, kSingleKey, opt);
869}