Csound Csound-dev Csound-tekno Search About

[CSOUND-DEV:5225] gen20

Date2004-08-18 09:41
FromTrond Lossius
Subject[CSOUND-DEV:5225] gen20
Hi!

I'm currently porting the gen20 window functions from fgens.c to Java in
order to use them in Max 4.5 inside the mxj object. While doing this I'm
wondering about a detail of the implementations:

In spite of the mathematical functions for the various windows being
symetrical, the implementation is not. E.g. Hamming, Hanning, etc. are
calculated this way:

    for (i = 0.0 ; i < TWOPI ; i += arg)
      *ft++ = (MYFLT)(xarg * (cf[0] - cf[1]*cos(i) +
                              cf[2]*cos(2.0*i) - cf[3]*cos(3.0*i)));

In this loop i will stop just short of TWOPI, and thus the last sample do
not equal the first.

The same goes for the Gaussian window:

    case 6:                     /* Gaussian */
        arg = 12.0/flen;
        for (i = -6.0 ; i < 0.0 ; i += arg)
          *ft++ = (MYFLT)(xarg * (pow( 2.71828, -(i*i)/2.0)));
        for (i = arg ; i < 6.0 ; i += arg)
          *ft++ = (MYFLT)(xarg * (pow( 2.71828, -(i*i)/2.0)));
        return;

The counter will be running from -6 (inclusive) to +6 (exclusive)

It probably doesn't make much of a difference for large buffers, but
still...

Best,
Trond

Date2004-08-18 10:41
FromTrond Lossius
Subject[CSOUND-DEV:5226] Re: gen20
Still working on the code for window functions. I've now discovered a bug in
the code for the sinc window. This is current code:

    case 9:                     /* Sinc */
        arg = TWOPI / flen;
        for (i = -PI ; i < 0 ; i += arg)
          *ft++ = (MYFLT)(xarg * sin(i) / i);
        *ft++ = FL(1.0);
        for (i = arg ; i < PI ; i += arg)
          *ft++ = (MYFLT)(xarg * sin(i) / i);
        return;

The reason why the *ft++ = FL(1.0); is used in between the two loops is that
the calculation ( sin(i) / i ) is illegal for i=0. From math calculus (using
L'Hopital's rule) it is known that

lim ( sin(i) / i) = 0  when  i -> 0

In gen20 amplitude xarg is implemented, and hence this value has to be
scaled:

*ft++ = FL(1.0 * xarg);

The current bug can be seen clearly in Kevin's illustration of the various
window functions: 

http://www.kevindumpscore.com/docs/csound-manual/miscwindows.html

In the sinc example, max amplitude = 0.75, but at i=0 there's an erraneous
spike hitting the value 1.

Corrected code should read:

    case 9:                     /* Sinc */
        arg = TWOPI / flen;
        for (i = -PI ; i < 0 ; i += arg)
          *ft++ = (MYFLT)(xarg * sin(i) / i);
        *ft++ = FL(xarg); // <------- This line is changed
        for (i = arg ; i < PI ; i += arg)
          *ft++ = (MYFLT)(xarg * sin(i) / i);
        return;


Best,
Trond



On 2004-08-18 10:41, "Trond Lossius"  wrote:

> Hi!
> 
> I'm currently porting the gen20 window functions from fgens.c to Java in order
> to use them in Max 4.5 inside the mxj object. While doing this I'm wondering
> about a detail of the implementations:
> 
> In spite of the mathematical functions for the various windows being
> symetrical, the implementation is not. E.g. Hamming, Hanning, etc. are
> calculated this way:
> 
>   for (i = 0.0 ; i < TWOPI ; i += arg)
>     *ft++ = (MYFLT)(xarg * (cf[0] - cf[1]*cos(i) +
>                             cf[2]*cos(2.0*i) - cf[3]*cos(3.0*i)));
> 
> In this loop i will stop just short of TWOPI, and thus the last sample do not
> equal the first.
> 
> The same goes for the Gaussian window:
> 
>   case 6:                     /* Gaussian */
>       arg = 12.0/flen;
>       for (i = -6.0 ; i < 0.0 ; i += arg)
>         *ft++ = (MYFLT)(xarg * (pow( 2.71828, -(i*i)/2.0)));
>       for (i = arg ; i < 6.0 ; i += arg)
>         *ft++ = (MYFLT)(xarg * (pow( 2.71828, -(i*i)/2.0)));
>       return;
> 
> The counter will be running from -6 (inclusive) to +6 (exclusive)
> 
> It probably doesn't make much of a difference for large buffers, but still...
> 
> Best,
> Trond

Date2004-08-18 11:35
Fromjpff@codemist.co.uk
Subject[CSOUND-DEV:5227] Re: gen20
The line 
*ft++ = FL(1.0 * xarg);
is not correct.  Should be 
*ft++ = FL(1.0) * xarg;
or just 
*ft++ = xarg;

==John ffitch