一、采样原理和方法

时域采样定理的要点是：

  1. 对模拟信号xa(t)x_a(t)以间隔T进行时域等间隔理想采样，形成的采样信号的频谱X(jΩ)X(j\Omega)是原模拟信号频谱Xa(jΩ)X_a(j\Omega)以采样角频率Ωs=2πT{\Omega}_s=\frac{2\pi}{T}为周期进行周期延拓。公式为：

X(jΩ)=FT[xa(t)]=1T∑n=−∞∞Xa(jΩ−jnΩs)X(j\Omega) =FT[x_a(t)]=\frac{1}{T}\sum_{n=-\infty}^{\infty}X_a(j\Omega-jn{\Omega}_s)

  1. 奈奎斯特采样定理：采样频率Ωs{\Omega}_s必须大于等于模拟信号最高频率的两倍以上，才能使采样信号的频谱不产生频率谱混叠。

上式说明理想采样信号的傅里叶变换可用相应的采样序列的傅里叶变换得到，只要将自变量ww用ΩT{\Omega}T替代即可。

结论：时域采样频谱周期延拓，频域采样时域周期延拓。

二、验证采样原理

2.1 时域采样定理的验证

给定模拟信号，xa(t)=Ae−αtsin(Ω0t)u(t)x_a(t)=Ae^{-{\alpha}t}sin({\Omega}_{0}t)u(t)，式中A=444.128,α=502π,Ω0=502πrad/sA=444.128,\alpha=50\sqrt{2}\pi,\Omega_0 = 50\sqrt{2}{\pi} rad/s，其幅频特性曲线为：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
A = 444.128;
a = 50*sqrt(2)*pi;
w0 = 50*sqrt(2)*pi;
F1 = 1000;
F2 = 300;
F3 = 200;

T1 = 1 / F1;
T2 = 1 / F2;
T3 = 1 / F3;
Tp = 64 / 1000;

n1 = 0: Tp*F1-1; 
n2 = 0: Tp*F2-1;
n3 = 0: Tp*F3-1;

N1 = Tp * F1;
N2 = Tp * F2;
N3 = Tp * F3;

y1 = A * exp(-a *n1 *T1).*sin(w0 *n1 *T1);
y2 = A * exp(-a *n2 *T2).*sin(w0 *n2 *T2);
y3 = A * exp(-a *n3 *T3).*sin(w0 *n3 *T3);

fft1 = fft(y1, length(n1));
fft2 = fft(y2, length(n2));
fft3 = fft(y3, length(n3));

% 画图

subplot(3,2,1);
stem(n1, y1, ".");
title("(a)Fs=1000Hz")
grid on;
subplot(3,2,2);
plot((0:length(fft1)-1)/Tp,abs(fft1)*T1);
title("(a)T*FT[xa(nT)],Fs=1000Hz");
xlabel("f(Hz)");
ylabel("幅度");

subplot(3,2,3);
stem(n2, y2,".");
title("(a)Fs=200Hz")
grid on;
subplot(3,2,4);
plot((0:length(fft2)-1)/Tp,abs(fft2)*T2);
title("(a)T*FT[xa(nT)],Fs=200Hz");
xlabel("f(Hz)");
ylabel("幅度");

subplot(3,2,5);
stem(n3, y3, ".");
title("(a)Fs=300Hz")
grid on;
subplot(3,2,6);
plot((0:N3-1)/Tp,abs(fft3)*T3);
title("(a)T*FT[xa(nT)],Fs=300Hz");
xlabel("f(Hz)");
ylabel("幅度");

练习：近似绘制X(n)=R4(n)X(n)=R_4(n)在(0,2π)(0,2\pi)上的幅频响应曲线(∣FT[x(n)]∣)(|FT[x(n)]|)：

1
2
3
4
5
6
7
8
9
10
11
x = [1 1 1 1];
N = 64;
xk = fft(x, N);
figure;
subplot(2, 1, 1);
stem(0:3, x, ".");
subplot(2,1,2);
k=0:N-1;
plot(2*k/N, abs(xk));
title("|FT[x(n)]|");
xlabel("\omega/\pi");

注意：在这段MATLAB代码中，plot函数中传入参数2*k/N的目的是为了绘制频谱图。在频谱图中，横轴表示频率，纵轴表示信号的幅度。频率的范围是从0到π，对应着k的范围从0到N-1。通过将k除以N，可以将频率范围映射到0到1之间。然后乘以2，可以将频率范围映射到0到2之间，这样就可以在横轴上正确显示频率的范围。因此，传入参数2*k/N的意义是对频率进行适当的缩放，以便正确显示频谱图。

三、频率采样理论的验证

给定信号如下：

x(n)={n+10≤n≤1327−n14≤n≤260其他x(n)=\begin{cases} n+1 & 0{\leq}n{\leq}13 \\ 27-n & 14 {\leq}n{\leq}26 \\ 0 & 其他 \end{cases} 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
%数字信号处理-实验5-频域域采样理论论证
clc;
clear;

%序列一共27长
M = 27;
N = 32;
n = 0:M-1;
xa = 1:14; xb = 13:-1:1; xn = [xa,xb];

%32点采样、还原
x32k = fft(xn,32);
x32n = ifft(x32k);
%16点采样、还原,隔一个取一个
x16k = x32k(1:2:N);
x16n = ifft(x16k,N/2);

subplot(3,2,2);
stem(n,xn,'.');
xlabel('n');ylabel('x(n)');title('(b)三角波序列x(n)');
axis([0,32,0,20]);
grid on;

%取2的整次幂
k = 0 : 1023;
wk = 2*k/1024;
X = fft(xn,1024);
subplot(3,2,1);
plot(wk,abs(X));
axis([0 1 0 200]);
xlabel('ω/π');ylabel('|X(e^j^ω)|');title('(a)FT[x(n)]');
grid on;

k = 0:N/2-1;
subplot(3,2,3);
stem(k,abs(x16k),'.');
xlabel('k'); ylabel('|X_1_6(k)|'); title('(c)16点频域采样');
axis([0,8,0,200]);
grid on;

subplot(3,2,4);
stem(k,x16n,'.');
xlabel('n'); ylabel('|X_1_6(n)|');title('(d)16点IDFT[X_1_6(k)]');
axis([0,32,0,20]);
grid on;

k = 0:N-1;
subplot(3,2,5);
stem(k,abs(x32k),'.');
xlabel('k'); ylabel('|X_3_2(k)|');title('(e)32点频域采样');
axis([0,16,0,200]);
grid on;
subplot(3,2,6);
stem(k,x32n,'.');
xlabel('n'); ylabel('|X_3_2(k)|');title('(f)32点IDFT[X_3_2(k)]');
axis([0,32,0,20]);
grid on;