Merge remote-tracking branch 'qatar/master'
[ffmpeg.git] / 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 FFmpeg.
7  *
8  * FFmpeg 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  * FFmpeg 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 FFmpeg; 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 #include "libavutil/mem.h"
35 #define NR_ALLPASS_BANDS20 30
36 #define NR_ALLPASS_BANDS34 50
37 #define PS_AP_LINKS 3
38 static float pd_re_smooth[8*8*8];
39 static float pd_im_smooth[8*8*8];
40 static float HA[46][8][4];
41 static float HB[46][8][4];
42 static DECLARE_ALIGNED(16, float, f20_0_8) [ 8][8][2];
43 static DECLARE_ALIGNED(16, float, f34_0_12)[12][8][2];
44 static DECLARE_ALIGNED(16, float, f34_1_8) [ 8][8][2];
45 static DECLARE_ALIGNED(16, float, f34_2_4) [ 4][8][2];
46 static DECLARE_ALIGNED(16, float, Q_fract_allpass)[2][50][3][2];
47 static DECLARE_ALIGNED(16, float, phi_fract)[2][50][2];
48
49 static const float g0_Q8[] = {
50     0.00746082949812f, 0.02270420949825f, 0.04546865930473f, 0.07266113929591f,
51     0.09885108575264f, 0.11793710567217f, 0.125f
52 };
53
54 static const float g0_Q12[] = {
55     0.04081179924692f, 0.03812810994926f, 0.05144908135699f, 0.06399831151592f,
56     0.07428313801106f, 0.08100347892914f, 0.08333333333333f
57 };
58
59 static const float g1_Q8[] = {
60     0.01565675600122f, 0.03752716391991f, 0.05417891378782f, 0.08417044116767f,
61     0.10307344158036f, 0.12222452249753f, 0.125f
62 };
63
64 static const float g2_Q4[] = {
65     -0.05908211155639f, -0.04871498374946f, 0.0f,   0.07778723915851f,
66      0.16486303567403f,  0.23279856662996f, 0.25f
67 };
68
69 static void make_filters_from_proto(float (*filter)[8][2], const float *proto, int bands)
70 {
71     int q, n;
72     for (q = 0; q < bands; q++) {
73         for (n = 0; n < 7; n++) {
74             double theta = 2 * M_PI * (q + 0.5) * (n - 6) / bands;
75             filter[q][n][0] = proto[n] *  cos(theta);
76             filter[q][n][1] = proto[n] * -sin(theta);
77         }
78     }
79 }
80
81 static void ps_tableinit(void)
82 {
83     static const float ipdopd_sin[] = { 0, M_SQRT1_2, 1,  M_SQRT1_2,  0, -M_SQRT1_2, -1, -M_SQRT1_2 };
84     static const float ipdopd_cos[] = { 1, M_SQRT1_2, 0, -M_SQRT1_2, -1, -M_SQRT1_2,  0,  M_SQRT1_2 };
85     int pd0, pd1, pd2;
86
87     static const float iid_par_dequant[] = {
88         //iid_par_dequant_default
89         0.05623413251903, 0.12589254117942, 0.19952623149689, 0.31622776601684,
90         0.44668359215096, 0.63095734448019, 0.79432823472428, 1,
91         1.25892541179417, 1.58489319246111, 2.23872113856834, 3.16227766016838,
92         5.01187233627272, 7.94328234724282, 17.7827941003892,
93         //iid_par_dequant_fine
94         0.00316227766017, 0.00562341325190, 0.01,             0.01778279410039,
95         0.03162277660168, 0.05623413251903, 0.07943282347243, 0.11220184543020,
96         0.15848931924611, 0.22387211385683, 0.31622776601684, 0.39810717055350,
97         0.50118723362727, 0.63095734448019, 0.79432823472428, 1,
98         1.25892541179417, 1.58489319246111, 1.99526231496888, 2.51188643150958,
99         3.16227766016838, 4.46683592150963, 6.30957344480193, 8.91250938133745,
100         12.5892541179417, 17.7827941003892, 31.6227766016838, 56.2341325190349,
101         100,              177.827941003892, 316.227766016837,
102     };
103     static const float icc_invq[] = {
104         1, 0.937,      0.84118,    0.60092,    0.36764,   0,      -0.589,    -1
105     };
106     static const float acos_icc_invq[] = {
107         0, 0.35685527, 0.57133466, 0.92614472, 1.1943263, M_PI/2, 2.2006171, M_PI
108     };
109     int iid, icc;
110
111     int k, m;
112     static const int8_t f_center_20[] = {
113         -3, -1, 1, 3, 5, 7, 10, 14, 18, 22,
114     };
115     static const int8_t f_center_34[] = {
116          2,  6, 10, 14, 18, 22, 26, 30,
117         34,-10, -6, -2, 51, 57, 15, 21,
118         27, 33, 39, 45, 54, 66, 78, 42,
119        102, 66, 78, 90,102,114,126, 90,
120     };
121     static const float fractional_delay_links[] = { 0.43f, 0.75f, 0.347f };
122     const float fractional_delay_gain = 0.39f;
123
124     for (pd0 = 0; pd0 < 8; pd0++) {
125         float pd0_re = ipdopd_cos[pd0];
126         float pd0_im = ipdopd_sin[pd0];
127         for (pd1 = 0; pd1 < 8; pd1++) {
128             float pd1_re = ipdopd_cos[pd1];
129             float pd1_im = ipdopd_sin[pd1];
130             for (pd2 = 0; pd2 < 8; pd2++) {
131                 float pd2_re = ipdopd_cos[pd2];
132                 float pd2_im = ipdopd_sin[pd2];
133                 float re_smooth = 0.25f * pd0_re + 0.5f * pd1_re + pd2_re;
134                 float im_smooth = 0.25f * pd0_im + 0.5f * pd1_im + pd2_im;
135                 float pd_mag = 1 / sqrt(im_smooth * im_smooth + re_smooth * re_smooth);
136                 pd_re_smooth[pd0*64+pd1*8+pd2] = re_smooth * pd_mag;
137                 pd_im_smooth[pd0*64+pd1*8+pd2] = im_smooth * pd_mag;
138             }
139         }
140     }
141
142     for (iid = 0; iid < 46; iid++) {
143         float c = iid_par_dequant[iid]; ///< Linear Inter-channel Intensity Difference
144         float c1 = (float)M_SQRT2 / sqrtf(1.0f + c*c);
145         float c2 = c * c1;
146         for (icc = 0; icc < 8; icc++) {
147             /*if (PS_BASELINE || ps->icc_mode < 3)*/ {
148                 float alpha = 0.5f * acos_icc_invq[icc];
149                 float beta  = alpha * (c1 - c2) * (float)M_SQRT1_2;
150                 HA[iid][icc][0] = c2 * cosf(beta + alpha);
151                 HA[iid][icc][1] = c1 * cosf(beta - alpha);
152                 HA[iid][icc][2] = c2 * sinf(beta + alpha);
153                 HA[iid][icc][3] = c1 * sinf(beta - alpha);
154             } /* else */ {
155                 float alpha, gamma, mu, rho;
156                 float alpha_c, alpha_s, gamma_c, gamma_s;
157                 rho = FFMAX(icc_invq[icc], 0.05f);
158                 alpha = 0.5f * atan2f(2.0f * c * rho, c*c - 1.0f);
159                 mu = c + 1.0f / c;
160                 mu = sqrtf(1 + (4 * rho * rho - 4)/(mu * mu));
161                 gamma = atanf(sqrtf((1.0f - mu)/(1.0f + mu)));
162                 if (alpha < 0) alpha += M_PI/2;
163                 alpha_c = cosf(alpha);
164                 alpha_s = sinf(alpha);
165                 gamma_c = cosf(gamma);
166                 gamma_s = sinf(gamma);
167                 HB[iid][icc][0] =  M_SQRT2 * alpha_c * gamma_c;
168                 HB[iid][icc][1] =  M_SQRT2 * alpha_s * gamma_c;
169                 HB[iid][icc][2] = -M_SQRT2 * alpha_s * gamma_s;
170                 HB[iid][icc][3] =  M_SQRT2 * alpha_c * gamma_s;
171             }
172         }
173     }
174
175     for (k = 0; k < NR_ALLPASS_BANDS20; k++) {
176         double f_center, theta;
177         if (k < FF_ARRAY_ELEMS(f_center_20))
178             f_center = f_center_20[k] * 0.125;
179         else
180             f_center = k - 6.5f;
181         for (m = 0; m < PS_AP_LINKS; m++) {
182             theta = -M_PI * fractional_delay_links[m] * f_center;
183             Q_fract_allpass[0][k][m][0] = cos(theta);
184             Q_fract_allpass[0][k][m][1] = sin(theta);
185         }
186         theta = -M_PI*fractional_delay_gain*f_center;
187         phi_fract[0][k][0] = cos(theta);
188         phi_fract[0][k][1] = sin(theta);
189     }
190     for (k = 0; k < NR_ALLPASS_BANDS34; k++) {
191         double f_center, theta;
192         if (k < FF_ARRAY_ELEMS(f_center_34))
193             f_center = f_center_34[k] / 24.;
194         else
195             f_center = k - 26.5f;
196         for (m = 0; m < PS_AP_LINKS; m++) {
197             theta = -M_PI * fractional_delay_links[m] * f_center;
198             Q_fract_allpass[1][k][m][0] = cos(theta);
199             Q_fract_allpass[1][k][m][1] = sin(theta);
200         }
201         theta = -M_PI*fractional_delay_gain*f_center;
202         phi_fract[1][k][0] = cos(theta);
203         phi_fract[1][k][1] = sin(theta);
204     }
205
206     make_filters_from_proto(f20_0_8,  g0_Q8,   8);
207     make_filters_from_proto(f34_0_12, g0_Q12, 12);
208     make_filters_from_proto(f34_1_8,  g1_Q8,   8);
209     make_filters_from_proto(f34_2_4,  g2_Q4,   4);
210 }
211 #endif /* CONFIG_HARDCODED_TABLES */
212
213 #endif /* AACPS_TABLEGEN_H */