PWGCF/Correlations/DPhi/AliAnalysisTaskDiHadronPID.cxx

   1 // ----------------------------------------------------------------------------\r
   2 // This class makes di-hadron correlations, with TOF and TPC signals for\r
   3 // the associated particles.\r
   4 //\r
   5 //   Last Update:\r
   6 //     - Added a number of setters/getters.\r
   7 //     - Added spectra in pt, eta and phi as a function of performed cuts.\r
   8 //     - Added DCA histogram.\r
   9 //     - Added pp functionality.\r
  10 //     - Added the option to make a DCA cut and an ITS cut.\r
  11 //     - Variable centrality.\r
  12 //     - Variable maximum p_T for triggers.\r
  13 //     - Removed di-hadron correlations with one PID signal.\r
  14 //\r
  15 // ----------------------------------------------------------------------------\r
  16 // Author: Misha Veldhoen (misha.veldhoen@cern.ch)\r
  17 // Last edit: Apr 12th 2012. (v 8.00)\r
  18 // ----------------------------------------------------------------------------\r
  19 \r
  20 #include <iostream>\r
  21 #include "TChain.h"\r
  22 #include "TTree.h"\r
  23 #include "TH1F.h"\r
  24 #include "TH2F.h"\r
  25 #include "TH3F.h"\r
  26 #include "TCanvas.h"\r
  27 #include "TFile.h"\r
  28 \r
  29 #include "AliAODTrack.h"\r
  30 #include "AliAODEvent.h"\r
  31 #include "AliAODInputHandler.h"\r
  32 #include "AliAODVertex.h"\r
  33 //#include "AliAODPid.h"\r
  34 \r
  35 #include "AliAnalysisManager.h"\r
  36 #include "AliInputEventHandler.h"\r
  37 #include "AliPIDResponse.h"\r
  38 #include "AliTPCPIDResponse.h"\r
  39 //#include "AliTOFPIDResponse.h"\r
  40 \r
  41 using namespace std;\r
  42 \r
  43 #include "AliAnalysisTaskDiHadronPID.h"\r
  44 \r
  45 ClassImp(AliAnalysisTaskDiHadronPID);\r
  46 \r
  47 //_____________________________________________________________________________\r
  48 AliAnalysisTaskDiHadronPID::AliAnalysisTaskDiHadronPID():\r
  49         AliAnalysisTaskSE(),\r
  50         fPIDResponse(0x0),\r
  51         fAODEvent(0x0),\r
  52         fAODHeader(0x0),\r
  53         fAODVertex(0x0),\r
  54         fAODTrack(0x0),\r
  55         fGlobalTracks(0x0),\r
  56         fCentrality(0x0),\r
  57         fVertexZ(0x0),\r
  58     fDCA(0x0),\r
  59     fDCAZoomed(0x0),\r
  60     fDCAZoomedTwice(0x0),\r
  61     fDCACut(0x0),\r
  62     fDCAZoomedCut(0x0),\r
  63     fDCAZoomedTwiceCut(0x0),\r
  64     fITSHits(0x0),\r
  65     fTrackCutsCount(0x0),\r
  66     fTrackCutsPt(0x0),\r
  67     fTrackCutsEta(0x0),\r
  68     fTrackCutsPhi(0x0),\r
  69     fEtaSpectrumTrig(0x0),\r
  70     fEtaSpectrumAssoc(0x0),\r
  71     fPhiSpectrumAssoc(0x0),\r
  72         fTPCnSigmaProton(0x0),\r
  73         fTOFnSigmaProton(0x0),\r
  74         fTPCnSigmaPion(0x0),\r
  75         fTOFnSigmaPion(0x0),\r
  76         fTPCnSigmaKaon(0x0),\r
  77         fTOFnSigmaKaon(0x0),\r
  78         fTPCSignal(0x0),\r
  79         fTOFSignal(0x0),\r
  80     fMixedEvents(0x0),\r
  81         fHistoList(0x0),\r
  82     fCalculateMixedEvents(kFALSE),\r
  83     fBeamType("PbPb"),\r
  84     fMaxEta(0.8),\r
  85     fMaxPlotEta(0.9),\r
  86     fMaxPt(10.),\r
  87     fNEtaBins(25),\r
  88     fNPhiBins(36),\r
  89     fVertexZMixedEvents(2.),\r
  90     fCentralityCutMax(0.),\r
  91     fCentralityCutMin(10.),\r
  92     fZoomed(kFALSE),\r
  93     fDoITSCut(kFALSE),\r
  94     fDoDCACut(kFALSE),\r
  95     fDemandNoMismatch(kFALSE),\r
  96     fVerbose(0),\r
  97     fPrintBufferSize(kFALSE),\r
  98         fTrigBufferIndex(0),    \r
  99         fTrigBufferSize(0),\r
 100         fTrigBufferMaxSize(1000)\r
 101 \r
 102 {\r
 103     \r
 104         //\r
 105         // Default Constructor.\r
 106         //\r
 107         \r
 108         // Trigger buffer.\r
 109         for(Int_t ii=0; ii<25000; ii++) {\r
 110                 for(Int_t jj=0; jj<4; jj++) {\r
 111                         fTrigBuffer[ii][jj]=0;\r
 112                 }                               \r
 113         }       \r
 114         \r
 115     // The identified di-hadron correlations.\r
 116         for (Int_t ii = 0; ii < 3; ii++) {\r
 117                 for (Int_t jj = 0; jj < 10; jj++) {\r
 118                         fDiHadronTPCTOF[ii][jj]=0x0;\r
 119                 }\r
 120         }\r
 121     \r
 122     // Track cut labels.\r
 123     for (Int_t ii=0; ii<8; ii++) {\r
 124         fTrackCutLabelNumbers[ii]=ii+1;\r
 125     }\r
 126 }\r
 127 \r
 128 //_____________________________________________________________________________\r
 129 AliAnalysisTaskDiHadronPID::AliAnalysisTaskDiHadronPID(const char *name):\r
 130         AliAnalysisTaskSE(name),\r
 131     fPIDResponse(0x0),\r
 132     fAODEvent(0x0),\r
 133     fAODHeader(0x0),\r
 134     fAODVertex(0x0),\r
 135     fAODTrack(0x0),\r
 136     fGlobalTracks(0x0),\r
 137     fCentrality(0x0),\r
 138     fVertexZ(0x0),\r
 139     fDCA(0x0),\r
 140     fDCAZoomed(0x0),\r
 141     fDCAZoomedTwice(0x0),\r
 142     fDCACut(0x0),\r
 143     fDCAZoomedCut(0x0),\r
 144     fDCAZoomedTwiceCut(0x0),\r
 145     fITSHits(0x0),\r
 146     fTrackCutsCount(0x0),\r
 147     fTrackCutsPt(0x0),\r
 148     fTrackCutsEta(0x0),\r
 149     fTrackCutsPhi(0x0),\r
 150     fEtaSpectrumTrig(0x0),\r
 151     fEtaSpectrumAssoc(0x0),\r
 152     fPhiSpectrumAssoc(0x0),\r
 153     fTPCnSigmaProton(0x0),\r
 154     fTOFnSigmaProton(0x0),\r
 155     fTPCnSigmaPion(0x0),\r
 156     fTOFnSigmaPion(0x0),\r
 157     fTPCnSigmaKaon(0x0),\r
 158     fTOFnSigmaKaon(0x0),\r
 159     fTPCSignal(0x0),\r
 160     fTOFSignal(0x0),\r
 161     fMixedEvents(0x0),\r
 162     fHistoList(0x0),\r
 163     fCalculateMixedEvents(kFALSE),\r
 164     fBeamType("PbPb"),\r
 165     fMaxEta(0.8),\r
 166     fMaxPlotEta(0.9),\r
 167     fMaxPt(10.),\r
 168     fNEtaBins(25),\r
 169     fNPhiBins(36),\r
 170     fVertexZMixedEvents(2.),\r
 171     fCentralityCutMax(0.),\r
 172     fCentralityCutMin(10.),\r
 173     fZoomed(kFALSE),\r
 174     fDoITSCut(kFALSE),\r
 175     fDoDCACut(kFALSE),\r
 176     fDemandNoMismatch(kFALSE),\r
 177     fVerbose(0),\r
 178     fPrintBufferSize(kFALSE),\r
 179     fTrigBufferIndex(0),        \r
 180     fTrigBufferSize(0),\r
 181         fTrigBufferMaxSize(1000)\r
 182 \r
 183 {\r
 184     \r
 185         //\r
 186         // Named Constructor.\r
 187         //\r
 188         \r
 189         DefineInput(0, TChain::Class());\r
 190         DefineOutput(1, TList::Class());\r
 191         \r
 192         // Trigger buffer.\r
 193         for(Int_t ii=0; ii<25000; ii++) {\r
 194                 for(Int_t jj=0; jj<4; jj++) {\r
 195                         fTrigBuffer[ii][jj]=0;\r
 196                 }                               \r
 197         }       \r
 198         \r
 199     // The identified di-hadron correlations.\r
 200         for (Int_t ii = 0; ii < 3; ii++) {\r
 201                 for (Int_t jj = 0; jj < 10; jj++) {\r
 202                         fDiHadronTPCTOF[ii][jj]=0x0;\r
 203                 }\r
 204         }\r
 205     \r
 206     // Track cut labels.\r
 207     for (Int_t ii=0; ii<8; ii++) {\r
 208         fTrackCutLabelNumbers[ii]=ii+1;\r
 209     }\r
 210 }\r
 211 \r
 212 //_____________________________________________________________________________\r
 213 AliAnalysisTaskDiHadronPID::~AliAnalysisTaskDiHadronPID() {\r
 214 \r
 215         //\r
 216         // Destructor.\r
 217         //\r
 218         \r
 219     if(fGlobalTracks) {\r
 220         delete fGlobalTracks;\r
 221         fGlobalTracks=0;\r
 222     }\r
 223 \r
 224 }\r
 225 \r
 226 //_____________________________________________________________________________\r
 227 void AliAnalysisTaskDiHadronPID::UserCreateOutputObjects() \r
 228 \r
 229 {\r
 230         //\r
 231         // Output objects.\r
 232         //\r
 233     \r
 234     // Print the settings of the analysis task.\r
 235     cout<<endl;\r
 236     cout<<"-----------------------------------------------"<<endl;\r
 237     cout<<" AliAnalysisTaskDiHadronPID Settings:"<<endl;\r
 238     cout<<"-----------------------------------------------"<<endl;\r
 239     cout<<"Verbose Level: "<<fVerbose<<endl;\r
 240     cout<<"Mixed Events Calculated: "<<fCalculateMixedEvents<<endl;\r
 241     cout<<"Beam Type: "<<fBeamType<<endl;\r
 242     cout<<Form("Max eta: %3.1f",fMaxEta)<<endl;\r
 243     cout<<Form("Max eta plotted: %3.1f",fMaxPlotEta)<<endl;\r
 244     cout<<Form("Max p_T for the triggers: %3.1f GeV/c.",fMaxPt)<<endl;\r
 245     cout<<"Nbins in eta and delta eta: "<<fNEtaBins<<endl;\r
 246     cout<<"Nbins in phi and delta phi: "<<fNPhiBins<<endl;\r
 247     cout<<Form("Events mixed if vertices differ %3.1f cm.",fVertexZMixedEvents)<<endl;\r
 248     cout<<Form("Centrality between %3.1f and %3.1f percent.",fCentralityCutMax,fCentralityCutMin)<<endl;\r
 249     cout<<"Tracks are cut if less than 2 SPD hits: "<<fDoITSCut<<endl;\r
 250     cout<<"Tracks cut if DCA is too big: "<<fDoDCACut<<endl;\r
 251     cout<<"Tracks cut if there is a TPC-TOF mismatch: "<<fDemandNoMismatch<<endl;\r
 252     cout<<"Maximum number of triggers stored: "<<fTrigBufferMaxSize<<endl;\r
 253     cout<<"-----------------------------------------------"<<endl;\r
 254     cout<<endl;\r
 255     \r
 256         // Obtain a pointer to the analysis manager.\r
 257         AliAnalysisManager* manager = AliAnalysisManager::GetAnalysisManager();\r
 258     if (!manager) {\r
 259         if (fVerbose>0) cout<<"AliAnalysisTaskDiHadronPID::UserCreateOutputObjects -> ERROR: Analysis manager not found."<<endl;\r
 260         return;\r
 261     }\r
 262     if (fVerbose>1) {\r
 263         cout<<"AliAnalysisTaskDiHadronPID::UserCreateOutputObjects -> Analysis manager found."<<endl;\r
 264     }\r
 265     \r
 266     // Obtain a pointer to the input handler.\r
 267     AliInputEventHandler* inputHandler = dynamic_cast<AliInputEventHandler*> (manager->GetInputEventHandler());\r
 268     if (!inputHandler) {\r
 269         if (fVerbose>0) cout<<"AliAnalysisTaskDiHadronPID::UserCreateOutputObjects -> ERROR: Input handler not found."<<endl;\r
 270         return;\r
 271     }\r
 272     if (fVerbose>1) {\r
 273         cout<<"AliAnalysisTaskDiHadronPID::UserCreateOutputObjects -> Input handler found."<<endl;\r
 274     }\r
 275     \r
 276         // Obtain a pointer to the PID response object. \r
 277         fPIDResponse = inputHandler->GetPIDResponse();\r
 278     if (!fPIDResponse) {\r
 279         if (fVerbose>0) cout<<"AliAnalysisTaskDiHadronPID::UserCreateOutputObjects -> ERROR: PID response object not found."<<endl;\r
 280         return;\r
 281     }\r
 282     if (fVerbose>1) {\r
 283         cout<<"AliAnalysisTaskDiHadronPID::UserCreateOutputObjects -> PID response object found."<<endl;\r
 284     }\r
 285 \r
 286         // Create the output of the task.       \r
 287         fHistoList = new TList();\r
 288         fHistoList->SetOwner(kTRUE); \r
 289         \r
 290         // Ranges in dPhi, dEta, and the number of bins for di-hadron correlations.\r
 291         Int_t binsHisto[4] = {fNPhiBins,fNEtaBins};\r
 292         Double_t minHisto[4] = {-TMath::Pi()/2.,-2*fMaxPlotEta};\r
 293         Double_t maxHisto[4] = {3.*TMath::Pi()/2,2*fMaxPlotEta};\r
 294         \r
 295     // --- EVENT QA PLOTS ---\r
 296     \r
 297         fCentrality = new TH1F("fCentrality","Centrality;Centrality;N_{evt}",100,0,100);\r
 298         fHistoList->Add(fCentrality);\r
 299         fVertexZ = new TH1F("fVertexZ","Vertex Z position;z (cm);N_{evt}",30,-15,15);\r
 300         fHistoList->Add(fVertexZ);\r
 301         \r
 302     // --- TRACK QA PLOTS ---\r
 303     \r
 304     fDCA = new TH2F("fDCA","DCA positions TPC only cuts;xy (cm);z (cm)",100,-5,5,100,-5,5);\r
 305     fHistoList->Add(fDCA);\r
 306     fDCAZoomed = new TH2F("fDCAZoomed","DCA positions TPC only cuts;xy (cm);z (cm)",100,-.5,.5,100,-.5,.5);\r
 307     fHistoList->Add(fDCAZoomed);\r
 308     fDCAZoomedTwice = new TH2F("fDCAZoomedTwice","DCA positions TPC only cuts;xy (cm);z (cm)",100,-.05,.05,100,-.05,.05);\r
 309     fHistoList->Add(fDCAZoomedTwice);\r
 310     fDCACut = new TH2F("fDCACut","DCA positions after DCA cut;xy (cm);z (cm)",100,-5,5,100,-5,5);\r
 311     fHistoList->Add(fDCACut);\r
 312     fDCAZoomedCut = new TH2F("fDCAZoomedCut","DCA positions after DCA cut;xy (cm);z (cm)",100,-.5,.5,100,-.5,.5);\r
 313     fHistoList->Add(fDCAZoomedCut);\r
 314     fDCAZoomedTwiceCut = new TH2F("fDCAZoomedTwiceCut","DCA positions after DCA cut;xy (cm);z (cm)",100,-.05,.05,100,-.05,.05);\r
 315     fHistoList->Add(fDCAZoomedTwiceCut);\r
 316     \r
 317     fITSHits = new TH1F("fITSHits","ITS hits",3,0,3);\r
 318     (fITSHits->GetXaxis())->SetBinLabel(1,"No SPD hit");\r
 319     (fITSHits->GetXaxis())->SetBinLabel(2,"One SPD hit");\r
 320     (fITSHits->GetXaxis())->SetBinLabel(3,"Two SPD hits");\r
 321     fHistoList->Add(fITSHits);\r
 322     \r
 323     // Determine the bins used in the track cut histograms.\r
 324     const Int_t ncuts = 5 + fDoITSCut + fDoDCACut + fDemandNoMismatch;\r
 325     \r
 326     if (!fDoITSCut) {\r
 327         fTrackCutLabelNumbers[3]--;\r
 328         fTrackCutLabelNumbers[4]--;\r
 329         fTrackCutLabelNumbers[5]--;\r
 330         fTrackCutLabelNumbers[6]--;\r
 331         fTrackCutLabelNumbers[7]--;\r
 332     }\r
 333     \r
 334     if (!fDoDCACut) {\r
 335         fTrackCutLabelNumbers[4]--;\r
 336         fTrackCutLabelNumbers[5]--;\r
 337         fTrackCutLabelNumbers[6]--;\r
 338         fTrackCutLabelNumbers[7]--;\r
 339     }\r
 340     \r
 341     fTrackCutsCount = new TH1F("fTrackCutsCount","Track Cuts;;Count",ncuts,0,ncuts);\r
 342     (fTrackCutsCount->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[0],"Std. TPC Only");\r
 343     (fTrackCutsCount->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[1],Form("#eta < %3.1f",fMaxEta));\r
 344     if (fDoITSCut) (fTrackCutsCount->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[2],"1-2 SPD hits");\r
 345     if (fDoDCACut) (fTrackCutsCount->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[3],"DCA cut");\r
 346     (fTrackCutsCount->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[4],"p_{T} < 5.0 GeV/c");\r
 347     (fTrackCutsCount->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[5],"TPCpid (p_{T} < 5.0 GeV/c)");\r
 348     (fTrackCutsCount->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[6],"TOFpid (p_{T} < 5.0 GeV/c)");\r
 349     if (fDemandNoMismatch) (fTrackCutsCount->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[7],"No Mismatch (p_{T} < 5.0 GeV/c)");\r
 350     fHistoList->Add(fTrackCutsCount);\r
 351     \r
 352     fTrackCutsPt = new TH2F("fTrackCutsPt","Track Cuts vs p_{T};;p_{T};Count",ncuts,0,ncuts,40,0,fMaxPt);\r
 353     (fTrackCutsPt->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[0],"Std. TPC Only");\r
 354     (fTrackCutsPt->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[1],Form("#eta < %3.1f",fMaxEta));\r
 355     if (fDoITSCut) (fTrackCutsPt->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[2],"1-2 SPD hits");\r
 356     if (fDoDCACut) (fTrackCutsPt->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[3],"DCA cut");\r
 357     (fTrackCutsPt->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[4],"p_{T} < 5.0 GeV/c");\r
 358     (fTrackCutsPt->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[5],"TPCpid (p_{T} < 5.0 GeV/c)");\r
 359     (fTrackCutsPt->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[6],"TOFpid (p_{T} < 5.0 GeV/c)");\r
 360     if (fDemandNoMismatch) (fTrackCutsPt->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[7],"No Mismatch (p_{T} < 5.0 GeV/c)");\r
 361     fHistoList->Add(fTrackCutsPt);\r
 362     \r
 363     fTrackCutsEta = new TH2F("fTrackCutsEta","Track Cuts vs #eta;;#eta;Count",ncuts,0,ncuts,4*binsHisto[1],-fMaxPlotEta,fMaxPlotEta);\r
 364     (fTrackCutsEta->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[0],"Std. TPC Only");\r
 365     (fTrackCutsEta->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[1],Form("#eta < %3.1f",fMaxEta));\r
 366     if (fDoITSCut) (fTrackCutsEta->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[2],"1-2 SPD hits");\r
 367     if (fDoDCACut) (fTrackCutsEta->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[3],"DCA cut");\r
 368     (fTrackCutsEta->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[4],"p_{T} < 5.0 GeV/c");\r
 369     (fTrackCutsEta->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[5],"TPCpid (p_{T} < 5.0 GeV/c)");\r
 370     (fTrackCutsEta->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[6],"TOFpid (p_{T} < 5.0 GeV/c)");\r
 371     if (fDemandNoMismatch) (fTrackCutsEta->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[7],"No Mismatch (p_{T} < 5.0 GeV/c)");\r
 372     fHistoList->Add(fTrackCutsEta);\r
 373     \r
 374     fTrackCutsPhi = new TH2F("fTrackCutsPhi","Track Cuts vs #phi;;#phi;Count",ncuts,0,ncuts,72,0,2.*TMath::Pi());\r
 375     (fTrackCutsPhi->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[0],"Std. TPC Only");\r
 376     (fTrackCutsPhi->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[1],Form("#eta < %3.1f",fMaxEta));\r
 377     if (fDoITSCut) (fTrackCutsPhi->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[2],"1-2 SPD hits");\r
 378     if (fDoDCACut) (fTrackCutsPhi->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[3],"DCA cut");\r
 379     (fTrackCutsPhi->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[4],"p_{T} < 5.0 GeV/c");\r
 380     (fTrackCutsPhi->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[5],"TPCpid (p_{T} < 5.0 GeV/c)");\r
 381     (fTrackCutsPhi->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[6],"TOFpid (p_{T} < 5.0 GeV/c)");\r
 382     if (fDemandNoMismatch) (fTrackCutsPhi->GetXaxis())->SetBinLabel(fTrackCutLabelNumbers[7],"No Mismatch (p_{T} < 5.0 GeV/c)");\r
 383     fHistoList->Add(fTrackCutsPhi);\r
 384     \r
 385     fEtaSpectrumTrig = new TH1F("fEtaSpectrumTrig","#eta Spectrum Triggers;#eta;Count",4*binsHisto[1],-fMaxPlotEta,fMaxPlotEta);\r
 386     fHistoList->Add(fEtaSpectrumTrig);\r
 387     fEtaSpectrumAssoc = new TH2F("fEtaSpectrumAssoc","#eta Spectrum Associateds;#eta;p_{T};Count",4*binsHisto[1],-fMaxPlotEta,fMaxPlotEta,10,0.,5.);\r
 388     fHistoList->Add(fEtaSpectrumAssoc); \r
 389     fPhiSpectrumAssoc = new TH2F("fPhiSpectrumAssoc","#phi Spectrum Associateds;#phi;p_{T};Count",72,0,2.*TMath::Pi(),10,0.,5.);\r
 390     fHistoList->Add(fPhiSpectrumAssoc);\r
 391         \r
 392     // --- QA PLOTS PID ---\r
 393     \r
 394         fTPCnSigmaProton = new TH2F("fTPCnSigmaProton","n#sigma plot for the TPC (Protons);p (GeV/c);n#sigma",100,0.,5.,100,-10.,10.);\r
 395         fHistoList->Add(fTPCnSigmaProton);\r
 396         fTOFnSigmaProton = new TH2F("fTOFnSigmaProton","n#sigma plot for the TOF (Protons);p (GeV/c);n#sigma",100,0.,5.,100,-10.,10.);\r
 397         fHistoList->Add(fTOFnSigmaProton);\r
 398         fTPCnSigmaPion = new TH2F("fTPCnSigmaPion","n#sigma plot for the TPC (Pions);p (GeV/c);n#sigma",100,0.,5.,100,-10.,10.);\r
 399         fHistoList->Add(fTPCnSigmaPion);\r
 400         fTOFnSigmaPion = new TH2F("fTOFnSigmaPion","n#sigma plot for the TOF (Pions);p (GeV/c);n#sigma",100,0.,5.,100,-10.,10.);\r
 401         fHistoList->Add(fTOFnSigmaPion);\r
 402         fTPCnSigmaKaon = new TH2F("fTPCnSigmaKaon","n#sigma plot for the TPC (Kaons);p (GeV/c);n#sigma",100,0.,5.,100,-10.,10.);\r
 403         fHistoList->Add(fTPCnSigmaKaon);\r
 404         fTOFnSigmaKaon = new TH2F("fTOFnSigmaKaon","n#sigma plot for the TOF (Kaons);p (GeV/c);n#sigma",100,0.,5.,100,-10.,10.);\r
 405         fHistoList->Add(fTOFnSigmaKaon);        \r
 406         fTPCSignal = new TH3F("fTPCSignal","TPC Signal;#eta;p_{T} GeV/c;dE/dx",25,-fMaxEta,fMaxEta,100,0.,5.,150,0.,300.);\r
 407         fHistoList->Add(fTPCSignal);\r
 408         fTOFSignal = new TH3F("fTOFSignal","TOF Signal;#eta;p_{T} GeV/c;t",25,-fMaxEta,fMaxEta,100,0.,5.,150,10000.,25000.);\r
 409         fHistoList->Add(fTOFSignal);\r
 410                         \r
 411     // --- DI-HADRON CORRELATIONS WITH TPC AND TOF SIGNALS ---\r
 412     \r
 413         // Di Hadron Correlations with two PID signals, for each p_T bin and particle species separately. (i.e. 30 histo's) \r
 414         // Axes: {Dphi, Deta, TPC signal, TOF signal}\r
 415     \r
 416     // Unzoomed pictures.\r
 417         Int_t binsTPCnoZoom[3][10] = {{100,100,100,100,100,100,100,100,100,100},\r
 418                                   {100,100,100,100,100,100,100,100,100,100},\r
 419                                   {100,100,100,100,100,100,100,100,100,100}};\r
 420                 \r
 421         Double_t minTPCnoZoom[3][10] = {{   -50., -50.,-30.,-30.,-30.,-30.,-30.,-30.,-30.,-30.},\r
 422                                     {   -100., -50.,-20.,-20.,-20.,-20.,-20.,-20.,-20.,-20.},\r
 423                                     {  -200.,-150.,-50.,-30.,-20.,-20.,-20.,-20.,-20.,-20.}};\r
 424         Double_t maxTPCnoZoom[3][10] = {{ 150., 100., 50., 30., 25., 25., 25., 25., 25., 25.},\r
 425                                     {  50., 100., 40., 30., 30., 30., 30., 30., 30., 30.},\r
 426                                     { 100.,  50., 50., 30., 30., 30., 30., 30., 30., 30.}};\r
 427         \r
 428         Int_t binsTOFnoZoom[3][10] = {{100,100,100,100,100,100,100,100,100,100},\r
 429                                   {100,100,100,100,100,100,100,100,100,100},\r
 430                                   {100,100,100,100,100,100,100,100,100,100}};\r
 431         \r
 432         Double_t minTOFnoZoom[3][10] = {{-2000.,-2000.,-1000., -500., -500., -500., -500.,      -500., -500., -500.},\r
 433                                     {-2000.,-2000.,-2000.,-1000.,-1000.,-1000.,-1000.,-1000.,-1000.,-1000.},\r
 434                                     {-2000.,-2000.,-6000.,-3000.,-2000.,-1500.,-1500.,-1500.,-1500.,-1500.}};\r
 435         Double_t maxTOFnoZoom[3][10] = {{ 2000., 2000., 5000., 3000., 2000., 1500., 1500., 1500., 1500., 1500.},\r
 436                                     { 2000., 2000., 5000., 2000., 1500., 1500., 1500., 1500., 1500., 1500.},\r
 437                                     { 2000., 2000., 2000., 1000., 1000., 1000., 1000., 1000., 1000., 1000.}};\r
 438     \r
 439     // Zoomed pictures.\r
 440     Int_t binsTPCZoom[3][10] = {{100,100,100,100,100,100,100,100,100,100},\r
 441                                 {100,100,100,100,100,100,100,100,100,100},\r
 442                                 {100,100,100,100,100,100,100,100,100,100}};\r
 443     \r
 444         Double_t minTPCZoom[3][10] = {{   -20., -20.,-20.,-20.,-20.,-20.,-20.,-20.,-20.,-20.},\r
 445                                   {      -20., -20.,-20.,-20.,-20.,-20.,-20.,-20.,-20.,-20.},\r
 446                                   {   -40., -20.,-20.,-20.,-20.,-20.,-20.,-20.,-20.,-20.}};\r
 447         Double_t maxTPCZoom[3][10] = {{  20.,  20., 20., 20., 20., 20., 20., 20., 20., 20.},\r
 448                                   {  20.,  20., 20., 20., 20., 20., 20., 20., 20., 20.},\r
 449                                   {  40.,  20., 20., 20., 20., 30., 30., 30., 30., 30.}};\r
 450         \r
 451         Int_t binsTOFZoom[3][10] = {{100,100,100,100,100,100,100,100,100,100},\r
 452                                 {100,100,100,100,100,100,100,100,100,100},\r
 453                                 {100,100,100,100,100,100,100,100,100,100}};\r
 454         \r
 455         Double_t minTOFZoom[3][10] = {{-1000.,-1000., -500., -500., -500., -400., -400., -400., -400., -400.},\r
 456                                   { -800., -800., -800., -800., -800., -600., -500., -500., -400., -400.},\r
 457                                   {-1000.,-1000.,-1000.,-1000., -800.,-1000.,-1000., -800., -700., -700.}};\r
 458         Double_t maxTOFZoom[3][10] = {{ 1000., 1000., 1000., 1000., 1000., 1000., 1000.,  900.,  800.,  700.},\r
 459                                   { 1000., 1000.,  500.,  500.,  500.,  900.,  700.,  700.,  600.,  500.},\r
 460                                   { 1000., 1000., 1000.,  500.,  500.,  500.,  500.,  500.,  500.,  500.}};\r
 461     \r
 462         TString basenameTPCTOF("fDiHadronTPCTOF");\r
 463         TString basetitle("Di-Hadron correlation");\r
 464         TString finalname, finaltitle;\r
 465         TString species[3] = {"Pion","Kaon","Proton"};\r
 466         TString ptbins[11] = {"0.0","0.5","1.0","1.5","2.0","2.5","3.0","3.5","4.0","4.5","5.0"};\r
 467     \r
 468         // Recall that AliPID::kPion = 2, AliPID::kKaon = 3, AliPID::kProton = 4.\r
 469         for (Int_t iSpecies = 0; iSpecies < 3; iSpecies++) {\r
 470                                 \r
 471                 for (Int_t iPtBin = 0; iPtBin < 10; iPtBin++) {\r
 472                         \r
 473             if (fZoomed) {\r
 474                 binsHisto[2] = binsTPCZoom[iSpecies][iPtBin];\r
 475                 binsHisto[3] = binsTOFZoom[iSpecies][iPtBin];\r
 476                 minHisto[2] = minTPCZoom[iSpecies][iPtBin];\r
 477                 minHisto[3] = minTOFZoom[iSpecies][iPtBin];\r
 478                 maxHisto[2] = maxTPCZoom[iSpecies][iPtBin];\r
 479                 maxHisto[3] = maxTOFZoom[iSpecies][iPtBin];\r
 480             } else {\r
 481                 binsHisto[2] = binsTPCnoZoom[iSpecies][iPtBin];\r
 482                 binsHisto[3] = binsTOFnoZoom[iSpecies][iPtBin];\r
 483                 minHisto[2] = minTPCnoZoom[iSpecies][iPtBin];\r
 484                 minHisto[3] = minTOFnoZoom[iSpecies][iPtBin];\r
 485                 maxHisto[2] = maxTPCnoZoom[iSpecies][iPtBin];\r
 486                 maxHisto[3] = maxTOFnoZoom[iSpecies][iPtBin];\r
 487             }\r
 488 \r
 489             // Make the di-hadron correlations with different pid signals.\r
 490             // TODO: Rewrite this with Form().\r
 491             finaltitle = basetitle;\r
 492                         (((((finaltitle += " (") += species[iSpecies]) += ") ") += ptbins[iPtBin]) += " < P_t < ") += ptbins[iPtBin+1];\r
 493             finaltitle+=";#Delta#phi;#Delta#eta;dE/dx;t (ms)";\r
 494             finalname = basenameTPCTOF;\r
 495                         (((finalname += "_") += species[iSpecies]) += "_") += iPtBin;\r
 496             \r
 497             fDiHadronTPCTOF[iSpecies][iPtBin] = new THnSparseF(finalname,finaltitle,4,binsHisto,minHisto,maxHisto);\r
 498             fHistoList->Add(fDiHadronTPCTOF[iSpecies][iPtBin]);\r
 499             \r
 500                 }\r
 501         }\r
 502     \r
 503     // --- MIXED EVENTS ---\r
 504     \r
 505     if (fCalculateMixedEvents) {\r
 506         fMixedEvents = new TH3F("fMixedEvents","Mixed Events;#Delta#phi;#Delta#eta;p_{T_assoc}",binsHisto[0],minHisto[0],maxHisto[0],binsHisto[1],minHisto[1],maxHisto[1],10,0.,5.);\r
 507         fMixedEvents->Sumw2();\r
 508         fHistoList->Add(fMixedEvents);\r
 509     }\r
 510     \r
 511         PostData(1, fHistoList);\r
 512         \r
 513 }\r
 514 \r
 515 //_____________________________________________________________________________\r
 516 Int_t AliAnalysisTaskDiHadronPID::ClassifyTrack(AliAODTrack *track)\r
 517 \r
 518 {\r
 519         //\r
 520     // This function both classifies tracks, and fills the track QA histograms.\r
 521     //\r
 522     // Classifies the track in:\r
 523     //  0 -> Not Useful,\r
 524     //  1 -> Associated,\r
 525     //  2 -> Trigger.\r
 526     //\r
 527     // IDEA: later we can do this with filterbits.\r
 528     //\r
 529     // The following track cuts are applied for triggers:\r
 530     //\r
 531     // 1a) pT > 5.0 GeV/c, pT < fPtMax,\r
 532     //  2) StandardTPCOnlyTrackCuts (filterbit 7 for AOD's),\r
 533     //  3) eta < fMaxEta,\r
 534     //  4) ITS track cut, a track is only selected if it has at least one hit in the SPD,\r
 535     //     that is, in one of the first two layers of the ITS. (can be switched on and off\r
 536     //     using (Bool_t)fDoITSCut),\r
 537     //  5) DCA track cut, of all tracks with at least one SPD hit, the DCA is constrained as\r
 538     //     a function of pT, in order to remove lots of the secondaries. (can be switched on and \r
 539     //     off by using (Bool_t)fDoDCACut),\r
 540     //\r
 541     // For associateds all the same track cuts are applied, except 1a is\r
 542     // replaced by 2b. Moreover the following cuts are applied too:\r
 543     //\r
 544     // 1b) pT < 5.0 GeV/c,\r
 545     //  6) TPCpid,\r
 546     //  7) TOFpid,\r
 547     //  8) no TPCTOF mismatch.\r
 548     //\r
 549     \r
 550     // Check if a track is supplied.\r
 551     if (!track) {\r
 552         if (fVerbose>0) cout<<"AliAnalysisTaskDiHadronPID::SelectTrack -> ERROR: No track found."<<endl;\r
 553         return 0;\r
 554     }\r
 555     \r
 556     // Get basic track information.\r
 557     Int_t classification=0;\r
 558     Double_t pt = track->Pt();\r
 559     Double_t eta = track->Eta();\r
 560     Double_t phi = track->Phi();\r
 561         \r
 562     // 1) pT cut: First separation between triggers and associateds.\r
 563     if (pt<5.0) classification=1;\r
 564     if ((pt>5.0)&&(pt<fMaxPt)) classification = 2;\r
 565     if (!classification) return 0;\r
 566     \r
 567     // 2) StandardTPCOnlyTrackCuts.\r
 568     if (!(track->TestFilterMask(1<<7))) {\r
 569                 //if (fVerbose>3) cout<<"Track Ignored: Did Not pass filterbit."<<endl;\r
 570                 return 0;\r
 571         }\r
 572        \r
 573         if (fVerbose>3) {\r
 574                 cout<<endl;\r
 575                 cout<<"pt: "<<pt<<" eta: "<<eta<<" phi: "<<phi<<endl;\r
 576         }\r
 577         \r
 578     fTrackCutsCount->Fill(fTrackCutLabelNumbers[0]-.5);\r
 579     fTrackCutsPt->Fill(fTrackCutLabelNumbers[0]-.5,pt);\r
 580     fTrackCutsEta->Fill(fTrackCutLabelNumbers[0]-.5,eta);\r
 581     fTrackCutsPhi->Fill(fTrackCutLabelNumbers[0]-.5,phi);\r
 582     \r
 583     // 3) eta cut.\r
 584     if (TMath::Abs(eta)>fMaxEta) {\r
 585                 if (fVerbose>3) cout<<"Track Ignored: Eta too large."<<endl;\r
 586                 return 0;\r
 587     }\r
 588         \r
 589     fTrackCutsCount->Fill(fTrackCutLabelNumbers[1]-.5);\r
 590     fTrackCutsPt->Fill(fTrackCutLabelNumbers[1]-.5,pt);\r
 591     fTrackCutsEta->Fill(fTrackCutLabelNumbers[1]-.5,eta);\r
 592     fTrackCutsPhi->Fill(fTrackCutLabelNumbers[1]-.5,phi);\r
 593     \r
 594     // Obtaining ITS information.\r
 595     AliAODTrack* globaltrack = GetGlobalTrack(track);\r
 596     Bool_t ITSLayerHit[6];\r
 597     for (Int_t iITSLayer=0; iITSLayer<6; iITSLayer++) {\r
 598         ITSLayerHit[iITSLayer] = globaltrack->HasPointOnITSLayer(iITSLayer);\r
 599     }\r
 600     Int_t SPDHits=ITSLayerHit[0]+ITSLayerHit[1];\r
 601     \r
 602         if (fVerbose>3) cout<<"SPD hits: "<<SPDHits<<endl;\r
 603         \r
 604     // Fill the ITS hist.\r
 605     fITSHits->Fill(SPDHits+0.5);\r
 606 \r
 607     // 4) ITS cut.\r
 608     if (fDoITSCut) {\r
 609         if (!SPDHits) {\r
 610                         if (fVerbose>3) cout<<"Track Ignored: Not enough SPD hits."<<endl;\r
 611                         return 0;\r
 612                 }\r
 613         fTrackCutsCount->Fill(fTrackCutLabelNumbers[2]-.5);\r
 614         fTrackCutsPt->Fill(fTrackCutLabelNumbers[2]-.5,pt);\r
 615         fTrackCutsEta->Fill(fTrackCutLabelNumbers[2]-.5,eta);\r
 616         fTrackCutsPhi->Fill(fTrackCutLabelNumbers[2]-.5,phi);\r
 617     }\r
 618         \r
 619     // Propagate the global track to the DCA.\r
 620     Double_t PosAtDCA[2] = {-999,-999};\r
 621     Double_t covar[3] = {-999,-999,-999};\r
 622     globaltrack->PropagateToDCA(fAODVertex,fAODEvent->GetMagneticField(),100.,PosAtDCA,covar);\r
 623         \r
 624     // Fill the DCA hist (before cut)\r
 625     fDCA->Fill(PosAtDCA[0],PosAtDCA[1]);\r
 626     fDCAZoomed->Fill(PosAtDCA[0],PosAtDCA[1]);\r
 627     fDCAZoomedTwice->Fill(PosAtDCA[0],PosAtDCA[1]);\r
 628     \r
 629     // 5) DCA cut (See R_AA paper).\r
 630     Double_t DCAcutvalue[2];\r
 631     DCAcutvalue[0] = 0.018 + 0.035*TMath::Power(pt,-1.01);\r
 632     DCAcutvalue[1] = 2.; \r
 633     if (fDoDCACut) {\r
 634         if (SPDHits&&((TMath::Abs(PosAtDCA[0])>DCAcutvalue[0])||(TMath::Abs(PosAtDCA[1])>DCAcutvalue[1]))) {\r
 635                         if (fVerbose>3) cout<<"Track Ignored: Enough SPD hits, but out of DCA range."<<endl;\r
 636                         return 0;\r
 637                 }\r
 638         fTrackCutsCount->Fill(fTrackCutLabelNumbers[3]-.5);\r
 639         fTrackCutsPt->Fill(fTrackCutLabelNumbers[3]-.5,pt);\r
 640         fTrackCutsEta->Fill(fTrackCutLabelNumbers[3]-.5,eta);\r
 641         fTrackCutsPhi->Fill(fTrackCutLabelNumbers[3]-.5,phi);\r
 642         \r
 643         // Fill the DCA hist (after cut)\r
 644         fDCACut->Fill(PosAtDCA[0],PosAtDCA[1]);\r
 645         fDCAZoomedCut->Fill(PosAtDCA[0],PosAtDCA[1]);\r
 646         fDCAZoomedTwiceCut->Fill(PosAtDCA[0],PosAtDCA[1]);\r
 647     }\r
 648     \r
 649     // Now all the common cuts have been performed. Tracks identified as a trigger will\r
 650     // be returned. Note that they will not appear in the histograms after track cut 5.\r
 651     if (classification==2) return 2;\r
 652     \r
 653     // Now we're left with only tracks with pT < 5.0 GeV/c. \r
 654     fTrackCutsCount->Fill(fTrackCutLabelNumbers[4]-.5);\r
 655     fTrackCutsPt->Fill(fTrackCutLabelNumbers[4]-.5,pt);\r
 656     fTrackCutsEta->Fill(fTrackCutLabelNumbers[4]-.5,eta);\r
 657     fTrackCutsPhi->Fill(fTrackCutLabelNumbers[4]-.5,phi);\r
 658     \r
 659     // Obtain the status of the track.\r
 660     ULong_t status = globaltrack->GetStatus();\r
 661     \r
 662     // 6) TPCpid\r
 663     if (!((status&AliAODTrack::kTPCpid)==AliAODTrack::kTPCpid)) {\r
 664                 if (fVerbose>3) cout<<"Track Ignored: No TPC pid."<<endl;\r
 665                 return 0;\r
 666         }\r
 667     \r
 668     fTrackCutsCount->Fill(fTrackCutLabelNumbers[5]-.5);\r
 669     fTrackCutsPt->Fill(fTrackCutLabelNumbers[5]-.5,pt);\r
 670     fTrackCutsEta->Fill(fTrackCutLabelNumbers[5]-.5,eta);\r
 671     fTrackCutsPhi->Fill(fTrackCutLabelNumbers[5]-.5,phi);\r
 672     \r
 673     // 7) TOFpid\r
 674     if (!((status&AliAODTrack::kTOFpid)==AliAODTrack::kTOFpid)) {\r
 675                 if (fVerbose>3) cout<<"Track Ignored: No TOF pid."<<endl;\r
 676                 return 0;\r
 677         }\r
 678     \r
 679     fTrackCutsCount->Fill(fTrackCutLabelNumbers[6]-.5);\r
 680     fTrackCutsPt->Fill(fTrackCutLabelNumbers[6]-.5,pt);\r
 681     fTrackCutsEta->Fill(fTrackCutLabelNumbers[6]-.5,eta);\r
 682     fTrackCutsPhi->Fill(fTrackCutLabelNumbers[6]-.5,phi);\r
 683     \r
 684     // 8) TPC, TOF mismatch.\r
 685     if (fDemandNoMismatch) {\r
 686         if ((status&AliAODTrack::kTOFmismatch)==AliAODTrack::kTOFmismatch) {\r
 687                         if (fVerbose>3) cout<<"Track Ignored: TOF mismatch found."<<endl;\r
 688                         return 0;\r
 689                 }\r
 690         fTrackCutsCount->Fill(fTrackCutLabelNumbers[7]-.5);\r
 691         fTrackCutsPt->Fill(fTrackCutLabelNumbers[7]-.5,pt);\r
 692         fTrackCutsEta->Fill(fTrackCutLabelNumbers[7]-.5,eta);\r
 693         fTrackCutsPhi->Fill(fTrackCutLabelNumbers[7]-.5,phi);\r
 694     }\r
 695     \r
 696     // All tracks which made it up to here are classified as associateds.\r
 697     \r
 698     // Fill the PID QA histograms.\r
 699     Double_t mom, nSigma;\r
 700     \r
 701     mom = globaltrack->GetTPCmomentum();\r
 702     nSigma = fPIDResponse->NumberOfSigmasTPC(globaltrack,AliPID::kProton);\r
 703     fTPCnSigmaProton->Fill(mom,nSigma);\r
 704     nSigma = fPIDResponse->NumberOfSigmasTPC(globaltrack,AliPID::kPion);\r
 705     fTPCnSigmaPion->Fill(mom,nSigma);\r
 706     nSigma = fPIDResponse->NumberOfSigmasTPC(globaltrack,AliPID::kKaon);\r
 707     fTPCnSigmaKaon->Fill(mom,nSigma);\r
 708     \r
 709     fTPCSignal->Fill(eta,pt,globaltrack->GetTPCsignal());\r
 710     \r
 711     mom =globaltrack->P();\r
 712     nSigma = fPIDResponse->NumberOfSigmasTOF(globaltrack,AliPID::kProton);\r
 713     fTOFnSigmaProton->Fill(mom,nSigma);                 \r
 714     nSigma = fPIDResponse->NumberOfSigmasTOF(globaltrack,AliPID::kPion);\r
 715     fTOFnSigmaPion->Fill(mom,nSigma);   \r
 716     nSigma = fPIDResponse->NumberOfSigmasTOF(globaltrack,AliPID::kKaon);\r
 717     fTOFnSigmaKaon->Fill(mom,nSigma);\r
 718     \r
 719     fTOFSignal->Fill(eta,pt,globaltrack->GetTOFsignal());\r
 720     \r
 721     // Return associated.\r
 722     return 1;\r
 723     \r
 724 } \r
 725     \r
 726 //____________________________________________________________________________\r
 727 Bool_t AliAnalysisTaskDiHadronPID::SelectEvent(AliAODVertex* vertex)\r
 728 \r
 729 {\r
 730         \r
 731         //\r
 732         // Event Selection.\r
 733         //\r
 734         \r
 735         Double_t primVtx[3];\r
 736         vertex->GetXYZ(primVtx);\r
 737         if (TMath::Sqrt(primVtx[0]*primVtx[0] + primVtx[1]*primVtx[1])>1. || TMath::Abs(primVtx[2])>10.) {\r
 738                 if (fVerbose>2) cout << "AliAnalysisTaskDiHadronPID::SelectEvent -> Vertex Out of Range." << endl;\r
 739         if (fVerbose>2) cout << "AliAnalysisTaskDiHadronPID::SelectEvent -> Event not selected." << endl;\r
 740                 return kFALSE;\r
 741         }\r
 742         if (fVerbose>2) cout << "AliAnalysisTaskDiHadronPID::SelectEvent -> Vertex is OK." << endl;\r
 743     \r
 744     // We also wish to make centrality cut, but only if run on Pb+Pb.\r
 745     if (fBeamType=="PbPb") {\r
 746         Double_t cent = fAODHeader->GetCentrality();\r
 747         if (cent>fCentralityCutMin||cent<fCentralityCutMax) {\r
 748             if (fVerbose>2) cout<<"AliAnalysisTaskDiHadronPID::SelectEvent -> Event did not pass centaltity cut."<<endl;\r
 749             return kFALSE;\r
 750         }\r
 751     }\r
 752 \r
 753     if (fVerbose>2) cout << "AliAnalysisTaskDiHadronPID::SelectEvent -> Event selected." << endl;\r
 754         return kTRUE;\r
 755         \r
 756 }\r
 757 \r
 758 //_____________________________________________________________________________\r
 759 void AliAnalysisTaskDiHadronPID::FillGlobalTracksArray() {\r
 760 \r
 761         // Initialize the mapping for corresponding PID tracks. (see\r
 762         // GetGlobalTrack(AliAODTrack* track)). \r
 763         //\r
 764         \r
 765         if (!fAODEvent) {\r
 766         if (fVerbose>0) cout << "AliAnalysisTaskDiHadronPID::FillGlobalTracksArray -> ERROR: fAODEvent not set." << endl;\r
 767                 return;\r
 768         }\r
 769         \r
 770         fGlobalTracks = new TObjArray();\r
 771         AliAODTrack* track = 0x0;\r
 772                 \r
 773         for (Int_t iTrack = 0; iTrack < fAODEvent->GetNumberOfTracks(); iTrack++) {\r
 774                 \r
 775                 track = fAODEvent->GetTrack(iTrack);\r
 776                 \r
 777                 // I.e., if it does NOT pass the filtermask.\r
 778                 if (!(track->TestFilterMask(1<<7))) {\r
 779             if (track->GetID()>-1) fGlobalTracks->AddAtAndExpand(track,track->GetID());\r
 780             //if (track->GetID()<1) cout<<"Track ID: "<<track->GetID()<<" Partner ID: "<<(-track->GetID()-1)<<endl;\r
 781                 }\r
 782         \r
 783         }\r
 784         \r
 785 }\r
 786 \r
 787 //_____________________________________________________________________________\r
 788 AliAODTrack* AliAnalysisTaskDiHadronPID::GetGlobalTrack(AliAODTrack* track) {\r
 789         \r
 790         //\r
 791         // Returns the "parner track" of track, which contains the pid information.\r
 792         //\r
 793         \r
 794         AliAODTrack* partner = 0x0;\r
 795             \r
 796     partner = (AliAODTrack*)(fGlobalTracks->At(-track->GetID()-1));\r
 797             \r
 798         if (!partner&&(fVerbose>3)) cout<<"AliAnalysisTaskDiHadronPID::GetGlobalTrack -> No Global Track Found!"<<endl;\r
 799         \r
 800         return partner;\r
 801         \r
 802 }\r
 803 \r
 804 //_____________________________________________________________________________\r
 805 Double_t AliAnalysisTaskDiHadronPID::PhiRange(Double_t DPhi)\r
 806 \r
 807 {\r
 808         //\r
 809         // Puts the argument in the range [-pi/2,3 pi/2].\r
 810         //\r
 811         \r
 812         if (DPhi < -TMath::Pi()/2) DPhi += 2*TMath::Pi();\r
 813         if (DPhi > 3*TMath::Pi()/2) DPhi -= 2*TMath::Pi();      \r
 814 \r
 815         return DPhi;\r
 816         \r
 817 }\r
 818 \r
 819 //_____________________________________________________________________________\r
 820 void AliAnalysisTaskDiHadronPID::UserExec(Option_t *)\r
 821 \r
 822 {\r
 823         //\r
 824         // UserExec.\r
 825         //\r
 826         \r
 827         // Input the event.\r
 828         fAODEvent = dynamic_cast<AliAODEvent*>(InputEvent());\r
 829         if (!fAODEvent) {\r
 830                 if (fVerbose>0) cout << "AliAnalysisTaskDiHadronPID::UserExec -> ERROR: No AliAODEvent pointer could be created." << endl;\r
 831                 return;\r
 832         }\r
 833         \r
 834         // Get the event header.\r
 835         fAODHeader = fAODEvent->GetHeader();\r
 836         if (!fAODHeader) {\r
 837                 if (fVerbose>0) cout << "AliAnalysisTaskDiHadronPID::UserExec -> ERROR: No AliAODHeader pointer could be created."<<endl;\r
 838                 return;\r
 839         }\r
 840         \r
 841         // Get the event vertex.\r
 842         fAODVertex = fAODEvent->GetPrimaryVertex();\r
 843         if (!fAODVertex) {\r
 844                 if (fVerbose>0) cout << "AliAnalysisTaskDiHadronPID::UserExec -> ERROR: No AliAODVertex pointer could be created." << endl;\r
 845                 return;\r
 846         }\r
 847 \r
 848     // See if the event passes the event selection.\r
 849     if (!SelectEvent(fAODVertex)) return;\r
 850     \r
 851         // Display basic event information.\r
 852         if ((fVerbose>2)) cout << endl;\r
 853         if ((fVerbose>2)&&(fBeamType=="PbPb")) cout << "Event centrality: " << fAODHeader->GetCentrality() << endl;\r
 854         if ((fVerbose>2)) cout << "Event Vertex Z: " << fAODVertex->GetZ() << endl;\r
 855         if ((fVerbose>2)) cout << "Event tracks in AOD: " << fAODEvent->GetNumberOfTracks() << endl;\r
 856         if ((fVerbose>2)) cout << endl;\r
 857 \r
 858     // Filling Event QA plots.\r
 859         if (fBeamType=="PbPb") fCentrality->Fill(fAODHeader->GetCentrality());\r
 860         fVertexZ->Fill(fAODVertex->GetZ());\r
 861     \r
 862         // Fill the TObjArray which holds Global tracks.\r
 863         FillGlobalTracksArray();\r
 864         \r
 865         // Create object arrays for triggers and associateds.\r
 866         TObjArray *triggers             = new TObjArray();\r
 867         TObjArray *associateds  = new TObjArray();\r
 868                     \r
 869     // In this loop the triggers and associateds will be identified, track QA and PID QA histograms will be filled.\r
 870         for (Int_t iTrack = 0; iTrack < fAODEvent->GetNumberOfTracks(); iTrack++) {\r
 871         \r
 872         // Obtain a pointer to the track.\r
 873                 fAODTrack = fAODEvent->GetTrack(iTrack);\r
 874         if (!fAODTrack&&(fVerbose>0)) {\r
 875             cout << "AliAnalysisTaskDiHadronPID::UserExec -> ERROR: Track object not found." << endl;\r
 876             continue;\r
 877         }\r
 878         \r
 879         // Find the track classification.\r
 880         Int_t tracktype = ClassifyTrack(fAODTrack);\r
 881                 \r
 882         if (tracktype==0) {\r
 883                         continue;\r
 884                 }\r
 885                 \r
 886         if (tracktype==1) {\r
 887                         if (fVerbose>3) cout<<"Track added to associated buffer."<<endl;\r
 888             associateds->AddLast(fAODTrack);\r
 889             fEtaSpectrumAssoc->Fill(fAODTrack->Eta(),fAODTrack->Pt());\r
 890             fPhiSpectrumAssoc->Fill(fAODTrack->Phi(),fAODTrack->Pt());\r
 891         }\r
 892         \r
 893         if (tracktype==2) {\r
 894                         if (fVerbose>3) cout<<"Track added to trigger buffer."<<endl;\r
 895             triggers->AddLast(fAODTrack);\r
 896             fEtaSpectrumTrig->Fill(fAODTrack->Eta());\r
 897         }\r
 898         \r
 899         \r
 900     }\r
 901     \r
 902     // In This Loop the di-hadron correlation will be made.\r
 903     Double_t histoFill[4];                  // {Dphi, Deta, TPC signal, TOF signal}\r
 904     AliAODTrack* currentTrigger = 0x0;\r
 905         AliAODTrack* currentAssociated = 0x0;\r
 906         AliAODTrack* currentAssociatedGlobal = 0x0;\r
 907     \r
 908     AliTPCPIDResponse& TPCPIDResponse = fPIDResponse->GetTPCResponse();\r
 909     \r
 910     for (Int_t iTrig = 0; iTrig < triggers->GetEntriesFast(); iTrig++){\r
 911         \r
 912                 currentTrigger = (AliAODTrack*)(triggers->At(iTrig));\r
 913                 \r
 914                 for (Int_t iAssoc = 0; iAssoc < associateds->GetEntriesFast(); iAssoc++) {\r
 915             \r
 916                         currentAssociated = (AliAODTrack*)(associateds->At(iAssoc));\r
 917                         currentAssociatedGlobal = GetGlobalTrack(currentAssociated);\r
 918             \r
 919                         Double_t pt = currentAssociated->Pt();\r
 920                         histoFill[0] = PhiRange(currentTrigger->Phi() - currentAssociated->Phi());\r
 921                         histoFill[1] = currentTrigger->Eta() - currentAssociated->Eta();\r
 922 \r
 923                         // Is there a caveat here when Pt = 5.00000000?\r
 924                         const Int_t ptbin = (Int_t)(2*pt);\r
 925                         \r
 926                         if (currentAssociatedGlobal) {\r
 927                                 \r
 928                 // Get TPC (expected) signals.\r
 929                                 Double_t TPCmom = currentAssociatedGlobal->GetTPCmomentum();\r
 930                 Double_t TPCsignal = currentAssociatedGlobal->GetTPCsignal();\r
 931                                 Double_t expectedTPCsignalPion = TPCPIDResponse.GetExpectedSignal(TPCmom,AliPID::kPion);\r
 932                 Double_t expectedTPCsignalKaon = TPCPIDResponse.GetExpectedSignal(TPCmom,AliPID::kKaon);\r
 933                                 Double_t expectedTPCsignalProton = TPCPIDResponse.GetExpectedSignal(TPCmom,AliPID::kProton);\r
 934 \r
 935                 // Get TOF (expected) signals.\r
 936                 Double_t TOFsignal = currentAssociatedGlobal->GetTOFsignal();\r
 937                 Double_t times[AliPID::kSPECIES];\r
 938                 currentAssociatedGlobal->GetIntegratedTimes(times);\r
 939                 Double_t expectedTOFsignalPion = times[AliPID::kPion];\r
 940                 Double_t expectedTOFsignalKaon = times[AliPID::kKaon];\r
 941                 Double_t expectedTOFsignalProton = times[AliPID::kProton]; \r
 942             \r
 943                 // Fill the histograms.\r
 944                 histoFill[2] = TPCsignal - expectedTPCsignalPion;\r
 945                 histoFill[3] = TOFsignal - expectedTOFsignalPion;\r
 946                 fDiHadronTPCTOF[0][ptbin]->Fill(histoFill);\r
 947                 \r
 948                 histoFill[2] = TPCsignal - expectedTPCsignalKaon;\r
 949                 histoFill[3] = TOFsignal - expectedTOFsignalKaon;\r
 950                 fDiHadronTPCTOF[1][ptbin]->Fill(histoFill);\r
 951                 \r
 952                 histoFill[2] = TPCsignal - expectedTPCsignalProton;\r
 953                 histoFill[3] = TOFsignal - expectedTOFsignalProton;\r
 954                 fDiHadronTPCTOF[2][ptbin]->Fill(histoFill);\r
 955                 \r
 956                         }\r
 957                 }\r
 958         }\r
 959 \r
 960     // In this loop we calculate the mixed events.\r
 961     if (fCalculateMixedEvents) {\r
 962         \r
 963         // Loop over the trigger buffer.\r
 964         if (fVerbose>3) cout << "AliAnalysisTaskDiHadronPID::UserExec -> Mixing the events with "<<fTrigBufferSize<<" triggers from the buffer." <<endl;\r
 965         if (fVerbose>3) cout << "AliAnalysisTaskDiHadronPID::UserExec -> Buffer size: "<<fTrigBufferIndex<<endl;\r
 966         \r
 967         for (Int_t iTrig=0;iTrig<fTrigBufferSize;iTrig++) {\r
 968             \r
 969             // Check if the trigger and the associated have a reconstructed\r
 970             // vertext no further than 2cm apart.\r
 971             \r
 972             // fTrigBuffer[i][0] = z\r
 973             // fTrigBuffer[i][1] = phi\r
 974             // fTrigBuffer[i][2] = eta\r
 975             // fTrigBuffer[i][3] = p_t\r
 976             \r
 977             if (TMath::Abs(fTrigBuffer[iTrig][0]-fAODVertex->GetZ())<fVertexZMixedEvents) {\r
 978                 \r
 979                 if (fVerbose>3) cout<<"AliAnalysisTaskDiHadronPID::UserExec -> Mixing with trigger Z: "<<fTrigBuffer[iTrig][0]<<", Pt: "<<fTrigBuffer[iTrig][3]<<endl;\r
 980                 \r
 981                 for (Int_t iAssoc = 0; iAssoc < associateds->GetEntriesFast(); iAssoc++) {\r
 982                     \r
 983                     currentAssociated = (AliAODTrack*)(associateds->At(iAssoc));\r
 984                     currentAssociatedGlobal = GetGlobalTrack(currentAssociated);\r
 985                     \r
 986                     if (currentAssociatedGlobal) {\r
 987                         \r
 988                         Double_t DPhi = PhiRange(fTrigBuffer[iTrig][1] - currentAssociated->Phi());\r
 989                         Double_t DEta = fTrigBuffer[iTrig][2] - currentAssociated->Eta();\r
 990                         Double_t pt = currentAssociated->Pt();\r
 991                         \r
 992                         fMixedEvents->Fill(DPhi,DEta,pt);\r
 993                         \r
 994                     }\r
 995                 }\r
 996             }\r
 997         }\r
 998     \r
 999         // Copy the triggers from the current event into the buffer.\r
1000         if (fAODVertex->GetZ()<10.) {\r
1001         \r
1002             if (fVerbose>3) cout<<"AliAnalysisTaskDiHadronPID::UserExec -> Copying "<<triggers->GetEntriesFast()<<" triggers with vertex z = "<<fAODVertex->GetZ()<<" to the buffer."<<endl;\r
1003         \r
1004             for (Int_t iTrig = 0; iTrig<triggers->GetEntriesFast(); iTrig++) {\r
1005             \r
1006                 currentTrigger = (AliAODTrack*)(triggers->At(iTrig));\r
1007                 if (fVerbose>3) cout<<"AliAnalysisTaskDiHadronPID::UserExec -> Trigger pt = "<<currentTrigger->Pt()<<endl;\r
1008             \r
1009                 fTrigBuffer[fTrigBufferIndex][0] = fAODVertex->GetZ();\r
1010                 fTrigBuffer[fTrigBufferIndex][1] = currentTrigger->Phi();\r
1011                 fTrigBuffer[fTrigBufferIndex][2] = currentTrigger->Eta();\r
1012                 fTrigBuffer[fTrigBufferIndex][3] = currentTrigger->Pt();\r
1013                 fTrigBufferIndex++;\r
1014                 if (fTrigBufferSize<fTrigBufferMaxSize) {fTrigBufferSize++;} // 250 triggers should be enough to get 10 times more data in mixed events.\r
1015                 if (fTrigBufferIndex==fTrigBufferMaxSize) {fTrigBufferIndex=0;}\r
1016             }\r
1017         }\r
1018     }        \r
1019         \r
1020     if (fPrintBufferSize) cout<<"AliAnalysisTaskDiHadronPID::UserExec -> Trigger buffer index: "<<fTrigBufferIndex<<", and size: "<<fTrigBufferSize<<endl;\r
1021     \r
1022         delete triggers;\r
1023         delete associateds;\r
1024                           \r
1025         PostData(1,fHistoList);\r
1026         \r
1027 }\r
1028 \r
1029 //_____________________________________________________________________________\r
1030 void AliAnalysisTaskDiHadronPID::Terminate(Option_t *)\r
1031 \r
1032 {\r
1033         //\r
1034         // Terminate.\r
1035     //\r
1036     \r
1037 }\r
1038 \r