模拟幅度调制相干解调系统抗噪声性能仿真分析

目录

• • 1. 引言
  • 2. 系统模型
  • • 2.1 常规调幅调制`(AM)`
    • 2.2 抑制载波双边带调制`(DSB-SC)`
    • 2.3 单边带调制`(SSB)`
    • • 2.3.1 SSB 调制
      • 2.3.2 SSB 解调
  • 3. 抗干扰性能理论分析
  • 4. 仿真实现与仿真结果
  • • 4.1 常规调幅调制解调仿真`(AM)`
    • 4.2 抑制载波双边带调制解调仿真`(DSB-SC)`
    • 4.3 上边带调制解调仿真`(USSB)`
    • 4.4 下边带调制解调仿真`(LSSB)`
  • 5. 小结
  • 6. 参考资料
  • 7. 附录

1. 引言

1.1 研究目的

• 理解通信系统模型, 掌握信号时域与频域特性分析方法, 能够分析模拟/数字基带/频带等各种通信系统对信号时域以及频域特性的变换关系;
• 理解加性高斯白噪声与频带受限信道, 理解匹配滤波器接收机与相干接收机的工作原理, 掌握带宽无限与频带受限信道条件下传输波形的设计方法;
• 理解模拟通信系统接收机输入与输出信噪比计算方法, 能够用信噪比对模拟通信系统性能进行分析;
• 能够根据系统模型实现链路级仿真, 掌握仿真参数设置原则, 分析信号的时域以及频域特性, 获得误码性能仿真结果;

1.2 研究方法

• 理论分析
• 软件仿真

1.3 主要内容

• 使用理论分析与软件仿真两种方法, 分析比较 AM, DSB-SC 与 SSC 三种系统的抗噪声性能
  • 对原理及模型进行必要的介绍
  • 给出必要的系统框图以及分析结果
  • 完成仿真报告

---

2. 系统模型

2.1 常规调幅调制(AM)

2.1.1 AM 调制

图2.1.1 AM调制器模型

时域表示:
$$\{\begin{array}{ll}{s}_{AM}(t)& =A_c[1+m(t)]\cos (2\pi f_{c}t)\\ {\beta }_{AM}& =max\{|m(t)|\}\end{array}\text{\hspace{1em}(2.1.1)}$$

• $m(t)$ 为基带信号
• $A_c$ 为载波信号的幅值
• $f_c$ 为载波信号的频率
• $s_{AM}(t)$ 为已调AM信号
• ${\beta }_{AM}$ 为调幅指数
  • ${\beta }_{AM}<1$ 时为欠调制
  • ${\beta }_{AM}=1$ 时为临界调制
  • ${\beta }_{AM}>1$ 时为过调制
  • 为避免包络检波时出错, 调制时应确保 ${\beta }_{AM}\le 1$

频域表示:
$$S_{AM}(f)=\frac{A_c}{2}[\delta (f-f_c)+\delta (f+f_c)]+\frac{A_c}{2}[M(f-f_c)+M(f+f_c)]\text{\hspace{1em}(2.1.2)}$$

• $S_{AM}(f)$ 为已调AM信号的傅里叶变换

---

2.1.2 AM 解调

图2.1.2 AM相干解调器模型

解调器输入的已调AM信号 $s(t)$:

时域表示:
$$s(t)=A_c[1+m(t)]\cos (2\pi f_{c}t)\text{\hspace{1em}(2.1.3)}$$

频域表示:
$$S(f)=\frac{A_c}{2}[\delta (f-f_c)+\delta (f+f_c)]+\frac{A_c}{2}[M(f-f_c)+M(f+f_c)]\text{\hspace{1em}(2.1.4)}$$

---

通过带通滤波器(BPF)并移频的AM信号 $s_d(t)$:

时域表示:
$$s_d(t)=\frac{A_c}{2}[1+m(t)]\cdot [1+\cos (2\cdot 2\pi f_{c}t)]\text{\hspace{1em}(2.1.5)}$$

频域表示:
$$S_d(f)=\frac{A_c}{2}[M(f)+\delta (f)]+\frac{A_c}{4}[M(f-2f_c)+M(f+2f_c)+\delta (f-2f_c)+\delta (f+2f_c)]\text{\hspace{1em}(2.1.6)}$$

---

通过低通滤波器(LPF)滤除高频成分的解调信号 $s_o(t)$:

时域表示:
$$s_o(t)=\frac{A_c}{2}+\frac{A_c}{2}m(t)\text{\hspace{1em}(2.1.7)}$$

频域表示:
$$S_o(f)=\frac{A_c}{2}\delta (f)+\frac{A_c}{2}M(f)\text{\hspace{1em}(2.1.8)}$$

所滤除的高频成分 $s_h(t)$:

时域表示:
$$s_h(t)=\frac{A_c}{2}[1+m(t)]\cdot \cos (2\cdot 2\pi f_{c}t)\text{\hspace{1em}(2.1.9)}$$

频域表示:
$$S_h(f)=\frac{A_c}{4}[M(f-2f_c)+M(f+2f_c)+\delta (f-2f_c)+\delta (f+2f_c)]\text{\hspace{1em}(2.1.10)}$$

---

2.2 抑制载波双边带调制(DSB-SC)

2.2.1 DSB_SC 调制

图2.2.1 DSB-SC调制器模型

时域表示:
$$s_{DSB}(t)=A_{c}m(t)\cos (2\pi f_{c}t)\text{\hspace{1em}(2.2.1)}$$

• $m(t)$ 为基带信号
• $A_c$ 为载波信号的幅值
• $f_c$ 为载波信号的频率
• $s_{DSB}(t)$ 为双边带(DSB)信号

频域表示:
$$S_{DSB}(f)=\frac{A_c}{2}[M(f-f_c)+M(f+f_c)]\text{\hspace{1em}(2.2.2)}$$

• $S_{DSB}(f)$ 为双边带(DSB)信号的傅里叶变换

---

2.2.2 DSB_SC 解调

图2.2.2 DSB-SC相干解调器模型

解调器输入的DSB信号 $s(t)$:

时域表示:
$$s(t)=A_{c}m(t)\cos (2\pi f_{c}t)\text{\hspace{1em}(2.2.3)}$$

频域表示:
$$S(f)=\frac{A_c}{2}[M(f-f_c)+M(f+f_c)]\text{\hspace{1em}(2.2.4)}$$

---

通过带通滤波器(BPF)的并移频的DSB信号 $s_d(t)$:

时域表示:
$$s_d(t)=\frac{A_c}{2}m(t)\cdot [1+\cos (2\cdot 2\pi f_{c}t)]\text{\hspace{1em}(2.2.5)}$$

频域表示:
$$S_d(f)=\frac{A_c}{2}M(f)+\frac{A_c}{4}[M(f-2f_c)+M(f+2f_c)]\text{\hspace{1em}(2.2.6)}$$

---

通过低通滤波器(LPF)滤除高频成分的解调信号 $s_o(t)$:

时域表示:
$$s_o(t)=\frac{A_c}{2}m(t)\text{\hspace{1em}(2.2.7)}$$

频域表示:
$$S_o(f)=\frac{A_c}{2}M(f)\text{\hspace{1em}(2.2.8)}$$

所滤除的高频成分 $s_h(t)$:

时域表示:
$$s_h(t)=\frac{A_c}{2}m(t)\cdot \cos (2\cdot 2\pi f_{c}t)\text{\hspace{1em}(2.2.9)}$$

频域表示:
$$S_h(f)=\frac{A_c}{4}[M(f-2f_c)+M(f+2f_c)]\text{\hspace{1em}(2.2.10)}$$

---

2.3 单边带调制(SSB)

2.3.1 SSB 调制

2.3.1.1 滤波法调制

图2.3.1 SSB调制器模型(滤波法)

时域表示:
$$s_{SSB}(t)=A_{c}m(t)\cos (2\pi f_{c}t)\ast h_{SSB}(t)\text{\hspace{1em}(2.3.1)}$$

• $m(t)$ 为基带信号

• $A_c$ 为载波信号的幅值

• $f_c$ 为载波信号的频率

• $s_{SSB}(t)$ 为已调SSB信号

• $h_{SSB}(t)$ 为滤波器冲激响应

  • 若为下边带调制(LSSB), 则滤波器为低通滤波器(LPF), 其上限截止频率为 $f_c$
  • 若为上边带调制(USSB), 则滤波器为高通滤波器(HPF), 其下限截止频率为 $f_c$

• $s_{SSB}(t)$ 为单边带(SSB)信号

频域表示:
$$S_{SSB}(f)=\frac{A_c}{2}[M(f-f_c)+M(f+f_c)]\cdot H_{SSB}\text{\hspace{1em}(2.3.2)}$$

• $M(f)$ 为基带信号的傅里叶变换
• $S_{SSB}(f)$ 为已调SSB信号的傅里叶变换
• $H_{SSB}(f)$ 为滤波器的频域表示
  • 下边带滤波器的频域表示:
    $$H_{LSSB}(f)=\frac{1}{2}[sgn(f+f_c)-sgn(f-f_c)]\text{\hspace{1em}(2.3.3)}$$
  • 上边带滤波器的频域表示:
    $$\begin{array}{rl}{H}_{USSB}(f)& =1-H_{LSSB}(f)\\ & =\frac{1}{2}[sgn(f-f_c)-sgn(f+f_c)]\end{array}\text{\hspace{1em}(2.3.4)}$$

---

2.3.1.2 相移法调制

图2.3.2 SSB调制器模型(相移法)

时域表示:
$$s_{LSSB}(t)=A_c[m(t)\cos (2\pi f_{c}t)+\hat{m}(t)\sin (2\pi f_{c}t)]\text{\hspace{1em}(2.3.5)}$$
$$\begin{array}{rl}{s}_{USSB}(t)& =s_{DSB}(t)-s_{LSSB}(t)\\ & =A_{c}m(t)\cos (2\pi f_{c}t)-A_c[m(t)\cos (2\pi f_{c}t)+\hat{m}(t)\sin (2\pi f_{c}t)]\\ & =A_c[m(t)\cos (2\pi f_{c}t)-\hat{m}(t)\sin (2\pi f_{c}t)]\end{array}\text{\hspace{1em}(2.3.6)}$$

• $m(t)$ 为基带信号
• $\hat{m}(t)$ 为基带信号的希尔伯特变换
  $$\hat{m}(t)=m(t)\ast \frac{1}{\pi t}\text{\hspace{1em}(2.3.7)}$$
• $A_c$ 为载波信号的幅值
• $f_c$ 为载波信号的频率
• $s_{DSB}(t)$ 为双边带(DSB)信号
• $s_{LSSB}(t)$ 为上边带(LSSB)信号
• $s_{USSB}(t)$ 为下边带(USSB)信号

频域表示:
$$S_{LSSB}(f)=\frac{A_c}{2}[M(f-f_c)+M(f+f_c)]+\frac{A_c}{2}[M(f+f_c)\cdot sgn(f+f_c)-M(f+f_c)\cdot sgn(f-f_c)]\text{\hspace{1em}(2.3.8)}$$
$$\begin{array}{rl}{S}_{USSB}& =S_{DSB}(f)-S_{LSSB}(f)\\ & =\frac{A_c}{2}[M(f-f_c)+M(f+f_c)]-\frac{A_c}{2}[M(f+f_c)\cdot sgn(f+f_c)-M(f+f_c)\cdot sgn(f-f_c)]\end{array}\text{\hspace{1em}(2.3.9)}$$

• $S_{DSB}(f)$ 为DSB信号的傅里叶变换
• $S_{LSSB}(f)$ 为LSSB信号的傅里叶变换
• $S_{USSB}(f)$ 为USSB信号的傅里叶变换

---

2.3.2 SSB 解调

图2.3.3 SSB相干解调器模型

解调器输入的SSB信号 $s(t)$:

时域表示:
$$s(t)=A_c[m(t)\cos (2\pi f_{c}t)\pm \hat{m}(t)\sin (2\pi f_{c}t)]\text{\hspace{1em}(2.3.10)}$$

频域表示:
$$S(f)=\frac{A_c}{2}[M(f-f_c)+M(f+f_c)]\pm \frac{A_c}{2}[M(f+f_c)\cdot sgn(f+f_c)-M(f+f_c)\cdot sgn(f-f_c)]\text{\hspace{1em}(2.3.11)}$$

---

通过带通滤波器(BPF)的并移频的SSB信号 $s_d(t)$:

时域表示:
$$s_d(t)=\frac{A_c}{2}m(t)+\frac{A_c}{2}[m(t)\cos (2\cdot 2\pi f_{c}t)\pm \hat{m}(t)\sin (2\cdot 2\pi f_{c}t)]\text{\hspace{1em}(2.3.12)}$$

频域表示:
$$S_d(f)=\frac{A_c}{2}M(f)+\frac{A_c}{4}[M(f-2f_c)+M(f+2f_c)]\pm \frac{A_c}{4}[M(f+2f_c)\cdot sgn(f+2f_c)-M(f+2f_c)\cdot sgn(f-2f_c)]\text{\hspace{1em}(2.2.13)}$$

---

通过低通滤波器(LPF)滤除高频成分的解调信号 $s_o(t)$:

时域表示:
$$s_o(t)=\frac{A_c}{2}m(t)\text{\hspace{1em}(2.2.14)}$$

频域表示:
$$S_o(f)=\frac{A_c}{2}M(f)\text{\hspace{1em}(2.2.15)}$$

所滤除的高频成分 $s_h(t)$:

时域表示:
$$s_h(t)=\frac{A_c}{2}[m(t)\cos (2\cdot 2\pi f_{c}t)\pm \hat{m}(t)\sin (2\cdot 2\pi f_{c}t)]\text{\hspace{1em}(2.2.16)}$$

频域表示:
$$S_h(f)=\frac{A_c}{4}[M(f-2f_c)+M(f+2f_c)]\pm \frac{A_c}{4}[M(f+2f_c)\cdot sgn(f+2f_c)-M(f+2f_c)\cdot sgn(f-2f_c)]\text{\hspace{1em}(2.2.17)}$$

---

3. 抗干扰性能理论分析

图3.1 相干解调器模型

---

3.1 加性高斯白噪声(AWGN)通过相干解调器

由带通AWGN的解析信号
$$\tilde{n}(t)=n_c(t)+j{n}_s(t)\text{\hspace{1em}(3.1.1)}$$

可得带通AWGN的同相-正交形式:
$$\begin{array}{rl}n(t)& =Re\{\tilde{n}(t)e^{j2\pi f_{c}t}\}\\ & =Re\{[n_c(t)+j{n}_s(t)]\cdot [\cos (2\pi f_{c}t)+j\sin (2\pi f_{c}t)]\}\\ & =n_c(t)\cos (2\pi f_{c}t)-n_s(t)\sin (2\pi f_{c}t)\end{array}\text{\hspace{1em}(3.1.2)}$$

带通AWGN通过乘法器移频后
$$\begin{array}{rl}{n}_d(t)& =n(t)\cdot A_c\cos (2\pi f_{c}t)\\ & =n_c(t){\cos }^2(2\pi f_{c}t)-n_s(t)\sin (2\pi f_{c}t)\cos (2\pi f_{c}t)\\ & =\frac{1}{2}{n}_c(t)+\frac{1}{2}[n_c(t)\cos (2\cdot 2\pi f_{c}t)-n_s(t)\sin (2\cdot 2\pi f_{c}t)]\end{array}\text{\hspace{1em}(3.1.3)}$$

相干解调器输出的AWGN成分:
$$n_o(t)=\frac{1}{2}{n}_c(t)\text{\hspace{1em}(3.1.4)}$$

为方便分析, 假设单边带信号带通滤波器的通带带宽为 $B$, 双边带信号带通滤波器的通带带宽为 $2B$, AWGN的单边功率谱密度 $P_{N单}(f)=n_0$, AWGN的双边功率谱密度 $P_{N双}=\frac{n_0}{2}$

双边带调制时输入相干解调器的带通AWGN的功率
$$N_{in双}=E[n^2(t)]=\frac{n_0}{2}\cdot 2\cdot 2B=2n_{0}B\text{\hspace{1em}(3.1.5)}$$

单边带调制时输入相干解调器的带通AWGN的功率
$$N_{in单}=E[n^2(t)]=\frac{n_0}{2}\cdot 2B=n_{0}B\text{\hspace{1em}(3.1.6)}$$

又因为
$$n(t)\sim n_c(t)\sim n_s(t)\sim N(m_0,{\sigma }_0^2)\text{\hspace{1em}(3.1.7)}$$

综合AWGN的特性后可得
$$\{\begin{array}{llccc}E[n(t)]& =E[n_c(t)]& =E[n_s(t)]& =m_0& =0\\ E[n^2(t)]& =E[n_c^2(t)]& =E[n_s^2(t)]& ={\sigma }_0^2& \end{array}\text{\hspace{1em}(3.1.8)}$$

再结合公式 (3.1.4) 可得输出相干解调器得带通AWGN功率

双边带调制时输出相干解调器的带通AWGN的功率:
$$N_{out双}=E[n_o^2(t)]=\frac{1}{4}E[n_c^2(t)]=\frac{1}{4}E[n^2(t)]=\frac{1}{4}{N}_{in双}=\frac{1}{2}{n}_{0}B\text{\hspace{1em}(3.1.9)}$$

单边带调制时输出相干解调器的带通AWGN的功率:
$$N_{out单}=E[n_o^2(t)]=\frac{1}{4}E[n_c^2(t)]=\frac{1}{4}E[n^2(t)]=\frac{1}{4}{N}_{in单}=\frac{1}{4}{n}_{0}B\text{\hspace{1em}(3.1.10)}$$

---

3.2 常规调幅调制解调(AM)

由 (2.1.3) 与归一化功率的定义
$$P=⟨w^2(t)⟩\text{\hspace{1em}(3.2.1)}$$

可得输入相干解调器的AM信号的功率

$$S_{in}=⟨s^2(t)⟩=\frac{A_c^2(1+P_m)}{2}\text{\hspace{1em}(3.2.2)}$$

其中 $P_m$ 为基带信号 $m(t)$ 的归一化平均功率
$$P_m=⟨m^2(t)⟩\text{\hspace{1em}(3.2.3)}$$

再结合公式 (2.1.7) 可得输出相干解调器的AM信号的功率

$$S_{out}=⟨s_o^2(t)⟩=\frac{A_c^2}{4}⟨(1+m(t){)}^2⟩=\frac{A_c^2\cdot (1+P_m)}{4}\text{\hspace{1em}(3.2.4)}$$

结合公式 (3.2.2) 与公式 (3.1.5) 可得到相干解调器输入信号的信噪比

$$SN{R}_{in}=\frac{S_{in}}{N_{in双}}=\frac{\frac{A_c^2}{2}(1+P_m)}{2n_{0}B}=\frac{A_c^2(1+P_m)}{4n_{0}B}\text{\hspace{1em}(3.2.5)}$$

结合公式 (3.2.4) 与公式 (3.1.10) 可得到相干解调器输出信号的信噪比

$$SN{R}_{out}=\frac{S_{out}}{N_{out}}=\frac{\frac{A_c^2\cdot (1+P_m)}{4}}{\frac{1}{2}{n}_{0}B}=\frac{A_c^2\cdot (1+P_m)}{2n_{0}B}\text{\hspace{1em}(3.2.6)}$$

结合公式 (3.2.5) 与公式 (3.2.6) 可得到相干解调器关于AM信号的信噪比增益

$$G_{AM}=\frac{SN{R}_{out}}{SN{R}_{in}}=\frac{2\cdot (1+P_m)}{1+P_m}=2\text{\hspace{1em}(3.2.7)}$$

---

3.3 抑制载波双边带调制解调(DSB-SC)

由 (2.2.1) 与归一化功率的定义
$$P=⟨w^2(t)⟩\text{\hspace{1em}(3.3.1)}$$

可得输入相干解调器的DSB-SC信号的功率

$$S_{in}=⟨s^2(t)⟩=\frac{A_c^2{P}_m}{2}\text{\hspace{1em}(3.3.2)}$$

其中 $P_m$ 为基带信号 $m(t)$ 的归一化平均功率
$$P_m=⟨m^2(t)⟩\text{\hspace{1em}(3.3.3)}$$

再结合公式 (2.2.7) 可得输出相干解调器的DSB-SC信号的功率

$$S_{out}=⟨s_o^2(t)⟩=\frac{A_c^2}{4}⟨m^2(t)⟩=\frac{A_c^2{P}_m}{4}\text{\hspace{1em}(3.3.4)}$$

结合公式 (3.3.2) 与公式 (3.1.5) 可得到相干解调器输入信号的信噪比

$$SN{R}_{in}=\frac{S_{in}}{N_{in}}=\frac{\frac{A_c^2}{2}{P}_m}{2n_{0}B}=\frac{A_c^2{P}_m}{4n_{0}B}\text{\hspace{1em}(3.3.5)}$$

结合公式 (3.3.4) 与公式 (3.1.10) 可得到相干解调器输出信号的信噪比

$$SN{R}_{out}=\frac{S_{out}}{N_{out}}=\frac{\frac{A_c^2{P}_m}{4}}{\frac{1}{2}{n}_{0}B}=\frac{A_c^2{P}_m}{2n_{0}B}\text{\hspace{1em}(3.3.6)}$$

结合公式 (3.3.5) 与公式 (3.3.6) 可得到相干解调器关于DSB-SC信号的信噪比增益

$$G_{DSB}=\frac{SN{R}_{out}}{SN{R}_{in}}=2\text{\hspace{1em}(3.3.7)}$$

---

3.4 单边带调制解调(SSB)

由 (2.3.10) 与归一化功率的定义
$$P=⟨w^2(t)⟩\text{\hspace{1em}(3.4.1)}$$

可得输入相干解调器的DSB-SC信号的功率

$$S_{in}=⟨s^2(t)⟩=A_c^2{P}_m\text{\hspace{1em}(3.4.2)}$$

其中 $P_m$ 为基带信号 $m(t)$ 的归一化平均功率
$$P_m=⟨m^2(t)⟩\text{\hspace{1em}(3.4.3)}$$

再结合公式 (2.3.14) 可得输出相干解调器的DSB-SC信号的功率

$$S_{out}=⟨s_o^2(t)⟩=\frac{A_c^2}{4}⟨m^2(t)⟩=\frac{A_c^2{P}_m}{4}\text{\hspace{1em}(3.4.4)}$$

结合公式 (3.4.2) 与公式 (3.1.6) 可得到相干解调器输入信号的信噪比

$$SN{R}_{in}=\frac{S_{in}}{N_{in}}=\frac{A_c^2{P}_m}{n_{0}B}=\frac{A_c^2{P}_m}{n_{0}B}\text{\hspace{1em}(3.4.5)}$$

结合公式 (3.4.4) 与公式 (3.1.11) 可得到相干解调器输出信号的信噪比

$$SN{R}_{out}=\frac{S_{out}}{N_{out}}=\frac{\frac{A_c^2{P}_m}{4}}{\frac{1}{4}{n}_{0}B}=\frac{A_c^2{P}_m}{n_{0}B}\text{\hspace{1em}(3.4.6)}$$

结合公式 (3.4.5) 与公式 (3.4.6) 可得到相干解调器关于DSB-SC信号的信噪比增益

$$G_{DSB}=\frac{SN{R}_{out}}{SN{R}_{in}}=1\text{\hspace{1em}(3.4.7)}$$

---

4. 仿真实现与仿真结果

4.1 常规调幅调制解调仿真(AM)

仿真程序见 附录7.2 AM调制解调仿真程序 (AM.py)

---

4.1.1 基带信号 $m(t)$, 载波信号 $c(t)$ 与 AM调制信号 $s(t)$ 时域波形

图4.1.1 基带信号, 载波信号与AM调制信号时域波形

---

4.1.2 基带信号 $m(t)$ 频域波形

图4.1.2 基带信号频域波形

---

4.1.3 载波信号 $c(t)$ 频域波形

图4.1.3 载波信号频域波形

---

4.1.4 AM调制信号 $s(t)$ 频域波形

图4.1.4 AM调制信号频域波形

---

4.1.5 AM调制信号 $s(t)$, 同步载波信号 $c(t)$ 与 移频信号 $s_d(t)$ 时域波形

图4.1.5 AM调制信号, 同步载波信号与移频信号时域波形

---

4.1.6 移频信号 $s_d(t)$ 频域波形

图4.1.6 移频信号频域波形

---

4.1.7 解调信号 $m_o(t)$ 时域波形

图4.1.7 解调信号时域波形

---

4.1.8 解调信号 $m_o(t)$ 频域波形

图4.1.8 解调信号频域波形

---

4.1.9 AWGN噪声 $n_i(t)$, 带通AWGN噪声 $n(t)$, 移频带通AWGN噪声 $n_d(t)$ 与输出AWGN噪声 $n_o(t)$ 时域波形

图4.1.9 相干解调器中噪声的时域波形

---

4.1.10 AWGN噪声 $n_i(t)$, 带通AWGN噪声 $n(t)$, 移频带通AWGN噪声 $n_d(t)$ 与输出AWGN噪声 $n_o(t)$ 功率谱密度波形

图4.1.10 相干解调器中噪声的功率谱

---

4.1.11 功率, 信噪比与信噪比增益计算与分析

表4.1 AM调制信号相干解调时功率与信噪比

信号特征	理论值	仿真值	相对误差
输入信号功率 $S_{in}$	5.625000e-01	5.641479e-01	2.92E-03
输入噪声功率 $N_{in}$	1.103856e-02	1.200238e-02	8.03E-02
输入信噪比 $SN{R}_{in}$	5.095772e+01	4.700301e+01	8.41E-02
输出信号功率 $S_{out}$	2.812500e-01	2.901258e-01	3.06E-02
输出噪声功率 $N_{out}$	2.759640e-03	3.321398e-03	1.69E-01
输出信噪比 $SN{R}_{out}$	1.019154e+02	8.735052e+01	1.67E-01
信噪比增益 $G$	2.000000e+00	1.858403e+00	7.62E-02

• 注: 通过BPF前信号的信噪比为6dB.

在误差允许范围内可以认为理论分析符合实际情况.

误差产生原因分析:

1. 相干解调器中的滤波器(BPF, LPF)并非理想滤波器, 因此有一部分通带外的频率成分通过, 其主要导致输入与输出噪声功率仿真值的功率仿真值偏大, 该误差可以通过提高滤波器的阶数减小.
2. 由于抽样信号的频率较低, 因此产生的AWGN并不完全符合随机分布, 且高度依赖随机种子, 不同的随机种子产生的带通AWGN功率有一个较大的散布范围, 其主要影响输入与输出的噪声功率, 间接影响输入与输出的信噪比.
3. 由于AM调制信号仅使用相干解调器解调的结果还保留了直流分量, 因此我使用了修正后的理论值计算公式 (3.2.2), (3.2.4), (3.2.5), (3.2.6) 和 (3.2.7), 现实中对AM调制信号的解调还会使用包络检波器获取其包络, 使用耦合电路过滤其直流分量, 因此在仿真中直流分量对仿真结果又一定的影响.

---

4.2 抑制载波双边带调制解调仿真(DSB-SC)

仿真程序见 附录7.3 DSB-SC调制解调仿真程序 (DSB-SC.py)

---

4.2.1 基带信号 $m(t)$, 载波信号 $c(t)$ 与 AM调制信号 $s(t)$ 时域波形

图4.2.1 基带信号, 载波信号与DSB-SC调制信号时域波形

---

4.2.2 基带信号 $m(t)$ 频域波形

图4.2.2 基带信号频域波形

---

4.2.3 载波信号 $c(t)$ 频域波形

图4.2.3 载波信号频域波形

---

4.2.4 DSB-SC调制信号 $s(t)$ 频域波形

图4.2.4 DSB-SC调制信号频域波形

---

4.2.5 DSB-SC调制信号 $s(t)$, 同步载波信号 $c(t)$ 与 移频信号 $s_d(t)$ 时域波形

图4.2.5 DSB-SC调制信号, 同步载波信号与移频信号时域波形

---

4.2.6 移频信号 $s_d(t)$ 频域波形

图4.2.6 移频信号频域波形

---

4.2.7 解调信号 $m_o(t)$ 时域波形

图4.2.7 解调信号时域波形

---

4.2.8 解调信号 $m_o(t)$ 频域波形

图4.2.8 解调信号频域波形

---

4.2.9 AWGN噪声 $n_i(t)$, 带通AWGN噪声 $n(t)$, 移频带通AWGN噪声 $n_d(t)$ 与输出AWGN噪声 $n_o(t)$ 时域波形

图4.2.9 相干解调器中噪声的时域波形

---

4.2.10 AWGN噪声 $n_i(t)$, 带通AWGN噪声 $n(t)$, 移频带通AWGN噪声 $n_d(t)$ 与输出AWGN噪声 $n_o(t)$ 功率谱密度波形

图4.2.10 相干解调器中噪声的功率谱

---

4.2.11 功率, 信噪比与信噪比增益计算与分析

表4.2 DSB-SC调制信号相干解调时功率与信噪比

信号特征	理论值	仿真值	相对误差
输入信号功率 $S_{in}$	2.500000e-01	2.507324e-01	2.92E-03
输入噪声功率 $N_{in}$	4.906028e-03	5.334388e-03	8.03E-02
输入信噪比 $SN{R}_{in}$	5.095772e+01	4.700304e+01	8.41E-02
输出信号功率 $S_{out}$	1.250000e-01	1.269980e-01	1.57E-02
输出噪声功率 $N_{out}$	1.226507e-03	1.476175e-03	1.69E-01
输出信噪比 $SN{R}_{out}$	1.019153e+02	8.603178e+01	1.85E-01
信噪比增益 $G$	2.000000e+00	1.830345e+00	9.27E-02

• 注: 通过BPF前信号的信噪比为6dB.

在误差允许范围内可以认为理论分析符合实际情况.

误差产生原因分析:

1. 相干解调器中的滤波器(BPF, LPF)并非理想滤波器, 因此有一部分通带外的频率成分通过, 其主要导致输入与输出噪声功率仿真值的功率仿真值偏大, 该误差可以通过提高滤波器的阶数减小.
2. 由于抽样信号的频率较低, 因此产生的AWGN并不完全符合随机分布, 且高度依赖随机种子, 不同的随机种子产生的带通AWGN功率有一个较大的散布范围, 其主要影响输入与输出的噪声功率, 间接影响输入与输出的信噪比.

---

4.3 上边带调制解调仿真(USSB)

仿真程序见 附录7.4 SSB调制解调仿真程序 (SSB.py)

---

4.3.1 基带信号 $m(t)$, 载波信号 $c(t)$ 与 AM调制信号 $s(t)$ 时域波形

图4.3.1 基带信号, 载波信号与USSB调制信号时域波形

---

4.3.2 基带信号 $m(t)$ 频域波形

图4.3.2 基带信号频域波形

---

4.3.3 载波信号 $c(t)$ 频域波形

图4.3.3 载波信号频域波形

---

4.3.4 USSB调制信号 $s(t)$ 频域波形

图4.3.4 USSB调制信号频域波形

---

4.3.5 USSB调制信号 $s(t)$, 同步载波信号 $c(t)$ 与 移频信号 $s_d(t)$ 时域波形

图4.3.5 USSB调制信号, 同步载波信号与移频信号时域波形

---

4.3.6 移频信号 $s_d(t)$ 频域波形

图4.3.6 移频信号频域波形

---

4.3.7 解调信号 $m_o(t)$ 时域波形

图4.3.7 解调信号时域波形

---

4.3.8 解调信号 $m_o(t)$ 频域波形

图4.3.8 解调信号频域波形

---

4.3.9 AWGN噪声 $n_i(t)$, 带通AWGN噪声 $n(t)$, 移频带通AWGN噪声 $n_d(t)$ 与输出AWGN噪声 $n_o(t)$ 时域波形

图4.3.9 相干解调器中噪声的时域波形

---

4.3.10 AWGN噪声 $n_i(t)$, 带通AWGN噪声 $n(t)$, 移频带通AWGN噪声 $n_d(t)$ 与输出AWGN噪声 $n_o(t)$ 功率谱密度波形

图4.3.10 相干解调器中噪声的功率谱

---

4.3.11 功率, 信噪比与信噪比增益计算与误差分析

表4.3 USSB调制信号相干解调时功率与信噪比

信号特征	理论值	仿真值	相对误差
输入信号功率 $S_{in}$	5.000000e-01	5.009397e-01	1.88E-03
输入噪声功率 $N_{in}$	4.906028e-03	4.968031e-03	1.25E-02
输入信噪比 $SN{R}_{in}$	1.019154e+02	1.008327e+02	1.07E-02
输出信号功率 $S_{out}$	1.250000e-01	1.265897e-01	1.26E-02
输出噪声功率 $N_{out}$	1.226507e-03	1.134983e-03	8.06E-02
输出信噪比 $SN{R}_{out}$	1.019153e+02	1.115344e+02	8.62E-02
信噪比增益 $G$	1.000000e+00	1.106134e+00	9.60E-02

• 注: 通过BPF前信号的信噪比为6dB.

在误差允许范围内可以认为理论分析符合实际情况.

误差产生原因分析:

1. 相干解调器中的滤波器(BPF, LPF)并非理想滤波器, 因此有一部分通带外的频率成分通过, 其主要导致输入与输出噪声功率仿真值的功率仿真值偏大, 该误差可以通过提高滤波器的阶数减小.
2. 由于抽样信号的频率较低, 因此产生的AWGN并不完全符合随机分布, 且高度依赖随机种子, 不同的随机种子产生的带通AWGN功率有一个较大的散布范围, 其主要影响输入与输出的噪声功率, 间接影响输入与输出的信噪比.
3. 还是由于产生的AWGN并不完全符合随机分布, 因此相同的AWGN通过USSB的带通滤波器与通过LSSB的带通滤波器后的功率并不相同, 且其同样高度依赖随机种子, 不同的随机种子对USSB与LSSB的噪声功率之差也有决定性影响.

---

4.4 下边带调制解调仿真(LSSB)

仿真程序见 附录7.4 SSB调制解调仿真程序 (SSB.py)

---

4.4.1 基带信号 $m(t)$, 载波信号 $c(t)$ 与 AM调制信号 $s(t)$ 时域波形

图4.4.1 基带信号, 载波信号与LSSB调制信号时域波形

---

4.4.2 基带信号 $m(t)$ 频域波形

图4.4.2 基带信号频域波形

---

4.4.3 载波信号 $c(t)$ 频域波形

图4.4.3 载波信号频域波形

---

4.4.4 LSSB调制信号 $s(t)$ 频域波形

图4.4.4 LSSB调制信号频域波形

---

4.4.5 LSSB调制信号 $s(t)$, 同步载波信号 $c(t)$ 与 移频信号 $s_d(t)$ 时域波形

图4.4.5 LSSB调制信号, 同步载波信号与移频信号时域波形

---

4.4.6 移频信号 $s_d(t)$ 频域波形

图4.4.6 移频信号频域波形

---

4.4.7 解调信号 $m_o(t)$ 时域波形

图4.4.7 解调信号时域波形

---

4.4.8 解调信号 $m_o(t)$ 频域波形

图4.4.8 解调信号频域波形

---

4.4.9 AWGN噪声 $n_i(t)$, 带通AWGN噪声 $n(t)$, 移频带通AWGN噪声 $n_d(t)$ 与输出AWGN噪声 $n_o(t)$ 时域波形

图4.4.9 相干解调器中噪声的时域波形

---

4.4.10 AWGN噪声 $n_i(t)$, 带通AWGN噪声 $n(t)$, 移频带通AWGN噪声 $n_d(t)$ 与输出AWGN噪声 $n_o(t)$ 功率谱密度波形

图4.4.10 相干解调器中噪声的功率谱

---

4.4.11 功率, 信噪比与信噪比增益计算与误差分析

表4.3 LSSB调制信号相干解调时功率与信噪比

信号特征	理论值	仿真值	相对误差
输入信号功率 $S_{in}$	5.000000e-01	5.010163e-01	2.03E-03
输入噪声功率 $N_{in}$	4.906028e-03	5.601459e-03	1.24E-01
输入信噪比 $SN{R}_{in}$	1.019154e+02	8.944388e+01	1.39E-01
输出信号功率 $S_{out}$	1.250000e-01	1.274379e-01	1.91E-02
输出噪声功率 $N_{out}$	1.226507e-03	1.353658e-03	9.39E-02
输出信噪比 $SN{R}_{out}$	1.019153e+02	9.414334e+01	8.26E-02
信噪比增益 $G$	1.000000e+00	1.052541e+00	4.99E-02

• 注: 通过BPF前信号的信噪比为6dB.

在误差允许范围内可以认为理论分析符合实际情况.

误差产生原因分析:

1. 相干解调器中的滤波器(BPF, LPF)并非理想滤波器, 因此有一部分通带外的频率成分通过, 其主要导致输入与输出噪声功率仿真值的功率仿真值偏大, 该误差可以通过提高滤波器的阶数减小.
2. 由于抽样信号的频率较低, 因此产生的AWGN并不完全符合随机分布, 且高度依赖随机种子, 不同的随机种子产生的带通AWGN功率有一个较大的散布范围, 其主要影响输入与输出的噪声功率, 间接影响输入与输出的信噪比.
3. 还是由于产生的AWGN并不完全符合随机分布, 因此相同的AWGN通过USSB的带通滤波器与通过LSSB的带通滤波器后的功率并不相同, 且其同样高度依赖随机种子, 不同的随机种子对USSB与LSSB的噪声功率之差也有决定性影响.

---

5. 小结

表5.1 不同调制信号的时域表达式

信号时域表达式	AM	DSB-SC	USSB	LSSB
基带信号 $m(t)$	$m(t)$	$m(t)$	$m(t)$	$m(t)$
调制信号 $s(t)$	$A_c[1+m(t)]\cos (2\pi f_{c}t)$	$A_{c}m(t)\cos (2\pi f_{c}t)$	$A_c[m(t)\cos (2\pi f_{c}t)-\hat{m}(t)\sin (2\pi f_{c}t)]$	$A_c[m(t)\cos (2\pi f_{c}t)+\hat{m}(t)\sin (2\pi f_{c}t)]$
移频后调制信号 $s_d(t)$	$\frac{A_c}{2}[1+m(t)]\cdot [1+\cos (2\cdot 2\pi f_{c}t)]$	$\frac{A_c}{2}m(t)\cdot [1+\cos (2\cdot 2\pi f_{c}t)]$	$\frac{A_c}{2}m(t)+\frac{A_c}{2}[m(t)\cos (2\cdot 2\pi f_{c}t)-\hat{m}(t)\sin (2\cdot 2\pi f_{c}t)]$	$\frac{A_c}{2}m(t)+\frac{A_c}{2}[m(t)\cos (2\cdot 2\pi f_{c}t)+\hat{m}(t)\sin (2\cdot 2\pi f_{c}t)]$
解调信号 $m_o(t)$	$\frac{A_c}{2}m(t)$	$\frac{A_c}{2}m(t)$	$\frac{A_c}{2}m(t)$	$\frac{A_c}{2}m(t)$

• 注: AM调制信号解调后未消除直流分量

表5.2 不同调制信号的频域表达式

信号频域表达式	AM	DSB-SC	USSB	LSSB
基带信号频域表达式 $M(f)$	$M(f)$	$M(f)$	$M(f)$	$M(f)$
调制信号频域表达式 $S(f)$	$\frac{A_c}{2}[\delta (f-f_c)+\delta (f+f_c)]+\frac{A_c}{2}[M(f-f_c)+M(f+f_c)]$	$\frac{A_c}{2}[M(f-f_c)+M(f+f_c)]$	$\frac{A_c}{2}[M(f-f_c)+M(f+f_c)]-\frac{A_c}{2}[M(f+f_c)\cdot sgn(f+f_c)-M(f+f_c)\cdot sgn(f-f_c)]$	$\frac{A_c}{2}[M(f-f_c)+M(f+f_c)]+\frac{A_c}{2}[M(f+f_c)\cdot sgn(f+f_c)-M(f+f_c)\cdot sgn(f-f_c)]$
移频后调制信号频域表达式 $S_d(f)$	$\frac{A_c}{2}[M(f)+\delta (f)]+\frac{A_c}{4}[M(f-2f_c)+M(f+2f_c)+\delta (f-2f_c)+\delta (f+2f_c)]$	$\frac{A_c}{2}M(f)+\frac{A_c}{4}[M(f-2f_c)+M(f+2f_c)]$	$\frac{A_c}{2}M(f)+\frac{A_c}{4}[M(f-2f_c)+M(f+2f_c)]-\frac{A_c}{4}[M(f+2f_c)\cdot sgn(f+2f_c)-M(f+2f_c)\cdot sgn(f-2f_c)]$	$\frac{A_c}{2}M(f)+\frac{A_c}{4}[M(f-2f_c)+M(f+2f_c)]+\frac{A_c}{4}[M(f+2f_c)\cdot sgn(f+2f_c)-M(f+2f_c)\cdot sgn(f-2f_c)]$
解调信号频域表达式 $M_o(f)$	$\frac{A_c}{2}\delta (f)+\frac{A_c}{2}M(f)$	$\frac{A_c}{2}M(f)$	$\frac{A_c}{2}M(f)$	$\frac{A_c}{2}M(f)$

• 注: AM调制信号解调后未消除直流分量

表5.3 相干解调器解调不同调制信号的抗噪能力对比

信号特征	AM	DSB-SC	SSB
带宽 $BW$	$2B$	$2B$	$B$
输入信号功率 $S_{in}$	$\frac{A_c^2(1+P_m)}{2}$	$\frac{A_c^2{P}_m}{2}$	$A_c^2{P}_m$
输入噪声功率 $N_{in}$	$2n_{0}B$	$2n_{0}B$	$n_{0}B$
输入信噪比 $SN{R}_{in}$	$\frac{A_c^2(1+P_m)}{4n_{0}B}$	$\frac{A_c^2{P}_m}{4n_{0}B}$	$\frac{A_c^2{P}_m}{n_{0}B}$
输出信号功率 $S_{out}$	$\frac{A_c^2\cdot (1+P_m)}{4}$	$\frac{A_c^2{P}_m}{4}$	$\frac{A_c^2{P}_m}{4}$
输出噪声功率 $N_{out}$	$\frac{1}{2}{n}_{0}B$	$\frac{1}{2}{n}_{0}B$	$\frac{1}{4}{n}_{0}B$
输出信噪比 $SN{R}_{out}$	$\frac{A_c^2\cdot (1+P_m)}{2n_{0}B}$	$\frac{A_c^2{P}_m}{2n_{0}B}$	$\frac{A_c^2{P}_m}{n_{0}B}$
信噪比增益 $G$	$2$	$2$	$1$

• 注: AM调制信号解调后未消除直流分量

---

6. 参考资料

[1] Markdown 公式编号.
[2] markdown让文字居中和带颜色.
[3] LaTeX公式手册(全网最全).
[4] python signal滤波器使用说明.
[5] Python基于scipy实现信号滤波功能.
[7] Python 加性高斯白噪声 AWGN.
[8] 高斯白噪声 python.
[9] 频域特征提取的Python实现（频谱、功率谱、倒频谱）.
[10] Python求均值，方差，标准差.
[11] python基础_格式化输出.
[12] 现代通信原理6.1 常规调幅调制(AM)与抑制载波双边带(DSB-SC)调制.
[13] 现代通信原理6.2：单边带(SSB)调制.
[14] 仿真作业3：噪声通过DSB-SC解调器.
[15] MATLAB通信仿真实例1：无噪声信道下DSB-SC调制解调器.
[16] 现代通信原理A.3：随机信号的功率谱估计.
[17] 现代通信原理4.3：白噪声.
[18] 模拟幅度调制相干解调系统抗噪声性能仿真分析.
[19] 模拟幅度调制相干解调系统抗噪声性能仿真分析.
[20] 现代通信原理 - 仿真2 - DSB-SC(双边带抑制载波)调制解调器的仿真.

---

7. 附录

7.1 实用程序 utils.py

# -*- coding: utf-8 -*-
# utils.py

# %% 导入库
import numpy as np
import matplotlib.pyplot as plt
# help(plt.psd)s


# %% 加性高斯白噪声生成函数
def wgn(x, snr, seed=7):
    '''
    加入高斯白噪声 Additive White Gaussian Noise
    :param x: 原始信号
    :param snr: 信噪比
    :param seed: 随机种子
    :return: 加入噪声后的信号
    '''
    np.random.seed(seed)  # 设置随机种子
    snr = 10 ** (snr / 10.0)
    xpower = np.sum(x ** 2) / len(x)
    npower = xpower / snr
    return np.random.randn(len(x)) * np.sqrt(npower)


# %% 波形绘制函数
def drawWave(
        figure_i,
        row_i,
        column_j,
        label_x,
        axis_xs,
        axis_ys,
        formats,
        labels,
        label_locs,
        limit_ys,
        linewidths):
    """
    param:
        figure_i: 图窗编号
        row_i: 图表的行数
        column_j: 图表的列数
        label_x: 横轴标签
        axis_xs: 横轴
        axis_ys: 纵轴
        formats: 曲线格式(颜色 风格 标记)
            颜色:
                'b' 蓝色(blue)(默认)
                'g' 绿色(green)
                'r' 红色(red)
                'y' 黄色(yellow)
                'k' 黑色(black)
                'w' 白色(white)
                'c' 青绿色(cyan)
                'm' 洋红色(magenta)
                '#xxxxxx'   RGB颜色
                '0.xxx' 灰度
            风格:
                '-' 实线(默认)
                ':' 虚线
                ' ' 无线条
                '--'    破折线
                '-.'    点划线
            标记:
                '.' 点(默认)
                ',' 极小点
                'o' 实心圆
                '^' 正三角
                'v' 倒三角
                '>' 右三角
                '<' 左三角
                '1' 下花三角
                '2' 上花三角
                '3' 左花三角
                '4' 右花三角
                's' 实心方形(square)
                'p' 实心五边形(pentagon)
                '*' 星形
                'h' 竖六边形(hexagon)
                'H' 横六边形
                '+' 十字
                'x' 叉形
                'D' 菱形(diamond)
                'd' 瘦菱形
                '|' 垂直线
        labels: 曲线标签内容
        label_locs: 曲线标签位置
            0 "best 最佳位置(默认)
            1 "upper right" 右上
            2 "upper left" 左上
            3 "lower right" 右下
            4 "lower left" 左下
        limit_ys: 坐标轴下限与上限
        linewidths: 线宽
    """

    fig = plt.figure(figure_i)
    chart_formate = row_i * 100 + column_j * 10
    for i in range(row_i * column_j):
        ax = fig.add_subplot(chart_formate + 1 + i)
        ax.plot(
            axis_xs[i],
            axis_ys[i],
            formats[i],
            label=labels[i],
            linewidth=linewidths[i]
        )  # 横轴, 纵轴, 标签
        ax.legend(loc=label_locs[i])  # 标签位置
        if limit_ys[i] is not None:
            ax.set_ylim(limit_ys[i])  # 纵轴刻度范围
        ax.grid()  # 显示网格线
    ax.set_xlabel(label_x)


# %% 功率密度谱绘制函数
def drawPSD(
        figure_i,
        row_i,
        column_j,
        label_x,
        axis_ys,
        F_samples,
        limit_ys):
    """
    param:
        figure_i: 图窗编号
        row_i: 图表的行数
        column_j: 图表的列数
        label_x: 横轴标签
        axis_ys: 纵轴
        Fs_s:
        limit_ys: 坐标轴下限与上限
    """

    fig = plt.figure(figure_i)
    chart_formate = row_i * 100 + column_j * 10
    for i in range(row_i * column_j):
        ax = fig.add_subplot(chart_formate + 1 + i)
        ax.psd(axis_ys[i], Fs=F_samples[i])
        if limit_ys[i] is not None:
            ax.set_ylim(limit_ys[i])  # 纵轴刻度范围
        ax.grid()  # 显示网格线
    ax.set_xlabel(label_x)


# %% 显示所有图窗
def show():
    plt.show()

---

7.2 AM调制解调仿真程序 AM.py

# -*- coding: utf-8 -*-
# AM.PY
# %% 配置导入
from scipy import signal
from scipy import fft
import numpy as np
import utils


# %% 参数设置
time_begin = 0  # 起始时间
time_end = 1  # 终止时间
points = 1024  # 栅网采样点数(采样频率)
f_sample = points / (time_end - time_begin)  # 采样频率
T = (time_end - time_begin) / (points - 1)  # 栅网采样时间间距
frequency = [10, 100]  # 频率
amplitude = [0.5, 1]  # 调幅系数/振幅
phase = [0, 0]  # 相位(cos(t)相位为0)
zoom = 4  # 频谱缩放倍数
display_interval = 1.25  # 显示区间的系数
tick_max = max(amplitude) * 2  # 波形纵轴刻度最大值
tick_min = -tick_max  # 波形纵轴刻度最小值

snr = 6  # AWGN生成时的信噪比, 单位(dB)
seed = 3  # AWGN发生器随机种子

filter_order = 8  # 滤波器的阶数
f_H_LPF = frequency[0] * 2  # 低通滤波器的上限截止频率
f_H_BPF_DSB = frequency[1] + frequency[0] * 2  # 带通滤波器的上限截止频率
f_L_BPF_DSB = frequency[1] - frequency[0] * 2  # 带通滤波器的下限截止频率


# %% 生成时域波形
t = np.linspace(time_begin, time_end, points)  # 生成时域栅网
waves = []
for i in range(len(frequency)):
    waves.append(
        signal.chirp(
            t,
            f0=frequency[i],
            t1=time_end,
            f1=frequency[i],
            phi=phase[i],
            method='linear'
        ) * amplitude[i]
    )
s = None
for wave in waves:
    if s is None:
        s = 1 + wave.copy()
    else:
        s *= wave
s_d = s * waves[1]


# %% 获得频域波形
f = np.linspace(0, 1 / T, points)  # 生成频域栅网
# 输入信号 M(f)
M = fft.fft(waves[0])  # 快速傅里叶变换
M_nor = M / points  # 归一化
M_mod = np.abs(M_nor)  # 获取幅度频谱
M_pha = np.angle(M)  # 获取相位频谱
# 移频信号 C(f)
C = fft.fft(waves[1])  # 快速傅里叶变换
C_nor = C / points  # 归一化
C_mod = np.abs(C_nor)  # 获取幅度频谱
C_pha = np.angle(C)  # 获取相位频谱
# 第一次移频 -> S(f)
S = fft.fft(s)  # 快速傅里叶变换
S_nor = S / points  # 归一化
S_mod = np.abs(S_nor)  # 获取幅度频谱
S_pha = np.angle(S)  # 获取相位频谱
# 第二次移频 -> S_d(f)
S_d = fft.fft(s_d)
S_d_nor = S_d / points  # 归一化
S_d_mod = np.abs(S_d_nor)  # 获取幅度频谱
S_d_pha = np.angle(S_d)  # 获取相位频谱


# %% 滤波器参数
"""
    滤波器构造函数(仅介绍Butterworth滤波器)
        scipy.signal.butter(N, Wn, btype='low', analog=False, output='ba')
    输入参数：
        N:滤波器的阶数
        Wn：归一化截止频率。计算公式Wn=2*截止频率/采样频率。（注意：根据采样定理，采样频率要大于两倍的信号本身最大的频率，才能还原信号。截止频率一定小于信号本身最大的频率，所以Wn一定在0和1之间）。当构造带通滤波器或者带阻滤波器时，Wn为长度为2的列表。
        btype : 滤波器类型{‘lowpass’, ‘highpass’, ‘bandpass’, ‘bandstop’},
        output : 输出类型{‘ba’, ‘zpk’, ‘sos’},
    输出参数：
        b，a: IIR滤波器的分子（b）和分母（a）多项式系数向量。output='ba'
        z,p,k: IIR滤波器传递函数的零点、极点和系统增益. output= 'zpk'
        sos: IIR滤波器的二阶截面表示。output= 'sos'
    函数的使用
        信号滤波中最常用的无非低通滤波、高通滤波和带通滤波。下面简单介绍这三种滤波的使用过程：
    (1).高通滤波
        # 这里假设采样频率为1000hz,信号本身最大的频率为500hz，要滤除10hz以下频率成分，即截至频率为10hz，则wn=2*10/1000=0.02
        # from scipy import signal
        # b, a = signal.butter(8, 0.02, 'highpass')
        # filtedData = signal.filtfilt(b, a, data)#data为要过滤的信号
    (2).低通滤波
        # 这里假设采样频率为1000hz,信号本身最大的频率为500hz，要滤除10hz以上频率成分，即截至频率为10hz，则wn=2*10/1000=0.02
        # from scipy import signal
        # b, a = signal.butter(8, 0.02, 'lowpass')
        # filtedData = signal.filtfilt(b, a, data)       #data为要过滤的信号
    (3).带通滤波
        # 这里假设采样频率为1000hz,信号本身最大的频率为500hz，要滤除10hz以下和400hz以上频率成分，即截至频率为10hz和400hz,则wn1=2*10/1000=0.02,wn2=2*400/1000=0.8。Wn=[0.02,0.8]
        # from scipy import signal
        # b, a = signal.butter(8, [0.02,0.8], 'bandpass')
        # filtedData = signal.filtfilt(b, a, data)   #data为要过滤的信号
"""
wn_LPF = 2 * f_H_LPF / f_sample
b_LPF, a_LPF = signal.butter(filter_order, wn_LPF, 'lowpass')  # 配置滤波器

wn_BPF_DSB = [2 * f_L_BPF_DSB / f_sample, 2 * f_H_BPF_DSB / f_sample]
b_BPF_DSB, a_BPF_DSB = signal.butter(filter_order, wn_BPF_DSB, 'bandpass')  # 配置滤波器


# %% 进行滤波, 获得 m_o(t) 的时域波形
"""
    滤波函数
        scipy.signal.filtfilt(b, a, x, axis=-1, padtype='odd', padlen=None, method='pad', irlen=None)
    输入参数：
        b: 滤波器的分子系数向量
        a: 滤波器的分母系数向量
        x: 要过滤的数据数组。（array型）
        axis: 指定要过滤的数据数组x的轴
        padtype: 必须是“奇数”、“偶数”、“常数”或“无”。这决定了用于过滤器应用的填充信号的扩展类型。{‘odd’, ‘even’, ‘constant’, None}
        padlen：在应用滤波器之前在轴两端延伸X的元素数目。此值必须小于要滤波元素个数 - 1。（int型或None）
        method：确定处理信号边缘的方法。当method为“pad”时，填充信号；填充类型padtype和padlen决定，irlen被忽略。当method为“gust”时，使用古斯塔夫森方法，而忽略padtype和padlen。{“pad” ，“gust”}
        irlen：当method为“gust”时，irlen指定滤波器的脉冲响应的长度。如果irlen是None，则脉冲响应的任何部分都被忽略。对于长信号，指定irlen可以显著改善滤波器的性能。（int型或None）
    输出参数：
        y: 滤波后的数据数组
"""
m_o = signal.filtfilt(b_LPF, a_LPF, s_d)  # 为要过滤的信号
t_m_original = np.linspace(time_begin, T * (len(m_o) - 1), len(m_o))  # m_o(t) 信号的栅网


# %% 获得 m_o(t) 的频域波形
M_o = fft.fft(m_o)
M_o_nor = M_o / points  # 归一化
M_o_mod = np.abs(M_o_nor)  # 获取幅度频谱
M_o_pha = np.angle(M_o)  # 获取相位频谱


# %% AWGN生成与处理
n_i = utils.wgn(s, snr, seed)  # 生成AWGN
n = signal.filtfilt(b_BPF_DSB, a_BPF_DSB, n_i)  # 获得带通AWGN
n_d = n * waves[1]  # 获得移频后的带通AWGN
n_o = signal.filtfilt(b_LPF, a_LPF, n_d)  # 获得通过低通滤波器的带通AWGN


# %% 绘制 m(t) c(t) s(t) 时域波形
utils.drawWave(
    figure_i=1,
    row_i=3,
    column_j=1,
    label_x=r"t / s",
    axis_xs=[t] * 3,
    axis_ys=waves + [s],
    formats=['b-,', 'g-,', 'r-,'],
    labels=[r"m(t) = 0.5cos(20πt)",
            r"c(t) = cos(200πt)",
            r"s(t) = [1 + m(t)]c(t)"],
    label_locs=[r"upper right"] * 3,
    limit_ys=[(tick_min, tick_max)] * 3,
    linewidths=[1] * 3
)


# %% 绘制 m(t) 频域波形
utils.drawWave(
    figure_i=2,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 + [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((M_mod[-1:0:-1], M_mod)),
             np.hstack((-M_pha[-1:0:-1], M_pha)),
             M_mod[0: points // zoom],
             M_pha[0: points // zoom]],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|M(f)|", r"θᴍ(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(M_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制 c(t) 频域波形
utils.drawWave(
    figure_i=3,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 + [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((C_mod[-1:0:-1], C_mod)),
             np.hstack((-C_pha[-1:0:-1], C_pha)),
             C_mod[0: points // zoom],
             C_pha[0: points // zoom]],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|C(f)|", r"θᴄ(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(C_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制 s(t) 频域波形
utils.drawWave(
    figure_i=4,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 + [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((S_mod[-1:0:-1], S_mod)),
             np.hstack((-S_pha[-1:0:-1], S_pha)),
             S_mod[0: points // zoom],
             S_pha[0: points // zoom]],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|S(f)|", r"θs(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(S_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制 s(t) c(t) s_d(t) 时域波形
utils.drawWave(
    figure_i=5,
    row_i=3,
    column_j=1,
    label_x=r"t / s",
    axis_xs=[t] * 3,
    axis_ys=[s, waves[1], s_d],
    formats=['b-,', 'g-,', 'r-,'],
    labels=[r"s(t) = m(t)c(t)",
            r"c(t) = cos(200πt)",
            r"s_d(t) = s(t)c(t)"],
    label_locs=[r"upper right"] * 3,
    limit_ys=[(tick_min, tick_max)] * 3,
    linewidths=[1] * 3
)


# %% 绘制 s_d(t) 频域波形
utils.drawWave(
    figure_i=6,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 + [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((S_d_mod[-1:0:-1], S_d_mod)),
             np.hstack((-S_d_pha[-1:0:-1], S_d_pha)),
             S_d_mod[0: points // zoom],
             S_d_pha[0: points // zoom]],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|Sᴅ(f)|", r"θsᴅ(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(S_d_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制 m_o(t) 时域波形
utils.drawWave(
    figure_i=9,
    row_i=2,
    column_j=1,
    label_x=r"t / s",
    axis_xs=[t_m_original, t],
    axis_ys=[m_o, waves[0]],
    formats=['b-,', 'g-,'],
    labels=[r"m_o(t) = s_d(t) * h(t)", r"m(t)"],
    label_locs=[r"upper right"] * 2,
    limit_ys=[(tick_min, tick_max)] * 2,
    linewidths=[1] * 2
)


# %% 绘制 m_o(t) 频域波形
utils.drawWave(
    figure_i=10,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 +
    [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((M_o_mod[-1:0:-1], M_o_mod)),
             np.hstack((-M_o_pha[-1:0:-1], M_o_pha)),
             M_o_mod[0: points // zoom],
             M_o_pha[0: points // zoom]],
    # axis_xs=[f_h] * 2,
    # axis_ys=[H_mod,
    #          H_pha],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|Mo(f)|", r"θᴍo(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(M_o_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制噪声 n_i(t), n(t), n_d(t), n_o(t) 时域波形
utils.drawWave(
    figure_i=11,
    row_i=4,
    column_j=1,
    label_x=r"t / s",
    axis_xs=[t] * 4,
    axis_ys=[n_i, n, n_d, n_o],
    formats=['b-,', 'g-,', 'r-,', 'k-,'],
    labels=[r"n_i(t)",
            r"n(t)",
            r"n_d(t)",
            r"n_o(t)"],
    label_locs=[r"upper right"] * 4,
    limit_ys=[None] * 4,
    linewidths=[1] * 4
)


# %% 绘制噪声 n_i(t), n(t), n(t), n_o(t) 功率谱密度波形
utils.drawPSD(
    figure_i=12,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_ys=[n_i, n, n_d, n_o],
    F_samples=[f_sample] * 4,
    limit_ys=[None] * 4)


# %% 求功率, 信噪比与信噪比增益
S_in = np.mean(s ** 2)
N_in = np.mean(n ** 2)
SNR_in = S_in / N_in

S_out = np.mean(m_o ** 2)
N_out = np.mean(n_o ** 2)
SNR_out = S_out / N_out

G = SNR_out / SNR_in

print("S_in = %e\nN_in = %e\nSNR_in = %e\nS_out = %e\nN_out = %e\nSNR_out = %e\nG = %e\n" % (S_in, N_in, SNR_in, S_out, N_out, SNR_out, G))


# %% 显示所有图窗
utils.show()

---

7.3 DSB-SC调制解调仿真程序 DSB-SC.py

# -*- coding: utf-8 -*-
# DSB-SC.PY
# %% 配置导入
from scipy import signal
from scipy import fft
import numpy as np
import utils


# %% 参数设置
time_begin = 0  # 起始时间
time_end = 1  # 终止时间
points = 1024  # 栅网采样点数(采样频率)
f_sample = points / (time_end - time_begin)  # 采样频率
T = (time_end - time_begin) / (points - 1)  # 栅网采样时间间距
frequency = [10, 100]  # 频率
amplitude = [1, 1]  # 振幅
phase = [0, 0]  # 相位(cos(t)相位为0)
zoom = 4  # 频谱缩放倍数
display_interval = 1.25  # 显示区间的系数
tick_max = max(amplitude) * 1.5  # 波形纵轴刻度最大值
tick_min = -tick_max  # 波形纵轴刻度最小值

snr = 6  # AWGN生成时的信噪比, 单位(dB)
seed = 3  # AWGN发生器随机种子

filter_order = 8  # 滤波器的阶数

f_H_LPF = frequency[0] * 2  # 低通滤波器的上限截止频率
f_H_BPF_DSB = frequency[1] + frequency[0] * 2  # 带通滤波器的上限截止频率
f_L_BPF_DSB = frequency[1] - frequency[0] * 2  # 带通滤波器的下限截止频率


# %% 生成时域波形
t = np.linspace(time_begin, time_end, points)  # 生成时域栅网
waves = []
for i in range(len(frequency)):
    waves.append(
        signal.chirp(
            t,
            f0=frequency[i],
            t1=time_end,
            f1=frequency[i],
            phi=phase[i],
            method='linear'
        ) * amplitude[i]
    )
s = None
for wave in waves:
    if s is None:
        s = wave.copy()
    else:
        s *= wave
s_d = waves[0] * waves[1] * waves[1]


# %% 获得频域波形
f = np.linspace(0, 1 / T, points)  # 生成频域栅网
# 输入信号 M(f)
M = fft.fft(waves[0])  # 快速傅里叶变换
M_nor = M / points  # 归一化
M_mod = np.abs(M_nor)  # 获取幅度频谱
M_pha = np.angle(M)  # 获取相位频谱
# 移频信号 C(f)
C = fft.fft(waves[1])  # 快速傅里叶变换
C_nor = C / points  # 归一化
C_mod = np.abs(C_nor)  # 获取幅度频谱
C_pha = np.angle(C)  # 获取相位频谱
# 第一次移频 -> S(f)
S = fft.fft(s)  # 快速傅里叶变换
S_nor = S / points  # 归一化
S_mod = np.abs(S_nor)  # 获取幅度频谱
S_pha = np.angle(S)  # 获取相位频谱
# 第二次移频 -> S_d(f)
S_d = fft.fft(s_d)
S_d_nor = S_d / points  # 归一化
S_d_mod = np.abs(S_d_nor)  # 获取幅度频谱
S_d_pha = np.angle(S_d)  # 获取相位频谱


# %% 滤波器参数
"""
    滤波器构造函数(仅介绍Butterworth滤波器)
        scipy.signal.butter(N, Wn, btype='low', analog=False, output='ba')
    输入参数：
        N:滤波器的阶数
        Wn：归一化截止频率。计算公式Wn=2*截止频率/采样频率。（注意：根据采样定理，采样频率要大于两倍的信号本身最大的频率，才能还原信号。截止频率一定小于信号本身最大的频率，所以Wn一定在0和1之间）。当构造带通滤波器或者带阻滤波器时，Wn为长度为2的列表。
        btype : 滤波器类型{‘lowpass’, ‘highpass’, ‘bandpass’, ‘bandstop’},
        output : 输出类型{‘ba’, ‘zpk’, ‘sos’},
    输出参数：
        b，a: IIR滤波器的分子（b）和分母（a）多项式系数向量。output='ba'
        z,p,k: IIR滤波器传递函数的零点、极点和系统增益. output= 'zpk'
        sos: IIR滤波器的二阶截面表示。output= 'sos'
    函数的使用
        信号滤波中最常用的无非低通滤波、高通滤波和带通滤波。下面简单介绍这三种滤波的使用过程：
    (1).高通滤波
        # 这里假设采样频率为1000hz,信号本身最大的频率为500hz，要滤除10hz以下频率成分，即截至频率为10hz，则wn=2*10/1000=0.02
        # from scipy import signal
        # b, a = signal.butter(8, 0.02, 'highpass')
        # filtedData = signal.filtfilt(b, a, data)#data为要过滤的信号
    (2).低通滤波
        # 这里假设采样频率为1000hz,信号本身最大的频率为500hz，要滤除10hz以上频率成分，即截至频率为10hz，则wn=2*10/1000=0.02
        # from scipy import signal
        # b, a = signal.butter(8, 0.02, 'lowpass')
        # filtedData = signal.filtfilt(b, a, data)       #data为要过滤的信号
    (3).带通滤波
        # 这里假设采样频率为1000hz,信号本身最大的频率为500hz，要滤除10hz以下和400hz以上频率成分，即截至频率为10hz和400hz,则wn1=2*10/1000=0.02,wn2=2*400/1000=0.8。Wn=[0.02,0.8]
        # from scipy import signal
        # b, a = signal.butter(8, [0.02,0.8], 'bandpass')
        # filtedData = signal.filtfilt(b, a, data)   #data为要过滤的信号
"""
wn_LPF = 2 * f_H_LPF / f_sample
b_LPF, a_LPF = signal.butter(filter_order, wn_LPF, 'lowpass')  # 配置滤波器

wn_BPF_DSB = [2 * f_L_BPF_DSB / f_sample, 2 * f_H_BPF_DSB / f_sample]
b_BPF_DSB, a_BPF_DSB = signal.butter(filter_order, wn_BPF_DSB, 'bandpass')  # 配置滤波器


# %% 进行滤波, 获得 m_o(t) 的时域波形
"""
    滤波函数
        scipy.signal.filtfilt(b, a, x, axis=-1, padtype='odd', padlen=None, method='pad', irlen=None)
    输入参数：
        b: 滤波器的分子系数向量
        a: 滤波器的分母系数向量
        x: 要过滤的数据数组。（array型）
        axis: 指定要过滤的数据数组x的轴
        padtype: 必须是“奇数”、“偶数”、“常数”或“无”。这决定了用于过滤器应用的填充信号的扩展类型。{‘odd’, ‘even’, ‘constant’, None}
        padlen：在应用滤波器之前在轴两端延伸X的元素数目。此值必须小于要滤波元素个数 - 1。（int型或None）
        method：确定处理信号边缘的方法。当method为“pad”时，填充信号；填充类型padtype和padlen决定，irlen被忽略。当method为“gust”时，使用古斯塔夫森方法，而忽略padtype和padlen。{“pad” ，“gust”}
        irlen：当method为“gust”时，irlen指定滤波器的脉冲响应的长度。如果irlen是None，则脉冲响应的任何部分都被忽略。对于长信号，指定irlen可以显著改善滤波器的性能。（int型或None）
    输出参数：
        y: 滤波后的数据数组
"""
m_o = signal.filtfilt(b_LPF, a_LPF, s_d)  # 为要过滤的信号
t_m_original = np.linspace(time_begin, T * (len(m_o) - 1), len(m_o))  # m_o(t) 信号的栅网


# %% 获得 m_o(t) 的频域波形
M_o = fft.fft(m_o)
M_o_nor = M_o / points  # 归一化
M_o_mod = np.abs(M_o_nor)  # 获取幅度频谱
M_o_pha = np.angle(M_o)  # 获取相位频谱


# %% AWGN生成与处理
n_i = utils.wgn(s, snr, seed)  # 生成AWGN
n = signal.filtfilt(b_BPF_DSB, a_BPF_DSB, n_i)  # 获得带通AWGN
n_d = n * waves[1]  # 获得移频后的带通AWGN
n_o = signal.filtfilt(b_LPF, a_LPF, n_d)  # 获得通过低通滤波器的带通AWGN


# %% 绘制 m(t) c(t) s(t) 时域波形
utils.drawWave(
    figure_i=1,
    row_i=3,
    column_j=1,
    label_x=r"t / s",
    axis_xs=[t] * 3,
    axis_ys=waves + [s],
    formats=['b-,', 'g-,', 'r-,'],
    labels=[r"m(t) = cos(20πt)",
            r"c(t) = cos(200πt)",
            r"s(t) = m(t)c(t)"],
    label_locs=[r"upper right"] * 3,
    limit_ys=[(tick_min, tick_max)] * 3,
    linewidths=[1] * 3
)


# %% 绘制 m(t) 频域波形
utils.drawWave(
    figure_i=2,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 + [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((M_mod[-1:0:-1], M_mod)),
             np.hstack((-M_pha[-1:0:-1], M_pha)),
             M_mod[0: points // zoom],
             M_pha[0: points // zoom]],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|M(f)|", r"θᴍ(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(M_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制 c(t) 频域波形
utils.drawWave(
    figure_i=3,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 + [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((C_mod[-1:0:-1], C_mod)),
             np.hstack((-C_pha[-1:0:-1], C_pha)),
             C_mod[0: points // zoom],
             C_pha[0: points // zoom]],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|C(f)|", r"θᴄ(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(C_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制 s(t) 频域波形
utils.drawWave(
    figure_i=4,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 + [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((S_mod[-1:0:-1], S_mod)),
             np.hstack((-S_pha[-1:0:-1], S_pha)),
             S_mod[0: points // zoom],
             S_pha[0: points // zoom]],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|S(f)|", r"θs(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(S_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制 s(t) c(t) s_d(t) 时域波形
utils.drawWave(
    figure_i=5,
    row_i=3,
    column_j=1,
    label_x=r"t / s",
    axis_xs=[t] * 3,
    axis_ys=[s, waves[1], s_d],
    formats=['b-,', 'g-,', 'r-,'],
    labels=[r"s(t) = m(t)c(t)",
            r"c(t) = cos(200πt)",
            r"s_d(t) = s(t)c(t)"],
    label_locs=[r"upper right"] * 3,
    limit_ys=[(tick_min, tick_max)] * 3,
    linewidths=[1] * 3
)


# %% 绘制 s_d(t) 频域波形
utils.drawWave(
    figure_i=6,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 + [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((S_d_mod[-1:0:-1], S_d_mod)),
             np.hstack((-S_d_pha[-1:0:-1], S_d_pha)),
             S_d_mod[0: points // zoom],
             S_d_pha[0: points // zoom]],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|Sᴅ(f)|", r"θsᴅ(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(S_d_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制 m_o(t) 时域波形
utils.drawWave(
    figure_i=9,
    row_i=2,
    column_j=1,
    label_x=r"t / s",
    axis_xs=[t_m_original, t],
    axis_ys=[m_o, waves[0]],
    formats=['b-,', 'g-,'],
    labels=[r"m_o(t) = s_d(t) * h(t)", r"m(t)"],
    label_locs=[r"upper right"] * 2,
    limit_ys=[(tick_min, tick_max)] * 2,
    linewidths=[1] * 2
)


# %% 绘制 m_o(t) 频域波形
utils.drawWave(
    figure_i=10,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 +
    [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((M_o_mod[-1:0:-1], M_o_mod)),
             np.hstack((-M_o_pha[-1:0:-1], M_o_pha)),
             M_o_mod[0: points // zoom],
             M_o_pha[0: points // zoom]],
    # axis_xs=[f_h] * 2,
    # axis_ys=[H_mod,
    #          H_pha],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|Mo(f)|", r"θᴍo(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(M_o_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制噪声 n_i(t), n(t), n(t), n_o(t) 时域波形
utils.drawWave(
    figure_i=11,
    row_i=4,
    column_j=1,
    label_x=r"t / s",
    axis_xs=[t] * 4,
    axis_ys=[n_i, n, n_d, n_o],
    formats=['b-,', 'g-,', 'r-,', 'k-,'],
    labels=[r"n_i(t)",
            r"n(t)",
            r"n_d(t)",
            r"n_o(t)"],
    label_locs=[r"upper right"] * 4,
    limit_ys=[None] * 4,
    linewidths=[1] * 4
)


# %% 绘制噪声 n_i(t), n(t), n(t), n_o(t) 功率谱密度波形
utils.drawPSD(
    figure_i=12,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_ys=[n_i, n, n_d, n_o],
    F_samples=[f_sample] * 4,
    limit_ys=[None] * 4)


# %% 求功率, 信噪比与信噪比增益
S_in = np.mean(s ** 2)
N_in = np.mean(n ** 2)
SNR_in = S_in / N_in

S_out = np.mean(m_o ** 2)
N_out = np.mean(n_o ** 2)
SNR_out = S_out / N_out

G = SNR_out / SNR_in

print("S_in = %e\nN_in = %e\nSNR_in = %e\nS_out = %e\nN_out = %e\nSNR_out = %e\nG = %e\n" % (S_in, N_in, SNR_in, S_out, N_out, SNR_out, G))


# %% 显示所有图窗
utils.show()

---

7.4 SSB调制解调仿真程序 SSB.py

# -*- coding: utf-8 -*-
# AM.PY
# %% 配置导入
from scipy import fftpack
from scipy import signal
from scipy import fft
import numpy as np
import utils


# %% 参数设置
time_begin = 0  # 起始时间
time_end = 1  # 终止时间
points = 1024  # 栅网采样点数(采样频率)
f_sample = points / (time_end - time_begin)  # 采样频率
T = (time_end - time_begin) / (points - 1)  # 栅网采样时间间距
frequency = [10, 100, 100]  # 频率
amplitude = [1, 1, 1]  # 振幅
phase = [0, 0, -90]  # 相位(cos(t)相位为0)
zoom = 4  # 频谱缩放倍数
display_interval = 1.25  # 显示区间的系数
tick_max = max(amplitude) * 1.5  # 波形纵轴刻度最大值
tick_min = -tick_max  # 波形纵轴刻度最小值

snr = 6  # AWGN生成时的信噪比, 单位(dB)
seed = 3  # AWGN发生器随机种子

USSB_FLAG = True  # 是否为上边带调制

filter_order = 8  # 滤波器的阶数
f_H_LPF = frequency[0] * 2  # 低通滤波器的上限截止频率

if USSB_FLAG:
    # 上边带 USSB
    f_H_BPF_SSB = frequency[1] + frequency[0] * 2  # 带通滤波器的上限截止频率
    f_L_BPF_SSB = frequency[1]  # 带通滤波器的下限截止频率
    sign = -1  # 正负符号
else:
    # 下边带 LSSB
    f_H_BPF_SSB = frequency[1]  # 带通滤波器的上限截止频率
    f_L_BPF_SSB = frequency[1] - frequency[0] * 2  # 带通滤波器的下限截止频率
    sign = 1  # 正负符号



# %% 生成时域波形
t = np.linspace(time_begin, time_end, points)  # 生成时域栅网
waves = []
for i in range(len(frequency)):
    waves.append(
        signal.chirp(
            t,
            f0=frequency[i],
            t1=time_end,
            f1=frequency[i],
            phi=phase[i],
            method='linear'
        ) * amplitude[i]
    )
s = waves[0] * waves[1] - sign * fftpack.hilbert(waves[0]) * waves[2]
# 该函数希尔伯特变换将相位向前移动90度, cos(t) -> -sin(t), 因此需要取负值
s_d = s * waves[1]


# %% 获得频域波形
f = np.linspace(0, 1 / T, points)  # 生成频域栅网
# 输入信号 M(f)
M = fft.fft(waves[0])  # 快速傅里叶变换
M_nor = M / points  # 归一化
M_mod = np.abs(M_nor)  # 获取幅度频谱
M_pha = np.angle(M)  # 获取相位频谱
# 移频信号 C(f)
C = fft.fft(waves[1])  # 快速傅里叶变换
C_nor = C / points  # 归一化
C_mod = np.abs(C_nor)  # 获取幅度频谱
C_pha = np.angle(C)  # 获取相位频谱
# 第一次移频 -> S(f)
S = fft.fft(s)  # 快速傅里叶变换
S_nor = S / points  # 归一化
S_mod = np.abs(S_nor)  # 获取幅度频谱
S_pha = np.angle(S)  # 获取相位频谱
# 第二次移频 -> S_d(f)
S_d = fft.fft(s_d)
S_d_nor = S_d / points  # 归一化
S_d_mod = np.abs(S_d_nor)  # 获取幅度频谱
S_d_pha = np.angle(S_d)  # 获取相位频谱


# %% 滤波器参数
"""
    滤波器构造函数(仅介绍Butterworth滤波器)
        scipy.signal.butter(N, Wn, btype='low', analog=False, output='ba')
    输入参数：
        N:滤波器的阶数
        Wn：归一化截止频率。计算公式Wn=2*截止频率/采样频率。（注意：根据采样定理，采样频率要大于两倍的信号本身最大的频率，才能还原信号。截止频率一定小于信号本身最大的频率，所以Wn一定在0和1之间）。当构造带通滤波器或者带阻滤波器时，Wn为长度为2的列表。
        btype : 滤波器类型{‘lowpass’, ‘highpass’, ‘bandpass’, ‘bandstop’},
        output : 输出类型{‘ba’, ‘zpk’, ‘sos’},
    输出参数：
        b，a: IIR滤波器的分子（b）和分母（a）多项式系数向量。output='ba'
        z,p,k: IIR滤波器传递函数的零点、极点和系统增益. output= 'zpk'
        sos: IIR滤波器的二阶截面表示。output= 'sos'
    函数的使用
        信号滤波中最常用的无非低通滤波、高通滤波和带通滤波。下面简单介绍这三种滤波的使用过程：
    (1).高通滤波
        # 这里假设采样频率为1000hz,信号本身最大的频率为500hz，要滤除10hz以下频率成分，即截至频率为10hz，则wn=2*10/1000=0.02
        # from scipy import signal
        # b, a = signal.butter(8, 0.02, 'highpass')
        # filtedData = signal.filtfilt(b, a, data)#data为要过滤的信号
    (2).低通滤波
        # 这里假设采样频率为1000hz,信号本身最大的频率为500hz，要滤除10hz以上频率成分，即截至频率为10hz，则wn=2*10/1000=0.02
        # from scipy import signal
        # b, a = signal.butter(8, 0.02, 'lowpass')
        # filtedData = signal.filtfilt(b, a, data)       #data为要过滤的信号
    (3).带通滤波
        # 这里假设采样频率为1000hz,信号本身最大的频率为500hz，要滤除10hz以下和400hz以上频率成分，即截至频率为10hz和400hz,则wn1=2*10/1000=0.02,wn2=2*400/1000=0.8。Wn=[0.02,0.8]
        # from scipy import signal
        # b, a = signal.butter(8, [0.02,0.8], 'bandpass')
        # filtedData = signal.filtfilt(b, a, data)   #data为要过滤的信号
"""
wn_LPF = 2 * f_H_LPF / f_sample
b_LPF, a_LPF = signal.butter(filter_order, wn_LPF, 'lowpass')  # 配置滤波器

wn_BPF_DSB = [2 * f_L_BPF_SSB / f_sample, 2 * f_H_BPF_SSB / f_sample]
b_BPF_DSB, a_BPF_DSB = signal.butter(filter_order, wn_BPF_DSB, 'bandpass')  # 配置滤波器


# %% 进行滤波, 获得 m_o(t) 的时域波形
"""
    滤波函数
        scipy.signal.filtfilt(b, a, x, axis=-1, padtype='odd', padlen=None, method='pad', irlen=None)
    输入参数：
        b: 滤波器的分子系数向量
        a: 滤波器的分母系数向量
        x: 要过滤的数据数组。（array型）
        axis: 指定要过滤的数据数组x的轴
        padtype: 必须是“奇数”、“偶数”、“常数”或“无”。这决定了用于过滤器应用的填充信号的扩展类型。{‘odd’, ‘even’, ‘constant’, None}
        padlen：在应用滤波器之前在轴两端延伸X的元素数目。此值必须小于要滤波元素个数 - 1。（int型或None）
        method：确定处理信号边缘的方法。当method为“pad”时，填充信号；填充类型padtype和padlen决定，irlen被忽略。当method为“gust”时，使用古斯塔夫森方法，而忽略padtype和padlen。{“pad” ，“gust”}
        irlen：当method为“gust”时，irlen指定滤波器的脉冲响应的长度。如果irlen是None，则脉冲响应的任何部分都被忽略。对于长信号，指定irlen可以显著改善滤波器的性能。（int型或None）
    输出参数：
        y: 滤波后的数据数组
"""
m_o = signal.filtfilt(b_LPF, a_LPF, s_d)  # 为要过滤的信号
t_m_original = np.linspace(time_begin, T * (len(m_o) - 1), len(m_o))  # m_o(t) 信号的栅网


# %% 获得 m_o(t) 的频域波形
M_o = fft.fft(m_o)
M_o_nor = M_o / points  # 归一化
M_o_mod = np.abs(M_o_nor)  # 获取幅度频谱
M_o_pha = np.angle(M_o)  # 获取相位频谱


# %% AWGN生成与处理
n_i = utils.wgn(s, snr, seed)  # 生成AWGN
n = signal.filtfilt(b_BPF_DSB, a_BPF_DSB, n_i)  # 获得带通AWGN
n_d = n * waves[1]  # 获得移频后的带通AWGN
n_o = signal.filtfilt(b_LPF, a_LPF, n_d)  # 获得通过低通滤波器的带通AWGN


# %% 绘制 m(t) c(t) s(t) 时域波形
utils.drawWave(
    figure_i=1,
    row_i=3,
    column_j=1,
    label_x=r"t / s",
    axis_xs=[t] * 3,
    axis_ys=waves[0:2] + [s],
    formats=['b-,', 'g-,', 'r-,'],
    labels=[r"m(t) = cos(20πt)",
            r"c(t) = cos(200πt)",
            r"s(t)"],
    label_locs=[r"upper right"] * 3,
    limit_ys=[(tick_min, tick_max)] * 3,
    linewidths=[1] * 3
)


# %% 绘制 m(t) 频域波形
utils.drawWave(
    figure_i=2,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 + [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((M_mod[-1:0:-1], M_mod)),
             np.hstack((-M_pha[-1:0:-1], M_pha)),
             M_mod[0: points // zoom],
             M_pha[0: points // zoom]],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|M(f)|", r"θᴍ(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(M_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制 c(t) 频域波形
utils.drawWave(
    figure_i=3,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 + [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((C_mod[-1:0:-1], C_mod)),
             np.hstack((-C_pha[-1:0:-1], C_pha)),
             C_mod[0: points // zoom],
             C_pha[0: points // zoom]],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|C(f)|", r"θᴄ(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(C_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制 s(t) 频域波形
utils.drawWave(
    figure_i=4,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 + [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((S_mod[-1:0:-1], S_mod)),
             np.hstack((-S_pha[-1:0:-1], S_pha)),
             S_mod[0: points // zoom],
             S_pha[0: points // zoom]],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|S(f)|", r"θs(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(S_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制 s(t) c(t) s_d(t) 时域波形
utils.drawWave(
    figure_i=5,
    row_i=3,
    column_j=1,
    label_x=r"t / s",
    axis_xs=[t] * 3,
    axis_ys=[s, waves[1], s_d],
    formats=['b-,', 'g-,', 'r-,'],
    labels=[r"s(t) = m(t)c(t)",
            r"c(t) = cos(200πt)",
            r"s_d(t) = s(t)c(t)"],
    label_locs=[r"upper right"] * 3,
    limit_ys=[(tick_min, tick_max)] * 3,
    linewidths=[1] * 3
)


# %% 绘制 s_d(t) 频域波形
utils.drawWave(
    figure_i=6,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 + [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((S_d_mod[-1:0:-1], S_d_mod)),
             np.hstack((-S_d_pha[-1:0:-1], S_d_pha)),
             S_d_mod[0: points // zoom],
             S_d_pha[0: points // zoom]],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|Sᴅ(f)|", r"θsᴅ(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(S_d_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制 m_o(t) 时域波形
utils.drawWave(
    figure_i=9,
    row_i=2,
    column_j=1,
    label_x=r"t / s",
    axis_xs=[t_m_original, t],
    axis_ys=[m_o, waves[0]],
    formats=['b-,', 'g-,'],
    labels=[r"m_o(t) = s_d(t) * h(t)", r"m(t)"],
    label_locs=[r"upper right"] * 2,
    limit_ys=[(tick_min, tick_max)] * 2,
    linewidths=[1] * 2
)


# %% 绘制 m_o(t) 频域波形
utils.drawWave(
    figure_i=10,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_xs=[np.hstack((-f[-1:0:-1], f))] * 2 +
    [f[0:points//zoom]] * 2,
    axis_ys=[np.hstack((M_o_mod[-1:0:-1], M_o_mod)),
             np.hstack((-M_o_pha[-1:0:-1], M_o_pha)),
             M_o_mod[0: points // zoom],
             M_o_pha[0: points // zoom]],
    # axis_xs=[f_h] * 2,
    # axis_ys=[H_mod,
    #          H_pha],
    formats=['b-,', 'g-,'] * 2,
    labels=[r"|Mo(f)|", r"θᴍo(f)"] * 2,
    label_locs=[r"upper right"] * 4,
    limit_ys=[(0, max(M_o_mod) * display_interval), None] * 2,
    linewidths=[1] * 4
)


# %% 绘制噪声 n_i(t), n(t), n(t), n_o(t) 时域波形
utils.drawWave(
    figure_i=11,
    row_i=4,
    column_j=1,
    label_x=r"t / s",
    axis_xs=[t] * 4,
    axis_ys=[n_i, n, n_d, n_o],
    formats=['b-,', 'g-,', 'r-,', 'k-,'],
    labels=[r"n_i(t)",
            r"n(t)",
            r"n_d(t)",
            r"n_o(t)"],
    label_locs=[r"upper right"] * 4,
    limit_ys=[None] * 4,
    linewidths=[1] * 4
)


# %% 绘制噪声 n_i(t), n(t), n(t), n_o(t) 功率谱密度波形
utils.drawPSD(
    figure_i=12,
    row_i=4,
    column_j=1,
    label_x=r"f / Hz",
    axis_ys=[n_i, n, n_d, n_o],
    F_samples=[f_sample * 2] * 4,
    limit_ys=[None] * 4)


# %% 求功率, 信噪比与信噪比增益
S_in = np.mean(s ** 2)
N_in = np.mean(n ** 2)
SNR_in = S_in / N_in

S_out = np.mean(m_o ** 2)
N_out = np.mean(n_o ** 2)
SNR_out = S_out / N_out

G = SNR_out / SNR_in

print("S_in = %e\nN_in = %e\nSNR_in = %e\nS_out = %e\nN_out = %e\nSNR_out = %e\nG = %e\n" % (S_in, N_in, SNR_in, S_out, N_out, SNR_out, G))


# %% 显示所有图窗
utils.show()

本文地址：https://blog.csdn.net/weixin_43157144/article/details/109616435