1// Copyright 2012 The Go Authors. All rights reserved.
2// Use of this source code is governed by a BSD-style
3// license that can be found in the LICENSE file.
4
5package jpeg
6
7import (
8	"bytes"
9	"fmt"
10	"math"
11	"math/rand"
12	"testing"
13)
14
15func benchmarkDCT(b *testing.B, f func(*block)) {
16	b.StopTimer()
17	blocks := make([]block, 0, b.N*len(testBlocks))
18	for i := 0; i < b.N; i++ {
19		blocks = append(blocks, testBlocks[:]...)
20	}
21	b.StartTimer()
22	for i := range blocks {
23		f(&blocks[i])
24	}
25}
26
27func BenchmarkFDCT(b *testing.B) {
28	benchmarkDCT(b, fdct)
29}
30
31func BenchmarkIDCT(b *testing.B) {
32	benchmarkDCT(b, idct)
33}
34
35func TestDCT(t *testing.T) {
36	blocks := make([]block, len(testBlocks))
37	copy(blocks, testBlocks[:])
38
39	// Append some randomly generated blocks of varying sparseness.
40	r := rand.New(rand.NewSource(123))
41	for i := 0; i < 100; i++ {
42		b := block{}
43		n := r.Int() % 64
44		for j := 0; j < n; j++ {
45			b[r.Int()%len(b)] = r.Int31() % 256
46		}
47		blocks = append(blocks, b)
48	}
49
50	// Check that the FDCT and IDCT functions are inverses, after a scale and
51	// level shift. Scaling reduces the rounding errors in the conversion from
52	// floats to ints.
53	for i, b := range blocks {
54		got, want := b, b
55		for j := range got {
56			got[j] = (got[j] - 128) * 8
57		}
58		slowFDCT(&got)
59		slowIDCT(&got)
60		for j := range got {
61			got[j] = got[j]/8 + 128
62		}
63		if differ(&got, &want) {
64			t.Errorf("i=%d: IDCT(FDCT)\nsrc\n%s\ngot\n%s\nwant\n%s\n", i, &b, &got, &want)
65		}
66	}
67
68	// Check that the optimized and slow FDCT implementations agree.
69	// The fdct function already does a scale and level shift.
70	for i, b := range blocks {
71		got, want := b, b
72		fdct(&got)
73		for j := range want {
74			want[j] = (want[j] - 128) * 8
75		}
76		slowFDCT(&want)
77		if differ(&got, &want) {
78			t.Errorf("i=%d: FDCT\nsrc\n%s\ngot\n%s\nwant\n%s\n", i, &b, &got, &want)
79		}
80	}
81
82	// Check that the optimized and slow IDCT implementations agree.
83	for i, b := range blocks {
84		got, want := b, b
85		idct(&got)
86		slowIDCT(&want)
87		if differ(&got, &want) {
88			t.Errorf("i=%d: IDCT\nsrc\n%s\ngot\n%s\nwant\n%s\n", i, &b, &got, &want)
89		}
90	}
91}
92
93// differ reports whether any pair-wise elements in b0 and b1 differ by 2 or
94// more. That tolerance is because there isn't a single definitive decoding of
95// a given JPEG image, even before the YCbCr to RGB conversion; implementations
96// can have different IDCT rounding errors.
97func differ(b0, b1 *block) bool {
98	for i := range b0 {
99		delta := b0[i] - b1[i]
100		if delta < -2 || +2 < delta {
101			return true
102		}
103	}
104	return false
105}
106
107// alpha returns 1 if i is 0 and returns √2 otherwise.
108func alpha(i int) float64 {
109	if i == 0 {
110		return 1
111	}
112	return math.Sqrt2
113}
114
115var cosines [32]float64 // cosines[k] = cos(π/2 * k/8)
116
117func init() {
118	for k := range cosines {
119		cosines[k] = math.Cos(math.Pi * float64(k) / 16)
120	}
121}
122
123// slowFDCT performs the 8*8 2-dimensional forward discrete cosine transform:
124//
125//	dst[u,v] = (1/8) * Σ_x Σ_y alpha(u) * alpha(v) * src[x,y] *
126//		cos((π/2) * (2*x + 1) * u / 8) *
127//		cos((π/2) * (2*y + 1) * v / 8)
128//
129// x and y are in pixel space, and u and v are in transform space.
130//
131// b acts as both dst and src.
132func slowFDCT(b *block) {
133	var dst [blockSize]float64
134	for v := 0; v < 8; v++ {
135		for u := 0; u < 8; u++ {
136			sum := 0.0
137			for y := 0; y < 8; y++ {
138				for x := 0; x < 8; x++ {
139					sum += alpha(u) * alpha(v) * float64(b[8*y+x]) *
140						cosines[((2*x+1)*u)%32] *
141						cosines[((2*y+1)*v)%32]
142				}
143			}
144			dst[8*v+u] = sum / 8
145		}
146	}
147	// Convert from float64 to int32.
148	for i := range dst {
149		b[i] = int32(dst[i] + 0.5)
150	}
151}
152
153// slowIDCT performs the 8*8 2-dimensional inverse discrete cosine transform:
154//
155//	dst[x,y] = (1/8) * Σ_u Σ_v alpha(u) * alpha(v) * src[u,v] *
156//		cos((π/2) * (2*x + 1) * u / 8) *
157//		cos((π/2) * (2*y + 1) * v / 8)
158//
159// x and y are in pixel space, and u and v are in transform space.
160//
161// b acts as both dst and src.
162func slowIDCT(b *block) {
163	var dst [blockSize]float64
164	for y := 0; y < 8; y++ {
165		for x := 0; x < 8; x++ {
166			sum := 0.0
167			for v := 0; v < 8; v++ {
168				for u := 0; u < 8; u++ {
169					sum += alpha(u) * alpha(v) * float64(b[8*v+u]) *
170						cosines[((2*x+1)*u)%32] *
171						cosines[((2*y+1)*v)%32]
172				}
173			}
174			dst[8*y+x] = sum / 8
175		}
176	}
177	// Convert from float64 to int32.
178	for i := range dst {
179		b[i] = int32(dst[i] + 0.5)
180	}
181}
182
183func (b *block) String() string {
184	s := bytes.NewBuffer(nil)
185	fmt.Fprintf(s, "{\n")
186	for y := 0; y < 8; y++ {
187		fmt.Fprintf(s, "\t")
188		for x := 0; x < 8; x++ {
189			fmt.Fprintf(s, "0x%04x, ", uint16(b[8*y+x]))
190		}
191		fmt.Fprintln(s)
192	}
193	fmt.Fprintf(s, "}")
194	return s.String()
195}
196
197// testBlocks are the first 10 pre-IDCT blocks from ../testdata/video-001.jpeg.
198var testBlocks = [10]block{
199	{
200		0x7f, 0xf6, 0x01, 0x07, 0xff, 0x00, 0x00, 0x00,
201		0xf5, 0x01, 0xfa, 0x01, 0xfe, 0x00, 0x01, 0x00,
202		0x05, 0x05, 0x01, 0x00, 0x00, 0x00, 0x00, 0x00,
203		0x01, 0xff, 0xf8, 0x00, 0x01, 0xff, 0x00, 0x00,
204		0x00, 0x01, 0x00, 0x01, 0x00, 0xff, 0xff, 0x00,
205		0xff, 0x0c, 0x00, 0x00, 0x00, 0x00, 0xff, 0x01,
206		0x00, 0x00, 0x00, 0x01, 0x00, 0x00, 0x00, 0x00,
207		0x01, 0x00, 0x00, 0x01, 0xff, 0x01, 0x00, 0xfe,
208	},
209	{
210		0x29, 0x07, 0x00, 0xfc, 0x01, 0x01, 0x00, 0x00,
211		0x07, 0x00, 0x03, 0x00, 0x01, 0x00, 0xff, 0xff,
212		0xff, 0xfd, 0xff, 0x00, 0x00, 0x00, 0x00, 0x00,
213		0x00, 0x00, 0x04, 0x00, 0xff, 0x01, 0x00, 0x00,
214		0x01, 0x00, 0x01, 0xff, 0x00, 0x00, 0x00, 0x00,
215		0x01, 0xfa, 0x01, 0x00, 0x01, 0x00, 0x01, 0xff,
216		0x00, 0x00, 0xff, 0x00, 0x00, 0x00, 0x00, 0x00,
217		0x00, 0x00, 0x00, 0xff, 0x00, 0xff, 0x00, 0x02,
218	},
219	{
220		0xc5, 0xfa, 0x01, 0x00, 0x00, 0x01, 0x00, 0xff,
221		0x02, 0xff, 0x01, 0x00, 0x01, 0x00, 0xff, 0x00,
222		0xff, 0xff, 0x00, 0xff, 0x01, 0x00, 0x00, 0x00,
223		0xff, 0x00, 0x01, 0x00, 0x00, 0x00, 0xff, 0x00,
224		0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0xff,
225		0x00, 0xff, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
226		0x01, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
227		0xff, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
228	},
229	{
230		0x86, 0x05, 0x00, 0x02, 0x00, 0x00, 0x01, 0x00,
231		0xf2, 0x06, 0x00, 0x00, 0x01, 0x02, 0x00, 0x00,
232		0xf6, 0xfa, 0xf9, 0x00, 0xff, 0x01, 0x00, 0x00,
233		0xf9, 0x00, 0x00, 0xff, 0x00, 0x00, 0x00, 0x00,
234		0x00, 0xff, 0x00, 0xff, 0xff, 0xff, 0x00, 0x00,
235		0xff, 0x00, 0x00, 0x01, 0x00, 0xff, 0x01, 0x00,
236		0x00, 0x00, 0x00, 0xff, 0x00, 0x00, 0x00, 0x01,
237		0x00, 0x01, 0xff, 0x01, 0x00, 0xff, 0x00, 0x00,
238	},
239	{
240		0x24, 0xfe, 0x00, 0xff, 0x00, 0xff, 0xff, 0x00,
241		0x08, 0xfd, 0x00, 0x01, 0x01, 0x00, 0x01, 0x00,
242		0x06, 0x03, 0x03, 0xff, 0x00, 0x00, 0x00, 0x00,
243		0x04, 0xff, 0x00, 0x00, 0x00, 0x00, 0x00, 0xff,
244		0x00, 0x00, 0x00, 0x00, 0x01, 0x00, 0x00, 0x01,
245		0x01, 0x00, 0x01, 0xff, 0x00, 0x01, 0x00, 0x00,
246		0x01, 0x01, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
247		0x01, 0x00, 0x01, 0x00, 0x00, 0x00, 0xff, 0x01,
248	},
249	{
250		0xcd, 0xff, 0x00, 0x00, 0x00, 0x00, 0x01, 0x01,
251		0x03, 0xff, 0x00, 0x00, 0x00, 0x00, 0x00, 0xff,
252		0x01, 0x01, 0x01, 0x01, 0x01, 0x00, 0x00, 0x00,
253		0x01, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
254		0x01, 0x00, 0x00, 0x00, 0x00, 0x01, 0x01, 0x00,
255		0x00, 0x00, 0x00, 0x01, 0x00, 0x00, 0x00, 0x00,
256		0x00, 0xff, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
257		0x00, 0x00, 0x00, 0x00, 0x00, 0x01, 0x00, 0xff,
258	},
259	{
260		0x81, 0xfe, 0x05, 0xff, 0x01, 0xff, 0x01, 0x00,
261		0xef, 0xf9, 0x00, 0xf9, 0x00, 0xff, 0x00, 0xff,
262		0x05, 0xf9, 0x00, 0xf8, 0x01, 0xff, 0x01, 0xff,
263		0x00, 0xff, 0x07, 0x00, 0x01, 0x00, 0x00, 0x00,
264		0x01, 0x00, 0x01, 0x01, 0x00, 0x00, 0x00, 0x00,
265		0x01, 0x00, 0x00, 0x00, 0xff, 0xff, 0x00, 0x01,
266		0xff, 0x01, 0x01, 0x00, 0xff, 0x00, 0x00, 0x00,
267		0x01, 0x01, 0x00, 0xff, 0x00, 0x00, 0x00, 0xff,
268	},
269	{
270		0x28, 0x00, 0xfe, 0x00, 0x00, 0x00, 0x00, 0x00,
271		0x0b, 0x02, 0x01, 0x03, 0x00, 0xff, 0x00, 0x01,
272		0xfe, 0x02, 0x01, 0x03, 0xff, 0x00, 0x00, 0x00,
273		0x01, 0x00, 0xfd, 0x00, 0x01, 0x00, 0xff, 0x00,
274		0x01, 0xff, 0x00, 0xff, 0x01, 0x00, 0x00, 0x00,
275		0x00, 0x00, 0x00, 0xff, 0x01, 0x01, 0x00, 0xff,
276		0x01, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00,
277		0xff, 0xff, 0x00, 0x00, 0x00, 0xff, 0x00, 0x01,
278	},
279	{
280		0xdf, 0xf9, 0xfe, 0x00, 0x03, 0x01, 0xff, 0xff,
281		0x04, 0x01, 0x00, 0x01, 0x00, 0x00, 0x00, 0x00,
282		0xff, 0x01, 0x01, 0x01, 0x00, 0x00, 0x00, 0x01,
283		0x00, 0x00, 0xfe, 0x01, 0x00, 0x00, 0x00, 0x00,
284		0x00, 0x00, 0xff, 0x01, 0x00, 0x00, 0x00, 0x01,
285		0xff, 0x00, 0x00, 0x00, 0x00, 0x01, 0x00, 0x00,
286		0x00, 0xff, 0x00, 0xff, 0x01, 0x00, 0x00, 0x01,
287		0xff, 0xff, 0x00, 0x00, 0x00, 0x01, 0x00, 0x00,
288	},
289	{
290		0x88, 0xfd, 0x00, 0x00, 0xff, 0x00, 0x01, 0xff,
291		0xe1, 0x06, 0x06, 0x01, 0xff, 0x00, 0x01, 0x00,
292		0x08, 0x00, 0xfa, 0x00, 0xff, 0xff, 0xff, 0xff,
293		0x08, 0x01, 0x00, 0xff, 0x01, 0xff, 0x00, 0x00,
294		0xf5, 0xff, 0x00, 0x01, 0xff, 0x01, 0x01, 0x00,
295		0xff, 0xff, 0x01, 0xff, 0x01, 0x00, 0x01, 0x00,
296		0x00, 0x01, 0x01, 0xff, 0x00, 0xff, 0x00, 0x01,
297		0x02, 0x00, 0x00, 0xff, 0xff, 0x00, 0xff, 0x00,
298	},
299}
300