]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/SPECTRA/Kinks/AliAnalysisKinkESDMC.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'master' of https://git.cern.ch/reps/AliRoot
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / Kinks / AliAnalysisKinkESDMC.cxx
CommitLineData
10eaad41 1/**************************************************************************
2 * Authors: Martha Spyropoulou-Stassinaki and the members
3 * of the Greek group at Physics Department of Athens University
4 * Paraskevi Ganoti, Anastasia Belogianni and Filimon Roukoutakis
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16//-----------------------------------------------------------------
17// AliAnalysisKinkESDMC class
18// Example of an analysis task for kink topology study
19// Kaons from kink topology are 'identified' in this code
20//-----------------------------------------------------------------
21
a0271296 22#include "TChain.h"
23#include "TTree.h"
10eaad41 24#include "TH1F.h"
25#include "TH2F.h"
a0271296 26#include "TH3F.h"
27#include "TH1D.h"
28#include "TH2D.h"
29#include "TParticle.h"
30#include <TVector3.h>
31#include "TF1.h"
32
33#include "AliAnalysisTask.h"
34#include "AliAnalysisManager.h"
35#include "AliTrackReference.h"
10eaad41 36
a0271296 37#include "AliVEvent.h"
10eaad41 38#include "AliESDEvent.h"
10eaad41 39#include "AliMCEvent.h"
a0271296 40#include "AliAnalysisKinkESDMC.h"
10eaad41 41#include "AliStack.h"
a0271296 42#include "AliESDpid.h"
10eaad41 43#include "AliESDkink.h"
a0271296 44#include "AliESDtrack.h"
45#include "AliESDtrackCuts.h"
46#include "AliPhysicsSelectionTask.h"
47#include "AliInputEventHandler.h"
48#include "AliESDInputHandler.h"
49#include "AliAnalysisManager.h"
50#include "AliVEvent.h"
51 #include "AliPIDResponse.h"
52///#include "AliTPCpidESD.h"
10eaad41 53
be1a7181 54ClassImp(AliAnalysisKinkESDMC)
a0271296 55
56
10eaad41 57//________________________________________________________________________
58AliAnalysisKinkESDMC::AliAnalysisKinkESDMC(const char *name)
92adf4f6 59 : AliAnalysisTaskSE(name), fHistPtESD(0),fHistPt(0),fHistQtAll(0),fHistQt1(0),fHistQt2(0)
a0271296 60 , fHistPtKaon(0),fHistPtKPDG(0),fHistEta(0),fHistEtaK(0),fptKMC(0),fMultiplMC(0),fESDMult(0),frad(0),
61 fradMC(0), fKinkKaon(0), fKinkKaonBg(0), fM1kaon(0), fgenPtEtR(0),fPtKink(0),
92adf4f6 62 fcodeH(0), fdcodeH(0), fAngMomK(0),fAngMomPi(0), fAngMomKC(0), fMultESDK(0), fMultMCK(0),
a0271296 63 fSignPtNcl(0), fSignPtEta(0), fSignPtEtaMC(0), fSignPtMC(0), fEtaNcl(0), fSignPt(0), fChi2NclTPC(0), fRatChi2Ncl(0),
64 fRadiusNcl(0), fTPCSgnlP(0), fTPCSgnlPa(0), fSignPtGen(0),
65 fRpr(0),fZpr(0), fdcatoVxXY(0), fMCEtaKaon(0),
66 fZvXv(0),fZvYv(0),fXvYv(0),fPtPrKink(0),fgenPtEtRP(0),fgenPtEtRN(0),fkinkKaonP(0),fkinkKaonN(0),
67 frapidESDK(0), frapidKMC(0), fPtKPlMC(0), fPtKMnMC(0),
68 fHistPtKaoP(0), fHistPtKaoN(0), fHiPtKPDGP(0), fHiPtKPDGN(0),fKinKBGP(0),fKinKBGN(0),
69 fQtKMu(0),fQtKPi(0),fQtKEl(0),fFakepipi(0), fFakeKPi(0),
70 fDCAkink(0), fDCAkinkBG(0), fPosiKink(0), fPosiKinkK(0), fPosiKinKXZ(0), fPosiKinKYZ(0), fPosiKinKBgZY(0),
71 fcode2(0), fcode4(0), fZkinkZDau(0),
72 fQtKMuMC(0), fQtKElMC(0), fQtKPiMC(0), fQtK3PiP(0),fQtK3PiM(0), fmaxAngMomKmu(0),
73 fPosiKinKBgZX(0), fPosiKinKBgXY(0), fMinvPi(0),fMinvKa(0),fMinvPr(0),
74 fTPCSgnlPtpc(0),
75 fTPCMomNSgnl(0), fMothKinkMomSgnl(0), fNSigmTPC(0), fTPCSgnlKinkDau(0),fcodeDau1(0),fcodeDau2(0), fMothKinkMomSgnlD(0),
fda107c4 76 fInvMassMuNuAll(0), fInvMassMuNuPt(0), fInvMassMuNuPtAll(0), fRatioCrossedRows(0), fRatioCrossedRowsKink(0), fRadiusPt(0), fRadiusPtcln(0),
5220a87d 77 fRadiusPtKaon(0), fRadiusPtPion(0), fRadiusPtFake(0), fPtCut1(0), fPtCut2(0), fPtCut3(0), fAngMomKKinks(0),
a0271296 78 flengthMCK(0), flifetiMCK(0), flifetim2(0), fLHelESDK(0),flifeInt(0), flifeYuri(0), flenYuri(0), flenTrRef(0),flifeSmall(0), flifetime(0),flifTiESDK(0),
79 flifeKink(), flenHelx(0), fradPtRapMC(0), fradPtRapDC(0), fradPtRapESD(0), fRadNclcln(0),
80 f1(0), f2(0),
5220a87d 81 fListOfHistos(0),fLowMulcut(-1),fUpMulcut(-1), fKinkRadUp(200),fKinkRadLow(130), fLowCluster(20), fLowQt(.12), fCutsMul(0),fMaxDCAtoVtxCut(0), fPIDResponse(0)
10eaad41 82
83{
84 // Constructor
85
86 // Define input and output slots here
a0271296 87
88 DefineOutput(1, TList::Class());
10eaad41 89}
90
91//________________________________________________________________________
92adf4f6 92void AliAnalysisKinkESDMC::UserCreateOutputObjects()
10eaad41 93{
a0271296 94
95 // Input slot #0 works with a TChain
96 // DefineInput(0, TChain::Class());liESDtrackCuts("Mul","Mul");
97 fCutsMul=new AliESDtrackCuts("Mul","Mul");
98 fCutsMul->SetMinNClustersTPC(70);
99 fCutsMul->SetMaxChi2PerClusterTPC(4);
100 fCutsMul->SetAcceptKinkDaughters(kFALSE);
101 fCutsMul->SetRequireTPCRefit(kTRUE);
102 // ITS
103 fCutsMul->SetRequireITSRefit(kTRUE);
104 fCutsMul->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
105 AliESDtrackCuts::kAny);
106 fCutsMul->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
107
108 fCutsMul->SetMaxDCAToVertexZ(2);
109 fCutsMul->SetDCAToVertex2D(kFALSE);
110 fCutsMul->SetRequireSigmaToVertex(kFALSE);
111
112 fCutsMul->SetEtaRange(-0.8,+0.8);
113 fCutsMul->SetPtRange(0.15, 1e10);
114//
115 fMaxDCAtoVtxCut=new AliESDtrackCuts("fMaxDCA", "fMaxDCA");
116 fMaxDCAtoVtxCut->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
117 fMaxDCAtoVtxCut->SetMaxChi2TPCConstrainedGlobal(36);
118
119 // Create histograms
120 // Called once
121
122 f1=new TF1("f1","((atan([0]*[1]*(1.0/(sqrt((x^2)*(1.0-([1]^2))-([0]^2)*([1]^2))))))*180.)/[2]",1.1,10.0);
123 f1->SetParameter(0,0.493677);
124 f1->SetParameter(1,0.9127037);
125 f1->SetParameter(2,TMath::Pi());
10eaad41 126
10eaad41 127
a0271296 128 f2=new TF1("f2","((atan([0]*[1]*(1.0/(sqrt((x^2)*(1.0-([1]^2))-([0]^2)*([1]^2))))))*180.)/[2]",0.1,10.0);
129 f2->SetParameter(0,0.13957018);
130 f2->SetParameter(1,0.2731374);
131 f2->SetParameter(2,TMath::Pi());
132//
133 Double_t gPt7K0[45] = {0., 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8,0.9,1.0,
134 1.1, 1.2, 1.3, 1.4, 1.5, 1.6,1.7,1.8,1.9, 2.0,
135 2.2, 2.4, 2.6, 2.8, 3.0, 3.3, 3.6, 3.9,
136 4.2, 4.6,5.0, 5.4, 5.9, 6.5, 7.0,7.5, 8.0,8.5, 9.2, 10., 11., 12., 13.5,15.0 }; // David K0
ab0877b8 137//
ab0877b8 138//! ! ! ! ! KINK FROM HERE --------------->
fa0ec057 139 Double_t gPt7Comb[48] = {
ab0877b8 1400.25,0.30,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0, 1.1, 1.2,
1411.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0, 2.1, 2.2, 2.3, 2.4, 2.5, 2.6, 2.7, 2.8, 2.9,
fa0ec057 1423.0, 3.2, 3.4, 3.6, 3.8, 4.0, 4.2, 4.4, 4.6, 4.8,5.0, 5.5, 6.0
ab0877b8 143 }; // 25/11/2013 from Francesco
a0271296 144/*
145 Double_t gPt7TOF[47] = { 0.2,0.25, 0.3,0.35, 0.4,0.45, 0.5,0.55, 0.6,0.65, 0.7,0.75, 0.8, 0.85, 0.9, 0.95, 1.0,
146 1.1, 1.2, 1.3, 1.4, 1.5, 1.6,1.7,1.8,1.9, 2.0,
147 2.1, 2.2, 2.3, 2.4, 2.5, 2.6, 2.7,2.8, 2.9, 3.0,
148 3.2, 3.4, 3.6, 3.8, 4.0, 4.2, 4.4, 4.6, 4.8,5.0 }; // Barbara TOF Kch
149*/
150 fHistPtESD = new TH1F("fHistPtESD", "P_{T} distribution",200, 0.0,10.0);
151 fHistPt = new TH1F("fHistPt", "P_{T} distribution",200, 0.0,10.0);
152 fHistQtAll = new TH1F("fHistQtAll", "Q_{T} distr All Kinks ",100, 0.0,.300);
153 fHistQt1= new TH1F("fHistQt1", "Q_{T} distribution",100, 0.0,.300);
154 fHistQt2= new TH1F("fHistQt2", "Q_{T} distribution",100, 0.0,.300);
155// fHistPtKaon = new TH1F("fHistPtKaon", "P_{T}Kaon distribution",200, 0.0,10.0);
ab0877b8 156 fHistPtKaon = new TH1F("fHistPtKaon", "P_{T}Kaon distribution",47,gPt7Comb );
157 fHistPtKPDG = new TH1F("fHistPtKPDG", "P_{T}Kaon distribution",47, gPt7Comb );
a0271296 158 fHistEta= new TH1F("fHistEta", "Eta distribution", 26,-1.3, 1.3);
159 fHistEtaK= new TH1F("fHistEtaK", "EtaK distribution", 26,-1.3, 1.3);
ab0877b8 160 fptKMC= new TH1F("fptKMC", "P_{T}Kaon generated",47, gPt7Comb );
a0271296 161 fMultiplMC= new TH1F("fMultiplMC", " charge particle multipl",100, 0.0,2500.);
162 fESDMult= new TH1F("fESDMult", "charge multipliESD",100, 0.0,100.);
163 frad= new TH1F("frad", "radius K ESD recon",100,0.,1000.);
164 fradMC= new TH1F("fradMC", "radius K generated",100,0.,1000.);
ab0877b8 165 fKinkKaon= new TH1F("fKinkKaon", "P_{T}Kaon kinks identi", 47, gPt7Comb );
166 fKinkKaonBg= new TH1F("fKinkKaonBg", "P_{T}Kaon kinks backgr",47 , gPt7Comb );
fda107c4 167 //fM1kaon= new TH1F("fM1kaon","Invar m(kaon) from kink->mu+netrino decay",180,0.1, 1.0);
168 fM1kaon= new TH1F("fM1kaon","Invar m(kaon) from kink->mu+netrino decay",600,0.1, 0.7);
ab0877b8 169 fgenPtEtR= new TH1F("fgenPtEtR", "P_{T}Kaon distribution", 47, gPt7Comb );
170 fPtKink= new TH1F("fPtKink", "P_{T}Kaon Kink bution",47, gPt7Comb );
a0271296 171 fcodeH = new TH2F("fcodeH", "code vrs dcode dist. kinks,K",100,0.,2500.,100,0.,2500.);
172 fdcodeH = new TH2F("fdcodeH", "code vrs dcode dist. kinks,K",100,0.,2500.,100,0.,2500.);
173 fAngMomK= new TH2F("fAngMomK","Decay angle vrs Mother Mom,K",100,0.0,10.0,80,0.,80.);
174 fAngMomPi= new TH2F("fAngMomPi","Decay angle vrs Mother Mom,Pi",100,0.0,10.0,80,0.,80.);
175 fAngMomKC= new TH2F("fAngMomKC","Decay angle vrs Mother Mom,K",100,0.0,10.0,80,0.,80.);
176 fMultESDK=new TH1F("fMultESDK", "charge multipliESD kaons",100, 0.,2500.);
177 fMultMCK=new TH1F("fMultMCK", "charge multipli MC kaons",100, 0.,2500.);
178 fSignPtNcl= new TH2F("fSignPtNcl","SignPt vrs Ncl,K",80,-4.,4.0,70,20.,160.);
179 fSignPtEta= new TH2F("fSignPtEta","SignPt vrs Eta,K",80,-4.0,4.0,30,-1.5,1.5);
180 fSignPtEtaMC= new TH2F("fSignPtEtaMC","SignPt vrs Eta,K",80,-4.0,4.0,30,-1.5,1.5);
181 fSignPtMC= new TH1F("fSignPtMC","SignPt ,K",100,-5.,5.0);
182 fEtaNcl= new TH2F("fEtaNcl","Eta vrs Ncl,K",26,-1.3,1.3,70,20.,160.);
183 fSignPt= new TH1F("fSignPt","SignPt ,K",100,-5.,5.0);
184 fChi2NclTPC= new TH2F("fChi2NclTPC","Chi2vrs nclust,K",100,0.,500., 70,20, 160);
185 fRatChi2Ncl= new TH1F("fRatChi2Ncl","Ratio chi2/nclusters in TPC,K",50,0.0,5.0);
186 fRadiusNcl = new TH2F("fRadiusNcl","KinkRadius Ncl,K",75,100.,250., 80,0, 160);
187 fTPCSgnlP = new TH2F("fTPCSgnlP","Kink TCP de/dx,K",300,0.,15.,150,0.,300);
188 fTPCSgnlPa= new TH2F("fTPCSgnlPa","Kink TCP de/dx,a",300,0.,15.,150,0.,300);
189 fSignPtGen= new TH1F("fSignPtGen","SignPtGen ,K",100,-5.0,5.0);
190 fRpr = new TH1D("fRpr", "rad distribution PID pr",100,-1.0,1.0);
191 fZpr = new TH1D("fZpr", "z distribution PID pr ",60,-15.,15.);
192 fdcatoVxXY = new TH1D("fdcatoVxXY", "dca distribution PID ",20,-1.,1.);
193 fMCEtaKaon = new TH1F("fMCEtaKaon"," Hist of Eta K -Kink Selecied",26,-1.3,1.3);
194 fZvXv= new TH2F("fZvXv","Xv-Zv main vtx",60,-0.5,0.5, 60, -15., 15.);
195 fZvYv= new TH2F("fZvYv","Yv-Zv main vtx",60,-0.5,0.5, 60, -15., 15.);
196 fXvYv= new TH2F("fXvYv","Xv-Yv main vtx", 60,-1.5,1.5, 60, -1.5, 1.5);
197 fPtPrKink=new TH1F("fPtPrKink","pt of ESD kaonKink tracks",300, 0.0,15.0);
ab0877b8 198 fgenPtEtRP= new TH1F("fgenPtEtRP", "P_{T}Kaon distribution positive", 47, gPt7Comb );
199 fgenPtEtRN= new TH1F("fgenPtEtRN", "P_{T}Kaon distribution negative", 47, gPt7Comb );
200 fkinkKaonP= new TH1F("fKinkKaonP", "P_{T}Kaon distribution positive", 47, gPt7Comb );
201 fkinkKaonN= new TH1F("fKinkKaonN", "P_{T}Kaon distribution negative", 47, gPt7Comb );
a0271296 202 frapidESDK= new TH1F("frapidESDK", "rapidity distribution", 26,-1.3, 1.3);
203 frapidKMC = new TH1F("frapidKMC ", "rapidity distri MC ",26,-1.3, 1.3);
ab0877b8 204 fPtKPlMC= new TH1F("fPtKPlMC", "P_{T}Kaon Pos generated", 47, gPt7Comb );
205 fPtKMnMC= new TH1F("fPtKMnMC", "P_{T}Kaon Minusgenerated",47 , gPt7Comb );
206 fHistPtKaoP= new TH1F("fHistPtKaoP", "P_{T}Kaon Pos ESD", 47, gPt7Comb );
207 fHistPtKaoN= new TH1F("fHistPtKaoN", "P_{T}Kaon Neg ESD", 47, gPt7Comb );
208 fHiPtKPDGP= new TH1F("fHiPtKPDGP", "P_{T}Kaon Pos ESD", 47, gPt7Comb );
209 fHiPtKPDGN= new TH1F("fHiPtKPDGN", "P_{T}Kaon neg ESD", 47, gPt7Comb );
210 fKinKBGP = new TH1F("fKinKBGP ", "P_{T}Kaon Pos ESD", 47, gPt7Comb );
211 fKinKBGN= new TH1F("fKinKBGN", "P_{T}Kaon neg ESD", 47, gPt7Comb );
a0271296 212 fQtKMu= new TH1F("fQtKMu", "Q_{T} distribution K to mu ",100, 0.0,.300);
213 fQtKPi= new TH1F("fQtKPi", "Q_{T} distribution K to pi",100, 0.0,.300);
214 fQtKEl= new TH1F("fQtKEl", "Q_{T} distribution K to elec",100, 0.0,.300);
ab0877b8 215 fFakepipi = new TH1F("fFakepipi", "P_{T}fake pipi ",47 , gPt7Comb );
216 fFakeKPi = new TH1F("fFakeKPi", "P_{T}fake Kpi ", 47, gPt7Comb );
a0271296 217 fDCAkink = new TH1F("fDCAkink", "DCA kink vetrex ",50, 0.0,1.0);
218 fDCAkinkBG = new TH1F("fDCAkinkBG", "DCA kink vetrex ",50, 0.0,1.0);
219 fPosiKink= new TH2F("fPosiKink", "Y vrx kink Vrex ",100, -300.0,300.0,100, -300, 300.);
220 fPosiKinkK= new TH2F("fPosiKinkK", "Y vrx kink VrexK ",100, -300.0,300.0,100, -300, 300.);
221 fPosiKinKXZ= new TH2F("fPosiKinKXZ", "Y vrx kink VrexK ",100, -300.0,300.0,100, -300, 300.);
222 fPosiKinKYZ= new TH2F("fPosiKinKYZ", "Y vrx kink VrexK ",100, -300.0,300.0,100, -300, 300.);
223 fPosiKinKBgZY= new TH2F("fPosiKinKBgZY", "Y vrx Z kink Vrexbg ",100, -300.0,300.0,100, -300, 300.);
224 fcode2 = new TH2F("fcode2", "code vrs dcode dist. kinks,K",100,0.,2500.,100,0.,2500.);
225 fcode4 = new TH2F("fcode4", "code vrs dcode dist. kinks,K",100,0.,2500.,100,0.,2500.);
226 fZkinkZDau = new TH2F("fZkinkZDau", "Y vrx kink VrexK ",100, -300.0,300.0,100, -300, 300.);
227 fQtKMuMC= new TH1F("fQtKMuMC", "Q_{T} distribution K to mu MC",100, 0.0,.300);
228 fQtKPiMC= new TH1F("fQtKPiMC", "Q_{T} distribution K to pi MC",100, 0.0,.300);
229 fQtKElMC= new TH1F("fQtKElMC", "Q_{T} distribution K to elec MC",100, 0.0,.300);
230 fQtK3PiP= new TH1F("fQtK3PiP", "Q_{T} distribution K to 3pi ",100, 0.0,.300);
231 fQtK3PiM= new TH1F("fQtK3PiM", "Q_{T} distribution K to 3pi ",100, 0.0,.300);
232 fmaxAngMomKmu= new TH2F("fmaxAngMomKmu","Decay angle vrs Mother Mom,Kmu",100,0.0,10.0,80,0.,80.);
233 fPosiKinKBgZX= new TH2F("fPosiKinKBgZX", "X vrx Z kink Vrexbg ",100, -20.0,20.0,100, 0., 300.);
234 fPosiKinKBgXY= new TH2F("fPosiKinKBgXY", "Y vrx X kink Vrexbg ",100, -300.0,300.0,100, -300, 300.);
235 fMinvPi= new TH1F("fMinvPi","Invar m(kaon) from kink-> decay",100,0.0, 1.2);
236 fMinvKa= new TH1F("fMinvKa","Invar m(kaon) from kink-> decay",100,0.0, 2.0);
237 fMinvPr= new TH1F("fMinvPr","Invar m(kaon) from kink-> decay",100,0.0, 1.2);
238 fTPCSgnlPtpc = new TH2F("fTPCSgnlPtpc","TPC signal de/dx Mom TPC,K ",100,0.0,8.0,100, 0., 250. );
239 fTPCMomNSgnl = new TH2F("fTPCMomNsgnl","TPC signal de/dx Mom TPC,K ",100,0.0,8.0,20 , -10., 10.);
240 fMothKinkMomSgnl = new TH2F("fMothKinkMomSgnl","TPC signal de/dx Mom TPC,Kink ",100,0.0,250.,100, 0., 250. );
241 fNSigmTPC = new TH1F("fNSigmTPC","TPC Nsigma de/dx TPC,K ", 30 , -7.5, 7.5);
242 fTPCSgnlKinkDau = new TH2F("fTPCSgnlKinkDau","TPC signal de/dx Mom,K",100,0.0,8.0,100,0.,250.);
243 fcodeDau1 = new TH2F("fcodeDau1", "code vrs dcode dist. kinks,K",100,0.,500.,100,0.,500.);
244 fcodeDau2 = new TH2F("fcodeDau2", "code vrs dcode dist. kinks,K",100,0.,500.,100,0.,500.);
245 fMothKinkMomSgnlD = new TH2F("fMothKinkMomSgnlD","TPC signal de/dx Mom TPC,Kink ",100,0.0,250.,100, 0., 250. );
fda107c4 246 //fInvMassMuNuAll = new TH1F("fInvMassMuNuAll","Invar from kink->mu+netrino decay",180,0.1, 1.0);
247 fInvMassMuNuAll = new TH1F("fInvMassMuNuAll","Invar from kink->mu+netrino decay",600,0.1, 0.7); // 22/8/2013
248 //fInvMassMuNuPt = new TH2F("fInvMassMuNuPt","Invar from kink->mu+netrino decay vs Pt",180,0.1, 1.0, 100, 0. , 10.);
249 fInvMassMuNuPt = new TH2F("fInvMassMuNuPt","Invar from kink->mu+netrino decay vs Pt",600,0.1, 0.7, 100, 0. , 10.); // 22/8/2013
250 fInvMassMuNuPtAll = new TH2F("fInvMassMuNuPtAll","Invar from kink->mu+netrino decay vs Pt",600,0.1, 0.7, 100, 0. , 10.); // 22/8/2013
a0271296 251 fRatioCrossedRows = new TH1F("fRatioCrossedRows","Ratio crossed rows in TPC",20,0.0,1.0);
252 fRatioCrossedRowsKink = new TH1F("fRatioCrossedRowsKink","Ratio crossed rows in TPC for kinks",20,0.0,1.0);
253 fRadiusPt =new TH2F("fRadiusPt","radius vs pt ",80, 90.,250.,100, 0.,10. );
254 fRadiusPtcln =new TH2F("fRadiusPtcln","radius vs pt clean ",80, 90.,250.,100, 0.,10. );
5220a87d 255 fRadiusPtKaon =new TH2F("fRadiusPtKaon","radius vs pt Kaon PDG ",80, 90.,250.,100, 0.,10. );
256 fRadiusPtPion =new TH2F("fRadiusPtPion","radius vs pt Pion PDG ",80, 90.,250.,100, 0.,10. );
257 fRadiusPtFake =new TH2F("fRadiusPtFake","radius vs pt Pion Fake ",80, 90.,250.,100, 0.,10. );
a0271296 258 fPtCut1 = new TH1F("fPtCut1", "P_{T}Kaon distribution",300, 0.0,15.0);
259 fPtCut2 = new TH1F("fPtCut2", "P_{T}Kaon distribution",300, 0.0,15.0);
260 fPtCut3 = new TH1F("fPtCut3", "P_{T}Kaon distribution",300, 0.0,15.0);
261 fAngMomKKinks = new TH2F("fAngMomKKinks","Decay angle vrs Mother Mom,Kinks",300,0.0,15.0,100,0.,100.);
262 flengthMCK=new TH1F("flengthMCK", "length of K MCref decay ",100,0.,1000.0);
263 flifetiMCK=new TH1F("flifetiMCK", "lifetime ref K Decay ",100,0.,1000.0);
264 flifetim2 =new TH1F("flifetim2", "lifetime ref K Decay ",100,0.,1000.0);
265 fLHelESDK =new TH1F("fLHelESDK", "lifetime ref K Decay ",100,0.,1000.0);
266 flifeInt =new TH1F("flifeInt", "lifetime ref K Decay ",100,0.,1000.0);
267 flifeYuri=new TH1F("flifeYuri","lifetime ref K Decay ",100,0.,1000.0);
268 flenYuri=new TH1F("flenYuri","lifetime ref K Decay ",100,0.,1000.0);
269 flenTrRef =new TH1F("flenTrRef","lifetime ref K Decay ",100,0.,1000.0);
270 flifeSmall=new TH1F("flifeSmall","lifetime ref K Decay ",100,0.,1000.0);
271 flifetime =new TH1F("flifetime","lifetime ref K Decay ",100,0.,1000.0);
272 flifTiESDK=new TH1F("flifTiESDK","lifetime ref K Decay ",100,0.,1000.0);
273 flifeKink =new TH1F("flifeKink", "lifetime ref K Decay ",100,0.,1000.0);
274 flenHelx =new TH1F("flenHelx", "lifetime ref K Decay ",100,0.,1000.0);
275 fradPtRapMC=new TH3F("fradPtRapMC","rad pt rap dat",28,100.,240., 20, 0., 5., 20, -1., 1. );
276 fradPtRapDC=new TH3F("fradPtRapDC","rad pt rap dat",28,100.,240., 20, 0., 5., 20, -1., 1. );
277 fradPtRapESD=new TH3F("fradPtRapESD","rad pt rap dat",28,100.,240., 20, 0., 5., 20, -1., 1. );
278 fRadNclcln = new TH2F("fRadNclcln","KinkRadius Ncl,K",75,100.,250., 80,0, 160);
279
280
281
282
283 fListOfHistos=new TList();
284 fListOfHistos->SetOwner();
10eaad41 285 fListOfHistos->Add(fHistPtESD);
286 fListOfHistos->Add(fHistPt);
287 fListOfHistos->Add(fHistQtAll);
288 fListOfHistos->Add(fHistQt1);
289 fListOfHistos->Add(fHistQt2);
290 fListOfHistos->Add(fHistPtKaon);
291 fListOfHistos->Add(fHistPtKPDG);
292 fListOfHistos->Add(fHistEta);
293 fListOfHistos->Add(fHistEtaK);
294 fListOfHistos->Add(fptKMC);
295 fListOfHistos->Add(fMultiplMC);
296 fListOfHistos->Add(fESDMult);
10eaad41 297 fListOfHistos->Add(frad);
a0271296 298 fListOfHistos->Add(fradMC);
10eaad41 299 fListOfHistos->Add(fKinkKaon);
300 fListOfHistos->Add(fKinkKaonBg);
301 fListOfHistos->Add(fM1kaon);
302 fListOfHistos->Add(fgenPtEtR);
303 fListOfHistos->Add(fPtKink);
10eaad41 304 fListOfHistos->Add(fcodeH);
305 fListOfHistos->Add(fdcodeH);
10eaad41 306 fListOfHistos->Add(fAngMomK);
307 fListOfHistos->Add(fAngMomPi);
92adf4f6 308 fListOfHistos->Add(fAngMomKC);
309 fListOfHistos->Add(fMultESDK);
310 fListOfHistos->Add(fMultMCK);
a0271296 311 fListOfHistos->Add(fSignPtNcl);
312 fListOfHistos->Add(fSignPtEta);
313 fListOfHistos->Add(fSignPtEtaMC);
314 fListOfHistos->Add(fSignPtMC);
315 fListOfHistos->Add(fEtaNcl);
316 fListOfHistos->Add(fSignPt);
317 fListOfHistos->Add(fChi2NclTPC);
318 fListOfHistos->Add(fRatChi2Ncl);
319 fListOfHistos->Add(fRadiusNcl);
320 fListOfHistos->Add(fTPCSgnlP);
321 fListOfHistos->Add(fTPCSgnlPa);
322 fListOfHistos->Add(fSignPtGen);
10eaad41 323 fListOfHistos->Add(fRpr);
324 fListOfHistos->Add(fZpr);
a0271296 325 fListOfHistos->Add(fdcatoVxXY);
326 fListOfHistos->Add(fMCEtaKaon);
10eaad41 327 fListOfHistos->Add(fZvXv);
328 fListOfHistos->Add(fZvYv);
329 fListOfHistos->Add(fXvYv);
330 fListOfHistos->Add(fPtPrKink);
a0271296 331 fListOfHistos->Add(fgenPtEtRP);
332 fListOfHistos->Add(fgenPtEtRN);
333 fListOfHistos->Add(fkinkKaonP);
334 fListOfHistos->Add(fkinkKaonN);
335 fListOfHistos->Add(frapidESDK);
336 fListOfHistos->Add(frapidKMC);
337 fListOfHistos->Add(fPtKPlMC);
338 fListOfHistos->Add(fPtKMnMC);
339 fListOfHistos->Add(fHistPtKaoP);
340 fListOfHistos->Add(fHistPtKaoN);
341 fListOfHistos->Add(fHiPtKPDGP);
342 fListOfHistos->Add(fHiPtKPDGN);
343 fListOfHistos->Add(fKinKBGP);
344 fListOfHistos->Add(fKinKBGN);
345 fListOfHistos->Add(fQtKMu);
346 fListOfHistos->Add(fQtKPi);
347 fListOfHistos->Add(fQtKEl);
348 fListOfHistos->Add(fFakepipi);
349 fListOfHistos->Add(fFakeKPi);
350 fListOfHistos->Add(fDCAkink);
351 fListOfHistos->Add(fDCAkinkBG);
352 fListOfHistos->Add(fPosiKink);
353 fListOfHistos->Add(fPosiKinkK);
354 fListOfHistos->Add(fPosiKinKXZ);
355 fListOfHistos->Add(fPosiKinKYZ);
356 fListOfHistos->Add(fPosiKinKBgZY);
357 fListOfHistos->Add(fcode2);
358 fListOfHistos->Add(fcode4);
359 fListOfHistos->Add(fZkinkZDau);
360 fListOfHistos->Add(fQtKMuMC);
361 fListOfHistos->Add(fQtKPiMC);
362 fListOfHistos->Add(fQtKElMC);
363 fListOfHistos->Add(fQtK3PiP);
364 fListOfHistos->Add(fQtK3PiM);
365 fListOfHistos->Add(fmaxAngMomKmu);
366 fListOfHistos->Add(fPosiKinKBgZX);
367 fListOfHistos->Add(fPosiKinKBgXY);
368 fListOfHistos->Add(fMinvPi);
369 fListOfHistos->Add(fMinvKa);
370 fListOfHistos->Add(fMinvPr);
371 fListOfHistos->Add(fTPCSgnlPtpc);
372 fListOfHistos->Add(fTPCMomNSgnl);
373 fListOfHistos->Add(fMothKinkMomSgnl);
374 fListOfHistos->Add(fNSigmTPC);
375 fListOfHistos->Add(fTPCSgnlKinkDau);
376 fListOfHistos->Add(fcodeDau1);
377 fListOfHistos->Add(fcodeDau2);
378 fListOfHistos->Add(fMothKinkMomSgnlD);
379 fListOfHistos->Add(fInvMassMuNuAll);
380 fListOfHistos->Add(fInvMassMuNuPt);
fda107c4 381 fListOfHistos->Add(fInvMassMuNuPtAll);
a0271296 382 fListOfHistos->Add(fRatioCrossedRows);
383 fListOfHistos->Add(fRatioCrossedRowsKink);
384 fListOfHistos->Add(fRadiusPt);
385 fListOfHistos->Add(fRadiusPtcln);
5220a87d 386 fListOfHistos->Add(fRadiusPtKaon);
387 fListOfHistos->Add(fRadiusPtPion);
388 fListOfHistos->Add(fRadiusPtFake);
a0271296 389 fListOfHistos->Add(fPtCut1);
390 fListOfHistos->Add(fPtCut2);
391 fListOfHistos->Add(fPtCut3);
392 fListOfHistos->Add(fAngMomKKinks);
393 fListOfHistos->Add(flengthMCK);
394 fListOfHistos->Add(flifetiMCK);
395 fListOfHistos->Add(flifetim2);
396 fListOfHistos->Add(fLHelESDK);
397 fListOfHistos->Add(flifeInt);
398 fListOfHistos->Add(flifeYuri);
399 fListOfHistos->Add(flenYuri);
400 fListOfHistos->Add(flenTrRef);
401 fListOfHistos->Add(flifeSmall);
402 fListOfHistos->Add(flifetime);
403 fListOfHistos->Add(flifTiESDK);
404 fListOfHistos->Add(flifeKink);
405 fListOfHistos->Add(flenHelx);
406 fListOfHistos->Add(fradPtRapMC);
407 fListOfHistos->Add(fradPtRapDC);
408 fListOfHistos->Add(fradPtRapESD);
409 fListOfHistos->Add(fRadNclcln);
410
411 PostData(1, fListOfHistos);
10eaad41 412}
413
414//________________________________________________________________________
92adf4f6 415void AliAnalysisKinkESDMC::UserExec(Option_t *)
10eaad41 416{
417 // Main loop
418 // Called for each event
419
420 // Process MC truth, therefore we receive the AliAnalysisManager and ask it for the AliMCEventHandler
421 // This handler can return the current MC event
92adf4f6 422
423 AliVEvent *event = InputEvent();
424 if (!event) {
425 Printf("ERROR: Could not retrieve event");
426 return;
427 }
428
429 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
430 if (!esd) {
431 Printf("ERROR: Could not retrieve esd");
432 return;
a0271296 433}
10eaad41 434
92adf4f6 435 AliMCEvent* mcEvent = MCEvent();
10eaad41 436 if (!mcEvent) {
437 Printf("ERROR: Could not retrieve MC event");
438 return;
439 }
a0271296 440 fMultMCK->Fill(mcEvent->GetNumberOfTracks() );
441//
442
443
444// multiplicity selection
445 Float_t refmultiplicity=fCutsMul->CountAcceptedTracks(esd);
446 if(fLowMulcut>-1)
447 {
448 if(refmultiplicity<fLowMulcut)
449 return;
450 }
451 if(fUpMulcut>-1)
452 {
453 if(refmultiplicity>fUpMulcut)
454 return;
455 }
456
457
458 fMultESDK->Fill(refmultiplicity);
459
460
461
10eaad41 462
7ccf0419 463 Printf("MC particles: %d", mcEvent->GetNumberOfTracks());
10eaad41 464
10eaad41 465 AliStack* stack=mcEvent->Stack();
466
467//primary tracks in MC
468 Int_t nPrim = stack->GetNprimary();
469//
a0271296 470
471 // loop over mc particles
472
473 // variable to count tracks
474 Int_t nMCTracks = 0;
475 Int_t nMCKinkKs = 0;
476
477// 15/2/12 OK for (Int_t iMc = 0; iMc < nPrim; iMc++)
478 //for (Int_t iMc = 0; iMc < stack->GetNtrack(); iMc++)
479 for (Int_t iMc = 0; iMc < mcEvent->GetNumberOfTracks(); iMc++)
480 {
481
482 TParticle* particle = stack->Particle(iMc);
483
484 if (!particle)
485 {
486 // AliDebug(AliLog::kError, Form("UNEXPECTED: particle with label %d not found in stack (mc loop).", iMc));
487 continue;
488 }
489// keep only primaries
490 if (!stack->IsPhysicalPrimary(iMc)) continue;
491
492 //
493
494 Double_t ptK = particle->Pt();
495
496 if( ptK <0.200) continue; // 12/2/2012
497//
498
499 Float_t charg=0;
500 Float_t code = particle->GetPdgCode();
501 Int_t mcProcess=-1011;
502//---------------------------------------kaon selection
503 if ((code==321)||(code==-321)){
504
505
506 Double_t etracKMC= TMath::Sqrt(particle->P() *particle->P() + 0.493677 *0.493677 );
507 Double_t rapidiKMC = 0.5 * (TMath::Log( (etracKMC +particle->Pz())/( etracKMC-particle->Pz() )) ) ;
508
509 if ( TMath::Abs( rapidiKMC) > 0.7) continue; //
510 frapidKMC ->Fill(rapidiKMC) ; //18/feb rapiddistr of PDG kink ESD kaons
511
512
513// maximum angle vrs momentum
514 // Double_t maxDecAnKmu=f1->Eval(particle->P(), 0.,0.,0.);
515 // fmaxAngMomKmu->Fill(particle->P() , maxDecAnKmu);
516
517 if(code > 0 ) charg =1;
518 if(code < 0 ) charg =-1;
519 Float_t chargPt= charg*ptK;
520
521 fptKMC->Fill(ptK);
522 fSignPtGen->Fill(chargPt);// kaon gensign pt
523 if (charg==1 ) fPtKPlMC->Fill( ptK );
524 if ( charg==-1) fPtKMnMC->Fill( ptK );
525// primary vertex
526 // Double_t mVx=particle->Vx();
5220a87d 527 // Double_t mVy=particle->Vy();
a0271296 528 Double_t mVz=particle->Vz();
529// 25/11/2012 ???????back 10/1/2013
530 TClonesArray* trArray=0;
531 TParticle* tempParticle=0;
532
533 TVector3 DecayMomentumK(0,0,0);
534 Float_t lengthKMC=0;
535 if (mcEvent->GetParticleAndTR(iMc, tempParticle, trArray) != -1) {
536 AliTrackReference* MCtrackReference = static_cast<AliTrackReference*>(trArray->Last());
537 lengthKMC = MCtrackReference->GetLength();
538
539// DecayMomentumK.SetXYZ(MCtrackReference->Px(), MCtrackReference->Py(), MCtrackReference->Pz());
540 }
541 flenTrRef ->Fill(lengthKMC);
542 flifetime ->Fill( (lengthKMC*0.493667 /particle->P())); // 19/7
543 if ((lengthKMC>100.)&& (lengthKMC<300.) ) flifeSmall->Fill( (lengthKMC*0.493667/particle->P() ) );
544
545 Int_t nMCKpi =0;
546 Int_t mcProc4 =0;
547 Int_t mcProc13=0;
548 Double_t radiusD=0;
549 // Double_t lengthK =0.;
550 Double_t LengthK =0.;
551 Double_t lenYuri =0.;
552 Double_t MCQt =0.;
553 Double_t MCQt3[2];
554 Int_t firstD=particle->GetFirstDaughter();
555 Int_t lastD=particle->GetLastDaughter();
556
557 if( (lastD<=0) || (firstD<=0) ) continue;
558
559 if ( lastD > mcEvent->GetNumberOfTracks() ) continue;
560 if (firstD > mcEvent->GetNumberOfTracks() ) continue;
561// 25/112012
562// TClonesArray* trArray=0;
563// TParticle* tempParticle=0;
564 // TVector3 DecayMomentumK(0,0,0);
565//loop on secondaries
566 for (Int_t k=firstD;k<=lastD;k++) {
567 if ( k > 0 ) {
568 TParticle* daughter1=stack->Particle(k); // 27/8
569 Float_t dcode = daughter1->GetPdgCode();
570
571// mother momentum trackrefs and QtMC // 17/9/2010, test Feb 2011
572 if (mcEvent->GetParticleAndTR(iMc, tempParticle, trArray) != -1) {
573 AliTrackReference* MCtrackReference = static_cast<AliTrackReference*>(trArray->Last());
574 DecayMomentumK.SetXYZ(MCtrackReference->Px(), MCtrackReference->Py(), MCtrackReference->Pz());
575 }
576 const TVector3 MCP3d(daughter1->Px(), daughter1->Py(), daughter1->Pz()); //MC daughter's momentum
577 MCQt = MCP3d.Perp(DecayMomentumK); //MC daughter's transverse momentum in mother's frame
578//
579 Double_t MCKinkAngle = TMath::ASin(MCQt/daughter1->P() );
580 Double_t MCKinkAngle2= TMath::RadToDeg() * MCKinkAngle; // in degrees
581 // fmaxAngMomKmu->Fill(particle->P() , MCKinkAngle2);// MC
582//
583 mcProcess=daughter1->GetUniqueID();
584 radiusD=daughter1->R();
585// secondary vertex
5220a87d 586 // Double_t hVx=daughter1->Vx();
587 // Double_t hVy=daughter1->Vy();
a0271296 588 Double_t hVz=daughter1->Vz();
589
590 LengthK = TMath::Sqrt( radiusD*radiusD + ( mVz-hVz) * (mVz-hVz) ); // 19/7/2010 mss
591
592// lengthK = TMath::Sqrt( (mVx -hVx)*( mVx -hVx) + ( mVy -hVy)* (mVy -hVy ) + ( mVz -hVz ) * (mVz -hVz) );
593 lenYuri = (TMath::Abs( mVz-hVz))* (TMath::Sqrt( 1.+ ( ptK*ptK)/ (particle->Pz() * particle->Pz()) )) ;
594
595 if(mcProcess==13) {
596 mcProc13++;
597
598 if(mcProc13==1) flifeInt ->Fill( (lengthKMC*0.493667 /particle->P())); // 19/7
599
600 // if( (charg==1)&&(mcProc13==1 )) fradIntKP->Fill(daughter1->R());
601
602 // if( ( charg ==-1)&&(mcProc13==1)) fradIntKM->Fill(daughter1->R());
603 }
604
605
606//
607 if (mcProcess==4) {
608
609 mcProc4++;
610 if ( mcProc4==1) {
611 // 10/1/13 if( (radiusD >120.)&&( radiusD< 210.) ) {
612 flifeYuri ->Fill( (lenYuri *0.493667 /particle->P())); // 19/7
613 flifetiMCK->Fill(LengthK);
614 flenYuri ->Fill(lenYuri);
615// 10/1/13 }
616
617 // flengthMCK->Fill(lengthK); //
618 // flifetime ->Fill( (lengthK*0.493667 /particle->P())); // 19/7
619 flengthMCK->Fill(lengthKMC); //
620 flifetime ->Fill( (lengthKMC*0.493667 /particle->P())); // 19/7
621 fradPtRapMC->Fill( radiusD, 1./ptK, rapidiKMC); // systematics 26/8
622 }
623
5220a87d 624
a0271296 625//
626 if ( ( ( code==321 )&&( dcode ==-13 ))||( ( code ==-321)&&(dcode== 13) ) || ( ( code==321 )&&( dcode ==-11 )) || ( (code ==-321)&&(dcode== 11))) {
627 flifetim2 ->Fill( (lengthKMC*0.493667 /particle->P()));
628 // 8/2/2013allgh radius if( (radiusD >120.)&&( radiusD< 210.) )
629 // 14/2/13 if( (radiusD >130.)&&( radiusD< 200.) )
630 if( (radiusD >fKinkRadLow )&&( radiusD< fKinkRadUp) )
631 fmaxAngMomKmu->Fill(particle->P() , MCKinkAngle2);// MC
632 }
633
634 if (( (TMath::Abs(code)==321 )&&(TMath::Abs(dcode) ==211 ))&& ( mcProc4<2)) flifetim2->Fill( lengthKMC *0.493667 /particle->P()) ;//19/7
635
5220a87d 636 // test feb 2013 if ((TMath::Abs(hVz)<0.5) || (TMath::Abs(hVz )>225)) continue;
a0271296 637/// inside radius region ----------------------------------------------
638 if(MCKinkAngle2 < 2.) continue; // as in ESD
639 // ====== 8/2/13 if (((daughter1->R())>120)&&((daughter1->R())<210)&& (MCQt>0.120) ){
886c7202 640 if (((daughter1->R())> fKinkRadLow )&&((daughter1->R())< fKinkRadUp )&& (MCQt>fLowQt) ){
a0271296 641
642 if ( ( code==321 )&&( dcode ==-13 )) {
643 fradPtRapDC->Fill( radiusD, 1./ptK, rapidiKMC); // systematics 26/8
644 fradMC->Fill(daughter1->R());
645 fQtKMuMC ->Fill(MCQt );//to muon
646 fgenPtEtR->Fill( ptK );//to muon
647 fgenPtEtRP->Fill( ptK );//to muon
648 fMCEtaKaon ->Fill(rapidiKMC );//to muon
649 fSignPtEtaMC->Fill(ptK,rapidiKMC );//to muon
650 fSignPtMC->Fill(ptK);//to muon
651 flifetiMCK->Fill( lenYuri*0.4933667/particle->P() );
652 // flifetiMCK->Fill( LengthK*0.4933667/ ptK
653 } // positive kaon
654 if ( ( code ==-321)&&(dcode== 13)){
655 fgenPtEtR->Fill( ptK );//to muon
656 fQtKMuMC ->Fill(MCQt );//to muon
657 fgenPtEtRN->Fill(ptK); //
658 fSignPtEtaMC->Fill(chargPt,rapidiKMC );//to muon
659 fMCEtaKaon ->Fill(rapidiKMC );//to muon
660 fradPtRapDC->Fill( radiusD, 1./ptK, rapidiKMC); // systematics 26/8
661 fradMC->Fill(daughter1->R());
662 fSignPtMC->Fill(chargPt);//to muon
663 flifetiMCK->Fill( lenYuri*0.4933667/particle->P() ) ;
664 // flifetiMCK->Fill( LengthK*0.4933667/ ptK ) ;
665
666 } // negative code
667 // if (( ( code==321 )&& ( dcode ==211 ))|| (( code == -321 )&& ( dcode ==-211))) fgenPtEtR->Fill( ptK );//to pion
668 if ( ( code==321 )&&( dcode ==-11 )) {
669 fQtKElMC ->Fill(MCQt );//to muon
670 fgenPtEtR->Fill( ptK );//to electron
671 fgenPtEtRP->Fill( ptK );//to muon
672 fMCEtaKaon ->Fill(rapidiKMC );//to electron
673 fradPtRapDC->Fill( radiusD, 1./ptK, rapidiKMC); // systematics 26/8
674 fradMC->Fill(daughter1->R());
675 fSignPtEtaMC->Fill(ptK,rapidiKMC );//to electron
676 fSignPtMC->Fill(ptK);
677 flifetiMCK->Fill( lenYuri*0.4933667/particle->P() );
678 // flifetiMCK->Fill( LengthK*0.4933667/ ptK );
679
680 } // positive kaon
681 if ( ( code ==-321)&&(dcode== 11)){
682 fgenPtEtR->Fill( ptK );//to electron
683 fQtKElMC ->Fill(MCQt );//to muon
684 fgenPtEtRN->Fill(ptK); //
685 fSignPtEtaMC->Fill(chargPt,rapidiKMC );//to electron
686 fMCEtaKaon ->Fill(rapidiKMC );//to electron
687 fradPtRapDC->Fill( radiusD, 1./ptK, rapidiKMC); // systematics 26/8
688 fradMC->Fill(daughter1->R());
689 fSignPtMC->Fill(chargPt);//to electron
690 flifetiMCK->Fill( lenYuri*0.4933667/particle->P() );
691 // flifetiMCK->Fill( LengthK*0.4933667/ ptK );
692
693 } // negative code
694
695 if (( ( code==321 )&& ( dcode ==211 ))|| (( code == -321 )&& ( dcode ==-211))) nMCKpi++ ;
696 if (( ( code==321 )&& ( dcode ==211 ))|| (( code == -321 )&& ( dcode ==-211))) {
697 if ( nMCKpi > 0) {
698 MCQt3[nMCKpi-1] = MCQt ;// k to pipipi
699}
700 }
701 nMCKinkKs++;
702 }
703
704 }// decay
705 } // positive k
706 }// daughters
707
708
709
710 if ( code > 0) {
711 if( nMCKpi == 1) fgenPtEtR->Fill(ptK); // k to pipi
712 if( nMCKpi == 1) fgenPtEtRP->Fill(ptK); // k to pipi
713 if( nMCKpi == 1) fSignPtEtaMC->Fill(ptK,rapidiKMC ); // k to pipi
714 if( nMCKpi == 1) fSignPtMC->Fill(ptK); // k to pipi
715 if( nMCKpi == 1) fMCEtaKaon->Fill(rapidiKMC ); // k to pipi
716 if(nMCKpi==1) fradMC->Fill(radiusD );
717 if(nMCKpi==1) fQtKPiMC->Fill( MCQt );
718 if (nMCKpi==1) fradPtRapDC->Fill( radiusD, 1./ptK, rapidiKMC); // systematics 26/8
719 if(nMCKpi==1) flifetiMCK->Fill( lenYuri*0.4933667/particle->P() );
720 // if(nMCKpi==1) flifetiMCK->Fill( LengthK*0.4933667/ ptK );
721 } //positive kaon
722 //
723 if ( code < 0) {
724
725 if( nMCKpi == 1) fgenPtEtR->Fill( ptK ); // k to pipi
726 if( nMCKpi == 1) fgenPtEtRN->Fill(ptK); // k to pipi
727 if( nMCKpi == 1) fSignPtEtaMC->Fill(chargPt,rapidiKMC ); // k to pipi
728 if( nMCKpi == 1) fSignPtMC->Fill(chargPt); // k to pipi
729 if( nMCKpi == 1) fMCEtaKaon->Fill(rapidiKMC ); // k to pipi
730 if(nMCKpi==1) fradMC->Fill(radiusD );
731 if(nMCKpi==1) fQtKPiMC->Fill( MCQt );
732 if( nMCKpi== 1) fradPtRapDC->Fill( radiusD, 1./ptK, rapidiKMC); // systematics 26/8
733 if(nMCKpi==1) flifetiMCK->Fill( lenYuri*0.4933667/particle->P() );
734 // if(nMCKpi==1) flifetiMCK->Fill( LengthK*0.4933667/ ptK );
735
736 } //negative K
737
738 } /// kaons loop
739
740
741
742 nMCTracks++;
743 }// end of mc particle
744
745// Phys sel 2012 EFF calculation
746Bool_t isSelected =
747((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected()&AliVEvent::kMB;
748
749 if ( isSelected ==kFALSE) return; // 24/6/11 apo MF
750//
751 fMultiplMC->Fill(nPrim);
752//=======================================================================================
753// main vertex selection
92adf4f6 754 const AliESDVertex *vertex=GetEventVertex(esd);
10eaad41 755 if(!vertex) return;
a0271296 756///
757// fMultiplMC->Fill(nPrim);
10eaad41 758//
a0271296 759/* / apo Alexander Feb 2012
760 AliESDpid *fESDpid = new AliESDpid();
761 if (!fESDpid) fESDpid =
762 ((AliESDInputHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->GetESDpid();
763*/ ///========================================================================
764// apo Eftihi
765 if(!fPIDResponse) {
766 AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
767 AliInputEventHandler* inputHandler =
768(AliInputEventHandler*)(man->GetInputEventHandler());
769 fPIDResponse = inputHandler->GetPIDResponse();
770 }
771
772
773 Double_t vpos[3];
774 vertex->GetXYZ(vpos);
10eaad41 775 fZpr->Fill(vpos[2]);
a0271296 776 if ( TMath::Abs( vpos[2]) > 10. ) return; /// it is applied on ESD and generation
10eaad41 777
778 Double_t vtrack[3], ptrack[3];
779
780
92adf4f6 781 // Int_t nESDTracK = 0;
10eaad41 782
92adf4f6 783 Int_t nGoodTracks = esd->GetNumberOfTracks();
a0271296 784//
10eaad41 785 fESDMult->Fill(nGoodTracks);
92adf4f6 786
10eaad41 787//
a0271296 788 Double_t nsigmall = 100.0;
789 Double_t nsigma = 100.0;
790 Double_t nsigmaPion =-100.0;
791 Double_t nsigmaPi=-100.0;
792 // Double_t dEdxDauMC = 0.0;
793
794//
795for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
796
797 AliESDtrack* trackD = esd->GetTrack(iTrack);
798 if (!trackD) {
799 Printf("ERROR: Could not receive track %d", iTrack);
800 continue;
801 }
802
5220a87d 803 // Int_t indexKinkDau=trackD->GetKinkIndex(0);
a0271296 804// daughter kink
805 nsigmaPion = (fPIDResponse->NumberOfSigmasTPC(trackD , AliPID::kPion));// 26/10 eftihis
806 // nsigmaPion= (fESDpid->NumberOfSigmasTPC(trackD,AliPID::kPion));
807 // 22/11/12 if((indexKinkDau >0)&& (nsigmaPion>1.2)) fTPCSgnlKinkDau->Fill(trackD->P(), (trackD->GetTPCsignal() ) ) ; // daughter kink
808 //if((indexKinkDau >0)) fTPCSgnlKinkDau->Fill(trackD->P(), (trackD->GetTPCsignal() ) ) ; // daughter kink
809// if((indexKinkDau >0)) dEdxDauMC = trackD->GetTPCsignal() ; // daughter kink
810 }
811
812// loop on ESD tracks
813
814//
92adf4f6 815 for (Int_t iTracks = 0; iTracks < esd->GetNumberOfTracks(); iTracks++) {
92adf4f6 816
817 AliESDtrack* track = esd->GetTrack(iTracks);
10eaad41 818 if (!track) {
7ccf0419 819 Printf("ERROR: Could not receive track %d", iTracks);
10eaad41 820 continue;
821 }
822
823
a0271296 824
825// sigmas
826 nsigmall = (fPIDResponse->NumberOfSigmasTPC(track, AliPID::kKaon));
827// nsigma = TMath::Abs(fESDpid->NumberOfSigmasTPC(track,AliPID::kKaon));
828 // nsigmall = (fESDpid->NumberOfSigmasTPC(track,AliPID::kKaon));
829 nsigma = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(track, AliPID::kKaon));
830
831//==================================
832 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
833 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
834 if (track->GetTPCNclsF()>0) {
835 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC / track->GetTPCNclsF();
836 fRatioCrossedRows->Fill(ratioCrossedRowsOverFindableClustersTPC);
837 }
838
10eaad41 839 fHistPt->Fill(track->Pt());
840
92adf4f6 841
a0271296 842 Double_t tpcNCl = track->GetTPCclusters(0); //
843 //Int_t tpcNCl = track->GetTPCclusters(0); //
844 Double_t tpcSign = track->GetSign(); //
845
846 Int_t label = track->GetLabel();
847
92adf4f6 848 label = TMath::Abs(label);
a0271296 849
850 if(label > mcEvent->GetNumberOfTracks()) continue; //
851
92adf4f6 852 TParticle * part = stack->Particle(label);
853 if (!part) continue;
854// loop only on Primary tracks
a0271296 855 if (label > nPrim) continue; /// primary tracks only ,21/3/10 EFF study
92adf4f6 856
a0271296 857 // pt cut
858 if ( (track->Pt())<.200)continue; //12/2/2012 -------------------------------------------
10eaad41 859
a0271296 860 UInt_t status=track->GetStatus();
10eaad41 861
a0271296 862 if((status&AliESDtrack::kITSrefit)==0) continue;
863 if((status&AliESDtrack::kTPCrefit)==0) continue; //6 feb
864 // if((track->GetTPCchi2()/track->GetTPCclusters(0))>3.8) continue;
865 //if((track->GetTPCchi2()/track->GetTPCclusters(0))>4.5) continue; // test 10/3/2012
866 if((track->GetTPCchi2()/track->GetTPCclusters(0))>4.0) continue; // test 10/3/2012
10eaad41 867
868 Double_t extCovPos[15];
869 track->GetExternalCovariance(extCovPos);
a0271296 870//
10eaad41 871
872
873 track->GetXYZ(vtrack);
874 fXvYv->Fill(vtrack[0],vtrack[1]);
92adf4f6 875 fZvYv->Fill(vtrack[0],vtrack[2]);
876 fZvXv->Fill(vtrack[1],vtrack[2]);
10eaad41 877
878// track momentum
879 track->GetPxPyPz(ptrack);
880
881 TVector3 trackMom(ptrack[0],ptrack[1],ptrack[2]);
882
a0271296 883 Double_t etracK= TMath::Sqrt(trackMom.Mag()*trackMom.Mag() + 0.493677 *0.493677 );
884 Double_t rapiditK = 0.5 * (TMath::Log( (etracK + ptrack[2] ) / ( etracK - ptrack[2]) )) ;
10eaad41 885 Double_t trackEta=trackMom.Eta();
a0271296 886 // Double_t trMoment= trackMom.Mag();
887 Double_t trackPt = track->Pt();
10eaad41 888
10eaad41 889
890 Float_t bpos[2];
891 Float_t bCovpos[3];
892 track->GetImpactParameters(bpos,bCovpos);
893
894 if (bCovpos[0]<=0 || bCovpos[2]<=0) {
7ccf0419 895 Printf("Estimated b resolution lower or equal zero!");
10eaad41 896 bCovpos[0]=0; bCovpos[2]=0;
897 }
898
899 Float_t dcaToVertexXYpos = bpos[0];
5220a87d 900// Float_t dcaToVertexZpos = bpos[1];
10eaad41 901
a0271296 902 //fRpr->Fill(dcaToVertexZpos);
903 fRpr->Fill(dcaToVertexXYpos);
904
905 //if((TMath::Abs(dcaToVertexXYpos) >0.3)||(TMath::Abs(dcaToVertexZpos)>2.5)) continue; // 22/7/11
906 // if((TMath::Abs(dcaToVertexXYpos) >0.24)||(TMath::Abs(dcaToVertexZpos)>2.5)) continue; // 12/2/13
907 if (!fMaxDCAtoVtxCut->AcceptTrack(track)) continue;
10eaad41 908
10eaad41 909
910 // cut on eta
a0271296 911 // if( (TMath::Abs(trackEta )) > 0.9 ) continue;
912 if( (TMath::Abs(rapiditK )) > 0.7 ) continue; //// rapid K, Feb 20
10eaad41 913 fHistPtESD->Fill(track->Pt());
914
915 // Add Kink analysis
916
917 Int_t indexKinkPos=track->GetKinkIndex(0);
a0271296 918// loop on mother kinks
10eaad41 919 if(indexKinkPos<0){
920 fPtKink->Fill(track->Pt()); /// pt from track
921
a0271296 922 fRatioCrossedRowsKink->Fill(ratioCrossedRowsOverFindableClustersTPC);
923
10eaad41 924 // select kink class
925
92adf4f6 926 AliESDkink *kink=esd->GetKink(TMath::Abs(indexKinkPos)-1);
10eaad41 927//
928
a0271296 929// DCA kink
930 Double_t Dist2 = kink->GetDistance();
931 // fDCAkink->Fill( Dist2 );
932//
10eaad41 933 Int_t eSDfLabel1=kink->GetLabel(0);
934 TParticle *particle1= stack->Particle(TMath::Abs(eSDfLabel1));
935 Int_t code1= particle1->GetPdgCode();
a0271296 936 // Int_t mcProcssMo= particle1->GetUniqueID();
10eaad41 937
938 Int_t eSDfLabeld=kink->GetLabel(1);
939 TParticle *particled= stack->Particle(TMath::Abs(eSDfLabeld));
940 Int_t dcode1= particled->GetPdgCode();
a0271296 941 Int_t mcProcssDa= particled->GetUniqueID();
942//
943// loop on MC daugthres for 3Pi 24/9/2010
944 Int_t nESDKpi =0;
945 if(mcProcssDa==4) {
946 Int_t firstDa=particle1->GetFirstDaughter();
947 Int_t lastDa=particle1->GetLastDaughter();
948
949 if( (lastDa<=0) || (firstDa<=0) ) continue;
950
951 if ( lastDa > mcEvent->GetNumberOfTracks() ) continue;
952 if (firstDa > mcEvent->GetNumberOfTracks() ) continue;
953//loop on secondaries
954 for (Int_t kk=firstDa;kk<=lastDa;kk++) {
955 if ( kk > 0 ) {
956 TParticle* daughter3=stack->Particle(kk); // 24/9
957 Float_t dcode3= daughter3->GetPdgCode();
958 if (( ( code1==321 )&& ( dcode3==211 ))|| (( code1 == -321 )&& ( dcode3==-211))) nESDKpi++ ;
959 }
960 }
961 }
962//---------------------------edw telos 9/2010
963 Double_t hVzdau=particled->Vz();
10eaad41 964
965 const TVector3 motherMfromKink(kink->GetMotherP());
966 const TVector3 daughterMKink(kink->GetDaughterP());
967 Float_t qT=kink->GetQt();
a0271296 968 // Float_t motherPt=motherMfromKink.Pt();
969// Kink mother momentum
5220a87d 970// Double_t trMomTPCKink=motherMfromKink.Mag();
a0271296 971// TPC mother momentun
5220a87d 972 // Double_t trMomTPC=track->GetTPCmomentum();
a0271296 973 // Float_t etaMother=motherMfromKink.Eta();
974
10eaad41 975
10eaad41 976 fHistQtAll->Fill(qT) ; // Qt distr
a0271296 977
978 const TVector3 vposKink(kink->GetPosition());
979 fPosiKink ->Fill( vposKink[0], vposKink[1] );
980
981 Double_t dxKink = vpos[0]-vposKink[0], dyKink=vpos[1]-vposKink[1], dzKink=vpos[2]-vposKink[2] ;
982 // Double_t dzKink=vpos[2]-vposKink[2] ; /// ??
983 Double_t lifeKink= TMath::Sqrt( dxKink*dxKink + dyKink*dyKink + dzKink*dzKink ) ;
984
985 Double_t tanLamda = track-> GetTgl() ;// ??
986 if (tanLamda ==0 ) continue;// ??
987 Double_t lenHelx = (TMath::Abs(dzKink ) ) *(TMath::Sqrt( 1. + tanLamda *tanLamda ) ) / ( TMath::Abs( tanLamda)) ;// ??
988
989
990
10eaad41 991
92adf4f6 992 Float_t kinkAngle=TMath::RadToDeg()*kink->GetAngle(2);
993
10eaad41 994// fake kinks, small Qt and small kink angle
a0271296 995 if(( kinkAngle<1.)) fHistQt1 ->Fill(qT) ; // Qt distr
996//
997 if( ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==321))) fFakeKPi->Fill(track->Pt());
998 if( ( (TMath::Abs(code1)==211)&&(TMath::Abs(dcode1)==211))) fFakepipi->Fill(track->Pt());
999//
10eaad41 1000
10eaad41 1001//remove the double taracks
a0271296 1002 if( (kinkAngle<2.) ) continue; // test 15/7 2010 , it removes 3 from 10000 good Kaons
1003 // BG ?????==============
1004 if ( TMath::Abs(vposKink[2]) > 225. ) continue ;
1005 if ( TMath::Abs(vposKink[2]) < 0.5 ) continue ;
7ccf0419 1006
a0271296 1007 fPtCut1 ->Fill(trackPt );
1008 fChi2NclTPC->Fill( (track->GetTPCchi2() ) , tpcNCl );
1009
1010 Float_t signPt= tpcSign*trackPt;
10eaad41 1011//
5220a87d 1012 if ( (TMath::Abs(code1)==211)&&(TMath::Abs(dcode1)==13))
1013 fRadiusPtPion->Fill( kink->GetR(), track->Pt()); //
10eaad41 1014
5220a87d 1015 if( ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==13))||
1016 ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==11)) ||
1017 ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==211)) ) {
1018 fRadiusPtKaon->Fill( kink->GetR(), track->Pt()); //
1019 }
1020//
a0271296 1021 // ======8/1/13 if((kink->GetR()>120.)&&(kink->GetR()<210.)&&(TMath::Abs(rapiditK)<0.7)&&(label<nPrim)) {
1022 if((kink->GetR()> fKinkRadLow )&&(kink->GetR()< fKinkRadUp )&&(TMath::Abs(rapiditK)<0.7)&&(label<nPrim)) {
1023 if( ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==13))) fQtKMu->Fill(qT);
1024 if ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==11)) fQtKEl->Fill(qT);
1025 if ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==211)) fQtKPi->Fill(qT);
1026 if (( nESDKpi>1) && ( ((code1)==321)&&((dcode1)==211)) ) fQtK3PiP->Fill(qT);
1027 if (( nESDKpi>1) && ( ((code1)==-321)&&((dcode1)==-211)) ) fQtK3PiM->Fill(qT);
92adf4f6 1028 if( ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==13))||
1029 ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==11)) ||
1030 ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==211)) ) {
5220a87d 1031 if(qT>fLowQt) fHistPtKPDG->Fill(track->Pt()); // ALL KAONS (pdg) inside ESD kink sample
1032 if(qT>fLowQt ) {
1033 if(code1>0.) fHiPtKPDGP->Fill(trackPt ); // 26/feb // ALL KAONS (pdg) inside ESD kink sample
1034 if(code1<0.) fHiPtKPDGN->Fill( trackPt ); // 26/feb // ALL KAONS (pdg) inside ESD kink sample
a0271296 1035 }
1036 fHistEta->Fill(trackEta) ; // Eta distr of PDG kink ESD kaons
1037 frapidESDK->Fill(rapiditK) ; //18/feb rapiddistr of PDG kink ESD kaons
5220a87d 1038 if( qT > fLowQt ) fHistQt2->Fill(qT); // PDG ESD kaons
1039 fRadiusPt->Fill( kink->GetR(), track->Pt()); //
92adf4f6 1040 }
1041 }
10eaad41 1042
a0271296 1043
1044 //Double_t maxDecAngKmu=f1->Eval(motherMfromKink.Mag(),0.,0.,0.);
1045 Double_t maxDecAngKmu=f1->Eval(track->P(), 0.,0.,0.);
1046 Double_t maxDecAngpimu=f2->Eval(track->P(), 0.,0.,0.);
10eaad41 1047// two dimensional plot
a0271296 1048 if(TMath::Abs(code1)==321) fAngMomK->Fill(track->P(), kinkAngle);
1049 if(TMath::Abs(code1)==211) fAngMomPi->Fill( track->P(), kinkAngle);
1050//_______
1051
10eaad41 1052//
1053// invariant mass of mother track decaying to mu
1054 Float_t energyDaughterMu=TMath::Sqrt(daughterMKink.Mag()*daughterMKink.Mag()+0.105658*0.105658);
a0271296 1055 Float_t energyDaughterPi=TMath::Sqrt(daughterMKink.Mag()*daughterMKink.Mag()+0.139570*0.139570);
1056 Float_t energyDaughterKa=TMath::Sqrt(daughterMKink.Mag()*daughterMKink.Mag()+0.493677*0.493677);
1057// Float_t energyDaughterPr=TMath::Sqrt(daughterMKink.Mag()*daughterMKink.Mag()+0.938658*0.938658);
10eaad41 1058 Float_t p1XM= motherMfromKink.Px();
1059 Float_t p1YM= motherMfromKink.Py();
1060 Float_t p1ZM= motherMfromKink.Pz();
1061 Float_t p2XM= daughterMKink.Px();
1062 Float_t p2YM= daughterMKink.Py();
1063 Float_t p2ZM= daughterMKink.Pz();
1064 Float_t p3Daughter=TMath::Sqrt(((p1XM-p2XM)*(p1XM-p2XM))+((p1YM-p2YM)*(p1YM-p2YM))+((p1ZM-p2ZM)*(p1ZM-p2ZM)));
1065 Double_t invariantMassKmu= TMath::Sqrt((energyDaughterMu+p3Daughter)*(energyDaughterMu+p3Daughter)-motherMfromKink.Mag()*motherMfromKink.Mag());
a0271296 1066 Double_t invariantMassKpi= TMath::Sqrt((energyDaughterPi+p3Daughter)*(energyDaughterPi+p3Daughter)-motherMfromKink.Mag()*motherMfromKink.Mag());
1067 Double_t invariantMassKK = TMath::Sqrt((energyDaughterKa+p3Daughter)*(energyDaughterKa+p3Daughter)-motherMfromKink.Mag()*motherMfromKink.Mag());
1068 fInvMassMuNuAll ->Fill(invariantMassKmu);
fda107c4 1069 fInvMassMuNuPtAll ->Fill(invariantMassKmu, trackPt);
5220a87d 1070 // 20/4 testRadiusPt->Fill( kink->GetR(), track->Pt()); //
10eaad41 1071
a0271296 1072
5220a87d 1073 if (qT> fLowQt ) fSignPtNcl->Fill( signPt , tpcNCl );
a0271296 1074
1075//
1076 // if((qT>0.12)&&((kink->GetR()>120.)&&(kink->GetR()<210.))&&(TMath::Abs(rapiditK )<0.7)) {
5220a87d 1077 if((qT> fLowQt )&&((kink->GetR()> fKinkRadLow )&&(kink->GetR()< fKinkRadUp ))&&(TMath::Abs(rapiditK )<0.7)) {
10eaad41 1078 fM1kaon->Fill(invariantMassKmu);
a0271296 1079 fMinvPi->Fill(invariantMassKpi);
1080 fMinvKa->Fill(invariantMassKK);
1081 fRadiusNcl->Fill( (kink->GetR()) ,(track->GetTPCclusters(0) ) ) ;
1082 }
1083
1084//
1085 // if ( tpcNCl<30 ) continue;
5220a87d 1086 // test for systematics , march 13 if ( tpcNCl<50. ) continue;
1087 if ( tpcNCl<fLowCluster ) continue;
1088 // if ( ( tpcNCl<20. )|| ( tpcNCl > 100.) ) continue;// test , den edwse kati shmantiko
1089
a0271296 1090 //if( ( ( track->GetTPCclusters(0) ) / (kink->GetR() ) ) > 0.63 ) continue;
1091 Double_t tpcNClHigh = -51.67+ (11./12.) *( kink->GetR() ) ;
5220a87d 1092 if (tpcNCl > tpcNClHigh ) fcodeH->Fill(TMath::Abs(code1), TMath::Abs(dcode1));
1093 if ( tpcNCl >tpcNClHigh) fZkinkZDau->Fill( vposKink[2],hVzdau );
1094 if (tpcNCl > tpcNClHigh) fRadiusPtFake->Fill( kink->GetR(), track->Pt()); //
a0271296 1095
1096 Double_t tpcNClMin = -85.5 + (65./95.) *( kink->GetR() ) ;
1097 // if ( tpcNClMin < tpcNCl ) continue;
5220a87d 1098 if ( tpcNCl > tpcNClHigh) continue;
a0271296 1099 if ( tpcNCl < tpcNClMin ) continue;
a0271296 1100 //
10eaad41 1101
1102// kaon selection from kinks
a0271296 1103
1104 //===== 8/2/13 if((kinkAngle>maxDecAngpimu)&&(qT>0.12)&&(qT<0.30)&&((kink->GetR()>120.)&&(kink->GetR()<210.))&&(TMath::Abs(rapiditK )<0.7)&&(invariantMassKmu<0.8)) {
5220a87d 1105 if((kinkAngle>maxDecAngpimu)&&(qT>fLowQt )&&(qT<0.30)&&((kink->GetR()> fKinkRadLow )&&(kink->GetR()< fKinkRadUp ))&&(TMath::Abs(rapiditK )<0.7)&&(invariantMassKmu<0.8)) {
a0271296 1106// 29092010 if((kinkAngle>maxDecAngpimu)&&(qT>0.120)&&(qT<0.25)&&((kink->GetR()>120.)&&(kink->GetR()<210.))&&(TMath::Abs(rapiditK )<0.7)&&(invariantMassKmu<0.6)) {
1107// if((kinkAngle>maxDecAngpimu)&&(qT>0.04)&&(qT<0.30)&&((kink->GetR()>133.)&&(kink->GetR()<179.))&&(TMath::Abs(rapiditK )<0.5)&&(invariantMassKmu<0.6)) {
1108
1109 fAngMomKKinks->Fill(track->P(), kinkAngle);
1110 fPtCut2 ->Fill(trackPt );
1111
1112 if ( (kinkAngle>maxDecAngKmu*0.98)&& ( track->P() > 1.2 ) ) continue; // maximum angle s
1113 if ( (kinkAngle<maxDecAngpimu*1.20) ) continue; // maximum angle s
1114
1115 fPtCut3 ->Fill(trackPt );
1116 //fTPCSgnlPa->Fill(track->P(),track->GetTPCsignal());
1117 fTPCSgnlPa->Fill(track->GetInnerParam()->GetP(),track->GetTPCsignal());
1118 // if( nsigma > 3.5 ) fcode2->Fill(TMath::Abs(code1), TMath::Abs(dcode1));
5220a87d 1119 if(nsigma > 3.5) continue; // 1/11/12
1120 // test 16/2/13 if(nsigma > 3.5) continue; // 15/2/13
a0271296 1121 // if(nsigma > 4.0) continue; // test 17/2/2011 4% or more ? bg?
1122//
1123 fTPCSgnlP->Fill(track->GetInnerParam()->GetP(), (track->GetTPCsignal() ) ) ;
1124
1125 fInvMassMuNuPt ->Fill(invariantMassKmu, trackPt);
1126// loop on kink daughters inside mother's loop
1127 Int_t ESDLabelM = 0. ;
1128 Int_t ESDLabelD = 0. ;
1129 Double_t dEdxDauMC = 0.0;
1130 Double_t raDAU=0.;
1131 Int_t Ikink =0;
1132 Int_t IRkink =0;
1133for (Int_t jTrack = 0; jTrack < esd->GetNumberOfTracks(); jTrack++) {
1134
1135 AliESDtrack* trackDau = esd->GetTrack(jTrack);
1136 if (!trackDau) {
1137 Printf("ERROR: Could not receive track %d", jTrack);
1138 continue;
1139 }
1140 Int_t indexKinkDAU =trackDau->GetKinkIndex(0);
1141 if (indexKinkDAU <0 ){
1142 AliESDkink *kinkDau=esd->GetKink(TMath::Abs(indexKinkDAU)-1);
1143 raDAU= kinkDau->GetR();
1144 ESDLabelM=kinkDau->GetLabel(0); // mothers's label
1145 ESDLabelM = TMath::Abs(ESDLabelM);
1146 ESDLabelD=kinkDau->GetLabel(1); // Daughter'slabel
1147 ESDLabelD = TMath::Abs(ESDLabelD);
1148 if ( kink->GetR() == kinkDau->GetR() ) IRkink++;
1149 if ( ESDLabelM == label ) Ikink++ ;
1150 }
1151 // if (( ESDLabelM== eSDfLabel1)) {
1152 if ( (Ikink >0) && (IRkink>0 ) ) {
1153// daughter kink
1154 //if(indexKinkDAU >0) nsigmaPi = (fPIDResponse->NumberOfSigmasTPC(trackDau,AliPID::kPion));// 26/10 eftihis
1155 if(indexKinkDAU >0) nsigmaPi = (fPIDResponse->NumberOfSigmasTPC(trackDau,AliPID::kKaon));// 26/10 eftihis
1156 // nsigmaPion= (fESDpid->NumberOfSigmasTPC(trackD,AliPID::kPion));
1157 // 22/11/12 if((indexKinkDau >0)&& (nsigmaPion>1.2)) fTPCSgnlKinkDau->Fill(trackD->P(), (trackD->GetTPCsignal() ) ) ; // daughter kink
1158 //if((indexKinkDau >0)) fTPCSgnlKinkDau->Fill(trackD->P(), (trackD->GetTPCsignal() ) ) ; // daughter kink
1159 if((indexKinkDAU >0)) dEdxDauMC = trackDau->GetTPCsignal() ; // daughter kink
1160 }
1161 }
1162// end internal loop for kink daughters
1163
1164 // fTPCSgnlP->Fill(trMomTPC , (track->GetTPCsignal() ) ) ;
1165 // fTPCSgnlPtpc->Fill(trMomTPC , (track->GetTPCsignal() ) ) ;
1166 //fMothKinkMomSgnl ->Fill(trMomTPCKink , (track->GetTPCsignal() ) ) ;
1167 // fMothKinkMomSgnl ->Fill( dEdxDauMC , (track->GetTPCsignal() ) ) ;
1168 // fTPCMomNSgnl->Fill(trMomTPC ,pidResponse->NumberOfSigmasTPC(track, AliPID::kKaon) );
1169 fNSigmTPC ->Fill(nsigmall );
1170// daughter selection
1171 //fTPCSgnlKinkDau->Fill( daughterMKink.Mag() , dEdxDauMC ) ; // daughter kink
1172 // if( nsigmaPion > 1.0 ) fcode4->Fill(TMath::Abs(code1), TMath::Abs(dcode1));
1173 // if( dEdxDauMC > 1.5 *(track->GetTPCsignal() ) ) fcode4->Fill(TMath::Abs(code1), TMath::Abs(dcode1));
1174//
1175 //fTPCMomNSgnl->Fill(trMomTPC ,nsigmall );
1176 fTPCMomNSgnl->Fill(track->GetInnerParam()->GetP() ,nsigmall );
1177//
7ccf0419 1178
a0271296 1179 fRadNclcln->Fill( (kink->GetR()) ,(track->GetTPCclusters(0) ) ) ;
1180 fRadiusPtcln->Fill( kink->GetR(), trackPt); //
10eaad41 1181
a0271296 1182 fAngMomKC->Fill(track->P() , kinkAngle);
1183 fHistPtKaon->Fill(track->Pt() ); //all PID kink-kaon
1184 if( code1>0.) fHistPtKaoP->Fill(track->Pt() ); //all PID kink-kaon
1185 if( code1<0.) fHistPtKaoN->Fill(track->Pt() ); //all PID kink-kaon
1186// systematics
1187 fradPtRapESD->Fill(kink->GetR(), 1./ track->Pt(), rapiditK);
1188//
1189 fHistEtaK->Fill(rapiditK );
1190 frad->Fill( kink->GetR() );
1191 flenHelx->Fill( lenHelx ); //??
1192 flifeKink ->Fill(lifeKink );//??
1193 fLHelESDK ->Fill( ( lenHelx /track->P() )*0.493677);// for all 'PID' kaons 31/7/11// ??
1194 flifTiESDK->Fill( ( lifeKink /track->P() )*0.493677); // ??
10eaad41 1195
10eaad41 1196
a0271296 1197
1198 fSignPtNcl->Fill( signPt , tpcNCl );
1199 fSignPtEta->Fill(signPt , rapiditK );
1200 fEtaNcl->Fill(rapiditK, tpcNCl );
1201 fSignPt->Fill( signPt );
1202 fRatChi2Ncl-> Fill((track->GetTPCchi2()/track->GetTPCclusters(0) )) ;
1203 fdcatoVxXY->Fill(dcaToVertexXYpos);
10eaad41 1204
92adf4f6 1205// kaons from k to mun and k to pipi and to e decay
1206 if( ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==13))||
1207 ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==11)) ||
1208 ( (TMath::Abs(code1)==321)&&(TMath::Abs(dcode1)==211)) ) {
10eaad41 1209
7ccf0419 1210 if ( label<=nPrim ) fPtPrKink->Fill(track->Pt());
a0271296 1211
1212 // flifetim2 ->Fill( ( lenHelx /track->P() )*0.493667); // to compare with fLHelESDK
1213 flifTiESDK->Fill( ( lifeKink /track->P() )*0.493667);
1214
1215
1216
1217 fKinkKaon->Fill(track->Pt());
1218 fDCAkink->Fill( Dist2 );
1219 fPosiKinkK->Fill( vposKink[0], vposKink[1] );
1220 fPosiKinKXZ->Fill( vposKink[0], vposKink[2] );
1221 fPosiKinKYZ->Fill( vposKink[1], vposKink[2] );
1222 if( code1>0.) fkinkKaonP->Fill(trackPt); // kPtPID kink-kaon
1223 if( code1<0.) fkinkKaonN->Fill(trackPt); // PID kink-kaon
1224// daughters
1225 if((((nsigmaPi) > 0.)&& ( dEdxDauMC > 70. ) )) fcode4->Fill(TMath::Abs(code1), TMath::Abs(dcode1));
1226 if ( TMath::Abs(dEdxDauMC - track->GetTPCsignal() ) < 2) fcode2->Fill( TMath::Abs(code1), TMath::Abs(dcode1));
1227 // fTPCSgnlPtpc->Fill(trMomTPC , (track->GetTPCsignal() ) ) ;
1228 fTPCSgnlPtpc ->Fill( daughterMKink.Mag() , dEdxDauMC ) ; // daughter kink
1229 fMothKinkMomSgnlD->Fill( dEdxDauMC , (track->GetTPCsignal() ) ) ;
92adf4f6 1230 }
10eaad41 1231 else {
a0271296 1232 fKinkKaonBg->Fill(track->Pt());
1233 fMothKinkMomSgnl ->Fill( dEdxDauMC , (track->GetTPCsignal() ) ) ;
1234// fTPCSgnlKinkDau->Fill( daughterMKink.Mag() , dEdxDauMC ) ; // daughter kink
1235 // if( dEdxDauMC > 1.5 *(track->GetTPCsignal() ) ) fcode4->Fill(TMath::Abs(code1), TMath::Abs(dcode1));
1236 // if(( (track->P())<1. )&& ( dEdxDauMC > 1.5 *(track->GetTPCsignal() ) )) fcodeDau1 ->Fill(TMath::Abs(code1), TMath::Abs(dcode1));
1237 if( (( nsigmaPi) > 0. ) && (( dEdxDauMC > 70 ) )) fcodeDau1 ->Fill(TMath::Abs(code1), TMath::Abs(dcode1));
1238 if ( TMath::Abs(dEdxDauMC - track->GetTPCsignal() ) < 2) fcodeDau2->Fill( TMath::Abs(code1), TMath::Abs(dcode1));
1239 fTPCSgnlKinkDau->Fill( daughterMKink.Mag() , dEdxDauMC ) ; // daughter kink
1240//
1241 fMinvPr->Fill(invariantMassKmu);
1242//
1243 fDCAkinkBG->Fill( Dist2 );
1244 fPosiKinKBgXY->Fill( vposKink[0], vposKink[1] );
1245 fPosiKinKBgZY->Fill( vposKink[2], vposKink[1] );
1246 fPosiKinKBgZX->Fill( vposKink[2], kink->GetR() ); // 31/7/11
1247 if( code1>0.) fKinKBGP ->Fill( trackPt ); //all PID kink-kaon
1248 if( code1<0.) fKinKBGN ->Fill( trackPt ); //all PID kink-kaonl
7ccf0419 1249 fdcodeH->Fill( TMath::Abs(code1), TMath::Abs(dcode1)); // put it here, 22/10/2009
5220a87d 1250 // if (eSDfLabel1==eSDfLabeld) fcodeH->Fill(TMath::Abs(code1), TMath::Abs(dcode1));
1251 // if (eSDfLabeld>nPrim ) fZkinkZDau->Fill( vposKink[2],hVzdau );
a0271296 1252
7ccf0419 1253 } // primary and all +BG
10eaad41 1254
1255 } // kink selection
1256
1257
1258 } //End Kink Information
1259
1260
1261 } //track loop
1262
a0271296 1263// } // vx 10 cm only on esd
10eaad41 1264 PostData(1, fListOfHistos);
1265
1266}
1267
1268//________________________________________________________________________
1269void AliAnalysisKinkESDMC::Terminate(Option_t *)
1270{
1271 // Draw result to the screen
1272 // Called once at the end of the query
1273
1274}
10eaad41 1275
1276const AliESDVertex* AliAnalysisKinkESDMC::GetEventVertex(const AliESDEvent* esd) const
1277 // Get the vertex from the ESD and returns it if the vertex is valid
1278
1279{
1280 // Get the vertex
1281
a0271296 1282// 10/4 const AliESDVertex* vertex = esd->GetPrimaryVertex();
1283 const AliESDVertex* vertex = esd->GetPrimaryVertexTracks();
10eaad41 1284
a0271296 1285 //if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>2)) return vertex;
1286 if((vertex->GetStatus()==kTRUE)) return vertex;
10eaad41 1287 else
1288 {
1289 vertex = esd->GetPrimaryVertexSPD();
a0271296 1290 if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>0)) return vertex;
1291 // if((vertex->GetStatus()==kTRUE)) return vertex;
10eaad41 1292 else
1293 return 0;
1294 }
f92b626a 1295}