1 /************************************************************************
2 *
3 * idct.c, inverse fast DCT for tmndecode (H.263 decoder)
4 * Copyright (C) 1995, 1996 Telenor R&D, Norway
5 *
6 * Contacts:
7 * Robert Danielsen <Robert.Danielsen@nta.no>
8 *
9 * Telenor Research and Development http://www.nta.no/brukere/DVC/
10 * P.O.Box 83 tel.: +47 63 84 84 00
11 * N-2007 Kjeller, Norway fax.: +47 63 81 00 76
12 *
13 * Copyright (C) 1997 University of BC, Canada
14 * Modified by: Michael Gallant <mikeg@ee.ubc.ca>
15 * Guy Cote <guyc@ee.ubc.ca>
16 * Berna Erol <bernae@ee.ubc.ca>
17 *
18 * Contacts:
19 * Michael Gallant <mikeg@ee.ubc.ca>
20 *
21 * UBC Image Processing Laboratory http://www.ee.ubc.ca/image
22 * 2356 Main Mall tel.: +1 604 822 4051
23 * Vancouver BC Canada V6T1Z4 fax.: +1 604 822 5949
24 *
25 ************************************************************************/
26
27 /* Disclaimer of Warranty
28 *
29 * These software programs are available to the user without any license fee
30 * or royalty on an "as is" basis. The University of British Columbia
31 * disclaims any and all warranties, whether express, implied, or
32 * statuary, including any implied warranties or merchantability or of
33 * fitness for a particular purpose. In no event shall the
34 * copyright-holder be liable for any incidental, punitive, or
35 * consequential damages of any kind whatsoever arising from the use of
36 * these programs.
37 *
38 * This disclaimer of warranty extends to the user of these programs and
39 * user's customers, employees, agents, transferees, successors, and
40 * assigns.
41 *
42 * The University of British Columbia does not represent or warrant that the
43 * programs furnished hereunder are free of infringement of any
44 * third-party patents.
45 *
46 * Commercial implementations of H.263, including shareware, are subject to
47 * royalty fees to patent holders. Many of these patents are general
48 * enough such that they are unavoidable regardless of implementation
49 * design.
50 *
51 */
52
53
54
55 /* based on mpeg2decode, (C) 1994, MPEG Software Simulation Group and
56 * mpeg2play, (C) 1994 Stefan Eckart <stefan@lis.e-technik.tu-muenchen.de>
57 *
58 */
59
60
61 /**********************************************************/
62 /* inverse two dimensional DCT, Chen-Wang algorithm */
63 /* (cf. IEEE ASSP-32, pp. 803-816, Aug. 1984) */
64 /* 32-bit integer arithmetic (8 bit coefficients) */
65 /* 11 mults, 29 adds per DCT */
66 /* sE, 18.8.91 */
67 /**********************************************************/
68 /* coefficients extended to 12 bit for IEEE1180-1990 */
69 /* compliance sE, 2.1.94 */
70 /**********************************************************/
71
72 /* this code assumes >> to be a two's-complement arithmetic */
73 /* right shift: (-2)>>1 == -1 , (-3)>>1 == -2 */
74
75 #include "config.h"
76
77 #define W1 2841 /* 2048*sqrt(2)*cos(1*pi/16) */
78 #define W2 2676 /* 2048*sqrt(2)*cos(2*pi/16) */
79 #define W3 2408 /* 2048*sqrt(2)*cos(3*pi/16) */
80 #define W5 1609 /* 2048*sqrt(2)*cos(5*pi/16) */
81 #define W6 1108 /* 2048*sqrt(2)*cos(6*pi/16) */
82 #define W7 565 /* 2048*sqrt(2)*cos(7*pi/16) */
83
84 /* global declarations */
85 void init_idct _ANSI_ARGS_ ((void));
86 void idct _ANSI_ARGS_ ((short *block));
87
88 /* private data */
89 static short iclip[1024]; /* clipping table */
90 static short *iclp;
91
92 /* private prototypes */
93 static void idctrow _ANSI_ARGS_ ((short *blk));
94 static void idctcol _ANSI_ARGS_ ((short *blk));
95
96 /* row (horizontal) IDCT
97 *
98 * 7 pi 1 dst[k] = sum c[l] * src[l] * cos( -- *
99 * ( k + - ) * l ) l=0 8 2
100 *
101 * where: c[0] = 128 c[1..7] = 128*sqrt(2) */
102
103 static void idctrow (short *blk)
104 {
105 int x0, x1, x2, x3, x4, x5, x6, x7, x8;
106
107 /* shortcut */
108 if (!((x1 = blk[4] << 11) | (x2 = blk[6]) | (x3 = blk[2]) |
109 (x4 = blk[1]) | (x5 = blk[7]) | (x6 = blk[5]) | (x7 = blk[3])))
110 {
111 blk[0] = blk[1] = blk[2] = blk[3] = blk[4] = blk[5] = blk[6] = blk[7] = blk[0] << 3;
112 return;
113 }
114 x0 = (blk[0] << 11) + 128; /* for proper rounding in the fourth stage */
115
116 /* first stage */
117 x8 = W7 * (x4 + x5);
118 x4 = x8 + (W1 - W7) * x4;
119 x5 = x8 - (W1 + W7) * x5;
120 x8 = W3 * (x6 + x7);
121 x6 = x8 - (W3 - W5) * x6;
122 x7 = x8 - (W3 + W5) * x7;
123
124 /* second stage */
125 x8 = x0 + x1;
126 x0 -= x1;
127 x1 = W6 * (x3 + x2);
128 x2 = x1 - (W2 + W6) * x2;
129 x3 = x1 + (W2 - W6) * x3;
130 x1 = x4 + x6;
131 x4 -= x6;
132 x6 = x5 + x7;
133 x5 -= x7;
134
135 /* third stage */
136 x7 = x8 + x3;
137 x8 -= x3;
138 x3 = x0 + x2;
139 x0 -= x2;
140 x2 = (181 * (x4 + x5) + 128) >> 8;
141 x4 = (181 * (x4 - x5) + 128) >> 8;
142
143 /* fourth stage */
144 blk[0] = (x7 + x1) >> 8;
145 blk[1] = (x3 + x2) >> 8;
146 blk[2] = (x0 + x4) >> 8;
147 blk[3] = (x8 + x6) >> 8;
148 blk[4] = (x8 - x6) >> 8;
149 blk[5] = (x0 - x4) >> 8;
150 blk[6] = (x3 - x2) >> 8;
151 blk[7] = (x7 - x1) >> 8;
152 }
153
154 /* column (vertical) IDCT
155 *
156 * 7 pi 1 dst[8*k] = sum c[l] * src[8*l] *
157 * cos( -- * ( k + - ) * l ) l=0 8 2
158 *
159 * where: c[0] = 1/1024 c[1..7] = (1/1024)*sqrt(2) */
160 static void idctcol (short *blk)
161 {
162 int x0, x1, x2, x3, x4, x5, x6, x7, x8;
163
164 /* shortcut */
165 if (!((x1 = (blk[8 * 4] << 8)) | (x2 = blk[8 * 6]) | (x3 = blk[8 * 2]) |
166 (x4 = blk[8 * 1]) | (x5 = blk[8 * 7]) | (x6 = blk[8 * 5]) | (x7 = blk[8 * 3])))
167 {
168 blk[8 * 0] = blk[8 * 1] = blk[8 * 2] = blk[8 * 3] = blk[8 * 4] = blk[8 * 5] = blk[8 * 6] = blk[8 * 7] =
169 iclp[(blk[8 * 0] + 32) >> 6];
170 return;
171 }
172 x0 = (blk[8 * 0] << 8) + 8192;
173
174 /* first stage */
175 x8 = W7 * (x4 + x5) + 4;
176 x4 = (x8 + (W1 - W7) * x4) >> 3;
177 x5 = (x8 - (W1 + W7) * x5) >> 3;
178 x8 = W3 * (x6 + x7) + 4;
179 x6 = (x8 - (W3 - W5) * x6) >> 3;
180 x7 = (x8 - (W3 + W5) * x7) >> 3;
181
182 /* second stage */
183 x8 = x0 + x1;
184 x0 -= x1;
185 x1 = W6 * (x3 + x2) + 4;
186 x2 = (x1 - (W2 + W6) * x2) >> 3;
187 x3 = (x1 + (W2 - W6) * x3) >> 3;
188 x1 = x4 + x6;
189 x4 -= x6;
190 x6 = x5 + x7;
191 x5 -= x7;
192
193 /* third stage */
194 x7 = x8 + x3;
195 x8 -= x3;
196 x3 = x0 + x2;
197 x0 -= x2;
198 x2 = (181 * (x4 + x5) + 128) >> 8;
199 x4 = (181 * (x4 - x5) + 128) >> 8;
200
201 /* fourth stage */
202 blk[8 * 0] = iclp[(x7 + x1) >> 14];
203 blk[8 * 1] = iclp[(x3 + x2) >> 14];
204 blk[8 * 2] = iclp[(x0 + x4) >> 14];
205 blk[8 * 3] = iclp[(x8 + x6) >> 14];
206 blk[8 * 4] = iclp[(x8 - x6) >> 14];
207 blk[8 * 5] = iclp[(x0 - x4) >> 14];
208 blk[8 * 6] = iclp[(x3 - x2) >> 14];
209 blk[8 * 7] = iclp[(x7 - x1) >> 14];
210 }
211
212 /* two dimensional inverse discrete cosine transform */
213 void idct (short *block)
214 {
215 int i;
216
217 for (i = 0; i < 8; i++)
218 idctrow (block + 8 * i);
219
220 for (i = 0; i < 8; i++)
221 idctcol (block + i);
222 }
223
224 void init_idct ()
225 {
226 int i;
227
228 iclp = iclip + 512;
229 for (i = -512; i < 512; i++)
230 iclp[i] = (i < -256) ? -256 : ((i > 255) ? 255 : i);
231 }
232
This page was automatically generated by the
LXR engine.
Visit the LXR main site for more
information.