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.92 2006/08/28 10:01:56 kharlov
22 * Effective C++ warnings fixed (Timur Pocheptsov)
24 * Revision 1.91 2006/04/29 20:25:30 hristov
25 * Decalibration is implemented (Yu.Kharlov)
27 * Revision 1.90 2006/04/22 10:30:17 hristov
28 * Add fEnergy to AliPHOSDigit and operate with EMC amplitude in energy units (Yu.Kharlov)
30 * Revision 1.89 2006/04/11 15:22:59 hristov
31 * run number in query set to -1: forces AliCDBManager to use its run number (A.Colla)
33 * Revision 1.88 2006/03/13 14:05:43 kharlov
34 * Calibration objects for EMC and CPV
36 * Revision 1.87 2005/08/24 15:33:49 kharlov
37 * Calibration data for raw digits
39 * Revision 1.86 2005/07/12 20:07:35 hristov
40 * Changes needed to run simulation and reconstrruction in the same AliRoot session
42 * Revision 1.85 2005/05/28 14:19:04 schutz
43 * Compilation warnings fixed by T.P.
47 //_________________________________________________________________________
48 //*-- Author : Dmitri Peressounko (SUBATECH & Kurchatov Institute)
49 //////////////////////////////////////////////////////////////////////////////
50 // This TTask performs digitization of Summable digits (in the PHOS case it is just
51 // the sum of contributions from all primary particles into a given cell).
52 // In addition it performs mixing of summable digits from different events.
53 // The name of the TTask is also the title of the branch that will contain
54 // the created SDigits
55 // The title of the TTAsk is the name of the file that contains the hits from
56 // which the SDigits are created
58 // For each event two branches are created in TreeD:
59 // "PHOS" - list of digits
60 // "AliPHOSDigitizer" - AliPHOSDigitizer with all parameters used in digitization
62 // Note, that one can set a title for new digits branch, and repeat digitization with
63 // another set of parameters.
66 // root[0] AliPHOSDigitizer * d = new AliPHOSDigitizer() ;
67 // root[1] d->ExecuteTask()
68 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
69 // //Digitizes SDigitis in all events found in file galice.root
71 // root[2] AliPHOSDigitizer * d1 = new AliPHOSDigitizer("galice1.root") ;
72 // // Will read sdigits from galice1.root
73 // root[3] d1->MixWith("galice2.root")
74 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
75 // // Reads another set of sdigits from galice2.root
76 // root[3] d1->MixWith("galice3.root")
77 // // Reads another set of sdigits from galice3.root
78 // root[4] d->ExecuteTask("deb timing")
79 // // Reads SDigits from files galice1.root, galice2.root ....
80 // // mixes them and stores produced Digits in file galice1.root
81 // // deb - prints number of produced digits
82 // // deb all - prints list of produced digits
83 // // timing - prints time used for digitization
86 // --- ROOT system ---
89 #include "TBenchmark.h"
92 // --- Standard library ---
94 // --- AliRoot header files ---
96 #include "AliRunDigitizer.h"
97 #include "AliPHOSDigit.h"
98 #include "AliPHOSGetter.h"
99 #include "AliPHOSDigitizer.h"
100 #include "AliPHOSSDigitizer.h"
101 #include "AliPHOSGeometry.h"
102 #include "AliPHOSTick.h"
104 ClassImp(AliPHOSDigitizer)
107 //____________________________________________________________________________
108 AliPHOSDigitizer::AliPHOSDigitizer() :
114 fInputFileNames(0x0),
118 fEMCDigitThreshold(0.f),
120 fCPVDigitThreshold(0.f),
121 fTimeResolution(0.f),
123 fTimeSignalLength(0.f),
125 fADCpedestalEmc(0.f),
128 fADCpedestalCpv(0.f),
130 fEventFolderName(""),
136 fManager = 0 ; // We work in the standalong mode
139 //____________________________________________________________________________
140 AliPHOSDigitizer::AliPHOSDigitizer(TString alirunFileName,
141 TString eventFolderName):
142 AliDigitizer("PHOS"+AliConfig::Instance()->GetDigitizerTaskName(), alirunFileName),
143 fDefaultInit(kFALSE),
147 fInputFileNames(0x0),
151 fEMCDigitThreshold(0.f),
153 fCPVDigitThreshold(0.f),
154 fTimeResolution(0.f),
156 fTimeSignalLength(0.f),
158 fADCpedestalEmc(0.f),
161 fADCpedestalCpv(0.f),
163 fEventFolderName(eventFolderName),
170 fDefaultInit = kFALSE ;
171 fManager = 0 ; // We work in the standalong mode
174 //____________________________________________________________________________
175 AliPHOSDigitizer::AliPHOSDigitizer(const AliPHOSDigitizer & d) :
177 fDefaultInit(d.fDefaultInit),
178 fDigitsInRun(d.fDigitsInRun),
181 fInputFileNames(0x0),//?
183 fEmcCrystals(d.fEmcCrystals),
184 fPinNoise(d.fPinNoise),
185 fEMCDigitThreshold(d.fEMCDigitThreshold),
186 fCPVNoise(d.fCPVNoise),
187 fCPVDigitThreshold(d.fCPVDigitThreshold),
188 fTimeResolution(d.fTimeResolution),
189 fTimeThreshold(d.fTimeThreshold),
190 fTimeSignalLength(d.fTimeSignalLength),
191 fADCchanelEmc(d.fADCchanelEmc),
192 fADCpedestalEmc(d.fADCpedestalEmc),
193 fNADCemc(d.fNADCemc),
194 fADCchanelCpv(d.fADCchanelCpv),
195 fADCpedestalCpv(d.fADCpedestalCpv),
196 fNADCcpv(d.fNADCcpv),
197 fEventFolderName(d.fEventFolderName),
198 fFirstEvent(d.fFirstEvent),
199 fLastEvent(d.fLastEvent)
202 SetName(d.GetName()) ;
203 SetTitle(d.GetTitle()) ;
206 //____________________________________________________________________________
207 AliPHOSDigitizer::AliPHOSDigitizer(AliRunDigitizer * rd) :
208 AliDigitizer(rd,"PHOS"+AliConfig::Instance()->GetDigitizerTaskName()),
209 fDefaultInit(kFALSE),
213 fInputFileNames(0x0),
217 fEMCDigitThreshold(0.f),
219 fCPVDigitThreshold(0.f),
220 fTimeResolution(0.f),
222 fTimeSignalLength(0.f),
224 fADCpedestalEmc(0.f),
227 fADCpedestalCpv(0.f),
229 fEventFolderName(fManager->GetInputFolderName(0)),
233 // ctor Init() is called by RunDigitizer
235 SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
237 fDefaultInit = kFALSE ;
240 //____________________________________________________________________________
241 AliPHOSDigitizer::~AliPHOSDigitizer()
243 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
245 // Clean Digitizer from the white board
246 gime->PhosLoader()->CleanDigitizer() ;
248 delete [] fInputFileNames ;
249 delete [] fEventNames ;
253 //____________________________________________________________________________
254 void AliPHOSDigitizer::Digitize(Int_t event)
257 // Makes the digitization of the collected summable digits.
258 // It first creates the array of all PHOS modules
259 // filled with noise (different for EMC, and CPV) and
260 // then adds contributions from SDigits.
261 // This design avoids scanning over the list of digits to add
262 // contribution to new SDigits only.
264 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
265 Int_t ReadEvent = event ;
267 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetCurrentEventNumber() ;
268 AliDebug(1,Form("Adding event %d from input stream 0 %s %s",
269 ReadEvent, GetTitle(), fEventFolderName.Data())) ;
270 gime->Event(ReadEvent, "S") ;
271 TClonesArray * digits = gime->Digits() ;
274 const AliPHOSGeometry *geom = gime->PHOSGeometry() ;
275 //Making digits with noise, first EMC
276 Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
281 nCPV = nEMC + geom->GetNumberOfCPVPadsZ() * geom->GetNumberOfCPVPadsPhi() * geom->GetNModules() ;
283 digits->Expand(nCPV) ;
285 // get first the sdigitizer from the tasks list
286 if ( !gime->SDigitizer() )
287 gime->LoadSDigitizer();
288 AliPHOSSDigitizer * sDigitizer = gime->SDigitizer();
291 AliFatal(Form("SDigitizer with name %s %s not found",
292 GetTitle(), fEventFolderName.Data() )) ;
294 //take all the inputs to add together and load the SDigits
295 TObjArray * sdigArray = new TObjArray(fInput) ;
296 sdigArray->AddAt(gime->SDigits(), 0) ;
298 for(i = 1 ; i < fInput ; i++){
299 TString tempo(fEventNames[i]) ;
301 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
303 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(i))->GetCurrentEventNumber() ;
304 AliInfo(Form("Adding event %d from input stream %d %s %s",
305 ReadEvent, i, fInputFileNames[i].Data(), tempo.Data())) ;
306 gime->Event(ReadEvent,"S");
307 sdigArray->AddAt(gime->SDigits(), i) ;
310 //Find the first crystal with signal
311 Int_t nextSig = 200000 ;
312 TClonesArray * sdigits ;
313 for(i = 0 ; i < fInput ; i++){
314 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
315 if ( !sdigits->GetEntriesFast() )
317 Int_t curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(0))->GetId() ;
318 if(curNext < nextSig)
322 TArrayI index(fInput) ;
323 index.Reset() ; //Set all indexes to zero
325 AliPHOSDigit * digit ;
326 AliPHOSDigit * curSDigit ;
328 TClonesArray * ticks = new TClonesArray("AliPHOSTick",1000) ;
330 //Put Noise contribution
331 for(absID = 1 ; absID <= nEMC ; absID++){
332 Float_t noise = gRandom->Gaus(0., fPinNoise) ;
333 // YVK: do not digitize amplitudes for EMC
334 // new((*digits)[absID-1]) AliPHOSDigit( -1, absID, sDigitizer->Digitize(noise), TimeOfNoise() ) ;
335 new((*digits)[absID-1]) AliPHOSDigit( -1, absID, noise, TimeOfNoise() ) ;
336 //look if we have to add signal?
337 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
340 //Add SDigits from all inputs
343 Float_t a = digit->GetEnergy() ;
344 Float_t b = TMath::Abs( a / fTimeSignalLength) ;
345 //Mark the beginning of the signal
346 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime(),0, b);
347 //Mark the end of the signal
348 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime()+fTimeSignalLength, -a, -b);
351 for(i = 0 ; i < fInput ; i++){
352 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
353 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
356 //May be several digits will contribute from the same input
357 while(curSDigit && curSDigit->GetId() == absID){
358 //Shift primary to separate primaries belonging different inputs
359 Int_t primaryoffset ;
361 primaryoffset = fManager->GetMask(i) ;
363 primaryoffset = 10000000*i ;
364 curSDigit->ShiftPrimary(primaryoffset) ;
366 a = curSDigit->GetEnergy() ;
367 b = a /fTimeSignalLength ;
368 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime(),0, b);
369 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b);
371 *digit += *curSDigit ; //add energies
374 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
375 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
381 //calculate and set time
382 Float_t time = FrontEdgeTime(ticks) ;
383 digit->SetTime(time) ;
385 //Find next signal module
387 for(i = 0 ; i < fInput ; i++){
388 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
389 Int_t curNext = nextSig ;
390 if(sdigits->GetEntriesFast() > index[i] ){
391 curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(index[i]))->GetId() ;
393 if(curNext < nextSig) nextSig = curNext ;
401 //Now CPV digits (different noise and no timing)
402 for(absID = nEMC+1; absID <= nCPV; absID++){
403 Float_t noise = gRandom->Gaus(0., fCPVNoise) ;
404 new((*digits)[absID-1]) AliPHOSDigit( -1,absID,sDigitizer->Digitize(noise), TimeOfNoise() ) ;
405 //look if we have to add signal?
407 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
408 //Add SDigits from all inputs
409 for(i = 0 ; i < fInput ; i++){
410 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
411 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
415 //May be several digits will contribute from the same input
416 while(curSDigit && curSDigit->GetId() == absID){
417 //Shift primary to separate primaries belonging different inputs
418 Int_t primaryoffset ;
420 primaryoffset = fManager->GetMask(i) ;
422 primaryoffset = 10000000*i ;
423 curSDigit->ShiftPrimary(primaryoffset) ;
426 *digit += *curSDigit ;
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 //Find next signal module
437 for(i = 0 ; i < fInput ; i++){
438 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
439 Int_t curNext = nextSig ;
440 if(sdigits->GetEntriesFast() > index[i] )
441 curNext = dynamic_cast<AliPHOSDigit *>( sdigits->At(index[i]) )->GetId() ;
442 if(curNext < nextSig) nextSig = curNext ;
448 delete sdigArray ; //We should not delete its contents
450 //remove digits below thresholds
451 for(i = 0 ; i < nEMC ; i++){
452 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
453 DecalibrateEMC(digit);
454 if(digit->GetEnergy() < fEMCDigitThreshold)
455 digits->RemoveAt(i) ;
457 digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ;
461 for(i = nEMC; i < nCPV ; i++)
462 // if( sDigitizer->Calibrate( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetAmp() ) < fCPVDigitThreshold )
463 if( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetEnergy() < fCPVDigitThreshold )
464 digits->RemoveAt(i) ;
468 Int_t ndigits = digits->GetEntriesFast() ;
469 digits->Expand(ndigits) ;
471 //Set indexes in list of digits and make true digitization of the energy
472 for (i = 0 ; i < ndigits ; i++) {
473 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
474 digit->SetIndexInList(i) ;
475 if(digit->GetId() > fEmcCrystals){ //digitize CPV only
476 digit->SetAmp(DigitizeCPV(digit->GetEnergy(),digit->GetId()) ) ;
481 //____________________________________________________________________________
482 void AliPHOSDigitizer::DecalibrateEMC(AliPHOSDigit *digit)
484 // Decalibrate EMC digit, i.e. change its energy by a factor read from CDB
486 AliPHOSGetter* gime = AliPHOSGetter::Instance();
488 if(!gime->CalibData()) {
489 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1);
490 gime->SetCalibData(cdb);
493 //Determine rel.position of the cell absolute ID
495 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
496 Int_t module=relId[0];
498 Int_t column=relId[3];
499 Float_t decalibration = gime->CalibData()->GetADCchannelEmc(module,column,row);
500 Float_t energy = digit->GetEnergy() * decalibration;
501 digit->SetEnergy(energy);
503 //____________________________________________________________________________
504 Int_t AliPHOSDigitizer::DigitizeCPV(Float_t charge, Int_t absId)
506 // Returns digitized value of the CPV charge in a pad absId
508 AliPHOSGetter* gime = AliPHOSGetter::Instance();
510 if(!gime->CalibData()) {
511 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1); // use AliCDBManager's run number
512 gime->SetCalibData(cdb);
515 //Determine rel.position of the cell absId
517 gime->PHOSGeometry()->AbsToRelNumbering(absId,relId);
518 Int_t module=relId[0];
520 Int_t column=relId[3];
524 if(absId > fEmcCrystals){ //digitize CPV only
526 //reading calibration data for cell absId.
527 //If no calibration DB found, accept default values.
529 if(gime->CalibData()) {
530 fADCpedestalCpv = gime->CalibData()->GetADCpedestalCpv(module,column,row);
531 fADCchanelCpv = gime->CalibData()->GetADCchannelCpv( module,column,row);
534 channel = (Int_t) TMath::Ceil((charge - fADCpedestalCpv)/fADCchanelCpv) ;
535 if(channel > fNADCcpv ) channel = fNADCcpv ;
540 //____________________________________________________________________________
541 void AliPHOSDigitizer::Exec(Option_t *option)
543 // Steering method to process digitization for events
544 // in the range from fFirstEvent to fLastEvent.
545 // This range is optionally set by SetEventRange().
546 // if fLastEvent=-1, then process events until the end.
547 // by default fLastEvent = fFirstEvent (process only one event)
549 if (!fInit) { // to prevent overwrite existing file
550 AliError(Form("Give a version name different from %s",
551 fEventFolderName.Data() )) ;
555 if (strstr(option,"print")) {
560 if(strstr(option,"tim"))
561 gBenchmark->Start("PHOSDigitizer");
563 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
565 // Post Digitizer to the white board
566 gime->PostDigitizer(this) ;
568 if (fLastEvent == -1)
569 fLastEvent = gime->MaxEvent() - 1 ;
571 fLastEvent = fFirstEvent ;
573 Int_t nEvents = fLastEvent - fFirstEvent + 1;
577 for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
579 gime->Event(ievent,"S") ;
581 Digitize(ievent) ; //Add prepared SDigits to digits and add the noise
585 if(strstr(option,"deb"))
588 //increment the total number of Digits per run
589 fDigitsInRun += gime->Digits()->GetEntriesFast() ;
592 gime->PhosLoader()->CleanDigitizer();
594 if(strstr(option,"tim")){
595 gBenchmark->Stop("PHOSDigitizer");
597 message = " took %f seconds for Digitizing %f seconds per event\n" ;
598 AliInfo(Form( message.Data(),
599 gBenchmark->GetCpuTime("PHOSDigitizer"),
600 gBenchmark->GetCpuTime("PHOSDigitizer")/nEvents ));
604 //____________________________________________________________________________
605 Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks) const
607 // Returns the shortest time among all time ticks
609 ticks->Sort() ; //Sort in accordance with times of ticks
611 AliPHOSTick * ctick = (AliPHOSTick *) it.Next() ;
612 Float_t time = ctick->CrossingTime(fTimeThreshold) ;
615 while((t=(AliPHOSTick*) it.Next())){
616 if(t->GetTime() < time) //This tick starts before crossing
621 time = ctick->CrossingTime(fTimeThreshold) ;
626 //____________________________________________________________________________
627 Bool_t AliPHOSDigitizer::Init()
629 // Makes all memory allocations
631 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
633 AliFatal(Form("Could not obtain the Getter object for file %s and event %s !",
634 GetTitle(), fEventFolderName.Data()));
638 const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
640 fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ;
642 TString opt("Digits") ;
643 if(gime->VersionExists(opt) ) {
644 AliError(Form("Give a version name different from %s",
645 fEventFolderName.Data() )) ;
650 fLastEvent = fFirstEvent ;
652 fInput = fManager->GetNinputs() ;
656 fInputFileNames = new TString[fInput] ;
657 fEventNames = new TString[fInput] ;
658 fInputFileNames[0] = GetTitle() ;
659 fEventNames[0] = fEventFolderName.Data() ;
661 for (index = 1 ; index < fInput ; index++) {
662 fInputFileNames[index] = dynamic_cast<AliStream*>(fManager->GetInputStream(index))->GetFileName(0);
663 TString tempo = fManager->GetInputFolderName(index) ;
664 fEventNames[index] = tempo.Remove(tempo.Length()-1) ; // strip of the stream number added by fManager
667 //to prevent cleaning of this object while GetEvent is called
668 gime->PhosLoader()->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
673 //____________________________________________________________________________
674 void AliPHOSDigitizer::InitParameters()
676 // Set initial parameters Digitizer
678 fPinNoise = 0.004 ; // [GeV]
679 fEMCDigitThreshold = 0.012 ; // [GeV]
680 fCPVNoise = 0.01; // [aux units]
681 fCPVDigitThreshold = 0.09 ; // [aux units]
682 fTimeResolution = 0.5e-9 ; // [sec]
683 fTimeSignalLength = 1.0e-9 ; // [sec]
685 fADCchanelEmc = 1.0; // Coefficient between real and measured energies in EMC
686 fADCpedestalEmc = 0. ; //
687 fNADCemc = (Int_t) TMath::Power(2,16) ; // number of channels in EMC ADC
689 fADCchanelCpv = 0.0012 ; // width of one ADC channel in CPV 'popugais'
690 fADCpedestalCpv = 0.012 ; //
691 fNADCcpv = (Int_t) TMath::Power(2,12); // number of channels in CPV ADC
693 // fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
694 fTimeThreshold = 0.001 ; // [GeV]
695 SetEventRange(0,-1) ;
699 //__________________________________________________________________
700 void AliPHOSDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
702 // Allows to produce digits by superimposing background and signal event.
703 // It is assumed, that headers file with SIGNAL events is opened in
705 // Sets the BACKGROUND event, with which the SIGNAL event is to be mixed
706 // Thus we avoid writing (changing) huge and expensive
707 // backgound files: all output will be writen into SIGNAL, i.e.
708 // opened in constructor file.
710 // One can open as many files to mix with as one needs.
711 // However only Sdigits with the same name (i.e. constructed with the same SDigitizer)
714 if( strcmp(fEventFolderName, "") == 0 )
718 Warning("MixWith", "Cannot use this method with AliRunDigitizer\n" ) ;
721 // looking for file which contains AliRun
722 if (gSystem->AccessPathName(alirunFileName)) {// file does not exist
723 AliError(Form("File %s does not exist!", alirunFileName.Data())) ;
726 // looking for the file which contains SDigits
727 AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
728 TString fileName( gime->GetSDigitsFileName() ) ;
729 if ( eventFolderName != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
730 fileName = fileName.ReplaceAll(".root", "") + "_" + eventFolderName + ".root" ;
731 if ( (gSystem->AccessPathName(fileName)) ) {
732 AliError(Form("The file %s does not exist!", fileName.Data())) ;
735 // need to increase the arrays
736 TString tempo = fInputFileNames[fInput-1] ;
737 delete [] fInputFileNames ;
738 fInputFileNames = new TString[fInput+1] ;
739 fInputFileNames[fInput-1] = tempo ;
741 tempo = fEventNames[fInput-1] ;
742 delete [] fEventNames ;
743 fEventNames = new TString[fInput+1] ;
744 fEventNames[fInput-1] = tempo ;
746 fInputFileNames[fInput] = alirunFileName ;
747 fEventNames[fInput] = eventFolderName ;
751 //__________________________________________________________________
752 void AliPHOSDigitizer::Print(const Option_t *)const
754 // Print Digitizer's parameters
755 AliInfo(Form("\n------------------- %s -------------", GetName() )) ;
756 if( strcmp(fEventFolderName.Data(), "") != 0 ){
757 printf(" Writing Digits to branch with title %s\n", fEventFolderName.Data()) ;
761 nStreams = GetNInputStreams() ;
766 for (index = 0 ; index < nStreams ; index++) {
767 TString tempo(fEventNames[index]) ;
769 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[index], tempo) ;
770 TString fileName( gime->GetSDigitsFileName() ) ;
771 if ( fEventNames[index] != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
772 fileName = fileName.ReplaceAll(".root", "") + "_" + fEventNames[index] + ".root" ;
773 printf ("Adding SDigits from %s %s\n", fInputFileNames[index].Data(), fileName.Data()) ;
775 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
776 printf("\nWriting digits to %s", gime->GetDigitsFileName().Data()) ;
778 printf("\nWith following parameters:\n") ;
779 printf(" Electronics noise in EMC (fPinNoise) = %f GeV\n", fPinNoise ) ;
780 printf(" Threshold in EMC (fEMCDigitThreshold) = %f GeV\n", fEMCDigitThreshold ) ;
781 printf(" Noise in CPV (fCPVNoise) = %f aux units\n", fCPVNoise ) ;
782 printf(" Threshold in CPV (fCPVDigitThreshold) = %f aux units\n",fCPVDigitThreshold ) ;
783 printf(" ---------------------------------------------------\n") ;
786 AliInfo(Form("AliPHOSDigitizer not initialized" )) ;
790 //__________________________________________________________________
791 void AliPHOSDigitizer::PrintDigits(Option_t * option)
793 // Print a table of digits
795 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
796 TClonesArray * digits = gime->Digits() ;
798 AliInfo(Form("%d", digits->GetEntriesFast())) ;
799 printf("\nevent %d", gAlice->GetEvNumber()) ;
800 printf("\n Number of entries in Digits list %d", digits->GetEntriesFast() ) ;
803 if(strstr(option,"all")||strstr(option,"EMC")){
805 AliPHOSDigit * digit;
806 printf("\nEMC digits (with primaries):\n") ;
807 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
808 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
810 for (index = 0 ; (index < digits->GetEntriesFast()) &&
811 (dynamic_cast<AliPHOSDigit *>(digits->At(index))->GetId() <= maxEmc) ; index++) {
812 digit = (AliPHOSDigit * ) digits->At(index) ;
813 if(digit->GetNprimary() == 0)
815 // printf("%6d %8d %6.5e %4d %2d :",
816 // digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ; // YVK
817 printf("%6d %.4f %6.5e %4d %2d :",
818 digit->GetId(), digit->GetEnergy(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
820 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
821 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
827 if(strstr(option,"all")||strstr(option,"CPV")){
829 //loop over CPV digits
830 AliPHOSDigit * digit;
831 printf("\nCPV digits (with primaries):\n") ;
832 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
833 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
835 for (index = 0 ; index < digits->GetEntriesFast(); index++) {
836 digit = (AliPHOSDigit * ) digits->At(index) ;
837 if(digit->GetId() > maxEmc){
838 printf("%6d %8d %4d %2d :",
839 digit->GetId(), digit->GetAmp(), digit->GetIndexInList(), digit->GetNprimary()) ;
841 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
842 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
851 //__________________________________________________________________
852 Float_t AliPHOSDigitizer::TimeOfNoise(void) const
853 { // Calculates the time signal generated by noise
854 //PH Info("TimeOfNoise", "Change me") ;
855 return gRandom->Rndm() * 1.28E-5;
858 //__________________________________________________________________
859 void AliPHOSDigitizer::Unload()
863 for(i = 1 ; i < fInput ; i++){
864 TString tempo(fEventNames[i]) ;
866 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
867 gime->PhosLoader()->UnloadSDigits() ;
870 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
871 gime->PhosLoader()->UnloadDigits() ;
874 //____________________________________________________________________________
875 void AliPHOSDigitizer::WriteDigits()
878 // Makes TreeD in the output file.
879 // Check if branch already exists:
880 // if yes, exit without writing: ROOT TTree does not support overwriting/updating of
881 // already existing branches.
882 // else creates branch with Digits, named "PHOS", title "...",
883 // and branch "AliPHOSDigitizer", with the same title to keep all the parameters
884 // and names of files, from which digits are made.
886 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
887 const TClonesArray * digits = gime->Digits() ;
888 TTree * treeD = gime->TreeD();
890 // -- create Digits branch
891 Int_t bufferSize = 32000 ;
892 TBranch * digitsBranch = treeD->Branch("PHOS","TClonesArray",&digits,bufferSize);
893 digitsBranch->SetTitle(fEventFolderName);
894 digitsBranch->Fill() ;
896 gime->WriteDigits("OVERWRITE");
897 gime->WriteDigitizer("OVERWRITE");