libs/vamp-plugins/TruePeak.cpp

   1 /*
   2  *  Copyright (C) 2013-2016 Robin Gareus <robin@gareus.org>
   3  *  Copyright (C) 2006-2012 Fons Adriaensen <fons@linuxaudio.org>
   4  *  Copyright (C) 2006 Chris Cannam
   5  *
   6  *  This program 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  *  This program 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 this program.  If not, see <http://www.gnu.org/licenses/>.
  18  */
  19
  20 #include <stdlib.h>
  21 #include <stdio.h>
  22 #include <string.h>
  23 #include <math.h>
  24
  25 #include "TruePeak.h"
  26
  27 namespace TruePeakMeter {
  28
  29 static double sinc (double x)
  30 {
  31         x = fabs (x);
  32         if (x < 1e-6) return 1.0;
  33         x *= M_PI;
  34         return sin (x) / x;
  35 }
  36
  37 static double wind (double x)
  38 {
  39         x = fabs (x);
  40         if (x >= 1.0) return 0.0f;
  41         x *= M_PI;
  42         return 0.384 + 0.500 * cos (x) + 0.116 * cos (2 * x);
  43 }
  44
  45 Resampler_table  *Resampler_table::_list = 0;
  46 Resampler_mutex   Resampler_table::_mutex;
  47
  48 Resampler_table::Resampler_table (double fr, unsigned int hl, unsigned int np)
  49         : _next (0)
  50         , _refc (0)
  51         , _fr (fr)
  52         , _hl (hl)
  53         , _np (np)
  54 {
  55         unsigned int  i, j;
  56         double        t;
  57         float        *p;
  58
  59         _ctab = new float [hl * (np + 1)];
  60         p = _ctab;
  61         for (j = 0; j <= np; j++)
  62         {
  63                 t = (double) j / (double) np;
  64                 for (i = 0; i < hl; i++)
  65                 {
  66                         p [hl - i - 1] = (float)(fr * sinc (t * fr) * wind (t / hl));
  67                         t += 1;
  68                 }
  69                 p += hl;
  70         }
  71 }
  72
  73 Resampler_table::~Resampler_table (void)
  74 {
  75         delete[] _ctab;
  76 }
  77
  78 Resampler_table *
  79 Resampler_table::create (double fr, unsigned int hl, unsigned int np)
  80 {
  81         Resampler_table *P;
  82
  83         _mutex.lock ();
  84         P = _list;
  85         while (P)
  86         {
  87                 if ((fr >= P->_fr * 0.999) && (fr <= P->_fr * 1.001) && (hl == P->_hl) && (np == P->_np))
  88                 {
  89                         P->_refc++;
  90                         _mutex.unlock ();
  91                         return P;
  92                 }
  93                 P = P->_next;
  94         }
  95         P = new Resampler_table (fr, hl, np);
  96         P->_refc = 1;
  97         P->_next = _list;
  98         _list = P;
  99         _mutex.unlock ();
 100         return P;
 101 }
 102
 103 void
 104 Resampler_table::destroy (Resampler_table *T)
 105 {
 106         Resampler_table *P, *Q;
 107
 108         _mutex.lock ();
 109         if (T)
 110         {
 111                 T->_refc--;
 112                 if (T->_refc == 0)
 113                 {
 114                         P = _list;
 115                         Q = 0;
 116                         while (P)
 117                         {
 118                                 if (P == T)
 119                                 {
 120                                         if (Q) Q->_next = T->_next;
 121                                         else      _list = T->_next;
 122                                         break;
 123                                 }
 124                                 Q = P;
 125                                 P = P->_next;
 126                         }
 127                         delete T;
 128                 }
 129         }
 130         _mutex.unlock ();
 131 }
 132
 133 static unsigned int
 134 gcd (unsigned int a, unsigned int b)
 135 {
 136         if (a == 0) return b;
 137         if (b == 0) return a;
 138         while (1)
 139         {
 140                 if (a > b)
 141                 {
 142                         a = a % b;
 143                         if (a == 0) return b;
 144                         if (a == 1) return 1;
 145                 }
 146                 else
 147                 {
 148                         b = b % a;
 149                         if (b == 0) return a;
 150                         if (b == 1) return 1;
 151                 }
 152         }
 153         return 1;
 154 }
 155
 156 Resampler::Resampler (void)
 157         : _table (0)
 158         , _nchan (0)
 159         , _buff  (0)
 160 {
 161         reset ();
 162 }
 163
 164 Resampler::~Resampler (void)
 165 {
 166         clear ();
 167 }
 168
 169 int
 170 Resampler::setup (unsigned int fs_inp,
 171                   unsigned int fs_out,
 172                   unsigned int nchan,
 173                   unsigned int hlen)
 174 {
 175         if ((hlen < 8) || (hlen > 96)) return 1;
 176         return setup (fs_inp, fs_out, nchan, hlen, 1.0 - 2.6 / hlen);
 177 }
 178
 179 int
 180 Resampler::setup (unsigned int fs_inp,
 181                   unsigned int fs_out,
 182                   unsigned int nchan,
 183                   unsigned int hlen,
 184                   double       frel)
 185 {
 186         unsigned int       g, h, k, n, s;
 187         double             r;
 188         float              *B = 0;
 189         Resampler_table    *T = 0;
 190
 191         k = s = 0;
 192         if (fs_inp && fs_out && nchan)
 193         {
 194                 r = (double) fs_out / (double) fs_inp;
 195                 g = gcd (fs_out, fs_inp);
 196                 n = fs_out / g;
 197                 s = fs_inp / g;
 198                 if ((16 * r >= 1) && (n <= 1000))
 199                 {
 200                         h = hlen;
 201                         k = 250;
 202                         if (r < 1)
 203                         {
 204                                 frel *= r;
 205                                 h = (unsigned int)(ceil (h / r));
 206                                 k = (unsigned int)(ceil (k / r));
 207                         }
 208                         T = Resampler_table::create (frel, h, n);
 209                         B = new float [nchan * (2 * h - 1 + k)];
 210                 }
 211         }
 212         clear ();
 213         if (T)
 214         {
 215                 _table = T;
 216                 _buff  = B;
 217                 _nchan = nchan;
 218                 _inmax = k;
 219                 _pstep = s;
 220                 return reset ();
 221         }
 222         else return 1;
 223 }
 224
 225 void
 226 Resampler::clear (void)
 227 {
 228         Resampler_table::destroy (_table);
 229         delete[] _buff;
 230         _buff  = 0;
 231         _table = 0;
 232         _nchan = 0;
 233         _inmax = 0;
 234         _pstep = 0;
 235         reset ();
 236 }
 237
 238 double
 239 Resampler::inpdist (void) const
 240 {
 241         if (!_table) return 0;
 242         return (int)(_table->_hl + 1 - _nread) - (double)_phase / _table->_np;
 243 }
 244
 245 int
 246 Resampler::inpsize (void) const
 247 {
 248         if (!_table) return 0;
 249         return 2 * _table->_hl;
 250 }
 251
 252 int
 253 Resampler::reset (void)
 254 {
 255         if (!_table) return 1;
 256
 257         inp_count = 0;
 258         out_count = 0;
 259         inp_data = 0;
 260         out_data = 0;
 261         _index = 0;
 262         _nread = 0;
 263         _nzero = 0;
 264         _phase = 0;
 265         if (_table)
 266         {
 267                 _nread = 2 * _table->_hl;
 268                 return 0;
 269         }
 270         return 1;
 271 }
 272
 273 int
 274 Resampler::process (void)
 275 {
 276         unsigned int   hl, ph, np, dp, in, nr, nz, i, n, c;
 277         float          *p1, *p2;
 278
 279         if (!_table) return 1;
 280
 281         hl = _table->_hl;
 282         np = _table->_np;
 283         dp = _pstep;
 284         in = _index;
 285         nr = _nread;
 286         ph = _phase;
 287         nz = _nzero;
 288         n = (2 * hl - nr) * _nchan;
 289         p1 = _buff + in * _nchan;
 290         p2 = p1 + n;
 291
 292         while (out_count)
 293         {
 294                 if (nr)
 295                 {
 296                         if (inp_count == 0) break;
 297                         if (inp_data)
 298                         {
 299                                 for (c = 0; c < _nchan; c++) p2 [c] = inp_data [c];
 300                                 inp_data += _nchan;
 301                                 nz = 0;
 302                         }
 303                         else
 304                         {
 305                                 for (c = 0; c < _nchan; c++) p2 [c] = 0;
 306                                 if (nz < 2 * hl) nz++;
 307                         }
 308                         nr--;
 309                         p2 += _nchan;
 310                         inp_count--;
 311                 }
 312                 else
 313                 {
 314                         if (out_data)
 315                         {
 316                                 if (nz < 2 * hl)
 317                                 {
 318                                         float *c1 = _table->_ctab + hl * ph;
 319                                         float *c2 = _table->_ctab + hl * (np - ph);
 320                                         for (c = 0; c < _nchan; c++)
 321                                         {
 322                                                 float *q1 = p1 + c;
 323                                                 float *q2 = p2 + c;
 324                                                 float s = 1e-20f;
 325                                                 for (i = 0; i < hl; i++)
 326                                                 {
 327                                                         q2 -= _nchan;
 328                                                         s += *q1 * c1 [i] + *q2 * c2 [i];
 329                                                         q1 += _nchan;
 330                                                 }
 331                                                 *out_data++ = s - 1e-20f;
 332                                         }
 333                                 }
 334                                 else
 335                                 {
 336                                         for (c = 0; c < _nchan; c++) *out_data++ = 0;
 337                                 }
 338                         }
 339                         out_count--;
 340
 341                         ph += dp;
 342                         if (ph >= np)
 343                         {
 344                                 nr = ph / np;
 345                                 ph -= nr * np;
 346                                 in += nr;
 347                                 p1 += nr * _nchan;;
 348                                 if (in >= _inmax)
 349                                 {
 350                                         n = (2 * hl - nr) * _nchan;
 351                                         memcpy (_buff, p1, n * sizeof (float));
 352                                         in = 0;
 353                                         p1 = _buff;
 354                                         p2 = p1 + n;
 355                                 }
 356                         }
 357                 }
 358         }
 359         _index = in;
 360         _nread = nr;
 361         _phase = ph;
 362         _nzero = nz;
 363
 364         return 0;
 365 }
 366
 367 TruePeakdsp::TruePeakdsp (void)
 368         : _m (0)
 369         , _p (0)
 370         , _res (true)
 371         , _buf (NULL)
 372 {
 373 }
 374
 375 TruePeakdsp::~TruePeakdsp (void)
 376 {
 377         free(_buf);
 378 }
 379
 380 void
 381 TruePeakdsp::process (float const *d, int n)
 382 {
 383         _src.inp_count = n;
 384         _src.inp_data = d;
 385         _src.out_count = n * 4;
 386         _src.out_data = _buf;
 387         _src.process ();
 388
 389         float x = 0;
 390         float v;
 391         float *b = _buf;
 392         while (n--) {
 393                 v = fabsf(*b++);
 394                 if (v > x) x = v;
 395                 v = fabsf(*b++);
 396                 if (v > x) x = v;
 397                 v = fabsf(*b++);
 398                 if (v > x) x = v;
 399                 v = fabsf(*b++);
 400                 if (v > x) x = v;
 401         }
 402
 403         if (_res) {
 404                 _m = x;
 405                 _res = false;
 406         } else if (x > _m) {
 407                 _m = x;
 408         }
 409
 410         if (_res_peak) {
 411                 _p = x;
 412                 _res_peak = false;
 413         } else if (x > _p) {
 414                 _p = x;
 415         }
 416 }
 417
 418 float
 419 TruePeakdsp::read (void)
 420 {
 421         _res = true;
 422         return _m;
 423 }
 424
 425 void
 426 TruePeakdsp::read (float &m, float &p)
 427 {
 428         _res = true;
 429         _res_peak = true;
 430         m = _m;
 431         p = _p;
 432 }
 433
 434 void
 435 TruePeakdsp::reset ()
 436 {
 437         _res = true;
 438         _m = 0;
 439         _p = 0;
 440 }
 441
 442 bool
 443 TruePeakdsp::init (float fsamp)
 444 {
 445         _src.setup(fsamp, fsamp * 4.0, 1, 24, 1.0);
 446         _buf = (float*) malloc(32768 * sizeof(float));
 447         if (!_buf) {
 448                 return false;
 449         }
 450
 451         /* q/d initialize */
 452         float zero[8192];
 453         for (int i = 0; i < 8192; ++i) {
 454                 zero[i]= 0.0;
 455         }
 456         _src.inp_count = 8192;
 457         _src.inp_data = zero;
 458         _src.out_count = 32768;
 459         _src.out_data = _buf;
 460         _src.process ();
 461         return true;
 462 }
 463
 464 }
 465
 466 ///////////////////////////////////////////////////////////////////////////////
 467
 468 using std::string;
 469 using std::vector;
 470 using std::cerr;
 471 using std::endl;
 472 using namespace TruePeakMeter;
 473
 474 VampTruePeak::VampTruePeak(float inputSampleRate)
 475     : Plugin(inputSampleRate)
 476     , m_blockSize(0)
 477     , m_rate (inputSampleRate)
 478 {
 479 }
 480
 481 VampTruePeak::~VampTruePeak()
 482 {
 483 }
 484
 485 string
 486 VampTruePeak::getIdentifier() const
 487 {
 488         return "dBTP";
 489 }
 490
 491 string
 492 VampTruePeak::getName() const
 493 {
 494         return "dBTP Meter";
 495 }
 496
 497 string
 498 VampTruePeak::getDescription() const
 499 {
 500         return "True Peak Meter (4x Oversampling)";
 501 }
 502
 503 string
 504 VampTruePeak::getMaker() const
 505 {
 506         return "Robin Gareus, Fons Adrianesen";
 507 }
 508
 509 int
 510 VampTruePeak::getPluginVersion() const
 511 {
 512         return 2;
 513 }
 514
 515 string
 516 VampTruePeak::getCopyright() const
 517 {
 518         return "GPL version 3 or later";
 519 }
 520
 521 bool
 522 VampTruePeak::initialise(size_t channels, size_t stepSize, size_t blockSize)
 523 {
 524         if (channels < getMinChannelCount() ||
 525                         channels > getMaxChannelCount()) {
 526                 return false;
 527         }
 528
 529         if (blockSize == 0 || blockSize > 8192) {
 530                 return false;
 531         }
 532
 533         if (!_meter.init (m_inputSampleRate)) {
 534                 return false;
 535         }
 536
 537         m_blockSize = blockSize;
 538
 539         return true;
 540 }
 541
 542 void
 543 VampTruePeak::reset()
 544 {
 545         _meter.reset ();
 546 }
 547
 548 VampTruePeak::OutputList
 549 VampTruePeak::getOutputDescriptors() const
 550 {
 551         OutputList list;
 552
 553         OutputDescriptor zc;
 554         zc.identifier = "level";
 555         zc.name = "TruePeak";
 556         zc.description = "TruePeak (4x Oversampling)";
 557         zc.unit = "dbTP";
 558         zc.hasFixedBinCount = true;
 559         zc.binCount = 0;
 560         zc.hasKnownExtents = false;
 561         zc.isQuantized = false;
 562         zc.sampleType = OutputDescriptor::OneSamplePerStep;
 563         list.push_back(zc);
 564
 565         zc.identifier = "peaks";
 566         zc.name = "TruePeakPeaks";
 567         zc.description = "Location of Peaks above -1dBTP";
 568         zc.unit = "sec";
 569         zc.hasFixedBinCount = true;
 570         zc.binCount = 0;
 571         zc.hasKnownExtents = false;
 572         zc.isQuantized = false;
 573         zc.sampleType = OutputDescriptor::OneSamplePerStep;
 574         list.push_back(zc);
 575
 576         return list;
 577 }
 578
 579 VampTruePeak::FeatureSet
 580 VampTruePeak::process(const float *const *inputBuffers,
 581                       Vamp::RealTime timestamp)
 582 {
 583         if (m_blockSize == 0) {
 584                 cerr << "ERROR: VampTruePeak::process: "
 585                         << "VampTruePeak has not been initialised"
 586                         << endl;
 587                 return FeatureSet();
 588         }
 589
 590         _meter.process (inputBuffers[0], m_blockSize);
 591
 592         // TODO optional (not rt safe)
 593         if (_meter.read () >= .89125 /* -1dBTP */) {
 594                 long f = Vamp::RealTime::realTime2Frame (timestamp, m_rate);
 595                 _above_m1.values.push_back ((float) f);
 596         }
 597
 598         return FeatureSet();
 599 }
 600
 601 VampTruePeak::FeatureSet
 602 VampTruePeak::getRemainingFeatures()
 603 {
 604         FeatureSet returnFeatures;
 605
 606         float m, p;
 607         _meter.read(m, p);
 608
 609         Feature dbtp;
 610         dbtp.hasTimestamp = false;
 611         dbtp.values.push_back(p);
 612         returnFeatures[0].push_back(dbtp);
 613
 614         _above_m1.hasTimestamp = false;
 615         returnFeatures[1].push_back(_above_m1);
 616
 617         return returnFeatures;
 618 }