使用矩阵运算替代 for 循环实现信号的DTFT

1. 问题的引入

  最近的数字信号处理实验中有这样的一道题

对 $x_a(t)=e^{-1000|t|}$信号进行采样分析，讨论信号的抽样的最低频率，验证奈奎斯特采样定理。要求对信号 $x_a(t)$ 分析，确定采样频率，求出采样后的序列 $x(n)$ 的频谱 ，比较不同采样频率情况下（满足和不满足采样定理） FT频谱 $x(j\Omega)$ 和 DTFT频谱$x(e^{j\omega})$ 关系。

实验思路如下。首先对信号进行截断，因为计算机只能分析有限长的序列，所以需要截断，截断的依据是当幅值衰减到小于某一阈值（如$e^{-5}$）时截断；其次对截断后的连续非周期信号进行傅里叶变换（FT），这里通过 MATLAB 的 Fourier 命令可以直接得到；得到原信号的FT之后根据频谱判断截止频率（最高频率）大概在什么位置，从而确定采样信号的频率（根据香农采样定理，采样频率是最高频率的二倍）；最后对采样后的信号进行DTFT。因为在 MATLAB 中没有直接求取 DTFT 的函数，所以需要自己编写。而本文的主要灵感也来自于自己编写DTFT的过程。

  假设我们现在已经得到了经过抽样之后的离散非周期序列 $x(n)$，现在要对其进行 DTFT 变换。由于 DTFT 的定义式如下

$$X(e^{j\omega})=\Sigma_{n=-\infty}^{\infty}x(n)e^{-j\omega n}\tag{1}\label{1}$$
所以如果从定义式下手，需要对上述的 $x(n)e^{-j\omega}$ 在使得 $x(n)$ 有意义上的区间内求和，也就是说是有限项累加求和的形式，这种形式的计算往往是通过 for 循环实现的。但是无论是在MATLAB 还是在 python 的 Numpy 计算库上，矩阵相乘的计算速度又要高于 for 循环的计算速度，所以将 for 循环改写成矩阵相乘的形式也是一件很应该做的事情，在机器学习中这个常常叫做“向量化”。

2. 使用矩阵运算实现DTFT

  假设$x(n)$ 是一个5 点的序列，即

$$\begin{bmatrix} x(0) & x(1) & x(2) & x(3) & x(4) \end{bmatrix}$$
根据公式 (1) 可知， $x(n)$ 的 DTFT 为
$$X(e^{j\omega}) =x(0)e^{-j\omega 0} +x(1)e^{-j\omega 1}+x(2)e^{-j\omega 2}+x(3)e^{-j\omega 3}+x(4)e^{-j\omega 4}$$
上式虽然是 5 项相加，但是实际上一个关于 $\omega$ 的函数。这 5 项求和的过程可以使用 for 进行求和得到，当然也可以通过矩阵相乘的形式获得

X(ejω)=[x(0)x(1)x(2)