FMD/AliFMDInput.cxx

   1 /**************************************************************************
   2  * Copyright(c) 2004, 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 /* $Id$ */
  16 /** @file    AliFMDInput.cxx
  17     @author  Christian Holm Christensen <cholm@nbi.dk>
  18     @date    Mon Mar 27 12:42:40 2006
  19     @brief   FMD utility classes for reading FMD data
  20 */
  21 //___________________________________________________________________
  22 //
  23 // The classes defined here, are utility classes for reading in data
  24 // for the FMD.  They are  put in a seperate library to not polute the
  25 // normal libraries.  The classes are intended to be used as base
  26 // classes for customized class that do some sort of analysis on the
  27 // various types of data produced by the FMD.
  28 //
  29 // Latest changes by Christian Holm Christensen
  30 //
  31 #include "AliFMDInput.h"        // ALIFMDHIT_H
  32 #include "AliLog.h"             // ALILOG_H
  33 #include "AliLoader.h"          // ALILOADER_H
  34 #include "AliRunLoader.h"       // ALIRUNLOADER_H
  35 #include "AliRun.h"             // ALIRUN_H
  36 #include "AliStack.h"           // ALISTACK_H
  37 #include "AliRawReaderFile.h"   // ALIRAWREADERFILE_H
  38 #include "AliFMD.h"             // ALIFMD_H
  39 #include "AliFMDHit.h"          // ALIFMDHIT_H
  40 #include "AliFMDDigit.h"        // ALIFMDDigit_H
  41 #include "AliFMDSDigit.h"       // ALIFMDDigit_H
  42 #include "AliFMDRecPoint.h"     // ALIFMDRECPOINT_H
  43 #include "AliFMDRawReader.h"    // ALIFMDRAWREADER_H
  44 #include <AliESD.h>
  45 #include <AliESDFMD.h>
  46 #include <AliCDBManager.h>
  47 #include <AliCDBEntry.h>
  48 #include <AliAlignObjAngles.h>
  49 #include <TTree.h>              // ROOT_TTree
  50 #include <TChain.h>             // ROOT_TChain
  51 #include <TParticle.h>          // ROOT_TParticle
  52 #include <TString.h>            // ROOT_TString
  53 #include <TDatabasePDG.h>       // ROOT_TDatabasePDG
  54 #include <TMath.h>              // ROOT_TMath
  55 #include <TGeoManager.h>        // ROOT_TGeoManager
  56 #include <TSystemDirectory.h>   // ROOT_TSystemDirectory
  57 #include <Riostream.h>          // ROOT_Riostream
  58
  59 //____________________________________________________________________
  60 ClassImp(AliFMDInput)
  61 #if 0
  62   ; // This is here to keep Emacs for indenting the next line
  63 #endif
  64
  65
  66 //____________________________________________________________________
  67 AliFMDInput::AliFMDInput()
  68   : fGAliceFile(""),
  69     fLoader(0),
  70     fRun(0),
  71     fStack(0),
  72     fFMDLoader(0),
  73     fReader(0),
  74     fFMD(0),
  75     fMainESD(0),
  76     fESD(0),
  77     fTreeE(0),
  78     fTreeH(0),
  79     fTreeD(0),
  80     fTreeS(0),
  81     fTreeR(0),
  82     fTreeA(0),
  83     fChainE(0),
  84     fArrayE(0),
  85     fArrayH(0),
  86     fArrayD(0),
  87     fArrayS(0),
  88     fArrayR(0),
  89     fArrayA(0),
  90     fGeoManager(0),
  91     fTreeMask(0),
  92     fIsInit(kFALSE)
  93 {
  94   // Constructor of an FMD input object.  Specify what data to read in
  95   // using the AddLoad member function.   Sub-classes should at a
  96   // minimum overload the member function Event.   A full job can be
  97   // executed using the member function Run.
  98 }
  99
 100
 101
 102 //____________________________________________________________________
 103 AliFMDInput::AliFMDInput(const char* gAliceFile)
 104   : fGAliceFile(gAliceFile),
 105     fLoader(0),
 106     fRun(0),
 107     fStack(0),
 108     fFMDLoader(0),
 109     fReader(0),
 110     fFMD(0),
 111     fMainESD(0),
 112     fESD(0),
 113     fTreeE(0),
 114     fTreeH(0),
 115     fTreeD(0),
 116     fTreeS(0),
 117     fTreeR(0),
 118     fTreeA(0),
 119     fChainE(0),
 120     fArrayE(0),
 121     fArrayH(0),
 122     fArrayD(0),
 123     fArrayS(0),
 124     fArrayR(0),
 125     fArrayA(0),
 126     fGeoManager(0),
 127     fTreeMask(0),
 128     fIsInit(kFALSE)
 129 {
 130   // Constructor of an FMD input object.  Specify what data to read in
 131   // using the AddLoad member function.   Sub-classes should at a
 132   // minimum overload the member function Event.   A full job can be
 133   // executed using the member function Run.
 134 }
 135
 136 //____________________________________________________________________
 137 Int_t
 138 AliFMDInput::NEvents() const
 139 {
 140   // Get number of events
 141   if (fTreeE) return fTreeE->GetEntries();
 142   return -1;
 143 }
 144
 145 //____________________________________________________________________
 146 Bool_t
 147 AliFMDInput::Init()
 148 {
 149   // Initialize the object.  Get the needed loaders, and such.
 150
 151   // Check if we have been initialized
 152   if (fIsInit) {
 153     AliWarning("Already initialized");
 154     return fIsInit;
 155   }
 156   if (fGAliceFile.IsNull()) fGAliceFile = "galice.root";
 157   // Get the loader
 158   fLoader = AliRunLoader::Open(fGAliceFile.Data(), "Alice", "read");
 159   if (!fLoader) {
 160     AliError(Form("Coulnd't read the file %s", fGAliceFile.Data()));
 161     return kFALSE;
 162   }
 163
 164   // Get the run
 165   if  (fLoader->LoadgAlice()) return kFALSE;
 166   fRun = fLoader->GetAliRun();
 167
 168   // Get the FMD
 169   fFMD = static_cast<AliFMD*>(fRun->GetDetector("FMD"));
 170   if (!fFMD) {
 171     AliError("Failed to get detector FMD from loader");
 172     return kFALSE;
 173   }
 174
 175   // Get the FMD loader
 176   fFMDLoader = fLoader->GetLoader("FMDLoader");
 177   if (!fFMDLoader) {
 178     AliError("Failed to get detector FMD loader from loader");
 179     return kFALSE;
 180   }
 181   if (fLoader->LoadHeader()) {
 182     AliError("Failed to get event header information from loader");
 183     return kFALSE;
 184   }
 185   fTreeE = fLoader->TreeE();
 186
 187   // Optionally, get the ESD files
 188   if (TESTBIT(fTreeMask, kESD)) {
 189     fChainE = new TChain("esdTree");
 190     TSystemDirectory dir(".",".");
 191     TList*           files = dir.GetListOfFiles();
 192     TSystemFile*     file = 0;
 193     if (!files) {
 194       AliError("No files");
 195       return kFALSE;
 196     }
 197     files->Sort();
 198     TIter            next(files);
 199     while ((file = static_cast<TSystemFile*>(next()))) {
 200       TString fname(file->GetName());
 201       if (fname.Contains("AliESDs")) fChainE->AddFile(fname.Data());
 202     }
 203     fChainE->SetBranchAddress("ESD", &fMainESD);
 204   }
 205
 206   if (TESTBIT(fTreeMask, kRaw)) {
 207     AliInfo("Getting FMD raw data digits");
 208     fArrayA = new TClonesArray("AliFMDDigit");
 209     fReader = new AliRawReaderFile(-1);
 210   }
 211
 212   // Optionally, get the geometry
 213   if (TESTBIT(fTreeMask, kGeometry)) {
 214     TString fname(fRun->GetGeometryFileName());
 215     if (fname.IsNull()) {
 216       Warning("Init", "No file name for the geometry from AliRun");
 217       fname = gSystem->DirName(fGAliceFile);
 218       fname.Append("/geometry.root");
 219     }
 220     fGeoManager = TGeoManager::Import(fname.Data());
 221     if (!fGeoManager) {
 222       Fatal("Init", "No geometry manager found");
 223       return kFALSE;
 224     }
 225     AliCDBManager* cdb   = AliCDBManager::Instance();
 226     AliCDBEntry*   align = cdb->Get("FMD/Align/Data");
 227     if (align) {
 228       AliInfo("Got alignment data from CDB");
 229       TClonesArray* array = dynamic_cast<TClonesArray*>(align->GetObject());
 230       if (!array) {
 231         AliWarning("Invalid align data from CDB");
 232       }
 233       else {
 234         Int_t nAlign = array->GetEntries();
 235         for (Int_t i = 0; i < nAlign; i++) {
 236           AliAlignObjAngles* a = static_cast<AliAlignObjAngles*>(array->At(i));
 237           if (!a->ApplyToGeometry()) {
 238             AliWarning(Form("Failed to apply alignment to %s",
 239                             a->GetSymName()));
 240           }
 241         }
 242       }
 243     }
 244   }
 245
 246
 247   fIsInit = kTRUE;
 248   return fIsInit;
 249 }
 250
 251 //____________________________________________________________________
 252 Bool_t
 253 AliFMDInput::Begin(Int_t event)
 254 {
 255   // Called at the begining of each event.  Per default, it gets the
 256   // data trees and gets pointers to the output arrays.   Users can
 257   // overload this, but should call this member function in the
 258   // overloaded member function of the derived class.
 259
 260   // Check if we have been initialized
 261   if (!fIsInit) {
 262     AliError("Not initialized");
 263     return fIsInit;
 264   }
 265   // Get the event
 266   if (fLoader->GetEvent(event)) return kFALSE;
 267   AliInfo(Form("Now in event %d/%d", event, NEvents()));
 268
 269   // Possibly load global kinematics information
 270   if (TESTBIT(fTreeMask, kKinematics)) {
 271     AliInfo("Getting kinematics");
 272     if (fLoader->LoadKinematics()) return kFALSE;
 273     fStack = fLoader->Stack();
 274   }
 275   // Possibly load FMD Hit information
 276   if (TESTBIT(fTreeMask, kHits)) {
 277     AliInfo("Getting FMD hits");
 278     if (fFMDLoader->LoadHits()) return kFALSE;
 279     fTreeH = fFMDLoader->TreeH();
 280     if (!fArrayH) fArrayH = fFMD->Hits();
 281   }
 282   // Possibly load FMD Digit information
 283   if (TESTBIT(fTreeMask, kDigits)) {
 284     AliInfo("Getting FMD digits");
 285     if (fFMDLoader->LoadDigits()) return kFALSE;
 286     fTreeD = fFMDLoader->TreeD();
 287     if (!fArrayD) fArrayD = fFMD->Digits();
 288   }
 289   // Possibly load FMD Sdigit information
 290   if (TESTBIT(fTreeMask, kSDigits)) {
 291     AliInfo("Getting FMD summable digits");
 292     if (fFMDLoader->LoadSDigits()) return kFALSE;
 293     fTreeS = fFMDLoader->TreeS();
 294     if (!fArrayS) fArrayS = fFMD->SDigits();
 295   }
 296   // Possibly load FMD RecPoints information
 297   if (TESTBIT(fTreeMask, kRecPoints)) {
 298     AliInfo("Getting FMD reconstructed points");
 299     if (fFMDLoader->LoadRecPoints()) return kFALSE;
 300     fTreeR = fFMDLoader->TreeR();
 301     if (!fArrayR) fArrayR = new TClonesArray("AliFMDRecPoint");
 302     fTreeR->SetBranchAddress("FMD",  &fArrayR);
 303   }  // Possibly load FMD ESD information
 304   if (TESTBIT(fTreeMask, kESD)) {
 305     AliInfo("Getting FMD event summary data");
 306     Int_t read = fChainE->GetEntry(event);
 307     if (read <= 0) return kFALSE;
 308     fESD = fMainESD->GetFMDData();
 309     if (!fESD) return kFALSE;
 310   }
 311   // Possibly load FMD Digit information
 312   if (TESTBIT(fTreeMask, kRaw)) {
 313     AliInfo("Getting FMD raw data digits");
 314     if (!fReader->NextEvent()) return kFALSE;
 315     AliFMDRawReader r(fReader, 0);
 316     fArrayA->Clear();
 317     r.ReadAdcs(fArrayA);
 318   }
 319
 320   return kTRUE;
 321 }
 322
 323
 324 //____________________________________________________________________
 325 Bool_t
 326 AliFMDInput::Event()
 327 {
 328   // Process one event.  The default implementation one or more of
 329   //
 330   //   -  ProcessHits       if the hits are loaded.
 331   //   -  ProcessDigits     if the digits are loaded.
 332   //   -  ProcessSDigits    if the sumbable digits are loaded.
 333   //   -  ProcessRecPoints  if the reconstructed points are loaded.
 334   //   -  ProcessESD        if the event summary data is loaded
 335   //
 336   if (TESTBIT(fTreeMask, kHits))
 337     if (!ProcessHits()) return kFALSE;
 338   if (TESTBIT(fTreeMask, kDigits))
 339     if (!ProcessDigits()) return kFALSE;
 340   if (TESTBIT(fTreeMask, kSDigits))
 341     if (!ProcessSDigits()) return kFALSE;
 342   if (TESTBIT(fTreeMask, kRaw))
 343     if (!ProcessRawDigits()) return kFALSE;
 344   if (TESTBIT(fTreeMask, kRecPoints))
 345     if (!ProcessRecPoints()) return kFALSE;
 346   if (TESTBIT(fTreeMask, kESD))
 347     if (!ProcessESD(fESD)) return kFALSE;
 348
 349   return kTRUE;
 350 }
 351
 352 //____________________________________________________________________
 353 Bool_t
 354 AliFMDInput::ProcessHits()
 355 {
 356   // Read the hit tree, and pass each hit to the member function
 357   // ProcessHit.
 358   if (!fTreeH) {
 359     AliError("No hit tree defined");
 360     return kFALSE;
 361   }
 362   Int_t nTracks = fTreeH->GetEntries();
 363   for (Int_t i = 0; i < nTracks; i++) {
 364     Int_t hitRead  = fTreeH->GetEntry(i);
 365     if (hitRead <= 0) continue;
 366     if (!fArrayH) {
 367       AliError("No hit array defined");
 368       return kFALSE;
 369     }
 370     Int_t nHit = fArrayH->GetEntries();
 371     if (nHit <= 0) continue;
 372     for (Int_t j = 0; j < nHit; j++) {
 373       AliFMDHit* hit = static_cast<AliFMDHit*>(fArrayH->At(j));
 374       if (!hit) continue;
 375       TParticle* track = 0;
 376       if (TESTBIT(fTreeMask, kKinematics) && fStack) {
 377         Int_t trackno = hit->Track();
 378         track = fStack->Particle(trackno);
 379       }
 380       if (!ProcessHit(hit, track)) return kFALSE;
 381     }
 382   }
 383   return kTRUE;
 384 }
 385
 386 //____________________________________________________________________
 387 Bool_t
 388 AliFMDInput::ProcessDigits()
 389 {
 390   // Read the digit tree, and pass each digit to the member function
 391   // ProcessDigit.
 392   Int_t nEv = fTreeD->GetEntries();
 393   for (Int_t i = 0; i < nEv; i++) {
 394     Int_t digitRead  = fTreeD->GetEntry(i);
 395     if (digitRead <= 0) continue;
 396     Int_t nDigit = fArrayD->GetEntries();
 397     if (nDigit <= 0) continue;
 398     for (Int_t j = 0; j < nDigit; j++) {
 399       AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayD->At(j));
 400       if (!digit) continue;
 401       if (!ProcessDigit(digit)) return kFALSE;
 402     }
 403   }
 404   return kTRUE;
 405 }
 406
 407 //____________________________________________________________________
 408 Bool_t
 409 AliFMDInput::ProcessSDigits()
 410 {
 411   // Read the summable digit tree, and pass each sumable digit to the
 412   // member function ProcessSdigit.
 413   Int_t nEv = fTreeD->GetEntries();
 414   for (Int_t i = 0; i < nEv; i++) {
 415     Int_t sdigitRead  = fTreeS->GetEntry(i);
 416     if (sdigitRead <= 0) continue;
 417     Int_t nSdigit = fArrayS->GetEntries();
 418     if (nSdigit <= 0) continue;
 419     for (Int_t j = 0; j < nSdigit; j++) {
 420       AliFMDSDigit* sdigit = static_cast<AliFMDSDigit*>(fArrayS->At(j));
 421       if (!sdigit) continue;
 422       if (!ProcessSDigit(sdigit)) return kFALSE;
 423     }
 424   }
 425   return kTRUE;
 426 }
 427
 428 //____________________________________________________________________
 429 Bool_t
 430 AliFMDInput::ProcessRawDigits()
 431 {
 432   // Read the digit tree, and pass each digit to the member function
 433   // ProcessDigit.
 434   Int_t nDigit = fArrayA->GetEntries();
 435   if (nDigit <= 0) return kTRUE;
 436   for (Int_t j = 0; j < nDigit; j++) {
 437     AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayA->At(j));
 438     if (!digit) continue;
 439     if (!ProcessRawDigit(digit)) return kFALSE;
 440   }
 441   return kTRUE;
 442 }
 443
 444 //____________________________________________________________________
 445 Bool_t
 446 AliFMDInput::ProcessRecPoints()
 447 {
 448   // Read the reconstrcted points tree, and pass each reconstruction
 449   // object (AliFMDRecPoint) to either ProcessRecPoint.
 450   Int_t nEv = fTreeR->GetEntries();
 451   for (Int_t i = 0; i < nEv; i++) {
 452     Int_t recRead  = fTreeR->GetEntry(i);
 453     if (recRead <= 0) continue;
 454     Int_t nRecPoint = fArrayR->GetEntries();
 455     for (Int_t j = 0; j < nRecPoint; j++) {
 456       AliFMDRecPoint* recPoint = static_cast<AliFMDRecPoint*>(fArrayR->At(j));
 457       if (!recPoint) continue;
 458       if (!ProcessRecPoint(recPoint)) return kFALSE;
 459     }
 460   }
 461   return kTRUE;
 462 }
 463
 464 //____________________________________________________________________
 465 Bool_t
 466 AliFMDInput::End()
 467 {
 468   // Called at the end of each event.  Per default, it unloads the
 469   // data trees and resets the pointers to the output arrays.   Users
 470   // can overload this, but should call this member function in the
 471   // overloaded member function of the derived class.
 472
 473   // Check if we have been initialized
 474   if (!fIsInit) {
 475     AliError("Not initialized");
 476     return fIsInit;
 477   }
 478   // Possibly unload global kinematics information
 479   if (TESTBIT(fTreeMask, kKinematics)) {
 480     fLoader->UnloadKinematics();
 481     // fTreeK = 0;
 482     fStack = 0;
 483   }
 484   // Possibly unload FMD Hit information
 485   if (TESTBIT(fTreeMask, kHits)) {
 486     fFMDLoader->UnloadHits();
 487     fTreeH = 0;
 488   }
 489   // Possibly unload FMD Digit information
 490   if (TESTBIT(fTreeMask, kDigits)) {
 491     fFMDLoader->UnloadDigits();
 492     fTreeD = 0;
 493   }
 494   // Possibly unload FMD Sdigit information
 495   if (TESTBIT(fTreeMask, kSDigits)) {
 496     fFMDLoader->UnloadSDigits();
 497     fTreeS = 0;
 498   }
 499   // Possibly unload FMD RecPoints information
 500   if (TESTBIT(fTreeMask, kRecPoints)) {
 501     fFMDLoader->UnloadRecPoints();
 502     fTreeR = 0;
 503   }
 504   AliInfo("Now out event");
 505   return kTRUE;
 506 }
 507
 508 //____________________________________________________________________
 509 Bool_t
 510 AliFMDInput::Run()
 511 {
 512   // Run over all events and files references in galice.root
 513
 514   Bool_t retval;
 515   if (!(retval = Init())) return retval;
 516
 517   Int_t nEvents = NEvents();
 518   for (Int_t event = 0; event < nEvents; event++) {
 519     if (!(retval = Begin(event))) break;
 520     if (!(retval = Event())) break;
 521     if (!(retval = End())) break;
 522   }
 523   if (!retval) return retval;
 524   retval = Finish();
 525   return retval;
 526 }
 527
 528
 529 //____________________________________________________________________
 530 //
 531 // EOF
 532 //