3

我在 matlab 中使用 pwelch 方法来计算一些风速测量的功率谱。因此,到目前为止,我已经编写了以下代码作为示例:

t = 10800; % number of seconds in 3 hours
t = 1:t; % generate time vector
fs = 1; % sampling frequency (seconds)
A = 2; % amplitude
P = 1000; % period (seconds), the time it takes for the signal to repeat itself
f1 = 1/P; % number of cycles per second (i.e. how often the signal repeats itself every second).
y = A*sin(2*pi*f1*t); % signal

fh = figure(1);
set(fh,'color','white','Units', 'Inches', 'Position', [0,0,6,6],...
    'PaperUnits', 'Inches', 'PaperSize', [6,6]);
[pxx, f] = pwelch(y,[],[],[],fs);
loglog(f,10*(pxx),'k','linewidth',1.2);
xlabel('log10(cycles per s)');
ylabel('Spectral Density (dB Hz^{-1})');

我无法包含该情节,因为我没有足够的声望点

这有意义吗?我正在为在情节右侧产生噪音的想法而苦苦挣扎。分解后的信号是一个没有噪声的正弦波,这个噪声是从哪里来的呢?y轴上的值为负的事实是否表明这些频率可以忽略不计?此外,如果风速以 m/s 为单位测量,那么在 y 轴上写单位的最佳方法是什么,这可以转换为对环境科学家更有意义的东西吗?

4

1 回答 1

3

你的结果很好。dB可能会令人困惑。

线性图将获得良好的视野,

Fs = 1000;                    % Sampling frequency
T = 1/Fs;                     % Sample time
L = 1000;                     % Length of signal
t = (0:L-1)*T;                % Time vector
y = sin(2 * pi * 50 * t);     % 50Hz signal

一种fft方法,

NFFT = 2^nextpow2(L); % Next power of 2 from length of y
Y = fft(y,NFFT)/L;
f = Fs/2*linspace(0,1,NFFT/2+1);
subplot(1,2,1);
plot(f,2*abs(Y(1:NFFT/2+1))) 
xlabel('Frequency (Hz)') 
ylabel('|Y(f)|')

pwelch方法,

subplot(1,2,2);
[pxx, freq] = pwelch(y,[],[],[],Fs);
plot(freq,10*(pxx),'k','linewidth',1.2);
xlabel('Frequency (Hz)'); 
ylabel('Spectral Density (Hz^{-1})');

在此处输入图像描述

如您所见,它们的峰值都在50Hz.

用于loglog两者,

在此处输入图像描述

所以“噪音”是存在的,也1e-6存在于fft其中,可以忽略不计。

对于您的第二个问题,我认为轴不会frequency再次改变。因为Fs您应该使用风速的采样频率,例如如果您有10一秒钟内的速度样本,那么您Fs就是10。图表中较高的频率意味着风速的更多变化,而较低的频率表示速度的变化较小。

于 2014-11-22T17:17:11.143 回答