libavcodec/aacps_tablegen.h

   1 /*
   2  * Header file for hardcoded Parametric Stereo tables
   3  *
   4  * Copyright (c) 2010 Alex Converse <alex.converse@gmail.com>
   5  *
   6  * This file is part of Libav.
   7  *
   8  * Libav is free software; you can redistribute it and/or
   9  * modify it under the terms of the GNU Lesser General Public
  10  * License as published by the Free Software Foundation; either
  11  * version 2.1 of the License, or (at your option) any later version.
  12  *
  13  * Libav is distributed in the hope that it will be useful,
  14  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  15  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  16  * Lesser General Public License for more details.
  17  *
  18  * You should have received a copy of the GNU Lesser General Public
  19  * License along with Libav; if not, write to the Free Software
  20  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
  21  */
  22
  23 #ifndef AACPS_TABLEGEN_H
  24 #define AACPS_TABLEGEN_H
  25
  26 #include <stdint.h>
  27
  28 #if CONFIG_HARDCODED_TABLES
  29 #define ps_tableinit()
  30 #include "libavcodec/aacps_tables.h"
  31 #else
  32 #include "libavutil/common.h"
  33 #include "libavutil/mathematics.h"
  34 #define NR_ALLPASS_BANDS20 30
  35 #define NR_ALLPASS_BANDS34 50
  36 #define PS_AP_LINKS 3
  37 static float pd_re_smooth[8*8*8];
  38 static float pd_im_smooth[8*8*8];
  39 static float HA[46][8][4];
  40 static float HB[46][8][4];
  41 static float f20_0_8 [ 8][7][2];
  42 static float f34_0_12[12][7][2];
  43 static float f34_1_8 [ 8][7][2];
  44 static float f34_2_4 [ 4][7][2];
  45 static float Q_fract_allpass[2][50][3][2];
  46 static float phi_fract[2][50][2];
  47
  48 static const float g0_Q8[] = {
  49     0.00746082949812f, 0.02270420949825f, 0.04546865930473f, 0.07266113929591f,
  50     0.09885108575264f, 0.11793710567217f, 0.125f
  51 };
  52
  53 static const float g0_Q12[] = {
  54     0.04081179924692f, 0.03812810994926f, 0.05144908135699f, 0.06399831151592f,
  55     0.07428313801106f, 0.08100347892914f, 0.08333333333333f
  56 };
  57
  58 static const float g1_Q8[] = {
  59     0.01565675600122f, 0.03752716391991f, 0.05417891378782f, 0.08417044116767f,
  60     0.10307344158036f, 0.12222452249753f, 0.125f
  61 };
  62
  63 static const float g2_Q4[] = {
  64     -0.05908211155639f, -0.04871498374946f, 0.0f,   0.07778723915851f,
  65      0.16486303567403f,  0.23279856662996f, 0.25f
  66 };
  67
  68 static void make_filters_from_proto(float (*filter)[7][2], const float *proto, int bands)
  69 {
  70     int q, n;
  71     for (q = 0; q < bands; q++) {
  72         for (n = 0; n < 7; n++) {
  73             double theta = 2 * M_PI * (q + 0.5) * (n - 6) / bands;
  74             filter[q][n][0] = proto[n] *  cos(theta);
  75             filter[q][n][1] = proto[n] * -sin(theta);
  76         }
  77     }
  78 }
  79
  80 static void ps_tableinit(void)
  81 {
  82     static const float ipdopd_sin[] = { 0, M_SQRT1_2, 1,  M_SQRT1_2,  0, -M_SQRT1_2, -1, -M_SQRT1_2 };
  83     static const float ipdopd_cos[] = { 1, M_SQRT1_2, 0, -M_SQRT1_2, -1, -M_SQRT1_2,  0,  M_SQRT1_2 };
  84     int pd0, pd1, pd2;
  85
  86     static const float iid_par_dequant[] = {
  87         //iid_par_dequant_default
  88         0.05623413251903, 0.12589254117942, 0.19952623149689, 0.31622776601684,
  89         0.44668359215096, 0.63095734448019, 0.79432823472428, 1,
  90         1.25892541179417, 1.58489319246111, 2.23872113856834, 3.16227766016838,
  91         5.01187233627272, 7.94328234724282, 17.7827941003892,
  92         //iid_par_dequant_fine
  93         0.00316227766017, 0.00562341325190, 0.01,             0.01778279410039,
  94         0.03162277660168, 0.05623413251903, 0.07943282347243, 0.11220184543020,
  95         0.15848931924611, 0.22387211385683, 0.31622776601684, 0.39810717055350,
  96         0.50118723362727, 0.63095734448019, 0.79432823472428, 1,
  97         1.25892541179417, 1.58489319246111, 1.99526231496888, 2.51188643150958,
  98         3.16227766016838, 4.46683592150963, 6.30957344480193, 8.91250938133745,
  99         12.5892541179417, 17.7827941003892, 31.6227766016838, 56.2341325190349,
 100         100,              177.827941003892, 316.227766016837,
 101     };
 102     static const float icc_invq[] = {
 103         1, 0.937,      0.84118,    0.60092,    0.36764,   0,      -0.589,    -1
 104     };
 105     static const float acos_icc_invq[] = {
 106         0, 0.35685527, 0.57133466, 0.92614472, 1.1943263, M_PI/2, 2.2006171, M_PI
 107     };
 108     int iid, icc;
 109
 110     int k, m;
 111     static const int8_t f_center_20[] = {
 112         -3, -1, 1, 3, 5, 7, 10, 14, 18, 22,
 113     };
 114     static const int8_t f_center_34[] = {
 115          2,  6, 10, 14, 18, 22, 26, 30,
 116         34,-10, -6, -2, 51, 57, 15, 21,
 117         27, 33, 39, 45, 54, 66, 78, 42,
 118        102, 66, 78, 90,102,114,126, 90,
 119     };
 120     static const float fractional_delay_links[] = { 0.43f, 0.75f, 0.347f };
 121     const float fractional_delay_gain = 0.39f;
 122
 123     for (pd0 = 0; pd0 < 8; pd0++) {
 124         float pd0_re = ipdopd_cos[pd0];
 125         float pd0_im = ipdopd_sin[pd0];
 126         for (pd1 = 0; pd1 < 8; pd1++) {
 127             float pd1_re = ipdopd_cos[pd1];
 128             float pd1_im = ipdopd_sin[pd1];
 129             for (pd2 = 0; pd2 < 8; pd2++) {
 130                 float pd2_re = ipdopd_cos[pd2];
 131                 float pd2_im = ipdopd_sin[pd2];
 132                 float re_smooth = 0.25f * pd0_re + 0.5f * pd1_re + pd2_re;
 133                 float im_smooth = 0.25f * pd0_im + 0.5f * pd1_im + pd2_im;
 134                 float pd_mag = 1 / sqrt(im_smooth * im_smooth + re_smooth * re_smooth);
 135                 pd_re_smooth[pd0*64+pd1*8+pd2] = re_smooth * pd_mag;
 136                 pd_im_smooth[pd0*64+pd1*8+pd2] = im_smooth * pd_mag;
 137             }
 138         }
 139     }
 140
 141     for (iid = 0; iid < 46; iid++) {
 142         float c = iid_par_dequant[iid]; ///< Linear Inter-channel Intensity Difference
 143         float c1 = (float)M_SQRT2 / sqrtf(1.0f + c*c);
 144         float c2 = c * c1;
 145         for (icc = 0; icc < 8; icc++) {
 146             /*if (PS_BASELINE || ps->icc_mode < 3)*/ {
 147                 float alpha = 0.5f * acos_icc_invq[icc];
 148                 float beta  = alpha * (c1 - c2) * (float)M_SQRT1_2;
 149                 HA[iid][icc][0] = c2 * cosf(beta + alpha);
 150                 HA[iid][icc][1] = c1 * cosf(beta - alpha);
 151                 HA[iid][icc][2] = c2 * sinf(beta + alpha);
 152                 HA[iid][icc][3] = c1 * sinf(beta - alpha);
 153             } /* else */ {
 154                 float alpha, gamma, mu, rho;
 155                 float alpha_c, alpha_s, gamma_c, gamma_s;
 156                 rho = FFMAX(icc_invq[icc], 0.05f);
 157                 alpha = 0.5f * atan2f(2.0f * c * rho, c*c - 1.0f);
 158                 mu = c + 1.0f / c;
 159                 mu = sqrtf(1 + (4 * rho * rho - 4)/(mu * mu));
 160                 gamma = atanf(sqrtf((1.0f - mu)/(1.0f + mu)));
 161                 if (alpha < 0) alpha += M_PI/2;
 162                 alpha_c = cosf(alpha);
 163                 alpha_s = sinf(alpha);
 164                 gamma_c = cosf(gamma);
 165                 gamma_s = sinf(gamma);
 166                 HB[iid][icc][0] =  M_SQRT2 * alpha_c * gamma_c;
 167                 HB[iid][icc][1] =  M_SQRT2 * alpha_s * gamma_c;
 168                 HB[iid][icc][2] = -M_SQRT2 * alpha_s * gamma_s;
 169                 HB[iid][icc][3] =  M_SQRT2 * alpha_c * gamma_s;
 170             }
 171         }
 172     }
 173
 174     for (k = 0; k < NR_ALLPASS_BANDS20; k++) {
 175         double f_center, theta;
 176         if (k < FF_ARRAY_ELEMS(f_center_20))
 177             f_center = f_center_20[k] * 0.125;
 178         else
 179             f_center = k - 6.5f;
 180         for (m = 0; m < PS_AP_LINKS; m++) {
 181             theta = -M_PI * fractional_delay_links[m] * f_center;
 182             Q_fract_allpass[0][k][m][0] = cos(theta);
 183             Q_fract_allpass[0][k][m][1] = sin(theta);
 184         }
 185         theta = -M_PI*fractional_delay_gain*f_center;
 186         phi_fract[0][k][0] = cos(theta);
 187         phi_fract[0][k][1] = sin(theta);
 188     }
 189     for (k = 0; k < NR_ALLPASS_BANDS34; k++) {
 190         double f_center, theta;
 191         if (k < FF_ARRAY_ELEMS(f_center_34))
 192             f_center = f_center_34[k] / 24.;
 193         else
 194             f_center = k - 26.5f;
 195         for (m = 0; m < PS_AP_LINKS; m++) {
 196             theta = -M_PI * fractional_delay_links[m] * f_center;
 197             Q_fract_allpass[1][k][m][0] = cos(theta);
 198             Q_fract_allpass[1][k][m][1] = sin(theta);
 199         }
 200         theta = -M_PI*fractional_delay_gain*f_center;
 201         phi_fract[1][k][0] = cos(theta);
 202         phi_fract[1][k][1] = sin(theta);
 203     }
 204
 205     make_filters_from_proto(f20_0_8,  g0_Q8,   8);
 206     make_filters_from_proto(f34_0_12, g0_Q12, 12);
 207     make_filters_from_proto(f34_1_8,  g1_Q8,   8);
 208     make_filters_from_proto(f34_2_4,  g2_Q4,   4);
 209 }
 210 #endif /* CONFIG_HARDCODED_TABLES */
 211
 212 #endif /* AACPS_TABLEGEN_H */