add yinfft method, make it default for pd plugin
authorPaul Brossier <piem@altern.org>
Tue, 21 Mar 2006 22:59:17 +0000 (22:59 +0000)
committerPaul Brossier <piem@altern.org>
Tue, 21 Mar 2006 22:59:17 +0000 (22:59 +0000)
add yinfft method, make it default for pd plugin

plugins/puredata/aubiopitch~.c
python/aubio/task/pitch.py
python/aubio/task/utils.py
python/aubiopitch
src/Makefile.am
src/pitchdetection.c
src/pitchdetection.h
src/pitchyinfft.c [new file with mode: 0644]
src/pitchyinfft.h [new file with mode: 0644]
swig/aubio.i

index 3a088cb..ffa4c13 100644 (file)
@@ -12,7 +12,7 @@
 
 char aubiopitch_version[] = "aubiopitch~ version 0.1";
 
-aubio_pitchdetection_type type_pitch = aubio_pitch_mcomb;
+aubio_pitchdetection_type type_pitch = aubio_pitch_yinfft;
 aubio_pitchdetection_mode mode_pitch = aubio_pitchm_freq;
 
 static t_class *aubiopitch_tilde_class;
@@ -80,6 +80,7 @@ static void *aubiopitch_tilde_new (void)
        //FIXME: get the real samplerate
        x->o = new_aubio_pitchdetection(x->bufsize*4, 
                     x->hopsize, 1, 44100., type_pitch, mode_pitch);
+       aubio_pitchdetection_set_yinthresh(x->o, 0.7);
        x->vec = (fvec_t *)new_fvec(x->hopsize,1);
 
        //floatinlet_new (&x->x_obj, &x->threshold);
index a4b2003..f2e8312 100644 (file)
@@ -7,13 +7,19 @@ class taskpitch(task):
        def __init__(self,input,params=None):
                task.__init__(self,input,params=params)
                self.shortlist = [0. for i in range(self.params.pitchsmooth)]
+               if self.params.pitchmode == 'yinfft':
+                       yinthresh = self.params.yinfftthresh
+               elif self.params.pitchmode == 'yin':
+                       yinthresh = self.params.yinthresh
+               else:
+                       yinthresh = 0.
                self.pitchdet  = pitchdetection(mode=get_pitch_mode(self.params.pitchmode),
                        bufsize=self.params.bufsize,
                        hopsize=self.params.hopsize,
                        channels=self.channels,
                        samplerate=self.srate,
                        omode=self.params.omode,
-                       yinthresh=self.params.yinthresh)
+                       yinthresh=yinthresh)
 
        def __call__(self):
                from aubio.median import short_find
@@ -67,7 +73,7 @@ class taskpitch(task):
                                time,pitch =[],[]
                                while(tasksil.readsize==tasksil.params.hopsize):
                                        tasksil()
-                                       time.append(tasksil.params.step*tasksil.frameread)
+                                       time.append(tasksil.params.step*(tasksil.frameread))
                                        if not tasksil.issilence:
                                                pitch.append(self.truth)
                                        else:
@@ -155,7 +161,7 @@ class taskpitch(task):
                        title=self.params.pitchmode))
 
                        
-       def plotplot(self,wplot,oplots,outplot=None,multiplot = 1, midi = 1):
+       def plotplot(self,wplot,oplots,outplot=None,multiplot = 0, midi = 1):
                from aubio.gnuplot import gnuplot_init, audio_to_array, make_audio_plot
                import re
                import Gnuplot
@@ -197,7 +203,7 @@ class taskpitch(task):
                        g('set yrange [100:%f]' % self.params.pitchmax) 
                else: 
                        g.ylabel('pitch (midi)')
-                       g('set yrange [%f:%f]' % (40, aubio_freqtomidi(self.params.pitchmax)))
+                       g('set yrange [%f:%f]' % (aubio_freqtomidi(self.params.pitchmin), aubio_freqtomidi(self.params.pitchmax)))
                g('set key right top')
                g('set noclip one') 
                g('set format x ""')
index 929bf5e..ef49718 100644 (file)
@@ -33,6 +33,8 @@ def get_pitch_mode(nvalue):
                 return aubio_pitch_fcomb
        elif nvalue == 'schmitt':
                 return aubio_pitch_schmitt
+       elif nvalue == 'yinfft':
+                return aubio_pitch_yinfft
        else:
                 import sys
                 print "error: unknown pitch detection function selected"
index 07cdcf7..7d3cbec 100755 (executable)
@@ -89,7 +89,7 @@ params.samplerate = float(sndfile(filename).samplerate())
 params.hopsize    = int(options.hopsize)
 params.bufsize    = int(options.bufsize)
 params.step       = params.samplerate/float(params.hopsize)
-params.threshold  = float(options.threshold)
+params.yinthresh  = float(options.threshold)
 params.silence    = float(options.silence)
 params.verbose    = options.verbose
 #mintol     = float(options.mintol)*step
@@ -100,7 +100,6 @@ else:             delay = 2./params.step
 if options.note:
         exit("not implemented yet")
 
-
 wplot,oplots = [],[]
 modes = options.mode.split(',')
 for i in range(len(modes)):
index 75e3b1c..1b68cbf 100644 (file)
@@ -17,6 +17,7 @@ pkginclude_HEADERS = aubio.h \
        pitchyin.h \
        pitchschmitt.h \
        pitchfcomb.h \
+       pitchyinfft.h \
        beattracking.h \
        filter.h
 nodist_pkginclude_HEADERS = config.h
@@ -56,6 +57,8 @@ libaubio_la_SOURCES = aubio.h \
        pitchschmitt.h \
        pitchfcomb.c \
        pitchfcomb.h \
+       pitchyinfft.c \
+       pitchyinfft.h \
        beattracking.c \
        beattracking.h \
        filter.c \
index 10d1224..29cd24d 100644 (file)
@@ -25,6 +25,7 @@
 #include "pitchyin.h"
 #include "pitchfcomb.h"
 #include "pitchschmitt.h"
+#include "pitchyinfft.h"
 #include "pitchdetection.h"
 
 typedef smpl_t (*aubio_pitchdetection_func_t)(aubio_pitchdetection_t *p, 
@@ -43,6 +44,7 @@ struct _aubio_pitchdetection_t {
        aubio_pitchmcomb_t * mcomb;
        aubio_pitchfcomb_t * fcomb;
        aubio_pitchschmitt_t * schmitt;
+       aubio_pitchyinfft_t * yinfft;
        aubio_filter_t * filter;
        /* for yin */
        fvec_t * buf;
@@ -91,7 +93,7 @@ aubio_pitchdetection_t * new_aubio_pitchdetection(uint_t bufsize,
                case aubio_pitch_mcomb:
                        p->pv       = new_aubio_pvoc(bufsize, hopsize, channels);
                        p->fftgrain = new_cvec(bufsize, channels);
-                       p->mcomb    = new_aubio_pitchmcomb(bufsize,channels,samplerate);
+                       p->mcomb    = new_aubio_pitchmcomb(bufsize,hopsize,channels,samplerate);
                        p->filter   = new_aubio_cdsgn_filter(samplerate);
                         p->callback = aubio_pitchdetection_mcomb;
                        break;
@@ -105,6 +107,12 @@ aubio_pitchdetection_t * new_aubio_pitchdetection(uint_t bufsize,
                         p->schmitt  = new_aubio_pitchschmitt(bufsize,samplerate);
                         p->callback = aubio_pitchdetection_schmitt;
                         break;
+               case aubio_pitch_yinfft:
+                       p->buf      = new_fvec(bufsize,channels);
+                        p->yinfft   = new_aubio_pitchyinfft(bufsize);
+                        p->callback = aubio_pitchdetection_yinfft;
+                       p->yinthres = 0.2;
+                        break;
                 default:
                         break;
        }
@@ -147,6 +155,10 @@ void del_aubio_pitchdetection(aubio_pitchdetection_t * p) {
                        del_fvec(p->buf);
                         del_aubio_pitchfcomb(p->fcomb);
                         break;
+                case aubio_pitch_yinfft:
+                       del_fvec(p->buf);
+                        del_aubio_pitchyinfft(p->yinfft);
+                        break;
                default:
                        break;
        }
@@ -205,6 +217,18 @@ smpl_t aubio_pitchdetection_yin(aubio_pitchdetection_t *p, fvec_t *ibuf) {
 }
 
 
+smpl_t aubio_pitchdetection_yinfft(aubio_pitchdetection_t *p, fvec_t *ibuf){
+       smpl_t pitch = 0.;
+       aubio_pitchdetection_slideblock(p,ibuf);
+       pitch = aubio_pitchyinfft_detect(p->yinfft,p->buf,p->yinthres);
+        if (pitch>0) {
+                pitch = p->srate/(pitch+0.);
+        } else {
+                pitch = 0.;
+        }
+       return pitch; 
+}
+
 smpl_t aubio_pitchdetection_fcomb(aubio_pitchdetection_t *p, fvec_t *ibuf){
         aubio_pitchdetection_slideblock(p,ibuf);
         return aubio_pitchfcomb_detect(p->fcomb,p->buf);
index 1238f19..d640345 100644 (file)
@@ -27,7 +27,8 @@ typedef enum {
         aubio_pitch_yin,
         aubio_pitch_mcomb,
         aubio_pitch_schmitt,
-        aubio_pitch_fcomb
+        aubio_pitch_fcomb,
+       aubio_pitch_yinfft
 } aubio_pitchdetection_type;
 
 typedef enum {
@@ -44,6 +45,7 @@ smpl_t aubio_pitchdetection_mcomb(aubio_pitchdetection_t *p, fvec_t * ibuf);
 smpl_t aubio_pitchdetection_yin(aubio_pitchdetection_t *p, fvec_t *ibuf);
 smpl_t aubio_pitchdetection_schmitt(aubio_pitchdetection_t *p, fvec_t *ibuf);
 smpl_t aubio_pitchdetection_fcomb(aubio_pitchdetection_t *p, fvec_t *ibuf);
+smpl_t aubio_pitchdetection_yinfft(aubio_pitchdetection_t *p, fvec_t *ibuf);
 
 void aubio_pitchdetection_set_yinthresh(aubio_pitchdetection_t *p, smpl_t thres);
 void del_aubio_pitchdetection(aubio_pitchdetection_t * p);
diff --git a/src/pitchyinfft.c b/src/pitchyinfft.c
new file mode 100644 (file)
index 0000000..45cf83d
--- /dev/null
@@ -0,0 +1,163 @@
+/*
+   Copyright (C) 2003 Paul Brossier
+
+   This program 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 2 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 General Public License for more details.
+
+   You should have received a copy of the GNU General Public License
+   along with this program; if not, write to the Free Software
+   Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+*/
+
+/* This algorithm was developped by A. de Cheveigne and H. Kawahara and
+ * published in:
+ * 
+ * de Cheveigné, A., Kawahara, H. (2002) "YIN, a fundamental frequency
+ * estimator for speech and music", J. Acoust. Soc. Am. 111, 1917-1930.  
+ *
+ * see http://recherche.ircam.fr/equipes/pcm/pub/people/cheveign.html
+ *
+ * This implementation is using an FFT to compute the square difference
+ * function, which allows spectral weighting 
+ * 
+ */
+
+#include "aubio_priv.h"
+#include "sample.h"
+#include "mathutils.h"
+#include "fft.h"
+#include "pitchyinfft.h"
+
+struct _aubio_pitchyinfft_t {
+        uint_t bufsize;
+        uint_t rate;
+       fvec_t * win;
+       fvec_t * winput;
+       cvec_t * res;
+       fvec_t * sqrmag;
+       fvec_t * weight;
+       cvec_t * fftout;
+       aubio_mfft_t * fft;
+       fvec_t * yinfft;
+};
+
+static const smpl_t freqs[] = {0., 20., 25., 31.5, 40., 50., 63., 80., 100.,
+       125., 160., 200., 250., 315., 400., 500., 630., 800., 1000., 1250.,
+       1600., 2000., 2500., 3150., 4000., 5000., 6300., 8000., 9000., 10000.,
+       12500., 15000., 20000.,  25100};
+
+static const smpl_t weight[] = {-75.8, -70.1, -60.8, -52.1, -44.2, -37.5,
+       -31.3, -25.6, -20.9, -16.5, -12.6, -9.6, -7.0, -4.7, -3.0, -1.8, -0.8,
+       -0.2, -0.0, 0.5, 1.6, 3.2, 5.4, 7.8, 8.1, 5.3, -2.4, -11.1, -12.8,
+       -12.2, -7.4, -17.8, -17.8, -17.8};
+
+aubio_pitchyinfft_t * new_aubio_pitchyinfft (uint_t bufsize)
+{
+  aubio_pitchyinfft_t * p = AUBIO_NEW(aubio_pitchyinfft_t);
+  p->bufsize      = bufsize;
+  p->winput       = new_fvec(bufsize,1);
+  p->fft          = new_aubio_mfft(bufsize, 1);
+  p->fftout       = new_cvec(bufsize,1);
+  p->sqrmag       = new_fvec(bufsize,1);
+  p->res          = new_cvec(bufsize,1);
+  p->yinfft       = new_fvec(bufsize/2+1,1);
+  p->win         = new_fvec(bufsize,1);
+  aubio_window(p->win->data[0], bufsize, aubio_win_hanningz);
+  p->weight      = new_fvec(bufsize/2+1,1);
+  {
+         uint_t i = 0, j = 1;
+         smpl_t freq = 0, a0 = 0, a1 = 0, f0 = 0, f1 = 0;
+         for (i=0; i<p->weight->length; i++) {
+                 freq = (smpl_t)i/(smpl_t)bufsize*(smpl_t)44100.;
+                 while (freq > freqs[j]) {
+                         j +=1;
+                 }
+                 a0 = weight[j-1];
+                 f0 = freqs[j-1];
+                 a1 = weight[j];
+                 f1 = freqs[j];
+                 if (f0 == f1) { // just in case 
+                         p->weight->data[0][i] = a0;
+                 } else if (f0 == 0) { // y = ax+b
+                         p->weight->data[0][i] = (a1-a0)/f1*freq + a0;
+                 } else {
+                         p->weight->data[0][i] = (a1-a0)/(f1-f0)*freq + 
+                                 (a0 - (a1 - a0)/(f1/f0 - 1.));
+                 }
+                 while (freq > freqs[j]) {
+                         j +=1;
+                 }
+                 //AUBIO_DBG("%f\n",p->weight->data[0][i]);
+                 p->weight->data[0][i] = DB2LIN(p->weight->data[0][i]);
+                 //p->weight->data[0][i] = SQRT(DB2LIN(p->weight->data[0][i]));
+         }
+  }
+  return p;
+}
+
+smpl_t aubio_pitchyinfft_detect(aubio_pitchyinfft_t * p, fvec_t * input, smpl_t tol) {
+  uint_t tau, l = 0;
+  uint_t halfperiod;
+  smpl_t tmp = 0, sum = 0;
+  cvec_t * res = (cvec_t *)p->res;
+  fvec_t * yin = (fvec_t *)p->yinfft;
+  for (l=0; l < input->length; l++){
+         p->winput->data[0][l] = p->win->data[0][l] * input->data[0][l];
+  }
+  aubio_mfft_do(p->fft,p->winput,p->fftout);
+  for (l=0; l < p->fftout->length; l++){
+         p->sqrmag->data[0][l] = SQR(p->fftout->norm[0][l]);
+         p->sqrmag->data[0][(p->fftout->length-1)*2-l] = 
+               SQR(p->fftout->norm[0][l]);
+         p->sqrmag->data[0][l] *= p->weight->data[0][l]; 
+         p->sqrmag->data[0][(p->fftout->length-1)*2-l] *=
+                p->weight->data[0][l];
+  }
+  for (l=0; l < p->sqrmag->length/2+1; l++) {
+         sum += p->sqrmag->data[0][l];
+  }
+  sum *= 2.;
+  aubio_mfft_do(p->fft,p->sqrmag,res);
+  yin->data[0][0] = 1.; 
+  for (tau=1; tau < yin->length; tau++) {
+         yin->data[0][tau] = sum -
+                 res->norm[0][tau+1]*COS(res->phas[0][tau+1]); 
+         tmp += yin->data[0][tau];
+         yin->data[0][tau] *= tau/tmp;
+  }
+  tau = vec_min_elem(yin); 
+  if (yin->data[0][tau] < tol) {
+         /* no interpolation */
+         //return tau+2;
+         /* 3 point quadratic interpolation */
+         //return vec_quadint_min(yin,tau,1);
+         /* additional check nlikely octave doubling in higher frequencies */
+         if (tau>35) {
+                 return vec_quadint_min(yin,tau,1)+1;
+         } else {
+                 /* should compare the minimum value of each interpolated peaks */
+                 halfperiod = FLOOR(tau/2+.5);
+                 if (yin->data[0][halfperiod] < tol)
+                         return vec_quadint_min(yin,halfperiod,1)+1;
+                 else
+                         return vec_quadint_min(yin,tau,1)+1;
+         }
+  } else
+         return 0;
+}
+
+void del_aubio_pitchyinfft(aubio_pitchyinfft_t *p){
+       del_fvec(p->win);
+       del_aubio_mfft(p->fft);
+       del_fvec(p->yinfft);
+       del_fvec(p->sqrmag);
+       del_cvec(p->res);
+       del_cvec(p->fftout);
+}
diff --git a/src/pitchyinfft.h b/src/pitchyinfft.h
new file mode 100644 (file)
index 0000000..11161fd
--- /dev/null
@@ -0,0 +1,45 @@
+/*
+   Copyright (C) 2003 Paul Brossier
+
+   This program 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 2 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 General Public License for more details.
+
+   You should have received a copy of the GNU General Public License
+   along with this program; if not, write to the Free Software
+   Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+*/
+
+/* This algorithm was developped by A. de Cheveigne and H. Kawahara and
+ * published in:
+ * 
+ * de Cheveigné, A., Kawahara, H. (2002) "YIN, a fundamental frequency
+ * estimator for speech and music", J. Acoust. Soc. Am. 111, 1917-1930.  
+ *
+ * see http://recherche.ircam.fr/equipes/pcm/pub/people/cheveign.html
+ */
+
+#ifndef PITCHYINFFT_H
+#define PITCHYINFFT_H
+
+#ifdef __cplusplus
+extern "C" {
+#endif
+
+typedef struct _aubio_pitchyinfft_t aubio_pitchyinfft_t;
+
+smpl_t aubio_pitchyinfft_detect (aubio_pitchyinfft_t *p, fvec_t * input, smpl_t tol);
+aubio_pitchyinfft_t * new_aubio_pitchyinfft (uint_t bufsize);
+void del_aubio_pitchyinfft (aubio_pitchyinfft_t *p);
+
+#ifdef __cplusplus
+}
+#endif
+
+#endif /*PITCHYINFFT_H*/ 
index c16a853..f8e314e 100644 (file)
@@ -190,7 +190,8 @@ typedef enum {
         aubio_pitch_yin,
         aubio_pitch_mcomb,
         aubio_pitch_schmitt,
-        aubio_pitch_fcomb
+        aubio_pitch_fcomb,
+       aubio_pitch_yinfft
 } aubio_pitchdetection_type;
 
 typedef enum {
@@ -217,7 +218,7 @@ aubio_pitchdetection_t * new_aubio_pitchdetection(uint_t bufsize,
 
 
 /* pitch mcomb */
-aubio_pitchmcomb_t * new_aubio_pitchmcomb(uint_t size, uint_t channels, uint_t samplerate);
+aubio_pitchmcomb_t * new_aubio_pitchmcomb(uint_t bufsize, uint_t hopsize, uint_t channels, uint_t samplerate);
 smpl_t aubio_pitchmcomb_detect(aubio_pitchmcomb_t * p, cvec_t * fftgrain);
 uint_t aubio_pitch_cands(aubio_pitchmcomb_t * p, cvec_t * fftgrain, smpl_t * cands);
 void del_aubio_pitchmcomb (aubio_pitchmcomb_t *p);
@@ -226,7 +227,7 @@ void del_aubio_pitchmcomb (aubio_pitchmcomb_t *p);
 void aubio_pitchyin_diff(fvec_t *input, fvec_t *yin);
 void aubio_pitchyin_getcum(fvec_t *yin);
 uint_t aubio_pitchyin_getpitch(fvec_t *yin);
-uint_t aubio_pitchyin_getpitchfast(fvec_t * input, fvec_t *yin, smpl_t tol);
+smpl_t aubio_pitchyin_getpitchfast(fvec_t * input, fvec_t *yin, smpl_t tol);
 
 /* pitch schmitt */
 aubio_pitchschmitt_t * new_aubio_pitchschmitt (uint_t size, uint_t samplerate);