三分之一倍频程分析

cd20140704

《matlab在振动信号处理中的应用》中关于三分之一倍频程的处理程序见下，对于红色部分有点不理解，请各位大神批评指导下。
对于上下限频率对应的序号nl和nu为何如此求解。后面的滤波的实现原因？麻烦请各位大神指点下，对于三分之一倍频程的认识还是不够。谢谢。相关附件我上传。

%三分之一倍频程处理
clear
clc
close all hidden
format long
%fni=input('三分之一倍频程处理－输入数据文件名','s');
%fid=fopen(fni,'r')
%sf=fscanf(fid,'%f',1);     %读入采样频率值
%fno=fscanf(fid,'%d',1);    %读入输出数据文件名
%x=fscanf(fid,'%f',[1 inf]); %读入输入数据存成列向量
%status=fclose(fid);
sf=500;fno='out6_4.mat';
load y
x=y;
%定义三分之一倍频程的中心频率
f=[1.00 1.25 1.60 2.00 2.50 3.15 4.00 5.00 6.30 8.00];
fc=[f,10*f,100*f,1000*f,10000*f];
%中心频率与下限频率的比值
oc6=2^(1/6);
nc=length(fc);
n=length(x);
nfft=2^nextpow2(n);
a=fft(x,nfft);
for j=1:nc
    fl=fc(j)/oc6;
    fu=fc(j)*oc6;
    nl=round(fl*nfft/sf+1);
    nu=round(fu*nfft/sf+1);
    if fu>sf/2
        m=j-1;break
    end
    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);
    yc(j)=sqrt(var(real(b(1:n))));
end
subplot(2,1,1);
t=0:1/sf:(n-1)/sf;
plot(t,x);
xlabel('时间(s)');
ylabel('加速度(g)');
grid on;
subplot(2,1,2);
plot(fc(1:m),yc(1:m));
xlabel('频率(Hz)');
ylabel('有效值');
grid on;
fid=fopen(fno,'w');
for k=1:m
   fprintf(fid,'%f%f\n',fc(k),yc(k));
end
status=fclose(fid);

simplebinbin

相当于1/3倍频程的频带滤波器，并计算各频带有效值
” yc(j)=sqrt(var(real(b(1:n))))“里面的b(1:n)应改为c(1:n)才正确