[求助]请教根据声压信号求A计权声压值的MATLAB程序
我自己写了一个，但是求出来跟声级计得到的信号差距还蛮大的有没有哪位能给我提供一个根据声压信号求A计权声压值的MATLAB程序啊
不胜感激！！！
回复 楼主 wangwenting 的帖子
程序如下，请教各位帮忙看看问题在哪儿谢啦！
=oct3bank(noi); %noi是从声级计取出来的声压的电压信号
plot(ff,p);
Cj=[-19.1 -16.1 -13.4 -10.9 -8.6 -6.6 -4.8 -3.2 -1.9 -0.8 0 0.6 1.0 1.2 1.3 1.2 1.0 0.5 -0.1 -1.1 -2.5];
Lj=zeros(1,21);
for n=1:21
Lj(n)=10^(0.1*(p(n)+Cj(n)));
end
Lp(i)=10*log10(sum(Lj));
后面是两个函数，是从网站上面下的程序改了一下
function = oct3bank(x);
pi = 3.14159265358979;
Fs = 44100; % Sampling Frequency
N = 3; % Order of analysis filters.
F = [ 25 31.5 40 50 63 80 100 125 160 200 250 315 400 500 630 800 1000 1250 1600 2000 2500 3150 4000 5000 6300 8000 10000]; % Preferred labeling freq.
ff = (1000).*((2^(1/3)).^[-16:10]); % Exact center freq.
P = zeros(1,27);
m = length(x);
% Design filters and compute RMS powers in 1/3-oct. bands
% 5000 Hz band to 1600 Hz band, direct implementation of filters.
for i = 27:-1:22
= oct3dsgn(ff(i),Fs,N);
y = filter(B,A,x);
P(i) = sum(y.^2)/m;
end
% 1250 Hz to 100 Hz, multirate filter implementation (see ).
= oct3dsgn(ff(24),Fs,N); % Upper 1/3-oct. band in last octave.
= oct3dsgn(ff(23),Fs,N); % Center 1/3-oct. band in last octave.
= oct3dsgn(ff(22),Fs,N); % Lower 1/3-oct. band in last octave.
for j = 6:-1:0
x = decimate(x,2);
m = length(x);
y = filter(Bu,Au,x);
P(j*3+3) = sum(y.^2)/m;
y = filter(Bc,Ac,x);
P(j*3+2) = sum(y.^2)/m;
y = filter(Bl,Al,x);
P(j*3+1) = sum(y.^2)/m;
end
% Convert to decibels.
Pref = 0.00002; % Reference level for dB scale.
idx = (P>0);
P(idx) = 20*log10(P(idx)/Pref);
P(~idx) = NaN*ones(sum(~idx),1);
% Generate the plot
if (nargout == 0)
bar(P);
ax = axis;
axis()
set(gca,'XTick',); % Label frequency axis on octaves.
set(gca,'XTickLabels',F(2:3:length(F)));% MATLAB 4.1c
%set(gca,'XTickLabel',F(2:3:length(F)));% MATLAB 5.1
xlabel('Frequency band '); ylabel('Power ');
title('One-third-octave spectrum')
% Set up output parameters
elseif (nargout == 1)
p = P;
elseif (nargout == 2)
p = P;
f = F;
end
function = oct3dsgn(Fc,Fs,N);
if (nargin > 3) | (nargin < 2)
error('Invalide number of arguments.');
end
if (nargin == 2)
N = 3;
end
if (Fc > 0.88*(Fs/2))
error('Design not possible. Check frequencies.');
end
% Design Butterworth 2Nth-order one-third-octave filter
% Note: BUTTER is based on a bilinear transformation, as suggested in .
pi = 3.14159265358979;
f1 = Fc/(2^(1/6));
f2 = Fc*(2^(1/6));
Qr = Fc/(f2-f1);
Qd = (pi/2/N)/(sin(pi/2/N))*Qr;
alpha = (1 + sqrt(1+4*Qd^2))/2/Qd;
W1 = Fc/(Fs/2)/alpha;
W2 = Fc/(Fs/2)*alpha;
= butter(N,);
页:
[1]