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