信息处理 - 矩形脉冲傅里叶变换与FFT的区别 - 吾爱随笔录

矩形脉冲傅里叶变换与FFT的区别

信息处理 matlab fft 傅里叶变换频谱脉冲

2022-02-19 09:02:07

我试图找到非周期信号的傅里叶变换与它们的 FFT 之间的联系。因此，从一个基本示例开始，让我们采用一个宽度为 0.1s、幅度为 1 的矩形脉冲偏移 0.05。使用对应关系，我可以计算预期的频谱： $X(f) = 0.1 \cdot sinc(0.1f) \cdot e^{j 2 \pi f \cdot 0.05}$

但是现在，当我使用以下 Matlab 代码生成信号时：

f_abt = 50e3;
x=0:1/f_abt:1;
y=zeros(1,length(x));
for ii=1:length(x)
    if x(ii)<=.1
        y(ii)=1;
    end
end

并计算它的频谱，结果取决于信号的长度。因此，当我从上面生成的信号（1s 持续时间）计算单边频谱时，我得到：

然后，当我将信号长度设置为 2s 时（其他一切不变）：

x=0:1/f_abt:2;

我得到以下光谱：

我想差异来自我使用的 FFT 算法。在进行 FFT 时，我通过 Nfft 对值进行归一化，因此当我改变信号长度时，我的幅度发生变化是完全有道理的。

我的问题是：我如何获得正确的频谱以及我如何知道它是正确的，例如，当我无法使用通信“手动”计算它时？我在寻找我的“真实”限时信号及其 FFT 和“理论”矩形脉冲之间的联系时遇到问题。

我用于计算单边光谱的代码：

function [f_xa, mag, phase] = calc_fft_f(ta, xa)
N_a = numel(xa);
fft_xa = fft(xa); 
P2_norm = fft_xa/(N_a);
if (mod(N_a,2))
    P1_norm_single = P2_norm(1:ceil(end/2));    
    P1_norm_single(2:end) = 2*P1_norm_single(2:end);
else 
    P1_norm_single = P2_norm(1:(end/2)+1);
    P1_norm_single(2:end-1) = 2*P1_norm_single(2:end-1);
end

mag = abs(P1_norm_single);
phase = rad2deg(angle(P1_norm_single));
Fsa = 1/(ta(2)-ta(1));
f_xa = Fsa*(0:(length(mag)-1))/N_a;

end

提前致谢！

1个回答

假设连续时间信号的相关部分 $x(t)$ 位于（或已移至）区间内 $[0,T]$ ，信号的采样版本的 DFT 以下列方式逼近连续时间傅里叶变换 (CTFT)：

\begin{aligned} X (f) & = \int_{- \infty}^{\infty} x (t) e^{- j 2 π f t} d t \\ \overset{truncation}{\approx} \int_{0}^{T} x (t) e^{- j 2 π f t} d t \\ (1) & \overset{sampling}{\approx} \sum_{n = 0}^{N - 1} x (n Δ t) e^{- j 2 π f n Δ t} Δ t \end{aligned}

$\begin{align}X(f)&=\int_{-\infty}^{\infty}x(t)e^{-j2\pi ft}dt\\&\stackrel{\textrm{truncation}}{\approx}\int_{0}^{T}x(t)e^{-j2\pi ft}dt\\&\stackrel{\textrm{sampling}}{\approx}\sum_{n=0}^{N-1}x(n\Delta t)e^{-j2\pi f n\Delta t}\Delta t\tag{1}\end{align}$

和 $T=N\Delta t$ . 从 $(1)$ 和 $\Delta t=T/N$ 与 $f=k/T$ , 的采样版本 $X(f)$ 可以近似为

\begin{matrix} (2) & X (\frac{k}{T}) \approx Δ t \sum_{n = 0}^{N - 1} x (n Δ t) e^{- j 2 π k n / N} = Δ t \cdot X_{d} [k] \end{matrix}

$X\left(\frac{k}{T}\right)\approx \Delta t \sum_{n=0}^{N-1}x(n\Delta t)e^{-j2\pi k n/N}=\Delta t \cdot X_d[k]\tag{2}$

在哪里 $X_d[k]$ 是长度 $N$ 的 DFT $x_d[n]=x(n\Delta t)$ .

请注意，对于时间受限的信号，截断误差可以为零，而对于完全带限的信号，采样误差可以为零。由于一个信号不能同时受到时间限制和频带限制，因此总是至少存在两个错误中的一个。在实践中，您通常必须处理这两种类型的错误。

以下 Matlab/Octave 代码显示了一个示例：

Fs = 1e3;   % sampling frequency
Ts = 1/Fs;
T1 = 0.1;
T2 = 2;
tgrid = 0:Ts:T2;
N = length(tgrid);
x = zeros(1,N);
x( find( tgrid <= T1 ) ) = 1;
fgrid = (0:N-1)*Fs/N;

% analytic continuous-time Fourier transform
X = T * sin( T*fgrid*pi ) ./ (T*fgrid*pi) .* exp( -1i*pi*fgrid*T );
X(1) = T;

% DFT approximation
X2 = fft(x,N) * Ts;

plot( fgrid,abs(X),fgrid,abs(X2),'r' )
axis([0,Fs/2,0,T]), grid on
legend('analytic','DFT')

其它你可能感兴趣的问题

上一篇对参数均衡器进行编程 - 防止削波下一篇MUSIC DF 算法中的非零特征值过多