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