other variable names, t_float to t_sample change, array access 64 bit

svn path=/trunk/externals/iemlib/; revision=16193

1 files changed, 27 insertions, 19 deletions

diff --git a/iemlib1/src/biquad_freq_resp.c b/iemlib1/src/biquad_freq_resp.c
index 4fe9d31..aaadfac 100644
--- a/iemlib1/src/biquad_freq_resp.c
+++ b/iemlib1/src/biquad_freq_resp.c
@@ -1,7 +1,7 @@
 /* For information on usage and redistribution, and for a DISCLAIMER OF ALL
 * WARRANTIES, see the file, "LICENSE.txt," in this distribution.
 
-iemlib1 written by Thomas Musil, Copyright (c) IEM KUG Graz Austria 2000 - 2006 */
+iemlib1 written by Thomas Musil, Copyright (c) IEM KUG Graz Austria 2000 - 2011 */
 
 
 #include "m_pd.h"
@@ -21,33 +21,39 @@ typedef struct _biquad_freq_resp
   t_float   b2;
   t_outlet  *x_out_re;
   t_outlet  *x_out_im;
+  t_outlet  *x_out_abs;
+  t_outlet  *x_out_arg;
 } t_biquad_freq_resp;
 
 static t_class *biquad_freq_resp_class;
 
 static void biquad_freq_resp_float(t_biquad_freq_resp *x, t_floatarg freq_phase)
 {
-  t_float re1, im1, re2, im2;
-  t_float c, s, a;
-  t_float pi = 4.0f * atan(1.0f);
+  t_float re1, im1, re2, im2, re, im;
+  t_float c, s, ra;
+  t_float pi = 4.0 * atan(1.0);
   
-  if(freq_phase < 0.0f)
-    freq_phase = 0.0f;
-  else if(freq_phase > 180.0f)
+  if(freq_phase < 0.0)
+    freq_phase = 0.0;
+  else if(freq_phase > 180.0)
     freq_phase = 180.0;
   freq_phase *= pi;
-  freq_phase /= 180.0f;
+  freq_phase /= 180.0;
   
   c = cos(freq_phase);
   s = sin(freq_phase);
   
   re1 = x->a0 + x->a1*c + x->a2*(c*c - s*s);
-  im1 = x->a1*s + x->a2*2.0f*(s*c);
-  re2 = 1.0f - x->b1*c - x->b2*(c*c - s*s);
-  im2 = -x->b1*s - x->b2*2.0f*(s*c);
-  a = re2*re2 + im2*im2;
-  outlet_float(x->x_out_im, (re1*im2 - re2*im1)/a);/* because z^-1 = e^-jwt, negative sign */
-  outlet_float(x->x_out_re, (re1*re2 + im1*im2)/a);
+  im1 = x->a1*s + x->a2*2.0*(s*c);
+  re2 = 1.0 - x->b1*c - x->b2*(c*c - s*s);
+  im2 = -x->b1*s - x->b2*2.0*(s*c);
+  ra = 1.0 / (re2*re2 + im2*im2);
+  im = (re1*im2 - re2*im1) * ra; /* because z^-1 = e^-jwt, negative sign */
+  re = (re1*re2 + im1*im2) * ra;
+  outlet_float(x->x_out_arg, 180.0*atan2(im, re)/pi);
+  outlet_float(x->x_out_abs, sqrt(re*re + im*im));
+  outlet_float(x->x_out_im, im);
+  outlet_float(x->x_out_re, re);
   
 }
 /* y/x = (a0 + a1*z-1 + a2*z-2)/(1 - b1*z-1 - b2*z-2);*/
@@ -69,11 +75,13 @@ static void *biquad_freq_resp_new(void)
   t_biquad_freq_resp *x = (t_biquad_freq_resp *)pd_new(biquad_freq_resp_class);
   x->x_out_re = outlet_new(&x->x_obj, &s_float);
   x->x_out_im = outlet_new(&x->x_obj, &s_float);
-  x->b1 = 0.0f;
-  x->b2 = 0.0f;
-  x->a0 = 0.0f;
-  x->a1 = 0.0f;
-  x->a2 = 0.0f;
+  x->x_out_abs = outlet_new(&x->x_obj, &s_float);
+  x->x_out_arg = outlet_new(&x->x_obj, &s_float);
+  x->b1 = 0.0;
+  x->b2 = 0.0;
+  x->a0 = 0.0;
+  x->a1 = 0.0;
+  x->a2 = 0.0;
   return (x);
 }