1 /**************************************************************************
2 * Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 /////////////////////////////////////////////////////////////
20 // Base class for cuts on AOD reconstructed heavy-flavour decay
22 // Author: A.Dainese, andrea.dainese@pd.infn.it
23 /////////////////////////////////////////////////////////////
24 #include <Riostream.h>
26 #include "AliVEvent.h"
27 #include "AliESDEvent.h"
28 #include "AliAODEvent.h"
29 #include "AliVVertex.h"
30 #include "AliESDVertex.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 "AliRDHFCuts.h"
42 #include "AliAnalysisManager.h"
43 #include "AliInputEventHandler.h"
44 #include "AliPIDResponse.h"
48 //--------------------------------------------------------------------------
49 AliRDHFCuts::AliRDHFCuts(const Char_t* name, const Char_t* title) :
50 AliAnalysisCuts(name,title),
55 fMinSPDMultiplicity(0),
57 fTriggerClass("CINT1"),
74 fRemoveDaughtersFromPrimary(kFALSE),
89 // Default Constructor
92 //--------------------------------------------------------------------------
93 AliRDHFCuts::AliRDHFCuts(const AliRDHFCuts &source) :
94 AliAnalysisCuts(source),
95 fMinVtxType(source.fMinVtxType),
96 fMinVtxContr(source.fMinVtxContr),
97 fMaxVtxRedChi2(source.fMaxVtxRedChi2),
98 fMaxVtxZ(source.fMaxVtxZ),
99 fMinSPDMultiplicity(source.fMinSPDMultiplicity),
100 fTriggerMask(source.fTriggerMask),
101 fTriggerClass(source.fTriggerClass),
103 fnPtBins(source.fnPtBins),
104 fnPtBinLimits(source.fnPtBinLimits),
106 fnVars(source.fnVars),
108 fnVarsForOpt(source.fnVarsForOpt),
110 fGlobalIndex(source.fGlobalIndex),
113 fUsePID(source.fUsePID),
114 fUseAOD049(source.fUseAOD049),
116 fWhyRejection(source.fWhyRejection),
117 fEvRejectionBits(source.fEvRejectionBits),
118 fRemoveDaughtersFromPrimary(source.fRemoveDaughtersFromPrimary),
119 fUseMCVertex(source.fUseMCVertex),
120 fOptPileup(source.fOptPileup),
121 fMinContrPileup(source.fMinContrPileup),
122 fMinDzPileup(source.fMinDzPileup),
123 fUseCentrality(source.fUseCentrality),
124 fMinCentrality(source.fMinCentrality),
125 fMaxCentrality(source.fMaxCentrality),
126 fFixRefs(source.fFixRefs),
127 fIsSelectedCuts(source.fIsSelectedCuts),
128 fIsSelectedPID(source.fIsSelectedPID),
129 fMinPtCand(source.fMinPtCand),
130 fMaxPtCand(source.fMaxPtCand)
135 cout<<"Copy constructor"<<endl;
136 if(source.GetTrackCuts()) AddTrackCuts(source.GetTrackCuts());
137 if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);
138 if(source.fVarNames) SetVarNames(source.fnVars,source.fVarNames,source.fIsUpperCut);
139 if(source.fCutsRD) SetCuts(source.fGlobalIndex,source.fCutsRD);
140 if(source.fVarsForOpt) SetVarsForOpt(source.fnVarsForOpt,source.fVarsForOpt);
141 if(source.fPidHF) SetPidHF(source.fPidHF);
145 //--------------------------------------------------------------------------
146 AliRDHFCuts &AliRDHFCuts::operator=(const AliRDHFCuts &source)
149 // assignment operator
151 if(&source == this) return *this;
153 AliAnalysisCuts::operator=(source);
155 fMinVtxType=source.fMinVtxType;
156 fMinVtxContr=source.fMinVtxContr;
157 fMaxVtxRedChi2=source.fMaxVtxRedChi2;
158 fMaxVtxZ=source.fMaxVtxZ;
159 fMinSPDMultiplicity=source.fMinSPDMultiplicity;
160 fTriggerMask=source.fTriggerMask;
161 fTriggerClass=source.fTriggerClass;
162 fnPtBins=source.fnPtBins;
163 fnVars=source.fnVars;
164 fGlobalIndex=source.fGlobalIndex;
165 fnVarsForOpt=source.fnVarsForOpt;
166 fUsePID=source.fUsePID;
167 fUseAOD049=source.fUseAOD049;
168 SetPidHF(source.GetPidHF());
169 fWhyRejection=source.fWhyRejection;
170 fEvRejectionBits=source.fEvRejectionBits;
171 fRemoveDaughtersFromPrimary=source.fRemoveDaughtersFromPrimary;
172 fUseMCVertex=source.fUseMCVertex;
173 fOptPileup=source.fOptPileup;
174 fMinContrPileup=source.fMinContrPileup;
175 fMinDzPileup=source.fMinDzPileup;
176 fUseCentrality=source.fUseCentrality;
177 fMinCentrality=source.fMinCentrality;
178 fMaxCentrality=source.fMaxCentrality;
179 fFixRefs=source.fFixRefs;
180 fIsSelectedCuts=source.fIsSelectedCuts;
181 fIsSelectedPID=source.fIsSelectedPID;
182 fMinPtCand=source.fMinPtCand;
183 fMaxPtCand=source.fMaxPtCand;
185 if(source.GetTrackCuts()) AddTrackCuts(source.GetTrackCuts());
186 if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);
187 if(source.fVarNames) SetVarNames(source.fnVars,source.fVarNames,source.fIsUpperCut);
188 if(source.fCutsRD) SetCuts(source.fGlobalIndex,source.fCutsRD);
189 if(source.fVarsForOpt) SetVarsForOpt(source.fnVarsForOpt,source.fVarsForOpt);
194 //--------------------------------------------------------------------------
195 AliRDHFCuts::~AliRDHFCuts() {
197 // Default Destructor
199 if(fPtBinLimits) {delete [] fPtBinLimits; fPtBinLimits=0;}
200 if(fVarNames) {delete [] fVarNames; fVarNames=0;}
201 if(fVarsForOpt) {delete [] fVarsForOpt; fVarsForOpt=0;}
206 if(fIsUpperCut) {delete [] fIsUpperCut; fIsUpperCut=0;}
212 //---------------------------------------------------------------------------
213 Int_t AliRDHFCuts::IsEventSelectedInCentrality(AliVEvent *event) {
215 // Centrality selection
217 if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid){
218 AliWarning("Centrality estimator not valid");
221 Float_t centvalue=GetCentrality((AliAODEvent*)event);
222 if (centvalue<-998.){//-999 if no centralityP
225 if (centvalue<fMinCentrality || centvalue>fMaxCentrality){
232 //---------------------------------------------------------------------------
233 Bool_t AliRDHFCuts::IsEventSelected(AliVEvent *event) {
237 //if(fTriggerMask && event->GetTriggerMask()!=fTriggerMask) return kFALSE;
245 TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)event)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
246 if(mcArray) {isMC=kTRUE;fUseAOD049=kFALSE;}
248 // settings for the TPC dE/dx BB parameterization
249 if(fPidHF && fPidHF->GetOldPid()) {
250 // pp, from LHC10d onwards
251 if((event->GetRunNumber()>121693 && event->GetRunNumber()<136851) ||
252 event->GetRunNumber()>139517) fPidHF->SetOnePad(kTRUE);
254 if(event->GetRunNumber()>=136851 && event->GetRunNumber()<=139517) fPidHF->SetPbPb(kTRUE);
256 if(isMC) fPidHF->SetMC(kTRUE);
257 if(isMC && (event->GetRunNumber()>=146686 && event->GetRunNumber()<=146860))
258 fPidHF->SetMClowenpp2011(kTRUE);
260 else if(fPidHF && !fPidHF->GetOldPid()) {
261 if(fPidHF->GetPidResponse()==0x0){
262 AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
263 AliInputEventHandler *inputHandler=(AliInputEventHandler*)mgr->GetInputEventHandler();
264 AliPIDResponse *pidResp=inputHandler->GetPIDResponse();
265 fPidHF->SetPidResponse(pidResp);
271 TString firedTriggerClasses=((AliAODEvent*)event)->GetFiredTriggerClasses();
272 // don't do for MC and for PbPb 2010 data
273 if(!isMC && (event->GetRunNumber()<136851 || event->GetRunNumber()>139517)) {
274 if(!firedTriggerClasses.Contains(fTriggerClass.Data())) {
276 fEvRejectionBits+=1<<kNotSelTrigger;
281 // TEMPORARY FIX FOR REFERENCES
282 // Fix references to daughter tracks
284 // AliAnalysisVertexingHF *fixer = new AliAnalysisVertexingHF();
285 // fixer->FixReferences((AliAODEvent*)event);
292 // vertex requirements
294 const AliVVertex *vertex = event->GetPrimaryVertex();
298 fEvRejectionBits+=1<<kNoVertex;
300 TString title=vertex->GetTitle();
301 if(title.Contains("Z") && fMinVtxType>1){
303 fEvRejectionBits+=1<<kNoVertex;
305 else if(title.Contains("3D") && fMinVtxType>2){
307 fEvRejectionBits+=1<<kNoVertex;
309 if(vertex->GetNContributors()<fMinVtxContr){
311 fEvRejectionBits+=1<<kTooFewVtxContrib;
313 if(TMath::Abs(vertex->GetZ())>fMaxVtxZ) {
314 fEvRejectionBits+=1<<kZVtxOutFid;
315 if(accept) fWhyRejection=6;
322 if(fOptPileup==kRejectPileupEvent){
323 Int_t cutc=(Int_t)fMinContrPileup;
324 Double_t cutz=(Double_t)fMinDzPileup;
325 if(event->IsPileupFromSPD(cutc,cutz,3.,2.,10.)) {
326 if(accept) fWhyRejection=1;
327 fEvRejectionBits+=1<<kPileupSPD;
332 // centrality selection
333 if (fUseCentrality!=kCentOff) {
334 Int_t rejection=IsEventSelectedInCentrality(event);
336 if(accept) fWhyRejection=rejection;
337 fEvRejectionBits+=1<<kOutsideCentrality;
344 //---------------------------------------------------------------------------
345 Bool_t AliRDHFCuts::AreDaughtersSelected(AliAODRecoDecayHF *d) const {
347 // Daughter tracks selection
349 if(!fTrackCuts) return kTRUE;
351 Int_t ndaughters = d->GetNDaughters();
352 AliAODVertex *vAOD = d->GetPrimaryVtx();
353 Double_t pos[3],cov[6];
355 vAOD->GetCovarianceMatrix(cov);
356 const AliESDVertex vESD(pos,cov,100.,100);
360 for(Int_t idg=0; idg<ndaughters; idg++) {
361 AliAODTrack *dgTrack = (AliAODTrack*)d->GetDaughter(idg);
362 if(!dgTrack) {retval = kFALSE; continue;}
363 //printf("charge %d\n",dgTrack->Charge());
364 if(dgTrack->Charge()==0) continue; // it's not a track, but a V0
366 if(!IsDaughterSelected(dgTrack,&vESD,fTrackCuts)) retval = kFALSE;
371 //---------------------------------------------------------------------------
372 Bool_t AliRDHFCuts::IsDaughterSelected(AliAODTrack *track,const AliESDVertex *primary,AliESDtrackCuts *cuts) const {
374 // Convert to ESDtrack, relate to vertex and check cuts
376 if(!cuts) return kTRUE;
380 // convert to ESD track here
381 AliESDtrack esdTrack(track);
382 // set the TPC cluster info
383 esdTrack.SetTPCClusterMap(track->GetTPCClusterMap());
384 esdTrack.SetTPCSharedMap(track->GetTPCSharedMap());
385 esdTrack.SetTPCPointsF(track->GetTPCNclsF());
386 // needed to calculate the impact parameters
387 esdTrack.RelateToVertex(primary,0.,3.);
389 if(!cuts->IsSelected(&esdTrack)) retval = kFALSE;
391 if(fOptPileup==kRejectTracksFromPileupVertex){
393 // we need either to have here the AOD Event,
394 // or to have the pileup vertex object
398 //---------------------------------------------------------------------------
399 void AliRDHFCuts::SetPtBins(Int_t nPtBinLimits,Float_t *ptBinLimits) {
403 delete [] fPtBinLimits;
405 printf("Changing the pt bins\n");
408 if(nPtBinLimits != fnPtBins+1){
409 cout<<"Warning: ptBinLimits dimention "<<nPtBinLimits<<" != nPtBins+1 ("<<fnPtBins+1<<")\nSetting nPtBins to "<<nPtBinLimits-1<<endl;
410 SetNPtBins(nPtBinLimits-1);
413 fnPtBinLimits = nPtBinLimits;
415 //cout<<"Changing also Global Index -> "<<fGlobalIndex<<endl;
416 fPtBinLimits = new Float_t[fnPtBinLimits];
417 for(Int_t ib=0; ib<nPtBinLimits; ib++) fPtBinLimits[ib]=ptBinLimits[ib];
421 //---------------------------------------------------------------------------
422 void AliRDHFCuts::SetVarNames(Int_t nVars,TString *varNames,Bool_t *isUpperCut){
423 // Set the variable names
428 //printf("Changing the variable names\n");
431 printf("Wrong number of variables: it has to be %d\n",fnVars);
435 fVarNames = new TString[nVars];
436 fIsUpperCut = new Bool_t[nVars];
437 for(Int_t iv=0; iv<nVars; iv++) {
438 fVarNames[iv] = varNames[iv];
439 fIsUpperCut[iv] = isUpperCut[iv];
444 //---------------------------------------------------------------------------
445 void AliRDHFCuts::SetVarsForOpt(Int_t nVars,Bool_t *forOpt) {
446 // Set the variables to be used for cuts optimization
449 delete [] fVarsForOpt;
451 //printf("Changing the variables for cut optimization\n");
454 if(nVars==0){//!=fnVars) {
455 printf("%d not accepted as number of variables: it has to be %d\n",nVars,fnVars);
460 fVarsForOpt = new Bool_t[fnVars];
461 for(Int_t iv=0; iv<fnVars; iv++) {
462 fVarsForOpt[iv]=forOpt[iv];
463 if(fVarsForOpt[iv]) fnVarsForOpt++;
469 //---------------------------------------------------------------------------
470 void AliRDHFCuts::SetUseCentrality(Int_t flag) {
472 // set centrality estimator
475 if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid) AliWarning("Centrality estimator not valid");
481 //---------------------------------------------------------------------------
482 void AliRDHFCuts::SetCuts(Int_t nVars,Int_t nPtBins,Float_t **cutsRD) {
487 printf("Wrong number of variables: it has to be %d\n",fnVars);
490 if(nPtBins!=fnPtBins) {
491 printf("Wrong number of pt bins: it has to be %d\n",fnPtBins);
495 if(!fCutsRD) fCutsRD = new Float_t[fGlobalIndex];
498 for(Int_t iv=0; iv<fnVars; iv++) {
500 for(Int_t ib=0; ib<fnPtBins; ib++) {
503 if(GetGlobalIndex(iv,ib)>=fGlobalIndex) {
504 cout<<"Overflow, exit..."<<endl;
508 fCutsRD[GetGlobalIndex(iv,ib)] = cutsRD[iv][ib];
514 //---------------------------------------------------------------------------
515 void AliRDHFCuts::SetCuts(Int_t glIndex,Float_t* cutsRDGlob){
519 if(glIndex != fGlobalIndex){
520 cout<<"Wrong array size: it has to be "<<fGlobalIndex<<endl;
523 if(!fCutsRD) fCutsRD = new Float_t[fGlobalIndex];
525 for(Int_t iGl=0;iGl<fGlobalIndex;iGl++){
526 fCutsRD[iGl] = cutsRDGlob[iGl];
530 //---------------------------------------------------------------------------
531 void AliRDHFCuts::PrintAll() const {
533 // print all cuts values
536 printf("Minimum vtx type %d\n",fMinVtxType);
537 printf("Minimum vtx contr %d\n",fMinVtxContr);
538 printf("Max vtx red chi2 %f\n",fMaxVtxRedChi2);
539 printf("Min SPD mult %d\n",fMinSPDMultiplicity);
540 printf("Use PID %d\n",(Int_t)fUsePID);
541 printf("Remove daughters from vtx %d\n",(Int_t)fRemoveDaughtersFromPrimary);
542 printf("Pileup rejection: %s\n",(fOptPileup > 0) ? "Yes" : "No");
543 if(fOptPileup==1) printf(" -- Reject pileup event");
544 if(fOptPileup==2) printf(" -- Reject tracks from pileup vtx");
545 if(fUseCentrality>0) {
546 TString estimator="";
547 if(fUseCentrality==1) estimator = "V0";
548 if(fUseCentrality==2) estimator = "Tracks";
549 if(fUseCentrality==3) estimator = "Tracklets";
550 if(fUseCentrality==4) estimator = "SPD clusters outer";
551 printf("Centrality class considered: %.1f-%.1f, estimated with %s",fMinCentrality,fMaxCentrality,estimator.Data());
555 cout<<"Array of variables"<<endl;
556 for(Int_t iv=0;iv<fnVars;iv++){
557 cout<<fVarNames[iv]<<"\t";
562 cout<<"Array of optimization"<<endl;
563 for(Int_t iv=0;iv<fnVars;iv++){
564 cout<<fVarsForOpt[iv]<<"\t";
569 cout<<"Array of upper/lower cut"<<endl;
570 for(Int_t iv=0;iv<fnVars;iv++){
571 cout<<fIsUpperCut[iv]<<"\t";
576 cout<<"Array of ptbin limits"<<endl;
577 for(Int_t ib=0;ib<fnPtBinLimits;ib++){
578 cout<<fPtBinLimits[ib]<<"\t";
583 cout<<"Matrix of cuts"<<endl;
584 for(Int_t iv=0;iv<fnVars;iv++){
585 for(Int_t ib=0;ib<fnPtBins;ib++){
586 cout<<"fCutsRD["<<iv<<"]["<<ib<<"] = "<<fCutsRD[GetGlobalIndex(iv,ib)]<<"\t";
594 //---------------------------------------------------------------------------
595 void AliRDHFCuts::GetCuts(Float_t**& cutsRD) const{
600 //cout<<"Give back a "<<fnVars<<"x"<<fnPtBins<<" matrix."<<endl;
605 //cout<<"Initialization..."<<endl;
606 cutsRD=new Float_t*[fnVars];
607 for(iv=0; iv<fnVars; iv++) {
608 cutsRD[iv] = new Float_t[fnPtBins];
612 for(Int_t iGlobal=0; iGlobal<fGlobalIndex; iGlobal++) {
613 GetVarPtIndex(iGlobal,iv,ib);
614 cutsRD[iv][ib] = fCutsRD[iGlobal];
620 //---------------------------------------------------------------------------
621 Int_t AliRDHFCuts::GetGlobalIndex(Int_t iVar,Int_t iPtBin) const{
623 // give the global index from variable and pt bin
625 return iPtBin*fnVars+iVar;
628 //---------------------------------------------------------------------------
629 void AliRDHFCuts::GetVarPtIndex(Int_t iGlob, Int_t& iVar, Int_t& iPtBin) const {
631 //give the index of the variable and of the pt bin from the global index
633 iPtBin=(Int_t)iGlob/fnVars;
639 //---------------------------------------------------------------------------
640 Int_t AliRDHFCuts::PtBin(Double_t pt) const {
642 //give the pt bin where the pt lies.
645 if(pt<fPtBinLimits[0])return ptbin;
646 for (Int_t i=0;i<fnPtBins;i++){
647 if(pt<fPtBinLimits[i+1]) {
654 //-------------------------------------------------------------------
655 Float_t AliRDHFCuts::GetCutValue(Int_t iVar,Int_t iPtBin) const {
657 // Give the value of cut set for the variable iVar and the pt bin iPtBin
660 cout<<"Cuts not iniziaisez yet"<<endl;
663 return fCutsRD[GetGlobalIndex(iVar,iPtBin)];
665 //-------------------------------------------------------------------
666 Float_t AliRDHFCuts::GetCentrality(AliAODEvent* aodEvent,AliRDHFCuts::ECentrality estimator) {
668 // Get centrality percentile
671 TClonesArray *mcArray = (TClonesArray*)((AliAODEvent*)aodEvent)->GetList()->FindObject(AliAODMCParticle::StdBranchName());
672 if(mcArray) {fUseAOD049=kFALSE;}
674 AliAODHeader *header=aodEvent->GetHeader();
675 AliCentrality *centrality=header->GetCentralityP();
677 Bool_t isSelRun=kFALSE;
678 Int_t selRun[5]={138364, 138826, 138828, 138836, 138871};
679 if(!centrality) return cent;
681 if (estimator==kCentV0M){
682 cent=(Float_t)(centrality->GetCentralityPercentile("V0M"));
684 Int_t quality = centrality->GetQuality();
686 cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");
688 Int_t runnum=aodEvent->GetRunNumber();
689 for(Int_t ir=0;ir<5;ir++){
690 if(runnum==selRun[ir]){
695 if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");
699 //temporary fix for AOD049 outliers
700 if(fUseAOD049&¢>=0){
702 AliAODVZERO* aodV0 = aodEvent->GetVZEROData();
703 v0+=aodV0->GetMTotV0A();
704 v0+=aodV0->GetMTotV0C();
705 if(cent==0&&v0<19500)return -1;//filtering issue
706 Float_t tkl = (Float_t)(aodEvent->GetTracklets()->GetNumberOfTracklets());
707 Float_t val= 1.30552 + 0.147931 * v0;
708 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};
709 if ( TMath::Abs(tkl-val) > 6.*tklSigma[(Int_t)cent] )return -1;//outlier
713 if (estimator==kCentTRK) {
714 cent=(Float_t)(centrality->GetCentralityPercentile("TRK"));
716 Int_t quality = centrality->GetQuality();
718 cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TRK");
720 Int_t runnum=aodEvent->GetRunNumber();
721 for(Int_t ir=0;ir<5;ir++){
722 if(runnum==selRun[ir]){
727 if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TRK");
732 if (estimator==kCentTKL){
733 cent=(Float_t)(centrality->GetCentralityPercentile("TKL"));
735 Int_t quality = centrality->GetQuality();
737 cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TKL");
739 Int_t runnum=aodEvent->GetRunNumber();
740 for(Int_t ir=0;ir<5;ir++){
741 if(runnum==selRun[ir]){
746 if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TKL");
751 if (estimator==kCentCL1){
752 cent=(Float_t)(centrality->GetCentralityPercentile("CL1"));
754 Int_t quality = centrality->GetQuality();
756 cent=(Float_t)centrality->GetCentralityPercentileUnchecked("CL1");
758 Int_t runnum=aodEvent->GetRunNumber();
759 for(Int_t ir=0;ir<5;ir++){
760 if(runnum==selRun[ir]){
765 if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("CL1");
770 AliWarning("Centrality estimator not valid");
779 //-------------------------------------------------------------------
780 Bool_t AliRDHFCuts::CompareCuts(const AliRDHFCuts *obj) const {
782 // Compare two cuts objects
785 Bool_t areEqual=kTRUE;
787 if(fMinVtxType!=obj->fMinVtxType) { printf("Minimum vtx type %d %d\n",fMinVtxType,obj->fMinVtxType); areEqual=kFALSE;}
789 if(fMinVtxContr!=obj->fMinVtxContr) { printf("Minimum vtx contr %d %d\n",fMinVtxContr,obj->fMinVtxContr); areEqual=kFALSE;}
791 if(TMath::Abs(fMaxVtxRedChi2-obj->fMaxVtxRedChi2)>1.e-10) { printf("Max vtx red chi2 %f %f\n",fMaxVtxRedChi2,obj->fMaxVtxRedChi2);areEqual=kFALSE;}
793 if(fMinSPDMultiplicity!=obj->fMinSPDMultiplicity) { printf("Min SPD mult %d\n %d",fMinSPDMultiplicity,obj->fMinSPDMultiplicity);areEqual=kFALSE;}
795 if(fUsePID!=obj->fUsePID) { printf("Use PID %d %d\n",(Int_t)fUsePID,(Int_t)obj->fUsePID); areEqual=kFALSE;}
797 if(fRemoveDaughtersFromPrimary!=obj->fRemoveDaughtersFromPrimary) {printf("Remove daughters from vtx %d %d\n",(Int_t)fRemoveDaughtersFromPrimary,(Int_t)obj->fRemoveDaughtersFromPrimary); areEqual=kFALSE;}
799 if(fTrackCuts->GetMinNClusterTPC()!=obj->fTrackCuts->GetMinNClusterTPC()) {printf("MinNClsTPC %d %d\n",fTrackCuts->GetMinNClusterTPC(),obj->fTrackCuts->GetMinNClusterTPC()); areEqual=kFALSE;}
801 if(fTrackCuts->GetMinNClustersITS()!=obj->fTrackCuts->GetMinNClustersITS()) {printf("MinNClsITS %d %d\n",fTrackCuts->GetMinNClustersITS(),obj->fTrackCuts->GetMinNClustersITS()); areEqual=kFALSE;}
803 if(TMath::Abs(fTrackCuts->GetMaxChi2PerClusterTPC()-obj->fTrackCuts->GetMaxChi2PerClusterTPC())>1.e-10) {printf("MaxChi2ClsTPC %f %f\n",fTrackCuts->GetMaxChi2PerClusterTPC(),obj->fTrackCuts->GetMaxChi2PerClusterTPC()); areEqual=kFALSE;}
805 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;}
809 for(Int_t iv=0;iv<fnVars;iv++) {
810 for(Int_t ib=0;ib<fnPtBins;ib++) {
811 if(TMath::Abs(fCutsRD[GetGlobalIndex(iv,ib)]-obj->fCutsRD[GetGlobalIndex(iv,ib)])>1.e-10) {
812 cout<<"fCutsRD["<<iv<<"]["<<ib<<"] = "<<fCutsRD[GetGlobalIndex(iv,ib)]<<" "<<obj->fCutsRD[GetGlobalIndex(iv,ib)]<<"\n";
821 //---------------------------------------------------------------------------
822 void AliRDHFCuts::MakeTable() const {
824 // print cuts values in table format
827 TString ptString = "pT range";
828 if(fVarNames && fPtBinLimits && fCutsRD){
829 TString firstLine(Form("* %-15s",ptString.Data()));
830 for (Int_t ivar=0; ivar<fnVars; ivar++){
831 firstLine+=Form("* %-15s ",fVarNames[ivar].Data());
836 Printf("%s",firstLine.Data());
838 for (Int_t ipt=0; ipt<fnPtBins; ipt++){
840 if (ipt==fnPtBins-1){
841 line=Form("* %5.1f < pt < inf ",fPtBinLimits[ipt]);
844 line=Form("* %5.1f < pt < %4.1f ",fPtBinLimits[ipt],fPtBinLimits[ipt+1]);
846 for (Int_t ivar=0; ivar<fnVars; ivar++){
847 line+=Form("* %-15f ",fCutsRD[GetGlobalIndex(ivar,ipt)]);
849 Printf("%s",line.Data());
857 //--------------------------------------------------------------------------
858 Bool_t AliRDHFCuts::RecalcOwnPrimaryVtx(AliAODRecoDecayHF *d,
859 AliAODEvent *aod) const
862 // Recalculate primary vertex without daughters
866 AliError("Can not remove daughters from vertex without AOD event");
870 AliAODVertex *recvtx=d->RemoveDaughtersFromPrimaryVtx(aod);
872 AliDebug(2,"Removal of daughter tracks failed");
877 //set recalculed primary vertex
878 d->SetOwnPrimaryVtx(recvtx);
883 //--------------------------------------------------------------------------
884 Bool_t AliRDHFCuts::SetMCPrimaryVtx(AliAODRecoDecayHF *d,AliAODEvent *aod) const
887 // Recalculate primary vertex without daughters
891 AliError("Can not get MC vertex without AOD event");
896 AliAODMCHeader *mcHeader =
897 (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
899 AliError("Can not get MC vertex without AODMCHeader event");
903 Double_t covmatrix[6]={0.,0.,0.,0.,0.,0.};
904 mcHeader->GetVertex(pos);
905 AliAODVertex *recvtx=new AliAODVertex(pos,covmatrix);
908 AliDebug(2,"Removal of daughter tracks failed");
912 //set recalculed primary vertex
913 d->SetOwnPrimaryVtx(recvtx);
915 d->RecalculateImpPars(recvtx,aod);
921 //--------------------------------------------------------------------------
922 void AliRDHFCuts::CleanOwnPrimaryVtx(AliAODRecoDecayHF *d,
924 AliAODVertex *origownvtx) const
927 // Clean-up own primary vertex if needed
930 if(fRemoveDaughtersFromPrimary || fUseMCVertex) {
931 d->UnsetOwnPrimaryVtx();
933 d->SetOwnPrimaryVtx(origownvtx);
934 delete origownvtx; origownvtx=NULL;
936 d->RecalculateImpPars(d->GetPrimaryVtx(),aod);
939 delete origownvtx; origownvtx=NULL;