FFTした後、周波数​と時間ごとのパワー強​度のデータが欲しい。

11 visualizzazioni (ultimi 30 giorni)
友一 安藤
友一 安藤 il 22 Mar 2022
Commentato: 友一 安藤 il 24 Mar 2022
FFTをして X軸:周波数 , Y軸 : パワーのプロットはできたのですが、
1秒後毎の、周波数毎のパワー強度の求めからを教えてください。 よろしくお願いします。
サンプリングレート: 250
サンプリング時間:  0.004sec
切り出すデータ時間: 1sec
オーバーラップ: 0
トータル時間: 6min
周波数: 0-125 Hz

Risposte (2)

Shunichi Kusano
Shunichi Kusano il 22 Mar 2022
短時間フーリエ変換のことかと思います。関数がありますので、任意のパラメータ設定してお使いください。

Hernia Baby
Hernia Baby il 22 Mar 2022
Modificato: Hernia Baby il 22 Mar 2022
二つアプローチがありますが、やってることはだいたい同じです。
■信号作成
サンプルに10Hzの信号を作ります
clc,clear;
sec = 6;
Fs = 1/4e-3;
t = (0:1/Fs:sec)';
f1 = 10;
x = sin(2*pi*f1.*t);
figure
hold on
plot(t,x,'Color',[.5 .5 .5])
plot(t,ones(size(t)).*(max(x)/sqrt(2)),'--r')
xlim([0 .2])
xlabel 'time'
stft
0.5秒分ごとに区切り、リーク誤差の軽減のためにハニングウィンドウを使いました
L = 125;
w = hann(L,"periodic");
[s,f,tt]=stft(x,Fs,'Window',w,...
'Overlaplength',0,...
'FrequencyRange','onesided');
窓の正規化をして図示します
figure
c = sum(w)/length(w);
waterfall(f,tt,(((abs(s)/(L/2))/c/sqrt(2)).^2)')
xlim([0 125])
ylim([0 sec])
xlabel('frequency')
ylabel('time')
やっていることはstftと同じです
違うのはブロックサイズの正規化はしなくていい点ですかね
窓の正規化は必要です
[p,fp,tp] = pspectrum(x,Fs,'spectrogram','TimeResolution',.5,'OverlapPercent',0,'Leakage',0.75);
figure
waterfall(fp,tp,p'/c^2/4)
xlim([0 125])
ylim([0 sec])
xlabel('frequency')
ylabel('time')
■結論
1secのデータを中央値である0.5secだと考えれば、
ブロックサイズをサンプリング周波数として取り出すことが可能です
今回だと見栄え上 0.5sec ずつで刻んだので、その刻みごとのデータがとれるというわけですね
  3 Commenti
Hernia Baby
Hernia Baby il 23 Mar 2022
ウィンドウサイズの設定かもしれません 私の書いたコードでいうLはどのようにしましたか?
友一 安藤
友一 安藤 il 24 Mar 2022
L は 時間ではなく、ポイント数にしていたため、 上手く表示出来なかったようです。
アドバイス、ありがとうございました。

Accedi per commentare.

Prodotti


Release

R2021b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!