]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG3/vertexingHF/AliAnalysisTaskSEDplus.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Obsolete macro removed
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliAnalysisTaskSEDplus.cxx
CommitLineData
d486095a 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
27de2dfb 16/* $Id$ */
17
993bfba1 18//*************************************************************************
19// Class AliAnalysisTaskSEDplus
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)
4afc48a2 23// Authors: Renu Bala, bala@to.infn.it
24// F. Prino, prino@to.infn.it
25// G. Ortona, ortona@to.infn.it
d486095a 26/////////////////////////////////////////////////////////////
27
28#include <TClonesArray.h>
29#include <TNtuple.h>
13808a30 30#include <TCanvas.h>
d486095a 31#include <TList.h>
1f4e9722 32#include <TString.h>
d486095a 33#include <TH1F.h>
4afc48a2 34#include <TDatabasePDG.h>
b557eb43 35
36#include "AliAnalysisManager.h"
4c230f6d 37#include "AliRDHFCutsDplustoKpipi.h"
b557eb43 38#include "AliAODHandler.h"
d486095a 39#include "AliAODEvent.h"
40#include "AliAODVertex.h"
41#include "AliAODTrack.h"
42#include "AliAODMCHeader.h"
43#include "AliAODMCParticle.h"
44#include "AliAODRecoDecayHF3Prong.h"
45#include "AliAnalysisVertexingHF.h"
46#include "AliAnalysisTaskSE.h"
47#include "AliAnalysisTaskSEDplus.h"
a96083b9 48#include "AliNormalizationCounter.h"
d486095a 49
50ClassImp(AliAnalysisTaskSEDplus)
51
52
53//________________________________________________________________________
54AliAnalysisTaskSEDplus::AliAnalysisTaskSEDplus():
55AliAnalysisTaskSE(),
595cc7e2 56 fOutput(0),
57 fHistNEvents(0),
58 fPtVsMass(0),
59 fPtVsMassTC(0),
60 fYVsPt(0),
61 fYVsPtTC(0),
62 fYVsPtSig(0),
63 fYVsPtSigTC(0),
64 fNtupleDplus(0),
65 fUpmasslimit(1.965),
66 fLowmasslimit(1.765),
67 fNPtBins(0),
68 fBinWidth(0.002),
69 fListCuts(0),
70 fRDCutsProduction(0),
71 fRDCutsAnalysis(0),
a96083b9 72 fCounter(0),
595cc7e2 73 fFillNtuple(kFALSE),
74 fReadMC(kFALSE),
a96083b9 75 fUseStrangeness(kFALSE),
5fc4893f 76 fUseBit(kTRUE),
595cc7e2 77 fDoLS(kFALSE)
d486095a 78{
82bb8d4b 79 // Default constructor
d486095a 80}
81
82//________________________________________________________________________
4c230f6d 83AliAnalysisTaskSEDplus::AliAnalysisTaskSEDplus(const char *name,AliRDHFCutsDplustoKpipi *dpluscutsana,AliRDHFCutsDplustoKpipi *dpluscutsprod,Bool_t fillNtuple):
d486095a 84AliAnalysisTaskSE(name),
4afc48a2 85fOutput(0),
86fHistNEvents(0),
595cc7e2 87fPtVsMass(0),
88fPtVsMassTC(0),
89fYVsPt(0),
90fYVsPtTC(0),
91fYVsPtSig(0),
92fYVsPtSigTC(0),
d486095a 93fNtupleDplus(0),
82bb8d4b 94fUpmasslimit(1.965),
95fLowmasslimit(1.765),
96fNPtBins(0),
993bfba1 97fBinWidth(0.002),
4c230f6d 98fListCuts(0),
99fRDCutsProduction(dpluscutsprod),
100fRDCutsAnalysis(dpluscutsana),
a96083b9 101fCounter(0),
1f4e9722 102fFillNtuple(fillNtuple),
4afc48a2 103fReadMC(kFALSE),
a96083b9 104fUseStrangeness(kFALSE),
5fc4893f 105fUseBit(kTRUE),
4c230f6d 106fDoLS(kFALSE)
d486095a 107{
4c230f6d 108 //
109 // Standrd constructor
110 //
111 //Double_t ptlim[5]={0.,2.,3.,5,9999999.};
112 //SetPtBinLimit(5, ptlim);
113 SetPtBinLimit(fRDCutsAnalysis->GetNPtBins()+1,fRDCutsAnalysis->GetPtBinLimits());
d486095a 114 // Default constructor
82bb8d4b 115 // Output slot #1 writes into a TList container
d486095a 116 DefineOutput(1,TList::Class()); //My private output
4c230f6d 117 // Output slot #2 writes cut to private output
118 // DefineOutput(2,AliRDHFCutsDplustoKpipi::Class());
119 DefineOutput(2,TList::Class());
a96083b9 120// Output slot #3 writes cut to private output
121 DefineOutput(3,AliNormalizationCounter::Class());
122
1f4e9722 123 if(fFillNtuple){
a96083b9 124 // Output slot #4 writes into a TNtuple container
125 DefineOutput(4,TNtuple::Class()); //My private output
1f4e9722 126 }
d486095a 127}
128
129//________________________________________________________________________
130AliAnalysisTaskSEDplus::~AliAnalysisTaskSEDplus()
131{
4c230f6d 132 //
d486095a 133 // Destructor
4c230f6d 134 //
d486095a 135 if (fOutput) {
136 delete fOutput;
137 fOutput = 0;
138 }
137c7bcd 139 if(fHistNEvents){
140 delete fHistNEvents;
141 fHistNEvents=0;
142 }
37fc80e6 143 for(Int_t i=0;i<3;i++){
144 if(fHistCentrality[i]){delete fHistCentrality[i]; fHistCentrality[i]=0;}
145 }
137c7bcd 146 for(Int_t i=0;i<3*fNPtBins;i++){
147 if(fMassHist[i]){ delete fMassHist[i]; fMassHist[i]=0;}
148 if(fCosPHist[i]){ delete fCosPHist[i]; fCosPHist[i]=0;}
149 if(fDLenHist[i]){ delete fDLenHist[i]; fDLenHist[i]=0;}
150 if(fSumd02Hist[i]){ delete fSumd02Hist[i]; fSumd02Hist[i]=0;}
151 if(fSigVertHist[i]){ delete fSigVertHist[i]; fSigVertHist[i]=0;}
152 if(fPtMaxHist[i]){ delete fPtMaxHist[i]; fPtMaxHist[i]=0;}
153 if(fDCAHist[i]){ delete fDCAHist[i]; fDCAHist[i]=0;}
154 if(fMassHistTC[i]){ delete fMassHistTC[i]; fMassHistTC[i]=0;}
bb69f127 155 if(fMassHistTCPlus[i]){ delete fMassHistTCPlus[i]; fMassHistTCPlus[i]=0;}
156 if(fMassHistTCMinus[i]){ delete fMassHistTCMinus[i]; fMassHistTCMinus[i]=0;}
4c230f6d 157
137c7bcd 158 if(fMassHistLS[i]){ delete fMassHistLS[i]; fMassHistLS[i]=0;}
159 if(fCosPHistLS[i]){ delete fCosPHistLS[i]; fCosPHistLS[i]=0;}
160 if(fDLenHistLS[i]){ delete fDLenHistLS[i]; fDLenHistLS[i]=0;}
161 if(fSumd02HistLS[i]){ delete fSumd02HistLS[i]; fSumd02HistLS[i]=0;}
162 if(fSigVertHistLS[i]){ delete fSigVertHistLS[i]; fSigVertHistLS[i]=0;}
163 if(fPtMaxHistLS[i]){ delete fPtMaxHistLS[i]; fPtMaxHistLS[i]=0;}
164 if(fDCAHistLS[i]){ delete fDCAHistLS[i]; fDCAHistLS[i]=0;}
165 if(fMassHistLSTC[i]){ delete fMassHistLSTC[i]; fMassHistLSTC[i]=0;}
166 }
167 if(fPtVsMass){
168 delete fPtVsMass;
169 fPtVsMass=0;
170 }
171 if(fPtVsMassTC){
172 delete fPtVsMassTC;
173 fPtVsMassTC=0;
174 }
175 if(fYVsPt){
176 delete fYVsPt;
177 fYVsPt=0;
178 }
179 if(fYVsPtTC){
180 delete fYVsPtTC;
181 fYVsPtTC=0;
182 }
183 if(fYVsPtSig){
184 delete fYVsPtSig;
185 fYVsPtSig=0;
186 }
187 if(fYVsPtSigTC){
188 delete fYVsPtSigTC;
189 fYVsPtSigTC=0;
190 }
191 if(fNtupleDplus){
192 delete fNtupleDplus;
193 fNtupleDplus=0;
194 }
4c230f6d 195 if (fListCuts) {
196 delete fListCuts;
197 fListCuts = 0;
d486095a 198 }
4c230f6d 199 if(fRDCutsProduction){
200 delete fRDCutsProduction;
201 fRDCutsProduction = 0;
202 }
137c7bcd 203 if(fRDCutsAnalysis){
4c230f6d 204 delete fRDCutsAnalysis;
205 fRDCutsAnalysis = 0;
206 }
207
a96083b9 208 if(fCounter){
209 delete fCounter;
210 fCounter = 0;
211 }
137c7bcd 212
213
d486095a 214}
82bb8d4b 215//_________________________________________________________________
216void AliAnalysisTaskSEDplus::SetMassLimits(Float_t range){
4c230f6d 217 // set invariant mass limits
993bfba1 218 Float_t bw=GetBinWidth();
82bb8d4b 219 fUpmasslimit = 1.865+range;
220 fLowmasslimit = 1.865-range;
993bfba1 221 SetBinWidth(bw);
82bb8d4b 222}
223//_________________________________________________________________
224void AliAnalysisTaskSEDplus::SetMassLimits(Float_t lowlimit, Float_t uplimit){
4c230f6d 225 // set invariant mass limits
82bb8d4b 226 if(uplimit>lowlimit)
227 {
993bfba1 228 Float_t bw=GetBinWidth();
82bb8d4b 229 fUpmasslimit = lowlimit;
230 fLowmasslimit = uplimit;
993bfba1 231 SetBinWidth(bw);
82bb8d4b 232 }
233}
82bb8d4b 234//________________________________________________________________________
4c230f6d 235void AliAnalysisTaskSEDplus::SetPtBinLimit(Int_t n, Float_t* lim){
82bb8d4b 236 // define pt bins for analysis
237 if(n>kMaxPtBins){
238 printf("Max. number of Pt bins = %d\n",kMaxPtBins);
239 fNPtBins=kMaxPtBins;
240 fArrayBinLimits[0]=0.;
241 fArrayBinLimits[1]=2.;
242 fArrayBinLimits[2]=3.;
243 fArrayBinLimits[3]=5.;
244 for(Int_t i=4; i<kMaxPtBins+1; i++) fArrayBinLimits[i]=99999999.;
245 }else{
246 fNPtBins=n-1;
247 fArrayBinLimits[0]=lim[0];
248 for(Int_t i=1; i<fNPtBins+1; i++)
249 if(lim[i]>fArrayBinLimits[i-1]){
250 fArrayBinLimits[i]=lim[i];
251 }
252 else {
253 fArrayBinLimits[i]=fArrayBinLimits[i-1];
254 }
255 for(Int_t i=fNPtBins; i<kMaxPtBins+1; i++) fArrayBinLimits[i]=99999999.;
256 }
257 if(fDebug > 1){
258 printf("Number of Pt bins = %d\n",fNPtBins);
e11ae259 259 for(Int_t i=0; i<fNPtBins; i++) printf(" Bin%d = %8.2f-%8.2f\n",i,fArrayBinLimits[i],fArrayBinLimits[i+1]);
82bb8d4b 260 }
261}
993bfba1 262//________________________________________________________________
263void AliAnalysisTaskSEDplus::SetBinWidth(Float_t w){
264 Float_t width=w;
265 Int_t nbins=(Int_t)((fUpmasslimit-fLowmasslimit)/width+0.5);
266 Int_t missingbins=4-nbins%4;
267 nbins=nbins+missingbins;
268 width=(fUpmasslimit-fLowmasslimit)/nbins;
269 if(missingbins!=0){
270 printf("AliAnalysisTaskSEDplus::SetBinWidth: W-bin width of %f will produce histograms not rebinnable by 4. New width set to %f\n",w,width);
271 }
272 else{
273 if(fDebug>1) printf("AliAnalysisTaskSEDplus::SetBinWidth: width set to %f\n",width);
274 }
275 fBinWidth=width;
276}
82bb8d4b 277//_________________________________________________________________
278Double_t AliAnalysisTaskSEDplus::GetPtBinLimit(Int_t ibin){
4c230f6d 279 // get pt bin limit
82bb8d4b 280 if(ibin>fNPtBins)return -1;
281 return fArrayBinLimits[ibin];
282}
993bfba1 283//_________________________________________________________________
284Int_t AliAnalysisTaskSEDplus::GetNBinsHistos(){
285 return (Int_t)((fUpmasslimit-fLowmasslimit)/fBinWidth+0.5);
286}
82bb8d4b 287//_________________________________________________________________
288void AliAnalysisTaskSEDplus::LSAnalysis(TClonesArray *arrayOppositeSign,TClonesArray *arrayLikeSign,AliAODEvent *aod,AliAODVertex *vtx1, Int_t nDplusOS){
4c230f6d 289 //
290 //
291 // Fill the Like Sign histograms
292 //
5fc4893f 293 printf("started LS\n");
294 Bool_t fPlotVariables=kTRUE;
295 //histograms for like sign
296 Int_t nbins=GetNBinsHistos();;
297 TH1F *histLSPlus = new TH1F("LSPlus","LSPlus",nbins,fLowmasslimit,fUpmasslimit);
298 TH1F *histLSMinus = new TH1F("LSMinus","LSMinus",nbins,fLowmasslimit,fUpmasslimit);
299
82bb8d4b 300 Int_t nPosTrks=0,nNegTrks=0;
82bb8d4b 301 Int_t nOStriplets = arrayOppositeSign->GetEntriesFast();
82bb8d4b 302 Int_t nDplusLS=0;
5fc4893f 303 Int_t nDminusLS=0;
82bb8d4b 304 Int_t nLikeSign = arrayLikeSign->GetEntriesFast();
5fc4893f 305 Int_t index=0;
306
307 // loop over like sign candidates
82bb8d4b 308 for(Int_t iLikeSign = 0; iLikeSign < nLikeSign; iLikeSign++) {
309 AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)arrayLikeSign->UncheckedAt(iLikeSign);
5fc4893f 310 if(fUseBit && !d->HasSelectionBit(AliRDHFCuts::kDplusCuts))continue;
82bb8d4b 311 Bool_t unsetvtx=kFALSE;
312 if(!d->GetOwnPrimaryVtx()) {
313 d->SetOwnPrimaryVtx(vtx1); // needed to compute all variables
314 unsetvtx=kTRUE;
315 }
5fc4893f 316
317 // if(fRDCutsProduction->IsSelected(d,AliRDHFCuts::kCandidate,aod)){
318 fRDCutsAnalysis->IsSelected(d,AliRDHFCuts::kCandidate,aod);
319 Int_t passTightCuts=fRDCutsAnalysis->GetIsSelectedCuts();
320 if(passTightCuts>0){
321
82bb8d4b 322 Double_t ptCand = d->Pt();
5fc4893f 323 Int_t iPtBin = fRDCutsAnalysis->PtBin(ptCand);
324 if(iPtBin<0)continue;
82bb8d4b 325
82bb8d4b 326 Int_t sign= d->GetCharge();
82bb8d4b 327 Float_t invMass = d->InvMassDplus();
5fc4893f 328
82bb8d4b 329 index=GetLSHistoIndex(iPtBin);
82bb8d4b 330 fMassHistLS[index+1]->Fill(invMass);
5fc4893f 331 if(sign>0){
332 histLSPlus->Fill(invMass);
333 nDplusLS++;
334 }else{
335 histLSMinus->Fill(invMass);
336 nDminusLS++;
337 }
338 if(fPlotVariables){
339 Double_t dlen=d->DecayLength();
340 Double_t cosp=d->CosPointingAngle();
341 Double_t sumD02=d->Getd0Prong(0)*d->Getd0Prong(0)+d->Getd0Prong(1)*d->Getd0Prong(1)+d->Getd0Prong(2)*d->Getd0Prong(2);
342 Double_t dca=d->GetDCA();
343 Double_t sigvert=d->GetSigmaVert();
344 Double_t ptmax=0;
345 for(Int_t i=0;i<3;i++){
346 if(d->PtProng(i)>ptmax)ptmax=d->PtProng(i);
347 }
348 fCosPHistLS[iPtBin]->Fill(cosp);
349 fDLenHistLS[iPtBin]->Fill(dlen);
350 fSumd02HistLS[iPtBin]->Fill(sumD02);
351 fSigVertHistLS[iPtBin]->Fill(sigvert);
352 fPtMaxHistLS[iPtBin]->Fill(ptmax);
353 fDCAHistLS[iPtBin]->Fill(dca);
82bb8d4b 354 }
355 }
356 if(unsetvtx) d->UnsetOwnPrimaryVtx();
357 }
5fc4893f 358 //wei:normalized to the number of combinations (production)
359 //wei2:normalized to the number of LS/OS (production)
360 //wei3:normalized to the number of LS/OS (analysis)
361 //wei4:normalized to the number of combinations (analysis)
362 Float_t wei2=0;
363 if(nLikeSign!=0)wei2 = (Float_t)nOStriplets/(Float_t)nLikeSign;
364 Float_t wei3=0;
365 if(nDplusLS!=0)wei3 = (Float_t)nDplusOS/(Float_t)nDplusLS;
366 Float_t weiplus=1.,weiminus=1.;
367 Float_t wei4plus=1.,wei4minus=1.;
368 //wei* should be automatically protected, since to get a triplet there must be at least 3 good tracks in the event
369 if(nPosTrks>2)weiplus=3.*(Float_t)nNegTrks/((Float_t)nPosTrks-2.);
370 if(nDplusLS>2)wei4plus=3.*(Float_t)nDminusLS/((Float_t)nDplusLS-2.);
371 if(nNegTrks>2)weiminus=3.*(Float_t)nPosTrks/((Float_t)nNegTrks-2.);
372 if(nDminusLS>2)wei4minus=3.*(Float_t)nDplusLS/((Float_t)nDminusLS-2.);
373
374 fMassHistLS[index]->Add(histLSPlus,weiplus);
375 fMassHistLS[index]->Add(histLSMinus,weiminus);
376 fMassHistLS[index+2]->Add(histLSPlus,wei2);
377 fMassHistLS[index+2]->Add(histLSMinus,wei2);
378 fMassHistLS[index+3]->Add(histLSPlus,wei3);
379 fMassHistLS[index+3]->Add(histLSMinus,wei3);
380 fMassHistLS[index+4]->Add(histLSPlus,wei4plus);
381 fMassHistLS[index+4]->Add(histLSMinus,wei4minus);
382
383 delete histLSPlus;histLSPlus=0;
384 delete histLSMinus;histLSMinus=0;
82bb8d4b 385
5fc4893f 386 if(fDebug>1) printf("LS analysis terminated\n");
82bb8d4b 387}
d486095a 388
4c230f6d 389//__________________________________________
390void AliAnalysisTaskSEDplus::Init(){
391 //
392 // Initialization
393 //
d486095a 394 if(fDebug > 1) printf("AnalysisTaskSEDplus::Init() \n");
4c230f6d 395
396 //PostData(2,fRDCutsloose);//we should then put those cuts in a tlist if we have more than 1
397 fListCuts=new TList();
90993728 398 AliRDHFCutsDplustoKpipi *production = new AliRDHFCutsDplustoKpipi(*fRDCutsProduction);
399 production->SetName("ProductionCuts");
400 AliRDHFCutsDplustoKpipi *analysis = new AliRDHFCutsDplustoKpipi(*fRDCutsAnalysis);
401 analysis->SetName("AnalysisCuts");
4c230f6d 402
403 fListCuts->Add(production);
404 fListCuts->Add(analysis);
405 PostData(2,fListCuts);
406
d486095a 407 return;
408}
409
410//________________________________________________________________________
411void AliAnalysisTaskSEDplus::UserCreateOutputObjects()
412{
413 // Create the output container
414 //
415 if(fDebug > 1) printf("AnalysisTaskSEDplus::UserCreateOutputObjects() \n");
416
417 // Several histograms are more conveniently managed in a TList
418 fOutput = new TList();
419 fOutput->SetOwner();
1f4e9722 420 fOutput->SetName("OutputHistos");
421
1f4e9722 422 TString hisname;
82bb8d4b 423 Int_t index=0;
424 Int_t indexLS=0;
993bfba1 425 Int_t nbins=GetNBinsHistos();
37fc80e6 426 fHistCentrality[0]=new TH1F("centrality","centrality",100,0.5,30000.5);
427 fHistCentrality[1]=new TH1F("centrality(selectedCent)","centrality(selectedCent)",100,0.5,30000.5);
428 fHistCentrality[2]=new TH1F("centrality(OutofCent)","centrality(OutofCent)",100,0.5,30000.5);
429 for(Int_t i=0;i<3;i++){
430 fHistCentrality[i]->Sumw2();
431 fOutput->Add(fHistCentrality[i]);
432 }
82bb8d4b 433 for(Int_t i=0;i<fNPtBins;i++){
1f4e9722 434
82bb8d4b 435 index=GetHistoIndex(i);
436 indexLS=GetLSHistoIndex(i);
ba9ae5b2 437
82bb8d4b 438 hisname.Form("hMassPt%d",i);
993bfba1 439 fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 440 fMassHist[index]->Sumw2();
ba9ae5b2 441 hisname.Form("hCosPAllPt%d",i);
993bfba1 442 fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.5,1.);
ba9ae5b2 443 fCosPHist[index]->Sumw2();
444 hisname.Form("hDLenAllPt%d",i);
993bfba1 445 fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.5);
ba9ae5b2 446 fDLenHist[index]->Sumw2();
447 hisname.Form("hSumd02AllPt%d",i);
993bfba1 448 fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,1.);
ba9ae5b2 449 fSumd02Hist[index]->Sumw2();
450 hisname.Form("hSigVertAllPt%d",i);
993bfba1 451 fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.1);
ba9ae5b2 452 fSigVertHist[index]->Sumw2();
453 hisname.Form("hPtMaxAllPt%d",i);
993bfba1 454 fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.5,5.);
ba9ae5b2 455 fPtMaxHist[index]->Sumw2();
456
457 hisname.Form("hDCAAllPt%d",i);
993bfba1 458 fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.1);
ba9ae5b2 459 fDCAHist[index]->Sumw2();
460
461
462
1f4e9722 463 hisname.Form("hMassPt%dTC",i);
993bfba1 464 fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 465 fMassHistTC[index]->Sumw2();
bb69f127 466 hisname.Form("hMassPt%dTCPlus",i);
467 fMassHistTCPlus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
468 fMassHistTCPlus[index]->Sumw2();
469 hisname.Form("hMassPt%dTCMinus",i);
470 fMassHistTCMinus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
471 fMassHistTCMinus[index]->Sumw2();
ba9ae5b2 472
473
474
475
476 hisname.Form("hCosPAllPt%dLS",i);
993bfba1 477 fCosPHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.5,1.);
ba9ae5b2 478 fCosPHistLS[index]->Sumw2();
479 hisname.Form("hDLenAllPt%dLS",i);
993bfba1 480 fDLenHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.5);
ba9ae5b2 481 fDLenHistLS[index]->Sumw2();
482 hisname.Form("hSumd02AllPt%dLS",i);
993bfba1 483 fSumd02HistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,1.);
ba9ae5b2 484 fSumd02HistLS[index]->Sumw2();
485 hisname.Form("hSigVertAllPt%dLS",i);
993bfba1 486 fSigVertHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.1);
ba9ae5b2 487 fSigVertHistLS[index]->Sumw2();
488 hisname.Form("hPtMaxAllPt%dLS",i);
993bfba1 489 fPtMaxHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.5,5.);
ba9ae5b2 490 fPtMaxHistLS[index]->Sumw2();
491
492 hisname.Form("hDCAAllPt%dLS",i);
993bfba1 493 fDCAHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.1);
ba9ae5b2 494 fDCAHistLS[index]->Sumw2();
495
82bb8d4b 496 hisname.Form("hLSPt%dLC",i);
993bfba1 497 fMassHistLS[indexLS] = new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 498 fMassHistLS[indexLS]->Sumw2();
499
ba9ae5b2 500 hisname.Form("hLSPt%dTC",i);
993bfba1 501 fMassHistLSTC[indexLS] = new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
ba9ae5b2 502 fMassHistLSTC[indexLS]->Sumw2();
503
504
505
82bb8d4b 506 index=GetSignalHistoIndex(i);
507 indexLS++;
508 hisname.Form("hSigPt%d",i);
993bfba1 509 fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 510 fMassHist[index]->Sumw2();
ba9ae5b2 511 hisname.Form("hCosPSigPt%d",i);
993bfba1 512 fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.5,1.);
ba9ae5b2 513 fCosPHist[index]->Sumw2();
514 hisname.Form("hDLenSigPt%d",i);
993bfba1 515 fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.5);
ba9ae5b2 516 fDLenHist[index]->Sumw2();
517 hisname.Form("hSumd02SigPt%d",i);
993bfba1 518 fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,1.);
ba9ae5b2 519 fSumd02Hist[index]->Sumw2();
520 hisname.Form("hSigVertSigPt%d",i);
993bfba1 521 fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.1);
ba9ae5b2 522 fSigVertHist[index]->Sumw2();
523 hisname.Form("hPtMaxSigPt%d",i);
993bfba1 524 fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.5,5.);
ba9ae5b2 525 fPtMaxHist[index]->Sumw2();
526
527 hisname.Form("hDCASigPt%d",i);
993bfba1 528 fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.1);
ba9ae5b2 529 fDCAHist[index]->Sumw2();
530
531
1f4e9722 532 hisname.Form("hSigPt%dTC",i);
993bfba1 533 fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 534 fMassHistTC[index]->Sumw2();
bb69f127 535 hisname.Form("hSigPt%dTCPlus",i);
536 fMassHistTCPlus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
537 fMassHistTCPlus[index]->Sumw2();
538 hisname.Form("hSigPt%dTCMinus",i);
539 fMassHistTCMinus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
540 fMassHistTCMinus[index]->Sumw2();
ba9ae5b2 541
82bb8d4b 542 hisname.Form("hLSPt%dLCnw",i);
993bfba1 543 fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 544 fMassHistLS[indexLS]->Sumw2();
545 hisname.Form("hLSPt%dTCnw",i);
993bfba1 546 fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 547 fMassHistLSTC[indexLS]->Sumw2();
548
ba9ae5b2 549
550
551 hisname.Form("hCosPSigPt%dLS",i);
993bfba1 552 fCosPHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.5,1.);
ba9ae5b2 553 fCosPHistLS[index]->Sumw2();
554 hisname.Form("hDLenSigPt%dLS",i);
993bfba1 555 fDLenHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.5);
ba9ae5b2 556 fDLenHistLS[index]->Sumw2();
557 hisname.Form("hSumd02SigPt%dLS",i);
993bfba1 558 fSumd02HistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,1.);
ba9ae5b2 559 fSumd02HistLS[index]->Sumw2();
560 hisname.Form("hSigVertSigPt%dLS",i);
993bfba1 561 fSigVertHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.1);
ba9ae5b2 562 fSigVertHistLS[index]->Sumw2();
563 hisname.Form("hPtMaxSigPt%dLS",i);
993bfba1 564 fPtMaxHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.5,5.);
ba9ae5b2 565 fPtMaxHistLS[index]->Sumw2();
566
567 hisname.Form("hDCASigPt%dLS",i);
993bfba1 568 fDCAHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.1);
ba9ae5b2 569 fDCAHistLS[index]->Sumw2();
570
571
572
82bb8d4b 573 index=GetBackgroundHistoIndex(i);
574 indexLS++;
575 hisname.Form("hBkgPt%d",i);
993bfba1 576 fMassHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 577 fMassHist[index]->Sumw2();
ba9ae5b2 578 hisname.Form("hCosPBkgPt%d",i);
993bfba1 579 fCosPHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.5,1.);
ba9ae5b2 580 fCosPHist[index]->Sumw2();
581 hisname.Form("hDLenBkgPt%d",i);
993bfba1 582 fDLenHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.5);
ba9ae5b2 583 fDLenHist[index]->Sumw2();
584 hisname.Form("hSumd02BkgPt%d",i);
993bfba1 585 fSumd02Hist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,1.);
ba9ae5b2 586 fSumd02Hist[index]->Sumw2();
587 hisname.Form("hSigVertBkgPt%d",i);
993bfba1 588 fSigVertHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.1);
ba9ae5b2 589 fSigVertHist[index]->Sumw2();
590 hisname.Form("hPtMaxBkgPt%d",i);
993bfba1 591 fPtMaxHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.5,5.);
ba9ae5b2 592 fPtMaxHist[index]->Sumw2();
593
594 hisname.Form("hDCABkgPt%d",i);
993bfba1 595 fDCAHist[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.1);
ba9ae5b2 596 fDCAHist[index]->Sumw2();
597
598
1f4e9722 599 hisname.Form("hBkgPt%dTC",i);
993bfba1 600 fMassHistTC[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 601 fMassHistTC[index]->Sumw2();
bb69f127 602 hisname.Form("hBkgPt%dTCPlus",i);
603 fMassHistTCPlus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
604 fMassHistTCPlus[index]->Sumw2();
605 hisname.Form("hBkgPt%dTCMinus",i);
606 fMassHistTCMinus[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
607 fMassHistTCMinus[index]->Sumw2();
ba9ae5b2 608
82bb8d4b 609 hisname.Form("hLSPt%dLCntrip",i);
993bfba1 610 fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 611 fMassHistLS[indexLS]->Sumw2();
612 hisname.Form("hLSPt%dTCntrip",i);
993bfba1 613 fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 614 fMassHistLSTC[indexLS]->Sumw2();
615
ba9ae5b2 616
617 hisname.Form("hCosPBkgPt%dLS",i);
993bfba1 618 fCosPHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.5,1.);
ba9ae5b2 619 fCosPHistLS[index]->Sumw2();
620 hisname.Form("hDLenBkgPt%dLS",i);
993bfba1 621 fDLenHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.5);
ba9ae5b2 622 fDLenHistLS[index]->Sumw2();
623 hisname.Form("hSumd02BkgPt%dLS",i);
993bfba1 624 fSumd02HistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,1.);
ba9ae5b2 625 fSumd02HistLS[index]->Sumw2();
626 hisname.Form("hSigVertBkgPt%dLS",i);
993bfba1 627 fSigVertHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.1);
ba9ae5b2 628 fSigVertHistLS[index]->Sumw2();
629 hisname.Form("hPtMaxBkgPt%dLS",i);
993bfba1 630 fPtMaxHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.5,5.);
ba9ae5b2 631 fPtMaxHistLS[index]->Sumw2();
632 hisname.Form("hDCABkgPt%dLS",i);
993bfba1 633 fDCAHistLS[index]=new TH1F(hisname.Data(),hisname.Data(),nbins,0.,0.1);
ba9ae5b2 634 fDCAHistLS[index]->Sumw2();
635
636
82bb8d4b 637 indexLS++;
638 hisname.Form("hLSPt%dLCntripsinglecut",i);
993bfba1 639 fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 640 fMassHistLS[indexLS]->Sumw2();
641 hisname.Form("hLSPt%dTCntripsinglecut",i);
993bfba1 642 fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 643 fMassHistLSTC[indexLS]->Sumw2();
644
645 indexLS++;
646 hisname.Form("hLSPt%dLCspc",i);
993bfba1 647 fMassHistLS[indexLS]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 648 fMassHistLS[indexLS]->Sumw2();
649 hisname.Form("hLSPt%dTCspc",i);
993bfba1 650 fMassHistLSTC[indexLS]=new TH1F(hisname.Data(),hisname.Data(),nbins,fLowmasslimit,fUpmasslimit);
82bb8d4b 651 fMassHistLSTC[indexLS]->Sumw2();
1f4e9722 652 }
8c34bd86 653
82bb8d4b 654 for(Int_t i=0; i<3*fNPtBins; i++){
655 fOutput->Add(fMassHist[i]);
ba9ae5b2 656 fOutput->Add(fCosPHist[i]);
657 fOutput->Add(fDLenHist[i]);
658 fOutput->Add(fSumd02Hist[i]);
659 fOutput->Add(fSigVertHist[i]);
660 fOutput->Add(fPtMaxHist[i]);
661 fOutput->Add(fDCAHist[i]);
82bb8d4b 662 fOutput->Add(fMassHistTC[i]);
bb69f127 663 fOutput->Add(fMassHistTCPlus[i]);
664 fOutput->Add(fMassHistTCMinus[i]);
82bb8d4b 665 }
ba9ae5b2 666 for(Int_t i=0; i<3*fNPtBins&&fDoLS; i++){
667 fOutput->Add(fCosPHistLS[i]);
668 fOutput->Add(fDLenHistLS[i]);
669 fOutput->Add(fSumd02HistLS[i]);
670 fOutput->Add(fSigVertHistLS[i]);
671 fOutput->Add(fPtMaxHistLS[i]);
672 fOutput->Add(fDCAHistLS[i]);
673 }
82bb8d4b 674 for(Int_t i=0; i<5*fNPtBins&&fDoLS; i++){
675 fOutput->Add(fMassHistLS[i]);
676 fOutput->Add(fMassHistLSTC[i]);
677 }
ba9ae5b2 678
26da4c42 679
5fc4893f 680 fHistNEvents = new TH1F("fHistNEvents", "number of events ",9,-0.5,8.5);
de1b82ff 681 fHistNEvents->GetXaxis()->SetBinLabel(1,"nEventsAnal");
26da4c42 682 fHistNEvents->GetXaxis()->SetBinLabel(2,"nEvents with good vertex");
de1b82ff 683 fHistNEvents->GetXaxis()->SetBinLabel(3,"nEvents with PbPb HM trigger");
684 fHistNEvents->GetXaxis()->SetBinLabel(4,"no. of Rejected pileup events");
685 fHistNEvents->GetXaxis()->SetBinLabel(5,"no. of candidate");
686 fHistNEvents->GetXaxis()->SetBinLabel(6,"no. of D+ after loose cuts");
687 fHistNEvents->GetXaxis()->SetBinLabel(7,"no. of D+ after tight cuts");
37fc80e6 688 fHistNEvents->GetXaxis()->SetBinLabel(8,"no. of out centrality events");
5fc4893f 689 fHistNEvents->GetXaxis()->SetBinLabel(9,"no. of cand wo bitmask");
37fc80e6 690
26da4c42 691 fHistNEvents->GetXaxis()->SetNdivisions(1,kFALSE);
692
ba9ae5b2 693 fHistNEvents->Sumw2();
694 fHistNEvents->SetMinimum(0);
695 fOutput->Add(fHistNEvents);
595cc7e2 696
4d797793 697 fPtVsMass=new TH2F("hPtVsMass","PtVsMass (prod. cuts)",nbins,fLowmasslimit,fUpmasslimit,200,0.,20.);
698 fPtVsMassTC=new TH2F("hPtVsMassTC","PtVsMass (analysis cuts)",nbins,fLowmasslimit,fUpmasslimit,200,0.,20.);
595cc7e2 699 fYVsPt=new TH2F("hYVsPt","YvsPt (prod. cuts)",40,0.,20.,80,-2.,2.);
700 fYVsPtTC=new TH2F("hYVsPtTC","YvsPt (analysis cuts)",40,0.,20.,80,-2.,2.);
701 fYVsPtSig=new TH2F("hYVsPtSig","YvsPt (MC, only sig., prod. cuts)",40,0.,20.,80,-2.,2.);
702 fYVsPtSigTC=new TH2F("hYVsPtSigTC","YvsPt (MC, only Sig, analysis cuts)",40,0.,20.,80,-2.,2.);
703
704 fOutput->Add(fPtVsMass);
705 fOutput->Add(fPtVsMassTC);
706 fOutput->Add(fYVsPt);
707 fOutput->Add(fYVsPtTC);
708 fOutput->Add(fYVsPtSig);
709 fOutput->Add(fYVsPtSigTC);
4afc48a2 710
a96083b9 711
712 //Counter for Normalization
e0fdfa52 713 TString normName="NormalizationCounter";
714 AliAnalysisDataContainer *cont = GetOutputSlot(3)->GetContainer();
715 if(cont)normName=(TString)cont->GetName();
716 fCounter = new AliNormalizationCounter(normName.Data());
de1b82ff 717 fCounter->SetRejectPileUp(fRDCutsProduction->GetOptPileUp());
a96083b9 718
1f4e9722 719 if(fFillNtuple){
e0fdfa52 720 OpenFile(4); // 4 is the slot number of the ntuple
ba9ae5b2 721
722 fNtupleDplus = new TNtuple("fNtupleDplus","D +","pdg:Px:Py:Pz:PtTrue:VxTrue:VyTrue:VzTrue:Ptpi:PtK:Ptpi2:PtRec:PointingAngle:DecLeng:VxRec:VyRec:VzRec:InvMass:sigvert:d0Pi:d0K:d0Pi2:dca:d0square");
723
1f4e9722 724 }
ba9ae5b2 725
d486095a 726 return;
727}
728
729//________________________________________________________________________
730void AliAnalysisTaskSEDplus::UserExec(Option_t */*option*/)
731{
732 // Execute analysis for current event:
733 // heavy flavor candidates association to MC truth
734
4c230f6d 735 AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
4afc48a2 736
a96083b9 737
738
b557eb43 739 TClonesArray *array3Prong = 0;
4afc48a2 740 TClonesArray *arrayLikeSign =0;
b557eb43 741 if(!aod && AODEvent() && IsStandardAOD()) {
742 // In case there is an AOD handler writing a standard AOD, use the AOD
743 // event in memory rather than the input (ESD) event.
744 aod = dynamic_cast<AliAODEvent*> (AODEvent());
4c230f6d 745 // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
746 // have to taken from the AOD event hold by the AliAODExtension
b557eb43 747 AliAODHandler* aodHandler = (AliAODHandler*)
748 ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
749 if(aodHandler->GetExtensions()) {
750 AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
751 AliAODEvent *aodFromExt = ext->GetAOD();
752 array3Prong=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
4afc48a2 753 arrayLikeSign=(TClonesArray*)aodFromExt->GetList()->FindObject("LikeSign3Prong");
b557eb43 754 }
dc222f77 755 } else if(aod) {
b557eb43 756 array3Prong=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
4afc48a2 757 arrayLikeSign=(TClonesArray*)aod->GetList()->FindObject("LikeSign3Prong");
b557eb43 758 }
8931c313 759
90993728 760 if(!aod || !array3Prong) {
d486095a 761 printf("AliAnalysisTaskSEDplus::UserExec: Charm3Prong branch not found!\n");
762 return;
763 }
3a3c90f4 764 if(!arrayLikeSign && fDoLS) {
4afc48a2 765 printf("AliAnalysisTaskSEDplus::UserExec: LikeSign3Prong branch not found!\n");
82bb8d4b 766 return;
4afc48a2 767 }
768
7c23877d 769
770 // fix for temporary bug in ESDfilter
771 // the AODs with null vertex pointer didn't pass the PhysSel
a96083b9 772 if(!aod->GetPrimaryVertex()||TMath::Abs(aod->GetMagneticField())<0.001) return;
773 fCounter->StoreEvent(aod,fReadMC);
de1b82ff 774 fHistNEvents->Fill(0); // count event
37fc80e6 775
776 Bool_t isEvSel=fRDCutsAnalysis->IsEventSelected(aod);
777 Float_t centrality=aod->GetNTracks();//fRDCutsAnalysis->GetCentrality(aod);
778 fHistCentrality[0]->Fill(centrality);
de1b82ff 779 // trigger class for PbPb C0SMH-B-NOPF-ALLNOTRD
780 TString trigclass=aod->GetFiredTriggerClasses();
781 if(trigclass.Contains("C0SMH-B-NOPF-ALLNOTRD")||trigclass.Contains("C0SMH-B-NOPF-ALL")) fHistNEvents->Fill(2);
782 if(fRDCutsAnalysis->GetWhyRejection()==1)fHistNEvents->Fill(3);
37fc80e6 783 if(fRDCutsAnalysis->GetWhyRejection()==2){fHistNEvents->Fill(7);fHistCentrality[2]->Fill(centrality);}
784
785 // Post the data already here
de1b82ff 786 PostData(1,fOutput);
37fc80e6 787 if(!isEvSel)return;
788
789 fHistCentrality[1]->Fill(centrality);
790 fHistNEvents->Fill(1);
791
4afc48a2 792 TClonesArray *arrayMC=0;
793 AliAODMCHeader *mcHeader=0;
d486095a 794
795 // AOD primary vertex
1f4e9722 796 AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
797 // vtx1->Print();
26da4c42 798 TString primTitle = vtx1->GetTitle();
799 //if(primTitle.Contains("VertexerTracks") && vtx1->GetNContributors()>0)fHistNEvents->Fill(2);
800
d486095a 801 // load MC particles
4afc48a2 802 if(fReadMC){
803
804 arrayMC = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
805 if(!arrayMC) {
806 printf("AliAnalysisTaskSEDplus::UserExec: MC particles branch not found!\n");
90993728 807 return;
4afc48a2 808 }
1f4e9722 809
d486095a 810 // load MC header
4afc48a2 811 mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
812 if(!mcHeader) {
82bb8d4b 813 printf("AliAnalysisTaskSEDplus::UserExec: MC header branch not found!\n");
814 return;
4afc48a2 815 }
d486095a 816 }
4afc48a2 817
1f4e9722 818 Int_t n3Prong = array3Prong->GetEntriesFast();
2bf2e62b 819 printf("Number of D+->Kpipi: %d and of tracks: %d\n",n3Prong,aod->GetNumberOfTracks());
4afc48a2 820
821
822 Int_t nOS=0;
82bb8d4b 823 Int_t index;
1f4e9722 824 Int_t pdgDgDplustoKpipi[3]={321,211,211};
4c230f6d 825 // Double_t cutsDplus[12]={0.2,0.4,0.4,0.,0.,0.01,0.06,0.02,0.,0.85,0.,10000000000.};//TO REMOVE
826 //Double_t *cutsDplus = new (Double_t*)fRDCuts->GetCuts();
a96083b9 827 Int_t nSelectedloose=0,nSelectedtight=0;
1f4e9722 828 for (Int_t i3Prong = 0; i3Prong < n3Prong; i3Prong++) {
829 AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)array3Prong->UncheckedAt(i3Prong);
de1b82ff 830 fHistNEvents->Fill(4);
5fc4893f 831 if(fUseBit && !d->HasSelectionBit(AliRDHFCuts::kDplusCuts)){
832 fHistNEvents->Fill(8);
833 continue;
834 }
1f4e9722 835 Bool_t unsetvtx=kFALSE;
836 if(!d->GetOwnPrimaryVtx()){
837 d->SetOwnPrimaryVtx(vtx1);
838 unsetvtx=kTRUE;
839 }
4afc48a2 840
1275bf81 841 if(fRDCutsProduction->IsSelected(d,AliRDHFCuts::kCandidate,aod)) {
a96083b9 842
2bf2e62b 843
26da4c42 844
a96083b9 845
82bb8d4b 846 Int_t iPtBin = -1;
1f4e9722 847 Double_t ptCand = d->Pt();
4c230f6d 848
82bb8d4b 849 for(Int_t ibin=0;ibin<fNPtBins&&iPtBin<0&&ptCand>fArrayBinLimits[0]&&ptCand<fArrayBinLimits[fNPtBins];ibin++){
850 if(ptCand<fArrayBinLimits[ibin+1])iPtBin=ibin;
1f4e9722 851 }
82bb8d4b 852
1275bf81 853 Int_t passTightCuts=fRDCutsAnalysis->IsSelected(d,AliRDHFCuts::kCandidate,aod);
a96083b9 854
2bf2e62b 855
4afc48a2 856 Int_t labDp=-1;
857 Float_t deltaPx=0.;
858 Float_t deltaPy=0.;
859 Float_t deltaPz=0.;
860 Float_t truePt=0.;
861 Float_t xDecay=0.;
862 Float_t yDecay=0.;
863 Float_t zDecay=0.;
864 Float_t pdgCode=-2;
865 if(fReadMC){
866 labDp = d->MatchToMC(411,arrayMC,3,pdgDgDplustoKpipi);
867 if(labDp>=0){
868 AliAODMCParticle *partDp = (AliAODMCParticle*)arrayMC->At(labDp);
869 AliAODMCParticle *dg0 = (AliAODMCParticle*)arrayMC->At(partDp->GetDaughter(0));
870 deltaPx=partDp->Px()-d->Px();
871 deltaPy=partDp->Py()-d->Py();
872 deltaPz=partDp->Pz()-d->Pz();
873 truePt=partDp->Pt();
874 xDecay=dg0->Xv();
875 yDecay=dg0->Yv();
876 zDecay=dg0->Zv();
877 pdgCode=TMath::Abs(partDp->GetPdgCode());
878 }else{
879 pdgCode=-1;
880 }
881 }
1f4e9722 882 Double_t invMass=d->InvMassDplus();
595cc7e2 883 Double_t rapid=d->YDplus();
884 fYVsPt->Fill(ptCand,rapid);
885 if(passTightCuts) fYVsPtTC->Fill(ptCand,rapid);
886 Bool_t isFidAcc=fRDCutsAnalysis->IsInFiducialAcceptance(ptCand,rapid);
887 if(isFidAcc){
888 fPtVsMass->Fill(invMass,ptCand);
889 if(passTightCuts) fPtVsMassTC->Fill(invMass,ptCand);
890 }
ba9ae5b2 891 Float_t tmp[24];
4afc48a2 892 if(fFillNtuple){
893 tmp[0]=pdgCode;
894 tmp[1]=deltaPx;
895 tmp[2]=deltaPy;
896 tmp[3]=deltaPz;
897 tmp[4]=truePt;
898 tmp[5]=xDecay;
899 tmp[6]=yDecay;
900 tmp[7]=zDecay;
901 tmp[8]=d->PtProng(0);
902 tmp[9]=d->PtProng(1);
903 tmp[10]=d->PtProng(2);
904 tmp[11]=d->Pt();
905 tmp[12]=d->CosPointingAngle();
906 tmp[13]=d->DecayLength();
907 tmp[14]=d->Xv();
908 tmp[15]=d->Yv();
909 tmp[16]=d->Zv();
910 tmp[17]=d->InvMassDplus();
911 tmp[18]=d->GetSigmaVert();
912 tmp[19]=d->Getd0Prong(0);
913 tmp[20]=d->Getd0Prong(1);
ba9ae5b2 914 tmp[21]=d->Getd0Prong(2);
915 tmp[22]=d->GetDCA();
916 tmp[23]=d->Prodd0d0();
4afc48a2 917 fNtupleDplus->Fill(tmp);
de1b82ff 918 PostData(4,fNtupleDplus);
4afc48a2 919 }
ba9ae5b2 920 Double_t dlen=d->DecayLength();
921 Double_t cosp=d->CosPointingAngle();
922 Double_t sumD02=d->Getd0Prong(0)*d->Getd0Prong(0)+d->Getd0Prong(1)*d->Getd0Prong(1)+d->Getd0Prong(2)*d->Getd0Prong(2);
0db12536 923 Double_t dca=d->GetDCA();
924 Double_t sigvert=d->GetSigmaVert();
925 Double_t ptmax=0;
ba9ae5b2 926 for(Int_t i=0;i<3;i++){
927 if(d->PtProng(i)>ptmax)ptmax=d->PtProng(i);
928 }
82bb8d4b 929 if(iPtBin>=0){
930
931 index=GetHistoIndex(iPtBin);
595cc7e2 932 if(isFidAcc){
de1b82ff 933 fHistNEvents->Fill(5);
a96083b9 934 nSelectedloose++;
595cc7e2 935 fMassHist[index]->Fill(invMass);
936 fCosPHist[index]->Fill(cosp);
937 fDLenHist[index]->Fill(dlen);
938 fSumd02Hist[index]->Fill(sumD02);
939 fSigVertHist[index]->Fill(sigvert);
940 fPtMaxHist[index]->Fill(ptmax);
941 fDCAHist[index]->Fill(dca);
942
de1b82ff 943 if(passTightCuts){ fHistNEvents->Fill(6);
a96083b9 944 nSelectedtight++;
595cc7e2 945 fMassHistTC[index]->Fill(invMass);
bb69f127 946 if(d->GetCharge()>0) fMassHistTCPlus[index]->Fill(invMass);
947 else if(d->GetCharge()<0) fMassHistTCMinus[index]->Fill(invMass);
595cc7e2 948 }
82bb8d4b 949 }
950
951 if(fReadMC){
952 if(labDp>=0) {
953 index=GetSignalHistoIndex(iPtBin);
595cc7e2 954 if(isFidAcc){
a96083b9 955 Float_t factor[3]={1.,1.,1.};
956 if(fUseStrangeness){
957 for(Int_t iprong=0;iprong<3;iprong++){
958 AliAODTrack *trad = (AliAODTrack*)d->GetDaughter(iprong);
959 Int_t labd= trad->GetLabel();
960 if(labd>=0){
961 AliAODMCParticle *dau = (AliAODMCParticle*)arrayMC->At(labd);
962 if(dau){
963 Int_t labm = dau->GetMother();
964 if(labm>=0){
965 AliAODMCParticle *mot = (AliAODMCParticle*)arrayMC->At(labm);
966 if(mot){
967 if(TMath::Abs(mot->GetPdgCode())==310 || TMath::Abs(mot->GetPdgCode())==130 || TMath::Abs(mot->GetPdgCode())==321){ //K0_S, K0_L, K^+-
968 if(d->PtProng(iprong)<=1)factor[iprong]=1./.7;
969 else factor[iprong]=1./.6;
970 // fNentries->Fill(12);
971 }
972 if(TMath::Abs(mot->GetPdgCode())==3122) { //Lambda
973 factor[iprong]=1./0.25;
974 // fNentries->Fill(13);
975 }//if 3122
976 }//if(mot)
977 }//if labm>0
978 }//if(dau)
979 }//if labd>=0
980 }//prong loop
981 }
982 Float_t fact=1.;for(Int_t k=0;k<3;k++)fact=fact*factor[k];
595cc7e2 983 fMassHist[index]->Fill(invMass);
a96083b9 984 fCosPHist[index]->Fill(cosp,fact);
985 fDLenHist[index]->Fill(dlen,fact);
986 Float_t sumd02s=d->Getd0Prong(0)*d->Getd0Prong(0)*factor[0]*factor[0]+d->Getd0Prong(1)*d->Getd0Prong(1)*factor[1]*factor[1]+d->Getd0Prong(2)*d->Getd0Prong(2)*factor[2]*factor[2];
987 fSumd02Hist[index]->Fill(sumd02s);
988 fSigVertHist[index]->Fill(sigvert,fact);
989 fPtMaxHist[index]->Fill(ptmax,fact);
990 fDCAHist[index]->Fill(dca,fact);
595cc7e2 991 if(passTightCuts){
992 fMassHistTC[index]->Fill(invMass);
bb69f127 993 if(d->GetCharge()>0) fMassHistTCPlus[index]->Fill(invMass);
994 else if(d->GetCharge()<0) fMassHistTCMinus[index]->Fill(invMass);
595cc7e2 995 }
996 }
997 fYVsPtSig->Fill(ptCand,rapid);
998 if(passTightCuts) fYVsPtSigTC->Fill(ptCand,rapid);
82bb8d4b 999 }else{
1000 index=GetBackgroundHistoIndex(iPtBin);
595cc7e2 1001 if(isFidAcc){
a96083b9 1002 Float_t factor[3]={1.,1.,1.};
1003 if(fUseStrangeness){
1004 for(Int_t iprong=0;iprong<3;iprong++){
1005 AliAODTrack *trad = (AliAODTrack*)d->GetDaughter(iprong);
1006 Int_t labd= trad->GetLabel();
1007 if(labd>=0){
1008 AliAODMCParticle *dau = (AliAODMCParticle*)arrayMC->At(labd);
1009 if(dau){
1010 Int_t labm = dau->GetMother();
1011 if(labm>=0){
1012 AliAODMCParticle *mot = (AliAODMCParticle*)arrayMC->At(labm);
1013 if(mot){
1014 if(TMath::Abs(mot->GetPdgCode())==310 || TMath::Abs(mot->GetPdgCode())==130 || TMath::Abs(mot->GetPdgCode())==321){ //K0_S, K0_L, K^+-
1015 if(d->PtProng(iprong)<=1)factor[iprong]=1./.7;
1016 else factor[iprong]=1./.6;
1017 // fNentries->Fill(12);
1018 }
1019 if(TMath::Abs(mot->GetPdgCode())==3122) { //Lambda
1020 factor[iprong]=1./0.25;
1021 // fNentries->Fill(13);
1022 }//if 3122
1023 }//if(mot)
1024 }//if labm>0
1025 }//if(dau)
1026 }//if labd>=0
1027 }//prong loop
1028 }
1029 Float_t fact=1.;for(Int_t k=0;k<3;k++)fact=fact*factor[k];
595cc7e2 1030 fMassHist[index]->Fill(invMass);
a96083b9 1031 fCosPHist[index]->Fill(cosp,fact);
1032 fDLenHist[index]->Fill(dlen,fact);
1033 Float_t sumd02s=d->Getd0Prong(0)*d->Getd0Prong(0)*factor[0]*factor[0]+d->Getd0Prong(1)*d->Getd0Prong(1)*factor[1]*factor[1]+d->Getd0Prong(2)*d->Getd0Prong(2)*factor[2]*factor[2];
1034 fSumd02Hist[index]->Fill(sumd02s);
1035 fSigVertHist[index]->Fill(sigvert,fact);
1036 fPtMaxHist[index]->Fill(ptmax,fact);
1037 fDCAHist[index]->Fill(dca,fact);
595cc7e2 1038 if(passTightCuts){
1039 fMassHistTC[index]->Fill(invMass);
bb69f127 1040 if(d->GetCharge()>0) fMassHistTCPlus[index]->Fill(invMass);
1041 else if(d->GetCharge()<0) fMassHistTCMinus[index]->Fill(invMass);
595cc7e2 1042 }
82bb8d4b 1043 }
4afc48a2 1044 }
2bf2e62b 1045
1046 }
1047
1048 }
1049
d486095a 1050 }
1f4e9722 1051 if(unsetvtx) d->UnsetOwnPrimaryVtx();
1052 }
a96083b9 1053 fCounter->StoreCandidates(aod,nSelectedloose,kTRUE);
1054 fCounter->StoreCandidates(aod,nSelectedtight,kFALSE);
4c230f6d 1055
4afc48a2 1056 //start LS analysis
1057 if(fDoLS && arrayLikeSign) LSAnalysis(array3Prong,arrayLikeSign,aod,vtx1,nOS);
82bb8d4b 1058
4afc48a2 1059 PostData(1,fOutput);
a96083b9 1060 PostData(3,fCounter);
d486095a 1061 return;
1062}
1063
4afc48a2 1064
1065
82bb8d4b 1066//________________________________________________________________________
d486095a 1067void AliAnalysisTaskSEDplus::Terminate(Option_t */*option*/)
1068{
1069 // Terminate analysis
1070 //
1071 if(fDebug > 1) printf("AnalysisTaskSEDplus: Terminate() \n");
82bb8d4b 1072
d486095a 1073 fOutput = dynamic_cast<TList*> (GetOutputData(1));
1074 if (!fOutput) {
1075 printf("ERROR: fOutput not available\n");
1076 return;
1077 }
595cc7e2 1078 fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));
82bb8d4b 1079
595cc7e2 1080 TString hisname;
1081 Int_t index=0;
ba9ae5b2 1082
595cc7e2 1083 for(Int_t i=0;i<fNPtBins;i++){
82bb8d4b 1084 index=GetHistoIndex(i);
ba9ae5b2 1085
37fc80e6 1086 hisname.Form("hMassPt%dTC",i);
82bb8d4b 1087 fMassHistTC[index]=dynamic_cast<TH1F*>(fOutput->FindObject(hisname.Data()));
37fc80e6 1088 }
ba9ae5b2 1089
13808a30 1090 TCanvas *c1=new TCanvas("c1","D+ invariant mass distribution",500,500);
1091 c1->cd();
1092 TH1F *hMassPt=(TH1F*)fOutput->FindObject("hMassPt3TC");
1093 hMassPt->SetLineColor(kBlue);
1094 hMassPt->SetXTitle("M[GeV/c^{2}]");
1095 hMassPt->Draw();
1096
595cc7e2 1097 return;
d486095a 1098}