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:
20 * $Log: AliPHOSDigitizer.cxx,v $
21 * Revision 1.104 2007/12/18 09:08:18 hristov
22 * Splitting of the QA maker into simulation and reconstruction dependent parts (Yves)
24 * Revision 1.103 2007/11/07 11:25:06 schutz
25 * Comment out the QA checking before starting digitization
27 * Revision 1.102 2007/10/19 18:04:29 schutz
28 * The standalone QA data maker is called from AliSimulation and AliReconstruction outside the event loop; i.e. re-reading the data. The QA data making in the event loop has been commented out.
30 * Revision 1.101 2007/10/14 21:08:10 schutz
31 * Introduced the checking of QA results from previous step before entering the event loop
33 * Revision 1.100 2007/10/10 09:05:10 schutz
34 * Changing name QualAss to QA
36 * Revision 1.99 2007/09/30 17:08:20 schutz
37 * Introducing the notion of QA data acquisition cycle (needed by online)
39 * Revision 1.98 2007/09/26 14:22:17 cvetan
40 * 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.
42 * Revision 1.97 2007/08/07 14:12:03 kharlov
43 * Quality assurance added (Yves Schutz)
45 * Revision 1.96 2007/04/28 10:43:36 policheh
46 * Dead channels simulation: digit energy sets to 0.
48 * Revision 1.95 2007/04/10 07:20:52 kharlov
49 * Decalibration should use the same CDB as calibration in AliPHOSClusterizerv1
51 * Revision 1.94 2007/02/01 10:34:47 hristov
52 * Removing warnings on Solaris x86
54 * Revision 1.93 2006/10/17 13:17:01 kharlov
55 * Replace AliInfo by AliDebug
57 * Revision 1.92 2006/08/28 10:01:56 kharlov
58 * Effective C++ warnings fixed (Timur Pocheptsov)
60 * Revision 1.91 2006/04/29 20:25:30 hristov
61 * Decalibration is implemented (Yu.Kharlov)
63 * Revision 1.90 2006/04/22 10:30:17 hristov
64 * Add fEnergy to AliPHOSDigit and operate with EMC amplitude in energy units (Yu.Kharlov)
66 * Revision 1.89 2006/04/11 15:22:59 hristov
67 * run number in query set to -1: forces AliCDBManager to use its run number (A.Colla)
69 * Revision 1.88 2006/03/13 14:05:43 kharlov
70 * Calibration objects for EMC and CPV
72 * Revision 1.87 2005/08/24 15:33:49 kharlov
73 * Calibration data for raw digits
75 * Revision 1.86 2005/07/12 20:07:35 hristov
76 * Changes needed to run simulation and reconstrruction in the same AliRoot session
78 * Revision 1.85 2005/05/28 14:19:04 schutz
79 * Compilation warnings fixed by T.P.
83 //_________________________________________________________________________
84 //*-- Author : Dmitri Peressounko (SUBATECH & Kurchatov Institute)
85 //////////////////////////////////////////////////////////////////////////////
86 // This TTask performs digitization of Summable digits (in the PHOS case it is just
87 // the sum of contributions from all primary particles into a given cell).
88 // In addition it performs mixing of summable digits from different events.
89 // The name of the TTask is also the title of the branch that will contain
90 // the created SDigits
91 // The title of the TTAsk is the name of the file that contains the hits from
92 // which the SDigits are created
94 // For each event two branches are created in TreeD:
95 // "PHOS" - list of digits
96 // "AliPHOSDigitizer" - AliPHOSDigitizer with all parameters used in digitization
98 // Note, that one can set a title for new digits branch, and repeat digitization with
99 // another set of parameters.
102 // root[0] AliPHOSDigitizer * d = new AliPHOSDigitizer() ;
103 // root[1] d->ExecuteTask()
104 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
105 // //Digitizes SDigitis in all events found in file galice.root
107 // root[2] AliPHOSDigitizer * d1 = new AliPHOSDigitizer("galice1.root") ;
108 // // Will read sdigits from galice1.root
109 // root[3] d1->MixWith("galice2.root")
110 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
111 // // Reads another set of sdigits from galice2.root
112 // root[3] d1->MixWith("galice3.root")
113 // // Reads another set of sdigits from galice3.root
114 // root[4] d->ExecuteTask("deb timing")
115 // // Reads SDigits from files galice1.root, galice2.root ....
116 // // mixes them and stores produced Digits in file galice1.root
117 // // deb - prints number of produced digits
118 // // deb all - prints list of produced digits
119 // // timing - prints time used for digitization
122 // --- ROOT system ---
125 #include "TBenchmark.h"
128 // --- Standard library ---
130 // --- AliRoot header files ---
132 #include "AliRunDigitizer.h"
133 #include "AliPHOSDigit.h"
134 #include "AliPHOSGetter.h"
135 #include "AliPHOSDigitizer.h"
136 #include "AliPHOSSDigitizer.h"
137 #include "AliPHOSGeometry.h"
138 #include "AliPHOSTick.h"
140 ClassImp(AliPHOSDigitizer)
143 //____________________________________________________________________________
144 AliPHOSDigitizer::AliPHOSDigitizer() :
150 fInputFileNames(0x0),
154 fEMCDigitThreshold(0.f),
156 fCPVDigitThreshold(0.f),
157 fTimeResolution(0.f),
159 fTimeSignalLength(0.f),
161 fADCpedestalEmc(0.f),
164 fADCpedestalCpv(0.f),
166 fEventFolderName(""),
174 fManager = 0 ; // We work in the standalong mode
177 //____________________________________________________________________________
178 AliPHOSDigitizer::AliPHOSDigitizer(TString alirunFileName,
179 TString eventFolderName):
180 AliDigitizer("PHOS"+AliConfig::Instance()->GetDigitizerTaskName(), alirunFileName),
181 fDefaultInit(kFALSE),
185 fInputFileNames(0x0),
189 fEMCDigitThreshold(0.f),
191 fCPVDigitThreshold(0.f),
192 fTimeResolution(0.f),
194 fTimeSignalLength(0.f),
196 fADCpedestalEmc(0.f),
199 fADCpedestalCpv(0.f),
201 fEventFolderName(eventFolderName),
210 fDefaultInit = kFALSE ;
211 fManager = 0 ; // We work in the standalong mode
212 // //Initialize the quality assurance data maker only once
213 // fQADM = new AliPHOSQADataMaker() ;
214 // //FIXME: get the run number
217 // GetQADataMaker()->Init(AliQA::kDIGITS, run, fgkCycles) ;
218 // GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
221 //____________________________________________________________________________
222 AliPHOSDigitizer::AliPHOSDigitizer(const AliPHOSDigitizer & d) :
224 fDefaultInit(d.fDefaultInit),
225 fDigitsInRun(d.fDigitsInRun),
228 fInputFileNames(0x0),//?
230 fEmcCrystals(d.fEmcCrystals),
231 fPinNoise(d.fPinNoise),
232 fEMCDigitThreshold(d.fEMCDigitThreshold),
233 fCPVNoise(d.fCPVNoise),
234 fCPVDigitThreshold(d.fCPVDigitThreshold),
235 fTimeResolution(d.fTimeResolution),
236 fTimeThreshold(d.fTimeThreshold),
237 fTimeSignalLength(d.fTimeSignalLength),
238 fADCchanelEmc(d.fADCchanelEmc),
239 fADCpedestalEmc(d.fADCpedestalEmc),
240 fNADCemc(d.fNADCemc),
241 fADCchanelCpv(d.fADCchanelCpv),
242 fADCpedestalCpv(d.fADCpedestalCpv),
243 fNADCcpv(d.fNADCcpv),
244 fEventFolderName(d.fEventFolderName),
245 fFirstEvent(d.fFirstEvent),
246 fLastEvent(d.fLastEvent),
252 SetName(d.GetName()) ;
253 SetTitle(d.GetTitle()) ;
254 //Initialize the quality assurance data maker only once
255 //FIXME: get the run number
258 // GetQADataMaker()->Init(AliQA::kDIGITS, run, fgkCycles) ;
259 // GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
262 //____________________________________________________________________________
263 AliPHOSDigitizer::AliPHOSDigitizer(AliRunDigitizer * rd) :
264 AliDigitizer(rd,"PHOS"+AliConfig::Instance()->GetDigitizerTaskName()),
265 fDefaultInit(kFALSE),
269 fInputFileNames(0x0),
273 fEMCDigitThreshold(0.f),
275 fCPVDigitThreshold(0.f),
276 fTimeResolution(0.f),
278 fTimeSignalLength(0.f),
280 fADCpedestalEmc(0.f),
283 fADCpedestalCpv(0.f),
285 fEventFolderName(fManager->GetInputFolderName(0)),
292 // ctor Init() is called by RunDigitizer
294 SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
296 fDefaultInit = kFALSE ;
297 //Initialize the quality assurance data maker only once
298 // fQADM = new AliPHOSQADataMaker() ;
299 // //FIXME: get the run number
302 // GetQADataMaker()->Init(AliQA::kDIGITS, run) ;
303 // GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
306 //____________________________________________________________________________
307 AliPHOSDigitizer::~AliPHOSDigitizer()
309 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
311 // Clean Digitizer from the white board
312 gime->PhosLoader()->CleanDigitizer() ;
314 delete [] fInputFileNames ;
315 delete [] fEventNames ;
321 //____________________________________________________________________________
322 void AliPHOSDigitizer::Digitize(Int_t event)
325 // Makes the digitization of the collected summable digits.
326 // It first creates the array of all PHOS modules
327 // filled with noise (different for EMC, and CPV) and
328 // then adds contributions from SDigits.
329 // This design avoids scanning over the list of digits to add
330 // contribution to new SDigits only.
332 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
333 Int_t ReadEvent = event ;
335 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetCurrentEventNumber() ;
336 AliDebug(1,Form("Adding event %d from input stream 0 %s %s",
337 ReadEvent, GetTitle(), fEventFolderName.Data())) ;
338 gime->Event(ReadEvent, "S") ;
339 TClonesArray * digits = gime->Digits() ;
342 const AliPHOSGeometry *geom = gime->PHOSGeometry() ;
343 //Making digits with noise, first EMC
344 Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
349 nCPV = nEMC + geom->GetNumberOfCPVPadsZ() * geom->GetNumberOfCPVPadsPhi() * geom->GetNModules() ;
351 digits->Expand(nCPV) ;
353 // get first the sdigitizer from the tasks list
354 if ( !gime->SDigitizer() )
355 gime->LoadSDigitizer();
356 AliPHOSSDigitizer * sDigitizer = gime->SDigitizer();
359 AliFatal(Form("SDigitizer with name %s %s not found",
360 GetTitle(), fEventFolderName.Data() )) ;
362 //take all the inputs to add together and load the SDigits
363 TObjArray * sdigArray = new TObjArray(fInput) ;
364 sdigArray->AddAt(gime->SDigits(), 0) ;
366 for(i = 1 ; i < fInput ; i++){
367 TString tempo(fEventNames[i]) ;
369 AliPHOSGetter * gime1 = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
371 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(i))->GetCurrentEventNumber() ;
372 AliInfo(Form("Adding event %d from input stream %d %s %s",
373 ReadEvent, i, fInputFileNames[i].Data(), tempo.Data())) ;
374 gime1->Event(ReadEvent,"S");
375 sdigArray->AddAt(gime1->SDigits(), i) ;
378 //Find the first crystal with signal
379 Int_t nextSig = 200000 ;
380 TClonesArray * sdigits ;
381 for(i = 0 ; i < fInput ; i++){
382 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
383 if ( !sdigits->GetEntriesFast() )
385 Int_t curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(0))->GetId() ;
386 if(curNext < nextSig)
390 TArrayI index(fInput) ;
391 index.Reset() ; //Set all indexes to zero
393 AliPHOSDigit * digit ;
394 AliPHOSDigit * curSDigit ;
396 TClonesArray * ticks = new TClonesArray("AliPHOSTick",1000) ;
398 //Put Noise contribution
399 for(absID = 1 ; absID <= nEMC ; absID++){
400 Float_t noise = gRandom->Gaus(0., fPinNoise) ;
401 // YVK: do not digitize amplitudes for EMC
402 // new((*digits)[absID-1]) AliPHOSDigit( -1, absID, sDigitizer->Digitize(noise), TimeOfNoise() ) ;
403 new((*digits)[absID-1]) AliPHOSDigit( -1, absID, noise, TimeOfNoise() ) ;
404 //look if we have to add signal?
405 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
408 //Add SDigits from all inputs
411 Float_t a = digit->GetEnergy() ;
412 Float_t b = TMath::Abs( a / fTimeSignalLength) ;
413 //Mark the beginning of the signal
414 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime(),0, b);
415 //Mark the end of the signal
416 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime()+fTimeSignalLength, -a, -b);
419 for(i = 0 ; i < fInput ; i++){
420 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
421 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) ;
434 a = curSDigit->GetEnergy() ;
435 b = a /fTimeSignalLength ;
436 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime(),0, b);
437 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b);
439 *digit += *curSDigit ; //add energies
442 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
443 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
449 //calculate and set time
450 Float_t time = FrontEdgeTime(ticks) ;
451 digit->SetTime(time) ;
453 //Find next signal module
455 for(i = 0 ; i < fInput ; i++){
456 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
457 Int_t curNext = nextSig ;
458 if(sdigits->GetEntriesFast() > index[i] ){
459 curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(index[i]))->GetId() ;
461 if(curNext < nextSig) nextSig = curNext ;
469 //Now CPV digits (different noise and no timing)
470 for(absID = nEMC+1; absID <= nCPV; absID++){
471 Float_t noise = gRandom->Gaus(0., fCPVNoise) ;
472 new((*digits)[absID-1]) AliPHOSDigit( -1,absID,sDigitizer->Digitize(noise), TimeOfNoise() ) ;
473 //look if we have to add signal?
475 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
476 //Add SDigits from all inputs
477 for(i = 0 ; i < fInput ; i++){
478 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
479 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
483 //May be several digits will contribute from the same input
484 while(curSDigit && curSDigit->GetId() == absID){
485 //Shift primary to separate primaries belonging different inputs
486 Int_t primaryoffset ;
488 primaryoffset = fManager->GetMask(i) ;
490 primaryoffset = 10000000*i ;
491 curSDigit->ShiftPrimary(primaryoffset) ;
494 *digit += *curSDigit ;
496 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
497 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i]) ) ;
503 //Find next signal module
505 for(i = 0 ; i < fInput ; i++){
506 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
507 Int_t curNext = nextSig ;
508 if(sdigits->GetEntriesFast() > index[i] )
509 curNext = dynamic_cast<AliPHOSDigit *>( sdigits->At(index[i]) )->GetId() ;
510 if(curNext < nextSig) nextSig = curNext ;
516 delete sdigArray ; //We should not delete its contents
518 //remove digits below thresholds
519 for(i = 0 ; i < nEMC ; i++){
520 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
521 //By default no decalibration should be applyed
522 // DecalibrateEMC(digit);
524 if(digit->GetEnergy() < fEMCDigitThreshold)
525 digits->RemoveAt(i) ;
527 digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ;
531 for(i = nEMC; i < nCPV ; i++)
532 // if( sDigitizer->Calibrate( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetAmp() ) < fCPVDigitThreshold )
533 if( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetEnergy() < fCPVDigitThreshold )
534 digits->RemoveAt(i) ;
538 Int_t ndigits = digits->GetEntriesFast() ;
539 digits->Expand(ndigits) ;
541 //Set indexes in list of digits and make true digitization of the energy
542 for (i = 0 ; i < ndigits ; i++) {
543 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
544 digit->SetIndexInList(i) ;
545 if(digit->GetId() > fEmcCrystals){ //digitize CPV only
546 digit->SetAmp(DigitizeCPV(digit->GetEnergy(),digit->GetId()) ) ;
552 //set amplitudes in bad channels to zero
553 if(!gime->CalibData()) {
554 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1);
555 gime->SetCalibData(cdb);
557 for(i = 0 ; i <digits->GetEntries(); i++){
558 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
559 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
560 if(relId[1] == 0) // Emc
561 if(gime->CalibData()->IsBadChannelEmc(relId[0],relId[3],relId[2])) digit->SetEnergy(0.);
566 //____________________________________________________________________________
567 void AliPHOSDigitizer::DecalibrateEMC(AliPHOSDigit *digit)
569 // Decalibrate EMC digit, i.e. change its energy by a factor read from CDB
571 AliPHOSGetter* gime = AliPHOSGetter::Instance();
573 if(!gime->CalibData()) {
574 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1);
575 gime->SetCalibData(cdb);
578 //Determine rel.position of the cell absolute ID
580 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
581 Int_t module=relId[0];
583 Int_t column=relId[3];
584 Float_t decalibration = gime->CalibData()->GetADCchannelEmc(module,column,row);
585 Float_t energy = digit->GetEnergy() / decalibration;
586 digit->SetEnergy(energy);
588 //____________________________________________________________________________
589 Int_t AliPHOSDigitizer::DigitizeCPV(Float_t charge, Int_t absId)
591 // Returns digitized value of the CPV charge in a pad absId
593 AliPHOSGetter* gime = AliPHOSGetter::Instance();
595 if(!gime->CalibData()) {
596 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1); // use AliCDBManager's run number
597 gime->SetCalibData(cdb);
600 //Determine rel.position of the cell absId
602 gime->PHOSGeometry()->AbsToRelNumbering(absId,relId);
603 Int_t module=relId[0];
605 Int_t column=relId[3];
609 if(absId > fEmcCrystals){ //digitize CPV only
611 //reading calibration data for cell absId.
612 //If no calibration DB found, accept default values.
614 if(gime->CalibData()) {
615 fADCpedestalCpv = gime->CalibData()->GetADCpedestalCpv(module,column,row);
616 fADCchanelCpv = gime->CalibData()->GetADCchannelCpv( module,column,row);
619 channel = (Int_t) TMath::Ceil((charge - fADCpedestalCpv)/fADCchanelCpv) ;
620 if(channel > fNADCcpv ) channel = fNADCcpv ;
625 //____________________________________________________________________________
626 void AliPHOSDigitizer::Exec(Option_t *option)
628 // Steering method to process digitization for events
629 // in the range from fFirstEvent to fLastEvent.
630 // This range is optionally set by SetEventRange().
631 // if fLastEvent=-1, then process events until the end.
632 // by default fLastEvent = fFirstEvent (process only one event)
634 if (!fInit) { // to prevent overwrite existing file
635 AliError(Form("Give a version name different from %s",
636 fEventFolderName.Data() )) ;
640 if (strstr(option,"print")) {
645 // // check the QA result for Hits and SDigits
646 // AliQA * qa = AliQA::Instance(AliQA::kPHOS) ;
647 // if ( qa->IsSet(AliQA::kPHOS, AliQA::kSIM, AliQA::kFATAL)) {
648 // AliFatal("QA status in Hits and/or SDIGITS was Fatal") ;
649 // } else if ( qa->IsSet(AliQA::kPHOS, AliQA::kSIM, AliQA::kERROR)) {
650 // AliError("QA status in Hits and/or SDIGITS was Error") ;
651 // } else if ( qa->IsSet(AliQA::kPHOS, AliQA::kSIM, AliQA::kWARNING) ) {
652 // AliWarning("QA status in Hits and/or SDIGITS was Warning") ;
653 // } else if ( qa->IsSet(AliQA::kPHOS, AliQA::kSIM, AliQA::kINFO) ) {
654 // AliInfo("QA status in Hits and/or SDIGITS was Info") ;
657 if(strstr(option,"tim"))
658 gBenchmark->Start("PHOSDigitizer");
660 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
662 // Post Digitizer to the white board
663 gime->PostDigitizer(this) ;
665 if (fLastEvent == -1)
666 fLastEvent = gime->MaxEvent() - 1 ;
668 fLastEvent = fFirstEvent ;
670 Int_t nEvents = fLastEvent - fFirstEvent + 1;
674 for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
676 gime->Event(ievent,"S") ;
678 Digitize(ievent) ; //Add prepared SDigits to digits and add the noise
680 // //makes the quality assurance data
681 // if (GetQADataMaker()->IsCycleDone() ) {
682 // GetQADataMaker()->EndOfCycle(AliQA::kDIGITS) ;
683 // GetQADataMaker()->StartOfCycle(AliQA::kDIGITS) ;
685 // GetQADataMaker()->Exec(AliQA::kDIGITS, gime->Digits()) ;
686 // GetQADataMaker()->Increment() ;
690 if(strstr(option,"deb"))
693 //increment the total number of Digits per run
694 fDigitsInRun += gime->Digits()->GetEntriesFast() ;
697 // //Write the quality assurance data only after the last event
698 // if ( fEventCounter == gime->MaxEvent() ) {
699 // GetQADataMaker()->EndOfCycle(AliQA::kDIGITS) ;
700 // GetQADataMaker()->Finish() ;
703 gime->PhosLoader()->CleanDigitizer();
705 if(strstr(option,"tim")){
706 gBenchmark->Stop("PHOSDigitizer");
708 message = " took %f seconds for Digitizing %f seconds per event\n" ;
709 AliInfo(Form( message.Data(),
710 gBenchmark->GetCpuTime("PHOSDigitizer"),
711 gBenchmark->GetCpuTime("PHOSDigitizer")/nEvents ));
715 //____________________________________________________________________________
716 Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks) const
718 // Returns the shortest time among all time ticks
720 ticks->Sort() ; //Sort in accordance with times of ticks
722 AliPHOSTick * ctick = (AliPHOSTick *) it.Next() ;
723 Float_t time = ctick->CrossingTime(fTimeThreshold) ;
726 while((t=(AliPHOSTick*) it.Next())){
727 if(t->GetTime() < time) //This tick starts before crossing
732 time = ctick->CrossingTime(fTimeThreshold) ;
737 //____________________________________________________________________________
738 Bool_t AliPHOSDigitizer::Init()
740 // Makes all memory allocations
742 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
744 AliFatal(Form("Could not obtain the Getter object for file %s and event %s !",
745 GetTitle(), fEventFolderName.Data()));
749 const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
751 fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ;
753 TString opt("Digits") ;
754 if(gime->VersionExists(opt) ) {
755 AliError(Form("Give a version name different from %s",
756 fEventFolderName.Data() )) ;
761 fLastEvent = fFirstEvent ;
763 fInput = fManager->GetNinputs() ;
767 fInputFileNames = new TString[fInput] ;
768 fEventNames = new TString[fInput] ;
769 fInputFileNames[0] = GetTitle() ;
770 fEventNames[0] = fEventFolderName.Data() ;
772 for (index = 1 ; index < fInput ; index++) {
773 fInputFileNames[index] = dynamic_cast<AliStream*>(fManager->GetInputStream(index))->GetFileName(0);
774 TString tempo = fManager->GetInputFolderName(index) ;
775 fEventNames[index] = tempo.Remove(tempo.Length()-1) ; // strip of the stream number added by fManager
778 //to prevent cleaning of this object while GetEvent is called
779 gime->PhosLoader()->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
784 //____________________________________________________________________________
785 void AliPHOSDigitizer::InitParameters()
787 // Set initial parameters Digitizer
789 fPinNoise = 0.004 ; // [GeV]
790 fEMCDigitThreshold = 0.012 ; // [GeV]
791 fCPVNoise = 0.01; // [aux units]
792 fCPVDigitThreshold = 0.09 ; // [aux units]
793 fTimeResolution = 0.5e-9 ; // [sec]
794 fTimeSignalLength = 1.0e-9 ; // [sec]
796 fADCchanelEmc = 1.0; // Coefficient between real and measured energies in EMC
797 fADCpedestalEmc = 0. ; //
798 fNADCemc = (Int_t) TMath::Power(2,16) ; // number of channels in EMC ADC
800 fADCchanelCpv = 0.0012 ; // width of one ADC channel in CPV 'popugais'
801 fADCpedestalCpv = 0.012 ; //
802 fNADCcpv = (Int_t) TMath::Power(2,12); // number of channels in CPV ADC
804 // fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
805 fTimeThreshold = 0.001 ; // [GeV]
806 SetEventRange(0,-1) ;
810 //__________________________________________________________________
811 void AliPHOSDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
813 // Allows to produce digits by superimposing background and signal event.
814 // It is assumed, that headers file with SIGNAL events is opened in
816 // Sets the BACKGROUND event, with which the SIGNAL event is to be mixed
817 // Thus we avoid writing (changing) huge and expensive
818 // backgound files: all output will be writen into SIGNAL, i.e.
819 // opened in constructor file.
821 // One can open as many files to mix with as one needs.
822 // However only Sdigits with the same name (i.e. constructed with the same SDigitizer)
825 if( strcmp(fEventFolderName, "") == 0 )
829 Warning("MixWith", "Cannot use this method with AliRunDigitizer\n" ) ;
832 // looking for file which contains AliRun
833 if (gSystem->AccessPathName(alirunFileName)) {// file does not exist
834 AliError(Form("File %s does not exist!", alirunFileName.Data())) ;
837 // looking for the file which contains SDigits
838 AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
839 TString fileName( gime->GetSDigitsFileName() ) ;
840 if ( eventFolderName != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
841 fileName = fileName.ReplaceAll(".root", "") + "_" + eventFolderName + ".root" ;
842 if ( (gSystem->AccessPathName(fileName)) ) {
843 AliError(Form("The file %s does not exist!", fileName.Data())) ;
846 // need to increase the arrays
847 TString tempo = fInputFileNames[fInput-1] ;
848 delete [] fInputFileNames ;
849 fInputFileNames = new TString[fInput+1] ;
850 fInputFileNames[fInput-1] = tempo ;
852 tempo = fEventNames[fInput-1] ;
853 delete [] fEventNames ;
854 fEventNames = new TString[fInput+1] ;
855 fEventNames[fInput-1] = tempo ;
857 fInputFileNames[fInput] = alirunFileName ;
858 fEventNames[fInput] = eventFolderName ;
862 //__________________________________________________________________
863 void AliPHOSDigitizer::Print(const Option_t *)const
865 // Print Digitizer's parameters
866 AliInfo(Form("\n------------------- %s -------------", GetName() )) ;
867 if( strcmp(fEventFolderName.Data(), "") != 0 ){
868 printf(" Writing Digits to branch with title %s\n", fEventFolderName.Data()) ;
872 nStreams = GetNInputStreams() ;
877 for (index = 0 ; index < nStreams ; index++) {
878 TString tempo(fEventNames[index]) ;
880 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[index], tempo) ;
881 TString fileName( gime->GetSDigitsFileName() ) ;
882 if ( fEventNames[index] != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
883 fileName = fileName.ReplaceAll(".root", "") + "_" + fEventNames[index] + ".root" ;
884 printf ("Adding SDigits from %s %s\n", fInputFileNames[index].Data(), fileName.Data()) ;
886 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
887 printf("\nWriting digits to %s", gime->GetDigitsFileName().Data()) ;
889 printf("\nWith following parameters:\n") ;
890 printf(" Electronics noise in EMC (fPinNoise) = %f GeV\n", fPinNoise ) ;
891 printf(" Threshold in EMC (fEMCDigitThreshold) = %f GeV\n", fEMCDigitThreshold ) ;
892 printf(" Noise in CPV (fCPVNoise) = %f aux units\n", fCPVNoise ) ;
893 printf(" Threshold in CPV (fCPVDigitThreshold) = %f aux units\n",fCPVDigitThreshold ) ;
894 printf(" ---------------------------------------------------\n") ;
897 AliInfo(Form("AliPHOSDigitizer not initialized" )) ;
901 //__________________________________________________________________
902 void AliPHOSDigitizer::PrintDigits(Option_t * option)
904 // Print a table of digits
906 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
907 TClonesArray * digits = gime->Digits() ;
909 AliInfo(Form("%d", digits->GetEntriesFast())) ;
910 printf("\nevent %d", gAlice->GetEvNumber()) ;
911 printf("\n Number of entries in Digits list %d", digits->GetEntriesFast() ) ;
914 if(strstr(option,"all")||strstr(option,"EMC")){
916 AliPHOSDigit * digit;
917 printf("\nEMC digits (with primaries):\n") ;
918 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
919 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
921 for (index = 0 ; (index < digits->GetEntriesFast()) &&
922 (dynamic_cast<AliPHOSDigit *>(digits->At(index))->GetId() <= maxEmc) ; index++) {
923 digit = (AliPHOSDigit * ) digits->At(index) ;
924 if(digit->GetNprimary() == 0)
926 // printf("%6d %8d %6.5e %4d %2d :",
927 // digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ; // YVK
928 printf("%6d %.4f %6.5e %4d %2d :",
929 digit->GetId(), digit->GetEnergy(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
931 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
932 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
938 if(strstr(option,"all")||strstr(option,"CPV")){
940 //loop over CPV digits
941 AliPHOSDigit * digit;
942 printf("\nCPV digits (with primaries):\n") ;
943 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
944 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
946 for (index = 0 ; index < digits->GetEntriesFast(); index++) {
947 digit = (AliPHOSDigit * ) digits->At(index) ;
948 if(digit->GetId() > maxEmc){
949 printf("%6d %8d %4d %2d :",
950 digit->GetId(), digit->GetAmp(), digit->GetIndexInList(), digit->GetNprimary()) ;
952 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
953 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
962 //__________________________________________________________________
963 Float_t AliPHOSDigitizer::TimeOfNoise(void) const
964 { // Calculates the time signal generated by noise
965 //PH Info("TimeOfNoise", "Change me") ;
966 return gRandom->Rndm() * 1.28E-5;
969 //__________________________________________________________________
970 void AliPHOSDigitizer::Unload()
974 for(i = 1 ; i < fInput ; i++){
975 TString tempo(fEventNames[i]) ;
977 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
978 gime->PhosLoader()->UnloadSDigits() ;
981 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
982 gime->PhosLoader()->UnloadDigits() ;
985 //____________________________________________________________________________
986 void AliPHOSDigitizer::WriteDigits()
989 // Makes TreeD in the output file.
990 // Check if branch already exists:
991 // if yes, exit without writing: ROOT TTree does not support overwriting/updating of
992 // already existing branches.
993 // else creates branch with Digits, named "PHOS", title "...",
994 // and branch "AliPHOSDigitizer", with the same title to keep all the parameters
995 // and names of files, from which digits are made.
997 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
998 const TClonesArray * digits = gime->Digits() ;
999 TTree * treeD = gime->TreeD();
1001 // -- create Digits branch
1002 Int_t bufferSize = 32000 ;
1003 TBranch * digitsBranch = treeD->Branch("PHOS","TClonesArray",&digits,bufferSize);
1004 digitsBranch->SetTitle(fEventFolderName);
1005 digitsBranch->Fill() ;
1007 gime->WriteDigits("OVERWRITE");
1008 gime->WriteDigitizer("OVERWRITE");