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