バンドパスフィルター後の音圧測定について

特定の周波数帯域のみの音圧を測定するために，バンドパスフィルターをかけ，音圧（RMS値）を算出しています．

以下に，コードを記します．

－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－

[y, fs] = audioread('demo.wav');

cal = 176.8; %絶対音圧に補正するためのキャリブレーション値

calin = 10.^(abs(cal)/20); %絶対音圧に補正するための式

data_row = y*calin; %絶対音圧に補正

bandpass(data_row, [80000 120000], fs); %作図

[p, f] = bandpass(data_row, [80000 120000], fs);

pdata_bp = abs(p).^2;

RMS_bp = 10*log10(sum(pdata_bp)/length(pdata_bp));

table(RMS_bp)

－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－－

〈ご質問〉

bandpass関数を用いた作図の音圧（60 dB程度）と，RMSの音圧（RMS_bp）（82.9 dB）の値が大きく異なってしまいます．

解析に用いた音データと，作図の結果を添付いたします．

なぜこのような違いが生じてしまうのか，ご存じの方がいらっしゃいましたら，ご教示いただけますと幸いです．

また，どこをどのように修正するべきか，ご教示いただきたく存じます．

よろしくお願いいたします．

Answer 1

[Hernia Baby]

見ているものが違うからです。

@Tomoさんが見ているのは合計値すなわちOverAllです。

これは全ての周波数帯におけるパワーを合算した数値です。

clear,clc;
load('demo.mat')

cal = 176.8;   %絶対音圧に補正するためのキャリブレーション値
calin = 10.^(abs(cal)/20);   %絶対音圧に補正するための式
data_row = y*calin;   %絶対音圧に補正

bandpass(data_row, [80000 120000], fs);    %作図

バンドパスをかけたものを見てみましょう。

[p, f] = bandpass(data_row, [80000 120000], fs);
% pdata_bp = abs(p).^2;
% RMS_bp = 10*log10(sum(pdata_bp)/length(pdata_bp));
RMS_bp = 10*log10(rms(p).^2);
table(RMS_bp)

ans = table
    RMS_bp
    ______

    82.901

ここで、RMS_bpはOverAllです。

----

周波数帯での結果をpspectrumで確認します。

[pp0,ff0] = pspectrum(data_row,fs,"Leakage",0.85);
[pp,ff] = pspectrum(p,fs,"Leakage",0.85);
figure
plot(ff0./1e3,10*log10(pp0),ff./1e3,10*log10(pp))
grid on
legend({"original","Filtered"})
ylim([0 120])

音圧が異なっていないことがわかりました。

[Hernia Baby]

※補足：コメントアウトした部分と書き換えた部分は同じ操作です。

[Tomo]

@Hernia Baby様

ご親切にご回答いただきまして，ありがとうございました．

Over allが原因であることは分かりました．

それでは，特定の周波数帯のみのRMSを計算するにはどのようにすればよいのでしょうか？

お手数をおかけいたしますが，何卒よろしくお願い申し上げます．

[Hernia Baby]

@Tomoさん、返信ありがとうございます。

今回された操作は80k~120kHzのバンドパスフィルタをかけた時系列データのRMS値をパワーで算出しています。なのでRMS_bpがこの周波数帯でのRMS値(の二乗をレベルにしたもの)になります。

区間を指定した Over All は Partial Over All (P.O.A)と呼ばれます。

単にこの周波数帯の実効値を求めたい場合はrmsをつけるだけでオッケーです。

参考：オーバオール値とは何ですか、またどのように算出するのか？

clear,clc;
load('demo.mat')
cal = 176.8;                % 絶対音圧に補正するためのキャリブレーション値
calin = 10.^(abs(cal)/20);  % 絶対音圧に補正するための式
data_row = y*calin;         % 絶対音圧に補正

[p, f] = bandpass(data_row, [80000 120000], fs);    % バンドパスフィルタ
RMS_bp = 10*log10(rms(p).^2)                        % P.O.A(レベル)

RMS_bp = 82.9013

figure
plot(p)
hold on
yline(rms(p),'k:')
legend({'フィルタ済み信号','実効値(RMS値)'})

[Tomo]

@Hernia Baby様

ご親切にご回答いただきまして，ありがとうございました．

理解いたしました．ありがとうございました．