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.95 2007/04/10 07:20:52 kharlov
22 * Decalibration should use the same CDB as calibration in AliPHOSClusterizerv1
24 * Revision 1.94 2007/02/01 10:34:47 hristov
25 * Removing warnings on Solaris x86
27 * Revision 1.93 2006/10/17 13:17:01 kharlov
28 * Replace AliInfo by AliDebug
30 * Revision 1.92 2006/08/28 10:01:56 kharlov
31 * Effective C++ warnings fixed (Timur Pocheptsov)
33 * Revision 1.91 2006/04/29 20:25:30 hristov
34 * Decalibration is implemented (Yu.Kharlov)
36 * Revision 1.90 2006/04/22 10:30:17 hristov
37 * Add fEnergy to AliPHOSDigit and operate with EMC amplitude in energy units (Yu.Kharlov)
39 * Revision 1.89 2006/04/11 15:22:59 hristov
40 * run number in query set to -1: forces AliCDBManager to use its run number (A.Colla)
42 * Revision 1.88 2006/03/13 14:05:43 kharlov
43 * Calibration objects for EMC and CPV
45 * Revision 1.87 2005/08/24 15:33:49 kharlov
46 * Calibration data for raw digits
48 * Revision 1.86 2005/07/12 20:07:35 hristov
49 * Changes needed to run simulation and reconstrruction in the same AliRoot session
51 * Revision 1.85 2005/05/28 14:19:04 schutz
52 * Compilation warnings fixed by T.P.
56 //_________________________________________________________________________
57 //*-- Author : Dmitri Peressounko (SUBATECH & Kurchatov Institute)
58 //////////////////////////////////////////////////////////////////////////////
59 // This TTask performs digitization of Summable digits (in the PHOS case it is just
60 // the sum of contributions from all primary particles into a given cell).
61 // In addition it performs mixing of summable digits from different events.
62 // The name of the TTask is also the title of the branch that will contain
63 // the created SDigits
64 // The title of the TTAsk is the name of the file that contains the hits from
65 // which the SDigits are created
67 // For each event two branches are created in TreeD:
68 // "PHOS" - list of digits
69 // "AliPHOSDigitizer" - AliPHOSDigitizer with all parameters used in digitization
71 // Note, that one can set a title for new digits branch, and repeat digitization with
72 // another set of parameters.
75 // root[0] AliPHOSDigitizer * d = new AliPHOSDigitizer() ;
76 // root[1] d->ExecuteTask()
77 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
78 // //Digitizes SDigitis in all events found in file galice.root
80 // root[2] AliPHOSDigitizer * d1 = new AliPHOSDigitizer("galice1.root") ;
81 // // Will read sdigits from galice1.root
82 // root[3] d1->MixWith("galice2.root")
83 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
84 // // Reads another set of sdigits from galice2.root
85 // root[3] d1->MixWith("galice3.root")
86 // // Reads another set of sdigits from galice3.root
87 // root[4] d->ExecuteTask("deb timing")
88 // // Reads SDigits from files galice1.root, galice2.root ....
89 // // mixes them and stores produced Digits in file galice1.root
90 // // deb - prints number of produced digits
91 // // deb all - prints list of produced digits
92 // // timing - prints time used for digitization
95 // --- ROOT system ---
98 #include "TBenchmark.h"
101 // --- Standard library ---
103 // --- AliRoot header files ---
105 #include "AliRunDigitizer.h"
106 #include "AliPHOSDigit.h"
107 #include "AliPHOSGetter.h"
108 #include "AliPHOSDigitizer.h"
109 #include "AliPHOSSDigitizer.h"
110 #include "AliPHOSGeometry.h"
111 #include "AliPHOSTick.h"
113 ClassImp(AliPHOSDigitizer)
116 //____________________________________________________________________________
117 AliPHOSDigitizer::AliPHOSDigitizer() :
123 fInputFileNames(0x0),
127 fEMCDigitThreshold(0.f),
129 fCPVDigitThreshold(0.f),
130 fTimeResolution(0.f),
132 fTimeSignalLength(0.f),
134 fADCpedestalEmc(0.f),
137 fADCpedestalCpv(0.f),
139 fEventFolderName(""),
145 fManager = 0 ; // We work in the standalong mode
148 //____________________________________________________________________________
149 AliPHOSDigitizer::AliPHOSDigitizer(TString alirunFileName,
150 TString eventFolderName):
151 AliDigitizer("PHOS"+AliConfig::Instance()->GetDigitizerTaskName(), alirunFileName),
152 fDefaultInit(kFALSE),
156 fInputFileNames(0x0),
160 fEMCDigitThreshold(0.f),
162 fCPVDigitThreshold(0.f),
163 fTimeResolution(0.f),
165 fTimeSignalLength(0.f),
167 fADCpedestalEmc(0.f),
170 fADCpedestalCpv(0.f),
172 fEventFolderName(eventFolderName),
179 fDefaultInit = kFALSE ;
180 fManager = 0 ; // We work in the standalong mode
183 //____________________________________________________________________________
184 AliPHOSDigitizer::AliPHOSDigitizer(const AliPHOSDigitizer & d) :
186 fDefaultInit(d.fDefaultInit),
187 fDigitsInRun(d.fDigitsInRun),
190 fInputFileNames(0x0),//?
192 fEmcCrystals(d.fEmcCrystals),
193 fPinNoise(d.fPinNoise),
194 fEMCDigitThreshold(d.fEMCDigitThreshold),
195 fCPVNoise(d.fCPVNoise),
196 fCPVDigitThreshold(d.fCPVDigitThreshold),
197 fTimeResolution(d.fTimeResolution),
198 fTimeThreshold(d.fTimeThreshold),
199 fTimeSignalLength(d.fTimeSignalLength),
200 fADCchanelEmc(d.fADCchanelEmc),
201 fADCpedestalEmc(d.fADCpedestalEmc),
202 fNADCemc(d.fNADCemc),
203 fADCchanelCpv(d.fADCchanelCpv),
204 fADCpedestalCpv(d.fADCpedestalCpv),
205 fNADCcpv(d.fNADCcpv),
206 fEventFolderName(d.fEventFolderName),
207 fFirstEvent(d.fFirstEvent),
208 fLastEvent(d.fLastEvent)
211 SetName(d.GetName()) ;
212 SetTitle(d.GetTitle()) ;
215 //____________________________________________________________________________
216 AliPHOSDigitizer::AliPHOSDigitizer(AliRunDigitizer * rd) :
217 AliDigitizer(rd,"PHOS"+AliConfig::Instance()->GetDigitizerTaskName()),
218 fDefaultInit(kFALSE),
222 fInputFileNames(0x0),
226 fEMCDigitThreshold(0.f),
228 fCPVDigitThreshold(0.f),
229 fTimeResolution(0.f),
231 fTimeSignalLength(0.f),
233 fADCpedestalEmc(0.f),
236 fADCpedestalCpv(0.f),
238 fEventFolderName(fManager->GetInputFolderName(0)),
242 // ctor Init() is called by RunDigitizer
244 SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
246 fDefaultInit = kFALSE ;
249 //____________________________________________________________________________
250 AliPHOSDigitizer::~AliPHOSDigitizer()
252 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
254 // Clean Digitizer from the white board
255 gime->PhosLoader()->CleanDigitizer() ;
257 delete [] fInputFileNames ;
258 delete [] fEventNames ;
262 //____________________________________________________________________________
263 void AliPHOSDigitizer::Digitize(Int_t event)
266 // Makes the digitization of the collected summable digits.
267 // It first creates the array of all PHOS modules
268 // filled with noise (different for EMC, and CPV) and
269 // then adds contributions from SDigits.
270 // This design avoids scanning over the list of digits to add
271 // contribution to new SDigits only.
273 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
274 Int_t ReadEvent = event ;
276 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetCurrentEventNumber() ;
277 AliDebug(1,Form("Adding event %d from input stream 0 %s %s",
278 ReadEvent, GetTitle(), fEventFolderName.Data())) ;
279 gime->Event(ReadEvent, "S") ;
280 TClonesArray * digits = gime->Digits() ;
283 const AliPHOSGeometry *geom = gime->PHOSGeometry() ;
284 //Making digits with noise, first EMC
285 Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
290 nCPV = nEMC + geom->GetNumberOfCPVPadsZ() * geom->GetNumberOfCPVPadsPhi() * geom->GetNModules() ;
292 digits->Expand(nCPV) ;
294 // get first the sdigitizer from the tasks list
295 if ( !gime->SDigitizer() )
296 gime->LoadSDigitizer();
297 AliPHOSSDigitizer * sDigitizer = gime->SDigitizer();
300 AliFatal(Form("SDigitizer with name %s %s not found",
301 GetTitle(), fEventFolderName.Data() )) ;
303 //take all the inputs to add together and load the SDigits
304 TObjArray * sdigArray = new TObjArray(fInput) ;
305 sdigArray->AddAt(gime->SDigits(), 0) ;
307 for(i = 1 ; i < fInput ; i++){
308 TString tempo(fEventNames[i]) ;
310 AliPHOSGetter * gime1 = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
312 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(i))->GetCurrentEventNumber() ;
313 AliInfo(Form("Adding event %d from input stream %d %s %s",
314 ReadEvent, i, fInputFileNames[i].Data(), tempo.Data())) ;
315 gime1->Event(ReadEvent,"S");
316 sdigArray->AddAt(gime1->SDigits(), i) ;
319 //Find the first crystal with signal
320 Int_t nextSig = 200000 ;
321 TClonesArray * sdigits ;
322 for(i = 0 ; i < fInput ; i++){
323 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
324 if ( !sdigits->GetEntriesFast() )
326 Int_t curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(0))->GetId() ;
327 if(curNext < nextSig)
331 TArrayI index(fInput) ;
332 index.Reset() ; //Set all indexes to zero
334 AliPHOSDigit * digit ;
335 AliPHOSDigit * curSDigit ;
337 TClonesArray * ticks = new TClonesArray("AliPHOSTick",1000) ;
339 //Put Noise contribution
340 for(absID = 1 ; absID <= nEMC ; absID++){
341 Float_t noise = gRandom->Gaus(0., fPinNoise) ;
342 // YVK: do not digitize amplitudes for EMC
343 // new((*digits)[absID-1]) AliPHOSDigit( -1, absID, sDigitizer->Digitize(noise), TimeOfNoise() ) ;
344 new((*digits)[absID-1]) AliPHOSDigit( -1, absID, noise, TimeOfNoise() ) ;
345 //look if we have to add signal?
346 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
349 //Add SDigits from all inputs
352 Float_t a = digit->GetEnergy() ;
353 Float_t b = TMath::Abs( a / fTimeSignalLength) ;
354 //Mark the beginning of the signal
355 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime(),0, b);
356 //Mark the end of the signal
357 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime()+fTimeSignalLength, -a, -b);
360 for(i = 0 ; i < fInput ; i++){
361 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
362 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
365 //May be several digits will contribute from the same input
366 while(curSDigit && curSDigit->GetId() == absID){
367 //Shift primary to separate primaries belonging different inputs
368 Int_t primaryoffset ;
370 primaryoffset = fManager->GetMask(i) ;
372 primaryoffset = 10000000*i ;
373 curSDigit->ShiftPrimary(primaryoffset) ;
375 a = curSDigit->GetEnergy() ;
376 b = a /fTimeSignalLength ;
377 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime(),0, b);
378 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b);
380 *digit += *curSDigit ; //add energies
383 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
384 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
390 //calculate and set time
391 Float_t time = FrontEdgeTime(ticks) ;
392 digit->SetTime(time) ;
394 //Find next signal module
396 for(i = 0 ; i < fInput ; i++){
397 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
398 Int_t curNext = nextSig ;
399 if(sdigits->GetEntriesFast() > index[i] ){
400 curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(index[i]))->GetId() ;
402 if(curNext < nextSig) nextSig = curNext ;
410 //Now CPV digits (different noise and no timing)
411 for(absID = nEMC+1; absID <= nCPV; absID++){
412 Float_t noise = gRandom->Gaus(0., fCPVNoise) ;
413 new((*digits)[absID-1]) AliPHOSDigit( -1,absID,sDigitizer->Digitize(noise), TimeOfNoise() ) ;
414 //look if we have to add signal?
416 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
417 //Add SDigits from all inputs
418 for(i = 0 ; i < fInput ; i++){
419 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
420 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
424 //May be several digits will contribute from the same input
425 while(curSDigit && curSDigit->GetId() == absID){
426 //Shift primary to separate primaries belonging different inputs
427 Int_t primaryoffset ;
429 primaryoffset = fManager->GetMask(i) ;
431 primaryoffset = 10000000*i ;
432 curSDigit->ShiftPrimary(primaryoffset) ;
435 *digit += *curSDigit ;
437 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
438 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i]) ) ;
444 //Find next signal module
446 for(i = 0 ; i < fInput ; i++){
447 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
448 Int_t curNext = nextSig ;
449 if(sdigits->GetEntriesFast() > index[i] )
450 curNext = dynamic_cast<AliPHOSDigit *>( sdigits->At(index[i]) )->GetId() ;
451 if(curNext < nextSig) nextSig = curNext ;
457 delete sdigArray ; //We should not delete its contents
459 //remove digits below thresholds
460 for(i = 0 ; i < nEMC ; i++){
461 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
462 DecalibrateEMC(digit);
463 if(digit->GetEnergy() < fEMCDigitThreshold)
464 digits->RemoveAt(i) ;
466 digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ;
470 for(i = nEMC; i < nCPV ; i++)
471 // if( sDigitizer->Calibrate( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetAmp() ) < fCPVDigitThreshold )
472 if( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetEnergy() < fCPVDigitThreshold )
473 digits->RemoveAt(i) ;
477 Int_t ndigits = digits->GetEntriesFast() ;
478 digits->Expand(ndigits) ;
480 //Set indexes in list of digits and make true digitization of the energy
481 for (i = 0 ; i < ndigits ; i++) {
482 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
483 digit->SetIndexInList(i) ;
484 if(digit->GetId() > fEmcCrystals){ //digitize CPV only
485 digit->SetAmp(DigitizeCPV(digit->GetEnergy(),digit->GetId()) ) ;
491 //set amplitudes in bad channels to zero
492 for(i = 0 ; i <digits->GetEntries(); i++){
493 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
494 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
495 if(relId[1] == 0) // Emc
496 if(gime->CalibData()->IsBadChannelEmc(relId[0],relId[3],relId[2])) digit->SetEnergy(0.);
501 //____________________________________________________________________________
502 void AliPHOSDigitizer::DecalibrateEMC(AliPHOSDigit *digit)
504 // Decalibrate EMC digit, i.e. change its energy by a factor read from CDB
506 AliPHOSGetter* gime = AliPHOSGetter::Instance();
508 if(!gime->CalibData()) {
509 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1);
510 gime->SetCalibData(cdb);
513 //Determine rel.position of the cell absolute ID
515 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
516 Int_t module=relId[0];
518 Int_t column=relId[3];
519 Float_t decalibration = gime->CalibData()->GetADCchannelEmc(module,column,row);
520 Float_t energy = digit->GetEnergy() / decalibration;
521 digit->SetEnergy(energy);
523 //____________________________________________________________________________
524 Int_t AliPHOSDigitizer::DigitizeCPV(Float_t charge, Int_t absId)
526 // Returns digitized value of the CPV charge in a pad absId
528 AliPHOSGetter* gime = AliPHOSGetter::Instance();
530 if(!gime->CalibData()) {
531 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1); // use AliCDBManager's run number
532 gime->SetCalibData(cdb);
535 //Determine rel.position of the cell absId
537 gime->PHOSGeometry()->AbsToRelNumbering(absId,relId);
538 Int_t module=relId[0];
540 Int_t column=relId[3];
544 if(absId > fEmcCrystals){ //digitize CPV only
546 //reading calibration data for cell absId.
547 //If no calibration DB found, accept default values.
549 if(gime->CalibData()) {
550 fADCpedestalCpv = gime->CalibData()->GetADCpedestalCpv(module,column,row);
551 fADCchanelCpv = gime->CalibData()->GetADCchannelCpv( module,column,row);
554 channel = (Int_t) TMath::Ceil((charge - fADCpedestalCpv)/fADCchanelCpv) ;
555 if(channel > fNADCcpv ) channel = fNADCcpv ;
560 //____________________________________________________________________________
561 void AliPHOSDigitizer::Exec(Option_t *option)
563 // Steering method to process digitization for events
564 // in the range from fFirstEvent to fLastEvent.
565 // This range is optionally set by SetEventRange().
566 // if fLastEvent=-1, then process events until the end.
567 // by default fLastEvent = fFirstEvent (process only one event)
569 if (!fInit) { // to prevent overwrite existing file
570 AliError(Form("Give a version name different from %s",
571 fEventFolderName.Data() )) ;
575 if (strstr(option,"print")) {
580 if(strstr(option,"tim"))
581 gBenchmark->Start("PHOSDigitizer");
583 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
585 // Post Digitizer to the white board
586 gime->PostDigitizer(this) ;
588 if (fLastEvent == -1)
589 fLastEvent = gime->MaxEvent() - 1 ;
591 fLastEvent = fFirstEvent ;
593 Int_t nEvents = fLastEvent - fFirstEvent + 1;
597 for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
599 gime->Event(ievent,"S") ;
601 Digitize(ievent) ; //Add prepared SDigits to digits and add the noise
605 if(strstr(option,"deb"))
608 //increment the total number of Digits per run
609 fDigitsInRun += gime->Digits()->GetEntriesFast() ;
612 gime->PhosLoader()->CleanDigitizer();
614 if(strstr(option,"tim")){
615 gBenchmark->Stop("PHOSDigitizer");
617 message = " took %f seconds for Digitizing %f seconds per event\n" ;
618 AliInfo(Form( message.Data(),
619 gBenchmark->GetCpuTime("PHOSDigitizer"),
620 gBenchmark->GetCpuTime("PHOSDigitizer")/nEvents ));
624 //____________________________________________________________________________
625 Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks) const
627 // Returns the shortest time among all time ticks
629 ticks->Sort() ; //Sort in accordance with times of ticks
631 AliPHOSTick * ctick = (AliPHOSTick *) it.Next() ;
632 Float_t time = ctick->CrossingTime(fTimeThreshold) ;
635 while((t=(AliPHOSTick*) it.Next())){
636 if(t->GetTime() < time) //This tick starts before crossing
641 time = ctick->CrossingTime(fTimeThreshold) ;
646 //____________________________________________________________________________
647 Bool_t AliPHOSDigitizer::Init()
649 // Makes all memory allocations
651 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
653 AliFatal(Form("Could not obtain the Getter object for file %s and event %s !",
654 GetTitle(), fEventFolderName.Data()));
658 const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
660 fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ;
662 TString opt("Digits") ;
663 if(gime->VersionExists(opt) ) {
664 AliError(Form("Give a version name different from %s",
665 fEventFolderName.Data() )) ;
670 fLastEvent = fFirstEvent ;
672 fInput = fManager->GetNinputs() ;
676 fInputFileNames = new TString[fInput] ;
677 fEventNames = new TString[fInput] ;
678 fInputFileNames[0] = GetTitle() ;
679 fEventNames[0] = fEventFolderName.Data() ;
681 for (index = 1 ; index < fInput ; index++) {
682 fInputFileNames[index] = dynamic_cast<AliStream*>(fManager->GetInputStream(index))->GetFileName(0);
683 TString tempo = fManager->GetInputFolderName(index) ;
684 fEventNames[index] = tempo.Remove(tempo.Length()-1) ; // strip of the stream number added by fManager
687 //to prevent cleaning of this object while GetEvent is called
688 gime->PhosLoader()->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
693 //____________________________________________________________________________
694 void AliPHOSDigitizer::InitParameters()
696 // Set initial parameters Digitizer
698 fPinNoise = 0.004 ; // [GeV]
699 fEMCDigitThreshold = 0.012 ; // [GeV]
700 fCPVNoise = 0.01; // [aux units]
701 fCPVDigitThreshold = 0.09 ; // [aux units]
702 fTimeResolution = 0.5e-9 ; // [sec]
703 fTimeSignalLength = 1.0e-9 ; // [sec]
705 fADCchanelEmc = 1.0; // Coefficient between real and measured energies in EMC
706 fADCpedestalEmc = 0. ; //
707 fNADCemc = (Int_t) TMath::Power(2,16) ; // number of channels in EMC ADC
709 fADCchanelCpv = 0.0012 ; // width of one ADC channel in CPV 'popugais'
710 fADCpedestalCpv = 0.012 ; //
711 fNADCcpv = (Int_t) TMath::Power(2,12); // number of channels in CPV ADC
713 // fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
714 fTimeThreshold = 0.001 ; // [GeV]
715 SetEventRange(0,-1) ;
719 //__________________________________________________________________
720 void AliPHOSDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
722 // Allows to produce digits by superimposing background and signal event.
723 // It is assumed, that headers file with SIGNAL events is opened in
725 // Sets the BACKGROUND event, with which the SIGNAL event is to be mixed
726 // Thus we avoid writing (changing) huge and expensive
727 // backgound files: all output will be writen into SIGNAL, i.e.
728 // opened in constructor file.
730 // One can open as many files to mix with as one needs.
731 // However only Sdigits with the same name (i.e. constructed with the same SDigitizer)
734 if( strcmp(fEventFolderName, "") == 0 )
738 Warning("MixWith", "Cannot use this method with AliRunDigitizer\n" ) ;
741 // looking for file which contains AliRun
742 if (gSystem->AccessPathName(alirunFileName)) {// file does not exist
743 AliError(Form("File %s does not exist!", alirunFileName.Data())) ;
746 // looking for the file which contains SDigits
747 AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
748 TString fileName( gime->GetSDigitsFileName() ) ;
749 if ( eventFolderName != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
750 fileName = fileName.ReplaceAll(".root", "") + "_" + eventFolderName + ".root" ;
751 if ( (gSystem->AccessPathName(fileName)) ) {
752 AliError(Form("The file %s does not exist!", fileName.Data())) ;
755 // need to increase the arrays
756 TString tempo = fInputFileNames[fInput-1] ;
757 delete [] fInputFileNames ;
758 fInputFileNames = new TString[fInput+1] ;
759 fInputFileNames[fInput-1] = tempo ;
761 tempo = fEventNames[fInput-1] ;
762 delete [] fEventNames ;
763 fEventNames = new TString[fInput+1] ;
764 fEventNames[fInput-1] = tempo ;
766 fInputFileNames[fInput] = alirunFileName ;
767 fEventNames[fInput] = eventFolderName ;
771 //__________________________________________________________________
772 void AliPHOSDigitizer::Print(const Option_t *)const
774 // Print Digitizer's parameters
775 AliInfo(Form("\n------------------- %s -------------", GetName() )) ;
776 if( strcmp(fEventFolderName.Data(), "") != 0 ){
777 printf(" Writing Digits to branch with title %s\n", fEventFolderName.Data()) ;
781 nStreams = GetNInputStreams() ;
786 for (index = 0 ; index < nStreams ; index++) {
787 TString tempo(fEventNames[index]) ;
789 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[index], tempo) ;
790 TString fileName( gime->GetSDigitsFileName() ) ;
791 if ( fEventNames[index] != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
792 fileName = fileName.ReplaceAll(".root", "") + "_" + fEventNames[index] + ".root" ;
793 printf ("Adding SDigits from %s %s\n", fInputFileNames[index].Data(), fileName.Data()) ;
795 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
796 printf("\nWriting digits to %s", gime->GetDigitsFileName().Data()) ;
798 printf("\nWith following parameters:\n") ;
799 printf(" Electronics noise in EMC (fPinNoise) = %f GeV\n", fPinNoise ) ;
800 printf(" Threshold in EMC (fEMCDigitThreshold) = %f GeV\n", fEMCDigitThreshold ) ;
801 printf(" Noise in CPV (fCPVNoise) = %f aux units\n", fCPVNoise ) ;
802 printf(" Threshold in CPV (fCPVDigitThreshold) = %f aux units\n",fCPVDigitThreshold ) ;
803 printf(" ---------------------------------------------------\n") ;
806 AliInfo(Form("AliPHOSDigitizer not initialized" )) ;
810 //__________________________________________________________________
811 void AliPHOSDigitizer::PrintDigits(Option_t * option)
813 // Print a table of digits
815 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
816 TClonesArray * digits = gime->Digits() ;
818 AliInfo(Form("%d", digits->GetEntriesFast())) ;
819 printf("\nevent %d", gAlice->GetEvNumber()) ;
820 printf("\n Number of entries in Digits list %d", digits->GetEntriesFast() ) ;
823 if(strstr(option,"all")||strstr(option,"EMC")){
825 AliPHOSDigit * digit;
826 printf("\nEMC digits (with primaries):\n") ;
827 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
828 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
830 for (index = 0 ; (index < digits->GetEntriesFast()) &&
831 (dynamic_cast<AliPHOSDigit *>(digits->At(index))->GetId() <= maxEmc) ; index++) {
832 digit = (AliPHOSDigit * ) digits->At(index) ;
833 if(digit->GetNprimary() == 0)
835 // printf("%6d %8d %6.5e %4d %2d :",
836 // digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ; // YVK
837 printf("%6d %.4f %6.5e %4d %2d :",
838 digit->GetId(), digit->GetEnergy(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
840 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
841 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
847 if(strstr(option,"all")||strstr(option,"CPV")){
849 //loop over CPV digits
850 AliPHOSDigit * digit;
851 printf("\nCPV digits (with primaries):\n") ;
852 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
853 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
855 for (index = 0 ; index < digits->GetEntriesFast(); index++) {
856 digit = (AliPHOSDigit * ) digits->At(index) ;
857 if(digit->GetId() > maxEmc){
858 printf("%6d %8d %4d %2d :",
859 digit->GetId(), digit->GetAmp(), digit->GetIndexInList(), digit->GetNprimary()) ;
861 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
862 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
871 //__________________________________________________________________
872 Float_t AliPHOSDigitizer::TimeOfNoise(void) const
873 { // Calculates the time signal generated by noise
874 //PH Info("TimeOfNoise", "Change me") ;
875 return gRandom->Rndm() * 1.28E-5;
878 //__________________________________________________________________
879 void AliPHOSDigitizer::Unload()
883 for(i = 1 ; i < fInput ; i++){
884 TString tempo(fEventNames[i]) ;
886 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
887 gime->PhosLoader()->UnloadSDigits() ;
890 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
891 gime->PhosLoader()->UnloadDigits() ;
894 //____________________________________________________________________________
895 void AliPHOSDigitizer::WriteDigits()
898 // Makes TreeD in the output file.
899 // Check if branch already exists:
900 // if yes, exit without writing: ROOT TTree does not support overwriting/updating of
901 // already existing branches.
902 // else creates branch with Digits, named "PHOS", title "...",
903 // and branch "AliPHOSDigitizer", with the same title to keep all the parameters
904 // and names of files, from which digits are made.
906 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
907 const TClonesArray * digits = gime->Digits() ;
908 TTree * treeD = gime->TreeD();
910 // -- create Digits branch
911 Int_t bufferSize = 32000 ;
912 TBranch * digitsBranch = treeD->Branch("PHOS","TClonesArray",&digits,bufferSize);
913 digitsBranch->SetTitle(fEventFolderName);
914 digitsBranch->Fill() ;
916 gime->WriteDigits("OVERWRITE");
917 gime->WriteDigitizer("OVERWRITE");