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 /* History of cvs commits:
21 * Revision 1.96 2007/04/28 10:43:36 policheh
22 * Dead channels simulation: digit energy sets to 0.
24 * Revision 1.95 2007/04/10 07:20:52 kharlov
25 * Decalibration should use the same CDB as calibration in AliPHOSClusterizerv1
27 * Revision 1.94 2007/02/01 10:34:47 hristov
28 * Removing warnings on Solaris x86
30 * Revision 1.93 2006/10/17 13:17:01 kharlov
31 * Replace AliInfo by AliDebug
33 * Revision 1.92 2006/08/28 10:01:56 kharlov
34 * Effective C++ warnings fixed (Timur Pocheptsov)
36 * Revision 1.91 2006/04/29 20:25:30 hristov
37 * Decalibration is implemented (Yu.Kharlov)
39 * Revision 1.90 2006/04/22 10:30:17 hristov
40 * Add fEnergy to AliPHOSDigit and operate with EMC amplitude in energy units (Yu.Kharlov)
42 * Revision 1.89 2006/04/11 15:22:59 hristov
43 * run number in query set to -1: forces AliCDBManager to use its run number (A.Colla)
45 * Revision 1.88 2006/03/13 14:05:43 kharlov
46 * Calibration objects for EMC and CPV
48 * Revision 1.87 2005/08/24 15:33:49 kharlov
49 * Calibration data for raw digits
51 * Revision 1.86 2005/07/12 20:07:35 hristov
52 * Changes needed to run simulation and reconstrruction in the same AliRoot session
54 * Revision 1.85 2005/05/28 14:19:04 schutz
55 * Compilation warnings fixed by T.P.
59 //_________________________________________________________________________
60 //*-- Author : Dmitri Peressounko (SUBATECH & Kurchatov Institute)
61 //////////////////////////////////////////////////////////////////////////////
62 // This TTask performs digitization of Summable digits (in the PHOS case it is just
63 // the sum of contributions from all primary particles into a given cell).
64 // In addition it performs mixing of summable digits from different events.
65 // The name of the TTask is also the title of the branch that will contain
66 // the created SDigits
67 // The title of the TTAsk is the name of the file that contains the hits from
68 // which the SDigits are created
70 // For each event two branches are created in TreeD:
71 // "PHOS" - list of digits
72 // "AliPHOSDigitizer" - AliPHOSDigitizer with all parameters used in digitization
74 // Note, that one can set a title for new digits branch, and repeat digitization with
75 // another set of parameters.
78 // root[0] AliPHOSDigitizer * d = new AliPHOSDigitizer() ;
79 // root[1] d->ExecuteTask()
80 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
81 // //Digitizes SDigitis in all events found in file galice.root
83 // root[2] AliPHOSDigitizer * d1 = new AliPHOSDigitizer("galice1.root") ;
84 // // Will read sdigits from galice1.root
85 // root[3] d1->MixWith("galice2.root")
86 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
87 // // Reads another set of sdigits from galice2.root
88 // root[3] d1->MixWith("galice3.root")
89 // // Reads another set of sdigits from galice3.root
90 // root[4] d->ExecuteTask("deb timing")
91 // // Reads SDigits from files galice1.root, galice2.root ....
92 // // mixes them and stores produced Digits in file galice1.root
93 // // deb - prints number of produced digits
94 // // deb all - prints list of produced digits
95 // // timing - prints time used for digitization
98 // --- ROOT system ---
101 #include "TBenchmark.h"
104 // --- Standard library ---
106 // --- AliRoot header files ---
108 #include "AliRunDigitizer.h"
109 #include "AliPHOSDigit.h"
110 #include "AliPHOSGetter.h"
111 #include "AliPHOSDigitizer.h"
112 #include "AliPHOSSDigitizer.h"
113 #include "AliPHOSGeometry.h"
114 #include "AliPHOSTick.h"
115 #include "AliPHOSQualAssDataMaker.h"
117 ClassImp(AliPHOSDigitizer)
120 //____________________________________________________________________________
121 AliPHOSDigitizer::AliPHOSDigitizer() :
127 fInputFileNames(0x0),
131 fEMCDigitThreshold(0.f),
133 fCPVDigitThreshold(0.f),
134 fTimeResolution(0.f),
136 fTimeSignalLength(0.f),
138 fADCpedestalEmc(0.f),
141 fADCpedestalCpv(0.f),
143 fEventFolderName(""),
151 fManager = 0 ; // We work in the standalong mode
154 //____________________________________________________________________________
155 AliPHOSDigitizer::AliPHOSDigitizer(TString alirunFileName,
156 TString eventFolderName):
157 AliDigitizer("PHOS"+AliConfig::Instance()->GetDigitizerTaskName(), alirunFileName),
158 fDefaultInit(kFALSE),
162 fInputFileNames(0x0),
166 fEMCDigitThreshold(0.f),
168 fCPVDigitThreshold(0.f),
169 fTimeResolution(0.f),
171 fTimeSignalLength(0.f),
173 fADCpedestalEmc(0.f),
176 fADCpedestalCpv(0.f),
178 fEventFolderName(eventFolderName),
187 fDefaultInit = kFALSE ;
188 fManager = 0 ; // We work in the standalong mode
189 //Initialize the quality assurance data maker only once
190 fQADM = new AliPHOSQualAssDataMaker() ;
191 GetQualAssDataMaker()->Init(AliQualAss::kDIGITS) ;
194 //____________________________________________________________________________
195 AliPHOSDigitizer::AliPHOSDigitizer(const AliPHOSDigitizer & d) :
197 fDefaultInit(d.fDefaultInit),
198 fDigitsInRun(d.fDigitsInRun),
201 fInputFileNames(0x0),//?
203 fEmcCrystals(d.fEmcCrystals),
204 fPinNoise(d.fPinNoise),
205 fEMCDigitThreshold(d.fEMCDigitThreshold),
206 fCPVNoise(d.fCPVNoise),
207 fCPVDigitThreshold(d.fCPVDigitThreshold),
208 fTimeResolution(d.fTimeResolution),
209 fTimeThreshold(d.fTimeThreshold),
210 fTimeSignalLength(d.fTimeSignalLength),
211 fADCchanelEmc(d.fADCchanelEmc),
212 fADCpedestalEmc(d.fADCpedestalEmc),
213 fNADCemc(d.fNADCemc),
214 fADCchanelCpv(d.fADCchanelCpv),
215 fADCpedestalCpv(d.fADCpedestalCpv),
216 fNADCcpv(d.fNADCcpv),
217 fEventFolderName(d.fEventFolderName),
218 fFirstEvent(d.fFirstEvent),
219 fLastEvent(d.fLastEvent),
225 SetName(d.GetName()) ;
226 SetTitle(d.GetTitle()) ;
227 //Initialize the quality assurance data maker only once
228 GetQualAssDataMaker()->Init(AliQualAss::kDIGITS) ;
231 //____________________________________________________________________________
232 AliPHOSDigitizer::AliPHOSDigitizer(AliRunDigitizer * rd) :
233 AliDigitizer(rd,"PHOS"+AliConfig::Instance()->GetDigitizerTaskName()),
234 fDefaultInit(kFALSE),
238 fInputFileNames(0x0),
242 fEMCDigitThreshold(0.f),
244 fCPVDigitThreshold(0.f),
245 fTimeResolution(0.f),
247 fTimeSignalLength(0.f),
249 fADCpedestalEmc(0.f),
252 fADCpedestalCpv(0.f),
254 fEventFolderName(fManager->GetInputFolderName(0)),
261 // ctor Init() is called by RunDigitizer
263 SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
265 fDefaultInit = kFALSE ;
266 //Initialize the quality assurance data maker only once
267 fQADM = new AliPHOSQualAssDataMaker() ;
268 GetQualAssDataMaker()->Init(AliQualAss::kDIGITS) ;
271 //____________________________________________________________________________
272 AliPHOSDigitizer::~AliPHOSDigitizer()
274 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
276 // Clean Digitizer from the white board
277 gime->PhosLoader()->CleanDigitizer() ;
279 delete [] fInputFileNames ;
280 delete [] fEventNames ;
286 //____________________________________________________________________________
287 void AliPHOSDigitizer::Digitize(Int_t event)
290 // Makes the digitization of the collected summable digits.
291 // It first creates the array of all PHOS modules
292 // filled with noise (different for EMC, and CPV) and
293 // then adds contributions from SDigits.
294 // This design avoids scanning over the list of digits to add
295 // contribution to new SDigits only.
297 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
298 Int_t ReadEvent = event ;
300 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetCurrentEventNumber() ;
301 AliDebug(1,Form("Adding event %d from input stream 0 %s %s",
302 ReadEvent, GetTitle(), fEventFolderName.Data())) ;
303 gime->Event(ReadEvent, "S") ;
304 TClonesArray * digits = gime->Digits() ;
307 const AliPHOSGeometry *geom = gime->PHOSGeometry() ;
308 //Making digits with noise, first EMC
309 Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
314 nCPV = nEMC + geom->GetNumberOfCPVPadsZ() * geom->GetNumberOfCPVPadsPhi() * geom->GetNModules() ;
316 digits->Expand(nCPV) ;
318 // get first the sdigitizer from the tasks list
319 if ( !gime->SDigitizer() )
320 gime->LoadSDigitizer();
321 AliPHOSSDigitizer * sDigitizer = gime->SDigitizer();
324 AliFatal(Form("SDigitizer with name %s %s not found",
325 GetTitle(), fEventFolderName.Data() )) ;
327 //take all the inputs to add together and load the SDigits
328 TObjArray * sdigArray = new TObjArray(fInput) ;
329 sdigArray->AddAt(gime->SDigits(), 0) ;
331 for(i = 1 ; i < fInput ; i++){
332 TString tempo(fEventNames[i]) ;
334 AliPHOSGetter * gime1 = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
336 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(i))->GetCurrentEventNumber() ;
337 AliInfo(Form("Adding event %d from input stream %d %s %s",
338 ReadEvent, i, fInputFileNames[i].Data(), tempo.Data())) ;
339 gime1->Event(ReadEvent,"S");
340 sdigArray->AddAt(gime1->SDigits(), i) ;
343 //Find the first crystal with signal
344 Int_t nextSig = 200000 ;
345 TClonesArray * sdigits ;
346 for(i = 0 ; i < fInput ; i++){
347 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
348 if ( !sdigits->GetEntriesFast() )
350 Int_t curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(0))->GetId() ;
351 if(curNext < nextSig)
355 TArrayI index(fInput) ;
356 index.Reset() ; //Set all indexes to zero
358 AliPHOSDigit * digit ;
359 AliPHOSDigit * curSDigit ;
361 TClonesArray * ticks = new TClonesArray("AliPHOSTick",1000) ;
363 //Put Noise contribution
364 for(absID = 1 ; absID <= nEMC ; absID++){
365 Float_t noise = gRandom->Gaus(0., fPinNoise) ;
366 // YVK: do not digitize amplitudes for EMC
367 // new((*digits)[absID-1]) AliPHOSDigit( -1, absID, sDigitizer->Digitize(noise), TimeOfNoise() ) ;
368 new((*digits)[absID-1]) AliPHOSDigit( -1, absID, noise, TimeOfNoise() ) ;
369 //look if we have to add signal?
370 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
373 //Add SDigits from all inputs
376 Float_t a = digit->GetEnergy() ;
377 Float_t b = TMath::Abs( a / fTimeSignalLength) ;
378 //Mark the beginning of the signal
379 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime(),0, b);
380 //Mark the end of the signal
381 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime()+fTimeSignalLength, -a, -b);
384 for(i = 0 ; i < fInput ; i++){
385 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
386 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
389 //May be several digits will contribute from the same input
390 while(curSDigit && curSDigit->GetId() == absID){
391 //Shift primary to separate primaries belonging different inputs
392 Int_t primaryoffset ;
394 primaryoffset = fManager->GetMask(i) ;
396 primaryoffset = 10000000*i ;
397 curSDigit->ShiftPrimary(primaryoffset) ;
399 a = curSDigit->GetEnergy() ;
400 b = a /fTimeSignalLength ;
401 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime(),0, b);
402 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b);
404 *digit += *curSDigit ; //add energies
407 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
408 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
414 //calculate and set time
415 Float_t time = FrontEdgeTime(ticks) ;
416 digit->SetTime(time) ;
418 //Find next signal module
420 for(i = 0 ; i < fInput ; i++){
421 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
422 Int_t curNext = nextSig ;
423 if(sdigits->GetEntriesFast() > index[i] ){
424 curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(index[i]))->GetId() ;
426 if(curNext < nextSig) nextSig = curNext ;
434 //Now CPV digits (different noise and no timing)
435 for(absID = nEMC+1; absID <= nCPV; absID++){
436 Float_t noise = gRandom->Gaus(0., fCPVNoise) ;
437 new((*digits)[absID-1]) AliPHOSDigit( -1,absID,sDigitizer->Digitize(noise), TimeOfNoise() ) ;
438 //look if we have to add signal?
440 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
441 //Add SDigits from all inputs
442 for(i = 0 ; i < fInput ; i++){
443 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
444 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
448 //May be several digits will contribute from the same input
449 while(curSDigit && curSDigit->GetId() == absID){
450 //Shift primary to separate primaries belonging different inputs
451 Int_t primaryoffset ;
453 primaryoffset = fManager->GetMask(i) ;
455 primaryoffset = 10000000*i ;
456 curSDigit->ShiftPrimary(primaryoffset) ;
459 *digit += *curSDigit ;
461 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
462 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i]) ) ;
468 //Find next signal module
470 for(i = 0 ; i < fInput ; i++){
471 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
472 Int_t curNext = nextSig ;
473 if(sdigits->GetEntriesFast() > index[i] )
474 curNext = dynamic_cast<AliPHOSDigit *>( sdigits->At(index[i]) )->GetId() ;
475 if(curNext < nextSig) nextSig = curNext ;
481 delete sdigArray ; //We should not delete its contents
483 //remove digits below thresholds
484 for(i = 0 ; i < nEMC ; i++){
485 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
486 DecalibrateEMC(digit);
487 if(digit->GetEnergy() < fEMCDigitThreshold)
488 digits->RemoveAt(i) ;
490 digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ;
494 for(i = nEMC; i < nCPV ; i++)
495 // if( sDigitizer->Calibrate( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetAmp() ) < fCPVDigitThreshold )
496 if( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetEnergy() < fCPVDigitThreshold )
497 digits->RemoveAt(i) ;
501 Int_t ndigits = digits->GetEntriesFast() ;
502 digits->Expand(ndigits) ;
504 //Set indexes in list of digits and make true digitization of the energy
505 for (i = 0 ; i < ndigits ; i++) {
506 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
507 digit->SetIndexInList(i) ;
508 if(digit->GetId() > fEmcCrystals){ //digitize CPV only
509 digit->SetAmp(DigitizeCPV(digit->GetEnergy(),digit->GetId()) ) ;
515 //set amplitudes in bad channels to zero
516 for(i = 0 ; i <digits->GetEntries(); i++){
517 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
518 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
519 if(relId[1] == 0) // Emc
520 if(gime->CalibData()->IsBadChannelEmc(relId[0],relId[3],relId[2])) digit->SetEnergy(0.);
525 //____________________________________________________________________________
526 void AliPHOSDigitizer::DecalibrateEMC(AliPHOSDigit *digit)
528 // Decalibrate EMC digit, i.e. change its energy by a factor read from CDB
530 AliPHOSGetter* gime = AliPHOSGetter::Instance();
532 if(!gime->CalibData()) {
533 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1);
534 gime->SetCalibData(cdb);
537 //Determine rel.position of the cell absolute ID
539 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
540 Int_t module=relId[0];
542 Int_t column=relId[3];
543 Float_t decalibration = gime->CalibData()->GetADCchannelEmc(module,column,row);
544 Float_t energy = digit->GetEnergy() / decalibration;
545 digit->SetEnergy(energy);
547 //____________________________________________________________________________
548 Int_t AliPHOSDigitizer::DigitizeCPV(Float_t charge, Int_t absId)
550 // Returns digitized value of the CPV charge in a pad absId
552 AliPHOSGetter* gime = AliPHOSGetter::Instance();
554 if(!gime->CalibData()) {
555 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1); // use AliCDBManager's run number
556 gime->SetCalibData(cdb);
559 //Determine rel.position of the cell absId
561 gime->PHOSGeometry()->AbsToRelNumbering(absId,relId);
562 Int_t module=relId[0];
564 Int_t column=relId[3];
568 if(absId > fEmcCrystals){ //digitize CPV only
570 //reading calibration data for cell absId.
571 //If no calibration DB found, accept default values.
573 if(gime->CalibData()) {
574 fADCpedestalCpv = gime->CalibData()->GetADCpedestalCpv(module,column,row);
575 fADCchanelCpv = gime->CalibData()->GetADCchannelCpv( module,column,row);
578 channel = (Int_t) TMath::Ceil((charge - fADCpedestalCpv)/fADCchanelCpv) ;
579 if(channel > fNADCcpv ) channel = fNADCcpv ;
584 //____________________________________________________________________________
585 void AliPHOSDigitizer::Exec(Option_t *option)
587 // Steering method to process digitization for events
588 // in the range from fFirstEvent to fLastEvent.
589 // This range is optionally set by SetEventRange().
590 // if fLastEvent=-1, then process events until the end.
591 // by default fLastEvent = fFirstEvent (process only one event)
593 if (!fInit) { // to prevent overwrite existing file
594 AliError(Form("Give a version name different from %s",
595 fEventFolderName.Data() )) ;
599 if (strstr(option,"print")) {
604 if(strstr(option,"tim"))
605 gBenchmark->Start("PHOSDigitizer");
607 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
609 // Post Digitizer to the white board
610 gime->PostDigitizer(this) ;
612 if (fLastEvent == -1)
613 fLastEvent = gime->MaxEvent() - 1 ;
615 fLastEvent = fFirstEvent ;
617 Int_t nEvents = fLastEvent - fFirstEvent + 1;
621 for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
623 gime->Event(ievent,"S") ;
625 Digitize(ievent) ; //Add prepared SDigits to digits and add the noise
627 //makes the quality assurance data
628 GetQualAssDataMaker()->SetData(gime->Digits()) ;
629 GetQualAssDataMaker()->Exec(AliQualAss::kDIGITS) ;
633 if(strstr(option,"deb"))
636 //increment the total number of Digits per run
637 fDigitsInRun += gime->Digits()->GetEntriesFast() ;
640 //Write the quality assurance data only after the last event
641 if ( fEventCounter == gime->MaxEvent() )
642 GetQualAssDataMaker()->Finish(AliQualAss::kDIGITS) ;
644 gime->PhosLoader()->CleanDigitizer();
646 if(strstr(option,"tim")){
647 gBenchmark->Stop("PHOSDigitizer");
649 message = " took %f seconds for Digitizing %f seconds per event\n" ;
650 AliInfo(Form( message.Data(),
651 gBenchmark->GetCpuTime("PHOSDigitizer"),
652 gBenchmark->GetCpuTime("PHOSDigitizer")/nEvents ));
656 //____________________________________________________________________________
657 Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks) const
659 // Returns the shortest time among all time ticks
661 ticks->Sort() ; //Sort in accordance with times of ticks
663 AliPHOSTick * ctick = (AliPHOSTick *) it.Next() ;
664 Float_t time = ctick->CrossingTime(fTimeThreshold) ;
667 while((t=(AliPHOSTick*) it.Next())){
668 if(t->GetTime() < time) //This tick starts before crossing
673 time = ctick->CrossingTime(fTimeThreshold) ;
678 //____________________________________________________________________________
679 Bool_t AliPHOSDigitizer::Init()
681 // Makes all memory allocations
683 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
685 AliFatal(Form("Could not obtain the Getter object for file %s and event %s !",
686 GetTitle(), fEventFolderName.Data()));
690 const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
692 fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ;
694 TString opt("Digits") ;
695 if(gime->VersionExists(opt) ) {
696 AliError(Form("Give a version name different from %s",
697 fEventFolderName.Data() )) ;
702 fLastEvent = fFirstEvent ;
704 fInput = fManager->GetNinputs() ;
708 fInputFileNames = new TString[fInput] ;
709 fEventNames = new TString[fInput] ;
710 fInputFileNames[0] = GetTitle() ;
711 fEventNames[0] = fEventFolderName.Data() ;
713 for (index = 1 ; index < fInput ; index++) {
714 fInputFileNames[index] = dynamic_cast<AliStream*>(fManager->GetInputStream(index))->GetFileName(0);
715 TString tempo = fManager->GetInputFolderName(index) ;
716 fEventNames[index] = tempo.Remove(tempo.Length()-1) ; // strip of the stream number added by fManager
719 //to prevent cleaning of this object while GetEvent is called
720 gime->PhosLoader()->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
725 //____________________________________________________________________________
726 void AliPHOSDigitizer::InitParameters()
728 // Set initial parameters Digitizer
730 fPinNoise = 0.004 ; // [GeV]
731 fEMCDigitThreshold = 0.012 ; // [GeV]
732 fCPVNoise = 0.01; // [aux units]
733 fCPVDigitThreshold = 0.09 ; // [aux units]
734 fTimeResolution = 0.5e-9 ; // [sec]
735 fTimeSignalLength = 1.0e-9 ; // [sec]
737 fADCchanelEmc = 1.0; // Coefficient between real and measured energies in EMC
738 fADCpedestalEmc = 0. ; //
739 fNADCemc = (Int_t) TMath::Power(2,16) ; // number of channels in EMC ADC
741 fADCchanelCpv = 0.0012 ; // width of one ADC channel in CPV 'popugais'
742 fADCpedestalCpv = 0.012 ; //
743 fNADCcpv = (Int_t) TMath::Power(2,12); // number of channels in CPV ADC
745 // fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
746 fTimeThreshold = 0.001 ; // [GeV]
747 SetEventRange(0,-1) ;
751 //__________________________________________________________________
752 void AliPHOSDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
754 // Allows to produce digits by superimposing background and signal event.
755 // It is assumed, that headers file with SIGNAL events is opened in
757 // Sets the BACKGROUND event, with which the SIGNAL event is to be mixed
758 // Thus we avoid writing (changing) huge and expensive
759 // backgound files: all output will be writen into SIGNAL, i.e.
760 // opened in constructor file.
762 // One can open as many files to mix with as one needs.
763 // However only Sdigits with the same name (i.e. constructed with the same SDigitizer)
766 if( strcmp(fEventFolderName, "") == 0 )
770 Warning("MixWith", "Cannot use this method with AliRunDigitizer\n" ) ;
773 // looking for file which contains AliRun
774 if (gSystem->AccessPathName(alirunFileName)) {// file does not exist
775 AliError(Form("File %s does not exist!", alirunFileName.Data())) ;
778 // looking for the file which contains SDigits
779 AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
780 TString fileName( gime->GetSDigitsFileName() ) ;
781 if ( eventFolderName != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
782 fileName = fileName.ReplaceAll(".root", "") + "_" + eventFolderName + ".root" ;
783 if ( (gSystem->AccessPathName(fileName)) ) {
784 AliError(Form("The file %s does not exist!", fileName.Data())) ;
787 // need to increase the arrays
788 TString tempo = fInputFileNames[fInput-1] ;
789 delete [] fInputFileNames ;
790 fInputFileNames = new TString[fInput+1] ;
791 fInputFileNames[fInput-1] = tempo ;
793 tempo = fEventNames[fInput-1] ;
794 delete [] fEventNames ;
795 fEventNames = new TString[fInput+1] ;
796 fEventNames[fInput-1] = tempo ;
798 fInputFileNames[fInput] = alirunFileName ;
799 fEventNames[fInput] = eventFolderName ;
803 //__________________________________________________________________
804 void AliPHOSDigitizer::Print(const Option_t *)const
806 // Print Digitizer's parameters
807 AliInfo(Form("\n------------------- %s -------------", GetName() )) ;
808 if( strcmp(fEventFolderName.Data(), "") != 0 ){
809 printf(" Writing Digits to branch with title %s\n", fEventFolderName.Data()) ;
813 nStreams = GetNInputStreams() ;
818 for (index = 0 ; index < nStreams ; index++) {
819 TString tempo(fEventNames[index]) ;
821 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[index], tempo) ;
822 TString fileName( gime->GetSDigitsFileName() ) ;
823 if ( fEventNames[index] != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
824 fileName = fileName.ReplaceAll(".root", "") + "_" + fEventNames[index] + ".root" ;
825 printf ("Adding SDigits from %s %s\n", fInputFileNames[index].Data(), fileName.Data()) ;
827 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
828 printf("\nWriting digits to %s", gime->GetDigitsFileName().Data()) ;
830 printf("\nWith following parameters:\n") ;
831 printf(" Electronics noise in EMC (fPinNoise) = %f GeV\n", fPinNoise ) ;
832 printf(" Threshold in EMC (fEMCDigitThreshold) = %f GeV\n", fEMCDigitThreshold ) ;
833 printf(" Noise in CPV (fCPVNoise) = %f aux units\n", fCPVNoise ) ;
834 printf(" Threshold in CPV (fCPVDigitThreshold) = %f aux units\n",fCPVDigitThreshold ) ;
835 printf(" ---------------------------------------------------\n") ;
838 AliInfo(Form("AliPHOSDigitizer not initialized" )) ;
842 //__________________________________________________________________
843 void AliPHOSDigitizer::PrintDigits(Option_t * option)
845 // Print a table of digits
847 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
848 TClonesArray * digits = gime->Digits() ;
850 AliInfo(Form("%d", digits->GetEntriesFast())) ;
851 printf("\nevent %d", gAlice->GetEvNumber()) ;
852 printf("\n Number of entries in Digits list %d", digits->GetEntriesFast() ) ;
855 if(strstr(option,"all")||strstr(option,"EMC")){
857 AliPHOSDigit * digit;
858 printf("\nEMC digits (with primaries):\n") ;
859 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
860 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
862 for (index = 0 ; (index < digits->GetEntriesFast()) &&
863 (dynamic_cast<AliPHOSDigit *>(digits->At(index))->GetId() <= maxEmc) ; index++) {
864 digit = (AliPHOSDigit * ) digits->At(index) ;
865 if(digit->GetNprimary() == 0)
867 // printf("%6d %8d %6.5e %4d %2d :",
868 // digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ; // YVK
869 printf("%6d %.4f %6.5e %4d %2d :",
870 digit->GetId(), digit->GetEnergy(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
872 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
873 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
879 if(strstr(option,"all")||strstr(option,"CPV")){
881 //loop over CPV digits
882 AliPHOSDigit * digit;
883 printf("\nCPV digits (with primaries):\n") ;
884 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
885 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
887 for (index = 0 ; index < digits->GetEntriesFast(); index++) {
888 digit = (AliPHOSDigit * ) digits->At(index) ;
889 if(digit->GetId() > maxEmc){
890 printf("%6d %8d %4d %2d :",
891 digit->GetId(), digit->GetAmp(), digit->GetIndexInList(), digit->GetNprimary()) ;
893 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
894 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
903 //__________________________________________________________________
904 Float_t AliPHOSDigitizer::TimeOfNoise(void) const
905 { // Calculates the time signal generated by noise
906 //PH Info("TimeOfNoise", "Change me") ;
907 return gRandom->Rndm() * 1.28E-5;
910 //__________________________________________________________________
911 void AliPHOSDigitizer::Unload()
915 for(i = 1 ; i < fInput ; i++){
916 TString tempo(fEventNames[i]) ;
918 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
919 gime->PhosLoader()->UnloadSDigits() ;
922 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
923 gime->PhosLoader()->UnloadDigits() ;
926 //____________________________________________________________________________
927 void AliPHOSDigitizer::WriteDigits()
930 // Makes TreeD in the output file.
931 // Check if branch already exists:
932 // if yes, exit without writing: ROOT TTree does not support overwriting/updating of
933 // already existing branches.
934 // else creates branch with Digits, named "PHOS", title "...",
935 // and branch "AliPHOSDigitizer", with the same title to keep all the parameters
936 // and names of files, from which digits are made.
938 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
939 const TClonesArray * digits = gime->Digits() ;
940 TTree * treeD = gime->TreeD();
942 // -- create Digits branch
943 Int_t bufferSize = 32000 ;
944 TBranch * digitsBranch = treeD->Branch("PHOS","TClonesArray",&digits,bufferSize);
945 digitsBranch->SetTitle(fEventFolderName);
946 digitsBranch->Fill() ;
948 gime->WriteDigits("OVERWRITE");
949 gime->WriteDigitizer("OVERWRITE");