4 * Copyright (C) 1994-1996, Thomas G. Lane.
5 * Modified 2003-2015 by Guido Vollbeding.
6 * This file is part of the Independent JPEG Group's software.
7 * For conditions of distribution and use, see the accompanying README file.
9 * This file contains a floating-point implementation of the
10 * forward DCT (Discrete Cosine Transform).
12 * This implementation should be more accurate than either of the integer
13 * DCT implementations. However, it may not give the same results on all
14 * machines because of differences in roundoff behavior. Speed will depend
15 * on the hardware's floating point capacity.
17 * A 2-D DCT can be done by 1-D DCT on each row followed by 1-D DCT
18 * on each column. Direct algorithms are also available, but they are
19 * much more complex and seem not to be any faster when reduced to code.
21 * This implementation is based on Arai, Agui, and Nakajima's algorithm for
22 * scaled DCT. Their original paper (Trans. IEICE E-71(11):1095) is in
23 * Japanese, but the algorithm is described in the Pennebaker & Mitchell
24 * JPEG textbook (see REFERENCES section in file README). The following code
25 * is based directly on figure 4-8 in P&M.
26 * While an 8-point DCT cannot be done in less than 11 multiplies, it is
27 * possible to arrange the computation so that many of the multiplies are
28 * simple scalings of the final outputs. These multiplies can then be
29 * folded into the multiplications or divisions by the JPEG quantization
30 * table entries. The AA&N method leaves only 5 multiplies and 29 adds
31 * to be done in the DCT itself.
32 * The primary disadvantage of this method is that with a fixed-point
33 * implementation, accuracy is lost due to imprecise representation of the
34 * scaled quantization values. However, that problem does not arise if
35 * we use floating point arithmetic.
38 #define JPEG_INTERNALS
41 #include "jdct.h" /* Private declarations for DCT subsystem */
43 #ifdef DCT_FLOAT_SUPPORTED
47 * This module is specialized to the case DCTSIZE = 8.
51 Sorry
, this code only copes with
8x8 DCTs
. /* deliberate syntax err */
56 * Perform the forward DCT on one block of samples.
58 * cK represents cos(K*pi/16).
62 jpeg_fdct_float (FAST_FLOAT
* data
, JSAMPARRAY sample_data
, JDIMENSION start_col
)
64 FAST_FLOAT tmp0
, tmp1
, tmp2
, tmp3
, tmp4
, tmp5
, tmp6
, tmp7
;
65 FAST_FLOAT tmp10
, tmp11
, tmp12
, tmp13
;
66 FAST_FLOAT z1
, z2
, z3
, z4
, z5
, z11
, z13
;
71 /* Pass 1: process rows. */
74 for (ctr
= 0; ctr
< DCTSIZE
; ctr
++) {
75 elemptr
= sample_data
[ctr
] + start_col
;
77 /* Load data into workspace */
78 tmp0
= (FAST_FLOAT
) (GETJSAMPLE(elemptr
[0]) + GETJSAMPLE(elemptr
[7]));
79 tmp7
= (FAST_FLOAT
) (GETJSAMPLE(elemptr
[0]) - GETJSAMPLE(elemptr
[7]));
80 tmp1
= (FAST_FLOAT
) (GETJSAMPLE(elemptr
[1]) + GETJSAMPLE(elemptr
[6]));
81 tmp6
= (FAST_FLOAT
) (GETJSAMPLE(elemptr
[1]) - GETJSAMPLE(elemptr
[6]));
82 tmp2
= (FAST_FLOAT
) (GETJSAMPLE(elemptr
[2]) + GETJSAMPLE(elemptr
[5]));
83 tmp5
= (FAST_FLOAT
) (GETJSAMPLE(elemptr
[2]) - GETJSAMPLE(elemptr
[5]));
84 tmp3
= (FAST_FLOAT
) (GETJSAMPLE(elemptr
[3]) + GETJSAMPLE(elemptr
[4]));
85 tmp4
= (FAST_FLOAT
) (GETJSAMPLE(elemptr
[3]) - GETJSAMPLE(elemptr
[4]));
89 tmp10
= tmp0
+ tmp3
; /* phase 2 */
94 /* Apply unsigned->signed conversion. */
95 dataptr
[0] = tmp10
+ tmp11
- 8 * CENTERJSAMPLE
; /* phase 3 */
96 dataptr
[4] = tmp10
- tmp11
;
98 z1
= (tmp12
+ tmp13
) * ((FAST_FLOAT
) 0.707106781); /* c4 */
99 dataptr
[2] = tmp13
+ z1
; /* phase 5 */
100 dataptr
[6] = tmp13
- z1
;
104 tmp10
= tmp4
+ tmp5
; /* phase 2 */
108 /* The rotator is modified from fig 4-8 to avoid extra negations. */
109 z5
= (tmp10
- tmp12
) * ((FAST_FLOAT
) 0.382683433); /* c6 */
110 z2
= ((FAST_FLOAT
) 0.541196100) * tmp10
+ z5
; /* c2-c6 */
111 z4
= ((FAST_FLOAT
) 1.306562965) * tmp12
+ z5
; /* c2+c6 */
112 z3
= tmp11
* ((FAST_FLOAT
) 0.707106781); /* c4 */
114 z11
= tmp7
+ z3
; /* phase 5 */
117 dataptr
[5] = z13
+ z2
; /* phase 6 */
118 dataptr
[3] = z13
- z2
;
119 dataptr
[1] = z11
+ z4
;
120 dataptr
[7] = z11
- z4
;
122 dataptr
+= DCTSIZE
; /* advance pointer to next row */
125 /* Pass 2: process columns. */
128 for (ctr
= DCTSIZE
-1; ctr
>= 0; ctr
--) {
129 tmp0
= dataptr
[DCTSIZE
*0] + dataptr
[DCTSIZE
*7];
130 tmp7
= dataptr
[DCTSIZE
*0] - dataptr
[DCTSIZE
*7];
131 tmp1
= dataptr
[DCTSIZE
*1] + dataptr
[DCTSIZE
*6];
132 tmp6
= dataptr
[DCTSIZE
*1] - dataptr
[DCTSIZE
*6];
133 tmp2
= dataptr
[DCTSIZE
*2] + dataptr
[DCTSIZE
*5];
134 tmp5
= dataptr
[DCTSIZE
*2] - dataptr
[DCTSIZE
*5];
135 tmp3
= dataptr
[DCTSIZE
*3] + dataptr
[DCTSIZE
*4];
136 tmp4
= dataptr
[DCTSIZE
*3] - dataptr
[DCTSIZE
*4];
140 tmp10
= tmp0
+ tmp3
; /* phase 2 */
145 dataptr
[DCTSIZE
*0] = tmp10
+ tmp11
; /* phase 3 */
146 dataptr
[DCTSIZE
*4] = tmp10
- tmp11
;
148 z1
= (tmp12
+ tmp13
) * ((FAST_FLOAT
) 0.707106781); /* c4 */
149 dataptr
[DCTSIZE
*2] = tmp13
+ z1
; /* phase 5 */
150 dataptr
[DCTSIZE
*6] = tmp13
- z1
;
154 tmp10
= tmp4
+ tmp5
; /* phase 2 */
158 /* The rotator is modified from fig 4-8 to avoid extra negations. */
159 z5
= (tmp10
- tmp12
) * ((FAST_FLOAT
) 0.382683433); /* c6 */
160 z2
= ((FAST_FLOAT
) 0.541196100) * tmp10
+ z5
; /* c2-c6 */
161 z4
= ((FAST_FLOAT
) 1.306562965) * tmp12
+ z5
; /* c2+c6 */
162 z3
= tmp11
* ((FAST_FLOAT
) 0.707106781); /* c4 */
164 z11
= tmp7
+ z3
; /* phase 5 */
167 dataptr
[DCTSIZE
*5] = z13
+ z2
; /* phase 6 */
168 dataptr
[DCTSIZE
*3] = z13
- z2
;
169 dataptr
[DCTSIZE
*1] = z11
+ z4
;
170 dataptr
[DCTSIZE
*7] = z11
- z4
;
172 dataptr
++; /* advance pointer to next column */
176 #endif /* DCT_FLOAT_SUPPORTED */