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