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 **************************************************************************/
18 //_________________________________________________________________________
20 //////////////////////////////////////////////////////////////////////////////
21 // Class performs digitization of Summable digits
23 // In addition it performs mixing of summable digits from different events.
25 // For each event two branches are created in TreeD:
26 // "EMCAL" - list of digits
27 // "AliEMCALDigitizer" - AliEMCALDigitizer with all parameters used in digitization
29 // Note, that one cset title for new digits branch, and repeat digitization with
30 // another set of parameters.
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
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 ////////////////////////////////////////////////////////////////////////////////////
53 //*-- Author: Sahal Yacoob (LBL)
54 // based on : AliEMCALDigitizer
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 //_________________________________________________________________________________
61 // --- ROOT system ---
65 #include <TBenchmark.h>
69 #include <TObjectTable.h>
73 // --- AliRoot header files ---
76 #include "AliRunDigitizer.h"
77 #include "AliRunLoader.h"
78 #include "AliCDBManager.h"
79 #include "AliCDBEntry.h"
80 #include "AliEMCALDigit.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"
89 ClassImp(AliEMCALDigitizer)
92 //____________________________________________________________________________
93 AliEMCALDigitizer::AliEMCALDigitizer()
94 : AliDigitizer("",""),
102 fMeanPhotonElectron(0),
108 fTimeSignalLength(0),
112 fEventFolderName(""),
116 fHists(0),fCalibData(0x0)
120 fManager = 0 ; // We work in the standalong mode
123 //____________________________________________________________________________
124 AliEMCALDigitizer::AliEMCALDigitizer(TString alirunFileName, TString eventFolderName)
125 : AliDigitizer("EMCAL"+AliConfig::Instance()->GetDigitizerTaskName(), alirunFileName),
126 fDefaultInit(kFALSE),
133 fMeanPhotonElectron(0),
139 fTimeSignalLength(0),
143 fEventFolderName(eventFolderName),
147 fHists(0),fCalibData(0x0)
152 fManager = 0 ; // We work in the standalong mode
155 //____________________________________________________________________________
156 AliEMCALDigitizer::AliEMCALDigitizer(const AliEMCALDigitizer & d)
157 : AliDigitizer(d.GetName(),d.GetTitle()),
158 fDefaultInit(d.fDefaultInit),
159 fDigitsInRun(d.fDigitsInRun),
162 fInputFileNames(d.fInputFileNames),
163 fEventNames(d.fEventNames),
164 fDigitThreshold(d.fDigitThreshold),
165 fMeanPhotonElectron(d.fMeanPhotonElectron),
166 fPedestal(d.fPedestal),
168 fPinNoise(d.fPinNoise),
169 fTimeResolution(d.fTimeResolution),
170 fTimeThreshold(d.fTimeThreshold),
171 fTimeSignalLength(d.fTimeSignalLength),
172 fADCchannelEC(d.fADCchannelEC),
173 fADCpedestalEC(d.fADCpedestalEC),
175 fEventFolderName(d.fEventFolderName),
176 fFirstEvent(d.fFirstEvent),
177 fLastEvent(d.fLastEvent),
178 fControlHists(d.fControlHists),
179 fHists(d.fHists),fCalibData(d.fCalibData)
184 //____________________________________________________________________________
185 AliEMCALDigitizer::AliEMCALDigitizer(AliRunDigitizer * rd)
186 : AliDigitizer(rd,"EMCAL"+AliConfig::Instance()->GetDigitizerTaskName()),
187 fDefaultInit(kFALSE),
194 fMeanPhotonElectron(0),
200 fTimeSignalLength(0),
208 fHists(0),fCalibData(0x0)
210 // ctor Init() is called by RunDigitizer
212 fEventFolderName = fManager->GetInputFolderName(0) ;
213 SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
217 //____________________________________________________________________________
218 AliEMCALDigitizer::~AliEMCALDigitizer()
221 if (AliRunLoader::GetRunLoader()) {
222 AliLoader *emcalLoader=0;
223 if ((emcalLoader = AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL")))
224 emcalLoader->CleanDigitizer();
227 AliDebug(1," no runloader present");
228 delete [] fInputFileNames ;
229 delete [] fEventNames ;
231 if(fHists) delete fHists;
234 //____________________________________________________________________________
235 void AliEMCALDigitizer::Digitize(Int_t event)
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 isTrd1Geom = -1; // -1 - mean undefined
245 static int nEMC=0; //max number of digits possible
247 AliRunLoader *rl = AliRunLoader::GetRunLoader();
248 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(rl->GetDetectorLoader("EMCAL"));
249 Int_t readEvent = event ;
250 // fManager is data member from AliDigitizer
252 readEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetCurrentEventNumber() ;
253 AliDebug(1,Form("Adding event %d from input stream 0 %s %s",
254 readEvent, GetTitle(), fEventFolderName.Data())) ;
255 rl->GetEvent(readEvent);
257 TClonesArray * digits = emcalLoader->Digits() ;
261 AliEMCALGeometry *geom = 0;
262 if (rl->GetAliRun()) {
263 AliEMCAL * emcal = (AliEMCAL*)rl->GetAliRun()->GetDetector("EMCAL");
264 geom = emcal->GetGeometry();
267 AliFatal("Could not get AliRun from runLoader");
270 AliDebug(1, Form(" get Geometry %s : %s ", geom->GetName(),geom->GetTitle()));
271 TString ng(geom->GetName());
273 if(ng.Contains("SHISH") && ng.Contains("TRD1")) isTrd1Geom = 1;
275 if(isTrd1Geom == 0) nEMC = geom->GetNPhi()*geom->GetNZ();
276 else nEMC = geom->GetNCells();
277 AliDebug(1,Form("nEMC %i (number cells in EMCAL) | %s | isTrd1Geom %i\n", nEMC, geom->GetName(), isTrd1Geom));
281 digits->Expand(nEMC) ;
283 // get first the sdigitizer from the tasks list (must have same name as the digitizer)
284 if ( !emcalLoader->SDigitizer() )
285 emcalLoader->LoadSDigitizer();
286 AliEMCALSDigitizer * sDigitizer = dynamic_cast<AliEMCALSDigitizer*>(emcalLoader->SDigitizer());
289 Fatal("Digitize", "SDigitizer with name %s %s not found", fEventFolderName.Data(), GetTitle() ) ;
291 //take all the inputs to add together and load the SDigits
292 TObjArray * sdigArray = new TObjArray(fInput) ;
293 sdigArray->AddAt(emcalLoader->SDigits(), 0) ;
296 for(i = 1 ; i < fInput ; i++){
297 TString tempo(fEventNames[i]) ;
300 AliRunLoader * rl2 = AliRunLoader::GetRunLoader(tempo) ;
303 rl2 = AliRunLoader::Open(fInputFileNames[i], tempo) ;
306 readEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(i))->GetCurrentEventNumber() ;
307 Info("Digitize", "Adding event %d from input stream %d %s %s", readEvent, i, fInputFileNames[i].Data(), tempo.Data()) ;
309 rl2->GetEvent(readEvent);
310 AliEMCALLoader *emcalLoader2 = dynamic_cast<AliEMCALLoader*>(rl2->GetDetectorLoader("EMCAL"));
311 sdigArray->AddAt(emcalLoader2->SDigits(), i) ;
314 //Find the first tower with signal
315 Int_t nextSig = nEMC + 1 ;
316 TClonesArray * sdigits ;
317 for(i = 0 ; i < fInput ; i++){
318 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
319 if ( !sdigits->GetEntriesFast() )
321 Int_t curNext = dynamic_cast<AliEMCALDigit *>(sdigits->At(0))->GetId() ;
322 if(curNext < nextSig)
324 AliDebug(1,Form("input %i : #sdigits %i \n",
325 i, sdigits->GetEntriesFast()));
327 AliDebug(1,Form("FIRST tower with signal %i \n", nextSig));
329 TArrayI index(fInput) ;
330 index.Reset() ; //Set all indexes to zero
332 AliEMCALDigit * digit ;
333 AliEMCALDigit * curSDigit ;
335 TClonesArray * ticks = new TClonesArray("AliEMCALTick",1000) ;
337 //Put Noise contribution
338 for(absID = 0; absID < nEMC; absID++){ // Nov 30, 2006 by PAI; was from 1 to nEMC
340 // amplitude set to zero, noise will be added later
341 new((*digits)[absID]) AliEMCALDigit( -1, -1, absID, 0, TimeOfNoise() ); // absID-1->absID
342 //look if we have to add signal?
343 digit = dynamic_cast<AliEMCALDigit *>(digits->At(absID)); // absID-1->absID
346 //Add SDigits from all inputs
349 Float_t a = digit->GetAmp() ;
350 Float_t b = TMath::Abs( a /fTimeSignalLength) ;
351 //Mark the beginning of the signal
352 new((*ticks)[contrib++]) AliEMCALTick(digit->GetTime(),0, b);
353 //Mark the end of the signal
354 new((*ticks)[contrib++]) AliEMCALTick(digit->GetTime()+fTimeSignalLength, -a, -b);
357 for(i = 0; i< fInput ; i++){ //loop over (possible) merge sources
358 if(dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
359 curSDigit = dynamic_cast<AliEMCALDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
362 //May be several digits will contribute from the same input
363 while(curSDigit && (curSDigit->GetId() == absID)){
364 //Shift primary to separate primaries belonging different inputs
365 Int_t primaryoffset ;
367 primaryoffset = fManager->GetMask(i) ;
370 curSDigit->ShiftPrimary(primaryoffset) ;
372 a = curSDigit->GetAmp() ;
373 b = a /fTimeSignalLength ;
374 new((*ticks)[contrib++]) AliEMCALTick(curSDigit->GetTime(),0, b);
375 new((*ticks)[contrib++]) AliEMCALTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b);
377 *digit = *digit + *curSDigit ; //add energies
380 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
381 curSDigit = dynamic_cast<AliEMCALDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
386 // add fluctuations for photo-electron creation
387 amp = sDigitizer->Calibrate(digit->GetAmp()) ; // GeV
388 amp *= static_cast<Float_t>(gRandom->Poisson(fMeanPhotonElectron)) / static_cast<Float_t>(fMeanPhotonElectron) ;
390 //calculate and set time
391 Float_t time = FrontEdgeTime(ticks) ;
392 digit->SetTime(time) ;
394 //Find next signal module
396 for(i = 0 ; i < fInput ; i++){
397 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
398 Int_t curNext = nextSig ;
399 if(sdigits->GetEntriesFast() > index[i] ){
400 curNext = dynamic_cast<AliEMCALDigit *>(sdigits->At(index[i]))->GetId() ;
402 if(curNext < nextSig) nextSig = curNext ;
406 amp += TMath::Abs(gRandom->Gaus(0., fPinNoise)) ;
407 digit->SetAmp(sDigitizer->Digitize(amp)) ;
408 AliDebug(10,Form(" absID %5i amp %f nextSig %5i\n",
409 absID, amp, nextSig));
410 } // for(absID = 1; absID <= nEMC; absID++)
415 delete sdigArray ; //We should not delete its contents
417 //remove digits below thresholds
418 for(i = 0 ; i < nEMC ; i++){
419 digit = dynamic_cast<AliEMCALDigit*>( digits->At(i) ) ;
420 Float_t threshold = fDigitThreshold ;
421 if(sDigitizer->Calibrate( digit->GetAmp() ) < threshold)
422 digits->RemoveAt(i) ;
424 digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ;
429 Int_t ndigits = digits->GetEntriesFast() ;
431 //Set indexes in list of digits and fill hists.
432 AliEMCALHistoUtilities::FillH1(fHists, 0, Double_t(ndigits));
433 Float_t energy=0., esum=0.;
434 for (i = 0 ; i < ndigits ; i++) {
435 digit = dynamic_cast<AliEMCALDigit *>( digits->At(i) ) ;
436 digit->SetIndexInList(i) ;
437 energy = sDigitizer->Calibrate(digit->GetAmp()) ;
439 digit->SetAmp(DigitizeEnergy(energy, digit->GetId()) ) ; // for what ??
440 AliEMCALHistoUtilities::FillH1(fHists, 2, double(digit->GetAmp()));
441 AliEMCALHistoUtilities::FillH1(fHists, 3, double(energy));
442 AliEMCALHistoUtilities::FillH1(fHists, 4, double(digit->GetId()));
443 // if(digit->GetId() == nEMC) {
444 // printf(" i %i \n", i );
449 AliEMCALHistoUtilities::FillH1(fHists, 1, esum);
452 // //_____________________________________________________________________
453 Int_t AliEMCALDigitizer::DigitizeEnergy(Float_t energy, Int_t AbsId)
455 // Returns digitized value of the energy in a cell absId
458 const AliEMCALGeometry * geom = AliEMCALGeometry::GetInstance();
461 AliFatal("Did not get geometry from EMCALLoader");
469 Int_t channel = -999;
471 Bool_t bCell = geom->GetCellIndex(AbsId, iSupMod, nModule, nIphi, nIeta) ;
473 Error("DigitizeEnergy","Wrong cell id number : AbsId %i ", AbsId) ;
474 geom->GetCellPhiEtaIndexInSModule(iSupMod,nModule,nIphi, nIeta,iphi,ieta);
477 fADCpedestalEC = fCalibData->GetADCpedestal(iSupMod,ieta,iphi);
478 fADCchannelEC = fCalibData->GetADCchannel(iSupMod,ieta,iphi);
481 channel = static_cast<Int_t> (TMath::Ceil( (energy + fADCpedestalEC)/fADCchannelEC )) ;
483 if(channel > fNADCEC )
489 //____________________________________________________________________________
490 void AliEMCALDigitizer::Exec(Option_t *option)
492 // Steering method to process digitization for events
493 // in the range from fFirstEvent to fLastEvent.
494 // This range is optionally set by SetEventRange().
495 // if fLastEvent=-1, then process events until the end.
496 // by default fLastEvent = fFirstEvent (process only one event)
498 if (!fInit) { // to prevent overwrite existing file
499 Error( "Exec", "Give a version name different from %s", fEventFolderName.Data() ) ;
503 if (strstr(option,"print")) {
509 if(strstr(option,"tim"))
510 gBenchmark->Start("EMCALDigitizer");
512 AliRunLoader *rl = AliRunLoader::GetRunLoader();
513 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(rl->GetDetectorLoader("EMCAL"));
515 // Post Digitizer to the white board
516 emcalLoader->PostDigitizer(this) ;
518 if (fLastEvent == -1) {
519 fLastEvent = rl->GetNumberOfEvents() - 1 ;
522 fLastEvent = fFirstEvent ; // what is this ??
524 Int_t nEvents = fLastEvent - fFirstEvent + 1;
527 rl->LoadSDigits("EMCAL");
528 for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
530 rl->GetEvent(ievent);
532 Digitize(ievent) ; //Add prepared SDigits to digits and add the noise
536 if(strstr(option,"deb"))
538 if(strstr(option,"table")) gObjectTable->Print();
540 //increment the total number of Digits per run
541 fDigitsInRun += emcalLoader->Digits()->GetEntriesFast() ;
544 emcalLoader->CleanDigitizer() ;
546 if(strstr(option,"tim")){
547 gBenchmark->Stop("EMCALDigitizer");
548 AliInfo(Form("Exec: took %f seconds for Digitizing %f seconds per event",
549 gBenchmark->GetCpuTime("EMCALDigitizer"), gBenchmark->GetCpuTime("EMCALDigitizer")/nEvents )) ;
553 //____________________________________________________________________________
554 Float_t AliEMCALDigitizer::FrontEdgeTime(TClonesArray * ticks)
556 // Returns the shortest time among all time ticks
558 ticks->Sort() ; //Sort in accordance with times of ticks
560 AliEMCALTick * ctick = (AliEMCALTick *) it.Next() ;
561 Float_t time = ctick->CrossingTime(fTimeThreshold) ;
564 while((t=(AliEMCALTick*) it.Next())){
565 if(t->GetTime() < time) //This tick starts before crossing
570 time = ctick->CrossingTime(fTimeThreshold) ;
575 //____________________________________________________________________________
576 Bool_t AliEMCALDigitizer::Init()
578 // Makes all memory allocations
580 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
582 if ( emcalLoader == 0 ) {
583 Fatal("Init", "Could not obtain the AliEMCALLoader");
588 fLastEvent = fFirstEvent ;
591 fInput = fManager->GetNinputs() ;
595 fInputFileNames = new TString[fInput] ;
596 fEventNames = new TString[fInput] ;
597 fInputFileNames[0] = GetTitle() ;
598 fEventNames[0] = fEventFolderName.Data() ;
600 for (index = 1 ; index < fInput ; index++) {
601 fInputFileNames[index] = dynamic_cast<AliStream*>(fManager->GetInputStream(index))->GetFileName(0);
602 TString tempo = fManager->GetInputFolderName(index) ;
603 fEventNames[index] = tempo.Remove(tempo.Length()-1) ; // strip of the stream number added bt fManager
606 //to prevent cleaning of this object while GetEvent is called
607 emcalLoader->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
610 //Calibration instance
611 fCalibData = emcalLoader->CalibData();
615 //____________________________________________________________________________
616 void AliEMCALDigitizer::InitParameters()
618 // Parameter initialization for digitizer
619 // Tune parameters - 24-nov-04; Apr 29, 2007
621 fMeanPhotonElectron = 3300; // electrons per GeV
622 fPinNoise = 0.010; // pin noise in GEV from analysis test beam data
623 if (fPinNoise == 0. )
624 Warning("InitParameters", "No noise added\n") ;
625 fDigitThreshold = fPinNoise * 3; // 3 * sigma
626 fTimeResolution = 0.3e-9 ; // 300 psc
627 fTimeSignalLength = 1.0e-9 ;
629 // These defaults are normally not used.
630 // Values are read from calibration database instead
631 fADCchannelEC = 0.0153; // Update 24 Apr 2007: 250./16/1024 - width of one ADC channel in GeV
632 fADCpedestalEC = 0.0 ; // GeV
633 fNADCEC = (Int_t) TMath::Power(2,16) ; // number of channels in Tower ADC - 65536
635 fTimeThreshold = 0.001*10000000 ; // Means 1 MeV in terms of SDigits amplitude ??
636 // hists. for control; no hists on default
641 //__________________________________________________________________
642 void AliEMCALDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
644 // Allows to produce digits by superimposing background and signal event.
645 // It is assumed, that headers file with SIGNAL events is opened in
647 // Sets the BACKGROUND event, with which the SIGNAL event is to be mixed
648 // Thus we avoid writing (changing) huge and expensive
649 // backgound files: all output will be writen into SIGNAL, i.e.
650 // opened in constructor file.
652 // One can open as many files to mix with as one needs.
653 // However only Sdigits with the same name (i.e. constructed with the same SDigitizer)
656 if( strcmp(GetName(), "") == 0 )
660 Error("MixWith", "Cannot use this method under AliRunDigitizer") ;
663 // looking for file which contains AliRun
664 if (gSystem->AccessPathName(alirunFileName)) {// file does not exist
665 Error("MixWith", "File %s does not exist!", alirunFileName.Data()) ;
668 // looking for the file which contains SDigits
669 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
670 TString fileName( emcalLoader->GetSDigitsFileName() ) ;
671 if ( eventFolderName != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
672 fileName = fileName.ReplaceAll(".root", "") + "_" + eventFolderName + ".root" ;
673 if ( (gSystem->AccessPathName(fileName)) ) {
674 Error("MixWith", "The file %s does not exist!", fileName.Data()) ;
677 // need to increase the arrays
678 // MvL: This code only works when fInput == 1, I think.
679 TString tempo = fInputFileNames[fInput-1] ;
680 delete [] fInputFileNames ;
681 fInputFileNames = new TString[fInput+1] ;
682 fInputFileNames[fInput-1] = tempo ;
684 tempo = fEventNames[fInput-1] ;
685 delete [] fEventNames ;
686 fEventNames = new TString[fInput+1] ;
687 fEventNames[fInput-1] = tempo ;
689 fInputFileNames[fInput] = alirunFileName ;
690 fEventNames[fInput] = eventFolderName ;
694 void AliEMCALDigitizer::Print1(Option_t * option)
695 { // 19-nov-04 - just for convinience
700 //__________________________________________________________________
701 void AliEMCALDigitizer::Print(Option_t*)const
703 // Print Digitizer's parameters
704 printf("Print: \n------------------- %s -------------", GetName() ) ;
705 if( strcmp(fEventFolderName.Data(), "") != 0 ){
706 printf(" Writing Digits to branch with title %s\n", fEventFolderName.Data()) ;
710 nStreams = GetNInputStreams() ;
718 for (index = 0 ; index < nStreams ; index++) {
719 TString tempo(fEventNames[index]) ;
721 if ((rl = AliRunLoader::GetRunLoader(tempo)) == 0)
722 rl = AliRunLoader::Open(fInputFileNames[index], tempo) ;
723 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(rl->GetDetectorLoader("EMCAL"));
724 TString fileName( emcalLoader->GetSDigitsFileName() ) ;
725 if ( fEventNames[index] != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
726 fileName = fileName.ReplaceAll(".root", "") + "_" + fEventNames[index] + ".root" ;
727 printf ("Adding SDigits from %s %s\n", fInputFileNames[index].Data(), fileName.Data()) ;
730 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
732 printf("\nWriting digits to %s", emcalLoader->GetDigitsFileName().Data()) ;
734 printf("\nWith following parameters:\n") ;
736 printf(" Electronics noise in EMC (fPinNoise) = %f\n", fPinNoise) ;
737 printf(" Threshold in EMC (fDigitThreshold) = %f\n", fDigitThreshold) ;
738 printf("---------------------------------------------------\n") ;
741 printf("Print: AliEMCALDigitizer not initialized") ;
744 //__________________________________________________________________
745 void AliEMCALDigitizer::PrintDigits(Option_t * option)
747 //utility method for printing digit information
749 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
750 TClonesArray * digits = emcalLoader->Digits() ;
751 TClonesArray * sdigits = emcalLoader->SDigits() ;
753 printf("\n #Digits: %d : sdigits %d ", digits->GetEntriesFast(), sdigits->GetEntriesFast()) ;
754 printf("\n event %d", emcalLoader->GetRunLoader()->GetEventNumber());
756 if(strstr(option,"all")){
758 AliEMCALDigit * digit;
759 printf("\nEMC digits (with primaries):\n") ;
760 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
762 for (index = 0 ; index < digits->GetEntries() ; index++) {
763 digit = dynamic_cast<AliEMCALDigit *>(digits->At(index)) ;
764 printf("\n%6d %8d %6.5e %4d %2d : ",
765 digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
767 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
768 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
775 //__________________________________________________________________
776 Float_t AliEMCALDigitizer::TimeOfNoise(void)
778 // Calculates the time signal generated by noise
779 //PH Info("TimeOfNoise", "Change me") ;
780 return gRandom->Rndm() * 1.28E-5;
783 //__________________________________________________________________
784 void AliEMCALDigitizer::Unload()
786 // Unloads the SDigits and Digits
790 for(i = 1 ; i < fInput ; i++){
791 TString tempo(fEventNames[i]) ;
793 if ((rl = AliRunLoader::GetRunLoader(tempo)))
794 rl->GetDetectorLoader("EMCAL")->UnloadSDigits() ;
796 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
797 emcalLoader->UnloadDigits() ;
800 //_________________________________________________________________________________________
801 void AliEMCALDigitizer::WriteDigits()
804 // Makes TreeD in the output file.
805 // Check if branch already exists:
806 // if yes, exit without writing: ROOT TTree does not support overwriting/updating of
807 // already existing branches.
808 // else creates branch with Digits, named "EMCAL", title "...",
809 // and branch "AliEMCALDigitizer", with the same title to keep all the parameters
810 // and names of files, from which digits are made.
812 AliEMCALLoader *emcalLoader = dynamic_cast<AliEMCALLoader*>(AliRunLoader::GetRunLoader()->GetDetectorLoader("EMCAL"));
814 const TClonesArray * digits = emcalLoader->Digits() ;
815 TTree * treeD = emcalLoader->TreeD();
817 emcalLoader->MakeDigitsContainer();
818 treeD = emcalLoader->TreeD();
821 // -- create Digits branch
822 Int_t bufferSize = 32000 ;
823 TBranch * digitsBranch = 0;
824 if ((digitsBranch = treeD->GetBranch("EMCAL"))) {
825 digitsBranch->SetAddress(&digits);
826 AliWarning("Digits Branch already exists. Not all digits will be visible");
829 treeD->Branch("EMCAL","TClonesArray",&digits,bufferSize);
830 //digitsBranch->SetTitle(fEventFolderName);
833 emcalLoader->WriteDigits("OVERWRITE");
834 emcalLoader->WriteDigitizer("OVERWRITE");
840 void AliEMCALDigitizer::Browse(TBrowser* b)
842 if(fHists) b->Add(fHists);
846 TList *AliEMCALDigitizer::BookControlHists(int var)
849 // histograms for monitoring digitizer performance
851 Info("BookControlHists"," started ");
853 const AliEMCALGeometry *geom = AliEMCALGeometry::GetInstance();
855 new TH1F("hDigiN", "#EMCAL digits with fAmp > fDigitThreshold",
856 fNADCEC+1, -0.5, Double_t(fNADCEC));
857 new TH1F("HDigiSumEnergy","Sum.EMCAL energy from digi", 1000, 0.0, 200.);
858 new TH1F("hDigiAmp", "EMCAL digital amplitude", fNADCEC+1, -0.5, Double_t(fNADCEC));
859 new TH1F("hDigiEnergy","EMCAL cell energy", 2000, 0.0, 200.);
860 new TH1F("hDigiAbsId","EMCAL absId cells with fAmp > fDigitThreshold ",
861 geom->GetNCells(), 0.5, Double_t(geom->GetNCells())+0.5);
864 fHists = AliEMCALHistoUtilities::MoveHistsToList("EmcalDigiControlHists", kFALSE);
865 fHists = 0; //huh? JLK 03-Mar-2006
869 void AliEMCALDigitizer::SaveHists(const char* name, Bool_t kSingleKey, const char* opt)
871 AliEMCALHistoUtilities::SaveListOfHists(fHists, name, kSingleKey, opt);