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 (fTreeD) {
 288       if (!fArrayD) fArrayD = fFMD->Digits();
 289     }
 290     else {
 291       fArrayD = 0;
 292       AliWarning(Form("Failed to load FMD Digits"));
 293     }
 294   }
 295   // Possibly load FMD Sdigit information
 296   if (TESTBIT(fTreeMask, kSDigits)) {
 297     AliInfo("Getting FMD summable digits");
 298     if (fFMDLoader->LoadSDigits()) return kFALSE;
 299     fTreeS = fFMDLoader->TreeS();
 300     if (!fArrayS) fArrayS = fFMD->SDigits();
 301   }
 302   // Possibly load FMD RecPoints information
 303   if (TESTBIT(fTreeMask, kRecPoints)) {
 304     AliInfo("Getting FMD reconstructed points");
 305     if (fFMDLoader->LoadRecPoints()) return kFALSE;
 306     fTreeR = fFMDLoader->TreeR();
 307     if (!fArrayR) fArrayR = new TClonesArray("AliFMDRecPoint");
 308     fTreeR->SetBranchAddress("FMD",  &fArrayR);
 309   }  // Possibly load FMD ESD information
 310   if (TESTBIT(fTreeMask, kESD)) {
 311     AliInfo("Getting FMD event summary data");
 312     Int_t read = fChainE->GetEntry(event);
 313     if (read <= 0) return kFALSE;
 314     fESD = fMainESD->GetFMDData();
 315     if (!fESD) return kFALSE;
 316     fESD->CheckNeedUShort(fChainE->GetFile());
 317   }
 318   // Possibly load FMD Digit information
 319   if (TESTBIT(fTreeMask, kRaw)) {
 320     AliInfo("Getting FMD raw data digits");
 321     if (!fReader->NextEvent()) return kFALSE;
 322     AliFMDRawReader r(fReader, 0);
 323     fArrayA->Clear();
 324     r.ReadAdcs(fArrayA);
 325   }
 326
 327   return kTRUE;
 328 }
 329
 330
 331 //____________________________________________________________________
 332 Bool_t
 333 AliFMDInput::Event()
 334 {
 335   // Process one event.  The default implementation one or more of
 336   //
 337   //   -  ProcessHits       if the hits are loaded.
 338   //   -  ProcessDigits     if the digits are loaded.
 339   //   -  ProcessSDigits    if the sumbable digits are loaded.
 340   //   -  ProcessRecPoints  if the reconstructed points are loaded.
 341   //   -  ProcessESD        if the event summary data is loaded
 342   //
 343   if (TESTBIT(fTreeMask, kHits))
 344     if (!ProcessHits()) return kFALSE;
 345   if (TESTBIT(fTreeMask, kDigits))
 346     if (!ProcessDigits()) return kFALSE;
 347   if (TESTBIT(fTreeMask, kSDigits))
 348     if (!ProcessSDigits()) return kFALSE;
 349   if (TESTBIT(fTreeMask, kRaw))
 350     if (!ProcessRawDigits()) return kFALSE;
 351   if (TESTBIT(fTreeMask, kRecPoints))
 352     if (!ProcessRecPoints()) return kFALSE;
 353   if (TESTBIT(fTreeMask, kESD))
 354     if (!ProcessESDs()) return kFALSE;
 355
 356   return kTRUE;
 357 }
 358
 359 //____________________________________________________________________
 360 Bool_t
 361 AliFMDInput::ProcessHits()
 362 {
 363   // Read the hit tree, and pass each hit to the member function
 364   // ProcessHit.
 365   if (!fTreeH) {
 366     AliError("No hit tree defined");
 367     return kFALSE;
 368   }
 369   Int_t nTracks = fTreeH->GetEntries();
 370   for (Int_t i = 0; i < nTracks; i++) {
 371     Int_t hitRead  = fTreeH->GetEntry(i);
 372     if (hitRead <= 0) continue;
 373     if (!fArrayH) {
 374       AliError("No hit array defined");
 375       return kFALSE;
 376     }
 377     Int_t nHit = fArrayH->GetEntries();
 378     if (nHit <= 0) continue;
 379     for (Int_t j = 0; j < nHit; j++) {
 380       AliFMDHit* hit = static_cast<AliFMDHit*>(fArrayH->At(j));
 381       if (!hit) continue;
 382       TParticle* track = 0;
 383       if (TESTBIT(fTreeMask, kKinematics) && fStack) {
 384         Int_t trackno = hit->Track();
 385         track = fStack->Particle(trackno);
 386       }
 387       if (!ProcessHit(hit, track)) return kFALSE;
 388     }
 389   }
 390   return kTRUE;
 391 }
 392
 393 //____________________________________________________________________
 394 Bool_t
 395 AliFMDInput::ProcessDigits()
 396 {
 397   // Read the digit tree, and pass each digit to the member function
 398   // ProcessDigit.
 399   Int_t nEv = fTreeD->GetEntries();
 400   for (Int_t i = 0; i < nEv; i++) {
 401     Int_t digitRead  = fTreeD->GetEntry(i);
 402     if (digitRead <= 0) continue;
 403     Int_t nDigit = fArrayD->GetEntries();
 404     if (nDigit <= 0) continue;
 405     for (Int_t j = 0; j < nDigit; j++) {
 406       AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayD->At(j));
 407       if (!digit) continue;
 408       if (!ProcessDigit(digit)) return kFALSE;
 409     }
 410   }
 411   return kTRUE;
 412 }
 413
 414 //____________________________________________________________________
 415 Bool_t
 416 AliFMDInput::ProcessSDigits()
 417 {
 418   // Read the summable digit tree, and pass each sumable digit to the
 419   // member function ProcessSdigit.
 420   Int_t nEv = fTreeD->GetEntries();
 421   for (Int_t i = 0; i < nEv; i++) {
 422     Int_t sdigitRead  = fTreeS->GetEntry(i);
 423     if (sdigitRead <= 0) continue;
 424     Int_t nSdigit = fArrayS->GetEntries();
 425     if (nSdigit <= 0) continue;
 426     for (Int_t j = 0; j < nSdigit; j++) {
 427       AliFMDSDigit* sdigit = static_cast<AliFMDSDigit*>(fArrayS->At(j));
 428       if (!sdigit) continue;
 429       if (!ProcessSDigit(sdigit)) return kFALSE;
 430     }
 431   }
 432   return kTRUE;
 433 }
 434
 435 //____________________________________________________________________
 436 Bool_t
 437 AliFMDInput::ProcessRawDigits()
 438 {
 439   // Read the digit tree, and pass each digit to the member function
 440   // ProcessDigit.
 441   Int_t nDigit = fArrayA->GetEntries();
 442   if (nDigit <= 0) return kTRUE;
 443   for (Int_t j = 0; j < nDigit; j++) {
 444     AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayA->At(j));
 445     if (!digit) continue;
 446     if (!ProcessRawDigit(digit)) return kFALSE;
 447   }
 448   return kTRUE;
 449 }
 450
 451 //____________________________________________________________________
 452 Bool_t
 453 AliFMDInput::ProcessRecPoints()
 454 {
 455   // Read the reconstrcted points tree, and pass each reconstruction
 456   // object (AliFMDRecPoint) to either ProcessRecPoint.
 457   Int_t nEv = fTreeR->GetEntries();
 458   for (Int_t i = 0; i < nEv; i++) {
 459     Int_t recRead  = fTreeR->GetEntry(i);
 460     if (recRead <= 0) continue;
 461     Int_t nRecPoint = fArrayR->GetEntries();
 462     for (Int_t j = 0; j < nRecPoint; j++) {
 463       AliFMDRecPoint* recPoint = static_cast<AliFMDRecPoint*>(fArrayR->At(j));
 464       if (!recPoint) continue;
 465       if (!ProcessRecPoint(recPoint)) return kFALSE;
 466     }
 467   }
 468   return kTRUE;
 469 }
 470
 471 //____________________________________________________________________
 472 Bool_t
 473 AliFMDInput::ProcessESDs()
 474 {
 475   // Process event summary data
 476   if (!fESD) return kFALSE;
 477   for (UShort_t det = 1; det <= 3; det++) {
 478     Char_t rings[] = { 'I', (det == 1 ? '\0' : 'O'), '\0' };
 479     for (Char_t* rng = rings; *rng != '\0'; rng++) {
 480       UShort_t nsec = (*rng == 'I' ?  20 :  40);
 481       UShort_t nstr = (*rng == 'I' ? 512 : 256);
 482       for (UShort_t sec = 0; sec < nsec; sec++) {
 483         for (UShort_t str = 0; str < nstr; str++) {
 484           Float_t eta  = fESD->Eta(det,*rng,sec,str);
 485           Float_t mult = fESD->Multiplicity(det,*rng,sec,str);
 486           if (!ProcessESD(det, *rng, sec, str, eta, mult)) continue;
 487         }
 488       }
 489     }
 490   }
 491   return kTRUE;
 492 }
 493
 494 //____________________________________________________________________
 495 Bool_t
 496 AliFMDInput::End()
 497 {
 498   // Called at the end of each event.  Per default, it unloads the
 499   // data trees and resets the pointers to the output arrays.   Users
 500   // can overload this, but should call this member function in the
 501   // overloaded member function of the derived class.
 502
 503   // Check if we have been initialized
 504   if (!fIsInit) {
 505     AliError("Not initialized");
 506     return fIsInit;
 507   }
 508   // Possibly unload global kinematics information
 509   if (TESTBIT(fTreeMask, kKinematics)) {
 510     fLoader->UnloadKinematics();
 511     // fTreeK = 0;
 512     fStack = 0;
 513   }
 514   // Possibly unload FMD Hit information
 515   if (TESTBIT(fTreeMask, kHits)) {
 516     fFMDLoader->UnloadHits();
 517     fTreeH = 0;
 518   }
 519   // Possibly unload FMD Digit information
 520   if (TESTBIT(fTreeMask, kDigits)) {
 521     fFMDLoader->UnloadDigits();
 522     fTreeD = 0;
 523   }
 524   // Possibly unload FMD Sdigit information
 525   if (TESTBIT(fTreeMask, kSDigits)) {
 526     fFMDLoader->UnloadSDigits();
 527     fTreeS = 0;
 528   }
 529   // Possibly unload FMD RecPoints information
 530   if (TESTBIT(fTreeMask, kRecPoints)) {
 531     fFMDLoader->UnloadRecPoints();
 532     fTreeR = 0;
 533   }
 534   AliInfo("Now out event");
 535   return kTRUE;
 536 }
 537
 538 //____________________________________________________________________
 539 Bool_t
 540 AliFMDInput::Run()
 541 {
 542   // Run over all events and files references in galice.root
 543
 544   Bool_t retval;
 545   if (!(retval = Init())) return retval;
 546
 547   Int_t nEvents = NEvents();
 548   for (Int_t event = 0; event < nEvents; event++) {
 549     if (!(retval = Begin(event))) break;
 550     if (!(retval = Event())) break;
 551     if (!(retval = End())) break;
 552   }
 553   if (!retval) return retval;
 554   retval = Finish();
 555   return retval;
 556 }
 557
 558
 559 //____________________________________________________________________
 560 //
 561 // EOF
 562 //