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.1 by buchmann, Mon Jan 30 14:46:25 2012 UTC vs.
Revision 1.34 by buchmann, Thu Feb 7 14:30:09 2013 UTC

# Line 17 | Line 17
17   #include <netdb.h>
18   #include <stdio.h>
19   #include <string.h>
20 + #include <assert.h>
21  
22   #include <TFile.h>
23   #include <TTree.h>
# Line 38 | 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 58 | Line 62 | namespace PlottingSetup {
62      vector<float> global_ratio_binning;
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 +
71   }
72  
73   bool dopng=false;
# Line 119 | Line 129 | string topdir(string child) {
129          return tempdirectory.substr(0,ichar);
130        }
131      }
132 +    return "TOPDIRFAILURE";
133   }
134  
135   template < typename CHAR_TYPE,
# Line 250 | 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;
277 +  dout << "\033[1;31m  __| | ___  __ _  __| |   ___ _ __   __| | " << 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 << "" << endl;
283 +  dout << " A totally fatal error has occurred: " << endl;
284 +  dout << " " << message << endl;
285 +  dout << endl;
286 +  dout << "If you do not know how to deal with this error message please call Customer Support Mumbai (or Zurich...) \033[0m" << endl;
287 +  dout << "   ... aborting now ... " << endl;
288 +  assert(0);
289 + }
290 +
291 +
292 + void SanityChecks() {
293 +  //this function gets called whenever samples are initialized (see ActiveSamples.C)
294 +  //whatever is necessary for the code to run should be put here.
295 +  
296 +  if(gROOT->GetVersionInt()<53200) {
297 +    DeadEnd("ROOT version outdated, you are using "+any2string(gROOT->GetVersionInt())+" but the earliest recommended version is 5.32.00");
298 +  }
299 +  
300 +  
301 + }
302 +
303   void initialize_log() {
304    if(!PlottingSetup::publicmode) {
305      dout << "____________________________________________________________" << endl;
# Line 278 | Line 328 | void initialize_log() {
328      dout << "        PUBLIC MODE     " << endl;
329    }
330    time_t rawtime;
281  struct tm * timeinfo;
331    time ( &rawtime );
332    dout << "       Analysis run on " << asctime (localtime ( &rawtime ));
333    dout << "____________________________________________________________" << endl;
# Line 289 | Line 338 | void extract_cbaf_dir(string curpath) {
338    int position=curpath.find("/Plotting");
339    if(position<0) position=curpath.find("/DistributedModelCalculations");
340    if(position<0) position=curpath.find("/various_assignments");
341 +  if(position<0) position=curpath.find("/exchange");
342    PlottingSetup::cbafbasedir=curpath.substr(0,position);
343   }
344  
# Line 296 | Line 346 | void set_directory(string basedir="") {
346    if(basedir.substr(0,1)=="/") basedir=basedir.substr(1,basedir.length()-1);
347    if(basedir.substr(basedir.length()-1,1)!="/") basedir+="/";
348    char currentpath[1024];
349 <  char *path = getcwd(currentpath,1024);
349 >  getcwd(currentpath,1024);
350    PlottingSetup::basedirectory=(string)currentpath+"/Plots/"+basedir;
351    ensure_directory_exists(PlottingSetup::basedirectory);
352    initialize_log();
# Line 321 | Line 371 | string extract_directory(string savethis
371   string extract_root_dir(string name) {
372      int position = -1;
373      if(name.substr(0,1)=="/") name=name.substr(1,name.length()-1);
374 <    for(int ipos=0;ipos<name.length();ipos++) {
374 >    for(int ipos=0;ipos<(int)name.length();ipos++) {
375          if(name.substr(ipos,1)=="/") position=ipos;
376      }
377      if(position==-1) return "";
# Line 331 | Line 381 | string extract_root_dir(string name) {
381   string extract_root_filename(string name) {
382      int position = -1;
383      if(name.substr(0,1)=="/") name=name.substr(1,name.length()-1);
384 <    for(int ipos=0;ipos<name.length();ipos++) {
384 >    for(int ipos=0;ipos<(int)name.length();ipos++) {
385          if(name.substr(ipos,1)=="/") position=ipos;
386      }
387      return name.substr(position+1,name.length()-position-1);
# Line 389 | 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 414 | Line 464 | void write_entry(string entry,int width,
464    int currwidth=entry.size();
465    while(currwidth<width) {
466      entry=" "+entry;
467 <    if(entry.size()<width) entry=entry+" ";
467 >    if((int)entry.size()<width) entry=entry+" ";
468      currwidth=entry.size();
469    }
470    bool do_special=false;
# Line 428 | Line 478 | void write_entry(string entry,int width,
478   void make_nice_table(vector<vector <string> > &entries) {
479    int ncolumns=write_first_line(entries);
480    int cellwidth=(int)(((float)(60-(ncolumns+1)))/(ncolumns));
481 <  for(int iline=0;iline<entries.size();iline++) {
481 >  for(int iline=0;iline<(int)entries.size();iline++) {
482      vector<string> currline = entries[iline];
483      dout << " |";
484 <    for(int ientry=0;ientry<currline.size();ientry++) {
484 >    for(int ientry=0;ientry<(int)currline.size();ientry++) {
485        write_entry(currline[ientry],cellwidth);
486      }
487      dout << endl;
# Line 443 | Line 493 | void make_nice_table(vector<vector <stri
493   void make_nice_jzb_table(vector<vector <string> > &entries) {
494    int ncolumns=write_first_line(entries);
495    int cellwidth=(int)(((float)(60-(ncolumns+1)))/(ncolumns));
496 <  for(int iline=0;iline<entries.size();iline++) {
496 >  for(int iline=0;iline<(int)entries.size();iline++) {
497      vector<string> currline = entries[iline];
498      dout << " |";
499 <    for(int ientry=0;ientry<currline.size();ientry++) {
499 >    for(int ientry=0;ientry<(int)currline.size();ientry++) {
500        write_entry(currline[ientry],cellwidth,iline,ientry);
501      }
502      dout << endl;
# Line 547 | Line 597 | void DrawPrelim(float writelumi=generalt
597    string prel=" Preliminary";
598    if(PlottingSetup::PaperMode) prel="";
599    //prelimtext << "CMS Preliminary 2011 , #sqrt{s}= 7 TeV, L= O(1) fb^{-1}"; //temporary replacement
600 +  if(PlottingSetup::is53reco) prel += " 53X";
601 +  string energy="7 TeV";
602 +  if(PlottingSetup::is2012) energy="8 TeV";
603    if(writelumi == 0) {
604 <    if(isMC) prelimtext << "CMS Simulation, #sqrt{s} = 7 TeV";
605 <    else prelimtext << "CMS" << prel << ", #sqrt{s} = 7 TeV";
604 >    if(isMC) prelimtext << "CMS Simulation, #sqrt{s} = " << energy;
605 >    else prelimtext << "CMS" << prel << ", #sqrt{s} = " << energy;
606    } else {
607 <    if(isMC) prelimtext << "CMS Simulation, #sqrt{s} = 7 TeV, L_{int} = " << std::setprecision(2) <<writelumi<<" "<<barn<<"^{-1}";
608 <    else prelimtext << "CMS" << prel << ", #sqrt{s} = 7 TeV, L_{int} = " << std::setprecision(2) <<writelumi<<" "<<barn<<"^{-1}";
607 >    if(isMC) prelimtext << "CMS Simulation, #sqrt{s} = " << energy << ", L_{int} = " << std::setprecision(3) <<writelumi<<" "<<barn<<"^{-1}";
608 >     else prelimtext << "CMS" << prel << ", #sqrt{s} = " << energy << ", L_{int} = " << std::setprecision(3) <<writelumi<<" "<<barn<<"^{-1}";
609    }
610    TPaveText *eventSelectionPaveText = new TPaveText(0.27, 0.93,0.77, 1.0,"blNDC");
611    eventSelectionPaveText->SetFillStyle(4000);
612    eventSelectionPaveText->SetBorderSize(0);
613    eventSelectionPaveText->SetFillColor(kWhite);
614    eventSelectionPaveText->SetTextFont(42);
615 +  //eventSelectionPaveText->SetTextFont(62); // paper style
616    eventSelectionPaveText->SetTextSize(0.042);
617    eventSelectionPaveText->AddText(prelimtext.str().c_str());
618    eventSelectionPaveText->Draw();
# Line 596 | Line 650 | TGraph* make_nice_ratio(int nbins,float
650      bincontent[i]=0;
651      errors[i]=0;
652    }
653 <  float currlimit=binning[0];
653 > //  float currlimit=binning[0];
654    int currtoplim=1;
655    for(int ibin=1;ibin<=histo->GetNbinsX();ibin++)
656    {
# Line 629 | Line 683 | float computeProductError(float a, float
683      return val;
684   }
685  
686 + TGraphAsymmErrors *SysAndStatRatio(TH1F *h1,TH1F *h2, TH1F *sys, int id, vector<float>binning, bool precise=false) {
687 +    int absJZBbinsNumber = binning.size()-1;
688 +    TGraphAsymmErrors* graph = new TGraphAsymmErrors(absJZBbinsNumber);
689 +    for(int i=0;i<absJZBbinsNumber;i++)
690 +        {
691 +                float xCenter=h1->GetBinCenter(i+1);
692 +                float xWidth=(h1->GetBinWidth(i+1))*0.5;
693 +                float nominatorError = h1->GetBinError(i+1);
694 +                float nominator=h1->GetBinContent(i+1);
695 +                float denominatorError=h2->GetBinError(i+1);
696 +                float denominator=h2->GetBinContent(i+1);
697 +                float errorN = 0;
698 +                float errorP = computeRatioError(nominator,nominatorError,denominator,denominatorError);
699 +                float syserr=sys->GetBinContent(i+1);
700 +                float errorsys=0;
701 +                if(id==1) // (is data)
702 +                {
703 +                        if(!precise) errorP = computeRatioError(nominator,statErrorP(nominator),denominator,statErrorP(denominator));
704 +                        else         errorP = computeRatioError(nominator,nominatorError,denominator,denominatorError);
705 +                        errorsys=computeRatioError(nominator,syserr,denominator,0);
706 +                        errorP = TMath::Sqrt(errorsys*errorsys+errorP*errorP);
707 +                        errorN = errorP; // symmetrize using statErrorP
708 +                } else {
709 +                        errorN = computeRatioError(nominator,nominatorError,denominator,denominatorError);
710 +                        errorsys=computeRatioError(nominator,syserr,denominator,0);
711 +                        errorN = TMath::Sqrt(errorsys*errorsys+errorP*errorP);
712 +                        errorP = errorN; // symmetrize using statErrorP
713 +                }
714 +                if(denominator!=0) {
715 +                        graph->SetPoint(i, xCenter, nominator/denominator);
716 +                        graph->SetPointError(i,xWidth,xWidth,errorN,errorP);
717 +                }
718 +                else {
719 +                        graph->SetPoint(i, xCenter, -999);
720 +                        graph->SetPointError(i,xWidth,xWidth,errorN,errorP);
721 +                }
722 +        }      
723 +        return graph;
724 + }
725 +    
726 +
727   TGraphAsymmErrors *histRatio(TH1F *h1,TH1F *h2, int id, vector<float>binning, bool precise=false)
728   {
729          int absJZBbinsNumber = binning.size()-1;
730          TGraphAsymmErrors* graph = new TGraphAsymmErrors(absJZBbinsNumber);
731          
732 <        for(unsigned int i=0;i<absJZBbinsNumber;i++)
732 >        for(int i=0;i<absJZBbinsNumber;i++)
733          {
734                  float xCenter=h1->GetBinCenter(i+1);
735                  float xWidth=(h1->GetBinWidth(i+1))*0.5;
# Line 746 | Line 841 | void splitupcut(string incut, vector<str
841      int paranthesis_open=0;
842      int substr_start=0;
843      string currchar="";
844 <    for (int ichar=0;ichar<incut.length();ichar++)
844 >    for (int ichar=0;ichar<(int)incut.length();ichar++)
845      {
846          currchar=incut.substr(ichar,1);
847   //        if(currchar=="(") paranthesis_open++;
# Line 759 | Line 854 | void splitupcut(string incut, vector<str
854      partvector.push_back(incut.substr(substr_start,incut.length()-substr_start));
855      if(Verbosity>1) {
856          dout << "[ splitupcut() ] : The cut vector now contains the following elements: "<< endl;
857 <        for (int ipart=0;ipart<partvector.size();ipart++)
857 >        for (int ipart=0;ipart<(int)partvector.size();ipart++)
858          {
859              dout << "     - " << partvector[ipart] << endl;
860          }
# Line 768 | Line 863 | void splitupcut(string incut, vector<str
863  
864   int atleastvalue(string expression, int &morethanlessthan) // takes in an expression such as ">2" or ">=3" and returns e.g. 3 (in both examples)
865   {
866 <  int retval=0;
866 > //  int retval=0;
867    if(expression.substr(0,1)==">"&&expression.substr(1,1)=="=") {
868   //    dout << "The expression " << expression << " is saying that we have at least " << atoi(expression.substr(2,1).c_str()) << " jets" << endl;
869      morethanlessthan=1;
# Line 794 | Line 889 | int atleastvalue(string expression, int
889      morethanlessthan=-1;
890      return 1+atoi(expression.substr(2,1).c_str());
891    }
892 +  
893 +  return -1234567;
894   }
895  
896   int do_jet_cut(string incut, int *nJets) {
# Line 829 | Line 926 | int do_jet_cut(string incut, int *nJets)
926      dout << "Is of type >" << endl;
927      return 0;
928    }
929 +  return -12345;
930   }
931  
932   string interpret_cut(string incut, bool &isJetCut, int *permittednJets)
# Line 862 | 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 894 | Line 992 | string interpret_cuts(vector<string> &cu
992      int finalJetCut=-1;
993      int permittednJets[20];
994      for(int ijet=0;ijet<20;ijet++) permittednJets[ijet]=1;
995 <    int morethanlessthan=0;//-1: less than, 0: exactly, 1: more than
996 <    for(int icut=0;icut<cutparts.size();icut++)
995 > //    int morethanlessthan=0;//-1: less than, 0: exactly, 1: more than
996 >    for(int icut=0;icut<(int)cutparts.size();icut++)
997      {
998          if(icut==0) nicecut<<interpret_cut(cutparts[icut],isJetCut,permittednJets);
999          else {
# Line 960 | Line 1058 | TH1F * makehistofromfunction(TF1 *f1,TH1
1058    histo->SetTitle(histoname.str().c_str());
1059    histo->SetName(histoname.str().c_str());
1060    int nbins=histo->GetNbinsX();
1061 <  float low=histo->GetBinLowEdge(1);
1062 <  float hi=histo->GetBinLowEdge(histo->GetNbinsX())+histo->GetBinWidth(histo->GetNbinsX());
1061 > //  float low=histo->GetBinLowEdge(1);
1062 > //  float hi=histo->GetBinLowEdge(histo->GetNbinsX())+histo->GetBinWidth(histo->GetNbinsX());
1063    
1064    for(int i=0;i<=nbins;i++) {
1065      histo->SetBinContent(i,(f1->Integral(histo->GetBinLowEdge(i),histo->GetBinLowEdge(i)+histo->GetBinWidth(i)))/histo->GetBinWidth(i));
# Line 1209 | Line 1307 | TGraphAsymmErrors* produce_ratio_graph(T
1307   }
1308  
1309  
1310 + TText* WriteSelection(int njets) {
1311 +  string sel="Loose selection";
1312 +  if(njets==3) sel="Tight selection";
1313 +  assert(njets==2||njets==3);
1314 +  TText *sele = new TText(0.97,0.135,sel.c_str());
1315 +  sele->SetNDC(true);
1316 +  sele->SetTextColor(TColor::GetColor("#848484"));
1317 +  sele->SetTextFont(42);
1318 +  sele->SetTextAlign(32);
1319 +  sele->SetTextSize(0.03);
1320 +  sele->SetTextAngle(270);
1321 +  return sele;
1322 + }
1323 +
1324   Double_t MarcosChi2TestX(const TH1* h1, const TH1* h2, Double_t &chi2, Int_t &ndf, Int_t &igood, Option_t *option)
1325   {
1326  
# Line 1231 | Line 1343 | Double_t MarcosChi2TestX(const TH1* h1,
1343     return prob;
1344   }
1345  
1346 < void save_with_ratio(TH1F *nominator, TH1F *denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false, bool extendrange=false, string yaxistitle="ratio") {
1346 > 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) {
1347    //this function saves the pad being passed as well as a new one including the ratio.
1348 <  CompleteSave(canvas,savemeas);
1349 <
1348 >  
1349 >  orig_canvas->cd();
1350 >  orig_canvas->Update();
1351 >  CompleteSave(orig_canvas,savemeas);
1352 >  TVirtualPad *canvas = (TVirtualPad*) orig_canvas->Clone("TempCanvas");//otherwise the Ratio_main_canvas will own our pad and destroy it upon deletion
1353 >  
1354    float bottommargin=gStyle->GetPadBottomMargin();
1355    float canvas_height=gStyle->GetCanvasDefH();
1356    float canvas_width=gStyle->GetCanvasDefW();
# Line 1245 | Line 1361 | void save_with_ratio(TH1F *nominator, TH
1361    
1362    float xstretchfactor=((1-ratiospace)*(1-gStyle->GetPadTopMargin()))/((1)*ratiospace);
1363    
1364 <  TCanvas *main_canvas = new TCanvas("main_canvas","main_canvas",(Int_t)canvas_width,(Int_t)(canvas_height*(1+ratiospace)));
1364 >  TCanvas *Ratio_main_canvas = new TCanvas("Ratio_main_canvas","Ratio_main_canvas",(Int_t)canvas_width,(Int_t)(canvas_height*(1+ratiospace)));
1365    TPad *mainpad = new TPad("mainpad","mainpad",0,1-(1.0/(1+ratiospace)),1,1);//top (main) pad
1366    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
1367    TPad *bottompad = new TPad("bottompad", "Ratio Pad",0,0,1,(1-(1-bottommargin)/(1+ratiospace))-0.015); //bottom pad
1368    
1369 <  main_canvas->Range(0,0,1,1);
1370 <  main_canvas->SetBorderSize(0);
1371 <  main_canvas->SetFrameFillColor(0);
1369 >  Ratio_main_canvas->Range(0,0,1,1);
1370 >  Ratio_main_canvas->SetBorderSize(0);
1371 >  Ratio_main_canvas->SetFrameFillColor(0);
1372    
1373    mainpad->Draw();
1374    mainpad->cd();
# Line 1263 | Line 1379 | void save_with_ratio(TH1F *nominator, TH
1379    canvas->Range(0,0,1,1);
1380    canvas->Draw("same");
1381    mainpad->Modified();
1382 <  main_canvas->cd();
1382 >  Ratio_main_canvas->cd();
1383    coverpad->Draw();
1384    coverpad->cd();
1385    coverpad->Range(0,0,1,1);
# Line 1271 | Line 1387 | void save_with_ratio(TH1F *nominator, TH
1387    coverpad->SetBorderSize(0);
1388    coverpad->SetFrameFillColor(0);
1389    coverpad->Modified();
1390 <  main_canvas->cd();
1390 >  Ratio_main_canvas->cd();
1391    bottompad->SetTopMargin(ratiotopmargin);
1392    bottompad->SetBottomMargin(ratiobottommargin);
1393    bottompad->Draw();
# Line 1281 | Line 1397 | void save_with_ratio(TH1F *nominator, TH
1397    TH1F *ratio = (TH1F*)nominator->Clone(GetNumericHistoName().c_str());
1398    ratio->Divide(denominator);
1399  
1284
1400    TGraphAsymmErrors *eratio;
1401 +  TH1F *SystDown;
1402 +  TH1F *SystUp;
1403 +  
1404    if(!do_bpred_ratio) eratio = produce_ratio_graph(ratio);
1405    else {
1406       bool using_data=false;
1407       if((int)savemeas.find("Data")>0) using_data=true;
1408       eratio = histRatio(nominator,denominator,using_data,PlottingSetup::global_ratio_binning,false);
1409 +     if(syshisto!=0) {
1410 +       SystDown = (TH1F*)syshisto->Clone("SystDown");
1411 +       SystUp = (TH1F*)syshisto->Clone("SystUp");
1412 +       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1413 +         SystDown->SetBinContent(i,1-syshisto->GetBinContent(i));
1414 +         SystUp->SetBinContent(i,1+syshisto->GetBinContent(i));
1415 +       }
1416 +       SystDown->SetLineColor(TColor::GetColor("#006DE1"));
1417 +       SystUp->SetLineColor(TColor::GetColor("#006DE1"));
1418 +     }
1419       for(int i=1;i<=ratio->GetNbinsX();i++) {
1420          ratio->SetBinContent(i,0);
1421          ratio->SetBinError(i,0);
1422       }
1423 +    
1424 +  }
1425 +  
1426 +  if(syshisto && !do_bpred_ratio) {
1427 +       SystDown = (TH1F*)syshisto->Clone("SystDown");
1428 +       SystUp = (TH1F*)syshisto->Clone("SystUp");
1429 +       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1430 +         SystDown->SetBinContent(i,1-syshisto->GetBinContent(i));
1431 +         SystUp->SetBinContent(i,1+syshisto->GetBinContent(i));
1432 +       }
1433 +       SystDown->SetLineColor(TColor::GetColor("#006DE1"));
1434 +       SystUp->SetLineColor(TColor::GetColor("#006DE1"));
1435    }
1436    eratio->SetFillColor(TColor::GetColor("#00ADE1"));
1437  
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 1313 | Line 1453 | void save_with_ratio(TH1F *nominator, TH
1453    ratio->SetFillColor(TColor::GetColor("#58D3F7"));
1454    ratio->SetMarkerSize(0);
1455    ratio->Draw("e2");
1456 <  eratio->Draw("2");
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 >  
1462 >  ratio_center->SetFillColor(TColor::GetColor("#006381"));
1463 >  
1464 >  if(syshisto!=0) {
1465 > //    sysratio->Draw("2");
1466 > //    eratio->Draw("2,same");
1467 >    eratio->Draw("2");
1468 >    SystDown->Draw("histo,same");
1469 >    SystUp->Draw("histo,same");
1470 >  } else {
1471 >    eratio->Draw("2");
1472 >  }
1473 >  ratio_center->Draw("2");
1474    ratio->Draw("same,axis");
1475    TLine *oneline = new TLine(ratio->GetXaxis()->GetBinLowEdge(1),1,ratio->GetXaxis()->GetBinLowEdge(ratio->GetNbinsX())+ratio->GetXaxis()->GetBinWidth(ratio->GetNbinsX()),1);
1476    oneline->SetLineStyle(2);
1477    oneline->SetLineColor(kBlue);
1478    oneline->Draw("same");
1479  
1480 <  main_canvas->cd();
1481 <  main_canvas->Modified();
1482 <  main_canvas->cd();
1483 <  main_canvas->SetSelected(main_canvas);
1480 >  Ratio_main_canvas->cd();
1481 >  Ratio_main_canvas->Modified();
1482 >  Ratio_main_canvas->cd();
1483 >  Ratio_main_canvas->SetSelected(Ratio_main_canvas);
1484    
1485 <  CompleteSave(main_canvas,savemeas+"_withratio");
1485 >  CompleteSave(Ratio_main_canvas,savemeas+"_withratio");
1486    bottompad->cd();
1487    
1488    Double_t chi2;
1489    Int_t ndf,igood;
1490 <  Double_t res=0.0;
1490 > //  Double_t res=0.0;
1491    Double_t chi2prob = MarcosChi2TestX(nominator,denominator,chi2,ndf,igood,"UW");
1492  
1493    stringstream Chi2text;
# Line 1345 | Line 1502 | void save_with_ratio(TH1F *nominator, TH
1502    chi2probbox->SetTextSize(0.08);
1503    chi2probbox->SetTextAlign(11);
1504    chi2probbox->Draw();
1505 <  CompleteSave(main_canvas,savemeas+"_withratio_and_Chi2");
1505 >  
1506 >  stringstream CompRatio;
1507 >  CompRatio << "Ratio: " << nominator->Integral() << " / " << denominator->Integral() << " \n " << nominator->Integral()/denominator->Integral();
1508 >  TText *CompleteRatio = write_text(0.98,0.04,CompRatio.str());
1509 >  CompleteRatio->SetTextSize(0.08);
1510 >  CompleteRatio->SetTextAlign(31);
1511 >  CompleteRatio->Draw();
1512 >  
1513 >  CompleteSave(Ratio_main_canvas,savemeas+"_withratio_and_Chi2");
1514  
1515   //  float KS   = nominator->KolmogorovTest(denominator);
1516   //  stringstream KStext;
1517   //  Chi2text << "KS = " << KS << endl;
1518   //cout << "Found : " << KStext.str() << endl;
1519  
1520 <
1521 <  delete main_canvas;
1520 >  delete eratio;
1521 >  delete ratio_center;
1522 >  delete ratio;
1523 >  delete Ratio_main_canvas;
1524   }
1525  
1526 + void save_with_ratio_and_sys_band(TH1F *nominator, TH1F *denominator, TVirtualPad *orig_canvas, string savemeas, bool do_bpred_ratio=false, bool extendrange=false, string yaxistitle="ratio",TH1F *syshisto=NULL) {
1527 +  //this function saves the pad being passed as well as a new one including the SysRatio.
1528 +  orig_canvas->cd();
1529 +  orig_canvas->Update();
1530 +  CompleteSave(orig_canvas,savemeas);
1531 +  
1532 +  TVirtualPad *canvas = (TVirtualPad*) orig_canvas->Clone("TempCanvas");//otherwise the Ratio_main_canvas will own our pad and destroy it upon deletion
1533  
1534 < TH1F* CollapseStack(THStack stack) {
1534 >  float bottommargin=gStyle->GetPadBottomMargin();
1535 >  float canvas_height=gStyle->GetCanvasDefH();
1536 >  float canvas_width=gStyle->GetCanvasDefW();
1537 >  float ratiospace=0.25;// space the ratio should take up (relative to original pad)
1538 >  
1539 >  float ratiobottommargin=0.3;
1540 >  float ratiotopmargin=0.1;
1541 >  
1542 >  float xstretchfactor=((1-ratiospace)*(1-gStyle->GetPadTopMargin()))/((1)*ratiospace);
1543 >  
1544 >  TCanvas *Ratio_main_canvas = new TCanvas("Ratio_main_canvas","Ratio_main_canvas",(Int_t)canvas_width,(Int_t)(canvas_height*(1+ratiospace)));
1545 >  TPad *mainpad = new TPad("mainpad","mainpad",0,1-(1.0/(1+ratiospace)),1,1);//top (main) pad
1546 >  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
1547 >  TPad *bottompad = new TPad("bottompad", "Ratio Pad",0,0,1,(1-(1-bottommargin)/(1+ratiospace))-0.015); //bottom pad
1548 >  
1549 >  Ratio_main_canvas->Range(0,0,1,1);
1550 >  Ratio_main_canvas->SetBorderSize(0);
1551 >  Ratio_main_canvas->SetFrameFillColor(0);
1552 >  
1553 >  mainpad->Draw();
1554 >  mainpad->cd();
1555 >  mainpad->Range(0,0,1,1);
1556 >  mainpad->SetFillColor(kWhite);
1557 >  mainpad->SetBorderSize(0);
1558 >  mainpad->SetFrameFillColor(0);
1559 >  canvas->Range(0,0,1,1);
1560 >  canvas->Draw("same");
1561 >  mainpad->Modified();
1562 >  Ratio_main_canvas->cd();
1563 >  coverpad->Draw();
1564 >  coverpad->cd();
1565 >  coverpad->Range(0,0,1,1);
1566 >  coverpad->SetFillColor(kWhite);
1567 >  coverpad->SetBorderSize(0);
1568 >  coverpad->SetFrameFillColor(0);
1569 >  coverpad->Modified();
1570 >  Ratio_main_canvas->cd();
1571 >  bottompad->SetTopMargin(ratiotopmargin);
1572 >  bottompad->SetBottomMargin(ratiobottommargin);
1573 >  bottompad->Draw();
1574 >  bottompad->cd();
1575 >  bottompad->Range(0,0,1,1);
1576 >  bottompad->SetFillColor(kWhite);
1577 >  TH1F *ratio = (TH1F*)nominator->Clone(GetNumericHistoName().c_str());
1578 >  ratio->Divide(denominator);
1579 >
1580 >
1581 >  TGraphAsymmErrors *eratio;
1582 >  TGraphAsymmErrors *SysEnvelope = new TGraphAsymmErrors();
1583 >  
1584 >  if(syshisto && !do_bpred_ratio) {
1585 >       for(int i=1;i<=syshisto->GetNbinsX();i++) {
1586 >         float dx=syshisto->GetBinWidth(i)/2.0;
1587 >         float dy=syshisto->GetBinContent(i);
1588 >         if(dy<0.01) dy=0.07;
1589 >         SysEnvelope->SetPoint(i-1,syshisto->GetBinCenter(i),1.0);
1590 >         SysEnvelope->SetPointError(i-1,dx,dx,dy,dy);
1591 >       }
1592 >       SysEnvelope->SetFillColor(TColor::GetColor("#006DE1"));
1593 >  }
1594 >
1595 >  ratio->SetTitle("");
1596 >  ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1597 >  if(do_bpred_ratio) ratio->GetYaxis()->SetRangeUser(0.0,2.0);
1598 >  if(extendrange) ratio->GetYaxis()->SetRangeUser(0.0,4.0);
1599 >  ratio->GetXaxis()->SetTitle(nominator->GetXaxis()->GetTitle());
1600 >  ratio->GetXaxis()->CenterTitle();
1601 >  ratio->GetYaxis()->SetTitle(yaxistitle.c_str());
1602 >  ratio->GetYaxis()->SetTitleOffset(0.4);
1603 >  ratio->GetYaxis()->CenterTitle();
1604 >  ratio->SetStats(0);
1605 >  ratio->GetXaxis()->SetLabelSize(xstretchfactor*ratio->GetXaxis()->GetLabelSize());
1606 >  ratio->GetYaxis()->SetLabelSize(xstretchfactor*ratio->GetYaxis()->GetLabelSize());
1607 >  ratio->GetXaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
1608 >  ratio->GetYaxis()->SetTitleSize(xstretchfactor*gStyle->GetTitleSize());
1609 >  ratio->GetYaxis()->SetNdivisions(502,false);
1610 >  ratio->SetFillColor(TColor::GetColor("#58D3F7"));
1611 >  ratio->Draw("e1");
1612 >  
1613 >  if(syshisto!=0) {
1614 >    SysEnvelope->SetFillColor(TColor::GetColor("#FE9A2E"));
1615 >    SysEnvelope->Draw("2,same");
1616 >    ratio->Draw("e1,same");
1617 >  } else {
1618 >    eratio->Draw("1");
1619 >  }
1620 >  ratio->Draw("same,axis");
1621 >  TLine *oneline = new TLine(ratio->GetXaxis()->GetBinLowEdge(1),1,ratio->GetXaxis()->GetBinLowEdge(ratio->GetNbinsX())+ratio->GetXaxis()->GetBinWidth(ratio->GetNbinsX()),1);
1622 >  oneline->SetLineStyle(2);
1623 >  oneline->SetLineColor(kBlue);
1624 >  oneline->Draw("same");
1625 >
1626 >  Ratio_main_canvas->cd();
1627 >  Ratio_main_canvas->Modified();
1628 >  Ratio_main_canvas->cd();
1629 >  Ratio_main_canvas->SetSelected(Ratio_main_canvas);
1630 >  
1631 >  CompleteSave(Ratio_main_canvas,savemeas+"_withSysRatio");
1632 >  bottompad->cd();
1633 >  
1634 >  Double_t chi2;
1635 >  Int_t ndf,igood;
1636 >  Double_t chi2prob = MarcosChi2TestX(nominator,denominator,chi2,ndf,igood,"UW");
1637 >
1638 >  stringstream Chi2text;
1639 >  Chi2text << "#chi^{2} / ndf: " << chi2 << " / " << ndf;
1640 >  stringstream Chi2probtext;
1641 >  Chi2probtext << "#chi^{2} prob: " << chi2prob;
1642 >  TText* chi2box = write_text(0.02,0.11,Chi2text.str());
1643 >  chi2box->SetTextSize(0.08);
1644 >  chi2box->SetTextAlign(11);
1645 >  chi2box->Draw();
1646 >  TText* chi2probbox = write_text(0.02,0.04,Chi2probtext.str());
1647 >  chi2probbox->SetTextSize(0.08);
1648 >  chi2probbox->SetTextAlign(11);
1649 >  chi2probbox->Draw();
1650 >  
1651 >  stringstream CompRatio;
1652 >  CompRatio << "Ratio: " << nominator->Integral() << " / " << denominator->Integral() << " \n " << nominator->Integral()/denominator->Integral();
1653 >  TText *CompleteRatio = write_text(0.98,0.04,CompRatio.str());
1654 >  CompleteRatio->SetTextSize(0.08);
1655 >  CompleteRatio->SetTextAlign(31);
1656 >  CompleteRatio->Draw();
1657 >  
1658 >  CompleteSave(Ratio_main_canvas,savemeas+"_withSysRatio_and_Chi2");
1659 >  delete Ratio_main_canvas;
1660 >  delete ratio;
1661 > }
1662 >
1663 > TH1F* CollapseStack(THStack stack,TString hname="CollapsedStack") {
1664     TH1F *bhist = ((TH1F*)((stack.GetHists())->At(0)));
1665 <   TH1F *basehisto = (TH1F*)bhist->Clone("base");
1665 >   TH1F *basehisto = (TH1F*)bhist->Clone(hname);
1666     TIter next(stack.GetHists());
1667     TH1F *h;
1668     int counter=0;
# Line 1371 | Line 1674 | TH1F* CollapseStack(THStack stack) {
1674     return basehisto;
1675   }
1676  
1677 < void save_with_ratio(TH1F *nominator, THStack denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false) {
1678 <  save_with_ratio(nominator, CollapseStack(denominator), canvas, savemeas, do_bpred_ratio);
1677 > void Save_With_Ratio(TH1F *nominator, THStack denominator, TVirtualPad *canvas, string savemeas, bool do_bpred_ratio=false) {
1678 >  TH1F *denominator_histo = (TH1F*) CollapseStack(denominator);
1679 >  Save_With_Ratio(nominator, denominator_histo, canvas, savemeas, do_bpred_ratio);
1680 >  delete denominator_histo;
1681   }
1682  
1683   void flag_this_change(string function, int line, int checked=0) {
# Line 1388 | Line 1693 | void flag_this_change(string function, i
1693    if(checked==2) write_info(function,peakmodificationwarning.str());
1694   }
1695  
1696 < void write_analysis_type(bool isonpeak) {
1696 > void write_analysis_type(bool isonpeak, bool dobtag) {
1697   //http://www.network-science.de/ascii/, ogre
1698    if(!PlottingSetup::publicmode) {
1699      if(isonpeak) {
1700 + //font: big    
1701          dout << "\033[1;34m" << endl;
1702 <        dout <<  "    //\\\\                                  _                          _           _     " << endl;
1703 <        dout <<  "   //  \\\\     ___  _ __  _ __   ___  __ _| | __     __ _ _ __   __ _| |_   _ ___(_)___" << endl;
1704 <        dout <<  "  //    \\\\   / _ \\| '_ \\| '_ \\ / _ \\/ _` | |/ /    / _` | '_ \\ / _` | | | | / __| / __|" << endl;
1705 <        dout <<  " //      \\\\ | (_) | | | | |_) |  __/ (_| |   <    | (_| | | | | (_| | | |_| \\__ \\ \\__ \\" << endl;
1706 <        dout <<  "//        \\\\ \\___/|_| |_| .__/ \\___|\\__,_|_|\\_\\    \\__,_|_| |_|\\__,_|_|\\__, |___/_|___/" << endl;
1707 <        dout <<  "                        |_|                                            |___/" << endl;
1702 >        dout <<  "                             _ __________                                                         " << endl;
1703 >        dout <<  "                            | |___  /  _ \\                                                       " << endl;
1704 >        dout <<  "       /\\           ___     | |  / /| |_) |                                                       " << endl;
1705 >        dout <<  "      /  \\         /  _|_   | | / / |  _ <                                                         " << endl;
1706 >        dout <<  "     /    \\        | (_| |__| |/ /__| |_) |                                                       " << endl;
1707 >        dout <<  "                   \\___|\\____//_____|____/                                                       " << endl;
1708      } else {
1709 + //font: big    
1710          dout << "\033[1;33m" << endl;
1711 <        dout <<  "                    __  __                  _                          _           _     " << endl;
1712 <        dout <<  "              ___  / _|/ _|_ __   ___  __ _| | __     __ _ _ __   __ _| |_   _ ___(_)___ " << endl;
1713 <        dout <<  "    /\\       / _ \\| |_| |_| '_ \\ / _ \\/ _` | |/ /    / _` | '_ \\ / _` | | | | / __| / __|" << endl;
1714 <        dout <<  " __/  \\__   | (_) |  _|  _| |_) |  __/ (_| |   <    | (_| | | | | (_| | | |_| \\__ \\ \\__ \\" << endl;
1715 <        dout <<  "             \\___/|_| |_| | .__/ \\___|\\__,_|_|\\_\\    \\__,_|_| |_|\\__,_|_|\\__, |___/_|___/" << endl;
1716 <        dout <<  "                          |_|                                            |___/           " << endl;
1711 >        dout <<  "                    _      _ __________                                                             " << endl;
1712 >        dout <<  "                   (_)    | |___  /  _ \\                                                            " << endl;
1713 >        dout <<  "       /\\           _     | |  / /| |_) |                                                           " << endl;
1714 >        dout <<  "    /\\/  \\         | |_   | | / / |  _ <                                                            " << endl;
1715 >        dout <<  "   /      \\        | | |__| |/ /__| |_) |                                                           " << endl;
1716 >        dout <<  "                   |_|\\____//_____|____/                                                            " << endl;
1717      }
1718 <    dout << "\033[0m" << endl;
1718 >  
1719 >    if(dobtag) {
1720 > //font: big    
1721 >        dout <<  "\033[1;32m                  \\     /                                                                " << endl;
1722 >        dout <<  "             \\    o ^ o    /                                                             " << endl;
1723 >        dout <<  "               \\ (     ) /                                                               " << endl;
1724 >        dout <<  "    ____________(%%%%%%%)____________        _          _ __________                     " << endl;
1725 >        dout <<  "   (     /   /  )%%%%%%%(  \\   \\     )      | |        | |___  /  _ \\                    " << endl;
1726 >        dout <<  "   (___/___/__/           \\__\\___\\___)      | |__      | |  / /| |_) |                   " << endl;
1727 >        dout <<  "      (     /  /(%%%%%%%)\\  \\     )         | '_ \\ _   | | / / |  _ <                    " << endl;
1728 >        dout <<  "       (__/___/ (%%%%%%%) \\___\\__)          | |_) | |__| |/ /__| |_) |                   " << endl;
1729 >        dout <<  "               /(       )\\                  |_.__/ \\____//_____|____/                    " << endl;
1730 >        dout <<  "             /   (%%%%%)   \\                                                             " << endl;
1731 >        dout <<  "                  (%%%)                                                                  " << endl;
1732 >        dout <<  "                    !                                                                    " << endl;
1733    }
1734 +  } else {
1735 +    //we're in public! don't advertise what we're up to :-)
1736 +       dout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << endl;
1737 +       dout << "                              Starting the analysis                     " << endl;
1738 +       dout << "                                     i:" << !isonpeak << " , b " << dobtag << endl;
1739 +       dout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" << endl;
1740 +  }
1741 +  dout << "\033[0m" << endl;
1742 +
1743   }
1744  
1745  
# Line 1418 | Line 1748 | vector<string> StringSplit(string str, s
1748          vector<string> results;
1749          
1750          int cutAt;
1751 <        while( (cutAt = str.find_first_of(delim)) != str.npos )
1751 >        while( (cutAt = str.find_first_of(delim)) != (int)str.npos )
1752          {
1753                  if(cutAt > 0)
1754                  {
# Line 1435 | Line 1765 | vector<string> StringSplit(string str, s
1765  
1766   void manually_set_jzb_cuts(vector<float> &jzb_cut,string jzbcut_string) {
1767          vector<string> jzbcutvector = StringSplit(jzbcut_string,",");
1768 <        for(int i=0;i<jzbcutvector.size();i++) {
1768 >        for(int i=0;i<(int)jzbcutvector.size();i++) {
1769                  jzb_cut.push_back(atoi(jzbcutvector[i].c_str()));
1770                  dout << "Added a JZB cut manually at " << atoi(jzbcutvector[i].c_str()) << endl;
1771          }
# Line 1449 | Line 1779 | void process_directory(TString directory
1779          dp = opendir (directory);
1780          if (dp != NULL)
1781          {
1782 <                while (ep = readdir (dp))
1782 >                while ((ep = readdir (dp)))
1783                  {
1784                          string filename=(string)ep->d_name;
1785 <                        if(filename.find(".root")!=-1)
1785 >                        if((int)filename.find(".root")!=-1)
1786                          {
1787                                  files.push_back(string(directory)+filename);
1788                          }
# Line 1463 | Line 1793 | void process_directory(TString directory
1793                  perror ("Couldn't open the directory");
1794   }
1795  
1796 + void ClearHisto(TH1* histo) {
1797 +    for(int ix=0;ix<=histo->GetNbinsX()+1;ix++) {//also reset over- and underflow
1798 +        histo->SetBinContent(ix,0);
1799 +        histo->SetBinContent(ix,0);
1800 +    }
1801 + }
1802 +
1803 + void AddSquared(TH1F *histo, TH1F *addthis, float weight) {
1804 +    for(int ix=1;ix<=addthis->GetNbinsX()+1;ix++) {//also include over- and underflow
1805 +        float bincenter=addthis->GetBinCenter(ix);
1806 +        int nBinHisto= histo->FindBin(bincenter);
1807 +        float old_content=histo->GetBinContent(nBinHisto);
1808 +        histo->SetBinContent(nBinHisto,old_content+((addthis->GetBinContent(ix))*(addthis->GetBinContent(ix))*weight));
1809 +    }
1810 + }
1811 +
1812 + void SQRT(TH1* histo) {
1813 +    for(int ix=0;ix<=histo->GetNbinsX()+1;ix++) {//also include over- and underflow
1814 +        histo->SetBinContent(ix,TMath::Sqrt(histo->GetBinContent(ix)));
1815 +    }
1816 + }
1817 +
1818   /*string GetCompleteHostname() {
1819    struct addrinfo hints, *info, *p;
1820    int gai_result;
# Line 1494 | Line 1846 | const char* concatenate (const char* a,
1846    return bla.str().c_str();
1847   }
1848  
1849 +
1850 + string ReplaceCharacter(string originalstring,string replacethis,string replacewiththis)
1851 + {
1852 +        int pos = originalstring.find(replacethis);
1853 +        while(pos != -1) {
1854 +          originalstring.replace(pos, replacewiththis.length(), replacewiththis);
1855 +          pos = originalstring.find(replacethis);
1856 +        }
1857 +        return originalstring;
1858 + }
1859 +
1860 + string removefunnystring(string name) {
1861 +  name=ReplaceCharacter(name,"[","_");
1862 +  name=ReplaceCharacter(name,"]","_");
1863 +  name=ReplaceCharacter(name,"{","_");
1864 +  name=ReplaceCharacter(name,"}","_");
1865 +  name=ReplaceCharacter(name,"(","_");
1866 +  name=ReplaceCharacter(name,")","_");
1867 +  name=ReplaceCharacter(name," ","_");
1868 +  name=ReplaceCharacter(name,".","_");
1869 +  name=ReplaceCharacter(name,",","_");
1870 +  name=ReplaceCharacter(name,";","_");
1871 +  name=ReplaceCharacter(name,":","_");
1872 +  name=ReplaceCharacter(name,"'","_");
1873 +  name=ReplaceCharacter(name,"$","_");
1874 +  name=ReplaceCharacter(name,"@","_");
1875 +  name=ReplaceCharacter(name,"#","_");
1876 +  name=ReplaceCharacter(name,"+","_");
1877 +  name=ReplaceCharacter(name,"-","_");
1878 +  name=ReplaceCharacter(name,"/","_");
1879 +  return name;
1880 + }
1881 +
1882 + TPolyLine* GetFitUncertaintyShape(TF1 *fit, float low, float high,float x0,float x1) {
1883 +  int nPoints=1000;
1884 +  double x[nPoints+1];
1885 +  double y[nPoints];
1886 +  
1887 +  float par1=fit->GetParameter(0);
1888 +  float dpar1=fit->GetParError(0);
1889 +  float par2=fit->GetParameter(1);
1890 +  float dpar2=fit->GetParError(1);
1891 +    
1892 +  float step=(x1-x0)/(nPoints/2.0 -1);
1893 +  for(int i=0;i<nPoints/2.0;i++) {
1894 +    x[i]=x0+i*step;
1895 +    y[i]=(par1+x[i]*par2)    -    sqrt(dpar1*dpar1+dpar2*dpar2*x[i]*x[i]);
1896 +    x[nPoints-1-i]=x[i];
1897 +    y[nPoints-1-i]=(par1+x[i]*par2)    +    sqrt(dpar1*dpar1+dpar2*dpar2*x[i]*x[i]);
1898 +    if(y[i]<low) y[i]=low;
1899 +    if(y[i]>high) y[i]=high;
1900 +    if(y[nPoints-1-i]<low) y[nPoints-1-i]=low;
1901 +    if(y[nPoints-1-i]>high) y[nPoints-1-i]=high;
1902 +  }
1903 +  
1904 +  x[nPoints]=x[0];
1905 +  y[nPoints]=y[0];
1906 +  
1907 +  TPolyLine *l = new TPolyLine(nPoints+1,x,y);
1908 +  l->SetFillColor(TColor::GetColor("#5858FA"));
1909 +  l->SetLineColor(TColor::GetColor("#5858FA"));
1910 +  l->SetLineWidth(1);
1911 +  return l;
1912 + }  
1913 +
1914 + //code for ReplaceAll copied from
1915 + //http://stackoverflow.com/questions/5343190/c-how-do-i-replace-all-instances-of-of-a-string-with-another-string
1916 + string ReplaceAll(string str, string from, string to) {
1917 +    size_t start_pos = 0;
1918 +    while((start_pos = str.find(from, start_pos)) != std::string::npos) {
1919 +         str.replace(start_pos, from.length(), to);
1920 +         start_pos += to.length(); // ...
1921 +    }
1922 +    return str;
1923 + }
1924 +
1925 + TPolyLine* GetFitUncertaintyShape(TH1 *histo, string fitname, float low, float high) {
1926 +  TF1 *fit = (TF1*)histo->GetFunction(fitname.c_str());
1927 +  if(!fit) {
1928 +    cout << "CANNOT PRODUCE FIT UNCERTAINTY AS THERE IS NO FIT CALLED " << fitname << " IN " << histo->GetName() << endl;
1929 +    return NULL;
1930 +  }
1931 +  float x0=histo->GetBinLowEdge(1);
1932 +  float x1=histo->GetBinLowEdge(histo->GetNbinsX())+histo->GetBinWidth(histo->GetNbinsX());
1933 +  return GetFitUncertaintyShape(fit, low, high,x0,x1);
1934 + }
1935 +
1936 + TPolyLine* GetFitUncertaintyShape(TGraphErrors *gr, string fitname, float low, float high, float minx=-999, float maxx=-999) {
1937 +  TF1 *fit = (TF1*)gr->GetFunction(fitname.c_str());
1938 +  if(!fit) {
1939 +    cout << "CANNOT PRODUCE FIT UNCERTAINTY AS THERE IS NO FIT CALLED " << fitname << " IN " << gr->GetName() << endl;
1940 +    return NULL;
1941 +  }
1942 +  double x,y;
1943 +  gr->GetPoint(0,x,y);
1944 +  float min=x;
1945 +  gr->GetPoint(gr->GetN()-1,x,y);
1946 +  float max=x;
1947 +  if(minx>-999 && maxx>-999 && minx<maxx) {
1948 +    max=maxx;
1949 +    min=minx;
1950 +  }
1951 +  return GetFitUncertaintyShape(fit, low, high,min,max);
1952 + }
1953 +
1954   stringstream all_bugs;
1955  
1956   void bug_tracker(string function, int line, string description) {
# Line 1530 | Line 1987 | void bug_tracker(string function, int li
1987  
1988   //TODO: Write a bug summary at the end.
1989  
1990 + float pSTAR(float mglu, float mlsp, float mchi2) {
1991 +        float mz=91.0;
1992 +        float res=((mchi2*mchi2)-(mlsp+mz)*(mlsp+mz));
1993 +        res*=((mchi2*mchi2)-(mlsp-mz)*(mlsp-mz));
1994 +        res=TMath::Sqrt(res)/(2*mchi2);
1995 +        return res;
1996 + }
1997 +
1998 + float generalizedpSTAR(float massmother, float massdaughter1, float massdaughter2) {
1999 +        float res=((massmother*massmother)-(massdaughter1+massdaughter2)*(massdaughter1+massdaughter2));
2000 +        res*=((massmother*massmother)-(massdaughter1-massdaughter2)*(massdaughter1-massdaughter2));
2001 +        res=TMath::Sqrt(res)/(2*massmother);
2002 +        return res;
2003 + }
2004 +
2005 + void FindMinMax(TGraphErrors *gr, float &min, float &max) {
2006 +  double x,y;
2007 +  gr->GetPoint(0,x,y);
2008 +  min=500*y;
2009 +  max=0.0002*y;
2010 +  for(int i=0;i<gr->GetN();i++) {
2011 +    gr->GetPoint(i,x,y);
2012 +    cout << "Point at " << x << "/" << y << " with error " << gr->GetErrorY(i) << endl;
2013 +    float potmin=y-1.5*gr->GetErrorY(i); // taking 1.5 times the error to allow for enough space to present the plot
2014 +    float potmax=y+1.5*gr->GetErrorY(i); // taking 1.5 times the error to allow for enough space to present the plot
2015 +    if(potmin<min) min=potmin;
2016 +    if(potmax>max) max=potmax;
2017 +  }
2018 + }
2019 +
2020 +
2021 + void CleanLegends() {
2022 +  for(int ihist=(int)PlottingSetup::FakeHistoHeap.size()-1;ihist>=0;ihist--) {
2023 +    PlottingSetup::FakeHistoHeap[ihist]->Delete();
2024 +    PlottingSetup::FakeHistoHeap.pop_back();
2025 +  }
2026 + }
2027 +
2028 + string DigitsAfterComma(float number, int digits) {
2029 +  float temp=number*pow(10.0,digits);
2030 +  temp=int(temp)/pow(10.0,digits);
2031 +  return any2string(temp);
2032 + }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines