测量锯齿的周期

Question

我有一系列二维测量结果（x轴上的时间），绘制成非光滑（但非常好）的锯齿波。在理想的世界中，数据点会形成完美的锯齿波（在任一端都有部分振幅数据点）。有没有使用OCTAVE/MATLAB来计算波的（平均）周期的方法？我尝试使用公式维基百科（Sawtooth_wave）锯齿：

P = mean(time.*pi./acot(tan(y./4))), -pi < y < +pi 

也试过：

P = mean(abs(time.*pi./acot(tan(y./4)))) 

，但它没有工作，或至少它给了我，我知道的答案了。

绘制数据的一个例子：

我也尝试下面的方法 - 应该工作 - 但它不给我我所知道的是接近正确的答案。可能是我的代码简单和错误。什么？

slopes = diff(y)./diff(x); % form vector of slopes for each two adjacent points 
for n = 1:length(diff(y)) % delete slope of any two points that form the 'cliff' 
    if abs(diff(y(n,1))) > pi 
    slopes(n,:) = []; 
    end 
    end 
P = median((2*pi)./slopes); % Amplitude is 2*pi 

Answer 1

旧帖子，但认为我会提供我的2美分的价值。我认为有两个合理的方式来做到这一点：

1. 进行傅立叶变换，计算基本
2. 执行阶段，周期，幅度的曲线拟合，并且偏移到理想的方波。

鉴于曲线拟合可能很困难，因为锯齿波的不连续性，所以我建议傅里叶变换。下面是一个独立的示例：

f_s = 10;    # Sampling freq. in Hz 
record_length = 1000; # length of recording in sec. 

% Create noisy saw-tooth wave, with known period and phase 
saw_period = 50; 
saw_phase = 10; 
t = (1/f_s):(1/f_s):record_length; 
saw_function = @(t) mod((t-saw_phase)*(2*pi/saw_period), 2*pi) - pi; 

noise_lvl = 2.0; 
saw_wave = saw_function(t) + noise_lvl*randn(size(t)); 
num_tsteps = length(t); 

% Plot time-series data 
figure(); 
plot(t, saw_wave, '*r', t, saw_function(t)); 
xlabel('Time [s]'); 
ylabel('Measurement'); 
legend('measurements', 'ideal'); 

% Perform fast-Fourier transform (and plot it) 
dft = fft(saw_wave); 
freq = 0:(f_s/length(saw_wave)):(f_s/2); 
dft = dft(1:(length(saw_wave)/2+1)); 

figure(); 
plot(freq, abs(dft)); 
xlabel('Freqency [Hz]'); 
ylabel('FFT of Measurement'); 

% Estimate fundamental frequency: 
[~, idx] = max(abs(dft)); 
peak_f = abs(freq(idx)); 
peak_period = 1/peak_f; 
disp(strcat('Estimated period [s]: ', num2str(peak_period))) 

其中输出一对图表，以及锯齿波的估计周期。你可以玩弄噪音的数量，并看到它正确地获得了50秒的时间，直到非常高的噪音水平。

Estimated period [s]: 50