在LZ给出的SOS matrix和Scale Factors参数中已包含有滤波器的信息。这里给出2种方法形成滤波器参数，这是一个18阶的IIR滤波器，一种是分解成9个2阶IIR滤波器串接，另一种是18阶IIR滤波器：

SOS=[1  -1.753478595924798  1  1  -1.715929423633677  0.989596471431086;...

1  -1.753478595924798  1  1  -1.770907559180061  0.990608353245765;...

1  -1.753478595924798  1  1  -1.702996452119821  0.970500292633477;...

1  -1.753478595924798  1  1  -1.752339882018184  0.973017861982348;...

1  -1.753478595924798  1  1  -1.696640367244520  0.955631016351785;...

1  -1.753478595924798  1  1  -1.733851156227253  0.958455470476445;...

1  -1.753478595924798  1  1  -1.697331544768989  0.946633764727194;...

1  -1.753478595924798  1  1  -1.717331241642967  0.948461822736893;...

1  -1.753478595924798  1  1  -1.704617601062175  0.944269642097509];

Sv=[0.994654917981474

0.994654917981474

0.985575452484598

0.985575452484598

0.978354972045005

0.978354972045005

0.973726850565062

0.973726850565062

0.972134821048755

1.000000000000000];

Tdb=zeros(1,501);

for k=1: 9

B(k,1:3)=SOS(k,1:3)*Sv(k);;

A(k,1:3)=SOS(k,4:6);

[db, mag, pha, w]=freqz_m1(B(k,:),A(k,:));

DB(k,:)=db;

Tdb=Tdb+db;

end

b0=Sv(10);

Tdb=Tdb+20*log10(b0);

figure

plot(w/pi*50,Tdb); grid; axis([0 30 -400 5]);

xlabel('Frequency (kHz)'); ylabel('Magnitude (dB)');

fprintf('B=\n');

for k=1 : 9

fprintf('%4d  %5.6f  %5.6f  %5.6f\n',k,B(k,:));

end

fprintf('\n\n');

fprintf('A=\n');

for k=1 : 9

fprintf('%4d  %5.6f  %5.6f  %5.6f\n',k,A(k,:));

end

fprintf('\n\n');

[b, a]=cas2dir(b0, B, A);

figure

[db, mag, pha, w]=freqz_m1(b,a);

plot(w/pi*50,db); grid; axis([0 30 -400 5]);

xlabel('Frequency (kHz)'); ylabel('Magnitude (dB)');

fprintf('b=\n');

fprintf('%5.6f  %5.6f  %5.6f  %5.6f  %5.6f\n',b)

fprintf('\n\n');

fprintf('a=\n');

fprintf('%5.6f  %5.6f  %5.6f  %5.6f  %5.6f\n',a)

fprintf('\n\n');

其中用到freqz_m1为

function [db, mag, pha, w]=freqz_m1(b,a);

%Modified version of freqz subroutine

[H,w]=freqz(b,a,1000,'whole');

H=(H(1:501))'; w=(w(1:501))';

mag=abs(H);

db=20*log10(mag+eps);

pha=angle(H);

9个2阶IIR滤波器串接和18阶IIR滤波器的响应图如下，它们的参数为：

B=

1  0.994655  -1.744106  0.994655

2  0.994655  -1.744106  0.994655

3  0.985575  -1.728185  0.985575

4  0.985575  -1.728185  0.985575

5  0.978355  -1.715525  0.978355

6  0.978355  -1.715525  0.978355

7  0.973727  -1.707409  0.973727

8  0.973727  -1.707409  0.973727

9  0.972135  -1.704618  0.972135

A=

1  1.000000  -1.715929  0.989596

2  1.000000  -1.770908  0.990608

3  1.000000  -1.702996  0.970500

4  1.000000  -1.752340  0.973018

5  1.000000  -1.696640  0.955631

6  1.000000  -1.733851  0.958455

7  1.000000  -1.697332  0.946634

8  1.000000  -1.717331  0.948462

9  1.000000  -1.704618  0.944270

b=

0.847849  -13.380160  101.477935  -491.012808  1697.382412

-4449.362152  9161.799827  -15148.109667  20383.171259  -22485.628984

20383.171259  -15148.109667  9161.799827  -4449.362152  1697.382412

-491.012808  101.477935  -13.380160  0.847849

a=

1.000000  -15.491945  115.341293  -547.875022  1859.309795

-4784.778030  9672.660185  -15701.246921  20742.817747  -22466.281798

19995.763845  -14590.650556  8664.764207  -4131.835186  1547.759375

-439.646743  89.223268  -11.552358  0.718847

从图上可以看出9个2阶IIR滤波器串接的响应曲线要比18阶IIR滤波器的响应曲线好，这是因为18阶IIR滤波器对系数的精度和计算的精度都有更高的要求，而在PC机中的双精度已不能满足其要求，而产生的误差造成的。

[本帖最后由 songzy41 于 2009-6-20 06:42 编辑]

qz33a.jpg

(32.72 KB, 下载次数: 5)

2009-6-20 06:40 上传

9个2阶IIR滤波器串接的响应曲线

qz33b.jpg

(33.04 KB, 下载次数: 0)

2009-6-20 06:40 上传

18阶IIR滤波器的响应曲线