Commit | Line | Data |
---|---|---|
2ba45a60 DM |
1 | /* |
2 | * (I)DCT Transforms | |
3 | * Copyright (c) 2009 Peter Ross <pross@xvid.org> | |
4 | * Copyright (c) 2010 Alex Converse <alex.converse@gmail.com> | |
5 | * Copyright (c) 2010 Vitor Sessak | |
6 | * | |
7 | * This file is part of FFmpeg. | |
8 | * | |
9 | * FFmpeg is free software; you can redistribute it and/or | |
10 | * modify it under the terms of the GNU Lesser General Public | |
11 | * License as published by the Free Software Foundation; either | |
12 | * version 2.1 of the License, or (at your option) any later version. | |
13 | * | |
14 | * FFmpeg is distributed in the hope that it will be useful, | |
15 | * but WITHOUT ANY WARRANTY; without even the implied warranty of | |
16 | * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU | |
17 | * Lesser General Public License for more details. | |
18 | * | |
19 | * You should have received a copy of the GNU Lesser General Public | |
20 | * License along with FFmpeg; if not, write to the Free Software | |
21 | * Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA | |
22 | */ | |
23 | ||
24 | /** | |
25 | * @file | |
26 | * (Inverse) Discrete Cosine Transforms. These are also known as the | |
27 | * type II and type III DCTs respectively. | |
28 | */ | |
29 | ||
30 | #include <math.h> | |
31 | #include <string.h> | |
32 | ||
33 | #include "libavutil/mathematics.h" | |
34 | #include "dct.h" | |
35 | #include "dct32.h" | |
36 | ||
37 | /* sin((M_PI * x / (2 * n)) */ | |
38 | #define SIN(s, n, x) (s->costab[(n) - (x)]) | |
39 | ||
40 | /* cos((M_PI * x / (2 * n)) */ | |
41 | #define COS(s, n, x) (s->costab[x]) | |
42 | ||
43 | static void dst_calc_I_c(DCTContext *ctx, FFTSample *data) | |
44 | { | |
45 | int n = 1 << ctx->nbits; | |
46 | int i; | |
47 | ||
48 | data[0] = 0; | |
49 | for (i = 1; i < n / 2; i++) { | |
50 | float tmp1 = data[i ]; | |
51 | float tmp2 = data[n - i]; | |
52 | float s = SIN(ctx, n, 2 * i); | |
53 | ||
54 | s *= tmp1 + tmp2; | |
55 | tmp1 = (tmp1 - tmp2) * 0.5f; | |
56 | data[i] = s + tmp1; | |
57 | data[n - i] = s - tmp1; | |
58 | } | |
59 | ||
60 | data[n / 2] *= 2; | |
61 | ctx->rdft.rdft_calc(&ctx->rdft, data); | |
62 | ||
63 | data[0] *= 0.5f; | |
64 | ||
65 | for (i = 1; i < n - 2; i += 2) { | |
66 | data[i + 1] += data[i - 1]; | |
67 | data[i] = -data[i + 2]; | |
68 | } | |
69 | ||
70 | data[n - 1] = 0; | |
71 | } | |
72 | ||
73 | static void dct_calc_I_c(DCTContext *ctx, FFTSample *data) | |
74 | { | |
75 | int n = 1 << ctx->nbits; | |
76 | int i; | |
77 | float next = -0.5f * (data[0] - data[n]); | |
78 | ||
79 | for (i = 0; i < n / 2; i++) { | |
80 | float tmp1 = data[i]; | |
81 | float tmp2 = data[n - i]; | |
82 | float s = SIN(ctx, n, 2 * i); | |
83 | float c = COS(ctx, n, 2 * i); | |
84 | ||
85 | c *= tmp1 - tmp2; | |
86 | s *= tmp1 - tmp2; | |
87 | ||
88 | next += c; | |
89 | ||
90 | tmp1 = (tmp1 + tmp2) * 0.5f; | |
91 | data[i] = tmp1 - s; | |
92 | data[n - i] = tmp1 + s; | |
93 | } | |
94 | ||
95 | ctx->rdft.rdft_calc(&ctx->rdft, data); | |
96 | data[n] = data[1]; | |
97 | data[1] = next; | |
98 | ||
99 | for (i = 3; i <= n; i += 2) | |
100 | data[i] = data[i - 2] - data[i]; | |
101 | } | |
102 | ||
103 | static void dct_calc_III_c(DCTContext *ctx, FFTSample *data) | |
104 | { | |
105 | int n = 1 << ctx->nbits; | |
106 | int i; | |
107 | ||
108 | float next = data[n - 1]; | |
109 | float inv_n = 1.0f / n; | |
110 | ||
111 | for (i = n - 2; i >= 2; i -= 2) { | |
112 | float val1 = data[i]; | |
113 | float val2 = data[i - 1] - data[i + 1]; | |
114 | float c = COS(ctx, n, i); | |
115 | float s = SIN(ctx, n, i); | |
116 | ||
117 | data[i] = c * val1 + s * val2; | |
118 | data[i + 1] = s * val1 - c * val2; | |
119 | } | |
120 | ||
121 | data[1] = 2 * next; | |
122 | ||
123 | ctx->rdft.rdft_calc(&ctx->rdft, data); | |
124 | ||
125 | for (i = 0; i < n / 2; i++) { | |
126 | float tmp1 = data[i] * inv_n; | |
127 | float tmp2 = data[n - i - 1] * inv_n; | |
128 | float csc = ctx->csc2[i] * (tmp1 - tmp2); | |
129 | ||
130 | tmp1 += tmp2; | |
131 | data[i] = tmp1 + csc; | |
132 | data[n - i - 1] = tmp1 - csc; | |
133 | } | |
134 | } | |
135 | ||
136 | static void dct_calc_II_c(DCTContext *ctx, FFTSample *data) | |
137 | { | |
138 | int n = 1 << ctx->nbits; | |
139 | int i; | |
140 | float next; | |
141 | ||
142 | for (i = 0; i < n / 2; i++) { | |
143 | float tmp1 = data[i]; | |
144 | float tmp2 = data[n - i - 1]; | |
145 | float s = SIN(ctx, n, 2 * i + 1); | |
146 | ||
147 | s *= tmp1 - tmp2; | |
148 | tmp1 = (tmp1 + tmp2) * 0.5f; | |
149 | ||
150 | data[i] = tmp1 + s; | |
151 | data[n-i-1] = tmp1 - s; | |
152 | } | |
153 | ||
154 | ctx->rdft.rdft_calc(&ctx->rdft, data); | |
155 | ||
156 | next = data[1] * 0.5; | |
157 | data[1] *= -1; | |
158 | ||
159 | for (i = n - 2; i >= 0; i -= 2) { | |
160 | float inr = data[i ]; | |
161 | float ini = data[i + 1]; | |
162 | float c = COS(ctx, n, i); | |
163 | float s = SIN(ctx, n, i); | |
164 | ||
165 | data[i] = c * inr + s * ini; | |
166 | data[i + 1] = next; | |
167 | ||
168 | next += s * inr - c * ini; | |
169 | } | |
170 | } | |
171 | ||
172 | static void dct32_func(DCTContext *ctx, FFTSample *data) | |
173 | { | |
174 | ctx->dct32(data, data); | |
175 | } | |
176 | ||
177 | av_cold int ff_dct_init(DCTContext *s, int nbits, enum DCTTransformType inverse) | |
178 | { | |
179 | int n = 1 << nbits; | |
180 | int i; | |
181 | ||
182 | memset(s, 0, sizeof(*s)); | |
183 | ||
184 | s->nbits = nbits; | |
185 | s->inverse = inverse; | |
186 | ||
187 | if (inverse == DCT_II && nbits == 5) { | |
188 | s->dct_calc = dct32_func; | |
189 | } else { | |
190 | ff_init_ff_cos_tabs(nbits + 2); | |
191 | ||
192 | s->costab = ff_cos_tabs[nbits + 2]; | |
193 | s->csc2 = av_malloc_array(n / 2, sizeof(FFTSample)); | |
194 | ||
195 | if (ff_rdft_init(&s->rdft, nbits, inverse == DCT_III) < 0) { | |
f6fa7814 | 196 | av_freep(&s->csc2); |
2ba45a60 DM |
197 | return -1; |
198 | } | |
199 | ||
200 | for (i = 0; i < n / 2; i++) | |
201 | s->csc2[i] = 0.5 / sin((M_PI / (2 * n) * (2 * i + 1))); | |
202 | ||
203 | switch (inverse) { | |
204 | case DCT_I : s->dct_calc = dct_calc_I_c; break; | |
205 | case DCT_II : s->dct_calc = dct_calc_II_c; break; | |
206 | case DCT_III: s->dct_calc = dct_calc_III_c; break; | |
207 | case DST_I : s->dct_calc = dst_calc_I_c; break; | |
208 | } | |
209 | } | |
210 | ||
211 | s->dct32 = ff_dct32_float; | |
212 | if (ARCH_X86) | |
213 | ff_dct_init_x86(s); | |
214 | ||
215 | return 0; | |
216 | } | |
217 | ||
218 | av_cold void ff_dct_end(DCTContext *s) | |
219 | { | |
220 | ff_rdft_end(&s->rdft); | |
f6fa7814 | 221 | av_freep(&s->csc2); |
2ba45a60 | 222 | } |