|
|
1/3倍频程程序
function [fc,yc] = octspectrum(yy,fs,fl,fh)
%yy信号;fs采样率;fl为分析频率下限,fh为分析频率上限
f=[1.00 1.25 1.60 2.00 2.50 3.15 4.00 5.00 6.30 8.00];
fc0=[f,10*f,10^2*f,10^3*f,10^4*f,10^5*f,10^6*f,10^7*f,10^8*f];
oc6=2^(1/6);
index = find(fc0 <= fh & fc0 >= fl);
fc = fc0(index);
nc=length(fc);
n=length(yy);
nfft=2^nextpow2(n);
a=fft(yy,nfft);
for j=1:nc,
fl=fc(j)/oc6;
fu=fc(j)*oc6;
nl=round(fl*nfft/fs+1);
nu=round(fu*nfft/fs+1);
b=zeros(1,nfft);
b(nl:nu)=a(nl:nu);
b(nfft-nu+1:nfft-nl+1)=a(nfft-nu+1:nfft-nl+1);
c=ifft(b,nfft);
%%%计算对应每个中心频率段的有效值,与幅度谱相差3分贝
yc(j)=sqrt(var(real(c(1:n))));
end
return |
|
发表于 2007-6-5 10:29
显身卡