src/pitch/pitchyinfft.c: comment out debug output
[aubio.git] / src / pitch / pitchyinfft.c
index 6e42cf0..b613f60 100644 (file)
@@ -36,38 +36,46 @@ struct _aubio_pitchyinfft_t
   aubio_fft_t *fft;   /**< fft object to compute square difference function */
   fvec_t *yinfft;     /**< Yin function */
   smpl_t tol;         /**< Yin tolerance */
-  smpl_t confidence;  /**< confidence */
+  uint_t peak_pos;    /**< currently selected peak pos*/
+  uint_t short_period; /** shortest period under which to check for octave error */
 };
 
-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 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., -1.
 };
 
-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
+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.60,  -7.00,  -4.70,  -3.00,  -1.80,  -0.80,  -0.20,  -0.00,   0.50,
+   1.60,   3.20,   5.40,   7.80,   8.10,   5.30,  -2.40,  -11.1,  -12.8,  -12.2,
+  -7.40,  -17.8,  -17.8,  -17.8
 };
 
 aubio_pitchyinfft_t *
-new_aubio_pitchyinfft (uint_t bufsize)
+new_aubio_pitchyinfft (uint_t samplerate, uint_t bufsize)
 {
+  uint_t i = 0, j = 1;
+  smpl_t freq = 0, a0 = 0, a1 = 0, f0 = 0, f1 = 0;
   aubio_pitchyinfft_t *p = AUBIO_NEW (aubio_pitchyinfft_t);
   p->winput = new_fvec (bufsize);
   p->fft = new_aubio_fft (bufsize);
+  if (!p->fft) goto beach;
   p->fftout = new_fvec (bufsize);
   p->sqrmag = new_fvec (bufsize);
   p->yinfft = new_fvec (bufsize / 2 + 1);
   p->tol = 0.85;
+  p->peak_pos = 0;
   p->win = new_aubio_window ("hanningz", bufsize);
   p->weight = new_fvec (bufsize / 2 + 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]) {
+    freq = (smpl_t) i / (smpl_t) bufsize *(smpl_t) samplerate;
+    while (freq > freqs[j] && freqs[j] > 0) {
+      //AUBIO_DBG("freq %3.5f > %3.5f \tsamplerate %d (Hz) \t"
+      //    "(weight length %d, bufsize %d) %d %d\n", freq, freqs[j],
+      //    samplerate, p->weight->length, bufsize, i, j);
       j += 1;
     }
     a0 = weight[j - 1];
@@ -89,11 +97,18 @@ new_aubio_pitchyinfft (uint_t bufsize)
     p->weight->data[i] = DB2LIN (p->weight->data[i]);
     //p->weight->data[i] = SQRT(DB2LIN(p->weight->data[i]));
   }
+  // check for octave errors above 1300 Hz
+  p->short_period = (uint_t)ROUND(samplerate / 1300.);
   return p;
+
+beach:
+  if (p->winput) del_fvec(p->winput);
+  AUBIO_FREE(p);
+  return NULL;
 }
 
 void
-aubio_pitchyinfft_do (aubio_pitchyinfft_t * p, fvec_t * input, fvec_t * output)
+aubio_pitchyinfft_do (aubio_pitchyinfft_t * p, const fvec_t * input, fvec_t * output)
 {
   uint_t tau, l;
   uint_t length = p->fftout->length;
@@ -102,9 +117,7 @@ aubio_pitchyinfft_do (aubio_pitchyinfft_t * p, fvec_t * input, fvec_t * output)
   fvec_t *yin = p->yinfft;
   smpl_t tmp = 0., sum = 0.;
   // window the input
-  for (l = 0; l < input->length; l++) {
-    p->winput->data[l] = p->win->data[l] * input->data[l];
-  }
+  fvec_weighted_copy(input, p->win, p->winput);
   // get the real / imag parts of its fft
   aubio_fft_do_complex (p->fft, p->winput, fftout);
   // get the squared magnitude spectrum, applying some weight
@@ -130,7 +143,11 @@ aubio_pitchyinfft_do (aubio_pitchyinfft_t * p, fvec_t * input, fvec_t * output)
     yin->data[tau] = sum - fftout->data[tau];
     // and the cumulative mean normalized difference function
     tmp += yin->data[tau];
-    yin->data[tau] *= tau / tmp;
+    if (tmp != 0) {
+      yin->data[tau] *= tau / tmp;
+    } else {
+      yin->data[tau] = 1.;
+    }
   }
   // find best candidates
   tau = fvec_min_elem (yin);
@@ -140,19 +157,21 @@ aubio_pitchyinfft_do (aubio_pitchyinfft_t * p, fvec_t * input, fvec_t * output)
     //return;
 
     // 3 point quadratic interpolation
-    //return fvec_quadint_min(yin,tau,1);
+    //return fvec_quadratic_peak_pos (yin,tau,1);
     /* additional check for (unlikely) octave doubling in higher frequencies */
-    if (tau > 35) {
-      output->data[0] = fvec_quadint (yin, tau);
+    if (tau > p->short_period) {
+      output->data[0] = fvec_quadratic_peak_pos (yin, tau);
     } else {
       /* should compare the minimum value of each interpolated peaks */
       halfperiod = FLOOR (tau / 2 + .5);
       if (yin->data[halfperiod] < p->tol)
-        output->data[0] = fvec_quadint (yin, halfperiod);
+        p->peak_pos = halfperiod;
       else
-        output->data[0] = fvec_quadint (yin, tau);
+        p->peak_pos = tau;
+      output->data[0] = fvec_quadratic_peak_pos (yin, p->peak_pos);
     }
   } else {
+    p->peak_pos = 0;
     output->data[0] = 0.;
   }
 }
@@ -172,8 +191,7 @@ del_aubio_pitchyinfft (aubio_pitchyinfft_t * p)
 
 smpl_t
 aubio_pitchyinfft_get_confidence (aubio_pitchyinfft_t * o) {
-  o->confidence = 1. - fvec_min (o->yinfft);
-  return o->confidence;
+  return 1. - o->yinfft->data[o->peak_pos];
 }
 
 uint_t