python频域三维图_深入浅出通信原理连载1-21（Python代码版）

深入浅出通信原理Python代码版

深入浅出通信原理是陈爱军的心血之作，于通信人家园连载，此处仅作python代码笔记训练所用

陈老师的连载从多项式乘法讲起，一步一步引出卷积、傅立叶级数展开、旋转向量、三维频谱、IQ调制、数字调制等一系列通信原理知识

连载1：从多项式乘法说起

\[(x+1)(x^2+2x+5)=x^3+3x^2+7x+5

\]

import sympy

x = sympy.Symbol('x')

sympy.expand((x+1)*(x*x+2*x+5))

这种计算方法总结起来就是：

反褶：一般多项式都是按x的降幂排列，这里将其中一个多项式的各项按x的升幂排列。

平移：将按x的升幂排列的多项式每次向右平移一个项。

相乘：垂直对齐的项分别相乘。

求和：相乘的各结果相加。

反褶、平移、相乘、求和－这就是通信原理中最常用的一个概念“卷积”的计算过程。

连载2：卷积的表达式

将多项式转换成矩阵形式

\[x+1 -> [1,1]\\

x^2+2x+5 -> [1,2,5]\\

x^3+3x^2+7x+5=[1,1]*[1,2,5]=[1,3,7,5]\\

其中卷积(*)的运算如下：\\

c(n)=a(n)*b(n)=\sum_{k=0}^n{a(k)b(n-k)},n=0~(n1+n2)\\

其中n1是a(n)的系数总个数，n2是b(n)的系数总个数

\]

连载3: Python计算卷积

import numpy as np

np.convolve([1,1],[1,2,5]) #上例,result:[1, 3, 7, 5]

# 杨辉三角，输入行数，输出对应行的值

def pascal_triangle(n):

if n == 0:

return 1

elif n == 1:

return [1,1]

else:

return np.convolve([1,1],pascal_triangle(n-1))

# test code

for i in range(0,7):

print(pascal_triangle(i))

# 杨辉三角，迭代器形式

def triangles():

nlist=[1]

while True:

yield nlist

nlist.append(0)

nlist = [nlist[i] + nlist[i-1] for i in range(len(nlist))]

# test code

tr = triangles()

for i in range(0,7):

print(next(tr))

连载4：将信号表示成多项式的形式

\[Fourier Transform:signal(x) = x_nx^n+\cdots+a_2x^2+a_1x+a_0\\

x^n=f(nw_0)\\

x = cosw_0t+jsinw_0t = e^{jw_0t} = f(w_0)\\

x^2 = (cosw_0t+jsinw_0t)^2 = cos^2w_0t-sin^2w_0t+j2sinw_0tcosw_0t=f(2w_0)

\]

连载5：欧拉公式证明

\[欧拉公式:cosw_0t+jsinw_0t=e^{jw_0t}

\]

import sympy

x = sympy.Symbol('x')

sympy.series(sympy.exp(sympy.I*x), x, 0, 10) #series perform Taylor expansion

sympy.expand(series(sin(x), x, 0, 10))

sympy.expand(series(cos(x), x, 0, 10))

利用泰勒展开，从输出结果易看出欧拉公式成立

理解复数

引自[连载23,利用欧拉公式理解复数)

连载6：利用卷积计算两个信号的乘积

若两个信号均表示成多项式形式，则可通过卷积计算信号相乘结果。

连载7： 信号进行傅里叶级数展开

周期T的函数f(t)展开成傅里叶级数(直流分量、余弦分量，正选分量多项式之和)：

\[f(t)=\frac{a_0}{2}+\sum_{k=1}^{\infty}a_kcoskw_0t+b_ksinkw_0t\\

w=\frac{2\pi}{T}\\

a_k=\frac{2}{T}\int_{-2/T}^{2/T}f(t)cosk\omega_0tdt(k=0,1,2,3\cdots)\\

b_k=\frac{2}{T}\int_{-2/T}^{2/T}f(t)sink\omega_0tdt(k=1,2,3\cdots)

\]

展开成复指数信号之和：

\[f(t)=\sum_{k=-\infty}^{+\infty}c_ke^{jkw_0t}\\

c_k=

\begin{cases}

\frac{a_0}{2} &,& k=0 \\

\frac{1}{2}(a_k-jb_k) & , & k=1,2,\cdots \\

\frac{1}{2}(a_{-k}+jb_{-k}) & , & k=-1,-2,\cdots

\end{cases}

\]

ck是复傅里叶系数

以角频率w为横轴，傅里叶系数ck为纵轴，画出频谱图

# Further improvement:https://blog.csdn.net/u012915522/article/details/91506120

import matplotlib.pyplot as plt

plt.subplot(3,1,1)

plt.stem([0,1,2,3],[6,5,1,0],use_line_collection=True)

plt.yticks([0,1,5,6])

# plt.axis('off')

plt.xticks([0,1,2,3],('0','w','2w','3w'))

plt.subplot(3,1,2)

plt.stem([0,1,2,3],[2,3,0,0],use_line_collection=True)

plt.xticks([0,1,2,3],('0','w','2w','3w'))

plt.yticks([0,1,2,3])

plt.subplot(3,1,3)

plt.stem([0,1,2,3],[12,28,17,3],use_line_collection=True)

plt.xticks([0,1,2,3],('0',r'$\omega$w','2w','3w'))

plt.yticks([3,17,12,28])

# plt.axis('off')

plt.tight_layout()

plt.show()

连载8：时域信号相乘相当于频域卷积

由连载6显然如此

时域

FFT

频域

f(t)

\(f(t)=\sum_{n=-\infty}^{+\infty}f[n]^{jkw_0t}\)

f[n]

g(t)

$$g(t)=\sum_{n=-\infty}^{+\infty}g[n]^{jkw_0t}$$

g[n]

y(t)=f(t)g(t)

\(y(t)=\sum_{n=-\infty}^{+\infty}y[n]^{jkw_0t}\)

y[n]=f[n]*g[n]

连载9：用余弦信号合成方波信号

plt.suptitle('逐步合成方波', fontproperties='stsong')

plt.subplot(3,1,1)

x = np.arange(0,20,0.1)

y1 = 0.5+0.637*np.cos(x)

plt.plot(x,y1)

plt.title(r'$y_1(x)=0.5+0.637*cos(x)$',fontsize=10)

plt.subplot(3,1,2)

y2 = 0.5+0.637*np.cos(x)-0.212*np.cos(3*x)

plt.plot(x,y2)

plt.title(r'$y_2(x)=0.5+0.637*cos(x)-0.212cos(3x)$',fontsize=10)

plt.subplot(3,1,3)

y3 = 0.5+0.637*np.cos(x)-0.212*np.cos(3*x)+0.127*np.cos(5*x)

plt.plot(x,y3)

plt.title(r'$y_3(x)=0.5+0.637*cos(x))-0.212cos(3x)+0.127cos(5x)$',fontsize=10)

plt.tight_layout()

plt.show()

连载12：正弦信号及余弦信号频谱图

绘制正弦及余弦信号的频谱图

# sin(w0t) & cos(w0t) spectrum

fig = plt.figure('Sine & Cosine Wave', (10,8))

ax = axisartist.Subplot(fig, 2,2,1)

# plt.subplot(4,1,1)

fig.add_axes(ax)

# ax.axis[:].set_visible(False) # 隐藏默认坐标轴

# ax.axis["x"] = ax.new_floating_axis(0, 0)

# ax.axis["y"] = ax.new_floating_axis(1, 0)# 新建坐标轴X-Y

# ax.axis["x"].set_axis_direction('top')

# ax.axis["y"].set_axis_direction('left')# 设置刻度标识方向

# ax.axis["x"].set_axisline_style("->", size = 20.0)

ax.axis["bottom"].set_axisline_style("-|>", size = 1.5)

ax.axis["left"].set_axisline_style("->", size = 1.5)

ax.axis["top"].set_visible(False)

ax.axis["right"].set_visible(False)

# ax.axis["y"].set_axisline_style("->", size = 20.0)#加上坐标轴箭头

ax.stem([-1,1],[0.5,0.5],use_line_collection=True)

ax.set_title('y = cos($\omega_0t$),$c_{-1}=1/2,c_1=1/2$',fontsize = 15, pad = 20)

ax.set_xticks([-1,0,1])

ax.set_xticklabels([r'$-\omega$','0',r'$\omega$'])

ax.set_xlabel('Angular frequency')

ax.set_yticks([0.5])

ax.set_ylabel('magnitude')

ax2 = axisartist.Subplot(fig, 2,2,3)

fig.add_axes(ax2)

ax2.axis["left"].set_axisline_style("->", size = 1.5)

ax2.axis["top"].set_visible(False)

ax2.axis["right"].set_visible(False)

ax2.stem([-1,1],[0,0],use_line_collection=True)

ax2.set_xticks([-1,0,1])

ax2.set_xticklabels([r'$-\omega$','0',r'$\omega$'])

ax2.set_xlabel('Angular frequency')

ax2.set_yticks([0])

ax2.set_ylabel('phase')

fig.subplots_adjust(hspace=10)

ax3 = axisartist.Subplot(fig, 2,2,2)

fig.add_axes(ax3)

ax3.axis["left"].set_axisline_style("->", size = 1.5)

ax3.axis["top"].set_visible(False)

ax3.axis["right"].set_visible(False)

ax3.stem([-1,1],[0.5,0.5],use_line_collection=True)

ax3.set_title('y = sin($\omega_0t$),$ c_{-1}=j/2,c_1=-j/2 $,fontsize = 15, pad = 20)

ax3.set_xticks([-1,0,1])

ax3.set_xticklabels([r'$-\omega$','0',r'$\omega$'])

ax3.set_xlabel('Angular frequency')

ax3.set_yticks([0.5])

ax3.set_ylabel('magnitude')

ax4 = axisartist.Subplot(fig, 2,2,4)

fig.add_axes(ax4)

ax4.axis["left"].set_axisline_style("->", size = 1.5)

ax4.axis["top"].set_visible(False)

ax4.axis["right"].set_visible(False)

ax4.stem([-1,1],[1.5708,-1.5708],use_line_collection=True)

ax4.set_xticks([-1,0,1])

ax4.set_xticklabels([r'$-\omega$','0',r'$\omega$'])

ax4.set_xlabel('Angular frequency')

ax4.set_yticks([np.pi/2,0,-np.pi/2])

ax4.set_yticklabels([r'$\pi/2$','0',r'$-\pi/$'])

ax4.set_ylabel('phase')

plt.tight_layout()

plt.show()

由该图可看出，实信号的频谱是共轭对称的，及\(c_k&c_{-k}\)是互为共轭：模相等，相位相反。

而且出现负频率

余正弦信号的三维频谱图

这个太难，以我的浅薄知识，用Python绘制不出来，故从陈老师处直接盗图，但是可以从正弦信号三维频谱看出即使复傅里叶系数是虚数，其原始信号仍有可能是实数

通过复傅立叶级数展开将实信号分解为了一系列的旋转向量之和(由此引出了复数，使得实信号的表达式中出现了复数)，但由于逆时针和顺时针旋转的向量成对出现，而且成对出现的旋转向量的初始相位关于实轴对称，旋转的角速度相同，旋转方向相反，所以这些旋转向量合成的结果最终还是一个实信号(只在实轴上有分量，虚轴上的分量相互抵消掉了)。

李萨育图形

使用互相成谐波频率关系的两个信号分别作为X和Y偏转信号送入示波器时，这两个信号分别在X轴、Y轴方向同时作用于电子束而描绘出稳定的图形，这些稳定的图形就叫“李萨育图形”

Lissajous curve, 用来展示两个相互垂直的简谐振动的合成。

公式：

\[\left\{

\begin{aligned}

x(\theta) & = a\cos(\theta) \\

y(\theta) & = b\sin(n\theta+\phi) \\

\end{aligned}

\right.

\]

def lissajous(q,fignum):

theta = np.arange(0, 2*np.pi, 0.01)

x = np.cos(theta)

y = np.sin(q*theta)

plt.subplot(1,fignum,q)

plt.axis('off')

plt.plot(x,y)

return

plt.figure(figsize=(10,2))

fignum=5

for i in range(fignum):

lissajous(i+1,fignum)

实信号和复信号的波形对比

直接采用陈老师的三维频谱图：

在下面两张图中：x轴(实轴)、y轴(虚轴)所在的平面是复平面，t轴(时间轴)垂直于复平面。

上图为实信号f(t)=cos(2πt)的波形图。

下图为复信号f(t)=cos(2πt)+jsin(2πt)的波形图。

对比这两张图，很容易得出：实信号在复平面上投影时只有实轴方向有分量，而复信号在复平面上投影时实轴和虚轴方向都有分量。

最后说一句，由于博客部分图片显示不出，特贴上个人备份GitHub地址