%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 利用polarplot函数绘制天线二维方向图
% 电流源
theta = -pi: 0.1: pi;
rho = abs(sin(theta));
polarplot(theta, rho); %利用polarplot函数绘制曲线
%偶极子天线
labmda = 1e-3;
k = 2*pi./labmda;
L = labmda./2; %这里计算半波振子的方向图
theta = 0: pi/100: 2*pi;
rho = (cos(k*L/2*cos(theta))-cos(k*L/2))./sin(theta); %偶极子的方向图函数
polarplot(theta, abs(rho), 'b'); %线颜色为蓝色
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 利用surf函数绘制半波振子三维方向图
%方法一：使用for循环
theta = 0: pi/100: pi;
phi = 0: pi/100: 2*pi;
for m = 1: 1: length(theta)
F(m) = cos(pi*cos(theta(m))/2)/sin(theta(m)); %半波振子方向图函数
for n = 1: 1: length(phi)
x(m, n) = F(m)*sin(theta(m))*cos(phi(n));
y(m, n) = F(m)*sin(theta(m))*sin(phi(n));
z(m, n) = F(m)*cos(theta(m));
end
end
surf(x, y, z);
%方法二：省略for循环，加快运行速度
theta = 0: pi/100: pi;
phi = 0: pi/100: 2*pi;
[theta, phi] = meshgrid(theta, phi); %定义网格
F = cos(pi*cos(theta)/2)./sin(theta); %半波振子方向图函数
x = F.*sin(theta).*cos(phi);
y = F.*sin(theta).*sin(phi);
z = F.*cos(theta);
surf(x, y, z );
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 4元端射天线阵因子（沿z轴均匀直线排列）的方向图绘制
%2D
theta = linspace(0, 2*pi); %仰角
phi = linspace(0, 2*pi); %方位角
lambda = 1; %设置波长（频率为0.3GHz）
k0 = 2*pi/lambda;
d = lambda/2;
N = 4;
F = sin(N*k0*d/2*(cos(theta)-1))./(N*sin(k0*d/2*(cos(theta)-1))); %N元端射阵因子
figure
polarplot(theta, F); %绘制二维方向图(H面，包含z轴的平面，theta=0°为端射方向）H面为xoy面
title('2D plot')
%3D方法一（F与phi无关，可将其绕z轴旋转一圈得到3D方向图）
z1 = (F.*cos(theta))'*ones(size(phi)); %绕z轴旋转得到曲面，再计算z值
x1 = (F.*sin(theta))'*cos(phi); %绕z轴旋转得到曲面，再计算x值
y1 = (F.*sin(theta))'*sin(phi); %绕z轴旋转得到曲面，再计算y值，'代表对向量的转置
figure
surf(x1, y1, z1); %绘制三维方向图
title('3D plot')
axis equal %设置各个轴等比例显示
% 3D方法二（利用meshgrid函数）
[theta, phi] = meshgrid(theta, phi); %创建网格
x1 = F.*sin(theta).*cos(phi);
y1 = F.*sin(theta).*sin(phi);
z1 = F.*cos(theta);
figure
surf(x1, y1, z1); %绘制三维方向图
title('3D plot')
axis equal %设置各个轴等比例显示
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 6*6反射阵的远区方向图
% 第一种相位排布方式
freq = 1e+9;
c = 3e+8;
lambda = c/freq;
k = 2*pi/lambda;
m = 6; n= 6; %反射阵的行数和列数
reflect_phi = repmat([pi 0], m, n/2); %将[1 0]重复m行n/2列，构成m*n的矩阵 （反射阵的相位分布情况）
reflect_amp = ones(size(reflect_phi)) ; %取reflect_phi矩阵的行数和列数，再生成相同行列数的全1矩阵（反射阵的幅度分布情况）
d = lambda; %单元间距为波长
M =500; N = 500; %方位角和仰角离散化的个数
theta = linspace(0, pi/2, M);
phi = linspace(0, pi*2, N);
[theta, phi] = meshgrid(theta, phi); %创建网格
F_total = zeros(M, N); %初始方向函数，各个方向均为0
for ii = 1: m
for jj = 1: n
F_total = F_total +reflect_amp(ii, jj).*exp(1i.*(reflect_phi(ii, jj)+k.*d.*sin(theta).*((ii-1/2).*cos(phi)+(jj-1/2).*sin(phi)))); %利用叠加原理，计算阵因子
end
end
Fx = abs(F_total).*sin(theta).*cos(phi);
Fy = abs(F_total).*sin(theta).*sin(phi);
Fz = abs(F_total).*cos(theta);
figure
surf(Fx, Fy, Fz, 'EdgeColor', 'none'); %绘制立体方向图
xlabel('Ex', 'FontSize', 12, 'FontWeight', 'b', 'color', 'r');
ylabel('Ey', 'FontSize', 12, 'FontWeight', 'b', 'color', 'b');
zlabel('Ez', 'FontSize', 12, 'FontWeight', 'b');
xx = sin(theta).*cos(phi);
yy = sin(theta).*sin(phi);
Fzz = abs(F_total);
figure(2);
contourf(yy, xx, Fzz, 'LineStyle','none'); %绘制二维等高线图
xlabel('Ey', 'FontSize', 12, 'FontWeight', 'b', 'color', 'r');
ylabel('Ex', 'FontSize', 12, 'FontWeight', 'b', 'color', 'b');
% 第二种相位排布方式
freq = 1e+9;
c = 3e+8;
lambda = c/freq;
k = 2*pi/lambda;
m = 6; n= 6;
reflect_phi = pi.*repmat(eye(2), m/2, n/2); %将2阶单位阵重复m/2行n/2列，构成m*n矩阵（反射阵的相位分布情况）
reflect_amp = ones(size(reflect_phi)) ; %取reflect_phi矩阵的行数和列数，再生成相同行列数的全1矩阵（反射阵的幅度分布情况）
d = lambda; %单元间距为波长
M =500; N = 500; %方位角和仰角离散化的个数
theta = linspace(0, pi/2, M);
phi = linspace(0, pi*2, N);
[theta, phi] = meshgrid(theta, phi); %创建网格
F_total = zeros(M, N); %初始方向函数，各个方向均为0
for ii = 1: m
for jj = 1: n
F_total = F_total +reflect_amp(ii, jj).*exp(1i.*(reflect_phi(ii, jj)+k.*d.*sin(theta).*((ii-1/2).*cos(phi)+(jj-1/2).*sin(phi)))); %利用叠加原理，计算阵因子
end
end
Fx = abs(F_total).*sin(theta).*cos(phi);
Fy = abs(F_total).*sin(theta).*sin(phi);
Fz = abs(F_total).*cos(theta);
figure
surf(Fx, Fy, Fz, 'EdgeColor', 'none'); %绘制立体方向图
xlabel('Ex', 'FontSize', 12, 'FontWeight', 'b', 'color', 'r');
ylabel('Ey', 'FontSize', 12, 'FontWeight', 'b', 'color', 'b');
zlabel('Ez', 'FontSize', 12, 'FontWeight', 'b');
xx = sin(theta).*cos(phi);
yy = sin(theta).*sin(phi);
Fzz = abs(F_total);
figure(2);
contourf(yy, xx, Fzz, 'LineStyle','none'); %绘制二维等高线图
xlabel('Ey', 'FontSize', 12, 'FontWeight', 'b', 'color', 'r');
ylabel('Ex', 'FontSize', 12, 'FontWeight', 'b', 'color', 'b');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 正弦波_行波传输动态演示
omega = 2*pi;
t = 0;
z = 0: 0.01: 15; %电磁波沿z轴的传输距离
k = 1; %设置波数（相移常数）为1
v = VideoWriter( '正弦波传输.avi');
open(v);
for i = 1: 300 %总帧数为300
y = sin(omega*t-k*z);
plot(z, y);
axis([0 15 -2 2]); %设置坐标轴范围，x轴为0~15，y轴为-2~2
hold on
frame = getframe(gcf); %捕捉当前图像作为一帧
pause(0.1); %波形显示0.1秒（matlab暂停执行0.1秒）
t = t+0.1;
hold off
writeVideo(v, frame);
end
close(v);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% 正弦波_驻波动态演示
omega = 2*pi;
t = 0;
z = 0: 0.01: 30; %电磁波沿z轴的传输距离
k = 1; %设置波数（相移常数）为1
v = VideoWriter( '驻波.avi');
open(v);
for i = 1: 300 %帧数
y1 = sin(omega*t-k*z);
y2 = sin(omega*t+k*z);
y = y1+y2;
plot(z, y1, 'b', z, y2, 'g', z, y, 'r');
axis([0 30 -2.5 2.5]);
hold on
frame = getframe(gcf); %捕捉当前图像作为一帧
pause(0.1); %波形显示0.1秒（matlab暂停执行0.1秒）
t = t+0.01;
hold off
writeVideo(v, frame);
end
close(v);
