PlotEfficiency.C

Go to the documentation of this file.

/**************************************************************************
 * This file is property of and copyright by the ALICE HLT Project        *
 * All rights reserved.                                                   *
 *                                                                        *
 * Primary Authors:                                                       *
 *   Artur Szostak <artursz@iafrica.com>                                  *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/



#if !defined(__CINT__) || defined(__MAKECINT__)
#include "Riostream.h"
#include "TH1.h"
#include "TH2.h"
#include "TFile.h"
#include "TList.h"
#include "TGraphAsymmErrors.h"
#include "TCanvas.h"
#endif


void EfficiencyIntegral(const char* prefix, TH1* pass, TH1* total, double& eff, double& error, bool includeOverflow = false)
{
  
  TH1D passhist("temppasshist", "", 1, 0, 1);
  TH1D totalhist("temptotalhist", "", 1, 0, 1);
  if (includeOverflow)
  {
    passhist.SetBinContent(1, pass->Integral(0, pass->GetNbinsX()+1));
    totalhist.SetBinContent(1, total->Integral(0, total->GetNbinsX()+1));
  }
  else
  {
    passhist.SetBinContent(1, pass->Integral(1, pass->GetNbinsX()));
    totalhist.SetBinContent(1, total->Integral(1, total->GetNbinsX()));
  }
  TGraphAsymmErrors* effgraph = new TGraphAsymmErrors(&passhist, &totalhist);
  eff = effgraph->GetY()[0];
  double effhigh = effgraph->GetEYhigh()[0];
  double efflow = effgraph->GetEYlow()[0];
  error = effhigh > efflow ? effhigh : efflow;
  cout << prefix << eff << " + " << effhigh << " - " << efflow << endl;
}


void PlotEfficiency(const char* histfilename = "hists.root", bool includeOverflow = true)
{
  
  TFile file(histfilename);
  
  TList* histolist = NULL;
  TIter next(file.GetListOfKeys());
  for (TObject* key = next(); key != NULL; key = next())
  {
    TList* list = dynamic_cast<TList*>( file.Get(key->GetName()) );
    if (list != NULL)
    {
      histolist = list;
      break;
    }
  }
  
  if (histolist == NULL)
  {
    cout << "ERROR: Could not find the histogram list in file '" << histfilename << "'" << endl;
    return;
  }
  
  TH2* findableTracksHist = dynamic_cast<TH2*>( histolist->FindObject("findableTracksHist") );
  TH2* foundTracksHistMC = dynamic_cast<TH2*>( histolist->FindObject("foundTracksHistMC") );
  TH2* fakeTracksHist = dynamic_cast<TH2*>( histolist->FindObject("fakeTracksHist") );
  
  if (findableTracksHist == NULL)
  {
    cout << "ERROR: Could not find the histogram findableTracksHist in TList '" << histolist->GetName() << "'" << endl;
    return;
  }
  if (foundTracksHistMC == NULL)
  {
    cout << "ERROR: Could not find the histogram foundTracksHistMC in TList '" << histolist->GetName() << "'" << endl;
    return;
  }
  if (fakeTracksHist == NULL)
  {
    cout << "ERROR: Could not find the histogram fakeTracksHist in TList '" << histolist->GetName() << "'" << endl;
    return;
  }
  
  TH1* projeffx = foundTracksHistMC->ProjectionX();
  TH1* projfindx = findableTracksHist->ProjectionX();
  TGraphAsymmErrors* effpt = new TGraphAsymmErrors(projeffx, projfindx);
  effpt->Draw("ap");
  effpt->GetHistogram()->SetTitle("Efficiency of track reconstruction.");
  effpt->GetHistogram()->SetXTitle(projeffx->GetXaxis()->GetTitle());
  effpt->GetHistogram()->SetYTitle("#epsilon");
  effpt->GetYaxis()->SetRangeUser(0, 1.1);
  
  new TCanvas;
  TH1* projfakex = fakeTracksHist->ProjectionX();
  TGraphAsymmErrors* fakept = new TGraphAsymmErrors(projfakex, projfindx);
  fakept->Draw("ap");
  fakept->GetHistogram()->SetTitle("Fake reconstructed tracks ratio.");
  fakept->GetHistogram()->SetXTitle(projfakex->GetXaxis()->GetTitle());
  fakept->GetHistogram()->SetYTitle("fake tracks / findable tracks");
  
  new TCanvas;
  TH1* projeffy = foundTracksHistMC->ProjectionY();
  TH1* projfindy = findableTracksHist->ProjectionY();
  TGraphAsymmErrors* effy = new TGraphAsymmErrors(projeffy, projfindy);
  effy->Draw("ap");
  effy->GetHistogram()->SetTitle("Efficiency of track reconstruction.");
  effy->GetHistogram()->SetXTitle(projeffy->GetXaxis()->GetTitle());
  effy->GetHistogram()->SetYTitle("#epsilon");
  effy->GetYaxis()->SetRangeUser(0, 1.1);
  
  new TCanvas;
  TH1* projfakey = fakeTracksHist->ProjectionY();
  TGraphAsymmErrors* fakey = new TGraphAsymmErrors(projfakey, projfindy);
  fakey->Draw("ap");
  fakey->GetHistogram()->SetTitle("Fake reconstructed tracks ratio.");
  fakey->GetHistogram()->SetXTitle(projfakey->GetXaxis()->GetTitle());
  fakey->GetHistogram()->SetYTitle("fake tracks / findable tracks");
  
  if (includeOverflow)
  {
    cout << "findable tracks = " << projfindx->Integral(0, projfindx->GetNbinsX()+1) << endl;
    cout << "found tracks    = " << projeffx->Integral(0, projeffx->GetNbinsX()+1) << endl;
    cout << "fake tracks     = " << projfakex->Integral(0, projfakex->GetNbinsX()+1) << endl;
  }
  else
  {
    cout << "findable tracks = " << projfindx->Integral(0, projfindx->GetNbinsX()+1) << endl;
    cout << "found tracks    = " << projeffx->Integral(0, projeffx->GetNbinsX()+1) << endl;
    cout << "fake tracks     = " << projfakex->Integral(0, projfakex->GetNbinsX()+1) << endl;
  }
  double eff, fake, deff, dfake;
  EfficiencyIntegral("Total efficiency = ", projeffx, projfindx, eff, deff, includeOverflow);
  EfficiencyIntegral("Total fake ratio = ", projfakex, projfindx, fake, dfake, includeOverflow);
  cout << "Sum = " << eff + fake << " +/- " << sqrt(deff*deff + dfake*dfake) << endl;
}