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.93 2006/10/17 13:17:01 kharlov
22 * Replace AliInfo by AliDebug
24 * Revision 1.92 2006/08/28 10:01:56 kharlov
25 * Effective C++ warnings fixed (Timur Pocheptsov)
27 * Revision 1.91 2006/04/29 20:25:30 hristov
28 * Decalibration is implemented (Yu.Kharlov)
30 * Revision 1.90 2006/04/22 10:30:17 hristov
31 * Add fEnergy to AliPHOSDigit and operate with EMC amplitude in energy units (Yu.Kharlov)
33 * Revision 1.89 2006/04/11 15:22:59 hristov
34 * run number in query set to -1: forces AliCDBManager to use its run number (A.Colla)
36 * Revision 1.88 2006/03/13 14:05:43 kharlov
37 * Calibration objects for EMC and CPV
39 * Revision 1.87 2005/08/24 15:33:49 kharlov
40 * Calibration data for raw digits
42 * Revision 1.86 2005/07/12 20:07:35 hristov
43 * Changes needed to run simulation and reconstrruction in the same AliRoot session
45 * Revision 1.85 2005/05/28 14:19:04 schutz
46 * Compilation warnings fixed by T.P.
50 //_________________________________________________________________________
51 //*-- Author : Dmitri Peressounko (SUBATECH & Kurchatov Institute)
52 //////////////////////////////////////////////////////////////////////////////
53 // This TTask performs digitization of Summable digits (in the PHOS case it is just
54 // the sum of contributions from all primary particles into a given cell).
55 // In addition it performs mixing of summable digits from different events.
56 // The name of the TTask is also the title of the branch that will contain
57 // the created SDigits
58 // The title of the TTAsk is the name of the file that contains the hits from
59 // which the SDigits are created
61 // For each event two branches are created in TreeD:
62 // "PHOS" - list of digits
63 // "AliPHOSDigitizer" - AliPHOSDigitizer with all parameters used in digitization
65 // Note, that one can set a title for new digits branch, and repeat digitization with
66 // another set of parameters.
69 // root[0] AliPHOSDigitizer * d = new AliPHOSDigitizer() ;
70 // root[1] d->ExecuteTask()
71 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
72 // //Digitizes SDigitis in all events found in file galice.root
74 // root[2] AliPHOSDigitizer * d1 = new AliPHOSDigitizer("galice1.root") ;
75 // // Will read sdigits from galice1.root
76 // root[3] d1->MixWith("galice2.root")
77 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
78 // // Reads another set of sdigits from galice2.root
79 // root[3] d1->MixWith("galice3.root")
80 // // Reads another set of sdigits from galice3.root
81 // root[4] d->ExecuteTask("deb timing")
82 // // Reads SDigits from files galice1.root, galice2.root ....
83 // // mixes them and stores produced Digits in file galice1.root
84 // // deb - prints number of produced digits
85 // // deb all - prints list of produced digits
86 // // timing - prints time used for digitization
89 // --- ROOT system ---
92 #include "TBenchmark.h"
95 // --- Standard library ---
97 // --- AliRoot header files ---
99 #include "AliRunDigitizer.h"
100 #include "AliPHOSDigit.h"
101 #include "AliPHOSGetter.h"
102 #include "AliPHOSDigitizer.h"
103 #include "AliPHOSSDigitizer.h"
104 #include "AliPHOSGeometry.h"
105 #include "AliPHOSTick.h"
107 ClassImp(AliPHOSDigitizer)
110 //____________________________________________________________________________
111 AliPHOSDigitizer::AliPHOSDigitizer() :
117 fInputFileNames(0x0),
121 fEMCDigitThreshold(0.f),
123 fCPVDigitThreshold(0.f),
124 fTimeResolution(0.f),
126 fTimeSignalLength(0.f),
128 fADCpedestalEmc(0.f),
131 fADCpedestalCpv(0.f),
133 fEventFolderName(""),
139 fManager = 0 ; // We work in the standalong mode
142 //____________________________________________________________________________
143 AliPHOSDigitizer::AliPHOSDigitizer(TString alirunFileName,
144 TString eventFolderName):
145 AliDigitizer("PHOS"+AliConfig::Instance()->GetDigitizerTaskName(), alirunFileName),
146 fDefaultInit(kFALSE),
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(eventFolderName),
173 fDefaultInit = kFALSE ;
174 fManager = 0 ; // We work in the standalong mode
177 //____________________________________________________________________________
178 AliPHOSDigitizer::AliPHOSDigitizer(const AliPHOSDigitizer & d) :
180 fDefaultInit(d.fDefaultInit),
181 fDigitsInRun(d.fDigitsInRun),
184 fInputFileNames(0x0),//?
186 fEmcCrystals(d.fEmcCrystals),
187 fPinNoise(d.fPinNoise),
188 fEMCDigitThreshold(d.fEMCDigitThreshold),
189 fCPVNoise(d.fCPVNoise),
190 fCPVDigitThreshold(d.fCPVDigitThreshold),
191 fTimeResolution(d.fTimeResolution),
192 fTimeThreshold(d.fTimeThreshold),
193 fTimeSignalLength(d.fTimeSignalLength),
194 fADCchanelEmc(d.fADCchanelEmc),
195 fADCpedestalEmc(d.fADCpedestalEmc),
196 fNADCemc(d.fNADCemc),
197 fADCchanelCpv(d.fADCchanelCpv),
198 fADCpedestalCpv(d.fADCpedestalCpv),
199 fNADCcpv(d.fNADCcpv),
200 fEventFolderName(d.fEventFolderName),
201 fFirstEvent(d.fFirstEvent),
202 fLastEvent(d.fLastEvent)
205 SetName(d.GetName()) ;
206 SetTitle(d.GetTitle()) ;
209 //____________________________________________________________________________
210 AliPHOSDigitizer::AliPHOSDigitizer(AliRunDigitizer * rd) :
211 AliDigitizer(rd,"PHOS"+AliConfig::Instance()->GetDigitizerTaskName()),
212 fDefaultInit(kFALSE),
216 fInputFileNames(0x0),
220 fEMCDigitThreshold(0.f),
222 fCPVDigitThreshold(0.f),
223 fTimeResolution(0.f),
225 fTimeSignalLength(0.f),
227 fADCpedestalEmc(0.f),
230 fADCpedestalCpv(0.f),
232 fEventFolderName(fManager->GetInputFolderName(0)),
236 // ctor Init() is called by RunDigitizer
238 SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
240 fDefaultInit = kFALSE ;
243 //____________________________________________________________________________
244 AliPHOSDigitizer::~AliPHOSDigitizer()
246 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
248 // Clean Digitizer from the white board
249 gime->PhosLoader()->CleanDigitizer() ;
251 delete [] fInputFileNames ;
252 delete [] fEventNames ;
256 //____________________________________________________________________________
257 void AliPHOSDigitizer::Digitize(Int_t event)
260 // Makes the digitization of the collected summable digits.
261 // It first creates the array of all PHOS modules
262 // filled with noise (different for EMC, and CPV) and
263 // then adds contributions from SDigits.
264 // This design avoids scanning over the list of digits to add
265 // contribution to new SDigits only.
267 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
268 Int_t ReadEvent = event ;
270 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetCurrentEventNumber() ;
271 AliDebug(1,Form("Adding event %d from input stream 0 %s %s",
272 ReadEvent, GetTitle(), fEventFolderName.Data())) ;
273 gime->Event(ReadEvent, "S") ;
274 TClonesArray * digits = gime->Digits() ;
277 const AliPHOSGeometry *geom = gime->PHOSGeometry() ;
278 //Making digits with noise, first EMC
279 Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
284 nCPV = nEMC + geom->GetNumberOfCPVPadsZ() * geom->GetNumberOfCPVPadsPhi() * geom->GetNModules() ;
286 digits->Expand(nCPV) ;
288 // get first the sdigitizer from the tasks list
289 if ( !gime->SDigitizer() )
290 gime->LoadSDigitizer();
291 AliPHOSSDigitizer * sDigitizer = gime->SDigitizer();
294 AliFatal(Form("SDigitizer with name %s %s not found",
295 GetTitle(), fEventFolderName.Data() )) ;
297 //take all the inputs to add together and load the SDigits
298 TObjArray * sdigArray = new TObjArray(fInput) ;
299 sdigArray->AddAt(gime->SDigits(), 0) ;
301 for(i = 1 ; i < fInput ; i++){
302 TString tempo(fEventNames[i]) ;
304 AliPHOSGetter * gime1 = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
306 ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(i))->GetCurrentEventNumber() ;
307 AliInfo(Form("Adding event %d from input stream %d %s %s",
308 ReadEvent, i, fInputFileNames[i].Data(), tempo.Data())) ;
309 gime1->Event(ReadEvent,"S");
310 sdigArray->AddAt(gime1->SDigits(), i) ;
313 //Find the first crystal with signal
314 Int_t nextSig = 200000 ;
315 TClonesArray * sdigits ;
316 for(i = 0 ; i < fInput ; i++){
317 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
318 if ( !sdigits->GetEntriesFast() )
320 Int_t curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(0))->GetId() ;
321 if(curNext < nextSig)
325 TArrayI index(fInput) ;
326 index.Reset() ; //Set all indexes to zero
328 AliPHOSDigit * digit ;
329 AliPHOSDigit * curSDigit ;
331 TClonesArray * ticks = new TClonesArray("AliPHOSTick",1000) ;
333 //Put Noise contribution
334 for(absID = 1 ; absID <= nEMC ; absID++){
335 Float_t noise = gRandom->Gaus(0., fPinNoise) ;
336 // YVK: do not digitize amplitudes for EMC
337 // new((*digits)[absID-1]) AliPHOSDigit( -1, absID, sDigitizer->Digitize(noise), TimeOfNoise() ) ;
338 new((*digits)[absID-1]) AliPHOSDigit( -1, absID, noise, TimeOfNoise() ) ;
339 //look if we have to add signal?
340 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
343 //Add SDigits from all inputs
346 Float_t a = digit->GetEnergy() ;
347 Float_t b = TMath::Abs( a / fTimeSignalLength) ;
348 //Mark the beginning of the signal
349 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime(),0, b);
350 //Mark the end of the signal
351 new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime()+fTimeSignalLength, -a, -b);
354 for(i = 0 ; i < fInput ; i++){
355 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
356 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
359 //May be several digits will contribute from the same input
360 while(curSDigit && curSDigit->GetId() == absID){
361 //Shift primary to separate primaries belonging different inputs
362 Int_t primaryoffset ;
364 primaryoffset = fManager->GetMask(i) ;
366 primaryoffset = 10000000*i ;
367 curSDigit->ShiftPrimary(primaryoffset) ;
369 a = curSDigit->GetEnergy() ;
370 b = a /fTimeSignalLength ;
371 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime(),0, b);
372 new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b);
374 *digit += *curSDigit ; //add energies
377 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
378 curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
384 //calculate and set time
385 Float_t time = FrontEdgeTime(ticks) ;
386 digit->SetTime(time) ;
388 //Find next signal module
390 for(i = 0 ; i < fInput ; i++){
391 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
392 Int_t curNext = nextSig ;
393 if(sdigits->GetEntriesFast() > index[i] ){
394 curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(index[i]))->GetId() ;
396 if(curNext < nextSig) nextSig = curNext ;
404 //Now CPV digits (different noise and no timing)
405 for(absID = nEMC+1; absID <= nCPV; absID++){
406 Float_t noise = gRandom->Gaus(0., fCPVNoise) ;
407 new((*digits)[absID-1]) AliPHOSDigit( -1,absID,sDigitizer->Digitize(noise), TimeOfNoise() ) ;
408 //look if we have to add signal?
410 digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
411 //Add SDigits from all inputs
412 for(i = 0 ; i < fInput ; i++){
413 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
414 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
418 //May be several digits will contribute from the same input
419 while(curSDigit && curSDigit->GetId() == absID){
420 //Shift primary to separate primaries belonging different inputs
421 Int_t primaryoffset ;
423 primaryoffset = fManager->GetMask(i) ;
425 primaryoffset = 10000000*i ;
426 curSDigit->ShiftPrimary(primaryoffset) ;
429 *digit += *curSDigit ;
431 if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
432 curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i]) ) ;
438 //Find next signal module
440 for(i = 0 ; i < fInput ; i++){
441 sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
442 Int_t curNext = nextSig ;
443 if(sdigits->GetEntriesFast() > index[i] )
444 curNext = dynamic_cast<AliPHOSDigit *>( sdigits->At(index[i]) )->GetId() ;
445 if(curNext < nextSig) nextSig = curNext ;
451 delete sdigArray ; //We should not delete its contents
453 //remove digits below thresholds
454 for(i = 0 ; i < nEMC ; i++){
455 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
456 DecalibrateEMC(digit);
457 if(digit->GetEnergy() < fEMCDigitThreshold)
458 digits->RemoveAt(i) ;
460 digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ;
464 for(i = nEMC; i < nCPV ; i++)
465 // if( sDigitizer->Calibrate( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetAmp() ) < fCPVDigitThreshold )
466 if( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetEnergy() < fCPVDigitThreshold )
467 digits->RemoveAt(i) ;
471 Int_t ndigits = digits->GetEntriesFast() ;
472 digits->Expand(ndigits) ;
474 //Set indexes in list of digits and make true digitization of the energy
475 for (i = 0 ; i < ndigits ; i++) {
476 digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
477 digit->SetIndexInList(i) ;
478 if(digit->GetId() > fEmcCrystals){ //digitize CPV only
479 digit->SetAmp(DigitizeCPV(digit->GetEnergy(),digit->GetId()) ) ;
484 //____________________________________________________________________________
485 void AliPHOSDigitizer::DecalibrateEMC(AliPHOSDigit *digit)
487 // Decalibrate EMC digit, i.e. change its energy by a factor read from CDB
489 AliPHOSGetter* gime = AliPHOSGetter::Instance();
491 if(!gime->CalibData()) {
492 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1);
493 gime->SetCalibData(cdb);
496 //Determine rel.position of the cell absolute ID
498 gime->PHOSGeometry()->AbsToRelNumbering(digit->GetId(),relId);
499 Int_t module=relId[0];
501 Int_t column=relId[3];
502 Float_t decalibration = gime->CalibData()->GetADCchannelEmc(module,column,row);
503 Float_t energy = digit->GetEnergy() * decalibration;
504 digit->SetEnergy(energy);
506 //____________________________________________________________________________
507 Int_t AliPHOSDigitizer::DigitizeCPV(Float_t charge, Int_t absId)
509 // Returns digitized value of the CPV charge in a pad absId
511 AliPHOSGetter* gime = AliPHOSGetter::Instance();
513 if(!gime->CalibData()) {
514 AliPHOSCalibData* cdb = new AliPHOSCalibData(-1); // use AliCDBManager's run number
515 gime->SetCalibData(cdb);
518 //Determine rel.position of the cell absId
520 gime->PHOSGeometry()->AbsToRelNumbering(absId,relId);
521 Int_t module=relId[0];
523 Int_t column=relId[3];
527 if(absId > fEmcCrystals){ //digitize CPV only
529 //reading calibration data for cell absId.
530 //If no calibration DB found, accept default values.
532 if(gime->CalibData()) {
533 fADCpedestalCpv = gime->CalibData()->GetADCpedestalCpv(module,column,row);
534 fADCchanelCpv = gime->CalibData()->GetADCchannelCpv( module,column,row);
537 channel = (Int_t) TMath::Ceil((charge - fADCpedestalCpv)/fADCchanelCpv) ;
538 if(channel > fNADCcpv ) channel = fNADCcpv ;
543 //____________________________________________________________________________
544 void AliPHOSDigitizer::Exec(Option_t *option)
546 // Steering method to process digitization for events
547 // in the range from fFirstEvent to fLastEvent.
548 // This range is optionally set by SetEventRange().
549 // if fLastEvent=-1, then process events until the end.
550 // by default fLastEvent = fFirstEvent (process only one event)
552 if (!fInit) { // to prevent overwrite existing file
553 AliError(Form("Give a version name different from %s",
554 fEventFolderName.Data() )) ;
558 if (strstr(option,"print")) {
563 if(strstr(option,"tim"))
564 gBenchmark->Start("PHOSDigitizer");
566 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
568 // Post Digitizer to the white board
569 gime->PostDigitizer(this) ;
571 if (fLastEvent == -1)
572 fLastEvent = gime->MaxEvent() - 1 ;
574 fLastEvent = fFirstEvent ;
576 Int_t nEvents = fLastEvent - fFirstEvent + 1;
580 for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
582 gime->Event(ievent,"S") ;
584 Digitize(ievent) ; //Add prepared SDigits to digits and add the noise
588 if(strstr(option,"deb"))
591 //increment the total number of Digits per run
592 fDigitsInRun += gime->Digits()->GetEntriesFast() ;
595 gime->PhosLoader()->CleanDigitizer();
597 if(strstr(option,"tim")){
598 gBenchmark->Stop("PHOSDigitizer");
600 message = " took %f seconds for Digitizing %f seconds per event\n" ;
601 AliInfo(Form( message.Data(),
602 gBenchmark->GetCpuTime("PHOSDigitizer"),
603 gBenchmark->GetCpuTime("PHOSDigitizer")/nEvents ));
607 //____________________________________________________________________________
608 Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks) const
610 // Returns the shortest time among all time ticks
612 ticks->Sort() ; //Sort in accordance with times of ticks
614 AliPHOSTick * ctick = (AliPHOSTick *) it.Next() ;
615 Float_t time = ctick->CrossingTime(fTimeThreshold) ;
618 while((t=(AliPHOSTick*) it.Next())){
619 if(t->GetTime() < time) //This tick starts before crossing
624 time = ctick->CrossingTime(fTimeThreshold) ;
629 //____________________________________________________________________________
630 Bool_t AliPHOSDigitizer::Init()
632 // Makes all memory allocations
634 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
636 AliFatal(Form("Could not obtain the Getter object for file %s and event %s !",
637 GetTitle(), fEventFolderName.Data()));
641 const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
643 fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ;
645 TString opt("Digits") ;
646 if(gime->VersionExists(opt) ) {
647 AliError(Form("Give a version name different from %s",
648 fEventFolderName.Data() )) ;
653 fLastEvent = fFirstEvent ;
655 fInput = fManager->GetNinputs() ;
659 fInputFileNames = new TString[fInput] ;
660 fEventNames = new TString[fInput] ;
661 fInputFileNames[0] = GetTitle() ;
662 fEventNames[0] = fEventFolderName.Data() ;
664 for (index = 1 ; index < fInput ; index++) {
665 fInputFileNames[index] = dynamic_cast<AliStream*>(fManager->GetInputStream(index))->GetFileName(0);
666 TString tempo = fManager->GetInputFolderName(index) ;
667 fEventNames[index] = tempo.Remove(tempo.Length()-1) ; // strip of the stream number added by fManager
670 //to prevent cleaning of this object while GetEvent is called
671 gime->PhosLoader()->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
676 //____________________________________________________________________________
677 void AliPHOSDigitizer::InitParameters()
679 // Set initial parameters Digitizer
681 fPinNoise = 0.004 ; // [GeV]
682 fEMCDigitThreshold = 0.012 ; // [GeV]
683 fCPVNoise = 0.01; // [aux units]
684 fCPVDigitThreshold = 0.09 ; // [aux units]
685 fTimeResolution = 0.5e-9 ; // [sec]
686 fTimeSignalLength = 1.0e-9 ; // [sec]
688 fADCchanelEmc = 1.0; // Coefficient between real and measured energies in EMC
689 fADCpedestalEmc = 0. ; //
690 fNADCemc = (Int_t) TMath::Power(2,16) ; // number of channels in EMC ADC
692 fADCchanelCpv = 0.0012 ; // width of one ADC channel in CPV 'popugais'
693 fADCpedestalCpv = 0.012 ; //
694 fNADCcpv = (Int_t) TMath::Power(2,12); // number of channels in CPV ADC
696 // fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
697 fTimeThreshold = 0.001 ; // [GeV]
698 SetEventRange(0,-1) ;
702 //__________________________________________________________________
703 void AliPHOSDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
705 // Allows to produce digits by superimposing background and signal event.
706 // It is assumed, that headers file with SIGNAL events is opened in
708 // Sets the BACKGROUND event, with which the SIGNAL event is to be mixed
709 // Thus we avoid writing (changing) huge and expensive
710 // backgound files: all output will be writen into SIGNAL, i.e.
711 // opened in constructor file.
713 // One can open as many files to mix with as one needs.
714 // However only Sdigits with the same name (i.e. constructed with the same SDigitizer)
717 if( strcmp(fEventFolderName, "") == 0 )
721 Warning("MixWith", "Cannot use this method with AliRunDigitizer\n" ) ;
724 // looking for file which contains AliRun
725 if (gSystem->AccessPathName(alirunFileName)) {// file does not exist
726 AliError(Form("File %s does not exist!", alirunFileName.Data())) ;
729 // looking for the file which contains SDigits
730 AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
731 TString fileName( gime->GetSDigitsFileName() ) ;
732 if ( eventFolderName != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
733 fileName = fileName.ReplaceAll(".root", "") + "_" + eventFolderName + ".root" ;
734 if ( (gSystem->AccessPathName(fileName)) ) {
735 AliError(Form("The file %s does not exist!", fileName.Data())) ;
738 // need to increase the arrays
739 TString tempo = fInputFileNames[fInput-1] ;
740 delete [] fInputFileNames ;
741 fInputFileNames = new TString[fInput+1] ;
742 fInputFileNames[fInput-1] = tempo ;
744 tempo = fEventNames[fInput-1] ;
745 delete [] fEventNames ;
746 fEventNames = new TString[fInput+1] ;
747 fEventNames[fInput-1] = tempo ;
749 fInputFileNames[fInput] = alirunFileName ;
750 fEventNames[fInput] = eventFolderName ;
754 //__________________________________________________________________
755 void AliPHOSDigitizer::Print(const Option_t *)const
757 // Print Digitizer's parameters
758 AliInfo(Form("\n------------------- %s -------------", GetName() )) ;
759 if( strcmp(fEventFolderName.Data(), "") != 0 ){
760 printf(" Writing Digits to branch with title %s\n", fEventFolderName.Data()) ;
764 nStreams = GetNInputStreams() ;
769 for (index = 0 ; index < nStreams ; index++) {
770 TString tempo(fEventNames[index]) ;
772 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[index], tempo) ;
773 TString fileName( gime->GetSDigitsFileName() ) ;
774 if ( fEventNames[index] != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
775 fileName = fileName.ReplaceAll(".root", "") + "_" + fEventNames[index] + ".root" ;
776 printf ("Adding SDigits from %s %s\n", fInputFileNames[index].Data(), fileName.Data()) ;
778 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
779 printf("\nWriting digits to %s", gime->GetDigitsFileName().Data()) ;
781 printf("\nWith following parameters:\n") ;
782 printf(" Electronics noise in EMC (fPinNoise) = %f GeV\n", fPinNoise ) ;
783 printf(" Threshold in EMC (fEMCDigitThreshold) = %f GeV\n", fEMCDigitThreshold ) ;
784 printf(" Noise in CPV (fCPVNoise) = %f aux units\n", fCPVNoise ) ;
785 printf(" Threshold in CPV (fCPVDigitThreshold) = %f aux units\n",fCPVDigitThreshold ) ;
786 printf(" ---------------------------------------------------\n") ;
789 AliInfo(Form("AliPHOSDigitizer not initialized" )) ;
793 //__________________________________________________________________
794 void AliPHOSDigitizer::PrintDigits(Option_t * option)
796 // Print a table of digits
798 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
799 TClonesArray * digits = gime->Digits() ;
801 AliInfo(Form("%d", digits->GetEntriesFast())) ;
802 printf("\nevent %d", gAlice->GetEvNumber()) ;
803 printf("\n Number of entries in Digits list %d", digits->GetEntriesFast() ) ;
806 if(strstr(option,"all")||strstr(option,"EMC")){
808 AliPHOSDigit * digit;
809 printf("\nEMC digits (with primaries):\n") ;
810 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
811 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
813 for (index = 0 ; (index < digits->GetEntriesFast()) &&
814 (dynamic_cast<AliPHOSDigit *>(digits->At(index))->GetId() <= maxEmc) ; index++) {
815 digit = (AliPHOSDigit * ) digits->At(index) ;
816 if(digit->GetNprimary() == 0)
818 // printf("%6d %8d %6.5e %4d %2d :",
819 // digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ; // YVK
820 printf("%6d %.4f %6.5e %4d %2d :",
821 digit->GetId(), digit->GetEnergy(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
823 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
824 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
830 if(strstr(option,"all")||strstr(option,"CPV")){
832 //loop over CPV digits
833 AliPHOSDigit * digit;
834 printf("\nCPV digits (with primaries):\n") ;
835 printf("\n Id Amplitude Time Index Nprim: Primaries list \n") ;
836 Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
838 for (index = 0 ; index < digits->GetEntriesFast(); index++) {
839 digit = (AliPHOSDigit * ) digits->At(index) ;
840 if(digit->GetId() > maxEmc){
841 printf("%6d %8d %4d %2d :",
842 digit->GetId(), digit->GetAmp(), digit->GetIndexInList(), digit->GetNprimary()) ;
844 for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
845 printf("%d ",digit->GetPrimary(iprimary+1) ) ;
854 //__________________________________________________________________
855 Float_t AliPHOSDigitizer::TimeOfNoise(void) const
856 { // Calculates the time signal generated by noise
857 //PH Info("TimeOfNoise", "Change me") ;
858 return gRandom->Rndm() * 1.28E-5;
861 //__________________________________________________________________
862 void AliPHOSDigitizer::Unload()
866 for(i = 1 ; i < fInput ; i++){
867 TString tempo(fEventNames[i]) ;
869 AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
870 gime->PhosLoader()->UnloadSDigits() ;
873 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
874 gime->PhosLoader()->UnloadDigits() ;
877 //____________________________________________________________________________
878 void AliPHOSDigitizer::WriteDigits()
881 // Makes TreeD in the output file.
882 // Check if branch already exists:
883 // if yes, exit without writing: ROOT TTree does not support overwriting/updating of
884 // already existing branches.
885 // else creates branch with Digits, named "PHOS", title "...",
886 // and branch "AliPHOSDigitizer", with the same title to keep all the parameters
887 // and names of files, from which digits are made.
889 AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
890 const TClonesArray * digits = gime->Digits() ;
891 TTree * treeD = gime->TreeD();
893 // -- create Digits branch
894 Int_t bufferSize = 32000 ;
895 TBranch * digitsBranch = treeD->Branch("PHOS","TClonesArray",&digits,bufferSize);
896 digitsBranch->SetTitle(fEventFolderName);
897 digitsBranch->Fill() ;
899 gime->WriteDigits("OVERWRITE");
900 gime->WriteDigitizer("OVERWRITE");