AliPhysics  6cf2591 (6cf2591)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskSEDvsEventShapes.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 /* $Id$ */
17 
18 //*************************************************************************
19 // Class AliAnalysisTaskSEDvsEventShapes
20 // AliAnalysisTaskSE for the D meson vs. Event shape analysis in different mutiplicity window
21 // Authors: Renu Bala, Manoj Bhanudas Jadhav
22 //*************************************************************************
23 
24 #include <TClonesArray.h>
25 #include <TCanvas.h>
26 #include <TList.h>
27 #include <TMath.h>
28 #include <TString.h>
29 #include <TDatabasePDG.h>
30 #include <TH1F.h>
31 #include <TH2F.h>
32 #include <TH3F.h>
33 #include <THnSparse.h>
34 #include <TProfile.h>
35 #include "AliAnalysisManager.h"
36 #include "AliRDHFCuts.h"
39 #include "AliRDHFCutsD0toKpi.h"
40 #include "AliAODHandler.h"
41 #include "AliAODEvent.h"
42 #include "AliAODVertex.h"
43 #include "AliAODTrack.h"
44 #include "AliAODRecoDecayHF.h"
45 #include "AliAODRecoCascadeHF.h"
46 #include "AliAnalysisVertexingHF.h"
47 #include "AliAnalysisTaskSE.h"
50 #include "AliVertexingHFUtils.h"
51 #include "AliAODVZERO.h"
52 #include "AliESDUtils.h"
54 
55 //________________________________________________________________________
58 fOutput(0),
59 fListCuts(0),
60 fOutputCounters(0),
61 fListProfiles(0),
62 fOutputEffCorr(0),
63 fHistNEvents(0),
64 fHistNtrVsZvtx(0),
65 fHistNtrCorrVsZvtx(0),
66 fHistNtrVsnTrackEvWithCand(0),
67 fHistNtrVsSo(0),
68 fHistNtrCorrVsSo(0),
69 fHistNtrVsSpheri(0),
70 fHistNtrCorrVsSpheri(0),
71 fHistNtrVsNchMC(0),
72 fHistNtrCorrVsNchMC(0),
73 fHistNtrVsNchMCPrimary(0),
74 fHistNtrCorrVsNchMCPrimary(0),
75 fHistNtrVsNchMCPhysicalPrimary(0),
76 fHistNtrCorrVsNchMCPhysicalPrimary(0),
77 fHistGenPrimaryParticlesInelGt0(0),
78 fHistNchMCVsNchMCPrimaryVsNchMCPhysicalPrimary(0),
79 fHistNtrCorrPSSel(0),
80 fHistNtrCorrEvSel(0),
81 fHistNtrCorrEvWithCand(0),
82 fHistNtrCorrEvWithD(0),
83 fHistnTrackvsEtavsPhi(0),
84 fHistnTrackvsEtavsPhiEvWithCand(0),
85 fHistTrueSovsMeasSo(0),
86 fHistTrueSovsMeasSoEvWithCand(0),
87 fSparseEvtShape(0),
88 fSparseEvtShapewithNoPid(0),
89 fSparseEvtShapePrompt(0),
90 fSparseEvtShapeFeeddown(0),
91 fSparseEvtShapeRecSphero(0),
92 fMCAccGenPrompt(0),
93 fMCAccGenFeeddown(0),
94 fMCRecoPrompt(0),
95 fMCRecoFeeddown(0),
96 fMCRecoBothPromptFD(0),
97 fMCAccGenPromptSpheri(0),
98 fMCAccGenFeeddownSpheri(0),
99 fMCRecoPromptSpheri(0),
100 fMCRecoFeeddownSpheri(0),
101 fMCRecoBothPromptFDSpheri(0),
102 fMCAccGenPromptEvSel(0),
103 fMCAccGenFeeddownEvSel(0),
104 fUpmasslimit(1.965),
105 fLowmasslimit(1.765),
106 fNMassBins(200),
107 fRDCutsAnalysis(0),
108 fCounterC(0),
109 fCounterU(0),
110 fCounterCandidates(0),
111 fDoImpPar(kFALSE),
112 fNImpParBins(400),
113 fLowerImpPar(-2000.),
114 fHigherImpPar(2000.),
115 fReadMC(kFALSE),
116 fMCOption(0),
117 fisPPbData(kFALSE),
118 fUseBit(kTRUE),
119 fSubtractTrackletsFromDau(kFALSE),
120 fCalculateSphericity(kFALSE),
121 fRecomputeSpherocity(kFALSE),
122 fRemoveD0fromDstar(kFALSE),
123 fUseNchWeight(0),
124 fHistoMCNch(0),
125 fHistoMeasNch(0),
126 fUsePtWeight(kFALSE),
127 fWeight(1.),
128 fRefMult(9.26),
129 fPdgMeson(411),
130 fMultiplicityEstimator(kNtrk10),
131 fMCPrimariesEstimator(kEta10),
132 fDoVZER0ParamVertexCorr(1),
133 fFillSoSparseChecks(0),
134 fUseQuarkTag(kTRUE),
135 fFillTrackHisto(kFALSE),
136 fEtaAccCut(0.9),
137 fPtAccCut(0.1),
138 fetaMin(-0.8),
139 fetaMax(0.8),
140 fptMin(0.15),
141 fptMax(10.),
142 fminMult(3),
143 ffiltbit1(256),
144 ffiltbit2(512),
145 fphiStepSizeDeg(0.1)
146 {
147  // Default constructor
148  for(Int_t i=0; i<5; i++) fHistMassPtImpPar[i]=0;
149  for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=0;
150 }
151 
152 //________________________________________________________________________
154 AliAnalysisTaskSE(name),
155 fOutput(0),
156 fListCuts(0),
157 fOutputCounters(0),
158 fListProfiles(0),
159 fOutputEffCorr(0),
160 fHistNEvents(0),
161 fHistNtrVsZvtx(0),
162 fHistNtrCorrVsZvtx(0),
163 fHistNtrVsnTrackEvWithCand(0),
164 fHistNtrVsSo(0),
165 fHistNtrCorrVsSo(0),
166 fHistNtrVsSpheri(0),
167 fHistNtrCorrVsSpheri(0),
168 fHistNtrVsNchMC(0),
169 fHistNtrCorrVsNchMC(0),
170 fHistNtrVsNchMCPrimary(0),
171 fHistNtrCorrVsNchMCPrimary(0),
172 fHistNtrVsNchMCPhysicalPrimary(0),
173 fHistNtrCorrVsNchMCPhysicalPrimary(0),
174 fHistGenPrimaryParticlesInelGt0(0),
175 fHistNchMCVsNchMCPrimaryVsNchMCPhysicalPrimary(0),
176 fHistNtrCorrPSSel(0),
177 fHistNtrCorrEvSel(0),
178 fHistNtrCorrEvWithCand(0),
179 fHistNtrCorrEvWithD(0),
180 fHistnTrackvsEtavsPhi(0),
181 fHistnTrackvsEtavsPhiEvWithCand(0),
182 fHistTrueSovsMeasSo(0),
183 fHistTrueSovsMeasSoEvWithCand(0),
184 fSparseEvtShape(0),
185 fSparseEvtShapewithNoPid(0),
186 fSparseEvtShapePrompt(0),
187 fSparseEvtShapeFeeddown(0),
188 fSparseEvtShapeRecSphero(0),
189 fMCAccGenPrompt(0),
190 fMCAccGenFeeddown(0),
191 fMCRecoPrompt(0),
192 fMCRecoFeeddown(0),
193 fMCRecoBothPromptFD(0),
194 fMCAccGenPromptSpheri(0),
195 fMCAccGenFeeddownSpheri(0),
196 fMCRecoPromptSpheri(0),
197 fMCRecoFeeddownSpheri(0),
198 fMCRecoBothPromptFDSpheri(0),
199 fMCAccGenPromptEvSel(0),
200 fMCAccGenFeeddownEvSel(0),
201 fUpmasslimit(1.965),
202 fLowmasslimit(1.765),
203 fNMassBins(200),
204 fRDCutsAnalysis(cuts),
205 fCounterC(0),
206 fCounterU(0),
207 fCounterCandidates(0),
208 fDoImpPar(kFALSE),
209 fNImpParBins(400),
210 fLowerImpPar(-2000.),
211 fHigherImpPar(2000.),
212 fReadMC(kFALSE),
213 fMCOption(0),
214 fisPPbData(switchPPb),
215 fUseBit(kTRUE),
216 fSubtractTrackletsFromDau(kFALSE),
217 fCalculateSphericity(kFALSE),
218 fRecomputeSpherocity(kFALSE),
219 fRemoveD0fromDstar(kFALSE),
220 fUseNchWeight(0),
221 fHistoMCNch(0),
222 fHistoMeasNch(0),
223 fUsePtWeight(kFALSE),
224 fWeight(1.),
225 fRefMult(9.26),
226 fPdgMeson(pdgMeson),
227 fMultiplicityEstimator(kNtrk10),
228 fMCPrimariesEstimator(kEta10),
229 fDoVZER0ParamVertexCorr(1),
230 fFillSoSparseChecks(0),
231 fUseQuarkTag(kTRUE),
232 fFillTrackHisto(kFALSE),
233 fEtaAccCut(0.9),
234 fPtAccCut(0.1),
235 fetaMin(-0.8),
236 fetaMax(0.8),
237 fptMin(0.15),
238 fptMax(10.),
239 fminMult(3),
240 ffiltbit1(256),
241 ffiltbit2(512),
242 fphiStepSizeDeg(0.1)
243 {
244  //
245  // Standard constructor
246  //
247 
248  for(Int_t i=0; i<5; i++) fHistMassPtImpPar[i]=0;
249  for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=0;
250  if(fPdgMeson==413){
251  fNMassBins=200;
252  SetMassLimits(0.12,0.2);
253  }else{
254  fNMassBins=200;
256  }
257  // Default constructor
258  // Otput slot #1 writes into a TList container
259  DefineOutput(1,TList::Class()); //My private output
260  // Output slot #2 writes cut to private output
261  DefineOutput(2,TList::Class());
262  // Output slot #3 writes cut to private output
263  DefineOutput(3,TList::Class());
264  // Output slot #4 writes cut to private output
265  DefineOutput(4,TList::Class());
266  // Output slot #5 writes cut to private output
267  DefineOutput(5,TList::Class());
268 }
269 //________________________________________________________________________
271 {
272  //
273  // Destructor
274  //
275  delete fOutput;
276  delete fHistNEvents;
277  delete fListCuts;
278  delete fListProfiles;
279  delete fOutputEffCorr;
280  delete fRDCutsAnalysis;
281  delete fCounterC;
282  delete fCounterU;
283  delete fCounterCandidates;
284 
285  for(Int_t i=0; i<4; i++) {
286  if (fMultEstimatorAvg[i]) delete fMultEstimatorAvg[i];
287  }
288 
289  for(Int_t i=0; i<5; i++){
290  delete fHistMassPtImpPar[i];
291  }
292  if(fHistoMCNch) delete fHistoMCNch;
293  if(fHistoMeasNch) delete fHistoMeasNch;
294 }
295 
296 //_________________________________________________________________
298  // set invariant mass limits
299  if(uplimit>lowlimit){
300  fLowmasslimit = lowlimit;
301  fUpmasslimit = uplimit;
302  }else{
303  AliError("Wrong mass limits: upper value should be larger than lower one");
304  }
305 }
306 //_________________________________________________________________
308  // set invariant mass limits
309  Double_t mass=TDatabasePDG::Instance()->GetParticle(TMath::Abs(pdg))->Mass();
310  SetMassLimits(mass-range,mass+range);
311 }
312 //________________________________________________________________________
314  //
315  // Initialization
316  //
317  printf("AnalysisTaskSEDvsMultiplicity_0::Init() \n");
318 
319  if(fUseNchWeight && !fReadMC){ AliFatal("Nch weights can only be used in MC mode"); return; }
320  if(fUsePtWeight && !fReadMC){ AliFatal("pT weights can only be used in MC mode"); return; }
321  if(fUsePtWeight && fUseNchWeight) { AliInfo("Beware, using at the same time pT and Nch weights, please check"); }
322  if(fUseNchWeight && !fHistoMCNch){ AliFatal("Nch weights can only be used without histogram"); return; }
323  if(fUseNchWeight==1 && !fHistoMeasNch) {//Nch weights
324  if(fisPPbData){ AliFatal("Nch weights can only be used with MC and data histogram in pPb"); return; }
325  else CreateMeasuredNchHisto();
326  }
327  if(fUseNchWeight==2 && !fHistoMeasNch){ AliFatal("Ntrk weights can only be used with MC and data histogram"); return; } //for pp, pPb Ntrk weights
328 
329  fListCuts=new TList();
330  fListCuts->SetOwner();
331  fListCuts->SetName("CutsList");
332 
333  if(fPdgMeson==411){
334  AliRDHFCutsDplustoKpipi* copycut=new AliRDHFCutsDplustoKpipi(*(static_cast<AliRDHFCutsDplustoKpipi*>(fRDCutsAnalysis)));
335  copycut->SetName("AnalysisCutsDplus");
336  fListCuts->Add(copycut);
337  }else if(fPdgMeson==421){
338  AliRDHFCutsD0toKpi* copycut=new AliRDHFCutsD0toKpi(*(static_cast<AliRDHFCutsD0toKpi*>(fRDCutsAnalysis)));
339  copycut->SetName("AnalysisCutsDzero");
340  fListCuts->Add(copycut);
341  }else if(fPdgMeson==413){
342  AliRDHFCutsDStartoKpipi* copycut=new AliRDHFCutsDStartoKpipi(*(static_cast<AliRDHFCutsDStartoKpipi*>(fRDCutsAnalysis)));
343  copycut->SetName("AnalysisCutsDStar");
344  fListCuts->Add(copycut);
345  }
348 
349  PostData(2,fListCuts);
350 
351  fListProfiles = new TList();
352  fListProfiles->SetOwner();
353  TString period[4];
354  Int_t nProfiles=4;
355  if (fisPPbData) {period[0]="LHC13b"; period[1]="LHC13c"; nProfiles = 2;}
356  else {period[0]="LHC10b"; period[1]="LHC10c"; period[2]="LHC10d"; period[3]="LHC10e"; nProfiles = 4;}
357 
358  for(Int_t i=0; i<nProfiles; i++){
359  if(fMultEstimatorAvg[i]){
360  TProfile* hprof=new TProfile(*fMultEstimatorAvg[i]);
361  hprof->SetName(Form("ProfileTrkVsZvtx%s\n",period[i].Data()));
362  fListProfiles->Add(hprof);
363  }
364  }
365 
366  PostData(4,fListProfiles);
367 
368  return;
369 }
370 
371 //________________________________________________________________________
373 {
374  // Create the output container
375  //
376  if(fDebug > 1) printf("AnalysisTaskSEDvsMultiplicity::UserCreateOutputObjects() \n");
377 
378  // Several histograms are more conveniently managed in a TList
379  fOutput = new TList();
380  fOutput->SetOwner();
381  fOutput->SetName("OutputHistos");
382 
383  fOutputEffCorr = new TList();
384  fOutputEffCorr->SetOwner();
385  fOutputEffCorr->SetName("OutputEffCorrHistos");
386 
387  Int_t nMultBins = 200;
388  Float_t firstMultBin = -0.5;
389  Float_t lastMultBin = 199.5;
390  Int_t nMultBinsNtrk = nMultBins;
391  Float_t lastMultBinNtrk = lastMultBin;
392  Int_t nMultBinsV0 = 400;
393  Float_t lastMultBinV0 = 799.5;
394  const char *estimatorName="tracklets";
395  if(fisPPbData) {
396  nMultBinsNtrk = 375;
397  lastMultBinNtrk = 374.5;
398  nMultBins = nMultBinsNtrk;
399  lastMultBin = lastMultBinNtrk;
400  }
402  nMultBins = nMultBinsV0;
403  lastMultBin = lastMultBinV0;
404  estimatorName = "vzero";
405  }
406 
407  fHistNtrCorrPSSel = new TH1F("hNtrCorrPSSel",Form("Corrected %s multiplicity for PS selected events; %s ; Entries",estimatorName,estimatorName),nMultBins,firstMultBin,lastMultBin);
408  fHistNtrCorrEvSel = new TH1F("hNtrCorrEvSel",Form("Corrected %s multiplicity for selected events; %s ; Entries",estimatorName,estimatorName),nMultBins,firstMultBin,lastMultBin);
409  fHistNtrCorrEvWithCand = new TH1F("hNtrCorrEvWithCand", Form("%s multiplicity for events with D candidates; %s ; Entries",estimatorName,estimatorName),nMultBins,firstMultBin,lastMultBin);// Total multiplicity
410  fHistNtrCorrEvWithD = new TH1F("hNtrCorrEvWithD", Form("%s multiplicity for events with D in mass region ; %s ; Entries",estimatorName,estimatorName),nMultBins,firstMultBin,lastMultBin); //
411 
412  fHistNtrVsZvtx = new TH2F("hNtrVsZvtx",Form("N%s vs VtxZ; VtxZ;N_{%s};",estimatorName,estimatorName),300,-15,15,nMultBins,firstMultBin,lastMultBin); //
413  fHistNtrCorrVsZvtx = new TH2F("hNtrCorrVsZvtx",Form("N%s vs VtxZ; VtxZ;N_{%s};",estimatorName,estimatorName),300,-15,15,nMultBins,firstMultBin,lastMultBin); //
414  fHistNtrVsnTrackEvWithCand = new TH2F("hNtrVsnTrackEvWithCand",Form("N%s vs nTracks; nTracks; N_{%s};",estimatorName,estimatorName),nMultBins,firstMultBin,lastMultBin,nMultBins,firstMultBin,lastMultBin); //
415 
416  TString histoNtrName;
417  TString histoNtrCorrName;
418  TString parNameNtr;
419 
420  histoNtrName = "hNtrVsSphero";
421  histoNtrCorrName = "hNtrCorrVsSphero";
422  parNameNtr = "Sphero";
423 
424  fHistNtrVsSo = new TH2F(histoNtrName.Data(),Form("N_{%s} vs %s; %s; N_{%s};",estimatorName,parNameNtr.Data(),parNameNtr.Data(),estimatorName), 20, 0., 1., nMultBins,firstMultBin,lastMultBin); //
425  fHistNtrCorrVsSo = new TH2F(histoNtrCorrName.Data(),Form("N_{%s} vs %s; %s; N_{%s};",estimatorName,parNameNtr.Data(),parNameNtr.Data(),estimatorName), 20, 0., 1., nMultBins, firstMultBin,lastMultBin); //
426 
427  fHistnTrackvsEtavsPhi = new TH3F("hnTrackvsEtavsPhi", "Eta vs Phi vs nTracks; #eta; #phi[rad]; nTracks;", 100, -1.5, 1.5, 100, 0., 6.28,nMultBins,firstMultBin,lastMultBin);
428  fHistnTrackvsEtavsPhiEvWithCand = new TH3F("hnTrackvsEtavsPhiEvWithCand", "Eta vs Phi vs nTracks; #eta; #phi[rad]; nTracks;", 100, -1.5, 1.5, 100, 0., 6.28,nMultBins,firstMultBin,lastMultBin);
429 
430  fHistTrueSovsMeasSo = new TH3F("hTrueSovsMeasSo", "trueSo vs measSo; S_{o} (true); S_{o} (meas); tracklets;", 100, 0., 1., 100, 0., 1., nMultBins, firstMultBin, lastMultBin);
431  fHistTrueSovsMeasSoEvWithCand = new TH3F("hTrueSovsMeasSoEvWithCand", "trueSo vs measSo; S_{o} (true); S_{o} (meas); tracklets;", 100, 0., 1., 100, 0., 1., nMultBins, firstMultBin, lastMultBin);
432 
433  TString histoNtrSphriName;
434  TString histoNtrCorrSphriName;
435  TString parNameNtrSphri;
436 
437  histoNtrSphriName = "hNtrVsSpheri";
438  histoNtrCorrSphriName = "hNtrCorrVsSpheri";
439  parNameNtrSphri = "Spheri";
440 
442  fHistNtrVsSpheri = new TH2F(histoNtrSphriName.Data(),Form("N_{%s} vs %s; %s; N_{%s};",estimatorName,parNameNtrSphri.Data(),parNameNtrSphri.Data(),estimatorName), 20, 0., 1., nMultBins,firstMultBin,lastMultBin); //
443  fHistNtrCorrVsSpheri = new TH2F(histoNtrCorrSphriName.Data(),Form("N_{%s} vs %s; %s; N_{%s};",estimatorName,parNameNtrSphri.Data(),parNameNtrSphri.Data(),estimatorName), 20, 0., 1., nMultBins, firstMultBin,lastMultBin); //
444  }
445 
446  fHistNtrVsNchMC = new TH2F("hNtrVsNchMC",Form("N%s vs NchMC; Nch;N_{%s};",estimatorName,estimatorName),nMultBins,firstMultBin,lastMultBin,nMultBins,firstMultBin,lastMultBin); //
447  fHistNtrCorrVsNchMC = new TH2F("hNtrCorrVsNchMC",Form("N%s vs Nch; Nch;N_{%s};",estimatorName,estimatorName),nMultBins,firstMultBin,lastMultBin,nMultBins,firstMultBin,lastMultBin); //
448 
449  fHistNtrVsNchMCPrimary = new TH2F("hNtrVsNchMCPrimary",Form("N%s vs Nch (Primary); Nch (Primary);N_{%s};",estimatorName,estimatorName),nMultBins,firstMultBin,lastMultBin,nMultBins,firstMultBin,lastMultBin); //
450  fHistNtrCorrVsNchMCPrimary = new TH2F("hNtrCorrVsNchMCPrimary",Form("N%s vs Nch (Primary); Nch(Primary) ;N_{%s};",estimatorName,estimatorName),nMultBins,firstMultBin,lastMultBin,nMultBins,firstMultBin,lastMultBin); //
451 
452  fHistNtrVsNchMCPhysicalPrimary = new TH2F("hNtrVsNchMCPhysicalPrimary",Form("N%s vs Nch (Physical Primary); Nch (Physical Primary);N_{%s};",estimatorName,estimatorName),nMultBins,firstMultBin,lastMultBin,nMultBins,firstMultBin,lastMultBin); //
453  fHistNtrCorrVsNchMCPhysicalPrimary = new TH2F("hNtrCorrVsMCPhysicalPrimary",Form("N%s vs Nch (Physical Primary); Nch (Physical Primary);N_{%s};",estimatorName,estimatorName),nMultBins,firstMultBin,lastMultBin,nMultBins,firstMultBin,lastMultBin); //
454 
455  fHistGenPrimaryParticlesInelGt0 = new TH1F("hGenPrimaryParticlesInelGt0","Multiplcity of generated charged particles ; Nparticles ; Entries",nMultBins,firstMultBin,lastMultBin);
456 
457  fHistNchMCVsNchMCPrimaryVsNchMCPhysicalPrimary = new TH3F("fHistNchMCVsNchMCPrimaryVsNchMCPhysicalPrimary", "MC: Nch (Physical Primary) vs Nch (Primary) vs Nch (Generated); Nch (Generated); Nch (Primary); Nch (Physical Primary)",nMultBins,firstMultBin,lastMultBin,nMultBins,firstMultBin,lastMultBin,nMultBins,firstMultBin,lastMultBin);
458 
463 
464  fOutput->Add(fHistNtrVsZvtx);
471 
472  fOutput->Add(fHistNtrVsSo);
477  }
478  fOutput->Add(fHistNtrVsNchMC);
486 
487  fHistNEvents = new TH1F("fHistNEvents", "number of events ",11,-0.5,10.5);
488  fHistNEvents->GetXaxis()->SetBinLabel(1,"nEvents total");
489  fHistNEvents->GetXaxis()->SetBinLabel(2,"nEvents with Z vertex");
490  fHistNEvents->GetXaxis()->SetBinLabel(3,"nEvents selected");
491  fHistNEvents->GetXaxis()->SetBinLabel(4,"Rejected due to trigger");
492  fHistNEvents->GetXaxis()->SetBinLabel(5,"Rejected due to phys sel");
493  fHistNEvents->GetXaxis()->SetBinLabel(6,"Rejected due to vertex cuts");
494  fHistNEvents->GetXaxis()->SetBinLabel(7,"Rejected due to pileup");
495  fHistNEvents->GetXaxis()->SetBinLabel(8,"Total no. of candidate");
496  fHistNEvents->GetXaxis()->SetBinLabel(9,"no. of cand wo bitmask");
497  fHistNEvents->GetXaxis()->SetBinLabel(10,"D after cuts (No PID)");
498  fHistNEvents->GetXaxis()->SetBinLabel(11,"D after cuts + PID)");
499  fHistNEvents->GetXaxis()->SetNdivisions(1,kFALSE);
500  fHistNEvents->SetMinimum(0);
501  fOutput->Add(fHistNEvents);
502 
503  // With flag fFillSoSparseChecks to fill THnSparse with MultUncorr and NoPid cases ( 0 = only Mult, 1 = Mult and multUncorr, 2 = NoPid and 3 is All)
504  Int_t nbinsSo[4]={48, fNMassBins, 20, nMultBins};
505  Double_t xminSo[4]={0., fLowmasslimit,0., firstMultBin};
506  Double_t xmaxSo[4]={24., fUpmasslimit, 1., lastMultBin};
507 
508  Int_t nbinsSoSpheri[5]={48, fNMassBins, 20, nMultBins, 20};
509  Double_t xminSoSpheri[5]={0., fLowmasslimit,0., firstMultBin, 0.};
510  Double_t xmaxSoSpheri[5]={24., fUpmasslimit, 1., lastMultBin, 1.};
511 
512  Int_t nbinsSowithMultUncorr[5]={48, fNMassBins, 20, nMultBins, nMultBins};
513  Double_t xminSowithMultUncorr[5]={0., fLowmasslimit,0., firstMultBin, firstMultBin};
514  Double_t xmaxSowithMultUncorr[5]={24., fUpmasslimit, 1., lastMultBin, lastMultBin};
515 
516  Int_t nbinsSoSpheriwithMultUncorr[6]={48, fNMassBins, 20, nMultBins, nMultBins, 20};
517  Double_t xminSoSpheriwithMultUncorr[6]={0., fLowmasslimit,0., firstMultBin, firstMultBin, 0.};
518  Double_t xmaxSoSpheriwithMultUncorr[6]={24., fUpmasslimit, 1., lastMultBin, lastMultBin, 1.};
519 
520  TString histoName = "hSparseEvtShape";
521  TString histoNameNoPid = "hSparseEvtShapewithNoPid";
522  TString parNameSo = "Spherocity";
523  TString parNameSpheri = "Sphericity";
524  TString histoNamePrompt = "hSparseEvtShapePrompt";
525  TString histoNameFeeddown = "hSparseEvtShapeFeeddown";
526  TString histoNameRecSphero = "hSparseEvtShapeRecSphero";
527 
528  if(fFillSoSparseChecks == 1 || fFillSoSparseChecks == 3){
529  if(fCalculateSphericity) fSparseEvtShape = new THnSparseD(histoName.Data(), Form("D candidates:; p_{T} [GeV/c]; InvMass [GeV/c^{2}]; %s; Multipicity; MultipicityUncorr; %s;", parNameSo.Data(), parNameSpheri.Data()), 6 , nbinsSoSpheriwithMultUncorr, xminSoSpheriwithMultUncorr, xmaxSoSpheriwithMultUncorr);
530  else fSparseEvtShape = new THnSparseD(histoName.Data(), Form("D candidates:; p_{T} [GeV/c]; InvMass [GeV/c^{2}]; %s; Multipicity; MultipicityUncorr;", parNameSo.Data()), 5 , nbinsSowithMultUncorr, xminSowithMultUncorr, xmaxSowithMultUncorr);
531  }
532  else{
533  if(fCalculateSphericity) fSparseEvtShape = new THnSparseD(histoName.Data(), Form("D candidates:; p_{T} [GeV/c]; InvMass [GeV/c^{2}]; %s; Multipicity; %s;", parNameSo.Data(), parNameSpheri.Data()), 5 , nbinsSoSpheri, xminSoSpheri, xmaxSoSpheri);
534  else fSparseEvtShape = new THnSparseD(histoName.Data(), Form("D candidates:; p_{T} [GeV/c]; InvMass [GeV/c^{2}]; %s; Multipicity;", parNameSo.Data()), 4 , nbinsSo, xminSo, xmaxSo);
535  }
536 
537  if(fFillSoSparseChecks == 2|| fFillSoSparseChecks == 3) {
538  if(fCalculateSphericity) fSparseEvtShapewithNoPid = new THnSparseD(histoNameNoPid.Data(), Form("D candidates with NoPID:; p_{T} [GeV/c]; InvMass [GeV/c^{2}]; %s; Multipicity; %s;", parNameSo.Data(), parNameSpheri.Data()), 5 , nbinsSoSpheri, xminSoSpheri, xmaxSoSpheri);
539  else fSparseEvtShapewithNoPid = new THnSparseD(histoNameNoPid.Data(), Form("D candidates with NoPID:; p_{T} [GeV/c]; InvMass [GeV/c^{2}]; %s; Multipicity;", parNameSo.Data()), 4 , nbinsSo, xminSo, xmaxSo);
540  }
541 
542  if(fRecomputeSpherocity) fSparseEvtShapeRecSphero = new THnSparseD(histoNameRecSphero.Data(), Form("D candidates:; p_{T} [GeV/c]; InvMass [GeV/c^{2}]; %s; Multipicity; RecSpherocity;", parNameSo.Data()), 5 , nbinsSoSpheri, xminSoSpheri, xmaxSoSpheri);
543 
544  fOutput->Add(fSparseEvtShape);
547 
548  Int_t nbinsPrompt[4]={48, nMultBins, 20, 100};
549  Int_t nbinsFeeddown[4]={48, nMultBins, 20, 100};
550  Double_t xminPrompt[4] = {0.,firstMultBin, 0., -1.};
551  Double_t xmaxPrompt[4] = {24.,lastMultBin, 1., 1.};
552  Double_t xminFeeddown[4] = {0.,firstMultBin, 0., -1.};
553  Double_t xmaxFeeddown[4] = {24.,lastMultBin, 1., 1.};
554 
555  Int_t nbinsRecSpheroPrompt[5]={48, nMultBins, 20, 100, 20};
556  Int_t nbinsRecSpheroFeeddown[5]={48, nMultBins, 20, 100, 20};
557  Double_t xminRecSpheroPrompt[5] = {0.,firstMultBin, 0., -1., 0.};
558  Double_t xmaxRecSpheroPrompt[5] = {24.,lastMultBin, 1., 1., 1.};
559  Double_t xminRecSpheroFeeddown[5] = {0.,firstMultBin, 0., -1., 0.};
560  Double_t xmaxRecSpheroFeeddown[5] = {24.,lastMultBin, 1., 1., 1.};
561 
562  if(fReadMC){
564  fSparseEvtShapePrompt = new THnSparseD(histoNamePrompt.Data(), Form("D candidates:; p_{T} [GeV/c]; InvMass [GeV/c^{2}]; %s; Multipicity; RecSpherocity;", parNameSo.Data()), 5 , nbinsSoSpheri, xminSoSpheri, xmaxSoSpheri);
565  fSparseEvtShapeFeeddown = new THnSparseD(histoNameFeeddown.Data(), Form("D candidates:; p_{T} [GeV/c]; InvMass [GeV/c^{2}]; %s; Multipicity; RecSpherocity;", parNameSo.Data()), 5 , nbinsSoSpheri, xminSoSpheri, xmaxSoSpheri);
566  }
567  else{
568  fSparseEvtShapePrompt = new THnSparseD(histoNamePrompt.Data(), Form("D candidates:; p_{T} [GeV/c]; InvMass [GeV/c^{2}]; %s; Multipicity;", parNameSo.Data()), 4 , nbinsSo, xminSo, xmaxSo);
569  fSparseEvtShapeFeeddown = new THnSparseD(histoNameFeeddown.Data(), Form("D candidates:; p_{T} [GeV/c]; InvMass [GeV/c^{2}]; %s; Multipicity;", parNameSo.Data()), 4 , nbinsSo, xminSo, xmaxSo);
570  }
571 
574 
575  //Prompt
577  fMCAccGenPrompt = new THnSparseD("hMCAccGenPrompt", "kStepMCAcceptance:; p_{T} [GeV/c]; Multipicity; Spherocity; y; RecSpherocity; - promptD",5,nbinsRecSpheroPrompt,xminRecSpheroPrompt,xmaxRecSpheroPrompt);
578  fMCRecoPrompt = new THnSparseD("hMCRecoPrompt","kStepRecoPID:; p_{T} [GeV/c]; Multipicity; Spherocity; y; RecSpherocity; - promptD",5,nbinsRecSpheroPrompt,xminRecSpheroPrompt,xmaxRecSpheroPrompt);
579  }
580  else{
581  fMCAccGenPrompt = new THnSparseD("hMCAccGenPrompt","kStepMCAcceptance:; p_{T} [GeV/c]; Multipicity; Spherocity; y; - promptD",4,nbinsPrompt,xminPrompt,xmaxPrompt);
582  fMCRecoPrompt = new THnSparseD("hMCRecoPrompt","kStepRecoPID:; p_{T} [GeV/c]; Multipicity; Spherocity; y; - promptD",4,nbinsPrompt,xminPrompt,xmaxPrompt);
583  }
585  fMCAccGenPromptSpheri = new THnSparseD("hMCAccGenPromptSpheri","kStepMCAcceptance:; p_{T} [GeV/c]; Multipicity; Sphericity; y; - promptD",4,nbinsPrompt,xminPrompt,xmaxPrompt);
586  fMCRecoPromptSpheri = new THnSparseD("hMCRecoPromptSpheri","kStepRecoPID:; p_{T} [GeV/c]; Multipicity; Sphericity; y; - promptD",4,nbinsPrompt,xminPrompt,xmaxPrompt);
587  }
588  fMCAccGenPromptEvSel = new THnSparseD("hMCAccGenPromptEvSel","kStepMCAcceptanceEvSel:; p_{T} [GeV/c]; Multipicity; Spherocity; y; - promptD",4,nbinsPrompt,xminPrompt,xmaxPrompt);
589 
590  //Feeddown
592  fMCAccGenFeeddown = new THnSparseD("hMCAccGenBFeeddown","kStepMCAcceptance:; p_{T} [GeV/c]; Multipicity; Spherocity; y; RecSpherocity; - DfromB",5,nbinsRecSpheroFeeddown,xminRecSpheroFeeddown,xmaxRecSpheroFeeddown);
593  fMCRecoFeeddown = new THnSparseD("hMCRecoFeeddown","kStepRecoPID:; p_{T} [GeV/c]; Multipicity; Spherocity; y; RecSpherocity; - DfromB",5,nbinsRecSpheroFeeddown,xminRecSpheroFeeddown,xmaxRecSpheroFeeddown);
594  }
595  else{
596  fMCAccGenFeeddown = new THnSparseD("hMCAccGenBFeeddown","kStepMCAcceptance:; p_{T} [GeV/c]; Multipicity; Spherocity; y; - DfromB",4,nbinsFeeddown,xminFeeddown,xmaxFeeddown);
597  fMCRecoFeeddown = new THnSparseD("hMCRecoFeeddown","kStepRecoPID:; p_{T} [GeV/c]; Multipicity; Spherocity; y; - DfromB",4,nbinsFeeddown,xminFeeddown,xmaxFeeddown);
598  }
600  fMCAccGenFeeddownSpheri = new THnSparseD("hMCAccGenBFeeddownSpheri","kStepMCAcceptance:; p_{T} [GeV/c]; Multipicity; Sphericity; y; - DfromB",4,nbinsFeeddown,xminFeeddown,xmaxFeeddown);
601  fMCRecoFeeddownSpheri = new THnSparseD("hMCRecoFeeddownSpheri","kStepRecoPID:; p_{T} [GeV/c]; Multipicity; Sphericity; y; - DfromB",4,nbinsFeeddown,xminFeeddown,xmaxFeeddown);
602  }
603  fMCAccGenFeeddownEvSel = new THnSparseD("hMCAccGenBFeeddownEvSel","kStepMCAcceptance:; p_{T} [GeV/c]; Multipicity; Spherocity; y; - DfromB",4,nbinsFeeddown,xminFeeddown,xmaxFeeddown);
604 
605  //BothPromptFeeddown
606  fMCRecoBothPromptFD = new THnSparseD("hMCRecoBothPromptFD","kStepRecoPID:; p_{T} [GeV/c]; Multipicity; Spherocity; y; - BothPromptFD",4,nbinsPrompt,xminPrompt,xmaxPrompt);
607 
608  if(fCalculateSphericity) fMCRecoBothPromptFDSpheri = new THnSparseD("hMCRecoBothPromptFDSpheri","kStepRecoPID:; p_{T} [GeV/c]; Multipicity; Sphericity; y; - BothPromptFD",4,nbinsPrompt,xminPrompt,xmaxPrompt);
609 
621  }
624  }
625 
627 
628  fCounterC = new AliNormalizationCounter("NormCounterCorrMult");
629  fCounterC->SetStudyMultiplicity(kTRUE,1.);
630  fCounterC->SetStudySpherocity(kTRUE,10.);
631  fCounterC->Init();
632 
633  fCounterU = new AliNormalizationCounter("NormCounterUnCorrMult");
634  fCounterU->SetStudyMultiplicity(kTRUE,1.);
635  fCounterU->SetStudySpherocity(kTRUE,10.);
636  fCounterU->Init();
637 
638  fCounterCandidates = new AliNormalizationCounter("NormCounterCorrMultCandidates");
642 
643  fOutputCounters = new TList();
644  fOutputCounters->SetOwner();
645  fOutputCounters->SetName("OutputCounters");
649 
650  PostData(1,fOutput);
651  PostData(2,fListCuts);
652  PostData(3,fOutputCounters);
653  PostData(4,fListProfiles);
654  if(fReadMC) PostData(5,fOutputEffCorr);
655 
656  return;
657 }
658 
659 //________________________________________________________________________
661 {
662  // Execute analysis for current event:
663  // heavy flavor candidates association to MC truth
664 
665  AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
666 
667  TClonesArray *arrayCand = 0;
668  TString arrayName="";
669  UInt_t pdgDau[3];
670  Int_t nDau=0;
671  Int_t selbit=0;
672  if(fPdgMeson==411){
673  arrayName="Charm3Prong";
674  pdgDau[0]=211; pdgDau[1]=321; pdgDau[2]=211;
675  nDau=3;
677  }else if(fPdgMeson==421){
678  arrayName="D0toKpi";
679  pdgDau[0]=211; pdgDau[1]=321; pdgDau[2]=0;
680  nDau=2;
682  }else if(fPdgMeson==413){
683  arrayName="Dstar";
684  pdgDau[0]=321; pdgDau[1]=211; pdgDau[2]=0; // Quoting here D0 daughters (D* ones on another variable later)
685  nDau=2;
687  }
688 
689  if(!aod && AODEvent() && IsStandardAOD()) {
690  // In case there is an AOD handler writing a standard AOD, use the AOD
691  // event in memory rather than the input (ESD) event.
692  aod = dynamic_cast<AliAODEvent*> (AODEvent());
693  // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
694  // have to taken from the AOD event hold by the AliAODExtension
695  AliAODHandler* aodHandler = (AliAODHandler*)
696  ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
697  if(aodHandler->GetExtensions()) {
698  AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
699  AliAODEvent *aodFromExt = ext->GetAOD();
700  arrayCand=(TClonesArray*)aodFromExt->GetList()->FindObject(arrayName.Data());
701  }
702  } else if(aod) {
703  arrayCand=(TClonesArray*)aod->GetList()->FindObject(arrayName.Data());
704  }
705 
706  if(!aod || !arrayCand) {
707  printf("AliAnalysisTaskSEDvsEventShapes::UserExec: Charm3Prong branch not found!\n");
708  return;
709  }
710 
711  if(fisPPbData && fReadMC){
712  Int_t runnumber = aod->GetRunNumber();
713  if(aod->GetTriggerMask()==0 &&
714  (runnumber>=195344 && runnumber<=195677)){
715  AliDebug(3,"Event rejected because of null trigger mask");
716  return;
717  }
718  }
719 
720  // fix for temporary bug in ESDfilter
721  // the AODs with null vertex pointer didn't pass the PhysSel
722  if(!aod->GetPrimaryVertex()||TMath::Abs(aod->GetMagneticField())<0.001) return;
723 
724  Int_t countTreta1=0, countTreta03=0, countTreta05=0, countTreta16=0;
725  AliAODTracklets* tracklets=aod->GetTracklets();
726  Int_t nTr=tracklets->GetNumberOfTracklets();
727  for(Int_t iTr=0; iTr<nTr; iTr++){
728  Double_t theta=tracklets->GetTheta(iTr);
729  Double_t eta=-TMath::Log(TMath::Tan(theta/2.));
730  if(eta>-0.3 && eta<0.3) countTreta03++;
731  if(eta>-0.5 && eta<0.5) countTreta05++;
732  if(eta>-1.0 && eta<1.0) countTreta1++;
733  if(eta>-1.6 && eta<1.6) countTreta16++;
734  }
735 
736  Int_t vzeroMult=0, vzeroMultA=0, vzeroMultC=0;
737  Int_t vzeroMultEq=0, vzeroMultAEq=0, vzeroMultCEq=0;
738  AliAODVZERO *vzeroAOD = (AliAODVZERO*)aod->GetVZEROData();
739  if(vzeroAOD) {
740  vzeroMultA = static_cast<Int_t>(vzeroAOD->GetMTotV0A());
741  vzeroMultC = static_cast<Int_t>(vzeroAOD->GetMTotV0C());
742  vzeroMult = vzeroMultA + vzeroMultC;
743  vzeroMultAEq = static_cast<Int_t>(AliVertexingHFUtils::GetVZEROAEqualizedMultiplicity(aod));
744  vzeroMultCEq = static_cast<Int_t>(AliVertexingHFUtils::GetVZEROCEqualizedMultiplicity(aod));
745  vzeroMultEq = vzeroMultAEq + vzeroMultCEq;
746  }
747 
748  Int_t countMult = countTreta1;
749  if(fMultiplicityEstimator==kNtrk03) { countMult = countTreta03; }
750  else if(fMultiplicityEstimator==kNtrk05) { countMult = countTreta05; }
751  else if(fMultiplicityEstimator==kNtrk10to16) { countMult = countTreta16 - countTreta1; }
752  else if(fMultiplicityEstimator==kVZERO) { countMult = vzeroMult; }
753  else if(fMultiplicityEstimator==kVZEROA) { countMult = vzeroMultA; }
754  else if(fMultiplicityEstimator==kVZEROEq) { countMult = vzeroMultEq; }
755  else if(fMultiplicityEstimator==kVZEROAEq) { countMult = vzeroMultAEq; }
756 
757  Double_t spherocity = -0.5;
758  Double_t sphericity = -0.5;
759  if(fCalculateSphericity){ //When kTRUE, it calculates Sphericity and THnSparse filled for sphericity
761  }
763 
764  Double_t St=1;
765  fCounterU->StoreEvent(aod,fRDCutsAnalysis,fReadMC,countMult,spherocity);
766  fHistNEvents->Fill(0); // count event
767 
768  Double_t countTreta1corr=countTreta1;
769  Double_t countCorr=countMult;
770  AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
771  // In case of VZERO multiplicity, consider the zvtx correction flag
772  // fDoVZER0ParamVertexCorr: 0= none, 1= usual d2h, 2=AliESDUtils
773  Bool_t isDataDrivenZvtxCorr=kTRUE;
774  Bool_t isVtxOk=kFALSE;
775  Int_t vzeroMultACorr=vzeroMultA, vzeroMultCCorr=vzeroMultC, vzeroMultCorr=vzeroMult;
776  Int_t vzeroMultAEqCorr=vzeroMultAEq, vzeroMultCEqCorr=vzeroMultCEq, vzeroMultEqCorr=vzeroMultEq;
777  if(vtx1){
778  if(vtx1->GetNContributors()>0){
779  fHistNEvents->Fill(1);
780  isVtxOk=kTRUE;
781  }
782  }
783  if(isVtxOk){
787  // do not correct
788  isDataDrivenZvtxCorr=kFALSE;
789  } else if (fDoVZER0ParamVertexCorr==2){
790  // use AliESDUtils correction
791  Float_t zvtx = vtx1->GetZ();
792  isDataDrivenZvtxCorr=kFALSE;
793  vzeroMultACorr = static_cast<Int_t>(AliESDUtils::GetCorrV0A(vzeroMultA,zvtx));
794  vzeroMultCCorr = static_cast<Int_t>(AliESDUtils::GetCorrV0C(vzeroMultC,zvtx));
795  vzeroMultCorr = vzeroMultACorr + vzeroMultCCorr;
796  vzeroMultAEqCorr = static_cast<Int_t>(AliESDUtils::GetCorrV0A(vzeroMultAEq,zvtx));
797  vzeroMultCEqCorr =static_cast<Int_t>( AliESDUtils::GetCorrV0C(vzeroMultCEq,zvtx));
798  vzeroMultEqCorr = vzeroMultAEqCorr + vzeroMultCEqCorr;
799  if(fMultiplicityEstimator==kVZERO) { countCorr = vzeroMultCorr; }
800  else if(fMultiplicityEstimator==kVZEROA) { countCorr = vzeroMultACorr; }
801  else if(fMultiplicityEstimator==kVZEROEq) { countCorr = vzeroMultEqCorr; }
802  else if(fMultiplicityEstimator==kVZEROAEq) { countCorr = vzeroMultAEqCorr; }
803  }
804  }
805  }
806  // Data driven multiplicity z-vertex correction
807  if(isVtxOk && isDataDrivenZvtxCorr){
808  TProfile* estimatorAvg = GetEstimatorHistogram(aod);
809  if(estimatorAvg){
810  countTreta1corr=static_cast<Int_t>(AliVertexingHFUtils::GetCorrectedNtracklets(estimatorAvg,countTreta1,vtx1->GetZ(),fRefMult));
811  countCorr=static_cast<Int_t>(AliVertexingHFUtils::GetCorrectedNtracklets(estimatorAvg,countMult,vtx1->GetZ(),fRefMult));
812  }
813  }
814 
815  fCounterC->StoreEvent(aod,fRDCutsAnalysis,fReadMC,countCorr,spherocity);
816 
817  Bool_t isEvSel=fRDCutsAnalysis->IsEventSelected(aod);
818 
819  if(fRDCutsAnalysis->GetWhyRejection()==5) fHistNEvents->Fill(3);
820  if(fRDCutsAnalysis->GetWhyRejection()==7) fHistNEvents->Fill(4);
821  if(fRDCutsAnalysis->GetWhyRejection()==6) fHistNEvents->Fill(5);
822  if(fRDCutsAnalysis->GetWhyRejection()==1) fHistNEvents->Fill(6);
823 
825 
826  if(!isEvPSRejected){
827  fHistNtrCorrPSSel->Fill(countCorr);
828  }
829 
830  TClonesArray *arrayMC=0;
831  AliAODMCHeader *mcHeader=0;
832 
833  fWeight=1.;
834  Double_t nchWeight=1.0;
835 
836  Int_t nChargedMCEta10=0, nChargedMCEta03=0, nChargedMCEta05=0, nChargedMCEta16=0, nChargedMCEtam37tm17=0, nChargedMCEta28t51=0;
837  Int_t nChargedMCPrimaryEta10=0, nChargedMCPrimaryEta03=0, nChargedMCPrimaryEta05=0, nChargedMCPrimaryEta16=0, nChargedMCPrimaryEtam37tm17=0, nChargedMCPrimaryEta28t51=0;
838  Int_t nChargedMCPhysicalPrimaryEta10=0, nChargedMCPhysicalPrimaryEta03=0, nChargedMCPhysicalPrimaryEta05=0, nChargedMCPhysicalPrimaryEta16=0, nChargedMCPhysicalPrimaryEtam37tm17=0, nChargedMCPhysicalPrimaryEta28t51=0;
839  Int_t nChargedMC=0, nChargedMCPrimary=0, nChargedMCPhysicalPrimary=0;
840 
841  Double_t genspherocity = -0.5;
842  // load MC particles and get weight on Nch
843  if(fReadMC){
844 
845  arrayMC = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
846  if(!arrayMC) {
847  printf("AliAnalysisTaskSEDvsEventShapes::UserExec: MC particles branch not found!\n");
848  return;
849  }
850  // load MC header
851  mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
852  if(!mcHeader) {
853  printf("AliAnalysisTaskSEDvsEventShapes::UserExec: MC header branch not found!\n");
854  return;
855  }
856 
857  for(Int_t i=0; i<arrayMC->GetEntriesFast(); i++){
858  AliAODMCParticle *part=(AliAODMCParticle*)arrayMC->UncheckedAt(i);
859  Int_t charge = part->Charge();
860  Double_t eta = part->Eta();
861  Bool_t isPrim = part->IsPrimary();
862  Bool_t isPhysPrim = part->IsPhysicalPrimary();
863  if(charge!=0) {
864  if(eta>-0.3 && eta< 0.3) {
865  nChargedMCEta03++;
866  if(isPrim) nChargedMCPrimaryEta03++;
867  if(isPhysPrim) nChargedMCPhysicalPrimaryEta03++;
868  }
869  if(eta>-0.5 && eta< 0.5) {
870  nChargedMCEta05++;
871  if(isPrim) nChargedMCPrimaryEta05++;
872  if(isPhysPrim) nChargedMCPhysicalPrimaryEta05++;
873  }
874  if(eta>-1.0 && eta< 1.0) {
875  nChargedMCEta10++;
876  if(isPrim) nChargedMCPrimaryEta10++;
877  if(isPhysPrim) nChargedMCPhysicalPrimaryEta10++;
878  }
879  if(eta>-1.6 && eta< 1.6) {
880  nChargedMCEta16++;
881  if(isPrim) nChargedMCPrimaryEta16++;
882  if(isPhysPrim) nChargedMCPhysicalPrimaryEta16++;
883  }
884  if(eta>-3.7 && eta<-1.7) {
885  nChargedMCEtam37tm17++;
886  if(isPrim) nChargedMCPrimaryEtam37tm17++;
887  if(isPhysPrim) nChargedMCPhysicalPrimaryEtam37tm17++;
888  }
889  if(eta> 2.8 && eta< 5.1) {
890  nChargedMCEta28t51++;
891  if(isPrim) nChargedMCPrimaryEta28t51++;
892  if(isPhysPrim) nChargedMCPhysicalPrimaryEta28t51++;
893  }
894  }
895  }
896 
897  nChargedMC=nChargedMCEta10;
898  nChargedMCPrimary=nChargedMCPrimaryEta10;
899  nChargedMCPhysicalPrimary=nChargedMCPhysicalPrimaryEta10;
900 
901  // Compute the Nch weights (reference is Ntracklets within |eta|<1.0)
902  if(fUseNchWeight>0){
903 
904  Double_t tmpweight = 1.0;
905  Double_t tmpXweight=nChargedMCPhysicalPrimary; // Nch weights
906  if(fUseNchWeight==2) tmpXweight=countMult; // Ntrk weights
907 
908  if(tmpXweight<=0) tmpweight = 0.0;
909  else{
910  Double_t pMeas = fHistoMeasNch->GetBinContent(fHistoMeasNch->FindBin(tmpXweight));
911  //printf(" pMeas=%2.2f and histo MCNch %s \n",pMeas,fHistoMCNch);
912  Double_t pMC = fHistoMCNch->GetBinContent(fHistoMCNch->FindBin(tmpXweight));
913  tmpweight = pMC>0 ? pMeas/pMC : 0.;
914  }
915  nchWeight *= tmpweight;
916  AliDebug(2,Form("Using Nch weights, Mult=%f Weight=%f\n",tmpXweight,nchWeight));
917  }
918 
919  FillMCGenAccHistos(aod, arrayMC, mcHeader, countCorr, spherocity, sphericity, isEvSel, nchWeight);//Fill 2 separate THnSparses, one for prompt andf one for feeddown
921 
922  }
923 
924  if(!isEvSel)return;
925  fHistNEvents->Fill(2);
926 
927  if(vtx1){
928  fHistNtrVsZvtx->Fill(vtx1->GetZ(),countMult);
929  fHistNtrCorrVsZvtx->Fill(vtx1->GetZ(),countCorr);
930  fHistNtrVsSo->Fill(spherocity,countMult);
931  fHistNtrCorrVsSo->Fill(spherocity,countCorr);
933  fHistNtrVsSpheri->Fill(sphericity,countMult);
934  fHistNtrCorrVsSpheri->Fill(sphericity,countCorr);
935  }
936  }
937 
938  if(fReadMC){
939  // Now recompute the variables in case another MC estimator is considered
941  nChargedMC = nChargedMCEta16 - nChargedMCEta10;
942  nChargedMCPrimary = nChargedMCPrimaryEta16 - nChargedMCPrimaryEta10;
943  nChargedMCPhysicalPrimary = nChargedMCPhysicalPrimaryEta16 - nChargedMCPhysicalPrimaryEta10;
944  } else if(fMCPrimariesEstimator==kEta05){
945  nChargedMC = nChargedMCEta05;
946  nChargedMCPrimary = nChargedMCPrimaryEta05;
947  nChargedMCPhysicalPrimary = nChargedMCPhysicalPrimaryEta05;
948  } else if(fMCPrimariesEstimator==kEta03){
949  nChargedMC = nChargedMCEta03;
950  nChargedMCPrimary = nChargedMCPrimaryEta03;
951  nChargedMCPhysicalPrimary = nChargedMCPhysicalPrimaryEta03;
952  } else if(fMCPrimariesEstimator==kEtaVZERO){
953  nChargedMC = nChargedMCEtam37tm17 + nChargedMCEta28t51;
954  nChargedMCPrimary = nChargedMCPrimaryEtam37tm17 + nChargedMCPrimaryEta28t51;
955  nChargedMCPhysicalPrimary = nChargedMCPhysicalPrimaryEtam37tm17 + nChargedMCPhysicalPrimaryEta28t51;
956  } else if(fMCPrimariesEstimator==kEtaVZEROA){
957  nChargedMC = nChargedMCEta28t51;
958  nChargedMCPrimary = nChargedMCPrimaryEta28t51;
959  nChargedMCPhysicalPrimary = nChargedMCPhysicalPrimaryEta28t51;
960  }
961  // Here fill the MC correlation plots
962  if(nChargedMCPhysicalPrimary>0){ // INEL>0 for |eta|<1
963  fHistGenPrimaryParticlesInelGt0->Fill(nChargedMCPhysicalPrimary,nchWeight);
964  }
965 
966  fHistNtrVsNchMC->Fill(nChargedMC,countMult,nchWeight);
967  fHistNtrCorrVsNchMC->Fill(nChargedMC,countCorr,nchWeight);
968 
969  fHistNtrVsNchMCPrimary->Fill(nChargedMCPrimary,countMult,nchWeight);
970  fHistNtrCorrVsNchMCPrimary->Fill(nChargedMCPrimary,countCorr,nchWeight);
971 
972  fHistNtrVsNchMCPhysicalPrimary->Fill(nChargedMCPhysicalPrimary,countMult,nchWeight);
973  fHistNtrCorrVsNchMCPhysicalPrimary->Fill(nChargedMCPhysicalPrimary,countCorr,nchWeight);
974 
975  fHistNchMCVsNchMCPrimaryVsNchMCPhysicalPrimary->Fill(nChargedMC,nChargedMCPrimary,nChargedMCPhysicalPrimary,nchWeight);
976  }
977 
978  Int_t nCand = arrayCand->GetEntriesFast();
979  Int_t nSelectedNoPID=0,nSelectedPID=0,nSelectedInMassPeak=0;
980  Double_t mD0PDG = TDatabasePDG::Instance()->GetParticle(421)->Mass();
981  Double_t mDplusPDG = TDatabasePDG::Instance()->GetParticle(411)->Mass();
982  Double_t mDstarPDG = TDatabasePDG::Instance()->GetParticle(413)->Mass();
983 
984  // pdg of daughters needed for D* too
985  UInt_t pdgDgDStartoD0pi[2]={421,211};
986 
987  Double_t aveMult=0.;
988  Double_t nSelCand=0.;
989  for (Int_t iCand = 0; iCand < nCand; iCand++) {
990  AliAODRecoDecayHF *d = (AliAODRecoDecayHF*)arrayCand->UncheckedAt(iCand);
991  AliAODRecoCascadeHF *dCascade = NULL;
992  if(fPdgMeson==413) dCascade = (AliAODRecoCascadeHF*)d;
993 
994  fHistNEvents->Fill(7);
995  if(fUseBit && !d->HasSelectionBit(selbit)){
996  fHistNEvents->Fill(8);
997  continue;
998  }
999 
1000  Double_t ptCand = d->Pt();
1001  Double_t rapid=d->Y(fPdgMeson);
1002  Bool_t isFidAcc=fRDCutsAnalysis->IsInFiducialAcceptance(ptCand,rapid);
1003  if(!isFidAcc) continue;
1004 
1005  Int_t passAllCuts=fRDCutsAnalysis->IsSelected(d,AliRDHFCuts::kAll,aod);
1006  Int_t passTopolCuts=fRDCutsAnalysis->GetIsSelectedCuts();
1007  if(passTopolCuts==0) continue;
1008  nSelectedNoPID++;
1009  fHistNEvents->Fill(9);
1010  if(passAllCuts){
1011  nSelectedPID++;
1012  fHistNEvents->Fill(10);
1013  }
1014  Double_t multForCand = countCorr;
1015 
1016  //Weight according to ptcand, using FONLL
1017  Double_t ptWeight = 1.;
1018  if(fUsePtWeight) ptWeight = GetPtWeight(ptCand);
1019  fWeight = nchWeight*ptWeight;
1020 
1022  // For the D* case, subtract only the D0 daughter tracks <=== FIXME !!
1023  AliAODRecoDecayHF2Prong* d0fromDstar = NULL;
1024  if(fPdgMeson==413) d0fromDstar = (AliAODRecoDecayHF2Prong*)dCascade->Get2Prong();
1025 
1026  for(Int_t iDau=0; iDau<nDau; iDau++){
1027  AliAODTrack *t = NULL;
1028  if(fPdgMeson==413){ t = (AliAODTrack*)d0fromDstar->GetDaughter(iDau); }
1029  else{ t = (AliAODTrack*)d->GetDaughter(iDau); }
1030  if(!t) continue;
1031  if(t->HasPointOnITSLayer(0) && t->HasPointOnITSLayer(1)){
1032  if(multForCand>0) multForCand-=1;
1033  }
1034  }
1035  }
1036 
1037  Bool_t isPrimary=kTRUE;
1038  Double_t trueImpParXY=9999.;
1039  Double_t impparXY=d->ImpParXY()*10000.;
1040  Double_t dlen=0.1; //FIXME
1041  Double_t mass[2];
1042  if(fPdgMeson==411){
1043  mass[0]=d->InvMass(nDau,pdgDau);
1044  mass[1]=-1.;
1045  if(TMath::Abs(mass[0]-mDplusPDG)<0.02) nSelectedInMassPeak++; //20 MeV for now... FIXME
1046  }else if(fPdgMeson==421){
1047  UInt_t pdgdaughtersD0[2]={211,321};//pi,K
1048  UInt_t pdgdaughtersD0bar[2]={321,211};//K,pi
1049  mass[0]=d->InvMass(2,pdgdaughtersD0);
1050  mass[1]=d->InvMass(2,pdgdaughtersD0bar);
1051  if(TMath::Abs(mass[0]-mD0PDG)<0.02 || TMath::Abs(mass[1]-mD0PDG)<0.02 ) nSelectedInMassPeak++; //20 MeV for now... FIXME
1052  }else if(fPdgMeson==413){
1053  // FIXME
1054  mass[0]=dCascade->DeltaInvMass();
1055  mass[1]=-1.;
1056  if(TMath::Abs(mass[0]-(mDstarPDG-mD0PDG))<0.0015) nSelectedInMassPeak++; //1 MeV for now... FIXME
1057  }
1058 
1059  Int_t labD0=-1;
1060  Double_t recSpherocity = -0.5;
1061 
1062  // subtract D-meson daughters from spherocity calculation !!
1064  Int_t totTrkToSkip = d->GetNDaughters();
1065  const Int_t nTrkToSkip = totTrkToSkip;
1066  Int_t idToSkip[nTrkToSkip];
1067  for(Int_t i=0; i<nTrkToSkip; i++) idToSkip[i]=-1;
1068 
1069  for(Int_t iDau=0; iDau<nTrkToSkip; iDau++){
1070  AliAODTrack *t = NULL;
1071  t = dynamic_cast<AliAODTrack*>(d->GetDaughter(iDau));
1072  if(!t) continue;
1073  idToSkip[iDau] = t->GetID();
1074  }
1076  }
1077 
1078  // remove D0 from Dstar at reconstruction !!
1079  if(fReadMC && fRemoveD0fromDstar){
1080  if(fPdgMeson==421){
1081  labD0 = d->MatchToMC(fPdgMeson,arrayMC,nDau,(Int_t*)pdgDau);
1082  if(labD0>=0){
1083  Bool_t keep=kTRUE;
1084  AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(arrayMC->At(labD0));
1085  Int_t motherD0 = mcMoth->GetMother();
1086  AliAODMCParticle* mcMothD0 = dynamic_cast<AliAODMCParticle*>(arrayMC->At(motherD0));
1087  if(!mcMothD0) continue;
1088  if(TMath::Abs(mcMothD0->GetPdgCode())==413) keep=kFALSE;
1089  if(!keep) continue;
1090  }
1091  }
1092  }
1093 
1094  Int_t labD=-1;
1095  Int_t Origin = 0;
1096 
1097  for(Int_t iHyp=0; iHyp<2; iHyp++){
1098  if(mass[iHyp]<0.) continue; // for D+ and D* we have 1 mass hypothesis
1099  Double_t invMass=mass[iHyp];
1100  Double_t arrayForSparse[5]={invMass,ptCand,impparXY,dlen,multForCand};
1101 
1102  if(fReadMC){
1103  if(fPdgMeson==413){
1104  labD = dCascade->MatchToMC(fPdgMeson,421,(Int_t*)pdgDgDStartoD0pi,(Int_t*)pdgDau,arrayMC);
1105  } else {
1106  labD = d->MatchToMC(fPdgMeson,arrayMC,nDau,(Int_t*)pdgDau);
1107  }
1108 
1109  Bool_t fillHisto=fDoImpPar;
1110  if(labD>=0){
1111  AliAODMCParticle *partD = (AliAODMCParticle*)arrayMC->At(labD);
1112  Int_t code=partD->GetPdgCode();
1113  Origin = AliVertexingHFUtils::CheckOrigin(arrayMC,partD,fUseQuarkTag);
1114  if(Origin==5) isPrimary=kFALSE;
1115  if(code<0 && iHyp==0) fillHisto=kFALSE;
1116  if(code>0 && iHyp==1) fillHisto=kFALSE;
1117  if(!isPrimary){
1118  if(fPdgMeson==411){
1119  trueImpParXY=AliVertexingHFUtils::GetTrueImpactParameterDplus(mcHeader,arrayMC,partD)*10000.;
1120  }else if(fPdgMeson==421){
1121  trueImpParXY=AliVertexingHFUtils::GetTrueImpactParameterDzero(mcHeader,arrayMC,partD)*10000.;
1122  }else if(fPdgMeson==413){
1123  trueImpParXY=0.;
1124  }
1125  Double_t arrayForSparseTrue[5]={invMass,ptCand,trueImpParXY,dlen,multForCand};
1126  if(fillHisto && passAllCuts){
1127  fHistMassPtImpPar[2]->Fill(arrayForSparse, fWeight);
1128  fHistMassPtImpPar[3]->Fill(arrayForSparseTrue, fWeight);
1129  }
1130  }else{
1131  if(fillHisto && passAllCuts) fHistMassPtImpPar[1]->Fill(arrayForSparse, fWeight);
1132  }
1133  }else{
1134  if(fillHisto && passAllCuts)fHistMassPtImpPar[4]->Fill(arrayForSparse, fWeight);
1135  }
1136  if(TMath::Abs(labD)==fPdgMeson && fMCOption==2) continue;
1137  if(TMath::Abs(labD)!=fPdgMeson && fMCOption==1) continue;
1138  }
1139  if(fPdgMeson==421){
1140  if(iHyp==0 && !(passTopolCuts&1)) continue; // candidate not passing as D0
1141  if(iHyp==1 && !(passTopolCuts&2)) continue; // candidate not passing as D0bar
1142  }
1143  if(fFillSoSparseChecks == 2 || fFillSoSparseChecks == 3){ //Filling THnSparse for Spherocity without PID
1145  Double_t arrayForSparseSoNoPid[5]={ptCand, invMass, spherocity, multForCand, sphericity};
1146  fSparseEvtShapewithNoPid->Fill(arrayForSparseSoNoPid, fWeight);
1147  }else{
1148  Double_t arrayForSparseSoNoPid[4]={ptCand, invMass, spherocity, multForCand};
1149  fSparseEvtShapewithNoPid->Fill(arrayForSparseSoNoPid, fWeight);
1150  }
1151  }
1152  if(fPdgMeson==421){
1153  if(iHyp==0 && !(passAllCuts&1)) continue; // candidate not passing as D0
1154  if(iHyp==1 && !(passAllCuts&2)) continue; // candidate not passing as D0bar
1155  }
1156  if(passAllCuts){
1157  aveMult+=multForCand;
1158  nSelCand+=1.;
1159 
1160  if(fFillSoSparseChecks == 1 || fFillSoSparseChecks == 3){
1162  Double_t arrayForSparseSowithMultUnncorr[6]={ptCand, invMass, spherocity, multForCand, (Double_t)countTreta1, sphericity};
1163  fSparseEvtShape->Fill(arrayForSparseSowithMultUnncorr, fWeight);
1164  }else{
1165  Double_t arrayForSparseSowithMultUnncorr[5]={ptCand, invMass, spherocity, multForCand, (Double_t)countTreta1};
1166  fSparseEvtShape->Fill(arrayForSparseSowithMultUnncorr, fWeight);
1167  }
1168  }
1169  else{
1171  Double_t arrayForSparseSo[5]={ptCand, invMass, spherocity, multForCand, sphericity};
1172  fSparseEvtShape->Fill(arrayForSparseSo, fWeight);
1173  }else{
1174  Double_t arrayForSparseSo[4]={ptCand, invMass, spherocity, multForCand};
1175  fSparseEvtShape->Fill(arrayForSparseSo, fWeight);
1176  }
1177  }
1178 
1180  Double_t arrayForSparseRecSphero[5]={ptCand, invMass, spherocity, multForCand, recSpherocity};
1181  fSparseEvtShapeRecSphero->Fill(arrayForSparseRecSphero, fWeight);
1182  }
1183 
1184  if(fReadMC){
1186  Double_t arrayForSparseSoPromptFD[5]={ptCand, invMass, spherocity, multForCand, recSpherocity};
1187 
1188  if(Origin==4) fSparseEvtShapePrompt->Fill(arrayForSparseSoPromptFD, fWeight);
1189  else if(Origin==5) fSparseEvtShapeFeeddown->Fill(arrayForSparseSoPromptFD, fWeight);
1190  }else{
1191  Double_t arrayForSparseSoPromptFD[4]={ptCand, invMass, spherocity, multForCand};
1192 
1193  if(Origin==4) fSparseEvtShapePrompt->Fill(arrayForSparseSoPromptFD, fWeight);
1194  else if(Origin==5) fSparseEvtShapeFeeddown->Fill(arrayForSparseSoPromptFD, fWeight);
1195  }
1196  }
1197  if(labD>=0){
1198  Bool_t keepCase=kTRUE;
1199  if(fPdgMeson==421){
1200  AliAODMCParticle *partD = (AliAODMCParticle*)arrayMC->At(labD);
1201  Int_t code=partD->GetPdgCode();
1202  if(code<0 && iHyp==0) keepCase=kFALSE;
1203  if(code>0 && iHyp==1) keepCase=kFALSE;
1204  }
1205  if(keepCase) FillMCMassHistos(arrayMC,labD, multForCand, spherocity, sphericity, recSpherocity, nchWeight);
1206  }
1207  if(fDoImpPar) fHistMassPtImpPar[0]->Fill(arrayForSparse, fWeight);
1208  }
1209  }
1210  }
1211  if(fSubtractTrackletsFromDau && nSelCand>0){
1212  aveMult/=nSelCand;
1213  fCounterCandidates->StoreEvent(aod,fRDCutsAnalysis,fReadMC,(Int_t)(aveMult+0.5001),spherocity);
1214  }else{
1215  fCounterCandidates->StoreEvent(aod,fRDCutsAnalysis,fReadMC,(Int_t)countCorr,spherocity);
1216  }
1217 
1218  fCounterCandidates->StoreCandidates(aod,nSelectedNoPID,kTRUE);
1219  fCounterCandidates->StoreCandidates(aod,nSelectedPID,kFALSE);
1220  fHistNtrCorrEvSel->Fill(countCorr,nchWeight);
1221  if(nSelectedPID>0) fHistNtrCorrEvWithCand->Fill(countCorr,nchWeight);
1222  if(nSelectedInMassPeak>0) fHistNtrCorrEvWithD->Fill(countCorr,nchWeight);
1223 
1224  if(fFillTrackHisto && fReadMC) FillTrackControlHisto(aod, countCorr, spherocity, genspherocity, nSelectedInMassPeak);
1225 
1226  PostData(1,fOutput);
1227  PostData(2,fListCuts);
1228  PostData(3,fOutputCounters);
1229  if(fReadMC) PostData(5,fOutputEffCorr);
1230 
1231  return;
1232 }
1233 
1234 //________________________________________________________________________
1236 {
1237  // Histos for impact paramter study
1238  // mass . pt , impact parameter , decay length , multiplicity
1239 
1240  Int_t nbins[5]={fNMassBins,200,fNImpParBins,50,100};
1241  Double_t xmin[5]={fLowmasslimit,0.,fLowerImpPar,0.,0.};
1242  Double_t xmax[5]={fUpmasslimit,20.,fHigherImpPar,1.,100.};
1243 
1244  fHistMassPtImpPar[0]=new THnSparseF("hMassPtImpParAll",
1245  "Mass vs. pt vs.imppar - All",
1246  5,nbins,xmin,xmax);
1247  fHistMassPtImpPar[1]=new THnSparseF("hMassPtImpParPrompt",
1248  "Mass vs. pt vs.imppar - promptD",
1249  5,nbins,xmin,xmax);
1250  fHistMassPtImpPar[2]=new THnSparseF("hMassPtImpParBfeed",
1251  "Mass vs. pt vs.imppar - DfromB",
1252  5,nbins,xmin,xmax);
1253  fHistMassPtImpPar[3]=new THnSparseF("hMassPtImpParTrueBfeed",
1254  "Mass vs. pt vs.true imppar -DfromB",
1255  5,nbins,xmin,xmax);
1256  fHistMassPtImpPar[4]=new THnSparseF("hMassPtImpParBkg",
1257  "Mass vs. pt vs.imppar - backgr.",
1258  5,nbins,xmin,xmax);
1259  for(Int_t i=0; i<5;i++){
1260  fOutput->Add(fHistMassPtImpPar[i]);
1261  }
1262 }
1263 
1264 //________________________________________________________________________
1266 {
1267  // Terminate analysis
1268  //
1269  if(fDebug > 1) printf("AnalysisTaskSEDvsMultiplicity: Terminate() \n");
1270 
1271  fOutput = dynamic_cast<TList*> (GetOutputData(1));
1272  if (!fOutput) {
1273  printf("ERROR: fOutput not available\n");
1274  return;
1275  }
1276 
1277  fHistNEvents = dynamic_cast<TH1F*>(fOutput->FindObject("fHistNEvents"));
1278  if(!fHistNEvents){
1279  printf("ERROR: fHistNEvents not available\n");
1280  return;
1281  }
1282  printf("Number of Analyzed Events = %d\n",(Int_t)fHistNEvents->GetBinContent(3));
1283 
1284  return;
1285 }
1286 
1287 //____________________________________________________________________________
1289 {
1290  // Get Estimator Histogram from period event->GetRunNumber();
1291  //
1292  // If you select SPD tracklets in |eta|<1 you should use type == 1
1293  //
1294 
1295  Int_t runNo = event->GetRunNumber();
1296  Int_t period = -1; // pp: 0-LHC10b, 1-LHC10c, 2-LHC10d, 3-LHC10e
1297  // pPb: 0-LHC13b, 1-LHC13c
1298  if (fisPPbData) {
1299  if (runNo>195343 && runNo<195484) period = 0;
1300  if (runNo>195528 && runNo<195678) period = 1;
1301  if (period < 0 || period > 1) return 0;
1302  }
1303  else {
1304  if(runNo>114930 && runNo<117223) period = 0;
1305  if(runNo>119158 && runNo<120830) period = 1;
1306  if(runNo>122373 && runNo<126438) period = 2;
1307  if(runNo>127711 && runNo<130851) period = 3;
1308  if(period<0 || period>3) return 0;
1309  }
1310 
1311  return fMultEstimatorAvg[period];
1312 }
1313 //__________________________________________________________________________________________________
1314 
1316 {
1317  //
1318  // calculating the weight to fill the container
1319  //
1320 
1321  // FNOLL central:
1322  // p0 = 1.63297e-01 --> 0.322643
1323  // p1 = 2.96275e+00
1324  // p2 = 2.30301e+00
1325  // p3 = 2.50000e+00
1326 
1327  // PYTHIA
1328  // p0 = 1.85906e-01 --> 0.36609
1329  // p1 = 1.94635e+00
1330  // p2 = 1.40463e+00
1331  // p3 = 2.50000e+00
1332 
1333  Double_t func1[4] = {0.322643,2.96275,2.30301,2.5};
1334  Double_t func2[4] = {0.36609,1.94635,1.40463,2.5};
1335 
1336  Double_t dndpt_func1 = dNdptFit(pt,func1);
1337  //if(fUseFlatPtWeight) dndpt_func1 = 1./30.;
1338  Double_t dndpt_func2 = dNdptFit(pt,func2);
1339  AliDebug(2,Form("pt = %f, FONLL = %f, Pythia = %f, ratio = %f",pt,dndpt_func1,dndpt_func2,dndpt_func1/dndpt_func2));
1340  return dndpt_func1/dndpt_func2;
1341 }
1342 
1343 //__________________________________________________________________________________________________
1345 {
1346  //
1347  // calculating dNdpt
1348  //
1349 
1350  Double_t denom = TMath::Power((pt/par[1]), par[3] );
1351  Double_t dNdpt = par[0]*pt/TMath::Power(1.+denom, par[2]);
1352 
1353  return dNdpt;
1354 }
1355 
1356 
1357 //__________________________________________________________________________________________________
1359 {
1360  // creates historgam with measured multiplcity distribution in pp 7 TeV collisions (from Eur. Phys. J. C (2010) 68: 345–354)
1361  //
1362  // for Nch > 70 the points were obtainedwith a double NBD distribution
1363  // TF1 *fit1 = new TF1("fit1","[0]*(TMath::Gamma(x+[1])/(TMath::Gamma(x+1)*TMath::Gamma([1])))*(TMath::Power(([2]/[1]),x))*(TMath::Power((1+([2]/[1])),-x-[1]))"); fit1->SetParameter(0,1.);// normalization constant
1364  // fit1->SetParameter(1,1.63); // k parameter
1365  // fit1->SetParameter(2,12.8); // mean multiplicity
1366  Double_t nchbins[82]={0.50,1.50,2.50,3.50,4.50,5.50,6.50,7.50,8.50,9.50,
1367  10.50,11.50,12.50,13.50,14.50,15.50,16.50,17.50,18.50,19.50,
1368  20.50,21.50,22.50,23.50,24.50,25.50,26.50,27.50,28.50,29.50,
1369  30.50,31.50,32.50,33.50,34.50,35.50,36.50,37.50,38.50,39.50,
1370  40.50,41.50,42.50,43.50,44.50,45.50,46.50,47.50,48.50,49.50,
1371  50.50,51.50,52.50,53.50,54.50,55.50,56.50,57.50,58.50,59.50,
1372  60.50,62.50,64.50,66.50,68.50,70.50,72.50,74.50,76.50,78.50,
1373  80.50,82.50,84.50,86.50,88.50,90.50,92.50,94.50,96.50,98.50,
1374  100.50,102.50};
1375  Double_t pch[81]={0.062011,0.072943,0.070771,0.067245,0.062834,0.057383,0.051499,0.04591,0.041109,0.036954,
1376  0.03359,0.030729,0.028539,0.026575,0.024653,0.0229,0.021325,0.019768,0.018561,0.017187,
1377  0.01604,0.014836,0.013726,0.012576,0.011481,0.010393,0.009502,0.008776,0.008024,0.007452,
1378  0.006851,0.006428,0.00594,0.005515,0.005102,0.00469,0.004162,0.003811,0.003389,0.003071,
1379  0.002708,0.002422,0.002184,0.001968,0.00186,0.00165,0.001577,0.001387,0.001254,0.001118,
1380  0.001037,0.000942,0.000823,0.000736,0.000654,0.000579,0.000512,0.00049,0.00045,0.000355,
1381  0.000296,0.000265,0.000193,0.00016,0.000126,0.0000851, 0.0000676,0.0000537,0.0000426, 0.0000338,
1382  0.0000268,0.0000213,0.0000166,0.0000133,0.0000106,0.00000837,0.00000662, 0.00000524,0.00000414, 0.00000327,
1383  0.00000258};
1384 
1385  if(fHistoMeasNch) delete fHistoMeasNch;
1386  fHistoMeasNch=new TH1F("hMeaseNch","",81,nchbins);
1387  for(Int_t i=0; i<81; i++){
1388  fHistoMeasNch->SetBinContent(i+1,pch[i]);
1389  fHistoMeasNch->SetBinError(i+1,0.);
1390  }
1391 }
1392 
1393 //________________________________________________________________________
1394 Bool_t AliAnalysisTaskSEDvsEventShapes::FillTrackControlHisto(AliAODEvent* aod, Int_t nSelTrkCorr, Double_t spherocity, Double_t genspherocity, Int_t nSelectedEvwithCand)
1395 {
1397 
1398  Int_t nTracks=aod->GetNumberOfTracks();
1399  Int_t nSelTracks=0;
1400 
1401  Double_t* etaArr=new Double_t[nTracks];
1402  Double_t* phiArr=new Double_t[nTracks];
1403  Double_t sumpt=0.;
1404 
1405  for(Int_t it=0; it<nTracks; it++) {
1406  AliAODTrack *tr=dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
1407  if(!tr) continue;
1408  Float_t eta = tr->Eta();
1409  Float_t pt = tr->Pt();
1410  Float_t phi = tr->Phi();
1411  if(eta<fetaMin || eta>fetaMax) continue;
1412  if(pt<fptMin || pt>fptMax) continue;
1413  Bool_t fb1 = tr->TestFilterBit(ffiltbit1);
1414  Bool_t fb2 = tr->TestFilterBit(ffiltbit2);
1415  if( !(fb1 || fb2) ) continue;
1416  if(ffiltbit1==0 || ffiltbit2==0){
1417  Int_t status=tr->GetStatus();
1418  if(!(status & AliAODTrack::kTPCrefit)) continue;
1419  }
1420  etaArr[nSelTracks]=eta;
1421  phiArr[nSelTracks]=phi;
1422  nSelTracks++;
1423  sumpt+=pt;
1424  }
1425 
1426  if(nSelTracks<fminMult) return kFALSE;
1427 
1428  if(nSelectedEvwithCand>0) fHistNtrVsnTrackEvWithCand->Fill(nSelTracks, nSelTrkCorr);
1429 
1430  for(Int_t iselt=0; iselt<nSelTracks; iselt++) {
1431  fHistnTrackvsEtavsPhi->Fill(etaArr[iselt], phiArr[iselt], nSelTracks);
1432  fHistTrueSovsMeasSo->Fill(genspherocity, spherocity, nSelTrkCorr);
1433  if(nSelectedEvwithCand>0){
1434  fHistnTrackvsEtavsPhiEvWithCand->Fill(etaArr[iselt], phiArr[iselt], nSelTracks);
1435  fHistTrueSovsMeasSoEvWithCand->Fill(genspherocity, spherocity, nSelTrkCorr);
1436  }
1437  }
1438 
1439  return kTRUE;
1440 
1441 }
1442 
1443 //__________________________________________________________________________________________________
1444 void AliAnalysisTaskSEDvsEventShapes::FillMCMassHistos(TClonesArray *arrayMC, Int_t labD, Double_t countMult, Double_t spherocity, Double_t sphericity, Double_t recSpherocity, Double_t nchWeight)
1445 {
1446  // Function to fill the true MC signal
1447 
1448  AliAODMCParticle *partD = (AliAODMCParticle*)arrayMC->At(labD);
1449 
1450  Double_t mass = partD->M();
1451  Double_t pt = partD->Pt();
1452  Double_t rapid = partD->Y();
1453 
1454  //Weight according to pt, using FONLL
1455  Double_t ptWeight = 1.;
1456  if(fUsePtWeight) ptWeight = GetPtWeight(pt);
1457  fWeight = nchWeight*ptWeight;
1458 
1459  Int_t orig=AliVertexingHFUtils::CheckOrigin(arrayMC,partD,fUseQuarkTag); // Prompt = 4, FeedDown = 5
1460 
1461  //for prompt
1462  if(orig == 4){
1463  //fill histo for prompt
1464  Double_t arrayMCRecoRecSpheroPrompt[5] = {pt, countMult, spherocity, rapid, recSpherocity};
1465  Double_t arrayMCRecoPrompt[4] = {pt, countMult, spherocity, rapid};
1466  if(fRecomputeSpherocity) fMCRecoPrompt->Fill(arrayMCRecoRecSpheroPrompt, fWeight);
1467  else fMCRecoPrompt->Fill(arrayMCRecoPrompt, fWeight);
1468  if(fCalculateSphericity) {
1469  Double_t arrayMCRecoPromptSpheri[4] = {pt, countMult, sphericity, rapid};
1470  fMCRecoPromptSpheri->Fill(arrayMCRecoPromptSpheri, fWeight);
1471  }
1472  }
1473  //for FD
1474  else if(orig == 5){
1475  //fill histo for FD
1476  Double_t arrayMCRecoRecSpheroFeeddown[5] = {pt, countMult, spherocity, rapid, recSpherocity};
1477  Double_t arrayMCRecoFeeddown[4] = {pt, countMult, spherocity, rapid};
1478  if(fRecomputeSpherocity) fMCRecoFeeddown->Fill(arrayMCRecoRecSpheroFeeddown, fWeight);
1479  else fMCRecoFeeddown->Fill(arrayMCRecoFeeddown, fWeight);
1481  Double_t arrayMCRecoFeeddownSpheri[4] = {pt, countMult, sphericity, rapid};
1482  fMCRecoFeeddownSpheri->Fill(arrayMCRecoFeeddownSpheri, fWeight);
1483  }
1484  }
1485 
1486  Double_t arrayMCReco[4] = {pt, countMult, spherocity, rapid};
1487  fMCRecoBothPromptFD->Fill(arrayMCReco, fWeight);
1489  Double_t arrayMCRecoSpheri[4] = {pt, countMult, sphericity, rapid};
1490  fMCRecoBothPromptFDSpheri->Fill(arrayMCRecoSpheri, fWeight);
1491  }
1492 
1493 }
1494 
1495 //__________________________________________________________________________________________________
1496 void AliAnalysisTaskSEDvsEventShapes::FillMCGenAccHistos(AliAODEvent* aod, TClonesArray *arrayMC, AliAODMCHeader *mcHeader, Double_t countMult, Double_t spherocity, Double_t sphericity, Bool_t isEvSel, Double_t nchWeight)
1497 {
1498 
1500 
1501  Int_t nProng=2;
1502  Int_t totPart = arrayMC->GetEntriesFast(); //number of particles
1503  Int_t totTracks = aod->GetNumberOfTracks(); // number of tracks
1504  Double_t recSpherocity = -0.5;
1505 
1506  const Int_t nPart = totPart;
1507  Int_t trkToSkip[nPart];
1508  for(Int_t i=0; i<nPart; i++) trkToSkip[i]=-1; //stores ID of tracks at i'th label of particle
1509 
1511  for(Int_t it=0; it<totTracks; it++){
1512  AliAODTrack *t=dynamic_cast<AliAODTrack*>(aod->GetTrack(it));
1513  if(!t) continue;
1514  if(!(t->TestFilterMask(BIT(4)))) continue;
1515  //if(!(t->HasPointOnITSLayer(0)) && !(t->HasPointOnITSLayer(1))) continue;
1516 
1517  Int_t lab=TMath::Abs(t->GetLabel());
1518  Int_t id=t->GetID();
1519  Double_t pt = t->Pt();
1520  if(id<0) continue;
1521  if(pt<0.3) continue;
1522  AliAODMCParticle* genDup = dynamic_cast<AliAODMCParticle*>(arrayMC->At(lab));
1523  Double_t xver = genDup->Xv();
1524  Double_t yver = genDup->Yv();
1525  Double_t rver = TMath::Sqrt(xver*xver + yver*yver);
1526  if(rver>3) continue;
1527  trkToSkip[lab] = id;
1528  }
1529  }
1530 
1531  if(fPdgMeson==421){
1532  nProng=2;
1533  }else if(fPdgMeson==411 || fPdgMeson==431){
1534  nProng=3;
1535  }
1536 
1537  Int_t nTrkToSkip = nProng;
1538  Int_t idToSkip[nTrkToSkip];
1539  for(Int_t i=0; i<nTrkToSkip; i++) idToSkip[i]=-1;
1540 
1541  Double_t zMCVertex = mcHeader->GetVtxZ(); //vertex MC
1542 
1543  for(Int_t iPart=0; iPart<totPart; iPart++){
1544  AliAODMCParticle* mcGenPart = dynamic_cast<AliAODMCParticle*>(arrayMC->At(iPart));
1545 
1546  if (TMath::Abs(mcGenPart->GetPdgCode()) == fPdgMeson){
1547  Int_t orig=AliVertexingHFUtils::CheckOrigin(arrayMC,mcGenPart,fUseQuarkTag);//Prompt = 4, FeedDown = 5
1548 
1549  Int_t deca = 0;
1550  Bool_t isGoodDecay=kFALSE;
1551  Int_t labDau[4]={-1,-1,-1,-1};
1552  Bool_t isInAcc = kFALSE;
1553  Bool_t isFidAcc = kFALSE;
1554 
1555  if(fPdgMeson==421){
1556  deca=AliVertexingHFUtils::CheckD0Decay(arrayMC,mcGenPart,labDau);
1557  if(mcGenPart->GetNDaughters()!=2) continue;
1558  if(deca==1) isGoodDecay=kTRUE;
1559  }else if(fPdgMeson==411){
1560  deca=AliVertexingHFUtils::CheckDplusDecay(arrayMC,mcGenPart,labDau);
1561  if(deca>0) isGoodDecay=kTRUE;
1562  }else if(fPdgMeson==431){
1563  deca=AliVertexingHFUtils::CheckDsDecay(arrayMC,mcGenPart,labDau);
1564  if(deca==1) isGoodDecay=kTRUE;
1565  }else if(fPdgMeson==413){
1566  deca=AliVertexingHFUtils::CheckDstarDecay(arrayMC,mcGenPart,labDau);
1567  if(deca==1) isGoodDecay=kTRUE;
1568  }
1569 
1570  if(labDau[0]==-1){
1571  continue; //protection against unfilled array of labels
1572  }
1573 
1574  if(fRecomputeSpherocity && isGoodDecay){
1575  for(Int_t iDau=0; iDau<nTrkToSkip; iDau++){
1576  Int_t indexDau = TMath::Abs(mcGenPart->GetDaughter(iDau)); //index of daughter i.e. label
1577  idToSkip[iDau] = trkToSkip[indexDau];
1578  }
1580  }
1581  if(fPdgMeson==421){
1582  //Removal of D0 from Dstar at Generation !!
1583  if(fRemoveD0fromDstar){
1584  Int_t mother = mcGenPart->GetMother();
1585  AliAODMCParticle* mcMoth = dynamic_cast<AliAODMCParticle*>(arrayMC->At(mother));
1586  if(!mcMoth) continue;
1587  if(TMath::Abs(mcMoth->GetPdgCode())==413) continue;
1588  }
1589  }
1590 
1591  Double_t pt = mcGenPart->Pt();
1592  Double_t rapid = mcGenPart->Y();
1593 
1594  //Weight according to pt, using FONLL
1595  Double_t ptWeight = 1.;
1596  if(fUsePtWeight) ptWeight = GetPtWeight(pt);
1597  fWeight = nchWeight*ptWeight;
1598 
1599  isFidAcc=fRDCutsAnalysis->IsInFiducialAcceptance(pt,rapid);
1600  isInAcc=CheckGenAcc(arrayMC,nProng,labDau);
1601 
1602  if(isGoodDecay && TMath::Abs(zMCVertex) < fRDCutsAnalysis->GetMaxVtxZ() && isFidAcc && isInAcc) {
1603  //for prompt
1604  if(orig == 4){
1605  //fill histo for prompt
1606  Double_t arrayMCGenRecSpheroPrompt[5] = {pt, countMult, spherocity, rapid, recSpherocity};
1607  Double_t arrayMCGenPrompt[4] = {pt, countMult, spherocity, rapid};
1608  if(fRecomputeSpherocity) fMCAccGenPrompt->Fill(arrayMCGenRecSpheroPrompt, fWeight);
1609  else fMCAccGenPrompt->Fill(arrayMCGenPrompt, fWeight);
1610  if(isEvSel) fMCAccGenPromptEvSel->Fill(arrayMCGenPrompt, fWeight);
1612  Double_t arrayMCGenPromptSpheri[4] = {pt, countMult, sphericity, rapid};
1613  fMCAccGenPromptSpheri->Fill(arrayMCGenPromptSpheri, fWeight);
1614  }
1615  }
1616  //for FD
1617  else if(orig == 5){
1618  //fill histo for FD
1619  Double_t arrayMCGenRecSpheroFeeddown[5] = {pt, countMult, spherocity, rapid, recSpherocity};
1620  Double_t arrayMCGenFeeddown[4] = {pt, countMult, spherocity, rapid};
1621  if(fRecomputeSpherocity) fMCAccGenFeeddown->Fill(arrayMCGenRecSpheroFeeddown, fWeight);
1622  else fMCAccGenFeeddown->Fill(arrayMCGenFeeddown, fWeight);
1623  if(isEvSel) fMCAccGenFeeddownEvSel->Fill(arrayMCGenFeeddown, fWeight);
1625  Double_t arrayMCGenFeeddownSpheri[4] = {pt, countMult, sphericity, rapid};
1626  fMCAccGenFeeddownSpheri->Fill(arrayMCGenFeeddownSpheri, fWeight);
1627  }
1628  }
1629  else
1630  continue;
1631  }
1632  }
1633  }
1634 }
1635 
1636 //_________________________________________________________________
1637 Bool_t AliAnalysisTaskSEDvsEventShapes::CheckGenAcc(TClonesArray* arrayMC, Int_t nProng, Int_t *labDau)
1638 {
1640  for (Int_t iProng = 0; iProng<nProng; iProng++){
1641  AliAODMCParticle* mcPartDaughter=dynamic_cast<AliAODMCParticle*>(arrayMC->At(labDau[iProng]));
1642  if(!mcPartDaughter) return kFALSE;
1643  Double_t eta = mcPartDaughter->Eta();
1644  Double_t pt = mcPartDaughter->Pt();
1645  if (TMath::Abs(eta) > fEtaAccCut || pt < fPtAccCut) return kFALSE;
1646  }
1647  return kTRUE;
1648 }
Double_t fetaMin
pt limits for acceptance step
Int_t charge
Bool_t CheckGenAcc(TClonesArray *arrayMC, Int_t nProng, Int_t *labDau)
Int_t pdg
TH2F * fHistNtrCorrVsSo
hist of ntracklets vs So
static Int_t CheckDplusDecay(AliStack *stack, Int_t label, Int_t *arrayDauLab)
THnSparseD * fMCRecoPromptSpheri
histo for StepMCGenAcc for D meson feeddown for Sphericity
double Double_t
Definition: External.C:58
Definition: External.C:260
void StoreCandidates(AliVEvent *, Int_t nCand=0, Bool_t flagFilter=kTRUE)
Definition: External.C:236
Double_t DeltaInvMass() const
TH2F * fHistNtrCorrVsSpheri
hist of ntracklets vs Spheri
AliNormalizationCounter * fCounterCandidates
Counter for normalization, uncorrected multiplicity.
Int_t GetIsSelectedCuts() const
Definition: AliRDHFCuts.h:343
static Int_t CheckDstarDecay(AliStack *stack, Int_t label, Int_t *arrayDauLab)
Bool_t HasSelectionBit(Int_t i) const
Int_t MatchToMC(Int_t pdgabs, Int_t pdgabs2prong, Int_t *pdgDg, Int_t *pdgDg2prong, TClonesArray *mcArray, Bool_t isV0=kFALSE) const
Double_t mass
THnSparseD * fMCAccGenFeeddown
histo for StepMCGenAcc for D meson prompt
Double_t ImpParXY() const
TH3F * fHistnTrackvsEtavsPhi
hist. of ntracklets for evnts with a candidate in D mass peak
THnSparseD * fMCRecoFeeddown
histo for StepMCReco for D meson feeddown
THnSparseD * fSparseEvtShapeFeeddown
THnSparse histograms for Prompt D0 vs. Spherocity.
TH2F * fHistNtrVsNchMCPhysicalPrimary
! hist of ntracklets vs Nch (Physical Primary)
TH2F * fHistNtrCorrVsNchMCPhysicalPrimary
! hist of ntracklets vs Nch (Physical Primary)
static Int_t CheckDsDecay(AliStack *stack, Int_t label, Int_t *arrayDauLab)
Int_t GetWhyRejection() const
Definition: AliRDHFCuts.h:300
THnSparseD * fMCAccGenFeeddownEvSel
histo for StepMCGenAcc for D meson prompt with Vertex selection (IsEvSel = kTRUE) ...
TH2F * fHistNtrVsNchMC
hist of ntracklets vs Spheri
static Int_t CheckOrigin(TClonesArray *arrayMC, AliAODMCParticle *mcPart, Bool_t searchUpToQuark=kTRUE)
Double_t GetMaxVtxZ() const
Definition: AliRDHFCuts.h:243
TH2F * fHistNtrCorrVsZvtx
hist of ntracklets vs Zvertex
THnSparseF * fHistMassPtImpPar[5]
histo for StepMCGenAcc for D meson feeddown with Vertex selection (IsEvSel = kTRUE) ...
TH2F * fHistNtrVsZvtx
hist. for No. of events
Bool_t fDoImpPar
Counter for normalization, corrected multiplicity for candidates.
TList * fListProfiles
list send on output slot 3
THnSparseD * fMCRecoPrompt
histo for StepMCGenAcc for D meson feeddown
static Double_t GetVZEROCEqualizedMultiplicity(AliAODEvent *ev)
THnSparseD * fSparseEvtShapewithNoPid
THnSparse histograms for Spherocity.
Class for cuts on AOD reconstructed D+->Kpipi.
TH2F * fHistNtrVsnTrackEvWithCand
hist of ntracklets vs Zvertex
void SetStudyMultiplicity(Bool_t flag, Float_t etaRange)
void FillMCGenAccHistos(AliAODEvent *aod, TClonesArray *arrayMC, AliAODMCHeader *mcHeader, Double_t countMult, Double_t spherocity, Double_t sphericity, Bool_t isEvSel, Double_t nchWeight)
TH2F * fHistNtrCorrVsNchMC
! hist of ntracklets vs Nch (Generated)
int Int_t
Definition: External.C:63
THnSparseD * fMCRecoFeeddownSpheri
histo for StepMCReco for D meson feeddown for Sphericity
unsigned int UInt_t
Definition: External.C:33
float Float_t
Definition: External.C:68
THnSparseD * fMCRecoBothPromptFDSpheri
histo for StepMCReco for D meson feeddown for Sphericity
TH3F * fHistTrueSovsMeasSo
! hist. of number of tracks passing track selection for spherocity calculation vs eta vs...
TH2F * fHistNtrVsNchMCPrimary
! hist of ntracklets vs Nch (Primary)
TH3F * fHistNchMCVsNchMCPrimaryVsNchMCPhysicalPrimary
! hist of Nch (generated) vs Nch (Primary) vs Nch (Physical Primary)
Bool_t fFillTrackHisto
flag for quark/hadron level identification of prompt and feeddown
AliNormalizationCounter * fCounterU
Counter for normalization, corrected multiplicity.
static Int_t CheckD0Decay(AliStack *stack, Int_t label, Int_t *arrayDauLab)
static Double_t GetGeneratedSpherocity(TClonesArray *arrayMC, Double_t etaMin=-0.8, Double_t etaMax=0.8, Double_t ptMin=0.15, Double_t ptMax=10., Int_t minMult=3, Double_t phiStepSizeDeg=0.1)
void SetStudySpherocity(Bool_t flag, Double_t nsteps=100.)
THnSparseD * fSparseEvtShapePrompt
THnSparse histograms for D0 vs. Spherocity.
TH3F * fHistTrueSovsMeasSoEvWithCand
! hist. of number of tracks passing track selection for spherocity calculation vs eta vs...
TH1F * fHistNtrCorrEvSel
hist. of ntracklets for physics selection only selected events
void SetMassLimits(Double_t lowlimit, Double_t uplimit)
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
TH1F * fHistNtrCorrEvWithD
hist. of ntracklets for evnts with a candidate
Bool_t IsEventRejectedDuePhysicsSelection() const
Definition: AliRDHFCuts.h:329
TProfile * GetEstimatorHistogram(const AliVEvent *event)
TH1F * fHistNtrCorrEvWithCand
hist. of ntracklets for selected events
Bool_t FillTrackControlHisto(AliAODEvent *aod, Int_t nSelTrkCorr, Double_t spherocity, Double_t genspherocity, Int_t nSelectedEvwithCand)
THnSparseD * fMCAccGenFeeddownSpheri
histo for StepMCGenAcc for D meson prompt for Sphericity
ClassImp(AliAnalysisTaskSEDvsEventShapes) AliAnalysisTaskSEDvsEventShapes
static Double_t GetCorrectedNtracklets(TProfile *estimatorAvg, Double_t uncorrectedNacc, Double_t vtxZ, Double_t refMult)
THnSparseD * fMCRecoBothPromptFD
histo for StepMCReco for D meson feeddown
Bool_t IsEventSelected(AliVEvent *event)
THnSparseD * fMCAccGenPromptSpheri
histo for StepMCReco for D meson Both Prompt Feeddown
void StoreEvent(AliVEvent *, AliRDHFCuts *, Bool_t mc=kFALSE, Int_t multiplicity=-9999, Double_t spherocity=-99.)
Double_t fUpmasslimit
histograms for impact paramter studies
Double_t dNdptFit(Float_t pt, Double_t *par)
THnSparseD * fMCAccGenPromptEvSel
histo for StepMCReco for D meson Both Prompt Feeddown for Sphericity
void FillMCMassHistos(TClonesArray *arrayMC, Int_t labD, Double_t countMult, Double_t spherocity, Double_t sphericity, Double_t recSpherocity, Double_t nchWeight)
TH3F * fHistnTrackvsEtavsPhiEvWithCand
! hist. of number of tracks passing track selection for spherocity calculation vs eta vs...
Double_t fPtAccCut
eta limits for acceptance step
Bool_t IsSelected(TObject *obj)
Definition: AliRDHFCuts.h:277
const char Option_t
Definition: External.C:48
TList * fListCuts
list send on output slot 1
THnSparseD * fSparseEvtShapeRecSphero
THnSparse histograms for feeddown D0 vs. Spherocity.
TH1F * fHistGenPrimaryParticlesInelGt0
!hist. of geenrated multiplcity
static Double_t GetVZEROAEqualizedMultiplicity(AliAODEvent *ev)
Utilities for V0 multiplicity checks.
static Double_t GetTrueImpactParameterDzero(AliAODMCHeader *mcHeader, TClonesArray *arrayMC, AliAODMCParticle *partDp)
Functions for computing true impact parameter of D meson.
const Int_t nbins
bool Bool_t
Definition: External.C:53
TH2F * fHistNtrCorrVsNchMCPrimary
! hist of ntracklets vs Nch (Primary)
AliAODRecoDecayHF2Prong * Get2Prong() const
THnSparseD * fMCAccGenPrompt
THnSparse histograms for Both Prompt and feeddown D0 vs. Spherocity.
static Double_t GetTrueImpactParameterDplus(AliAODMCHeader *mcHeader, TClonesArray *arrayMC, AliAODMCParticle *partDp)
Double_t fEtaAccCut
flag for filling track control histograms
virtual Bool_t IsInFiducialAcceptance(Double_t, Double_t) const
Definition: AliRDHFCuts.h:297
static Double_t GetSphericity(AliAODEvent *aod, Double_t etaMin=-0.8, Double_t etaMax=0.8, Double_t ptMin=0.15, Double_t ptMax=10., Int_t filtbit1=256, Int_t filtbit2=512, Int_t minMult=3)
TH2F * fHistNtrVsSpheri
hist of ntracklets vs So
TH1F * fHistNEvents
list send on output slot 5
static Double_t GetSpherocity(AliAODEvent *aod, Double_t etaMin=-0.8, Double_t etaMax=0.8, Double_t ptMin=0.15, Double_t ptMax=10., Int_t filtbit1=256, Int_t filtbit2=512, Int_t minMult=3, Double_t phiStepSizeDeg=0.1, Int_t nTrksToSkip=0, Int_t *idToSkip=0x0)
Functions for event shape variables.