4

我有一组周期性的（但不是正弦的）数据。我在一个向量中有一组时间值，在第二个向量中有一组幅度。我想快速估计函数的周期。有什么建议么？

具体来说，这是我当前的代码。我想根据向量 t 来近似向量 x(:,2) 的周期。最终，我想对许多初始条件执行此操作并计算每个条件的周期并绘制结果。

function xdot = f (x,t)
         xdot(1) =x(2);
         xdot(2) =-sin(x(1));
endfunction

x0=[1;1.75];     #eventually, I'd like to try lots of values for x0(2)
t = linspace (0, 50, 200);


x = lsode ("f", x0, t)

plot(x(:,1),x(:,2));

谢谢！

约翰

4

2 回答 2

6

看一下自动相关功能。

来自维基百科

自相关是信号与其自身的互相关。非正式地，观察之间的相似性是它们之间时间间隔的函数。它是一种数学工具，用于查找重复模式，例如是否存在已被噪声掩埋的周期性信号，或识别其谐波频率所暗示的信号中缺失的基频。它通常用于信号处理中，用于分析函数或一系列值，例如时域信号。

Paul Bourke 描述了如何基于快速傅立叶变换（链接）有效地计算自相关函数。

于 2010-04-04T08:43:36.010 回答
2

离散傅里叶变换可以为您提供周期性。更长的时间窗口为您提供更高的频率分辨率，因此我将您的t定义更改为t = linspace(0, 500, 2000). 时域 http://img402.imageshack.us/img402/8775/timedomain.png（这是情节的链接，在托管网站上看起来更好）。你可以这样做：

h = hann(length(x), 'periodic'); %# use a Hann window to reduce leakage
y = fft(x .* [h h]); %# window each time signal and calculate FFT
df = 1/t(end); %# if t is in seconds, df is in Hz
ym = abs(y(1:(length(y)/2), :)); %# we just want amplitude of 0..pi frequency components
semilogy(((1:length(ym))-1)*df, ym);

频域 http://img406.imageshack.us/img406/2696/freqdomain.png 绘图链接。

查看图表，第一个峰值在 0.06 Hz 左右，对应于图中所示的 16 秒周期plot(t,x)。

不过，这在计算上并不是那么快。FFT 是 N*log(N) 次操作。

于 2010-04-04T08:18:53.413 回答