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