src/mathutils.h

   1 /*
   2   Copyright (C) 2003-2015 Paul Brossier <piem@aubio.org>
   3
   4   This file is part of aubio.
   5
   6   aubio is free software: you can redistribute it and/or modify
   7   it under the terms of the GNU General Public License as published by
   8   the Free Software Foundation, either version 3 of the License, or
   9   (at your option) any later version.
  10
  11   aubio is distributed in the hope that it will be useful,
  12   but WITHOUT ANY WARRANTY; without even the implied warranty of
  13   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  14   GNU General Public License for more details.
  15
  16   You should have received a copy of the GNU General Public License
  17   along with aubio.  If not, see <http://www.gnu.org/licenses/>.
  18
  19 */
  20
  21 /** \file
  22
  23   Various math functions
  24
  25   \example test-mathutils.c
  26   \example test-mathutils-window.c
  27
  28  */
  29
  30 #ifndef AUBIO_MATHUTILS_H
  31 #define AUBIO_MATHUTILS_H
  32
  33 #include "fvec.h"
  34 #include "musicutils.h"
  35
  36 #ifdef __cplusplus
  37 extern "C" {
  38 #endif
  39
  40 /** compute the mean of a vector
  41
  42   \param s vector to compute mean from
  43   \return the mean of `v`
  44
  45 */
  46 smpl_t fvec_mean (fvec_t * s);
  47
  48 /** find the max of a vector
  49
  50   \param s vector to get the max from
  51
  52   \return the value of the minimum of v
  53
  54 */
  55 smpl_t fvec_max (fvec_t * s);
  56
  57 /** find the min of a vector
  58
  59   \param s vector to get the min from
  60
  61   \return the value of the maximum of v
  62
  63 */
  64 smpl_t fvec_min (fvec_t * s);
  65
  66 /** find the index of the min of a vector
  67
  68   \param s vector to get the index from
  69
  70   \return the index of the minimum element of v
  71
  72 */
  73 uint_t fvec_min_elem (fvec_t * s);
  74
  75 /** find the index of the max of a vector
  76
  77   \param s vector to get the index from
  78
  79   \return the index of the maximum element of v
  80
  81 */
  82 uint_t fvec_max_elem (fvec_t * s);
  83
  84 /** swap the left and right halves of a vector
  85
  86   This function swaps the left part of the signal with the right part of the
  87 signal. Therefore
  88
  89   \f$ a[0], a[1], ..., a[\frac{N}{2}], a[\frac{N}{2}+1], ..., a[N-1], a[N] \f$
  90
  91   becomes
  92
  93   \f$ a[\frac{N}{2}+1], ..., a[N-1], a[N], a[0], a[1], ..., a[\frac{N}{2}] \f$
  94
  95   This operation, known as 'fftshift' in the Matlab Signal Processing Toolbox,
  96 can be used before computing the FFT to simplify the phase relationship of the
  97 resulting spectrum. See Amalia de Götzen's paper referred to above.
  98
  99 */
 100 void fvec_shift (fvec_t * v);
 101
 102 /** swap the left and right halves of a vector
 103
 104   This function swaps the left part of the signal with the right part of the
 105 signal. Therefore
 106
 107   \f$ a[0], a[1], ..., a[\frac{N}{2}], a[\frac{N}{2}+1], ..., a[N-1], a[N] \f$
 108
 109   becomes
 110
 111   \f$ a[\frac{N}{2}+1], ..., a[N-1], a[N], a[0], a[1], ..., a[\frac{N}{2}] \f$
 112
 113   This operation, known as 'ifftshift' in the Matlab Signal Processing Toolbox,
 114 can be used after computing the inverse FFT to simplify the phase relationship
 115 of the resulting spectrum. See Amalia de Götzen's paper referred to above.
 116
 117 */
 118 void fvec_ishift (fvec_t * v);
 119
 120 /** compute the sum of all elements of a vector
 121
 122   \param v vector to compute the sum of
 123
 124   \return the sum of v
 125
 126 */
 127 smpl_t fvec_sum (fvec_t * v);
 128
 129 /** compute the High Frequency Content of a vector
 130
 131   The High Frequency Content is defined as \f$ \sum_0^{N-1} (k+1) v[k] \f$.
 132
 133   \param v vector to get the energy from
 134
 135   \return the HFC of v
 136
 137 */
 138 smpl_t fvec_local_hfc (fvec_t * v);
 139
 140 /** computes the p-norm of a vector
 141
 142   Computes the p-norm of a vector for \f$ p = \alpha \f$
 143
 144   \f$ L^p = ||x||_p = (|x_1|^p + |x_2|^p + ... + |x_n|^p ) ^ \frac{1}{p} \f$
 145
 146   If p = 1, the result is the Manhattan distance.
 147
 148   If p = 2, the result is the Euclidean distance.
 149
 150   As p tends towards large values, \f$ L^p \f$ tends towards the maximum of the
 151 input vector.
 152
 153   References:
 154
 155     - <a href="http://en.wikipedia.org/wiki/Lp_space">\f$L^p\f$ space</a> on
 156   Wikipedia
 157
 158   \param v vector to compute norm from
 159   \param p order of the computed norm
 160
 161   \return the p-norm of v
 162
 163 */
 164 smpl_t fvec_alpha_norm (fvec_t * v, smpl_t p);
 165
 166 /**  alpha normalisation
 167
 168   This function divides all elements of a vector by the p-norm as computed by
 169 fvec_alpha_norm().
 170
 171   \param v vector to compute norm from
 172   \param p order of the computed norm
 173
 174 */
 175 void fvec_alpha_normalise (fvec_t * v, smpl_t p);
 176
 177 /** add a constant to each elements of a vector
 178
 179   \param v vector to add constant to
 180   \param c constant to add to v
 181
 182 */
 183 void fvec_add (fvec_t * v, smpl_t c);
 184
 185 /** remove the minimum value of the vector to each elements
 186
 187   \param v vector to remove minimum from
 188
 189 */
 190 void fvec_min_removal (fvec_t * v);
 191
 192 /** compute moving median threshold of a vector
 193
 194   This function computes the moving median threshold value of at the given
 195 position of a vector, taking the median among post elements before and up to
 196 pre elements after pos.
 197
 198   \param v input vector
 199   \param tmp temporary vector of length post+1+pre
 200   \param post length of causal part to take before pos
 201   \param pre length of anti-causal part to take after pos
 202   \param pos index to compute threshold for
 203
 204   \return moving median threshold value
 205
 206 */
 207 smpl_t fvec_moving_thres (fvec_t * v, fvec_t * tmp, uint_t post, uint_t pre,
 208     uint_t pos);
 209
 210 /** apply adaptive threshold to a vector
 211
 212   For each points at position p of an input vector, this function remove the
 213 moving median threshold computed at p.
 214
 215   \param v input vector
 216   \param tmp temporary vector of length post+1+pre
 217   \param post length of causal part to take before pos
 218   \param pre length of anti-causal part to take after pos
 219
 220 */
 221 void fvec_adapt_thres (fvec_t * v, fvec_t * tmp, uint_t post, uint_t pre);
 222
 223 /** returns the median of a vector
 224
 225   The QuickSelect routine is based on the algorithm described in "Numerical
 226 recipes in C", Second Edition, Cambridge University Press, 1992, Section 8.5,
 227 ISBN 0-521-43108-5
 228
 229   This implementation of the QuickSelect routine is based on Nicolas
 230 Devillard's implementation, available at http://ndevilla.free.fr/median/median/
 231 and in the Public Domain.
 232
 233   \param v vector to get median from
 234
 235   \return the median of v
 236
 237 */
 238 smpl_t fvec_median (fvec_t * v);
 239
 240 /** finds exact peak index by quadratic interpolation
 241
 242   See [Quadratic Interpolation of Spectral
 243   Peaks](https://ccrma.stanford.edu/~jos/sasp/Quadratic_Peak_Interpolation.html),
 244   by Julius O. Smith III
 245
 246   \f$ p_{frac} = \frac{1}{2} \frac {x[p-1] - x[p+1]} {x[p-1] - 2 x[p] + x[p+1]} \in [ -.5, .5] \f$
 247
 248   \param x vector to get the interpolated peak position from
 249   \param p index of the peak in vector `x`
 250   \return \f$ p + p_{frac} \f$ exact peak position of interpolated maximum or minimum
 251
 252 */
 253 smpl_t fvec_quadratic_peak_pos (const fvec_t * x, uint_t p);
 254
 255 /** finds magnitude of peak by quadratic interpolation
 256
 257   See [Quadratic Interpolation of Spectral
 258   Peaks](https://ccrma.stanford.edu/~jos/sasp/Quadratic_Peak_Interpolation.html),
 259   by Julius O. Smith III
 260
 261   \param x vector to get the magnitude of the interpolated peak position from
 262   \param p index of the peak in vector `x`
 263   \return magnitude of interpolated peak
 264
 265 */
 266 smpl_t fvec_quadratic_peak_mag (fvec_t * x, smpl_t p);
 267
 268 /** Quadratic interpolation using Lagrange polynomial.
 269
 270   Inspired from ``Comparison of interpolation algorithms in real-time sound
 271 processing'', Vladimir Arnost,
 272
 273   \param s0,s1,s2 are 3 consecutive samples of a curve
 274   \param pf is the floating point index [0;2]
 275
 276   \return \f$ s0 + (pf/2.)*((pf-3.)*s0-2.*(pf-2.)*s1+(pf-1.)*s2); \f$
 277
 278 */
 279 smpl_t aubio_quadfrac (smpl_t s0, smpl_t s1, smpl_t s2, smpl_t pf);
 280
 281 /** return 1 if v[p] is a peak and positive, 0 otherwise
 282
 283   This function returns 1 if a peak is found at index p in the vector v. The
 284 peak is defined as follows:
 285
 286   - v[p] is positive
 287   - v[p-1] < v[p]
 288   - v[p] > v[p+1]
 289
 290   \param v input vector
 291   \param p position of supposed for peak
 292
 293   \return 1 if a peak is found, 0 otherwise
 294
 295 */
 296 uint_t fvec_peakpick (const fvec_t * v, uint_t p);
 297
 298 /** return 1 if a is a power of 2, 0 otherwise */
 299 uint_t aubio_is_power_of_two(uint_t a);
 300
 301 /** return the next power of power of 2 greater than a */
 302 uint_t aubio_next_power_of_two(uint_t a);
 303
 304 /** compute normalised autocorrelation function
 305
 306   \param input vector to compute autocorrelation from
 307   \param output vector to store autocorrelation function to
 308
 309 */
 310 void aubio_autocorr (const fvec_t * input, fvec_t * output);
 311
 312 #ifdef __cplusplus
 313 }
 314 #endif
 315
 316 #endif /* AUBIO_MATHUTILS_H */