フーリエ変換のsin,cos成分について

48 ビュー (過去 30 日間)
ryu
ryu 2022 年 3 月 16 日
コメント済み: ryu 2022 年 3 月 24 日
初学者です。初歩的な質問失礼いたします。
フーリエ変換のさいにsin波とcos波に分解されると思いますが,fft関数を使用した後に、それぞれの成分ごとのスペクトルパワーを取得したいです。方法を教えてくださいますでしょうか。
また加速度計などの解析で使用されるharmonicRatioは偶関数成分(cos)と奇関数成分(sin)の比率と考えておりますが,単純にパワーの割り算でよろしいのでしょうか。
お詳しい方がいらっしゃいましたらお願いいたします。

回答 (2 件)

Hernia Baby
Hernia Baby 2022 年 3 月 17 日
時系列信号でものを言います。
フーリエ変換は時間領域を周波数領域に写像変換します。
例でみてみましょう。10Hzと30Hzの正弦波を合成します。
fft を使う
dt = 0.01;
L = 512;
t = 0:dt:dt*(L-1);
Fs = 1/dt; %サンプリング周波数
y = .5*sin(2*pi*10*t) + 2*sin(2*pi*30*t);
figure
plot(t,y)
xlim([0 t(end)])
xlabel '時間[sec]'
ylabel '信号'
FFTをかけるとどうなるか?複素数に変換されます。
cy = fft(y);
n = 100;
cys = cy([n end-n+2]);
tt = linspace(0,2*pi,100);
% 図示
figure
plot(sin(tt),cos(tt),'--')
axis([-2,2,-2,2])
hold on
plot(cys./(abs(cys)),'o:')
xlabel 'Real'
ylabel 'Imag'
xline(0)
yline(0)
hold off
上記のように複素共役の状態で出されます。
パワースペクトルとはフーリエスペクトルを とすると です。
ここで注意すべきは実効値の二乗でパワースペクトルは算出されるということです。
つまり振幅を A とすると となります。
正規化して以下に示します。
f = Fs*(0:(L/2))/L;
f = f(1:end-1);
P = abs(cy(1:ceil(length(cy)/2)))./(length(y)/2);
P = P.^2;
% 図示
figure
stem(f,P)
hold on
idx = knnsearch(f',[10;30]);
scatter(f(idx),P(idx),'o')
name = round(P(idx),1)';
text(f(idx),P(idx)+0.1,cellstr(num2str(name)));
ylim([0 2.5])
xlabel '周波数[Hz]'
ylabel 'P(f)/Hz'
リーク誤差が発生していますが、上記のように求めることができました。
----
次はもっと簡単な方法を紹介します
pspectrum を使う
ぶっちゃけ計算は1行で終わります。めちゃくちゃすごいです。
[p,f] = pspectrum(y,Fs);
% 図示
figure
plot(f,p,'--');
xlabel '周波数[Hz]'
ylabel 'P1(f)/Hz'
ylim([0 2.5])
hold on
idx = knnsearch(f,[10;30]);
scatter(f(idx),p(idx),'o')
name = round(p(idx),1);
text(f(idx),p(idx)+0.1,cellstr(num2str(name)));
hold off
■Harmonic Ratio(HR)について
  だと思います。
 こちらの 式 (2) をご参照ください。
 とはいえ、harmonicRatio は相関係数のような形をしており、少し違う気がしますが…。
  2 件のコメント
Hernia Baby
Hernia Baby 2022 年 3 月 17 日
前半の長々としたものは、僕のQiita記事からとってます。
Hernia Baby
Hernia Baby 2022 年 3 月 18 日
■上記の例では具体的にsin,cos成分はどれになるのでしょうか
 cos成分とsin成分は実部虚部です。
 パワースペクトルは絶対値をどうこうするので r の部分に相当します。
■この時左右対象なのでHRはほぼ1になるという解釈でよろしいでしょうか
 基本周波数に対して奇数次成分と偶数次成分の比較から歪度をみるんじゃなかったでしたっけ?
 あまり詳しくはないのではっきりとした回答はできませんが…
 基本周波数 の考えでいくと 0 だと思っています。

サインインしてコメントする。


ryu
ryu 2022 年 3 月 17 日
大変丁寧な回答ありがとうございます。
pspectrumは知りませんでした。
上記の例では具体的にsin,cos成分はどれになるのでしょうか
複素共役の図で言うとx軸より左がsin,右がsinという解釈でよろしいでしょうか
またこの時左右対象なのでHRはほぼ1になるという解釈でよろしいでしょうか
  2 件のコメント
Hernia Baby
Hernia Baby 2022 年 3 月 18 日
自分のコメントに記載しました
ryu
ryu 2022 年 3 月 24 日
大変丁寧に説明して頂きありがとうございました!

サインインしてコメントする。

カテゴリ

Help Center および File ExchangeDiscrete Fourier and Cosine Transforms についてさらに検索

Community Treasure Hunt

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

Start Hunting!