From: albinali on 25 May 2010 12:21 Hi, I am trying to implement the filter function in Matlab and it is described as "The filter function is implemented as a direct form II transposed structure: y(n) = b(1)*x(n) + b(2)*x(n-1) + ... + b(nb+1)*x(n-nb) - a(2)*y(n-1) - ... - a(na+1)*y(n-na) where n-1 is the filter order, which handles both FIR and IIR filters, na is the feedback filter order, and nb is the feedforward filter order. I wrote the following code to implement it but I am still getting different results. My filtered output seems to be monotonically increasing. I would appreciate it if someone can take a look at my code and let me know what I did wrong. Thanks! Fahd static double[] b = new double[5] {1.0,0.0,-2.0,0.0,1.0}; static double[] a = new double[5] { 1.00, -3.9946673, 5.9840174, -3.9840329, 0.9946827 }; static int NZEROS = 4; static int NPOLES = 4; static int ORDER = 4; static double[] xv = new double[NZEROS + 1]; static double[] yv = new double[NPOLES + 1]; static double[] Filter(double[] data) { double[] filtered = new double[data.Length]; for (int i = 0; (i < data.Length); i++) { xv[0] = xv[1]; xv[1] = xv[2]; xv[2] = xv[3]; xv[3] = xv[4]; xv[4] = data[i]; yv[0] = yv[1]; yv[1] = yv[2]; yv[2] = yv[3]; yv[3] = yv[4]; yv[4] = b[0] * xv[4] + b[1] * xv[3] + b[2] * xv[2] + b[3] * xv[1] + b[4] * xv[0] - a[1] * yv[3] - a[2] * yv[2] - a[3] * yv[1] - a[4] * yv[0]; filtered[i] = yv[4]; } return filtered; }
From: Tim Wescott on 25 May 2010 13:01 On 05/25/2010 09:21 AM, albinali wrote: > Hi, > > I am trying to implement the filter function in Matlab and it is > described as "The filter function is implemented as a direct form II > transposed structure: > > y(n) = b(1)*x(n) + b(2)*x(n-1) + ... + b(nb+1)*x(n-nb) > - a(2)*y(n-1) - ... - a(na+1)*y(n-na) > > where n-1 is the filter order, which handles both FIR and IIR filters, > na is the feedback filter order, and nb is the feedforward filter > order. > > I wrote the following code to implement it but I am still getting > different results. My filtered output seems to be monotonically > increasing. I would appreciate it if someone can take a look at my > code and let me know what I did wrong. Thanks! > > Fahd > > static double[] b = new double[5] {1.0,0.0,-2.0,0.0,1.0}; > static double[] a = new double[5] { 1.00, -3.9946673, > 5.9840174, -3.9840329, 0.9946827 }; First, the unstable pole at z = 1.016 may have something to do with the monotonic (and no doubt exponential) increase. Second, as polynomial order goes up the roots become exceedingly sensitive to variations in the coefficients. This creates all sorts of numerical havoc with direct-form IIR filters. Thus, any time you implement an IIR filter you really, _really_ want to factor it as far down as you can, with second-order filters for the complex pole pairs, with first-order filters for any real poles. Then implement your complete filter as a cascade of the individual filters. Second and a half, your unstable pole may be an unintended consequence of this sensitivity, and the fact that you're expressing the filter as a 4-th order polynomial. HTH -- Tim Wescott Control system and signal processing consulting www.wescottdesign.com
From: Mikolaj on 25 May 2010 13:26 Code looks fine. Your filter is unstable so your results are fine.
|
Pages: 1 Prev: Correlation error derivation? Next: Soft decision decoding of Repetition code |