c22fb594838d0c9ef44122f1aee512ab25b4f171
[u/mrichter/AliRoot.git] / FMD / AliFMDInput.cxx
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
16 /* $Id$ */
17
18 //___________________________________________________________________
19 //
20 // The classes defined here, are utility classes for reading in data
21 // for the FMD.  They are  put in a seperate library to not polute the
22 // normal libraries.  The classes are intended to be used as base
23 // classes for customized class that do some sort of analysis on the
24 // various types of data produced by the FMD. 
25 //
26 // Latest changes by Christian Holm Christensen
27 //
28 #include "AliFMDInput.h"        // ALIFMDHIT_H
29 #include "AliLog.h"             // ALILOG_H
30 #include "AliLoader.h"          // ALILOADER_H
31 #include "AliRunLoader.h"       // ALIRUNLOADER_H
32 #include "AliRun.h"             // ALIRUN_H
33 #include "AliStack.h"           // ALISTACK_H
34 #include "AliFMD.h"             // ALIFMD_H
35 #include "AliFMDHit.h"          // ALIFMDHIT_H
36 #include "AliFMDDigit.h"        // ALIFMDDigit_H
37 #include "AliFMDRecPoint.h"     // ALIFMDRECPOINT_H
38 #include <AliESD.h>
39 #include <AliESDFMD.h>
40 #include <TTree.h>              // ROOT_TTree
41 #include <TChain.h>             // ROOT_TChain
42 #include <TParticle.h>          // ROOT_TParticle
43 #include <TString.h>            // ROOT_TString
44 #include <TDatabasePDG.h>       // ROOT_TDatabasePDG
45 #include <TMath.h>              // ROOT_TMath
46 #include <TGeoManager.h>        // ROOT_TGeoManager 
47 #include <TSystemDirectory.h>   // ROOT_TSystemDirectory
48 #include <Riostream.h>          // ROOT_Riostream
49
50 //____________________________________________________________________
51 ClassImp(AliFMDInput)
52 #if 0
53   ; // This is here to keep Emacs for indenting the next line
54 #endif
55
56
57 //____________________________________________________________________
58 AliFMDInput::AliFMDInput()
59   : fGAliceFile(""), 
60     fLoader(0),
61     fRun(0), 
62     fStack(0),
63     fFMDLoader(0), 
64     fFMD(0),
65     fMainESD(0),
66     fESD(0),
67     fTreeE(0),
68     fTreeH(0),
69     fTreeD(0),
70     fTreeS(0),
71     fTreeR(0), 
72     fChainE(0),
73     fArrayE(0),
74     fArrayH(0),
75     fArrayD(0),
76     fArrayS(0), 
77     fArrayR(0), 
78     fTreeMask(0), 
79     fIsInit(kFALSE)
80 {
81   // Constructor of an FMD input object.  Specify what data to read in
82   // using the AddLoad member function.   Sub-classes should at a
83   // minimum overload the member function Event.   A full job can be
84   // executed using the member function Run. 
85 }
86
87   
88
89 //____________________________________________________________________
90 AliFMDInput::AliFMDInput(const char* gAliceFile)
91   : fGAliceFile(gAliceFile),
92     fLoader(0),
93     fRun(0), 
94     fStack(0),
95     fFMDLoader(0), 
96     fFMD(0),
97     fMainESD(0),
98     fESD(0),
99     fTreeE(0),
100     fTreeH(0),
101     fTreeD(0),
102     fTreeS(0),
103     fTreeR(0), 
104     fChainE(0),
105     fArrayE(0),
106     fArrayH(0),
107     fArrayD(0),
108     fArrayS(0), 
109     fArrayR(0), 
110     fTreeMask(0), 
111     fIsInit(kFALSE)
112 {
113   // Constructor of an FMD input object.  Specify what data to read in
114   // using the AddLoad member function.   Sub-classes should at a
115   // minimum overload the member function Event.   A full job can be
116   // executed using the member function Run. 
117 }
118
119 //____________________________________________________________________
120 Int_t
121 AliFMDInput::NEvents() const 
122 {
123   // Get number of events
124   if (fTreeE) return fTreeE->GetEntries();
125   return -1;
126 }
127
128 //____________________________________________________________________
129 Bool_t
130 AliFMDInput::Init()
131 {
132   // Initialize the object.  Get the needed loaders, and such. 
133
134   // Check if we have been initialized
135   if (fIsInit) { 
136     AliWarning("Already initialized");
137     return fIsInit;
138   }
139   if (fGAliceFile.IsNull()) fGAliceFile = "galice.root";
140   // Get the loader
141   fLoader = AliRunLoader::Open(fGAliceFile.Data(), "Alice", "read");
142   if (!fLoader) {
143     AliError(Form("Coulnd't read the file %s", fGAliceFile.Data()));
144     return kFALSE;
145   }
146   
147   // Get the run 
148   if  (fLoader->LoadgAlice()) return kFALSE;
149   fRun = fLoader->GetAliRun();
150   
151   // Get the FMD 
152   fFMD = static_cast<AliFMD*>(fRun->GetDetector("FMD"));
153   if (!fFMD) {
154     AliError("Failed to get detector FMD from loader");
155     return kFALSE;
156   }
157   
158   // Get the FMD loader
159   fFMDLoader = fLoader->GetLoader("FMDLoader");
160   if (!fFMDLoader) {
161     AliError("Failed to get detector FMD loader from loader");
162     return kFALSE;
163   }
164   if (fLoader->LoadHeader()) { 
165     AliError("Failed to get event header information from loader");
166     return kFALSE;
167   }
168   fTreeE = fLoader->TreeE();
169
170   // Optionally, get the ESD files
171   if (TESTBIT(fTreeMask, kESD)) {
172     fChainE = new TChain("esdTree");
173     TSystemDirectory dir;
174     TList*           files = dir.GetListOfFiles();
175     TSystemFile*     file;
176     files->Sort();
177     TIter            next(files);
178     while ((file = static_cast<TSystemFile*>(next()))) {
179       TString fname(file->GetName());
180       if (fname.Contains("AliESDs")) fChainE->AddFile(fname.Data());
181     }
182     fChainE->SetBranchAddress("ESD", &fMainESD);
183   }
184     
185   // Optionally, get the geometry 
186   if (TESTBIT(fTreeMask, kGeometry)) {
187     TString fname(fRun->GetGeometryFileName());
188     if (fname.IsNull()) {
189       Warning("Init", "No file name for the geometry from AliRun");
190       fname = gSystem->DirName(fGAliceFile);
191       fname.Append("/geometry.root");
192     }
193     fGeoManager = TGeoManager::Import(fname.Data());
194     if (!fGeoManager) {
195       Fatal("Init", "No geometry manager found");
196       return kFALSE;
197     }
198   }
199
200   
201   fIsInit = kTRUE;
202   return fIsInit;
203 }
204
205 //____________________________________________________________________
206 Bool_t
207 AliFMDInput::Begin(Int_t event)
208 {
209   // Called at the begining of each event.  Per default, it gets the
210   // data trees and gets pointers to the output arrays.   Users can
211   // overload this, but should call this member function in the
212   // overloaded member function of the derived class. 
213
214   // Check if we have been initialized
215   if (!fIsInit) { 
216     AliError("Not initialized");
217     return fIsInit;
218   }
219   // Get the event 
220   if (fLoader->GetEvent(event)) return kFALSE;
221   AliInfo(Form("Now in event %d/%d", event, NEvents()));
222
223   // Possibly load global kinematics information 
224   if (TESTBIT(fTreeMask, kKinematics)) {
225     AliInfo("Getting kinematics");
226     if (fLoader->LoadKinematics()) return kFALSE;
227     fStack = fLoader->Stack();
228   }
229   // Possibly load FMD Hit information 
230   if (TESTBIT(fTreeMask, kHits)) {
231     AliInfo("Getting FMD hits");
232     if (fFMDLoader->LoadHits()) return kFALSE;
233     fTreeH = fFMDLoader->TreeH();
234     if (!fArrayH) fArrayH = fFMD->Hits(); 
235   }
236   // Possibly load FMD Digit information 
237   if (TESTBIT(fTreeMask, kDigits)) {
238     AliInfo("Getting FMD digits");
239     if (fFMDLoader->LoadDigits()) return kFALSE;
240     fTreeD = fFMDLoader->TreeD();
241     if (!fArrayD) fArrayD = fFMD->Digits();
242   }
243   // Possibly load FMD Sdigit information 
244   if (TESTBIT(fTreeMask, kSDigits)) {
245     AliInfo("Getting FMD summable digits");
246     if (fFMDLoader->LoadSDigits()) return kFALSE;
247     fTreeS = fFMDLoader->TreeS();
248     if (!fArrayS) fArrayS = fFMD->SDigits();
249   }
250   // Possibly load FMD RecPoints information 
251   if (TESTBIT(fTreeMask, kRecPoints)) {
252     AliInfo("Getting FMD reconstructed points");
253     if (fFMDLoader->LoadRecPoints()) return kFALSE;
254     fTreeR = fFMDLoader->TreeR();
255     if (!fArrayR) fArrayR = new TClonesArray("AliFMDRecPoint");
256     fTreeR->SetBranchAddress("FMD",  &fArrayR);
257   }
258
259   // Possibly load FMD ESD information 
260   if (TESTBIT(fTreeMask, kESD)) {
261     AliInfo("Getting FMD event summary data");
262     Int_t read = fChainE->GetEntry(event);
263     if (read <= 0) return kFALSE;
264     fESD = fMainESD->GetFMDData();
265     if (!fESD) return kFALSE;
266   }
267   
268   return kTRUE;
269 }
270
271
272 //____________________________________________________________________
273 Bool_t 
274 AliFMDInput::Event()
275 {
276   // Process one event.  The default implementation one or more of 
277   //
278   //   -  ProcessHits       if the hits are loaded. 
279   //   -  ProcessDigits     if the digits are loaded. 
280   //   -  ProcessSDigits    if the sumbable digits are loaded. 
281   //   -  ProcessRecPoints  if the reconstructed points are loaded. 
282   //   -  ProcessESD        if the event summary data is loaded
283   // 
284   if (TESTBIT(fTreeMask, kHits)) 
285     if (!ProcessHits()) return kFALSE; 
286   if (TESTBIT(fTreeMask, kDigits)) 
287     if (!ProcessDigits()) return kFALSE;
288   if (TESTBIT(fTreeMask, kSDigits)) 
289     if (!ProcessSDigits()) return kFALSE;
290   if (TESTBIT(fTreeMask, kRecPoints)) 
291     if (!ProcessRecPoints()) return kFALSE;
292   if (TESTBIT(fTreeMask, kESD))
293     if (!ProcessESD(fESD)) return kFALSE;
294   
295   return kTRUE;
296 }
297
298 //____________________________________________________________________
299 Bool_t 
300 AliFMDInput::ProcessHits()
301 {
302   if (!fTreeH) {
303     AliError("No hit tree defined");
304     return kFALSE;
305   }
306   Int_t nTracks = fTreeH->GetEntries();
307   for (Int_t i = 0; i < nTracks; i++) {
308     Int_t hitRead  = fTreeH->GetEntry(i);
309     if (hitRead <= 0) continue;
310     if (!fArrayH) {
311       AliError("No hit array defined");
312       return kFALSE;
313     }
314     Int_t nHit = fArrayH->GetEntries();
315     if (nHit <= 0) continue;
316     for (Int_t j = 0; j < nHit; j++) {
317       AliFMDHit* hit = static_cast<AliFMDHit*>(fArrayH->At(j));
318       if (!hit) continue;
319       TParticle* track = 0;
320       if (TESTBIT(fTreeMask, kKinematics) && fStack) {
321         Int_t trackno = hit->Track();
322         track = fStack->Particle(trackno);
323       }
324       if (!ProcessHit(hit, track)) return kFALSE;
325     }    
326   }
327   return kTRUE;
328 }
329
330 //____________________________________________________________________
331 Bool_t 
332 AliFMDInput::ProcessDigits()
333 {
334   // Read the digit tree, and pass each digit to the member function
335   // ProcessDigit.
336   Int_t nEv = fTreeD->GetEntries();
337   for (Int_t i = 0; i < nEv; i++) {
338     Int_t digitRead  = fTreeD->GetEntry(i);
339     if (digitRead <= 0) continue;
340     Int_t nDigit = fArrayD->GetEntries();
341     if (nDigit <= 0) continue;
342     for (Int_t j = 0; j < nDigit; j++) {
343       AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayD->At(j));
344       if (!digit) continue;
345       if (!ProcessDigit(digit)) return kFALSE;
346     }    
347   }
348   return kTRUE;
349 }
350
351 //____________________________________________________________________
352 Bool_t 
353 AliFMDInput::ProcessSDigits()
354 {
355   // Read the summable digit tree, and pass each sumable digit to the
356   // member function ProcessSdigit.
357   Int_t nEv = fTreeD->GetEntries();
358   for (Int_t i = 0; i < nEv; i++) {
359     Int_t sdigitRead  = fTreeS->GetEntry(i);
360     if (sdigitRead <= 0) continue;
361     Int_t nSdigit = fArrayS->GetEntries();
362     if (nSdigit <= 0) continue;
363     for (Int_t j = 0; j < nSdigit; j++) {
364       AliFMDSDigit* sdigit = static_cast<AliFMDSDigit*>(fArrayS->At(j));
365       if (!sdigit) continue;
366       if (!ProcessSDigit(sdigit)) return kFALSE;
367     }    
368   }
369   return kTRUE;
370 }
371
372 //____________________________________________________________________
373 Bool_t 
374 AliFMDInput::ProcessRecPoints()
375 {
376   // Read the reconstrcted points tree, and pass each reconstruction
377   // object (AliFMDRecPoint) to either ProcessRecPoint.
378   Int_t nEv = fTreeR->GetEntries();
379   for (Int_t i = 0; i < nEv; i++) {
380     Int_t recRead  = fTreeR->GetEntry(i);
381     if (recRead <= 0) continue;
382     Int_t nRecPoint = fArrayR->GetEntries();
383     for (Int_t j = 0; j < nRecPoint; j++) {
384       AliFMDRecPoint* recPoint = static_cast<AliFMDRecPoint*>(fArrayR->At(j));
385       if (!recPoint) continue;
386       if (!ProcessRecPoint(recPoint)) return kFALSE;
387     }    
388   }
389   return kTRUE;
390 }
391
392 //____________________________________________________________________
393 Bool_t
394 AliFMDInput::End()
395 {
396   // Called at the end of each event.  Per default, it unloads the
397   // data trees and resets the pointers to the output arrays.   Users
398   // can overload this, but should call this member function in the
399   // overloaded member function of the derived class. 
400
401   // Check if we have been initialized
402   if (!fIsInit) { 
403     AliError("Not initialized");
404     return fIsInit;
405   }
406   // Possibly unload global kinematics information 
407   if (TESTBIT(fTreeMask, kKinematics)) {
408     fLoader->UnloadKinematics();
409     // fTreeK = 0;
410     fStack = 0;
411   }
412   // Possibly unload FMD Hit information 
413   if (TESTBIT(fTreeMask, kHits)) {
414     fFMDLoader->UnloadHits();
415     fTreeH = 0;
416   }
417   // Possibly unload FMD Digit information 
418   if (TESTBIT(fTreeMask, kDigits)) {
419     fFMDLoader->UnloadDigits();
420     fTreeD = 0;
421   }
422   // Possibly unload FMD Sdigit information 
423   if (TESTBIT(fTreeMask, kSDigits)) {
424     fFMDLoader->UnloadSDigits();
425     fTreeS = 0;
426   }
427   // Possibly unload FMD RecPoints information 
428   if (TESTBIT(fTreeMask, kRecPoints)) {
429     fFMDLoader->UnloadRecPoints();
430     fTreeR = 0;
431   }
432   AliInfo("Now out event");
433   return kTRUE;
434 }
435
436 //____________________________________________________________________
437 Bool_t
438 AliFMDInput::Run()
439 {
440   // Run over all events and files references in galice.root 
441
442   Bool_t retval;
443   if (!(retval = Init())) return retval;
444
445   Int_t nEvents = NEvents();
446   for (Int_t event = 0; event < nEvents; event++) {
447     if (!(retval = Begin(event))) break;
448     if (!(retval = Event())) break;
449     if (!(retval = End())) break;
450   }
451   if (!retval) return retval;
452   retval = Finish();
453   return retval;
454 }
455
456 //====================================================================
457 ClassImp(AliFMDInputHits)
458 #if 0
459   ;
460 #endif
461
462
463 //====================================================================
464 ClassImp(AliFMDInputDigits)
465 #if 0
466   ;
467 #endif
468
469 //====================================================================
470 ClassImp(AliFMDInputSDigits)
471 #if 0
472   ;
473 #endif
474
475 //====================================================================
476 ClassImp(AliFMDInputRecPoints)
477 #if 0
478   ;
479 #endif
480
481 //____________________________________________________________________
482 //
483 // EOF
484 //