[vlc-devel] [PATCH 1/2] Introducing FFT windowing routines for spectrum-based visualizations

Jean-Baptiste Kempf jb at videolan.org
Sat Feb 8 15:40:28 CET 2014


Applied, both.

Btw, it would be nice to have a LGPL FFT in visual, instead of the GPL
one...

On 07 Feb, Ronald Wright wrote :
> This patch introduces FFT windowing routines to reduce spectral leakage in
> spectrum-based visualizations.  It supports four types of windows that are
> commonly used in spectrum analyzers (besides the pre-existing rectangular
> window):  Hann, flat top, Blackman-Harris, and Kaiser.
> ---
>  modules/visualization/visual/window.c         | 217 ++++++++++++++++++++++++++
>  modules/visualization/visual/window.h         |  66 ++++++++
>  modules/visualization/visual/window_presets.h |  37 +++++
>  3 files changed, 320 insertions(+)
>  create mode 100644 modules/visualization/visual/window.c
>  create mode 100644 modules/visualization/visual/window.h
>  create mode 100644 modules/visualization/visual/window_presets.h
> 
> diff --git a/modules/visualization/visual/window.c b/modules/visualization/visual/window.c
> new file mode 100644
> index 0000000..f30ffc8
> --- /dev/null
> +++ b/modules/visualization/visual/window.c
> @@ -0,0 +1,217 @@
> +/*****************************************************************************
> + * window.c : Implementation of FFT window routines
> + *****************************************************************************
> + * Copyright (C) 2014 Ronald Wright
> + * $Id$
> + *
> + * Author: Ronald Wright <logiconcepts819 at gmail.com>
> + *
> + * This program is free software; you can redistribute it and/or modify it
> + * under the terms of the GNU Lesser General Public License as published by
> + * the Free Software Foundation; either version 2.1 of the License, or
> + * (at your option) any later version.
> + *
> + * This program is distributed in the hope that it will be useful,
> + * but WITHOUT ANY WARRANTY; without even the implied warranty of
> + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
> + * GNU Lesser General Public License for more details.
> + *
> + * You should have received a copy of the GNU Lesser General Public License
> + * along with this program; if not, write to the Free Software Foundation,
> + * Inc., 51 Franklin Street, Fifth Floor, Boston MA 02110-1301, USA.
> + *****************************************************************************/
> +
> +#ifdef HAVE_CONFIG_H
> +# include "config.h"
> +#endif
> +
> +#include <math.h>
> +#include "window.h"
> +#include "window_presets.h"
> +
> +/* Flat top window coefficients */
> +#define FT_A0 1.000f
> +#define FT_A1 1.930f
> +#define FT_A2 1.290f
> +#define FT_A3 0.388f
> +#define FT_A4 0.028f
> +
> +/* Blackman-Harris window coefficients */
> +#define BH_A0 0.35875f
> +#define BH_A1 0.48829f
> +#define BH_A2 0.14128f
> +#define BH_A3 0.01168f
> +
> +/*
> + * The modified Bessel function I0(x).  See Chapter 6 of the "Numerical Recipes
> + * in C: The Art of Scientific Computing" book at
> + * http://www.aip.de/groups/soe/local/numres/bookcpdf/c6-6.pdf
> + */
> +static float bessi0( float x )
> +{
> +    float ax, ans;
> +    double y; /* Accumulate polynomials in double precision. */
> +    if( ( ax = fabsf( x ) ) < 3.75f ) /* Polynomial fit. */
> +    {
> +        y = x / 3.75;
> +        y *= y;
> +        ans = 1.0 + y * ( 3.5156229 + y * ( 3.0899424 + y * ( 1.2067492
> +                  + y * ( 0.2659732 + y * ( 0.360768e-1
> +                  + y * 0.45813e-2 ) ) ) ) );
> +    }
> +    else
> +    {
> +        y = 3.75 / ax;
> +        ans = ( exp( ax ) / sqrt( ax ) ) * ( 0.39894228 + y * ( 0.1328592e-1
> +            + y * ( 0.225319e-2 + y * ( -0.157565e-2 + y * ( 0.916281e-2
> +            + y * ( -0.2057706e-1 + y * ( 0.2635537e-1 + y * ( -0.1647633e-1
> +            + y * 0.392377e-2 ) ) ) ) ) ) ) );
> +    }
> +    return ans;
> +}
> +
> +/*
> + * Obtain the window type from the window type variable.
> + */
> +void window_get_param( vlc_object_t * p_aout, window_param * p_param )
> +{
> +    /* Fetch Kaiser parameter */
> +    p_param->f_kaiser_alpha = var_InheritFloat( p_aout, "effect-kaiser-param" );
> +
> +    /* Fetch window type */
> +    char * psz_preset = var_InheritString( p_aout, "effect-fft-window" );
> +    if( !psz_preset )
> +    {
> +        goto no_preset;
> +    }
> +
> +    for( int i = 0; i < NB_WINDOWS; i++ )
> +    {
> +        if( !strcasecmp( psz_preset, window_list[i] ) )
> +        {
> +            p_param->wind_type = i;
> +            return;
> +        }
> +    }
> +
> +no_preset:
> +    msg_Warn( p_aout, "No matching window preset found; using rectangular "
> +                      "window (i.e. no window)" );
> +    p_param->wind_type = NONE;
> +}
> +
> +/*
> + * Initialization routine - sets up a lookup table for scaling a sample of data
> + * by window data.  If the lookup table is successfully allocated, its memory
> + * location and its specified size are stored at the specified memory location
> + * of the internal context.
> + * Returns true if initialization succeeded and returns false otherwise.
> + * The internal context should be freed when it is finished with, by
> + * window_close().
> + */
> +bool window_init( int i_buffer_size, window_param * p_param,
> +                  window_context * p_ctx )
> +{
> +    float * pf_table = NULL;
> +    window_type wind_type = p_param->wind_type;
> +
> +    if( wind_type != HANN && wind_type != FLATTOP
> +                          && wind_type != BLACKMANHARRIS
> +                          && wind_type != KAISER )
> +    {
> +        /* Assume a rectangular window (i.e. no window) */
> +        i_buffer_size = 0;
> +        goto exit;
> +    }
> +
> +    pf_table = malloc( i_buffer_size * sizeof( *pf_table ) );
> +    if( !pf_table )
> +    {
> +        /* Memory allocation failed */
> +        return false;
> +    }
> +
> +    int i_buffer_size_minus_1 = i_buffer_size - 1;
> +    switch( wind_type )
> +    {
> +    case HANN:
> +        /* Hann window */
> +        for( int i = 0; i < i_buffer_size; i++ )
> +        {
> +            float f_val = (float) i / (float) i_buffer_size_minus_1;
> +            pf_table[i] = 0.5f - 0.5f * cosf( 2.0f * (float) M_PI * f_val );
> +        }
> +        break;
> +    case FLATTOP:
> +        /* Flat top window */
> +        for( int i = 0; i < i_buffer_size; i++ )
> +        {
> +            float f_val = (float) i / (float) i_buffer_size_minus_1;
> +            pf_table[i] = FT_A0
> +                        - FT_A1 * cosf( 2.0f * (float) M_PI * f_val )
> +                        + FT_A2 * cosf( 4.0f * (float) M_PI * f_val )
> +                        - FT_A3 * cosf( 6.0f * (float) M_PI * f_val )
> +                        + FT_A4 * cosf( 8.0f * (float) M_PI * f_val );
> +        }
> +        break;
> +    case BLACKMANHARRIS:
> +        /* Blackman-Harris window */
> +        for( int i = 0; i < i_buffer_size; i++ )
> +        {
> +            float f_val = (float) i / (float) i_buffer_size_minus_1;
> +            pf_table[i] = BH_A0
> +                        - BH_A1 * cosf( 2.0f * (float) M_PI * f_val )
> +                        + BH_A2 * cosf( 4.0f * (float) M_PI * f_val )
> +                        - BH_A3 * cosf( 6.0f * (float) M_PI * f_val );
> +        }
> +        break;
> +    case KAISER:
> +    {
> +        /* Kaiser window */
> +        float f_pialph = (float) M_PI * p_param->f_kaiser_alpha;
> +        float f_bessi0_pialph = bessi0( f_pialph );
> +        for( int i = 0; i < i_buffer_size; i++ )
> +        {
> +            float f_val = (float) i / (float) i_buffer_size_minus_1;
> +            float f_term_to_square = 2.0f * f_val - 1.0f;
> +            float f_sqd_term = f_term_to_square * f_term_to_square;
> +            float f_sqr_term = sqrtf( 1.0f - f_sqd_term );
> +            pf_table[i] = bessi0( f_pialph * f_sqr_term ) / f_bessi0_pialph;
> +        }
> +        break;
> +    }
> +    default:
> +        /* We should not reach here */
> +        break;
> +    }
> +
> +exit:
> +    p_ctx->pf_window_table = pf_table;
> +    p_ctx->i_buffer_size = i_buffer_size;
> +    return true;
> +}
> +
> +/*
> + * Perform an in-place scaling of the input buffer by the window data
> + * referenced from the specified context.
> + */
> +void window_scale_in_place( int16_t * p_buffer, window_context * p_ctx )
> +{
> +    for( int i = 0; i < p_ctx->i_buffer_size; i++ )
> +    {
> +        p_buffer[i] *= p_ctx->pf_window_table[i];
> +    }
> +}
> +
> +/*
> + * Free the context.
> + */
> +void window_close( window_context * p_ctx )
> +{
> +    if( p_ctx->pf_window_table )
> +    {
> +        free( p_ctx->pf_window_table );
> +        p_ctx->pf_window_table = NULL;
> +        p_ctx->i_buffer_size = 0;
> +    }
> +}
> diff --git a/modules/visualization/visual/window.h b/modules/visualization/visual/window.h
> new file mode 100644
> index 0000000..d8403aa
> --- /dev/null
> +++ b/modules/visualization/visual/window.h
> @@ -0,0 +1,66 @@
> +/*****************************************************************************
> + * window.h : Header for FFT window routines
> + *****************************************************************************
> + * Copyright (C) 2014 Ronald Wright
> + * $Id$
> + *
> + * Author: Ronald Wright <logiconcepts819 at gmail.com>
> + *
> + * This program is free software; you can redistribute it and/or modify it
> + * under the terms of the GNU Lesser General Public License as published by
> + * the Free Software Foundation; either version 2.1 of the License, or
> + * (at your option) any later version.
> + *
> + * This program is distributed in the hope that it will be useful,
> + * but WITHOUT ANY WARRANTY; without even the implied warranty of
> + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
> + * GNU Lesser General Public License for more details.
> + *
> + * You should have received a copy of the GNU Lesser General Public License
> + * along with this program; if not, write to the Free Software Foundation,
> + * Inc., 51 Franklin Street, Fifth Floor, Boston MA 02110-1301, USA.
> + *****************************************************************************/
> +
> +#ifndef _WINDOW_H_
> +#define _WINDOW_H_
> +
> +#include <vlc_common.h>
> +
> +/* Window type enum */
> +enum _enum_window_type { NONE, HANN, FLATTOP, BLACKMANHARRIS, KAISER };
> +
> +/* Window context structure */
> +struct _struct_window_context {
> +
> +    /* Storage for window function values */
> +    float * pf_window_table;
> +
> +    /* Buffer size for the window */
> +    int i_buffer_size;
> +};
> +
> +typedef enum _enum_window_type window_type;
> +
> +/* Window parameter structure */
> +struct _struct_window_param {
> +
> +    /* Window type */
> +    window_type wind_type;
> +
> +    /* Kaiser window parameter */
> +    float f_kaiser_alpha;
> +};
> +
> +/* Prototypes for the window function */
> +typedef struct _struct_window_context window_context;
> +typedef struct _struct_window_param window_param;
> +void window_get_param( vlc_object_t * p_aout, window_param * p_param );
> +bool window_init( int i_buffer_size, window_param * p_param,
> +                  window_context * p_ctx );
> +void window_scale_in_place( int16_t * p_buffer, window_context * p_ctx );
> +void window_close( window_context * p_ctx );
> +
> +/* Macro for defining a new window context */
> +#define DEFINE_WIND_CONTEXT(name) window_context name = {NULL, 0}
> +
> +#endif /* _WINDOW_H_ */
> diff --git a/modules/visualization/visual/window_presets.h b/modules/visualization/visual/window_presets.h
> new file mode 100644
> index 0000000..09e8d3a
> --- /dev/null
> +++ b/modules/visualization/visual/window_presets.h
> @@ -0,0 +1,37 @@
> +/*****************************************************************************
> + * window_presets.h : Variable names and strings of FFT window presets
> + *****************************************************************************
> + * Copyright (C) 2014 Ronald Wright
> + * $Id$
> + *
> + * Author: Ronald Wright <logiconcepts819 at gmail.com>
> + *
> + * This program is free software; you can redistribute it and/or modify it
> + * under the terms of the GNU Lesser General Public License as published by
> + * the Free Software Foundation; either version 2.1 of the License, or
> + * (at your option) any later version.
> + *
> + * This program is distributed in the hope that it will be useful,
> + * but WITHOUT ANY WARRANTY; without even the implied warranty of
> + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
> + * GNU Lesser General Public License for more details.
> + *
> + * You should have received a copy of the GNU Lesser General Public License
> + * along with this program; if not, write to the Free Software Foundation,
> + * Inc., 51 Franklin Street, Fifth Floor, Boston MA 02110-1301, USA.
> + *****************************************************************************/
> +
> +#ifndef _WINDOW_PRESETS_
> +#define _WINDOW_PRESETS_
> +
> +/* Window functions supported by VLC. These are the typical window types used
> + * in spectrum analyzers. */
> +#define NB_WINDOWS 5
> +static const char * const window_list[NB_WINDOWS] = {
> +    "none", "hann", "flattop", "blackmanharris", "kaiser",
> +};
> +static const char * const window_list_text[NB_WINDOWS] = {
> +    N_("None"), N_("Hann"), N_("Flat Top"), N_("Blackman-Harris"), N_("Kaiser"),
> +};
> +
> +#endif /* _WINDOW_PRESETS_ */
> -- 
> 1.8.3.2
> 
> _______________________________________________
> vlc-devel mailing list
> To unsubscribe or modify your subscription options:
> https://mailman.videolan.org/listinfo/vlc-devel

-- 
With my kindest regards,

-- 
Jean-Baptiste Kempf
http://www.jbkempf.com/ - +33 672 704 734
Sent from my Electronic Device



More information about the vlc-devel mailing list