Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
hist047_Graphics_candle_decay.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_hist
3/// \notebook
4/// Candle Decay, illustrate a time development of a certain value.
5///
6/// \macro_image (tcanvas_js)
7/// \macro_code
8///
9/// \author Georg Troska
10
12{
13 // create a new 2D histogram with x-axis title probability density and y-axis title time
14 // it is of a type TH2I (2D histogram with integer values)
15 auto *hist = new TH2I("hist", "Decay; probability density; time", 1000, 0, 1000, 20, 0, 20);
16 TRandom rand; // create a random number generator outside the loop to avoid reseeding
17
18 for (int iBin = 0; iBin < 19; iBin++) {
19 for (int j = 0; j < 1000000; j++) {
20 // generate a random number from a 2D Gaussian distribution
21 // a simulation of a time development of a certain value
22 float myRand = rand.Gaus(350 + iBin * 8, 20 + 2 * iBin);
23 hist->Fill(myRand, iBin);
24 }
25 }
26 hist->SetBarWidth(3);
27 hist->SetFillStyle(0);
28 hist->SetFillColor(kGray);
29 hist->SetLineColor(kBlue);
30 // create a new canvas
31 auto *can = new TCanvas("can", "Candle Decay", 800, 600);
32 can->Divide(2, 1); // divide the canvas into 2 pads (2 rows, 1 column)
33
34 can->cd(1);
35 hist->Draw("violiny(112000000)");
36 can->cd(2);
37 auto *hist2 = static_cast<TH2I *>(hist->Clone("hist2"));
38 hist2->SetBarWidth(0.8);
39 // There are six predefined candle-plot representations: (X can be replaced by Y)
40 // "CANDLEX1": Standard candle (whiskers cover the whole distribution)
41 // "CANDLEX2": Standard candle with better whisker definition + outliers. It is a good compromise
42 // "CANDLEX3": Like candle2 but with a mean as a circle. It is easier to distinguish mean and median
43 // "CANDLEX4": Like candle3 but showing the uncertainty of the median as well (notched candle plots). For bigger
44 // datasets per candle "CANDLEX5": Like candle2 but showing all data points. For very small datasets "CANDLEX6":
45 // Like candle2 but showing all datapoints scattered. For huge datasets
46 hist2->DrawCopy("candley2");
47}
@ kGray
Definition Rtypes.h:65
@ kBlue
Definition Rtypes.h:66
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
The Canvas class.
Definition TCanvas.h:23
2-D histogram with an int per channel (see TH1 documentation)
Definition TH2.h:227
This is the base class for the ROOT Random number generators.
Definition TRandom.h:27