src/spectral/phasevoc.c

   1 /*
   2   Copyright (C) 2003-2014 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 #include "aubio_priv.h"
  22 #include "fvec.h"
  23 #include "cvec.h"
  24 #include "mathutils.h"
  25 #include "spectral/fft.h"
  26 #include "spectral/phasevoc.h"
  27
  28 /** phasevocoder internal object */
  29 struct _aubio_pvoc_t {
  30   uint_t win_s;       /** grain length */
  31   uint_t hop_s;       /** overlap step */
  32   aubio_fft_t * fft;  /** fft object */
  33   fvec_t * data;      /** current input grain, [win_s] frames */
  34   fvec_t * dataold;   /** memory of past grain, [win_s-hop_s] frames */
  35   fvec_t * synth;     /** current output grain, [win_s] frames */
  36   fvec_t * synthold;  /** memory of past grain, [win_s-hop_s] frames */
  37   fvec_t * w;         /** grain window [win_s] */
  38   uint_t start;       /** where to start additive synthesis */
  39   uint_t end;         /** where to end it */
  40   smpl_t scale;       /** scaling factor for synthesis */
  41   uint_t end_datasize;  /** size of memory to end */
  42   uint_t hop_datasize;  /** size of memory to hop_s */
  43 };
  44
  45
  46 /** returns data and dataold slided by hop_s */
  47 static void aubio_pvoc_swapbuffers(aubio_pvoc_t *pv, const fvec_t *new);
  48
  49 /** do additive synthesis from 'old' and 'cur' */
  50 static void aubio_pvoc_addsynth(aubio_pvoc_t *pv, fvec_t * synthnew);
  51
  52 void aubio_pvoc_do(aubio_pvoc_t *pv, const fvec_t * datanew, cvec_t *fftgrain) {
  53   /* slide  */
  54   aubio_pvoc_swapbuffers(pv, datanew);
  55   /* windowing */
  56   fvec_weight(pv->data, pv->w);
  57   /* shift */
  58   fvec_shift(pv->data);
  59   /* calculate fft */
  60   aubio_fft_do (pv->fft,pv->data,fftgrain);
  61 }
  62
  63 void aubio_pvoc_rdo(aubio_pvoc_t *pv,cvec_t * fftgrain, fvec_t * synthnew) {
  64   /* calculate rfft */
  65   aubio_fft_rdo(pv->fft,fftgrain,pv->synth);
  66   /* unshift */
  67   fvec_ishift(pv->synth);
  68   /* windowing */
  69   // if overlap = 50%, do not apply window (identity)
  70   if (pv->hop_s * 2 < pv->win_s) {
  71     fvec_weight(pv->synth, pv->w);
  72   }
  73   /* additive synthesis */
  74   aubio_pvoc_addsynth(pv, synthnew);
  75 }
  76
  77 aubio_pvoc_t * new_aubio_pvoc (uint_t win_s, uint_t hop_s) {
  78   aubio_pvoc_t * pv = AUBIO_NEW(aubio_pvoc_t);
  79
  80   /* if (win_s < 2*hop_s) {
  81     AUBIO_WRN("Hop size bigger than half the window size!\n");
  82   } */
  83
  84   if ((sint_t)hop_s < 1) {
  85     AUBIO_ERR("pvoc: got hop_size %d, but can not be < 1\n", hop_s);
  86     goto beach;
  87   } else if ((sint_t)win_s < 2) {
  88     AUBIO_ERR("pvoc: got buffer_size %d, but can not be < 2\n", win_s);
  89     goto beach;
  90   } else if (win_s < hop_s) {
  91     AUBIO_ERR("pvoc: hop size (%d) is larger than win size (%d)\n", win_s, hop_s);
  92     goto beach;
  93   }
  94
  95   pv->fft      = new_aubio_fft (win_s);
  96   if (pv->fft == NULL) {
  97     goto beach;
  98   }
  99
 100   /* remember old */
 101   pv->data     = new_fvec (win_s);
 102   pv->synth    = new_fvec (win_s);
 103
 104   /* new input output */
 105   if (win_s > hop_s) {
 106     pv->dataold  = new_fvec  (win_s-hop_s);
 107     pv->synthold = new_fvec (win_s-hop_s);
 108   } else {
 109     pv->dataold  = new_fvec  (1);
 110     pv->synthold = new_fvec (1);
 111   }
 112   pv->w        = new_aubio_window ("hanningz", win_s);
 113
 114   pv->hop_s    = hop_s;
 115   pv->win_s    = win_s;
 116
 117   /* more than 50% overlap, overlap anyway */
 118   if (win_s < 2 * hop_s) pv->start = 0;
 119   /* less than 50% overlap, reset latest grain trail */
 120   else pv->start = win_s - hop_s - hop_s;
 121
 122   if (win_s > hop_s) pv->end = win_s - hop_s;
 123   else pv->end = 0;
 124
 125   pv->end_datasize = pv->end * sizeof(smpl_t);
 126   pv->hop_datasize = pv->hop_s * sizeof(smpl_t);
 127
 128   // for reconstruction with 75% overlap
 129   if (win_s == hop_s * 4) {
 130     pv->scale = 2./3.;
 131   } else if (win_s == hop_s * 8) {
 132     pv->scale = 1./3.;
 133   } else if (win_s == hop_s * 2) {
 134     pv->scale = 1.;
 135   } else {
 136     pv->scale = .5;
 137   }
 138
 139   return pv;
 140
 141 beach:
 142   AUBIO_FREE (pv);
 143   return NULL;
 144 }
 145
 146 void del_aubio_pvoc(aubio_pvoc_t *pv) {
 147   del_fvec(pv->data);
 148   del_fvec(pv->synth);
 149   del_fvec(pv->dataold);
 150   del_fvec(pv->synthold);
 151   del_fvec(pv->w);
 152   del_aubio_fft(pv->fft);
 153   AUBIO_FREE(pv);
 154 }
 155
 156 static void aubio_pvoc_swapbuffers(aubio_pvoc_t *pv, const fvec_t *new)
 157 {
 158   /* some convenience pointers */
 159   smpl_t * data = pv->data->data;
 160   smpl_t * dataold = pv->dataold->data;
 161   smpl_t * datanew = new->data;
 162 #ifndef HAVE_MEMCPY_HACKS
 163   uint_t i;
 164   for (i = 0; i < pv->end; i++)
 165     data[i] = dataold[i];
 166   for (i = 0; i < pv->hop_s; i++)
 167     data[pv->end + i] = datanew[i];
 168   for (i = 0; i < pv->end; i++)
 169     dataold[i] = data[i + pv->hop_s];
 170 #else
 171   memcpy(data, dataold, pv->end_datasize);
 172   data += pv->end;
 173   memcpy(data, datanew, pv->hop_datasize);
 174   data -= pv->end;
 175   data += pv->hop_s;
 176   memcpy(dataold, data, pv->end_datasize);
 177 #endif
 178 }
 179
 180 static void aubio_pvoc_addsynth(aubio_pvoc_t *pv, fvec_t *synth_new)
 181 {
 182   uint_t i;
 183   /* some convenience pointers */
 184   smpl_t * synth    = pv->synth->data;
 185   smpl_t * synthold = pv->synthold->data;
 186   smpl_t * synthnew = synth_new->data;
 187
 188   /* put new result in synthnew */
 189   for (i = 0; i < pv->hop_s; i++)
 190     synthnew[i] = synth[i] * pv->scale;
 191
 192   /* no overlap, nothing else to do */
 193   if (pv->end == 0) return;
 194
 195   /* add new synth to old one */
 196   for (i = 0; i < pv->hop_s; i++)
 197     synthnew[i] += synthold[i];
 198
 199   /* shift synthold */
 200   for (i = 0; i < pv->start; i++)
 201     synthold[i] = synthold[i + pv->hop_s];
 202
 203   /* erase last frame in synthold */
 204   for (i = pv->start; i < pv->end; i++)
 205     synthold[i] = 0.;
 206
 207   /* additive synth */
 208   for (i = 0; i < pv->end; i++)
 209     synthold[i] += synth[i + pv->hop_s] * pv->scale;
 210 }