周波数特性を知りたい

clear;
[p1,fs]=audioread('40kph_RN.wav');
[p2,fs]=audioread('2000rpm_Knock.wav');
fs=44100;
p2(448512,[1 2]) = [0,0];
smix=p1+p2;
soundsc(smix,fs)

このような音を二つ重ねた時の周波数特性を表示する方法を教えていただきたいです.基本的なことなのですが,久しぶりにMATLABを触ったため忘れてしまいました.

[Kaede]

このとき縦軸はmag2dbを使ってdBに変換したいです.

Answer 1

[Kaede]

clear;
[p1,fs]=audioread('40kph_RN.wav');
[p2,fs]=audioread('2000rpm_Knock.wav');
fs=44100;
N=1024;
p2(448512,[1 2]) = [0,0];
smix=p1+p2;
%soundsc(smix,fs)
p1_fft=fft(p1,N);
p2_fft=fft(p2,N);
A=zeros(1,N/2+1);
B=zeros(1,N/2+1);

frequency=zeros(1,N/2+1);
for k=1:N/2+1;
    A(k)=abs(p1_fft(k));
    B(k)=abs(p2_fft(k));
    frequency(k)=(k-1)*fs/N;
end
dB_A=mag2db(A);
dB_B=mag2db(B);
plot(frequency,dB_A);

hold on
plot(frequency,dB_B);
hold off

自分でも行ってみましたが表示されたグラフ的に違うような気がします.

[takemoto]

合成した後の音の周波数特性を求めるのであれば、"smix"に対してfft演算を行うべきではないでしょうか？

例えば、以下でどうでしょう。

Y = fft(smix,N);

P = abs(Y/N);

frequency = 0 : fs/N : fs - fs/N;

plot(frequency,10*log10(P)), xlim([0,fs/2]), grid

[Kaede]

コメントありがとうございます.

合成した音ということは"smix"をグラフにすると一本のグラフになってしまいますよね？

であれば合成したあとは音をきくだけでよくて周波数特性は合成前の音をそれぞれ表示したいです.

ただデータ数を1024と指定せずに全部のデータに対して行いたいのですがそのような際はどうしたらよいですか.

[takemoto]

別のスレッドで既に解決されている様ですが、fft関数の第二引数を省略することで、入力の全データに対して変換処理を行いますので、この例ですと、

p1_fft=fft(p1,N); ==> p1_fft=fft(p1);

p2_fft=fft(p2,N); ==> p2_fft=fft(p2);

とすればよさそうです。

[Kaede]

コメントありがとうございます.

すでに解決いたしましたがわかりやすく書いていただきありがとうございます.