AliPhysics  9b6b435 (9b6b435)
AliCFTaskVertexingHF.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 #include <TCanvas.h>
33 #include <TDatabasePDG.h>
34 #include <TProfile.h>
35 #include <TH1I.h>
36 #include <TStyle.h>
37 #include <TFile.h>
38 #include <TF1.h>
39 
40 #include "AliCFTaskVertexingHF.h"
41 #include "AliMCEvent.h"
42 #include "AliCFManager.h"
43 #include "AliCFContainer.h"
44 #include "AliLog.h"
45 #include "AliInputEventHandler.h"
46 #include "AliAnalysisManager.h"
47 #include "AliAODHandler.h"
48 #include "AliAODEvent.h"
49 #include "AliAODRecoDecay.h"
50 #include "AliAODRecoDecayHF.h"
54 #include "AliAODRecoCascadeHF.h"
55 #include "AliAODMCParticle.h"
56 #include "AliAODMCHeader.h"
57 #include "AliESDtrack.h"
58 #include "TChain.h"
59 #include "THnSparse.h"
60 #include "TH2D.h"
61 #include "AliESDtrackCuts.h"
62 #include "AliRDHFCuts.h"
63 #include "AliRDHFCutsD0toKpi.h"
66 #include "AliRDHFCutsDstoKKpi.h"
67 #include "AliRDHFCutsLctopKpi.h"
68 #include "AliRDHFCutsD0toKpipipi.h"
69 #include "AliRDHFCutsLctoV0.h"
70 #include "AliCFVertexingHF2Prong.h"
71 #include "AliCFVertexingHF3Prong.h"
74 #include "AliCFVertexingHF.h"
75 #include "AliVertexingHFUtils.h"
76 #include "AliAnalysisDataSlot.h"
77 #include "AliAnalysisDataContainer.h"
78 #include "AliAnalysisVertexingHF.h"
79 #include "AliPIDResponse.h"
80 
81 //__________________________________________________________________________
84  fCFManager(0x0),
85  fHistEventsProcessed(0x0),
86  fCorrelation(0x0),
87  fListProfiles(0),
88  fCountMC(0),
89  fCountGenLimAcc(0),
90  fCountGenLimAccNoAcc(0),
91  fCountAcc(0),
92  fCountVertex(0),
93  fCountRefit(0),
94  fCountReco(0),
95  fCountRecoAcc(0),
96  fCountRecoITSClusters(0),
97  fCountRecoPPR(0),
98  fCountRecoPID(0),
99  fEvents(0),
100  fDecayChannel(0),
101  fFillFromGenerated(kFALSE),
102  fOriginDselection(0),
103  fAcceptanceUnf(kTRUE),
104  fCuts(0),
105  fUseWeight(kFALSE),
106  fWeight(1.),
107  fUseFlatPtWeight(kFALSE),
108  fUseZWeight(kFALSE),
109  fUseNchWeight(kFALSE),
110  fUseTrackletsWeight(kFALSE),
111  fUseMultRatioAsWeight(kFALSE),
112  fNvar(0),
113  fPartName(""),
114  fDauNames(""),
115  fSign(2),
116  fCentralitySelection(kTRUE),
117  fFakeSelection(0),
118  fRejectIfNoQuark(kTRUE),
119  fUseMCVertex(kFALSE),
120  fDsOption(1),
121  fGenDsOption(3),
122  fConfiguration(kCheetah), // by default, setting the fast configuration
123  fFuncWeight(0x0),
124  fHistoPtWeight(0x0),
125  fHistoMeasNch(0x0),
126  fHistoMCNch(0x0),
127  fResonantDecay(0),
128  fLctoV0bachelorOption(1),
129  fGenLctoV0bachelorOption(0),
130  fUseSelectionBit(kTRUE),
131  fPDGcode(0),
132  fMultiplicityEstimator(kNtrk10),
133  fRefMult(9.26),
134  fZvtxCorrectedNtrkEstimator(kFALSE),
135  fIsPPData(kFALSE),
136  fIsPPbData(kFALSE),
137  fUseAdditionalCuts(kFALSE),
138  fUseCutsForTMVA(kFALSE),
139  fUseCascadeTaskForLctoV0bachelor(kFALSE),
140  fFillMinimumSteps(kFALSE),
141  fCutOnMomConservation(0.00001)
142 {
143  //
144  //Default ctor
145  //
146  for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=0;
147 }
148 //___________________________________________________________________________
150  AliAnalysisTaskSE(name),
151  fCFManager(0x0),
153  fCorrelation(0x0),
154  fListProfiles(0),
155  fCountMC(0),
156  fCountGenLimAcc(0),
158  fCountAcc(0),
159  fCountVertex(0),
160  fCountRefit(0),
161  fCountReco(0),
162  fCountRecoAcc(0),
164  fCountRecoPPR(0),
165  fCountRecoPID(0),
166  fEvents(0),
167  fDecayChannel(0),
168  fFillFromGenerated(kFALSE),
170  fAcceptanceUnf(kTRUE),
171  fCuts(cuts),
172  fUseWeight(kFALSE),
173  fWeight(1.),
174  fUseFlatPtWeight(kFALSE),
175  fUseZWeight(kFALSE),
176  fUseNchWeight(kFALSE),
177  fUseTrackletsWeight(kFALSE),
178  fUseMultRatioAsWeight(kFALSE),
179  fNvar(0),
180  fPartName(""),
181  fDauNames(""),
182  fSign(2),
183  fCentralitySelection(kTRUE),
184  fFakeSelection(0),
185  fRejectIfNoQuark(kTRUE),
186  fUseMCVertex(kFALSE),
187  fDsOption(1),
188  fGenDsOption(3),
189  fConfiguration(kCheetah), // by default, setting the fast configuration
190  fFuncWeight(func),
191  fHistoPtWeight(0x0),
192  fHistoMeasNch(0x0),
193  fHistoMCNch(0x0),
194  fResonantDecay(0),
197  fUseSelectionBit(kTRUE),
198  fPDGcode(0),
200  fRefMult(9.26),
202  fIsPPData(kFALSE),
203  fIsPPbData(kFALSE),
204  fUseAdditionalCuts(kFALSE),
205  fUseCutsForTMVA(kFALSE),
207  fFillMinimumSteps(kFALSE),
208  fCutOnMomConservation(0.00001)
209 {
210  //
211  // Constructor. Initialization of Inputs and Outputs
212  //
213  /*
214  DefineInput(0) and DefineOutput(0)
215  are taken care of by AliAnalysisTaskSE constructor
216  */
217  DefineOutput(1,TH1I::Class());
218  DefineOutput(2,AliCFContainer::Class());
219  DefineOutput(3,THnSparseD::Class());
220  DefineOutput(4,AliRDHFCuts::Class());
221  for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=0;
222  DefineOutput(5,TList::Class()); // slot #5 keeps the zvtx Ntrakclets correction profiles
223 
224  fCuts->PrintAll();
225 }
226 
227 //___________________________________________________________________________
229 {
230  //
231  // Assignment operator
232  //
233  if (this!=&c) {
234  AliAnalysisTaskSE::operator=(c) ;
237  fCuts = c.fCuts;
242  for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=c.fMultEstimatorAvg[i];
243  }
244  return *this;
245 }
246 
247 //___________________________________________________________________________
254  fCountMC(c.fCountMC),
257  fCountAcc(c.fCountAcc),
265  fEvents(c.fEvents),
270  fCuts(c.fCuts),
272  fWeight(c.fWeight),
278  fNvar(c.fNvar),
279  fPartName(c.fPartName),
280  fDauNames(c.fDauNames),
281  fSign(c.fSign),
286  fDsOption(c.fDsOption),
297  fPDGcode(c.fPDGcode),
299  fRefMult(c.fRefMult),
301  fIsPPData(c.fIsPPData),
308 {
309  //
310  // Copy Constructor
311  //
312  for(Int_t i=0; i<4; i++) fMultEstimatorAvg[i]=c.fMultEstimatorAvg[i];
313 }
314 
315 //___________________________________________________________________________
317 {
318  //
319  //destructor
320  //
321  if (fCFManager) delete fCFManager ;
323  if (fCorrelation) delete fCorrelation ;
324  if (fListProfiles) delete fListProfiles;
325  if (fCuts) delete fCuts;
326  if (fFuncWeight) delete fFuncWeight;
327  if (fHistoPtWeight) delete fHistoPtWeight;
328  if (fHistoMeasNch) delete fHistoMeasNch;
329  if (fHistoMCNch) delete fHistoMCNch;
330  for(Int_t i=0; i<4; i++) { if(fMultEstimatorAvg[i]) delete fMultEstimatorAvg[i]; }
331 }
332 
333 //_________________________________________________________________________-
335 {
336  //
337  // Initialization
338  //
339 
340  if (fDebug>1) printf("AliCFTaskVertexingHF::Init()");
341  if(fUseWeight && fUseZWeight) { AliFatal("Can not use at the same time pt and z-vtx weights, please choose"); return; }
342  if(fUseWeight && fUseNchWeight) { AliInfo("Beware, using at the same time pt and Nch weights, please check"); }
343  if(fUseNchWeight && !fHistoMCNch) { AliFatal("Need to pass the MC Nch distribution to use Nch weights"); return; }
345 
346  AliRDHFCuts *copyfCuts = 0x0;
347  if (!fCuts){
348  AliFatal("No cuts defined - Exiting...");
349  return;
350  }
351 
352  switch (fDecayChannel){
353  case 2:{
354  fPDGcode = 421;
355  copyfCuts = new AliRDHFCutsD0toKpi(*(static_cast<AliRDHFCutsD0toKpi*>(fCuts)));
356  switch (fConfiguration) {
357  case kSnail: // slow configuration: all variables in
358  fNvar = 16;
359  break;
360  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
361  fNvar = 8;
362  break;
363  case kFalcon:// super fast configuration: only pt_candidate, y, centrality
364  fNvar = 4;
365  break;
366  case kESE:// configuration with variables for ESE analysis (pt,y,centrality,mult,q2)
367  fNvar = 6;
368  break;
369  }
370  fPartName="D0";
371  fDauNames="K+pi";
372  break;
373  }
374  case 21:{
375  fPDGcode = 413;
376  copyfCuts = new AliRDHFCutsDStartoKpipi(*(static_cast<AliRDHFCutsDStartoKpipi*>(fCuts)));
377  switch (fConfiguration) {
378  case kSnail: // slow configuration: all variables in
379  fNvar = 16;
380  break;
381  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
382  fNvar = 8;
383  break;
384  case kFalcon:// super fast configuration: only pt_candidate, y, centrality
385  fNvar = 4;
386  break;
387  case kESE:// configuration with variables for ESE analysis (pt,y,centrality,mult,q2)
388  fNvar = 6;
389  break;
390  }
391  fPartName="Dstar";
392  fDauNames="K+pi+pi";
393  break;
394  }
395  case 22:{
396  fPDGcode = 4122;
397  copyfCuts = new AliRDHFCutsLctoV0(*(static_cast<AliRDHFCutsLctoV0*>(fCuts)));
398  switch (fConfiguration) {
399  case kSnail: // slow configuration: all variables in
400  fNvar = 16;
401  break;
402  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
403  fNvar = 8;
404  break;
405  case kFalcon:// super fast configuration: only pt_candidate, y, centrality
406  fNvar = 4;
407  break;
408  }
409  fPartName="Lambdac";
410  fDauNames="V0+bachelor";
411  break;
412  }
413  case 31:{
414  fPDGcode = 411;
415  copyfCuts = new AliRDHFCutsDplustoKpipi(*(static_cast<AliRDHFCutsDplustoKpipi*>(fCuts)));
416  switch (fConfiguration) {
417  case kSnail: // slow configuration: all variables in
418  fNvar = 14;
419  break;
420  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
421  fNvar = 8;
422  break;
423  case kFalcon:// super fast configuration: only pt_candidate, y, centrality
424  fNvar = 4;
425  break;
426  case kESE:// configuration with variables for ESE analysis (pt,y,centrality,mult,q2)
427  fNvar = 6;
428  break;
429  }
430  fPartName="Dplus";
431  fDauNames="K+pi+pi";
432  break;
433  }
434  case 32:{
435  fPDGcode = 4122;
436  copyfCuts = new AliRDHFCutsLctopKpi(*(static_cast<AliRDHFCutsLctopKpi*>(fCuts)));
437  switch (fConfiguration) {
438  case kSnail: // slow configuration: all variables in
439  fNvar = 18;
440  break;
441  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
442  fNvar = 8;
443  break;
444  case kFalcon:// super fast configuration: only pt_candidate, y, centrality
445  fNvar = 4;
446  break;
447  }
448  fPartName="Lambdac";
449  fDauNames="p+K+pi";
450  break;
451  }
452  case 33:{
453  fPDGcode = 431;
454  copyfCuts = new AliRDHFCutsDstoKKpi(*(static_cast<AliRDHFCutsDstoKKpi*>(fCuts)));
455  switch (fConfiguration) {
456  case kSnail: // slow configuration: all variables in
457  fNvar = 14;
458  break;
459  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
460  fNvar = 8;
461  break;
462  case kFalcon:// super fast configuration: only pt_candidate, y, centrality
463  fNvar = 4;
464  break;
465  }
466  fPartName="Ds";
467  fDauNames="K+K+pi";
468  break;
469  }
470  case 4:{
471  fPDGcode = 421;
472  copyfCuts = new AliRDHFCutsD0toKpipipi(*(static_cast<AliRDHFCutsD0toKpipipi*>(fCuts)));
473  switch (fConfiguration) {
474  case kSnail: // slow configuration: all variables in
475  fNvar = 16;
476  break;
477  case kCheetah:// fast configuration: only pt_candidate, y, phi, ct, fake, z_vtx, centrality, multiplicity will be filled
478  fNvar = 8;
479  break;
480  case kFalcon:// super fast configuration: only pt_candidate, y, centrality
481  fNvar = 4;
482  break;
483  }
484  fPartName="D0";
485  fDauNames="K+pi+pi+pi";
486  break;
487  }
488  default:
489  AliFatal("The decay channel MUST be defined according to AliCFVertexing::DecayChannel - Exiting...");
490  break;
491  }
492 
493  const char* nameoutput=GetOutputSlot(4)->GetContainer()->GetName();
494  if (copyfCuts){
495  copyfCuts->SetName(nameoutput);
496 
497  //Post the data
498  PostData(4, copyfCuts);
499  }
500  else{
501  AliFatal("Failing initializing AliRDHFCuts object - Exiting...");
502  }
503 
504  fListProfiles = new TList();
505  fListProfiles->SetOwner();
506  TString period[4];
507  Int_t nProfiles=4;
508 
509  if (fIsPPbData) { //if pPb, use only two estimator histos
510  period[0] = "LHC13b"; period[1] = "LHC13c";
511  nProfiles = 2;
512  } else { // else assume pp (four histos for LHC10)
513  period[0] = "LHC10b"; period[1] = "LHC10c"; period[2] = "LHC10d"; period[3] = "LHC10e";
514  nProfiles = 4;
515  }
516 
517  for(Int_t i=0; i<nProfiles; i++){
518  if(fMultEstimatorAvg[i]){
519  TProfile* hprof=new TProfile(*fMultEstimatorAvg[i]);
520  hprof->SetName(Form("ProfileTrkVsZvtx%s\n",period[i].Data()));
521  fListProfiles->Add(hprof);
522  }
523  }
524 
525  // Save also the weight functions or histograms
530 
531  PostData(5,fListProfiles);
532 
533  return;
534 }
535 
536 //_________________________________________________
538 {
539  //
540  // Main loop function
541  //
542 
543  PostData(1,fHistEventsProcessed) ;
544  PostData(2,fCFManager->GetParticleContainer()) ;
545  PostData(3,fCorrelation) ;
546 
547  AliDebug(3,Form("*** Processing event %d\n", fEvents));
548 
549  if (fFillFromGenerated){
550  AliWarning("Flag to fill container with generated value ON ---> dca, d0pi, d0K, d0xd0, cosPointingAngle will be set as dummy!");
551  }
552 
553  if (!fInputEvent) {
554  Error("UserExec","NO EVENT FOUND!");
555  return;
556  }
557 
558  AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);
559 
560  TClonesArray *arrayBranch=0;
561 
562  if(!aodEvent && AODEvent() && IsStandardAOD()) {
563  // In case there is an AOD handler writing a standard AOD, use the AOD
564  // event in memory rather than the input (ESD) event.
565  aodEvent = dynamic_cast<AliAODEvent*> (AODEvent());
566  // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
567  // have to taken from the AOD event hold by the AliAODExtension
568  AliAODHandler* aodHandler = (AliAODHandler*)
569  ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
570  if(aodHandler->GetExtensions()) {
571  AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
572  AliAODEvent *aodFromExt = ext->GetAOD();
573 
574  switch (fDecayChannel){
575  case 2:{
576  arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi");
577  break;
578  }
579  case 21:{
580  arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");
581  break;
582  }
583  case 22:{
584  arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("CascadesHF");
585  break;
586  }
587  case 31:
588  case 32:
589  case 33:{
590  arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
591  break;
592  }
593  case 4:{
594  arrayBranch=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm4Prong");
595  break;
596  }
597  default:
598  break;
599  }
600  }
601  }
602  else {
603  switch (fDecayChannel){
604  case 2:{
605  arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("D0toKpi");
606  break;
607  }
608  case 21:{
609  arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("Dstar");
610  break;
611  }
612  case 22:{
613  arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("CascadesHF");
614  break;
615  }
616  case 31:
617  case 32:
618  case 33:{
619  arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("Charm3Prong");
620  break;
621  }
622  case 4:{
623  arrayBranch=(TClonesArray*)aodEvent->GetList()->FindObject("Charm4Prong");
624  break;
625  }
626  default:
627  break;
628  }
629  }
630 
631  AliAODVertex *aodVtx = (AliAODVertex*)aodEvent->GetPrimaryVertex();
632  if (!aodVtx) {
633  AliDebug(3, "The event was skipped due to missing vertex");
634  return;
635  }
636 
637  if (!arrayBranch) {
638  AliError("Could not find array of HF vertices");
639  return;
640  }
641 
642  fEvents++;
643 
644  fCFManager->SetRecEventInfo(aodEvent);
645  fCFManager->SetMCEventInfo(aodEvent);
646 
647  //******** DEFINE number of variables of the container***** for now set at 13, in the future in the config macro.
648 
649  //loop on the MC event
650 
651  TClonesArray* mcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));
652  if (!mcArray) {
653  AliError("Could not find Monte-Carlo in AOD");
654  return;
655  }
656  Int_t icountMC = 0;
657  Int_t icountGenLimAcc = 0;
658  Int_t icountGenLimAccNoAcc = 0;
659  Int_t icountAcc = 0;
660  Int_t icountReco = 0;
661  Int_t icountVertex = 0;
662  Int_t icountRefit = 0;
663  Int_t icountRecoAcc = 0;
664  Int_t icountRecoITSClusters = 0;
665  Int_t icountRecoPPR = 0;
666  Int_t icountRecoPID = 0;
667  Int_t cquarks = 0;
668 
669  AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(aodEvent->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
670  if (!mcHeader) {
671  AliError("Could not find MC Header in AOD");
672  return;
673  }
674 
675  fHistEventsProcessed->Fill(0.5);
676 
677  Double_t* containerInput = new Double_t[fNvar];
678  Double_t* containerInputMC = new Double_t[fNvar];
679 
680 
681  AliCFVertexingHF* cfVtxHF=0x0;
682  switch (fDecayChannel){
683  case 2:{
684  cfVtxHF = new AliCFVertexingHF2Prong(mcArray, fOriginDselection);
685  break;
686  }
687  case 21:{
688  cfVtxHF = new AliCFVertexingHFCascade(mcArray, fOriginDselection);
689  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGcascade(413);
690  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGbachelor(211);
691  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaugh(421);
692  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaughForMC(421);
693  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaughPositive(211);
694  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaughNegative(321);
695  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPrimaryVertex(aodVtx);
696  break;
697  }
698  case 22:{
699  // Lc -> K0S+proton
701  cfVtxHF = new AliCFVertexingHFCascade(mcArray, fOriginDselection);
702  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGcascade(4122);
703  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGbachelor(2212);
704  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaugh(310);
705  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaughForMC(311);
706  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaughPositive(211);
707  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPDGneutrDaughNegative(211);
708  ((AliCFVertexingHFCascade*)cfVtxHF)->SetPrimaryVertex(aodVtx);
709  ((AliCFVertexingHFCascade*)cfVtxHF)->SetCutOnMomConservation(fCutOnMomConservation);
710  if (fUseAdditionalCuts) ((AliCFVertexingHFCascade*)cfVtxHF)->SetUseCutsForTMVA(fUseCutsForTMVA);
711  }
712  else {
714  }
715  break;
716  }
717  case 31:
718  // case 32:
719  case 33:{
720  cfVtxHF = new AliCFVertexingHF3Prong(mcArray, fOriginDselection, fDecayChannel);
721  if(fDecayChannel==33){
722  ((AliCFVertexingHF3Prong*)cfVtxHF)->SetGeneratedDsOption(fGenDsOption);
723  }
724  break;
725  }
726  case 32:{
728  }
729  case 4:{
730  //cfVtxHF = new AliCFVertexingHF4Prong(mcArray, originDselection); // not there yet
731  break;
732  }
733  default:
734  break;
735  }
736  if (!cfVtxHF){
737  AliError("No AliCFVertexingHF initialized");
738  delete[] containerInput;
739  delete[] containerInputMC;
740  return;
741  }
742 
743  Double_t zPrimVertex = aodVtx ->GetZ();
744  Double_t zMCVertex = mcHeader->GetVtxZ();
745  Int_t runnumber = aodEvent->GetRunNumber();
746 
747  // Multiplicity definition with tracklets
748  Double_t nTracklets = 0;
749  Int_t nTrackletsEta10 = static_cast<Int_t>(AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aodEvent,-1.,1.));
750  Int_t nTrackletsEta16 = static_cast<Int_t>(AliVertexingHFUtils::GetNumberOfTrackletsInEtaRange(aodEvent,-1.6,1.6));
751  nTracklets = (Double_t)nTrackletsEta10;
752  if(fMultiplicityEstimator==kNtrk10to16) { nTracklets = (Double_t)(nTrackletsEta16 - nTrackletsEta10); }
753 
754  // Apply the Ntracklets z-vtx data driven correction
756  TProfile* estimatorAvg = GetEstimatorHistogram(aodEvent);
757  if(estimatorAvg) {
758  Int_t nTrackletsEta10Corr = static_cast<Int_t>(AliVertexingHFUtils::GetCorrectedNtracklets(estimatorAvg,nTrackletsEta10,zPrimVertex,fRefMult));
759  Int_t nTrackletsEta16Corr = static_cast<Int_t>(AliVertexingHFUtils::GetCorrectedNtracklets(estimatorAvg,nTrackletsEta16,zPrimVertex,fRefMult));
760  nTracklets = (Double_t)nTrackletsEta10Corr;
761  if(fMultiplicityEstimator==kNtrk10to16) { nTracklets = (Double_t)(nTrackletsEta16Corr - nTrackletsEta10Corr); }
762  }
763  }
764 
765 
766  fWeight=1.;
767  if(fUseZWeight) fWeight *= GetZWeight(zMCVertex,runnumber);
768  if(fUseNchWeight){
769  Int_t nChargedMCPhysicalPrimary=AliVertexingHFUtils::GetGeneratedPhysicalPrimariesInEtaRange(mcArray,-1.0,1.0);
770  if(!fUseTrackletsWeight) fWeight *= GetNchWeight(nChargedMCPhysicalPrimary);
771  else fWeight *= GetNchWeight(static_cast<Int_t>(nTracklets));
772  AliDebug(2,Form("Using Nch weights, Mult=%d Weight=%f\n",nChargedMCPhysicalPrimary,fWeight));
773  }
774  Double_t eventWeight=fWeight;
775 
776  if (TMath::Abs(zMCVertex) > fCuts->GetMaxVtxZ()){
777  AliDebug(3,Form("z coordinate of MC vertex = %f, it was required to be within [-%f, +%f], skipping event", zMCVertex, fCuts->GetMaxVtxZ(), fCuts->GetMaxVtxZ()));
778  delete[] containerInput;
779  delete[] containerInputMC;
780  delete cfVtxHF;
781  return;
782  }
783 
784  if(aodEvent->GetTriggerMask()==0 &&
785  (runnumber>=195344 && runnumber<=195677)){
786  AliDebug(3,"Event rejected because of null trigger mask");
787  delete[] containerInput;
788  delete[] containerInputMC;
789  delete cfVtxHF;
790  return;
791  }
792 
793  AliESDtrackCuts** trackCuts = new AliESDtrackCuts*[cfVtxHF->GetNProngs()];
794  if (fDecayChannel == 21){
795  // for the D*, setting the third element of the array of the track cuts to those for the soft pion
796  for (Int_t iProng = 0; iProng<cfVtxHF->GetNProngs()-1; iProng++){
797  trackCuts[iProng]=fCuts->GetTrackCuts();
798  }
799  trackCuts[2] = fCuts->GetTrackCutsSoftPi();
800  }
801  else if (fDecayChannel == 22) {
802  // for the Lc->V0+bachelor, setting the second and third elements of the array of the track cuts to those for the V0 daughters
803  trackCuts[0]=fCuts->GetTrackCuts();
804  trackCuts[1]=fCuts->GetTrackCutsV0daughters();
805  trackCuts[2]=fCuts->GetTrackCutsV0daughters();
806  }
807  else {
808  for (Int_t iProng = 0; iProng<cfVtxHF->GetNProngs(); iProng++){
809  trackCuts[iProng]=fCuts->GetTrackCuts();
810  }
811  }
812 
813  //General settings: vertex, feed down and fill reco container with generated values.
814  cfVtxHF->SetRecoPrimVertex(zPrimVertex);
815  cfVtxHF->SetMCPrimaryVertex(zMCVertex);
817  cfVtxHF->SetNVar(fNvar);
821 
822  // switch-off the trigger class selection (doesn't work for MC)
823  fCuts->SetTriggerClass("");
824 
825  // MC vertex, to be used, in case, for pp
827 
828  if (fCentralitySelection){ // keep only the requested centrality
829 
830  if(fCuts->IsEventSelectedInCentrality(aodEvent)!=0) {
831  delete[] containerInput;
832  delete[] containerInputMC;
833  delete [] trackCuts;
834  delete cfVtxHF;
835  return;
836  }
837  } else { // keep all centralities
838  fCuts->SetMinCentrality(0.);
839  fCuts->SetMaxCentrality(100.);
840  }
841 
842  Float_t centValue = 0.;
843  if(!fIsPPData) centValue = fCuts->GetCentrality(aodEvent);
844  cfVtxHF->SetCentralityValue(centValue);
845 
846  // multiplicity estimator with VZERO
847  Double_t vzeroMult=0;
848  AliAODVZERO *vzeroAOD = (AliAODVZERO*)aodEvent->GetVZEROData();
849  if(vzeroAOD) vzeroMult = vzeroAOD->GetMTotV0A() + vzeroAOD->GetMTotV0C();
850 
851  Double_t multiplicity = nTracklets; // set to the Ntracklet estimator
852  if(fMultiplicityEstimator==kVZERO) { multiplicity = vzeroMult; }
853 
854  cfVtxHF->SetMultiplicity(multiplicity);
855 
856  Double_t q2=0;
857  if(fConfiguration==kESE) {
858  //set q2 in case of kESE configuration
859  q2=ComputeTPCq2(aodEvent,mcHeader,-0.8,0.8,0.2,5.);
860  cfVtxHF->Setq2Value(q2);
861 
862  //set track array in case of kESE configuration
863  cfVtxHF->SetTrackArray(aodEvent->GetTracks());
864  }
865 
866  // printf("Multiplicity estimator %d, value %2.2f\n",fMultiplicityEstimator,multiplicity);
867 
868  for (Int_t iPart=0; iPart<mcArray->GetEntriesFast(); iPart++) {
869  AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(mcArray->At(iPart));
870  if (!mcPart){
871  AliError("Failed casting particle from MC array!, Skipping particle");
872  continue;
873  }
874 
875  //counting c quarks
876  cquarks += cfVtxHF->MCcquarkCounting(mcPart);
877 
878  // check the MC-level cuts, must be the desidered particle
879  if (!fCFManager->CheckParticleCuts(0, mcPart)) {
880  AliDebug(2,"Check the MC-level cuts - not desidered particle");
881  continue; // 0 stands for MC level
882  }
883  else {
884  AliDebug(3, Form("\n\n---> COOL! we found a particle (particle %d)!!! with PDG code = %d \n\n", iPart, mcPart->GetPdgCode()));
885  }
886  cfVtxHF->SetMCCandidateParam(iPart);
887 
888 
889  if (!(cfVtxHF->SetLabelArray())){
890  AliDebug(2,Form("Impossible to set the label array for particle %d (decaychannel = %d)", iPart, fDecayChannel));
891  continue;
892  }
893 
894  //check the candiate family at MC level
895  if (!(cfVtxHF->CheckMCPartFamily(mcPart, mcArray))) {
896  AliDebug(2,Form("Check on the family wrong for particle %d!!! (decaychannel = %d)", iPart, fDecayChannel));
897  continue;
898  }
899  else{
900  AliDebug(2,Form("Check on the family OK for particle %d!!! (decaychannel = %d)", iPart, fDecayChannel));
901  }
902 
903  //Fill the MC container
904  Bool_t mcContainerFilled = cfVtxHF -> FillMCContainer(containerInputMC);
905  AliDebug(2, Form("particle = %d mcContainerFilled = %d", iPart, mcContainerFilled));
906  if (mcContainerFilled) {
907  if (fUseWeight){
908  if (fHistoPtWeight) { // using an histogram as weight function
909  AliDebug(2,"Using Histogram as Pt weight function");
910  fWeight = eventWeight*GetPtWeightFromHistogram(containerInputMC[0]);
911  }
912  else if (fFuncWeight){ // using user-defined function
913  AliDebug(2,"Using function");
914  fWeight = eventWeight*fFuncWeight->Eval(containerInputMC[0]);
915  }
916  else{ // using FONLL
917  AliDebug(2,"Using FONLL");
918  fWeight = eventWeight*GetWeight(containerInputMC[0]);
919  }
920  AliDebug(2,Form("pt = %f, weight = %f",containerInputMC[0], fWeight));
921  }
922  if (!fCuts->IsInFiducialAcceptance(containerInputMC[0],containerInputMC[1])) {
923  AliDebug(3, Form("Not in limited acceptance, containerInputMC[0] = %f, containerInputMC[1] = %f", containerInputMC[0], containerInputMC[1]));
924  continue;
925  }
926  else{
927  AliDebug(3, Form("YES!! in limited acceptance, containerInputMC[0] = %f, containerInputMC[1] = %f", containerInputMC[0],containerInputMC[1]));
928  }
929 
930  //MC Limited Acceptance
931  if (TMath::Abs(containerInputMC[1]) < 0.5) {
932  fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepGeneratedLimAcc, fWeight);
933  icountGenLimAcc++;
934  AliDebug(3,"MC Lim Acc container filled\n");
935  if (fCFManager->GetParticleContainer()->GetNStep() == kStepGenLimAccNoAcc+1) {
936  if (!(cfVtxHF-> MCAcceptanceStep())) {
937  if(!fFillMinimumSteps)
938  fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepGenLimAccNoAcc, fWeight);
939  icountGenLimAccNoAcc++;
940  AliDebug(3,"MC Lim Acc No Acc container filled\n");
941  }
942  }
943  }
944 
945  //MC
946  if(!fFillMinimumSteps)
947  fCFManager->GetParticleContainer()->Fill(containerInputMC, kStepGenerated, fWeight);
948  icountMC++;
949  AliDebug(3,"MC container filled \n");
950 
951  // MC in acceptance
952  // check the MC-Acceptance level cuts
953  // since standard CF functions are not applicable, using Kine Cuts on daughters
954  Bool_t mcAccepStep = cfVtxHF-> MCAcceptanceStep();
955  if (mcAccepStep){
956  fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepAcceptance, fWeight);
957  AliDebug(3,"MC acceptance cut passed\n");
958  icountAcc++;
959 
960  //MC Vertex step
961  if (fCuts->IsEventSelected(aodEvent)){
962  // filling the container if the vertex is ok
963  if(!fFillMinimumSteps)
964  fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepVertex, fWeight) ;
965  AliDebug(3,"Vertex cut passed and container filled\n");
966  icountVertex++;
967 
968  //mc Refit requirement
969  Bool_t mcRefitStep = cfVtxHF->MCRefitStep(aodEvent, &trackCuts[0]);
970  if (mcRefitStep){
971  if(!fFillMinimumSteps)
972  fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepRefit, fWeight);
973  AliDebug(3,"MC Refit cut passed and container filled\n");
974  icountRefit++;
975  }
976  else{
977  AliDebug(3,"MC Refit cut not passed\n");
978  continue;
979  }
980  }
981  else{
982  AliDebug (3, "MC vertex step not passed\n");
983  continue;
984  }
985  }
986  else{
987  AliDebug (3, "MC in acceptance step not passed\n");
988  continue;
989  }
990  }
991  else {
992  AliDebug (3, "MC container not filled\n");
993  }
994  }
995 
996  if (cquarks<2) AliDebug(2,Form("Event with %d c-quarks", cquarks));
997  AliDebug(2,Form("Found %i MC particles that are %s!!",icountMC,fPartName.Data()));
998  AliDebug(2,Form("Found %i MC particles that are %s and satisfy Acc cuts!!",icountAcc,fPartName.Data()));
999  AliDebug(2,Form("Found %i MC particles that are %s and satisfy Vertex cuts!!",icountVertex,fPartName.Data()));
1000  AliDebug(2,Form("Found %i MC particles that are %s and satisfy Refit cuts!!",icountRefit,fPartName.Data()));
1001 
1002  // Now go to rec level
1003  fCountMC += icountMC;
1004  fCountGenLimAcc += icountGenLimAcc;
1005  fCountGenLimAccNoAcc += icountGenLimAccNoAcc;
1006  fCountAcc += icountAcc;
1007  fCountVertex+= icountVertex;
1008  fCountRefit+= icountRefit;
1009 
1010  AliDebug(2,Form("Found %d vertices for decay channel %d",arrayBranch->GetEntriesFast(),fDecayChannel));
1012  for(Int_t iCandid = 0; iCandid<arrayBranch->GetEntriesFast();iCandid++){
1013  fHistEventsProcessed->Fill(2.5);
1014  AliAODRecoDecayHF* charmCandidate=0x0;
1015  switch (fDecayChannel){
1016  case 2:{
1017  charmCandidate = (AliAODRecoDecayHF2Prong*)arrayBranch->At(iCandid);
1018  if(charmCandidate->GetIsFilled()!=0) fHistEventsProcessed->Fill(3.5);
1019  if(!vHF->FillRecoCand(aodEvent,(AliAODRecoDecayHF2Prong*)charmCandidate)){
1020  fHistEventsProcessed->Fill(5.5);
1021  continue;
1022  }else{
1023  fHistEventsProcessed->Fill(4.5);
1024  }
1025  break;
1026  }
1027  case 21:{
1028  charmCandidate = (AliAODRecoCascadeHF*)arrayBranch->At(iCandid);
1029  if(charmCandidate->GetIsFilled()!=0) fHistEventsProcessed->Fill(3.5);
1030  if(!vHF->FillRecoCasc(aodEvent,((AliAODRecoCascadeHF*)charmCandidate),kTRUE)){
1031  fHistEventsProcessed->Fill(5.5);
1032  continue; //DStar
1033  }else{
1034  fHistEventsProcessed->Fill(4.5);
1035  }
1036  break;
1037  }
1038  case 22:{
1039  charmCandidate = (AliAODRecoCascadeHF*)arrayBranch->At(iCandid);
1040  if(charmCandidate->GetIsFilled()!=0) fHistEventsProcessed->Fill(3.5);
1041  if(!vHF->FillRecoCasc(aodEvent,((AliAODRecoCascadeHF*)charmCandidate),kFALSE)){
1042  fHistEventsProcessed->Fill(5.5);
1043  continue; //Cascade
1044  }else{
1045  fHistEventsProcessed->Fill(4.5);
1046  }
1047  break;
1048  }
1049  case 31:
1050  case 32:
1051  case 33:{
1052  charmCandidate = (AliAODRecoDecayHF3Prong*)arrayBranch->At(iCandid);
1053  if(charmCandidate->GetIsFilled()!=0) fHistEventsProcessed->Fill(3.5);
1054  if(!vHF->FillRecoCand(aodEvent,(AliAODRecoDecayHF3Prong*)charmCandidate)){
1055  fHistEventsProcessed->Fill(5.5);
1056  continue;
1057  }else{
1058  fHistEventsProcessed->Fill(4.5);
1059  }
1060  break;
1061  }
1062  case 4:{
1063  charmCandidate = (AliAODRecoDecayHF4Prong*)arrayBranch->At(iCandid);
1064  break;
1065  }
1066  default:
1067  break;
1068  }
1069 
1070  Bool_t unsetvtx=kFALSE;
1071  if(!charmCandidate->GetOwnPrimaryVtx()) {
1072  charmCandidate->SetOwnPrimaryVtx(aodVtx); // needed to compute all variables
1073  unsetvtx=kTRUE;
1074  }
1075 
1076  Bool_t signAssociation = cfVtxHF->SetRecoCandidateParam((AliAODRecoDecayHF*)charmCandidate);
1077  if (!signAssociation){
1078  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1079  continue;
1080  }
1081 
1082  Int_t isPartOrAntipart = cfVtxHF->CheckReflexion(fSign);
1083  if (isPartOrAntipart == 0){
1084  AliDebug(2, Form("The candidate pdg code doesn't match the requirement set in the task (fSign = %d)",fSign));
1085  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1086  continue;
1087  }
1088 
1089  AliDebug(3,Form("iCandid=%d - signAssociation=%d, isPartOrAntipart=%d",iCandid, signAssociation, isPartOrAntipart));
1090 
1091  Bool_t recoContFilled = cfVtxHF->FillRecoContainer(containerInput);
1092  AliDebug(3, Form("CF task: RecoContFilled for candidate %d is %d", iCandid, (Int_t)recoContFilled));
1093  if (recoContFilled){
1094 
1095  // weight according to pt
1096  if (fUseWeight){
1097  if (fHistoPtWeight) {
1098  AliDebug(2,"Using Histogram as Pt weight function");
1099  fWeight = eventWeight*GetPtWeightFromHistogram(containerInput[0]);
1100  }
1101  else if (fFuncWeight){ // using user-defined function
1102  AliDebug(2, "Using function");
1103  fWeight = eventWeight*fFuncWeight->Eval(containerInput[0]);
1104  }
1105  else{ // using FONLL
1106  AliDebug(2, "Using FONLL");
1107  fWeight = eventWeight*GetWeight(containerInput[0]);
1108  }
1109  AliDebug(2, Form("pt = %f, weight = %f",containerInput[0], fWeight));
1110  }
1111 
1112  if (!fCuts->IsInFiducialAcceptance(containerInput[0],containerInput[1])){
1113  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1114  continue;
1115  }
1116  //Reco Step
1117  Bool_t recoStep = cfVtxHF->RecoStep();
1118  if (recoStep) AliDebug(2, Form("particle = %d --> CF task: Reco step for candidate %d is %d", iCandid, iCandid, (Int_t)recoStep));
1119  Bool_t vtxCheck = fCuts->IsEventSelected(aodEvent);
1120 
1121 
1122  // Selection on the filtering bit
1123  Bool_t isBitSelected = kTRUE;
1124  if(fDecayChannel==2) {
1125  if(fUseSelectionBit && !charmCandidate->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts)) isBitSelected = kFALSE;
1126  }else if(fDecayChannel==31){
1127  if(fUseSelectionBit && !charmCandidate->HasSelectionBit(AliRDHFCuts::kDplusCuts)) isBitSelected = kFALSE;
1128  }else if(fDecayChannel==33){
1129  if(fUseSelectionBit && !charmCandidate->HasSelectionBit(AliRDHFCuts::kDsCuts)) isBitSelected = kFALSE;
1130  }else if(fDecayChannel==32){
1131  if(fUseSelectionBit && !charmCandidate->HasSelectionBit(AliRDHFCuts::kLcCuts)) isBitSelected = kFALSE;
1132  }else if(fDecayChannel==22){
1133  if(!((dynamic_cast<AliAODRecoCascadeHF*>(charmCandidate))->CheckCascadeFlags())) isBitSelected = kFALSE; // select only Lc among cascade candidates
1134  }
1135  if(!isBitSelected){
1136  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1137  continue;
1138  }
1139 
1140 
1141  if (recoStep && recoContFilled && vtxCheck){
1142  if(!fFillMinimumSteps)
1143  fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructed, fWeight) ;
1144  icountReco++;
1145  AliDebug(3,"Reco step passed and container filled\n");
1146 
1147  //Reco in the acceptance -- take care of UNFOLDING!!!!
1148  Bool_t recoAcceptanceStep = cfVtxHF->RecoAcceptStep(&trackCuts[0]);
1149  if (recoAcceptanceStep) {
1150  if(!fFillMinimumSteps)
1151  fCFManager->GetParticleContainer()->Fill(containerInput,kStepRecoAcceptance, fWeight) ;
1152  icountRecoAcc++;
1153  AliDebug(3,"Reco acceptance cut passed and container filled\n");
1154 
1155  if(fAcceptanceUnf){
1156  Double_t fill[4]; //fill response matrix
1157  Bool_t bUnfolding = cfVtxHF -> FillUnfoldingMatrix(fPDGcode,fill);
1158  if (bUnfolding) fCorrelation->Fill(fill);
1159  }
1160 
1161  //Number of ITS cluster requirements
1162  Int_t recoITSnCluster = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kTracks, aodEvent);
1163  if (recoITSnCluster){
1164  if(!fFillMinimumSteps)
1165  fCFManager->GetParticleContainer()->Fill(containerInput,kStepRecoITSClusters, fWeight) ;
1166  icountRecoITSClusters++;
1167  AliDebug(3,"Reco n ITS cluster cut passed and container filled\n");
1168 
1169  Bool_t iscutsusingpid = fCuts->GetIsUsePID();
1170  Int_t recoAnalysisCuts = -1, recoPidSelection = -1;
1171  fCuts->SetUsePID(kFALSE);
1172  recoAnalysisCuts = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kCandidate, aodEvent);
1173 
1174  if (fDecayChannel==33){ // Ds case, where more possibilities are considered
1175  Bool_t keepDs=ProcessDs(recoAnalysisCuts);
1176  if(keepDs) recoAnalysisCuts=3;
1177  }
1178  else if (fDecayChannel==22){ // Lc->V0+bachelor case, where more possibilities are considered
1179  Bool_t keepLctoV0bachelor = ProcessLctoV0Bachelor(recoAnalysisCuts);
1180  if (keepLctoV0bachelor) recoAnalysisCuts=3;
1181  AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
1182  AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler());
1183  AliPIDResponse* pidResponse = inputHandler->GetPIDResponse();
1184  if (fUseAdditionalCuts){
1185  if (!((AliCFVertexingHFCascade*)cfVtxHF)->CheckAdditionalCuts(pidResponse)) recoAnalysisCuts = -1;
1186  }
1187  }
1188 
1189  fCuts->SetUsePID(iscutsusingpid); //restoring usage of the PID from the cuts object
1190  Bool_t tempAn=(recoAnalysisCuts == 3 || recoAnalysisCuts == isPartOrAntipart);
1191  if (fDecayChannel == 22) tempAn = (recoAnalysisCuts == 3);
1192  if (fDecayChannel == 32) tempAn=(recoAnalysisCuts >0 || recoAnalysisCuts == isPartOrAntipart);
1193 
1194  if (tempAn){
1195  fCFManager->GetParticleContainer()->Fill(containerInput, kStepRecoPPR, fWeight);
1196  icountRecoPPR++;
1197  AliDebug(3,"Reco Analysis cuts passed and container filled \n");
1198  //pid selection
1199  //recoPidSelection = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kPID);
1200  //if((fCuts->CombineSelectionLevels(3,recoAnalysisCuts,recoPidSelection)==isPartOrAntipart)||(fCuts->CombineSelectionLevels(3,recoAnalysisCuts,recoPidSelection)==3)){
1201  recoPidSelection = fCuts->IsSelected(charmCandidate, AliRDHFCuts::kCandidate, aodEvent);
1202 
1203  if (fDecayChannel==33){ // Ds case, where more possibilities are considered
1204  Bool_t keepDs=ProcessDs(recoPidSelection);
1205  if(keepDs) recoPidSelection=3;
1206  } else if (fDecayChannel==22){ // Lc->V0+bachelor case, where more possibilities are considered
1207  Bool_t keepLctoV0bachelor=ProcessLctoV0Bachelor(recoPidSelection);
1208  if (keepLctoV0bachelor) recoPidSelection=3;
1209  }
1210 
1211  Bool_t tempPid=(recoPidSelection == 3 || recoPidSelection == isPartOrAntipart);
1212  if (fDecayChannel == 22) tempPid = (recoPidSelection == 3);
1213  if (fDecayChannel == 32) tempPid=(recoPidSelection >0 || recoPidSelection == isPartOrAntipart);
1214 
1215  if (tempPid){
1216  Double_t weigPID = 1.;
1217  if (fDecayChannel == 2){ // D0 with Bayesian PID using weights
1218  if(((AliRDHFCutsD0toKpi*)fCuts)->GetCombPID() && (((AliRDHFCutsD0toKpi*)fCuts)->GetBayesianStrategy() == AliRDHFCutsD0toKpi::kBayesWeight || ((AliRDHFCutsD0toKpi*)fCuts)->GetBayesianStrategy() == AliRDHFCutsD0toKpi::kBayesWeightNoFilter)){
1219  if (isPartOrAntipart == 1){
1220  weigPID = ((AliRDHFCutsD0toKpi*)fCuts)->GetWeightsNegative()[AliPID::kKaon] * ((AliRDHFCutsD0toKpi*)fCuts)->GetWeightsPositive()[AliPID::kPion];
1221  }else if (isPartOrAntipart == 2){
1222  weigPID = ((AliRDHFCutsD0toKpi*)fCuts)->GetWeightsPositive()[AliPID::kKaon] * ((AliRDHFCutsD0toKpi*)fCuts)->GetWeightsNegative()[AliPID::kPion];
1223  }
1224  if ((weigPID < 0) || (weigPID > 1)) weigPID = 0.;
1225  }
1226  }else if (fDecayChannel == 33){ // Ds with Bayesian PID using weights
1228  Int_t labDau0=((AliAODTrack*)charmCandidate->GetDaughter(0))->GetLabel();
1229  AliAODMCParticle* firstDau=(AliAODMCParticle*)mcArray->UncheckedAt(TMath::Abs(labDau0));
1230  if(firstDau){
1231  Int_t pdgCode0=TMath::Abs(firstDau->GetPdgCode());
1232  if(pdgCode0==321){
1233  weigPID=((AliRDHFCutsDstoKKpi*)fCuts)->GetWeightForKKpi();
1234  }else if(pdgCode0==211){
1235  weigPID=((AliRDHFCutsDstoKKpi*)fCuts)->GetWeightForpiKK();
1236  }
1237  if ((weigPID < 0) || (weigPID > 1)) weigPID = 0.;
1238  }else{
1239  weigPID=0.;
1240  }
1241  }
1242  }
1243  fCFManager->GetParticleContainer()->Fill(containerInput, kStepRecoPID, fWeight*weigPID);
1244  icountRecoPID++;
1245  AliDebug(3,"Reco PID cuts passed and container filled \n");
1246  if(!fAcceptanceUnf){
1247  Double_t fill[4]; //fill response matrix
1248  Bool_t bUnfolding = cfVtxHF -> FillUnfoldingMatrix(fPDGcode,fill);
1249  if (bUnfolding) fCorrelation->Fill(fill);
1250  }
1251  }
1252  else {
1253  AliDebug(3, "Analysis Cuts step not passed \n");
1254  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1255  continue;
1256  }
1257  }
1258  else {
1259  AliDebug(3, "PID selection not passed \n");
1260  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1261  continue;
1262  }
1263  }
1264  else{
1265  AliDebug(3, "Number of ITS cluster step not passed\n");
1266  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1267  continue;
1268  }
1269  }
1270  else{
1271  AliDebug(3, "Reco acceptance step not passed\n");
1272  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1273  continue;
1274  }
1275  }
1276  else {
1277  AliDebug(3, "Reco step not passed\n");
1278  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1279  continue;
1280  }
1281  }
1282 
1283  if(unsetvtx) charmCandidate->UnsetOwnPrimaryVtx();
1284  } // end loop on candidate
1285  delete vHF;
1286  fCountReco+= icountReco;
1287  fCountRecoAcc+= icountRecoAcc;
1288  fCountRecoITSClusters+= icountRecoITSClusters;
1289  fCountRecoPPR+= icountRecoPPR;
1290  fCountRecoPID+= icountRecoPID;
1291 
1292  fHistEventsProcessed->Fill(1.5);
1293 
1294  delete[] containerInput;
1295  delete[] containerInputMC;
1296  delete cfVtxHF;
1297  if (trackCuts){
1298  // for (Int_t i=0; i<cfVtxHF->GetNProngs(); i++){
1299  // delete [] trackCuts[i];
1300  // }
1301  delete [] trackCuts;
1302  }
1303 
1304 
1305 }
1306 
1307 //___________________________________________________________________________
1309 {
1310  // The Terminate() function is the last function to be called during
1311  // a query. It always runs on the client, it can be used to present
1312  // the results graphically or save the results to file.
1313 
1314  AliAnalysisTaskSE::Terminate();
1315 
1316  AliInfo(Form("Found %i MC particles that are %s in MC, in %d events",fCountMC,fPartName.Data(),fEvents));
1317  AliInfo(Form("Found %i MC particles that are %s in MC in limited acceptance, in %d events",fCountGenLimAcc,fPartName.Data(),fEvents));
1318  AliInfo(Form("Found %i MC particles that are %s in MC in limited acceptance and don't satisfy Acc cuts, in %d events",fCountGenLimAccNoAcc,fPartName.Data(),fEvents));
1319  AliInfo(Form("Found %i MC particles that are %s in MC and satisfy Acc cuts, in %d events",fCountAcc,fPartName.Data(),fEvents));
1320  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));
1321  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));
1322  AliInfo(Form("Found %i reco %s that are decaying in %s, in %d events",fCountReco,fPartName.Data(),fDauNames.Data(),fEvents));
1323  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));
1324  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));
1325  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));
1326  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));
1327 
1328  // draw some example plots....
1329  AliCFContainer *cont= dynamic_cast<AliCFContainer*> (GetOutputData(2));
1330  if(!cont) {
1331  printf("CONTAINER NOT FOUND\n");
1332  return;
1333  }
1334  // projecting the containers to obtain histograms
1335  // first argument = variable, second argument = step
1336 
1337  TH1D** h = new TH1D*[3];
1338  Int_t nvarToPlot = 0;
1339  if (fConfiguration == kSnail){
1340  //h = new TH1D[3][12];
1341  for (Int_t ih = 0; ih<3; ih++){
1342  if(fDecayChannel==22){
1343  nvarToPlot = 16;
1344  } else {
1345  nvarToPlot = 12;
1346  }
1347  h[ih] = new TH1D[nvarToPlot];
1348  }
1349  for(Int_t iC=1;iC<nvarToPlot; iC++){
1350  // MC-level
1351  h[0][iC] = *(cont->ShowProjection(iC,0));
1352  // MC-Acceptance level
1353  h[1][iC] = *(cont->ShowProjection(iC,1));
1354  // Reco-level
1355  h[2][iC] = *(cont->ShowProjection(iC,4));
1356  }
1357  }
1358  else if(fConfiguration == kCheetah){
1359  //h = new TH1D[3][12];
1360  nvarToPlot = 8;
1361  for (Int_t ih = 0; ih<3; ih++){
1362  h[ih] = new TH1D[nvarToPlot];
1363  }
1364  for(Int_t iC=0;iC<nvarToPlot; iC++){
1365  // MC-level
1366  h[0][iC] = *(cont->ShowProjection(iC,0));
1367  // MC-Acceptance level
1368  h[1][iC] = *(cont->ShowProjection(iC,1));
1369  // Reco-level
1370  h[2][iC] = *(cont->ShowProjection(iC,4));
1371  }
1372  }
1373  else if(fConfiguration == kFalcon){
1374  //h = new TH1D[3][12];
1375  nvarToPlot = 4;
1376  for (Int_t ih = 0; ih<3; ih++){
1377  h[ih] = new TH1D[nvarToPlot];
1378  }
1379  for(Int_t iC=0;iC<nvarToPlot; iC++){
1380  // MC-level
1381  h[0][iC] = *(cont->ShowProjection(iC,0));
1382  // MC-Acceptance level
1383  h[1][iC] = *(cont->ShowProjection(iC,1));
1384  // Reco-level
1385  h[2][iC] = *(cont->ShowProjection(iC,4));
1386  }
1387  }
1388  else if(fConfiguration == kESE){
1389  //h = new TH1D[3][12];
1390  nvarToPlot = 6;
1391  for (Int_t ih = 0; ih<3; ih++){
1392  h[ih] = new TH1D[nvarToPlot];
1393  }
1394  for(Int_t iC=0;iC<nvarToPlot; iC++){
1395  // MC-level
1396  h[0][iC] = *(cont->ShowProjection(iC,0));
1397  // MC-Acceptance level
1398  h[1][iC] = *(cont->ShowProjection(iC,1));
1399  // Reco-level
1400  h[2][iC] = *(cont->ShowProjection(iC,4));
1401  }
1402  }
1403  TString* titles;
1404  //Int_t nvarToPlot = 0;
1405  if (fConfiguration == kSnail){
1406  if(fDecayChannel==31){
1407  //nvarToPlot = 12;
1408  titles = new TString[nvarToPlot];
1409  titles[0]="pT_Dplus (GeV/c)";
1410  titles[1]="rapidity";
1411  titles[2]="phi (rad)";
1412  titles[3]="cT (#mum)";
1413  titles[4]="cosPointingAngle";
1414  titles[5]="pT_1 (GeV/c)";
1415  titles[6]="pT_2 (GeV/c)";
1416  titles[7]="pT_3 (GeV/c)";
1417  titles[8]="d0_1 (#mum)";
1418  titles[9]="d0_2 (#mum)";
1419  titles[10]="d0_3 (#mum)";
1420  titles[11]="zVertex (cm)";
1421  } else if (fDecayChannel==22) {
1422  //nvarToPlot = 16;
1423  titles = new TString[nvarToPlot];
1424  titles[0]="p_{T}(#Lambda_{c}) [GeV/c]";
1425  titles[1]="y(#Lambda_{c})";
1426  titles[2]="#varphi(#Lambda_{c}) [rad]";
1427  titles[3]="onTheFlyStatusV0";
1428  titles[4]="z_{vtx} [cm]";
1429  titles[5]="centrality";
1430  titles[6]="fake";
1431  titles[7]="multiplicity";
1432  //titles[8]="pT(bachelor) [GeV/c]";
1433  titles[8]="p(bachelor) [GeV/c]";
1434  titles[9]="p_{T}(V0) [GeV/c]";
1435  titles[10]="y(V0)";
1436  titles[11]="#varphi(V0) [rad]";
1437  titles[12]="m_{inv}(#pi^{+}#pi^{+}) [GeV/c^{2}]";
1438  titles[13]="dcaV0 (nSigma)";
1439  titles[14]="cosine pointing angle (V0)";
1440  titles[15]="cosine pointing angle (#Lambda_{c})";
1441  //titles[16]="c#tauV0 (#mum)";
1442  //titles[17]="c#tau (#mum)";
1443  } else {
1444  //nvarToPlot = 12;
1445  titles = new TString[nvarToPlot];
1446  titles[0]="pT_D0 (GeV/c)";
1447  titles[1]="rapidity";
1448  titles[2]="cosThetaStar";
1449  titles[3]="pT_pi (GeV/c)";
1450  titles[4]="pT_K (Gev/c)";
1451  titles[5]="cT (#mum)";
1452  titles[6]="dca (#mum)";
1453  titles[7]="d0_pi (#mum)";
1454  titles[8]="d0_K (#mum)";
1455  titles[9]="d0xd0 (#mum^2)";
1456  titles[10]="cosPointingAngle";
1457  titles[11]="phi (rad)";
1458  }
1459  }
1460  else if(fConfiguration == kCheetah){
1461  //nvarToPlot = 8;
1462  titles = new TString[nvarToPlot];
1463  if (fDecayChannel==22) {
1464  titles[0]="p_{T}(#Lambda_{c}) [GeV/c]";
1465  titles[1]="y(#Lambda_{c})";
1466  titles[2]="#varphi(#Lambda_{c}) [rad]";
1467  titles[3]="onTheFlyStatusV0";
1468  titles[4]="z_{vtx} [cm]";
1469  titles[5]="centrality";
1470  titles[6]="fake";
1471  titles[7]="multiplicity";
1472  } else {
1473  titles[0]="pT_candidate (GeV/c)";
1474  titles[1]="rapidity";
1475  titles[2]="cT (#mum)";
1476  titles[3]="phi";
1477  titles[4]="z_{vtx}";
1478  titles[5]="centrality";
1479  titles[6]="fake";
1480  titles[7]="multiplicity";
1481  }
1482  }
1483  else if(fConfiguration == kFalcon){
1484  //nvarToPlot = 4;
1485  titles = new TString[nvarToPlot];
1486  if (fDecayChannel==22) {
1487  titles[0]="p_{T}(#Lambda_{c}) [GeV/c]";
1488  titles[1]="y(#Lambda_{c})";
1489  titles[2]="centrality";
1490  titles[3]="multiplicity";
1491  } else {
1492  titles[0]="pT_candidate (GeV/c)";
1493  titles[1]="rapidity";
1494  titles[2]="centrality";
1495  titles[3]="multiplicity";
1496  }
1497  }
1498  else if(fConfiguration == kESE){
1499  //nvarToPlot = 4;
1500  titles = new TString[nvarToPlot];
1501  titles[0]="pT_candidate (GeV/c)";
1502  titles[1]="rapidity";
1503  titles[2]="centrality";
1504  titles[3]="multiplicity";
1505  titles[4]="N_{tracks} (R<0.4)";
1506  titles[5]="q_{2}";
1507  }
1508 
1509  Int_t markers[16]={20,24,21,25,27,28,
1510  20,24,21,25,27,28,
1511  20,24,21,25};
1512  Int_t colors[3]={2,8,4};
1513  for(Int_t iC=0;iC<nvarToPlot; iC++){
1514  for(Int_t iStep=0;iStep<3;iStep++){
1515  h[iStep][iC].SetTitle(titles[iC].Data());
1516  h[iStep][iC].GetXaxis()->SetTitle(titles[iC].Data());
1517  Double_t maxh=h[iStep][iC].GetMaximum();
1518  h[iStep][iC].GetYaxis()->SetRangeUser(0,maxh*1.2);
1519  h[iStep][iC].SetMarkerStyle(markers[iC]);
1520  h[iStep][iC].SetMarkerColor(colors[iStep]);
1521  }
1522  }
1523 
1524  gStyle->SetCanvasColor(0);
1525  gStyle->SetFrameFillColor(0);
1526  gStyle->SetTitleFillColor(0);
1527  gStyle->SetStatColor(0);
1528 
1529  // drawing in 2 separate canvas for a matter of clearity
1530  TCanvas * c1 =new TCanvas(Form("c1New_%d",fDecayChannel),"Vars 0, 1, 2, 3",1100,1200);
1531  c1->Divide(3,4);
1532  Int_t iPad=1;
1533  for(Int_t iVar=0; iVar<4; iVar++){
1534  if(iVar>=fNvar) continue;
1535  c1->cd(iPad++);
1536  h[0][iVar].DrawCopy("p");
1537  c1->cd(iPad++);
1538  h[1][iVar].DrawCopy("p");
1539  c1->cd(iPad++);
1540  h[2][iVar].DrawCopy("p");
1541  }
1542 
1544  TCanvas * c2 =new TCanvas(Form("c2New_%d",fDecayChannel),"Vars 4, 5, 6, 7",1100,1200);
1545  c2->Divide(3,4);
1546  iPad=1;
1547  for(Int_t iVar=4; iVar<8; iVar++){
1548  c2->cd(iPad++);
1549  h[0][iVar].DrawCopy("p");
1550  c2->cd(iPad++);
1551  h[1][iVar].DrawCopy("p");
1552  c2->cd(iPad++);
1553  h[2][iVar].DrawCopy("p");
1554  }
1555  }
1556 
1557  if (fConfiguration == kSnail){
1558  TCanvas * c3 =new TCanvas(Form("c3New_%d",fDecayChannel),"Vars 8, 9, 10, 11",1100,1200);
1559  c3->Divide(3,4);
1560  iPad=1;
1561  for(Int_t iVar=8; iVar<12; iVar++){
1562  c3->cd(iPad++);
1563  h[0][iVar].DrawCopy("p");
1564  c3->cd(iPad++);
1565  h[1][iVar].DrawCopy("p");
1566  c3->cd(iPad++);
1567  h[2][iVar].DrawCopy("p");
1568  }
1569  if (fDecayChannel==22) {
1570  TCanvas * c4 =new TCanvas(Form("c4New_%d",fDecayChannel),"Vars 12, 13, 14, 15",1100,1200);
1571  c4->Divide(3,4);
1572  iPad=1;
1573  for(Int_t iVar=12; iVar<16; iVar++){
1574  c4->cd(iPad++);
1575  h[0][iVar].DrawCopy("p");
1576  c4->cd(iPad++);
1577  h[1][iVar].DrawCopy("p");
1578  c4->cd(iPad++);
1579  h[2][iVar].DrawCopy("p");
1580  }
1581  }
1582  }
1583 
1584  /*
1585  THnSparseD* hcorr = dynamic_cast<THnSparseD*> (GetOutputData(3));
1586 
1587  TH2D* corr1 = hcorr->Projection(0,2);
1588  TH2D* corr2 = hcorr->Projection(1,3);
1589 
1590  TCanvas * c7 =new TCanvas(Form("c7New_%d",fDecayChannel),"",800,400);
1591  c7->Divide(2,1);
1592  c7->cd(1);
1593  corr1->DrawCopy("text");
1594  c7->cd(2);
1595  corr2->DrawCopy("text");
1596  */
1597  TFile* file_projection = new TFile("CFtaskHFprojectionNew.root","RECREATE");
1598 
1599  // corr1->Write();
1600  // corr2->Write();
1601 
1602  for(Int_t iC=0;iC<nvarToPlot; iC++){
1603  for(Int_t iStep=0;iStep<3;iStep++){
1604  h[iStep][iC].Write(Form("Step%d_%s",iStep,titles[iC].Data()));
1605  }
1606  }
1607  file_projection->Close();
1608  for (Int_t ih = 0; ih<3; ih++) delete [] h[ih];
1609  delete [] h;
1610  delete [] titles;
1611 
1612 }
1613 
1614 //___________________________________________________________________________
1616 {
1617  //HERE ONE CAN CREATE OUTPUT OBJECTS, IN PARTICULAR IF THE OBJECT PARAMETERS DON'T NEED
1618  //TO BE SET BEFORE THE EXECUTION OF THE TASK
1619  //
1620  Info("UserCreateOutputObjects","CreateOutputObjects of task %s\n", GetName());
1621 
1622  //slot #1
1623  OpenFile(1);
1624  const char* nameoutput=GetOutputSlot(1)->GetContainer()->GetName();
1625  fHistEventsProcessed = new TH1I(nameoutput,"",6,0,6) ;
1626  fHistEventsProcessed->GetXaxis()->SetBinLabel(1,"Events processed (all)");
1627  fHistEventsProcessed->GetXaxis()->SetBinLabel(2,"Events analyzed (after selection)");
1628  fHistEventsProcessed->GetXaxis()->SetBinLabel(3,"Candidates processed (all)");
1629  fHistEventsProcessed->GetXaxis()->SetBinLabel(4,"Candidates already filled");
1630  fHistEventsProcessed->GetXaxis()->SetBinLabel(5,"Candidates OK in FillRecoCand");
1631  fHistEventsProcessed->GetXaxis()->SetBinLabel(6,"Candidates failing in FillRecoCand");
1632 
1633  PostData(1,fHistEventsProcessed) ;
1634  PostData(2,fCFManager->GetParticleContainer()) ;
1635  PostData(3,fCorrelation) ;
1636 
1637 }
1638 
1639 
1640 //_________________________________________________________________________
1642  // ad-hoc weight function from ratio of
1643  // pt spectra from FONLL 2.76 TeV and
1644  // pt spectra from MC production LHC12a17a (PYTHIA Perugia0 with pthard bins)
1645  if(fFuncWeight) delete fFuncWeight;
1646  fFuncWeight=new TF1("funcWeight","[0]+[1]*TMath::Exp(-[2]*x)",0.,50.);
1647  fFuncWeight->SetParameter(0,4.63891e-02);
1648  fFuncWeight->SetParameter(1,1.51674e+01);
1649  fFuncWeight->SetParameter(2,4.09941e-01);
1650  fUseWeight=kTRUE;
1651 }
1652 //_________________________________________________________________________
1654  // ad-hoc weight function from ratio of
1655  // D0 pt spectra in PbPb 2011 0-10% centrality and
1656  // pt spectra from MC production LHC12a17a (PYTHIA Perugia0 with pthard bins)
1657  if(fFuncWeight) delete fFuncWeight;
1658  fFuncWeight=new TF1("funcWeight","[0]+[1]/TMath::Power(x,[2])",0.05,50.);
1659  fFuncWeight->SetParameter(0,1.43116e-02);
1660  fFuncWeight->SetParameter(1,4.37758e+02);
1661  fFuncWeight->SetParameter(2,3.08583);
1662  fUseWeight=kTRUE;
1663 }
1664 
1665 //_________________________________________________________________________
1667  // weight function from the ratio of the LHC12a17b MC
1668  // and FONLL calculations for pp data
1669  if(fFuncWeight) delete fFuncWeight;
1670  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.);
1671  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);
1672  fUseWeight=kTRUE;
1673 }
1674 
1675 //_________________________________________________________________________
1677  // weight function from the ratio of the LHC12a17b MC
1678  // and FONLL calculations for pp data
1679  // corrected by the BAMPS Raa calculation for 30-50% CC
1680  if(fFuncWeight) delete fFuncWeight;
1681  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.);
1682  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);
1683  fUseWeight=kTRUE;
1684 }
1685 
1686 //_________________________________________________________________________
1688  // weight function from the ratio of the LHC16i2a MC
1689  // 1.5-14 GeV/c using data and 1-1.5, 14-50 GeV/c using FONLL calculations
1690 
1691  if(fHistoPtWeight) delete fHistoPtWeight;
1692  fHistoPtWeight = new TH1F("histoWeight","histoWeight",500,0.,50.);
1693  Float_t binc[500]={ 1.695705, 1.743693, 1.790289, 1.835410, 1.878981, 1.920938, 1.961223, 1.999787, 2.036589, 2.071597, 2.104784, 2.136132, 2.165629, 2.193270, 2.219057, 2.174545, 2.064698, 1.959489, 1.858770, 1.762396, 1.670224, 1.582115, 1.497931, 1.417541, 1.340814, 1.267622, 1.197842, 1.131352, 1.068033, 1.007770, 0.950450, 0.895963, 0.844202, 0.795062, 0.748441, 0.704241, 0.662363, 0.622715, 0.585204, 0.549742, 0.516242, 0.484620, 0.454795, 0.426686, 0.400217, 0.375314, 0.351903, 0.329915, 0.309281, 0.289936, 0.271816, 0.254860, 0.239007, 0.224201, 0.210386, 0.197508, 0.185516, 0.174360, 0.163992, 0.154366, 0.145438, 0.137166, 0.129508, 0.122426, 0.115882, 0.109840, 0.104266, 0.099128, 0.094395, 0.090036, 0.086023, 0.082331, 0.078933, 0.075805, 0.072925, 0.070271, 0.067823, 0.065562, 0.063471, 0.061532, 0.059730, 0.058051, 0.056481, 0.055007, 0.053619, 0.052306, 0.051059, 0.049867, 0.048725, 0.047624, 0.046558, 0.045522, 0.044511, 0.043521, 0.042548, 0.041590, 0.040643, 0.039706, 0.038778, 0.037857, 0.036944, 0.036039, 0.035141, 0.034251, 0.033370, 0.032500, 0.031641, 0.030796, 0.029966, 0.029153, 0.028359, 0.027587, 0.026837, 0.026113, 0.025416, 0.024748, 0.024111, 0.023507, 0.022937, 0.022402, 0.021904, 0.021443, 0.021020, 0.020634, 0.020286, 0.019974, 0.019698, 0.019455, 0.019244, 0.019062, 0.018905, 0.018770, 0.018652, 0.018545, 0.018444, 0.018342, 0.018231, 0.018102, 0.017947, 0.017755, 0.017536, 0.017327, 0.017120, 0.016915, 0.016713, 0.016514, 0.016317, 0.016122, 0.015929, 0.015739, 0.015551, 0.015366, 0.015182, 0.015001, 0.014822, 0.014645, 0.014470, 0.014297, 0.014127, 0.013958, 0.013791, 0.013627, 0.013464, 0.013303, 0.013145, 0.012988, 0.012833, 0.012679, 0.012528, 0.012378, 0.012231, 0.012085, 0.011940, 0.011798, 0.011657, 0.011518, 0.011380, 0.011244, 0.011110, 0.010978, 0.010846, 0.010717, 0.010589, 0.010463, 0.010338, 0.010214, 0.010092, 0.009972, 0.009853, 0.009735, 0.009619, 0.009504, 0.009391, 0.009279, 0.009168, 0.009058, 0.008950, 0.008843, 0.008738, 0.008633, 0.008530, 0.008429, 0.008328, 0.008229, 0.008130, 0.008033, 0.007937, 0.007843, 0.007749, 0.007656, 0.007565, 0.007475, 0.007385, 0.007297, 0.007210, 0.007124, 0.007039, 0.006955, 0.006872, 0.006790, 0.006709, 0.006629, 0.006550, 0.006471, 0.006394, 0.006318, 0.006242, 0.006168, 0.006094, 0.006022, 0.005950, 0.005879, 0.005808, 0.005739, 0.005671, 0.005603, 0.005536, 0.005470, 0.005405, 0.005340, 0.005276, 0.005213, 0.005151, 0.005090, 0.005029, 0.004969, 0.004909, 0.004851, 0.004793, 0.004736, 0.004679, 0.004623, 0.004568, 0.004514, 0.004460, 0.004406, 0.004354, 0.004302, 0.004251, 0.004200, 0.004150, 0.004100, 0.004051, 0.004003, 0.003955, 0.003908, 0.003861, 0.003815, 0.003770, 0.003725, 0.003680, 0.003636, 0.003593, 0.003550, 0.003507, 0.003466, 0.003424, 0.003383, 0.003343, 0.003303, 0.003264, 0.003225, 0.003186, 0.003148, 0.003110, 0.003073, 0.003037, 0.003000, 0.002965, 0.002929, 0.002894, 0.002860, 0.002826, 0.002792, 0.002758, 0.002726, 0.002693, 0.002661, 0.002629, 0.002598, 0.002567, 0.002536, 0.002506, 0.002476, 0.002446, 0.002417, 0.002388, 0.002360, 0.002332, 0.002304, 0.002276, 0.002249, 0.002222, 0.002196, 0.002169, 0.002144, 0.002118, 0.002093, 0.002068, 0.002043, 0.002019, 0.001995, 0.001971, 0.001947, 0.001924, 0.001901, 0.001878, 0.001856, 0.001834, 0.001812, 0.001790, 0.001769, 0.001748, 0.001727, 0.001706, 0.001686, 0.001666, 0.001646, 0.001626, 0.001607, 0.001588, 0.001569, 0.001550, 0.001531, 0.001513, 0.001495, 0.001477, 0.001460, 0.001442, 0.001425, 0.001408, 0.001391, 0.001374, 0.001358, 0.001342, 0.001326, 0.001310, 0.001294, 0.001279, 0.001264, 0.001249, 0.001234, 0.001219, 0.001204, 0.001190, 0.001176, 0.001162, 0.001148, 0.001134, 0.001121, 0.001107, 0.001094, 0.001081, 0.001068, 0.001055, 0.001043, 0.001030, 0.001018, 0.001006, 0.000994, 0.000982, 0.000970, 0.000959, 0.000947, 0.000936, 0.000925, 0.000914, 0.000903, 0.000892, 0.000881, 0.000871, 0.000860, 0.000850, 0.000840, 0.000830, 0.000820, 0.000810, 0.000801, 0.000791, 0.000782, 0.000772, 0.000763, 0.000754, 0.000745, 0.000736, 0.000727, 0.000719, 0.000710, 0.000702, 0.000693, 0.000685, 0.000677, 0.000669, 0.000661, 0.000653, 0.000645, 0.000637, 0.000630, 0.000622, 0.000615, 0.000607, 0.000600, 0.000593, 0.000586, 0.000579, 0.000572, 0.000565, 0.000558, 0.000552, 0.000545, 0.000539, 0.000532, 0.000526, 0.000520, 0.000513, 0.000507, 0.000501, 0.000495, 0.000489, 0.000483, 0.000478, 0.000472, 0.000466, 0.000461, 0.000455, 0.000450, 0.000444, 0.000439, 0.000434, 0.000429, 0.000424, 0.000419, 0.000414, 0.000409, 0.000404, 0.000399, 0.000394, 0.000389, 0.000385, 0.000380, 0.000376, 0.000371, 0.000367, 0.000362, 0.000358, 0.000354, 0.000350, 0.000345, 0.000341, 0.000337, 0.000333, 0.000329, 0.000325, 0.000321, 0.000318, 0.000314, 0.000310, 0.000306, 0.000303, 0.000299, 0.000295, 0.000292, 0.000288, 0.000285, 0.000282, 0.000278, 0.000275, 0.000272, 0.000268, 0.000265, 0.000262, 0.000259, 0.000256, 0.000253, 0.000250, 0.000247, 0.000244, 0.000241, 0.000238, 0.000235};
1694  for(Int_t i=0; i<500; i++){
1695  fHistoPtWeight->SetBinContent(i+1,binc[i]);
1696  }
1697  //SetWeightHistogram();
1698  fUseWeight=kTRUE;
1699 }
1700 //_________________________________________________________________________
1702  // weight function from the ratio of the LHC16i2a MC
1703  // using FONLL*Raa from Xu.Cao.Bass(LBT model)
1704 
1705  if(fHistoPtWeight) delete fHistoPtWeight;
1706  fHistoPtWeight = new TH1F("histoWeight","histoWeight",500,0.,50.);
1707  Float_t binc[500]={ -3.491009, -2.151267, -1.169817, -0.462113, 0.038852, 0.385901, 0.620406, 0.774498, 0.872869, 0.934235, 0.972536, 0.997905, 1.017465, 1.035981, 1.056389, 1.080233, 1.108019, 1.139503, 1.173928, 1.210210, 1.247090, 1.283254, 1.317421, 1.348413, 1.375204, 1.396948, 1.413000, 1.422920, 1.426469, 1.423596, 1.414424, 1.399224, 1.378394, 1.352429, 1.321897, 1.287417, 1.249627, 1.209171, 1.166675, 1.122735, 1.077901, 1.032675, 0.987496, 0.942748, 0.898751, 0.855766, 0.813999, 0.773607, 0.734700, 0.697350, 0.661596, 0.627448, 0.594895, 0.563911, 0.534455, 0.506477, 0.479923, 0.454732, 0.430845, 0.408201, 0.386738, 0.366399, 0.347127, 0.328866, 0.311564, 0.295172, 0.279641, 0.264928, 0.250988, 0.237781, 0.248506, 0.240198, 0.232287, 0.224750, 0.217563, 0.210707, 0.204160, 0.197906, 0.191927, 0.186209, 0.180735, 0.175494, 0.170472, 0.165657, 0.161039, 0.156606, 0.152350, 0.148261, 0.144331, 0.140552, 0.136916, 0.133416, 0.130046, 0.126800, 0.123671, 0.120654, 0.117744, 0.114937, 0.112226, 0.109608, 0.107080, 0.104636, 0.102273, 0.099988, 0.097778, 0.095638, 0.093567, 0.091561, 0.089618, 0.087735, 0.085910, 0.084141, 0.082424, 0.080759, 0.079143, 0.077574, 0.076051, 0.074571, 0.073134, 0.071737, 0.070379, 0.069059, 0.067775, 0.066526, 0.065311, 0.064128, 0.062977, 0.061856, 0.060764, 0.059700, 0.058664, 0.057654, 0.056670, 0.055710, 0.054775, 0.053862, 0.052972, 0.052103, 0.051255, 0.050428, 0.049620, 0.048831, 0.048061, 0.047308, 0.046573, 0.045855, 0.045152, 0.044466, 0.043795, 0.043139, 0.042498, 0.041870, 0.041256, 0.040656, 0.040068, 0.039493, 0.038930, 0.038379, 0.037839, 0.037310, 0.036793, 0.036286, 0.035789, 0.035302, 0.034825, 0.034358, 0.033900, 0.033451, 0.033010, 0.032578, 0.032155, 0.031740, 0.031332, 0.030932, 0.030540, 0.030155, 0.029778, 0.029407, 0.029043, 0.028686, 0.028335, 0.027991, 0.027653, 0.027320, 0.026994, 0.026674, 0.026359, 0.026049, 0.025745, 0.025447, 0.025153, 0.024864, 0.024581, 0.024302, 0.024027, 0.023758, 0.023492, 0.023232, 0.022975, 0.022723, 0.022474, 0.022230, 0.021990, 0.021754, 0.021521, 0.021292, 0.021066, 0.020845, 0.020626, 0.020411, 0.020199, 0.019991, 0.019786, 0.019584, 0.019385, 0.019189, 0.018995, 0.018805, 0.018618, 0.018433, 0.018251, 0.018072, 0.017895, 0.017721, 0.017549, 0.017380, 0.017214, 0.017049, 0.016887, 0.016727, 0.016570, 0.016415, 0.016262, 0.016111, 0.015962, 0.015815, 0.015670, 0.015527, 0.015386, 0.015247, 0.015110, 0.014974, 0.014841, 0.014709, 0.014579, 0.014450, 0.014324, 0.014199, 0.014075, 0.013953, 0.013833, 0.013714, 0.013597, 0.013481, 0.013367, 0.013254, 0.013143, 0.013033, 0.012924, 0.012817, 0.012711, 0.012606, 0.012503, 0.012401, 0.012300, 0.012200, 0.012102, 0.012004, 0.011908, 0.011813, 0.011719, 0.011626, 0.011535, 0.011444, 0.011355, 0.011266, 0.011179, 0.011092, 0.011007, 0.010922, 0.010839, 0.010756, 0.010675, 0.010594, 0.010514, 0.010435, 0.010357, 0.010280, 0.010204, 0.010128, 0.010053, 0.009980, 0.009907, 0.009834, 0.009763, 0.009692, 0.009622, 0.009553, 0.009485, 0.009417, 0.009350, 0.009284, 0.009218, 0.009153, 0.009089, 0.009026, 0.008963, 0.008901, 0.008839, 0.008778, 0.008718, 0.008658, 0.008599, 0.008540, 0.008483, 0.008425, 0.008368, 0.008312, 0.008257, 0.008202, 0.008147, 0.008093, 0.008040, 0.007987, 0.007934, 0.007883, 0.007831, 0.007780, 0.007730, 0.007680, 0.007631, 0.007582, 0.007533, 0.007485, 0.007438, 0.007390, 0.007344, 0.007298, 0.007252, 0.007206, 0.007161, 0.007117, 0.007073, 0.007029, 0.006986, 0.006943, 0.006900, 0.006858, 0.006817, 0.006775, 0.006734, 0.006694, 0.006653, 0.006613, 0.006574, 0.006535, 0.006496, 0.006457, 0.006419, 0.006381, 0.006344, 0.006307, 0.006270, 0.006233, 0.006197, 0.006161, 0.006126, 0.006090, 0.006055, 0.006021, 0.005986, 0.005952, 0.005918, 0.005885, 0.005852, 0.005819, 0.005786, 0.005754, 0.005722, 0.005690, 0.005658, 0.005627, 0.005596, 0.005565, 0.005534, 0.005504, 0.005474, 0.005444, 0.005415, 0.005385, 0.005356, 0.005327, 0.005299, 0.005270, 0.005242, 0.005214, 0.005186, 0.005159, 0.005132, 0.005104, 0.005078, 0.005051, 0.005025, 0.004998, 0.004972, 0.004946, 0.004921, 0.004895, 0.004870, 0.004845, 0.004820, 0.004796, 0.004771, 0.004747, 0.004723, 0.004699, 0.004675, 0.004651, 0.004628, 0.004605, 0.004582, 0.004559, 0.004536, 0.004514, 0.004491, 0.004469, 0.004447, 0.004425, 0.004404, 0.004382, 0.004361, 0.004339, 0.004318, 0.004297, 0.004277, 0.004256, 0.004235, 0.004215, 0.004195, 0.004175, 0.004155, 0.004135, 0.004116, 0.004096, 0.004077, 0.004058, 0.004039, 0.004020, 0.004001, 0.003982, 0.003964, 0.003945, 0.003927, 0.003909, 0.003891, 0.003873, 0.003855, 0.003837, 0.003820, 0.003802, 0.003785, 0.003768, 0.003751, 0.003734, 0.003717, 0.003700, 0.003684, 0.003667, 0.003651, 0.003634, 0.003618, 0.003602, 0.003586, 0.003570, 0.003555, 0.003539, 0.003523, 0.003508, 0.003493, 0.003477, 0.003462, 0.003447, 0.003432, 0.003417, 0.003403, 0.003388, 0.003373, 0.003359, 0.003345, 0.003330, 0.003316, 0.003302, 0.003288, 0.003274, 0.003260, 0.003246, 0.003233, 0.003219, 0.003206, 0.003192};
1708  for(Int_t i=0; i<500; i++){
1709  fHistoPtWeight->SetBinContent(i+1,binc[i]);
1710  }
1711  //SetWeightHistogram();
1712  fUseWeight=kTRUE;
1713 }
1714 //_________________________________________________________________________
1716  // weight function from the ratio of the LHC16i2a+b+c MC
1717  // and FONLL calculations for pp data
1718 
1719  if(fHistoPtWeight) delete fHistoPtWeight;
1720  fHistoPtWeight = new TH1F("histoWeight","histoWeight",400,0.,40.);
1721  Float_t binc[400]={1.118416, 1.003458, 0.935514, 0.907222, 0.904359, 0.913668, 0.933906, 0.963898, 0.996388, 1.031708, 1.066404, 1.099683, 1.125805, 1.145181, 1.165910, 1.181905, 1.193425, 1.203891, 1.204726, 1.209411, 1.209943, 1.204763, 1.205291, 1.198912, 1.197390, 1.182005, 1.184194, 1.175994, 1.167881, 1.158348, 1.147190, 1.139833, 1.126940, 1.123322, 1.108389, 1.102199, 1.089464, 1.075874, 1.061964, 1.051429, 1.038113, 1.026668, 1.011441, 0.998567, 0.987658, 0.972434, 0.950068, 0.940758, 0.916880, 0.911931, 0.894512, 0.878691, 0.860589, 0.848025, 0.830774, 0.819399, 0.801134, 0.775276, 0.766382, 0.750495, 0.736935, 0.717529, 0.702637, 0.689152, 0.671334, 0.652030, 0.635696, 0.621365, 0.608362, 0.599019, 0.576024, 0.562136, 0.550938, 0.533587, 0.516410, 0.509744, 0.501655, 0.487402, 0.476469, 0.463762, 0.445979, 0.438088, 0.422214, 0.417467, 0.404357, 0.391450, 0.379996, 0.371201, 0.361497, 0.352912, 0.343189, 0.329183, 0.327662, 0.310783, 0.304525, 0.301007, 0.293306, 0.278332, 0.274419, 0.267361, 0.261459, 0.255514, 0.249293, 0.241129, 0.237600, 0.231343, 0.221982, 0.216872, 0.211094, 0.206954, 0.202333, 0.196572, 0.193274, 0.188240, 0.181817, 0.178364, 0.173614, 0.167135, 0.166055, 0.163423, 0.156557, 0.155821, 0.151985, 0.144909, 0.145062, 0.139720, 0.138873, 0.131892, 0.129969, 0.126509, 0.126978, 0.120451, 0.117661, 0.116300, 0.115604, 0.112215, 0.109237, 0.107720, 0.106419, 0.102050, 0.102777, 0.097406, 0.098447, 0.095964, 0.093868, 0.092430, 0.089329, 0.088249, 0.085881, 0.084417, 0.085498, 0.082444, 0.079151, 0.079565, 0.077811, 0.077293, 0.075218, 0.072445, 0.073054, 0.071545, 0.070279, 0.068046, 0.067854, 0.068092, 0.065378, 0.064405, 0.062060, 0.063391, 0.061718, 0.059616, 0.058913, 0.058895, 0.058311, 0.056320, 0.056527, 0.055349, 0.053701, 0.054735, 0.052264, 0.051277, 0.051554, 0.050545, 0.048995, 0.049507, 0.048466, 0.048156, 0.046809, 0.047600, 0.046078, 0.044801, 0.044113, 0.043700, 0.043530, 0.043396, 0.042556, 0.041048, 0.041657, 0.040394, 0.041314, 0.040720, 0.039656, 0.038478, 0.039276, 0.038777, 0.037730, 0.036918, 0.036466, 0.035827, 0.035285, 0.035963, 0.034371, 0.034757, 0.033205, 0.033666, 0.033266, 0.032583, 0.033570, 0.032102, 0.032107, 0.031464, 0.032160, 0.030091, 0.030564, 0.029464, 0.029613, 0.029626, 0.029512, 0.029324, 0.028607, 0.027628, 0.027251, 0.027072, 0.027077, 0.026724, 0.026961, 0.026303, 0.026237, 0.025454, 0.025133, 0.025365, 0.026014, 0.024807, 0.023901, 0.023459, 0.023405, 0.023654, 0.023981, 0.023675, 0.022493, 0.022781, 0.021801, 0.021704, 0.022372, 0.021189, 0.020681, 0.020779, 0.021324, 0.020558, 0.020901, 0.020586, 0.020808, 0.019276, 0.019516, 0.019706, 0.018935, 0.018632, 0.018516, 0.019187, 0.018916, 0.018039, 0.018208, 0.018045, 0.017628, 0.017916, 0.017711, 0.017838, 0.017222, 0.016565, 0.015733, 0.016264, 0.015826, 0.016090, 0.016622, 0.015802, 0.016621, 0.015441, 0.015309, 0.014860, 0.014935, 0.014968, 0.014443, 0.014485, 0.015136, 0.014078, 0.014414, 0.013908, 0.014071, 0.014078, 0.013766, 0.013436, 0.013507, 0.013480, 0.013224, 0.013041, 0.013935, 0.012885, 0.012453, 0.012528, 0.012492, 0.012225, 0.012542, 0.012706, 0.012136, 0.011902, 0.011560, 0.011448, 0.011861, 0.011271, 0.011831, 0.011159, 0.011171, 0.010966, 0.011311, 0.011002, 0.011130, 0.010995, 0.010450, 0.010663, 0.010678, 0.010492, 0.009861, 0.010507, 0.009916, 0.010121, 0.010029, 0.010046, 0.009370, 0.009647, 0.010104, 0.009282, 0.009830, 0.009403, 0.009148, 0.009172, 0.008893, 0.009158, 0.009019, 0.008780, 0.008579, 0.009063, 0.008634, 0.008988, 0.008265, 0.008581, 0.008575, 0.008690, 0.008181, 0.008352, 0.008150, 0.008430, 0.008256, 0.008119, 0.008453, 0.008447, 0.008021, 0.007938, 0.008025, 0.007718, 0.008127, 0.007651, 0.007590, 0.007316, 0.007839, 0.007504, 0.007341, 0.007527, 0.007263, 0.007668, 0.007306, 0.007271, 0.006910, 0.007257, 0.007260, 0.006810, 0.006967, 0.006887, 0.006867, 0.007202, 0.006829, 0.006370, 0.006710, 0.006417, 0.006361, 0.006800, 0.006410, 0.006323, 0.006790, 0.006322, 0.006673, 0.006547};
1722  for(Int_t i=0; i<400; i++){
1723  fHistoPtWeight->SetBinContent(i+1,binc[i]);
1724  }
1725  //SetWeightHistogram();
1726  fUseWeight=kTRUE;
1727 }
1728 
1729 //_________________________________________________________________________
1731  // weight function from the ratio of the LHC16i2a+b+c MC
1732  // and FONLL calculations for pp data
1733  // corrected by the BAMPS Raa calculation for 30-50% CC
1734  if(fHistoPtWeight) delete fHistoPtWeight;
1735  fHistoPtWeight = new TH1F("histoWeight","histoWeight",400,0.,40.);
1736  Float_t binc[400]={2.166180, 1.866117, 1.667595, 1.547176, 1.486661, 1.457891, 1.426949, 1.399055, 1.383278, 1.349383, 1.317009, 1.282321, 1.234257, 1.181136, 1.136655, 1.087523, 1.037912, 0.993256, 0.944746, 0.900948, 0.865869, 0.827193, 0.794424, 0.757723, 0.733020, 0.700164, 0.682189, 0.659872, 0.637918, 0.615749, 0.593020, 0.574402, 0.556158, 0.542663, 0.525494, 0.516038, 0.503629, 0.490980, 0.479143, 0.469005, 0.457749, 0.447668, 0.436803, 0.427073, 0.418282, 0.407867, 0.395093, 0.387861, 0.374742, 0.369462, 0.360146, 0.351991, 0.342990, 0.336259, 0.327730, 0.322382, 0.314602, 0.303874, 0.299820, 0.293049, 0.287539, 0.280329, 0.274866, 0.269939, 0.263299, 0.256057, 0.249215, 0.242170, 0.235704, 0.230709, 0.220529, 0.213921, 0.208394, 0.202424, 0.196700, 0.194943, 0.192620, 0.187894, 0.184411, 0.180204, 0.172915, 0.169077, 0.162201, 0.159636, 0.153904, 0.148296, 0.143282, 0.139306, 0.135561, 0.132342, 0.128696, 0.123444, 0.122873, 0.116544, 0.114197, 0.112878, 0.110018, 0.104547, 0.103222, 0.100707, 0.098622, 0.096513, 0.094295, 0.091334, 0.090122, 0.087870, 0.084894, 0.083729, 0.082265, 0.081404, 0.080323, 0.078750, 0.078132, 0.076781, 0.074823, 0.074050, 0.072614, 0.070093, 0.069828, 0.068907, 0.066189, 0.066054, 0.064600, 0.061757, 0.061986, 0.059862, 0.059656, 0.056807, 0.055956, 0.054386, 0.054507, 0.051629, 0.050358, 0.049702, 0.049331, 0.047814, 0.046476, 0.045762, 0.045142, 0.043224, 0.043484, 0.041282, 0.041794, 0.040809, 0.039985, 0.039439, 0.038181, 0.037782, 0.036831, 0.036264, 0.036790, 0.035535, 0.034173, 0.034409, 0.033659, 0.033308, 0.032290, 0.030981, 0.031121, 0.030361, 0.029708, 0.028653, 0.028461, 0.028449, 0.027208, 0.026697, 0.025623, 0.026069, 0.025279, 0.024332, 0.024341, 0.024629, 0.024677, 0.024117, 0.024490, 0.024257, 0.023804, 0.024537, 0.023692, 0.023502, 0.023888, 0.023673, 0.023193, 0.023684, 0.023429, 0.023521, 0.023014, 0.023346, 0.022544, 0.021866, 0.021477, 0.021224, 0.021089, 0.020972, 0.020515, 0.019739, 0.019982, 0.019328, 0.019719, 0.019387, 0.018833, 0.018227, 0.018558, 0.018276, 0.017738, 0.017460, 0.017365, 0.017178, 0.017033, 0.017478, 0.016817, 0.017119, 0.016463, 0.016802, 0.016711, 0.016475, 0.017083, 0.016441, 0.016548, 0.016320, 0.016786, 0.015804, 0.016153, 0.015668, 0.015843, 0.015810, 0.015651, 0.015454, 0.014981, 0.014376, 0.014089, 0.013906, 0.013818, 0.013549, 0.013580, 0.013160, 0.013040, 0.012566, 0.012324, 0.012353, 0.012582, 0.011915, 0.011401, 0.011112, 0.011008, 0.011046, 0.011119, 0.010954, 0.010439, 0.010604, 0.010179, 0.010163, 0.010507, 0.009981, 0.009771, 0.009846, 0.010134, 0.009798, 0.009991, 0.009869, 0.010005, 0.009295, 0.009438, 0.009557, 0.009210, 0.009088, 0.009057, 0.009412, 0.009306, 0.008899, 0.009009, 0.008952, 0.008764, 0.008926, 0.008842, 0.008924, 0.008634, 0.008322, 0.007920, 0.008205, 0.008000, 0.008151, 0.008438, 0.008037, 0.008472, 0.007886, 0.007835, 0.007621, 0.007675, 0.007707, 0.007452, 0.007489, 0.007841, 0.007308, 0.007497, 0.007248, 0.007348, 0.007367, 0.007227, 0.007097, 0.007179, 0.007209, 0.007115, 0.007059, 0.007588, 0.007058, 0.006862, 0.006945, 0.006965, 0.006856, 0.007075, 0.007209, 0.006925, 0.006830, 0.006672, 0.006645, 0.006923, 0.006615, 0.006982, 0.006622, 0.006666, 0.006579, 0.006823, 0.006673, 0.006786, 0.006740, 0.006440, 0.006606, 0.006650, 0.006568, 0.006206, 0.006646, 0.006305, 0.006468, 0.006442, 0.006486, 0.006080, 0.006291, 0.006622, 0.006113, 0.006506, 0.006254, 0.006114, 0.006161, 0.006002, 0.006211, 0.006146, 0.006012, 0.005902, 0.006264, 0.005996, 0.006271, 0.005793, 0.006043, 0.006067, 0.006177, 0.005842, 0.005991, 0.005872, 0.006102, 0.006003, 0.005930, 0.006201, 0.006224, 0.005937, 0.005901, 0.005992, 0.005788, 0.006121, 0.005787, 0.005766, 0.005582, 0.006006, 0.005774, 0.005672, 0.005841, 0.005660, 0.006000, 0.005741, 0.005737, 0.005475, 0.005773, 0.005799, 0.005462, 0.005610, 0.005569, 0.005574, 0.005871, 0.005589, 0.005234, 0.005535, 0.005314, 0.005288, 0.005676, 0.005371, 0.005319, 0.005734, 0.005360, 0.005679, 0.005593};
1737  for(Int_t i=0; i<400; i++){
1738  fHistoPtWeight->SetBinContent(i+1,binc[i]);
1739  }
1740  fUseWeight=kTRUE;
1741 }
1742 
1743 //_________________________________________________________________________
1745  // weight function from the ratio of the LHC16i2a+b+c MC
1746  // and FONLL calculations for pp data
1747  // corrected by the TAMU Raa calculation for 0-10% CC (not available in 30-50% CC)
1748  if(fHistoPtWeight) delete fHistoPtWeight;
1749  fHistoPtWeight = new TH1F("histoWeight","histoWeight",400,0.,40.);
1750  Float_t binc[400]={1.179906, 1.091249, 1.047774, 1.045579, 1.071679, 1.112413, 1.167414, 1.236240, 1.310301, 1.390289, 1.471711, 1.553389, 1.626886, 1.692115, 1.760647, 1.813658, 1.850817, 1.886699, 1.907671, 1.934832, 1.955433, 1.966727, 1.987262, 1.996316, 2.013326, 1.973926, 1.931144, 1.871654, 1.812942, 1.752718, 1.690846, 1.635303, 1.572611, 1.523510, 1.459790, 1.402510, 1.331908, 1.261575, 1.192241, 1.127915, 1.061798, 0.998830, 0.933514, 0.871774, 0.812936, 0.762844, 0.719340, 0.686587, 0.644108, 0.615714, 0.579512, 0.545254, 0.510508, 0.479884, 0.447423, 0.426154, 0.408934, 0.388264, 0.376424, 0.361389, 0.347757, 0.331685, 0.318029, 0.305285, 0.290922, 0.278523, 0.269807, 0.262025, 0.254878, 0.249325, 0.238179, 0.230899, 0.224792, 0.216253, 0.207879, 0.204465, 0.201153, 0.195373, 0.190926, 0.185773, 0.178589, 0.175371, 0.168959, 0.167004, 0.161705, 0.156809, 0.152788, 0.149806, 0.146429, 0.143478, 0.140037, 0.134813, 0.134679, 0.128205, 0.126078, 0.125038, 0.122214, 0.116329, 0.115044, 0.112427, 0.110279, 0.108098, 0.105784, 0.102628, 0.101429, 0.099101, 0.095464, 0.093631, 0.091491, 0.090045, 0.088374, 0.086188, 0.085067, 0.083168, 0.080636, 0.079414, 0.077610, 0.075013, 0.074825, 0.073932, 0.071106, 0.071050, 0.069574, 0.066593, 0.066924, 0.064876, 0.065064, 0.062345, 0.061980, 0.060859, 0.061616, 0.058952, 0.058079, 0.057894, 0.058031, 0.056604, 0.055180, 0.054490, 0.053909, 0.051768, 0.052210, 0.049552, 0.050152, 0.048955, 0.047953, 0.047224, 0.045588, 0.044985, 0.043728, 0.042934, 0.043434, 0.041834, 0.040118, 0.040281, 0.039348, 0.038987, 0.037793, 0.036258, 0.036420, 0.035528, 0.034761, 0.033524, 0.033296, 0.033280, 0.031825, 0.031351, 0.030329, 0.031103, 0.030401, 0.029481, 0.029247, 0.029352, 0.029174, 0.028286, 0.028500, 0.028017, 0.027293, 0.027932, 0.026779, 0.026379, 0.026628, 0.026211, 0.025508, 0.025877, 0.025433, 0.025328, 0.024636, 0.025069, 0.024282, 0.023625, 0.023278, 0.023074, 0.023000, 0.022943, 0.022514, 0.021767, 0.022180, 0.021594, 0.022175, 0.021944, 0.021456, 0.020901, 0.021419, 0.021230, 0.020738, 0.020322, 0.020055, 0.019686, 0.019371, 0.019725, 0.018835, 0.019029, 0.018163, 0.018398, 0.018163, 0.017719, 0.018126, 0.017208, 0.017086, 0.016622, 0.016865, 0.015663, 0.015791, 0.015108, 0.015069, 0.015033, 0.015006, 0.014940, 0.014604, 0.014133, 0.013968, 0.013904, 0.013934, 0.013780, 0.013930, 0.013727, 0.013940, 0.013763, 0.013826, 0.014192, 0.014801, 0.014347, 0.014048, 0.014009, 0.014197, 0.014571, 0.014999, 0.015030, 0.014491, 0.014891, 0.014456, 0.014596, 0.015256, 0.014648, 0.014492, 0.014756, 0.015344, 0.014986, 0.015433, 0.015394, 0.015756, 0.014778, 0.015145, 0.015478, 0.015051, 0.014986, 0.015067, 0.015793, 0.015748, 0.015188, 0.015502, 0.015533, 0.015340, 0.015759, 0.015745, 0.016026, 0.015635, 0.015194, 0.014579, 0.015225, 0.014963, 0.015365, 0.016030, 0.015387, 0.016341, 0.015327, 0.015340, 0.015030, 0.015246, 0.015420, 0.015015, 0.015195, 0.016021, 0.015034, 0.015528, 0.015114, 0.015423, 0.015564, 0.015348, 0.015107, 0.015314, 0.015411, 0.015243, 0.015154, 0.016324, 0.015215, 0.014823, 0.015030, 0.015104, 0.014896, 0.015400, 0.015721, 0.015131, 0.014951, 0.014630, 0.014597, 0.015235, 0.014583, 0.015418, 0.014648, 0.014769, 0.014601, 0.015167, 0.014857, 0.015134, 0.015053, 0.014405, 0.014800, 0.014921, 0.014760, 0.013966, 0.014979, 0.014230, 0.014620, 0.014581, 0.014701, 0.013799, 0.014299, 0.015071, 0.013931, 0.014846, 0.014290, 0.013988, 0.014113, 0.013767, 0.014263, 0.014131, 0.013840, 0.013604, 0.014456, 0.013853, 0.014505, 0.013416, 0.014010, 0.014081, 0.014352, 0.013589, 0.013952, 0.013690, 0.014241, 0.014024, 0.013868, 0.014517, 0.014587, 0.013927, 0.013857, 0.014084, 0.013619, 0.014417, 0.013644, 0.013607, 0.013185, 0.014200, 0.013665, 0.013437, 0.013849, 0.013431, 0.014252, 0.013648, 0.013652, 0.013039, 0.013761, 0.013836, 0.013043, 0.013408, 0.013319, 0.013344, 0.014065, 0.013400, 0.012560, 0.013294, 0.012773, 0.012721, 0.013663, 0.012939, 0.012823, 0.013835, 0.012942, 0.013723, 0.013525};
1751  for(Int_t i=0; i<400; i++){
1752  fHistoPtWeight->SetBinContent(i+1,binc[i]);
1753  }
1754  fUseWeight=kTRUE;
1755 }
1756 
1757 //_________________________________________________________________________
1759  // weight function from the ratio of the LHC16i2a MC
1760  // 2.5-14 GeV/c using data and 0-2.5, 14-50 GeV/c using FONLL calculations
1761 
1762  if(fHistoPtWeight) delete fHistoPtWeight;
1763  fHistoPtWeight = new TH1F("histoWeight","histoWeight",500,0.,50.);
1764  Float_t binc[500]={ 1.977386, 1.092380, 0.982069, 0.990772, 1.031869, 1.075578, 1.121814, 1.160738, 1.194926, 1.223790, 1.248790, 1.265487,
1765  1.280118, 1.292664, 1.302578, 1.308750, 1.316914, 1.321336, 1.325862, 1.335865, 1.330973, 1.338854, 1.347197, 1.355202,
1766  1.345868, 1.306835, 1.217713, 1.150661, 1.082328, 1.016336, 0.961361, 0.909377, 0.859229, 0.808615, 0.770923, 0.726400,
1767  0.690444, 0.656639, 0.626977, 0.601282, 0.565640, 0.541549, 0.517613, 0.492610, 0.467891, 0.450866, 0.427611, 0.405115,
1768  0.393432, 0.367097, 0.355256, 0.344176, 0.327230, 0.308159, 0.299667, 0.288000, 0.278949, 0.263341, 0.250445, 0.243064,
1769  0.232859, 0.221800, 0.214465, 0.207126, 0.197129, 0.188575, 0.186439, 0.179715, 0.170882, 0.164423, 0.158682, 0.154019,
1770  0.150277, 0.144272, 0.137206, 0.135669, 0.128931, 0.125974, 0.121932, 0.119566, 0.113806, 0.112343, 0.107325, 0.106804,
1771  0.104220, 0.098339, 0.095757, 0.095992, 0.092098, 0.089896, 0.087719, 0.086209, 0.085378, 0.082050, 0.078781, 0.078904,
1772  0.077688, 0.077198, 0.075539, 0.071257, 0.071342, 0.070415, 0.069445, 0.067272, 0.064635, 0.065670, 0.063817, 0.062011,
1773  0.060226, 0.058354, 0.058203, 0.057580, 0.055477, 0.055900, 0.054993, 0.053452, 0.052578, 0.051701, 0.050199, 0.049844,
1774  0.047829, 0.047402, 0.046873, 0.044984, 0.045633, 0.044068, 0.043317, 0.042371, 0.042091, 0.042391, 0.041007, 0.039471,
1775  0.038575, 0.038065, 0.039097, 0.037241, 0.035288, 0.035440, 0.035264, 0.033790, 0.033902, 0.033853, 0.032486, 0.031943,
1776  0.031140, 0.030939, 0.030288, 0.029462, 0.029969, 0.028097, 0.028703, 0.028627, 0.027646, 0.027564, 0.026505, 0.026356,
1777  0.025599, 0.026124, 0.025300, 0.024600, 0.024092, 0.024186, 0.023315, 0.023520, 0.023087, 0.022226, 0.022345, 0.021455,
1778  0.020492, 0.021519, 0.020416, 0.020084, 0.020114, 0.019717, 0.019508, 0.018334, 0.018877, 0.018658, 0.018210, 0.017702,
1779  0.018155, 0.017342, 0.016846, 0.017155, 0.016306, 0.016247, 0.016075, 0.015419, 0.015835, 0.016264, 0.015473, 0.015212,
1780  0.014994, 0.014338, 0.015062, 0.014039, 0.014740, 0.014412, 0.014111, 0.013556, 0.013372, 0.012951, 0.012977, 0.013086,
1781  0.012727, 0.012118, 0.012176, 0.012056, 0.012350, 0.011744, 0.011840, 0.011740, 0.012028, 0.011028, 0.010922, 0.011132,
1782  0.011492, 0.010939, 0.011151, 0.010446, 0.010507, 0.009898, 0.010788, 0.010553, 0.010305, 0.010137, 0.009395, 0.009406,
1783  0.009338, 0.009635, 0.009477, 0.009497, 0.008915, 0.009116, 0.009080, 0.009294, 0.008765, 0.008222, 0.008483, 0.008666,
1784  0.008659, 0.008378, 0.008100, 0.008121, 0.008041, 0.007868, 0.008032, 0.008100, 0.008079, 0.007581, 0.007421, 0.007118,
1785  0.007416, 0.007828, 0.007152, 0.007088, 0.007122, 0.007133, 0.006934, 0.007103, 0.006555, 0.007113, 0.006829, 0.006883,
1786  0.006288, 0.006291, 0.006312, 0.006090, 0.006457, 0.006410, 0.006449, 0.005728, 0.006082, 0.006032, 0.005872, 0.005989,
1787  0.005934, 0.005962, 0.005698, 0.005851, 0.005615, 0.005442, 0.005838, 0.005589, 0.005422, 0.005286, 0.005219, 0.005153,
1788  0.005222, 0.005520, 0.004971, 0.005049, 0.005134, 0.004912, 0.004870, 0.004990, 0.005039, 0.004782, 0.004498, 0.004962,
1789  0.004755, 0.004995, 0.004437, 0.004686, 0.004753, 0.004703, 0.004486, 0.004614, 0.004152, 0.004107, 0.004131, 0.004349,
1790  0.004254, 0.004112, 0.004103, 0.004045, 0.004095, 0.004148, 0.004417, 0.004104, 0.003765, 0.003868, 0.004080, 0.003937,
1791  0.004091, 0.003892, 0.003849, 0.003681, 0.003834, 0.003530, 0.003934, 0.003552, 0.003714, 0.003282, 0.003898, 0.003909,
1792  0.003685, 0.003800, 0.003469, 0.003311, 0.003483, 0.003242, 0.003083, 0.003145, 0.003158, 0.003073, 0.003331, 0.003388,
1793  0.003156, 0.002935, 0.003256, 0.003261, 0.002972, 0.002966, 0.003312, 0.003253, 0.003122, 0.003075, 0.002978, 0.003029,
1794  0.002883, 0.002793, 0.003039, 0.002905, 0.002511, 0.002779, 0.002964, 0.002755, 0.002747, 0.002863, 0.002777, 0.002514,
1795  0.002775, 0.002845, 0.002427, 0.002540, 0.002406, 0.002490, 0.002592, 0.002457, 0.002560, 0.002223, 0.002582, 0.002725,
1796  0.002456, 0.002258, 0.002419, 0.002630, 0.002263, 0.002158, 0.002518, 0.002364, 0.002034, 0.002077, 0.002132, 0.002602,
1797  0.002309, 0.002357, 0.002333, 0.002194, 0.002213, 0.002270, 0.002202, 0.002308, 0.002042, 0.002417, 0.002220, 0.002359,
1798  0.002305, 0.002082, 0.002036, 0.002075, 0.002176, 0.002291, 0.002174, 0.002057, 0.001983, 0.001972, 0.001997, 0.001875,
1799  0.001958, 0.001838, 0.001895, 0.001868, 0.001996, 0.001807, 0.001765, 0.002055, 0.001933, 0.001797, 0.001908, 0.001816,
1800  0.001790, 0.001774, 0.001731, 0.002141, 0.001511, 0.001603, 0.002065, 0.001799, 0.001736, 0.001954, 0.001836, 0.001729,
1801  0.001817, 0.001671, 0.001766, 0.001652, 0.001706, 0.001519, 0.001649, 0.001518, 0.001505, 0.001552, 0.001531, 0.001730,
1802  0.001717, 0.001442, 0.001626, 0.001668, 0.001591, 0.001509, 0.001313, 0.001555, 0.001393, 0.001587, 0.001363, 0.001360,
1803  0.001785, 0.001410, 0.001363, 0.001278, 0.001427, 0.001407, 0.001102, 0.001358, 0.001408, 0.001509, 0.001399, 0.001411,
1804  0.001391, 0.001284, 0.001641, 0.001205, 0.001172, 0.001277, 0.001163, 0.001321, 0.001243, 0.001284, 0.001353, 0.001219,
1805  0.001338, 0.001204, 0.001206, 0.001250, 0.001127, 0.001366, 0.001240, 0.001124 };
1806 
1807  for(Int_t i=0; i<500; i++){
1808  fHistoPtWeight->SetBinContent(i+1,binc[i]);
1809  }
1810  fUseWeight=kTRUE;
1811 }
1812 
1813 //_________________________________________________________________________
1815  //
1816  // weight function from the ratio of the D0 0-80% over lhc16i2abc
1817  // This is for Lc/D0 analysis, so the functional form is derived using only 3-16 GeV/c
1818  //
1819  if(fFuncWeight) delete fFuncWeight;
1820  fFuncWeight=new TF1("funcWeight","([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]*x+[6])",3.,16.);
1821  fFuncWeight->SetParameters(-93.179,2.15711,6.45459,1.33679,0.373436,-0.070921,2.52372);
1822  fUseWeight=kTRUE;
1823 }
1824 
1825 //_________________________________________________________________________
1827  //
1828  // weight function from the ratio of the D0 0-80% over lhc16i2abc
1829  // Modified taking into account Lc/D0 ratio given in PRC 79, 044905 (2009)
1830  // This is for Lc/D0 analysis, so the functional form is derived using only 3-16 GeV/c
1831  //
1832  if(fFuncWeight) delete fFuncWeight;
1833  fFuncWeight=new TF1("funcWeight","(([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]*x+[6]))*[7]*TMath::Gaus(x,[8],[9])",3.,16.);
1834  fFuncWeight->SetParameters(-93.179,2.15711,6.45459,1.33679,0.373436,-0.070921,2.52372,1.02369,1.41492,3.09519);
1835  fUseWeight=kTRUE;
1836 }
1837 
1838 //_________________________________________________________________________
1840  //
1841  // weight function from the ratio of the D0 0-80% over lhc16i2abc
1842  // Modified taking into account Lc/D0 ratio given in PLB 663, 55-60 (2008)
1843  // This is for Lc/D0 analysis, so the functional form is derived using only 3-16 GeV/c
1844  //
1845  if(fFuncWeight) delete fFuncWeight;
1846  fFuncWeight=new TF1("funcWeight","(([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]*x+[6]))*[7]*TMath::Gaus(x,[8],[9])",3.,16.);
1847  fFuncWeight->SetParameters(-93.179,2.15711,6.45459,1.33679,0.373436,-0.070921,2.52372,0.9,5.0,2.9);
1848  fUseWeight=kTRUE;
1849 }
1850 
1851 //_________________________________________________________________________
1853  // weight function from the ratio of the LHC13d3 MC
1854  // and FONLL calculations for pp data
1855  if(fFuncWeight) delete fFuncWeight;
1856  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.);
1857  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);
1858  fUseWeight=kTRUE;
1859 }
1860 
1861 //_________________________________________________________________________
1863  // weight function from the ratio of the LHC10f6a MC
1864  // and FONLL calculations for pp data
1865  if(fFuncWeight) delete fFuncWeight;
1866  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.);
1867  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);
1868  fUseWeight=kTRUE;
1869 }
1870 
1871 //_________________________________________________________________________
1873  // weight function from the ratio of the LHC12a12 MC
1874  // and FONLL calculations for pp data
1875  if(fFuncWeight) delete fFuncWeight;
1876  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.);
1877  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);
1878  fUseWeight=kTRUE;
1879 }
1880 
1881 //_________________________________________________________________________
1883  // weight function from the ratio of the LHC12a12bis MC
1884  // and FONLL calculations for pp data
1885  if(fFuncWeight) delete fFuncWeight;
1886  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.);
1887  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);
1888  fUseWeight=kTRUE;
1889 }
1890 
1891 //_________________________________________________________________________
1893  // weight function from the ratio of the LHC13e2fix MC
1894  // and FONLL calculations for pp data
1895  if(fFuncWeight) delete fFuncWeight;
1896  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.);
1897  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);
1898  fUseWeight=kTRUE;
1899 }
1900 
1901 //________________________________________________________________
1903  // weight function from the ratio of the LHC10f6a MC
1904  // and FONLL calculations for pp data
1905  if(fFuncWeight) delete fFuncWeight;
1906  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.);
1907  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);
1908  fUseWeight=kTRUE;
1909 }
1910 
1911 //________________________________________________________________
1913  // weight function from the ratio of the LHC10f6a MC
1914  // and FONLL calculations for pp data
1915  if(fFuncWeight) delete fFuncWeight;
1916  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.);
1917  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);
1918  fUseWeight=kTRUE;
1919 }
1920 
1921 //_________________________________________________________________________
1923  // weight function from the ratio of the LHC13d3 MC
1924  // and FONLL calculations for pPb data
1925  // using Lc simulated pt distribution
1926  if(fFuncWeight) delete fFuncWeight;
1927  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,20.);
1928  fFuncWeight->SetParameters(5.94428e+01,1.63585e+01,9.65555e+00,6.71944e+00,8.88338e-02,2.40477e+00,-4.88649e-02,-6.78599e-01,-2.10951e-01);
1929  fUseWeight=kTRUE;
1930 }
1931 //_________________________________________________________________________
1933  // weight function from the ratio of the LHC16i6a MC
1934  // and FONLL calculations for pp data
1935  // using D0 simulated pt distribution
1936  if(fFuncWeight) delete fFuncWeight;
1937  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.2,20.);
1938  fFuncWeight->SetParameters(2.70821e+03,2.98122e+00,8.67776e+01,1.16611e+00,1.17276e-02,5.41670e+00,6.01099e-02,-2.04524e+00,6.69208e-02);
1939  fUseWeight=kTRUE;
1940 }
1941 
1942 //_________________________________________________________________________
1944  // weight function from the ratio of the LHC11b2 MC
1945  // and FONLL calculations for pp data
1946  // using Lc simulated pt distribution
1947  if(fFuncWeight) delete fFuncWeight;
1948  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)",1.,20.);
1949  fFuncWeight->SetParameters(2.11879e+02,3.73290e+00,2.01235e+01,1.41508e+00,1.06268e-01,1.86285e+00,-4.52956e-02,-9.90631e-01,-1.31615e+00);
1950  fUseWeight=kTRUE;
1951 }
1952 
1953 //_________________________________________________________________________
1955  // weight function from the ratio of the LHC10f7a MC
1956  // and FONLL calculations for pp data
1957  // using Lc simulated pt distribution
1958  if(fFuncWeight) delete fFuncWeight;
1959  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)",1.,20.);
1960  fFuncWeight->SetParameters(2.84268e+02,2.18850e+01,2.36298e+01,7.46144e+00,1.69747e-01,1.66993e+00,-5.54726e-02,-1.53869e+00,-1.18404e+00);
1961  fUseWeight=kTRUE;
1962 }
1963 
1964 
1965 //_________________________________________________________________________
1967  // weight function from the ratio of the LHC15l2a2 MC
1968  // and FONLL calculations for pp data
1969  // using D0 simulated pt distribution
1970  if(fFuncWeight) delete fFuncWeight;
1971  fFuncWeight=new TF1("funcWeight","(x<28.5)*(([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x))+(x>=28.5)*0.140023",0.15,40);
1972  fFuncWeight->SetParameters(4.66092e+01, 4.27321e+01, 4.46858e+01, 1.67788e+00, -2.88457e-02, 4.40656e+00, -7.31064e-01, 2.96431e+00, -2.79976e-01);
1973  fUseWeight=kTRUE;
1974 }
1975 
1976 
1977 //_________________________________________________________________________
1979  // weight function from the ratio of the LHC17c3a1-LHC17c3a2 MC
1980  // and FONLL calculations for pp data
1981  // using D0 simulated pt distribution
1982  if(fFuncWeight) delete fFuncWeight;
1983  fFuncWeight=new TF1("funcWeight","(x<10.71)*(([0]*x)/TMath::Power([2],(1+TMath::Power([3],x/[1])))+[4]*TMath::Exp([5]+[6]*x)+[7]*TMath::Exp([8]*x))+(x>=10.71)*2.192",0.15,40);
1984  fFuncWeight->SetParameters(-8.28294e+04,1.04367e+00,1.43041e+02,1.40041e+00,-7.23198e-03,4.71629e+00,-3.15545e+00,2.36495e+00,-7.08928e-03);
1985  fUseWeight=kTRUE;
1986 }
1987 
1988 //_________________________________________________________________________
1990  // Weight function from the ratio of the D0 spectrum in LHC18a4a2 (fast+cent) MC
1991  // and FONLL calculations for pp data at 5 TeV
1992 
1993  if(fHistoPtWeight) delete fHistoPtWeight;
1994  fHistoPtWeight = new TH1F("histoWeight","histoWeight",500,0.,50.);
1995  fHistoPtWeight->Sumw2();
1996  Float_t binc[500]={0.7836, 0.7203, 0.6840, 0.6749, 0.6840, 0.7038, 0.7337, 0.7691, 0.8098, 0.8520, 0.8917, 0.9331, 0.9673, 0.9985, 1.0297, 1.0569, 1.0781, 1.0993, 1.1166, 1.1285, 1.1425, 1.1520, 1.1598, 1.1685, 1.1762, 1.1791, 1.1858, 1.1899, 1.1944, 1.1986, 1.1968, 1.2041, 1.2072, 1.2031, 1.2068, 1.2090, 1.2079, 1.2089, 1.2081, 1.2047, 1.2088, 1.2093, 1.2085, 1.2105, 1.2099, 1.2125, 1.2108, 1.2090, 1.2071, 1.2066, 1.2122, 1.2126, 1.2131, 1.2136, 1.2141, 1.2145, 1.2150, 1.2155, 1.2160, 1.2165, 1.2169, 1.2174, 1.2179, 1.2184, 1.2188, 1.2193, 1.2198, 1.2203, 1.2207, 1.2212, 1.2217, 1.2222, 1.2227, 1.2231, 1.2236, 1.2241, 1.2246, 1.2250, 1.2255, 1.2260, 1.2265, 1.2269, 1.2274, 1.2279, 1.2284, 1.2289, 1.2293, 1.2298, 1.2303, 1.2308, 1.2312, 1.2317, 1.2322, 1.2327, 1.2331, 1.2336, 1.2341, 1.2346, 1.2351, 1.2355, 1.2360, 1.2365, 1.2370, 1.2374, 1.2379, 1.2384, 1.2389, 1.2393, 1.2398, 1.2403, 1.2408, 1.2413, 1.2417, 1.2422, 1.2427, 1.2432, 1.2436, 1.2441, 1.2446, 1.2451, 1.2455, 1.2460, 1.2465, 1.2470, 1.2475, 1.2479, 1.2484, 1.2489, 1.2494, 1.2498, 1.2503, 1.2508, 1.2513, 1.2517, 1.2522, 1.2527, 1.2532, 1.2537, 1.2541, 1.2546, 1.2551, 1.2556, 1.2560, 1.2565, 1.2570, 1.2575, 1.2579, 1.2584, 1.2589, 1.2594, 1.2599, 1.2603, 1.2608, 1.2613, 1.2618, 1.2622, 1.2627, 1.2632, 1.2637, 1.2641, 1.2646, 1.2651, 1.2656, 1.2661, 1.2665, 1.2670, 1.2675, 1.2680, 1.2684, 1.2689, 1.2694, 1.2699, 1.2703, 1.2708, 1.2713, 1.2718, 1.2723, 1.2727, 1.2732, 1.2737, 1.2742, 1.2746, 1.2751, 1.2756, 1.2761, 1.2765, 1.2770, 1.2775, 1.2780, 1.2785, 1.2789, 1.2794, 1.2799, 1.2804, 1.2808, 1.2813, 1.2818, 1.2823, 1.2827, 1.2832, 1.2837, 1.2842, 1.2847, 1.2851, 1.2856, 1.2861, 1.2866, 1.2870, 1.2875, 1.2880, 1.2885, 1.2889, 1.2894, 1.2899, 1.2904, 1.2909, 1.2913, 1.2918, 1.2923, 1.2928, 1.2932, 1.2937, 1.2942, 1.2947, 1.2951, 1.2956, 1.2961, 1.2966, 1.2971, 1.2975, 1.2980, 1.2985, 1.2990, 1.2994, 1.2999, 1.3004, 1.3009, 1.3013, 1.3018, 1.3023, 1.3028, 1.3033, 1.3037, 1.3042, 1.3047, 1.3052, 1.3056, 1.3061, 1.3066, 1.3071, 1.3075, 1.3080, 1.3085, 1.3090, 1.3095, 1.3099, 1.3104, 1.3109, 1.3114, 1.3118, 1.3123, 1.3128, 1.3133, 1.3137, 1.3142, 1.3147, 1.3152, 1.3157, 1.3161, 1.3166, 1.3171, 1.3176, 1.3180, 1.3185, 1.3190, 1.3195, 1.3199, 1.3204, 1.3209, 1.3214, 1.3219, 1.3223, 1.3228, 1.3233, 1.3238, 1.3242, 1.3247, 1.3252, 1.3257, 1.3262, 1.3266, 1.3271, 1.3276, 1.3281, 1.3285, 1.3290, 1.3295, 1.3300, 1.3304, 1.3309, 1.3314, 1.3319, 1.3324, 1.3328, 1.3333, 1.3338, 1.3343, 1.3347, 1.3352, 1.3357, 1.3362, 1.3366, 1.3371, 1.3376, 1.3381, 1.3386, 1.3390, 1.3395, 1.3400, 1.3405, 1.3409, 1.3414, 1.3419, 1.3424, 1.3428, 1.3433, 1.3438, 1.3443, 1.3448, 1.3452, 1.3457, 1.3462, 1.3467, 1.3471, 1.3476, 1.3481, 1.3486, 1.3490, 1.3495, 1.3500, 1.3505, 1.3510, 1.3514, 1.3519, 1.3524, 1.3529, 1.3533, 1.3538, 1.3543, 1.3548, 1.3552, 1.3557, 1.3562, 1.3567, 1.3572, 1.3576, 1.3581, 1.3586, 1.3591, 1.3595, 1.3600, 1.3605, 1.3610, 1.3614, 1.3619, 1.3624, 1.3629, 1.3634, 1.3638, 1.3643, 1.3648, 1.3653, 1.3657, 1.3662, 1.3667, 1.3672, 1.3676, 1.3681, 1.3686, 1.3691, 1.3696, 1.3700, 1.3705, 1.3710, 1.3715, 1.3719, 1.3724, 1.3729, 1.3734, 1.3738, 1.3743, 1.3748, 1.3753, 1.3758, 1.3762, 1.3767, 1.3772, 1.3777, 1.3781, 1.3786, 1.3791, 1.3796, 1.3800, 1.3805, 1.3810, 1.3815, 1.3820, 1.3824, 1.3829, 1.3834, 1.3839, 1.3843, 1.3848, 1.3853, 1.3858, 1.3862, 1.3867, 1.3872, 1.3877, 1.3882, 1.3886, 1.3891, 1.3896, 1.3901, 1.3905, 1.3910, 1.3915, 1.3920, 1.3924, 1.3929, 1.3934, 1.3939, 1.3944, 1.3948, 1.3953, 1.3958, 1.3963, 1.3967, 1.3972, 1.3977, 1.3982, 1.3986, 1.3991, 1.3996, 1.4001, 1.4006, 1.4010, 1.4015, 1.4020, 1.4025, 1.4029, 1.4034, 1.4039, 1.4044, 1.4048, 1.4053, 1.4058, 1.4063, 1.4068, 1.4072, 1.4077, 1.4082, 1.4087, 1.4091, 1.4096, 1.4101, 1.4106, 1.4110, 1.4115, 1.4120, 1.4125, 1.4130, 1.4134, 1.4139, 1.4144, 1.4149, 1.4153, 1.4158, 1.4163, 1.4168, 1.4172, 1.4177, 1.4182, 1.4187, 1.4192, 1.4196, 1.4201, 1.4206, 1.4211, 1.4215, 1.4220, 1.4225, 1.4230, 1.4234, 1.4239, 1.4244, 1.4249, 1.4254, 1.4258, 1.4263};
1997  for(Int_t i=0; i<500; i++){
1998  fHistoPtWeight->SetBinContent(i+1,binc[i]);
1999  }
2000  fUseWeight=kTRUE;
2001 }
2002 
2003 //_________________________________________________________________________
2005 {
2006  //
2007  // calculating the weight to fill the container
2008  //
2009 
2010  // FNOLL central:
2011  // p0 = 1.63297e-01 --> 0.322643
2012  // p1 = 2.96275e+00
2013  // p2 = 2.30301e+00
2014  // p3 = 2.50000e+00
2015 
2016  // PYTHIA
2017  // p0 = 1.85906e-01 --> 0.36609
2018  // p1 = 1.94635e+00
2019  // p2 = 1.40463e+00
2020  // p3 = 2.50000e+00
2021 
2022  Double_t func1[4] = {0.322643,2.96275,2.30301,2.5};
2023  Double_t func2[4] = {0.36609,1.94635,1.40463,2.5};
2024 
2025  Double_t dndpt_func1 = dNdptFit(pt,func1);
2026  if(fUseFlatPtWeight) dndpt_func1 = 1./30.;
2027  Double_t dndpt_func2 = dNdptFit(pt,func2);
2028  AliDebug(2,Form("pt = %f, FONLL = %f, Pythia = %f, ratio = %f",pt,dndpt_func1,dndpt_func2,dndpt_func1/dndpt_func2));
2029  return dndpt_func1/dndpt_func2;
2030 }
2031 
2032 //__________________________________________________________________________________________________
2034 {
2035  //
2036  // calculating dNdpt
2037  //
2038 
2039  Double_t denom = TMath::Power((pt/par[1]), par[3] );
2040  Double_t dNdpt = par[0]*pt/TMath::Power(1.+denom, par[2]);
2041 
2042  return dNdpt;
2043 }
2044 
2045 //_________________________________________________________________________
2047 {
2048  //
2049  // Using an histogram as weight function
2050  // weight = 0 in the range outside the function
2051  //
2052  Double_t weight = 0.0;
2053  Int_t histoNbins = fHistoPtWeight->GetNbinsX();
2054  Int_t bin2 = fHistoPtWeight->FindBin(pt);
2055  if( (bin2>0) && (bin2<=histoNbins) ) {
2056  Int_t bin1=bin2-1;
2057  Int_t bin3=bin2+1;
2058  if(bin2==1) bin1=bin2+2;
2059  if(bin2==histoNbins) bin3=bin2-2;
2060  Float_t x_1=fHistoPtWeight->GetXaxis()->GetBinCenter(bin1);
2061  Float_t x_2=fHistoPtWeight->GetXaxis()->GetBinCenter(bin2);
2062  Float_t x_3=fHistoPtWeight->GetXaxis()->GetBinCenter(bin3);
2063  Float_t y_1=fHistoPtWeight->GetBinContent(bin1);
2064  Float_t y_2=fHistoPtWeight->GetBinContent(bin2);
2065  Float_t y_3=fHistoPtWeight->GetBinContent(bin3);
2066  Double_t a=( (y_3-y_2)*(x_1-x_2) - (y_1-y_2)*(x_3-x_2) )/( (x_3*x_3-x_2*x_2)*(x_1-x_2) - (x_1*x_1-x_2*x_2)*(x_3-x_2) );
2067  Double_t b=((y_1-y_2)-a*(x_1*x_1-x_2*x_2))/(x_1-x_2);
2068  Double_t c=y_3-a*(x_3*x_3)-b*x_3;
2069  weight = a*pt*pt+b*pt+c;
2070  }
2071  return weight;
2072 }
2073 
2074 //__________________________________________________________________________________________________
2076  //
2077  // calculating the z-vtx weight for the given run range
2078  //
2079 
2080  if(runnumber>146824 || runnumber<146803) return 1.0;
2081 
2082  Double_t func1[3] = {1.0, -0.5, 6.5 };
2083  Double_t func2[3] = {1.0, -0.5, 5.5 };
2084 
2085  Double_t dzFunc1 = DodzFit(z,func1);
2086  Double_t dzFunc2 = DodzFit(z,func2);
2087 
2088  return dzFunc1/dzFunc2;
2089 }
2090 
2091 //__________________________________________________________________________________________________
2093 
2094  //
2095  // Gaussian z-vtx shape
2096  //
2097  //gaussian = [0]/TMath::Sqrt(2.*TMath::Pi())/[2]*exp[-(x-[1])*(x-[1])/(2*[2]*[2])]
2098 
2099  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]);
2100 
2101  return value;
2102 }
2103 //__________________________________________________________________________________________________
2105  //
2106  // calculates the Nch weight using the measured and generateed Nch distributions
2107  //
2108  if(nch<=0) return 0.;
2109  if(!fHistoMeasNch || !fHistoMCNch) { AliError("Input histos to evaluate Nch weights missing"); return 0.; }
2110  Double_t pMeas=fHistoMeasNch->GetBinContent(fHistoMeasNch->FindBin(nch));
2111  Double_t pMC=fHistoMCNch->GetBinContent(fHistoMCNch->FindBin(nch));
2112  Double_t weight = pMC>0 ? pMeas/pMC : 0.;
2113  if(fUseMultRatioAsWeight) weight = pMC;
2114  return weight;
2115 }
2116 //__________________________________________________________________________________________________
2118  // creates historgam with measured multiplcity distribution in pp 7 TeV collisions (from Eur. Phys. J. C (2010) 68: 345–354)
2119  //
2120  // for Nch > 70 the points were obtained with a double NBD distribution fit
2121  // 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
2122  // fit1->SetParameter(1,1.63); // k parameter
2123  // fit1->SetParameter(2,12.8); // mean multiplicity
2124  //
2125  Double_t nchbins[82]={0.50,1.50,2.50,3.50,4.50,5.50,6.50,7.50,8.50,9.50,
2126  10.50,11.50,12.50,13.50,14.50,15.50,16.50,17.50,18.50,19.50,
2127  20.50,21.50,22.50,23.50,24.50,25.50,26.50,27.50,28.50,29.50,
2128  30.50,31.50,32.50,33.50,34.50,35.50,36.50,37.50,38.50,39.50,
2129  40.50,41.50,42.50,43.50,44.50,45.50,46.50,47.50,48.50,49.50,
2130  50.50,51.50,52.50,53.50,54.50,55.50,56.50,57.50,58.50,59.50,
2131  60.50,62.50,64.50,66.50,68.50,70.50,72.50,74.50,76.50,78.50,
2132  80.50,82.50,84.50,86.50,88.50,90.50,92.50,94.50,96.50,98.50,
2133  100.50,102.50};
2134  Double_t pch[81]={0.062011,0.072943,0.070771,0.067245,0.062834,0.057383,0.051499,0.04591,0.041109,0.036954,
2135  0.03359,0.030729,0.028539,0.026575,0.024653,0.0229,0.021325,0.019768,0.018561,0.017187,
2136  0.01604,0.014836,0.013726,0.012576,0.011481,0.010393,0.009502,0.008776,0.008024,0.007452,
2137  0.006851,0.006428,0.00594,0.005515,0.005102,0.00469,0.004162,0.003811,0.003389,0.003071,
2138  0.002708,0.002422,0.002184,0.001968,0.00186,0.00165,0.001577,0.001387,0.001254,0.001118,
2139  0.001037,0.000942,0.000823,0.000736,0.000654,0.000579,0.000512,0.00049,0.00045,0.000355,
2140  0.000296,0.000265,0.000193,0.00016,0.000126,0.0000851, 0.0000676,0.0000537,0.0000426, 0.0000338,
2141  0.0000268,0.0000213,0.0000166,0.0000133,0.0000106,0.00000837,0.00000662, 0.00000524,0.00000414, 0.00000327,
2142  0.00000258};
2143 
2144  if(fHistoMeasNch) delete fHistoMeasNch;
2145  fHistoMeasNch=new TH1F("hMeaseNch","",81,nchbins);
2146  for(Int_t i=0; i<81; i++){
2147  fHistoMeasNch->SetBinContent(i+1,pch[i]);
2148  fHistoMeasNch->SetBinError(i+1,0.);
2149  }
2150 }
2151 
2152 //__________________________________________________________________________________________________
2154  // processes output of Ds is selected
2155  Bool_t keep=kFALSE;
2156  if(recoAnalysisCuts > 0){
2157  Int_t isKKpi=recoAnalysisCuts&1;
2158  Int_t ispiKK=recoAnalysisCuts&2;
2159  Int_t isPhiKKpi=recoAnalysisCuts&4;
2160  Int_t isPhipiKK=recoAnalysisCuts&8;
2161  Int_t isK0starKKpi=recoAnalysisCuts&16;
2162  Int_t isK0starpiKK=recoAnalysisCuts&32;
2163  if(fDsOption==1){
2164  if(isKKpi && isPhiKKpi) keep=kTRUE;
2165  if(ispiKK && isPhipiKK) keep=kTRUE;
2166  }
2167  else if(fDsOption==2){
2168  if(isKKpi && isK0starKKpi) keep=kTRUE;
2169  if(ispiKK && isK0starpiKK) keep=kTRUE;
2170  }
2171  else if(fDsOption==3)keep=kTRUE;
2172  }
2173  return keep;
2174 }
2175 //__________________________________________________________________________________________________
2177 
2178  // processes output of Lc->V0+bnachelor is selected
2179 
2180  Bool_t keep=kFALSE;
2181 
2182  if (recoAnalysisCuts > 0){
2183 
2184  Int_t isK0Sp = recoAnalysisCuts&1;
2185  Int_t isLambdaBarpi = recoAnalysisCuts&2;
2186  Int_t isLambdapi = recoAnalysisCuts&4;
2187 
2188  if(fLctoV0bachelorOption == 1){
2189  if(isK0Sp) keep=kTRUE;
2190  }
2191  else if(fLctoV0bachelorOption == 2){
2192  if(isLambdaBarpi) keep=kTRUE;
2193  }
2194  else if(fLctoV0bachelorOption == 4){
2195  if(isLambdapi) keep=kTRUE;
2196  }
2197  else if(fLctoV0bachelorOption == 7) {
2198  if (isK0Sp || isLambdaBarpi || isLambdapi) keep=kTRUE;
2199  }
2200  }
2201  return keep;
2202 }
2203 
2204 //____________________________________________________________________________
2205 TProfile* AliCFTaskVertexingHF::GetEstimatorHistogram(const AliVEvent* event){
2206  // Get Estimator Histogram from period event->GetRunNumber();
2207  //
2208  // If you select SPD tracklets in |eta|<1 you should use type == 1
2209  //
2210 
2211  Int_t runNo = event->GetRunNumber();
2212  Int_t period = -1; // pp: 0-LHC10b, 1-LHC10c, 2-LHC10d, 3-LHC10e
2213  // pPb: 0-LHC13b, 1-LHC13c
2214 
2215  if (fIsPPbData) { // setting run numbers for LHC13 if pPb
2216  if (runNo>195343 && runNo<195484) period = 0;
2217  if (runNo>195528 && runNo<195678) period = 1;
2218  if (period<0 || period>1) return 0;
2219  } else { //else assume pp 2010
2220  if(runNo>114930 && runNo<117223) period = 0;
2221  if(runNo>119158 && runNo<120830) period = 1;
2222  if(runNo>122373 && runNo<126438) period = 2;
2223  if(runNo>127711 && runNo<130841) period = 3;
2224  if(period<0 || period>3) return 0;
2225  }
2226 
2227  return fMultEstimatorAvg[period];
2228 }
2229 
2230 //________________________________________________________________________
2233 
2234  Int_t nTracks=aod->GetNumberOfTracks();
2235  Double_t nHarmonic=2.;
2236  Double_t q2Vec[2] = {0.,0.};
2237  Int_t multQvec=0;
2238 
2239  for(Int_t it=0; it<nTracks; it++){
2240  AliAODTrack* track=(AliAODTrack*)aod->GetTrack(it);
2241  if(!track) continue;
2242  TString genname = AliVertexingHFUtils::GetGenerator(track->GetLabel(),mcHeader);
2243  if(genname.Contains("Hijing")) {
2244  if(track->TestFilterBit(BIT(8))||track->TestFilterBit(BIT(9))) {
2245  Double_t pt=track->Pt();
2246  Double_t eta=track->Eta();
2247  Double_t phi=track->Phi();
2248  Double_t qx=TMath::Cos(nHarmonic*phi);
2249  Double_t qy=TMath::Sin(nHarmonic*phi);
2250  if(eta<etamax && eta>etamin && pt>ptmin && pt<ptmax) {
2251  q2Vec[0]+=qx;
2252  q2Vec[1]+=qy;
2253  multQvec++;
2254  }
2255  }
2256  }
2257  }
2258 
2259  Double_t q2 = 0.;
2260  if(multQvec>0) q2 = TMath::Sqrt(q2Vec[0]*q2Vec[0]+q2Vec[1]*q2Vec[1])/TMath::Sqrt(multQvec);
2261 
2262  return q2;
2263 }
2264 
void SetTrackArray(TClonesArray *trkarray)
fast configuration, only a subset of variables
void SetCentralityValue(Float_t centValue)
virtual AliESDtrackCuts * GetTrackCutsV0daughters() const
Definition: AliRDHFCuts.h:265
void SetFillFromGenerated(Bool_t flag)
Double_t GetPtWeightFromHistogram(Float_t pt)
Bool_t fUseAdditionalCuts
flag for pPb data (used for multiplicity corrections)
Bool_t fUseZWeight
flag to decide to use a flat pt shape
AliCFTaskVertexingHF & operator=(const AliCFTaskVertexingHF &c)
Bool_t fFillFromGenerated
decay channel to configure the task
Double_t GetWeight(Float_t pt)
double Double_t
Definition: External.C:58
virtual AliESDtrackCuts * GetTrackCutsSoftPi() const
Definition: AliRDHFCuts.h:264
Int_t fCountAcc
MC particle found in limited acceptance that doesn&#39;t satisfy acceptance cuts.
void SetUseMCVertex()
Definition: AliRDHFCuts.h:376
Int_t fCountRefit
Reco particle found that satisfy vertex constrained.
Int_t MCcquarkCounting(AliAODMCParticle *mcPart) const
Class for HF corrections as a function of many variables and steps For D* and other cascades...
Int_t IsEventSelectedInCentrality(AliVEvent *event)
Int_t fCountReco
Reco particle found that satisfy kTPCrefit and kITSrefit.
Bool_t HasSelectionBit(Int_t i) const
UInt_t fPDGcode
flag to use selection bit
Class for HF corrections as a function of many variables and step.
Bool_t FillRecoContainer(Double_t *containerInput)
Double_t ComputeTPCq2(AliAODEvent *aod, AliAODMCHeader *mcHeader, Double_t etamin, Double_t etamax, Double_t ptmin, Double_t ptmax) const
Int_t fCountGenLimAcc
MC particle found.
Bool_t fIsPPbData
flag for pp data (not checking centrality)
void SetPtWeightsFromD0Cent080dataModOhoverLHC16i2abc()
TH1F * fHistoPtWeight
user-defined function to be used to calculate weights
Bool_t fUseSelectionBit
Lc->V0+bachelor decay option (generation level)
char Char_t
Definition: External.C:18
AliCFTaskVertexingHF()
multiplicity estimators
Int_t fCountRecoITSClusters
Reco particle found that satisfy cuts in requested acceptance.
static TString GetGenerator(Int_t label, AliAODMCHeader *header)
Int_t fNvar
flag to use directly the ratio of the distributions (fHistoMCNch) instead of computing it ...
void SetNVar(Int_t nVar)
void SetRejectCandidateIfNotFromQuark(Bool_t opt)
TCanvas * c
Definition: TestFitELoss.C:172
virtual Bool_t SetLabelArray()
Int_t fGenLctoV0bachelorOption
Lc->V0+bachelor decay option (selection level)
Double_t GetZWeight(Float_t z, Int_t runnumber)
void SetPtWeightsFromFONLL276andBAMPSoverLHC12a17b()
Bool_t FillRecoCand(AliVEvent *event, AliAODRecoDecayHF3Prong *rd3)
Int_t fCountGenLimAccNoAcc
MC particle found in limited acceptance.
Double_t dNdptFit(Float_t pt, Double_t *par)
Int_t GetNProngs() const
TString fPartName
number of variables for the container
Double_t GetMaxVtxZ() const
Definition: AliRDHFCuts.h:258
Bool_t FillRecoCasc(AliVEvent *event, AliAODRecoCascadeHF *rc, Bool_t isDStar, Bool_t recoSecVtx=kFALSE)
Double_t fRefMult
TProfile with mult vas. Z per period.
Bool_t ProcessDs(Int_t returnCodeDs) const
Float_t fCutOnMomConservation
Skip filling the unneed steps for most of the analyses to save disk space.
Bool_t fZvtxCorrectedNtrkEstimator
refrence multiplcity (period b)
TList * fListProfiles
response matrix for unfolding
const Double_t etamin
Class for cuts on AOD reconstructed D+->Kpipi.
Bool_t fAcceptanceUnf
flag to select D0 origins. 0 Only from charm 1 only from beauty 2 both from charm and beauty ...
static Int_t GetNumberOfTrackletsInEtaRange(AliAODEvent *ev, Double_t mineta, Double_t maxeta)
Bool_t fUseTrackletsWeight
flag to decide whether to use Ncharged weights != 1 when filling the container or not ...
int Int_t
Definition: External.C:63
Int_t fCountVertex
MC particle found that satisfy acceptance cuts.
Definition: External.C:204
void SetMinCentrality(Float_t minCentrality=0.)
Definition: AliRDHFCuts.h:52
Int_t GetIsFilled() const
float Float_t
Definition: External.C:68
Double_t DodzFit(Float_t z, Double_t *par)
Bool_t fUseCascadeTaskForLctoV0bachelor
these are the pre-selection cuts for the TMVA
const Double_t ptmax
Double_t GetNchWeight(Int_t nch)
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
Definition: External.C:212
AliAODVertex * GetOwnPrimaryVtx() const
Char_t fSign
daughter in fin state
AliESDtrackCuts * GetTrackCuts() const
Definition: AliRDHFCuts.h:263
slow configuration, all variables
const Double_t ptmin
Float_t GetCentrality(AliAODEvent *aodEvent)
Definition: AliRDHFCuts.h:259
AliRDHFCuts * fCuts
flag for unfolding before or after cuts.
Int_t fFakeSelection
flag to switch off the centrality selection
void SetFakeSelection(Int_t fakeSel)
TF1 * fFuncWeight
configuration (slow / fast) of the CF –> different variables will be allocated (all / reduced number...
UInt_t fResonantDecay
histogram with Nch distribution from MC production
Bool_t fCentralitySelection
flag to decide wheter to keep D0 only (0), D0bar only (1), or both D0 and D0bar (2) ...
Int_t fCountRecoAcc
Reco particle found that satisfy cuts.
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)
Int_t fGenDsOption
Ds decay option (selection level)
void SetMaxCentrality(Float_t maxCentrality=100.)
Definition: AliRDHFCuts.h:53
Int_t colors[nPtBins]
void SetOwnPrimaryVtx(const AliAODVertex *vtx)
void SetRecoPrimVertex(Double_t zPrimVertex)
Bool_t fIsPPData
flag to use the z-vtx corrected (if not use uncorrected) multiplicity estimator
void SetMCCandidateParam(Int_t label)
Bool_t CheckMCPartFamily(AliAODMCParticle *, TClonesArray *) const
virtual Bool_t SetRecoCandidateParam(AliAODRecoDecayHF *)
static Double_t GetCorrectedNtracklets(TProfile *estimatorAvg, Double_t uncorrectedNacc, Double_t vtxZ, Double_t refMult)
Int_t fConfiguration
Ds decay option (generation level)
void SetPtWeightsFromFONLL5andDplusdataoverLHC16i2a()
Int_t CheckReflexion(Char_t isSign)
Bool_t IsEventSelected(AliVEvent *event)
Int_t fLctoV0bachelorOption
resonant deacy channel to be used if the CF should be run on resonant channels only ...
Bool_t fRejectIfNoQuark
selection flag for fakes tracks
void SetUsePID(Bool_t flag=kTRUE)
Definition: AliRDHFCuts.h:210
void SetPtWeightsFromD0Cent080dataModMartinezoverLHC16i2abc()
TH1I * fHistEventsProcessed
pointer to the CF manager
void Setq2Value(Double_t q2)
virtual void PrintAll() const
Bool_t fFillMinimumSteps
flag to define which task to use for Lc –> K0S+p
TString fDauNames
D meson name.
Bool_t RecoAcceptStep(AliESDtrackCuts **trackCuts) const
Int_t fCountRecoPID
Reco particle found that satisfy cuts in PPR.
const Double_t etamax
Bool_t fUseMultRatioAsWeight
flag to decide whether to use Ncharged weights != 1 when filling the container or not ...
static Int_t GetGeneratedPhysicalPrimariesInEtaRange(TClonesArray *arrayMC, Double_t mineta, Double_t maxeta)
Bool_t IsSelected(TObject *obj)
Definition: AliRDHFCuts.h:293
Int_t fEvents
Reco PID step.
Bool_t MCRefitStep(AliAODEvent *aodEvent, AliESDtrackCuts **trackCuts) const
Bool_t fUseCutsForTMVA
flag to use additional cuts needed for Lc –> K0S + p, TMVA
Int_t fDecayChannel
n. of events
Bool_t GetIsUsePID() const
Definition: AliRDHFCuts.h:271
const char Option_t
Definition: External.C:48
void UserExec(Option_t *option)
Int_t fCountRecoPPR
Reco particle found that satisfy cuts in n. of ITS clusters.
Bool_t fUseMCVertex
flag to remove events not geenrated with PYTHIA
bool Bool_t
Definition: External.C:53
void SetTriggerClass(TString trclass0, TString trclass1="")
Definition: AliRDHFCuts.h:198
Int_t fDsOption
flag to use MC vertex (useful when runnign in pp)
Bool_t fUseNchWeight
flag to decide whether to use z-vtx weights != 1 when filling the container or not ...
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:314
TProfile * GetEstimatorHistogram(const AliVEvent *event)
void SetConfiguration(Int_t configuration)
void SetMultiplicity(Double_t multiplicity)
void SetMCPrimaryVertex(Double_t zMCVertex)
Int_t fMultiplicityEstimator
PDG code.
void UserCreateOutputObjects()
ANALYSIS FRAMEWORK STUFF to loop on data and fill output objects.
TH1F * fHistoMCNch
histogram with measured Nch distribution (pp 7 TeV)
TProfile * fMultEstimatorAvg[4]
Definition of the multiplicity estimator: kNtrk10=0, kNtrk10to16=1, kVZERO=2.
Double_t fWeight
flag to decide whether to use pt-weights != 1 when filling the container or not
TH1F * fHistoMeasNch
user-defined histogram to calculate the Pt weights
Bool_t ProcessLctoV0Bachelor(Int_t returnCodeDs) const
Class for HF corrections as a function of many variables and step.
TList * OpenFile(const char *fname)
Definition: DrawAnaELoss.C:65
super fast configuration, only (pt,y,centrality)