1 /**************************************************************************
2 * Copyright(c) 2007-2009, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 ///////////////////////////////////////////////////////////////////
20 // Analysis task to produce Ds candidates mass spectra //
21 // Origin: F.Prino, Torino, prino@to.infn.it //
23 ///////////////////////////////////////////////////////////////////
25 #include <TClonesArray.h>
31 #include <TDatabasePDG.h>
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"
42 #include "AliRDHFCutsDstoKKpi.h"
43 #include "AliAnalysisTaskSE.h"
44 #include "AliAnalysisTaskSEDs.h"
45 #include "AliNormalizationCounter.h"
47 ClassImp(AliAnalysisTaskSEDs)
50 //________________________________________________________________________
51 AliAnalysisTaskSEDs::AliAnalysisTaskSEDs():
70 // Default constructor
73 //________________________________________________________________________
74 AliAnalysisTaskSEDs::AliAnalysisTaskSEDs(const char *name, AliRDHFCutsDstoKKpi* productioncuts, AliRDHFCutsDstoKKpi* analysiscuts):
75 AliAnalysisTaskSE(name),
90 fProdCuts(productioncuts),
91 fAnalysisCuts(analysiscuts)
93 // Default constructor
94 // Output slot #1 writes into a TList container
95 Int_t nptbins=fAnalysisCuts->GetNPtBins();
96 Float_t *ptlim=fAnalysisCuts->GetPtBinLimits();
97 SetPtBins(nptbins,ptlim);
99 DefineOutput(1,TList::Class()); //My private output
101 DefineOutput(2,TList::Class());
103 DefineOutput(3,AliNormalizationCounter::Class());
106 //________________________________________________________________________
107 void AliAnalysisTaskSEDs::SetPtBins(Int_t n, Float_t* lim){
108 // define pt bins for analysis
110 printf("Max. number of Pt bins = %d\n",kMaxPtBins);
117 for(Int_t i=5; i<kMaxPtBins+1; i++) fPtLimits[i]=99999999.;
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.;
124 printf("Number of Pt bins = %d\n",fNPtBins);
125 for(Int_t i=0; i<fNPtBins; i++) printf(" Bin%d = %8.2f-%8.2f\n",i,fPtLimits[i],fPtLimits[i+1]);
128 //________________________________________________________________________
129 AliAnalysisTaskSEDs::~AliAnalysisTaskSEDs()
147 for(Int_t i=0;i<4*fNPtBins;i++){
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;}
205 delete fAnalysisCuts;
210 //________________________________________________________________________
211 void AliAnalysisTaskSEDs::Init()
215 if(fDebug > 1) printf("AnalysisTaskSEDs::Init() \n");
217 fListCuts=new TList();
218 AliRDHFCutsDstoKKpi *production = new AliRDHFCutsDstoKKpi(*fProdCuts);
219 production->SetName("ProductionCuts");
220 AliRDHFCutsDstoKKpi *analysis = new AliRDHFCutsDstoKKpi(*fAnalysisCuts);
221 analysis->SetName("AnalysisCuts");
223 fListCuts->Add(production);
224 fListCuts->Add(analysis);
225 PostData(2,fListCuts);
229 //________________________________________________________________________
230 void AliAnalysisTaskSEDs::UserCreateOutputObjects()
232 // Create the output container
234 if(fDebug > 1) printf("AnalysisTaskSEDs::UserCreateOutputObjects() \n");
236 // Several histograms are more conveniently managed in a TList
237 fOutput = new TList();
239 fOutput->SetName("OutputHistos");
241 Double_t massDs=TDatabasePDG::Instance()->GetParticle(431)->Mass();
242 Int_t nInvMassBins=(Int_t)(fMassRange/fMassBinSize+0.5);
243 if(nInvMassBins%2==1) nInvMassBins++;
244 // Double_t minMass=massDs-1.0*nInvMassBins*fMassBinSize;
245 Double_t minMass=massDs-0.5*nInvMassBins*fMassBinSize;
246 // Double_t maxMass=massDs+1.0*nInvMassBins*fMassBinSize;
247 Double_t maxMass=massDs+0.5*nInvMassBins*fMassBinSize;
252 for(Int_t iType=0; iType<4; iType++){
253 for(Int_t i=0;i<fNPtBins;i++){
256 index=GetHistoIndex(i);
259 index=GetSignalHistoIndex(i);
262 index=GetBackgroundHistoIndex(i);
265 index=GetReflSignalHistoIndex(i);
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();
327 for(Int_t i=0; i<4*fNPtBins; i++){
328 fOutput->Add(fMassHist[i]);
329 fOutput->Add(fMassHistCuts[i]);
330 fOutput->Add(fMassHistPhi[i]);
331 fOutput->Add(fMassHistCutsPhi[i]);
332 fOutput->Add(fMassHistK0st[i]);
333 fOutput->Add(fMassHistCutsK0st[i]);
334 fOutput->Add(fCosPHist[i]);
335 fOutput->Add(fDLenHist[i]);
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]);
344 fOutput->Add(fDalitz[i]);
345 fOutput->Add(fDalitzPhi[i]);
346 fOutput->Add(fDalitzK0st[i]);
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);
352 fChanHist[3] = new TH1F("hChanReflSig", "KKpi and piKK candidates",4,-0.5,3.5);
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);
356 fChanHistCuts[3] = new TH1F("hChanReflSigCuts", "KKpi and piKK candidates",4,-0.5,3.5);
357 for(Int_t i=0;i<4;i++){
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]);
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);
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.);
378 fOutput->Add(fPtVsMass);
379 fOutput->Add(fPtVsMassAC);
380 fOutput->Add(fYVsPt);
381 fOutput->Add(fYVsPtAC);
382 fOutput->Add(fYVsPtSig);
383 fOutput->Add(fYVsPtSigAC);
385 //Counter for Normalization
386 fCounter = new AliNormalizationCounter("NormalizationCounter");
391 //________________________________________________________________________
392 void AliAnalysisTaskSEDs::UserExec(Option_t */*option*/)
394 // Ds selection for current event, fill mass histos and selecetion variable histo
395 // separate signal and backgound if fReadMC is activated
397 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
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");
415 array3Prong=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
418 if(!aod || !array3Prong) {
419 printf("AliAnalysisTaskSEDs::UserExec: Charm3Prong branch not found!\n");
424 // fix for temporary bug in ESDfilter
425 // the AODs with null vertex pointer didn't pass the PhysSel
426 if(!aod->GetPrimaryVertex() || TMath::Abs(aod->GetMagneticField())<0.001) return;
429 fHistNEvents->Fill(0); // count event
430 // Post the data already here
433 fCounter->StoreEvent(aod,fReadMC);
435 TClonesArray *arrayMC=0;
436 AliAODMCHeader *mcHeader=0;
438 // AOD primary vertex
439 AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
445 arrayMC = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
447 printf("AliAnalysisTaskSEDs::UserExec: MC particles branch not found!\n");
452 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
454 printf("AliAnalysisTaskSEDs::UserExec: MC header branch not found!\n");
459 Int_t n3Prong = array3Prong->GetEntriesFast();
460 if(fDebug>1) printf("Number of Ds->KKpi: %d\n",n3Prong);
463 Int_t pdgDstoKKpi[3]={321,321,211};
464 Int_t nSelectedloose=0;
465 Int_t nSelectedtight=0;
467 for (Int_t i3Prong = 0; i3Prong < n3Prong; i3Prong++) {
468 AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)array3Prong->UncheckedAt(i3Prong);
471 Bool_t unsetvtx=kFALSE;
472 if(!d->GetOwnPrimaryVtx()){
473 d->SetOwnPrimaryVtx(vtx1);
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;
481 Int_t isPhi=retCodeProductionCuts&4;
482 Int_t isK0star=retCodeProductionCuts&8;
483 Double_t ptCand = d->Pt();
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;
488 Int_t isPhiAC=retCodeAnalysisCuts&4;
489 Int_t isK0starAC=retCodeAnalysisCuts&8;
493 labDs = d->MatchToMC(431,arrayMC,3,pdgDstoKKpi);
496 Double_t dlen=d->DecayLength();
497 Double_t cosp=d->CosPointingAngle();
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();
505 for(Int_t i=0;i<3;i++){
506 if(d->PtProng(i)>ptmax)ptmax=d->PtProng(i);
509 Double_t rapid=d->YDs();
510 fYVsPt->Fill(ptCand,rapid);
511 if(retCodeAnalysisCuts) fYVsPtAC->Fill(ptCand,rapid);
513 Bool_t isFidAcc=fAnalysisCuts->IsInFiducialAcceptance(ptCand,rapid);
517 if(retCodeAnalysisCuts>0)nSelectedtight++;
520 Int_t index=GetHistoIndex(iPtBin);
525 if(isKKpiAC) typeAC+=1;
526 if(ispiKKAC) typeAC+=2;
527 fCosPHist[index]->Fill(cosp);
528 fDLenHist[index]->Fill(dlen);
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);
534 fChanHist[0]->Fill(type);
535 fPtProng0Hist[index]->Fill(pt0);
536 fPtProng1Hist[index]->Fill(pt1);
537 fPtProng2Hist[index]->Fill(pt2);
539 if(retCodeAnalysisCuts>0) fChanHistCuts[0]->Fill(typeAC);
542 index=GetSignalHistoIndex(iPtBin);
543 fChanHist[1]->Fill(type);
544 if(retCodeAnalysisCuts>0) fChanHistCuts[1]->Fill(typeAC);
546 index=GetBackgroundHistoIndex(iPtBin);
547 fChanHist[2]->Fill(type);
548 if(retCodeAnalysisCuts>0) fChanHistCuts[2]->Fill(typeAC);
552 index=GetHistoIndex(iPtBin);
553 Double_t invMass=d->InvMassDsKKpi();
554 Double_t massKK=d->InvMass2Prongs(0,1,321,321);
555 Double_t massKp=d->InvMass2Prongs(1,2,321,211);
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);
573 Int_t labDau0=((AliAODTrack*)d->GetDaughter(0))->GetLabel();
574 AliAODMCParticle* p=(AliAODMCParticle*)arrayMC->UncheckedAt(TMath::Abs(labDau0));
575 Int_t pdgCode0=TMath::Abs(p->GetPdgCode());
578 index=GetSignalHistoIndex(iPtBin);
579 fYVsPtSig->Fill(ptCand,rapid);
580 if(retCodeAnalysisCuts>0 && isKKpiAC) fYVsPtSigAC->Fill(ptCand,rapid);
582 index=GetReflSignalHistoIndex(iPtBin);
585 index=GetBackgroundHistoIndex(iPtBin);
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);
613 index=GetHistoIndex(iPtBin);
614 Double_t invMass=d->InvMassDspiKK();
615 Double_t massKK=d->InvMass2Prongs(1,2,321,321);
616 Double_t massKp=d->InvMass2Prongs(0,1,211,321);
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);
635 Int_t labDau0=((AliAODTrack*)d->GetDaughter(0))->GetLabel();
636 AliAODMCParticle* p=(AliAODMCParticle*)arrayMC->UncheckedAt(TMath::Abs(labDau0));
637 Int_t pdgCode0=TMath::Abs(p->GetPdgCode());
640 index=GetSignalHistoIndex(iPtBin);
641 fYVsPtSig->Fill(ptCand,rapid);
642 if(retCodeAnalysisCuts>0 && isKKpiAC) fYVsPtSigAC->Fill(ptCand,rapid);
644 index=GetReflSignalHistoIndex(iPtBin);
647 index=GetBackgroundHistoIndex(iPtBin);
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);
676 if(unsetvtx) d->UnsetOwnPrimaryVtx();
679 fCounter->StoreCandidates(aod,nSelectedloose,kTRUE);
680 fCounter->StoreCandidates(aod,nSelectedtight,kFALSE);
683 PostData(3,fCounter);
688 //_________________________________________________________________
690 void AliAnalysisTaskSEDs::Terminate(Option_t */*option*/)
692 // Terminate analysis
694 if(fDebug > 1) printf("AnalysisTaskSEDs: Terminate() \n");
695 fOutput = dynamic_cast<TList*> (GetOutputData(1));
697 printf("ERROR: fOutput not available\n");
700 fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("hNEvents"));
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"));
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"));
710 fChanHist[3] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanReflSig"));
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"));
714 fChanHistCuts[3] = dynamic_cast<TH1F*>(fOutput->FindObject("hChanReflSigCuts"));
716 fCounter = dynamic_cast<AliNormalizationCounter*>(GetOutputData(3));
721 for(Int_t iType=0; iType<4; iType++){
722 for(Int_t i=0;i<fNPtBins;i++){
725 index=GetHistoIndex(i);
728 index=GetSignalHistoIndex(i);
731 index=GetBackgroundHistoIndex(i);
734 index=GetReflSignalHistoIndex(i);
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()));