Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
LegendreAssoc.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_math
3/// \notebook
4/// Example describing the usage of different kinds of Associate Legendre Polynomials.
5///
6/// To execute the macro type in:
7///
8/// ~~~{.cpp}
9/// root[0] .x LegendreAssoc.C
10/// ~~~
11///
12/// It draws common graphs for first 5
13/// Associate Legendre Polynomials
14/// and Spherical Associate Legendre Polynomials.
15/// Their integrals on the range [-1, 1] are calculated.
16///
17/// \macro_image
18/// \macro_output
19/// \macro_code
20///
21/// \author Magdalena Slawinska
22
23
24#include "TMath.h"
25#include "TF1.h"
26#include "TCanvas.h"
27
28#include <Riostream.h>
29#include "TLegend.h"
30#include "TLegendEntry.h"
31
32#include "Math/IFunction.h"
33#include <cmath>
34#include "TSystem.h"
35
36void LegendreAssoc()
37{
38 std::cout <<"Drawing associate Legendre Polynomials.." << std::endl;
39 TCanvas *Canvas = new TCanvas("DistCanvas", "Associate Legendre polynomials", 10, 10, 800, 500);
40 Canvas->Divide(2,1);
41 TLegend *leg1 = new TLegend(0.5, 0.7, 0.8, 0.89);
42 TLegend *leg2 = new TLegend(0.5, 0.7, 0.8, 0.89);
43
44 //-------------------------------------------
45 //drawing the set of Legendre functions
46 TF1* L[5];
47
48 L[0]= new TF1("L_0", "ROOT::Math::assoc_legendre(1, 0,x)", -1, 1);
49 L[1]= new TF1("L_1", "ROOT::Math::assoc_legendre(1, 1,x)", -1, 1);
50 L[2]= new TF1("L_2", "ROOT::Math::assoc_legendre(2, 0,x)", -1, 1);
51 L[3]= new TF1("L_3", "ROOT::Math::assoc_legendre(2, 1,x)", -1, 1);
52 L[4]= new TF1("L_4", "ROOT::Math::assoc_legendre(2, 2,x)", -1, 1);
53
54 TF1* SL[5];
55
56 SL[0]= new TF1("SL_0", "ROOT::Math::sph_legendre(1, 0,x)", -TMath::Pi(), TMath::Pi());
57 SL[1]= new TF1("SL_1", "ROOT::Math::sph_legendre(1, 1,x)", -TMath::Pi(), TMath::Pi());
58 SL[2]= new TF1("SL_2", "ROOT::Math::sph_legendre(2, 0,x)", -TMath::Pi(), TMath::Pi());
59 SL[3]= new TF1("SL_3", "ROOT::Math::sph_legendre(2, 1,x)", -TMath::Pi(), TMath::Pi());
60 SL[4]= new TF1("SL_4", "ROOT::Math::sph_legendre(2, 2,x)", -TMath::Pi(), TMath::Pi() );
61
62 Canvas->cd(1);
63 gPad->SetGrid();
64 gPad->SetFillColor(kWhite);
65 L[0]->SetMaximum(3);
66 L[0]->SetMinimum(-2);
67 L[0]->SetTitle("Associate Legendre Polynomials");
68 for (int nu = 0; nu < 5; nu++) {
69 L[nu]->SetLineStyle(kSolid);
70 L[nu]->SetLineWidth(2);
71 L[nu]->SetLineColor(nu+1);
72 }
73
74 leg1->AddEntry(L[0]->DrawCopy(), " P^{1}_{0}(x)", "l");
75 leg1->AddEntry(L[1]->DrawCopy("same"), " P^{1}_{1}(x)", "l");
76 leg1->AddEntry(L[2]->DrawCopy("same"), " P^{2}_{0}(x)", "l");
77 leg1->AddEntry(L[3]->DrawCopy("same"), " P^{2}_{1}(x)", "l");
78 leg1->AddEntry(L[4]->DrawCopy("same"), " P^{2}_{2}(x)", "l");
79 leg1->Draw();
80
81 Canvas->cd(2);
82 gPad->SetGrid();
83 gPad->SetFillColor(kWhite);
84 SL[0]->SetMaximum(1);
85 SL[0]->SetMinimum(-1);
86 SL[0]->SetTitle("Spherical Legendre Polynomials");
87 for (int nu = 0; nu < 5; nu++) {
88 SL[nu]->SetLineStyle(kSolid);
89 SL[nu]->SetLineWidth(2);
90 SL[nu]->SetLineColor(nu+1);
91 }
92
93 leg2->AddEntry(SL[0]->DrawCopy(), " P^{1}_{0}(x)", "l");
94 leg2->AddEntry(SL[1]->DrawCopy("same"), " P^{1}_{1}(x)", "l");
95 leg2->AddEntry(SL[2]->DrawCopy("same"), " P^{2}_{0}(x)", "l");
96 leg2->AddEntry(SL[3]->DrawCopy("same"), " P^{2}_{1}(x)", "l");
97 leg2->AddEntry(SL[4]->DrawCopy("same"), " P^{2}_{2}(x)", "l");
98 leg2->Draw();
99
100
101 //integration
102
103 std::cout << "Calculating integrals of Associate Legendre Polynomials on [-1, 1]" << std::endl;
104 double integral[5];
105 for (int nu = 0; nu < 5; nu++) {
106 integral[nu] = L[nu]->Integral(-1.0, 1.0);
107 std::cout <<"Integral [-1,1] for Associated Legendre Polynomial of Degree " << nu << "\t = \t" << integral[nu] << std::endl;
108 }
109}
@ kWhite
Definition Rtypes.h:66
@ kSolid
Definition TAttLine.h:52
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
#define gPad
The Canvas class.
Definition TCanvas.h:23
TVirtualPad * cd(Int_t subpadnumber=0) override
Set current canvas & pad.
Definition TCanvas.cxx:719
1-Dim function class
Definition TF1.h:182
This class displays a legend box (TPaveText) containing several legend entries.
Definition TLegend.h:23
void Divide(Int_t nx=1, Int_t ny=1, Float_t xmargin=0.01, Float_t ymargin=0.01, Int_t color=0) override
Automatic pad generation by division.
Definition TPad.cxx:1260
virtual void SetGrid(Int_t valuex=1, Int_t valuey=1)=0
RooArgList L(Args_t &&... args)
Definition RooArgList.h:156
constexpr Double_t Pi()
Definition TMath.h:38