ref: 4daa7ec629cb0273e68aff1fcd4eb9fc9a52e93a
parent: b236757ac4ae678d7f762e1d7cabb18ca22486fc
author: Paul Brossier <piem@piem.org>
date: Wed Aug 8 13:23:57 EDT 2018
src/spectral/dct.c: add dct factory
--- /dev/null
+++ b/src/spectral/dct.c
@@ -1,0 +1,169 @@
+/*
+ Copyright (C) 2018 Paul Brossier <piem@aubio.org>
+
+ This file is part of aubio.
+
+ aubio is free software: you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation, either version 3 of the License, or
+ (at your option) any later version.
+
+ aubio 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 General Public License for more details.
+
+ You should have received a copy of the GNU General Public License
+ along with aubio. If not, see <http://www.gnu.org/licenses/>.
+
+*/
+
+/** \file
+
+ Discrete Cosine Transform
+
+ Functions aubio_dct_do() and aubio_dct_rdo() are equivalent to MATLAB/Octave
+ dct() and idct() functions, as well as scipy.fftpack.dct(x, norm='ortho') and
+ scipy.fftpack.idct(x, norm='ortho')
+
+ \example spectral/test-dct.c
+
+*/
+
+#include "aubio_priv.h"
+#include "fvec.h"
+#include "spectral/dct.h"
+
+// function pointers prototypes
+typedef void (*aubio_dct_do_t)(aubio_dct_t * s, const fvec_t * input, fvec_t * output);
+typedef void (*aubio_dct_rdo_t)(aubio_dct_t * s, const fvec_t * input, fvec_t * output);
+typedef void (*del_aubio_dct_t)(aubio_dct_t * s);
+
+#if defined(HAVE_ACCELERATE)
+typedef struct _aubio_dct_opt_t aubio_dct_accelerate_t;
+extern aubio_dct_accelerate_t * new_aubio_dct_accelerate (uint_t size);
+extern void aubio_dct_accelerate_do(aubio_dct_accelerate_t *s, const fvec_t *input, fvec_t *output);
+extern void aubio_dct_accelerate_rdo(aubio_dct_accelerate_t *s, const fvec_t *input, fvec_t *output);
+extern void del_aubio_dct_accelerate (aubio_dct_accelerate_t *s);
+#elif defined(HAVE_FFTW3)
+typedef struct _aubio_dct_opt_t aubio_dct_fftw_t;
+extern aubio_dct_fftw_t * new_aubio_dct_fftw (uint_t size);
+extern void aubio_dct_fftw_do(aubio_dct_fftw_t *s, const fvec_t *input, fvec_t *output);
+extern void aubio_dct_fftw_rdo(aubio_dct_fftw_t *s, const fvec_t *input, fvec_t *output);
+extern void del_aubio_dct_fftw (aubio_dct_fftw_t *s);
+#elif defined(HAVE_INTEL_IPP)
+typedef struct _aubio_dct_ipp_t aubio_dct_ipp_t;
+extern aubio_dct_ipp_t * new_aubio_dct_ipp (uint_t size);
+extern void aubio_dct_ipp_do(aubio_dct_ipp_t *s, const fvec_t *input, fvec_t *output);
+extern void aubio_dct_ipp_rdo(aubio_dct_ipp_t *s, const fvec_t *input, fvec_t *output);
+extern void del_aubio_dct_ipp (aubio_dct_ipp_t *s);
+#else
+typedef struct _aubio_dct_ooura_t aubio_dct_ooura_t;
+extern aubio_dct_ooura_t * new_aubio_dct_ooura (uint_t size);
+extern void aubio_dct_ooura_do(aubio_dct_ooura_t *s, const fvec_t *input, fvec_t *output);
+extern void aubio_dct_ooura_rdo(aubio_dct_ooura_t *s, const fvec_t *input, fvec_t *output);
+extern void del_aubio_dct_ooura (aubio_dct_ooura_t *s);
+#endif
+
+// plain mode
+typedef struct _aubio_dct_plain_t aubio_dct_plain_t;
+extern aubio_dct_plain_t * new_aubio_dct_plain (uint_t size);
+extern void aubio_dct_plain_do(aubio_dct_plain_t *s, const fvec_t *input, fvec_t *output);
+extern void aubio_dct_plain_rdo(aubio_dct_plain_t *s, const fvec_t *input, fvec_t *output);
+extern void del_aubio_dct_plain (aubio_dct_plain_t *s);
+
+struct _aubio_dct_t {+ void *dct;
+ aubio_dct_do_t dct_do;
+ aubio_dct_rdo_t dct_rdo;
+ del_aubio_dct_t del_dct;
+};
+
+aubio_dct_t* new_aubio_dct (uint_t size) {+ aubio_dct_t * s = AUBIO_NEW(aubio_dct_t);
+ if ((sint_t)size <= 0) goto beach;
+#if defined(HAVE_ACCELERATE)
+ // TODO check
+ // vDSP supports sizes = f * 2 ** n, where n >= 4 and f in [1, 3, 5, 15]
+ // see https://developer.apple.com/documentation/accelerate/1449930-vdsp_dct_createsetup
+ if (aubio_is_power_of_two(size/16) != 1
+ || (size/16 != 3 && size/16 != 5 && size/16 != 15)) {+ goto plain;
+ }
+ s->dct = (void *)new_aubio_dct_accelerate (size);
+ if (s->dct) {+ s->dct_do = (aubio_dct_do_t)aubio_dct_accelerate_do;
+ s->dct_rdo = (aubio_dct_rdo_t)aubio_dct_accelerate_rdo;
+ s->del_dct = (del_aubio_dct_t)del_aubio_dct_accelerate;
+ return s;
+ }
+#elif defined(HAVE_FFTW3)
+ // fftw supports any positive integer size
+ s->dct = (void *)new_aubio_dct_fftw (size);
+ if (s->dct) {+ s->dct_do = (aubio_dct_do_t)aubio_dct_fftw_do;
+ s->dct_rdo = (aubio_dct_rdo_t)aubio_dct_fftw_rdo;
+ s->del_dct = (del_aubio_dct_t)del_aubio_dct_fftw;
+ return s;
+ } else {+ AUBIO_WRN("dct: unexcepected error while creating dct_fftw with size %d",+ size);
+ goto plain;
+ }
+#elif defined(HAVE_INTEL_IPP)
+ // unclear from the docs, but intel ipp seems to support any size
+ s->dct = (void *)new_aubio_dct_ipp (size);
+ if (s->dct) {+ s->dct_do = (aubio_dct_do_t)aubio_dct_ipp_do;
+ s->dct_rdo = (aubio_dct_rdo_t)aubio_dct_ipp_rdo;
+ s->del_dct = (del_aubio_dct_t)del_aubio_dct_ipp;
+ return s;
+ } else {+ AUBIO_WRN("dct: unexcepected error while creating dct_ipp with size %d",+ size);
+ goto plain;
+ }
+#else
+ // ooura does not support sizes < 16 or non-power of 2
+ if (aubio_is_power_of_two(size) != 1 || size < 16) {+ goto plain;
+ }
+ s->dct = (void *)new_aubio_dct_ooura (size);
+ if (s->dct) {+ s->dct_do = (aubio_dct_do_t)aubio_dct_ooura_do;
+ s->dct_rdo = (aubio_dct_rdo_t)aubio_dct_ooura_rdo;
+ s->del_dct = (del_aubio_dct_t)del_aubio_dct_ooura;
+ return s;
+ }
+#endif
+ // falling back to plain mode
+ AUBIO_WRN("dct: d no optimised implementation could be created for size %d",+ size);
+plain:
+ s->dct = (void *)new_aubio_dct_plain (size);
+ if (s->dct) {+ s->dct_do = (aubio_dct_do_t)aubio_dct_plain_do;
+ s->dct_rdo = (aubio_dct_rdo_t)aubio_dct_plain_rdo;
+ s->del_dct = (del_aubio_dct_t)del_aubio_dct_plain;
+ return s;
+ } else {+ goto beach;
+ }
+beach:
+ AUBIO_ERROR("dct: failed creating with size %d, should be > 0", size);+ AUBIO_FREE (s);
+ return NULL;
+}
+
+void del_aubio_dct(aubio_dct_t *s) {+ if (s->dct && s->del_dct) s->del_dct (s->dct);
+ AUBIO_FREE (s);
+}
+
+void aubio_dct_do(aubio_dct_t *s, const fvec_t *input, fvec_t *output) {+ s->dct_do ((void *)s->dct, input, output);
+}
+
+void aubio_dct_rdo(aubio_dct_t *s, const fvec_t *input, fvec_t *output) {+ s->dct_rdo ((void *)s->dct, input, output);
+}