]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/SPECTRA/Kinks/AliAnalysisKinkESDat.cxx
updated kink task for dcaXY
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / Kinks / AliAnalysisKinkESDat.cxx
CommitLineData
828a6e0c 1/**************************************************************************
fcd98e58 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 * The method is applied in pp and Pb-Pb real data. *
828a6e0c 6 * *
828a6e0c 7 **************************************************************************/
8
9//-----------------------------------------------------------------
10// AliAnalysisKinkESDat class
11// Example of an analysis task for kink topology study
12// Kaons from kink topology are 'identified' in this code
13//-----------------------------------------------------------------
14
15#include "TChain.h"
16#include "TTree.h"
17#include "TH1F.h"
18#include "TH2F.h"
19#include "TH3F.h"
20#include "TH1D.h"
21#include "TH2D.h"
22#include "TParticle.h"
23#include <TVector3.h>
24#include "TF1.h"
25
26#include "AliAnalysisTask.h"
27#include "AliAnalysisManager.h"
28
29#include "AliVEvent.h"
30#include "AliESDEvent.h"
31#include "AliMCEvent.h"
32#include "AliAnalysisKinkESDat.h"
33#include "AliStack.h"
34#include "AliESDpid.h"
35#include "AliPID.h"
36#include "AliESDkink.h"
37#include "AliESDtrack.h"
38#include "AliPhysicsSelectionTask.h"
39#include "AliInputEventHandler.h"
ad1036b2 40#include "AliESDInputHandler.h"
828a6e0c 41#include "AliAnalysisManager.h"
42#include "AliVEvent.h"
828a6e0c 43#include "AliESDtrackCuts.h"
fcd98e58 44#include "AliPIDResponse.h"
8e57b720 45/////////#include "AliTENDERSupplies.h"
828a6e0c 46ClassImp(AliAnalysisKinkESDat)
47
48
49//________________________________________________________________________
50AliAnalysisKinkESDat::AliAnalysisKinkESDat(const char *name)
51 : AliAnalysisTaskSE(name), fHistPtESD(0),fHistPt(0),fHistQtAll(0),fHistQt1(0),fHistQt2(0)
52 , fHistPtKaon(0),fHistPtKPDG(0),fHistEta(0),fHistEtaK(0),fptKMC(0),fMultiplMC(0),fESDMult(0),fgenpt(0),frad(0),
ad1036b2 53 fKinkKaon(0),fKinKRbn(0), fKinkKaonBg(0), fM1kaon(0), fPtKink(0), fptKink(0),
54 fAngMomK(0),fAngMomPi(0), fAngMomKC(0), fMultESDK(0), fMultMCK(0),
828a6e0c 55 fSignPtNcl(0), fSignPtEta(0), fEtaNcl(0), fSignPt(0), fChi2NclTPC(0), fRatChi2Ncl(0), fRadiusNcl(0), fTPCSgnlP(0),
56 fTPCSgnlPa(0), fRpr(0),fZpr(0), fdcatoVxXY(0), fnSigmToVx(0), fKinkMothDau(0),
ad1036b2 57 fZvXv(0),fZvYv(0), fXvYv(0), fHistPtKaoP(0), fHistPtKaoN(0),frapiKESD(0), flifetime(), fradLK(0),
828a6e0c 58 fradPtRpDt(0), fInvMuNuAll(0), fQtInvM(0),
ad1036b2 59 fDCAkink(0), fPosiKink(0), fPosiKinkK(0),fPosiKinKXZ(0), fPosiKinKYZ(0), fPosiKinKBg(0), fQtMothP(0), fTPCSgnlPtpc(0),
60 fTPCMomNSgnl(0), fMothKinkMomSgnl(0), fNSigmTPC(0), fTPCSgnlKinkDau(0), fPtKinkPos(0), fPtKinkNeg(0), fRadNclCln(0),
8e57b720 61 fRatioCrossedRows(0), fRatioCrossedRowsKink(0),fRadiusPt(0), fRadiusPtcln(0), fInvMassMuNuPt(0), fPtCut1(0), fPtCut2(0), fPtCut3(0),
62 fAngMomKKinks(0),
828a6e0c 63 f1(0), f2(0),
8e57b720 64 fListOfHistos(0),fLowMulcut(-1),fUpMulcut(-1),fCutsMul(0), fMaxDCAtoVtxCut(0), fPIDResponse(0)
828a6e0c 65{
66 // Constructor
67
68 // Define input and output slots here
69 // Input slot #0 works with a TChain
70 // DefineInput(0, TChain::Class());
71//----------------------Marek multiplicity bins
72 fCutsMul=new AliESDtrackCuts("Mul","Mul");
73 fCutsMul->SetMinNClustersTPC(70);
74 fCutsMul->SetMaxChi2PerClusterTPC(4);
75 fCutsMul->SetAcceptKinkDaughters(kFALSE);
76 fCutsMul->SetRequireTPCRefit(kTRUE);
77 // ITS
78 fCutsMul->SetRequireITSRefit(kTRUE);
79 fCutsMul->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
80 AliESDtrackCuts::kAny);
81 fCutsMul->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
82
83 fCutsMul->SetMaxDCAToVertexZ(2);
84 fCutsMul->SetDCAToVertex2D(kFALSE);
85 fCutsMul->SetRequireSigmaToVertex(kFALSE);
86
87 fCutsMul->SetEtaRange(-0.8,+0.8);
88 fCutsMul->SetPtRange(0.15, 1e10);
89
8e57b720 90 fMaxDCAtoVtxCut=new AliESDtrackCuts("fMaxDCA", "fMaxDCA");
91 fMaxDCAtoVtxCut->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
92 fMaxDCAtoVtxCut->SetMaxChi2TPCConstrainedGlobal(36);
93// esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
94 // esdTrackCuts->SetMaxChi2TPCConstrainedGlobal(36);
95
96
828a6e0c 97
98 DefineOutput(1, TList::Class());
99}
100
101//________________________________________________________________________
102void AliAnalysisKinkESDat::UserCreateOutputObjects()
103{
104 // Create histograms
105 // Called once
106
107 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);
108 f1->SetParameter(0,0.493677);
109 f1->SetParameter(1,0.9127037);
110 f1->SetParameter(2,TMath::Pi());
111
112
113 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);
114 f2->SetParameter(0,0.13957018);
115 f2->SetParameter(1,0.2731374);
116 f2->SetParameter(2,TMath::Pi());
117 //Open file 1= CAF
118 //OpenFile(1);
ad1036b2 119 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,
120 1.1, 1.2, 1.3, 1.4, 1.5, 1.6,1.7,1.8,1.9, 2.0,
121 2.2, 2.4, 2.6, 2.8, 3.0, 3.3, 3.6, 3.9,
fcd98e58 122 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
ad1036b2 123
124 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,
828a6e0c 125 1.1, 1.2, 1.3, 1.4, 1.5, 1.6,1.7,1.8,1.9, 2.0,
126 2.1, 2.2, 2.3, 2.4, 2.5, 2.6, 2.7,2.8, 2.9, 3.0,
ad1036b2 127 3.2, 3.4, 3.6, 3.8, 4.0, 4.2, 4.4, 4.6, 4.8,5.0 }; // Barbara TOF Kch
828a6e0c 128
ad1036b2 129 fHistPtESD = new TH1F("fHistPtESD", "P_{T} distribution",300, 0.0,15.0);
828a6e0c 130 fHistPtESD->GetXaxis()->SetTitle("P_{T} (GeV/c)");
131 fHistPtESD->GetYaxis()->SetTitle("dN/dP_{T} (c/GeV)");
132 fHistPtESD->SetMarkerStyle(kFullCircle);
ad1036b2 133 fHistPt = new TH1F("fHistPt", "P_{T} distribution",300, 0.0,15.0);
828a6e0c 134 fHistQtAll = new TH1F("fHistQtAll", "Q_{T} distr All Kinks ",100, 0.0,.300);
135 fHistQt1= new TH1F("fHistQt1", "Q_{T} distribution",100, 0.0,.300);
136 fHistQt2= new TH1F("fHistQt2", "Q_{T} distribution",100, 0.0,.300);
ad1036b2 137 fHistPtKaon = new TH1F("fHistPtKaon", "P_{T}Kaon distribution",300, 0.0,15.0);
138 fHistPtKPDG = new TH1F("fHistPtKPDG", "P_{T}Kaon distribution",300, 0.0,15.0);
828a6e0c 139 fHistEta= new TH1F("fHistEta", "Eta distribution", 26,-1.3, 1.3);
140 fHistEtaK= new TH1F("fHistEtaK", "EtaK distribution", 26,-1.3, 1.3);
ad1036b2 141 fptKMC= new TH1F("fptKMC", "P_{T}Kaon generated",300, 0.0,15.0);
828a6e0c 142 fMultiplMC= new TH1F("fMultiplMC", "charge multiplicity MC",100, 0.0,300.0);
143 fESDMult= new TH1F("fESDMult", "charge multipliESD",100, 0.0,300.0);
ad1036b2 144 fgenpt= new TH1F("fgenpt", "genpt K distribution",300, 0.0,15.0);
828a6e0c 145 //frad= new TH1F("frad", "radius K generated",100, 50., 250.0);
146 frad= new TH1F("frad", "radius K generated",100, 0.,1000.0);
ad1036b2 147 fKinkKaon= new TH1F("fKinkKaon", "P_{T}Kaon kinks identi",300, 0.0,15.0);
148 fKinKRbn= new TH1F("fKinKRbn", "p_{t}Kaon kinks identi[GeV/c],Entries",46,gPt7TOF);
149 fKinkKaonBg= new TH1F("fKinkKaonBg", "P_{T}Kaon kinks backgr",300, 0.0,15.0);
8e57b720 150 fM1kaon= new TH1F("fM1kaon","Invar m(kaon) from kink->mu+netrino decay",180,0.10, 1.0);
ad1036b2 151 fPtKink= new TH1F("fPtKink", "P_{T}Kaon Kink distribution, counts",44, gPt7K0);
152 fptKink= new TH1F("fptKink", "P_{T}Kaon Kink bution",300, 0.0,15.0);
828a6e0c 153 fAngMomK= new TH2F("fAngMomK","Decay angle vrs Mother Mom,K",100,0.0,5.0,80,0.,80.);
154 fAngMomPi= new TH2F("fAngMomPi","Decay angle vrs Mother Mom,Pi",100,0.0,5.0,80,0.,80.);
155 fAngMomKC= new TH2F("fAngMomKC","Decay angle vrs Mother Mom,K",100,0.0,5.0,80,0.,80.);
156 fMultESDK=new TH1F("fMultESDK", "charge multipliESD kaons",100, 0.0,100.0);
157 fMultMCK=new TH1F("fMultMCK", "charge multipli MC kaons",100, 0.0,100.0);
158 fSignPtNcl= new TH2F("fSignPtNcl","SignPt vrs Ncl,K",80,-4.,4.0,70,20.,160.);
159 fSignPtEta= new TH2F("fSignPtEta","SignPt vrs Eta,K",80,-4.0,4.0,30,-1.5,1.5);
160 fEtaNcl= new TH2F("fEtaNcl","Eta vrs nclust,K",30,-1.5,1.5, 70,20, 160);
828a6e0c 161 fSignPt= new TH1F("fSignPt","SignPt ,K",80,-4.0,4.0);
162 fChi2NclTPC= new TH2F("fChi2NclTPC","Chi2vrs nclust,K",100,0.,500., 70,20, 160);
163 fRatChi2Ncl= new TH1F("fRatChi2Ncl","Ratio chi2/nclusters in TPC,K",50,0.0,5.0);
ad1036b2 164 fRadiusNcl = new TH2F("fRadiusNcl","kink radius vrs Nclust,K",75,100.,250., 80,0, 160);
8e57b720 165 fTPCSgnlP = new TH2F("fTPCSgnlP","TPC signal de/dx Mom,K",1000,0.0,20.0,150,0.,300.);
166 fTPCSgnlPa= new TH2F("fTPCSgnlPa","TPC signal de/dx Mom,K",1000,0.0,20.,150, 0.,300.);
828a6e0c 167 fRpr = new TH1D("fRpr", "rad distribution PID pr",100,-10.0, 10.0);
168 fZpr = new TH1D("fZpr", "z distribution PID pr ",80,-20.,20.);
169 fdcatoVxXY = new TH1D("fdcatoVxXY", "dca distribution PID ",20,-1.,1.);
170 fnSigmToVx = new TH1D("fnSigmToVx", "dca distribution PID ",80,0.,8.);
171 fKinkMothDau= new TH2F("fKinkMothDau","TPC kink Moth Daugh ,K",50,0.0,2.5,50, 0., 2.5);
172 fZvXv= new TH2F("fZvXv","Xv-Zv main vtx",60,-0.5,0.5,60, -15., 15.0);
173 fZvYv= new TH2F("fZvYv","Yv-Zv main vtx",60,-0.5,0.5, 60, -15., 15.);
174 fXvYv= new TH2F("fXvYv","Xv-Yv main vtx", 60,-1.5,1.5, 60, -1.5, 1.5);
ad1036b2 175 fHistPtKaoP = new TH1F("fHistPtKaoP", "P_{T}KaonP distribution",300, 0.0,15.0);
ad1036b2 176 fHistPtKaoN = new TH1F("fHistPtKaoN", "P_{T}KaonN distribution",300, 0.0,15.0);
828a6e0c 177 frapiKESD=new TH1F("frapiKESD","rapid Kdistribution", 26,-1.3, 1.3);
178 flifetime= new TH1F("flifetime", "ct study of K-kinks",100,0.,1000.);
179 fradLK= new TH1F("fradLK", "Length of K generated",100,0.,1000.);
180 fradPtRpDt=new TH3F("fradPtRpDt","rad pt rap dat",28,100.,240., 20, 0., 5., 20, -1., 1. );
8e57b720 181 fInvMuNuAll= new TH1F("fInvMuNuAll", " Inv Mass MuNu all kink",180,0.1,1.0);
828a6e0c 182 fQtInvM= new TH2F("fQtInvM", "Q_{T} Versus Inv MuNu ",80, 0., 0.80 , 100 , 0., 0.300);
183 fDCAkink = new TH1F("fDCAkink ", "DCA kink vetrex ",50, 0.0,1.0);
828a6e0c 184 fPosiKink= new TH2F("fPosiKink", "Y vrx kink Vrex ",100, -300.0,300.0,100, -300, 300.);
185 fPosiKinkK= new TH2F("fPosiKinkK", "Y vrx kink VrexK ",100, -300.0,300.0,100, -300, 300.);
186 fPosiKinKXZ= new TH2F("fPosiKinKXZ", "Y vrx kink VrexK ",100, -300.0,300.0,100, -300, 300.);
187 fPosiKinKYZ= new TH2F("fPosiKinKYZ", "Y vrx kink VrexK ",100, -300.0,300.0,100, -300, 300.);
188 fPosiKinKBg= new TH2F("fPosiKinKBg", "z vrx kink rad ",100, -300.0,300.0,100,100., 300.);
189 fQtMothP = new TH2F("fQtMothP", " Qt vrs Mother P", 100, 0., 5.0,100, 0.,0.300);
ad1036b2 190 fTPCSgnlPtpc = new TH2F("fTPCSgnlPtpc","TPC signal de/dx Mom TPC,K ",300,0.0,15.0,100, 0., 250. );
191 fTPCMomNSgnl = new TH2F("fTPCMomNsgnl","TPC signal de/dx Mom TPC,K ",300,0.0,15.0,20 , -10., 10.);
8e57b720 192 fMothKinkMomSgnl = new TH2F("fMothKinkMomSgnl","TPC signal de/dx Mom TPC,Kink ",100,0.0,250.0,100, 0., 250. );
ad1036b2 193 fNSigmTPC = new TH1F("fNSigmTPC","TPC Nsigma de/dx TPC,K ", 30 , -7.5, 7.5);
8e57b720 194 fTPCSgnlKinkDau = new TH2F("fTPCSgnlKinkDau","TPC signal de/dx Mom,K",500,0.0,10.0,100,0.,250.);
ad1036b2 195 fPtKinkPos= new TH1F("fPtKinkPos", "Pos P_{T}Kaon Kink distribution, counts",44, gPt7K0);
196 fPtKinkNeg= new TH1F("fPtKinkNeg", "Neg P_{T}Kaon Kink distribution, counts",44, gPt7K0);
197 fRadNclCln = new TH2F("fRadNclCln","kink radius vrs Nclust,K Clean ",75,100.,250., 80,0, 160);
8e57b720 198 fRatioCrossedRows = new TH1F("fRatioCrossedRows","Ratio crossed rows in TPC",20,0.0,1.0);
199 fRatioCrossedRowsKink = new TH1F("fRatioCrossedRowsKink","Ratio crossed rows in TPC for kinks",20,0.0,1.0);
200 fRadiusPt =new TH2F("fRadiusPt","radius vs pt ",80, 90.,250.,100, 0.,10. );
201 fRadiusPtcln =new TH2F("fRadiusPtcln","radius vs pt clean ",80, 90.,250.,100, 0.,10. );
202 fInvMassMuNuPt =new TH2F("fInvMassMuNuPt","Invariant mass-munu vs pt ",180, 0.10, 1.00, 100, 0.0, 10.0 );
203 fPtCut1 = new TH1F("fPtCut1", "P_{T}Kaon distribution",300, 0.0,15.0);
204 fPtCut2 = new TH1F("fPtCut2", "P_{T}Kaon distribution",300, 0.0,15.0);
205 fPtCut3 = new TH1F("fPtCut3", "P_{T}Kaon distribution",300, 0.0,15.0);
206 fAngMomKKinks = new TH2F("fAngMomKKinks","Decay angle vrs Mother Mom,Kinks",300,0.0,15.0,100,0.,100.);
828a6e0c 207
208 fListOfHistos=new TList();
209
210 fListOfHistos->Add(fHistPtESD);
211 fListOfHistos->Add(fHistPt);
212 fListOfHistos->Add(fHistQtAll);
213 fListOfHistos->Add(fHistQt1);
214 fListOfHistos->Add(fHistQt2);
215 fListOfHistos->Add(fHistPtKaon);
216 fListOfHistos->Add(fHistPtKPDG);
217 fListOfHistos->Add(fHistEta);
218 fListOfHistos->Add(fHistEtaK);
219 fListOfHistos->Add(fptKMC);
220 fListOfHistos->Add(fMultiplMC);
221 fListOfHistos->Add(fESDMult);
222 fListOfHistos->Add(fgenpt);
223 fListOfHistos->Add(frad);
224 fListOfHistos->Add(fKinkKaon);
225 fListOfHistos->Add(fKinKRbn);
226 fListOfHistos->Add(fKinkKaonBg);
227 fListOfHistos->Add(fM1kaon);
828a6e0c 228 fListOfHistos->Add(fPtKink);
229 fListOfHistos->Add(fptKink);
828a6e0c 230 fListOfHistos->Add(fAngMomK);
231 fListOfHistos->Add(fAngMomPi);
232 fListOfHistos->Add(fAngMomKC);
233 fListOfHistos->Add(fMultESDK);
234 fListOfHistos->Add(fMultMCK);
235 fListOfHistos->Add(fSignPtNcl);
236 fListOfHistos->Add(fSignPtEta);
237 fListOfHistos->Add(fEtaNcl);
238 fListOfHistos->Add(fSignPt);
239 fListOfHistos->Add(fChi2NclTPC);
240 fListOfHistos->Add(fRatChi2Ncl);
241 fListOfHistos->Add(fRadiusNcl);
242 fListOfHistos->Add(fTPCSgnlP);
243 fListOfHistos->Add(fTPCSgnlPa);
244 fListOfHistos->Add(fRpr);
245 fListOfHistos->Add(fZpr);
246 fListOfHistos->Add(fdcatoVxXY);
247 fListOfHistos->Add(fnSigmToVx);
248 fListOfHistos->Add(fKinkMothDau);
249 fListOfHistos->Add(fZvXv);
250 fListOfHistos->Add(fZvYv);
251 fListOfHistos->Add(fXvYv);
828a6e0c 252 fListOfHistos->Add(fHistPtKaoP);
253 fListOfHistos->Add(fHistPtKaoN);
254 fListOfHistos->Add(frapiKESD);
255 fListOfHistos->Add(flifetime);
256 fListOfHistos->Add(fradLK);
257 fListOfHistos->Add(fradPtRpDt);
258 fListOfHistos->Add(fInvMuNuAll);
259 fListOfHistos->Add(fQtInvM);
260 fListOfHistos->Add(fDCAkink);
261 fListOfHistos->Add(fPosiKink);
262 fListOfHistos->Add(fPosiKinkK);
263 fListOfHistos->Add(fPosiKinKXZ);
264 fListOfHistos->Add(fPosiKinKYZ);
265 fListOfHistos->Add(fPosiKinKBg);
266 fListOfHistos->Add(fQtMothP);
ad1036b2 267 fListOfHistos->Add(fTPCSgnlPtpc);
268 fListOfHistos->Add(fTPCMomNSgnl);
269 fListOfHistos->Add(fMothKinkMomSgnl);
270 fListOfHistos->Add(fNSigmTPC);
271 fListOfHistos->Add(fTPCSgnlKinkDau);
272 fListOfHistos->Add(fPtKinkPos);
273 fListOfHistos->Add(fPtKinkNeg);
274 fListOfHistos->Add(fRadNclCln);
8e57b720 275 fListOfHistos->Add(fRatioCrossedRows);
276 fListOfHistos->Add(fRatioCrossedRowsKink);
277 fListOfHistos->Add(fRadiusPt);
278 fListOfHistos->Add(fRadiusPtcln);
279 fListOfHistos->Add(fInvMassMuNuPt);
280 fListOfHistos->Add(fPtCut1);
281 fListOfHistos->Add(fPtCut2);
282 fListOfHistos->Add(fPtCut3);
283 fListOfHistos->Add(fAngMomKKinks);
828a6e0c 284
285}
828a6e0c 286//________________________________________________________________________
287void AliAnalysisKinkESDat::UserExec(Option_t *)
288{
289 // Main loop
290 // Called for each event
291
292 // Process MC truth, therefore we receive the AliAnalysisManager and ask it for the AliMCEventHandler
293 // This handler can return the current MC event
294
295 AliVEvent *event = InputEvent();
296 if (!event) {
297 Printf("ERROR: Could not retrieve event");
298 return;
299 }
300
301 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
302 if (!esd) {
303 Printf("ERROR: Could not retrieve esd");
304 return;
305 }
ad1036b2 306// Number ESD tracks
307 Int_t nESDTracks = esd->GetNumberOfTracks();
308 fMultiplMC->Fill(nESDTracks);
309//
828a6e0c 310//==================check of Physics selection?
311 Bool_t isSelected =
312((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected()&AliVEvent::kMB;
313
314 if ( isSelected ==kFALSE) return; // 24/6/11 apo MF
315//
828a6e0c 316 fMultMCK->Fill(nESDTracks);
ad1036b2 317//
318//===============Marek's multiplicity
828a6e0c 319 Float_t refmultiplicity=fCutsMul->CountAcceptedTracks(esd);
320 if(fLowMulcut>-1)
321 {
322 if(refmultiplicity<fLowMulcut)
323 return;
324 }
325 if(fUpMulcut>-1)
326 {
327 if(refmultiplicity>fUpMulcut)
328 return;
329 }
330
331
332
333 fMultESDK->Fill(refmultiplicity);
334
335//
828a6e0c 336
337
338
339 const AliESDVertex *vertex=GetEventVertex(esd); // 22/8
ad1036b2 340 if(!vertex) return;
828a6e0c 341//
342 Double_t vpos[3];
343 vertex->GetXYZ(vpos);
344 fZpr->Fill(vpos[2]);
ad1036b2 345 if (TMath::Abs( vpos[2] ) > 10. ) return;
828a6e0c 346
347
348
349 Double_t vtrack[3], ptrack[3];
350
351
352 Int_t nESDTracK = 0;
353 Int_t nESDTrKink = 0;
354
355 Int_t nGoodTracks = esd->GetNumberOfTracks();
356 fESDMult->Fill(nGoodTracks);
357
ad1036b2 358 Double_t nsigmall = 100.0;
828a6e0c 359 Double_t nsigma = 100.0;
ad1036b2 360 Double_t nsigmaPion =-100.0;
8e57b720 361 Double_t nsigmaDau =-100.0;
362 Double_t dEdxKinkDau =0.0;
363 Double_t KinkDauCl =0.0;
fcd98e58 364// apo Eftihi
365 if(!fPIDResponse) {
366 AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
367 AliInputEventHandler* inputHandler =
368(AliInputEventHandler*)(man->GetInputEventHandler());
369 fPIDResponse = inputHandler->GetPIDResponse();
370 }
371
ad1036b2 372// loop on kink daughters
373 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
374
375 AliESDtrack* trackD = esd->GetTrack(iTrack);
376 if (!trackD) {
377 Printf("ERROR: Could not receive track %d", iTrack);
378 continue;
379 }
8e57b720 380//
ad1036b2 381 Int_t indexKinkDau=trackD->GetKinkIndex(0);
382// daughter kink
383// AliESDkink *kink=esd->GetKink(TMath::Abs(indexKinkDau)-1);
8e57b720 384 if ( indexKinkDau > 0 ) {
385 Int_t labelD = trackD->GetLabel();
386 labelD = TMath::Abs(labelD);
387 nsigmaPion = (fPIDResponse->NumberOfSigmasTPC(trackD , AliPID::kPion));// 26/10 eftihis
388 dEdxKinkDau = (trackD->GetTPCsignal() ) ; // daughter kink dEdx
389 KinkDauCl=(trackD->GetTPCclusters(0) ) ;
390 }
fcd98e58 391//if((indexKinkDau >0)&& (nsigmaPion>1.2)) fTPCSgnlKinkDau->Fill(trackD->P(), (trackD->GetTPCsignal() ) ) ; // daughter kink
8e57b720 392// Ayto mexri 26/11/2012 if(indexKinkDau >0) fTPCSgnlKinkDau->Fill(trackD->P(), (trackD->GetTPCsignal() ) ) ; // daughter kink
ad1036b2 393 }
fcd98e58 394
828a6e0c 395// track loop
ad1036b2 396//
828a6e0c 397 for (Int_t iTracks = 0; iTracks < esd->GetNumberOfTracks(); iTracks++) {
398
399 AliESDtrack* track = esd->GetTrack(iTracks);
400 if (!track) {
401 Printf("ERROR: Could not receive track %d", iTracks);
402 continue;
403 }
404
405 fHistPt->Fill(track->Pt());
406
407
408 // sigmas
fcd98e58 409 nsigmall = (fPIDResponse->NumberOfSigmasTPC(track, AliPID::kKaon));
ad1036b2 410 // nsigmaPion= (fESDpid->NumberOfSigmasTPC(track,AliPID::kPion));
fcd98e58 411 nsigma = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(track, AliPID::kKaon));
828a6e0c 412
8e57b720 413//=======================new
414 Float_t nCrossedRowsTPC = track->GetTPCClusterInfo(2,1);
415 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;
416 if (track->GetTPCNclsF()>0) {
417 ratioCrossedRowsOverFindableClustersTPC = nCrossedRowsTPC / track->GetTPCNclsF();
418 fRatioCrossedRows->Fill(ratioCrossedRowsOverFindableClustersTPC);
419 }
420//_______
828a6e0c 421
fcd98e58 422 Int_t indexKinkPos=track->GetKinkIndex(0); // kink index
828a6e0c 423
424 Int_t tpcNCl = track->GetTPCclusters(0);
425 Double_t tpcSign = track->GetSign();
426
427 Int_t label = track->GetLabel();
428 label = TMath::Abs(label);
429
430
431 UInt_t status=track->GetStatus();
432
433 if((status&AliESDtrack::kITSrefit)==0) continue;
434 if((status&AliESDtrack::kTPCrefit)==0) continue;
ad1036b2 435 if((track->GetTPCchi2()/track->GetTPCclusters(0))>4.0) continue;
828a6e0c 436
437 Double_t extCovPos[15];
438 track->GetExternalCovariance(extCovPos);
828a6e0c 439
440
441 track->GetXYZ(vtrack);
442 fXvYv->Fill(vtrack[0],vtrack[1]);
443 fZvYv->Fill(vtrack[0],vtrack[2]);
444 fZvXv->Fill(vtrack[1],vtrack[2]);
445
446// track momentum, rapidity calculation
447 track->GetPxPyPz(ptrack);
448
449 TVector3 trackMom(ptrack[0],ptrack[1],ptrack[2]);
450
ad1036b2 451// K-rapidity calcualtion
828a6e0c 452 Double_t etracK= TMath::Sqrt(trackMom.Mag()*trackMom.Mag() + 0.493677 *0.493677 );
453 Double_t rapiditK = 0.5 * (TMath::Log( (etracK + ptrack[2] ) / ( etracK - ptrack[2]) )) ;
454
455 Double_t trackEta=trackMom.Eta();
828a6e0c 456 Double_t trackPt = track->Pt();
457
458
459
460 Float_t bpos[2];
461 Float_t bCovpos[3];
462 track->GetImpactParameters(bpos,bCovpos);
463
464 if (bCovpos[0]<=0 || bCovpos[2]<=0) {
465 Printf("Estimated b resolution lower or equal zero!");
466 bCovpos[0]=0; bCovpos[2]=0;
467 }
468
469 Float_t dcaToVertexXYpos = bpos[0];
470 Float_t dcaToVertexZpos = bpos[1];
471
472 fRpr->Fill(dcaToVertexZpos);
473
ad1036b2 474 if((TMath::Abs(dcaToVertexXYpos)>0.3)||(TMath::Abs(dcaToVertexZpos)>2.5)) nESDTrKink++; // count of second 23Jul11
8e57b720 475 //if((TMath::Abs(dcaToVertexXYpos)>0.3)||(TMath::Abs(dcaToVertexZpos)>2.5))
476// continue; // allagi 23Jul11
477
478 // apo Filhmona if (!fMaxDCAtoVtxCut->AcceptTrack(track)) continue;
479// Float_t MaxDCAxy = fMaxDCAtoVtxCut->GetMaxDCAToVertexXYPtDep(track );
480 // if (MaxDCAxy > 2.4 ) continue ;
481
482 fdcatoVxXY->Fill(dcaToVertexXYpos);
828a6e0c 483//
828a6e0c 484
485// track Mult. after selection
486 nESDTracK++;
487 //
488//=========================================
489 fHistPtESD->Fill(track->Pt());
490
491 // Add Kink analysis =============================
492
ad1036b2 493// daughter kink
494//if(indexKinkPos >0)fTPCSgnlKinkDau->Fill(track->P(), (track->GetTPCsignal() ) ) ; // daughter kink
495
828a6e0c 496// loop on kinks
497 if(indexKinkPos<0){ ////mother kink
828a6e0c 498
fcd98e58 499 fptKMC ->Fill( track->Pt() ); // Pt from tracks , all kinks
500
8e57b720 501 fRatioCrossedRowsKink->Fill(ratioCrossedRowsOverFindableClustersTPC);
828a6e0c 502 // select kink class
503
504 AliESDkink *kink=esd->GetKink(TMath::Abs(indexKinkPos)-1);
505//
506 // DCA kink
507 Double_t Dist2 = kink->GetDistance();
ad1036b2 508 // fDCAkink->Fill( Dist2 );
509 // if (Dist2 > 0.2) continue; // systematics 11/8/11
828a6e0c 510//
511
ad1036b2 512// TPC mother momentum
828a6e0c 513
514 const TVector3 vposKink(kink->GetPosition());
515 fPosiKink ->Fill( vposKink[0], vposKink[1] );
828a6e0c 516 Double_t dzKink=vpos[2]-vposKink[2];
828a6e0c 517//
8e57b720 518// lifitime
828a6e0c 519 Double_t tanLamda = track->GetTgl(); // 25/6/2010
520
521 Double_t lifeKink= (TMath::Abs( dzKink ))*( TMath::Sqrt(1.+ tanLamda*tanLamda) ) / (TMath::Abs( tanLamda)) ;
8e57b720 522//
828a6e0c 523 const TVector3 motherMfromKink(kink->GetMotherP());
524 const TVector3 daughterMKink(kink->GetDaughterP());
525
526 Float_t qT=kink->GetQt();
527 Float_t motherPt=motherMfromKink.Pt();
ad1036b2 528// Kink mother momentum
529 Double_t trMomTPCKink=motherMfromKink.Mag();
530// TPC mother momentun
531 Double_t trMomTPC=track->GetTPCmomentum();
8e57b720 532 // fTPCSgnlKinkDau->Fill( daughterMKink.Mag(), dEdxKinkDau ) ; // daughter kink
533 //
828a6e0c 534 fHistQtAll->Fill(qT) ; // Qt distr
535
536 fptKink->Fill(motherMfromKink.Pt()); /// pt from kink
537
538 Float_t kinkAngle=TMath::RadToDeg()*kink->GetAngle(2);
539
fcd98e58 540// rapiditya nd pt selection
828a6e0c 541 if( (TMath::Abs(rapiditK )) > 0.7 ) continue;
ad1036b2 542 if ( (track->Pt())<.200)continue; // new Feb 2012
828a6e0c 543
544 fQtMothP->Fill( track->P(), qT);
545
fcd98e58 546 if ( qT> 0.04) fHistQt1 ->Fill(qT) ; // Qt distr
828a6e0c 547
548
549
fcd98e58 550 fHistEta->Fill(trackEta) ; // Eta distr
551 fHistQt2->Fill(qT); //
8e57b720 552 fKinkKaonBg->Fill(motherPt);
828a6e0c 553
554// maximum decay angle at a given mother momentum
555 //Double_t maxDecAngKmu=f1->Eval(motherMfromKink.Mag(),0.,0.,0.);
556 Double_t maxDecAngKmu=f1->Eval(track->P() ,0.,0.,0.);
557 Double_t maxDecAngpimu=f2->Eval( track->P(), 0.,0.,0.);
fcd98e58 558
559// fake kinks are removed
828a6e0c 560 if( (kinkAngle<2.) ) continue;
561
562
563 // BG ?????==============
564 if ( TMath::Abs(vposKink[2]) > 225. ) continue ;
565 if ( TMath::Abs(vposKink[2]) < 0.5 ) continue ;
828a6e0c 566//
8e57b720 567 fPtCut1 ->Fill(trackPt );
fcd98e58 568 fAngMomPi->Fill( track->P(), kinkAngle);
828a6e0c 569//
570// invariant mass of mother track decaying to mu
571 Float_t energyDaughterMu=TMath::Sqrt(daughterMKink.Mag()*daughterMKink.Mag()+0.105658*0.105658);
572 Float_t p1XM= motherMfromKink.Px();
573 Float_t p1YM= motherMfromKink.Py();
574 Float_t p1ZM= motherMfromKink.Pz();
575 Float_t p2XM= daughterMKink.Px();
576 Float_t p2YM= daughterMKink.Py();
577 Float_t p2ZM= daughterMKink.Pz();
578 Float_t p3Daughter=TMath::Sqrt(((p1XM-p2XM)*(p1XM-p2XM))+((p1YM-p2YM)*(p1YM-p2YM))+((p1ZM-p2ZM)*(p1ZM-p2ZM)));
579 Double_t invariantMassKmu= TMath::Sqrt((energyDaughterMu+p3Daughter)*(energyDaughterMu+p3Daughter)-motherMfromKink.Mag()*motherMfromKink.Mag());
580 fQtInvM -> Fill ( invariantMassKmu, qT);
581 fInvMuNuAll->Fill(invariantMassKmu);
fcd98e58 582//
8e57b720 583 fRadiusPt->Fill( kink->GetR(), trackPt); //
fcd98e58 584 // radius and Minv selection
585 if( ( kink->GetR()> 120 ) && ( kink->GetR() < 210 ) ) {
586 // for systematics if( ( kink->GetR()> 130 ) && ( kink->GetR() < 200 ) ) {
828a6e0c 587 if (qT>0.12) fAngMomKC->Fill(track->P(), kinkAngle);
588 if ( qT>0.12) fM1kaon->Fill(invariantMassKmu);
589 if ( qT > 0.12)
590 fRadiusNcl->Fill( (kink->GetR()) ,(track->GetTPCclusters(0) ) ) ;
591 }
fcd98e58 592// tails cleaning
593 if( ( tpcNCl<20) ) continue; // test 27 feb 2012 ,, OK
594// cleaning BG in tails
ad1036b2 595 Int_t tpcNClHigh = -51.67+ (11./12.) *( kink->GetR() ) ;
596 if ( tpcNCl > tpcNClHigh) continue;
597
ad1036b2 598 Int_t tpcNClMin = -85.5 + (65./95.) *( kink->GetR() ) ;
ad1036b2 599 if ( tpcNCl < tpcNClMin ) continue;
828a6e0c 600
8e57b720 601 // back, 20/1/2013 if (ratioCrossedRowsOverFindableClustersTPC< 0.5) continue;// test 14/1/2013
828a6e0c 602//
fcd98e58 603 fHistPtKPDG->Fill(track->Pt()); // ALL K-candidates until now
ad1036b2 604 // if((kinkAngle>maxDecAngpimu)&&(qT>0.12)&&(qT<0.30)&&((kink->GetR()>=120.)&&(kink->GetR()<=210.))&&(TMath::Abs(rapiditK)<0.7)&&(invariantMassKmu<0.6)){
fcd98e58 605 if((kinkAngle>maxDecAngpimu)&&(qT>0.12)&&(qT<0.30)&&((kink->GetR()>=120.)&&(kink->GetR()<=210.))&&(TMath::Abs(rapiditK)<0.7)&&(invariantMassKmu<0.8)){
606 // systematics if((kinkAngle>maxDecAngpimu)&&(qT>0.12)&&(qT<0.30)&&((kink->GetR()>=130.)&&(kink->GetR()<=200.))&&(TMath::Abs(rapiditK)<0.7)&&(invariantMassKmu<0.8)){
607//
8e57b720 608 fAngMomKKinks->Fill(track->P(), kinkAngle);
609 fPtCut2 ->Fill(trackPt );
fcd98e58 610// maximum angles selection with some error cut
828a6e0c 611 if( (kinkAngle<maxDecAngpimu*1.2) ) continue;
612 if ( (kinkAngle>maxDecAngKmu*.98) && ( track->P() >1.2 )) continue; ///5/5/2010
613
8e57b720 614 fPtCut3 ->Fill(trackPt );
fcd98e58 615// here the kaons selected by the decay features
ad1036b2 616 fTPCSgnlPa->Fill( track->GetInnerParam()->GetP() ,(track->GetTPCsignal() ) ) ;
fcd98e58 617//
8e57b720 618 // NO dEdx cut test 9/2/13 if ( nsigma > 3.5) continue;
fcd98e58 619 // system if ( nsigma > 4.0) continue; // gia systamatic error
ad1036b2 620//
fcd98e58 621// next plots for the identified kaons by the kink analysis
622
623 fHistPtKaon->Fill(track->Pt()); //
624 if(tpcSign >0.) fHistPtKaoP->Fill( track->Pt() ) ; //
625 if ( tpcSign <0.) fHistPtKaoN->Fill( track->Pt() ) ; //
ad1036b2 626 fTPCSgnlP->Fill(track->GetInnerParam()->GetP(), (track->GetTPCsignal() ) ) ;
627 fRadNclCln->Fill( (kink->GetR()) ,(track->GetTPCclusters(0) ) ) ;
8e57b720 628 fRadiusPtcln->Fill( kink->GetR(), trackPt); //
629 fInvMassMuNuPt ->Fill(invariantMassKmu, trackPt);
ad1036b2 630
8e57b720 631 // fTPCSgnlPtpc->Fill(trMomTPC , (track->GetTPCsignal() ) ) ;
632 //fMothKinkMomSgnl ->Fill(trMomTPCKink , (track->GetTPCsignal() ) ) ;
633 fMothKinkMomSgnl ->Fill( dEdxKinkDau , (track->GetTPCsignal() ) ) ;
634//
635 fTPCSgnlKinkDau->Fill( daughterMKink.Mag(), dEdxKinkDau ) ; // daughter kink
ad1036b2 636//
637 fTPCMomNSgnl->Fill(trMomTPC ,nsigmall );
638 fNSigmTPC ->Fill(nsigmall );
639//
828a6e0c 640 frad->Fill(kink->GetR()); // kink
641 fradLK->Fill(lifeKink ); // kink
642 fHistEtaK->Fill(trackEta);
643 frapiKESD ->Fill(rapiditK); // rapidityof kaons
644 fPosiKinKBg->Fill( vposKink[2], kink->GetR() );
645
646 Float_t signPt= tpcSign*trackPt;
647 fSignPtNcl->Fill( signPt , tpcNCl ); /// 28/4/2010
648 fSignPtEta->Fill( signPt , rapiditK );
649 fEtaNcl->Fill( rapiditK, tpcNCl );
650 fSignPt->Fill( signPt );
651 fChi2NclTPC->Fill( (track->GetTPCchi2() ) , tpcNCl );
652 fRatChi2Ncl-> Fill ( (track->GetTPCchi2()/track->GetTPCclusters(0) )) ;
8e57b720 653 // fdcatoVxXY->Fill(dcaToVertexXYpos);
654 //if( dEdxKinkDau> 1.5* (track->GetTPCsignal() ) ) fKinkMothDau->Fill(track->P(),daughterMKink.Mag());
655 // if((dEdxKinkDau> 80. ) && (dEdxKinkDau > 4.*nsigmaPion) ) fKinkMothDau->Fill(track->P(),daughterMKink.Mag());
656 if (nsigmaPion> 3. ) fTPCSgnlPtpc->Fill( daughterMKink.Mag(), dEdxKinkDau ) ;
657 //if (TMath::Abs(dEdxKinkDau -(track->GetTPCsignal() )> 10. )) fTPCSgnlPtpc->Fill( daughterMKink.Mag(), dEdxKinkDau ) ;
828a6e0c 658 flifetime->Fill(( lifeKink*.493667 ) /track->P() ) ;
659 fKinkKaon->Fill(track->Pt());
ad1036b2 660 fDCAkink->Fill( Dist2 );
661
662 fPtKink->Fill(track->Pt()); /// K0 bins
663 if(tpcSign >0.) fPtKinkPos ->Fill( track->Pt() ) ; //K-plus bins K0
664 if(tpcSign <0.) fPtKinkNeg ->Fill( track->Pt() ) ; //K-minus bins K0
665 fKinKRbn->Fill(track->Pt()); // TOF
828a6e0c 666
828a6e0c 667 fradPtRpDt->Fill( kink->GetR(), 1./track->Pt(), rapiditK);
668 fAngMomK->Fill( track->P(), kinkAngle);
669 fPosiKinkK->Fill( vposKink[0], vposKink[1] );
670 fPosiKinKXZ->Fill( vposKink[2], vposKink[0] );
671 fPosiKinKYZ->Fill( vposKink[2], vposKink[1] );
672
673 } // kink selection
674
675
676 } //End Kink Information
677
678
679 } //track loop
ad1036b2 680// } //daughter loop
828a6e0c 681
682 // fMultiplMC->Fill(nESDTracK );
683
684 PostData(1, fListOfHistos);
685
686}
687
688//________________________________________________________________________
689void AliAnalysisKinkESDat::Terminate(Option_t *)
690{
691 // Draw result to the screen
692 // Called once at the end of the query
693
694}
695
696//____________________________________________________________________//
697
698
699const AliESDVertex* AliAnalysisKinkESDat::GetEventVertex(const AliESDEvent* esd) const
700 // Get the vertex from the ESD and returns it if the vertex is valid
701
702{
703 // Get the vertex
704
828a6e0c 705 const AliESDVertex* vertex = esd->GetPrimaryVertexTracks();
706
828a6e0c 707 if((vertex->GetStatus()==kTRUE)) return vertex;
708 else
709 {
710 vertex = esd->GetPrimaryVertexSPD();
711 if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>0)) return vertex;
828a6e0c 712 else
713 return 0;
714 }
ad1036b2 715}