[vlc-devel] [PATCH 1/6] Add monotone cubic interpolation support
Rémi Denis-Courmont
remi at remlab.net
Sat Nov 28 10:38:24 CET 2015
On Saturday 28 November 2015 00:31:01 Ronald Wright wrote:
> ---
> include/vlc_interpolation.h | 129 +++++++++++++++++++++++++++
> po/POTFILES.in | 4 +
> src/Makefile.am | 4 +
> src/libvlccore.sym | 6 ++
> src/misc/interpolation.c | 212
> ++++++++++++++++++++++++++++++++++++++++++++ 5 files changed, 355
> insertions(+)
> create mode 100644 include/vlc_interpolation.h
> create mode 100644 src/misc/interpolation.c
>
> diff --git a/include/vlc_interpolation.h b/include/vlc_interpolation.h
> new file mode 100644
> index 0000000..18aa353
> --- /dev/null
> +++ b/include/vlc_interpolation.h
> @@ -0,0 +1,129 @@
> +/**************************************************************************
> *** + * vlc_interpolation.h: Interpolation utility
> +
> ***************************************************************************
> ** + * Copyright (C) 2015 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 VLC_INTERPOLATION_H
> +#define VLC_INTERPOLATION_H
> +
> +/**
> + * \file
> + * This file is the interface definition for a modification of
> Fritsch-Carlson + * monotone cubic spline interpolation (implementation in
> + * src/misc/interpolation.c)
> + */
Probably too long for Doxygen. Use multiple paragraphs if needed.
> +
> +/**************************************************************************
> *** + * Documentation
> +
> ***************************************************************************
> **/ +/*
> + **** Background
> + *
> + * This implements a modification of Fritsch-Carlson monotone cubic spline
> + * interpolation. Fritsch-Carlson interpolation creates an interpolant as
> + * follows:
> + *
> + * 1. The slopes between successive data points are computed.
> + * 2. The tangents at each data point is computed as the average of the
> + * secants.
> + * 3. To preserve monotonicity, slopes are set to zero wherever two
> successive + * points are equal.
> + * 4. To prevent overshoot, two successive slopes, normalized by the
> tangent + * at the leftmost point, are restricted so that the sum of
> their squares do + * not exceed a certain threshold (9 in this
> algorithm).
> + *
> + * The interpolation is modified such that all sample points left of the
> spline + * have values equal to the left-most data point and all sample
> points right of + * the spline have values equal to the right-most data
> point. Hence, the + * resulting line extending outward from the left-most
> data point has zero + * slope, and the resulting line extending outward
> from the right-most data + * point has zero slope.
> + *
> + **** Example usage
> + *
> + * #define NUM_POINTS 4
> + * ...
> + * my_function()
> + * {
> + * ...
> + * float points[NUM_POINTS] = { -0.5, 1.0, 2.5, 5.0 };
> + * float values[NUM_POINTS] = { 8.0, -5.0, 0.0, 3.0 };
> + * ...
> + * vlc_interpolant_t itp;
> + * if( vlc_create_interpolant( &itp, points, values, NUM_POINTS ) ==
> VLC_SUCCESS ) + * {
> + * ...
> + * float y = vlc_interpolate( itp, 0.5 );
> + * ...
> + * vlc_destroy_interpolant( itp );
> + * ...
> + * }
> + * ...
> + * }
> + * */
Consider Doxygen formatting for \file or \defgroup.
> +
> +/**************************************************************************
> *** + * Interpolant type
> +
> ***************************************************************************
> **/ +
> +/**
> + * The interpolant data type. This is a pointer to a private structure
> defined
> + * in misc/interpolation.c.
> + */
Missing new lines. Doxygen won´t format well.
> +typedef struct _vlc_interpolant_t * vlc_interpolant_t;
No pointer type please.
> +
> +/**************************************************************************
> *** + * Interpolation routines
> +
> ***************************************************************************
> **/ +
> +/**
> + * Create a new interpolant from the given data points and data point
> values.
> + *
> + * \param p_interpolant The object that shall hold interpolant data
> + * \param pf_points The set of data points (or data x-values) to use
> + * \param pf_values The set of data point values (or data y-values) to use
> + * \param i_point_count The number of points and values
> + * \return VLC_SUCCESS for success, VLC_ENOMEM for memory failure, or
> + * VLC_EBADVAR if there are not enough data points to create the
> spline
\retval
> + */
> +VLC_API int vlc_create_interpolant( vlc_interpolant_t * p_interpolant,
> + const float * pf_points,
> + const float * pf_values,
> + int i_point_count );
vlc_interpolant_create
> +
> +/**
> + * Sample the spline at the given sample point.
> + *
> + * \param interpolant The object that holds interpolant data
> + * \param f_sample_point The point at which the spline is to be sampled
> + * \return The value of the spline at the given sample point
> + */
> +VLC_API float vlc_interpolate( const vlc_interpolant_t interpolant,
> + float f_sample_point );
> +
> +/**
> + * Destroy an interpolant that is no longer in use.
> + *
> + * \param interpolant The object that holds interpolant data
> + */
> +VLC_API void vlc_destroy_interpolant( vlc_interpolant_t interpolant );
vlc_interpolant_destroy
> +
> +#endif
> diff --git a/po/POTFILES.in b/po/POTFILES.in
> index d3d3bb0..3a2dec2 100644
> --- a/po/POTFILES.in
> +++ b/po/POTFILES.in
> @@ -16,6 +16,7 @@ include/vlc_configuration.h
> include/vlc_demux.h
> include/vlc/deprecated.h
> include/vlc_epg.h
> +include/vlc_eqz_util.h
Missing file. This breaks incremental builds. Same below.
> include/vlc_es.h
> include/vlc_es_out.h
> include/vlc_events.h
> @@ -25,6 +26,7 @@ include/vlc_gcrypt.h
> include/vlc_httpd.h
> include/vlc_image.h
> include/vlc_input.h
> +include/vlc_interpolation.h
> include/vlc_intf_strings.h
> include/vlc_iso_lang.h
> include/vlc_keys.h
> @@ -105,11 +107,13 @@ src/libvlc.h
> src/libvlc-module.c
> src/misc/block.c
> src/misc/cpu.c
> +src/misc/eqz_util.c
> src/misc/error.c
> src/misc/es_format.c
> src/misc/events.c
> src/misc/filter_chain.c
> src/misc/image.c
> +src/misc/interpolation.c
> src/misc/md5.c
> src/misc/messages.c
> src/misc/mtime.c
> diff --git a/src/Makefile.am b/src/Makefile.am
> index 2520d67..b60257d 100644
> --- a/src/Makefile.am
> +++ b/src/Makefile.am
> @@ -42,6 +42,7 @@ pluginsinclude_HEADERS = \
> ../include/vlc_dialog.h \
> ../include/vlc_demux.h \
> ../include/vlc_epg.h \
> + ../include/vlc_eqz_util.h \
> ../include/vlc_es.h \
> ../include/vlc_es_out.h \
> ../include/vlc_events.h \
> @@ -57,6 +58,7 @@ pluginsinclude_HEADERS = \
> ../include/vlc_input.h \
> ../include/vlc_input_item.h \
> ../include/vlc_interface.h \
> + ../include/vlc_interpolation.h \
> ../include/vlc_keys.h \
> ../include/vlc_main.h \
> ../include/vlc_md5.h \
> @@ -434,7 +436,9 @@ SOURCES_libvlc_common = \
> misc/fifo.c \
> misc/fourcc.c \
> misc/fourcc_list.h \
> + misc/eqz_util.c \
> misc/es_format.c \
> + misc/interpolation.c \
> misc/picture.c \
> misc/picture.h \
> misc/picture_fifo.c \
> diff --git a/src/libvlccore.sym b/src/libvlccore.sym
> index a9e2636..b4faaa0 100644
> --- a/src/libvlccore.sym
> +++ b/src/libvlccore.sym
> @@ -526,8 +526,13 @@ vlc_sem_destroy
> vlc_sem_post
> vlc_sem_wait
> vlc_control_cancel
> +vlc_create_interpolant
> +vlc_destroy_interpolant
> vlc_GetCPUCount
> vlc_CPU
> +vlc_eqz_util_band_parser_destroy
> +vlc_eqz_util_band_parser_init
> +vlc_eqz_util_get_amp_array
> vlc_error
> vlc_event_attach
> vlc_event_detach
> @@ -554,6 +559,7 @@ vlc_ngettext
> vlc_iconv
> vlc_iconv_close
> vlc_iconv_open
> +vlc_interpolate
> vlc_poll_i11e
> vlc_read_i11e
> vlc_readv_i11e
> diff --git a/src/misc/interpolation.c b/src/misc/interpolation.c
> new file mode 100644
> index 0000000..4bd0a29
> --- /dev/null
> +++ b/src/misc/interpolation.c
> @@ -0,0 +1,212 @@
> +/**************************************************************************
> *** + * vlc_interpolation.h: Interpolation utility
> +
> ***************************************************************************
> ** + * Copyright (C) 2015 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. +
> ***************************************************************************
> **/ +
> +/**************************************************************************
> *** + * Preamble
> +
> ***************************************************************************
> **/ +
> +#ifdef HAVE_CONFIG_H
> +# include "config.h"
> +#endif
> +
> +#include <vlc_common.h>
> +#include <vlc_interpolation.h>
> +
> +/**************************************************************************
> *** + * Documentation : Read vlc_interpolation.h
> +
> ***************************************************************************
> **/ +
> +/**************************************************************************
> *** + * Private types.
> +
> ***************************************************************************
> **/ +
> +struct _vlc_interpolant_t
> +{
> + float * pf_points;
> + float * pf_values;
> +
> + float * pf_c1s;
> + float * pf_c2s;
> + float * pf_c3s;
> +
> + int i_point_count;
> +};
> +
> +/**************************************************************************
> *** + * Implementation of interpolation routines.
> +
> ***************************************************************************
> **/ +
> +int vlc_create_interpolant( vlc_interpolant_t * p_interpolant,
> + const float * pf_points, const float *
> pf_values,
> + int i_point_count )
> +{
> + /* This subroutine assumes that each sample point in the sequence
> + * referenced by pf_points is unique and that the sequence is sorted.
> */
> +
> + /* We should have at least two data points for the computation */
> + if( i_point_count < 2 )
> + return VLC_EBADVAR;
> +
> + int i_segment_count = i_point_count - 1;
> +
> + /* Allocate memory for pf_dys, pf_dxs, and pf_ms, which will all be
> arrays
> + * of length i_segment_count. */
> +
> + float * pf_dys = malloc( 3 * i_segment_count * sizeof( *pf_dys ) );
> + if( unlikely( !pf_dys ) )
> + return VLC_ENOMEM;
> +
> + float * pf_dxs = pf_dys + i_segment_count;
> + float * pf_ms = pf_dxs + i_segment_count;
> +
> + /* Allocate memory for the interpolant struct and its data members.
> pf_c1s,
> + * pf_points_, and pf_values_ will be arrays of length
> i_point_count, and
> + * pf_c2s and pf_c3s will be arrays of length
> i_segment_count. */
> +
> + vlc_interpolant_t interpolant;
> + float * pf_c1s, * pf_c2s, * pf_c3s, * pf_points_, * pf_values_;
> +
> + int i_ret_val = VLC_ENOMEM;
> + int i_data_len = 3 * i_point_count + 2 * i_segment_count;
> + int i_total_size = i_data_len * sizeof( *pf_c1s ) + sizeof(
*interpolant );
> + if( likely( ( interpolant = malloc( i_total_size ) )
> != NULL ) )
No need to merge assignment and condition here.
> + {
> + i_ret_val = VLC_SUCCESS;
> +
> + pf_c1s = ( float * ) ( interpolant + 1 );
Append a flexible array at the bottom of the structure if you want to do that.
Otherwise, the compiler won´t align the structure size.
> + pf_c2s = pf_c1s + i_point_count;
> + pf_c3s = pf_c2s + i_segment_count;
> + pf_points_ = pf_c3s + i_segment_count;
> + pf_values_ = pf_points_ + i_point_count;
> +
> + /* Get consecutive differences and slopes */
> + for( int i = 0; i < i_segment_count; i++ )
> + {
> + float f_dx = pf_points[i + 1] - pf_points[i];
> + float f_dy = pf_values[i + 1] - pf_values[i];
> + pf_dxs[i] = f_dx;
> + pf_dys[i] = f_dy;
> + pf_ms[i] = f_dy / f_dx;
> + }
> +
> + /* Get degree-1 coefficients */
> + pf_c1s[0] = pf_ms[0];
> + for( int i = 0; i < i_segment_count - 1; i++ )
> + {
> + float f_m = pf_ms[i];
> + float f_mNext = pf_ms[i + 1];
> + float f_c1;
> + if( f_m * f_mNext <= 0 )
> + {
> + f_c1 = 0.0f;
> + }
> + else
> + {
> + float f_dx = pf_dxs[i];
> + float f_dxNext = pf_dxs[i + 1];
> + float f_common = f_dx + f_dxNext;
> + f_c1 = 3.0f * f_common / ( ( f_common + f_dxNext ) / f_m
> + + ( f_common + f_dx ) / f_mNext );
> + }
> + pf_c1s[i + 1] = f_c1;
> + }
> + pf_c1s[i_segment_count] = pf_ms[i_segment_count - 1];
> +
> + /* Get degree-2 and degree-3 coefficients */
> + for( int i = 0; i < i_segment_count; i++ )
> + {
> + float f_c1 = pf_c1s[i];
> + float f_m = pf_ms[i];
> + float f_invDx = 1.0f / pf_dxs[i];
> + float f_common = f_c1 + pf_c1s[i + 1] - f_m - f_m;
> + pf_c2s[i] = ( f_m - f_c1 - f_common ) * f_invDx;
> + pf_c3s[i] = f_common * f_invDx * f_invDx;
> + }
> +
> + /* Set up structure */
> + memcpy( pf_points_, pf_points, i_point_count * sizeof( *pf_points )
> ); + memcpy( pf_values_, pf_values, i_point_count * sizeof(
> *pf_values ) );
> + interpolant->pf_c1s = pf_c1s;
> + interpolant->pf_c2s = pf_c2s;
> + interpolant->pf_c3s = pf_c3s;
> + interpolant->pf_points = pf_points_;
> + interpolant->pf_values = pf_values_;
> + interpolant->i_point_count = i_point_count;
> + *p_interpolant = interpolant;
> + }
> + free( pf_dys );
> + return i_ret_val;
> +}
> +
> +float vlc_interpolate( const vlc_interpolant_t interpolant,
> + float f_sample_point )
> +{
> + const float * pf_points = interpolant->pf_points;
> + const float * pf_values = interpolant->pf_values;
> +
> + /* Here, the implementation is modified to give the closest endpoints
> for
> + * sample points that lie outside the spline region */
> + int i = interpolant->i_point_count - 1;
> + if( f_sample_point >= pf_points[i] )
> + {
> + return pf_values[i];
> + }
> + if( f_sample_point <= pf_points[0] )
> + {
> + return pf_values[0];
> + }
> +
> + /* Use binary search to find the matching data value or the region
> where
> + * the matching value is to be interpolated */
> + int i_low = 0;
> + int i_high = i - 1;
> + while( i_low <= i_high )
> + {
> + int i_mid = ( i_low + i_high ) >> 1;
> + float f_point_here = pf_points[i_mid];
> + if( f_point_here < f_sample_point )
> + {
> + i_low = i_mid + 1;
> + }
> + else if( f_point_here > f_sample_point )
> + {
> + i_high = i_mid - 1;
> + }
> + else
> + {
> + return pf_values[i_mid];
> + }
> + }
> + i = i_high > 0 ? i_high : 0;
> +
> + /* Interpolate */
> + float f_diff = f_sample_point - pf_points[i];
> + return pf_values[i] + f_diff * ( interpolant->pf_c1s[i]
> + + f_diff * ( interpolant->pf_c2s[i]
> + + f_diff * interpolant->pf_c3s[i] ) );
> +}
> +
> +void vlc_destroy_interpolant( vlc_interpolant_t interpolant )
> +{
> + free( interpolant );
> +}
--
Rémi Denis-Courmont
http://www.remlab.net/
More information about the vlc-devel
mailing list