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 "AliRunLoader.h"
  64 #include "AliStack.h"
  65 #include "AliPHOSRawStream.h"
  66 #include "AliRawReaderFile.h"
  67 #include "AliLog.h"
  68 #include "AliCDBLocal.h"
  69 #include "AliCDBStorage.h"
  70 #include "AliCDBManager.h"
  71
  72 ClassImp(AliPHOSGetter)
  73
  74 AliPHOSGetter   * AliPHOSGetter::fgObjGetter  = 0 ;
  75 AliPHOSLoader   * AliPHOSGetter::fgPhosLoader = 0;
  76 AliPHOSCalibData* AliPHOSGetter::fgCalibData  = 0;
  77 Int_t AliPHOSGetter::fgDebug = 0;
  78
  79 //  TFile * AliPHOSGetter::fgFile = 0 ;
  80
  81 //____________________________________________________________________________
  82 AliPHOSGetter::AliPHOSGetter(const char* headerFile, const char* version, Option_t * openingOption)
  83 {
  84   // ctor only called by Instance()
  85
  86   AliRunLoader* rl = AliRunLoader::GetRunLoader(version) ;
  87   if (!rl) {
  88     rl = AliRunLoader::Open(headerFile, version, openingOption);
  89     if (!rl) {
  90       Fatal("AliPHOSGetter", "Could not find the Run Loader for %s - %s",headerFile, version) ;
  91       return ;
  92     }
  93     if (rl->GetAliRun() == 0x0) {
  94       rl->LoadgAlice();
  95       gAlice = rl->GetAliRun(); // should be removed
  96       rl->LoadHeader();
  97     }
  98   }
  99   fgPhosLoader = dynamic_cast<AliPHOSLoader*>(rl->GetLoader("PHOSLoader"));
 100   if ( !fgPhosLoader )
 101     Error("AliPHOSGetter", "Could not find PHOSLoader") ;
 102   else
 103     fgPhosLoader->SetTitle(version);
 104
 105   // initialize data members
 106   SetDebug(0) ;
 107   fBTE = 0 ;
 108   fPrimaries = 0 ;
 109   fLoadingStatus = "" ;
 110
 111   fESDFileName = rl->GetFileName()  ; // this should be the galice.root file
 112   fESDFileName.ReplaceAll("galice.root", "AliESDs.root") ;
 113   fESDFile = 0 ;
 114   fESD = 0 ;
 115   fESDTree = 0 ;
 116   fRawDigits = kFALSE ;
 117
 118 }
 119
 120 //____________________________________________________________________________
 121 AliPHOSGetter::~AliPHOSGetter()
 122 {
 123   // dtor
 124   if(fgPhosLoader){
 125     delete fgPhosLoader ;
 126     fgPhosLoader = 0 ;
 127   }
 128   if(fBTE){
 129     delete fBTE ;
 130     fBTE = 0 ;
 131   }
 132   if(fPrimaries){
 133     fPrimaries->Delete() ;
 134     delete fPrimaries ;
 135   }
 136   if (fESD)
 137     delete fESD ;
 138   if (fESDTree)
 139     delete fESDTree ;
 140
 141   fgObjGetter = 0;
 142 }
 143
 144 //____________________________________________________________________________
 145 void AliPHOSGetter::Reset()
 146 {
 147   // resets things in case the getter is called consecutively with different files
 148   // the PHOS Loader is already deleted by the Run Loader
 149
 150   if (fPrimaries) {
 151     fPrimaries->Delete() ;
 152     delete fPrimaries ;
 153   }
 154   fgPhosLoader = 0;
 155   fgObjGetter = 0;
 156 }
 157
 158 //____________________________________________________________________________
 159 AliPHOSClusterizer * AliPHOSGetter::Clusterizer()
 160 {
 161   // Returns pointer to the Clusterizer task
 162   AliPHOSClusterizer * rv ;
 163   rv =  dynamic_cast<AliPHOSClusterizer *>(PhosLoader()->Reconstructioner()) ;
 164   if (!rv) {
 165     Event(0, "R") ;
 166     rv =  dynamic_cast<AliPHOSClusterizer*>(PhosLoader()->Reconstructioner()) ;
 167   }
 168   return rv ;
 169 }
 170
 171 //____________________________________________________________________________
 172 TObjArray * AliPHOSGetter::CpvRecPoints() const
 173 {
 174   // asks the Loader to return the CPV RecPoints container
 175
 176   TObjArray * rv = 0 ;
 177
 178   rv = PhosLoader()->CpvRecPoints() ;
 179   if (!rv) {
 180     PhosLoader()->MakeRecPointsArray() ;
 181     rv = PhosLoader()->CpvRecPoints() ;
 182   }
 183   return rv ;
 184 }
 185
 186 //____________________________________________________________________________
 187 TClonesArray * AliPHOSGetter::Digits() const
 188 {
 189   // asks the Loader to return the Digits container
 190
 191   TClonesArray * rv = 0 ;
 192   rv = PhosLoader()->Digits() ;
 193
 194   if( !rv ) {
 195     PhosLoader()->MakeDigitsArray() ;
 196     rv = PhosLoader()->Digits() ;
 197   }
 198   return rv ;
 199 }
 200
 201 //____________________________________________________________________________
 202 AliPHOSDigitizer * AliPHOSGetter::Digitizer()
 203 {
 204   // Returns pointer to the Digitizer task
 205   AliPHOSDigitizer * rv ;
 206   rv =  dynamic_cast<AliPHOSDigitizer *>(PhosLoader()->Digitizer()) ;
 207   if (!rv) {
 208     Event(0, "D") ;
 209     rv =  dynamic_cast<AliPHOSDigitizer *>(PhosLoader()->Digitizer()) ;
 210   }
 211   return rv ;
 212 }
 213
 214
 215 //____________________________________________________________________________
 216 TObjArray * AliPHOSGetter::EmcRecPoints() const
 217 {
 218   // asks the Loader to return the EMC RecPoints container
 219
 220   TObjArray * rv = 0 ;
 221
 222   rv = PhosLoader()->EmcRecPoints() ;
 223   if (!rv) {
 224     PhosLoader()->MakeRecPointsArray() ;
 225     rv = PhosLoader()->EmcRecPoints() ;
 226   }
 227   return rv ;
 228 }
 229
 230 //____________________________________________________________________________
 231 TClonesArray * AliPHOSGetter::TrackSegments() const
 232 {
 233   // asks the Loader to return the TrackSegments container
 234
 235   TClonesArray * rv = 0 ;
 236
 237   rv = PhosLoader()->TrackSegments() ;
 238   if (!rv) {
 239     PhosLoader()->MakeTrackSegmentsArray() ;
 240     rv = PhosLoader()->TrackSegments() ;
 241   }
 242   return rv ;
 243 }
 244
 245 //____________________________________________________________________________
 246 AliPHOSTrackSegmentMaker * AliPHOSGetter::TrackSegmentMaker()
 247 {
 248   // Returns pointer to the TrackSegmentMaker task
 249   AliPHOSTrackSegmentMaker * rv ;
 250   rv =  dynamic_cast<AliPHOSTrackSegmentMaker *>(PhosLoader()->TrackSegmentMaker()) ;
 251   if (!rv) {
 252     Event(0, "T") ;
 253     rv =  dynamic_cast<AliPHOSTrackSegmentMaker *>(PhosLoader()->TrackSegmentMaker()) ;
 254   }
 255   return rv ;
 256 }
 257
 258 //____________________________________________________________________________
 259 TClonesArray * AliPHOSGetter::RecParticles() const
 260 {
 261   // asks the Loader to return the TrackSegments container
 262
 263   TClonesArray * rv = 0 ;
 264
 265   rv = PhosLoader()->RecParticles() ;
 266   if (!rv) {
 267     PhosLoader()->MakeRecParticlesArray() ;
 268     rv = PhosLoader()->RecParticles() ;
 269   }
 270   return rv ;
 271 }
 272 //____________________________________________________________________________
 273 void AliPHOSGetter::Event(Int_t event, const char* opt)
 274 {
 275   // Reads the content of all Tree's S, D and R
 276
 277 //   if ( event >= MaxEvent() ) {
 278 //     Error("Event", "%d not found in TreeE !", event) ;
 279 //     return ;
 280 //   }
 281
 282   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 283
 284 //   // checks if we are dealing with test-beam data
 285 //   TBranch * btb = rl->TreeE()->GetBranch("AliPHOSBeamTestEvent") ;
 286 //   if(btb){
 287 //     if(!fBTE)
 288 //       fBTE = new AliPHOSBeamTestEvent() ;
 289 //     btb->SetAddress(&fBTE) ;
 290 //     btb->GetEntry(event) ;
 291 //   }
 292 //   else{
 293 //     if(fBTE){
 294 //       delete fBTE ;
 295 //       fBTE = 0 ;
 296 //     }
 297 //   }
 298
 299   // Loads the type of object(s) requested
 300
 301   rl->GetEvent(event) ;
 302
 303   if(strstr(opt,"X") || (strcmp(opt,"")==0)){
 304     ReadPrimaries() ;
 305   }
 306
 307   if(strstr(opt,"H")  ){
 308     ReadTreeH();
 309   }
 310
 311   if(strstr(opt,"S")  ){
 312     ReadTreeS() ;
 313   }
 314
 315   if(strstr(opt,"D") ){
 316     ReadTreeD() ;
 317   }
 318
 319   if(strstr(opt,"R") ){
 320     ReadTreeR() ;
 321   }
 322
 323   if( strstr(opt,"T") ){
 324     ReadTreeT() ;
 325   }
 326
 327   if( strstr(opt,"P") ){
 328     ReadTreeP() ;
 329   }
 330
 331   if( strstr(opt,"E") ){
 332     ReadTreeE(event) ;
 333   }
 334
 335 }
 336
 337
 338 //____________________________________________________________________________
 339 void AliPHOSGetter::Event(AliRawReader *rawReader, const char* opt)
 340 {
 341   // Reads the raw event from rawReader
 342
 343   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 344
 345   if( strstr(opt,"W")  ){
 346     rawReader->NextEvent();
 347     Int_t iEvent = rl->GetEventNumber();
 348     rl->GetEvent(iEvent);
 349     ReadRaw(rawReader) ;
 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() const
 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) const
 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, "QRO", "", 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, "QRO", "", 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   gLowGain ->SetTitle("Low gain channel");
 686   gHighGain->SetTitle("High gain channel");
 687
 688   while ( in.Next() ) { // PHOS entries loop
 689     if ( (in.IsNewRow() || in.IsNewColumn() || in.IsNewModule()) ) {
 690       if (!first) {
 691         FitRaw(lowGainFlag, gLowGain, gHighGain, signalF, energy, time) ;
 692         amp = (Int_t)energy;
 693         if (energy > 0) {
 694           new((*digits)[idigit]) AliPHOSDigit( -1, id, (Float_t)energy, time) ;
 695           idigit++ ;
 696         }
 697         Int_t index ;
 698         for (index = 0; index < PHOS()->GetRawFormatTimeBins(); index++) {
 699           gLowGain ->SetPoint(index,
 700                               index * PHOS()->GetRawFormatTimeMax() /
 701                               PHOS()->GetRawFormatTimeBins(), 0) ;
 702           gHighGain->SetPoint(index,
 703                               index * PHOS()->GetRawFormatTimeMax() /
 704                               PHOS()->GetRawFormatTimeBins(), 0) ;
 705         }
 706       }
 707       first = kFALSE ;
 708       relId[0] = in.GetModule() ;
 709       if ( relId[0] >= PHOS()->GetRawFormatLowGainOffset() ) {
 710         relId[0] -=  PHOS()->GetRawFormatLowGainOffset() ;
 711         lowGainFlag = kTRUE ;
 712       } else
 713         lowGainFlag = kFALSE ;
 714       relId[1] = 0 ;
 715       relId[2] = in.GetRow() ;
 716       relId[3] = in.GetColumn() ;
 717       PHOSGeometry()->RelToAbsNumbering(relId, id) ;
 718     }
 719     if (lowGainFlag)
 720       gLowGain->SetPoint(in.GetTime(),
 721                          in.GetTime()* PHOS()->GetRawFormatTimeMax() /
 722                          PHOS()->GetRawFormatTimeBins(),
 723                          in.GetSignal()) ;
 724     else
 725       gHighGain->SetPoint(in.GetTime(),
 726                           in.GetTime() * PHOS()->GetRawFormatTimeMax() /
 727                           PHOS()->GetRawFormatTimeBins(),
 728                           in.GetSignal() ) ;
 729
 730   } // PHOS entries loop
 731
 732   FitRaw(lowGainFlag, gLowGain, gHighGain, signalF, energy, time) ;
 733   amp = (Int_t)energy;
 734   if (energy > 0 ) {
 735     new((*digits)[idigit]) AliPHOSDigit( -1, id, (Float_t)energy, time) ;
 736     idigit++ ;
 737   }
 738   digits->Sort() ;
 739
 740   delete signalF ;
 741   delete gLowGain;
 742   delete gHighGain ;
 743
 744   return digits->GetEntriesFast() ;
 745 }
 746
 747 //____________________________________________________________________________
 748 Int_t AliPHOSGetter::ReadTreeD()
 749 {
 750   // Read the Digits
 751
 752   PhosLoader()->CleanDigits() ;
 753   PhosLoader()->LoadDigits("UPDATE") ;
 754   PhosLoader()->LoadDigitizer("UPDATE") ;
 755   return Digits()->GetEntries() ;
 756 }
 757
 758 //____________________________________________________________________________
 759 Int_t AliPHOSGetter::ReadTreeH()
 760 {
 761   // Read the Hits
 762   PhosLoader()->CleanHits() ;
 763   // gets TreeH from the root file (PHOS.Hit.root)
 764   //if ( !IsLoaded("H") ) {
 765     PhosLoader()->LoadHits("UPDATE") ;
 766   //  SetLoaded("H") ;
 767   //}
 768   return Hits()->GetEntries() ;
 769 }
 770
 771 //____________________________________________________________________________
 772 Int_t AliPHOSGetter::ReadTreeR()
 773 {
 774   // Read the RecPoints
 775
 776   PhosLoader()->CleanRecPoints() ;
 777   // gets TreeR from the root file (PHOS.RecPoints.root)
 778   //if ( !IsLoaded("R") ) {
 779     PhosLoader()->LoadRecPoints("UPDATE") ;
 780     PhosLoader()->LoadClusterizer("UPDATE") ;
 781     //  SetLoaded("R") ;
 782     //}
 783
 784   return EmcRecPoints()->GetEntries() ;
 785 }
 786
 787 //____________________________________________________________________________
 788 Int_t AliPHOSGetter::ReadTreeT()
 789 {
 790   // Read the TrackSegments
 791
 792   PhosLoader()->CleanTracks() ;
 793   // gets TreeT from the root file (PHOS.TrackSegments.root)
 794   //if ( !IsLoaded("T") ) {
 795     PhosLoader()->LoadTracks("UPDATE") ;
 796     PhosLoader()->LoadTrackSegmentMaker("UPDATE") ;
 797     //    SetLoaded("T") ;
 798     //}
 799
 800   return TrackSegments()->GetEntries() ;
 801 }
 802 //____________________________________________________________________________
 803 Int_t AliPHOSGetter::ReadTreeP()
 804 {
 805   // Read the RecParticles
 806
 807   PhosLoader()->CleanRecParticles() ;
 808
 809   // gets TreeT from the root file (PHOS.TrackSegments.root)
 810   //  if ( !IsLoaded("P") ) {
 811     PhosLoader()->LoadRecParticles("UPDATE") ;
 812     PhosLoader()->LoadPID("UPDATE") ;
 813     //  SetLoaded("P") ;
 814     //}
 815
 816   return RecParticles()->GetEntries() ;
 817 }
 818 //____________________________________________________________________________
 819 Int_t AliPHOSGetter::ReadTreeS()
 820 {
 821   // Read the SDigits
 822
 823   PhosLoader()->CleanSDigits() ;
 824   // gets TreeS from the root file (PHOS.SDigits.root)
 825   //if ( !IsLoaded("S") ) {
 826     PhosLoader()->LoadSDigits("READ") ;
 827     PhosLoader()->LoadSDigitizer("READ") ;
 828     //  SetLoaded("S") ;
 829     //}
 830
 831   return SDigits()->GetEntries() ;
 832 }
 833
 834 //____________________________________________________________________________
 835 Int_t AliPHOSGetter::ReadTreeE(Int_t event)
 836 {
 837   // Read the ESD
 838
 839   // gets esdTree from the root file (AliESDs.root)
 840   if (!fESDFile)
 841     if ( !OpenESDFile() )
 842       return -1 ;
 843
 844   fESDTree = static_cast<TTree*>(fESDFile->Get("esdTree")) ;
 845   fESD = new AliESD;
 846    if (!fESDTree) {
 847
 848      Error("ReadTreeE", "no ESD tree found");
 849      return -1;
 850    }
 851    fESDTree->SetBranchAddress("ESD", &fESD);
 852    fESDTree->GetEvent(event);
 853
 854    return event ;
 855 }
 856
 857 //____________________________________________________________________________
 858 TClonesArray * AliPHOSGetter::SDigits() const
 859 {
 860   // asks the Loader to return the Digits container
 861
 862   TClonesArray * rv = 0 ;
 863
 864   rv = PhosLoader()->SDigits() ;
 865   if (!rv) {
 866     PhosLoader()->MakeSDigitsArray() ;
 867     rv = PhosLoader()->SDigits() ;
 868   }
 869   return rv ;
 870 }
 871
 872 //____________________________________________________________________________
 873 AliPHOSSDigitizer * AliPHOSGetter::SDigitizer()
 874 {
 875   // Returns pointer to the SDigitizer task
 876   AliPHOSSDigitizer * rv ;
 877   rv =  dynamic_cast<AliPHOSSDigitizer *>(PhosLoader()->SDigitizer()) ;
 878   if (!rv) {
 879     Event(0, "S") ;
 880     rv =  dynamic_cast<AliPHOSSDigitizer *>(PhosLoader()->SDigitizer()) ;
 881   }
 882   return rv ;
 883 }
 884
 885 //____________________________________________________________________________
 886 TParticle * AliPHOSGetter::Secondary(const TParticle* p, Int_t index) const
 887 {
 888   // Return first (index=1) or second (index=2) secondary particle of primary particle p
 889
 890   if(index <= 0)
 891     return 0 ;
 892   if(index > 2)
 893     return 0 ;
 894
 895   if(p) {
 896   Int_t daughterIndex = p->GetDaughter(index-1) ;
 897   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 898   return  rl->GetAliRun()->GetMCApp()->Particle(daughterIndex) ;
 899   }
 900   else
 901     return 0 ;
 902 }
 903
 904 //____________________________________________________________________________
 905 void AliPHOSGetter::Track(Int_t itrack)
 906 {
 907   // Read the first entry of PHOS branch in hit tree gAlice->TreeH()
 908
 909  AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
 910
 911   if( !TreeH() ) // load treeH the first time
 912     rl->LoadHits() ;
 913
 914   // first time create the container
 915   TClonesArray * hits = Hits() ;
 916   if ( hits )
 917     hits->Clear() ;
 918
 919   TBranch * phosbranch = dynamic_cast<TBranch*>(TreeH()->GetBranch("PHOS")) ;
 920   phosbranch->SetAddress(&hits) ;
 921   phosbranch->GetEntry(itrack) ;
 922 }
 923
 924 //____________________________________________________________________________
 925 TTree * AliPHOSGetter::TreeD() const
 926 {
 927   // Returns pointer to the Digits Tree
 928   TTree * rv = 0 ;
 929   rv = PhosLoader()->TreeD() ;
 930   if ( !rv ) {
 931     PhosLoader()->MakeTree("D");
 932     rv = PhosLoader()->TreeD() ;
 933   }
 934
 935   return rv ;
 936 }
 937
 938 //____________________________________________________________________________
 939 TTree * AliPHOSGetter::TreeH() const
 940 {
 941   // Returns pointer to the Hits Tree
 942   TTree * rv = 0 ;
 943   rv = PhosLoader()->TreeH() ;
 944   if ( !rv ) {
 945     PhosLoader()->MakeTree("H");
 946     rv = PhosLoader()->TreeH() ;
 947   }
 948
 949   return rv ;
 950 }
 951
 952 //____________________________________________________________________________
 953 TTree * AliPHOSGetter::TreeR() const
 954 {
 955   // Returns pointer to the RecPoints Tree
 956   TTree * rv = 0 ;
 957   rv = PhosLoader()->TreeR() ;
 958   if ( !rv ) {
 959     PhosLoader()->MakeTree("R");
 960     rv = PhosLoader()->TreeR() ;
 961   }
 962
 963   return rv ;
 964 }
 965
 966 //____________________________________________________________________________
 967 TTree * AliPHOSGetter::TreeT() const
 968 {
 969   // Returns pointer to the TrackSegments Tree
 970   TTree * rv = 0 ;
 971   rv = PhosLoader()->TreeT() ;
 972   if ( !rv ) {
 973     PhosLoader()->MakeTree("T");
 974     rv = PhosLoader()->TreeT() ;
 975   }
 976
 977   return rv ;
 978 }
 979 //____________________________________________________________________________
 980 TTree * AliPHOSGetter::TreeP() const
 981 {
 982   // Returns pointer to the RecParticles  Tree
 983   TTree * rv = 0 ;
 984   rv = PhosLoader()->TreeP() ;
 985   if ( !rv ) {
 986     PhosLoader()->MakeTree("P");
 987     rv = PhosLoader()->TreeP() ;
 988   }
 989
 990   return rv ;
 991 }
 992
 993 //____________________________________________________________________________
 994 TTree * AliPHOSGetter::TreeS() const
 995 {
 996  // Returns pointer to the SDigits Tree
 997   TTree * rv = 0 ;
 998   rv = PhosLoader()->TreeS() ;
 999   if ( !rv ) {
1000     PhosLoader()->MakeTree("S");
1001     rv = PhosLoader()->TreeS() ;
1002   }
1003
1004   return rv ;
1005 }
1006
1007 //____________________________________________________________________________
1008 Bool_t AliPHOSGetter::VersionExists(TString & opt) const
1009 {
1010   // checks if the version with the present name already exists in the same directory
1011
1012   Bool_t rv = kFALSE ;
1013
1014   AliRunLoader * rl = AliRunLoader::GetRunLoader(PhosLoader()->GetTitle());
1015   TString version( rl->GetEventFolder()->GetName() ) ;
1016
1017   opt.ToLower() ;
1018
1019   if ( opt == "sdigits") {
1020     // add the version name to the root file name
1021     TString fileName( PhosLoader()->GetSDigitsFileName() ) ;
1022     if (version != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
1023       fileName = fileName.ReplaceAll(".root", "") + "_" + version + ".root" ;
1024     if ( !(gSystem->AccessPathName(fileName)) ) {
1025       Warning("VersionExists", "The file %s already exists", fileName.Data()) ;
1026       rv = kTRUE ;
1027     }
1028     PhosLoader()->SetSDigitsFileName(fileName) ;
1029   }
1030
1031   if ( opt == "digits") {
1032     // add the version name to the root file name
1033     TString fileName( PhosLoader()->GetDigitsFileName() ) ;
1034     if (version != AliConfig::GetDefaultEventFolderName()) // only if not the default folder name
1035       fileName = fileName.ReplaceAll(".root", "") + "_" + version + ".root" ;
1036     if ( !(gSystem->AccessPathName(fileName)) ) {
1037       Warning("VersionExists", "The file %s already exists", fileName.Data()) ;
1038       rv = kTRUE ;
1039     }
1040   }
1041
1042   return rv ;
1043
1044 }
1045
1046 //____________________________________________________________________________
1047 UShort_t AliPHOSGetter::EventPattern(void) const
1048 {
1049   // Return the pattern (trigger bit register) of the beam-test event
1050   if(fBTE)
1051     return fBTE->GetPattern() ;
1052   else
1053     return 0 ;
1054 }
1055 //____________________________________________________________________________
1056 Float_t AliPHOSGetter::BeamEnergy(void) const
1057 {
1058   // Return the beam energy of the beam-test event
1059   if(fBTE)
1060     return fBTE->GetBeamEnergy() ;
1061   else
1062     return 0 ;
1063 }
1064 //____________________________________________________________________________
1065
1066 AliPHOSCalibData* AliPHOSGetter::CalibData()
1067 {
1068   // Check if the instance of AliPHOSCalibData exists, and return it
1069
1070   return fgCalibData;
1071 }