/*
  * Copyright (c) 2011 Jan Kokemüller
  *
  * This file is part of FFmpeg.
  *
  * FFmpeg is free software; you can redistribute it and/or
  * modify it under the terms of the GNU Lesser General Public
  * License as published by the Free Software Foundation; either
  * version 2.1 of the License, or (at your option) any later version.
  *
  * FFmpeg 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
  * Lesser General Public License for more details.
  *
  * You should have received a copy of the GNU Lesser General Public
  * License along with FFmpeg; if not, write to the Free Software
  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
  *
  * This file is based on libebur128 which is available at
  * https://github.com/jiixyj/libebur128/
  *
  * Libebur128 has the following copyright:
  *
  * Permission is hereby granted, free of charge, to any person obtaining a copy
  * of this software and associated documentation files (the "Software"), to deal
  * in the Software without restriction, including without limitation the rights
  * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
  * copies of the Software, and to permit persons to whom the Software is
  * furnished to do so, subject to the following conditions:
  *
  * The above copyright notice and this permission notice shall be included in
  * all copies or substantial portions of the Software.
  *
  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
  * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
  * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
  * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
  * THE SOFTWARE.
 */
 
 #include "ebur128.h"
 
 #include <float.h>
 #include <limits.h>
 #include <math.h>               /* You may have to define _USE_MATH_DEFINES if you use MSVC */
 
 #include "libavutil/common.h"
 #include "libavutil/mem.h"
 #include "libavutil/thread.h"
 
 #define CHECK_ERROR(condition, errorcode, goto_point)                          \
     if ((condition)) {                                                         \
         errcode = (errorcode);                                                 \
         goto goto_point;                                                       \
     }
 
 #define ALMOST_ZERO 0.000001
 
 #define RELATIVE_GATE         (-10.0)
 #define RELATIVE_GATE_FACTOR  pow(10.0, RELATIVE_GATE / 10.0)
 #define MINUS_20DB            pow(10.0, -20.0 / 10.0)
 
 struct FFEBUR128StateInternal {
     /** Filtered audio data (used as ring buffer). */
     double *audio_data;
     /** Size of audio_data array. */
     size_t audio_data_frames;
     /** Current index for audio_data. */
     size_t audio_data_index;
     /** How many frames are needed for a gating block. Will correspond to 400ms
      *  of audio at initialization, and 100ms after the first block (75% overlap
      *  as specified in the 2011 revision of BS1770). */
     unsigned long needed_frames;
     /** The channel map. Has as many elements as there are channels. */
     int *channel_map;
     /** How many samples fit in 100ms (rounded). */
     unsigned long samples_in_100ms;
     /** BS.1770 filter coefficients (nominator). */
     double b[5];
     /** BS.1770 filter coefficients (denominator). */
     double a[5];
     /** BS.1770 filter state. */
     double v[5][5];
     /** Histograms, used to calculate LRA. */
     unsigned long *block_energy_histogram;
     unsigned long *short_term_block_energy_histogram;
     /** Keeps track of when a new short term block is needed. */
     size_t short_term_frame_counter;
     /** Maximum sample peak, one per channel */
     double *sample_peak;
     /** The maximum window duration in ms. */
     unsigned long window;
     /** Data pointer array for interleaved data */
     void **data_ptrs;
 };
 
 static AVOnce histogram_init = AV_ONCE_INIT;
 static DECLARE_ALIGNED(32, double, histogram_energies)[1000];
 static DECLARE_ALIGNED(32, double, histogram_energy_boundaries)[1001];
 
 static void ebur128_init_filter(FFEBUR128State * st)
 {
     int i, j;
 
     double f0 = 1681.974450955533;
     double G = 3.999843853973347;
     double Q = 0.7071752369554196;
 
     double K = tan(M_PI * f0 / (double) st->samplerate);
     double Vh = pow(10.0, G / 20.0);
     double Vb = pow(Vh, 0.4996667741545416);
 
     double pb[3] = { 0.0, 0.0, 0.0 };
     double pa[3] = { 1.0, 0.0, 0.0 };
     double rb[3] = { 1.0, -2.0, 1.0 };
     double ra[3] = { 1.0, 0.0, 0.0 };
 
     double a0 = 1.0 + K / Q + K * K;
     pb[0] = (Vh + Vb * K / Q + K * K) / a0;
     pb[1] = 2.0 * (K * K - Vh) / a0;
     pb[2] = (Vh - Vb * K / Q + K * K) / a0;
     pa[1] = 2.0 * (K * K - 1.0) / a0;
     pa[2] = (1.0 - K / Q + K * K) / a0;
 
     f0 = 38.13547087602444;
     Q = 0.5003270373238773;
     K = tan(M_PI * f0 / (double) st->samplerate);
 
     ra[1] = 2.0 * (K * K - 1.0) / (1.0 + K / Q + K * K);
     ra[2] = (1.0 - K / Q + K * K) / (1.0 + K / Q + K * K);
 
     st->d->b[0] = pb[0] * rb[0];
     st->d->b[1] = pb[0] * rb[1] + pb[1] * rb[0];
     st->d->b[2] = pb[0] * rb[2] + pb[1] * rb[1] + pb[2] * rb[0];
     st->d->b[3] = pb[1] * rb[2] + pb[2] * rb[1];
     st->d->b[4] = pb[2] * rb[2];
 
     st->d->a[0] = pa[0] * ra[0];
     st->d->a[1] = pa[0] * ra[1] + pa[1] * ra[0];
     st->d->a[2] = pa[0] * ra[2] + pa[1] * ra[1] + pa[2] * ra[0];
     st->d->a[3] = pa[1] * ra[2] + pa[2] * ra[1];
     st->d->a[4] = pa[2] * ra[2];
 
     for (i = 0; i < 5; ++i) {
         for (j = 0; j < 5; ++j) {
             st->d->v[i][j] = 0.0;
         }
     }
 }
 
 static int ebur128_init_channel_map(FFEBUR128State * st)
 {
     size_t i;
     st->d->channel_map =

         (int *) av_malloc_array(st->channels, sizeof(*st->d->channel_map));

     if (!st->d->channel_map)
         return AVERROR(ENOMEM);
     if (st->channels == 4) {
         st->d->channel_map[0] = FF_EBUR128_LEFT;
         st->d->channel_map[1] = FF_EBUR128_RIGHT;
         st->d->channel_map[2] = FF_EBUR128_LEFT_SURROUND;
         st->d->channel_map[3] = FF_EBUR128_RIGHT_SURROUND;
     } else if (st->channels == 5) {
         st->d->channel_map[0] = FF_EBUR128_LEFT;
         st->d->channel_map[1] = FF_EBUR128_RIGHT;
         st->d->channel_map[2] = FF_EBUR128_CENTER;
         st->d->channel_map[3] = FF_EBUR128_LEFT_SURROUND;
         st->d->channel_map[4] = FF_EBUR128_RIGHT_SURROUND;
     } else {
         for (i = 0; i < st->channels; ++i) {
             switch (i) {
             case 0:
                 st->d->channel_map[i] = FF_EBUR128_LEFT;
                 break;
             case 1:
                 st->d->channel_map[i] = FF_EBUR128_RIGHT;
                 break;
             case 2:
                 st->d->channel_map[i] = FF_EBUR128_CENTER;
                 break;
             case 3:
                 st->d->channel_map[i] = FF_EBUR128_UNUSED;
                 break;
             case 4:
                 st->d->channel_map[i] = FF_EBUR128_LEFT_SURROUND;
                 break;
             case 5:
                 st->d->channel_map[i] = FF_EBUR128_RIGHT_SURROUND;
                 break;
             default:
                 st->d->channel_map[i] = FF_EBUR128_UNUSED;
                 break;
             }
         }
     }
     return 0;
 }
 
 static inline void init_histogram(void)
 {
     int i;
     /* initialize static constants */
     histogram_energy_boundaries[0] = pow(10.0, (-70.0 + 0.691) / 10.0);
     for (i = 0; i < 1000; ++i) {
         histogram_energies[i] =
             pow(10.0, ((double) i / 10.0 - 69.95 + 0.691) / 10.0);
     }
     for (i = 1; i < 1001; ++i) {
         histogram_energy_boundaries[i] =
             pow(10.0, ((double) i / 10.0 - 70.0 + 0.691) / 10.0);
     }
 }
 
 FFEBUR128State *ff_ebur128_init(unsigned int channels,
                                 unsigned long samplerate,
                                 unsigned long window, int mode)
 {
     int errcode;
     FFEBUR128State *st;
 

     st = (FFEBUR128State *) av_malloc(sizeof(*st));

     CHECK_ERROR(!st, 0, exit)
     st->d = (struct FFEBUR128StateInternal *)

         av_malloc(sizeof(*st->d));

     CHECK_ERROR(!st->d, 0, free_state)
     st->channels = channels;
     errcode = ebur128_init_channel_map(st);
     CHECK_ERROR(errcode, 0, free_internal)
 
     st->d->sample_peak =

         (double *) av_mallocz_array(channels, sizeof(*st->d->sample_peak));

     CHECK_ERROR(!st->d->sample_peak, 0, free_channel_map)
 
     st->samplerate = samplerate;
     st->d->samples_in_100ms = (st->samplerate + 5) / 10;
     st->mode = mode;
     if ((mode & FF_EBUR128_MODE_S) == FF_EBUR128_MODE_S) {
         st->d->window = FFMAX(window, 3000);
     } else if ((mode & FF_EBUR128_MODE_M) == FF_EBUR128_MODE_M) {
         st->d->window = FFMAX(window, 400);
     } else {
         goto free_sample_peak;
     }
     st->d->audio_data_frames = st->samplerate * st->d->window / 1000;
     if (st->d->audio_data_frames % st->d->samples_in_100ms) {
         /* round up to multiple of samples_in_100ms */
         st->d->audio_data_frames = st->d->audio_data_frames
             + st->d->samples_in_100ms
             - (st->d->audio_data_frames % st->d->samples_in_100ms);
     }
     st->d->audio_data =
         (double *) av_mallocz_array(st->d->audio_data_frames,

                                     st->channels * sizeof(*st->d->audio_data));

     CHECK_ERROR(!st->d->audio_data, 0, free_sample_peak)
 
     ebur128_init_filter(st);
 
     st->d->block_energy_histogram =

         av_mallocz(1000 * sizeof(*st->d->block_energy_histogram));

     CHECK_ERROR(!st->d->block_energy_histogram, 0, free_audio_data)
     st->d->short_term_block_energy_histogram =

         av_mallocz(1000 * sizeof(*st->d->short_term_block_energy_histogram));

     CHECK_ERROR(!st->d->short_term_block_energy_histogram, 0,
                 free_block_energy_histogram)
     st->d->short_term_frame_counter = 0;
 
     /* the first block needs 400ms of audio data */
     st->d->needed_frames = st->d->samples_in_100ms * 4;
     /* start at the beginning of the buffer */
     st->d->audio_data_index = 0;
 
     if (ff_thread_once(&histogram_init, &init_histogram) != 0)
         goto free_short_term_block_energy_histogram;
 

     st->d->data_ptrs = av_malloc_array(channels, sizeof(*st->d->data_ptrs));

     CHECK_ERROR(!st->d->data_ptrs, 0,
                 free_short_term_block_energy_histogram);
 
     return st;
 
 free_short_term_block_energy_histogram:
     av_free(st->d->short_term_block_energy_histogram);
 free_block_energy_histogram:
     av_free(st->d->block_energy_histogram);
 free_audio_data:
     av_free(st->d->audio_data);
 free_sample_peak:
     av_free(st->d->sample_peak);
 free_channel_map:
     av_free(st->d->channel_map);
 free_internal:
     av_free(st->d);
 free_state:
     av_free(st);
 exit:
     return NULL;
 }
 
 void ff_ebur128_destroy(FFEBUR128State ** st)
 {
     av_free((*st)->d->block_energy_histogram);
     av_free((*st)->d->short_term_block_energy_histogram);
     av_free((*st)->d->audio_data);
     av_free((*st)->d->channel_map);
     av_free((*st)->d->sample_peak);
     av_free((*st)->d->data_ptrs);
     av_free((*st)->d);
     av_free(*st);
     *st = NULL;
 }
 

 #define EBUR128_FILTER(type, scaling_factor)                                       \

 static void ebur128_filter_##type(FFEBUR128State* st, const type** srcs,           \
                                   size_t src_index, size_t frames,                 \
                                   int stride) {                                    \
     double* audio_data = st->d->audio_data + st->d->audio_data_index;              \
     size_t i, c;                                                                   \
                                                                                    \
     if ((st->mode & FF_EBUR128_MODE_SAMPLE_PEAK) == FF_EBUR128_MODE_SAMPLE_PEAK) { \
         for (c = 0; c < st->channels; ++c) {                                       \
             double max = 0.0;                                                      \
             for (i = 0; i < frames; ++i) {                                         \
                 type v = srcs[c][src_index + i * stride];                          \
                 if (v > max) {                                                     \
                     max =        v;                                                \
                 } else if (-v > max) {                                             \
                     max = -1.0 * v;                                                \
                 }                                                                  \
             }                                                                      \
             max /= scaling_factor;                                                 \
             if (max > st->d->sample_peak[c]) st->d->sample_peak[c] = max;          \
         }                                                                          \
     }                                                                              \
     for (c = 0; c < st->channels; ++c) {                                           \
         int ci = st->d->channel_map[c] - 1;                                        \
         if (ci < 0) continue;                                                      \
         else if (ci == FF_EBUR128_DUAL_MONO - 1) ci = 0; /*dual mono */            \
         for (i = 0; i < frames; ++i) {                                             \
             st->d->v[ci][0] = (double) (srcs[c][src_index + i * stride] / scaling_factor) \
                          - st->d->a[1] * st->d->v[ci][1]                           \
                          - st->d->a[2] * st->d->v[ci][2]                           \
                          - st->d->a[3] * st->d->v[ci][3]                           \
                          - st->d->a[4] * st->d->v[ci][4];                          \
             audio_data[i * st->channels + c] =                                     \
                            st->d->b[0] * st->d->v[ci][0]                           \
                          + st->d->b[1] * st->d->v[ci][1]                           \
                          + st->d->b[2] * st->d->v[ci][2]                           \
                          + st->d->b[3] * st->d->v[ci][3]                           \
                          + st->d->b[4] * st->d->v[ci][4];                          \
             st->d->v[ci][4] = st->d->v[ci][3];                                     \
             st->d->v[ci][3] = st->d->v[ci][2];                                     \
             st->d->v[ci][2] = st->d->v[ci][1];                                     \
             st->d->v[ci][1] = st->d->v[ci][0];                                     \
         }                                                                          \
         st->d->v[ci][4] = fabs(st->d->v[ci][4]) < DBL_MIN ? 0.0 : st->d->v[ci][4]; \
         st->d->v[ci][3] = fabs(st->d->v[ci][3]) < DBL_MIN ? 0.0 : st->d->v[ci][3]; \
         st->d->v[ci][2] = fabs(st->d->v[ci][2]) < DBL_MIN ? 0.0 : st->d->v[ci][2]; \
         st->d->v[ci][1] = fabs(st->d->v[ci][1]) < DBL_MIN ? 0.0 : st->d->v[ci][1]; \
     }                                                                              \
 }

 EBUR128_FILTER(short, -((double)SHRT_MIN))
 EBUR128_FILTER(int, -((double)INT_MIN))
 EBUR128_FILTER(float,  1.0)
 EBUR128_FILTER(double, 1.0)

 
 static double ebur128_energy_to_loudness(double energy)
 {

     return 10 * log10(energy) - 0.691;

 }
 
 static size_t find_histogram_index(double energy)
 {
     size_t index_min = 0;
     size_t index_max = 1000;
     size_t index_mid;
 
     do {
         index_mid = (index_min + index_max) / 2;
         if (energy >= histogram_energy_boundaries[index_mid]) {
             index_min = index_mid;
         } else {
             index_max = index_mid;
         }
     } while (index_max - index_min != 1);
 
     return index_min;
 }
 
 static void ebur128_calc_gating_block(FFEBUR128State * st,
                                       size_t frames_per_block,
                                       double *optional_output)
 {
     size_t i, c;
     double sum = 0.0;
     double channel_sum;
     for (c = 0; c < st->channels; ++c) {
         if (st->d->channel_map[c] == FF_EBUR128_UNUSED)
             continue;
         channel_sum = 0.0;
         if (st->d->audio_data_index < frames_per_block * st->channels) {
             for (i = 0; i < st->d->audio_data_index / st->channels; ++i) {
                 channel_sum += st->d->audio_data[i * st->channels + c] *
                     st->d->audio_data[i * st->channels + c];
             }
             for (i = st->d->audio_data_frames -
                  (frames_per_block -
                   st->d->audio_data_index / st->channels);
                  i < st->d->audio_data_frames; ++i) {
                 channel_sum += st->d->audio_data[i * st->channels + c] *
                     st->d->audio_data[i * st->channels + c];
             }
         } else {
             for (i =
                  st->d->audio_data_index / st->channels - frames_per_block;
                  i < st->d->audio_data_index / st->channels; ++i) {
                 channel_sum +=
                     st->d->audio_data[i * st->channels +
                                       c] * st->d->audio_data[i *
                                                              st->channels +
                                                              c];
             }
         }
         if (st->d->channel_map[c] == FF_EBUR128_Mp110 ||
             st->d->channel_map[c] == FF_EBUR128_Mm110 ||
             st->d->channel_map[c] == FF_EBUR128_Mp060 ||
             st->d->channel_map[c] == FF_EBUR128_Mm060 ||
             st->d->channel_map[c] == FF_EBUR128_Mp090 ||
             st->d->channel_map[c] == FF_EBUR128_Mm090) {
             channel_sum *= 1.41;
         } else if (st->d->channel_map[c] == FF_EBUR128_DUAL_MONO) {
             channel_sum *= 2.0;
         }
         sum += channel_sum;
     }
     sum /= (double) frames_per_block;
     if (optional_output) {
         *optional_output = sum;
     } else if (sum >= histogram_energy_boundaries[0]) {
         ++st->d->block_energy_histogram[find_histogram_index(sum)];
     }
 }
 
 int ff_ebur128_set_channel(FFEBUR128State * st,
                            unsigned int channel_number, int value)
 {
     if (channel_number >= st->channels) {
         return 1;
     }
     if (value == FF_EBUR128_DUAL_MONO &&
         (st->channels != 1 || channel_number != 0)) {
         return 1;
     }
     st->d->channel_map[channel_number] = value;
     return 0;
 }
 
 static int ebur128_energy_shortterm(FFEBUR128State * st, double *out);
 #define FF_EBUR128_ADD_FRAMES_PLANAR(type)                                             \
 void ff_ebur128_add_frames_planar_##type(FFEBUR128State* st, const type** srcs,        \
                                  size_t frames, int stride) {                          \
     size_t src_index = 0;                                                              \
     while (frames > 0) {                                                               \
         if (frames >= st->d->needed_frames) {                                          \
             ebur128_filter_##type(st, srcs, src_index, st->d->needed_frames, stride);  \
             src_index += st->d->needed_frames * stride;                                \
             frames -= st->d->needed_frames;                                            \
             st->d->audio_data_index += st->d->needed_frames * st->channels;            \
             /* calculate the new gating block */                                       \
             if ((st->mode & FF_EBUR128_MODE_I) == FF_EBUR128_MODE_I) {                 \
                 ebur128_calc_gating_block(st, st->d->samples_in_100ms * 4, NULL);      \
             }                                                                          \
             if ((st->mode & FF_EBUR128_MODE_LRA) == FF_EBUR128_MODE_LRA) {             \
                 st->d->short_term_frame_counter += st->d->needed_frames;               \
                 if (st->d->short_term_frame_counter == st->d->samples_in_100ms * 30) { \
                     double st_energy;                                                  \
                     ebur128_energy_shortterm(st, &st_energy);                          \
                     if (st_energy >= histogram_energy_boundaries[0]) {                 \
                         ++st->d->short_term_block_energy_histogram[                    \
                                                     find_histogram_index(st_energy)];  \
                     }                                                                  \
                     st->d->short_term_frame_counter = st->d->samples_in_100ms * 20;    \
                 }                                                                      \
             }                                                                          \
             /* 100ms are needed for all blocks besides the first one */                \
             st->d->needed_frames = st->d->samples_in_100ms;                            \
             /* reset audio_data_index when buffer full */                              \
             if (st->d->audio_data_index == st->d->audio_data_frames * st->channels) {  \
                 st->d->audio_data_index = 0;                                           \
             }                                                                          \
         } else {                                                                       \
             ebur128_filter_##type(st, srcs, src_index, frames, stride);                \
             st->d->audio_data_index += frames * st->channels;                          \
             if ((st->mode & FF_EBUR128_MODE_LRA) == FF_EBUR128_MODE_LRA) {             \
                 st->d->short_term_frame_counter += frames;                             \
             }                                                                          \
             st->d->needed_frames -= frames;                                            \
             frames = 0;                                                                \
         }                                                                              \
     }                                                                                  \
 }
 FF_EBUR128_ADD_FRAMES_PLANAR(short)
 FF_EBUR128_ADD_FRAMES_PLANAR(int)
 FF_EBUR128_ADD_FRAMES_PLANAR(float)
 FF_EBUR128_ADD_FRAMES_PLANAR(double)
 #define FF_EBUR128_ADD_FRAMES(type)                                            \
 void ff_ebur128_add_frames_##type(FFEBUR128State* st, const type* src,         \
                                     size_t frames) {                           \
   int i;                                                                       \
   const type **buf = (const type**)st->d->data_ptrs;                           \
   for (i = 0; i < st->channels; i++)                                           \
     buf[i] = src + i;                                                          \
   ff_ebur128_add_frames_planar_##type(st, buf, frames, st->channels);          \
 }
 FF_EBUR128_ADD_FRAMES(short)
 FF_EBUR128_ADD_FRAMES(int)
 FF_EBUR128_ADD_FRAMES(float)
 FF_EBUR128_ADD_FRAMES(double)
 

 static int ebur128_calc_relative_threshold(FFEBUR128State **sts, size_t size,

                                            double *relative_threshold)
 {

     size_t i, j;
     int above_thresh_counter = 0;

     *relative_threshold = 0.0;
 

     for (i = 0; i < size; i++) {
         unsigned long *block_energy_histogram = sts[i]->d->block_energy_histogram;
         for (j = 0; j < 1000; ++j) {
             *relative_threshold += block_energy_histogram[j] * histogram_energies[j];
             above_thresh_counter += block_energy_histogram[j];
         }

     }
 

     if (above_thresh_counter != 0) {
         *relative_threshold /= (double)above_thresh_counter;

         *relative_threshold *= RELATIVE_GATE_FACTOR;
     }
 

     return above_thresh_counter;

 }
 
 static int ebur128_gated_loudness(FFEBUR128State ** sts, size_t size,
                                   double *out)
 {
     double gated_loudness = 0.0;
     double relative_threshold;
     size_t above_thresh_counter;
     size_t i, j, start_index;
 

     for (i = 0; i < size; i++)
         if ((sts[i]->mode & FF_EBUR128_MODE_I) != FF_EBUR128_MODE_I)

             return AVERROR(EINVAL);
 

     if (!ebur128_calc_relative_threshold(sts, size, &relative_threshold)) {

         *out = -HUGE_VAL;
         return 0;
     }
 
     above_thresh_counter = 0;
     if (relative_threshold < histogram_energy_boundaries[0]) {
         start_index = 0;
     } else {
         start_index = find_histogram_index(relative_threshold);
         if (relative_threshold > histogram_energies[start_index]) {
             ++start_index;
         }
     }
     for (i = 0; i < size; i++) {
         for (j = start_index; j < 1000; ++j) {
             gated_loudness += sts[i]->d->block_energy_histogram[j] *
                 histogram_energies[j];
             above_thresh_counter += sts[i]->d->block_energy_histogram[j];
         }
     }
     if (!above_thresh_counter) {
         *out = -HUGE_VAL;
         return 0;
     }
     gated_loudness /= (double) above_thresh_counter;
     *out = ebur128_energy_to_loudness(gated_loudness);
     return 0;
 }
 
 int ff_ebur128_relative_threshold(FFEBUR128State * st, double *out)
 {
     double relative_threshold;
 

     if ((st->mode & FF_EBUR128_MODE_I) != FF_EBUR128_MODE_I)

         return AVERROR(EINVAL);
 

     if (!ebur128_calc_relative_threshold(&st, 1, &relative_threshold)) {

         *out = -70.0;
         return 0;
     }
 
     *out = ebur128_energy_to_loudness(relative_threshold);
     return 0;
 }
 
 int ff_ebur128_loudness_global(FFEBUR128State * st, double *out)
 {
     return ebur128_gated_loudness(&st, 1, out);
 }
 
 int ff_ebur128_loudness_global_multiple(FFEBUR128State ** sts, size_t size,
                                         double *out)
 {
     return ebur128_gated_loudness(sts, size, out);
 }
 
 static int ebur128_energy_in_interval(FFEBUR128State * st,
                                       size_t interval_frames, double *out)
 {
     if (interval_frames > st->d->audio_data_frames) {
         return AVERROR(EINVAL);
     }
     ebur128_calc_gating_block(st, interval_frames, out);
     return 0;
 }
 
 static int ebur128_energy_shortterm(FFEBUR128State * st, double *out)
 {
     return ebur128_energy_in_interval(st, st->d->samples_in_100ms * 30,
                                       out);
 }
 
 int ff_ebur128_loudness_momentary(FFEBUR128State * st, double *out)
 {
     double energy;
     int error = ebur128_energy_in_interval(st, st->d->samples_in_100ms * 4,
                                            &energy);
     if (error) {
         return error;
     } else if (energy <= 0.0) {
         *out = -HUGE_VAL;
         return 0;
     }
     *out = ebur128_energy_to_loudness(energy);
     return 0;
 }
 
 int ff_ebur128_loudness_shortterm(FFEBUR128State * st, double *out)
 {
     double energy;
     int error = ebur128_energy_shortterm(st, &energy);
     if (error) {
         return error;
     } else if (energy <= 0.0) {
         *out = -HUGE_VAL;
         return 0;
     }
     *out = ebur128_energy_to_loudness(energy);
     return 0;
 }
 
 int ff_ebur128_loudness_window(FFEBUR128State * st,
                                unsigned long window, double *out)
 {
     double energy;
     size_t interval_frames = st->samplerate * window / 1000;
     int error = ebur128_energy_in_interval(st, interval_frames, &energy);
     if (error) {
         return error;
     } else if (energy <= 0.0) {
         *out = -HUGE_VAL;
         return 0;
     }
     *out = ebur128_energy_to_loudness(energy);
     return 0;
 }
 
 /* EBU - TECH 3342 */
 int ff_ebur128_loudness_range_multiple(FFEBUR128State ** sts, size_t size,
                                        double *out)
 {
     size_t i, j;
     size_t stl_size;
     double stl_power, stl_integrated;
     /* High and low percentile energy */
     double h_en, l_en;
     unsigned long hist[1000] = { 0 };
     size_t percentile_low, percentile_high;
     size_t index;
 
     for (i = 0; i < size; ++i) {
         if (sts[i]) {
             if ((sts[i]->mode & FF_EBUR128_MODE_LRA) !=
                 FF_EBUR128_MODE_LRA) {
                 return AVERROR(EINVAL);
             }
         }
     }
 
     stl_size = 0;
     stl_power = 0.0;
     for (i = 0; i < size; ++i) {
         if (!sts[i])
             continue;
         for (j = 0; j < 1000; ++j) {
             hist[j] += sts[i]->d->short_term_block_energy_histogram[j];
             stl_size += sts[i]->d->short_term_block_energy_histogram[j];
             stl_power += sts[i]->d->short_term_block_energy_histogram[j]
                 * histogram_energies[j];
         }
     }
     if (!stl_size) {
         *out = 0.0;
         return 0;
     }
 
     stl_power /= stl_size;
     stl_integrated = MINUS_20DB * stl_power;
 
     if (stl_integrated < histogram_energy_boundaries[0]) {
         index = 0;
     } else {
         index = find_histogram_index(stl_integrated);
         if (stl_integrated > histogram_energies[index]) {
             ++index;
         }
     }
     stl_size = 0;
     for (j = index; j < 1000; ++j) {
         stl_size += hist[j];
     }
     if (!stl_size) {
         *out = 0.0;
         return 0;
     }
 
     percentile_low = (size_t) ((stl_size - 1) * 0.1 + 0.5);
     percentile_high = (size_t) ((stl_size - 1) * 0.95 + 0.5);
 
     stl_size = 0;
     j = index;
     while (stl_size <= percentile_low) {
         stl_size += hist[j++];
     }
     l_en = histogram_energies[j - 1];
     while (stl_size <= percentile_high) {
         stl_size += hist[j++];
     }
     h_en = histogram_energies[j - 1];
     *out =
         ebur128_energy_to_loudness(h_en) -
         ebur128_energy_to_loudness(l_en);
     return 0;
 }
 
 int ff_ebur128_loudness_range(FFEBUR128State * st, double *out)
 {
     return ff_ebur128_loudness_range_multiple(&st, 1, out);
 }
 
 int ff_ebur128_sample_peak(FFEBUR128State * st,
                            unsigned int channel_number, double *out)
 {
     if ((st->mode & FF_EBUR128_MODE_SAMPLE_PEAK) !=
         FF_EBUR128_MODE_SAMPLE_PEAK) {
         return AVERROR(EINVAL);
     } else if (channel_number >= st->channels) {
         return AVERROR(EINVAL);
     }
     *out = st->d->sample_peak[channel_number];
     return 0;
 }