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.17 by fronga, Fri Aug 31 07:59:34 2012 UTC vs.
Revision 1.28 by buchmann, Wed Jan 16 17:11:23 2013 UTC

# Line 40 | Line 40
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 62 | Line 64 | namespace PlottingSetup {
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=false;
67 >    bool is53reco=true;
68 >    bool openBox = true;
69 >    vector<TH1F*> FakeHistoHeap;
70 >
71   }
72  
73   bool dopng=false;
# Line 256 | Line 261 | void ensure_directory_exists(string this
261    }
262   }
263  
264 + void DeleteStack(THStack mcstack) {
265 +  TH1F* h;
266 +  TIter nextHist(mcstack.GetHists());
267 +  while ( h = (TH1F*)nextHist() ) {
268 +    h->Delete();
269 +  }
270 +  delete h;
271 + }
272 +
273 +
274   void DeadEnd(string message) {
275    dout << endl;
276    dout << "\033[1;31m     _                _                  _  " << endl;
# Line 424 | Line 439 | void CompleteSave(TVirtualPad *can, stri
439  
440  
441   void setlumi(float l) {
442 < generaltoolboxlumi=l;
442 >  generaltoolboxlumi=l;
443   }
444  
445   int write_first_line(vector<vector<string> > &entries) {
# Line 945 | Line 960 | string interpret_cut(string incut, bool
960    
961      if(incut=="id1==0") return "ee";
962      if(incut=="id1==1") return "#mu#mu";
963 <    if(incut=="abs(mll-91.2)<20") return "|m_{l^{+}l^{-}}-m_{Z}|<20";
963 >    if(incut=="abs(mll-91)<10") return "|m_{l^{+}l^{-}}-m_{Z}|<10";
964      if(incut=="pfJetGoodID[0]") return "";
965      if(incut=="pfJetGoodID[1]") return "";
966      if((int)incut.find("pfJetGoodNum")>-1) {
# Line 1314 | Line 1329 | Double_t MarcosChi2TestX(const TH1* h1,
1329     return prob;
1330   }
1331  
1317 /*void save_with_difference(TH1F *obs, vector<TH1F*> predictions, vector<float> prediction_weights, vector<float> prediction_uncerts, TVirtualPad *canvas, string savemeas) {
1318  if(predictions.size()<1) {
1319    write_warning(__FUNCTION__,"Cannot make a pred-obs prediction plot because the prediction vector is empty. aborting!");
1320    return;
1321  }
1322  TH1F *comppred = predictions[0]->Clone("comppred");
1323  comppred->Scale(prediction_weights[0]);
1324  for(int i=1;i<predictions.size();i++) comppred->Add(predictions[i],prediction_weights[i]);
1325  
1326  TH1F *statpred = (TH1F*)comppred->Clone("statpred");
1327  for(int i=1;i<=statpred->GetNbinsX();i++) statpred->SetBinContent(0);
1328  
1329  TH1F *syspred = predictions[0];
1330  for(int i=1;i<=syspred->GetNbinsX();i++) syspred->SetBinContent(i,prediction_weights[0]*prediction_weights[0]*prediction_uncerts[0]*prediction_uncerts[0]*syspred->GetBinContent(i)*syspred->GetBinContent(i));
1331  for(int i=1;i<predictions.size();i++) {
1332    for(int ibin=1;ibin<=predictions[i]->GetNbinsX();ibin++) {
1333      float newentry=syspred->GetBinContent(ibin);
1334      newentry+=(predictions[i]->GetBinContent(ibin))*(predictions[i]->GetBinContent(ibin))*prediction_weights[i]*prediction_weights[i]*prediction_uncerts[i]*prediction_uncerts[i];
1335      syspred->SetBinContent(ibin,newentry);
1336    }
1337  }
1338  for(int ibin=1;ibin<=syspred->GetNbinsX();ibin++) {
1339    syspred->SetBinError(ibin,syspred->GetBinContent(i)+statpred->GetBinError());
1340    syspred->SetBinContent(ibin,0);
1341  }
1342  
1343  
1344 }
1345 */
1332   void save_with_ratio(TH1F *nominator, TH1F *denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false, bool extendrange=false, string yaxistitle="ratio",TH1F *syshisto=NULL) {
1333    //this function saves the pad being passed as well as a new one including the ratio.
1334    CompleteSave(canvas,savemeas);
# Line 1393 | Line 1379 | void save_with_ratio(TH1F *nominator, TH
1379    TH1F *ratio = (TH1F*)nominator->Clone(GetNumericHistoName().c_str());
1380    ratio->Divide(denominator);
1381  
1396
1382    TGraphAsymmErrors *eratio;
1383    TH1F *SystDown;
1384    TH1F *SystUp;
# Line 1419 | Line 1404 | void save_with_ratio(TH1F *nominator, TH
1404       }
1405      
1406    }
1422  eratio->SetFillColor(TColor::GetColor("#00ADE1"));
1407    
1408 +  if(syshisto && !do_bpred_ratio) {
1409 +       SystDown = (TH1F*)syshisto->Clone("SystDown");
1410 +       SystUp = (TH1F*)syshisto->Clone("SystUp");
1411 +       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1412 +         SystDown->SetBinContent(i,1-syshisto->GetBinContent(i));
1413 +         SystUp->SetBinContent(i,1+syshisto->GetBinContent(i));
1414 +       }
1415 +       SystDown->SetLineColor(TColor::GetColor("#006DE1"));
1416 +       SystUp->SetLineColor(TColor::GetColor("#006DE1"));
1417 +  }
1418 +  eratio->SetFillColor(TColor::GetColor("#00ADE1"));
1419  
1420    ratio->SetTitle("");
1421 <  ratio->GetYaxis()->SetRangeUser(0.5,1.5);
1421 >  ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1422    if(do_bpred_ratio) ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1423    if(extendrange) ratio->GetYaxis()->SetRangeUser(0.0,4.0);
1424    ratio->GetXaxis()->SetTitle(nominator->GetXaxis()->GetTitle());
# Line 1440 | Line 1435 | void save_with_ratio(TH1F *nominator, TH
1435    ratio->SetFillColor(TColor::GetColor("#58D3F7"));
1436    ratio->SetMarkerSize(0);
1437    ratio->Draw("e2");
1438 +  
1439 +  TGraphAsymmErrors *ratio_center = (TGraphAsymmErrors*)eratio->Clone("ratio_center");
1440 +  for(int i=0;i<ratio_center->GetN();i++) {
1441 +    ratio_center->SetPointError(i,ratio_center->GetErrorXlow(i),ratio_center->GetErrorXhigh(i),0.005,0.005);
1442 +  }
1443 +  
1444 +  ratio_center->SetFillColor(TColor::GetColor("#006381"));
1445 +  
1446    if(syshisto!=0) {
1447   //    sysratio->Draw("2");
1448   //    eratio->Draw("2,same");
# Line 1449 | Line 1452 | void save_with_ratio(TH1F *nominator, TH
1452    } else {
1453      eratio->Draw("2");
1454    }
1455 +  ratio_center->Draw("2");
1456    ratio->Draw("same,axis");
1457    TLine *oneline = new TLine(ratio->GetXaxis()->GetBinLowEdge(1),1,ratio->GetXaxis()->GetBinLowEdge(ratio->GetNbinsX())+ratio->GetXaxis()->GetBinWidth(ratio->GetNbinsX()),1);
1458    oneline->SetLineStyle(2);
# Line 1495 | Line 1499 | void save_with_ratio(TH1F *nominator, TH
1499   //  Chi2text << "KS = " << KS << endl;
1500   //cout << "Found : " << KStext.str() << endl;
1501  
1498
1502    delete main_canvas;
1503   }
1504  
1505 + void save_with_ratio_and_sys_band(TH1F *nominator, TH1F *denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false, bool extendrange=false, string yaxistitle="ratio",TH1F *syshisto=NULL) {
1506 +  //this function saves the pad being passed as well as a new one including the SysRatio.
1507 +  CompleteSave(canvas,savemeas);
1508 +
1509 +  float bottommargin=gStyle->GetPadBottomMargin();
1510 +  float canvas_height=gStyle->GetCanvasDefH();
1511 +  float canvas_width=gStyle->GetCanvasDefW();
1512 +  float ratiospace=0.25;// space the ratio should take up (relative to original pad)
1513 +  
1514 +  float ratiobottommargin=0.3;
1515 +  float ratiotopmargin=0.1;
1516 +  
1517 +  float xstretchfactor=((1-ratiospace)*(1-gStyle->GetPadTopMargin()))/((1)*ratiospace);
1518 +  
1519 +  TCanvas *main_canvas = new TCanvas("main_canvas","main_canvas",(Int_t)canvas_width,(Int_t)(canvas_height*(1+ratiospace)));
1520 +  TPad *mainpad = new TPad("mainpad","mainpad",0,1-(1.0/(1+ratiospace)),1,1);//top (main) pad
1521 +  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
1522 +  TPad *bottompad = new TPad("bottompad", "Ratio Pad",0,0,1,(1-(1-bottommargin)/(1+ratiospace))-0.015); //bottom pad
1523 +  
1524 +  main_canvas->Range(0,0,1,1);
1525 +  main_canvas->SetBorderSize(0);
1526 +  main_canvas->SetFrameFillColor(0);
1527 +  
1528 +  mainpad->Draw();
1529 +  mainpad->cd();
1530 +  mainpad->Range(0,0,1,1);
1531 +  mainpad->SetFillColor(kWhite);
1532 +  mainpad->SetBorderSize(0);
1533 +  mainpad->SetFrameFillColor(0);
1534 +  canvas->Range(0,0,1,1);
1535 +  canvas->Draw("same");
1536 +  mainpad->Modified();
1537 +  main_canvas->cd();
1538 +  coverpad->Draw();
1539 +  coverpad->cd();
1540 +  coverpad->Range(0,0,1,1);
1541 +  coverpad->SetFillColor(kWhite);
1542 +  coverpad->SetBorderSize(0);
1543 +  coverpad->SetFrameFillColor(0);
1544 +  coverpad->Modified();
1545 +  main_canvas->cd();
1546 +  bottompad->SetTopMargin(ratiotopmargin);
1547 +  bottompad->SetBottomMargin(ratiobottommargin);
1548 +  bottompad->Draw();
1549 +  bottompad->cd();
1550 +  bottompad->Range(0,0,1,1);
1551 +  bottompad->SetFillColor(kWhite);
1552 +  TH1F *ratio = (TH1F*)nominator->Clone(GetNumericHistoName().c_str());
1553 +  ratio->Divide(denominator);
1554 +
1555  
1556 < TH1F* CollapseStack(THStack stack) {
1556 >  TGraphAsymmErrors *eratio;
1557 >  TGraphAsymmErrors *SysEnvelope = new TGraphAsymmErrors();
1558 >  
1559 >  if(syshisto && !do_bpred_ratio) {
1560 >       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1561 >         float dx=syshisto->GetBinWidth(i)/2.0;
1562 >         float dy=syshisto->GetBinContent(i);
1563 >         if(dy<0.01) dy=0.07;
1564 >         SysEnvelope->SetPoint(i-1,syshisto->GetBinCenter(i),1.0);
1565 >         SysEnvelope->SetPointError(i-1,dx,dx,dy,dy);
1566 >       }
1567 >       SysEnvelope->SetFillColor(TColor::GetColor("#006DE1"));
1568 >  }
1569 >
1570 >  ratio->SetTitle("");
1571 >  ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1572 >  if(do_bpred_ratio) ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1573 >  if(extendrange) ratio->GetYaxis()->SetRangeUser(0.0,4.0);
1574 >  ratio->GetXaxis()->SetTitle(nominator->GetXaxis()->GetTitle());
1575 >  ratio->GetXaxis()->CenterTitle();
1576 >  ratio->GetYaxis()->SetTitle(yaxistitle.c_str());
1577 >  ratio->GetYaxis()->SetTitleOffset(0.4);
1578 >  ratio->GetYaxis()->CenterTitle();
1579 >  ratio->SetStats(0);
1580 >  ratio->GetXaxis()->SetLabelSize(xstretchfactor*ratio->GetXaxis()->GetLabelSize());
1581 >  ratio->GetYaxis()->SetLabelSize(xstretchfactor*ratio->GetYaxis()->GetLabelSize());
1582 >  ratio->GetXaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
1583 >  ratio->GetYaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
1584 >  ratio->GetYaxis()->SetNdivisions(502,false);
1585 >  ratio->SetFillColor(TColor::GetColor("#58D3F7"));
1586 >  ratio->Draw("e1");
1587 >  
1588 >  if(syshisto!=0) {
1589 >    SysEnvelope->SetFillColor(TColor::GetColor("#FE9A2E"));
1590 >    SysEnvelope->Draw("2,same");
1591 >    ratio->Draw("e1,same");
1592 >  } else {
1593 >    eratio->Draw("1");
1594 >  }
1595 >  ratio->Draw("same,axis");
1596 >  TLine *oneline = new TLine(ratio->GetXaxis()->GetBinLowEdge(1),1,ratio->GetXaxis()->GetBinLowEdge(ratio->GetNbinsX())+ratio->GetXaxis()->GetBinWidth(ratio->GetNbinsX()),1);
1597 >  oneline->SetLineStyle(2);
1598 >  oneline->SetLineColor(kBlue);
1599 >  oneline->Draw("same");
1600 >
1601 >  main_canvas->cd();
1602 >  main_canvas->Modified();
1603 >  main_canvas->cd();
1604 >  main_canvas->SetSelected(main_canvas);
1605 >  
1606 >  CompleteSave(main_canvas,savemeas+"_withSysRatio");
1607 >  bottompad->cd();
1608 >  
1609 >  Double_t chi2;
1610 >  Int_t ndf,igood;
1611 >  Double_t chi2prob = MarcosChi2TestX(nominator,denominator,chi2,ndf,igood,"UW");
1612 >
1613 >  stringstream Chi2text;
1614 >  Chi2text << "#chi^{2} / ndf: " << chi2 << " / " << ndf;
1615 >  stringstream Chi2probtext;
1616 >  Chi2probtext << "#chi^{2} prob: " << chi2prob;
1617 >  TText* chi2box = write_text(0.02,0.11,Chi2text.str());
1618 >  chi2box->SetTextSize(0.08);
1619 >  chi2box->SetTextAlign(11);
1620 >  chi2box->Draw();
1621 >  TText* chi2probbox = write_text(0.02,0.04,Chi2probtext.str());
1622 >  chi2probbox->SetTextSize(0.08);
1623 >  chi2probbox->SetTextAlign(11);
1624 >  chi2probbox->Draw();
1625 >  
1626 >  stringstream CompRatio;
1627 >  CompRatio << "Ratio: " << nominator->Integral() << " / " << denominator->Integral() << " \n " << nominator->Integral()/denominator->Integral();
1628 >  TText *CompleteRatio = write_text(0.98,0.04,CompRatio.str());
1629 >  CompleteRatio->SetTextSize(0.08);
1630 >  CompleteRatio->SetTextAlign(31);
1631 >  CompleteRatio->Draw();
1632 >  
1633 >  CompleteSave(main_canvas,savemeas+"_withSysRatio_and_Chi2");
1634 >  delete main_canvas;
1635 > }
1636 >
1637 > TH1F* CollapseStack(THStack stack,TString hname="base") {
1638     TH1F *bhist = ((TH1F*)((stack.GetHists())->At(0)));
1639 <   TH1F *basehisto = (TH1F*)bhist->Clone("base");
1639 >   TH1F *basehisto = (TH1F*)bhist->Clone(hname);
1640     TIter next(stack.GetHists());
1641     TH1F *h;
1642     int counter=0;
# Line 1700 | Line 1834 | string removefunnystring(string name) {
1834    name=ReplaceCharacter(name,"]","_");
1835    name=ReplaceCharacter(name,"{","_");
1836    name=ReplaceCharacter(name,"}","_");
1837 +  name=ReplaceCharacter(name,"(","_");
1838 +  name=ReplaceCharacter(name,")","_");
1839 +  name=ReplaceCharacter(name," ","_");
1840    name=ReplaceCharacter(name,".","_");
1841    name=ReplaceCharacter(name,",","_");
1842    name=ReplaceCharacter(name,";","_");
# Line 1714 | Line 1851 | string removefunnystring(string name) {
1851    return name;
1852   }
1853  
1854 + TPolyLine* GetFitUncertaintyShape(TF1 *fit, float low, float high,float x0,float x1) {
1855 +  int nPoints=1000;
1856 +  double x[nPoints+1];
1857 +  double y[nPoints];
1858 +  
1859 +  float par1=fit->GetParameter(0);
1860 +  float dpar1=fit->GetParError(0);
1861 +  float par2=fit->GetParameter(1);
1862 +  float dpar2=fit->GetParError(1);
1863 +    
1864 +  float step=(x1-x0)/(nPoints/2.0 -1);
1865 +  for(int i=0;i<nPoints/2.0;i++) {
1866 +    x[i]=x0+i*step;
1867 +    y[i]=(par1+x[i]*par2)    -    sqrt(dpar1*dpar1+dpar2*dpar2*x[i]*x[i]);
1868 +    x[nPoints-1-i]=x[i];
1869 +    y[nPoints-1-i]=(par1+x[i]*par2)    +    sqrt(dpar1*dpar1+dpar2*dpar2*x[i]*x[i]);
1870 +    if(y[i]<low) y[i]=low;
1871 +    if(y[i]>high) y[i]=high;
1872 +    if(y[nPoints-1-i]<low) y[nPoints-1-i]=low;
1873 +    if(y[nPoints-1-i]>high) y[nPoints-1-i]=high;
1874 +  }
1875 +  
1876 +  x[nPoints]=x[0];
1877 +  y[nPoints]=y[0];
1878 +  
1879 +  TPolyLine *l = new TPolyLine(nPoints+1,x,y);
1880 +  l->SetFillColor(TColor::GetColor("#5858FA"));
1881 +  l->SetLineColor(TColor::GetColor("#5858FA"));
1882 +  l->SetLineWidth(1);
1883 +  return l;
1884 + }  
1885 +
1886 + //code for ReplaceAll copied from
1887 + //http://stackoverflow.com/questions/5343190/c-how-do-i-replace-all-instances-of-of-a-string-with-another-string
1888 + string ReplaceAll(string str, string from, string to) {
1889 +    size_t start_pos = 0;
1890 +    while((start_pos = str.find(from, start_pos)) != std::string::npos) {
1891 +         str.replace(start_pos, from.length(), to);
1892 +         start_pos += to.length(); // ...
1893 +    }
1894 +    return str;
1895 + }
1896 +
1897 + TPolyLine* GetFitUncertaintyShape(TH1 *histo, string fitname, float low, float high) {
1898 +  TF1 *fit = (TF1*)histo->GetFunction(fitname.c_str());
1899 +  if(!fit) {
1900 +    cout << "CANNOT PRODUCE FIT UNCERTAINTY AS THERE IS NO FIT CALLED " << fitname << " IN " << histo->GetName() << endl;
1901 +    return NULL;
1902 +  }
1903 +  float x0=histo->GetBinLowEdge(1);
1904 +  float x1=histo->GetBinLowEdge(histo->GetNbinsX())+histo->GetBinWidth(histo->GetNbinsX());
1905 +  return GetFitUncertaintyShape(fit, low, high,x0,x1);
1906 + }
1907 +
1908 + TPolyLine* GetFitUncertaintyShape(TGraphErrors *gr, string fitname, float low, float high, float minx=-999, float maxx=-999) {
1909 +  TF1 *fit = (TF1*)gr->GetFunction(fitname.c_str());
1910 +  if(!fit) {
1911 +    cout << "CANNOT PRODUCE FIT UNCERTAINTY AS THERE IS NO FIT CALLED " << fitname << " IN " << gr->GetName() << endl;
1912 +    return NULL;
1913 +  }
1914 +  double x,y;
1915 +  gr->GetPoint(0,x,y);
1916 +  float min=x;
1917 +  gr->GetPoint(gr->GetN()-1,x,y);
1918 +  float max=x;
1919 +  if(minx>-999 && maxx>-999 && minx<maxx) {
1920 +    max=maxx;
1921 +    min=minx;
1922 +  }
1923 +  return GetFitUncertaintyShape(fit, low, high,min,max);
1924 + }
1925 +
1926   stringstream all_bugs;
1927  
1928   void bug_tracker(string function, int line, string description) {
# Line 1751 | Line 1960 | void bug_tracker(string function, int li
1960   //TODO: Write a bug summary at the end.
1961  
1962   float pSTAR(float mglu, float mlsp, float mchi2) {
1963 <        float mz=91.2;
1963 >        float mz=91.0;
1964          float res=((mchi2*mchi2)-(mlsp+mz)*(mlsp+mz));
1965          res*=((mchi2*mchi2)-(mlsp-mz)*(mlsp-mz));
1966          res=TMath::Sqrt(res)/(2*mchi2);
# Line 1764 | Line 1973 | float generalizedpSTAR(float massmother,
1973          res=TMath::Sqrt(res)/(2*massmother);
1974          return res;
1975   }
1976 +
1977 + void FindMinMax(TGraphErrors *gr, float &min, float &max) {
1978 +  double x,y;
1979 +  gr->GetPoint(0,x,y);
1980 +  min=500*y;
1981 +  max=0.0002*y;
1982 +  for(int i=0;i<gr->GetN();i++) {
1983 +    gr->GetPoint(i,x,y);
1984 +    cout << "Point at " << x << "/" << y << " with error " << gr->GetErrorY(i) << endl;
1985 +    float potmin=y-1.5*gr->GetErrorY(i); // taking 1.5 times the error to allow for enough space to present the plot
1986 +    float potmax=y+1.5*gr->GetErrorY(i); // taking 1.5 times the error to allow for enough space to present the plot
1987 +    if(potmin<min) min=potmin;
1988 +    if(potmax>max) max=potmax;
1989 +  }
1990 + }
1991 +
1992 +
1993 + void CleanLegends() {
1994 +  cout << "----------------------- before deleting anything -----------------------" << endl;
1995 +  gROOT->ls();
1996 +  cout << "---------------------- /before deleting anything -----------------------" << endl;
1997 +  
1998 +  for(int ihist=(int)PlottingSetup::FakeHistoHeap.size()-1;ihist>=0;ihist--) {
1999 +    cout << "Deleting item " << ihist << endl;
2000 +    cout << "    " << PlottingSetup::FakeHistoHeap[ihist]->GetName() << " , " << PlottingSetup::FakeHistoHeap[ihist]->GetTitle() << endl;
2001 +    PlottingSetup::FakeHistoHeap[ihist]->Delete();
2002 +    PlottingSetup::FakeHistoHeap.pop_back();
2003 +  }
2004 + }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines