1

我有一系列 2D 测量值（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
4

1 回答 1

1

旧帖子，但我想我会提供两美分的价值。我认为有两种合理的方法可以做到这一点：

  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
于 2017-03-21T07:26:58.710 回答