PWGLF/SPECTRA/Kinks/AliAnalysisPionKinksMCESD.cxx

   1 /**************************************************************************
   2  * Authors: Eftychios Cheiladakis, for his Master Thesis                  *
   3  * at the  Physics Department of Athens University                        *
   4  * under the supervision of  Prof. Martha Spyropoulou-Stassinaki          *
   5  *                                                                        *
   6  * Permission to use, copy, modify and distribute this software and its   *
   7  * documentation strictly for non-commercial purposes is hereby granted   *
   8  * without fee, provided that the above copyright notice appears in all   *
   9  * copies and that both the copyright notice and this permission notice   *
  10  * appear in the supporting documentation. The authors make no claims     *
  11  * about the suitability of this software for any purpose. It is          *
  12  * provided "as is" without express or implied warranty.                  *
  13  **************************************************************************/
  14
  15 //-----------------------------------------------------------------
  16 //                 AliAnalysisPionKinksESDMC class
  17 //       Example of an analysis task for kink topology study
  18 //      pions from kink topology are 'identified' in this code
  19 //-----------------------------------------------------------------
  20
  21 #include "TCanvas.h"
  22 #include "TVector3.h"
  23 #include "TLorentzVector.h"
  24 #include "TMath.h"
  25 #include "TF1.h"
  26 #include "TH1.h"
  27 #include "TH2.h"
  28 #include "TH3.h"
  29 #include "TList.h"
  30 #include "TParticle.h"
  31
  32 #include "AliMCEvent.h"
  33 #include "AliVParticle.h"
  34 #include "AliMCParticle.h"
  35 #include "AliESDEvent.h"
  36 #include "AliESDkink.h"
  37 #include "AliESDpid.h"
  38 #include "AliPID.h"
  39 #include "AliStack.h"
  40
  41 #include "AliAnalysisTask.h"
  42 #include "AliInputEventHandler.h"
  43 #include "AliESDInputHandler.h"
  44 #include "AliPIDResponse.h"
  45 #include "AliAnalysisManager.h"
  46 #include "AliMCEventHandler.h"
  47 #include "AliESDtrackCuts.h"
  48
  49 #include "AliAnalysisPionKinksMCESD.h"
  50
  51 ClassImp(AliAnalysisPionKinksMCESD)
  52
  53 //________________________________________________________________________
  54 AliAnalysisPionKinksMCESD::AliAnalysisPionKinksMCESD(const char *name)
  55 :AliAnalysisTaskSE(name),
  56 fMaxKinkAngKmu(0),
  57 fMaxKinkAngPimu(0), //functions
  58 hMCMult(0),
  59 hMCMultPrim(0),
  60 hMCPtAll(0),
  61 hMCEtaAll(0),
  62 hMCPtPrim(0),
  63 hMCEtaPrim(0),
  64 hMCPt(0),
  65 hMCEta(0),
  66 hMCPdg(0),
  67 hMCMultPiPlus(0),
  68 hMCPtPiPlus(0),
  69 hMCEtaPiPlus(0),
  70 hMCRapidityPiPlus(0),
  71 hMCNDaughtersPlus(0),
  72 hMCRadPiDauPlus(0),
  73 hMCKinkPosZPlus(0),
  74 hMCUIDPiDauPlus(0),
  75 hMCPdgPiNonDecayedPlus(0),
  76 hMCPdgPiDauPlus(0),
  77 hMCQtPlus(0),
  78 hMCKinkAnglePlus(0),
  79 hMCPKinkAngPlus(0),
  80 hMCPdgCodemdPlus(0),
  81 hMCPtmdPlus(0),
  82 hMCPtPimuonPlus(0),
  83 hMCEtaPimuonPlus(0),
  84 hMCRapidityPimuonPlus(0),
  85 hMCQtPimuonPlus(0),
  86 hMCPKinkAngPimuonPlus(0),
  87 hMCPtPiotherPlus(0),
  88 hMCEtaPiotherPlus(0),
  89 hMCRapidityPiotherPlus(0),
  90 hMCQtPiotherPlus(0),
  91 hMCPKinkAngPiotherPlus(0),
  92 hMCMultPiMinus(0),
  93 hMCPtPiMinus(0),
  94 hMCEtaPiMinus(0),
  95 hMCRapidityPiMinus(0),
  96 hMCNDaughtersMinus(0),
  97 hMCRadPiDauMinus(0),
  98 hMCKinkPosZMinus(0),
  99 hMCUIDPiDauMinus(0),
 100 hMCPdgPiNonDecayedMinus(0),
 101 hMCPdgPiDauMinus(0),
 102 hMCQtMinus(0),
 103 hMCKinkAngleMinus(0),
 104 hMCPKinkAngMinus(0),
 105 hMCPdgCodemdMinus(0),
 106 hMCPtmdMinus(0),
 107 hMCPtPimuonMinus(0),
 108 hMCEtaPimuonMinus(0),
 109 hMCRapidityPimuonMinus(0),
 110 hMCQtPimuonMinus(0),
 111 hMCPKinkAngPimuonMinus(0),
 112 hMCPtPiotherMinus(0),
 113 hMCEtaPiotherMinus(0),
 114 hMCRapidityPiotherMinus(0),
 115 hMCQtPiotherMinus(0),
 116 hMCPKinkAngPiotherMinus(0),//MC histograms
 117 hMult(0),
 118 hAcceptedMult(0),
 119 hMultPS(0),
 120 hvtx(0),
 121 hvtxy(0),
 122 hvtyz(0),
 123 hvtxz(0),
 124 hMultPSV(0),
 125 hPtAll(0),
 126 hEtaAll(0),
 127 hTrackPos(0),
 128 hTrackPosxy(0),
 129 hTrackPosyz(0),
 130 hTrackPosxz(0),
 131 //hTPCchi2clusters(0),
 132 //hdcaToVertexXY(0),
 133 //hdcaToVertexZ(0),
 134 hMultPrim(0),
 135 hPtPrim(0),
 136 hEtaPrim(0),
 137 hPrimTrackPos(0),
 138 hPrimTrackPosxy(0),
 139 hPrimTrackPosyz(0),
 140 hPrimTrackPosxz(0),
 141 hPt(0),
 142 hEta(0),
 143 //hRapidity(0),
 144 hPtKink(0),
 145 hEtaKink(0),
 146 hRapidityKink(0),
 147 hPmP(0),
 148 hKinkPosRTPCclusters1(0),
 149 hKinkPosRTPCclusters2(0),
 150 hQt(0),
 151 hKinkAngle(0),
 152 hDCAkink(0),
 153 hPmKinkAng(0),
 154 hKinkPosXY(0),
 155 hKinkPosZY(0),
 156 hKinkPosZR(0),
 157 hKinkPosR(0),
 158 hKinkPosZ(0),
 159 hKinkPosZMCKinkPosZ(0),
 160 hPdgCodemd(0),
 161 hPmd(0),
 162 hMinvPimu(0),
 163 hUIDKinkDau(0),
 164 hdEdx(0),
 165 hPtKinkFake(0),
 166 hEtaKinkFake(0),
 167 hRapidityKinkFake(0),
 168 hPmPFake(0),
 169 hKinkPosRTPCclusters1Fake(0),
 170 hKinkPosRTPCclusters2Fake(0),
 171 hQtFake(0),
 172 hKinkAngleFake(0),
 173 hDCAkinkFake(0),
 174 hPmKinkAngFake(0),
 175 hKinkPosXYFake(0),
 176 hKinkPosZYFake(0),
 177 hKinkPosZRFake(0),
 178 hKinkPosRFake(0),
 179 hKinkPosZFake(0),
 180 hKinkPosZMCKinkPosZFake(0),
 181 hPdgCodemdFake(0),
 182 hPmdFake(0),
 183 hMinvPimuFake(0),
 184 hUIDKinkDauFake(0),
 185 hdEdxFake(0),
 186 hPtPosRSelected(0),
 187 hPdgCodemdZRejected(0),
 188 hPtZSelected(0),
 189 hPdgCodemdAngRejected(0),
 190 hPtAngSelected(0),
 191 hPdgCodemdPmRejected(0),
 192 hPtPmSelected(0),
 193 hPdgCodemdQtLowRejected(0),
 194 hPtGoodKink(0),
 195 hEtaGoodKink(0),
 196 hRapidityGoodKink(0),
 197 hQtGoodKink(0),
 198 hPmGoodKinkAng(0),
 199 hPdgCodemdGoodKink(0),
 200 hPmdGoodKink(0),
 201 hUIDGoodKinkDau(0),
 202 hdEdxGoodKink(0),
 203 hUIDPiDauPlus(0),
 204 hMultPiPlus(0),
 205 hPtPiPlus(0),
 206 hEtaPiPlus(0),
 207 hRapidityPiPlus(0),
 208 hQtPiPlus(0),
 209 hKinkAnglePiPlus(0),
 210 hPmKinkAngPiPlus(0),
 211 hKinkPosXYPiPlus(0),
 212 hKinkPosZRPiPlus(0),
 213 hKinkPosRPiPlus(0),
 214 hDCAkinkPiPlus(0),
 215 hPdgCodemdPiPlus(0),
 216 hPmdPiPlus(0),
 217 hdEdxPiPlus(0),
 218 hQtPimuPlus(0),
 219 hPmKinkAngPimuPlus(0),
 220 hQtPiotherPlus(0),
 221 hPmKinkAngPiotherPlus(0),
 222 hPdgCodemdPiotherPlus(0),
 223 hUIDPiDauMinus(0),
 224 hMultPiMinus(0),
 225 hPtPiMinus(0),
 226 hEtaPiMinus(0),
 227 hRapidityPiMinus(0),
 228 hQtPiMinus(0),
 229 hKinkAnglePiMinus(0),
 230 hPmKinkAngPiMinus(0),
 231 hKinkPosXYPiMinus(0),
 232 hKinkPosZRPiMinus(0),
 233 hKinkPosRPiMinus(0),
 234 hDCAkinkPiMinus(0),
 235 hPdgCodemdPiMinus(0),
 236 hPmdPiMinus(0),
 237 hdEdxPiMinus(0),
 238 hQtPimuMinus(0),
 239 hPmKinkAngPimuMinus(0),
 240 hQtPiotherMinus(0),
 241 hPmKinkAngPiotherMinus(0),
 242 hPdgCodemdPiotherMinus(0),
 243 hPdgCodemdQtRejected(0),
 244 hPtQtSelected(0),
 245 hPdgCodemdMaxAngRejected(0),
 246 hPtMaxAngSelected(0),
 247 hPdgCodemdRTPCclustersRejected(0),
 248 hPtRTPCclustersSelected(0),
 249 hRTPCclustersRTPCclustersSelected(0),
 250 hPdgCodemdMinvRejected(0),
 251 hPtSelected(0),
 252 hEtaSelected(0),
 253 hRapiditySelected(0),
 254 hQtSelected(0),
 255 hKinkAngleSelected(0),
 256 hDCAkinkSelected(0),
 257 hPmKinkAngSelected(0),
 258 hKinkPosXYSelected(0),
 259 hKinkPosZRSelected(0),
 260 hKinkPosRSelected(0),
 261 hPdgCodemdSelected(0),
 262 hPmdSelected(0),
 263 hMinvPimuSelected(0),
 264 hUIDKinkDauSelected(0),
 265 hdEdxSelected(0),
 266 hPtSelectedFake(0),
 267 hEtaSelectedFake(0),
 268 hRapiditySelectedFake(0),
 269 hQtSelectedFake(0),
 270 hKinkAngleSelectedFake(0),
 271 hDCAkinkSelectedFake(0),
 272 hPmKinkAngSelectedFake(0),
 273 hKinkPosXYSelectedFake(0),
 274 hKinkPosZRSelectedFake(0),
 275 hKinkPosRSelectedFake(0),
 276 hPmdSelectedFake(0),
 277 hMinvPimuSelectedFake(0),
 278 hdEdxSelectedFake(0),
 279 hPdgCodemddEdxRejected(0),
 280 hPtPiSelected(0),
 281 hEtaPiSelected(0),
 282 hRapidityPiSelected(0),
 283 hQtPiSelected(0),
 284 hKinkAnglePiSelected(0),
 285 hDCAkinkPiSelected(0),
 286 hPmKinkAngPiSelected(0),
 287 hKinkPosRTPCclusters1PiSelected(0),
 288 hKinkPosRTPCclusters2PiSelected(0),
 289 hKinkPosXYPiSelected(0),
 290 hKinkPosZRPiSelected(0),
 291 hKinkPosRPiSelected(0),
 292 hKinkPosZPiSelected(0),
 293 hPmPPiSelected(0),
 294 hPdgCodemdPiSelected(0),
 295 hPmdPiSelected(0),
 296 hMinvPimuPiSelected(0),
 297 hUIDKinkDauPiSelected(0),
 298 hdEdxPiSelected(0),
 299 hPtPiSelectedPlus(0),
 300 hEtaPiSelectedPlus(0),
 301 hRapidityPiSelectedPlus(0),
 302 hQtPiSelectedPlus(0),
 303 hKinkAnglePiSelectedPlus(0),
 304 hDCAkinkPiSelectedPlus(0),
 305 hPmKinkAngPiSelectedPlus(0),
 306 hKinkPosXYPiSelectedPlus(0),
 307 hKinkPosZRPiSelectedPlus(0),
 308 hPdgCodemdPiSelectedPlus(0),
 309 hPmdPiSelectedPlus(0),
 310 hMinvPimuPiSelectedPlus(0),
 311 hUIDPiDaumuSelectedPlus(0),
 312 hdEdxPiSelectedPlus(0),
 313 hPtPrimPiKinksPlus(0),
 314 hEtaPrimPiKinksPlus(0),
 315 hRapidityPrimPiKinksPlus(0),
 316 hPtSecondPiKinksPlus(0),
 317 hEtaSecondPiKinksPlus(0),
 318 hRapiditySecondPiKinksPlus(0),
 319 hPtNonPiKinksPlus(0),
 320 hEtaNonPiKinksPlus(0),
 321 hRapidityNonPiKinksPlus(0),
 322 hPdgCodemdNonPiKinksPlus(0),
 323 hPtPiSelectedMinus(0),
 324 hEtaPiSelectedMinus(0),
 325 hRapidityPiSelectedMinus(0),
 326 hQtPiSelectedMinus(0),
 327 hKinkAnglePiSelectedMinus(0),
 328 hDCAkinkPiSelectedMinus(0),
 329 hPmKinkAngPiSelectedMinus(0),
 330 hKinkPosXYPiSelectedMinus(0),
 331 hKinkPosZRPiSelectedMinus(0),
 332 hPdgCodemdPiSelectedMinus(0),
 333 hPmdPiSelectedMinus(0),
 334 hMinvPimuPiSelectedMinus(0),
 335 hUIDPiDaumuSelectedMinus(0),
 336 hdEdxPiSelectedMinus(0),
 337 hPtPrimPiKinksMinus(0),
 338 hEtaPrimPiKinksMinus(0),
 339 hRapidityPrimPiKinksMinus(0),
 340 hPtSecondPiKinksMinus(0),
 341 hEtaSecondPiKinksMinus(0),
 342 hRapiditySecondPiKinksMinus(0),
 343 hPtNonPiKinksMinus(0),
 344 hEtaNonPiKinksMinus(0),
 345 hRapidityNonPiKinksMinus(0),
 346 hPdgCodemdNonPiKinksMinus(0),// reconstruction histograms
 347 fListOfHistos(0),
 348 fLowMulcut(-1), fUpMulcut(-1),
 349 cLowPt(0), cRapidityLim(0),
 350 cLowR(0), cUpR(0),
 351 cLowZ(0), cUpZ(0),
 352 cLowKinkAngle(0),
 353 cLowQt(0), cUpQt(0),
 354 cLowInvMass(0), cUpInvMass(0),
 355 cSigmaCut(0),
 356 cPdgKaon(321), cPdgPion(211), cPdgMuon(13), cPdgElectron(11),
 357 cKaonMass(0), cPionMass(0), cMuonMass(0), cElectronMass(0),
 358 nBinsMult(0), hLowMult(0), hUpMult(0),
 359 nBinsPt(0), hLowPt(0), hUpPt(0),
 360 nBinsEta(0), hLowEta(0), hUpEta(0),
 361 nBinsQt(0), hLowQt(0), hUpQt(0),
 362 nBinsPdg(0), hLowPdg(0), hUpPdg(0),
 363 nBinsPdg2(0), hLowPdg2(0), hUpPdg2(0),
 364 nBinsUID(0), hLowUID(0), hUpUID(0),
 365 nBinsR(0), hLowR(0), hUpR(0),
 366 nBinsZ(0), hLowZ(0), hUpZ(0),
 367 nBinsXY(0), hLowXY(0), hUpXY(0),
 368 nBinsAngle(0), hLowAngle(0), hUpAngle(0),
 369 nBinsZV(0), hLowZV(0), hUpZV(0),
 370 nBinsXYV(0), hLowXYV(0), hUpXYV(0),
 371 nBinsInvMass(0), hLowInvMass(0), hUpInvMass(0),
 372 nBinsdEdx(0), hLowdEdx(0), hUpdEdx(0), fPIDResponse(0),
 373 fMaxDCAtoVtxCut(0), fTrackCuts(0)
 374
 375 {
 376 //Multiplicity bins
 377 fMaxDCAtoVtxCut=new AliESDtrackCuts("fMaxDCAtoVtxCut","fMaxDCAtoVtxCut");
 378 fMaxDCAtoVtxCut->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
 379 fMaxDCAtoVtxCut->SetMaxChi2TPCConstrainedGlobal(36);
 380
 381 fTrackCuts = new AliESDtrackCuts("Multiplicity bins","Multiplicity bins");
 382 fTrackCuts->SetMinNClustersTPC(70);
 383 fTrackCuts->SetMaxChi2PerClusterTPC(4);
 384 fTrackCuts->SetAcceptKinkDaughters(kFALSE);
 385 fTrackCuts->SetRequireTPCRefit(kTRUE);
 386 fTrackCuts->SetRequireITSRefit(kTRUE);
 387 fTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny);
 388 fTrackCuts->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
 389 fTrackCuts->SetMaxDCAToVertexZ(2);
 390 fTrackCuts->SetDCAToVertex2D(kFALSE);
 391 fTrackCuts->SetRequireSigmaToVertex(kFALSE);
 392 fTrackCuts->SetEtaRange(-0.8,+0.8);
 393 fTrackCuts->SetPtRange(0.15, 1e10);
 394
 395 //DefineOutput(0, TList::Class());
 396 DefineOutput(1, TList::Class());
 397 }
 398
 399 //________________________________________________________________________
 400 void AliAnalysisPionKinksMCESD::UserCreateOutputObjects() {
 401 fListOfHistos=new TList();
 402
 403 //maximum kink angle for kaons to muons
 404 fMaxKinkAngKmu=new TF1("fMaxKinkAngKmu","((atan([0]*[1]*(1.0/(sqrt((x^2)*(1.0-([1]^2))-([0]^2)*([1]^2))))))*180.)/[2]",1.1,10.0);
 405 fMaxKinkAngKmu->SetParameter(0,cKaonMass);
 406 fMaxKinkAngKmu->SetParameter(1,0.9127037);
 407 fMaxKinkAngKmu->SetParameter(2,TMath::Pi());
 408
 409 //maximum kink angle for pions to muons
 410 fMaxKinkAngPimu=new TF1("fMaxKinkAngPimu","((atan([0]*[1]*(1.0/(sqrt((x^2)*(1.0-([1]^2))-([0]^2)*([1]^2))))))*180.)/[2]",0.1,10.0);
 411 fMaxKinkAngPimu->SetParameter(0,cPionMass);
 412 fMaxKinkAngPimu->SetParameter(1,0.2731374);
 413 fMaxKinkAngPimu->SetParameter(2,TMath::Pi());
 414
 415 //Create histograms
 416 TH1::SetDefaultSumw2();
 417 TH2::SetDefaultSumw2();
 418
 419 //MC histograms
 420 hMCMult = new TH1F("hMCMult", "MC multiplicity; Number of tracks; Number of events", 100, 0.0, 2000);
 421 hMCMultPrim = new TH1F("hMCMultPrim", "MC primary tracks multiplicity; Number of tracks; Number of events", nBinsMult, hLowMult, hUpMult);
 422 hMCPtAll = new TH1F("hMCPtAll", "Transverse momentum of all MC tracks; p_{T} (GeV/c); dN/dp_{T}", nBinsPt, hLowPt, hUpPt);
 423 hMCEtaAll = new TH1F("hMCEtaAll", "Pseudorapidity of all MC tracks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 424 hMCPtPrim = new TH1F("hMCPtPrim", "Transverse momentum of primary MC tracks; p_{T} (GeV/c); dN/dp_{T}",nBinsPt, hLowPt, hUpPt);
 425 hMCEtaPrim = new TH1F("hMCEtaPrim", "Pseudorapidity of primary MC tracks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 426 hMCPt = new TH1F("hMCPt", "Transverse momentum of selected MC tracks; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 427 hMCEta = new TH1F("hMCEta", "Pseudorapidity of selected MC tracks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 428 hMCPdg = new TH1F("hMCPdg", "Pdg code of selected MC tracks; Pdg code; Number of particles", nBinsPdg, hLowPdg, hUpPdg);
 429 hMCMultPiPlus = new TH1F("hMCMultPiPlus", "MC pion multiplicity; Number of pion tracks; Number of events", nBinsMult, hLowMult, hUpMult);
 430 hMCPtPiPlus = new TH1F("hMCPtPiPlus", "Transverse momentum of selected MC pions; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 431 hMCEtaPiPlus = new TH1F("hMCEtaPiPlus", "Pseudorapidity of selected MC piaons; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 432 hMCRapidityPiPlus = new TH1F("hMCRapidityPiPlus", "Pseudorapidity of selected MC pions; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 433 hMCNDaughtersPlus = new TH1F("hMCNDaughtersPlus", "Number of daughters; number of daughers; number of pions", 10,0,10);
 434 hMCRadPiDauPlus = new TH1F("hMCRadPiDauPlus", "Radius of MC daughter generation position; R (cm); dN/dR", nBinsR, hLowR, hUpR);
 435 hMCKinkPosZPlus = new TH1F("hMCKinkPosZPlus", "z position of MC daughter generation vertex; z (cm); dN/dz", 100, 0.0, 500.0);
 436 hMCUIDPiDauPlus = new TH1F("hMCUIDPiDauPlus", "UID (method of production) of MC pion daughters; UID (method of production); Number of particles", 21, 0.0, 20);
 437 hMCPdgPiNonDecayedPlus = new TH1F("hMCPdgPiNonDecayedPlus", "Pdg code of MC pion non-decayed products; Pdg code; Number of particles", nBinsPdg, hLowPdg, hUpPdg);
 438 hMCPdgPiDauPlus = new TH1F("hMCPdgPiDauPlus", "Pdg code of MC pion daughters; Pdg code; Number of particles", nBinsPdg, hLowPdg, hUpPdg);
 439 hMCQtPlus = new TH1F("hMCQtPlus", "Daughter's transverse momentum in mother's frame for MC pion daughters; q_{T} (GeV/c); dN/dq_{T}", nBinsQt, hLowQt, hUpQt);
 440 hMCKinkAnglePlus = new TH1F("hMCKinkAnglePlus", "MC angle between pion mother's and daughter's momentum; #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 441 hMCPKinkAngPlus = new TH2F("hMCPKinkAngPlus", "MC mother's P vs kink angle; P (GeV/c); #theta (#circ); dN/d#theta",   nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 442 hMCPdgCodemdPlus = new TH2F("hMCPdgCodemdPlus", "MC mother vs daughter pdg code; Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 443 hMCPtmdPlus = new TH2F("hMCPtmdPlus", "MC mother vs daughter transverse momentum; Mother's p_{T} (GeV/c); Daughter's p_{T} (GeV/c)",   nBinsPt, hLowPt, hUpPt,   nBinsPt, hLowPt, hUpPt);
 444 hMCPtPimuonPlus = new TH1F("hMCPtPimuonPlus", "Transverse momentum of selected pions MC decaying to muons; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 445 hMCEtaPimuonPlus = new TH1F("hMCEtaPimuonPlus", "Pseudorapidity; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 446 hMCRapidityPimuonPlus = new TH1F("hMCRapidityPimuonPlus", "Rapidity of selected MC tracks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 447 hMCQtPimuonPlus = new TH1F("hMCQtPimuonPlus", "MC daughter's (muon) transverse momentum in mother's (pion) frame; q_{T} (GeV/c); dN/dq_{T}", nBinsQt, hLowQt, hUpQt);
 448 hMCPKinkAngPimuonPlus = new TH1F("hMCPKinkAngPimuonPlus", "MC angle between pion mother's and muon daughter's momentum; #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 449 hMCPtPiotherPlus = new TH1F("hMCPtPiotherPlus", "Transverse momentum of selected pions MC not decaying to muons; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 450 hMCEtaPiotherPlus = new TH1F("hMCEtaPiotherPlus", "Pseudorapidity; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 451 hMCRapidityPiotherPlus = new TH1F("hMCRapidityPiotherPlus", "Rapidity of selected MC tracks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 452 hMCQtPiotherPlus = new TH1F("hMCQtPiotherPlus", "MC daughter's (muon) transverse momentum in mother's (pion) frame; q_{T} (GeV/c); dN/dq_{T}", nBinsQt, hLowQt, hUpQt);
 453 hMCPKinkAngPiotherPlus = new TH1F("hMCPKinkAngPiotherPlus", "MC angle between pion mother's and non-muon daughter's momentum; #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 454 hMCMultPiMinus = new TH1F("hMCMultPiMinus", "MC pion multiplicity; Number of pion tracks; Number of events", nBinsMult, hLowMult, hUpMult);
 455 hMCPtPiMinus = new TH1F("hMCPtPiMinus", "Transverse momentum of selected MC pions; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 456 hMCEtaPiMinus = new TH1F("hMCEtaPiMinus", "Pseudorapidity of selected MC pions; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 457 hMCRapidityPiMinus = new TH1F("hMCRapidityPiMinus", "Pseudorapidity of selected MC pions; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 458 hMCNDaughtersMinus = new TH1F("hMCNDaughtersMinus", "Number of daughters; number of daughers; number of pions", 10,0,10);
 459 hMCRadPiDauMinus = new TH1F("hMCRadPiDauMinus", "Radius of MC daughter generation position; R (cm); dN/dR", nBinsR, hLowR, hUpR);
 460 hMCKinkPosZMinus = new TH1F("hMCKinkPosZMinus", "z position of MC daughter generation vertex; z (cm); dN/dz", 100, 0.0, 500.0);
 461 hMCUIDPiDauMinus = new TH1F("hMCUIDPiDauMinus", "UID (method of production) of MC pion daughters; UID (method of production); Number of particles", 21, 0.0, 20);
 462 hMCPdgPiNonDecayedMinus = new TH1F("hMCPdgPiNonDecayedMinus", "Pdg code of MC pion non-decayed products; Pdg code; Number of particles", nBinsPdg, hLowPdg, hUpPdg);
 463 hMCPdgPiDauMinus = new TH1F("hMCPdgPiDauMinus", "Pdg code of MC pion daughters; Pdg code; Number of particles", nBinsPdg, hLowPdg, hUpPdg);
 464 hMCQtMinus = new TH1F("hMCQtMinus", "Daughter's transverse momentum in mother's frame for MC pion daughters; q_{T} (GeV/c); dN/dq_{T}", nBinsQt, hLowQt, hUpQt);
 465 hMCKinkAngleMinus = new TH1F("hMCKinkAngleMinus", "MC angle between pion mother's and daughter's momentum; #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 466 hMCPKinkAngMinus = new TH2F("hMCPKinkAngMinus", "MC mother's P vs kink angle; P (GeV/c); #theta (#circ); dN/d#theta",   nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 467 hMCPdgCodemdMinus = new TH2F("hMCPdgCodemdMinus", "MC mother vs daughter pdg code; Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 468 hMCPtmdMinus = new TH2F("hMCPtmdMinus", "MC mother vs daughter transverse momentum; Mother's p_{T} (GeV/c); Daughter's p_{T} (GeV/c)",   nBinsPt, hLowPt, hUpPt,   nBinsPt, hLowPt, hUpPt);
 469 hMCPtPimuonMinus = new TH1F("hMCPtPimuonMinus", "Transverse momentum of selected pions MC decaying to muons; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 470 hMCEtaPimuonMinus = new TH1F("hMCEtaPimuonMinus", "Pseudorapidity; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 471 hMCRapidityPimuonMinus = new TH1F("hMCRapidityPimuonMinus", "Rapidity of selected MC tracks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 472 hMCQtPimuonMinus = new TH1F("hMCQtPimuonMinus", "MC daughter's (muon) transverse momentum in mother's (pion) frame; q_{T} (GeV/c); dN/dq_{T}", nBinsQt, hLowQt, hUpQt);
 473 hMCPKinkAngPimuonMinus = new TH1F("hMCPKinkAngPimuonMinus", "MC angle betweenpion mother's and muon daughter's momentum; #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 474 hMCPtPiotherMinus = new TH1F("hMCPtPiotherMinus", "Transverse momentum of selected pions MC not decaying to muons; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 475 hMCEtaPiotherMinus = new TH1F("hMCEtaPiotherMinus", "Pseudorapidity; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 476 hMCRapidityPiotherMinus = new TH1F("hMCRapidityPiotherMinus", "Rapidity of selected MC tracks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 477 hMCQtPiotherMinus = new TH1F("hMCQtPiotherMinus", "MC daughter's (muon) transverse momentum in mother's (pion) frame; q_{T} (GeV/c); dN/dq_{T}", nBinsQt, hLowQt, hUpQt);
 478 hMCPKinkAngPiotherMinus = new TH1F("hMCPKinkAngPiotherMinus", "MC angle betweenpion mother's and non-muon daughter's momentum; #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 479
 480
 481 //Reconstruction histograms
 482 hMult = new TH1F("hMult", "Multiplicity (unbiased); Number of tracks; Number of events", nBinsMult, hLowMult, hUpMult);
 483 hAcceptedMult = new TH1F("hAcceptedMult", "Multiplicity (biased); Number of tracks; Number of events", nBinsMult, hLowMult, hUpMult);
 484 hMultPS = new TH1F("hMultPS", "Multiplicity after physics selection; Number of tracks; Number of events", nBinsMult, hLowMult, hUpMult);
 485 hvtx = new TH3F("hvtx", "Reconstructed primary vertex position; x axis; y axis; z axis", nBinsXYV, hLowXYV, hUpXYV, nBinsXYV, hLowXYV, hUpXYV, nBinsZV, hLowZV, hUpZV);
 486 hvtxy = new TH2F("hvtxy", "Reconstructed primary vertex position in x-y plane; x axis; y axis", nBinsXYV, hLowXYV, hUpXYV, nBinsXYV, hLowXYV, hUpXYV);
 487 hvtyz = new TH2F("hvtyz", "Reconstructed primary vertex position in y-z plane; y axis; z axis", nBinsXYV, hLowXYV, hUpXYV, nBinsZV, hLowZV, hUpZV);
 488 hvtxz = new TH2F("hvtxz", "Reconstructed primary vertex position in x-z plane; x axis; z axis", nBinsXYV, hLowXYV, hUpXYV, nBinsZV, hLowZV, hUpZV);
 489 hMultPSV = new TH1F("hMultPSV", "Multiplicity after physics selection & vertex cut; Number of tracks; Number of events", nBinsMult, hLowMult, hUpMult);
 490 hPtAll = new TH1F("hPtAll", "Transverse momentum of all tracks; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 491 hEtaAll = new TH1F("hEtaAll", "Pseudorapidity of all tracks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 492 hTrackPos = new TH3F("hTrackPos", "Generetion position of all reconstructed tracks", nBinsXYV, hLowXYV, hUpXYV, nBinsXYV, hLowXYV, hUpXYV, nBinsZV, hLowZV, hUpZV);
 493 hTrackPosxy = new TH2F("hTrackPosxy", "Generetion position of all reconstructed tracks in x-y plane; x axis; y axis", nBinsXYV, hLowXYV, hUpXYV, nBinsXYV, hLowXYV, hUpXYV);
 494 hTrackPosyz = new TH2F("hTrackPosyz", "Generetion position of all reconstructed tracks in y-z plane; y axis; z axis", nBinsXYV, hLowXYV, hUpXYV, nBinsZV, hLowZV, hUpZV);
 495 hTrackPosxz = new TH2F("hTrackPosxz", "Generetion position of all reconstructed tracks in x-z plane; x axis; z axis", nBinsXYV, hLowXYV, hUpXYV, nBinsZV, hLowZV, hUpZV);
 496 //hTPCchi2clusters = new TH1F("hTPCchi2clusters", "#chi^{2}/Number of TPC clusters; #chi^{2}; TPC clusters)", 100, 0.0, 2.0);//
 497 //hdcaToVertexXY = new TH1F("hdcaToVertexXY", "Track to vertex impact parameter in x-y plane; DCA_{z} (cm); dN/d(DCA_{z})", 100, 0.0, 2.0);//
 498 //hdcaToVertexZ = new TH1F("hdcaToVertexZ", "Track to vertex impact parameter in z axis; DCA_{z} (cm); dN/d(DCA_{z})", 100, 0.0, 2.0);//
 499 hMultPrim = new TH1F("hMultPrim", "ESD primary - supposed tracks multiplicity; Number of tracks; Number of events", nBinsMult, hLowMult, hUpMult);
 500 hPtPrim = new TH1F("hPtPrim", "Transverse momentum of primary - supposed ESD tracks; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 501 hEtaPrim = new TH1F("hEtaPrim", "Pseudorapidity of primary - supposed ESD tracks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 502 hPrimTrackPos = new TH3F("hPrimTrackPos", "Generetion position of selected tracks (DCA and quality cuts)", nBinsXYV, hLowXYV, hUpXYV, nBinsXYV, hLowXYV, hUpXYV, nBinsZV, hLowZV, hUpZV);
 503 hPrimTrackPosxy = new TH2F("hPrimTrackPosxy", "Generetion position of selected tracks in x-y plane (DCA and quality cuts); x axis; y axis", nBinsXYV, hLowXYV, hUpXYV, nBinsXYV, hLowXYV, hUpXYV);
 504 hPrimTrackPosyz = new TH2F("hPrimTrackPosyz", "Generetion position of selected tracks in y-z plane (DCA and quality cuts); y axis; z axis", nBinsXYV, hLowXYV, hUpXYV, nBinsZV, hLowZV, hUpZV);
 505 hPrimTrackPosxz = new TH2F("hPrimTrackPosxz", "Generetion position of selected tracks in x-z plane (DCA and quality cuts); x axis; z axis", nBinsXYV, hLowXYV, hUpXYV, nBinsZV, hLowZV, hUpZV);
 506 hPt = new TH1F("hPt", "Transverse momentum of selected ESD tracks; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 507 hEta = new TH1F("hEta", "Pseudorapidity of selected ESD tracks; n; Number of tracks dN/dn", nBinsEta, hLowEta, hUpEta);
 508 //hRapidity = new TH1F("hRapidity", "Rapidity of selected ESD tracks; n; Number of tracks dN/dn", nBinsEta, hLowEta, hUpEta);
 509 hPtKink = new TH1F("hPtKink", "Mother's transverse momentum for all ESD kinks; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 510 hEtaKink = new TH1F("hEtaKink", "Mother's pseudorapidity for all ESD kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 511 hRapidityKink = new TH1F("hRapidityKink", "Mother's rapidity for all ESD kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 512 hPmP = new TH2F("hPmP", "Mother's momentum as calculated by kink and by track; P_{kink}; P_{track}", nBinsPt, hLowPt, hUpPt, nBinsPt, hLowPt, hUpPt);
 513 hKinkPosRTPCclusters1 = new TH1F("hKinkPosRTPCclusters1", " ;kinkposR; tpc clusters",100,0,1);
 514 hKinkPosRTPCclusters2 = new TH2F("hKinkPosRTPCclusters2", " ;kinkposR; tpc clusters",nBinsR, hLowR, hUpR,100,0,200 );
 515 hQt = new TH1F("hQt", "Daughter's transverse momentum in mother's frame for all ESD kinks; q_{T} (GeV/c); dN/dq_{T}", nBinsQt, hLowQt, hUpQt);
 516 hKinkAngle = new TH1F("hKinkAngle", "Kink angle for all ESD kinks; #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 517 hDCAkink = new TH1F("hDCAkink", "DCA between the two kink tracks; DCA(cm); Number of kinks", 100, 0.0, 2.0);
 518 hPmKinkAng = new TH2F("hPmKinkAng", "k, p_(GeV/c);  #theta (#circ); d^{2}N/dpd#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 519 hKinkPosXY = new TH2F("hKinkPosXY", "X-Y Position of all kinks; X (cm); Y (cm)", nBinsXY, hLowXY, hUpXY, nBinsXY, hLowXY, hUpXY);
 520 hKinkPosZY = new TH2F("hKinkPosZY", "z vs y position of kinks (DCA, quality, p_{T} & y cuts); z (cm); y (cm)", nBinsZ, hLowZ, hUpZ, nBinsXY, hLowXY, hUpXY);
 521 hKinkPosZR = new TH2F("hKinkPosZR", "Z vs radius of all kinks position; Z (cm); R (cm)", nBinsZ, hLowZ, hUpZ, nBinsR, hLowR, hUpR);
 522 hKinkPosR = new TH1F("hKinkPosR", "Position radius of all ESD kinks; R (cm); dN/dR", nBinsR, hLowR, hUpR);
 523 hKinkPosZ = new TH1F("hKinkPosZ", "z position of kinks (DCA, quality, p_{T} & y cuts); z (cm); dN/dz", nBinsZ, -1, 1);
 524 hKinkPosZMCKinkPosZ = new TH2F("hKinkPosZMCKinkPosZ", "Reconstructed vs generated z position of kinks (DCA, quality, p_{T} & y cuts); z (cm); dN/dz", nBinsZ, hLowZ, hUpZ, nBinsZ, hLowZ, hUpZ);
 525 hPdgCodemd = new TH2F("hPdgCodemd", "Mother vs daughter pdg code for all ESD kinks; Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 526 hPmd = new TH2F("hPmd", "ESD mother vs daughter momentum magnitude; Mother's P (GeV/c); Daughter's P (GeV/c)", nBinsPt, hLowPt, hUpPt, nBinsPt, hLowPt, hUpPt);
 527 hMinvPimu = new TH1F("hMinvPimu", "Invariant mass of ESD pions decaying to muons; m (GeV/c^{2}); dN/dm", nBinsInvMass, hLowInvMass, hUpInvMass);
 528 hUIDKinkDau = new TH1F("hUIDKinkDau", "UID (method of production) of all ESD pion kink daughters (perfect PID); UID (method of production); Number of particles", 21, 0.0, 20);
 529 hdEdx = new TH2F("hdEdx", "dE/dx vs mother's momentum; p (GeV/c); dE/dx (GeV/cm)",   nBinsPt, hLowPt, hUpPt, 100, 0.0, 300.0);
 530 hPtKinkFake = new TH1F("hPtKinkFake", "Mother's transverse momentum for all fake kinks; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 531 hEtaKinkFake = new TH1F("hEtaKinkFake", "Mother's pseudorapidity for all fake kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 532 hRapidityKinkFake = new TH1F("hRapidityKinkFake", "Mother's rapidity for all fake kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 533 hPmPFake = new TH2F("hPmPFake", "Mother's momentum magnitude calculated from kink vs calculated from track for all fake kinks; p_{kink} (GeV/c); p_{track} (GeV/c)",   nBinsPt, hLowPt, hUpPt,   nBinsPt, hLowPt, hUpPt);
 534 hKinkPosRTPCclusters1Fake = new TH1F("hKinkPosRTPCclusters1Fake", " ;kinkposR; tpc clusters",100,0,1);
 535 hKinkPosRTPCclusters2Fake = new TH2F("hKinkPosRTPCclusters2Fake", " ;kinkposR; tpc clusters",nBinsR, hLowR, hUpR,100,0,200 );
 536 hQtFake = new TH1F("hQtFake", "Daughter's transverse momentum in mother's frame for all fake kinks; q_{T} (GeV/c); dN/dq_{T}", nBinsQt, hLowQt, hUpQt);
 537 hKinkAngleFake = new TH1F("hKinkAngleFake", "Kink angle for all fake kinks; #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 538 hDCAkinkFake = new TH1F("hDCAkinkFake", "DCA between the two fake kink tracks; DCA(cm); Number of kinks", 100, 0.0, 2.0);
 539 hPmKinkAngFake = new TH2F("hPmKinkAngFake", "k, p_(GeV/c);  #theta (#circ); d^{2}N/dpd#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 540 hKinkPosXYFake = new TH2F("hKinkPosXYFake", "X-Y Position of all fake kinks; X (cm); Y (cm)", nBinsXY, hLowXY, hUpXY, nBinsXY, hLowXY, hUpXY);
 541 hKinkPosZYFake = new TH2F("hKinkPosZYFake", "z vs y position of fake kinks (DCA, quality, p_{T} & y cuts); z (cm); y (cm)", nBinsZ, hLowZ, hUpZ, nBinsXY, hLowXY, hUpXY);
 542 hKinkPosZRFake = new TH2F("hKinkPosZRFake", "Z vs radius of all fake kinks position; Z (cm); R (cm)", nBinsZ, hLowZ, hUpZ, nBinsR, hLowR, hUpR);
 543 hKinkPosRFake = new TH1F("hKinkPosRFake", "Position radius of all fake kinks; R (cm); dN/dR", nBinsR, hLowR, hUpR);
 544 hKinkPosZFake = new TH1F("hKinkPosZFake", "z position of fake kinks (DCA, quality, p_{T} & y cuts); z (cm); dN/dz", nBinsZ, -1, 1);
 545 hKinkPosZMCKinkPosZFake = new TH2F("hKinkPosZMCKinkPosZFake", "Reconstructed vs generated z position of fake kinks (DCA, quality, p_{T} & y cuts); z (cm); dN/dz", nBinsZ, hLowZ, hUpZ, nBinsZ, hLowZ, hUpZ);
 546 hPdgCodemdFake = new TH2F("hPdgCodemdFake", "Mother vs daughter pdg code for fakes; Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 547 hPmdFake = new TH2F("hPmdFake", "ESD mother vs daughter momentum magnitude (fake kinks); Mother's P (GeV/c); Daughter's P (GeV/c)", nBinsPt, hLowPt, hUpPt, nBinsPt, hLowPt, hUpPt);
 548 hMinvPimuFake = new TH1F("hMinvPimuFake", "Invariant mass of ESD pions decaying to muons (fake kinks); m (GeV/c^{2}); dN/dm", nBinsInvMass, hLowInvMass, hUpInvMass);
 549 hUIDKinkDauFake = new TH1F("hUIDKinkDauFake", "UID (method of production) of all ESD pion kink daughters (perfect PID); UID (method of production); Number of particles", 21, 0.0, 20);
 550 hdEdxFake = new TH2F("hdEdxFake", "dE/dx vs mother's momentum; p (GeV/c); dE/dx (GeV/cm)",   nBinsPt, hLowPt, hUpPt, 100, 0.0, 300.0);
 551 hPtPosRSelected = new TH1F("hPtPosRSelected", "Mother's transverse momentum for selected kinks (DCA, quality, p_{T}, y & R cuts); R (cm); dN/dR", nBinsPt, hLowPt, hUpPt);
 552 hPdgCodemdZRejected = new TH2F("hPdgCodemdZRejected", "Mother's vs daughter's pdg code for rejected kinks (DCA, quality, p_{T}, y, R & z cuts); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 553 hPtZSelected = new TH1F("hPtZSelected",  "Mother's transverse momentum for selected kinks (DCA, quality, p_{T}, y, R & z cuts); p_{T} (GeV/c); dN/dp_{T}", nBinsPt, hLowPt, hUpPt);
 554 hPdgCodemdAngRejected = new TH2F("hPdgCodemdAngRejected", "Mother's vs daughter's pdg code for rejected kinks (DCA, quality, p_{T}, y, R, z & #theta cuts); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 555 hPtAngSelected = new TH1F("hPtAngSelected", "Mother's transverse momentum for selected kinks (DCA, quality, p_{T}, y, R, z & #theta cuts); p_{T} (GeV/c); dN/dp_{T}", nBinsPt, hLowPt, hUpPt);
 556 hPdgCodemdPmRejected = new TH2F("hPdgCodemdPmRejected", "Mother's vs daughter's pdg code for rejected kinks (DCA, quality, p_{T}, y, R, z, #theta & p_{track}/p{kink} cuts); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 557 hPtPmSelected = new TH1F("hPtPmSelected", "Mother's transverse momentum for selected kinks (DCA, quality, p_{T}, y, R, z, #theta & p_{track}/p{kink} cuts); p_{T} (GeV/c); dN/dp_{T}", nBinsPt, hLowPt, hUpPt);
 558 hPdgCodemdQtLowRejected = new TH2F("hPdgCodemdQtLowRejected", "Mother's vs daughter's pdg code for rejected kinks (DCA, quality, p_{T}, y, R, z, #theta, p_{track}/p{kink} & low q_{T} cuts); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 559 hPtGoodKink = new TH1F("hPtGoodKink", "Mother's transverse momentum for real ESD kinks (realistic PID); p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 560 hEtaGoodKink = new TH1F("hEtaGoodKink", "Mother's pseudorapidity for real ESD kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 561 hRapidityGoodKink = new TH1F("hRapidityGoodKink", "Mother's rapidity for real ESD kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 562 hQtGoodKink = new TH1F("hQtGoodKink", "Daughter's transverse momentum in mother's frame for real ESD kinks (realistic PID); q_{T} (GeV/c); dN/dq_{T}", 200, 0.0, 0.3);
 563 hPmGoodKinkAng = new TH2F("hPmGoodKinkAng", "Mother's momentum magnitude vs kink angle for real ESD kinks (realistic PID); p_{T} (GeV/c);  #theta (#circ); d^{2}N/dp_{T}d#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 564 hPdgCodemdGoodKink = new TH2F("hPdgCodemdGoodKink", "Mother vs daughter pdg code for real ESD kinks; Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 565 hPmdGoodKink = new TH2F("hPmdGoodKink", "Mother vs daughter momentum magnitude for real ESD kinks (realistic PID); Mother's P; Daughter's P", nBinsPt, hLowPt, hUpPt, nBinsPt, hLowPt, hUpPt);
 566 hUIDGoodKinkDau = new TH1F("hUIDGoodKinkDau", "UID (method of production); UID (method of production); Number of particles", 21, 0.0, 20);
 567 hdEdxGoodKink = new TH2F("hdEdxGoodKink", "dE/dx vs mother's momentum; p (GeV/c); dE/dx (GeV/cm)",   nBinsPt, hLowPt, hUpPt, 100, 0.0, 300.0);
 568 hUIDPiDauPlus = new TH1F("hUIDPiDauPlus", "UID (method of production) of all ESD pion kink daughters (perfect PID); UID (method of production); Number of particles", 21, 0.0, 20);
 569 hMultPiPlus = new TH1F("hMultPiPlus", "ESD pion multiplicity in selected TPC area; Number of pions; Number of events", 100, 0.0, 200.0);
 570 hPtPiPlus = new TH1F("hPtPiPlus", "Mother's transverse momentum for ESD pion kinks (perfect PID); p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 571 hEtaPiPlus = new TH1F("hEtaPiPlus", "Mother's pseudorapidity for ESD pion kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 572 hRapidityPiPlus = new TH1F("hRapidityPiPlus", "Mother's rapidity for ESD pion kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 573 hQtPiPlus = new TH1F("hQtPiPlus", "Daughter's transverse momentum in mother's frame for ESD pion kinks (perfect PID); q_{T} (GeV/c); dN/dq_{T}", 100, 0.0, 0.1);
 574 hKinkAnglePiPlus = new TH1F("hKinkAnglePiPlus", "Kink angle of ESD pion kinks (perfect PID); #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 575 hPmKinkAngPiPlus = new TH2F("hPmKinkAngPiPlus", "Mother's momentum magnitude vs kink angle for ESD pion kinks (perfect PID); p_{T} (GeV/c); #theta (#circ); d^{2}N/dp_{T}d#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 576 hKinkPosXYPiPlus = new TH2F("hKinkPosXYPiPlus", "X-Y Position ; X (cm); Y (cm)", nBinsXY, hLowXY, hUpXY, nBinsXY, hLowXY, hUpXY);
 577 hKinkPosZRPiPlus = new TH2F("hKinkPosZRPiPlus", "Z vs radius of all fake kinks position; Z (cm); R (cm)", nBinsZ, hLowZ, hUpZ, nBinsR, hLowR, hUpR);
 578 hKinkPosRPiPlus = new TH1F("hKinkPosRPiPlus", "Position radius of all fake kinks; R (cm); dN/dR", nBinsR, hLowR, hUpR);
 579 hDCAkinkPiPlus = new TH1F("hDCAkinkPiPlus", "DCA; DCA(cm); Number of kinks", 100, 0.0, 2.0);
 580 hPdgCodemdPiPlus = new TH2F("hPdgCodemdPiPlus", "Mother vs daughter pdg code for pi kinks; Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 581 hPmdPiPlus = new TH2F("hPmdPiPlus", "Mother vs daughter momentum magnitude for pi kinks (realistic PID); Mother's P; Daughter's P", nBinsPt, hLowPt, hUpPt, nBinsPt, hLowPt, hUpPt);
 582 hdEdxPiPlus = new TH2F("hdEdxPiPlus", "dE/dx vs mother's momentum; p (GeV/c); dE/dx (GeV/cm)",   nBinsPt, hLowPt, hUpPt, 100, 0.0, 300.0);
 583 hQtPimuPlus = new TH1F("hQtPimuPlus", "Daughter's transverse momentum in mother's frame for ESD pion to muon kinks (perfect PID); q_{T} (GeV/c); dN/dq_{T}", 100, 0.0, 0.1);
 584 hPmKinkAngPimuPlus = new TH2F("hPmKinkAngPimuPlus", "Mother's momentum magnitude vs kink angle for ESD pion to muon kinks (perfect PID); p_{T} (GeV/c);  #theta (#circ); d^{2}N/dp_{T}d#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 585 hQtPiotherPlus = new TH1F("hQtPiotherPlus", "Daughter's transverse momentum in mother's frame for ESD pion to other (not muon) kinks (perfect PID); q_{T} (GeV/c); dN/dq_{T}", 100, 0.0, 0.1);
 586 hPmKinkAngPiotherPlus = new TH2F("hPmKinkAngPiotherPlus", "Mother's momentum magnitude vs kink angle for ESD pion to other (not muon) kinks (perfect PID); p_{T} (GeV/c);  #theta (#circ); d^{2}N/dp_{T}d#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 587 hPdgCodemdPiotherPlus = new TH2F("hPdgCodemdPiotherPlus", "Mother vs daughter pdg code for ESD pion to other (not muon) kinks (perfect PID); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 588 hUIDPiDauMinus = new TH1F("hUIDPiDauMinus", "UID (method of production) of all ESD pion kink daughters (perfect PID); UID (method of production); Number of particles", 21, 0.0, 20);
 589 hMultPiMinus = new TH1F("hMultPiMinus", "ESD pion multiplicity in selected TPC area; Number of pions; Number of events", 100, 0.0, 200.0);
 590 hPtPiMinus = new TH1F("hPtPiMinus", "Mother's transverse momentum for ESD pion kinks (perfect PID); p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 591 hEtaPiMinus = new TH1F("hEtaPiMinus", "Mother's pseudorapidity for ESD pion kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 592 hRapidityPiMinus = new TH1F("hRapidityPiMinus", "Mother's rapidity for ESD pion kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 593 hQtPiMinus = new TH1F("hQtPiMinus", "Daughter's transverse momentum in mother's frame for ESD pion kinks (perfect PID); q_{T} (GeV/c); dN/dq_{T}", 100, 0.0, 0.1);
 594 hKinkAnglePiMinus = new TH1F("hKinkAnglePiMinus", "Kink angle of ESD pion kinks (perfect PID); #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 595 hPmKinkAngPiMinus = new TH2F("hPmKinkAngPiMinus", "Mother's momentum magnitude vs kink angle for ESD pion kinks (perfect PID); p_{T} (GeV/c); #theta (#circ); d^{2}N/dp_{T}d#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 596 hKinkPosXYPiMinus = new TH2F("hKinkPosXYPiMinus", "X-Y Position ; X (cm); Y (cm)", nBinsXY, hLowXY, hUpXY, nBinsXY, hLowXY, hUpXY);
 597 hKinkPosZRPiMinus = new TH2F("hKinkPosZRPiMinus", "Z vs radius of all fake kinks position; Z (cm); R (cm)", nBinsZ, hLowZ, hUpZ, nBinsR, hLowR, hUpR);
 598 hKinkPosRPiMinus = new TH1F("hKinkPosRPiMinus", "Position radius of all fake kinks; R (cm); dN/dR", nBinsR, hLowR, hUpR);
 599 hDCAkinkPiMinus = new TH1F("hDCAkinkPiMinus", "DCA; DCA(cm); Number of kinks", 100, 0.0, 2.0);
 600 hPdgCodemdPiMinus = new TH2F("hPdgCodemdPiMinus", "Mother vs daughter pdg code for pi kinks; Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 601 hPmdPiMinus = new TH2F("hPmdPiMinus", "Mother vs daughter momentum magnitude for pi kinks (realistic PID); Mother's P; Daughter's P", nBinsPt, hLowPt, hUpPt, nBinsPt, hLowPt, hUpPt);
 602 hdEdxPiMinus = new TH2F("hdEdxPiMinus", "dE/dx vs mother's momentum; p (GeV/c); dE/dx (GeV/cm)",   nBinsPt, hLowPt, hUpPt, 100, 0.0, 300.0);
 603 hQtPimuMinus = new TH1F("hQtPimuMinus", "Daughter's transverse momentum in mother's frame for ESD pion to muon kinks (perfect PID); q_{T} (GeV/c); dN/dq_{T}", 100, 0.0, 0.1);
 604 hPmKinkAngPimuMinus = new TH2F("hPmKinkAngPimuMinus", "Mother's momentum magnitude vs kink angle for ESD pion to muon kinks (perfect PID); p_{T} (GeV/c);  #theta (#circ); d^{2}N/dp_{T}d#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 605 hQtPiotherMinus = new TH1F("hQtPiotherMinus", "Daughter's transverse momentum in mother's frame for ESD pion to other (not muon) kinks (perfect PID); q_{T} (GeV/c); dN/dq_{T}", 100, 0.0, 0.1);
 606 hPmKinkAngPiotherMinus = new TH2F("hPmKinkAngPiotherMinus", "Mother's momentum magnitude vs kink angle for ESD pion to other (not muon) kinks (perfect PID); p_{T} (GeV/c);  #theta (#circ); d^{2}N/dp_{T}d#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 607 hPdgCodemdPiotherMinus = new TH2F("hPdgCodemdPiotherMinus", "Mother vs daughter pdg code for ESD pion to other (not muon) kinks (perfect PID); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 608 hPdgCodemdQtRejected = new TH2F("hPdgCodemdQtRejected", "Mother's vs daughter's pdg code for rejected kinks (DCA, quality, p_{T}, y, R, z, #theta, p_{track}/p{kink} & q^{T} cuts); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 609 hPtQtSelected = new TH1F("hPtQtSelected", "Mother's transverse momentum for selected kinks (DCA, quality, p_{T}, y, R, z, #theta, p_{track}/p{kink} & q^{T} cuts); p_{T} (GeV/c); dN/dp_{T}", nBinsPt, hLowPt, hUpPt);
 610 hPdgCodemdMaxAngRejected = new TH2F("hPdgCodemdMaxAngRejected", "Mother's vs daughter's pdg code for rejected kinks (DCA, quality, p_{T}, y, R, z, #theta, p_{track}/p{kink}, q^{T} & #theta_{max} cuts); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 611 hPtMaxAngSelected = new TH1F("hPtMaxAngSelected", "Mother's transverse momentum for selected kinks (DCA, quality, p_{T}, y, R, z, #theta, p_{track}/p{kink}, q^{T} & #theta_{max} cuts); p_{T} (GeV/c); dN/dp_{T}", nBinsPt, hLowPt, hUpPt);
 612 hPdgCodemdRTPCclustersRejected = new TH2F("hPdgCodemdRTPCclustersRejected", "Mother's vs daughter's pdg code for rejected kinks (DCA, quality, p_{T}, y, R, z, #theta, p_{track}/p{kink}, q^{T}, #theta_{max} & R/TPC clusters cuts); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 613 hPtRTPCclustersSelected = new TH1F("hPtRTPCclustersSelected", "Mother's transverse momentum for selected kinks (DCA, quality, p_{T}, y, R, z, #theta, p_{track}/p{kink}, q^{T}, #theta_{max} & R/TPC clusters cuts); p_{T} (GeV/c); dN/dp_{T}", nBinsPt, hLowPt, hUpPt);
 614 hRTPCclustersRTPCclustersSelected = new TH2F("hRTPCclustersRTPCclustersSelected", "Number of TPC clusters vs radius for selected kinks (DCA, quality, p_{T}, y, R, z, #theta, p_{track}/p{kink}, q^{T}, #theta_{max} & R/TPC clusters cuts) R (cm); number of TPC clusters", nBinsR, hLowR, hUpR, 100, 0, 200);
 615 hPdgCodemdMinvRejected = new TH2F("hPdgCodemdMinvRejected", "Mother's vs daughter's pdg code for rejected kinks (DCA, quality, p_{T}, y, R, z, #theta, p_{track}/p{kink}, q^{T}, #theta_{max}, R/TPC clusters & m_{inv} cuts); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 616 hPtSelected = new TH1F("hPtSelected", "Mother's transverse momentum for selected kinks (all cuts except dE/dx); p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 617 hEtaSelected = new TH1F("hEtaSelected", "Mother's pseudorapidity for selected kinks (all cuts except dE/dx); n; dN/dn", nBinsEta, hLowEta, hUpEta);
 618 hRapiditySelected = new TH1F("hRapiditySelected", "Mother's rapidity for selected kinks (all cuts except dE/dx); n; dN/dn", nBinsEta, hLowEta, hUpEta);
 619 hQtSelected = new TH1F("hQtSelected", "Daughter's transverse momentum in mother's frame for selected kinks (all cuts except dE/dx); q_{T} (GeV/c); dN/dq_{T}", 200, 0.0, 0.3);
 620 hKinkAngleSelected = new TH1F("hKinkAngleSelected", "Kink angle of selected kinks (all cuts except dE/dx); #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 621 hDCAkinkSelected = new TH1F("hDCAkinkSelected", "DCA; DCA(cm); Number of kinks", 100, 0.0, 2.0);
 622 hPmKinkAngSelected = new TH2F("hPmKinkAngSelected", "Mother's momentum magnitude vs kink angle for selected kinks (all cuts except dE/dx); p_{T} (GeV/c);  #theta (#circ); d^{2}N/dp_{T}d#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 623 hKinkPosXYSelected = new TH2F("hKinkPosXYSelected", "X-Y Position ; X (cm); Y (cm)", nBinsXY, hLowXY, hUpXY, nBinsXY, hLowXY, hUpXY);
 624 hKinkPosZRSelected = new TH2F("hKinkPosZRSelected", "Z vs radius of selected kinks (all cuts except dE/dx); Z (cm); R (cm)", nBinsZ, hLowZ, hUpZ, nBinsR, hLowR, hUpR);
 625 hKinkPosRSelected = new TH1F("hKinkPosRSelected", "Position radius of selected kinks (all cuts except dE/dx); R (cm); dN/dR", nBinsR, hLowR, hUpR);
 626 hPdgCodemdSelected = new TH2F("hPdgCodemdSelected", "Mother vs daughter pdg code for selected kinks (all cuts except dE/dx); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 627 hPmdSelected = new TH2F("hPmdSelected", "Mother vs daughter momentum magnitude for selected kinks (all cuts except dE/dx); Mother's P; Daughter's P", nBinsPt, hLowPt, hUpPt, nBinsPt, hLowPt, hUpPt);
 628 hMinvPimuSelected = new TH1F("hMinvPimuSelected", "Invariant mass for #pi^{#pm} decaying to #mu^{#pm} for selected kinks (all cuts except dE/dx); m (GeV/c^{2}); dN/dm", nBinsInvMass, hLowInvMass, hUpInvMass);
 629 hUIDKinkDauSelected = new TH1F("hUIDKinkDauSelected", "UID (method of production); UID (method of production); Number of particles", 100, 0.0, 20);
 630 hdEdxSelected = new TH2F("hdEdxSelected", "dE/dx vs mother's momentum for selected kinks (all cuts except dE/dx); p (GeV/c); dE/dx (GeV/cm)",   nBinsPt, hLowPt, hUpPt, 100, 0.0, 300.0);
 631 hPtSelectedFake = new TH1F("hPtSelectedFake", "Mother's transverse momentum for fake kinks (all cuts except dE/dx); p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 632 hEtaSelectedFake = new TH1F("hEtaSelectedFake", "Mother's pseudorapidity for fake kinks (all cuts except dE/dx); n; dN/dn", nBinsEta, hLowEta, hUpEta);
 633 hRapiditySelectedFake = new TH1F("hRapiditySelectedFake", "Mother's rapidity ffor fake kinks (all cuts except dE/dx); n; dN/dn", nBinsEta, hLowEta, hUpEta);
 634 hQtSelectedFake = new TH1F("hQtSelectedFake", "Daughter's transverse momentum in mother's frame ffor fake kinks (all cuts except dE/dx); q_{T} (GeV/c); dN/dq_{T}", 200, 0.0, 0.3);
 635 hKinkAngleSelectedFake = new TH1F("hKinkAngleSelectedFake", "Kink angle of fake kinks (all cuts except dE/dx); #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 636 hDCAkinkSelectedFake = new TH1F("hDCAkinkSelectedFake", "Kink DCA for fake kinks (all cuts except dE/dx); DCA; Number of kinks", 100, 0.0, 2.0);
 637 hPmKinkAngSelectedFake = new TH2F("hPmKinkAngSelectedFake", "Mother's momentum magnitude vs kink angle for fake kinks (all cuts except dE/dx);  #theta (#circ); d^{2}N/dp_{T}d#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 638 hKinkPosXYSelectedFake = new TH2F("hKinkPosXYSelectedFake", "X-Y Position for fake kinks (all cuts except dE/dx); X (cm); Y (cm)", nBinsXY, hLowXY, hUpXY, nBinsXY, hLowXY, hUpXY);
 639 hKinkPosZRSelectedFake = new TH2F("hKinkPosZRSelectedFake", "Z vs radius of fake kinks (all cuts except dE/dx); Z (cm); R (cm)", nBinsZ, hLowZ, hUpZ, nBinsR, hLowR, hUpR);
 640 hKinkPosRSelectedFake = new TH1F("hKinkPosRSelectedFake", "Position radius of fake kinks (all cuts except dE/dx); R (cm); dN/dR", nBinsR, hLowR, hUpR);
 641 hPmdSelectedFake = new TH2F("hPmdSelectedFake", "Mother vs daughter momentum magnitude for fake kinks (all cuts except dE/dx); Mother's P; Daughter's P", nBinsPt, hLowPt, hUpPt, nBinsPt, hLowPt, hUpPt);
 642 hMinvPimuSelectedFake = new TH1F("hMinvPimuSelectedFake", "Invariant mass for #pi^{#pm} decaying to #mu^{#pm} for selected kinks (all cuts) for fake kinks (all cuts except dE/dx); m (GeV/c^{2}); dN/dm", nBinsInvMass, hLowInvMass, hUpInvMass);
 643 hdEdxSelectedFake = new TH2F("hdEdxSelectedFake", "dE/dx vs mother's momentum for fake kinks (all cuts except dE/dx); p (GeV/c); dE/dx (GeV/cm)",   nBinsPt, hLowPt, hUpPt, 100, 0.0, 300.0);
 644 hPdgCodemddEdxRejected = new TH2F("hPdgCodemddEdxRejected", "Mother's vs daughter's pdg code for rejected kinks (all cuts); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 645 hPtPiSelected = new TH1F("hPtPiSelected", "Mother's transverse momentum for selected kinks (all cuts); p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 646 hEtaPiSelected = new TH1F("hEtaPiSelected", "Mother's pseudorapidity for selected kinks (all cuts); n; dN/dn", nBinsEta, hLowEta, hUpEta);
 647 hRapidityPiSelected = new TH1F("hRapidityPiSelected", "Mother's rapidity for selected kinks (all cuts); n; dN/dn", nBinsEta, hLowEta, hUpEta);
 648 hQtPiSelected = new TH1F("hQtPiSelected", "Daughter's transverse momentum in mother's frame for selected kinks (all cuts); q_{T} (GeV/c); dN/dq_{T}", 200, 0.0, 0.3);
 649 hKinkAnglePiSelected = new TH1F("hKinkAnglePiSelected", "Kink angle of selected kinks (all cuts); #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 650 hDCAkinkPiSelected = new TH1F("hDCAkinkPiSelected", "DCA; DCA(cm); Number of kinks", 100, 0.0, 2.0);
 651 hPmKinkAngPiSelected = new TH2F("hPmKinkAngPiSelected", "Mother's momentum magnitude vs kink angle for selected kinks (all cuts); p_{T} (GeV/c);  #theta (#circ); d^{2}N/dp_{T}d#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 652 hKinkPosRTPCclusters1PiSelected = new TH1F("hKinkPosRTPCclusters1PiSelected", " ;kinkposR; tpc clusters",100,0,1);
 653 hKinkPosRTPCclusters2PiSelected = new TH2F("hKinkPosRTPCclusters2PiSelected", " ;kinkposR; tpc clusters",nBinsR, hLowR, hUpR,100,0,200 );
 654 hKinkPosXYPiSelected = new TH2F("hKinkPosXYPiSelected", "X-Y Position ; X (cm); Y (cm)", nBinsXY, hLowXY, hUpXY, nBinsXY, hLowXY, hUpXY);
 655 hKinkPosZRPiSelected = new TH2F("hKinkPosZRPiSelected", "Z vs radius of selected kinks (all cuts); Z (cm); R (cm)", nBinsZ, hLowZ, hUpZ, nBinsR, hLowR, hUpR);
 656 hKinkPosRPiSelected = new TH1F("hKinkPosRPiSelected", "Position radius of selected kinks (all cuts); R (cm); dN/dR", nBinsR, hLowR, hUpR);
 657 hKinkPosZPiSelected = new TH1F("hKinkPosZPiSelected", "z position of selected kinks (all cuts); R (cm); dN/dR", nBinsZ, hLowZ, hUpZ);
 658 hPmPPiSelected = new TH2F("hPmPPiSelected", "Mother's momentum as calculated by kink and by track of selected kinks (all cuts); R (cm); dN/dR", nBinsPt, hLowPt, hUpPt, nBinsPt, hLowPt, hUpPt);
 659 hPdgCodemdPiSelected = new TH2F("hPdgCodemdPiSelected", "Mother vs daughter pdg code for selected kinks (all cuts); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 660 hPmdPiSelected = new TH2F("hPmdPiSelected", "Mother vs daughter momentum magnitude for selected kinks (all cuts); Mother's P; Daughter's P", nBinsPt, hLowPt, hUpPt, nBinsPt, hLowPt, hUpPt);
 661 hMinvPimuPiSelected = new TH1F("hMinvPimuPiSelected", "Invariant mass for #pi^{#pm} decaying to #mu^{#pm} for selected kinks (all cuts); m (GeV/c^{2}); dN/dm", nBinsInvMass, hLowInvMass, hUpInvMass);
 662 hUIDKinkDauPiSelected = new TH1F("hUIDKinkDauPiSelected", "UID (method of production); UID (method of production); Number of particles", 100, 0.0, 20);
 663 hdEdxPiSelected = new TH2F("hdEdxPiSelected", "dE/dx vs mother's momentum for selected kinks (all cuts); p (GeV/c); dE/dx (GeV/cm)",   nBinsPt, hLowPt, hUpPt, 100, 0.0, 300.0);
 664
 665 hPtPiSelectedPlus = new TH1F("hPtPiSelectedPlus", "Mother's transverse momentum for selected positive kinks (all cuts); p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 666 hEtaPiSelectedPlus = new TH1F("hEtaPiSelectedPlus", "Mother's pseudorapidity for selected positive kinks (all cuts); n; dN/dn", nBinsEta, hLowEta, hUpEta);
 667 hRapidityPiSelectedPlus = new TH1F("hRapidityPiSelectedPlus", "Mother's rapidity for selected positive kinks (all cuts); n; dN/dn", nBinsEta, hLowEta, hUpEta);
 668 hQtPiSelectedPlus = new TH1F("hQtPiSelectedPlus", "Daughter's transverse momentum in mother's frame for selected positive kinks (all cuts); q_{T} (GeV/c); dN/dq_{T}", 200, 0.0, 0.3);
 669 hKinkAnglePiSelectedPlus = new TH1F("hKinkAnglePiSelectedPlus", "Kink angle for selected positive kinks (all cuts); #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 670 hDCAkinkPiSelectedPlus = new TH1F("hDCAkinkPiSelectedPlus", "DCA for selected positive kinks (all cuts); DCA(cm); Number of kinks", 100, 0.0, 2.0);
 671 hPmKinkAngPiSelectedPlus = new TH2F("hPmKinkAngPiSelectedPlus", "Mother's momentum magnitude vs kink angle for selected positive kinks (all cuts); p_{T} (GeV/c);  #theta (#circ); d^{2}N/dp_{T}d#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 672 hKinkPosXYPiSelectedPlus = new TH2F("hKinkPosXYPiSelectedPlus", "X-Y Position of selected positive kinks (all cuts); X (cm); Y (cm)", nBinsXY, hLowXY, hUpXY, nBinsXY, hLowXY, hUpXY);
 673 hKinkPosZRPiSelectedPlus = new TH2F("hKinkPosZRPiSelectedPlus", "Z vs radius of selected positive kinks (all cuts); Z (cm); R (cm)", nBinsZ, hLowZ, hUpZ,nBinsR, hLowR, hUpR);
 674 hPdgCodemdPiSelectedPlus = new TH2F("hPdgCodemdPiSelectedPlus", "Mother vs daughter pdg code for selected positive kinks (all cuts); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 675 hPmdPiSelectedPlus = new TH2F("hPmdPiSelectedPlus", "Mother vs daughter momentum magnitude for selected positive kinks (all cuts); Mother's P; Daughter's P", nBinsPt, hLowPt, hUpPt, nBinsPt, hLowPt, hUpPt);
 676 hMinvPimuPiSelectedPlus = new TH1F("hMinvPimuPiSelectedPlus", "Invariant mass for #pi^{+} decaying to #mu^{+} for selected kinks (all cuts); m (GeV/c^{2}); dN/dm",nBinsInvMass, hLowInvMass, hUpInvMass);
 677 hUIDPiDaumuSelectedPlus = new TH1F("hUIDPiDaumuSelectedPlus", "UID (method of production) of selected positive kinks (all cuts); UID (method of production); Number of particles", 100, 0.0, 20);
 678 hdEdxPiSelectedPlus = new TH2F("hdEdxPiSelectedPlus", "dE/dx vs mother's momentum for selected positive kinks (all cuts); p (GeV/c); dE/dx (GeV/cm)",   nBinsPt, hLowPt, hUpPt, 100, 0.0, 300.0);
 679 hPtPrimPiKinksPlus = new TH1F("hPtPrimPiKinksPlus", "Mother's transverse momentum for selected real #pi^{+} primary kinks; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 680 hEtaPrimPiKinksPlus = new TH1F("hEtaPrimPiKinksPlus", "Mother's pseudorapidity for selected real #pi^{+} primary kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 681 hRapidityPrimPiKinksPlus = new TH1F("hRapidityPrimPiKinksPlus", "Mother's rapidity for selected real #pi^{+} primary kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 682 hPtSecondPiKinksPlus = new TH1F("hPtSecondPiKinksPlus", "Mother's transverse momentum for selected real #pi^{+} secondary kinks; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 683 hEtaSecondPiKinksPlus = new TH1F("hEtaSecondPiKinksPlus", "Mother's pseudorapidity for selected real #pi^{+} secondary kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 684 hRapiditySecondPiKinksPlus = new TH1F("hRapiditySecondPiKinksPlus", "Mother's rapidity for selected real #pi^{+} secondary kinks; y; dN/dy", nBinsEta, hLowEta, hUpEta);
 685 hPtNonPiKinksPlus = new TH1F("hPtNonPiKinksPlus", "Mother's transverse momentum for selected non #pi^{+} kinks; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 686 hEtaNonPiKinksPlus = new TH1F("hEtaNonPiKinksPlus", "Mother's pseudorapidity for selected non #pi^{+} kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 687 hRapidityNonPiKinksPlus = new TH1F("hRapidityNonPiKinksPlus", "Mother's rapidity for selected non #pi^{+} kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 688 hPdgCodemdNonPiKinksPlus = new TH2F("hPdgCodemdNonPiKinksPlus", "Mother vs daughter pdg code for selected non #pi^{+} kinks; Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 689
 690 hPtPiSelectedMinus = new TH1F("hPtPiSelectedMinus", "Mother's transverse momentum for selected negative kinks (all cuts); p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 691 hEtaPiSelectedMinus = new TH1F("hEtaPiSelectedMinus", "Mother's pseudorapidity for selected negative kinks (all cuts); n; dN/dn", nBinsEta, hLowEta, hUpEta);
 692 hRapidityPiSelectedMinus = new TH1F("hRapidityPiSelectedMinus", "Mother's rapidity for selected negative kinks (all cuts); n; dN/dn", nBinsEta, hLowEta, hUpEta);
 693 hQtPiSelectedMinus = new TH1F("hQtPiSelectedMinus", "Daughter's transverse momentum in mother's frame for selected negative kinks (all cuts); q_{T} (GeV/c); dN/dq_{T}", 200, 0.0, 0.3);
 694 hKinkAnglePiSelectedMinus = new TH1F("hKinkAnglePiSelectedMinus", "Kink angle for selected negative kinks (all cuts); #theta (#circ); dN/d#theta", nBinsAngle, hLowAngle, hUpAngle);
 695 hDCAkinkPiSelectedMinus = new TH1F("hDCAkinkPiSelectedMinus", "DCA for selected negative kinks (all cuts); DCA(cm); Number of kinks", 100, 0.0, 2.0);
 696 hPmKinkAngPiSelectedMinus = new TH2F("hPmKinkAngPiSelectedMinus", "Mother's momentum magnitude vs kink angle for selected negative kinks (all cuts); p_{T} (GeV/c);  #theta (#circ); d^{2}N/dp_{T}d#theta", nBinsPt, hLowPt, hUpPt, nBinsAngle, hLowAngle, hUpAngle);
 697 hKinkPosXYPiSelectedMinus = new TH2F("hKinkPosXYPiSelectedMinus", "X-Y Position for selected negative kinks (all cuts); X (cm); Y (cm)", nBinsXY, hLowXY, hUpXY, nBinsXY, hLowXY, hUpXY);
 698 hKinkPosZRPiSelectedMinus = new TH2F("hKinkPosZRPiSelectedMinus", "Z vs radius for selected negative kinks (all cuts); Z (cm); R (cm)",nBinsZ, hLowZ, hUpZ,nBinsR, hLowR, hUpR);
 699 hPdgCodemdPiSelectedMinus = new TH2F("hPdgCodemdPiSelectedMinus", "Mother vs daughter pdg code for selected negative kinks (all cuts); Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 700 hPmdPiSelectedMinus = new TH2F("hPmdPiSelectedMinus", "Mother vs daughter momentum magnitude for selected negative kinks (all cuts); Mother's P; Daughter's P", nBinsPt, hLowPt, hUpPt, nBinsPt, hLowPt, hUpPt);
 701 hMinvPimuPiSelectedMinus = new TH1F("hMinvPimuPiSelectedMinus", "Invariant mass for #pi^{-} decaying to #mu^{-} for selected kinks (all cuts); m (GeV/c^{2}); dN/dm", nBinsInvMass, hLowInvMass, hUpInvMass);
 702 hUIDPiDaumuSelectedMinus = new TH1F("hUIDPiDaumuSelectedMinus", "UID (method of production) of selected negative kinks (all cuts); UID (method of production); Number of particles", 100, 0.0, 20);
 703 hdEdxPiSelectedMinus = new TH2F("hdEdxPiSelectedMinus", "dE/dx vs mother's momentum for selected negative kinks (all cuts); p (GeV/c); dE/dx (GeV/cm)",   nBinsPt, hLowPt, hUpPt, 100, 0.0, 300.0);
 704 hPtPrimPiKinksMinus = new TH1F("hPtPrimPiKinksMinus", "Mother's transverse momentum for selected real #pi^{-} primary kinks; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 705 hEtaPrimPiKinksMinus = new TH1F("hEtaPrimPiKinksMinus", "Mother's pseudorapidity for selected real #pi^{-} primary kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 706 hRapidityPrimPiKinksMinus = new TH1F("hRapidityPrimPiKinksMinus", "Mother's rapidity for selected real #pi^{-} primary kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 707 hPtSecondPiKinksMinus = new TH1F("hPtSecondPiKinksMinus", "Mother's transverse momentum for selected real #pi^{-} secondary kinks; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 708 hEtaSecondPiKinksMinus = new TH1F("hEtaSecondPiKinksMinus", "Mother's pseudorapidity for selected real #pi^{-} secondary kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 709 hRapiditySecondPiKinksMinus = new TH1F("hRapiditySecondPiKinksMinus", "Mother's rapidity for selected real #pi^{-} secondary kinks; y; dN/dy", nBinsEta, hLowEta, hUpEta);
 710 hPtNonPiKinksMinus = new TH1F("hPtNonPiKinksMinus", "Mother's transverse momentum for selected non #pi^{-} kinks; p_{T} (GeV/c); dN/dp_{T}",   nBinsPt, hLowPt, hUpPt);
 711 hEtaNonPiKinksMinus = new TH1F("hEtaNonPiKinksMinus", "Mother's pseudorapidity for selected non #pi^{-} kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 712 hRapidityNonPiKinksMinus = new TH1F("hRapidityNonPiKinksMinus", "Mother's rapidity for selected non #pi^{-} kinks; n; dN/dn", nBinsEta, hLowEta, hUpEta);
 713 hPdgCodemdNonPiKinksMinus = new TH2F("hPdgCodemdNonPiKinksMinus", "Mother vs daughter pdg code for selected non #pi^{-} kinks; Mother's pdg code; Daughter's pdg code", nBinsPdg, hLowPdg, hUpPdg, nBinsPdg, hLowPdg, hUpPdg);
 714
 715
 716
 717 fListOfHistos->Add(hMCMult);
 718 fListOfHistos->Add(hMCMultPrim);
 719 fListOfHistos->Add(hMCPtAll);
 720 fListOfHistos->Add(hMCEtaAll);
 721 fListOfHistos->Add(hMCPtPrim);
 722 fListOfHistos->Add(hMCEtaPrim);
 723 fListOfHistos->Add(hMCPt);
 724 fListOfHistos->Add(hMCEta);
 725 fListOfHistos->Add(hMCPdg);
 726 fListOfHistos->Add(hMCMultPiPlus);
 727 fListOfHistos->Add(hMCPtPiPlus);
 728 fListOfHistos->Add(hMCEtaPiPlus);
 729 fListOfHistos->Add(hMCRapidityPiPlus);
 730 fListOfHistos->Add(hMCNDaughtersPlus);
 731 fListOfHistos->Add(hMCRadPiDauPlus);
 732 fListOfHistos->Add(hMCKinkPosZPlus);
 733 fListOfHistos->Add(hMCUIDPiDauPlus);
 734 fListOfHistos->Add(hMCPdgPiNonDecayedPlus);
 735 fListOfHistos->Add(hMCPdgPiDauPlus);
 736 fListOfHistos->Add(hMCQtPlus);
 737 fListOfHistos->Add(hMCKinkAnglePlus);
 738 fListOfHistos->Add(hMCPKinkAngPlus);
 739 fListOfHistos->Add(hMCPdgCodemdPlus);
 740 fListOfHistos->Add(hMCPtmdPlus);
 741 fListOfHistos->Add(hMCPtPimuonPlus);
 742 fListOfHistos->Add(hMCEtaPimuonPlus);
 743 fListOfHistos->Add(hMCRapidityPimuonPlus);
 744 fListOfHistos->Add(hMCQtPimuonPlus);
 745 fListOfHistos->Add(hMCPKinkAngPimuonPlus);
 746 fListOfHistos->Add(hMCPtPiotherPlus);
 747 fListOfHistos->Add(hMCEtaPiotherPlus);
 748 fListOfHistos->Add(hMCRapidityPiotherPlus);
 749 fListOfHistos->Add(hMCQtPiotherPlus);
 750 fListOfHistos->Add(hMCPKinkAngPiotherPlus);
 751 fListOfHistos->Add(hMCMultPiMinus);
 752 fListOfHistos->Add(hMCPtPiMinus);
 753 fListOfHistos->Add(hMCEtaPiMinus);
 754 fListOfHistos->Add(hMCRapidityPiMinus);
 755 fListOfHistos->Add(hMCNDaughtersMinus);
 756 fListOfHistos->Add(hMCRadPiDauMinus);
 757 fListOfHistos->Add(hMCKinkPosZMinus);
 758 fListOfHistos->Add(hMCUIDPiDauMinus);
 759 fListOfHistos->Add(hMCPdgPiNonDecayedMinus);
 760 fListOfHistos->Add(hMCPdgPiDauMinus);
 761 fListOfHistos->Add(hMCQtMinus);
 762 fListOfHistos->Add(hMCKinkAngleMinus);
 763 fListOfHistos->Add(hMCPKinkAngMinus);
 764 fListOfHistos->Add(hMCPdgCodemdMinus);
 765 fListOfHistos->Add(hMCPtmdMinus);
 766 fListOfHistos->Add(hMCPtPimuonMinus);
 767 fListOfHistos->Add(hMCEtaPimuonMinus);
 768 fListOfHistos->Add(hMCRapidityPimuonMinus);
 769 fListOfHistos->Add(hMCQtPimuonMinus);
 770 fListOfHistos->Add(hMCPKinkAngPimuonMinus);
 771 fListOfHistos->Add(hMCPtPiotherMinus);
 772 fListOfHistos->Add(hMCEtaPiotherMinus);
 773 fListOfHistos->Add(hMCRapidityPiotherMinus);
 774 fListOfHistos->Add(hMCQtPiotherMinus);
 775 fListOfHistos->Add(hMCPKinkAngPiotherMinus);
 776
 777 fListOfHistos->Add(hMult);
 778 fListOfHistos->Add(hAcceptedMult);
 779 fListOfHistos->Add(hMultPS);
 780 fListOfHistos->Add(hvtx);
 781 fListOfHistos->Add(hvtxy);
 782 fListOfHistos->Add(hvtyz);
 783 fListOfHistos->Add(hvtxz);
 784 fListOfHistos->Add(hMultPSV);
 785 fListOfHistos->Add(hPtAll);
 786 fListOfHistos->Add(hEtaAll);
 787 fListOfHistos->Add(hTrackPos);
 788 fListOfHistos->Add(hTrackPosxy);
 789 fListOfHistos->Add(hTrackPosyz);
 790 fListOfHistos->Add(hTrackPosxz);
 791 //fListOfHistos->Add(hTPCchi2clusters);
 792 //fListOfHistos->Add(hdcaToVertexXY);
 793 //fListOfHistos->Add(hdcaToVertexZ);
 794 fListOfHistos->Add(hMultPrim);
 795 fListOfHistos->Add(hPtPrim);
 796 fListOfHistos->Add(hEtaPrim);
 797 fListOfHistos->Add(hPrimTrackPos);
 798 fListOfHistos->Add(hPrimTrackPosxy);
 799 fListOfHistos->Add(hPrimTrackPosyz);
 800 fListOfHistos->Add(hPrimTrackPosxz);
 801 fListOfHistos->Add(hPt);
 802 fListOfHistos->Add(hEta);
 803 //fListOfHistos->Add(hRapidity);
 804 fListOfHistos->Add(hPtKink);
 805 fListOfHistos->Add(hEtaKink);
 806 fListOfHistos->Add(hRapidityKink);
 807 fListOfHistos->Add(hPmP);
 808 fListOfHistos->Add(hKinkPosRTPCclusters1);
 809 fListOfHistos->Add(hKinkPosRTPCclusters2);
 810 fListOfHistos->Add(hQt);
 811 fListOfHistos->Add(hKinkAngle);
 812 fListOfHistos->Add(hDCAkink);
 813 fListOfHistos->Add(hPmKinkAng);
 814 fListOfHistos->Add(hKinkPosXY);
 815 fListOfHistos->Add(hKinkPosZY);
 816 fListOfHistos->Add(hKinkPosZR);
 817 fListOfHistos->Add(hKinkPosR);
 818 fListOfHistos->Add(hKinkPosZ);
 819 fListOfHistos->Add(hKinkPosZMCKinkPosZ);
 820 fListOfHistos->Add(hPdgCodemd);
 821 fListOfHistos->Add(hPmd);
 822 fListOfHistos->Add(hMinvPimu);
 823 fListOfHistos->Add(hUIDKinkDau);
 824 fListOfHistos->Add(hdEdx);
 825 fListOfHistos->Add(hPtKinkFake);
 826 fListOfHistos->Add(hEtaKinkFake);
 827 fListOfHistos->Add(hRapidityKinkFake);
 828 fListOfHistos->Add(hPmPFake);
 829 fListOfHistos->Add(hKinkPosRTPCclusters1Fake);
 830 fListOfHistos->Add(hKinkPosRTPCclusters2Fake);
 831 fListOfHistos->Add(hQtFake);
 832 fListOfHistos->Add(hKinkAngleFake);
 833 fListOfHistos->Add(hDCAkinkFake);
 834 fListOfHistos->Add(hPmKinkAngFake);
 835 fListOfHistos->Add(hKinkPosXYFake);
 836 fListOfHistos->Add(hKinkPosZYFake);
 837 fListOfHistos->Add(hKinkPosZRFake);
 838 fListOfHistos->Add(hKinkPosRFake);
 839 fListOfHistos->Add(hKinkPosZFake);
 840 fListOfHistos->Add(hKinkPosZMCKinkPosZFake);
 841 fListOfHistos->Add(hPdgCodemdFake);
 842 fListOfHistos->Add(hPmdFake);
 843 fListOfHistos->Add(hMinvPimuFake);
 844 fListOfHistos->Add(hUIDKinkDauFake);
 845 fListOfHistos->Add(hdEdxFake);
 846 fListOfHistos->Add(hPtPosRSelected);
 847 fListOfHistos->Add(hPdgCodemdZRejected);
 848 fListOfHistos->Add(hPtZSelected);
 849 fListOfHistos->Add(hPdgCodemdAngRejected);
 850 fListOfHistos->Add(hPtAngSelected);
 851 fListOfHistos->Add(hPdgCodemdPmRejected);
 852 fListOfHistos->Add(hPtPmSelected);
 853 fListOfHistos->Add(hPdgCodemdQtLowRejected);
 854 fListOfHistos->Add(hPtGoodKink);
 855 fListOfHistos->Add(hEtaGoodKink);
 856 fListOfHistos->Add(hRapidityGoodKink);
 857 fListOfHistos->Add(hQtGoodKink);
 858 fListOfHistos->Add(hPmGoodKinkAng);
 859 fListOfHistos->Add(hPdgCodemdGoodKink);
 860 fListOfHistos->Add(hPmdGoodKink);
 861 fListOfHistos->Add(hUIDGoodKinkDau);
 862 fListOfHistos->Add(hdEdxGoodKink);
 863 fListOfHistos->Add(hUIDPiDauPlus);
 864 fListOfHistos->Add(hMultPiPlus);
 865 fListOfHistos->Add(hPtPiPlus);
 866 fListOfHistos->Add(hEtaPiPlus);
 867 fListOfHistos->Add(hRapidityPiPlus);
 868 fListOfHistos->Add(hQtPiPlus);
 869 fListOfHistos->Add(hKinkAnglePiPlus);
 870 fListOfHistos->Add(hPmKinkAngPiPlus);
 871 fListOfHistos->Add(hKinkPosXYPiPlus);
 872 fListOfHistos->Add(hKinkPosZRPiPlus);
 873 fListOfHistos->Add(hKinkPosRPiPlus);
 874 fListOfHistos->Add(hDCAkinkPiPlus);
 875 fListOfHistos->Add(hPdgCodemdPiPlus);
 876 fListOfHistos->Add(hPmdPiPlus);
 877 fListOfHistos->Add(hdEdxPiPlus);
 878 fListOfHistos->Add(hQtPimuPlus);
 879 fListOfHistos->Add(hPmKinkAngPimuPlus);
 880 fListOfHistos->Add(hQtPiotherPlus);
 881 fListOfHistos->Add(hPmKinkAngPiotherPlus);
 882 fListOfHistos->Add(hPdgCodemdPiotherPlus);
 883 fListOfHistos->Add(hUIDPiDauMinus);
 884 fListOfHistos->Add(hMultPiMinus);
 885 fListOfHistos->Add(hPtPiMinus);
 886 fListOfHistos->Add(hEtaPiMinus);
 887 fListOfHistos->Add(hRapidityPiMinus);
 888 fListOfHistos->Add(hQtPiMinus);
 889 fListOfHistos->Add(hKinkAnglePiMinus);
 890 fListOfHistos->Add(hPmKinkAngPiMinus);
 891 fListOfHistos->Add(hKinkPosXYPiMinus);
 892 fListOfHistos->Add(hKinkPosZRPiMinus);
 893 fListOfHistos->Add(hKinkPosRPiMinus);
 894 fListOfHistos->Add(hDCAkinkPiMinus);
 895 fListOfHistos->Add(hPdgCodemdPiMinus);
 896 fListOfHistos->Add(hPmdPiMinus);
 897 fListOfHistos->Add(hdEdxPiMinus);
 898 fListOfHistos->Add(hQtPimuMinus);
 899 fListOfHistos->Add(hPmKinkAngPimuMinus);
 900 fListOfHistos->Add(hQtPiotherMinus);
 901 fListOfHistos->Add(hPmKinkAngPiotherMinus);
 902 fListOfHistos->Add(hPdgCodemdPiotherMinus);
 903 fListOfHistos->Add(hPdgCodemdQtRejected);
 904 fListOfHistos->Add(hPtQtSelected);
 905 fListOfHistos->Add(hPdgCodemdMaxAngRejected);
 906 fListOfHistos->Add(hPtMaxAngSelected);
 907 fListOfHistos->Add(hPdgCodemdRTPCclustersRejected);
 908 fListOfHistos->Add(hPtRTPCclustersSelected);
 909 fListOfHistos->Add(hRTPCclustersRTPCclustersSelected);
 910 fListOfHistos->Add(hPdgCodemdMinvRejected);
 911 fListOfHistos->Add(hPtSelected);
 912 fListOfHistos->Add(hEtaSelected);
 913 fListOfHistos->Add(hRapiditySelected);
 914 fListOfHistos->Add(hQtSelected);
 915 fListOfHistos->Add(hKinkAngleSelected);
 916 fListOfHistos->Add(hDCAkinkSelected);
 917 fListOfHistos->Add(hPmKinkAngSelected);
 918 fListOfHistos->Add(hKinkPosXYSelected);
 919 fListOfHistos->Add(hKinkPosZRSelected);
 920 fListOfHistos->Add(hKinkPosRSelected);
 921 fListOfHistos->Add(hPdgCodemdSelected);
 922 fListOfHistos->Add(hPmdSelected);
 923 fListOfHistos->Add(hMinvPimuSelected);
 924 fListOfHistos->Add(hUIDKinkDauSelected);
 925 fListOfHistos->Add(hdEdxSelected);
 926 fListOfHistos->Add(hPtSelectedFake);
 927 fListOfHistos->Add(hEtaSelectedFake);
 928 fListOfHistos->Add(hRapiditySelectedFake);
 929 fListOfHistos->Add(hQtSelectedFake);
 930 fListOfHistos->Add(hKinkAngleSelectedFake);
 931 fListOfHistos->Add(hDCAkinkSelectedFake);
 932 fListOfHistos->Add(hPmKinkAngSelectedFake);
 933 fListOfHistos->Add(hKinkPosXYSelectedFake);
 934 fListOfHistos->Add(hKinkPosZRSelectedFake);
 935 fListOfHistos->Add(hKinkPosRSelectedFake);
 936 fListOfHistos->Add(hPmdSelectedFake);
 937 fListOfHistos->Add(hMinvPimuSelectedFake);
 938 fListOfHistos->Add(hdEdxSelectedFake);
 939 fListOfHistos->Add(hPdgCodemddEdxRejected);
 940 fListOfHistos->Add(hPtPiSelected);
 941 fListOfHistos->Add(hEtaPiSelected);
 942 fListOfHistos->Add(hRapidityPiSelected);
 943 fListOfHistos->Add(hQtPiSelected);
 944 fListOfHistos->Add(hKinkAnglePiSelected);
 945 fListOfHistos->Add(hDCAkinkPiSelected);
 946 fListOfHistos->Add(hPmKinkAngPiSelected);
 947 fListOfHistos->Add(hKinkPosRTPCclusters1PiSelected);
 948 fListOfHistos->Add(hKinkPosRTPCclusters2PiSelected);
 949 fListOfHistos->Add(hKinkPosXYPiSelected);
 950 fListOfHistos->Add(hKinkPosZRPiSelected);
 951 fListOfHistos->Add(hKinkPosRPiSelected);
 952 fListOfHistos->Add(hKinkPosZPiSelected);
 953 fListOfHistos->Add(hPmPPiSelected);
 954 fListOfHistos->Add(hPdgCodemdPiSelected);
 955 fListOfHistos->Add(hPmdPiSelected);
 956 fListOfHistos->Add(hMinvPimuPiSelected);
 957 fListOfHistos->Add(hUIDKinkDauPiSelected);
 958 fListOfHistos->Add(hdEdxPiSelected);
 959
 960 fListOfHistos->Add(hPtPiSelectedPlus);
 961 fListOfHistos->Add(hEtaPiSelectedPlus);
 962 fListOfHistos->Add(hRapidityPiSelectedPlus);
 963 fListOfHistos->Add(hQtPiSelectedPlus);
 964 fListOfHistos->Add(hKinkAnglePiSelectedPlus);
 965 fListOfHistos->Add(hDCAkinkPiSelectedPlus);
 966 fListOfHistos->Add(hPmKinkAngPiSelectedPlus);
 967 fListOfHistos->Add(hKinkPosXYPiSelectedPlus);
 968 fListOfHistos->Add(hKinkPosZRPiSelectedPlus);
 969 fListOfHistos->Add(hPdgCodemdPiSelectedPlus);
 970 fListOfHistos->Add(hPmdPiSelectedPlus);
 971 fListOfHistos->Add(hMinvPimuPiSelectedPlus);
 972 fListOfHistos->Add(hUIDPiDaumuSelectedPlus);
 973 fListOfHistos->Add(hdEdxPiSelectedPlus);
 974 fListOfHistos->Add(hPtPrimPiKinksPlus);
 975 fListOfHistos->Add(hEtaPrimPiKinksPlus);
 976 fListOfHistos->Add(hRapidityPrimPiKinksPlus);
 977 fListOfHistos->Add(hPtSecondPiKinksPlus);
 978 fListOfHistos->Add(hEtaSecondPiKinksPlus);
 979 fListOfHistos->Add(hRapiditySecondPiKinksPlus);
 980 fListOfHistos->Add(hPtNonPiKinksPlus);
 981 fListOfHistos->Add(hEtaNonPiKinksPlus);
 982 fListOfHistos->Add(hRapidityNonPiKinksPlus);
 983 fListOfHistos->Add(hPdgCodemdNonPiKinksPlus);
 984
 985 fListOfHistos->Add(hPtPiSelectedMinus);
 986 fListOfHistos->Add(hEtaPiSelectedMinus);
 987 fListOfHistos->Add(hRapidityPiSelectedMinus);
 988 fListOfHistos->Add(hQtPiSelectedMinus);
 989 fListOfHistos->Add(hKinkAnglePiSelectedMinus);
 990 fListOfHistos->Add(hDCAkinkPiSelectedMinus);
 991 fListOfHistos->Add(hPmKinkAngPiSelectedMinus);
 992 fListOfHistos->Add(hKinkPosXYPiSelectedMinus);
 993 fListOfHistos->Add(hKinkPosZRPiSelectedMinus);
 994 fListOfHistos->Add(hPdgCodemdPiSelectedMinus);
 995 fListOfHistos->Add(hPmdPiSelectedMinus);
 996 fListOfHistos->Add(hMinvPimuPiSelectedMinus);
 997 fListOfHistos->Add(hUIDPiDaumuSelectedMinus);
 998 fListOfHistos->Add(hdEdxPiSelectedMinus);
 999 fListOfHistos->Add(hPtPrimPiKinksMinus);
1000 fListOfHistos->Add(hEtaPrimPiKinksMinus);
1001 fListOfHistos->Add(hRapidityPrimPiKinksMinus);
1002 fListOfHistos->Add(hPtSecondPiKinksMinus);
1003 fListOfHistos->Add(hEtaSecondPiKinksMinus);
1004 fListOfHistos->Add(hRapiditySecondPiKinksMinus);
1005 fListOfHistos->Add(hPtNonPiKinksMinus);
1006 fListOfHistos->Add(hEtaNonPiKinksMinus);
1007 fListOfHistos->Add(hRapidityNonPiKinksMinus);
1008 fListOfHistos->Add(hPdgCodemdNonPiKinksMinus);
1009
1010 fListOfHistos->SetOwner(kTRUE);
1011 PostData(1, fListOfHistos);
1012 }
1013
1014 //________________________________________________________________________
1015 void AliAnalysisPionKinksMCESD::UserExec(Option_t *) {
1016         AliVEvent *event = InputEvent();
1017         if (!event) {
1018                 Printf("ERROR: Could not retrieve event");
1019                 return;
1020         }
1021
1022         AliMCEvent* mcEvent = MCEvent();
1023         if (!mcEvent) {
1024                 Printf("ERROR: Could not retrieve MC event");
1025                 return;
1026         }
1027         AliESDEvent* esdEvent = dynamic_cast<AliESDEvent*>(event);
1028         if (!esdEvent) {
1029                 Printf("ERROR: Could not retrieve esd");
1030                 return;
1031         }
1032
1033
1034 //------------------------------ MC Multiplicity unbiased------------------------//
1035         AliStack* MCstack = mcEvent->Stack(); //stack of MC events
1036         Int_t MCNTracks = mcEvent->GetNumberOfTracks();
1037 //      Int_t MCNTracks2 = MCstack->GetNtrack();
1038         hMCMult->Fill(MCNTracks);
1039         Int_t MCNPrimTracks = MCstack->GetNprimary();
1040         //Printf("There are %d number of primary tracks in this Monte Carlo event", MCNPrimTracks);
1041         hMCMultPrim->Fill(MCNPrimTracks);
1042
1043 //------------------------ Data Multiplicity unbiased --------------------------//
1044         Int_t NTracks = esdEvent->GetNumberOfTracks(); //number of ESD tracks
1045         hMult->Fill(NTracks);
1046
1047 //----------------------------- Accepted Multiplicity -----------------------------//
1048         Float_t NAcceptedTracks = fTrackCuts->CountAcceptedTracks(esdEvent);
1049         if(fLowMulcut>-1) {
1050                 if(NAcceptedTracks<fLowMulcut) return;
1051         }
1052         if(fUpMulcut>-1) {
1053                 if(NAcceptedTracks>fUpMulcut) return;
1054         }
1055         hAcceptedMult->Fill(NAcceptedTracks); //to check if the multiplicity limits are ok
1056
1057
1058 //-------------------------------- MC data Analysis -----------------------------//
1059         for (Int_t iMC = 0; iMC < MCNTracks; iMC++) { //loop on all accepted MC tracks
1060                 TParticle* MCparticle = MCstack->Particle(iMC);
1061                 if (!MCparticle) {
1062                         Printf("ERROR: Could not receive MC particle %d", iMC);
1063                         continue;
1064                 }
1065
1066                 Double_t MCPt = MCparticle->Pt();
1067                 Double_t MCP = MCparticle->P();
1068                 Double_t MCEta = MCparticle->Eta();
1069                 hMCPtAll->Fill(MCPt);
1070                 hMCEtaAll->Fill(MCEta);
1071
1072                 if (!MCstack->IsPhysicalPrimary(iMC)) continue;
1073                 hMCPtPrim->Fill(MCPt);
1074                 hMCEtaPrim->Fill(MCEta);
1075
1076                 if (MCPt<cLowPt) continue;
1077
1078                 hMCPt->Fill(MCPt);
1079                 hMCEta->Fill(MCEta);
1080
1081                 Int_t MCPdg = MCparticle->GetPdgCode();
1082                 hMCPdg->Fill(MCPdg);
1083
1084                 if (MCPdg==cPdgPion) { //positive pion selection
1085                 Double_t MCRapidity = GetMCRapidity(MCparticle);
1086                 if (TMath::Abs(MCRapidity)>cRapidityLim) continue;
1087                         hMCMultPiPlus->Fill(iMC); //MC pion multiplicity
1088                         hMCPtPiPlus->Fill(MCPt);
1089                         hMCEtaPiPlus->Fill(MCEta);
1090                         hMCRapidityPiPlus->Fill(MCRapidity);
1091
1092                         Int_t MCNDaughtersPlus = MCparticle->GetNDaughters();
1093                         hMCNDaughtersPlus->Fill(MCNDaughtersPlus);
1094                         Int_t FirstDPlus = MCparticle->GetFirstDaughter(); //first daughter's label
1095                         Int_t LastDPlus = MCparticle->GetLastDaughter(); //last daughter's label
1096                         if ((FirstDPlus>MCNTracks)||(LastDPlus>MCNTracks)) continue;
1097
1098
1099                         for (Int_t iMCdPlus=FirstDPlus; iMCdPlus<=LastDPlus; iMCdPlus++) { //loop on pion daughters
1100                                 if (iMCdPlus<0) continue; //debug
1101                                 TParticle* MCdaughterPlus = MCstack->Particle(iMCdPlus);
1102
1103                                 Double_t RdPlus = MCdaughterPlus->R(); //position radius of daughter's generation vertex
1104                                 hMCRadPiDauPlus->Fill(RdPlus);
1105
1106                                 Double_t MCKinkPosZPlus = MCdaughterPlus->Vz();
1107                                 hMCKinkPosZPlus->Fill(MCKinkPosZPlus);
1108
1109                                 if ((TMath::Abs(MCKinkPosZPlus)<cLowZ) || (TMath::Abs(MCKinkPosZPlus)>cUpZ)) continue;
1110
1111                                 if ((RdPlus<cLowR) || (RdPlus>cUpR)) continue; //selection of daughters that where generated in TPC
1112
1113                                 Int_t dcodePlus = MCdaughterPlus->GetPdgCode();
1114
1115                                 UInt_t MCDauUIDPlus = MCdaughterPlus->GetUniqueID(); //daughter's unique id (= method of production)
1116                                 hMCUIDPiDauPlus->Fill(MCDauUIDPlus);
1117
1118                                 if (MCDauUIDPlus!=4) { //selection of decayed daughters
1119                                         hMCPdgPiNonDecayedPlus->Fill(dcodePlus);
1120                                         continue;
1121                                 }
1122
1123                                 hMCPdgPiDauPlus->Fill(TMath::Abs(dcodePlus));
1124
1125                                 Double_t MCPtdPlus = MCdaughterPlus->Pt(); //daughter's transverse momentum in lab frame
1126                                 Double_t MCQtPlus = MCPQt(mcEvent, iMC, MCdaughterPlus); //daughter's transverse momentum in mother's frame (Qt)
1127                                 Double_t MCKinkAnglePlus2 = fuMCKinkAngle(mcEvent, iMC, MCdaughterPlus, kTRUE); //kink angle in degrees
1128                                 Double_t MCMaxKinkAngPimuPlus=fMaxKinkAngPimu->Eval(MCP,0.,0.,0.); //maximum decay angle in lab for pion decaying to muon
1129
1130                                 // if  (MCKinkAnglePlus2>(MCMaxKinkAngPimuPlus*1.1)) continue;
1131                                 if  (MCKinkAnglePlus2>(MCMaxKinkAngPimuPlus)) continue;
1132
1133                                 hMCQtPlus->Fill(MCQtPlus);
1134                                 hMCKinkAnglePlus->Fill(MCKinkAnglePlus2);
1135                                 hMCPKinkAngPlus->Fill(MCP,MCKinkAnglePlus2);
1136                                 hMCPdgCodemdPlus->Fill(TMath::Abs(MCPdg),TMath::Abs(dcodePlus));
1137                                 hMCPtmdPlus->Fill(MCPt, MCPtdPlus);
1138
1139                                 if ((MCQtPlus<cLowQt) || (MCQtPlus>cUpQt) || (MCKinkAnglePlus2<cLowKinkAngle)) continue;
1140
1141                                 if (dcodePlus==-cPdgMuon) { //muon daughters selection
1142                                         hMCPtPimuonPlus->Fill(MCPt);
1143                                         hMCEtaPimuonPlus->Fill(MCEta);
1144                                         hMCRapidityPimuonPlus->Fill(MCRapidity);
1145                                         hMCQtPimuonPlus->Fill(MCQtPlus);
1146                                         hMCPKinkAngPimuonPlus->Fill(MCKinkAnglePlus2);
1147                                 } else {hMCPtPiotherPlus->Fill(MCPt);
1148                                         hMCEtaPiotherPlus->Fill(MCEta);
1149                                         hMCRapidityPiotherPlus->Fill(MCRapidity);
1150                                         hMCQtPiotherPlus->Fill(MCQtPlus);
1151                                         hMCPKinkAngPiotherPlus->Fill(MCKinkAnglePlus2);
1152                                 }
1153                         } //end of pion daughters' loop
1154                 } else if (MCPdg==-cPdgPion) { //negative pion selection
1155                 Double_t MCRapidity = GetMCRapidity(MCparticle);
1156                 if (TMath::Abs(MCRapidity)>cRapidityLim) continue;
1157                         hMCMultPiMinus->Fill(iMC); //MC pion multiplicity
1158                         hMCPtPiMinus->Fill(MCPt);
1159                         hMCEtaPiMinus->Fill(MCEta);
1160                         hMCRapidityPiMinus->Fill(MCRapidity);
1161
1162                         Int_t MCNDaughtersMinus = MCparticle->GetNDaughters();
1163                         hMCNDaughtersMinus->Fill(MCNDaughtersMinus);
1164                         Int_t FirstDMinus = MCparticle->GetFirstDaughter(); //first daughter's label
1165                         Int_t LastDMinus = MCparticle->GetLastDaughter(); //last daughter's label
1166                         if ((FirstDMinus>MCNTracks)||(LastDMinus>MCNTracks)) continue;
1167
1168                         for (Int_t iMCdMinus=FirstDMinus; iMCdMinus<=LastDMinus; iMCdMinus++) { //loop on pion daughters
1169                                 if (iMCdMinus<0) continue; //debug
1170                                 TParticle* MCdaughterMinus = MCstack->Particle(iMCdMinus);
1171
1172                                 Double_t RdMinus = MCdaughterMinus->R(); //position radius of daughter's generation vertex
1173                                 hMCRadPiDauMinus->Fill(RdMinus);
1174
1175                                 Double_t MCKinkPosZMinus = MCdaughterMinus->Vz();
1176                                 hMCKinkPosZMinus->Fill(MCKinkPosZMinus);
1177
1178                                 if ((TMath::Abs(MCKinkPosZMinus)<cLowZ) || (TMath::Abs(MCKinkPosZMinus)>cUpZ)) continue;
1179
1180                                 if ((RdMinus<cLowR) || (RdMinus>cUpR)) continue; //selection of daughters that where generated in TPC
1181
1182                                 Int_t dcodeMinus = MCdaughterMinus->GetPdgCode();
1183
1184                                 UInt_t MCDauUIDMinus = MCdaughterMinus->GetUniqueID(); //daughter's unique id (= method of production)
1185                                 hMCUIDPiDauMinus->Fill(MCDauUIDMinus);
1186
1187                                 if (MCDauUIDMinus!=4) { //selection of decayed daughters
1188                                         hMCPdgPiNonDecayedMinus->Fill(dcodeMinus);
1189                                         continue;
1190                                 }
1191
1192                                 hMCPdgPiDauMinus->Fill(TMath::Abs(dcodeMinus));
1193
1194
1195                                 Double_t MCPtdMinus = MCdaughterMinus->Pt(); //daughter's transverse momentum in lab frame
1196                                 Double_t MCQtMinus = MCPQt(mcEvent, iMC, MCdaughterMinus); //daughter's transverse momentum in mother's frame (Qt)
1197                                 Double_t MCKinkAngleMinus2 = fuMCKinkAngle(mcEvent, iMC, MCdaughterMinus, kTRUE); //kink angle in degrees
1198                                 Double_t MCMaxKinkAngPimuMinus=fMaxKinkAngPimu->Eval(MCP,0.,0.,0.); //maximum decay angle in lab for pion decaying to muon
1199
1200                                 if  (MCKinkAngleMinus2>(MCMaxKinkAngPimuMinus*1.1)) continue;
1201
1202                                 hMCQtMinus->Fill(MCQtMinus);
1203                                 hMCKinkAngleMinus->Fill(MCKinkAngleMinus2);
1204                                 hMCPKinkAngMinus->Fill(MCP,MCKinkAngleMinus2);
1205                                 hMCPdgCodemdMinus->Fill(TMath::Abs(MCPdg),TMath::Abs(dcodeMinus));
1206                                 hMCPtmdMinus->Fill(MCPt, MCPtdMinus);
1207
1208                                 if ((MCQtMinus<cLowQt) || (MCQtMinus>cUpQt) || (MCKinkAngleMinus2<cLowKinkAngle)) continue;
1209
1210                                 if (dcodeMinus==cPdgMuon) { //muon daughters selection
1211                                         hMCPtPimuonMinus->Fill(MCPt);
1212                                         hMCEtaPimuonMinus->Fill(MCEta);
1213                                         hMCRapidityPimuonMinus->Fill(MCRapidity);
1214                                         hMCQtPimuonMinus->Fill(MCQtMinus);
1215                                         hMCPKinkAngPimuonMinus->Fill(MCKinkAngleMinus2);
1216                                 } else {hMCPtPiotherMinus->Fill(MCPt);///////////////////////////
1217                                         hMCEtaPiotherMinus->Fill(MCEta);
1218                                         hMCRapidityPiotherMinus->Fill(MCRapidity);
1219                                         hMCQtPiotherMinus->Fill(MCQtMinus);
1220                                         hMCPKinkAngPiotherMinus->Fill(MCKinkAngleMinus2);
1221                                 }
1222                         } //end of pion daughters' loop
1223                 }
1224
1225         } // end of MC tracks loop
1226
1227
1228
1229 //----------------------------- Physics selection ------------------------------//
1230         Bool_t IsSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected()&AliVEvent::kMB;
1231         if ( IsSelected ==kFALSE) return;
1232
1233 //--------------- Data Multiplicity after Physics selection --------------------//
1234         hMultPS->Fill(NTracks);
1235
1236 //------------------------------------ Vertex ----------------------------------//
1237         const AliESDVertex* vtx = GetEventVertex(esdEvent); //ESD primary vertex
1238         if (!vtx) return;
1239
1240         Double_t vtxpos[3]; //vertex position vector
1241         vtx->GetXYZ(vtxpos);
1242         hvtx->Fill(vtxpos[0], vtxpos[1], vtxpos[2]); //ESD primary vertex position (x-y-z)
1243         hvtxy->Fill(vtxpos[0], vtxpos[1]);
1244         hvtyz->Fill(vtxpos[1], vtxpos[2]);
1245         hvtxz->Fill(vtxpos[0], vtxpos[2]);
1246
1247         if (TMath::Abs(vtxpos[2])>10.) return;
1248
1249 //-------------------- Data Multiplicity after vertex cut -----------------------//
1250         hMultPSV->Fill(NTracks);
1251
1252 //------------------------------- dE/dx parameters -----------------------------//
1253   if(!fPIDResponse) {
1254     AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
1255     AliInputEventHandler* inputHandler = (AliInputEventHandler*)(man->GetInputEventHandler());
1256     fPIDResponse = inputHandler->GetPIDResponse();
1257   }
1258
1259 //------------------------- Reconstructed data Analysis -----------------------//
1260         for (Int_t iData = 0; iData<NTracks; iData++) { //loop on all ESD tracks
1261                 AliESDtrack *ESDtrack = esdEvent->GetTrack(iData);
1262                 if (!ESDtrack) {
1263                         Printf("ERROR: Could not receive ESD track %d", iData);
1264                         continue;
1265                 }
1266
1267                 Double_t Pt = ESDtrack->Pt();
1268                 Double_t Eta = ESDtrack->Eta();
1269
1270                 hPtAll->Fill(Pt);
1271                 hEtaAll->Fill(Eta);
1272                 //hRapidityAll->Fill(Rapidity(ESDtrack));
1273
1274                 Double_t TrackPos[3]; //starting position of ESD tracks (x-y-z)
1275                 ESDtrack->GetXYZ(TrackPos);
1276                 hTrackPos->Fill(TrackPos[0],TrackPos[1],TrackPos[2]);
1277                 hTrackPosxy->Fill(TrackPos[0], TrackPos[1]);
1278                 hTrackPosyz->Fill(TrackPos[1], TrackPos[2]);
1279                 hTrackPosxz->Fill(TrackPos[0], TrackPos[2]);
1280
1281                 if ((IsPrimaryTrack(ESDtrack) == kFALSE) || (IsGoodTrack(ESDtrack) == kFALSE)) continue; //reject bad & secondary tracks
1282                 hMultPrim->Fill(NTracks);
1283                 hPtPrim->Fill(Pt);
1284                 hEtaPrim->Fill(Eta);
1285                 //hRapidityPrim->Fill(Rapidity(ESDtrack));
1286
1287                 hPrimTrackPos->Fill(TrackPos[0],TrackPos[1],TrackPos[2]); //starting position of primary - supposed ESD tracks (x-y-z)
1288                 hPrimTrackPosxy->Fill(TrackPos[0], TrackPos[1]);
1289                 hPrimTrackPosyz->Fill(TrackPos[1], TrackPos[2]);
1290                 hPrimTrackPosxz->Fill(TrackPos[0], TrackPos[2]);
1291
1292                 if (Pt<cLowPt) continue;
1293
1294                 hPt->Fill(Pt);
1295                 hEta->Fill(Eta);
1296                 //hRapidity->Fill(fuRapidity(ESDtrack));
1297
1298                 Int_t KinkIndex = ESDtrack->GetKinkIndex(0); //kink index (1st component is negative if the track is a kink candidate)
1299                 if (KinkIndex>=0) continue; //kink selection
1300
1301                 Double_t Rapidity = fuRapidity(ESDtrack);
1302                 if (TMath::Abs(Rapidity)>cRapidityLim) continue;
1303
1304                 Double_t dEdx = ESDtrack->GetTPCsignal();
1305                 Double_t NTPCclusters = ESDtrack->GetTPCclusters(0);
1306
1307                 AliESDkink *kink = esdEvent->GetKink(TMath::Abs(KinkIndex)-1);
1308
1309                 Int_t Label1 = kink->GetLabel(0); //mother's track MC label (first component is mother)
1310                 Int_t Label2 = kink->GetLabel(1); //daughter's track MC label (second component is daughter)
1311                 if (Label1>MCNTracks) continue;
1312                 if (Label2>MCNTracks) continue;
1313
1314                 const TVector3 KinkPos(kink->GetPosition());
1315                 Double_t KinkPosR = kink->GetR(); //kink's position radius
1316                 Double_t KinkDistance = kink->GetDistance();
1317                 Double_t Qt = kink->GetQt(); //daughter's transverse momentum in mother's frame
1318                 Double_t KinkAngle =TMath::RadToDeg()*kink->GetAngle(2); //kink angle in mother frame (3rd component is angle in rads)
1319
1320                 const TVector3 PMother(kink->GetMotherP()); //ESD mother's momentum
1321                 Double_t Pmx = PMother.Px();
1322                 Double_t Pmy = PMother.Py();
1323                 Double_t Pmz = PMother.Pz();
1324                 Double_t Pm = PMother.Mag(); //ESD mother's momentum magnitude
1325                 Double_t PTrack[3];
1326                 ESDtrack->GetPxPyPz(PTrack);
1327                 TVector3 Pvector3(PTrack[0], PTrack[1], PTrack[2]);
1328                 Double_t P = Pvector3.Mag();
1329                 Double_t PTPC = ESDtrack->GetInnerParam()->GetP();
1330
1331                 const TVector3 PDaughter(kink->GetDaughterP()); //ESD daughter's momentum
1332                 Double_t Pdx = PDaughter.Px();
1333                 Double_t Pdy = PDaughter.Py();
1334                 Double_t Pdz = PDaughter.Pz();
1335                 Double_t Pd = PDaughter.Mag(); //ESD daugter's momentum magnitude
1336                 Double_t Edmu = TMath::Sqrt(TMath::Power(Pd,2)+TMath::Power(cMuonMass,2)); //ESD muon daughter's energy
1337
1338                 Double_t DP = TMath::Sqrt((Pmx-Pdx)*(Pmx-Pdx)+(Pmy-Pdy)*(Pmy-Pdy)+(Pmz-Pdz)*(Pmz-Pdz)); //transferred momentum magnitude
1339                 Double_t MinvPimu = TMath::Sqrt((Edmu+DP)*(Edmu+DP)-TMath::Power(Pm,2)); //pion mother's invariant mass when decaying to muon
1340                 Double_t MaxKinkAngPimu=fMaxKinkAngPimu->Eval(P,0.,0.,0.); //maximum decay angle in lab frame for a pion decaying to muon
1341
1342                 TParticle *particle1 = MCstack->Particle(TMath::Abs(kink->GetLabel(0))); //mother MC particle object
1343                 TParticle *particle2 = MCstack->Particle(TMath::Abs(kink->GetLabel(1))); //daughter MC particle object
1344
1345                 Int_t code1 = particle1->GetPdgCode(); //mother's pdg code obtained from MC mother object
1346                 Int_t code2 = particle2->GetPdgCode(); //daughter's pdg code obtained from MC daughter object
1347
1348                 UInt_t DauUID = particle2->GetUniqueID(); //ESD daughter's unique id (method of production)
1349
1350                 Double_t MCKinkPosZ = particle2->Vz();
1351
1352                 hPtKink->Fill(Pt);
1353                 hEtaKink->Fill(Eta);
1354                 hRapidityKink->Fill(Rapidity);
1355                 hPmP->Fill(Pm,P);
1356                 hKinkPosRTPCclusters1->Fill(NTPCclusters/KinkPosR);
1357                 hKinkPosRTPCclusters2->Fill(KinkPosR, NTPCclusters);
1358                 hQt->Fill(Qt);
1359                 hKinkAngle->Fill(KinkAngle);
1360                 hDCAkink->Fill(KinkDistance);
1361                 hPmKinkAng->Fill(P, KinkAngle);
1362                 hKinkPosXY->Fill(KinkPos[0],KinkPos[1]);
1363                 hKinkPosZY->Fill(KinkPos[2],KinkPos[1]);
1364                 hKinkPosZR->Fill(KinkPos[2],KinkPosR);
1365                 hKinkPosR->Fill(KinkPosR);
1366                 hKinkPosZ->Fill(KinkPos[2]);
1367                 hKinkPosZMCKinkPosZ->Fill(KinkPos[2], MCKinkPosZ);
1368                 hPdgCodemd->Fill(TMath::Abs(code1),TMath::Abs(code2));
1369                 hPmd->Fill(Pm,Pd);
1370                 hMinvPimu->Fill(MinvPimu);
1371                 hUIDKinkDau->Fill(DauUID);
1372                 hdEdx->Fill(PTPC, dEdx);
1373
1374                 if ( (IsRealPionKink(kink, MCstack) == kFALSE) || (Label1>MCNPrimTracks) ) {
1375                         hPtKinkFake->Fill(Pt);
1376                         hEtaKinkFake->Fill(Eta);
1377                         hRapidityKinkFake->Fill(Rapidity);
1378                         hPmPFake->Fill(Pm,P);
1379                         hKinkPosRTPCclusters1Fake->Fill(NTPCclusters/KinkPosR);
1380                         hKinkPosRTPCclusters2Fake->Fill(KinkPosR,NTPCclusters);
1381                         hQtFake->Fill(Qt);
1382                         hKinkAngleFake->Fill(KinkAngle);
1383                         hDCAkinkFake->Fill(KinkDistance);
1384                         hPmKinkAngFake->Fill(P, KinkAngle);
1385                         hKinkPosXYFake->Fill(KinkPos[0],KinkPos[1]);
1386                         hKinkPosZYFake->Fill(KinkPos[2],KinkPos[1]);
1387                         hKinkPosZRFake->Fill(KinkPos[2],KinkPosR);
1388                         hKinkPosRFake->Fill(KinkPosR);
1389                         hKinkPosZFake->Fill(KinkPos[2]);
1390                         hKinkPosZMCKinkPosZFake->Fill(KinkPos[2], MCKinkPosZ);
1391                         hPdgCodemdFake->Fill(TMath::Abs(code1),TMath::Abs(code2));
1392                         hPmdFake->Fill(Pm,Pd);
1393                         hMinvPimuFake->Fill(MinvPimu);
1394                         hUIDKinkDauFake->Fill(DauUID);
1395                         hdEdxFake->Fill(PTPC, dEdx);
1396                 }
1397
1398                 if ((KinkPosR<cLowR)||(KinkPosR>cUpR)) continue; //selection of kinks that are detected in the main region of TPC
1399
1400                 hPtPosRSelected->Fill(Pt);
1401
1402                 if ((TMath::Abs(KinkPos[2])<0.5) || (TMath::Abs(KinkPos[2])>225)){
1403                         hPdgCodemdZRejected->Fill(TMath::Abs(code1),TMath::Abs(code2));
1404                         continue;
1405                 }
1406                 hPtZSelected->Fill(Pt);
1407
1408
1409                 if  (KinkAngle<cLowKinkAngle) { //reject fake kinks (with very small kink angle)
1410                         hPdgCodemdAngRejected->Fill(TMath::Abs(code1), TMath::Abs(code2));
1411                         continue;
1412                 }
1413                 hPtAngSelected->Fill(Pt);
1414
1415                 if ((Pm/P<0.7) || (1.3<Pm/P)) {
1416                         hPdgCodemdPmRejected->Fill(TMath::Abs(code1),TMath::Abs(code2));
1417                         continue;
1418                 }//good tracks check mother track and kink if near
1419                 hPtPmSelected->Fill(Pt);
1420
1421                 if (Qt<cLowQt) { //pion kinks Qt cuts
1422                         hPdgCodemdQtLowRejected->Fill(TMath::Abs(code1),TMath::Abs(code2));
1423                         continue;
1424                 }
1425
1426                 //Good Kinks
1427                 hPtGoodKink->Fill(Pt);
1428                 hEtaGoodKink->Fill(Eta);
1429                 hRapidityGoodKink->Fill(Rapidity);
1430                 hQtGoodKink->Fill(Qt);
1431                 hPmGoodKinkAng->Fill(P, KinkAngle);
1432                 hPdgCodemdGoodKink->Fill(TMath::Abs(code1),TMath::Abs(code2));
1433                 hPmdGoodKink->Fill(Pm,Pd);
1434                 hUIDGoodKinkDau->Fill(DauUID);
1435                 hdEdxGoodKink->Fill(PTPC, dEdx);
1436
1437
1438
1439                 //---------------------- perfect PID from pdg code -------------------------//
1440
1441                 if (code1==cPdgPion) { //selection of pion plus kinks
1442                         hUIDPiDauPlus->Fill(DauUID);
1443                         if ((Label1<=MCNPrimTracks) && (DauUID==4)) { //select primary pions that decay
1444                                 hMultPiPlus->Fill(NTracks);
1445                                 hPtPiPlus->Fill(Pt);
1446                                 hEtaPiPlus->Fill(Eta);
1447                                 hRapidityPiPlus->Fill(Rapidity);
1448                                 hQtPiPlus->Fill(Qt);
1449                                 hKinkAnglePiPlus->Fill(KinkAngle);
1450                                 hPmKinkAngPiPlus->Fill(P, KinkAngle);
1451                                 hKinkPosXYPiPlus->Fill(KinkPos[0],KinkPos[1]);
1452                                 hKinkPosZRPiPlus->Fill(KinkPos[2],KinkPosR);
1453                                 hKinkPosRPiPlus->Fill(KinkPosR);
1454                                 hDCAkinkPiPlus->Fill(KinkDistance);
1455                                 hPdgCodemdPiPlus->Fill(TMath::Abs(code1),TMath::Abs(code2));
1456                                 hPmdPiPlus->Fill(Pm,Pd);
1457                                 hdEdxPiPlus->Fill(PTPC, dEdx);
1458
1459                                 if (code2==-cPdgMuon) { //selection of muon daughters
1460                                         hQtPimuPlus->Fill(Qt);
1461                                         hPmKinkAngPimuPlus->Fill(P, KinkAngle);
1462                                 } else { //other daughters?
1463                                         hQtPiotherPlus->Fill(Qt);
1464                                         hPmKinkAngPiotherPlus->Fill(P, KinkAngle);
1465                                         hPdgCodemdPiotherPlus->Fill(TMath::Abs(code1), TMath::Abs(code2));
1466                                 }
1467                         }
1468                 } else if (code1==-cPdgPion) {
1469                         hUIDPiDauMinus->Fill(DauUID);
1470                         if ((Label1<=MCNPrimTracks) && (DauUID==4)) { //select primary pions that decay
1471                                 hMultPiMinus->Fill(NTracks);
1472                                 hPtPiMinus->Fill(Pt);
1473                                 hEtaPiMinus->Fill(Eta);
1474                                 hRapidityPiMinus->Fill(Rapidity);
1475                                 hQtPiMinus->Fill(Qt);
1476                                 hKinkAnglePiMinus->Fill(KinkAngle);
1477                                 hPmKinkAngPiMinus->Fill(P, KinkAngle);
1478                                 hKinkPosXYPiMinus->Fill(KinkPos[0],KinkPos[1]);
1479                                 hKinkPosZRPiMinus->Fill(KinkPos[2],KinkPosR);
1480                                 hKinkPosRPiMinus->Fill(KinkPosR);
1481                                 hDCAkinkPiMinus->Fill(KinkDistance);
1482                                 hPdgCodemdPiMinus->Fill(TMath::Abs(code1),TMath::Abs(code2));
1483                                 hPmdPiMinus->Fill(Pm,Pd);
1484                                 hdEdxPiMinus->Fill(PTPC, dEdx);
1485
1486                                 if (code2==cPdgMuon) { //selection of muon daughters
1487                                         hQtPimuMinus->Fill(Qt);
1488                                         hPmKinkAngPimuMinus->Fill(P, KinkAngle);
1489                                 } else { //other daughters?
1490                                         hQtPiotherMinus->Fill(Qt);
1491                                         hPmKinkAngPiotherMinus->Fill(P, KinkAngle);
1492                                         hPdgCodemdPiotherMinus->Fill(TMath::Abs(code1), TMath::Abs(code2));
1493                                 }
1494                         }
1495                 }// end of pion selection
1496
1497
1498                 //------------------------ realistic PID from physical criteria --------------------//
1499
1500
1501
1502                 if (Qt>cUpQt){ //pion kinks Qt cuts
1503                         hPdgCodemdQtRejected->Fill(TMath::Abs(code1),TMath::Abs(code2));
1504                         continue;
1505                 }
1506                 hPtQtSelected->Fill(Pt);
1507
1508                 if  (KinkAngle>(MaxKinkAngPimu*1.1)) { //pion kink angle cuts
1509                         hPdgCodemdMaxAngRejected->Fill(TMath::Abs(code1),TMath::Abs(code2));
1510                         continue;
1511                 }
1512                 hPtMaxAngSelected->Fill(Pt);
1513
1514 //              if ( ((NTPCclusters/KinkPosR)>0.63) || ((NTPCclusters/KinkPosR)<0.20) ) { //good TPC tracks selection
1515                 Double_t tpcNClHigh = -51.67+ (11./12.)*KinkPosR;
1516                 Double_t tpcNClMin  = -85.5 + (65./95.)*KinkPosR;
1517                 if ( (NTPCclusters>tpcNClHigh) || (NTPCclusters<tpcNClMin) ){ //good TPC tracks selection
1518                         hPdgCodemdRTPCclustersRejected->Fill(TMath::Abs(code1),TMath::Abs(code2));
1519                         continue;
1520                 }
1521                 hPtRTPCclustersSelected->Fill(Pt);
1522                 hRTPCclustersRTPCclustersSelected->Fill(KinkPosR,NTPCclusters);
1523
1524                 if ((MinvPimu>cUpInvMass) || (MinvPimu<cLowInvMass)) { //pion kinks invariant mass cuts
1525                         hPdgCodemdMinvRejected->Fill(TMath::Abs(code1),TMath::Abs(code2));
1526                         continue;
1527                 }
1528
1529
1530 /////////////////selected without dE/dx cut
1531                 hPtSelected->Fill(Pt);
1532                 hEtaSelected->Fill(Eta);
1533                 hRapiditySelected->Fill(Rapidity);
1534                 hQtSelected->Fill(Qt);
1535                 hKinkAngleSelected->Fill(KinkAngle);
1536                 hDCAkinkSelected->Fill(KinkDistance);
1537                 hPmKinkAngSelected->Fill(P, KinkAngle);
1538                 hKinkPosXYSelected->Fill(KinkPos[0],KinkPos[1]);
1539                 hKinkPosZRSelected->Fill(KinkPos[2],KinkPosR);
1540                 hKinkPosRSelected->Fill(KinkPosR);
1541                 hPdgCodemdSelected->Fill(TMath::Abs(code1),TMath::Abs(code2));
1542                 hPmdSelected->Fill(Pm,Pd);
1543                 hMinvPimuSelected->Fill(MinvPimu);
1544                 hUIDKinkDauSelected->Fill(DauUID);
1545                 hdEdxSelected->Fill(PTPC, dEdx);
1546                 if ( (IsRealPionKink(kink, MCstack) == kFALSE) || (Label1>MCNPrimTracks) ) {
1547                         hPtSelectedFake->Fill(Pt);
1548                         hEtaSelectedFake->Fill(Eta);
1549                         hRapiditySelectedFake->Fill(Rapidity);
1550                         hQtSelectedFake->Fill(Qt);
1551                         hKinkAngleSelectedFake->Fill(KinkAngle);
1552                         hDCAkinkSelectedFake->Fill(KinkDistance);
1553                         hPmKinkAngSelectedFake->Fill(P, KinkAngle);
1554                         hKinkPosXYSelectedFake->Fill(KinkPos[0],KinkPos[1]);
1555                         hKinkPosZRSelectedFake->Fill(KinkPos[2],KinkPosR);
1556                         hKinkPosRSelectedFake->Fill(KinkPosR);
1557                         hPmdSelectedFake->Fill(Pm,Pd);
1558                         hMinvPimuSelectedFake->Fill(MinvPimu);
1559                         hdEdxSelectedFake->Fill(PTPC, dEdx);
1560                 }
1561
1562
1563                 Double_t NSigmaTPC = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(ESDtrack, AliPID::kPion));
1564                 if (NSigmaTPC>cSigmaCut)  { // dE/dx cut
1565                         hPdgCodemddEdxRejected->Fill(TMath::Abs(code1),TMath::Abs(code2));
1566                         continue;
1567                 }
1568
1569                 Double_t Sign = ESDtrack->GetSign();
1570
1571 ///////////////// dEdx selected
1572                 hPtPiSelected->Fill(Pt);
1573                 hEtaPiSelected->Fill(Eta);
1574                 hRapidityPiSelected->Fill(Rapidity);
1575                 hQtPiSelected->Fill(Qt);
1576                 hKinkAnglePiSelected->Fill(KinkAngle);
1577                 hDCAkinkPiSelected->Fill(KinkDistance);
1578                 hPmKinkAngPiSelected->Fill(P, KinkAngle);
1579                 hKinkPosRTPCclusters1PiSelected->Fill(NTPCclusters/KinkPosR);
1580                 hKinkPosRTPCclusters2PiSelected->Fill(KinkPosR, NTPCclusters);
1581                 hKinkPosXYPiSelected->Fill(KinkPos[0],KinkPos[1]);
1582                 hKinkPosZRPiSelected->Fill(KinkPos[2],KinkPosR);
1583                 hKinkPosRPiSelected->Fill(KinkPosR);
1584                 hKinkPosZPiSelected->Fill(KinkPos[2]);
1585                 hPmPPiSelected->Fill(Pm,P);
1586                 hPdgCodemdPiSelected->Fill(TMath::Abs(code1),TMath::Abs(code2));
1587                 hPmdPiSelected->Fill(Pm,Pd);
1588                 hMinvPimuPiSelected->Fill(MinvPimu);
1589                 hUIDKinkDauPiSelected->Fill(DauUID);
1590                 hdEdxPiSelected->Fill(PTPC, dEdx);
1591
1592                 if (Sign>0) {
1593                         hPtPiSelectedPlus->Fill(Pt);
1594                         hEtaPiSelectedPlus->Fill(Eta);
1595                         hRapidityPiSelectedPlus->Fill(Rapidity);
1596                         hQtPiSelectedPlus->Fill(Qt);
1597                         hKinkAnglePiSelectedPlus->Fill(KinkAngle);
1598                         hDCAkinkPiSelectedPlus->Fill(KinkDistance);
1599                         hPmKinkAngPiSelectedPlus->Fill(P, KinkAngle);
1600                         hKinkPosXYPiSelectedPlus->Fill(KinkPos[0],KinkPos[1]);
1601                         hKinkPosZRPiSelectedPlus->Fill(KinkPos[2],KinkPosR);
1602                         hPdgCodemdPiSelectedPlus->Fill(TMath::Abs(code1),TMath::Abs(code2));
1603                         hPmdPiSelectedPlus->Fill(Pm,Pd);
1604                         hMinvPimuPiSelectedPlus->Fill(MinvPimu);
1605                         hUIDPiDaumuSelectedPlus->Fill(DauUID);
1606                         hdEdxPiSelectedPlus->Fill(PTPC, dEdx);
1607
1608 //------------------- contamination of results of realistic PID ----------------//
1609                         if (IsRealPionKink(kink, MCstack) == kTRUE) { //real pion kinks
1610                                 if (Label1<=MCNPrimTracks) { //primaries
1611                                         hPtPrimPiKinksPlus->Fill(Pt);
1612                                         hEtaPrimPiKinksPlus->Fill(Eta);
1613                                         hRapidityPrimPiKinksPlus->Fill(Rapidity);
1614                                 } else { //secondary real kinks
1615                                         hPtSecondPiKinksPlus->Fill(Pt);
1616                                         hEtaSecondPiKinksPlus->Fill(Eta);
1617                                         hRapiditySecondPiKinksPlus->Fill(Rapidity);
1618                                 }
1619                         } else { //fakes
1620                                 hPtNonPiKinksPlus->Fill(Pt);
1621                                 hEtaNonPiKinksPlus->Fill(Eta);
1622                                 hRapidityNonPiKinksPlus->Fill(Rapidity);
1623                                 hPdgCodemdNonPiKinksPlus->Fill(TMath::Abs(code1),TMath::Abs(code2));
1624                 }
1625
1626                 } else if (Sign<0) {
1627                         hPtPiSelectedMinus->Fill(Pt);
1628                         hEtaPiSelectedMinus->Fill(Eta);
1629                         hRapidityPiSelectedMinus->Fill(Rapidity);
1630                         hQtPiSelectedMinus->Fill(Qt);
1631                         hKinkAnglePiSelectedMinus->Fill(KinkAngle);
1632                         hDCAkinkPiSelectedMinus->Fill(KinkDistance);
1633                         hPmKinkAngPiSelectedMinus->Fill(P, KinkAngle);
1634                         hKinkPosXYPiSelectedMinus->Fill(KinkPos[0],KinkPos[1]);
1635                         hKinkPosZRPiSelectedMinus->Fill(KinkPos[2],KinkPosR);
1636                         hPdgCodemdPiSelectedMinus->Fill(TMath::Abs(code1),TMath::Abs(code2));
1637                         hPmdPiSelectedMinus->Fill(Pm,Pd);
1638                         hMinvPimuPiSelectedMinus->Fill(MinvPimu);
1639                         hUIDPiDaumuSelectedMinus->Fill(DauUID);
1640                         hdEdxPiSelectedMinus->Fill(PTPC, dEdx);
1641
1642 //------------------- contamination of results of realistic PID ----------------//
1643                         if (IsRealPionKink(kink, MCstack) == kTRUE) { //real pion kinks
1644                                 if (Label1<=MCNPrimTracks) { //primaries
1645                                         hPtPrimPiKinksMinus->Fill(Pt);
1646                                         hEtaPrimPiKinksMinus->Fill(Eta);
1647                                         hRapidityPrimPiKinksMinus->Fill(Rapidity);
1648                                 } else { //secondary real kinks
1649                                         hPtSecondPiKinksMinus->Fill(Pt);
1650                                         hEtaSecondPiKinksMinus->Fill(Eta);
1651                                         hRapiditySecondPiKinksMinus->Fill(Rapidity);
1652                                 }
1653                         } else { //fakes
1654                                 hPtNonPiKinksMinus->Fill(Pt);
1655                                 hEtaNonPiKinksMinus->Fill(Eta);
1656                                 hRapidityNonPiKinksMinus->Fill(Rapidity);
1657                                 hPdgCodemdNonPiKinksMinus->Fill(TMath::Abs(code1),TMath::Abs(code2));
1658                 }
1659                 }
1660     } //end of ESD track loop
1661         PostData(1, fListOfHistos);
1662 }
1663
1664
1665 //________________________________________________________________________
1666 void AliAnalysisPionKinksMCESD::Terminate(Option_t *) {
1667 }
1668
1669 //_________________________________________________________________________
1670 const AliESDVertex* AliAnalysisPionKinksMCESD::GetEventVertex(AliESDEvent* esd) { //Gets ESD vertex and returns it if it is valid
1671         const AliESDVertex* vertex = esd->GetPrimaryVertexTracks();
1672         if(vertex->GetStatus()==kTRUE) return vertex;
1673         else {
1674                 vertex = esd->GetPrimaryVertexSPD();
1675                 if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>0)) return vertex;
1676                 else return 0;
1677         }
1678 }
1679
1680 //________________________________________________________________________
1681 Double_t AliAnalysisPionKinksMCESD::Energy(AliESDtrack* track) const { //calculates the energy for a pion track
1682         Double_t TrackMom[3];
1683         track->GetPxPyPz(TrackMom);
1684         TVector3 P(TrackMom[0], TrackMom[1], TrackMom[2]);
1685         //Double_t EnergyK = TMath::Sqrt(P.Mag()*P.Mag()+0.493677*0.493677); //kaon's energy
1686         Double_t EnergyPi = TMath::Sqrt(P.Mag()*P.Mag()+0.13957018*0.13957018); //pion's energy
1687         return EnergyPi;
1688 }
1689
1690 //________________________________________________________________________
1691 Double_t AliAnalysisPionKinksMCESD::GetMCRapidity(TParticle* particle) const { //calculates the rapidity of a track
1692         TVector3 P(particle->Px(),particle->Py(),particle->Pz());
1693         Double_t MCRapidity = 0.5*(TMath::Log(((particle->Energy())+P[2])/((particle->Energy())-P[2])));
1694         return MCRapidity;
1695 }
1696
1697 //________________________________________________________________________
1698 Double_t AliAnalysisPionKinksMCESD::fuRapidity(AliESDtrack* track) const { //calculates the rapidity for a track
1699         Double_t TrackMom[3];
1700         track->GetPxPyPz(TrackMom);
1701         TVector3 P(TrackMom[0], TrackMom[1], TrackMom[2]);
1702         Double_t RapidityPi = 0.5*(TMath::Log((Energy(track)+P[2])/(Energy(track)-P[2])));
1703         return RapidityPi;
1704 }
1705
1706 //________________________________________________________________________
1707 Double_t AliAnalysisPionKinksMCESD::MCPQt(AliMCEvent* mcEvent, Int_t iMC, TParticle* MCdaughter) const {
1708         TVector3 DecayMomentumPi(0,0,0); //mother's momentum when it decays
1709         TClonesArray* trArray=0;
1710         TParticle* tempParticle=0;
1711         if (mcEvent->GetParticleAndTR(iMC, tempParticle, trArray) != -1) {
1712                 AliTrackReference* MCtrackReference = static_cast<AliTrackReference*>(trArray->Last());
1713                 DecayMomentumPi.SetXYZ(MCtrackReference->Px(), MCtrackReference->Py(), MCtrackReference->Pz());
1714         } else return 0;
1715
1716         const TVector3 MCP3d(MCdaughter->Px(), MCdaughter->Py(), MCdaughter->Pz()); //daughter's momentum when producedframe
1717         Double_t MCQt = MCP3d.Perp(DecayMomentumPi); //daughter's transverse momentum in mother's frame (Qt)
1718         return MCQt;
1719 }
1720
1721 //________________________________________________________________________
1722 Double_t AliAnalysisPionKinksMCESD::fuMCKinkAngle(AliMCEvent* mcEvent, Int_t iMC, TParticle* MCdaughter, Bool_t degrees) const {
1723         Double_t MCQt=MCPQt(mcEvent, iMC, MCdaughter);
1724         Double_t MCKinkAngle = TMath::ASin(MCQt/(MCdaughter->P())); //kink angle in rads
1725         if (degrees==kFALSE) {
1726                 return MCKinkAngle;
1727      } else {
1728                 MCKinkAngle = TMath::RadToDeg()*MCKinkAngle; //kink angle in degrees
1729                 return MCKinkAngle;
1730         }
1731 }
1732
1733 //________________________________________________________________________
1734 Bool_t AliAnalysisPionKinksMCESD::IsGoodTrack(AliESDtrack* ESDtrack) const { //Checks if a track is acceptable as good
1735         UInt_t status = ESDtrack->GetStatus();
1736         if ((status&AliESDtrack::kITSrefit)==0) return kFALSE; //cut tracks that cannot be reconstructed by ITS only
1737         if ((status&AliESDtrack::kTPCrefit)==0) return kFALSE; //cut tracks that cannot be reconstructed by TPC only
1738
1739         Double_t NTPCclusters = ESDtrack->GetTPCclusters(0);
1740         Double_t TPCchi2 = ESDtrack->GetTPCchi2();
1741         if (NTPCclusters<20) return kFALSE;
1742         //if (NTPCclusters<30) return kFALSE; // cut sta 30 gia syst
1743
1744         //hTPCchi2clusters->Fill(TPCchi2/NTPCclusters);//histo to check???????
1745         if (TPCchi2/NTPCclusters>3.8) return kFALSE; //cut tracks of bad quality fit with the contributing clusters
1746
1747         /*Double_t ExtCov[15]; //external covariances matrix
1748         ESDtrack->GetExternalCovariance(ExtCov);
1749         if(ExtCov[0]>2) return kFALSE; //sigma(y^2)
1750         if(ExtCov[2]>2) return kFALSE; //sigma(z^2)
1751         if(ExtCov[5]>0.5) return kFALSE; //sigma(sinphi^2)
1752         if(ExtCov[9]>0.5) return kFALSE; //sigma(tanlamda^2)
1753         if(ExtCov[14]>2) return kFALSE; //sigma(1/Pt^2)*/
1754
1755         return kTRUE;
1756 }
1757
1758 //________________________________________________________________________
1759 Bool_t AliAnalysisPionKinksMCESD::IsPrimaryTrack(AliESDtrack* ESDtrack) const { //Checks if a track is acceptable as a primary
1760         Float_t ImpParam[2]; //DCA to vertex, 0->in x-y plane, 1->in z
1761         Float_t ImpParamCov[3]; //covariances of DCA to vertex
1762         ESDtrack->GetImpactParameters(ImpParam,ImpParamCov);
1763         //hdcaToVertexXY->Fill(ImpParam[0]);
1764         //hdcaToVertexZ->Fill(ImpParam[1]);
1765         if (ImpParamCov[0]<=0 || ImpParamCov[2]<=0) {
1766                 AliDebug (1, "Estimated DCA covariance lower or equal zero!");
1767                 ImpParamCov[0]=0; ImpParamCov[2]=0;
1768         }
1769
1770         //if((TMath::Abs(ImpParam[0])>0.3) || (TMath::Abs((ImpParam[1])>2.5))) return kFALSE; //absolute DCA cut
1771         if (!fMaxDCAtoVtxCut->AcceptTrack(ESDtrack))    return kFALSE;
1772         else return kTRUE;
1773 }
1774
1775
1776 //________________________________________________________________________
1777 Bool_t AliAnalysisPionKinksMCESD::IsRealPionKink(AliESDkink* kink, AliStack* MCstack) const { //checks if a reconstructed kink is a pion kink
1778         Int_t Label1 = kink->GetLabel(0); //mother's track MC label (first component is mother)
1779         Int_t Label2 = kink->GetLabel(1); //daughter's track MC label (second component is daughter)
1780         //if (Label1>MCNTracks) return kFALSE;
1781         //if (Label2>MCNTracks) return kFALSE;
1782
1783         TParticle *particle1 = MCstack->Particle(TMath::Abs(Label1)); //mother MC particle object
1784         TParticle *particle2 = MCstack->Particle(TMath::Abs(Label2)); //daughter MC particle object
1785
1786         Int_t code1 = particle1->GetPdgCode(); //mother's pdg code obtained from MC mother object
1787         Int_t code2 = particle2->GetPdgCode(); //daughter's pdg code obtained from MC daughter object
1788
1789         UInt_t UID = particle2->GetUniqueID(); //ESD daughter's unique id (method of production)
1790         if (UID!=4) return kFALSE;
1791
1792         if( (code1==cPdgPion) && (code2==-cPdgMuon) ){
1793                 return kTRUE;
1794         } else if ( (code1==-cPdgPion) && (code2==cPdgMuon) ){
1795                 return kTRUE;
1796         } else {
1797                 return kFALSE;
1798         }
1799 }
1800