X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=FMD%2FAliFMDInput.cxx;h=49a60a433495c27c241d8bb4635cc8c9b99b2738;hb=ec6e635345a030dc256a85bfb212137830296ae7;hp=a8facd7819a7af66cc93290c549a346044dd0e0a;hpb=b5ee4425c93874c0ce306795eec31f8c0a87308d;p=u%2Fmrichter%2FAliRoot.git diff --git a/FMD/AliFMDInput.cxx b/FMD/AliFMDInput.cxx index a8facd7819a..49a60a43349 100644 --- a/FMD/AliFMDInput.cxx +++ b/FMD/AliFMDInput.cxx @@ -29,23 +29,29 @@ // Latest changes by Christian Holm Christensen // #include "AliFMDInput.h" // ALIFMDHIT_H -#include "AliLog.h" // ALILOG_H +#include "AliFMDDebug.h" // ALIFMDDEBUG_H ALILOG_H #include "AliLoader.h" // ALILOADER_H #include "AliRunLoader.h" // ALIRUNLOADER_H #include "AliRun.h" // ALIRUN_H #include "AliStack.h" // ALISTACK_H #include "AliRawReaderFile.h" // ALIRAWREADERFILE_H +#include "AliRawReaderRoot.h" // ALIRAWREADERROOT_H +#include "AliRawReaderDate.h" // ALIRAWREADERDATE_H +#include "AliRawEventHeaderBase.h" #include "AliFMD.h" // ALIFMD_H #include "AliFMDHit.h" // ALIFMDHIT_H #include "AliFMDDigit.h" // ALIFMDDigit_H #include "AliFMDSDigit.h" // ALIFMDDigit_H #include "AliFMDRecPoint.h" // ALIFMDRECPOINT_H #include "AliFMDRawReader.h" // ALIFMDRAWREADER_H +#include "AliFMDGeometry.h" #include #include +#include #include #include -#include +#include +#include #include // ROOT_TTree #include // ROOT_TChain #include // ROOT_TParticle @@ -55,6 +61,9 @@ #include // ROOT_TGeoManager #include // ROOT_TSystemDirectory #include // ROOT_Riostream +#include // ROOT_TFile +#include +#include //____________________________________________________________________ ClassImp(AliFMDInput) @@ -62,20 +71,36 @@ ClassImp(AliFMDInput) ; // This is here to keep Emacs for indenting the next line #endif +//____________________________________________________________________ +const AliFMDInput::ETrees AliFMDInput::fgkAllLoads[] = { kHits, + kKinematics, + kDigits, + kSDigits, + kHeader, + kRecPoints, + kESD, + kRaw, + kGeometry, + kTrackRefs, + kRawCalib, + kUser }; //____________________________________________________________________ AliFMDInput::AliFMDInput() - : fGAliceFile(""), + : TNamed("AliFMDInput", "Input handler for various FMD data"), + fGAliceFile(""), fLoader(0), fRun(0), fStack(0), fFMDLoader(0), fReader(0), + fFMDReader(0), fFMD(0), - fMainESD(0), fESD(0), + fESDEvent(0), fTreeE(0), fTreeH(0), + fTreeTR(0), fTreeD(0), fTreeS(0), fTreeR(0), @@ -83,14 +108,21 @@ AliFMDInput::AliFMDInput() fChainE(0), fArrayE(0), fArrayH(0), + fArrayTR(0), fArrayD(0), fArrayS(0), fArrayR(0), fArrayA(0), + fHeader(0), fGeoManager(0), fTreeMask(0), - fIsInit(kFALSE) + fRawFile(""), + fInputDir("."), + fIsInit(kFALSE), + fEventCount(0), + fNEvents(-1) { + // Constructor of an FMD input object. Specify what data to read in // using the AddLoad member function. Sub-classes should at a // minimum overload the member function Event. A full job can be @@ -101,17 +133,20 @@ AliFMDInput::AliFMDInput() //____________________________________________________________________ AliFMDInput::AliFMDInput(const char* gAliceFile) - : fGAliceFile(gAliceFile), + : TNamed("AliFMDInput", "Input handler for various FMD data"), + fGAliceFile(gAliceFile), fLoader(0), fRun(0), fStack(0), fFMDLoader(0), fReader(0), + fFMDReader(0), fFMD(0), - fMainESD(0), fESD(0), + fESDEvent(0), fTreeE(0), fTreeH(0), + fTreeTR(0), fTreeD(0), fTreeS(0), fTreeR(0), @@ -119,25 +154,142 @@ AliFMDInput::AliFMDInput(const char* gAliceFile) fChainE(0), fArrayE(0), fArrayH(0), + fArrayTR(0), fArrayD(0), fArrayS(0), fArrayR(0), fArrayA(0), + fHeader(0), fGeoManager(0), fTreeMask(0), - fIsInit(kFALSE) + fRawFile(""), + fInputDir("."), + fIsInit(kFALSE), + fEventCount(0), + fNEvents(-1) { + // Constructor of an FMD input object. Specify what data to read in // using the AddLoad member function. Sub-classes should at a // minimum overload the member function Event. A full job can be // executed using the member function Run. } +//____________________________________________________________________ +void +AliFMDInput::SetLoads(UInt_t mask) +{ + for (UInt_t i = 0; i < sizeof(mask); i++) { + if (!(mask & (1 << i))) continue; + const ETrees *ptype = fgkAllLoads; + do { + ETrees type = *ptype; + if (i != UInt_t(type)) continue; + AddLoad(type); + break; + } while (*ptype++ != kUser); + } +} + +//____________________________________________________________________ +void +AliFMDInput::SetLoads(const char* what) +{ + TString l(what); + TObjArray* ll = l.Tokenize(", "); + TIter next(ll); + TObject* os = 0; + while ((os = next())) { + ETrees type = ParseLoad(os->GetName()); + AddLoad(type); + } + delete ll; +} + + +//____________________________________________________________________ +AliFMDInput::ETrees +AliFMDInput::ParseLoad(const char* what) +{ + TString opt(what); + opt.ToLower(); + const ETrees* ptype = fgkAllLoads; + do { + ETrees type = *ptype; + if (opt.Contains(TreeName(type,true), TString::kIgnoreCase)) + return type; + } while (*ptype++ != kUser); + return kUser; +} +//____________________________________________________________________ +const char* +AliFMDInput::LoadedString(Bool_t dataOnly) const +{ + static TString ret; + if (!ret.IsNull()) return ret.Data(); + + const ETrees* ptype = fgkAllLoads; + do { + ETrees type = *ptype; + if (dataOnly && + (type == kKinematics || + type == kHeader || + type == kGeometry || + type == kTrackRefs)) continue; + if (!IsLoaded(*ptype)) continue; + + if (!ret.IsNull()) ret.Append(","); + ret.Append(TreeName(type)); + } while (*ptype++ != kUser); + return ret.Data(); +} + +//____________________________________________________________________ +const char* +AliFMDInput::TreeName(ETrees tree, Bool_t shortest) +{ + if (shortest) { + switch (tree) { + case kHits: return "hit"; + case kKinematics: return "kin"; + case kDigits: return "dig"; + case kSDigits: return "sdig"; + case kHeader: return "hea"; + case kRecPoints: return "recp"; + case kESD: return "esd"; + case kRaw: return "raw"; + case kGeometry: return "geo"; + case kTrackRefs: return "trackr"; + case kRawCalib: return "rawc"; + case kUser: return "user"; + } + return 0; + } + switch (tree) { + case kHits: return "Hits"; + case kKinematics: return "Kinematics"; + case kDigits: return "Digits"; + case kSDigits: return "SDigits"; + case kHeader: return "Header"; + case kRecPoints: return "RecPoints"; + case kESD: return "ESD"; + case kRaw: return "Raw"; + case kGeometry: return "Geometry"; + case kTrackRefs: return "TrackRefs"; + case kRawCalib: return "RawCalib"; + case kUser: return "User"; + } + return 0; +} + //____________________________________________________________________ Int_t AliFMDInput::NEvents() const { // Get number of events + if (IsLoaded(kRaw) || + IsLoaded(kRawCalib)) return fReader->GetNumberOfEvents(); + if (fChainE) return fChainE->GetEntriesFast(); if (fTreeE) return fTreeE->GetEntries(); return -1; } @@ -153,76 +305,103 @@ AliFMDInput::Init() AliWarning("Already initialized"); return fIsInit; } - if (fGAliceFile.IsNull()) fGAliceFile = "galice.root"; - // Get the loader - fLoader = AliRunLoader::Open(fGAliceFile.Data(), "Alice", "read"); - if (!fLoader) { - AliError(Form("Coulnd't read the file %s", fGAliceFile.Data())); - return kFALSE; - } + TString what; + const ETrees* ptype = fgkAllLoads; + do { + ETrees type = *ptype; + what.Append(Form("\n\t%-20s: %s", TreeName(type), + IsLoaded(type) ? "yes" : "no")); + } while (*ptype++ != kUser); + Info("Init","Initialising w/mask 0x%04x%s", fTreeMask, what.Data()); // Get the run - if (fLoader->LoadgAlice()) return kFALSE; - fRun = fLoader->GetAliRun(); - - // Get the FMD - fFMD = static_cast(fRun->GetDetector("FMD")); - if (!fFMD) { - AliError("Failed to get detector FMD from loader"); - return kFALSE; - } - - // Get the FMD loader - fFMDLoader = fLoader->GetLoader("FMDLoader"); - if (!fFMDLoader) { - AliError("Failed to get detector FMD loader from loader"); - return kFALSE; - } - if (fLoader->LoadHeader()) { - AliError("Failed to get event header information from loader"); - return kFALSE; + if (IsLoaded(kDigits) || + IsLoaded(kSDigits) || + IsLoaded(kKinematics) || + IsLoaded(kTrackRefs) || + IsLoaded(kHeader)) { + if (!gSystem->FindFile(".:/", fGAliceFile)) { + AliWarning(Form("Cannot find file %s in .:/", fGAliceFile.Data())); + } + else { + fLoader = AliRunLoader::Open(fGAliceFile.Data(), "Alice", "read"); + if (!fLoader) { + AliError(Form("Coulnd't read the file %s", fGAliceFile.Data())); + return kFALSE; + } + AliInfo(Form("Opened GAlice file %s", fGAliceFile.Data())); + + if (fLoader->LoadgAlice()) return kFALSE; + + fRun = fLoader->GetAliRun(); + + // Get the FMD + fFMD = static_cast(fRun->GetDetector("FMD")); + if (!fFMD) { + AliError("Failed to get detector FMD from loader"); + return kFALSE; + } + + // Get the FMD loader + fFMDLoader = fLoader->GetLoader("FMDLoader"); + if (!fFMDLoader) { + AliError("Failed to get detector FMD loader from loader"); + return kFALSE; + } + if (fLoader->LoadHeader()) { + AliError("Failed to get event header information from loader"); + return kFALSE; + } + fTreeE = fLoader->TreeE(); + } } - fTreeE = fLoader->TreeE(); // Optionally, get the ESD files - if (TESTBIT(fTreeMask, kESD)) { - fChainE = new TChain("esdTree"); - TSystemDirectory dir(".","."); - TList* files = dir.GetListOfFiles(); - TSystemFile* file = 0; - if (!files) { - AliError("No files"); - return kFALSE; - } - files->Sort(); - TIter next(files); - while ((file = static_cast(next()))) { - TString fname(file->GetName()); - if (fname.Contains("AliESDs")) fChainE->AddFile(fname.Data()); - } - fChainE->SetBranchAddress("ESD", &fMainESD); + if (IsLoaded(kESD)) { + fChainE = MakeChain("ESD", fInputDir, true); + fESDEvent = new AliESDEvent(); + fESDEvent->ReadFromTree(fChainE); + // fChainE->SetBranchAddress("ESD", &fMainESD); + } - if (TESTBIT(fTreeMask, kRaw)) { + if (IsLoaded(kRaw) || + IsLoaded(kRawCalib)) { AliInfo("Getting FMD raw data digits"); fArrayA = new TClonesArray("AliFMDDigit"); - fReader = new AliRawReaderFile(-1); +#if 0 + if (!fRawFile.IsNull() && fRawFile.EndsWith(".root")) + fReader = new AliRawReaderRoot(fRawFile.Data()); + else if (!fRawFile.IsNull() && fRawFile.EndsWith(".raw")) + fReader = new AliRawReaderDate(fRawFile.Data()); + else + fReader = new AliRawReaderFile(-1); +#else + if(!fRawFile.IsNull()) + fReader = AliRawReader::Create(fRawFile.Data()); + else + fReader = new AliRawReaderFile(-1); +#endif + fFMDReader = new AliFMDRawReader(fReader, 0); } // Optionally, get the geometry - if (TESTBIT(fTreeMask, kGeometry)) { - TString fname(fRun->GetGeometryFileName()); - if (fname.IsNull()) { - Warning("Init", "No file name for the geometry from AliRun"); + if (IsLoaded(kGeometry)) { + TString fname; + if (fRun) { fname = gSystem->DirName(fGAliceFile); fname.Append("/geometry.root"); } - fGeoManager = TGeoManager::Import(fname.Data()); - if (!fGeoManager) { - Fatal("Init", "No geometry manager found"); - return kFALSE; - } + if (!gSystem->AccessPathName(fname.Data())) + fname = ""; AliCDBManager* cdb = AliCDBManager::Instance(); + if (!cdb->IsDefaultStorageSet()) { + cdb->SetDefaultStorage("local://$ALICE_ROOT/OCDB"); + cdb->SetRun(0); + } + + AliGeomManager::LoadGeometry(fname.IsNull() ? 0 : fname.Data()); + AliCDBEntry* align = cdb->Get("FMD/Align/Data"); if (align) { AliInfo("Got alignment data from CDB"); @@ -233,17 +412,20 @@ AliFMDInput::Init() else { Int_t nAlign = array->GetEntries(); for (Int_t i = 0; i < nAlign; i++) { - AliAlignObjAngles* a = static_cast(array->At(i)); + AliAlignObjParams* a = static_cast(array->At(i)); if (!a->ApplyToGeometry()) { AliWarning(Form("Failed to apply alignment to %s", - a->GetVolPath())); + a->GetSymName())); } } } } + AliFMDGeometry* geom = AliFMDGeometry::Instance(); + geom->Init(); + geom->InitTransformations(); } - + fEventCount = 0; fIsInit = kTRUE; return fIsInit; } @@ -262,61 +444,109 @@ AliFMDInput::Begin(Int_t event) AliError("Not initialized"); return fIsInit; } + // Get the event - if (fLoader->GetEvent(event)) return kFALSE; - AliInfo(Form("Now in event %d/%d", event, NEvents())); + if (fLoader && fLoader->GetEvent(event)) return kFALSE; // Possibly load global kinematics information - if (TESTBIT(fTreeMask, kKinematics)) { - AliInfo("Getting kinematics"); - if (fLoader->LoadKinematics()) return kFALSE; + if (IsLoaded(kKinematics)) { + // AliInfo("Getting kinematics"); + if (fLoader->LoadKinematics("READ")) return kFALSE; fStack = fLoader->Stack(); } + // Possibly load FMD Hit information - if (TESTBIT(fTreeMask, kHits)) { - AliInfo("Getting FMD hits"); - if (fFMDLoader->LoadHits()) return kFALSE; + if (IsLoaded(kHits)) { + // AliInfo("Getting FMD hits"); + if (!fFMDLoader || fFMDLoader->LoadHits("READ")) return kFALSE; fTreeH = fFMDLoader->TreeH(); if (!fArrayH) fArrayH = fFMD->Hits(); } + + // Possibly load FMD TrackReference information + if (IsLoaded(kTrackRefs)) { + // AliInfo("Getting FMD hits"); + if (!fLoader || fLoader->LoadTrackRefs("READ")) return kFALSE; + fTreeTR = fLoader->TreeTR(); + if (!fArrayTR) fArrayTR = new TClonesArray("AliTrackReference"); + fTreeTR->SetBranchAddress("TrackReferences", &fArrayTR); + } + + // Possibly load heaedr information + if (IsLoaded(kHeader)) { + // AliInfo("Getting FMD hits"); + if (!fLoader /* || fLoader->LoadHeader()*/) return kFALSE; + fHeader = fLoader->GetHeader(); + } + // Possibly load FMD Digit information - if (TESTBIT(fTreeMask, kDigits)) { - AliInfo("Getting FMD digits"); - if (fFMDLoader->LoadDigits()) return kFALSE; + if (IsLoaded(kDigits)) { + // AliInfo("Getting FMD digits"); + if (!fFMDLoader || fFMDLoader->LoadDigits("READ")) return kFALSE; fTreeD = fFMDLoader->TreeD(); - if (!fArrayD) fArrayD = fFMD->Digits(); + if (fTreeD) { + if (!fArrayD) fArrayD = fFMD->Digits(); + } + else { + fArrayD = 0; + AliWarning(Form("Failed to load FMD Digits")); + } } + // Possibly load FMD Sdigit information - if (TESTBIT(fTreeMask, kSDigits)) { - AliInfo("Getting FMD summable digits"); - if (fFMDLoader->LoadSDigits()) return kFALSE; + if (IsLoaded(kSDigits)) { + // AliInfo("Getting FMD summable digits"); + if (!fFMDLoader || fFMDLoader->LoadSDigits("READ")) { + AliWarning("Failed to load SDigits!"); + return kFALSE; + } fTreeS = fFMDLoader->TreeS(); if (!fArrayS) fArrayS = fFMD->SDigits(); } + // Possibly load FMD RecPoints information - if (TESTBIT(fTreeMask, kRecPoints)) { - AliInfo("Getting FMD reconstructed points"); - if (fFMDLoader->LoadRecPoints()) return kFALSE; + if (IsLoaded(kRecPoints)) { + // AliInfo("Getting FMD reconstructed points"); + if (!fFMDLoader || fFMDLoader->LoadRecPoints("READ")) return kFALSE; fTreeR = fFMDLoader->TreeR(); if (!fArrayR) fArrayR = new TClonesArray("AliFMDRecPoint"); fTreeR->SetBranchAddress("FMD", &fArrayR); - } // Possibly load FMD ESD information - if (TESTBIT(fTreeMask, kESD)) { - AliInfo("Getting FMD event summary data"); + } + + // Possibly load FMD ESD information + if (IsLoaded(kESD)) { + // AliInfo("Getting FMD event summary data"); Int_t read = fChainE->GetEntry(event); if (read <= 0) return kFALSE; - fESD = fMainESD->GetFMDData(); + fESD = fESDEvent->GetFMDData(); if (!fESD) return kFALSE; } + // Possibly load FMD Digit information - if (TESTBIT(fTreeMask, kRaw)) { - AliInfo("Getting FMD raw data digits"); - if (!fReader->NextEvent()) return kFALSE; - AliFMDRawReader r(fReader, 0); + if (IsLoaded(kRaw) || IsLoaded(kRawCalib)) { + Bool_t mon = fRawFile.Contains("mem://"); + // AliInfo("Getting FMD raw data digits"); + if (mon) std::cout << "Waiting for event ..." << std::flush; + do { + if (!fReader->NextEvent()) { + if (mon) { + gSystem->Sleep(3); + continue; + } + return kFALSE; + } + UInt_t eventType = fReader->GetType(); + if(eventType == AliRawEventHeaderBase::kPhysicsEvent || + eventType == AliRawEventHeaderBase::kCalibrationEvent) + break; + } while (true); + if (mon) std::cout << "got it" << std::endl; + // AliFMDRawReader r(fReader, 0); fArrayA->Clear(); - r.ReadAdcs(fArrayA); + fFMDReader->ReadAdcs(fArrayA); + AliFMDDebug(1, ("Got a total of %d digits", fArrayA->GetEntriesFast())); } - + fEventCount++; return kTRUE; } @@ -333,18 +563,18 @@ AliFMDInput::Event() // - ProcessRecPoints if the reconstructed points are loaded. // - ProcessESD if the event summary data is loaded // - if (TESTBIT(fTreeMask, kHits)) - if (!ProcessHits()) return kFALSE; - if (TESTBIT(fTreeMask, kDigits)) - if (!ProcessDigits()) return kFALSE; - if (TESTBIT(fTreeMask, kSDigits)) - if (!ProcessSDigits()) return kFALSE; - if (TESTBIT(fTreeMask, kRaw)) - if (!ProcessRawDigits()) return kFALSE; - if (TESTBIT(fTreeMask, kRecPoints)) - if (!ProcessRecPoints()) return kFALSE; - if (TESTBIT(fTreeMask, kESD)) - if (!ProcessESD(fESD)) return kFALSE; + if (IsLoaded(kHits)) if (!ProcessHits()) return kFALSE; + if (IsLoaded(kTrackRefs))if (!ProcessTrackRefs()) return kFALSE; + if (IsLoaded(kKinematics) && + IsLoaded(kHits)) if (!ProcessTracks()) return kFALSE; + if (IsLoaded(kKinematics))if (!ProcessStack()) return kFALSE; + if (IsLoaded(kSDigits)) if (!ProcessSDigits()) return kFALSE; + if (IsLoaded(kDigits)) if (!ProcessDigits()) return kFALSE; + if (IsLoaded(kRaw)) if (!ProcessRawDigits()) return kFALSE; + if (IsLoaded(kRawCalib)) if (!ProcessRawCalibDigits())return kFALSE; + if (IsLoaded(kRecPoints))if (!ProcessRecPoints()) return kFALSE; + if (IsLoaded(kESD)) if (!ProcessESDs()) return kFALSE; + if (IsLoaded(kUser)) if (!ProcessUsers()) return kFALSE; return kTRUE; } @@ -359,21 +589,25 @@ AliFMDInput::ProcessHits() AliError("No hit tree defined"); return kFALSE; } + if (!fArrayH) { + AliError("No hit array defined"); + return kFALSE; + } + Int_t nTracks = fTreeH->GetEntries(); for (Int_t i = 0; i < nTracks; i++) { Int_t hitRead = fTreeH->GetEntry(i); if (hitRead <= 0) continue; - if (!fArrayH) { - AliError("No hit array defined"); - return kFALSE; - } + Int_t nHit = fArrayH->GetEntries(); if (nHit <= 0) continue; + for (Int_t j = 0; j < nHit; j++) { AliFMDHit* hit = static_cast(fArrayH->At(j)); if (!hit) continue; + TParticle* track = 0; - if (TESTBIT(fTreeMask, kKinematics) && fStack) { + if (IsLoaded(kKinematics) && fStack) { Int_t trackno = hit->Track(); track = fStack->Particle(trackno); } @@ -382,18 +616,125 @@ AliFMDInput::ProcessHits() } return kTRUE; } +//____________________________________________________________________ +Bool_t +AliFMDInput::ProcessTrackRefs() +{ + // Read the reconstrcted points tree, and pass each reconstruction + // object (AliFMDRecPoint) to either ProcessRecPoint. + if (!fTreeTR) { + AliError("No track reference tree defined"); + return kFALSE; + } + if (!fArrayTR) { + AliError("No track reference array defined"); + return kFALSE; + } + + Int_t nEv = fTreeTR->GetEntries(); + for (Int_t i = 0; i < nEv; i++) { + Int_t trRead = fTreeTR->GetEntry(i); + if (trRead <= 0) continue; + Int_t nTrackRefs = fArrayTR->GetEntries(); + for (Int_t j = 0; j < nTrackRefs; j++) { + AliTrackReference* trackRef = + static_cast(fArrayTR->At(j)); + if (!trackRef) continue; + // if (trackRef->DetectorId() != AliTrackReference::kFMD) continue; + TParticle* track = 0; + if (IsLoaded(kKinematics) && fStack) { + Int_t trackno = trackRef->GetTrack(); + track = fStack->Particle(trackno); + } + if (!ProcessTrackRef(trackRef,track)) return kFALSE; + } + } + return kTRUE; +} +//____________________________________________________________________ +Bool_t +AliFMDInput::ProcessTracks() +{ + // Read the hit tree, and pass each hit to the member function + // ProcessTrack. + if (!fStack) { + AliError("No track tree defined"); + return kFALSE; + } + if (!fTreeH) { + AliError("No hit tree defined"); + return kFALSE; + } + if (!fArrayH) { + AliError("No hit array defined"); + return kFALSE; + } + + // Int_t nTracks = fStack->GetNtrack(); + Int_t nTracks = fTreeH->GetEntries(); + for (Int_t i = 0; i < nTracks; i++) { + Int_t trackno = nTracks - i - 1; + TParticle* track = fStack->Particle(trackno); + if (!track) continue; + + // Get the hits for this track. + Int_t hitRead = fTreeH->GetEntry(i); + Int_t nHit = fArrayH->GetEntries(); + if (nHit == 0 || hitRead <= 0) { + // Let user code see the track, even if there's no hits. + if (!ProcessTrack(trackno, track, 0)) return kFALSE; + continue; + } + + // Loop over the hits corresponding to this track. + for (Int_t j = 0; j < nHit; j++) { + AliFMDHit* hit = static_cast(fArrayH->At(j)); + if (!ProcessTrack(trackno, track, hit)) return kFALSE; + } + } + return kTRUE; +} +//____________________________________________________________________ +Bool_t +AliFMDInput::ProcessStack() +{ + // Read the hit tree, and pass each hit to the member function + // ProcessTrack. + if (!fStack) { + AliError("No track tree defined"); + return kFALSE; + } + Int_t nTracks = fStack->GetNtrack(); + for (Int_t i = 0; i < nTracks; i++) { + Int_t trackno = nTracks - i - 1; + TParticle* track = fStack->Particle(trackno); + if (!track) continue; + if (!ProcessParticle(trackno, track)) return kFALSE; + } + return kTRUE; +} //____________________________________________________________________ Bool_t AliFMDInput::ProcessDigits() { // Read the digit tree, and pass each digit to the member function // ProcessDigit. + if (!fTreeD) { + AliError("No digit tree defined"); + return kFALSE; + } + if (!fArrayD) { + AliError("No digit array defined"); + return kFALSE; + } + Int_t nEv = fTreeD->GetEntries(); for (Int_t i = 0; i < nEv; i++) { Int_t digitRead = fTreeD->GetEntry(i); if (digitRead <= 0) continue; Int_t nDigit = fArrayD->GetEntries(); + AliFMDDebug(0, ("Got %5d digits for this event", nDigit)); if (nDigit <= 0) continue; for (Int_t j = 0; j < nDigit; j++) { AliFMDDigit* digit = static_cast(fArrayD->At(j)); @@ -410,11 +751,25 @@ AliFMDInput::ProcessSDigits() { // Read the summable digit tree, and pass each sumable digit to the // member function ProcessSdigit. - Int_t nEv = fTreeD->GetEntries(); + if (!fTreeS) { + AliWarning("No sdigit tree defined"); + return kTRUE; // Empty SDigits is fine + } + if (!fArrayS) { + AliWarning("No sdigit array defined"); + return kTRUE; // Empty SDigits is fine + } + + Int_t nEv = fTreeS->GetEntries(); for (Int_t i = 0; i < nEv; i++) { Int_t sdigitRead = fTreeS->GetEntry(i); - if (sdigitRead <= 0) continue; - Int_t nSdigit = fArrayS->GetEntries(); + if (sdigitRead <= 0) { + AliInfo(Form("Read nothing from tree")); + continue; + } + Int_t nSdigit = fArrayS->GetEntriesFast(); + AliFMDDebug(0, ("Got %5d digits for this event", nSdigit)); + AliInfo(Form("Got %5d digits for this event", nSdigit)); if (nSdigit <= 0) continue; for (Int_t j = 0; j < nSdigit; j++) { AliFMDSDigit* sdigit = static_cast(fArrayS->At(j)); @@ -431,22 +786,61 @@ AliFMDInput::ProcessRawDigits() { // Read the digit tree, and pass each digit to the member function // ProcessDigit. + if (!fArrayA) { + AliError("No raw digit array defined"); + return kFALSE; + } + Int_t nDigit = fArrayA->GetEntries(); if (nDigit <= 0) return kTRUE; for (Int_t j = 0; j < nDigit; j++) { AliFMDDigit* digit = static_cast(fArrayA->At(j)); if (!digit) continue; + if (AliLog::GetDebugLevel("FMD","") >= 40 && j < 30) + digit->Print(); if (!ProcessRawDigit(digit)) return kFALSE; } return kTRUE; } +//____________________________________________________________________ +Bool_t +AliFMDInput::ProcessRawCalibDigits() +{ + // Read the digit tree, and pass each digit to the member function + // ProcessDigit. + if (!fArrayA) { + AliError("No raw digit array defined"); + return kFALSE; + } + + Int_t nDigit = fArrayA->GetEntries(); + if (nDigit <= 0) return kTRUE; + for (Int_t j = 0; j < nDigit; j++) { + AliFMDDigit* digit = static_cast(fArrayA->At(j)); + if (!digit) continue; + if (AliLog::GetDebugLevel("FMD","") >= 40 && j < 30) + digit->Print(); + if (!ProcessRawCalibDigit(digit)) return kFALSE; + } + return kTRUE; +} + //____________________________________________________________________ Bool_t AliFMDInput::ProcessRecPoints() { // Read the reconstrcted points tree, and pass each reconstruction // object (AliFMDRecPoint) to either ProcessRecPoint. + if (!fTreeR) { + AliError("No recpoint tree defined"); + return kFALSE; + } + if (!fArrayR) { + AliError("No recpoints array defined"); + return kFALSE; + } + Int_t nEv = fTreeR->GetEntries(); for (Int_t i = 0; i < nEv; i++) { Int_t recRead = fTreeR->GetEntry(i); @@ -461,6 +855,52 @@ AliFMDInput::ProcessRecPoints() return kTRUE; } +//____________________________________________________________________ +Bool_t +AliFMDInput::ProcessESDs() +{ + // Process event summary data + if (!fESD) return kFALSE; + for (UShort_t det = 1; det <= 3; det++) { + Char_t rings[] = { 'I', (det == 1 ? '\0' : 'O'), '\0' }; + for (Char_t* rng = rings; *rng != '\0'; rng++) { + UShort_t nsec = (*rng == 'I' ? 20 : 40); + UShort_t nstr = (*rng == 'I' ? 512 : 256); + for (UShort_t sec = 0; sec < nsec; sec++) { + for (UShort_t str = 0; str < nstr; str++) { + Float_t eta = fESD->Eta(det,*rng,sec,str); + Float_t mult = fESD->Multiplicity(det,*rng,sec,str); + if (!fESD->IsAngleCorrected()) + mult *= TMath::Abs(TMath::Cos(2.*TMath::ATan(TMath::Exp(-eta)))); + if (!ProcessESD(det, *rng, sec, str, eta, mult)) continue; + } + } + } + } + return kTRUE; +} + +//____________________________________________________________________ +Bool_t +AliFMDInput::ProcessUsers() +{ + // Process event summary data + for (UShort_t det = 1; det <= 3; det++) { + Char_t rings[] = { 'I', (det == 1 ? '\0' : 'O'), '\0' }; + for (Char_t* rng = rings; *rng != '\0'; rng++) { + UShort_t nsec = (*rng == 'I' ? 20 : 40); + UShort_t nstr = (*rng == 'I' ? 512 : 256); + for (UShort_t sec = 0; sec < nsec; sec++) { + for (UShort_t str = 0; str < nstr; str++) { + Float_t v = GetSignal(det,*rng,sec,str); + if (!ProcessUser(det, *rng, sec, str, v)) continue; + } + } + } + } + return kTRUE; +} + //____________________________________________________________________ Bool_t AliFMDInput::End() @@ -476,55 +916,164 @@ AliFMDInput::End() return fIsInit; } // Possibly unload global kinematics information - if (TESTBIT(fTreeMask, kKinematics)) { + if (IsLoaded(kKinematics)) { fLoader->UnloadKinematics(); // fTreeK = 0; fStack = 0; } // Possibly unload FMD Hit information - if (TESTBIT(fTreeMask, kHits)) { + if (IsLoaded(kHits)) { fFMDLoader->UnloadHits(); fTreeH = 0; } // Possibly unload FMD Digit information - if (TESTBIT(fTreeMask, kDigits)) { + if (IsLoaded(kDigits)) { fFMDLoader->UnloadDigits(); fTreeD = 0; } // Possibly unload FMD Sdigit information - if (TESTBIT(fTreeMask, kSDigits)) { + if (IsLoaded(kSDigits)) { fFMDLoader->UnloadSDigits(); fTreeS = 0; } // Possibly unload FMD RecPoints information - if (TESTBIT(fTreeMask, kRecPoints)) { + if (IsLoaded(kRecPoints)) { fFMDLoader->UnloadRecPoints(); fTreeR = 0; } - AliInfo("Now out event"); + // AliInfo("Now out event"); return kTRUE; } //____________________________________________________________________ Bool_t -AliFMDInput::Run() +AliFMDInput::Run(UInt_t maxEvents) { // Run over all events and files references in galice.root Bool_t retval; if (!(retval = Init())) return retval; - Int_t nEvents = NEvents(); - for (Int_t event = 0; event < nEvents; event++) { + fNEvents = NEvents(); + if (fNEvents < 0) fNEvents = maxEvents; + else if (maxEvents > 0) fNEvents = TMath::Min(fNEvents,Int_t(maxEvents)); + + Int_t event = 0; + for (; fNEvents < 0 || event < fNEvents; event++) { + printf("\rEvent %8d/%8d ...", event, fNEvents); if (!(retval = Begin(event))) break; if (!(retval = Event())) break; if (!(retval = End())) break; } + printf("Looped over %8d events\n", event+1); if (!retval) return retval; retval = Finish(); return retval; } +//__________________________________________________________________ +TArrayF +AliFMDInput::MakeLogScale(Int_t n, Double_t min, Double_t max) +{ + // Service function to define a logarithmic axis. + // Parameters: + // n Number of bins + // min Minimum of axis + // max Maximum of axis + TArrayF bins(n+1); + bins[0] = min; + if (n <= 20) { + for (Int_t i = 1; i < n+1; i++) bins[i] = bins[i-1] + (max-min)/n; + return bins; + } + Float_t dp = n / TMath::Log10(max / min); + Float_t pmin = TMath::Log10(min); + for (Int_t i = 1; i < n+1; i++) { + Float_t p = pmin + i / dp; + bins[i] = TMath::Power(10, p); + } + return bins; +} + +//____________________________________________________________________ +void +AliFMDInput::ScanDirectory(TSystemDirectory* dir, + const TString& olddir, + TChain* chain, + const char* pattern, bool recursive) +{ + // Get list of files, and go back to old working directory + TString oldDir(gSystem->WorkingDirectory()); + TList* files = dir->GetListOfFiles(); + gSystem->ChangeDirectory(oldDir); + + // Sort list of files and check if we should add it + if (!files) return; + files->Sort(); + TIter next(files); + TSystemFile* file = 0; + while ((file = static_cast(next()))) { + TString name(file->GetName()); + + // Ignore special links + if (name == "." || name == "..") continue; + + // Check if this is a directory + if (file->IsDirectory()) { + if (recursive) + ScanDirectory(static_cast(file), + olddir, chain, + pattern,recursive); + continue; + } + + // If this is not a root file, ignore + if (!name.EndsWith(".root")) continue; + + // If this file does not contain the pattern, ignore + if (!name.Contains(pattern)) continue; + if (name.Contains("friends")) continue; + + // Get the path + TString data(Form("%s/%s", file->GetTitle(), name.Data())); + + TFile* test = TFile::Open(data.Data(), "READ"); + if (!test || test->IsZombie()) { + ::Warning("ScanDirectory", "Failed to open file %s", data.Data()); + continue; + } + test->Close(); + chain->Add(data); + } +} + +//____________________________________________________________________ +TChain* +AliFMDInput::MakeChain(const char* what, const char* datadir, bool recursive) +{ + TString w(what); + w.ToUpper(); + const char* treeName = 0; + const char* pattern = 0; + if (w.Contains("ESD")) { treeName = "esdTree"; pattern = "AliESD"; } + else if (w.Contains("MC")) { treeName = "TE"; pattern = "galice"; } + else { + ::Error("MakeChain", "Unknown mode '%s' (not one of ESD, or MC)", what); + return 0; + } + + // --- Our data chain ---------------------------------------------- + TChain* chain = new TChain(treeName); + + // --- Get list of ESDs -------------------------------------------- + // Open source directory, and make sure we go back to were we were + TString oldDir(gSystem->WorkingDirectory()); + TSystemDirectory d(datadir, datadir); + ScanDirectory(&d, oldDir, chain, pattern, recursive); + + return chain; +} + //____________________________________________________________________ //