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.99 2007/09/30 17:08:20 schutz
22 * Introducing the notion of QA data acquisition cycle (needed by online)
24 * Revision 1.98 2007/09/26 14:22:17 cvetan
25 * 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.
27 * Revision 1.97 2007/08/07 14:12:03 kharlov
28 * Quality assurance added (Yves Schutz)
30 * Revision 1.96 2007/04/28 10:43:36 policheh
31 * Dead channels simulation: digit energy sets to 0.
33 * Revision 1.95 2007/04/10 07:20:52 kharlov
34 * Decalibration should use the same CDB as calibration in AliPHOSClusterizerv1
36 * Revision 1.94 2007/02/01 10:34:47 hristov
37 * Removing warnings on Solaris x86
39 * Revision 1.93 2006/10/17 13:17:01 kharlov
40 * Replace AliInfo by AliDebug
42 * Revision 1.92 2006/08/28 10:01:56 kharlov
43 * Effective C++ warnings fixed (Timur Pocheptsov)
45 * Revision 1.91 2006/04/29 20:25:30 hristov
46 * Decalibration is implemented (Yu.Kharlov)
48 * Revision 1.90 2006/04/22 10:30:17 hristov
49 * Add fEnergy to AliPHOSDigit and operate with EMC amplitude in energy units (Yu.Kharlov)
51 * Revision 1.89 2006/04/11 15:22:59 hristov
52 * run number in query set to -1: forces AliCDBManager to use its run number (A.Colla)
54 * Revision 1.88 2006/03/13 14:05:43 kharlov
55 * Calibration objects for EMC and CPV
57 * Revision 1.87 2005/08/24 15:33:49 kharlov
58 * Calibration data for raw digits
60 * Revision 1.86 2005/07/12 20:07:35 hristov
61 * Changes needed to run simulation and reconstrruction in the same AliRoot session
63 * Revision 1.85 2005/05/28 14:19:04 schutz
64 * Compilation warnings fixed by T.P.
68 //_________________________________________________________________________
69 //*-- Author : Dmitri Peressounko (SUBATECH & Kurchatov Institute)
70 //////////////////////////////////////////////////////////////////////////////
71 // This TTask performs digitization of Summable digits (in the PHOS case it is just
72 // the sum of contributions from all primary particles into a given cell).
73 // In addition it performs mixing of summable digits from different events.
74 // The name of the TTask is also the title of the branch that will contain
75 // the created SDigits
76 // The title of the TTAsk is the name of the file that contains the hits from
77 // which the SDigits are created
79 // For each event two branches are created in TreeD:
80 // "PHOS" - list of digits
81 // "AliPHOSDigitizer" - AliPHOSDigitizer with all parameters used in digitization
83 // Note, that one can set a title for new digits branch, and repeat digitization with
84 // another set of parameters.
87 // root[0] AliPHOSDigitizer * d = new AliPHOSDigitizer() ;
88 // root[1] d->ExecuteTask()
89 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
90 // //Digitizes SDigitis in all events found in file galice.root
92 // root[2] AliPHOSDigitizer * d1 = new AliPHOSDigitizer("galice1.root") ;
93 // // Will read sdigits from galice1.root
94 // root[3] d1->MixWith("galice2.root")
95 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
96 // // Reads another set of sdigits from galice2.root
97 // root[3] d1->MixWith("galice3.root")
98 // // Reads another set of sdigits from galice3.root
99 // root[4] d->ExecuteTask("deb timing")
100 // // Reads SDigits from files galice1.root, galice2.root ....
101 // // mixes them and stores produced Digits in file galice1.root
102 // // deb - prints number of produced digits
103 // // deb all - prints list of produced digits
104 // // timing - prints time used for digitization
107 // --- ROOT system ---
110 #include "TBenchmark.h"
113 // --- Standard library ---
115 // --- AliRoot header files ---
117 #include "AliRunDigitizer.h"
118 #include "AliPHOSDigit.h"
119 #include "AliPHOSGetter.h"
120 #include "AliPHOSDigitizer.h"
121 #include "AliPHOSSDigitizer.h"
122 #include "AliPHOSGeometry.h"
123 #include "AliPHOSTick.h"
124 #include "AliPHOSQADataMaker.h"
126 ClassImp(AliPHOSDigitizer)
129 //____________________________________________________________________________
130 AliPHOSDigitizer::AliPHOSDigitizer() :
136 fInputFileNames(0x0),
140 fEMCDigitThreshold(0.f),
142 fCPVDigitThreshold(0.f),
143 fTimeResolution(0.f),
145 fTimeSignalLength(0.f),
147 fADCpedestalEmc(0.f),
150 fADCpedestalCpv(0.f),
152 fEventFolderName(""),
160 fManager = 0 ; // We work in the standalong mode
163 //____________________________________________________________________________
164 AliPHOSDigitizer::AliPHOSDigitizer(TString alirunFileName,
165 TString eventFolderName):
166 AliDigitizer("PHOS"+AliConfig::Instance()->GetDigitizerTaskName(), alirunFileName),
167 fDefaultInit(kFALSE),
171 fInputFileNames(0x0),
175 fEMCDigitThreshold(0.f),
177 fCPVDigitThreshold(0.f),
178 fTimeResolution(0.f),
180 fTimeSignalLength(0.f),
182 fADCpedestalEmc(0.f),
185 fADCpedestalCpv(0.f),
187 fEventFolderName(eventFolderName),
196 fDefaultInit = kFALSE ;
197 fManager = 0 ; // We work in the standalong mode
198 //Initialize the quality assurance data maker only once
199 fQADM = new AliPHOSQADataMaker() ;
200 //FIXME: get the run number
203 GetQADataMaker()->Init(AliQA::kDIGITS, run, fgkCycles) ;
204 GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
207 //____________________________________________________________________________
208 AliPHOSDigitizer::AliPHOSDigitizer(const AliPHOSDigitizer & d) :
210 fDefaultInit(d.fDefaultInit),
211 fDigitsInRun(d.fDigitsInRun),
214 fInputFileNames(0x0),//?
216 fEmcCrystals(d.fEmcCrystals),
217 fPinNoise(d.fPinNoise),
218 fEMCDigitThreshold(d.fEMCDigitThreshold),
219 fCPVNoise(d.fCPVNoise),
220 fCPVDigitThreshold(d.fCPVDigitThreshold),
221 fTimeResolution(d.fTimeResolution),
222 fTimeThreshold(d.fTimeThreshold),
223 fTimeSignalLength(d.fTimeSignalLength),
224 fADCchanelEmc(d.fADCchanelEmc),
225 fADCpedestalEmc(d.fADCpedestalEmc),
226 fNADCemc(d.fNADCemc),
227 fADCchanelCpv(d.fADCchanelCpv),
228 fADCpedestalCpv(d.fADCpedestalCpv),
229 fNADCcpv(d.fNADCcpv),
230 fEventFolderName(d.fEventFolderName),
231 fFirstEvent(d.fFirstEvent),
232 fLastEvent(d.fLastEvent),
238 SetName(d.GetName()) ;
239 SetTitle(d.GetTitle()) ;
240 //Initialize the quality assurance data maker only once
241 //FIXME: get the run number
244 GetQADataMaker()->Init(AliQA::kDIGITS, run, fgkCycles) ;
245 GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
248 //____________________________________________________________________________
249 AliPHOSDigitizer::AliPHOSDigitizer(AliRunDigitizer * rd) :
250 AliDigitizer(rd,"PHOS"+AliConfig::Instance()->GetDigitizerTaskName()),
251 fDefaultInit(kFALSE),
255 fInputFileNames(0x0),
259 fEMCDigitThreshold(0.f),
261 fCPVDigitThreshold(0.f),
262 fTimeResolution(0.f),
264 fTimeSignalLength(0.f),
266 fADCpedestalEmc(0.f),
269 fADCpedestalCpv(0.f),
271 fEventFolderName(fManager->GetInputFolderName(0)),
278 // ctor Init() is called by RunDigitizer
280 SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
282 fDefaultInit = kFALSE ;
283 //Initialize the quality assurance data maker only once
284 fQADM = new AliPHOSQADataMaker() ;
285 //FIXME: get the run number
288 GetQADataMaker()->Init(AliQA::kDIGITS, run) ;
289 GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
292 //____________________________________________________________________________
293 AliPHOSDigitizer::~AliPHOSDigitizer()
295 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
297 // Clean Digitizer from the white board
298 gime->PhosLoader()->CleanDigitizer() ;
300 delete [] fInputFileNames ;
301 delete [] fEventNames ;
307 //____________________________________________________________________________
308 void AliPHOSDigitizer::Digitize(Int_t event)
311 // Makes the digitization of the collected summable digits.
312 // It first creates the array of all PHOS modules
313 // filled with noise (different for EMC, and CPV) and
314 // then adds contributions from SDigits.
315 // This design avoids scanning over the list of digits to add
316 // contribution to new SDigits only.
318 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
319 Int_t ReadEvent = event ;
321 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetCurrentEventNumber() ;
322 AliDebug(1,Form("Adding event %d from input stream 0 %s %s",
323 ReadEvent, GetTitle(), fEventFolderName.Data())) ;
324 gime->Event(ReadEvent, "S") ;
325 TClonesArray * digits = gime->Digits() ;
328 const AliPHOSGeometry *geom = gime->PHOSGeometry() ;
329 //Making digits with noise, first EMC
330 Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
335 nCPV = nEMC + geom->GetNumberOfCPVPadsZ() * geom->GetNumberOfCPVPadsPhi() * geom->GetNModules() ;
337 digits->Expand(nCPV) ;
339 // get first the sdigitizer from the tasks list
340 if ( !gime->SDigitizer() )
341 gime->LoadSDigitizer();
342 AliPHOSSDigitizer * sDigitizer = gime->SDigitizer();
345 AliFatal(Form("SDigitizer with name %s %s not found",
346 GetTitle(), fEventFolderName.Data() )) ;
348 //take all the inputs to add together and load the SDigits
349 TObjArray * sdigArray = new TObjArray(fInput) ;
350 sdigArray->AddAt(gime->SDigits(), 0) ;
352 for(i = 1 ; i < fInput ; i++){
353 TString tempo(fEventNames[i]) ;
355 AliPHOSGetter * gime1 = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
357 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(i))->GetCurrentEventNumber() ;
358 AliInfo(Form("Adding event %d from input stream %d %s %s",
359 ReadEvent, i, fInputFileNames[i].Data(), tempo.Data())) ;
360 gime1->Event(ReadEvent,"S");
361 sdigArray->AddAt(gime1->SDigits(), i) ;
364 //Find the first crystal with signal
365 Int_t nextSig = 200000 ;
366 TClonesArray * sdigits ;
367 for(i = 0 ; i < fInput ; i++){
368 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
369 if ( !sdigits->GetEntriesFast() )
371 Int_t curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(0))->GetId() ;
372 if(curNext < nextSig)
376 TArrayI index(fInput) ;
377 index.Reset() ; //Set all indexes to zero
379 AliPHOSDigit * digit ;
380 AliPHOSDigit * curSDigit ;
382 TClonesArray * ticks = new TClonesArray("AliPHOSTick",1000) ;
384 //Put Noise contribution
385 for(absID = 1 ; absID <= nEMC ; absID++){
386 Float_t noise = gRandom->Gaus(0., fPinNoise) ;
387 // YVK: do not digitize amplitudes for EMC
388 // new((*digits)[absID-1]) AliPHOSDigit( -1, absID, sDigitizer->Digitize(noise), TimeOfNoise() ) ;
389 new((*digits)[absID-1]) AliPHOSDigit( -1, absID, noise, TimeOfNoise() ) ;
390 //look if we have to add signal?
391 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
394 //Add SDigits from all inputs
397 Float_t a = digit->GetEnergy() ;
398 Float_t b = TMath::Abs( a / fTimeSignalLength) ;
399 //Mark the beginning of the signal
400 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime(),0, b);
401 //Mark the end of the signal
402 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime()+fTimeSignalLength, -a, -b);
405 for(i = 0 ; i < fInput ; i++){
406 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
407 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
410 //May be several digits will contribute from the same input
411 while(curSDigit && curSDigit->GetId() == absID){
412 //Shift primary to separate primaries belonging different inputs
413 Int_t primaryoffset ;
415 primaryoffset = fManager->GetMask(i) ;
417 primaryoffset = 10000000*i ;
418 curSDigit->ShiftPrimary(primaryoffset) ;
420 a = curSDigit->GetEnergy() ;
421 b = a /fTimeSignalLength ;
422 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime(),0, b);
423 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b);
425 *digit += *curSDigit ; //add energies
428 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
429 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
435 //calculate and set time
436 Float_t time = FrontEdgeTime(ticks) ;
437 digit->SetTime(time) ;
439 //Find next signal module
441 for(i = 0 ; i < fInput ; i++){
442 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
443 Int_t curNext = nextSig ;
444 if(sdigits->GetEntriesFast() > index[i] ){
445 curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(index[i]))->GetId() ;
447 if(curNext < nextSig) nextSig = curNext ;
455 //Now CPV digits (different noise and no timing)
456 for(absID = nEMC+1; absID <= nCPV; absID++){
457 Float_t noise = gRandom->Gaus(0., fCPVNoise) ;
458 new((*digits)[absID-1]) AliPHOSDigit( -1,absID,sDigitizer->Digitize(noise), TimeOfNoise() ) ;
459 //look if we have to add signal?
461 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
462 //Add SDigits from all inputs
463 for(i = 0 ; i < fInput ; i++){
464 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
465 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
469 //May be several digits will contribute from the same input
470 while(curSDigit && curSDigit->GetId() == absID){
471 //Shift primary to separate primaries belonging different inputs
472 Int_t primaryoffset ;
474 primaryoffset = fManager->GetMask(i) ;
476 primaryoffset = 10000000*i ;
477 curSDigit->ShiftPrimary(primaryoffset) ;
480 *digit += *curSDigit ;
482 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
483 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i]) ) ;
489 //Find next signal module
491 for(i = 0 ; i < fInput ; i++){
492 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
493 Int_t curNext = nextSig ;
494 if(sdigits->GetEntriesFast() > index[i] )
495 curNext = dynamic_cast<AliPHOSDigit *>( sdigits->At(index[i]) )->GetId() ;
496 if(curNext < nextSig) nextSig = curNext ;
502 delete sdigArray ; //We should not delete its contents
504 //remove digits below thresholds
505 for(i = 0 ; i < nEMC ; i++){
506 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
507 DecalibrateEMC(digit);
508 if(digit->GetEnergy() < fEMCDigitThreshold)
509 digits->RemoveAt(i) ;
511 digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ;
515 for(i = nEMC; i < nCPV ; i++)
516 // if( sDigitizer->Calibrate( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetAmp() ) < fCPVDigitThreshold )
517 if( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetEnergy() < fCPVDigitThreshold )
518 digits->RemoveAt(i) ;
522 Int_t ndigits = digits->GetEntriesFast() ;
523 digits->Expand(ndigits) ;
525 //Set indexes in list of digits and make true digitization of the energy
526 for (i = 0 ; i < ndigits ; i++) {
527 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
528 digit->SetIndexInList(i) ;
529 if(digit->GetId() > fEmcCrystals){ //digitize CPV only
530 digit->SetAmp(DigitizeCPV(digit->GetEnergy(),digit->GetId()) ) ;
536 //set amplitudes in bad channels to zero
537 for(i = 0 ; i <digits->GetEntries(); i++){
538 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
539 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
540 if(relId[1] == 0) // Emc
541 if(gime->CalibData()->IsBadChannelEmc(relId[0],relId[3],relId[2])) digit->SetEnergy(0.);
546 //____________________________________________________________________________
547 void AliPHOSDigitizer::DecalibrateEMC(AliPHOSDigit *digit)
549 // Decalibrate EMC digit, i.e. change its energy by a factor read from CDB
551 AliPHOSGetter* gime = AliPHOSGetter::Instance();
553 if(!gime->CalibData()) {
554 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1);
555 gime->SetCalibData(cdb);
558 //Determine rel.position of the cell absolute ID
560 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
561 Int_t module=relId[0];
563 Int_t column=relId[3];
564 Float_t decalibration = gime->CalibData()->GetADCchannelEmc(module,column,row);
565 Float_t energy = digit->GetEnergy() / decalibration;
566 digit->SetEnergy(energy);
568 //____________________________________________________________________________
569 Int_t AliPHOSDigitizer::DigitizeCPV(Float_t charge, Int_t absId)
571 // Returns digitized value of the CPV charge in a pad absId
573 AliPHOSGetter* gime = AliPHOSGetter::Instance();
575 if(!gime->CalibData()) {
576 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1); // use AliCDBManager's run number
577 gime->SetCalibData(cdb);
580 //Determine rel.position of the cell absId
582 gime->PHOSGeometry()->AbsToRelNumbering(absId,relId);
583 Int_t module=relId[0];
585 Int_t column=relId[3];
589 if(absId > fEmcCrystals){ //digitize CPV only
591 //reading calibration data for cell absId.
592 //If no calibration DB found, accept default values.
594 if(gime->CalibData()) {
595 fADCpedestalCpv = gime->CalibData()->GetADCpedestalCpv(module,column,row);
596 fADCchanelCpv = gime->CalibData()->GetADCchannelCpv( module,column,row);
599 channel = (Int_t) TMath::Ceil((charge - fADCpedestalCpv)/fADCchanelCpv) ;
600 if(channel > fNADCcpv ) channel = fNADCcpv ;
605 //____________________________________________________________________________
606 void AliPHOSDigitizer::Exec(Option_t *option)
608 // Steering method to process digitization for events
609 // in the range from fFirstEvent to fLastEvent.
610 // This range is optionally set by SetEventRange().
611 // if fLastEvent=-1, then process events until the end.
612 // by default fLastEvent = fFirstEvent (process only one event)
614 if (!fInit) { // to prevent overwrite existing file
615 AliError(Form("Give a version name different from %s",
616 fEventFolderName.Data() )) ;
620 if (strstr(option,"print")) {
625 if(strstr(option,"tim"))
626 gBenchmark->Start("PHOSDigitizer");
628 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
630 // Post Digitizer to the white board
631 gime->PostDigitizer(this) ;
633 if (fLastEvent == -1)
634 fLastEvent = gime->MaxEvent() - 1 ;
636 fLastEvent = fFirstEvent ;
638 Int_t nEvents = fLastEvent - fFirstEvent + 1;
642 for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
644 gime->Event(ievent,"S") ;
646 Digitize(ievent) ; //Add prepared SDigits to digits and add the noise
648 //makes the quality assurance data
649 if (GetQADataMaker()->IsCycleDone() ) {
650 GetQADataMaker()->EndOfCycle(AliQA::kDIGITS) ;
651 GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
653 GetQADataMaker()->Exec(AliQA::kDIGITS, gime->Digits()) ;
654 GetQADataMaker()->Increment() ;
658 if(strstr(option,"deb"))
661 //increment the total number of Digits per run
662 fDigitsInRun += gime->Digits()->GetEntriesFast() ;
665 //Write the quality assurance data only after the last event
666 if ( fEventCounter == gime->MaxEvent() ) {
667 GetQADataMaker()->EndOfCycle(AliQA::kDIGITS) ;
668 GetQADataMaker()->Finish(AliQA::kDIGITS) ;
671 gime->PhosLoader()->CleanDigitizer();
673 if(strstr(option,"tim")){
674 gBenchmark->Stop("PHOSDigitizer");
676 message = " took %f seconds for Digitizing %f seconds per event\n" ;
677 AliInfo(Form( message.Data(),
678 gBenchmark->GetCpuTime("PHOSDigitizer"),
679 gBenchmark->GetCpuTime("PHOSDigitizer")/nEvents ));
683 //____________________________________________________________________________
684 Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks) const
686 // Returns the shortest time among all time ticks
688 ticks->Sort() ; //Sort in accordance with times of ticks
690 AliPHOSTick * ctick = (AliPHOSTick *) it.Next() ;
691 Float_t time = ctick->CrossingTime(fTimeThreshold) ;
694 while((t=(AliPHOSTick*) it.Next())){
695 if(t->GetTime() < time) //This tick starts before crossing
700 time = ctick->CrossingTime(fTimeThreshold) ;
705 //____________________________________________________________________________
706 Bool_t AliPHOSDigitizer::Init()
708 // Makes all memory allocations
710 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
712 AliFatal(Form("Could not obtain the Getter object for file %s and event %s !",
713 GetTitle(), fEventFolderName.Data()));
717 const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
719 fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ;
721 TString opt("Digits") ;
722 if(gime->VersionExists(opt) ) {
723 AliError(Form("Give a version name different from %s",
724 fEventFolderName.Data() )) ;
729 fLastEvent = fFirstEvent ;
731 fInput = fManager->GetNinputs() ;
735 fInputFileNames = new TString[fInput] ;
736 fEventNames = new TString[fInput] ;
737 fInputFileNames[0] = GetTitle() ;
738 fEventNames[0] = fEventFolderName.Data() ;
740 for (index = 1 ; index < fInput ; index++) {
741 fInputFileNames[index] = dynamic_cast<AliStream*>(fManager->GetInputStream(index))->GetFileName(0);
742 TString tempo = fManager->GetInputFolderName(index) ;
743 fEventNames[index] = tempo.Remove(tempo.Length()-1) ; // strip of the stream number added by fManager
746 //to prevent cleaning of this object while GetEvent is called
747 gime->PhosLoader()->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
752 //____________________________________________________________________________
753 void AliPHOSDigitizer::InitParameters()
755 // Set initial parameters Digitizer
757 fPinNoise = 0.004 ; // [GeV]
758 fEMCDigitThreshold = 0.012 ; // [GeV]
759 fCPVNoise = 0.01; // [aux units]
760 fCPVDigitThreshold = 0.09 ; // [aux units]
761 fTimeResolution = 0.5e-9 ; // [sec]
762 fTimeSignalLength = 1.0e-9 ; // [sec]
764 fADCchanelEmc = 1.0; // Coefficient between real and measured energies in EMC
765 fADCpedestalEmc = 0. ; //
766 fNADCemc = (Int_t) TMath::Power(2,16) ; // number of channels in EMC ADC
768 fADCchanelCpv = 0.0012 ; // width of one ADC channel in CPV 'popugais'
769 fADCpedestalCpv = 0.012 ; //
770 fNADCcpv = (Int_t) TMath::Power(2,12); // number of channels in CPV ADC
772 // fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
773 fTimeThreshold = 0.001 ; // [GeV]
774 SetEventRange(0,-1) ;
778 //__________________________________________________________________
779 void AliPHOSDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
781 // Allows to produce digits by superimposing background and signal event.
782 // It is assumed, that headers file with SIGNAL events is opened in
784 // Sets the BACKGROUND event, with which the SIGNAL event is to be mixed
785 // Thus we avoid writing (changing) huge and expensive
786 // backgound files: all output will be writen into SIGNAL, i.e.
787 // opened in constructor file.
789 // One can open as many files to mix with as one needs.
790 // However only Sdigits with the same name (i.e. constructed with the same SDigitizer)
793 if( strcmp(fEventFolderName, "") == 0 )
797 Warning("MixWith", "Cannot use this method with AliRunDigitizer\n" ) ;
800 // looking for file which contains AliRun
801 if (gSystem->AccessPathName(alirunFileName)) {// file does not exist
802 AliError(Form("File %s does not exist!", alirunFileName.Data())) ;
805 // looking for the file which contains SDigits
806 AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
807 TString fileName( gime->GetSDigitsFileName() ) ;
808 if ( eventFolderName != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
809 fileName = fileName.ReplaceAll(".root", "") + "_" + eventFolderName + ".root" ;
810 if ( (gSystem->AccessPathName(fileName)) ) {
811 AliError(Form("The file %s does not exist!", fileName.Data())) ;
814 // need to increase the arrays
815 TString tempo = fInputFileNames[fInput-1] ;
816 delete [] fInputFileNames ;
817 fInputFileNames = new TString[fInput+1] ;
818 fInputFileNames[fInput-1] = tempo ;
820 tempo = fEventNames[fInput-1] ;
821 delete [] fEventNames ;
822 fEventNames = new TString[fInput+1] ;
823 fEventNames[fInput-1] = tempo ;
825 fInputFileNames[fInput] = alirunFileName ;
826 fEventNames[fInput] = eventFolderName ;
830 //__________________________________________________________________
831 void AliPHOSDigitizer::Print(const Option_t *)const
833 // Print Digitizer's parameters
834 AliInfo(Form("\n------------------- %s -------------", GetName() )) ;
835 if( strcmp(fEventFolderName.Data(), "") != 0 ){
836 printf(" Writing Digits to branch with title %s\n", fEventFolderName.Data()) ;
840 nStreams = GetNInputStreams() ;
845 for (index = 0 ; index < nStreams ; index++) {
846 TString tempo(fEventNames[index]) ;
848 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[index], tempo) ;
849 TString fileName( gime->GetSDigitsFileName() ) ;
850 if ( fEventNames[index] != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
851 fileName = fileName.ReplaceAll(".root", "") + "_" + fEventNames[index] + ".root" ;
852 printf ("Adding SDigits from %s %s\n", fInputFileNames[index].Data(), fileName.Data()) ;
854 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
855 printf("\nWriting digits to %s", gime->GetDigitsFileName().Data()) ;
857 printf("\nWith following parameters:\n") ;
858 printf(" Electronics noise in EMC (fPinNoise) = %f GeV\n", fPinNoise ) ;
859 printf(" Threshold in EMC (fEMCDigitThreshold) = %f GeV\n", fEMCDigitThreshold ) ;
860 printf(" Noise in CPV (fCPVNoise) = %f aux units\n", fCPVNoise ) ;
861 printf(" Threshold in CPV (fCPVDigitThreshold) = %f aux units\n",fCPVDigitThreshold ) ;
862 printf(" ---------------------------------------------------\n") ;
865 AliInfo(Form("AliPHOSDigitizer not initialized" )) ;
869 //__________________________________________________________________
870 void AliPHOSDigitizer::PrintDigits(Option_t * option)
872 // Print a table of digits
874 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
875 TClonesArray * digits = gime->Digits() ;
877 AliInfo(Form("%d", digits->GetEntriesFast())) ;
878 printf("\nevent %d", gAlice->GetEvNumber()) ;
879 printf("\n Number of entries in Digits list %d", digits->GetEntriesFast() ) ;
882 if(strstr(option,"all")||strstr(option,"EMC")){
884 AliPHOSDigit * digit;
885 printf("\nEMC digits (with primaries):\n") ;
886 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
887 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
889 for (index = 0 ; (index < digits->GetEntriesFast()) &&
890 (dynamic_cast<AliPHOSDigit *>(digits->At(index))->GetId() <= maxEmc) ; index++) {
891 digit = (AliPHOSDigit * ) digits->At(index) ;
892 if(digit->GetNprimary() == 0)
894 // printf("%6d %8d %6.5e %4d %2d :",
895 // digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ; // YVK
896 printf("%6d %.4f %6.5e %4d %2d :",
897 digit->GetId(), digit->GetEnergy(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
899 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
900 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
906 if(strstr(option,"all")||strstr(option,"CPV")){
908 //loop over CPV digits
909 AliPHOSDigit * digit;
910 printf("\nCPV digits (with primaries):\n") ;
911 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
912 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
914 for (index = 0 ; index < digits->GetEntriesFast(); index++) {
915 digit = (AliPHOSDigit * ) digits->At(index) ;
916 if(digit->GetId() > maxEmc){
917 printf("%6d %8d %4d %2d :",
918 digit->GetId(), digit->GetAmp(), digit->GetIndexInList(), digit->GetNprimary()) ;
920 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
921 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
930 //__________________________________________________________________
931 Float_t AliPHOSDigitizer::TimeOfNoise(void) const
932 { // Calculates the time signal generated by noise
933 //PH Info("TimeOfNoise", "Change me") ;
934 return gRandom->Rndm() * 1.28E-5;
937 //__________________________________________________________________
938 void AliPHOSDigitizer::Unload()
942 for(i = 1 ; i < fInput ; i++){
943 TString tempo(fEventNames[i]) ;
945 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
946 gime->PhosLoader()->UnloadSDigits() ;
949 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
950 gime->PhosLoader()->UnloadDigits() ;
953 //____________________________________________________________________________
954 void AliPHOSDigitizer::WriteDigits()
957 // Makes TreeD in the output file.
958 // Check if branch already exists:
959 // if yes, exit without writing: ROOT TTree does not support overwriting/updating of
960 // already existing branches.
961 // else creates branch with Digits, named "PHOS", title "...",
962 // and branch "AliPHOSDigitizer", with the same title to keep all the parameters
963 // and names of files, from which digits are made.
965 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
966 const TClonesArray * digits = gime->Digits() ;
967 TTree * treeD = gime->TreeD();
969 // -- create Digits branch
970 Int_t bufferSize = 32000 ;
971 TBranch * digitsBranch = treeD->Branch("PHOS","TClonesArray",&digits,bufferSize);
972 digitsBranch->SetTitle(fEventFolderName);
973 digitsBranch->Fill() ;
975 gime->WriteDigits("OVERWRITE");
976 gime->WriteDigitizer("OVERWRITE");