1 /* See COPYING file for copyright and license details. */
2
3 #include "ebur128.h"
4
5 #include <float.h>
6 #include <limits.h>
7 #include <math.h> /* You may have to define _USE_MATH_DEFINES if you use MSVC */
8 #include <stdio.h>
9 #include <stdlib.h>
10
11 /* This can be replaced by any BSD-like queue implementation. */
12 #include <sys/queue.h>
13
14 #define CHECK_ERROR(condition, errorcode, goto_point) \
15 if ((condition)) { \
16 errcode = (errorcode); \
17 goto goto_point; \
18 }
19
20 STAILQ_HEAD(ebur128_double_queue, ebur128_dq_entry);
21 struct ebur128_dq_entry {
22 double z;
23 STAILQ_ENTRY(ebur128_dq_entry) entries;
24 };
25
26 #define ALMOST_ZERO 0.000001
27
28 typedef struct { /* Data structure for polyphase FIR interpolator */
29 unsigned int factor; /* Interpolation factor of the interpolator */
30 unsigned int taps; /* Taps (prefer odd to increase zero coeffs) */
31 unsigned int channels; /* Number of channels */
32 unsigned int delay; /* Size of delay buffer */
33 struct {
34 unsigned int count; /* Number of coefficients in this subfilter */
35 unsigned int* index; /* Delay index of corresponding filter coeff */
36 double* coeff; /* List of subfilter coefficients */
37 }* filter; /* List of subfilters (one for each factor) */
38 float** z; /* List of delay buffers (one for each channel) */
39 unsigned int zi; /* Current delay buffer index */
40 } interpolator;
41
42 struct ebur128_state_internal {
43 /** Filtered audio data (used as ring buffer). */
44 double* audio_data;
45 /** Size of audio_data array. */
46 size_t audio_data_frames;
47 /** Current index for audio_data. */
48 size_t audio_data_index;
49 /** How many frames are needed for a gating block. Will correspond to 400ms
50 * of audio at initialization, and 100ms after the first block (75% overlap
51 * as specified in the 2011 revision of BS1770). */
52 unsigned long needed_frames;
53 /** The channel map. Has as many elements as there are channels. */
54 int* channel_map;
55 /** How many samples fit in 100ms (rounded). */
56 unsigned long samples_in_100ms;
57 /** BS.1770 filter coefficients (nominator). */
58 double b[5];
59 /** BS.1770 filter coefficients (denominator). */
60 double a[5];
61 /** BS.1770 filter state. */
62 double v[5][5];
63 /** Linked list of block energies. */
64 struct ebur128_double_queue block_list;
65 unsigned long block_list_max;
66 unsigned long block_list_size;
67 /** Linked list of 3s-block energies, used to calculate LRA. */
68 struct ebur128_double_queue short_term_block_list;
69 unsigned long st_block_list_max;
70 unsigned long st_block_list_size;
71 int use_histogram;
72 unsigned long *block_energy_histogram;
73 unsigned long *short_term_block_energy_histogram;
74 /** Keeps track of when a new short term block is needed. */
75 size_t short_term_frame_counter;
76 /** Maximum sample peak, one per channel */
77 double* sample_peak;
78 double* prev_sample_peak;
79 /** Maximum true peak, one per channel */
80 double* true_peak;
81 double* prev_true_peak;
82 interpolator* interp;
83 float* resampler_buffer_input;
84 size_t resampler_buffer_input_frames;
85 float* resampler_buffer_output;
86 size_t resampler_buffer_output_frames;
87 /** The maximum window duration in ms. */
88 unsigned long window;
89 unsigned long history;
90 };
91
92 static double relative_gate = -10.0;
93
94 /* Those will be calculated when initializing the library */
95 static double relative_gate_factor;
96 static double minus_twenty_decibels;
97 static double histogram_energies[1000];
98 static double histogram_energy_boundaries[1001];
99
interp_create(unsigned int taps,unsigned int factor,unsigned int channels)100 static interpolator* interp_create(unsigned int taps, unsigned int factor, unsigned int channels) {
101 interpolator* interp = calloc(1, sizeof(interpolator));
102 unsigned int j = 0;
103
104 interp->taps = taps;
105 interp->factor = factor;
106 interp->channels = channels;
107 interp->delay = (interp->taps + interp->factor - 1) / interp->factor;
108
109 /* Initialize the filter memory
110 * One subfilter per interpolation factor. */
111 interp->filter = calloc(interp->factor, sizeof(*interp->filter));
112 for (j = 0; j < interp->factor; j++) {
113 interp->filter[j].index = calloc(interp->delay, sizeof(unsigned int));
114 interp->filter[j].coeff = calloc(interp->delay, sizeof(double));
115 }
116 /* One delay buffer per channel. */
117 interp->z = calloc(interp->channels, sizeof(float*));
118 for (j = 0; j < interp->channels; j++) {
119 interp->z[j] = calloc( interp->delay, sizeof(float) );
120 }
121
122 /* Calculate the filter coefficients */
123 for (j = 0; j < interp->taps; j++) {
124 /* Calculate sinc */
125 double m = (double)j - (double)(interp->taps - 1) / 2.0;
126 double c = 1.0;
127 if (fabs(m) > ALMOST_ZERO) {
128 c = sin(m * M_PI / interp->factor) / (m * M_PI / interp->factor);
129 }
130 /* Apply Hanning window */
131 c *= 0.5 * (1 - cos(2 * M_PI * j / (interp->taps - 1)));
132
133 if (fabs(c) > ALMOST_ZERO) { /* Ignore any zero coeffs. */
134 /* Put the coefficient into the correct subfilter */
135 unsigned int f = j % interp->factor;
136 unsigned int t = interp->filter[f].count++;
137 interp->filter[f].coeff[t] = c;
138 interp->filter[f].index[t] = j / interp->factor;
139 }
140 }
141 return interp;
142 }
143
interp_destroy(interpolator * interp)144 static void interp_destroy(interpolator* interp) {
145 unsigned int j = 0;
146 if (!interp) {
147 return;
148 }
149 for (j = 0; j < interp->factor; j++) {
150 free(interp->filter[j].index);
151 free(interp->filter[j].coeff);
152 }
153 free(interp->filter);
154 for (j = 0; j < interp->channels; j++) {
155 free(interp->z[j]);
156 }
157 free(interp->z);
158 free(interp);
159 }
160
interp_process(interpolator * interp,size_t frames,float * in,float * out)161 static size_t interp_process(interpolator* interp, size_t frames, float* in, float* out) {
162 size_t frame = 0;
163 unsigned int chan = 0;
164 unsigned int f = 0;
165 unsigned int t = 0;
166 unsigned int out_stride = interp->channels * interp->factor;
167 float* outp = 0;
168 double acc = 0;
169 double c = 0;
170
171 for (frame = 0; frame < frames; frame++) {
172 for (chan = 0; chan < interp->channels; chan++) {
173 /* Add sample to delay buffer */
174 interp->z[chan][interp->zi] = *in++;
175 /* Apply coefficients */
176 outp = out + chan;
177 for (f = 0; f < interp->factor; f++) {
178 acc = 0.0;
179 for (t = 0; t < interp->filter[f].count; t++) {
180 int i = (int)interp->zi - (int)interp->filter[f].index[t];
181 if (i < 0) {
182 i += interp->delay;
183 }
184 c = interp->filter[f].coeff[t];
185 acc += interp->z[chan][i] * c;
186 }
187 *outp = (float)acc;
188 outp += interp->channels;
189 }
190 }
191 out += out_stride;
192 interp->zi++;
193 if (interp->zi == interp->delay) {
194 interp->zi = 0;
195 }
196 }
197
198 return frames * interp->factor;
199 }
200
ebur128_init_filter(ebur128_state * st)201 static void ebur128_init_filter(ebur128_state* st) {
202 int i, j;
203
204 double f0 = 1681.974450955533;
205 double G = 3.999843853973347;
206 double Q = 0.7071752369554196;
207
208 double K = tan(M_PI * f0 / (double) st->samplerate);
209 double Vh = pow(10.0, G / 20.0);
210 double Vb = pow(Vh, 0.4996667741545416);
211
212 double pb[3] = {0.0, 0.0, 0.0};
213 double pa[3] = {1.0, 0.0, 0.0};
214 double rb[3] = {1.0, -2.0, 1.0};
215 double ra[3] = {1.0, 0.0, 0.0};
216
217 double a0 = 1.0 + K / Q + K * K ;
218 pb[0] = (Vh + Vb * K / Q + K * K) / a0;
219 pb[1] = 2.0 * (K * K - Vh) / a0;
220 pb[2] = (Vh - Vb * K / Q + K * K) / a0;
221 pa[1] = 2.0 * (K * K - 1.0) / a0;
222 pa[2] = (1.0 - K / Q + K * K) / a0;
223
224 /* fprintf(stderr, "%.14f %.14f %.14f %.14f %.14f\n",
225 b1[0], b1[1], b1[2], a1[1], a1[2]); */
226
227 f0 = 38.13547087602444;
228 Q = 0.5003270373238773;
229 K = tan(M_PI * f0 / (double) st->samplerate);
230
231 ra[1] = 2.0 * (K * K - 1.0) / (1.0 + K / Q + K * K);
232 ra[2] = (1.0 - K / Q + K * K) / (1.0 + K / Q + K * K);
233
234 /* fprintf(stderr, "%.14f %.14f\n", a2[1], a2[2]); */
235
236 st->d->b[0] = pb[0] * rb[0];
237 st->d->b[1] = pb[0] * rb[1] + pb[1] * rb[0];
238 st->d->b[2] = pb[0] * rb[2] + pb[1] * rb[1] + pb[2] * rb[0];
239 st->d->b[3] = pb[1] * rb[2] + pb[2] * rb[1];
240 st->d->b[4] = pb[2] * rb[2];
241
242 st->d->a[0] = pa[0] * ra[0];
243 st->d->a[1] = pa[0] * ra[1] + pa[1] * ra[0];
244 st->d->a[2] = pa[0] * ra[2] + pa[1] * ra[1] + pa[2] * ra[0];
245 st->d->a[3] = pa[1] * ra[2] + pa[2] * ra[1];
246 st->d->a[4] = pa[2] * ra[2];
247
248 for (i = 0; i < 5; ++i) {
249 for (j = 0; j < 5; ++j) {
250 st->d->v[i][j] = 0.0;
251 }
252 }
253 }
254
ebur128_init_channel_map(ebur128_state * st)255 static int ebur128_init_channel_map(ebur128_state* st) {
256 size_t i;
257 st->d->channel_map = (int*) malloc(st->channels * sizeof(int));
258 if (!st->d->channel_map) {
259 return EBUR128_ERROR_NOMEM;
260 }
261 if (st->channels == 4) {
262 st->d->channel_map[0] = EBUR128_LEFT;
263 st->d->channel_map[1] = EBUR128_RIGHT;
264 st->d->channel_map[2] = EBUR128_LEFT_SURROUND;
265 st->d->channel_map[3] = EBUR128_RIGHT_SURROUND;
266 } else if (st->channels == 5) {
267 st->d->channel_map[0] = EBUR128_LEFT;
268 st->d->channel_map[1] = EBUR128_RIGHT;
269 st->d->channel_map[2] = EBUR128_CENTER;
270 st->d->channel_map[3] = EBUR128_LEFT_SURROUND;
271 st->d->channel_map[4] = EBUR128_RIGHT_SURROUND;
272 } else {
273 for (i = 0; i < st->channels; ++i) {
274 switch (i) {
275 case 0: st->d->channel_map[i] = EBUR128_LEFT; break;
276 case 1: st->d->channel_map[i] = EBUR128_RIGHT; break;
277 case 2: st->d->channel_map[i] = EBUR128_CENTER; break;
278 case 3: st->d->channel_map[i] = EBUR128_UNUSED; break;
279 case 4: st->d->channel_map[i] = EBUR128_LEFT_SURROUND; break;
280 case 5: st->d->channel_map[i] = EBUR128_RIGHT_SURROUND; break;
281 default: st->d->channel_map[i] = EBUR128_UNUSED; break;
282 }
283 }
284 }
285 return EBUR128_SUCCESS;
286 }
287
ebur128_init_resampler(ebur128_state * st)288 static int ebur128_init_resampler(ebur128_state* st) {
289 int errcode = EBUR128_SUCCESS;
290
291 if (st->samplerate < 96000) {
292 st->d->interp = interp_create(49, 4, st->channels);
293 CHECK_ERROR(!st->d->interp, EBUR128_ERROR_NOMEM, exit)
294 } else if (st->samplerate < 192000) {
295 st->d->interp = interp_create(49, 2, st->channels);
296 CHECK_ERROR(!st->d->interp, EBUR128_ERROR_NOMEM, exit)
297 } else {
298 st->d->resampler_buffer_input = NULL;
299 st->d->resampler_buffer_output = NULL;
300 st->d->interp = NULL;
301 goto exit;
302 }
303
304 st->d->resampler_buffer_input_frames = st->d->samples_in_100ms * 4;
305 st->d->resampler_buffer_input = malloc(st->d->resampler_buffer_input_frames *
306 st->channels *
307 sizeof(float));
308 CHECK_ERROR(!st->d->resampler_buffer_input, EBUR128_ERROR_NOMEM, free_interp)
309
310 st->d->resampler_buffer_output_frames =
311 st->d->resampler_buffer_input_frames *
312 st->d->interp->factor;
313 st->d->resampler_buffer_output = malloc
314 (st->d->resampler_buffer_output_frames *
315 st->channels *
316 sizeof(float));
317 CHECK_ERROR(!st->d->resampler_buffer_output, EBUR128_ERROR_NOMEM, free_input)
318
319 return errcode;
320
321 free_interp:
322 interp_destroy(st->d->interp);
323 st->d->interp = NULL;
324 free_input:
325 free(st->d->resampler_buffer_input);
326 st->d->resampler_buffer_input = NULL;
327 exit:
328 return errcode;
329 }
330
ebur128_destroy_resampler(ebur128_state * st)331 static void ebur128_destroy_resampler(ebur128_state* st) {
332 free(st->d->resampler_buffer_input);
333 st->d->resampler_buffer_input = NULL;
334 free(st->d->resampler_buffer_output);
335 st->d->resampler_buffer_output = NULL;
336 interp_destroy(st->d->interp);
337 st->d->interp = NULL;
338 }
339
ebur128_get_version(int * major,int * minor,int * patch)340 void ebur128_get_version(int* major, int* minor, int* patch) {
341 *major = EBUR128_VERSION_MAJOR;
342 *minor = EBUR128_VERSION_MINOR;
343 *patch = EBUR128_VERSION_PATCH;
344 }
345
ebur128_init(unsigned int channels,unsigned long samplerate,int mode)346 ebur128_state* ebur128_init(unsigned int channels,
347 unsigned long samplerate,
348 int mode) {
349 int result;
350 int errcode;
351 ebur128_state* st;
352 unsigned int i;
353 size_t j;
354
355 if (channels == 0 || samplerate < 5) {
356 return NULL;
357 }
358
359 st = (ebur128_state*) malloc(sizeof(ebur128_state));
360 CHECK_ERROR(!st, 0, exit)
361 st->d = (struct ebur128_state_internal*)
362 malloc(sizeof(struct ebur128_state_internal));
363 CHECK_ERROR(!st->d, 0, free_state)
364 st->channels = channels;
365 errcode = ebur128_init_channel_map(st);
366 CHECK_ERROR(errcode, 0, free_internal)
367
368 st->d->sample_peak = (double*) malloc(channels * sizeof(double));
369 CHECK_ERROR(!st->d->sample_peak, 0, free_channel_map)
370 st->d->prev_sample_peak = (double*) malloc(channels * sizeof(double));
371 CHECK_ERROR(!st->d->prev_sample_peak, 0, free_sample_peak)
372 st->d->true_peak = (double*) malloc(channels * sizeof(double));
373 CHECK_ERROR(!st->d->true_peak, 0, free_prev_sample_peak)
374 st->d->prev_true_peak = (double*) malloc(channels * sizeof(double));
375 CHECK_ERROR(!st->d->prev_true_peak, 0, free_true_peak)
376 for (i = 0; i < channels; ++i) {
377 st->d->sample_peak[i] = 0.0;
378 st->d->prev_sample_peak[i] = 0.0;
379 st->d->true_peak[i] = 0.0;
380 st->d->prev_true_peak[i] = 0.0;
381 }
382
383 st->d->use_histogram = mode & EBUR128_MODE_HISTOGRAM ? 1 : 0;
384 st->d->history = ULONG_MAX;
385 st->samplerate = samplerate;
386 st->d->samples_in_100ms = (st->samplerate + 5) / 10;
387 st->mode = mode;
388 if ((mode & EBUR128_MODE_S) == EBUR128_MODE_S) {
389 st->d->window = 3000;
390 } else if ((mode & EBUR128_MODE_M) == EBUR128_MODE_M) {
391 st->d->window = 400;
392 } else {
393 goto free_prev_true_peak;
394 }
395 st->d->audio_data_frames = st->samplerate * st->d->window / 1000;
396 if (st->d->audio_data_frames % st->d->samples_in_100ms) {
397 /* round up to multiple of samples_in_100ms */
398 st->d->audio_data_frames = st->d->audio_data_frames
399 + st->d->samples_in_100ms
400 - (st->d->audio_data_frames % st->d->samples_in_100ms);
401 }
402 st->d->audio_data = (double*) malloc(st->d->audio_data_frames *
403 st->channels *
404 sizeof(double));
405 CHECK_ERROR(!st->d->audio_data, 0, free_true_peak)
406 for (j = 0; j < st->d->audio_data_frames * st->channels; ++j) {
407 st->d->audio_data[j] = 0.0;
408 }
409
410 ebur128_init_filter(st);
411
412 if (st->d->use_histogram) {
413 st->d->block_energy_histogram = malloc(1000 * sizeof(unsigned long));
414 CHECK_ERROR(!st->d->block_energy_histogram, 0, free_audio_data)
415 for (i = 0; i < 1000; ++i) {
416 st->d->block_energy_histogram[i] = 0;
417 }
418 } else {
419 st->d->block_energy_histogram = NULL;
420 }
421 if (st->d->use_histogram) {
422 st->d->short_term_block_energy_histogram = malloc(1000 * sizeof(unsigned long));
423 CHECK_ERROR(!st->d->short_term_block_energy_histogram, 0, free_block_energy_histogram)
424 for (i = 0; i < 1000; ++i) {
425 st->d->short_term_block_energy_histogram[i] = 0;
426 }
427 } else {
428 st->d->short_term_block_energy_histogram = NULL;
429 }
430 STAILQ_INIT(&st->d->block_list);
431 st->d->block_list_size = 0;
432 st->d->block_list_max = st->d->history / 100;
433 STAILQ_INIT(&st->d->short_term_block_list);
434 st->d->st_block_list_size = 0;
435 st->d->st_block_list_max = st->d->history / 3000;
436 st->d->short_term_frame_counter = 0;
437
438 result = ebur128_init_resampler(st);
439 CHECK_ERROR(result, 0, free_short_term_block_energy_histogram)
440
441 /* the first block needs 400ms of audio data */
442 st->d->needed_frames = st->d->samples_in_100ms * 4;
443 /* start at the beginning of the buffer */
444 st->d->audio_data_index = 0;
445
446 /* initialize static constants */
447 relative_gate_factor = pow(10.0, relative_gate / 10.0);
448 minus_twenty_decibels = pow(10.0, -20.0 / 10.0);
449 histogram_energy_boundaries[0] = pow(10.0, (-70.0 + 0.691) / 10.0);
450 if (st->d->use_histogram) {
451 for (i = 0; i < 1000; ++i) {
452 histogram_energies[i] = pow(10.0, ((double) i / 10.0 - 69.95 + 0.691) / 10.0);
453 }
454 for (i = 1; i < 1001; ++i) {
455 histogram_energy_boundaries[i] = pow(10.0, ((double) i / 10.0 - 70.0 + 0.691) / 10.0);
456 }
457 }
458
459 return st;
460
461 free_short_term_block_energy_histogram:
462 free(st->d->short_term_block_energy_histogram);
463 free_block_energy_histogram:
464 free(st->d->block_energy_histogram);
465 free_audio_data:
466 free(st->d->audio_data);
467 free_prev_true_peak:
468 free(st->d->prev_true_peak);
469 free_true_peak:
470 free(st->d->true_peak);
471 free_prev_sample_peak:
472 free(st->d->prev_sample_peak);
473 free_sample_peak:
474 free(st->d->sample_peak);
475 free_channel_map:
476 free(st->d->channel_map);
477 free_internal:
478 free(st->d);
479 free_state:
480 free(st);
481 exit:
482 return NULL;
483 }
484
ebur128_destroy(ebur128_state ** st)485 void ebur128_destroy(ebur128_state** st) {
486 struct ebur128_dq_entry* entry;
487 free((*st)->d->block_energy_histogram);
488 free((*st)->d->short_term_block_energy_histogram);
489 free((*st)->d->audio_data);
490 free((*st)->d->channel_map);
491 free((*st)->d->sample_peak);
492 free((*st)->d->prev_sample_peak);
493 free((*st)->d->true_peak);
494 free((*st)->d->prev_true_peak);
495 while (!STAILQ_EMPTY(&(*st)->d->block_list)) {
496 entry = STAILQ_FIRST(&(*st)->d->block_list);
497 STAILQ_REMOVE_HEAD(&(*st)->d->block_list, entries);
498 free(entry);
499 }
500 while (!STAILQ_EMPTY(&(*st)->d->short_term_block_list)) {
501 entry = STAILQ_FIRST(&(*st)->d->short_term_block_list);
502 STAILQ_REMOVE_HEAD(&(*st)->d->short_term_block_list, entries);
503 free(entry);
504 }
505 ebur128_destroy_resampler(*st);
506 free((*st)->d);
507 free(*st);
508 *st = NULL;
509 }
510
ebur128_check_true_peak(ebur128_state * st,size_t frames)511 static void ebur128_check_true_peak(ebur128_state* st, size_t frames) {
512 size_t c, i, frames_out;
513
514 frames_out = interp_process(st->d->interp, frames,
515 st->d->resampler_buffer_input,
516 st->d->resampler_buffer_output);
517
518 for (i = 0; i < frames_out; ++i) {
519 for (c = 0; c < st->channels; ++c) {
520 float val = st->d->resampler_buffer_output[i * st->channels + c];
521
522 if (val > st->d->prev_true_peak[c]) {
523 st->d->prev_true_peak[c] = val;
524 } else if (-val > st->d->prev_true_peak[c]) {
525 st->d->prev_true_peak[c] = -val;
526 }
527 }
528 }
529 }
530
531 #if defined(__SSE2_MATH__) || defined(_M_X64) || _M_IX86_FP >= 2
532 #include <xmmintrin.h>
533 #define TURN_ON_FTZ \
534 unsigned int mxcsr = _mm_getcsr(); \
535 _mm_setcsr(mxcsr | _MM_FLUSH_ZERO_ON);
536 #define TURN_OFF_FTZ _mm_setcsr(mxcsr);
537 #define FLUSH_MANUALLY
538 #else
539 #warning "manual FTZ is being used, please enable SSE2 (-msse2 -mfpmath=sse)"
540 #define TURN_ON_FTZ
541 #define TURN_OFF_FTZ
542 #define FLUSH_MANUALLY \
543 st->d->v[ci][4] = fabs(st->d->v[ci][4]) < DBL_MIN ? 0.0 : st->d->v[ci][4]; \
544 st->d->v[ci][3] = fabs(st->d->v[ci][3]) < DBL_MIN ? 0.0 : st->d->v[ci][3]; \
545 st->d->v[ci][2] = fabs(st->d->v[ci][2]) < DBL_MIN ? 0.0 : st->d->v[ci][2]; \
546 st->d->v[ci][1] = fabs(st->d->v[ci][1]) < DBL_MIN ? 0.0 : st->d->v[ci][1];
547 #endif
548
549 #define EBUR128_FILTER(type, min_scale, max_scale) \
550 static void ebur128_filter_##type(ebur128_state* st, const type* src, \
551 size_t frames) { \
552 static double scaling_factor = \
553 -((double) (min_scale)) > (double) (max_scale) ? \
554 -((double) (min_scale)) : (double) (max_scale); \
555 double* audio_data = st->d->audio_data + st->d->audio_data_index; \
556 size_t i, c; \
557 \
558 TURN_ON_FTZ \
559 \
560 if ((st->mode & EBUR128_MODE_SAMPLE_PEAK) == EBUR128_MODE_SAMPLE_PEAK) { \
561 for (c = 0; c < st->channels; ++c) { \
562 double max = 0.0; \
563 for (i = 0; i < frames; ++i) { \
564 if (src[i * st->channels + c] > max) { \
565 max = src[i * st->channels + c]; \
566 } else if (-src[i * st->channels + c] > max) { \
567 max = -1.0 * src[i * st->channels + c]; \
568 } \
569 } \
570 max /= scaling_factor; \
571 if (max > st->d->prev_sample_peak[c]) st->d->prev_sample_peak[c] = max; \
572 } \
573 } \
574 if ((st->mode & EBUR128_MODE_TRUE_PEAK) == EBUR128_MODE_TRUE_PEAK && \
575 st->d->interp) { \
576 for (c = 0; c < st->channels; ++c) { \
577 for (i = 0; i < frames; ++i) { \
578 st->d->resampler_buffer_input[i * st->channels + c] = \
579 (float) (src[i * st->channels + c] / scaling_factor); \
580 } \
581 } \
582 ebur128_check_true_peak(st, frames); \
583 } \
584 for (c = 0; c < st->channels; ++c) { \
585 int ci = st->d->channel_map[c] - 1; \
586 if (ci < 0) continue; \
587 else if (ci == EBUR128_DUAL_MONO - 1) ci = 0; /*dual mono */ \
588 for (i = 0; i < frames; ++i) { \
589 st->d->v[ci][0] = (double) (src[i * st->channels + c] / scaling_factor) \
590 - st->d->a[1] * st->d->v[ci][1] \
591 - st->d->a[2] * st->d->v[ci][2] \
592 - st->d->a[3] * st->d->v[ci][3] \
593 - st->d->a[4] * st->d->v[ci][4]; \
594 audio_data[i * st->channels + c] = \
595 st->d->b[0] * st->d->v[ci][0] \
596 + st->d->b[1] * st->d->v[ci][1] \
597 + st->d->b[2] * st->d->v[ci][2] \
598 + st->d->b[3] * st->d->v[ci][3] \
599 + st->d->b[4] * st->d->v[ci][4]; \
600 st->d->v[ci][4] = st->d->v[ci][3]; \
601 st->d->v[ci][3] = st->d->v[ci][2]; \
602 st->d->v[ci][2] = st->d->v[ci][1]; \
603 st->d->v[ci][1] = st->d->v[ci][0]; \
604 } \
605 FLUSH_MANUALLY \
606 } \
607 TURN_OFF_FTZ \
608 }
EBUR128_FILTER(short,SHRT_MIN,SHRT_MAX)609 EBUR128_FILTER(short, SHRT_MIN, SHRT_MAX)
610 EBUR128_FILTER(int, INT_MIN, INT_MAX)
611 EBUR128_FILTER(float, -1.0f, 1.0f)
612 EBUR128_FILTER(double, -1.0, 1.0)
613
614 static double ebur128_energy_to_loudness(double energy) {
615 return 10 * (log(energy) / log(10.0)) - 0.691;
616 }
617
find_histogram_index(double energy)618 static size_t find_histogram_index(double energy) {
619 size_t index_min = 0;
620 size_t index_max = 1000;
621 size_t index_mid;
622
623 do {
624 index_mid = (index_min + index_max) / 2;
625 if (energy >= histogram_energy_boundaries[index_mid]) {
626 index_min = index_mid;
627 } else {
628 index_max = index_mid;
629 }
630 } while (index_max - index_min != 1);
631
632 return index_min;
633 }
634
ebur128_calc_gating_block(ebur128_state * st,size_t frames_per_block,double * optional_output)635 static int ebur128_calc_gating_block(ebur128_state* st, size_t frames_per_block,
636 double* optional_output) {
637 size_t i, c;
638 double sum = 0.0;
639 double channel_sum;
640 for (c = 0; c < st->channels; ++c) {
641 if (st->d->channel_map[c] == EBUR128_UNUSED) {
642 continue;
643 }
644 channel_sum = 0.0;
645 if (st->d->audio_data_index < frames_per_block * st->channels) {
646 for (i = 0; i < st->d->audio_data_index / st->channels; ++i) {
647 channel_sum += st->d->audio_data[i * st->channels + c] *
648 st->d->audio_data[i * st->channels + c];
649 }
650 for (i = st->d->audio_data_frames -
651 (frames_per_block -
652 st->d->audio_data_index / st->channels);
653 i < st->d->audio_data_frames; ++i) {
654 channel_sum += st->d->audio_data[i * st->channels + c] *
655 st->d->audio_data[i * st->channels + c];
656 }
657 } else {
658 for (i = st->d->audio_data_index / st->channels - frames_per_block;
659 i < st->d->audio_data_index / st->channels;
660 ++i) {
661 channel_sum += st->d->audio_data[i * st->channels + c] *
662 st->d->audio_data[i * st->channels + c];
663 }
664 }
665 if (st->d->channel_map[c] == EBUR128_Mp110 ||
666 st->d->channel_map[c] == EBUR128_Mm110 ||
667 st->d->channel_map[c] == EBUR128_Mp060 ||
668 st->d->channel_map[c] == EBUR128_Mm060 ||
669 st->d->channel_map[c] == EBUR128_Mp090 ||
670 st->d->channel_map[c] == EBUR128_Mm090) {
671 channel_sum *= 1.41;
672 } else if (st->d->channel_map[c] == EBUR128_DUAL_MONO) {
673 channel_sum *= 2.0;
674 }
675 sum += channel_sum;
676 }
677 sum /= (double) frames_per_block;
678 if (optional_output) {
679 *optional_output = sum;
680 return EBUR128_SUCCESS;
681 } else if (sum >= histogram_energy_boundaries[0]) {
682 if (st->d->use_histogram) {
683 ++st->d->block_energy_histogram[find_histogram_index(sum)];
684 } else {
685 struct ebur128_dq_entry* block;
686 if (st->d->block_list_size == st->d->block_list_max) {
687 block = STAILQ_FIRST(&st->d->block_list);
688 STAILQ_REMOVE_HEAD(&st->d->block_list, entries);
689 } else {
690 block = (struct ebur128_dq_entry*) malloc(sizeof(struct ebur128_dq_entry));
691 if (!block) {
692 return EBUR128_ERROR_NOMEM;
693 }
694 st->d->block_list_size++;
695 }
696 block->z = sum;
697 STAILQ_INSERT_TAIL(&st->d->block_list, block, entries);
698 }
699 return EBUR128_SUCCESS;
700 } else {
701 return EBUR128_SUCCESS;
702 }
703 }
704
ebur128_set_channel(ebur128_state * st,unsigned int channel_number,int value)705 int ebur128_set_channel(ebur128_state* st,
706 unsigned int channel_number,
707 int value) {
708 if (channel_number >= st->channels) {
709 return 1;
710 }
711 if (value == EBUR128_DUAL_MONO &&
712 (st->channels != 1 || channel_number != 0)) {
713 fprintf(stderr, "EBUR128_DUAL_MONO only works with mono files!\n");
714 return 1;
715 }
716 st->d->channel_map[channel_number] = value;
717 return 0;
718 }
719
ebur128_change_parameters(ebur128_state * st,unsigned int channels,unsigned long samplerate)720 int ebur128_change_parameters(ebur128_state* st,
721 unsigned int channels,
722 unsigned long samplerate) {
723 int errcode = EBUR128_SUCCESS;
724 size_t j;
725
726 if (channels == 0 || samplerate < 5) {
727 return EBUR128_ERROR_NOMEM;
728 }
729
730 if (channels == st->channels &&
731 samplerate == st->samplerate) {
732 return EBUR128_ERROR_NO_CHANGE;
733 }
734
735 free(st->d->audio_data);
736 st->d->audio_data = NULL;
737
738 if (channels != st->channels) {
739 unsigned int i;
740
741 free(st->d->channel_map); st->d->channel_map = NULL;
742 free(st->d->sample_peak); st->d->sample_peak = NULL;
743 free(st->d->prev_sample_peak); st->d->prev_sample_peak = NULL;
744 free(st->d->true_peak); st->d->true_peak = NULL;
745 free(st->d->prev_true_peak); st->d->prev_true_peak = NULL;
746 st->channels = channels;
747
748 errcode = ebur128_init_channel_map(st);
749 CHECK_ERROR(errcode, EBUR128_ERROR_NOMEM, exit)
750
751 st->d->sample_peak = (double*) malloc(channels * sizeof(double));
752 CHECK_ERROR(!st->d->sample_peak, EBUR128_ERROR_NOMEM, exit)
753 st->d->prev_sample_peak = (double*) malloc(channels * sizeof(double));
754 CHECK_ERROR(!st->d->prev_sample_peak, EBUR128_ERROR_NOMEM, exit)
755 st->d->true_peak = (double*) malloc(channels * sizeof(double));
756 CHECK_ERROR(!st->d->true_peak, EBUR128_ERROR_NOMEM, exit)
757 st->d->prev_true_peak = (double*) malloc(channels * sizeof(double));
758 CHECK_ERROR(!st->d->prev_true_peak, EBUR128_ERROR_NOMEM, exit)
759 for (i = 0; i < channels; ++i) {
760 st->d->sample_peak[i] = 0.0;
761 st->d->prev_sample_peak[i] = 0.0;
762 st->d->true_peak[i] = 0.0;
763 st->d->prev_true_peak[i] = 0.0;
764 }
765 }
766 if (samplerate != st->samplerate) {
767 st->samplerate = samplerate;
768 st->d->samples_in_100ms = (st->samplerate + 5) / 10;
769 ebur128_init_filter(st);
770 }
771 st->d->audio_data_frames = st->samplerate * st->d->window / 1000;
772 if (st->d->audio_data_frames % st->d->samples_in_100ms) {
773 /* round up to multiple of samples_in_100ms */
774 st->d->audio_data_frames = st->d->audio_data_frames
775 + st->d->samples_in_100ms
776 - (st->d->audio_data_frames % st->d->samples_in_100ms);
777 }
778 st->d->audio_data = (double*) malloc(st->d->audio_data_frames *
779 st->channels *
780 sizeof(double));
781 CHECK_ERROR(!st->d->audio_data, EBUR128_ERROR_NOMEM, exit)
782 for (j = 0; j < st->d->audio_data_frames * st->channels; ++j) {
783 st->d->audio_data[j] = 0.0;
784 }
785
786 ebur128_destroy_resampler(st);
787 errcode = ebur128_init_resampler(st);
788 CHECK_ERROR(errcode, EBUR128_ERROR_NOMEM, exit)
789
790 /* the first block needs 400ms of audio data */
791 st->d->needed_frames = st->d->samples_in_100ms * 4;
792 /* start at the beginning of the buffer */
793 st->d->audio_data_index = 0;
794 /* reset short term frame counter */
795 st->d->short_term_frame_counter = 0;
796
797 exit:
798 return errcode;
799 }
800
ebur128_set_max_window(ebur128_state * st,unsigned long window)801 int ebur128_set_max_window(ebur128_state* st, unsigned long window)
802 {
803 int errcode = EBUR128_SUCCESS;
804 size_t j;
805
806 if ((st->mode & EBUR128_MODE_S) == EBUR128_MODE_S && window < 3000) {
807 window = 3000;
808 } else if ((st->mode & EBUR128_MODE_M) == EBUR128_MODE_M && window < 400) {
809 window = 400;
810 }
811 if (window == st->d->window) {
812 return EBUR128_ERROR_NO_CHANGE;
813 }
814
815 st->d->window = window;
816 free(st->d->audio_data);
817 st->d->audio_data = NULL;
818 st->d->audio_data_frames = st->samplerate * st->d->window / 1000;
819 if (st->d->audio_data_frames % st->d->samples_in_100ms) {
820 /* round up to multiple of samples_in_100ms */
821 st->d->audio_data_frames = st->d->audio_data_frames
822 + st->d->samples_in_100ms
823 - (st->d->audio_data_frames % st->d->samples_in_100ms);
824 }
825 st->d->audio_data = (double*) malloc(st->d->audio_data_frames *
826 st->channels *
827 sizeof(double));
828 CHECK_ERROR(!st->d->audio_data, EBUR128_ERROR_NOMEM, exit)
829 for (j = 0; j < st->d->audio_data_frames * st->channels; ++j) {
830 st->d->audio_data[j] = 0.0;
831 }
832
833 /* the first block needs 400ms of audio data */
834 st->d->needed_frames = st->d->samples_in_100ms * 4;
835 /* start at the beginning of the buffer */
836 st->d->audio_data_index = 0;
837 /* reset short term frame counter */
838 st->d->short_term_frame_counter = 0;
839
840 exit:
841 return errcode;
842 }
843
ebur128_set_max_history(ebur128_state * st,unsigned long history)844 int ebur128_set_max_history(ebur128_state* st, unsigned long history)
845 {
846 if ((st->mode & EBUR128_MODE_LRA) == EBUR128_MODE_LRA && history < 3000) {
847 history = 3000;
848 } else if ((st->mode & EBUR128_MODE_M) == EBUR128_MODE_M && history < 400) {
849 history = 400;
850 }
851 if (history == st->d->history) {
852 return EBUR128_ERROR_NO_CHANGE;
853 }
854 st->d->history = history;
855 st->d->block_list_max = st->d->history / 100;
856 st->d->st_block_list_max = st->d->history / 3000;
857 while (st->d->block_list_size > st->d->block_list_max) {
858 struct ebur128_dq_entry* block = STAILQ_FIRST(&st->d->block_list);
859 STAILQ_REMOVE_HEAD(&st->d->block_list, entries);
860 free(block);
861 st->d->block_list_size--;
862 }
863 while (st->d->st_block_list_size > st->d->st_block_list_max) {
864 struct ebur128_dq_entry* block = STAILQ_FIRST(&st->d->short_term_block_list);
865 STAILQ_REMOVE_HEAD(&st->d->short_term_block_list, entries);
866 free(block);
867 st->d->st_block_list_size--;
868 }
869 return EBUR128_SUCCESS;
870 }
871
872 static int ebur128_energy_shortterm(ebur128_state* st, double* out);
873 #define EBUR128_ADD_FRAMES(type) \
874 int ebur128_add_frames_##type(ebur128_state* st, \
875 const type* src, size_t frames) { \
876 size_t src_index = 0; \
877 unsigned int c = 0; \
878 for (c = 0; c < st->channels; c++) { \
879 st->d->prev_sample_peak[c] = 0.0; \
880 st->d->prev_true_peak[c] = 0.0; \
881 } \
882 while (frames > 0) { \
883 if (frames >= st->d->needed_frames) { \
884 ebur128_filter_##type(st, src + src_index, st->d->needed_frames); \
885 src_index += st->d->needed_frames * st->channels; \
886 frames -= st->d->needed_frames; \
887 st->d->audio_data_index += st->d->needed_frames * st->channels; \
888 /* calculate the new gating block */ \
889 if ((st->mode & EBUR128_MODE_I) == EBUR128_MODE_I) { \
890 if (ebur128_calc_gating_block(st, st->d->samples_in_100ms * 4, NULL)) {\
891 return EBUR128_ERROR_NOMEM; \
892 } \
893 } \
894 if ((st->mode & EBUR128_MODE_LRA) == EBUR128_MODE_LRA) { \
895 st->d->short_term_frame_counter += st->d->needed_frames; \
896 if (st->d->short_term_frame_counter == st->d->samples_in_100ms * 30) { \
897 struct ebur128_dq_entry* block; \
898 double st_energy; \
899 if (ebur128_energy_shortterm(st, &st_energy) == EBUR128_SUCCESS && \
900 st_energy >= histogram_energy_boundaries[0]) { \
901 if (st->d->use_histogram) { \
902 ++st->d->short_term_block_energy_histogram[ \
903 find_histogram_index(st_energy)];\
904 } else { \
905 if (st->d->st_block_list_size == st->d->st_block_list_max) { \
906 block = STAILQ_FIRST(&st->d->short_term_block_list); \
907 STAILQ_REMOVE_HEAD(&st->d->short_term_block_list, entries); \
908 } else { \
909 block = (struct ebur128_dq_entry*) \
910 malloc(sizeof(struct ebur128_dq_entry)); \
911 if (!block) return EBUR128_ERROR_NOMEM; \
912 st->d->st_block_list_size++; \
913 } \
914 block->z = st_energy; \
915 STAILQ_INSERT_TAIL(&st->d->short_term_block_list, \
916 block, entries); \
917 } \
918 } \
919 st->d->short_term_frame_counter = st->d->samples_in_100ms * 20; \
920 } \
921 } \
922 /* 100ms are needed for all blocks besides the first one */ \
923 st->d->needed_frames = st->d->samples_in_100ms; \
924 /* reset audio_data_index when buffer full */ \
925 if (st->d->audio_data_index == st->d->audio_data_frames * st->channels) {\
926 st->d->audio_data_index = 0; \
927 } \
928 } else { \
929 ebur128_filter_##type(st, src + src_index, frames); \
930 st->d->audio_data_index += frames * st->channels; \
931 if ((st->mode & EBUR128_MODE_LRA) == EBUR128_MODE_LRA) { \
932 st->d->short_term_frame_counter += frames; \
933 } \
934 st->d->needed_frames -= frames; \
935 frames = 0; \
936 } \
937 } \
938 for (c = 0; c < st->channels; c++) { \
939 if (st->d->prev_sample_peak[c] > st->d->sample_peak[c]) { \
940 st->d->sample_peak[c] = st->d->prev_sample_peak[c]; \
941 } \
942 if (st->d->prev_true_peak[c] > st->d->true_peak[c]) { \
943 st->d->true_peak[c] = st->d->prev_true_peak[c]; \
944 } \
945 } \
946 return EBUR128_SUCCESS; \
947 }
948 EBUR128_ADD_FRAMES(short)
EBUR128_ADD_FRAMES(int)949 EBUR128_ADD_FRAMES(int)
950 EBUR128_ADD_FRAMES(float)
951 EBUR128_ADD_FRAMES(double)
952
953 static int ebur128_calc_relative_threshold(ebur128_state* st,
954 size_t* above_thresh_counter,
955 double* relative_threshold) {
956 struct ebur128_dq_entry* it;
957 size_t i;
958
959 if (st->d->use_histogram) {
960 for (i = 0; i < 1000; ++i) {
961 *relative_threshold += st->d->block_energy_histogram[i] *
962 histogram_energies[i];
963 *above_thresh_counter += st->d->block_energy_histogram[i];
964 }
965 } else {
966 STAILQ_FOREACH(it, &st->d->block_list, entries) {
967 ++*above_thresh_counter;
968 *relative_threshold += it->z;
969 }
970 }
971
972 return EBUR128_SUCCESS;
973 }
974
ebur128_gated_loudness(ebur128_state ** sts,size_t size,double * out)975 static int ebur128_gated_loudness(ebur128_state** sts, size_t size,
976 double* out) {
977 struct ebur128_dq_entry* it;
978 double gated_loudness = 0.0;
979 double relative_threshold = 0.0;
980 size_t above_thresh_counter = 0;
981 size_t i, j, start_index;
982
983 for (i = 0; i < size; i++) {
984 if (sts[i] && (sts[i]->mode & EBUR128_MODE_I) != EBUR128_MODE_I) {
985 return EBUR128_ERROR_INVALID_MODE;
986 }
987 }
988
989 for (i = 0; i < size; i++) {
990 if (!sts[i]) {
991 continue;
992 }
993 ebur128_calc_relative_threshold(sts[i], &above_thresh_counter, &relative_threshold);
994 }
995 if (!above_thresh_counter) {
996 *out = -HUGE_VAL;
997 return EBUR128_SUCCESS;
998 }
999
1000 relative_threshold /= (double)above_thresh_counter;
1001 relative_threshold *= relative_gate_factor;
1002
1003 above_thresh_counter = 0;
1004 if (relative_threshold < histogram_energy_boundaries[0]) {
1005 start_index = 0;
1006 } else {
1007 start_index = find_histogram_index(relative_threshold);
1008 if (relative_threshold > histogram_energies[start_index]) {
1009 ++start_index;
1010 }
1011 }
1012 for (i = 0; i < size; i++) {
1013 if (!sts[i]) {
1014 continue;
1015 }
1016 if (sts[i]->d->use_histogram) {
1017 for (j = start_index; j < 1000; ++j) {
1018 gated_loudness += sts[i]->d->block_energy_histogram[j] *
1019 histogram_energies[j];
1020 above_thresh_counter += sts[i]->d->block_energy_histogram[j];
1021 }
1022 } else {
1023 STAILQ_FOREACH(it, &sts[i]->d->block_list, entries) {
1024 if (it->z >= relative_threshold) {
1025 ++above_thresh_counter;
1026 gated_loudness += it->z;
1027 }
1028 }
1029 }
1030 }
1031 if (!above_thresh_counter) {
1032 *out = -HUGE_VAL;
1033 return EBUR128_SUCCESS;
1034 }
1035 gated_loudness /= (double) above_thresh_counter;
1036 *out = ebur128_energy_to_loudness(gated_loudness);
1037 return EBUR128_SUCCESS;
1038 }
1039
ebur128_relative_threshold(ebur128_state * st,double * out)1040 int ebur128_relative_threshold(ebur128_state* st, double* out) {
1041 double relative_threshold = 0.0;
1042 size_t above_thresh_counter = 0;
1043
1044 if ((st->mode & EBUR128_MODE_I) != EBUR128_MODE_I) {
1045 return EBUR128_ERROR_INVALID_MODE;
1046 }
1047
1048 ebur128_calc_relative_threshold(st, &above_thresh_counter, &relative_threshold);
1049
1050 if (!above_thresh_counter) {
1051 *out = -70.0;
1052 return EBUR128_SUCCESS;
1053 }
1054
1055 relative_threshold /= (double)above_thresh_counter;
1056 relative_threshold *= relative_gate_factor;
1057
1058 *out = ebur128_energy_to_loudness(relative_threshold);
1059 return EBUR128_SUCCESS;
1060 }
1061
ebur128_loudness_global(ebur128_state * st,double * out)1062 int ebur128_loudness_global(ebur128_state* st, double* out) {
1063 return ebur128_gated_loudness(&st, 1, out);
1064 }
1065
ebur128_loudness_global_multiple(ebur128_state ** sts,size_t size,double * out)1066 int ebur128_loudness_global_multiple(ebur128_state** sts, size_t size,
1067 double* out) {
1068 return ebur128_gated_loudness(sts, size, out);
1069 }
1070
ebur128_energy_in_interval(ebur128_state * st,size_t interval_frames,double * out)1071 static int ebur128_energy_in_interval(ebur128_state* st,
1072 size_t interval_frames,
1073 double* out) {
1074 if (interval_frames > st->d->audio_data_frames) {
1075 return EBUR128_ERROR_INVALID_MODE;
1076 }
1077 ebur128_calc_gating_block(st, interval_frames, out);
1078 return EBUR128_SUCCESS;
1079 }
1080
ebur128_energy_shortterm(ebur128_state * st,double * out)1081 static int ebur128_energy_shortterm(ebur128_state* st, double* out) {
1082 return ebur128_energy_in_interval(st, st->d->samples_in_100ms * 30, out);
1083 }
1084
ebur128_loudness_momentary(ebur128_state * st,double * out)1085 int ebur128_loudness_momentary(ebur128_state* st, double* out) {
1086 double energy;
1087 int error = ebur128_energy_in_interval(st, st->d->samples_in_100ms * 4,
1088 &energy);
1089 if (error) {
1090 return error;
1091 } else if (energy <= 0.0) {
1092 *out = -HUGE_VAL;
1093 return EBUR128_SUCCESS;
1094 }
1095 *out = ebur128_energy_to_loudness(energy);
1096 return EBUR128_SUCCESS;
1097 }
1098
ebur128_loudness_shortterm(ebur128_state * st,double * out)1099 int ebur128_loudness_shortterm(ebur128_state* st, double* out) {
1100 double energy;
1101 int error = ebur128_energy_shortterm(st, &energy);
1102 if (error) {
1103 return error;
1104 } else if (energy <= 0.0) {
1105 *out = -HUGE_VAL;
1106 return EBUR128_SUCCESS;
1107 }
1108 *out = ebur128_energy_to_loudness(energy);
1109 return EBUR128_SUCCESS;
1110 }
1111
ebur128_loudness_window(ebur128_state * st,unsigned long window,double * out)1112 int ebur128_loudness_window(ebur128_state* st,
1113 unsigned long window,
1114 double* out) {
1115 double energy;
1116 size_t interval_frames = st->samplerate * window / 1000;
1117 int error = ebur128_energy_in_interval(st, interval_frames, &energy);
1118 if (error) {
1119 return error;
1120 } else if (energy <= 0.0) {
1121 *out = -HUGE_VAL;
1122 return EBUR128_SUCCESS;
1123 }
1124 *out = ebur128_energy_to_loudness(energy);
1125 return EBUR128_SUCCESS;
1126 }
1127
ebur128_double_cmp(const void * p1,const void * p2)1128 static int ebur128_double_cmp(const void *p1, const void *p2) {
1129 const double* d1 = (const double*) p1;
1130 const double* d2 = (const double*) p2;
1131 return (*d1 > *d2) - (*d1 < *d2);
1132 }
1133
1134 /* EBU - TECH 3342 */
ebur128_loudness_range_multiple(ebur128_state ** sts,size_t size,double * out)1135 int ebur128_loudness_range_multiple(ebur128_state** sts, size_t size,
1136 double* out) {
1137 size_t i, j;
1138 struct ebur128_dq_entry* it;
1139 double* stl_vector;
1140 size_t stl_size;
1141 double* stl_relgated;
1142 size_t stl_relgated_size;
1143 double stl_power, stl_integrated;
1144 /* High and low percentile energy */
1145 double h_en, l_en;
1146 int use_histogram = 0;
1147
1148 for (i = 0; i < size; ++i) {
1149 if (sts[i]) {
1150 if ((sts[i]->mode & EBUR128_MODE_LRA) != EBUR128_MODE_LRA) {
1151 return EBUR128_ERROR_INVALID_MODE;
1152 }
1153 if (i == 0 && sts[i]->mode & EBUR128_MODE_HISTOGRAM) {
1154 use_histogram = 1;
1155 } else if (use_histogram != !!(sts[i]->mode & EBUR128_MODE_HISTOGRAM)) {
1156 return EBUR128_ERROR_INVALID_MODE;
1157 }
1158 }
1159 }
1160
1161 if (use_histogram) {
1162 unsigned long hist[1000] = { 0 };
1163 size_t percentile_low, percentile_high;
1164 size_t index;
1165
1166 stl_size = 0;
1167 stl_power = 0.0;
1168 for (i = 0; i < size; ++i) {
1169 if (!sts[i]) {
1170 continue;
1171 }
1172 for (j = 0; j < 1000; ++j) {
1173 hist[j] += sts[i]->d->short_term_block_energy_histogram[j];
1174 stl_size += sts[i]->d->short_term_block_energy_histogram[j];
1175 stl_power += sts[i]->d->short_term_block_energy_histogram[j]
1176 * histogram_energies[j];
1177 }
1178 }
1179 if (!stl_size) {
1180 *out = 0.0;
1181 return EBUR128_SUCCESS;
1182 }
1183
1184 stl_power /= stl_size;
1185 stl_integrated = minus_twenty_decibels * stl_power;
1186
1187 if (stl_integrated < histogram_energy_boundaries[0]) {
1188 index = 0;
1189 } else {
1190 index = find_histogram_index(stl_integrated);
1191 if (stl_integrated > histogram_energies[index]) {
1192 ++index;
1193 }
1194 }
1195 stl_size = 0;
1196 for (j = index; j < 1000; ++j) {
1197 stl_size += hist[j];
1198 }
1199 if (!stl_size) {
1200 *out = 0.0;
1201 return EBUR128_SUCCESS;
1202 }
1203
1204 percentile_low = (size_t) ((stl_size - 1) * 0.1 + 0.5);
1205 percentile_high = (size_t) ((stl_size - 1) * 0.95 + 0.5);
1206
1207 stl_size = 0;
1208 j = index;
1209 while (stl_size <= percentile_low) {
1210 stl_size += hist[j++];
1211 }
1212 l_en = histogram_energies[j - 1];
1213 while (stl_size <= percentile_high) {
1214 stl_size += hist[j++];
1215 }
1216 h_en = histogram_energies[j - 1];
1217 *out = ebur128_energy_to_loudness(h_en) - ebur128_energy_to_loudness(l_en);
1218 return EBUR128_SUCCESS;
1219
1220 } else {
1221 stl_size = 0;
1222 for (i = 0; i < size; ++i) {
1223 if (!sts[i]) {
1224 continue;
1225 }
1226 STAILQ_FOREACH(it, &sts[i]->d->short_term_block_list, entries) {
1227 ++stl_size;
1228 }
1229 }
1230 if (!stl_size) {
1231 *out = 0.0;
1232 return EBUR128_SUCCESS;
1233 }
1234 stl_vector = (double*) malloc(stl_size * sizeof(double));
1235 if (!stl_vector) {
1236 return EBUR128_ERROR_NOMEM;
1237 }
1238
1239 j = 0;
1240 for (i = 0; i < size; ++i) {
1241 if (!sts[i]) {
1242 continue;
1243 }
1244 STAILQ_FOREACH(it, &sts[i]->d->short_term_block_list, entries) {
1245 stl_vector[j] = it->z;
1246 ++j;
1247 }
1248 }
1249 qsort(stl_vector, stl_size, sizeof(double), ebur128_double_cmp);
1250 stl_power = 0.0;
1251 for (i = 0; i < stl_size; ++i) {
1252 stl_power += stl_vector[i];
1253 }
1254 stl_power /= (double) stl_size;
1255 stl_integrated = minus_twenty_decibels * stl_power;
1256
1257 stl_relgated = stl_vector;
1258 stl_relgated_size = stl_size;
1259 while (stl_relgated_size > 0 && *stl_relgated < stl_integrated) {
1260 ++stl_relgated;
1261 --stl_relgated_size;
1262 }
1263
1264 if (stl_relgated_size) {
1265 h_en = stl_relgated[(size_t) ((stl_relgated_size - 1) * 0.95 + 0.5)];
1266 l_en = stl_relgated[(size_t) ((stl_relgated_size - 1) * 0.1 + 0.5)];
1267 free(stl_vector);
1268 *out = ebur128_energy_to_loudness(h_en) - ebur128_energy_to_loudness(l_en);
1269 return EBUR128_SUCCESS;
1270 } else {
1271 free(stl_vector);
1272 *out = 0.0;
1273 return EBUR128_SUCCESS;
1274 }
1275 }
1276 }
1277
ebur128_loudness_range(ebur128_state * st,double * out)1278 int ebur128_loudness_range(ebur128_state* st, double* out) {
1279 return ebur128_loudness_range_multiple(&st, 1, out);
1280 }
1281
ebur128_sample_peak(ebur128_state * st,unsigned int channel_number,double * out)1282 int ebur128_sample_peak(ebur128_state* st,
1283 unsigned int channel_number,
1284 double* out) {
1285 if ((st->mode & EBUR128_MODE_SAMPLE_PEAK) != EBUR128_MODE_SAMPLE_PEAK) {
1286 return EBUR128_ERROR_INVALID_MODE;
1287 } else if (channel_number >= st->channels) {
1288 return EBUR128_ERROR_INVALID_CHANNEL_INDEX;
1289 }
1290 *out = st->d->sample_peak[channel_number];
1291 return EBUR128_SUCCESS;
1292 }
1293
ebur128_prev_sample_peak(ebur128_state * st,unsigned int channel_number,double * out)1294 int ebur128_prev_sample_peak(ebur128_state* st,
1295 unsigned int channel_number,
1296 double* out) {
1297 if ((st->mode & EBUR128_MODE_SAMPLE_PEAK) != EBUR128_MODE_SAMPLE_PEAK) {
1298 return EBUR128_ERROR_INVALID_MODE;
1299 } else if (channel_number >= st->channels) {
1300 return EBUR128_ERROR_INVALID_CHANNEL_INDEX;
1301 }
1302 *out = st->d->prev_sample_peak[channel_number];
1303 return EBUR128_SUCCESS;
1304 }
1305
ebur128_true_peak(ebur128_state * st,unsigned int channel_number,double * out)1306 int ebur128_true_peak(ebur128_state* st,
1307 unsigned int channel_number,
1308 double* out) {
1309 if ((st->mode & EBUR128_MODE_TRUE_PEAK) != EBUR128_MODE_TRUE_PEAK) {
1310 return EBUR128_ERROR_INVALID_MODE;
1311 } else if (channel_number >= st->channels) {
1312 return EBUR128_ERROR_INVALID_CHANNEL_INDEX;
1313 }
1314 *out = st->d->true_peak[channel_number] > st->d->sample_peak[channel_number]
1315 ? st->d->true_peak[channel_number]
1316 : st->d->sample_peak[channel_number];
1317 return EBUR128_SUCCESS;
1318 }
1319
ebur128_prev_true_peak(ebur128_state * st,unsigned int channel_number,double * out)1320 int ebur128_prev_true_peak(ebur128_state* st,
1321 unsigned int channel_number,
1322 double* out) {
1323 if ((st->mode & EBUR128_MODE_TRUE_PEAK) != EBUR128_MODE_TRUE_PEAK) {
1324 return EBUR128_ERROR_INVALID_MODE;
1325 } else if (channel_number >= st->channels) {
1326 return EBUR128_ERROR_INVALID_CHANNEL_INDEX;
1327 }
1328 *out = st->d->prev_true_peak[channel_number]
1329 > st->d->prev_sample_peak[channel_number]
1330 ? st->d->prev_true_peak[channel_number]
1331 : st->d->prev_sample_peak[channel_number];
1332 return EBUR128_SUCCESS;
1333 }
1334