Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
hist052_Graphics_candle_plot_whiskers.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_hist
3/// \notebook
4/// Example of candle plot showing the whiskers definition.
5///
6/// \macro_image (tcanvas_js)
7/// \macro_output
8/// \macro_code
9///
10/// \date June 2022
11/// \author Georg Troska
12
14{
15 auto c1 = new TCanvas("c1", "Candle Presets", 700, 800);
16 c1->Divide(1, 2);
17
18 auto rng = new TRandom();
19 auto h1 = new TH2I("h1", "Gaus", 100, -5, 5, 1, 0, 1);
20 auto h2 = new TH1I("h2", "Gaus", 100, -5, 5);
21
22 h1->GetXaxis()->SetTitle("Standard deviation #sigma");
23 h2->GetXaxis()->SetTitle("Standard deviation #sigma");
24 h2->GetYaxis()->SetTitle("dN/d#sigma");
25
26 float myRand;
27 for (int i = 0; i < 100000; i++) {
28 myRand = rng->Gaus(0, 1);
29 h1->Fill(myRand, 0);
30 h2->Fill(myRand);
31 }
32
33 Double_t *q = new Double_t[3];
34 Double_t *p = new Double_t[3];
35 q[0] = 0.;
36 q[1] = 0.;
37 q[2] = 0.;
38 p[0] = 0.25;
39 p[1] = 0.5;
40 p[2] = 0.75;
41
42 h2->GetQuantiles(3, q, p);
43 cout << "Q1 (-25%): " << q[0] << " Median: " << q[1] << " Q3 (+25%): " << q[2] << endl;
44 double iqr = q[2] - q[0];
45 auto mygaus_1_middle = new TF1("mygaus_1_middle", "gaus", q[0], q[2]);
46 auto mygaus_1_left = new TF1("mygaus_1_left", "gaus", q[0] - 1.5 * iqr, q[0]);
47 mygaus_1_left->SetLineColor(kGreen);
48 auto mygaus_1_right = new TF1("mygaus_1_right", "gaus", q[2], q[2] + 1.5 * iqr);
49 mygaus_1_right->SetLineColor(kGreen);
50 c1->cd(1);
51 h1->SetLineWidth(3);
52 h1->SetFillStyle(0);
53 h1->Draw("candley2 scat");
54
55 c1->cd(2);
56 h2->Draw("");
57 h2->Fit("mygaus_1_left", "R");
58 mygaus_1_left->Draw("same");
59 auto l3 = new TLine(q[0] - 1.5 * iqr, 0, q[0] - 1.5 * iqr, mygaus_1_left->Eval(q[0] - 1.5 * iqr));
60 l3->SetLineColor(kGreen);
61 l3->SetLineWidth(2);
62 l3->Draw("");
63 auto l1 = new TLine(q[0], 0, q[0], mygaus_1_left->Eval(q[0]));
64 l1->SetLineWidth(2);
65 l1->SetLineColor(kGreen);
66 l1->Draw("");
67
68 h2->Fit("mygaus_1_right", "R", "");
69 mygaus_1_right->Draw("same");
70 auto l4 = new TLine(q[2] + 1.5 * iqr, 0, q[2] + 1.5 * iqr, mygaus_1_left->Eval(q[2] + 1.5 * iqr));
71 l4->SetLineColor(kGreen);
72 l4->SetLineWidth(2);
73 l4->Draw("");
74 auto l5 = new TLine(q[2], 0, q[2], mygaus_1_right->Eval(q[2]));
75 l5->SetLineWidth(2);
76 l5->SetLineColor(kGreen);
77 l5->Draw("");
78
79 h2->Fit("mygaus_1_middle", "R");
80 mygaus_1_middle->Draw("same");
81
82 // In principal one could calculate these values by h2->Integral() as well
83 TText t;
84 t.SetTextFont(42);
85 t.DrawText(0, mygaus_1_middle->Eval(0) / 2, "50%");
86 t.DrawText(-1.5, mygaus_1_middle->Eval(-1.5) / 2, "24.65%");
87 t.DrawText(+1, mygaus_1_middle->Eval(+1.5) / 2, "24.65%");
88 t.DrawText(q[0] - 1.5 * iqr, 1000, Form("%.3f", q[0] - 1.5 * iqr))->SetTextAngle(90);
89 t.DrawText(q[2] + 1.5 * iqr, 1000, Form("%.3f", q[2] + 1.5 * iqr))->SetTextAngle(90);
90 t.DrawText(q[0], 1000, Form("%.3f", q[0]))->SetTextAngle(90);
91 t.DrawText(q[2], 1000, Form("%.3f", q[2]))->SetTextAngle(90);
92}
double Double_t
Definition RtypesCore.h:59
@ kGreen
Definition Rtypes.h:66
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
winID h TVirtualViewer3D TVirtualGLPainter p
float * q
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2489
virtual void SetFillStyle(Style_t fstyle)
Set the fill area style.
Definition TAttFill.h:40
virtual void SetLineWidth(Width_t lwidth)
Set the line width.
Definition TAttLine.h:45
virtual void SetTextFont(Font_t tfont=62)
Set the text font.
Definition TAttText.h:48
The Canvas class.
Definition TCanvas.h:23
1-Dim function class
Definition TF1.h:233
1-D histogram with an int per channel (see TH1 documentation)
Definition TH1.h:563
TAxis * GetXaxis()
Definition TH1.h:340
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition TH1.cxx:3315
void Draw(Option_t *option="") override
Draw this histogram with options.
Definition TH1.cxx:3037
2-D histogram with an int per channel (see TH1 documentation)
Definition TH2.h:227
Use the TLine constructor to create a simple line.
Definition TLine.h:22
virtual void SetTitle(const char *title="")
Set the title of the TNamed.
Definition TNamed.cxx:164
This is the base class for the ROOT Random number generators.
Definition TRandom.h:27
Base class for several text objects.
Definition TText.h:22
virtual TText * DrawText(Double_t x, Double_t y, const char *text)
Draw this text with new coordinates.
Definition TText.cxx:176
return c1
Definition legend1.C:41
TH1F * h1
Definition legend1.C:5