信号処理 演習課題のMatlab

3 views (last 30 days)
Ru
Ru on 27 Jan 2022
Edited: Hernia Baby on 27 Jan 2022
講義中に出た演習問題の解き方を教えて頂きたいです。
この問題と似たものが期末課題として後日でるので、参考にさせて頂きます。

Accepted Answer

Hernia Baby
Hernia Baby on 27 Jan 2022
Edited: Hernia Baby on 27 Jan 2022
すみませんが、以下のものでわからなければ追加でコメントください。課題等は直接解答を書くことができません…あくまでヒントの提示になってしまいます。
  3 Comments
Hernia Baby
Hernia Baby on 27 Jan 2022
Edited: Hernia Baby on 27 Jan 2022
コメントありがとうございます。
あくまで例を書いていきます。
fft関数 とほぼ同じ流れです。
----------------------------------------------
時間の設定
まず時間の刻み幅が dt =0.01としましょう。
そしてサンプリングブロックが512点とします。
ここでわかるのは 秒分のデータをとってきてることです。
時間軸に直すと以下のようになります。
dt = 0.01;
L = 512;
t = 0:dt:dt*(L-1)
t = 1×512
0 0.0100 0.0200 0.0300 0.0400 0.0500 0.0600 0.0700 0.0800 0.0900 0.1000 0.1100 0.1200 0.1300 0.1400 0.1500 0.1600 0.1700 0.1800 0.1900 0.2000 0.2100 0.2200 0.2300 0.2400 0.2500 0.2600 0.2700 0.2800 0.2900
これは0からdt刻みでL個分数値をとってきていることになります。
※L-1は0スタートのためです。
時系列データ
yは正弦波で表現しましょう。
ですから、前述の t を使うことになります。
周波数 10Hz, 30Hz の正弦波を合成します。
Fs = 1/dt; %サンプリング周波数
y = .5*sin(2*pi*10*t) + 2*sin(2*pi*30*t);
plot(t,y)
xlim([0 t(end)])
xlabel '時間[sec]'
ylabel '信号'
FFT
時間軸の世界から周波数軸の世界へ写像変換します。
wikipedia の動画がめちゃくちゃわかりやすいです。
fft関数 は複素数(信号のスペクトル)で出力されます。
図示してちょっと見てみますか。
stem(abs(fft(y)))
xlim([0 L])
なんか対称になってるのがわかりますかね?
これが ナイキストの標本化定理 ってやつとかかわってきます。
右半分はエイリアシングといってストロボ効果が例によくあげられます。
つまり、サンプリング点(サンプリング周波数)の半分未満しか見れません。
これは複素共役が理由です
cy = fft(y);
n = 100;
cys = cy([n end-n+2])
cys =
2.7465 + 1.3219i 2.7465 - 1.3219i
% グラフ化
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'
hold off
正規化
「じゃあ半分の点だけ取ったらオッケーか!やったー!」とはなりません。
写像変換をしたときはだいたい歪みます。
それを0~1のサイズにして整えるのをざっくりと正規化と呼びます。
横軸(周波数軸)と縦軸(スペクトル)を整える必要があります。
まずは横軸から
これは点数分で割ってサンプリング周波数でかければいいですね?
f = Fs*(0:(L/2))/L;
f = f(1:end-1);
縦軸はどうでしょうか?
絶対値だけ題材にします。パーシバルの定理 がここで出てきます。
今回、片側スペクトルをとるので以下のようになります。
P = abs(cy(1:ceil(length(cy)/2)))./(length(y)/2);
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 1.7])
xlabel '周波数[Hz]'
ylabel 'P1(f)'
最後に
何故ヒントにするかというと、課題の答えは質問者のためになりません。
そして見る人が見たら何大学の何先生の課題か分かってしまいます。
大学共通テストでもありました通り、Cheatには厳罰が下されます。
なのでヒントしか差し上げることができないのです。
今後、質問する際は
課題を解くうえで、何が自分にとって問題なのかを明記してください。
きっと優しいユーザーたちが素敵な答えを書いてくれます。
それでは素敵なキャンパスライフを!
Ru
Ru on 27 Jan 2022
すごくご丁寧にありがとうございます!!
分からないことだらけですが1つずつやってみます!
本当にありがとうございました!

Sign in to comment.

More Answers (0)

Categories

Find more on プログラミング in Help Center and File Exchange

Tags

Community Treasure Hunt

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

Start Hunting!