PHOS/AliPHOSGetter.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 //_________________________________________________________________________
  19 //  A singleton. This class should be used in the analysis stage to get
  20 //  reconstructed objects: Digits, RecPoints, TrackSegments and RecParticles,
  21 //  instead of directly reading them from galice.root file. This container
  22 //  ensures, that one reads Digits, made of these particular digits, RecPoints,
  23 //  made of these particular RecPoints, TrackSegments and RecParticles.
  24 //  This becomes non trivial if there are several identical branches, produced with
  25 //  different set of parameters.
  26 //
  27 //  An example of how to use (see also class AliPHOSAnalyser):
  28 //  AliPHOSGetter * gime = AliPHOSGetter::GetInstance("galice.root","test") ;
  29 //  for(Int_t irecp = 0; irecp < gime->NRecParticles() ; irecp++)
  30 //     AliPHOSRecParticle * part = gime->RecParticle(1) ;
  31 //     ................
  32 //  gime->Event(event) ;    // reads new event from galice.root
  33 //
  34 //*-- Author: Yves Schutz (SUBATECH) & Dmitri Peressounko (RRC KI & SUBATECH)
  35 //*--         Completely redesigned by Dmitri Peressounko March 2001
  36 //
  37 //*-- YS June 2001 : renamed the original AliPHOSIndexToObject and make
  38 //*--         systematic usage of TFolders without changing the interface
  39 //////////////////////////////////////////////////////////////////////////////
  40
  41 // --- ROOT system ---
  42
  43 #include <TFile.h>
  44 #include <TROOT.h>
  45 #include <TSystem.h>
  46 #include <TParticle.h>
  47 #include <TF1.h>
  48 #include <TGraph.h>
  49 //#include <TCanvas.h>
  50 //#include <TFrame.h>
  51
  52 // --- Standard library ---
  53
  54 // --- AliRoot header files ---
  55 #include "AliESD.h"
  56 #include "AliHeader.h"
  57 #include "AliMC.h"
  58 #include "AliPHOS.h"
  59 #include "AliPHOSBeamTestEvent.h"
  60 #include "AliPHOSGetter.h"
  61 #include "AliPHOSLoader.h"
  62 #include "AliRunLoader.h"
  63 #include "AliStack.h"
  64 #include "AliPHOSRawStream.h"
  65 #include "AliRawReaderFile.h"
  66 #include "AliLog.h"
  67 #include "AliCDBLocal.h"
  68 #include "AliCDBStorage.h"
  69 #include "AliCDBManager.h"
  70
  71 ClassImp(AliPHOSGetter)
  72
  73 AliPHOSGetter * AliPHOSGetter::fgObjGetter = 0 ;
  74 AliPHOSLoader * AliPHOSGetter::fgPhosLoader = 0;
  75 Int_t AliPHOSGetter::fgDebug = 0;
  76 AliPHOSCalibData* AliPHOSGetter::fCalibData = 0;
  77
  78 //  TFile * AliPHOSGetter::fgFile = 0 ;
  79
  80 //____________________________________________________________________________
  81 AliPHOSGetter::AliPHOSGetter(const char* headerFile, const char* version, Option_t * openingOption)
  82 {
  83   // ctor only called by Instance()
  84
  85   AliRunLoader* rl = AliRunLoader::GetRunLoader(version) ;
  86   if (!rl) {
  87     rl = AliRunLoader::Open(headerFile, version, openingOption);
  88     if (!rl) {
  89       Fatal("AliPHOSGetter", "Could not find the Run Loader for %s - %s",headerFile, version) ;
  90       return ;
  91     }
  92     if (rl->GetAliRun() == 0x0) {
  93       rl->LoadgAlice();
  94       gAlice = rl->GetAliRun(); // should be removed
  95       rl->LoadHeader();
  96     }
  97   }
  98   fgPhosLoader = dynamic_cast<AliPHOSLoader*>(rl->GetLoader("PHOSLoader"));
  99   if ( !fgPhosLoader )
 100     Error("AliPHOSGetter", "Could not find PHOSLoader") ;
 101   else
 102     fgPhosLoader->SetTitle(version);
 103
 104   // initialize data members
 105   SetDebug(0) ;
 106   fBTE = 0 ;
 107   fPrimaries = 0 ;
 108   fLoadingStatus = "" ;
 109
 110   fESDFileName = rl->GetFileName()  ; // this should be the galice.root file
 111   fESDFileName.ReplaceAll("galice.root", "AliESDs.root") ;
 112   fESDFile = 0 ;
 113   fESD = 0 ;
 114   fESDTree = 0 ;
 115   fRawDigits = kFALSE ;
 116
 117 }
 118
 119 //____________________________________________________________________________
 120 AliPHOSGetter::~AliPHOSGetter()
 121 {
 122   // dtor
 123   if(fgPhosLoader){
 124     delete fgPhosLoader ;
 125     fgPhosLoader = 0 ;
 126   }
 127   if(fBTE){
 128     delete fBTE ;
 129     fBTE = 0 ;
 130   }
 131   if(fPrimaries){
 132     fPrimaries->Delete() ;
 133     delete fPrimaries ;
 134   }
 135   if (fESD)
 136     delete fESD ;
 137   if (fESDTree)
 138     delete fESDTree ;
 139
 140   fgObjGetter = 0;
 141 }
 142
 143 //____________________________________________________________________________
 144 void AliPHOSGetter::Reset()
 145 {
 146   // resets things in case the getter is called consecutively with different files
 147   // the PHOS Loader is already deleted by the Run Loader
 148
 149   if (fPrimaries) {
 150     fPrimaries->Delete() ;
 151     delete fPrimaries ;
 152   }
 153   fgPhosLoader = 0;
 154   fgObjGetter = 0;
 155 }
 156
 157 //____________________________________________________________________________
 158 AliPHOSClusterizer * AliPHOSGetter::Clusterizer()
 159 {
 160   // Returns pointer to the Clusterizer task
 161   AliPHOSClusterizer * rv ;
 162   rv =  dynamic_cast<AliPHOSClusterizer *>(PhosLoader()->Reconstructioner()) ;
 163   if (!rv) {
 164     Event(0, "R") ;
 165     rv =  dynamic_cast<AliPHOSClusterizer*>(PhosLoader()->Reconstructioner()) ;
 166   }
 167   return rv ;
 168 }
 169
 170 //____________________________________________________________________________
 171 TObjArray * AliPHOSGetter::CpvRecPoints()
 172 {
 173   // asks the Loader to return the CPV RecPoints container
 174
 175   TObjArray * rv = 0 ;
 176
 177   rv = PhosLoader()->CpvRecPoints() ;
 178   if (!rv) {
 179     PhosLoader()->MakeRecPointsArray() ;
 180     rv = PhosLoader()->CpvRecPoints() ;
 181   }
 182   return rv ;
 183 }
 184
 185 //____________________________________________________________________________
 186 TClonesArray * AliPHOSGetter::Digits()
 187 {
 188   // asks the Loader to return the Digits container
 189
 190   TClonesArray * rv = 0 ;
 191   rv = PhosLoader()->Digits() ;
 192
 193   if( !rv ) {
 194     PhosLoader()->MakeDigitsArray() ;
 195     rv = PhosLoader()->Digits() ;
 196   }
 197   return rv ;
 198 }
 199
 200 //____________________________________________________________________________
 201 AliPHOSDigitizer * AliPHOSGetter::Digitizer()
 202 {
 203   // Returns pointer to the Digitizer task
 204   AliPHOSDigitizer * rv ;
 205   rv =  dynamic_cast<AliPHOSDigitizer *>(PhosLoader()->Digitizer()) ;
 206   if (!rv) {
 207     Event(0, "D") ;
 208     rv =  dynamic_cast<AliPHOSDigitizer *>(PhosLoader()->Digitizer()) ;
 209   }
 210   return rv ;
 211 }
 212
 213
 214 //____________________________________________________________________________
 215 TObjArray * AliPHOSGetter::EmcRecPoints()
 216 {
 217   // asks the Loader to return the EMC RecPoints container
 218
 219   TObjArray * rv = 0 ;
 220
 221   rv = PhosLoader()->EmcRecPoints() ;
 222   if (!rv) {
 223     PhosLoader()->MakeRecPointsArray() ;
 224     rv = PhosLoader()->EmcRecPoints() ;
 225   }
 226   return rv ;
 227 }
 228
 229 //____________________________________________________________________________
 230 TClonesArray * AliPHOSGetter::TrackSegments()
 231 {
 232   // asks the Loader to return the TrackSegments container
 233
 234   TClonesArray * rv = 0 ;
 235
 236   rv = PhosLoader()->TrackSegments() ;
 237   if (!rv) {
 238     PhosLoader()->MakeTrackSegmentsArray() ;
 239     rv = PhosLoader()->TrackSegments() ;
 240   }
 241   return rv ;
 242 }
 243
 244 //____________________________________________________________________________
 245 AliPHOSTrackSegmentMaker * AliPHOSGetter::TrackSegmentMaker()
 246 {
 247   // Returns pointer to the TrackSegmentMaker task
 248   AliPHOSTrackSegmentMaker * rv ;
 249   rv =  dynamic_cast<AliPHOSTrackSegmentMaker *>(PhosLoader()->TrackSegmentMaker()) ;
 250   if (!rv) {
 251     Event(0, "T") ;
 252     rv =  dynamic_cast<AliPHOSTrackSegmentMaker *>(PhosLoader()->TrackSegmentMaker()) ;
 253   }
 254   return rv ;
 255 }
 256
 257 //____________________________________________________________________________
 258 TClonesArray * AliPHOSGetter::RecParticles()
 259 {
 260   // asks the Loader to return the TrackSegments container
 261
 262   TClonesArray * rv = 0 ;
 263
 264   rv = PhosLoader()->RecParticles() ;
 265   if (!rv) {
 266     PhosLoader()->MakeRecParticlesArray() ;
 267     rv = PhosLoader()->RecParticles() ;
 268   }
 269   return rv ;
 270 }
 271 //____________________________________________________________________________
 272 void AliPHOSGetter::Event(Int_t event, const char* opt)
 273 {
 274   // Reads the content of all Tree's S, D and R
 275
 276 //   if ( event >= MaxEvent() ) {
 277 //     Error("Event", "%d not found in TreeE !", event) ;
 278 //     return ;
 279 //   }
 280
 281   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 282
 283 //   // checks if we are dealing with test-beam data
 284 //   TBranch * btb = rl->TreeE()->GetBranch("AliPHOSBeamTestEvent") ;
 285 //   if(btb){
 286 //     if(!fBTE)
 287 //       fBTE = new AliPHOSBeamTestEvent() ;
 288 //     btb->SetAddress(&fBTE) ;
 289 //     btb->GetEntry(event) ;
 290 //   }
 291 //   else{
 292 //     if(fBTE){
 293 //       delete fBTE ;
 294 //       fBTE = 0 ;
 295 //     }
 296 //   }
 297
 298   // Loads the type of object(s) requested
 299
 300   rl->GetEvent(event) ;
 301
 302   if(strstr(opt,"X") || (strcmp(opt,"")==0)){
 303     ReadPrimaries() ;
 304   }
 305
 306   if(strstr(opt,"H")  ){
 307     ReadTreeH();
 308   }
 309
 310   if(strstr(opt,"S")  ){
 311     ReadTreeS() ;
 312   }
 313
 314   if(strstr(opt,"D") ){
 315     ReadTreeD() ;
 316   }
 317
 318   if(strstr(opt,"R") ){
 319     ReadTreeR() ;
 320   }
 321
 322   if( strstr(opt,"T") ){
 323     ReadTreeT() ;
 324   }
 325
 326   if( strstr(opt,"P") ){
 327     ReadTreeP() ;
 328   }
 329
 330   if( strstr(opt,"E") ){
 331     ReadTreeE(event) ;
 332   }
 333
 334 }
 335
 336
 337 //____________________________________________________________________________
 338 void AliPHOSGetter::Event(AliRawReader *rawReader, const char* opt)
 339 {
 340   // Reads the raw event from rawReader
 341
 342   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 343
 344   if( strstr(opt,"W")  ){
 345     rawReader->NextEvent();
 346     Int_t iEvent = rl->GetEventNumber();
 347     rl->GetEvent(iEvent);
 348     ReadRaw(rawReader) ;
 349     rl->SetEventNumber(++iEvent);
 350   }
 351
 352 }
 353
 354
 355 //____________________________________________________________________________
 356 Int_t AliPHOSGetter::EventNumber() const
 357   {
 358   // return the current event number
 359   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 360   return static_cast<Int_t>(rl->GetEventNumber()) ;
 361 }
 362
 363 //____________________________________________________________________________
 364   TClonesArray * AliPHOSGetter::Hits()
 365 {
 366   // asks the loader to return  the Hits container
 367
 368   TClonesArray * rv = 0 ;
 369
 370   rv = PhosLoader()->Hits() ;
 371   if ( !rv ) {
 372     PhosLoader()->LoadHits("read");
 373     rv = PhosLoader()->Hits() ;
 374   }
 375   return rv ;
 376 }
 377
 378 //____________________________________________________________________________
 379 AliPHOSGetter * AliPHOSGetter::Instance(const char* alirunFileName, const char* version, Option_t * openingOption)
 380 {
 381   // Creates and returns the pointer of the unique instance
 382   // Must be called only when the environment has changed
 383
 384   if(!fgObjGetter){ // first time the getter is called
 385     fgObjGetter = new AliPHOSGetter(alirunFileName, version, openingOption) ;
 386   }
 387   else { // the getter has been called previously
 388     AliRunLoader * rl = AliRunLoader::GetRunLoader(fgPhosLoader->GetTitle());
 389     if ( rl->GetFileName() == alirunFileName ) {// the alirunFile has the same name
 390       // check if the file is already open
 391       TFile * galiceFile = dynamic_cast<TFile *>(gROOT->FindObject(rl->GetFileName()) ) ;
 392
 393       if ( !galiceFile )
 394         fgObjGetter = new AliPHOSGetter(alirunFileName, version, openingOption) ;
 395
 396       else {  // the file is already open check the version name
 397         TString currentVersionName = rl->GetEventFolder()->GetName() ;
 398         TString newVersionName(version) ;
 399         if (currentVersionName == newVersionName)
 400           if(fgDebug)
 401             ::Warning( "Instance", "Files with version %s already open", currentVersionName.Data() ) ;
 402         else {
 403           fgObjGetter = new AliPHOSGetter(alirunFileName, version, openingOption) ;
 404         }
 405       }
 406     }
 407     else {
 408       AliRunLoader * rl = AliRunLoader::GetRunLoader(fgPhosLoader->GetTitle()) ;
 409       if ( strstr(version, AliConfig::GetDefaultEventFolderName()) ) // false in case of merging
 410         delete rl ;
 411       fgObjGetter = new AliPHOSGetter(alirunFileName, version, openingOption) ;
 412     }
 413   }
 414   if (!fgObjGetter)
 415     ::Error("AliPHOSGetter::Instance", "Failed to create the PHOS Getter object") ;
 416   else
 417     if (fgDebug)
 418       Print() ;
 419
 420   return fgObjGetter ;
 421 }
 422
 423 //____________________________________________________________________________
 424 AliPHOSGetter *  AliPHOSGetter::Instance()
 425 {
 426   // Returns the pointer of the unique instance already defined
 427
 428   if(!fgObjGetter && fgDebug)
 429      ::Warning("AliPHOSGetter::Instance", "Getter not initialized") ;
 430
 431    return fgObjGetter ;
 432
 433 }
 434
 435 //____________________________________________________________________________
 436 Int_t AliPHOSGetter::MaxEvent() const
 437 {
 438   // returns the number of events in the run (from TE)
 439
 440   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 441   return static_cast<Int_t>(rl->GetNumberOfEvents()) ;
 442 }
 443
 444 //____________________________________________________________________________
 445 TParticle * AliPHOSGetter::Primary(Int_t index) const
 446 {
 447   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 448   return rl->Stack()->Particle(index) ;
 449 }
 450
 451 //____________________________________________________________________________
 452 AliPHOS * AliPHOSGetter:: PHOS() const
 453 {
 454   // returns the PHOS object
 455   AliPHOSLoader *    loader = 0;
 456   static AliPHOSLoader * oldloader = 0;
 457   static AliPHOS * phos = 0;
 458
 459   loader = PhosLoader();
 460
 461   if ( loader != oldloader) {
 462     phos = dynamic_cast<AliPHOS*>(loader->GetModulesFolder()->FindObject("PHOS")) ;
 463     oldloader = loader;
 464   }
 465   if (!phos)
 466     if (fgDebug)
 467       Warning("PHOS", "PHOS module not found in module folders: %s", PhosLoader()->GetModulesFolder()->GetName() ) ;
 468   return phos ;
 469 }
 470
 471
 472
 473 //____________________________________________________________________________
 474 AliPHOSPID * AliPHOSGetter::PID()
 475 {
 476   // Returns pointer to the PID task
 477   AliPHOSPID * rv ;
 478   rv =  dynamic_cast<AliPHOSPID *>(PhosLoader()->PIDTask()) ;
 479   if (!rv) {
 480     Event(0, "P") ;
 481     rv =  dynamic_cast<AliPHOSPID *>(PhosLoader()->PIDTask()) ;
 482   }
 483   return rv ;
 484 }
 485
 486 //____________________________________________________________________________
 487 AliPHOSGeometry * AliPHOSGetter::PHOSGeometry() const
 488 {
 489   // Returns PHOS geometry
 490
 491   AliPHOSGeometry * rv = 0 ;
 492   if (PHOS() )
 493     rv =  PHOS()->GetGeometry() ;
 494   return rv ;
 495 }
 496
 497 //____________________________________________________________________________
 498 TClonesArray * AliPHOSGetter::Primaries()
 499 {
 500   // creates the Primaries container if needed
 501   if ( !fPrimaries ) {
 502     if (fgDebug)
 503       Info("Primaries", "Creating a new TClonesArray for primaries") ;
 504     fPrimaries = new TClonesArray("TParticle", 1000) ;
 505   }
 506   return fPrimaries ;
 507 }
 508
 509 //____________________________________________________________________________
 510 void  AliPHOSGetter::Print()
 511 {
 512   // Print usefull information about the getter
 513
 514   AliRunLoader * rl = AliRunLoader::GetRunLoader(fgPhosLoader->GetTitle());
 515   ::Info( "Print", "gAlice file is %s -- version name is %s", (rl->GetFileName()).Data(), rl->GetEventFolder()->GetName() ) ;
 516 }
 517
 518 //____________________________________________________________________________
 519 void AliPHOSGetter::ReadPrimaries()
 520 {
 521   // Read Primaries from Kinematics.root
 522
 523   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 524
 525   // gets kine tree from the root file (Kinematics.root)
 526   if ( ! rl->TreeK() ) { // load treeK the first time
 527     rl->LoadKinematics() ;
 528   }
 529
 530   fNPrimaries = (rl->GetHeader())->GetNtrack();
 531   if (fgDebug)
 532     Info( "ReadTreeK", "Found %d particles in event # %d", fNPrimaries, EventNumber() ) ;
 533
 534
 535   // first time creates the container
 536   if ( Primaries() )
 537     fPrimaries->Clear() ;
 538
 539   Int_t index = 0 ;
 540   for (index = 0 ; index < fNPrimaries; index++) {
 541     new ((*fPrimaries)[index]) TParticle(*(Primary(index)));
 542   }
 543 }
 544
 545 //____________________________________________________________________________
 546 Bool_t AliPHOSGetter::OpenESDFile()
 547 {
 548   //Open the ESD file
 549   Bool_t rv = kTRUE ;
 550   if (!fESDFile) {
 551     fESDFile = TFile::Open(fESDFileName) ;
 552     if (!fESDFile )
 553       return kFALSE ;
 554   }
 555   else if (fESDFile->IsOpen()) {
 556     fESDFile->Close() ;
 557     fESDFile = TFile::Open(fESDFileName) ;
 558   }
 559   if (!fESDFile->IsOpen())
 560     rv = kFALSE ;
 561   return rv ;
 562 }
 563
 564 //____________________________________________________________________________
 565 void AliPHOSGetter::FitRaw(Bool_t lowGainFlag, TGraph * gLowGain, TGraph * gHighGain, TF1* signalF, Double_t & energy, Double_t & time)
 566 {
 567   // Fits the raw signal time distribution
 568
 569   const Int_t kNoiseThreshold = 0 ;
 570   Double_t timezero1 = 0., timezero2 = 0., timemax = 0. ;
 571   Double_t signal = 0., signalmax = 0. ;
 572   time   = 0. ;
 573   energy = 0. ;
 574
 575   if (lowGainFlag) {
 576     timezero1 = timezero2 = signalmax = timemax = 0. ;
 577     signalF->FixParameter(0, PHOS()->GetRawFormatLowCharge()) ;
 578     signalF->FixParameter(1, PHOS()->GetRawFormatLowGain()) ;
 579     Int_t index ;
 580     for (index = 0; index < PHOS()->GetRawFormatTimeBins(); index++) {
 581       gLowGain->GetPoint(index, time, signal) ;
 582       if (signal > kNoiseThreshold && timezero1 == 0.)
 583         timezero1 = time ;
 584       if (signal <= kNoiseThreshold && timezero1 > 0. && timezero2 == 0.)
 585         timezero2 = time ;
 586       if (signal > signalmax) {
 587         signalmax = signal ;
 588         timemax   = time ;
 589       }
 590     }
 591     signalmax /= PHOS()->RawResponseFunctionMax(PHOS()->GetRawFormatLowCharge(),
 592                                                 PHOS()->GetRawFormatLowGain()) ;
 593     if ( timezero1 + PHOS()->GetRawFormatTimePeak() < PHOS()->GetRawFormatTimeMax() * 0.4 ) { // else its noise
 594       signalF->SetParameter(2, signalmax) ;
 595       signalF->SetParameter(3, timezero1) ;
 596       gLowGain->Fit(signalF, "QRON", "", 0., timezero2); //, "QRON") ;
 597       energy = signalF->GetParameter(2) ;
 598       time   = signalF->GetMaximumX() - PHOS()->GetRawFormatTimePeak() - PHOS()->GetRawFormatTimeTrigger() ;
 599     }
 600   } else {
 601     timezero1 = timezero2 = signalmax = timemax = 0. ;
 602     signalF->FixParameter(0, PHOS()->GetRawFormatHighCharge()) ;
 603     signalF->FixParameter(1, PHOS()->GetRawFormatHighGain()) ;
 604     Int_t index ;
 605     for (index = 0; index < PHOS()->GetRawFormatTimeBins(); index++) {
 606       gHighGain->GetPoint(index, time, signal) ;
 607       if (signal > kNoiseThreshold && timezero1 == 0.)
 608         timezero1 = time ;
 609       if (signal <= kNoiseThreshold && timezero1 > 0. && timezero2 == 0.)
 610         timezero2 = time ;
 611       if (signal > signalmax) {
 612         signalmax = signal ;
 613         timemax   = time ;
 614       }
 615     }
 616     signalmax /= PHOS()->RawResponseFunctionMax(PHOS()->GetRawFormatHighCharge(),
 617                                                 PHOS()->GetRawFormatHighGain()) ;;
 618     if ( timezero1 + PHOS()->GetRawFormatTimePeak() < PHOS()->GetRawFormatTimeMax() * 0.4 ) { // else its noise
 619       signalF->SetParameter(2, signalmax) ;
 620       signalF->SetParameter(3, timezero1) ;
 621       gHighGain->Fit(signalF, "QRON", "", 0., timezero2) ;
 622       energy = signalF->GetParameter(2) ;
 623       time   = signalF->GetMaximumX() - PHOS()->GetRawFormatTimePeak() - PHOS()->GetRawFormatTimeTrigger() ;
 624     }
 625   }
 626   if (time == 0) energy = 0 ;
 627 }
 628
 629 //____________________________________________________________________________
 630 Int_t AliPHOSGetter::CalibrateRaw(Double_t energy, Int_t *relId)
 631 {
 632   // Convert energy into digitized amplitude for a cell relId
 633   // It is a user responsilibity to open CDB and set
 634   // AliPHOSCalibData object by the following operators:
 635   //
 636   // AliCDBLocal *loc = new AliCDBLocal("deCalibDB");
 637   // AliPHOSCalibData* clb = (AliPHOSCalibData*)AliCDBStorage::Instance()
 638   //    ->Get(path_to_calibdata,run_number);
 639   // AliPHOSGetter* gime = AliPHOSGetter::Instance("galice.root");
 640   // gime->SetCalibData(clb);
 641
 642   if (CalibData() == 0)
 643     Warning("CalibrateRaw","Calibration DB is not initiated!");
 644
 645   Int_t   module = relId[0];
 646   Int_t   column = relId[3];
 647   Int_t   row    = relId[2];
 648
 649   Float_t gainFactor = 0.0015; // width of one Emc ADC channel in GeV
 650   Float_t pedestal   = 0.005;  // Emc pedestals
 651
 652   if(CalibData()) {
 653     gainFactor = CalibData()->GetADCchannelEmc (module,column,row);
 654     pedestal   = CalibData()->GetADCpedestalEmc(module,column,row);
 655   }
 656
 657   Int_t   amp = static_cast<Int_t>( (energy - pedestal) / gainFactor + 0.5 ) ;
 658   return amp;
 659 }
 660 //____________________________________________________________________________
 661 Int_t AliPHOSGetter::ReadRaw(AliRawReader *rawReader)
 662 {
 663   // reads the raw format data, converts it into digits format and store digits in Digits()
 664   // container.
 665
 666   AliPHOSRawStream in(rawReader);
 667
 668   Bool_t first = kTRUE ;
 669
 670   TF1 * signalF = new TF1("signal", AliPHOS::RawResponseFunction, 0, PHOS()->GetRawFormatTimeMax(), 4);
 671   signalF->SetParNames("Charge", "Gain", "Amplitude", "TimeZero") ;
 672
 673   Int_t relId[4], id =0;
 674   Bool_t lowGainFlag = kFALSE ;
 675
 676   TClonesArray * digits = Digits() ;
 677   digits->Clear() ;
 678   Int_t idigit = 0 ;
 679   Int_t amp    = 0 ;
 680   Double_t energy = 0. ;
 681   Double_t time   = 0. ;
 682
 683   TGraph * gLowGain = new TGraph(PHOS()->GetRawFormatTimeBins()) ;
 684   TGraph * gHighGain= new TGraph(PHOS()->GetRawFormatTimeBins()) ;
 685
 686   while ( in.Next() ) { // PHOS entries loop
 687     if ( (in.IsNewRow() || in.IsNewColumn() || in.IsNewModule()) ) {
 688       relId[0] = in.GetModule() ;
 689       if ( relId[0] >= PHOS()->GetRawFormatLowGainOffset() ) {
 690         relId[0] -=  PHOS()->GetRawFormatLowGainOffset() ;
 691         lowGainFlag = kTRUE ;
 692       } else
 693         lowGainFlag = kFALSE ;
 694       relId[1] = 0 ;
 695       relId[2] = in.GetRow() ;
 696       relId[3] = in.GetColumn() ;
 697       PHOSGeometry()->RelToAbsNumbering(relId, id) ;
 698       if (!first) {
 699         FitRaw(lowGainFlag, gLowGain, gHighGain, signalF, energy, time) ;
 700         amp = CalibrateRaw(energy,relId);
 701         if (amp > 0) {
 702           new((*digits)[idigit]) AliPHOSDigit( -1, id, amp, time) ;
 703           idigit++ ;
 704         }
 705         Int_t index ;
 706         for (index = 0; index < PHOS()->GetRawFormatTimeBins(); index++) {
 707           gLowGain->SetPoint(index, index * PHOS()->GetRawFormatTimeMax() / PHOS()->GetRawFormatTimeBins(), 0) ;
 708           gHighGain->SetPoint(index, index * PHOS()->GetRawFormatTimeMax() / PHOS()->GetRawFormatTimeBins(), 0) ;
 709         }
 710       }
 711       first = kFALSE ;
 712     }
 713     if (lowGainFlag)
 714       gLowGain->SetPoint(in.GetTime(),
 715                      in.GetTime()* PHOS()->GetRawFormatTimeMax() / PHOS()->GetRawFormatTimeBins(),
 716                      in.GetSignal()) ;
 717     else
 718       gHighGain->SetPoint(in.GetTime(),
 719                          in.GetTime() * PHOS()->GetRawFormatTimeMax() / PHOS()->GetRawFormatTimeBins(),
 720                          in.GetSignal() ) ;
 721
 722   } // PHOS entries loop
 723
 724   FitRaw(lowGainFlag, gLowGain, gHighGain, signalF, energy, time) ;
 725   amp = CalibrateRaw(energy,relId);
 726   if (amp > 0 ) {
 727     new((*digits)[idigit]) AliPHOSDigit( -1, id, amp, time) ;
 728     idigit++ ;
 729   }
 730   digits->Sort() ;
 731
 732   delete signalF ;
 733   delete gLowGain;
 734   delete gHighGain ;
 735
 736 //   AliInfo(Form("Event %d, digits: %d\n", EventNumber(),digits->GetEntries())); //bvp
 737   return digits->GetEntriesFast() ;
 738 }
 739
 740 //   TClonesArray * digits = Digits() ;
 741 //   digits->Clear() ;
 742 //   Int_t idigit = 0 ;
 743
 744 //   while ( in.Next() ) { // PHOS entries loop
 745
 746 //     Int_t amp = in.GetSignal() ;
 747 //     Double_t time = in.GetTime() ;
 748 //     Int_t relId[4], id ;
 749
 750 //     relId[0] = in.GetModule() ;
 751 //     if ( relId[0] >= PHOS()->GetRawFormatLowGainOffset() ) {
 752 //       relId[0] -=  PHOS()->GetRawFormatLowGainOffset() ;
 753 //     }
 754 //     relId[1] = 0 ;
 755 //     relId[2] = in.GetRow() ;
 756 //     relId[3] = in.GetColumn() ;
 757 //     PHOSGeometry()->RelToAbsNumbering(relId, id) ;
 758
 759 //     if (amp > 0 && id >=0 ) {
 760 //       new((*digits)[idigit]) AliPHOSDigit( -1, id, amp, time) ;
 761 //       idigit++ ;
 762 //     }
 763
 764 //   } // PHOS entries loop
 765
 766 //   digits->Sort() ;
 767
 768 //   return digits->GetEntriesFast() ;
 769 // }
 770 //____________________________________________________________________________
 771 Int_t AliPHOSGetter::ReadTreeD()
 772 {
 773   // Read the Digits
 774
 775   PhosLoader()->CleanDigits() ;
 776   PhosLoader()->LoadDigits("UPDATE") ;
 777   PhosLoader()->LoadDigitizer("UPDATE") ;
 778   return Digits()->GetEntries() ;
 779 }
 780
 781 //____________________________________________________________________________
 782 Int_t AliPHOSGetter::ReadTreeH()
 783 {
 784   // Read the Hits
 785   PhosLoader()->CleanHits() ;
 786   // gets TreeH from the root file (PHOS.Hit.root)
 787   //if ( !IsLoaded("H") ) {
 788     PhosLoader()->LoadHits("UPDATE") ;
 789   //  SetLoaded("H") ;
 790   //}
 791   return Hits()->GetEntries() ;
 792 }
 793
 794 //____________________________________________________________________________
 795 Int_t AliPHOSGetter::ReadTreeR()
 796 {
 797   // Read the RecPoints
 798
 799   PhosLoader()->CleanRecPoints() ;
 800   // gets TreeR from the root file (PHOS.RecPoints.root)
 801   //if ( !IsLoaded("R") ) {
 802     PhosLoader()->LoadRecPoints("UPDATE") ;
 803     PhosLoader()->LoadClusterizer("UPDATE") ;
 804     //  SetLoaded("R") ;
 805     //}
 806
 807   return EmcRecPoints()->GetEntries() ;
 808 }
 809
 810 //____________________________________________________________________________
 811 Int_t AliPHOSGetter::ReadTreeT()
 812 {
 813   // Read the TrackSegments
 814
 815   PhosLoader()->CleanTracks() ;
 816   // gets TreeT from the root file (PHOS.TrackSegments.root)
 817   //if ( !IsLoaded("T") ) {
 818     PhosLoader()->LoadTracks("UPDATE") ;
 819     PhosLoader()->LoadTrackSegmentMaker("UPDATE") ;
 820     //    SetLoaded("T") ;
 821     //}
 822
 823   return TrackSegments()->GetEntries() ;
 824 }
 825 //____________________________________________________________________________
 826 Int_t AliPHOSGetter::ReadTreeP()
 827 {
 828   // Read the RecParticles
 829
 830   PhosLoader()->CleanRecParticles() ;
 831
 832   // gets TreeT from the root file (PHOS.TrackSegments.root)
 833   //  if ( !IsLoaded("P") ) {
 834     PhosLoader()->LoadRecParticles("UPDATE") ;
 835     PhosLoader()->LoadPID("UPDATE") ;
 836     //  SetLoaded("P") ;
 837     //}
 838
 839   return RecParticles()->GetEntries() ;
 840 }
 841 //____________________________________________________________________________
 842 Int_t AliPHOSGetter::ReadTreeS()
 843 {
 844   // Read the SDigits
 845
 846   PhosLoader()->CleanSDigits() ;
 847   // gets TreeS from the root file (PHOS.SDigits.root)
 848   //if ( !IsLoaded("S") ) {
 849     PhosLoader()->LoadSDigits("READ") ;
 850     PhosLoader()->LoadSDigitizer("READ") ;
 851     //  SetLoaded("S") ;
 852     //}
 853
 854   return SDigits()->GetEntries() ;
 855 }
 856
 857 //____________________________________________________________________________
 858 Int_t AliPHOSGetter::ReadTreeE(Int_t event)
 859 {
 860   // Read the ESD
 861
 862   // gets esdTree from the root file (AliESDs.root)
 863   if (!fESDFile)
 864     if ( !OpenESDFile() )
 865       return -1 ;
 866
 867   fESDTree = static_cast<TTree*>(fESDFile->Get("esdTree")) ;
 868   fESD = new AliESD;
 869    if (!fESDTree) {
 870
 871      Error("ReadTreeE", "no ESD tree found");
 872      return -1;
 873    }
 874    fESDTree->SetBranchAddress("ESD", &fESD);
 875    fESDTree->GetEvent(event);
 876
 877    return event ;
 878 }
 879
 880 //____________________________________________________________________________
 881 TClonesArray * AliPHOSGetter::SDigits()
 882 {
 883   // asks the Loader to return the Digits container
 884
 885   TClonesArray * rv = 0 ;
 886
 887   rv = PhosLoader()->SDigits() ;
 888   if (!rv) {
 889     PhosLoader()->MakeSDigitsArray() ;
 890     rv = PhosLoader()->SDigits() ;
 891   }
 892   return rv ;
 893 }
 894
 895 //____________________________________________________________________________
 896 AliPHOSSDigitizer * AliPHOSGetter::SDigitizer()
 897 {
 898   // Returns pointer to the SDigitizer task
 899   AliPHOSSDigitizer * rv ;
 900   rv =  dynamic_cast<AliPHOSSDigitizer *>(PhosLoader()->SDigitizer()) ;
 901   if (!rv) {
 902     Event(0, "S") ;
 903     rv =  dynamic_cast<AliPHOSSDigitizer *>(PhosLoader()->SDigitizer()) ;
 904   }
 905   return rv ;
 906 }
 907
 908 //____________________________________________________________________________
 909 TParticle * AliPHOSGetter::Secondary(const TParticle* p, Int_t index) const
 910 {
 911   // Return first (index=1) or second (index=2) secondary particle of primary particle p
 912
 913   if(index <= 0)
 914     return 0 ;
 915   if(index > 2)
 916     return 0 ;
 917
 918   if(p) {
 919   Int_t daughterIndex = p->GetDaughter(index-1) ;
 920   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 921   return  rl->GetAliRun()->GetMCApp()->Particle(daughterIndex) ;
 922   }
 923   else
 924     return 0 ;
 925 }
 926
 927 //____________________________________________________________________________
 928 void AliPHOSGetter::Track(Int_t itrack)
 929 {
 930   // Read the first entry of PHOS branch in hit tree gAlice->TreeH()
 931
 932  AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 933
 934   if( !TreeH() ) // load treeH the first time
 935     rl->LoadHits() ;
 936
 937   // first time create the container
 938   TClonesArray * hits = Hits() ;
 939   if ( hits )
 940     hits->Clear() ;
 941
 942   TBranch * phosbranch = dynamic_cast<TBranch*>(TreeH()->GetBranch("PHOS")) ;
 943   phosbranch->SetAddress(&hits) ;
 944   phosbranch->GetEntry(itrack) ;
 945 }
 946
 947 //____________________________________________________________________________
 948 TTree * AliPHOSGetter::TreeD() const
 949 {
 950   // Returns pointer to the Digits Tree
 951   TTree * rv = 0 ;
 952   rv = PhosLoader()->TreeD() ;
 953   if ( !rv ) {
 954     PhosLoader()->MakeTree("D");
 955     rv = PhosLoader()->TreeD() ;
 956   }
 957
 958   return rv ;
 959 }
 960
 961 //____________________________________________________________________________
 962 TTree * AliPHOSGetter::TreeH() const
 963 {
 964   // Returns pointer to the Hits Tree
 965   TTree * rv = 0 ;
 966   rv = PhosLoader()->TreeH() ;
 967   if ( !rv ) {
 968     PhosLoader()->MakeTree("H");
 969     rv = PhosLoader()->TreeH() ;
 970   }
 971
 972   return rv ;
 973 }
 974
 975 //____________________________________________________________________________
 976 TTree * AliPHOSGetter::TreeR() const
 977 {
 978   // Returns pointer to the RecPoints Tree
 979   TTree * rv = 0 ;
 980   rv = PhosLoader()->TreeR() ;
 981   if ( !rv ) {
 982     PhosLoader()->MakeTree("R");
 983     rv = PhosLoader()->TreeR() ;
 984   }
 985
 986   return rv ;
 987 }
 988
 989 //____________________________________________________________________________
 990 TTree * AliPHOSGetter::TreeT() const
 991 {
 992   // Returns pointer to the TrackSegments Tree
 993   TTree * rv = 0 ;
 994   rv = PhosLoader()->TreeT() ;
 995   if ( !rv ) {
 996     PhosLoader()->MakeTree("T");
 997     rv = PhosLoader()->TreeT() ;
 998   }
 999
1000   return rv ;
1001 }
1002 //____________________________________________________________________________
1003 TTree * AliPHOSGetter::TreeP() const
1004 {
1005   // Returns pointer to the RecParticles  Tree
1006   TTree * rv = 0 ;
1007   rv = PhosLoader()->TreeP() ;
1008   if ( !rv ) {
1009     PhosLoader()->MakeTree("P");
1010     rv = PhosLoader()->TreeP() ;
1011   }
1012
1013   return rv ;
1014 }
1015
1016 //____________________________________________________________________________
1017 TTree * AliPHOSGetter::TreeS() const
1018 {
1019  // Returns pointer to the SDigits Tree
1020   TTree * rv = 0 ;
1021   rv = PhosLoader()->TreeS() ;
1022   if ( !rv ) {
1023     PhosLoader()->MakeTree("S");
1024     rv = PhosLoader()->TreeS() ;
1025   }
1026
1027   return rv ;
1028 }
1029
1030 //____________________________________________________________________________
1031 Bool_t AliPHOSGetter::VersionExists(TString & opt) const
1032 {
1033   // checks if the version with the present name already exists in the same directory
1034
1035   Bool_t rv = kFALSE ;
1036
1037   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
1038   TString version( rl->GetEventFolder()->GetName() ) ;
1039
1040   opt.ToLower() ;
1041
1042   if ( opt == "sdigits") {
1043     // add the version name to the root file name
1044     TString fileName( PhosLoader()->GetSDigitsFileName() ) ;
1045     if (version != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
1046       fileName = fileName.ReplaceAll(".root", "") + "_" + version + ".root" ;
1047     if ( !(gSystem->AccessPathName(fileName)) ) {
1048       Warning("VersionExists", "The file %s already exists", fileName.Data()) ;
1049       rv = kTRUE ;
1050     }
1051     PhosLoader()->SetSDigitsFileName(fileName) ;
1052   }
1053
1054   if ( opt == "digits") {
1055     // add the version name to the root file name
1056     TString fileName( PhosLoader()->GetDigitsFileName() ) ;
1057     if (version != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
1058       fileName = fileName.ReplaceAll(".root", "") + "_" + version + ".root" ;
1059     if ( !(gSystem->AccessPathName(fileName)) ) {
1060       Warning("VersionExists", "The file %s already exists", fileName.Data()) ;
1061       rv = kTRUE ;
1062     }
1063   }
1064
1065   return rv ;
1066
1067 }
1068
1069 //____________________________________________________________________________
1070 UShort_t AliPHOSGetter::EventPattern(void) const
1071 {
1072   // Return the pattern (trigger bit register) of the beam-test event
1073   if(fBTE)
1074     return fBTE->GetPattern() ;
1075   else
1076     return 0 ;
1077 }
1078 //____________________________________________________________________________
1079 Float_t AliPHOSGetter::BeamEnergy(void) const
1080 {
1081   // Return the beam energy of the beam-test event
1082   if(fBTE)
1083     return fBTE->GetBeamEnergy() ;
1084   else
1085     return 0 ;
1086 }
1087 //____________________________________________________________________________
1088
1089 AliPHOSCalibData* AliPHOSGetter::CalibData()
1090 {
1091
1092   if( !(AliCDBManager::Instance()->IsDefaultStorageSet()) ) {
1093     fCalibData=0x0;
1094     Warning("CalibData","Calibration DB is not initiated!");
1095     return fCalibData;
1096   }
1097
1098   return fCalibData;
1099
1100
1101 }