ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/Development/Plotting/Modules/GeneralToolBox.C
(Generate patch)

Comparing UserCode/cbrown/Development/Plotting/Modules/GeneralToolBox.C (file contents):
Revision 1.11 by buchmann, Wed Jun 6 08:24:13 2012 UTC vs.
Revision 1.43 by buchmann, Sun May 12 20:22:53 2013 UTC

# Line 39 | Line 39
39   #include <TPaveText.h>
40   #include <TRandom.h>
41   #include <TGraphErrors.h>
42 + #include <TROOT.h>
43 + #include <TPolyLine.h>
44 +
45 + #include "ValueClass.C"
46   #ifndef Verbosity
47   #define Verbosity 0
48   #endif
# Line 61 | Line 65 | namespace PlottingSetup {
65      int PaperMode=0; // PaperMode=true will suppress "Preliminary" in DrawPrelim()
66      int Approved=0; // Approved=true will only plot approved plots
67      bool is2012=true;
68 +    bool is53reco=true;
69 +    bool openBox = true;
70 +    vector<TH1F*> FakeHistoHeap;
71 +    bool DrawMetSignalRegionMllLines=false; // whether to draw the lines in mll plots for the MET signal region
72 +    float ConsideredZWidth=20;
73   }
74  
75   bool dopng=false;
# Line 87 | Line 96 | void write_error(string funcname, string
96   void write_info(string funcname, string text);
97   string get_directory();
98   bool Contains(string wholestring, string findme);
99 + TH1F* CollapseStack(THStack stack,TString hname);
100 + TH1F* CollapseStack(THStack* stack,TString hname);
101   //-------------------------------------------------------------------------------------
102  
103   template<typename U>
# Line 254 | Line 265 | void ensure_directory_exists(string this
265    }
266   }
267  
268 + void DeleteStack(THStack mcstack) {
269 +  TH1F* h;
270 +  TIter nextHist(mcstack.GetHists());
271 +  while ( h = (TH1F*)nextHist() ) {
272 +    h->Delete();
273 +  }
274 +  delete h;
275 + }
276 +
277 +
278 + void DeadEnd(string message) {
279 +  dout << endl;
280 +  dout << "\033[1;31m     _                _                  _  " << endl;
281 +  dout << "\033[1;31m  __| | ___  __ _  __| |   ___ _ __   __| | " << endl;
282 +  dout << "\033[1;31m / _` |/ _ \\/ _` |/ _` |  / _ \\ '_ \\ / _` | " << endl;
283 +  dout << "\033[1;31m| (_| |  __/ (_| | (_| | |  __/ | | | (_| | " << endl;
284 +  dout << "\033[1;31m \\__,_|\\___|\\__,_|\\__,_|  \\___|_| |_|\\__,_| " << endl;
285 +  dout << "\033[1;31m                                            " << endl;
286 +  dout << "" << endl;
287 +  dout << " A totally fatal error has occurred: " << endl;
288 +  dout << " " << message << endl;
289 +  dout << endl;
290 +  dout << "If you do not know how to deal with this error message please call Customer Support Mumbai (or Zurich...) \033[0m" << endl;
291 +  dout << "   ... aborting now ... " << endl;
292 +  assert(0);
293 + }
294 +
295 +
296 + void SanityChecks() {
297 +  //this function gets called whenever samples are initialized (see ActiveSamples.C)
298 +  //whatever is necessary for the code to run should be put here.
299 +  
300 +  if(gROOT->GetVersionInt()<53200) {
301 +    DeadEnd("ROOT version outdated, you are using "+any2string(gROOT->GetVersionInt())+" but the earliest recommended version is 5.32.00");
302 +  }
303 +  
304 +  
305 + }
306 +
307   void initialize_log() {
308    if(!PlottingSetup::publicmode) {
309      dout << "____________________________________________________________" << endl;
# Line 393 | Line 443 | void CompleteSave(TVirtualPad *can, stri
443  
444  
445   void setlumi(float l) {
446 < generaltoolboxlumi=l;
446 >  generaltoolboxlumi=l;
447   }
448  
449   int write_first_line(vector<vector<string> > &entries) {
# Line 551 | Line 601 | void DrawPrelim(float writelumi=generalt
601    string prel=" Preliminary";
602    if(PlottingSetup::PaperMode) prel="";
603    //prelimtext << "CMS Preliminary 2011 , #sqrt{s}= 7 TeV, L= O(1) fb^{-1}"; //temporary replacement
604 + //  if(PlottingSetup::is53reco) prel += " 53X";
605    string energy="7 TeV";
606    if(PlottingSetup::is2012) energy="8 TeV";
607    if(writelumi == 0) {
608      if(isMC) prelimtext << "CMS Simulation, #sqrt{s} = " << energy;
609      else prelimtext << "CMS" << prel << ", #sqrt{s} = " << energy;
610    } else {
611 <    if(isMC) prelimtext << "CMS Simulation, #sqrt{s} = " << energy << ", L_{int} = " << std::setprecision(3) <<writelumi<<" "<<barn<<"^{-1}";
612 <     else prelimtext << "CMS" << prel << ", #sqrt{s} = " << energy << ", L_{int} = " << std::setprecision(3) <<writelumi<<" "<<barn<<"^{-1}";
611 >    if(isMC) prelimtext << "CMS Simulation, #sqrt{s} = " << energy << ", #scale[0.6]{#int}L dt = " << std::setprecision(3) <<writelumi<<" "<<barn<<"^{-1}";
612 >     else prelimtext << "CMS" << prel << ", #sqrt{s} = " << energy << ", #scale[0.6]{#int}L dt = " << std::setprecision(3) <<writelumi<<" "<<barn<<"^{-1}";
613    }
614    TPaveText *eventSelectionPaveText = new TPaveText(0.27, 0.93,0.77, 1.0,"blNDC");
615    eventSelectionPaveText->SetFillStyle(4000);
# Line 636 | Line 687 | float computeProductError(float a, float
687      return val;
688   }
689  
690 + TGraphAsymmErrors *SysAndStatRatio(TH1F *h1,TH1F *h2, TH1F *sys, int id, vector<float>binning, bool precise=false) {
691 +    int absJZBbinsNumber = binning.size()-1;
692 +    TGraphAsymmErrors* graph = new TGraphAsymmErrors(absJZBbinsNumber);
693 +    for(int i=0;i<absJZBbinsNumber;i++)
694 +        {
695 +                float xCenter=h1->GetBinCenter(i+1);
696 +                float xWidth=(h1->GetBinWidth(i+1))*0.5;
697 +                float nominatorError = h1->GetBinError(i+1);
698 +                float nominator=h1->GetBinContent(i+1);
699 +                float denominatorError=h2->GetBinError(i+1);
700 +                float denominator=h2->GetBinContent(i+1);
701 +                float errorN = 0;
702 +                float errorP = computeRatioError(nominator,nominatorError,denominator,denominatorError);
703 +                float syserr=sys->GetBinContent(i+1);
704 +                float errorsys=0;
705 +                if(id==1) // (is data)
706 +                {
707 +                        if(!precise) errorP = computeRatioError(nominator,statErrorP(nominator),denominator,statErrorP(denominator));
708 +                        else         errorP = computeRatioError(nominator,nominatorError,denominator,denominatorError);
709 +                        errorsys=computeRatioError(nominator,syserr,denominator,0);
710 +                        errorP = TMath::Sqrt(errorsys*errorsys+errorP*errorP);
711 +                        errorN = errorP; // symmetrize using statErrorP
712 +                } else {
713 +                        errorN = computeRatioError(nominator,nominatorError,denominator,denominatorError);
714 +                        errorsys=computeRatioError(nominator,syserr,denominator,0);
715 +                        errorN = TMath::Sqrt(errorsys*errorsys+errorP*errorP);
716 +                        errorP = errorN; // symmetrize using statErrorP
717 +                }
718 +                if(denominator!=0) {
719 +                        graph->SetPoint(i, xCenter, nominator/denominator);
720 +                        graph->SetPointError(i,xWidth,xWidth,errorN,errorP);
721 +                }
722 +                else {
723 +                        graph->SetPoint(i, xCenter, -999);
724 +                        graph->SetPointError(i,xWidth,xWidth,errorN,errorP);
725 +                }
726 +        }      
727 +        return graph;
728 + }
729 +    
730 +
731   TGraphAsymmErrors *histRatio(TH1F *h1,TH1F *h2, int id, vector<float>binning, bool precise=false)
732   {
733          int absJZBbinsNumber = binning.size()-1;
# Line 872 | Line 964 | string interpret_cut(string incut, bool
964    
965      if(incut=="id1==0") return "ee";
966      if(incut=="id1==1") return "#mu#mu";
967 <    if(incut=="abs(mll-91.2)<20") return "|m_{l^{+}l^{-}}-m_{Z}|<20";
967 >    if(incut=="abs(mll-91)<10") return "|m_{l^{+}l^{-}}-m_{Z}|<10";
968      if(incut=="pfJetGoodID[0]") return "";
969      if(incut=="pfJetGoodID[1]") return "";
970      if((int)incut.find("pfJetGoodNum")>-1) {
# Line 1148 | Line 1240 | bool Contains(string wholestring, string
1240    else return false;
1241   }
1242  
1243 +
1244   //////////////////////////////////////////////////////////////////////////////
1245   //
1246   // http://stackoverflow.com/questions/669438/how-to-get-memory-usage-at-run-time-in-c
# Line 1160 | Line 1253 | bool Contains(string wholestring, string
1253   /* usage:
1254   double vm2, rss2;
1255   process_mem_usage(vm2, rss2);
1256 < cout << "Memory usage: VM: " << vm << "; RSS: " << rss << endl;
1256 > cout << "Memory usage: VM: " << vm2 << "; RSS: " << rss2 << endl;
1257   */
1258  
1259   void process_mem_usage(double& vm_usage, double& resident_set)
# Line 1219 | Line 1312 | TGraphAsymmErrors* produce_ratio_graph(T
1312   }
1313  
1314  
1315 + TLatex* WriteSelection(int njets) {
1316 +  string sel="Loose selection";
1317 +  if(njets==3) sel="Low E_{T}^{miss} selection";
1318 +  assert(njets==2||njets==3);
1319 +  TLatex *sele = new TLatex(0.97,0.135,sel.c_str());
1320 +  sele->SetNDC(true);
1321 +  sele->SetTextColor(TColor::GetColor("#848484"));
1322 +  sele->SetTextFont(42);
1323 +  sele->SetTextAlign(32);
1324 +  sele->SetTextSize(0.03);
1325 +  sele->SetTextAngle(270);
1326 +  return sele;
1327 + }
1328 +
1329   Double_t MarcosChi2TestX(const TH1* h1, const TH1* h2, Double_t &chi2, Int_t &ndf, Int_t &igood, Option_t *option)
1330   {
1331  
# Line 1241 | Line 1348 | Double_t MarcosChi2TestX(const TH1* h1,
1348     return prob;
1349   }
1350  
1351 < void save_with_ratio(TH1F *nominator, TH1F *denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false, bool extendrange=false, string yaxistitle="ratio") {
1351 > void Save_With_Ratio(TH1F *nominator, TH1F *denominator, TVirtualPad *orig_canvas, string savemeas, bool do_bpred_ratio=false, bool extendrange=false, string yaxistitle="ratio",TH1F *syshisto=NULL) {
1352    //this function saves the pad being passed as well as a new one including the ratio.
1353 <  CompleteSave(canvas,savemeas);
1354 <
1353 >  
1354 >  orig_canvas->cd();
1355 >  orig_canvas->Update();
1356 >  CompleteSave(orig_canvas,savemeas);
1357 >  TVirtualPad *canvas = (TVirtualPad*) orig_canvas->Clone("TempCanvas");//otherwise the Ratio_main_canvas will own our pad and destroy it upon deletion
1358 >  
1359    float bottommargin=gStyle->GetPadBottomMargin();
1360    float canvas_height=gStyle->GetCanvasDefH();
1361    float canvas_width=gStyle->GetCanvasDefW();
# Line 1255 | Line 1366 | void save_with_ratio(TH1F *nominator, TH
1366    
1367    float xstretchfactor=((1-ratiospace)*(1-gStyle->GetPadTopMargin()))/((1)*ratiospace);
1368    
1369 <  TCanvas *main_canvas = new TCanvas("main_canvas","main_canvas",(Int_t)canvas_width,(Int_t)(canvas_height*(1+ratiospace)));
1369 >  TCanvas *Ratio_main_canvas = new TCanvas("Ratio_main_canvas","Ratio_main_canvas",(Int_t)canvas_width,(Int_t)(canvas_height*(1+ratiospace)));
1370    TPad *mainpad = new TPad("mainpad","mainpad",0,1-(1.0/(1+ratiospace)),1,1);//top (main) pad
1371    TPad *coverpad = new TPad("coverpad","coverpad",gStyle->GetPadLeftMargin()-0.008,1-(1.0/(1+ratiospace))-0.04,1,1-(1.0/(1+ratiospace))+0.103);//pad covering up the x scale
1372    TPad *bottompad = new TPad("bottompad", "Ratio Pad",0,0,1,(1-(1-bottommargin)/(1+ratiospace))-0.015); //bottom pad
1373    
1374 <  main_canvas->Range(0,0,1,1);
1375 <  main_canvas->SetBorderSize(0);
1376 <  main_canvas->SetFrameFillColor(0);
1266 <  
1374 >  Ratio_main_canvas->Range(0,0,1,1);
1375 >  Ratio_main_canvas->SetBorderSize(0);
1376 >  Ratio_main_canvas->SetFrameFillColor(0);
1377    mainpad->Draw();
1378    mainpad->cd();
1379    mainpad->Range(0,0,1,1);
# Line 1273 | Line 1383 | void save_with_ratio(TH1F *nominator, TH
1383    canvas->Range(0,0,1,1);
1384    canvas->Draw("same");
1385    mainpad->Modified();
1386 <  main_canvas->cd();
1386 >  Ratio_main_canvas->cd();
1387    coverpad->Draw();
1388    coverpad->cd();
1389    coverpad->Range(0,0,1,1);
# Line 1281 | Line 1391 | void save_with_ratio(TH1F *nominator, TH
1391    coverpad->SetBorderSize(0);
1392    coverpad->SetFrameFillColor(0);
1393    coverpad->Modified();
1394 <  main_canvas->cd();
1394 >  Ratio_main_canvas->cd();
1395    bottompad->SetTopMargin(ratiotopmargin);
1396    bottompad->SetBottomMargin(ratiobottommargin);
1397    bottompad->Draw();
# Line 1290 | Line 1400 | void save_with_ratio(TH1F *nominator, TH
1400    bottompad->SetFillColor(kWhite);
1401    TH1F *ratio = (TH1F*)nominator->Clone(GetNumericHistoName().c_str());
1402    ratio->Divide(denominator);
1293
1294
1403    TGraphAsymmErrors *eratio;
1404 +  TH1F *SystDown;
1405 +  TH1F *SystUp;
1406 +  
1407    if(!do_bpred_ratio) eratio = produce_ratio_graph(ratio);
1408    else {
1409       bool using_data=false;
1410       if((int)savemeas.find("Data")>0) using_data=true;
1411       eratio = histRatio(nominator,denominator,using_data,PlottingSetup::global_ratio_binning,false);
1412 +     if(syshisto!=0) {
1413 +       SystDown = (TH1F*)syshisto->Clone("SystDown");
1414 +       SystUp = (TH1F*)syshisto->Clone("SystUp");
1415 +       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1416 +         SystDown->SetBinContent(i,1-syshisto->GetBinContent(i));
1417 +         SystUp->SetBinContent(i,1+syshisto->GetBinContent(i));
1418 +       }
1419 +       SystDown->SetLineColor(TColor::GetColor("#006DE1"));
1420 +       SystUp->SetLineColor(TColor::GetColor("#006DE1"));
1421 +     }
1422       for(int i=1;i<=ratio->GetNbinsX();i++) {
1423          ratio->SetBinContent(i,0);
1424          ratio->SetBinError(i,0);
1425       }
1426 +    
1427 +  }
1428 +  
1429 +  if(syshisto && !do_bpred_ratio) {
1430 +       SystDown = (TH1F*)syshisto->Clone("SystDown");
1431 +       SystUp = (TH1F*)syshisto->Clone("SystUp");
1432 +       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1433 +         SystDown->SetBinContent(i,1-syshisto->GetBinContent(i));
1434 +         SystUp->SetBinContent(i,1+syshisto->GetBinContent(i));
1435 +       }
1436 +       SystDown->SetLineColor(TColor::GetColor("#006DE1"));
1437 +       SystUp->SetLineColor(TColor::GetColor("#006DE1"));
1438    }
1439    eratio->SetFillColor(TColor::GetColor("#00ADE1"));
1307
1440    ratio->SetTitle("");
1441 <  ratio->GetYaxis()->SetRangeUser(0.5,1.5);
1441 >  ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1442    if(do_bpred_ratio) ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1443    if(extendrange) ratio->GetYaxis()->SetRangeUser(0.0,4.0);
1444    ratio->GetXaxis()->SetTitle(nominator->GetXaxis()->GetTitle());
# Line 1322 | Line 1454 | void save_with_ratio(TH1F *nominator, TH
1454    ratio->GetYaxis()->SetNdivisions(502,false);
1455    ratio->SetFillColor(TColor::GetColor("#58D3F7"));
1456    ratio->SetMarkerSize(0);
1457 <  ratio->Draw("e2");
1458 <  eratio->Draw("2");
1457 >  ratio->Draw("e20");
1458 >  
1459 >  TGraphAsymmErrors *ratio_center = (TGraphAsymmErrors*)eratio->Clone("ratio_center");
1460 >  for(int i=0;i<ratio_center->GetN();i++) {
1461 >    ratio_center->SetPointError(i,ratio_center->GetErrorXlow(i),ratio_center->GetErrorXhigh(i),0.005,0.005);
1462 >  }
1463 >  ratio_center->SetFillColor(TColor::GetColor("#006381"));
1464 >  
1465 >  if(syshisto!=0) {
1466 > //    sysratio->Draw("2");
1467 > //    eratio->Draw("2,same");
1468 >    eratio->Draw("20");
1469 >    SystDown->Draw("histo,same");
1470 >    SystUp->Draw("histo,same");
1471 >  } else {
1472 >    eratio->Draw("20");
1473 >  }
1474 >  ratio_center->Draw("20");
1475    ratio->Draw("same,axis");
1476    TLine *oneline = new TLine(ratio->GetXaxis()->GetBinLowEdge(1),1,ratio->GetXaxis()->GetBinLowEdge(ratio->GetNbinsX())+ratio->GetXaxis()->GetBinWidth(ratio->GetNbinsX()),1);
1477    oneline->SetLineStyle(2);
1478    oneline->SetLineColor(kBlue);
1479    oneline->Draw("same");
1480 <
1481 <  main_canvas->cd();
1482 <  main_canvas->Modified();
1483 <  main_canvas->cd();
1484 <  main_canvas->SetSelected(main_canvas);
1480 >  if(PlottingSetup::DrawMetSignalRegionMllLines) {
1481 >    cout << "Drawing extra lines in ratio this time around ... " << endl;
1482 >    float RatioYMax=2.0;
1483 >    if(extendrange) RatioYMax=4.0;
1484 >    
1485 >    TLine *SRline   = new TLine(70,0,70,RatioYMax);
1486 >    TLine *ZLowLine = new TLine(91.2-PlottingSetup::ConsideredZWidth,0,91.2-PlottingSetup::ConsideredZWidth,RatioYMax);
1487 >    TLine *ZHiLine  = new TLine(91.2+PlottingSetup::ConsideredZWidth,0,91.2+PlottingSetup::ConsideredZWidth,RatioYMax);
1488 >  
1489 >    SRline->SetLineStyle(2);
1490 >    ZLowLine->SetLineStyle(2);
1491 >    ZHiLine->SetLineStyle(2);
1492 >    SRline->SetLineColor(kGray+2);
1493 >    ZLowLine->SetLineColor(kGray+2);
1494 >    ZHiLine->SetLineColor(kGray+2);
1495 >    SRline->Draw();
1496 >    ZLowLine->Draw();
1497 >    ZHiLine->Draw();
1498 >  }
1499 >      
1500 >  Ratio_main_canvas->cd();
1501 >  Ratio_main_canvas->Modified();
1502 >  Ratio_main_canvas->cd();
1503 >  Ratio_main_canvas->SetSelected(Ratio_main_canvas);
1504    
1505 <  CompleteSave(main_canvas,savemeas+"_withratio");
1505 >  CompleteSave(Ratio_main_canvas,savemeas+"_withratio");
1506    bottompad->cd();
1507    
1508    Double_t chi2;
# Line 1363 | Line 1530 | void save_with_ratio(TH1F *nominator, TH
1530    CompleteRatio->SetTextAlign(31);
1531    CompleteRatio->Draw();
1532    
1533 <  CompleteSave(main_canvas,savemeas+"_withratio_and_Chi2");
1533 >  CompleteSave(Ratio_main_canvas,savemeas+"_withratio_and_Chi2");
1534  
1535   //  float KS   = nominator->KolmogorovTest(denominator);
1536   //  stringstream KStext;
1537   //  Chi2text << "KS = " << KS << endl;
1538   //cout << "Found : " << KStext.str() << endl;
1539  
1540 +  delete eratio;
1541 +  delete ratio_center;
1542 +  delete ratio;
1543 +  delete Ratio_main_canvas;
1544 + }
1545  
1546 <  delete main_canvas;
1546 > void Save_With_Ratio_And_Line(TH1F *nominator, TH1F *denominator, TVirtualPad *orig_canvas, string savemeas, bool do_bpred_ratio=false, bool extendrange=false, string yaxistitle="ratio",TH1F *syshisto=NULL) {
1547 >  PlottingSetup::ConsideredZWidth=10.0;
1548 >  PlottingSetup::DrawMetSignalRegionMllLines=true;
1549 >  Save_With_Ratio(nominator, denominator, orig_canvas, savemeas, do_bpred_ratio, extendrange, yaxistitle,syshisto);
1550 >  PlottingSetup::DrawMetSignalRegionMllLines=false;
1551   }
1552 +  
1553 + void Save_With_Ratio_And_Line(TH1F *nominator, THStack denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false) {
1554 +  TH1F *denominator_histo = (TH1F*) CollapseStack(denominator,"TemporaryStack");
1555 +  Save_With_Ratio_And_Line(nominator, denominator_histo, canvas, savemeas, do_bpred_ratio);
1556 +  delete denominator_histo;
1557 + }
1558 +
1559 + void save_with_ratio_and_sys_band(float ConsideredZWidth, TH1F *nominator, TH1F *denominator, TVirtualPad *orig_canvas, string savemeas, bool do_bpred_ratio=false, bool extendrange=false, string yaxistitle="ratio",TH1F *syshisto=NULL) {
1560 +  //this function saves the pad being passed as well as a new one including the SysRatio.
1561 +  orig_canvas->cd();
1562 +  orig_canvas->Update();
1563 +  CompleteSave(orig_canvas,savemeas);
1564 +  
1565 +  TVirtualPad *canvas = (TVirtualPad*) orig_canvas->Clone("TempCanvas");//otherwise the Ratio_main_canvas will own our pad and destroy it upon deletion
1566 +
1567 +  float bottommargin=gStyle->GetPadBottomMargin();
1568 +  float canvas_height=gStyle->GetCanvasDefH();
1569 +  float canvas_width=gStyle->GetCanvasDefW();
1570 +  float ratiospace=0.25;// space the ratio should take up (relative to original pad)
1571 +  
1572 +  float ratiobottommargin=0.3;
1573 +  float ratiotopmargin=0.1;
1574 +  
1575 +  float xstretchfactor=((1-ratiospace)*(1-gStyle->GetPadTopMargin()))/((1)*ratiospace);
1576 +  
1577 +  TCanvas *Ratio_main_canvas = new TCanvas("Ratio_main_canvas","Ratio_main_canvas",(Int_t)canvas_width,(Int_t)(canvas_height*(1+ratiospace)));
1578 +  TPad *mainpad = new TPad("mainpad","mainpad",0,1-(1.0/(1+ratiospace)),1,1);//top (main) pad
1579 +  TPad *coverpad = new TPad("coverpad","coverpad",gStyle->GetPadLeftMargin()-0.008,1-(1.0/(1+ratiospace))-0.04,1,1-(1.0/(1+ratiospace))+0.103);//pad covering up the x scale
1580 +  TPad *bottompad = new TPad("bottompad", "Ratio Pad",0,0,1,(1-(1-bottommargin)/(1+ratiospace))-0.015); //bottom pad
1581 +  
1582 +  Ratio_main_canvas->Range(0,0,1,1);
1583 +  Ratio_main_canvas->SetBorderSize(0);
1584 +  Ratio_main_canvas->SetFrameFillColor(0);
1585 +  
1586 +  mainpad->Draw();
1587 +  mainpad->cd();
1588 +  mainpad->Range(0,0,1,1);
1589 +  mainpad->SetFillColor(kWhite);
1590 +  mainpad->SetBorderSize(0);
1591 +  mainpad->SetFrameFillColor(0);
1592 +  canvas->Range(0,0,1,1);
1593 +  canvas->Draw("same");
1594 +  mainpad->Modified();
1595 +  Ratio_main_canvas->cd();
1596 +  coverpad->Draw();
1597 +  coverpad->cd();
1598 +  coverpad->Range(0,0,1,1);
1599 +  coverpad->SetFillColor(kWhite);
1600 +  coverpad->SetBorderSize(0);
1601 +  coverpad->SetFrameFillColor(0);
1602 +  coverpad->Modified();
1603 +  Ratio_main_canvas->cd();
1604 +  bottompad->SetTopMargin(ratiotopmargin);
1605 +  bottompad->SetBottomMargin(ratiobottommargin);
1606 +  bottompad->Draw();
1607 +  bottompad->cd();
1608 +  bottompad->Range(0,0,1,1);
1609 +  bottompad->SetFillColor(kWhite);
1610 +  TH1F *ratio = (TH1F*)nominator->Clone(GetNumericHistoName().c_str());
1611 +  ratio->Divide(denominator);
1612 +
1613 +
1614 +  TGraphAsymmErrors *eratio;
1615 +  TGraphAsymmErrors *SysEnvelope = new TGraphAsymmErrors();
1616 +  
1617 +  if(syshisto && !do_bpred_ratio) {
1618 +       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1619 +         float dx=syshisto->GetBinWidth(i)/2.0;
1620 +         float dy=syshisto->GetBinContent(i);
1621 +         if(dy<0.01) dy=0.07;
1622 +         SysEnvelope->SetPoint(i-1,syshisto->GetBinCenter(i),1.0);
1623 +         SysEnvelope->SetPointError(i-1,dx,dx,dy,dy);
1624 +       }
1625 +       SysEnvelope->SetFillColor(TColor::GetColor("#006DE1"));
1626 +  }
1627 +
1628 +  ratio->SetTitle("");
1629 +  ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1630 +  if(do_bpred_ratio) ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1631 +  if(extendrange) ratio->GetYaxis()->SetRangeUser(0.0,4.0);
1632 +  ratio->GetXaxis()->SetTitle(nominator->GetXaxis()->GetTitle());
1633 +  ratio->GetXaxis()->CenterTitle();
1634 +  ratio->GetYaxis()->SetTitle(yaxistitle.c_str());
1635 +  ratio->GetYaxis()->SetTitleOffset(0.4);
1636 +  ratio->GetYaxis()->CenterTitle();
1637 +  ratio->SetStats(0);
1638 +  ratio->GetXaxis()->SetLabelSize(xstretchfactor*ratio->GetXaxis()->GetLabelSize());
1639 +  ratio->GetYaxis()->SetLabelSize(xstretchfactor*ratio->GetYaxis()->GetLabelSize());
1640 +  ratio->GetXaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
1641 +  ratio->GetYaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
1642 +  ratio->GetYaxis()->SetNdivisions(502,false);
1643 +  ratio->SetFillColor(TColor::GetColor("#58D3F7"));
1644 +  ratio->Draw("e0");
1645 +  
1646 +  if(syshisto!=0) {
1647 +    SysEnvelope->SetFillColor(TColor::GetColor("#FE9A2E"));
1648 +    SysEnvelope->Draw("2,0,same");
1649 +    ratio->Draw("e0,same");
1650 +  } else {
1651 +    eratio->Draw("0");
1652 +  }
1653 +  ratio->Draw("same,axis");
1654 +  TLine *oneline = new TLine(ratio->GetXaxis()->GetBinLowEdge(1),1,ratio->GetXaxis()->GetBinLowEdge(ratio->GetNbinsX())+ratio->GetXaxis()->GetBinWidth(ratio->GetNbinsX()),1);
1655 +  oneline->SetLineStyle(2);
1656 +  oneline->SetLineColor(kBlue);
1657 +  oneline->Draw("same");
1658 +  
1659 +  
1660 +  
1661 +  if(ConsideredZWidth>0.1) {
1662 +    cout << "Drawing extra lines in ratio this time around ... " << endl;
1663 +    float RatioYMax=2.0;
1664 +    if(extendrange) RatioYMax=4.0;
1665 +    
1666 +    TLine *SRline   = new TLine(70,0,70,RatioYMax);
1667 +    TLine *ZLowLine = new TLine(91.2-ConsideredZWidth,0,91.2-ConsideredZWidth,RatioYMax);
1668 +    TLine *ZHiLine  = new TLine(91.2+ConsideredZWidth,0,91.2+ConsideredZWidth,RatioYMax);
1669 +  
1670 +    SRline->SetLineStyle(2);
1671 +    ZLowLine->SetLineStyle(2);
1672 +    ZHiLine->SetLineStyle(2);
1673 +    SRline->SetLineColor(kGray+2);
1674 +    ZLowLine->SetLineColor(kGray+2);
1675 +    ZHiLine->SetLineColor(kGray+2);
1676 +    SRline->Draw();
1677 +    ZLowLine->Draw();
1678 +    ZHiLine->Draw();
1679 +  }
1680 +
1681 +  Ratio_main_canvas->cd();
1682 +  Ratio_main_canvas->Modified();
1683 +  Ratio_main_canvas->cd();
1684 +  Ratio_main_canvas->SetSelected(Ratio_main_canvas);
1685 +  
1686 +  CompleteSave(Ratio_main_canvas,savemeas+"_withSysRatio");
1687 +  bottompad->cd();
1688 +  
1689 +  Double_t chi2;
1690 +  Int_t ndf,igood;
1691 +  Double_t chi2prob = MarcosChi2TestX(nominator,denominator,chi2,ndf,igood,"UW");
1692  
1693 +  stringstream Chi2text;
1694 +  Chi2text << "#chi^{2} / ndf: " << chi2 << " / " << ndf;
1695 +  stringstream Chi2probtext;
1696 +  Chi2probtext << "#chi^{2} prob: " << chi2prob;
1697 +  TText* chi2box = write_text(0.02,0.11,Chi2text.str());
1698 +  chi2box->SetTextSize(0.08);
1699 +  chi2box->SetTextAlign(11);
1700 +  chi2box->Draw();
1701 +  TText* chi2probbox = write_text(0.02,0.04,Chi2probtext.str());
1702 +  chi2probbox->SetTextSize(0.08);
1703 +  chi2probbox->SetTextAlign(11);
1704 +  chi2probbox->Draw();
1705 +  
1706 +  stringstream CompRatio;
1707 +  CompRatio << "Ratio: " << nominator->Integral() << " / " << denominator->Integral() << " \n " << nominator->Integral()/denominator->Integral();
1708 +  TText *CompleteRatio = write_text(0.98,0.04,CompRatio.str());
1709 +  CompleteRatio->SetTextSize(0.08);
1710 +  CompleteRatio->SetTextAlign(31);
1711 +  CompleteRatio->Draw();
1712 +  
1713 +  CompleteSave(Ratio_main_canvas,savemeas+"_withSysRatio_and_Chi2");
1714 +  delete Ratio_main_canvas;
1715 +  delete ratio;
1716 + }
1717  
1718 < TH1F* CollapseStack(THStack stack) {
1718 > TH1F* CollapseStack(THStack stack,TString hname="CollapsedStack") {
1719     TH1F *bhist = ((TH1F*)((stack.GetHists())->At(0)));
1720 <   TH1F *basehisto = (TH1F*)bhist->Clone("base");
1720 >   TH1F *basehisto = (TH1F*)bhist->Clone(hname);
1721     TIter next(stack.GetHists());
1722     TH1F *h;
1723     int counter=0;
# Line 1389 | Line 1729 | TH1F* CollapseStack(THStack stack) {
1729     return basehisto;
1730   }
1731  
1732 < void save_with_ratio(TH1F *nominator, THStack denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false) {
1733 <  save_with_ratio(nominator, CollapseStack(denominator), canvas, savemeas, do_bpred_ratio);
1732 > TH1F* CollapseStack(THStack *stack,TString hname="CollapsedStack") {
1733 >   TH1F *bhist = ((TH1F*)((stack->GetHists())->At(0)));
1734 >   TH1F *basehisto = (TH1F*)bhist->Clone(hname);
1735 >   TIter next(stack->GetHists());
1736 >   TH1F *h;
1737 >   int counter=0;
1738 >   while ((h=(TH1F*)next())) {
1739 >     counter++;
1740 >     if(counter==1) continue;
1741 >     basehisto->Add(h);
1742 >   }
1743 >   return basehisto;
1744 > }
1745 >
1746 > void Save_With_Ratio(TH1F *nominator, THStack denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false) {
1747 >  TH1F *denominator_histo = (TH1F*) CollapseStack(denominator);
1748 >  Save_With_Ratio(nominator, denominator_histo, canvas, savemeas, do_bpred_ratio);
1749 >  delete denominator_histo;
1750   }
1751  
1752   void flag_this_change(string function, int line, int checked=0) {
# Line 1506 | Line 1862 | void process_directory(TString directory
1862                  perror ("Couldn't open the directory");
1863   }
1864  
1865 + void ClearHisto(TH1* histo) {
1866 +    for(int ix=0;ix<=histo->GetNbinsX()+1;ix++) {//also reset over- and underflow
1867 +        histo->SetBinContent(ix,0);
1868 +        histo->SetBinContent(ix,0);
1869 +    }
1870 + }
1871 +
1872 + void AddSquared(TH1F *histo, TH1F *addthis, float weight) {
1873 +    for(int ix=1;ix<=addthis->GetNbinsX()+1;ix++) {//also include over- and underflow
1874 +        float bincenter=addthis->GetBinCenter(ix);
1875 +        int nBinHisto= histo->FindBin(bincenter);
1876 +        float old_content=histo->GetBinContent(nBinHisto);
1877 +        histo->SetBinContent(nBinHisto,old_content+((addthis->GetBinContent(ix))*(addthis->GetBinContent(ix))*weight));
1878 +    }
1879 + }
1880 +
1881 + void SQRT(TH1* histo) {
1882 +    for(int ix=0;ix<=histo->GetNbinsX()+1;ix++) {//also include over- and underflow
1883 +        histo->SetBinContent(ix,TMath::Sqrt(histo->GetBinContent(ix)));
1884 +    }
1885 + }
1886 +
1887   /*string GetCompleteHostname() {
1888    struct addrinfo hints, *info, *p;
1889    int gai_result;
# Line 1553 | Line 1931 | string removefunnystring(string name) {
1931    name=ReplaceCharacter(name,"]","_");
1932    name=ReplaceCharacter(name,"{","_");
1933    name=ReplaceCharacter(name,"}","_");
1934 +  name=ReplaceCharacter(name,"(","_");
1935 +  name=ReplaceCharacter(name,")","_");
1936 +  name=ReplaceCharacter(name," ","_");
1937    name=ReplaceCharacter(name,".","_");
1938    name=ReplaceCharacter(name,",","_");
1939    name=ReplaceCharacter(name,";","_");
# Line 1567 | Line 1948 | string removefunnystring(string name) {
1948    return name;
1949   }
1950  
1951 + TPolyLine* GetFitUncertaintyShape(TF1 *fit, float low, float high,float x0,float x1) {
1952 +  int nPoints=1000;
1953 +  double x[nPoints+1];
1954 +  double y[nPoints];
1955 +  
1956 +  float par1=fit->GetParameter(0);
1957 +  float dpar1=fit->GetParError(0);
1958 +  float par2=fit->GetParameter(1);
1959 +  float dpar2=fit->GetParError(1);
1960 +    
1961 +  float step=(x1-x0)/(nPoints/2.0 -1);
1962 +  for(int i=0;i<nPoints/2.0;i++) {
1963 +    x[i]=x0+i*step;
1964 +    y[i]=(par1+x[i]*par2)    -    sqrt(dpar1*dpar1+dpar2*dpar2*x[i]*x[i]);
1965 +    x[nPoints-1-i]=x[i];
1966 +    y[nPoints-1-i]=(par1+x[i]*par2)    +    sqrt(dpar1*dpar1+dpar2*dpar2*x[i]*x[i]);
1967 +    if(y[i]<low) y[i]=low;
1968 +    if(y[i]>high) y[i]=high;
1969 +    if(y[nPoints-1-i]<low) y[nPoints-1-i]=low;
1970 +    if(y[nPoints-1-i]>high) y[nPoints-1-i]=high;
1971 +  }
1972 +  
1973 +  x[nPoints]=x[0];
1974 +  y[nPoints]=y[0];
1975 +  
1976 +  TPolyLine *l = new TPolyLine(nPoints+1,x,y);
1977 +  l->SetFillColor(TColor::GetColor("#A2A2FA"));
1978 +  l->SetLineColor(TColor::GetColor("#A2A2FA"));
1979 +  l->SetLineWidth(1);
1980 +  return l;
1981 + }  
1982 +
1983 + //code for ReplaceAll copied from
1984 + //http://stackoverflow.com/questions/5343190/c-how-do-i-replace-all-instances-of-of-a-string-with-another-string
1985 + string ReplaceAll(string str, string from, string to) {
1986 +    size_t start_pos = 0;
1987 +    while((start_pos = str.find(from, start_pos)) != std::string::npos) {
1988 +         str.replace(start_pos, from.length(), to);
1989 +         start_pos += to.length(); // ...
1990 +    }
1991 +    return str;
1992 + }
1993 +
1994 + TPolyLine* GetFitUncertaintyShape(TH1 *histo, string fitname, float low, float high) {
1995 +  TF1 *fit = (TF1*)histo->GetFunction(fitname.c_str());
1996 +  if(!fit) {
1997 +    cout << "CANNOT PRODUCE FIT UNCERTAINTY AS THERE IS NO FIT CALLED " << fitname << " IN " << histo->GetName() << endl;
1998 +    return NULL;
1999 +  }
2000 +  float x0=histo->GetBinLowEdge(1);
2001 +  float x1=histo->GetBinLowEdge(histo->GetNbinsX())+histo->GetBinWidth(histo->GetNbinsX());
2002 +  return GetFitUncertaintyShape(fit, low, high,x0,x1);
2003 + }
2004 +
2005 + TPolyLine* GetFitUncertaintyShape(TGraphErrors *gr, string fitname, float low, float high, float minx=-999, float maxx=-999) {
2006 +  TF1 *fit = (TF1*)gr->GetFunction(fitname.c_str());
2007 +  if(!fit) {
2008 +    cout << "CANNOT PRODUCE FIT UNCERTAINTY AS THERE IS NO FIT CALLED " << fitname << " IN " << gr->GetName() << endl;
2009 +    return NULL;
2010 +  }
2011 +  double x,y;
2012 +  gr->GetPoint(0,x,y);
2013 +  float min=x;
2014 +  gr->GetPoint(gr->GetN()-1,x,y);
2015 +  float max=x;
2016 +  if(minx>-999 && maxx>-999 && minx<maxx) {
2017 +    max=maxx;
2018 +    min=minx;
2019 +  }
2020 +  return GetFitUncertaintyShape(fit, low, high,min,max);
2021 + }
2022 +
2023   stringstream all_bugs;
2024  
2025   void bug_tracker(string function, int line, string description) {
# Line 1604 | Line 2057 | void bug_tracker(string function, int li
2057   //TODO: Write a bug summary at the end.
2058  
2059   float pSTAR(float mglu, float mlsp, float mchi2) {
2060 <        float mz=91.2;
2060 >        float mz=91.0;
2061          float res=((mchi2*mchi2)-(mlsp+mz)*(mlsp+mz));
2062          res*=((mchi2*mchi2)-(mlsp-mz)*(mlsp-mz));
2063          res=TMath::Sqrt(res)/(2*mchi2);
# Line 1617 | Line 2070 | float generalizedpSTAR(float massmother,
2070          res=TMath::Sqrt(res)/(2*massmother);
2071          return res;
2072   }
2073 +
2074 + void FindMinMax(TGraphErrors *gr, float &min, float &max) {
2075 +  double x,y;
2076 +  gr->GetPoint(0,x,y);
2077 +  min=500*y;
2078 +  max=0.0002*y;
2079 +  for(int i=0;i<gr->GetN();i++) {
2080 +    gr->GetPoint(i,x,y);
2081 +    cout << "Point at " << x << "/" << y << " with error " << gr->GetErrorY(i) << endl;
2082 +    float potmin=y-1.5*gr->GetErrorY(i); // taking 1.5 times the error to allow for enough space to present the plot
2083 +    float potmax=y+1.5*gr->GetErrorY(i); // taking 1.5 times the error to allow for enough space to present the plot
2084 +    if(potmin<min) min=potmin;
2085 +    if(potmax>max) max=potmax;
2086 +  }
2087 + }
2088 +
2089 +
2090 + void CleanLegends() {
2091 +  for(int ihist=(int)PlottingSetup::FakeHistoHeap.size()-1;ihist>=0;ihist--) {
2092 +    PlottingSetup::FakeHistoHeap[ihist]->Delete();
2093 +    PlottingSetup::FakeHistoHeap.pop_back();
2094 +  }
2095 + }
2096 +
2097 + string DigitsAfterComma(float number, int digits) {
2098 +  float temp=number*pow(10.0,digits);
2099 +  temp=int(temp)/pow(10.0,digits);
2100 +  return any2string(temp);
2101 + }
2102 +
2103 +
2104 + Value GetYield(TH1F *histo, float min, float max) {
2105 +  float res=0.0;
2106 +  float reserr=0.0;
2107 +  for(int i=1;i<=histo->GetNbinsX();i++) {
2108 +      if(histo->GetBinLowEdge(i)+histo->GetBinWidth(i)<=min) continue;// considered in the next bin
2109 +      if(histo->GetBinLowEdge(i)>=max) continue;//above the threshold, out!
2110 + //      cout << "         added yield for bin [" << histo->GetBinLowEdge(i) << "," << histo->GetBinLowEdge(i)+histo->GetBinWidth(i) << "] which is " << histo->GetBinContent(i) << "   in GetYield" << endl;
2111 +      res+=histo->GetBinContent(i);
2112 +      reserr=sqrt(reserr*reserr+histo->GetBinError(i)*histo->GetBinError(i));
2113 +  }
2114 +  return Value(res,reserr);
2115 + }
2116 +
2117 + void ProduceYields(float min, float max, TH1F *data, THStack *stack) {
2118 +  dout << "   *************** <MC YIELDS> ********* " << endl;
2119 +  dout << "   Considering " << min << " < mll < " << max << "           " << endl;
2120 +  dout << "   Data : " << GetYield(data,min,max)  << endl;
2121 +  TIter nextSF(stack->GetHists());
2122 +  TH1F* h;
2123 +  while ( h = (TH1F*)nextSF() ) {
2124 +    dout << "   " << h->GetName() << " : " << GetYield(h,min,max)  << endl;
2125 +  }
2126 +  dout << "   *************** </MC YIELDS> ********* " << endl;
2127 + }
2128 +
2129 +  
2130 + void WriteYield(THStack *mc, float low, float high) {
2131 +  TH1F *h;
2132 +  TIter NextHisto(mc->GetHists());
2133 +  while ( h = (TH1F*)NextHisto() ) {
2134 +    float CurrBKG=0;
2135 +    float CurrBKGErr=0;
2136 +    for(int i=1;i<=h->GetNbinsX()+1;i++) {
2137 +      if(h->GetBinLowEdge(i)+h->GetBinWidth(i)<=low) continue;// considered in the next bin
2138 +      if(h->GetBinLowEdge(i)>=high) continue;//above the threshold, out!
2139 +      CurrBKG+=h->GetBinContent(i);
2140 +      CurrBKGErr=sqrt(CurrBKGErr*CurrBKGErr+h->GetBinError(i)*h->GetBinError(i));
2141 +    }
2142 +    dout << "      " << h->GetName() << " : " << CurrBKG << " +/- " << CurrBKGErr << endl;
2143 +  }
2144 + }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines