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.101 2007/10/14 21:08:10 schutz
22 * Introduced the checking of QA results from previous step before entering the event loop
24 * Revision 1.100 2007/10/10 09:05:10 schutz
25 * Changing name QualAss to QA
27 * Revision 1.99 2007/09/30 17:08:20 schutz
28 * Introducing the notion of QA data acquisition cycle (needed by online)
30 * Revision 1.98 2007/09/26 14:22:17 cvetan
31 * 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.
33 * Revision 1.97 2007/08/07 14:12:03 kharlov
34 * Quality assurance added (Yves Schutz)
36 * Revision 1.96 2007/04/28 10:43:36 policheh
37 * Dead channels simulation: digit energy sets to 0.
39 * Revision 1.95 2007/04/10 07:20:52 kharlov
40 * Decalibration should use the same CDB as calibration in AliPHOSClusterizerv1
42 * Revision 1.94 2007/02/01 10:34:47 hristov
43 * Removing warnings on Solaris x86
45 * Revision 1.93 2006/10/17 13:17:01 kharlov
46 * Replace AliInfo by AliDebug
48 * Revision 1.92 2006/08/28 10:01:56 kharlov
49 * Effective C++ warnings fixed (Timur Pocheptsov)
51 * Revision 1.91 2006/04/29 20:25:30 hristov
52 * Decalibration is implemented (Yu.Kharlov)
54 * Revision 1.90 2006/04/22 10:30:17 hristov
55 * Add fEnergy to AliPHOSDigit and operate with EMC amplitude in energy units (Yu.Kharlov)
57 * Revision 1.89 2006/04/11 15:22:59 hristov
58 * run number in query set to -1: forces AliCDBManager to use its run number (A.Colla)
60 * Revision 1.88 2006/03/13 14:05:43 kharlov
61 * Calibration objects for EMC and CPV
63 * Revision 1.87 2005/08/24 15:33:49 kharlov
64 * Calibration data for raw digits
66 * Revision 1.86 2005/07/12 20:07:35 hristov
67 * Changes needed to run simulation and reconstrruction in the same AliRoot session
69 * Revision 1.85 2005/05/28 14:19:04 schutz
70 * Compilation warnings fixed by T.P.
74 //_________________________________________________________________________
75 //*-- Author : Dmitri Peressounko (SUBATECH & Kurchatov Institute)
76 //////////////////////////////////////////////////////////////////////////////
77 // This TTask performs digitization of Summable digits (in the PHOS case it is just
78 // the sum of contributions from all primary particles into a given cell).
79 // In addition it performs mixing of summable digits from different events.
80 // The name of the TTask is also the title of the branch that will contain
81 // the created SDigits
82 // The title of the TTAsk is the name of the file that contains the hits from
83 // which the SDigits are created
85 // For each event two branches are created in TreeD:
86 // "PHOS" - list of digits
87 // "AliPHOSDigitizer" - AliPHOSDigitizer with all parameters used in digitization
89 // Note, that one can set a title for new digits branch, and repeat digitization with
90 // another set of parameters.
93 // root[0] AliPHOSDigitizer * d = new AliPHOSDigitizer() ;
94 // root[1] d->ExecuteTask()
95 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
96 // //Digitizes SDigitis in all events found in file galice.root
98 // root[2] AliPHOSDigitizer * d1 = new AliPHOSDigitizer("galice1.root") ;
99 // // Will read sdigits from galice1.root
100 // root[3] d1->MixWith("galice2.root")
101 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
102 // // Reads another set of sdigits from galice2.root
103 // root[3] d1->MixWith("galice3.root")
104 // // Reads another set of sdigits from galice3.root
105 // root[4] d->ExecuteTask("deb timing")
106 // // Reads SDigits from files galice1.root, galice2.root ....
107 // // mixes them and stores produced Digits in file galice1.root
108 // // deb - prints number of produced digits
109 // // deb all - prints list of produced digits
110 // // timing - prints time used for digitization
113 // --- ROOT system ---
116 #include "TBenchmark.h"
119 // --- Standard library ---
121 // --- AliRoot header files ---
123 #include "AliRunDigitizer.h"
124 #include "AliPHOSDigit.h"
125 #include "AliPHOSGetter.h"
126 #include "AliPHOSDigitizer.h"
127 #include "AliPHOSSDigitizer.h"
128 #include "AliPHOSGeometry.h"
129 #include "AliPHOSTick.h"
130 #include "AliPHOSQADataMaker.h"
132 ClassImp(AliPHOSDigitizer)
135 //____________________________________________________________________________
136 AliPHOSDigitizer::AliPHOSDigitizer() :
142 fInputFileNames(0x0),
146 fEMCDigitThreshold(0.f),
148 fCPVDigitThreshold(0.f),
149 fTimeResolution(0.f),
151 fTimeSignalLength(0.f),
153 fADCpedestalEmc(0.f),
156 fADCpedestalCpv(0.f),
158 fEventFolderName(""),
166 fManager = 0 ; // We work in the standalong mode
169 //____________________________________________________________________________
170 AliPHOSDigitizer::AliPHOSDigitizer(TString alirunFileName,
171 TString eventFolderName):
172 AliDigitizer("PHOS"+AliConfig::Instance()->GetDigitizerTaskName(), alirunFileName),
173 fDefaultInit(kFALSE),
177 fInputFileNames(0x0),
181 fEMCDigitThreshold(0.f),
183 fCPVDigitThreshold(0.f),
184 fTimeResolution(0.f),
186 fTimeSignalLength(0.f),
188 fADCpedestalEmc(0.f),
191 fADCpedestalCpv(0.f),
193 fEventFolderName(eventFolderName),
202 fDefaultInit = kFALSE ;
203 fManager = 0 ; // We work in the standalong mode
204 // //Initialize the quality assurance data maker only once
205 // fQADM = new AliPHOSQADataMaker() ;
206 // //FIXME: get the run number
209 // GetQADataMaker()->Init(AliQA::kDIGITS, run, fgkCycles) ;
210 // GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
213 //____________________________________________________________________________
214 AliPHOSDigitizer::AliPHOSDigitizer(const AliPHOSDigitizer & d) :
216 fDefaultInit(d.fDefaultInit),
217 fDigitsInRun(d.fDigitsInRun),
220 fInputFileNames(0x0),//?
222 fEmcCrystals(d.fEmcCrystals),
223 fPinNoise(d.fPinNoise),
224 fEMCDigitThreshold(d.fEMCDigitThreshold),
225 fCPVNoise(d.fCPVNoise),
226 fCPVDigitThreshold(d.fCPVDigitThreshold),
227 fTimeResolution(d.fTimeResolution),
228 fTimeThreshold(d.fTimeThreshold),
229 fTimeSignalLength(d.fTimeSignalLength),
230 fADCchanelEmc(d.fADCchanelEmc),
231 fADCpedestalEmc(d.fADCpedestalEmc),
232 fNADCemc(d.fNADCemc),
233 fADCchanelCpv(d.fADCchanelCpv),
234 fADCpedestalCpv(d.fADCpedestalCpv),
235 fNADCcpv(d.fNADCcpv),
236 fEventFolderName(d.fEventFolderName),
237 fFirstEvent(d.fFirstEvent),
238 fLastEvent(d.fLastEvent),
244 SetName(d.GetName()) ;
245 SetTitle(d.GetTitle()) ;
246 //Initialize the quality assurance data maker only once
247 //FIXME: get the run number
250 // GetQADataMaker()->Init(AliQA::kDIGITS, run, fgkCycles) ;
251 // GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
254 //____________________________________________________________________________
255 AliPHOSDigitizer::AliPHOSDigitizer(AliRunDigitizer * rd) :
256 AliDigitizer(rd,"PHOS"+AliConfig::Instance()->GetDigitizerTaskName()),
257 fDefaultInit(kFALSE),
261 fInputFileNames(0x0),
265 fEMCDigitThreshold(0.f),
267 fCPVDigitThreshold(0.f),
268 fTimeResolution(0.f),
270 fTimeSignalLength(0.f),
272 fADCpedestalEmc(0.f),
275 fADCpedestalCpv(0.f),
277 fEventFolderName(fManager->GetInputFolderName(0)),
284 // ctor Init() is called by RunDigitizer
286 SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
288 fDefaultInit = kFALSE ;
289 //Initialize the quality assurance data maker only once
290 // fQADM = new AliPHOSQADataMaker() ;
291 // //FIXME: get the run number
294 // GetQADataMaker()->Init(AliQA::kDIGITS, run) ;
295 // GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
298 //____________________________________________________________________________
299 AliPHOSDigitizer::~AliPHOSDigitizer()
301 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
303 // Clean Digitizer from the white board
304 gime->PhosLoader()->CleanDigitizer() ;
306 delete [] fInputFileNames ;
307 delete [] fEventNames ;
313 //____________________________________________________________________________
314 void AliPHOSDigitizer::Digitize(Int_t event)
317 // Makes the digitization of the collected summable digits.
318 // It first creates the array of all PHOS modules
319 // filled with noise (different for EMC, and CPV) and
320 // then adds contributions from SDigits.
321 // This design avoids scanning over the list of digits to add
322 // contribution to new SDigits only.
324 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
325 Int_t ReadEvent = event ;
327 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetCurrentEventNumber() ;
328 AliDebug(1,Form("Adding event %d from input stream 0 %s %s",
329 ReadEvent, GetTitle(), fEventFolderName.Data())) ;
330 gime->Event(ReadEvent, "S") ;
331 TClonesArray * digits = gime->Digits() ;
334 const AliPHOSGeometry *geom = gime->PHOSGeometry() ;
335 //Making digits with noise, first EMC
336 Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
341 nCPV = nEMC + geom->GetNumberOfCPVPadsZ() * geom->GetNumberOfCPVPadsPhi() * geom->GetNModules() ;
343 digits->Expand(nCPV) ;
345 // get first the sdigitizer from the tasks list
346 if ( !gime->SDigitizer() )
347 gime->LoadSDigitizer();
348 AliPHOSSDigitizer * sDigitizer = gime->SDigitizer();
351 AliFatal(Form("SDigitizer with name %s %s not found",
352 GetTitle(), fEventFolderName.Data() )) ;
354 //take all the inputs to add together and load the SDigits
355 TObjArray * sdigArray = new TObjArray(fInput) ;
356 sdigArray->AddAt(gime->SDigits(), 0) ;
358 for(i = 1 ; i < fInput ; i++){
359 TString tempo(fEventNames[i]) ;
361 AliPHOSGetter * gime1 = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
363 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(i))->GetCurrentEventNumber() ;
364 AliInfo(Form("Adding event %d from input stream %d %s %s",
365 ReadEvent, i, fInputFileNames[i].Data(), tempo.Data())) ;
366 gime1->Event(ReadEvent,"S");
367 sdigArray->AddAt(gime1->SDigits(), i) ;
370 //Find the first crystal with signal
371 Int_t nextSig = 200000 ;
372 TClonesArray * sdigits ;
373 for(i = 0 ; i < fInput ; i++){
374 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
375 if ( !sdigits->GetEntriesFast() )
377 Int_t curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(0))->GetId() ;
378 if(curNext < nextSig)
382 TArrayI index(fInput) ;
383 index.Reset() ; //Set all indexes to zero
385 AliPHOSDigit * digit ;
386 AliPHOSDigit * curSDigit ;
388 TClonesArray * ticks = new TClonesArray("AliPHOSTick",1000) ;
390 //Put Noise contribution
391 for(absID = 1 ; absID <= nEMC ; absID++){
392 Float_t noise = gRandom->Gaus(0., fPinNoise) ;
393 // YVK: do not digitize amplitudes for EMC
394 // new((*digits)[absID-1]) AliPHOSDigit( -1, absID, sDigitizer->Digitize(noise), TimeOfNoise() ) ;
395 new((*digits)[absID-1]) AliPHOSDigit( -1, absID, noise, TimeOfNoise() ) ;
396 //look if we have to add signal?
397 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
400 //Add SDigits from all inputs
403 Float_t a = digit->GetEnergy() ;
404 Float_t b = TMath::Abs( a / fTimeSignalLength) ;
405 //Mark the beginning of the signal
406 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime(),0, b);
407 //Mark the end of the signal
408 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime()+fTimeSignalLength, -a, -b);
411 for(i = 0 ; i < fInput ; i++){
412 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
413 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
416 //May be several digits will contribute from the same input
417 while(curSDigit && curSDigit->GetId() == absID){
418 //Shift primary to separate primaries belonging different inputs
419 Int_t primaryoffset ;
421 primaryoffset = fManager->GetMask(i) ;
423 primaryoffset = 10000000*i ;
424 curSDigit->ShiftPrimary(primaryoffset) ;
426 a = curSDigit->GetEnergy() ;
427 b = a /fTimeSignalLength ;
428 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime(),0, b);
429 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b);
431 *digit += *curSDigit ; //add energies
434 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
435 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
441 //calculate and set time
442 Float_t time = FrontEdgeTime(ticks) ;
443 digit->SetTime(time) ;
445 //Find next signal module
447 for(i = 0 ; i < fInput ; i++){
448 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
449 Int_t curNext = nextSig ;
450 if(sdigits->GetEntriesFast() > index[i] ){
451 curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(index[i]))->GetId() ;
453 if(curNext < nextSig) nextSig = curNext ;
461 //Now CPV digits (different noise and no timing)
462 for(absID = nEMC+1; absID <= nCPV; absID++){
463 Float_t noise = gRandom->Gaus(0., fCPVNoise) ;
464 new((*digits)[absID-1]) AliPHOSDigit( -1,absID,sDigitizer->Digitize(noise), TimeOfNoise() ) ;
465 //look if we have to add signal?
467 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
468 //Add SDigits from all inputs
469 for(i = 0 ; i < fInput ; i++){
470 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
471 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
475 //May be several digits will contribute from the same input
476 while(curSDigit && curSDigit->GetId() == absID){
477 //Shift primary to separate primaries belonging different inputs
478 Int_t primaryoffset ;
480 primaryoffset = fManager->GetMask(i) ;
482 primaryoffset = 10000000*i ;
483 curSDigit->ShiftPrimary(primaryoffset) ;
486 *digit += *curSDigit ;
488 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
489 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i]) ) ;
495 //Find next signal module
497 for(i = 0 ; i < fInput ; i++){
498 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
499 Int_t curNext = nextSig ;
500 if(sdigits->GetEntriesFast() > index[i] )
501 curNext = dynamic_cast<AliPHOSDigit *>( sdigits->At(index[i]) )->GetId() ;
502 if(curNext < nextSig) nextSig = curNext ;
508 delete sdigArray ; //We should not delete its contents
510 //remove digits below thresholds
511 for(i = 0 ; i < nEMC ; i++){
512 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
513 DecalibrateEMC(digit);
514 if(digit->GetEnergy() < fEMCDigitThreshold)
515 digits->RemoveAt(i) ;
517 digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ;
521 for(i = nEMC; i < nCPV ; i++)
522 // if( sDigitizer->Calibrate( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetAmp() ) < fCPVDigitThreshold )
523 if( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetEnergy() < fCPVDigitThreshold )
524 digits->RemoveAt(i) ;
528 Int_t ndigits = digits->GetEntriesFast() ;
529 digits->Expand(ndigits) ;
531 //Set indexes in list of digits and make true digitization of the energy
532 for (i = 0 ; i < ndigits ; i++) {
533 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
534 digit->SetIndexInList(i) ;
535 if(digit->GetId() > fEmcCrystals){ //digitize CPV only
536 digit->SetAmp(DigitizeCPV(digit->GetEnergy(),digit->GetId()) ) ;
542 //set amplitudes in bad channels to zero
543 for(i = 0 ; i <digits->GetEntries(); i++){
544 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
545 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
546 if(relId[1] == 0) // Emc
547 if(gime->CalibData()->IsBadChannelEmc(relId[0],relId[3],relId[2])) digit->SetEnergy(0.);
552 //____________________________________________________________________________
553 void AliPHOSDigitizer::DecalibrateEMC(AliPHOSDigit *digit)
555 // Decalibrate EMC digit, i.e. change its energy by a factor read from CDB
557 AliPHOSGetter* gime = AliPHOSGetter::Instance();
559 if(!gime->CalibData()) {
560 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1);
561 gime->SetCalibData(cdb);
564 //Determine rel.position of the cell absolute ID
566 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
567 Int_t module=relId[0];
569 Int_t column=relId[3];
570 Float_t decalibration = gime->CalibData()->GetADCchannelEmc(module,column,row);
571 Float_t energy = digit->GetEnergy() / decalibration;
572 digit->SetEnergy(energy);
574 //____________________________________________________________________________
575 Int_t AliPHOSDigitizer::DigitizeCPV(Float_t charge, Int_t absId)
577 // Returns digitized value of the CPV charge in a pad absId
579 AliPHOSGetter* gime = AliPHOSGetter::Instance();
581 if(!gime->CalibData()) {
582 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1); // use AliCDBManager's run number
583 gime->SetCalibData(cdb);
586 //Determine rel.position of the cell absId
588 gime->PHOSGeometry()->AbsToRelNumbering(absId,relId);
589 Int_t module=relId[0];
591 Int_t column=relId[3];
595 if(absId > fEmcCrystals){ //digitize CPV only
597 //reading calibration data for cell absId.
598 //If no calibration DB found, accept default values.
600 if(gime->CalibData()) {
601 fADCpedestalCpv = gime->CalibData()->GetADCpedestalCpv(module,column,row);
602 fADCchanelCpv = gime->CalibData()->GetADCchannelCpv( module,column,row);
605 channel = (Int_t) TMath::Ceil((charge - fADCpedestalCpv)/fADCchanelCpv) ;
606 if(channel > fNADCcpv ) channel = fNADCcpv ;
611 //____________________________________________________________________________
612 void AliPHOSDigitizer::Exec(Option_t *option)
614 // Steering method to process digitization for events
615 // in the range from fFirstEvent to fLastEvent.
616 // This range is optionally set by SetEventRange().
617 // if fLastEvent=-1, then process events until the end.
618 // by default fLastEvent = fFirstEvent (process only one event)
620 if (!fInit) { // to prevent overwrite existing file
621 AliError(Form("Give a version name different from %s",
622 fEventFolderName.Data() )) ;
626 if (strstr(option,"print")) {
631 // check the QA result for Hits and SDigits
632 AliQA * qa = AliQA::Instance(AliQA::kPHOS) ;
633 if ( qa->IsSet(AliQA::kPHOS, AliQA::kSIM, AliQA::kFATAL)) {
634 AliFatal("QA status in Hits and/or SDIGITS was Fatal") ;
635 } else if ( qa->IsSet(AliQA::kPHOS, AliQA::kSIM, AliQA::kERROR)) {
636 AliError("QA status in Hits and/or SDIGITS was Error") ;
637 } else if ( qa->IsSet(AliQA::kPHOS, AliQA::kSIM, AliQA::kWARNING) ) {
638 AliWarning("QA status in Hits and/or SDIGITS was Warning") ;
639 } else if ( qa->IsSet(AliQA::kPHOS, AliQA::kSIM, AliQA::kINFO) ) {
640 AliInfo("QA status in Hits and/or SDIGITS was Info") ;
643 if(strstr(option,"tim"))
644 gBenchmark->Start("PHOSDigitizer");
646 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
648 // Post Digitizer to the white board
649 gime->PostDigitizer(this) ;
651 if (fLastEvent == -1)
652 fLastEvent = gime->MaxEvent() - 1 ;
654 fLastEvent = fFirstEvent ;
656 Int_t nEvents = fLastEvent - fFirstEvent + 1;
660 for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
662 gime->Event(ievent,"S") ;
664 Digitize(ievent) ; //Add prepared SDigits to digits and add the noise
666 // //makes the quality assurance data
667 // if (GetQADataMaker()->IsCycleDone() ) {
668 // GetQADataMaker()->EndOfCycle(AliQA::kDIGITS) ;
669 // GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
671 // GetQADataMaker()->Exec(AliQA::kDIGITS, gime->Digits()) ;
672 // GetQADataMaker()->Increment() ;
676 if(strstr(option,"deb"))
679 //increment the total number of Digits per run
680 fDigitsInRun += gime->Digits()->GetEntriesFast() ;
683 // //Write the quality assurance data only after the last event
684 // if ( fEventCounter == gime->MaxEvent() ) {
685 // GetQADataMaker()->EndOfCycle(AliQA::kDIGITS) ;
686 // GetQADataMaker()->Finish() ;
689 gime->PhosLoader()->CleanDigitizer();
691 if(strstr(option,"tim")){
692 gBenchmark->Stop("PHOSDigitizer");
694 message = " took %f seconds for Digitizing %f seconds per event\n" ;
695 AliInfo(Form( message.Data(),
696 gBenchmark->GetCpuTime("PHOSDigitizer"),
697 gBenchmark->GetCpuTime("PHOSDigitizer")/nEvents ));
701 //____________________________________________________________________________
702 Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks) const
704 // Returns the shortest time among all time ticks
706 ticks->Sort() ; //Sort in accordance with times of ticks
708 AliPHOSTick * ctick = (AliPHOSTick *) it.Next() ;
709 Float_t time = ctick->CrossingTime(fTimeThreshold) ;
712 while((t=(AliPHOSTick*) it.Next())){
713 if(t->GetTime() < time) //This tick starts before crossing
718 time = ctick->CrossingTime(fTimeThreshold) ;
723 //____________________________________________________________________________
724 Bool_t AliPHOSDigitizer::Init()
726 // Makes all memory allocations
728 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
730 AliFatal(Form("Could not obtain the Getter object for file %s and event %s !",
731 GetTitle(), fEventFolderName.Data()));
735 const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
737 fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ;
739 TString opt("Digits") ;
740 if(gime->VersionExists(opt) ) {
741 AliError(Form("Give a version name different from %s",
742 fEventFolderName.Data() )) ;
747 fLastEvent = fFirstEvent ;
749 fInput = fManager->GetNinputs() ;
753 fInputFileNames = new TString[fInput] ;
754 fEventNames = new TString[fInput] ;
755 fInputFileNames[0] = GetTitle() ;
756 fEventNames[0] = fEventFolderName.Data() ;
758 for (index = 1 ; index < fInput ; index++) {
759 fInputFileNames[index] = dynamic_cast<AliStream*>(fManager->GetInputStream(index))->GetFileName(0);
760 TString tempo = fManager->GetInputFolderName(index) ;
761 fEventNames[index] = tempo.Remove(tempo.Length()-1) ; // strip of the stream number added by fManager
764 //to prevent cleaning of this object while GetEvent is called
765 gime->PhosLoader()->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
770 //____________________________________________________________________________
771 void AliPHOSDigitizer::InitParameters()
773 // Set initial parameters Digitizer
775 fPinNoise = 0.004 ; // [GeV]
776 fEMCDigitThreshold = 0.012 ; // [GeV]
777 fCPVNoise = 0.01; // [aux units]
778 fCPVDigitThreshold = 0.09 ; // [aux units]
779 fTimeResolution = 0.5e-9 ; // [sec]
780 fTimeSignalLength = 1.0e-9 ; // [sec]
782 fADCchanelEmc = 1.0; // Coefficient between real and measured energies in EMC
783 fADCpedestalEmc = 0. ; //
784 fNADCemc = (Int_t) TMath::Power(2,16) ; // number of channels in EMC ADC
786 fADCchanelCpv = 0.0012 ; // width of one ADC channel in CPV 'popugais'
787 fADCpedestalCpv = 0.012 ; //
788 fNADCcpv = (Int_t) TMath::Power(2,12); // number of channels in CPV ADC
790 // fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
791 fTimeThreshold = 0.001 ; // [GeV]
792 SetEventRange(0,-1) ;
796 //__________________________________________________________________
797 void AliPHOSDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
799 // Allows to produce digits by superimposing background and signal event.
800 // It is assumed, that headers file with SIGNAL events is opened in
802 // Sets the BACKGROUND event, with which the SIGNAL event is to be mixed
803 // Thus we avoid writing (changing) huge and expensive
804 // backgound files: all output will be writen into SIGNAL, i.e.
805 // opened in constructor file.
807 // One can open as many files to mix with as one needs.
808 // However only Sdigits with the same name (i.e. constructed with the same SDigitizer)
811 if( strcmp(fEventFolderName, "") == 0 )
815 Warning("MixWith", "Cannot use this method with AliRunDigitizer\n" ) ;
818 // looking for file which contains AliRun
819 if (gSystem->AccessPathName(alirunFileName)) {// file does not exist
820 AliError(Form("File %s does not exist!", alirunFileName.Data())) ;
823 // looking for the file which contains SDigits
824 AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
825 TString fileName( gime->GetSDigitsFileName() ) ;
826 if ( eventFolderName != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
827 fileName = fileName.ReplaceAll(".root", "") + "_" + eventFolderName + ".root" ;
828 if ( (gSystem->AccessPathName(fileName)) ) {
829 AliError(Form("The file %s does not exist!", fileName.Data())) ;
832 // need to increase the arrays
833 TString tempo = fInputFileNames[fInput-1] ;
834 delete [] fInputFileNames ;
835 fInputFileNames = new TString[fInput+1] ;
836 fInputFileNames[fInput-1] = tempo ;
838 tempo = fEventNames[fInput-1] ;
839 delete [] fEventNames ;
840 fEventNames = new TString[fInput+1] ;
841 fEventNames[fInput-1] = tempo ;
843 fInputFileNames[fInput] = alirunFileName ;
844 fEventNames[fInput] = eventFolderName ;
848 //__________________________________________________________________
849 void AliPHOSDigitizer::Print(const Option_t *)const
851 // Print Digitizer's parameters
852 AliInfo(Form("\n------------------- %s -------------", GetName() )) ;
853 if( strcmp(fEventFolderName.Data(), "") != 0 ){
854 printf(" Writing Digits to branch with title %s\n", fEventFolderName.Data()) ;
858 nStreams = GetNInputStreams() ;
863 for (index = 0 ; index < nStreams ; index++) {
864 TString tempo(fEventNames[index]) ;
866 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[index], tempo) ;
867 TString fileName( gime->GetSDigitsFileName() ) ;
868 if ( fEventNames[index] != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
869 fileName = fileName.ReplaceAll(".root", "") + "_" + fEventNames[index] + ".root" ;
870 printf ("Adding SDigits from %s %s\n", fInputFileNames[index].Data(), fileName.Data()) ;
872 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
873 printf("\nWriting digits to %s", gime->GetDigitsFileName().Data()) ;
875 printf("\nWith following parameters:\n") ;
876 printf(" Electronics noise in EMC (fPinNoise) = %f GeV\n", fPinNoise ) ;
877 printf(" Threshold in EMC (fEMCDigitThreshold) = %f GeV\n", fEMCDigitThreshold ) ;
878 printf(" Noise in CPV (fCPVNoise) = %f aux units\n", fCPVNoise ) ;
879 printf(" Threshold in CPV (fCPVDigitThreshold) = %f aux units\n",fCPVDigitThreshold ) ;
880 printf(" ---------------------------------------------------\n") ;
883 AliInfo(Form("AliPHOSDigitizer not initialized" )) ;
887 //__________________________________________________________________
888 void AliPHOSDigitizer::PrintDigits(Option_t * option)
890 // Print a table of digits
892 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
893 TClonesArray * digits = gime->Digits() ;
895 AliInfo(Form("%d", digits->GetEntriesFast())) ;
896 printf("\nevent %d", gAlice->GetEvNumber()) ;
897 printf("\n Number of entries in Digits list %d", digits->GetEntriesFast() ) ;
900 if(strstr(option,"all")||strstr(option,"EMC")){
902 AliPHOSDigit * digit;
903 printf("\nEMC digits (with primaries):\n") ;
904 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
905 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
907 for (index = 0 ; (index < digits->GetEntriesFast()) &&
908 (dynamic_cast<AliPHOSDigit *>(digits->At(index))->GetId() <= maxEmc) ; index++) {
909 digit = (AliPHOSDigit * ) digits->At(index) ;
910 if(digit->GetNprimary() == 0)
912 // printf("%6d %8d %6.5e %4d %2d :",
913 // digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ; // YVK
914 printf("%6d %.4f %6.5e %4d %2d :",
915 digit->GetId(), digit->GetEnergy(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
917 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
918 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
924 if(strstr(option,"all")||strstr(option,"CPV")){
926 //loop over CPV digits
927 AliPHOSDigit * digit;
928 printf("\nCPV digits (with primaries):\n") ;
929 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
930 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
932 for (index = 0 ; index < digits->GetEntriesFast(); index++) {
933 digit = (AliPHOSDigit * ) digits->At(index) ;
934 if(digit->GetId() > maxEmc){
935 printf("%6d %8d %4d %2d :",
936 digit->GetId(), digit->GetAmp(), digit->GetIndexInList(), digit->GetNprimary()) ;
938 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
939 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
948 //__________________________________________________________________
949 Float_t AliPHOSDigitizer::TimeOfNoise(void) const
950 { // Calculates the time signal generated by noise
951 //PH Info("TimeOfNoise", "Change me") ;
952 return gRandom->Rndm() * 1.28E-5;
955 //__________________________________________________________________
956 void AliPHOSDigitizer::Unload()
960 for(i = 1 ; i < fInput ; i++){
961 TString tempo(fEventNames[i]) ;
963 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
964 gime->PhosLoader()->UnloadSDigits() ;
967 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
968 gime->PhosLoader()->UnloadDigits() ;
971 //____________________________________________________________________________
972 void AliPHOSDigitizer::WriteDigits()
975 // Makes TreeD in the output file.
976 // Check if branch already exists:
977 // if yes, exit without writing: ROOT TTree does not support overwriting/updating of
978 // already existing branches.
979 // else creates branch with Digits, named "PHOS", title "...",
980 // and branch "AliPHOSDigitizer", with the same title to keep all the parameters
981 // and names of files, from which digits are made.
983 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
984 const TClonesArray * digits = gime->Digits() ;
985 TTree * treeD = gime->TreeD();
987 // -- create Digits branch
988 Int_t bufferSize = 32000 ;
989 TBranch * digitsBranch = treeD->Branch("PHOS","TClonesArray",&digits,bufferSize);
990 digitsBranch->SetTitle(fEventFolderName);
991 digitsBranch->Fill() ;
993 gime->WriteDigits("OVERWRITE");
994 gime->WriteDigitizer("OVERWRITE");