1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
19 #include <TFileMerger.h>
20 #include <TPluginManager.h>
25 #include "AliESDEvent.h"
26 #include "AliHeader.h"
28 #include "AliModule.h"
30 #include "AliQADataMaker.h"
31 #include "AliQADataMakerSteer.h"
32 #include "AliRawReaderDate.h"
33 #include "AliRawReaderFile.h"
34 #include "AliRawReaderRoot.h"
36 #include "AliRunLoader.h"
38 ClassImp(AliQADataMakerSteer)
40 //_____________________________________________________________________________
41 AliQADataMakerSteer::AliQADataMakerSteer(const char* gAliceFilename, const char * name, const char * title) :
48 fGAliceFileName(gAliceFilename),
52 fRawReaderDelete(kTRUE),
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 ;
65 //_____________________________________________________________________________
66 AliQADataMakerSteer::AliQADataMakerSteer(const AliQADataMakerSteer & qas) :
69 fDetectors(qas.fDetectors),
73 fGAliceFileName(qas.fGAliceFileName),
74 fRunNumber(qas.fRunNumber),
75 fNumberOfEvents(qas.fNumberOfEvents),
77 fRawReaderDelete(kTRUE),
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] ;
88 //_____________________________________________________________________________
89 AliQADataMakerSteer & AliQADataMakerSteer::operator = (const AliQADataMakerSteer & qas)
91 // assignment operator
92 this->~AliQADataMakerSteer() ;
93 new(this) AliQADataMakerSteer(qas) ;
97 //_____________________________________________________________________________
98 AliQADataMakerSteer::~AliQADataMakerSteer()
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 ;
112 if (fRawReaderDelete) {
119 //_____________________________________________________________________________
120 Bool_t AliQADataMakerSteer::DoIt(const AliQA::TASKINDEX taskIndex)
122 // Runs all the QA data Maker for every detector
125 // Fill QA data in event loop
126 for (UInt_t iEvent = 0 ; iEvent < fNumberOfEvents ; iEvent++) {
128 AliDebug(1, Form("processing event %d", iEvent));
129 if ( taskIndex == AliQA::kRAWS ) {
130 if ( !fRawReader->NextEvent() )
132 } else if ( taskIndex == AliQA::kESDS ) {
133 fESDTree->GetEntry(iEvent) ;
135 fRunLoader->GetEvent(iEvent);
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++) {
144 AliModule* det = static_cast<AliModule*>(detArray->FindObject(AliQA::GetDetName(iDet))) ;
145 if (!det || !det->IsActive())
148 if (!IsSelected(AliQA::GetDetName(iDet)))
150 AliQADataMaker * qadm = GetQADataMaker(iDet) ;
154 if ( qadm->IsCycleDone() ) {
155 qadm->EndOfCycle(AliQA::kRAWS) ;
156 qadm->StartOfCycle(AliQA::kRAWS) ;
161 qadm->Exec(taskIndex, fRawReader) ;
164 GetLoader(iDet)->LoadHits() ;
165 data = GetLoader(iDet)->TreeH() ;
167 AliWarning(Form(" Hit Tree not found for %s", AliQA::GetDetName(iDet))) ;
169 qadm->Exec(taskIndex, data) ;
172 case AliQA::kSDIGITS :
173 GetLoader(iDet)->LoadSDigits() ;
174 data = GetLoader(iDet)->TreeS() ;
176 AliWarning(Form(" SDigit Tree not found for %s", AliQA::GetDetName(iDet))) ;
178 qadm->Exec(taskIndex, data) ;
181 case AliQA::kDIGITS :
182 GetLoader(iDet)->LoadDigits() ;
183 data = GetLoader(iDet)->TreeD() ;
185 AliWarning(Form(" Digit Tree not found for %s", AliQA::GetDetName(iDet))) ;
187 qadm->Exec(taskIndex, data) ;
190 case AliQA::kRECPOINTS :
191 GetLoader(iDet)->LoadRecPoints() ;
192 data = GetLoader(iDet)->TreeR() ;
194 AliWarning(Form("RecPoints not found for %s", AliQA::GetDetName(iDet))) ;
196 qadm->Exec(taskIndex, data) ;
199 case AliQA::kTRACKSEGMENTS :
201 case AliQA::kRECPARTICLES :
204 qadm->Exec(taskIndex, fESD) ;
211 // Save QA data for all detectors
212 rv = Finish(taskIndex) ;
216 //_____________________________________________________________________________
217 AliLoader * AliQADataMakerSteer::GetLoader(Int_t iDet)
219 // get the loader for a detector
221 TString detName = AliQA::GetDetName(iDet) ;
222 fLoader[iDet] = fRunLoader->GetLoader(detName + "Loader");
224 return fLoader[iDet] ;
226 // load the QA data maker object
227 TPluginManager* pluginManager = gROOT->GetPluginManager() ;
228 TString loaderName = "Ali" + detName + "Loader" ;
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 + "()") ;
240 pluginManager->AddHandler("AliLoader", detName, loaderName, detName, loaderName + "()") ;
242 pluginHandler = pluginManager->FindHandler("AliLoader", detName) ;
244 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
245 loader = (AliLoader *) pluginHandler->ExecPlugin(0) ;
248 fLoader[iDet] = loader ;
252 //_____________________________________________________________________________
253 AliQADataMaker * AliQADataMakerSteer::GetQADataMaker(Int_t iDet)
255 // get the quality assurance data maker for a detector
257 if (fQADataMaker[iDet])
258 return fQADataMaker[iDet] ;
260 AliQADataMaker * qadm = NULL ;
263 // load the QA data maker object
264 TPluginManager* pluginManager = gROOT->GetPluginManager() ;
265 TString detName = AliQA::GetDetName(iDet) ;
266 TString qadmName = "Ali" + detName + "QADataMaker" ;
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 + "()") ;
277 pluginManager->AddHandler("AliQADataMaker", detName, qadmName, detName, qadmName + "()") ;
279 pluginHandler = pluginManager->FindHandler("AliQADataMaker", detName) ;
281 if (pluginHandler && (pluginHandler->LoadPlugin() == 0)) {
282 qadm = (AliQADataMaker *) pluginHandler->ExecPlugin(0) ;
285 fQADataMaker[iDet] = qadm ;
290 //_____________________________________________________________________________
291 Bool_t AliQADataMakerSteer::Init(const AliQA::TASKINDEX taskIndex, const char * input )
293 // Initialize the event source and QA data makers
295 if (taskIndex == AliQA::kRAWS) {
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);
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() ;
323 AliError("AliESDs.root not found") ;
327 if ( !InitRunLoader() ) {
328 AliError("Run Loader not found") ;
330 //if (fRunLoader->GetHeader())
331 // fRunNumber = fRunLoader->GetHeader()->GetRun() ;
334 fNumberOfEvents = fRunLoader->GetNumberOfEvents() ;
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) ;
342 AliWarning(Form("AliQADataMaker not found for %s", AliQA::GetDetName(iDet))) ;
344 AliInfo(Form("Data Maker found for %s", qadm->GetName())) ;
345 qadm->Init(taskIndex, fRunNumber, GetQACycles(iDet)) ;
346 qadm->StartOfCycle(taskIndex, fCycleSame) ;
354 //_____________________________________________________________________________
355 Bool_t AliQADataMakerSteer::IsSelected(const char * det)
357 // check whether detName is contained in detectors
358 // if yes, it is removed from detectors
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 ")) {
370 // search for the given detector
372 if ((fDetectors.CompareTo(detName) == 0) ||
373 fDetectors.BeginsWith(detName+" ") ||
374 fDetectors.EndsWith(" "+detName) ||
375 fDetectors.Contains(" "+detName+" ")) {
376 // fDetectors.ReplaceAll(detName, "");
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);
391 //_____________________________________________________________________________
392 Bool_t AliQADataMakerSteer::InitRunLoader()
394 // get or create the run loader
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)))
406 TString detName = AliQA::GetDetName(iDet) ;
407 if (detName == "HLT")
409 if (libs.Contains("lib" + detName + "base.so"))
411 gSystem->Load("lib" + detName + "base.so");
413 fRunLoader = AliRunLoader::Open(fGAliceFileName.Data());
415 AliError(Form("no run loader found in file %s", fGAliceFileName.Data()));
418 fRunLoader->CdGAFile();
419 if (fRunLoader->LoadgAlice() == 0) {
420 gAlice = fRunLoader->GetAliRun();
424 AliError(Form("no gAlice object found in file %s", fGAliceFileName.Data()));
428 } else { // galice.root does not exist
429 AliError(Form("the file %s does not exist", fGAliceFileName.Data()));
436 //_____________________________________________________________________________
437 Bool_t AliQADataMakerSteer::Finish(const AliQA::TASKINDEX taskIndex)
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) ;
444 qadm->EndOfCycle(taskIndex) ;
451 //_____________________________________________________________________________
452 Bool_t AliQADataMakerSteer::Merge()
454 // Merge all the cycles from all detectors in one single file per run
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] ;
464 in >> file[index++] ;
468 Int_t previousRun = -1 ;
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 ;
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)) ;
486 for (Int_t irun = 0 ; irun < runIndex ; irun++) {
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++) {
493 sprintf(pattern, "QA.%d.", runIndex-1) ;
494 TString tmp(file[ifile]) ;
495 if (tmp.Contains(pattern))
496 merger.AddFile(tmp) ;
504 //_____________________________________________________________________________
505 void AliQADataMakerSteer::Reset()
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 ;
519 if (fRawReaderDelete) {
524 fCycleSame = kFALSE ;
528 fNumberOfEvents = 0 ;
531 //_____________________________________________________________________________
532 Bool_t AliQADataMakerSteer::Run(const char * detectors, AliRawReader * rawReader)
534 //Runs all the QA data Maker for Raws only
535 fRawReader = rawReader ;
536 fRawReaderDelete = kFALSE ;
538 fDetectors = detectors ;
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) ;
545 AliWarning(Form("AliQADataMaker not found for %s", AliQA::GetDetName(iDet))) ;
547 AliInfo(Form("Data Maker found for %s", qadm->GetName())) ;
548 qadm->Init(AliQA::kRAWS, fRunNumber, GetQACycles(iDet)) ;
549 qadm->StartOfCycle(AliQA::kRAWS, fCycleSame) ;
555 return DoIt(AliQA::kRAWS) ;
558 //_____________________________________________________________________________
559 Bool_t AliQADataMakerSteer::Run(const char * detectors, const AliQA::TASKINDEX taskIndex, const char * fileName )
561 // Runs all the QA data Maker for every detector
564 fDetectors = detectors ;
566 if ( !Init(taskIndex, fileName) )
569 rv = DoIt(taskIndex) ;