CbmRoot
CbmDrawHist.cxx
Go to the documentation of this file.
1 
6 #include "CbmDrawHist.h"
7 
8 #include "CbmUtils.h" // for NumberToString
9 
10 #include <TAxis.h> // for TAxis
11 #include <TF1.h> // for TF1
12 #include <TGraph.h> // for TGraph
13 #include <TGraph2D.h> // for TGraph2D
14 #include <TH1.h> // for TH1, TH1D
15 #include <TH2.h> // for TH2, TH2D
16 #include <TH3.h> // for TH3
17 #include <TLatex.h> // for TLatex
18 #include <TLegend.h> // for TLegend
19 #include <TStyle.h> // for TStyle, gStyle
20 #include <TVirtualPad.h> // for TVirtualPad, gPad
21 
22 #include <algorithm> // for max, min
23 #include <cassert> // for assert
24 #include <iostream> // for string, operator<<, basic_ostream, strin...
25 #include <limits> // for numeric_limits
26 #include <string> // for operator+, allocator, basic_string, char...
27 
28 using std::string;
29 using std::stringstream;
30 using std::vector;
31 
32 /* Set default styles for histograms. */
34  gStyle->SetOptStat("rm");
35  gStyle->SetOptFit(1);
36  gStyle->SetOptTitle(0);
37 
38  gStyle->SetCanvasColor(kWhite);
39  gStyle->SetFrameFillColor(kWhite);
40  gStyle->SetFrameBorderMode(0);
41  gStyle->SetPadColor(kWhite);
42  gStyle->SetStatColor(kWhite);
43  gStyle->SetTitleFillColor(kWhite);
44  gStyle->SetPalette(55, 0);
45  // gStyle->SetPalette(1);
46 }
47 
48 /* Draw TH1 histogram.*/
49 void DrawH1(TH1* hist,
50  HistScale logx,
51  HistScale logy,
52  const string& drawOpt,
53  Int_t color,
54  Int_t lineWidth,
55  Int_t lineStyle,
56  Int_t markerSize,
57  Int_t markerStyle) {
58  Double_t textSize = CbmDrawingOptions::TextSize();
59  hist->SetLineColor(color);
60  hist->SetLineWidth(lineWidth);
61  hist->SetLineStyle(lineStyle);
62  hist->SetMarkerColor(color);
63  hist->SetMarkerSize(markerSize);
64  hist->SetMarkerStyle(markerStyle);
65  if (logx == kLog) { gPad->SetLogx(); }
66  if (logy == kLog) { gPad->SetLogy(); }
67  hist->GetXaxis()->SetLabelSize(textSize);
68  hist->GetXaxis()->SetNdivisions(505, kTRUE);
69  hist->GetYaxis()->SetLabelSize(textSize);
70  hist->GetXaxis()->SetTitleSize(textSize);
71  hist->GetYaxis()->SetTitleSize(textSize);
72  hist->GetXaxis()->SetTitleOffset(1.0);
73  hist->GetYaxis()->SetTitleOffset(1.3);
74  gPad->SetLeftMargin(0.17);
75  gPad->SetBottomMargin(0.15);
76  gPad->SetTopMargin(0.12);
77  gPad->SetTicks(1, 1);
78  hist->Draw(drawOpt.c_str());
79  gPad->SetGrid(true, true);
80  hist->SetStats(false);
81 }
82 
83 /* Draw TH2 histogram.*/
84 void DrawH2(TH2* hist,
85  HistScale logx,
86  HistScale logy,
87  HistScale logz,
88  const string& drawOpt) {
89  Double_t textSize = CbmDrawingOptions::TextSize();
90  if (logx == kLog) { gPad->SetLogx(); }
91  if (logy == kLog) { gPad->SetLogy(); }
92  if (logz == kLog) { gPad->SetLogz(); }
93  hist->GetXaxis()->SetLabelSize(textSize);
94  hist->GetXaxis()->SetNdivisions(505, kTRUE);
95  hist->GetYaxis()->SetLabelSize(textSize);
96  hist->GetYaxis()->SetNdivisions(505, kTRUE);
97  hist->GetZaxis()->SetLabelSize(textSize);
98  // hist->GetZaxis()->SetNdivisions(505, kTRUE);
99  hist->GetXaxis()->SetTitleSize(textSize);
100  hist->GetYaxis()->SetTitleSize(textSize);
101  hist->GetZaxis()->SetTitleSize(textSize);
102  hist->GetXaxis()->SetTitleOffset(1.0);
103  hist->GetYaxis()->SetTitleOffset(1.3);
104  hist->GetZaxis()->SetTitleOffset(1.5);
105  gPad->SetLeftMargin(0.17);
106  gPad->SetRightMargin(0.30);
107  gPad->SetBottomMargin(0.15);
108  gPad->SetTicks(1, 1);
109  hist->Draw(drawOpt.c_str());
110  gPad->SetGrid(true, true);
111  hist->SetStats(false);
112 }
113 
114 /* Draw several TH1 histograms. */
115 void DrawH1(const vector<TH1*>& histos,
116  const vector<string>& histLabels,
117  HistScale logx,
118  HistScale logy,
119  Bool_t drawLegend,
120  Double_t x1,
121  Double_t y1,
122  Double_t x2,
123  Double_t y2,
124  const string& drawOpt) {
125  assert(histos.size() != 0 && histLabels.size() == histos.size());
127  UInt_t nofHistos = histos.size();
128  TLegend* legend = new TLegend(x1, y1, x2, y2);
129  legend->SetFillColor(kWhite);
130  for (UInt_t iHist = 0; iHist < nofHistos; iHist++) {
131  TH1* hist = histos[iHist];
132  string opt = (iHist == 0) ? drawOpt
133  : (iHist == nofHistos - 1) ? "SAME" + drawOpt
134  : "SAME" + drawOpt;
135  DrawH1(hist,
136  logx,
137  logy,
138  opt,
144  max = std::max(max, hist->GetMaximum());
145  legend->AddEntry(hist, histLabels[iHist].c_str(), "lp");
146  }
147  histos[0]->SetMaximum(max * 1.17);
148  if (drawLegend) { legend->Draw(); }
149 }
150 
151 void DrawGraph(TGraph* graph,
152  HistScale logx,
153  HistScale logy,
154  const string& drawOpt,
155  Int_t color,
156  Int_t lineWidth,
157  Int_t lineStyle,
158  Int_t markerSize,
159  Int_t markerStyle) {
160  Double_t textSize = CbmDrawingOptions::TextSize();
161  graph->SetLineColor(color);
162  graph->SetLineWidth(lineWidth);
163  graph->SetLineStyle(lineStyle);
164  graph->SetMarkerColor(color);
165  graph->SetMarkerSize(markerSize);
166  graph->SetMarkerStyle(markerStyle);
167  if (drawOpt.find("A") != string::npos) {
168  if (logx == kLog) { gPad->SetLogx(); }
169  if (logy == kLog) { gPad->SetLogy(); }
170  graph->GetXaxis()->SetLabelSize(textSize);
171  graph->GetXaxis()->SetNdivisions(505, kTRUE);
172  graph->GetYaxis()->SetLabelSize(textSize);
173  graph->GetXaxis()->SetTitleSize(textSize);
174  graph->GetYaxis()->SetTitleSize(textSize);
175  graph->GetXaxis()->SetTitleOffset(1.0);
176  graph->GetYaxis()->SetTitleOffset(1.3);
177  }
178  gPad->SetLeftMargin(0.17);
179  gPad->SetBottomMargin(0.15);
180  graph->Draw(drawOpt.c_str());
181  gPad->SetGrid(true, true);
182 }
183 
184 
185 /* Draw several TGraphs. */
186 void DrawGraph(const vector<TGraph*>& graphs,
187  const vector<string>& graphLabels,
188  HistScale logx,
189  HistScale logy,
190  Bool_t drawLegend,
191  Double_t x1,
192  Double_t y1,
193  Double_t x2,
194  Double_t y2) {
195  assert(graphs.size() != 0 && graphs.size() == graphLabels.size());
196 
199  TLegend* legend = new TLegend(x1, y1, x2, y2);
200  legend->SetFillColor(kWhite);
201  UInt_t nofGraphs = graphs.size();
202  for (UInt_t iGraph = 0; iGraph < nofGraphs; iGraph++) {
203  TGraph* graph = graphs[iGraph];
204  string opt = (iGraph == 0) ? "ACP" : "CP";
205  DrawGraph(graph,
206  logx,
207  logy,
208  opt,
209  CbmDrawingOptions::Color(iGraph),
214  max = std::max(graph->GetYaxis()->GetXmax(), max);
215  min = std::min(graph->GetYaxis()->GetXmin(), min);
216  legend->AddEntry(graph, graphLabels[iGraph].c_str(), "lp");
217  }
218  graphs[0]->SetMaximum(max);
219  graphs[0]->SetMinimum(min);
220  if (drawLegend) { legend->Draw(); }
221 }
222 
223 /* Draws 2D graph.*/
224 void DrawGraph2D(TGraph2D* graph,
225  HistScale logx,
226  HistScale logy,
227  HistScale logz,
228  const string& drawOpt) {
229  Double_t textSize = CbmDrawingOptions::TextSize();
230  if (logx == kLog) { gPad->SetLogx(); }
231  if (logy == kLog) { gPad->SetLogy(); }
232  if (logz == kLog) { gPad->SetLogz(); }
233  graph->GetXaxis()->SetLabelSize(textSize);
234  graph->GetXaxis()->SetNdivisions(505, kTRUE);
235  graph->GetYaxis()->SetLabelSize(textSize);
236  graph->GetYaxis()->SetNdivisions(505, kTRUE);
237  graph->GetZaxis()->SetLabelSize(textSize);
238  // graph->GetZaxis()->SetNdivisions(505, kTRUE);
239  graph->GetXaxis()->SetTitleSize(textSize);
240  graph->GetYaxis()->SetTitleSize(textSize);
241  graph->GetZaxis()->SetTitleSize(textSize);
242  graph->GetXaxis()->SetTitleOffset(1.0);
243  graph->GetYaxis()->SetTitleOffset(1.3);
244  graph->GetZaxis()->SetTitleOffset(1.7);
245  gPad->SetLeftMargin(0.17);
246  gPad->SetRightMargin(0.2);
247  gPad->SetBottomMargin(0.15);
248  gPad->SetTicks(1, 1);
249  graph->Draw(drawOpt.c_str());
250  gPad->SetGrid(true, true);
251 }
252 
253 void DrawTextOnPad(const string& text,
254  Double_t x1,
255  Double_t y1,
256  Double_t x2,
257  Double_t y2) {
258  TLegend* leg = new TLegend(x1, y1, x2, y2);
259  leg->AddEntry(new TH2D(), text.c_str(), "");
260  leg->SetFillColor(kWhite);
261  leg->SetFillStyle(0);
262  leg->SetBorderSize(0);
263  leg->Draw();
264 }
265 
266 void DrawH1andFitGauss(TH1* hist,
267  Bool_t drawResults,
268  Bool_t doScale,
269  Double_t userRangeMin,
270  Double_t userRangeMax) {
271  if (hist == nullptr) return;
272 
273  hist->GetYaxis()->SetTitle("Yield (a.u.)");
274  DrawH1(hist, kLinear, kLinear, "hist");
275  if (doScale) hist->Scale(1. / hist->Integral());
276  if (!(userRangeMin == 0. && userRangeMax == 0.))
277  hist->GetXaxis()->SetRangeUser(userRangeMin, userRangeMax);
278  hist->Fit("gaus", "Q");
279  TF1* func = hist->GetFunction("gaus");
280  if (func == nullptr) return;
281  func->SetLineColor(kBlack);
282  func->Draw("same");
283 
284  if (drawResults) {
285  double m = func->GetParameter(1);
286  double s = func->GetParameter(2);
287  string txt1 = Cbm::NumberToString<Double_t>(m, 2) + " / "
288  + Cbm::NumberToString<Double_t>(s, 2);
289  TLatex text;
290  text.SetTextAlign(21);
291  text.SetTextSize(0.06);
292  text.DrawTextNDC(0.5, 0.92, txt1.c_str());
293  }
294 }
295 
296 void DrawH2WithProfile(TH2* hist,
297  Bool_t doGaussFit,
298  Bool_t drawOnlyMean,
299  const string& drawOpt2D,
300  Int_t profileColor,
301  Int_t profileLineWidth) {
302  if (hist == nullptr) return;
303 
304  // TProfile does not allow to fit individual slice with gauss
305  /* TProfile* profX = (TProfile*)hist->ProfileX("_pfx", 1, -1, "s")->Clone();
306  DrawH2(hist, kLinear, kLinear, kLinear, drawOpt);
307  if (!drawOnlyMean) {
308  DrawH1(profX, kLinear, kLinear, "same", profileColor, profileLineWidth, 1, 1, kOpenCircle);
309  } else {
310  DrawH1(profX, kLinear, kLinear, "same hist p", profileColor, profileLineWidth, 1, 1, kFullCircle);
311  }*/
312 
313  TH1D* hMean =
314  (TH1D*) hist->ProjectionX((string(hist->GetName()) + "_mean").c_str())
315  ->Clone();
316  string yTitle = (doGaussFit) ? "Mean and sigma. " : "Mean and RMS. ";
317  hMean->GetYaxis()->SetTitle(
318  (yTitle + string(hist->GetYaxis()->GetTitle())).c_str());
319  for (Int_t i = 1; i <= hist->GetXaxis()->GetNbins(); i++) {
320  stringstream ss;
321  ss << string(hist->GetName()) << "_py" << i;
322  TH1D* pr = hist->ProjectionY(ss.str().c_str(), i, i);
323  if (hMean == nullptr || pr == nullptr) continue;
324 
325  Double_t m = 0., s = 0.;
326  if (doGaussFit) {
327  pr->Fit("gaus", "QO");
328  TF1* func = pr->GetFunction("gaus");
329  if (func != nullptr) {
330  m = func->GetParameter(1);
331  s = func->GetParameter(2);
332  }
333  } else {
334  m = pr->GetMean();
335  s = pr->GetRMS();
336  }
337 
338  hMean->SetBinContent(i, m);
339  if (!drawOnlyMean) {
340  hMean->SetBinError(i, s);
341  } else {
342  hMean->SetBinError(i, 0.);
343  }
344  }
345  DrawH2(hist, kLinear, kLinear, kLinear, drawOpt2D);
346  if (!drawOnlyMean) {
347  DrawH1(hMean,
348  kLinear,
349  kLinear,
350  "same",
351  profileColor,
352  profileLineWidth,
353  1,
354  1,
355  kOpenCircle);
356  } else {
357  DrawH1(hMean,
358  kLinear,
359  kLinear,
360  "same hist p",
361  profileColor,
362  profileLineWidth,
363  1,
364  1,
365  kFullCircle);
366  }
367 }
368 
369 TH2D* DrawH3Profile(TH3* h,
370  Bool_t drawMean,
371  Bool_t doGaussFit,
372  Double_t zMin,
373  Double_t zMax) {
374  Int_t nBinsX = h->GetNbinsX();
375  Int_t nBinsY = h->GetNbinsY();
376  TH2D* h2 = (TH2D*) h->Project3D("yx")->Clone();
377 
378  for (Int_t x = 1; x <= nBinsX; x++) {
379  for (Int_t y = 1; y <= nBinsY; y++) {
380  stringstream ss;
381  ss << h->GetName() << "_z_" << x << "_" << y;
382  TH1D* hz = h->ProjectionZ(ss.str().c_str(), x, x, y, y);
383  Double_t ms = 0.;
384  if (doGaussFit) {
385  hz->Fit("gaus", "QO");
386  TF1* func = hz->GetFunction("gaus");
387  if (func != nullptr) {
388  ms = (drawMean) ? func->GetParameter(1) : func->GetParameter(2);
389  }
390  } else {
391  ms = (drawMean) ? hz->GetMean() : hz->GetRMS();
392  }
393  h2->SetBinContent(x, y, ms);
394  }
395  }
396 
397  string zAxisTitle = string(h->GetZaxis()->GetTitle());
398  string sigmaRms = (doGaussFit) ? "Sigma." : "RMS.";
399  zAxisTitle = (drawMean) ? "Mean." + zAxisTitle : sigmaRms + zAxisTitle;
400 
401  h2->GetZaxis()->SetTitle(zAxisTitle.c_str());
402  if (zMin < zMax) h2->GetZaxis()->SetRangeUser(zMin, zMax);
403  DrawH2(h2);
404 
405  return h2;
406 }
DrawH3Profile
TH2D * DrawH3Profile(TH3 *h, Bool_t drawMean, Bool_t doGaussFit, Double_t zMin, Double_t zMax)
Definition: CbmDrawHist.cxx:369
HistScale
HistScale
Define linear or logarithmic scale for drawing.
Definition: CbmDrawHist.h:77
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
DrawGraph2D
void DrawGraph2D(TGraph2D *graph, HistScale logx, HistScale logy, HistScale logz, const string &drawOpt)
Definition: CbmDrawHist.cxx:224
CbmDrawHist.h
Helper functions for drawing 1D and 2D histograms and graphs.
DrawGraph
void DrawGraph(TGraph *graph, HistScale logx, HistScale logy, const string &drawOpt, Int_t color, Int_t lineWidth, Int_t lineStyle, Int_t markerSize, Int_t markerStyle)
Definition: CbmDrawHist.cxx:151
min
friend F32vec4 min(const F32vec4 &a, const F32vec4 &b)
Definition: L1/vectors/P4_F32vec4.h:33
DrawTextOnPad
void DrawTextOnPad(const string &text, Double_t x1, Double_t y1, Double_t x2, Double_t y2)
Definition: CbmDrawHist.cxx:253
DrawH1
void DrawH1(TH1 *hist, HistScale logx, HistScale logy, const string &drawOpt, Int_t color, Int_t lineWidth, Int_t lineStyle, Int_t markerSize, Int_t markerStyle)
Definition: CbmDrawHist.cxx:49
CbmDrawingOptions::MarkerSize
static Int_t MarkerSize()
Definition: CbmDrawHist.h:56
CbmDrawingOptions::MarkerStyle
static Int_t MarkerStyle(Int_t markerIndex)
Definition: CbmDrawHist.h:58
h
Data class with information on a STS local track.
kLinear
@ kLinear
Definition: CbmDrawHist.h:79
DrawH2WithProfile
void DrawH2WithProfile(TH2 *hist, Bool_t doGaussFit, Bool_t drawOnlyMean, const string &drawOpt2D, Int_t profileColor, Int_t profileLineWidth)
Definition: CbmDrawHist.cxx:296
CbmUtils.h
CbmDrawingOptions::LineWidth
static Int_t LineWidth()
Definition: CbmDrawHist.h:52
CbmDrawingOptions::LineStyle
static Int_t LineStyle(Int_t lineStyleIndex)
Definition: CbmDrawHist.h:54
x
Double_t x
Definition: CbmMvdSensorDigiToHitTask.cxx:68
CbmDrawingOptions::TextSize
static Double_t TextSize()
Definition: CbmDrawHist.h:70
m
__m128 m
Definition: L1/vectors/P4_F32vec4.h:26
y
Double_t y
Definition: CbmMvdSensorDigiToHitTask.cxx:68
DrawH2
void DrawH2(TH2 *hist, HistScale logx, HistScale logy, HistScale logz, const string &drawOpt)
Definition: CbmDrawHist.cxx:84
SetDefaultDrawStyle
void SetDefaultDrawStyle()
Definition: CbmDrawHist.cxx:33
max
friend F32vec4 max(const F32vec4 &a, const F32vec4 &b)
Definition: L1/vectors/P4_F32vec4.h:36
kLog
@ kLog
Definition: CbmDrawHist.h:78
DrawH1andFitGauss
void DrawH1andFitGauss(TH1 *hist, Bool_t drawResults, Bool_t doScale, Double_t userRangeMin, Double_t userRangeMax)
Definition: CbmDrawHist.cxx:266
CbmDrawingOptions::Color
static Int_t Color(Int_t colorIndex)
Definition: CbmDrawHist.h:32