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.16 by fronga, Thu Aug 23 15:53:06 2012 UTC vs.
Revision 1.27 by buchmann, Wed Dec 12 14:27:54 2012 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   }
70  
71   bool dopng=false;
# Line 582 | Line 585 | void DrawPrelim(float writelumi=generalt
585    string prel=" Preliminary";
586    if(PlottingSetup::PaperMode) prel="";
587    //prelimtext << "CMS Preliminary 2011 , #sqrt{s}= 7 TeV, L= O(1) fb^{-1}"; //temporary replacement
588 +  if(PlottingSetup::is53reco) prel += " 53X";
589    string energy="7 TeV";
590    if(PlottingSetup::is2012) energy="8 TeV";
591    if(writelumi == 0) {
# Line 944 | Line 948 | string interpret_cut(string incut, bool
948    
949      if(incut=="id1==0") return "ee";
950      if(incut=="id1==1") return "#mu#mu";
951 <    if(incut=="abs(mll-91.2)<20") return "|m_{l^{+}l^{-}}-m_{Z}|<20";
951 >    if(incut=="abs(mll-91)<10") return "|m_{l^{+}l^{-}}-m_{Z}|<10";
952      if(incut=="pfJetGoodID[0]") return "";
953      if(incut=="pfJetGoodID[1]") return "";
954      if((int)incut.find("pfJetGoodNum")>-1) {
# Line 1313 | Line 1317 | Double_t MarcosChi2TestX(const TH1* h1,
1317     return prob;
1318   }
1319  
1316 /*void save_with_difference(TH1F *obs, vector<TH1F*> predictions, vector<float> prediction_weights, vector<float> prediction_uncerts, TVirtualPad *canvas, string savemeas) {
1317  if(predictions.size()<1) {
1318    write_warning(__FUNCTION__,"Cannot make a pred-obs prediction plot because the prediction vector is empty. aborting!");
1319    return;
1320  }
1321  TH1F *comppred = predictions[0]->Clone("comppred");
1322  comppred->Scale(prediction_weights[0]);
1323  for(int i=1;i<predictions.size();i++) comppred->Add(predictions[i],prediction_weights[i]);
1324  
1325  TH1F *statpred = (TH1F*)comppred->Clone("statpred");
1326  for(int i=1;i<=statpred->GetNbinsX();i++) statpred->SetBinContent(0);
1327  
1328  TH1F *syspred = predictions[0];
1329  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));
1330  for(int i=1;i<predictions.size();i++) {
1331    for(int ibin=1;ibin<=predictions[i]->GetNbinsX();ibin++) {
1332      float newentry=syspred->GetBinContent(ibin);
1333      newentry+=(predictions[i]->GetBinContent(ibin))*(predictions[i]->GetBinContent(ibin))*prediction_weights[i]*prediction_weights[i]*prediction_uncerts[i]*prediction_uncerts[i];
1334      syspred->SetBinContent(ibin,newentry);
1335    }
1336  }
1337  for(int ibin=1;ibin<=syspred->GetNbinsX();ibin++) {
1338    syspred->SetBinError(ibin,syspred->GetBinContent(i)+statpred->GetBinError());
1339    syspred->SetBinContent(ibin,0);
1340  }
1341  
1342  
1343 }
1344 */
1320   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) {
1321    //this function saves the pad being passed as well as a new one including the ratio.
1322    CompleteSave(canvas,savemeas);
# Line 1418 | Line 1393 | void save_with_ratio(TH1F *nominator, TH
1393       }
1394      
1395    }
1421  eratio->SetFillColor(TColor::GetColor("#00ADE1"));
1396    
1397 +  if(syshisto && !do_bpred_ratio) {
1398 +       SystDown = (TH1F*)syshisto->Clone("SystDown");
1399 +       SystUp = (TH1F*)syshisto->Clone("SystUp");
1400 +       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1401 +         SystDown->SetBinContent(i,1-syshisto->GetBinContent(i));
1402 +         SystUp->SetBinContent(i,1+syshisto->GetBinContent(i));
1403 +       }
1404 +       SystDown->SetLineColor(TColor::GetColor("#006DE1"));
1405 +       SystUp->SetLineColor(TColor::GetColor("#006DE1"));
1406 +  }
1407 +  eratio->SetFillColor(TColor::GetColor("#00ADE1"));
1408  
1409    ratio->SetTitle("");
1410 <  ratio->GetYaxis()->SetRangeUser(0.5,1.5);
1410 >  ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1411    if(do_bpred_ratio) ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1412    if(extendrange) ratio->GetYaxis()->SetRangeUser(0.0,4.0);
1413    ratio->GetXaxis()->SetTitle(nominator->GetXaxis()->GetTitle());
# Line 1498 | Line 1483 | void save_with_ratio(TH1F *nominator, TH
1483    delete main_canvas;
1484   }
1485  
1486 + 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) {
1487 +  //this function saves the pad being passed as well as a new one including the SysRatio.
1488 +  CompleteSave(canvas,savemeas);
1489  
1490 < TH1F* CollapseStack(THStack stack) {
1490 >  float bottommargin=gStyle->GetPadBottomMargin();
1491 >  float canvas_height=gStyle->GetCanvasDefH();
1492 >  float canvas_width=gStyle->GetCanvasDefW();
1493 >  float ratiospace=0.25;// space the ratio should take up (relative to original pad)
1494 >  
1495 >  float ratiobottommargin=0.3;
1496 >  float ratiotopmargin=0.1;
1497 >  
1498 >  float xstretchfactor=((1-ratiospace)*(1-gStyle->GetPadTopMargin()))/((1)*ratiospace);
1499 >  
1500 >  TCanvas *main_canvas = new TCanvas("main_canvas","main_canvas",(Int_t)canvas_width,(Int_t)(canvas_height*(1+ratiospace)));
1501 >  TPad *mainpad = new TPad("mainpad","mainpad",0,1-(1.0/(1+ratiospace)),1,1);//top (main) pad
1502 >  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
1503 >  TPad *bottompad = new TPad("bottompad", "Ratio Pad",0,0,1,(1-(1-bottommargin)/(1+ratiospace))-0.015); //bottom pad
1504 >  
1505 >  main_canvas->Range(0,0,1,1);
1506 >  main_canvas->SetBorderSize(0);
1507 >  main_canvas->SetFrameFillColor(0);
1508 >  
1509 >  mainpad->Draw();
1510 >  mainpad->cd();
1511 >  mainpad->Range(0,0,1,1);
1512 >  mainpad->SetFillColor(kWhite);
1513 >  mainpad->SetBorderSize(0);
1514 >  mainpad->SetFrameFillColor(0);
1515 >  canvas->Range(0,0,1,1);
1516 >  canvas->Draw("same");
1517 >  mainpad->Modified();
1518 >  main_canvas->cd();
1519 >  coverpad->Draw();
1520 >  coverpad->cd();
1521 >  coverpad->Range(0,0,1,1);
1522 >  coverpad->SetFillColor(kWhite);
1523 >  coverpad->SetBorderSize(0);
1524 >  coverpad->SetFrameFillColor(0);
1525 >  coverpad->Modified();
1526 >  main_canvas->cd();
1527 >  bottompad->SetTopMargin(ratiotopmargin);
1528 >  bottompad->SetBottomMargin(ratiobottommargin);
1529 >  bottompad->Draw();
1530 >  bottompad->cd();
1531 >  bottompad->Range(0,0,1,1);
1532 >  bottompad->SetFillColor(kWhite);
1533 >  TH1F *ratio = (TH1F*)nominator->Clone(GetNumericHistoName().c_str());
1534 >  ratio->Divide(denominator);
1535 >
1536 >
1537 >  TGraphAsymmErrors *eratio;
1538 >  TGraphAsymmErrors *SysEnvelope = new TGraphAsymmErrors();
1539 >  
1540 >  if(syshisto && !do_bpred_ratio) {
1541 >       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1542 >         float dx=syshisto->GetBinWidth(i)/2.0;
1543 >         float dy=syshisto->GetBinContent(i);
1544 >         if(dy<0.01) dy=0.07;
1545 >         SysEnvelope->SetPoint(i-1,syshisto->GetBinCenter(i),1.0);
1546 >         SysEnvelope->SetPointError(i-1,dx,dx,dy,dy);
1547 >       }
1548 >       SysEnvelope->SetFillColor(TColor::GetColor("#006DE1"));
1549 >  }
1550 >
1551 >  ratio->SetTitle("");
1552 >  ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1553 >  if(do_bpred_ratio) ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1554 >  if(extendrange) ratio->GetYaxis()->SetRangeUser(0.0,4.0);
1555 >  ratio->GetXaxis()->SetTitle(nominator->GetXaxis()->GetTitle());
1556 >  ratio->GetXaxis()->CenterTitle();
1557 >  ratio->GetYaxis()->SetTitle(yaxistitle.c_str());
1558 >  ratio->GetYaxis()->SetTitleOffset(0.4);
1559 >  ratio->GetYaxis()->CenterTitle();
1560 >  ratio->SetStats(0);
1561 >  ratio->GetXaxis()->SetLabelSize(xstretchfactor*ratio->GetXaxis()->GetLabelSize());
1562 >  ratio->GetYaxis()->SetLabelSize(xstretchfactor*ratio->GetYaxis()->GetLabelSize());
1563 >  ratio->GetXaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
1564 >  ratio->GetYaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
1565 >  ratio->GetYaxis()->SetNdivisions(502,false);
1566 >  ratio->SetFillColor(TColor::GetColor("#58D3F7"));
1567 >  ratio->Draw("e1");
1568 >  
1569 >  if(syshisto!=0) {
1570 >    SysEnvelope->SetFillColor(TColor::GetColor("#FE9A2E"));
1571 >    SysEnvelope->Draw("2,same");
1572 >    ratio->Draw("e1,same");
1573 >  } else {
1574 >    eratio->Draw("1");
1575 >  }
1576 >  ratio->Draw("same,axis");
1577 >  TLine *oneline = new TLine(ratio->GetXaxis()->GetBinLowEdge(1),1,ratio->GetXaxis()->GetBinLowEdge(ratio->GetNbinsX())+ratio->GetXaxis()->GetBinWidth(ratio->GetNbinsX()),1);
1578 >  oneline->SetLineStyle(2);
1579 >  oneline->SetLineColor(kBlue);
1580 >  oneline->Draw("same");
1581 >
1582 >  main_canvas->cd();
1583 >  main_canvas->Modified();
1584 >  main_canvas->cd();
1585 >  main_canvas->SetSelected(main_canvas);
1586 >  
1587 >  CompleteSave(main_canvas,savemeas+"_withSysRatio");
1588 >  bottompad->cd();
1589 >  
1590 >  Double_t chi2;
1591 >  Int_t ndf,igood;
1592 >  Double_t chi2prob = MarcosChi2TestX(nominator,denominator,chi2,ndf,igood,"UW");
1593 >
1594 >  stringstream Chi2text;
1595 >  Chi2text << "#chi^{2} / ndf: " << chi2 << " / " << ndf;
1596 >  stringstream Chi2probtext;
1597 >  Chi2probtext << "#chi^{2} prob: " << chi2prob;
1598 >  TText* chi2box = write_text(0.02,0.11,Chi2text.str());
1599 >  chi2box->SetTextSize(0.08);
1600 >  chi2box->SetTextAlign(11);
1601 >  chi2box->Draw();
1602 >  TText* chi2probbox = write_text(0.02,0.04,Chi2probtext.str());
1603 >  chi2probbox->SetTextSize(0.08);
1604 >  chi2probbox->SetTextAlign(11);
1605 >  chi2probbox->Draw();
1606 >  
1607 >  stringstream CompRatio;
1608 >  CompRatio << "Ratio: " << nominator->Integral() << " / " << denominator->Integral() << " \n " << nominator->Integral()/denominator->Integral();
1609 >  TText *CompleteRatio = write_text(0.98,0.04,CompRatio.str());
1610 >  CompleteRatio->SetTextSize(0.08);
1611 >  CompleteRatio->SetTextAlign(31);
1612 >  CompleteRatio->Draw();
1613 >  
1614 >  CompleteSave(main_canvas,savemeas+"_withSysRatio_and_Chi2");
1615 >  delete main_canvas;
1616 > }
1617 >
1618 > TH1F* CollapseStack(THStack stack,TString hname="base") {
1619     TH1F *bhist = ((TH1F*)((stack.GetHists())->At(0)));
1620 <   TH1F *basehisto = (TH1F*)bhist->Clone("base");
1620 >   TH1F *basehisto = (TH1F*)bhist->Clone(hname);
1621     TIter next(stack.GetHists());
1622     TH1F *h;
1623     int counter=0;
# Line 1713 | Line 1829 | string removefunnystring(string name) {
1829    return name;
1830   }
1831  
1832 + TPolyLine* GetFitUncertaintyShape(TF1 *fit, float low, float high,float x0,float x1) {
1833 +  int nPoints=1000;
1834 +  double x[nPoints+1];
1835 +  double y[nPoints];
1836 +  
1837 +  float par1=fit->GetParameter(0);
1838 +  float dpar1=fit->GetParError(0);
1839 +  float par2=fit->GetParameter(1);
1840 +  float dpar2=fit->GetParError(1);
1841 +    
1842 +  float step=(x1-x0)/(nPoints/2.0 -1);
1843 +  for(int i=0;i<nPoints/2.0;i++) {
1844 +    x[i]=x0+i*step;
1845 +    y[i]=(par1+x[i]*par2)    -    sqrt(dpar1*dpar1+dpar2*dpar2*x[i]*x[i]);
1846 +    x[nPoints-1-i]=x[i];
1847 +    y[nPoints-1-i]=(par1+x[i]*par2)    +    sqrt(dpar1*dpar1+dpar2*dpar2*x[i]*x[i]);
1848 +    if(y[i]<low) y[i]=low;
1849 +    if(y[i]>high) y[i]=high;
1850 +    if(y[nPoints-1-i]<low) y[nPoints-1-i]=low;
1851 +    if(y[nPoints-1-i]>high) y[nPoints-1-i]=high;
1852 +  }
1853 +  
1854 +  x[nPoints]=x[0];
1855 +  y[nPoints]=y[0];
1856 +  
1857 +  TPolyLine *l = new TPolyLine(nPoints+1,x,y);
1858 +  l->SetFillColor(TColor::GetColor("#5858FA"));
1859 +  l->SetLineColor(TColor::GetColor("#5858FA"));
1860 +  l->SetLineWidth(1);
1861 +  return l;
1862 + }  
1863 +
1864 + //code for ReplaceAll copied from
1865 + //http://stackoverflow.com/questions/5343190/c-how-do-i-replace-all-instances-of-of-a-string-with-another-string
1866 + string ReplaceAll(string str, string from, string to) {
1867 +    size_t start_pos = 0;
1868 +    while((start_pos = str.find(from, start_pos)) != std::string::npos) {
1869 +         str.replace(start_pos, from.length(), to);
1870 +         start_pos += to.length(); // ...
1871 +    }
1872 +    return str;
1873 + }
1874 +
1875 + TPolyLine* GetFitUncertaintyShape(TH1 *histo, string fitname, float low, float high) {
1876 +  TF1 *fit = (TF1*)histo->GetFunction(fitname.c_str());
1877 +  if(!fit) {
1878 +    cout << "CANNOT PRODUCE FIT UNCERTAINTY AS THERE IS NO FIT CALLED " << fitname << " IN " << histo->GetName() << endl;
1879 +    return NULL;
1880 +  }
1881 +  float x0=histo->GetBinLowEdge(1);
1882 +  float x1=histo->GetBinLowEdge(histo->GetNbinsX())+histo->GetBinWidth(histo->GetNbinsX());
1883 +  return GetFitUncertaintyShape(fit, low, high,x0,x1);
1884 + }
1885 +
1886 + TPolyLine* GetFitUncertaintyShape(TGraphErrors *gr, string fitname, float low, float high, float minx=-999, float maxx=-999) {
1887 +  TF1 *fit = (TF1*)gr->GetFunction(fitname.c_str());
1888 +  if(!fit) {
1889 +    cout << "CANNOT PRODUCE FIT UNCERTAINTY AS THERE IS NO FIT CALLED " << fitname << " IN " << gr->GetName() << endl;
1890 +    return NULL;
1891 +  }
1892 +  double x,y;
1893 +  gr->GetPoint(0,x,y);
1894 +  float min=x;
1895 +  gr->GetPoint(gr->GetN()-1,x,y);
1896 +  float max=x;
1897 +  if(minx>-999 && maxx>-999 && minx<maxx) {
1898 +    max=maxx;
1899 +    min=minx;
1900 +  }
1901 +  return GetFitUncertaintyShape(fit, low, high,min,max);
1902 + }
1903 +
1904   stringstream all_bugs;
1905  
1906   void bug_tracker(string function, int line, string description) {
# Line 1750 | Line 1938 | void bug_tracker(string function, int li
1938   //TODO: Write a bug summary at the end.
1939  
1940   float pSTAR(float mglu, float mlsp, float mchi2) {
1941 <        float mz=91.2;
1941 >        float mz=91.0;
1942          float res=((mchi2*mchi2)-(mlsp+mz)*(mlsp+mz));
1943          res*=((mchi2*mchi2)-(mlsp-mz)*(mlsp-mz));
1944          res=TMath::Sqrt(res)/(2*mchi2);
# Line 1763 | Line 1951 | float generalizedpSTAR(float massmother,
1951          res=TMath::Sqrt(res)/(2*massmother);
1952          return res;
1953   }
1954 +
1955 + void FindMinMax(TGraphErrors *gr, float &min, float &max) {
1956 +  double x,y;
1957 +  gr->GetPoint(0,x,y);
1958 +  min=500*y;
1959 +  max=0.0002*y;
1960 +  for(int i=0;i<gr->GetN();i++) {
1961 +    gr->GetPoint(i,x,y);
1962 +    cout << "Point at " << x << "/" << y << " with error " << gr->GetErrorY(i) << endl;
1963 +    float potmin=y-1.5*gr->GetErrorY(i); // taking 1.5 times the error to allow for enough space to present the plot
1964 +    float potmax=y+1.5*gr->GetErrorY(i); // taking 1.5 times the error to allow for enough space to present the plot
1965 +    if(potmin<min) min=potmin;
1966 +    if(potmax>max) max=potmax;
1967 +  }
1968 + }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines