This repository has been archived by the owner on Mar 7, 2024. It is now read-only.
-
Notifications
You must be signed in to change notification settings - Fork 1
/
Guo_Spherical_harmonic.m
72 lines (57 loc) · 1.98 KB
/
Guo_Spherical_harmonic.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Spherical Harmonic Function Plot %
% 2021/6 PHY104 Modern Physics %
% Written by Guo Yiming %
% PHY2009481 %
% %
%%%%%%%%%%%%%%%%%%% wait for input ell,ml %%%%%%%%%%%%%%%%%%%%%%%
ell=input('(ell)=');
m=input('(m_l)=');
%%%%%%%%%%%%%%%%%%% Check if the put is valid %%%%%%%%%%%%%%%%%%%%%%
if mod(ell,1)==0
elseif mod(ell,1)~=0 disp('ell not Integer[ FAIL ]'), return
elseif ell==0
end
if ell/abs(ell)==1
elseif ell==0
elseif ell/abs(ell)~=0 disp('ell is Negative[ FAIL ]'), return
end
if mod(m,1)==0
elseif mod(m,1)~=0 disp('m not integer[ FAIL ]'), return
end
if ell >= abs(m)
elseif abs(m) >= ell disp('wtf m>ell[ FAIL ]'), return
end
theta = [0 : 50]' *pi/50;
phi =[-50 : 50] *pi/50;
absm=abs(m);
%%%%%%%%%%%%%%% Processing Data %%%%%%%%%%%%%%%%%%
% Create the surface of Re Ylm(theta,phi)/r^(l+1) = 1;
% Solves for r on a theta, phi grid
gamma = 1/(ell+1);
all = legendre(ell, cos(theta));
if (ell == 0) all=all'; end
Y = all(absm+1, :)' * cos(absm*phi);
[A,B]=size(Y);
r = (abs(Y) .^ gamma); %
%%%%%%%%%%%%%%%% Ploting %%%%%%%%%%%%%%%%%%%%
%Spherical2C
X = r.* (sin(theta)*cos(phi)) ;
Y = r.* (sin(theta)*sin(phi));
Z = r.* (cos(theta)*ones(size(phi)));
C = r.^2;
figure('Name','Spherical Harmonic Function','NumberTitle','off')
surf(X, Y, Z,C,'FaceColor','texturemap','FaceAlpha',0.5,'FaceLighting','flat') %
shading interp
h=camlight('left');
grid on
axis on
axis equal
axis vis3d
ti=['l=' int2str(ell) ', m_l=' int2str(m)];
title(ti);
% Rotate graph
%for i=1:60;
% camorbit(1,0,'data',[0 0 1])
% drawnow
%end