281efaad7cf026dd07076b5682f69349bcf5a5a2
[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 /* $Id$ */
16 /** @file    AliFMDInput.cxx
17     @author  Christian Holm Christensen <cholm@nbi.dk>
18     @date    Mon Mar 27 12:42:40 2006
19     @brief   FMD utility classes for reading FMD data
20 */
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 "AliLog.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 "AliFMD.h"             // ALIFMD_H
39 #include "AliFMDHit.h"          // ALIFMDHIT_H
40 #include "AliFMDDigit.h"        // ALIFMDDigit_H
41 #include "AliFMDSDigit.h"       // ALIFMDDigit_H
42 #include "AliFMDRecPoint.h"     // ALIFMDRECPOINT_H
43 #include "AliFMDRawReader.h"    // ALIFMDRAWREADER_H
44 #include <AliESD.h>
45 #include <AliESDFMD.h>
46 #include <AliCDBManager.h>
47 #include <AliCDBEntry.h>
48 #include <AliAlignObjAngles.h>
49 #include <TTree.h>              // ROOT_TTree
50 #include <TChain.h>             // ROOT_TChain
51 #include <TParticle.h>          // ROOT_TParticle
52 #include <TString.h>            // ROOT_TString
53 #include <TDatabasePDG.h>       // ROOT_TDatabasePDG
54 #include <TMath.h>              // ROOT_TMath
55 #include <TGeoManager.h>        // ROOT_TGeoManager 
56 #include <TSystemDirectory.h>   // ROOT_TSystemDirectory
57 #include <Riostream.h>          // ROOT_Riostream
58
59 //____________________________________________________________________
60 ClassImp(AliFMDInput)
61 #if 0
62   ; // This is here to keep Emacs for indenting the next line
63 #endif
64
65
66 //____________________________________________________________________
67 AliFMDInput::AliFMDInput()
68   : fGAliceFile(""), 
69     fLoader(0),
70     fRun(0), 
71     fStack(0),
72     fFMDLoader(0), 
73     fReader(0),
74     fFMD(0),
75     fMainESD(0),
76     fESD(0),
77     fTreeE(0),
78     fTreeH(0),
79     fTreeD(0),
80     fTreeS(0),
81     fTreeR(0), 
82     fTreeA(0), 
83     fChainE(0),
84     fArrayE(0),
85     fArrayH(0),
86     fArrayD(0),
87     fArrayS(0), 
88     fArrayR(0), 
89     fArrayA(0), 
90     fGeoManager(0),
91     fTreeMask(0), 
92     fIsInit(kFALSE)
93 {
94   // Constructor of an FMD input object.  Specify what data to read in
95   // using the AddLoad member function.   Sub-classes should at a
96   // minimum overload the member function Event.   A full job can be
97   // executed using the member function Run. 
98 }
99
100   
101
102 //____________________________________________________________________
103 AliFMDInput::AliFMDInput(const char* gAliceFile)
104   : fGAliceFile(gAliceFile),
105     fLoader(0),
106     fRun(0), 
107     fStack(0),
108     fFMDLoader(0), 
109     fReader(0),
110     fFMD(0),
111     fMainESD(0),
112     fESD(0),
113     fTreeE(0),
114     fTreeH(0),
115     fTreeD(0),
116     fTreeS(0),
117     fTreeR(0), 
118     fTreeA(0), 
119     fChainE(0),
120     fArrayE(0),
121     fArrayH(0),
122     fArrayD(0),
123     fArrayS(0), 
124     fArrayR(0), 
125     fArrayA(0), 
126     fGeoManager(0),
127     fTreeMask(0), 
128     fIsInit(kFALSE)
129 {
130   // Constructor of an FMD input object.  Specify what data to read in
131   // using the AddLoad member function.   Sub-classes should at a
132   // minimum overload the member function Event.   A full job can be
133   // executed using the member function Run. 
134 }
135
136 //____________________________________________________________________
137 Int_t
138 AliFMDInput::NEvents() const 
139 {
140   // Get number of events
141   if (fTreeE) return fTreeE->GetEntries();
142   return -1;
143 }
144
145 //____________________________________________________________________
146 Bool_t
147 AliFMDInput::Init()
148 {
149   // Initialize the object.  Get the needed loaders, and such. 
150
151   // Check if we have been initialized
152   if (fIsInit) { 
153     AliWarning("Already initialized");
154     return fIsInit;
155   }
156   if (fGAliceFile.IsNull()) fGAliceFile = "galice.root";
157   // Get the loader
158   fLoader = AliRunLoader::Open(fGAliceFile.Data(), "Alice", "read");
159   if (!fLoader) {
160     AliError(Form("Coulnd't read the file %s", fGAliceFile.Data()));
161     return kFALSE;
162   }
163   
164   // Get the run 
165   if  (fLoader->LoadgAlice()) return kFALSE;
166   fRun = fLoader->GetAliRun();
167   
168   // Get the FMD 
169   fFMD = static_cast<AliFMD*>(fRun->GetDetector("FMD"));
170   if (!fFMD) {
171     AliError("Failed to get detector FMD from loader");
172     return kFALSE;
173   }
174   
175   // Get the FMD loader
176   fFMDLoader = fLoader->GetLoader("FMDLoader");
177   if (!fFMDLoader) {
178     AliError("Failed to get detector FMD loader from loader");
179     return kFALSE;
180   }
181   if (fLoader->LoadHeader()) { 
182     AliError("Failed to get event header information from loader");
183     return kFALSE;
184   }
185   fTreeE = fLoader->TreeE();
186
187   // Optionally, get the ESD files
188   if (TESTBIT(fTreeMask, kESD)) {
189     fChainE = new TChain("esdTree");
190     TSystemDirectory dir(".",".");
191     TList*           files = dir.GetListOfFiles();
192     TSystemFile*     file = 0;
193     if (!files) {
194       AliError("No files");
195       return kFALSE;
196     }
197     files->Sort();
198     TIter            next(files);
199     while ((file = static_cast<TSystemFile*>(next()))) {
200       TString fname(file->GetName());
201       if (fname.Contains("AliESDs")) fChainE->AddFile(fname.Data());
202     }
203     fChainE->SetBranchAddress("ESD", &fMainESD);
204   }
205     
206   if (TESTBIT(fTreeMask, kRaw)) {
207     AliInfo("Getting FMD raw data digits");
208     fArrayA = new TClonesArray("AliFMDDigit");
209     fReader = new AliRawReaderFile(-1);
210   }
211   
212   // Optionally, get the geometry 
213   if (TESTBIT(fTreeMask, kGeometry)) {
214     TString fname(fRun->GetGeometryFileName());
215     if (fname.IsNull()) {
216       Warning("Init", "No file name for the geometry from AliRun");
217       fname = gSystem->DirName(fGAliceFile);
218       fname.Append("/geometry.root");
219     }
220     fGeoManager = TGeoManager::Import(fname.Data());
221     if (!fGeoManager) {
222       Fatal("Init", "No geometry manager found");
223       return kFALSE;
224     }
225     AliCDBManager* cdb   = AliCDBManager::Instance();
226     AliCDBEntry*   align = cdb->Get("FMD/Align/Data");
227     if (align) {
228       AliInfo("Got alignment data from CDB");
229       TClonesArray* array = dynamic_cast<TClonesArray*>(align->GetObject());
230       if (!array) {
231         AliWarning("Invalid align data from CDB");
232       }
233       else {
234         Int_t nAlign = array->GetEntries();
235         for (Int_t i = 0; i < nAlign; i++) {
236           AliAlignObjAngles* a = static_cast<AliAlignObjAngles*>(array->At(i));
237           if (!a->ApplyToGeometry()) {
238             AliWarning(Form("Failed to apply alignment to %s", 
239                             a->GetSymName()));
240           }
241         }
242       }
243     }
244   }
245
246   
247   fIsInit = kTRUE;
248   return fIsInit;
249 }
250
251 //____________________________________________________________________
252 Bool_t
253 AliFMDInput::Begin(Int_t event)
254 {
255   // Called at the begining of each event.  Per default, it gets the
256   // data trees and gets pointers to the output arrays.   Users can
257   // overload this, but should call this member function in the
258   // overloaded member function of the derived class. 
259
260   // Check if we have been initialized
261   if (!fIsInit) { 
262     AliError("Not initialized");
263     return fIsInit;
264   }
265   // Get the event 
266   if (fLoader->GetEvent(event)) return kFALSE;
267   AliInfo(Form("Now in event %d/%d", event, NEvents()));
268
269   // Possibly load global kinematics information 
270   if (TESTBIT(fTreeMask, kKinematics)) {
271     AliInfo("Getting kinematics");
272     if (fLoader->LoadKinematics()) return kFALSE;
273     fStack = fLoader->Stack();
274   }
275   // Possibly load FMD Hit information 
276   if (TESTBIT(fTreeMask, kHits)) {
277     AliInfo("Getting FMD hits");
278     if (fFMDLoader->LoadHits()) return kFALSE;
279     fTreeH = fFMDLoader->TreeH();
280     if (!fArrayH) fArrayH = fFMD->Hits(); 
281   }
282   // Possibly load FMD Digit information 
283   if (TESTBIT(fTreeMask, kDigits)) {
284     AliInfo("Getting FMD digits");
285     if (fFMDLoader->LoadDigits()) return kFALSE;
286     fTreeD = fFMDLoader->TreeD();
287     if (fTreeD) {
288       if (!fArrayD) fArrayD = fFMD->Digits();
289     }
290     else {
291       fArrayD = 0;
292       AliWarning(Form("Failed to load FMD Digits"));
293     } 
294   }
295   // Possibly load FMD Sdigit information 
296   if (TESTBIT(fTreeMask, kSDigits)) {
297     AliInfo("Getting FMD summable digits");
298     if (fFMDLoader->LoadSDigits()) return kFALSE;
299     fTreeS = fFMDLoader->TreeS();
300     if (!fArrayS) fArrayS = fFMD->SDigits();
301   }
302   // Possibly load FMD RecPoints information 
303   if (TESTBIT(fTreeMask, kRecPoints)) {
304     AliInfo("Getting FMD reconstructed points");
305     if (fFMDLoader->LoadRecPoints()) return kFALSE;
306     fTreeR = fFMDLoader->TreeR();
307     if (!fArrayR) fArrayR = new TClonesArray("AliFMDRecPoint");
308     fTreeR->SetBranchAddress("FMD",  &fArrayR);
309   }  // Possibly load FMD ESD information 
310   if (TESTBIT(fTreeMask, kESD)) {
311     AliInfo("Getting FMD event summary data");
312     Int_t read = fChainE->GetEntry(event);
313     if (read <= 0) return kFALSE;
314     fESD = fMainESD->GetFMDData();
315     if (!fESD) return kFALSE;
316     fESD->CheckNeedUShort(fChainE->GetFile());
317   }
318   // Possibly load FMD Digit information 
319   if (TESTBIT(fTreeMask, kRaw)) {
320     AliInfo("Getting FMD raw data digits");
321     if (!fReader->NextEvent()) return kFALSE;
322     AliFMDRawReader r(fReader, 0);
323     fArrayA->Clear();
324     r.ReadAdcs(fArrayA);
325   }
326   
327   return kTRUE;
328 }
329
330
331 //____________________________________________________________________
332 Bool_t 
333 AliFMDInput::Event()
334 {
335   // Process one event.  The default implementation one or more of 
336   //
337   //   -  ProcessHits       if the hits are loaded. 
338   //   -  ProcessDigits     if the digits are loaded. 
339   //   -  ProcessSDigits    if the sumbable digits are loaded. 
340   //   -  ProcessRecPoints  if the reconstructed points are loaded. 
341   //   -  ProcessESD        if the event summary data is loaded
342   // 
343   if (TESTBIT(fTreeMask, kHits)) 
344     if (!ProcessHits()) return kFALSE; 
345   if (TESTBIT(fTreeMask, kDigits)) 
346     if (!ProcessDigits()) return kFALSE;
347   if (TESTBIT(fTreeMask, kSDigits)) 
348     if (!ProcessSDigits()) return kFALSE;
349   if (TESTBIT(fTreeMask, kRaw)) 
350     if (!ProcessRawDigits()) return kFALSE;
351   if (TESTBIT(fTreeMask, kRecPoints)) 
352     if (!ProcessRecPoints()) return kFALSE;
353   if (TESTBIT(fTreeMask, kESD))
354     if (!ProcessESDs()) return kFALSE;
355   
356   return kTRUE;
357 }
358
359 //____________________________________________________________________
360 Bool_t 
361 AliFMDInput::ProcessHits()
362 {
363   // Read the hit tree, and pass each hit to the member function
364   // ProcessHit.
365   if (!fTreeH) {
366     AliError("No hit tree defined");
367     return kFALSE;
368   }
369   Int_t nTracks = fTreeH->GetEntries();
370   for (Int_t i = 0; i < nTracks; i++) {
371     Int_t hitRead  = fTreeH->GetEntry(i);
372     if (hitRead <= 0) continue;
373     if (!fArrayH) {
374       AliError("No hit array defined");
375       return kFALSE;
376     }
377     Int_t nHit = fArrayH->GetEntries();
378     if (nHit <= 0) continue;
379     for (Int_t j = 0; j < nHit; j++) {
380       AliFMDHit* hit = static_cast<AliFMDHit*>(fArrayH->At(j));
381       if (!hit) continue;
382       TParticle* track = 0;
383       if (TESTBIT(fTreeMask, kKinematics) && fStack) {
384         Int_t trackno = hit->Track();
385         track = fStack->Particle(trackno);
386       }
387       if (!ProcessHit(hit, track)) return kFALSE;
388     }    
389   }
390   return kTRUE;
391 }
392
393 //____________________________________________________________________
394 Bool_t 
395 AliFMDInput::ProcessDigits()
396 {
397   // Read the digit tree, and pass each digit to the member function
398   // ProcessDigit.
399   Int_t nEv = fTreeD->GetEntries();
400   for (Int_t i = 0; i < nEv; i++) {
401     Int_t digitRead  = fTreeD->GetEntry(i);
402     if (digitRead <= 0) continue;
403     Int_t nDigit = fArrayD->GetEntries();
404     if (nDigit <= 0) continue;
405     for (Int_t j = 0; j < nDigit; j++) {
406       AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayD->At(j));
407       if (!digit) continue;
408       if (!ProcessDigit(digit)) return kFALSE;
409     }    
410   }
411   return kTRUE;
412 }
413
414 //____________________________________________________________________
415 Bool_t 
416 AliFMDInput::ProcessSDigits()
417 {
418   // Read the summable digit tree, and pass each sumable digit to the
419   // member function ProcessSdigit.
420   Int_t nEv = fTreeD->GetEntries();
421   for (Int_t i = 0; i < nEv; i++) {
422     Int_t sdigitRead  = fTreeS->GetEntry(i);
423     if (sdigitRead <= 0) continue;
424     Int_t nSdigit = fArrayS->GetEntries();
425     if (nSdigit <= 0) continue;
426     for (Int_t j = 0; j < nSdigit; j++) {
427       AliFMDSDigit* sdigit = static_cast<AliFMDSDigit*>(fArrayS->At(j));
428       if (!sdigit) continue;
429       if (!ProcessSDigit(sdigit)) return kFALSE;
430     }    
431   }
432   return kTRUE;
433 }
434
435 //____________________________________________________________________
436 Bool_t 
437 AliFMDInput::ProcessRawDigits()
438 {
439   // Read the digit tree, and pass each digit to the member function
440   // ProcessDigit.
441   Int_t nDigit = fArrayA->GetEntries();
442   if (nDigit <= 0) return kTRUE;
443   for (Int_t j = 0; j < nDigit; j++) {
444     AliFMDDigit* digit = static_cast<AliFMDDigit*>(fArrayA->At(j));
445     if (!digit) continue;
446     if (!ProcessRawDigit(digit)) return kFALSE;
447   }    
448   return kTRUE;
449 }
450
451 //____________________________________________________________________
452 Bool_t 
453 AliFMDInput::ProcessRecPoints()
454 {
455   // Read the reconstrcted points tree, and pass each reconstruction
456   // object (AliFMDRecPoint) to either ProcessRecPoint.
457   Int_t nEv = fTreeR->GetEntries();
458   for (Int_t i = 0; i < nEv; i++) {
459     Int_t recRead  = fTreeR->GetEntry(i);
460     if (recRead <= 0) continue;
461     Int_t nRecPoint = fArrayR->GetEntries();
462     for (Int_t j = 0; j < nRecPoint; j++) {
463       AliFMDRecPoint* recPoint = static_cast<AliFMDRecPoint*>(fArrayR->At(j));
464       if (!recPoint) continue;
465       if (!ProcessRecPoint(recPoint)) return kFALSE;
466     }    
467   }
468   return kTRUE;
469 }
470
471 //____________________________________________________________________
472 Bool_t 
473 AliFMDInput::ProcessESDs()
474 {
475   // Process event summary data
476   if (!fESD) return kFALSE;
477   for (UShort_t det = 1; det <= 3; det++) {
478     Char_t rings[] = { 'I', (det == 1 ? '\0' : 'O'), '\0' };
479     for (Char_t* rng = rings; *rng != '\0'; rng++) {
480       UShort_t nsec = (*rng == 'I' ?  20 :  40);
481       UShort_t nstr = (*rng == 'I' ? 512 : 256);
482       for (UShort_t sec = 0; sec < nsec; sec++) {
483         for (UShort_t str = 0; str < nstr; str++) {
484           Float_t eta  = fESD->Eta(det,*rng,sec,str);
485           Float_t mult = fESD->Multiplicity(det,*rng,sec,str);
486           if (!ProcessESD(det, *rng, sec, str, eta, mult)) continue;
487         }
488       }
489     }
490   }
491   return kTRUE;
492 }
493
494 //____________________________________________________________________
495 Bool_t
496 AliFMDInput::End()
497 {
498   // Called at the end of each event.  Per default, it unloads the
499   // data trees and resets the pointers to the output arrays.   Users
500   // can overload this, but should call this member function in the
501   // overloaded member function of the derived class. 
502
503   // Check if we have been initialized
504   if (!fIsInit) { 
505     AliError("Not initialized");
506     return fIsInit;
507   }
508   // Possibly unload global kinematics information 
509   if (TESTBIT(fTreeMask, kKinematics)) {
510     fLoader->UnloadKinematics();
511     // fTreeK = 0;
512     fStack = 0;
513   }
514   // Possibly unload FMD Hit information 
515   if (TESTBIT(fTreeMask, kHits)) {
516     fFMDLoader->UnloadHits();
517     fTreeH = 0;
518   }
519   // Possibly unload FMD Digit information 
520   if (TESTBIT(fTreeMask, kDigits)) {
521     fFMDLoader->UnloadDigits();
522     fTreeD = 0;
523   }
524   // Possibly unload FMD Sdigit information 
525   if (TESTBIT(fTreeMask, kSDigits)) {
526     fFMDLoader->UnloadSDigits();
527     fTreeS = 0;
528   }
529   // Possibly unload FMD RecPoints information 
530   if (TESTBIT(fTreeMask, kRecPoints)) {
531     fFMDLoader->UnloadRecPoints();
532     fTreeR = 0;
533   }
534   AliInfo("Now out event");
535   return kTRUE;
536 }
537
538 //____________________________________________________________________
539 Bool_t
540 AliFMDInput::Run()
541 {
542   // Run over all events and files references in galice.root 
543
544   Bool_t retval;
545   if (!(retval = Init())) return retval;
546
547   Int_t nEvents = NEvents();
548   for (Int_t event = 0; event < nEvents; event++) {
549     if (!(retval = Begin(event))) break;
550     if (!(retval = Event())) break;
551     if (!(retval = End())) break;
552   }
553   if (!retval) return retval;
554   retval = Finish();
555   return retval;
556 }
557
558
559 //____________________________________________________________________
560 //
561 // EOF
562 //