Mercurial > libavcodec.hg
annotate fft-test.c @ 2497:69adfbbdcdeb libavcodec
- samples from mplayer ftp in the "adv" profile seem to have profile=2,
which isn't the advanced one; and indeed, using adv. profile parser fails.
Using normal parser works, and that's what is done
- attempt at taking care of stride for NORM2 bitplane decoding
- duplication of much code from msmpeg4.c; this code isn't yet used, but
goes down as far as the block layer (mainly Transform Type stuff, the
remains are wild editing without checking). Unusable yet, and lacks the AC
decoding (but a step further in bitstream parsing)
patch by anonymous
author | michael |
---|---|
date | Fri, 04 Feb 2005 02:20:38 +0000 |
parents | dd63cb7e5080 |
children | 786ccf72ccd5 |
rev | line source |
---|---|
1106 | 1 /** |
2 * @file fft-test.c | |
3 * FFT and MDCT tests. | |
4 */ | |
5 | |
781 | 6 #include "dsputil.h" |
7 #include <math.h> | |
980 | 8 #include <unistd.h> |
781 | 9 #include <sys/time.h> |
10 | |
11 int mm_flags; | |
12 | |
13 /* reference fft */ | |
14 | |
15 #define MUL16(a,b) ((a) * (b)) | |
16 | |
17 #define CMAC(pre, pim, are, aim, bre, bim) \ | |
18 {\ | |
19 pre += (MUL16(are, bre) - MUL16(aim, bim));\ | |
20 pim += (MUL16(are, bim) + MUL16(bre, aim));\ | |
21 } | |
22 | |
23 FFTComplex *exptab; | |
24 | |
25 void fft_ref_init(int nbits, int inverse) | |
26 { | |
27 int n, i; | |
28 float c1, s1, alpha; | |
29 | |
30 n = 1 << nbits; | |
31 exptab = av_malloc((n / 2) * sizeof(FFTComplex)); | |
32 | |
33 for(i=0;i<(n/2);i++) { | |
34 alpha = 2 * M_PI * (float)i / (float)n; | |
35 c1 = cos(alpha); | |
36 s1 = sin(alpha); | |
37 if (!inverse) | |
38 s1 = -s1; | |
39 exptab[i].re = c1; | |
40 exptab[i].im = s1; | |
41 } | |
42 } | |
43 | |
44 void fft_ref(FFTComplex *tabr, FFTComplex *tab, int nbits) | |
45 { | |
46 int n, i, j, k, n2; | |
47 float tmp_re, tmp_im, s, c; | |
48 FFTComplex *q; | |
49 | |
50 n = 1 << nbits; | |
51 n2 = n >> 1; | |
52 for(i=0;i<n;i++) { | |
53 tmp_re = 0; | |
54 tmp_im = 0; | |
55 q = tab; | |
56 for(j=0;j<n;j++) { | |
57 k = (i * j) & (n - 1); | |
58 if (k >= n2) { | |
59 c = -exptab[k - n2].re; | |
60 s = -exptab[k - n2].im; | |
61 } else { | |
62 c = exptab[k].re; | |
63 s = exptab[k].im; | |
64 } | |
65 CMAC(tmp_re, tmp_im, c, s, q->re, q->im); | |
66 q++; | |
67 } | |
68 tabr[i].re = tmp_re; | |
69 tabr[i].im = tmp_im; | |
70 } | |
71 } | |
72 | |
73 void imdct_ref(float *out, float *in, int n) | |
74 { | |
75 int k, i, a; | |
76 float sum, f; | |
77 | |
78 for(i=0;i<n;i++) { | |
79 sum = 0; | |
80 for(k=0;k<n/2;k++) { | |
81 a = (2 * i + 1 + (n / 2)) * (2 * k + 1); | |
82 f = cos(M_PI * a / (double)(2 * n)); | |
83 sum += f * in[k]; | |
84 } | |
85 out[i] = -sum; | |
86 } | |
87 } | |
88 | |
89 /* NOTE: no normalisation by 1 / N is done */ | |
90 void mdct_ref(float *output, float *input, int n) | |
91 { | |
92 int k, i; | |
93 float a, s; | |
94 | |
95 /* do it by hand */ | |
96 for(k=0;k<n/2;k++) { | |
97 s = 0; | |
98 for(i=0;i<n;i++) { | |
99 a = (2*M_PI*(2*i+1+n/2)*(2*k+1) / (4 * n)); | |
100 s += input[i] * cos(a); | |
101 } | |
102 output[k] = s; | |
103 } | |
104 } | |
105 | |
106 | |
107 float frandom(void) | |
108 { | |
109 return (float)((random() & 0xffff) - 32768) / 32768.0; | |
110 } | |
111 | |
1064 | 112 int64_t gettime(void) |
781 | 113 { |
114 struct timeval tv; | |
115 gettimeofday(&tv,NULL); | |
1064 | 116 return (int64_t)tv.tv_sec * 1000000 + tv.tv_usec; |
781 | 117 } |
118 | |
119 void check_diff(float *tab1, float *tab2, int n) | |
120 { | |
121 int i; | |
122 | |
123 for(i=0;i<n;i++) { | |
124 if (fabsf(tab1[i] - tab2[i]) >= 1e-3) { | |
125 printf("ERROR %d: %f %f\n", | |
126 i, tab1[i], tab2[i]); | |
127 } | |
128 } | |
129 } | |
130 | |
131 | |
132 void help(void) | |
133 { | |
134 printf("usage: fft-test [-h] [-s] [-i] [-n b]\n" | |
135 "-h print this help\n" | |
136 "-s speed test\n" | |
137 "-m (I)MDCT test\n" | |
138 "-i inverse transform test\n" | |
139 "-n b set the transform size to 2^b\n" | |
140 ); | |
141 exit(1); | |
142 } | |
143 | |
144 | |
145 | |
146 int main(int argc, char **argv) | |
147 { | |
148 FFTComplex *tab, *tab1, *tab_ref; | |
149 FFTSample *tabtmp, *tab2; | |
150 int it, i, c; | |
151 int do_speed = 0; | |
152 int do_mdct = 0; | |
153 int do_inverse = 0; | |
154 FFTContext s1, *s = &s1; | |
155 MDCTContext m1, *m = &m1; | |
156 int fft_nbits, fft_size; | |
157 | |
158 mm_flags = 0; | |
159 fft_nbits = 9; | |
160 for(;;) { | |
161 c = getopt(argc, argv, "hsimn:"); | |
162 if (c == -1) | |
163 break; | |
164 switch(c) { | |
165 case 'h': | |
166 help(); | |
167 break; | |
168 case 's': | |
169 do_speed = 1; | |
170 break; | |
171 case 'i': | |
172 do_inverse = 1; | |
173 break; | |
174 case 'm': | |
175 do_mdct = 1; | |
176 break; | |
177 case 'n': | |
178 fft_nbits = atoi(optarg); | |
179 break; | |
180 } | |
181 } | |
182 | |
183 fft_size = 1 << fft_nbits; | |
184 tab = av_malloc(fft_size * sizeof(FFTComplex)); | |
185 tab1 = av_malloc(fft_size * sizeof(FFTComplex)); | |
186 tab_ref = av_malloc(fft_size * sizeof(FFTComplex)); | |
187 tabtmp = av_malloc(fft_size / 2 * sizeof(FFTSample)); | |
188 tab2 = av_malloc(fft_size * sizeof(FFTSample)); | |
189 | |
190 if (do_mdct) { | |
191 if (do_inverse) | |
192 printf("IMDCT"); | |
193 else | |
194 printf("MDCT"); | |
969 | 195 ff_mdct_init(m, fft_nbits, do_inverse); |
781 | 196 } else { |
197 if (do_inverse) | |
198 printf("IFFT"); | |
199 else | |
200 printf("FFT"); | |
1879
dd63cb7e5080
fft_*() renamed into ff_fft_*() patch by (Gildas Bazin <gbazin at altern dot org>)
michael
parents:
1106
diff
changeset
|
201 ff_fft_init(s, fft_nbits, do_inverse); |
781 | 202 fft_ref_init(fft_nbits, do_inverse); |
203 } | |
204 printf(" %d test\n", fft_size); | |
205 | |
206 /* generate random data */ | |
207 | |
208 for(i=0;i<fft_size;i++) { | |
209 tab1[i].re = frandom(); | |
210 tab1[i].im = frandom(); | |
211 } | |
212 | |
213 /* checking result */ | |
214 printf("Checking...\n"); | |
215 | |
216 if (do_mdct) { | |
217 if (do_inverse) { | |
218 imdct_ref((float *)tab_ref, (float *)tab1, fft_size); | |
969 | 219 ff_imdct_calc(m, tab2, (float *)tab1, tabtmp); |
781 | 220 check_diff((float *)tab_ref, tab2, fft_size); |
221 } else { | |
222 mdct_ref((float *)tab_ref, (float *)tab1, fft_size); | |
223 | |
969 | 224 ff_mdct_calc(m, tab2, (float *)tab1, tabtmp); |
781 | 225 |
226 check_diff((float *)tab_ref, tab2, fft_size / 2); | |
227 } | |
228 } else { | |
229 memcpy(tab, tab1, fft_size * sizeof(FFTComplex)); | |
1879
dd63cb7e5080
fft_*() renamed into ff_fft_*() patch by (Gildas Bazin <gbazin at altern dot org>)
michael
parents:
1106
diff
changeset
|
230 ff_fft_permute(s, tab); |
dd63cb7e5080
fft_*() renamed into ff_fft_*() patch by (Gildas Bazin <gbazin at altern dot org>)
michael
parents:
1106
diff
changeset
|
231 ff_fft_calc(s, tab); |
781 | 232 |
233 fft_ref(tab_ref, tab1, fft_nbits); | |
234 check_diff((float *)tab_ref, (float *)tab, fft_size * 2); | |
235 } | |
236 | |
237 /* do a speed test */ | |
238 | |
239 if (do_speed) { | |
1064 | 240 int64_t time_start, duration; |
781 | 241 int nb_its; |
242 | |
243 printf("Speed test...\n"); | |
244 /* we measure during about 1 seconds */ | |
245 nb_its = 1; | |
246 for(;;) { | |
247 time_start = gettime(); | |
248 for(it=0;it<nb_its;it++) { | |
249 if (do_mdct) { | |
250 if (do_inverse) { | |
969 | 251 ff_imdct_calc(m, (float *)tab, (float *)tab1, tabtmp); |
781 | 252 } else { |
969 | 253 ff_mdct_calc(m, (float *)tab, (float *)tab1, tabtmp); |
781 | 254 } |
255 } else { | |
256 memcpy(tab, tab1, fft_size * sizeof(FFTComplex)); | |
1879
dd63cb7e5080
fft_*() renamed into ff_fft_*() patch by (Gildas Bazin <gbazin at altern dot org>)
michael
parents:
1106
diff
changeset
|
257 ff_fft_calc(s, tab); |
781 | 258 } |
259 } | |
260 duration = gettime() - time_start; | |
261 if (duration >= 1000000) | |
262 break; | |
263 nb_its *= 2; | |
264 } | |
265 printf("time: %0.1f us/transform [total time=%0.2f s its=%d]\n", | |
266 (double)duration / nb_its, | |
267 (double)duration / 1000000.0, | |
268 nb_its); | |
269 } | |
270 | |
271 if (do_mdct) { | |
969 | 272 ff_mdct_end(m); |
781 | 273 } else { |
1879
dd63cb7e5080
fft_*() renamed into ff_fft_*() patch by (Gildas Bazin <gbazin at altern dot org>)
michael
parents:
1106
diff
changeset
|
274 ff_fft_end(s); |
781 | 275 } |
276 return 0; | |
277 } |