Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
Bessel.py
Go to the documentation of this file.
1## \file
2## \ingroup tutorial_math
3## \notebook
4## Show the different kinds of Bessel functions available in ROOT
5## To execute the macro type in:
6##
7## ~~~{.cpp}
8## root[0] .x Bessel.C
9## ~~~
10##
11## It will create one canvas with the representation
12## of the cylindrical and spherical Bessel functions
13## regular and modified
14##
15## Based on Bessel.C by Magdalena Slawinska
16##
17## \macro_image
18## \macro_code
19##
20## \author Juan Fernando Jaramillo Botero
21
22from ROOT import TCanvas, TF1, gSystem, gPad, TLegend, TPaveLabel, kBlack
23
24
25gSystem.Load("libMathMore")
26
27DistCanvas = TCanvas("DistCanvas", "Bessel functions example", 10, 10, 800, 600)
28DistCanvas.SetFillColor(17)
29DistCanvas.Divide(2, 2)
30DistCanvas.cd(1)
31gPad.SetGrid()
32gPad.SetFrameFillColor(19)
33leg = TLegend(0.75, 0.7, 0.89, 0.89)
34
35# Drawing the set of Bessel J functions
36#
37# n is the number of functions in each pad
38n = 5
39JBessel = []
40for nu in range(n):
41 jbessel = TF1("J_0", "ROOT::Math::cyl_bessel_j([0],x)", 0, 10)
42 jbessel.SetParameters(nu, 0.0)
43 jbessel.SetTitle("")
44 jbessel.SetLineStyle(1)
45 jbessel.SetLineWidth(3)
46 jbessel.SetLineColor(nu + 1)
47 JBessel.append(jbessel)
48
49# Setting x axis for JBessel
50xaxis = JBessel[0].GetXaxis()
51xaxis.SetTitle("x")
52xaxis.SetTitleSize(0.06)
53xaxis.SetTitleOffset(.7)
54
55# setting the title in a label style
56p1 = TPaveLabel(.0, .90, .0 + .50, .90 + .10, "Bessel J functions", "NDC")
57p1.SetFillColor(0)
58p1.SetTextFont(22)
59p1.SetTextColor(kBlack)
60
61# setting the legend
62leg.AddEntry(JBessel[0].DrawCopy(), " J_0(x)", "l")
63leg.AddEntry(JBessel[1].DrawCopy("same"), " J_1(x)", "l")
64leg.AddEntry(JBessel[2].DrawCopy("same"), " J_2(x)", "l")
65leg.AddEntry(JBessel[3].DrawCopy("same"), " J_3(x)", "l")
66leg.AddEntry(JBessel[4].DrawCopy("same"), " J_4(x)", "l")
67
68leg.Draw()
69p1.Draw()
70
71# Set canvas 2
72DistCanvas.cd(2)
73gPad.SetGrid()
74gPad.SetFrameFillColor(19)
75leg2 = TLegend(0.75, 0.7, 0.89, 0.89)
76
77# Drawing Bessel k
78KBessel = []
79for nu in range(n):
80 kbessel = TF1("J_0", "ROOT::Math::cyl_bessel_k([0],x)", 0, 10)
81 kbessel.SetParameters(nu, 0.0)
82 kbessel.SetTitle("Bessel K functions")
83 kbessel.SetLineStyle(1)
84 kbessel.SetLineWidth(3)
85 kbessel.SetLineColor(nu+1)
86 KBessel.append(kbessel)
87kxaxis = KBessel[0].GetXaxis()
88kxaxis.SetTitle("x")
89kxaxis.SetTitleSize(0.06)
90kxaxis.SetTitleOffset(.7)
91
92# setting title
93p2 = TPaveLabel(.0, .90, .0 + .50, .90 + .10, "Bessel K functions", "NDC")
94p2.SetFillColor(0)
95p2.SetTextFont(22)
96p2.SetTextColor(kBlack)
97
98# setting legend
99leg2.AddEntry(KBessel[0].DrawCopy(), " K_0(x)", "l")
100leg2.AddEntry(KBessel[1].DrawCopy("same"), " K_1(x)", "l")
101leg2.AddEntry(KBessel[2].DrawCopy("same"), " K_2(x)", "l")
102leg2.AddEntry(KBessel[3].DrawCopy("same"), " K_3(x)", "l")
103leg2.AddEntry(KBessel[4].DrawCopy("same"), " K_4(x)", "l")
104leg2.Draw()
105p2.Draw()
106
107# Set canvas 3
108DistCanvas.cd(3)
109gPad.SetGrid()
110gPad.SetFrameFillColor(19)
111leg3 = TLegend(0.75, 0.7, 0.89, 0.89)
112
113# Drawing Bessel i
114iBessel = []
115for nu in range(5):
116 ibessel = TF1("J_0", "ROOT::Math::cyl_bessel_i([0],x)", 0, 10)
117 ibessel.SetParameters(nu, 0.0)
118 ibessel.SetTitle("Bessel I functions")
119 ibessel.SetLineStyle(1)
120 ibessel.SetLineWidth(3)
121 ibessel.SetLineColor(nu + 1)
122 iBessel.append(ibessel)
123
124iaxis = iBessel[0].GetXaxis()
125iaxis.SetTitle("x")
126iaxis.SetTitleSize(0.06)
127iaxis.SetTitleOffset(.7)
128
129# setting title
130p3 = TPaveLabel(.0, .90, .0 + .50, .90 + .10, "Bessel I functions", "NDC")
131p3.SetFillColor(0)
132p3.SetTextFont(22)
133p3.SetTextColor(kBlack)
134
135# setting legend
136leg3.AddEntry(iBessel[0].DrawCopy(), " I_0", "l")
137leg3.AddEntry(iBessel[1].DrawCopy("same"), " I_1(x)", "l")
138leg3.AddEntry(iBessel[2].DrawCopy("same"), " I_2(x)", "l")
139leg3.AddEntry(iBessel[3].DrawCopy("same"), " I_3(x)", "l")
140leg3.AddEntry(iBessel[4].DrawCopy("same"), " I_4(x)", "l")
141leg3.Draw()
142p3.Draw()
143
144# Set canvas 4
145DistCanvas.cd(4)
146gPad.SetGrid()
147gPad.SetFrameFillColor(19)
148leg4 = TLegend(0.75, 0.7, 0.89, 0.89)
149
150# Drawing sph_bessel
151jBessel = []
152for nu in range(5):
153 jbessel = TF1("J_0", "ROOT::Math::sph_bessel([0],x)", 0, 10)
154 jbessel.SetParameters(nu, 0.0)
155 jbessel.SetTitle("Bessel j functions")
156 jbessel.SetLineStyle(1)
157 jbessel.SetLineWidth(3)
158 jbessel.SetLineColor(nu+1)
159 jBessel.append(jbessel)
160jaxis = jBessel[0].GetXaxis()
161jaxis.SetTitle("x")
162jaxis.SetTitleSize(0.06)
163jaxis.SetTitleOffset(.7)
164
165# setting title
166p4 = TPaveLabel(.0, .90, .0 + .50, .90 + .10, "Bessel j functions", "NDC")
167p4.SetFillColor(0)
168p4.SetTextFont(22)
169p4.SetTextColor(kBlack)
170
171# setting legend
172leg4.AddEntry(jBessel[0].DrawCopy(), " j_0(x)", "l")
173leg4.AddEntry(jBessel[1].DrawCopy("same"), " j_1(x)", "l")
174leg4.AddEntry(jBessel[2].DrawCopy("same"), " j_2(x)", "l")
175leg4.AddEntry(jBessel[3].DrawCopy("same"), " j_3(x)", "l")
176leg4.AddEntry(jBessel[4].DrawCopy("same"), " j_4(x)", "l")
177
178leg4.Draw()
179p4.Draw()
180
181DistCanvas.cd()
The Canvas class.
Definition TCanvas.h:23
1-Dim function class
Definition TF1.h:213
This class displays a legend box (TPaveText) containing several legend entries.
Definition TLegend.h:23
A Pave (see TPave) with a text centered in the Pave.
Definition TPaveLabel.h:20