]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - PWG3/vertexingHF/AliRDHFCuts.cxx
Updates for proper treatment of fQuality flag in centrality selection and for usage...
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliRDHFCuts.cxx
... / ...
CommitLineData
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
18/////////////////////////////////////////////////////////////
19//
20// Base class for cuts on AOD reconstructed heavy-flavour decay
21//
22// Author: A.Dainese, andrea.dainese@pd.infn.it
23/////////////////////////////////////////////////////////////
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 "AliRDHFCuts.h"
40
41ClassImp(AliRDHFCuts)
42
43//--------------------------------------------------------------------------
44AliRDHFCuts::AliRDHFCuts(const Char_t* name, const Char_t* title) :
45AliAnalysisCuts(name,title),
46fMinVtxType(3),
47fMinVtxContr(1),
48fMaxVtxRedChi2(1e6),
49fMaxVtxZ(1e6),
50fMinSPDMultiplicity(0),
51fTriggerMask(0),
52fTriggerClass("CINT1"),
53fTrackCuts(0),
54fnPtBins(1),
55fnPtBinLimits(1),
56fPtBinLimits(0),
57fnVars(1),
58fVarNames(0),
59fnVarsForOpt(0),
60fVarsForOpt(0),
61fGlobalIndex(1),
62fCutsRD(0),
63fIsUpperCut(0),
64fUsePID(kFALSE),
65fPidHF(0),
66fWhyRejection(0),
67fRemoveDaughtersFromPrimary(kFALSE),
68fOptPileup(0),
69fMinContrPileup(3),
70fMinDzPileup(0.6),
71fUseCentrality(0),
72fMinCentrality(0.),
73fMaxCentrality(100.),
74fFixRefs(kFALSE),
75fIsSelectedCuts(0),
76fIsSelectedPID(0),
77fMinPtCand(-1.),
78fMaxPtCand(100000.)
79{
80 //
81 // Default Constructor
82 //
83}
84//--------------------------------------------------------------------------
85AliRDHFCuts::AliRDHFCuts(const AliRDHFCuts &source) :
86 AliAnalysisCuts(source),
87 fMinVtxType(source.fMinVtxType),
88 fMinVtxContr(source.fMinVtxContr),
89 fMaxVtxRedChi2(source.fMaxVtxRedChi2),
90 fMaxVtxZ(source.fMaxVtxZ),
91 fMinSPDMultiplicity(source.fMinSPDMultiplicity),
92 fTriggerMask(source.fTriggerMask),
93 fTriggerClass(source.fTriggerClass),
94 fTrackCuts(0),
95 fnPtBins(source.fnPtBins),
96 fnPtBinLimits(source.fnPtBinLimits),
97 fPtBinLimits(0),
98 fnVars(source.fnVars),
99 fVarNames(0),
100 fnVarsForOpt(source.fnVarsForOpt),
101 fVarsForOpt(0),
102 fGlobalIndex(source.fGlobalIndex),
103 fCutsRD(0),
104 fIsUpperCut(0),
105 fUsePID(source.fUsePID),
106 fPidHF(0),
107 fWhyRejection(source.fWhyRejection),
108 fRemoveDaughtersFromPrimary(source.fRemoveDaughtersFromPrimary),
109 fOptPileup(source.fOptPileup),
110 fMinContrPileup(source.fMinContrPileup),
111 fMinDzPileup(source.fMinDzPileup),
112 fUseCentrality(source.fUseCentrality),
113 fMinCentrality(source.fMinCentrality),
114 fMaxCentrality(source.fMaxCentrality),
115 fFixRefs(source.fFixRefs),
116 fIsSelectedCuts(source.fIsSelectedCuts),
117 fIsSelectedPID(source.fIsSelectedPID),
118 fMinPtCand(source.fMinPtCand),
119 fMaxPtCand(source.fMaxPtCand)
120{
121 //
122 // Copy constructor
123 //
124 cout<<"Copy constructor"<<endl;
125 if(source.GetTrackCuts()) AddTrackCuts(source.GetTrackCuts());
126 if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);
127 if(source.fVarNames) SetVarNames(source.fnVars,source.fVarNames,source.fIsUpperCut);
128 if(source.fCutsRD) SetCuts(source.fGlobalIndex,source.fCutsRD);
129 if(source.fVarsForOpt) SetVarsForOpt(source.fnVarsForOpt,source.fVarsForOpt);
130 if(source.fPidHF) SetPidHF(source.fPidHF);
131 PrintAll();
132
133}
134//--------------------------------------------------------------------------
135AliRDHFCuts &AliRDHFCuts::operator=(const AliRDHFCuts &source)
136{
137 //
138 // assignment operator
139 //
140 if(&source == this) return *this;
141
142 AliAnalysisCuts::operator=(source);
143
144 fMinVtxType=source.fMinVtxType;
145 fMinVtxContr=source.fMinVtxContr;
146 fMaxVtxRedChi2=source.fMaxVtxRedChi2;
147 fMaxVtxZ=source.fMaxVtxZ;
148 fMinSPDMultiplicity=source.fMinSPDMultiplicity;
149 fTriggerMask=source.fTriggerMask;
150 fTriggerClass=source.fTriggerClass;
151 fnPtBins=source.fnPtBins;
152 fnVars=source.fnVars;
153 fGlobalIndex=source.fGlobalIndex;
154 fnVarsForOpt=source.fnVarsForOpt;
155 fUsePID=source.fUsePID;
156 SetPidHF(source.GetPidHF());
157 fWhyRejection=source.fWhyRejection;
158 fRemoveDaughtersFromPrimary=source.fRemoveDaughtersFromPrimary;
159 fOptPileup=source.fOptPileup;
160 fMinContrPileup=source.fMinContrPileup;
161 fMinDzPileup=source.fMinDzPileup;
162 fUseCentrality=source.fUseCentrality;
163 fMinCentrality=source.fMinCentrality;
164 fMaxCentrality=source.fMaxCentrality;
165 fFixRefs=source.fFixRefs;
166 fIsSelectedCuts=source.fIsSelectedCuts;
167 fIsSelectedPID=source.fIsSelectedPID;
168 fMinPtCand=source.fMinPtCand;
169 fMaxPtCand=source.fMaxPtCand;
170
171 if(source.GetTrackCuts()) AddTrackCuts(source.GetTrackCuts());
172 if(source.fPtBinLimits) SetPtBins(source.fnPtBinLimits,source.fPtBinLimits);
173 if(source.fVarNames) SetVarNames(source.fnVars,source.fVarNames,source.fIsUpperCut);
174 if(source.fCutsRD) SetCuts(source.fGlobalIndex,source.fCutsRD);
175 if(source.fVarsForOpt) SetVarsForOpt(source.fnVarsForOpt,source.fVarsForOpt);
176 PrintAll();
177
178 return *this;
179}
180//--------------------------------------------------------------------------
181AliRDHFCuts::~AliRDHFCuts() {
182 //
183 // Default Destructor
184 //
185 if(fPtBinLimits) {delete [] fPtBinLimits; fPtBinLimits=0;}
186 if(fVarNames) {delete [] fVarNames; fVarNames=0;}
187 if(fVarsForOpt) {delete [] fVarsForOpt; fVarsForOpt=0;}
188 if(fCutsRD) {
189 delete [] fCutsRD;
190 fCutsRD=0;
191 }
192 if(fIsUpperCut) {delete [] fIsUpperCut; fIsUpperCut=0;}
193 if(fPidHF){
194 delete fPidHF;
195 fPidHF=0;
196 }
197}
198//---------------------------------------------------------------------------
199Int_t AliRDHFCuts::IsEventSelectedInCentrality(AliVEvent *event) {
200 //
201 // Centrality selection
202 //
203
204 if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid){
205 AliWarning("Centrality estimator not valid");
206 return 3;
207 }else{
208 Float_t centvalue=GetCentrality((AliAODEvent*)event);
209 if (centvalue<-998.){//-999 if no centralityP
210 return 0;
211 }else{
212 if (centvalue<fMinCentrality || centvalue>fMaxCentrality){
213 return 2;
214 }
215 }
216 }
217 return 0;
218}
219//---------------------------------------------------------------------------
220Bool_t AliRDHFCuts::IsEventSelected(AliVEvent *event) {
221 //
222 // Event selection
223 //
224 //if(fTriggerMask && event->GetTriggerMask()!=fTriggerMask) return kFALSE;
225
226 fWhyRejection=0;
227
228 // trigger class
229 TString firedTriggerClasses=((AliAODEvent*)event)->GetFiredTriggerClasses();
230 // don't do for PbPb 2010 data
231 if(event->GetRunNumber()<136851 || event->GetRunNumber()>139517) {
232 if(!firedTriggerClasses.Contains(fTriggerClass.Data())) {
233 fWhyRejection=5;
234 return kFALSE;
235 }
236 }
237
238 // TEMPORARY FIX FOR REFERENCES
239 // Fix references to daughter tracks
240 if(fFixRefs) {
241 AliAnalysisVertexingHF *fixer = new AliAnalysisVertexingHF();
242 fixer->FixReferences((AliAODEvent*)event);
243 delete fixer;
244 }
245 //
246
247
248
249 // multiplicity cuts no implemented yet
250
251 const AliVVertex *vertex = event->GetPrimaryVertex();
252
253 if(!vertex) return kFALSE;
254
255 TString title=vertex->GetTitle();
256 if(title.Contains("Z") && fMinVtxType>1) return kFALSE;
257 if(title.Contains("3D") && fMinVtxType>2) return kFALSE;
258
259 if(vertex->GetNContributors()<fMinVtxContr) return kFALSE;
260
261 if(TMath::Abs(vertex->GetZ())>fMaxVtxZ) return kFALSE;
262
263 // switch to settings for 1-pad cls in TPC
264 if(fPidHF) {
265 if(event->GetRunNumber()>121693 && event->GetRunNumber()<136851)
266 fPidHF->SetOnePad(kTRUE);
267 if(event->GetRunNumber()>=136851 && event->GetRunNumber()<=139517)
268 fPidHF->SetPbPb(kTRUE);
269 }
270
271 if(fOptPileup==kRejectPileupEvent){
272 Int_t cutc=(Int_t)fMinContrPileup;
273 Double_t cutz=(Double_t)fMinDzPileup;
274 if(event->IsPileupFromSPD(cutc,cutz,3.,2.,10.)) {
275 fWhyRejection=1;
276 return kFALSE;
277 }
278 }
279
280 //centrality selection
281 if (fUseCentrality!=kCentOff) {
282 Int_t rejection=IsEventSelectedInCentrality(event);
283 if(rejection>1){
284 fWhyRejection=rejection;
285 return kFALSE;
286 }
287 }
288
289 return kTRUE;
290}
291//---------------------------------------------------------------------------
292Bool_t AliRDHFCuts::AreDaughtersSelected(AliAODRecoDecayHF *d) const {
293 //
294 // Daughter tracks selection
295 //
296 if(!fTrackCuts) return kTRUE;
297
298 Int_t ndaughters = d->GetNDaughters();
299 AliAODVertex *vAOD = d->GetPrimaryVtx();
300 Double_t pos[3],cov[6];
301 vAOD->GetXYZ(pos);
302 vAOD->GetCovarianceMatrix(cov);
303 const AliESDVertex vESD(pos,cov,100.,100);
304
305 Bool_t retval=kTRUE;
306
307 for(Int_t idg=0; idg<ndaughters; idg++) {
308 AliAODTrack *dgTrack = (AliAODTrack*)d->GetDaughter(idg);
309 if(!dgTrack) {retval = kFALSE; continue;}
310 //printf("charge %d\n",dgTrack->Charge());
311 if(dgTrack->Charge()==0) continue; // it's not a track, but a V0
312
313 if(!IsDaughterSelected(dgTrack,&vESD,fTrackCuts)) retval = kFALSE;
314 }
315
316 return retval;
317}
318//---------------------------------------------------------------------------
319Bool_t AliRDHFCuts::IsDaughterSelected(AliAODTrack *track,const AliESDVertex *primary,AliESDtrackCuts *cuts) const {
320 //
321 // Convert to ESDtrack, relate to vertex and check cuts
322 //
323 if(!cuts) return kTRUE;
324
325 Bool_t retval=kTRUE;
326
327 // convert to ESD track here
328 AliESDtrack esdTrack(track);
329 // set the TPC cluster info
330 esdTrack.SetTPCClusterMap(track->GetTPCClusterMap());
331 esdTrack.SetTPCSharedMap(track->GetTPCSharedMap());
332 esdTrack.SetTPCPointsF(track->GetTPCNclsF());
333 // needed to calculate the impact parameters
334 esdTrack.RelateToVertex(primary,0.,3.);
335
336 if(!cuts->IsSelected(&esdTrack)) retval = kFALSE;
337
338 if(fOptPileup==kRejectTracksFromPileupVertex){
339 // to be implemented
340 // we need either to have here the AOD Event,
341 // or to have the pileup vertex object
342 }
343 return retval;
344}
345//---------------------------------------------------------------------------
346void AliRDHFCuts::SetPtBins(Int_t nPtBinLimits,Float_t *ptBinLimits) {
347 // Set the pt bins
348
349 if(fPtBinLimits) {
350 delete [] fPtBinLimits;
351 fPtBinLimits = NULL;
352 printf("Changing the pt bins\n");
353 }
354
355 if(nPtBinLimits != fnPtBins+1){
356 cout<<"Warning: ptBinLimits dimention "<<nPtBinLimits<<" != nPtBins+1 ("<<fnPtBins+1<<")\nSetting nPtBins to "<<nPtBinLimits-1<<endl;
357 SetNPtBins(nPtBinLimits-1);
358 }
359
360 fnPtBinLimits = nPtBinLimits;
361 SetGlobalIndex();
362 //cout<<"Changing also Global Index -> "<<fGlobalIndex<<endl;
363 fPtBinLimits = new Float_t[fnPtBinLimits];
364 for(Int_t ib=0; ib<nPtBinLimits; ib++) fPtBinLimits[ib]=ptBinLimits[ib];
365
366 return;
367}
368//---------------------------------------------------------------------------
369void AliRDHFCuts::SetVarNames(Int_t nVars,TString *varNames,Bool_t *isUpperCut){
370 // Set the variable names
371
372 if(fVarNames) {
373 delete [] fVarNames;
374 fVarNames = NULL;
375 //printf("Changing the variable names\n");
376 }
377 if(nVars!=fnVars){
378 printf("Wrong number of variables: it has to be %d\n",fnVars);
379 return;
380 }
381 //fnVars=nVars;
382 fVarNames = new TString[nVars];
383 fIsUpperCut = new Bool_t[nVars];
384 for(Int_t iv=0; iv<nVars; iv++) {
385 fVarNames[iv] = varNames[iv];
386 fIsUpperCut[iv] = isUpperCut[iv];
387 }
388
389 return;
390}
391//---------------------------------------------------------------------------
392void AliRDHFCuts::SetVarsForOpt(Int_t nVars,Bool_t *forOpt) {
393 // Set the variables to be used for cuts optimization
394
395 if(fVarsForOpt) {
396 delete [] fVarsForOpt;
397 fVarsForOpt = NULL;
398 //printf("Changing the variables for cut optimization\n");
399 }
400
401 if(nVars==0){//!=fnVars) {
402 printf("%d not accepted as number of variables: it has to be %d\n",nVars,fnVars);
403 return;
404 }
405
406 fnVarsForOpt = 0;
407 fVarsForOpt = new Bool_t[fnVars];
408 for(Int_t iv=0; iv<fnVars; iv++) {
409 fVarsForOpt[iv]=forOpt[iv];
410 if(fVarsForOpt[iv]) fnVarsForOpt++;
411 }
412
413 return;
414}
415
416//---------------------------------------------------------------------------
417void AliRDHFCuts::SetUseCentrality(Int_t flag) {
418 //
419 // set centrality estimator
420 //
421 fUseCentrality=flag;
422 if(fUseCentrality<kCentOff||fUseCentrality>=kCentInvalid) AliWarning("Centrality estimator not valid");
423
424 return;
425}
426
427
428//---------------------------------------------------------------------------
429void AliRDHFCuts::SetCuts(Int_t nVars,Int_t nPtBins,Float_t **cutsRD) {
430 //
431 // store the cuts
432 //
433 if(nVars!=fnVars) {
434 printf("Wrong number of variables: it has to be %d\n",fnVars);
435 return;
436 }
437 if(nPtBins!=fnPtBins) {
438 printf("Wrong number of pt bins: it has to be %d\n",fnPtBins);
439 return;
440 }
441
442 if(!fCutsRD) fCutsRD = new Float_t[fGlobalIndex];
443
444
445 for(Int_t iv=0; iv<fnVars; iv++) {
446
447 for(Int_t ib=0; ib<fnPtBins; ib++) {
448
449 //check
450 if(GetGlobalIndex(iv,ib)>=fGlobalIndex) {
451 cout<<"Overflow, exit..."<<endl;
452 return;
453 }
454
455 fCutsRD[GetGlobalIndex(iv,ib)] = cutsRD[iv][ib];
456
457 }
458 }
459 return;
460}
461//---------------------------------------------------------------------------
462void AliRDHFCuts::SetCuts(Int_t glIndex,Float_t* cutsRDGlob){
463 //
464 // store the cuts
465 //
466 if(glIndex != fGlobalIndex){
467 cout<<"Wrong array size: it has to be "<<fGlobalIndex<<endl;
468 return;
469 }
470 if(!fCutsRD) fCutsRD = new Float_t[fGlobalIndex];
471
472 for(Int_t iGl=0;iGl<fGlobalIndex;iGl++){
473 fCutsRD[iGl] = cutsRDGlob[iGl];
474 }
475 return;
476}
477//---------------------------------------------------------------------------
478void AliRDHFCuts::PrintAll() const {
479 //
480 // print all cuts values
481 //
482
483 printf("Minimum vtx type %d\n",fMinVtxType);
484 printf("Minimum vtx contr %d\n",fMinVtxContr);
485 printf("Max vtx red chi2 %f\n",fMaxVtxRedChi2);
486 printf("Min SPD mult %d\n",fMinSPDMultiplicity);
487 printf("Use PID %d\n",(Int_t)fUsePID);
488 printf("Remove daughters from vtx %d\n",(Int_t)fRemoveDaughtersFromPrimary);
489 printf("Pileup rejection: %s\n",(fOptPileup > 0) ? "Yes" : "No");
490 if(fOptPileup==1) printf(" -- Reject pileup event");
491 if(fOptPileup==2) printf(" -- Reject tracks from pileup vtx");
492 if(fUseCentrality>0) {
493 TString estimator="";
494 if(fUseCentrality==1) estimator = "V0";
495 if(fUseCentrality==2) estimator = "Tracks";
496 if(fUseCentrality==3) estimator = "Tracklets";
497 if(fUseCentrality==4) estimator = "SPD clusters outer";
498 printf("Centrality class considered: %.1f-%.1f, estimated with %s",fMinCentrality,fMaxCentrality,estimator.Data());
499 }
500
501 if(fVarNames){
502 cout<<"Array of variables"<<endl;
503 for(Int_t iv=0;iv<fnVars;iv++){
504 cout<<fVarNames[iv]<<"\t";
505 }
506 cout<<endl;
507 }
508 if(fVarsForOpt){
509 cout<<"Array of optimization"<<endl;
510 for(Int_t iv=0;iv<fnVars;iv++){
511 cout<<fVarsForOpt[iv]<<"\t";
512 }
513 cout<<endl;
514 }
515 if(fIsUpperCut){
516 cout<<"Array of upper/lower cut"<<endl;
517 for(Int_t iv=0;iv<fnVars;iv++){
518 cout<<fIsUpperCut[iv]<<"\t";
519 }
520 cout<<endl;
521 }
522 if(fPtBinLimits){
523 cout<<"Array of ptbin limits"<<endl;
524 for(Int_t ib=0;ib<fnPtBinLimits;ib++){
525 cout<<fPtBinLimits[ib]<<"\t";
526 }
527 cout<<endl;
528 }
529 if(fCutsRD){
530 cout<<"Matrix of cuts"<<endl;
531 for(Int_t iv=0;iv<fnVars;iv++){
532 for(Int_t ib=0;ib<fnPtBins;ib++){
533 cout<<"fCutsRD["<<iv<<"]["<<ib<<"] = "<<fCutsRD[GetGlobalIndex(iv,ib)]<<"\t";
534 }
535 cout<<endl;
536 }
537 cout<<endl;
538 }
539 return;
540}
541//---------------------------------------------------------------------------
542void AliRDHFCuts::GetCuts(Float_t**& cutsRD) const{
543 //
544 // get the cuts
545 //
546
547 //cout<<"Give back a "<<fnVars<<"x"<<fnPtBins<<" matrix."<<endl;
548
549
550 Int_t iv,ib;
551 if(!cutsRD) {
552 //cout<<"Initialization..."<<endl;
553 cutsRD=new Float_t*[fnVars];
554 for(iv=0; iv<fnVars; iv++) {
555 cutsRD[iv] = new Float_t[fnPtBins];
556 }
557 }
558
559 for(Int_t iGlobal=0; iGlobal<fGlobalIndex; iGlobal++) {
560 GetVarPtIndex(iGlobal,iv,ib);
561 cutsRD[iv][ib] = fCutsRD[iGlobal];
562 }
563
564 return;
565}
566
567//---------------------------------------------------------------------------
568Int_t AliRDHFCuts::GetGlobalIndex(Int_t iVar,Int_t iPtBin) const{
569 //
570 // give the global index from variable and pt bin
571 //
572 return iPtBin*fnVars+iVar;
573}
574
575//---------------------------------------------------------------------------
576void AliRDHFCuts::GetVarPtIndex(Int_t iGlob, Int_t& iVar, Int_t& iPtBin) const {
577 //
578 //give the index of the variable and of the pt bin from the global index
579 //
580 iPtBin=(Int_t)iGlob/fnVars;
581 iVar=iGlob%fnVars;
582
583 return;
584}
585
586//---------------------------------------------------------------------------
587Int_t AliRDHFCuts::PtBin(Double_t pt) const {
588 //
589 //give the pt bin where the pt lies.
590 //
591 Int_t ptbin=-1;
592 if(pt<fPtBinLimits[0])return ptbin;
593 for (Int_t i=0;i<fnPtBins;i++){
594 if(pt<fPtBinLimits[i+1]) {
595 ptbin=i;
596 break;
597 }
598 }
599 return ptbin;
600}
601//-------------------------------------------------------------------
602Float_t AliRDHFCuts::GetCutValue(Int_t iVar,Int_t iPtBin) const {
603 //
604 // Give the value of cut set for the variable iVar and the pt bin iPtBin
605 //
606 if(!fCutsRD){
607 cout<<"Cuts not iniziaisez yet"<<endl;
608 return 0;
609 }
610 return fCutsRD[GetGlobalIndex(iVar,iPtBin)];
611}
612//-------------------------------------------------------------------
613Float_t AliRDHFCuts::GetCentrality(AliAODEvent* aodEvent,AliRDHFCuts::ECentrality estimator) const {
614 //
615 // Get centrality percentile
616 //
617 AliAODHeader *header=aodEvent->GetHeader();
618 AliCentrality *centrality=header->GetCentralityP();
619 Float_t cent=-999.;
620 Bool_t isSelRun=kFALSE;
621 Int_t selRun[5]={138364, 138826, 138828, 138836, 138871};
622 if(!centrality) return cent;
623 else{
624 if (estimator==kCentV0M){
625 cent=(Float_t)(centrality->GetCentralityPercentile("V0M"));
626 if(cent<0){
627 Int_t quality = centrality->GetQuality();
628 if(quality<=1){
629 cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");
630 }else{
631 Int_t runnum=aodEvent->GetRunNumber();
632 for(Int_t ir=0;ir<5;ir++){
633 if(runnum==selRun[ir]){
634 isSelRun=kTRUE;
635 break;
636 }
637 }
638 if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("V0M");
639 }
640 }
641 }
642 else {
643 if (estimator==kCentTRK) {
644 cent=(Float_t)(centrality->GetCentralityPercentile("TRK"));
645 if(cent<0){
646 Int_t quality = centrality->GetQuality();
647 if(quality<=1){
648 cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TRK");
649 }else{
650 Int_t runnum=aodEvent->GetRunNumber();
651 for(Int_t ir=0;ir<5;ir++){
652 if(runnum==selRun[ir]){
653 isSelRun=kTRUE;
654 break;
655 }
656 }
657 if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TRK");
658 }
659 }
660 }
661 else{
662 if (estimator==kCentTKL){
663 cent=(Float_t)(centrality->GetCentralityPercentile("TKL"));
664 if(cent<0){
665 Int_t quality = centrality->GetQuality();
666 if(quality<=1){
667 cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TKL");
668 }else{
669 Int_t runnum=aodEvent->GetRunNumber();
670 for(Int_t ir=0;ir<5;ir++){
671 if(runnum==selRun[ir]){
672 isSelRun=kTRUE;
673 break;
674 }
675 }
676 if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("TKL");
677 }
678 }
679 }
680 else{
681 if (estimator==kCentCL1){
682 cent=(Float_t)(centrality->GetCentralityPercentile("CL1"));
683 if(cent<0){
684 Int_t quality = centrality->GetQuality();
685 if(quality<=1){
686 cent=(Float_t)centrality->GetCentralityPercentileUnchecked("CL1");
687 }else{
688 Int_t runnum=aodEvent->GetRunNumber();
689 for(Int_t ir=0;ir<5;ir++){
690 if(runnum==selRun[ir]){
691 isSelRun=kTRUE;
692 break;
693 }
694 }
695 if((quality==8||quality==9)&&isSelRun)cent=(Float_t)centrality->GetCentralityPercentileUnchecked("CL1");
696 }
697 }
698 }
699 else {
700 AliWarning("Centrality estimator not valid");
701
702 }
703 }
704 }
705 }
706 }
707 return cent;
708}
709//-------------------------------------------------------------------
710Bool_t AliRDHFCuts::CompareCuts(const AliRDHFCuts *obj) const {
711 //
712 // Compare two cuts objects
713 //
714
715 Bool_t areEqual=kTRUE;
716
717 if(fMinVtxType!=obj->fMinVtxType) { printf("Minimum vtx type %d %d\n",fMinVtxType,obj->fMinVtxType); areEqual=kFALSE;}
718
719 if(fMinVtxContr!=obj->fMinVtxContr) { printf("Minimum vtx contr %d %d\n",fMinVtxContr,obj->fMinVtxContr); areEqual=kFALSE;}
720
721 if(TMath::Abs(fMaxVtxRedChi2-obj->fMaxVtxRedChi2)>1.e-10) { printf("Max vtx red chi2 %f %f\n",fMaxVtxRedChi2,obj->fMaxVtxRedChi2);areEqual=kFALSE;}
722
723 if(fMinSPDMultiplicity!=obj->fMinSPDMultiplicity) { printf("Min SPD mult %d\n %d",fMinSPDMultiplicity,obj->fMinSPDMultiplicity);areEqual=kFALSE;}
724
725 if(fUsePID!=obj->fUsePID) { printf("Use PID %d %d\n",(Int_t)fUsePID,(Int_t)obj->fUsePID); areEqual=kFALSE;}
726
727 if(fRemoveDaughtersFromPrimary!=obj->fRemoveDaughtersFromPrimary) {printf("Remove daughters from vtx %d %d\n",(Int_t)fRemoveDaughtersFromPrimary,(Int_t)obj->fRemoveDaughtersFromPrimary); areEqual=kFALSE;}
728 if(fTrackCuts){
729 if(fTrackCuts->GetMinNClusterTPC()!=obj->fTrackCuts->GetMinNClusterTPC()) {printf("MinNClsTPC %d %d\n",fTrackCuts->GetMinNClusterTPC(),obj->fTrackCuts->GetMinNClusterTPC()); areEqual=kFALSE;}
730
731 if(fTrackCuts->GetMinNClustersITS()!=obj->fTrackCuts->GetMinNClustersITS()) {printf("MinNClsITS %d %d\n",fTrackCuts->GetMinNClustersITS(),obj->fTrackCuts->GetMinNClustersITS()); areEqual=kFALSE;}
732
733 if(TMath::Abs(fTrackCuts->GetMaxChi2PerClusterTPC()-obj->fTrackCuts->GetMaxChi2PerClusterTPC())>1.e-10) {printf("MaxChi2ClsTPC %f %f\n",fTrackCuts->GetMaxChi2PerClusterTPC(),obj->fTrackCuts->GetMaxChi2PerClusterTPC()); areEqual=kFALSE;}
734
735 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;}
736 }
737
738 if(fCutsRD) {
739 for(Int_t iv=0;iv<fnVars;iv++) {
740 for(Int_t ib=0;ib<fnPtBins;ib++) {
741 if(TMath::Abs(fCutsRD[GetGlobalIndex(iv,ib)]-obj->fCutsRD[GetGlobalIndex(iv,ib)])>1.e-10) {
742 cout<<"fCutsRD["<<iv<<"]["<<ib<<"] = "<<fCutsRD[GetGlobalIndex(iv,ib)]<<" "<<obj->fCutsRD[GetGlobalIndex(iv,ib)]<<"\n";
743 areEqual=kFALSE;
744 }
745 }
746 }
747 }
748
749 return areEqual;
750}
751//---------------------------------------------------------------------------
752void AliRDHFCuts::MakeTable() const {
753 //
754 // print cuts values in table format
755 //
756
757 TString ptString = "pT range";
758 if(fVarNames && fPtBinLimits && fCutsRD){
759 TString firstLine(Form("* %-15s",ptString.Data()));
760 for (Int_t ivar=0; ivar<fnVars; ivar++){
761 firstLine+=Form("* %-15s ",fVarNames[ivar].Data());
762 if (ivar == fnVars){
763 firstLine+="*\n";
764 }
765 }
766 Printf("%s",firstLine.Data());
767
768 for (Int_t ipt=0; ipt<fnPtBins; ipt++){
769 TString line;
770 if (ipt==fnPtBins-1){
771 line=Form("* %5.1f < pt < inf ",fPtBinLimits[ipt]);
772 }
773 else{
774 line=Form("* %5.1f < pt < %4.1f ",fPtBinLimits[ipt],fPtBinLimits[ipt+1]);
775 }
776 for (Int_t ivar=0; ivar<fnVars; ivar++){
777 line+=Form("* %-15f ",fCutsRD[GetGlobalIndex(ivar,ipt)]);
778 }
779 Printf("%s",line.Data());
780 }
781
782 }
783
784
785 return;
786}
787//--------------------------------------------------------------------------
788Bool_t AliRDHFCuts::RecalcOwnPrimaryVtx(AliAODRecoDecayHF *d,AliAODEvent *aod,
789 AliAODVertex *origownvtx,AliAODVertex *recvtx) const
790{
791 //
792 // Recalculate primary vertex without daughters
793 //
794
795 if(!aod) {
796 AliError("Can not remove daughters from vertex without AOD event");
797 return kFALSE;
798 }
799 if(d->GetOwnPrimaryVtx()) origownvtx=new AliAODVertex(*d->GetOwnPrimaryVtx());
800 recvtx=d->RemoveDaughtersFromPrimaryVtx(aod);
801 if(!recvtx){
802 AliDebug(2,"Removal of daughter tracks failed");
803 if(origownvtx){
804 delete origownvtx;
805 origownvtx=NULL;
806 }
807 return kFALSE;
808 }
809 //set recalculed primary vertex
810 d->SetOwnPrimaryVtx(recvtx);
811 delete recvtx; recvtx=NULL;
812
813 return kTRUE;
814}
815//--------------------------------------------------------------------------
816void AliRDHFCuts::CleanOwnPrimaryVtx(AliAODRecoDecayHF *d,AliAODVertex *origownvtx) const
817{
818 //
819 // Clean-up own primary vertex if needed
820 //
821
822 if(origownvtx) {
823 d->SetOwnPrimaryVtx(origownvtx);
824 delete origownvtx;
825 origownvtx=NULL;
826 } else if(fRemoveDaughtersFromPrimary) {
827 d->UnsetOwnPrimaryVtx();
828 AliDebug(3,"delete new vertex\n");
829 }
830
831 return;
832}
833
834
835