CORRFW/test/AliCFSingleTrackTask.C

   1 //DEFINITION OF A FEW CONSTANTS
   2 const Double_t ymin  = -1.0 ;
   3 const Double_t ymax  =  1.0 ;
   4 const Double_t ptmin =  0.0 ;
   5 const Double_t ptmax =  8.0 ;
   6 const Int_t    mintrackrefsTPC = 2 ;
   7 const Int_t    mintrackrefsITS = 3 ;
   8 const Int_t    charge  = 1 ;
   9 const Int_t    PDG = 2212;
  10 const Int_t    minclustersTPC = 50 ;
  11 //----------------------------------------------------
  12
  13 Bool_t AliCFSingleTrackTask(
  14                             const Bool_t useGrid = 1,
  15                             const Bool_t readAOD = 0,
  16                             const Bool_t readTPCTracks = 0,
  17                             const char * kTagXMLFile="wn.xml" // XML file containing tags
  18                             )
  19 {
  20
  21   TBenchmark benchmark;
  22   benchmark.Start("AliSingleTrackTask");
  23
  24   AliLog::SetGlobalDebugLevel(0);
  25
  26   Load() ; //load the required libraries
  27
  28   TChain * analysisChain ;
  29
  30   if (useGrid) { //data located on AliEn
  31     TGrid::Connect("alien://") ;    //  Create an AliRunTagCuts and an AliEventTagCuts Object
  32                                     //  and impose some selection criteria
  33     AliRunTagCuts      *runCuts   = new AliRunTagCuts();
  34     AliEventTagCuts    *eventCuts = new AliEventTagCuts();
  35     AliLHCTagCuts      *lhcCuts   = new AliLHCTagCuts();
  36     AliDetectorTagCuts *detCuts   = new AliDetectorTagCuts();
  37     eventCuts->SetMultiplicityRange(0,2000);
  38
  39     //  Create an AliTagAnalysis Object and chain the tags
  40     AliTagAnalysis   *tagAna = new AliTagAnalysis();
  41     if (readAOD) tagAna->SetType("AOD");  //for aliroot > v4-05
  42     else         tagAna->SetType("ESD");  //for aliroot > v4-05
  43     TAlienCollection *coll   = TAlienCollection::Open(kTagXMLFile);
  44     TGridResult      *tagResult = coll->GetGridResult("",0,0);
  45     tagResult->Print();
  46     tagAna->ChainGridTags(tagResult);
  47
  48     //  Create a new esd chain and assign the chain that is returned by querying the tags
  49     analysisChain = tagAna->QueryTags(runCuts,lhcCuts,detCuts,eventCuts);
  50   }
  51
  52   else {// local data
  53     //here put your input data path
  54     printf("\n\nRunning on local file, please check the path\n\n");
  55
  56     if (readAOD) {
  57       analysisChain = new TChain("aodTree");
  58       analysisChain->Add("your_data_path/001/AliAOD.root");
  59       analysisChain->Add("your_data_path/002/AliAOD.root");
  60     }
  61     else {
  62       analysisChain = new TChain("esdTree");
  63       analysisChain->Add("your_data_path/001/AliESDs.root");
  64       analysisChain->Add("your_data_path/002/AliESDs.root");
  65     }
  66   }
  67
  68
  69   Info("AliCFSingleTrackTask",Form("CHAIN HAS %d ENTRIES",(Int_t)analysisChain->GetEntries()));
  70
  71   //CONTAINER DEFINITION
  72   Info("AliCFSingleTrackTask","SETUP CONTAINER");
  73   //the sensitive variables (2 in this example), their indices
  74   UInt_t ipt = 0;
  75   UInt_t iy  = 1;
  76   //Setting up the container grid...
  77   UInt_t nstep = 4 ; //number of selection steps MC
  78   const Int_t nvar   = 2 ; //number of variables on the grid:pt,y
  79   const Int_t nbin1  = 8 ; //bins in pt
  80   const Int_t nbin2  = 8 ; //bins in y
  81
  82   //arrays for the number of bins in each dimension
  83   Int_t iBin[nvar];
  84   iBin[0]=nbin1;
  85   iBin[1]=nbin2;
  86
  87   //arrays for lower bounds :
  88   Double_t *binLim1=new Double_t[nbin1+1];
  89   Double_t *binLim2=new Double_t[nbin2+1];
  90
  91   //values for bin lower bounds
  92   for(Int_t i=0; i<=nbin1; i++) binLim1[i]=(Double_t)ptmin + (ptmax-ptmin)/nbin1*(Double_t)i ;
  93   for(Int_t i=0; i<=nbin2; i++) binLim2[i]=(Double_t)ymin  + (ymax-ymin)  /nbin2*(Double_t)i ;
  94   //one "container" for MC
  95   AliCFContainer* container = new AliCFContainer("container","container for tracks",nstep,nvar,iBin);
  96   //setting the bin limits
  97   container -> SetBinLimits(ipt,binLim1);
  98   container -> SetBinLimits(iy,binLim2);
  99   container -> SetVarTitle(ipt,"pt");
 100   container -> SetVarTitle(iy, "y");
 101   container -> SetStepTitle(0, "generated");
 102   container -> SetStepTitle(1, "in acceptance");
 103   container -> SetStepTitle(2, "reconstructed");
 104   container -> SetStepTitle(3, "after PID");
 105
 106   // SET TLIST FOR QA HISTOS
 107   TList* qaList = new TList();
 108
 109   //CREATE THE  CUTS -----------------------------------------------
 110
 111   //Event-level cuts:
 112   AliCFEventRecCuts* evtRecCuts = new AliCFEventRecCuts("evtRecCuts","Rec-event cuts");
 113 //   evtRecCuts->SetUseTPCVertex();
 114 //   evtRecCuts->SetRequireVtxCuts(kTRUE);
 115 //   evtRecCuts->SetVertexNContributors(-2,5);
 116   evtRecCuts->SetQAOn(qaList);
 117
 118   // Gen-Level kinematic cuts
 119   AliCFTrackKineCuts *mcKineCuts = new AliCFTrackKineCuts("mcKineCuts","MC-level kinematic cuts");
 120   mcKineCuts->SetPtRange(ptmin,ptmax);
 121   mcKineCuts->SetRapidityRange(ymin,ymax);
 122   mcKineCuts->SetChargeMC(charge);
 123   mcKineCuts->SetQAOn(qaList);
 124
 125   //Particle-Level cuts:
 126   AliCFParticleGenCuts* mcGenCuts = new AliCFParticleGenCuts("mcGenCuts","MC particle generation cuts");
 127   mcGenCuts->SetRequireIsPrimary();
 128   mcGenCuts->SetRequirePdgCode(PDG,/*absolute=*/kTRUE);
 129   mcGenCuts->SetQAOn(qaList);
 130
 131   //Acceptance Cuts
 132   AliCFAcceptanceCuts *mcAccCuts = new AliCFAcceptanceCuts("mcAccCuts","MC acceptance cuts");
 133   mcAccCuts->SetMinNHitITS(mintrackrefsITS);
 134   mcAccCuts->SetMinNHitTPC(mintrackrefsTPC);
 135   mcAccCuts->SetQAOn(qaList);
 136
 137   // Rec-Level kinematic cuts
 138   AliCFTrackKineCuts *recKineCuts = new AliCFTrackKineCuts("recKineCuts","rec-level kine cuts");
 139   recKineCuts->SetPtRange(ptmin,ptmax);
 140   recKineCuts->SetRapidityRange(ymin,ymax);
 141   recKineCuts->SetChargeRec(charge);
 142   recKineCuts->SetQAOn(qaList);
 143
 144   AliCFTrackQualityCuts *recQualityCuts = new AliCFTrackQualityCuts("recQualityCuts","rec-level quality cuts");
 145   if (!readAOD)       {
 146 //     recQualityCuts->SetMinNClusterTRD(0);
 147 //     recQualityCuts->SetMaxChi2PerClusterTRD(10.);
 148   }
 149   recQualityCuts->SetStatus(AliESDtrack::kTPCrefit);
 150   recQualityCuts->SetQAOn(qaList);
 151
 152   AliCFTrackIsPrimaryCuts *recIsPrimaryCuts = new AliCFTrackIsPrimaryCuts("recIsPrimaryCuts","rec-level isPrimary cuts");
 153   if (readAOD) recIsPrimaryCuts->SetAODType(AliAODTrack::kPrimary);
 154   else         recIsPrimaryCuts->SetMaxNSigmaToVertex(3);
 155   recIsPrimaryCuts->SetQAOn(qaList);
 156
 157   AliCFTrackCutPid* cutPID = new AliCFTrackCutPid("cutPID","ESD_PID") ;
 158   int n_species = AliPID::kSPECIES ;
 159   Double_t* prior = new Double_t[n_species];
 160
 161   prior[0] = 0.0244519 ;
 162   prior[1] = 0.0143988 ;
 163   prior[2] = 0.805747  ;
 164   prior[3] = 0.0928785 ;
 165   prior[4] = 0.0625243 ;
 166
 167   cutPID->SetPriors(prior);
 168   cutPID->SetProbabilityCut(0.0);
 169   if (readTPCTracks) cutPID->SetDetectors("TPC");
 170   else               cutPID->SetDetectors("ALL");
 171   if (readAOD) cutPID->SetAODmode(kTRUE );
 172   else         cutPID->SetAODmode(kFALSE);
 173   switch(TMath::Abs(PDG)) {
 174   case 11   : cutPID->SetParticleType(AliPID::kElectron, kTRUE); break;
 175   case 13   : cutPID->SetParticleType(AliPID::kMuon    , kTRUE); break;
 176   case 211  : cutPID->SetParticleType(AliPID::kPion    , kTRUE); break;
 177   case 321  : cutPID->SetParticleType(AliPID::kKaon    , kTRUE); break;
 178   case 2212 : cutPID->SetParticleType(AliPID::kProton  , kTRUE); break;
 179   default   : printf("UNDEFINED PID\n"); break;
 180   }
 181   cutPID->SetQAOn(qaList);
 182
 183   printf("CREATE EVENT LEVEL CUTS\n");
 184   TObjArray* evtList = new TObjArray(0) ;
 185 //   evtList->AddLast(evtRecCuts);
 186
 187   printf("CREATE MC KINE CUTS\n");
 188   TObjArray* mcList = new TObjArray(0) ;
 189   mcList->AddLast(mcKineCuts);
 190   mcList->AddLast(mcGenCuts);
 191
 192   printf("CREATE ACCEPTANCE CUTS\n");
 193   TObjArray* accList = new TObjArray(0) ;
 194   accList->AddLast(mcAccCuts);
 195
 196   printf("CREATE RECONSTRUCTION CUTS\n");
 197   TObjArray* recList = new TObjArray(0) ;
 198   recList->AddLast(recKineCuts);
 199   recList->AddLast(recQualityCuts);
 200   recList->AddLast(recIsPrimaryCuts);
 201
 202   printf("CREATE PID CUTS\n");
 203   TObjArray* fPIDCutList = new TObjArray(0) ;
 204   fPIDCutList->AddLast(cutPID);
 205
 206   //CREATE THE INTERFACE TO CORRECTION FRAMEWORK USED IN THE TASK
 207   printf("CREATE INTERFACE AND CUTS\n");
 208   AliCFManager* man = new AliCFManager() ;
 209
 210   man->SetNStepEvent(1);
 211   man->SetEventCutsList(0,evtList);
 212
 213   man->SetParticleContainer(container);
 214   man->SetParticleCutsList(0,mcList);
 215   man->SetParticleCutsList(1,accList);
 216   man->SetParticleCutsList(2,recList);
 217   man->SetParticleCutsList(3,fPIDCutList);
 218
 219
 220   //CREATE THE TASK
 221   printf("CREATE TASK\n");
 222   // create the task
 223   AliCFSingleTrackTask *task = new AliCFSingleTrackTask("AliSingleTrackTask");
 224   task->SetCFManager(man); //here is set the CF manager
 225   task->SetQAList(qaList);
 226   if (readAOD)       task->SetReadAODData() ;
 227   if (readTPCTracks) task->SetReadTPCTracks();
 228
 229   //SETUP THE ANALYSIS MANAGER TO READ INPUT CHAIN AND WRITE DESIRED OUTPUTS
 230   printf("CREATE ANALYSIS MANAGER\n");
 231   // Make the analysis manager
 232   AliAnalysisManager *mgr = new AliAnalysisManager("TestManager");
 233
 234   if (useGrid) mgr->SetAnalysisType(AliAnalysisManager::kGridAnalysis);
 235   else mgr->SetAnalysisType(AliAnalysisManager::kLocalAnalysis);
 236
 237
 238   AliMCEventHandler*  mcHandler = new AliMCEventHandler();
 239   mgr->SetMCtruthEventHandler(mcHandler);
 240
 241   AliInputEventHandler* dataHandler ;
 242
 243   if   (readAOD) dataHandler = new AliAODInputHandler();
 244   else           dataHandler = new AliESDInputHandler();
 245   mgr->SetInputEventHandler(dataHandler);
 246
 247   // Create and connect containers for input/output
 248
 249   //------ input data ------
 250   AliAnalysisDataContainer *cinput0  = mgr->CreateContainer("cchain0",TChain::Class(),AliAnalysisManager::kInputContainer);
 251
 252   // ----- output data -----
 253
 254   //slot 0 : default output tree (by default handled by AliAnalysisTaskSE)
 255   AliAnalysisDataContainer *coutput0 = mgr->CreateContainer("ctree0", TTree::Class(),AliAnalysisManager::kOutputContainer,"output.root");
 256
 257   //now comes user's output objects :
 258
 259   // output TH1I for event counting
 260   AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("chist0", TH1I::Class(),AliAnalysisManager::kOutputContainer,"output.root");
 261   // output Correction Framework Container (for acceptance & efficiency calculations)
 262   AliAnalysisDataContainer *coutput2 = mgr->CreateContainer("ccontainer0", AliCFContainer::Class(),AliAnalysisManager::kOutputContainer,"output.root");
 263   // output QA histograms
 264   AliAnalysisDataContainer *coutput3 = mgr->CreateContainer("clist0", TList::Class(),AliAnalysisManager::kOutputContainer,"output.root");
 265
 266   cinput0->SetData(analysisChain);
 267
 268   mgr->AddTask(task);
 269   mgr->ConnectInput(task,0,mgr->GetCommonInputContainer());
 270   mgr->ConnectOutput(task,0,coutput0);
 271   mgr->ConnectOutput(task,1,coutput1);
 272   mgr->ConnectOutput(task,2,coutput2);
 273   mgr->ConnectOutput(task,3,coutput3);
 274
 275   printf("READY TO RUN\n");
 276   //RUN !!!
 277   if (mgr->InitAnalysis()) {
 278     mgr->PrintStatus();
 279     mgr->StartAnalysis("local",analysisChain);
 280   }
 281
 282   benchmark.Stop("AliSingleTrackTask");
 283   benchmark.Show("AliSingleTrackTask");
 284
 285   return kTRUE ;
 286 }
 287
 288 void Load() {
 289
 290   //load the required aliroot libraries
 291   gSystem->Load("libANALYSIS") ;
 292   gSystem->Load("libANALYSISalice") ;
 293   gSystem->Load("libCORRFW.so") ;
 294
 295   //compile online the task class
 296   gSystem->SetIncludePath("-I. -I$ALICE_ROOT/include -I$ROOTSYS/include");
 297   gROOT->LoadMacro("./AliCFSingleTrackTask.cxx++g");
 298 }