]> git.uio.no Git - u/mrichter/AliRoot.git/blob - STEER/AliQADataMakerSteer.cxx
Added list of active detectors
[u/mrichter/AliRoot.git] / STEER / AliQADataMakerSteer.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, 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 #include <TFile.h>
19 #include <TFileMerger.h>
20 #include <TPluginManager.h>
21 #include <TROOT.h>
22 #include <TString.h>
23 #include <TSystem.h>
24
25 #include "AliESDEvent.h"
26 #include "AliHeader.h"
27 #include "AliLog.h"
28 #include "AliModule.h"
29 #include "AliQA.h"
30 #include "AliQADataMaker.h"
31 #include "AliQADataMakerSteer.h" 
32 #include "AliRawReaderDate.h"
33 #include "AliRawReaderFile.h"
34 #include "AliRawReaderRoot.h"
35 #include "AliRun.h"
36 #include "AliRunLoader.h"
37
38 ClassImp(AliQADataMakerSteer) 
39
40 //_____________________________________________________________________________
41 AliQADataMakerSteer::AliQADataMakerSteer(const char* gAliceFilename, const char * name, const char * title) :
42         TNamed(name, title), 
43         fCycleSame(kFALSE),
44     fDetectors("ALL"), 
45         fESD(NULL), 
46         fESDTree(NULL),
47         fFirst(kTRUE),  
48         fGAliceFileName(gAliceFilename), 
49         fRunNumber(0), 
50         fNumberOfEvents(0), 
51         fRawReader(NULL), 
52         fRawReaderDelete(kTRUE), 
53         fRunLoader(NULL)  
54 {
55         // default ctor
56         for (UInt_t iDet = 0; iDet < fgkNDetectors; iDet++) {
57                 if (IsSelected(AliQA::GetDetName(iDet))) {
58                         fLoader[iDet]      = NULL ;
59                         fQADataMaker[iDet] = NULL ;
60                         fQACycles[iDet]    = 999999 ;
61                 }
62         }       
63 }
64
65 //_____________________________________________________________________________
66 AliQADataMakerSteer::AliQADataMakerSteer(const AliQADataMakerSteer & qas) : 
67         TNamed(qas), 
68         fCycleSame(kFALSE),
69     fDetectors(qas.fDetectors), 
70         fESD(NULL), 
71         fESDTree(NULL), 
72         fFirst(qas.fFirst),  
73         fGAliceFileName(qas.fGAliceFileName), 
74         fRunNumber(qas.fRunNumber), 
75         fNumberOfEvents(qas.fNumberOfEvents), 
76         fRawReader(NULL), 
77         fRawReaderDelete(kTRUE), 
78         fRunLoader(NULL)  
79 {
80         // cpy ctor
81         for (UInt_t iDet = 0; iDet < fgkNDetectors; iDet++) {
82                 fLoader[iDet]      = qas.fLoader[iDet] ;
83                 fQADataMaker[iDet] = qas.fQADataMaker[iDet] ;
84                 fQACycles[iDet]    = qas.fQACycles[iDet] ;      
85         }
86 }
87
88 //_____________________________________________________________________________
89 AliQADataMakerSteer & AliQADataMakerSteer::operator = (const AliQADataMakerSteer & qas) 
90 {
91         // assignment operator
92   this->~AliQADataMakerSteer() ;
93   new(this) AliQADataMakerSteer(qas) ;
94   return *this ;
95 }
96
97 //_____________________________________________________________________________
98 AliQADataMakerSteer::~AliQADataMakerSteer() 
99 {
100         // dtor
101   for (UInt_t iDet = 0; iDet < fgkNDetectors; iDet++) {
102           if (IsSelected(AliQA::GetDetName(iDet))) {
103                   fLoader[iDet] = NULL;
104                   if (fQADataMaker[iDet]) {
105                           (fQADataMaker[iDet])->Finish() ; 
106                           delete fQADataMaker[iDet] ;
107                           fQADataMaker[iDet] = NULL ;
108                   }
109           }
110   }
111
112   if (fRawReaderDelete) { 
113         fRunLoader = NULL ;
114         delete fRawReader ;
115         fRawReader = NULL ;
116   }
117 }
118
119 //_____________________________________________________________________________
120 Bool_t AliQADataMakerSteer::DoIt(const AliQA::TASKINDEX taskIndex)
121 {
122         // Runs all the QA data Maker for every detector
123         Bool_t rv = kFALSE ;
124         
125     // Fill QA data in event loop 
126         for (UInt_t iEvent = 0 ; iEvent < fNumberOfEvents ; iEvent++) {
127                 // Get the event
128                 AliDebug(1, Form("processing event %d", iEvent));
129                 if ( taskIndex == AliQA::kRAWS ) {
130                         if ( !fRawReader->NextEvent() )
131                                 break ;
132                 } else if ( taskIndex == AliQA::kESDS ) {
133                         fESDTree->GetEntry(iEvent) ;
134                 } else {
135                         fRunLoader->GetEvent(iEvent);
136                 }
137                 // loop over detectors
138                 TObjArray* detArray = NULL ; 
139                 if (fRunLoader) // check if RunLoader exists 
140                         if ( fRunLoader->GetAliRun() ) // check if AliRun exists in gAlice.root
141                                 detArray = fRunLoader->GetAliRun()->Detectors() ;
142                 for (UInt_t iDet = 0 ; iDet < fgkNDetectors ; iDet++) {
143                         if (detArray) {
144                                 AliModule* det = static_cast<AliModule*>(detArray->FindObject(AliQA::GetDetName(iDet))) ;
145                                 if (!det || !det->IsActive()) 
146                                         continue ;
147                         }
148                         if (!IsSelected(AliQA::GetDetName(iDet)))
149                                 continue ;
150                         AliQADataMaker * qadm = GetQADataMaker(iDet) ;
151                         if (!qadm) {
152                                 rv = kFALSE ;
153                         } else {
154                                 if ( qadm->IsCycleDone() ) {
155                                         qadm->EndOfCycle(AliQA::kRAWS) ;
156                                         qadm->StartOfCycle(AliQA::kRAWS) ;
157                                 }
158                                 TTree * data ; 
159                                 switch (taskIndex) {
160                                         case AliQA::kRAWS :
161                                                 qadm->Exec(taskIndex, fRawReader) ; 
162                                                 break ; 
163                                         case AliQA::kHITS :
164                                                 GetLoader(iDet)->LoadHits() ; 
165                                                 data = GetLoader(iDet)->TreeH() ; 
166                                                 if ( ! data ) {
167                                                         AliWarning(Form(" Hit Tree not found for  %s", AliQA::GetDetName(iDet))) ; 
168                                                 } else {
169                                                         qadm->Exec(taskIndex, data) ;
170                                                 } 
171                                                 break ;
172                                                 case AliQA::kSDIGITS :
173                                                 GetLoader(iDet)->LoadSDigits() ; 
174                                                 data = GetLoader(iDet)->TreeS() ; 
175                                                 if ( ! data ) {
176                                                         AliWarning(Form(" SDigit Tree not found for  %s", AliQA::GetDetName(iDet))) ; 
177                                                 } else {
178                                                         qadm->Exec(taskIndex, data) ; 
179                                                 }
180                                                 break; 
181                                                 case AliQA::kDIGITS :
182                                                 GetLoader(iDet)->LoadDigits() ; 
183                                                 data = GetLoader(iDet)->TreeD() ; 
184                                                 if ( ! data ) {
185                                                         AliWarning(Form(" Digit Tree not found for  %s", AliQA::GetDetName(iDet))) ; 
186                                                 } else {
187                                                         qadm->Exec(taskIndex, data) ;
188                                                 }
189                                                 break; 
190                                                 case AliQA::kRECPOINTS :
191                                                 GetLoader(iDet)->LoadRecPoints() ; 
192                                                 data = GetLoader(iDet)->TreeR() ; 
193                                                 if (!data) {
194                                                         AliWarning(Form("RecPoints not found for %s", AliQA::GetDetName(iDet))) ; 
195                                                 } else {
196                                                         qadm->Exec(taskIndex, data) ; 
197                                                 }
198                                                 break; 
199                                                 case AliQA::kTRACKSEGMENTS :
200                                                 break; 
201                                                 case AliQA::kRECPARTICLES :
202                                                 break; 
203                                                 case AliQA::kESDS :
204                                                 qadm->Exec(taskIndex, fESD) ;
205                                                 break; 
206                                 } //task switch
207                                 qadm->Increment() ; 
208                         } //data maker exist
209                 } // detector loop
210         } // event loop 
211         // Save QA data for all detectors
212         rv = Finish(taskIndex) ;
213         return rv ; 
214 }
215
216 //_____________________________________________________________________________
217 AliLoader * AliQADataMakerSteer::GetLoader(Int_t iDet)
218 {
219         // get the loader for a detector
220         
221         TString detName = AliQA::GetDetName(iDet) ;
222     fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
223         if (fLoader[iDet]) 
224                 return fLoader[iDet] ;
225         
226         // load the QA data maker object
227         TPluginManager* pluginManager = gROOT->GetPluginManager() ;
228         TString loaderName = "Ali" + detName + "Loader" ;
229
230         AliLoader * loader = NULL ;
231         // first check if a plugin is defined for the quality assurance data maker
232         TPluginHandler* pluginHandler = pluginManager->FindHandler("AliLoader", detName) ;
233         // if not, add a plugin for it
234         if (!pluginHandler) {
235                 AliDebug(1, Form("defining plugin for %s", loaderName.Data())) ;
236                 TString libs = gSystem->GetLibraries() ;
237                 if (libs.Contains("lib" + detName + "base.so") || (gSystem->Load("lib" + detName + "base.so") >= 0)) {
238                         pluginManager->AddHandler("AliQADataMaker", detName, loaderName, detName + "loader", loaderName + "()") ;
239                 } else {
240                         pluginManager->AddHandler("AliLoader", detName, loaderName, detName, loaderName + "()") ;
241                 }
242                 pluginHandler = pluginManager->FindHandler("AliLoader", detName) ;
243         }
244         if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
245                 loader = (AliLoader *) pluginHandler->ExecPlugin(0) ;
246         }
247         if (loader) 
248                 fLoader[iDet] = loader ;
249         return loader ;
250 }
251
252 //_____________________________________________________________________________
253 AliQADataMaker * AliQADataMakerSteer::GetQADataMaker(Int_t iDet)
254 {
255         // get the quality assurance data maker for a detector
256         
257         if (fQADataMaker[iDet]) 
258                 return fQADataMaker[iDet] ;
259         
260         AliQADataMaker * qadm = NULL ;
261         
262         if (fFirst) {
263                 // load the QA data maker object
264                 TPluginManager* pluginManager = gROOT->GetPluginManager() ;
265                 TString detName = AliQA::GetDetName(iDet) ;
266                 TString qadmName = "Ali" + detName + "QADataMaker" ;
267
268                 // first check if a plugin is defined for the quality assurance data maker
269                 TPluginHandler* pluginHandler = pluginManager->FindHandler("AliQADataMaker", detName) ;
270                 // if not, add a plugin for it
271                 if (!pluginHandler) {
272                         AliDebug(1, Form("defining plugin for %s", qadmName.Data())) ;
273                         TString libs = gSystem->GetLibraries() ;
274                         if (libs.Contains("lib" + detName + "base.so") || (gSystem->Load("lib" + detName + "base.so") >= 0)) {
275                                 pluginManager->AddHandler("AliQADataMaker", detName, qadmName, detName + "qadm", qadmName + "()") ;
276                         } else {
277                                 pluginManager->AddHandler("AliQADataMaker", detName, qadmName, detName, qadmName + "()") ;
278                         }
279                         pluginHandler = pluginManager->FindHandler("AliQADataMaker", detName) ;
280                 }
281                 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
282                         qadm = (AliQADataMaker *) pluginHandler->ExecPlugin(0) ;
283                 }
284                 if (qadm) 
285                         fQADataMaker[iDet] = qadm ;
286         }
287         return qadm ;
288 }
289
290 //_____________________________________________________________________________
291 Bool_t AliQADataMakerSteer::Init(const AliQA::TASKINDEX taskIndex, const  char * input )
292 {
293         // Initialize the event source and QA data makers
294         
295         if (taskIndex == AliQA::kRAWS) { 
296                 if (!fRawReader) {
297                         TString fileName(input);
298                         if (fileName.EndsWith("/")) {
299                                 fRawReader = new AliRawReaderFile(fileName);
300                         } else if (fileName.EndsWith(".root")) {
301                                 fRawReader = new AliRawReaderRoot(fileName);
302                         } else if (!fileName.IsNull()) {
303                                 fRawReader = new AliRawReaderDate(fileName);
304                                 fRawReader->SelectEvents(7);
305                         }
306                 }
307             if ( ! fRawReader ) 
308                         return kFALSE ; 
309                 fRawReader->NextEvent() ; 
310                 fRunNumber = fRawReader->GetRunNumber() ; 
311                 fRawReader->RewindEvents();
312                 fNumberOfEvents = 999999 ;
313         } else if (taskIndex == AliQA::kESDS) {
314                 if (!gSystem->AccessPathName("AliESDs.root")) { // AliESDs.root exists
315                         TFile * esdFile = TFile::Open("AliESDs.root") ;
316                         fESDTree = dynamic_cast<TTree *> (esdFile->Get("esdTree")) ; 
317                         fESD     = new AliESDEvent() ;
318                         fESD->ReadFromTree(fESDTree) ;
319                         fESDTree->GetEntry(0) ; 
320                         fRunNumber = fESD->GetRunNumber() ; 
321                         fNumberOfEvents = fESDTree->GetEntries() ;
322                 } else {
323                         AliError("AliESDs.root not found") ; 
324                         return kFALSE ; 
325                 }                       
326         } else {
327                 if ( !InitRunLoader() ) {
328                         AliError("Run Loader not found") ; 
329                 } else {
330                         //if (fRunLoader->GetHeader()) 
331 //                              fRunNumber      = fRunLoader->GetHeader()->GetRun() ;
332 //                      else
333                                 fRunNumber      = 0 ; 
334                         fNumberOfEvents = fRunLoader->GetNumberOfEvents() ;
335                 }
336         }
337         // Initialize all QA data makers for all detectors
338         for (UInt_t iDet = 0; iDet < fgkNDetectors ; iDet++) {
339                 if (IsSelected(AliQA::GetDetName(iDet))) {
340                         AliQADataMaker * qadm = GetQADataMaker(iDet) ;
341                         if (!qadm) {
342                                 AliWarning(Form("AliQADataMaker not found for %s", AliQA::GetDetName(iDet))) ; 
343                         } else {
344                                 AliInfo(Form("Data Maker found for %s", qadm->GetName())) ; 
345                                 qadm->Init(taskIndex, fRunNumber, GetQACycles(iDet)) ;
346                                 qadm->StartOfCycle(taskIndex, fCycleSame) ;
347                         }
348                 }
349         } 
350         fFirst = kFALSE ;
351         return kTRUE ; 
352 }
353
354 //_____________________________________________________________________________
355 Bool_t AliQADataMakerSteer::IsSelected(const char * det) 
356 {
357         // check whether detName is contained in detectors
358         // if yes, it is removed from detectors
359         
360         const TString detName(det) ;
361         // check if all detectors are selected
362         if ((fDetectors.CompareTo("ALL") == 0) ||
363                 fDetectors.BeginsWith("ALL ") ||
364                 fDetectors.EndsWith(" ALL") ||
365                 fDetectors.Contains(" ALL ")) {
366                 fDetectors = "ALL";
367                 return kTRUE;
368         }
369         
370         // search for the given detector
371         Bool_t rv = kFALSE;
372         if ((fDetectors.CompareTo(detName) == 0) ||
373                 fDetectors.BeginsWith(detName+" ") ||
374                 fDetectors.EndsWith(" "+detName) ||
375                 fDetectors.Contains(" "+detName+" ")) {
376 //              fDetectors.ReplaceAll(detName, "");
377                 rv = kTRUE;
378         }
379         
380         // clean up the detectors string
381 //      while (fDetectors.Contains("  ")) 
382 //              fDetectors.ReplaceAll("  ", " ");
383 //      while (fDetectors.BeginsWith(" ")) 
384 //              fDetectors.Remove(0, 1);
385 //      while (fDetectors.EndsWith(" ")) 
386 //              fDetectors.Remove(fDetectors.Length()-1, 1);
387         
388         return rv ;
389 }
390
391 //_____________________________________________________________________________
392 Bool_t AliQADataMakerSteer::InitRunLoader()
393 {
394         // get or create the run loader
395         if (fRunLoader) {
396                 fCycleSame = kTRUE ; 
397                 return kTRUE ;
398         } 
399                 
400         if (!gSystem->AccessPathName(fGAliceFileName.Data())) { // galice.root exists
401     // load all base libraries to get the loader classes
402                 TString libs = gSystem->GetLibraries() ;
403                 for (UInt_t iDet = 0; iDet < fgkNDetectors; iDet++) {
404                         if (!IsSelected(AliQA::GetDetName(iDet))) 
405                                 continue ; 
406                         TString detName = AliQA::GetDetName(iDet) ;
407                         if (detName == "HLT") 
408                                 continue;
409                         if (libs.Contains("lib" + detName + "base.so")) 
410                                 continue;
411                         gSystem->Load("lib" + detName + "base.so");
412                 }
413                 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
414                 if (!fRunLoader) {
415                         AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
416                         return kFALSE;
417                 }
418                 fRunLoader->CdGAFile();
419                 if (fRunLoader->LoadgAlice() == 0) {
420                         gAlice = fRunLoader->GetAliRun();
421                 }
422
423                 if (!gAlice) {
424                         AliError(Form("no gAlice object found in file %s", fGAliceFileName.Data()));
425                         return kFALSE;
426                 }
427
428         } else {               // galice.root does not exist
429                 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
430                 return kFALSE;
431     }
432
433   return kTRUE;
434 }
435
436 //_____________________________________________________________________________
437 Bool_t AliQADataMakerSteer::Finish(const AliQA::TASKINDEX taskIndex) 
438 {
439         // write output to file for all detectors
440         for (UInt_t iDet = 0; iDet < fgkNDetectors ; iDet++) {
441                 if (IsSelected(AliQA::GetDetName(iDet))) {
442                         AliQADataMaker * qadm = GetQADataMaker(iDet) ;
443                         if (qadm) {
444                                 qadm->EndOfCycle(taskIndex) ; 
445                         }
446                 }
447         }
448         return kTRUE ; 
449 }
450
451 //_____________________________________________________________________________
452 Bool_t AliQADataMakerSteer::Merge() 
453 {
454         // Merge all the cycles from all detectors in one single file per run
455         
456         gROOT->ProcessLine(".! ls *QA*.*.*.root > tempo.txt") ; 
457         ifstream in("tempo.txt") ; 
458         const Int_t runMax = 10 ;  
459         TString file[AliQA::kNDET*runMax] ;
460         Int_t run[AliQA::kNDET*runMax] ;
461         
462         Int_t index = 0 ; 
463         while ( 1 ) {
464                 in >> file[index++] ; 
465                 if ( !in.good() ) 
466                         break ; 
467         }
468         Int_t previousRun = -1 ;
469         Int_t runIndex = 0 ;  
470         for (Int_t ifile = 0 ; ifile < index-1 ; ifile++) {
471                 TString tmp(file[ifile]) ; 
472                 tmp.ReplaceAll(".root", "") ; 
473                 TString det = tmp(0, tmp.Index(".")) ; 
474                 tmp.Remove(0, tmp.Index(".QA.")+4) ; 
475                 TString ttmp = tmp(0, tmp.Index(".")) ; 
476                 Int_t newRun = ttmp.Atoi() ;
477                 if (newRun != previousRun) {
478                         run[runIndex] = newRun ; 
479                         previousRun = newRun ; 
480                         runIndex++ ; 
481                 }
482                 ttmp = tmp(tmp.Index("."), tmp.Length()) ; 
483                 Int_t cycle = ttmp.Atoi() ;  
484                 AliInfo(Form("%s : det = %s run = %d cycle = %d \n", file[ifile].Data(), det.Data(), newRun, cycle)) ; 
485         }
486         for (Int_t irun = 0 ; irun < runIndex ; irun++) {
487                 TFileMerger merger ; 
488                 char outFileName[runMax] ; 
489                 sprintf(outFileName, "Merged.QA.%d.root", runIndex-1) ; 
490                 merger.OutputFile(outFileName) ; 
491                 for (Int_t ifile = 0 ; ifile < index-1 ; ifile++) {
492                         char pattern[100] ; 
493                         sprintf(pattern, "QA.%d.", runIndex-1) ; 
494                         TString tmp(file[ifile]) ; 
495                         if (tmp.Contains(pattern))
496                                 merger.AddFile(tmp) ; 
497                 }
498                 merger.Merge() ; 
499         }
500         
501         return kTRUE ; 
502 }
503
504 //_____________________________________________________________________________
505 void AliQADataMakerSteer::Reset()
506 {
507         // Reset the default data members
508         for (UInt_t iDet = 0; iDet < fgkNDetectors; iDet++) {
509                 if (IsSelected(AliQA::GetDetName(iDet))) {
510                         fLoader[iDet] = NULL;
511                         if (fQADataMaker[iDet]) {
512                                 (fQADataMaker[iDet])->Reset() ; 
513                                 //delete fQADataMaker[iDet] ;
514                                 //fQADataMaker[iDet] = NULL ;
515                         }
516                 }
517         }
518
519         if (fRawReaderDelete) { 
520                 delete fRawReader ;
521                 fRawReader      = NULL ;
522         }
523
524         fCycleSame      = kFALSE ; 
525         fESD            = NULL ; 
526         fESDTree        = NULL ; 
527         fFirst          = kTRUE ;   
528         fNumberOfEvents = 0 ;  
529 }
530
531 //_____________________________________________________________________________
532 Bool_t AliQADataMakerSteer::Run(const char * detectors, AliRawReader * rawReader) 
533 {
534         //Runs all the QA data Maker for Raws only
535         fRawReader       = rawReader ;          
536         fRawReaderDelete = kFALSE ; 
537         fCycleSame       = kTRUE ; 
538         fDetectors       = detectors ; 
539
540         // Initialize all QA data makers for all detectors
541         for (UInt_t iDet = 0; iDet < fgkNDetectors ; iDet++) {
542                 if (IsSelected(AliQA::GetDetName(iDet))) {
543                         AliQADataMaker * qadm = GetQADataMaker(iDet) ;
544                         if (!qadm) {
545                                 AliWarning(Form("AliQADataMaker not found for %s", AliQA::GetDetName(iDet))) ; 
546                         } else {
547                                 AliInfo(Form("Data Maker found for %s", qadm->GetName())) ; 
548                                 qadm->Init(AliQA::kRAWS, fRunNumber, GetQACycles(iDet)) ;
549                                 qadm->StartOfCycle(AliQA::kRAWS, fCycleSame) ;
550                         }
551                 }
552         } 
553         fFirst = kFALSE ;
554                 
555         return DoIt(AliQA::kRAWS) ; 
556 }
557
558 //_____________________________________________________________________________
559 Bool_t AliQADataMakerSteer::Run(const char * detectors, const AliQA::TASKINDEX taskIndex, const  char * fileName )
560 {
561         // Runs all the QA data Maker for every detector
562
563         Bool_t rv  = kFALSE ;
564         fDetectors = detectors ; 
565         
566         if ( !Init(taskIndex, fileName) ) 
567                 return kFALSE ; 
568
569         rv = DoIt(taskIndex) ;
570         
571         return rv ;   
572
573 }