DrawUtils.h
Go to the documentation of this file.
1 
2 #include "TH1.h"
3 #include "TH2.h"
4 #include "TMarker.h"
5 #include "TProfile.h"
6 #include "TProfile2D.h"
7 #include "TLine.h"
8 #include "TLegend.h"
9 #include "TCanvas.h"
10 #include "TGraph.h"
11 #include "TGraphErrors.h"
12 #include "TText.h"
13 
14 
15 typedef struct Plot Plot;
16 struct Plot {
19  TH1* h;
22  int color;
23  Plot(){
24  h=0;
25  }
26  Plot(TH1* h_){
27  h=h_;
28  if(h)
29  type=h->ClassName();
30  }
31  Plot(TH1* h_, std::string name_,
32  std::string legLabel_){
33  h=h_;
34  name=name_;
35  legLabel=legLabel_;
36  if(h)
37  type=h->ClassName();
38  }
39  Plot(std::string file_,
40  std::string name_,
41  std::string legLabel_){
42  file=file_;
43  name=name_;
44  legLabel=legLabel_;
45  TFile * f = TFile::Open(file.c_str());
46  f->GetObject(name.c_str(),h);
47  if(!h){
48  std::cout << "did not find " << name
49  << " in " << f->GetName()
50  << std::endl;
51  std::abort();
52  } else {
53  h->SetDirectory(0);
54  type=h->ClassName();
55  if(h->GetEntries()==0){
56  std::cout << "Zero entries in " << name
57  << " from " << f->GetName()
58  << std::endl;
59  //std::abort(); // occasionally okay to have zero entries if not ratio denom
60  }
61  }
62  f->Close();
63  }
64 };
65 
66 
67 class DrawUtils{
68 
69 public:
70 
71  DrawUtils();
72  ~DrawUtils();
73 
74  // Stacks plot_vec colored with "plc pmc" onto canv.
75  // If hDenom is specified, it is drawn first, grey and dotted.
76  // If ratio_title is provided, the ratios of
77  // plot_vec over hDenom are drawn with the same colors.
78  // Instead of drawing a legend on canv, potentially obscuring info,
79  // draw a legend on canvleg for the user to place manually
80  // Plot is a class with the TH1*, label, etc..
81  void DrawPlotStack( TCanvas *&canv, TCanvas *&canvleg,
82  const std::vector< Plot > &plot_vec,
83  const Plot &plot_denom,
84  std::string ratio_title=""); // default empty
85  void DrawPlotStack( TCanvas *&canv, TCanvas *&canvleg,
86  const std::vector< Plot > &plot_vec );
87 
88  void MakeRatio( TPad *&pad,
89  const std::vector<Plot> &plot_vec,
90  const Plot &plot_denom,
91  std::string ratio_name );
92  void SetRatioErrors(TGraphErrors* &gRat, TH1* hNum, TH1* hDenom);
93 
94  void Draw2DPlot( TCanvas *&canv, TH2* h2 );
95 
96  void MakeLegend( TCanvas *&cleg, const std::vector< Plot > &plots, const Plot &denom);
97 
98  void GetBestYRange( float & ymin, float & ymax,
99  std::vector< TProfile * > profiles,
100  float maxpad,
101  int rebin = 0,
102  float outlierCut = 0.);
103  void GetBestYRange( float & ymin, float & ymax,
104  const std::vector< Plot > &plots,
105  float maxpad, int rebin = 0, float outlierCut = 0.);
106 
107  float GetLastYMin(){return fLastYMin;}
108  float GetLastYMax(){return fLastYMax;}
109 
110 private:
111 
112  int Cnum(){return fCnum; fCnum++;}
113  int Pnum(){return fPnum; fPnum++;}
114 
115  int fCnum; //unique canvas number
116  int fPnum; //unique plot nummber
117 
118  float fLastYMin;
119  float fLastYMax;
120 
121 };
122 
123 
124 //^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^
126  : fCnum(0)
127  , fPnum(0)
128 {
129 }
130 //^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^
132 {
133 }
134 
135 
136 //^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^
137 void DrawUtils::DrawPlotStack( TCanvas *&c, TCanvas *&cleg,
138  const std::vector< Plot > &plot_vec )
139 {
140  DrawPlotStack( c, cleg, plot_vec, Plot(0,"null","null plot"), "" );
141 }
142 
143 
144 //^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^
145 void DrawUtils::DrawPlotStack( TCanvas *&c, TCanvas *&cleg,
146  const std::vector< Plot > &plot_vec,
147  const Plot &denom,
148  std::string ratio_title )
149 {
150  // Make legend seperately, placing on plots automatically
151  // potentially covers up something important
152  MakeLegend(cleg,plot_vec, denom);
153  c->Clear(); c->cd();
154  gStyle->SetOptStat(0);
155 
156  // temporary until debugged for TH1F
157  if(denom.type!="TProfile") return;
158 
159  std::string plotname;
160  bool isProf;
161  if(plot_vec.size()==0){
162  std::cout << "DrawUtils::Draw Requires at least one plot" << std::endl;
163  return;
164  } else {
165  plotname = plot_vec[0].name;
166  isProf = plot_vec[0].type=="TProfile";
167  if( denom.h &&
168  plot_vec[0].type != denom.type ){
169  std::cout << "Denom type: "<< plot_vec[0].type
170  <<" not the same as plot vector: "
171  << denom.type << std::endl;
172  return;
173  }
174  }
175 
176 
177  // If ratio left blank, still draw denom_plot as
178  // grey in the main plot, but don't draw ratio
179  bool drawRatio = (ratio_title!="");
180 
181 
182  // Prepare pads for optional ratio plots
183  TPad *p1 = new TPad("p1","p1", 0, 0, 1, 1);
184  TPad *p2 = new TPad("p2","p2", 0, 0, 1, 1);
185  p1->SetBottomMargin(0.45);
186  p2->SetTopMargin(0.55);
187  p1->SetFillStyle(0); p2->SetFillStyle(0);
188  if(drawRatio) p1->cd();
189  else c->cd();
190 
191 
192  // AutoZooming on Y
193  float ymin(0), ymax(0);
194  int rebin(0);
195  float outlierCut(0.);
196  if( plotname.find("_time_") != std::string::npos ){
197  rebin = 0; outlierCut = 10.; }
198  // TODO: Make GetBestYRange work with rebin, change line above to use it
199  std::vector< Plot > all_plot_vec = plot_vec;
200  all_plot_vec.push_back(denom);
201  if( isProf ) GetBestYRange( ymin, ymax,
202  all_plot_vec,
203  0.1, rebin, outlierCut );
204  //bool autozoom = false;
205  bool autozoom = ( ymin != std::numeric_limits<float>::max() &&
207 
208 
209  // Should Eventually remove name-specific tweaks like this
210  if( plotname.find("_time_") != std::string::npos ){
211  ymin = ymin - 0.2*(ymax - ymin);
212  fLastYMin = ymin;
213  }
214 
215  // 1d plots: Overlay all samples for this plot
216  // If "denom" sample is specified, draw it first
217  // Grey dotted line for denom sample (MC or control sample or...)
218  bool firstdraw=true;
219  if(denom.h){
220  if(drawRatio) denom.h->GetXaxis()->SetLabelSize(0.);
221  denom.h->SetLineWidth(1);
222  denom.h->SetLineColor(kGray+2);
223  denom.h->SetMarkerColor(kGray+2);
224  if(autozoom){
225  denom.h->SetMinimum(ymin);
226  denom.h->SetMaximum(ymax);
227  }
228  denom.h->Draw( Form( "%s", isProf ? "" : "hist" ));
229  firstdraw=false;
230  }
231 
232 
233  // Now the rest of the plots
234  for( auto& plot: plot_vec ){
235  if( !plot.h ){
236  std::cout << "Plot " << plot.name
237  << " is null, skipping" << std::endl;
238  return;
239  }
240  if( plot.h==denom.h ) continue; // already drew
241 
242  if( isProf ){
243  plot.h->SetLineWidth(2);
244  //if(rebin) plot.h->Rebin(rebin);
245  if(autozoom){
246  plot.h->SetMinimum(ymin);
247  plot.h->SetMaximum(ymax);
248  }
249  } else {
250  plot.h->Scale(1/plot.h->Integral());
251  plot.h->GetYaxis()->SetTitle("Arb. Units");
252  }
253  plot.h->SetStats(0);
254 
255  //if(!isProf) PrintUnderOverFlow(h,plotname);
256 
257  if(firstdraw){
258  if(drawRatio) plot.h->GetXaxis()->SetLabelSize(0.);
259  plot.h->Draw(Form("%s plc pmc", isProf ? "" : "hist"));
260  firstdraw=false;
261  } else {
262  plot.h->Draw(Form("%s same plc pmc", isProf ? "" : "hist"));
263  }
264 
265  } // plot_vec loop
266 
267 
268  if(drawRatio){
269  MakeRatio( p2, plot_vec, denom, ratio_title );
270  c->cd();
271  p1->Draw();
272  p2->Draw();
273  }
274 
275 
276 } // DrawPlotStack
277 
278 
279 //^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^
280 void DrawUtils::Draw2DPlot( TCanvas *&canv,
281  TH2* h2 )
282 {
283  std::string plotname;
284  if(!h2->InheritsFrom(TH2::Class())){
285  std::cout << "not a TH2: " << plotname << std::endl; return; }
286 
287  canv->cd();
288  h2->SetStats(0);
289  gStyle->SetPalette(kBird);
290  h2->Draw("colz");
291  //gROOT->SetStyle("novaStyle");
292 }
293 
294 
295 
296 
297 //^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^
298 void DrawUtils::MakeRatio( TPad *&pad,
299  const std::vector<Plot> &plot_vec,
300  const Plot &plot_denom,
301  std::string ratio_title )
302 {
303  // Build vector of ratios
304  std::vector<TProfile*> ratio_vec;
305  if(plot_denom.h->InheritsFrom(TProfile::Class()))
306  for(auto& pl : plot_vec)
307  ratio_vec.push_back((TProfile*)pl.h->Clone(Form("%s_%i",pl.h->GetName(),Pnum())));
308  //if(plot_denom.h->InheritsFrom(TH1F::Class()))
309  // for(auto& pl : plot_vec)
310  // ratio_vec.push_back((TH1F*)pl.h->Clone(Form("%s_%i",pl.h->GetName(),Pnum())));
311  if( plot_denom.type == plot_vec[0].type ){
312  for(auto& rat : ratio_vec)
313  rat->Divide(plot_denom.h);
314  } else {
315  std::cout << "Denom different type, cannot take ratio" << std::endl;
316  return;
317  }
318 
319  pad->cd();
320  TH1 *hAxes = (TH1*)plot_vec[0].h->Clone();
321  hAxes->Draw("AXIS");
322  hAxes->GetYaxis()->SetTitle(ratio_title.c_str());
323 
324  for( unsigned int i=0; i<ratio_vec.size(); i++ ){
325 
326  TGraphErrors *gRat = new TGraphErrors(ratio_vec[i]);
327  gRat->SetLineColor (plot_vec[i].h->GetLineColor());
328  gRat->SetMarkerColor(plot_vec[i].h->GetLineColor());
329  SetRatioErrors(gRat, plot_vec[i].h, plot_denom.h);
330 
331  if(i==0){
332  gRat->Draw("P PLC PMC");
333  double x[2] = {gRat->GetXaxis()->GetXmin(), gRat->GetXaxis()->GetXmax()};
334  double y[2] = {1.0, 1.0};
335  TGraph* l_Xaxis = new TGraph(2, x, y);
336  l_Xaxis->SetLineColor(kGray);
337  l_Xaxis->SetLineWidth(3);
338  l_Xaxis->SetLineStyle(2);
339  l_Xaxis->Draw("L SAME");
340  }
341  else {
342  gRat->Draw("P PLC PMC");
343  }
344  } // ratio_vec loop
345 
346  float ymin, ymax;
347  GetBestYRange( ymin, ymax, ratio_vec, 0.1, 0, 10);
348  hAxes->GetYaxis()->CenterTitle();
349  hAxes->GetYaxis()->ChangeLabel();
350  hAxes->GetYaxis()->SetRangeUser(ymin, ymax);
351 
352  return;
353 } // MakeRatio
354 
355 
356 //^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^
357 void DrawUtils::SetRatioErrors(TGraphErrors* &gRat, TH1* hNum, TH1* hDenom)
358 {
359  for(int bin = 1; bin <= hNum->GetNbinsX(); bin++)
360  {
361  double num, denom, errnum, errdenom;
362  num = hNum->GetBinContent(bin); errnum = hNum->GetBinError(bin);
363  denom = hDenom->GetBinContent(bin); errdenom = hDenom->GetBinError(bin);
364 
365  if(num!=0. && denom!=0.){
366  gRat->SetPointError(bin-1, 0., (num/denom)*std::sqrt((errnum/num)*(errnum/num)+(errdenom/denom)*(errdenom/denom)));
367  }
368  else{
369  gRat->SetPointError(bin-1, 0., 0.);
370  }
371  }
372 
373  return;
374 } // SetRatioErrors
375 
376 
377 
378 
379 //^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^
380 void DrawUtils::MakeLegend( TCanvas *&cleg, const std::vector< Plot > &plots, const Plot &denom)
381 {
382  cleg->cd();
383  TLegend *leg = new TLegend(0, 0, 1, 1);
384  if(plots.size()+1>2 && plots.size()+1<7) leg->SetNColumns(3);
385  if(plots.size()+1>=7) leg->SetNColumns(5);
386 
387  leg->AddEntry(denom.h, denom.legLabel.c_str(), "L");
388  for(auto &p : plots)
389  leg->AddEntry(p.h, p.legLabel.c_str(), "L");
390 
391  leg->Draw();
392  return;
393 } // MakeLegend
394 
395 
396 
397 //^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^
398 void DrawUtils::GetBestYRange( float & ymin, float & ymax,
399  const std::vector< Plot > &plots,
400  float maxpad, int rebin, float outlierCut )
401 {
402  if(plots[0].type!="TProfile"){
403  std::cout << "GetBestYRange only works for profiles, skip"
404  << plots[0].name
405  << std::endl;
406 
407  return;
408  }
409  std::vector< TProfile* > profiles;
410  for( auto& pl : plots )
411  profiles.push_back((TProfile*)pl.h);
412 
413  GetBestYRange( ymin, ymax, profiles, maxpad, rebin, outlierCut );
414 }
415 
416 
417 
418 //^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^~^
419 void DrawUtils::GetBestYRange( float & ymin, float & ymax,
420  std::vector< TProfile* > profiles,
421  float maxpad,
422  int rebin,
423  float outlierCut )
424 {
427  TProfile * h;
428  for( auto& p: profiles ){
429  h = (TProfile *)p->Clone(Form("%s_%i",p->GetName(),Pnum()));
430  //if(rebin) h->Rebin(rebin);
431 
432  for( int b = 1; b < h->GetXaxis()->GetNbins()-1; b++ ){
433  // TProfile returns content=0 for bins with no entries.
434  // Number of entries isn't publicly accessible, but the
435  // bin error is always 0 when there are no entries.
436  if( h->GetBinError(b) == 0 ) continue;
437 
438  if( h->GetBinContent(b) < ymin ) ymin = h->GetBinContent(b);
439  if( ymax < h->GetBinContent(b) ) ymax = h->GetBinContent(b);
440  }
441  }
442 
443 
444  if(outlierCut!=0.){
445 
446  // Get RMS of all bins from all profiles combined by filling a TH1
447  TH1F * binContDist = new TH1F(Form("%s_bin_contents_%i",
448  profiles[0]->GetName(), Pnum()),"",
449  1000,ymin,ymax);
450  for( auto& p: profiles ){
451  h = (TProfile *)p->Clone(Form("%s_%i",p->GetName(),Pnum()));
452  //if(rebin) h->Rebin(rebin);
453 
454  for( int b = 1; b < h->GetXaxis()->GetNbins()-1; b++ ){
455  if( h->GetBinError(b) == 0 ) continue;
456  binContDist->Fill(h->GetBinContent(b));
457  }
458  }
459  float combinedMean = binContDist->GetMean();
460  float combinedRMS = binContDist->GetRMS();
461  //std::cout << "Combined Bins: mean = " << combinedMean
462  // << ", rms = " << combinedRMS << std::endl;
463 
464  // Last loop cutting outliers using combined RMS
467  for( auto& p: profiles ){
468  h = (TProfile *)p->Clone(Form("%s_%i",p->GetName(),Pnum()));
469  //if(rebin) h->Rebin(rebin);
470  for( int b = 1; b < h->GetXaxis()->GetNbins()-1; b++ ){
471  if( h->GetBinError(b) == 0 ) continue;
472 
473  if( outlierCut*combinedRMS < std::abs(h->GetBinContent(b) - combinedMean) ) continue;
474 
475  if( h->GetBinContent(b) < ymin ) ymin = h->GetBinContent(b);
476  if( ymax < h->GetBinContent(b) ) ymax = h->GetBinContent(b);
477  }
478  }
479 
480  }
481 
482  if( ymin == std::numeric_limits<float>::max() ||
484  std::cout << "Min/Max Y Failed in " << profiles[0]->GetName()
485  << ": min = " << ymin << ", max = " << ymax
486  << "\n The plot may be empty."
487  << std::endl;
488  return;
489  }
490 
491  float ymid = (ymax+ymin)/2;
492  float fracRange = (ymax-ymin)/ymid;
493 
494  float pad = ( (maxpad - 0.1) / (0.01 - 0.2) ) * fracRange;
495  if( pad > maxpad ) pad = maxpad; // with ceiling
496  if( pad < 0.1 ) pad = 0.1; // and floor
497 
498  ymax = ymax + pad*ymid;
499  ymin = ymin - pad*ymid;
500 
501  fLastYMax = ymax;
502  fLastYMin = ymin;
503 
504  if(ymax < ymin){
505  std::cout << "Min/Max Y Failed in " << profiles[0]->GetName()
506  << ": min = " << ymin << ", max = " << ymax
507  << std::endl;
508  return;
509  }
510 
511 } // GetBestYRange
512 
Plot(std::string file_, std::string name_, std::string legLabel_)
Definition: DrawUtils.h:39
std::string legLabel
Definition: DrawUtils.h:21
Plot(TH1 *h_, std::string name_, std::string legLabel_)
Definition: DrawUtils.h:31
void MakeRatio(TPad *&pad, const std::vector< Plot > &plot_vec, const Plot &plot_denom, std::string ratio_name)
Definition: DrawUtils.h:298
struct Plot Plot
Definition: DrawUtils.h:15
float fLastYMin
Definition: DrawUtils.h:118
float GetLastYMax()
Definition: DrawUtils.h:108
const char * p
Definition: xmltok.h:285
T sqrt(T number)
Definition: d0nt_math.hpp:156
int fPnum
Definition: DrawUtils.h:116
Plot(TH1 *h_)
Definition: DrawUtils.h:26
TCanvas * canv
~DrawUtils()
Definition: DrawUtils.h:131
int Pnum()
Definition: DrawUtils.h:113
float abs(float number)
Definition: d0nt_math.hpp:39
int fCnum
Definition: DrawUtils.h:115
TH1 * h
Definition: DrawUtils.h:19
int color
Definition: DrawUtils.h:22
Double_t ymax
Definition: plot.C:25
TH1D * MakeRatio(TH1D *num, TH1D *denom, int Col, std::string FType)
Definition: DoThePlots.C:533
std::string GetName(int i)
float fLastYMax
Definition: DrawUtils.h:119
====================================================================== ///
Definition: CutFlow_Data.C:28
void GetBestYRange(float &ymin, float &ymax, std::vector< TProfile * > profiles, float maxpad, int rebin=0, float outlierCut=0.)
Definition: DrawUtils.h:419
void DrawPlotStack(TCanvas *&canv, TCanvas *&canvleg, const std::vector< Plot > &plot_vec, const Plot &plot_denom, std::string ratio_title="")
Definition: DrawUtils.h:145
const std::vector< Plot > plots
int Cnum()
Definition: DrawUtils.h:112
void Draw2DPlot(TCanvas *&canv, TH2 *h2)
Definition: DrawUtils.h:280
std::string file
Definition: DrawUtils.h:17
float bin[41]
Definition: plottest35.C:14
float GetLastYMin()
Definition: DrawUtils.h:107
OStream cout
Definition: OStream.cxx:6
TH1F * h2
Definition: plot.C:45
static float min(const float a, const float b, const float c)
Definition: absgeo.cxx:45
int num
Definition: f2_nu.C:119
void SetRatioErrors(TGraphErrors *&gRat, TH1 *hNum, TH1 *hDenom)
Definition: DrawUtils.h:357
const hit & b
Definition: hits.cxx:21
void MakeLegend(TCanvas *&cleg, const std::vector< Plot > &plots, const Plot &denom)
Definition: DrawUtils.h:380
Double_t ymin
Definition: plot.C:24
Plot()
Definition: DrawUtils.h:23
T max(sqlite3 *const db, std::string const &table_name, std::string const &column_name)
Definition: statistics.h:68
TLegend * MakeLegend(TString beam, TString defname, TH2 *hpred, TH2 *hcosm, TPolyLine *line, TGraph *gr)
std::string type
Definition: DrawUtils.h:20
std::string name
Definition: DrawUtils.h:18
void plot(std::string label, std::map< std::string, std::map< std::string, Spectrum * >> &plots, bool log)
enum BeamMode string