2

我在将频谱转换为时间序列时遇到了一个小问题。我已经阅读了很多文章，我认为我正在应用正确的程序，但我没有得到正确的结果。你能帮忙找出错误吗？

我有一个时间序列，例如： 在此处输入图像描述

当我计算频谱时，我会： %number of points nPoints=length(timeSeries);

%time interval
dt=time(2)-time(1);

%Fast Fourier transform
p=abs(fft(timeSeries))./(nPoints/2);

%power of positive frequencies
spectrum=p(1:(nPoints/2)).^2;

%frequency
dfFFT=1/tDur;
frequency=(1:nPoints)*dfFFT;
frequency=frequency(1:(nPoints)/2);

%plot spectrum
semilogy(frequency,spectrum); grid on;
xlabel('Frequency [Hz]');
ylabel('Power Spectrum [N*m]^2/[Hz]');
title('SPD load signal');

我得到：

在此处输入图像描述

我认为频谱计算得很好。但是，现在我需要返回并从该频谱中获取时间序列，我会这样做：

df=frequency(2)-frequency(1);
ap = sqrt(2.*spectrum*df)';
%random number form -pi to pi    
epsilon=-pi + 2*pi*rand(1,length(ap));
%transform to time series
randomSeries=length(time).*real(ifft(pad(ap.*exp(epsilon.*i.*2.*pi),length(time))));
%Add the mean value
randomSeries=randomSeries+mean(timeSeries);

但是，情节看起来像：

在此处输入图像描述

它比原始系列低一个数量级。有什么推荐吗？

4

1 回答 1

4

这里（至少）发生了两件事。首先是您丢弃信息，然后用随机数代替该信息。

实数序列的 FFT 是由实部和虚部组成的复数序列。将这些数字转换为极坐标形式可以为您提供幅度和相位角。您使用 捕获幅度部分p=aps(fft(...))，但没有捕获相位角（这将涉及atan2(...)）。然后，您正在编造随机数（epsilon=...) 并在您重建时间序列时使用它们替换原始数字。此外，由于实数序列的 FFT 具有特定的对称性，用随机数代替相位角会破坏这种对称性，这意味着 IFFT 通常不再是实数序列，而是复数序列 - 再次，您只查看 IFFT 的真实部分，因此您再次丢弃信息。如果这是一个音频信号，结果可能听起来有点像原始信号（或者它们可能完全不同），但波形肯定不匹配......

第二个问题是，在许多实现中，ifft(fft(...))将通过信号中的点数来缩放结果。有几种不同的方法可以避免这种情况，产生不同的结果，但有时在不同的情况下更有吸引力，这取决于你想要做什么。您可以在fft()执行之前缩放结果，也可以在最后ifft()缩放ifft()结果，或者在某些情况下，我什至看到两者都被缩放了一个因子sqrt(N)- 执行两次的最终结果是将最终结果缩放N，但它的效率有点低，因为你做了两次缩放......

于 2014-01-22T15:14:08.023 回答