AliPhysics  88b7ad0 (88b7ad0)
AliRDHFCuts.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 /* $Id$ */
17 
19 //
20 // Base class for cuts on AOD reconstructed heavy-flavour decay
21 //
22 // Author: A.Dainese, andrea.dainese@pd.infn.it
24 #include <Riostream.h>
25 
26 #include "AliVEvent.h"
27 #include "AliESDEvent.h"
28 #include "AliAODEvent.h"
29 #include "AliVVertex.h"
30 #include "AliESDVertex.h"
31 #include "AliLog.h"
32 #include "AliAODVertex.h"
33 #include "AliESDtrack.h"
34 #include "AliAODTrack.h"
35 #include "AliESDtrackCuts.h"
36 #include "AliCentrality.h"
37 #include "AliAODRecoDecayHF.h"
38 #include "AliAnalysisVertexingHF.h"
39 #include "AliAODMCHeader.h"
40 #include "AliAODMCParticle.h"
41 #include "AliVertexerTracks.h"
42 #include "AliRDHFCuts.h"
43 #include "AliAnalysisManager.h"
44 #include "AliAODHandler.h"
45 #include "AliInputEventHandler.h"
46 #include "AliPIDResponse.h"
47 #include "AliAnalysisUtils.h"
48 #include "AliMultSelection.h"
49 #include "AliAODVZERO.h"
50 #include "TRandom.h"
51 #include <TF1.h>
52 #include <TFile.h>
53 #include <TKey.h>
54 
55 using std::cout;
56 using std::endl;
57 
59 ClassImp(AliRDHFCuts);
61 
62 
63 //--------------------------------------------------------------------------
65 AliAnalysisCuts(name,title),
66 fMinVtxType(3),
67 fMinVtxContr(1),
68 fMaxVtxRedChi2(1e6),
69 fMaxVtxZ(10.),
70 fMinSPDMultiplicity(0),
71 fMaxVtxChi2PileupMV(5.),
72 fMinWDzPileupMV(15.),
73 fRejectPlpFromDiffBCMV(kFALSE),
74 fTriggerMask(AliVEvent::kAnyINT),
75 fUseOnlyOneTrigger(kFALSE),
76 fTrackCuts(0),
77 fnPtBins(1),
78 fnPtBinLimits(1),
79 fPtBinLimits(0),
80 fnVars(1),
81 fVarNames(0),
82 fnVarsForOpt(0),
83 fVarsForOpt(0),
84 fGlobalIndex(1),
85 fCutsRD(0),
86 fIsUpperCut(0),
87 fUsePID(kFALSE),
88 fUseAOD049(kFALSE),
89 fPidHF(0),
90 fWhyRejection(0),
91 fEvRejectionBits(0),
92 fRemoveDaughtersFromPrimary(kFALSE),
93 fUseMCVertex(kFALSE),
94 fUsePhysicsSelection(kTRUE),
95 fOptPileup(0),
96 fMinContrPileup(3),
97 fMinDzPileup(0.6),
98 fUseMultDepPileupCut(kFALSE),
99 fUseCentrality(0),
100 fMinCentrality(0.),
101 fMaxCentrality(100.),
102 fMultSelectionObjectName("MultSelection"),
103 fFixRefs(kFALSE),
104 fIsSelectedCuts(0),
105 fIsSelectedPID(0),
106 fMinPtCand(-1.),
107 fMaxPtCand(100000.),
108 fMaxRapidityCand(-999.),
109 fKeepSignalMC(kFALSE),
110 fIsCandTrackSPDFirst(kFALSE),
111 fMaxPtCandTrackSPDFirst(0.),
112 fApplySPDDeadPbPb2011(kFALSE),
113 fApplySPDMisalignedPP2012(kFALSE),
114 fMaxDiffTRKV0Centr(-1.),
115 fRemoveTrackletOutliers(kFALSE),
116 fCutOnzVertexSPD(0),
117 fKinkReject(kFALSE),
118 fUseTrackSelectionWithFilterBits(kTRUE),
119 fUseCentrFlatteningInMC(kFALSE),
120 fHistCentrDistr(0x0),
121 fCutRatioClsOverCrossRowsTPC(0),
122 fCutRatioSignalNOverCrossRowsTPC(0),
123 fCutMinCrossedRowsTPCPtDep(""),
124 f1CutMinNCrossedRowsTPCPtDep(0x0),
125 fUseCutGeoNcrNcl(kFALSE),
126 fDeadZoneWidth(3.),
127 fCutGeoNcrNclLength(130.),
128 fCutGeoNcrNclGeom1Pt(1.5),
129 fCutGeoNcrNclFractionNcr(0.85),
130 fCutGeoNcrNclFractionNcl(0.7),
131 fUseV0ANDSelectionOffline(kFALSE),
132 fUseTPCtrackCutsOnThisDaughter(kTRUE)
133 {
134  //
135  // Default Constructor
136  //
137  fTriggerClass[0]="CINT1"; fTriggerClass[1]="";
138 }
139 //--------------------------------------------------------------------------
141  AliAnalysisCuts(source),
142  fMinVtxType(source.fMinVtxType),
143  fMinVtxContr(source.fMinVtxContr),
145  fMaxVtxZ(source.fMaxVtxZ),
150  fTriggerMask(source.fTriggerMask),
152  fTriggerClass(),
153  fTrackCuts(0),
154  fnPtBins(source.fnPtBins),
156  fPtBinLimits(0),
157  fnVars(source.fnVars),
158  fVarNames(0),
159  fnVarsForOpt(source.fnVarsForOpt),
160  fVarsForOpt(0),
161  fGlobalIndex(source.fGlobalIndex),
162  fCutsRD(0),
163  fIsUpperCut(0),
164  fUsePID(source.fUsePID),
165  fUseAOD049(source.fUseAOD049),
166  fPidHF(0),
170  fUseMCVertex(source.fUseMCVertex),
172  fOptPileup(source.fOptPileup),
174  fMinDzPileup(source.fMinDzPileup),
180  fFixRefs(source.fFixRefs),
183  fMinPtCand(source.fMinPtCand),
184  fMaxPtCand(source.fMaxPtCand),
194  fKinkReject(source.fKinkReject),
197  fHistCentrDistr(0x0),
210 {
211  //
212  // Copy constructor
213  //
214  cout<<"Copy constructor"<<endl;
215  fTriggerClass[0] = source.fTriggerClass[0];
216  fTriggerClass[1] = source.fTriggerClass[1];
217  if(source.GetTrackCuts()) AddTrackCuts(source.GetTrackCuts());
218  if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);
219  if(source.fVarNames) SetVarNames(source.fnVars,source.fVarNames,source.fIsUpperCut);
220  if(source.fCutsRD) SetCuts(source.fGlobalIndex,source.fCutsRD);
221  if(source.fVarsForOpt) SetVarsForOpt(source.fnVarsForOpt,source.fVarsForOpt);
222  if(source.fPidHF) SetPidHF(source.fPidHF);
223  if(source.fHistCentrDistr) fHistCentrDistr=(TH1F*)(source.fHistCentrDistr->Clone());
226  PrintAll();
227 
228 }
229 //--------------------------------------------------------------------------
231 {
232  //
233  // assignment operator
234  //
235  if(&source == this) return *this;
236 
237  AliAnalysisCuts::operator=(source);
238 
239  fMinVtxType=source.fMinVtxType;
240  fMinVtxContr=source.fMinVtxContr;
242  fMaxVtxZ=source.fMaxVtxZ;
247  fTriggerMask=source.fTriggerMask;
249  fTriggerClass[0]=source.fTriggerClass[0];
250  fTriggerClass[1]=source.fTriggerClass[1];
251  fnPtBins=source.fnPtBins;
253  fnVars=source.fnVars;
254  fGlobalIndex=source.fGlobalIndex;
255  fnVarsForOpt=source.fnVarsForOpt;
256  fUsePID=source.fUsePID;
257  fUseAOD049=source.fUseAOD049;
258  if(fPidHF) delete fPidHF;
259  fPidHF=new AliAODPidHF(*(source.GetPidHF()));
263  fUseMCVertex=source.fUseMCVertex;
265  fOptPileup=source.fOptPileup;
267  fMinDzPileup=source.fMinDzPileup;
273  fFixRefs=source.fFixRefs;
276  fMinPtCand=source.fMinPtCand;
277  fMaxPtCand=source.fMaxPtCand;
287  fKinkReject=source.fKinkReject;
289  if(fHistCentrDistr) delete fHistCentrDistr;
291  if(source.fHistCentrDistr)fHistCentrDistr=(TH1F*)(source.fHistCentrDistr->Clone());
292 
293  if(source.GetTrackCuts()) {delete fTrackCuts; fTrackCuts=new AliESDtrackCuts(*(source.GetTrackCuts()));}
294  if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);
295  if(source.fVarNames) SetVarNames(source.fnVars,source.fVarNames,source.fIsUpperCut);
296  if(source.fCutsRD) SetCuts(source.fGlobalIndex,source.fCutsRD);
297  if(source.fVarsForOpt) SetVarsForOpt(source.fnVarsForOpt,source.fVarsForOpt);
298 
304 
313 
314  PrintAll();
315 
316  return *this;
317 }
318 //--------------------------------------------------------------------------
320  //
321  // Default Destructor
322  //
323  if(fTrackCuts) { delete fTrackCuts; fTrackCuts=0; }
324  if(fPtBinLimits) {delete [] fPtBinLimits; fPtBinLimits=0;}
325  if(fVarNames) {delete [] fVarNames; fVarNames=0;}
326  if(fVarsForOpt) {delete [] fVarsForOpt; fVarsForOpt=0;}
327  if(fCutsRD) {
328  delete [] fCutsRD;
329  fCutsRD=0;
330  }
331  if(fIsUpperCut) {delete [] fIsUpperCut; fIsUpperCut=0;}
332  if(fPidHF){
333  delete fPidHF;
334  fPidHF=0;
335  }
337 
341  }
342 
343 }
344 //---------------------------------------------------------------------------
346  //
347  // Centrality selection
348  //
349  if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid){
350  AliWarning("Centrality estimator not valid");
351  return 3;
352  }else{
353  Float_t centvalue=GetCentrality((AliAODEvent*)event);
354  if (centvalue<-998.){//-999 if no centralityP
355  return 3;
357  return 3;
358  }else{
359  if (centvalue<fMinCentrality || centvalue>fMaxCentrality){
360  return 2;
361  }
362  // selection to flatten the centrality distribution
363  if(fHistCentrDistr){
364  if(!IsEventSelectedForCentrFlattening(centvalue))return 4;
365  }
366  }
367  }
368  return 0;
369 }
370 
371 
372 //-------------------------------------------------
373 void AliRDHFCuts::SetHistoForCentralityFlattening(TH1F *h,Double_t minCentr,Double_t maxCentr,Double_t centrRef,Int_t switchTRand){
374  // set the histo for centrality flattening
375  // the centrality is flatten in the range minCentr,maxCentr
376  // if centrRef is zero, the minimum in h within (minCentr,maxCentr) defines the reference
377  // positive, the value of h(centrRef) defines the reference (-> the centrality distribution might be not flat in the whole desired range)
378  // negative, h(bin with max in range)*centrRef is used to define the reference (-> defines the maximum loss of events, also in this case the distribution might be not flat)
379  // switchTRand is used to set the unerflow bin of the histo: if it is < -1 in the analysis the random event selection will be done on using TRandom
380 
381  if(maxCentr<minCentr){
382  AliWarning("AliRDHFCuts::Wrong centralities values while setting the histogram for centrality flattening");
383  }
384 
386  fHistCentrDistr=(TH1F*)h->Clone("hCentralityFlat");
387  fHistCentrDistr->SetTitle("Reference histo for centrality flattening");
388  Int_t minbin=fHistCentrDistr->FindBin(minCentr*1.00001); // fast if fix bin width
389  Int_t maxbin=fHistCentrDistr->FindBin(maxCentr*0.9999);
390  fHistCentrDistr->GetXaxis()->SetRange(minbin,maxbin);
391  Double_t ref=0.,bincont=0.,binrefwidth=1.;
392  Int_t binref=0;
393  if(TMath::Abs(centrRef)<0.0001){
394  binref=fHistCentrDistr->GetMinimumBin();
395  binrefwidth=fHistCentrDistr->GetBinWidth(binref);
396  ref=fHistCentrDistr->GetBinContent(binref)/binrefwidth;
397  }
398  else if(centrRef>0.){
399  binref=h->FindBin(centrRef);
400  if(binref<1||binref>h->GetNbinsX()){
401  AliWarning("AliRDHFCuts::Wrong centrality reference value while setting the histogram for centrality flattening");
402  }
403  binrefwidth=fHistCentrDistr->GetBinWidth(binref);
404  ref=fHistCentrDistr->GetBinContent(binref)/binrefwidth;
405  }
406  else{
407  if(centrRef<-1) AliWarning("AliRDHFCuts: with this centrality reference no flattening will be applied");
408  binref=fHistCentrDistr->GetMaximumBin();
409  binrefwidth=fHistCentrDistr->GetBinWidth(binref);
410  ref=fHistCentrDistr->GetMaximum()*TMath::Abs(centrRef)/binrefwidth;
411  }
412 
413  for(Int_t j=1;j<=h->GetNbinsX();j++){// Now set the "probabilities"
414  if(h->GetBinLowEdge(j)*1.0001>=minCentr&&h->GetBinLowEdge(j+1)*0.9999<=maxCentr){
415  bincont=h->GetBinContent(j);
416  fHistCentrDistr->SetBinContent(j,ref/bincont*h->GetBinWidth(j));
417  fHistCentrDistr->SetBinError(j,h->GetBinError(j)*ref/bincont);
418  }
419  else{
420  h->SetBinContent(j,1.1);// prob > 1 to assure that events will not be rejected
421  }
422  }
423 
424  fHistCentrDistr->SetBinContent(0,switchTRand);
425  return;
426 
427 }
428 
429 //-------------------------------------------------
431  //
432  // Random event selection, based on fHistCentrDistr, to flatten the centrality distribution
433  // Can be faster if it was required that fHistCentrDistr covers
434  // exactly the desired centrality range (e.g. part of the lines below should be done during the
435  // setting of the histo) and TH1::SetMinimum called
436  //
437 
438  if(!fHistCentrDistr) return kTRUE;
439  // Int_t maxbin=fHistCentrDistr->FindBin(fMaxCentrality*0.9999);
440  // if(maxbin>fHistCentrDistr->GetNbinsX()){
441  // AliWarning("AliRDHFCuts: The maximum centrality exceeds the x-axis limit of the histogram for centrality flattening");
442  // }
443 
444  Int_t bin=fHistCentrDistr->FindBin(centvalue); // Fast if the histo has a fix bin
445  Double_t bincont=fHistCentrDistr->GetBinContent(bin);
446  Double_t centDigits=centvalue-(Int_t)(centvalue*100.)/100.;// this is to extract a random number between 0 and 0.01
447 
448  if(fHistCentrDistr->GetBinContent(0)<-0.9999){
449  if(gRandom->Uniform(1.)<bincont)return kTRUE;
450  return kFALSE;
451  }
452 
453  if(centDigits*100.<bincont)return kTRUE;
454  return kFALSE;
455 
456 }
457 //---------------------------------------------------------------------------
458 void AliRDHFCuts::SetupPID(AliVEvent *event) {
459  // Set the PID response object in the AliAODPidHF
460  // in case of old PID sets the TPC dE/dx BB parameterization
461 
462  if(fPidHF){
463  if(fPidHF->GetPidResponse()==0x0){
464  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
465  AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
466  AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
467  fPidHF->SetPidResponse(pidResp);
468  }
470  if(fPidHF->GetOldPid()) {
471 
472  Bool_t isMC=kFALSE;
473  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)event)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
474  if(mcArray) {isMC=kTRUE;fUseAOD049=kFALSE;}
475 
476  // pp, from LHC10d onwards
477  if((event->GetRunNumber()>121693 && event->GetRunNumber()<136851) ||
478  event->GetRunNumber()>139517) fPidHF->SetOnePad(kTRUE);
479  // pp, 2011 low energy run
480  if((event->GetRunNumber()>=146686 && event->GetRunNumber()<=146860)){
481  fPidHF->SetppLowEn2011(kTRUE);
482  fPidHF->SetOnePad(kFALSE);
483  }
484  // PbPb LHC10h
485  if(event->GetRunNumber()>=136851 && event->GetRunNumber()<=139517) fPidHF->SetPbPb(kTRUE);
486  // MC
487  if(isMC) fPidHF->SetMC(kTRUE);
488  if(isMC && (event->GetRunNumber()>=146686 && event->GetRunNumber()<=146860))
489  fPidHF->SetMClowenpp2011(kTRUE);
491  }else{
492  // check that AliPIDResponse object was properly set in case of using OADB
493  if(fPidHF->GetPidResponse()==0x0) AliFatal("AliPIDResponse object not set");
494  }
495  }
496 }
497 //---------------------------------------------------------------------------
499  //
500  // Event selection
501  //
502  //if(fTriggerMask && event->GetTriggerMask()!=fTriggerMask) return kFALSE;
503 
504 
505 
506  fWhyRejection=0;
508  Bool_t accept=kTRUE;
509 
510 
511  // check if it's MC
512  Bool_t isMC=kFALSE;
513  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)event)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
514  if(mcArray) {isMC=kTRUE;fUseAOD049=kFALSE;}
515 
516 
517  SetupPID(event);
518 
519  // trigger class
520  TString firedTriggerClasses=((AliAODEvent*)event)->GetFiredTriggerClasses();
521  // don't do for MC and for PbPb 2010 data
522  if(!isMC && (event->GetRunNumber()<136851 || event->GetRunNumber()>139517)) {
523  if(!firedTriggerClasses.Contains(fTriggerClass[0].Data()) &&
524  (fTriggerClass[1].CompareTo("")==0 || !firedTriggerClasses.Contains(fTriggerClass[1].Data())) ) {
525  fWhyRejection=5;
527  accept=kFALSE;
528  }
529  }
530 
531  // TEMPORARY FIX FOR GetEvent
532  Int_t nTracks=((AliAODEvent*)event)->GetNumberOfTracks();
533  for(Int_t itr=0; itr<nTracks; itr++){
534  AliAODTrack* tr=(AliAODTrack*)((AliAODEvent*)event)->GetTrack(itr);
535  tr->SetAODEvent((AliAODEvent*)event);
536  }
537 
538  // TEMPORARY FIX FOR REFERENCES
539  // Fix references to daughter tracks
540  // if(fFixRefs) {
541  // AliAnalysisVertexingHF *fixer = new AliAnalysisVertexingHF();
542  // fixer->FixReferences((AliAODEvent*)event);
543  // delete fixer;
544  // }
545  //
546 
547 
548  // physics selection requirements
550  Bool_t isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & fTriggerMask);
551  if(!isSelected) {
552  if(accept) fWhyRejection=7;
554  accept=kFALSE;
555  }else{
556  if(fUseOnlyOneTrigger){
557  if(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected()!=fTriggerMask){
558  if(accept) fWhyRejection=7;
560  accept=kFALSE;
561  }
562  }
564  AliAODVZERO* v0data=(AliAODVZERO*)((AliAODEvent*)event)->GetVZEROData();
565  Int_t tv0a=v0data->GetV0ADecision();
566  Int_t tv0c=v0data->GetV0CDecision();
567  if(!(tv0a==1 && tv0c==1)){
568  if(accept) fWhyRejection=7;
570  accept=kFALSE;
571  }
572  }
573  }
574  }
575 
576  // centrality selection
577  if (fUseCentrality!=kCentOff) {
578  Int_t rejection=IsEventSelectedInCentrality(event);
579  Bool_t okCent=kFALSE;
580  if(rejection==0) okCent=kTRUE;
581  if(isMC && rejection==4 && !fUseCentrFlatteningInMC) okCent=kTRUE;
582  if(!okCent){
583  if(accept) fWhyRejection=rejection;
586  accept=kFALSE;
587  }
588 
589  }
590 
591  // PbPb2011 outliers in tracklets vs. VZERO and centTRK vs. centV0
592  if(event->GetRunNumber()>=167693 && event->GetRunNumber()<=170593){
594  Double_t v0cent=GetCentrality((AliAODEvent*)event,kCentV0M);
595  Double_t ntracklets=((AliAODEvent*)event)->GetTracklets()->GetNumberOfTracklets();
596  Double_t cutval=60.-0.08*ntracklets+1./50000.*ntracklets*ntracklets;
597  if(ntracklets<1000. && v0cent<cutval){
598  if(accept) fWhyRejection=2;
600  accept=kFALSE;
601  }
602  }
603  if(fMaxDiffTRKV0Centr>0.){
604  Double_t v0cent=GetCentrality((AliAODEvent*)event,kCentV0M);
605  Double_t trkcent=GetCentrality((AliAODEvent*)event,kCentTRK);
606  if(TMath::Abs(trkcent-v0cent)>fMaxDiffTRKV0Centr){
607  if(accept) fWhyRejection=1;
609  accept=kFALSE;
610  }
611  }
612  }
613 
614  // vertex requirements
615 
616  const AliVVertex *vertex = event->GetPrimaryVertex();
617 
618  if(!vertex){
619  accept=kFALSE;
621  }else{
622  TString title=vertex->GetTitle();
623  if(title.Contains("Z") && fMinVtxType>1){
624  accept=kFALSE;
626  }
627  else if(title.Contains("3D") && fMinVtxType>2){
628  accept=kFALSE;
630  }
631  if(vertex->GetNContributors()<fMinVtxContr){
632  accept=kFALSE;
634  }
635  if(TMath::Abs(vertex->GetZ())>fMaxVtxZ) {
637  if(accept) fWhyRejection=6;
638  accept=kFALSE;
639  }
640  }
641 
642  if(fCutOnzVertexSPD>0){
643  const AliVVertex *vSPD = ((AliAODEvent*)event)->GetPrimaryVertexSPD();
644  if(!vSPD || (vSPD && vSPD->GetNContributors()<fMinVtxContr)){
645  accept=kFALSE;
647  }else{
648  if(fCutOnzVertexSPD==1 && TMath::Abs(vSPD->GetZ())>12.) {
649  // protection for events with bad reconstructed track vertex (introduced for 2011 Pb-Pb)
651  if(accept) fWhyRejection=6;
652  accept=kFALSE;
653  }
654  if(fCutOnzVertexSPD>=2 && vertex){
655  Double_t dz = vSPD->GetZ()-vertex->GetZ();
656  // cut on absolute distance between track and SPD vertex (introduced for 2011 Pb-Pb)
657  if(TMath::Abs(dz)>0.5) {
659  if(accept) fWhyRejection=0;
660  accept=kFALSE;
661  }
662  if(accept && fCutOnzVertexSPD==3){
663  // cut on nsigma distance between track and SPD vertex (for 2015 Pb-Pb)
664  double covTrc[6],covSPD[6];
665  vertex->GetCovarianceMatrix(covTrc);
666  vSPD->GetCovarianceMatrix(covSPD);
667  double errTot = TMath::Sqrt(covTrc[5]+covSPD[5]);
668  double errTrc = TMath::Sqrt(covTrc[5]);
669  double nsigTot = TMath::Abs(dz)/errTot, nsigTrc = TMath::Abs(dz)/errTrc;
670  if (TMath::Abs(dz)>0.2 || nsigTot>10 || nsigTrc>20){
671  // reject, bad reconstructed track vertex
673  if(accept) fWhyRejection=0;
674  accept=kFALSE;
675  }
676  }
677  }
678  }
679  }
680 
681  // pile-up rejection
683  Bool_t isPileup=kFALSE;
685  isPileup=event->IsPileupFromSPDInMultBins();
686  }else{
689  isPileup=event->IsPileupFromSPD(cutc,cutz,3.,2.,10.);
690  }
691  if(isPileup){
692  if(accept) fWhyRejection=1;
694  accept=kFALSE;
695  }
696  }
698  AliAnalysisUtils utils;
699  utils.SetMinPlpContribMV(fMinContrPileup); // min. multiplicity of the pile-up vertex to consider
700  utils.SetMaxPlpChi2MV(fMaxVtxChi2PileupMV); // max chi2 per contributor of the pile-up vertex to consider.
701  utils.SetMinWDistMV(fMinWDzPileupMV); // minimum weighted distance in Z between 2 vertices (i.e. (zv1-zv2)/sqrt(sigZv1^2+sigZv2^2) )
702  utils.SetCheckPlpFromDifferentBCMV(fRejectPlpFromDiffBCMV); // vertex with |BCID|>2 will trigger pile-up (based on TOF)
703  Bool_t isPUMV = utils.IsPileUpMV(event);
704  if(isPUMV) {
705  if(accept) fWhyRejection=1;
707  accept=kFALSE;
708  }
709  }
710 
711 
712  // Correcting PP2012 flag to remoce tracks crossing SPD misaligned staves for periods 12def
713  if(fApplySPDMisalignedPP2012 && !(event->GetRunNumber()>=195681 && event->GetRunNumber()<=197388)) fApplySPDMisalignedPP2012=false;
714 
715  return accept;
716 }
717 //---------------------------------------------------------------------------
719  //
720  // Daughter tracks selection
721  //
722  if(!fTrackCuts) return kTRUE;
723 
724  Int_t ndaughters = d->GetNDaughters();
725  AliAODVertex *vAOD = d->GetPrimaryVtx();
726  Double_t pos[3],cov[6];
727  vAOD->GetXYZ(pos);
728  vAOD->GetCovarianceMatrix(cov);
729  const AliESDVertex vESD(pos,cov,100.,100);
730 
731  Bool_t retval=kTRUE;
732 
733  for(Int_t idg=0; idg<ndaughters; idg++) {
734  AliAODTrack *dgTrack = (AliAODTrack*)d->GetDaughter(idg);
735  if(!dgTrack) {retval = kFALSE; continue;}
736  //printf("charge %d\n",dgTrack->Charge());
737  if(dgTrack->Charge()==0) continue; // it's not a track, but a V0
738 
740  { if(!dgTrack->HasPointOnITSLayer(0)) { retval = kFALSE; continue; } }
741 
742  if(!IsDaughterSelected(dgTrack,&vESD,fTrackCuts,aod)) retval = kFALSE;
743  }
744 
745  return retval;
746 }
747 //---------------------------------------------------------------------------
749  //
750  // Check if AOD and deltaAOD files are composed of the same events:
751  // mismatches were observed in the merged AODs of LHC15o
752  //
753  // When AOD+deltaAOD are produced from ESD, mismatches can be found looking at:
754  // - the AOD trees in AliAOD.root and AliAOD.VertexingHF.root have different number of entries
755  // - the titles of the TProcessID objects do not match
756  // When deltaAOD are produced from AOD, mismatches can be found looking at:
757  // - the AOD trees in AliAOD.root and AliAOD.VertexingHF.root have different number of entries
758  //
759  // Return values:
760  // -1: AOD and deltaAOD trees have different number of entries
761  // 0: AOD and deltaAOD trees have same number of entries + the titles of the TProcessID objects do not match
762  // 1: AOD and deltaAOD trees have same number of entries + the titles of the TProcessID objects match
763  Bool_t okTProcessNames = kTRUE;
764  AliAODHandler* aodHandler = (AliAODHandler*)((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
765  TTree *treeAOD = aodHandler->GetTree();
766  TTree *treeDeltaAOD = treeAOD->GetFriend("aodTree");
767  if(!treeDeltaAOD || !treeAOD) return -1;
768  if(treeDeltaAOD && treeAOD){
769  if(treeAOD->GetEntries()!=treeDeltaAOD->GetEntries()){
770  printf("AliRDHFCuts::CheckMatchingAODdeltaAODevents: Difference in number of entries in main and friend tree, skipping event\n");
771  return -1;
772  }
773  TFile *mfile = treeAOD->GetCurrentFile();
774  TFile *dfile = treeDeltaAOD->GetCurrentFile();
775  TList* lm=mfile->GetListOfKeys();
776  TList* ld=dfile->GetListOfKeys();
777  Int_t nentm=lm->GetEntries();
778  for(Int_t jm=0; jm<nentm; jm++){
779  TKey* o=(TKey*)lm->At(jm);
780  TString clnam=o->GetClassName();
781  if(clnam=="TProcessID"){
782  TString pname=o->GetName();
783  TString ptit=o->GetTitle();
784  if(pname.Contains("ProcessID")){
785  TObject* od=(TObject*)ld->FindObject(pname.Data());
786  if(od){
787  TString ptit2=od->GetTitle();
788  if(ptit2!=ptit){
789  printf("AliRDHFCuts::CheckMatchingAODdeltaAODevents: mismatch in %s: AOD: %s -- deltaAOD: %s\n",pname.Data(),ptit.Data(),ptit2.Data());
790  okTProcessNames = kFALSE;
791  }
792  }
793  }
794  }
795  }
796  }
797 
798  if (okTProcessNames) return 1;
799  else return 0;
800 }
801 //---------------------------------------------------------------------------
803  //
804  // Check the correctness of the string syntax
805  //
806  Bool_t retval=kTRUE;
807 
808  if(!rows.Contains("pt")) {
809  if(print) AliError("string must contain \"pt\"");
810  retval= kFALSE;
811  }
812  return retval;
813 }
814 //---------------------------------------------------------------------------
816  //
817  //Create the TFormula from TString for TPC crossed rows pT dependent cut
818  //
819 
820 
821  // setting data member that describes the TPC crossed rows pT dependent cut
823 
824  // creating TFormula from TString
827  // resetting TFormula
829  }
830  if(!CheckPtDepCrossedRows(rows,kTRUE))return;
831 
832  TString tmp(rows);
833  tmp.ReplaceAll("pt","x");
834  f1CutMinNCrossedRowsTPCPtDep = new TFormula("f1CutMinNCrossedRowsTPCPtDep",tmp.Data());
835 
836 
837 }
838 //---------------------------------------------------------------------------
839 Bool_t AliRDHFCuts::IsDaughterSelected(AliAODTrack *track,const AliESDVertex *primary,AliESDtrackCuts *cuts, const AliAODEvent* aod) const{
840  //
841  // Convert to ESDtrack, relate to vertex and check cuts
842  //
843  if(!cuts) return kTRUE;
844 
845  if(cuts->GetFlagCutTOFdistance()) cuts->SetFlagCutTOFdistance(kFALSE);
846 
847 
848  // convert to ESD track here
849  AliESDtrack esdTrack(track);
850  // set the TPC cluster info
851  esdTrack.SetTPCClusterMap(track->GetTPCClusterMap());
852  esdTrack.SetTPCSharedMap(track->GetTPCSharedMap());
853  esdTrack.SetTPCPointsF(track->GetTPCNclsF());
854  // needed to calculate the impact parameters
855  esdTrack.RelateToVertex(primary,0.,3.);
856 
857  //applying ESDtrackCut
858  if(!cuts->IsSelected(&esdTrack)) return kFALSE;
859 
860  //appliyng kink rejection
861  if(fKinkReject){
862  AliAODVertex *maybeKink=track->GetProdVertex();
863  if(maybeKink->GetType()==AliAODVertex::kKink) return kFALSE;
864  }
865 
866  //appliyng TPC crossed rows pT dependent cut
868  Float_t nCrossedRowsTPC = esdTrack.GetTPCCrossedRows();
869  if(nCrossedRowsTPC<f1CutMinNCrossedRowsTPCPtDep->Eval(esdTrack.Pt())) return kFALSE;
870  }
871 
872  //appliyng NTPCcls/NTPCcrossedRows cut
874  Float_t nCrossedRowsTPC = esdTrack.GetTPCCrossedRows();
875  Float_t nClustersTPC = esdTrack.GetTPCNcls();
876  if(nCrossedRowsTPC!=0){
877  Float_t ratio = nClustersTPC/nCrossedRowsTPC;
878  if(ratio<fCutRatioClsOverCrossRowsTPC) return kFALSE;
879  }
880  else return kFALSE;
881  }
882 
883  //appliyng TPCsignalN/NTPCcrossedRows cut
885  Float_t nCrossedRowsTPC = esdTrack.GetTPCCrossedRows();
886  Float_t nTPCsignal = esdTrack.GetTPCsignalN();
887  if(nCrossedRowsTPC!=0){
888  Float_t ratio = nTPCsignal/nCrossedRowsTPC;
889  if(ratio<fCutRatioSignalNOverCrossRowsTPC) return kFALSE;
890  }
891  else return kFALSE;
892  }
893 
894  // geometrical cut (note uses track at vertex instead of at TPC inner wall)
896  Float_t nCrossedRowsTPC = esdTrack.GetTPCCrossedRows();
897  Float_t lengthInActiveZoneTPC=esdTrack.GetLengthInActiveZone(0,fDeadZoneWidth,220.,aod->GetMagneticField());
898  Double_t cutGeoNcrNclLength=fCutGeoNcrNclLength-TMath::Power(TMath::Abs(esdTrack.GetSigned1Pt()),fCutGeoNcrNclGeom1Pt);
899  Bool_t isOK=kTRUE;
900  if (lengthInActiveZoneTPC<cutGeoNcrNclLength) isOK=kFALSE;
901  if (nCrossedRowsTPC<fCutGeoNcrNclFractionNcr*cutGeoNcrNclLength) isOK=kFALSE;
902  if (esdTrack.GetTPCncls()<fCutGeoNcrNclFractionNcl*cutGeoNcrNclLength) isOK=kFALSE;
903  if(!isOK) return kFALSE;
904  }
905 
906 
908  // to be implemented
909  // we need either to have here the AOD Event,
910  // or to have the pileup vertex object
911  }
912 
914  Bool_t deadSPDLay1PbPb2011[20][4]={
915  {kTRUE,kTRUE,kTRUE,kTRUE},
916  {kTRUE,kTRUE,kTRUE,kTRUE},
917  {kTRUE,kTRUE,kTRUE,kTRUE},
918  {kTRUE,kTRUE,kTRUE,kTRUE},
919  {kTRUE,kTRUE,kTRUE,kTRUE},
920  {kFALSE,kFALSE,kTRUE,kTRUE},
921  {kTRUE,kTRUE,kFALSE,kFALSE},
922  {kTRUE,kTRUE,kTRUE,kTRUE},
923  {kFALSE,kFALSE,kTRUE,kTRUE},
924  {kTRUE,kTRUE,kTRUE,kTRUE},
925  {kTRUE,kTRUE,kFALSE,kFALSE},
926  {kTRUE,kTRUE,kTRUE,kTRUE},
927  {kFALSE,kFALSE,kFALSE,kFALSE},
928  {kFALSE,kFALSE,kTRUE,kTRUE},
929  {kFALSE,kFALSE,kFALSE,kFALSE},
930  {kFALSE,kFALSE,kFALSE,kFALSE},
931  {kTRUE,kTRUE,kTRUE,kTRUE},
932  {kTRUE,kTRUE,kFALSE,kFALSE},
933  {kFALSE,kFALSE,kFALSE,kFALSE},
934  {kFALSE,kFALSE,kFALSE,kFALSE}
935  };
936  Bool_t deadSPDLay2PbPb2011[40][4]={
937  {kTRUE,kTRUE,kTRUE,kTRUE},
938  {kTRUE,kTRUE,kTRUE,kTRUE},
939  {kTRUE,kTRUE,kTRUE,kTRUE},
940  {kTRUE,kTRUE,kTRUE,kTRUE},
941  {kTRUE,kTRUE,kTRUE,kTRUE},
942  {kTRUE,kTRUE,kTRUE,kTRUE},
943  {kTRUE,kTRUE,kTRUE,kTRUE},
944  {kTRUE,kTRUE,kTRUE,kTRUE},
945  {kTRUE,kTRUE,kTRUE,kTRUE},
946  {kTRUE,kTRUE,kTRUE,kTRUE},
947  {kTRUE,kTRUE,kTRUE,kTRUE},
948  {kTRUE,kTRUE,kTRUE,kTRUE},
949  {kFALSE,kFALSE,kFALSE,kFALSE},
950  {kFALSE,kFALSE,kTRUE,kTRUE},
951  {kTRUE,kTRUE,kTRUE,kTRUE},
952  {kTRUE,kTRUE,kTRUE,kTRUE},
953  {kTRUE,kTRUE,kFALSE,kFALSE},
954  {kTRUE,kTRUE,kTRUE,kTRUE},
955  {kTRUE,kTRUE,kTRUE,kTRUE},
956  {kTRUE,kTRUE,kTRUE,kTRUE},
957  {kFALSE,kFALSE,kFALSE,kFALSE},
958  {kFALSE,kFALSE,kFALSE,kFALSE},
959  {kTRUE,kTRUE,kTRUE,kTRUE},
960  {kTRUE,kTRUE,kTRUE,kTRUE},
961  {kFALSE,kFALSE,kFALSE,kFALSE},
962  {kFALSE,kFALSE,kFALSE,kFALSE},
963  {kTRUE,kTRUE,kTRUE,kTRUE},
964  {kTRUE,kTRUE,kTRUE,kTRUE},
965  {kFALSE,kFALSE,kFALSE,kFALSE},
966  {kFALSE,kFALSE,kFALSE,kFALSE},
967  {kFALSE,kFALSE,kFALSE,kFALSE},
968  {kFALSE,kFALSE,kFALSE,kFALSE},
969  {kTRUE,kTRUE,kTRUE,kTRUE},
970  {kTRUE,kTRUE,kTRUE,kTRUE},
971  {kTRUE,kTRUE,kFALSE,kFALSE},
972  {kTRUE,kTRUE,kTRUE,kTRUE},
973  {kFALSE,kFALSE,kFALSE,kFALSE},
974  {kFALSE,kFALSE,kFALSE,kFALSE},
975  {kFALSE,kFALSE,kFALSE,kFALSE},
976  {kFALSE,kFALSE,kFALSE,kFALSE}
977  };
978  Double_t xyz1[3],xyz2[3];
979  esdTrack.GetXYZAt(3.9,0.,xyz1);
980  esdTrack.GetXYZAt(7.6,0.,xyz2);
981  Double_t phi1=TMath::ATan2(xyz1[1],xyz1[0]);
982  if(phi1<0) phi1+=2*TMath::Pi();
983  Int_t lad1=(Int_t)(phi1/(2.*TMath::Pi()/20.));
984  Double_t phi2=TMath::ATan2(xyz2[1],xyz2[0]);
985  if(phi2<0) phi2+=2*TMath::Pi();
986  Int_t lad2=(Int_t)(phi2/(2.*TMath::Pi()/40.));
987  Int_t mod1=(Int_t)((xyz1[2]+14)/7.);
988  Int_t mod2=(Int_t)((xyz2[2]+14)/7.);
989  Bool_t lay1ok=kFALSE;
990  if(mod1>=0 && mod1<4 && lad1<20){
991  lay1ok=deadSPDLay1PbPb2011[lad1][mod1];
992  }
993  Bool_t lay2ok=kFALSE;
994  if(mod2>=0 && mod2<4 && lad2<40){
995  lay2ok=deadSPDLay2PbPb2011[lad2][mod2];
996  }
997  if(!lay1ok && !lay2ok) return kFALSE;
998  }
999 
1001  // Cut tracks crossing the SPD at 5.6<phi<2pi
1002  Double_t xyz1[3],xyz2[3];
1003  esdTrack.GetXYZAt(3.9,0.,xyz1);
1004  esdTrack.GetXYZAt(7.6,0.,xyz2);
1005  Double_t phi1=TMath::ATan2(xyz1[1],xyz1[0]);
1006  if(phi1<0) phi1+=2*TMath::Pi();
1007  Double_t phi2=TMath::ATan2(xyz2[1],xyz2[0]);
1008  if(phi2<0) phi2+=2*TMath::Pi();
1009  Bool_t lay1ok=kTRUE;
1010  if(phi1>5.6 && phi1<2.*TMath::Pi()) lay1ok=kFALSE;
1011  Bool_t lay2ok=kTRUE;
1012  if(phi2>5.6 && phi2<2.*TMath::Pi()) lay2ok=kFALSE;
1013  if(!lay1ok || !lay2ok) return kFALSE;
1014  }
1015 
1016  return kTRUE;
1017 }
1018 //---------------------------------------------------------------------------
1019 void AliRDHFCuts::SetPtBins(Int_t nPtBinLimits,Float_t *ptBinLimits) {
1020  // Set the pt bins
1021 
1022  if(fPtBinLimits) {
1023  delete [] fPtBinLimits;
1024  fPtBinLimits = NULL;
1025  printf("Changing the pt bins\n");
1026  }
1027 
1028  if(nPtBinLimits != fnPtBins+1){
1029  cout<<"Warning: ptBinLimits dimention "<<nPtBinLimits<<" != nPtBins+1 ("<<fnPtBins+1<<")\nSetting nPtBins to "<<nPtBinLimits-1<<endl;
1030  SetNPtBins(nPtBinLimits-1);
1031  }
1032 
1033  fnPtBinLimits = nPtBinLimits;
1034  SetGlobalIndex();
1035  //cout<<"Changing also Global Index -> "<<fGlobalIndex<<endl;
1037  for(Int_t ib=0; ib<nPtBinLimits; ib++) fPtBinLimits[ib]=ptBinLimits[ib];
1038 
1039  return;
1040 }
1041 //---------------------------------------------------------------------------
1042 void AliRDHFCuts::SetVarNames(Int_t nVars,TString *varNames,Bool_t *isUpperCut){
1043  // Set the variable names
1044 
1045  if(fVarNames) {
1046  delete [] fVarNames;
1047  fVarNames = NULL;
1048  //printf("Changing the variable names\n");
1049  }
1050  if(nVars!=fnVars){
1051  printf("Wrong number of variables: it has to be %d\n",fnVars);
1052  return;
1053  }
1054  //fnVars=nVars;
1055  fVarNames = new TString[nVars];
1056  fIsUpperCut = new Bool_t[nVars];
1057  for(Int_t iv=0; iv<nVars; iv++) {
1058  fVarNames[iv] = varNames[iv];
1059  fIsUpperCut[iv] = isUpperCut[iv];
1060  }
1061 
1062  return;
1063 }
1064 //---------------------------------------------------------------------------
1066  // Set the variables to be used for cuts optimization
1067 
1068  if(fVarsForOpt) {
1069  delete [] fVarsForOpt;
1070  fVarsForOpt = NULL;
1071  //printf("Changing the variables for cut optimization\n");
1072  }
1073 
1074  if(nVars==0){
1075  printf("%d not accepted as number of variables: it has to be %d\n",nVars,fnVars);
1076  return;
1077  }
1078 
1079  fnVarsForOpt = 0;
1080  fVarsForOpt = new Bool_t[fnVars];
1081  for(Int_t iv=0; iv<fnVars; iv++) {
1082  fVarsForOpt[iv]=forOpt[iv];
1083  if(fVarsForOpt[iv]) fnVarsForOpt++;
1084  }
1085 
1086  return;
1087 }
1088 
1089 //---------------------------------------------------------------------------
1091  //
1092  // set centrality estimator
1093  //
1094  fUseCentrality=flag;
1095  if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid) AliWarning("Centrality estimator not valid");
1096 
1097  return;
1098 }
1099 
1100 
1101 //---------------------------------------------------------------------------
1103  //
1104  // store the cuts
1105  //
1106  if(nVars!=fnVars) {
1107  printf("Wrong number of variables: it has to be %d\n",fnVars);
1108  AliFatal("exiting");
1109  }
1110  if(nPtBins!=fnPtBins) {
1111  printf("Wrong number of pt bins: it has to be %d\n",fnPtBins);
1112  AliFatal("exiting");
1113  }
1114 
1115  if(!fCutsRD) fCutsRD = new Float_t[fGlobalIndex];
1116 
1117 
1118  for(Int_t iv=0; iv<fnVars; iv++) {
1119 
1120  for(Int_t ib=0; ib<fnPtBins; ib++) {
1121 
1122  //check
1123  if(GetGlobalIndex(iv,ib)>=fGlobalIndex) {
1124  cout<<"Overflow, exit..."<<endl;
1125  return;
1126  }
1127 
1128  fCutsRD[GetGlobalIndex(iv,ib)] = cutsRD[iv][ib];
1129 
1130  }
1131  }
1132  return;
1133 }
1134 //---------------------------------------------------------------------------
1135 void AliRDHFCuts::SetCuts(Int_t glIndex,Float_t* cutsRDGlob){
1136  //
1137  // store the cuts
1138  //
1139  if(glIndex != fGlobalIndex){
1140  cout<<"Wrong array size: it has to be "<<fGlobalIndex<<endl;
1141  AliFatal("exiting");
1142  }
1143  if(!fCutsRD) fCutsRD = new Float_t[fGlobalIndex];
1144 
1145  for(Int_t iGl=0;iGl<fGlobalIndex;iGl++){
1146  fCutsRD[iGl] = cutsRDGlob[iGl];
1147  }
1148  return;
1149 }
1150 //---------------------------------------------------------------------------
1152  //
1153  // print all cuts values
1154  //
1155 
1156  printf("Minimum vtx type %d\n",fMinVtxType);
1157  printf("Minimum vtx contr %d\n",fMinVtxContr);
1158  printf("Max vtx red chi2 %f\n",fMaxVtxRedChi2);
1159  printf("Min SPD mult %d\n",fMinSPDMultiplicity);
1160  printf("Use PID %d OldPid=%d\n",(Int_t)fUsePID,fPidHF ? fPidHF->GetOldPid() : -1);
1161  printf("Remove daughters from vtx %d\n",(Int_t)fRemoveDaughtersFromPrimary);
1162  printf("Physics selection: %s\n",fUsePhysicsSelection ? "Yes" : "No");
1163  printf("Pileup rejection: %s\n",(fOptPileup > 0) ? "Yes" : "No");
1164  if(fOptPileup==1) printf(" -- Reject pileup event");
1165  if(fOptPileup==2) printf(" -- Reject tracks from pileup vtx");
1166  if(fUseCentrality>0) {
1167  TString estimator="";
1168  if(fUseCentrality==kCentV0M) estimator = "V0";
1169  if(fUseCentrality==kCentTRK) estimator = "Tracks";
1170  if(fUseCentrality==kCentTKL) estimator = "Tracklets";
1171  if(fUseCentrality==kCentCL1) estimator = "SPD clusters outer";
1172  if(fUseCentrality==kCentZNA) estimator = "ZNA";
1173  if(fUseCentrality==kCentZPA) estimator = "ZPA";
1174  if(fUseCentrality==kCentV0A) estimator = "V0A";
1175  if(fUseCentrality==kCentCL0) estimator = "SPD clusters inner";
1176  printf("Centrality class considered: %.1f-%.1f, estimated with %s\n",fMinCentrality,fMaxCentrality,estimator.Data());
1177  }
1178  if(fIsCandTrackSPDFirst) printf("Check for candidates with pt < %2.2f, that daughters fullfill kFirst criteria\n",fMaxPtCandTrackSPDFirst);
1179 
1180  if(fCutRatioClsOverCrossRowsTPC) printf("N TPC Clusters > %f N TPC Crossed Rows\n", fCutRatioClsOverCrossRowsTPC);
1181  if(fCutRatioSignalNOverCrossRowsTPC) printf("N TPC Points for dE/dx > %f N TPC Crossed Rows\n", fCutRatioSignalNOverCrossRowsTPC);
1182  if(f1CutMinNCrossedRowsTPCPtDep) printf("N TPC Crossed Rows pT-dependent cut: %s\n", fCutMinCrossedRowsTPCPtDep.Data());
1183 
1184  if(fVarNames){
1185  cout<<"Array of variables"<<endl;
1186  for(Int_t iv=0;iv<fnVars;iv++){
1187  cout<<fVarNames[iv]<<"\t";
1188  }
1189  cout<<endl;
1190  }
1191  if(fVarsForOpt){
1192  cout<<"Array of optimization"<<endl;
1193  for(Int_t iv=0;iv<fnVars;iv++){
1194  cout<<fVarsForOpt[iv]<<"\t";
1195  }
1196  cout<<endl;
1197  }
1198  if(fIsUpperCut){
1199  cout<<"Array of upper/lower cut"<<endl;
1200  for(Int_t iv=0;iv<fnVars;iv++){
1201  cout<<fIsUpperCut[iv]<<"\t";
1202  }
1203  cout<<endl;
1204  }
1205  if(fPtBinLimits){
1206  cout<<"Array of ptbin limits"<<endl;
1207  for(Int_t ib=0;ib<fnPtBinLimits;ib++){
1208  cout<<fPtBinLimits[ib]<<"\t";
1209  }
1210  cout<<endl;
1211  }
1212  if(fCutsRD){
1213  cout<<"Matrix of cuts"<<endl;
1214  for(Int_t iv=0;iv<fnVars;iv++){
1215  for(Int_t ib=0;ib<fnPtBins;ib++){
1216  cout<<"fCutsRD["<<iv<<"]["<<ib<<"] = "<<fCutsRD[GetGlobalIndex(iv,ib)]<<"\t";
1217  }
1218  cout<<endl;
1219  }
1220  cout<<endl;
1221  }
1222  if(fPidHF) fPidHF->PrintAll();
1223  return;
1224 }
1225 
1226 //--------------------------------------------------------------------------
1228  // print the trigger selection
1229 
1230  printf("Selected trigger classes: %s %s\n",fTriggerClass[0].Data(),fTriggerClass[1].Data());
1231 
1232  cout<<" Trigger selection pattern: ";
1233  if( fTriggerMask & AliVEvent::kAny ) cout<<" kAny ";
1234  if( fTriggerMask & AliVEvent::kAnyINT ) cout<<" kAnyINT ";
1235  if( fTriggerMask & AliVEvent::kINT7 ) cout<<" kINT7 ";
1236  if( fTriggerMask & AliVEvent::kMB ) cout<<" kMB ";
1237  if( fTriggerMask & AliVEvent::kCINT5 ) cout<<" kCINT5 ";
1238  if( fTriggerMask & AliVEvent::kCentral ) cout<<" kCentral ";
1239  if( fTriggerMask & AliVEvent::kSemiCentral ) cout<<" kSemiCentral ";
1240  if( fTriggerMask & AliVEvent::kEMCEGA ) cout<<" kEMCEGA ";
1241  if( fTriggerMask & AliVEvent::kHighMult ) cout<<" kHighMult ";
1242  if( fTriggerMask & AliVEvent::kFastOnly ) cout<<" kFastOnly ";
1243  cout << endl<< endl;
1244 
1245 }
1246 
1247 //---------------------------------------------------------------------------
1248 void AliRDHFCuts::GetCuts(Float_t**& cutsRD) const{
1249  //
1250  // get the cuts
1251  //
1252 
1253  //cout<<"Give back a "<<fnVars<<"x"<<fnPtBins<<" matrix."<<endl;
1254 
1255 
1256  Int_t iv,ib;
1257  if(!cutsRD) {
1258  //cout<<"Initialization..."<<endl;
1259  cutsRD=new Float_t*[fnVars];
1260  for(iv=0; iv<fnVars; iv++) {
1261  cutsRD[iv] = new Float_t[fnPtBins];
1262  }
1263  }
1264 
1265  for(Int_t iGlobal=0; iGlobal<fGlobalIndex; iGlobal++) {
1266  GetVarPtIndex(iGlobal,iv,ib);
1267  cutsRD[iv][ib] = fCutsRD[iGlobal];
1268  }
1269 
1270  return;
1271 }
1272 
1273 //---------------------------------------------------------------------------
1275  //
1276  // give the global index from variable and pt bin
1277  //
1278  return iPtBin*fnVars+iVar;
1279 }
1280 
1281 //---------------------------------------------------------------------------
1282 void AliRDHFCuts::GetVarPtIndex(Int_t iGlob, Int_t& iVar, Int_t& iPtBin) const {
1283  //
1284  //give the index of the variable and of the pt bin from the global index
1285  //
1286  iPtBin=(Int_t)iGlob/fnVars;
1287  iVar=iGlob%fnVars;
1288 
1289  return;
1290 }
1291 
1292 //---------------------------------------------------------------------------
1294  //
1295  //give the pt bin where the pt lies.
1296  //
1297  Int_t ptbin=-1;
1298  if(pt<fPtBinLimits[0])return ptbin;
1299  for (Int_t i=0;i<fnPtBins;i++){
1300  if(pt<fPtBinLimits[i+1]) {
1301  ptbin=i;
1302  break;
1303  }
1304  }
1305  return ptbin;
1306 }
1307 //-------------------------------------------------------------------
1309  //
1310  // Give the value of cut set for the variable iVar and the pt bin iPtBin
1311  //
1312  if(!fCutsRD){
1313  cout<<"Cuts not iniziaisez yet"<<endl;
1314  return 0;
1315  }
1316  return fCutsRD[GetGlobalIndex(iVar,iPtBin)];
1317 }
1318 
1319 //-------------------------------------------------------------------
1321 
1322  if(aodEvent->GetRunNumber()<244824)return GetCentralityOldFramework(aodEvent,estimator);
1323  Double_t cent=-999;
1324 
1325  if(estimator==kCentOff) return -999;
1326 
1327  AliMultSelection *multSelection = (AliMultSelection*)aodEvent->FindListObject(fMultSelectionObjectName);
1328  if(!multSelection){
1329  AliWarning("AliMultSelection could not be found in the aod event list of objects");
1330  return cent;
1331  }
1332 
1333  // Compare centrality with the centrality at AOD filtering and on-the-fly
1334  if( fMultSelectionObjectName.CompareTo("MultSelection")!=0 ){
1335  AliMultSelection *defaultmultSelection = (AliMultSelection*)aodEvent->FindListObject("MultSelection");
1336  if(!defaultmultSelection){
1337  AliWarning("AliMultSelection default method could not be found in the aod event list of objects");
1338  return cent;
1339  }
1340  Float_t defaultCent = defaultmultSelection->GetMultiplicityPercentile("V0M");
1341  Float_t newCent = multSelection->GetMultiplicityPercentile("V0M");
1342  if( defaultCent<20. && newCent>20.) fEvRejectionBits+=1<<kMismatchOldNewCentrality;
1343  else if (defaultCent>20. && newCent<20.) fEvRejectionBits+=1<<kMismatchOldNewCentrality;
1344  }
1345 
1346  if(estimator==kCentV0M){
1347  cent=multSelection->GetMultiplicityPercentile("V0M");
1348  }else if(estimator==kCentV0A){
1349  cent=multSelection->GetMultiplicityPercentile("V0A");
1350  }else if(estimator==kCentZNA){
1351  cent=multSelection->GetMultiplicityPercentile("ZNA");
1352  }else if(estimator==kCentCL1){
1353  cent=multSelection->GetMultiplicityPercentile("CL1");
1354  }else if(estimator==kCentCL0){
1355  cent=multSelection->GetMultiplicityPercentile("CL0");
1356  }else {
1357  AliWarning(Form("CENTRALITY ESTIMATE WITH ESTIMATOR %d NOT YET IMPLEMENTED FOR NEW FRAMEWORK",(Int_t)estimator));
1358  return cent;
1359  }
1360  Int_t qual = multSelection->GetEvSelCode();
1361  if(qual == 199 ) cent=-999;
1362  return cent;
1363 }
1364 //-------------------------------------------------------------------
1366  //
1367  // Get centrality percentile
1368  //
1369 
1370  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)aodEvent)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
1371  if(mcArray) {fUseAOD049=kFALSE;}
1372 
1373  AliAODHeader *header=dynamic_cast<AliAODHeader*>(aodEvent->GetHeader());
1374  if(!header) AliFatal("Not a standard AOD");
1375  AliCentrality *centrality=header->GetCentralityP();
1376  Float_t cent=-999.;
1377  Bool_t isSelRun=kFALSE;
1378  Int_t selRun[5]={138364, 138826, 138828, 138836, 138871};
1379  if(!centrality) return cent;
1380  else{
1381  if (estimator==kCentV0M){
1382  cent=(Float_t)(centrality->GetCentralityPercentile("V0M"));
1383  if(cent<0){
1384  Int_t quality = centrality->GetQuality();
1385  if(quality<=1){ // fQuality==1 means rejected by zVertex cut that we apply a part and we want to keep separate (Giacomo)
1386  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");
1387  }else{
1388  Int_t runnum=aodEvent->GetRunNumber();
1389  for(Int_t ir=0;ir<5;ir++){
1390  if(runnum==selRun[ir]){
1391  isSelRun=kTRUE;
1392  break;
1393  }
1394  }
1395  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");
1396  }
1397  }
1398 
1399  //temporary fix for AOD049 outliers
1400  if(fUseAOD049&&cent>=0){
1401  Float_t v0=0;
1402  AliAODVZERO* aodV0 = aodEvent->GetVZEROData();
1403  v0+=aodV0->GetMTotV0A();
1404  v0+=aodV0->GetMTotV0C();
1405  if(cent==0&&v0<19500)return -1;//filtering issue
1406  Float_t tkl = (Float_t)(aodEvent->GetTracklets()->GetNumberOfTracklets());
1407  Float_t val= 1.30552 + 0.147931 * v0;
1408  Float_t tklSigma[101]={176.644, 156.401, 153.789, 153.015, 142.476, 137.951, 136.127, 129.852, 127.436, 124.86, 120.788, 115.611, 113.172, 110.496, 109.127, 104.421, 102.479, 99.9766, 97.5152, 94.0654, 92.4602, 89.3364, 87.1342, 83.3497, 82.6216, 81.1084, 78.0793, 76.1234, 72.9434, 72.1334, 68.0056, 68.2755, 66.0376, 62.9666, 62.4274, 59.65, 58.3776, 56.6361, 54.5184, 53.4224, 51.932, 50.8922, 48.2848, 47.912, 46.5717, 43.4114, 43.2083, 41.3065, 40.1863, 38.5255, 37.2851, 37.5396, 34.4949, 33.8366, 31.8043, 31.7412, 30.8392, 30.0274, 28.8793, 27.6398, 26.6488, 25.0183, 25.1489, 24.4185, 22.9107, 21.2002, 21.6977, 20.1242, 20.4963, 19.0235, 19.298, 17.4103, 16.868, 15.2939, 15.2939, 16.0295, 14.186, 14.186, 15.2173, 12.9504, 12.9504, 12.9504, 15.264, 12.3674, 12.3674, 12.3674, 12.3674, 12.3674, 18.3811, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544, 13.7544};
1409  if ( TMath::Abs(tkl-val) > 6.*tklSigma[(Int_t)cent] )return -1;//outlier
1410  }
1411  }
1412  else {
1413  if (estimator==kCentTRK) {
1414  cent=(Float_t)(centrality->GetCentralityPercentile("TRK"));
1415  if(cent<0){
1416  Int_t quality = centrality->GetQuality();
1417  if(quality<=1){
1418  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TRK");
1419  }else{
1420  Int_t runnum=aodEvent->GetRunNumber();
1421  for(Int_t ir=0;ir<5;ir++){
1422  if(runnum==selRun[ir]){
1423  isSelRun=kTRUE;
1424  break;
1425  }
1426  }
1427  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TRK");
1428  }
1429  }
1430  }
1431  else{
1432  if (estimator==kCentTKL){
1433  cent=(Float_t)(centrality->GetCentralityPercentile("TKL"));
1434  if(cent<0){
1435  Int_t quality = centrality->GetQuality();
1436  if(quality<=1){
1437  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TKL");
1438  }else{
1439  Int_t runnum=aodEvent->GetRunNumber();
1440  for(Int_t ir=0;ir<5;ir++){
1441  if(runnum==selRun[ir]){
1442  isSelRun=kTRUE;
1443  break;
1444  }
1445  }
1446  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TKL");
1447  }
1448  }
1449  }
1450  else{
1451  if (estimator==kCentCL1){
1452  cent=(Float_t)(centrality->GetCentralityPercentile("CL1"));
1453  if(cent<0){
1454  Int_t quality = centrality->GetQuality();
1455  if(quality<=1){
1456  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("CL1");
1457  }else{
1458  Int_t runnum=aodEvent->GetRunNumber();
1459  for(Int_t ir=0;ir<5;ir++){
1460  if(runnum==selRun[ir]){
1461  isSelRun=kTRUE;
1462  break;
1463  }
1464  }
1465  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("CL1");
1466  }
1467  }
1468  }
1469  else{
1470  if (estimator==kCentZNA){
1471  cent=(Float_t)(centrality->GetCentralityPercentile("ZNA"));
1472  if(cent<0){
1473  Int_t quality = centrality->GetQuality();
1474  if(quality<=1){
1475  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("ZNA");
1476  }else{
1477  Int_t runnum=aodEvent->GetRunNumber();
1478  for(Int_t ir=0;ir<5;ir++){
1479  if(runnum==selRun[ir]){
1480  isSelRun=kTRUE;
1481  break;
1482  }
1483  }
1484  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("ZNA");
1485  }
1486  }
1487  }
1488  else{
1489  if (estimator==kCentZPA){
1490  cent=(Float_t)(centrality->GetCentralityPercentile("ZPA"));
1491  if(cent<0){
1492  Int_t quality = centrality->GetQuality();
1493  if(quality<=1){
1494  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("ZPA");
1495  }else{
1496  Int_t runnum=aodEvent->GetRunNumber();
1497  for(Int_t ir=0;ir<5;ir++){
1498  if(runnum==selRun[ir]){
1499  isSelRun=kTRUE;
1500  break;
1501  }
1502  }
1503  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("ZPA");
1504  }
1505  }
1506  }
1507  else{
1508  if (estimator==kCentV0A){
1509  cent=(Float_t)(centrality->GetCentralityPercentile("V0A"));
1510  if(cent<0){
1511  Int_t quality = centrality->GetQuality();
1512  if(quality<=1){
1513  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0A");
1514  }else{
1515  Int_t runnum=aodEvent->GetRunNumber();
1516  for(Int_t ir=0;ir<5;ir++){
1517  if(runnum==selRun[ir]){
1518  isSelRun=kTRUE;
1519  break;
1520  }
1521  }
1522  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0A");
1523  }
1524  }
1525  }
1526  else {
1527  AliWarning("Centrality estimator not valid");
1528 
1529  }
1530  }
1531  }
1532  }
1533  }
1534  }
1535  }
1536  }
1537  return cent;
1538 }
1539 //-------------------------------------------------------------------
1541  //
1542  // Compare two cuts objects
1543  //
1544 
1545  Bool_t areEqual=kTRUE;
1546 
1547  if(fMinVtxType!=obj->fMinVtxType) { printf("Minimum vtx type %d %d\n",fMinVtxType,obj->fMinVtxType); areEqual=kFALSE;}
1548 
1549  if(fMinVtxContr!=obj->fMinVtxContr) { printf("Minimum vtx contr %d %d\n",fMinVtxContr,obj->fMinVtxContr); areEqual=kFALSE;}
1550 
1551  if(TMath::Abs(fMaxVtxRedChi2-obj->fMaxVtxRedChi2)>1.e-10) { printf("Max vtx red chi2 %f %f\n",fMaxVtxRedChi2,obj->fMaxVtxRedChi2);areEqual=kFALSE;}
1552 
1553  if(fMinSPDMultiplicity!=obj->fMinSPDMultiplicity) { printf("Min SPD mult %d\n %d",fMinSPDMultiplicity,obj->fMinSPDMultiplicity);areEqual=kFALSE;}
1554 
1555  if(fUsePID!=obj->fUsePID) { printf("Use PID %d %d\n",(Int_t)fUsePID,(Int_t)obj->fUsePID); areEqual=kFALSE;}
1556 
1557  if(fRemoveDaughtersFromPrimary!=obj->fRemoveDaughtersFromPrimary) {printf("Remove daughters from vtx %d %d\n",(Int_t)fRemoveDaughtersFromPrimary,(Int_t)obj->fRemoveDaughtersFromPrimary); areEqual=kFALSE;}
1558  if(fTrackCuts){
1559  if(fTrackCuts->GetMinNClusterTPC()!=obj->fTrackCuts->GetMinNClusterTPC()) {printf("MinNClsTPC %d %d\n",fTrackCuts->GetMinNClusterTPC(),obj->fTrackCuts->GetMinNClusterTPC()); areEqual=kFALSE;}
1560 
1561  if(fTrackCuts->GetMinNClustersITS()!=obj->fTrackCuts->GetMinNClustersITS()) {printf("MinNClsITS %d %d\n",fTrackCuts->GetMinNClustersITS(),obj->fTrackCuts->GetMinNClustersITS()); areEqual=kFALSE;}
1562 
1563  if(TMath::Abs(fTrackCuts->GetMaxChi2PerClusterTPC()-obj->fTrackCuts->GetMaxChi2PerClusterTPC())>1.e-10) {printf("MaxChi2ClsTPC %f %f\n",fTrackCuts->GetMaxChi2PerClusterTPC(),obj->fTrackCuts->GetMaxChi2PerClusterTPC()); areEqual=kFALSE;}
1564 
1565  if(fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)!=obj->fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)) {printf("ClusterReq SPD %d %d\n",fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD),obj->fTrackCuts->GetClusterRequirementITS(AliESDtrackCuts::kSPD)); areEqual=kFALSE;}
1566  }
1567 
1568  if(fCutsRD) {
1569  for(Int_t iv=0;iv<fnVars;iv++) {
1570  for(Int_t ib=0;ib<fnPtBins;ib++) {
1571  if(TMath::Abs(fCutsRD[GetGlobalIndex(iv,ib)]-obj->fCutsRD[GetGlobalIndex(iv,ib)])>1.e-10) {
1572  cout<<"fCutsRD["<<iv<<"]["<<ib<<"] = "<<fCutsRD[GetGlobalIndex(iv,ib)]<<" "<<obj->fCutsRD[GetGlobalIndex(iv,ib)]<<"\n";
1573  areEqual=kFALSE;
1574  }
1575  }
1576  }
1577  }
1578 
1579  return areEqual;
1580 }
1581 //---------------------------------------------------------------------------
1583  //
1584  // print cuts values in table format
1585  //
1586 
1587  TString ptString = "pT range";
1588  if(fVarNames && fPtBinLimits && fCutsRD){
1589  TString firstLine(Form("* %-15s",ptString.Data()));
1590  for (Int_t ivar=0; ivar<fnVars; ivar++){
1591  firstLine+=Form("* %-15s ",fVarNames[ivar].Data());
1592  if (ivar == fnVars){
1593  firstLine+="*\n";
1594  }
1595  }
1596  Printf("%s",firstLine.Data());
1597 
1598  for (Int_t ipt=0; ipt<fnPtBins; ipt++){
1599  TString line;
1600  if (ipt==fnPtBins-1){
1601  line=Form("* %5.1f < pt < inf ",fPtBinLimits[ipt]);
1602  }
1603  else{
1604  line=Form("* %5.1f < pt < %4.1f ",fPtBinLimits[ipt],fPtBinLimits[ipt+1]);
1605  }
1606  for (Int_t ivar=0; ivar<fnVars; ivar++){
1607  line+=Form("* %-15f ",fCutsRD[GetGlobalIndex(ivar,ipt)]);
1608  }
1609  Printf("%s",line.Data());
1610  }
1611 
1612  }
1613 
1614 
1615  return;
1616 }
1617 //--------------------------------------------------------------------------
1619  AliAODEvent *aod) const
1620 {
1621  //
1622  // Recalculate primary vertex without daughters
1623  //
1624 
1625  if(!aod) {
1626  AliError("Can not remove daughters from vertex without AOD event");
1627  return 0;
1628  }
1629 
1630  AliAODVertex *recvtx=d->RemoveDaughtersFromPrimaryVtx(aod);
1631  if(!recvtx){
1632  AliDebug(2,"Removal of daughter tracks failed");
1633  return kFALSE;
1634  }
1635 
1636 
1637  //set recalculed primary vertex
1638  d->SetOwnPrimaryVtx(recvtx);
1639  delete recvtx;
1640 
1641  return kTRUE;
1642 }
1643 //--------------------------------------------------------------------------
1645 {
1646  //
1647  // Recalculate primary vertex without daughters
1648  //
1649 
1650  if(!aod) {
1651  AliError("Can not get MC vertex without AOD event");
1652  return kFALSE;
1653  }
1654 
1655  // load MC header
1656  AliAODMCHeader *mcHeader =
1657  (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
1658  if(!mcHeader) {
1659  AliError("Can not get MC vertex without AODMCHeader event");
1660  return kFALSE;
1661  }
1662  Double_t pos[3];
1663  Double_t covmatrix[6]={0.,0.,0.,0.,0.,0.};
1664  mcHeader->GetVertex(pos);
1665  AliAODVertex *recvtx=new AliAODVertex(pos,covmatrix);
1666 
1667  if(!recvtx){
1668  AliDebug(2,"Removal of daughter tracks failed");
1669  return kFALSE;
1670  }
1671 
1672  //set recalculed primary vertex
1673  d->SetOwnPrimaryVtx(recvtx);
1674 
1675  d->RecalculateImpPars(recvtx,aod);
1676 
1677  delete recvtx;
1678 
1679  return kTRUE;
1680 }
1681 //--------------------------------------------------------------------------
1683  AliAODEvent *aod,
1684  AliAODVertex *origownvtx) const
1685 {
1686  //
1687  // Clean-up own primary vertex if needed
1688  //
1689 
1691  d->UnsetOwnPrimaryVtx();
1692  if(origownvtx) {
1693  d->SetOwnPrimaryVtx(origownvtx);
1694  delete origownvtx; origownvtx=NULL;
1695  }
1696  d->RecalculateImpPars(d->GetPrimaryVtx(),aod);
1697  } else {
1698  if(origownvtx) {
1699  delete origownvtx; origownvtx=NULL;
1700  }
1701  }
1702  return;
1703 }
1704 //--------------------------------------------------------------------------
1705 Bool_t AliRDHFCuts::IsSignalMC(AliAODRecoDecay *d,AliAODEvent *aod,Int_t pdg) const
1706 {
1707  //
1708  // Checks if this candidate is matched to MC signal
1709  //
1710 
1711  if(!aod) return kFALSE;
1712 
1713  // get the MC array
1714  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)aod)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
1715 
1716  if(!mcArray) return kFALSE;
1717 
1718  // try to match
1719  Int_t label = d->MatchToMC(pdg,mcArray);
1720 
1721  if(label>=0) {
1722  //printf("MATCH!\n");
1723  return kTRUE;
1724  }
1725 
1726  return kFALSE;
1727 }
1728 
1729 
1730 //--------------------------------------------------------------------------
1732  // recompute event primary vertex from AOD tracks
1733 
1734  AliVertexerTracks *vertexer = new AliVertexerTracks(event->GetMagneticField());
1735  vertexer->SetITSMode();
1736  vertexer->SetMinClusters(3);
1737 
1738  AliAODVertex* pvtx=event->GetPrimaryVertex();
1739  if(strstr(pvtx->GetTitle(),"VertexerTracksWithConstraint")) {
1740  Float_t diamondcovxy[3];
1741  event->GetDiamondCovXY(diamondcovxy);
1742  Double_t pos[3]={event->GetDiamondX(),event->GetDiamondY(),0.};
1743  Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.*10.};
1744  AliESDVertex *diamond = new AliESDVertex(pos,cov,1.,1);
1745  vertexer->SetVtxStart(diamond);
1746  delete diamond; diamond=NULL;
1747  }
1748 
1749  AliESDVertex* vertexESD = (AliESDVertex*)vertexer->FindPrimaryVertex(event);
1750  if(!vertexESD) return kFALSE;
1751  if(vertexESD->GetNContributors()<=0) {
1752  //AliDebug(2,"vertexing failed");
1753  delete vertexESD; vertexESD=NULL;
1754  return kFALSE;
1755  }
1756  delete vertexer; vertexer=NULL;
1757 
1758  // convert to AliAODVertex
1759  Double_t pos[3],cov[6],chi2perNDF;
1760  vertexESD->GetXYZ(pos); // position
1761  vertexESD->GetCovMatrix(cov); //covariance matrix
1762  chi2perNDF = vertexESD->GetChi2toNDF();
1763  delete vertexESD; vertexESD=NULL;
1764 
1765  pvtx->SetPosition(pos[0],pos[1],pos[2]);
1766  pvtx->SetChi2perNDF(chi2perNDF);
1767  pvtx->SetCovMatrix(cov);
1768 
1769  return kTRUE;
1770 }
void SetNPtBins(Int_t nptBins)
Definition: AliRDHFCuts.h:393
Int_t fIsSelectedCuts
fix the daughter track references
Definition: AliRDHFCuts.h:445
Int_t pdg
Bool_t IsEventSelectedForCentrFlattening(Float_t centvalue)
Float_t fMinDzPileup
min. n. of tracklets in pileup vertex
Definition: AliRDHFCuts.h:434
Double_t fMaxPtCandTrackSPDFirst
flag to select the track kFirst criteria for pt < ptlimit
Definition: AliRDHFCuts.h:452
Int_t fWhyRejection
PID for heavy flavours manager.
Definition: AliRDHFCuts.h:427
double Double_t
Definition: External.C:58
Int_t fCutOnzVertexSPD
flag to apply cut on tracklets vs. centrality for 2011 data
Definition: AliRDHFCuts.h:457
void SetupPID(AliVEvent *event)
Bool_t IsSignalMC(AliAODRecoDecay *d, AliAODEvent *aod, Int_t pdg) const
const char * title
Definition: MakeQAPdf.C:27
Double_t fCutGeoNcrNclFractionNcr
3rd parameter of GeoNcrNcl cut
Definition: AliRDHFCuts.h:470
Bool_t fUseMCVertex
flag to switch on the removal of duaghters from the primary vertex computation
Definition: AliRDHFCuts.h:430
Bool_t fUseOnlyOneTrigger
trigger mask
Definition: AliRDHFCuts.h:409
Int_t IsEventSelectedInCentrality(AliVEvent *event)
Bool_t * fIsUpperCut
Definition: AliRDHFCuts.h:423
Bool_t SetMCPrimaryVtx(AliAODRecoDecayHF *d, AliAODEvent *aod) const
Bool_t fRemoveDaughtersFromPrimary
Definition: AliRDHFCuts.h:429
void SetHistoForCentralityFlattening(TH1F *h, Double_t minCentr, Double_t maxCentr, Double_t centrRef=0., Int_t switchTRand=0)
Bool_t GetUseCombined()
Definition: AliAODPidHF.h:165
Bool_t fUseMultDepPileupCut
min deltaz between main and pileup vertices
Definition: AliRDHFCuts.h:435
void SetUseCentrality(Int_t flag=1)
Float_t * fPtBinLimits
"number of limits", that is fnPtBins+1
Definition: AliRDHFCuts.h:416
Int_t fMinVtxContr
0: not cut; 1: SPDZ; 2: SPD3D; 3: Tracks
Definition: AliRDHFCuts.h:401
centrality
char Char_t
Definition: External.C:18
static Int_t CheckMatchingAODdeltaAODevents()
Float_t fMaxVtxRedChi2
minimum vertex contributors
Definition: AliRDHFCuts.h:402
void SetMClowenpp2011(Bool_t mc)
Definition: AliAODPidHF.h:102
const Float_t * GetCuts() const
Definition: AliRDHFCuts.h:254
void SetMinCrossedRowsTPCPtDep(const char *rows="")
Float_t fCutRatioClsOverCrossRowsTPC
histogram with reference centrality distribution for centrality distribution flattening ...
Definition: AliRDHFCuts.h:462
Double_t fDeadZoneWidth
flag for enabling/disabling geometrical cut on TPC track
Definition: AliRDHFCuts.h:467
Int_t fMinVtxType
cuts on the event
Definition: AliRDHFCuts.h:400
void SetppLowEn2011(Bool_t opt)
Definition: AliAODPidHF.h:104
Double_t fMaxRapidityCand
minimum pt of the candidate
Definition: AliRDHFCuts.h:449
Double_t fCutGeoNcrNclFractionNcl
4th parameter of GeoNcrNcl cut
Definition: AliRDHFCuts.h:471
void SetGlobalIndex()
Definition: AliRDHFCuts.h:201
AliRDHFCuts & operator=(const AliRDHFCuts &source)
Float_t GetCutValue(Int_t iVar, Int_t iPtBin) const
Int_t fUseCentrality
flag to use a multiplicity dependent pileup selection
Definition: AliRDHFCuts.h:436
Bool_t fUsePID
Definition: AliRDHFCuts.h:424
TRandom * gRandom
Bool_t fRejectPlpFromDiffBCMV
minimum weighted distance in Z between 2 vertices (multi-vertexer)
Definition: AliRDHFCuts.h:407
void PrintAll() const
void SetPidHF(AliAODPidHF *pidObj)
see enum below
Definition: AliRDHFCuts.h:215
Bool_t GetOldPid()
Definition: AliAODPidHF.h:157
TList * GetList(const TDirectory *d, const char *what)
Definition: ExtractData.C:172
AliAODPidHF * GetPidHF() const
Definition: AliRDHFCuts.h:247
void RecalculateImpPars(AliAODVertex *vtxAODNew, AliAODEvent *aod)
Bool_t fRemoveTrackletOutliers
Max. difference between TRK and V0 centrality (remove TPC pileup for PbPb 2011)
Definition: AliRDHFCuts.h:456
Bool_t fUseV0ANDSelectionOffline
5th parameter of GeoNcrNcl cut
Definition: AliRDHFCuts.h:472
AliRDHFCuts(const Char_t *name="RDHFCuts", const Char_t *title="")
Definition: AliRDHFCuts.cxx:64
const Int_t nPtBins
ULong64_t fTriggerMask
flag to reject pileup from different BC (multi-vertexer)
Definition: AliRDHFCuts.h:408
Bool_t fUseCutGeoNcrNcl
pT-dep cut in TPC minimum n crossed rows
Definition: AliRDHFCuts.h:466
int Int_t
Definition: External.C:63
Bool_t fUseTrackSelectionWithFilterBits
flag to reject kink daughters
Definition: AliRDHFCuts.h:459
void SetCuts(Int_t nVars, Int_t nPtBins, Float_t **cutsRD)
Int_t fnVarsForOpt
Definition: AliRDHFCuts.h:419
Float_t fMinWDzPileupMV
max chi2 per contributor of the pile-up vertex to consider (multi-vertexer).
Definition: AliRDHFCuts.h:406
Double_t fCutGeoNcrNclLength
1st parameter of GeoNcrNcl cut
Definition: AliRDHFCuts.h:468
float Float_t
Definition: External.C:68
virtual ~AliRDHFCuts()
void SetPbPb(Bool_t pbpb)
Definition: AliAODPidHF.h:105
TString * fVarNames
number of cut vars for candidates
Definition: AliRDHFCuts.h:418
Double_t fMaxPtCand
minimum pt of the candidate
Definition: AliRDHFCuts.h:448
AliESDtrackCuts * fTrackCuts
quality cuts on the daughter tracks
Definition: AliRDHFCuts.h:412
void MakeTable() const
Bool_t fFixRefs
name of the AliMultSelection object to be considered
Definition: AliRDHFCuts.h:444
Bool_t fKeepSignalMC
max rapidity of candidate (if !=-999 overrides IsInFiducialAcceptance)
Definition: AliRDHFCuts.h:450
Int_t fOptPileup
use Physics selection criteria
Definition: AliRDHFCuts.h:432
AliESDtrackCuts * GetTrackCuts() const
Definition: AliRDHFCuts.h:262
void PrintTrigger() const
Float_t fMaxVtxZ
maximum chi2/ndf
Definition: AliRDHFCuts.h:403
Float_t GetCentrality(AliAODEvent *aodEvent)
Definition: AliRDHFCuts.h:258
Int_t fIsSelectedPID
outcome of cuts selection
Definition: AliRDHFCuts.h:446
Float_t fMinCentrality
Definition: AliRDHFCuts.h:441
AliAODVertex * RemoveDaughtersFromPrimaryVtx(AliAODEvent *aod)
void SetUpCombinedPID()
AliPIDResponse * GetPidResponse() const
Definition: AliAODPidHF.h:160
UInt_t fEvRejectionBits
used to code the step at which candidate was rejected
Definition: AliRDHFCuts.h:428
Double_t fMaxDiffTRKV0Centr
flag to apply cut on tracks crossing SPD misaligned modules for PP2012 data
Definition: AliRDHFCuts.h:455
Float_t GetCentralityOldFramework(AliAODEvent *aodEvent, AliRDHFCuts::ECentrality estimator)
Int_t fnPtBinLimits
number of pt bins for cuts
Definition: AliRDHFCuts.h:415
Float_t * fCutsRD
fnVars*fnPtBins
Definition: AliRDHFCuts.h:422
Bool_t CheckPtDepCrossedRows(TString rows, Bool_t print=kFALSE) const
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
void SetOwnPrimaryVtx(const AliAODVertex *vtx)
Bool_t fUseCentrFlatteningInMC
flag to enable/disable the check on filter bits
Definition: AliRDHFCuts.h:460
void SetVarsForOpt(Int_t nVars, Bool_t *forOpt)
void SetVarNames(Int_t nVars, TString *varNames, Bool_t *isUpperCut)
void GetVarPtIndex(Int_t iGlob, Int_t &iVar, Int_t &iPtBin) const
void SetBetheBloch()
Bool_t * fVarsForOpt
number of cut vars to be optimized for candidates
Definition: AliRDHFCuts.h:420
Bool_t isMC
Bool_t CompareCuts(const AliRDHFCuts *obj) const
Bool_t AreDaughtersSelected(AliAODRecoDecayHF *rd, const AliAODEvent *aod=0x0) const
Bool_t fApplySPDMisalignedPP2012
flag to apply SPD dead module map of PbPb2011
Definition: AliRDHFCuts.h:454
Bool_t IsEventSelected(AliVEvent *event)
Int_t fMinSPDMultiplicity
maximum |z| of primary vertex
Definition: AliRDHFCuts.h:404
void SetMC(Bool_t mc)
Definition: AliAODPidHF.h:101
virtual void PrintAll() const
Bool_t fUseAOD049
enable PID usage (off by default)
Definition: AliRDHFCuts.h:425
Bool_t fApplySPDDeadPbPb2011
maximum pt of the candidate for which to check if the daughters fulfill kFirst criteria ...
Definition: AliRDHFCuts.h:453
void CleanOwnPrimaryVtx(AliAODRecoDecayHF *d, AliAODEvent *aod, AliAODVertex *origownvtx) const
Float_t fMaxVtxChi2PileupMV
SPD multiplicity.
Definition: AliRDHFCuts.h:405
Bool_t fIsCandTrackSPDFirst
IsSelected returns always kTRUE for MC signal.
Definition: AliRDHFCuts.h:451
TH1F * fHistCentrDistr
flag for enabling/diabling centrality flattening in MC
Definition: AliRDHFCuts.h:461
void SetPtBins(Int_t nPtBinLimits, Float_t *ptBinLimits)
Bool_t fKinkReject
cut on zSPD vertex to remove outliers in centrality vs. tracklets (0=no cut, 1= cut at 12 cm...
Definition: AliRDHFCuts.h:458
TString fMultSelectionObjectName
maximum centrality for selected events
Definition: AliRDHFCuts.h:443
Int_t fMinContrPileup
option for pielup selection
Definition: AliRDHFCuts.h:433
Bool_t fUseTPCtrackCutsOnThisDaughter
flag to apply V0AND selection offline
Definition: AliRDHFCuts.h:473
Bool_t fUsePhysicsSelection
use MC primary vertex
Definition: AliRDHFCuts.h:431
TFormula * f1CutMinNCrossedRowsTPCPtDep
pT-dep cut in TPC minimum n crossed rows
Definition: AliRDHFCuts.h:465
AliAODVertex * GetPrimaryVtx() const
void AddTrackCuts(const AliESDtrackCuts *cuts)
Definition: AliRDHFCuts.h:207
TString fTriggerClass[2]
flag to select one trigger only
Definition: AliRDHFCuts.h:410
bool Bool_t
Definition: External.C:53
Int_t fnPtBins
cuts on the candidate
Definition: AliRDHFCuts.h:414
Double_t fCutGeoNcrNclGeom1Pt
2nd parameter of GeoNcrNcl cut
Definition: AliRDHFCuts.h:469
AliAODPidHF * fPidHF
enable AOD049 centrality cleanup
Definition: AliRDHFCuts.h:426
Float_t fCutRatioSignalNOverCrossRowsTPC
min. value ratio NTPCClusters/NTPCCrossedRows, cut if !=0
Definition: AliRDHFCuts.h:463
Int_t fGlobalIndex
Definition: AliRDHFCuts.h:421
Bool_t RecalcOwnPrimaryVtx(AliAODRecoDecayHF *d, AliAODEvent *aod) const
Int_t PtBin(Double_t pt) const
Int_t GetGlobalIndex(Int_t iVar, Int_t iPtBin) const
TString fCutMinCrossedRowsTPCPtDep
min. value ratio TPCPointsUsedForPID/NTPCCrossedRows, cut if !=0
Definition: AliRDHFCuts.h:464
Double_t fMinPtCand
outcome of PID selection
Definition: AliRDHFCuts.h:447
Bool_t IsDaughterSelected(AliAODTrack *track, const AliESDVertex *primary, AliESDtrackCuts *cuts, const AliAODEvent *aod=0x0) const
void SetOnePad(Bool_t onepad)
Definition: AliAODPidHF.h:103
Float_t fMaxCentrality
minimum centrality for selected events
Definition: AliRDHFCuts.h:442
void SetPidResponse(AliPIDResponse *pidResp)
Definition: AliAODPidHF.h:111
Bool_t RecomputePrimaryVertex(AliAODEvent *event) const
Int_t fnVars
Definition: AliRDHFCuts.h:417