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