FFmpeg
rdft.c
Go to the documentation of this file.
1 /*
2  * (I)RDFT transforms
3  * Copyright (c) 2009 Alex Converse <alex dot converse at gmail dot com>
4  *
5  * This file is part of FFmpeg.
6  *
7  * FFmpeg is free software; you can redistribute it and/or
8  * modify it under the terms of the GNU Lesser General Public
9  * License as published by the Free Software Foundation; either
10  * version 2.1 of the License, or (at your option) any later version.
11  *
12  * FFmpeg is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
15  * Lesser General Public License for more details.
16  *
17  * You should have received a copy of the GNU Lesser General Public
18  * License along with FFmpeg; if not, write to the Free Software
19  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
20  */
21 #include <stdlib.h>
22 #include <math.h>
23 #include "libavutil/mathematics.h"
24 #include "rdft.h"
25 
26 /**
27  * @file
28  * (Inverse) Real Discrete Fourier Transforms.
29  */
30 
31 /** Map one real FFT into two parallel real even and odd FFTs. Then interleave
32  * the two real FFTs into one complex FFT. Unmangle the results.
33  * ref: http://www.engineeringproductivitytools.com/stuff/T0001/PT10.HTM
34  */
36 {
37  int i, i1, i2;
38  FFTComplex ev, od, odsum;
39  const int n = 1 << s->nbits;
40  const float k1 = 0.5;
41  const float k2 = 0.5 - s->inverse;
42  const FFTSample *tcos = s->tcos;
43  const FFTSample *tsin = s->tsin;
44 
45  if (!s->inverse) {
46  s->fft.fft_permute(&s->fft, (FFTComplex*)data);
47  s->fft.fft_calc(&s->fft, (FFTComplex*)data);
48  }
49  /* i=0 is a special case because of packing, the DC term is real, so we
50  are going to throw the N/2 term (also real) in with it. */
51  ev.re = data[0];
52  data[0] = ev.re+data[1];
53  data[1] = ev.re-data[1];
54 
55 #define RDFT_UNMANGLE(sign0, sign1) \
56  for (i = 1; i < (n>>2); i++) { \
57  i1 = 2*i; \
58  i2 = n-i1; \
59  /* Separate even and odd FFTs */ \
60  ev.re = k1*(data[i1 ]+data[i2 ]); \
61  od.im = k2*(data[i2 ]-data[i1 ]); \
62  ev.im = k1*(data[i1+1]-data[i2+1]); \
63  od.re = k2*(data[i1+1]+data[i2+1]); \
64  /* Apply twiddle factors to the odd FFT and add to the even FFT */ \
65  odsum.re = od.re*tcos[i] sign0 od.im*tsin[i]; \
66  odsum.im = od.im*tcos[i] sign1 od.re*tsin[i]; \
67  data[i1 ] = ev.re + odsum.re; \
68  data[i1+1] = ev.im + odsum.im; \
69  data[i2 ] = ev.re - odsum.re; \
70  data[i2+1] = odsum.im - ev.im; \
71  }
72 
73  if (s->negative_sin) {
74  RDFT_UNMANGLE(+,-)
75  } else {
76  RDFT_UNMANGLE(-,+)
77  }
78 
79  data[2*i+1]=s->sign_convention*data[2*i+1];
80  if (s->inverse) {
81  data[0] *= k1;
82  data[1] *= k1;
83  s->fft.fft_permute(&s->fft, (FFTComplex*)data);
84  s->fft.fft_calc(&s->fft, (FFTComplex*)data);
85  }
86 }
87 
88 av_cold int ff_rdft_init(RDFTContext *s, int nbits, enum RDFTransformType trans)
89 {
90  int n = 1 << nbits;
91  int ret;
92 
93  s->nbits = nbits;
94  s->inverse = trans == IDFT_C2R || trans == DFT_C2R;
95  s->sign_convention = trans == IDFT_R2C || trans == DFT_C2R ? 1 : -1;
96  s->negative_sin = trans == DFT_C2R || trans == DFT_R2C;
97 
98  if (nbits < 4 || nbits > 16)
99  return AVERROR(EINVAL);
100 
101  if ((ret = ff_fft_init(&s->fft, nbits-1, trans == IDFT_C2R || trans == IDFT_R2C)) < 0)
102  return ret;
103 
104  ff_init_ff_cos_tabs(nbits);
105  s->tcos = ff_cos_tabs[nbits];
106  s->tsin = ff_cos_tabs[nbits] + (n >> 2);
107  s->rdft_calc = rdft_calc_c;
108 
109  if (ARCH_ARM) ff_rdft_init_arm(s);
110 
111  return 0;
112 }
113 
115 {
116  ff_fft_end(&s->fft);
117 }
ff_fft_init
#define ff_fft_init
Definition: fft.h:149
ff_init_ff_cos_tabs
#define ff_init_ff_cos_tabs
Definition: fft.h:141
AVERROR
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
rdft.h
DFT_C2R
@ DFT_C2R
Definition: avfft.h:75
ff_rdft_init_arm
av_cold void ff_rdft_init_arm(RDFTContext *s)
Definition: rdft_init_arm.c:27
data
const char data[16]
Definition: mxf.c:91
mathematics.h
IDFT_C2R
@ IDFT_C2R
Definition: avfft.h:73
ff_rdft_end
av_cold void ff_rdft_end(RDFTContext *s)
Definition: rdft.c:114
rdft_calc_c
static void rdft_calc_c(RDFTContext *s, FFTSample *data)
Map one real FFT into two parallel real even and odd FFTs.
Definition: rdft.c:35
av_cold
#define av_cold
Definition: attributes.h:90
s
#define s(width, name)
Definition: cbs_vp9.c:257
ff_fft_end
#define ff_fft_end
Definition: fft.h:150
IDFT_R2C
@ IDFT_R2C
Definition: avfft.h:74
DFT_R2C
@ DFT_R2C
Definition: avfft.h:72
FFTSample
float FFTSample
Definition: avfft.h:35
RDFT_UNMANGLE
#define RDFT_UNMANGLE(sign0, sign1)
FFTComplex::re
FFTSample re
Definition: avfft.h:38
ff_rdft_init
av_cold int ff_rdft_init(RDFTContext *s, int nbits, enum RDFTransformType trans)
Set up a real FFT.
Definition: rdft.c:88
i
#define i(width, name, range_min, range_max)
Definition: cbs_h2645.c:269
RDFTContext
Definition: rdft.h:28
RDFTransformType
RDFTransformType
Definition: avfft.h:71
ret
ret
Definition: filter_design.txt:187
FFTComplex
Definition: avfft.h:37