FFmpeg
dct.c
Go to the documentation of this file.
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 
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 
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 
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 
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 
173 {
174  ctx->dct32(data, data);
175 }
176 
178 {
179  int n = 1 << nbits;
180  int i;
181  int ret;
182 
183  memset(s, 0, sizeof(*s));
184 
185  s->nbits = nbits;
186  s->inverse = inverse;
187 
188  if (inverse == DCT_II && nbits == 5) {
189  s->dct_calc = dct32_func;
190  } else {
191  ff_init_ff_cos_tabs(nbits + 2);
192 
193  s->costab = ff_cos_tabs[nbits + 2];
194  s->csc2 = av_malloc_array(n / 2, sizeof(FFTSample));
195  if (!s->csc2)
196  return AVERROR(ENOMEM);
197 
198  if ((ret = ff_rdft_init(&s->rdft, nbits, inverse == DCT_III)) < 0) {
199  av_freep(&s->csc2);
200  return ret;
201  }
202 
203  for (i = 0; i < n / 2; i++)
204  s->csc2[i] = 0.5 / sin((M_PI / (2 * n) * (2 * i + 1)));
205 
206  switch (inverse) {
207  case DCT_I : s->dct_calc = dct_calc_I_c; break;
208  case DCT_II : s->dct_calc = dct_calc_II_c; break;
209  case DCT_III: s->dct_calc = dct_calc_III_c; break;
210  case DST_I : s->dct_calc = dst_calc_I_c; break;
211  }
212  }
213 
214  s->dct32 = ff_dct32_float;
215  if (ARCH_X86)
216  ff_dct_init_x86(s);
217 
218  return 0;
219 }
220 
222 {
223  ff_rdft_end(&s->rdft);
224  av_freep(&s->csc2);
225 }
av_cold void ff_rdft_end(RDFTContext *s)
Definition: rdft.c:114
static void dct32_func(DCTContext *ctx, FFTSample *data)
Definition: dct.c:172
ptrdiff_t const GLvoid * data
Definition: opengl_enc.c:100
const float * costab
Definition: dct.h:36
Definition: avfft.h:95
static void dct_calc_I_c(DCTContext *ctx, FFTSample *data)
Definition: dct.c:73
RDFTContext rdft
Definition: dct.h:35
void ff_dct_init_x86(DCTContext *s)
Definition: dct_init.c:29
#define av_cold
Definition: attributes.h:82
#define f(width, name)
Definition: cbs_vp9.c:255
Undefined Behavior In the C some operations are like signed integer dereferencing freed accessing outside allocated Undefined Behavior must not occur in a C it is not safe even if the output of undefined operations is unused The unsafety may seem nit picking but Optimizing compilers have in fact optimized code on the assumption that no undefined Behavior occurs Optimizing code based on wrong assumptions can and has in some cases lead to effects beyond the output of computations The signed integer overflow problem in speed critical code Code which is highly optimized and works with signed integers sometimes has the problem that often the output of the computation does not c
Definition: undefined.txt:32
int nbits
Definition: dct.h:33
#define i(width, name, range_min, range_max)
Definition: cbs_h2645.c:259
DCTTransformType
Definition: avfft.h:93
float FFTSample
Definition: avfft.h:35
Definition: avfft.h:97
static void dct_calc_II_c(DCTContext *ctx, FFTSample *data)
Definition: dct.c:136
void(* dct_calc)(struct DCTContext *s, FFTSample *data)
Definition: dct.h:38
AVFormatContext * ctx
Definition: movenc.c:48
Definition: dct.h:32
#define COS(s, n, x)
Definition: dct.c:41
#define s(width, name)
Definition: cbs_vp9.c:257
int n
Definition: avisynth_c.h:760
void(* rdft_calc)(struct RDFTContext *s, FFTSample *z)
Definition: rdft.h:38
void ff_dct32_float(float *dst, const float *src)
static void dct_calc_III_c(DCTContext *ctx, FFTSample *data)
Definition: dct.c:103
int inverse
Definition: dct.h:34
av_cold int ff_dct_init(DCTContext *s, int nbits, enum DCTTransformType inverse)
Set up DCT.
Definition: dct.c:177
void(* dct32)(FFTSample *out, const FFTSample *in)
Definition: dct.h:39
#define ff_init_ff_cos_tabs
Definition: fft.h:141
Definition: avfft.h:94
FFTSample * csc2
Definition: dct.h:37
static void dst_calc_I_c(DCTContext *ctx, FFTSample *data)
Definition: dct.c:43
#define SIN(s, n, x)
Definition: dct.c:38
av_cold void ff_dct_end(DCTContext *s)
Definition: dct.c:221
#define av_freep(p)
#define M_PI
Definition: mathematics.h:52
static uint32_t inverse(uint32_t v)
find multiplicative inverse modulo 2 ^ 32
Definition: asfcrypt.c:35
#define av_malloc_array(a, b)
Definition: avfft.h:96
Filter the word “frame” indicates either a video frame or a group of audio as stored in an AVFrame structure Format for each input and each output the list of supported formats For video that means pixel format For audio that means channel sample they are references to shared objects When the negotiation mechanism computes the intersection of the formats supported at each end of a all references to both lists are replaced with a reference to the intersection And when a single format is eventually chosen for a link amongst the remaining all references to the list are updated That means that if a filter requires that its input and output have the same format amongst a supported all it has to do is use a reference to the same list of formats query_formats can leave some formats unset and return AVERROR(EAGAIN) to cause the negotiation mechanism toagain later.That can be used by filters with complex requirements to use the format negotiated on one link to set the formats supported on another.Frame references ownership and permissions
av_cold int ff_rdft_init(RDFTContext *s, int nbits, enum RDFTransformType trans)
Set up a real FFT.
Definition: rdft.c:88