信号频域处理-＞恢复到时域

时域采样 vs 频域采样

摘自信号处理里面频域采样是怎么一回事呢？

频域采样仿真

仿真代码及仿真结果如下，

clc;clear all;close all;

GHZ = 1e9;
MHZ = 1e6;
KHZ = 1e3;
HZ  = 1;

FS  = 100*MHZ;
TS  = 1/FS;
N   = 1024;
t   = (0:N-1)*TS;
AMP = 100;
f   = 10*MHZ;
sig = AMP*cos(2*pi*f*t);

y1 = fft(sig);
y2 = y1(1:4:end);

s_y1 = ifft(y1);
s_y2 = ifft(y2,'symmetric');

figure;
subplot(211);
hold on;
plot(sig,'r');
plot(s_y1,'k');
hold off;
legend('原始时域信号','恢复时域信号');
title('原始信号FFT，再IFFT恢复');

subplot(212);
hold on;
plot(sig,'r');
plot(s_y2(1:end),'g');
plot(s_y2(2:end),'k');
hold off;
legend('原始时域信号','恢复时域信号1','恢复时域信号2');
title('原始信号FFT，频域4倍抽取，再IFFT恢复');

仿真过程很粗暴，直接对FFT结果进行抽取采样，再恢复到时域。仿真结果与上面提到的频域采样是对应的上的，即：对FFT结果直接抽取，并不等价于对时域信号进行抽取，而是在时域上对原信号进行周期延拓，IFFT恢复出的时域信号的数据率并没有改变。

Figure1的结果自不必说，Figure2的结果有点要注意，直接恢复出来的时域波形相较原始波形有一点延迟（绿色波形），将这个延迟偏移量去掉后就与原始波形完全重合了（黑色波形）。
（由于单次FFT中并不含有时间信息，而对频谱直接抽取，在不混叠的情况下，频谱信息并没有丢失，再逆回去的时域并不一定能和原始数据对应上）

---

clc;clear all;close all;

GHZ  = 1e9;
MHZ  = 1e6;
KHZ  = 1e3;
HZ   = 1;

FS   = 100*MHZ;
TS   = 1/FS;
N    = 1024;
t    = (0:N-1)*TS;
f1   = 5*MHZ;
f2   = 3*MHZ;
sig1 = 50*cos(2*pi*f1*t);
% sig2 = 80*cos(2*pi*f2*t);
sig2 = 80*square(2*pi*30*t,50);

sig_use1  = [sig1(1:512) sig2(1:512)];
sig_use2  = sig1 + sig2;

y1   = fft(sig_use1);
y2   = y1(1:2:end);
s_y2 = ifft(y2);
figure;
subplot(211);
hold on;
plot(sig_use1,'r');
plot(sig_use2,'b');
plot(s_y2(1:end),'k');
hold off;
legend('sig use1(orig)','sig use2','恢复时域信号');
title('原始信号FFT，频域抽取，再IFFT恢复');

y1   = fft(sig_use2);
y2   = y1(1:2:end);
s_y2 = ifft(y2);
subplot(212);
hold on;
plot(sig_use1,'r');
plot(sig_use2,'b');
plot(s_y2(1:end),'k');
hold off;
legend('sig use1','sig use2(orig)','恢复时域信号');
title('原始信号FFT，频域抽取，再IFFT恢复');

频域上进行信号抽取、信号插值

时域信号抽取后的频谱分析

clc;clear all;close all;

GHZ = 1e9;
MHZ = 1e6;
KHZ = 1e3;
HZ  = 1;

FS  = 100*MHZ;
TS  = 1/FS;
N   = 1024;
t   = (0:N-1)*TS;
AMP = 100;
f   = 10*MHZ;

dec_ratio = 4;
sig    = AMP*cos(2*pi*f*t);
sig_D  = sig(1:dec_ratio:end);

y1 = abs(fft(sig));
y2 = abs(fft(sig_D))*dec_ratio;
figure;
hold on;
plot(y1,'r');
plot(y2,'k');
hold off;
legend('原始信号频谱','4倍抽取后信号频谱');
title('未经过fftshift处理,抽取频谱进行幅度补偿');

y1 = abs(fftshift(fft(sig)));
y2 = abs(fftshift(fft(sig_D)))*dec_ratio;
zeros_pad = zeros(1,(N-N/dec_ratio)/2);
y2_pad = [zeros_pad y2 zeros_pad];
figure;
hold on;
plot(y1,'r');
plot(y2_pad,'k');
hold off;
legend('原始信号频谱','4倍抽取后信号频谱');
title('经过fftshift处理,抽取频谱进行幅度补偿、两边补0');

在时域进行抽取后，FFT变换后，信号幅度有变化，为了便于分析，画图时进行了补偿；另外，可以发现，不进行fftshift时，原始信号频谱与抽取信号的频谱，有相当一部分是重叠的；在进行fftshift后，对抽取信号频谱进行两端补零，就和原始信号频谱基本一致了。

频域实现抽取

% compute spectrum
z = fftshift(fft(y))/numel(y);

% half of total number of reduced samples
n = numel(z)*(m-1)/2/m; 

% drop n number of samples from beginning and end of spectrum 
z(1:n) = [];
z(end-n+1:end) = [];

% Compute inverse FFT by forcing conjugate symmetric
yp = ifft(ifftshift(z)*numel(z),'symmetric');

频域实现插值

% Calculate number of padding zeros
padlen = length(y)*(m - 1);

% Compute number of half points in FFT 
zlen = ceil((length(y)+1)/2);

% Compute FFT
z = fft(y);

% Construct a new spectrum (row vector) by centering zeros
zp = [z(1:zlen) zeros(1, padlen) z(zlen+1:end)];

if  ~mod(length(y), 2); % even number
    zp(zlen) = zp(zlen)/2; zp(zlen+padlen) = zp(zlen);
end

% Compute inverse FFT and scale by m
yp = real(ifft(zp))*m;