AliRoot Core  edcc906 (edcc906)
AliFMDInput.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 2004, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 /* $Id$ */
21 //___________________________________________________________________
22 //
23 // The classes defined here, are utility classes for reading in data
24 // for the FMD. They are put in a seperate library to not polute the
25 // normal libraries. The classes are intended to be used as base
26 // classes for customized class that do some sort of analysis on the
27 // various types of data produced by the FMD.
28 //
29 // Latest changes by Christian Holm Christensen
30 //
31 #include "AliFMDInput.h" // ALIFMDHIT_H
32 #include "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
47 #include "AliFMDGeometry.h"
48 #include <AliESD.h>
49 #include <AliESDFMD.h>
50 #include <AliESDEvent.h>
51 #include <AliCDBManager.h>
52 #include <AliCDBEntry.h>
53 #include <AliAlignObjParams.h>
54 #include <AliTrackReference.h>
55 #include <TTree.h> // ROOT_TTree
56 #include <TChain.h> // ROOT_TChain
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
62 #include <TSystemDirectory.h> // ROOT_TSystemDirectory
63 #include <Riostream.h> // ROOT_Riostream
64 #include <TFile.h> // ROOT_TFile
65 #include <TStreamerInfo.h>
66 #include <TArrayF.h>
67 
68 //____________________________________________________________________
69 ClassImp(AliFMDInput)
70 #if 0
71  ; // This is here to keep Emacs for indenting the next line
72 #endif
73 
74 //____________________________________________________________________
76  kKinematics,
77  kDigits,
78  kSDigits,
79  kHeader,
80  kRecPoints,
81  kESD,
82  kRaw,
83  kGeometry,
84  kTrackRefs,
85  kRawCalib,
86  kUser };
87 
88 //____________________________________________________________________
90  : TNamed("AliFMDInput", "Input handler for various FMD data"),
91  fGAliceFile(""),
92  fLoader(0),
93  fRun(0),
94  fStack(0),
95  fFMDLoader(0),
96  fReader(0),
97  fFMDReader(0),
98  fFMD(0),
99  fESD(0),
100  fESDEvent(0),
101  fTreeE(0),
102  fTreeH(0),
103  fTreeTR(0),
104  fTreeD(0),
105  fTreeS(0),
106  fTreeR(0),
107  fTreeA(0),
108  fChainE(0),
109  fArrayE(0),
110  fArrayH(0),
111  fArrayTR(0),
112  fArrayD(0),
113  fArrayS(0),
114  fArrayR(0),
115  fArrayA(0),
116  fHeader(0),
117  fGeoManager(0),
118  fTreeMask(0),
119  fRawFile(""),
120  fInputDir("."),
121  fIsInit(kFALSE),
122  fEventCount(0),
123  fNEvents(-1)
124 {
125 
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 //____________________________________________________________________
135 AliFMDInput::AliFMDInput(const char* gAliceFile)
136  : TNamed("AliFMDInput", "Input handler for various FMD data"),
137  fGAliceFile(gAliceFile),
138  fLoader(0),
139  fRun(0),
140  fStack(0),
141  fFMDLoader(0),
142  fReader(0),
143  fFMDReader(0),
144  fFMD(0),
145  fESD(0),
146  fESDEvent(0),
147  fTreeE(0),
148  fTreeH(0),
149  fTreeTR(0),
150  fTreeD(0),
151  fTreeS(0),
152  fTreeR(0),
153  fTreeA(0),
154  fChainE(0),
155  fArrayE(0),
156  fArrayH(0),
157  fArrayTR(0),
158  fArrayD(0),
159  fArrayS(0),
160  fArrayR(0),
161  fArrayA(0),
162  fHeader(0),
163  fGeoManager(0),
164  fTreeMask(0),
165  fRawFile(""),
166  fInputDir("."),
167  fIsInit(kFALSE),
168  fEventCount(0),
169  fNEvents(-1)
170 {
171 
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 //____________________________________________________________________
179 void
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 //____________________________________________________________________
195 void
196 AliFMDInput::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  delete ll;
207 }
208 
209 
210 //____________________________________________________________________
212 AliFMDInput::ParseLoad(const char* what)
213 {
214  TString opt(what);
215  opt.ToLower();
216  const ETrees* ptype = fgkAllLoads;
217  do {
218  ETrees type = *ptype;
219  if (opt.Contains(TreeName(type,true), TString::kIgnoreCase))
220  return type;
221  } while (*ptype++ != kUser);
222  return kUser;
223 }
224 //____________________________________________________________________
225 const char*
226 AliFMDInput::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 //____________________________________________________________________
248 const char*
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  }
282  return 0;
283 }
284 
285 //____________________________________________________________________
286 Int_t
288 {
289  // Get number of events
290  if (IsLoaded(kRaw) ||
291  IsLoaded(kRawCalib)) return fReader->GetNumberOfEvents();
292  if (fChainE) return fChainE->GetEntriesFast();
293  if (fTreeE) return fTreeE->GetEntries();
294  return -1;
295 }
296 
297 //____________________________________________________________________
298 Bool_t
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  }
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"));
314  } while (*ptype++ != kUser);
315 
316  Info("Init","Initialising w/mask 0x%04x%s", fTreeMask, what.Data());
317  // Get the run
318  if (IsLoaded(kDigits) ||
319  IsLoaded(kSDigits) ||
320  IsLoaded(kKinematics) ||
321  IsLoaded(kTrackRefs) ||
322  IsLoaded(kHeader)) {
323  if (!gSystem->FindFile(".:/", fGAliceFile)) {
324  AliWarning(Form("Cannot find file %s in .:/", fGAliceFile.Data()));
325  }
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();
356  }
357  }
358 
359  // Optionally, get the ESD files
360  if (IsLoaded(kESD)) {
361  fChainE = MakeChain("ESD", fInputDir, true);
362  fESDEvent = new AliESDEvent();
364  // fChainE->SetBranchAddress("ESD", &fMainESD);
365 
366  }
367 
368  if (IsLoaded(kRaw) ||
369  IsLoaded(kRawCalib)) {
370  AliInfo("Getting FMD raw data digits");
371  fArrayA = new TClonesArray("AliFMDDigit");
372 #if 0
373  if (!fRawFile.IsNull() && fRawFile.EndsWith(".root"))
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);
379 #else
380  if(!fRawFile.IsNull())
381  fReader = AliRawReader::Create(fRawFile.Data());
382  else
383  fReader = new AliRawReaderFile(-1);
384 #endif
386  }
387 
388  // Optionally, get the geometry
389  if (IsLoaded(kGeometry)) {
390  TString fname;
391  if (fRun) {
392  fname = gSystem->DirName(fGAliceFile);
393  fname.Append("/geometry.root");
394  }
395  if (!gSystem->AccessPathName(fname.Data()))
396  fname = "";
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 
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++) {
415  AliAlignObjParams* a = static_cast<AliAlignObjParams*>(array->At(i));
416  if (!a->ApplyToGeometry()) {
417  AliWarning(Form("Failed to apply alignment to %s",
418  a->GetSymName()));
419  }
420  }
421  }
422  }
424  geom->Init();
425  geom->InitTransformations();
426  }
427 
428  fEventCount = 0;
429  fIsInit = kTRUE;
430  return fIsInit;
431 }
432 
433 //____________________________________________________________________
434 Bool_t
435 AliFMDInput::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  }
447 
448  // Get the event
449  if (fLoader && fLoader->GetEvent(event)) return kFALSE;
450 
451  // Possibly load global kinematics information
452  if (IsLoaded(kKinematics)) {
453  // AliInfo("Getting kinematics");
454  if (fLoader->LoadKinematics("READ")) return kFALSE;
455  fStack = fLoader->Stack();
456  }
457 
458  // Possibly load FMD Hit information
459  if (IsLoaded(kHits)) {
460  // AliInfo("Getting FMD hits");
461  if (!fFMDLoader || fFMDLoader->LoadHits("READ")) return kFALSE;
462  fTreeH = fFMDLoader->TreeH();
463  if (!fArrayH) fArrayH = fFMD->Hits();
464  }
465 
466  // Possibly load FMD TrackReference information
467  if (IsLoaded(kTrackRefs)) {
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 
475  // Possibly load heaedr information
476  if (IsLoaded(kHeader)) {
477  // AliInfo("Getting FMD hits");
478  if (!fLoader /* || fLoader->LoadHeader()*/) return kFALSE;
480  }
481 
482  // Possibly load FMD Digit information
483  if (IsLoaded(kDigits)) {
484  // AliInfo("Getting FMD digits");
485  if (!fFMDLoader || fFMDLoader->LoadDigits("READ")) return kFALSE;
486  fTreeD = fFMDLoader->TreeD();
487  if (fTreeD) {
488  if (!fArrayD) fArrayD = fFMD->Digits();
489  }
490  else {
491  fArrayD = 0;
492  AliWarning(Form("Failed to load FMD Digits"));
493  }
494  }
495 
496  // Possibly load FMD Sdigit information
497  if (IsLoaded(kSDigits)) {
498  // AliInfo("Getting FMD summable digits");
499  if (!fFMDLoader || fFMDLoader->LoadSDigits("READ")) {
500  AliWarning("Failed to load SDigits!");
501  return kFALSE;
502  }
503  fTreeS = fFMDLoader->TreeS();
504  if (!fArrayS) fArrayS = fFMD->SDigits();
505  }
506 
507  // Possibly load FMD RecPoints information
508  if (IsLoaded(kRecPoints)) {
509  // AliInfo("Getting FMD reconstructed points");
510  if (!fFMDLoader || fFMDLoader->LoadRecPoints("READ")) return kFALSE;
511  fTreeR = fFMDLoader->TreeR();
512  if (!fArrayR) fArrayR = new TClonesArray("AliFMDRecPoint");
513  fTreeR->SetBranchAddress("FMD", &fArrayR);
514  }
515 
516  // Possibly load FMD ESD information
517  if (IsLoaded(kESD)) {
518  // AliInfo("Getting FMD event summary data");
519  Int_t read = fChainE->GetEntry(event);
520  if (read <= 0) return kFALSE;
522  if (!fESD) return kFALSE;
523  }
524 
525  // Possibly load FMD Digit information
526  if (IsLoaded(kRaw) || IsLoaded(kRawCalib)) {
527  Bool_t mon = fRawFile.Contains("mem://");
528  // AliInfo("Getting FMD raw data digits");
529  if (mon) std::cout << "Waiting for event ..." << std::flush;
530  do {
531  if (!fReader->NextEvent()) {
532  if (mon) {
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);
543  if (mon) std::cout << "got it" << std::endl;
544  // AliFMDRawReader r(fReader, 0);
545  fArrayA->Clear();
547  AliFMDDebug(1, ("Got a total of %d digits", fArrayA->GetEntriesFast()));
548  }
549  fEventCount++;
550  return kTRUE;
551 }
552 
553 
554 //____________________________________________________________________
555 Bool_t
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  //
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;
578 
579  return kTRUE;
580 }
581 
582 //____________________________________________________________________
583 Bool_t
585 {
586  // Read the hit tree, and pass each hit to the member function
587  // ProcessHit.
588  if (!fTreeH) {
589  AliError("No hit tree defined");
590  return kFALSE;
591  }
592  if (!fArrayH) {
593  AliError("No hit array defined");
594  return kFALSE;
595  }
596 
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;
601 
602  Int_t nHit = fArrayH->GetEntries();
603  if (nHit <= 0) continue;
604 
605  for (Int_t j = 0; j < nHit; j++) {
606  AliFMDHit* hit = static_cast<AliFMDHit*>(fArrayH->At(j));
607  if (!hit) continue;
608 
609  TParticle* track = 0;
610  if (IsLoaded(kKinematics) && fStack) {
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 }
619 //____________________________________________________________________
620 Bool_t
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  }
633 
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++) {
640  AliTrackReference* trackRef =
641  static_cast<AliTrackReference*>(fArrayTR->At(j));
642  if (!trackRef) continue;
643  // if (trackRef->DetectorId() != AliTrackReference::kFMD) continue;
644  TParticle* track = 0;
645  if (IsLoaded(kKinematics) && fStack) {
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 }
654 //____________________________________________________________________
655 Bool_t
657 {
658  // Read the hit tree, and pass each hit to the member function
659  // ProcessTrack.
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  }
668  if (!fArrayH) {
669  AliError("No hit array defined");
670  return kFALSE;
671  }
672 
673  // Int_t nTracks = fStack->GetNtrack();
674  Int_t nTracks = fTreeH->GetEntries();
675  for (Int_t i = 0; i < nTracks; i++) {
676  Int_t trackno = nTracks - i - 1;
677  TParticle* track = fStack->Particle(trackno);
678  if (!track) continue;
679 
680  // Get the hits for this track.
681  Int_t hitRead = fTreeH->GetEntry(i);
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;
687  }
688 
689  // Loop over the hits corresponding to this track.
690  for (Int_t j = 0; j < nHit; j++) {
691  AliFMDHit* hit = static_cast<AliFMDHit*>(fArrayH->At(j));
692  if (!ProcessTrack(trackno, track, hit)) return kFALSE;
693  }
694  }
695  return kTRUE;
696 }
697 //____________________________________________________________________
698 Bool_t
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;
712 
713  if (!ProcessParticle(trackno, track)) return kFALSE;
714  }
715  return kTRUE;
716 }
717 //____________________________________________________________________
718 Bool_t
720 {
721  // Read the digit tree, and pass each digit to the member function
722  // ProcessDigit.
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 
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();
737  AliFMDDebug(0, ("Got %5d digits for this event", nDigit));
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 //____________________________________________________________________
749 Bool_t
751 {
752  // Read the summable digit tree, and pass each sumable digit to the
753  // member function ProcessSdigit.
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();
764  for (Int_t i = 0; i < nEv; i++) {
765  Int_t sdigitRead = fTreeS->GetEntry(i);
766  if (sdigitRead <= 0) {
767  AliInfo(Form("Read nothing from tree"));
768  continue;
769  }
770  Int_t nSdigit = fArrayS->GetEntriesFast();
771  AliFMDDebug(0, ("Got %5d digits for this event", nSdigit));
772  AliInfo(Form("Got %5d digits for this event", nSdigit));
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 
783 //____________________________________________________________________
784 Bool_t
786 {
787  // Read the digit tree, and pass each digit to the member function
788  // ProcessDigit.
789  if (!fArrayA) {
790  AliError("No raw digit array defined");
791  return kFALSE;
792  }
793 
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;
799  if (AliLog::GetDebugLevel("FMD","") >= 40 && j < 30)
800  digit->Print();
801  if (!ProcessRawDigit(digit)) return kFALSE;
802  }
803  return kTRUE;
804 }
805 
806 //____________________________________________________________________
807 Bool_t
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 
829 //____________________________________________________________________
830 Bool_t
832 {
833  // Read the reconstrcted points tree, and pass each reconstruction
834  // object (AliFMDRecPoint) to either ProcessRecPoint.
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 
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  }
854  }
855  return kTRUE;
856 }
857 
858 //____________________________________________________________________
859 Bool_t
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);
873  if (!fESD->IsAngleCorrected())
874  mult *= TMath::Abs(TMath::Cos(2.*TMath::ATan(TMath::Exp(-eta))));
875  if (!ProcessESD(det, *rng, sec, str, eta, mult)) continue;
876  }
877  }
878  }
879  }
880  return kTRUE;
881 }
882 
883 //____________________________________________________________________
884 Bool_t
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 
904 //____________________________________________________________________
905 Bool_t
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
919  if (IsLoaded(kKinematics)) {
921  // fTreeK = 0;
922  fStack = 0;
923  }
924  // Possibly unload FMD Hit information
925  if (IsLoaded(kHits)) {
927  fTreeH = 0;
928  }
929  // Possibly unload FMD Digit information
930  if (IsLoaded(kDigits)) {
932  fTreeD = 0;
933  }
934  // Possibly unload FMD Sdigit information
935  if (IsLoaded(kSDigits)) {
937  fTreeS = 0;
938  }
939  // Possibly unload FMD RecPoints information
940  if (IsLoaded(kRecPoints)) {
942  fTreeR = 0;
943  }
944  // AliInfo("Now out event");
945  return kTRUE;
946 }
947 
948 //____________________________________________________________________
949 Bool_t
950 AliFMDInput::Run(UInt_t maxEvents)
951 {
952  // Run over all events and files references in galice.root
953 
954  Bool_t retval;
955  if (!(retval = Init())) return retval;
956 
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);
964  if (!(retval = Begin(event))) break;
965  if (!(retval = Event())) break;
966  if (!(retval = End())) break;
967  }
968  printf("Looped over %8d events\n", event+1);
969  if (!retval) return retval;
970  retval = Finish();
971  return retval;
972 }
973 
974 //__________________________________________________________________
975 TArrayF
976 AliFMDInput::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 
998 //____________________________________________________________________
999 void
1000 AliFMDInput::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 //____________________________________________________________________
1051 TChain*
1052 AliFMDInput::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 }
1076 
1077 
1078 //____________________________________________________________________
1079 //
1080 // EOF
1081 //
class for digits
Definition: AliFMDDigit.h:28
virtual Bool_t ProcessSDigit(AliFMDSDigit *sdigit)
Definition: AliFMDInput.h:457
Reconstructed FMD points. It contains the pseudo-inclusive multiplicity.
Int_t LoadgAlice()
printf("Chi2/npoints = %f\n", TMath::Sqrt(chi2/npoints))
virtual Bool_t ProcessTracks()
TFile * Open(const char *filename, Long64_t &nevents)
Geometry mananger for the FMD.
void Print(Option_t *opt="") const
void UnloadRecPoints() const
Definition: AliLoader.h:132
TClonesArray * Hits() const
Definition: AliDetector.h:35
TClonesArray * Digits() const
Definition: AliDetector.h:34
virtual const char * LoadedString(Bool_t dataOnly=false) const
Hit in the FMD.
virtual Float_t GetSignal(UShort_t d, Char_t r, UShort_t s, UShort_t t)
Definition: AliFMDInput.h:465
virtual Bool_t Finish()
Definition: AliFMDInput.h:192
virtual Bool_t ProcessTrack(Int_t i, TParticle *p, AliFMDHit *h)
Definition: AliFMDInput.h:453
virtual Bool_t ProcessRecPoints()
#define TObjArray
AliTPCcalibAlign align
Definition: CalibAlign.C:43
virtual void InitTransformations(Bool_t force=kFALSE)
TTree * fTreeR
Definition: AliFMDInput.h:428
virtual Bool_t ProcessHits()
virtual Bool_t IsLoaded(ETrees tree) const
Definition: AliFMDInput.h:141
Digits for the FMD.
virtual void SetLoads(UInt_t mask)
FMD utility classes for reading FMD data.
TGeoManager * fGeoManager
Definition: AliFMDInput.h:439
static const char * TreeName(ETrees tree, bool shortest=false)
Class to read ALTRO formated raw data from an AliRawReader object.
static const ETrees fgkAllLoads[kUser+1]
Definition: AliFMDInput.h:446
TTree * fTreeD
Definition: AliFMDInput.h:426
TTree * TreeD() const
Definition: AliLoader.h:87
Float_t p[]
Definition: kNNTest.C:133
virtual Bool_t Begin(Int_t event)
AliLoader * GetLoader(const char *detname) const
Bool_t ApplyToGeometry(Bool_t ovlpcheck=kFALSE)
AliLoader * fFMDLoader
Definition: AliFMDInput.h:417
virtual Bool_t ProcessHit(AliFMDHit *h, TParticle *p)
Definition: AliFMDInput.h:450
Bool_t IsAngleCorrected() const
Definition: AliESDFMD.h:219
virtual Bool_t Event()
AliRunLoader * fLoader
Definition: AliFMDInput.h:414
Forward Multiplicity Detector based on Silicon wafers. This class is the driver for especially simula...
Definition: AliFMD.h:306
TTree * TreeS() const
Definition: AliLoader.h:85
TChain * chain
TTree * fTreeE
Definition: AliFMDInput.h:423
static ETrees ParseLoad(const char *what)
AliTPCfastTrack * track
Int_t LoadRecPoints(Option_t *opt="")
Definition: AliLoader.h:113
TChain * fChainE
Definition: AliFMDInput.h:430
TTree * tree
TString fInputDir
Definition: AliFMDInput.h:442
AliHeader * fHeader
Definition: AliFMDInput.h:438
AliFMDRawReader * fFMDReader
Definition: AliFMDInput.h:419
TClonesArray * fArrayR
Definition: AliFMDInput.h:436
Pseudo reconstructed charged particle multiplicity.
#define AliWarning(message)
Definition: AliLog.h:541
TTree * fTreeH
Definition: AliFMDInput.h:424
virtual Bool_t ProcessDigit(AliFMDDigit *digit)
Definition: AliFMDInput.h:456
Int_t Track() const
Definition: AliHit.h:24
TObjArray * array
Definition: AnalyzeLaser.C:12
AliCDBEntry * Get(const AliCDBId &query, Bool_t forceCaching=kFALSE)
static AliRunLoader * Open(const char *filename="galice.root", const char *eventfoldername=AliConfig::GetDefaultEventFolderName(), Option_t *option="READ")
TTree * TreeR() const
Definition: AliLoader.h:89
virtual Bool_t ProcessUsers()
virtual Bool_t ProcessStack()
virtual Int_t GetTrack() const
virtual Bool_t Init()
AliESDFMD * GetFMDData() const
Definition: AliESDEvent.h:252
Int_t LoadSDigits(Option_t *opt="")
Definition: AliLoader.h:101
virtual Bool_t ProcessRawDigit(AliFMDDigit *digit)
Definition: AliFMDInput.h:458
virtual Int_t GetNtrack() const
Definition: AliStack.h:134
virtual Bool_t ProcessRawCalibDigit(AliFMDDigit *digit)
Definition: AliFMDInput.h:459
void UnloadDigits() const
Definition: AliLoader.h:131
virtual Bool_t ProcessRawCalibDigits()
virtual Int_t NEvents() const
Int_t fTreeMask
Definition: AliFMDInput.h:440
AliRun * fRun
Definition: AliFMDInput.h:415
Singleton object of FMD geometry descriptions and parameters. This class is a singleton that handles ...
AliHeader * GetHeader() const
#define AliInfo(message)
Definition: AliLog.h:484
virtual TClonesArray * SDigits()
Definition: AliFMD.h:427
Float_t Multiplicity(UShort_t detector, Char_t ring, UShort_t sector, UShort_t strip) const
Definition: AliESDFMD.cxx:185
Declaration of AliFMD detector driver.
virtual Bool_t ProcessESD(UShort_t d, Char_t r, UShort_t s, UShort_t t, Float_t eta, Float_t mult)
Definition: AliFMDInput.h:461
TTree * fTreeS
Definition: AliFMDInput.h:427
AliESDEvent * fESDEvent
Definition: AliFMDInput.h:422
TString fRawFile
Definition: AliFMDInput.h:441
static TArrayF MakeLogScale(Int_t n, Double_t min, Double_t max)
Float_t Eta(UShort_t detector, Char_t ring, UShort_t sector, UShort_t strip) const
Definition: AliESDFMD.cxx:198
TClonesArray * fArrayA
Definition: AliFMDInput.h:437
TClonesArray * fArrayE
Definition: AliFMDInput.h:431
void UnloadKinematics()
Bool_t fIsInit
Definition: AliFMDInput.h:443
virtual Bool_t ProcessDigits()
virtual Bool_t ProcessTrackRefs()
#define AliFMDDebug(N, A)
Definition: AliFMDDebug.h:39
TTree * fTreeTR
Definition: AliFMDInput.h:425
virtual Bool_t ReadAdcs(TClonesArray *array)
virtual void AddLoad(ETrees tree)
Definition: AliFMDInput.h:134
Int_t LoadHeader()
virtual Bool_t ProcessSDigits()
static void ScanDirectory(TSystemDirectory *dir, const TString &olddir, TChain *chain, const char *pattern, bool recursive)
AliRawReader * fReader
Definition: AliFMDInput.h:418
TClonesArray * fArrayS
Definition: AliFMDInput.h:435
TString fGAliceFile
Definition: AliFMDInput.h:413
virtual void Init()
void SetRun(Int_t run)
Class to read raw data.
Definition: AliCDBEntry.h:18
AliStack * fStack
Definition: AliFMDInput.h:416
TTree * TreeTR() const
TTree * TreeE() const
Int_t LoadHits(Option_t *opt="")
Definition: AliLoader.h:96
TClonesArray * fArrayTR
Definition: AliFMDInput.h:433
AliFMD * fFMD
Definition: AliFMDInput.h:420
void SetDefaultStorage(const char *dbString)
Digits for the FMD.
virtual Bool_t ProcessRecPoint(AliFMDRecPoint *point)
Definition: AliFMDInput.h:460
Int_t GetEvent(Int_t evno)
void test()
Definition: interpolTest.C:100
TClonesArray * fArrayD
Definition: AliFMDInput.h:434
Base class for reading in various FMD data. The class loops over all found events. For each event the specified data is read in. The class then loops over all elements of the read data, and process these with user defined code.
Definition: AliFMDInput.h:106
TParticle * Particle(Int_t id, Bool_t useInEmbedding=kFALSE)
Definition: AliStack.cxx:688
virtual Bool_t End()
virtual Bool_t ProcessESDs()
Int_t fNEvents
Definition: AliFMDInput.h:445
void ReadFromTree(TTree *tree, Option_t *opt="")
void UnloadHits() const
Definition: AliLoader.h:129
TTree * TreeH() const
Definition: AliLoader.h:83
virtual Bool_t Run(UInt_t maxEvents=0)
AliESDFMD * fESD
Definition: AliFMDInput.h:421
class for summable digits
Definition: AliFMDSDigit.h:27
Int_t LoadTrackRefs(Option_t *option="READ")
static TChain * MakeChain(const char *what, const char *datadir, bool recursive=false)
void UnloadSDigits() const
Definition: AliLoader.h:130
AliFMDhit is the hit class for the FMD. Hits are the information that comes from a Monte Carlo at eac...
Definition: AliFMDHit.h:30
virtual Bool_t ProcessRawDigits()
TTree * fTreeA
Definition: AliFMDInput.h:429
const char * GetSymName() const
Definition: AliAlignObj.h:63
static Int_t GetDebugLevel(const char *module, const char *className)
Definition: AliLog.cxx:843
AliDetector * GetDetector(const char *name) const
Definition: AliRun.cxx:200
Int_t fEventCount
Definition: AliFMDInput.h:444
#define AliError(message)
Definition: AliLog.h:591
virtual Bool_t ProcessTrackRef(AliTrackReference *trackRef, TParticle *track)
Definition: AliFMDInput.h:451
Int_t LoadDigits(Option_t *opt="")
Definition: AliLoader.h:106
static AliCDBManager * Instance(TMap *entryCache=NULL, Int_t run=-1)
TEveGeoShape * geom
Definition: tpc_tracks.C:10
static AliFMDGeometry * Instance()
Bool_t IsDefaultStorageSet() const
Definition: AliCDBManager.h:60
char * fname
static void LoadGeometry(const char *geomFileName=NULL)
virtual Bool_t ProcessParticle(Int_t i, TParticle *p)
Definition: AliFMDInput.h:455
virtual Bool_t ProcessUser(UShort_t d, Char_t r, UShort_t s, UShort_t t, Float_t v)
Definition: AliFMDInput.h:463
Int_t LoadKinematics(Option_t *option="READ")
AliRun * GetAliRun() const
AliStack * Stack() const
Definition: AliRunLoader.h:95
TClonesArray * fArrayH
Definition: AliFMDInput.h:432