Hi everyone
I am new to matlab and I am working on designing of a cascaded IIR filter
(biquads with 5 coeff each) design. My sampling rate is 48KHz, cut off is 6KHz.
I am currently using fdatool. But I am not sure about its stability. I also
tried using the commands butter and then tf2sos.
Some one kindly suggest how stable these filters can be? Which is a better
method. Kindly provide your suggestions.
Thanking You
Megha Daga
IIR filter design
Started by ●May 18, 2006
Reply by ●May 22, 20062006-05-22
Dear Amit
Can you pl explain what u mean by normalization criteria. I was normalizing by the highest coefficient an dlater on i realised ther was overflow in som ecase, so i normalized by highest value (represnted by 2 bits) by 3. Thats what i did. Kinldy correct me.
The input is going to be a speech signal. I just chose a random input vector. I tested with other vectors too (which had positive and negative values), but the result is same.
Actually amit I am using matlab to design the filter and finally I am going to use these coefficients on a filter which I am implementing on a TI DSP C5416. I tested the coefficients on the TI DSP. i am not getting correct result.
Following are the calculations:
Fs = 48000;
>> HalfFs = Fs/2;
>> Wp = 6000/HalfFs;
>> Ws = 8600/HalfFs;
>> Rp = 3;
>> Rs = 20;
>> [n,Wn] = buttord(Wp,Ws,Rp,Rs)
n
6
Wn
0.2587
>> [b,a] = butter(n,Wn);
>> [sos,g] = tf2sos(b,a)
sos
1.0000 2.0085 1.0085 1.0000 -0.8084 0.1756
1.0000 2.0000 1.0000 1.0000 -0.9087 0.3215
1.0000 1.9915 0.9915 1.0000 -1.1577 0.6836
g
0.0012
x
-1233
424
-933
-933
-1233
-1984
1009
-225
1766
-140
-333
1416
102
-1217
704
1383
-1966
>> R = filter(b,a,x)
R
-1.5358
-13.1021
-52.7019
-135.9341
-260.9256
-414.7120
-591.4483
-783.6995
-942.4701
-962.3276
-742.9542
-291.2615
241.9823
647.3840
786.9343
655.2580
360.7719
Now calculations with normalized vector. I am using a factor of 3 as thats the highest number which can be with 2 bits.
>> sos = sos/3
sos
0.3333 0.6695 0.3362 0.3333 -0.2695 0.0585
0.3333 0.6667 0.3333 0.3333 -0.3029 0.1072
0.3333 0.6638 0.3305 0.3333 -0.3859 0.2279
>> sos = sos*32768
sos
1.0e+004 *
1.0923 2.1938 1.1016 1.0923 -0.8829 0.1918
1.0923 2.1845 1.0923 1.0923 -0.9926 0.3512
1.0923 2.1753 1.0830 1.0923 -1.2645 0.7467
[b,a] = sos2tf(sosnew,g)
R = filter(b,a,x)
R
-1.5358
-13.1021
-52.7019
-135.9341
-260.9256
-414.7120
-591.4483
-783.6995
-942.4701
-962.3276
-742.9542
-291.2615
241.9823
647.3840
786.9343
655.2580
360.7719
So its the same as before.
Now I need to check with the output I am getting at the DSP. In DSP for the computation they multiply the a1 coefficient of the biquad by a factor of 2. Now sos is nothing but:
sos = [b01, b11, b21, 1 , a11, a21
b02, b12, b22, 1, a12, a22
: : : : : :
b0L, b1L, b2L, 1, a1L, a2L]
for biquads. In DSP I need to pass 5 coefficients for each biquads, and hence for each biquad I am passing b0, b1, b2 , a1 and a2. Now as the DSP does calculation with a1 multiplied by 2, i am simply multiplying the a1 in sos by 2 and using that for filter calculations. I have checked this method and its working with other calculations.
sos(1,5) = sos(1,5)*2
sos(2,5) = sos(2,5)*2
sos(3,5) = sos(3,5)*2
sos
1.0e+004 *
1.0923 2.1938 1.1016 1.0923 -1.7659 0.1918
1.0923 2.1845 1.0923 1.0923 -1.9851 0.3512
1.0923 2.1753 1.0830 1.0923 -2.5291 0.7467
>> [b,a] = sos2tf(sos,g)
b
1.0e+010 *
0.1623 0.9739 2.4348 3.2464 2.4348 0.9739 0.1623
a
1.0e+013 *
0.1303 -0.7492 1.5729 -1.4518 0.5665 -0.0918 0.0050
>> R = filter(b,a,x)
R
1.0e+007 *
-0.0000
-0.0000
-0.0000
-0.0000
-0.0001
-0.0004
-0.0012
-0.0029
-0.0070
-0.0162
-0.0362
-0.0790
-0.1693
-0.3573
-0.7445
-1.5361
-3.1435
And if I run that on the DSP, I am getting following output:
-46
-347
-1216
-2740
-4785
-7423
-10612
-13227
-13343
-10014
-4438
1009
4733
6482
6491
5308
Hence i am not getting the correct output. But if sos is in following form
sos
1.0e+004 *
1.0923 2.1938 1.1016 3.2768 -1.7659 0.1918
1.0923 2.1845 1.0923 3.2768 -1.9851 0.3512
1.0923 2.1753 1.0830 3.2768 -2.5291 0.7467
instead of prev one. that is if the 4 col is 1 as in original sos matrix and if g=1, I get same values as I am getting from DSP. That is the output is same.
I hope I am clear. Hence I am wondering am I doing normalization wrong or am I passing wrong values?
Kindly tell me where am I going wrong.
Kindly reply
Thanks
Megha
megha daga wrote:
Amit Shaw wrote:
Hi Megha,
Its a low pass filter and check your normalization criretia (see below).
Also I do not understand why did you select that specific test vector
(all positive values).
Your filter will be BIBO stable (balanced input / balanced outout).
Though I donot find any issue with quantization but you need to check the
intermediates for the practical input signal for saturations.
Where are you using this filter, i mean what is the inout signal expected??
Regards,
Amit
---------------------------------
From: megha daga [mailto:m...@yahoo.com]
Sent: Monday, May 22, 2006 6:22 AM
To: Amit Shaw
Subject: RE: [matlab] IIR filter design
Hi amit
I did normalize and quantize the sos to Q15. But the thing I am not getting any difference in the filtered output of non quantized and quantized coefficients. I checked the pole zero plot in both the case and its the same.
Following are my calculations:
Fs = 48000;
HalfFs = Fs/2;
Wp = 6000/HalfFs;
Ws = 8600/HalfFs;
Rp = 3;
Rs = 20;
[n,Wn] = buttord(Wp,Ws,Rp,Rs)
n
6
Wn
0.2587
[b,a] = butter(n,Wn);
b
0.0012 0.0075 0.0187 0.0249 0.0187 0.0075 0.0012
a
1.0000 -2.8748 3.9031 -3.0192 1.3840 -0.3521 0.0386
sos
1.0000 2.0085 1.0085 1.0000 -0.8084 0.1756
1.0000 2.0000 1.0000 1.0000 -0.9087 0.3215
1.0000 1.9915 0.9915 1.0000 -1.1577 0.6836
g
0.0012
With above values and i/p as
X = [100, 200, 150, 345, 23, 678, 23, 567, 890, 46]
my output came :
R = filter(b,a,X)
R
Columns 1 through 8
0.1246 1.3546 6.9579 22.8703 54.7673 103.0378 160.8801 218.3556
Columns 9 through 10
268.7207 312.9030
Now I divided whole sos by 2.0085 (normalize) and multiplied by 32768 (for
> Q15) and then cal a and b:
Amit >> What is you normalization criteria? Do you want to normalizer the
i/p vs o/p energy to be same.i.e. |H(z)|^2 = 1 or |H(e(jw))|^2 = 1.
In that case norrmalization cosntant will be for coefficients will be
Norm = sum(coeffs of numerator) / sum(coeffs of denominator), you can prove that..
YOu need to check whether SOS normalization is fine!!
sosnew=sos/2.0085
sosnew
0.4979 1.0000 0.5021 0.4979 -0.4025 0.0874
0.4979 0.9958 0.4979 0.4979 -0.4524 0.1601
0.4979 0.9915 0.4937 0.4979 -0.5764 0.3404
sosnew = sosnew*32768
sosnew
1.0e+004 *
1.6315 3.2768 1.6454 1.6315 -1.3188 0.2864
1.6315 3.2629 1.6315 1.6315 -1.4826 0.5246
1.6315 3.2491 1.6177 1.6315 -1.8888 1.1153
[b,a] = sos2tf(sosnew,g)
b
1.0e+011 *
0.0541 0.3245 0.8113 1.0818 0.8113 0.3245 0.0541
a
1.0e+013 *
0.4342 -1.2484 1.6949 -1.3111 0.6010 -0.1529 0.0168
and calculated filter on same data with these values:
R = filter(b,a,X)
R
Columns 1 through 8
0.1246 1.3546 6.9579 22.8703 54.7673 103.0378 160.8801 218.3556
Columns 9 through 10
268.7207 312.9030
So now you can see that, both the times R value is same. I dont understand why the values are not changing. Kindly provide your opinion.
Kinldy reply.
Thanking You
Megha Daga
Amit Shaw wrote:
Hi Megha,
I cannot comment on the data being same at input and output!! May need more inputs to
comment on that.
Can you give some input on the structure of the filter. Because this is important while
doing quantization. Though I have never done is by myself for a filter with coefficients
more than 2, but general trend is that quantization will change the location of the poles
and zeros depending on the Q-format used. You may want to plot the pole zero for the
quantized coefficients. Also it is important that data in the filter path donot get
saturated/underflow so a proper scaling if the stages if required will help improve the
SQNR characteristics of the filter.
Also please update me on your finding, because this way I am also learning with you:)
Btw, what r u working for this, I mean what is the end application? Are u doing some
internship project?
Regards,
Amit
---------------------------------
From: megha daga [mailto:m...@yahoo.com]
Sent: Friday, May 19, 2006 11:36 PM
To: Amit Shaw
Subject: RE: [matlab] IIR filter design
Dear Amit
Thanks fro the reply. Actually I am working on A TI device C5416 and fro that I am designing this filter. Fir that device I need to convert the coefficients into Q15 format (quantize). you Said that I need to be careful while quantization. How can I check the stability of my filter after quantizing. i checked the output data of IIR before quantization and after quantization and its surprisingly the same.
Kindly provide your opinion.
Thanks
megha
Amit Shaw wrote: Hi megha,
You can plot the pole-zero and see whether all the poles are inside
the unit circle for stability.
However in general the IIR filter with fixed point are more prone to
stability issues than the floating point case. So while deciding
the quantization you have to be careful.
Regards,
Amit
---------------------------------
From: m... [mailto:m...] On Behalf Of megha daga
Sent: Wednesday, May 17, 2006 10:14 PM
To: m...
Subject: [matlab] IIR filter design
Hi everyone
I am new to matlab and I am working on designing of a cascaded IIR filter (biquads with 5 coeff each) design. My sampling rate is 48KHz, cut off is 6KHz. I am currently using fdatool. But I am not sure about its stability. I also tried using the commands butter and then tf2sos.
Some one kindly suggest how stable these filters can be? Which is a better method. Kindly provide your suggestions.
Thanking You
Megha Daga
Can you pl explain what u mean by normalization criteria. I was normalizing by the highest coefficient an dlater on i realised ther was overflow in som ecase, so i normalized by highest value (represnted by 2 bits) by 3. Thats what i did. Kinldy correct me.
The input is going to be a speech signal. I just chose a random input vector. I tested with other vectors too (which had positive and negative values), but the result is same.
Actually amit I am using matlab to design the filter and finally I am going to use these coefficients on a filter which I am implementing on a TI DSP C5416. I tested the coefficients on the TI DSP. i am not getting correct result.
Following are the calculations:
Fs = 48000;
>> HalfFs = Fs/2;
>> Wp = 6000/HalfFs;
>> Ws = 8600/HalfFs;
>> Rp = 3;
>> Rs = 20;
>> [n,Wn] = buttord(Wp,Ws,Rp,Rs)
n
6
Wn
0.2587
>> [b,a] = butter(n,Wn);
>> [sos,g] = tf2sos(b,a)
sos
1.0000 2.0085 1.0085 1.0000 -0.8084 0.1756
1.0000 2.0000 1.0000 1.0000 -0.9087 0.3215
1.0000 1.9915 0.9915 1.0000 -1.1577 0.6836
g
0.0012
x
-1233
424
-933
-933
-1233
-1984
1009
-225
1766
-140
-333
1416
102
-1217
704
1383
-1966
>> R = filter(b,a,x)
R
-1.5358
-13.1021
-52.7019
-135.9341
-260.9256
-414.7120
-591.4483
-783.6995
-942.4701
-962.3276
-742.9542
-291.2615
241.9823
647.3840
786.9343
655.2580
360.7719
Now calculations with normalized vector. I am using a factor of 3 as thats the highest number which can be with 2 bits.
>> sos = sos/3
sos
0.3333 0.6695 0.3362 0.3333 -0.2695 0.0585
0.3333 0.6667 0.3333 0.3333 -0.3029 0.1072
0.3333 0.6638 0.3305 0.3333 -0.3859 0.2279
>> sos = sos*32768
sos
1.0e+004 *
1.0923 2.1938 1.1016 1.0923 -0.8829 0.1918
1.0923 2.1845 1.0923 1.0923 -0.9926 0.3512
1.0923 2.1753 1.0830 1.0923 -1.2645 0.7467
[b,a] = sos2tf(sosnew,g)
R = filter(b,a,x)
R
-1.5358
-13.1021
-52.7019
-135.9341
-260.9256
-414.7120
-591.4483
-783.6995
-942.4701
-962.3276
-742.9542
-291.2615
241.9823
647.3840
786.9343
655.2580
360.7719
So its the same as before.
Now I need to check with the output I am getting at the DSP. In DSP for the computation they multiply the a1 coefficient of the biquad by a factor of 2. Now sos is nothing but:
sos = [b01, b11, b21, 1 , a11, a21
b02, b12, b22, 1, a12, a22
: : : : : :
b0L, b1L, b2L, 1, a1L, a2L]
for biquads. In DSP I need to pass 5 coefficients for each biquads, and hence for each biquad I am passing b0, b1, b2 , a1 and a2. Now as the DSP does calculation with a1 multiplied by 2, i am simply multiplying the a1 in sos by 2 and using that for filter calculations. I have checked this method and its working with other calculations.
sos(1,5) = sos(1,5)*2
sos(2,5) = sos(2,5)*2
sos(3,5) = sos(3,5)*2
sos
1.0e+004 *
1.0923 2.1938 1.1016 1.0923 -1.7659 0.1918
1.0923 2.1845 1.0923 1.0923 -1.9851 0.3512
1.0923 2.1753 1.0830 1.0923 -2.5291 0.7467
>> [b,a] = sos2tf(sos,g)
b
1.0e+010 *
0.1623 0.9739 2.4348 3.2464 2.4348 0.9739 0.1623
a
1.0e+013 *
0.1303 -0.7492 1.5729 -1.4518 0.5665 -0.0918 0.0050
>> R = filter(b,a,x)
R
1.0e+007 *
-0.0000
-0.0000
-0.0000
-0.0000
-0.0001
-0.0004
-0.0012
-0.0029
-0.0070
-0.0162
-0.0362
-0.0790
-0.1693
-0.3573
-0.7445
-1.5361
-3.1435
And if I run that on the DSP, I am getting following output:
-46
-347
-1216
-2740
-4785
-7423
-10612
-13227
-13343
-10014
-4438
1009
4733
6482
6491
5308
Hence i am not getting the correct output. But if sos is in following form
sos
1.0e+004 *
1.0923 2.1938 1.1016 3.2768 -1.7659 0.1918
1.0923 2.1845 1.0923 3.2768 -1.9851 0.3512
1.0923 2.1753 1.0830 3.2768 -2.5291 0.7467
instead of prev one. that is if the 4 col is 1 as in original sos matrix and if g=1, I get same values as I am getting from DSP. That is the output is same.
I hope I am clear. Hence I am wondering am I doing normalization wrong or am I passing wrong values?
Kindly tell me where am I going wrong.
Kindly reply
Thanks
Megha
megha daga wrote:
Amit Shaw wrote:
Hi Megha,
Its a low pass filter and check your normalization criretia (see below).
Also I do not understand why did you select that specific test vector
(all positive values).
Your filter will be BIBO stable (balanced input / balanced outout).
Though I donot find any issue with quantization but you need to check the
intermediates for the practical input signal for saturations.
Where are you using this filter, i mean what is the inout signal expected??
Regards,
Amit
---------------------------------
From: megha daga [mailto:m...@yahoo.com]
Sent: Monday, May 22, 2006 6:22 AM
To: Amit Shaw
Subject: RE: [matlab] IIR filter design
Hi amit
I did normalize and quantize the sos to Q15. But the thing I am not getting any difference in the filtered output of non quantized and quantized coefficients. I checked the pole zero plot in both the case and its the same.
Following are my calculations:
Fs = 48000;
HalfFs = Fs/2;
Wp = 6000/HalfFs;
Ws = 8600/HalfFs;
Rp = 3;
Rs = 20;
[n,Wn] = buttord(Wp,Ws,Rp,Rs)
n
6
Wn
0.2587
[b,a] = butter(n,Wn);
b
0.0012 0.0075 0.0187 0.0249 0.0187 0.0075 0.0012
a
1.0000 -2.8748 3.9031 -3.0192 1.3840 -0.3521 0.0386
sos
1.0000 2.0085 1.0085 1.0000 -0.8084 0.1756
1.0000 2.0000 1.0000 1.0000 -0.9087 0.3215
1.0000 1.9915 0.9915 1.0000 -1.1577 0.6836
g
0.0012
With above values and i/p as
X = [100, 200, 150, 345, 23, 678, 23, 567, 890, 46]
my output came :
R = filter(b,a,X)
R
Columns 1 through 8
0.1246 1.3546 6.9579 22.8703 54.7673 103.0378 160.8801 218.3556
Columns 9 through 10
268.7207 312.9030
Now I divided whole sos by 2.0085 (normalize) and multiplied by 32768 (for
> Q15) and then cal a and b:
Amit >> What is you normalization criteria? Do you want to normalizer the
i/p vs o/p energy to be same.i.e. |H(z)|^2 = 1 or |H(e(jw))|^2 = 1.
In that case norrmalization cosntant will be for coefficients will be
Norm = sum(coeffs of numerator) / sum(coeffs of denominator), you can prove that..
YOu need to check whether SOS normalization is fine!!
sosnew=sos/2.0085
sosnew
0.4979 1.0000 0.5021 0.4979 -0.4025 0.0874
0.4979 0.9958 0.4979 0.4979 -0.4524 0.1601
0.4979 0.9915 0.4937 0.4979 -0.5764 0.3404
sosnew = sosnew*32768
sosnew
1.0e+004 *
1.6315 3.2768 1.6454 1.6315 -1.3188 0.2864
1.6315 3.2629 1.6315 1.6315 -1.4826 0.5246
1.6315 3.2491 1.6177 1.6315 -1.8888 1.1153
[b,a] = sos2tf(sosnew,g)
b
1.0e+011 *
0.0541 0.3245 0.8113 1.0818 0.8113 0.3245 0.0541
a
1.0e+013 *
0.4342 -1.2484 1.6949 -1.3111 0.6010 -0.1529 0.0168
and calculated filter on same data with these values:
R = filter(b,a,X)
R
Columns 1 through 8
0.1246 1.3546 6.9579 22.8703 54.7673 103.0378 160.8801 218.3556
Columns 9 through 10
268.7207 312.9030
So now you can see that, both the times R value is same. I dont understand why the values are not changing. Kindly provide your opinion.
Kinldy reply.
Thanking You
Megha Daga
Amit Shaw wrote:
Hi Megha,
I cannot comment on the data being same at input and output!! May need more inputs to
comment on that.
Can you give some input on the structure of the filter. Because this is important while
doing quantization. Though I have never done is by myself for a filter with coefficients
more than 2, but general trend is that quantization will change the location of the poles
and zeros depending on the Q-format used. You may want to plot the pole zero for the
quantized coefficients. Also it is important that data in the filter path donot get
saturated/underflow so a proper scaling if the stages if required will help improve the
SQNR characteristics of the filter.
Also please update me on your finding, because this way I am also learning with you:)
Btw, what r u working for this, I mean what is the end application? Are u doing some
internship project?
Regards,
Amit
---------------------------------
From: megha daga [mailto:m...@yahoo.com]
Sent: Friday, May 19, 2006 11:36 PM
To: Amit Shaw
Subject: RE: [matlab] IIR filter design
Dear Amit
Thanks fro the reply. Actually I am working on A TI device C5416 and fro that I am designing this filter. Fir that device I need to convert the coefficients into Q15 format (quantize). you Said that I need to be careful while quantization. How can I check the stability of my filter after quantizing. i checked the output data of IIR before quantization and after quantization and its surprisingly the same.
Kindly provide your opinion.
Thanks
megha
Amit Shaw wrote: Hi megha,
You can plot the pole-zero and see whether all the poles are inside
the unit circle for stability.
However in general the IIR filter with fixed point are more prone to
stability issues than the floating point case. So while deciding
the quantization you have to be careful.
Regards,
Amit
---------------------------------
From: m... [mailto:m...] On Behalf Of megha daga
Sent: Wednesday, May 17, 2006 10:14 PM
To: m...
Subject: [matlab] IIR filter design
Hi everyone
I am new to matlab and I am working on designing of a cascaded IIR filter (biquads with 5 coeff each) design. My sampling rate is 48KHz, cut off is 6KHz. I am currently using fdatool. But I am not sure about its stability. I also tried using the commands butter and then tf2sos.
Some one kindly suggest how stable these filters can be? Which is a better method. Kindly provide your suggestions.
Thanking You
Megha Daga
