% SSB调制解调过程
%% 基本参数
clear all; % 清除变量
close all; % 关闭所有窗口图像
fm = 100; % 基带信号频率
T = 2; % 信号时长
fs = 20000; % 采样频率 奈奎斯特采样定理为最大频率的两倍，这里取20倍为了绘制更多的细节，让时域信号更平滑
dt=1/fs; % 时间采样间隔，采样频率的倒数
N=T/dt; % 采样点个数，总时长除以采样间隔
t=[0:N-1]*dt; % 采样点的时间序列，作为横坐标
%% ******************调制信号时域波形******************
Am=1; % 基带信号幅度
mt=Am*cos(2*pi*fm*t); % 基带信号
figure(1); % 绘制第一幅图
subplot(221); % 窗口分割，将一幅图分割成2*2的
plot(t,mt,'Linewidth',2); % 时间t为横坐标，基带信号mt为纵坐标绘图，线宽为2
xlabel('t/时间'); % 横坐标标注
ylabel('幅度'); % 纵坐标标注
title('基带信号'); % 图标题标注
axis([0,0.1,-1.1,1.1]); % 横纵坐标范围设置
line([0,0.1],[0,0],'color','b','Linewidth',2);% 绘制一条从(0,0)到(0.1,0)的蓝色实线，线宽为2
%% ******************调制信号频域波形******************
subplot(222);
[mf,msf]=T2F(t,mt); % 傅里叶变换，得到纵坐标频谱和横坐标频率
plot(mf,abs(msf),'Linewidth',2); % 画出基带信号频谱，线宽为2
title('基带信号的频谱'); % 图标题标注
xlabel('f/Hz'); % 横坐标标注
ylabel('幅度/H(f)'); % 纵坐标标注
axis([-150 150 -inf inf]); % 横纵坐标范围设置
%% ******************载波信号时域波形******************
subplot(223);
fc=1000; % 载波频率
zaibo=cos(2*pi*fc*t); % 载波时域信号
plot(t,zaibo,'r','Linewidth',2); % 时间t为横坐标，载波信号zaibo为纵坐标绘图，线宽为2，红色
xlabel('t/时间'); % 横坐标标注
ylabel('幅度'); % 纵坐标标注
title('载波信号'); % 图标题标注
axis([0,0.01,-1.1,1.1]); % 横纵坐标范围设置
line([0,0.01],[0,0],'color','b','Linewidth',2);% 绘制一条从(0,0)到(0.01,0)的蓝色实线，线宽为2
%% ******************载波信号频域波形******************
subplot(224);
[mf1,msf1]=T2F(t,zaibo); % 傅里叶变换，得到纵坐标频谱和横坐标频率
plot(mf1,abs(msf1),'r','Linewidth',2); % 载波信号频谱，线宽为2,红色
title('载波信号的频谱'); % 图标题标注
xlabel('f/Hz'); % 横坐标标注
ylabel('幅度/H(f)'); % 纵坐标标注
axis([-1200 1200 -inf inf]); % 横纵坐标范围设置
%% ******************DSB波信号时域波形******************
dsb=mt.*zaibo; % 画出DSB信号波形
figure(2); % 绘制第二幅图
subplot(211); % 窗口分割，将一幅图分割成2*1的
plot(t,dsb,'Linewidth',2); % 画出AM信号波形，线宽为2
title('DSB调制信号'); % 图标题标注
xlabel('t/时间'); % 横坐标标注
ylabel('幅度'); % 纵坐标标注
axis([0,0.02,-3.1,3.1]); % 横纵坐标范围设置
line([0,0.02],[0,0],'color','b','Linewidth',2);% 绘制一条从(0,0)到(0.02,0)的蓝色实线，线宽为2
%% ******************DSB波信号频域波形******************
[mf2,msf2]=T2F(t,dsb); % 傅里叶变换，得到纵坐标频谱和横坐标频率
subplot(212);
plot(mf2,abs(msf2),'Linewidth',2); % 画出DSB信号频谱
title('DSB波信号的频谱'); % 图标题标注
xlabel('f/Hz'); % 横坐标标注
ylabel('幅度/H(f)'); % 纵坐标标注
axis([-1500 1500 -inf inf]); % 横纵坐标范围设置
%% ******************VSB波信号时域波形******************
[samsf1]=lpf_VSB(mf2,msf2,fc); % VSB滤波
[t,vsb]=F2T(mf2,samsf1); %逆傅里叶变换
%% 加噪声
SNR=5; %信噪比为 dB
vsb=awgn(vsb,SNR,'measured');
figure(3); % 绘制第三幅图
subplot(211); % 窗口分割，将一幅图分割成2*1的
plot(t,2*vsb,'Linewidth',2); % 绘制VSB波形
title('SSB信号波形'); % 图标题标注
xlabel('t/时间'); % 横坐标标注
ylabel('幅度'); % 纵坐标标注
axis([0 0.01 -2 2]); % 横纵坐标范围设置
line([0,0.01],[0,0],'color','b','Linewidth',2);% 绘制一条从(0,0)到(0.1,0)的蓝色实线，线宽为2
%% ******************SSB波信号频域波形******************
[mf3,msf3]=T2F(t,vsb); % 傅里叶变换，得到纵坐标频谱和横坐标频率
subplot(212);
plot(mf3,abs(msf3),'Linewidth',2); % 画出SSB信号频谱
title('VSB波信号的频谱'); % 图标题标注
xlabel('f/Hz'); % 横坐标标注
ylabel('幅度/H(f)'); % 纵坐标标注
axis([-1500 1500 -inf inf]); % 横纵坐标范围设置
%% ******************相干解调******************
%% ******************已调信号与载波信号相乘******************
st=vsb.*zaibo; % 已调信号与载波信号相乘
figure(4); % 绘制第四幅图
subplot(211); % 窗口分割，将一幅图分割成2*1的
plot(t,st,'Linewidth',2); % 时间t为横坐标，相乘信号st为纵坐标绘图，线宽为2
title('已调信号与载波信号相乘'); % 图标题标注
xlabel('t/时间'); % 横坐标标注
ylabel('幅度'); % 纵坐标标注
axis([0 0.04 -1.5,1.5]); % 横纵坐标范围设置
line([0,0.04],[0,0],'color','b','Linewidth',2);% 绘制一条从(0,0)到(0.04,0)的蓝色实线，线宽为2
%% ******************已调信号与载波信号相乘的频谱******************
[f,sf]=T2F(t,st); % 傅里叶变换，得到纵坐标频谱和横坐标频率
subplot(212); % 窗口分割，将一幅图分割成2*1的
plot(f,sf,'Linewidth',2); % 绘制相乘信号st的频谱，线宽为2
title('已调信号与载波信号相乘的频谱');% 图标题标注
xlabel('f/Hz'); % 横坐标标注
ylabel('幅度/H(f)'); % 纵坐标标注
axis([-2200 2200 -inf inf]); % 横纵坐标范围设置
%% ******************解调出来的信号******************
[t,st]=lpf(f,sf,2*fm); % 频域低通滤波
figure(5); % 绘制第五幅图
subplot(211); % 窗口分割，将一幅图分割成2*1的
plot(t,4*st,'Linewidth',2); % 绘制解调波形
title('经过低通滤波的相干解调信号波形');% 图标题标注
xlabel('t/时间'); % 横坐标标注
ylabel('幅度'); % 纵坐标标注
axis([0 0.1 -2 2]); % 横纵坐标范围设置
line([0,0.1],[0,0],'color','b','Linewidth',2);% 绘制一条从(0,0)到(0.1,0)的蓝色实线，线宽为2
%% ******************原信号******************
subplot(212);
plot(t,mt,'r-','Linewidth',2); % 绘制原始基带信号
title('原始基带信号'); % 图标题标注
xlabel('t/时间'); % 横坐标标注
ylabel('幅度'); % 纵坐标标注
axis([0 0.1 -1.1 1.1]); % 横纵坐标范围设置
line([0,0.1],[0,0],'color','b','Linewidth',2);% 绘制一条从(0,0)到(0.1,0)的蓝色实线，线宽为2
function [f,sf]= T2F(t,st) % FFT
% dt = t(2)-t(1);
T=t(end); % 输入信号的时间最大值为T
df = 1/T; % dt=1/fs; 时间采样间隔，采样频率的倒数;
% N=T/dt; 采样点个数，总时长除以采样间隔
% 两式联合推导 df = 1/T
N = length(st); % 输入信号时间的长度为采样点数
f=-N/2*df : df : N/2 * df-df; % 频率分布
sf = fft(st); % 做FFT
sf = T/N * fftshift(sf); % 最后输出，将0-fs频谱搬移到-fs/2-fs/2频谱
function[t,st]=F2T(f,Sf) % IFFT
df=f(2)-f(1); % 频率间隔
fmax=(f(end)-f(1)+df); % 最大频率减最低频率加上频率间隔为带宽
dt=1/fmax; % 采样间隔
N=length(f); % 采样点数
t=[0:N-1] * dt; % 时间分布
Sf=fftshift(Sf); % 将0-fs频谱搬移到-fs/2-fs/2频谱
st=fmax * ifft(Sf); % 做IFFT
st=real(st); % 取实部
function[yf]=lpf_VSB(f,sf,B)
df=f(2)-f(1);
fN=length(f);
ym=zeros(1,fN);
xm=floor(B/df);
ym(-xm+floor(fN/2))=0.5;
ym(xm-1+floor(fN/2))=0.5;
for i=1:floor(xm/4)
ym(-xm+floor(fN/2)-i)=0.4;
ym(-xm+floor(fN/2)+i)=0.6;
ym(xm-1+floor(fN/2)-i)=0.6;
ym(xm-1+floor(fN/2)+i)=0.4;
end
yf=ym.*sf;
end
