复正弦信号、实数信号的DFT

本文将介绍复正弦信号和实数信号的DFT性质，内容包括：

复正弦信号的定义以及DFT性质
实数信号的定义以及性质
代码实例

复正弦信号

首先说明，复正弦信号只是一种数学定义上的信号，在实际生活中没有发现过这种信号

我们先来看一种较为简单且比较有规律的复正弦信号，它的定义如下：

x1[n]=ej2πk0n/Nn=0,1,2,...,N-1x_1[n] = e^{j2\pi k_0 n / N} \\ \text{n=0,1,2,...,N-1} x1[n]=ej2πk0n/Nn=0,1,2,...,N-1

其中，eee是自然常熟，N表示信号的个数。

至于k0k_0k0，它是很关键的一个参数，我们可以将它理解为x1[n]x_1[n]x1[n]信号的周期个数。通常我们只能对一小段信号做DFT，因此需要先对信号截断，如果截断的信号刚好是整数周期，即k0k_0k0是一个整数，那么信号不发生泄露；相反，如果截断信号不是整数周期，即k0k_0k0是一个浮点数，那么发现信号泄露。

关于信号泄露，请参看这篇文章什么是泄露

复正弦信号的DFT

将复正弦信号的定义直接套上DFT的公式可以推出其DFT形式：
X1[k]=∑n=0N−1x1[n]e−j2πkn/N=∑n=0N−1ej2πk0n/Ne−j2πkn/N=∑n=0N−1ej2π(k0−k)n/N(1)\begin{array}{l} X_1[k] &= \sum_{n=0}^{N-1} x_1[n] e^{-j2\pi k n / N} \\ &= \sum_{n=0}^{N-1} e^{j2\pi k_0 n / N}e^{-j2\pi k n / N} \\ &= \sum_{n=0}^{N-1} e^{j2\pi(k_0 - k)n / N}\\ \end{array}\tag{1} X1[k]=∑n=0N−1x1[n]e−j2πkn/N=∑n=0N−1ej2πk0n/Ne−j2πkn/N=∑n=0N−1ej2π(k0−k)n/N(1)

当 k0=kk_0 = kk0=k 时，ej2π(k0−k)n/N=e0=1e^{j2\pi(k_0 - k)n / N} = e^{0} = 1ej2π(k0−k)n/N=e0=1，因此
X1[k]=∑n=0N−11=N,if k0=kX_1[k] = \sum_{n=0}^{N-1} 1 = N, \text{if } k_0=kX1[k]=n=0∑N−11=N,if k0=k

当 k0≠kk_0 \neq kk0=k 时，利用等比数列求和公式可得
X1[k]=1−ej2π(k0−k)n1−ej2π(k0−k)n/N=0\begin{array}{l} X_1[k] &= \frac{1-e^{j2\pi(k_0 - k)n}}{1-e^{j2\pi(k_0 - k)n/N}} \\ &= 0 \end{array} X1[k]=1−ej2π(k0−k)n/N1−ej2π(k0−k)n=0

即
X1[k]={N,if k0=k0,if k0≠kX_1[k] = \begin{cases} N, & \text{if $k_0 = k$ } \\ 0, & \text{if $k_0 \neq k$} \end{cases} X1[k]={N,0,if k0=k if k0=k

这里仔细说明下X1[0]=0X_1[0] = 0X1[0]=0是怎么来的

k0k_0k0和kkk都是整数，那他们的差k0−kk_0-kk0−k也是整数
当N是整数时，ej2πN=1e^{j2\pi N} = 1ej2πN=1，因此1−ej2π(k0−k)n=01-e^{j2\pi(k_0 - k)n} = 01−ej2π(k0−k)n=0
根据1,2，推出X1[0]=0X_1[0] = 0X1[0]=0

等比数列求和的推导，可以参考这个公式:
∑k=0Nzk={1−zN+11−z,if z≠1N,if z=1\sum_{k=0}^{N} z^k = \begin{cases} \frac{1-z^{N+1}}{1-z}, & \text{if $z\neq 1$ } \\ N, & \text{if $z = 1$} \end{cases} k=0∑Nzk={1−z1−zN+1,N,if z=1 if z=1

复正弦信号DFT实例

我们将通过代码来展示复正弦信号的性质，首先生成一段复正弦信号，然后进行DFT，画出结果

import 需要的包

import numpy as np
import matplotlib.pyplot as pltfrom scipy.fftpack import fft%matplotlib inline

定义生成复正弦信号的函数

def generate_complex_signal(num_sample, k0):'''generate a complex signalnum_sample : 信号的个数，即公式中的Nk0 : 周期个数returnsx : 复正弦信号'''n = np.arange(num_sample)x = np.exp(1j*2*np.pi*k0*n/num_sample)return x

生成信号，进行DFT，画出结果

当k0是整数时，DFT后的结果与预期的一致

# DFT and plot the results
num_sample = 100
k0 = 20
x = generate_complex_signal(num_sample, k0)X = fft(x)plt.figure(figsize=(10,6))
plt.subplot(211)
plt.plot(np.real(x))plt.subplot(212)
plt.plot(np.real(X))
plt.show()

实数信号

我们生活中遇到的信号都是实数信号，它的定义如下：
x2[n]=A0cos⁡(2πk0n/N)=A02(ej2πk0n/N+e−j2πk0n/N)x_2[n] = A_0 \cos (2\pi k_0 n/N) = \frac{A_0}{2}(e^{j2\pi k_0 n/N} + e^{-j2\pi k_0 n/N}) x2[n]=A0cos(2πk0n/N)=2A0(ej2πk0n/N+e−j2πk0n/N)

额…虽然生活中的信号都是实数信号，但是上面的式子是对实数信号的超高度浓缩版本，是最简单的版本。实际中的信号模型更为复杂。

实数信号的DFT

将实数信号带入DFT等式中

X2[k]=∑n=0N−1x2[n]e−j2πkn/N=∑n=0N−1A02(ej2πk0n/N+e−j2πk0n/N)e−j2πkn/N=∑n=0N−1A02e−j2π(k−k0)n/N+∑n=0N−1A02e−j2π(k+k0)n/N(2)\begin{array}{l} X_2[k] &= \sum_{n=0}^{N-1} x_2[n] e^{-j2\pi k n / N} \\ &= \sum_{n=0}^{N-1} \frac{A_0}{2}(e^{j2\pi k_0 n/N} + e^{-j2\pi k_0 n/N}) e^{-j2\pi k n / N} \\ &= \sum_{n=0}^{N-1} \frac{A_0}{2}e^{-j2\pi(k-k_0) n / N} + \sum_{n=0}^{N-1} \frac{A_0}{2} e^{-j2\pi(k+k_0) n / N} \\ \end{array}\tag{2} X2[k]=∑n=0N−1x2[n]e−j2πkn/N=∑n=0N−12A0(ej2πk0n/N+e−j2πk0n/N)e−j2πkn/N=∑n=0N−12A0e−j2π(k−k0)n/N+∑n=0N−12A0e−j2π(k+k0)n/N(2)

当k=k0k=k_0k=k0或者k=−k0k=-k_0k=−k0时，(2)式中一项为0，另一项为NA02N\frac{A_0}{2}N2A0

否则(2)式为0

即
X2[k]={NA02,if k=k0,−k00,if k≠k0X_2[k] = \begin{cases} N\frac{A_0}{2}, & \text{if $k = k_0, -k_0$ } \\ 0, & \text{if $k \neq k_0$} \end{cases} X2[k]={N2A0,0,if k=k0,−k0 if k=k0

实数信号DFT实例

我们将通过代码来展示实数信号的性质，首先生成一段复正弦信号，然后进行DFT，画出结果

def generate_real_signal(num_sample, A, k0):'''generate real signalnum_sample : 信号的个数，即公式中的NA  : 振幅k0 : 周期个数returnsx : 实数信号'''hN = num_sample//2n = np.arange(-hN, hN)x = A * np.cos( 2*np.pi*k0*n/num_sample )return x

通过打印出结果，我们发现确实有2个峰值，但是并不是k0k_0k0和−k0-k_0−k0

其实只要做一次交换就可以，具体请看代码Y的实现

# DFT and plot the results
num_sample = 100
k0 = 20
A = 0.8
x = generate_real_signal(num_sample, A, k0)X = fft(fftbuffer)plt.figure(figsize=)
plt.subplot(311)
plt.plot(np.real(x))plt.subplot(312)
plt.plot(np.real(X))hM1 = (num_sample+1)//2
hM2 = num_sample//2
Y = np.zeros(num_sample)
Y[:hM1] = X[-hM1:]
Y[-hM2:] = X[:hM1]
plt.subplot(313)
x_axis = np.arange(-hM1, hM2)
plt.plot(x_axis, np.real(Y))
plt.show()