ra4_draw  4bd0201e3d922d42bd545d4b045ed44db33454a4
plot_limit.cxx
Go to the documentation of this file.
1 #include <fstream>
2 #include <sstream>
3 #include <string>
4 #include <algorithm>
5 
6 #include <unistd.h>
7 #include <getopt.h>
8 
9 #include "TCanvas.h"
10 #include "TGraph.h"
11 #include "TGraphAsymmErrors.h"
12 #include "TH1D.h"
13 #include "TStyle.h"
14 #include "TLegend.h"
15 #include "TLatex.h"
16 #include "TLine.h"
17 
18 #include "core/styles.hpp"
19 #include "core/utilities.hpp"
20 #include "core/plot_opt.hpp"
21 
22 using namespace std;
23 
24 namespace{
25  TString lumi = "40";
26  TString filename = "txt/limits/limits_TChiHH_lumi"+lumi+".txt";
27  TString model = "TChiHH";
28 }
29 
30 void GetOptions(int argc, char *argv[]);
31 
32 int main(int argc, char *argv[]){
33  GetOptions(argc, argv);
34 
36  PlotOpt opts("txt/plot_styles.txt", "Ratio");
37  setPlotStyle(opts);
38  gStyle->SetGridStyle(3);
39 
40  if(filename == "") ERROR("No input file provided");
41  ifstream infile(filename);
42 
43  vector<double> vmx, vmy, vxsec, vexsec, vobs, vobsup, vobsdown;
44  vector<double> vexp, vup, vdown, v2up, v2down, vsigobs, vsigexp, zeroes, ones;
45  double maxy=-99., miny=1e99;
46 
47  string line_s;
48  while(getline(infile, line_s)){
49  istringstream iss(line_s);
50  double pmx, pmy, pxsec, pexsec, pobs, pobsup, pobsdown, pexp, pup, pdown, p2up, p2down, sigobs, sigexp;
51  iss >> pmx >> pmy >> pxsec >> pexsec >> pobs >> pobsup >> pobsdown
52  >> pexp >> pup >> pdown >> p2up >> p2down >> sigobs >> sigexp;
53  if(pmx<200) continue;
54  vmx.push_back(pmx);
55  vmy.push_back(pmy);
56  vxsec.push_back(pxsec);
57  vexsec.push_back(pexsec);
58  vobs.push_back(pobs);
59  vobsup.push_back(pobsup);
60  vobsdown.push_back(pobsdown);
61  vexp.push_back(pexp);
62  vup.push_back(pup-pexp);
63  vdown.push_back(pexp-pdown);
64  v2up.push_back(p2up-pexp);
65  v2down.push_back(pexp-p2down);
66  vsigobs.push_back(sigobs);
67  vsigexp.push_back(sigexp);
68  zeroes.push_back(0);
69  ones.push_back(1);
70  if(miny > min(vobs.back(), 1.)) miny = min(vobs.back(), 1.);
71  if(maxy < max(vobs.back(), 1.)) maxy = max(vobs.back(), 1.);
72  }
73  infile.close();
74 
75  if(vmx.size() <= 0) ERROR("Need at least 1 model to draw limits");
76  if(vmx.size() != vmy.size()
77  || vmx.size() != vxsec.size()
78  || vmx.size() != vexsec.size()
79  || vmx.size() != vobs.size()
80  || vmx.size() != vobsup.size()
81  || vmx.size() != vobsdown.size()
82  || vmx.size() != vexp.size()
83  || vmx.size() != vup.size()
84  || vmx.size() != vdown.size()
85  || vmx.size() != v2up.size()
86  || vmx.size() != v2down.size()
87  || vmx.size() != vsigobs.size()
88  || vmx.size() != vsigexp.size()) ERROR("Error parsing text file. Model point not fully specified");
89 
90  // Sorting vectors
91  vector<size_t> perm = SortPermutation(vmx);
92  vmx = ApplyPermutation(vmx , perm);
93  vmy = ApplyPermutation(vmy , perm);
94  vxsec = ApplyPermutation(vxsec , perm);
95  vexsec = ApplyPermutation(vexsec , perm);
96  vobs = ApplyPermutation(vobs , perm);
97  vobsup = ApplyPermutation(vobsup , perm);
98  vobsdown = ApplyPermutation(vobsdown , perm);
99  vexp = ApplyPermutation(vexp , perm);
100  vup = ApplyPermutation(vup , perm);
101  vdown = ApplyPermutation(vdown , perm);
102  v2up = ApplyPermutation(v2up , perm);
103  v2down = ApplyPermutation(v2down , perm);
104  vsigobs = ApplyPermutation(vsigobs , perm);
105  vsigexp = ApplyPermutation(vsigexp , perm);
106 
107  TCanvas can;
108  //can.SetGrid();
109  can.SetFillStyle(4000);
110  TString chi1n = "#lower[-0.12]{#tilde{#chi}}#lower[0.2]{#scale[0.85]{^{0}}}#kern[-1.3]{#scale[0.85]{_{1}}}";
111  TString chi2n = "#lower[-0.12]{#tilde{#chi}}#lower[0.2]{#scale[0.85]{^{0}}}#kern[-1.]{#scale[0.85]{_{2}}}";
112  TString chi1pm= "#lower[-0.12]{#tilde{#chi}}#lower[0.2]{#scale[0.85]{^{#pm}}}#kern[-1.3]{#scale[0.85]{_{1}}}";
113  TString chii= "#lower[-0.12]{#tilde{#chi}}#lower[0.2]{#scale[0.85]{^{0,#pm}}}#kern[-3.]{#scale[0.85]{_{i}}}";
114  TString chij= "#lower[-0.12]{#tilde{#chi}}#lower[0.2]{#scale[0.85]{^{0,#mp}}}#kern[-3.]{#scale[0.85]{_{j}}}";
115  TString mass_ = "m#kern[0.1]{#lower[-0.12]{_{";
116  float minh=200, maxh=1000;
117  TH1D histo("histo", "", 18, minh, maxh);
118  histo.SetMinimum(0);
119  histo.SetMaximum(4.5);
120  histo.GetYaxis()->CenterTitle(true);
121  histo.GetXaxis()->SetLabelOffset(0.02);
122  histo.SetXTitle("Higgsino mass "+mass_+chi1n+"}}} [GeV]");
123  histo.SetYTitle("#sigma_{excl}^{95% CL}/#sigma_{theor}");
124  histo.Draw();
125 
126  TLine line;
127  line.SetLineColor(4); line.SetLineStyle(2); line.SetLineWidth(4);
128  TLatex cmslabel;
129 
130  cmslabel.SetNDC(kTRUE);
131 
132  TGraphAsymmErrors grexp2(vmx.size(), &(vmx[0]), &(vexp[0]), &(zeroes[0]), &(zeroes[0]), &(v2down[0]), &(v2up[0]));
133  grexp2.SetLineColor(1); grexp2.SetFillColor(5); grexp2.SetLineWidth(3); grexp2.SetLineStyle(2);
134  grexp2.Draw("e3 same");
135  TGraphAsymmErrors grexp1(vmx.size(), &(vmx[0]), &(vexp[0]), &(zeroes[0]), &(zeroes[0]), &(vdown[0]), &(vup[0]));
136  grexp1.SetLineColor(1); grexp1.SetFillColor(3); grexp1.SetLineWidth(3); grexp1.SetLineStyle(2);
137  grexp1.Draw("e3 same");
138  TGraph grexp(vmx.size(), &(vmx[0]), &(vexp[0]));
139  grexp.SetLineWidth(3); grexp.SetLineStyle(2);
140  grexp.Draw("same");
141  TGraph grobs(vmx.size(), &(vmx[0]), &(vobs[0]));
142  grobs.SetLineWidth(1);
143  grobs.Draw("same");
144 
146  TString cmsPrel = "#font[62]{CMS} #scale[0.8]{#font[52]{Simulation}}";
147  TString lumiEner = "#font[42]{"+lumi+" fb^{-1} (13 TeV)}"; lumiEner.ReplaceAll("p",".");
148  TString ppChiChi = "pp #rightarrow "+chii+"#kern[0.6]{"+chij+"} #rightarrow hh#tilde{G}#tilde{G}";
149  TString mChis = mass_+chi2n+"}}} #approx "+mass_+chi1pm+"}}} #approx "+mass_+chi1n+"}}}, "+mass_+"#tilde{G}}}} = 1 GeV";
150  cmslabel.SetTextAlign(11); cmslabel.SetTextSize(0.06);
151  cmslabel.DrawLatex(opts.LeftMargin()+0.005, 1-opts.TopMargin()+0.015, cmsPrel);
152  cmslabel.SetTextAlign(31); cmslabel.SetTextSize(0.056);
153  cmslabel.DrawLatex(1-opts.RightMargin()-0.005, 1-opts.TopMargin()+0.015, lumiEner);
154  line.DrawLine(minh, 1, maxh, 1);
156  cmslabel.SetTextAlign(31); cmslabel.SetTextSize(0.045);
157  cmslabel.SetTextFont(132);
158  cmslabel.DrawLatex(1-opts.RightMargin()-0.023, opts.BottomMargin()+0.09, ppChiChi);
159  cmslabel.DrawLatex(1-opts.RightMargin()-0.023, opts.BottomMargin()+0.04, mChis);
160 
161 
162  double legX(0.54), legY(1-opts.TopMargin()-0.04), legSingle = 0.05;
163  double legW = 0.26, legH = legSingle*4;
164  TLegend leg(legX-legW, legY-legH, legX, legY);
165  leg.SetTextSize(0.04); leg.SetFillColor(0);
166  leg.SetFillStyle(0); leg.SetBorderSize(0);
167  leg.AddEntry(&line, "NLO+NLL", "l");
168  leg.AddEntry(&grobs, "\"Observed\"", "l");
169  leg.AddEntry(&grexp1, "Expected #pm #sigma");
170  leg.AddEntry(&grexp2, "Expected #pm 2#sigma");
171  leg.Draw();
172 
173  histo.Draw("axis same");
174  can.SaveAs("plots/higgsino_limits_lumi"+lumi+".pdf");
175 
176  // for(size_t i = 0; i < vxsec.size(); ++i)
177  // cout<<vmx[i]<<" -> "<<vexp[i]<<"+"<<vup[i]<<"++"<<v2up[i]<<" -"<<vdown[i]<<"--"<<v2down[i]<<endl;
178 
182  maxy=-99.; miny=1e99;
183  for(size_t i = 0; i < vxsec.size(); ++i){
184  vxsec[i] *= 1000; // Converting it to fb
185  vexsec[i] *= vxsec[i];
186  vobs[i] *= vxsec[i];
187  vobsup[i] *= vxsec[i];
188  vobsdown[i]*= vxsec[i];
189  vexp[i] *= vxsec[i];
190  vup[i] *= vxsec[i];
191  vdown[i] *= vxsec[i];
192  v2up[i] *= vxsec[i];
193  v2down[i] *= vxsec[i];
194  if(miny > min(vexp[i]-v2down[i], vxsec[i])) miny = min(vexp[i]-v2down[i], vxsec[i]);
195  if(maxy < max(vexp[i]+v2up[i], vxsec[i])) maxy = max(vexp[i]+v2up[i], vxsec[i]);
196  //cout<<vmx[i]<<" -> "<<vobs[i]<<endl;
197  }
198 
199  histo.GetXaxis()->SetLabelOffset(0.01);
200  histo.SetMinimum(miny/2.);
201  histo.SetMaximum(maxy*1.2);
202  histo.SetYTitle("#sigma_{excl}^{95% CL} #times BF(hh #rightarrow bbbb) [fb]");
203  histo.Draw();
204  TGraphAsymmErrors gexp2(vmx.size(), &(vmx[0]), &(vexp[0]), &(zeroes[0]), &(zeroes[0]), &(v2down[0]), &(v2up[0]));
205  gexp2.SetLineColor(1); gexp2.SetFillColor(5); gexp2.SetLineWidth(3); gexp2.SetLineStyle(2);
206  gexp2.Draw("e3 same");
207  TGraphAsymmErrors gexp1(vmx.size(), &(vmx[0]), &(vexp[0]), &(zeroes[0]), &(zeroes[0]), &(vdown[0]), &(vup[0]));
208  gexp1.SetLineColor(1); gexp1.SetFillColor(3); gexp1.SetLineWidth(3); gexp1.SetLineStyle(2);
209  gexp1.Draw("e3 same");
210  TGraph gexp(vmx.size(), &(vmx[0]), &(vexp[0]));
211  gexp.SetLineWidth(3); gexp.SetLineStyle(2);
212  gexp.Draw("same");
213  TGraph gobs(vmx.size(), &(vmx[0]), &(vobs[0]));
214  gobs.SetLineWidth(1);
215  gobs.Draw("same");
216  TGraph gxsec(vmx.size(), &(vmx[0]), &(vxsec[0]));
217  gxsec.SetLineWidth(4); gxsec.SetLineColor(4); gxsec.SetLineStyle(2);
218  gxsec.Draw("same");
219 
220  can.SetLogy(true);
221 
222  legX = 1-opts.RightMargin()-0.02;
223  leg.SetX1NDC(legX-legW); leg.SetX2NDC(legX);
224  leg.Draw();
225 
227  cmslabel.SetTextAlign(11); cmslabel.SetTextSize(0.06);
228  cmslabel.DrawLatex(opts.LeftMargin()+0.005, 1-opts.TopMargin()+0.015, cmsPrel);
229  cmslabel.SetTextAlign(31); cmslabel.SetTextSize(0.056);
230  cmslabel.DrawLatex(1-opts.RightMargin()-0.005, 1-opts.TopMargin()+0.015, lumiEner);
231 
233  cmslabel.SetTextAlign(11); cmslabel.SetTextSize(0.045);
234  cmslabel.SetTextFont(132);
235  cmslabel.DrawLatex(opts.LeftMargin()+0.03, opts.BottomMargin()+0.09, ppChiChi);
236  cmslabel.DrawLatex(opts.LeftMargin()+0.03, opts.BottomMargin()+0.04, mChis);
237 
238 
239  histo.Draw("axis same");
240  can.SaveAs("plots/higgsino_limits_fb_lumi"+lumi+".pdf");
241 
245  can.SetLogy(false);
246  histo.GetXaxis()->SetLabelOffset(0.02);
247  histo.SetMinimum(0);
248  histo.SetMaximum(3.5);
249  if(lumi=="40") histo.SetMaximum(4.5);
250  histo.SetYTitle(" Expected discovery significance [#sigma]");
251  histo.Draw();
252 
253  TGraph gsig(vmx.size(), &(vmx[0]), &(vsigexp[0]));
254  gsig.SetLineWidth(3); gsig.SetLineColor(4);
255  gsig.Draw("same");
257  cmslabel.SetTextAlign(11); cmslabel.SetTextSize(0.06);
258  cmslabel.DrawLatex(opts.LeftMargin()+0.005, 1-opts.TopMargin()+0.015, cmsPrel);
259  cmslabel.SetTextAlign(31); cmslabel.SetTextSize(0.056);
260  cmslabel.DrawLatex(1-opts.RightMargin()-0.005, 1-opts.TopMargin()+0.015, lumiEner);
262  cmslabel.SetTextAlign(33); cmslabel.SetTextSize(0.045);
263  cmslabel.SetTextFont(132);
264  cmslabel.DrawLatex(1-opts.RightMargin()-0.023, 1-opts.TopMargin()-0.025, ppChiChi);
265  cmslabel.DrawLatex(1-opts.RightMargin()-0.023, 1-opts.TopMargin()-0.09, mChis);
266 
267  can.SaveAs("plots/higgsino_significance_lumi"+lumi+".pdf");
268 
269 
270 }
271 
272 void GetOptions(int argc, char *argv[]){
273  while(true){
274  static struct option long_options[] = {
275  {"model", required_argument, 0, 'm'},
276  {"file", required_argument, 0, 'f'},
277  {0, 0, 0, 0}
278  };
279 
280  char opt = -1;
281  int option_index;
282  opt = getopt_long(argc, argv, "f:m:", long_options, &option_index);
283  if( opt == -1) break;
284 
285  string optname;
286  switch(opt){
287  case 'm':
288  model = optarg;
289  break;
290  case 'f':
291  filename = optarg;
292  break;
293  case 0:
294  optname = long_options[option_index].name;
295  if(optname == ""){
296  printf("Bad option! Found option name %s\n", optname.c_str());
297  }else{
298  printf("Bad option! Found option name %s\n", optname.c_str());
299  }
300  break;
301  default:
302  printf("Bad option! getopt_long returned character code 0%o\n", opt);
303  break;
304  }
305  }
306 }
PlotOpt & TopMargin(double top)
Definition: plot_opt.cpp:260
void setPlotStyle(PlotOpt opts)
Definition: styles.cpp:7
PlotOpt & LeftMargin(double left)
Definition: plot_opt.cpp:233
STL namespace.
std::vector< T > ApplyPermutation(const std::vector< T > &vec, const std::vector< std::size_t > &perm)
Definition: utilities.hpp:47
#define ERROR(x)
Definition: utilities.hpp:17
int main(int argc, char *argv[])
Definition: plot_limit.cxx:32
PlotOpt & BottomMargin(double bottom)
Definition: plot_opt.cpp:251
void GetOptions(int argc, char *argv[])
Definition: plot_limit.cxx:272
std::vector< std::size_t > SortPermutation(const std::vector< T > &vec)
Definition: utilities.hpp:39
PlotOpt & RightMargin(double right)
Definition: plot_opt.cpp:242