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