Search All of the Math Forum:
Views expressed in these public forums are not endorsed by
Drexel University or The Math Forum.


Wang
Posts:
4
Registered:
8/13/10


Help for the IIR Direct Form II Transposed Implementation
Posted:
Aug 13, 2010 5:50 PM


I've done with a implementation with a IIR Direct Form II Transposed filter, both in MFile code and Fortran 77 code. The algorithm is supposed to be completely the same as shown in the Matlab "Signal Processing Tool User's Guide". which is as follows
The MFile is as follows:  % Function: TDFII(b, a, x). % Direct realizations of rational transfer functions. % Input parameters: % typ: 1 for direct realization, 2 for transposed % b, a: numerator and denominator polynomials % x: input sequence. % Output: % y: output sequence. p = length(a)1; q = length(b)1; pq = max(p,q); a = a(2:p+1); u = zeros(1,pq); % u: the internal state for i = 1:length(x), y(i) = b(1)*x(i)+u(1); u = [u(2:pq),0]; u(1:q) = u(1:q) + b(2:q+1)*x(i); u(1:p) = u(1:p)  a*y(i); end
The Fortran 77 code is as follows:  SUBROUTINE Filter( x, y, SL ) INTEGER SL ! Sample length REAL x(SL), y(SL) ! x and y vector is input and output respectively ! Initialize the input PARAMETER IIR_NSEC = 23 DOUBLE PRECISION in(SL), out(SL) DOUBLE PRECISION CF_A(IIR_NSEC), CF_B(IIR_NSEC), z(IIR_NSEC1) INTEGER m, i, j, k, kk SAVE z
! ! Initialization of filter coefficients ! ! Initializing the filter coeffcients, which should be constant while filtering SAVE CF_A, CF_B DATA CF_A / Here is the coefficeint generated by FDATools /
DATA CF_B / Here is the coefficeint generated by FDATools /
! ! Initialize the input and output vector for filtering ! DO 5, i = 1, SL in(i) = DBLE( x(i) ) out(i) = DBLE( 0 ) 5 CONTINUE
! ! filtering the signal vector !
! Reset past status of delay output DO 10, j = 1, IIR_NSEC1 z(j) = 0.0 10 CONTINUE
! ! Start filtering the kth sample ! Omit CF_A(1) since it has been normalized to 1 DO 30, m = 1, SL out(m) = CF_B(1) * in(m) + z(1) ! Update the past status of delay output DO 20, kk=1, IIR_NSEC2 z(kk) = z(kk+1) 20 CONTINUE z(IIR_NSEC1) = 0.0 ! Calculate the delay outputs z(k) DO 15, k = 1, IIR_NSEC1 z(k) = z(k) + CF_B(k+1)*in(m)  CF_A(k+1)*out(m) 15 CONTINUE 30 CONTINUE
DO 50, im = 1, SL y(SL) = REAL( out(SL) ) 50 CONTINUE
END
I use the IDENTICAL coefficient both in M code and Fortran code, but the M code works fine in Matlab and the filtered signal becomes divergence after only about 19th loop.
I’ve tried the Fortran complier of Intel Visual Fortran 10 and Compaq Digital Fortran 6, almost the same error occurred during runtime.
Is it because the precision of data? I did found some difference between the results come from the same expression in Mcode and Fortran code.
Please help me about the solution !!! This problem becomes a hurdle to my project! So somewhat rush!



