AliPhysics  vAN-20150810 (e477c3b)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliCFTaskVertexingHFCutVarFDSub.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 //-----------------------------------------------------------------------
17 // Class for HF corrections as a function of many variables
18 // 6 Steps introduced: MC, MC Acc, Reco, Reco Acc, Reco Acc + ITS Cl,
19 // Reco Acc + ITS Cl + PPR cuts
20 // 12 variables used: pt, y, cosThetaStar, ptPi, ptK, ct,
21 // dca, d0Pi, d0K, d0Pixd0K, cosPointingAngle, phi
22 //
23 //-----------------------------------------------------------------------
24 // Author : C. Zampolli, CERN
25 // D. Caffarri, Univ & INFN Padova caffarri@pd.infn.it
26 //-----------------------------------------------------------------------
27 //-----------------------------------------------------------------------
28 // Base class for HF Unfolding (pt and eta)
29 // correlation matrix filled at Acceptance and PPR level
30 // Author: A.Grelli , Utrecht - agrelli@uu.nl
31 //-----------------------------------------------------------------------
32 //-----------------------------------------------------------------------
33 // Modified by for analysis of D0 candidates using the cut-variation
34 // feed-down subtraction method
35 // by
36 // Andrea Rossi <Andrea.Rossi@cern.ch>
37 // Felix Reidt <Felix.Reidt@cern.ch>
38 //-----------------------------------------------------------------------
40 
41 #include "TCanvas.h"
42 #include "TProfile.h"
43 #include "TH1I.h"
44 #include "TH3F.h"
45 #include "TH1F.h"
46 #include "TStyle.h"
47 #include "TFile.h"
48 #include "TF1.h"
49 
50 #include "AliCFManager.h"
51 #include "AliCFContainer.h"
52 #include "AliLog.h"
53 #include "AliInputEventHandler.h"
54 #include "AliAnalysisManager.h"
55 #include "AliAODHandler.h"
56 #include "AliAODEvent.h"
57 #include "AliAODRecoDecayHF.h"
61 #include "AliAODRecoCascadeHF.h"
62 #include "AliAODMCParticle.h"
63 #include "AliAODMCHeader.h"
64 #include "AliESDtrack.h"
65 
66 #include "THnSparse.h"
67 #include "AliESDtrackCuts.h"
68 #include "AliRDHFCuts.h"
69 #include "AliRDHFCutsD0toKpi.h"
72 #include "AliRDHFCutsDstoKKpi.h"
73 #include "AliRDHFCutsLctopKpi.h"
74 #include "AliRDHFCutsD0toKpipipi.h"
75 #include "AliRDHFCutsLctoV0.h"
76 #include "AliCFVertexingHF2Prong.h"
77 #include "AliCFVertexingHF3Prong.h"
80 #include "AliCFVertexingHF.h"
81 #include "AliVertexingHFUtils.h"
82 #include "AliAnalysisDataSlot.h"
83 #include "AliAnalysisDataContainer.h"
84 #include "AliAnalysisTaskSE.h"
85 #include "AliPIDResponse.h"
86 #include "AliHFsubtractBFDcuts.h"
87 
88 //__________________________________________________________________________
90  AliAnalysisTaskSE(),
91  fCFManager(0x0),
92  fHistEventsProcessed(0x0),
93  fCorrelation(0x0),
94  fListProfiles(0),
95  fCountMC(0),
96  fCountAcc(0),
97  fCountVertex(0),
98  fCountRefit(0),
99  fCountReco(0),
100  fCountRecoAcc(0),
101  fCountRecoITSClusters(0),
102  fCountRecoPPR(0),
103  fCountRecoPID(0),
104  fEvents(0),
105  fDecayChannel(0),
106  fFillFromGenerated(kFALSE),
107  fOriginDselection(0),
108  fAcceptanceUnf(kTRUE),
109  fCuts(0),
110  fUseWeight(kFALSE),
111  fWeight(1.),
112  fUseFlatPtWeight(kFALSE),
113  fUseZWeight(kFALSE),
114  fUseNchWeight(kFALSE),
115  fUseTrackletsWeight(kFALSE),
116  fUseMultRatioAsWeight(kFALSE),
117  fNvar(0),
118  fPartName(""),
119  fDauNames(""),
120  fSign(2),
121  fCentralitySelection(kTRUE),
122  fFakeSelection(0),
123  fRejectIfNoQuark(kTRUE),
124  fUseMCVertex(kFALSE),
125  fDsOption(1),
126  fGenDsOption(3),
127  fConfiguration(kCheetah), // by default, setting the fast configuration
128  fFuncWeight(0x0),
129  fHistoPtWeight(0x0),
130  fHistoMeasNch(0x0),
131  fHistoMCNch(0x0),
132  fResonantDecay(0),
133  fLctoV0bachelorOption(1),
134  fGenLctoV0bachelorOption(0),
135  fUseSelectionBit(kTRUE),
136  fPDGcode(0),
137  fMultiplicityEstimator(kNtrk10),
138  fRefMult(9.26),
139  fZvtxCorrectedNtrkEstimator(kFALSE),
140  fIsPPData(kFALSE),
141  fIsPPbData(kFALSE),
142  fUseAdditionalCuts(kFALSE),
143  fUseCutsForTMVA(kFALSE),
144  fUseCascadeTaskForLctoV0bachelor(kFALSE),
145  fCutOnMomConservation(0.00001),
146  fObjSpr(0x0),
147  fhSparseCutVar(0x0),
148  fhPtCutVar(0x0),
149  fhBptCutVar(0x0),
150  fListBdecays(0x0),
151  fQAHists(0x0)
152 {
153  //
154  //Default ctor
155  //
156  for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=0;
157 }
158 //___________________________________________________________________________
160  AliAnalysisTaskSE(name),
161  fCFManager(0x0),
162  fHistEventsProcessed(0x0),
163  fCorrelation(0x0),
164  fListProfiles(0),
165  fCountMC(0),
166  fCountAcc(0),
167  fCountVertex(0),
168  fCountRefit(0),
169  fCountReco(0),
170  fCountRecoAcc(0),
171  fCountRecoITSClusters(0),
172  fCountRecoPPR(0),
173  fCountRecoPID(0),
174  fEvents(0),
175  fDecayChannel(0),
176  fFillFromGenerated(kFALSE),
177  fOriginDselection(0),
178  fAcceptanceUnf(kTRUE),
179  fCuts(cuts),
180  fUseWeight(kFALSE),
181  fWeight(1.),
182  fUseFlatPtWeight(kFALSE),
183  fUseZWeight(kFALSE),
184  fUseNchWeight(kFALSE),
185  fUseTrackletsWeight(kFALSE),
186  fUseMultRatioAsWeight(kFALSE),
187  fNvar(0),
188  fPartName(""),
189  fDauNames(""),
190  fSign(2),
191  fCentralitySelection(kTRUE),
192  fFakeSelection(0),
193  fRejectIfNoQuark(kTRUE),
194  fUseMCVertex(kFALSE),
195  fDsOption(1),
196  fGenDsOption(3),
197  fConfiguration(kCheetah), // by default, setting the fast configuration
198  fFuncWeight(func),
199  fHistoPtWeight(0x0),
200  fHistoMeasNch(0x0),
201  fHistoMCNch(0x0),
202  fResonantDecay(0),
203  fLctoV0bachelorOption(1),
204  fGenLctoV0bachelorOption(0),
205  fUseSelectionBit(kTRUE),
206  fPDGcode(0),
207  fMultiplicityEstimator(kNtrk10),
208  fRefMult(9.26),
209  fZvtxCorrectedNtrkEstimator(kFALSE),
210  fIsPPData(kFALSE),
211  fIsPPbData(kFALSE),
212  fUseAdditionalCuts(kFALSE),
213  fUseCutsForTMVA(kFALSE),
214  fUseCascadeTaskForLctoV0bachelor(kFALSE),
215  fCutOnMomConservation(0.00001),
216  fObjSpr(0x0),
217  fhSparseCutVar(0x0),
218  fhPtCutVar(0x0),
219  fhBptCutVar(0x0),
220  fListBdecays(0x0),
221  fQAHists(0x0)
222 {
223  //
224  // Constructor. Initialization of Inputs and Outputs
225  //
226  /*
227  DefineInput(0) and DefineOutput(0)
228  are taken care of by AliAnalysisTaskSE constructor
229  */
230  DefineOutput(1,TH1I::Class());
231  DefineOutput(2,AliCFContainer::Class());
232  DefineOutput(3,THnSparseD::Class());
233  DefineOutput(4,AliRDHFCuts::Class());
234  for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=0;
235  DefineOutput(5,TList::Class()); // slot #5 keeps the zvtx Ntrakclets correction profiles
236  DefineOutput(6,THnSparseF::Class());
237  DefineOutput(7,TH3F::Class());
238  DefineOutput(8,TH1F::Class());
239  DefineOutput(9,TList::Class());
240  DefineOutput(10,TList::Class());
241 
242  fCuts->PrintAll();
243 }
244 
245 //___________________________________________________________________________
247 {
248  //
249  // Assignment operator
250  //
251  if (this!=&c) {
252  AliAnalysisTaskSE::operator=(c) ;
255  fCuts = c.fCuts;
260  for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=c.fMultEstimatorAvg[i];
261  fObjSpr=c.fObjSpr;
265  fListBdecays=c.fListBdecays; // FIXME: TList copy contructor not implemented
266  fQAHists=c.fQAHists; // FIXME: TList copy contructor not implemented
267  }
268  return *this;
269 }
270 
271 //___________________________________________________________________________
273  AliAnalysisTaskSE(c),
274  fCFManager(c.fCFManager),
275  fHistEventsProcessed(c.fHistEventsProcessed),
276  fCorrelation(c.fCorrelation),
277  fListProfiles(c.fListProfiles),
278  fCountMC(c.fCountMC),
279  fCountAcc(c.fCountAcc),
280  fCountVertex(c.fCountVertex),
281  fCountRefit(c.fCountRefit),
282  fCountReco(c.fCountReco),
283  fCountRecoAcc(c.fCountRecoAcc),
284  fCountRecoITSClusters(c.fCountRecoITSClusters),
285  fCountRecoPPR(c.fCountRecoPPR),
286  fCountRecoPID(c.fCountRecoPID),
287  fEvents(c.fEvents),
288  fDecayChannel(c.fDecayChannel),
289  fFillFromGenerated(c.fFillFromGenerated),
290  fOriginDselection(c.fOriginDselection),
291  fAcceptanceUnf(c.fAcceptanceUnf),
292  fCuts(c.fCuts),
293  fUseWeight(c.fUseWeight),
294  fWeight(c.fWeight),
295  fUseFlatPtWeight(c.fUseFlatPtWeight),
296  fUseZWeight(c.fUseZWeight),
297  fUseNchWeight(c.fUseNchWeight),
298  fUseTrackletsWeight(c.fUseTrackletsWeight),
299  fUseMultRatioAsWeight(c.fUseMultRatioAsWeight),
300  fNvar(c.fNvar),
301  fPartName(c.fPartName),
302  fDauNames(c.fDauNames),
303  fSign(c.fSign),
304  fCentralitySelection(c.fCentralitySelection),
305  fFakeSelection(c.fFakeSelection),
306  fRejectIfNoQuark(c.fRejectIfNoQuark),
307  fUseMCVertex(c.fUseMCVertex),
308  fDsOption(c.fDsOption),
309  fGenDsOption(c.fGenDsOption),
310  fConfiguration(c.fConfiguration),
311  fFuncWeight(c.fFuncWeight),
312  fHistoPtWeight(c.fHistoPtWeight),
313  fHistoMeasNch(c.fHistoMeasNch),
314  fHistoMCNch(c.fHistoMCNch),
315  fResonantDecay(c.fResonantDecay),
316  fLctoV0bachelorOption(c.fLctoV0bachelorOption),
317  fGenLctoV0bachelorOption(c.fGenLctoV0bachelorOption),
318  fUseSelectionBit(c.fUseSelectionBit),
319  fPDGcode(c.fPDGcode),
320  fMultiplicityEstimator(c.fMultiplicityEstimator),
321  fRefMult(c.fRefMult),
322  fZvtxCorrectedNtrkEstimator(c.fZvtxCorrectedNtrkEstimator),
323  fIsPPData(c.fIsPPData),
324  fIsPPbData(c.fIsPPbData),
325  fUseAdditionalCuts(c.fUseAdditionalCuts),
326  fUseCutsForTMVA(c.fUseCutsForTMVA),
327  fUseCascadeTaskForLctoV0bachelor(c.fUseCascadeTaskForLctoV0bachelor),
328  fCutOnMomConservation(c.fCutOnMomConservation),
329  fObjSpr(c.fObjSpr),
330  fhSparseCutVar(c.fhSparseCutVar),
331  fhPtCutVar(c.fhPtCutVar),
332  fhBptCutVar(c.fhBptCutVar),
333  fListBdecays(c.fListBdecays), // FIXME: TList copy contructor not implemented
334  fQAHists(c.fQAHists) // FIXME: TList copy contructor not implemented
335 {
336  //
337  // Copy Constructor
338  //
339  for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=c.fMultEstimatorAvg[i];
340 }
341 
342 //___________________________________________________________________________
344 {
345  //
346  //destructor
347  //
348  if (fCFManager) { delete fCFManager; fCFManager=0x0; }
350  if (fCorrelation) { delete fCorrelation; fCorrelation=0x0; }
351  if (fListProfiles) { delete fListProfiles; fListProfiles=0x0; }
352  if (fCuts) { delete fCuts; fCuts=0x0; }
353  if (fFuncWeight) { delete fFuncWeight; fFuncWeight=0x0; }
354  if (fHistoPtWeight) { delete fHistoPtWeight; fHistoPtWeight=0x0; }
355  if (fHistoMeasNch) { delete fHistoMeasNch; fHistoMeasNch=0x0; }
356  if (fHistoMCNch) { delete fHistoMCNch; fHistoMCNch=0x0; }
357  if (fObjSpr) { delete fObjSpr; fObjSpr=0x0; }
358  if (fhSparseCutVar) { delete fhSparseCutVar; fhSparseCutVar=0x0; }
359  if (fhPtCutVar) { delete fhPtCutVar; fhPtCutVar=0x0; }
360  if (fhBptCutVar) { delete fhBptCutVar; fhBptCutVar=0x0; }
361  if (fListBdecays) { delete fListBdecays; fListBdecays=0x0; }
362  if (fQAHists) { delete fQAHists; fQAHists=0x0; }
363  for(Int_t i=0; i<4; i++) {
364  if(fMultEstimatorAvg[i]) { delete fMultEstimatorAvg[i]; fMultEstimatorAvg[i]=0x0; }
365  }
366 }
367 
368 //__________________________________________________________________________
370 {
371  //
372  // Initialization
373  //
374 
375  if (fDebug>1) printf("AliCFTaskVertexingHFCutVarFDSub::Init()");
376  if(fUseWeight && fUseZWeight) { AliFatal("Can not use at the same time pt and z-vtx weights, please choose"); return; }
377  if(fUseWeight && fUseNchWeight) { AliInfo("Beware, using at the same time pt and Nch weights, please check"); }
378  if(fUseNchWeight && !fHistoMCNch) { AliFatal("Need to pass the MC Nch distribution to use Nch weights"); return; }
380 
381  AliRDHFCuts *copyfCuts = 0x0;
382  if (!fCuts){
383  AliFatal("No cuts defined - Exiting...");
384  return;
385  }
386 
387  switch (fDecayChannel){
388  case 2:{
389  fPDGcode = 421;
390  copyfCuts = new AliRDHFCutsD0toKpi(*(static_cast<AliRDHFCutsD0toKpi*>(fCuts)));
391  switch (fConfiguration) {
392  case kSnail: // slow configuration: all variables in
393  fNvar = 16;
394  break;
395  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
396  fNvar = 8;
397  break;
398  }
399  fPartName="D0";
400  fDauNames="K+pi";
401  break;
402  }
403  case 21:{
404  fPDGcode = 413;
405  copyfCuts = new AliRDHFCutsDStartoKpipi(*(static_cast<AliRDHFCutsDStartoKpipi*>(fCuts)));
406  switch (fConfiguration) {
407  case kSnail: // slow configuration: all variables in
408  fNvar = 16;
409  break;
410  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
411  fNvar = 8;
412  break;
413  }
414  fPartName="Dstar";
415  fDauNames="K+pi+pi";
416  break;
417  }
418  case 22:{
419  fPDGcode = 4122;
420  copyfCuts = new AliRDHFCutsLctoV0(*(static_cast<AliRDHFCutsLctoV0*>(fCuts)));
421  switch (fConfiguration) {
422  case kSnail: // slow configuration: all variables in
423  fNvar = 16;
424  break;
425  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
426  fNvar = 8;
427  break;
428  }
429  fPartName="Lambdac";
430  fDauNames="V0+bachelor";
431  break;
432  }
433  case 31:{
434  fPDGcode = 411;
435  copyfCuts = new AliRDHFCutsDplustoKpipi(*(static_cast<AliRDHFCutsDplustoKpipi*>(fCuts)));
436  switch (fConfiguration) {
437  case kSnail: // slow configuration: all variables in
438  fNvar = 14;
439  break;
440  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
441  fNvar = 8;
442  break;
443  }
444  fPartName="Dplus";
445  fDauNames="K+pi+pi";
446  break;
447  }
448  case 32:{
449  fPDGcode = 4122;
450  copyfCuts = new AliRDHFCutsLctopKpi(*(static_cast<AliRDHFCutsLctopKpi*>(fCuts)));
451  switch (fConfiguration) {
452  case kSnail: // slow configuration: all variables in
453  fNvar = 18;
454  break;
455  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
456  fNvar = 8;
457  break;
458  }
459  fPartName="Lambdac";
460  fDauNames="p+K+pi";
461  break;
462  }
463  case 33:{
464  fPDGcode = 431;
465  copyfCuts = new AliRDHFCutsDstoKKpi(*(static_cast<AliRDHFCutsDstoKKpi*>(fCuts)));
466  switch (fConfiguration) {
467  case kSnail: // slow configuration: all variables in
468  fNvar = 14;
469  break;
470  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
471  fNvar = 8;
472  break;
473  }
474  fPartName="Ds";
475  fDauNames="K+K+pi";
476  break;
477  }
478  case 4:{
479  fPDGcode = 421;
480  copyfCuts = new AliRDHFCutsD0toKpipipi(*(static_cast<AliRDHFCutsD0toKpipipi*>(fCuts)));
481  switch (fConfiguration) {
482  case kSnail: // slow configuration: all variables in
483  fNvar = 16;
484  break;
485  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
486  fNvar = 8;
487  break;
488  }
489  fPartName="D0";
490  fDauNames="K+pi+pi+pi";
491  break;
492  }
493  default:
494  AliFatal("The decay channel MUST be defined according to AliCFVertexing::DecayChannel - Exiting...");
495  break;
496  }
497 
498  const char* nameoutput=GetOutputSlot(4)->GetContainer()->GetName();
499  if (copyfCuts){
500  copyfCuts->SetName(nameoutput);
501 
502  //Post the data
503  PostData(4, copyfCuts);
504  }
505  else{
506  AliFatal("Failing initializing AliRDHFCuts object - Exiting...");
507  }
508 
509  fListProfiles = new TList();
510  fListProfiles->SetOwner();
511  TString period[4];
512  Int_t nProfiles=4;
513 
514  if (fIsPPbData) { //if pPb, use only two estimator histos
515  period[0] = "LHC13b"; period[1] = "LHC13c";
516  nProfiles = 2;
517  } else { // else assume pp (four histos for LHC10)
518  period[0] = "LHC10b"; period[1] = "LHC10c"; period[2] = "LHC10d"; period[3] = "LHC10e";
519  nProfiles = 4;
520  }
521 
522  for(Int_t i=0; i<nProfiles; i++){
523  if(fMultEstimatorAvg[i]){
524  TProfile* hprof=new TProfile(*fMultEstimatorAvg[i]);
525  hprof->SetName(Form("ProfileTrkVsZvtx%s\n",period[i].Data()));
526  fListProfiles->Add(hprof);
527  }
528  }
529 
530  // Save also the weight functions or histograms
535 
536  PostData(5,fListProfiles);
537 
538  return;
539 }
540 
541 //_________________________________________________
543 {
544  //
545  // Main loop function
546  //
547 
548  PostData(1,fHistEventsProcessed) ;
549  PostData(2,fCFManager->GetParticleContainer()) ;
550  PostData(3,fCorrelation);
551 
552  AliDebug(3,Form("*** Processing event %d\n", fEvents));
553 
554  if (fFillFromGenerated){
555  AliWarning("Flag to fill container with generated value ON ---> dca, d0pi, d0K, d0xd0, cosPointingAngle will be set as dummy!");
556  }
557 
558  if (!fInputEvent) {
559  Error("UserExec","NO EVENT FOUND!");
560  return;
561  }
562 
563  AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);
564 
565  TClonesArray *arrayBranch=0;
566 
567  if(!aodEvent && AODEvent() && IsStandardAOD()) {
568  // In case there is an AOD handler writing a standard AOD, use the AOD
569  // event in memory rather than the input (ESD) event.
570  aodEvent = dynamic_cast<AliAODEvent*> (AODEvent());
571  // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
572  // have to taken from the AOD event hold by the AliAODExtension
573  AliAODHandler* aodHandler = (AliAODHandler*)
574  ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
575  if(aodHandler->GetExtensions()) {
576  AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
577  AliAODEvent *aodFromExt = ext->GetAOD();
578 
579  switch (fDecayChannel){
580  case 2:{
581  arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi");
582  break;
583  }
584  case 21:{
585  arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");
586  break;
587  }
588  case 22:{
589  arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("CascadesHF");
590  break;
591  }
592  case 31:
593  case 32:
594  case 33:{
595  arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
596  break;
597  }
598  case 4:{
599  arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm4Prong");
600  break;
601  }
602  default:
603  break;
604  }
605  }
606  }
607  else {
608  switch (fDecayChannel){
609  case 2:{
610  arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("D0toKpi");
611  break;
612  }
613  case 21:{
614  arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("Dstar");
615  break;
616  }
617  case 22:{
618  arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("CascadesHF");
619  break;
620  }
621  case 31:
622  case 32:
623  case 33:{
624  arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("Charm3Prong");
625  break;
626  }
627  case 4:{
628  arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("Charm4Prong");
629  break;
630  }
631  default:
632  break;
633  }
634  }
635 
636  AliAODVertex *aodVtx = (AliAODVertex*)aodEvent->GetPrimaryVertex();
637  if (!aodVtx) {
638  AliDebug(3, "The event was skipped due to missing vertex");
639  return;
640  }
641 
642  if (!arrayBranch) {
643  AliError("Could not find array of HF vertices");
644  return;
645  }
646 
647  fEvents++;
648 
649  fCFManager->SetRecEventInfo(aodEvent);
650  fCFManager->SetMCEventInfo(aodEvent);
651 
652  //******** DEFINE number of variables of the container***** for now set at 13, in the future in the config macro.
653 
654  //loop on the MC event
655 
656  TClonesArray* mcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));
657  if (!mcArray) {
658  AliError("Could not find Monte-Carlo in AOD");
659  return;
660  }
661  Int_t icountMC = 0;
662  Int_t icountAcc = 0;
663  Int_t icountReco = 0;
664  Int_t icountVertex = 0;
665  Int_t icountRefit = 0;
666  Int_t icountRecoAcc = 0;
667  Int_t icountRecoITSClusters = 0;
668  Int_t icountRecoPPR = 0;
669  Int_t icountRecoPID = 0;
670  Int_t cquarks = 0;
671 
672  AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
673  if (!mcHeader) {
674  AliError("Could not find MC Header in AOD");
675  return;
676  }
677 
678  fHistEventsProcessed->Fill(0.5);
679 
680  Double_t* containerInput = new Double_t[fNvar];
681  Double_t* containerInputMC = new Double_t[fNvar];
682 
683 
684  AliCFVertexingHF* cfVtxHF=0x0;
685  switch (fDecayChannel){
686  case 2:{
687  cfVtxHF = new AliCFVertexingHF2Prong(mcArray, fOriginDselection);
688  break;
689  }
690  case 21:{
691  cfVtxHF = new AliCFVertexingHFCascade(mcArray, fOriginDselection);
692  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGcascade(413);
693  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGbachelor(211);
694  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaugh(421);
695  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaughForMC(421);
696  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaughPositive(211);
697  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaughNegative(321);
698  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPrimaryVertex(aodVtx);
699  break;
700  }
701  case 22:{
702  // Lc -> K0S+proton
704  cfVtxHF = new AliCFVertexingHFCascade(mcArray, fOriginDselection);
705  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGcascade(4122);
706  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGbachelor(2212);
707  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaugh(310);
708  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaughForMC(311);
709  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaughPositive(211);
710  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaughNegative(211);
711  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPrimaryVertex(aodVtx);
712  ((AliCFVertexingHFCascade*)cfVtxHF)->SetCutOnMomConservation(fCutOnMomConservation);
713  if (fUseAdditionalCuts) ((AliCFVertexingHFCascade*)cfVtxHF)->SetUseCutsForTMVA(fUseCutsForTMVA);
714  }
715  else {
717  }
718  break;
719  }
720  case 31:
721  // case 32:
722  case 33:{
723  cfVtxHF = new AliCFVertexingHF3Prong(mcArray, fOriginDselection, fDecayChannel);
724  if(fDecayChannel==33){
725  ((AliCFVertexingHF3Prong*)cfVtxHF)->SetGeneratedDsOption(fGenDsOption);
726  }
727  break;
728  }
729  case 32:{
731  }
732  case 4:{
733  //cfVtxHF = new AliCFVertexingHF4Prong(mcArray, originDselection); // not there yet
734  break;
735  }
736  default:
737  break;
738  }
739  if (!cfVtxHF){
740  AliError("No AliCFVertexingHF initialized");
741  delete[] containerInput;
742  delete[] containerInputMC;
743  return;
744  }
745 
746  Double_t zPrimVertex = aodVtx ->GetZ();
747  Double_t zMCVertex = mcHeader->GetVtxZ();
748  Int_t runnumber = aodEvent->GetRunNumber();
749 
750  // Multiplicity definition with tracklets
751  Double_t nTracklets = 0;
752  Int_t nTrackletsEta10 = static_cast<Int_t>(AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aodEvent,-1.,1.));
753  Int_t nTrackletsEta16 = static_cast<Int_t>(AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aodEvent,-1.6,1.6));
754  nTracklets = (Double_t)nTrackletsEta10;
755  if(fMultiplicityEstimator==kNtrk10to16) { nTracklets = (Double_t)(nTrackletsEta16 - nTrackletsEta10); }
756 
757  // Apply the Ntracklets z-vtx data driven correction
759  TProfile* estimatorAvg = GetEstimatorHistogram(aodEvent);
760  if(estimatorAvg) {
761  Int_t nTrackletsEta10Corr = static_cast<Int_t>(AliVertexingHFUtils::GetCorrectedNtracklets(estimatorAvg,nTrackletsEta10,zPrimVertex,fRefMult));
762  Int_t nTrackletsEta16Corr = static_cast<Int_t>(AliVertexingHFUtils::GetCorrectedNtracklets(estimatorAvg,nTrackletsEta16,zPrimVertex,fRefMult));
763  nTracklets = (Double_t)nTrackletsEta10Corr;
764  if(fMultiplicityEstimator==kNtrk10to16) { nTracklets = (Double_t)(nTrackletsEta16Corr - nTrackletsEta10Corr); }
765  }
766  }
767 
768 
769  fWeight=1.;
770  if(fUseZWeight) fWeight *= GetZWeight(zMCVertex,runnumber);
771  if(fUseNchWeight){
772  Int_t nChargedMCPhysicalPrimary=AliVertexingHFUtils::GetGeneratedPhysicalPrimariesInEtaRange(mcArray,-1.0,1.0);
773  if(!fUseTrackletsWeight) fWeight *= GetNchWeight(nChargedMCPhysicalPrimary);
774  else fWeight *= GetNchWeight(static_cast<Int_t>(nTracklets));
775  AliDebug(2,Form("Using Nch weights, Mult=%d Weight=%f\n",nChargedMCPhysicalPrimary,fWeight));
776  }
777  Double_t eventWeight=fWeight;
778 
779  if (TMath::Abs(zMCVertex) > fCuts->GetMaxVtxZ()){
780  AliDebug(3,Form("z coordinate of MC vertex = %f, it was required to be within [-%f, +%f], skipping event", zMCVertex, fCuts->GetMaxVtxZ(), fCuts->GetMaxVtxZ()));
781  delete[] containerInput;
782  delete[] containerInputMC;
783  delete cfVtxHF;
784  return;
785  }
786 
787  if(aodEvent->GetTriggerMask()==0 &&
788  (runnumber>=195344 && runnumber<=195677)){
789  AliDebug(3,"Event rejected because of null trigger mask");
790  delete[] containerInput;
791  delete[] containerInputMC;
792  delete cfVtxHF;
793  return;
794  }
795 
796  AliESDtrackCuts** trackCuts = new AliESDtrackCuts*[cfVtxHF->GetNProngs()];
797  if (fDecayChannel == 21){
798  // for the D*, setting the third element of the array of the track cuts to those for the soft pion
799  for (Int_t iProng = 0; iProng<cfVtxHF->GetNProngs()-1; iProng++){
800  trackCuts[iProng]=fCuts->GetTrackCuts();
801  }
802  trackCuts[2] = fCuts->GetTrackCutsSoftPi();
803  }
804  else if (fDecayChannel == 22) {
805  // for the Lc->V0+bachelor, setting the second and third elements of the array of the track cuts to those for the V0 daughters
806  trackCuts[0]=fCuts->GetTrackCuts();
807  trackCuts[1]=fCuts->GetTrackCutsV0daughters();
808  trackCuts[2]=fCuts->GetTrackCutsV0daughters();
809  }
810  else {
811  for (Int_t iProng = 0; iProng<cfVtxHF->GetNProngs(); iProng++){
812  trackCuts[iProng]=fCuts->GetTrackCuts();
813  }
814  }
815 
816  //General settings: vertex, feed down and fill reco container with generated values.
817  cfVtxHF->SetRecoPrimVertex(zPrimVertex);
818  cfVtxHF->SetMCPrimaryVertex(zMCVertex);
820  cfVtxHF->SetNVar(fNvar);
824 
825  // switch-off the trigger class selection (doesn't work for MC)
826  fCuts->SetTriggerClass("");
827 
828  // MC vertex, to be used, in case, for pp
830 
831  if (fCentralitySelection){ // keep only the requested centrality
832 
833  if(fCuts->IsEventSelectedInCentrality(aodEvent)!=0) {
834  delete[] containerInput;
835  delete[] containerInputMC;
836  delete [] trackCuts;
837  delete cfVtxHF;
838  return;
839  }
840  } else { // keep all centralities
841  fCuts->SetMinCentrality(0.);
842  fCuts->SetMaxCentrality(100.);
843  }
844 
845  Float_t centValue = 0.;
846  if(!fIsPPData) centValue = fCuts->GetCentrality(aodEvent);
847  cfVtxHF->SetCentralityValue(centValue);
848 
849  // multiplicity estimator with VZERO
850  Double_t vzeroMult=0;
851  AliAODVZERO *vzeroAOD = (AliAODVZERO*)aodEvent->GetVZEROData();
852  if(vzeroAOD) vzeroMult = vzeroAOD->GetMTotV0A() + vzeroAOD->GetMTotV0C();
853 
854  Double_t multiplicity = nTracklets; // set to the Ntracklet estimator
855  if(fMultiplicityEstimator==kVZERO) { multiplicity = vzeroMult; }
856 
857  cfVtxHF->SetMultiplicity(multiplicity);
858 
859  // printf("Multiplicity estimator %d, value %2.2f\n",fMultiplicityEstimator,multiplicity);
860 
861  for (Int_t iPart=0; iPart<mcArray->GetEntriesFast(); iPart++) {
862  AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(mcArray->At(iPart));
863  if (!mcPart){
864  AliError("Failed casting particle from MC array!, Skipping particle");
865  continue;
866  }
867 
868  // Obtain B0 pt spectrum
869  if ((((mcPart->GetPdgCode()%1000)/100==5||(mcPart->GetPdgCode()%1000)/100==-5)) &&
870  (mcPart->Y()<1.)&&(mcPart->Y()>-1.)) {
871  fhBptCutVar->Fill(mcPart->Pt());
872  }
873 
874  //counting c quarks
875  cquarks += cfVtxHF->MCcquarkCounting(mcPart);
876 
877  // check the MC-level cuts, must be the desidered particle
878  if (!fCFManager->CheckParticleCuts(0, mcPart)) {
879  AliDebug(2,"Check the MC-level cuts - not desidered particle");
880  continue; // 0 stands for MC level
881  }
882  else {
883  AliDebug(3, Form("\n\n---> COOL! we found a particle (particle %d)!!! with PDG code = %d \n\n", iPart, mcPart->GetPdgCode()));
884  }
885  cfVtxHF->SetMCCandidateParam(iPart);
886 
887 
888  if (!(cfVtxHF->SetLabelArray())){
889  AliDebug(2,Form("Impossible to set the label array for particle %d (decaychannel = %d)", iPart, fDecayChannel));
890  continue;
891  }
892 
893  //check the candiate family at MC level
894  if (!(cfVtxHF->CheckMCPartFamily(mcPart, mcArray))) {
895  AliDebug(2,Form("Check on the family wrong for particle %d!!! (decaychannel = %d)", iPart, fDecayChannel));
896  continue;
897  }
898  else{
899  AliDebug(2,Form("Check on the family OK for particle %d!!! (decaychannel = %d)", iPart, fDecayChannel));
900  }
901 
902  //Fill the MC container
903  Bool_t mcContainerFilled = cfVtxHF -> FillMCContainer(containerInputMC);
904  AliDebug(2, Form("particle = %d mcContainerFilled = %d", iPart, mcContainerFilled));
905  if (mcContainerFilled) {
906  if (fUseWeight){
907  if (fHistoPtWeight) { // using an histogram as weight function
908  AliDebug(2,"Using Histogram as Pt weight function");
909  fWeight = eventWeight*GetPtWeightFromHistogram(containerInputMC[0]);
910  }
911  else if (fFuncWeight){ // using user-defined function
912  AliDebug(2,"Using function");
913  fWeight = eventWeight*fFuncWeight->Eval(containerInputMC[0]);
914  }
915  else{ // using FONLL
916  AliDebug(2,"Using FONLL");
917  fWeight = eventWeight*GetWeight(containerInputMC[0]);
918  }
919  AliDebug(2,Form("pt = %f, weight = %f",containerInputMC[0], fWeight));
920  }
921  if (!fCuts->IsInFiducialAcceptance(containerInputMC[0],containerInputMC[1])) {
922  AliDebug(3, Form("Not in limited acceptance, containerInputMC[0] = %f, containerInputMC[1] = %f", containerInputMC[0], containerInputMC[1]));
923  continue;
924  }
925  else{
926  AliDebug(3, Form("YES!! in limited acceptance, containerInputMC[0] = %f, containerInputMC[1] = %f", containerInputMC[0],containerInputMC[1]));
927  }
928 
929  //MC Limited Acceptance
930  if (TMath::Abs(containerInputMC[1]) < 0.5) {
931  fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepGeneratedLimAcc, fWeight);
932  AliDebug(3,"MC Lim Acc container filled\n");
933  }
934 
935  //MC
936  fCFManager->GetParticleContainer()->Fill(containerInputMC, kStepGenerated, fWeight);
937  icountMC++;
938  AliDebug(3,"MC container filled \n");
939 
940  // MC in acceptance
941  // check the MC-Acceptance level cuts
942  // since standard CF functions are not applicable, using Kine Cuts on daughters
943  Bool_t mcAccepStep = cfVtxHF-> MCAcceptanceStep();
944  if (mcAccepStep){
945  fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepAcceptance, fWeight);
946  AliDebug(3,"MC acceptance cut passed\n");
947  icountAcc++;
948  if(fhPtCutVar){fObjSpr->FillGenStep(mcPart, mcPart->Pt(), fWeight, mcArray);}
949 
950  //MC Vertex step
951  if (fCuts->IsEventSelected(aodEvent)){
952  // filling the container if the vertex is ok
953  fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepVertex, fWeight) ;
954  AliDebug(3,"Vertex cut passed and container filled\n");
955  icountVertex++;
956 
957  //mc Refit requirement
958  Bool_t mcRefitStep = cfVtxHF->MCRefitStep(aodEvent, &trackCuts[0]);
959  if (mcRefitStep){
960  fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepRefit, fWeight);
961  AliDebug(3,"MC Refit cut passed and container filled\n");
962  icountRefit++;
963  }
964  else{
965  AliDebug(3,"MC Refit cut not passed\n");
966  continue;
967  }
968  }
969  else{
970  AliDebug (3, "MC vertex step not passed\n");
971  continue;
972  }
973  }
974  else{
975  AliDebug (3, "MC in acceptance step not passed\n");
976  continue;
977  }
978  }
979  else {
980  AliDebug (3, "MC container not filled\n");
981  }
982  }
983 
984  if (cquarks<2) AliDebug(2,Form("Event with %d c-quarks", cquarks));
985  AliDebug(2,Form("Found %i MC particles that are %s!!",icountMC,fPartName.Data()));
986  AliDebug(2,Form("Found %i MC particles that are %s and satisfy Acc cuts!!",icountAcc,fPartName.Data()));
987  AliDebug(2,Form("Found %i MC particles that are %s and satisfy Vertex cuts!!",icountVertex,fPartName.Data()));
988  AliDebug(2,Form("Found %i MC particles that are %s and satisfy Refit cuts!!",icountRefit,fPartName.Data()));
989 
990  // Now go to rec level
991  fCountMC += icountMC;
992  fCountAcc += icountAcc;
993  fCountVertex+= icountVertex;
994  fCountRefit+= icountRefit;
995 
996  AliDebug(2,Form("Found %d vertices for decay channel %d",arrayBranch->GetEntriesFast(),fDecayChannel));
997 
998  for(Int_t iCandid = 0; iCandid<arrayBranch->GetEntriesFast();iCandid++){
999  AliAODRecoDecayHF* charmCandidate=0x0;
1000  switch (fDecayChannel){
1001  case 2:{
1002  charmCandidate = (AliAODRecoDecayHF2Prong*)arrayBranch->At(iCandid);
1003  break;
1004  }
1005  case 21:
1006  case 22:{
1007  charmCandidate = (AliAODRecoCascadeHF*)arrayBranch->At(iCandid);
1008  break;
1009  }
1010  case 31:
1011  case 32:
1012  case 33:{
1013  charmCandidate = (AliAODRecoDecayHF3Prong*)arrayBranch->At(iCandid);
1014  break;
1015  }
1016  case 4:{
1017  charmCandidate = (AliAODRecoDecayHF4Prong*)arrayBranch->At(iCandid);
1018  break;
1019  }
1020  default:
1021  break;
1022  }
1023 
1024  Bool_t unsetvtx=kFALSE;
1025  if(!charmCandidate->GetOwnPrimaryVtx()) {
1026  charmCandidate->SetOwnPrimaryVtx(aodVtx); // needed to compute all variables
1027  unsetvtx=kTRUE;
1028  }
1029 
1030  Bool_t signAssociation = cfVtxHF->SetRecoCandidateParam((AliAODRecoDecayHF*)charmCandidate);
1031  if (!signAssociation){
1032  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1033  continue;
1034  }
1035 
1036  Int_t isPartOrAntipart = cfVtxHF->CheckReflexion(fSign);
1037  if (isPartOrAntipart == 0){
1038  AliDebug(2, Form("The candidate pdg code doesn't match the requirement set in the task (fSign = %d)",fSign));
1039  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1040  continue;
1041  }
1042 
1043  AliDebug(3,Form("iCandid=%d - signAssociation=%d, isPartOrAntipart=%d",iCandid, signAssociation, isPartOrAntipart));
1044 
1045  Bool_t recoContFilled = cfVtxHF->FillRecoContainer(containerInput);
1046  AliDebug(3, Form("CF task: RecoContFilled for candidate %d is %d", iCandid, (Int_t)recoContFilled));
1047  if (recoContFilled){
1048 
1049  // weight according to pt
1050  if (fUseWeight){
1051  if (fHistoPtWeight) {
1052  AliDebug(2,"Using Histogram as Pt weight function");
1053  fWeight = eventWeight*GetPtWeightFromHistogram(containerInput[0]);
1054  }
1055  else if (fFuncWeight){ // using user-defined function
1056  AliDebug(2, "Using function");
1057  fWeight = eventWeight*fFuncWeight->Eval(containerInput[0]);
1058  }
1059  else{ // using FONLL
1060  AliDebug(2, "Using FONLL");
1061  fWeight = eventWeight*GetWeight(containerInput[0]);
1062  }
1063  AliDebug(2, Form("pt = %f, weight = %f",containerInput[0], fWeight));
1064  }
1065 
1066  if (!fCuts->IsInFiducialAcceptance(containerInput[0],containerInput[1])){
1067  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1068  continue;
1069  }
1070  //Reco Step
1071  Bool_t recoStep = cfVtxHF->RecoStep();
1072  if (recoStep) AliDebug(2, Form("particle = %d --> CF task: Reco step for candidate %d is %d", iCandid, iCandid, (Int_t)recoStep));
1073  Bool_t vtxCheck = fCuts->IsEventSelected(aodEvent);
1074 
1075 
1076  // Selection on the filtering bit
1077  Bool_t isBitSelected = kTRUE;
1078  if(fDecayChannel==2) {
1079  if(fUseSelectionBit && !charmCandidate->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts)) isBitSelected = kFALSE;
1080  }else if(fDecayChannel==31){
1081  if(fUseSelectionBit && !charmCandidate->HasSelectionBit(AliRDHFCuts::kDplusCuts)) isBitSelected = kFALSE;
1082  }else if(fDecayChannel==33){
1083  if(fUseSelectionBit && !charmCandidate->HasSelectionBit(AliRDHFCuts::kDsCuts)) isBitSelected = kFALSE;
1084  }else if(fDecayChannel==32){
1085  if(fUseSelectionBit && !charmCandidate->HasSelectionBit(AliRDHFCuts::kLcCuts)) isBitSelected = kFALSE;
1086  }
1087  if(!isBitSelected){
1088  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1089  continue;
1090  }
1091 
1092 
1093  if (recoStep && recoContFilled && vtxCheck){
1094  fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructed, fWeight) ;
1095  icountReco++;
1096  AliDebug(3,"Reco step passed and container filled\n");
1097 
1098  //Reco in the acceptance -- take care of UNFOLDING!!!!
1099  Bool_t recoAcceptanceStep = cfVtxHF->RecoAcceptStep(&trackCuts[0]);
1100  if (recoAcceptanceStep) {
1101  fCFManager->GetParticleContainer()->Fill(containerInput,kStepRecoAcceptance, fWeight) ;
1102  icountRecoAcc++;
1103  AliDebug(3,"Reco acceptance cut passed and container filled\n");
1104 
1105  if(fAcceptanceUnf){
1106  Double_t fill[4]; //fill response matrix
1107  Bool_t bUnfolding = cfVtxHF -> FillUnfoldingMatrix(fPDGcode,fill);
1108  if (bUnfolding) fCorrelation->Fill(fill);
1109  }
1110 
1111  //Number of ITS cluster requirements
1112  Int_t recoITSnCluster = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kTracks);
1113  if (recoITSnCluster){
1114  fCFManager->GetParticleContainer()->Fill(containerInput,kStepRecoITSClusters, fWeight) ;
1115  icountRecoITSClusters++;
1116  AliDebug(3,"Reco n ITS cluster cut passed and container filled\n");
1117 
1118  Bool_t iscutsusingpid = fCuts->GetIsUsePID();
1119  Int_t recoAnalysisCuts = -1, recoPidSelection = -1;
1120  fCuts->SetUsePID(kFALSE);
1121  recoAnalysisCuts = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kCandidate, aodEvent);
1122 
1123  if (fDecayChannel==33){ // Ds case, where more possibilities are considered
1124  Bool_t keepDs=ProcessDs(recoAnalysisCuts);
1125  if(keepDs) recoAnalysisCuts=3;
1126  }
1127  else if (fDecayChannel==22){ // Lc->V0+bachelor case, where more possibilities are considered
1128  Bool_t keepLctoV0bachelor = ProcessLctoV0Bachelor(recoAnalysisCuts);
1129  if (keepLctoV0bachelor) recoAnalysisCuts=3;
1130  AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
1131  AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
1132  AliPIDResponse* pidResponse = inputHandler->GetPIDResponse();
1133  if (fUseAdditionalCuts){
1134  if (!((AliCFVertexingHFCascade*)cfVtxHF)->CheckAdditionalCuts(pidResponse)) recoAnalysisCuts = -1;
1135  }
1136  }
1137 
1138  fCuts->SetUsePID(iscutsusingpid); //restoring usage of the PID from the cuts object
1139  Bool_t tempAn=(recoAnalysisCuts == 3 || recoAnalysisCuts == isPartOrAntipart);
1140  if (fDecayChannel == 22) tempAn = (recoAnalysisCuts == 3);
1141  if (fDecayChannel == 32) tempAn=(recoAnalysisCuts >0 || recoAnalysisCuts == isPartOrAntipart);
1142 
1143  if (tempAn){
1144  fCFManager->GetParticleContainer()->Fill(containerInput, kStepRecoPPR, fWeight);
1145  icountRecoPPR++;
1146  AliDebug(3,"Reco Analysis cuts passed and container filled \n");
1147  //pid selection
1148  //recoPidSelection = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kPID);
1149  //if((fCuts->CombineSelectionLevels(3,recoAnalysisCuts,recoPidSelection)==isPartOrAntipart)||(fCuts->CombineSelectionLevels(3,recoAnalysisCuts,recoPidSelection)==3)){
1150  recoPidSelection = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kCandidate, aodEvent);
1151 
1152  if (fDecayChannel==33){ // Ds case, where more possibilities are considered
1153  Bool_t keepDs=ProcessDs(recoPidSelection);
1154  if(keepDs) recoPidSelection=3;
1155  } else if (fDecayChannel==22){ // Lc->V0+bachelor case, where more possibilities are considered
1156  Bool_t keepLctoV0bachelor=ProcessLctoV0Bachelor(recoPidSelection);
1157  if (keepLctoV0bachelor) recoPidSelection=3;
1158  }
1159 
1160  Bool_t tempPid=(recoPidSelection == 3 || recoPidSelection == isPartOrAntipart);
1161  if (fDecayChannel == 22) tempPid = (recoPidSelection == 3);
1162  if (fDecayChannel == 32) tempPid=(recoPidSelection >0 || recoPidSelection == isPartOrAntipart);
1163 
1164  if (tempPid){
1165  Double_t weigPID = 1.;
1166  if (fDecayChannel == 2){ // D0 with Bayesian PID using weights
1167  if(((AliRDHFCutsD0toKpi*)fCuts)->GetCombPID() && (((AliRDHFCutsD0toKpi*)fCuts)->GetBayesianStrategy() == AliRDHFCutsD0toKpi::kBayesWeight || ((AliRDHFCutsD0toKpi*)fCuts)->GetBayesianStrategy() == AliRDHFCutsD0toKpi::kBayesWeightNoFilter)){
1168  if (isPartOrAntipart == 1){
1169  weigPID = ((AliRDHFCutsD0toKpi*)fCuts)->GetWeightsNegative()[AliPID::kKaon] * ((AliRDHFCutsD0toKpi*)fCuts)->GetWeightsPositive()[AliPID::kPion];
1170  }else if (isPartOrAntipart == 2){
1171  weigPID = ((AliRDHFCutsD0toKpi*)fCuts)->GetWeightsPositive()[AliPID::kKaon] * ((AliRDHFCutsD0toKpi*)fCuts)->GetWeightsNegative()[AliPID::kPion];
1172  }
1173  if ((weigPID < 0) || (weigPID > 1)) weigPID = 0.;
1174  }
1175  }else if (fDecayChannel == 33){ // Ds with Bayesian PID using weights
1177  Int_t labDau0=((AliAODTrack*)charmCandidate->GetDaughter(0))->GetLabel();
1178  AliAODMCParticle* firstDau=(AliAODMCParticle*)mcArray->UncheckedAt(TMath::Abs(labDau0));
1179  if(firstDau){
1180  Int_t pdgCode0=TMath::Abs(firstDau->GetPdgCode());
1181  if(pdgCode0==321){
1182  weigPID=((AliRDHFCutsDstoKKpi*)fCuts)->GetWeightForKKpi();
1183  }else if(pdgCode0==211){
1184  weigPID=((AliRDHFCutsDstoKKpi*)fCuts)->GetWeightForpiKK();
1185  }
1186  if ((weigPID < 0) || (weigPID > 1)) weigPID = 0.;
1187  }else{
1188  weigPID=0.;
1189  }
1190  }
1191  }
1192  fCFManager->GetParticleContainer()->Fill(containerInput, kStepRecoPID, fWeight*weigPID);
1193 
1194  AliAODRecoDecayHF2Prong *d0toKpi = (AliAODRecoDecayHF2Prong*)charmCandidate;
1195  fObjSpr->SetFillMC(kTRUE);
1196  fObjSpr->FillSparses(d0toKpi, recoAnalysisCuts, d0toKpi->Pt(),d0toKpi->InvMassD0(), d0toKpi->InvMassD0bar(), fWeight*weigPID, mcArray, aodEvent);
1197 
1198  icountRecoPID++;
1199  AliDebug(3,"Reco PID cuts passed and container filled \n");
1200  if(!fAcceptanceUnf){
1201  Double_t fill[4]; //fill response matrix
1202  Bool_t bUnfolding = cfVtxHF -> FillUnfoldingMatrix(fPDGcode,fill);
1203  if (bUnfolding) fCorrelation->Fill(fill);
1204  }
1205  }
1206  else {
1207  AliDebug(3, "Analysis Cuts step not passed \n");
1208  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1209  continue;
1210  }
1211  }
1212  else {
1213  AliDebug(3, "PID selection not passed \n");
1214  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1215  continue;
1216  }
1217  }
1218  else{
1219  AliDebug(3, "Number of ITS cluster step not passed\n");
1220  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1221  continue;
1222  }
1223  }
1224  else{
1225  AliDebug(3, "Reco acceptance step not passed\n");
1226  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1227  continue;
1228  }
1229  }
1230  else {
1231  AliDebug(3, "Reco step not passed\n");
1232  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1233  continue;
1234  }
1235  }
1236 
1237  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1238  } // end loop on candidate
1239 
1240  fCountReco+= icountReco;
1241  fCountRecoAcc+= icountRecoAcc;
1242  fCountRecoITSClusters+= icountRecoITSClusters;
1243  fCountRecoPPR+= icountRecoPPR;
1244  fCountRecoPID+= icountRecoPID;
1245 
1246  fHistEventsProcessed->Fill(1.5);
1247 
1248  delete[] containerInput;
1249  delete[] containerInputMC;
1250  delete cfVtxHF;
1251  if (trackCuts){
1252  // for (Int_t i=0; i<cfVtxHF->GetNProngs(); i++){
1253  // delete [] trackCuts[i];
1254  // }
1255  delete [] trackCuts;
1256  }
1257 
1258 
1259 }
1260 
1261 //___________________________________________________________________________
1263 {
1264  // The Terminate() function is the last function to be called during
1265  // a query. It always runs on the client, it can be used to present
1266  // the results graphically or save the results to file.
1267 
1268  AliAnalysisTaskSE::Terminate();
1269 
1270  AliInfo(Form("Found %i MC particles that are %s in MC, in %d events",fCountMC,fPartName.Data(),fEvents));
1271  AliInfo(Form("Found %i MC particles that are %s in MC and satisfy Acc cuts, in %d events",fCountAcc,fPartName.Data(),fEvents));
1272  AliInfo(Form("Found %i MC particles that are %s in MC and satisfy Acc cuts, and satisfy Vertex requirement in %d events",fCountVertex,fPartName.Data(),fEvents));
1273  AliInfo(Form("Found %i MC particles that are %s in MC and satisfy Acc cuts, and satisfy ITS+TPC refit requirementin %d events",fCountRefit,fPartName.Data(),fEvents));
1274  AliInfo(Form("Found %i reco %s that are decaying in %s, in %d events",fCountReco,fPartName.Data(),fDauNames.Data(),fEvents));
1275  AliInfo(Form("Among the above, found %i reco %s that are decaying in %s and are in the requested acceptance, in %d events",fCountRecoAcc,fPartName.Data(),fDauNames.Data(),fEvents));
1276  AliInfo(Form("Among the above, found %i reco %s that are decaying in %s and have at least 5 clusters in ITS, in %d events",fCountRecoITSClusters,fPartName.Data(),fDauNames.Data(),fEvents));
1277  AliInfo(Form("Among the above, found %i reco %s that are decaying in %s and satisfy PPR cuts, in %d events",fCountRecoPPR,fPartName.Data(),fDauNames.Data(),fEvents));
1278  AliInfo(Form("Among the above, found %i reco %s that are decaying in %s and satisfy PPR+PID cuts, in %d events",fCountRecoPID,fPartName.Data(),fDauNames.Data(),fEvents));
1279 
1280  // draw some example plots....
1281  AliCFContainer *cont= dynamic_cast<AliCFContainer*> (GetOutputData(2));
1282  if(!cont) {
1283  printf("CONTAINER NOT FOUND\n");
1284  return;
1285  }
1286  // projecting the containers to obtain histograms
1287  // first argument = variable, second argument = step
1288 
1289  TH1D** h = new TH1D*[3];
1290  Int_t nvarToPlot = 0;
1291  if (fConfiguration == kSnail){
1292  //h = new TH1D[3][12];
1293  for (Int_t ih = 0; ih<3; ih++){
1294  if(fDecayChannel==22){
1295  nvarToPlot = 16;
1296  } else {
1297  nvarToPlot = 12;
1298  }
1299  h[ih] = new TH1D[nvarToPlot];
1300  }
1301  for(Int_t iC=1;iC<nvarToPlot; iC++){
1302  // MC-level
1303  h[0][iC] = *(cont->ShowProjection(iC,0));
1304  // MC-Acceptance level
1305  h[1][iC] = *(cont->ShowProjection(iC,1));
1306  // Reco-level
1307  h[2][iC] = *(cont->ShowProjection(iC,4));
1308  }
1309  }
1310  else {
1311  //h = new TH1D[3][12];
1312  nvarToPlot = 8;
1313  for (Int_t ih = 0; ih<3; ih++){
1314  h[ih] = new TH1D[nvarToPlot];
1315  }
1316  for(Int_t iC=0;iC<nvarToPlot; iC++){
1317  // MC-level
1318  h[0][iC] = *(cont->ShowProjection(iC,0));
1319  // MC-Acceptance level
1320  h[1][iC] = *(cont->ShowProjection(iC,1));
1321  // Reco-level
1322  h[2][iC] = *(cont->ShowProjection(iC,4));
1323  }
1324  }
1325  TString* titles;
1326  //Int_t nvarToPlot = 0;
1327  if (fConfiguration == kSnail){
1328  if(fDecayChannel==31){
1329  //nvarToPlot = 12;
1330  titles = new TString[nvarToPlot];
1331  titles[0]="pT_Dplus (GeV/c)";
1332  titles[1]="rapidity";
1333  titles[2]="phi (rad)";
1334  titles[3]="cT (#mum)";
1335  titles[4]="cosPointingAngle";
1336  titles[5]="pT_1 (GeV/c)";
1337  titles[6]="pT_2 (GeV/c)";
1338  titles[7]="pT_3 (GeV/c)";
1339  titles[8]="d0_1 (#mum)";
1340  titles[9]="d0_2 (#mum)";
1341  titles[10]="d0_3 (#mum)";
1342  titles[11]="zVertex (cm)";
1343  } else if (fDecayChannel==22) {
1344  //nvarToPlot = 16;
1345  titles = new TString[nvarToPlot];
1346  titles[0]="p_{T}(#Lambda_{c}) [GeV/c]";
1347  titles[1]="y(#Lambda_{c})";
1348  titles[2]="#varphi(#Lambda_{c}) [rad]";
1349  titles[3]="onTheFlyStatusV0";
1350  titles[4]="z_{vtx} [cm]";
1351  titles[5]="centrality";
1352  titles[6]="fake";
1353  titles[7]="multiplicity";
1354  //titles[8]="pT(bachelor) [GeV/c]";
1355  titles[8]="p(bachelor) [GeV/c]";
1356  titles[9]="p_{T}(V0) [GeV/c]";
1357  titles[10]="y(V0)";
1358  titles[11]="#varphi(V0) [rad]";
1359  titles[12]="m_{inv}(#pi^{+}#pi^{+}) [GeV/c^{2}]";
1360  titles[13]="dcaV0 (nSigma)";
1361  titles[14]="cosine pointing angle (V0)";
1362  titles[15]="cosine pointing angle (#Lambda_{c})";
1363  //titles[16]="c#tauV0 (#mum)";
1364  //titles[17]="c#tau (#mum)";
1365  } else {
1366  //nvarToPlot = 12;
1367  titles = new TString[nvarToPlot];
1368  titles[0]="pT_D0 (GeV/c)";
1369  titles[1]="rapidity";
1370  titles[2]="cosThetaStar";
1371  titles[3]="pT_pi (GeV/c)";
1372  titles[4]="pT_K (Gev/c)";
1373  titles[5]="cT (#mum)";
1374  titles[6]="dca (#mum)";
1375  titles[7]="d0_pi (#mum)";
1376  titles[8]="d0_K (#mum)";
1377  titles[9]="d0xd0 (#mum^2)";
1378  titles[10]="cosPointingAngle";
1379  titles[11]="phi (rad)";
1380  }
1381  }
1382  else {
1383  //nvarToPlot = 8;
1384  titles = new TString[nvarToPlot];
1385  if (fDecayChannel==22) {
1386  titles[0]="p_{T}(#Lambda_{c}) [GeV/c]";
1387  titles[1]="y(#Lambda_{c})";
1388  titles[2]="#varphi(#Lambda_{c}) [rad]";
1389  titles[3]="onTheFlyStatusV0";
1390  titles[4]="z_{vtx} [cm]";
1391  titles[5]="centrality";
1392  titles[6]="fake";
1393  titles[7]="multiplicity";
1394  } else {
1395  titles[0]="pT_candidate (GeV/c)";
1396  titles[1]="rapidity";
1397  titles[2]="cT (#mum)";
1398  titles[3]="phi";
1399  titles[4]="z_{vtx}";
1400  titles[5]="centrality";
1401  titles[6]="fake";
1402  titles[7]="multiplicity";
1403  }
1404  }
1405 
1406  Int_t markers[16]={20,24,21,25,27,28,
1407  20,24,21,25,27,28,
1408  20,24,21,25};
1409  Int_t colors[3]={2,8,4};
1410  for(Int_t iC=0;iC<nvarToPlot; iC++){
1411  for(Int_t iStep=0;iStep<3;iStep++){
1412  h[iStep][iC].SetTitle(titles[iC].Data());
1413  h[iStep][iC].GetXaxis()->SetTitle(titles[iC].Data());
1414  Double_t maxh=h[iStep][iC].GetMaximum();
1415  h[iStep][iC].GetYaxis()->SetRangeUser(0,maxh*1.2);
1416  h[iStep][iC].SetMarkerStyle(markers[iC]);
1417  h[iStep][iC].SetMarkerColor(colors[iStep]);
1418  }
1419  }
1420 
1421  gStyle->SetCanvasColor(0);
1422  gStyle->SetFrameFillColor(0);
1423  gStyle->SetTitleFillColor(0);
1424  gStyle->SetStatColor(0);
1425 
1426  // drawing in 2 separate canvas for a matter of clearity
1427  TCanvas * c1 =new TCanvas(Form("c1New_%d",fDecayChannel),"Vars 0, 1, 2, 3",1100,1200);
1428  c1->Divide(3,4);
1429  Int_t iPad=1;
1430  for(Int_t iVar=0; iVar<4; iVar++){
1431  c1->cd(iPad++);
1432  h[0][iVar].DrawCopy("p");
1433  c1->cd(iPad++);
1434  h[1][iVar].DrawCopy("p");
1435  c1->cd(iPad++);
1436  h[2][iVar].DrawCopy("p");
1437  }
1438 
1439  TCanvas * c2 =new TCanvas(Form("c2New_%d",fDecayChannel),"Vars 4, 5, 6, 7",1100,1200);
1440  c2->Divide(3,4);
1441  iPad=1;
1442  for(Int_t iVar=4; iVar<8; iVar++){
1443  c2->cd(iPad++);
1444  h[0][iVar].DrawCopy("p");
1445  c2->cd(iPad++);
1446  h[1][iVar].DrawCopy("p");
1447  c2->cd(iPad++);
1448  h[2][iVar].DrawCopy("p");
1449  }
1450 
1451  if (fConfiguration == kSnail){
1452  TCanvas * c3 =new TCanvas(Form("c3New_%d",fDecayChannel),"Vars 8, 9, 10, 11",1100,1200);
1453  c3->Divide(3,4);
1454  iPad=1;
1455  for(Int_t iVar=8; iVar<12; iVar++){
1456  c3->cd(iPad++);
1457  h[0][iVar].DrawCopy("p");
1458  c3->cd(iPad++);
1459  h[1][iVar].DrawCopy("p");
1460  c3->cd(iPad++);
1461  h[2][iVar].DrawCopy("p");
1462  }
1463  if (fDecayChannel==22) {
1464  TCanvas * c4 =new TCanvas(Form("c4New_%d",fDecayChannel),"Vars 12, 13, 14, 15",1100,1200);
1465  c4->Divide(3,4);
1466  iPad=1;
1467  for(Int_t iVar=12; iVar<16; iVar++){
1468  c4->cd(iPad++);
1469  h[0][iVar].DrawCopy("p");
1470  c4->cd(iPad++);
1471  h[1][iVar].DrawCopy("p");
1472  c4->cd(iPad++);
1473  h[2][iVar].DrawCopy("p");
1474  }
1475  }
1476  }
1477 
1478  /*
1479  THnSparseD* hcorr = dynamic_cast<THnSparseD*> (GetOutputData(3));
1480 
1481  TH2D* corr1 = hcorr->Projection(0,2);
1482  TH2D* corr2 = hcorr->Projection(1,3);
1483 
1484  TCanvas * c7 =new TCanvas(Form("c7New_%d",fDecayChannel),"",800,400);
1485  c7->Divide(2,1);
1486  c7->cd(1);
1487  corr1->DrawCopy("text");
1488  c7->cd(2);
1489  corr2->DrawCopy("text");
1490  */
1491  TFile* file_projection = new TFile("CFtaskHFprojectionNew.root","RECREATE");
1492 
1493  // corr1->Write();
1494  // corr2->Write();
1495 
1496  for(Int_t iC=0;iC<nvarToPlot; iC++){
1497  for(Int_t iStep=0;iStep<3;iStep++){
1498  h[iStep][iC].Write(Form("Step%d_%s",iStep,titles[iC].Data()));
1499  }
1500  }
1501  file_projection->Close();
1502  for (Int_t ih = 0; ih<3; ih++) delete [] h[ih];
1503  delete [] h;
1504  delete [] titles;
1505 
1506 }
1507 
1508 //___________________________________________________________________________
1510 {
1511  //HERE ONE CAN CREATE OUTPUT OBJECTS, IN PARTICULAR IF THE OBJECT PARAMETERS DON'T NEED
1512  //TO BE SET BEFORE THE EXECUTION OF THE TASK
1513  //
1514  Info("UserCreateOutputObjects","CreateOutputObjects of task %s\n", GetName());
1515 
1516  //slot #1
1517  OpenFile(1);
1518  const char* nameoutput=GetOutputSlot(1)->GetContainer()->GetName();
1519  fHistEventsProcessed = new TH1I(nameoutput,"",2,0,2) ;
1520  fHistEventsProcessed->GetXaxis()->SetBinLabel(1,"Events processed (all)");
1521  fHistEventsProcessed->GetXaxis()->SetBinLabel(2,"Events analyzed (after selection)");
1522 
1523  if(!fObjSpr){
1524  fObjSpr=new AliHFsubtractBFDcuts("feedDownCuts","feedDownCuts");
1525  fObjSpr->InitHistos();
1530  }
1531 
1532  fhBptCutVar = new TH1F("hBptCutVar", "B meson #it{p}_{T} spectrum;#it{p}_{T};Counts (a.u.)",201,0.,50.25);
1533 
1534  PostData( 1,fHistEventsProcessed) ;
1535  PostData( 2,fCFManager->GetParticleContainer()) ;
1536  PostData( 3,fCorrelation) ;
1537  PostData( 6,fhSparseCutVar);
1538  PostData( 7,fhPtCutVar);
1539  PostData( 8,fhBptCutVar);
1540  PostData( 9,fListBdecays);
1541  PostData(10,fQAHists);
1542 }
1543 
1544 
1545 //_________________________________________________________________________
1547  // ad-hoc weight function from ratio of
1548  // pt spectra from FONLL 2.76 TeV and
1549  // pt spectra from MC production LHC12a17a (PYTHIA Perugia0 with pthard bins)
1550  if(fFuncWeight) delete fFuncWeight;
1551  fFuncWeight=new TF1("funcWeight","[0]+[1]*TMath::Exp(-[2]*x)",0.,50.);
1552  fFuncWeight->SetParameter(0,4.63891e-02);
1553  fFuncWeight->SetParameter(1,1.51674e+01);
1554  fFuncWeight->SetParameter(2,4.09941e-01);
1555  fUseWeight=kTRUE;
1556 }
1557 //_________________________________________________________________________
1559  // ad-hoc weight function from ratio of
1560  // D0 pt spectra in PbPb 2011 0-10% centrality and
1561  // pt spectra from MC production LHC12a17a (PYTHIA Perugia0 with pthard bins)
1562  if(fFuncWeight) delete fFuncWeight;
1563  fFuncWeight=new TF1("funcWeight","[0]+[1]/TMath::Power(x,[2])",0.05,50.);
1564  fFuncWeight->SetParameter(0,1.43116e-02);
1565  fFuncWeight->SetParameter(1,4.37758e+02);
1566  fFuncWeight->SetParameter(2,3.08583);
1567  fUseWeight=kTRUE;
1568 }
1569 
1570 //_________________________________________________________________________
1572  // weight function from the ratio of the LHC12a17b MC
1573  // and FONLL calculations for pp data
1574  if(fFuncWeight) delete fFuncWeight;
1575  fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x)",0.15,50.);
1576  fFuncWeight->SetParameters(1.92381e+01, 5.05055e+00, 1.05314e+01, 2.5, 1.88214e-03, 3.44871e+00, -9.74325e-02, 1.97671e+00, -3.21278e-01);
1577  fUseWeight=kTRUE;
1578 }
1579 
1580 //_________________________________________________________________________
1582  // weight function from the ratio of the LHC12a17b MC
1583  // and FONLL calculations for pp data
1584  // corrected by the BAMPS Raa calculation for 30-50% CC
1585  if(fFuncWeight) delete fFuncWeight;
1586  fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x)",0.15,50.);
1587  fFuncWeight->SetParameters(6.10443e+00, 1.53487e+00, 1.99474e+00, 2.5, 5.51172e-03, 5.86590e+00, -5.46963e-01, 9.41201e-02, -1.64323e-01);
1588  fUseWeight=kTRUE;
1589 }
1590 
1591 //_________________________________________________________________________
1593  // weight function from the ratio of the LHC13d3 MC
1594  // and FONLL calculations for pp data
1595  if(fFuncWeight) delete fFuncWeight;
1596  fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x)",0.15,30.);
1597  fFuncWeight->SetParameters(2.94999e+00,3.47032e+00,2.81278e+00,2.5,1.93370e-02,3.86865e+00,-1.54113e-01,8.86944e-02,2.56267e-02);
1598  fUseWeight=kTRUE;
1599 }
1600 
1601 //_________________________________________________________________________
1603  // weight function from the ratio of the LHC10f6a MC
1604  // and FONLL calculations for pp data
1605  if(fFuncWeight) delete fFuncWeight;
1606  fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x)",0.15,40.);
1607  fFuncWeight->SetParameters(2.41522e+01,4.92146e+00,6.72495e+00,2.5,6.15361e-03,4.78995e+00,-4.29135e-01,3.99421e-01,-1.57220e-02);
1608  fUseWeight=kTRUE;
1609 }
1610 
1611 //_________________________________________________________________________
1613  // weight function from the ratio of the LHC12a12 MC
1614  // and FONLL calculations for pp data
1615  if(fFuncWeight) delete fFuncWeight;
1616  fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x+[9])",0.15,50.);
1617  fFuncWeight->SetParameters(1.31497e+01,3.30503e+00,3.45594e+00,2.5,2.28642e-02,1.42372e+00,2.32892e-04,5.21986e-02,-2.14236e-01,3.86200e+00);
1618  fUseWeight=kTRUE;
1619 }
1620 
1621 //_________________________________________________________________________
1623  // weight function from the ratio of the LHC12a12bis MC
1624  // and FONLL calculations for pp data
1625  if(fFuncWeight) delete fFuncWeight;
1626  fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x+[9])",0.15,50.);
1627  fFuncWeight->SetParameters(6.54519e+00,2.74007e+00,2.48325e+00,2.5,1.61113e-01,-5.32546e-01,-3.75916e-04,2.38189e-01,-2.17561e-01,2.35975e+00);
1628  fUseWeight=kTRUE;
1629 }
1630 
1631 //_________________________________________________________________________
1633  // weight function from the ratio of the LHC13e2fix MC
1634  // and FONLL calculations for pp data
1635  if(fFuncWeight) delete fFuncWeight;
1636  fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x+[9])",0.15,50.);
1637  fFuncWeight->SetParameters(1.85862e+01,2.48171e+00,3.39356e+00,2.5,1.70426e-02,2.28453e+00,-4.57749e-02,5.84585e-02,-3.19719e-01,4.16789e+00);
1638  fUseWeight=kTRUE;
1639 }
1640 
1641 //________________________________________________________________
1643  // weight function from the ratio of the LHC10f6a MC
1644  // and FONLL calculations for pp data
1645  if(fFuncWeight) delete fFuncWeight;
1646  fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x)",0.15,40.);
1647  fFuncWeight->SetParameters(2.77730e+01,4.78942e+00,7.45378e+00,2.5,9.86255e-02,2.30120e+00,-4.16435e-01,3.43770e-01,-2.29380e-02);
1648  fUseWeight=kTRUE;
1649 }
1650 
1651 //________________________________________________________________
1653  // weight function from the ratio of the LHC10f6a MC
1654  // and FONLL calculations for pp data
1655  if(fFuncWeight) delete fFuncWeight;
1656  fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x+[9])",0.15,40.);
1657  fFuncWeight->SetParameters(1.34412e+01,3.20068e+00,5.14481e+00,2.5,7.59405e-04,7.51821e+00,-3.93811e-01,2.16849e-02,-3.37768e-02,2.40308e+00);
1658  fUseWeight=kTRUE;
1659 }
1660 
1661 //_________________________________________________________________________
1663 {
1664  //
1665  // calculating the weight to fill the container
1666  //
1667 
1668  // FNOLL central:
1669  // p0 = 1.63297e-01 --> 0.322643
1670  // p1 = 2.96275e+00
1671  // p2 = 2.30301e+00
1672  // p3 = 2.50000e+00
1673 
1674  // PYTHIA
1675  // p0 = 1.85906e-01 --> 0.36609
1676  // p1 = 1.94635e+00
1677  // p2 = 1.40463e+00
1678  // p3 = 2.50000e+00
1679 
1680  Double_t func1[4] = {0.322643,2.96275,2.30301,2.5};
1681  Double_t func2[4] = {0.36609,1.94635,1.40463,2.5};
1682 
1683  Double_t dndpt_func1 = dNdptFit(pt,func1);
1684  if(fUseFlatPtWeight) dndpt_func1 = 1./30.;
1685  Double_t dndpt_func2 = dNdptFit(pt,func2);
1686  AliDebug(2,Form("pt = %f, FONLL = %f, Pythia = %f, ratio = %f",pt,dndpt_func1,dndpt_func2,dndpt_func1/dndpt_func2));
1687  return dndpt_func1/dndpt_func2;
1688 }
1689 
1690 //__________________________________________________________________________________________________
1691 Double_t AliCFTaskVertexingHFCutVarFDSub::dNdptFit(Float_t pt, Double_t* par)
1692 {
1693  //
1694  // calculating dNdpt
1695  //
1696 
1697  Double_t denom = TMath::Power((pt/par[1]), par[3] );
1698  Double_t dNdpt = par[0]*pt/TMath::Power(1.+denom, par[2]);
1699 
1700  return dNdpt;
1701 }
1702 
1703 //_________________________________________________________________________
1705 {
1706  //
1707  // Using an histogram as weight function
1708  // weight = 0 in the range outside the function
1709  //
1710  Double_t weight = 0.0;
1711  Int_t histoNbins = fHistoPtWeight->GetNbinsX();
1712  Int_t histobin = fHistoPtWeight->FindBin(pt);
1713  if( (histobin>0) && (histobin<=histoNbins) ) {
1714  weight = fHistoPtWeight->GetBinContent(histobin);
1715  }
1716 
1717  return weight;
1718 }
1719 
1720 //__________________________________________________________________________________________________
1721 Double_t AliCFTaskVertexingHFCutVarFDSub::GetZWeight(Float_t z, Int_t runnumber){
1722  //
1723  // calculating the z-vtx weight for the given run range
1724  //
1725 
1726  if(runnumber>146824 || runnumber<146803) return 1.0;
1727 
1728  Double_t func1[3] = {1.0, -0.5, 6.5 };
1729  Double_t func2[3] = {1.0, -0.5, 5.5 };
1730 
1731  Double_t dzFunc1 = DodzFit(z,func1);
1732  Double_t dzFunc2 = DodzFit(z,func2);
1733 
1734  return dzFunc1/dzFunc2;
1735 }
1736 
1737 //__________________________________________________________________________________________________
1738 Double_t AliCFTaskVertexingHFCutVarFDSub::DodzFit(Float_t z, Double_t* par) {
1739 
1740  //
1741  // Gaussian z-vtx shape
1742  //
1743  //gaussian = [0]/TMath::Sqrt(2.*TMath::Pi())/[2]*exp[-(x-[1])*(x-[1])/(2*[2]*[2])]
1744 
1745  Double_t value = par[0]/TMath::Sqrt(2.*TMath::Pi())/par[2]*TMath::Exp(-(z-par[1])*(z-par[1])/2./par[2]/par[2]);
1746 
1747  return value;
1748 }
1749 //__________________________________________________________________________________________________
1751  //
1752  // calculates the Nch weight using the measured and generateed Nch distributions
1753  //
1754  if(nch<=0) return 0.;
1755  if(!fHistoMeasNch || !fHistoMCNch) { AliError("Input histos to evaluate Nch weights missing"); return 0.; }
1756  Double_t pMeas=fHistoMeasNch->GetBinContent(fHistoMeasNch->FindBin(nch));
1757  Double_t pMC=fHistoMCNch->GetBinContent(fHistoMCNch->FindBin(nch));
1758  Double_t weight = pMC>0 ? pMeas/pMC : 0.;
1759  if(fUseMultRatioAsWeight) weight = pMC;
1760  return weight;
1761 }
1762 //__________________________________________________________________________________________________
1764  // creates historgam with measured multiplcity distribution in pp 7 TeV collisions (from Eur. Phys. J. C (2010) 68: 345–354)
1765  //
1766  // for Nch > 70 the points were obtained with a double NBD distribution fit
1767  // 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
1768  // fit1->SetParameter(1,1.63); // k parameter
1769  // fit1->SetParameter(2,12.8); // mean multiplicity
1770  //
1771  Double_t nchbins[82]={0.50,1.50,2.50,3.50,4.50,5.50,6.50,7.50,8.50,9.50,
1772  10.50,11.50,12.50,13.50,14.50,15.50,16.50,17.50,18.50,19.50,
1773  20.50,21.50,22.50,23.50,24.50,25.50,26.50,27.50,28.50,29.50,
1774  30.50,31.50,32.50,33.50,34.50,35.50,36.50,37.50,38.50,39.50,
1775  40.50,41.50,42.50,43.50,44.50,45.50,46.50,47.50,48.50,49.50,
1776  50.50,51.50,52.50,53.50,54.50,55.50,56.50,57.50,58.50,59.50,
1777  60.50,62.50,64.50,66.50,68.50,70.50,72.50,74.50,76.50,78.50,
1778  80.50,82.50,84.50,86.50,88.50,90.50,92.50,94.50,96.50,98.50,
1779  100.50,102.50};
1780  Double_t pch[81]={0.062011,0.072943,0.070771,0.067245,0.062834,0.057383,0.051499,0.04591,0.041109,0.036954,
1781  0.03359,0.030729,0.028539,0.026575,0.024653,0.0229,0.021325,0.019768,0.018561,0.017187,
1782  0.01604,0.014836,0.013726,0.012576,0.011481,0.010393,0.009502,0.008776,0.008024,0.007452,
1783  0.006851,0.006428,0.00594,0.005515,0.005102,0.00469,0.004162,0.003811,0.003389,0.003071,
1784  0.002708,0.002422,0.002184,0.001968,0.00186,0.00165,0.001577,0.001387,0.001254,0.001118,
1785  0.001037,0.000942,0.000823,0.000736,0.000654,0.000579,0.000512,0.00049,0.00045,0.000355,
1786  0.000296,0.000265,0.000193,0.00016,0.000126,0.0000851, 0.0000676,0.0000537,0.0000426, 0.0000338,
1787  0.0000268,0.0000213,0.0000166,0.0000133,0.0000106,0.00000837,0.00000662, 0.00000524,0.00000414, 0.00000327,
1788  0.00000258};
1789 
1790  if(fHistoMeasNch) delete fHistoMeasNch;
1791  fHistoMeasNch=new TH1F("hMeaseNch","",81,nchbins);
1792  for(Int_t i=0; i<81; i++){
1793  fHistoMeasNch->SetBinContent(i+1,pch[i]);
1794  fHistoMeasNch->SetBinError(i+1,0.);
1795  }
1796 }
1797 
1798 //__________________________________________________________________________________________________
1799 Bool_t AliCFTaskVertexingHFCutVarFDSub::ProcessDs(Int_t recoAnalysisCuts) const{
1800  // processes output of Ds is selected
1801  Bool_t keep=kFALSE;
1802  if(recoAnalysisCuts > 0){
1803  Int_t isKKpi=recoAnalysisCuts&1;
1804  Int_t ispiKK=recoAnalysisCuts&2;
1805  Int_t isPhiKKpi=recoAnalysisCuts&4;
1806  Int_t isPhipiKK=recoAnalysisCuts&8;
1807  Int_t isK0starKKpi=recoAnalysisCuts&16;
1808  Int_t isK0starpiKK=recoAnalysisCuts&32;
1809  if(fDsOption==1){
1810  if(isKKpi && isPhiKKpi) keep=kTRUE;
1811  if(ispiKK && isPhipiKK) keep=kTRUE;
1812  }
1813  else if(fDsOption==2){
1814  if(isKKpi && isK0starKKpi) keep=kTRUE;
1815  if(ispiKK && isK0starpiKK) keep=kTRUE;
1816  }
1817  else if(fDsOption==3)keep=kTRUE;
1818  }
1819  return keep;
1820 }
1821 //__________________________________________________________________________________________________
1822 Bool_t AliCFTaskVertexingHFCutVarFDSub::ProcessLctoV0Bachelor(Int_t recoAnalysisCuts) const{
1823  // processes output of Lc->V0+bnachelor is selected
1824  Bool_t keep=kFALSE;
1825  if(recoAnalysisCuts > 0){
1826 
1827  Int_t isK0Sp = recoAnalysisCuts&1;
1828  Int_t isLambdaBarpi = recoAnalysisCuts&2;
1829  Int_t isLambdapi = recoAnalysisCuts&4;
1830 
1831  if(fLctoV0bachelorOption==1){
1832  if(isK0Sp) keep=kTRUE;
1833  }
1834  else if(fLctoV0bachelorOption==2){
1835  if(isLambdaBarpi) keep=kTRUE;
1836  }
1837  else if(fLctoV0bachelorOption==4){
1838  if(isLambdapi) keep=kTRUE;
1839  }
1840  else if(fLctoV0bachelorOption==7) {
1841  if (isK0Sp || isLambdaBarpi || isLambdapi) keep=kTRUE;
1842  }
1843  }
1844  return keep;
1845 }
1846 
1847 //____________________________________________________________________________
1849  // Get Estimator Histogram from period event->GetRunNumber();
1850  //
1851  // If you select SPD tracklets in |eta|<1 you should use type == 1
1852  //
1853 
1854  Int_t runNo = event->GetRunNumber();
1855  Int_t period = -1; // pp: 0-LHC10b, 1-LHC10c, 2-LHC10d, 3-LHC10e
1856  // pPb: 0-LHC13b, 1-LHC13c
1857 
1858  if (fIsPPbData) { // setting run numbers for LHC13 if pPb
1859  if (runNo>195343 && runNo<195484) period = 0;
1860  if (runNo>195528 && runNo<195678) period = 1;
1861  if (period<0 || period>1) return 0;
1862  } else { //else assume pp 2010
1863  if(runNo>114930 && runNo<117223) period = 0;
1864  if(runNo>119158 && runNo<120830) period = 1;
1865  if(runNo>122373 && runNo<126438) period = 2;
1866  if(runNo>127711 && runNo<130841) period = 3;
1867  if(period<0 || period>3) return 0;
1868  }
1869 
1870  return fMultEstimatorAvg[period];
1871 }
Int_t fNvar
flag to use directly the ratio of the distributions (fHistoMCNch) instead of computing it ...
void SetCentralityValue(Float_t centValue)
Bool_t fAcceptanceUnf
flag to select D0 origins. 0 Only from charm 1 only from beauty 2 both from charm and beauty ...
virtual AliESDtrackCuts * GetTrackCutsV0daughters() const
Definition: AliRDHFCuts.h:247
void SetFillFromGenerated(Bool_t flag)
Bool_t fRejectIfNoQuark
selection flag for fakes tracks
Int_t fCountRecoAcc
Reco particle found that satisfy cuts.
virtual AliESDtrackCuts * GetTrackCutsSoftPi() const
Definition: AliRDHFCuts.h:246
void SetUseMCVertex()
Definition: AliRDHFCuts.h:336
Bool_t fFillFromGenerated
decay channel to configure the task
Int_t MCcquarkCounting(AliAODMCParticle *mcPart) const
UShort_t fOriginDselection
flag to indicate whether data container should be filled with generated values also for reconstructed...
Bool_t fUseFlatPtWeight
weight used to fill the container
Class for HF corrections as a function of many variables and steps For D* and other cascades...
Int_t IsEventSelectedInCentrality(AliVEvent *event)
Bool_t HasSelectionBit(Int_t i) const
Class for HF corrections as a function of many variables and step.
Bool_t fIsPPData
flag to use the z-vtx corrected (if not use uncorrected) multiplicity estimator
AliHFsubtractBFDcuts * fObjSpr
cut on momentum conservation
Bool_t FillRecoContainer(Double_t *containerInput)
Double_t fWeight
flag to decide whether to use pt-weights != 1 when filling the container or not
Int_t fConfiguration
Ds decay option (generation level)
THnSparseF * GetSparseMC() const
void SetNVar(Int_t nVar)
Int_t fCountRecoITSClusters
Reco particle found that satisfy cuts in requested acceptance.
Bool_t fIsPPbData
flag for pp data (not checking centrality)
Double_t dNdptFit(Float_t pt, Double_t *par)
void SetRejectCandidateIfNotFromQuark(Bool_t opt)
virtual Bool_t SetLabelArray()
void FillSparses(AliAODRecoDecayHF2Prong *dzeroPart, Int_t isSelected, Double_t pt=-1, Double_t massD0=-1, Double_t massD0bar=-1, Double_t weight=1., TClonesArray *mcArray=0x0, AliAODEvent *aodEvent=0x0)
Bool_t fCentralitySelection
flag to decide wheter to keep D0 only (0), D0bar only (1), or both D0 and D0bar (2) ...
UInt_t fPDGcode
flag to use selection bit
Bool_t fUseMCVertex
flag to remove events not geenrated with PYTHIA
Int_t GetNProngs() const
TH1F * fHistoMCNch
histogram with measured Nch distribution (pp 7 TeV)
Double_t GetMaxVtxZ() const
Definition: AliRDHFCuts.h:241
void SetFillMC(Bool_t fillMC=kTRUE)
Int_t fCountMC
list of profile histos for z-vtx correction
AliCFTaskVertexingHFCutVarFDSub()
multiplicity estimators
TH1F * fHistoMeasNch
user-defined histogram to calculate the Pt weights
TH3F * GetHistoPtMCgen() const
TH1F * fHistoPtWeight
user-defined function to be used to calculate weights
Class for cuts on AOD reconstructed D+->Kpipi.
TF1 * fFuncWeight
configuration (slow / fast) of the CF –> different variables will be allocated (all / reduced number)...
TH1F * fhBptCutVar
Generator level histogram D0 pt, Nprongs of the decay, Mother pt.
static Int_t GetNumberOfTrackletsInEtaRange(AliAODEvent *ev, Double_t mineta, Double_t maxeta)
Int_t fCountReco
Reco particle found that satisfy kTPCrefit and kITSrefit.
void FillGenStep(AliAODMCParticle *dzeroMC, Double_t pt=-1, Double_t weight=1., TClonesArray *mcArray=0x0)
Bool_t ProcessDs(Int_t returnCodeDs) const
void SetMinCentrality(Float_t minCentrality=0.)
Definition: AliRDHFCuts.h:51
Bool_t fUseAdditionalCuts
flag for pPb data (used for multiplicity corrections)
Bool_t fUseMultRatioAsWeight
flag to decide whether to use Ncharged weights != 1 when filling the container or not ...
Int_t fGenLctoV0bachelorOption
Lc->V0+bachelor decay option (selection level)
Bool_t fUseTrackletsWeight
flag to decide whether to use Ncharged weights != 1 when filling the container or not ...
Double_t DodzFit(Float_t z, Double_t *par)
AliAODVertex * GetOwnPrimaryVtx() const
AliESDtrackCuts * GetTrackCuts() const
Definition: AliRDHFCuts.h:245
Float_t GetCentrality(AliAODEvent *aodEvent)
Definition: AliRDHFCuts.h:242
Int_t fLctoV0bachelorOption
resonant deacy channel to be used if the CF should be run on resonant channels only ...
void SetFakeSelection(Int_t fakeSel)
Bool_t fUseSelectionBit
Lc->V0+bachelor decay option (generation level)
Bool_t fZvtxCorrectedNtrkEstimator
refrence multiplcity (period b)
Int_t fCountRefit
Reco particle found that satisfy vertex constrained.
Double_t fRefMult
TProfile with mult vas. Z per period.
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)
void SetMaxCentrality(Float_t maxCentrality=100.)
Definition: AliRDHFCuts.h:52
void SetOwnPrimaryVtx(const AliAODVertex *vtx)
void SetRecoPrimVertex(Double_t zPrimVertex)
Double_t GetZWeight(Float_t z, Int_t runnumber)
TList * fQAHists
List with B hadron decays (unique entries)
void SetMCCandidateParam(Int_t label)
TH1I * fHistEventsProcessed
pointer to the CF manager
Int_t fCountRecoPPR
Reco particle found that satisfy cuts in n. of ITS clusters.
Bool_t CheckMCPartFamily(AliAODMCParticle *, TClonesArray *) const
THnSparseF * fhSparseCutVar
object for cut variation study
Bool_t fUseCutsForTMVA
flag to use additional cuts needed for Lc –> K0S + p, TMVA
virtual Bool_t SetRecoCandidateParam(AliAODRecoDecayHF *)
static Double_t GetCorrectedNtracklets(TProfile *estimatorAvg, Double_t uncorrectedNacc, Double_t vtxZ, Double_t refMult)
Int_t CheckReflexion(Char_t isSign)
Bool_t IsEventSelected(AliVEvent *event)
Bool_t fUseZWeight
flag to decide to use a flat pt shape
void SetUsePID(Bool_t flag=kTRUE)
Definition: AliRDHFCuts.h:204
virtual void PrintAll() const
AliCFTaskVertexingHFCutVarFDSub & operator=(const AliCFTaskVertexingHFCutVarFDSub &c)
Bool_t RecoAcceptStep(AliESDtrackCuts **trackCuts) const
TList * fListProfiles
response matrix for unfolding
static Int_t GetGeneratedPhysicalPrimariesInEtaRange(TClonesArray *arrayMC, Double_t mineta, Double_t maxeta)
Bool_t IsSelected(TObject *obj)
Definition: AliRDHFCuts.h:273
TProfile * fMultEstimatorAvg[4]
Definition of the multiplicity estimator: kNtrk10=0, kNtrk10to16=1, kVZERO=2.
Bool_t MCRefitStep(AliAODEvent *aodEvent, AliESDtrackCuts **trackCuts) const
Int_t fFakeSelection
flag to switch off the centrality selection
Bool_t GetIsUsePID() const
Definition: AliRDHFCuts.h:253
UInt_t fResonantDecay
histogram with Nch distribution from MC production
Float_t fCutOnMomConservation
flag to define which task to use for Lc –> K0S+p
AliRDHFCuts * fCuts
flag for unfolding before or after cuts.
void SetTriggerClass(TString trclass0, TString trclass1="")
Definition: AliRDHFCuts.h:192
Class to compute variables for correction framework // for 3-body decays of D mesons (D+...
virtual Bool_t IsInFiducialAcceptance(Double_t, Double_t) const
Definition: AliRDHFCuts.h:290
TProfile * GetEstimatorHistogram(const AliVEvent *event)
void SetConfiguration(Int_t configuration)
void SetMultiplicity(Double_t multiplicity)
Int_t fGenDsOption
Ds decay option (selection level)
void SetMCPrimaryVertex(Double_t zMCVertex)
Int_t fCountRecoPID
Reco particle found that satisfy cuts in PPR.
Int_t fDsOption
flag to use MC vertex (useful when runnign in pp)
void UserCreateOutputObjects()
ANALYSIS FRAMEWORK STUFF to loop on data and fill output objects.
Bool_t ProcessLctoV0Bachelor(Int_t returnCodeDs) const
TString fPartName
number of variables for the container
Class for HF corrections as a function of many variables and step.
Int_t fCountVertex
MC particle found that satisfy acceptance cuts.