Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
gr012_polar.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_graphs
3## \notebook
4## \preview Create and draw a polar graph. See the [TGraphPolar documentation](https://root.cern/doc/master/classTGraphPolar.html)
5##
6## Since TGraphPolar is a TGraphErrors, it is painted with
7## [TGraphPainter](https://root.cern/doc/master/classTGraphPainter.html) options.
8##
9## With GetPolargram we retrieve the polar axis to format it; see the
10## [TGraphPolargram documentation](https://root.cern/doc/master/classTGraphPolargram.html)
11##
12## \macro_image
13## \macro_code
14## \author Olivier Couet, Jamie Gooding
15
16# Illustrates how to use TGraphPolar
17
18import math
19
20import numpy as np
21import ROOT
22
23CPol = ROOT.TCanvas("CPol", "TGraphPolar Examples", 1200, 600)
24CPol.Divide(2, 1)
25
26# Left-side pad. Two graphs without errors
27CPol.cd(1)
28xmin = 0
29xmax = math.pi * 2
30
31x = np.array([])
32y = np.array([])
33xval1 = np.array([])
34yval1 = np.array([])
35
36# Graph 1 to be drawn with line and fill
37fplot = ROOT.TF1("fplot", "cos(2*x)*cos(20*x)", xmin, xmax)
38for ipt in range(1000):
39 x = np.append(x, ipt * (xmax - xmin) / 1000 + xmin)
40 y = np.append(y, fplot.Eval(x[ipt]))
41
42grP = ROOT.TGraphPolar(1000, x, y)
47grP.Draw("AFL")
48
49# Graph 2 to be drawn superposed over graph 1, with curve and polymarker
50for ipt in range(20):
51 xval1 = np.append(xval1, x[math.floor(1000 / 20 * ipt)])
52 yval1 = np.append(yval1, y[math.floor(1000 / 20 * ipt)])
53
54grP1 = ROOT.TGraphPolar(20, xval1, yval1)
59grP1.Draw("CP")
60
61# To format the polar axis, we retrieve the TGraphPolargram.
62# First update the canvas, otherwise GetPolargram returns 0
66 grP1.GetPolargram().SetRangePolar(-math.pi, math.pi)
67 grP1.GetPolargram().SetNdivPolar(703)
68 grP1.GetPolargram().SetToRadian() # tell ROOT that the x and xval1 are in radians
69
70
71# Right-side pad. One graph with errors
72CPol.cd(2)
73x2 = np.array([])
74y2 = np.array([])
75ex = np.array([])
76ey = np.array([])
77for ipt in range(30):
78 x2 = np.append(x2, x[math.floor(1000 / 30 * ipt)])
79 y2 = np.append(y2, 1.2 + 0.4 * math.sin(math.pi * 2 * ipt / 30))
80 ex = np.append(ex, 0.2 + 0.1 * math.cos(2 * math.pi / 30 * ipt))
81 ey = np.append(ey, 0.2)
82
83# Grah to be drawn with polymarker and errors
84grPE = ROOT.TGraphPolar(30, x2, y2, ex, ey)
90grPE.Draw("EP")
91
92# To format the polar axis, we retrieve the TGraphPolargram.
93# First update the canvas, otherwise GetPolargram returns 0
97 grPE.GetPolargram().SetTwoPi()
98 grPE.GetPolargram().SetToRadian() # tell ROOT that the x2 values are in radians
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t SetTextSize
pt SetTextColor(4)