1 /**************************************************************************
2 * Copyright(c) 2004, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
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
21 //___________________________________________________________________
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.
29 // Latest changes by Christian Holm Christensen
31 #include "AliFMDInput.h" // ALIFMDHIT_H
32 #include "AliFMDDebug.h" // ALIFMDDEBUG_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 "AliRawReaderRoot.h" // ALIRAWREADERROOT_H
39 #include "AliRawReaderDate.h" // ALIRAWREADERDATE_H
40 #include "AliRawEventHeaderBase.h"
41 #include "AliFMD.h" // ALIFMD_H
42 #include "AliFMDHit.h" // ALIFMDHIT_H
43 #include "AliFMDDigit.h" // ALIFMDDigit_H
44 #include "AliFMDSDigit.h" // ALIFMDDigit_H
45 #include "AliFMDRecPoint.h" // ALIFMDRECPOINT_H
46 #include "AliFMDRawReader.h" // ALIFMDRAWREADER_H
48 #include <AliESDFMD.h>
49 #include <AliESDEvent.h>
50 #include <AliCDBManager.h>
51 #include <AliCDBEntry.h>
52 #include <AliAlignObjParams.h>
53 #include <AliTrackReference.h>
54 #include <TTree.h> // ROOT_TTree
55 #include <TChain.h> // ROOT_TChain
56 #include <TParticle.h> // ROOT_TParticle
57 #include <TString.h> // ROOT_TString
58 #include <TDatabasePDG.h> // ROOT_TDatabasePDG
59 #include <TMath.h> // ROOT_TMath
60 #include <TGeoManager.h> // ROOT_TGeoManager
61 #include <TSystemDirectory.h> // ROOT_TSystemDirectory
62 #include <Riostream.h> // ROOT_Riostream
63 #include <TFile.h> // ROOT_TFile
64 #include <TStreamerInfo.h>
67 //____________________________________________________________________
70 ; // This is here to keep Emacs for indenting the next line
74 //____________________________________________________________________
75 AliFMDInput::AliFMDInput()
76 : TNamed("AliFMDInput", "Input handler for various FMD data"),
110 // Constructor of an FMD input object. Specify what data to read in
111 // using the AddLoad member function. Sub-classes should at a
112 // minimum overload the member function Event. A full job can be
113 // executed using the member function Run.
118 //____________________________________________________________________
119 AliFMDInput::AliFMDInput(const char* gAliceFile)
120 : TNamed("AliFMDInput", "Input handler for various FMD data"),
121 fGAliceFile(gAliceFile),
154 // Constructor of an FMD input object. Specify what data to read in
155 // using the AddLoad member function. Sub-classes should at a
156 // minimum overload the member function Event. A full job can be
157 // executed using the member function Run.
160 //____________________________________________________________________
162 AliFMDInput::NEvents() const
164 // Get number of events
165 if (TESTBIT(fTreeMask, kRaw) ||
166 TESTBIT(fTreeMask, kRawCalib)) return fReader->GetNumberOfEvents();
167 if (fChainE) return fChainE->GetEntriesFast();
168 if (fTreeE) return fTreeE->GetEntries();
172 //____________________________________________________________________
176 // Initialize the object. Get the needed loaders, and such.
178 // Check if we have been initialized
180 AliWarning("Already initialized");
183 Info("Init","Initialising w/mask 0x%04x\n"
197 TESTBIT(fTreeMask, kHits),
198 TESTBIT(fTreeMask, kKinematics),
199 TESTBIT(fTreeMask, kDigits),
200 TESTBIT(fTreeMask, kSDigits),
201 TESTBIT(fTreeMask, kHeader),
202 TESTBIT(fTreeMask, kRecPoints),
203 TESTBIT(fTreeMask, kESD),
204 TESTBIT(fTreeMask, kRaw),
205 TESTBIT(fTreeMask, kRawCalib),
206 TESTBIT(fTreeMask, kGeometry),
207 TESTBIT(fTreeMask, kTracks),
208 TESTBIT(fTreeMask, kTrackRefs));
210 if (TESTBIT(fTreeMask, kDigits) ||
211 TESTBIT(fTreeMask, kSDigits) ||
212 TESTBIT(fTreeMask, kKinematics) ||
213 TESTBIT(fTreeMask, kTrackRefs) ||
214 TESTBIT(fTreeMask, kTracks) ||
215 TESTBIT(fTreeMask, kHeader)) {
216 if (!gSystem->FindFile(".:/", fGAliceFile)) {
217 AliWarning(Form("Cannot find file %s in .:/", fGAliceFile.Data()));
220 fLoader = AliRunLoader::Open(fGAliceFile.Data(), "Alice", "read");
222 AliError(Form("Coulnd't read the file %s", fGAliceFile.Data()));
225 AliInfo(Form("Opened GAlice file %s", fGAliceFile.Data()));
227 if (fLoader->LoadgAlice()) return kFALSE;
229 fRun = fLoader->GetAliRun();
232 fFMD = static_cast<AliFMD*>(fRun->GetDetector("FMD"));
234 AliError("Failed to get detector FMD from loader");
238 // Get the FMD loader
239 fFMDLoader = fLoader->GetLoader("FMDLoader");
241 AliError("Failed to get detector FMD loader from loader");
244 if (fLoader->LoadHeader()) {
245 AliError("Failed to get event header information from loader");
248 fTreeE = fLoader->TreeE();
252 // Optionally, get the ESD files
253 if (TESTBIT(fTreeMask, kESD)) {
254 fChainE = new TChain("esdTree");
255 TSystemDirectory dir(".",".");
256 TList* files = dir.GetListOfFiles();
257 TSystemFile* file = 0;
259 AliError("No files");
264 while ((file = static_cast<TSystemFile*>(next()))) {
265 TString fname(file->GetName());
266 if (fname.Contains("AliESDs")) fChainE->AddFile(fname.Data());
268 fESDEvent = new AliESDEvent();
269 fESDEvent->ReadFromTree(fChainE);
270 // fChainE->SetBranchAddress("ESD", &fMainESD);
274 if (TESTBIT(fTreeMask, kRaw) ||
275 TESTBIT(fTreeMask, kRawCalib)) {
276 AliInfo("Getting FMD raw data digits");
277 fArrayA = new TClonesArray("AliFMDDigit");
279 if (!fRawFile.IsNull() && fRawFile.EndsWith(".root"))
280 fReader = new AliRawReaderRoot(fRawFile.Data());
281 else if (!fRawFile.IsNull() && fRawFile.EndsWith(".raw"))
282 fReader = new AliRawReaderDate(fRawFile.Data());
284 fReader = new AliRawReaderFile(-1);
286 if(!fRawFile.IsNull())
287 fReader = AliRawReader::Create(fRawFile.Data());
289 fReader = new AliRawReaderFile(-1);
291 fFMDReader = new AliFMDRawReader(fReader, 0);
294 // Optionally, get the geometry
295 if (TESTBIT(fTreeMask, kGeometry)) {
297 TString fname(""/*fRun->GetGeometryFileName()*/);
298 // if (fname.IsNull()) {
299 Warning("Init", "No file name for the geometry from AliRun");
300 fname = gSystem->DirName(fGAliceFile);
301 fname.Append("/geometry.root");
303 fGeoManager = TGeoManager::Import(fname.Data());
305 Fatal("Init", "No geometry manager found");
310 AliGeomManager::LoadGeometry();
312 AliCDBManager* cdb = AliCDBManager::Instance();
313 AliCDBEntry* align = cdb->Get("FMD/Align/Data");
315 AliInfo("Got alignment data from CDB");
316 TClonesArray* array = dynamic_cast<TClonesArray*>(align->GetObject());
318 AliWarning("Invalid align data from CDB");
321 Int_t nAlign = array->GetEntries();
322 for (Int_t i = 0; i < nAlign; i++) {
323 AliAlignObjParams* a = static_cast<AliAlignObjParams*>(array->At(i));
324 if (!a->ApplyToGeometry()) {
325 AliWarning(Form("Failed to apply alignment to %s",
338 //____________________________________________________________________
340 AliFMDInput::Begin(Int_t event)
342 // Called at the begining of each event. Per default, it gets the
343 // data trees and gets pointers to the output arrays. Users can
344 // overload this, but should call this member function in the
345 // overloaded member function of the derived class.
347 // Check if we have been initialized
349 AliError("Not initialized");
354 if (fLoader && fLoader->GetEvent(event)) return kFALSE;
355 AliInfo(Form("Now in event %8d/%8d", event, NEvents()));
357 // Possibly load global kinematics information
358 if (TESTBIT(fTreeMask, kKinematics) || TESTBIT(fTreeMask, kTracks)) {
359 // AliInfo("Getting kinematics");
360 if (fLoader->LoadKinematics("READ")) return kFALSE;
361 fStack = fLoader->Stack();
364 // Possibly load FMD Hit information
365 if (TESTBIT(fTreeMask, kHits) || TESTBIT(fTreeMask, kTracks)) {
366 // AliInfo("Getting FMD hits");
367 if (!fFMDLoader || fFMDLoader->LoadHits("READ")) return kFALSE;
368 fTreeH = fFMDLoader->TreeH();
369 if (!fArrayH) fArrayH = fFMD->Hits();
372 // Possibly load FMD TrackReference information
373 if (TESTBIT(fTreeMask, kTrackRefs) || TESTBIT(fTreeMask, kTracks)) {
374 // AliInfo("Getting FMD hits");
375 if (!fLoader || fLoader->LoadTrackRefs("READ")) return kFALSE;
376 fTreeTR = fLoader->TreeTR();
377 if (!fArrayTR) fArrayTR = new TClonesArray("AliTrackReference");
378 fTreeTR->SetBranchAddress("TrackReferences", &fArrayTR);
381 // Possibly load heaedr information
382 if (TESTBIT(fTreeMask, kHeader)) {
383 // AliInfo("Getting FMD hits");
384 if (!fLoader /* || fLoader->LoadHeader()*/) return kFALSE;
385 fHeader = fLoader->GetHeader();
388 // Possibly load FMD Digit information
389 if (TESTBIT(fTreeMask, kDigits)) {
390 // AliInfo("Getting FMD digits");
391 if (!fFMDLoader || fFMDLoader->LoadDigits("READ")) return kFALSE;
392 fTreeD = fFMDLoader->TreeD();
394 if (!fArrayD) fArrayD = fFMD->Digits();
398 AliWarning(Form("Failed to load FMD Digits"));
402 // Possibly load FMD Sdigit information
403 if (TESTBIT(fTreeMask, kSDigits)) {
404 // AliInfo("Getting FMD summable digits");
405 if (!fFMDLoader || fFMDLoader->LoadSDigits("READ")) {
406 AliWarning("Failed to load SDigits!");
409 fTreeS = fFMDLoader->TreeS();
410 if (!fArrayS) fArrayS = fFMD->SDigits();
413 // Possibly load FMD RecPoints information
414 if (TESTBIT(fTreeMask, kRecPoints)) {
415 // AliInfo("Getting FMD reconstructed points");
416 if (!fFMDLoader || fFMDLoader->LoadRecPoints("READ")) return kFALSE;
417 fTreeR = fFMDLoader->TreeR();
418 if (!fArrayR) fArrayR = new TClonesArray("AliFMDRecPoint");
419 fTreeR->SetBranchAddress("FMD", &fArrayR);
422 // Possibly load FMD ESD information
423 if (TESTBIT(fTreeMask, kESD)) {
424 // AliInfo("Getting FMD event summary data");
425 Int_t read = fChainE->GetEntry(event);
426 if (read <= 0) return kFALSE;
427 fESD = fESDEvent->GetFMDData();
428 if (!fESD) return kFALSE;
430 TFile* f = fChainE->GetFile();
432 TObject* o = f->GetStreamerInfoList()->FindObject("AliFMDMap");
434 TStreamerInfo* info = static_cast<TStreamerInfo*>(o);
435 std::cout << "AliFMDMap class version read is "
436 << info->GetClassVersion() << std::endl;
439 // fESD->CheckNeedUShort(fChainE->GetFile());
443 // Possibly load FMD Digit information
444 if (TESTBIT(fTreeMask, kRaw) || TESTBIT(fTreeMask, kRawCalib)) {
445 Bool_t mon = fRawFile.Contains("mem://");
446 // AliInfo("Getting FMD raw data digits");
447 if (mon) std::cout << "Waiting for event ..." << std::flush;
449 if (!fReader->NextEvent()) {
456 UInt_t eventType = fReader->GetType();
457 if(eventType == AliRawEventHeaderBase::kPhysicsEvent ||
458 eventType == AliRawEventHeaderBase::kCalibrationEvent)
461 if (mon) std::cout << "got it" << std::endl;
462 // AliFMDRawReader r(fReader, 0);
464 fFMDReader->ReadAdcs(fArrayA);
465 AliFMDDebug(1, ("Got a total of %d digits", fArrayA->GetEntriesFast()));
472 //____________________________________________________________________
476 // Process one event. The default implementation one or more of
478 // - ProcessHits if the hits are loaded.
479 // - ProcessDigits if the digits are loaded.
480 // - ProcessSDigits if the sumbable digits are loaded.
481 // - ProcessRecPoints if the reconstructed points are loaded.
482 // - ProcessESD if the event summary data is loaded
484 if (TESTBIT(fTreeMask, kHits))
485 if (!ProcessHits()) return kFALSE;
486 if (TESTBIT(fTreeMask, kTrackRefs))
487 if (!ProcessTrackRefs()) return kFALSE;
488 if (TESTBIT(fTreeMask, kTracks))
489 if (!ProcessTracks()) return kFALSE;
490 if (TESTBIT(fTreeMask, kSDigits))
491 if (!ProcessSDigits()) return kFALSE;
492 if (TESTBIT(fTreeMask, kDigits))
493 if (!ProcessDigits()) return kFALSE;
494 if (TESTBIT(fTreeMask, kRaw))
495 if (!ProcessRawDigits()) return kFALSE;
496 if (TESTBIT(fTreeMask, kRawCalib))
497 if (!ProcessRawCalibDigits()) return kFALSE;
498 if (TESTBIT(fTreeMask, kRecPoints))
499 if (!ProcessRecPoints()) return kFALSE;
500 if (TESTBIT(fTreeMask, kESD))
501 if (!ProcessESDs()) return kFALSE;
502 if (TESTBIT(fTreeMask, kUser))
503 if (!ProcessUsers()) return kFALSE;
508 //____________________________________________________________________
510 AliFMDInput::ProcessHits()
512 // Read the hit tree, and pass each hit to the member function
515 AliError("No hit tree defined");
519 AliError("No hit array defined");
523 Int_t nTracks = fTreeH->GetEntries();
524 for (Int_t i = 0; i < nTracks; i++) {
525 Int_t hitRead = fTreeH->GetEntry(i);
526 if (hitRead <= 0) continue;
528 Int_t nHit = fArrayH->GetEntries();
529 if (nHit <= 0) continue;
531 for (Int_t j = 0; j < nHit; j++) {
532 AliFMDHit* hit = static_cast<AliFMDHit*>(fArrayH->At(j));
535 TParticle* track = 0;
536 if (TESTBIT(fTreeMask, kKinematics) && fStack) {
537 Int_t trackno = hit->Track();
538 track = fStack->Particle(trackno);
540 if (!ProcessHit(hit, track)) return kFALSE;
545 //____________________________________________________________________
547 AliFMDInput::ProcessTrackRefs()
549 // Read the reconstrcted points tree, and pass each reconstruction
550 // object (AliFMDRecPoint) to either ProcessRecPoint.
552 AliError("No track reference tree defined");
556 AliError("No track reference array defined");
560 Int_t nEv = fTreeTR->GetEntries();
561 for (Int_t i = 0; i < nEv; i++) {
562 Int_t trRead = fTreeTR->GetEntry(i);
563 if (trRead <= 0) continue;
564 Int_t nTrackRefs = fArrayTR->GetEntries();
565 for (Int_t j = 0; j < nTrackRefs; j++) {
566 AliTrackReference* trackRef = static_cast<AliTrackReference*>(fArrayTR->At(j));
567 if (!trackRef) continue;
568 if (trackRef->DetectorId() != AliTrackReference::kFMD) continue;
569 TParticle* track = 0;
570 if (TESTBIT(fTreeMask, kKinematics) && fStack) {
571 Int_t trackno = trackRef->GetTrack();
572 track = fStack->Particle(trackno);
574 if (!ProcessTrackRef(trackRef,track)) return kFALSE;
579 //____________________________________________________________________
581 AliFMDInput::ProcessTracks()
583 // Read the hit tree, and pass each hit to the member function
586 AliError("No track tree defined");
590 AliError("No hit tree defined");
594 AliError("No hit array defined");
598 // Int_t nTracks = fStack->GetNtrack();
599 Int_t nTracks = fTreeH->GetEntries();
600 for (Int_t i = 0; i < nTracks; i++) {
601 Int_t trackno = nTracks - i - 1;
602 TParticle* track = fStack->Particle(trackno);
603 if (!track) continue;
605 // Get the hits for this track.
606 Int_t hitRead = fTreeH->GetEntry(i);
607 Int_t nHit = fArrayH->GetEntries();
608 if (nHit == 0 || hitRead <= 0) {
609 // Let user code see the track, even if there's no hits.
610 if (!ProcessTrack(trackno, track, 0)) return kFALSE;
614 // Loop over the hits corresponding to this track.
615 for (Int_t j = 0; j < nHit; j++) {
616 AliFMDHit* hit = static_cast<AliFMDHit*>(fArrayH->At(j));
617 if (!ProcessTrack(trackno, track, hit)) return kFALSE;
623 //____________________________________________________________________
625 AliFMDInput::ProcessDigits()
627 // Read the digit tree, and pass each digit to the member function
630 AliError("No digit tree defined");
634 AliError("No digit array defined");
638 Int_t nEv = fTreeD->GetEntries();
639 for (Int_t i = 0; i < nEv; i++) {
640 Int_t digitRead = fTreeD->GetEntry(i);
641 if (digitRead <= 0) continue;
642 Int_t nDigit = fArrayD->GetEntries();
643 AliFMDDebug(0, ("Got %5d digits for this event", nDigit));
644 if (nDigit <= 0) continue;
645 for (Int_t j = 0; j < nDigit; j++) {
646 AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayD->At(j));
647 if (!digit) continue;
648 if (!ProcessDigit(digit)) return kFALSE;
654 //____________________________________________________________________
656 AliFMDInput::ProcessSDigits()
658 // Read the summable digit tree, and pass each sumable digit to the
659 // member function ProcessSdigit.
661 AliWarning("No sdigit tree defined");
662 return kTRUE; // Empty SDigits is fine
665 AliWarning("No sdigit array defined");
666 return kTRUE; // Empty SDigits is fine
669 Int_t nEv = fTreeS->GetEntries();
670 for (Int_t i = 0; i < nEv; i++) {
671 Int_t sdigitRead = fTreeS->GetEntry(i);
672 if (sdigitRead <= 0) {
673 AliInfo(Form("Read nothing from tree"));
676 Int_t nSdigit = fArrayS->GetEntriesFast();
677 AliFMDDebug(0, ("Got %5d digits for this event", nSdigit));
678 AliInfo(Form("Got %5d digits for this event", nSdigit));
679 if (nSdigit <= 0) continue;
680 for (Int_t j = 0; j < nSdigit; j++) {
681 AliFMDSDigit* sdigit = static_cast<AliFMDSDigit*>(fArrayS->At(j));
682 if (!sdigit) continue;
683 if (!ProcessSDigit(sdigit)) return kFALSE;
689 //____________________________________________________________________
691 AliFMDInput::ProcessRawDigits()
693 // Read the digit tree, and pass each digit to the member function
696 AliError("No raw digit array defined");
700 Int_t nDigit = fArrayA->GetEntries();
701 if (nDigit <= 0) return kTRUE;
702 for (Int_t j = 0; j < nDigit; j++) {
703 AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayA->At(j));
704 if (!digit) continue;
705 if (AliLog::GetDebugLevel("FMD","") >= 40 && j < 30)
707 if (!ProcessRawDigit(digit)) return kFALSE;
712 //____________________________________________________________________
714 AliFMDInput::ProcessRawCalibDigits()
716 // Read the digit tree, and pass each digit to the member function
719 AliError("No raw digit array defined");
723 Int_t nDigit = fArrayA->GetEntries();
724 if (nDigit <= 0) return kTRUE;
725 for (Int_t j = 0; j < nDigit; j++) {
726 AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayA->At(j));
727 if (!digit) continue;
728 if (AliLog::GetDebugLevel("FMD","") >= 40 && j < 30)
730 if (!ProcessRawCalibDigit(digit)) return kFALSE;
735 //____________________________________________________________________
737 AliFMDInput::ProcessRecPoints()
739 // Read the reconstrcted points tree, and pass each reconstruction
740 // object (AliFMDRecPoint) to either ProcessRecPoint.
742 AliError("No recpoint tree defined");
746 AliError("No recpoints array defined");
750 Int_t nEv = fTreeR->GetEntries();
751 for (Int_t i = 0; i < nEv; i++) {
752 Int_t recRead = fTreeR->GetEntry(i);
753 if (recRead <= 0) continue;
754 Int_t nRecPoint = fArrayR->GetEntries();
755 for (Int_t j = 0; j < nRecPoint; j++) {
756 AliFMDRecPoint* recPoint = static_cast<AliFMDRecPoint*>(fArrayR->At(j));
757 if (!recPoint) continue;
758 if (!ProcessRecPoint(recPoint)) return kFALSE;
764 //____________________________________________________________________
766 AliFMDInput::ProcessESDs()
768 // Process event summary data
769 if (!fESD) return kFALSE;
770 for (UShort_t det = 1; det <= 3; det++) {
771 Char_t rings[] = { 'I', (det == 1 ? '\0' : 'O'), '\0' };
772 for (Char_t* rng = rings; *rng != '\0'; rng++) {
773 UShort_t nsec = (*rng == 'I' ? 20 : 40);
774 UShort_t nstr = (*rng == 'I' ? 512 : 256);
775 for (UShort_t sec = 0; sec < nsec; sec++) {
776 for (UShort_t str = 0; str < nstr; str++) {
777 Float_t eta = fESD->Eta(det,*rng,sec,str);
778 Float_t mult = fESD->Multiplicity(det,*rng,sec,str);
779 if (!fESD->IsAngleCorrected())
780 mult *= TMath::Abs(TMath::Cos(2.*TMath::ATan(TMath::Exp(-eta))));
781 if (!ProcessESD(det, *rng, sec, str, eta, mult)) continue;
789 //____________________________________________________________________
791 AliFMDInput::ProcessUsers()
793 // Process event summary data
794 for (UShort_t det = 1; det <= 3; det++) {
795 Char_t rings[] = { 'I', (det == 1 ? '\0' : 'O'), '\0' };
796 for (Char_t* rng = rings; *rng != '\0'; rng++) {
797 UShort_t nsec = (*rng == 'I' ? 20 : 40);
798 UShort_t nstr = (*rng == 'I' ? 512 : 256);
799 for (UShort_t sec = 0; sec < nsec; sec++) {
800 for (UShort_t str = 0; str < nstr; str++) {
801 Float_t v = GetSignal(det,*rng,sec,str);
802 if (!ProcessUser(det, *rng, sec, str, v)) continue;
810 //____________________________________________________________________
814 // Called at the end of each event. Per default, it unloads the
815 // data trees and resets the pointers to the output arrays. Users
816 // can overload this, but should call this member function in the
817 // overloaded member function of the derived class.
819 // Check if we have been initialized
821 AliError("Not initialized");
824 // Possibly unload global kinematics information
825 if (TESTBIT(fTreeMask, kKinematics) || TESTBIT(fTreeMask, kTracks)) {
826 fLoader->UnloadKinematics();
830 // Possibly unload FMD Hit information
831 if (TESTBIT(fTreeMask, kHits) || TESTBIT(fTreeMask, kTracks)) {
832 fFMDLoader->UnloadHits();
835 // Possibly unload FMD Digit information
836 if (TESTBIT(fTreeMask, kDigits)) {
837 fFMDLoader->UnloadDigits();
840 // Possibly unload FMD Sdigit information
841 if (TESTBIT(fTreeMask, kSDigits)) {
842 fFMDLoader->UnloadSDigits();
845 // Possibly unload FMD RecPoints information
846 if (TESTBIT(fTreeMask, kRecPoints)) {
847 fFMDLoader->UnloadRecPoints();
850 // AliInfo("Now out event");
854 //____________________________________________________________________
858 // Run over all events and files references in galice.root
861 if (!(retval = Init())) return retval;
863 Int_t nEvents = NEvents();
864 for (Int_t event = 0; nEvents < 0 || event < nEvents; event++) {
865 if (!(retval = Begin(event))) break;
866 if (!(retval = Event())) break;
867 if (!(retval = End())) break;
869 if (!retval) return retval;
874 //__________________________________________________________________
876 AliFMDInput::MakeLogScale(Int_t n, Double_t min, Double_t max)
878 // Service function to define a logarithmic axis.
881 // min Minimum of axis
882 // max Maximum of axis
886 for (Int_t i = 1; i < n+1; i++) bins[i] = bins[i-1] + (max-min)/n;
889 Float_t dp = n / TMath::Log10(max / min);
890 Float_t pmin = TMath::Log10(min);
891 for (Int_t i = 1; i < n+1; i++) {
892 Float_t p = pmin + i / dp;
893 bins[i] = TMath::Power(10, p);
900 //____________________________________________________________________