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