]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG3/vertexingHF/AliAnalysisTaskSEDs.cxx
Major dielectron framework update; includes "alignment" to updates in
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliAnalysisTaskSEDs.cxx
CommitLineData
25c94fc8 1/**************************************************************************
2 * Copyright(c) 2007-2009, 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
c2a0a9d5 16/* $Id: $ */
25c94fc8 17
18///////////////////////////////////////////////////////////////////
19// //
20// Analysis task to produce Ds candidates mass spectra //
21// Origin: F.Prino, Torino, prino@to.infn.it //
22// //
23///////////////////////////////////////////////////////////////////
24
25#include <TClonesArray.h>
26#include <TNtuple.h>
27#include <TList.h>
28#include <TString.h>
29#include <TH1F.h>
30#include <TMath.h>
31#include <TDatabasePDG.h>
32
33#include "AliAnalysisManager.h"
34#include "AliAODHandler.h"
35#include "AliAODEvent.h"
36#include "AliAODVertex.h"
37#include "AliAODTrack.h"
38#include "AliAODMCHeader.h"
39#include "AliAODMCParticle.h"
40#include "AliAODRecoDecayHF3Prong.h"
41#include "AliAnalysisVertexingHF.h"
4c230f6d 42#include "AliRDHFCutsDstoKKpi.h"
25c94fc8 43#include "AliAnalysisTaskSE.h"
44#include "AliAnalysisTaskSEDs.h"
f93cbbb4 45#include "AliNormalizationCounter.h"
25c94fc8 46
47ClassImp(AliAnalysisTaskSEDs)
48
49
50//________________________________________________________________________
51AliAnalysisTaskSEDs::AliAnalysisTaskSEDs():
52 AliAnalysisTaskSE(),
53 fOutput(0),
54 fHistNEvents(0),
3787d195 55 fPtVsMass(0),
56 fPtVsMassAC(0),
57 fYVsPt(0),
58 fYVsPtAC(0),
59 fYVsPtSig(0),
60 fYVsPtSigAC(0),
25c94fc8 61 fReadMC(kFALSE),
62 fNPtBins(0),
4c230f6d 63 fListCuts(0),
25c94fc8 64 fMassRange(0.2),
a0bb2d06 65 fMassBinSize(0.002),
f93cbbb4 66 fCounter(0),
4c230f6d 67 fProdCuts(0),
68 fAnalysisCuts(0)
25c94fc8 69{
70 // Default constructor
71}
72
73//________________________________________________________________________
4c230f6d 74AliAnalysisTaskSEDs::AliAnalysisTaskSEDs(const char *name, AliRDHFCutsDstoKKpi* productioncuts, AliRDHFCutsDstoKKpi* analysiscuts):
25c94fc8 75 AliAnalysisTaskSE(name),
76 fOutput(0),
77 fHistNEvents(0),
3787d195 78 fPtVsMass(0),
79 fPtVsMassAC(0),
80 fYVsPt(0),
81 fYVsPtAC(0),
82 fYVsPtSig(0),
83 fYVsPtSigAC(0),
25c94fc8 84 fReadMC(kFALSE),
85 fNPtBins(0),
4c230f6d 86 fListCuts(0),
25c94fc8 87 fMassRange(0.2),
a0bb2d06 88 fMassBinSize(0.002),
f93cbbb4 89 fCounter(0),
4c230f6d 90 fProdCuts(productioncuts),
91 fAnalysisCuts(analysiscuts)
25c94fc8 92{
93 // Default constructor
94 // Output slot #1 writes into a TList container
4c230f6d 95 Int_t nptbins=fAnalysisCuts->GetNPtBins();
96 Float_t *ptlim=fAnalysisCuts->GetPtBinLimits();
97 SetPtBins(nptbins,ptlim);
25c94fc8 98
99 DefineOutput(1,TList::Class()); //My private output
100
4c230f6d 101 DefineOutput(2,TList::Class());
f93cbbb4 102
103 DefineOutput(3,AliNormalizationCounter::Class());
25c94fc8 104}
105
106//________________________________________________________________________
4c230f6d 107void AliAnalysisTaskSEDs::SetPtBins(Int_t n, Float_t* lim){
25c94fc8 108 // define pt bins for analysis
109 if(n>kMaxPtBins){
110 printf("Max. number of Pt bins = %d\n",kMaxPtBins);
111 fNPtBins=kMaxPtBins;
112 fPtLimits[0]=0.;
113 fPtLimits[1]=1.;
114 fPtLimits[2]=3.;
115 fPtLimits[3]=5.;
116 fPtLimits[4]=10.;
117 for(Int_t i=5; i<kMaxPtBins+1; i++) fPtLimits[i]=99999999.;
118 }else{
119 fNPtBins=n;
120 for(Int_t i=0; i<fNPtBins+1; i++) fPtLimits[i]=lim[i];
121 for(Int_t i=fNPtBins+1; i<kMaxPtBins+1; i++) fPtLimits[i]=99999999.;
122 }
123 if(fDebug > 1){
124 printf("Number of Pt bins = %d\n",fNPtBins);
125 for(Int_t i=0; i<fNPtBins+1; i++) printf(" Bin%d = %8.2f-%8.2f\n",i,fPtLimits[i],fPtLimits[i+1]);
126 }
127}
128//________________________________________________________________________
129AliAnalysisTaskSEDs::~AliAnalysisTaskSEDs()
130{
131 // Destructor
132 if (fOutput) {
133 delete fOutput;
134 fOutput = 0;
135 }
3787d195 136
137 if(fHistNEvents){
138 delete fHistNEvents;
139 fHistNEvents=0;
140 }
141
4c230f6d 142 if (fListCuts) {
143 delete fListCuts;
144 fListCuts = 0;
145 }
146
3787d195 147 for(Int_t i=0;i<4*fNPtBins;i++){
148
149 if(fMassHist[i]){ delete fMassHist[i]; fMassHist[i]=0;}
150 if(fMassHistCuts[i]){ delete fMassHistCuts[i]; fMassHistCuts[i]=0;}
151 if(fMassHistPhi[i]){ delete fMassHistPhi[i]; fMassHistPhi[i]=0;}
152 if(fMassHistCutsPhi[i]){ delete fMassHistCutsPhi[i]; fMassHistCutsPhi[i]=0;}
153 if(fMassHistK0st[i]){ delete fMassHistK0st[i]; fMassHistK0st[i]=0;}
154 if(fMassHistCutsK0st[i]){ delete fMassHistCutsK0st[i]; fMassHistCutsK0st[i]=0;}
155 if(fCosPHist[i]){ delete fCosPHist[i]; fCosPHist[i]=0;}
156 if(fDLenHist[i]){ delete fDLenHist[i]; fDLenHist[i]=0;}
157 if(fSumd02Hist[i]){ delete fSumd02Hist[i]; fSumd02Hist[i]=0;}
158 if(fSigVertHist[i]){ delete fSigVertHist[i]; fSigVertHist[i]=0;}
159 if(fPtMaxHist[i]){ delete fPtMaxHist[i]; fPtMaxHist[i]=0;}
160 if(fDCAHist[i]){ delete fDCAHist[i]; fDCAHist[i]=0;}
161 if(fPtProng0Hist[i]){ delete fPtProng0Hist[i]; fPtProng0Hist[i]=0;}
162 if(fPtProng1Hist[i]){ delete fPtProng1Hist[i]; fPtProng1Hist[i]=0;}
163 if(fPtProng2Hist[i]){ delete fPtProng2Hist[i]; fPtProng2Hist[i]=0;}
164 if(fDalitz[i]){ delete fDalitz[i]; fDalitz[i]=0;}
165 if(fDalitzPhi[i]){ delete fDalitzPhi[i]; fDalitzPhi[i]=0;}
166 if(fDalitzK0st[i]){ delete fDalitzK0st[i]; fDalitzK0st[i]=0;}
167
168 }
169
170 if(fPtVsMass){
171 delete fPtVsMass;
172 fPtVsMass=0;
173 }
174 if(fPtVsMassAC){
175 delete fPtVsMassAC;
176 fPtVsMassAC=0;
177 }
178 if(fYVsPt){
179 delete fYVsPt;
180 fYVsPt=0;
181 }
182 if(fYVsPtAC){
183 delete fYVsPtAC;
184 fYVsPtAC=0;
185 }
186 if(fYVsPtSig){
187 delete fYVsPtSig;
188 fYVsPtSig=0;
189 }
190 if(fYVsPtSigAC){
191 delete fYVsPtSigAC;
192 fYVsPtSigAC=0;
193 }
194
f93cbbb4 195 if(fCounter){
196 delete fCounter;
197 fCounter = 0;
198 }
199
4c230f6d 200 if (fProdCuts) {
201 delete fProdCuts;
202 fProdCuts = 0;
203 }
204 if (fAnalysisCuts) {
205 delete fAnalysisCuts;
206 fAnalysisCuts = 0;
25c94fc8 207 }
208}
209
210//________________________________________________________________________
211void AliAnalysisTaskSEDs::Init()
212{
213 // Initialization
214
215 if(fDebug > 1) printf("AnalysisTaskSEDs::Init() \n");
216
4c230f6d 217 fListCuts=new TList();
218 AliRDHFCutsDstoKKpi *production = new AliRDHFCutsDstoKKpi();
219 production=fProdCuts;
220 AliRDHFCutsDstoKKpi *analysis = new AliRDHFCutsDstoKKpi();
221 analysis=fAnalysisCuts;
222
223 fListCuts->Add(production);
224 fListCuts->Add(analysis);
225 PostData(2,fListCuts);
25c94fc8 226 return;
227}
228
229//________________________________________________________________________
230void AliAnalysisTaskSEDs::UserCreateOutputObjects()
231{
232 // Create the output container
233 //
234 if(fDebug > 1) printf("AnalysisTaskSEDs::UserCreateOutputObjects() \n");
235
236 // Several histograms are more conveniently managed in a TList
237 fOutput = new TList();
238 fOutput->SetOwner();
239 fOutput->SetName("OutputHistos");
240
241 Double_t massDs=TDatabasePDG::Instance()->GetParticle(431)->Mass();
a0bb2d06 242 Int_t nInvMassBins=(Int_t)(fMassRange/fMassBinSize+0.5);
243 if(nInvMassBins%2==1) nInvMassBins++;
3787d195 244 // Double_t minMass=massDs-1.0*nInvMassBins*fMassBinSize;
a0bb2d06 245 Double_t minMass=massDs-0.5*nInvMassBins*fMassBinSize;
3787d195 246 // Double_t maxMass=massDs+1.0*nInvMassBins*fMassBinSize;
a0bb2d06 247 Double_t maxMass=massDs+0.5*nInvMassBins*fMassBinSize;
248
25c94fc8 249 TString hisname;
86ede1c0 250 TString htype;
25c94fc8 251 Int_t index;
86ede1c0 252 for(Int_t iType=0; iType<4; iType++){
253 for(Int_t i=0;i<fNPtBins;i++){
254 if(iType==0){
255 htype="All";
256 index=GetHistoIndex(i);
257 }else if(iType==1){
258 htype="Sig";
259 index=GetSignalHistoIndex(i);
260 }else if(iType==2){
261 htype="Bkg";
262 index=GetBackgroundHistoIndex(i);
263 }else{
264 htype="ReflSig";
265 index=GetReflSignalHistoIndex(i);
266 }
267 hisname.Form("hMass%sPt%d",htype.Data(),i);
268 fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),nInvMassBins,minMass,maxMass);
269 fMassHist[index]->Sumw2();
270 hisname.Form("hMass%sPt%dCuts",htype.Data(),i);
271 fMassHistCuts[index]=new TH1F(hisname.Data(),hisname.Data(),nInvMassBins,minMass,maxMass);
272 fMassHistCuts[index]->Sumw2();
273 hisname.Form("hMass%sPt%dphi",htype.Data(),i);
274 fMassHistPhi[index]=new TH1F(hisname.Data(),hisname.Data(),100,minMass,maxMass);
275 fMassHistPhi[index]->Sumw2();
276 hisname.Form("hMass%sPt%dphiCuts",htype.Data(),i);
277 fMassHistCutsPhi[index]=new TH1F(hisname.Data(),hisname.Data(),100,minMass,maxMass);
278 fMassHistCutsPhi[index]->Sumw2();
279 hisname.Form("hMass%sPt%dk0st",htype.Data(),i);
280 fMassHistK0st[index]=new TH1F(hisname.Data(),hisname.Data(),100,minMass,maxMass);
281 fMassHistK0st[index]->Sumw2();
282 hisname.Form("hMass%sPt%dk0stCuts",htype.Data(),i);
283 fMassHistCutsK0st[index]=new TH1F(hisname.Data(),hisname.Data(),100,minMass,maxMass);
284 fMassHistCutsK0st[index]->Sumw2();
285 hisname.Form("hCosP%sPt%d",htype.Data(),i);
286 fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,1.);
287 fCosPHist[index]->Sumw2();
288 hisname.Form("hDLen%sPt%d",htype.Data(),i);
289 fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.5);
290 fDLenHist[index]->Sumw2();
291 hisname.Form("hSumd02%sPt%d",htype.Data(),i);
292 fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,1.);
293 fSumd02Hist[index]->Sumw2();
294 hisname.Form("hSigVert%sPt%d",htype.Data(),i);
295 fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
296 fSigVertHist[index]->Sumw2();
297 hisname.Form("hPtMax%sPt%d",htype.Data(),i);
298 fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,20.);
299 fPtMaxHist[index]->Sumw2();
300 hisname.Form("hPtCand%sPt%d",htype.Data(),i);
301 fPtCandHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.5,20.);
302 fPtCandHist[index]->Sumw2();
303 hisname.Form("hDCA%sPt%d",htype.Data(),i);
304 fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.,0.1);
305 fDCAHist[index]->Sumw2();
306 hisname.Form("hPtProng0%sPt%d",htype.Data(),i);
307 fPtProng0Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.0,20.);
308 fPtProng0Hist[index]->Sumw2();
309 hisname.Form("hPtProng1%sPt%d",htype.Data(),i);
310 fPtProng1Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.0,20.);
311 fPtProng1Hist[index]->Sumw2();
312 hisname.Form("hPtProng2%sPt%d",htype.Data(),i);
313 fPtProng2Hist[index]=new TH1F(hisname.Data(),hisname.Data(),100,0.0,20.);
314 fPtProng2Hist[index]->Sumw2();
315 hisname.Form("hDalitz%sPt%d",htype.Data(),i);
316 fDalitz[index]=new TH2F(hisname.Data(),hisname.Data(),100,0.,2.,100,0.,2.);
317 fDalitz[index]->Sumw2();
318 hisname.Form("hDalitz%sPt%dphi",htype.Data(),i);
319 fDalitzPhi[index]=new TH2F(hisname.Data(),hisname.Data(),100,0.,2.,100,0.,2.);
320 fDalitzPhi[index]->Sumw2();
321 hisname.Form("hDalitz%sPt%dk0st",htype.Data(),i);
322 fDalitzK0st[index]=new TH2F(hisname.Data(),hisname.Data(),100,0.,2.,100,0.,2.);
323 fDalitzK0st[index]->Sumw2();
324 }
25c94fc8 325 }
326
c2a0a9d5 327 for(Int_t i=0; i<4*fNPtBins; i++){
25c94fc8 328 fOutput->Add(fMassHist[i]);
329 fOutput->Add(fMassHistCuts[i]);
86ede1c0 330 fOutput->Add(fMassHistPhi[i]);
331 fOutput->Add(fMassHistCutsPhi[i]);
332 fOutput->Add(fMassHistK0st[i]);
333 fOutput->Add(fMassHistCutsK0st[i]);
25c94fc8 334 fOutput->Add(fCosPHist[i]);
335 fOutput->Add(fDLenHist[i]);
c2a0a9d5 336 fOutput->Add(fSumd02Hist[i]);
337 fOutput->Add(fSigVertHist[i]);
338 fOutput->Add(fPtMaxHist[i]);
339 fOutput->Add(fPtCandHist[i]);
340 fOutput->Add(fDCAHist[i]);
341 fOutput->Add(fPtProng0Hist[i]);
342 fOutput->Add(fPtProng1Hist[i]);
343 fOutput->Add(fPtProng2Hist[i]);
25c94fc8 344 fOutput->Add(fDalitz[i]);
86ede1c0 345 fOutput->Add(fDalitzPhi[i]);
346 fOutput->Add(fDalitzK0st[i]);
25c94fc8 347 }
348
349 fChanHist[0] = new TH1F("hChanAll", "KKpi and piKK candidates",4,-0.5,3.5);
350 fChanHist[1] = new TH1F("hChanSig", "KKpi and piKK candidates",4,-0.5,3.5);
351 fChanHist[2] = new TH1F("hChanBkg", "KKpi and piKK candidates",4,-0.5,3.5);
c2a0a9d5 352 fChanHist[3] = new TH1F("hChanReflSig", "KKpi and piKK candidates",4,-0.5,3.5);
25c94fc8 353 fChanHistCuts[0] = new TH1F("hChanAllCuts", "KKpi and piKK candidates",4,-0.5,3.5);
354 fChanHistCuts[1] = new TH1F("hChanSigCuts", "KKpi and piKK candidates",4,-0.5,3.5);
355 fChanHistCuts[2] = new TH1F("hChanBkgCuts", "KKpi and piKK candidates",4,-0.5,3.5);
c2a0a9d5 356 fChanHistCuts[3] = new TH1F("hChanReflSigCuts", "KKpi and piKK candidates",4,-0.5,3.5);
357 for(Int_t i=0;i<4;i++){
25c94fc8 358 fChanHist[i]->Sumw2();
359 fChanHist[i]->SetMinimum(0);
360 fChanHistCuts[i]->Sumw2();
361 fChanHistCuts[i]->SetMinimum(0);
362 fOutput->Add(fChanHist[i]);
363 fOutput->Add(fChanHistCuts[i]);
364 }
365
366 fHistNEvents = new TH1F("hNEvents", "Number of processed events",3,-1.5,1.5);
367 fHistNEvents->Sumw2();
368 fHistNEvents->SetMinimum(0);
369 fOutput->Add(fHistNEvents);
370
3787d195 371 fPtVsMass=new TH2F("hPtVsMass","PtVsMass (prod. cuts)",nInvMassBins,minMass,maxMass,40,0.,20.);
372 fPtVsMassAC=new TH2F("hPtVsMassAC","PtVsMass (analysis cuts)",nInvMassBins,minMass,maxMass,40,0.,20.);
373 fYVsPt=new TH2F("hYVsPt","YvsPt (prod. cuts)",40,0.,20.,80,-2.,2.);
374 fYVsPtAC=new TH2F("hYVsPtAC","YvsPt (analysis cuts)",40,0.,20.,80,-2.,2.);
375 fYVsPtSig=new TH2F("hYVsPtSig","YvsPt (MC, only sig., prod. cuts)",40,0.,20.,80,-2.,2.);
376 fYVsPtSigAC=new TH2F("hYVsPtSigAC","YvsPt (MC, only Sig, analysis cuts)",40,0.,20.,80,-2.,2.);
377
378 fOutput->Add(fPtVsMass);
379 fOutput->Add(fPtVsMassAC);
380 fOutput->Add(fYVsPt);
381 fOutput->Add(fYVsPtAC);
382 fOutput->Add(fYVsPtSig);
383 fOutput->Add(fYVsPtSigAC);
25c94fc8 384
f93cbbb4 385 //Counter for Normalization
386 fCounter = new AliNormalizationCounter("NormalizationCounter");
387
25c94fc8 388 return;
389}
390
391//________________________________________________________________________
392void AliAnalysisTaskSEDs::UserExec(Option_t */*option*/)
393{
394 // Ds selection for current event, fill mass histos and selecetion variable histo
395 // separate signal and backgound if fReadMC is activated
396
397 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
25c94fc8 398
399
400 TClonesArray *array3Prong = 0;
401 if(!aod && AODEvent() && IsStandardAOD()) {
402 // In case there is an AOD handler writing a standard AOD, use the AOD
403 // event in memory rather than the input (ESD) event.
404 aod = dynamic_cast<AliAODEvent*> (AODEvent());
405 // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
406 // have to taken from the AOD event hold by the AliAODExtension
407 AliAODHandler* aodHandler = (AliAODHandler*)
408 ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
409 if(aodHandler->GetExtensions()) {
410 AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
411 AliAODEvent *aodFromExt = ext->GetAOD();
412 array3Prong=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
413 }
414 } else {
415 array3Prong=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
416 }
417
418 if(!array3Prong) {
419 printf("AliAnalysisTaskSEDs::UserExec: Charm3Prong branch not found!\n");
420 return;
421 }
422
7c23877d 423
424 // fix for temporary bug in ESDfilter
425 // the AODs with null vertex pointer didn't pass the PhysSel
5806c290 426 if(!aod->GetPrimaryVertex() || TMath::Abs(aod->GetMagneticField())<0.001) return;
7c23877d 427
428
429 fHistNEvents->Fill(0); // count event
430 // Post the data already here
431 PostData(1,fOutput);
f93cbbb4 432
433 fCounter->StoreEvent(aod,fReadMC);
434
25c94fc8 435 TClonesArray *arrayMC=0;
436 AliAODMCHeader *mcHeader=0;
437
438 // AOD primary vertex
439 AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
440 // vtx1->Print();
441
442 // load MC particles
443 if(fReadMC){
444
445 arrayMC = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
446 if(!arrayMC) {
447 printf("AliAnalysisTaskSEDs::UserExec: MC particles branch not found!\n");
448 return;
449 }
450
451 // load MC header
452 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
453 if(!mcHeader) {
454 printf("AliAnalysisTaskSEDs::UserExec: MC header branch not found!\n");
455 return;
456 }
457 }
458
459 Int_t n3Prong = array3Prong->GetEntriesFast();
460 if(fDebug>1) printf("Number of Ds->KKpi: %d\n",n3Prong);
461
462
463 Int_t pdgDstoKKpi[3]={321,321,211};
f93cbbb4 464 Int_t nSelectedloose=0;
465 Int_t nSelectedtight=0;
466
25c94fc8 467 for (Int_t i3Prong = 0; i3Prong < n3Prong; i3Prong++) {
468 AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)array3Prong->UncheckedAt(i3Prong);
469
470
471 Bool_t unsetvtx=kFALSE;
472 if(!d->GetOwnPrimaryVtx()){
473 d->SetOwnPrimaryVtx(vtx1);
474 unsetvtx=kTRUE;
475 }
4c230f6d 476
477 Int_t retCodeProductionCuts=fProdCuts->IsSelected(d,AliRDHFCuts::kCandidate);
478 if(retCodeProductionCuts>0){
479 Int_t isKKpi=retCodeProductionCuts&1;
480 Int_t ispiKK=retCodeProductionCuts&2;
86ede1c0 481 Int_t isPhi=retCodeProductionCuts&4;
482 Int_t isK0star=retCodeProductionCuts&8;
25c94fc8 483 Double_t ptCand = d->Pt();
4c230f6d 484 Int_t iPtBin=TMath::BinarySearch(fNPtBins,fPtLimits,(Float_t)ptCand);
485 Int_t retCodeAnalysisCuts=fAnalysisCuts->IsSelected(d,AliRDHFCuts::kCandidate);
486 Int_t isKKpiAC=retCodeAnalysisCuts&1;
487 Int_t ispiKKAC=retCodeAnalysisCuts&2;
86ede1c0 488 Int_t isPhiAC=retCodeAnalysisCuts&4;
489 Int_t isK0starAC=retCodeAnalysisCuts&8;
4c230f6d 490
25c94fc8 491 Int_t labDs=-1;
492 if(fReadMC){
493 labDs = d->MatchToMC(431,arrayMC,3,pdgDstoKKpi);
494 }
495
496 Double_t dlen=d->DecayLength();
497 Double_t cosp=d->CosPointingAngle();
c2a0a9d5 498 Double_t pt0=d->PtProng(0);
499 Double_t pt1=d->PtProng(1);
500 Double_t pt2=d->PtProng(2);
501 Double_t sigvert=d->GetSigmaVert();
502 Double_t sumD02=d->Getd0Prong(0)*d->Getd0Prong(0)+d->Getd0Prong(1)*d->Getd0Prong(1)+d->Getd0Prong(2)*d->Getd0Prong(2);
503 Double_t dca=d->GetDCA();
504 Double_t ptmax=0;
505 for(Int_t i=0;i<3;i++){
506 if(d->PtProng(i)>ptmax)ptmax=d->PtProng(i);
507 }
3787d195 508
509 Double_t rapid=d->YDs();
510 fYVsPt->Fill(ptCand,rapid);
511 if(retCodeAnalysisCuts) fYVsPtAC->Fill(ptCand,rapid);
512
513 Bool_t isFidAcc=fAnalysisCuts->IsInFiducialAcceptance(ptCand,rapid);
c2a0a9d5 514
f93cbbb4 515 if(isFidAcc){
516 nSelectedloose++;
517 if(retCodeAnalysisCuts>0)nSelectedtight++;
518 }
519
25c94fc8 520 Int_t index=GetHistoIndex(iPtBin);
521 Int_t type=0;
522 if(isKKpi) type+=1;
523 if(ispiKK) type+=2;
4c230f6d 524 Int_t typeAC=0;
525 if(isKKpiAC) typeAC+=1;
526 if(ispiKKAC) typeAC+=2;
25c94fc8 527 fCosPHist[index]->Fill(cosp);
528 fDLenHist[index]->Fill(dlen);
c2a0a9d5 529 fSigVertHist[index]->Fill(sigvert);
530 fSumd02Hist[index]->Fill(sumD02);
531 fPtMaxHist[index]->Fill(ptmax);
532 fPtCandHist[index]->Fill(ptCand);
533 fDCAHist[index]->Fill(dca);
25c94fc8 534 fChanHist[0]->Fill(type);
c2a0a9d5 535 fPtProng0Hist[index]->Fill(pt0);
536 fPtProng1Hist[index]->Fill(pt1);
537 fPtProng2Hist[index]->Fill(pt2);
538
4c230f6d 539 if(retCodeAnalysisCuts>0) fChanHistCuts[0]->Fill(typeAC);
25c94fc8 540 if(fReadMC){
541 if(labDs>=0) {
542 index=GetSignalHistoIndex(iPtBin);
25c94fc8 543 fChanHist[1]->Fill(type);
4c230f6d 544 if(retCodeAnalysisCuts>0) fChanHistCuts[1]->Fill(typeAC);
25c94fc8 545 }else{
546 index=GetBackgroundHistoIndex(iPtBin);
25c94fc8 547 fChanHist[2]->Fill(type);
4c230f6d 548 if(retCodeAnalysisCuts>0) fChanHistCuts[2]->Fill(typeAC);
25c94fc8 549 }
550 }
551 if(isKKpi){
552 index=GetHistoIndex(iPtBin);
553 Double_t invMass=d->InvMassDsKKpi();
86ede1c0 554 Double_t massKK=d->InvMass2Prongs(0,1,321,321);
555 Double_t massKp=d->InvMass2Prongs(1,2,321,211);
3787d195 556 if(isFidAcc){
557 fMassHist[index]->Fill(invMass);
558 fPtVsMass->Fill(invMass,ptCand);
559 if(isPhi) fMassHistPhi[index]->Fill(invMass);
560 if(isK0star) fMassHistK0st[index]->Fill(invMass);
561 fDalitz[index]->Fill(massKK,massKp);
562 if(isPhi) fDalitzPhi[index]->Fill(massKK,massKp);
563 if(isK0star) fDalitzK0st[index]->Fill(massKK,massKp);
564 if(retCodeAnalysisCuts>0 && isKKpiAC){
565 fMassHistCuts[index]->Fill(invMass);
566 fPtVsMassAC->Fill(invMass,ptCand);
567 if(isPhiAC) fMassHistCutsPhi[index]->Fill(invMass);
568 if(isK0starAC) fMassHistCutsK0st[index]->Fill(invMass);
569 }
86ede1c0 570 }
25c94fc8 571 if(fReadMC){
3787d195 572 if(labDs>=0){
c2a0a9d5 573 Int_t labDau0=((AliAODTrack*)d->GetDaughter(0))->GetLabel();
574 AliAODMCParticle* p=(AliAODMCParticle*)arrayMC->UncheckedAt(labDau0);
575 Int_t pdgCode0=TMath::Abs(p->GetPdgCode());
3787d195 576 //if(labDs>=0){
c2a0a9d5 577 if(pdgCode0==321) {
578 index=GetSignalHistoIndex(iPtBin);
3787d195 579 fYVsPtSig->Fill(ptCand,rapid);
580 if(retCodeAnalysisCuts>0 && isKKpiAC) fYVsPtSigAC->Fill(ptCand,rapid);
c2a0a9d5 581 }else{
582 index=GetReflSignalHistoIndex(iPtBin);
583 }
25c94fc8 584 }else{
585 index=GetBackgroundHistoIndex(iPtBin);
25c94fc8 586 }
3787d195 587 if(isFidAcc){
588 fMassHist[index]->Fill(invMass);
589 if(isPhi) fMassHistPhi[index]->Fill(invMass);
590 if(isK0star) fMassHistK0st[index]->Fill(invMass);
591 fCosPHist[index]->Fill(cosp);
592 fDLenHist[index]->Fill(dlen);
593 fSigVertHist[index]->Fill(sigvert);
594 fSumd02Hist[index]->Fill(sumD02);
595 fPtMaxHist[index]->Fill(ptmax);
596 fPtCandHist[index]->Fill(ptCand);
597 fDCAHist[index]->Fill(dca);
598 fPtProng0Hist[index]->Fill(pt0);
599 fPtProng1Hist[index]->Fill(pt1);
600 fPtProng2Hist[index]->Fill(pt2);
601 fDalitz[index]->Fill(massKK,massKp);
602 if(isPhi) fDalitzPhi[index]->Fill(massKK,massKp);
603 if(isK0star) fDalitzK0st[index]->Fill(massKK,massKp);
604 if(retCodeAnalysisCuts>0 && isKKpiAC){
605 fMassHistCuts[index]->Fill(invMass);
606 if(isPhiAC) fMassHistCutsPhi[index]->Fill(invMass);
607 if(isK0starAC) fMassHistCutsK0st[index]->Fill(invMass);
608 }
86ede1c0 609 }
25c94fc8 610 }
611 }
612 if(ispiKK){
613 index=GetHistoIndex(iPtBin);
614 Double_t invMass=d->InvMassDspiKK();
86ede1c0 615 Double_t massKK=d->InvMass2Prongs(1,2,321,321);
616 Double_t massKp=d->InvMass2Prongs(0,1,211,321);
3787d195 617
618
619 if(isFidAcc){
620 fMassHist[index]->Fill(invMass);
621 if(isPhi) fMassHistPhi[index]->Fill(invMass);
622 if(isK0star) fMassHistK0st[index]->Fill(invMass);
623 fDalitz[index]->Fill(massKK,massKp);
624 if(isPhi) fDalitzPhi[index]->Fill(massKK,massKp);
625 if(isK0star) fDalitzK0st[index]->Fill(massKK,massKp);
626 if(retCodeAnalysisCuts>0 && ispiKKAC){
627 fMassHistCuts[index]->Fill(invMass);
628 fPtVsMassAC->Fill(invMass,ptCand);
629 if(isPhiAC) fMassHistCutsPhi[index]->Fill(invMass);
630 if(isK0starAC) fMassHistCutsK0st[index]->Fill(invMass);
631 }
86ede1c0 632 }
25c94fc8 633 if(fReadMC){
3787d195 634 if(labDs>=0) {
c2a0a9d5 635 Int_t labDau0=((AliAODTrack*)d->GetDaughter(0))->GetLabel();
636 AliAODMCParticle* p=(AliAODMCParticle*)arrayMC->UncheckedAt(labDau0);
637 Int_t pdgCode0=TMath::Abs(p->GetPdgCode());
3787d195 638 // if(labDs>=0) {
c2a0a9d5 639 if(pdgCode0==211) {
640 index=GetSignalHistoIndex(iPtBin);
3787d195 641 fYVsPtSig->Fill(ptCand,rapid);
642 if(retCodeAnalysisCuts>0 && isKKpiAC) fYVsPtSigAC->Fill(ptCand,rapid);
c2a0a9d5 643 }else{
644 index=GetReflSignalHistoIndex(iPtBin);
645 }
25c94fc8 646 }else{
647 index=GetBackgroundHistoIndex(iPtBin);
25c94fc8 648 }
3787d195 649 if(isFidAcc){
650 fMassHist[index]->Fill(invMass);
651 fPtVsMass->Fill(invMass,ptCand);
652 if(isPhi) fMassHistPhi[index]->Fill(invMass);
653 if(isK0star) fMassHistK0st[index]->Fill(invMass);
654 fCosPHist[index]->Fill(cosp);
655 fDLenHist[index]->Fill(dlen);
656 fSigVertHist[index]->Fill(sigvert);
657 fSumd02Hist[index]->Fill(sumD02);
658 fPtMaxHist[index]->Fill(ptmax);
659 fPtCandHist[index]->Fill(ptCand);
660 fDCAHist[index]->Fill(dca);
661 fPtProng0Hist[index]->Fill(pt0);
662 fPtProng1Hist[index]->Fill(pt1);
663 fPtProng2Hist[index]->Fill(pt2);
664 fDalitz[index]->Fill(massKK,massKp);
665 if(isPhi) fDalitzPhi[index]->Fill(massKK,massKp);
666 if(isK0star) fDalitzK0st[index]->Fill(massKK,massKp);
667 if(retCodeAnalysisCuts>0 && ispiKKAC){
668 fMassHistCuts[index]->Fill(invMass);
669 if(isPhiAC) fMassHistCutsPhi[index]->Fill(invMass);
670 if(isK0starAC) fMassHistCutsK0st[index]->Fill(invMass);
671 }
86ede1c0 672 }
25c94fc8 673 }
674 }
675 }
676 if(unsetvtx) d->UnsetOwnPrimaryVtx();
677 }
678
f93cbbb4 679 fCounter->StoreCandidates(aod,nSelectedloose,kTRUE);
680 fCounter->StoreCandidates(aod,nSelectedtight,kFALSE);
681
682 PostData(1,fOutput);
683 PostData(3,fCounter);
684
25c94fc8 685 return;
686}
687
688//_________________________________________________________________
689
690void AliAnalysisTaskSEDs::Terminate(Option_t */*option*/)
691{
692 // Terminate analysis
693 //
694 if(fDebug > 1) printf("AnalysisTaskSEDs: Terminate() \n");
695 fOutput = dynamic_cast<TList*> (GetOutputData(1));
696 if (!fOutput) {
697 printf("ERROR: fOutput not available\n");
698 return;
699 }
700 fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("hNEvents"));
3787d195 701 fYVsPt = dynamic_cast<TH2F*>(fOutput->FindObject("hYVsPt"));
702 fYVsPtAC = dynamic_cast<TH2F*>(fOutput->FindObject("hYVsPtAC"));
703 fYVsPtSig = dynamic_cast<TH2F*>(fOutput->FindObject("hYVsPtSig"));
704 fYVsPtSigAC = dynamic_cast<TH2F*>(fOutput->FindObject("hYVsPtSigAC"));
705 fPtVsMass = dynamic_cast<TH2F*>(fOutput->FindObject("hPtVsMass"));
706 fPtVsMassAC = dynamic_cast<TH2F*>(fOutput->FindObject("hPtVsMassAC"));
25c94fc8 707 fChanHist[0] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanAll"));
708 fChanHist[1] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanSig"));
709 fChanHist[2] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanBkg"));
c2a0a9d5 710 fChanHist[3] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanReflSig"));
25c94fc8 711 fChanHistCuts[0] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanAllCuts"));
712 fChanHistCuts[1] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanSigCuts"));
713 fChanHistCuts[2] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanBkgCuts"));
c2a0a9d5 714 fChanHistCuts[3] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanReflSigCuts"));
25c94fc8 715
f93cbbb4 716 fCounter = dynamic_cast<AliNormalizationCounter*>(GetOutputData(3));
25c94fc8 717
718 TString hisname;
86ede1c0 719 TString htype;
25c94fc8 720 Int_t index;
86ede1c0 721 for(Int_t iType=0; iType<4; iType++){
722 for(Int_t i=0;i<fNPtBins;i++){
723 if(iType==0){
724 htype="All";
725 index=GetHistoIndex(i);
726 }else if(iType==1){
727 htype="Sig";
728 index=GetSignalHistoIndex(i);
729 }else if(iType==2){
730 htype="Bkg";
731 index=GetBackgroundHistoIndex(i);
732 }else{
733 htype="ReflSig";
734 index=GetReflSignalHistoIndex(i);
735 }
736 hisname.Form("hMass%sPt%d",htype.Data(),i);
737 fMassHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
738 hisname.Form("hMass%sPt%dCuts",htype.Data(),i);
739 fMassHistCuts[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
740 hisname.Form("hMass%sPt%dphi",htype.Data(),i);
741 fMassHistPhi[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
742 hisname.Form("hMass%sPt%dphiCuts",htype.Data(),i);
743 fMassHistCutsPhi[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
744 hisname.Form("hMass%sPt%dk0st",htype.Data(),i);
745 fMassHistK0st[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
746 hisname.Form("hMass%sPt%dk0stCuts",htype.Data(),i);
747 fMassHistCutsK0st[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
748 hisname.Form("hCosP%sPt%d",htype.Data(),i);
749 fCosPHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
750 hisname.Form("hDLen%sPt%d",htype.Data(),i);
751 fDLenHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
752 hisname.Form("hSumd02%sPt%d",htype.Data(),i);
753 fSumd02Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
754 hisname.Form("hSigVert%sPt%d",htype.Data(),i);
755 fSigVertHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
756 hisname.Form("hPtMax%sPt%d",htype.Data(),i);
757 fPtMaxHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
758 hisname.Form("hPtCand%sPt%d",htype.Data(),i);
759 fPtCandHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
760 hisname.Form("hDCA%sPt%d",htype.Data(),i);
761 fDCAHist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
762 hisname.Form("hPtProng0%sPt%d",htype.Data(),i);
763 fPtProng0Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
764 hisname.Form("hPtProng1%sPt%d",htype.Data(),i);
765 fPtProng1Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
766 hisname.Form("hPtProng2%sPt%d",htype.Data(),i);
767 fPtProng2Hist[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
768 hisname.Form("hDalitz%sPt%d",htype.Data(),i);
769 fDalitz[index]=dynamic_cast<TH2F*>(fOutput->FindObject(hisname.Data()));
770 hisname.Form("hDalitz%sPt%dphi",htype.Data(),i);
771 fDalitzPhi[index]=dynamic_cast<TH2F*>(fOutput->FindObject(hisname.Data()));
772 hisname.Form("hDalitz%sPt%dk0st",htype.Data(),i);
773 fDalitzK0st[index]=dynamic_cast<TH2F*>(fOutput->FindObject(hisname.Data()));
774 }
25c94fc8 775 }
776 return;
777}
f93cbbb4 778