AliPhysics  9b6b435 (9b6b435)
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(3),
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 fApplyZcutOnSPDvtx(kFALSE),
134 fUsePreselect(0)
135 {
136  //
137  // Default Constructor
138  //
139  fTriggerClass[0]="CINT1"; fTriggerClass[1]="";
140 }
141 //--------------------------------------------------------------------------
143  AliAnalysisCuts(source),
144  fMinVtxType(source.fMinVtxType),
145  fMinVtxContr(source.fMinVtxContr),
147  fMaxVtxZ(source.fMaxVtxZ),
152  fTriggerMask(source.fTriggerMask),
154  fTriggerClass(),
155  fTrackCuts(0),
156  fnPtBins(source.fnPtBins),
158  fPtBinLimits(0),
159  fnVars(source.fnVars),
160  fVarNames(0),
161  fnVarsForOpt(source.fnVarsForOpt),
162  fVarsForOpt(0),
163  fGlobalIndex(source.fGlobalIndex),
164  fCutsRD(0),
165  fIsUpperCut(0),
166  fUsePID(source.fUsePID),
167  fUseAOD049(source.fUseAOD049),
168  fPidHF(0),
172  fUseMCVertex(source.fUseMCVertex),
174  fOptPileup(source.fOptPileup),
176  fMinDzPileup(source.fMinDzPileup),
182  fFixRefs(source.fFixRefs),
185  fMinPtCand(source.fMinPtCand),
186  fMaxPtCand(source.fMaxPtCand),
196  fKinkReject(source.fKinkReject),
199  fHistCentrDistr(0x0),
214 {
215  //
216  // Copy constructor
217  //
218  cout<<"Copy constructor"<<endl;
219  fTriggerClass[0] = source.fTriggerClass[0];
220  fTriggerClass[1] = source.fTriggerClass[1];
221  if(source.GetTrackCuts()) AddTrackCuts(source.GetTrackCuts());
222  if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);
223  if(source.fVarNames) SetVarNames(source.fnVars,source.fVarNames,source.fIsUpperCut);
224  if(source.fCutsRD) SetCuts(source.fGlobalIndex,source.fCutsRD);
225  if(source.fVarsForOpt) SetVarsForOpt(source.fnVarsForOpt,source.fVarsForOpt);
226  if(source.fPidHF) SetPidHF(source.fPidHF);
227  if(source.fHistCentrDistr) fHistCentrDistr=(TH1F*)(source.fHistCentrDistr->Clone());
230  PrintAll();
231 
232 }
233 //--------------------------------------------------------------------------
235 {
236  //
237  // assignment operator
238  //
239  if(&source == this) return *this;
240 
241  AliAnalysisCuts::operator=(source);
242 
243  fMinVtxType=source.fMinVtxType;
244  fMinVtxContr=source.fMinVtxContr;
246  fMaxVtxZ=source.fMaxVtxZ;
251  fTriggerMask=source.fTriggerMask;
253  fTriggerClass[0]=source.fTriggerClass[0];
254  fTriggerClass[1]=source.fTriggerClass[1];
255  fnPtBins=source.fnPtBins;
257  fnVars=source.fnVars;
258  fGlobalIndex=source.fGlobalIndex;
259  fnVarsForOpt=source.fnVarsForOpt;
260  fUsePID=source.fUsePID;
261  fUseAOD049=source.fUseAOD049;
262  if(fPidHF) delete fPidHF;
263  fPidHF=new AliAODPidHF(*(source.GetPidHF()));
267  fUseMCVertex=source.fUseMCVertex;
269  fOptPileup=source.fOptPileup;
271  fMinDzPileup=source.fMinDzPileup;
277  fFixRefs=source.fFixRefs;
280  fMinPtCand=source.fMinPtCand;
281  fMaxPtCand=source.fMaxPtCand;
291  fKinkReject=source.fKinkReject;
293  if(fHistCentrDistr) delete fHistCentrDistr;
295  if(source.fHistCentrDistr)fHistCentrDistr=(TH1F*)(source.fHistCentrDistr->Clone());
296 
297  if(source.GetTrackCuts()) {delete fTrackCuts; fTrackCuts=new AliESDtrackCuts(*(source.GetTrackCuts()));}
298  if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);
299  if(source.fVarNames) SetVarNames(source.fnVars,source.fVarNames,source.fIsUpperCut);
300  if(source.fCutsRD) SetCuts(source.fGlobalIndex,source.fCutsRD);
301  if(source.fVarsForOpt) SetVarsForOpt(source.fnVarsForOpt,source.fVarsForOpt);
302 
308 
318 
319  PrintAll();
320 
321  return *this;
322 }
323 //--------------------------------------------------------------------------
325  //
326  // Default Destructor
327  //
328  if(fTrackCuts) { delete fTrackCuts; fTrackCuts=0; }
329  if(fPtBinLimits) {delete [] fPtBinLimits; fPtBinLimits=0;}
330  if(fVarNames) {delete [] fVarNames; fVarNames=0;}
331  if(fVarsForOpt) {delete [] fVarsForOpt; fVarsForOpt=0;}
332  if(fCutsRD) {
333  delete [] fCutsRD;
334  fCutsRD=0;
335  }
336  if(fIsUpperCut) {delete [] fIsUpperCut; fIsUpperCut=0;}
337  if(fPidHF){
338  delete fPidHF;
339  fPidHF=0;
340  }
342 
346  }
347 
348 }
349 //---------------------------------------------------------------------------
351  //
352  // Centrality selection
353  //
354  if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid){
355  AliWarning("Centrality estimator not valid");
356  return 3;
357  }else{
358  Float_t centvalue=GetCentrality((AliAODEvent*)event);
359  if (centvalue<-998.){//-999 if no centralityP
360  return 3;
362  return 3;
363  }else{
364  if (centvalue<fMinCentrality || centvalue>fMaxCentrality){
365  return 2;
366  }
367  // selection to flatten the centrality distribution
368  if(fHistCentrDistr){
369  if(!IsEventSelectedForCentrFlattening(centvalue))return 4;
370  }
371  }
372  }
373  return 0;
374 }
375 
376 
377 //-------------------------------------------------
378 void AliRDHFCuts::SetHistoForCentralityFlattening(TH1F *h,Double_t minCentr,Double_t maxCentr,Double_t centrRef,Int_t switchTRand){
379  // set the histo for centrality flattening
380  // the centrality is flatten in the range minCentr,maxCentr
381  // if centrRef is zero, the minimum in h within (minCentr,maxCentr) defines the reference
382  // positive, the value of h(centrRef) defines the reference (-> the centrality distribution might be not flat in the whole desired range)
383  // 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)
384  // 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
385 
386  if(maxCentr<minCentr){
387  AliWarning("AliRDHFCuts::Wrong centralities values while setting the histogram for centrality flattening");
388  }
389 
391  fHistCentrDistr=(TH1F*)h->Clone("hCentralityFlat");
392  fHistCentrDistr->SetTitle("Reference histo for centrality flattening");
393  Int_t minbin=fHistCentrDistr->FindBin(minCentr*1.00001); // fast if fix bin width
394  Int_t maxbin=fHistCentrDistr->FindBin(maxCentr*0.9999);
395  fHistCentrDistr->GetXaxis()->SetRange(minbin,maxbin);
396  Double_t ref=0.,bincont=0.,binrefwidth=1.;
397  Int_t binref=0;
398  if(TMath::Abs(centrRef)<0.0001){
399  binref=fHistCentrDistr->GetMinimumBin();
400  binrefwidth=fHistCentrDistr->GetBinWidth(binref);
401  ref=fHistCentrDistr->GetBinContent(binref)/binrefwidth;
402  }
403  else if(centrRef>0.){
404  binref=h->FindBin(centrRef);
405  if(binref<1||binref>h->GetNbinsX()){
406  AliWarning("AliRDHFCuts::Wrong centrality reference value while setting the histogram for centrality flattening");
407  }
408  binrefwidth=fHistCentrDistr->GetBinWidth(binref);
409  ref=fHistCentrDistr->GetBinContent(binref)/binrefwidth;
410  }
411  else{
412  if(centrRef<-1) AliWarning("AliRDHFCuts: with this centrality reference no flattening will be applied");
413  binref=fHistCentrDistr->GetMaximumBin();
414  binrefwidth=fHistCentrDistr->GetBinWidth(binref);
415  ref=fHistCentrDistr->GetMaximum()*TMath::Abs(centrRef)/binrefwidth;
416  }
417 
418  for(Int_t j=1;j<=h->GetNbinsX();j++){// Now set the "probabilities"
419  if(h->GetBinLowEdge(j)*1.0001>=minCentr&&h->GetBinLowEdge(j+1)*0.9999<=maxCentr){
420  bincont=h->GetBinContent(j);
421  fHistCentrDistr->SetBinContent(j,ref/bincont*h->GetBinWidth(j));
422  fHistCentrDistr->SetBinError(j,h->GetBinError(j)*ref/bincont);
423  }
424  else{
425  h->SetBinContent(j,1.1);// prob > 1 to assure that events will not be rejected
426  }
427  }
428 
429  fHistCentrDistr->SetBinContent(0,switchTRand);
430  return;
431 
432 }
433 
434 //-------------------------------------------------
436  //
437  // Random event selection, based on fHistCentrDistr, to flatten the centrality distribution
438  // Can be faster if it was required that fHistCentrDistr covers
439  // exactly the desired centrality range (e.g. part of the lines below should be done during the
440  // setting of the histo) and TH1::SetMinimum called
441  //
442 
443  if(!fHistCentrDistr) return kTRUE;
444  // Int_t maxbin=fHistCentrDistr->FindBin(fMaxCentrality*0.9999);
445  // if(maxbin>fHistCentrDistr->GetNbinsX()){
446  // AliWarning("AliRDHFCuts: The maximum centrality exceeds the x-axis limit of the histogram for centrality flattening");
447  // }
448 
449  Int_t bin=fHistCentrDistr->FindBin(centvalue); // Fast if the histo has a fix bin
450  Double_t bincont=fHistCentrDistr->GetBinContent(bin);
451  Double_t centDigits=centvalue-(Int_t)(centvalue*100.)/100.;// this is to extract a random number between 0 and 0.01
452 
453  if(fHistCentrDistr->GetBinContent(0)<-0.9999){
454  if(gRandom->Uniform(1.)<bincont)return kTRUE;
455  return kFALSE;
456  }
457 
458  if(centDigits*100.<bincont)return kTRUE;
459  return kFALSE;
460 
461 }
462 //---------------------------------------------------------------------------
463 void AliRDHFCuts::SetupPID(AliVEvent *event) {
464  // Set the PID response object in the AliAODPidHF
465  // in case of old PID sets the TPC dE/dx BB parameterization
466 
467  if(fPidHF){
468  if(fPidHF->GetPidResponse()==0x0){
469  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
470  AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
471  AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
472  fPidHF->SetPidResponse(pidResp);
473  }
475  if(fPidHF->GetOldPid()) {
476 
477  Bool_t isMC=kFALSE;
478  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)event)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
479  if(mcArray) {isMC=kTRUE;fUseAOD049=kFALSE;}
480 
481  // pp, from LHC10d onwards
482  if((event->GetRunNumber()>121693 && event->GetRunNumber()<136851) ||
483  event->GetRunNumber()>139517) fPidHF->SetOnePad(kTRUE);
484  // pp, 2011 low energy run
485  if((event->GetRunNumber()>=146686 && event->GetRunNumber()<=146860)){
486  fPidHF->SetppLowEn2011(kTRUE);
487  fPidHF->SetOnePad(kFALSE);
488  }
489  // PbPb LHC10h
490  if(event->GetRunNumber()>=136851 && event->GetRunNumber()<=139517) fPidHF->SetPbPb(kTRUE);
491  // MC
492  if(isMC) fPidHF->SetMC(kTRUE);
493  if(isMC && (event->GetRunNumber()>=146686 && event->GetRunNumber()<=146860))
494  fPidHF->SetMClowenpp2011(kTRUE);
496  }else{
497  // check that AliPIDResponse object was properly set in case of using OADB
498  if(fPidHF->GetPidResponse()==0x0) AliFatal("AliPIDResponse object not set");
499  }
500  }
501 }
502 //---------------------------------------------------------------------------
504  //
505  // Event selection
506  //
507  //if(fTriggerMask && event->GetTriggerMask()!=fTriggerMask) return kFALSE;
508 
509 
510 
511  fWhyRejection=0;
513  Bool_t accept=kTRUE;
514 
515 
516  // check if it's MC
517  Bool_t isMC=kFALSE;
518  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)event)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
519  if(mcArray) {isMC=kTRUE;fUseAOD049=kFALSE;}
520 
521 
522  SetupPID(event);
523 
524  // trigger class
525  TString firedTriggerClasses=((AliAODEvent*)event)->GetFiredTriggerClasses();
526  // don't do for MC and for PbPb 2010 data
527  if(!isMC && (event->GetRunNumber()<136851 || event->GetRunNumber()>139517)) {
528  if(!firedTriggerClasses.Contains(fTriggerClass[0].Data()) &&
529  (fTriggerClass[1].CompareTo("")==0 || !firedTriggerClasses.Contains(fTriggerClass[1].Data())) ) {
530  fWhyRejection=5;
532  accept=kFALSE;
533  }
534  }
535 
536  // TEMPORARY FIX FOR GetEvent
537  Int_t nTracks=((AliAODEvent*)event)->GetNumberOfTracks();
538  for(Int_t itr=0; itr<nTracks; itr++){
539  AliAODTrack* tr=(AliAODTrack*)((AliAODEvent*)event)->GetTrack(itr);
540  tr->SetAODEvent((AliAODEvent*)event);
541  }
542 
543  // TEMPORARY FIX FOR REFERENCES
544  // Fix references to daughter tracks
545  // if(fFixRefs) {
546  // AliAnalysisVertexingHF *fixer = new AliAnalysisVertexingHF();
547  // fixer->FixReferences((AliAODEvent*)event);
548  // delete fixer;
549  // }
550  //
551 
552 
553  // physics selection requirements
555  Bool_t isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & fTriggerMask);
556  if(!isSelected) {
557  if(accept) fWhyRejection=7;
559  accept=kFALSE;
560  }else{
561  if(fUseOnlyOneTrigger){
562  if(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected()!=fTriggerMask){
563  if(accept) fWhyRejection=7;
565  accept=kFALSE;
566  }
567  }
569  AliAODVZERO* v0data=(AliAODVZERO*)((AliAODEvent*)event)->GetVZEROData();
570  Int_t tv0a=v0data->GetV0ADecision();
571  Int_t tv0c=v0data->GetV0CDecision();
572  if(!(tv0a==1 && tv0c==1)){
573  if(accept) fWhyRejection=7;
575  accept=kFALSE;
576  }
577  }
578  }
579  }
580 
581  // centrality selection
582  if (fUseCentrality!=kCentOff) {
583  Int_t rejection=IsEventSelectedInCentrality(event);
584  Bool_t okCent=kFALSE;
585  if(rejection==0) okCent=kTRUE;
586  if(isMC && rejection==4 && !fUseCentrFlatteningInMC) okCent=kTRUE;
587  if(!okCent){
588  if(accept) fWhyRejection=rejection;
591  accept=kFALSE;
592  }
593 
594  }
595 
596  // PbPb2011 outliers in tracklets vs. VZERO and centTRK vs. centV0
597  if(event->GetRunNumber()>=167693 && event->GetRunNumber()<=170593){
599  Double_t v0cent=GetCentrality((AliAODEvent*)event,kCentV0M);
600  Double_t ntracklets=((AliAODEvent*)event)->GetTracklets()->GetNumberOfTracklets();
601  Double_t cutval=60.-0.08*ntracklets+1./50000.*ntracklets*ntracklets;
602  if(ntracklets<1000. && v0cent<cutval){
603  if(accept) fWhyRejection=2;
605  accept=kFALSE;
606  }
607  }
608  if(fMaxDiffTRKV0Centr>0.){
609  Double_t v0cent=GetCentrality((AliAODEvent*)event,kCentV0M);
610  Double_t trkcent=GetCentrality((AliAODEvent*)event,kCentTRK);
611  if(TMath::Abs(trkcent-v0cent)>fMaxDiffTRKV0Centr){
612  if(accept) fWhyRejection=1;
614  accept=kFALSE;
615  }
616  }
617  }
618 
619  // vertex requirements
620 
621  const AliVVertex *vertex = event->GetPrimaryVertex();
622 
623  if(!vertex){
624  accept=kFALSE;
626  }else{
627  TString title=vertex->GetTitle();
628  if(title.Contains("Z") && fMinVtxType>1){
629  accept=kFALSE;
631  }
632  else if(title.Contains("3D") && fMinVtxType>2){
633  accept=kFALSE;
635  }
636  if(vertex->GetNContributors()<fMinVtxContr){
637  accept=kFALSE;
639  }
640  }
641 
642  const AliVVertex *vSPD = ((AliAODEvent*)event)->GetPrimaryVertexSPD();
643  if(fCutOnzVertexSPD>0){
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  Double_t zvert = -999;
682  if(!fApplyZcutOnSPDvtx && vertex) zvert=vertex->GetZ();
683  else if(fApplyZcutOnSPDvtx) {
684  if(!vSPD || (vSPD && vSPD->GetNContributors()<1)){
685  accept=kFALSE;
687  }
688  else zvert = vSPD->GetZ();
689  }
690 
691  if(TMath::Abs(zvert)>fMaxVtxZ) {
693  if(accept) fWhyRejection=6;
694  accept=kFALSE;
695  }
696 
697  // pile-up rejection
699  Bool_t isPileup=kFALSE;
701  isPileup=event->IsPileupFromSPDInMultBins();
702  }else{
705  isPileup=event->IsPileupFromSPD(cutc,cutz,3.,2.,10.);
706  }
707  if(isPileup){
708  if(accept) fWhyRejection=1;
710  accept=kFALSE;
711  }
712  }
714  AliAnalysisUtils utils;
715  utils.SetMinPlpContribMV(fMinContrPileup); // min. multiplicity of the pile-up vertex to consider
716  utils.SetMaxPlpChi2MV(fMaxVtxChi2PileupMV); // max chi2 per contributor of the pile-up vertex to consider.
717  utils.SetMinWDistMV(fMinWDzPileupMV); // minimum weighted distance in Z between 2 vertices (i.e. (zv1-zv2)/sqrt(sigZv1^2+sigZv2^2) )
718  utils.SetCheckPlpFromDifferentBCMV(fRejectPlpFromDiffBCMV); // vertex with |BCID|>2 will trigger pile-up (based on TOF)
719  Bool_t isPUMV = utils.IsPileUpMV(event);
720  if(isPUMV) {
721  if(accept) fWhyRejection=1;
723  accept=kFALSE;
724  }
725  }
726 
727 
728  // Correcting PP2012 flag to remoce tracks crossing SPD misaligned staves for periods 12def
729  if(fApplySPDMisalignedPP2012 && !(event->GetRunNumber()>=195681 && event->GetRunNumber()<=197388)) fApplySPDMisalignedPP2012=false;
730 
731  return accept;
732 }
733 //---------------------------------------------------------------------------
735  //
736  // Daughter tracks selection
737  //
738  if(!fTrackCuts) return kTRUE;
739 
740  Int_t ndaughters = d->GetNDaughters();
741  AliAODVertex *vAOD = d->GetPrimaryVtx();
742  Double_t pos[3],cov[6];
743  vAOD->GetXYZ(pos);
744  vAOD->GetCovarianceMatrix(cov);
745  const AliESDVertex vESD(pos,cov,100.,100);
746 
747  Bool_t retval=kTRUE;
748 
749  for(Int_t idg=0; idg<ndaughters; idg++) {
750  AliAODTrack *dgTrack = (AliAODTrack*)d->GetDaughter(idg);
751  if(!dgTrack) {retval = kFALSE; continue;}
752  //printf("charge %d\n",dgTrack->Charge());
753  if(dgTrack->Charge()==0) continue; // it's not a track, but a V0
754 
756  { if(!dgTrack->HasPointOnITSLayer(0)) { retval = kFALSE; continue; } }
757 
758  if(!IsDaughterSelected(dgTrack,&vESD,fTrackCuts,aod)) retval = kFALSE;
759  }
760 
761  return retval;
762 }
763 //---------------------------------------------------------------------------
765  //
766  // Check if AOD and deltaAOD files are composed of the same events:
767  // mismatches were observed in the merged AODs of LHC15o
768  //
769  // When AOD+deltaAOD are produced from ESD, mismatches can be found looking at:
770  // - the AOD trees in AliAOD.root and AliAOD.VertexingHF.root have different number of entries
771  // - the titles of the TProcessID objects do not match
772  // When deltaAOD are produced from AOD, mismatches can be found looking at:
773  // - the AOD trees in AliAOD.root and AliAOD.VertexingHF.root have different number of entries
774  //
775  // Return values:
776  // -1: AOD and deltaAOD trees have different number of entries
777  // 0: AOD and deltaAOD trees have same number of entries + the titles of the TProcessID objects do not match
778  // 1: AOD and deltaAOD trees have same number of entries + the titles of the TProcessID objects match
779  Bool_t okTProcessNames = kTRUE;
780  AliAODHandler* aodHandler = (AliAODHandler*)((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
781  TTree *treeAOD = aodHandler->GetTree();
782  TTree *treeDeltaAOD = treeAOD->GetFriend("aodTree");
783  if(!treeDeltaAOD || !treeAOD) return -1;
784  if(treeDeltaAOD && treeAOD){
785  if(treeAOD->GetEntries()!=treeDeltaAOD->GetEntries()){
786  printf("AliRDHFCuts::CheckMatchingAODdeltaAODevents: Difference in number of entries in main and friend tree, skipping event\n");
787  return -1;
788  }
789  TFile *mfile = treeAOD->GetCurrentFile();
790  TFile *dfile = treeDeltaAOD->GetCurrentFile();
791  TList* lm=mfile->GetListOfKeys();
792  TList* ld=dfile->GetListOfKeys();
793  Int_t nentm=lm->GetEntries();
794  for(Int_t jm=0; jm<nentm; jm++){
795  TKey* o=(TKey*)lm->At(jm);
796  TString clnam=o->GetClassName();
797  if(clnam=="TProcessID"){
798  TString pname=o->GetName();
799  TString ptit=o->GetTitle();
800  if(pname.Contains("ProcessID")){
801  TObject* od=(TObject*)ld->FindObject(pname.Data());
802  if(od){
803  TString ptit2=od->GetTitle();
804  if(ptit2!=ptit){
805  printf("AliRDHFCuts::CheckMatchingAODdeltaAODevents: mismatch in %s: AOD: %s -- deltaAOD: %s\n",pname.Data(),ptit.Data(),ptit2.Data());
806  okTProcessNames = kFALSE;
807  }
808  }
809  }
810  }
811  }
812  }
813 
814  if (okTProcessNames) return 1;
815  else return 0;
816 }
817 //---------------------------------------------------------------------------
819  //
820  // Check the correctness of the string syntax
821  //
822  Bool_t retval=kTRUE;
823 
824  if(!rows.Contains("pt")) {
825  if(print) AliError("string must contain \"pt\"");
826  retval= kFALSE;
827  }
828  return retval;
829 }
830 //---------------------------------------------------------------------------
832  //
833  //Create the TFormula from TString for TPC crossed rows pT dependent cut
834  //
835 
836 
837  // setting data member that describes the TPC crossed rows pT dependent cut
839 
840  // creating TFormula from TString
843  // resetting TFormula
845  }
846  if(!CheckPtDepCrossedRows(rows,kTRUE))return;
847 
848  TString tmp(rows);
849  tmp.ReplaceAll("pt","x");
850  f1CutMinNCrossedRowsTPCPtDep = new TFormula("f1CutMinNCrossedRowsTPCPtDep",tmp.Data());
851 
852 
853 }
854 //---------------------------------------------------------------------------
855 Bool_t AliRDHFCuts::IsDaughterSelected(AliAODTrack *track,const AliESDVertex *primary,AliESDtrackCuts *cuts, const AliAODEvent* aod) const{
856  //
857  // Convert to ESDtrack, relate to vertex and check cuts
858  //
859  if(!cuts) return kTRUE;
860 
861  if(cuts->GetFlagCutTOFdistance()) cuts->SetFlagCutTOFdistance(kFALSE);
862 
863 
864  // convert to ESD track here
865  AliESDtrack esdTrack(track);
866  // set the TPC cluster info
867  esdTrack.SetTPCClusterMap(track->GetTPCClusterMap());
868  esdTrack.SetTPCSharedMap(track->GetTPCSharedMap());
869  esdTrack.SetTPCPointsF(track->GetTPCNclsF());
870  // needed to calculate the impact parameters
871  esdTrack.RelateToVertex(primary,0.,3.);
872 
873  //applying ESDtrackCut
874  if(!cuts->IsSelected(&esdTrack)) return kFALSE;
875 
876  //appliyng kink rejection
877  if(fKinkReject){
878  AliAODVertex *maybeKink=track->GetProdVertex();
879  if(maybeKink->GetType()==AliAODVertex::kKink) return kFALSE;
880  }
881 
882  //appliyng TPC crossed rows pT dependent cut
884  Float_t nCrossedRowsTPC = esdTrack.GetTPCCrossedRows();
885  if(nCrossedRowsTPC<f1CutMinNCrossedRowsTPCPtDep->Eval(esdTrack.Pt())) return kFALSE;
886  }
887 
888  //appliyng NTPCcls/NTPCcrossedRows cut
890  Float_t nCrossedRowsTPC = esdTrack.GetTPCCrossedRows();
891  Float_t nClustersTPC = esdTrack.GetTPCNcls();
892  if(nCrossedRowsTPC!=0){
893  Float_t ratio = nClustersTPC/nCrossedRowsTPC;
894  if(ratio<fCutRatioClsOverCrossRowsTPC) return kFALSE;
895  }
896  else return kFALSE;
897  }
898 
899  //appliyng TPCsignalN/NTPCcrossedRows cut
901  Float_t nCrossedRowsTPC = esdTrack.GetTPCCrossedRows();
902  Float_t nTPCsignal = esdTrack.GetTPCsignalN();
903  if(nCrossedRowsTPC!=0){
904  Float_t ratio = nTPCsignal/nCrossedRowsTPC;
905  if(ratio<fCutRatioSignalNOverCrossRowsTPC) return kFALSE;
906  }
907  else return kFALSE;
908  }
909 
910  // geometrical cut (note uses track at vertex instead of at TPC inner wall)
912  Float_t nCrossedRowsTPC = esdTrack.GetTPCCrossedRows();
913  Float_t lengthInActiveZoneTPC=esdTrack.GetLengthInActiveZone(0,fDeadZoneWidth,220.,aod->GetMagneticField());
914  Double_t cutGeoNcrNclLength=fCutGeoNcrNclLength-TMath::Power(TMath::Abs(esdTrack.GetSigned1Pt()),fCutGeoNcrNclGeom1Pt);
915  Bool_t isOK=kTRUE;
916  if (lengthInActiveZoneTPC<cutGeoNcrNclLength) isOK=kFALSE;
917  if (nCrossedRowsTPC<fCutGeoNcrNclFractionNcr*cutGeoNcrNclLength) isOK=kFALSE;
918  if (esdTrack.GetTPCncls()<fCutGeoNcrNclFractionNcl*cutGeoNcrNclLength) isOK=kFALSE;
919  if(!isOK) return kFALSE;
920  }
921 
922 
924  // to be implemented
925  // we need either to have here the AOD Event,
926  // or to have the pileup vertex object
927  }
928 
930  Bool_t deadSPDLay1PbPb2011[20][4]={
931  {kTRUE,kTRUE,kTRUE,kTRUE},
932  {kTRUE,kTRUE,kTRUE,kTRUE},
933  {kTRUE,kTRUE,kTRUE,kTRUE},
934  {kTRUE,kTRUE,kTRUE,kTRUE},
935  {kTRUE,kTRUE,kTRUE,kTRUE},
936  {kFALSE,kFALSE,kTRUE,kTRUE},
937  {kTRUE,kTRUE,kFALSE,kFALSE},
938  {kTRUE,kTRUE,kTRUE,kTRUE},
939  {kFALSE,kFALSE,kTRUE,kTRUE},
940  {kTRUE,kTRUE,kTRUE,kTRUE},
941  {kTRUE,kTRUE,kFALSE,kFALSE},
942  {kTRUE,kTRUE,kTRUE,kTRUE},
943  {kFALSE,kFALSE,kFALSE,kFALSE},
944  {kFALSE,kFALSE,kTRUE,kTRUE},
945  {kFALSE,kFALSE,kFALSE,kFALSE},
946  {kFALSE,kFALSE,kFALSE,kFALSE},
947  {kTRUE,kTRUE,kTRUE,kTRUE},
948  {kTRUE,kTRUE,kFALSE,kFALSE},
949  {kFALSE,kFALSE,kFALSE,kFALSE},
950  {kFALSE,kFALSE,kFALSE,kFALSE}
951  };
952  Bool_t deadSPDLay2PbPb2011[40][4]={
953  {kTRUE,kTRUE,kTRUE,kTRUE},
954  {kTRUE,kTRUE,kTRUE,kTRUE},
955  {kTRUE,kTRUE,kTRUE,kTRUE},
956  {kTRUE,kTRUE,kTRUE,kTRUE},
957  {kTRUE,kTRUE,kTRUE,kTRUE},
958  {kTRUE,kTRUE,kTRUE,kTRUE},
959  {kTRUE,kTRUE,kTRUE,kTRUE},
960  {kTRUE,kTRUE,kTRUE,kTRUE},
961  {kTRUE,kTRUE,kTRUE,kTRUE},
962  {kTRUE,kTRUE,kTRUE,kTRUE},
963  {kTRUE,kTRUE,kTRUE,kTRUE},
964  {kTRUE,kTRUE,kTRUE,kTRUE},
965  {kFALSE,kFALSE,kFALSE,kFALSE},
966  {kFALSE,kFALSE,kTRUE,kTRUE},
967  {kTRUE,kTRUE,kTRUE,kTRUE},
968  {kTRUE,kTRUE,kTRUE,kTRUE},
969  {kTRUE,kTRUE,kFALSE,kFALSE},
970  {kTRUE,kTRUE,kTRUE,kTRUE},
971  {kTRUE,kTRUE,kTRUE,kTRUE},
972  {kTRUE,kTRUE,kTRUE,kTRUE},
973  {kFALSE,kFALSE,kFALSE,kFALSE},
974  {kFALSE,kFALSE,kFALSE,kFALSE},
975  {kTRUE,kTRUE,kTRUE,kTRUE},
976  {kTRUE,kTRUE,kTRUE,kTRUE},
977  {kFALSE,kFALSE,kFALSE,kFALSE},
978  {kFALSE,kFALSE,kFALSE,kFALSE},
979  {kTRUE,kTRUE,kTRUE,kTRUE},
980  {kTRUE,kTRUE,kTRUE,kTRUE},
981  {kFALSE,kFALSE,kFALSE,kFALSE},
982  {kFALSE,kFALSE,kFALSE,kFALSE},
983  {kFALSE,kFALSE,kFALSE,kFALSE},
984  {kFALSE,kFALSE,kFALSE,kFALSE},
985  {kTRUE,kTRUE,kTRUE,kTRUE},
986  {kTRUE,kTRUE,kTRUE,kTRUE},
987  {kTRUE,kTRUE,kFALSE,kFALSE},
988  {kTRUE,kTRUE,kTRUE,kTRUE},
989  {kFALSE,kFALSE,kFALSE,kFALSE},
990  {kFALSE,kFALSE,kFALSE,kFALSE},
991  {kFALSE,kFALSE,kFALSE,kFALSE},
992  {kFALSE,kFALSE,kFALSE,kFALSE}
993  };
994  Double_t xyz1[3],xyz2[3];
995  esdTrack.GetXYZAt(3.9,0.,xyz1);
996  esdTrack.GetXYZAt(7.6,0.,xyz2);
997  Double_t phi1=TMath::ATan2(xyz1[1],xyz1[0]);
998  if(phi1<0) phi1+=2*TMath::Pi();
999  Int_t lad1=(Int_t)(phi1/(2.*TMath::Pi()/20.));
1000  Double_t phi2=TMath::ATan2(xyz2[1],xyz2[0]);
1001  if(phi2<0) phi2+=2*TMath::Pi();
1002  Int_t lad2=(Int_t)(phi2/(2.*TMath::Pi()/40.));
1003  Int_t mod1=TMath::Floor((xyz1[2]+14)/7.);
1004  Int_t mod2=TMath::Floor((xyz2[2]+14)/7.);
1005  Bool_t lay1ok=kFALSE;
1006  if(mod1>=0 && mod1<4 && lad1<20){
1007  lay1ok=deadSPDLay1PbPb2011[lad1][mod1];
1008  }
1009  Bool_t lay2ok=kFALSE;
1010  if(mod2>=0 && mod2<4 && lad2<40){
1011  lay2ok=deadSPDLay2PbPb2011[lad2][mod2];
1012  }
1013  if(!lay1ok && !lay2ok) return kFALSE;
1014  }
1015 
1017  // Cut tracks crossing the SPD at 5.6<phi<2pi
1018  Double_t xyz1[3],xyz2[3];
1019  esdTrack.GetXYZAt(3.9,0.,xyz1);
1020  esdTrack.GetXYZAt(7.6,0.,xyz2);
1021  Double_t phi1=TMath::ATan2(xyz1[1],xyz1[0]);
1022  if(phi1<0) phi1+=2*TMath::Pi();
1023  Double_t phi2=TMath::ATan2(xyz2[1],xyz2[0]);
1024  if(phi2<0) phi2+=2*TMath::Pi();
1025  Bool_t lay1ok=kTRUE;
1026  if(phi1>5.6 && phi1<2.*TMath::Pi()) lay1ok=kFALSE;
1027  Bool_t lay2ok=kTRUE;
1028  if(phi2>5.6 && phi2<2.*TMath::Pi()) lay2ok=kFALSE;
1029  if(!lay1ok || !lay2ok) return kFALSE;
1030  }
1031 
1032  return kTRUE;
1033 }
1034 //---------------------------------------------------------------------------
1035 void AliRDHFCuts::SetPtBins(Int_t nPtBinLimits,Float_t *ptBinLimits) {
1036  // Set the pt bins
1037 
1038  if(fPtBinLimits) {
1039  delete [] fPtBinLimits;
1040  fPtBinLimits = NULL;
1041  printf("Changing the pt bins\n");
1042  }
1043 
1044  if(nPtBinLimits != fnPtBins+1){
1045  cout<<"Warning: ptBinLimits dimention "<<nPtBinLimits<<" != nPtBins+1 ("<<fnPtBins+1<<")\nSetting nPtBins to "<<nPtBinLimits-1<<endl;
1046  SetNPtBins(nPtBinLimits-1);
1047  }
1048 
1049  fnPtBinLimits = nPtBinLimits;
1050  SetGlobalIndex();
1051  //cout<<"Changing also Global Index -> "<<fGlobalIndex<<endl;
1053  for(Int_t ib=0; ib<nPtBinLimits; ib++) fPtBinLimits[ib]=ptBinLimits[ib];
1054 
1055  return;
1056 }
1057 //---------------------------------------------------------------------------
1058 void AliRDHFCuts::SetVarNames(Int_t nVars,TString *varNames,Bool_t *isUpperCut){
1059  // Set the variable names
1060 
1061  if(fVarNames) {
1062  delete [] fVarNames;
1063  fVarNames = NULL;
1064  //printf("Changing the variable names\n");
1065  }
1066  if(nVars!=fnVars){
1067  printf("Wrong number of variables: it has to be %d\n",fnVars);
1068  return;
1069  }
1070  //fnVars=nVars;
1071  fVarNames = new TString[nVars];
1072  fIsUpperCut = new Bool_t[nVars];
1073  for(Int_t iv=0; iv<nVars; iv++) {
1074  fVarNames[iv] = varNames[iv];
1075  fIsUpperCut[iv] = isUpperCut[iv];
1076  }
1077 
1078  return;
1079 }
1080 //---------------------------------------------------------------------------
1082  // Set the variables to be used for cuts optimization
1083 
1084  if(fVarsForOpt) {
1085  delete [] fVarsForOpt;
1086  fVarsForOpt = NULL;
1087  //printf("Changing the variables for cut optimization\n");
1088  }
1089 
1090  if(nVars==0){
1091  printf("%d not accepted as number of variables: it has to be %d\n",nVars,fnVars);
1092  return;
1093  }
1094 
1095  fnVarsForOpt = 0;
1096  fVarsForOpt = new Bool_t[fnVars];
1097  for(Int_t iv=0; iv<fnVars; iv++) {
1098  fVarsForOpt[iv]=forOpt[iv];
1099  if(fVarsForOpt[iv]) fnVarsForOpt++;
1100  }
1101 
1102  return;
1103 }
1104 
1105 //---------------------------------------------------------------------------
1107  //
1108  // set centrality estimator
1109  //
1110  fUseCentrality=flag;
1111  if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid) AliWarning("Centrality estimator not valid");
1112 
1113  return;
1114 }
1115 
1116 
1117 //---------------------------------------------------------------------------
1119  //
1120  // store the cuts
1121  //
1122  if(nVars!=fnVars) {
1123  printf("Wrong number of variables: it has to be %d\n",fnVars);
1124  AliFatal("exiting");
1125  }
1126  if(nPtBins!=fnPtBins) {
1127  printf("Wrong number of pt bins: it has to be %d\n",fnPtBins);
1128  AliFatal("exiting");
1129  }
1130 
1131  if(!fCutsRD) fCutsRD = new Float_t[fGlobalIndex];
1132 
1133 
1134  for(Int_t iv=0; iv<fnVars; iv++) {
1135 
1136  for(Int_t ib=0; ib<fnPtBins; ib++) {
1137 
1138  //check
1139  if(GetGlobalIndex(iv,ib)>=fGlobalIndex) {
1140  cout<<"Overflow, exit..."<<endl;
1141  return;
1142  }
1143 
1144  fCutsRD[GetGlobalIndex(iv,ib)] = cutsRD[iv][ib];
1145 
1146  }
1147  }
1148  return;
1149 }
1150 //---------------------------------------------------------------------------
1151 void AliRDHFCuts::SetCuts(Int_t glIndex,Float_t* cutsRDGlob){
1152  //
1153  // store the cuts
1154  //
1155  if(glIndex != fGlobalIndex){
1156  cout<<"Wrong array size: it has to be "<<fGlobalIndex<<endl;
1157  AliFatal("exiting");
1158  }
1159  if(!fCutsRD) fCutsRD = new Float_t[fGlobalIndex];
1160 
1161  for(Int_t iGl=0;iGl<fGlobalIndex;iGl++){
1162  fCutsRD[iGl] = cutsRDGlob[iGl];
1163  }
1164  return;
1165 }
1166 //---------------------------------------------------------------------------
1168  //
1169  // print all cuts values
1170  //
1171 
1172  printf("Minimum vtx type %d\n",fMinVtxType);
1173  printf("Minimum vtx contr %d\n",fMinVtxContr);
1174  printf("Max vtx red chi2 %f\n",fMaxVtxRedChi2);
1175  printf("Min SPD mult %d\n",fMinSPDMultiplicity);
1176  printf("Use PID %d OldPid=%d\n",(Int_t)fUsePID,fPidHF ? fPidHF->GetOldPid() : -1);
1177  printf("Remove daughters from vtx %d\n",(Int_t)fRemoveDaughtersFromPrimary);
1178  printf("Physics selection: %s\n",fUsePhysicsSelection ? "Yes" : "No");
1179  printf("Pileup rejection: %s\n",(fOptPileup > 0) ? "Yes" : "No");
1180  if(fOptPileup==1) printf(" -- Reject pileup event");
1181  if(fOptPileup==2) printf(" -- Reject tracks from pileup vtx");
1182  if(fUseCentrality>0) {
1183  TString estimator="";
1184  if(fUseCentrality==kCentV0M) estimator = "V0";
1185  if(fUseCentrality==kCentTRK) estimator = "Tracks";
1186  if(fUseCentrality==kCentTKL) estimator = "Tracklets";
1187  if(fUseCentrality==kCentCL1) estimator = "SPD clusters outer";
1188  if(fUseCentrality==kCentZNA) estimator = "ZNA";
1189  if(fUseCentrality==kCentZPA) estimator = "ZPA";
1190  if(fUseCentrality==kCentV0A) estimator = "V0A";
1191  if(fUseCentrality==kCentCL0) estimator = "SPD clusters inner";
1192  printf("Centrality class considered: %.1f-%.1f, estimated with %s\n",fMinCentrality,fMaxCentrality,estimator.Data());
1193  }
1194  if(fIsCandTrackSPDFirst) printf("Check for candidates with pt < %2.2f, that daughters fullfill kFirst criteria\n",fMaxPtCandTrackSPDFirst);
1195 
1196  if(fCutRatioClsOverCrossRowsTPC) printf("N TPC Clusters > %f N TPC Crossed Rows\n", fCutRatioClsOverCrossRowsTPC);
1197  if(fCutRatioSignalNOverCrossRowsTPC) printf("N TPC Points for dE/dx > %f N TPC Crossed Rows\n", fCutRatioSignalNOverCrossRowsTPC);
1198  if(f1CutMinNCrossedRowsTPCPtDep) printf("N TPC Crossed Rows pT-dependent cut: %s\n", fCutMinCrossedRowsTPCPtDep.Data());
1199 
1200  if(fVarNames){
1201  cout<<"Array of variables"<<endl;
1202  for(Int_t iv=0;iv<fnVars;iv++){
1203  cout<<fVarNames[iv]<<"\t";
1204  }
1205  cout<<endl;
1206  }
1207  if(fVarsForOpt){
1208  cout<<"Array of optimization"<<endl;
1209  for(Int_t iv=0;iv<fnVars;iv++){
1210  cout<<fVarsForOpt[iv]<<"\t";
1211  }
1212  cout<<endl;
1213  }
1214  if(fIsUpperCut){
1215  cout<<"Array of upper/lower cut"<<endl;
1216  for(Int_t iv=0;iv<fnVars;iv++){
1217  cout<<fIsUpperCut[iv]<<"\t";
1218  }
1219  cout<<endl;
1220  }
1221  if(fPtBinLimits){
1222  cout<<"Array of ptbin limits"<<endl;
1223  for(Int_t ib=0;ib<fnPtBinLimits;ib++){
1224  cout<<fPtBinLimits[ib]<<"\t";
1225  }
1226  cout<<endl;
1227  }
1228  if(fCutsRD){
1229  cout<<"Matrix of cuts"<<endl;
1230  for(Int_t iv=0;iv<fnVars;iv++){
1231  for(Int_t ib=0;ib<fnPtBins;ib++){
1232  cout<<"fCutsRD["<<iv<<"]["<<ib<<"] = "<<fCutsRD[GetGlobalIndex(iv,ib)]<<"\t";
1233  }
1234  cout<<endl;
1235  }
1236  cout<<endl;
1237  }
1238  printf("fUsePreselect=%d \n",fUsePreselect);
1239  if(fPidHF) fPidHF->PrintAll();
1240  return;
1241 }
1242 
1243 //--------------------------------------------------------------------------
1245  // print the trigger selection
1246 
1247  printf("Selected trigger classes: %s %s\n",fTriggerClass[0].Data(),fTriggerClass[1].Data());
1248 
1249  cout<<" Trigger selection pattern: ";
1250  if( fTriggerMask & AliVEvent::kAny ) cout<<" kAny ";
1251  if( fTriggerMask & AliVEvent::kAnyINT ) cout<<" kAnyINT ";
1252  if( fTriggerMask & AliVEvent::kINT7 ) cout<<" kINT7 ";
1253  if( fTriggerMask & AliVEvent::kMB ) cout<<" kMB ";
1254  if( fTriggerMask & AliVEvent::kCINT5 ) cout<<" kCINT5 ";
1255  if( fTriggerMask & AliVEvent::kCentral ) cout<<" kCentral ";
1256  if( fTriggerMask & AliVEvent::kSemiCentral ) cout<<" kSemiCentral ";
1257  if( fTriggerMask & AliVEvent::kEMCEGA ) cout<<" kEMCEGA ";
1258  if( fTriggerMask & AliVEvent::kHighMult ) cout<<" kHighMult ";
1259  if( fTriggerMask & AliVEvent::kFastOnly ) cout<<" kFastOnly ";
1260  cout << endl<< endl;
1261 
1262 }
1263 
1264 //---------------------------------------------------------------------------
1265 void AliRDHFCuts::GetCuts(Float_t**& cutsRD) const{
1266  //
1267  // get the cuts
1268  //
1269 
1270  //cout<<"Give back a "<<fnVars<<"x"<<fnPtBins<<" matrix."<<endl;
1271 
1272 
1273  Int_t iv,ib;
1274  if(!cutsRD) {
1275  //cout<<"Initialization..."<<endl;
1276  cutsRD=new Float_t*[fnVars];
1277  for(iv=0; iv<fnVars; iv++) {
1278  cutsRD[iv] = new Float_t[fnPtBins];
1279  }
1280  }
1281 
1282  for(Int_t iGlobal=0; iGlobal<fGlobalIndex; iGlobal++) {
1283  GetVarPtIndex(iGlobal,iv,ib);
1284  cutsRD[iv][ib] = fCutsRD[iGlobal];
1285  }
1286 
1287  return;
1288 }
1289 
1290 //---------------------------------------------------------------------------
1292  //
1293  // give the global index from variable and pt bin
1294  //
1295  return iPtBin*fnVars+iVar;
1296 }
1297 
1298 //---------------------------------------------------------------------------
1299 void AliRDHFCuts::GetVarPtIndex(Int_t iGlob, Int_t& iVar, Int_t& iPtBin) const {
1300  //
1301  //give the index of the variable and of the pt bin from the global index
1302  //
1303  iPtBin=(Int_t)iGlob/fnVars;
1304  iVar=iGlob%fnVars;
1305 
1306  return;
1307 }
1308 
1309 //---------------------------------------------------------------------------
1311  //
1312  //give the pt bin where the pt lies.
1313  //
1314  Int_t ptbin=-1;
1315  if(pt<fPtBinLimits[0])return ptbin;
1316  for (Int_t i=0;i<fnPtBins;i++){
1317  if(pt<fPtBinLimits[i+1]) {
1318  ptbin=i;
1319  break;
1320  }
1321  }
1322  return ptbin;
1323 }
1324 //-------------------------------------------------------------------
1326  //
1327  // Give the value of cut set for the variable iVar and the pt bin iPtBin
1328  //
1329  if(!fCutsRD){
1330  cout<<"Cuts not iniziaisez yet"<<endl;
1331  return 0;
1332  }
1333  return fCutsRD[GetGlobalIndex(iVar,iPtBin)];
1334 }
1335 
1336 //-------------------------------------------------------------------
1338 
1339  if(aodEvent->GetRunNumber()<244824)return GetCentralityOldFramework(aodEvent,estimator);
1340  Double_t cent=-999;
1341 
1342  if(estimator==kCentOff) return -999;
1343 
1344  AliMultSelection *multSelection = (AliMultSelection*)aodEvent->FindListObject(fMultSelectionObjectName);
1345  if(!multSelection){
1346  AliWarning("AliMultSelection could not be found in the aod event list of objects");
1347  return cent;
1348  }
1349 
1350  // Compare centrality with the centrality at AOD filtering and on-the-fly
1351  if( fMultSelectionObjectName.CompareTo("MultSelection")!=0 ){
1352  AliMultSelection *defaultmultSelection = (AliMultSelection*)aodEvent->FindListObject("MultSelection");
1353  if(!defaultmultSelection){
1354  AliWarning("AliMultSelection default method could not be found in the aod event list of objects");
1355  return cent;
1356  }
1357  Float_t defaultCent = defaultmultSelection->GetMultiplicityPercentile("V0M");
1358  Float_t newCent = multSelection->GetMultiplicityPercentile("V0M");
1359  if( defaultCent<20. && newCent>20.) fEvRejectionBits+=1<<kMismatchOldNewCentrality;
1360  else if (defaultCent>20. && newCent<20.) fEvRejectionBits+=1<<kMismatchOldNewCentrality;
1361  }
1362 
1363  if(estimator==kCentV0M){
1364  cent=multSelection->GetMultiplicityPercentile("V0M");
1365  }else if(estimator==kCentV0A){
1366  cent=multSelection->GetMultiplicityPercentile("V0A");
1367  }else if(estimator==kCentZNA){
1368  cent=multSelection->GetMultiplicityPercentile("ZNA");
1369  }else if(estimator==kCentCL1){
1370  cent=multSelection->GetMultiplicityPercentile("CL1");
1371  }else if(estimator==kCentCL0){
1372  cent=multSelection->GetMultiplicityPercentile("CL0");
1373  }else {
1374  AliWarning(Form("CENTRALITY ESTIMATE WITH ESTIMATOR %d NOT YET IMPLEMENTED FOR NEW FRAMEWORK",(Int_t)estimator));
1375  return cent;
1376  }
1377  Int_t qual = multSelection->GetEvSelCode();
1378  if(qual == 199 ) cent=-999;
1379  return cent;
1380 }
1381 //-------------------------------------------------------------------
1383  //
1384  // Get centrality percentile
1385  //
1386 
1387  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)aodEvent)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
1388  if(mcArray) {fUseAOD049=kFALSE;}
1389 
1390  AliAODHeader *header=dynamic_cast<AliAODHeader*>(aodEvent->GetHeader());
1391  if(!header) AliFatal("Not a standard AOD");
1392  AliCentrality *centrality=header->GetCentralityP();
1393  Float_t cent=-999.;
1394  Bool_t isSelRun=kFALSE;
1395  Int_t selRun[5]={138364, 138826, 138828, 138836, 138871};
1396  if(!centrality) return cent;
1397  else{
1398  if (estimator==kCentV0M){
1399  cent=(Float_t)(centrality->GetCentralityPercentile("V0M"));
1400  if(cent<0){
1401  Int_t quality = centrality->GetQuality();
1402  if(quality<=1){ // fQuality==1 means rejected by zVertex cut that we apply a part and we want to keep separate (Giacomo)
1403  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");
1404  }else{
1405  Int_t runnum=aodEvent->GetRunNumber();
1406  for(Int_t ir=0;ir<5;ir++){
1407  if(runnum==selRun[ir]){
1408  isSelRun=kTRUE;
1409  break;
1410  }
1411  }
1412  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");
1413  }
1414  }
1415 
1416  //temporary fix for AOD049 outliers
1417  if(fUseAOD049&&cent>=0){
1418  Float_t v0=0;
1419  AliAODVZERO* aodV0 = aodEvent->GetVZEROData();
1420  v0+=aodV0->GetMTotV0A();
1421  v0+=aodV0->GetMTotV0C();
1422  if(cent==0&&v0<19500)return -1;//filtering issue
1423  Float_t tkl = (Float_t)(aodEvent->GetTracklets()->GetNumberOfTracklets());
1424  Float_t val= 1.30552 + 0.147931 * v0;
1425  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};
1426  if ( TMath::Abs(tkl-val) > 6.*tklSigma[(Int_t)cent] )return -1;//outlier
1427  }
1428  }
1429  else {
1430  if (estimator==kCentTRK) {
1431  cent=(Float_t)(centrality->GetCentralityPercentile("TRK"));
1432  if(cent<0){
1433  Int_t quality = centrality->GetQuality();
1434  if(quality<=1){
1435  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TRK");
1436  }else{
1437  Int_t runnum=aodEvent->GetRunNumber();
1438  for(Int_t ir=0;ir<5;ir++){
1439  if(runnum==selRun[ir]){
1440  isSelRun=kTRUE;
1441  break;
1442  }
1443  }
1444  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TRK");
1445  }
1446  }
1447  }
1448  else{
1449  if (estimator==kCentTKL){
1450  cent=(Float_t)(centrality->GetCentralityPercentile("TKL"));
1451  if(cent<0){
1452  Int_t quality = centrality->GetQuality();
1453  if(quality<=1){
1454  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TKL");
1455  }else{
1456  Int_t runnum=aodEvent->GetRunNumber();
1457  for(Int_t ir=0;ir<5;ir++){
1458  if(runnum==selRun[ir]){
1459  isSelRun=kTRUE;
1460  break;
1461  }
1462  }
1463  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TKL");
1464  }
1465  }
1466  }
1467  else{
1468  if (estimator==kCentCL1){
1469  cent=(Float_t)(centrality->GetCentralityPercentile("CL1"));
1470  if(cent<0){
1471  Int_t quality = centrality->GetQuality();
1472  if(quality<=1){
1473  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("CL1");
1474  }else{
1475  Int_t runnum=aodEvent->GetRunNumber();
1476  for(Int_t ir=0;ir<5;ir++){
1477  if(runnum==selRun[ir]){
1478  isSelRun=kTRUE;
1479  break;
1480  }
1481  }
1482  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("CL1");
1483  }
1484  }
1485  }
1486  else{
1487  if (estimator==kCentZNA){
1488  cent=(Float_t)(centrality->GetCentralityPercentile("ZNA"));
1489  if(cent<0){
1490  Int_t quality = centrality->GetQuality();
1491  if(quality<=1){
1492  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("ZNA");
1493  }else{
1494  Int_t runnum=aodEvent->GetRunNumber();
1495  for(Int_t ir=0;ir<5;ir++){
1496  if(runnum==selRun[ir]){
1497  isSelRun=kTRUE;
1498  break;
1499  }
1500  }
1501  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("ZNA");
1502  }
1503  }
1504  }
1505  else{
1506  if (estimator==kCentZPA){
1507  cent=(Float_t)(centrality->GetCentralityPercentile("ZPA"));
1508  if(cent<0){
1509  Int_t quality = centrality->GetQuality();
1510  if(quality<=1){
1511  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("ZPA");
1512  }else{
1513  Int_t runnum=aodEvent->GetRunNumber();
1514  for(Int_t ir=0;ir<5;ir++){
1515  if(runnum==selRun[ir]){
1516  isSelRun=kTRUE;
1517  break;
1518  }
1519  }
1520  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("ZPA");
1521  }
1522  }
1523  }
1524  else{
1525  if (estimator==kCentV0A){
1526  cent=(Float_t)(centrality->GetCentralityPercentile("V0A"));
1527  if(cent<0){
1528  Int_t quality = centrality->GetQuality();
1529  if(quality<=1){
1530  cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0A");
1531  }else{
1532  Int_t runnum=aodEvent->GetRunNumber();
1533  for(Int_t ir=0;ir<5;ir++){
1534  if(runnum==selRun[ir]){
1535  isSelRun=kTRUE;
1536  break;
1537  }
1538  }
1539  if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0A");
1540  }
1541  }
1542  }
1543  else {
1544  AliWarning("Centrality estimator not valid");
1545 
1546  }
1547  }
1548  }
1549  }
1550  }
1551  }
1552  }
1553  }
1554  return cent;
1555 }
1556 //-------------------------------------------------------------------
1558  //
1559  // Compare two cuts objects
1560  //
1561 
1562  Bool_t areEqual=kTRUE;
1563 
1564  if(fMinVtxType!=obj->fMinVtxType) { printf("Minimum vtx type %d %d\n",fMinVtxType,obj->fMinVtxType); areEqual=kFALSE;}
1565 
1566  if(fMinVtxContr!=obj->fMinVtxContr) { printf("Minimum vtx contr %d %d\n",fMinVtxContr,obj->fMinVtxContr); areEqual=kFALSE;}
1567 
1568  if(TMath::Abs(fMaxVtxRedChi2-obj->fMaxVtxRedChi2)>1.e-10) { printf("Max vtx red chi2 %f %f\n",fMaxVtxRedChi2,obj->fMaxVtxRedChi2);areEqual=kFALSE;}
1569 
1570  if(fMinSPDMultiplicity!=obj->fMinSPDMultiplicity) { printf("Min SPD mult %d\n %d",fMinSPDMultiplicity,obj->fMinSPDMultiplicity);areEqual=kFALSE;}
1571 
1572  if(fUsePID!=obj->fUsePID) { printf("Use PID %d %d\n",(Int_t)fUsePID,(Int_t)obj->fUsePID); areEqual=kFALSE;}
1573 
1574  if(fRemoveDaughtersFromPrimary!=obj->fRemoveDaughtersFromPrimary) {printf("Remove daughters from vtx %d %d\n",(Int_t)fRemoveDaughtersFromPrimary,(Int_t)obj->fRemoveDaughtersFromPrimary); areEqual=kFALSE;}
1575  if(fTrackCuts){
1576  if(fTrackCuts->GetMinNClusterTPC()!=obj->fTrackCuts->GetMinNClusterTPC()) {printf("MinNClsTPC %d %d\n",fTrackCuts->GetMinNClusterTPC(),obj->fTrackCuts->GetMinNClusterTPC()); areEqual=kFALSE;}
1577 
1578  if(fTrackCuts->GetMinNClustersITS()!=obj->fTrackCuts->GetMinNClustersITS()) {printf("MinNClsITS %d %d\n",fTrackCuts->GetMinNClustersITS(),obj->fTrackCuts->GetMinNClustersITS()); areEqual=kFALSE;}
1579 
1580  if(TMath::Abs(fTrackCuts->GetMaxChi2PerClusterTPC()-obj->fTrackCuts->GetMaxChi2PerClusterTPC())>1.e-10) {printf("MaxChi2ClsTPC %f %f\n",fTrackCuts->GetMaxChi2PerClusterTPC(),obj->fTrackCuts->GetMaxChi2PerClusterTPC()); areEqual=kFALSE;}
1581 
1582  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;}
1583  }
1584 
1585  if(fUsePreselect!=obj->fUsePreselect){printf("fUsePreselect: %d %d\n",fUsePreselect,obj->fUsePreselect);areEqual=kFALSE;}
1586 
1587  if(fCutsRD) {
1588  for(Int_t iv=0;iv<fnVars;iv++) {
1589  for(Int_t ib=0;ib<fnPtBins;ib++) {
1590  if(TMath::Abs(fCutsRD[GetGlobalIndex(iv,ib)]-obj->fCutsRD[GetGlobalIndex(iv,ib)])>1.e-10) {
1591  cout<<"fCutsRD["<<iv<<"]["<<ib<<"] = "<<fCutsRD[GetGlobalIndex(iv,ib)]<<" "<<obj->fCutsRD[GetGlobalIndex(iv,ib)]<<"\n";
1592  areEqual=kFALSE;
1593  }
1594  }
1595  }
1596  }
1597 
1598  return areEqual;
1599 }
1600 //---------------------------------------------------------------------------
1602  //
1603  // print cuts values in table format
1604  //
1605 
1606  TString ptString = "pT range";
1607  if(fVarNames && fPtBinLimits && fCutsRD){
1608  TString firstLine(Form("* %-15s",ptString.Data()));
1609  for (Int_t ivar=0; ivar<fnVars; ivar++){
1610  firstLine+=Form("* %-15s ",fVarNames[ivar].Data());
1611  if (ivar == fnVars){
1612  firstLine+="*\n";
1613  }
1614  }
1615  Printf("%s",firstLine.Data());
1616 
1617  for (Int_t ipt=0; ipt<fnPtBins; ipt++){
1618  TString line;
1619  if (ipt==fnPtBins-1){
1620  line=Form("* %5.1f < pt < inf ",fPtBinLimits[ipt]);
1621  }
1622  else{
1623  line=Form("* %5.1f < pt < %4.1f ",fPtBinLimits[ipt],fPtBinLimits[ipt+1]);
1624  }
1625  for (Int_t ivar=0; ivar<fnVars; ivar++){
1626  line+=Form("* %-15f ",fCutsRD[GetGlobalIndex(ivar,ipt)]);
1627  }
1628  Printf("%s",line.Data());
1629  }
1630 
1631  }
1632 
1633 
1634  return;
1635 }
1636 //--------------------------------------------------------------------------
1638  AliAODEvent *aod) const
1639 {
1640  //
1641  // Recalculate primary vertex without daughters
1642  //
1643 
1644  if(!aod) {
1645  AliError("Can not remove daughters from vertex without AOD event");
1646  return 0;
1647  }
1648 
1649  AliAODVertex *recvtx=d->RemoveDaughtersFromPrimaryVtx(aod);
1650  if(!recvtx){
1651  AliDebug(2,"Removal of daughter tracks failed");
1652  return kFALSE;
1653  }
1654 
1655 
1656  //set recalculed primary vertex
1657  d->SetOwnPrimaryVtx(recvtx);
1658  delete recvtx;
1659 
1660  return kTRUE;
1661 }
1662 //--------------------------------------------------------------------------
1664 {
1665  //
1666  // Recalculate primary vertex without daughters
1667  //
1668 
1669  if(!aod) {
1670  AliError("Can not get MC vertex without AOD event");
1671  return kFALSE;
1672  }
1673 
1674  // load MC header
1675  AliAODMCHeader *mcHeader =
1676  (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
1677  if(!mcHeader) {
1678  AliError("Can not get MC vertex without AODMCHeader event");
1679  return kFALSE;
1680  }
1681  Double_t pos[3];
1682  Double_t covmatrix[6]={0.,0.,0.,0.,0.,0.};
1683  mcHeader->GetVertex(pos);
1684  AliAODVertex *recvtx=new AliAODVertex(pos,covmatrix);
1685 
1686  if(!recvtx){
1687  AliDebug(2,"Removal of daughter tracks failed");
1688  return kFALSE;
1689  }
1690 
1691  //set recalculed primary vertex
1692  d->SetOwnPrimaryVtx(recvtx);
1693 
1694  d->RecalculateImpPars(recvtx,aod);
1695 
1696  delete recvtx;
1697 
1698  return kTRUE;
1699 }
1700 //--------------------------------------------------------------------------
1702  AliAODEvent *aod,
1703  AliAODVertex *origownvtx) const
1704 {
1705  //
1706  // Clean-up own primary vertex if needed
1707  //
1708 
1710  d->UnsetOwnPrimaryVtx();
1711  if(origownvtx) {
1712  d->SetOwnPrimaryVtx(origownvtx);
1713  delete origownvtx; origownvtx=NULL;
1714  }
1715  d->RecalculateImpPars(d->GetPrimaryVtx(),aod);
1716  } else {
1717  if(origownvtx) {
1718  delete origownvtx; origownvtx=NULL;
1719  }
1720  }
1721  return;
1722 }
1723 //--------------------------------------------------------------------------
1724 Bool_t AliRDHFCuts::IsSignalMC(AliAODRecoDecay *d,AliAODEvent *aod,Int_t pdg) const
1725 {
1726  //
1727  // Checks if this candidate is matched to MC signal
1728  //
1729 
1730  if(!aod) return kFALSE;
1731 
1732  // get the MC array
1733  TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)aod)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
1734 
1735  if(!mcArray) return kFALSE;
1736 
1737  // try to match
1738  Int_t label = d->MatchToMC(pdg,mcArray);
1739 
1740  if(label>=0) {
1741  //printf("MATCH!\n");
1742  return kTRUE;
1743  }
1744 
1745  return kFALSE;
1746 }
1747 
1748 
1749 //--------------------------------------------------------------------------
1751  // recompute event primary vertex from AOD tracks
1752 
1753  AliVertexerTracks *vertexer = new AliVertexerTracks(event->GetMagneticField());
1754  vertexer->SetITSMode();
1755  vertexer->SetMinClusters(3);
1756 
1757  AliAODVertex* pvtx=event->GetPrimaryVertex();
1758  if(strstr(pvtx->GetTitle(),"VertexerTracksWithConstraint")) {
1759  Float_t diamondcovxy[3];
1760  event->GetDiamondCovXY(diamondcovxy);
1761  Double_t pos[3]={event->GetDiamondX(),event->GetDiamondY(),0.};
1762  Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.*10.};
1763  AliESDVertex *diamond = new AliESDVertex(pos,cov,1.,1);
1764  vertexer->SetVtxStart(diamond);
1765  delete diamond; diamond=NULL;
1766  }
1767 
1768  AliESDVertex* vertexESD = (AliESDVertex*)vertexer->FindPrimaryVertex(event);
1769  if(!vertexESD) return kFALSE;
1770  if(vertexESD->GetNContributors()<=0) {
1771  //AliDebug(2,"vertexing failed");
1772  delete vertexESD; vertexESD=NULL;
1773  return kFALSE;
1774  }
1775  delete vertexer; vertexer=NULL;
1776 
1777  // convert to AliAODVertex
1778  Double_t pos[3],cov[6],chi2perNDF;
1779  vertexESD->GetXYZ(pos); // position
1780  vertexESD->GetCovMatrix(cov); //covariance matrix
1781  chi2perNDF = vertexESD->GetChi2toNDF();
1782  delete vertexESD; vertexESD=NULL;
1783 
1784  pvtx->SetPosition(pos[0],pos[1],pos[2]);
1785  pvtx->SetChi2perNDF(chi2perNDF);
1786  pvtx->SetCovMatrix(cov);
1787 
1788  return kTRUE;
1789 }
void SetNPtBins(Int_t nptBins)
Definition: AliRDHFCuts.h:406
Int_t fIsSelectedCuts
fix the daughter track references
Definition: AliRDHFCuts.h:458
Int_t pdg
Bool_t IsEventSelectedForCentrFlattening(Float_t centvalue)
Float_t fMinDzPileup
min. n. of tracklets in pileup vertex
Definition: AliRDHFCuts.h:447
Double_t fMaxPtCandTrackSPDFirst
flag to select the track kFirst criteria for pt < ptlimit
Definition: AliRDHFCuts.h:465
Int_t fWhyRejection
PID for heavy flavours manager.
Definition: AliRDHFCuts.h:440
double Double_t
Definition: External.C:58
Int_t fCutOnzVertexSPD
flag to apply cut on tracklets vs. centrality for 2011 data
Definition: AliRDHFCuts.h:470
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:483
Bool_t fUseMCVertex
flag to switch on the removal of duaghters from the primary vertex computation
Definition: AliRDHFCuts.h:443
Bool_t fUseOnlyOneTrigger
trigger mask
Definition: AliRDHFCuts.h:422
Int_t IsEventSelectedInCentrality(AliVEvent *event)
Bool_t * fIsUpperCut
Definition: AliRDHFCuts.h:436
Bool_t SetMCPrimaryVtx(AliAODRecoDecayHF *d, AliAODEvent *aod) const
Bool_t fRemoveDaughtersFromPrimary
Definition: AliRDHFCuts.h:442
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:448
void SetUseCentrality(Int_t flag=1)
Float_t * fPtBinLimits
"number of limits", that is fnPtBins+1
Definition: AliRDHFCuts.h:429
Int_t fMinVtxContr
0: not cut; 1: SPDZ; 2: SPD3D; 3: Tracks
Definition: AliRDHFCuts.h:414
centrality
char Char_t
Definition: External.C:18
static Int_t CheckMatchingAODdeltaAODevents()
Float_t fMaxVtxRedChi2
minimum vertex contributors
Definition: AliRDHFCuts.h:415
void SetMClowenpp2011(Bool_t mc)
Definition: AliAODPidHF.h:102
const Float_t * GetCuts() const
Definition: AliRDHFCuts.h:255
void SetMinCrossedRowsTPCPtDep(const char *rows="")
Float_t fCutRatioClsOverCrossRowsTPC
histogram with reference centrality distribution for centrality distribution flattening ...
Definition: AliRDHFCuts.h:475
Double_t fDeadZoneWidth
flag for enabling/disabling geometrical cut on TPC track
Definition: AliRDHFCuts.h:480
Int_t fMinVtxType
cuts on the event
Definition: AliRDHFCuts.h:413
void SetppLowEn2011(Bool_t opt)
Definition: AliAODPidHF.h:104
Double_t fMaxRapidityCand
minimum pt of the candidate
Definition: AliRDHFCuts.h:462
Double_t fCutGeoNcrNclFractionNcl
4th parameter of GeoNcrNcl cut
Definition: AliRDHFCuts.h:484
void SetGlobalIndex()
Definition: AliRDHFCuts.h:202
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:449
Bool_t fUsePID
Definition: AliRDHFCuts.h:437
TRandom * gRandom
Bool_t fRejectPlpFromDiffBCMV
minimum weighted distance in Z between 2 vertices (multi-vertexer)
Definition: AliRDHFCuts.h:420
void PrintAll() const
void SetPidHF(AliAODPidHF *pidObj)
see enum below
Definition: AliRDHFCuts.h:216
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:248
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:469
Bool_t fUseV0ANDSelectionOffline
5th parameter of GeoNcrNcl cut
Definition: AliRDHFCuts.h:485
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:421
Bool_t fUseCutGeoNcrNcl
pT-dep cut in TPC minimum n crossed rows
Definition: AliRDHFCuts.h:479
int Int_t
Definition: External.C:63
Bool_t fUseTrackSelectionWithFilterBits
flag to reject kink daughters
Definition: AliRDHFCuts.h:472
void SetCuts(Int_t nVars, Int_t nPtBins, Float_t **cutsRD)
Int_t fnVarsForOpt
Definition: AliRDHFCuts.h:432
Float_t fMinWDzPileupMV
max chi2 per contributor of the pile-up vertex to consider (multi-vertexer).
Definition: AliRDHFCuts.h:419
Double_t fCutGeoNcrNclLength
1st parameter of GeoNcrNcl cut
Definition: AliRDHFCuts.h:481
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:431
Double_t fMaxPtCand
minimum pt of the candidate
Definition: AliRDHFCuts.h:461
AliESDtrackCuts * fTrackCuts
quality cuts on the daughter tracks
Definition: AliRDHFCuts.h:425
void MakeTable() const
Bool_t fFixRefs
name of the AliMultSelection object to be considered
Definition: AliRDHFCuts.h:457
Bool_t fKeepSignalMC
max rapidity of candidate (if !=-999 overrides IsInFiducialAcceptance)
Definition: AliRDHFCuts.h:463
Int_t fOptPileup
use Physics selection criteria
Definition: AliRDHFCuts.h:445
AliESDtrackCuts * GetTrackCuts() const
Definition: AliRDHFCuts.h:263
void PrintTrigger() const
Float_t fMaxVtxZ
maximum chi2/ndf
Definition: AliRDHFCuts.h:416
Float_t GetCentrality(AliAODEvent *aodEvent)
Definition: AliRDHFCuts.h:259
Int_t fIsSelectedPID
outcome of cuts selection
Definition: AliRDHFCuts.h:459
Float_t fMinCentrality
Definition: AliRDHFCuts.h:454
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:441
Double_t fMaxDiffTRKV0Centr
flag to apply cut on tracks crossing SPD misaligned modules for PP2012 data
Definition: AliRDHFCuts.h:468
Float_t GetCentralityOldFramework(AliAODEvent *aodEvent, AliRDHFCuts::ECentrality estimator)
Int_t fnPtBinLimits
number of pt bins for cuts
Definition: AliRDHFCuts.h:428
Float_t * fCutsRD
fnVars*fnPtBins
Definition: AliRDHFCuts.h:435
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)
Int_t fUsePreselect
Definition: AliRDHFCuts.h:488
Bool_t fUseCentrFlatteningInMC
flag to enable/disable the check on filter bits
Definition: AliRDHFCuts.h:473
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:433
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:467
Bool_t IsEventSelected(AliVEvent *event)
Int_t fMinSPDMultiplicity
maximum |z| of primary vertex
Definition: AliRDHFCuts.h:417
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:438
Bool_t fApplySPDDeadPbPb2011
maximum pt of the candidate for which to check if the daughters fulfill kFirst criteria ...
Definition: AliRDHFCuts.h:466
void CleanOwnPrimaryVtx(AliAODRecoDecayHF *d, AliAODEvent *aod, AliAODVertex *origownvtx) const
Float_t fMaxVtxChi2PileupMV
SPD multiplicity.
Definition: AliRDHFCuts.h:418
Bool_t fIsCandTrackSPDFirst
IsSelected returns always kTRUE for MC signal.
Definition: AliRDHFCuts.h:464
TH1F * fHistCentrDistr
flag for enabling/diabling centrality flattening in MC
Definition: AliRDHFCuts.h:474
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:471
TString fMultSelectionObjectName
maximum centrality for selected events
Definition: AliRDHFCuts.h:456
Int_t fMinContrPileup
option for pielup selection
Definition: AliRDHFCuts.h:446
Bool_t fUseTPCtrackCutsOnThisDaughter
flag to apply V0AND selection offline
Definition: AliRDHFCuts.h:486
Bool_t fUsePhysicsSelection
use MC primary vertex
Definition: AliRDHFCuts.h:444
TFormula * f1CutMinNCrossedRowsTPCPtDep
pT-dep cut in TPC minimum n crossed rows
Definition: AliRDHFCuts.h:478
AliAODVertex * GetPrimaryVtx() const
void AddTrackCuts(const AliESDtrackCuts *cuts)
Definition: AliRDHFCuts.h:208
TString fTriggerClass[2]
flag to select one trigger only
Definition: AliRDHFCuts.h:423
bool Bool_t
Definition: External.C:53
Int_t fnPtBins
cuts on the candidate
Definition: AliRDHFCuts.h:427
Double_t fCutGeoNcrNclGeom1Pt
2nd parameter of GeoNcrNcl cut
Definition: AliRDHFCuts.h:482
Bool_t fApplyZcutOnSPDvtx
flag to apply TPC track quality cuts on specific D-meson daughter (used for different strategies for ...
Definition: AliRDHFCuts.h:487
AliAODPidHF * fPidHF
enable AOD049 centrality cleanup
Definition: AliRDHFCuts.h:439
Float_t fCutRatioSignalNOverCrossRowsTPC
min. value ratio NTPCClusters/NTPCCrossedRows, cut if !=0
Definition: AliRDHFCuts.h:476
Int_t fGlobalIndex
Definition: AliRDHFCuts.h:434
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:477
Double_t fMinPtCand
outcome of PID selection
Definition: AliRDHFCuts.h:460
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:455
void SetPidResponse(AliPIDResponse *pidResp)
Definition: AliAODPidHF.h:111
Bool_t RecomputePrimaryVertex(AliAODEvent *event) const
Int_t fnVars
Definition: AliRDHFCuts.h:430