]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG3/vertexingHF/AliCFVertexingHF.cxx
Modifications in the filtering task: 1) add switch for 3 prong LS candidates; 2)...
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliCFVertexingHF.cxx
index 3c88fae77b7994e562c407ae3f16871c70c5b02b..4c01793707f7b3a87a3583a5ad820eea9cc9594b 100644 (file)
@@ -13,6 +13,8 @@
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
 
+// $Id$
+
 //-----------------------------------------------------------------------
 // Class for HF corrections as a function of many variables and step 
 // Author : C. Zampolli, CERN
@@ -27,6 +29,7 @@
 #include "AliAODRecoDecayHF2Prong.h"
 #include "AliAODRecoDecayHF3Prong.h"
 #include "AliAODRecoDecayHF4Prong.h"
+#include "AliAODRecoCascadeHF.h"
 #include "AliAODMCHeader.h"
 #include "AliAODEvent.h"
 #include "AliLog.h"
@@ -50,7 +53,13 @@ AliCFVertexingHF::AliCFVertexingHF() :
        fKeepDfromBOnly(kFALSE),
        fmcLabel(0),
        fProngs(-1),
-       fLabelArray(0x0)
+       fLabelArray(0x0), 
+       fCentValue(0.),
+       fPtAccCut(0x0),
+       fEtaAccCut(0x0),
+       fFakeSelection(0),
+       fFake(1.), // setting to MC value
+       fRejectIfNoQuark(kFALSE)
 {
        //
        // constructor
@@ -77,7 +86,13 @@ AliCFVertexingHF::AliCFVertexingHF(TClonesArray *mcArray, UShort_t originDselect
        fKeepDfromBOnly(kFALSE),
        fmcLabel(0),
        fProngs(-1),
-       fLabelArray(0x0)
+       fLabelArray(0x0),
+       fCentValue(0.),
+       fPtAccCut(0x0),
+       fEtaAccCut(0x0),
+       fFakeSelection(0),
+       fFake(1.), // setting to MC value
+       fRejectIfNoQuark(kFALSE)
 {
        //
        // constructor with mcArray
@@ -98,8 +113,16 @@ AliCFVertexingHF::~AliCFVertexingHF()
        if (fRecoCandidate) fRecoCandidate = 0x0;
        if (fmcPartCandidate) fmcPartCandidate = 0x0;
        if (fLabelArray){
-               delete [] fLabelArray;
-               fLabelArray = 0x0;
+               delete [] fLabelArray;
+               fLabelArray = 0x0;
+       }       
+       if (fPtAccCut){
+               delete [] fPtAccCut;
+               fPtAccCut = 0x0;
+       }       
+       if (fEtaAccCut){
+               delete [] fEtaAccCut;
+               fEtaAccCut = 0x0;
        }       
 }
 
@@ -125,9 +148,19 @@ AliCFVertexingHF& AliCFVertexingHF::operator=(const AliCFVertexingHF& c)
                fKeepDfromBOnly = c.fKeepDfromBOnly;
                fmcLabel = c.fmcLabel;
                fProngs=c.fProngs;
+               fCentValue=c.fCentValue;
+               fFakeSelection=c.fFakeSelection;
+               fFake=c.fFake;
+               fRejectIfNoQuark=c.fRejectIfNoQuark;
                if (fProngs > 0){
                        fLabelArray = new Int_t[fProngs];
-                       for(Int_t iP=0; iP<fProngs; iP++)fLabelArray[iP]=c.fLabelArray[iP];
+                        fPtAccCut = new Float_t[fProngs];
+                        fEtaAccCut = new Float_t[fProngs];
+                       for(Int_t iP=0; iP<fProngs; iP++){
+                               fLabelArray[iP]=c.fLabelArray[iP];
+                               fPtAccCut[iP]=c.fPtAccCut[iP];
+                               fEtaAccCut[iP]=c.fEtaAccCut[iP];
+                       }
                }
        }
        
@@ -150,14 +183,24 @@ AliCFVertexingHF::AliCFVertexingHF(const AliCFVertexingHF &c) :
        fKeepDfromBOnly (c.fKeepDfromBOnly),
        fmcLabel(c.fmcLabel),
        fProngs(c.fProngs),
-       fLabelArray(0x0)
+       fLabelArray(0x0),
+       fCentValue(c.fCentValue),
+       fPtAccCut(0x0),
+       fEtaAccCut(0x0),
+       fFakeSelection(c.fFakeSelection),
+       fFake(c.fFake),
+       fRejectIfNoQuark(c.fRejectIfNoQuark)
 {  
        //
        //copy constructor
        //
        if (fProngs > 0){
                fLabelArray = new Int_t[fProngs];
+                fPtAccCut = new Float_t[fProngs];
+                fEtaAccCut = new Float_t[fProngs];
                if (c.fLabelArray) memcpy(fLabelArray,c.fLabelArray,fProngs*sizeof(Int_t));
+               if (c.fPtAccCut) memcpy(fPtAccCut,c.fPtAccCut,fProngs*sizeof(Int_t));
+               if (c.fEtaAccCut) memcpy(fEtaAccCut,c.fEtaAccCut,fProngs*sizeof(Int_t));
        }
 }
 
@@ -234,6 +277,11 @@ Bool_t AliCFVertexingHF::CheckMCPartFamily(AliAODMCParticle */*mcPart*/, TClones
        //
 
        Int_t pdgGranma = CheckOrigin();
+
+       if (pdgGranma == -99999){
+               AliDebug(2,"This particle does not have a quark in his genealogy\n");
+               return kFALSE;
+       }
        if (pdgGranma == -9999){
                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");   
                return kFALSE;
@@ -245,8 +293,8 @@ Bool_t AliCFVertexingHF::CheckMCPartFamily(AliAODMCParticle */*mcPart*/, TClones
        }       
        
        if (!CheckMCDaughters()) {
-               AliDebug(3, "CheckMCDaughters false");
-               return kFALSE;
+         AliDebug(3, "CheckMCDaughters false");
+         return kFALSE;
        }
        if (!CheckMCChannelDecay()) {
                AliDebug(3,"CheckMCChannelDecay false");
@@ -267,6 +315,8 @@ Int_t AliCFVertexingHF::CheckOrigin() const
        mother = fmcPartCandidate->GetMother();
        Int_t istep = 0;
        Int_t abspdgGranma =0;
+       Bool_t isFromB=kFALSE;
+       Bool_t isQuarkFound=kFALSE;
        while (mother >0 ){
                istep++;
                AliDebug(2,Form("mother at step %d = %d", istep, mother));
@@ -275,22 +325,23 @@ Int_t AliCFVertexingHF::CheckOrigin() const
                        pdgGranma = mcGranma->GetPdgCode();
                        AliDebug(2,Form("Pdg mother at step %d = %d", istep, pdgGranma));
                        abspdgGranma = TMath::Abs(pdgGranma);
-                       if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)) {
-                               if (!fKeepDfromB) return -9999; //skip particle if come from a B meson.
-                               
-                               else{
-                                       break;
-                               }
+                       if ((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000)){
+                         isFromB=kTRUE;
                        }
+                       if(abspdgGranma==4 || abspdgGranma==5) isQuarkFound=kTRUE;
                        mother = mcGranma->GetMother();
-               }
-               else {
+               }else{
                        AliError("Failed casting the mother particle!");
                        break;
                }
        }
-       if (!((abspdgGranma > 500 && abspdgGranma < 600) || (abspdgGranma > 5000 && abspdgGranma < 6000))){
-               if (fKeepDfromBOnly) return -999;
+
+       if(fRejectIfNoQuark && !isQuarkFound) return -99999;
+       if(isFromB){
+         if (!fKeepDfromB) return -9999; //skip particle if come from a B meson.
+       }
+       else{
+         if (fKeepDfromBOnly) return -999;
        }
        return pdgGranma;
 }
@@ -362,8 +413,8 @@ Bool_t AliCFVertexingHF::FillRecoContainer(Double_t *containerInput)
        Bool_t recoContainerFilled = kFALSE;
        Double_t* vectorValues = new Double_t[fNVar];
        Double_t* vectorReco = new Double_t[fNVar];  
-       
        for (Int_t iVar = 0; iVar<fNVar; iVar++) {
+
                vectorValues[iVar]= 9999.;
                vectorReco[iVar]=9999.;
        }
@@ -382,7 +433,7 @@ Bool_t AliCFVertexingHF::FillRecoContainer(Double_t *containerInput)
                
                if (GetRecoValuesFromCandidate(&vectorReco[0])){
                        for (Int_t iVar = 0; iVar<fNVar; iVar++){
-                               containerInput[iVar] = vectorReco[iVar];
+                         containerInput[iVar] = vectorReco[iVar];
                        }
                        recoContainerFilled = kTRUE;            
                }
@@ -426,8 +477,9 @@ Bool_t AliCFVertexingHF::MCAcceptanceStep() const
                Double_t pt = mcPartDaughter->Pt();
                
                //set values of eta and pt in the constructor.
-               if (TMath::Abs(eta) > 0.9 || pt < 0.1){
-                       AliDebug(3,"At least one daughter has eta>0.9 or pt < 0.1 \n"); 
+               //              if (TMath::Abs(eta) > 0.9 || pt < 0.1){
+               if (TMath::Abs(eta) > fEtaAccCut[iProng] || pt < fPtAccCut[iProng]){
+                       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])); 
                        return bMCAccStep;
                }
        }  
@@ -436,7 +488,7 @@ Bool_t AliCFVertexingHF::MCAcceptanceStep() const
        
 }
  //_____________________________________________________
-Bool_t AliCFVertexingHF::MCRefitStep(AliAODEvent *aodEvent, AliESDtrackCuts *trackCuts) const
+Bool_t AliCFVertexingHF::MCRefitStep(AliAODEvent *aodEvent, AliESDtrackCuts **trackCuts) const
 {              
        //
        // check on the kTPCrefit and kITSrefit conditions of the daughters
@@ -461,52 +513,57 @@ Bool_t AliCFVertexingHF::MCRefitStep(AliAODEvent *aodEvent, AliESDtrackCuts *tra
                temp[ilabel] = fLabelArray[ilabel];
        }
 
-       if (trackCuts->GetRequireTPCRefit() || trackCuts->GetRequireITSRefit()){
+       //      if (trackCuts->GetRequireTPCRefit() || trackCuts->GetRequireITSRefit()){
                
-               for(Int_t iaod =0; iaod<aodEvent->GetNumberOfTracks(); iaod++){
-                       AliAODTrack *track = (AliAODTrack*)aodEvent->GetTrack(iaod);
-                       if(track->GetStatus()&AliESDtrack::kITSpureSA) continue;
-                       Bool_t foundTrack = kFALSE;
-                       for (Int_t ilabel = 0; ilabel<fProngs; ilabel++){
-                               AliDebug(3,Form("fLabelArray[%d] = %d, track->GetLabel() = %d",ilabel,fLabelArray[ilabel],track->GetLabel()));
-                               if (track->GetLabel() == temp[ilabel]) {
-                                       foundTrack = kTRUE;
-                                       temp[ilabel] = 0;
-                                       break;
+       for(Int_t iaod =0; iaod<aodEvent->GetNumberOfTracks(); iaod++){
+               AliAODTrack *track = (AliAODTrack*)aodEvent->GetTrack(iaod);
+               if(track->GetStatus()&AliESDtrack::kITSpureSA) continue;
+               Bool_t foundTrack = kFALSE;
+               Int_t prongindex;
+               for (Int_t ilabel = 0; ilabel<fProngs; ilabel++){
+                       AliDebug(3,Form("fLabelArray[%d] = %d, track->GetLabel() = %d",ilabel,fLabelArray[ilabel],TMath::Abs(track->GetLabel())));
+                       if ((track->GetLabel()<0)&&(fFakeSelection==1)) continue;
+                       if ((track->GetLabel()>0)&&(fFakeSelection==2)) continue;
+
+                       if (TMath::Abs(track->GetLabel()) == temp[ilabel]) {
+                               foundTrack = kTRUE;
+                               temp[ilabel] = 0;
+                               prongindex=ilabel;
+                               break;
+                       }
+               }
+               if (foundTrack){
+                       foundDaughters++;
+                       AliDebug(4,Form("daughter %d \n",foundDaughters));
+                       if(trackCuts[prongindex]->GetRequireTPCRefit()){
+                               if(track->GetStatus()&AliESDtrack::kTPCrefit) {
+                                       bRefitStep = kTRUE;
+                               }
+                               else {
+                                       AliDebug(3, "Refit cut not passed , missing TPC refit\n");
+                                       delete [] temp;
+                                       temp = 0x0;
+                                       return kFALSE;
                                }
                        }
-                       if (foundTrack){
-                               foundDaughters++;
-                               AliDebug(4,Form("daughter %d \n",foundDaughters));
-                               if(trackCuts->GetRequireTPCRefit()){
-                                       if(track->GetStatus()&AliESDtrack::kTPCrefit) {
-                                               bRefitStep = kTRUE;
-                                       }
-                                       else {
-                                               AliDebug(3, "Refit cut not passed , missing TPC refit\n");
-                                               delete [] temp;
-                                               temp = 0x0;
-                                               return kFALSE;
-                                       }
+                       
+                       if (trackCuts[prongindex]->GetRequireITSRefit()) {
+                               if(track->GetStatus()&AliESDtrack::kITSrefit){
+                                       bRefitStep = kTRUE;
                                }
-                               
-                               if (trackCuts->GetRequireITSRefit()) {
-                                       if(track->GetStatus()&AliESDtrack::kITSrefit){
-                                               bRefitStep = kTRUE;
-                                       }
-                                       else {
-                                               AliDebug(3, "Refit cut not passed , missing ITS refit\n");
-                                               delete [] temp;
-                                               temp = 0x0;
-                                               return kFALSE;
-                                       }
+                               else {
+                                       AliDebug(3, "Refit cut not passed , missing ITS refit\n");
+                                       delete [] temp;
+                                       temp = 0x0;
+                                       return kFALSE;
                                }
-                       }       
-                       if (foundDaughters == fProngs){                         
-                               break;
-                       }                       
-               }    
-       }                                               
+                       }
+               }       
+               if (foundDaughters == fProngs){                         
+                       break;
+               }                       
+       }    
+       //}                                             
        delete [] temp;
        temp = 0x0;
        if (foundDaughters== fProngs)  return bRefitStep;
@@ -538,6 +595,10 @@ Bool_t AliCFVertexingHF::RecoStep()
        
        Int_t pdgGranma = CheckOrigin();
        
+       if (pdgGranma == -99999){
+               AliDebug(2,"This particle does not have a quark in his genealogy\n");
+               return bRecoStep;
+       }
        if (pdgGranma == -9999){
                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");
                return bRecoStep;
@@ -581,7 +642,7 @@ Double_t AliCFVertexingHF::GetPtProng(Int_t iProng) const
 
 //____________________________________________________________________
 
-Bool_t AliCFVertexingHF::RecoAcceptStep(AliESDtrackCuts *trackCuts) const
+Bool_t AliCFVertexingHF::RecoAcceptStep(AliESDtrackCuts **trackCuts) const
 {
        //
        // reco Acceptance step
@@ -590,13 +651,13 @@ Bool_t AliCFVertexingHF::RecoAcceptStep(AliESDtrackCuts *trackCuts) const
        Bool_t bRecoAccStep = kFALSE;
        
        Float_t etaCutMin, ptCutMin, etaCutMax, ptCutMax;
-       trackCuts->GetEtaRange(etaCutMin, etaCutMax);
-       trackCuts->GetPtRange(ptCutMin, ptCutMax);
        
        Float_t etaProng=0., ptProng=0.; 
        
        for (Int_t iProng =0; iProng<fProngs; iProng++){
                
+               trackCuts[iProng]->GetEtaRange(etaCutMin, etaCutMax);
+               trackCuts[iProng]->GetPtRange(ptCutMin, ptCutMax);
                etaProng = GetEtaProng(iProng);
                ptProng = GetPtProng(iProng);
                
@@ -715,6 +776,12 @@ Bool_t AliCFVertexingHF::SetLabelArray()
                for(Int_t iDau=0; iDau<fProngs-1; iDau++){
                        Int_t iLabelDau = labelFirstDau+iDau;
                        AliAODMCParticle* part = dynamic_cast<AliAODMCParticle*>(fmcArray->At(iLabelDau));
+                        if ( ! part ) {
+                          AliError("Wrong particle type in fmcArray");
+                         delete [] fLabelArray; 
+                         fLabelArray = 0x0; 
+                          return bLabelArray;
+                        }  
                        Int_t pdgCode=TMath::Abs(part->GetPdgCode());
                        if(pdgCode==211 || pdgCode==321 || pdgCode==2212){
                                if (part) {
@@ -764,6 +831,74 @@ Bool_t AliCFVertexingHF::SetLabelArray()
                fLabelArray = 0x0;  
                return bLabelArray;
        }
+       SetAccCut();   // setting the pt and eta acceptance cuts
        bLabelArray = kTRUE;
        return bLabelArray;
 }
+
+//___________________________________________________________
+
+void AliCFVertexingHF::SetPtAccCut(Float_t* ptAccCut)
+{
+       //
+       // setting the pt cut to be used in the Acceptance steps (MC+Reco)
+       //
+
+       if (fProngs>0){
+               for (Int_t iP=0; iP<fProngs; iP++){
+                       fPtAccCut[iP]=ptAccCut[iP];
+               }
+       }
+       return;
+}              
+
+
+
+//___________________________________________________________
+
+void AliCFVertexingHF::SetEtaAccCut(Float_t* etaAccCut)
+{
+       //
+       // setting the eta cut to be used in the Acceptance steps (MC+Reco)
+       //
+
+       if (fProngs>0){
+               for (Int_t iP=0; iP<fProngs; iP++){
+                       fEtaAccCut[iP]=etaAccCut[iP];
+               }
+       }
+       return;
+}      
+//___________________________________________________________
+
+void AliCFVertexingHF::SetAccCut(Float_t* ptAccCut, Float_t* etaAccCut)
+{
+       //
+       // setting the pt and eta cut to be used in the Acceptance steps (MC+Reco)
+       //
+
+       if (fProngs>0){
+               for (Int_t iP=0; iP<fProngs; iP++){
+                       fPtAccCut[iP]=ptAccCut[iP];
+                       fEtaAccCut[iP]=etaAccCut[iP];
+               }
+       }
+       return;
+}              
+
+//___________________________________________________________
+
+void AliCFVertexingHF::SetAccCut()
+{
+       //
+       // setting the pt and eta cut to be used in the Acceptance steps (MC+Reco)
+       //
+
+       if (fProngs>0){
+               for (Int_t iP=0; iP<fProngs; iP++){
+                       fPtAccCut[iP]=0.1;
+                       fEtaAccCut[iP]=0.9;
+               }
+       }
+       return;
+}