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.100 2007/10/10 09:05:10 schutz
22 * Changing name QualAss to QA
24 * Revision 1.99 2007/09/30 17:08:20 schutz
25 * Introducing the notion of QA data acquisition cycle (needed by online)
27 * Revision 1.98 2007/09/26 14:22:17 cvetan
28 * Important changes to the reconstructor classes. Complete elimination of the run-loaders, which are now steered only from AliReconstruction. Removal of the corresponding Reconstruct() and FillESD() methods.
30 * Revision 1.97 2007/08/07 14:12:03 kharlov
31 * Quality assurance added (Yves Schutz)
33 * Revision 1.96 2007/04/28 10:43:36 policheh
34 * Dead channels simulation: digit energy sets to 0.
36 * Revision 1.95 2007/04/10 07:20:52 kharlov
37 * Decalibration should use the same CDB as calibration in AliPHOSClusterizerv1
39 * Revision 1.94 2007/02/01 10:34:47 hristov
40 * Removing warnings on Solaris x86
42 * Revision 1.93 2006/10/17 13:17:01 kharlov
43 * Replace AliInfo by AliDebug
45 * Revision 1.92 2006/08/28 10:01:56 kharlov
46 * Effective C++ warnings fixed (Timur Pocheptsov)
48 * Revision 1.91 2006/04/29 20:25:30 hristov
49 * Decalibration is implemented (Yu.Kharlov)
51 * Revision 1.90 2006/04/22 10:30:17 hristov
52 * Add fEnergy to AliPHOSDigit and operate with EMC amplitude in energy units (Yu.Kharlov)
54 * Revision 1.89 2006/04/11 15:22:59 hristov
55 * run number in query set to -1: forces AliCDBManager to use its run number (A.Colla)
57 * Revision 1.88 2006/03/13 14:05:43 kharlov
58 * Calibration objects for EMC and CPV
60 * Revision 1.87 2005/08/24 15:33:49 kharlov
61 * Calibration data for raw digits
63 * Revision 1.86 2005/07/12 20:07:35 hristov
64 * Changes needed to run simulation and reconstrruction in the same AliRoot session
66 * Revision 1.85 2005/05/28 14:19:04 schutz
67 * Compilation warnings fixed by T.P.
71 //_________________________________________________________________________
72 //*-- Author : Dmitri Peressounko (SUBATECH & Kurchatov Institute)
73 //////////////////////////////////////////////////////////////////////////////
74 // This TTask performs digitization of Summable digits (in the PHOS case it is just
75 // the sum of contributions from all primary particles into a given cell).
76 // In addition it performs mixing of summable digits from different events.
77 // The name of the TTask is also the title of the branch that will contain
78 // the created SDigits
79 // The title of the TTAsk is the name of the file that contains the hits from
80 // which the SDigits are created
82 // For each event two branches are created in TreeD:
83 // "PHOS" - list of digits
84 // "AliPHOSDigitizer" - AliPHOSDigitizer with all parameters used in digitization
86 // Note, that one can set a title for new digits branch, and repeat digitization with
87 // another set of parameters.
90 // root[0] AliPHOSDigitizer * d = new AliPHOSDigitizer() ;
91 // root[1] d->ExecuteTask()
92 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
93 // //Digitizes SDigitis in all events found in file galice.root
95 // root[2] AliPHOSDigitizer * d1 = new AliPHOSDigitizer("galice1.root") ;
96 // // Will read sdigits from galice1.root
97 // root[3] d1->MixWith("galice2.root")
98 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
99 // // Reads another set of sdigits from galice2.root
100 // root[3] d1->MixWith("galice3.root")
101 // // Reads another set of sdigits from galice3.root
102 // root[4] d->ExecuteTask("deb timing")
103 // // Reads SDigits from files galice1.root, galice2.root ....
104 // // mixes them and stores produced Digits in file galice1.root
105 // // deb - prints number of produced digits
106 // // deb all - prints list of produced digits
107 // // timing - prints time used for digitization
110 // --- ROOT system ---
113 #include "TBenchmark.h"
116 // --- Standard library ---
118 // --- AliRoot header files ---
120 #include "AliRunDigitizer.h"
121 #include "AliPHOSDigit.h"
122 #include "AliPHOSGetter.h"
123 #include "AliPHOSDigitizer.h"
124 #include "AliPHOSSDigitizer.h"
125 #include "AliPHOSGeometry.h"
126 #include "AliPHOSTick.h"
127 #include "AliPHOSQADataMaker.h"
129 ClassImp(AliPHOSDigitizer)
132 //____________________________________________________________________________
133 AliPHOSDigitizer::AliPHOSDigitizer() :
139 fInputFileNames(0x0),
143 fEMCDigitThreshold(0.f),
145 fCPVDigitThreshold(0.f),
146 fTimeResolution(0.f),
148 fTimeSignalLength(0.f),
150 fADCpedestalEmc(0.f),
153 fADCpedestalCpv(0.f),
155 fEventFolderName(""),
163 fManager = 0 ; // We work in the standalong mode
166 //____________________________________________________________________________
167 AliPHOSDigitizer::AliPHOSDigitizer(TString alirunFileName,
168 TString eventFolderName):
169 AliDigitizer("PHOS"+AliConfig::Instance()->GetDigitizerTaskName(), alirunFileName),
170 fDefaultInit(kFALSE),
174 fInputFileNames(0x0),
178 fEMCDigitThreshold(0.f),
180 fCPVDigitThreshold(0.f),
181 fTimeResolution(0.f),
183 fTimeSignalLength(0.f),
185 fADCpedestalEmc(0.f),
188 fADCpedestalCpv(0.f),
190 fEventFolderName(eventFolderName),
199 fDefaultInit = kFALSE ;
200 fManager = 0 ; // We work in the standalong mode
201 //Initialize the quality assurance data maker only once
202 fQADM = new AliPHOSQADataMaker() ;
203 //FIXME: get the run number
206 GetQADataMaker()->Init(AliQA::kDIGITS, run, fgkCycles) ;
207 GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
210 //____________________________________________________________________________
211 AliPHOSDigitizer::AliPHOSDigitizer(const AliPHOSDigitizer & d) :
213 fDefaultInit(d.fDefaultInit),
214 fDigitsInRun(d.fDigitsInRun),
217 fInputFileNames(0x0),//?
219 fEmcCrystals(d.fEmcCrystals),
220 fPinNoise(d.fPinNoise),
221 fEMCDigitThreshold(d.fEMCDigitThreshold),
222 fCPVNoise(d.fCPVNoise),
223 fCPVDigitThreshold(d.fCPVDigitThreshold),
224 fTimeResolution(d.fTimeResolution),
225 fTimeThreshold(d.fTimeThreshold),
226 fTimeSignalLength(d.fTimeSignalLength),
227 fADCchanelEmc(d.fADCchanelEmc),
228 fADCpedestalEmc(d.fADCpedestalEmc),
229 fNADCemc(d.fNADCemc),
230 fADCchanelCpv(d.fADCchanelCpv),
231 fADCpedestalCpv(d.fADCpedestalCpv),
232 fNADCcpv(d.fNADCcpv),
233 fEventFolderName(d.fEventFolderName),
234 fFirstEvent(d.fFirstEvent),
235 fLastEvent(d.fLastEvent),
241 SetName(d.GetName()) ;
242 SetTitle(d.GetTitle()) ;
243 //Initialize the quality assurance data maker only once
244 //FIXME: get the run number
247 GetQADataMaker()->Init(AliQA::kDIGITS, run, fgkCycles) ;
248 GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
251 //____________________________________________________________________________
252 AliPHOSDigitizer::AliPHOSDigitizer(AliRunDigitizer * rd) :
253 AliDigitizer(rd,"PHOS"+AliConfig::Instance()->GetDigitizerTaskName()),
254 fDefaultInit(kFALSE),
258 fInputFileNames(0x0),
262 fEMCDigitThreshold(0.f),
264 fCPVDigitThreshold(0.f),
265 fTimeResolution(0.f),
267 fTimeSignalLength(0.f),
269 fADCpedestalEmc(0.f),
272 fADCpedestalCpv(0.f),
274 fEventFolderName(fManager->GetInputFolderName(0)),
281 // ctor Init() is called by RunDigitizer
283 SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
285 fDefaultInit = kFALSE ;
286 //Initialize the quality assurance data maker only once
287 fQADM = new AliPHOSQADataMaker() ;
288 //FIXME: get the run number
291 GetQADataMaker()->Init(AliQA::kDIGITS, run) ;
292 GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
295 //____________________________________________________________________________
296 AliPHOSDigitizer::~AliPHOSDigitizer()
298 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
300 // Clean Digitizer from the white board
301 gime->PhosLoader()->CleanDigitizer() ;
303 delete [] fInputFileNames ;
304 delete [] fEventNames ;
310 //____________________________________________________________________________
311 void AliPHOSDigitizer::Digitize(Int_t event)
314 // Makes the digitization of the collected summable digits.
315 // It first creates the array of all PHOS modules
316 // filled with noise (different for EMC, and CPV) and
317 // then adds contributions from SDigits.
318 // This design avoids scanning over the list of digits to add
319 // contribution to new SDigits only.
321 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
322 Int_t ReadEvent = event ;
324 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetCurrentEventNumber() ;
325 AliDebug(1,Form("Adding event %d from input stream 0 %s %s",
326 ReadEvent, GetTitle(), fEventFolderName.Data())) ;
327 gime->Event(ReadEvent, "S") ;
328 TClonesArray * digits = gime->Digits() ;
331 const AliPHOSGeometry *geom = gime->PHOSGeometry() ;
332 //Making digits with noise, first EMC
333 Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
338 nCPV = nEMC + geom->GetNumberOfCPVPadsZ() * geom->GetNumberOfCPVPadsPhi() * geom->GetNModules() ;
340 digits->Expand(nCPV) ;
342 // get first the sdigitizer from the tasks list
343 if ( !gime->SDigitizer() )
344 gime->LoadSDigitizer();
345 AliPHOSSDigitizer * sDigitizer = gime->SDigitizer();
348 AliFatal(Form("SDigitizer with name %s %s not found",
349 GetTitle(), fEventFolderName.Data() )) ;
351 //take all the inputs to add together and load the SDigits
352 TObjArray * sdigArray = new TObjArray(fInput) ;
353 sdigArray->AddAt(gime->SDigits(), 0) ;
355 for(i = 1 ; i < fInput ; i++){
356 TString tempo(fEventNames[i]) ;
358 AliPHOSGetter * gime1 = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
360 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(i))->GetCurrentEventNumber() ;
361 AliInfo(Form("Adding event %d from input stream %d %s %s",
362 ReadEvent, i, fInputFileNames[i].Data(), tempo.Data())) ;
363 gime1->Event(ReadEvent,"S");
364 sdigArray->AddAt(gime1->SDigits(), i) ;
367 //Find the first crystal with signal
368 Int_t nextSig = 200000 ;
369 TClonesArray * sdigits ;
370 for(i = 0 ; i < fInput ; i++){
371 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
372 if ( !sdigits->GetEntriesFast() )
374 Int_t curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(0))->GetId() ;
375 if(curNext < nextSig)
379 TArrayI index(fInput) ;
380 index.Reset() ; //Set all indexes to zero
382 AliPHOSDigit * digit ;
383 AliPHOSDigit * curSDigit ;
385 TClonesArray * ticks = new TClonesArray("AliPHOSTick",1000) ;
387 //Put Noise contribution
388 for(absID = 1 ; absID <= nEMC ; absID++){
389 Float_t noise = gRandom->Gaus(0., fPinNoise) ;
390 // YVK: do not digitize amplitudes for EMC
391 // new((*digits)[absID-1]) AliPHOSDigit( -1, absID, sDigitizer->Digitize(noise), TimeOfNoise() ) ;
392 new((*digits)[absID-1]) AliPHOSDigit( -1, absID, noise, TimeOfNoise() ) ;
393 //look if we have to add signal?
394 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
397 //Add SDigits from all inputs
400 Float_t a = digit->GetEnergy() ;
401 Float_t b = TMath::Abs( a / fTimeSignalLength) ;
402 //Mark the beginning of the signal
403 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime(),0, b);
404 //Mark the end of the signal
405 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime()+fTimeSignalLength, -a, -b);
408 for(i = 0 ; i < fInput ; i++){
409 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
410 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
413 //May be several digits will contribute from the same input
414 while(curSDigit && curSDigit->GetId() == absID){
415 //Shift primary to separate primaries belonging different inputs
416 Int_t primaryoffset ;
418 primaryoffset = fManager->GetMask(i) ;
420 primaryoffset = 10000000*i ;
421 curSDigit->ShiftPrimary(primaryoffset) ;
423 a = curSDigit->GetEnergy() ;
424 b = a /fTimeSignalLength ;
425 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime(),0, b);
426 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b);
428 *digit += *curSDigit ; //add energies
431 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
432 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
438 //calculate and set time
439 Float_t time = FrontEdgeTime(ticks) ;
440 digit->SetTime(time) ;
442 //Find next signal module
444 for(i = 0 ; i < fInput ; i++){
445 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
446 Int_t curNext = nextSig ;
447 if(sdigits->GetEntriesFast() > index[i] ){
448 curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(index[i]))->GetId() ;
450 if(curNext < nextSig) nextSig = curNext ;
458 //Now CPV digits (different noise and no timing)
459 for(absID = nEMC+1; absID <= nCPV; absID++){
460 Float_t noise = gRandom->Gaus(0., fCPVNoise) ;
461 new((*digits)[absID-1]) AliPHOSDigit( -1,absID,sDigitizer->Digitize(noise), TimeOfNoise() ) ;
462 //look if we have to add signal?
464 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
465 //Add SDigits from all inputs
466 for(i = 0 ; i < fInput ; i++){
467 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
468 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
472 //May be several digits will contribute from the same input
473 while(curSDigit && curSDigit->GetId() == absID){
474 //Shift primary to separate primaries belonging different inputs
475 Int_t primaryoffset ;
477 primaryoffset = fManager->GetMask(i) ;
479 primaryoffset = 10000000*i ;
480 curSDigit->ShiftPrimary(primaryoffset) ;
483 *digit += *curSDigit ;
485 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
486 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i]) ) ;
492 //Find next signal module
494 for(i = 0 ; i < fInput ; i++){
495 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
496 Int_t curNext = nextSig ;
497 if(sdigits->GetEntriesFast() > index[i] )
498 curNext = dynamic_cast<AliPHOSDigit *>( sdigits->At(index[i]) )->GetId() ;
499 if(curNext < nextSig) nextSig = curNext ;
505 delete sdigArray ; //We should not delete its contents
507 //remove digits below thresholds
508 for(i = 0 ; i < nEMC ; i++){
509 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
510 DecalibrateEMC(digit);
511 if(digit->GetEnergy() < fEMCDigitThreshold)
512 digits->RemoveAt(i) ;
514 digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ;
518 for(i = nEMC; i < nCPV ; i++)
519 // if( sDigitizer->Calibrate( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetAmp() ) < fCPVDigitThreshold )
520 if( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetEnergy() < fCPVDigitThreshold )
521 digits->RemoveAt(i) ;
525 Int_t ndigits = digits->GetEntriesFast() ;
526 digits->Expand(ndigits) ;
528 //Set indexes in list of digits and make true digitization of the energy
529 for (i = 0 ; i < ndigits ; i++) {
530 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
531 digit->SetIndexInList(i) ;
532 if(digit->GetId() > fEmcCrystals){ //digitize CPV only
533 digit->SetAmp(DigitizeCPV(digit->GetEnergy(),digit->GetId()) ) ;
539 //set amplitudes in bad channels to zero
540 for(i = 0 ; i <digits->GetEntries(); i++){
541 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
542 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
543 if(relId[1] == 0) // Emc
544 if(gime->CalibData()->IsBadChannelEmc(relId[0],relId[3],relId[2])) digit->SetEnergy(0.);
549 //____________________________________________________________________________
550 void AliPHOSDigitizer::DecalibrateEMC(AliPHOSDigit *digit)
552 // Decalibrate EMC digit, i.e. change its energy by a factor read from CDB
554 AliPHOSGetter* gime = AliPHOSGetter::Instance();
556 if(!gime->CalibData()) {
557 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1);
558 gime->SetCalibData(cdb);
561 //Determine rel.position of the cell absolute ID
563 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
564 Int_t module=relId[0];
566 Int_t column=relId[3];
567 Float_t decalibration = gime->CalibData()->GetADCchannelEmc(module,column,row);
568 Float_t energy = digit->GetEnergy() / decalibration;
569 digit->SetEnergy(energy);
571 //____________________________________________________________________________
572 Int_t AliPHOSDigitizer::DigitizeCPV(Float_t charge, Int_t absId)
574 // Returns digitized value of the CPV charge in a pad absId
576 AliPHOSGetter* gime = AliPHOSGetter::Instance();
578 if(!gime->CalibData()) {
579 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1); // use AliCDBManager's run number
580 gime->SetCalibData(cdb);
583 //Determine rel.position of the cell absId
585 gime->PHOSGeometry()->AbsToRelNumbering(absId,relId);
586 Int_t module=relId[0];
588 Int_t column=relId[3];
592 if(absId > fEmcCrystals){ //digitize CPV only
594 //reading calibration data for cell absId.
595 //If no calibration DB found, accept default values.
597 if(gime->CalibData()) {
598 fADCpedestalCpv = gime->CalibData()->GetADCpedestalCpv(module,column,row);
599 fADCchanelCpv = gime->CalibData()->GetADCchannelCpv( module,column,row);
602 channel = (Int_t) TMath::Ceil((charge - fADCpedestalCpv)/fADCchanelCpv) ;
603 if(channel > fNADCcpv ) channel = fNADCcpv ;
608 //____________________________________________________________________________
609 void AliPHOSDigitizer::Exec(Option_t *option)
611 // Steering method to process digitization for events
612 // in the range from fFirstEvent to fLastEvent.
613 // This range is optionally set by SetEventRange().
614 // if fLastEvent=-1, then process events until the end.
615 // by default fLastEvent = fFirstEvent (process only one event)
617 if (!fInit) { // to prevent overwrite existing file
618 AliError(Form("Give a version name different from %s",
619 fEventFolderName.Data() )) ;
623 if (strstr(option,"print")) {
628 // check the QA result for Hits and SDigits
629 AliQA * qa = AliQA::Instance(AliQA::kPHOS) ;
630 if ( qa->IsSet(AliQA::kPHOS, AliQA::kSIM, AliQA::kFATAL)) {
631 AliFatal("QA status in Hits and/or SDIGITS was Fatal") ;
632 } else if ( qa->IsSet(AliQA::kPHOS, AliQA::kSIM, AliQA::kERROR)) {
633 AliError("QA status in Hits and/or SDIGITS was Error") ;
634 } else if ( qa->IsSet(AliQA::kPHOS, AliQA::kSIM, AliQA::kWARNING) ) {
635 AliWarning("QA status in Hits and/or SDIGITS was Warning") ;
636 } else if ( qa->IsSet(AliQA::kPHOS, AliQA::kSIM, AliQA::kINFO) ) {
637 AliInfo("QA status in Hits and/or SDIGITS was Info") ;
640 if(strstr(option,"tim"))
641 gBenchmark->Start("PHOSDigitizer");
643 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
645 // Post Digitizer to the white board
646 gime->PostDigitizer(this) ;
648 if (fLastEvent == -1)
649 fLastEvent = gime->MaxEvent() - 1 ;
651 fLastEvent = fFirstEvent ;
653 Int_t nEvents = fLastEvent - fFirstEvent + 1;
657 for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
659 gime->Event(ievent,"S") ;
661 Digitize(ievent) ; //Add prepared SDigits to digits and add the noise
663 //makes the quality assurance data
664 if (GetQADataMaker()->IsCycleDone() ) {
665 GetQADataMaker()->EndOfCycle(AliQA::kDIGITS) ;
666 GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
668 GetQADataMaker()->Exec(AliQA::kDIGITS, gime->Digits()) ;
669 GetQADataMaker()->Increment() ;
673 if(strstr(option,"deb"))
676 //increment the total number of Digits per run
677 fDigitsInRun += gime->Digits()->GetEntriesFast() ;
680 //Write the quality assurance data only after the last event
681 if ( fEventCounter == gime->MaxEvent() ) {
682 GetQADataMaker()->EndOfCycle(AliQA::kDIGITS) ;
683 GetQADataMaker()->Finish(AliQA::kDIGITS) ;
686 gime->PhosLoader()->CleanDigitizer();
688 if(strstr(option,"tim")){
689 gBenchmark->Stop("PHOSDigitizer");
691 message = " took %f seconds for Digitizing %f seconds per event\n" ;
692 AliInfo(Form( message.Data(),
693 gBenchmark->GetCpuTime("PHOSDigitizer"),
694 gBenchmark->GetCpuTime("PHOSDigitizer")/nEvents ));
698 //____________________________________________________________________________
699 Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks) const
701 // Returns the shortest time among all time ticks
703 ticks->Sort() ; //Sort in accordance with times of ticks
705 AliPHOSTick * ctick = (AliPHOSTick *) it.Next() ;
706 Float_t time = ctick->CrossingTime(fTimeThreshold) ;
709 while((t=(AliPHOSTick*) it.Next())){
710 if(t->GetTime() < time) //This tick starts before crossing
715 time = ctick->CrossingTime(fTimeThreshold) ;
720 //____________________________________________________________________________
721 Bool_t AliPHOSDigitizer::Init()
723 // Makes all memory allocations
725 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
727 AliFatal(Form("Could not obtain the Getter object for file %s and event %s !",
728 GetTitle(), fEventFolderName.Data()));
732 const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
734 fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ;
736 TString opt("Digits") ;
737 if(gime->VersionExists(opt) ) {
738 AliError(Form("Give a version name different from %s",
739 fEventFolderName.Data() )) ;
744 fLastEvent = fFirstEvent ;
746 fInput = fManager->GetNinputs() ;
750 fInputFileNames = new TString[fInput] ;
751 fEventNames = new TString[fInput] ;
752 fInputFileNames[0] = GetTitle() ;
753 fEventNames[0] = fEventFolderName.Data() ;
755 for (index = 1 ; index < fInput ; index++) {
756 fInputFileNames[index] = dynamic_cast<AliStream*>(fManager->GetInputStream(index))->GetFileName(0);
757 TString tempo = fManager->GetInputFolderName(index) ;
758 fEventNames[index] = tempo.Remove(tempo.Length()-1) ; // strip of the stream number added by fManager
761 //to prevent cleaning of this object while GetEvent is called
762 gime->PhosLoader()->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
767 //____________________________________________________________________________
768 void AliPHOSDigitizer::InitParameters()
770 // Set initial parameters Digitizer
772 fPinNoise = 0.004 ; // [GeV]
773 fEMCDigitThreshold = 0.012 ; // [GeV]
774 fCPVNoise = 0.01; // [aux units]
775 fCPVDigitThreshold = 0.09 ; // [aux units]
776 fTimeResolution = 0.5e-9 ; // [sec]
777 fTimeSignalLength = 1.0e-9 ; // [sec]
779 fADCchanelEmc = 1.0; // Coefficient between real and measured energies in EMC
780 fADCpedestalEmc = 0. ; //
781 fNADCemc = (Int_t) TMath::Power(2,16) ; // number of channels in EMC ADC
783 fADCchanelCpv = 0.0012 ; // width of one ADC channel in CPV 'popugais'
784 fADCpedestalCpv = 0.012 ; //
785 fNADCcpv = (Int_t) TMath::Power(2,12); // number of channels in CPV ADC
787 // fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
788 fTimeThreshold = 0.001 ; // [GeV]
789 SetEventRange(0,-1) ;
793 //__________________________________________________________________
794 void AliPHOSDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
796 // Allows to produce digits by superimposing background and signal event.
797 // It is assumed, that headers file with SIGNAL events is opened in
799 // Sets the BACKGROUND event, with which the SIGNAL event is to be mixed
800 // Thus we avoid writing (changing) huge and expensive
801 // backgound files: all output will be writen into SIGNAL, i.e.
802 // opened in constructor file.
804 // One can open as many files to mix with as one needs.
805 // However only Sdigits with the same name (i.e. constructed with the same SDigitizer)
808 if( strcmp(fEventFolderName, "") == 0 )
812 Warning("MixWith", "Cannot use this method with AliRunDigitizer\n" ) ;
815 // looking for file which contains AliRun
816 if (gSystem->AccessPathName(alirunFileName)) {// file does not exist
817 AliError(Form("File %s does not exist!", alirunFileName.Data())) ;
820 // looking for the file which contains SDigits
821 AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
822 TString fileName( gime->GetSDigitsFileName() ) ;
823 if ( eventFolderName != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
824 fileName = fileName.ReplaceAll(".root", "") + "_" + eventFolderName + ".root" ;
825 if ( (gSystem->AccessPathName(fileName)) ) {
826 AliError(Form("The file %s does not exist!", fileName.Data())) ;
829 // need to increase the arrays
830 TString tempo = fInputFileNames[fInput-1] ;
831 delete [] fInputFileNames ;
832 fInputFileNames = new TString[fInput+1] ;
833 fInputFileNames[fInput-1] = tempo ;
835 tempo = fEventNames[fInput-1] ;
836 delete [] fEventNames ;
837 fEventNames = new TString[fInput+1] ;
838 fEventNames[fInput-1] = tempo ;
840 fInputFileNames[fInput] = alirunFileName ;
841 fEventNames[fInput] = eventFolderName ;
845 //__________________________________________________________________
846 void AliPHOSDigitizer::Print(const Option_t *)const
848 // Print Digitizer's parameters
849 AliInfo(Form("\n------------------- %s -------------", GetName() )) ;
850 if( strcmp(fEventFolderName.Data(), "") != 0 ){
851 printf(" Writing Digits to branch with title %s\n", fEventFolderName.Data()) ;
855 nStreams = GetNInputStreams() ;
860 for (index = 0 ; index < nStreams ; index++) {
861 TString tempo(fEventNames[index]) ;
863 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[index], tempo) ;
864 TString fileName( gime->GetSDigitsFileName() ) ;
865 if ( fEventNames[index] != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
866 fileName = fileName.ReplaceAll(".root", "") + "_" + fEventNames[index] + ".root" ;
867 printf ("Adding SDigits from %s %s\n", fInputFileNames[index].Data(), fileName.Data()) ;
869 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
870 printf("\nWriting digits to %s", gime->GetDigitsFileName().Data()) ;
872 printf("\nWith following parameters:\n") ;
873 printf(" Electronics noise in EMC (fPinNoise) = %f GeV\n", fPinNoise ) ;
874 printf(" Threshold in EMC (fEMCDigitThreshold) = %f GeV\n", fEMCDigitThreshold ) ;
875 printf(" Noise in CPV (fCPVNoise) = %f aux units\n", fCPVNoise ) ;
876 printf(" Threshold in CPV (fCPVDigitThreshold) = %f aux units\n",fCPVDigitThreshold ) ;
877 printf(" ---------------------------------------------------\n") ;
880 AliInfo(Form("AliPHOSDigitizer not initialized" )) ;
884 //__________________________________________________________________
885 void AliPHOSDigitizer::PrintDigits(Option_t * option)
887 // Print a table of digits
889 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
890 TClonesArray * digits = gime->Digits() ;
892 AliInfo(Form("%d", digits->GetEntriesFast())) ;
893 printf("\nevent %d", gAlice->GetEvNumber()) ;
894 printf("\n Number of entries in Digits list %d", digits->GetEntriesFast() ) ;
897 if(strstr(option,"all")||strstr(option,"EMC")){
899 AliPHOSDigit * digit;
900 printf("\nEMC digits (with primaries):\n") ;
901 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
902 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
904 for (index = 0 ; (index < digits->GetEntriesFast()) &&
905 (dynamic_cast<AliPHOSDigit *>(digits->At(index))->GetId() <= maxEmc) ; index++) {
906 digit = (AliPHOSDigit * ) digits->At(index) ;
907 if(digit->GetNprimary() == 0)
909 // printf("%6d %8d %6.5e %4d %2d :",
910 // digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ; // YVK
911 printf("%6d %.4f %6.5e %4d %2d :",
912 digit->GetId(), digit->GetEnergy(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
914 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
915 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
921 if(strstr(option,"all")||strstr(option,"CPV")){
923 //loop over CPV digits
924 AliPHOSDigit * digit;
925 printf("\nCPV digits (with primaries):\n") ;
926 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
927 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
929 for (index = 0 ; index < digits->GetEntriesFast(); index++) {
930 digit = (AliPHOSDigit * ) digits->At(index) ;
931 if(digit->GetId() > maxEmc){
932 printf("%6d %8d %4d %2d :",
933 digit->GetId(), digit->GetAmp(), digit->GetIndexInList(), digit->GetNprimary()) ;
935 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
936 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
945 //__________________________________________________________________
946 Float_t AliPHOSDigitizer::TimeOfNoise(void) const
947 { // Calculates the time signal generated by noise
948 //PH Info("TimeOfNoise", "Change me") ;
949 return gRandom->Rndm() * 1.28E-5;
952 //__________________________________________________________________
953 void AliPHOSDigitizer::Unload()
957 for(i = 1 ; i < fInput ; i++){
958 TString tempo(fEventNames[i]) ;
960 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
961 gime->PhosLoader()->UnloadSDigits() ;
964 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
965 gime->PhosLoader()->UnloadDigits() ;
968 //____________________________________________________________________________
969 void AliPHOSDigitizer::WriteDigits()
972 // Makes TreeD in the output file.
973 // Check if branch already exists:
974 // if yes, exit without writing: ROOT TTree does not support overwriting/updating of
975 // already existing branches.
976 // else creates branch with Digits, named "PHOS", title "...",
977 // and branch "AliPHOSDigitizer", with the same title to keep all the parameters
978 // and names of files, from which digits are made.
980 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
981 const TClonesArray * digits = gime->Digits() ;
982 TTree * treeD = gime->TreeD();
984 // -- create Digits branch
985 Int_t bufferSize = 32000 ;
986 TBranch * digitsBranch = treeD->Branch("PHOS","TClonesArray",&digits,bufferSize);
987 digitsBranch->SetTitle(fEventFolderName);
988 digitsBranch->Fill() ;
990 gime->WriteDigits("OVERWRITE");
991 gime->WriteDigitizer("OVERWRITE");