PHOS/AliPHOSDigitizer.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   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  **************************************************************************/
  15
  16 /* $Id$ */
  17
  18 /* History of cvs commits:
  19  *
  20  * $Log$
  21  * Revision 1.87  2005/08/24 15:33:49  kharlov
  22  * Calibration data for raw digits
  23  *
  24  * Revision 1.86  2005/07/12 20:07:35  hristov
  25  * Changes needed to run simulation and reconstrruction in the same AliRoot session
  26  *
  27  * Revision 1.85  2005/05/28 14:19:04  schutz
  28  * Compilation warnings fixed by T.P.
  29  *
  30  */
  31
  32 //_________________________________________________________________________
  33 //*-- Author :  Dmitri Peressounko (SUBATECH & Kurchatov Institute)
  34 //////////////////////////////////////////////////////////////////////////////
  35 // This TTask performs digitization of Summable digits (in the PHOS case it is just
  36 // the sum of contributions from all primary particles into a given cell).
  37 // In addition it performs mixing of summable digits from different events.
  38 // The name of the TTask is also the title of the branch that will contain
  39 // the created SDigits
  40 // The title of the TTAsk is the name of the file that contains the hits from
  41 // which the SDigits are created
  42 //
  43 // For each event two branches are created in TreeD:
  44 //   "PHOS" - list of digits
  45 //   "AliPHOSDigitizer" - AliPHOSDigitizer with all parameters used in digitization
  46 //
  47 // Note, that one can set a title for new digits branch, and repeat digitization with
  48 // another set of parameters.
  49 //
  50 // Use case:
  51 // root[0] AliPHOSDigitizer * d = new AliPHOSDigitizer() ;
  52 // root[1] d->ExecuteTask()
  53 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
  54 //                       //Digitizes SDigitis in all events found in file galice.root
  55 //
  56 // root[2] AliPHOSDigitizer * d1 = new AliPHOSDigitizer("galice1.root") ;
  57 //                       // Will read sdigits from galice1.root
  58 // root[3] d1->MixWith("galice2.root")
  59 // Warning in <TDatabasePDG::TDatabasePDG>: object already instantiated
  60 //                       // Reads another set of sdigits from galice2.root
  61 // root[3] d1->MixWith("galice3.root")
  62 //                       // Reads another set of sdigits from galice3.root
  63 // root[4] d->ExecuteTask("deb timing")
  64 //                       // Reads SDigits from files galice1.root, galice2.root ....
  65 //                       // mixes them and stores produced Digits in file galice1.root
  66 //                       // deb - prints number of produced digits
  67 //                       // deb all - prints list of produced digits
  68 //                       // timing  - prints time used for digitization
  69 //
  70
  71 // --- ROOT system ---
  72 #include "TTree.h"
  73 #include "TSystem.h"
  74 #include "TBenchmark.h"
  75 #include "TRandom.h"
  76
  77 // --- Standard library ---
  78
  79 // --- AliRoot header files ---
  80 #include "AliLog.h"
  81 #include "AliRunDigitizer.h"
  82 #include "AliPHOSDigit.h"
  83 #include "AliPHOSGetter.h"
  84 #include "AliPHOSDigitizer.h"
  85 #include "AliPHOSSDigitizer.h"
  86 #include "AliPHOSGeometry.h"
  87 #include "AliPHOSTick.h"
  88
  89 ClassImp(AliPHOSDigitizer)
  90
  91
  92 //____________________________________________________________________________
  93   AliPHOSDigitizer::AliPHOSDigitizer():AliDigitizer("",""),
  94                                        fInput(0),
  95                                        fInputFileNames(0x0),
  96                                        fEventNames(0x0)
  97 {
  98   // ctor
  99   InitParameters() ;
 100   fDefaultInit = kTRUE ;
 101   fManager = 0 ;                     // We work in the standalong mode
 102   fEventFolderName = "" ;
 103 }
 104
 105 //____________________________________________________________________________
 106 AliPHOSDigitizer::AliPHOSDigitizer(TString alirunFileName,
 107                                    TString eventFolderName):
 108   AliDigitizer("PHOS"+AliConfig::Instance()->GetDigitizerTaskName(),
 109                alirunFileName),
 110   fInputFileNames(0), fEventNames(0), fEventFolderName(eventFolderName)
 111 {
 112   // ctor
 113   InitParameters() ;
 114   Init() ;
 115   fDefaultInit = kFALSE ;
 116   fManager = 0 ;                     // We work in the standalong mode
 117 }
 118
 119 //____________________________________________________________________________
 120 AliPHOSDigitizer::AliPHOSDigitizer(const AliPHOSDigitizer & d)
 121   : AliDigitizer(d)
 122 {
 123   // copyy ctor
 124
 125   SetName(d.GetName()) ;
 126   SetTitle(d.GetTitle()) ;
 127   fPinNoise = d.fPinNoise ;
 128   fEMCDigitThreshold = d.fEMCDigitThreshold ;
 129   fCPVNoise          = d.fCPVNoise ;
 130   fCPVDigitThreshold = d.fCPVDigitThreshold ;
 131   fTimeResolution    = d.fTimeResolution ;
 132   fTimeThreshold     = d.fTimeThreshold ;
 133   fTimeSignalLength  = d.fTimeSignalLength ;
 134   fADCchanelEmc      = d.fADCchanelEmc ;
 135   fADCpedestalEmc    = d.fADCpedestalEmc ;
 136   fNADCemc           = d.fNADCemc ;
 137   fADCchanelCpv      = d.fADCchanelCpv ;
 138   fADCpedestalCpv    = d.fADCpedestalCpv ;
 139   fNADCcpv           = d.fNADCcpv ;
 140   fEventFolderName   = d.fEventFolderName;
 141 }
 142
 143 //____________________________________________________________________________
 144 AliPHOSDigitizer::AliPHOSDigitizer(AliRunDigitizer * rd):
 145  AliDigitizer(rd,"PHOS"+AliConfig::Instance()->GetDigitizerTaskName()),
 146  fEventFolderName(0)
 147 {
 148   // ctor Init() is called by RunDigitizer
 149   fManager = rd ;
 150   fEventFolderName = fManager->GetInputFolderName(0) ;
 151   SetTitle(dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetFileName(0));
 152   InitParameters() ;
 153   fDefaultInit = kFALSE ;
 154 }
 155
 156 //____________________________________________________________________________
 157   AliPHOSDigitizer::~AliPHOSDigitizer()
 158 {
 159   AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
 160
 161   // Clean Digitizer from the white board
 162   gime->PhosLoader()->CleanDigitizer() ;
 163   // dtor
 164   delete [] fInputFileNames ;
 165   delete [] fEventNames ;
 166
 167 }
 168
 169 //____________________________________________________________________________
 170 void AliPHOSDigitizer::Digitize(Int_t event)
 171 {
 172
 173   // Makes the digitization of the collected summable digits.
 174   //  It first creates the array of all PHOS modules
 175   //  filled with noise (different for EMC, and CPV) and
 176   //  then adds contributions from SDigits.
 177   // This design avoids scanning over the list of digits to add
 178   // contribution to new SDigits only.
 179
 180   AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
 181   Int_t ReadEvent = event ;
 182   if (fManager)
 183     ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(0))->GetCurrentEventNumber() ;
 184   AliInfo(Form("Adding event %d from input stream 0 %s %s",
 185                ReadEvent, GetTitle(), fEventFolderName.Data())) ;
 186   gime->Event(ReadEvent, "S") ;
 187   TClonesArray * digits = gime->Digits() ;
 188   digits->Clear() ;
 189
 190   const AliPHOSGeometry *geom = gime->PHOSGeometry() ;
 191   //Making digits with noise, first EMC
 192   Int_t nEMC = geom->GetNModules()*geom->GetNPhi()*geom->GetNZ();
 193
 194   Int_t nCPV ;
 195   Int_t absID ;
 196
 197   nCPV = nEMC + geom->GetNumberOfCPVPadsZ() * geom->GetNumberOfCPVPadsPhi() * geom->GetNModules() ;
 198
 199   digits->Expand(nCPV) ;
 200
 201   // get first the sdigitizer from the tasks list
 202   if ( !gime->SDigitizer() )
 203     gime->LoadSDigitizer();
 204   AliPHOSSDigitizer * sDigitizer = gime->SDigitizer();
 205
 206   if ( !sDigitizer )
 207     AliFatal(Form("SDigitizer with name %s %s not found",
 208                   GetTitle(), fEventFolderName.Data() )) ;
 209
 210   //take all the inputs to add together and load the SDigits
 211   TObjArray * sdigArray = new TObjArray(fInput) ;
 212   sdigArray->AddAt(gime->SDigits(), 0) ;
 213   Int_t i ;
 214   for(i = 1 ; i < fInput ; i++){
 215     TString tempo(fEventNames[i]) ;
 216     tempo += i ;
 217     AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
 218     if (fManager)
 219       ReadEvent = dynamic_cast<AliStream*>(fManager->GetInputStream(i))->GetCurrentEventNumber() ;
 220     AliInfo(Form("Adding event %d from input stream %d %s %s",
 221                  ReadEvent, i, fInputFileNames[i].Data(), tempo.Data())) ;
 222     gime->Event(ReadEvent,"S");
 223     sdigArray->AddAt(gime->SDigits(), i) ;
 224   }
 225
 226   //Find the first crystall with signal
 227   Int_t nextSig = 200000 ;
 228   TClonesArray * sdigits ;
 229   for(i = 0 ; i < fInput ; i++){
 230     sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
 231     if ( !sdigits->GetEntriesFast() )
 232       continue ;
 233     Int_t curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(0))->GetId() ;
 234     if(curNext < nextSig)
 235       nextSig = curNext ;
 236   }
 237
 238   TArrayI index(fInput) ;
 239   index.Reset() ;  //Set all indexes to zero
 240
 241   AliPHOSDigit * digit ;
 242   AliPHOSDigit * curSDigit ;
 243
 244   TClonesArray * ticks = new TClonesArray("AliPHOSTick",1000) ;
 245
 246   //Put Noise contribution
 247   for(absID = 1 ; absID <= nEMC ; absID++){
 248     Float_t noise = gRandom->Gaus(0., fPinNoise) ;
 249     new((*digits)[absID-1]) AliPHOSDigit( -1, absID, sDigitizer->Digitize(noise), TimeOfNoise() ) ;
 250     //look if we have to add signal?
 251     digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
 252
 253     if(absID==nextSig){
 254       //Add SDigits from all inputs
 255       ticks->Clear() ;
 256       Int_t contrib = 0 ;
 257       Float_t a = digit->GetAmp() ;
 258       Float_t b = TMath::Abs( a / fTimeSignalLength) ;
 259       //Mark the beginning of the signal
 260       new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime(),0, b);
 261       //Mark the end of the ignal
 262       new((*ticks)[contrib++]) AliPHOSTick(digit->GetTime()+fTimeSignalLength, -a, -b);
 263
 264       //loop over inputs
 265       for(i = 0 ; i < fInput ; i++){
 266         if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
 267           curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
 268         else
 269           curSDigit = 0 ;
 270         //May be several digits will contribute from the same input
 271         while(curSDigit && curSDigit->GetId() == absID){
 272           //Shift primary to separate primaries belonging different inputs
 273           Int_t primaryoffset ;
 274           if(fManager)
 275             primaryoffset = fManager->GetMask(i) ;
 276           else
 277             primaryoffset = 10000000*i ;
 278           curSDigit->ShiftPrimary(primaryoffset) ;
 279
 280           a = curSDigit->GetAmp() ;
 281           b = a /fTimeSignalLength ;
 282           new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime(),0, b);
 283           new((*ticks)[contrib++]) AliPHOSTick(curSDigit->GetTime()+fTimeSignalLength, -a, -b);
 284
 285           *digit = *digit + *curSDigit ;  //add energies
 286
 287           index[i]++ ;
 288           if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
 289             curSDigit = dynamic_cast<AliPHOSDigit*>(dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
 290           else
 291             curSDigit = 0 ;
 292         }
 293       }
 294
 295       //calculate and set time
 296       Float_t time = FrontEdgeTime(ticks) ;
 297       digit->SetTime(time) ;
 298
 299       //Find next signal module
 300       nextSig = 200000 ;
 301       for(i = 0 ; i < fInput ; i++){
 302         sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
 303         Int_t curNext = nextSig ;
 304         if(sdigits->GetEntriesFast() > index[i] ){
 305           curNext = dynamic_cast<AliPHOSDigit *>(sdigits->At(index[i]))->GetId() ;
 306         }
 307         if(curNext < nextSig) nextSig = curNext ;
 308       }
 309     }
 310   }
 311
 312   ticks->Delete() ;
 313   delete ticks ;
 314
 315   //Now CPV digits (different noise and no timing)
 316   for(absID = nEMC+1; absID <= nCPV; absID++){
 317     Float_t noise = gRandom->Gaus(0., fCPVNoise) ;
 318     new((*digits)[absID-1]) AliPHOSDigit( -1,absID,sDigitizer->Digitize(noise), TimeOfNoise() ) ;
 319     //look if we have to add signal?
 320     if(absID==nextSig){
 321       digit = dynamic_cast<AliPHOSDigit *>(digits->At(absID-1)) ;
 322       //Add SDigits from all inputs
 323       for(i = 0 ; i < fInput ; i++){
 324         if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
 325           curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i])) ;
 326         else
 327           curSDigit = 0 ;
 328
 329         //May be several digits will contribute from the same input
 330         while(curSDigit && curSDigit->GetId() == absID){
 331           //Shift primary to separate primaries belonging different inputs
 332           Int_t primaryoffset ;
 333           if(fManager)
 334             primaryoffset = fManager->GetMask(i) ;
 335           else
 336             primaryoffset = 10000000*i ;
 337           curSDigit->ShiftPrimary(primaryoffset) ;
 338
 339           //add energies
 340           *digit = *digit + *curSDigit ;
 341           index[i]++ ;
 342           if( dynamic_cast<TClonesArray *>(sdigArray->At(i))->GetEntriesFast() > index[i] )
 343             curSDigit = dynamic_cast<AliPHOSDigit*>( dynamic_cast<TClonesArray *>(sdigArray->At(i))->At(index[i]) ) ;
 344           else
 345             curSDigit = 0 ;
 346         }
 347       }
 348
 349       //Find next signal module
 350       nextSig = 200000 ;
 351       for(i = 0 ; i < fInput ; i++){
 352         sdigits = dynamic_cast<TClonesArray *>(sdigArray->At(i)) ;
 353         Int_t curNext = nextSig ;
 354         if(sdigits->GetEntriesFast() > index[i] )
 355           curNext = dynamic_cast<AliPHOSDigit *>( sdigits->At(index[i]) )->GetId() ;
 356         if(curNext < nextSig) nextSig = curNext ;
 357       }
 358
 359     }
 360   }
 361
 362   delete sdigArray ; //We should not delete its contents
 363
 364   //remove digits below thresholds
 365   for(i = 0 ; i < nEMC ; i++){
 366     digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
 367     if(sDigitizer->Calibrate( digit->GetAmp() ) < fEMCDigitThreshold)
 368       digits->RemoveAt(i) ;
 369     else
 370       digit->SetTime(gRandom->Gaus(digit->GetTime(),fTimeResolution) ) ;
 371   }
 372
 373
 374   for(i = nEMC; i < nCPV ; i++)
 375     if( sDigitizer->Calibrate( dynamic_cast<AliPHOSDigit*>(digits->At(i))->GetAmp() ) < fCPVDigitThreshold )
 376       digits->RemoveAt(i) ;
 377
 378   digits->Compress() ;
 379
 380   Int_t ndigits = digits->GetEntriesFast() ;
 381   digits->Expand(ndigits) ;
 382
 383   //Set indexes in list of digits and make true digitization of the energy
 384   for (i = 0 ; i < ndigits ; i++) {
 385     digit = dynamic_cast<AliPHOSDigit*>( digits->At(i) ) ;
 386     digit->SetIndexInList(i) ;
 387     Float_t energy = sDigitizer->Calibrate(digit->GetAmp()) ;
 388     digit->SetAmp(DigitizeEnergy(energy,digit->GetId()) ) ;
 389   }
 390 }
 391
 392 //____________________________________________________________________________
 393 Int_t AliPHOSDigitizer::DigitizeEnergy(Float_t energy, Int_t absId)
 394 {
 395   // Returns digitized value of the energy in a cell absId
 396
 397   AliPHOSGetter* gime = AliPHOSGetter::Instance();
 398
 399   if(!gime->CalibData()) {
 400     AliPHOSCalibData* cdb = new AliPHOSCalibData(gAlice->GetRunNumber());
 401     gime->SetCalibData(cdb);
 402   }
 403
 404   //Determine rel.position of the cell absId
 405   Int_t relId[4];
 406   gime->PHOSGeometry()->AbsToRelNumbering(absId,relId);
 407   Int_t module=relId[0];
 408   Int_t raw=relId[2];
 409   Int_t column=relId[3];
 410
 411   Int_t chanel ;
 412
 413   if(absId <= fEmcCrystals){ //digitize as EMC
 414
 415     //reading calibration data for cell absId.
 416     //If no calibration DB found, accept default values.
 417
 418     if(gime->CalibData()) {
 419       fADCpedestalEmc = gime->CalibData()->GetADCpedestalEmc(module,column,raw);
 420       fADCchanelEmc = gime->CalibData()->GetADCchannelEmc(module,column,raw);
 421     }
 422 //         printf("\t\tabsId %d ==>>module=%d column=%d raw=%d ped=%.4f gain=%.4f\n",
 423 //         absId,module,column,raw,fADCpedestalEmc,fADCchanelEmc);
 424
 425     chanel = (Int_t) TMath::Ceil((energy - fADCpedestalEmc)/fADCchanelEmc) ;
 426     if(chanel > fNADCemc ) chanel =  fNADCemc ;
 427   }
 428   else{ //Digitize as CPV
 429     if(gime->CalibData()) {
 430       fADCpedestalCpv = gime->CalibData()->GetADCpedestalCpv(module,column,raw);
 431       fADCchanelCpv = gime->CalibData()->GetADCchannelCpv(module,column,raw);
 432     }
 433
 434     chanel = (Int_t) TMath::Ceil((energy - fADCpedestalCpv)/fADCchanelCpv) ;
 435     if(chanel > fNADCcpv ) chanel =  fNADCcpv ;
 436   }
 437   return chanel ;
 438 }
 439
 440 //____________________________________________________________________________
 441 void AliPHOSDigitizer::Exec(Option_t *option)
 442 {
 443   // Steering method to process digitization for events
 444   // in the range from fFirstEvent to fLastEvent.
 445   // This range is optionally set by SetEventRange().
 446   // if fLastEvent=-1, then process events until the end.
 447   // by default fLastEvent = fFirstEvent (process only one event)
 448
 449   if (!fInit) { // to prevent overwrite existing file
 450     AliError(Form("Give a version name different from %s",
 451                   fEventFolderName.Data() )) ;
 452     return ;
 453   }
 454
 455   if (strstr(option,"print")) {
 456     Print();
 457     return ;
 458   }
 459
 460   if(strstr(option,"tim"))
 461     gBenchmark->Start("PHOSDigitizer");
 462
 463   AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
 464
 465   // Post Digitizer to the white board
 466   gime->PostDigitizer(this) ;
 467
 468   if (fLastEvent == -1)
 469     fLastEvent = gime->MaxEvent() - 1 ;
 470   else if (fManager)
 471     fLastEvent = fFirstEvent ;
 472
 473   Int_t nEvents   = fLastEvent - fFirstEvent + 1;
 474
 475   Int_t ievent ;
 476
 477   for (ievent = fFirstEvent; ievent <= fLastEvent; ievent++) {
 478
 479     gime->Event(ievent,"S") ;
 480
 481     Digitize(ievent) ; //Add prepared SDigits to digits and add the noise
 482
 483     WriteDigits() ;
 484
 485     if(strstr(option,"deb"))
 486       PrintDigits(option);
 487
 488     //increment the total number of Digits per run
 489     fDigitsInRun += gime->Digits()->GetEntriesFast() ;
 490  }
 491
 492   gime->PhosLoader()->CleanDigitizer();
 493
 494   if(strstr(option,"tim")){
 495     gBenchmark->Stop("PHOSDigitizer");
 496     TString message ;
 497     message = "  took %f seconds for Digitizing %f seconds per event\n" ;
 498     AliInfo(Form( message.Data(),
 499          gBenchmark->GetCpuTime("PHOSDigitizer"),
 500          gBenchmark->GetCpuTime("PHOSDigitizer")/nEvents ));
 501   }
 502 }
 503
 504 //____________________________________________________________________________
 505 Float_t AliPHOSDigitizer::FrontEdgeTime(TClonesArray * ticks) const
 506 {
 507   // Returns the shortest time among all time ticks
 508
 509   ticks->Sort() ; //Sort in accordance with times of ticks
 510   TIter it(ticks) ;
 511   AliPHOSTick * ctick = (AliPHOSTick *) it.Next() ;
 512   Float_t time = ctick->CrossingTime(fTimeThreshold) ;
 513
 514   AliPHOSTick * t ;
 515   while((t=(AliPHOSTick*) it.Next())){
 516     if(t->GetTime() < time)  //This tick starts before crossing
 517       *ctick+=*t ;
 518     else
 519       return time ;
 520
 521     time = ctick->CrossingTime(fTimeThreshold) ;
 522   }
 523   return time ;
 524 }
 525
 526 //____________________________________________________________________________
 527 Bool_t AliPHOSDigitizer::Init()
 528 {
 529   // Makes all memory allocations
 530   fInit = kTRUE ;
 531   AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
 532   if ( gime == 0 ) {
 533     AliFatal(Form("Could not obtain the Getter object for file %s and event %s !",
 534                   GetTitle(), fEventFolderName.Data()));
 535     return kFALSE;
 536   }
 537
 538   const AliPHOSGeometry * geom = gime->PHOSGeometry() ;
 539
 540   fEmcCrystals = geom->GetNModules() * geom->GetNCristalsInModule() ;
 541
 542   TString opt("Digits") ;
 543   if(gime->VersionExists(opt) ) {
 544     AliError(Form("Give a version name different from %s",
 545                   fEventFolderName.Data() )) ;
 546     fInit = kFALSE ;
 547   }
 548
 549   fFirstEvent = 0 ;
 550   fLastEvent = fFirstEvent ;
 551   if (fManager)
 552     fInput = fManager->GetNinputs() ;
 553   else
 554     fInput           = 1 ;
 555
 556   fInputFileNames  = new TString[fInput] ;
 557   fEventNames      = new TString[fInput] ;
 558   fInputFileNames[0] = GetTitle() ;
 559   fEventNames[0]     = fEventFolderName.Data() ;
 560   Int_t index ;
 561   for (index = 1 ; index < fInput ; index++) {
 562     fInputFileNames[index] = dynamic_cast<AliStream*>(fManager->GetInputStream(index))->GetFileName(0);
 563     TString tempo = fManager->GetInputFolderName(index) ;
 564     fEventNames[index] = tempo.Remove(tempo.Length()-1) ; // strip of the stream number added by fManager
 565   }
 566
 567   //to prevent cleaning of this object while GetEvent is called
 568   gime->PhosLoader()->GetDigitsDataLoader()->GetBaseTaskLoader()->SetDoNotReload(kTRUE);
 569
 570   return fInit ;
 571 }
 572
 573 //____________________________________________________________________________
 574 void AliPHOSDigitizer::InitParameters()
 575 {
 576   // Set initial parameters Digitizer
 577
 578   fPinNoise           = 0.004 ;
 579   fEMCDigitThreshold  = 0.012 ;
 580   fCPVNoise           = 0.01;
 581   fCPVDigitThreshold  = 0.09 ;
 582   fTimeResolution     = 0.5e-9 ;
 583   fTimeSignalLength   = 1.0e-9 ;
 584   fDigitsInRun  = 0 ;
 585   fADCchanelEmc = 0.0015;        // width of one ADC channel in GeV
 586   fADCpedestalEmc = 0.005 ;      //
 587   fNADCemc = (Int_t) TMath::Power(2,16) ;  // number of channels in EMC ADC
 588
 589   fADCchanelCpv = 0.0012 ;          // width of one ADC channel in CPV 'popugais'
 590   fADCpedestalCpv = 0.012 ;         //
 591   fNADCcpv = (Int_t) TMath::Power(2,12);      // number of channels in CPV ADC
 592
 593   fTimeThreshold = 0.001*10000000 ; //Means 1 MeV in terms of SDigits amplitude
 594   SetEventRange(0,-1) ;
 595
 596 }
 597
 598 //__________________________________________________________________
 599 void AliPHOSDigitizer::MixWith(TString alirunFileName, TString eventFolderName)
 600 {
 601   // Allows to produce digits by superimposing background and signal event.
 602   // It is assumed, that headers file with SIGNAL events is opened in
 603   // the constructor.
 604   // Sets the BACKGROUND event, with which the SIGNAL event is to be mixed
 605   // Thus we avoid writing (changing) huge and expensive
 606   // backgound files: all output will be writen into SIGNAL, i.e.
 607   // opened in constructor file.
 608   //
 609   // One can open as many files to mix with as one needs.
 610   // However only Sdigits with the same name (i.e. constructed with the same SDigitizer)
 611   // can be mixed.
 612
 613   if( strcmp(fEventFolderName, "") == 0 )
 614     Init() ;
 615
 616   if(fManager){
 617     Warning("MixWith", "Cannot use this method with AliRunDigitizer\n" ) ;
 618     return ;
 619   }
 620   // looking for file which contains AliRun
 621   if (gSystem->AccessPathName(alirunFileName)) {// file does not exist
 622     AliError(Form("File %s does not exist!", alirunFileName.Data())) ;
 623     return ;
 624   }
 625   // looking for the file which contains SDigits
 626   AliPHOSGetter * gime = AliPHOSGetter::Instance() ;
 627   TString fileName( gime->GetSDigitsFileName() ) ;
 628     if ( eventFolderName != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
 629       fileName = fileName.ReplaceAll(".root", "") + "_" + eventFolderName + ".root" ;
 630     if ( (gSystem->AccessPathName(fileName)) ) {
 631       AliError(Form("The file %s does not exist!", fileName.Data())) ;
 632       return ;
 633     }
 634     // need to increase the arrays
 635     TString tempo = fInputFileNames[fInput-1] ;
 636     delete [] fInputFileNames ;
 637     fInputFileNames = new TString[fInput+1] ;
 638     fInputFileNames[fInput-1] = tempo ;
 639
 640     tempo = fEventNames[fInput-1] ;
 641     delete [] fEventNames ;
 642     fEventNames = new TString[fInput+1] ;
 643     fEventNames[fInput-1] = tempo ;
 644
 645     fInputFileNames[fInput] = alirunFileName ;
 646     fEventNames[fInput]     = eventFolderName ;
 647     fInput++ ;
 648 }
 649
 650 //__________________________________________________________________
 651 void AliPHOSDigitizer::Print(const Option_t *)const
 652 {
 653   // Print Digitizer's parameters
 654   AliInfo(Form("\n------------------- %s -------------", GetName() )) ;
 655   if( strcmp(fEventFolderName.Data(), "") != 0 ){
 656     printf(" Writing Digits to branch with title  %s\n", fEventFolderName.Data()) ;
 657
 658     Int_t nStreams ;
 659     if (fManager)
 660       nStreams =  GetNInputStreams() ;
 661     else
 662       nStreams = fInput ;
 663
 664     Int_t index = 0 ;
 665     for (index = 0 ; index < nStreams ; index++) {
 666       TString tempo(fEventNames[index]) ;
 667       tempo += index ;
 668       AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[index], tempo) ;
 669       TString fileName( gime->GetSDigitsFileName() ) ;
 670       if ( fEventNames[index] != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
 671         fileName = fileName.ReplaceAll(".root", "") + "_" + fEventNames[index]  + ".root" ;
 672       printf ("Adding SDigits from %s %s\n", fInputFileNames[index].Data(), fileName.Data()) ;
 673     }
 674     AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
 675     printf("\nWriting digits to %s", gime->GetDigitsFileName().Data()) ;
 676
 677     printf("\nWith following parameters:\n") ;
 678     printf("  Electronics noise in EMC (fPinNoise) = %f\n", fPinNoise ) ;
 679     printf("  Threshold  in EMC  (fEMCDigitThreshold) = %f\n", fEMCDigitThreshold ) ;
 680     printf("  Noise in CPV (fCPVNoise) = %f\n", fCPVNoise ) ;
 681     printf("  Threshold in CPV (fCPVDigitThreshold) = %f\n",fCPVDigitThreshold ) ;
 682     printf(" ---------------------------------------------------\n") ;
 683   }
 684   else
 685     AliInfo(Form("AliPHOSDigitizer not initialized" )) ;
 686
 687 }
 688
 689 //__________________________________________________________________
 690  void AliPHOSDigitizer::PrintDigits(Option_t * option)
 691 {
 692   // Print a table of digits
 693
 694   AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
 695   TClonesArray * digits = gime->Digits() ;
 696
 697   AliInfo(Form("%d", digits->GetEntriesFast())) ;
 698   printf("\nevent %d", gAlice->GetEvNumber()) ;
 699   printf("\n       Number of entries in Digits list %d", digits->GetEntriesFast() )  ;
 700
 701
 702   if(strstr(option,"all")||strstr(option,"EMC")){
 703     //loop over digits
 704     AliPHOSDigit * digit;
 705     printf("\nEMC digits (with primaries):\n")  ;
 706     printf("\n   Id  Amplitude    Time          Index Nprim: Primaries list \n") ;
 707     Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
 708     Int_t index ;
 709     for (index = 0 ; (index < digits->GetEntriesFast()) &&
 710            (dynamic_cast<AliPHOSDigit *>(digits->At(index))->GetId() <= maxEmc) ; index++) {
 711       digit = (AliPHOSDigit * )  digits->At(index) ;
 712       if(digit->GetNprimary() == 0)
 713         continue;
 714       printf("%6d  %8d    %6.5e %4d      %2d :",
 715               digit->GetId(), digit->GetAmp(), digit->GetTime(), digit->GetIndexInList(), digit->GetNprimary()) ;
 716       Int_t iprimary;
 717       for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
 718         printf("%d ",digit->GetPrimary(iprimary+1) ) ;
 719       }
 720       printf("\n") ;
 721     }
 722   }
 723
 724   if(strstr(option,"all")||strstr(option,"CPV")){
 725
 726     //loop over CPV digits
 727     AliPHOSDigit * digit;
 728     printf("\nCPV digits (with primaries):\n")  ;
 729     printf("\n   Id  Amplitude    Time          Index Nprim: Primaries list \n") ;
 730     Int_t maxEmc = gime->PHOSGeometry()->GetNModules()*gime->PHOSGeometry()->GetNCristalsInModule() ;
 731     Int_t index ;
 732     for (index = 0 ; index < digits->GetEntriesFast(); index++) {
 733       digit = (AliPHOSDigit * )  digits->At(index) ;
 734       if(digit->GetId() > maxEmc){
 735         printf("%6d  %8d    %4d      %2d :",
 736                 digit->GetId(), digit->GetAmp(), digit->GetIndexInList(), digit->GetNprimary()) ;
 737         Int_t iprimary;
 738         for (iprimary=0; iprimary<digit->GetNprimary(); iprimary++) {
 739           printf("%d ",digit->GetPrimary(iprimary+1) ) ;
 740         }
 741         printf("\n") ;
 742       }
 743     }
 744   }
 745
 746 }
 747
 748 //__________________________________________________________________
 749 Float_t AliPHOSDigitizer::TimeOfNoise(void) const
 750 {  // Calculates the time signal generated by noise
 751   //PH  Info("TimeOfNoise", "Change me") ;
 752   return gRandom->Rndm() * 1.28E-5;
 753 }
 754
 755 //__________________________________________________________________
 756 void AliPHOSDigitizer::Unload()
 757 {
 758
 759   Int_t i ;
 760   for(i = 1 ; i < fInput ; i++){
 761     TString tempo(fEventNames[i]) ;
 762     tempo += i ;
 763     AliPHOSGetter * gime = AliPHOSGetter::Instance(fInputFileNames[i], tempo) ;
 764     gime->PhosLoader()->UnloadSDigits() ;
 765   }
 766
 767   AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle(), fEventFolderName) ;
 768   gime->PhosLoader()->UnloadDigits() ;
 769 }
 770
 771 //____________________________________________________________________________
 772 void AliPHOSDigitizer::WriteDigits()
 773 {
 774
 775   // Makes TreeD in the output file.
 776   // Check if branch already exists:
 777   //   if yes, exit without writing: ROOT TTree does not support overwriting/updating of
 778   //      already existing branches.
 779   //   else creates branch with Digits, named "PHOS", title "...",
 780   //      and branch "AliPHOSDigitizer", with the same title to keep all the parameters
 781   //      and names of files, from which digits are made.
 782
 783   AliPHOSGetter * gime = AliPHOSGetter::Instance(GetTitle()) ;
 784   const TClonesArray * digits = gime->Digits() ;
 785   TTree * treeD = gime->TreeD();
 786
 787   // -- create Digits branch
 788   Int_t bufferSize = 32000 ;
 789   TBranch * digitsBranch = treeD->Branch("PHOS","TClonesArray",&digits,bufferSize);
 790   digitsBranch->SetTitle(fEventFolderName);
 791   digitsBranch->Fill() ;
 792
 793   gime->WriteDigits("OVERWRITE");
 794   gime->WriteDigitizer("OVERWRITE");
 795
 796   Unload() ;
 797
 798 }