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