Main Archive Specials Wiki | FAQ Links Submit Forum

 Direct form IIType : genericReferences : Posted by FuzzpilzNotes : I've noticed there's no code for direct form II filters in general here, though probably many of the filter examples use it. I haven't looked at them all to verify that, but there certainly doesn't seem to be a snippet describing this. This is a simple direct form II implementation of a k-pole, k-zero filter. It's a little faster than (a naive, real-time implementation of) direct form I, as well as more numerically accurate.Code : Direct form I pseudocode: y[n] = a[0]*x[n] + a[1]*x[n-1] + .. + a[k]*x[n-k] - b[1]*y[n-1] - .. - b[k]*y[n-k]; Simple equivalent direct form II pseudocode: y[n] = a[0]*x[n] + d[0]; d[0] = a[1]*x[n] - b[1]*y[n] + d[1]; d[1] = a[2]*x[n] - b[2]*y[n] + d[2]; . . d[k-2] = a[k-1]*x[n] - b[k-1]*y[n] + d[k-1]; d[k-1] = a[k]*x[n] - b[k]*y[n]; For example, a biquad: out = a0*in + a1*h0 + a2*h1 - b1*h2 - b2*h3; h1 = h0; h0 = in; h3 = h2; h2 = out; becomes out = a0*in + d0; d0 = a1*in - b1*out + d1; d1 = a2*in - b2*out;

 CommentsAdded on : 21/02/07 by scoofy[ AT ]inf[ DOT ]elte[ DOT ]huComment : I think the per sample denormal number elimination on x87 FPUs is more difficult, since you need to check for denormals at 3 places instead of one (if I'm right).              Added on : 11/06/07 by gtjennings[ AT ]gmail[ DOT ]comComment : Are the constants (a and b) wrong here.  Don't they need to be switched?  If you look at like wikipedia that's the case and it makes more since.  I'm trying to implement a low pass filter at 25mhz passband edge.  I'm getting alot of fluctuation in my output more that expect.  Any suggestions? int main(int argc, char *argv[]) { double b[3] = {1,2,1}; double a1[3] = {1,-1.9995181705254206,0.99952100328066507}; //double a1[3] = {1,-1.9252217796690612,0.95315661147483732}; double a2[3] = {1,-1.9985996261556458,0.99860245760957123}; double a3[3] = {1,-1.9977949691405856,0.99779779945453828}; double a4[3] = {1,-1.9971690447494761,0.99717187417666975}; double a5[3] = {1,-1.9967721889631873,0.9967750178281477}; double a6[2] = {1, -0.99831813425055116}; double d[3] = {0}; double y[5][3] = {0}; double out[2] = {0}; double x[3]={0}, x1,x2,in; double i=0; char wait; while(i<10000000) { x1 = sin(2*10000*3.14159265*i); x2 = sin(2*10000*3.14159265*i-3.14159265); in = x1 * x2; x[0] = in * 7.0818881108085789e-7; y[0][0] = b[0]*x[0] + b[1]*x[1] + b[2]*x[2] - a1[1]*y[0][1] - a1[2]*y[0][1]; y[0][1] = y[0][0]; x[2] = x[1]; x[1] = x[0]; ////////////////////////////////////////////////////////////////////////////////// y[0][0] = y[0][0] * 7.0786348128153693e-7; y[1][0] = b[0]*y[0][0] + b[1]*y[0][1] + b[2]*y[0][2] - a2[1]*y[1][1] - a2[2]*y[1][1]; y[0][2] = y[0][1]; y[0][1] = y[0][0]; y[1][1] = y[1][0]; ////////////////////////////////////////////////////////////////////////////////// y[1][0] = y[1][0] * 7.0757848807506174e-7; y[2][0] = b[0]*y[1][0] + b[1]*y[1][1] + b[2]*y[1][2] - a3[1]*y[2][1] - a3[2]*y[2][1]; y[1][2] = y[1][1]; y[1][1] = y[1][0]; y[2][1] = y[2][0]; ////////////////////////////////////////////////////////////////////////////////// y[2][0] = y[2][0] * 7.0735679834155469e-7; y[3][0] = b[0]*y[2][0] + b[1]*y[2][1] + b[2]*y[2][2] - a4[1]*y[3][1] - a4[2]*y[3][1]; y[2][2] = y[2][1]; y[2][1] = y[2][0]; y[3][1] = y[3][0]; ////////////////////////////////////////////////////////////////////////////////// y[3][0] = y[3][0] * 7.0721624006526327e-007; y[4][0] = b[0]*y[3][0] + b[1]*y[3][1] + b[2]*y[3][2] - a5[1]*y[4][1] - a5[2]*y[4][1]; y[3][2] = y[3][1]; y[3][1] = y[3][0]; y[4][1] = y[4][0]; ////////////////////////////////////////////////////////////////////////////////// /*y[4][0] = y[4][0]* 0.000840932874724457; out[0] = 1*y[4][0] + 1*y[4][1] - a6[1]*out[1]; y[4][1] = y[4][0]; out[1] = out[0];*/ cout<

 Add your own comment Comments are displayed in fixed width, no HTML code allowed! Email: Comment: Are you human?

Site created and maintained by Bram
Graphic design by line.out | Server sponsered by fxpansion