方法一:%A计权声压级频谱分析
clc;
clear;
close all;
y=wavread('abc.wav');
fs=51200;%采样频率
p0=2e-5;%参考声压
f=[1.00 1.25 1.600 2.00 2.50 3.15 4.00 5.00 6.30 8.0]; %基准中心频率
f1=[20.00 25.0 31.5 40.0 50.0 63.0 80];
fc=[f1,100*f,1000*f,10000*f]; %%%%%%%%%中心频率%%%%%%%%
%20-16000Hz A声级计权值
cf=[-50.5,-44.7,-39.4,-34.6,-30.2,-26.2,-22.5,-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,-4.3,-6.6];
t1=1;
t2=2;
x=y(t1*fs+1:t2*fs);%截取需要处理的数据段
n=length(x);
t=(0:1/fs:(n-1)/fs);
subplot(221);
plot(t,x);%瞬时声压时程图
w=hanning(n); %汉宁窗
xx=1.633*x.*w; %加汉宁窗(恢复系数为1.633)
nfft=2^nextpow2(n);
%nextpow2(n)-取最接近的较大2次幂
a = fft(xx,nfft);
f = fs/2*linspace(0,1,nfft/2);
w=2*abs(a(1:nfft/2)/n);
subplot(222);
plot(f,w);%绘制频谱图
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%1/3倍频程计算
oc6=2^(1/6);
nc=length(cf);
%下面这个求1/3倍频程的程序是按照振动振级计算那个来的
for j=1:nc
fl=fc(j)/oc6;
fu=fc(j)*oc6;
nl=round(fl*nfft/fs+1);
nu&#