PWG4/PartCorrDep/AliAnaBtag.cxx

   1 /**************************************************************************\r
   2  * Copyright(c) 1998-2010, ALICE Experiment at CERN, All rights reserved. *\r
   3  *                                                                        *\r
   4 * Author: The ALICE Off-line Project.                                    *\r
   5  * Contributors are mentioned in the code where appropriate.              *\r
   6  *                                                                        *\r
   7  * Permission to use, copy, modify and distribute this software and its   *\r
   8  * documentation strictly for non-commercial purposes hereby granted      *\r
   9  * without fee, provided that the above copyright notice appears in all   *\r
  10  * copies and that both the copyright notice and this permission notice   *\r
  11  * appear in the supporting documentation. The authors make no claims     *\r
  12  * about the suitability of this software for any purpose. It is          *\r
  13  * provided "as is" without express or implied warranty.                  *\r
  14  **************************************************************************/\r
  15 /* $Id: $ */\r
  16 \r
  17 //_________________________________________________________________________\r
  18 //\r
  19 // Class for the electron identification and B-tagging.\r
  20 // Clusters from EMCAL matched to tracks\r
  21 // and kept in the AOD. Few histograms produced.\r
  22 //\r
  23 // -- Author: T.R.P.Aronsson (Yale) J.L. Klay (Cal Poly), M. Heinz (Yale)\r
  24 //////////////////////////////////////////////////////////////////////////////\r
  25   \r
  26 // --- ROOT system --- \r
  27 #include <TH2F.h>\r
  28 #include <TH3F.h>\r
  29 #include <TParticle.h>\r
  30 #include <TClonesArray.h>\r
  31 //#include <TObjString.h>\r
  32 //#include <Riostream.h>\r
  33 \r
  34 // --- Analysis system --- \r
  35 #include "AliAnaBtag.h" \r
  36 #include "AliCaloTrackReader.h"\r
  37 #include "AliMCAnalysisUtils.h"\r
  38 #include "AliAODCaloCluster.h"\r
  39 #include "AliFiducialCut.h"\r
  40 #include "AliAODTrack.h"\r
  41 #include "AliAODPid.h"\r
  42 #include "AliCaloPID.h"\r
  43 #include "AliAODMCParticle.h"\r
  44 #include "AliStack.h"\r
  45 #include "AliExternalTrackParam.h"\r
  46 #include "AliESDv0.h"\r
  47 #include "AliESDtrack.h"\r
  48 #include "AliAODJet.h"\r
  49 #include "AliAODEvent.h"\r
  50 #include "AliGenPythiaEventHeader.h"\r
  51 //#include "iostream.h"\r
  52 ClassImp(AliAnaBtag)\r
  53   \r
  54 //____________________________________________________________________________\r
  55 AliAnaBtag::AliAnaBtag() \r
  56 : AliAnaPartCorrBaseClass(),fCalorimeter(""),\r
  57   fpOverEmin(0.),fpOverEmax(0.),fResidualCut(0.),fMinClusEne(0.),\r
  58   fDrCut(0.),fPairDcaCut(0.),fDecayLenCut(0.),fImpactCut(0.),\r
  59   fAssocPtCut(0.),fMassCut(0.),fSdcaCut(0.),fITSCut(0),\r
  60   fNTagTrkCut(0),fIPSigCut(0.),fJetEtaCut(0.3),fJetPhiMin(1.8),fJetPhiMax(2.9),\r
  61   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
  62 {\r
  63   //default ctor\r
  64   \r
  65   //Initialize parameters\r
  66   InitParameters();\r
  67 \r
  68 }\r
  69 \r
  70 //____________________________________________________________________________\r
  71 AliAnaBtag::~AliAnaBtag() \r
  72 {\r
  73   //dtor\r
  74 \r
  75 }\r
  76 \r
  77 \r
  78 //________________________________________________________________________\r
  79 TList *  AliAnaBtag::GetCreateOutputObjects()\r
  80 {  \r
  81   // Create histograms to be saved in output file and \r
  82   // store them in outputContainer\r
  83   TList * outputContainer = new TList() ; \r
  84   outputContainer->SetName("ElectronHistos") ; \r
  85 \r
  86 \r
  87         fhEmcalElectrons = new TH1F("fhEmcalElectrons","",400,0,400);\r
  88     outputContainer->Add(fhEmcalElectrons);\r
  89     \r
  90     fhTRDElectrons = new TH1F("fhTRDElectrons","",400,0,400);\r
  91     outputContainer->Add(fhTRDElectrons);\r
  92 \r
  93     fhTPCElectrons = new TH1F("fhTPCElectrons","",400,0,400);\r
  94     outputContainer->Add(fhTPCElectrons);\r
  95 \r
  96     fhDVM1 = new TH1F("fhDVM1","",400,0,400);\r
  97     outputContainer->Add(fhDVM1);\r
  98 \r
  99     fhDVM2 = new TH1F("fhDVM2","",400,0,400);\r
 100     outputContainer->Add(fhDVM2);\r
 101 \r
 102 \r
 103     fhJets = new TH2F("fhJets","",400,0,400,20,0,20);\r
 104     outputContainer->Add(fhJets);\r
 105 \r
 106     fhJetsAllEtaPhi = new TH2F("fhJetsAllEtaPhi","",100,-2,2,100,-2,8);\r
 107     outputContainer->Add(fhJetsAllEtaPhi);\r
 108 \r
 109     fhJetsLeadingBElectronEtaPhi = new TH2F("fhJetsLeadingBElectronEtaPhi","",100,-5,5,200,-10,10);\r
 110     outputContainer->Add(fhJetsLeadingBElectronEtaPhi);\r
 111 \r
 112     fhDVM1EtaPhi = new TH2F("fhDVM1EtaPhi","",100,-2,2,100,-2,8);\r
 113     outputContainer->Add(fhDVM1EtaPhi);\r
 114 \r
 115     fhBJetElectronDetaDphi = new TH2F("fhBJetElectronDetaDphi","",100,-5,5,200,-10,10);\r
 116     outputContainer->Add(fhBJetElectronDetaDphi);\r
 117 \r
 118     fhClusterEnergy = new TH1F("fhClusterEnergy","",100,0,10);\r
 119     outputContainer->Add(fhClusterEnergy);\r
 120 \r
 121     fhTestalle = new TH1F("fhTestalle","",400,0,400);\r
 122     outputContainer->Add(fhTestalle);\r
 123 \r
 124     fhResidual = new TH1F("fhResidual","",500,0,5);\r
 125     outputContainer->Add(fhResidual);\r
 126 \r
 127     fhElectrons = new TH2F("fhElectrons","",200,0,100,20,0,20);\r
 128     outputContainer->Add(fhElectrons);\r
 129 \r
 130     fhTracks = new TH2F("fhTracks","",200,0,100,20,0,20);\r
 131     outputContainer->Add(fhTracks);\r
 132 \r
 133 \r
 134   \r
 135   return outputContainer ;\r
 136   \r
 137 }\r
 138 \r
 139 //____________________________________________________________________________\r
 140 void AliAnaBtag::Init()\r
 141 {\r
 142 \r
 143   //do some initialization\r
 144   if(fCalorimeter == "PHOS") {\r
 145     printf("AliAnaBtag::Init() - !!STOP: You want to use PHOS in analysis but this is not (yet) supported!!\n!!Check the configuration file!!\n");\r
 146     fCalorimeter = "EMCAL";\r
 147   }\r
 148   if(fCalorimeter == "EMCAL" && !GetReader()->IsEMCALSwitchedOn()){\r
 149     printf("AliAnaBtag::Init() - !!STOP: You want to use EMCAL in analysis but it is not read!!\n!!Check the configuration file!!\n");\r
 150     abort();\r
 151   }\r
 152 \r
 153 }\r
 154 \r
 155 \r
 156 //____________________________________________________________________________\r
 157 void AliAnaBtag::InitParameters()\r
 158 {\r
 159   \r
 160   //Initialize the parameters of the analysis.\r
 161   SetOutputAODClassName("AliAODPWG4Particle");\r
 162   SetOutputAODName("PWG4Particle");\r
 163 \r
 164   AddToHistogramsName("Btag_");\r
 165 \r
 166   fCalorimeter = "EMCAL" ;\r
 167   fpOverEmin = 0.5;\r
 168   fpOverEmax = 1.2;\r
 169   fResidualCut = 0.02;\r
 170   fMinClusEne = 4.0;\r
 171   //DVM B-tagging\r
 172   fDrCut       = 1.0; \r
 173   fPairDcaCut  = 0.02;\r
 174   fDecayLenCut = 1.0;\r
 175   fImpactCut   = 0.5;\r
 176   fAssocPtCut  = 1.0;\r
 177   fMassCut     = 1.5;\r
 178   fSdcaCut     = 0.1;\r
 179   fITSCut      = 4;\r
 180   //IPSig B-tagging\r
 181   fNTagTrkCut  = 2;\r
 182   fIPSigCut    = 3.0;\r
 183   //Jet fiducial cuts\r
 184   fJetEtaCut = 0.3;\r
 185   fJetPhiMin = 1.8;\r
 186   fJetPhiMax = 2.9;\r
 187 }\r
 188 \r
 189 //__________________________________________________________________\r
 190 void  AliAnaBtag::MakeAnalysisFillAOD() \r
 191 {\r
 192   //This reads in tracks, extrapolates to EMCAL, does p/E selectrons, identifies electron candidates\r
 193   //After candidates are obtained, btagging and saving into AOD.\r
 194 \r
 195   Double_t bfield = 0.;\r
 196   if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) bfield = GetReader()->GetBField();\r
 197 \r
 198   //Select the calorimeter of the electron\r
 199   if(fCalorimeter != "EMCAL") {\r
 200     printf("This class not yet implemented for PHOS\n");\r
 201     abort();\r
 202   }\r
 203   \r
 204   TObjArray *cl = GetAODEMCAL();\r
 205   \r
 206   if(!GetAODCTS() || GetAODCTS()->GetEntriesFast() == 0) return ;\r
 207   Int_t ntracks = GetAODCTS()->GetEntriesFast();\r
 208   if(GetDebug() > 0)\r
 209     printf("AliAnaBtag::MakeAnalysisFillAOD() - In CTS aod entries %d\n", ntracks);\r
 210 \r
 211   Int_t iCluster = -999;\r
 212   Int_t ntot = cl->GetEntriesFast();\r
 213 \r
 214   //CLUSTER STUFF\r
 215   if(1){\r
 216     for(Int_t iclus = 0; iclus < ntot; iclus++) {\r
 217       AliAODCaloCluster * clus = (AliAODCaloCluster*) (cl->At(iclus));\r
 218       if(!clus) continue;\r
 219       fhClusterEnergy->Fill(clus->E());\r
 220     }\r
 221   }\r
 222 \r
 223 \r
 224   for (Int_t itrk =  0; itrk <  ntracks; itrk++) {////////////// track loop\r
 225     iCluster = -999; //start with no match\r
 226     AliAODTrack * track = (AliAODTrack*) (GetAODCTS()->At(itrk)) ;\r
 227 \r
 228     Double_t imp[2] = {-999.,-999.}; Double_t cov[3] = {-999.,-999.,-999.};\r
 229     Bool_t dcaOkay = GetDCA(track,imp,cov);  //homegrown dca calculation until AOD is fixed\r
 230     if(!dcaOkay) printf("AliAnaBtag::Problem computing DCA to primary vertex for track %d.  Skipping it...\n",itrk);\r
 231 \r
 232     fhTracks->Fill(track->Pt(),1);\r
 233 \r
 234     AliAODPid* pid = (AliAODPid*) track->GetDetPid();\r
 235     if(pid == 0) {\r
 236       if(GetDebug() > 0) printf("AliAnaBtag::MakeAnalysisFillAOD() - No PID object - skipping track %d",itrk);\r
 237       continue;\r
 238     } \r
 239     else {\r
 240       Double_t emcpos[3];\r
 241       pid->GetEMCALPosition(emcpos);\r
 242       Double_t emcmom[3];\r
 243       pid->GetEMCALMomentum(emcmom);\r
 244       \r
 245       TVector3 pos(emcpos[0],emcpos[1],emcpos[2]);\r
 246       TVector3 mom(emcmom[0],emcmom[1],emcmom[2]);\r
 247       Double_t tphi = pos.Phi();\r
 248       Double_t teta = pos.Eta();\r
 249       Double_t tmom = mom.Mag();\r
 250 \r
 251       Bool_t in = kFALSE;\r
 252       if(mom.Phi()*180./TMath::Pi() > 80. && mom.Phi()*180./TMath::Pi() < 190. &&\r
 253          mom.Eta() > -0.7 && mom.Eta() < 0.7) in = kTRUE; //Checks momentum?\r
 254       //Also check the track\r
 255       if(track->Phi()*180./TMath::Pi() > 80. && track->Phi()*180./TMath::Pi() < 190. &&\r
 256          track->Eta() > -0.7 && track->Eta() < 0.7) in = kTRUE;\r
 257 \r
 258       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
 259 \r
 260       if(mom.Pt() > GetMinPt() && in) {\r
 261         fhTracks->Fill(track->Pt(),3);\r
 262         Double_t dEdx = pid->GetTPCsignal();\r
 263                 \r
 264         //Int_t ntot = cl->GetEntriesFast();\r
 265         Double_t res = 999.;\r
 266         Double_t pOverE = -999.;\r
 267         \r
 268         Int_t pidProb = track->GetMostProbablePID();\r
 269         Bool_t tpcEle = kFALSE; if(dEdx > 70.) tpcEle = kTRUE;\r
 270         Bool_t trkEle = kFALSE; if(pidProb == AliAODTrack::kElectron) trkEle = kTRUE;\r
 271         Bool_t trkChgHad = kFALSE; if(pidProb == AliAODTrack::kPion || pidProb == AliAODTrack::kKaon || pidProb == AliAODTrack::kProton) trkChgHad = kTRUE;\r
 272 \r
 273 \r
 274         //Track Matching!\r
 275         Bool_t emcEle = kFALSE;      \r
 276         double minRes=100.;\r
 277         Double_t minR  = 99;\r
 278         Double_t minPe =-1;\r
 279         Double_t minEp =-1;\r
 280         Double_t minMult = -1;\r
 281         Double_t minPt   = -1;\r
 282 \r
 283         for(Int_t iclus = 0; iclus < ntot; iclus++) {\r
 284           AliAODCaloCluster * clus = (AliAODCaloCluster*) (cl->At(iclus));\r
 285           if(!clus) continue;\r
 286 \r
 287           //As of 11-Oct-2009\r
 288           //only select "good" clusters   \r
 289 //        if (clus->GetNCells()       < 2    ) continue;\r
 290 //           if (clus->GetNCells()       > 30   ) continue;\r
 291 //           if (clus->E()               < fMinClusEne ) continue;\r
 292 //           if (clus->GetDispersion()   > 1    ) continue;\r
 293 //           if (clus->GetM20()          > 0.4  ) continue;\r
 294 //           if (clus->GetM02()          > 0.4  ) continue;\r
 295 //           if (clus->GetM20()          < 0.03 ) continue;\r
 296 //           if (clus->GetM02()          < 0.03 ) continue;\r
 297 // new optimized from ben. 2010May\r
 298           if (clus->GetNCells()       < 2    ) continue;\r
 299           if (clus->GetNCells()       > 35   ) continue;\r
 300           if (clus->E()               < 0 ) continue;\r
 301           if (clus->GetDispersion()   > 1.08    ) continue;\r
 302           if (clus->GetM20()          > 0.42  ) continue;\r
 303           if (clus->GetM02()          > 0.4  ) continue;\r
 304           if (clus->GetM20()          < 0 ) continue;\r
 305           if (clus->GetM02()          < 0.06 ) continue;\r
 306           \r
 307           \r
 308           Float_t x[3];\r
 309           clus->GetPosition(x);\r
 310           TVector3 cluspos(x[0],x[1],x[2]);\r
 311           Double_t deta = teta - cluspos.Eta();\r
 312           Double_t dphi = tphi - cluspos.Phi();\r
 313           if(dphi > TMath::Pi()) dphi -= 2*TMath::Pi();\r
 314           if(dphi < -TMath::Pi()) dphi += 2*TMath::Pi();\r
 315 \r
 316           res = sqrt(dphi*dphi + deta*deta);\r
 317           \r
 318           if(res<minRes){\r
 319             minRes=res;\r
 320             \r
 321           }\r
 322           \r
 323           if(res < 0.0275) { //   if(res < fResidualCut) { //Optimized from Ben\r
 324             iCluster = iclus;\r
 325             \r
 326             Double_t energy = clus->E(); \r
 327             if(energy > 0) pOverE = tmom/energy;\r
 328             \r
 329             if (res< minR) {\r
 330               minR  = res;\r
 331               minPe = pOverE;\r
 332               minEp = energy/tmom;\r
 333               minMult = clus->GetNCells() ;\r
 334               minPt = track->Pt();\r
 335             }\r
 336 \r
 337 \r
 338           } else {\r
 339               //unmatched\r
 340           }//res cut\r
 341         }//calo cluster loop\r
 342 \r
 343         fhResidual->Fill(minRes);\r
 344 \r
 345         if(minPe > 0.9 && minPe < 1.08) emcEle = kTRUE;//       if(minPe > fpOverEmin && minPe < fpOverEmax) emcEle = kTRUE;\r
 346         \r
 347 \r
 348         \r
 349         if(emcEle)\r
 350           fhEmcalElectrons->Fill(track->Pt());\r
 351         if(trkEle)\r
 352           fhTRDElectrons->Fill(track->Pt());\r
 353         if(tpcEle)\r
 354           fhTPCElectrons->Fill(track->Pt());\r
 355         \r
 356         \r
 357         //Take all emcal electrons, but the others only if pT < 10 GeV\r
 358         if(emcEle && (tpcEle || trkEle) ) {\r
 359           fhTestalle->Fill(track->Pt());\r
 360           //B-tagging\r
 361           if(GetDebug() > 1) printf("Found Electron - do b-tagging\n");\r
 362           Int_t dvmbtag = GetDVMBtag(track);\r
 363           \r
 364           \r
 365           if(dvmbtag>0){\r
 366             fhDVM1->Fill(track->Pt());\r
 367             fhDVM1EtaPhi->Fill(track->Eta(),track->Phi());        \r
 368             \r
 369           }\r
 370           if(dvmbtag>1)\r
 371             fhDVM2->Fill(track->Pt());\r
 372           \r
 373 \r
 374           //Create particle to save in AOD\r
 375           Double_t eMass = 0.511/1000; //mass in GeV\r
 376           Double_t eleE = sqrt(track->P()*track->P() + eMass*eMass);\r
 377           AliAODPWG4Particle tr = AliAODPWG4Particle(track->Px(),track->Py(),track->Pz(),eleE);\r
 378           tr.SetLabel(track->GetLabel());\r
 379           tr.SetCaloLabel(iCluster,-1); //sets the indices of the original caloclusters\r
 380           tr.SetTrackLabel(track->GetID(),-1); //sets the indices of the original tracks tr.SetTrackLabel(track->GetID(),-1) instead of itrk;\r
 381 \r
 382           if(emcEle) {//PID determined by EMCAL\r
 383             tr.SetDetector(fCalorimeter);\r
 384           } else {\r
 385             tr.SetDetector("CTS"); //PID determined by CTS\r
 386           }\r
 387 \r
 388           if(GetReader()->GetAODCTSNormalInputEntries() <= itrk) tr.SetInputFileIndex(1);\r
 389           //Make this preserve sign of particle\r
 390           if(track->Charge() < 0) tr.SetPdg(11); //electron is 11\r
 391           else  tr.SetPdg(-11); //positron is -11\r
 392 \r
 393           tr.SetBtag(dvmbtag);\r
 394           \r
 395 \r
 396           //Check origin of the candiates\r
 397           if(IsDataMC()){\r
 398             tr.SetTag(GetMCAnalysisUtils()->CheckOrigin(tr.GetLabel(),GetReader(),tr.GetInputFileIndex()));\r
 399             if(GetDebug() > 0) printf("AliAnaBtag::MakeAnalysisFillAOD() - Origin of candidate %d\n",tr.GetTag());\r
 400           }//Work with stack also   \r
 401           \r
 402           AddAODParticle(tr);\r
 403           \r
 404           if(GetDebug() > 1) printf("AliAnaBtag::MakeAnalysisFillAOD() - Electron selection cuts passed: pT %3.2f, pdg %d\n",tr.Pt(),tr.GetPdg());      \r
 405         }//electron\r
 406       }//pt, fiducial selection\r
 407     }//pid check\r
 408   }//track loop                         \r
 409   \r
 410 \r
 411 \r
 412 \r
 413   if(GetDebug() > 1) printf("AliAnaBtag::MakeAnalysisFillAOD()  End fill AODs \n");  \r
 414   \r
 415 }\r
 416 \r
 417 //__________________________________________________________________\r
 418 void  AliAnaBtag::MakeAnalysisFillHistograms() \r
 419 {\r
 420   //Do analysis and fill histograms\r
 421 \r
 422   AliStack * stack = 0x0;\r
 423 //   TParticle * primary = 0x0;\r
 424 \r
 425 \r
 426 \r
 427 \r
 428   if(IsDataMC()) {\r
 429     if(GetReader()->ReadStack()){\r
 430       stack =  GetMCStack() ;      \r
 431       if(!stack)\r
 432         printf("AliAnaBtag::MakeAnalysisFillHistograms() *** no stack ***: \n");\r
 433       \r
 434     }\r
 435   }// is data and MC\r
 436 \r
 437   ////////////////////////////////////\r
 438   //Loop over jets and check for b-tag\r
 439   ////////////////////////////////////\r
 440   double maxjetEta=-4.;\r
 441   double maxjetPhi=-4.;\r
 442 \r
 443   Int_t njets = (GetReader()->GetOutputEvent())->GetNJets();\r
 444   if(njets > 0) {\r
 445     if(GetDebug() > 0) printf("AliAnaBtag::MakeAnalysisFillHistograms() - Jet AOD branch has %d jets.  Performing b-jet tag analysis\n",njets);\r
 446 \r
 447     for(Int_t ijet = 0; ijet < njets ; ijet++) {\r
 448       AliAODJet * jet = (AliAODJet*)(GetReader()->GetOutputEvent())->GetJet(ijet) ;\r
 449 \r
 450       if(ijet==0){\r
 451         maxjetEta=jet->Eta();\r
 452         maxjetPhi=jet->Phi();\r
 453       }\r
 454 \r
 455       fhJets->Fill(jet->Pt(),1);\r
 456       fhJetsAllEtaPhi->Fill(jet->Eta(),jet->Phi());\r
 457 \r
 458       if(jet->Pt() < 0.) continue; //This has to be adjusted depending on pp or AA!\r
 459       fhJets->Fill(jet->Pt(),3); //All jets after pt cut\r
 460 \r
 461       //Geometric EMCAL cut\r
 462       if(TMath::Abs(jet->Eta()) > fJetEtaCut) continue;\r
 463       if(jet->Phi() < fJetPhiMin || jet->Phi() > fJetPhiMax) continue;\r
 464       fhJets->Fill(jet->Pt(),4); //All jets after geometric cut\r
 465 \r
 466       Bool_t leadJet  = kFALSE;\r
 467       if (ijet==0){ \r
 468         fhJets->Fill(jet->Pt(),5); //Leading jets\r
 469         leadJet= kTRUE;\r
 470       }\r
 471 \r
 472 \r
 473       Bool_t dvmJet = kFALSE;  \r
 474       TRefArray* rt = jet->GetRefTracks();\r
 475       Int_t ntrk = rt->GetEntries();\r
 476 \r
 477       for(Int_t itrk = 0; itrk < ntrk; itrk++) {\r
 478         AliAODTrack* jetTrack = (AliAODTrack*)jet->GetTrack(itrk);\r
 479         Bool_t isDVM = CheckIfBjet(jetTrack);\r
 480         if(isDVM) dvmJet = kTRUE;\r
 481       }\r
 482 \r
 483       if(dvmJet)\r
 484         fhJets->Fill(jet->Pt(),6);\r
 485 \r
 486 \r
 487 \r
 488       if(IsDataMC()) {\r
 489         //determine tagging efficiency & mis-tagging rate\r
 490         //using b-quarks from stack\r
 491         Bool_t isTrueBjet = IsMcBJet(jet->Eta(), jet->Phi());\r
 492         Bool_t isTrueDjet = IsMcDJet(jet->Eta(), jet->Phi());\r
 493         if (isTrueBjet && GetDebug() > 0) printf("== True Bjet==\n");\r
 494         if (isTrueDjet && GetDebug() > 0) printf("== True Charm-jet==\n");\r
 495         if (dvmJet && GetDebug() > 0)     printf("== found DVM jet==\n");\r
 496       }\r
 497 \r
 498     } //jet loop\r
 499   } //jets exist\r
 500   \r
 501 \r
 502 \r
 503 \r
 504   //Electron loop, read back electrons, fill histos\r
 505   Int_t naod = GetOutputAODBranch()->GetEntriesFast();\r
 506   if(GetDebug() > 0) printf("AliAnaBtag::MakeAnalysisFillHistograms() - aod branch entries %d\n", naod);\r
 507   \r
 508   for(Int_t iaod = 0; iaod < naod ; iaod++){\r
 509     AliAODPWG4Particle* ele =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));\r
 510     Int_t pdg = ele->GetPdg();\r
 511 \r
 512 \r
 513     if(TMath::Abs(pdg) != AliCaloPID::kElectron) continue; //not necessary..\r
 514 \r
 515     //MC tag of this electron\r
 516     //    Int_t mctag = ele->GetTag();\r
 517 \r
 518 \r
 519     fhElectrons->Fill(ele->Pt(),1); //All electrons\r
 520     Bool_t photonic = kFALSE;\r
 521     photonic = PhotonicV0(ele->GetTrackLabel(0)); //check against V0s\r
 522     if(!photonic) fhElectrons->Fill(ele->Pt(),3); //nonphotonic electrons\r
 523     if(photonic) fhElectrons->Fill(ele->Pt(),4);  //photonic electrons\r
 524 \r
 525     //Fill electron histograms \r
 526     Float_t phiele = ele->Phi();\r
 527     Float_t etaele = ele->Eta();\r
 528 \r
 529 \r
 530     if(ele->GetBtag()>0){ // removed bit tag shit\r
 531       fhElectrons->Fill(ele->Pt(),5);\r
 532       if(!photonic) fhElectrons->Fill(ele->Pt(),6);\r
 533       if(photonic) fhElectrons->Fill(ele->Pt(),7);\r
 534       fhJetsLeadingBElectronEtaPhi->Fill(maxjetEta,maxjetPhi); \r
 535       double deta=etaele-maxjetEta;\r
 536       double dphi=phiele-maxjetPhi;\r
 537       \r
 538       //double r = sqrt((dphi*dphi)+(deta*deta));\r
 539       fhBJetElectronDetaDphi->Fill(deta,dphi);\r
 540       \r
 541     }\r
 542     \r
 543   }//electron aod loop\r
 544 \r
 545 }\r
 546 \r
 547 //__________________________________________________________________\r
 548 Int_t AliAnaBtag::GetDVMBtag(AliAODTrack * tr )\r
 549 {\r
 550   //This method uses the Displaced Vertex between electron-hadron\r
 551   //pairs and the primary vertex to determine whether an electron is\r
 552   //likely from a B hadron.\r
 553 \r
 554   Int_t ncls1 = 0;\r
 555   for(Int_t l = 0; l < 6; l++) if(TESTBIT(tr->GetITSClusterMap(),l)) ncls1++;\r
 556 \r
 557 \r
 558   if (ncls1 < fITSCut) return 0;\r
 559 \r
 560   Double_t imp[2] = {-999.,-999.}; Double_t cov[3] = {-999.,-999.,-999.};\r
 561   Bool_t dcaOkay = GetDCA(tr,imp,cov);  //homegrown dca calculation until AOD is fixed                  \r
 562   if(!dcaOkay) {\r
 563     printf("AliAnaBtag::Problem computing DCA to primary vertex for track %d",tr->GetID());\r
 564     return 0;\r
 565   }\r
 566 \r
 567   if (TMath::Abs(imp[0])   > fImpactCut ) return 0;\r
 568   if (TMath::Abs(imp[1])   > fImpactCut ) return 0;\r
 569 \r
 570   Int_t nvtx1 = 0;\r
 571   Int_t nvtx2 = 0;\r
 572   Int_t nvtx3 = 0;\r
 573 \r
 574   for (Int_t k2 =0; k2 < GetAODCTS()->GetEntriesFast() ; k2++) {\r
 575     //loop over assoc\r
 576     AliAODTrack* track2 = (AliAODTrack*)GetAODCTS()->At(k2);\r
 577     Int_t id1 = tr->GetID();\r
 578     Int_t id2 = track2->GetID();\r
 579     if(id1 == id2) continue;\r
 580 \r
 581     Int_t ncls2 = 0;\r
 582     for(Int_t l = 0; l < 6; l++) if(TESTBIT(track2->GetITSClusterMap(),l)) ncls2++;\r
 583     if (ncls2 < fITSCut) continue;\r
 584 \r
 585     if(tr->Pt()<6.&&track2->Pt() < 0.4) continue;\r
 586     if(tr->Pt()>6.&&tr->Pt()<10.&&track2->Pt() < 0.2) continue;\r
 587     if(tr->Pt()>10.&&track2->Pt() < 0.6) continue;\r
 588 \r
 589 \r
 590     Double_t dphi = tr->Phi() - track2->Phi();\r
 591     if(dphi > TMath::Pi()) dphi -= 2*TMath::Pi();\r
 592     if(dphi < -TMath::Pi()) dphi += 2*TMath::Pi();\r
 593     Double_t deta = tr->Eta() - track2->Eta();\r
 594     Double_t dr = sqrt(deta*deta + dphi*dphi);\r
 595 \r
 596     if(dr > fDrCut) continue;\r
 597     \r
 598     if(tr->Pt()<6.){\r
 599       Double_t sDca = ComputeSignDca(tr, track2, 1.4,0.025);\r
 600         if(sDca > 0.06) nvtx2++;\r
 601     } \r
 602     if(tr->Pt()>6.&&tr->Pt()<10.){\r
 603       Double_t sDca = ComputeSignDca(tr, track2, 1.7,0.012);\r
 604         if(sDca > 0.03) nvtx2++;\r
 605     } \r
 606     if(tr->Pt()>10.){\r
 607       Double_t sDca = ComputeSignDca(tr, track2, 1.5,0.14);\r
 608         if(sDca > 0.04) nvtx2++;\r
 609     } \r
 610 \r
 611 \r
 612   } //loop over hadrons\r
 613 \r
 614   if(GetDebug() > 0) {\r
 615     if (nvtx1>0) printf("result1 of btagging: %d \n",nvtx1);\r
 616     if (nvtx2>0) printf("result2 of btagging: %d \n",nvtx2);\r
 617     if (nvtx3>0) printf("result3 of btagging: %d \n",nvtx3);\r
 618   }\r
 619 \r
 620 \r
 621 \r
 622 \r
 623   return nvtx2;\r
 624 \r
 625 }\r
 626 \r
 627 //__________________________________________________________________\r
 628 Double_t AliAnaBtag::ComputeSignDca(AliAODTrack *tr, AliAODTrack *tr2 , float masscut, double pdcacut)\r
 629 {\r
 630   //Compute the signed dca between two tracks\r
 631   //and return the result\r
 632 \r
 633   Double_t signDca=-999.;\r
 634   if(GetDebug() > 2 ) printf(">>ComputeSdca:: track1 %d, track2 %d, masscut %f \n", tr->GetLabel(), tr2->GetLabel(), masscut);\r
 635 \r
 636   //=====Now calculate DCA between both tracks=======  \r
 637   Double_t massE = 0.000511;\r
 638   Double_t massK = 0.493677;\r
 639 \r
 640   Double_t vertex[3] = {-999.,-999.,-999}; //vertex\r
 641   if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) {\r
 642     GetReader()->GetVertex(vertex); //If only one file, get the vertex from there\r
 643     //FIXME:  Add a check for whether file 2 is PYTHIA or HIJING\r
 644     //If PYTHIA, then set the vertex from file 2, if not, use the\r
 645     //vertex from file 1\r
 646     if(GetReader()->GetSecondInputAODTree()) GetReader()->GetSecondInputAODVertex(vertex);\r
 647   }\r
 648   \r
 649   TVector3 primV(vertex[0],vertex[1],vertex[2]) ;\r
 650 \r
 651   if(GetDebug() > 5) printf(">>ComputeSdca:: primary vertex = %2.2f,%2.2f,%2.2f \n",vertex[0],vertex[1],vertex[2]) ;\r
 652 \r
 653   AliExternalTrackParam *param1 = new AliExternalTrackParam(tr);\r
 654   AliExternalTrackParam *param2 = new AliExternalTrackParam(tr2);\r
 655 \r
 656   Double_t bfield[3];\r
 657   param1->GetBxByBz(bfield);\r
 658   bfield[0]=0;\r
 659   bfield[1]=0;\r
 660   bfield[2]=5.;\r
 661   Double_t bz = 5.; // = param1->GetBz();\r
 662   //cout<<"In ComputeSignDCA, bfield[3] and bz: "<<bfield[0]<<" "<<bfield[1]<<" "<<bfield[2]<<" "<<bz<<endl;\r
 663   Double_t xplane1 = 0.; Double_t xplane2 = 0.;\r
 664   Double_t pairdca = param1->GetDCA(param2,bz,xplane1,xplane2);\r
 665 \r
 666   param1->PropagateToBxByBz(xplane1,bfield);\r
 667   param2->PropagateToBxByBz(xplane2,bfield);\r
 668 \r
 669   Int_t id1 = 0, id2 = 0;\r
 670   AliESDv0 bvertex(*param1,id1,*param2,id2);\r
 671   Double_t vx,vy,vz;\r
 672   bvertex.GetXYZ(vx,vy,vz);\r
 673 \r
 674   Double_t emom[3];\r
 675   Double_t hmom[3];\r
 676   param1->PxPyPz(emom);\r
 677   param2->PxPyPz(hmom);\r
 678   TVector3 emomAtB(emom[0],emom[1],emom[2]);\r
 679   TVector3 hmomAtB(hmom[0],hmom[1],hmom[2]);\r
 680   TVector3 secvtxpt(vx,vy,vz);\r
 681   TVector3 decayvector(0,0,0);\r
 682   decayvector = secvtxpt - primV; //decay vector from PrimVtx\r
 683   Double_t decaylength = decayvector.Mag();\r
 684 \r
 685   //printf("\t JLK pairDCA = %2.2f\n",pairdca);\r
 686 \r
 687   if(GetDebug() > 0) {\r
 688     printf(">>ComputeSdca:: mom1=%f, mom2=%f \n", emomAtB.Perp(), hmomAtB.Perp() );\r
 689     printf(">>ComputeSdca:: pairDCA=%f, length=%f \n", pairdca,decaylength );\r
 690   }\r
 691 \r
 692 \r
 693 \r
 694   if (emomAtB.Mag()>0 && pairdca < pdcacut && decaylength < fDecayLenCut ) {\r
 695     TVector3 sumMom = emomAtB+hmomAtB;\r
 696     Double_t ener1 = sqrt(pow(emomAtB.Mag(),2) + massE*massE);\r
 697     Double_t ener2 = sqrt(pow(hmomAtB.Mag(),2) + massK*massK);\r
 698     Double_t ener3 = sqrt(pow(hmomAtB.Mag(),2) + massE*massE);\r
 699     Double_t mass = sqrt(pow((ener1+ener2),2) - pow(sumMom.Mag(),2));\r
 700     Double_t massPhot = sqrt(pow((ener1+ener3),2) - pow(sumMom.Mag(),2));\r
 701     Double_t sDca = decayvector.Dot(emomAtB)/emomAtB.Mag();\r
 702 \r
 703 //     if (masscut<1.1) fhDVMBtagQA2->Fill(sDca, mass);\r
 704 \r
 705     if (mass > masscut && massPhot > 0.1) signDca = sDca;\r
 706     \r
 707     if(GetDebug() > 0) printf("\t>>ComputeSdca:: mass=%f \n", mass);\r
 708     if(GetDebug() > 0) printf("\t>>ComputeSdca:: sec vtx-signdca :%f\n",signDca);\r
 709   }\r
 710 \r
 711   //clean up\r
 712   delete param1;\r
 713   delete param2;\r
 714 \r
 715   return signDca;\r
 716 }\r
 717 \r
 718 \r
 719 //__________________________________________________________________\r
 720 Bool_t AliAnaBtag::PhotonicV0(Int_t id) \r
 721 {\r
 722   //This method checks to see whether a track that has been flagged as\r
 723   //an electron was determined to match to a V0 candidate with\r
 724   //invariant mass consistent with photon conversion\r
 725 \r
 726   Bool_t itIS = kFALSE;\r
 727 \r
 728   Double_t massEta = 0.547;\r
 729   Double_t massRho0 = 0.770;\r
 730   Double_t massOmega = 0.782;\r
 731   Double_t massPhi = 1.020;\r
 732   \r
 733   //---Get V0s---\r
 734   AliAODEvent *aod = (AliAODEvent*) GetReader()->GetInputEvent();\r
 735   int nv0s = aod->GetNumberOfV0s();\r
 736   for (Int_t iV0 = 0; iV0 < nv0s; iV0++) {\r
 737     AliAODv0 *v0 = aod->GetV0(iV0);\r
 738     if (!v0) continue;\r
 739     double radius = v0->RadiusV0();\r
 740     double mass = v0->InvMass2Prongs(0,1,11,11);\r
 741     if(GetDebug() > 0) {\r
 742       printf("## PhotonicV0() :: v0: %d, radius: %f \n", iV0 , radius );\r
 743       printf("## PhotonicV0() :: neg-id: %d, pos-id: %d, THIS id: %d\n", v0->GetNegID(), v0->GetPosID(), id);\r
 744       printf("## PhotonicV0() :: Minv(e,e): %f \n", v0->InvMass2Prongs(0,1,11,11) );\r
 745     }\r
 746     if (mass < 0.100 ||\r
 747         (mass > massEta-0.05 || mass < massEta+0.05) ||\r
 748         (mass > massRho0-0.05 || mass < massRho0+0.05) ||\r
 749         (mass > massOmega-0.05 || mass < massOmega+0.05) ||\r
 750         (mass > massPhi-0.05 || mass < massPhi+0.05)) {\r
 751       if ( id == v0->GetNegID() || id == v0->GetPosID()) {\r
 752         itIS=kTRUE;\r
 753         if(GetDebug() > 0) printf("## PhotonicV0() :: It's a conversion electron!!! \n" );\r
 754       }\r
 755     } }\r
 756   return itIS;\r
 757 \r
 758 }\r
 759 \r
 760 //__________________________________________________________________\r
 761 Bool_t AliAnaBtag::GetDCA(const AliAODTrack* track,Double_t impPar[2], Double_t cov[3]) \r
 762 {\r
 763   //Use the Event vertex and AOD track information to get\r
 764   //a real impact parameter for the track\r
 765 \r
 766 \r
 767   Double_t maxD = 100000.; //max transverse IP\r
 768   if(GetReader()->GetDataType() != AliCaloTrackReader::kMC) {\r
 769     AliVEvent* ve = (AliVEvent*)GetReader()->GetInputEvent();\r
 770     AliVVertex *vv = (AliVVertex*)ve->GetPrimaryVertex();\r
 771     AliESDtrack esdTrack(track);\r
 772     Double_t bfield[3];\r
 773     esdTrack.GetBxByBz(bfield);\r
 774     bfield[0]=0;\r
 775     bfield[1]=0;\r
 776     bfield[2]=5.;\r
 777     \r
 778     Bool_t gotit = esdTrack.PropagateToDCABxByBz(vv,bfield,maxD,impPar,cov);\r
 779     return gotit;\r
 780   }\r
 781 \r
 782   return kFALSE;\r
 783 \r
 784 }\r
 785 //__________________________________________________________________\r
 786 Bool_t AliAnaBtag::CheckIfBjet(const AliAODTrack* track)\r
 787 {\r
 788   Bool_t bjet = kFALSE;\r
 789   Int_t trackId = track->GetID(); //get the index in the reader\r
 790   Int_t naod = GetOutputAODBranch()->GetEntriesFast();\r
 791   if(GetDebug() > 3) printf("AliAnaBtag::CheckIfBjet() - aod branch entries %d\n", naod);\r
 792 \r
 793   for(Int_t iaod = 0; iaod < naod ; iaod++){\r
 794     AliAODPWG4Particle* ele =  (AliAODPWG4Particle*) (GetOutputAODBranch()->At(iaod));\r
 795     Int_t electronlabel = ele->GetTrackLabel(0);\r
 796     if(electronlabel != trackId) continue;  //skip to the next one if they don't match\r
 797     if(ele->GetBtag()>0)\r
 798       bjet = kTRUE;\r
 799   } \r
 800 \r
 801   return bjet;\r
 802 } \r
 803 \r
 804 \r
 805 \r
 806 //__________________________________________________________________\r
 807 AliAODMCParticle* AliAnaBtag::GetMCParticle(Int_t ipart) \r
 808 {\r
 809   //Get the MC particle at position ipart\r
 810 \r
 811   AliAODMCParticle* aodprimary = 0x0;\r
 812   TClonesArray * mcparticles0 = 0x0;\r
 813   TClonesArray * mcparticles1 = 0x0;\r
 814 \r
 815   if(GetReader()->ReadAODMCParticles()){\r
 816     //Get the list of MC particles                                                                                                                           \r
 817     mcparticles0 = GetReader()->GetAODMCParticles(0);\r
 818     if(!mcparticles0 && GetDebug() > 0) {\r
 819       printf("AliAnaBtag::MakeAnalysisFillHistograms() -  Standard MCParticles not available!\n");\r
 820     }\r
 821     if(GetReader()->GetSecondInputAODTree()){\r
 822       mcparticles1 = GetReader()->GetAODMCParticles(1);\r
 823       if(!mcparticles1 && GetDebug() > 0) {\r
 824         printf("AliAnaBtag::MakeAnalysisFillHistograms() -  Second input MCParticles not available!\n");\r
 825       }\r
 826     }\r
 827 \r
 828     Int_t npart0 = mcparticles0->GetEntriesFast();\r
 829     Int_t npart1 = 0;\r
 830     if(mcparticles1) npart1 = mcparticles1->GetEntriesFast();\r
 831     if(ipart < npart0) aodprimary = (AliAODMCParticle*)mcparticles0->At(ipart);\r
 832     else aodprimary = (AliAODMCParticle*)mcparticles1->At(ipart-npart0);\r
 833     if(!aodprimary) {\r
 834       printf("AliAnaBtag::GetMCParticle() *** no primary ***:  label %d \n", ipart);\r
 835       return 0x0;\r
 836     }\r
 837 \r
 838   } else {\r
 839     printf("AliAnaBtag::GetMCParticle() - Asked for AliAODMCParticle but we have a stack reader.\n");\r
 840   }\r
 841   return aodprimary;\r
 842 \r
 843 }\r
 844 \r
 845 //__________________________________________________________________\r
 846 Bool_t  AliAnaBtag::IsMcBJet(Double_t jeta, Double_t jphi)\r
 847 {\r
 848   //Check the jet eta,phi against that of the b-quark\r
 849   //to decide whether it is an MC B-jet\r
 850   Bool_t bjet=kFALSE;\r
 851 \r
 852   AliStack* stack = 0x0;\r
 853   \r
 854   for(Int_t ipart = 0; ipart < 100; ipart++) {\r
 855 \r
 856     Double_t pphi = -999.;\r
 857     Double_t peta = -999.;\r
 858     Int_t pdg = 0;\r
 859     if(GetReader()->ReadStack()) {\r
 860       stack = GetMCStack();\r
 861       if(!stack) {\r
 862         printf("AliAnaBtag::IsMCBJet() *** no stack ***: \n");\r
 863         return kFALSE;\r
 864       }\r
 865       TParticle* primary = stack->Particle(ipart);\r
 866       if (!primary) continue;\r
 867       pdg = primary->GetPdgCode();\r
 868       pphi = primary->Phi();\r
 869       peta = primary->Eta();\r
 870     } else if(GetReader()->ReadAODMCParticles()) {\r
 871       AliAODMCParticle* aodprimary = GetMCParticle(ipart);\r
 872       if(!aodprimary) continue;\r
 873       pdg = aodprimary->GetPdgCode();\r
 874       pphi = aodprimary->Phi();\r
 875       peta = aodprimary->Eta();\r
 876     }\r
 877     if ( TMath::Abs(pdg) != 5) continue;\r
 878     \r
 879     //      printf("MTH: IsMcBJet : %d, pdg=%d : pt=%f \n", ipart, pdgcode, primary->Pt());\r
 880     Double_t dphi = jphi - pphi;\r
 881     Double_t deta = jeta - peta;\r
 882     Double_t dr = sqrt(deta*deta + dphi*dphi);\r
 883     \r
 884     if (dr < 0.2) {\r
 885       bjet=kTRUE;\r
 886       //printf("MTH: **** found matching MC-Bjet: PDG=%d, pt=%f,dr=%f \n", pdgcode, primary->Pt(),dr );\r
 887       break;\r
 888     }\r
 889   }\r
 890   return bjet;\r
 891 \r
 892 }\r
 893 \r
 894 //__________________________________________________________________\r
 895 Bool_t  AliAnaBtag::IsMcDJet(Double_t jeta, Double_t jphi)\r
 896 {\r
 897   //Check if this jet is a charm jet\r
 898   Bool_t cjet=kFALSE;\r
 899 \r
 900   AliStack* stack = 0x0;\r
 901 \r
 902   for(Int_t ipart = 0; ipart < 100; ipart++) {\r
 903     \r
 904     Double_t pphi = -999.;\r
 905     Double_t peta = -999.;\r
 906     Int_t pdg = 0;\r
 907     if(GetReader()->ReadStack()) {\r
 908       stack = GetMCStack();\r
 909       if(!stack) {\r
 910         printf("AliAnaBtag::IsMCDJet() *** no stack ***: \n");\r
 911         return kFALSE;\r
 912       }\r
 913       TParticle* primary = stack->Particle(ipart);\r
 914       if (!primary) continue;\r
 915       pdg = primary->GetPdgCode();\r
 916       pphi = primary->Phi();\r
 917       peta = primary->Eta();\r
 918     } else if(GetReader()->ReadAODMCParticles()) {\r
 919       AliAODMCParticle* aodprimary = GetMCParticle(ipart);\r
 920       if(!aodprimary) continue;\r
 921       pdg = aodprimary->GetPdgCode();\r
 922       pphi = aodprimary->Phi();\r
 923       peta = aodprimary->Eta();\r
 924     }\r
 925 \r
 926     if ( TMath::Abs(pdg) != 4) continue;\r
 927 \r
 928     Double_t dphi = jphi - pphi;\r
 929     Double_t deta = jeta - peta;\r
 930     Double_t dr = sqrt(deta*deta + dphi*dphi);\r
 931     \r
 932     if (dr < 0.2) {\r
 933       cjet=kTRUE;\r
 934       break;\r
 935     }\r
 936   }\r
 937 \r
 938   return cjet;\r
 939 \r
 940 }\r
 941 \r
 942 //__________________________________________________________________\r
 943 void AliAnaBtag::Print(const Option_t * opt) const\r
 944 {\r
 945   //Print some relevant parameters set for the analysis\r
 946   \r
 947   if(! opt)\r
 948     return;\r
 949   \r
 950   printf("**** Print %s %s ****\n", GetName(), GetTitle() ) ;\r
 951   AliAnaPartCorrBaseClass::Print(" ");\r
 952 \r
 953   printf("Calorimeter            =     %s\n", fCalorimeter.Data()) ;\r
 954   printf("pOverE range           =     %f - %f\n",fpOverEmin,fpOverEmax);\r
 955   printf("residual cut           =     %f\n",fResidualCut);\r
 956   printf("---DVM Btagging\n");\r
 957   printf("max IP-cut (e,h)       =     %f\n",fImpactCut);\r
 958   printf("min ITS-hits           =     %d\n",fITSCut);\r
 959   printf("max dR (e,h)           =     %f\n",fDrCut);\r
 960   printf("max pairDCA            =     %f\n",fPairDcaCut);\r
 961   printf("max decaylength        =     %f\n",fDecayLenCut);\r
 962   printf("min Associated Pt      =     %f\n",fAssocPtCut);\r
 963   printf("---IPSig Btagging\n");\r
 964   printf("min tag track          =     %d\n",fNTagTrkCut);\r
 965   printf("min IP significance    =     %f\n",fIPSigCut);\r
 966   printf("    \n") ;\r
 967         \r
 968 } \r
 969 \r
 970 \r
 971 //__________________________________________________________________\r
 972 void  AliAnaBtag::Terminate(TList* outputList)\r
 973 {\r
 974  \r
 975   //Do some plots to end\r
 976   //Recover histograms from output histograms list, needed for\r
 977   //distributed analysis.                \r
 978   //ReadHistograms(outputList);\r
 979 \r
 980   printf(" AliAnaBtag::Terminate()  *** %s Report: %d outputs\n", GetName(), outputList->GetEntries()) ;\r
 981 \r
 982 }\r
 983 \r