Some fixes
[u/mrichter/AliRoot.git] / FMD / AliFMDInput.cxx
CommitLineData
a1f80595 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 **************************************************************************/
a1f80595 15/* $Id$ */
c2fc1258 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*/
a1f80595 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
69893a66 32#include "AliFMDDebug.h" // ALIFMDDEBUG_H ALILOG_H
8f6ee336 33#include "AliLoader.h" // ALILOADER_H
34#include "AliRunLoader.h" // ALIRUNLOADER_H
35#include "AliRun.h" // ALIRUN_H
36#include "AliStack.h" // ALISTACK_H
d760ea03 37#include "AliRawReaderFile.h" // ALIRAWREADERFILE_H
039842fe 38#include "AliRawReaderRoot.h" // ALIRAWREADERROOT_H
39#include "AliRawReaderDate.h" // ALIRAWREADERDATE_H
e1a9aea4 40#include "AliRawEventHeaderBase.h"
8f6ee336 41#include "AliFMD.h" // ALIFMD_H
a1f80595 42#include "AliFMDHit.h" // ALIFMDHIT_H
8f6ee336 43#include "AliFMDDigit.h" // ALIFMDDigit_H
02a27b50 44#include "AliFMDSDigit.h" // ALIFMDDigit_H
bf000c32 45#include "AliFMDRecPoint.h" // ALIFMDRECPOINT_H
d760ea03 46#include "AliFMDRawReader.h" // ALIFMDRAWREADER_H
5632c898 47#include "AliFMDGeometry.h"
bf000c32 48#include <AliESD.h>
49#include <AliESDFMD.h>
df137876 50#include <AliESDEvent.h>
1e8f773e 51#include <AliCDBManager.h>
52#include <AliCDBEntry.h>
90dbf5fb 53#include <AliAlignObjParams.h>
08d168d9 54#include <AliTrackReference.h>
8f6ee336 55#include <TTree.h> // ROOT_TTree
bf000c32 56#include <TChain.h> // ROOT_TChain
8f6ee336 57#include <TParticle.h> // ROOT_TParticle
58#include <TString.h> // ROOT_TString
59#include <TDatabasePDG.h> // ROOT_TDatabasePDG
60#include <TMath.h> // ROOT_TMath
61#include <TGeoManager.h> // ROOT_TGeoManager
bf000c32 62#include <TSystemDirectory.h> // ROOT_TSystemDirectory
8f6ee336 63#include <Riostream.h> // ROOT_Riostream
97b4001e 64#include <TFile.h> // ROOT_TFile
a1e17193 65#include <TStreamerInfo.h>
7003fbd0 66#include <TArrayF.h>
67
a1f80595 68//____________________________________________________________________
69ClassImp(AliFMDInput)
70#if 0
71 ; // This is here to keep Emacs for indenting the next line
72#endif
73
f5b098de 74//____________________________________________________________________
75const AliFMDInput::ETrees AliFMDInput::fgkAllLoads[] = { kHits,
76 kKinematics,
77 kDigits,
78 kSDigits,
79 kHeader,
80 kRecPoints,
81 kESD,
82 kRaw,
83 kGeometry,
84 kTrackRefs,
85 kRawCalib,
86 kUser };
a1f80595 87
88//____________________________________________________________________
89AliFMDInput::AliFMDInput()
42f1b2f5 90 : TNamed("AliFMDInput", "Input handler for various FMD data"),
91 fGAliceFile(""),
a1f80595 92 fLoader(0),
93 fRun(0),
94 fStack(0),
95 fFMDLoader(0),
d760ea03 96 fReader(0),
5cf05dbb 97 fFMDReader(0),
a1f80595 98 fFMD(0),
bf000c32 99 fESD(0),
df137876 100 fESDEvent(0),
a1f80595 101 fTreeE(0),
102 fTreeH(0),
08d168d9 103 fTreeTR(0),
a1f80595 104 fTreeD(0),
105 fTreeS(0),
106 fTreeR(0),
b5ee4425 107 fTreeA(0),
bf000c32 108 fChainE(0),
a1f80595 109 fArrayE(0),
110 fArrayH(0),
faf80567 111 fArrayTR(0),
a1f80595 112 fArrayD(0),
113 fArrayS(0),
bf000c32 114 fArrayR(0),
d760ea03 115 fArrayA(0),
17e542eb 116 fHeader(0),
b5ee4425 117 fGeoManager(0),
a1f80595 118 fTreeMask(0),
039842fe 119 fRawFile(""),
f5b098de 120 fInputDir("."),
17e542eb 121 fIsInit(kFALSE),
2180cab3 122 fEventCount(0),
123 fNEvents(-1)
a1f80595 124{
df137876 125
a1f80595 126 // Constructor of an FMD input object. Specify what data to read in
127 // using the AddLoad member function. Sub-classes should at a
128 // minimum overload the member function Event. A full job can be
129 // executed using the member function Run.
130}
131
132
133
134//____________________________________________________________________
135AliFMDInput::AliFMDInput(const char* gAliceFile)
42f1b2f5 136 : TNamed("AliFMDInput", "Input handler for various FMD data"),
137 fGAliceFile(gAliceFile),
a1f80595 138 fLoader(0),
139 fRun(0),
140 fStack(0),
141 fFMDLoader(0),
d760ea03 142 fReader(0),
5cf05dbb 143 fFMDReader(0),
a1f80595 144 fFMD(0),
bf000c32 145 fESD(0),
df137876 146 fESDEvent(0),
a1f80595 147 fTreeE(0),
148 fTreeH(0),
08d168d9 149 fTreeTR(0),
a1f80595 150 fTreeD(0),
151 fTreeS(0),
152 fTreeR(0),
b5ee4425 153 fTreeA(0),
bf000c32 154 fChainE(0),
a1f80595 155 fArrayE(0),
156 fArrayH(0),
faf80567 157 fArrayTR(0),
a1f80595 158 fArrayD(0),
159 fArrayS(0),
bf000c32 160 fArrayR(0),
d760ea03 161 fArrayA(0),
17e542eb 162 fHeader(0),
b5ee4425 163 fGeoManager(0),
a1f80595 164 fTreeMask(0),
039842fe 165 fRawFile(""),
f5b098de 166 fInputDir("."),
17e542eb 167 fIsInit(kFALSE),
2180cab3 168 fEventCount(0),
169 fNEvents(-1)
a1f80595 170{
df137876 171
a1f80595 172 // Constructor of an FMD input object. Specify what data to read in
173 // using the AddLoad member function. Sub-classes should at a
174 // minimum overload the member function Event. A full job can be
175 // executed using the member function Run.
176}
177
178//____________________________________________________________________
f5b098de 179void
180AliFMDInput::SetLoads(UInt_t mask)
181{
182 for (UInt_t i = 0; i < sizeof(mask); i++) {
183 if (!(mask & (1 << i))) continue;
184 const ETrees *ptype = fgkAllLoads;
185 do {
186 ETrees type = *ptype;
187 if (i != UInt_t(type)) continue;
188 AddLoad(type);
189 break;
190 } while (*ptype != kUser);
191 }
192}
193
194//____________________________________________________________________
195void
196AliFMDInput::SetLoads(const char* what)
197{
198 TString l(what);
199 TObjArray* ll = l.Tokenize(", ");
200 TIter next(ll);
201 TObject* os = 0;
202 while ((os = next())) {
203 ETrees type = ParseLoad(os->GetName());
204 AddLoad(type);
205 }
206}
207
208
209//____________________________________________________________________
210AliFMDInput::ETrees
372d25a5 211AliFMDInput::ParseLoad(const char* what)
212{
213 TString opt(what);
214 opt.ToLower();
f5b098de 215 const ETrees* ptype = fgkAllLoads;
216 do {
217 ETrees type = *ptype;
218 if (opt.Contains(TreeName(type,true), TString::kIgnoreCase))
219 return type;
220 } while (*ptype++);
221 return kUser;
222}
223//____________________________________________________________________
224const char*
225AliFMDInput::LoadedString(Bool_t dataOnly) const
226{
227 static TString ret;
228 if (!ret.IsNull()) return ret.Data();
229
230 const ETrees* ptype = fgkAllLoads;
231 do {
232 ETrees type = *ptype;
233 if (dataOnly &&
234 (type == kKinematics ||
235 type == kHeader ||
236 type == kGeometry ||
237 type == kTrackRefs)) continue;
238 if (!IsLoaded(*ptype)) continue;
239
240 if (!ret.IsNull()) ret.Append(",");
241 ret.Append(TreeName(type));
242 } while (*ptype++ != kUser);
243 return ret.Data();
244}
245
246//____________________________________________________________________
247const char*
248AliFMDInput::TreeName(ETrees tree, Bool_t shortest)
249{
250 if (shortest) {
251 switch (tree) {
252 case kHits: return "hit";
253 case kKinematics: return "kin";
254 case kDigits: return "dig";
255 case kSDigits: return "sdig";
256 case kHeader: return "hea";
257 case kRecPoints: return "recp";
258 case kESD: return "esd";
259 case kRaw: return "raw";
260 case kGeometry: return "geo";
261 case kTrackRefs: return "trackr";
262 case kRawCalib: return "rawc";
263 case kUser: return "user";
264 }
265 return 0;
266 }
267 switch (tree) {
268 case kHits: return "Hits";
269 case kKinematics: return "Kinematics";
270 case kDigits: return "Digits";
271 case kSDigits: return "SDigits";
272 case kHeader: return "Header";
273 case kRecPoints: return "RecPoints";
274 case kESD: return "ESD";
275 case kRaw: return "Raw";
276 case kGeometry: return "Geometry";
277 case kTrackRefs: return "TrackRefs";
278 case kRawCalib: return "RawCalib";
279 case kUser: return "User";
280 }
372d25a5 281 return 0;
282}
f5b098de 283
372d25a5 284//____________________________________________________________________
a1f80595 285Int_t
286AliFMDInput::NEvents() const
287{
288 // Get number of events
f5b098de 289 if (IsLoaded(kRaw) ||
290 IsLoaded(kRawCalib)) return fReader->GetNumberOfEvents();
faf80567 291 if (fChainE) return fChainE->GetEntriesFast();
a1f80595 292 if (fTreeE) return fTreeE->GetEntries();
293 return -1;
294}
295
296//____________________________________________________________________
297Bool_t
298AliFMDInput::Init()
299{
300 // Initialize the object. Get the needed loaders, and such.
301
302 // Check if we have been initialized
303 if (fIsInit) {
304 AliWarning("Already initialized");
305 return fIsInit;
306 }
f5b098de 307 TString what;
308 const ETrees* ptype = fgkAllLoads;
309 do {
310 ETrees type = *ptype;
311 what.Append(Form("\n\t%-20s: %s", TreeName(type),
312 IsLoaded(type) ? "yes" : "no"));
313 } while (*ptype++);
314
315 Info("Init","Initialising w/mask 0x%04x%s", fTreeMask, what.Data());
a1f80595 316 // Get the run
f5b098de 317 if (IsLoaded(kDigits) ||
318 IsLoaded(kSDigits) ||
319 IsLoaded(kKinematics) ||
320 IsLoaded(kTrackRefs) ||
321 IsLoaded(kHeader)) {
faf80567 322 if (!gSystem->FindFile(".:/", fGAliceFile)) {
323 AliWarning(Form("Cannot find file %s in .:/", fGAliceFile.Data()));
59194467 324 }
faf80567 325 else {
326 fLoader = AliRunLoader::Open(fGAliceFile.Data(), "Alice", "read");
327 if (!fLoader) {
328 AliError(Form("Coulnd't read the file %s", fGAliceFile.Data()));
329 return kFALSE;
330 }
331 AliInfo(Form("Opened GAlice file %s", fGAliceFile.Data()));
332
333 if (fLoader->LoadgAlice()) return kFALSE;
334
335 fRun = fLoader->GetAliRun();
336
337 // Get the FMD
338 fFMD = static_cast<AliFMD*>(fRun->GetDetector("FMD"));
339 if (!fFMD) {
340 AliError("Failed to get detector FMD from loader");
341 return kFALSE;
342 }
343
344 // Get the FMD loader
345 fFMDLoader = fLoader->GetLoader("FMDLoader");
346 if (!fFMDLoader) {
347 AliError("Failed to get detector FMD loader from loader");
348 return kFALSE;
349 }
350 if (fLoader->LoadHeader()) {
351 AliError("Failed to get event header information from loader");
352 return kFALSE;
353 }
354 fTreeE = fLoader->TreeE();
59194467 355 }
a1f80595 356 }
8f6ee336 357
bf000c32 358 // Optionally, get the ESD files
f5b098de 359 if (IsLoaded(kESD)) {
360 fChainE = MakeChain("ESD", fInputDir, true);
df137876 361 fESDEvent = new AliESDEvent();
362 fESDEvent->ReadFromTree(fChainE);
363 // fChainE->SetBranchAddress("ESD", &fMainESD);
364
bf000c32 365 }
366
f5b098de 367 if (IsLoaded(kRaw) ||
368 IsLoaded(kRawCalib)) {
d760ea03 369 AliInfo("Getting FMD raw data digits");
370 fArrayA = new TClonesArray("AliFMDDigit");
0b1b2076 371#if 0
372 if (!fRawFile.IsNull() && fRawFile.EndsWith(".root"))
039842fe 373 fReader = new AliRawReaderRoot(fRawFile.Data());
374 else if (!fRawFile.IsNull() && fRawFile.EndsWith(".raw"))
375 fReader = new AliRawReaderDate(fRawFile.Data());
376 else
377 fReader = new AliRawReaderFile(-1);
0b1b2076 378#else
379 if(!fRawFile.IsNull())
380 fReader = AliRawReader::Create(fRawFile.Data());
381 else
382 fReader = new AliRawReaderFile(-1);
383#endif
5cf05dbb 384 fFMDReader = new AliFMDRawReader(fReader, 0);
d760ea03 385 }
386
8f6ee336 387 // Optionally, get the geometry
f5b098de 388 if (IsLoaded(kGeometry)) {
5632c898 389 TString fname;
59194467 390 if (fRun) {
5632c898 391 fname = gSystem->DirName(fGAliceFile);
392 fname.Append("/geometry.root");
8f6ee336 393 }
5632c898 394 if (!gSystem->AccessPathName(fname.Data()))
395 fname = "";
1e8f773e 396 AliCDBManager* cdb = AliCDBManager::Instance();
5632c898 397 if (!cdb->IsDefaultStorageSet()) {
398 cdb->SetDefaultStorage("local://$ALICE_ROOT/OCDB");
399 cdb->SetRun(0);
400 }
401
402 AliGeomManager::LoadGeometry(fname.IsNull() ? 0 : fname.Data());
403
1e8f773e 404 AliCDBEntry* align = cdb->Get("FMD/Align/Data");
405 if (align) {
406 AliInfo("Got alignment data from CDB");
407 TClonesArray* array = dynamic_cast<TClonesArray*>(align->GetObject());
408 if (!array) {
409 AliWarning("Invalid align data from CDB");
410 }
411 else {
412 Int_t nAlign = array->GetEntries();
413 for (Int_t i = 0; i < nAlign; i++) {
90dbf5fb 414 AliAlignObjParams* a = static_cast<AliAlignObjParams*>(array->At(i));
1e8f773e 415 if (!a->ApplyToGeometry()) {
416 AliWarning(Form("Failed to apply alignment to %s",
b760c02e 417 a->GetSymName()));
1e8f773e 418 }
419 }
420 }
421 }
5632c898 422 AliFMDGeometry* geom = AliFMDGeometry::Instance();
423 geom->Init();
424 geom->InitTransformations();
8f6ee336 425 }
bf000c32 426
69893a66 427 fEventCount = 0;
a1f80595 428 fIsInit = kTRUE;
429 return fIsInit;
430}
431
432//____________________________________________________________________
433Bool_t
434AliFMDInput::Begin(Int_t event)
435{
436 // Called at the begining of each event. Per default, it gets the
437 // data trees and gets pointers to the output arrays. Users can
438 // overload this, but should call this member function in the
439 // overloaded member function of the derived class.
440
441 // Check if we have been initialized
442 if (!fIsInit) {
443 AliError("Not initialized");
444 return fIsInit;
445 }
7003fbd0 446
a1f80595 447 // Get the event
59194467 448 if (fLoader && fLoader->GetEvent(event)) return kFALSE;
a1f80595 449
450 // Possibly load global kinematics information
f5b098de 451 if (IsLoaded(kKinematics)) {
69893a66 452 // AliInfo("Getting kinematics");
42f1b2f5 453 if (fLoader->LoadKinematics("READ")) return kFALSE;
a1f80595 454 fStack = fLoader->Stack();
455 }
7003fbd0 456
a1f80595 457 // Possibly load FMD Hit information
f5b098de 458 if (IsLoaded(kHits)) {
69893a66 459 // AliInfo("Getting FMD hits");
42f1b2f5 460 if (!fFMDLoader || fFMDLoader->LoadHits("READ")) return kFALSE;
a1f80595 461 fTreeH = fFMDLoader->TreeH();
462 if (!fArrayH) fArrayH = fFMD->Hits();
463 }
08d168d9 464
465 // Possibly load FMD TrackReference information
f5b098de 466 if (IsLoaded(kTrackRefs)) {
08d168d9 467 // AliInfo("Getting FMD hits");
468 if (!fLoader || fLoader->LoadTrackRefs("READ")) return kFALSE;
469 fTreeTR = fLoader->TreeTR();
470 if (!fArrayTR) fArrayTR = new TClonesArray("AliTrackReference");
471 fTreeTR->SetBranchAddress("TrackReferences", &fArrayTR);
472 }
473
9b98d361 474 // Possibly load heaedr information
f5b098de 475 if (IsLoaded(kHeader)) {
9b98d361 476 // AliInfo("Getting FMD hits");
477 if (!fLoader /* || fLoader->LoadHeader()*/) return kFALSE;
478 fHeader = fLoader->GetHeader();
479 }
480
a1f80595 481 // Possibly load FMD Digit information
f5b098de 482 if (IsLoaded(kDigits)) {
69893a66 483 // AliInfo("Getting FMD digits");
42f1b2f5 484 if (!fFMDLoader || fFMDLoader->LoadDigits("READ")) return kFALSE;
a1f80595 485 fTreeD = fFMDLoader->TreeD();
a9579262 486 if (fTreeD) {
487 if (!fArrayD) fArrayD = fFMD->Digits();
488 }
489 else {
490 fArrayD = 0;
491 AliWarning(Form("Failed to load FMD Digits"));
492 }
a1f80595 493 }
7003fbd0 494
a1f80595 495 // Possibly load FMD Sdigit information
f5b098de 496 if (IsLoaded(kSDigits)) {
69893a66 497 // AliInfo("Getting FMD summable digits");
faf80567 498 if (!fFMDLoader || fFMDLoader->LoadSDigits("READ")) {
499 AliWarning("Failed to load SDigits!");
500 return kFALSE;
501 }
a1f80595 502 fTreeS = fFMDLoader->TreeS();
503 if (!fArrayS) fArrayS = fFMD->SDigits();
504 }
7003fbd0 505
a1f80595 506 // Possibly load FMD RecPoints information
f5b098de 507 if (IsLoaded(kRecPoints)) {
69893a66 508 // AliInfo("Getting FMD reconstructed points");
42f1b2f5 509 if (!fFMDLoader || fFMDLoader->LoadRecPoints("READ")) return kFALSE;
a1f80595 510 fTreeR = fFMDLoader->TreeR();
bf000c32 511 if (!fArrayR) fArrayR = new TClonesArray("AliFMDRecPoint");
512 fTreeR->SetBranchAddress("FMD", &fArrayR);
7003fbd0 513 }
514
515 // Possibly load FMD ESD information
f5b098de 516 if (IsLoaded(kESD)) {
69893a66 517 // AliInfo("Getting FMD event summary data");
bf000c32 518 Int_t read = fChainE->GetEntry(event);
519 if (read <= 0) return kFALSE;
df137876 520 fESD = fESDEvent->GetFMDData();
bf000c32 521 if (!fESD) return kFALSE;
522 }
7003fbd0 523
d760ea03 524 // Possibly load FMD Digit information
f5b098de 525 if (IsLoaded(kRaw) || IsLoaded(kRawCalib)) {
506dc39d 526 Bool_t mon = fRawFile.Contains("mem://");
5cf05dbb 527 // AliInfo("Getting FMD raw data digits");
506dc39d 528 if (mon) std::cout << "Waiting for event ..." << std::flush;
e1a9aea4 529 do {
530 if (!fReader->NextEvent()) {
506dc39d 531 if (mon) {
e1a9aea4 532 gSystem->Sleep(3);
533 continue;
534 }
535 return kFALSE;
536 }
537 UInt_t eventType = fReader->GetType();
538 if(eventType == AliRawEventHeaderBase::kPhysicsEvent ||
539 eventType == AliRawEventHeaderBase::kCalibrationEvent)
540 break;
541 } while (true);
506dc39d 542 if (mon) std::cout << "got it" << std::endl;
5cf05dbb 543 // AliFMDRawReader r(fReader, 0);
d760ea03 544 fArrayA->Clear();
5cf05dbb 545 fFMDReader->ReadAdcs(fArrayA);
59194467 546 AliFMDDebug(1, ("Got a total of %d digits", fArrayA->GetEntriesFast()));
d760ea03 547 }
69893a66 548 fEventCount++;
bf000c32 549 return kTRUE;
550}
551
552
553//____________________________________________________________________
554Bool_t
555AliFMDInput::Event()
556{
557 // Process one event. The default implementation one or more of
558 //
559 // - ProcessHits if the hits are loaded.
560 // - ProcessDigits if the digits are loaded.
561 // - ProcessSDigits if the sumbable digits are loaded.
562 // - ProcessRecPoints if the reconstructed points are loaded.
563 // - ProcessESD if the event summary data is loaded
564 //
f5b098de 565 if (IsLoaded(kHits)) if (!ProcessHits()) return kFALSE;
566 if (IsLoaded(kTrackRefs))if (!ProcessTrackRefs()) return kFALSE;
567 if (IsLoaded(kKinematics) &&
568 IsLoaded(kHits)) if (!ProcessTracks()) return kFALSE;
569 if (IsLoaded(kKinematics))if (!ProcessStack()) return kFALSE;
570 if (IsLoaded(kSDigits)) if (!ProcessSDigits()) return kFALSE;
571 if (IsLoaded(kDigits)) if (!ProcessDigits()) return kFALSE;
572 if (IsLoaded(kRaw)) if (!ProcessRawDigits()) return kFALSE;
573 if (IsLoaded(kRawCalib)) if (!ProcessRawCalibDigits())return kFALSE;
574 if (IsLoaded(kRecPoints))if (!ProcessRecPoints()) return kFALSE;
575 if (IsLoaded(kESD)) if (!ProcessESDs()) return kFALSE;
576 if (IsLoaded(kUser)) if (!ProcessUsers()) return kFALSE;
bf000c32 577
578 return kTRUE;
579}
580
581//____________________________________________________________________
582Bool_t
583AliFMDInput::ProcessHits()
584{
02a27b50 585 // Read the hit tree, and pass each hit to the member function
586 // ProcessHit.
bf000c32 587 if (!fTreeH) {
588 AliError("No hit tree defined");
589 return kFALSE;
590 }
ddaa8027 591 if (!fArrayH) {
592 AliError("No hit array defined");
593 return kFALSE;
594 }
595
bf000c32 596 Int_t nTracks = fTreeH->GetEntries();
597 for (Int_t i = 0; i < nTracks; i++) {
598 Int_t hitRead = fTreeH->GetEntry(i);
599 if (hitRead <= 0) continue;
ddaa8027 600
bf000c32 601 Int_t nHit = fArrayH->GetEntries();
602 if (nHit <= 0) continue;
ddaa8027 603
bf000c32 604 for (Int_t j = 0; j < nHit; j++) {
605 AliFMDHit* hit = static_cast<AliFMDHit*>(fArrayH->At(j));
606 if (!hit) continue;
ddaa8027 607
bf000c32 608 TParticle* track = 0;
f5b098de 609 if (IsLoaded(kKinematics) && fStack) {
bf000c32 610 Int_t trackno = hit->Track();
611 track = fStack->Particle(trackno);
612 }
613 if (!ProcessHit(hit, track)) return kFALSE;
614 }
615 }
616 return kTRUE;
617}
08d168d9 618//____________________________________________________________________
619Bool_t
620AliFMDInput::ProcessTrackRefs()
621{
622 // Read the reconstrcted points tree, and pass each reconstruction
623 // object (AliFMDRecPoint) to either ProcessRecPoint.
624 if (!fTreeTR) {
625 AliError("No track reference tree defined");
626 return kFALSE;
627 }
628 if (!fArrayTR) {
629 AliError("No track reference array defined");
630 return kFALSE;
631 }
bf000c32 632
08d168d9 633 Int_t nEv = fTreeTR->GetEntries();
634 for (Int_t i = 0; i < nEv; i++) {
635 Int_t trRead = fTreeTR->GetEntry(i);
636 if (trRead <= 0) continue;
637 Int_t nTrackRefs = fArrayTR->GetEntries();
638 for (Int_t j = 0; j < nTrackRefs; j++) {
5632c898 639 AliTrackReference* trackRef =
640 static_cast<AliTrackReference*>(fArrayTR->At(j));
08d168d9 641 if (!trackRef) continue;
5632c898 642 // if (trackRef->DetectorId() != AliTrackReference::kFMD) continue;
08d168d9 643 TParticle* track = 0;
f5b098de 644 if (IsLoaded(kKinematics) && fStack) {
08d168d9 645 Int_t trackno = trackRef->GetTrack();
646 track = fStack->Particle(trackno);
647 }
648 if (!ProcessTrackRef(trackRef,track)) return kFALSE;
649 }
650 }
651 return kTRUE;
652}
bf000c32 653//____________________________________________________________________
654Bool_t
f95a63c4 655AliFMDInput::ProcessTracks()
656{
657 // Read the hit tree, and pass each hit to the member function
ddaa8027 658 // ProcessTrack.
f95a63c4 659 if (!fStack) {
660 AliError("No track tree defined");
661 return kFALSE;
662 }
663 if (!fTreeH) {
664 AliError("No hit tree defined");
665 return kFALSE;
666 }
ddaa8027 667 if (!fArrayH) {
668 AliError("No hit array defined");
669 return kFALSE;
670 }
671
672 // Int_t nTracks = fStack->GetNtrack();
f95a63c4 673 Int_t nTracks = fTreeH->GetEntries();
674 for (Int_t i = 0; i < nTracks; i++) {
ddaa8027 675 Int_t trackno = nTracks - i - 1;
676 TParticle* track = fStack->Particle(trackno);
f95a63c4 677 if (!track) continue;
ddaa8027 678
679 // Get the hits for this track.
f95a63c4 680 Int_t hitRead = fTreeH->GetEntry(i);
ddaa8027 681 Int_t nHit = fArrayH->GetEntries();
682 if (nHit == 0 || hitRead <= 0) {
683 // Let user code see the track, even if there's no hits.
684 if (!ProcessTrack(trackno, track, 0)) return kFALSE;
685 continue;
f95a63c4 686 }
f95a63c4 687
ddaa8027 688 // Loop over the hits corresponding to this track.
f95a63c4 689 for (Int_t j = 0; j < nHit; j++) {
690 AliFMDHit* hit = static_cast<AliFMDHit*>(fArrayH->At(j));
ddaa8027 691 if (!ProcessTrack(trackno, track, hit)) return kFALSE;
692 }
f95a63c4 693 }
694 return kTRUE;
695}
5632c898 696//____________________________________________________________________
697Bool_t
698AliFMDInput::ProcessStack()
699{
700 // Read the hit tree, and pass each hit to the member function
701 // ProcessTrack.
702 if (!fStack) {
703 AliError("No track tree defined");
704 return kFALSE;
705 }
706 Int_t nTracks = fStack->GetNtrack();
707 for (Int_t i = 0; i < nTracks; i++) {
708 Int_t trackno = nTracks - i - 1;
709 TParticle* track = fStack->Particle(trackno);
710 if (!track) continue;
f95a63c4 711
5632c898 712 if (!ProcessParticle(trackno, track)) return kFALSE;
713 }
714 return kTRUE;
715}
f95a63c4 716//____________________________________________________________________
717Bool_t
bf000c32 718AliFMDInput::ProcessDigits()
719{
720 // Read the digit tree, and pass each digit to the member function
721 // ProcessDigit.
42f1b2f5 722 if (!fTreeD) {
723 AliError("No digit tree defined");
724 return kFALSE;
725 }
726 if (!fArrayD) {
727 AliError("No digit array defined");
728 return kFALSE;
729 }
730
bf000c32 731 Int_t nEv = fTreeD->GetEntries();
732 for (Int_t i = 0; i < nEv; i++) {
733 Int_t digitRead = fTreeD->GetEntry(i);
734 if (digitRead <= 0) continue;
735 Int_t nDigit = fArrayD->GetEntries();
42f1b2f5 736 AliFMDDebug(0, ("Got %5d digits for this event", nDigit));
bf000c32 737 if (nDigit <= 0) continue;
738 for (Int_t j = 0; j < nDigit; j++) {
739 AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayD->At(j));
740 if (!digit) continue;
741 if (!ProcessDigit(digit)) return kFALSE;
742 }
743 }
744 return kTRUE;
745}
746
747//____________________________________________________________________
748Bool_t
749AliFMDInput::ProcessSDigits()
750{
751 // Read the summable digit tree, and pass each sumable digit to the
752 // member function ProcessSdigit.
42f1b2f5 753 if (!fTreeS) {
754 AliWarning("No sdigit tree defined");
755 return kTRUE; // Empty SDigits is fine
756 }
757 if (!fArrayS) {
758 AliWarning("No sdigit array defined");
759 return kTRUE; // Empty SDigits is fine
760 }
761
762 Int_t nEv = fTreeS->GetEntries();
bf000c32 763 for (Int_t i = 0; i < nEv; i++) {
764 Int_t sdigitRead = fTreeS->GetEntry(i);
faf80567 765 if (sdigitRead <= 0) {
766 AliInfo(Form("Read nothing from tree"));
767 continue;
768 }
769 Int_t nSdigit = fArrayS->GetEntriesFast();
42f1b2f5 770 AliFMDDebug(0, ("Got %5d digits for this event", nSdigit));
faf80567 771 AliInfo(Form("Got %5d digits for this event", nSdigit));
bf000c32 772 if (nSdigit <= 0) continue;
773 for (Int_t j = 0; j < nSdigit; j++) {
774 AliFMDSDigit* sdigit = static_cast<AliFMDSDigit*>(fArrayS->At(j));
775 if (!sdigit) continue;
776 if (!ProcessSDigit(sdigit)) return kFALSE;
777 }
778 }
779 return kTRUE;
780}
781
782//____________________________________________________________________
783Bool_t
d760ea03 784AliFMDInput::ProcessRawDigits()
785{
786 // Read the digit tree, and pass each digit to the member function
787 // ProcessDigit.
42f1b2f5 788 if (!fArrayA) {
789 AliError("No raw digit array defined");
790 return kFALSE;
791 }
792
d760ea03 793 Int_t nDigit = fArrayA->GetEntries();
794 if (nDigit <= 0) return kTRUE;
795 for (Int_t j = 0; j < nDigit; j++) {
796 AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayA->At(j));
797 if (!digit) continue;
59194467 798 if (AliLog::GetDebugLevel("FMD","") >= 40 && j < 30)
799 digit->Print();
d760ea03 800 if (!ProcessRawDigit(digit)) return kFALSE;
801 }
802 return kTRUE;
803}
804
805//____________________________________________________________________
806Bool_t
e064ab4a 807AliFMDInput::ProcessRawCalibDigits()
808{
809 // Read the digit tree, and pass each digit to the member function
810 // ProcessDigit.
811 if (!fArrayA) {
812 AliError("No raw digit array defined");
813 return kFALSE;
814 }
815
816 Int_t nDigit = fArrayA->GetEntries();
817 if (nDigit <= 0) return kTRUE;
818 for (Int_t j = 0; j < nDigit; j++) {
819 AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayA->At(j));
820 if (!digit) continue;
821 if (AliLog::GetDebugLevel("FMD","") >= 40 && j < 30)
822 digit->Print();
823 if (!ProcessRawCalibDigit(digit)) return kFALSE;
824 }
825 return kTRUE;
826}
827
828//____________________________________________________________________
829Bool_t
bf000c32 830AliFMDInput::ProcessRecPoints()
831{
832 // Read the reconstrcted points tree, and pass each reconstruction
833 // object (AliFMDRecPoint) to either ProcessRecPoint.
42f1b2f5 834 if (!fTreeR) {
835 AliError("No recpoint tree defined");
836 return kFALSE;
837 }
838 if (!fArrayR) {
839 AliError("No recpoints array defined");
840 return kFALSE;
841 }
842
bf000c32 843 Int_t nEv = fTreeR->GetEntries();
844 for (Int_t i = 0; i < nEv; i++) {
845 Int_t recRead = fTreeR->GetEntry(i);
846 if (recRead <= 0) continue;
847 Int_t nRecPoint = fArrayR->GetEntries();
848 for (Int_t j = 0; j < nRecPoint; j++) {
849 AliFMDRecPoint* recPoint = static_cast<AliFMDRecPoint*>(fArrayR->At(j));
850 if (!recPoint) continue;
851 if (!ProcessRecPoint(recPoint)) return kFALSE;
852 }
a1f80595 853 }
854 return kTRUE;
855}
856
857//____________________________________________________________________
a9579262 858Bool_t
859AliFMDInput::ProcessESDs()
860{
861 // Process event summary data
862 if (!fESD) return kFALSE;
863 for (UShort_t det = 1; det <= 3; det++) {
864 Char_t rings[] = { 'I', (det == 1 ? '\0' : 'O'), '\0' };
865 for (Char_t* rng = rings; *rng != '\0'; rng++) {
866 UShort_t nsec = (*rng == 'I' ? 20 : 40);
867 UShort_t nstr = (*rng == 'I' ? 512 : 256);
868 for (UShort_t sec = 0; sec < nsec; sec++) {
869 for (UShort_t str = 0; str < nstr; str++) {
870 Float_t eta = fESD->Eta(det,*rng,sec,str);
871 Float_t mult = fESD->Multiplicity(det,*rng,sec,str);
97b4001e 872 if (!fESD->IsAngleCorrected())
873 mult *= TMath::Abs(TMath::Cos(2.*TMath::ATan(TMath::Exp(-eta))));
a9579262 874 if (!ProcessESD(det, *rng, sec, str, eta, mult)) continue;
875 }
876 }
877 }
878 }
879 return kTRUE;
880}
881
882//____________________________________________________________________
506dc39d 883Bool_t
884AliFMDInput::ProcessUsers()
885{
886 // Process event summary data
887 for (UShort_t det = 1; det <= 3; det++) {
888 Char_t rings[] = { 'I', (det == 1 ? '\0' : 'O'), '\0' };
889 for (Char_t* rng = rings; *rng != '\0'; rng++) {
890 UShort_t nsec = (*rng == 'I' ? 20 : 40);
891 UShort_t nstr = (*rng == 'I' ? 512 : 256);
892 for (UShort_t sec = 0; sec < nsec; sec++) {
893 for (UShort_t str = 0; str < nstr; str++) {
894 Float_t v = GetSignal(det,*rng,sec,str);
895 if (!ProcessUser(det, *rng, sec, str, v)) continue;
896 }
897 }
898 }
899 }
900 return kTRUE;
901}
902
903//____________________________________________________________________
a1f80595 904Bool_t
905AliFMDInput::End()
906{
907 // Called at the end of each event. Per default, it unloads the
908 // data trees and resets the pointers to the output arrays. Users
909 // can overload this, but should call this member function in the
910 // overloaded member function of the derived class.
911
912 // Check if we have been initialized
913 if (!fIsInit) {
914 AliError("Not initialized");
915 return fIsInit;
916 }
917 // Possibly unload global kinematics information
f5b098de 918 if (IsLoaded(kKinematics)) {
a1f80595 919 fLoader->UnloadKinematics();
920 // fTreeK = 0;
921 fStack = 0;
922 }
923 // Possibly unload FMD Hit information
f5b098de 924 if (IsLoaded(kHits)) {
a1f80595 925 fFMDLoader->UnloadHits();
926 fTreeH = 0;
927 }
928 // Possibly unload FMD Digit information
f5b098de 929 if (IsLoaded(kDigits)) {
a1f80595 930 fFMDLoader->UnloadDigits();
931 fTreeD = 0;
932 }
933 // Possibly unload FMD Sdigit information
f5b098de 934 if (IsLoaded(kSDigits)) {
a1f80595 935 fFMDLoader->UnloadSDigits();
936 fTreeS = 0;
937 }
938 // Possibly unload FMD RecPoints information
f5b098de 939 if (IsLoaded(kRecPoints)) {
a1f80595 940 fFMDLoader->UnloadRecPoints();
941 fTreeR = 0;
942 }
f48d9b11 943 // AliInfo("Now out event");
a1f80595 944 return kTRUE;
945}
946
947//____________________________________________________________________
948Bool_t
2180cab3 949AliFMDInput::Run(UInt_t maxEvents)
a1f80595 950{
951 // Run over all events and files references in galice.root
952
953 Bool_t retval;
954 if (!(retval = Init())) return retval;
955
2180cab3 956 fNEvents = NEvents();
957 if (fNEvents < 0) fNEvents = maxEvents;
958 else if (maxEvents > 0) fNEvents = TMath::Min(fNEvents,Int_t(maxEvents));
959
960 Int_t event = 0;
961 for (; fNEvents < 0 || event < fNEvents; event++) {
962 printf("\rEvent %8d/%8d ...", event, fNEvents);
a1f80595 963 if (!(retval = Begin(event))) break;
964 if (!(retval = Event())) break;
965 if (!(retval = End())) break;
966 }
2180cab3 967 printf("Looped over %8d events\n", event+1);
a1f80595 968 if (!retval) return retval;
969 retval = Finish();
970 return retval;
971}
972
69893a66 973//__________________________________________________________________
974TArrayF
975AliFMDInput::MakeLogScale(Int_t n, Double_t min, Double_t max)
976{
977 // Service function to define a logarithmic axis.
978 // Parameters:
979 // n Number of bins
980 // min Minimum of axis
981 // max Maximum of axis
982 TArrayF bins(n+1);
983 bins[0] = min;
984 if (n <= 20) {
985 for (Int_t i = 1; i < n+1; i++) bins[i] = bins[i-1] + (max-min)/n;
986 return bins;
987 }
988 Float_t dp = n / TMath::Log10(max / min);
989 Float_t pmin = TMath::Log10(min);
990 for (Int_t i = 1; i < n+1; i++) {
991 Float_t p = pmin + i / dp;
992 bins[i] = TMath::Power(10, p);
993 }
994 return bins;
995}
996
f5b098de 997//____________________________________________________________________
998void
999AliFMDInput::ScanDirectory(TSystemDirectory* dir,
1000 const TString& olddir,
1001 TChain* chain,
1002 const char* pattern, bool recursive)
1003{
1004 // Get list of files, and go back to old working directory
1005 TString oldDir(gSystem->WorkingDirectory());
1006 TList* files = dir->GetListOfFiles();
1007 gSystem->ChangeDirectory(oldDir);
1008
1009 // Sort list of files and check if we should add it
1010 if (!files) return;
1011 files->Sort();
1012 TIter next(files);
1013 TSystemFile* file = 0;
1014 while ((file = static_cast<TSystemFile*>(next()))) {
1015 TString name(file->GetName());
1016
1017 // Ignore special links
1018 if (name == "." || name == "..") continue;
1019
1020 // Check if this is a directory
1021 if (file->IsDirectory()) {
1022 if (recursive)
1023 ScanDirectory(static_cast<TSystemDirectory*>(file),
1024 olddir, chain,
1025 pattern,recursive);
1026 continue;
1027 }
1028
1029 // If this is not a root file, ignore
1030 if (!name.EndsWith(".root")) continue;
1031
1032 // If this file does not contain the pattern, ignore
1033 if (!name.Contains(pattern)) continue;
1034 if (name.Contains("friends")) continue;
1035
1036 // Get the path
1037 TString data(Form("%s/%s", file->GetTitle(), name.Data()));
1038
1039 TFile* test = TFile::Open(data.Data(), "READ");
1040 if (!test || test->IsZombie()) {
1041 ::Warning("ScanDirectory", "Failed to open file %s", data.Data());
1042 continue;
1043 }
1044 test->Close();
1045 chain->Add(data);
1046 }
1047}
1048
1049//____________________________________________________________________
1050TChain*
1051AliFMDInput::MakeChain(const char* what, const char* datadir, bool recursive)
1052{
1053 TString w(what);
1054 w.ToUpper();
1055 const char* treeName = 0;
1056 const char* pattern = 0;
1057 if (w.Contains("ESD")) { treeName = "esdTree"; pattern = "AliESD"; }
1058 else if (w.Contains("MC")) { treeName = "TE"; pattern = "galice"; }
1059 else {
1060 ::Error("MakeChain", "Unknown mode '%s' (not one of ESD, or MC)", what);
1061 return 0;
1062 }
1063
1064 // --- Our data chain ----------------------------------------------
1065 TChain* chain = new TChain(treeName);
1066
1067 // --- Get list of ESDs --------------------------------------------
1068 // Open source directory, and make sure we go back to were we were
1069 TString oldDir(gSystem->WorkingDirectory());
1070 TSystemDirectory d(datadir, datadir);
1071 ScanDirectory(&d, oldDir, chain, pattern, recursive);
1072
1073 return chain;
1074}
69893a66 1075
a1f80595 1076
a1f80595 1077//____________________________________________________________________
1078//
1079// EOF
1080//