PatExamples/test/patElectron_eidEfficiency.C

0001 void patElectron_eidEfficiency()
0002 {
0003   // define proper canvas style
0004   setNiceStyle();
0005   gStyle->SetOptStat(0);
0006
0007   // open file
0008   TFile* file = new TFile("analyzePatElectron.root");
0009
0010   // get probe histograms
0011   // TH1F* proPt_  = file->Get("tightElectronID/pt");
0012   // TH1F* proEta_ = file->Get("tightElectronID/eta");
0013   // TH1F* proPhi_ = file->Get("tightElectronID/phi");
0014
0015   // get probe histograms
0016   TH1F* proPt_  = file->Get("looseElectronID/pt");
0017   TH1F* proEta_ = file->Get("looseElectronID/eta");
0018   TH1F* proPhi_ = file->Get("looseElectronID/phi");
0019
0020   // get reference histograms
0021   TH1F* refPt_  = file->Get("plainElectronID/pt");
0022   TH1F* refEta_ = file->Get("plainElectronID/eta");
0023   TH1F* refPhi_ = file->Get("plainElectronID/phi");
0024
0025   // draw canvas with electron kinematics
0026   TCanvas* canv0 = new TCanvas("canv0", "electron kinemtics", 600, 300);
0027   canv0->Divide(2,1);
0028   canv0->cd(1);
0029   TH1F* kinPt_=proPt_->Clone();
0030   kinPt_->SetFillStyle(3005.);
0031   kinPt_->SetFillColor(4.);
0032   setHistStyle(kinPt_);
0033   kinPt_ ->DrawCopy();
0034
0035   canv0->cd(2);
0036   TH1F* kinEta_=proEta_->Clone();
0037   kinEta_->SetFillStyle(3005.);
0038   kinEta_->SetFillColor(4.);
0039   setHistStyle(kinEta_);
0040   kinEta_->DrawCopy();
0041
0042   // draw canvas with electronID efficiencies
0043   TCanvas* canv1 = new TCanvas("canv1", "electron ID efficiency", 600, 600);
0044   canv1->Divide(2,2);
0045   if(correlatedError(proPt_,  refPt_ )==0){
0046     canv1->cd(1);
0047     canv1->GetPad(1)->SetGridx(1);
0048     canv1->GetPad(1)->SetGridy(1);
0049     proPt_ ->SetMaximum(1.3);
0050     setHistStyle(proPt_);
0051     proPt_ ->DrawCopy();
0052   }
0053   if(correlatedError(proEta_, refEta_)==0){
0054     canv1->cd(2);
0055     canv1->GetPad(2)->SetGridx(1);
0056     canv1->GetPad(2)->SetGridy(1);
0057     proEta_ ->SetMaximum(1.3);
0058     setHistStyle(proEta_);
0059     proEta_->DrawCopy();
0060   }
0061   if(correlatedError(proPhi_, refPhi_)==0){
0062     canv1->cd(3);
0063     canv1->GetPad(3)->SetGridx(1);
0064     canv1->GetPad(3)->SetGridy(1);
0065     proPhi_ ->SetMaximum(1.3);
0066     setHistStyle(proPhi_);
0067     proPhi_->DrawCopy();
0068   }
0069 }
0070
0071 int correlatedError(TH1F* nominator, TH1F* denominator)
0072 {
0073   // --------------------------------------------------
0074   // get error of correlated ratio for two histograms
0075   // for gaussian distributed quantities the formular
0076   //
0077   //  * de=e*Sqrt((dr/r)^2+(1-2e)*(dv/d)^2)
0078   //
0079   // turns automatically into
0080   //
0081   //  * de=Sqrt(e*(1-e)/r)
0082   // --------------------------------------------------
0083   if(nominator->GetNbinsX()!=denominator->GetNbinsX()){
0084     // these histogramsd do not correspond to each other
0085     return -1;
0086   }
0087
0088   // loop over the denominator bins
0089   for(Int_t i=0; i<denominator->GetNbinsX(); ++i){
0090     float dval =  nominator->GetBinError(i+1);
0091     float val  =  nominator->GetBinContent(i+1);
0092     float dref =  denominator->GetBinError(i+1);
0093     float ref  =  denominator->GetBinContent(i+1);
0094
0095     float err;
0096     if(val<=0){
0097       // val should never be smaller than 0
0098       err=0; continue;
0099     }
0100     if(ref==0){
0101       // ref should never be 0
0102       err=0; continue;
0103     }
0104
0105     if(val/ref<1)
0106       err=(val/ref)*TMath::Sqrt(TMath::Abs((dref*dref)/(ref*ref)+(1.-2.*(val/ref))*(dval*dval)/(val*val)));
0107     else
0108       err=(ref/val)*TMath::Sqrt(TMath::Abs((dval*dval)/(val*val)+(1.-2.*(ref/val))*(dref*dref)/(ref*ref)));
0109
0110     // set proper values and errors
0111     nominator->SetBinContent(i+1, val/ref);
0112     nominator->SetBinError(i+1, err);
0113   }
0114   return 0;
0115 }
0116
0117 void setAxisStyle(TH1* hist) {
0118   // --------------------------------------------------
0119   // define proper axsis style for a given histogram
0120   // --------------------------------------------------
0121     hist->GetXaxis()->SetTitleSize( 0.06);
0122     hist->GetXaxis()->SetTitleColor( 1);
0123     hist->GetXaxis()->SetTitleOffset( 0.8);
0124     hist->GetXaxis()->SetTitleFont( 62);
0125     hist->GetXaxis()->SetLabelSize( 0.05);
0126     hist->GetXaxis()->SetLabelFont( 62);
0127     hist->GetXaxis()->CenterTitle();
0128     hist->GetXaxis()->SetNdivisions( 505);
0129
0130     hist->GetYaxis()->SetTitleSize( 0.07);
0131     hist->GetYaxis()->SetTitleColor( 1);
0132     hist->GetYaxis()->SetTitleOffset( 0.5);
0133     hist->GetYaxis()->SetTitleFont( 62);
0134     hist->GetYaxis()->SetLabelSize( 0.05);
0135     hist->GetYaxis()->SetLabelFont( 62);
0136 }
0137
0138 void setHistStyle(TH1F* hist)
0139 {
0140   // --------------------------------------------------
0141   // define proper histogram style
0142   // --------------------------------------------------
0143   setAxisStyle(hist);
0144   hist->GetXaxis()->SetTitle(hist->GetTitle());
0145   hist->SetTitle();
0146   hist->SetLineColor(4.);
0147   hist->SetLineWidth(3.);
0148   hist->SetMarkerSize(0.75);
0149   hist->SetMarkerColor(4.);
0150   hist->SetMarkerStyle(20.);
0151 }
0152
0153 void setNiceStyle()
0154 {
0155   // --------------------------------------------------
0156   // define proper canvas style
0157   // --------------------------------------------------
0158   TStyle *MyStyle = new TStyle ("MyStyle", "My style for nicer plots");
0159
0160   Float_t xoff = MyStyle->GetLabelOffset("X"),
0161           yoff = MyStyle->GetLabelOffset("Y"),
0162           zoff = MyStyle->GetLabelOffset("Z");
0163
0164   MyStyle->SetCanvasBorderMode ( 0 );
0165   MyStyle->SetPadBorderMode    ( 0 );
0166   MyStyle->SetPadColor         ( 0 );
0167   MyStyle->SetCanvasColor      ( 0 );
0168   MyStyle->SetTitleColor       ( 0 );
0169   MyStyle->SetStatColor        ( 0 );
0170   MyStyle->SetTitleBorderSize  ( 0 );
0171   MyStyle->SetTitleFillColor   ( 0 );
0172   MyStyle->SetTitleH        ( 0.07 );
0173   MyStyle->SetTitleW        ( 1.00 );
0174   MyStyle->SetTitleFont     (  132 );
0175
0176   MyStyle->SetLabelOffset (1.5*xoff, "X");
0177   MyStyle->SetLabelOffset (1.5*yoff, "Y");
0178   MyStyle->SetLabelOffset (1.5*zoff, "Z");
0179
0180   MyStyle->SetTitleOffset (0.9,      "X");
0181   MyStyle->SetTitleOffset (0.9,      "Y");
0182   MyStyle->SetTitleOffset (0.9,      "Z");
0183
0184   MyStyle->SetTitleSize   (0.045,    "X");
0185   MyStyle->SetTitleSize   (0.045,    "Y");
0186   MyStyle->SetTitleSize   (0.045,    "Z");
0187
0188   MyStyle->SetLabelFont   (132,      "X");
0189   MyStyle->SetLabelFont   (132,      "Y");
0190   MyStyle->SetLabelFont   (132,      "Z");
0191
0192   MyStyle->SetPalette(1);
0193
0194   MyStyle->cd();
0195 }