Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
hist102_TH2_contour_list.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_hist
3/// \notebook
4/// Getting Contours From TH2D.
5///
6/// #### Image produced by `.x ContourList.C`
7/// The contours values are drawn next to each contour.
8/// \macro_image
9///
10/// #### Output produced by `.x ContourList.C`
11/// It shows that 6 contours and 12 graphs were found.
12/// \macro_output
13///
14/// #### `ContourList.C`
15/// \macro_code
16///
17/// \date November 2022
18/// \authors Josh de Bever (CSI Medical Physics Group, The University of Western Ontario, London, Ontario, Canada),
19/// Olivier Couet
20
22
24{
25
26 const Double_t PI = TMath::Pi();
27
28 TCanvas *c = new TCanvas("c", "Contour List", 0, 0, 600, 600);
29 c->SetRightMargin(0.15);
30 c->SetTopMargin(0.15);
31
32 Int_t i, j;
33
34 Int_t nZsamples = 80;
35 Int_t nPhiSamples = 80;
36
37 Double_t HofZwavelength = 4.0; // 4 meters
39 Double_t dPhi = 2 * PI / (Double_t)(nPhiSamples - 1);
40
45
46 // Discretized Z and Phi Values
47 for (i = 0; i < nZsamples; i++) {
48 z[i] = (i)*dZ - HofZwavelength / 2.0;
49 HofZ[i] = SawTooth(z[i], HofZwavelength);
50 }
51
52 for (Int_t i = 0; i < nPhiSamples; i++) {
53 phi[i] = (i)*dPhi;
54 FofPhi[i] = sin(phi[i]);
55 }
56
57 // Create Histogram
59 new TH2D("HstreamFn",
60 "#splitline{Histogram with negative and positive contents. Six contours are defined.}{It is plotted "
61 "with options CONT LIST to retrieve the contours points in TGraphs}",
62 nZsamples, z[0], z[nZsamples - 1], nPhiSamples, phi[0], phi[nPhiSamples - 1]);
63
64 // Load Histogram Data
65 for (Int_t i = 0; i < nZsamples; i++) {
66 for (Int_t j = 0; j < nPhiSamples; j++) {
67 HistStreamFn->SetBinContent(i, j, HofZ[i] * FofPhi[j]);
68 }
69 }
70
72 gStyle->SetTitleW(0.99);
73 gStyle->SetTitleH(0.08);
74
76 contours[0] = -0.7;
77 contours[1] = -0.5;
78 contours[2] = -0.1;
79 contours[3] = 0.1;
80 contours[4] = 0.4;
81 contours[5] = 0.8;
82
83 HistStreamFn->SetContour(6, contours);
84
85 // Draw contours as filled regions, and Save points
86 HistStreamFn->Draw("CONT Z LIST");
87 c->Update(); // Needed to force the plotting and retrieve the contours in TGraphs
88
89 // Get Contours
90 TObjArray *conts = (TObjArray *)gROOT->GetListOfSpecials()->FindObject("contours");
91
92 if (!conts) {
93 printf("*** No Contours Were Extracted!\n");
94 return nullptr;
95 }
96
97 TList *contLevel = nullptr;
98 TGraph *curv = nullptr;
99 TGraph *gc = nullptr;
100
101 Int_t nGraphs = 0;
102 Int_t TotalConts = conts->GetSize();
103
104 printf("TotalConts = %d\n", TotalConts);
105
106 for (i = 0; i < TotalConts; i++) {
107 contLevel = (TList *)conts->At(i);
108 printf("Contour %d has %d Graphs\n", i, contLevel->GetSize());
109 nGraphs += contLevel->GetSize();
110 }
111
112 nGraphs = 0;
113
114 TCanvas *c1 = new TCanvas("c1", "Contour List", 610, 0, 600, 600);
115 c1->SetTopMargin(0.15);
116 TH2F *hr = new TH2F(
117 "hr",
118 "#splitline{Negative contours are returned first (highest to lowest). Positive contours are returned "
119 "from}{lowest to highest. On this plot Negative contours are drawn in red and positive contours in blue.}",
120 2, -2, 2, 2, 0, 6.5);
121
122 hr->Draw();
124 TLatex l;
125 l.SetTextSize(0.03);
126 char val[20];
127
128 for (i = 0; i < TotalConts; i++) {
129 contLevel = (TList *)conts->At(i);
130 if (i < 3)
131 zval0 = contours[2 - i];
132 else
133 zval0 = contours[i];
134 printf("Z-Level Passed in as: Z = %f\n", zval0);
135
136 // Get first graph from list on curves on this level
137 curv = (TGraph *)contLevel->First();
138 for (j = 0; j < contLevel->GetSize(); j++) {
139 curv->GetPoint(0, xval0, yval0);
140 if (zval0 < 0)
141 curv->SetLineColor(kRed);
142 if (zval0 > 0)
143 curv->SetLineColor(kBlue);
144 nGraphs++;
145 printf("\tGraph: %d -- %d Elements\n", nGraphs, curv->GetN());
146
147 // Draw clones of the graphs to avoid deletions in case the 1st
148 // pad is redrawn.
149 gc = (TGraph *)curv->Clone();
150 gc->Draw("C");
151
152 l.DrawLatex(xval0, yval0, Form("%g", zval0));
153 curv = (TGraph *)contLevel->After(curv); // Get Next graph
154 }
155 }
156 c1->Update();
157 printf("\n\n\tExtracted %d Contours and %d Graphs \n", TotalConts, nGraphs);
158 gStyle->SetTitleW(0.);
159 gStyle->SetTitleH(0.);
160 return c1;
161}
162
164{
165
166 // This function is specific to a sawtooth function with period
167 // WaveLen, symmetric about x = 0, and with amplitude = 1. Each segment
168 // is 1/4 of the wavelength.
169 //
170 // |
171 // /\ |
172 // / \ |
173 // / \ |
174 // / \|
175 // /--------\--------/------------
176 // |\ /
177 // | \ /
178 // | \ /
179 // | \/
180 //
181
182 Double_t y;
183 if ((x < -WaveLen / 2) || (x > WaveLen / 2))
184 y = -99999999; // Error X out of bounds
185 if (x <= -WaveLen / 4) {
186 y = x + 2.0;
187 } else if ((x > -WaveLen / 4) && (x <= WaveLen / 4)) {
188 y = -x;
189 } else if ((x > WaveLen / 4) && (x <= WaveLen / 2)) {
190 y = x - 2.0;
191 }
192 return y;
193}
#define c(i)
Definition RSha256.hxx:101
int Int_t
Definition RtypesCore.h:45
double Double_t
Definition RtypesCore.h:59
@ kRed
Definition Rtypes.h:66
@ kBlue
Definition Rtypes.h:66
#define PI
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void gc
#define gROOT
Definition TROOT.h:406
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2489
R__EXTERN TStyle * gStyle
Definition TStyle.h:442
Array of doubles (64 bits per element).
Definition TArrayD.h:27
The Canvas class.
Definition TCanvas.h:23
A TGraph is an object made of two arrays X and Y with npoints each.
Definition TGraph.h:41
2-D histogram with a double per channel (see TH1 documentation)
Definition TH2.h:356
2-D histogram with a float per channel (see TH1 documentation)
Definition TH2.h:307
To draw Mathematical Formula.
Definition TLatex.h:18
A doubly linked list.
Definition TList.h:38
An array of TObjects.
Definition TObjArray.h:31
void SetOptStat(Int_t stat=1)
The type of information printed in the histogram statistics box can be selected via the parameter mod...
Definition TStyle.cxx:1642
void SetTitleW(Float_t w=0)
Definition TStyle.h:419
void SetTitleH(Float_t h=0)
Definition TStyle.h:420
Double_t y[n]
Definition legend1.C:17
return c1
Definition legend1.C:41
Double_t x[n]
Definition legend1.C:17
constexpr Double_t Pi()
Definition TMath.h:37
TLine l
Definition textangle.C:4