【问题标题】:period of sawtooth from measurements测量的锯齿周期
【发布时间】:2017-03-21 15:26:58
【问题描述】:

我有一系列 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

【问题讨论】:

  • 计算过零次数?
  • 理论上是的,但是 y 中的度量通常没有那么精确。
  • 永远不会在 x 轴上
  • 我认为这可能与询问如何计算锯齿的平均斜率相同。
  • 如果您知道锯齿的连续部分的平均斜率并知道幅度,您可以从中得出它,正如您所说。

标签: matlab octave period measures


【解决方案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

【讨论】:

  • 哇,谢谢@kabdulla 的回复。直到最近我才自己解决了这个问题。我所做的是取 M 数据(y 轴)的绝对值,将其更改为三角波,然后在 link 的答案中执行该过程。
  • 不用担心。我怀疑您正在使用的库也在进行 dft 以确定基频。只要您的信号没有零(y)偏移,您首先采用绝对值的方法应该没问题。如果存在零偏移量,您可能无法获得正确的周期。
猜你喜欢
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2014-02-18
  • 2010-11-11
  • 1970-01-01
  • 1970-01-01
  • 2011-09-21
  • 1970-01-01
相关资源
最近更新 更多