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