Changes from Tomas:
authorgconesab <gconesab@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 5 Oct 2010 17:24:58 +0000 (17:24 +0000)
committergconesab <gconesab@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 5 Oct 2010 17:24:58 +0000 (17:24 +0000)
AliAnaElectron: Use positve track labels. Removed condition that the momentum vector should point to the emcal. Should be if only the track itself points to the emcal.
B-field problem solved.

parmaker: in PWG4/macros/electrons is not working properly. Changed it
and added new script to make it work properly. New file is
makeparfilesproperly

makeparfilesproperly:  new bash script to make parfiles in PWG4/macros/electrons

AliAnaBtag: cleaned it up a bit and added btagging to tracks not only associated to the
geometry of the EMCAL.

PWG4/PartCorrDep/AliAnaBtag.cxx
PWG4/PartCorrDep/AliAnaBtag.h
PWG4/PartCorrDep/AliAnaElectron.cxx
PWG4/macros/ConfigAnalysisBtag.C [moved from PWG4/macros/ConfigAnalysisBTag.C with 57% similarity]
PWG4/macros/electrons/makeparfilesproperly [new file with mode: 0644]
PWG4/macros/electrons/parmaker

index dac5e61..f039247 100644 (file)
 //_________________________________________________________________________\r
 //\r
 // Class for the electron identification and B-tagging.\r
-// Clusters from EMCAL matched to tracks\r
-// and kept in the AOD. Few histograms produced.\r
+// Clusters from EMCAL matched to tracks to id electrons.\r
+// Btagger is run on all electrons, then jets are tagged as well.\r
+// \r
 //\r
-// -- Author: T.R.P.Aronsson (Yale) J.L. Klay (Cal Poly), M. Heinz (Yale)\r
+// -- Author: T.R.P.Aronsson (Yale), M. Heinz (Yale)\r
 //////////////////////////////////////////////////////////////////////////////\r
   \r
-// --- ROOT system --- \r
 #include <TH2F.h>\r
 #include <TH3F.h>\r
 #include <TParticle.h>\r
+#include <TNtuple.h>\r
 #include <TClonesArray.h>\r
-//#include <TObjString.h>\r
-//#include <Riostream.h>\r
 \r
-// --- Analysis system --- \r
 #include "AliAnaBtag.h" \r
 #include "AliCaloTrackReader.h"\r
 #include "AliMCAnalysisUtils.h"\r
 #include "AliAODJet.h"\r
 #include "AliAODEvent.h"\r
 #include "AliGenPythiaEventHeader.h"\r
-//#include "iostream.h"\r
 ClassImp(AliAnaBtag)\r
   \r
 //____________________________________________________________________________\r
 AliAnaBtag::AliAnaBtag() \r
-: AliAnaPartCorrBaseClass(),fCalorimeter(""),\r
-  fpOverEmin(0.),fpOverEmax(0.),fResidualCut(0.),fMinClusEne(0.),\r
-  fDrCut(0.),fPairDcaCut(0.),fDecayLenCut(0.),fImpactCut(0.),\r
+: AliAnaPartCorrBaseClass(),\r
+  fWriteNtuple(0),electrons(0),pairs(0),events(0),fEventNumber(0),fNElec(0),fNElecEv(0),fNPair(0),fDrCut(0.),fPairDcaCut(0.),fDecayLenCut(0.),fImpactCut(0.),\r
   fAssocPtCut(0.),fMassCut(0.),fSdcaCut(0.),fITSCut(0),\r
   fNTagTrkCut(0),fIPSigCut(0.),fJetEtaCut(0.3),fJetPhiMin(1.8),fJetPhiMax(2.9),\r
-  fhEmcalElectrons(0),fhTRDElectrons(0),fhTPCElectrons(0),fhDVM1(0),fhDVM2(0),fhJets(0),fhJetsAllEtaPhi(0),fhJetsLeadingBElectronEtaPhi(0),fhDVM1EtaPhi(0),fhBJetElectronDetaDphi(0),fhClusterEnergy(0),fhTestalle(0),fhResidual(0),fhElectrons(0),fhTracks(0)\r
+  fhEmcalElectrons(0),fhTRDElectrons(0),fhTPCElectrons(0),fhEmcalMCE(0),fhTRDMCE(0),fhTPCMCE(0),fhEmcalMCEFake(0),fhTRDMCEFake(0),fhTPCMCEFake(0),fhEmcalMCP(0),fhTRDMCP(0),fhTPCMCP(0),fhEmcalMCK(0),fhTRDMCK(0),fhTPCMCK(0),fhSpecies(0),fhDVM1(0),fhDVM2(0),fhNVTX(0),fhNVTXMC(0),fhJets(0),fhJetsAllEtaPhi(0),fhJetsLeadingBElectronEtaPhi(0),fhDVM1EtaPhi(0),fhBJetElectronDetaDphi(0),fhClusterMap(0),fhClusterEnergy(0),fhTestalle(0),fhResidual(0),fhPairPt(0),fhElectrons(0),fhTracks(0)\r
 {\r
   //default ctor\r
-  \r
   //Initialize parameters\r
   InitParameters();\r
 \r
@@ -83,91 +78,135 @@ TList *  AliAnaBtag::GetCreateOutputObjects()
   TList * outputContainer = new TList() ; \r
   outputContainer->SetName("ElectronHistos") ; \r
 \r
+  if(fWriteNtuple){\r
+    electrons = new TNtuple("electrons","Electron Ntuple","electronnumber:eventnumber:electronineventnumber:pairs:start:stop:tpc:trd:emc:mcpdg:parentbit:btag:pt");\r
+    outputContainer->Add(electrons);\r
 \r
-       fhEmcalElectrons = new TH1F("fhEmcalElectrons","",400,0,400);\r
-    outputContainer->Add(fhEmcalElectrons);\r
-    \r
-    fhTRDElectrons = new TH1F("fhTRDElectrons","",400,0,400);\r
-    outputContainer->Add(fhTRDElectrons);\r
+    pairs = new TNtuple("pairs","Pair Ntuple","pairnumber:electronnumber:eventnumber:pdca:sdca:minv:pth:massphoton:decaylength");\r
+    outputContainer->Add(pairs);\r
 \r
-    fhTPCElectrons = new TH1F("fhTPCElectrons","",400,0,400);\r
-    outputContainer->Add(fhTPCElectrons);\r
+    events = new TNtuple("events","Event NTuple","event");\r
+    outputContainer->Add(events);\r
+  }\r
+  \r
+  fhEmcalElectrons = new TH1F("fhEmcalElectrons","",400,0,400);\r
+  outputContainer->Add(fhEmcalElectrons);\r
+  \r
+  fhTRDElectrons = new TH1F("fhTRDElectrons","",400,0,400);\r
+  outputContainer->Add(fhTRDElectrons);\r
+  \r
+  fhTPCElectrons = new TH1F("fhTPCElectrons","",400,0,400);\r
+  outputContainer->Add(fhTPCElectrons);\r
+  \r
 \r
+  if(IsDataMC()){\r
+    fhEmcalMCE = new TH1F("fhEmcalMCE","",400,0,400);\r
+    outputContainer->Add(fhEmcalMCE);\r
+    \r
+    fhTRDMCE = new TH1F("fhTRDMCE","",400,0,400);\r
+    outputContainer->Add(fhTRDMCE);\r
+    \r
+    fhTPCMCE = new TH1F("fhTPCMCE","",400,0,400);\r
+    outputContainer->Add(fhTPCMCE);\r
+    \r
+    fhEmcalMCEFake = new TH1F("fhEmcalMCEFake","",400,0,400);\r
+    outputContainer->Add(fhEmcalMCEFake);\r
+    \r
+    fhTRDMCEFake = new TH1F("fhTRDMCEFake","",400,0,400);\r
+    outputContainer->Add(fhTRDMCEFake);\r
+    \r
+    fhTPCMCEFake = new TH1F("fhTPCMCEFake","",400,0,400);\r
+    outputContainer->Add(fhTPCMCEFake);\r
+    \r
+    fhEmcalMCP = new TH1F("fhEmcalMCP","",400,0,400);\r
+    outputContainer->Add(fhEmcalMCP);\r
+    \r
+    fhTRDMCP = new TH1F("fhTRDMCP","",400,0,400);\r
+    outputContainer->Add(fhTRDMCP);\r
+    \r
+    fhTPCMCP = new TH1F("fhTPCMCP","",400,0,400);\r
+    outputContainer->Add(fhTPCMCP);\r
+    \r
+    fhEmcalMCK = new TH1F("fhEmcalMCK","",400,0,400);\r
+    outputContainer->Add(fhEmcalMCK);\r
+    \r
+    fhTRDMCK = new TH1F("fhTRDMCK","",400,0,400);\r
+    outputContainer->Add(fhTRDMCK);\r
+    \r
+    fhTPCMCK = new TH1F("fhTPCMCK","",400,0,400);\r
+    outputContainer->Add(fhTPCMCK);\r
+    \r
+    fhSpecies  = new TH1F("fhSpecies","",1000,0,1000);\r
+    outputContainer->Add(fhSpecies);\r
+    \r
     fhDVM1 = new TH1F("fhDVM1","",400,0,400);\r
     outputContainer->Add(fhDVM1);\r
-\r
+    \r
     fhDVM2 = new TH1F("fhDVM2","",400,0,400);\r
     outputContainer->Add(fhDVM2);\r
+    \r
+    fhNVTX = new TH1F("fhNVTX","",20,0,20);\r
+    outputContainer->Add(fhNVTX);\r
+    \r
+    fhNVTXMC = new TH1F("fhNVTXMC","",20,0,20);\r
+    outputContainer->Add(fhNVTXMC);\r
+  }\r
+  \r
 \r
+  fhJets = new TH2F("fhJets","",400,0,400,20,0,20);\r
+  outputContainer->Add(fhJets);\r
+  \r
+  fhJetsAllEtaPhi = new TH2F("fhJetsAllEtaPhi","",100,-2,2,100,-2,8);\r
+  outputContainer->Add(fhJetsAllEtaPhi);\r
+  \r
+  fhJetsLeadingBElectronEtaPhi = new TH2F("fhJetsLeadingBElectronEtaPhi","",100,-5,5,200,-10,10);\r
+  outputContainer->Add(fhJetsLeadingBElectronEtaPhi);\r
+  \r
+  fhDVM1EtaPhi = new TH2F("fhDVM1EtaPhi","",100,-2,2,100,-2,8);\r
+  outputContainer->Add(fhDVM1EtaPhi);\r
+  \r
+  fhBJetElectronDetaDphi = new TH2F("fhBJetElectronDetaDphi","",100,-5,5,200,-10,10);\r
+  outputContainer->Add(fhBJetElectronDetaDphi);\r
+  \r
+  fhClusterMap = new TH2F("fhClusterMap","",100,-2,2,100,-2,8);\r
+  outputContainer->Add(fhClusterMap);\r
+  \r
+  fhClusterEnergy = new TH1F("fhClusterEnergy","",100,0,10);\r
+  outputContainer->Add(fhClusterEnergy);\r
+  \r
+  fhTestalle = new TH1F("fhTestalle","",400,0,400);\r
+  outputContainer->Add(fhTestalle);\r
+  \r
+  fhResidual = new TH1F("fhResidual","",500,0,5);\r
+  outputContainer->Add(fhResidual);\r
 \r
-    fhJets = new TH2F("fhJets","",400,0,400,20,0,20);\r
-    outputContainer->Add(fhJets);\r
-\r
-    fhJetsAllEtaPhi = new TH2F("fhJetsAllEtaPhi","",100,-2,2,100,-2,8);\r
-    outputContainer->Add(fhJetsAllEtaPhi);\r
-\r
-    fhJetsLeadingBElectronEtaPhi = new TH2F("fhJetsLeadingBElectronEtaPhi","",100,-5,5,200,-10,10);\r
-    outputContainer->Add(fhJetsLeadingBElectronEtaPhi);\r
-\r
-    fhDVM1EtaPhi = new TH2F("fhDVM1EtaPhi","",100,-2,2,100,-2,8);\r
-    outputContainer->Add(fhDVM1EtaPhi);\r
-\r
-    fhBJetElectronDetaDphi = new TH2F("fhBJetElectronDetaDphi","",100,-5,5,200,-10,10);\r
-    outputContainer->Add(fhBJetElectronDetaDphi);\r
-\r
-    fhClusterEnergy = new TH1F("fhClusterEnergy","",100,0,10);\r
-    outputContainer->Add(fhClusterEnergy);\r
-\r
-    fhTestalle = new TH1F("fhTestalle","",400,0,400);\r
-    outputContainer->Add(fhTestalle);\r
-\r
-    fhResidual = new TH1F("fhResidual","",500,0,5);\r
-    outputContainer->Add(fhResidual);\r
-\r
-    fhElectrons = new TH2F("fhElectrons","",200,0,100,20,0,20);\r
-    outputContainer->Add(fhElectrons);\r
-\r
-    fhTracks = new TH2F("fhTracks","",200,0,100,20,0,20);\r
-    outputContainer->Add(fhTracks);\r
-\r
-\r
+  fhPairPt = new TH1F("fhPairPt","",400,0,400);\r
+  outputContainer->Add(fhPairPt);\r
   \r
-  return outputContainer ;\r
+  fhElectrons = new TH2F("fhElectrons","",200,0,100,20,0,20);\r
+  outputContainer->Add(fhElectrons);\r
   \r
+  fhTracks = new TH2F("fhTracks","",200,0,100,20,0,20);\r
+  outputContainer->Add(fhTracks);\r
+    \r
+  return outputContainer ; \r
 }\r
 \r
 //____________________________________________________________________________\r
 void AliAnaBtag::Init()\r
 {\r
-\r
   //do some initialization\r
-  if(fCalorimeter == "PHOS") {\r
-    printf("AliAnaBtag::Init() - !!STOP: You want to use PHOS in analysis but this is not (yet) supported!!\n!!Check the configuration file!!\n");\r
-    fCalorimeter = "EMCAL";\r
-  }\r
-  if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn()){\r
-    printf("AliAnaBtag::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!!\n!!Check the configuration file!!\n");\r
-    abort();\r
-  }\r
-\r
+    printf("Rubbish init step AliAnaBtag::Init()");\r
 }\r
 \r
 \r
 //____________________________________________________________________________\r
 void AliAnaBtag::InitParameters()\r
-{\r
-  \r
+{ \r
   //Initialize the parameters of the analysis.\r
   SetOutputAODClassName("AliAODPWG4Particle");\r
   SetOutputAODName("PWG4Particle");\r
-\r
-  AddToHistogramsName("Btag_");\r
-\r
-  fCalorimeter = "EMCAL" ;\r
-  fpOverEmin = 0.5;\r
-  fpOverEmax = 1.2;\r
-  fResidualCut = 0.02;\r
-  fMinClusEne = 4.0;\r
+  //AddToHistogramsName("");\r
   //DVM B-tagging\r
   fDrCut       = 1.0; \r
   fPairDcaCut  = 0.02;\r
@@ -189,18 +228,17 @@ void AliAnaBtag::InitParameters()
 //__________________________________________________________________\r
 void  AliAnaBtag::MakeAnalysisFillAOD() \r
 {\r
+  fEventNumber++;\r
+  fNElecEv=0;\r
+  if(fWriteNtuple)\r
+    events->Fill(fEventNumber);\r
+\r
   //This reads in tracks, extrapolates to EMCAL, does p/E selectrons, identifies electron candidates\r
   //After candidates are obtained, btagging and saving into AOD.\r
-\r
-  Double_t bfield = 0.;\r
-  if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) bfield = GetReader()->GetBField();\r
-\r
-  //Select the calorimeter of the electron\r
-  if(fCalorimeter != "EMCAL") {\r
-    printf("This class not yet implemented for PHOS\n");\r
-    abort();\r
-  }\r
-  \r
+  AliStack *stack =0x0;\r
+  //Double_t bfield = 0.;\r
+  if(GetDebug()>0) printf("AliAnaBtag::MakeAnalysisFillAOD() - Write ntuple flag is %d \n",fWriteNtuple);\r
+  //if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) bfield = GetReader()->GetBField();\r
   TObjArray *cl = GetAODEMCAL();\r
   \r
   if(!GetAODCTS() || GetAODCTS()->GetEntriesFast() == 0) return ;\r
@@ -211,26 +249,36 @@ void  AliAnaBtag::MakeAnalysisFillAOD()
   Int_t iCluster = -999;\r
   Int_t ntot = cl->GetEntriesFast();\r
 \r
-  //CLUSTER STUFF\r
-  if(1){\r
-    for(Int_t iclus = 0; iclus < ntot; iclus++) {\r
-      AliVCluster * clus = (AliVCluster*) (cl->At(iclus));\r
-      if(!clus) continue;\r
-      fhClusterEnergy->Fill(clus->E());\r
-    }\r
-  }\r
+  //CLUSTER STUFF \r
+  for(Int_t iclus = 0; iclus < ntot; iclus++) {\r
+    AliVCluster * clus = (AliVCluster*) (cl->At(iclus));\r
+    if(!clus) continue;\r
+    fhClusterEnergy->Fill(clus->E());\r
+    Float_t xclus[3];\r
+    clus->GetPosition(xclus);\r
+    TVector3 cluspos(xclus[0],xclus[1],xclus[2]);\r
 \r
+    fhClusterMap->Fill(cluspos.Eta(),cluspos.Phi());\r
+  }\r
+  \r
+  \r
 \r
   for (Int_t itrk =  0; itrk <  ntracks; itrk++) {////////////// track loop\r
     iCluster = -999; //start with no match\r
     AliAODTrack * track = (AliAODTrack*) (GetAODCTS()->At(itrk)) ;\r
-\r
+    if(track->GetLabel()<0){\r
+      if(GetDebug()>0)\r
+       printf("Negative track label, aborting!\n");\r
+      continue;\r
+    }\r
     Double_t imp[2] = {-999.,-999.}; Double_t cov[3] = {-999.,-999.,-999.};\r
     Bool_t dcaOkay = GetDCA(track,imp,cov);  //homegrown dca calculation until AOD is fixed\r
-    if(!dcaOkay) printf("AliAnaBtag::Problem computing DCA to primary vertex for track %d.  Skipping it...\n",itrk);\r
-\r
+    if(!dcaOkay&&GetDebug()>0) printf("AliAnaBtag::FillAOD - Problem computing DCA to primary vertex for track %d.  Skipping it...\n",itrk);\r
     fhTracks->Fill(track->Pt(),1);\r
 \r
+    if(track->Pt()<0)\r
+      continue;\r
+\r
     AliAODPid* pid = (AliAODPid*) track->GetDetPid();\r
     if(pid == 0) {\r
       if(GetDebug() > 0) printf("AliAnaBtag::MakeAnalysisFillAOD() - No PID object - skipping track %d",itrk);\r
@@ -249,31 +297,27 @@ void  AliAnaBtag::MakeAnalysisFillAOD()
       Double_t tmom = mom.Mag();\r
 \r
       Bool_t in = kFALSE;\r
-      if(mom.Phi()*180./TMath::Pi() > 80. && mom.Phi()*180./TMath::Pi() < 190. &&\r
-        mom.Eta() > -0.7 && mom.Eta() < 0.7) in = kTRUE; //Checks momentum?\r
-      //Also check the track\r
       if(track->Phi()*180./TMath::Pi() > 80. && track->Phi()*180./TMath::Pi() < 190. &&\r
         track->Eta() > -0.7 && track->Eta() < 0.7) in = kTRUE;\r
 \r
-      if(GetDebug() > 1) printf("AliAnaBtag::MakeAnalysisFillAOD() - Track(Extrap) pt %2.2f(%2.2f), phi %2.2f(%2.2f), eta %2.2f(%2.2f) in fiducial cut %d\n",track->Pt(), mom.Pt(), track->Phi(), mom.Phi(), track->Eta(),mom.Eta(), in);\r
 \r
-      if(mom.Pt() > GetMinPt() && in) {\r
+      Double_t dEdx = pid->GetTPCsignal();\r
+      Int_t pidProb = track->GetMostProbablePID();\r
+      Bool_t tpcEle = kFALSE; if(dEdx > 70.) tpcEle = kTRUE;\r
+      Bool_t trkEle = kFALSE; if(pidProb == AliAODTrack::kElectron) trkEle = kTRUE;\r
+      Bool_t emcEle = kFALSE;   \r
+\r
+\r
+\r
+\r
+      ////////////////////////////////////////////////EMCAL//////////////////////\r
+      if(mom.Pt() > 1.0 && in) {\r
        fhTracks->Fill(track->Pt(),3);\r
-       Double_t dEdx = pid->GetTPCsignal();\r
-               \r
-       //Int_t ntot = cl->GetEntriesFast();\r
        Double_t res = 999.;\r
        Double_t pOverE = -999.;\r
        \r
-       Int_t pidProb = track->GetMostProbablePID();\r
-       Bool_t tpcEle = kFALSE; if(dEdx > 70.) tpcEle = kTRUE;\r
-       Bool_t trkEle = kFALSE; if(pidProb == AliAODTrack::kElectron) trkEle = kTRUE;\r
-       Bool_t trkChgHad = kFALSE; if(pidProb == AliAODTrack::kPion || pidProb == AliAODTrack::kKaon || pidProb == AliAODTrack::kProton) trkChgHad = kTRUE;\r
-\r
-\r
-       //Track Matching!\r
-       Bool_t emcEle = kFALSE;      \r
-       double minRes=100.;\r
+       //Track Matching parameters\r
+       double minRes=100.;\r
        Double_t minR  = 99;\r
         Double_t minPe =-1;\r
         Double_t minEp =-1;\r
@@ -284,17 +328,7 @@ void  AliAnaBtag::MakeAnalysisFillAOD()
          AliVCluster * clus = (AliVCluster*) (cl->At(iclus));\r
          if(!clus) continue;\r
 \r
-         //As of 11-Oct-2009\r
-         //only select "good" clusters   \r
-//       if (clus->GetNCells()       < 2    ) continue;\r
-//           if (clus->GetNCells()       > 30   ) continue;\r
-//           if (clus->E()               < fMinClusEne ) continue;\r
-//           if (clus->GetDispersion()   > 1    ) continue;\r
-//           if (clus->GetM20()          > 0.4  ) continue;\r
-//           if (clus->GetM02()          > 0.4  ) continue;\r
-//           if (clus->GetM20()          < 0.03 ) continue;\r
-//           if (clus->GetM02()          < 0.03 ) continue;\r
-// new optimized from ben. 2010May\r
+         // new optimized from ben. 2010May\r
          if (clus->GetNCells()       < 2    ) continue;\r
           if (clus->GetNCells()       > 35   ) continue;\r
           if (clus->E()               < 0 ) continue;\r
@@ -314,18 +348,12 @@ void  AliAnaBtag::MakeAnalysisFillAOD()
          if(dphi < -TMath::Pi()) dphi += 2*TMath::Pi();\r
 \r
          res = sqrt(dphi*dphi + deta*deta);\r
-         \r
-         if(res<minRes){\r
-           minRes=res;\r
-           \r
-         }\r
-         \r
-         if(res < 0.0275) { //   if(res < fResidualCut) { //Optimized from Ben\r
+         if(res<minRes)  minRes=res;       \r
+\r
+         if(res < 0.0275) { // { //Optimized from Ben\r
            iCluster = iclus;\r
-           \r
            Double_t energy = clus->E(); \r
            if(energy > 0) pOverE = tmom/energy;\r
-           \r
            if (res< minR) {\r
               minR  = res;\r
               minPe = pOverE;\r
@@ -333,83 +361,132 @@ void  AliAnaBtag::MakeAnalysisFillAOD()
               minMult = clus->GetNCells() ;\r
               minPt = track->Pt();\r
             }\r
-\r
-\r
          } else {\r
              //unmatched\r
          }//res cut\r
-       }//calo cluster loop\r
 \r
+       }//calo cluster loop\r
        fhResidual->Fill(minRes);\r
 \r
        if(minPe > 0.9 && minPe < 1.08) emcEle = kTRUE;//       if(minPe > fpOverEmin && minPe < fpOverEmax) emcEle = kTRUE;\r
-       \r
+       \r
+      }//pt, fiducial selection = EMCAL ////////////////////END EMCAL/////////////////////////\r
+\r
+\r
+\r
 \r
+      ////////////////////////////////////////////////////Electrons/////////////////////     \r
+      if(emcEle ||tpcEle || trkEle) { //Obsolete (kinda...)\r
+       fhTestalle->Fill(track->Pt());\r
+\r
+       //B-tagging\r
+       if(GetDebug() > 1) printf("Found Electron - do b-tagging\n");\r
+       Int_t pairs1=0,start=0,stop=0;\r
+       Int_t dvmbtag = GetDVMBtag(track,pairs1,start,stop); //add: get back #pairs, start stop in pair-Ntuple.\r
+       if(GetDebug() > 0) printf("AliAnaBtag::MakeAnalysisFillAOD -  Analyze, got back result from dvm: pair counts: pairs: %d, start %d, stop %d. \n",pairs1,start,stop);\r
+       fhNVTX->Fill(dvmbtag);\r
+\r
+       if(dvmbtag>0){\r
+         fhDVM1->Fill(track->Pt());\r
+         fhDVM1EtaPhi->Fill(track->Eta(),track->Phi());            \r
+       }\r
+       if(dvmbtag>1)\r
+         fhDVM2->Fill(track->Pt());\r
        \r
-       if(emcEle)\r
+       //Create particle to save in AOD///////////// Purpose of this is the AODJets needs to check against this.\r
+       Double_t eMass = 0.511/1000; //mass in GeV\r
+       Double_t eleE = sqrt(track->P()*track->P() + eMass*eMass);\r
+       AliAODPWG4Particle tr = AliAODPWG4Particle(track->Px(),track->Py(),track->Pz(),eleE);\r
+       tr.SetLabel(track->GetLabel());\r
+       tr.SetCaloLabel(iCluster,-1); //sets the indices of the original caloclusters\r
+       tr.SetTrackLabel(track->GetID(),-1); //sets the indices of the original tracks tr.SetTrackLabel(track->GetID(),-1) instead of itrk;\r
+       tr.SetBtag(dvmbtag);\r
+       if(track->Charge() < 0) tr.SetPdg(11); //electron is 11\r
+       else  tr.SetPdg(-11); //positron is -11 \r
+\r
+       //Set detector flags\r
+       Int_t emcflag=0,tpcflag=0,trdflag=0;\r
+       if(emcEle){\r
          fhEmcalElectrons->Fill(track->Pt());\r
-       if(trkEle)\r
+         emcflag=1;}\r
+       if(trkEle){\r
          fhTRDElectrons->Fill(track->Pt());\r
-       if(tpcEle)\r
+         trdflag=1;}\r
+       if(tpcEle){\r
          fhTPCElectrons->Fill(track->Pt());\r
-       \r
-       \r
-       //Take all emcal electrons, but the others only if pT < 10 GeV\r
-       if(emcEle && (tpcEle || trkEle) ) {\r
-         fhTestalle->Fill(track->Pt());\r
-         //B-tagging\r
-         if(GetDebug() > 1) printf("Found Electron - do b-tagging\n");\r
-         Int_t dvmbtag = GetDVMBtag(track);\r
-         \r
-         \r
-         if(dvmbtag>0){\r
-           fhDVM1->Fill(track->Pt());\r
-            fhDVM1EtaPhi->Fill(track->Eta(),track->Phi());       \r
-           \r
-          }\r
-         if(dvmbtag>1)\r
-           fhDVM2->Fill(track->Pt());\r
-         \r
-\r
-         //Create particle to save in AOD\r
-         Double_t eMass = 0.511/1000; //mass in GeV\r
-         Double_t eleE = sqrt(track->P()*track->P() + eMass*eMass);\r
-         AliAODPWG4Particle tr = AliAODPWG4Particle(track->Px(),track->Py(),track->Pz(),eleE);\r
-         tr.SetLabel(track->GetLabel());\r
-         tr.SetCaloLabel(iCluster,-1); //sets the indices of the original caloclusters\r
-         tr.SetTrackLabel(track->GetID(),-1); //sets the indices of the original tracks tr.SetTrackLabel(track->GetID(),-1) instead of itrk;\r
-\r
-         if(emcEle) {//PID determined by EMCAL\r
-           tr.SetDetector(fCalorimeter);\r
-         } else {\r
-           tr.SetDetector("CTS"); //PID determined by CTS\r
+         tpcflag=1;}\r
+\r
+       if(emcEle) {//PID determined by EMCAL\r
+         tr.SetDetector("EMCAL");\r
+       } else {\r
+         tr.SetDetector("CTS"); //PID determined by CTS\r
+       }\r
+\r
+\r
+       /////////////////////////MC stuff////////////////////////////////////////////////\r
+       if(IsDataMC()){\r
+         stack=GetMCStack();\r
+         if(!stack) printf("AliAnaBtag::MakeAnalysisFillHistograms() - Crap, no stack: \n");\r
+\r
+         //Is it really an electron?\r
+         TParticle *partX = stack->Particle(TMath::Abs(track->GetLabel()));\r
+         int pdg = 0;\r
+         pdg = partX->GetPdgCode();\r
+         fhSpecies->Fill(TMath::Abs(pdg));\r
+         if(TMath::Abs(pdg)==11){ //Check MC electrons\r
+           if(emcEle)\r
+             fhEmcalMCE->Fill(track->Pt());\r
+           if(trkEle)\r
+             fhTRDMCE->Fill(track->Pt());\r
+           if(tpcEle)\r
+             fhTPCMCE->Fill(track->Pt());\r
+         }else{ //Fake histos!\r
+           if(emcEle)\r
+             fhEmcalMCEFake->Fill(track->Pt());\r
+           if(trkEle)\r
+             fhTRDMCEFake->Fill(track->Pt());\r
+           if(tpcEle)\r
+             fhTPCMCEFake->Fill(track->Pt());\r
+         }\r
+         if(TMath::Abs(pdg)==211){ //Check MC pions\r
+           if(emcEle)\r
+             fhEmcalMCP->Fill(track->Pt());\r
+           if(trkEle)\r
+             fhTRDMCP->Fill(track->Pt());\r
+           if(tpcEle)\r
+             fhTPCMCP->Fill(track->Pt());\r
+         }\r
+         if(TMath::Abs(pdg)==321){ //Check MC Kaons\r
+           if(emcEle)\r
+             fhEmcalMCK->Fill(track->Pt());\r
+           if(trkEle)\r
+             fhTRDMCK->Fill(track->Pt());\r
+           if(tpcEle)\r
+             fhTPCMCK->Fill(track->Pt());\r
          }\r
+         //Take care of where it came from (parent bit)\r
+         tr.SetTag(GetMCAnalysisUtils()->CheckOrigin(tr.GetLabel(),GetReader(),tr.GetInputFileIndex())); //Gets a tag bit which contains info about super grandfather particle. Use (tag&(1<<11)), 11 for direct b, and 9 for B->C\r
 \r
-         //if(GetReader()->GetAODCTSNormalInputEntries() <= itrk) tr.SetInputFileIndex(1);\r
-         //Make this preserve sign of particle\r
-         if(track->Charge() < 0) tr.SetPdg(11); //electron is 11\r
-         else  tr.SetPdg(-11); //positron is -11\r
+         if(tr.GetTag()&(1<<9)||tr.GetTag()&(1<<11)) //MC particle from b-decay\r
+           fhNVTXMC->Fill(dvmbtag);\r
 \r
-         tr.SetBtag(dvmbtag);\r
-         \r
+         if(fWriteNtuple){\r
+           fNElec++;\r
+           fNElecEv++;\r
+           electrons->Fill(fNElec,fEventNumber,fNElecEv,pairs1,start,stop,tpcflag,trdflag,emcflag,pdg,tr.GetTag(),tr.GetBtag(),tr.Pt());\r
 \r
-         //Check origin of the candiates\r
-         if(IsDataMC()){\r
-           tr.SetTag(GetMCAnalysisUtils()->CheckOrigin(tr.GetLabel(),GetReader(),tr.GetInputFileIndex()));\r
-           if(GetDebug() > 0) printf("AliAnaBtag::MakeAnalysisFillAOD() - Origin of candidate %d\n",tr.GetTag());\r
-         }//Work with stack also   \r
-         \r
-         AddAODParticle(tr);\r
+\r
+         }\r
          \r
-         if(GetDebug() > 1) printf("AliAnaBtag::MakeAnalysisFillAOD() - Electron selection cuts passed: pT %3.2f, pdg %d\n",tr.Pt(),tr.GetPdg());      \r
-       }//electron\r
-      }//pt, fiducial selection\r
-    }//pid check\r
+         if(GetDebug() > 0) \r
+           printf("AliAnaBtag::MakeAnalysisFillAOD() - Origin of candidate (parent bit) %d\n",tr.GetTag());\r
+       }//MonteCarlo MC done\r
+       \r
+       AddAODParticle(tr);             \r
+      }//electron\r
+    \r
+    } //pid check\r
   }//track loop                         \r
-  \r
-\r
-\r
-\r
   if(GetDebug() > 1) printf("AliAnaBtag::MakeAnalysisFillAOD()  End fill AODs \n");  \r
   \r
 }\r
@@ -418,25 +495,17 @@ void  AliAnaBtag::MakeAnalysisFillAOD()
 void  AliAnaBtag::MakeAnalysisFillHistograms() \r
 {\r
   //Do analysis and fill histograms\r
-  \r
   AliStack * stack = 0x0;\r
   //   TParticle * primary = 0x0;\r
-  \r
-  \r
-  \r
-  \r
   if(IsDataMC()) {\r
     if(GetReader()->ReadStack()){\r
       stack =  GetMCStack() ;      \r
       if(!stack)\r
-        printf("AliAnaBtag::MakeAnalysisFillHistograms() *** no stack ***: \n");\r
-      \r
+        printf("AliAnaBtag::MakeAnalysisFillHistograms() *** no stack ***: \n");   \r
     }\r
   }// is data and MC\r
   \r
-  ////////////////////////////////////\r
-  //Loop over jets and check for b-tag\r
-  ////////////////////////////////////\r
+\r
   double maxjetEta=-4.;\r
   double maxjetPhi=-4.;\r
   \r
@@ -445,63 +514,76 @@ void  AliAnaBtag::MakeAnalysisFillHistograms()
   if(njets > 0) {\r
     if(GetDebug() > 0) printf("AliAnaBtag::MakeAnalysisFillHistograms() - Jet AOD branch has %d jets.  Performing b-jet tag analysis\n",njets);\r
     \r
+    ///////////////////////////////////Jet loop//////////////////////////////////////////////\r
     for(Int_t ijet = 0; ijet < njets ; ijet++) {\r
       AliAODJet * jet = (AliAODJet*)(GetReader()->GetOutputEvent())->GetJet(ijet) ;\r
-      \r
-      if(ijet==0){\r
-        maxjetEta=jet->Eta();\r
-        maxjetPhi=jet->Phi();\r
-      }\r
-      \r
-      fhJets->Fill(jet->Pt(),1);\r
+      fhJets->Fill(jet->Pt(),1);                                                                                ////////////////FILL\r
       fhJetsAllEtaPhi->Fill(jet->Eta(),jet->Phi());\r
-      \r
       if(jet->Pt() < 0.) continue; //This has to be adjusted depending on pp or AA!\r
-      fhJets->Fill(jet->Pt(),3); //All jets after pt cut\r
       \r
       //Geometric EMCAL cut\r
-      if(TMath::Abs(jet->Eta()) > fJetEtaCut) continue;\r
-      if(jet->Phi() < fJetPhiMin || jet->Phi() > fJetPhiMax) continue;\r
-      fhJets->Fill(jet->Pt(),4); //All jets after geometric cut\r
+      if(TMath::Abs(jet->Eta()) > 0.3) continue;\r
+      if(jet->Phi() < 1.8 || jet->Phi() > 2.9) continue; //This is BAD FIXME \r
+      fhJets->Fill(jet->Pt(),4); //All jets after geometric cut                                                 ////////////////FILL\r
       \r
       Bool_t leadJet  = kFALSE;\r
-      if (ijet==0){ \r
-        fhJets->Fill(jet->Pt(),5); //Leading jets\r
-        leadJet= kTRUE;\r
-      }\r
+      if (ijet==0) {leadJet= kTRUE; fhJets->Fill(jet->Pt(),5);} //Leading jets                                   ////////////////FILL\r
       \r
-      \r
-      Bool_t dvmJet = kFALSE;  \r
+      /////////////////////////Track loop in Jet////////////////////////////////////\r
+      Bool_t dvmJet = kFALSE; \r
+      Bool_t dvmMCJet = kFALSE; \r
       TRefArray* rt = jet->GetRefTracks();\r
       Int_t ntrk = rt->GetEntries();\r
       \r
       for(Int_t itrk = 0; itrk < ntrk; itrk++) {\r
        AliAODTrack* jetTrack = (AliAODTrack*)jet->GetTrack(itrk);\r
-        Bool_t isDVM = CheckIfBjet(jetTrack);\r
-        if(isDVM) dvmJet = kTRUE;\r
-      }\r
-      \r
-      if(dvmJet)\r
-        fhJets->Fill(jet->Pt(),6);\r
-      \r
+       \r
+       Int_t trackId = jetTrack->GetID(); //get the index in the reader\r
+       Int_t naod = GetOutputAODBranch()->GetEntriesFast();\r
+       if(GetDebug() > 3) printf("AliAnaBtag::CheckIfBjet() - aod branch entries %d\n", naod);\r
+       \r
+       for(Int_t iaod = 0; iaod < naod ; iaod++){\r
+         AliAODPWG4Particle* ele =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));\r
+         Int_t electronlabel = ele->GetTrackLabel(0);\r
+         if(electronlabel != trackId) continue;  //skip to the next one if they don't match\r
+         if(ele->GetBtag()>0){\r
+           dvmJet = kTRUE;\r
+           if(ele->GetTag()&(1<<9) || ele->GetTag()&(1<<11) )\r
+             dvmMCJet = kTRUE;\r
+         }\r
+       } //Electron check\r
+       \r
+      }//Track loop of jet tracks\r
       \r
+      if(dvmJet) fhJets->Fill(jet->Pt(),6);                                                                      ////////////////FILL\r
+      //MC stuff\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());\r
-        Bool_t isTrueDjet = IsMcDJet(jet->Eta(), jet->Phi());\r
-        if (isTrueBjet && GetDebug() > 0) printf("== True Bjet==\n");\r
-        if (isTrueDjet && GetDebug() > 0) printf("== True Charm-jet==\n");\r
-        if (dvmJet && GetDebug() > 0)     printf("== found DVM jet==\n");\r
-      }\r
-      \r
-    } //jet loop\r
+        //Bool_t isTrueDjet = IsMcDJet(jet->Eta(), jet->Phi());\r
+       if(dvmJet){\r
+         if(dvmMCJet){\r
+           fhJets->Fill(jet->Pt(),8);   //True                                                                   ////////////////FILL\r
+         }else{\r
+           fhJets->Fill(jet->Pt(),9);   //False                                                                 ////////////////FILL\r
+         }\r
+         if(isTrueBjet){ \r
+           fhJets->Fill(jet->Pt(),10);     //True                                                                 ////////////////FILL\r
+         }else{\r
+           fhJets->Fill(jet->Pt(),11);     //False                                                                ////////////////FILL\r
+         }\r
+       }\r
+       if(isTrueBjet){ \r
+         fhJets->Fill(jet->Pt(),12);     //True                                                                 ////////////////FILL\r
+       }else{\r
+         fhJets->Fill(jet->Pt(),13);     //False                                                                ////////////////FILL\r
+       }\r
+\r
+      }//MC stuff\r
+\r
+    }//jet loop\r
   } //jets exist\r
   \r
-  \r
-  \r
-  \r
-  //Electron loop, read back electrons, fill histos\r
+  //Electron loop, read back electrons, fill histos; mostly photonic shit.\r
   Int_t naod = GetOutputAODBranch()->GetEntriesFast();\r
   if(GetDebug() > 0) printf("AliAnaBtag::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);\r
   \r
@@ -545,33 +627,35 @@ void  AliAnaBtag::MakeAnalysisFillHistograms()
 }\r
 \r
 //__________________________________________________________________\r
-Int_t AliAnaBtag::GetDVMBtag(AliAODTrack * tr )\r
+Int_t AliAnaBtag::GetDVMBtag(AliAODTrack * tr, Int_t &pair, Int_t &start, Int_t &stop)\r
 {\r
   //This method uses the Displaced Vertex between electron-hadron\r
   //pairs and the primary vertex to determine whether an electron is\r
   //likely from a B hadron.\r
-\r
+  Int_t pairstart = fNPair;\r
+  Int_t pairstop = 0;\r
+  Int_t pairn = 0;\r
   Int_t ncls1 = 0;\r
   for(Int_t l = 0; l < 6; l++) if(TESTBIT(tr->GetITSClusterMap(),l)) ncls1++;\r
 \r
-\r
   if (ncls1 < fITSCut) return 0;\r
 \r
   Double_t imp[2] = {-999.,-999.}; Double_t cov[3] = {-999.,-999.,-999.};\r
   Bool_t dcaOkay = GetDCA(tr,imp,cov);  //homegrown dca calculation until AOD is fixed                  \r
   if(!dcaOkay) {\r
-    printf("AliAnaBtag::Problem computing DCA to primary vertex for track %d",tr->GetID());\r
+    printf("AliAnaBtag::GetDVMBtag - Problem computing DCA to primary vertex for track %d",tr->GetID());\r
     return 0;\r
   }\r
 \r
   if (TMath::Abs(imp[0])   > fImpactCut ) return 0;\r
   if (TMath::Abs(imp[1])   > fImpactCut ) return 0;\r
 \r
-//  Int_t nvtx1 = 0;\r
-  Int_t nvtx2 = 0;\r
-//  Int_t nvtx3 = 0;\r
+  Int_t nvtx = 0;\r
+//   Int_t nvtx2 = 0;\r
+//   Int_t nvtx3 = 0;\r
 \r
   for (Int_t k2 =0; k2 < GetAODCTS()->GetEntriesFast() ; k2++) {\r
+\r
     //loop over assoc\r
     AliAODTrack* track2 = (AliAODTrack*)GetAODCTS()->At(k2);\r
     Int_t id1 = tr->GetID();\r
@@ -582,9 +666,13 @@ Int_t AliAnaBtag::GetDVMBtag(AliAODTrack * tr )
     for(Int_t l = 0; l < 6; l++) if(TESTBIT(track2->GetITSClusterMap(),l)) ncls2++;\r
     if (ncls2 < fITSCut) continue;\r
 \r
-    if(tr->Pt()<6.&&track2->Pt() < 0.4) continue;\r
-    if(tr->Pt()>6.&&tr->Pt()<10.&&track2->Pt() < 0.2) continue;\r
-    if(tr->Pt()>10.&&track2->Pt() < 0.6) continue;\r
+\r
+    Double_t sdca=0,pdca=0,minv=0,pth=0,massphoton=0,decaylength=0;\r
+\r
+    sdca = ComputeSignDca(tr,track2,minv,pdca,massphoton,decaylength);\r
+    pth=track2->Pt();\r
+\r
+\r
 \r
 \r
     Double_t dphi = tr->Phi() - track2->Phi();\r
@@ -594,35 +682,37 @@ Int_t AliAnaBtag::GetDVMBtag(AliAODTrack * tr )
     Double_t dr = sqrt(deta*deta + dphi*dphi);\r
 \r
     if(dr > fDrCut) continue;\r
-    \r
-    if(tr->Pt()<6.){\r
-      Double_t sDca = ComputeSignDca(tr, track2, 1.4,0.025);\r
-       if(sDca > 0.06) nvtx2++;\r
-    } \r
-    if(tr->Pt()>6.&&tr->Pt()<10.){\r
-      Double_t sDca = ComputeSignDca(tr, track2, 1.7,0.012);\r
-        if(sDca > 0.03) nvtx2++;\r
-    } \r
-    if(tr->Pt()>10.){\r
-      Double_t sDca = ComputeSignDca(tr, track2, 1.5,0.14);\r
-        if(sDca > 0.04) nvtx2++;\r
-    } \r
+    fNPair++;\r
+    pairn++;\r
+    if(GetDebug() > 0) \r
+      printf("pairs: %d !!!!!!!!!!! \n",fNPair);\r
+    if(fWriteNtuple){\r
+      pairs->Fill(fNPair,fNElec,fEventNumber,pdca,sdca,minv,pth,massphoton,decaylength);\r
+    }\r
+    pairstop=fNPair;\r
+    //pairs->Fill(1,1,1,1,1,1,1);\r
+    fhPairPt->Fill(pth);\r
+    if(decaylength>1.0) continue;\r
+    if(tr->Pt()<6. && pth>0.4 && minv>1.4 && pdca<0.025 && sdca>0.06 && massphoton>0.1) nvtx++; \r
+    if(tr->Pt()>6.&&tr->Pt()<10. && pth>0.2 && minv>1.7 && pdca<0.012 && sdca>0.06 && massphoton>0.1) nvtx++;\r
+    if(tr->Pt()>10.&& pth>0.6 && minv>1.5 && pdca<0.14 && sdca>0.04 && massphoton>0.1) nvtx++;\r
+\r
 \r
 \r
   } //loop over hadrons\r
 \r
-  if(GetDebug() > 0) {\r
-//    if (nvtx1>0) printf("result1 of btagging: %d \n",nvtx1);\r
-    if (nvtx2>0) printf("result2 of btagging: %d \n",nvtx2);\r
-//    if (nvtx3>0) printf("result3 of btagging: %d \n",nvtx3);\r
-  }\r
+  if(GetDebug() > 0) printf("AliAnaBtag::GetDVMBtag - result of btagging: %d \n",nvtx);\r
 \r
-  return nvtx2;\r
 \r
+  pair=pairn;\r
+  start=pairstart+1;\r
+  stop=pairstop;\r
+  if(GetDebug() > 0) printf("End of DVM, pair counts: pairs: %d, start %d, stop %d. \n",pair,start,stop);\r
+  return nvtx;\r
 }\r
 \r
 //__________________________________________________________________\r
-Double_t AliAnaBtag::ComputeSignDca(AliAODTrack *tr, AliAODTrack *tr2 , float masscut, double pdcacut)\r
+Double_t AliAnaBtag::ComputeSignDca(AliAODTrack *tr, AliAODTrack *tr2 , Double_t &masscut, Double_t &pdcacut, Double_t &massphoton, Double_t &decay)\r
 {\r
   //Compute the signed dca between two tracks\r
   //and return the result\r
@@ -633,19 +723,14 @@ Double_t AliAnaBtag::ComputeSignDca(AliAODTrack *tr, AliAODTrack *tr2 , float ma
   //=====Now calculate DCA between both tracks=======  \r
   Double_t massE = 0.000511;\r
   Double_t massK = 0.493677;\r
-\r
   Double_t vertex[3] = {-999.,-999.,-999}; //vertex\r
+\r
   if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) {\r
-    GetVertex(vertex); //If only one file, get the vertex from there\r
-    //FIXME:  Add a check for whether file 2 is PYTHIA or HIJING\r
-    //If PYTHIA, then set the vertex from file 2, if not, use the\r
-    //vertex from file 1\r
-    //if(GetReader()->GetSecondInputAODTree()) GetReader()->GetSecondInputAODVertex(vertex);\r
+    GetVertex(vertex); //If only one file, get the vertex from there  \r
   }\r
   \r
-  TVector3 primV(vertex[0],vertex[1],vertex[2]) ;\r
-\r
-  if(GetDebug() > 5) printf(">>ComputeSdca:: primary vertex = %2.2f,%2.2f,%2.2f \n",vertex[0],vertex[1],vertex[2]) ;\r
+  TVector3 primV(vertex[0],vertex[1],vertex[2]) ; \r
+  if(GetDebug()>0) printf(">>ComputeSdca:: primary vertex = %2.2f,%2.2f,%2.2f \n",vertex[0],vertex[1],vertex[2]) ;\r
 \r
   AliExternalTrackParam *param1 = new AliExternalTrackParam(tr);\r
   AliExternalTrackParam *param2 = new AliExternalTrackParam(tr2);\r
@@ -656,7 +741,6 @@ Double_t AliAnaBtag::ComputeSignDca(AliAODTrack *tr, AliAODTrack *tr2 , float ma
   bfield[1]=0;\r
   bfield[2]=5.;\r
   Double_t bz = 5.; // = param1->GetBz();\r
-  //cout<<"In ComputeSignDCA, bfield[3] and bz: "<<bfield[0]<<" "<<bfield[1]<<" "<<bfield[2]<<" "<<bz<<endl;\r
   Double_t xplane1 = 0.; Double_t xplane2 = 0.;\r
   Double_t pairdca = param1->GetDCA(param2,bz,xplane1,xplane2);\r
 \r
@@ -678,17 +762,14 @@ Double_t AliAnaBtag::ComputeSignDca(AliAODTrack *tr, AliAODTrack *tr2 , float ma
   TVector3 decayvector(0,0,0);\r
   decayvector = secvtxpt - primV; //decay vector from PrimVtx\r
   Double_t decaylength = decayvector.Mag();\r
-\r
-  //printf("\t JLK pairDCA = %2.2f\n",pairdca);\r
+  decay=decaylength;\r
 \r
   if(GetDebug() > 0) {\r
     printf(">>ComputeSdca:: mom1=%f, mom2=%f \n", emomAtB.Perp(), hmomAtB.Perp() );\r
     printf(">>ComputeSdca:: pairDCA=%f, length=%f \n", pairdca,decaylength );\r
   }\r
 \r
-\r
-\r
-  if (emomAtB.Mag()>0 && pairdca < pdcacut && decaylength < fDecayLenCut ) {\r
+  if (emomAtB.Mag()>0 /*&& decaylength < fDecayLenCut*/ ) {\r
     TVector3 sumMom = emomAtB+hmomAtB;\r
     Double_t ener1 = sqrt(pow(emomAtB.Mag(),2) + massE*massE);\r
     Double_t ener2 = sqrt(pow(hmomAtB.Mag(),2) + massK*massK);\r
@@ -696,19 +777,16 @@ Double_t AliAnaBtag::ComputeSignDca(AliAODTrack *tr, AliAODTrack *tr2 , float ma
     Double_t mass = sqrt(pow((ener1+ener2),2) - pow(sumMom.Mag(),2));\r
     Double_t massPhot = sqrt(pow((ener1+ener3),2) - pow(sumMom.Mag(),2));\r
     Double_t sDca = decayvector.Dot(emomAtB)/emomAtB.Mag();\r
-\r
-//     if (masscut<1.1) fhDVMBtagQA2->Fill(sDca, mass);\r
-\r
-    if (mass > masscut && massPhot > 0.1) signDca = sDca;\r
-    \r
+    pdcacut=pairdca;//\r
+    masscut=mass; // \r
+    massphoton=massPhot;//                                                                     Send it back!\r
+    signDca = sDca;\r
     if(GetDebug() > 0) printf("\t>>ComputeSdca:: mass=%f \n", mass);\r
     if(GetDebug() > 0) printf("\t>>ComputeSdca:: sec vtx-signdca :%f\n",signDca);\r
   }\r
-\r
   //clean up\r
   delete param1;\r
   delete param2;\r
-\r
   return signDca;\r
 }\r
 \r
@@ -719,9 +797,7 @@ Bool_t AliAnaBtag::PhotonicV0(Int_t id)
   //This method checks to see whether a track that has been flagged as\r
   //an electron was determined to match to a V0 candidate with\r
   //invariant mass consistent with photon conversion\r
-\r
   Bool_t itIS = kFALSE;\r
-\r
   Double_t massEta = 0.547;\r
   Double_t massRho0 = 0.770;\r
   Double_t massOmega = 0.782;\r
@@ -751,7 +827,6 @@ Bool_t AliAnaBtag::PhotonicV0(Int_t id)
       }\r
     } }\r
   return itIS;\r
-\r
 }\r
 \r
 //__________________________________________________________________\r
@@ -759,8 +834,6 @@ Bool_t AliAnaBtag::GetDCA(const AliAODTrack* track,Double_t impPar[2], Double_t
 {\r
   //Use the Event vertex and AOD track information to get\r
   //a real impact parameter for the track\r
-\r
-\r
   Double_t maxD = 100000.; //max transverse IP\r
   if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) {\r
     AliVEvent* ve = (AliVEvent*)GetReader()->GetInputEvent();\r
@@ -775,9 +848,7 @@ Bool_t AliAnaBtag::GetDCA(const AliAODTrack* track,Double_t impPar[2], Double_t
     Bool_t gotit = esdTrack.PropagateToDCABxByBz(vv,bfield,maxD,impPar,cov);\r
     return gotit;\r
   }\r
-\r
   return kFALSE;\r
-\r
 }\r
 //__________________________________________________________________\r
 Bool_t AliAnaBtag::CheckIfBjet(const AliAODTrack* track)\r
@@ -798,46 +869,36 @@ Bool_t AliAnaBtag::CheckIfBjet(const AliAODTrack* track)
   return bjet;\r
 } \r
 \r
-\r
-\r
 //__________________________________________________________________\r
 AliAODMCParticle* AliAnaBtag::GetMCParticle(Int_t ipart) \r
 {\r
   //Get the MC particle at position ipart\r
-  \r
+\r
   AliAODMCParticle* aodprimary = 0x0;\r
   TClonesArray * mcparticles0 = 0x0;\r
-  //TClonesArray * mcparticles1 = 0x0;\r
-  \r
+  TClonesArray * mcparticles1 = 0x0;\r
+\r
   if(GetReader()->ReadAODMCParticles()){\r
     //Get the list of MC particles                                                                                                                           \r
     mcparticles0 = GetReader()->GetAODMCParticles(0);\r
-    if(!mcparticles0) {\r
-      if(GetDebug() > 0)printf("AliAnaBtag::MakeAnalysisFillHistograms() -  Standard MCParticles not available!\n");\r
+    if(!mcparticles0 && GetDebug() > 0) {\r
+      printf("AliAnaBtag::MakeAnalysisFillHistograms() -  Standard MCParticles not available!\n");\r
     }\r
-    //    if(GetReader()->GetSecondInputAODTree()){\r
-    //      mcparticles1 = GetReader()->GetAODMCParticles(1);\r
-    //      if(!mcparticles1 && GetDebug() > 0) {\r
-    // printf("AliAnaBtag::MakeAnalysisFillHistograms() -  Second input MCParticles not available!\n");\r
-    //      }\r
-    //    }\r
-    \r
-    else{\r
-      Int_t npart0 = mcparticles0->GetEntriesFast();\r
-      //Int_t npart1 = 0;\r
-      //if(mcparticles1) npart1 = mcparticles1->GetEntriesFast();\r
-      if(ipart < npart0) aodprimary = (AliAODMCParticle*)mcparticles0->At(ipart);\r
-      //else aodprimary = (AliAODMCParticle*)mcparticles1->At(ipart-npart0);\r
-      if(!aodprimary) {\r
-        printf("AliAnaBtag::GetMCParticle() *** no primary ***:  label %d \n", ipart);\r
-        return 0x0;\r
-      }\r
+\r
+    Int_t npart0 = mcparticles0->GetEntriesFast();\r
+    Int_t npart1 = 0;\r
+    if(mcparticles1) npart1 = mcparticles1->GetEntriesFast();\r
+    if(ipart < npart0) aodprimary = (AliAODMCParticle*)mcparticles0->At(ipart);\r
+    else aodprimary = (AliAODMCParticle*)mcparticles1->At(ipart-npart0);\r
+    if(!aodprimary) {\r
+      printf("AliAnaBtag::GetMCParticle() *** no primary ***:  label %d \n", ipart);\r
+      return 0x0;\r
     }\r
+\r
   } else {\r
     printf("AliAnaBtag::GetMCParticle() - Asked for AliAODMCParticle but we have a stack reader.\n");\r
   }\r
   return aodprimary;\r
-  \r
 }\r
 \r
 //__________________________________________________________________\r
@@ -846,7 +907,6 @@ Bool_t  AliAnaBtag::IsMcBJet(Double_t jeta, Double_t jphi)
   //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
   AliStack* stack = 0x0;\r
   \r
   for(Int_t ipart = 0; ipart < 100; ipart++) {\r
@@ -873,20 +933,16 @@ Bool_t  AliAnaBtag::IsMcBJet(Double_t jeta, Double_t jphi)
       peta = aodprimary->Eta();\r
     }\r
     if ( TMath::Abs(pdg) != 5) continue;\r
-    \r
-    //      printf("MTH: IsMcBJet : %d, pdg=%d : pt=%f \n", ipart, pdgcode, primary->Pt());\r
     Double_t dphi = jphi - pphi;\r
     Double_t deta = jeta - peta;\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
@@ -922,60 +978,21 @@ Bool_t  AliAnaBtag::IsMcDJet(Double_t jeta, Double_t jphi)
     }\r
 \r
     if ( TMath::Abs(pdg) != 4) continue;\r
-\r
     Double_t dphi = jphi - pphi;\r
     Double_t deta = jeta - peta;\r
     Double_t dr = sqrt(deta*deta + dphi*dphi);\r
-    \r
+   \r
     if (dr < 0.2) {\r
       cjet=kTRUE;\r
       break;\r
     }\r
   }\r
-\r
   return cjet;\r
-\r
 }\r
 \r
-//__________________________________________________________________\r
-void AliAnaBtag::Print(const Option_t * opt) const\r
-{\r
-  //Print some relevant parameters set for the analysis\r
-  \r
-  if(! opt)\r
-    return;\r
-  \r
-  printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;\r
-  AliAnaPartCorrBaseClass::Print(" ");\r
-\r
-  printf("Calorimeter            =     %s\n", fCalorimeter.Data()) ;\r
-  printf("pOverE range           =     %f - %f\n",fpOverEmin,fpOverEmax);\r
-  printf("residual cut           =     %f\n",fResidualCut);\r
-  printf("---DVM Btagging\n");\r
-  printf("max IP-cut (e,h)       =     %f\n",fImpactCut);\r
-  printf("min ITS-hits           =     %d\n",fITSCut);\r
-  printf("max dR (e,h)           =     %f\n",fDrCut);\r
-  printf("max pairDCA            =     %f\n",fPairDcaCut);\r
-  printf("max decaylength        =     %f\n",fDecayLenCut);\r
-  printf("min Associated Pt      =     %f\n",fAssocPtCut);\r
-  printf("---IPSig Btagging\n");\r
-  printf("min tag track          =     %d\n",fNTagTrkCut);\r
-  printf("min IP significance    =     %f\n",fIPSigCut);\r
-  printf("    \n") ;\r
-       \r
-} \r
-\r
 \r
 //__________________________________________________________________\r
 void  AliAnaBtag::Terminate(TList* outputList)\r
 {\r
\r
-  //Do some plots to end\r
-  //Recover histograms from output histograms list, needed for\r
-  //distributed analysis.                \r
-  //ReadHistograms(outputList);\r
-\r
-  printf(" AliAnaBtag::Terminate()  *** %s Report: %d outputs\n", GetName(), outputList->GetEntries()) ;\r
-\r
+  printf(" AliAnaBtag::Terminate()  *** %s Report: %d outputs/histograms \n", GetName(), outputList->GetEntries()) ;\r
 }\r
-\r
index 6035bb1..4101e83 100644 (file)
@@ -11,7 +11,7 @@
 // and kept in the AOD. Few histograms produced.\r
 //\r
 \r
-//-- Author T.R.P-R.Aronsson, J.Klay\r
+//-- Author T.R.P-R.Aronsson\r
 \r
 // --- ROOT system ---\r
 class TH2F ;\r
@@ -46,7 +46,7 @@ public:
   void MakeAnalysisFillHistograms() ; \r
   \r
   //B-tagging\r
-  Int_t GetDVMBtag(AliAODTrack * tr);//Main tagger\r
+  Int_t GetDVMBtag(AliAODTrack * tr, Int_t &pairs, Int_t &start, Int_t &stop);//Main tagger\r
 \r
   //Temporary local method to get DCA\r
   Bool_t GetDCA(const AliAODTrack* tr,Double_t imp[2], Double_t cov[3]);\r
@@ -57,11 +57,7 @@ public:
   Bool_t IsMcBJet(Double_t x, Double_t y);\r
   Bool_t IsMcDJet(Double_t x, Double_t y);\r
 \r
-  void Print(const Option_t * opt)const;\r
-  \r
-  TString GetCalorimeter()   const {return fCalorimeter ; }\r
-  Double_t GetpOverEmin()   const {return fpOverEmin ; }\r
-  Double_t GetpOverEmax()   const {return fpOverEmax ; }\r
+\r
 \r
   Double_t GetDrCut() const { return fDrCut; }\r
   Double_t GetPairDcaCut() const { return fPairDcaCut; }\r
@@ -73,13 +69,9 @@ public:
   Int_t    GetITSCut() const { return fITSCut; }\r
   Int_t    GetNTagTrackCut() const { return fNTagTrkCut; }\r
   Double_t GetIPSigCut() const { return fIPSigCut; }\r
-  Double_t GetMinClusEne() const { return fMinClusEne; }\r
 \r
-  void SetCalorimeter(TString det)    {fCalorimeter = det ; }\r
-  void SetpOverEmin(Double_t min)     {fpOverEmin = min ; }\r
-  void SetpOverEmax(Double_t max)     {fpOverEmax = max ; }\r
-  void SetResidualCut(Double_t cut)     {fResidualCut = cut ; }\r
 \r
+  void SetWriteNtuple(Int_t w) { fWriteNtuple = w; }\r
   void SetDrCut(Double_t dr)  { fDrCut = dr; }\r
   void SetPairDcaCut(Double_t pdca) { fPairDcaCut = pdca; }\r
   void SetDecayLenCut(Double_t dlen) { fDecayLenCut = dlen; }\r
@@ -90,25 +82,28 @@ public:
   void SetITSCut(Int_t its) { fITSCut = its; }\r
   void SetNTagTrackCut(Int_t ntr) { fNTagTrkCut = ntr; }\r
   void SetIPSigCut(Double_t ips) { fIPSigCut = ips; }\r
-  void SetMinClusEne(Double_t ene) { fMinClusEne = ene; }\r
+\r
 \r
   void InitParameters();\r
   void Terminate(TList * outputList);\r
          \r
   private:\r
   //For DVM B-tag method\r
-  Double_t ComputeSignDca(AliAODTrack *track, AliAODTrack *track2, float cut1, double pdcacut);\r
+  Double_t ComputeSignDca(AliAODTrack *tr, AliAODTrack *tr2 , Double_t &masscut, Double_t &pdcacut, Double_t &massphoton, Double_t &decay);\r
 \r
-  //Int_t GetMCSource(Int_t mctag);\r
   AliAODMCParticle* GetMCParticle(Int_t part);\r
 \r
 \r
   private:\r
-  TString  fCalorimeter;  //! Which detector? EMCAL or PHOS\r
-  Double_t fpOverEmin;    //! Minimum p/E value for Electrons\r
-  Double_t fpOverEmax;    //! Maximum p/E value for Electrons\r
-  Double_t fResidualCut;  //! Track-cluster matching distance\r
-  Double_t fMinClusEne;   //! Min clus energy for matching\r
+  //NTuples!\r
+  Int_t fWriteNtuple;    //Will always be no, but might be set to yes in config file.\r
+  TNtuple * electrons;   //Electrons\r
+  TNtuple * pairs;       //Pairs to the electrons\r
+  TNtuple * events;\r
+  Int_t fEventNumber;    // For Ntuple to label events (starts at 0)\r
+  Int_t fNElec;\r
+  Int_t fNElecEv;\r
+  Int_t fNPair;\r
 \r
   //DVM B-tagging\r
   Double_t fDrCut;       //max dR\r
@@ -133,16 +128,34 @@ public:
   TH1F * fhEmcalElectrons;              //All electrons, as id:ed by EMCAL\r
   TH1F * fhTRDElectrons;                //Electrons from TRD\r
   TH1F * fhTPCElectrons;                //Electrons from TPC\r
+  TH1F * fhEmcalMCE;                    //All electrons, as id:ed by EMCAL MC\r
+  TH1F * fhTRDMCE;                      //Electrons from TRD MC\r
+  TH1F * fhTPCMCE;                      //Electrons from TPC MC\r
+  TH1F * fhEmcalMCEFake;                //All electrons, as id:ed by EMCAL MC, fake\r
+  TH1F * fhTRDMCEFake;                  //Electrons from TRD MC, fake\r
+  TH1F * fhTPCMCEFake;                  //Electrons from TPC MC, fake\r
+  TH1F * fhEmcalMCP;                    //Pions, as id:ed by EMCAL\r
+  TH1F * fhTRDMCP;                      //Pions from TRD\r
+  TH1F * fhTPCMCP;                      //Pions from TPC\r
+  TH1F * fhEmcalMCK;                    //Kaons from EMCAL\r
+  TH1F * fhTRDMCK;                      //Kaons from TRD\r
+  TH1F * fhTPCMCK;                      //Kaons from TPC\r
+  TH1F * fhSpecies;                     //PDG of id:ed electrons\r
   TH1F * fhDVM1;                        //initial b-tag dvm1\r
   TH1F * fhDVM2;                        //initial b-tag dvm2 \r
+  TH1F * fhNVTX;                        //NVtx of all btags \r
+  TH1F * fhNVTXMC;                        //NVtx of MC btags \r
   TH2F * fhJets;                        //All jets 2d\r
   TH2F * fhJetsAllEtaPhi;               //Eta phi for all jets\r
   TH2F * fhJetsLeadingBElectronEtaPhi;  //deta dphi btw leading jet and bele\r
   TH2F * fhDVM1EtaPhi;                  //eta phi for b-electrons\r
   TH2F * fhBJetElectronDetaDphi;        //eta phi for jet with b-electrons\r
+  TH2F * fhClusterMap;                  //2D eta-phi of EMCAL clusters\r
   TH1F * fhClusterEnergy;               //cluster E of EMCAL\r
   TH1F * fhTestalle;                                   //Temp histo for EMCAL cluster energy\r
   TH1F * fhResidual;                    //Residuals from trackmatching\r
+  TH1F * fhPairPt;                      //Pairs\r
+\r
 \r
   //Analysis of electrons\r
   TH2F * fhElectrons;\r
index 51f54ab..196d052 100755 (executable)
@@ -696,6 +696,12 @@ void  AliAnaElectron::MakeAnalysisFillAOD()
   for (Int_t itrk =  0; itrk <  ntracks; itrk++) {////////////// track loop\r
     iCluster = -999; //start with no match\r
     AliAODTrack * track = (AliAODTrack*) (GetAODCTS()->At(itrk)) ;\r
+    //Added negative track condition. Seems that about 3-4% of the tracks have negative label. Causes crashes sometimes.\r
+    if(track->GetLabel()<0){\r
+      printf("Negative track label! Not sure what it  means, abort track. \n");\r
+      continue;\r
+    }\r
+\r
     if (TMath::Abs(track->Eta())< 0.5) refmult++;\r
     Double_t imp[2] = {-999.,-999.}; Double_t cov[3] = {-999.,-999.,-999.};\r
     Bool_t dcaOkay = GetDCA(track,imp,cov);  //homegrown dca calculation until AOD is fixed\r
@@ -728,8 +734,9 @@ void  AliAnaElectron::MakeAnalysisFillAOD()
       \r
       //TLorentzVector mom2(mom,0.);\r
       Bool_t in = kFALSE;\r
-      if(mom.Phi()*180./TMath::Pi() > 80. && mom.Phi()*180./TMath::Pi() < 190. &&\r
-        mom.Eta() > -0.7 && mom.Eta() < 0.7) in = kTRUE;\r
+      //Removed two lines below. Seems wrong to demand that the momentum vector should point to the EMCAL.\r
+      //if(mom.Phi()*180./TMath::Pi() > 80. && mom.Phi()*180./TMath::Pi() < 190. &&\r
+      // mom.Eta() > -0.7 && mom.Eta() < 0.7) in = kTRUE;\r
       //Also check the track\r
       if(track->Phi()*180./TMath::Pi() > 80. && track->Phi()*180./TMath::Pi() < 190. &&\r
         track->Eta() > -0.7 && track->Eta() < 0.7) in = kTRUE;\r
@@ -885,7 +892,6 @@ void  AliAnaElectron::MakeAnalysisFillAOD()
          //B-tagging\r
          if(GetDebug() > 1) printf("Found Electron - do b-tagging\n");\r
          Int_t dvmbtag = GetDVMBtag(track); bt += dvmbtag;\r
-\r
          fh2EledEdx->Fill(track->P(),dEdx);\r
          \r
          Double_t eMass = 0.511/1000; //mass in GeV\r
@@ -1317,7 +1323,6 @@ Int_t AliAnaElectron::GetDVMBtag(AliAODTrack * tr )
 \r
   fhDVMBtagQA3->Fill(ncls1);\r
   if (ncls1 < fITSCut) return 0;\r
-\r
   Double_t imp[2] = {-999.,-999.}; Double_t cov[3] = {-999.,-999.,-999.};\r
   Bool_t dcaOkay = GetDCA(tr,imp,cov);  //homegrown dca calculation until AOD is fixed                  \r
   if(!dcaOkay) {\r
@@ -1354,7 +1359,6 @@ Int_t AliAnaElectron::GetDVMBtag(AliAODTrack * tr )
     Double_t dr = sqrt(deta*deta + dphi*dphi);\r
 \r
     if(dr > fDrCut) continue;\r
-    \r
     Double_t sDca1 = ComputeSignDca(tr, track2, 1.0);\r
     if (sDca1 > fSdcaCut) nvtx1++;\r
     Double_t sDca2 = ComputeSignDca(tr, track2, 1.5);\r
@@ -1408,10 +1412,14 @@ Double_t AliAnaElectron::ComputeSignDca(AliAODTrack *tr, AliAODTrack *tr2 , floa
   AliExternalTrackParam *param1 = new AliExternalTrackParam(tr);\r
   AliExternalTrackParam *param2 = new AliExternalTrackParam(tr2);\r
 \r
+  //Replaced functions to get the B-field. They do not work.\r
   Double_t bfield[3];\r
   param1->GetBxByBz(bfield);\r
-  Double_t bz = param1->GetBz();\r
-\r
+  bfield[0]=0;\r
+  bfield[1]=0;\r
+  bfield[2]=5.0;\r
+  Double_t bz = 5.0; //param1->GetBz();\r
+  printf("ZOMG1 bfield is: %f, %f",bz,bfield[2]);\r
   Double_t xplane1 = 0.; Double_t xplane2 = 0.;\r
   Double_t pairdca = param1->GetDCA(param2,bz,xplane1,xplane2);\r
 \r
@@ -1434,7 +1442,7 @@ Double_t AliAnaElectron::ComputeSignDca(AliAODTrack *tr, AliAODTrack *tr2 , floa
   decayvector = secvtxpt - primV; //decay vector from PrimVtx\r
   Double_t decaylength = decayvector.Mag();\r
 \r
-  printf("\t JLK pairDCA = %2.2f\n",pairdca);\r
+  //printf("\t JLK pairDCA = %2.2f\n",pairdca);\r
 \r
   if(GetDebug() > 0) {\r
     printf(">>ComputeSdca:: mom1=%f, mom2=%f \n", emomAtB.Perp(), hmomAtB.Perp() );\r
@@ -1500,8 +1508,13 @@ Double_t AliAnaElectron::GetIPSignificance(AliAODTrack *tr, Double_t jetPhi)
   AliVTrack* vTrack = (AliVTrack*)vEvent->GetTrack(trackIndex);\r
   if(!vTrack) return -99;\r
   AliESDtrack esdTrack(vTrack);\r
+  //Replaced function to get B-field, it does not work\r
   Double_t bfield[3];\r
-  esdTrack.GetBxByBz(bfield);\r
+  //esdTrack.GetBxByBz(bfield);\r
+  bfield[0]=0;\r
+  bfield[1]=0;\r
+  bfield[2]=5.;\r
+  //printf("ZOMG2 ESD bfield is: %f",bfield[2]);\r
   if(!esdTrack.PropagateToDCABxByBz(vv, bfield, maxD, impPar, ipCov)) return -100;\r
   if(ipCov[0]<0) return -101;\r
 \r
@@ -1513,7 +1526,7 @@ Double_t AliAnaElectron::GetIPSignificance(AliAODTrack *tr, Double_t jetPhi)
   Double_t cosTheta = TMath::Cos(jetPhi - phiIP);\r
   Double_t sign = cosTheta/TMath::Abs(cosTheta);\r
   significance = TMath::Abs(impPar[0])/TMath::Sqrt(ipCov[0])*sign;\r
-  printf("\t JLK significance = %2.2f\n",significance);\r
+  //printf("\t JLK significance = %2.2f\n",significance);\r
   //ip = fabs(impPar[0]);\r
   fhIPSigBtagQA2->Fill(significance);\r
   return significance;\r
@@ -1687,9 +1700,14 @@ Bool_t AliAnaElectron::GetDCA(const AliAODTrack* track,Double_t impPar[2], Doubl
     AliVVertex *vv = (AliVVertex*)ve->GetPrimaryVertex();\r
     AliESDtrack esdTrack(track);\r
     Double_t bfield[3];\r
-    esdTrack.GetBxByBz(bfield);\r
+    //This GetBxByBz function does not work properly.\r
+    //esdTrack.GetBxByBz(bfield);\r
+    bfield[0]=0;\r
+    bfield[1]=0;\r
+    bfield[2]=5.;\r
+    //printf("ZOMG2 ESD bfield is: %f",bfield[2]);\r
     Bool_t gotit = esdTrack.PropagateToDCABxByBz(vv,bfield,maxD,impPar,cov);\r
-    printf("\t JLK impPar = %2.2f\n",impPar[0]);\r
+    //printf("\t JLK impPar = %2.2f\n",impPar[0]);\r
     return gotit;\r
   }\r
 \r
similarity index 57%
rename from PWG4/macros/ConfigAnalysisBTag.C
rename to PWG4/macros/ConfigAnalysisBtag.C
index 8d177f2..f7ac734 100644 (file)
@@ -1,9 +1,7 @@
 //------------------------------------\r
-// Configuration macro example:\r
+// Config file for B-tagging\r
 //\r
-// Configure EMCal electron analysis\r
-//\r
-// Modified by: K. Read\r
+// Author T. Aronsson\r
 //\r
 //------------------------------------\r
 \r
@@ -13,7 +11,7 @@ AliAnaPartCorrMaker*  ConfigAnalysis()
   // Configuration goes here\r
   // \r
   printf("======================== \n");\r
-  printf("ConfigAnalysisElectron() \n");\r
+  printf("==Preforming Btag.cxx=== \n");\r
   printf("======================== \n");\r
   Bool_t kInputIsESD = kTRUE;     //uncomment for input ESD\r
 //Bool_t kInputIsESD = kFALSE;    //uncomment for input AODs\r
@@ -21,32 +19,7 @@ AliAnaPartCorrMaker*  ConfigAnalysis()
 //Bool_t kFollowsFilter = kFALSE; //uncomment if no ESD filter task\r
   Bool_t kMC = kTRUE; //set to kFALSE for data\r
 \r
-  //enum for the different electron cut sets\r
-  //defined for dR and p/E\r
-  //kTight2 is the default standard cuts\r
-  enum kCutTypes {kTight2, kLooseTight, kTightLoose, kLoose2};\r
-  Int_t kCutSet = kTight2;\r
-  Double_t pOverEmin = 0.8;  //tight\r
-  Double_t pOverEmax = 1.2;  //tight\r
-  Double_t dRmax     = 0.02; //tight\r
-  if (gSystem->Getenv("ELECUTSET")){\r
-    kCutSet = atoi(gSystem->Getenv("ELECUTSET"));\r
-  }\r
-  if(kCutSet == kLooseTight) {\r
-    pOverEmin = 0.6;  //loose\r
-    pOverEmax = 1.4;  //loose\r
-    dRmax     = 0.02; //tight\r
-  }\r
-  if(kCutSet == kTightLoose) {\r
-    pOverEmin = 0.8;  //tight\r
-    pOverEmax = 1.2;  //tight\r
-    dRmax     = 0.05; //loose\r
-  }\r
-  if(kCutSet == kLoose2) {\r
-    pOverEmin = 0.6;  //loose\r
-    pOverEmax = 1.4;  //loose\r
-    dRmax     = 0.05; //loose\r
-  }    \r
+   \r
 \r
   //Alternatively, select input via anaInputData environment variable.\r
   if (gSystem->Getenv("anaInputData")){\r
@@ -80,13 +53,12 @@ AliAnaPartCorrMaker*  ConfigAnalysis()
   if(kInputIsESD && !kFollowsFilter)AliCaloTrackESDReader *reader = new AliCaloTrackESDReader();\r
   else           AliCaloTrackAODReader *reader = new AliCaloTrackAODReader();\r
   reader->SetDebug(-1);//10 for lots of messages\r
-\r
+  reader->SwitchOnWriteDeltaAOD();\r
   //Switch on or off the detectors information that you want\r
   reader->SwitchOnEMCAL();\r
-  reader->SwitchOnCTS();\r
-  //reader->SwitchOffEMCALCells();     \r
+  reader->SwitchOnCTS();       \r
   reader->SwitchOffPHOS();\r
-  //reader->SwitchOffPHOSCells();      \r
+       \r
 \r
   //Kine\r
   if(kMC && !kInputIsESD){\r
@@ -102,25 +74,12 @@ AliAnaPartCorrMaker*  ConfigAnalysis()
   reader->SetEMCALPtMin(0.0); //new\r
   if(kFollowsFilter)reader->SetTrackStatus(0);  //to prevent automatic TPC and ITS refit\r
 \r
-  //In case of generating jet events (with PYTHIA), if pt hard bin is small\r
-  //reject events with large difference between ptHard and triggered jet\r
-  //reader->SetPtHardAndJetPtComparison(kTRUE);\r
+\r
 \r
   reader->SetFiducialCut(fidCut);\r
 \r
   if(!kInputIsESD){\r
-    // Analysis with tracks, select only tracks with\r
-    // following bits\r
-\r
-    //     //We want tracks fitted in the detectors:\r
-    //     ULong_t status=AliAODTrack::kTPCrefit;\r
-    //     status|=AliAODTrack::kITSrefit;\r
-   \r
-    //     We want tracks whose PID bit is set:\r
-    //     ULong_t status =AliAODTrack::kITSpid;\r
-    //     status|=AliAODTrack::kTPCpid;       \r
 \r
-    // reader->SetTrackStatus(status);\r
   }\r
 \r
   reader->Print("");\r
@@ -136,36 +95,36 @@ AliAnaPartCorrMaker*  ConfigAnalysis()
 \r
   fidCut2->Print("");\r
 \r
+\r
+\r
+\r
+\r
   //---------------------------------------------------------------------\r
   // Analysis algorithm\r
   //---------------------------------------------------------------------\r
 \r
-  AliAnaBtag *anaelectron = new AliAnaBtag();\r
-  anaelectron->SetDebug(-1); //10 for lots of messages\r
-  anaelectron->SetCalorimeter("EMCAL");\r
+  AliAnaBtag *btag = new AliAnaBtag();\r
+\r
+ //Base class\r
+  btag->SetDebug(-1); //10 for lots of messages\r
+  //btag->SetWriteNtuple(1); //Can be used to write out NTuples for local analysis (1000 times faster than AOD analysis), default is off.\r
   if(kMC){\r
-    anaelectron->SwitchOffDataMC();\r
-    anaelectron->SetMinPt(1.);\r
+    btag->SwitchOnDataMC();\r
+    btag->SetMinPt(1.);\r
   }\r
-  anaelectron->SetOutputAODName("Electrons");\r
-  anaelectron->SetOutputAODClassName("AliAODPWG4Particle");\r
-  anaelectron->SetWriteNtuple(kFALSE);\r
-  //Determine which cuts to use based on enum\r
-  anaelectron->SetpOverEmin(pOverEmin);\r
-  anaelectron->SetpOverEmax(pOverEmax);\r
-  anaelectron->SetResidualCut(dRmax);\r
-  //Set Histrograms bins and ranges\r
-  anaelectron->SetHistoPtRangeAndNBins(0, 100, 100) ;\r
-  anaelectron->SetHistoPhiRangeAndNBins(0, TMath::TwoPi(), 100) ;\r
-  anaelectron->SetHistoEtaRangeAndNBins(-0.7, 0.7, 100) ;\r
-  anaelectron->Print("");\r
+  btag->SetOutputAODName("Electrons");\r
+  btag->SetOutputAODClassName("AliAODPWG4Particle");\r
+\r
+\r
+\r
+\r
 \r
   //---------------------------------------------------------------------\r
   // Set  analysis algorithm and reader\r
   //---------------------------------------------------------------------\r
   maker = new AliAnaPartCorrMaker();\r
   maker->SetReader(reader);//pointer to reader\r
-  maker->AddAnalysis(anaelectron,0);\r
+  maker->AddAnalysis(btag,0);\r
   maker->SetAnaDebug(-1)  ;\r
   maker->SwitchOnHistogramsMaker()  ;\r
   maker->SwitchOnAODsMaker()  ;\r
diff --git a/PWG4/macros/electrons/makeparfilesproperly b/PWG4/macros/electrons/makeparfilesproperly
new file mode 100644 (file)
index 0000000..d6f99f9
--- /dev/null
@@ -0,0 +1,28 @@
+#! /bin/bash
+
+echo "Making parfiles..."
+
+source /home/aron/root/setup.sh
+
+/home/aron/root/alice/AliRoot/PWG4/macros/electrons/./parmaker all remote
+
+tar xvf EMCALUtils.par
+
+echo "...copying"
+#cp /home/aron/root/alice/AliRoot/EMCAL/AliCaloPeakFinderConstants.h EMCALUtils/.
+cp -v /home/aron/root/alice/AliRoot/EMCAL/AliCaloConstants.h EMCALUtils/.
+
+echo "...removing"
+rm -v EMCALUtils.par
+
+echo "...remaking EMCALUtils.par"
+/home/aron/root/alice/AliRoot/PWG4/macros/electrons/./parmaker EMCALUtils
+
+echo "...removing EMCALUtils dir"
+rm -rfv EMCALUtils
+
+echo "...making directory GRID and copying important stuff there"
+mkdir GRID
+
+cp -v *.par parmaker anaJete.C run.C ConfigAnalysisElectron.C ConfigJetAnalysisFastJet.C anaJete.jdl validate.sh mergeout.jdl GRID/.
+
index bf04e85..077ae2e 100755 (executable)
@@ -25,17 +25,17 @@ case $2 in
 
     case $1 in
       "all")
-         parmaker ANALYSIS remote
-         parmaker ANALYSISalice remote
-         parmaker AOD remote
-         parmaker EMCALUtils remote
-         parmaker ESD remote
-         parmaker FASTJETAN remote
-         parmaker JETAN remote
-         parmaker PHOSUtils remote
-         parmaker PWG4PartCorrBase remote
-         parmaker PWG4PartCorrDep remote
-         parmaker STEERBase remote
+         ./parmaker ANALYSIS remote
+         ./parmaker ANALYSISalice remote
+         ./parmaker AOD remote
+         ./parmaker EMCALUtils remote
+         ./parmaker ESD remote
+         ./parmaker FASTJETAN remote
+         ./parmaker JETAN remote
+         ./parmaker PHOSUtils remote
+         ./parmaker PWG4PartCorrBase remote
+         ./parmaker PWG4PartCorrDep remote
+         ./parmaker STEERBase remote
          echo " "
          echo "finished making complete set of par files"
          exit