fixing casts
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliCFVertexingHF.cxx
CommitLineData
379592af 1/**************************************************************************
2 * Copyright(c) 1998-2009, 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
661d6c29 16// $Id$
17
379592af 18//-----------------------------------------------------------------------
19// Class for HF corrections as a function of many variables and step
20// Author : C. Zampolli, CERN
21// D. Caffarri, Univ & INFN Padova caffarri@pd.infn.it
22// Base class for HF Unfolding - agrelli@uu.nl
23//-----------------------------------------------------------------------
24
25#include "TParticle.h"
26#include "TClonesArray.h"
27#include "AliAODMCParticle.h"
28#include "AliAODRecoDecayHF.h"
29#include "AliAODRecoDecayHF2Prong.h"
30#include "AliAODRecoDecayHF3Prong.h"
31#include "AliAODRecoDecayHF4Prong.h"
367e9aa3 32#include "AliAODRecoCascadeHF.h"
379592af 33#include "AliAODMCHeader.h"
34#include "AliAODEvent.h"
35#include "AliLog.h"
36#include "AliESDtrackCuts.h"
37#include "AliESDtrack.h"
ec5288c3 38#include "AliCFTaskVertexingHF.h"
379592af 39
40#include "AliCFVertexingHF.h"
41
42//___________________________________________________________
43AliCFVertexingHF::AliCFVertexingHF() :
44 fmcArray(0x0),
45 fRecoCandidate(0),
46 fmcPartCandidate(0x0),
47 fNDaughters(0),
48 fNVar(0),
49 fzPrimVertex(0),
50 fzMCVertex(0),
51 fFillFromGenerated(0),
52 fOriginDselection(0),
53 fKeepDfromB(kFALSE),
54 fKeepDfromBOnly(kFALSE),
55 fmcLabel(0),
6e2e4f50 56 fProngs(-1),
b7af2639 57 fLabelArray(0x0),
2bf2e62b 58 fCentValue(0.),
59 fPtAccCut(0x0),
fbec9fa9 60 fEtaAccCut(0x0),
1f780958 61 fFakeSelection(0),
0ada222f 62 fFake(1.), // setting to MC value
ec5288c3 63 fRejectIfNoQuark(kFALSE),
64 fMultiplicity(0.),
65 fConfiguration(AliCFTaskVertexingHF::kCheetah) // by default, setting the fast configuration
379592af 66{
f2703bd2 67 //
68 // constructor
69 //
70
6e2e4f50 71
379592af 72 return;
73}
74
75
76
77//_____________________________________________________
78AliCFVertexingHF::AliCFVertexingHF(TClonesArray *mcArray, UShort_t originDselection) :
79 fmcArray(mcArray),
80 fRecoCandidate(0),
81 fmcPartCandidate(0x0),
82 fNDaughters(0),
83 fNVar(0),
84 fzPrimVertex(0),
85 fzMCVertex(0),
86 fFillFromGenerated(0),
87 fOriginDselection(0),
88 fKeepDfromB(kFALSE),
89 fKeepDfromBOnly(kFALSE),
90 fmcLabel(0),
6e2e4f50 91 fProngs(-1),
b7af2639 92 fLabelArray(0x0),
2bf2e62b 93 fCentValue(0.),
94 fPtAccCut(0x0),
fbec9fa9 95 fEtaAccCut(0x0),
1f780958 96 fFakeSelection(0),
0ada222f 97 fFake(1.), // setting to MC value
ec5288c3 98 fRejectIfNoQuark(kFALSE),
99 fMultiplicity(0.),
100 fConfiguration(AliCFTaskVertexingHF::kCheetah) // by default, setting the fast configuration
379592af 101{
f2703bd2 102 //
103 // constructor with mcArray
104 //
105
106 SetDselection(originDselection);
107 return;
379592af 108}
109
379592af 110//_______________________________________________________
111AliCFVertexingHF::~AliCFVertexingHF()
112{
f2703bd2 113 //
114 // destructor
115 //
379592af 116
f2703bd2 117 if (fmcArray) fmcArray = 0x0;
118 if (fRecoCandidate) fRecoCandidate = 0x0;
119 if (fmcPartCandidate) fmcPartCandidate = 0x0;
6e2e4f50 120 if (fLabelArray){
b7af2639 121 delete [] fLabelArray;
122 fLabelArray = 0x0;
6e2e4f50 123 }
2bf2e62b 124 if (fPtAccCut){
125 delete [] fPtAccCut;
126 fPtAccCut = 0x0;
127 }
128 if (fEtaAccCut){
129 delete [] fEtaAccCut;
130 fEtaAccCut = 0x0;
131 }
f2703bd2 132}
379592af 133
134//_____________________________________________________
f2703bd2 135AliCFVertexingHF& AliCFVertexingHF::operator=(const AliCFVertexingHF& c)
136{
137 //
138 // assigment operator
139 //
140
379592af 141 if (this!= &c){
142 TObject::operator=(c);
59fee8a6 143 delete fmcArray;
d95be5e9 144 fmcArray = new TClonesArray(*(c.fmcArray));
59fee8a6 145 delete fRecoCandidate;
d95be5e9 146 fRecoCandidate = new AliAODRecoDecayHF(*(c.fRecoCandidate));
59fee8a6 147 delete fmcPartCandidate;
d95be5e9 148 fmcPartCandidate = new AliAODMCParticle(*(c.fmcPartCandidate));
379592af 149 fNDaughters = c.fNDaughters;
150 fNVar = c.fNVar;
151 fzPrimVertex = c.fzPrimVertex;
152 fzMCVertex = c.fzMCVertex;
153 fFillFromGenerated = c.fFillFromGenerated;
154 fOriginDselection = c.fOriginDselection;
155 fKeepDfromB = c.fKeepDfromB;
156 fKeepDfromBOnly = c.fKeepDfromBOnly;
157 fmcLabel = c.fmcLabel;
6e2e4f50 158 fProngs=c.fProngs;
b7af2639 159 fCentValue=c.fCentValue;
fbec9fa9 160 fFakeSelection=c.fFakeSelection;
1f780958 161 fFake=c.fFake;
0ada222f 162 fRejectIfNoQuark=c.fRejectIfNoQuark;
6e2e4f50 163 if (fProngs > 0){
59fee8a6 164 delete [] fLabelArray;
165 delete [] fPtAccCut;
166 delete [] fEtaAccCut;
6e2e4f50 167 fLabelArray = new Int_t[fProngs];
661d6c29 168 fPtAccCut = new Float_t[fProngs];
169 fEtaAccCut = new Float_t[fProngs];
2bf2e62b 170 for(Int_t iP=0; iP<fProngs; iP++){
171 fLabelArray[iP]=c.fLabelArray[iP];
172 fPtAccCut[iP]=c.fPtAccCut[iP];
173 fEtaAccCut[iP]=c.fEtaAccCut[iP];
174 }
6e2e4f50 175 }
ec5288c3 176 fMultiplicity=c.fMultiplicity;
177 fConfiguration=c.fConfiguration;
379592af 178 }
179
180 return *this;
181}
182
183//____________________________________________________
184AliCFVertexingHF::AliCFVertexingHF(const AliCFVertexingHF &c) :
185 TObject(c),
cf1c5f6e 186 fmcArray(0),
187 fRecoCandidate(0),
188 fmcPartCandidate(0),
379592af 189 fNDaughters(c.fNDaughters),
190 fNVar(c.fNVar),
191 fzPrimVertex(c.fzPrimVertex),
192 fzMCVertex(c.fzMCVertex),
193 fFillFromGenerated(c.fFillFromGenerated),
194 fOriginDselection (c.fOriginDselection),
195 fKeepDfromB (c.fKeepDfromB),
196 fKeepDfromBOnly (c.fKeepDfromBOnly),
197 fmcLabel(c.fmcLabel),
6e2e4f50 198 fProngs(c.fProngs),
b7af2639 199 fLabelArray(0x0),
2bf2e62b 200 fCentValue(c.fCentValue),
201 fPtAccCut(0x0),
fbec9fa9 202 fEtaAccCut(0x0),
1f780958 203 fFakeSelection(c.fFakeSelection),
0ada222f 204 fFake(c.fFake),
ec5288c3 205 fRejectIfNoQuark(c.fRejectIfNoQuark),
206 fMultiplicity(c.fMultiplicity),
207 fConfiguration(c.fConfiguration)
f2703bd2 208{
75f0bf8c 209 //
210 //copy constructor
211 //
212 delete fmcArray;
d95be5e9 213 fmcArray = new TClonesArray(*(c.fmcArray));
75f0bf8c 214 delete fRecoCandidate;
d95be5e9 215 fRecoCandidate = new AliAODRecoDecayHF(*(c.fRecoCandidate));
75f0bf8c 216 delete fmcPartCandidate;
d95be5e9 217 fmcPartCandidate = new AliAODMCParticle(*(c.fmcPartCandidate));
75f0bf8c 218 if (fProngs > 0){
219 delete [] fLabelArray;
220 delete [] fPtAccCut;
221 delete [] fEtaAccCut;
222 fLabelArray = new Int_t[fProngs];
223 fPtAccCut = new Float_t[fProngs];
224 fEtaAccCut = new Float_t[fProngs];
225 if (c.fLabelArray) memcpy(fLabelArray,c.fLabelArray,fProngs*sizeof(Int_t));
226 if (c.fPtAccCut) memcpy(fPtAccCut,c.fPtAccCut,fProngs*sizeof(Int_t));
227 if (c.fEtaAccCut) memcpy(fEtaAccCut,c.fEtaAccCut,fProngs*sizeof(Int_t));
228 }
379592af 229}
230
231//___________________________________________________________
f2703bd2 232void AliCFVertexingHF::SetDselection(UShort_t originDselection)
233{
234 // setting the way the D0 will be selected
235 // 0 --> only from c quarks
236 // 1 --> only from b quarks
237 // 2 --> from both c quarks and b quarks
238
239 fOriginDselection = originDselection;
240
241 if (fOriginDselection == 0) {
242 fKeepDfromB = kFALSE;
243 fKeepDfromBOnly = kFALSE;
244 }
245
246 if (fOriginDselection == 1) {
247 fKeepDfromB = kTRUE;
248 fKeepDfromBOnly = kTRUE;
249 }
250
251 if (fOriginDselection == 2) {
252 fKeepDfromB = kTRUE;
253 fKeepDfromBOnly = kFALSE;
254 }
255
256 return;
379592af 257}
258
259//______________________________________________________
f2703bd2 260void AliCFVertexingHF::SetMCCandidateParam(Int_t label)
261{
262 //
263 // setting the parameters (candidate and n. daughters)
264 //
379592af 265
8a983038 266 fmcPartCandidate = dynamic_cast <AliAODMCParticle*> (fmcArray->At(label));
267 if (fmcPartCandidate){
268 fNDaughters = fmcPartCandidate->GetNDaughters();
269 }
270 else {
271 AliError(Form("Dynamic cast failed, fNdaughters will remain set to %d",fNDaughters));
272 }
f2703bd2 273 return;
379592af 274}
275
379592af 276//____________________________________________________________
f2703bd2 277Int_t AliCFVertexingHF::MCcquarkCounting(AliAODMCParticle* mcPart) const
278{
279 //
280 // counting the c-quarks
281 //
282
283 Int_t cquarks = 0;
8a983038 284 if (mcPart) {
285 if (mcPart->GetPdgCode() == 4) cquarks++;
286 if (mcPart->GetPdgCode() == -4) cquarks++;
287 }
288 else {
f2703bd2 289 AliWarning("Particle not found in tree, skipping\n");
290 return cquarks;
291 }
292
293 return cquarks;
379592af 294}
295
379592af 296//________________________________________________________
f2703bd2 297Bool_t AliCFVertexingHF::CheckMCPartFamily(AliAODMCParticle */*mcPart*/, TClonesArray */*mcArray*/) const
298{
299 //
300 //checking the family
301 //
302
303 Int_t pdgGranma = CheckOrigin();
0ada222f 304
4539d15b 305 AliDebug(3, Form("pdgGranma = %d", pdgGranma));
306
0ada222f 307 if (pdgGranma == -99999){
4539d15b 308 AliDebug(2, "This particle does not have a quark in his genealogy\n");
309 return kFALSE;
0ada222f 310 }
f2703bd2 311 if (pdgGranma == -9999){
4539d15b 312 AliDebug(2,"This particle come from a B decay channel but according to the settings of the task, we keep only the prompt charm particles\n");
313 return kFALSE;
f2703bd2 314 }
315
316 if (pdgGranma == -999){
4539d15b 317 AliDebug(2,"This particle come from a prompt charm particles but according to the settings of the task, we want only the ones coming from B\n");
318 return kFALSE;
f2703bd2 319 }
320
321 if (!CheckMCDaughters()) {
b7af2639 322 AliDebug(3, "CheckMCDaughters false");
323 return kFALSE;
f2703bd2 324 }
325 if (!CheckMCChannelDecay()) {
4539d15b 326 AliDebug(3, "CheckMCChannelDecay false");
327 return kFALSE;
f2703bd2 328 }
329 return kTRUE;
379592af 330}
331
332//_________________________________________________________________________________________________
f2703bd2 333Int_t AliCFVertexingHF::CheckOrigin() const
334{
335 //
336 // checking whether the mother of the particles come from a charm or a bottom quark
337 //
338
339 Int_t pdgGranma = 0;
340 Int_t mother = 0;
341 mother = fmcPartCandidate->GetMother();
342 Int_t istep = 0;
f036863c 343 Int_t abspdgGranma =0;
0ada222f 344 Bool_t isFromB=kFALSE;
345 Bool_t isQuarkFound=kFALSE;
ff12b981 346 while (mother >=0 ){
f2703bd2 347 istep++;
348 AliDebug(2,Form("mother at step %d = %d", istep, mother));
349 AliAODMCParticle* mcGranma = dynamic_cast<AliAODMCParticle*>(fmcArray->At(mother));
8a983038 350 if (mcGranma){
351 pdgGranma = mcGranma->GetPdgCode();
352 AliDebug(2,Form("Pdg mother at step %d = %d", istep, pdgGranma));
353 abspdgGranma = TMath::Abs(pdgGranma);
0ada222f 354 if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
355 isFromB=kTRUE;
f2703bd2 356 }
0ada222f 357 if(abspdgGranma==4 || abspdgGranma==5) isQuarkFound=kTRUE;
8a983038 358 mother = mcGranma->GetMother();
4539d15b 359 AliDebug(3, Form("mother = %d", mother));
0ada222f 360 }else{
8a983038 361 AliError("Failed casting the mother particle!");
362 break;
f2703bd2 363 }
f2703bd2 364 }
0ada222f 365
4539d15b 366 if(fRejectIfNoQuark && !isQuarkFound) {
367 return -99999;
368 }
0ada222f 369 if(isFromB){
4539d15b 370 if (!fKeepDfromB) {
371 return -9999; //skip particle if come from a B meson.
372 }
0ada222f 373 }
374 else{
4539d15b 375 if (fKeepDfromBOnly) {
376 return -999;
377 }
f036863c 378 }
379592af 379 return pdgGranma;
379592af 380}
381
379592af 382//___________________________________________
f2703bd2 383Bool_t AliCFVertexingHF::CheckMCDaughters()const
384{
385 //
386 // checking the daughters
387 // at MC level
379592af 388
f2703bd2 389 AliAODMCParticle *mcPartDaughter;
390 Bool_t checkDaughters = kFALSE;
391
392 Int_t label0 = fmcPartCandidate->GetDaughter(0);
393 Int_t label1 = fmcPartCandidate->GetDaughter(1);
4539d15b 394 AliDebug(3,Form("label0 = %d, label1 = %d", label0, label1));
04b11915 395 if (label1<=0 || label0 <= 0){
4539d15b 396 AliDebug(2, Form("The MC particle doesn't have correct daughters, skipping!!"));
397 return checkDaughters;
f2703bd2 398 }
379592af 399
6e2e4f50 400 if (fLabelArray == 0x0) {
4539d15b 401 return checkDaughters;
6e2e4f50 402 }
4539d15b 403
f2703bd2 404 for (Int_t iProng = 0; iProng<fProngs; iProng++){
4539d15b 405 AliDebug(3, Form("fLabelArray[%d] = %d", iProng, fLabelArray[iProng]));
406 }
407 AliDebug(3, Form("Entries in MC array = %d (fast = %d)", fmcArray->GetEntries(), fmcArray->GetEntriesFast()));
408 for (Int_t iProng = 0; iProng<fProngs; iProng++){
409 AliDebug(3, Form("fLabelArray[%d] = %d", iProng, fLabelArray[iProng]));
410 mcPartDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(fLabelArray[iProng]));
411 if (!mcPartDaughter) {
412 AliWarning("At least one Daughter Particle not found in tree, skipping");
413 return checkDaughters;
414 }
f2703bd2 415 }
416
417 checkDaughters = kTRUE;
418 return checkDaughters;
419}
379592af 420
421//______________________________________________________
422Bool_t AliCFVertexingHF::FillMCContainer(Double_t *containerInputMC)
423{
f2703bd2 424 //
425 // fill the container for Generator level selection
426 //
379592af 427
f2703bd2 428 Bool_t mcContainerFilled = kFALSE;
429
430 Double_t* vectorMC = new Double_t[fNVar];
431 for (Int_t iVar = 0; iVar<fNVar; iVar++) vectorMC[iVar]= 9999.;
432
433 if (GetGeneratedValuesFromMCParticle(&vectorMC[0])){
434 for (Int_t iVar = 0; iVar<fNVar; iVar++){
435 containerInputMC[iVar] = vectorMC[iVar];
436 }
437 mcContainerFilled = kTRUE;
438 }
4539d15b 439 else {
440 AliDebug(3, "We could not fill the array for the container");
441 }
f2703bd2 442 delete [] vectorMC;
443 vectorMC = 0x0;
444 return mcContainerFilled;
379592af 445}
446
447//______________________________________________________
f2703bd2 448Bool_t AliCFVertexingHF::FillRecoContainer(Double_t *containerInput)
449{
450 //
451 // fill the container for Reconstrucred level selection
452 //
453
454 Bool_t recoContainerFilled = kFALSE;
455 Double_t* vectorValues = new Double_t[fNVar];
456 Double_t* vectorReco = new Double_t[fNVar];
f2703bd2 457 for (Int_t iVar = 0; iVar<fNVar; iVar++) {
b7af2639 458
f2703bd2 459 vectorValues[iVar]= 9999.;
460 vectorReco[iVar]=9999.;
461 }
462
463 if(fFillFromGenerated){
464 //filled with MC values
465 if (GetGeneratedValuesFromMCParticle(&vectorValues[0])){
466 for (Int_t iVar = 0; iVar<fNVar; iVar++){
467 containerInput[iVar] = vectorValues[iVar];
468 }
469 recoContainerFilled = kTRUE;
379592af 470 }
471 }
472 else{
f2703bd2 473 //filled with Reco values
379592af 474
475 if (GetRecoValuesFromCandidate(&vectorReco[0])){
476 for (Int_t iVar = 0; iVar<fNVar; iVar++){
b7af2639 477 containerInput[iVar] = vectorReco[iVar];
379592af 478 }
479 recoContainerFilled = kTRUE;
480 }
481 }
f2703bd2 482
483 delete [] vectorValues;
484 delete [] vectorReco;
485 vectorValues = 0x0;
486 vectorReco = 0x0;
5cd139bc 487
f2703bd2 488 return recoContainerFilled;
379592af 489}
490
491//_____________________________________________________
492Bool_t AliCFVertexingHF::MCAcceptanceStep() const
493{
f2703bd2 494 //
495 // checking the MC acceptance step
496 //
497
498 Bool_t bMCAccStep = kFALSE;
499
500 AliAODMCParticle *mcPartDaughter;
501 Int_t label0 = fmcPartCandidate->GetDaughter(0);
502 Int_t label1 = fmcPartCandidate->GetDaughter(1);
04b11915 503 if (label1<=0 || label0 <= 0){
f2703bd2 504 AliDebug(2, Form("The MC particle doesn't have correct daughters, skipping!!"));
505 return bMCAccStep;
506 }
507
6e2e4f50 508 if (fLabelArray == 0x0) {
509 return bMCAccStep;
510 }
511
f2703bd2 512 for (Int_t iProng = 0; iProng<fProngs; iProng++){
6e2e4f50 513 mcPartDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(fLabelArray[iProng]));
f2703bd2 514 if (!mcPartDaughter) {
515 AliWarning("At least one Daughter Particle not found in tree, skipping");
516 return bMCAccStep;
517 }
518 Double_t eta = mcPartDaughter->Eta();
519 Double_t pt = mcPartDaughter->Pt();
520
521 //set values of eta and pt in the constructor.
2bf2e62b 522 // if (TMath::Abs(eta) > 0.9 || pt < 0.1){
523 if (TMath::Abs(eta) > fEtaAccCut[iProng] || pt < fPtAccCut[iProng]){
524 AliDebug(3,Form("At least one daughter has eta or pt outside the required range (|eta| = %f, pt = %f, should be |eta| < %f, pt > %f \n", TMath::Abs(eta), pt, fEtaAccCut[iProng], fPtAccCut[iProng]));
f2703bd2 525 return bMCAccStep;
526 }
527 }
528 bMCAccStep = kTRUE;
529 return bMCAccStep;
530
531}
f2703bd2 532 //_____________________________________________________
2bf2e62b 533Bool_t AliCFVertexingHF::MCRefitStep(AliAODEvent *aodEvent, AliESDtrackCuts **trackCuts) const
379592af 534{
f2703bd2 535 //
536 // check on the kTPCrefit and kITSrefit conditions of the daughters
537 //
379592af 538 Bool_t bRefitStep = kFALSE;
f2703bd2 539
379592af 540 Int_t label0 = fmcPartCandidate->GetDaughter(0);
541 Int_t label1 = fmcPartCandidate->GetDaughter(1);
f2703bd2 542
04b11915 543 if (label1<=0 || label0 <= 0){
f2703bd2 544 AliDebug(2, Form("The MC particle doesn't have correct daughters, skipping!!"));
545 return bRefitStep;
379592af 546 }
f2703bd2 547
6e2e4f50 548 if (fLabelArray == 0x0) {
549 return bRefitStep;
550 }
f2703bd2 551
379592af 552 Int_t foundDaughters = 0;
6e2e4f50 553 Int_t* temp = new Int_t[fProngs];
554 for (Int_t ilabel = 0; ilabel<fProngs; ilabel++){
555 temp[ilabel] = fLabelArray[ilabel];
556 }
557
2bf2e62b 558 // if (trackCuts->GetRequireTPCRefit() || trackCuts->GetRequireITSRefit()){
f2703bd2 559
2bf2e62b 560 for(Int_t iaod =0; iaod<aodEvent->GetNumberOfTracks(); iaod++){
561 AliAODTrack *track = (AliAODTrack*)aodEvent->GetTrack(iaod);
562 if(track->GetStatus()&AliESDtrack::kITSpureSA) continue;
563 Bool_t foundTrack = kFALSE;
564 Int_t prongindex;
565 for (Int_t ilabel = 0; ilabel<fProngs; ilabel++){
566 AliDebug(3,Form("fLabelArray[%d] = %d, track->GetLabel() = %d",ilabel,fLabelArray[ilabel],TMath::Abs(track->GetLabel())));
fbec9fa9 567 if ((track->GetLabel()<0)&&(fFakeSelection==1)) continue;
568 if ((track->GetLabel()>0)&&(fFakeSelection==2)) continue;
569
2bf2e62b 570 if (TMath::Abs(track->GetLabel()) == temp[ilabel]) {
571 foundTrack = kTRUE;
572 temp[ilabel] = 0;
573 prongindex=ilabel;
574 break;
6e2e4f50 575 }
2bf2e62b 576 }
577 if (foundTrack){
578 foundDaughters++;
579 AliDebug(4,Form("daughter %d \n",foundDaughters));
580 if(trackCuts[prongindex]->GetRequireTPCRefit()){
581 if(track->GetStatus()&AliESDtrack::kTPCrefit) {
582 bRefitStep = kTRUE;
379592af 583 }
2bf2e62b 584 else {
585 AliDebug(3, "Refit cut not passed , missing TPC refit\n");
586 delete [] temp;
587 temp = 0x0;
588 return kFALSE;
379592af 589 }
2bf2e62b 590 }
591
592 if (trackCuts[prongindex]->GetRequireITSRefit()) {
593 if(track->GetStatus()&AliESDtrack::kITSrefit){
594 bRefitStep = kTRUE;
595 }
596 else {
597 AliDebug(3, "Refit cut not passed , missing ITS refit\n");
598 delete [] temp;
599 temp = 0x0;
600 return kFALSE;
601 }
602 }
603 }
604 if (foundDaughters == fProngs){
605 break;
606 }
607 }
608 //}
6e2e4f50 609 delete [] temp;
610 temp = 0x0;
f2703bd2 611 if (foundDaughters== fProngs) return bRefitStep;
612 else return kFALSE;
379592af 613}
614
615//____________________________________________________________________________
616
617Bool_t AliCFVertexingHF::RecoStep()
379592af 618{
f2703bd2 619 //
620 //check also vertex and ITS Refit and TPC Refit
621 //
379592af 622
f2703bd2 623 Bool_t bRecoStep = kFALSE;
624 Int_t mcLabel = GetMCLabel();
625
626 if (mcLabel == -1) {
627 AliDebug(2,"No MC particle found");
628 return bRecoStep;
629 }
630 else{
631 fmcPartCandidate = (AliAODMCParticle*)fmcArray->At(mcLabel);
632 if (!fmcPartCandidate){
633 AliWarning("Could not find associated MC in AOD MC tree");
634 return bRecoStep;
635 }
636 }
637
638 Int_t pdgGranma = CheckOrigin();
639
0ada222f 640 if (pdgGranma == -99999){
4539d15b 641 AliDebug(2,"This particle does not have a quark in his genealogy\n");
0ada222f 642 return bRecoStep;
643 }
f2703bd2 644 if (pdgGranma == -9999){
4539d15b 645 AliDebug(2,"This particle come from a B decay channel but according to the settings of the task, we keep only prompt charm particles\n");
f2703bd2 646 return bRecoStep;
647 }
648
649 if (pdgGranma == -999){
650 AliDebug(2,"This particle come from a prompt charm particle but according to the settings of the task, we want only the ones coming from B\n");
651 return bRecoStep;
652 }
379592af 653
f2703bd2 654 bRecoStep=kTRUE;
655 return bRecoStep;
379592af 656}
657//____________________________________________
f2703bd2 658Double_t AliCFVertexingHF::GetEtaProng(Int_t iProng) const
659{
660 //
661 // getting eta of the prong
662 //
663
664 if (fRecoCandidate){
665 Double_t etaProng = fRecoCandidate->EtaProng(iProng);
666 return etaProng;
667 }
668 return 999999;
379592af 669}
670//______________________________________________________
f2703bd2 671Double_t AliCFVertexingHF::GetPtProng(Int_t iProng) const
672{
673 //
674 // getting pt of the prong
675 //
379592af 676
f2703bd2 677 if (fRecoCandidate){
678 Double_t ptProng = fRecoCandidate->PtProng(iProng);
679 return ptProng;
680 }
681 return 999999;
682
379592af 683}
684
685//____________________________________________________________________
686
2bf2e62b 687Bool_t AliCFVertexingHF::RecoAcceptStep(AliESDtrackCuts **trackCuts) const
379592af 688{
f2703bd2 689 //
690 // reco Acceptance step
691 //
692
693 Bool_t bRecoAccStep = kFALSE;
694
e56e56b4 695 Float_t etaCutMin=0, ptCutMin=0, etaCutMax=0, ptCutMax=0;
f2703bd2 696
697 Float_t etaProng=0., ptProng=0.;
698
699 for (Int_t iProng =0; iProng<fProngs; iProng++){
700
2bf2e62b 701 trackCuts[iProng]->GetEtaRange(etaCutMin, etaCutMax);
702 trackCuts[iProng]->GetPtRange(ptCutMin, ptCutMax);
f2703bd2 703 etaProng = GetEtaProng(iProng);
704 ptProng = GetPtProng(iProng);
705
706 Bool_t acceptanceProng = (etaProng>etaCutMin && etaProng<etaCutMax && ptProng>ptCutMin && ptProng<ptCutMax);
707 if (!acceptanceProng) {
708 AliDebug(2,"At least one reconstructed prong isn't in the acceptance\n");
709 return bRecoAccStep;
710 }
711 }
712
713 bRecoAccStep=kTRUE;
714 return bRecoAccStep;
379592af 715}
716//___________________________________________________________
717
04b11915 718Bool_t AliCFVertexingHF::FillUnfoldingMatrix(UInt_t pdg, Double_t fill[4]) const
f2703bd2 719{
720 //
721 // filling the unfolding matrix
722 //
723
f2703bd2 724 if(fmcPartCandidate){
725
726 fill[0] = GetPtCand();
04b11915 727 fill[1] = GetYCand(pdg);
f2703bd2 728
729 fill[2] = fmcPartCandidate->Pt();
730 fill[3] = fmcPartCandidate->Y();
731
732 return kTRUE;
733 }
734
f2703bd2 735 return kFALSE;
379592af 736}
f2703bd2 737//___________________________________________________________
379592af 738
3ee5eb83 739Int_t AliCFVertexingHF::CheckReflexion(Char_t isSign)
f2703bd2 740{
741 //
742 // check for reflexion (particle/antiparticle)
743 //
379592af 744
f2703bd2 745 Int_t mcLabel = GetMCLabel();
746
747 if (mcLabel == -1) {
748 AliDebug(2,"No MC particle found");
749 return 0;
750 }
751 else{
752 fmcPartCandidate = (AliAODMCParticle*)fmcArray->At(mcLabel);
753 if (!fmcPartCandidate){
754 AliWarning("Could not find associated MC in AOD MC tree");
755 return 0;
756 }
757 }
758
3ee5eb83 759 if(fmcPartCandidate->GetPdgCode()>0) {
760 if (isSign == 1){ // I ask for antiparticle only
761 AliDebug(2,"candidate is particle, I ask for antiparticle only");
762 return 0;
763 }
764 return 1; // particle
765 }
766 else if(fmcPartCandidate->GetPdgCode()<0) {
767 if (isSign == 0){ // I ask for particle only
768 AliDebug(2,"candidate is antiparticle, I ask for particle only");
769 return 0;
770 }
771 return 2; // antiparticle
772 }
f2703bd2 773 else return 0; // ....shouldn't be...
3ee5eb83 774
f2703bd2 775}
6e2e4f50 776//___________________________________________________________
777
778Bool_t AliCFVertexingHF::SetLabelArray()
779{
780 //
781 // setting the label arrays
782 //
783
784 Bool_t bLabelArray = kFALSE;
785
786 fLabelArray = new Int_t[fProngs];
787
788 AliAODMCParticle *mcPartDaughter;
789 Int_t label0 = fmcPartCandidate->GetDaughter(0);
790 Int_t label1 = fmcPartCandidate->GetDaughter(1);
bee31217 791 AliDebug(2, Form("label0 = %d, label1 = %d", label0, label1));
4539d15b 792 AliAODMCParticle* tmp0 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(label0));
793 AliAODMCParticle* tmp1 = dynamic_cast<AliAODMCParticle*>(fmcArray->At(label1));
794
795 AliDebug(2, Form("label0 = %d (pdg = %d), label1 = %d (pdg = %d)", label0, tmp0->GetPdgCode(), label1, tmp1->GetPdgCode()));
04b11915 796 if (label1<=0 || label0 <= 0){
4539d15b 797 AliDebug(2, Form("The MC particle doesn't have correct daughters, skipping!!"));
6e2e4f50 798 delete [] fLabelArray;
799 fLabelArray = 0x0;
800 return bLabelArray;
801 }
802
803 if (label1 - label0 == fProngs-1){
804 for (Int_t iProng = 0; iProng<fProngs; iProng++){
805 mcPartDaughter = dynamic_cast<AliAODMCParticle*>(fmcArray->At(label0+iProng));
8a983038 806 if (mcPartDaughter){
807 fLabelArray[iProng] = mcPartDaughter->GetLabel();
808 }
809 else{
810 AliError("Failed casting the daughter particle, returning a NULL label array");
2bf2e62b 811 delete [] fLabelArray;
812 fLabelArray = 0x0;
8a983038 813 return bLabelArray;
814 }
6e2e4f50 815 }
816
817 }
818 // resonant decay channel
819 else if (label1 - label0 == fProngs-2 && fProngs > 2){
4539d15b 820 AliDebug(3, "In the resonance decay channel");
6e2e4f50 821 Int_t labelFirstDau = fmcPartCandidate->GetDaughter(0);
822 Int_t foundDaughters = 0;
823 for(Int_t iDau=0; iDau<fProngs-1; iDau++){
824 Int_t iLabelDau = labelFirstDau+iDau;
825 AliAODMCParticle* part = dynamic_cast<AliAODMCParticle*>(fmcArray->At(iLabelDau));
661d6c29 826 if ( ! part ) {
827 AliError("Wrong particle type in fmcArray");
828 delete [] fLabelArray;
829 fLabelArray = 0x0;
830 return bLabelArray;
831 }
6e2e4f50 832 Int_t pdgCode=TMath::Abs(part->GetPdgCode());
4539d15b 833 AliDebug(3, Form("Prong %d had pdg = %d", iDau, pdgCode));
6e2e4f50 834 if(pdgCode==211 || pdgCode==321 || pdgCode==2212){
8a983038 835 if (part) {
836 fLabelArray[foundDaughters] = part->GetLabel();
4539d15b 837 AliDebug(3, Form("part found at %d has label = %d", iLabelDau, part->GetLabel()));
bee31217 838 AliDebug(3, Form("fLabelArray[%d] = %d", foundDaughters, fLabelArray[foundDaughters]));
8a983038 839 foundDaughters++;
840 }
841 else{
842 AliError("Error while casting particle! returning a NULL array");
2bf2e62b 843 delete [] fLabelArray;
844 fLabelArray = 0x0;
8a983038 845 return bLabelArray;
846 }
6e2e4f50 847 }
5cd139bc 848 // added K0S case - Start
849 else if (pdgCode==311) {
bee31217 850 AliDebug(3, Form("K0S case, foundDaughters = %d", foundDaughters));
5cd139bc 851 if (part->GetNDaughters()!=1) {
852 delete [] fLabelArray;
853 fLabelArray = 0x0;
854 return bLabelArray;
855 }
856 Int_t labelK0Dau = part->GetDaughter(0);
857 AliAODMCParticle* partK0S = dynamic_cast<AliAODMCParticle*>(fmcArray->At(labelK0Dau));
ea6f40e6 858 if(!partK0S){
859 AliError("Error while casting particle! returning a NULL array");
860 delete [] fLabelArray;
861 fLabelArray = 0x0;
862 return bLabelArray;
863 }
4539d15b 864 Int_t nDauRes = partK0S->GetNDaughters();
865 AliDebug(3, Form("nDauRes = %d", nDauRes));
866 if (nDauRes!=2 || partK0S->GetPdgCode() != 310) {
867 AliDebug(2, "No K0S on no 2-body decay");
5cd139bc 868 delete [] fLabelArray;
869 fLabelArray = 0x0;
870 return bLabelArray;
871 }
49aeca2b 872 Int_t labelFirstDauRes = partK0S->GetDaughter(0);
bee31217 873 AliDebug(2, Form("Found K0S (%d)", labelK0Dau));
5cd139bc 874 for(Int_t iDauRes=0; iDauRes<nDauRes; iDauRes++){
875 Int_t iLabelDauRes = labelFirstDauRes+iDauRes;
876 AliAODMCParticle* dauRes = dynamic_cast<AliAODMCParticle*>(fmcArray->At(iLabelDauRes));
bee31217 877 AliDebug(3, Form("daughter = %d, pointer = %p, with label = %d", iLabelDauRes, dauRes, dauRes->GetLabel()));
5cd139bc 878 if (dauRes){
4539d15b 879 AliDebug(3, Form("PDG code = %d", dauRes->GetPdgCode()));
49aeca2b 880 if (TMath::Abs(dauRes->GetPdgCode())!=211) {
b69956f1 881 AliDebug(2,"K0S doesn't decay in 2 charged pions!");
49aeca2b 882 delete [] fLabelArray;
883 fLabelArray = 0x0;
884 return bLabelArray;
885 }
886 else {
bee31217 887 fLabelArray[foundDaughters] = iLabelDauRes; // N.B.: do not use dauRes->GetLabel()!!!! it is wrong!!!
888 AliDebug(3, Form("Setting fLabelArray[%d] = %d (before it was %d)", foundDaughters, iLabelDauRes, dauRes->GetLabel()));
49aeca2b 889 foundDaughters++;
890 }
5cd139bc 891 }
892 else {
893 AliError("Error while casting resonant daughter! returning a NULL array");
894 delete [] fLabelArray;
895 fLabelArray = 0x0;
896 return bLabelArray;
897 }
898 }
899 }
900 // added K0S case - End
6e2e4f50 901 else{
4539d15b 902 Int_t nDauRes=part->GetNDaughters();
903 AliDebug(3, Form("nDauRes = %d", nDauRes));
904 if(nDauRes!=2) {
905 AliDebug(3, Form("nDauRes = %d, different from 2", nDauRes));
906 Int_t labelFirstDauResTest = part->GetDaughter(0);
907 for(Int_t iDauRes=0; iDauRes<nDauRes; iDauRes++){
908 Int_t iLabelDauResTest = labelFirstDauResTest+iDauRes;
909 AliAODMCParticle* dauRes = dynamic_cast<AliAODMCParticle*>(fmcArray->At(iLabelDauResTest));
910 if (dauRes){
911 AliDebug(3, Form("pdg of daugh %d = %d", iDauRes, dauRes->GetPdgCode()));
912 }
913 }
914 delete [] fLabelArray;
915 fLabelArray = 0x0;
916 return bLabelArray;
917 }
918 Int_t labelFirstDauRes = part->GetDaughter(0);
919 for(Int_t iDauRes=0; iDauRes<nDauRes; iDauRes++){
920 Int_t iLabelDauRes = labelFirstDauRes+iDauRes;
921 AliAODMCParticle* dauRes = dynamic_cast<AliAODMCParticle*>(fmcArray->At(iLabelDauRes));
922 if (dauRes){
923 fLabelArray[foundDaughters] = dauRes->GetLabel();
924 foundDaughters++;
925 }
926 else{
927 AliError("Error while casting resonant daughter! returning a NULL array");
928 delete [] fLabelArray;
929 fLabelArray = 0x0;
930 return bLabelArray;
931 }
932 }
6e2e4f50 933 }
934 }
935 if (foundDaughters != fProngs){
4539d15b 936 AliDebug(3, Form("foundDaughters (%d) != fProngs (%d)", foundDaughters, fProngs));
937 delete [] fLabelArray;
938 fLabelArray = 0x0;
939 return bLabelArray;
6e2e4f50 940 }
941 }
942 // wrong correspondance label <--> prongs
943 else{
4539d15b 944 delete [] fLabelArray;
945 fLabelArray = 0x0;
946 return bLabelArray;
947 }
948 AliDebug(3, "Setting AccCuts");
949 SetAccCut(); // setting the pt and eta acceptance cuts
950 bLabelArray = kTRUE;
951 return bLabelArray;
6e2e4f50 952}
2bf2e62b 953
954//___________________________________________________________
955
956void AliCFVertexingHF::SetPtAccCut(Float_t* ptAccCut)
957{
958 //
959 // setting the pt cut to be used in the Acceptance steps (MC+Reco)
960 //
961
962 if (fProngs>0){
963 for (Int_t iP=0; iP<fProngs; iP++){
964 fPtAccCut[iP]=ptAccCut[iP];
965 }
966 }
967 return;
968}
969
970
971
972//___________________________________________________________
973
974void AliCFVertexingHF::SetEtaAccCut(Float_t* etaAccCut)
975{
976 //
977 // setting the eta cut to be used in the Acceptance steps (MC+Reco)
978 //
979
980 if (fProngs>0){
981 for (Int_t iP=0; iP<fProngs; iP++){
982 fEtaAccCut[iP]=etaAccCut[iP];
983 }
984 }
985 return;
986}
987//___________________________________________________________
988
989void AliCFVertexingHF::SetAccCut(Float_t* ptAccCut, Float_t* etaAccCut)
990{
991 //
992 // setting the pt and eta cut to be used in the Acceptance steps (MC+Reco)
993 //
994
995 if (fProngs>0){
996 for (Int_t iP=0; iP<fProngs; iP++){
997 fPtAccCut[iP]=ptAccCut[iP];
998 fEtaAccCut[iP]=etaAccCut[iP];
999 }
1000 }
1001 return;
1002}
1003
1004//___________________________________________________________
1005
1006void AliCFVertexingHF::SetAccCut()
1007{
1008 //
1009 // setting the pt and eta cut to be used in the Acceptance steps (MC+Reco)
1010 //
1011
1012 if (fProngs>0){
1013 for (Int_t iP=0; iP<fProngs; iP++){
1014 fPtAccCut[iP]=0.1;
1015 fEtaAccCut[iP]=0.9;
1016 }
1017 }
1018 return;
1019}