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