]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWGHF/vertexingHF/AliAnalysisTaskSELc2V0bachelorTMVA.cxx
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliAnalysisTaskSELc2V0bachelorTMVA.cxx
old mode 100644 (file)
new mode 100755 (executable)
index e4ce4cb..000a1fd
@@ -186,7 +186,12 @@ AliAnalysisTaskSE(),
   fCutKFDeviationFromVtx(999999.),
   fCutKFDeviationFromVtxV0(0.),
   fCurrentEvent(-1),
-  fBField(0)
+  fBField(0),
+  fKeepingOnlyPYTHIABkg(kFALSE),
+  fHistoMCLcK0SpGen(0x0),
+  fHistoMCLcK0SpGenAcc(0x0),
+  fHistoMCLcK0SpGenLimAcc(0x0),
+  fTriggerMask(0)
 {
   //
   // Default ctor
@@ -305,8 +310,12 @@ AliAnalysisTaskSELc2V0bachelorTMVA::AliAnalysisTaskSELc2V0bachelorTMVA(const Cha
   fCutKFDeviationFromVtx(999999.),
   fCutKFDeviationFromVtxV0(0.),
   fCurrentEvent(-1),
-  fBField(0)
-
+  fBField(0),
+  fKeepingOnlyPYTHIABkg(kFALSE),
+  fHistoMCLcK0SpGen(0x0),
+  fHistoMCLcK0SpGenAcc(0x0),
+  fHistoMCLcK0SpGenLimAcc(0x0),
+  fTriggerMask(0)
 {
   //
   // Constructor. Initialization of Inputs and Outputs
@@ -393,7 +402,7 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::Init() {
   fListCuts->Add(new AliRDHFCutsLctoV0(*fAnalCuts));
   PostData(3,fListCuts);
 
-  if (fUseMCInfo && fKeepingOnlyHIJINGBkg) fUtils = new AliVertexingHFUtils();
+  if (fUseMCInfo && (fKeepingOnlyHIJINGBkg || fKeepingOnlyPYTHIABkg)) fUtils = new AliVertexingHFUtils();
 
   return;
 }
@@ -405,7 +414,7 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::Terminate(Option_t*)
   // a query. It always runs on the client, it can be used to present
   // the results graphically or save the results to file.
   
-  //AliInfo("Terminate","");
+  AliInfo("Terminate");
   AliAnalysisTaskSE::Terminate();
   
   fOutput = dynamic_cast<TList*> (GetOutputData(1));
@@ -413,6 +422,27 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::Terminate(Option_t*)
     AliError("fOutput not available");
     return;
   }
+
+  
+  //AliDebug(2, Form("At MC level, %f Lc --> K0S + p were found", fHistoMCLcK0SpGen->GetEntries()));
+  //AliDebug(2, Form("At MC level, %f Lc --> K0S + p were found in the acceptance", fHistoMCLcK0SpGenAcc->GetEntries()));
+  //AliDebug(2, Form("At Reco level, %lld Lc --> K0S + p were found", fVariablesTreeSgn->GetEntries()));
+  if(fHistoMCLcK0SpGen) {
+    AliInfo(Form("At MC level, %f Lc --> K0S + p were found", fHistoMCLcK0SpGen->GetEntries()));
+  } else {
+    AliInfo("fHistoMCLcK0SpGen not available");
+  }
+  if(fHistoMCLcK0SpGenAcc) {
+    AliInfo(Form("At MC level, %f Lc --> K0S + p were found in the acceptance", fHistoMCLcK0SpGenAcc->GetEntries()));
+  } else {
+    AliInfo("fHistoMCLcK0SpGenAcc not available");
+  }
+  if(fVariablesTreeSgn) {
+    AliInfo(Form("At Reco level, %lld Lc --> K0S + p were found", fVariablesTreeSgn->GetEntries()));
+  } else {
+    AliInfo("fVariablesTreeSgn not available");
+  }    
+  
   fOutputKF = dynamic_cast<TList*> (GetOutputData(6));
   if (!fOutputKF) {     
     AliError("fOutputKF not available");
@@ -467,8 +497,8 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::UserCreateOutputObjects() {
   fCandidateVariableNames[24]="LcEta";
   fCandidateVariableNames[25]="V0positiveEta";
   fCandidateVariableNames[26]="V0negativeEta";
-  fCandidateVariableNames[27]="combinedPionProb";
-  fCandidateVariableNames[28]="combinedKaonProb";
+  fCandidateVariableNames[27]="TPCProtonProb";
+  fCandidateVariableNames[28]="TOFProtonProb";
   fCandidateVariableNames[29]="bachelorEta";
   fCandidateVariableNames[30]="LcP";
   fCandidateVariableNames[31]="bachelorP";
@@ -594,8 +624,8 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::UserCreateOutputObjects() {
     fHistoLcpKpiBeforeCuts->GetXaxis()->SetBinLabel(ibin, labelBeforeCuts[ibin-1].Data());
   }
 
-  fHistoBackground = new TH1F("fHistoBackground", "fHistoBackground", 2, -0.5, 1.5);
-  TString labelBkg[2] = {"Injected", "Non-injected"};
+  fHistoBackground = new TH1F("fHistoBackground", "fHistoBackground", 4, -0.5, 3.5);
+  TString labelBkg[4] = {"Injected", "Non-injected", "Non-PYTHIA", "PYTHIA"};
   for (Int_t ibin = 1; ibin <= fHistoBackground->GetNbinsX(); ibin++){
     fHistoBackground->GetXaxis()->SetBinLabel(ibin, labelBkg[ibin-1].Data());
   }
@@ -603,6 +633,12 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::UserCreateOutputObjects() {
   //fOutput->Add(fVariablesTreeSgn);
   //fOutput->Add(fVariablesTreeBkg);
 
+  const Float_t ptbins[15] = {0., 1., 2., 3., 4., 5., 6., 7., 8., 9., 10., 12., 17., 25., 35.};
+
+  fHistoMCLcK0SpGen = new TH1F("fHistoMCLcK0SpGen", "fHistoMCLcK0SpGen", 14, ptbins);
+  fHistoMCLcK0SpGenAcc = new TH1F("fHistoMCLcK0SpGenAcc", "fHistoMCLcK0SpGenAcc", 14, ptbins);
+  fHistoMCLcK0SpGenLimAcc = new TH1F("fHistoMCLcK0SpGenLimAcc", "fHistoMCLcK0SpGenLimAcc", 14, ptbins);
+
   fOutput->Add(fHistoEvents);
   fOutput->Add(fHistoLc);
   fOutput->Add(fHistoLcOnTheFly);
@@ -612,6 +648,9 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::UserCreateOutputObjects() {
   fOutput->Add(fHistoCodesBkg);
   fOutput->Add(fHistoLcpKpiBeforeCuts);
   fOutput->Add(fHistoBackground);
+  fOutput->Add(fHistoMCLcK0SpGen);
+  fOutput->Add(fHistoMCLcK0SpGenAcc);
+  fOutput->Add(fHistoMCLcK0SpGenLimAcc);
 
   PostData(1, fOutput);
   PostData(4, fVariablesTreeSgn);
@@ -822,7 +861,7 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::UserExec(Option_t *)
   }
 
   fCurrentEvent++;
-  Printf("Processing event = %d", fCurrentEvent);
+  AliDebug(2, Form("Processing event = %d", fCurrentEvent));
   AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);
   TClonesArray *arrayLctopKos=0;
 
@@ -852,34 +891,23 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::UserExec(Option_t *)
   
   if ( !fUseMCInfo && fIspA) {
     fAnalCuts->SetTriggerClass("");
-    fAnalCuts->SetTriggerMask(AliVEvent::kINT7);
+    fAnalCuts->SetTriggerMask(fTriggerMask);
   }
-
-  fCounter->StoreEvent(aodEvent,fAnalCuts,fUseMCInfo);
-
-  // AOD primary vertex
-  fVtx1 = (AliAODVertex*)aodEvent->GetPrimaryVertex();
-  if (!fVtx1) return;
-  if (fVtx1->GetNContributors()<1) return;
-
-  fIsEventSelected = fAnalCuts->IsEventSelected(aodEvent);
-
-  if ( !fIsEventSelected ) {
-    fHistoEvents->Fill(0);
-    return; // don't take into account not selected events 
+  
+  Int_t runnumber = aodEvent->GetRunNumber();
+  if (aodEvent->GetTriggerMask() == 0 && (runnumber >= 195344 && runnumber <= 195677)){
+    AliDebug(3,"Event rejected because of null trigger mask");
+    return;
   }
-  fHistoEvents->Fill(1);
-
-  // Setting magnetic field for KF vertexing
-  fBField = aodEvent->GetMagneticField();
-  AliKFParticle::SetField(fBField);
-
+  
+  fCounter->StoreEvent(aodEvent,fAnalCuts,fUseMCInfo);
+  
   // mc analysis 
   TClonesArray *mcArray = 0;
   AliAODMCHeader *mcHeader=0;
 
   if (fUseMCInfo) {
-    // MC array need for maching
+    // MC array need for matching
     mcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));
     if (!mcArray) {
       AliError("Could not find Monte-Carlo in AOD");
@@ -891,7 +919,34 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::UserExec(Option_t *)
       AliError("AliAnalysisTaskSELc2V0bachelorTMVA::UserExec: MC header branch not found!\n");
       return;
     }
+
+    Double_t zMCVertex = mcHeader->GetVtxZ();
+    if (TMath::Abs(zMCVertex) > fAnalCuts->GetMaxVtxZ()){
+      AliDebug(3,Form("z coordinate of MC vertex = %f, it was required to be within [-%f, +%f], skipping event", zMCVertex, fAnalCuts->GetMaxVtxZ(), fAnalCuts->GetMaxVtxZ()));
+      AliInfo(Form("z coordinate of MC vertex = %f, it was required to be within [-%f, +%f], skipping event", zMCVertex, fAnalCuts->GetMaxVtxZ(), fAnalCuts->GetMaxVtxZ()));
+      return;
+    }
+    
+    //Printf("Filling MC histo");
+    FillMCHisto(mcArray);
   }
+  
+  // AOD primary vertex
+  fVtx1 = (AliAODVertex*)aodEvent->GetPrimaryVertex();
+  if (!fVtx1) return;
+  if (fVtx1->GetNContributors()<1) return;
+
+  fIsEventSelected = fAnalCuts->IsEventSelected(aodEvent);
+
+  if ( !fIsEventSelected ) {
+    fHistoEvents->Fill(0);
+    return; // don't take into account not selected events 
+  }
+  fHistoEvents->Fill(1);
+
+  // Setting magnetic field for KF vertexing
+  fBField = aodEvent->GetMagneticField();
+  AliKFParticle::SetField(fBField);
 
   Int_t nSelectedAnal = 0;
   if (fIsK0sAnalysis) {
@@ -907,6 +962,112 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::UserExec(Option_t *)
   PostData(5, fVariablesTreeBkg);
   PostData(6, fOutputKF);
 
+}
+//-------------------------------------------------------------------------------
+void AliAnalysisTaskSELc2V0bachelorTMVA::FillMCHisto(TClonesArray *mcArray){
+
+  // method to fill MC histo: how many Lc --> K0S + p are there at MC level
+  for (Int_t iPart=0; iPart<mcArray->GetEntriesFast(); iPart++) { 
+    AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(mcArray->At(iPart));
+    if (!mcPart){
+      AliError("Failed casting particle from MC array!, Skipping particle");
+      continue;
+    }
+    Int_t pdg = mcPart->GetPdgCode();
+    if (TMath::Abs(pdg) != 4122){
+      AliDebug(2, Form("MC particle %d is not a Lc: its pdg code is %d", iPart, pdg));
+      continue;
+    }
+    AliDebug(2, Form("Step 0 ok: MC particle %d is a Lc: its pdg code is %d", iPart, pdg));
+    Int_t labeldaugh0 = mcPart->GetDaughter(0);
+    Int_t labeldaugh1 = mcPart->GetDaughter(1);
+    if (labeldaugh0 <= 0 || labeldaugh1 <= 0){
+      AliDebug(2, Form("The MC particle doesn't have correct daughters, skipping!!"));
+      continue;
+    }
+    else if (labeldaugh1 - labeldaugh0 == 1){
+      AliDebug(2, Form("Step 1 ok: The MC particle has correct daughters!!"));
+      AliAODMCParticle* daugh0 = dynamic_cast<AliAODMCParticle*>(mcArray->At(labeldaugh0));
+      AliAODMCParticle* daugh1 = dynamic_cast<AliAODMCParticle*>(mcArray->At(labeldaugh1));
+      if(!daugh0 || !daugh1){
+       AliDebug(2,"Particle daughters not properly retrieved!");
+       return;
+      }
+      Int_t pdgCodeDaugh0 = TMath::Abs(daugh0->GetPdgCode());
+      Int_t pdgCodeDaugh1 = TMath::Abs(daugh1->GetPdgCode());
+      AliAODMCParticle* bachelorMC = daugh0;
+      AliAODMCParticle* v0MC = daugh1;
+      AliDebug(2, Form("pdgCodeDaugh0 = %d, pdgCodeDaugh1 = %d", pdgCodeDaugh0, pdgCodeDaugh1));
+      if ((pdgCodeDaugh0 == 311 && pdgCodeDaugh1 == 2212) || (pdgCodeDaugh0 == 2212 && pdgCodeDaugh1 == 311)){ 
+       // we are in the case of Lc --> K0 + p; now we have to check if the K0 decays in K0S, and if this goes in pi+pi-
+       /// first, we set the bachelor and the v0: above we assumed first proton and second V0, but we could have to change it:
+       if (pdgCodeDaugh0 == 311 && pdgCodeDaugh1 == 2212) {
+         bachelorMC = daugh1;
+         v0MC = daugh0;
+       }
+       AliDebug(2, Form("Number of Daughters of v0 = %d", v0MC->GetNDaughters()));
+       if (v0MC->GetNDaughters() != 1) { 
+         AliDebug(2, "The K0 does not decay in 1 body only! Impossible... Continuing...");
+         continue;
+       }
+       else { // So far: Lc --> K0 + p, K0 with 1 daughter 
+         AliDebug(2, "Step 2 ok: The K0 does decay in 1 body only! ");
+         Int_t labelK0daugh = v0MC->GetDaughter(0);
+         AliAODMCParticle* partK0S = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelK0daugh));
+         if(!partK0S){
+           AliError("Error while casting particle! returning a NULL array");
+           continue;
+         }
+         else { // So far: Lc --> K0 + p, K0 with 1 daughter that we can access
+           if (partK0S->GetNDaughters() != 2 || TMath::Abs(partK0S->GetPdgCode() != 310)){
+             AliDebug(2, "The K0 daughter is not a K0S or does not decay in 2 bodies");
+             continue;
+           }
+           else { // So far: Lc --> K0 + p, K0 --> K0S, K0S in 2 bodies
+             AliDebug(2, "Step 3 ok: The K0 daughter is a K0S and does decay in 2 bodies");
+             Int_t labelK0Sdaugh0 = partK0S->GetDaughter(0);
+             Int_t labelK0Sdaugh1 = partK0S->GetDaughter(1);
+             AliAODMCParticle* daughK0S0 = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelK0Sdaugh0));
+             AliAODMCParticle* daughK0S1 = dynamic_cast<AliAODMCParticle*>(mcArray->At(labelK0Sdaugh1));
+             if (!daughK0S0 || ! daughK0S1){
+               AliDebug(2, "Could not access K0S daughters, continuing...");
+               continue;
+             }
+             else { // So far: Lc --> K0 + p, K0 --> K0S, K0S in 2 bodies that we can access
+               AliDebug(2, "Step 4 ok: Could access K0S daughters, continuing...");
+               Int_t pdgK0Sdaugh0 = daughK0S0->GetPdgCode();
+               Int_t pdgK0Sdaugh1 = daughK0S1->GetPdgCode();
+               if (TMath::Abs(pdgK0Sdaugh0) != 211 || TMath::Abs(pdgK0Sdaugh1) != 211){
+                 AliDebug(2, "The K0S does not decay in pi+pi-, continuing");
+                 //AliInfo("The K0S does not decay in pi+pi-, continuing");
+               }
+               else { // Full chain: Lc --> K0 + p, K0 --> K0S, K0S --> pi+pi-
+                 if (fAnalCuts->IsInFiducialAcceptance(mcPart->Pt(), mcPart->Y())) {
+                   AliDebug(2, Form("----> Filling histo with pt = %f", mcPart->Pt()));
+                   if(TMath::Abs(mcPart->Y()) < 0.5) fHistoMCLcK0SpGenLimAcc->Fill(mcPart->Pt());
+                   //AliInfo(Form("\nparticle = %d, Filling MC Gen histo\n", iPart));
+                   fHistoMCLcK0SpGen->Fill(mcPart->Pt());
+                   if(!(TMath::Abs(bachelorMC->Eta()) > 0.9 || bachelorMC->Pt() < 0.1 ||
+                        TMath::Abs(daughK0S0->Eta()) > 0.9 || daughK0S0->Pt() < 0.1 ||
+                        TMath::Abs(daughK0S1->Eta()) > 0.9 || daughK0S1->Pt() < 0.1)) {
+                     fHistoMCLcK0SpGenAcc->Fill(mcPart->Pt());
+                   }
+                 }
+                 else {
+                   AliDebug(2, "not in fiducial acceptance! Skipping");
+                   continue;
+                 }
+               }
+             }
+           }
+         }
+       }
+      }
+    }
+  } // closing loop over mcArray
+
+  return; 
+
 }
 
 //-------------------------------------------------------------------------------
@@ -974,7 +1135,7 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::MakeAnalysisForLc2prK0S(TClonesArray *a
       // find associated MC particle for Lc -> p+K0 and K0S->pi+pi
       fmcLabelLc = lcK0spr->MatchToMC(pdgCand, pdgDgLctoV0bachelor[1], pdgDgLctoV0bachelor, pdgDgV0toDaughters, mcArray, kTRUE);
       if (fmcLabelLc>=0) {
-       AliDebug(2,Form(" ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~cascade number %d (total cascade number = %d)", iLctopK0s, nCascades));
+       AliDebug(2, Form("----> cascade number %d (total cascade number = %d) is a Lc!", iLctopK0s, nCascades));
 
        AliAODMCParticle *partLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(fmcLabelLc));
        if(partLc){
@@ -1021,7 +1182,7 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::MakeAnalysisForLc2prK0S(TClonesArray *a
 
     AliAODTrack * v0Pos = dynamic_cast<AliAODTrack*>(lcK0spr->Getv0PositiveTrack());
     AliAODTrack * v0Neg = dynamic_cast<AliAODTrack*>(lcK0spr->Getv0NegativeTrack());
-    if (!v0Neg || !v0Neg) {
+    if (!v0Neg || !v0Pos) {
       AliDebug(2,Form("V0 by cascade %d has no V0positive of V0negative object",iLctopK0s));
       continue;
     }
@@ -1078,10 +1239,46 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::MakeAnalysisForLc2prK0S(TClonesArray *a
            continue; 
          }
        }
+       else if (fKeepingOnlyPYTHIABkg){
+         // we have decided to fill the background only when the candidate has the daugthers that all come from HIJING underlying event!
+         AliAODTrack *bachelor = (AliAODTrack*)lcK0spr->GetBachelor();
+         AliAODTrack *v0pos = (AliAODTrack*)lcK0spr->Getv0PositiveTrack();
+         AliAODTrack *v0neg = (AliAODTrack*)lcK0spr->Getv0NegativeTrack();
+         if (!bachelor || !v0pos || !v0neg) {
+           AliDebug(2, "Cannot retrieve one of the tracks while checking origin, continuing");
+           continue;
+         }
+         else {
+           Int_t labelbachelor = TMath::Abs(bachelor->GetLabel());
+           Int_t labelv0pos = TMath::Abs(v0pos->GetLabel());
+           Int_t labelv0neg = TMath::Abs(v0neg->GetLabel());
+           AliAODMCParticle* MCbachelor =  (AliAODMCParticle*)mcArray->At(labelbachelor);
+           AliAODMCParticle* MCv0pos =  (AliAODMCParticle*)mcArray->At(labelv0pos);
+           AliAODMCParticle* MCv0neg =  (AliAODMCParticle*)mcArray->At(labelv0neg);
+           if (!MCbachelor || !MCv0pos || !MCv0neg) {
+             AliDebug(2, "Cannot retrieve MC particle for one of the tracks while checking origin, continuing");
+             continue;
+           }
+           else {
+             Int_t isBachelorFromPythia = fUtils->CheckOrigin(mcArray, MCbachelor, kTRUE);
+             Int_t isv0posFromPythia = fUtils->CheckOrigin(mcArray, MCv0pos, kTRUE);
+             Int_t isv0negFromPythia = fUtils->CheckOrigin(mcArray, MCv0neg, kTRUE);
+             if (isBachelorFromPythia != 0 && isv0posFromPythia != 0 && isv0negFromPythia != 0){
+               AliDebug(2, "The candidate is from PYTHIA (i.e. all daughters originate from a quark), keeping it to fill background");
+               fHistoBackground->Fill(2);
+             }
+             else {
+               AliDebug(2, "The candidate is NOT from PYTHIA, we skip it when filling background");
+               fHistoBackground->Fill(3);
+               continue; 
+             }
+           }
+         }
+       }
       }
     }
 
-    FillLc2pK0Sspectrum(lcK0spr, isLc, nSelectedAnal, cutsAnal, mcArray);
+    FillLc2pK0Sspectrum(lcK0spr, isLc, nSelectedAnal, cutsAnal, mcArray, iLctopK0s);
     
   }
   
@@ -1093,7 +1290,7 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::FillLc2pK0Sspectrum(AliAODRecoCascadeHF
                                                             Int_t isLc,
                                                             Int_t &nSelectedAnal,
                                                             AliRDHFCutsLctoV0 *cutsAnal,
-                                                            TClonesArray *mcArray){
+                                                            TClonesArray *mcArray, Int_t iLctopK0s){
   //
   // Fill histos for Lc -> K0S+proton
   //
@@ -1151,8 +1348,8 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::FillLc2pK0Sspectrum(AliAODRecoCascadeHF
   Int_t isInV0window = (((cutsAnal->IsSelectedSingleCut(part, AliRDHFCuts::kCandidate, 2)) & (AliRDHFCutsLctoV0::kLcToK0Spr)) == (AliRDHFCutsLctoV0::kLcToK0Spr)); // cut on V0 invMass
 
   if (isInV0window == 0) {
-    AliDebug(2, "No: The candidate has NOT passed the mass cuts!");
-    if (isLc) Printf("SIGNAL candidate rejected");
+    AliDebug(2, "No: The candidate has NOT passed the V0 window cuts!");
+    if (isLc) Printf("SIGNAL candidate rejected: V0 window cuts");
     return;           
   }
   else AliDebug(2, "Yes: The candidate has passed the mass cuts!");
@@ -1161,12 +1358,13 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::FillLc2pK0Sspectrum(AliAODRecoCascadeHF
 
   if (!isInCascadeWindow) {
     AliDebug(2, "No: The candidate has NOT passed the cascade window cuts!");
-    if (isLc) Printf("SIGNAL candidate rejected");
+    if (isLc) Printf("SIGNAL candidate rejected: cascade window cuts");
     return;
   }
   else AliDebug(2, "Yes: The candidate has passed the cascade window cuts!");
 
   Bool_t isCandidateSelectedCuts = (((cutsAnal->IsSelected(part, AliRDHFCuts::kCandidate)) & (AliRDHFCutsLctoV0::kLcToK0Spr)) == (AliRDHFCutsLctoV0::kLcToK0Spr)); // kinematic/topological cuts
+  AliDebug(2, Form("recoAnalysisCuts = %d", cutsAnal->IsSelected(part, AliRDHFCuts::kCandidate) & (AliRDHFCutsLctoV0::kLcToK0Spr)));
   if (!isCandidateSelectedCuts){
     AliDebug(2, "No: Analysis cuts kCandidate level NOT passed");
     if (isLc) Printf("SIGNAL candidate rejected");
@@ -1183,18 +1381,49 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::FillLc2pK0Sspectrum(AliAODRecoCascadeHF
   }
 
   //Bool_t isBachelorID = (((cutsAnal->IsSelected(part,AliRDHFCuts::kPID))&(AliRDHFCutsLctoV0::kLcToK0Spr))==(AliRDHFCutsLctoV0::kLcToK0Spr)); // ID x bachelor
-  Double_t probTPCTOF[AliPID::kSPECIES]={0.};
+  Double_t probTPCTOF[AliPID::kSPECIES]={-1.};
        
   UInt_t detUsed = fPIDCombined->ComputeProbabilities(bachelor, fPIDResponse, probTPCTOF);
-  Double_t probProton = 0.;
-  Double_t probPion = 0.;
-  Double_t probKaon = 0.;
+  AliDebug(2, Form("detUsed (TPCTOF case) = %d", detUsed));
+  Double_t probProton = -1.;
+  //  Double_t probPion = -1.;
+  //  Double_t probKaon = -1.;
   if (detUsed == (UInt_t)fPIDCombined->GetDetectorMask() ) {
+    AliDebug(2, Form("We have found the detector mask for TOF + TPC: probProton will be set to %f", probTPCTOF[AliPID::kProton]));
     probProton = probTPCTOF[AliPID::kProton];
-    probPion = probTPCTOF[AliPID::kPion];
-    probKaon = probTPCTOF[AliPID::kKaon];
-  }
-
+    // probPion = probTPCTOF[AliPID::kPion];
+    // probKaon = probTPCTOF[AliPID::kKaon];
+  }
+  else { // if you don't have both TOF and TPC, try only TPC
+    fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC);
+    AliDebug(2, "We did not find the detector mask for TOF + TPC, let's see only TPC");
+    detUsed = fPIDCombined->ComputeProbabilities(bachelor, fPIDResponse, probTPCTOF);
+    AliDebug(2,Form(" detUsed (TPC case) = %d", detUsed));
+    if (detUsed == (UInt_t)fPIDCombined->GetDetectorMask()) {
+      probProton = probTPCTOF[AliPID::kProton];
+      // probPion = probTPCTOF[AliPID::kPion];
+      // probKaon = probTPCTOF[AliPID::kKaon];
+      AliDebug(2, Form("TPC only worked: probProton will be set to %f", probTPCTOF[AliPID::kProton]));
+    }
+    else {
+      AliDebug(2, "Only TPC did not work...");
+    }
+    // resetting mask to ask for both TPC+TOF
+    fPIDCombined->SetDetectorMask(AliPIDResponse::kDetTPC+AliPIDResponse::kDetTOF);
+  }
+  AliDebug(2, Form("probProton = %f", probProton));
+
+  // now we get the TPC and TOF single PID probabilities (only for Proton, or the tree will explode :) )
+  Double_t probProtonTPC = -1.;
+  Double_t probProtonTOF = -1.;
+  Double_t pidTPC[AliPID::kSPECIES]={-1.};
+  Double_t pidTOF[AliPID::kSPECIES]={-1.};
+  Int_t respTPC = fPIDResponse->ComputePIDProbability(AliPIDResponse::kDetTPC, bachelor, AliPID::kSPECIES, pidTPC);
+  Int_t respTOF = fPIDResponse->ComputePIDProbability(AliPIDResponse::kDetTOF, bachelor, AliPID::kSPECIES, pidTOF);
+  if (respTPC == AliPIDResponse::kDetPidOk) probProtonTPC = pidTPC[AliPID::kProton];
+  if (respTOF == AliPIDResponse::kDetPidOk) probProtonTOF = pidTOF[AliPID::kProton];
+
+  // checking V0 status (on-the-fly vs offline)
   if ( !( !onFlyV0 || (onFlyV0 && fUseOnTheFlyV0) ) ) {
     AliDebug(2, "On-the-fly discarded");
     return;            
@@ -1346,8 +1575,8 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::FillLc2pK0Sspectrum(AliAODRecoCascadeHF
     fCandidateVariables[24] = part->Eta();
     fCandidateVariables[25] = v0pos->Eta();
     fCandidateVariables[26] = v0neg->Eta();
-    fCandidateVariables[27] = probPion;
-    fCandidateVariables[28] = probKaon;
+    fCandidateVariables[27] = probProtonTPC;
+    fCandidateVariables[28] = probProtonTOF;
     fCandidateVariables[29] = bachelor->Eta();
 
     fCandidateVariables[30] = part->P();
@@ -1474,7 +1703,7 @@ void AliAnalysisTaskSELc2V0bachelorTMVA::FillLc2pK0Sspectrum(AliAODRecoCascadeHF
 
     if (fUseMCInfo) {
       if (isLc){
-       AliDebug(2, "Filling Sgn");
+       AliDebug(2, Form("Reco particle %d --> Filling Sgn", iLctopK0s));
        fVariablesTreeSgn->Fill();
        fHistoCodesSgn->Fill(bachCode, k0SCode);          
       }
@@ -1568,9 +1797,9 @@ Int_t AliAnalysisTaskSELc2V0bachelorTMVA::CallKFVertexing(AliAODRecoCascadeHF *c
     }  
   }
   
-  Double_t xn, xp, dca;
+  Double_t xn=0., xp=0.;//, dca;
   AliDebug(2, Form("bField = %f, esdv0Daugh1 = %p, esdv0Daugh2 = %p", fBField, esdv0Daugh1, esdv0Daugh2));
-  dca = esdv0Daugh1->GetDCA(esdv0Daugh2, fBField, xn, xp);
+  //  dca = esdv0Daugh1->GetDCA(esdv0Daugh2, fBField, xn, xp);
     
   AliExternalTrackParam tr1(*esdv0Daugh1);
   AliExternalTrackParam tr2(*esdv0Daugh2);
@@ -1648,11 +1877,13 @@ Int_t AliAnalysisTaskSELc2V0bachelorTMVA::CallKFVertexing(AliAODRecoCascadeHF *c
     if (!tmpdaughv02 && labelsv0daugh[1] > 0){
       AliDebug(2, "Could not access MC info for second daughter of V0, continuing");
     }
-    Double_t xPionMC = tmpdaughv01->Xv(); //Production vertex of Pion --> Where K0S decays
-    Double_t yPionMC = tmpdaughv01->Yv();
-    Double_t zPionMC = tmpdaughv01->Zv();
-    //Printf("Got MC vtx for Pion");
-    Printf("Vertices: MC:  x = %f, y = %f, z = %f", xPionMC, yPionMC, zPionMC);
+    if(tmpdaughv01){
+      Double_t xPionMC = tmpdaughv01->Xv(); //Production vertex of Pion --> Where K0S decays
+      Double_t yPionMC = tmpdaughv01->Yv();
+      Double_t zPionMC = tmpdaughv01->Zv();
+      //Printf("Got MC vtx for Pion");
+      Printf("Vertices: MC:  x = %f, y = %f, z = %f", xPionMC, yPionMC, zPionMC);
+    }
   }
   else {
     Printf("Not a true V0");
@@ -1715,7 +1946,7 @@ Int_t AliAnalysisTaskSELc2V0bachelorTMVA::CallKFVertexing(AliAODRecoCascadeHF *c
   
   if(ftopoConstraint){
     if(ntcheck>0) {
-      //* Add Lc to primary vertex to improve the primary vertex resolution
+      // Add Lc to primary vertex to improve the primary vertex resolution
       primVtxCopy += Lc;
       Lc.SetProductionVertex(primVtxCopy);
     }
@@ -2001,7 +2232,7 @@ Int_t AliAnalysisTaskSELc2V0bachelorTMVA::CallKFVertexing(AliAODRecoCascadeHF *c
            AliDebug(3, "Could not access MC info for second daughter of Lc");
          }
          else { // we can access safely the K0S mother in the MC
-           if( daughv01Lc->GetMother() ==  daughv02Lc->GetMother() && daughv01Lc->GetMother()>=0 ){  // This is a true cascade! bachelor and V0 come from the same mother
+           if( daughv01Lc && (daughv01Lc->GetMother() ==  daughv02Lc->GetMother()) && (daughv01Lc->GetMother()>=0) ){  // This is a true cascade! bachelor and V0 come from the same mother
              //Printf("Lc: The mother has label %d", daughv01Lc->GetMother());
              AliAODMCParticle *motherLc = dynamic_cast<AliAODMCParticle*>(mcArray->At(daughv01Lc->GetMother()));
              Int_t pdgMum = 0, pdgBach = 0, pdgV0 = 0;