updates to check for MC bjet

diff --git a/PWG4/PartCorrDep/AliAnaElectron.cxx b/PWG4/PartCorrDep/AliAnaElectron.cxx
index f13d0e32003159fe61c3861d39d4deca946163ee..ced4e473cc132592ffaa6b4cb8420956b71ebb80 100755 (executable)
--- a/PWG4/PartCorrDep/AliAnaElectron.cxx
+++ b/PWG4/PartCorrDep/AliAnaElectron.cxx
@@ -86,8 +86,8 @@ AliAnaElectron::AliAnaElectron()
   fhPhiConversion(0),fhEtaConversion(0),\r
   //for comparisons with tracking detectors\r
   fhPtHadron(0),fhPtNPEleTPC(0),fhPtNPEleTPCTRD(0),fhPtNPEleTTE(0),\r
-  //for computing efficiency of IPSig tag\r
-  fhBJetPt1x4(0),fhBJetPt2x3(0),fhBJetPt3x2(0),\r
+  //for computing efficiency of B-jet tags\r
+  fhBJetPt1x4(0),fhBJetPt2x3(0),fhBJetPt3x2(0),fhDVMJet(0),\r
   //MC rate histograms/ntuple\r
   fMCEleNtuple(0),fhMCBJetElePt(0),fhPtMCHadron(0),fhPtMCElectron(0)\r
 {\r
@@ -140,8 +140,9 @@ AliAnaElectron::AliAnaElectron(const AliAnaElectron & g)
     //for comparisons with tracking detectors\r
     fhPtHadron(g.fhPtHadron),fhPtNPEleTPC(g.fhPtNPEleTPC),\r
     fhPtNPEleTPCTRD(g.fhPtNPEleTPCTRD),fhPtNPEleTTE(g.fhPtNPEleTTE),\r
-    //for computing efficiency of IPSig tag\r
-    fhBJetPt1x4(g.fhBJetPt1x4),fhBJetPt2x3(g.fhBJetPt2x3),fhBJetPt3x2(g.fhBJetPt3x2),\r
+    //for computing efficiency of B-jet tags\r
+    fhBJetPt1x4(g.fhBJetPt1x4),fhBJetPt2x3(g.fhBJetPt2x3),\r
+    fhBJetPt3x2(g.fhBJetPt3x2),fhDVMJet(g.fhDVMJet),\r
     //MC rate histograms/ntuple\r
     fMCEleNtuple(g.fMCEleNtuple),fhMCBJetElePt(g.fhMCBJetElePt),\r
     fhPtMCHadron(g.fhPtMCHadron),fhPtMCElectron(g.fhPtMCElectron)\r
@@ -236,8 +237,9 @@ AliAnaElectron & AliAnaElectron::operator = (const AliAnaElectron & g)
   //for comparisons with tracking detectors\r
   fhPtHadron = g.fhPtHadron; fhPtNPEleTPC = g.fhPtNPEleTPC; \r
   fhPtNPEleTPCTRD = g.fhPtNPEleTPCTRD; fhPtNPEleTTE = g.fhPtNPEleTTE; \r
-  //for computing efficiency of IPSig tag\r
-  fhBJetPt1x4 = g.fhBJetPt1x4; fhBJetPt2x3 = g.fhBJetPt2x3; fhBJetPt3x2 = g.fhBJetPt3x2; \r
+  //for computing efficiency of B-jet tags\r
+  fhBJetPt1x4 = g.fhBJetPt1x4; fhBJetPt2x3 = g.fhBJetPt2x3; \r
+  fhBJetPt3x2 = g.fhBJetPt3x2; fhDVMJet = g.fhDVMJet;\r
   //MC rate histograms/ntuple\r
   fMCEleNtuple = g.fMCEleNtuple; fhMCBJetElePt = g.fhMCBJetElePt; \r
   fhPtMCHadron = g.fhPtMCHadron; fhPtMCElectron = g.fhPtMCElectron; \r
@@ -420,10 +422,12 @@ TList *  AliAnaElectron::GetCreateOutputObjects()
     fhBJetPt1x4 = new TH1F("hBJetPt1x4","tagged B-jet pT (1 track, ipSignif>4);p_{T}",1000,0.,100.);\r
     fhBJetPt2x3 = new TH1F("hBJetPt2x3","tagged B-jet pT (2 track, ipSignif>3);p_{T}",1000,0.,100.);\r
     fhBJetPt3x2 = new TH1F("hBJetPt3x2","tagged B-jet pT (3 track, ipSignif>2);p_{T}",1000,0.,100.);\r
+    fhDVMJet = new TH2F("hDVM_algo","# DVM jets passing vs Mc-Bjet",10,0,10,300,0.,300.);\r
 \r
     outputContainer->Add(fhBJetPt1x4);\r
     outputContainer->Add(fhBJetPt2x3);\r
     outputContainer->Add(fhBJetPt3x2);\r
+    outputContainer->Add(fhDVMJet);\r
 \r
     //MC Only histograms\r
     \r
@@ -589,6 +593,14 @@ void  AliAnaElectron::MakeAnalysisFillAOD()
     if(TMath::Abs(track->Eta())< 0.5 && TMath::Abs(imp[0])<1.0 && TMath::Abs(imp[1])<1.0) refmult2++;\r
     fhImpactXY->Fill(imp[0]);\r
 \r
+    //JLK CHECK\r
+    AliESDtrack esdTrack(track);\r
+    Double_t tpcpid[AliPID::kSPECIES];\r
+    esdTrack.GetTPCpid(tpcpid);\r
+    Double_t eProb = tpcpid[AliPID::kElectron];\r
+    printf("<%d> ESD eProb = %2.2f\n",itrk,eProb);\r
+\r
+\r
     AliAODPid* pid = (AliAODPid*) track->GetDetPid();\r
     if(pid == 0) {\r
       if(GetDebug() > 0) printf("AliAnaElectron::MakeAnalysisFillAOD() - No PID object - skipping track %d",itrk);\r
@@ -808,7 +820,7 @@ void  AliAnaElectron::MakeAnalysisFillHistograms()
   Int_t njets = (GetReader()->GetOutputEvent())->GetNJets();\r
   if(njets > 0) {\r
     if(GetDebug() > 0) printf("AliAnaElectron::MakeAnalysisFillHistograms() - Jet AOD branch has %d jets.  Performing b-jet tag analysis\n",njets);\r
-    \r
+\r
     for(Int_t ijet = 0; ijet < njets ; ijet++) {\r
       AliAODJet * jet = (AliAODJet*)(GetReader()->GetOutputEvent())->GetJet(ijet) ;\r
       if(GetDebug() > 3) {\r
@@ -849,6 +861,25 @@ void  AliAnaElectron::MakeAnalysisFillHistograms()
 \r
       if(trackCounter[0] > fNTagTrkCut) ipsigJet = kTRUE;\r
 \r
+      if(IsDataMC()) {\r
+       //determine tagging efficiency & mis-tagging rate\r
+       //using b-quarks from stack\r
+       Bool_t isTrueBjet = IsMcBJet(jet->Eta(), jet->Phi() ,stack);\r
+       if (isTrueBjet && GetDebug() > 0) printf("== True Bjet==\n");\r
+       if (dvmJet && GetDebug() > 0)     printf("== found DVM jet==\n");\r
+\r
+       if(isTrueBjet && dvmJet) fhDVMJet->Fill(0.,jet->Pt()); // good tagged\r
+       if(isTrueBjet && !dvmJet) fhDVMJet->Fill(1.,jet->Pt()); // missed tagged\r
+       if(!isTrueBjet && dvmJet) fhDVMJet->Fill(2.,jet->Pt());  // fake tagged\r
+       if(!isTrueBjet && !dvmJet) fhDVMJet->Fill(3.,jet->Pt());  // others\r
+\r
+       if(isTrueBjet) {\r
+         if(trackCounter[1]>0) fhBJetPt1x4->Fill(jet->Pt());\r
+         if(trackCounter[2]>1) fhBJetPt2x3->Fill(jet->Pt());\r
+         if(trackCounter[3]>2) fhBJetPt3x2->Fill(jet->Pt());\r
+       }\r
+      }\r
+\r
       //Fill bjet histograms here\r
       if(!(eJet || ipsigJet || dvmJet)) fhJetType->Fill(0.,jet->Pt()); //none\r
       if(eJet && !(ipsigJet || dvmJet)) fhJetType->Fill(1.,jet->Pt()); //only ejet\r
@@ -956,7 +987,7 @@ void  AliAnaElectron::MakeAnalysisFillHistograms()
   if(IsDataMC()) {\r
 \r
     //MC Jets\r
-    TVector3 jetVect[4];\r
+    TVector3 bjetVect[4];\r
     Int_t nPythiaGenJets = 0;\r
     AliGenPythiaEventHeader*  pythiaGenHeader = (AliGenPythiaEventHeader*)GetReader()->GetGenEventHeader();\r
     if(pythiaGenHeader){\r
@@ -969,8 +1000,11 @@ void  AliAnaElectron::MakeAnalysisFillHistograms()
        pythiaGenHeader->TriggerJet(ip,p);\r
        TVector3 tempVect(p[0],p[1],p[2]);\r
        if ( TMath::Abs(tempVect.Eta())>fJetEtaCut || tempVect.Phi() < fJetPhiMin || tempVect.Phi() > fJetPhiMax) continue;\r
-       jetVect[iCount].SetXYZ(p[0], p[1], p[2]);\r
-       iCount++;\r
+       //Only store it if it has a b-quark within dR < 0.2 of jet axis ?\r
+       if(IsMcBJet(tempVect.Eta(),tempVect.Phi(),stack)) {\r
+         bjetVect[iCount].SetXYZ(p[0], p[1], p[2]);\r
+         iCount++;\r
+       }\r
       }\r
     }\r
         \r
@@ -1000,13 +1034,13 @@ void  AliAnaElectron::MakeAnalysisFillHistograms()
        mctag = GetMCAnalysisUtils()->CheckOrigin(ipart,GetReader(),input);\r
 \r
        if(GetMCSource(mctag)==1) { //bottom electron\r
-         //See if it is within dR < 0.4 of a jet\r
+         //See if it is within dR < 0.4 of a bjet\r
          for(Int_t ij = 0; ij < nPythiaGenJets; ij++) {\r
-           Double_t deta = primary->Eta() - jetVect[ij].Eta();\r
-           Double_t dphi = primary->Phi() - jetVect[ij].Phi();\r
+           Double_t deta = primary->Eta() - bjetVect[ij].Eta();\r
+           Double_t dphi = primary->Phi() - bjetVect[ij].Phi();\r
            Double_t dR = TMath::Sqrt(deta*deta + dphi*dphi);\r
            if(dR < 0.4) {\r
-             fhMCBJetElePt->Fill(primary->Pt(),jetVect[ij].Pt());\r
+             fhMCBJetElePt->Fill(primary->Pt(),bjetVect[ij].Pt());\r
            }\r
          }\r
        }\r
@@ -1521,6 +1555,37 @@ Int_t AliAnaElectron::GetMCSource(Int_t tag)
 \r
 }\r
 \r
+//__________________________________________________________________\r
+Bool_t  AliAnaElectron::IsMcBJet(Double_t eta, Double_t phi, AliStack* stack)\r
+{\r
+  //Check the jet eta,phi against that of the b-quark\r
+  //to decide whether it is an MC B-jet\r
+  Bool_t bjet=kFALSE;\r
+\r
+  //      printf("MTH: McStack ,nparticles=%d \n", stack->GetNtrack() );\r
+  \r
+  for(Int_t ipart = 0; ipart < 100; ipart++) {\r
+\r
+    TParticle* primary = stack->Particle(ipart);\r
+    if (!primary) continue;\r
+    Int_t pdgcode = primary->GetPdgCode();\r
+    if ( TMath::Abs(pdgcode) != 5) continue;\r
+    \r
+    //      printf("MTH: IsMcBJet : %d, pdg=%d : pt=%f \n", ipart, pdgcode, primary->Pt());\r
+    Double_t dphi = phi - primary->Phi();\r
+    Double_t deta = eta - primary->Eta();\r
+    Double_t dr = sqrt(deta*deta + dphi*dphi);\r
+    \r
+    if (dr < 0.2) {\r
+      bjet=kTRUE;\r
+      //printf("MTH: **** found matching MC-Bjet: PDG=%d, pt=%f,dr=%f \n", pdgcode, primary->Pt(),dr );\r
+      break;\r
+    }\r
+  }\r
+  return bjet;\r
+\r
+}\r
+\r
 //__________________________________________________________________\r
 void AliAnaElectron::Print(const Option_t * opt) const\r
 {\r

diff --git a/PWG4/PartCorrDep/AliAnaElectron.h b/PWG4/PartCorrDep/AliAnaElectron.h
index 15ad84b0890e5f6dba6713ad78234f386f2aa8c9..773813dad5072aa04b7024a01f3adad7d3594467 100755 (executable)
--- a/PWG4/PartCorrDep/AliAnaElectron.h
+++ b/PWG4/PartCorrDep/AliAnaElectron.h
@@ -55,6 +55,7 @@ public:
   //check if track has been flagged as a non-photonic or DVM electron\r
   //used with the jet tracks to tag bjets\r
   Bool_t CheckTrack(const AliAODTrack* track,const char* type);  \r
+  Bool_t IsMcBJet(Double_t x, Double_t y, AliStack* st);\r
 \r
   void Print(const Option_t * opt)const;\r
   \r
@@ -226,6 +227,8 @@ public:
   TH1F * fhBJetPt2x3;    //! IPSig B-tagging : result for (2 track, ipSignif>3)\r
   TH1F * fhBJetPt3x2;    //! IPSig B-tagging : result for (3 track, ipSignif>2)\r
 \r
+  TH2F* fhDVMJet;        //! DVM jet algo check\r
+\r
   ////////////////////////////\r
   //MC Only Rate histograms\r
 \r
@@ -235,7 +238,7 @@ public:
   TH1F* fhPtMCHadron;    //! Pt distribution of MC charged hadrons (pi,k,p) in EMCAL acceptance\r
   TH2F* fhPtMCElectron;  //! Pt distribution of MC electrons from various sources in EMCAL\r
 \r
-  ClassDef(AliAnaElectron,7)\r
+  ClassDef(AliAnaElectron,8)\r
 \r
 } ;\r
  \r