externals/mrpeach/sqosc~ sqosc~.c,1.1,1.2 - Pd-cvs

29 Mar 2007


      Update of /cvsroot/pure-data/externals/mrpeach/sqosc~
In directory sc8-pr-cvs1.sourceforge.net:/tmp/cvs-serv28684
Modified Files:
    sqosc~.c 
Log Message:
Added finite call to eliminate NaNs
Cleaned up comments
Index: sqosc~.c
===================================================================
RCS file: /cvsroot/pure-data/externals/mrpeach/sqosc~/sqosc~.c,v
retrieving revision 1.1
retrieving revision 1.2
diff -C2 -d -r1.1 -r1.2
*** sqosc~.c	13 Nov 2006 17:54:11 -0000	1.1
--- sqosc~.c	29 Mar 2007 02:41:10 -0000	1.2
***************
*** 1,21 ****
  /* sqosc.c Martin Peach 20060613 based on d_osc.c */
  /* 20060707 using x-x^3/3 to smooth the ramp */
  /* Copyright (c) 1997-1999 Miller Puckette.
  * For information on usage and redistribution, and for a DISCLAIMER OF ALL
  * WARRANTIES, see the file, "LICENSE.txt," in this distribution.  */
! 
! /* sinusoidal oscillator and table lookup; see also tabosc4~ in d_array.c.
! */
#include "m_pd.h"
! #include "math.h"
! #include <stdio.h> /* for file io */
#define UNITBIT32 1572864.  /* 3*2^19; bit 32 has place value 1 */
!     /* machine-dependent definitions.  These ifdefs really
!     should have been by CPU type and not by operating system! */
  #ifdef IRIX
!     /* big-endian.  Most significant byte is at low address in memory */
  #define HIOFFSET 0    /* word offset to find MSB */
  #define LOWOFFSET 1    /* word offset to find LSB */
--- 1,20 ----
  /* sqosc.c Martin Peach 20060613 based on d_osc.c */
  /* 20060707 using x-x^3/3 to smooth the ramp */
+ /* 20070328 added call to finite() to eliminate possible NaNs */
+ /* d_osc.c is: */
  /* Copyright (c) 1997-1999 Miller Puckette.
  * For information on usage and redistribution, and for a DISCLAIMER OF ALL
  * WARRANTIES, see the file, "LICENSE.txt," in this distribution.  */
! /* sinusoidal oscillator and table lookup; see also tabosc4~ in d_array.c. */
#include "m_pd.h"
! #include <math.h>
#define UNITBIT32 1572864.  /* 3*2^19; bit 32 has place value 1 */
! /* machine-dependent definitions.  These ifdefs really
! * should have been by CPU type and not by operating system! */
  #ifdef IRIX
! /* big-endian.  Most significant byte is at low address in memory */
  #define HIOFFSET 0    /* word offset to find MSB */
  #define LOWOFFSET 1    /* word offset to find LSB */
***************
*** 24,31 ****
#ifdef MSW
!     /* little-endian; most significant byte is at highest address */
  #define HIOFFSET 1
  #define LOWOFFSET 0
  #define int32 long
  #endif
--- 23,32 ----
#ifdef MSW
! /* little-endian; most significant byte is at highest address */
  #define HIOFFSET 1
  #define LOWOFFSET 0
  #define int32 long
+ #include <float.h> /* for _finite */
+ #define finite _finite
  #endif
***************
*** 84,89 ****
      double    x_dpw; /* pulse width in this pulse */
      int       x_pulse_ended; /* nonzero if pulse has finished */
- //    FILE      *x_logfp;
- //    int       x_logcount;
  } t_sqosc;
--- 85,88 ----
***************
*** 101,108 ****
      float  *fp, phase, phsinc = (2. * 3.14159) / SQOSCTABSIZE; 
      union  tabfudge tf;
-     //FILE   *cosfp;
if (sqosc_table) return;
-     //cosfp = fopen("sqosctable.txt", "wb");
      sqosc_table = (float *)getbytes(sizeof(float) * (SQOSCTABSIZE+1));
      for (i = SQOSCTABSIZE + 1, fp = sqosc_table, phase = 0; i--;
--- 100,105 ----
***************
*** 110,116 ****
      {
          *fp = cos(phase);
-         //fprintf(cosfp, "%f: %f\n", phase, *fp);
      }
-     //fclose(cosfp);
      /* here we check at startup whether the byte alignment
          is as we declared it.  If not, the code has to be
--- 107,111 ----
***************
*** 147,152 ****
      x->x_dpw = x->x_pw; /* pulse width in this pulse */
      x->x_pulse_ended = 1; /* nonzero if pulse has finished */
- //    x->x_logfp = fopen("sqosclog.txt", "wb");
- //    x->x_logcount = 0;
      return (x);
  }
--- 142,145 ----
***************
*** 243,250 ****
      lastin = *in++; /* latest frequency */
      if (lastin < 0) lastin = -lastin;/* negative frequency is the same as positive here */
!     if (lastin > x->x_bw) lastin = x->x_bw;// limit frequency to bandwidth
      slewindex = x->x_slew*lastin;
      dphase += lastin * conv; /* new phase is old phase + (frequency * table period) */
!     //addr = tab + (tf.tf_i[HIOFFSET] & (SQOSCTABSIZE-1)); /* point to the current sample in the table */
      index = tf.tf_i[HIOFFSET] & (SQOSCTABSIZE-1);
      tf.tf_i[HIOFFSET] = normhipart; /* zero the non-fractional part of the phase */
--- 236,243 ----
      lastin = *in++; /* latest frequency */
      if (lastin < 0) lastin = -lastin;/* negative frequency is the same as positive here */
!     if (lastin > x->x_bw) lastin = x->x_bw;/* limit frequency to bandwidth */
      slewindex = x->x_slew*lastin;
      dphase += lastin * conv; /* new phase is old phase + (frequency * table period) */
!     /*addr = tab + (tf.tf_i[HIOFFSET] & (SQOSCTABSIZE-1)); */ /* point to the current sample in the table */
      index = tf.tf_i[HIOFFSET] & (SQOSCTABSIZE-1);
      tf.tf_i[HIOFFSET] = normhipart; /* zero the non-fractional part of the phase */
***************
*** 253,257 ****
      {
          tf.tf_d = dphase;
!         //f1 = addr[0]; /* first sample */
          if (index <= slewindex)
          { /* rising phase */
--- 246,250 ----
      {
          tf.tf_d = dphase;
!         /*f1 = addr[0]; */ /* first sample */
          if (index <= slewindex)
          { /* rising phase */
***************
*** 263,291 ****
                  x->x_pulse_ended = 0;
              }
!             //findex = (index/(x->x_slew*lastin))*HALFSQOSCTABSIZE;
!             //addr = tab + HALFSQOSCTABSIZE + (int)findex; 
!             f1 = 1.0-2.0*(slewindex-index)/slewindex;// a ramp from -1 to +1 // addr[0];
!             f1 = f1 - pow(f1, 3.0)*onethird;// smooth the ramp
! //            if (x->x_logcount < 1000)
! //            {
! //                fprintf(x->x_logfp, "rise index %d slewindex %f f1 %f frac %f\n", index, slewindex, f1, frac);
! //                ++x->x_logcount;
! //                if (x->x_logcount >= 1000) fclose(x->x_logfp);
! //            }
          }
          else if (index < x->x_dpw) f1 = twothirds; /* risen */
          else if (index <= slewindex+x->x_dpw)
          { /* falling phase */
! //            findex = ((index-HALFSQOSCTABSIZE)/(x->x_slew*lastin))*HALFSQOSCTABSIZE;
! //            addr = tab + (int)findex; 
!             f1 = -1.0+2.0*(slewindex-index+x->x_dpw)/slewindex;// a ramp from +1 to -1 // addr[0];
!             f1 = f1 - pow(f1, 3.0)*onethird;// smooth the ramp
              x->x_pulse_ended = 1;
- //            if (x->x_logcount < 1000)
- //            {
- //                fprintf(x->x_logfp, "fall index %d slewindex %f f1 %f frac %f\n", index, slewindex, f1, frac);
- //                ++x->x_logcount;
- //                if (x->x_logcount >= 1000) fclose(x->x_logfp);
- //            }
          }
          else
--- 256,270 ----
                  x->x_pulse_ended = 0;
              }
!             /*findex = (index/(x->x_slew*lastin))*HALFSQOSCTABSIZE;*/
!             /*addr = tab + HALFSQOSCTABSIZE + (int)findex; */
!             f1 = 1.0-2.0*(slewindex-index)/slewindex;/* a ramp from -1 to +1 */ /* addr[0];*/
!             f1 = f1 - pow(f1, 3.0)*onethird;/* smooth the ramp */
          }
          else if (index < x->x_dpw) f1 = twothirds; /* risen */
          else if (index <= slewindex+x->x_dpw)
          { /* falling phase */
!             f1 = -1.0+2.0*(slewindex-index+x->x_dpw)/slewindex;/* a ramp from +1 to -1 */ /* addr[0];*/
!             f1 = f1 - pow(f1, 3.0)*onethird;/* smooth the ramp */
              x->x_pulse_ended = 1;
          }
          else
***************
*** 295,324 ****
          lastin = *in++;
          if (lastin < 0) lastin = -lastin;/* negative frequency is the same as positive here */
!         if (lastin > x->x_bw) lastin = x->x_bw;// limit frequency to bandwidth
          slewindex = x->x_slew*lastin;
          dphase += lastin * conv; /* next phase */
!         //f2 = addr[1]; /* second sample */
          if (index+1 <= slewindex)
          {
!             f2 = 1.0-2.0*(slewindex-index-1)/slewindex;// addr[1];
              f2 = f2 - pow(f2, 3.0)*onethird;
- //            if (x->x_logcount < 1000)
- //            {
- //                fprintf(x->x_logfp, "rise index %d slewindex %f f2 %f frac %f\n", index+1, slewindex, f2, frac);
- //                ++x->x_logcount;
- //                if (x->x_logcount >= 1000) fclose(x->x_logfp);
- //            }
          }
          else if (index+1 < x->x_dpw) f2 = twothirds;
          else if (index+1 <= slewindex+x->x_dpw)
          {
!             f2 = -1.0+2.0*(slewindex-index-1+x->x_dpw)/slewindex;// addr[1];
              f2 = f2 - pow(f2, 3.0)*onethird;
- //            if (x->x_logcount < 1000)
- //            {
- //                fprintf(x->x_logfp, "fall index %d slewindex %f f2 %f frac %f\n", index+1, slewindex, f2, frac);
- //                ++x->x_logcount;
- //                if (x->x_logcount >= 1000) fclose(x->x_logfp);
- //            }
          }
          else f2 = -twothirds;
--- 274,291 ----
          lastin = *in++;
          if (lastin < 0) lastin = -lastin;/* negative frequency is the same as positive here */
!         if (lastin > x->x_bw) lastin = x->x_bw;/* limit frequency to bandwidth */
          slewindex = x->x_slew*lastin;
          dphase += lastin * conv; /* next phase */
!         /*f2 = addr[1]; */ /* second sample */
          if (index+1 <= slewindex)
          {
!             f2 = 1.0-2.0*(slewindex-index-1)/slewindex;/* addr[1]; */
              f2 = f2 - pow(f2, 3.0)*onethird;
          }
          else if (index+1 < x->x_dpw) f2 = twothirds;
          else if (index+1 <= slewindex+x->x_dpw)
          {
!             f2 = -1.0+2.0*(slewindex-index-1+x->x_dpw)/slewindex;/* addr[1]; */
              f2 = f2 - pow(f2, 3.0)*onethird;
          }
          else f2 = -twothirds;
***************
*** 326,336 ****
          sample = f1 + frac * (f2 - f1); /* output first sample plus fraction of second sample (linear interpolation) */
          *out++ = sample;
! //        if (x->x_logcount < 1000)
! //        {
! //            fprintf(x->x_logfp, "index %ld f1 %f f2 %f frac %f out %f\n", index, f1, f2, frac, sample);
! //            ++x->x_logcount;
! //            if (x->x_logcount >= 1000) fclose(x->x_logfp);
! //        }
!         //addr = tab + (tf.tf_i[HIOFFSET] & (SQOSCTABSIZE-1)); /* point to the next sample */
          index = tf.tf_i[HIOFFSET] & (SQOSCTABSIZE-1);
          tf.tf_i[HIOFFSET] = normhipart; /* zero the non-fractional part */
--- 293,297 ----
          sample = f1 + frac * (f2 - f1); /* output first sample plus fraction of second sample (linear interpolation) */
          *out++ = sample;
!         /* addr = tab + (tf.tf_i[HIOFFSET] & (SQOSCTABSIZE-1)); */ /* point to the next sample */
          index = tf.tf_i[HIOFFSET] & (SQOSCTABSIZE-1);
          tf.tf_i[HIOFFSET] = normhipart; /* zero the non-fractional part */
***************
*** 338,342 ****
          frac = tf.tf_d - UNITBIT32; /* get next fractional part */
      }
!     //f1 = addr[0];
      if (index <= slewindex)
      {
--- 299,303 ----
          frac = tf.tf_d - UNITBIT32; /* get next fractional part */
      }
!     /* f1 = addr[0]; */
      if (index <= slewindex)
      {
***************
*** 348,376 ****
              x->x_pulse_ended = 0;
          }
!         //findex = (index/(x->x_slew*lastin))*HALFSQOSCTABSIZE;
!         //addr = tab + HALFSQOSCTABSIZE + (int)findex; 
!         f1 = 1.0-2.0*(slewindex-index)/slewindex;// addr[0];
          f1 = f1 - pow(f1, 3.0)*onethird;
- //        if (x->x_logcount < 1000)
- //        {
- //            fprintf(x->x_logfp, "rise2 index %d slewindex %f f1 %f frac %f\n", index, slewindex, f1, frac);
- //            ++x->x_logcount;
- //            if (x->x_logcount >= 1000) fclose(x->x_logfp);
- //        }
      }
      else if (index < x->x_dpw) f1 = twothirds; /* risen */
      else if (index <= slewindex+x->x_dpw)
      { /* falling phase */
! //        findex = ((index-HALFSQOSCTABSIZE)/(x->x_slew*lastin))*HALFSQOSCTABSIZE;
! //        addr = tab + (int)findex; 
!         f1 = -1.0+2.0*(slewindex-index+x->x_dpw)/slewindex;// addr[0];
          f1 = f1 - pow(f1, 3.0)*onethird;
          x->x_pulse_ended = 1;
- //        if (x->x_logcount < 1000)
- //        {
- //            fprintf(x->x_logfp, "fall2 index %d slewindex %f f1 %f frac %f\n", index, slewindex, f1, frac);
- //            ++x->x_logcount;
- ///            if (x->x_logcount >= 1000) fclose(x->x_logfp);
- //        }
      }
      else
--- 309,325 ----
              x->x_pulse_ended = 0;
          }
!         /* findex = (index/(x->x_slew*lastin))*HALFSQOSCTABSIZE; */
!         /* addr = tab + HALFSQOSCTABSIZE + (int)findex; */
!         f1 = 1.0-2.0*(slewindex-index)/slewindex;/* addr[0]; */
          f1 = f1 - pow(f1, 3.0)*onethird;
      }
      else if (index < x->x_dpw) f1 = twothirds; /* risen */
      else if (index <= slewindex+x->x_dpw)
      { /* falling phase */
! /*        findex = ((index-HALFSQOSCTABSIZE)/(x->x_slew*lastin))*HALFSQOSCTABSIZE;*/
! /*        addr = tab + (int)findex; */
!         f1 = -1.0+2.0*(slewindex-index+x->x_dpw)/slewindex;/* addr[0]; */
          f1 = f1 - pow(f1, 3.0)*onethird;
          x->x_pulse_ended = 1;
      }
      else
***************
*** 378,414 ****
          f1 = -twothirds;
      }
!     //f2 = addr[1]; /* second sample */
      if (index+1 <= slewindex)
      {
!         f2 = 1.0-2.0*(slewindex-index-1)/slewindex;// addr[1];
          f2 = f2 - pow(f2, 3.0)*onethird;
- //        if (x->x_logcount < 1000)
- //        {
- //            fprintf(x->x_logfp, "rise2 index %d slewindex %f f2 %f frac %f\n", index+1, slewindex, f2, frac);
- //            ++x->x_logcount;
- //            if (x->x_logcount >= 1000) fclose(x->x_logfp);
- //        }
      }
      else if (index+1 < x->x_dpw) f2 = twothirds;
      else if (index+1 <= slewindex+x->x_dpw)
      {
!         f2 = -1.0+2.0*(slewindex-index-1+x->x_dpw)/slewindex;// addr[1];
          f2 = f2 - pow(f2, 3.0)*onethird;
- //        if (x->x_logcount < 1000)
- //        {
- //            fprintf(x->x_logfp, "fall2 index %d slewindex %f f2 %f frac %f\n", index+1, slewindex, f2, frac);
- //            ++x->x_logcount;
- //            if (x->x_logcount >= 1000) fclose(x->x_logfp);
- //        }
      }
      else f2 = -twothirds;
      sample = f1 + frac * (f2 - f1); /* the final sample */
!     *out++ = sample;
! //    if (x->x_logcount < 1000)
! //    {
! //        fprintf(x->x_logfp, "*index %ld f1 %f f2 %f frac %f out %f\n", index, f1, f2, frac, sample);
! //        ++x->x_logcount;
! //        if (x->x_logcount >= 1000) fclose(x->x_logfp);
! //    }
tf.tf_d = UNITBIT32 * SQOSCTABSIZE; /* this just changes the exponent if the table size is a power of 2 */
--- 327,346 ----
          f1 = -twothirds;
      }
!     /* f2 = addr[1]; */ /* second sample */
      if (index+1 <= slewindex)
      {
!         f2 = 1.0-2.0*(slewindex-index-1)/slewindex;/* addr[1];*/
          f2 = f2 - pow(f2, 3.0)*onethird;
      }
      else if (index+1 < x->x_dpw) f2 = twothirds;
      else if (index+1 <= slewindex+x->x_dpw)
      {
!         f2 = -1.0+2.0*(slewindex-index-1+x->x_dpw)/slewindex;/* addr[1];*/
          f2 = f2 - pow(f2, 3.0)*onethird;
      }
      else f2 = -twothirds;
      sample = f1 + frac * (f2 - f1); /* the final sample */
!     if (_finite(sample))*out++ = sample;
!     else *out++ = 0.0;
tf.tf_d = UNITBIT32 * SQOSCTABSIZE; /* this just changes the exponent if the table size is a power of 2 */
***************
*** 421,434 ****
  static void sqosc_dsp(t_sqosc *x, t_signal **sp)
  {
- //    static int once = 0;
x->x_conv = SQOSCTABSIZE/sp[0]->s_sr;
  /* conv = table period = (samples/cycle)/(samples/sec) = sec/cycle = 0.011610sec for 512/44100 */
- //    if (once == 0)
- //    {
- //        ++once;
- //        post ("x->x_slew = %f, x->x_bw = %f, sp[0]->s_sr = %f", x->x_slew, x->x_bw, sp[0]->s_sr);
- //    }
      dsp_add(sqosc_perform, 4, x, sp[0]->s_vec, sp[1]->s_vec, sp[0]->s_n);
  }
--- 353,360 ----
***************
*** 442,446 ****
  {
      if ((pw <= 0)||(pw >= 1)) return;
-         //post("sqosc: pulse width must be greater than 0 and less than 1");// this is an annoying message...
      x->x_pw = pw * SQOSCTABSIZE;
  }
--- 368,371 ----