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