4 * Copyright (C) 1994-1996, Thomas G. Lane.
5 * This file is part of the Independent JPEG Group's software.
6 * For conditions of distribution and use, see the accompanying README file.
8 * This file contains a floating-point implementation of the
9 * forward DCT (Discrete Cosine Transform).
11 * This implementation should be more accurate than either of the integer
12 * DCT implementations. However, it may not give the same results on all
13 * machines because of differences in roundoff behavior. Speed will depend
14 * on the hardware's floating point capacity.
16 * A 2-D DCT can be done by 1-D DCT on each row followed by 1-D DCT
17 * on each column. Direct algorithms are also available, but they are
18 * much more complex and seem not to be any faster when reduced to code.
20 * This implementation is based on Arai, Agui, and Nakajima's algorithm for
21 * scaled DCT. Their original paper (Trans. IEICE E-71(11):1095) is in
22 * Japanese, but the algorithm is described in the Pennebaker & Mitchell
23 * JPEG textbook (see REFERENCES section in file README). The following code
24 * is based directly on figure 4-8 in P&M.
25 * While an 8-point DCT cannot be done in less than 11 multiplies, it is
26 * possible to arrange the computation so that many of the multiplies are
27 * simple scalings of the final outputs. These multiplies can then be
28 * folded into the multiplications or divisions by the JPEG quantization
29 * table entries. The AA&N method leaves only 5 multiplies and 29 adds
30 * to be done in the DCT itself.
31 * The primary disadvantage of this method is that with a fixed-point
32 * implementation, accuracy is lost due to imprecise representation of the
33 * scaled quantization values. However, that problem does not arise if
34 * we use floating point arithmetic.
37 #define JPEG_INTERNALS
40 #include "jdct.h" /* Private declarations for DCT subsystem */
42 #ifdef DCT_FLOAT_SUPPORTED
46 * This module is specialized to the case DCTSIZE = 8.
50 Sorry
, this code only copes with
8x8 DCTs
. /* deliberate syntax err */
55 * Perform the forward DCT on one block of samples.
59 jpeg_fdct_float (FAST_FLOAT
* data
)
61 FAST_FLOAT tmp0
, tmp1
, tmp2
, tmp3
, tmp4
, tmp5
, tmp6
, tmp7
;
62 FAST_FLOAT tmp10
, tmp11
, tmp12
, tmp13
;
63 FAST_FLOAT z1
, z2
, z3
, z4
, z5
, z11
, z13
;
67 /* Pass 1: process rows. */
70 for (ctr
= DCTSIZE
-1; ctr
>= 0; ctr
--) {
71 tmp0
= dataptr
[0] + dataptr
[7];
72 tmp7
= dataptr
[0] - dataptr
[7];
73 tmp1
= dataptr
[1] + dataptr
[6];
74 tmp6
= dataptr
[1] - dataptr
[6];
75 tmp2
= dataptr
[2] + dataptr
[5];
76 tmp5
= dataptr
[2] - dataptr
[5];
77 tmp3
= dataptr
[3] + dataptr
[4];
78 tmp4
= dataptr
[3] - dataptr
[4];
82 tmp10
= tmp0
+ tmp3
; /* phase 2 */
87 dataptr
[0] = tmp10
+ tmp11
; /* phase 3 */
88 dataptr
[4] = tmp10
- tmp11
;
90 z1
= (tmp12
+ tmp13
) * ((FAST_FLOAT
) 0.707106781); /* c4 */
91 dataptr
[2] = tmp13
+ z1
; /* phase 5 */
92 dataptr
[6] = tmp13
- z1
;
96 tmp10
= tmp4
+ tmp5
; /* phase 2 */
100 /* The rotator is modified from fig 4-8 to avoid extra negations. */
101 z5
= (tmp10
- tmp12
) * ((FAST_FLOAT
) 0.382683433); /* c6 */
102 z2
= ((FAST_FLOAT
) 0.541196100) * tmp10
+ z5
; /* c2-c6 */
103 z4
= ((FAST_FLOAT
) 1.306562965) * tmp12
+ z5
; /* c2+c6 */
104 z3
= tmp11
* ((FAST_FLOAT
) 0.707106781); /* c4 */
106 z11
= tmp7
+ z3
; /* phase 5 */
109 dataptr
[5] = z13
+ z2
; /* phase 6 */
110 dataptr
[3] = z13
- z2
;
111 dataptr
[1] = z11
+ z4
;
112 dataptr
[7] = z11
- z4
;
114 dataptr
+= DCTSIZE
; /* advance pointer to next row */
117 /* Pass 2: process columns. */
120 for (ctr
= DCTSIZE
-1; ctr
>= 0; ctr
--) {
121 tmp0
= dataptr
[DCTSIZE
*0] + dataptr
[DCTSIZE
*7];
122 tmp7
= dataptr
[DCTSIZE
*0] - dataptr
[DCTSIZE
*7];
123 tmp1
= dataptr
[DCTSIZE
*1] + dataptr
[DCTSIZE
*6];
124 tmp6
= dataptr
[DCTSIZE
*1] - dataptr
[DCTSIZE
*6];
125 tmp2
= dataptr
[DCTSIZE
*2] + dataptr
[DCTSIZE
*5];
126 tmp5
= dataptr
[DCTSIZE
*2] - dataptr
[DCTSIZE
*5];
127 tmp3
= dataptr
[DCTSIZE
*3] + dataptr
[DCTSIZE
*4];
128 tmp4
= dataptr
[DCTSIZE
*3] - dataptr
[DCTSIZE
*4];
132 tmp10
= tmp0
+ tmp3
; /* phase 2 */
137 dataptr
[DCTSIZE
*0] = tmp10
+ tmp11
; /* phase 3 */
138 dataptr
[DCTSIZE
*4] = tmp10
- tmp11
;
140 z1
= (tmp12
+ tmp13
) * ((FAST_FLOAT
) 0.707106781); /* c4 */
141 dataptr
[DCTSIZE
*2] = tmp13
+ z1
; /* phase 5 */
142 dataptr
[DCTSIZE
*6] = tmp13
- z1
;
146 tmp10
= tmp4
+ tmp5
; /* phase 2 */
150 /* The rotator is modified from fig 4-8 to avoid extra negations. */
151 z5
= (tmp10
- tmp12
) * ((FAST_FLOAT
) 0.382683433); /* c6 */
152 z2
= ((FAST_FLOAT
) 0.541196100) * tmp10
+ z5
; /* c2-c6 */
153 z4
= ((FAST_FLOAT
) 1.306562965) * tmp12
+ z5
; /* c2+c6 */
154 z3
= tmp11
* ((FAST_FLOAT
) 0.707106781); /* c4 */
156 z11
= tmp7
+ z3
; /* phase 5 */
159 dataptr
[DCTSIZE
*5] = z13
+ z2
; /* phase 6 */
160 dataptr
[DCTSIZE
*3] = z13
- z2
;
161 dataptr
[DCTSIZE
*1] = z11
+ z4
;
162 dataptr
[DCTSIZE
*7] = z11
- z4
;
164 dataptr
++; /* advance pointer to next column */
168 #endif /* DCT_FLOAT_SUPPORTED */