]>
Commit | Line | Data |
---|---|---|
9bc79c29 | 1 | #include <TCanvas.h> |
2 | #include <TTree.h> | |
3 | #include <TFile.h> | |
4 | #include <TH1F.h> | |
5 | #include <TH2F.h> | |
6 | #include <TH3F.h> | |
7 | #include <TPDGCode.h> | |
8 | #include <TDatabasePDG.h> | |
9 | #include <TClonesArray.h> | |
10 | #include <TROOT.h> | |
11 | ||
12 | #include "AliAODMCHeader.h" | |
13 | #include "AliAODMCParticle.h" | |
14 | ||
15 | #include "AliAODEvent.h" | |
16 | #include "AliAODv0.h" | |
17 | #include "AliAODcascade.h" | |
18 | ||
19 | #include "AliCentrality.h" | |
20 | ||
21 | #include "AliPID.h" | |
22 | #include "AliPIDResponse.h" | |
23 | #include "AliAODPid.h" | |
24 | ||
25 | #include "AliInputEventHandler.h" | |
26 | #include "AliAnalysisManager.h" | |
27 | ||
28 | #include "AliAnalysisTaskCTauPbPbaod.h" | |
29 | ||
2942f542 | 30 | //extern TROOT *gROOT; |
9bc79c29 | 31 | |
32 | ClassImp(AliAnalysisTaskCTauPbPbaod) | |
33 | ||
9bc79c29 | 34 | |
35 | // | |
36 | // This is a little task for checking the c*tau of the strange particles | |
37 | // | |
38 | ||
39 | AliAnalysisTaskCTauPbPbaod::AliAnalysisTaskCTauPbPbaod(const char *name) : | |
40 | AliAnalysisTaskSE(name), | |
41 | fIsMC(kFALSE), | |
42 | fCMin(0.), | |
43 | fCMax(90.), | |
ec62001b | 44 | fCPA(0.9975), |
9ffbd580 | 45 | fDCA(1.0), |
8356c5a0 | 46 | fTPCcr(70.), |
47 | fTPCcrfd(0.8), | |
e7ca8d73 | 48 | fDCApv(0.1), |
96d611c0 | 49 | fRmin(0.9), |
50 | fRmax(100.), | |
51 | ||
9bc79c29 | 52 | fOutput(0), |
53 | fMult(0), | |
ec62001b | 54 | fCosPA(0), |
9ffbd580 | 55 | fDtrDCA(0), |
8356c5a0 | 56 | fTPCrows(0), |
57 | fTPCratio(0), | |
e7ca8d73 | 58 | fPrimDCA(0), |
96d611c0 | 59 | fR(0), |
9bc79c29 | 60 | fdEdx(0), |
61 | fdEdxPid(0), | |
62 | ||
63 | fK0sM(0), | |
64 | fK0sSi(0), | |
65 | fK0sMC(0), | |
66 | fK0sAs(0), | |
67 | ||
68 | fLambdaM(0), | |
69 | fLambdaSi(0), | |
70 | fLambdaMC(0), | |
71 | fLambdaAs(0), | |
72 | ||
7d14e4a8 | 73 | fLambdaEff(0), |
74 | fLambdaPt(0), | |
75 | ||
5bd3c9e4 | 76 | fLambdaBarM(0), |
77 | fLambdaBarSi(0), | |
78 | fLambdaBarMC(0), | |
79 | fLambdaBarAs(0), | |
80 | ||
7d14e4a8 | 81 | fLambdaBarEff(0), |
82 | fLambdaBarPt(0), | |
9bc79c29 | 83 | |
84 | fLambdaFromXi(0), | |
85 | fXiM(0), | |
5bd3c9e4 | 86 | fXiSiP(0), |
87 | ||
88 | fLambdaBarFromXiBar(0), | |
89 | fXiBarM(0), | |
90 | fXiBarSiP(0) | |
9bc79c29 | 91 | { |
92 | // Constructor. Initialization of pointers | |
93 | DefineOutput(1, TList::Class()); | |
94 | } | |
95 | ||
96 | void AliAnalysisTaskCTauPbPbaod::UserCreateOutputObjects() | |
97 | { | |
7d14e4a8 | 98 | Int_t nbins=100; // number of bins |
99 | Double_t ltMax=100.; | |
100 | Double_t ptMax=10.; | |
101 | ||
9bc79c29 | 102 | fOutput = new TList(); |
103 | fOutput->SetOwner(); | |
104 | ||
105 | ||
106 | fMult=new TH1F("fMult","Multiplicity",1100,0.,3300); | |
107 | fMult->GetXaxis()->SetTitle("N tracks"); | |
108 | fOutput->Add(fMult); | |
109 | ||
ec62001b | 110 | fCosPA=new TH1F("fCosPA","Cos(PA) distribution",50,0.9975,1.0005); |
111 | fCosPA->GetXaxis()->SetTitle("Cos(PA)"); | |
112 | fOutput->Add(fCosPA); | |
113 | ||
9ffbd580 | 114 | fDtrDCA=new TH1F("fDtrDCA","DCA between V0 daughters",50,0.0,1.5); |
115 | fDtrDCA->GetXaxis()->SetTitle("DCA (rel. u.)"); | |
116 | fOutput->Add(fDtrDCA); | |
117 | ||
8356c5a0 | 118 | fTPCrows=new TH1F("fTPCrows","TPC crossed pad rows",180,0.,180.); |
119 | fTPCrows->GetXaxis()->SetTitle("TPC crossed pad rows"); | |
120 | fOutput->Add(fTPCrows); | |
121 | ||
122 | fTPCratio=new TH1F("fTPCratio","TPC crossed/findable pad rows",50,0.0,1.5); | |
123 | fTPCratio->GetXaxis()->SetTitle("TPC crossed/findable pad rows"); | |
124 | fOutput->Add(fTPCratio); | |
125 | ||
e7ca8d73 | 126 | fPrimDCA=new TH1F("fPrimDCA","DCA wrt the primary vertex",50,0.0,1.5); |
127 | fPrimDCA->GetXaxis()->SetTitle("DCA wrt the PV (cm)"); | |
128 | fOutput->Add(fPrimDCA); | |
129 | ||
96d611c0 | 130 | fR=new TH1F("fR","Radius of the V0 vertices",101,0.0,102); |
131 | fR->GetXaxis()->SetTitle("R (cm)"); | |
132 | fOutput->Add(fR); | |
133 | ||
9bc79c29 | 134 | fdEdx=new TH2F("fdEdx","dE/dx",50,0.2,3,50,0.,6.); |
135 | fOutput->Add(fdEdx); | |
136 | ||
137 | fdEdxPid=new TH2F("fdEdxPid","dE/dx with PID",50,0.2,3,50,0.,6.); | |
138 | fOutput->Add(fdEdxPid); | |
139 | ||
140 | fK0sM = | |
7d14e4a8 | 141 | new TH2F("fK0sM", "Mass for K^{0}_{s}", nbins/2,0.448,0.548,nbins,0.,ptMax); |
9bc79c29 | 142 | fK0sM->GetXaxis()->SetTitle("Mass (GeV/c)"); |
143 | fOutput->Add(fK0sM); | |
144 | ||
145 | fK0sSi = | |
146 | new TH2F("fK0sSi","L_{T} vs p_{T} for K^{0}_{s}, side-band subtracted", | |
7d14e4a8 | 147 | nbins,0.,ptMax,nbins,0.,ltMax); |
9bc79c29 | 148 | fK0sSi->GetXaxis()->SetTitle("p_{T} (GeV/c)"); |
149 | fK0sSi->GetYaxis()->SetTitle("L_{T} (cm)"); | |
150 | fOutput->Add(fK0sSi); | |
151 | ||
152 | fK0sMC = | |
153 | new TH2F("fK0sMC","L_{T} vs p_{T} for K^{0}_{s}, from MC stack", | |
7d14e4a8 | 154 | nbins,0.,ptMax,nbins,0.,ltMax); |
9bc79c29 | 155 | fK0sMC->GetXaxis()->SetTitle("p_{T} (GeV/c)"); |
156 | fK0sMC->GetYaxis()->SetTitle("L_{T} (cm)"); | |
157 | fOutput->Add(fK0sMC); | |
158 | ||
159 | fK0sAs = | |
160 | new TH2F("fK0sAs", "L_{T} vs p_{T} for K^{0}_{s}, associated", | |
7d14e4a8 | 161 | nbins,0.,ptMax,nbins,0.,ltMax); |
9bc79c29 | 162 | fK0sAs->GetXaxis()->SetTitle("p_{T} (GeV/c)"); |
163 | fK0sAs->GetYaxis()->SetTitle("L_{T} (cm)"); | |
164 | fOutput->Add(fK0sAs); | |
165 | ||
166 | //---------------------- | |
167 | ||
168 | fLambdaM = | |
7d14e4a8 | 169 | new TH2F("fLambdaM","Mass for \\Lambda", nbins, 1.065, 1.165,nbins,0.,ptMax); |
9bc79c29 | 170 | fLambdaM->GetXaxis()->SetTitle("Mass (GeV/c)"); |
171 | fOutput->Add(fLambdaM); | |
172 | ||
173 | fLambdaSi = | |
174 | new TH2F("fLambdaSi","L_{T} vs p_{T} for \\Lambda, side-band subtructed", | |
7d14e4a8 | 175 | nbins,0.,ptMax,nbins,0.,ltMax); |
9bc79c29 | 176 | fLambdaSi->GetXaxis()->SetTitle("p_{T} (GeV/c)"); |
177 | fLambdaSi->GetYaxis()->SetTitle("L_{T} (cm)"); | |
178 | fOutput->Add(fLambdaSi); | |
179 | ||
180 | fLambdaMC = | |
181 | new TH2F("fLambdaMC","L_{T} vs p_{T} for \\Lambda, from MC stack", | |
7d14e4a8 | 182 | nbins,0.,ptMax,nbins,0.,ltMax); |
9bc79c29 | 183 | fLambdaMC->GetXaxis()->SetTitle("p_{T} (GeV/c)"); |
184 | fLambdaMC->GetYaxis()->SetTitle("L_{T} (cm)"); | |
185 | fOutput->Add(fLambdaMC); | |
186 | ||
187 | fLambdaAs = | |
188 | new TH2F("fLambdaAs","L_{T} vs p_{T} for \\Lambda, associated", | |
7d14e4a8 | 189 | nbins,0.,ptMax,nbins,0.,ltMax); |
9bc79c29 | 190 | fLambdaAs->GetXaxis()->SetTitle("p_{T} (GeV/c)"); |
191 | fLambdaAs->GetYaxis()->SetTitle("L_{T} (cm)"); | |
192 | fOutput->Add(fLambdaAs); | |
193 | ||
7d14e4a8 | 194 | //---------------------- |
195 | ||
196 | fLambdaEff=fLambdaAs->ProjectionX(); | |
197 | fLambdaEff->SetName("fLambdaEff"); | |
198 | fLambdaEff->SetTitle("Efficiency for #Lambda"); | |
199 | fOutput->Add(fLambdaEff); | |
200 | ||
201 | fLambdaPt=fLambdaAs->ProjectionX(); | |
202 | fLambdaPt->SetName("fLambdaPt"); | |
203 | fLambdaPt->SetTitle("Raw #Lambda pT spectrum"); | |
204 | fOutput->Add(fLambdaPt); | |
205 | ||
206 | //---------------------- | |
5bd3c9e4 | 207 | |
208 | fLambdaBarM = | |
7d14e4a8 | 209 | new TH2F("fLambdaBarM","Mass for anti-\\Lambda", nbins, 1.065, 1.165,nbins,0.,ptMax); |
5bd3c9e4 | 210 | fLambdaBarM->GetXaxis()->SetTitle("Mass (GeV/c)"); |
211 | fOutput->Add(fLambdaBarM); | |
212 | ||
213 | fLambdaBarSi = | |
214 | new TH2F("fLambdaBarSi","L_{T} vs p_{T} for anti-\\Lambda, side-band subtructed", | |
7d14e4a8 | 215 | nbins,0.,ptMax,nbins,0.,ltMax); |
5bd3c9e4 | 216 | fLambdaBarSi->GetXaxis()->SetTitle("p_{T} (GeV/c)"); |
217 | fLambdaBarSi->GetYaxis()->SetTitle("L_{T} (cm)"); | |
218 | fOutput->Add(fLambdaBarSi); | |
219 | ||
220 | fLambdaBarMC = | |
221 | new TH2F("fLambdaBarMC","L_{T} vs p_{T} for anti-\\Lambda, from MC stack", | |
7d14e4a8 | 222 | nbins,0.,ptMax,nbins,0.,ltMax); |
5bd3c9e4 | 223 | fLambdaBarMC->GetXaxis()->SetTitle("p_{T} (GeV/c)"); |
224 | fLambdaBarMC->GetYaxis()->SetTitle("L_{T} (cm)"); | |
225 | fOutput->Add(fLambdaBarMC); | |
226 | ||
227 | fLambdaBarAs = | |
228 | new TH2F("fLambdaBarAs","L_{T} vs p_{T} for anti-\\Lambda, associated", | |
7d14e4a8 | 229 | nbins,0.,ptMax,nbins,0.,ltMax); |
5bd3c9e4 | 230 | fLambdaBarAs->GetXaxis()->SetTitle("p_{T} (GeV/c)"); |
231 | fLambdaBarAs->GetYaxis()->SetTitle("L_{T} (cm)"); | |
232 | fOutput->Add(fLambdaBarAs); | |
233 | ||
234 | ||
7d14e4a8 | 235 | //---------------------- |
5bd3c9e4 | 236 | |
7d14e4a8 | 237 | fLambdaBarEff=fLambdaBarAs->ProjectionX(); |
238 | fLambdaBarEff->SetName("fLambdaBarEff"); | |
239 | fLambdaBarEff->SetTitle("Efficiency for anti-#Lambda"); | |
240 | fOutput->Add(fLambdaBarEff); | |
9bc79c29 | 241 | |
7d14e4a8 | 242 | fLambdaBarPt=fLambdaBarAs->ProjectionX(); |
243 | fLambdaBarPt->SetName("fLambdaBarPt"); | |
244 | fLambdaBarPt->SetTitle("Raw anti-#Lambda pT spectrum"); | |
245 | fOutput->Add(fLambdaBarPt); | |
9bc79c29 | 246 | |
247 | //---------------------- | |
248 | ||
5bd3c9e4 | 249 | fLambdaFromXi=new TH3F("fLambdaFromXi","L_{T} vs p_{T} vs p_{T} of \\Xi for \\Lambda from \\Xi", |
7d14e4a8 | 250 | nbins,0.,ptMax,nbins,0.,ltMax,33,0.,ptMax+2); |
9bc79c29 | 251 | fOutput->Add(fLambdaFromXi); |
252 | ||
253 | fXiM = | |
7d14e4a8 | 254 | new TH2F("fXiM", "\\Xi mass distribution", 50, 1.271, 1.371,33,0.,ptMax+2); |
9bc79c29 | 255 | fOutput->Add(fXiM); |
256 | ||
257 | fXiSiP = new TH1F("fXiSiP", "Pt for \\Xi, side-band subracted", | |
7d14e4a8 | 258 | 33,0.,ptMax+2); |
9bc79c29 | 259 | fOutput->Add(fXiSiP); |
260 | ||
261 | ||
5bd3c9e4 | 262 | fLambdaBarFromXiBar=new TH3F("fLambdaBarFromXiBar","L_{T} vs p_{T} vs p_{T} of anti-\\Xi for anti-\\Lambda from anti-\\Xi", |
7d14e4a8 | 263 | nbins,0.,ptMax,nbins,0.,ltMax,33,0.,ptMax+2); |
6105beb1 | 264 | fOutput->Add(fLambdaBarFromXiBar); |
5bd3c9e4 | 265 | |
266 | fXiBarM = | |
7d14e4a8 | 267 | new TH2F("fXiBarM", "anti-\\Xi mass distribution", 50, 1.271, 1.371,33,0.,ptMax+2); |
5bd3c9e4 | 268 | fOutput->Add(fXiBarM); |
269 | ||
270 | fXiBarSiP = new TH1F("fXiBarSiP", "Pt for anti-\\Xi, side-band subracted", | |
7d14e4a8 | 271 | 33,0.,ptMax+2); |
5bd3c9e4 | 272 | fOutput->Add(fXiBarSiP); |
273 | ||
274 | ||
9bc79c29 | 275 | PostData(1, fOutput); |
276 | } | |
277 | ||
8356c5a0 | 278 | Bool_t AliAnalysisTaskCTauPbPbaod::AcceptTrack(const AliAODTrack *t) { |
9bc79c29 | 279 | if (!t->IsOn(AliAODTrack::kTPCrefit)) return kFALSE; |
280 | //if (t->GetKinkIndex(0)>0) return kFALSE; | |
281 | ||
282 | Float_t nCrossedRowsTPC = t->GetTPCClusterInfo(2,1); | |
8356c5a0 | 283 | if (nCrossedRowsTPC < fTPCcr) return kFALSE; |
9bc79c29 | 284 | Int_t findable=t->GetTPCNclsF(); |
285 | if (findable <= 0) return kFALSE; | |
8356c5a0 | 286 | if (nCrossedRowsTPC/findable < fTPCcrfd) return kFALSE; |
9bc79c29 | 287 | |
ff909752 | 288 | if (TMath::Abs(t->Eta()) > 0.8) return kFALSE; |
289 | ||
8356c5a0 | 290 | fTPCrows->Fill(nCrossedRowsTPC); |
291 | fTPCratio->Fill(nCrossedRowsTPC/findable); | |
292 | ||
9bc79c29 | 293 | return kTRUE; |
294 | } | |
295 | ||
ec62001b | 296 | Bool_t |
297 | AliAnalysisTaskCTauPbPbaod::AcceptV0(const AliAODv0 *v0,const AliAODEvent *aod) | |
298 | { | |
7d14e4a8 | 299 | if (v0->GetOnFlyStatus()) return kFALSE; |
9bc79c29 | 300 | |
ec62001b | 301 | Double_t cpa=v0->CosPointingAngle(aod->GetPrimaryVertex()); |
302 | if (cpa < fCPA) return kFALSE; | |
303 | ||
9ffbd580 | 304 | Double_t dca=v0->DcaV0Daughters(); |
305 | if (dca > fDCA) return kFALSE; | |
306 | ||
7d14e4a8 | 307 | const AliAODTrack *ntrack=(AliAODTrack *)v0->GetDaughter(1); |
308 | if (!AcceptTrack(ntrack)) return kFALSE; | |
9bc79c29 | 309 | |
7d14e4a8 | 310 | const AliAODTrack *ptrack=(AliAODTrack *)v0->GetDaughter(0); |
311 | if (!AcceptTrack(ptrack)) return kFALSE; | |
9bc79c29 | 312 | |
e7ca8d73 | 313 | Float_t xyn=v0->DcaNegToPrimVertex(); |
314 | if (TMath::Abs(xyn)<fDCApv) return kFALSE; | |
315 | Float_t xyp=v0->DcaPosToPrimVertex(); | |
316 | if (TMath::Abs(xyp)<fDCApv) return kFALSE; | |
9bc79c29 | 317 | |
7d14e4a8 | 318 | Double_t xyz[3]; v0->GetSecondaryVtx(xyz); |
9bc79c29 | 319 | Double_t r2=xyz[0]*xyz[0] + xyz[1]*xyz[1]; |
96d611c0 | 320 | if (r2<fRmin*fRmin) return kFALSE; |
321 | if (r2>fRmax*fRmax) return kFALSE; | |
9bc79c29 | 322 | |
8356c5a0 | 323 | fCosPA->Fill(cpa); |
324 | fDtrDCA->Fill(dca); | |
e7ca8d73 | 325 | fPrimDCA->Fill(xyn); fPrimDCA->Fill(xyp); |
96d611c0 | 326 | fR->Fill(TMath::Sqrt(r2)); |
8356c5a0 | 327 | |
9bc79c29 | 328 | return kTRUE; |
329 | } | |
330 | ||
8356c5a0 | 331 | Bool_t AliAnalysisTaskCTauPbPbaod::AcceptCascade(const AliAODcascade *cs, const AliAODEvent *aod) { |
9bc79c29 | 332 | |
9bc79c29 | 333 | AliAODVertex *xi = cs->GetDecayVertexXi(); |
334 | const AliAODTrack *bach=(AliAODTrack *)xi->GetDaughter(0); | |
335 | if (!AcceptTrack(bach)) return kFALSE; | |
336 | ||
337 | Double_t xy=cs->DcaBachToPrimVertex(); | |
338 | if (TMath::Abs(xy) < 0.03) return kFALSE; | |
339 | ||
340 | const AliAODVertex *vtx=aod->GetPrimaryVertex(); | |
341 | Double_t xv=vtx->GetX(), yv=vtx->GetY(), zv=vtx->GetZ(); | |
342 | Double_t cpa=cs->CosPointingAngleXi(xv,yv,zv); | |
343 | if (cpa<0.999) return kFALSE; | |
344 | ||
345 | if (cs->DcaXiDaughters() > 0.3) return kFALSE; | |
346 | ||
347 | return kTRUE; | |
348 | } | |
349 | ||
1c102c9a | 350 | static Bool_t AcceptPID(const AliPIDResponse *pidResponse, |
351 | const AliAODTrack *ptrack, const TClonesArray *stack) { | |
352 | ||
f6357875 | 353 | const AliAODPid *pid=ptrack->GetDetPid(); |
354 | if (!pid) return kTRUE; | |
355 | if (pid->GetTPCmomentum() > 1.) return kTRUE; | |
1c102c9a | 356 | |
357 | if (stack) { | |
358 | // MC PID | |
1c102c9a | 359 | Int_t plab=TMath::Abs(ptrack->GetLabel()); |
7d14e4a8 | 360 | AliAODMCParticle *pp=(AliAODMCParticle*)(*stack)[plab]; |
361 | if (!pp) return kTRUE; | |
362 | if (pp->Charge() > 0) { | |
363 | if (pp->GetPdgCode() == kProton) return kTRUE; | |
364 | } else { | |
365 | if (pp->GetPdgCode() == kProtonBar) return kTRUE; | |
366 | } | |
1c102c9a | 367 | } else { |
368 | // Real PID | |
f6357875 | 369 | Double_t nsig=pidResponse->NumberOfSigmasTPC(ptrack,AliPID::kProton); |
370 | if (TMath::Abs(nsig) < 3.) return kTRUE; | |
1c102c9a | 371 | } |
372 | ||
f6357875 | 373 | return kFALSE; |
1c102c9a | 374 | } |
375 | ||
7d14e4a8 | 376 | static AliAODMCParticle* |
377 | AssociateV0(const AliAODTrack *ptrack, const AliAODTrack *ntrack, | |
378 | const TClonesArray *stack, AliAODMCParticle *&mcp) { | |
379 | // | |
380 | // Try to associate the V0 with the daughters ptrack and ntrack | |
381 | // with the Monte Carlo | |
382 | // | |
383 | if (!stack) return 0; | |
384 | ||
385 | Int_t nlab=TMath::Abs(ntrack->GetLabel()); | |
386 | AliAODMCParticle *n=(AliAODMCParticle*)(*stack)[nlab]; | |
387 | if (!n) return 0; | |
388 | ||
389 | Int_t plab=TMath::Abs(ptrack->GetLabel()); | |
390 | AliAODMCParticle *p=(AliAODMCParticle*)(*stack)[plab]; | |
391 | if (!p) return 0; | |
392 | ||
393 | Int_t imp=p->GetMother(); //V0 particle, the mother of the pos. track | |
394 | AliAODMCParticle *p0=(AliAODMCParticle*)(*stack)[imp]; | |
395 | if (!p0) return 0; | |
396 | ||
397 | Int_t imn=n->GetMother(); | |
398 | if (imp != imn) { // Check decays of the daughters | |
399 | return 0; // Fixme | |
400 | } | |
401 | ||
402 | Int_t code=p0->GetPdgCode(); | |
403 | if (code != kK0Short) | |
404 | if (code != kLambda0) | |
405 | if (code != kLambda0Bar) return 0; | |
406 | ||
407 | mcp=p; | |
408 | ||
409 | return p0; | |
410 | } | |
411 | ||
412 | ||
9bc79c29 | 413 | void AliAnalysisTaskCTauPbPbaod::UserExec(Option_t *) |
414 | { | |
7d14e4a8 | 415 | const Double_t yMax=0.5; |
416 | const Double_t pMin=0.0; | |
417 | const Double_t lMax=0.001; | |
9bc79c29 | 418 | |
419 | AliAODEvent *aod=(AliAODEvent *)InputEvent(); | |
420 | ||
421 | if (!aod) { | |
422 | Printf("ERROR: aod not available"); | |
423 | return; | |
424 | } | |
425 | ||
5bd3c9e4 | 426 | // Vertex selection |
427 | const AliAODVertex *vtx=aod->GetPrimaryVertex(); | |
428 | if (vtx->GetNContributors()<3) return; | |
429 | Double_t xv=vtx->GetX(), yv=vtx->GetY(), zv=vtx->GetZ(); | |
430 | ||
431 | if (TMath::Abs(zv) > 10.) return ; | |
432 | ||
9bc79c29 | 433 | |
434 | // Physics selection | |
435 | AliAnalysisManager *mgr= AliAnalysisManager::GetAnalysisManager(); | |
436 | AliInputEventHandler *hdr=(AliInputEventHandler*)mgr->GetInputEventHandler(); | |
437 | UInt_t maskIsSelected = hdr->IsEventSelected(); | |
438 | Bool_t isSelected = (maskIsSelected & AliVEvent::kMB); | |
439 | if (!isSelected) return; | |
440 | ||
5bd3c9e4 | 441 | fMult->Fill(-100); //event counter |
442 | ||
9bc79c29 | 443 | // Centrality selection |
444 | AliCentrality *cent=aod->GetCentrality(); | |
445 | if (!cent->IsEventInCentralityClass(fCMin,fCMax,"V0M")) return; | |
446 | ||
9bc79c29 | 447 | AliPIDResponse *pidResponse = hdr->GetPIDResponse(); |
448 | ||
5bd3c9e4 | 449 | |
9bc79c29 | 450 | |
451 | //+++++++ MC | |
452 | TClonesArray *stack = 0x0; | |
453 | Double_t mcXv=0., mcYv=0., mcZv=0.; | |
454 | ||
455 | if (fIsMC) { | |
456 | TList *lst = aod->GetList(); | |
457 | stack = (TClonesArray*)lst->FindObject(AliAODMCParticle::StdBranchName()); | |
458 | if (!stack) { | |
459 | Printf("ERROR: stack not available"); | |
460 | return; | |
461 | } | |
462 | AliAODMCHeader * | |
463 | mcHdr=(AliAODMCHeader*)lst->FindObject(AliAODMCHeader::StdBranchName()); | |
464 | ||
465 | mcXv=mcHdr->GetVtxX(); mcYv=mcHdr->GetVtxY(); mcZv=mcHdr->GetVtxZ(); | |
466 | ||
7d14e4a8 | 467 | Int_t ntrk=stack->GetEntriesFast(); |
468 | while (ntrk--) { | |
469 | AliAODMCParticle *p0=(AliAODMCParticle*)stack->UncheckedAt(ntrk); | |
9bc79c29 | 470 | Int_t code=p0->GetPdgCode(); |
471 | if (code != kK0Short) | |
5bd3c9e4 | 472 | if (code != kLambda0) |
473 | if (code != kLambda0Bar) continue; | |
9bc79c29 | 474 | |
475 | Int_t plab=p0->GetDaughter(0), nlab=p0->GetDaughter(1); | |
476 | if (nlab==plab) continue; | |
7d14e4a8 | 477 | AliAODMCParticle *part=(AliAODMCParticle*)(*stack)[plab]; |
9bc79c29 | 478 | if (!part) continue; |
479 | Double_t charge=part->Charge(); | |
480 | if (charge == 0.) continue; | |
481 | ||
482 | Double_t pt=p0->Pt(); | |
483 | if (pt<pMin) continue; | |
484 | if (TMath::Abs(p0->Y())>yMax) continue; | |
485 | ||
486 | Double_t x=p0->Xv(), y=p0->Yv(), z=p0->Zv(); | |
52345783 | 487 | Double_t dxmc=mcXv-x, dymc=mcYv-y, dzmc=mcZv-z; |
488 | Double_t l=TMath::Sqrt(dxmc*dxmc + dymc*dymc + dzmc*dzmc); | |
9bc79c29 | 489 | |
7d14e4a8 | 490 | if (l > lMax) continue; // secondary V0 |
9bc79c29 | 491 | |
492 | x=part->Xv(); y=part->Yv(); | |
52345783 | 493 | dxmc=mcXv-x; dymc=mcYv-y; |
494 | Double_t lt=TMath::Sqrt(dxmc*dxmc + dymc*dymc); | |
9bc79c29 | 495 | |
7d14e4a8 | 496 | switch (code) { |
497 | case kK0Short: | |
9bc79c29 | 498 | fK0sMC->Fill(pt,lt); |
7d14e4a8 | 499 | break; |
500 | case kLambda0: | |
9bc79c29 | 501 | fLambdaMC->Fill(pt,lt); |
7d14e4a8 | 502 | break; |
503 | case kLambda0Bar: | |
5bd3c9e4 | 504 | fLambdaBarMC->Fill(pt,lt); |
7d14e4a8 | 505 | break; |
506 | default: break; | |
5bd3c9e4 | 507 | } |
9bc79c29 | 508 | } |
509 | } | |
7d14e4a8 | 510 | //+++++++ |
9bc79c29 | 511 | |
512 | ||
7d14e4a8 | 513 | Int_t ntrk1=aod->GetNumberOfTracks(); |
9bc79c29 | 514 | Int_t mult=0; |
7d14e4a8 | 515 | for (Int_t i=0; i<ntrk1; i++) { |
9bc79c29 | 516 | AliAODTrack *t=aod->GetTrack(i); |
517 | if (t->IsMuonTrack()) continue; | |
518 | if (!t->IsOn(AliAODTrack::kTPCrefit)) continue; | |
519 | ||
520 | Double_t xyz[3]; | |
521 | if (t->GetPosition(xyz)) continue; | |
522 | if (TMath::Abs(xyz[0])>3.) continue; | |
523 | if (TMath::Abs(xyz[1])>3.) continue; | |
524 | ||
525 | Double_t pt=t->Pt(),pz=t->Pz(); | |
526 | if (TMath::Abs(pz/pt)>0.8) continue; | |
527 | ||
528 | mult++; | |
529 | ||
530 | const AliAODPid *pid=t->GetDetPid(); | |
531 | if (!pid) continue; | |
532 | ||
533 | Double_t p=pid->GetTPCmomentum(); | |
534 | Double_t dedx=pid->GetTPCsignal()/47.; | |
535 | fdEdx->Fill(p,dedx,1); | |
536 | ||
1c102c9a | 537 | Double_t nsig=pidResponse->NumberOfSigmasTPC(t,AliPID::kProton); |
9bc79c29 | 538 | if (TMath::Abs(nsig) < 3.) fdEdxPid->Fill(p,dedx,1); |
539 | ||
540 | } | |
541 | fMult->Fill(mult); | |
542 | ||
543 | ||
544 | Int_t nv0 = aod->GetNumberOfV0s(); | |
545 | while (nv0--) { | |
546 | AliAODv0 *v0=aod->GetV0(nv0); | |
547 | ||
7d14e4a8 | 548 | Double_t pt=TMath::Sqrt(v0->Pt2V0()); |
549 | if (pt < pMin) continue; | |
550 | ||
9bc79c29 | 551 | if (!AcceptV0(v0,aod)) continue; |
552 | ||
553 | const AliAODTrack *ntrack=(AliAODTrack *)v0->GetDaughter(1); | |
554 | const AliAODTrack *ptrack=(AliAODTrack *)v0->GetDaughter(0); | |
555 | ||
556 | Double_t xyz[3]; v0->GetSecondaryVtx(xyz); | |
52345783 | 557 | Double_t dxv=xyz[0]-xv, dyv=xyz[1]-yv; |
558 | Double_t lt=TMath::Sqrt(dxv*dxv + dyv*dyv); | |
9bc79c29 | 559 | |
7d14e4a8 | 560 | if (lt/pt > 3*7.89/1.1157) continue; |
1c102c9a | 561 | |
7d14e4a8 | 562 | //--- V0 switches |
563 | Bool_t isK0s=kTRUE; | |
564 | Bool_t isLambda=kTRUE; | |
565 | Bool_t isLambdaBar=kTRUE; | |
9bc79c29 | 566 | |
7d14e4a8 | 567 | if (0.4977*lt/pt > 3*2.68) isK0s=kFALSE; |
568 | if (1.1157*lt/pt > 3*7.89) isLambdaBar=isLambda=kFALSE; | |
1c102c9a | 569 | |
7d14e4a8 | 570 | if (!AcceptPID(pidResponse, ptrack, stack)) isLambda=kFALSE; |
571 | if (!AcceptPID(pidResponse, ntrack, stack)) isLambdaBar=kFALSE; | |
9bc79c29 | 572 | |
7d14e4a8 | 573 | Double_t yK0s=TMath::Abs(v0->RapK0Short()); |
574 | Double_t yLam=TMath::Abs(v0->RapLambda()); | |
575 | if (yK0s > yMax) isK0s=kFALSE; | |
576 | if (yLam > yMax) isLambda=isLambdaBar=kFALSE; | |
577 | //--- | |
9bc79c29 | 578 | |
579 | Double_t mass=0., m=0., s=0.; | |
7d14e4a8 | 580 | if (isK0s) { |
9bc79c29 | 581 | mass=v0->MassK0Short(); |
582 | fK0sM->Fill(mass,pt); | |
583 | ||
584 | m=TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass(); | |
585 | s=0.0044 + (0.008-0.0044)/(10-1)*(pt - 1.); | |
586 | if (TMath::Abs(m-mass) < 3*s) { | |
587 | fK0sSi->Fill(pt,lt); | |
7d14e4a8 | 588 | } else { |
589 | isK0s=kFALSE; | |
9bc79c29 | 590 | } |
591 | if (TMath::Abs(m-mass + 4.5*s) < 1.5*s) { | |
592 | fK0sSi->Fill(pt,lt,-1); | |
593 | } | |
594 | if (TMath::Abs(m-mass - 4.5*s) < 1.5*s) { | |
595 | fK0sSi->Fill(pt,lt,-1); | |
596 | } | |
597 | } | |
598 | ||
7d14e4a8 | 599 | if (isLambda) { |
9bc79c29 | 600 | mass=v0->MassLambda(); |
601 | fLambdaM->Fill(mass,pt); | |
602 | ||
603 | m=TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(); | |
604 | //s=0.0027 + (0.004-0.0027)/(10-1)*(pt-1); | |
605 | //s=0.0015 + (0.002-0.0015)/(2.6-1)*(pt-1); | |
606 | s=0.0023 + (0.004-0.0023)/(6-1)*(pt-1); | |
607 | if (TMath::Abs(m-mass) < 3*s) { | |
608 | fLambdaSi->Fill(pt,lt); | |
7d14e4a8 | 609 | } else { |
610 | isLambda=kFALSE; | |
611 | } | |
9bc79c29 | 612 | if (TMath::Abs(m-mass + 4.5*s) < 1.5*s) { |
613 | fLambdaSi->Fill(pt,lt,-1); | |
9bc79c29 | 614 | } |
615 | if (TMath::Abs(m-mass - 4.5*s) < 1.5*s) { | |
616 | fLambdaSi->Fill(pt,lt,-1); | |
9bc79c29 | 617 | } |
618 | } | |
5bd3c9e4 | 619 | |
7d14e4a8 | 620 | if (isLambdaBar) { |
5bd3c9e4 | 621 | mass=v0->MassAntiLambda(); |
622 | fLambdaBarM->Fill(mass,pt); | |
623 | ||
624 | m=TDatabasePDG::Instance()->GetParticle(kLambda0Bar)->Mass(); | |
625 | //s=0.0027 + (0.004-0.0027)/(10-1)*(pt-1); | |
626 | //s=0.0015 + (0.002-0.0015)/(2.6-1)*(pt-1); | |
627 | s=0.0023 + (0.004-0.0023)/(6-1)*(pt-1); | |
628 | if (TMath::Abs(m-mass) < 3*s) { | |
629 | fLambdaBarSi->Fill(pt,lt); | |
7d14e4a8 | 630 | } else { |
631 | isLambdaBar=kFALSE; | |
5bd3c9e4 | 632 | } |
633 | if (TMath::Abs(m-mass + 4.5*s) < 1.5*s) { | |
634 | fLambdaBarSi->Fill(pt,lt,-1); | |
5bd3c9e4 | 635 | } |
636 | if (TMath::Abs(m-mass - 4.5*s) < 1.5*s) { | |
637 | fLambdaBarSi->Fill(pt,lt,-1); | |
5bd3c9e4 | 638 | } |
639 | } | |
7d14e4a8 | 640 | |
641 | if (!fIsMC) continue; | |
642 | ||
643 | //++++++ MC | |
644 | if (!isK0s) | |
645 | if (!isLambda) | |
646 | if (!isLambdaBar) continue;//check MC only for the accepted V0s | |
647 | ||
648 | AliAODMCParticle *mcp=0; | |
649 | AliAODMCParticle *mc0=AssociateV0(ptrack,ntrack,stack,mcp); | |
650 | if (!mc0) continue; | |
651 | ||
652 | Double_t ptAs=mc0->Pt(); | |
653 | if (ptAs < pMin) continue; | |
654 | Double_t yAs=mc0->Y(); | |
655 | if (TMath::Abs(yAs) > yMax) continue; | |
656 | ||
657 | Int_t code=mc0->GetPdgCode(); | |
658 | ||
659 | Double_t dx = mcXv - mcp->Xv(), dy = mcYv - mcp->Yv(); | |
660 | Double_t ltAs=TMath::Sqrt(dx*dx + dy*dy); | |
661 | ||
662 | Double_t dz=mcZv - mc0->Zv(); dy=mcYv - mc0->Yv(); dx=mcXv - mc0->Xv(); | |
663 | Double_t l = TMath::Sqrt(dx*dx + dy*dy + dz*dz); | |
664 | if (l<lMax) { // Primary V0s | |
665 | switch (code) { | |
666 | case kK0Short: | |
667 | if (isK0s) fK0sAs->Fill(ptAs,ltAs); | |
668 | break; | |
669 | case kLambda0: | |
670 | if (isLambda) fLambdaAs->Fill(ptAs,ltAs); | |
671 | break; | |
672 | case kLambda0Bar: | |
673 | if (isLambdaBar) fLambdaBarAs->Fill(ptAs,ltAs); | |
674 | break; | |
675 | default: break; | |
676 | } | |
677 | } else { | |
678 | if (code==kK0Short) continue; | |
679 | ||
680 | Int_t nx=mc0->GetMother(); | |
681 | AliAODMCParticle *xi=(AliAODMCParticle*)(*stack)[nx]; | |
682 | if (!xi) continue; | |
683 | Int_t xcode=xi->GetPdgCode(); | |
684 | ||
685 | switch (code) { | |
686 | case kLambda0: | |
687 | if (isLambda) | |
688 | if ((xcode==kXiMinus) || (xcode==3322)) | |
689 | fLambdaFromXi->Fill(ptAs,ltAs,xi->Pt()); | |
690 | break; | |
691 | case kLambda0Bar: | |
692 | if (isLambdaBar) | |
693 | if ((xcode==kXiPlusBar)||(xcode==-3322)) | |
694 | fLambdaBarFromXiBar->Fill(ptAs,ltAs,xi->Pt()); | |
695 | break; | |
696 | default: break; | |
697 | } | |
698 | } | |
699 | //++++++ | |
700 | ||
9bc79c29 | 701 | } |
702 | ||
703 | Int_t ncs=aod->GetNumberOfCascades(); | |
704 | for (Int_t i=0; i<ncs; i++) { | |
705 | AliAODcascade *cs=aod->GetCascade(i); | |
706 | ||
7d14e4a8 | 707 | Double_t pt=TMath::Sqrt(cs->Pt2Xi()); |
708 | if (pt < pMin) continue; | |
9bc79c29 | 709 | if (TMath::Abs(cs->RapXi()) > yMax) continue; |
710 | if (!AcceptCascade(cs,aod)) continue; | |
711 | ||
712 | const AliAODv0 *v0=(AliAODv0*)cs; | |
713 | if (v0->RapLambda() > yMax) continue; | |
714 | if (!AcceptV0(v0,aod)) continue; | |
715 | ||
7d14e4a8 | 716 | //--- Cascade switches |
717 | Bool_t isXiMinus=kTRUE; | |
718 | Bool_t isXiPlusBar=kTRUE; | |
9bc79c29 | 719 | |
1c102c9a | 720 | const AliAODTrack *ptrack=(AliAODTrack *)v0->GetDaughter(0); |
7d14e4a8 | 721 | if (!AcceptPID(pidResponse, ptrack, stack)) isXiMinus=kFALSE; |
1c102c9a | 722 | |
5bd3c9e4 | 723 | const AliAODTrack *ntrack=(AliAODTrack *)v0->GetDaughter(1); |
7d14e4a8 | 724 | if (!AcceptPID(pidResponse, ntrack, stack)) isXiPlusBar=kFALSE; |
5bd3c9e4 | 725 | |
9bc79c29 | 726 | Int_t charge=cs->ChargeXi(); |
7d14e4a8 | 727 | if (charge > 0) isXiMinus=kFALSE; |
728 | if (charge < 0) isXiPlusBar=kFALSE; | |
729 | //--- | |
730 | ||
731 | if (isXiMinus) { | |
9bc79c29 | 732 | Double_t mass=cs->MassXi(); |
733 | fXiM->Fill(mass,pt); | |
734 | Double_t m=TDatabasePDG::Instance()->GetParticle(kXiMinus)->Mass(); | |
735 | //Double_t s=0.0037; | |
736 | Double_t s=0.002 + (0.0032-0.002)/(6-1.5)*(pt-1.5); | |
737 | if (TMath::Abs(m-mass) < 3*s) { | |
738 | fXiSiP->Fill(pt); | |
7d14e4a8 | 739 | } else { |
740 | isXiMinus=kFALSE; | |
9bc79c29 | 741 | } |
742 | if (TMath::Abs(m-mass + 4.5*s) < 1.5*s) { | |
743 | fXiSiP->Fill(pt,-1); | |
744 | } | |
745 | if (TMath::Abs(m-mass - 4.5*s) < 1.5*s) { | |
746 | fXiSiP->Fill(pt,-1); | |
747 | } | |
748 | } | |
7d14e4a8 | 749 | |
750 | if (isXiPlusBar) { | |
5bd3c9e4 | 751 | Double_t mass=cs->MassXi(); |
752 | fXiBarM->Fill(mass,pt); | |
753 | Double_t m=TDatabasePDG::Instance()->GetParticle(kXiPlusBar)->Mass(); | |
754 | //Double_t s=0.0037; | |
755 | Double_t s=0.002 + (0.0032-0.002)/(6-1.5)*(pt-1.5); | |
756 | if (TMath::Abs(m-mass) < 3*s) { | |
757 | fXiBarSiP->Fill(pt); | |
7d14e4a8 | 758 | } else { |
759 | isXiPlusBar=kFALSE; | |
5bd3c9e4 | 760 | } |
761 | if (TMath::Abs(m-mass + 4.5*s) < 1.5*s) { | |
762 | fXiBarSiP->Fill(pt,-1); | |
763 | } | |
764 | if (TMath::Abs(m-mass - 4.5*s) < 1.5*s) { | |
765 | fXiBarSiP->Fill(pt,-1); | |
766 | } | |
767 | } | |
7d14e4a8 | 768 | |
769 | if (!fIsMC) continue; | |
770 | ||
771 | //++++++ MC | |
772 | if (!isXiMinus) | |
773 | if (!isXiPlusBar) continue;//check MC only for the accepted cascades | |
774 | // Here is the future association with MC | |
775 | ||
9bc79c29 | 776 | } |
777 | ||
778 | } | |
779 | ||
780 | void AliAnalysisTaskCTauPbPbaod::Terminate(Option_t *) | |
781 | { | |
782 | // The Terminate() function is the last function to be called during | |
783 | // a query. It always runs on the client, it can be used to present | |
784 | // the results graphically or save the results to file. | |
785 | ||
786 | fOutput=(TList*)GetOutputData(1); | |
787 | if (!fOutput) { | |
788 | Printf("ERROR: fOutput not available"); | |
789 | return; | |
790 | } | |
791 | ||
792 | fMult = dynamic_cast<TH1F*>(fOutput->FindObject("fMult")) ; | |
793 | if (!fMult) { | |
794 | Printf("ERROR: fMult not available"); | |
795 | return; | |
796 | } | |
797 | ||
798 | fdEdx = dynamic_cast<TH2F*>(fOutput->FindObject("fdEdx")) ; | |
799 | if (!fdEdx) { | |
800 | Printf("ERROR: fdEdx not available"); | |
801 | return; | |
802 | } | |
803 | ||
804 | fdEdxPid = dynamic_cast<TH2F*>(fOutput->FindObject("fdEdxPid")) ; | |
805 | if (!fdEdxPid) { | |
806 | Printf("ERROR: fdEdxPid not available"); | |
807 | return; | |
808 | } | |
809 | ||
810 | ||
811 | fK0sMC = dynamic_cast<TH2F*>(fOutput->FindObject("fK0sMC")) ; | |
812 | if (!fK0sMC) { | |
813 | Printf("ERROR: fK0sMC not available"); | |
814 | return; | |
815 | } | |
816 | TH1D *k0sMcPx=fK0sMC->ProjectionX(); k0sMcPx->Sumw2(); | |
817 | fK0sAs = dynamic_cast<TH2F*>(fOutput->FindObject("fK0sAs")) ; | |
818 | if (!fK0sAs) { | |
819 | Printf("ERROR: fK0sAs not available"); | |
820 | return; | |
821 | } | |
822 | TH1D *k0sAsPx=fK0sAs->ProjectionX(); | |
823 | k0sAsPx->Sumw2(); //k0sAsPx->Scale(0.69); | |
824 | ||
825 | ||
826 | ||
7d14e4a8 | 827 | // Lambda histograms |
9bc79c29 | 828 | fLambdaFromXi = dynamic_cast<TH3F*>(fOutput->FindObject("fLambdaFromXi")) ; |
829 | if (!fLambdaFromXi) { | |
830 | Printf("ERROR: fLambdaFromXi not available"); | |
831 | return; | |
832 | } | |
833 | TH1D *lambdaFromXiPx=fLambdaFromXi->ProjectionX(); lambdaFromXiPx->Sumw2(); | |
834 | ||
9bc79c29 | 835 | fLambdaMC = dynamic_cast<TH2F*>(fOutput->FindObject("fLambdaMC")) ; |
836 | if (!fLambdaMC) { | |
837 | Printf("ERROR: fLambdaMC not available"); | |
838 | return; | |
839 | } | |
840 | TH1D *lambdaMcPx=fLambdaMC->ProjectionX(); lambdaMcPx->Sumw2(); | |
841 | ||
842 | fLambdaAs = dynamic_cast<TH2F*>(fOutput->FindObject("fLambdaAs")) ; | |
843 | if (!fLambdaAs) { | |
844 | Printf("ERROR: fLambdaAs not available"); | |
845 | return; | |
846 | } | |
847 | TH1D *lambdaAsPx=fLambdaAs->ProjectionX(); | |
848 | lambdaAsPx->Sumw2(); //lambdaAsPx->Scale(0.64); | |
849 | ||
850 | fLambdaSi = dynamic_cast<TH2F*>(fOutput->FindObject("fLambdaSi")) ; | |
851 | if (!fLambdaSi) { | |
852 | Printf("ERROR: fLambdaSi not available"); | |
853 | return; | |
854 | } | |
855 | TH1D *lambdaSiPx=fLambdaSi->ProjectionX(); | |
856 | lambdaSiPx->SetName("fLambdaPt"); | |
857 | lambdaSiPx->Sumw2(); | |
858 | ||
859 | fLambdaEff = dynamic_cast<TH1D*>(fOutput->FindObject("fLambdaEff")) ; | |
860 | if (!fLambdaEff) { | |
861 | Printf("ERROR: fLambdaEff not available"); | |
862 | return; | |
863 | } | |
864 | fLambdaPt = dynamic_cast<TH1D*>(fOutput->FindObject("fLambdaPt")) ; | |
865 | if (!fLambdaPt) { | |
866 | Printf("ERROR: fLambdaPt not available"); | |
867 | return; | |
868 | } | |
869 | ||
870 | ||
7d14e4a8 | 871 | // anti-Lambda histograms |
872 | fLambdaBarFromXiBar = | |
873 | dynamic_cast<TH3F*>(fOutput->FindObject("fLambdaBarFromXiBar")) ; | |
874 | if (!fLambdaBarFromXiBar) { | |
875 | Printf("ERROR: fLambdaBarFromXiBar not available"); | |
876 | return; | |
877 | } | |
878 | TH1D *lambdaBarFromXiBarPx= | |
879 | fLambdaBarFromXiBar->ProjectionX("Bar"); lambdaBarFromXiBarPx->Sumw2(); | |
880 | ||
881 | fLambdaBarMC = dynamic_cast<TH2F*>(fOutput->FindObject("fLambdaBarMC")) ; | |
882 | if (!fLambdaBarMC) { | |
883 | Printf("ERROR: fLambdaBarMC not available"); | |
884 | return; | |
885 | } | |
886 | TH1D *lambdaBarMcPx=fLambdaBarMC->ProjectionX(); lambdaBarMcPx->Sumw2(); | |
887 | ||
888 | fLambdaBarAs = dynamic_cast<TH2F*>(fOutput->FindObject("fLambdaBarAs")) ; | |
889 | if (!fLambdaBarAs) { | |
890 | Printf("ERROR: fLambdaBarAs not available"); | |
891 | return; | |
892 | } | |
893 | TH1D *lambdaBarAsPx=fLambdaBarAs->ProjectionX(); | |
894 | lambdaBarAsPx->Sumw2(); //lambdaBarAsPx->Scale(0.64); | |
895 | ||
896 | fLambdaBarSi = dynamic_cast<TH2F*>(fOutput->FindObject("fLambdaBarSi")) ; | |
897 | if (!fLambdaBarSi) { | |
898 | Printf("ERROR: fLambdaBarSi not available"); | |
899 | return; | |
900 | } | |
901 | TH1D *lambdaBarSiPx=fLambdaBarSi->ProjectionX(); | |
902 | lambdaBarSiPx->SetName("fLambdaBarPt"); | |
903 | lambdaBarSiPx->Sumw2(); | |
904 | ||
905 | fLambdaBarEff = dynamic_cast<TH1D*>(fOutput->FindObject("fLambdaBarEff")) ; | |
906 | if (!fLambdaBarEff) { | |
907 | Printf("ERROR: fLambdaBarEff not available"); | |
908 | return; | |
909 | } | |
910 | fLambdaBarPt = dynamic_cast<TH1D*>(fOutput->FindObject("fLambdaBarPt")) ; | |
911 | if (!fLambdaBarPt) { | |
912 | Printf("ERROR: fLambdaBarPt not available"); | |
913 | return; | |
914 | } | |
915 | ||
916 | ||
9bc79c29 | 917 | if (!gROOT->IsBatch()) { |
918 | ||
919 | TCanvas *c1 = new TCanvas("c1","Mulitplicity"); | |
920 | c1->SetLogy(); | |
921 | fMult->DrawCopy() ; | |
922 | ||
923 | new TCanvas("c2","dE/dx"); | |
924 | fdEdx->DrawCopy() ; | |
925 | ||
926 | new TCanvas("c3","dE/dx with PID"); | |
927 | fdEdxPid->DrawCopy() ; | |
928 | ||
929 | if (fIsMC) { | |
930 | /* | |
931 | TH1D effK(*k0sAsPx); effK.SetTitle("Efficiency for K0s"); | |
932 | effK.Divide(k0sAsPx,k0sMcPx,1,1,"b"); | |
933 | new TCanvas("c4","Efficiency for K0s"); | |
934 | effK.DrawCopy("E") ; | |
935 | */ | |
936 | ||
7d14e4a8 | 937 | //+++ Lambdas |
9bc79c29 | 938 | fLambdaEff->Divide(lambdaAsPx,lambdaMcPx,1,1,"b"); |
939 | new TCanvas("c5","Efficiency for #Lambda"); | |
940 | fLambdaEff->DrawCopy("E") ; | |
941 | ||
942 | lambdaSiPx->Add(lambdaFromXiPx,-1); | |
943 | lambdaSiPx->Divide(fLambdaEff); | |
944 | ||
945 | new TCanvas("c6","Corrected #Lambda pt"); | |
946 | lambdaSiPx->SetTitle("Corrected #Lambda pt"); | |
947 | *fLambdaPt = *lambdaSiPx; | |
948 | fLambdaPt->SetLineColor(2); | |
949 | fLambdaPt->DrawCopy("E"); | |
950 | ||
951 | lambdaMcPx->DrawCopy("same"); | |
952 | ||
7d14e4a8 | 953 | |
954 | //+++ anti-Lambdas | |
955 | fLambdaBarEff->Divide(lambdaBarAsPx,lambdaBarMcPx,1,1,"b"); | |
956 | new TCanvas("c7","Efficiency for anti-#Lambda"); | |
957 | fLambdaBarEff->DrawCopy("E") ; | |
958 | ||
959 | lambdaBarSiPx->Add(lambdaBarFromXiBarPx,-1); | |
960 | lambdaBarSiPx->Divide(fLambdaBarEff); | |
961 | ||
962 | new TCanvas("c8","Corrected anti-#Lambda pt"); | |
963 | lambdaBarSiPx->SetTitle("Corrected anti-#Lambda pt"); | |
964 | *fLambdaBarPt = *lambdaBarSiPx; | |
965 | fLambdaBarPt->SetLineColor(2); | |
966 | fLambdaBarPt->DrawCopy("E"); | |
967 | ||
968 | lambdaBarMcPx->DrawCopy("same"); | |
9bc79c29 | 969 | } else { |
970 | new TCanvas("c6","Raw #Lambda pt"); | |
971 | lambdaSiPx->SetTitle("Raw #Lambda pt"); | |
972 | *fLambdaPt = *lambdaSiPx; | |
973 | fLambdaPt->SetLineColor(2); | |
974 | fLambdaPt->DrawCopy("E"); | |
7d14e4a8 | 975 | |
976 | ||
977 | new TCanvas("c7","Raw anti-#Lambda pt"); | |
978 | lambdaBarSiPx->SetTitle("Raw anti-#Lambda pt"); | |
979 | *fLambdaBarPt = *lambdaBarSiPx; | |
980 | fLambdaBarPt->SetLineColor(2); | |
981 | fLambdaBarPt->DrawCopy("E"); | |
9bc79c29 | 982 | } |
983 | } | |
984 | } |