]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/correlationHF/AliAnalysisTaskSEDplusCorrelations.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Adding Eulogios analysis
[u/mrichter/AliRoot.git] / PWGHF / correlationHF / AliAnalysisTaskSEDplusCorrelations.cxx
CommitLineData
35151011 1/**************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/* $Id$ */
17
18//*************************************************************************
19// Class AliAnalysisTaskSEDplusCorrelations
20// AliAnalysisTaskSE for the D+ candidates Invariant Mass Histogram and
21//comparison of heavy-flavour decay candidates
22// to MC truth (kinematics stored in the AOD)
23// Authors: Sadhana Dash (correlation)
24/////////////////////////////////////////////////////////////
25
26#include <TClonesArray.h>
27#include <TNtuple.h>
28#include <TCanvas.h>
29#include <TList.h>
30#include <TString.h>
31#include <TDatabasePDG.h>
32
33#include <AliAnalysisDataSlot.h>
34#include <AliAnalysisDataContainer.h>
35#include "AliAnalysisManager.h"
36#include "AliRDHFCutsDplustoKpipi.h"
37#include "AliAODHandler.h"
38#include "AliAODEvent.h"
39#include "AliAODPidHF.h"
40#include "AliAODVertex.h"
41#include "AliAODTrack.h"
42#include "AliAODRecoDecayHF3Prong.h"
43#include "AliAnalysisVertexingHF.h"
44#include "AliAnalysisTaskSE.h"
45#include "AliAnalysisTaskSEDplusCorrelations.h"
46#include "AliNormalizationCounter.h"
47#include "AliVParticle.h"
48#include "AliHFAssociatedTrackCuts.h"
49#include "AliReducedParticle.h"
50#include "AliHFCorrelator.h"
51
52
53using std::cout;
54using std::endl;
55
56ClassImp(AliAnalysisTaskSEDplusCorrelations)
57
58
59//________________________________________________________________________
60AliAnalysisTaskSEDplusCorrelations::AliAnalysisTaskSEDplusCorrelations():
61AliAnalysisTaskSE(),
62 fOutput(0),
63 fCorrelator(0),
64 fSelect(0),
65 fDisplacement(0),
66 fHistNEvents(0),
67 fMCSources(0),
68 fK0Origin(0),
69 fKaonOrigin(0),
70 fInvMassK0S(0),
71 fEventMix(0),
72 fPtVsMass(0),
73 fPtVsMassTC(0),
74 fYVsPt(0),
75 fYVsPtTC(0),
76 fYVsPtSig(0),
77 fYVsPtSigTC(0),
78 fUpmasslimit(1.965),
79 fLowmasslimit(1.765),
80 fNPtBins(0),
81 fBinWidth(0.002),
82 fListCuts(0),
83 fListCutsAsso(0),
84 fRDCutsAnalysis(0),
85 fCuts(0),
86 fCounter(0),
87 fReadMC(kFALSE),
88 fUseBit(kTRUE),
89 fMixing(kFALSE),
90 fSystem(kFALSE)
91{
92 // Default constructor
347d41de 93
35151011 94 for(Int_t i=0;i<kMaxPtBins+1;i++){
95 fArrayBinLimits[i]=0;
96 }
97
98}
99
100//________________________________________________________________________
101AliAnalysisTaskSEDplusCorrelations::AliAnalysisTaskSEDplusCorrelations(const char *name,AliRDHFCutsDplustoKpipi *dpluscutsana, AliHFAssociatedTrackCuts *asstrkcuts):
102 AliAnalysisTaskSE(name),
103 fOutput(0),
104 fCorrelator(0),
105 fSelect(0),
106 fDisplacement(0),
107 fHistNEvents(0),
108 fMCSources(0),
109 fK0Origin(0),
110 fKaonOrigin(0),
111 fInvMassK0S(0),
112 fEventMix(0),
113 fPtVsMass(0),
114 fPtVsMassTC(0),
115 fYVsPt(0),
116 fYVsPtTC(0),
117 fYVsPtSig(0),
118 fYVsPtSigTC(0),
119 fUpmasslimit(1.965),
120 fLowmasslimit(1.765),
121 fNPtBins(0),
122 fBinWidth(0.002),
123 fListCuts(0),
124 fListCutsAsso(0),
125 fRDCutsAnalysis(dpluscutsana),
126 fCuts(asstrkcuts),
127 fCounter(0),
128 fReadMC(kFALSE),
129 fUseBit(kTRUE),
130 fMixing(kFALSE),
131 fSystem(kFALSE)
132{
133 //
134 // Standrd constructor
135 //
136 fNPtBins=fRDCutsAnalysis->GetNPtBins();
137
138
35151011 139
140 for(Int_t i=0;i<kMaxPtBins+1;i++){
141 fArrayBinLimits[i]=0;
142 }
143
144
145 // Default constructor
146 // Output slot #1 writes into a TList container
147 DefineOutput(1,TList::Class()); //My private output
148 // Output slot #2 writes cut to private output
149 // DefineOutput(2,AliRDHFCutsDplusCorrelationstoKpipi::Class());
150 DefineOutput(2,TList::Class());
151 // Output slot #3 writes cut to private output
152 DefineOutput(3,AliNormalizationCounter::Class());
153 DefineOutput(4,AliHFAssociatedTrackCuts::Class());
154
155
156}
157
158//________________________________________________________________________
159AliAnalysisTaskSEDplusCorrelations::~AliAnalysisTaskSEDplusCorrelations()
160{
161 //
162 // Destructor
163 //
164 delete fOutput;
165 delete fListCuts;
166 delete fRDCutsAnalysis;
167 delete fCuts;
168 delete fCounter;
169 delete fCorrelator;
170
171
172}
173//_________________________________________________________________
174void AliAnalysisTaskSEDplusCorrelations::SetMassLimits(Float_t range){
175 // set invariant mass limits
176 Float_t bw=GetBinWidth();
177 fUpmasslimit = 1.865+range;
178 fLowmasslimit = 1.865-range;
179 SetBinWidth(bw);
180}
181//_________________________________________________________________
182void AliAnalysisTaskSEDplusCorrelations::SetMassLimits(Float_t lowlimit, Float_t uplimit){
183 // set invariant mass limits
184 if(uplimit>lowlimit)
185 {
186 Float_t bw=GetBinWidth();
187 fUpmasslimit = lowlimit;
188 fLowmasslimit = uplimit;
189 SetBinWidth(bw);
190 }
191}
192//________________________________________________________________
193void AliAnalysisTaskSEDplusCorrelations::SetBinWidth(Float_t w){
194 // Define width of mass bins
195 Float_t width=w;
196 Int_t nbins=(Int_t)((fUpmasslimit-fLowmasslimit)/width+0.5);
197 Int_t missingbins=4-nbins%4;
198 nbins=nbins+missingbins;
199 width=(fUpmasslimit-fLowmasslimit)/nbins;
200 if(missingbins!=0){
201 printf("AliAnalysisTaskSEDplusCorrelations::SetBinWidth: W-bin width of %f will produce histograms not rebinnable by 4. New width set to %f\n",w,width);
202 }
203 else{
204 if(fDebug>1) printf("AliAnalysisTaskSEDplusCorrelations::SetBinWidth: width set to %f\n",width);
205 }
206 fBinWidth=width;
207}
208//_________________________________________________________________
209Int_t AliAnalysisTaskSEDplusCorrelations::GetNBinsHistos(){
210 // Compute number of mass bins
211 return (Int_t)((fUpmasslimit-fLowmasslimit)/fBinWidth+0.5);
212}
213
214
215//__________________________________________
216void AliAnalysisTaskSEDplusCorrelations::Init(){
217 //
218 // Initialization
219 //
220 if(fDebug > 1) printf("AnalysisTaskSEDplusCorrelations::Init() \n");
221
222 //PostData(2,fRDCutsloose);//we should then put those cuts in a tlist if we have more than 1
223 fListCuts=new TList();
224 // fListCutsAsso=new TList();
225
226 AliRDHFCutsDplustoKpipi *analysis = new AliRDHFCutsDplustoKpipi(*fRDCutsAnalysis);
227 analysis->SetName("AnalysisCuts");
228
229 // AliHFAssociatedTrackCuts *trkcuts = new AliHFAssociatedTrackCuts(*fCuts);
230 //trkcuts->SetName("Assotrkcuts");
231
232 fListCuts->Add(analysis);
233 //fListCuts->Add(trkcuts);
234
235
236
237 PostData(2,fListCuts);
238 PostData(4,fCuts);
239
240 return;
241}
242
243//________________________________________________________________________
244void AliAnalysisTaskSEDplusCorrelations::UserCreateOutputObjects()
245{
246 // Create the output container
247 //
248 if(fDebug > 1) printf("AnalysisTaskSEDplusCorrelations::UserCreateOutputObjects() \n");
249 // correlator creation and definition
250
251 Double_t Pi = TMath::Pi();
252 fCorrelator = new AliHFCorrelator("Correlator",fCuts,fSystem); // fCuts is the hadron cut object, fSystem to switch between pp or PbPb
253 fCorrelator->SetDeltaPhiInterval((-0.5-1./32)*Pi,(1.5-1./32)*Pi); // set correct phi interval
254 //fCorrelator->SetDeltaPhiInterval(-1.57,4.71);
255 fCorrelator->SetEventMixing(fMixing); //set kFALSE/kTRUE for mixing Off/On
256 fCorrelator->SetAssociatedParticleType(fSelect); // set 1/2/3 for hadron/kaons/kzeros
257 fCorrelator->SetApplyDisplacementCut(fDisplacement); //set kFALSE/kTRUE for using the displacement cut
258 fCorrelator->SetUseMC(fReadMC);
34538691 259 fCorrelator->SetPIDmode(2);
35151011 260
261 Bool_t pooldef = fCorrelator->DefineEventPool();
262
263 if(!pooldef) AliInfo("Warning:: Event pool not defined properly");
264
265
266 // Several histograms are more conveniently managed in a TList
267 fOutput = new TList();
268 fOutput->SetOwner();
269 fOutput->SetName("OutputHistos");
270
271 TString hisname;
272 Int_t index=0;
273 Int_t nbins=GetNBinsHistos();
8d35b368 274
275
276 Int_t nbinsphi = 32;
277 Double_t philow = -0.5*Pi - Pi/32; // shift the bin by half the width so that at 0 is it the bin center
278 Double_t phiup = 1.5*Pi - Pi/32;
279
280 Int_t nbinseta = 16;
281 Double_t etalow = -1.6; // shift the bin by half the width so that at 0 is it the bin center
282 Double_t etaup = +1.6;
283
284
35151011 285
286 for(Int_t i=0;i<fNPtBins;i++){
287
288 index=GetHistoIndex(i);
289
290
291 hisname.Form("hMassVsdPhiHad%d",i);
8d35b368 292 fMassVsdPhiHistHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
35151011 293 fMassVsdPhiHistHad[index]->Sumw2();
294
295 hisname.Form("hMassVsdEtaHad%d",i);
8d35b368 296 fMassVsdEtaHistHad[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
35151011 297 fMassVsdEtaHistHad[index]->Sumw2();
298
299 hisname.Form("hMassVsdPhiKaon%d",i);
8d35b368 300 fMassVsdPhiHistKaon[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
35151011 301 fMassVsdPhiHistKaon[index]->Sumw2();
302
303 hisname.Form("hMassVsdEtaKaon%d",i);
8d35b368 304 fMassVsdEtaHistKaon[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
35151011 305 fMassVsdEtaHistKaon[index]->Sumw2();
306
307 hisname.Form("hMassVsdPhiK0%d",i);
8d35b368 308 fMassVsdPhiHistKshort[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
35151011 309 fMassVsdPhiHistKshort[index]->Sumw2();
310
311 hisname.Form("hMassVsdEtaK0%d",i);
8d35b368 312 fMassVsdEtaHistKshort[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
35151011 313 fMassVsdEtaHistKshort[index]->Sumw2();
314
315 hisname.Form("hMassVsdPhiLeadHad%d",i);
8d35b368 316 fMassVsdPhiHistLeadHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
35151011 317 fMassVsdPhiHistLeadHad[index]->Sumw2();
318
319 hisname.Form("hMassVsdEtaLeadHad%d",i);
8d35b368 320 fMassVsdEtaHistLeadHad[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
35151011 321 fMassVsdEtaHistLeadHad[index]->Sumw2();
322
323 hisname.Form("hMassPtK0S%d",i);
324 fMassHistK0S[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.3,0.8);
325 fMassHistK0S[index]->Sumw2();
326
327 hisname.Form("hLeadPt%d",i);
328 fLeadPt[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.0,50);
329 fLeadPt[index]->Sumw2();
330
331
332 hisname.Form("hMassPt%d",i);
333 fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
334 fMassHist[index]->Sumw2();
335
336
337 hisname.Form("hMassPt%dTC",i);
338 fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
339 fMassHistTC[index]->Sumw2();
340
341 hisname.Form("hMassPt%dTCPlus",i);
342 fMassHistTCPlus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
343 fMassHistTCPlus[index]->Sumw2();
344
345 hisname.Form("hMassPt%dTCMinus",i);
346 fMassHistTCMinus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
347 fMassHistTCMinus[index]->Sumw2();
348
349
350
351 index=GetSignalHistoIndex(i);
352
353 hisname.Form("hMassVsdPhiHadSig%d",i);
8d35b368 354 fMassVsdPhiHistHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
35151011 355 fMassVsdPhiHistHad[index]->Sumw2();
356
357 hisname.Form("hMassVsdEtaHadSig%d",i);
8d35b368 358 fMassVsdEtaHistHad[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
35151011 359 fMassVsdEtaHistHad[index]->Sumw2();
360
361 hisname.Form("hMassVsdPhiKaonSig%d",i);
8d35b368 362 fMassVsdPhiHistKaon[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
35151011 363 fMassVsdPhiHistKaon[index]->Sumw2();
364
365 hisname.Form("hMassVsdEtaKaonSig%d",i);
8d35b368 366 fMassVsdEtaHistKaon[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
35151011 367 fMassVsdEtaHistKaon[index]->Sumw2();
368
369 hisname.Form("hMassVsdPhiK0Sig%d",i);
8d35b368 370 fMassVsdPhiHistKshort[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
35151011 371 fMassVsdPhiHistKshort[index]->Sumw2();
372
373 hisname.Form("hMassVsdEtaK0Sig%d",i);
8d35b368 374 fMassVsdEtaHistKshort[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
35151011 375 fMassVsdEtaHistKshort[index]->Sumw2();
376
377 hisname.Form("hMassVsdPhiLeadHadSig%d",i);
8d35b368 378 fMassVsdPhiHistLeadHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
35151011 379 fMassVsdPhiHistLeadHad[index]->Sumw2();
380
381 hisname.Form("hMassVsdEtaLeadHadSig%d",i);
8d35b368 382 fMassVsdEtaHistLeadHad[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
35151011 383 fMassVsdEtaHistLeadHad[index]->Sumw2();
384
385
386 hisname.Form("hSigPtK0S%d",i);
387 fMassHistK0S[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.3,0.8);
388 fMassHistK0S[index]->Sumw2();
389
390 hisname.Form("hSigLeadPt%d",i);
391 fLeadPt[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.0,50);
392 fLeadPt[index]->Sumw2();
393
394 hisname.Form("hSigPt%d",i);
395 fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
396 fMassHist[index]->Sumw2();
397
398 hisname.Form("hSigPt%dTC",i);
399 fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
400 fMassHistTC[index]->Sumw2();
401 hisname.Form("hSigPt%dTCPlus",i);
402 fMassHistTCPlus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
403 fMassHistTCPlus[index]->Sumw2();
404 hisname.Form("hSigPt%dTCMinus",i);
405 fMassHistTCMinus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
406 fMassHistTCMinus[index]->Sumw2();
407
408
409
410 index=GetBackgroundHistoIndex(i);
411
412 hisname.Form("hMassVsdPhiBkgHad%d",i);
8d35b368 413 fMassVsdPhiHistHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
35151011 414 fMassVsdPhiHistHad[index]->Sumw2();
415
416 hisname.Form("hMassVsdEtaBkgHad%d",i);
8d35b368 417 fMassVsdEtaHistHad[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
35151011 418 fMassVsdEtaHistHad[index]->Sumw2();
419
420 hisname.Form("hMassVsdPhiBkgKaon%d",i);
8d35b368 421 fMassVsdPhiHistKaon[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
35151011 422 fMassVsdPhiHistKaon[index]->Sumw2();
423
424 hisname.Form("hMassVsdEtaBkgKaon%d",i);
8d35b368 425 fMassVsdEtaHistKaon[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
35151011 426 fMassVsdEtaHistKaon[index]->Sumw2();
427
428 hisname.Form("hMassVsdPhiBkgKshort%d",i);
8d35b368 429 fMassVsdPhiHistKshort[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
35151011 430 fMassVsdPhiHistKshort[index]->Sumw2();
431
432
433 hisname.Form("hMassVsdPhiBkgKshort%d",i);
8d35b368 434 fMassVsdEtaHistKshort[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
35151011 435 fMassVsdEtaHistKshort[index]->Sumw2();
436
437 hisname.Form("hMassVsdPhiBkgLeadHad%d",i);
8d35b368 438 fMassVsdPhiHistLeadHad[index]=new TH2F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup);
35151011 439 fMassVsdPhiHistLeadHad[index]->Sumw2();
440
441 hisname.Form("hMassVsdPhiBkgKshort%d",i);
8d35b368 442 fMassVsdEtaHistLeadHad[index]=new TH3D(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit,nbinsphi,philow,phiup,nbinseta,etalow,etaup);
35151011 443 fMassVsdEtaHistLeadHad[index]->Sumw2();
444
445 hisname.Form("hBkgPtK0S%d",i);
446 fMassHistK0S[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.3,0.8);
447 fMassHistK0S[index]->Sumw2();
448
449 hisname.Form("hLeadBkgPt%d",i);
450 fLeadPt[index]=new TH1F(hisname.Data(),hisname.Data(),500,0.0,50);
451 fLeadPt[index]->Sumw2();
452
453 hisname.Form("hBkgPt%dTC",i);
454 fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
455 fMassHistTC[index]->Sumw2();
456
457 hisname.Form("hBkgPt%dTCPlus",i);
458 fMassHistTCPlus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
459 fMassHistTCPlus[index]->Sumw2();
460
461 hisname.Form("hBkgPt%dTCMinus",i);
462 fMassHistTCMinus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
463 fMassHistTCMinus[index]->Sumw2();
464 }
465
466
467 for(Int_t i=0; i<3*fNPtBins; i++){
468 fOutput->Add(fMassVsdPhiHistHad[i]);
469 fOutput->Add(fMassVsdEtaHistHad[i]);
470 fOutput->Add(fMassVsdPhiHistKaon[i]);
471 fOutput->Add(fMassVsdEtaHistKaon[i]);
472 fOutput->Add(fMassVsdPhiHistKshort[i]);
473 fOutput->Add(fMassVsdEtaHistKshort[i]);
474 fOutput->Add(fMassVsdPhiHistLeadHad[i]);
475 fOutput->Add(fMassVsdEtaHistLeadHad[i]);
35151011 476 fOutput->Add(fMassHistK0S[i]);
477 fOutput->Add(fLeadPt[i]);
478 fOutput->Add(fMassHist[i]);
479 fOutput->Add(fMassHistTC[i]);
480 fOutput->Add(fMassHistTCPlus[i]);
481 fOutput->Add(fMassHistTCMinus[i]);
482
483
484
485 }
486 fInvMassK0S = new TH2F("K0S","K0S", 500,0.3,0.8,500,0,50);
487 fInvMassK0S->GetXaxis()->SetTitle("Invariant Mass (#pi #pi) (GeV/c^{2})");
488 fInvMassK0S->GetYaxis()->SetTitle("K0S pt (GeV/c)");
489 fOutput->Add(fInvMassK0S);
490
491 fMCSources = new TH1F("MCSources","Origin of associated particles in MC", 10, -0.5, 9.5);
492 fMCSources->GetXaxis()->SetBinLabel(1,"All ");
493 fMCSources->GetXaxis()->SetBinLabel(2," from Heavy flavour");
494 fMCSources->GetXaxis()->SetBinLabel(3," from c->D");
495 fMCSources->GetXaxis()->SetBinLabel(4," from b->D");
496 fMCSources->GetXaxis()->SetBinLabel(5," from b->B");
497 if(fReadMC) fOutput->Add(fMCSources);
498
499 fK0Origin = new TH1F("K0Origin","Origin of K0 ", 10, -0.5, 5.5);
500 fK0Origin->GetXaxis()->SetBinLabel(1,"All K0s");
501 fK0Origin->GetXaxis()->SetBinLabel(2,"K0s from Heavy flavour");
502 fK0Origin->GetXaxis()->SetBinLabel(3,"K0s from Charm");
503 fK0Origin->GetXaxis()->SetBinLabel(4,"K0s from Beauty");
504 if(fReadMC) fOutput->Add(fK0Origin);
505
506 fKaonOrigin = new TH1F("K0Origin","Origin of Kaon ", 10, -0.5, 5.5);
507 fKaonOrigin->GetXaxis()->SetBinLabel(1,"All Kaons");
508 fKaonOrigin->GetXaxis()->SetBinLabel(2,"Kaons from Heavy flavour");
509 fKaonOrigin->GetXaxis()->SetBinLabel(3,"Kaons from Charm");
510 fKaonOrigin->GetXaxis()->SetBinLabel(4,"Kaons from Beauty");
511 if(fReadMC) fOutput->Add(fKaonOrigin);
512
513
514
515 fEventMix = new TH2F("EventMixingCheck","EventMixingCheck",5,-0.5,4.5,7,-0.5,6.5);
516 if(fMixing)fOutput->Add(fEventMix);
517
518
519 fHistNEvents = new TH1F("fHistNEvents", "number of events ",11,-0.5,10.5);
520 fHistNEvents->GetXaxis()->SetBinLabel(1,"nEventsAnal");
521 fHistNEvents->GetXaxis()->SetBinLabel(2,"nEvents accepted");
522 fHistNEvents->GetXaxis()->SetBinLabel(3,"Rejected due to trigger");
523 fHistNEvents->GetXaxis()->SetBinLabel(4,"Rejected pileup events");
524 fHistNEvents->GetXaxis()->SetBinLabel(5,"Rejected due to centrality");
525 fHistNEvents->GetXaxis()->SetBinLabel(6,"Rejected due to vtxz");
526 fHistNEvents->GetXaxis()->SetBinLabel(7,"Rejected due to Physics Sel");
527 fHistNEvents->GetXaxis()->SetBinLabel(8,"Total no. of candidate");
528 fHistNEvents->GetXaxis()->SetBinLabel(9,"no. of cand wo bitmask");
529 fHistNEvents->GetXaxis()->SetBinLabel(10,"D+ after loose cuts");
530 fHistNEvents->GetXaxis()->SetBinLabel(11,"D+ after tight cuts");
531
532 fHistNEvents->GetXaxis()->SetNdivisions(1,kFALSE);
533 fHistNEvents->Sumw2();
534 fHistNEvents->SetMinimum(0);
535 fOutput->Add(fHistNEvents);
536
537 fPtVsMass=new TH2F("hPtVsMass","PtVsMass (prod. cuts)",nbins,fLowmasslimit,fUpmasslimit,200,0.,20.);
538 fPtVsMassTC=new TH2F("hPtVsMassTC","PtVsMass (analysis cuts)",nbins,fLowmasslimit,fUpmasslimit,200,0.,20.);
539 fYVsPt=new TH2F("hYVsPt","YvsPt (prod. cuts)",40,0.,20.,80,-2.,2.);
540 fYVsPtTC=new TH2F("hYVsPtTC","YvsPt (analysis cuts)",40,0.,20.,80,-2.,2.);
541 fYVsPtSig=new TH2F("hYVsPtSig","YvsPt (MC, only sig., prod. cuts)",40,0.,20.,80,-2.,2.);
542 fYVsPtSigTC=new TH2F("hYVsPtSigTC","YvsPt (MC, only Sig, analysis cuts)",40,0.,20.,80,-2.,2.);
543
544 fOutput->Add(fPtVsMass);
545 fOutput->Add(fPtVsMassTC);
546 fOutput->Add(fYVsPt);
547 fOutput->Add(fYVsPtTC);
548 fOutput->Add(fYVsPtSig);
549 fOutput->Add(fYVsPtSigTC);
550
551
552 // Counter for Normalization
553 TString normName="NormalizationCounter";
554 AliAnalysisDataContainer *cont = GetOutputSlot(3)->GetContainer();
555 if(cont)normName=(TString)cont->GetName();
556 fCounter = new AliNormalizationCounter(normName.Data());
557 fCounter->Init();
558
559
560
561 PostData(1,fOutput);
562 PostData(3,fCounter);
563 return;
564}
565
566//________________________________________________________________________
567void AliAnalysisTaskSEDplusCorrelations::UserExec(Option_t */*option*/)
568{
569 // Do the analysis
570 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
571 TClonesArray *array3Prong = 0;
572
573 if(!fMixing){
574 if(fSelect==1) cout << "TASK::Correlation with hadrons on SE "<< endl;
575 if(fSelect==2) cout << "TASK::Correlation with kaons on SE "<< endl;
576 if(fSelect==3) cout << "TASK::Correlation with kzeros on SE "<< endl;
577 }
578 if(fMixing){
579 if(fSelect==1) cout << "TASK::Correlation with hadrons on ME "<< endl;
580 if(fSelect==2) cout << "TASK::Correlation with kaons on ME "<< endl;
581 if(fSelect==3) cout << "TASK::Correlation with kzeros on ME "<< endl;
582 }
583
584
585 if(!aod && AODEvent() && IsStandardAOD()) {
586
587 aod = dynamic_cast<AliAODEvent*> (AODEvent());
588 AliAODHandler* aodHandler = (AliAODHandler*)
589 ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
590 if(aodHandler->GetExtensions()) {
591 AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
592 AliAODEvent *aodFromExt = ext->GetAOD();
593 array3Prong=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
594 }
595 } else if(aod) {
596 array3Prong=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
597 }
598
599 if(!aod || !array3Prong) {
600 printf("AliAnalysisTaskSEDplusCorrelations::UserExec: Charm3Prong branch not found!\n");
601 return;
602 }
603
604
605 // the AODs with null vertex pointer didn't pass the PhysSel
606 if(!aod->GetPrimaryVertex()||TMath::Abs(aod->GetMagneticField())<0.001) return;
607 fCounter->StoreEvent(aod,fRDCutsAnalysis,fReadMC);
608 fHistNEvents->Fill(0); // count event
609
610
611 Bool_t isEvSel=fRDCutsAnalysis->IsEventSelected(aod);
612
613
614
615 PostData(1,fOutput);
616 if(!isEvSel)return;
617 fHistNEvents->Fill(1);
618
34538691 619 // set PIDResponse for associated tracks
620 fCorrelator->SetPidAssociated();
35151011 621
622 TClonesArray *arrayMC=0;
623 AliAODMCHeader *mcHeader=0;
624 // AOD primary vertex
625 AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
626 // vtx1->Print();
627 TString primTitle = vtx1->GetTitle();
628 //if(primTitle.Contains("VertexerTracks") && vtx1->GetNContributors()>0)fHistNEvents->Fill(2);
629
630 // load MC particles
631 if(fReadMC){
632
633 arrayMC = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
634 if(!arrayMC) {
635 printf("AliAnalysisTaskSEDplusCorrelations::UserExec: MC particles branch not found!\n");
636 return;
637 }
638
639 // load MC header
640 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
641 if(!mcHeader) {
642 printf("AliAnalysisTaskSEDplusCorrelations::UserExec: MC header branch not found!\n");
643 return;
644 }
645 }
646
647 //HFCorrelators initialization (for this event)
648
649 fCorrelator->SetAODEvent(aod); // set the event to be processedfCorrelator->
650 Bool_t correlatorON = fCorrelator->Initialize(); //define the pool for mixing
651 if(!correlatorON) {
652 AliInfo("AliHFCorrelator didn't initialize the pool correctly or processed a bad event");
653 return;
654 }
655 if(fReadMC) fCorrelator->SetMCArray(arrayMC);
656
657
658
659 Int_t n3Prong = array3Prong->GetEntriesFast();
8d35b368 660 printf("Number of D+->Kpipi: %d and of tracks: %d\n",n3Prong,aod->GetNumberOfTracks());
35151011 661 Int_t nOS=0;
662 Int_t index;
663 Int_t pdgDgDplustoKpipi[3]={321,211,211};
664 Int_t nSelectedloose=0,nSelectedtight=0;
665
666
667
668 for (Int_t i3Prong = 0; i3Prong < n3Prong; i3Prong++) {
669 AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)array3Prong->UncheckedAt(i3Prong);
670 fHistNEvents->Fill(7);
671
672 if(fUseBit && !d->HasSelectionBit(AliRDHFCuts::kDplusCuts)){
673 fHistNEvents->Fill(8);
674 continue;
675 }
676
677 Int_t passTightCuts=fRDCutsAnalysis->IsSelected(d,AliRDHFCuts::kAll,aod);
678
8d35b368 679 if(!fRDCutsAnalysis->GetIsSelectedCuts()) continue;
35151011 680
681 Bool_t unsetvtx=kFALSE;
682 if(!d->GetOwnPrimaryVtx()){
683 d->SetOwnPrimaryVtx(vtx1);
684 unsetvtx=kTRUE;
685 }
686
687
688 Double_t ptCand = d->Pt();
689 Int_t iPtBin = fRDCutsAnalysis->PtBin(ptCand);
690 Bool_t recVtx=kFALSE;
691 AliAODVertex *origownvtx=0x0;
692 if(fRDCutsAnalysis->GetIsPrimaryWithoutDaughters()){
693 if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());
694 if(fRDCutsAnalysis->RecalcOwnPrimaryVtx(d,aod))recVtx=kTRUE;
695 else fRDCutsAnalysis->CleanOwnPrimaryVtx(d,aod,origownvtx);
696 }
697
698
699
700
701 Int_t labDp=-1;
702 Bool_t isDplus = kFALSE;
703
704 if(fReadMC){
705 labDp = d->MatchToMC(411,arrayMC,3,pdgDgDplustoKpipi);
706 if(labDp>=0){
707 isDplus = kTRUE;
708 }
709 }
710
711
712 Double_t invMass=d->InvMassDplus();
713 Double_t rapid=d->YDplus();
714 fYVsPt->Fill(ptCand,rapid);
715 if(passTightCuts>0.) {fYVsPtTC->Fill(ptCand,rapid);nOS++;}
716 //printf("****************: %d and of tracks: %d\n",n3Prong,nOS);
717 Bool_t isFidAcc=fRDCutsAnalysis->IsInFiducialAcceptance(ptCand,rapid);
718 if(isFidAcc){
719 fPtVsMass->Fill(invMass,ptCand);
720 if(passTightCuts>0) fPtVsMassTC->Fill(invMass,ptCand);
721 }
722
723
724 if(iPtBin>=0){
725 index=GetHistoIndex(iPtBin);
726 if(isFidAcc){
727 fHistNEvents->Fill(9);
728 nSelectedloose++;
729 fMassHist[index]->Fill(invMass);
730
731 // loop for tight cuts
34538691 732 if(passTightCuts>0){
35151011 733 fHistNEvents->Fill(10);
734 nSelectedtight++;
735 fMassHistTC[index]->Fill(invMass);
736
737 //Dplus info
738
739 Double_t phiDplus = fCorrelator->SetCorrectPhiRange(d->Phi());
740 fCorrelator->SetTriggerParticleProperties(d->Pt(),phiDplus,d->Eta());
741
742 Int_t nIDs[3] = {-9999999};
743 nIDs[0] = ((AliAODTrack*)d->GetDaughter(0))->GetID();
744 nIDs[1] = ((AliAODTrack*)d->GetDaughter(1))->GetID();
745 nIDs[2] = ((AliAODTrack*)d->GetDaughter(2))->GetID();
746
747 Double_t ptlead = 0;
748 Double_t philead = 0;
8d35b368 749 Double_t etalead = 0;
35151011 750 Double_t refpt = 0;
751
752
753
754 Bool_t execPool = fCorrelator->ProcessEventPool();
755
756 // printf("*************: %d\n",execPool);
757 if(fMixing && !execPool) {
758 AliInfo("Mixed event analysis: pool is not ready");
759 continue;
760 }
761 Int_t NofEventsinPool = 1;
762 if(fMixing) {
763 NofEventsinPool = fCorrelator->GetNofEventsInPool();
764 }
765
766
767
768 for (Int_t jMix =0; jMix < NofEventsinPool; jMix++){// loop on events in the pool; if it is SE analysis, stops at one
769
770 Bool_t analyzetracks = fCorrelator->ProcessAssociatedTracks(jMix);
771 if(!analyzetracks) {
772 AliInfo("AliHFCorrelator::Cannot process the track array");
773 continue;
774 }
775
776 //start the track loop
777
34538691 778 // Int_t NofTracks = fCorrelator->GetNofTracks();
779
780 //cout<<"*******"<<NofTracks<<endl;
781
35151011 782 for (Int_t iTrack = 0;iTrack<fCorrelator->GetNofTracks();iTrack++) {
783 Bool_t runcorrelation = fCorrelator->Correlate(iTrack);
784
785 if(!runcorrelation) continue;
786 Double_t DeltaPhi = fCorrelator->GetDeltaPhi();
787 Double_t DeltaEta = fCorrelator->GetDeltaEta();
788
789 AliReducedParticle* redpart = fCorrelator->GetAssociatedParticle();
790 Double_t phiHad=redpart->Phi();
791 Double_t ptHad=redpart->Pt();
8d35b368 792 Double_t etaHad=redpart->Eta();
35151011 793 Int_t label = redpart->GetLabel();
794 Int_t trackid = redpart->GetID();
795 phiHad = fCorrelator->SetCorrectPhiRange(phiHad);
796
797
798 // discard the dplus daughters
799 if (!fMixing){
800 if( trackid == nIDs[0] || trackid == nIDs[1] || trackid == nIDs[2]) continue;
801 }
802 // discard the negative id tracks
803 if(trackid < 0) continue;
804
805
806 FillCorrelations(d,DeltaPhi,DeltaEta,index,fSelect);
807
808 // For leading particle
809
810 if (ptHad > refpt) {
811 refpt = ptHad; ptlead = ptHad;
812 philead = d->Phi() - phiHad;
8d35b368 813 etalead = d->Eta() - etaHad;
35151011 814 if (philead < (-1)*TMath::Pi()/2) philead += 2*TMath::Pi();
815 if (philead > 3*TMath::Pi()/2) philead -= 2*TMath::Pi();
816
817 }
818
819 // montecarlo
820
821 if(fReadMC && isDplus) {
822
823 index=GetSignalHistoIndex(iPtBin);
824
347d41de 825 Bool_t* partSource = fCuts->IsMCpartFromHF(label,arrayMC); // check source of associated particle (hadron/kaon/K0)
826 FillMCCorrelations(d,DeltaPhi,DeltaEta,index,partSource,fSelect);
827 delete partSource;
35151011 828
829
830 } // readMC
831
832 }//count good tracks
833
834 // For leading particle
835 fMassVsdPhiHistLeadHad[index]->Fill(invMass,philead);
8d35b368 836 fMassVsdEtaHistLeadHad[index]->Fill(invMass,philead,etalead);
837
35151011 838 fLeadPt[index]->Fill(ptlead);
839
840 if(fReadMC && isDplus) {
841 index=GetSignalHistoIndex(iPtBin);
8d35b368 842 fMassVsdPhiHistLeadHad[index]->Fill(invMass,philead);
843 fMassVsdEtaHistLeadHad[index]->Fill(invMass,philead,etalead);
35151011 844 fLeadPt[index]->Fill(ptlead);
845
846 }
847
848 }//jmix
849
34538691 850 }// tc
35151011 851 }//fid
852
853
854
855
856
857 }
858 if(recVtx)fRDCutsAnalysis->CleanOwnPrimaryVtx(d,aod,origownvtx);
859
860 if(unsetvtx) d->UnsetOwnPrimaryVtx();
861
862 }
863
864 if(fMixing){
865 Bool_t updated = fCorrelator->PoolUpdate();
866
867 if(!updated) AliInfo("Pool was not updated");
868 }
869 fCounter->StoreCandidates(aod,nSelectedloose,kTRUE);
870 fCounter->StoreCandidates(aod,nSelectedtight,kFALSE);
871
872
873 PostData(1,fOutput);
874 PostData(2,fListCuts);
875 PostData(3,fCounter);
876 return;
877}
878
879//________________________________________________________________________
880void AliAnalysisTaskSEDplusCorrelations::Terminate(Option_t */*option*/)
881{
882 // Terminate analysis
883 //
884 if(fDebug > 1) printf("AnalysisTaskSEDplusCorrelations: Terminate() \n");
885
886 fOutput = dynamic_cast<TList*> (GetOutputData(1));
887 if (!fOutput) {
888 printf("ERROR: fOutput not available\n");
889 return;
890 }
891 fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));
892
893 TString hisname;
894 Int_t index=0;
895
896 for(Int_t i=0;i<fNPtBins;i++){
897 index=GetHistoIndex(i);
898
899 hisname.Form("hMassPt%dTC",i);
900 fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
901 }
902
903 TCanvas *c1=new TCanvas("c1","D+ invariant mass distribution",500,500);
904 c1->cd();
905 TH1F *hMassPt=(TH1F*)fOutput->FindObject("hMassPt3TC");
906 hMassPt->SetLineColor(kBlue);
907 hMassPt->SetXTitle("M[GeV/c^{2}]");
908 hMassPt->Draw();
909
910 return;
911}
912
913
914//________________________________________________________________________
915void AliAnalysisTaskSEDplusCorrelations::FillCorrelations(AliAODRecoDecayHF3Prong* d, Double_t deltaPhi, Double_t deltaEta, Int_t ind, Int_t sel) const{
916 //Filling histogams
917
918 Double_t invMass=d->InvMassDplus();
919
920
921 if(sel==1){
922 fMassVsdPhiHistHad[ind]->Fill(invMass,deltaPhi);
8d35b368 923 fMassVsdEtaHistHad[ind]->Fill(invMass,deltaPhi,deltaEta);
35151011 924 }
925 if(sel==2){
926 fMassVsdPhiHistKaon[ind]->Fill(invMass,deltaPhi);
8d35b368 927 fMassVsdEtaHistKaon[ind]->Fill(invMass,deltaPhi,deltaEta);
35151011 928 }
929 if(sel==3){
930 fMassVsdPhiHistKshort[ind]->Fill(invMass,deltaPhi);
8d35b368 931 fMassVsdEtaHistKshort[ind]->Fill(invMass,deltaPhi,deltaEta);
35151011 932 }
933
934 return;
935}
936
937
938
939//________________________________________________________________________
53454b81 940void AliAnalysisTaskSEDplusCorrelations::FillMCCorrelations(AliAODRecoDecayHF3Prong* d, Double_t deltaPhi, Double_t deltaEta, Int_t ind,Bool_t* mcSource, Int_t sel) const{
35151011 941 // Filling histos with MC information
942
943 Double_t invMass=d->InvMassDplus();
944
945
946 if(sel==1){
947 fMassVsdPhiHistHad[ind]->Fill(invMass,deltaPhi);
8d35b368 948 fMassVsdEtaHistHad[ind]->Fill(invMass,deltaPhi,deltaEta);
35151011 949
950
951 fMCSources->Fill(0);
952
53454b81 953 if(mcSource[2]&&mcSource[0]){ // is from charm ->D
35151011 954 fMCSources->Fill(1);
955 fMCSources->Fill(2);
956 }
53454b81 957 if(mcSource[2]&&mcSource[1]){ // is from beauty -> D
35151011 958 fMCSources->Fill(1);
959 fMCSources->Fill(3);
960 }
53454b81 961 if(mcSource[3]&&mcSource[1]){ // is from beauty ->B
35151011 962 fMCSources->Fill(1);
963 fMCSources->Fill(4);
964 }
965 }
966
967 if(sel==2){
968 fMassVsdPhiHistKaon[ind]->Fill(invMass,deltaPhi);
8d35b368 969 fMassVsdEtaHistKaon[ind]->Fill(invMass,deltaPhi,deltaEta);
35151011 970 fKaonOrigin->Fill(0);
53454b81 971 if(mcSource[2]&&mcSource[0]){ // is from charm ->D
35151011 972 fKaonOrigin->Fill(1);
973 fKaonOrigin->Fill(2);
974 }
53454b81 975 if(mcSource[2]&&mcSource[1]){ // is from beauty -> D
35151011 976 fKaonOrigin->Fill(1);
977 fKaonOrigin->Fill(3);
978 }
53454b81 979 if(mcSource[3]&&mcSource[1]){ // is from beauty ->B
35151011 980 fKaonOrigin->Fill(1);
981 fKaonOrigin->Fill(4);
982 }
983 }
984 if(sel==3){
985 fMassVsdPhiHistKshort[ind]->Fill(invMass,deltaPhi);
8d35b368 986 fMassVsdEtaHistKshort[ind]->Fill(invMass,deltaPhi,deltaEta);
35151011 987 fK0Origin->Fill(0);
53454b81 988 if(mcSource[2]&&mcSource[0]){ // is from charm ->D
35151011 989 fK0Origin->Fill(1);
990 fK0Origin->Fill(2);
991 }
53454b81 992 if(mcSource[2]&&mcSource[1]){ // is from beauty -> D
35151011 993 fK0Origin->Fill(1);
994 fK0Origin->Fill(3);
995 }
53454b81 996 if(mcSource[3]&&mcSource[1]){ // is from beauty ->B
35151011 997 fK0Origin->Fill(1);
998 fK0Origin->Fill(4);
999 }
1000
1001 }
1002
1003 return;
1004}
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014