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 * The method is applied in pp and Pb-Pb real data. *
7 **************************************************************************/
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 //-----------------------------------------------------------------
22 #include "TParticle.h"
26 #include "AliAnalysisTask.h"
27 #include "AliAnalysisManager.h"
29 #include "AliVEvent.h"
30 #include "AliESDEvent.h"
31 #include "AliMCEvent.h"
32 #include "AliAnalysisKinkESDat.h"
34 #include "AliESDpid.h"
36 #include "AliESDkink.h"
37 #include "AliESDtrack.h"
38 #include "AliPhysicsSelectionTask.h"
39 #include "AliInputEventHandler.h"
40 #include "AliESDInputHandler.h"
41 #include "AliAnalysisManager.h"
42 #include "AliVEvent.h"
43 #include "AliESDtrackCuts.h"
44 #include "AliPIDResponse.h"
45 /////////#include "AliTENDERSupplies.h"
46 ClassImp(AliAnalysisKinkESDat)
49 //________________________________________________________________________
50 AliAnalysisKinkESDat::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),
53 fKinkKaon(0),fKinKRbn(0), fKinkKaonBg(0), fM1kaon(0), fPtKink(0), fptKink(0),
54 fAngMomK(0),fAngMomPi(0), fAngMomKC(0), fMultESDK(0), fMultMCK(0),
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),
57 fZvXv(0),fZvYv(0), fXvYv(0), fHistPtKaoP(0), fHistPtKaoN(0),frapiKESD(0), flifetime(), fradLK(0),
58 fradPtRpDt(0), fInvMuNuAll(0), fQtInvM(0),
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),
61 fRatioCrossedRows(0), fRatioCrossedRowsKink(0),fRadiusPt(0), fRadiusPtcln(0), fInvMassMuNuPt(0), fPtCut1(0), fPtCut2(0), fPtCut3(0),
64 fListOfHistos(0),fLowMulcut(-1),fUpMulcut(-1), fKinkRadUp(200), fKinkRadLow(130), fCutsMul(0), fMaxDCAtoVtxCut(0), fPIDResponse(0)
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);
78 fCutsMul->SetRequireITSRefit(kTRUE);
79 fCutsMul->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
80 AliESDtrackCuts::kAny);
81 fCutsMul->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01");
83 fCutsMul->SetMaxDCAToVertexZ(2);
84 fCutsMul->SetDCAToVertex2D(kFALSE);
85 fCutsMul->SetRequireSigmaToVertex(kFALSE);
87 fCutsMul->SetEtaRange(-0.8,+0.8);
88 fCutsMul->SetPtRange(0.15, 1e10);
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);
98 DefineOutput(1, TList::Class());
101 //________________________________________________________________________
102 void AliAnalysisKinkESDat::UserCreateOutputObjects()
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());
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());
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,
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
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,
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,
127 3.2, 3.4, 3.6, 3.8, 4.0, 4.2, 4.4, 4.6, 4.8,5.0 }; // Barbara TOF Kch
129 fHistPtESD = new TH1F("fHistPtESD", "P_{T} distribution",300, 0.0,15.0);
130 fHistPtESD->GetXaxis()->SetTitle("P_{T} (GeV/c)");
131 fHistPtESD->GetYaxis()->SetTitle("dN/dP_{T} (c/GeV)");
132 fHistPtESD->SetMarkerStyle(kFullCircle);
133 fHistPt = new TH1F("fHistPt", "P_{T} distribution",300, 0.0,15.0);
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);
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);
139 fHistEta= new TH1F("fHistEta", "Eta distribution", 26,-1.3, 1.3);
140 fHistEtaK= new TH1F("fHistEtaK", "EtaK distribution", 26,-1.3, 1.3);
141 fptKMC= new TH1F("fptKMC", "P_{T}Kaon generated",300, 0.0,15.0);
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);
144 fgenpt= new TH1F("fgenpt", "genpt K distribution",300, 0.0,15.0);
145 //frad= new TH1F("frad", "radius K generated",100, 50., 250.0);
146 frad= new TH1F("frad", "radius K generated",100, 0.,1000.0);
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);
150 fM1kaon= new TH1F("fM1kaon","Invar m(kaon) from kink->mu+netrino decay",180,0.10, 1.0);
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);
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);
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);
164 fRadiusNcl = new TH2F("fRadiusNcl","kink radius vrs Nclust,K",75,100.,250., 80,0, 160);
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.);
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);
175 fHistPtKaoP = new TH1F("fHistPtKaoP", "P_{T}KaonP distribution",300, 0.0,15.0);
176 fHistPtKaoN = new TH1F("fHistPtKaoN", "P_{T}KaonN distribution",300, 0.0,15.0);
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. );
181 fInvMuNuAll= new TH1F("fInvMuNuAll", " Inv Mass MuNu all kink",180,0.1,1.0);
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);
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);
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.);
192 fMothKinkMomSgnl = new TH2F("fMothKinkMomSgnl","TPC signal de/dx Mom TPC,Kink ",100,0.0,250.0,100, 0., 250. );
193 fNSigmTPC = new TH1F("fNSigmTPC","TPC Nsigma de/dx TPC,K ", 30 , -7.5, 7.5);
194 fTPCSgnlKinkDau = new TH2F("fTPCSgnlKinkDau","TPC signal de/dx Mom,K",500,0.0,10.0,100,0.,250.);
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);
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.);
208 fListOfHistos=new TList();
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);
228 fListOfHistos->Add(fPtKink);
229 fListOfHistos->Add(fptKink);
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);
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);
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);
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);
286 //________________________________________________________________________
287 void AliAnalysisKinkESDat::UserExec(Option_t *)
290 // Called for each event
292 // Process MC truth, therefore we receive the AliAnalysisManager and ask it for the AliMCEventHandler
293 // This handler can return the current MC event
295 AliVEvent *event = InputEvent();
297 Printf("ERROR: Could not retrieve event");
301 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
303 Printf("ERROR: Could not retrieve esd");
307 Int_t nESDTracks = esd->GetNumberOfTracks();
308 fMultiplMC->Fill(nESDTracks);
310 //==================check of Physics selection?
312 ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected()&AliVEvent::kMB;
314 if ( isSelected ==kFALSE) return; // 24/6/11 apo MF
316 fMultMCK->Fill(nESDTracks);
318 //===============Marek's multiplicity
319 Float_t refmultiplicity=fCutsMul->CountAcceptedTracks(esd);
322 if(refmultiplicity<fLowMulcut)
327 if(refmultiplicity>fUpMulcut)
333 fMultESDK->Fill(refmultiplicity);
339 const AliESDVertex *vertex=GetEventVertex(esd); // 22/8
343 vertex->GetXYZ(vpos);
345 if (TMath::Abs( vpos[2] ) > 10. ) return;
349 Double_t vtrack[3], ptrack[3];
353 Int_t nESDTrKink = 0;
355 Int_t nGoodTracks = esd->GetNumberOfTracks();
356 fESDMult->Fill(nGoodTracks);
358 Double_t nsigmall = 100.0;
359 Double_t nsigma = 100.0;
360 Double_t nsigmaPion =-100.0;
361 Double_t nsigmaDau =-100.0;
362 Double_t dEdxKinkDau =0.0;
363 Double_t KinkDauCl =0.0;
366 AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
367 AliInputEventHandler* inputHandler =
368 (AliInputEventHandler*)(man->GetInputEventHandler());
369 fPIDResponse = inputHandler->GetPIDResponse();
372 // loop on kink daughters
373 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
375 AliESDtrack* trackD = esd->GetTrack(iTrack);
377 Printf("ERROR: Could not receive track %d", iTrack);
381 Int_t indexKinkDau=trackD->GetKinkIndex(0);
383 // AliESDkink *kink=esd->GetKink(TMath::Abs(indexKinkDau)-1);
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) ) ;
391 //if((indexKinkDau >0)&& (nsigmaPion>1.2)) fTPCSgnlKinkDau->Fill(trackD->P(), (trackD->GetTPCsignal() ) ) ; // daughter kink
392 // Ayto mexri 26/11/2012 if(indexKinkDau >0) fTPCSgnlKinkDau->Fill(trackD->P(), (trackD->GetTPCsignal() ) ) ; // daughter kink
397 for (Int_t iTracks = 0; iTracks < esd->GetNumberOfTracks(); iTracks++) {
399 AliESDtrack* track = esd->GetTrack(iTracks);
401 Printf("ERROR: Could not receive track %d", iTracks);
405 fHistPt->Fill(track->Pt());
409 nsigmall = (fPIDResponse->NumberOfSigmasTPC(track, AliPID::kKaon));
410 // nsigmaPion= (fESDpid->NumberOfSigmasTPC(track,AliPID::kPion));
411 nsigma = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(track, AliPID::kKaon));
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);
422 Int_t indexKinkPos=track->GetKinkIndex(0); // kink index
424 Int_t tpcNCl = track->GetTPCclusters(0);
425 Double_t tpcSign = track->GetSign();
427 Int_t label = track->GetLabel();
428 label = TMath::Abs(label);
431 UInt_t status=track->GetStatus();
433 if((status&AliESDtrack::kITSrefit)==0) continue;
434 if((status&AliESDtrack::kTPCrefit)==0) continue;
435 if((track->GetTPCchi2()/track->GetTPCclusters(0))>4.0) continue;
437 Double_t extCovPos[15];
438 track->GetExternalCovariance(extCovPos);
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]);
446 // track momentum, rapidity calculation
447 track->GetPxPyPz(ptrack);
449 TVector3 trackMom(ptrack[0],ptrack[1],ptrack[2]);
451 // K-rapidity calcualtion
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]) )) ;
455 Double_t trackEta=trackMom.Eta();
456 Double_t trackPt = track->Pt();
462 track->GetImpactParameters(bpos,bCovpos);
464 if (bCovpos[0]<=0 || bCovpos[2]<=0) {
465 Printf("Estimated b resolution lower or equal zero!");
466 bCovpos[0]=0; bCovpos[2]=0;
469 Float_t dcaToVertexXYpos = bpos[0];
470 Float_t dcaToVertexZpos = bpos[1];
472 fRpr->Fill(dcaToVertexZpos);
474 // 14/2/13 /================/ if((TMath::Abs(dcaToVertexXYpos)>0.3)||(TMath::Abs(dcaToVertexZpos)>2.5)) nESDTrKink++; // count of second 23Jul11
475 //if((TMath::Abs(dcaToVertexXYpos)>0.3)||(TMath::Abs(dcaToVertexZpos)>2.5))
476 // continue; // allagi 23Jul11
478 if (!fMaxDCAtoVtxCut->AcceptTrack(track)) continue;
479 // Float_t MaxDCAxy = fMaxDCAtoVtxCut->GetMaxDCAToVertexXYPtDep(track );
480 // if (MaxDCAxy > 2.4 ) continue ;
482 fdcatoVxXY->Fill(dcaToVertexXYpos);
485 // track Mult. after selection
488 //=========================================
489 fHistPtESD->Fill(track->Pt());
491 // Add Kink analysis =============================
494 //if(indexKinkPos >0)fTPCSgnlKinkDau->Fill(track->P(), (track->GetTPCsignal() ) ) ; // daughter kink
497 if(indexKinkPos<0){ ////mother kink
499 fptKMC ->Fill( track->Pt() ); // Pt from tracks , all kinks
501 fRatioCrossedRowsKink->Fill(ratioCrossedRowsOverFindableClustersTPC);
504 AliESDkink *kink=esd->GetKink(TMath::Abs(indexKinkPos)-1);
507 Double_t Dist2 = kink->GetDistance();
508 // fDCAkink->Fill( Dist2 );
509 // if (Dist2 > 0.2) continue; // systematics 11/8/11
512 // TPC mother momentum
514 const TVector3 vposKink(kink->GetPosition());
515 fPosiKink ->Fill( vposKink[0], vposKink[1] );
516 Double_t dzKink=vpos[2]-vposKink[2];
519 Double_t tanLamda = track->GetTgl(); // 25/6/2010
521 Double_t lifeKink= (TMath::Abs( dzKink ))*( TMath::Sqrt(1.+ tanLamda*tanLamda) ) / (TMath::Abs( tanLamda)) ;
523 const TVector3 motherMfromKink(kink->GetMotherP());
524 const TVector3 daughterMKink(kink->GetDaughterP());
526 Float_t qT=kink->GetQt();
527 Float_t motherPt=motherMfromKink.Pt();
528 // Kink mother momentum
529 Double_t trMomTPCKink=motherMfromKink.Mag();
530 // TPC mother momentun
531 Double_t trMomTPC=track->GetTPCmomentum();
532 // fTPCSgnlKinkDau->Fill( daughterMKink.Mag(), dEdxKinkDau ) ; // daughter kink
534 fHistQtAll->Fill(qT) ; // Qt distr
536 fptKink->Fill(motherMfromKink.Pt()); /// pt from kink
538 Float_t kinkAngle=TMath::RadToDeg()*kink->GetAngle(2);
540 // rapiditya nd pt selection
541 if( (TMath::Abs(rapiditK )) > 0.7 ) continue;
542 if ( (track->Pt())<.200)continue; // new Feb 2012
544 fQtMothP->Fill( track->P(), qT);
546 if ( qT> 0.04) fHistQt1 ->Fill(qT) ; // Qt distr
550 fHistEta->Fill(trackEta) ; // Eta distr
551 fHistQt2->Fill(qT); //
552 fKinkKaonBg->Fill(motherPt);
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.);
559 // fake kinks are removed
560 if( (kinkAngle<2.) ) continue;
563 // BG ?????==============
564 if ( TMath::Abs(vposKink[2]) > 225. ) continue ;
565 if ( TMath::Abs(vposKink[2]) < 0.5 ) continue ;
567 fPtCut1 ->Fill(trackPt );
568 fAngMomPi->Fill( track->P(), kinkAngle);
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);
583 fRadiusPt->Fill( kink->GetR(), trackPt); //
584 // radius and Minv selection
585 //if( ( kink->GetR()> 120 ) && ( kink->GetR() < 210 ) ) {
586 if( ( kink->GetR()> fKinkRadLow ) && ( kink->GetR() <fKinkRadUp ) ) {
587 // for systematics if( ( kink->GetR()> 130 ) && ( kink->GetR() < 200 ) ) {
588 if (qT>0.12) fAngMomKC->Fill(track->P(), kinkAngle);
589 if ( qT>0.12) fM1kaon->Fill(invariantMassKmu);
591 fRadiusNcl->Fill( (kink->GetR()) ,(track->GetTPCclusters(0) ) ) ;
594 if( ( tpcNCl<20) ) continue; // test 27 feb 2012 ,, OK
595 // cleaning BG in tails
596 Int_t tpcNClHigh = -51.67+ (11./12.) *( kink->GetR() ) ;
597 if ( tpcNCl > tpcNClHigh) continue;
599 Int_t tpcNClMin = -85.5 + (65./95.) *( kink->GetR() ) ;
600 if ( tpcNCl < tpcNClMin ) continue;
602 // back, 20/1/2013 if (ratioCrossedRowsOverFindableClustersTPC< 0.5) continue;// test 14/1/2013
604 fHistPtKPDG->Fill(track->Pt()); // ALL K-candidates until now
605 // if((kinkAngle>maxDecAngpimu)&&(qT>0.12)&&(qT<0.30)&&((kink->GetR()>=120.)&&(kink->GetR()<=210.))&&(TMath::Abs(rapiditK)<0.7)&&(invariantMassKmu<0.6)){
606 //if((kinkAngle>maxDecAngpimu)&&(qT>0.12)&&(qT<0.30)&&((kink->GetR()>=120.)&&(kink->GetR()<=210.))&&(TMath::Abs(rapiditK)<0.7)&&(invariantMassKmu<0.8)){
607 if((kinkAngle>maxDecAngpimu)&&(qT>0.12)&&(qT<0.30)&&((kink->GetR()>= fKinkRadLow )&&(kink->GetR()<= fKinkRadUp ))&&(TMath::Abs(rapiditK)<0.7)&&(invariantMassKmu<0.8)){
608 // systematics if((kinkAngle>maxDecAngpimu)&&(qT>0.12)&&(qT<0.30)&&((kink->GetR()>=130.)&&(kink->GetR()<=200.))&&(TMath::Abs(rapiditK)<0.7)&&(invariantMassKmu<0.8)){
610 fAngMomKKinks->Fill(track->P(), kinkAngle);
611 fPtCut2 ->Fill(trackPt );
612 // maximum angles selection with some error cut
613 if( (kinkAngle<maxDecAngpimu*1.2) ) continue;
614 if ( (kinkAngle>maxDecAngKmu*.98) && ( track->P() >1.2 )) continue; ///5/5/2010
616 fPtCut3 ->Fill(trackPt );
617 // here the kaons selected by the decay features
618 fTPCSgnlPa->Fill( track->GetInnerParam()->GetP() ,(track->GetTPCsignal() ) ) ;
620 // NO dEdx cut test 9/2/13 if ( nsigma > 3.5) continue;
621 // system if ( nsigma > 4.0) continue; // gia systamatic error
623 // next plots for the identified kaons by the kink analysis
625 fHistPtKaon->Fill(track->Pt()); //
626 if(tpcSign >0.) fHistPtKaoP->Fill( track->Pt() ) ; //
627 if ( tpcSign <0.) fHistPtKaoN->Fill( track->Pt() ) ; //
628 fTPCSgnlP->Fill(track->GetInnerParam()->GetP(), (track->GetTPCsignal() ) ) ;
629 fRadNclCln->Fill( (kink->GetR()) ,(track->GetTPCclusters(0) ) ) ;
630 fRadiusPtcln->Fill( kink->GetR(), trackPt); //
631 fInvMassMuNuPt ->Fill(invariantMassKmu, trackPt);
633 // fTPCSgnlPtpc->Fill(trMomTPC , (track->GetTPCsignal() ) ) ;
634 //fMothKinkMomSgnl ->Fill(trMomTPCKink , (track->GetTPCsignal() ) ) ;
635 fMothKinkMomSgnl ->Fill( dEdxKinkDau , (track->GetTPCsignal() ) ) ;
637 fTPCSgnlKinkDau->Fill( daughterMKink.Mag(), dEdxKinkDau ) ; // daughter kink
639 fTPCMomNSgnl->Fill(trMomTPC ,nsigmall );
640 fNSigmTPC ->Fill(nsigmall );
642 frad->Fill(kink->GetR()); // kink
643 fradLK->Fill(lifeKink ); // kink
644 fHistEtaK->Fill(trackEta);
645 frapiKESD ->Fill(rapiditK); // rapidityof kaons
646 fPosiKinKBg->Fill( vposKink[2], kink->GetR() );
648 Float_t signPt= tpcSign*trackPt;
649 fSignPtNcl->Fill( signPt , tpcNCl ); /// 28/4/2010
650 fSignPtEta->Fill( signPt , rapiditK );
651 fEtaNcl->Fill( rapiditK, tpcNCl );
652 fSignPt->Fill( signPt );
653 fChi2NclTPC->Fill( (track->GetTPCchi2() ) , tpcNCl );
654 fRatChi2Ncl-> Fill ( (track->GetTPCchi2()/track->GetTPCclusters(0) )) ;
655 // fdcatoVxXY->Fill(dcaToVertexXYpos);
656 //if( dEdxKinkDau> 1.5* (track->GetTPCsignal() ) ) fKinkMothDau->Fill(track->P(),daughterMKink.Mag());
657 // if((dEdxKinkDau> 80. ) && (dEdxKinkDau > 4.*nsigmaPion) ) fKinkMothDau->Fill(track->P(),daughterMKink.Mag());
658 if (nsigmaPion> 3. ) fTPCSgnlPtpc->Fill( daughterMKink.Mag(), dEdxKinkDau ) ;
659 //if (TMath::Abs(dEdxKinkDau -(track->GetTPCsignal() )> 10. )) fTPCSgnlPtpc->Fill( daughterMKink.Mag(), dEdxKinkDau ) ;
660 flifetime->Fill(( lifeKink*.493667 ) /track->P() ) ;
661 fKinkKaon->Fill(track->Pt());
662 fDCAkink->Fill( Dist2 );
664 fPtKink->Fill(track->Pt()); /// K0 bins
665 if(tpcSign >0.) fPtKinkPos ->Fill( track->Pt() ) ; //K-plus bins K0
666 if(tpcSign <0.) fPtKinkNeg ->Fill( track->Pt() ) ; //K-minus bins K0
667 fKinKRbn->Fill(track->Pt()); // TOF
669 fradPtRpDt->Fill( kink->GetR(), 1./track->Pt(), rapiditK);
670 fAngMomK->Fill( track->P(), kinkAngle);
671 fPosiKinkK->Fill( vposKink[0], vposKink[1] );
672 fPosiKinKXZ->Fill( vposKink[2], vposKink[0] );
673 fPosiKinKYZ->Fill( vposKink[2], vposKink[1] );
678 } //End Kink Information
684 // fMultiplMC->Fill(nESDTracK );
686 PostData(1, fListOfHistos);
690 //________________________________________________________________________
691 void AliAnalysisKinkESDat::Terminate(Option_t *)
693 // Draw result to the screen
694 // Called once at the end of the query
698 //____________________________________________________________________//
701 const AliESDVertex* AliAnalysisKinkESDat::GetEventVertex(const AliESDEvent* esd) const
702 // Get the vertex from the ESD and returns it if the vertex is valid
707 const AliESDVertex* vertex = esd->GetPrimaryVertexTracks();
709 if((vertex->GetStatus()==kTRUE)) return vertex;
712 vertex = esd->GetPrimaryVertexSPD();
713 if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>0)) return vertex;