TRD/qaRec/run.C

   1 // Steer TRD QA train for Reconstruction (Clusterizer, Tracking and PID).
   2 //
   3 // Usage:
   4 //   run.C(tasks, files)
   5 //   tasks : "ALL" or one/more of the following:
   6 //     "EFF"  : TRD Tracking Efficiency
   7 //     "EFFC" : TRD Tracking Efficiency Combined (barrel + stand alone) - only in case of simulations
   8 //     "RES"  : TRD tracking Resolution
   9 //     "CLRES": clusters Resolution
  10 //     "CAL"  : TRD calibration
  11 //     "ALGN" : TRD alignment
  12 //     "PID"  : TRD PID - pion efficiency
  13 //     "PIDR" : TRD PID - reference data
  14 //     "DET"  : Basic TRD Detector checks
  15 //     "NOFR" : Data set does not have AliESDfriends.root
  16 //     "NOMC" : Data set does not have Monte Carlo Informations (real data), so all tasks which rely
  17 //              on MC information are switched off
  18 //
  19 // In compiled mode :
  20 // Don't forget to load first the libraries
  21 // gSystem->Load("libMemStat.so")
  22 // gSystem->Load("libMemStatGui.so")
  23 // gSystem->Load("libANALYSIS.so")
  24 // gSystem->Load("libTRDqaRec.so")
  25 // gSystem->Load("libNetx.so") ;
  26 // gSystem->Load("libRAliEn.so");
  27 //
  28 // Authors:
  29 //   Alex Bercuci (A.Bercuci@gsi.de)
  30 //   Markus Fasel (m.Fasel@gsi.de)
  31
  32 #ifndef __CINT__
  33 #include <Riostream.h>
  34
  35 #include "TStopwatch.h"
  36 #include "TMemStat.h"
  37 #include "TMemStatViewerGUI.h"
  38
  39 #include "TROOT.h"
  40 #include "TSystem.h"
  41 #include "TError.h"
  42 #include "TChain.h"
  43 #include "TGrid.h"
  44 #include "TAlienCollection.h"
  45 #include "TGridCollection.h"
  46 #include "TGridResult.h"
  47
  48 #include "AliMagFCheb.h"
  49 #include "AliTracker.h"
  50 #include "AliLog.h"
  51 #include "AliCDBManager.h"
  52 #include "AliGeomManager.h"
  53 #include "AliAnalysisManager.h"
  54 #include "AliAnalysisDataContainer.h"
  55 #include "AliMCEventHandler.h"
  56 #include "AliESDInputHandler.h"
  57
  58 #include "TRD/AliTRDtrackerV1.h"
  59 #include "TRD/AliTRDcalibDB.h"
  60 #include "TRD/qaRec/AliTRDtrackInfo/AliTRDeventInfo.h"
  61 #include "TRD/qaRec/AliTRDtrackInfoGen.h"
  62 #include "TRD/qaRec/AliTRDtrackingEfficiency.h"
  63 #include "TRD/qaRec/AliTRDtrackingEfficiencyCombined.h"
  64 #include "TRD/qaRec/AliTRDtrackingResolution.h"
  65 #include "TRD/qaRec/AliTRDcalibration.h"
  66 #include "TRD/qaRec/AliTRDalignmentTask.h"
  67 #include "TRD/qaRec/AliTRDpidChecker.h"
  68 #include "TRD/qaRec/AliTRDpidRefMaker.h"
  69 #include "TRD/qaRec/AliTRDcheckDetector.h"
  70 #include "TRD/qaRec/AliTRDclusterResolution.h"
  71 #endif
  72
  73 #include "run.h"
  74
  75 Bool_t MEM = kFALSE;
  76
  77 TChain* MakeChainLST(const char* filename = 0x0);
  78 TChain* MakeChainXML(const char* filename = 0x0);
  79 void run(Char_t *tasks="ALL", const Char_t *files=0x0)
  80 {
  81   TMemStat *mem = 0x0;
  82   if(MEM){
  83     gSystem->Load("libMemStat.so");
  84     gSystem->Load("libMemStatGui.so");
  85     mem = new TMemStat("new, gnubuildin");
  86     mem->AddStamp("Start");
  87   }
  88
  89   TStopwatch timer;
  90   timer.Start();
  91
  92   if(gSystem->Load("libANALYSIS.so")<0) return;
  93   if(gSystem->Load("libTRDqaRec.so")<0) return;
  94
  95   // DB INITIALIZATION
  96   // initialize OCDB manager
  97   AliCDBManager *cdbManager = AliCDBManager::Instance();
  98   cdbManager->SetDefaultStorage("local://$ALICE_ROOT");
  99   cdbManager->SetRun(0);
 100   cdbManager->SetCacheFlag(kFALSE);
 101
 102   // initialize magnetic field. TODO We should use the GRP !
 103   AliMagFCheb *field = 0x0;
 104   field = new AliMagFCheb("Maps","Maps", 2, 1., 10., AliMagFCheb::k5kG, kTRUE,"$(ALICE_ROOT)/data/maps/mfchebKGI_sym.root");
 105   //field = new AliMagFCheb("Maps","Maps", 2, 0., 10., AliMagFCheb::k2kG);
 106   AliTracker::SetFieldMap(field, kTRUE);
 107
 108   // initialize TRD settings
 109   AliTRDcalibDB *cal = AliTRDcalibDB::Instance();
 110   AliTRDtrackerV1::SetNTimeBins(cal->GetNumberOfTimeBins());
 111   AliGeomManager::LoadGeometry();
 112
 113   Bool_t fHasMCdata = kTRUE;
 114   Bool_t fHasFriends = kTRUE;
 115   TObjArray *tasksArray = TString(tasks).Tokenize(" ");
 116
 117   Int_t fSteerTask = 0; SETBIT(fSteerTask, kInfoGen);
 118   for(Int_t isel = 0; isel < tasksArray->GetEntriesFast(); isel++){
 119     TString s = (dynamic_cast<TObjString *>(tasksArray->UncheckedAt(isel)))->String();
 120     if(s.CompareTo("ALL") == 0){
 121       for(Int_t itask = 1; itask < NQATASKS; itask++) SETBIT(fSteerTask, itask);
 122       continue;
 123     } else if(s.CompareTo("NOFR") == 0){
 124       fHasFriends = kFALSE;
 125     } else if(s.CompareTo("NOMC") == 0){
 126       fHasMCdata = kFALSE;
 127     } else {
 128       Bool_t foundOpt = kFALSE;
 129       for(Int_t itask = 1; itask < NTRDTASKS; itask++){
 130         if(s.CompareTo(fgkTRDtaskOpt[itask]) != 0) continue;
 131         SETBIT(fSteerTask, itask);
 132         foundOpt = kTRUE;
 133         break;
 134       }
 135       if(!foundOpt) Info("run.C", Form("Task %s not implemented (yet).", s.Data()));
 136     }
 137   }
 138   // extra rules for calibration tasks
 139   if(TSTBIT(fSteerTask, kClErrParam)) SETBIT(fSteerTask, kResolution);
 140   if(TSTBIT(fSteerTask, kPIDRefMaker)) SETBIT(fSteerTask, kPIDChecker);
 141   if(TSTBIT(fSteerTask, kAlignment)) SETBIT(fSteerTask, kResolution);
 142
 143   // define task list pointers;
 144   AliTRDrecoTask *taskPtr[NTRDTASKS], *task = 0x0;
 145   memset(taskPtr, 0, NTRDTASKS*sizeof(AliAnalysisTask*));
 146
 147   //____________________________________________//
 148   // DEFINE DATA CHAIN
 149   TChain *chain = 0x0;
 150   if(!files) chain = MakeChainLST();
 151   else{
 152     TString fn(files);
 153     if(fn.EndsWith("xml")) chain = MakeChainXML(files);
 154     else chain = MakeChainLST(files);
 155   }
 156   if(!chain) return;
 157
 158   //chain->SetBranchStatus("*", 0);
 159   chain->SetBranchStatus("*FMD*",0);
 160   chain->SetBranchStatus("*Calo*",0);
 161   chain->SetBranchStatus("Tracks", 1);
 162   chain->SetBranchStatus("ESDfriend*",1);
 163   chain->Lookup();
 164   chain->GetListOfFiles()->Print();
 165   printf("\n ----> CHAIN HAS %d ENTRIES <----\n\n", (Int_t)chain->GetEntries());
 166
 167   AliLog::SetGlobalLogLevel(AliLog::kError);
 168
 169   //____________________________________________
 170   // Make the analysis manager
 171   AliAnalysisManager *mgr = new AliAnalysisManager("TRD QA Reconstruction Manager");
 172   //mgr->SetSpecialOutputLocation(source); // To Be Changed
 173   AliVEventHandler *esdH = 0x0, *mcH = 0x0;
 174   mgr->SetInputEventHandler(esdH = new AliESDInputHandler);
 175   if(fHasMCdata) mgr->SetMCtruthEventHandler(mcH = new AliMCEventHandler());
 176   //mgr->SetDebugLevel(10);
 177
 178   //____________________________________________
 179   // TRD track summary generator
 180   mgr->AddTask(task = new AliTRDtrackInfoGen());
 181   taskPtr[(Int_t)kInfoGen] = task;
 182   task->SetDebugLevel(1);
 183   task->SetMCdata(fHasMCdata);
 184   // Create containers for input/output
 185   AliAnalysisDataContainer *cinput1 = mgr->CreateContainer("data", TChain::Class(), AliAnalysisManager::kInputContainer);
 186   AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("trackInfo", TObjArray::Class(), AliAnalysisManager::kExchangeContainer);
 187   AliAnalysisDataContainer *coutput1a = mgr->CreateContainer("eventInfo", AliTRDeventInfo::Class(), AliAnalysisManager::kExchangeContainer);
 188   mgr->ConnectInput( task, 0, cinput1);
 189   mgr->ConnectOutput(task, 0, coutput1);
 190   mgr->ConnectOutput(task, 1, coutput1a);
 191
 192   //____________________________________________
 193   // TRD detector checker
 194         if(TSTBIT(fSteerTask, kCheckDetector)){
 195     mgr->AddTask(task = new AliTRDcheckDetector());
 196     taskPtr[(Int_t)kCheckDetector] = task;
 197     task->SetDebugLevel(4);
 198     task->SetMCdata(fHasMCdata);
 199
 200     // Create containers for input/output
 201     mgr->ConnectInput( task, 0, coutput1);
 202     mgr->ConnectInput( task, 1, coutput1a);
 203     mgr->ConnectOutput(task, 0, mgr->CreateContainer(task->GetName(), TObjArray::Class(), AliAnalysisManager::kOutputContainer, Form("TRD.Task%s.root", task->GetName())));
 204   }
 205
 206   //____________________________________________
 207   // TRD barrel tracking efficiency
 208   if(fHasMCdata && TSTBIT(fSteerTask, kTrackingEff)){
 209     mgr->AddTask(task = new AliTRDtrackingEfficiency());
 210     taskPtr[(Int_t)kTrackingEff] = task;
 211     task->SetDebugLevel(0);
 212
 213     //Create containers for input/output
 214     mgr->ConnectInput( task, 0, coutput1);
 215     mgr->ConnectOutput(task, 0, mgr->CreateContainer(task->GetName(), TObjArray::Class(), AliAnalysisManager::kOutputContainer, Form("TRD.Task%s.root", task->GetName())));
 216   }
 217
 218   //____________________________________________
 219   // TRD combined tracking efficiency
 220   if(fHasMCdata && TSTBIT(fSteerTask, kTrackingEffMC)){
 221     mgr->AddTask(task = new AliTRDtrackingEfficiencyCombined());
 222     taskPtr[(Int_t)kTrackingEffMC] = task;
 223     task->SetDebugLevel(0);
 224
 225     // Create containers for input/output
 226     mgr->ConnectInput( task, 0, coutput1);
 227     mgr->ConnectOutput(task, 0, mgr->CreateContainer(task->GetName(), TObjArray::Class(), AliAnalysisManager::kOutputContainer, Form("TRD.Task%s.root", task->GetName())));
 228   }
 229
 230   //____________________________________________
 231   // TRD tracking resolution
 232   if(TSTBIT(fSteerTask, kResolution)){
 233     mgr->AddTask(task = new AliTRDtrackingResolution());
 234     taskPtr[(Int_t)kResolution] = task;
 235     task->SetMCdata(fHasMCdata);
 236     task->SetPostProcess(kFALSE);
 237     task->SetDebugLevel(1);
 238
 239     // Create containers for input/output
 240     mgr->ConnectInput( task, 0, coutput1);
 241     mgr->ConnectOutput(task, 0, mgr->CreateContainer(task->GetName(), TObjArray::Class(), AliAnalysisManager::kOutputContainer, Form("TRD.Task%s.root", task->GetName())));
 242
 243     // Create output containers for calibration tasks
 244     const Int_t nc = 4;
 245     const Char_t *cn[nc] = {"Cl", "Trklt", "MC_Cl", "MC_Trklt"};
 246     AliAnalysisDataContainer *co[nc];
 247     for(Int_t ic = 0; ic<nc; ic++){
 248       co[ic] = mgr->CreateContainer(Form("%s%s", task->GetName(), cn[ic]), TObjArray::Class(), AliAnalysisManager::kExchangeContainer);
 249       mgr->ConnectOutput(task, 1+ic, co[ic]);
 250     }
 251
 252     // test reconstruction calibration plugin
 253     if(TSTBIT(fSteerTask, kClErrParam)){
 254       mgr->AddTask(task = new AliTRDclusterResolution());
 255       taskPtr[(Int_t)kClErrParam] = task;
 256       ((AliTRDclusterResolution*)task)->SetExB();
 257       mgr->ConnectInput(task, 0, co[0]);
 258       mgr->ConnectOutput(task, 0, mgr->CreateContainer(task->GetName(), TObjArray::Class(), AliAnalysisManager::kOutputContainer, Form("TRD.Task%s.root", task->GetName())));
 259
 260       mgr->AddTask(task = new AliTRDclusterResolution("ClErrParamMC"));
 261       taskPtr[(Int_t)kClErrParam+1] = task;
 262       ((AliTRDclusterResolution*)task)->SetExB();
 263       mgr->ConnectInput(task, 0, co[2]);
 264       mgr->ConnectOutput(task, 0, mgr->CreateContainer(task->GetName(), TObjArray::Class(), AliAnalysisManager::kOutputContainer, Form("TRD.Task%s.root", task->GetName())));
 265     }
 266   }
 267
 268   //____________________________________________
 269   // TRD calibration
 270   if(TSTBIT(fSteerTask, kCalibration)){
 271     mgr->AddTask(task = new AliTRDcalibration());
 272     taskPtr[(Int_t)kCalibration] = task;
 273     ((AliTRDcalibration*)task)->SetLow(0);
 274     ((AliTRDcalibration*)task)->SetHigh(30);
 275     ((AliTRDcalibration*)task)->SetFillZero(kFALSE);
 276     task->SetDebugLevel(0);
 277
 278     // Create containers for input/output
 279     mgr->ConnectInput(task, 0, coutput1);
 280     mgr->ConnectOutput(task, 0, mgr->CreateContainer(task->GetName(), TObjArray::Class(), AliAnalysisManager::kOutputContainer, Form("TRD.Task%s.root", task->GetName())));
 281   }
 282
 283   //____________________________________________
 284   // TRD alignment
 285   if(TSTBIT(fSteerTask, kAlignment)){
 286     mgr->AddTask(task = new AliTRDalignmentTask());
 287     taskPtr[(Int_t)kAlignment] = task;
 288     task->SetDebugLevel(0);
 289
 290     // Create containers for input/output
 291     mgr->ConnectInput(task, 0, coutput1);
 292     mgr->ConnectOutput(task, 0, mgr->CreateContainer(Form("h%s", task->GetName()), TObjArray::Class(), AliAnalysisManager::kExchangeContainer));
 293
 294     mgr->ConnectOutput(task, 1, mgr->CreateContainer(task->GetName(), TTree::Class(), AliAnalysisManager::kOutputContainer, Form("TRD.Task%s.root", task->GetName())));
 295   }
 296
 297   //____________________________________________
 298   // TRD PID
 299   if(TSTBIT(fSteerTask, kPIDChecker)){
 300     mgr->AddTask(task = new AliTRDpidChecker());
 301     taskPtr[(Int_t)kPIDChecker] = task;
 302     task->SetDebugLevel(0);
 303     task->SetMCdata(fHasMCdata);
 304
 305     // Create containers for input/output
 306     mgr->ConnectInput( task, 0, coutput1);
 307     mgr->ConnectOutput(task, 0, mgr->CreateContainer(task->GetName(), TObjArray::Class(), AliAnalysisManager::kOutputContainer, Form("TRD.Task%s.root", task->GetName())));
 308
 309     //____________________________________________
 310     // TRD pid reference
 311     if(TSTBIT(fSteerTask, kPIDRefMaker)){
 312       mgr->AddTask(task = new AliTRDpidRefMaker());
 313       taskPtr[(Int_t)kPIDRefMaker] = task;
 314       task->SetDebugLevel(0);
 315       task->SetMCdata(fHasMCdata);
 316
 317       // Create containers for input/output
 318       mgr->ConnectInput( task, 0, coutput1);
 319       mgr->ConnectOutput(task, 0, mgr->CreateContainer(task->GetName(), TObjArray::Class(), AliAnalysisManager::kOutputContainer, Form("TRD.Task%s.root", task->GetName())));
 320       mgr->ConnectOutput(task, 1, mgr->CreateContainer(Form("%sNN", task->GetName()), TTree::Class(), AliAnalysisManager::kOutputContainer, Form("TRD.Task%sNN.root", task->GetName())));
 321       mgr->ConnectOutput(task, 2, mgr->CreateContainer(Form("%sLQ", task->GetName()), TTree::Class(), AliAnalysisManager::kOutputContainer, Form("TRD.Task%sLQ.root", task->GetName())));
 322     }
 323   }
 324
 325
 326   if (!mgr->InitAnalysis()) return;
 327   printf("\n\tRUNNING TRAIN FOR TASKS:\n");
 328   for(Int_t itask = 1; itask < NTRDTASKS; itask++){
 329     if(TSTBIT(fSteerTask, itask)) printf("\t   %s [%s]\n",  taskPtr[itask]->GetName(), taskPtr[itask]->GetTitle());
 330   }
 331   printf("\n\n");
 332   //mgr->PrintStatus();
 333
 334   mgr->StartAnalysis("local",chain);
 335
 336   timer.Stop();
 337   timer.Print();
 338
 339   cal->Terminate();
 340   delete field;
 341   delete cdbManager;
 342   for(Int_t it=NTRDTASKS; it--; ){
 343     if(taskPtr[it]){
 344       printf("Cleaning %s [%s] ...\n", fgkTRDtaskClassName[it], taskPtr[it]->GetTitle());
 345       delete taskPtr[it];
 346     }
 347   }
 348   if(mcH) delete mcH;
 349   delete esdH;
 350   delete mgr;
 351   delete chain;
 352   if(MEM) delete mem;
 353   if(MEM) TMemStatViewerGUI::ShowGUI();
 354 }
 355
 356 //____________________________________________
 357 TChain* MakeChainLST(const char* filename)
 358 {
 359   // Create the chain
 360   TChain* chain = new TChain("esdTree");
 361
 362   if(!filename){
 363     chain->Add(Form("%s/AliESDs.root", gSystem->pwd()));
 364     return chain;
 365   }
 366
 367
 368   // read ESD files from the input list.
 369   ifstream in;
 370   in.open(filename);
 371   TString esdfile;
 372   while(in.good()) {
 373     in >> esdfile;
 374     if (!esdfile.Contains("root")) continue; // protection
 375     chain->Add(esdfile.Data());
 376   }
 377
 378   in.close();
 379
 380   return chain;
 381 }
 382
 383 //____________________________________________
 384 TChain* MakeChainXML(const char* xmlfile)
 385 {
 386   if (!TFile::Open(xmlfile)) {
 387     Error("MakeChainXML", Form("No file %s was found", xmlfile));
 388     return 0x0;
 389   }
 390
 391   if(gSystem->Load("libNetx.so")<0) return 0x0;
 392   if(gSystem->Load("libRAliEn.so")<0) return 0x0;
 393   TGrid::Connect("alien://") ;
 394
 395   TGridCollection *collection = (TGridCollection*) TAlienCollection::Open(xmlfile);
 396   if (!collection) {
 397     Error("MakeChainXML", Form("No collection found in %s", xmlfile)) ;
 398     return 0x0;
 399   }
 400   //collection->CheckIfOnline();
 401
 402   TGridResult* result = collection->GetGridResult("",0 ,0);
 403   if(!result->GetEntries()){
 404     Error("MakeChainXML", Form("No entries found in %s", xmlfile)) ;
 405     return 0x0;
 406   }
 407   // Makes the ESD chain
 408   TChain* chain = new TChain("esdTree");
 409   for (Int_t idx = 0; idx < result->GetEntries(); idx++) {
 410     chain->Add(result->GetKey(idx, "turl"));
 411   }
 412   return chain;
 413 }