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.8 by buchmann, Mon May 14 15:13:04 2012 UTC vs.
Revision 1.41 by buchmann, Fri Apr 26 07:00:24 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   #ifndef Verbosity
46   #define Verbosity 0
47   #endif
# Line 60 | Line 63 | namespace PlottingSetup {
63      int publicmode=0;
64      int PaperMode=0; // PaperMode=true will suppress "Preliminary" in DrawPrelim()
65      int Approved=0; // Approved=true will only plot approved plots
66 +    bool is2012=true;
67 +    bool is53reco=true;
68 +    bool openBox = true;
69 +    vector<TH1F*> FakeHistoHeap;
70 +    bool DrawMetSignalRegionMllLines=false; // whether to draw the lines in mll plots for the MET signal region
71 +    float ConsideredZWidth=20;
72   }
73  
74   bool dopng=false;
# Line 86 | Line 95 | void write_error(string funcname, string
95   void write_info(string funcname, string text);
96   string get_directory();
97   bool Contains(string wholestring, string findme);
98 + TH1F* CollapseStack(THStack stack,TString hname);
99   //-------------------------------------------------------------------------------------
100  
101   template<typename U>
# Line 253 | Line 263 | void ensure_directory_exists(string this
263    }
264   }
265  
266 + void DeleteStack(THStack mcstack) {
267 +  TH1F* h;
268 +  TIter nextHist(mcstack.GetHists());
269 +  while ( h = (TH1F*)nextHist() ) {
270 +    h->Delete();
271 +  }
272 +  delete h;
273 + }
274 +
275 +
276 + void DeadEnd(string message) {
277 +  dout << endl;
278 +  dout << "\033[1;31m     _                _                  _  " << endl;
279 +  dout << "\033[1;31m  __| | ___  __ _  __| |   ___ _ __   __| | " << 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 << "" << endl;
285 +  dout << " A totally fatal error has occurred: " << endl;
286 +  dout << " " << message << endl;
287 +  dout << endl;
288 +  dout << "If you do not know how to deal with this error message please call Customer Support Mumbai (or Zurich...) \033[0m" << endl;
289 +  dout << "   ... aborting now ... " << endl;
290 +  assert(0);
291 + }
292 +
293 +
294 + void SanityChecks() {
295 +  //this function gets called whenever samples are initialized (see ActiveSamples.C)
296 +  //whatever is necessary for the code to run should be put here.
297 +  
298 +  if(gROOT->GetVersionInt()<53200) {
299 +    DeadEnd("ROOT version outdated, you are using "+any2string(gROOT->GetVersionInt())+" but the earliest recommended version is 5.32.00");
300 +  }
301 +  
302 +  
303 + }
304 +
305   void initialize_log() {
306    if(!PlottingSetup::publicmode) {
307      dout << "____________________________________________________________" << endl;
# Line 392 | Line 441 | void CompleteSave(TVirtualPad *can, stri
441  
442  
443   void setlumi(float l) {
444 < generaltoolboxlumi=l;
444 >  generaltoolboxlumi=l;
445   }
446  
447   int write_first_line(vector<vector<string> > &entries) {
# Line 550 | Line 599 | void DrawPrelim(float writelumi=generalt
599    string prel=" Preliminary";
600    if(PlottingSetup::PaperMode) prel="";
601    //prelimtext << "CMS Preliminary 2011 , #sqrt{s}= 7 TeV, L= O(1) fb^{-1}"; //temporary replacement
602 + //  if(PlottingSetup::is53reco) prel += " 53X";
603 +  string energy="7 TeV";
604 +  if(PlottingSetup::is2012) energy="8 TeV";
605    if(writelumi == 0) {
606 <    if(isMC) prelimtext << "CMS Simulation, #sqrt{s} = 7 TeV";
607 <    else prelimtext << "CMS" << prel << ", #sqrt{s} = 7 TeV";
606 >    if(isMC) prelimtext << "CMS Simulation, #sqrt{s} = " << energy;
607 >    else prelimtext << "CMS" << prel << ", #sqrt{s} = " << energy;
608    } else {
609 <    if(isMC) prelimtext << "CMS Simulation, #sqrt{s} = 7 TeV, L_{int} = " << std::setprecision(3) <<writelumi<<" "<<barn<<"^{-1}";
610 <     else prelimtext << "CMS" << prel << ", #sqrt{s} = 7 TeV, L_{int} = " << std::setprecision(3) <<writelumi<<" "<<barn<<"^{-1}";
609 >    if(isMC) prelimtext << "CMS Simulation, #sqrt{s} = " << energy << ", #scale[0.6]{#int}L dt = " << std::setprecision(3) <<writelumi<<" "<<barn<<"^{-1}";
610 >     else prelimtext << "CMS" << prel << ", #sqrt{s} = " << energy << ", #scale[0.6]{#int}L dt = " << std::setprecision(3) <<writelumi<<" "<<barn<<"^{-1}";
611    }
612    TPaveText *eventSelectionPaveText = new TPaveText(0.27, 0.93,0.77, 1.0,"blNDC");
613    eventSelectionPaveText->SetFillStyle(4000);
# Line 633 | Line 685 | float computeProductError(float a, float
685      return val;
686   }
687  
688 + TGraphAsymmErrors *SysAndStatRatio(TH1F *h1,TH1F *h2, TH1F *sys, int id, vector<float>binning, bool precise=false) {
689 +    int absJZBbinsNumber = binning.size()-1;
690 +    TGraphAsymmErrors* graph = new TGraphAsymmErrors(absJZBbinsNumber);
691 +    for(int i=0;i<absJZBbinsNumber;i++)
692 +        {
693 +                float xCenter=h1->GetBinCenter(i+1);
694 +                float xWidth=(h1->GetBinWidth(i+1))*0.5;
695 +                float nominatorError = h1->GetBinError(i+1);
696 +                float nominator=h1->GetBinContent(i+1);
697 +                float denominatorError=h2->GetBinError(i+1);
698 +                float denominator=h2->GetBinContent(i+1);
699 +                float errorN = 0;
700 +                float errorP = computeRatioError(nominator,nominatorError,denominator,denominatorError);
701 +                float syserr=sys->GetBinContent(i+1);
702 +                float errorsys=0;
703 +                if(id==1) // (is data)
704 +                {
705 +                        if(!precise) errorP = computeRatioError(nominator,statErrorP(nominator),denominator,statErrorP(denominator));
706 +                        else         errorP = computeRatioError(nominator,nominatorError,denominator,denominatorError);
707 +                        errorsys=computeRatioError(nominator,syserr,denominator,0);
708 +                        errorP = TMath::Sqrt(errorsys*errorsys+errorP*errorP);
709 +                        errorN = errorP; // symmetrize using statErrorP
710 +                } else {
711 +                        errorN = computeRatioError(nominator,nominatorError,denominator,denominatorError);
712 +                        errorsys=computeRatioError(nominator,syserr,denominator,0);
713 +                        errorN = TMath::Sqrt(errorsys*errorsys+errorP*errorP);
714 +                        errorP = errorN; // symmetrize using statErrorP
715 +                }
716 +                if(denominator!=0) {
717 +                        graph->SetPoint(i, xCenter, nominator/denominator);
718 +                        graph->SetPointError(i,xWidth,xWidth,errorN,errorP);
719 +                }
720 +                else {
721 +                        graph->SetPoint(i, xCenter, -999);
722 +                        graph->SetPointError(i,xWidth,xWidth,errorN,errorP);
723 +                }
724 +        }      
725 +        return graph;
726 + }
727 +    
728 +
729   TGraphAsymmErrors *histRatio(TH1F *h1,TH1F *h2, int id, vector<float>binning, bool precise=false)
730   {
731          int absJZBbinsNumber = binning.size()-1;
# Line 869 | Line 962 | string interpret_cut(string incut, bool
962    
963      if(incut=="id1==0") return "ee";
964      if(incut=="id1==1") return "#mu#mu";
965 <    if(incut=="abs(mll-91.2)<20") return "|m_{l^{+}l^{-}}-m_{Z}|<20";
965 >    if(incut=="abs(mll-91)<10") return "|m_{l^{+}l^{-}}-m_{Z}|<10";
966      if(incut=="pfJetGoodID[0]") return "";
967      if(incut=="pfJetGoodID[1]") return "";
968      if((int)incut.find("pfJetGoodNum")>-1) {
# Line 1145 | Line 1238 | bool Contains(string wholestring, string
1238    else return false;
1239   }
1240  
1241 +
1242   //////////////////////////////////////////////////////////////////////////////
1243   //
1244   // http://stackoverflow.com/questions/669438/how-to-get-memory-usage-at-run-time-in-c
# Line 1157 | Line 1251 | bool Contains(string wholestring, string
1251   /* usage:
1252   double vm2, rss2;
1253   process_mem_usage(vm2, rss2);
1254 < cout << "Memory usage: VM: " << vm << "; RSS: " << rss << endl;
1254 > cout << "Memory usage: VM: " << vm2 << "; RSS: " << rss2 << endl;
1255   */
1256  
1257   void process_mem_usage(double& vm_usage, double& resident_set)
# Line 1216 | Line 1310 | TGraphAsymmErrors* produce_ratio_graph(T
1310   }
1311  
1312  
1313 + TLatex* WriteSelection(int njets) {
1314 +  string sel="Loose selection";
1315 +  if(njets==3) sel="Low E_{T}^{miss} selection";
1316 +  assert(njets==2||njets==3);
1317 +  TLatex *sele = new TLatex(0.97,0.135,sel.c_str());
1318 +  sele->SetNDC(true);
1319 +  sele->SetTextColor(TColor::GetColor("#848484"));
1320 +  sele->SetTextFont(42);
1321 +  sele->SetTextAlign(32);
1322 +  sele->SetTextSize(0.03);
1323 +  sele->SetTextAngle(270);
1324 +  return sele;
1325 + }
1326 +
1327   Double_t MarcosChi2TestX(const TH1* h1, const TH1* h2, Double_t &chi2, Int_t &ndf, Int_t &igood, Option_t *option)
1328   {
1329  
# Line 1238 | Line 1346 | Double_t MarcosChi2TestX(const TH1* h1,
1346     return prob;
1347   }
1348  
1349 < void save_with_ratio(TH1F *nominator, TH1F *denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false, bool extendrange=false, string yaxistitle="ratio") {
1349 > 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) {
1350    //this function saves the pad being passed as well as a new one including the ratio.
1351 <  CompleteSave(canvas,savemeas);
1352 <
1351 >  
1352 >  orig_canvas->cd();
1353 >  orig_canvas->Update();
1354 >  CompleteSave(orig_canvas,savemeas);
1355 >  TVirtualPad *canvas = (TVirtualPad*) orig_canvas->Clone("TempCanvas");//otherwise the Ratio_main_canvas will own our pad and destroy it upon deletion
1356 >  
1357    float bottommargin=gStyle->GetPadBottomMargin();
1358    float canvas_height=gStyle->GetCanvasDefH();
1359    float canvas_width=gStyle->GetCanvasDefW();
# Line 1252 | Line 1364 | void save_with_ratio(TH1F *nominator, TH
1364    
1365    float xstretchfactor=((1-ratiospace)*(1-gStyle->GetPadTopMargin()))/((1)*ratiospace);
1366    
1367 <  TCanvas *main_canvas = new TCanvas("main_canvas","main_canvas",(Int_t)canvas_width,(Int_t)(canvas_height*(1+ratiospace)));
1367 >  TCanvas *Ratio_main_canvas = new TCanvas("Ratio_main_canvas","Ratio_main_canvas",(Int_t)canvas_width,(Int_t)(canvas_height*(1+ratiospace)));
1368    TPad *mainpad = new TPad("mainpad","mainpad",0,1-(1.0/(1+ratiospace)),1,1);//top (main) pad
1369    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
1370    TPad *bottompad = new TPad("bottompad", "Ratio Pad",0,0,1,(1-(1-bottommargin)/(1+ratiospace))-0.015); //bottom pad
1371    
1372 <  main_canvas->Range(0,0,1,1);
1373 <  main_canvas->SetBorderSize(0);
1374 <  main_canvas->SetFrameFillColor(0);
1263 <  
1372 >  Ratio_main_canvas->Range(0,0,1,1);
1373 >  Ratio_main_canvas->SetBorderSize(0);
1374 >  Ratio_main_canvas->SetFrameFillColor(0);
1375    mainpad->Draw();
1376    mainpad->cd();
1377    mainpad->Range(0,0,1,1);
# Line 1270 | Line 1381 | void save_with_ratio(TH1F *nominator, TH
1381    canvas->Range(0,0,1,1);
1382    canvas->Draw("same");
1383    mainpad->Modified();
1384 <  main_canvas->cd();
1384 >  Ratio_main_canvas->cd();
1385    coverpad->Draw();
1386    coverpad->cd();
1387    coverpad->Range(0,0,1,1);
# Line 1278 | Line 1389 | void save_with_ratio(TH1F *nominator, TH
1389    coverpad->SetBorderSize(0);
1390    coverpad->SetFrameFillColor(0);
1391    coverpad->Modified();
1392 <  main_canvas->cd();
1392 >  Ratio_main_canvas->cd();
1393    bottompad->SetTopMargin(ratiotopmargin);
1394    bottompad->SetBottomMargin(ratiobottommargin);
1395    bottompad->Draw();
# Line 1287 | Line 1398 | void save_with_ratio(TH1F *nominator, TH
1398    bottompad->SetFillColor(kWhite);
1399    TH1F *ratio = (TH1F*)nominator->Clone(GetNumericHistoName().c_str());
1400    ratio->Divide(denominator);
1290
1291
1401    TGraphAsymmErrors *eratio;
1402 +  TH1F *SystDown;
1403 +  TH1F *SystUp;
1404 +  
1405    if(!do_bpred_ratio) eratio = produce_ratio_graph(ratio);
1406    else {
1407       bool using_data=false;
1408       if((int)savemeas.find("Data")>0) using_data=true;
1409       eratio = histRatio(nominator,denominator,using_data,PlottingSetup::global_ratio_binning,false);
1410 +     if(syshisto!=0) {
1411 +       SystDown = (TH1F*)syshisto->Clone("SystDown");
1412 +       SystUp = (TH1F*)syshisto->Clone("SystUp");
1413 +       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1414 +         SystDown->SetBinContent(i,1-syshisto->GetBinContent(i));
1415 +         SystUp->SetBinContent(i,1+syshisto->GetBinContent(i));
1416 +       }
1417 +       SystDown->SetLineColor(TColor::GetColor("#006DE1"));
1418 +       SystUp->SetLineColor(TColor::GetColor("#006DE1"));
1419 +     }
1420       for(int i=1;i<=ratio->GetNbinsX();i++) {
1421          ratio->SetBinContent(i,0);
1422          ratio->SetBinError(i,0);
1423       }
1424 +    
1425 +  }
1426 +  
1427 +  if(syshisto && !do_bpred_ratio) {
1428 +       SystDown = (TH1F*)syshisto->Clone("SystDown");
1429 +       SystUp = (TH1F*)syshisto->Clone("SystUp");
1430 +       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1431 +         SystDown->SetBinContent(i,1-syshisto->GetBinContent(i));
1432 +         SystUp->SetBinContent(i,1+syshisto->GetBinContent(i));
1433 +       }
1434 +       SystDown->SetLineColor(TColor::GetColor("#006DE1"));
1435 +       SystUp->SetLineColor(TColor::GetColor("#006DE1"));
1436    }
1437    eratio->SetFillColor(TColor::GetColor("#00ADE1"));
1304
1438    ratio->SetTitle("");
1439 <  ratio->GetYaxis()->SetRangeUser(0.5,1.5);
1439 >  ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1440    if(do_bpred_ratio) ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1441    if(extendrange) ratio->GetYaxis()->SetRangeUser(0.0,4.0);
1442    ratio->GetXaxis()->SetTitle(nominator->GetXaxis()->GetTitle());
# Line 1319 | Line 1452 | void save_with_ratio(TH1F *nominator, TH
1452    ratio->GetYaxis()->SetNdivisions(502,false);
1453    ratio->SetFillColor(TColor::GetColor("#58D3F7"));
1454    ratio->SetMarkerSize(0);
1455 <  ratio->Draw("e2");
1456 <  eratio->Draw("2");
1455 >  ratio->Draw("e20");
1456 >  
1457 >  TGraphAsymmErrors *ratio_center = (TGraphAsymmErrors*)eratio->Clone("ratio_center");
1458 >  for(int i=0;i<ratio_center->GetN();i++) {
1459 >    ratio_center->SetPointError(i,ratio_center->GetErrorXlow(i),ratio_center->GetErrorXhigh(i),0.005,0.005);
1460 >  }
1461 >  ratio_center->SetFillColor(TColor::GetColor("#006381"));
1462 >  
1463 >  if(syshisto!=0) {
1464 > //    sysratio->Draw("2");
1465 > //    eratio->Draw("2,same");
1466 >    eratio->Draw("20");
1467 >    SystDown->Draw("histo,same");
1468 >    SystUp->Draw("histo,same");
1469 >  } else {
1470 >    eratio->Draw("20");
1471 >  }
1472 >  ratio_center->Draw("20");
1473    ratio->Draw("same,axis");
1474    TLine *oneline = new TLine(ratio->GetXaxis()->GetBinLowEdge(1),1,ratio->GetXaxis()->GetBinLowEdge(ratio->GetNbinsX())+ratio->GetXaxis()->GetBinWidth(ratio->GetNbinsX()),1);
1475    oneline->SetLineStyle(2);
1476    oneline->SetLineColor(kBlue);
1477    oneline->Draw("same");
1478 <
1479 <  main_canvas->cd();
1480 <  main_canvas->Modified();
1481 <  main_canvas->cd();
1482 <  main_canvas->SetSelected(main_canvas);
1478 >  if(PlottingSetup::DrawMetSignalRegionMllLines) {
1479 >    cout << "Drawing extra lines in ratio this time around ... " << endl;
1480 >    float RatioYMax=2.0;
1481 >    if(extendrange) RatioYMax=4.0;
1482 >    
1483 >    TLine *SRline   = new TLine(70,0,70,RatioYMax);
1484 >    TLine *ZLowLine = new TLine(91.2-PlottingSetup::ConsideredZWidth,0,91.2-PlottingSetup::ConsideredZWidth,RatioYMax);
1485 >    TLine *ZHiLine  = new TLine(91.2+PlottingSetup::ConsideredZWidth,0,91.2+PlottingSetup::ConsideredZWidth,RatioYMax);
1486 >  
1487 >    SRline->SetLineStyle(2);
1488 >    ZLowLine->SetLineStyle(2);
1489 >    ZHiLine->SetLineStyle(2);
1490 >    SRline->SetLineColor(kGray+2);
1491 >    ZLowLine->SetLineColor(kGray+2);
1492 >    ZHiLine->SetLineColor(kGray+2);
1493 >    SRline->Draw();
1494 >    ZLowLine->Draw();
1495 >    ZHiLine->Draw();
1496 >  }
1497 >      
1498 >  Ratio_main_canvas->cd();
1499 >  Ratio_main_canvas->Modified();
1500 >  Ratio_main_canvas->cd();
1501 >  Ratio_main_canvas->SetSelected(Ratio_main_canvas);
1502    
1503 <  CompleteSave(main_canvas,savemeas+"_withratio");
1503 >  CompleteSave(Ratio_main_canvas,savemeas+"_withratio");
1504    bottompad->cd();
1505    
1506    Double_t chi2;
# Line 1352 | Line 1520 | void save_with_ratio(TH1F *nominator, TH
1520    chi2probbox->SetTextSize(0.08);
1521    chi2probbox->SetTextAlign(11);
1522    chi2probbox->Draw();
1523 <  CompleteSave(main_canvas,savemeas+"_withratio_and_Chi2");
1523 >  
1524 >  stringstream CompRatio;
1525 >  CompRatio << "Ratio: " << nominator->Integral() << " / " << denominator->Integral() << " \n " << nominator->Integral()/denominator->Integral();
1526 >  TText *CompleteRatio = write_text(0.98,0.04,CompRatio.str());
1527 >  CompleteRatio->SetTextSize(0.08);
1528 >  CompleteRatio->SetTextAlign(31);
1529 >  CompleteRatio->Draw();
1530 >  
1531 >  CompleteSave(Ratio_main_canvas,savemeas+"_withratio_and_Chi2");
1532  
1533   //  float KS   = nominator->KolmogorovTest(denominator);
1534   //  stringstream KStext;
1535   //  Chi2text << "KS = " << KS << endl;
1536   //cout << "Found : " << KStext.str() << endl;
1537  
1538 +  delete eratio;
1539 +  delete ratio_center;
1540 +  delete ratio;
1541 +  delete Ratio_main_canvas;
1542 + }
1543  
1544 <  delete main_canvas;
1544 > 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) {
1545 >  PlottingSetup::ConsideredZWidth=10.0;
1546 >  PlottingSetup::DrawMetSignalRegionMllLines=true;
1547 >  Save_With_Ratio(nominator, denominator, orig_canvas, savemeas, do_bpred_ratio, extendrange, yaxistitle,syshisto);
1548 >  PlottingSetup::DrawMetSignalRegionMllLines=false;
1549 > }
1550 >  
1551 > void Save_With_Ratio_And_Line(TH1F *nominator, THStack denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false) {
1552 >  TH1F *denominator_histo = (TH1F*) CollapseStack(denominator,"TemporaryStack");
1553 >  Save_With_Ratio_And_Line(nominator, denominator_histo, canvas, savemeas, do_bpred_ratio);
1554 >  delete denominator_histo;
1555   }
1556  
1557 + 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) {
1558 +  //this function saves the pad being passed as well as a new one including the SysRatio.
1559 +  orig_canvas->cd();
1560 +  orig_canvas->Update();
1561 +  CompleteSave(orig_canvas,savemeas);
1562 +  
1563 +  TVirtualPad *canvas = (TVirtualPad*) orig_canvas->Clone("TempCanvas");//otherwise the Ratio_main_canvas will own our pad and destroy it upon deletion
1564 +
1565 +  float bottommargin=gStyle->GetPadBottomMargin();
1566 +  float canvas_height=gStyle->GetCanvasDefH();
1567 +  float canvas_width=gStyle->GetCanvasDefW();
1568 +  float ratiospace=0.25;// space the ratio should take up (relative to original pad)
1569 +  
1570 +  float ratiobottommargin=0.3;
1571 +  float ratiotopmargin=0.1;
1572 +  
1573 +  float xstretchfactor=((1-ratiospace)*(1-gStyle->GetPadTopMargin()))/((1)*ratiospace);
1574 +  
1575 +  TCanvas *Ratio_main_canvas = new TCanvas("Ratio_main_canvas","Ratio_main_canvas",(Int_t)canvas_width,(Int_t)(canvas_height*(1+ratiospace)));
1576 +  TPad *mainpad = new TPad("mainpad","mainpad",0,1-(1.0/(1+ratiospace)),1,1);//top (main) pad
1577 +  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
1578 +  TPad *bottompad = new TPad("bottompad", "Ratio Pad",0,0,1,(1-(1-bottommargin)/(1+ratiospace))-0.015); //bottom pad
1579 +  
1580 +  Ratio_main_canvas->Range(0,0,1,1);
1581 +  Ratio_main_canvas->SetBorderSize(0);
1582 +  Ratio_main_canvas->SetFrameFillColor(0);
1583 +  
1584 +  mainpad->Draw();
1585 +  mainpad->cd();
1586 +  mainpad->Range(0,0,1,1);
1587 +  mainpad->SetFillColor(kWhite);
1588 +  mainpad->SetBorderSize(0);
1589 +  mainpad->SetFrameFillColor(0);
1590 +  canvas->Range(0,0,1,1);
1591 +  canvas->Draw("same");
1592 +  mainpad->Modified();
1593 +  Ratio_main_canvas->cd();
1594 +  coverpad->Draw();
1595 +  coverpad->cd();
1596 +  coverpad->Range(0,0,1,1);
1597 +  coverpad->SetFillColor(kWhite);
1598 +  coverpad->SetBorderSize(0);
1599 +  coverpad->SetFrameFillColor(0);
1600 +  coverpad->Modified();
1601 +  Ratio_main_canvas->cd();
1602 +  bottompad->SetTopMargin(ratiotopmargin);
1603 +  bottompad->SetBottomMargin(ratiobottommargin);
1604 +  bottompad->Draw();
1605 +  bottompad->cd();
1606 +  bottompad->Range(0,0,1,1);
1607 +  bottompad->SetFillColor(kWhite);
1608 +  TH1F *ratio = (TH1F*)nominator->Clone(GetNumericHistoName().c_str());
1609 +  ratio->Divide(denominator);
1610 +
1611 +
1612 +  TGraphAsymmErrors *eratio;
1613 +  TGraphAsymmErrors *SysEnvelope = new TGraphAsymmErrors();
1614 +  
1615 +  if(syshisto && !do_bpred_ratio) {
1616 +       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1617 +         float dx=syshisto->GetBinWidth(i)/2.0;
1618 +         float dy=syshisto->GetBinContent(i);
1619 +         if(dy<0.01) dy=0.07;
1620 +         SysEnvelope->SetPoint(i-1,syshisto->GetBinCenter(i),1.0);
1621 +         SysEnvelope->SetPointError(i-1,dx,dx,dy,dy);
1622 +       }
1623 +       SysEnvelope->SetFillColor(TColor::GetColor("#006DE1"));
1624 +  }
1625  
1626 < TH1F* CollapseStack(THStack stack) {
1626 >  ratio->SetTitle("");
1627 >  ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1628 >  if(do_bpred_ratio) ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1629 >  if(extendrange) ratio->GetYaxis()->SetRangeUser(0.0,4.0);
1630 >  ratio->GetXaxis()->SetTitle(nominator->GetXaxis()->GetTitle());
1631 >  ratio->GetXaxis()->CenterTitle();
1632 >  ratio->GetYaxis()->SetTitle(yaxistitle.c_str());
1633 >  ratio->GetYaxis()->SetTitleOffset(0.4);
1634 >  ratio->GetYaxis()->CenterTitle();
1635 >  ratio->SetStats(0);
1636 >  ratio->GetXaxis()->SetLabelSize(xstretchfactor*ratio->GetXaxis()->GetLabelSize());
1637 >  ratio->GetYaxis()->SetLabelSize(xstretchfactor*ratio->GetYaxis()->GetLabelSize());
1638 >  ratio->GetXaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
1639 >  ratio->GetYaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
1640 >  ratio->GetYaxis()->SetNdivisions(502,false);
1641 >  ratio->SetFillColor(TColor::GetColor("#58D3F7"));
1642 >  ratio->Draw("e0");
1643 >  
1644 >  if(syshisto!=0) {
1645 >    SysEnvelope->SetFillColor(TColor::GetColor("#FE9A2E"));
1646 >    SysEnvelope->Draw("2,0,same");
1647 >    ratio->Draw("e0,same");
1648 >  } else {
1649 >    eratio->Draw("0");
1650 >  }
1651 >  ratio->Draw("same,axis");
1652 >  TLine *oneline = new TLine(ratio->GetXaxis()->GetBinLowEdge(1),1,ratio->GetXaxis()->GetBinLowEdge(ratio->GetNbinsX())+ratio->GetXaxis()->GetBinWidth(ratio->GetNbinsX()),1);
1653 >  oneline->SetLineStyle(2);
1654 >  oneline->SetLineColor(kBlue);
1655 >  oneline->Draw("same");
1656 >  
1657 >  
1658 >  
1659 >  if(ConsideredZWidth>0.1) {
1660 >    cout << "Drawing extra lines in ratio this time around ... " << endl;
1661 >    float RatioYMax=2.0;
1662 >    if(extendrange) RatioYMax=4.0;
1663 >    
1664 >    TLine *SRline   = new TLine(70,0,70,RatioYMax);
1665 >    TLine *ZLowLine = new TLine(91.2-ConsideredZWidth,0,91.2-ConsideredZWidth,RatioYMax);
1666 >    TLine *ZHiLine  = new TLine(91.2+ConsideredZWidth,0,91.2+ConsideredZWidth,RatioYMax);
1667 >  
1668 >    SRline->SetLineStyle(2);
1669 >    ZLowLine->SetLineStyle(2);
1670 >    ZHiLine->SetLineStyle(2);
1671 >    SRline->SetLineColor(kGray+2);
1672 >    ZLowLine->SetLineColor(kGray+2);
1673 >    ZHiLine->SetLineColor(kGray+2);
1674 >    SRline->Draw();
1675 >    ZLowLine->Draw();
1676 >    ZHiLine->Draw();
1677 >  }
1678 >
1679 >  Ratio_main_canvas->cd();
1680 >  Ratio_main_canvas->Modified();
1681 >  Ratio_main_canvas->cd();
1682 >  Ratio_main_canvas->SetSelected(Ratio_main_canvas);
1683 >  
1684 >  CompleteSave(Ratio_main_canvas,savemeas+"_withSysRatio");
1685 >  bottompad->cd();
1686 >  
1687 >  Double_t chi2;
1688 >  Int_t ndf,igood;
1689 >  Double_t chi2prob = MarcosChi2TestX(nominator,denominator,chi2,ndf,igood,"UW");
1690 >
1691 >  stringstream Chi2text;
1692 >  Chi2text << "#chi^{2} / ndf: " << chi2 << " / " << ndf;
1693 >  stringstream Chi2probtext;
1694 >  Chi2probtext << "#chi^{2} prob: " << chi2prob;
1695 >  TText* chi2box = write_text(0.02,0.11,Chi2text.str());
1696 >  chi2box->SetTextSize(0.08);
1697 >  chi2box->SetTextAlign(11);
1698 >  chi2box->Draw();
1699 >  TText* chi2probbox = write_text(0.02,0.04,Chi2probtext.str());
1700 >  chi2probbox->SetTextSize(0.08);
1701 >  chi2probbox->SetTextAlign(11);
1702 >  chi2probbox->Draw();
1703 >  
1704 >  stringstream CompRatio;
1705 >  CompRatio << "Ratio: " << nominator->Integral() << " / " << denominator->Integral() << " \n " << nominator->Integral()/denominator->Integral();
1706 >  TText *CompleteRatio = write_text(0.98,0.04,CompRatio.str());
1707 >  CompleteRatio->SetTextSize(0.08);
1708 >  CompleteRatio->SetTextAlign(31);
1709 >  CompleteRatio->Draw();
1710 >  
1711 >  CompleteSave(Ratio_main_canvas,savemeas+"_withSysRatio_and_Chi2");
1712 >  delete Ratio_main_canvas;
1713 >  delete ratio;
1714 > }
1715 >
1716 > TH1F* CollapseStack(THStack stack,TString hname="CollapsedStack") {
1717     TH1F *bhist = ((TH1F*)((stack.GetHists())->At(0)));
1718 <   TH1F *basehisto = (TH1F*)bhist->Clone("base");
1718 >   TH1F *basehisto = (TH1F*)bhist->Clone(hname);
1719     TIter next(stack.GetHists());
1720     TH1F *h;
1721     int counter=0;
# Line 1378 | Line 1727 | TH1F* CollapseStack(THStack stack) {
1727     return basehisto;
1728   }
1729  
1730 < void save_with_ratio(TH1F *nominator, THStack denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false) {
1731 <  save_with_ratio(nominator, CollapseStack(denominator), canvas, savemeas, do_bpred_ratio);
1730 > void Save_With_Ratio(TH1F *nominator, THStack denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false) {
1731 >  TH1F *denominator_histo = (TH1F*) CollapseStack(denominator);
1732 >  Save_With_Ratio(nominator, denominator_histo, canvas, savemeas, do_bpred_ratio);
1733 >  delete denominator_histo;
1734   }
1735  
1736   void flag_this_change(string function, int line, int checked=0) {
# Line 1437 | Line 1788 | void write_analysis_type(bool isonpeak,
1788      //we're in public! don't advertise what we're up to :-)
1789         dout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << endl;
1790         dout << "                              Starting the analysis                     " << endl;
1791 <       dout << "                                     i:" << isonpeak << " , b " << dobtag << endl;
1791 >       dout << "                                     i:" << !isonpeak << " , b " << dobtag << endl;
1792         dout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << endl;
1793    }
1794    dout << "\033[0m" << endl;
# Line 1495 | Line 1846 | void process_directory(TString directory
1846                  perror ("Couldn't open the directory");
1847   }
1848  
1849 + void ClearHisto(TH1* histo) {
1850 +    for(int ix=0;ix<=histo->GetNbinsX()+1;ix++) {//also reset over- and underflow
1851 +        histo->SetBinContent(ix,0);
1852 +        histo->SetBinContent(ix,0);
1853 +    }
1854 + }
1855 +
1856 + void AddSquared(TH1F *histo, TH1F *addthis, float weight) {
1857 +    for(int ix=1;ix<=addthis->GetNbinsX()+1;ix++) {//also include over- and underflow
1858 +        float bincenter=addthis->GetBinCenter(ix);
1859 +        int nBinHisto= histo->FindBin(bincenter);
1860 +        float old_content=histo->GetBinContent(nBinHisto);
1861 +        histo->SetBinContent(nBinHisto,old_content+((addthis->GetBinContent(ix))*(addthis->GetBinContent(ix))*weight));
1862 +    }
1863 + }
1864 +
1865 + void SQRT(TH1* histo) {
1866 +    for(int ix=0;ix<=histo->GetNbinsX()+1;ix++) {//also include over- and underflow
1867 +        histo->SetBinContent(ix,TMath::Sqrt(histo->GetBinContent(ix)));
1868 +    }
1869 + }
1870 +
1871   /*string GetCompleteHostname() {
1872    struct addrinfo hints, *info, *p;
1873    int gai_result;
# Line 1542 | Line 1915 | string removefunnystring(string name) {
1915    name=ReplaceCharacter(name,"]","_");
1916    name=ReplaceCharacter(name,"{","_");
1917    name=ReplaceCharacter(name,"}","_");
1918 +  name=ReplaceCharacter(name,"(","_");
1919 +  name=ReplaceCharacter(name,")","_");
1920 +  name=ReplaceCharacter(name," ","_");
1921    name=ReplaceCharacter(name,".","_");
1922    name=ReplaceCharacter(name,",","_");
1923    name=ReplaceCharacter(name,";","_");
# Line 1556 | Line 1932 | string removefunnystring(string name) {
1932    return name;
1933   }
1934  
1935 + TPolyLine* GetFitUncertaintyShape(TF1 *fit, float low, float high,float x0,float x1) {
1936 +  int nPoints=1000;
1937 +  double x[nPoints+1];
1938 +  double y[nPoints];
1939 +  
1940 +  float par1=fit->GetParameter(0);
1941 +  float dpar1=fit->GetParError(0);
1942 +  float par2=fit->GetParameter(1);
1943 +  float dpar2=fit->GetParError(1);
1944 +    
1945 +  float step=(x1-x0)/(nPoints/2.0 -1);
1946 +  for(int i=0;i<nPoints/2.0;i++) {
1947 +    x[i]=x0+i*step;
1948 +    y[i]=(par1+x[i]*par2)    -    sqrt(dpar1*dpar1+dpar2*dpar2*x[i]*x[i]);
1949 +    x[nPoints-1-i]=x[i];
1950 +    y[nPoints-1-i]=(par1+x[i]*par2)    +    sqrt(dpar1*dpar1+dpar2*dpar2*x[i]*x[i]);
1951 +    if(y[i]<low) y[i]=low;
1952 +    if(y[i]>high) y[i]=high;
1953 +    if(y[nPoints-1-i]<low) y[nPoints-1-i]=low;
1954 +    if(y[nPoints-1-i]>high) y[nPoints-1-i]=high;
1955 +  }
1956 +  
1957 +  x[nPoints]=x[0];
1958 +  y[nPoints]=y[0];
1959 +  
1960 +  TPolyLine *l = new TPolyLine(nPoints+1,x,y);
1961 +  l->SetFillColor(TColor::GetColor("#A2A2FA"));
1962 +  l->SetLineColor(TColor::GetColor("#A2A2FA"));
1963 +  l->SetLineWidth(1);
1964 +  return l;
1965 + }  
1966 +
1967 + //code for ReplaceAll copied from
1968 + //http://stackoverflow.com/questions/5343190/c-how-do-i-replace-all-instances-of-of-a-string-with-another-string
1969 + string ReplaceAll(string str, string from, string to) {
1970 +    size_t start_pos = 0;
1971 +    while((start_pos = str.find(from, start_pos)) != std::string::npos) {
1972 +         str.replace(start_pos, from.length(), to);
1973 +         start_pos += to.length(); // ...
1974 +    }
1975 +    return str;
1976 + }
1977 +
1978 + TPolyLine* GetFitUncertaintyShape(TH1 *histo, string fitname, float low, float high) {
1979 +  TF1 *fit = (TF1*)histo->GetFunction(fitname.c_str());
1980 +  if(!fit) {
1981 +    cout << "CANNOT PRODUCE FIT UNCERTAINTY AS THERE IS NO FIT CALLED " << fitname << " IN " << histo->GetName() << endl;
1982 +    return NULL;
1983 +  }
1984 +  float x0=histo->GetBinLowEdge(1);
1985 +  float x1=histo->GetBinLowEdge(histo->GetNbinsX())+histo->GetBinWidth(histo->GetNbinsX());
1986 +  return GetFitUncertaintyShape(fit, low, high,x0,x1);
1987 + }
1988 +
1989 + TPolyLine* GetFitUncertaintyShape(TGraphErrors *gr, string fitname, float low, float high, float minx=-999, float maxx=-999) {
1990 +  TF1 *fit = (TF1*)gr->GetFunction(fitname.c_str());
1991 +  if(!fit) {
1992 +    cout << "CANNOT PRODUCE FIT UNCERTAINTY AS THERE IS NO FIT CALLED " << fitname << " IN " << gr->GetName() << endl;
1993 +    return NULL;
1994 +  }
1995 +  double x,y;
1996 +  gr->GetPoint(0,x,y);
1997 +  float min=x;
1998 +  gr->GetPoint(gr->GetN()-1,x,y);
1999 +  float max=x;
2000 +  if(minx>-999 && maxx>-999 && minx<maxx) {
2001 +    max=maxx;
2002 +    min=minx;
2003 +  }
2004 +  return GetFitUncertaintyShape(fit, low, high,min,max);
2005 + }
2006 +
2007   stringstream all_bugs;
2008  
2009   void bug_tracker(string function, int line, string description) {
# Line 1593 | Line 2041 | void bug_tracker(string function, int li
2041   //TODO: Write a bug summary at the end.
2042  
2043   float pSTAR(float mglu, float mlsp, float mchi2) {
2044 <        float mz=91.2;
2044 >        float mz=91.0;
2045          float res=((mchi2*mchi2)-(mlsp+mz)*(mlsp+mz));
2046          res*=((mchi2*mchi2)-(mlsp-mz)*(mlsp-mz));
2047          res=TMath::Sqrt(res)/(2*mchi2);
# Line 1606 | Line 2054 | float generalizedpSTAR(float massmother,
2054          res=TMath::Sqrt(res)/(2*massmother);
2055          return res;
2056   }
2057 +
2058 + void FindMinMax(TGraphErrors *gr, float &min, float &max) {
2059 +  double x,y;
2060 +  gr->GetPoint(0,x,y);
2061 +  min=500*y;
2062 +  max=0.0002*y;
2063 +  for(int i=0;i<gr->GetN();i++) {
2064 +    gr->GetPoint(i,x,y);
2065 +    cout << "Point at " << x << "/" << y << " with error " << gr->GetErrorY(i) << endl;
2066 +    float potmin=y-1.5*gr->GetErrorY(i); // taking 1.5 times the error to allow for enough space to present the plot
2067 +    float potmax=y+1.5*gr->GetErrorY(i); // taking 1.5 times the error to allow for enough space to present the plot
2068 +    if(potmin<min) min=potmin;
2069 +    if(potmax>max) max=potmax;
2070 +  }
2071 + }
2072 +
2073 +
2074 + void CleanLegends() {
2075 +  for(int ihist=(int)PlottingSetup::FakeHistoHeap.size()-1;ihist>=0;ihist--) {
2076 +    PlottingSetup::FakeHistoHeap[ihist]->Delete();
2077 +    PlottingSetup::FakeHistoHeap.pop_back();
2078 +  }
2079 + }
2080 +
2081 + string DigitsAfterComma(float number, int digits) {
2082 +  float temp=number*pow(10.0,digits);
2083 +  temp=int(temp)/pow(10.0,digits);
2084 +  return any2string(temp);
2085 + }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines