SchieferD/jetcalib/mcfakerate_x.cpp

////////////////////////////////////////////////////////////////////////////////
//
// mcfakerate_x
// ------------
//
//            06/28/2007 Philipp Schieferdecker <philipp.schieferdecker@cern.ch>
////////////////////////////////////////////////////////////////////////////////


#include "utils/cmdline.h"

#include <TROOT.h>
#include <TStyle.h>
#include <TRint.h>
#include <TFile.h>
#include <TTree.h>
#include <TH1F.h>
#include <TCanvas.h>
#include <TLegend.h>
#include <TMath.h>

#include <iostream>
#include <iomanip>
#include <cmath>
#include <string>
#include <vector>


using namespace std;


////////////////////////////////////////////////////////////////////////////////
// declare global functions
////////////////////////////////////////////////////////////////////////////////

//______________________________________________________________________________
void   setBinomialErrors(TH1F* hEff,const TH1F* hFake, const TH1F* hAll);
string getLegLabel(const string& jetalg);

////////////////////////////////////////////////////////////////////////////////
// main
////////////////////////////////////////////////////////////////////////////////

//______________________________________________________________________________
int main(int argc,char**argv)
{
  cmdline cl;
  cl.parse(argc,argv);

  vector<string> input    = cl.get_vector<string>("input");
  string         treename = cl.get_value<string> ("treename",   "t");
  float          drmax    = cl.get_value<float>  ("drmax",      0.3);
  int            netbins  = cl.get_value<int>    ("netbins",     50);
  int            netabins = cl.get_value<int>    ("netabins",    25);
  int            nphibins = cl.get_value<int>    ("nphibins",    25);
  int            nemfbins = cl.get_value<int>    ("nemfbins",    15);
  float          etmin    = cl.get_value<float>  ("etmin",      0.0);
  float          etmax    = cl.get_value<float>  ("etmax",    200.0);
  float          etamin   = cl.get_value<float>  ("etamin",    -5.0);
  float          etamax   = cl.get_value<float>  ("etamax",     5.0);
  float          phimin   = cl.get_value<float>  ("phimin",    -3.2);
  float          phimax   = cl.get_value<float>  ("phimax",     3.2);
  float          emfmin   = cl.get_value<float>  ("emfmin",     0.0);
  float          emfmax   = cl.get_value<float>  ("emfmax",     1.0);
  short          applyjes = cl.get_value<short>  ("applyjes",     0);
  
  if (!cl.check()) return 0;
  cl.print();

  
  // jetalgs
  vector<string> jetalgs;
  for (vector<string>::const_iterator it=input.begin();it!=input.end();++it) {
    string jetalg = *it;
    string::size_type pos;
    while ((pos=jetalg.find('/'))!=string::npos) jetalg = jetalg.substr(pos+1);
    jetalg = jetalg.substr(0,jetalg.find(".root"));
    jetalgs.push_back(jetalg);
  }
  
  // loop over all input samples
  vector<TH1F*> etHistos;
  vector<TH1F*> etaHistos;
  vector<TH1F*> phiHistos;
  vector<TH1F*> emfHistos;
  for (unsigned int i=0;i<input.size();++i) {
    TFile* f = new TFile(input[i].c_str(),"READ");
    TTree* t = (TTree*)f->Get(treename.c_str());
    cout<<jetalgs[i]<<" sample has "<<t->GetEntries()<<" entries."<<endl;


    string nameEtAll = "EtAll_"+jetalgs[i];
    string nameEtFake = "EtFake_"+jetalgs[i];
    string nameEtEff = "EtEff_"+jetalgs[i];
    
    TH1F* hEtAll  = new TH1F(nameEtAll.c_str(),nameEtAll.c_str(),
                             netbins,etmin,etmax);
    TH1F* hEtFake = new TH1F(nameEtFake.c_str(),nameEtFake.c_str(),
                             netbins,etmin,etmax);
    TH1F* hEtEff  = new TH1F(nameEtEff.c_str(),nameEtEff.c_str(),
                             netbins,etmin,etmax);
    
    hEtAll->Sumw2();
    hEtFake->Sumw2();
    hEtEff->Sumw2();
    
    string nameEtaAll = "EtaAll_"+jetalgs[i];
    string nameEtaFake = "EtaFake_"+jetalgs[i];
    string nameEtaEff = "EtaEff_"+jetalgs[i];
    
    TH1F* hEtaAll  = new TH1F(nameEtaAll.c_str(),nameEtaAll.c_str(),
                              netabins,etamin,etamax);
    TH1F* hEtaFake = new TH1F(nameEtaFake.c_str(),nameEtaFake.c_str(),
                              netabins,etamin,etamax);
    TH1F* hEtaEff  = new TH1F(nameEtaEff.c_str(),nameEtaEff.c_str(),
                              netabins,etamin,etamax);
    
    hEtaAll->Sumw2();
    hEtaFake->Sumw2();
    hEtaEff->Sumw2();
    
    string namePhiAll = "PhiAll_"+jetalgs[i];
    string namePhiFake = "PhiFake_"+jetalgs[i];
    string namePhiEff = "PhiEff_"+jetalgs[i];
    
    TH1F* hPhiAll = new TH1F(namePhiAll.c_str(),namePhiAll.c_str(),
                             nphibins,phimin,phimax);
    TH1F* hPhiFake = new TH1F(namePhiFake.c_str(),namePhiFake.c_str(),
                              nphibins,phimin,phimax);
    TH1F* hPhiEff = new TH1F(namePhiEff.c_str(),namePhiEff.c_str(),
                             nphibins,phimin,phimax);
    
    hPhiAll->Sumw2();
    hPhiFake->Sumw2();
    hPhiEff->Sumw2();
    
    string nameEmfAll = "EmfAll_"+jetalgs[i];
    string nameEmfFake = "EmfFake_"+jetalgs[i];
    string nameEmfEff = "EmfEff_"+jetalgs[i];
    
    TH1F* hEmfAll = new TH1F(nameEmfAll.c_str(),nameEmfAll.c_str(),
                             nemfbins,emfmin,emfmax);
    TH1F* hEmfFake = new TH1F(nameEmfFake.c_str(),nameEmfFake.c_str(),
                              nemfbins,emfmin,emfmax);
    TH1F* hEmfEff = new TH1F(nameEmfEff.c_str(),nameEmfEff.c_str(),
                             nemfbins,emfmin,emfmax);
    
    hEmfAll->Sumw2();
    hEmfFake->Sumw2();
    hEmfEff->Sumw2();
    
    float weight;       t->SetBranchAddress("weight",&weight);
    char  njt;          t->SetBranchAddress("njt",   &njt);
    float jtet[100];    t->SetBranchAddress("jtet",   jtet);
    float jteta[100];   t->SetBranchAddress("jteta",  jteta);
    float jtphi[100];   t->SetBranchAddress("jtphi",  jtphi);
    float jtemf[100];   t->SetBranchAddress("jtemf",  jtemf);
    float jtgendr[100]; t->SetBranchAddress("jtgendr",jtgendr);
    float jtjes[100][3];
    if (applyjes>0&&t->FindBranch("jtjes")) t->SetBranchAddress("jtjes",jtjes);

    int nevts=t->GetEntries();
    
    for (int ievt=0;ievt<nevts;ievt++) {

      t->GetEntry(ievt);

      for (int ijt=0;ijt<njt;ijt++) {
        
        float dr  = jtgendr[ijt];
        float et  = (applyjes==0) ? jtet[ijt]:jtet[ijt]*jtjes[ijt][applyjes-1];
        float eta = jteta[ijt];
        float phi = jtphi[ijt];
        float emf = jtemf[ijt]; if (TMath::IsNaN(emf)) emf = 0.0;
        
        if (et <etmin ||et >etmax)  continue;
        if (eta<etamin||eta>etamax) continue;
        if (phi<phimin||phi>phimax) continue;
        if (emf<emfmin||emf>emfmax) continue;
        
        hEtAll ->Fill(et, weight);
        hEtaAll->Fill(eta,weight);
        hPhiAll->Fill(phi,weight);
        hEmfAll->Fill(emf,weight);
        
        if (dr<0.0||dr>drmax) {
          hEtFake ->Fill(et, weight);
          hEtaFake->Fill(eta,weight);
          hPhiFake->Fill(phi,weight);
          hEmfFake->Fill(emf,weight);
        }
      }
    }
    
    hEtEff->Divide(hEtFake,hEtAll,1.0,1.0,"B");
    setBinomialErrors(hEtEff,hEtFake,hEtAll);
    etHistos.push_back(hEtEff);
    
    hEtaEff->Divide(hEtaFake,hEtaAll,1.0,1.0,"B");
    setBinomialErrors(hEtaEff,hEtaFake,hEtaAll);
    etaHistos.push_back(hEtaEff);
 
    hPhiEff->Divide(hPhiFake,hPhiAll,1.0,1.0,"B");
    setBinomialErrors(hPhiEff,hPhiFake,hPhiAll);
    phiHistos.push_back(hPhiEff);

    hEmfEff->Divide(hEmfFake,hEmfAll,1.0,1.0,"B");
    setBinomialErrors(hEmfEff,hEmfFake,hEmfAll);
    emfHistos.push_back(hEmfEff);
  }
  
  
  // make plots
  argc=1;
  TRint* app = new TRint(argv[0],&argc,argv);
  
  gStyle->SetOptStat(0);
  
  Color_t color;
  
  // emf canvas
  TCanvas* cEmfEff  = new TCanvas("FakeVsEmf", "FakeVsEmf",10,420,400,400);
  TLegend* lEmfEff  = new TLegend(0.65,0.85-emfHistos.size()*0.05,0.86,0.85);
  cEmfEff->cd();
  lEmfEff->SetFillColor(10);
  
  color = kBlack;
  for (unsigned int i=0;i<emfHistos.size();i++) {
    TH1F* h = emfHistos[i];
    lEmfEff->AddEntry(h,getLegLabel(jetalgs[i]).c_str(),"lp");
    h->SetLineWidth(2);
    h->SetMarkerColor(color);
    h->SetLineColor(color++);
    h->SetMarkerStyle(kFullCircle);
    h->SetMarkerSize(0.7);
    if (i==0) {
      h->SetTitle("#epsilon_{FAKE} vs emf");
      h->SetXTitle("emf");
      h->SetYTitle("#epsilon_{FAKE}");
      h->Draw("E");
      h->SetMinimum(0.0);
      h->SetMaximum(1.05);
    }
    else h->Draw("ESAME");
  }
  lEmfEff->Draw();

  // phi canvas
  TCanvas* cPhiEff  = new TCanvas("FakeVsPhi", "FakeVsPhi",830,10,400,400);
  TLegend* lPhiEff  = new TLegend(0.65,0.85-phiHistos.size()*0.05,0.86,0.85);
  cPhiEff->cd();
  lPhiEff->SetFillColor(10);
  
  color = kBlack;
  for (unsigned int i=0;i<phiHistos.size();i++) {
    TH1F* h = phiHistos[i];
    lPhiEff->AddEntry(h,getLegLabel(jetalgs[i]).c_str(),"lp");
    h->SetLineWidth(2);
    h->SetMarkerColor(color);
    h->SetLineColor(color++);
    h->SetMarkerStyle(kFullCircle);
    h->SetMarkerSize(0.7);
    if (i==0) {
      h->SetTitle("#epsilon_{FAKE} vs #phi");
      h->SetXTitle("#phi");
      h->SetYTitle("#epsilon_{FAKE}");
      h->Draw("E");
      h->SetMinimum(0.0);
      h->SetMaximum(1.05);
    }
    else h->Draw("ESAME");
  }
  lPhiEff->Draw();

  // eta canvas
  TCanvas* cEtaEff  = new TCanvas("FakeVsEta", "FakeVsEta",420,10,400,400);
  TLegend* lEtaEff  = new TLegend(0.65,0.85-etaHistos.size()*0.05,0.86,0.85);
  cEtaEff->cd();
  lEtaEff->SetFillColor(10);
  
  color = kBlack;
  for (unsigned int i=0;i<etaHistos.size();i++) {
    TH1F* h = etaHistos[i];
    lEtaEff->AddEntry(h,getLegLabel(jetalgs[i]).c_str(),"lp");
    h->SetLineWidth(2);
    h->SetMarkerColor(color);
    h->SetLineColor(color++);
    h->SetMarkerStyle(kFullCircle);
    h->SetMarkerSize(0.7);
    if (i==0) {
      h->SetTitle("#epsilon_{FAKE} vs #eta");
      h->SetXTitle("#eta");
      h->SetYTitle("#epsilon_{FAKE}");
      h->Draw("E");
      h->SetMinimum(0.0);
      h->SetMaximum(1.05);
    }
    else h->Draw("ESAME");
    
  }
  lEtaEff->Draw();

  // et canvas
  TCanvas* cEtEff  = new TCanvas("FakeVsEt", "FakeVsEt",10,10,400,400);
  TLegend* lEtEff  = new TLegend(0.65,0.85-etHistos.size()*0.05,0.86,0.85);
  cEtEff->cd();
  lEtEff->SetFillColor(10);
  
  color = kBlack;
  for (unsigned int i=0;i<etHistos.size();i++) {
    TH1F* h = etHistos[i];
    lEtEff->AddEntry(h,getLegLabel(jetalgs[i]).c_str(),"lp");
    h->SetLineWidth(2);
    h->SetMarkerColor(color);
    h->SetLineColor(color++);
    h->SetMarkerStyle(kFullCircle);
    h->SetMarkerSize(0.7);
    if (i==0) {
      h->SetTitle("#epsilon_{FAKE} vs E_{T}");
      h->SetXTitle("E_{T} [GeV]");
      h->SetYTitle("#epsilon_{FAKE}");
      h->Draw("E");
      h->SetMinimum(0.0);
      h->SetMaximum(1.05);
    }
    else h->Draw("ESAME");
  }
  lEtEff->Draw();
  
  app->Run();
  
  return 0;
}


////////////////////////////////////////////////////////////////////////////////
// implementation of global function
////////////////////////////////////////////////////////////////////////////////

//______________________________________________________________________________
void setBinomialErrors(TH1F* hEff,const TH1F* hFake, const TH1F* hAll)
{
  for (int i=1;i<hEff->GetNbinsX();i++) {
    float nfake = hFake->GetBinContent(i);
    float nall  = hAll ->GetBinContent(i);
    float eeff  = (nall>0.0) ? std::sqrt(nfake/(nall*nall)*(1-nfake/nall)):0.0;
    hEff->SetBinError(i,eeff);
  }
}


//______________________________________________________________________________
string getLegLabel(const string& jetalg)
{
  string result;
  string tmp(jetalg);
  
  if (tmp.find("fastjet")==0) {
    result = "fast k_{T}, D=";
    tmp=tmp.substr(7);
  }
  else if (tmp.find("icone")==0) {
    result = "icone, R=";
    tmp = tmp.substr(5);
  }
  else if (tmp.find("mcone")==0) {
    result = "mcone, R=";
    tmp = tmp.substr(5);
  }
  else return result;

  stringstream ssparam;
  ssparam<<tmp;
  float param;
  ssparam>>param;
  param = 0.1*param;
  stringstream ssparam2;
  ssparam2<<param;
  result += ssparam2.str();

  return result;
}
Revision:	1.4
Committed:	Tue Sep 4 13:31:20 2007 UTC (17 years, 7 months ago) by schiefer
Branch:	MAIN
CVS Tags:	V00-00-00
Changes since 1.3:	+45 -9 lines
Log Message:	* empty log message *
#	User	Rev	Content
1	schiefer	1.1	////////////////////////////////////////////////////////////////////////////////
2			//
3			// mcfakerate_x
4			// ------------
5			//
6			// 06/28/2007 Philipp Schieferdecker <philipp.schieferdecker@cern.ch>
7			////////////////////////////////////////////////////////////////////////////////
8
9
10			#include "utils/cmdline.h"
11
12			#include <TROOT.h>
13			#include <TStyle.h>
14			#include <TRint.h>
15			#include <TFile.h>
16			#include <TTree.h>
17			#include <TH1F.h>
18			#include <TCanvas.h>
19			#include <TLegend.h>
20	schiefer	1.3	#include <TMath.h>
21	schiefer	1.1
22			#include <iostream>
23			#include <iomanip>
24			#include <cmath>
25			#include <string>
26			#include <vector>
27
28
29			using namespace std;
30
31
32			////////////////////////////////////////////////////////////////////////////////
33			// declare global functions
34			////////////////////////////////////////////////////////////////////////////////
35
36			//______________________________________________________________________________
37	schiefer	1.4	void setBinomialErrors(TH1F* hEff,const TH1F* hFake, const TH1F* hAll);
38			string getLegLabel(const string& jetalg);
39	schiefer	1.1
40			////////////////////////////////////////////////////////////////////////////////
41			// main
42			////////////////////////////////////////////////////////////////////////////////
43
44			//______________________________________________________________________________
45			int main(int argc,char**argv)
46			{
47			cmdline cl;
48			cl.parse(argc,argv);
49
50			vector<string> input = cl.get_vector<string>("input");
51			string treename = cl.get_value<string> ("treename", "t");
52	schiefer	1.3	float drmax = cl.get_value<float> ("drmax", 0.3);
53			int netbins = cl.get_value<int> ("netbins", 50);
54	schiefer	1.1	int netabins = cl.get_value<int> ("netabins", 25);
55			int nphibins = cl.get_value<int> ("nphibins", 25);
56	schiefer	1.3	int nemfbins = cl.get_value<int> ("nemfbins", 15);
57			float etmin = cl.get_value<float> ("etmin", 0.0);
58			float etmax = cl.get_value<float> ("etmax", 200.0);
59			float etamin = cl.get_value<float> ("etamin", -5.0);
60			float etamax = cl.get_value<float> ("etamax", 5.0);
61			float phimin = cl.get_value<float> ("phimin", -3.2);
62			float phimax = cl.get_value<float> ("phimax", 3.2);
63			float emfmin = cl.get_value<float> ("emfmin", 0.0);
64			float emfmax = cl.get_value<float> ("emfmax", 1.0);
65			short applyjes = cl.get_value<short> ("applyjes", 0);
66	schiefer	1.1
67			if (!cl.check()) return 0;
68			cl.print();
69
70
71			// jetalgs
72			vector<string> jetalgs;
73			for (vector<string>::const_iterator it=input.begin();it!=input.end();++it) {
74			string jetalg = *it;
75			string::size_type pos;
76			while ((pos=jetalg.find('/'))!=string::npos) jetalg = jetalg.substr(pos+1);
77			jetalg = jetalg.substr(0,jetalg.find(".root"));
78			jetalgs.push_back(jetalg);
79			}
80
81			// loop over all input samples
82	schiefer	1.3	vector<TH1F*> etHistos;
83	schiefer	1.1	vector<TH1F*> etaHistos;
84			vector<TH1F*> phiHistos;
85	schiefer	1.3	vector<TH1F*> emfHistos;
86	schiefer	1.1	for (unsigned int i=0;i<input.size();++i) {
87			TFile* f = new TFile(input[i].c_str(),"READ");
88			TTree* t = (TTree*)f->Get(treename.c_str());
89			cout<<jetalgs[i]<<" sample has "<<t->GetEntries()<<" entries."<<endl;
90
91
92	schiefer	1.3	string nameEtAll = "EtAll_"+jetalgs[i];
93			string nameEtFake = "EtFake_"+jetalgs[i];
94			string nameEtEff = "EtEff_"+jetalgs[i];
95
96			TH1F* hEtAll = new TH1F(nameEtAll.c_str(),nameEtAll.c_str(),
97			netbins,etmin,etmax);
98			TH1F* hEtFake = new TH1F(nameEtFake.c_str(),nameEtFake.c_str(),
99			netbins,etmin,etmax);
100			TH1F* hEtEff = new TH1F(nameEtEff.c_str(),nameEtEff.c_str(),
101			netbins,etmin,etmax);
102
103			hEtAll->Sumw2();
104			hEtFake->Sumw2();
105			hEtEff->Sumw2();
106	schiefer	1.1
107			string nameEtaAll = "EtaAll_"+jetalgs[i];
108			string nameEtaFake = "EtaFake_"+jetalgs[i];
109			string nameEtaEff = "EtaEff_"+jetalgs[i];
110
111			TH1F* hEtaAll = new TH1F(nameEtaAll.c_str(),nameEtaAll.c_str(),
112			netabins,etamin,etamax);
113			TH1F* hEtaFake = new TH1F(nameEtaFake.c_str(),nameEtaFake.c_str(),
114			netabins,etamin,etamax);
115			TH1F* hEtaEff = new TH1F(nameEtaEff.c_str(),nameEtaEff.c_str(),
116			netabins,etamin,etamax);
117
118			hEtaAll->Sumw2();
119			hEtaFake->Sumw2();
120			hEtaEff->Sumw2();
121
122			string namePhiAll = "PhiAll_"+jetalgs[i];
123			string namePhiFake = "PhiFake_"+jetalgs[i];
124			string namePhiEff = "PhiEff_"+jetalgs[i];
125
126			TH1F* hPhiAll = new TH1F(namePhiAll.c_str(),namePhiAll.c_str(),
127			nphibins,phimin,phimax);
128			TH1F* hPhiFake = new TH1F(namePhiFake.c_str(),namePhiFake.c_str(),
129			nphibins,phimin,phimax);
130			TH1F* hPhiEff = new TH1F(namePhiEff.c_str(),namePhiEff.c_str(),
131			nphibins,phimin,phimax);
132
133			hPhiAll->Sumw2();
134			hPhiFake->Sumw2();
135			hPhiEff->Sumw2();
136
137	schiefer	1.3	string nameEmfAll = "EmfAll_"+jetalgs[i];
138			string nameEmfFake = "EmfFake_"+jetalgs[i];
139			string nameEmfEff = "EmfEff_"+jetalgs[i];
140
141			TH1F* hEmfAll = new TH1F(nameEmfAll.c_str(),nameEmfAll.c_str(),
142			nemfbins,emfmin,emfmax);
143			TH1F* hEmfFake = new TH1F(nameEmfFake.c_str(),nameEmfFake.c_str(),
144			nemfbins,emfmin,emfmax);
145			TH1F* hEmfEff = new TH1F(nameEmfEff.c_str(),nameEmfEff.c_str(),
146			nemfbins,emfmin,emfmax);
147
148			hEmfAll->Sumw2();
149			hEmfFake->Sumw2();
150			hEmfEff->Sumw2();
151
152	schiefer	1.1	float weight; t->SetBranchAddress("weight",&weight);
153			char njt; t->SetBranchAddress("njt", &njt);
154	schiefer	1.3	float jtet[100]; t->SetBranchAddress("jtet", jtet);
155	schiefer	1.1	float jteta[100]; t->SetBranchAddress("jteta", jteta);
156			float jtphi[100]; t->SetBranchAddress("jtphi", jtphi);
157	schiefer	1.3	float jtemf[100]; t->SetBranchAddress("jtemf", jtemf);
158	schiefer	1.1	float jtgendr[100]; t->SetBranchAddress("jtgendr",jtgendr);
159	schiefer	1.3	float jtjes[100][3];
160			if (applyjes>0&&t->FindBranch("jtjes")) t->SetBranchAddress("jtjes",jtjes);
161
162	schiefer	1.1	int nevts=t->GetEntries();
163
164			for (int ievt=0;ievt<nevts;ievt++) {
165	schiefer	1.4
166	schiefer	1.1	t->GetEntry(ievt);
167	schiefer	1.4
168	schiefer	1.1	for (int ijt=0;ijt<njt;ijt++) {
169
170	schiefer	1.3	float dr = jtgendr[ijt];
171			float et = (applyjes==0) ? jtet[ijt]:jtet[ijt]*jtjes[ijt][applyjes-1];
172			float eta = jteta[ijt];
173			float phi = jtphi[ijt];
174			float emf = jtemf[ijt]; if (TMath::IsNaN(emf)) emf = 0.0;
175
176			if (et <etmin \|\|et >etmax) continue;
177			if (eta<etamin\|\|eta>etamax) continue;
178			if (phi<phimin\|\|phi>phimax) continue;
179			if (emf<emfmin\|\|emf>emfmax) continue;
180
181			hEtAll ->Fill(et, weight);
182			hEtaAll->Fill(eta,weight);
183			hPhiAll->Fill(phi,weight);
184			hEmfAll->Fill(emf,weight);
185
186			if (dr<0.0\|\|dr>drmax) {
187			hEtFake ->Fill(et, weight);
188			hEtaFake->Fill(eta,weight);
189			hPhiFake->Fill(phi,weight);
190			hEmfFake->Fill(emf,weight);
191	schiefer	1.1	}
192			}
193			}
194
195	schiefer	1.3	hEtEff->Divide(hEtFake,hEtAll,1.0,1.0,"B");
196			setBinomialErrors(hEtEff,hEtFake,hEtAll);
197			etHistos.push_back(hEtEff);
198	schiefer	1.1
199			hEtaEff->Divide(hEtaFake,hEtaAll,1.0,1.0,"B");
200			setBinomialErrors(hEtaEff,hEtaFake,hEtaAll);
201			etaHistos.push_back(hEtaEff);
202
203			hPhiEff->Divide(hPhiFake,hPhiAll,1.0,1.0,"B");
204			setBinomialErrors(hPhiEff,hPhiFake,hPhiAll);
205			phiHistos.push_back(hPhiEff);
206	schiefer	1.3
207			hEmfEff->Divide(hEmfFake,hEmfAll,1.0,1.0,"B");
208			setBinomialErrors(hEmfEff,hEmfFake,hEmfAll);
209			emfHistos.push_back(hEmfEff);
210	schiefer	1.1	}
211
212
213			// make plots
214			argc=1;
215			TRint* app = new TRint(argv[0],&argc,argv);
216
217			gStyle->SetOptStat(0);
218
219			Color_t color;
220
221	schiefer	1.3	// emf canvas
222			TCanvas* cEmfEff = new TCanvas("FakeVsEmf", "FakeVsEmf",10,420,400,400);
223			TLegend* lEmfEff = new TLegend(0.65,0.85-emfHistos.size()*0.05,0.86,0.85);
224			cEmfEff->cd();
225			lEmfEff->SetFillColor(10);
226
227			color = kBlack;
228			for (unsigned int i=0;i<emfHistos.size();i++) {
229			TH1F* h = emfHistos[i];
230	schiefer	1.4	lEmfEff->AddEntry(h,getLegLabel(jetalgs[i]).c_str(),"lp");
231	schiefer	1.3	h->SetLineWidth(2);
232			h->SetMarkerColor(color);
233			h->SetLineColor(color++);
234			h->SetMarkerStyle(kFullCircle);
235			h->SetMarkerSize(0.7);
236			if (i==0) {
237			h->SetTitle("#epsilon_{FAKE} vs emf");
238	schiefer	1.4	h->SetXTitle("emf");
239	schiefer	1.3	h->SetYTitle("#epsilon_{FAKE}");
240			h->Draw("E");
241			h->SetMinimum(0.0);
242			h->SetMaximum(1.05);
243			}
244			else h->Draw("ESAME");
245			}
246			lEmfEff->Draw();
247
248	schiefer	1.1	// phi canvas
249			TCanvas* cPhiEff = new TCanvas("FakeVsPhi", "FakeVsPhi",830,10,400,400);
250			TLegend* lPhiEff = new TLegend(0.65,0.85-phiHistos.size()*0.05,0.86,0.85);
251			cPhiEff->cd();
252			lPhiEff->SetFillColor(10);
253
254			color = kBlack;
255			for (unsigned int i=0;i<phiHistos.size();i++) {
256			TH1F* h = phiHistos[i];
257	schiefer	1.4	lPhiEff->AddEntry(h,getLegLabel(jetalgs[i]).c_str(),"lp");
258	schiefer	1.3	h->SetLineWidth(2);
259	schiefer	1.1	h->SetMarkerColor(color);
260			h->SetLineColor(color++);
261			h->SetMarkerStyle(kFullCircle);
262	schiefer	1.3	h->SetMarkerSize(0.7);
263	schiefer	1.1	if (i==0) {
264	schiefer	1.3	h->SetTitle("#epsilon_{FAKE} vs #phi");
265	schiefer	1.4	h->SetXTitle("#phi");
266	schiefer	1.3	h->SetYTitle("#epsilon_{FAKE}");
267	schiefer	1.1	h->Draw("E");
268			h->SetMinimum(0.0);
269			h->SetMaximum(1.05);
270			}
271			else h->Draw("ESAME");
272			}
273			lPhiEff->Draw();
274
275			// eta canvas
276	schiefer	1.2	TCanvas* cEtaEff = new TCanvas("FakeVsEta", "FakeVsEta",420,10,400,400);
277	schiefer	1.1	TLegend* lEtaEff = new TLegend(0.65,0.85-etaHistos.size()*0.05,0.86,0.85);
278			cEtaEff->cd();
279			lEtaEff->SetFillColor(10);
280
281			color = kBlack;
282			for (unsigned int i=0;i<etaHistos.size();i++) {
283			TH1F* h = etaHistos[i];
284	schiefer	1.4	lEtaEff->AddEntry(h,getLegLabel(jetalgs[i]).c_str(),"lp");
285	schiefer	1.3	h->SetLineWidth(2);
286	schiefer	1.1	h->SetMarkerColor(color);
287			h->SetLineColor(color++);
288			h->SetMarkerStyle(kFullCircle);
289	schiefer	1.3	h->SetMarkerSize(0.7);
290	schiefer	1.1	if (i==0) {
291	schiefer	1.3	h->SetTitle("#epsilon_{FAKE} vs #eta");
292	schiefer	1.4	h->SetXTitle("#eta");
293	schiefer	1.3	h->SetYTitle("#epsilon_{FAKE}");
294	schiefer	1.1	h->Draw("E");
295			h->SetMinimum(0.0);
296			h->SetMaximum(1.05);
297			}
298			else h->Draw("ESAME");
299
300			}
301			lEtaEff->Draw();
302
303	schiefer	1.3	// et canvas
304			TCanvas* cEtEff = new TCanvas("FakeVsEt", "FakeVsEt",10,10,400,400);
305			TLegend* lEtEff = new TLegend(0.65,0.85-etHistos.size()*0.05,0.86,0.85);
306			cEtEff->cd();
307			lEtEff->SetFillColor(10);
308	schiefer	1.1
309			color = kBlack;
310	schiefer	1.3	for (unsigned int i=0;i<etHistos.size();i++) {
311			TH1F* h = etHistos[i];
312	schiefer	1.4	lEtEff->AddEntry(h,getLegLabel(jetalgs[i]).c_str(),"lp");
313	schiefer	1.3	h->SetLineWidth(2);
314	schiefer	1.1	h->SetMarkerColor(color);
315			h->SetLineColor(color++);
316			h->SetMarkerStyle(kFullCircle);
317	schiefer	1.3	h->SetMarkerSize(0.7);
318	schiefer	1.1	if (i==0) {
319	schiefer	1.3	h->SetTitle("#epsilon_{FAKE} vs E_{T}");
320	schiefer	1.4	h->SetXTitle("E_{T} [GeV]");
321	schiefer	1.3	h->SetYTitle("#epsilon_{FAKE}");
322	schiefer	1.1	h->Draw("E");
323			h->SetMinimum(0.0);
324			h->SetMaximum(1.05);
325			}
326			else h->Draw("ESAME");
327			}
328	schiefer	1.3	lEtEff->Draw();
329	schiefer	1.1
330			app->Run();
331
332			return 0;
333			}
334
335
336			////////////////////////////////////////////////////////////////////////////////
337			// implementation of global function
338			////////////////////////////////////////////////////////////////////////////////
339
340			//______________________________________________________________________________
341			void setBinomialErrors(TH1F* hEff,const TH1F* hFake, const TH1F* hAll)
342			{
343			for (int i=1;i<hEff->GetNbinsX();i++) {
344	schiefer	1.3	float nfake = hFake->GetBinContent(i);
345			float nall = hAll ->GetBinContent(i);
346			float eeff = (nall>0.0) ? std::sqrt(nfake/(nallnall)(1-nfake/nall)):0.0;
347	schiefer	1.1	hEff->SetBinError(i,eeff);
348			}
349			}
350	schiefer	1.4
351
352			//______________________________________________________________________________
353			string getLegLabel(const string& jetalg)
354			{
355			string result;
356			string tmp(jetalg);
357
358			if (tmp.find("fastjet")==0) {
359			result = "fast k_{T}, D=";
360			tmp=tmp.substr(7);
361			}
362			else if (tmp.find("icone")==0) {
363			result = "icone, R=";
364			tmp = tmp.substr(5);
365			}
366			else if (tmp.find("mcone")==0) {
367			result = "mcone, R=";
368			tmp = tmp.substr(5);
369			}
370			else return result;
371
372			stringstream ssparam;
373			ssparam<<tmp;
374			float param;
375			ssparam>>param;
376			param = 0.1*param;
377			stringstream ssparam2;
378			ssparam2<<param;
379			result += ssparam2.str();
380
381			return result;
382			}