1 /*
2  * jidctint.c
3  *
4  * This file was part of the Independent JPEG Group's software.
5  * Copyright (C) 1991-1998, Thomas G. Lane.
6  * Modification developed 2002-2009 by Guido Vollbeding.
7  * libjpeg-turbo Modifications:
8  * Copyright (C) 2015, D. R. Commander
9  * For conditions of distribution and use, see the accompanying README file.
10  *
11  * This file contains a slow-but-accurate integer implementation of the
12  * inverse DCT (Discrete Cosine Transform).  In the IJG code, this routine
13  * must also perform dequantization of the input coefficients.
14  *
15  * A 2-D IDCT can be done by 1-D IDCT on each column followed by 1-D IDCT
16  * on each row (or vice versa, but it's more convenient to emit a row at
17  * a time).  Direct algorithms are also available, but they are much more
18  * complex and seem not to be any faster when reduced to code.
19  *
20  * This implementation is based on an algorithm described in
21  *   C. Loeffler, A. Ligtenberg and G. Moschytz, "Practical Fast 1-D DCT
22  *   Algorithms with 11 Multiplications", Proc. Int'l. Conf. on Acoustics,
23  *   Speech, and Signal Processing 1989 (ICASSP '89), pp. 988-991.
24  * The primary algorithm described there uses 11 multiplies and 29 adds.
25  * We use their alternate method with 12 multiplies and 32 adds.
26  * The advantage of this method is that no data path contains more than one
27  * multiplication; this allows a very simple and accurate implementation in
28  * scaled fixed-point arithmetic, with a minimal number of shifts.
29  *
30  * We also provide IDCT routines with various output sample block sizes for
31  * direct resolution reduction or enlargement without additional resampling:
32  * NxN (N=1...16) pixels for one 8x8 input DCT block.
33  *
34  * For N<8 we simply take the corresponding low-frequency coefficients of
35  * the 8x8 input DCT block and apply an NxN point IDCT on the sub-block
36  * to yield the downscaled outputs.
37  * This can be seen as direct low-pass downsampling from the DCT domain
38  * point of view rather than the usual spatial domain point of view,
39  * yielding significant computational savings and results at least
40  * as good as common bilinear (averaging) spatial downsampling.
41  *
42  * For N>8 we apply a partial NxN IDCT on the 8 input coefficients as
43  * lower frequencies and higher frequencies assumed to be zero.
44  * It turns out that the computational effort is similar to the 8x8 IDCT
45  * regarding the output size.
46  * Furthermore, the scaling and descaling is the same for all IDCT sizes.
47  *
48  * CAUTION: We rely on the FIX() macro except for the N=1,2,4,8 cases
49  * since there would be too many additional constants to pre-calculate.
50  */
51 
52 #define JPEG_INTERNALS
53 #include "jinclude.h"
54 #include "jpeglib.h"
55 #include "jdct.h"               /* Private declarations for DCT subsystem */
56 
57 #ifdef DCT_ISLOW_SUPPORTED
58 
59 
60 /*
61  * This module is specialized to the case DCTSIZE = 8.
62  */
63 
64 #if DCTSIZE != 8
65   Sorry, this code only copes with 8x8 DCT blocks. /* deliberate syntax err */
66 #endif
67 
68 
69 /*
70  * The poop on this scaling stuff is as follows:
71  *
72  * Each 1-D IDCT step produces outputs which are a factor of sqrt(N)
73  * larger than the true IDCT outputs.  The final outputs are therefore
74  * a factor of N larger than desired; since N=8 this can be cured by
75  * a simple right shift at the end of the algorithm.  The advantage of
76  * this arrangement is that we save two multiplications per 1-D IDCT,
77  * because the y0 and y4 inputs need not be divided by sqrt(N).
78  *
79  * We have to do addition and subtraction of the integer inputs, which
80  * is no problem, and multiplication by fractional constants, which is
81  * a problem to do in integer arithmetic.  We multiply all the constants
82  * by CONST_SCALE and convert them to integer constants (thus retaining
83  * CONST_BITS bits of precision in the constants).  After doing a
84  * multiplication we have to divide the product by CONST_SCALE, with proper
85  * rounding, to produce the correct output.  This division can be done
86  * cheaply as a right shift of CONST_BITS bits.  We postpone shifting
87  * as long as possible so that partial sums can be added together with
88  * full fractional precision.
89  *
90  * The outputs of the first pass are scaled up by PASS1_BITS bits so that
91  * they are represented to better-than-integral precision.  These outputs
92  * require BITS_IN_JSAMPLE + PASS1_BITS + 3 bits; this fits in a 16-bit word
93  * with the recommended scaling.  (To scale up 12-bit sample data further, an
94  * intermediate INT32 array would be needed.)
95  *
96  * To avoid overflow of the 32-bit intermediate results in pass 2, we must
97  * have BITS_IN_JSAMPLE + CONST_BITS + PASS1_BITS <= 26.  Error analysis
98  * shows that the values given below are the most effective.
99  */
100 
101 #if BITS_IN_JSAMPLE == 8
102 #define CONST_BITS  13
103 #define PASS1_BITS  2
104 #else
105 #define CONST_BITS  13
106 #define PASS1_BITS  1           /* lose a little precision to avoid overflow */
107 #endif
108 
109 /* Some C compilers fail to reduce "FIX(constant)" at compile time, thus
110  * causing a lot of useless floating-point operations at run time.
111  * To get around this we use the following pre-calculated constants.
112  * If you change CONST_BITS you may want to add appropriate values.
113  * (With a reasonable C compiler, you can just rely on the FIX() macro...)
114  */
115 
116 #if CONST_BITS == 13
117 #define FIX_0_298631336  ((INT32)  2446)        /* FIX(0.298631336) */
118 #define FIX_0_390180644  ((INT32)  3196)        /* FIX(0.390180644) */
119 #define FIX_0_541196100  ((INT32)  4433)        /* FIX(0.541196100) */
120 #define FIX_0_765366865  ((INT32)  6270)        /* FIX(0.765366865) */
121 #define FIX_0_899976223  ((INT32)  7373)        /* FIX(0.899976223) */
122 #define FIX_1_175875602  ((INT32)  9633)        /* FIX(1.175875602) */
123 #define FIX_1_501321110  ((INT32)  12299)       /* FIX(1.501321110) */
124 #define FIX_1_847759065  ((INT32)  15137)       /* FIX(1.847759065) */
125 #define FIX_1_961570560  ((INT32)  16069)       /* FIX(1.961570560) */
126 #define FIX_2_053119869  ((INT32)  16819)       /* FIX(2.053119869) */
127 #define FIX_2_562915447  ((INT32)  20995)       /* FIX(2.562915447) */
128 #define FIX_3_072711026  ((INT32)  25172)       /* FIX(3.072711026) */
129 #else
130 #define FIX_0_298631336  FIX(0.298631336)
131 #define FIX_0_390180644  FIX(0.390180644)
132 #define FIX_0_541196100  FIX(0.541196100)
133 #define FIX_0_765366865  FIX(0.765366865)
134 #define FIX_0_899976223  FIX(0.899976223)
135 #define FIX_1_175875602  FIX(1.175875602)
136 #define FIX_1_501321110  FIX(1.501321110)
137 #define FIX_1_847759065  FIX(1.847759065)
138 #define FIX_1_961570560  FIX(1.961570560)
139 #define FIX_2_053119869  FIX(2.053119869)
140 #define FIX_2_562915447  FIX(2.562915447)
141 #define FIX_3_072711026  FIX(3.072711026)
142 #endif
143 
144 
145 /* Multiply an INT32 variable by an INT32 constant to yield an INT32 result.
146  * For 8-bit samples with the recommended scaling, all the variable
147  * and constant values involved are no more than 16 bits wide, so a
148  * 16x16->32 bit multiply can be used instead of a full 32x32 multiply.
149  * For 12-bit samples, a full 32-bit multiplication will be needed.
150  */
151 
152 #if BITS_IN_JSAMPLE == 8
153 #define MULTIPLY(var,const)  MULTIPLY16C16(var,const)
154 #else
155 #define MULTIPLY(var,const)  ((var) * (const))
156 #endif
157 
158 
159 /* Dequantize a coefficient by multiplying it by the multiplier-table
160  * entry; produce an int result.  In this module, both inputs and result
161  * are 16 bits or less, so either int or short multiply will work.
162  */
163 
164 #define DEQUANTIZE(coef,quantval)  (((ISLOW_MULT_TYPE) (coef)) * (quantval))
165 
166 
167 /*
168  * Perform dequantization and inverse DCT on one block of coefficients.
169  */
170 
171 GLOBAL(void)
172 jpeg_idct_islow (j_decompress_ptr cinfo, jpeg_component_info * compptr,
173                  JCOEFPTR coef_block,
174                  JSAMPARRAY output_buf, JDIMENSION output_col)
175 {
176   INT32 tmp0, tmp1, tmp2, tmp3;
177   INT32 tmp10, tmp11, tmp12, tmp13;
178   INT32 z1, z2, z3, z4, z5;
179   JCOEFPTR inptr;
180   ISLOW_MULT_TYPE * quantptr;
181   int * wsptr;
182   JSAMPROW outptr;
183   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
184   int ctr;
185   int workspace[DCTSIZE2];      /* buffers data between passes */
186   SHIFT_TEMPS
187 
188   /* Pass 1: process columns from input, store into work array. */
189   /* Note results are scaled up by sqrt(8) compared to a true IDCT; */
190   /* furthermore, we scale the results by 2**PASS1_BITS. */
191 
192   inptr = coef_block;
193   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
194   wsptr = workspace;
195   for (ctr = DCTSIZE; ctr > 0; ctr--) {
196     /* Due to quantization, we will usually find that many of the input
197      * coefficients are zero, especially the AC terms.  We can exploit this
198      * by short-circuiting the IDCT calculation for any column in which all
199      * the AC terms are zero.  In that case each output is equal to the
200      * DC coefficient (with scale factor as needed).
201      * With typical images and quantization tables, half or more of the
202      * column DCT calculations can be simplified this way.
203      */
204 
205     if (inptr[DCTSIZE*1] == 0 && inptr[DCTSIZE*2] == 0 &&
206         inptr[DCTSIZE*3] == 0 && inptr[DCTSIZE*4] == 0 &&
207         inptr[DCTSIZE*5] == 0 && inptr[DCTSIZE*6] == 0 &&
208         inptr[DCTSIZE*7] == 0) {
209       /* AC terms all zero */
210       int dcval = LEFT_SHIFT(DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]),
211                              PASS1_BITS);
212 
213       wsptr[DCTSIZE*0] = dcval;
214       wsptr[DCTSIZE*1] = dcval;
215       wsptr[DCTSIZE*2] = dcval;
216       wsptr[DCTSIZE*3] = dcval;
217       wsptr[DCTSIZE*4] = dcval;
218       wsptr[DCTSIZE*5] = dcval;
219       wsptr[DCTSIZE*6] = dcval;
220       wsptr[DCTSIZE*7] = dcval;
221 
222       inptr++;                  /* advance pointers to next column */
223       quantptr++;
224       wsptr++;
225       continue;
226     }
227 
228     /* Even part: reverse the even part of the forward DCT. */
229     /* The rotator is sqrt(2)*c(-6). */
230 
231     z2 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
232     z3 = DEQUANTIZE(inptr[DCTSIZE*6], quantptr[DCTSIZE*6]);
233 
234     z1 = MULTIPLY(z2 + z3, FIX_0_541196100);
235     tmp2 = z1 + MULTIPLY(z3, - FIX_1_847759065);
236     tmp3 = z1 + MULTIPLY(z2, FIX_0_765366865);
237 
238     z2 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
239     z3 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
240 
241     tmp0 = LEFT_SHIFT(z2 + z3, CONST_BITS);
242     tmp1 = LEFT_SHIFT(z2 - z3, CONST_BITS);
243 
244     tmp10 = tmp0 + tmp3;
245     tmp13 = tmp0 - tmp3;
246     tmp11 = tmp1 + tmp2;
247     tmp12 = tmp1 - tmp2;
248 
249     /* Odd part per figure 8; the matrix is unitary and hence its
250      * transpose is its inverse.  i0..i3 are y7,y5,y3,y1 respectively.
251      */
252 
253     tmp0 = DEQUANTIZE(inptr[DCTSIZE*7], quantptr[DCTSIZE*7]);
254     tmp1 = DEQUANTIZE(inptr[DCTSIZE*5], quantptr[DCTSIZE*5]);
255     tmp2 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
256     tmp3 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
257 
258     z1 = tmp0 + tmp3;
259     z2 = tmp1 + tmp2;
260     z3 = tmp0 + tmp2;
261     z4 = tmp1 + tmp3;
262     z5 = MULTIPLY(z3 + z4, FIX_1_175875602); /* sqrt(2) * c3 */
263 
264     tmp0 = MULTIPLY(tmp0, FIX_0_298631336); /* sqrt(2) * (-c1+c3+c5-c7) */
265     tmp1 = MULTIPLY(tmp1, FIX_2_053119869); /* sqrt(2) * ( c1+c3-c5+c7) */
266     tmp2 = MULTIPLY(tmp2, FIX_3_072711026); /* sqrt(2) * ( c1+c3+c5-c7) */
267     tmp3 = MULTIPLY(tmp3, FIX_1_501321110); /* sqrt(2) * ( c1+c3-c5-c7) */
268     z1 = MULTIPLY(z1, - FIX_0_899976223); /* sqrt(2) * (c7-c3) */
269     z2 = MULTIPLY(z2, - FIX_2_562915447); /* sqrt(2) * (-c1-c3) */
270     z3 = MULTIPLY(z3, - FIX_1_961570560); /* sqrt(2) * (-c3-c5) */
271     z4 = MULTIPLY(z4, - FIX_0_390180644); /* sqrt(2) * (c5-c3) */
272 
273     z3 += z5;
274     z4 += z5;
275 
276     tmp0 += z1 + z3;
277     tmp1 += z2 + z4;
278     tmp2 += z2 + z3;
279     tmp3 += z1 + z4;
280 
281     /* Final output stage: inputs are tmp10..tmp13, tmp0..tmp3 */
282 
283     wsptr[DCTSIZE*0] = (int) DESCALE(tmp10 + tmp3, CONST_BITS-PASS1_BITS);
284     wsptr[DCTSIZE*7] = (int) DESCALE(tmp10 - tmp3, CONST_BITS-PASS1_BITS);
285     wsptr[DCTSIZE*1] = (int) DESCALE(tmp11 + tmp2, CONST_BITS-PASS1_BITS);
286     wsptr[DCTSIZE*6] = (int) DESCALE(tmp11 - tmp2, CONST_BITS-PASS1_BITS);
287     wsptr[DCTSIZE*2] = (int) DESCALE(tmp12 + tmp1, CONST_BITS-PASS1_BITS);
288     wsptr[DCTSIZE*5] = (int) DESCALE(tmp12 - tmp1, CONST_BITS-PASS1_BITS);
289     wsptr[DCTSIZE*3] = (int) DESCALE(tmp13 + tmp0, CONST_BITS-PASS1_BITS);
290     wsptr[DCTSIZE*4] = (int) DESCALE(tmp13 - tmp0, CONST_BITS-PASS1_BITS);
291 
292     inptr++;                    /* advance pointers to next column */
293     quantptr++;
294     wsptr++;
295   }
296 
297   /* Pass 2: process rows from work array, store into output array. */
298   /* Note that we must descale the results by a factor of 8 == 2**3, */
299   /* and also undo the PASS1_BITS scaling. */
300 
301   wsptr = workspace;
302   for (ctr = 0; ctr < DCTSIZE; ctr++) {
303     outptr = output_buf[ctr] + output_col;
304     /* Rows of zeroes can be exploited in the same way as we did with columns.
305      * However, the column calculation has created many nonzero AC terms, so
306      * the simplification applies less often (typically 5% to 10% of the time).
307      * On machines with very fast multiplication, it's possible that the
308      * test takes more time than it's worth.  In that case this section
309      * may be commented out.
310      */
311 
312 #ifndef NO_ZERO_ROW_TEST
313     if (wsptr[1] == 0 && wsptr[2] == 0 && wsptr[3] == 0 && wsptr[4] == 0 &&
314         wsptr[5] == 0 && wsptr[6] == 0 && wsptr[7] == 0) {
315       /* AC terms all zero */
316       JSAMPLE dcval = range_limit[(int) DESCALE((INT32) wsptr[0], PASS1_BITS+3)
317                                   & RANGE_MASK];
318 
319       outptr[0] = dcval;
320       outptr[1] = dcval;
321       outptr[2] = dcval;
322       outptr[3] = dcval;
323       outptr[4] = dcval;
324       outptr[5] = dcval;
325       outptr[6] = dcval;
326       outptr[7] = dcval;
327 
328       wsptr += DCTSIZE;         /* advance pointer to next row */
329       continue;
330     }
331 #endif
332 
333     /* Even part: reverse the even part of the forward DCT. */
334     /* The rotator is sqrt(2)*c(-6). */
335 
336     z2 = (INT32) wsptr[2];
337     z3 = (INT32) wsptr[6];
338 
339     z1 = MULTIPLY(z2 + z3, FIX_0_541196100);
340     tmp2 = z1 + MULTIPLY(z3, - FIX_1_847759065);
341     tmp3 = z1 + MULTIPLY(z2, FIX_0_765366865);
342 
343     tmp0 = LEFT_SHIFT((INT32) wsptr[0] + (INT32) wsptr[4], CONST_BITS);
344     tmp1 = LEFT_SHIFT((INT32) wsptr[0] - (INT32) wsptr[4], CONST_BITS);
345 
346     tmp10 = tmp0 + tmp3;
347     tmp13 = tmp0 - tmp3;
348     tmp11 = tmp1 + tmp2;
349     tmp12 = tmp1 - tmp2;
350 
351     /* Odd part per figure 8; the matrix is unitary and hence its
352      * transpose is its inverse.  i0..i3 are y7,y5,y3,y1 respectively.
353      */
354 
355     tmp0 = (INT32) wsptr[7];
356     tmp1 = (INT32) wsptr[5];
357     tmp2 = (INT32) wsptr[3];
358     tmp3 = (INT32) wsptr[1];
359 
360     z1 = tmp0 + tmp3;
361     z2 = tmp1 + tmp2;
362     z3 = tmp0 + tmp2;
363     z4 = tmp1 + tmp3;
364     z5 = MULTIPLY(z3 + z4, FIX_1_175875602); /* sqrt(2) * c3 */
365 
366     tmp0 = MULTIPLY(tmp0, FIX_0_298631336); /* sqrt(2) * (-c1+c3+c5-c7) */
367     tmp1 = MULTIPLY(tmp1, FIX_2_053119869); /* sqrt(2) * ( c1+c3-c5+c7) */
368     tmp2 = MULTIPLY(tmp2, FIX_3_072711026); /* sqrt(2) * ( c1+c3+c5-c7) */
369     tmp3 = MULTIPLY(tmp3, FIX_1_501321110); /* sqrt(2) * ( c1+c3-c5-c7) */
370     z1 = MULTIPLY(z1, - FIX_0_899976223); /* sqrt(2) * (c7-c3) */
371     z2 = MULTIPLY(z2, - FIX_2_562915447); /* sqrt(2) * (-c1-c3) */
372     z3 = MULTIPLY(z3, - FIX_1_961570560); /* sqrt(2) * (-c3-c5) */
373     z4 = MULTIPLY(z4, - FIX_0_390180644); /* sqrt(2) * (c5-c3) */
374 
375     z3 += z5;
376     z4 += z5;
377 
378     tmp0 += z1 + z3;
379     tmp1 += z2 + z4;
380     tmp2 += z2 + z3;
381     tmp3 += z1 + z4;
382 
383     /* Final output stage: inputs are tmp10..tmp13, tmp0..tmp3 */
384 
385     outptr[0] = range_limit[(int) DESCALE(tmp10 + tmp3,
386                                           CONST_BITS+PASS1_BITS+3)
387                             & RANGE_MASK];
388     outptr[7] = range_limit[(int) DESCALE(tmp10 - tmp3,
389                                           CONST_BITS+PASS1_BITS+3)
390                             & RANGE_MASK];
391     outptr[1] = range_limit[(int) DESCALE(tmp11 + tmp2,
392                                           CONST_BITS+PASS1_BITS+3)
393                             & RANGE_MASK];
394     outptr[6] = range_limit[(int) DESCALE(tmp11 - tmp2,
395                                           CONST_BITS+PASS1_BITS+3)
396                             & RANGE_MASK];
397     outptr[2] = range_limit[(int) DESCALE(tmp12 + tmp1,
398                                           CONST_BITS+PASS1_BITS+3)
399                             & RANGE_MASK];
400     outptr[5] = range_limit[(int) DESCALE(tmp12 - tmp1,
401                                           CONST_BITS+PASS1_BITS+3)
402                             & RANGE_MASK];
403     outptr[3] = range_limit[(int) DESCALE(tmp13 + tmp0,
404                                           CONST_BITS+PASS1_BITS+3)
405                             & RANGE_MASK];
406     outptr[4] = range_limit[(int) DESCALE(tmp13 - tmp0,
407                                           CONST_BITS+PASS1_BITS+3)
408                             & RANGE_MASK];
409 
410     wsptr += DCTSIZE;           /* advance pointer to next row */
411   }
412 }
413 
414 #ifdef IDCT_SCALING_SUPPORTED
415 
416 
417 /*
418  * Perform dequantization and inverse DCT on one block of coefficients,
419  * producing a 7x7 output block.
420  *
421  * Optimized algorithm with 12 multiplications in the 1-D kernel.
422  * cK represents sqrt(2) * cos(K*pi/14).
423  */
424 
425 GLOBAL(void)
jpeg_idct_7x7(j_decompress_ptr cinfo,jpeg_component_info * compptr,JCOEFPTR coef_block,JSAMPARRAY output_buf,JDIMENSION output_col)426 jpeg_idct_7x7 (j_decompress_ptr cinfo, jpeg_component_info * compptr,
427                JCOEFPTR coef_block,
428                JSAMPARRAY output_buf, JDIMENSION output_col)
429 {
430   INT32 tmp0, tmp1, tmp2, tmp10, tmp11, tmp12, tmp13;
431   INT32 z1, z2, z3;
432   JCOEFPTR inptr;
433   ISLOW_MULT_TYPE * quantptr;
434   int * wsptr;
435   JSAMPROW outptr;
436   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
437   int ctr;
438   int workspace[7*7];   /* buffers data between passes */
439   SHIFT_TEMPS
440 
441   /* Pass 1: process columns from input, store into work array. */
442 
443   inptr = coef_block;
444   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
445   wsptr = workspace;
446   for (ctr = 0; ctr < 7; ctr++, inptr++, quantptr++, wsptr++) {
447     /* Even part */
448 
449     tmp13 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
450     tmp13 = LEFT_SHIFT(tmp13, CONST_BITS);
451     /* Add fudge factor here for final descale. */
452     tmp13 += ONE << (CONST_BITS-PASS1_BITS-1);
453 
454     z1 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
455     z2 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
456     z3 = DEQUANTIZE(inptr[DCTSIZE*6], quantptr[DCTSIZE*6]);
457 
458     tmp10 = MULTIPLY(z2 - z3, FIX(0.881747734));     /* c4 */
459     tmp12 = MULTIPLY(z1 - z2, FIX(0.314692123));     /* c6 */
460     tmp11 = tmp10 + tmp12 + tmp13 - MULTIPLY(z2, FIX(1.841218003)); /* c2+c4-c6 */
461     tmp0 = z1 + z3;
462     z2 -= tmp0;
463     tmp0 = MULTIPLY(tmp0, FIX(1.274162392)) + tmp13; /* c2 */
464     tmp10 += tmp0 - MULTIPLY(z3, FIX(0.077722536));  /* c2-c4-c6 */
465     tmp12 += tmp0 - MULTIPLY(z1, FIX(2.470602249));  /* c2+c4+c6 */
466     tmp13 += MULTIPLY(z2, FIX(1.414213562));         /* c0 */
467 
468     /* Odd part */
469 
470     z1 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
471     z2 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
472     z3 = DEQUANTIZE(inptr[DCTSIZE*5], quantptr[DCTSIZE*5]);
473 
474     tmp1 = MULTIPLY(z1 + z2, FIX(0.935414347));      /* (c3+c1-c5)/2 */
475     tmp2 = MULTIPLY(z1 - z2, FIX(0.170262339));      /* (c3+c5-c1)/2 */
476     tmp0 = tmp1 - tmp2;
477     tmp1 += tmp2;
478     tmp2 = MULTIPLY(z2 + z3, - FIX(1.378756276));    /* -c1 */
479     tmp1 += tmp2;
480     z2 = MULTIPLY(z1 + z3, FIX(0.613604268));        /* c5 */
481     tmp0 += z2;
482     tmp2 += z2 + MULTIPLY(z3, FIX(1.870828693));     /* c3+c1-c5 */
483 
484     /* Final output stage */
485 
486     wsptr[7*0] = (int) RIGHT_SHIFT(tmp10 + tmp0, CONST_BITS-PASS1_BITS);
487     wsptr[7*6] = (int) RIGHT_SHIFT(tmp10 - tmp0, CONST_BITS-PASS1_BITS);
488     wsptr[7*1] = (int) RIGHT_SHIFT(tmp11 + tmp1, CONST_BITS-PASS1_BITS);
489     wsptr[7*5] = (int) RIGHT_SHIFT(tmp11 - tmp1, CONST_BITS-PASS1_BITS);
490     wsptr[7*2] = (int) RIGHT_SHIFT(tmp12 + tmp2, CONST_BITS-PASS1_BITS);
491     wsptr[7*4] = (int) RIGHT_SHIFT(tmp12 - tmp2, CONST_BITS-PASS1_BITS);
492     wsptr[7*3] = (int) RIGHT_SHIFT(tmp13, CONST_BITS-PASS1_BITS);
493   }
494 
495   /* Pass 2: process 7 rows from work array, store into output array. */
496 
497   wsptr = workspace;
498   for (ctr = 0; ctr < 7; ctr++) {
499     outptr = output_buf[ctr] + output_col;
500 
501     /* Even part */
502 
503     /* Add fudge factor here for final descale. */
504     tmp13 = (INT32) wsptr[0] + (ONE << (PASS1_BITS+2));
505     tmp13 = LEFT_SHIFT(tmp13, CONST_BITS);
506 
507     z1 = (INT32) wsptr[2];
508     z2 = (INT32) wsptr[4];
509     z3 = (INT32) wsptr[6];
510 
511     tmp10 = MULTIPLY(z2 - z3, FIX(0.881747734));     /* c4 */
512     tmp12 = MULTIPLY(z1 - z2, FIX(0.314692123));     /* c6 */
513     tmp11 = tmp10 + tmp12 + tmp13 - MULTIPLY(z2, FIX(1.841218003)); /* c2+c4-c6 */
514     tmp0 = z1 + z3;
515     z2 -= tmp0;
516     tmp0 = MULTIPLY(tmp0, FIX(1.274162392)) + tmp13; /* c2 */
517     tmp10 += tmp0 - MULTIPLY(z3, FIX(0.077722536));  /* c2-c4-c6 */
518     tmp12 += tmp0 - MULTIPLY(z1, FIX(2.470602249));  /* c2+c4+c6 */
519     tmp13 += MULTIPLY(z2, FIX(1.414213562));         /* c0 */
520 
521     /* Odd part */
522 
523     z1 = (INT32) wsptr[1];
524     z2 = (INT32) wsptr[3];
525     z3 = (INT32) wsptr[5];
526 
527     tmp1 = MULTIPLY(z1 + z2, FIX(0.935414347));      /* (c3+c1-c5)/2 */
528     tmp2 = MULTIPLY(z1 - z2, FIX(0.170262339));      /* (c3+c5-c1)/2 */
529     tmp0 = tmp1 - tmp2;
530     tmp1 += tmp2;
531     tmp2 = MULTIPLY(z2 + z3, - FIX(1.378756276));    /* -c1 */
532     tmp1 += tmp2;
533     z2 = MULTIPLY(z1 + z3, FIX(0.613604268));        /* c5 */
534     tmp0 += z2;
535     tmp2 += z2 + MULTIPLY(z3, FIX(1.870828693));     /* c3+c1-c5 */
536 
537     /* Final output stage */
538 
539     outptr[0] = range_limit[(int) RIGHT_SHIFT(tmp10 + tmp0,
540                                               CONST_BITS+PASS1_BITS+3)
541                             & RANGE_MASK];
542     outptr[6] = range_limit[(int) RIGHT_SHIFT(tmp10 - tmp0,
543                                               CONST_BITS+PASS1_BITS+3)
544                             & RANGE_MASK];
545     outptr[1] = range_limit[(int) RIGHT_SHIFT(tmp11 + tmp1,
546                                               CONST_BITS+PASS1_BITS+3)
547                             & RANGE_MASK];
548     outptr[5] = range_limit[(int) RIGHT_SHIFT(tmp11 - tmp1,
549                                               CONST_BITS+PASS1_BITS+3)
550                             & RANGE_MASK];
551     outptr[2] = range_limit[(int) RIGHT_SHIFT(tmp12 + tmp2,
552                                               CONST_BITS+PASS1_BITS+3)
553                             & RANGE_MASK];
554     outptr[4] = range_limit[(int) RIGHT_SHIFT(tmp12 - tmp2,
555                                               CONST_BITS+PASS1_BITS+3)
556                             & RANGE_MASK];
557     outptr[3] = range_limit[(int) RIGHT_SHIFT(tmp13,
558                                               CONST_BITS+PASS1_BITS+3)
559                             & RANGE_MASK];
560 
561     wsptr += 7;         /* advance pointer to next row */
562   }
563 }
564 
565 
566 /*
567  * Perform dequantization and inverse DCT on one block of coefficients,
568  * producing a reduced-size 6x6 output block.
569  *
570  * Optimized algorithm with 3 multiplications in the 1-D kernel.
571  * cK represents sqrt(2) * cos(K*pi/12).
572  */
573 
574 GLOBAL(void)
jpeg_idct_6x6(j_decompress_ptr cinfo,jpeg_component_info * compptr,JCOEFPTR coef_block,JSAMPARRAY output_buf,JDIMENSION output_col)575 jpeg_idct_6x6 (j_decompress_ptr cinfo, jpeg_component_info * compptr,
576                JCOEFPTR coef_block,
577                JSAMPARRAY output_buf, JDIMENSION output_col)
578 {
579   INT32 tmp0, tmp1, tmp2, tmp10, tmp11, tmp12;
580   INT32 z1, z2, z3;
581   JCOEFPTR inptr;
582   ISLOW_MULT_TYPE * quantptr;
583   int * wsptr;
584   JSAMPROW outptr;
585   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
586   int ctr;
587   int workspace[6*6];   /* buffers data between passes */
588   SHIFT_TEMPS
589 
590   /* Pass 1: process columns from input, store into work array. */
591 
592   inptr = coef_block;
593   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
594   wsptr = workspace;
595   for (ctr = 0; ctr < 6; ctr++, inptr++, quantptr++, wsptr++) {
596     /* Even part */
597 
598     tmp0 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
599     tmp0 = LEFT_SHIFT(tmp0, CONST_BITS);
600     /* Add fudge factor here for final descale. */
601     tmp0 += ONE << (CONST_BITS-PASS1_BITS-1);
602     tmp2 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
603     tmp10 = MULTIPLY(tmp2, FIX(0.707106781));   /* c4 */
604     tmp1 = tmp0 + tmp10;
605     tmp11 = RIGHT_SHIFT(tmp0 - tmp10 - tmp10, CONST_BITS-PASS1_BITS);
606     tmp10 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
607     tmp0 = MULTIPLY(tmp10, FIX(1.224744871));   /* c2 */
608     tmp10 = tmp1 + tmp0;
609     tmp12 = tmp1 - tmp0;
610 
611     /* Odd part */
612 
613     z1 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
614     z2 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
615     z3 = DEQUANTIZE(inptr[DCTSIZE*5], quantptr[DCTSIZE*5]);
616     tmp1 = MULTIPLY(z1 + z3, FIX(0.366025404)); /* c5 */
617     tmp0 = tmp1 + LEFT_SHIFT(z1 + z2, CONST_BITS);
618     tmp2 = tmp1 + LEFT_SHIFT(z3 - z2, CONST_BITS);
619     tmp1 = LEFT_SHIFT(z1 - z2 - z3, PASS1_BITS);
620 
621     /* Final output stage */
622 
623     wsptr[6*0] = (int) RIGHT_SHIFT(tmp10 + tmp0, CONST_BITS-PASS1_BITS);
624     wsptr[6*5] = (int) RIGHT_SHIFT(tmp10 - tmp0, CONST_BITS-PASS1_BITS);
625     wsptr[6*1] = (int) (tmp11 + tmp1);
626     wsptr[6*4] = (int) (tmp11 - tmp1);
627     wsptr[6*2] = (int) RIGHT_SHIFT(tmp12 + tmp2, CONST_BITS-PASS1_BITS);
628     wsptr[6*3] = (int) RIGHT_SHIFT(tmp12 - tmp2, CONST_BITS-PASS1_BITS);
629   }
630 
631   /* Pass 2: process 6 rows from work array, store into output array. */
632 
633   wsptr = workspace;
634   for (ctr = 0; ctr < 6; ctr++) {
635     outptr = output_buf[ctr] + output_col;
636 
637     /* Even part */
638 
639     /* Add fudge factor here for final descale. */
640     tmp0 = (INT32) wsptr[0] + (ONE << (PASS1_BITS+2));
641     tmp0 = LEFT_SHIFT(tmp0, CONST_BITS);
642     tmp2 = (INT32) wsptr[4];
643     tmp10 = MULTIPLY(tmp2, FIX(0.707106781));   /* c4 */
644     tmp1 = tmp0 + tmp10;
645     tmp11 = tmp0 - tmp10 - tmp10;
646     tmp10 = (INT32) wsptr[2];
647     tmp0 = MULTIPLY(tmp10, FIX(1.224744871));   /* c2 */
648     tmp10 = tmp1 + tmp0;
649     tmp12 = tmp1 - tmp0;
650 
651     /* Odd part */
652 
653     z1 = (INT32) wsptr[1];
654     z2 = (INT32) wsptr[3];
655     z3 = (INT32) wsptr[5];
656     tmp1 = MULTIPLY(z1 + z3, FIX(0.366025404)); /* c5 */
657     tmp0 = tmp1 + LEFT_SHIFT(z1 + z2, CONST_BITS);
658     tmp2 = tmp1 + LEFT_SHIFT(z3 - z2, CONST_BITS);
659     tmp1 = LEFT_SHIFT(z1 - z2 - z3, CONST_BITS);
660 
661     /* Final output stage */
662 
663     outptr[0] = range_limit[(int) RIGHT_SHIFT(tmp10 + tmp0,
664                                               CONST_BITS+PASS1_BITS+3)
665                             & RANGE_MASK];
666     outptr[5] = range_limit[(int) RIGHT_SHIFT(tmp10 - tmp0,
667                                               CONST_BITS+PASS1_BITS+3)
668                             & RANGE_MASK];
669     outptr[1] = range_limit[(int) RIGHT_SHIFT(tmp11 + tmp1,
670                                               CONST_BITS+PASS1_BITS+3)
671                             & RANGE_MASK];
672     outptr[4] = range_limit[(int) RIGHT_SHIFT(tmp11 - tmp1,
673                                               CONST_BITS+PASS1_BITS+3)
674                             & RANGE_MASK];
675     outptr[2] = range_limit[(int) RIGHT_SHIFT(tmp12 + tmp2,
676                                               CONST_BITS+PASS1_BITS+3)
677                             & RANGE_MASK];
678     outptr[3] = range_limit[(int) RIGHT_SHIFT(tmp12 - tmp2,
679                                               CONST_BITS+PASS1_BITS+3)
680                             & RANGE_MASK];
681 
682     wsptr += 6;         /* advance pointer to next row */
683   }
684 }
685 
686 
687 /*
688  * Perform dequantization and inverse DCT on one block of coefficients,
689  * producing a reduced-size 5x5 output block.
690  *
691  * Optimized algorithm with 5 multiplications in the 1-D kernel.
692  * cK represents sqrt(2) * cos(K*pi/10).
693  */
694 
695 GLOBAL(void)
jpeg_idct_5x5(j_decompress_ptr cinfo,jpeg_component_info * compptr,JCOEFPTR coef_block,JSAMPARRAY output_buf,JDIMENSION output_col)696 jpeg_idct_5x5 (j_decompress_ptr cinfo, jpeg_component_info * compptr,
697                JCOEFPTR coef_block,
698                JSAMPARRAY output_buf, JDIMENSION output_col)
699 {
700   INT32 tmp0, tmp1, tmp10, tmp11, tmp12;
701   INT32 z1, z2, z3;
702   JCOEFPTR inptr;
703   ISLOW_MULT_TYPE * quantptr;
704   int * wsptr;
705   JSAMPROW outptr;
706   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
707   int ctr;
708   int workspace[5*5];   /* buffers data between passes */
709   SHIFT_TEMPS
710 
711   /* Pass 1: process columns from input, store into work array. */
712 
713   inptr = coef_block;
714   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
715   wsptr = workspace;
716   for (ctr = 0; ctr < 5; ctr++, inptr++, quantptr++, wsptr++) {
717     /* Even part */
718 
719     tmp12 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
720     tmp12 = LEFT_SHIFT(tmp12, CONST_BITS);
721     /* Add fudge factor here for final descale. */
722     tmp12 += ONE << (CONST_BITS-PASS1_BITS-1);
723     tmp0 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
724     tmp1 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
725     z1 = MULTIPLY(tmp0 + tmp1, FIX(0.790569415)); /* (c2+c4)/2 */
726     z2 = MULTIPLY(tmp0 - tmp1, FIX(0.353553391)); /* (c2-c4)/2 */
727     z3 = tmp12 + z2;
728     tmp10 = z3 + z1;
729     tmp11 = z3 - z1;
730     tmp12 -= LEFT_SHIFT(z2, 2);
731 
732     /* Odd part */
733 
734     z2 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
735     z3 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
736 
737     z1 = MULTIPLY(z2 + z3, FIX(0.831253876));     /* c3 */
738     tmp0 = z1 + MULTIPLY(z2, FIX(0.513743148));   /* c1-c3 */
739     tmp1 = z1 - MULTIPLY(z3, FIX(2.176250899));   /* c1+c3 */
740 
741     /* Final output stage */
742 
743     wsptr[5*0] = (int) RIGHT_SHIFT(tmp10 + tmp0, CONST_BITS-PASS1_BITS);
744     wsptr[5*4] = (int) RIGHT_SHIFT(tmp10 - tmp0, CONST_BITS-PASS1_BITS);
745     wsptr[5*1] = (int) RIGHT_SHIFT(tmp11 + tmp1, CONST_BITS-PASS1_BITS);
746     wsptr[5*3] = (int) RIGHT_SHIFT(tmp11 - tmp1, CONST_BITS-PASS1_BITS);
747     wsptr[5*2] = (int) RIGHT_SHIFT(tmp12, CONST_BITS-PASS1_BITS);
748   }
749 
750   /* Pass 2: process 5 rows from work array, store into output array. */
751 
752   wsptr = workspace;
753   for (ctr = 0; ctr < 5; ctr++) {
754     outptr = output_buf[ctr] + output_col;
755 
756     /* Even part */
757 
758     /* Add fudge factor here for final descale. */
759     tmp12 = (INT32) wsptr[0] + (ONE << (PASS1_BITS+2));
760     tmp12 = LEFT_SHIFT(tmp12, CONST_BITS);
761     tmp0 = (INT32) wsptr[2];
762     tmp1 = (INT32) wsptr[4];
763     z1 = MULTIPLY(tmp0 + tmp1, FIX(0.790569415)); /* (c2+c4)/2 */
764     z2 = MULTIPLY(tmp0 - tmp1, FIX(0.353553391)); /* (c2-c4)/2 */
765     z3 = tmp12 + z2;
766     tmp10 = z3 + z1;
767     tmp11 = z3 - z1;
768     tmp12 -= LEFT_SHIFT(z2, 2);
769 
770     /* Odd part */
771 
772     z2 = (INT32) wsptr[1];
773     z3 = (INT32) wsptr[3];
774 
775     z1 = MULTIPLY(z2 + z3, FIX(0.831253876));     /* c3 */
776     tmp0 = z1 + MULTIPLY(z2, FIX(0.513743148));   /* c1-c3 */
777     tmp1 = z1 - MULTIPLY(z3, FIX(2.176250899));   /* c1+c3 */
778 
779     /* Final output stage */
780 
781     outptr[0] = range_limit[(int) RIGHT_SHIFT(tmp10 + tmp0,
782                                               CONST_BITS+PASS1_BITS+3)
783                             & RANGE_MASK];
784     outptr[4] = range_limit[(int) RIGHT_SHIFT(tmp10 - tmp0,
785                                               CONST_BITS+PASS1_BITS+3)
786                             & RANGE_MASK];
787     outptr[1] = range_limit[(int) RIGHT_SHIFT(tmp11 + tmp1,
788                                               CONST_BITS+PASS1_BITS+3)
789                             & RANGE_MASK];
790     outptr[3] = range_limit[(int) RIGHT_SHIFT(tmp11 - tmp1,
791                                               CONST_BITS+PASS1_BITS+3)
792                             & RANGE_MASK];
793     outptr[2] = range_limit[(int) RIGHT_SHIFT(tmp12,
794                                               CONST_BITS+PASS1_BITS+3)
795                             & RANGE_MASK];
796 
797     wsptr += 5;         /* advance pointer to next row */
798   }
799 }
800 
801 
802 /*
803  * Perform dequantization and inverse DCT on one block of coefficients,
804  * producing a reduced-size 3x3 output block.
805  *
806  * Optimized algorithm with 2 multiplications in the 1-D kernel.
807  * cK represents sqrt(2) * cos(K*pi/6).
808  */
809 
810 GLOBAL(void)
jpeg_idct_3x3(j_decompress_ptr cinfo,jpeg_component_info * compptr,JCOEFPTR coef_block,JSAMPARRAY output_buf,JDIMENSION output_col)811 jpeg_idct_3x3 (j_decompress_ptr cinfo, jpeg_component_info * compptr,
812                JCOEFPTR coef_block,
813                JSAMPARRAY output_buf, JDIMENSION output_col)
814 {
815   INT32 tmp0, tmp2, tmp10, tmp12;
816   JCOEFPTR inptr;
817   ISLOW_MULT_TYPE * quantptr;
818   int * wsptr;
819   JSAMPROW outptr;
820   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
821   int ctr;
822   int workspace[3*3];   /* buffers data between passes */
823   SHIFT_TEMPS
824 
825   /* Pass 1: process columns from input, store into work array. */
826 
827   inptr = coef_block;
828   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
829   wsptr = workspace;
830   for (ctr = 0; ctr < 3; ctr++, inptr++, quantptr++, wsptr++) {
831     /* Even part */
832 
833     tmp0 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
834     tmp0 = LEFT_SHIFT(tmp0, CONST_BITS);
835     /* Add fudge factor here for final descale. */
836     tmp0 += ONE << (CONST_BITS-PASS1_BITS-1);
837     tmp2 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
838     tmp12 = MULTIPLY(tmp2, FIX(0.707106781)); /* c2 */
839     tmp10 = tmp0 + tmp12;
840     tmp2 = tmp0 - tmp12 - tmp12;
841 
842     /* Odd part */
843 
844     tmp12 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
845     tmp0 = MULTIPLY(tmp12, FIX(1.224744871)); /* c1 */
846 
847     /* Final output stage */
848 
849     wsptr[3*0] = (int) RIGHT_SHIFT(tmp10 + tmp0, CONST_BITS-PASS1_BITS);
850     wsptr[3*2] = (int) RIGHT_SHIFT(tmp10 - tmp0, CONST_BITS-PASS1_BITS);
851     wsptr[3*1] = (int) RIGHT_SHIFT(tmp2, CONST_BITS-PASS1_BITS);
852   }
853 
854   /* Pass 2: process 3 rows from work array, store into output array. */
855 
856   wsptr = workspace;
857   for (ctr = 0; ctr < 3; ctr++) {
858     outptr = output_buf[ctr] + output_col;
859 
860     /* Even part */
861 
862     /* Add fudge factor here for final descale. */
863     tmp0 = (INT32) wsptr[0] + (ONE << (PASS1_BITS+2));
864     tmp0 = LEFT_SHIFT(tmp0, CONST_BITS);
865     tmp2 = (INT32) wsptr[2];
866     tmp12 = MULTIPLY(tmp2, FIX(0.707106781)); /* c2 */
867     tmp10 = tmp0 + tmp12;
868     tmp2 = tmp0 - tmp12 - tmp12;
869 
870     /* Odd part */
871 
872     tmp12 = (INT32) wsptr[1];
873     tmp0 = MULTIPLY(tmp12, FIX(1.224744871)); /* c1 */
874 
875     /* Final output stage */
876 
877     outptr[0] = range_limit[(int) RIGHT_SHIFT(tmp10 + tmp0,
878                                               CONST_BITS+PASS1_BITS+3)
879                             & RANGE_MASK];
880     outptr[2] = range_limit[(int) RIGHT_SHIFT(tmp10 - tmp0,
881                                               CONST_BITS+PASS1_BITS+3)
882                             & RANGE_MASK];
883     outptr[1] = range_limit[(int) RIGHT_SHIFT(tmp2,
884                                               CONST_BITS+PASS1_BITS+3)
885                             & RANGE_MASK];
886 
887     wsptr += 3;         /* advance pointer to next row */
888   }
889 }
890 
891 
892 /*
893  * Perform dequantization and inverse DCT on one block of coefficients,
894  * producing a 9x9 output block.
895  *
896  * Optimized algorithm with 10 multiplications in the 1-D kernel.
897  * cK represents sqrt(2) * cos(K*pi/18).
898  */
899 
900 GLOBAL(void)
jpeg_idct_9x9(j_decompress_ptr cinfo,jpeg_component_info * compptr,JCOEFPTR coef_block,JSAMPARRAY output_buf,JDIMENSION output_col)901 jpeg_idct_9x9 (j_decompress_ptr cinfo, jpeg_component_info * compptr,
902                JCOEFPTR coef_block,
903                JSAMPARRAY output_buf, JDIMENSION output_col)
904 {
905   INT32 tmp0, tmp1, tmp2, tmp3, tmp10, tmp11, tmp12, tmp13, tmp14;
906   INT32 z1, z2, z3, z4;
907   JCOEFPTR inptr;
908   ISLOW_MULT_TYPE * quantptr;
909   int * wsptr;
910   JSAMPROW outptr;
911   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
912   int ctr;
913   int workspace[8*9];   /* buffers data between passes */
914   SHIFT_TEMPS
915 
916   /* Pass 1: process columns from input, store into work array. */
917 
918   inptr = coef_block;
919   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
920   wsptr = workspace;
921   for (ctr = 0; ctr < 8; ctr++, inptr++, quantptr++, wsptr++) {
922     /* Even part */
923 
924     tmp0 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
925     tmp0 = LEFT_SHIFT(tmp0, CONST_BITS);
926     /* Add fudge factor here for final descale. */
927     tmp0 += ONE << (CONST_BITS-PASS1_BITS-1);
928 
929     z1 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
930     z2 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
931     z3 = DEQUANTIZE(inptr[DCTSIZE*6], quantptr[DCTSIZE*6]);
932 
933     tmp3 = MULTIPLY(z3, FIX(0.707106781));      /* c6 */
934     tmp1 = tmp0 + tmp3;
935     tmp2 = tmp0 - tmp3 - tmp3;
936 
937     tmp0 = MULTIPLY(z1 - z2, FIX(0.707106781)); /* c6 */
938     tmp11 = tmp2 + tmp0;
939     tmp14 = tmp2 - tmp0 - tmp0;
940 
941     tmp0 = MULTIPLY(z1 + z2, FIX(1.328926049)); /* c2 */
942     tmp2 = MULTIPLY(z1, FIX(1.083350441));      /* c4 */
943     tmp3 = MULTIPLY(z2, FIX(0.245575608));      /* c8 */
944 
945     tmp10 = tmp1 + tmp0 - tmp3;
946     tmp12 = tmp1 - tmp0 + tmp2;
947     tmp13 = tmp1 - tmp2 + tmp3;
948 
949     /* Odd part */
950 
951     z1 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
952     z2 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
953     z3 = DEQUANTIZE(inptr[DCTSIZE*5], quantptr[DCTSIZE*5]);
954     z4 = DEQUANTIZE(inptr[DCTSIZE*7], quantptr[DCTSIZE*7]);
955 
956     z2 = MULTIPLY(z2, - FIX(1.224744871));           /* -c3 */
957 
958     tmp2 = MULTIPLY(z1 + z3, FIX(0.909038955));      /* c5 */
959     tmp3 = MULTIPLY(z1 + z4, FIX(0.483689525));      /* c7 */
960     tmp0 = tmp2 + tmp3 - z2;
961     tmp1 = MULTIPLY(z3 - z4, FIX(1.392728481));      /* c1 */
962     tmp2 += z2 - tmp1;
963     tmp3 += z2 + tmp1;
964     tmp1 = MULTIPLY(z1 - z3 - z4, FIX(1.224744871)); /* c3 */
965 
966     /* Final output stage */
967 
968     wsptr[8*0] = (int) RIGHT_SHIFT(tmp10 + tmp0, CONST_BITS-PASS1_BITS);
969     wsptr[8*8] = (int) RIGHT_SHIFT(tmp10 - tmp0, CONST_BITS-PASS1_BITS);
970     wsptr[8*1] = (int) RIGHT_SHIFT(tmp11 + tmp1, CONST_BITS-PASS1_BITS);
971     wsptr[8*7] = (int) RIGHT_SHIFT(tmp11 - tmp1, CONST_BITS-PASS1_BITS);
972     wsptr[8*2] = (int) RIGHT_SHIFT(tmp12 + tmp2, CONST_BITS-PASS1_BITS);
973     wsptr[8*6] = (int) RIGHT_SHIFT(tmp12 - tmp2, CONST_BITS-PASS1_BITS);
974     wsptr[8*3] = (int) RIGHT_SHIFT(tmp13 + tmp3, CONST_BITS-PASS1_BITS);
975     wsptr[8*5] = (int) RIGHT_SHIFT(tmp13 - tmp3, CONST_BITS-PASS1_BITS);
976     wsptr[8*4] = (int) RIGHT_SHIFT(tmp14, CONST_BITS-PASS1_BITS);
977   }
978 
979   /* Pass 2: process 9 rows from work array, store into output array. */
980 
981   wsptr = workspace;
982   for (ctr = 0; ctr < 9; ctr++) {
983     outptr = output_buf[ctr] + output_col;
984 
985     /* Even part */
986 
987     /* Add fudge factor here for final descale. */
988     tmp0 = (INT32) wsptr[0] + (ONE << (PASS1_BITS+2));
989     tmp0 = LEFT_SHIFT(tmp0, CONST_BITS);
990 
991     z1 = (INT32) wsptr[2];
992     z2 = (INT32) wsptr[4];
993     z3 = (INT32) wsptr[6];
994 
995     tmp3 = MULTIPLY(z3, FIX(0.707106781));      /* c6 */
996     tmp1 = tmp0 + tmp3;
997     tmp2 = tmp0 - tmp3 - tmp3;
998 
999     tmp0 = MULTIPLY(z1 - z2, FIX(0.707106781)); /* c6 */
1000     tmp11 = tmp2 + tmp0;
1001     tmp14 = tmp2 - tmp0 - tmp0;
1002 
1003     tmp0 = MULTIPLY(z1 + z2, FIX(1.328926049)); /* c2 */
1004     tmp2 = MULTIPLY(z1, FIX(1.083350441));      /* c4 */
1005     tmp3 = MULTIPLY(z2, FIX(0.245575608));      /* c8 */
1006 
1007     tmp10 = tmp1 + tmp0 - tmp3;
1008     tmp12 = tmp1 - tmp0 + tmp2;
1009     tmp13 = tmp1 - tmp2 + tmp3;
1010 
1011     /* Odd part */
1012 
1013     z1 = (INT32) wsptr[1];
1014     z2 = (INT32) wsptr[3];
1015     z3 = (INT32) wsptr[5];
1016     z4 = (INT32) wsptr[7];
1017 
1018     z2 = MULTIPLY(z2, - FIX(1.224744871));           /* -c3 */
1019 
1020     tmp2 = MULTIPLY(z1 + z3, FIX(0.909038955));      /* c5 */
1021     tmp3 = MULTIPLY(z1 + z4, FIX(0.483689525));      /* c7 */
1022     tmp0 = tmp2 + tmp3 - z2;
1023     tmp1 = MULTIPLY(z3 - z4, FIX(1.392728481));      /* c1 */
1024     tmp2 += z2 - tmp1;
1025     tmp3 += z2 + tmp1;
1026     tmp1 = MULTIPLY(z1 - z3 - z4, FIX(1.224744871)); /* c3 */
1027 
1028     /* Final output stage */
1029 
1030     outptr[0] = range_limit[(int) RIGHT_SHIFT(tmp10 + tmp0,
1031                                               CONST_BITS+PASS1_BITS+3)
1032                             & RANGE_MASK];
1033     outptr[8] = range_limit[(int) RIGHT_SHIFT(tmp10 - tmp0,
1034                                               CONST_BITS+PASS1_BITS+3)
1035                             & RANGE_MASK];
1036     outptr[1] = range_limit[(int) RIGHT_SHIFT(tmp11 + tmp1,
1037                                               CONST_BITS+PASS1_BITS+3)
1038                             & RANGE_MASK];
1039     outptr[7] = range_limit[(int) RIGHT_SHIFT(tmp11 - tmp1,
1040                                               CONST_BITS+PASS1_BITS+3)
1041                             & RANGE_MASK];
1042     outptr[2] = range_limit[(int) RIGHT_SHIFT(tmp12 + tmp2,
1043                                               CONST_BITS+PASS1_BITS+3)
1044                             & RANGE_MASK];
1045     outptr[6] = range_limit[(int) RIGHT_SHIFT(tmp12 - tmp2,
1046                                               CONST_BITS+PASS1_BITS+3)
1047                             & RANGE_MASK];
1048     outptr[3] = range_limit[(int) RIGHT_SHIFT(tmp13 + tmp3,
1049                                               CONST_BITS+PASS1_BITS+3)
1050                             & RANGE_MASK];
1051     outptr[5] = range_limit[(int) RIGHT_SHIFT(tmp13 - tmp3,
1052                                               CONST_BITS+PASS1_BITS+3)
1053                             & RANGE_MASK];
1054     outptr[4] = range_limit[(int) RIGHT_SHIFT(tmp14,
1055                                               CONST_BITS+PASS1_BITS+3)
1056                             & RANGE_MASK];
1057 
1058     wsptr += 8;         /* advance pointer to next row */
1059   }
1060 }
1061 
1062 
1063 /*
1064  * Perform dequantization and inverse DCT on one block of coefficients,
1065  * producing a 10x10 output block.
1066  *
1067  * Optimized algorithm with 12 multiplications in the 1-D kernel.
1068  * cK represents sqrt(2) * cos(K*pi/20).
1069  */
1070 
1071 GLOBAL(void)
jpeg_idct_10x10(j_decompress_ptr cinfo,jpeg_component_info * compptr,JCOEFPTR coef_block,JSAMPARRAY output_buf,JDIMENSION output_col)1072 jpeg_idct_10x10 (j_decompress_ptr cinfo, jpeg_component_info * compptr,
1073                  JCOEFPTR coef_block,
1074                  JSAMPARRAY output_buf, JDIMENSION output_col)
1075 {
1076   INT32 tmp10, tmp11, tmp12, tmp13, tmp14;
1077   INT32 tmp20, tmp21, tmp22, tmp23, tmp24;
1078   INT32 z1, z2, z3, z4, z5;
1079   JCOEFPTR inptr;
1080   ISLOW_MULT_TYPE * quantptr;
1081   int * wsptr;
1082   JSAMPROW outptr;
1083   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
1084   int ctr;
1085   int workspace[8*10];  /* buffers data between passes */
1086   SHIFT_TEMPS
1087 
1088   /* Pass 1: process columns from input, store into work array. */
1089 
1090   inptr = coef_block;
1091   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
1092   wsptr = workspace;
1093   for (ctr = 0; ctr < 8; ctr++, inptr++, quantptr++, wsptr++) {
1094     /* Even part */
1095 
1096     z3 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
1097     z3 = LEFT_SHIFT(z3, CONST_BITS);
1098     /* Add fudge factor here for final descale. */
1099     z3 += ONE << (CONST_BITS-PASS1_BITS-1);
1100     z4 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
1101     z1 = MULTIPLY(z4, FIX(1.144122806));         /* c4 */
1102     z2 = MULTIPLY(z4, FIX(0.437016024));         /* c8 */
1103     tmp10 = z3 + z1;
1104     tmp11 = z3 - z2;
1105 
1106     tmp22 = RIGHT_SHIFT(z3 - LEFT_SHIFT(z1 - z2, 1),
1107                         CONST_BITS-PASS1_BITS);  /* c0 = (c4-c8)*2 */
1108 
1109     z2 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
1110     z3 = DEQUANTIZE(inptr[DCTSIZE*6], quantptr[DCTSIZE*6]);
1111 
1112     z1 = MULTIPLY(z2 + z3, FIX(0.831253876));    /* c6 */
1113     tmp12 = z1 + MULTIPLY(z2, FIX(0.513743148)); /* c2-c6 */
1114     tmp13 = z1 - MULTIPLY(z3, FIX(2.176250899)); /* c2+c6 */
1115 
1116     tmp20 = tmp10 + tmp12;
1117     tmp24 = tmp10 - tmp12;
1118     tmp21 = tmp11 + tmp13;
1119     tmp23 = tmp11 - tmp13;
1120 
1121     /* Odd part */
1122 
1123     z1 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
1124     z2 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
1125     z3 = DEQUANTIZE(inptr[DCTSIZE*5], quantptr[DCTSIZE*5]);
1126     z4 = DEQUANTIZE(inptr[DCTSIZE*7], quantptr[DCTSIZE*7]);
1127 
1128     tmp11 = z2 + z4;
1129     tmp13 = z2 - z4;
1130 
1131     tmp12 = MULTIPLY(tmp13, FIX(0.309016994));        /* (c3-c7)/2 */
1132     z5 = LEFT_SHIFT(z3, CONST_BITS);
1133 
1134     z2 = MULTIPLY(tmp11, FIX(0.951056516));           /* (c3+c7)/2 */
1135     z4 = z5 + tmp12;
1136 
1137     tmp10 = MULTIPLY(z1, FIX(1.396802247)) + z2 + z4; /* c1 */
1138     tmp14 = MULTIPLY(z1, FIX(0.221231742)) - z2 + z4; /* c9 */
1139 
1140     z2 = MULTIPLY(tmp11, FIX(0.587785252));           /* (c1-c9)/2 */
1141     z4 = z5 - tmp12 - LEFT_SHIFT(tmp13, CONST_BITS - 1);
1142 
1143     tmp12 = LEFT_SHIFT(z1 - tmp13 - z3, PASS1_BITS);
1144 
1145     tmp11 = MULTIPLY(z1, FIX(1.260073511)) - z2 - z4; /* c3 */
1146     tmp13 = MULTIPLY(z1, FIX(0.642039522)) - z2 + z4; /* c7 */
1147 
1148     /* Final output stage */
1149 
1150     wsptr[8*0] = (int) RIGHT_SHIFT(tmp20 + tmp10, CONST_BITS-PASS1_BITS);
1151     wsptr[8*9] = (int) RIGHT_SHIFT(tmp20 - tmp10, CONST_BITS-PASS1_BITS);
1152     wsptr[8*1] = (int) RIGHT_SHIFT(tmp21 + tmp11, CONST_BITS-PASS1_BITS);
1153     wsptr[8*8] = (int) RIGHT_SHIFT(tmp21 - tmp11, CONST_BITS-PASS1_BITS);
1154     wsptr[8*2] = (int) (tmp22 + tmp12);
1155     wsptr[8*7] = (int) (tmp22 - tmp12);
1156     wsptr[8*3] = (int) RIGHT_SHIFT(tmp23 + tmp13, CONST_BITS-PASS1_BITS);
1157     wsptr[8*6] = (int) RIGHT_SHIFT(tmp23 - tmp13, CONST_BITS-PASS1_BITS);
1158     wsptr[8*4] = (int) RIGHT_SHIFT(tmp24 + tmp14, CONST_BITS-PASS1_BITS);
1159     wsptr[8*5] = (int) RIGHT_SHIFT(tmp24 - tmp14, CONST_BITS-PASS1_BITS);
1160   }
1161 
1162   /* Pass 2: process 10 rows from work array, store into output array. */
1163 
1164   wsptr = workspace;
1165   for (ctr = 0; ctr < 10; ctr++) {
1166     outptr = output_buf[ctr] + output_col;
1167 
1168     /* Even part */
1169 
1170     /* Add fudge factor here for final descale. */
1171     z3 = (INT32) wsptr[0] + (ONE << (PASS1_BITS+2));
1172     z3 = LEFT_SHIFT(z3, CONST_BITS);
1173     z4 = (INT32) wsptr[4];
1174     z1 = MULTIPLY(z4, FIX(1.144122806));         /* c4 */
1175     z2 = MULTIPLY(z4, FIX(0.437016024));         /* c8 */
1176     tmp10 = z3 + z1;
1177     tmp11 = z3 - z2;
1178 
1179     tmp22 = z3 - LEFT_SHIFT(z1 - z2, 1);         /* c0 = (c4-c8)*2 */
1180 
1181     z2 = (INT32) wsptr[2];
1182     z3 = (INT32) wsptr[6];
1183 
1184     z1 = MULTIPLY(z2 + z3, FIX(0.831253876));    /* c6 */
1185     tmp12 = z1 + MULTIPLY(z2, FIX(0.513743148)); /* c2-c6 */
1186     tmp13 = z1 - MULTIPLY(z3, FIX(2.176250899)); /* c2+c6 */
1187 
1188     tmp20 = tmp10 + tmp12;
1189     tmp24 = tmp10 - tmp12;
1190     tmp21 = tmp11 + tmp13;
1191     tmp23 = tmp11 - tmp13;
1192 
1193     /* Odd part */
1194 
1195     z1 = (INT32) wsptr[1];
1196     z2 = (INT32) wsptr[3];
1197     z3 = (INT32) wsptr[5];
1198     z3 = LEFT_SHIFT(z3, CONST_BITS);
1199     z4 = (INT32) wsptr[7];
1200 
1201     tmp11 = z2 + z4;
1202     tmp13 = z2 - z4;
1203 
1204     tmp12 = MULTIPLY(tmp13, FIX(0.309016994));        /* (c3-c7)/2 */
1205 
1206     z2 = MULTIPLY(tmp11, FIX(0.951056516));           /* (c3+c7)/2 */
1207     z4 = z3 + tmp12;
1208 
1209     tmp10 = MULTIPLY(z1, FIX(1.396802247)) + z2 + z4; /* c1 */
1210     tmp14 = MULTIPLY(z1, FIX(0.221231742)) - z2 + z4; /* c9 */
1211 
1212     z2 = MULTIPLY(tmp11, FIX(0.587785252));           /* (c1-c9)/2 */
1213     z4 = z3 - tmp12 - LEFT_SHIFT(tmp13, CONST_BITS - 1);
1214 
1215     tmp12 = LEFT_SHIFT(z1 - tmp13, CONST_BITS) - z3;
1216 
1217     tmp11 = MULTIPLY(z1, FIX(1.260073511)) - z2 - z4; /* c3 */
1218     tmp13 = MULTIPLY(z1, FIX(0.642039522)) - z2 + z4; /* c7 */
1219 
1220     /* Final output stage */
1221 
1222     outptr[0] = range_limit[(int) RIGHT_SHIFT(tmp20 + tmp10,
1223                                               CONST_BITS+PASS1_BITS+3)
1224                             & RANGE_MASK];
1225     outptr[9] = range_limit[(int) RIGHT_SHIFT(tmp20 - tmp10,
1226                                               CONST_BITS+PASS1_BITS+3)
1227                             & RANGE_MASK];
1228     outptr[1] = range_limit[(int) RIGHT_SHIFT(tmp21 + tmp11,
1229                                               CONST_BITS+PASS1_BITS+3)
1230                             & RANGE_MASK];
1231     outptr[8] = range_limit[(int) RIGHT_SHIFT(tmp21 - tmp11,
1232                                               CONST_BITS+PASS1_BITS+3)
1233                             & RANGE_MASK];
1234     outptr[2] = range_limit[(int) RIGHT_SHIFT(tmp22 + tmp12,
1235                                               CONST_BITS+PASS1_BITS+3)
1236                             & RANGE_MASK];
1237     outptr[7] = range_limit[(int) RIGHT_SHIFT(tmp22 - tmp12,
1238                                               CONST_BITS+PASS1_BITS+3)
1239                             & RANGE_MASK];
1240     outptr[3] = range_limit[(int) RIGHT_SHIFT(tmp23 + tmp13,
1241                                               CONST_BITS+PASS1_BITS+3)
1242                             & RANGE_MASK];
1243     outptr[6] = range_limit[(int) RIGHT_SHIFT(tmp23 - tmp13,
1244                                               CONST_BITS+PASS1_BITS+3)
1245                             & RANGE_MASK];
1246     outptr[4] = range_limit[(int) RIGHT_SHIFT(tmp24 + tmp14,
1247                                               CONST_BITS+PASS1_BITS+3)
1248                             & RANGE_MASK];
1249     outptr[5] = range_limit[(int) RIGHT_SHIFT(tmp24 - tmp14,
1250                                               CONST_BITS+PASS1_BITS+3)
1251                             & RANGE_MASK];
1252 
1253     wsptr += 8;         /* advance pointer to next row */
1254   }
1255 }
1256 
1257 
1258 /*
1259  * Perform dequantization and inverse DCT on one block of coefficients,
1260  * producing a 11x11 output block.
1261  *
1262  * Optimized algorithm with 24 multiplications in the 1-D kernel.
1263  * cK represents sqrt(2) * cos(K*pi/22).
1264  */
1265 
1266 GLOBAL(void)
jpeg_idct_11x11(j_decompress_ptr cinfo,jpeg_component_info * compptr,JCOEFPTR coef_block,JSAMPARRAY output_buf,JDIMENSION output_col)1267 jpeg_idct_11x11 (j_decompress_ptr cinfo, jpeg_component_info * compptr,
1268                  JCOEFPTR coef_block,
1269                  JSAMPARRAY output_buf, JDIMENSION output_col)
1270 {
1271   INT32 tmp10, tmp11, tmp12, tmp13, tmp14;
1272   INT32 tmp20, tmp21, tmp22, tmp23, tmp24, tmp25;
1273   INT32 z1, z2, z3, z4;
1274   JCOEFPTR inptr;
1275   ISLOW_MULT_TYPE * quantptr;
1276   int * wsptr;
1277   JSAMPROW outptr;
1278   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
1279   int ctr;
1280   int workspace[8*11];  /* buffers data between passes */
1281   SHIFT_TEMPS
1282 
1283   /* Pass 1: process columns from input, store into work array. */
1284 
1285   inptr = coef_block;
1286   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
1287   wsptr = workspace;
1288   for (ctr = 0; ctr < 8; ctr++, inptr++, quantptr++, wsptr++) {
1289     /* Even part */
1290 
1291     tmp10 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
1292     tmp10 = LEFT_SHIFT(tmp10, CONST_BITS);
1293     /* Add fudge factor here for final descale. */
1294     tmp10 += ONE << (CONST_BITS-PASS1_BITS-1);
1295 
1296     z1 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
1297     z2 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
1298     z3 = DEQUANTIZE(inptr[DCTSIZE*6], quantptr[DCTSIZE*6]);
1299 
1300     tmp20 = MULTIPLY(z2 - z3, FIX(2.546640132));     /* c2+c4 */
1301     tmp23 = MULTIPLY(z2 - z1, FIX(0.430815045));     /* c2-c6 */
1302     z4 = z1 + z3;
1303     tmp24 = MULTIPLY(z4, - FIX(1.155664402));        /* -(c2-c10) */
1304     z4 -= z2;
1305     tmp25 = tmp10 + MULTIPLY(z4, FIX(1.356927976));  /* c2 */
1306     tmp21 = tmp20 + tmp23 + tmp25 -
1307             MULTIPLY(z2, FIX(1.821790775));          /* c2+c4+c10-c6 */
1308     tmp20 += tmp25 + MULTIPLY(z3, FIX(2.115825087)); /* c4+c6 */
1309     tmp23 += tmp25 - MULTIPLY(z1, FIX(1.513598477)); /* c6+c8 */
1310     tmp24 += tmp25;
1311     tmp22 = tmp24 - MULTIPLY(z3, FIX(0.788749120));  /* c8+c10 */
1312     tmp24 += MULTIPLY(z2, FIX(1.944413522)) -        /* c2+c8 */
1313              MULTIPLY(z1, FIX(1.390975730));         /* c4+c10 */
1314     tmp25 = tmp10 - MULTIPLY(z4, FIX(1.414213562));  /* c0 */
1315 
1316     /* Odd part */
1317 
1318     z1 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
1319     z2 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
1320     z3 = DEQUANTIZE(inptr[DCTSIZE*5], quantptr[DCTSIZE*5]);
1321     z4 = DEQUANTIZE(inptr[DCTSIZE*7], quantptr[DCTSIZE*7]);
1322 
1323     tmp11 = z1 + z2;
1324     tmp14 = MULTIPLY(tmp11 + z3 + z4, FIX(0.398430003)); /* c9 */
1325     tmp11 = MULTIPLY(tmp11, FIX(0.887983902));           /* c3-c9 */
1326     tmp12 = MULTIPLY(z1 + z3, FIX(0.670361295));         /* c5-c9 */
1327     tmp13 = tmp14 + MULTIPLY(z1 + z4, FIX(0.366151574)); /* c7-c9 */
1328     tmp10 = tmp11 + tmp12 + tmp13 -
1329             MULTIPLY(z1, FIX(0.923107866));              /* c7+c5+c3-c1-2*c9 */
1330     z1    = tmp14 - MULTIPLY(z2 + z3, FIX(1.163011579)); /* c7+c9 */
1331     tmp11 += z1 + MULTIPLY(z2, FIX(2.073276588));        /* c1+c7+3*c9-c3 */
1332     tmp12 += z1 - MULTIPLY(z3, FIX(1.192193623));        /* c3+c5-c7-c9 */
1333     z1    = MULTIPLY(z2 + z4, - FIX(1.798248910));       /* -(c1+c9) */
1334     tmp11 += z1;
1335     tmp13 += z1 + MULTIPLY(z4, FIX(2.102458632));        /* c1+c5+c9-c7 */
1336     tmp14 += MULTIPLY(z2, - FIX(1.467221301)) +          /* -(c5+c9) */
1337              MULTIPLY(z3, FIX(1.001388905)) -            /* c1-c9 */
1338              MULTIPLY(z4, FIX(1.684843907));             /* c3+c9 */
1339 
1340     /* Final output stage */
1341 
1342     wsptr[8*0]  = (int) RIGHT_SHIFT(tmp20 + tmp10, CONST_BITS-PASS1_BITS);
1343     wsptr[8*10] = (int) RIGHT_SHIFT(tmp20 - tmp10, CONST_BITS-PASS1_BITS);
1344     wsptr[8*1]  = (int) RIGHT_SHIFT(tmp21 + tmp11, CONST_BITS-PASS1_BITS);
1345     wsptr[8*9]  = (int) RIGHT_SHIFT(tmp21 - tmp11, CONST_BITS-PASS1_BITS);
1346     wsptr[8*2]  = (int) RIGHT_SHIFT(tmp22 + tmp12, CONST_BITS-PASS1_BITS);
1347     wsptr[8*8]  = (int) RIGHT_SHIFT(tmp22 - tmp12, CONST_BITS-PASS1_BITS);
1348     wsptr[8*3]  = (int) RIGHT_SHIFT(tmp23 + tmp13, CONST_BITS-PASS1_BITS);
1349     wsptr[8*7]  = (int) RIGHT_SHIFT(tmp23 - tmp13, CONST_BITS-PASS1_BITS);
1350     wsptr[8*4]  = (int) RIGHT_SHIFT(tmp24 + tmp14, CONST_BITS-PASS1_BITS);
1351     wsptr[8*6]  = (int) RIGHT_SHIFT(tmp24 - tmp14, CONST_BITS-PASS1_BITS);
1352     wsptr[8*5]  = (int) RIGHT_SHIFT(tmp25, CONST_BITS-PASS1_BITS);
1353   }
1354 
1355   /* Pass 2: process 11 rows from work array, store into output array. */
1356 
1357   wsptr = workspace;
1358   for (ctr = 0; ctr < 11; ctr++) {
1359     outptr = output_buf[ctr] + output_col;
1360 
1361     /* Even part */
1362 
1363     /* Add fudge factor here for final descale. */
1364     tmp10 = (INT32) wsptr[0] + (ONE << (PASS1_BITS+2));
1365     tmp10 = LEFT_SHIFT(tmp10, CONST_BITS);
1366 
1367     z1 = (INT32) wsptr[2];
1368     z2 = (INT32) wsptr[4];
1369     z3 = (INT32) wsptr[6];
1370 
1371     tmp20 = MULTIPLY(z2 - z3, FIX(2.546640132));     /* c2+c4 */
1372     tmp23 = MULTIPLY(z2 - z1, FIX(0.430815045));     /* c2-c6 */
1373     z4 = z1 + z3;
1374     tmp24 = MULTIPLY(z4, - FIX(1.155664402));        /* -(c2-c10) */
1375     z4 -= z2;
1376     tmp25 = tmp10 + MULTIPLY(z4, FIX(1.356927976));  /* c2 */
1377     tmp21 = tmp20 + tmp23 + tmp25 -
1378             MULTIPLY(z2, FIX(1.821790775));          /* c2+c4+c10-c6 */
1379     tmp20 += tmp25 + MULTIPLY(z3, FIX(2.115825087)); /* c4+c6 */
1380     tmp23 += tmp25 - MULTIPLY(z1, FIX(1.513598477)); /* c6+c8 */
1381     tmp24 += tmp25;
1382     tmp22 = tmp24 - MULTIPLY(z3, FIX(0.788749120));  /* c8+c10 */
1383     tmp24 += MULTIPLY(z2, FIX(1.944413522)) -        /* c2+c8 */
1384              MULTIPLY(z1, FIX(1.390975730));         /* c4+c10 */
1385     tmp25 = tmp10 - MULTIPLY(z4, FIX(1.414213562));  /* c0 */
1386 
1387     /* Odd part */
1388 
1389     z1 = (INT32) wsptr[1];
1390     z2 = (INT32) wsptr[3];
1391     z3 = (INT32) wsptr[5];
1392     z4 = (INT32) wsptr[7];
1393 
1394     tmp11 = z1 + z2;
1395     tmp14 = MULTIPLY(tmp11 + z3 + z4, FIX(0.398430003)); /* c9 */
1396     tmp11 = MULTIPLY(tmp11, FIX(0.887983902));           /* c3-c9 */
1397     tmp12 = MULTIPLY(z1 + z3, FIX(0.670361295));         /* c5-c9 */
1398     tmp13 = tmp14 + MULTIPLY(z1 + z4, FIX(0.366151574)); /* c7-c9 */
1399     tmp10 = tmp11 + tmp12 + tmp13 -
1400             MULTIPLY(z1, FIX(0.923107866));              /* c7+c5+c3-c1-2*c9 */
1401     z1    = tmp14 - MULTIPLY(z2 + z3, FIX(1.163011579)); /* c7+c9 */
1402     tmp11 += z1 + MULTIPLY(z2, FIX(2.073276588));        /* c1+c7+3*c9-c3 */
1403     tmp12 += z1 - MULTIPLY(z3, FIX(1.192193623));        /* c3+c5-c7-c9 */
1404     z1    = MULTIPLY(z2 + z4, - FIX(1.798248910));       /* -(c1+c9) */
1405     tmp11 += z1;
1406     tmp13 += z1 + MULTIPLY(z4, FIX(2.102458632));        /* c1+c5+c9-c7 */
1407     tmp14 += MULTIPLY(z2, - FIX(1.467221301)) +          /* -(c5+c9) */
1408              MULTIPLY(z3, FIX(1.001388905)) -            /* c1-c9 */
1409              MULTIPLY(z4, FIX(1.684843907));             /* c3+c9 */
1410 
1411     /* Final output stage */
1412 
1413     outptr[0]  = range_limit[(int) RIGHT_SHIFT(tmp20 + tmp10,
1414                                                CONST_BITS+PASS1_BITS+3)
1415                              & RANGE_MASK];
1416     outptr[10] = range_limit[(int) RIGHT_SHIFT(tmp20 - tmp10,
1417                                                CONST_BITS+PASS1_BITS+3)
1418                              & RANGE_MASK];
1419     outptr[1]  = range_limit[(int) RIGHT_SHIFT(tmp21 + tmp11,
1420                                                CONST_BITS+PASS1_BITS+3)
1421                              & RANGE_MASK];
1422     outptr[9]  = range_limit[(int) RIGHT_SHIFT(tmp21 - tmp11,
1423                                                CONST_BITS+PASS1_BITS+3)
1424                              & RANGE_MASK];
1425     outptr[2]  = range_limit[(int) RIGHT_SHIFT(tmp22 + tmp12,
1426                                                CONST_BITS+PASS1_BITS+3)
1427                              & RANGE_MASK];
1428     outptr[8]  = range_limit[(int) RIGHT_SHIFT(tmp22 - tmp12,
1429                                                CONST_BITS+PASS1_BITS+3)
1430                              & RANGE_MASK];
1431     outptr[3]  = range_limit[(int) RIGHT_SHIFT(tmp23 + tmp13,
1432                                                CONST_BITS+PASS1_BITS+3)
1433                              & RANGE_MASK];
1434     outptr[7]  = range_limit[(int) RIGHT_SHIFT(tmp23 - tmp13,
1435                                                CONST_BITS+PASS1_BITS+3)
1436                              & RANGE_MASK];
1437     outptr[4]  = range_limit[(int) RIGHT_SHIFT(tmp24 + tmp14,
1438                                                CONST_BITS+PASS1_BITS+3)
1439                              & RANGE_MASK];
1440     outptr[6]  = range_limit[(int) RIGHT_SHIFT(tmp24 - tmp14,
1441                                                CONST_BITS+PASS1_BITS+3)
1442                              & RANGE_MASK];
1443     outptr[5]  = range_limit[(int) RIGHT_SHIFT(tmp25,
1444                                                CONST_BITS+PASS1_BITS+3)
1445                              & RANGE_MASK];
1446 
1447     wsptr += 8;         /* advance pointer to next row */
1448   }
1449 }
1450 
1451 
1452 /*
1453  * Perform dequantization and inverse DCT on one block of coefficients,
1454  * producing a 12x12 output block.
1455  *
1456  * Optimized algorithm with 15 multiplications in the 1-D kernel.
1457  * cK represents sqrt(2) * cos(K*pi/24).
1458  */
1459 
1460 GLOBAL(void)
jpeg_idct_12x12(j_decompress_ptr cinfo,jpeg_component_info * compptr,JCOEFPTR coef_block,JSAMPARRAY output_buf,JDIMENSION output_col)1461 jpeg_idct_12x12 (j_decompress_ptr cinfo, jpeg_component_info * compptr,
1462                  JCOEFPTR coef_block,
1463                  JSAMPARRAY output_buf, JDIMENSION output_col)
1464 {
1465   INT32 tmp10, tmp11, tmp12, tmp13, tmp14, tmp15;
1466   INT32 tmp20, tmp21, tmp22, tmp23, tmp24, tmp25;
1467   INT32 z1, z2, z3, z4;
1468   JCOEFPTR inptr;
1469   ISLOW_MULT_TYPE * quantptr;
1470   int * wsptr;
1471   JSAMPROW outptr;
1472   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
1473   int ctr;
1474   int workspace[8*12];  /* buffers data between passes */
1475   SHIFT_TEMPS
1476 
1477   /* Pass 1: process columns from input, store into work array. */
1478 
1479   inptr = coef_block;
1480   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
1481   wsptr = workspace;
1482   for (ctr = 0; ctr < 8; ctr++, inptr++, quantptr++, wsptr++) {
1483     /* Even part */
1484 
1485     z3 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
1486     z3 = LEFT_SHIFT(z3, CONST_BITS);
1487     /* Add fudge factor here for final descale. */
1488     z3 += ONE << (CONST_BITS-PASS1_BITS-1);
1489 
1490     z4 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
1491     z4 = MULTIPLY(z4, FIX(1.224744871)); /* c4 */
1492 
1493     tmp10 = z3 + z4;
1494     tmp11 = z3 - z4;
1495 
1496     z1 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
1497     z4 = MULTIPLY(z1, FIX(1.366025404)); /* c2 */
1498     z1 = LEFT_SHIFT(z1, CONST_BITS);
1499     z2 = DEQUANTIZE(inptr[DCTSIZE*6], quantptr[DCTSIZE*6]);
1500     z2 = LEFT_SHIFT(z2, CONST_BITS);
1501 
1502     tmp12 = z1 - z2;
1503 
1504     tmp21 = z3 + tmp12;
1505     tmp24 = z3 - tmp12;
1506 
1507     tmp12 = z4 + z2;
1508 
1509     tmp20 = tmp10 + tmp12;
1510     tmp25 = tmp10 - tmp12;
1511 
1512     tmp12 = z4 - z1 - z2;
1513 
1514     tmp22 = tmp11 + tmp12;
1515     tmp23 = tmp11 - tmp12;
1516 
1517     /* Odd part */
1518 
1519     z1 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
1520     z2 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
1521     z3 = DEQUANTIZE(inptr[DCTSIZE*5], quantptr[DCTSIZE*5]);
1522     z4 = DEQUANTIZE(inptr[DCTSIZE*7], quantptr[DCTSIZE*7]);
1523 
1524     tmp11 = MULTIPLY(z2, FIX(1.306562965));                  /* c3 */
1525     tmp14 = MULTIPLY(z2, - FIX_0_541196100);                 /* -c9 */
1526 
1527     tmp10 = z1 + z3;
1528     tmp15 = MULTIPLY(tmp10 + z4, FIX(0.860918669));          /* c7 */
1529     tmp12 = tmp15 + MULTIPLY(tmp10, FIX(0.261052384));       /* c5-c7 */
1530     tmp10 = tmp12 + tmp11 + MULTIPLY(z1, FIX(0.280143716));  /* c1-c5 */
1531     tmp13 = MULTIPLY(z3 + z4, - FIX(1.045510580));           /* -(c7+c11) */
1532     tmp12 += tmp13 + tmp14 - MULTIPLY(z3, FIX(1.478575242)); /* c1+c5-c7-c11 */
1533     tmp13 += tmp15 - tmp11 + MULTIPLY(z4, FIX(1.586706681)); /* c1+c11 */
1534     tmp15 += tmp14 - MULTIPLY(z1, FIX(0.676326758)) -        /* c7-c11 */
1535              MULTIPLY(z4, FIX(1.982889723));                 /* c5+c7 */
1536 
1537     z1 -= z4;
1538     z2 -= z3;
1539     z3 = MULTIPLY(z1 + z2, FIX_0_541196100);                 /* c9 */
1540     tmp11 = z3 + MULTIPLY(z1, FIX_0_765366865);              /* c3-c9 */
1541     tmp14 = z3 - MULTIPLY(z2, FIX_1_847759065);              /* c3+c9 */
1542 
1543     /* Final output stage */
1544 
1545     wsptr[8*0]  = (int) RIGHT_SHIFT(tmp20 + tmp10, CONST_BITS-PASS1_BITS);
1546     wsptr[8*11] = (int) RIGHT_SHIFT(tmp20 - tmp10, CONST_BITS-PASS1_BITS);
1547     wsptr[8*1]  = (int) RIGHT_SHIFT(tmp21 + tmp11, CONST_BITS-PASS1_BITS);
1548     wsptr[8*10] = (int) RIGHT_SHIFT(tmp21 - tmp11, CONST_BITS-PASS1_BITS);
1549     wsptr[8*2]  = (int) RIGHT_SHIFT(tmp22 + tmp12, CONST_BITS-PASS1_BITS);
1550     wsptr[8*9]  = (int) RIGHT_SHIFT(tmp22 - tmp12, CONST_BITS-PASS1_BITS);
1551     wsptr[8*3]  = (int) RIGHT_SHIFT(tmp23 + tmp13, CONST_BITS-PASS1_BITS);
1552     wsptr[8*8]  = (int) RIGHT_SHIFT(tmp23 - tmp13, CONST_BITS-PASS1_BITS);
1553     wsptr[8*4]  = (int) RIGHT_SHIFT(tmp24 + tmp14, CONST_BITS-PASS1_BITS);
1554     wsptr[8*7]  = (int) RIGHT_SHIFT(tmp24 - tmp14, CONST_BITS-PASS1_BITS);
1555     wsptr[8*5]  = (int) RIGHT_SHIFT(tmp25 + tmp15, CONST_BITS-PASS1_BITS);
1556     wsptr[8*6]  = (int) RIGHT_SHIFT(tmp25 - tmp15, CONST_BITS-PASS1_BITS);
1557   }
1558 
1559   /* Pass 2: process 12 rows from work array, store into output array. */
1560 
1561   wsptr = workspace;
1562   for (ctr = 0; ctr < 12; ctr++) {
1563     outptr = output_buf[ctr] + output_col;
1564 
1565     /* Even part */
1566 
1567     /* Add fudge factor here for final descale. */
1568     z3 = (INT32) wsptr[0] + (ONE << (PASS1_BITS+2));
1569     z3 = LEFT_SHIFT(z3, CONST_BITS);
1570 
1571     z4 = (INT32) wsptr[4];
1572     z4 = MULTIPLY(z4, FIX(1.224744871)); /* c4 */
1573 
1574     tmp10 = z3 + z4;
1575     tmp11 = z3 - z4;
1576 
1577     z1 = (INT32) wsptr[2];
1578     z4 = MULTIPLY(z1, FIX(1.366025404)); /* c2 */
1579     z1 = LEFT_SHIFT(z1, CONST_BITS);
1580     z2 = (INT32) wsptr[6];
1581     z2 = LEFT_SHIFT(z2, CONST_BITS);
1582 
1583     tmp12 = z1 - z2;
1584 
1585     tmp21 = z3 + tmp12;
1586     tmp24 = z3 - tmp12;
1587 
1588     tmp12 = z4 + z2;
1589 
1590     tmp20 = tmp10 + tmp12;
1591     tmp25 = tmp10 - tmp12;
1592 
1593     tmp12 = z4 - z1 - z2;
1594 
1595     tmp22 = tmp11 + tmp12;
1596     tmp23 = tmp11 - tmp12;
1597 
1598     /* Odd part */
1599 
1600     z1 = (INT32) wsptr[1];
1601     z2 = (INT32) wsptr[3];
1602     z3 = (INT32) wsptr[5];
1603     z4 = (INT32) wsptr[7];
1604 
1605     tmp11 = MULTIPLY(z2, FIX(1.306562965));                  /* c3 */
1606     tmp14 = MULTIPLY(z2, - FIX_0_541196100);                 /* -c9 */
1607 
1608     tmp10 = z1 + z3;
1609     tmp15 = MULTIPLY(tmp10 + z4, FIX(0.860918669));          /* c7 */
1610     tmp12 = tmp15 + MULTIPLY(tmp10, FIX(0.261052384));       /* c5-c7 */
1611     tmp10 = tmp12 + tmp11 + MULTIPLY(z1, FIX(0.280143716));  /* c1-c5 */
1612     tmp13 = MULTIPLY(z3 + z4, - FIX(1.045510580));           /* -(c7+c11) */
1613     tmp12 += tmp13 + tmp14 - MULTIPLY(z3, FIX(1.478575242)); /* c1+c5-c7-c11 */
1614     tmp13 += tmp15 - tmp11 + MULTIPLY(z4, FIX(1.586706681)); /* c1+c11 */
1615     tmp15 += tmp14 - MULTIPLY(z1, FIX(0.676326758)) -        /* c7-c11 */
1616              MULTIPLY(z4, FIX(1.982889723));                 /* c5+c7 */
1617 
1618     z1 -= z4;
1619     z2 -= z3;
1620     z3 = MULTIPLY(z1 + z2, FIX_0_541196100);                 /* c9 */
1621     tmp11 = z3 + MULTIPLY(z1, FIX_0_765366865);              /* c3-c9 */
1622     tmp14 = z3 - MULTIPLY(z2, FIX_1_847759065);              /* c3+c9 */
1623 
1624     /* Final output stage */
1625 
1626     outptr[0]  = range_limit[(int) RIGHT_SHIFT(tmp20 + tmp10,
1627                                                CONST_BITS+PASS1_BITS+3)
1628                              & RANGE_MASK];
1629     outptr[11] = range_limit[(int) RIGHT_SHIFT(tmp20 - tmp10,
1630                                                CONST_BITS+PASS1_BITS+3)
1631                              & RANGE_MASK];
1632     outptr[1]  = range_limit[(int) RIGHT_SHIFT(tmp21 + tmp11,
1633                                                CONST_BITS+PASS1_BITS+3)
1634                              & RANGE_MASK];
1635     outptr[10] = range_limit[(int) RIGHT_SHIFT(tmp21 - tmp11,
1636                                                CONST_BITS+PASS1_BITS+3)
1637                              & RANGE_MASK];
1638     outptr[2]  = range_limit[(int) RIGHT_SHIFT(tmp22 + tmp12,
1639                                                CONST_BITS+PASS1_BITS+3)
1640                              & RANGE_MASK];
1641     outptr[9]  = range_limit[(int) RIGHT_SHIFT(tmp22 - tmp12,
1642                                                CONST_BITS+PASS1_BITS+3)
1643                              & RANGE_MASK];
1644     outptr[3]  = range_limit[(int) RIGHT_SHIFT(tmp23 + tmp13,
1645                                                CONST_BITS+PASS1_BITS+3)
1646                              & RANGE_MASK];
1647     outptr[8]  = range_limit[(int) RIGHT_SHIFT(tmp23 - tmp13,
1648                                                CONST_BITS+PASS1_BITS+3)
1649                              & RANGE_MASK];
1650     outptr[4]  = range_limit[(int) RIGHT_SHIFT(tmp24 + tmp14,
1651                                                CONST_BITS+PASS1_BITS+3)
1652                              & RANGE_MASK];
1653     outptr[7]  = range_limit[(int) RIGHT_SHIFT(tmp24 - tmp14,
1654                                                CONST_BITS+PASS1_BITS+3)
1655                              & RANGE_MASK];
1656     outptr[5]  = range_limit[(int) RIGHT_SHIFT(tmp25 + tmp15,
1657                                                CONST_BITS+PASS1_BITS+3)
1658                              & RANGE_MASK];
1659     outptr[6]  = range_limit[(int) RIGHT_SHIFT(tmp25 - tmp15,
1660                                                CONST_BITS+PASS1_BITS+3)
1661                              & RANGE_MASK];
1662 
1663     wsptr += 8;         /* advance pointer to next row */
1664   }
1665 }
1666 
1667 
1668 /*
1669  * Perform dequantization and inverse DCT on one block of coefficients,
1670  * producing a 13x13 output block.
1671  *
1672  * Optimized algorithm with 29 multiplications in the 1-D kernel.
1673  * cK represents sqrt(2) * cos(K*pi/26).
1674  */
1675 
1676 GLOBAL(void)
jpeg_idct_13x13(j_decompress_ptr cinfo,jpeg_component_info * compptr,JCOEFPTR coef_block,JSAMPARRAY output_buf,JDIMENSION output_col)1677 jpeg_idct_13x13 (j_decompress_ptr cinfo, jpeg_component_info * compptr,
1678                  JCOEFPTR coef_block,
1679                  JSAMPARRAY output_buf, JDIMENSION output_col)
1680 {
1681   INT32 tmp10, tmp11, tmp12, tmp13, tmp14, tmp15;
1682   INT32 tmp20, tmp21, tmp22, tmp23, tmp24, tmp25, tmp26;
1683   INT32 z1, z2, z3, z4;
1684   JCOEFPTR inptr;
1685   ISLOW_MULT_TYPE * quantptr;
1686   int * wsptr;
1687   JSAMPROW outptr;
1688   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
1689   int ctr;
1690   int workspace[8*13];  /* buffers data between passes */
1691   SHIFT_TEMPS
1692 
1693   /* Pass 1: process columns from input, store into work array. */
1694 
1695   inptr = coef_block;
1696   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
1697   wsptr = workspace;
1698   for (ctr = 0; ctr < 8; ctr++, inptr++, quantptr++, wsptr++) {
1699     /* Even part */
1700 
1701     z1 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
1702     z1 = LEFT_SHIFT(z1, CONST_BITS);
1703     /* Add fudge factor here for final descale. */
1704     z1 += ONE << (CONST_BITS-PASS1_BITS-1);
1705 
1706     z2 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
1707     z3 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
1708     z4 = DEQUANTIZE(inptr[DCTSIZE*6], quantptr[DCTSIZE*6]);
1709 
1710     tmp10 = z3 + z4;
1711     tmp11 = z3 - z4;
1712 
1713     tmp12 = MULTIPLY(tmp10, FIX(1.155388986));                /* (c4+c6)/2 */
1714     tmp13 = MULTIPLY(tmp11, FIX(0.096834934)) + z1;           /* (c4-c6)/2 */
1715 
1716     tmp20 = MULTIPLY(z2, FIX(1.373119086)) + tmp12 + tmp13;   /* c2 */
1717     tmp22 = MULTIPLY(z2, FIX(0.501487041)) - tmp12 + tmp13;   /* c10 */
1718 
1719     tmp12 = MULTIPLY(tmp10, FIX(0.316450131));                /* (c8-c12)/2 */
1720     tmp13 = MULTIPLY(tmp11, FIX(0.486914739)) + z1;           /* (c8+c12)/2 */
1721 
1722     tmp21 = MULTIPLY(z2, FIX(1.058554052)) - tmp12 + tmp13;   /* c6 */
1723     tmp25 = MULTIPLY(z2, - FIX(1.252223920)) + tmp12 + tmp13; /* c4 */
1724 
1725     tmp12 = MULTIPLY(tmp10, FIX(0.435816023));                /* (c2-c10)/2 */
1726     tmp13 = MULTIPLY(tmp11, FIX(0.937303064)) - z1;           /* (c2+c10)/2 */
1727 
1728     tmp23 = MULTIPLY(z2, - FIX(0.170464608)) - tmp12 - tmp13; /* c12 */
1729     tmp24 = MULTIPLY(z2, - FIX(0.803364869)) + tmp12 - tmp13; /* c8 */
1730 
1731     tmp26 = MULTIPLY(tmp11 - z2, FIX(1.414213562)) + z1;      /* c0 */
1732 
1733     /* Odd part */
1734 
1735     z1 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
1736     z2 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
1737     z3 = DEQUANTIZE(inptr[DCTSIZE*5], quantptr[DCTSIZE*5]);
1738     z4 = DEQUANTIZE(inptr[DCTSIZE*7], quantptr[DCTSIZE*7]);
1739 
1740     tmp11 = MULTIPLY(z1 + z2, FIX(1.322312651));     /* c3 */
1741     tmp12 = MULTIPLY(z1 + z3, FIX(1.163874945));     /* c5 */
1742     tmp15 = z1 + z4;
1743     tmp13 = MULTIPLY(tmp15, FIX(0.937797057));       /* c7 */
1744     tmp10 = tmp11 + tmp12 + tmp13 -
1745             MULTIPLY(z1, FIX(2.020082300));          /* c7+c5+c3-c1 */
1746     tmp14 = MULTIPLY(z2 + z3, - FIX(0.338443458));   /* -c11 */
1747     tmp11 += tmp14 + MULTIPLY(z2, FIX(0.837223564)); /* c5+c9+c11-c3 */
1748     tmp12 += tmp14 - MULTIPLY(z3, FIX(1.572116027)); /* c1+c5-c9-c11 */
1749     tmp14 = MULTIPLY(z2 + z4, - FIX(1.163874945));   /* -c5 */
1750     tmp11 += tmp14;
1751     tmp13 += tmp14 + MULTIPLY(z4, FIX(2.205608352)); /* c3+c5+c9-c7 */
1752     tmp14 = MULTIPLY(z3 + z4, - FIX(0.657217813));   /* -c9 */
1753     tmp12 += tmp14;
1754     tmp13 += tmp14;
1755     tmp15 = MULTIPLY(tmp15, FIX(0.338443458));       /* c11 */
1756     tmp14 = tmp15 + MULTIPLY(z1, FIX(0.318774355)) - /* c9-c11 */
1757             MULTIPLY(z2, FIX(0.466105296));          /* c1-c7 */
1758     z1    = MULTIPLY(z3 - z2, FIX(0.937797057));     /* c7 */
1759     tmp14 += z1;
1760     tmp15 += z1 + MULTIPLY(z3, FIX(0.384515595)) -   /* c3-c7 */
1761              MULTIPLY(z4, FIX(1.742345811));         /* c1+c11 */
1762 
1763     /* Final output stage */
1764 
1765     wsptr[8*0]  = (int) RIGHT_SHIFT(tmp20 + tmp10, CONST_BITS-PASS1_BITS);
1766     wsptr[8*12] = (int) RIGHT_SHIFT(tmp20 - tmp10, CONST_BITS-PASS1_BITS);
1767     wsptr[8*1]  = (int) RIGHT_SHIFT(tmp21 + tmp11, CONST_BITS-PASS1_BITS);
1768     wsptr[8*11] = (int) RIGHT_SHIFT(tmp21 - tmp11, CONST_BITS-PASS1_BITS);
1769     wsptr[8*2]  = (int) RIGHT_SHIFT(tmp22 + tmp12, CONST_BITS-PASS1_BITS);
1770     wsptr[8*10] = (int) RIGHT_SHIFT(tmp22 - tmp12, CONST_BITS-PASS1_BITS);
1771     wsptr[8*3]  = (int) RIGHT_SHIFT(tmp23 + tmp13, CONST_BITS-PASS1_BITS);
1772     wsptr[8*9]  = (int) RIGHT_SHIFT(tmp23 - tmp13, CONST_BITS-PASS1_BITS);
1773     wsptr[8*4]  = (int) RIGHT_SHIFT(tmp24 + tmp14, CONST_BITS-PASS1_BITS);
1774     wsptr[8*8]  = (int) RIGHT_SHIFT(tmp24 - tmp14, CONST_BITS-PASS1_BITS);
1775     wsptr[8*5]  = (int) RIGHT_SHIFT(tmp25 + tmp15, CONST_BITS-PASS1_BITS);
1776     wsptr[8*7]  = (int) RIGHT_SHIFT(tmp25 - tmp15, CONST_BITS-PASS1_BITS);
1777     wsptr[8*6]  = (int) RIGHT_SHIFT(tmp26, CONST_BITS-PASS1_BITS);
1778   }
1779 
1780   /* Pass 2: process 13 rows from work array, store into output array. */
1781 
1782   wsptr = workspace;
1783   for (ctr = 0; ctr < 13; ctr++) {
1784     outptr = output_buf[ctr] + output_col;
1785 
1786     /* Even part */
1787 
1788     /* Add fudge factor here for final descale. */
1789     z1 = (INT32) wsptr[0] + (ONE << (PASS1_BITS+2));
1790     z1 = LEFT_SHIFT(z1, CONST_BITS);
1791 
1792     z2 = (INT32) wsptr[2];
1793     z3 = (INT32) wsptr[4];
1794     z4 = (INT32) wsptr[6];
1795 
1796     tmp10 = z3 + z4;
1797     tmp11 = z3 - z4;
1798 
1799     tmp12 = MULTIPLY(tmp10, FIX(1.155388986));                /* (c4+c6)/2 */
1800     tmp13 = MULTIPLY(tmp11, FIX(0.096834934)) + z1;           /* (c4-c6)/2 */
1801 
1802     tmp20 = MULTIPLY(z2, FIX(1.373119086)) + tmp12 + tmp13;   /* c2 */
1803     tmp22 = MULTIPLY(z2, FIX(0.501487041)) - tmp12 + tmp13;   /* c10 */
1804 
1805     tmp12 = MULTIPLY(tmp10, FIX(0.316450131));                /* (c8-c12)/2 */
1806     tmp13 = MULTIPLY(tmp11, FIX(0.486914739)) + z1;           /* (c8+c12)/2 */
1807 
1808     tmp21 = MULTIPLY(z2, FIX(1.058554052)) - tmp12 + tmp13;   /* c6 */
1809     tmp25 = MULTIPLY(z2, - FIX(1.252223920)) + tmp12 + tmp13; /* c4 */
1810 
1811     tmp12 = MULTIPLY(tmp10, FIX(0.435816023));                /* (c2-c10)/2 */
1812     tmp13 = MULTIPLY(tmp11, FIX(0.937303064)) - z1;           /* (c2+c10)/2 */
1813 
1814     tmp23 = MULTIPLY(z2, - FIX(0.170464608)) - tmp12 - tmp13; /* c12 */
1815     tmp24 = MULTIPLY(z2, - FIX(0.803364869)) + tmp12 - tmp13; /* c8 */
1816 
1817     tmp26 = MULTIPLY(tmp11 - z2, FIX(1.414213562)) + z1;      /* c0 */
1818 
1819     /* Odd part */
1820 
1821     z1 = (INT32) wsptr[1];
1822     z2 = (INT32) wsptr[3];
1823     z3 = (INT32) wsptr[5];
1824     z4 = (INT32) wsptr[7];
1825 
1826     tmp11 = MULTIPLY(z1 + z2, FIX(1.322312651));     /* c3 */
1827     tmp12 = MULTIPLY(z1 + z3, FIX(1.163874945));     /* c5 */
1828     tmp15 = z1 + z4;
1829     tmp13 = MULTIPLY(tmp15, FIX(0.937797057));       /* c7 */
1830     tmp10 = tmp11 + tmp12 + tmp13 -
1831             MULTIPLY(z1, FIX(2.020082300));          /* c7+c5+c3-c1 */
1832     tmp14 = MULTIPLY(z2 + z3, - FIX(0.338443458));   /* -c11 */
1833     tmp11 += tmp14 + MULTIPLY(z2, FIX(0.837223564)); /* c5+c9+c11-c3 */
1834     tmp12 += tmp14 - MULTIPLY(z3, FIX(1.572116027)); /* c1+c5-c9-c11 */
1835     tmp14 = MULTIPLY(z2 + z4, - FIX(1.163874945));   /* -c5 */
1836     tmp11 += tmp14;
1837     tmp13 += tmp14 + MULTIPLY(z4, FIX(2.205608352)); /* c3+c5+c9-c7 */
1838     tmp14 = MULTIPLY(z3 + z4, - FIX(0.657217813));   /* -c9 */
1839     tmp12 += tmp14;
1840     tmp13 += tmp14;
1841     tmp15 = MULTIPLY(tmp15, FIX(0.338443458));       /* c11 */
1842     tmp14 = tmp15 + MULTIPLY(z1, FIX(0.318774355)) - /* c9-c11 */
1843             MULTIPLY(z2, FIX(0.466105296));          /* c1-c7 */
1844     z1    = MULTIPLY(z3 - z2, FIX(0.937797057));     /* c7 */
1845     tmp14 += z1;
1846     tmp15 += z1 + MULTIPLY(z3, FIX(0.384515595)) -   /* c3-c7 */
1847              MULTIPLY(z4, FIX(1.742345811));         /* c1+c11 */
1848 
1849     /* Final output stage */
1850 
1851     outptr[0]  = range_limit[(int) RIGHT_SHIFT(tmp20 + tmp10,
1852                                                CONST_BITS+PASS1_BITS+3)
1853                              & RANGE_MASK];
1854     outptr[12] = range_limit[(int) RIGHT_SHIFT(tmp20 - tmp10,
1855                                                CONST_BITS+PASS1_BITS+3)
1856                              & RANGE_MASK];
1857     outptr[1]  = range_limit[(int) RIGHT_SHIFT(tmp21 + tmp11,
1858                                                CONST_BITS+PASS1_BITS+3)
1859                              & RANGE_MASK];
1860     outptr[11] = range_limit[(int) RIGHT_SHIFT(tmp21 - tmp11,
1861                                                CONST_BITS+PASS1_BITS+3)
1862                              & RANGE_MASK];
1863     outptr[2]  = range_limit[(int) RIGHT_SHIFT(tmp22 + tmp12,
1864                                                CONST_BITS+PASS1_BITS+3)
1865                              & RANGE_MASK];
1866     outptr[10] = range_limit[(int) RIGHT_SHIFT(tmp22 - tmp12,
1867                                                CONST_BITS+PASS1_BITS+3)
1868                              & RANGE_MASK];
1869     outptr[3]  = range_limit[(int) RIGHT_SHIFT(tmp23 + tmp13,
1870                                                CONST_BITS+PASS1_BITS+3)
1871                              & RANGE_MASK];
1872     outptr[9]  = range_limit[(int) RIGHT_SHIFT(tmp23 - tmp13,
1873                                                CONST_BITS+PASS1_BITS+3)
1874                              & RANGE_MASK];
1875     outptr[4]  = range_limit[(int) RIGHT_SHIFT(tmp24 + tmp14,
1876                                                CONST_BITS+PASS1_BITS+3)
1877                              & RANGE_MASK];
1878     outptr[8]  = range_limit[(int) RIGHT_SHIFT(tmp24 - tmp14,
1879                                                CONST_BITS+PASS1_BITS+3)
1880                              & RANGE_MASK];
1881     outptr[5]  = range_limit[(int) RIGHT_SHIFT(tmp25 + tmp15,
1882                                                CONST_BITS+PASS1_BITS+3)
1883                              & RANGE_MASK];
1884     outptr[7]  = range_limit[(int) RIGHT_SHIFT(tmp25 - tmp15,
1885                                                CONST_BITS+PASS1_BITS+3)
1886                              & RANGE_MASK];
1887     outptr[6]  = range_limit[(int) RIGHT_SHIFT(tmp26,
1888                                                CONST_BITS+PASS1_BITS+3)
1889                              & RANGE_MASK];
1890 
1891     wsptr += 8;         /* advance pointer to next row */
1892   }
1893 }
1894 
1895 
1896 /*
1897  * Perform dequantization and inverse DCT on one block of coefficients,
1898  * producing a 14x14 output block.
1899  *
1900  * Optimized algorithm with 20 multiplications in the 1-D kernel.
1901  * cK represents sqrt(2) * cos(K*pi/28).
1902  */
1903 
1904 GLOBAL(void)
jpeg_idct_14x14(j_decompress_ptr cinfo,jpeg_component_info * compptr,JCOEFPTR coef_block,JSAMPARRAY output_buf,JDIMENSION output_col)1905 jpeg_idct_14x14 (j_decompress_ptr cinfo, jpeg_component_info * compptr,
1906                  JCOEFPTR coef_block,
1907                  JSAMPARRAY output_buf, JDIMENSION output_col)
1908 {
1909   INT32 tmp10, tmp11, tmp12, tmp13, tmp14, tmp15, tmp16;
1910   INT32 tmp20, tmp21, tmp22, tmp23, tmp24, tmp25, tmp26;
1911   INT32 z1, z2, z3, z4;
1912   JCOEFPTR inptr;
1913   ISLOW_MULT_TYPE * quantptr;
1914   int * wsptr;
1915   JSAMPROW outptr;
1916   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
1917   int ctr;
1918   int workspace[8*14];  /* buffers data between passes */
1919   SHIFT_TEMPS
1920 
1921   /* Pass 1: process columns from input, store into work array. */
1922 
1923   inptr = coef_block;
1924   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
1925   wsptr = workspace;
1926   for (ctr = 0; ctr < 8; ctr++, inptr++, quantptr++, wsptr++) {
1927     /* Even part */
1928 
1929     z1 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
1930     z1 = LEFT_SHIFT(z1, CONST_BITS);
1931     /* Add fudge factor here for final descale. */
1932     z1 += ONE << (CONST_BITS-PASS1_BITS-1);
1933     z4 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
1934     z2 = MULTIPLY(z4, FIX(1.274162392));         /* c4 */
1935     z3 = MULTIPLY(z4, FIX(0.314692123));         /* c12 */
1936     z4 = MULTIPLY(z4, FIX(0.881747734));         /* c8 */
1937 
1938     tmp10 = z1 + z2;
1939     tmp11 = z1 + z3;
1940     tmp12 = z1 - z4;
1941 
1942     tmp23 = RIGHT_SHIFT(z1 - LEFT_SHIFT(z2 + z3 - z4, 1),
1943                         CONST_BITS-PASS1_BITS);  /* c0 = (c4+c12-c8)*2 */
1944 
1945     z1 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
1946     z2 = DEQUANTIZE(inptr[DCTSIZE*6], quantptr[DCTSIZE*6]);
1947 
1948     z3 = MULTIPLY(z1 + z2, FIX(1.105676686));    /* c6 */
1949 
1950     tmp13 = z3 + MULTIPLY(z1, FIX(0.273079590)); /* c2-c6 */
1951     tmp14 = z3 - MULTIPLY(z2, FIX(1.719280954)); /* c6+c10 */
1952     tmp15 = MULTIPLY(z1, FIX(0.613604268)) -     /* c10 */
1953             MULTIPLY(z2, FIX(1.378756276));      /* c2 */
1954 
1955     tmp20 = tmp10 + tmp13;
1956     tmp26 = tmp10 - tmp13;
1957     tmp21 = tmp11 + tmp14;
1958     tmp25 = tmp11 - tmp14;
1959     tmp22 = tmp12 + tmp15;
1960     tmp24 = tmp12 - tmp15;
1961 
1962     /* Odd part */
1963 
1964     z1 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
1965     z2 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
1966     z3 = DEQUANTIZE(inptr[DCTSIZE*5], quantptr[DCTSIZE*5]);
1967     z4 = DEQUANTIZE(inptr[DCTSIZE*7], quantptr[DCTSIZE*7]);
1968     tmp13 = LEFT_SHIFT(z4, CONST_BITS);
1969 
1970     tmp14 = z1 + z3;
1971     tmp11 = MULTIPLY(z1 + z2, FIX(1.334852607));           /* c3 */
1972     tmp12 = MULTIPLY(tmp14, FIX(1.197448846));             /* c5 */
1973     tmp10 = tmp11 + tmp12 + tmp13 - MULTIPLY(z1, FIX(1.126980169)); /* c3+c5-c1 */
1974     tmp14 = MULTIPLY(tmp14, FIX(0.752406978));             /* c9 */
1975     tmp16 = tmp14 - MULTIPLY(z1, FIX(1.061150426));        /* c9+c11-c13 */
1976     z1    -= z2;
1977     tmp15 = MULTIPLY(z1, FIX(0.467085129)) - tmp13;        /* c11 */
1978     tmp16 += tmp15;
1979     z1    += z4;
1980     z4    = MULTIPLY(z2 + z3, - FIX(0.158341681)) - tmp13; /* -c13 */
1981     tmp11 += z4 - MULTIPLY(z2, FIX(0.424103948));          /* c3-c9-c13 */
1982     tmp12 += z4 - MULTIPLY(z3, FIX(2.373959773));          /* c3+c5-c13 */
1983     z4    = MULTIPLY(z3 - z2, FIX(1.405321284));           /* c1 */
1984     tmp14 += z4 + tmp13 - MULTIPLY(z3, FIX(1.6906431334)); /* c1+c9-c11 */
1985     tmp15 += z4 + MULTIPLY(z2, FIX(0.674957567));          /* c1+c11-c5 */
1986 
1987     tmp13 = LEFT_SHIFT(z1 - z3, PASS1_BITS);
1988 
1989     /* Final output stage */
1990 
1991     wsptr[8*0]  = (int) RIGHT_SHIFT(tmp20 + tmp10, CONST_BITS-PASS1_BITS);
1992     wsptr[8*13] = (int) RIGHT_SHIFT(tmp20 - tmp10, CONST_BITS-PASS1_BITS);
1993     wsptr[8*1]  = (int) RIGHT_SHIFT(tmp21 + tmp11, CONST_BITS-PASS1_BITS);
1994     wsptr[8*12] = (int) RIGHT_SHIFT(tmp21 - tmp11, CONST_BITS-PASS1_BITS);
1995     wsptr[8*2]  = (int) RIGHT_SHIFT(tmp22 + tmp12, CONST_BITS-PASS1_BITS);
1996     wsptr[8*11] = (int) RIGHT_SHIFT(tmp22 - tmp12, CONST_BITS-PASS1_BITS);
1997     wsptr[8*3]  = (int) (tmp23 + tmp13);
1998     wsptr[8*10] = (int) (tmp23 - tmp13);
1999     wsptr[8*4]  = (int) RIGHT_SHIFT(tmp24 + tmp14, CONST_BITS-PASS1_BITS);
2000     wsptr[8*9]  = (int) RIGHT_SHIFT(tmp24 - tmp14, CONST_BITS-PASS1_BITS);
2001     wsptr[8*5]  = (int) RIGHT_SHIFT(tmp25 + tmp15, CONST_BITS-PASS1_BITS);
2002     wsptr[8*8]  = (int) RIGHT_SHIFT(tmp25 - tmp15, CONST_BITS-PASS1_BITS);
2003     wsptr[8*6]  = (int) RIGHT_SHIFT(tmp26 + tmp16, CONST_BITS-PASS1_BITS);
2004     wsptr[8*7]  = (int) RIGHT_SHIFT(tmp26 - tmp16, CONST_BITS-PASS1_BITS);
2005   }
2006 
2007   /* Pass 2: process 14 rows from work array, store into output array. */
2008 
2009   wsptr = workspace;
2010   for (ctr = 0; ctr < 14; ctr++) {
2011     outptr = output_buf[ctr] + output_col;
2012 
2013     /* Even part */
2014 
2015     /* Add fudge factor here for final descale. */
2016     z1 = (INT32) wsptr[0] + (ONE << (PASS1_BITS+2));
2017     z1 = LEFT_SHIFT(z1, CONST_BITS);
2018     z4 = (INT32) wsptr[4];
2019     z2 = MULTIPLY(z4, FIX(1.274162392));         /* c4 */
2020     z3 = MULTIPLY(z4, FIX(0.314692123));         /* c12 */
2021     z4 = MULTIPLY(z4, FIX(0.881747734));         /* c8 */
2022 
2023     tmp10 = z1 + z2;
2024     tmp11 = z1 + z3;
2025     tmp12 = z1 - z4;
2026 
2027     tmp23 = z1 - LEFT_SHIFT(z2 + z3 - z4, 1);    /* c0 = (c4+c12-c8)*2 */
2028 
2029     z1 = (INT32) wsptr[2];
2030     z2 = (INT32) wsptr[6];
2031 
2032     z3 = MULTIPLY(z1 + z2, FIX(1.105676686));    /* c6 */
2033 
2034     tmp13 = z3 + MULTIPLY(z1, FIX(0.273079590)); /* c2-c6 */
2035     tmp14 = z3 - MULTIPLY(z2, FIX(1.719280954)); /* c6+c10 */
2036     tmp15 = MULTIPLY(z1, FIX(0.613604268)) -     /* c10 */
2037             MULTIPLY(z2, FIX(1.378756276));      /* c2 */
2038 
2039     tmp20 = tmp10 + tmp13;
2040     tmp26 = tmp10 - tmp13;
2041     tmp21 = tmp11 + tmp14;
2042     tmp25 = tmp11 - tmp14;
2043     tmp22 = tmp12 + tmp15;
2044     tmp24 = tmp12 - tmp15;
2045 
2046     /* Odd part */
2047 
2048     z1 = (INT32) wsptr[1];
2049     z2 = (INT32) wsptr[3];
2050     z3 = (INT32) wsptr[5];
2051     z4 = (INT32) wsptr[7];
2052     z4 = LEFT_SHIFT(z4, CONST_BITS);
2053 
2054     tmp14 = z1 + z3;
2055     tmp11 = MULTIPLY(z1 + z2, FIX(1.334852607));           /* c3 */
2056     tmp12 = MULTIPLY(tmp14, FIX(1.197448846));             /* c5 */
2057     tmp10 = tmp11 + tmp12 + z4 - MULTIPLY(z1, FIX(1.126980169)); /* c3+c5-c1 */
2058     tmp14 = MULTIPLY(tmp14, FIX(0.752406978));             /* c9 */
2059     tmp16 = tmp14 - MULTIPLY(z1, FIX(1.061150426));        /* c9+c11-c13 */
2060     z1    -= z2;
2061     tmp15 = MULTIPLY(z1, FIX(0.467085129)) - z4;           /* c11 */
2062     tmp16 += tmp15;
2063     tmp13 = MULTIPLY(z2 + z3, - FIX(0.158341681)) - z4;    /* -c13 */
2064     tmp11 += tmp13 - MULTIPLY(z2, FIX(0.424103948));       /* c3-c9-c13 */
2065     tmp12 += tmp13 - MULTIPLY(z3, FIX(2.373959773));       /* c3+c5-c13 */
2066     tmp13 = MULTIPLY(z3 - z2, FIX(1.405321284));           /* c1 */
2067     tmp14 += tmp13 + z4 - MULTIPLY(z3, FIX(1.6906431334)); /* c1+c9-c11 */
2068     tmp15 += tmp13 + MULTIPLY(z2, FIX(0.674957567));       /* c1+c11-c5 */
2069 
2070     tmp13 = LEFT_SHIFT(z1 - z3, CONST_BITS) + z4;
2071 
2072     /* Final output stage */
2073 
2074     outptr[0]  = range_limit[(int) RIGHT_SHIFT(tmp20 + tmp10,
2075                                                CONST_BITS+PASS1_BITS+3)
2076                              & RANGE_MASK];
2077     outptr[13] = range_limit[(int) RIGHT_SHIFT(tmp20 - tmp10,
2078                                                CONST_BITS+PASS1_BITS+3)
2079                              & RANGE_MASK];
2080     outptr[1]  = range_limit[(int) RIGHT_SHIFT(tmp21 + tmp11,
2081                                                CONST_BITS+PASS1_BITS+3)
2082                              & RANGE_MASK];
2083     outptr[12] = range_limit[(int) RIGHT_SHIFT(tmp21 - tmp11,
2084                                                CONST_BITS+PASS1_BITS+3)
2085                              & RANGE_MASK];
2086     outptr[2]  = range_limit[(int) RIGHT_SHIFT(tmp22 + tmp12,
2087                                                CONST_BITS+PASS1_BITS+3)
2088                              & RANGE_MASK];
2089     outptr[11] = range_limit[(int) RIGHT_SHIFT(tmp22 - tmp12,
2090                                                CONST_BITS+PASS1_BITS+3)
2091                              & RANGE_MASK];
2092     outptr[3]  = range_limit[(int) RIGHT_SHIFT(tmp23 + tmp13,
2093                                                CONST_BITS+PASS1_BITS+3)
2094                              & RANGE_MASK];
2095     outptr[10] = range_limit[(int) RIGHT_SHIFT(tmp23 - tmp13,
2096                                                CONST_BITS+PASS1_BITS+3)
2097                              & RANGE_MASK];
2098     outptr[4]  = range_limit[(int) RIGHT_SHIFT(tmp24 + tmp14,
2099                                                CONST_BITS+PASS1_BITS+3)
2100                              & RANGE_MASK];
2101     outptr[9]  = range_limit[(int) RIGHT_SHIFT(tmp24 - tmp14,
2102                                                CONST_BITS+PASS1_BITS+3)
2103                              & RANGE_MASK];
2104     outptr[5]  = range_limit[(int) RIGHT_SHIFT(tmp25 + tmp15,
2105                                                CONST_BITS+PASS1_BITS+3)
2106                              & RANGE_MASK];
2107     outptr[8]  = range_limit[(int) RIGHT_SHIFT(tmp25 - tmp15,
2108                                                CONST_BITS+PASS1_BITS+3)
2109                              & RANGE_MASK];
2110     outptr[6]  = range_limit[(int) RIGHT_SHIFT(tmp26 + tmp16,
2111                                                CONST_BITS+PASS1_BITS+3)
2112                              & RANGE_MASK];
2113     outptr[7]  = range_limit[(int) RIGHT_SHIFT(tmp26 - tmp16,
2114                                                CONST_BITS+PASS1_BITS+3)
2115                              & RANGE_MASK];
2116 
2117     wsptr += 8;         /* advance pointer to next row */
2118   }
2119 }
2120 
2121 
2122 /*
2123  * Perform dequantization and inverse DCT on one block of coefficients,
2124  * producing a 15x15 output block.
2125  *
2126  * Optimized algorithm with 22 multiplications in the 1-D kernel.
2127  * cK represents sqrt(2) * cos(K*pi/30).
2128  */
2129 
2130 GLOBAL(void)
jpeg_idct_15x15(j_decompress_ptr cinfo,jpeg_component_info * compptr,JCOEFPTR coef_block,JSAMPARRAY output_buf,JDIMENSION output_col)2131 jpeg_idct_15x15 (j_decompress_ptr cinfo, jpeg_component_info * compptr,
2132                  JCOEFPTR coef_block,
2133                  JSAMPARRAY output_buf, JDIMENSION output_col)
2134 {
2135   INT32 tmp10, tmp11, tmp12, tmp13, tmp14, tmp15, tmp16;
2136   INT32 tmp20, tmp21, tmp22, tmp23, tmp24, tmp25, tmp26, tmp27;
2137   INT32 z1, z2, z3, z4;
2138   JCOEFPTR inptr;
2139   ISLOW_MULT_TYPE * quantptr;
2140   int * wsptr;
2141   JSAMPROW outptr;
2142   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
2143   int ctr;
2144   int workspace[8*15];  /* buffers data between passes */
2145   SHIFT_TEMPS
2146 
2147   /* Pass 1: process columns from input, store into work array. */
2148 
2149   inptr = coef_block;
2150   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
2151   wsptr = workspace;
2152   for (ctr = 0; ctr < 8; ctr++, inptr++, quantptr++, wsptr++) {
2153     /* Even part */
2154 
2155     z1 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
2156     z1 = LEFT_SHIFT(z1, CONST_BITS);
2157     /* Add fudge factor here for final descale. */
2158     z1 += ONE << (CONST_BITS-PASS1_BITS-1);
2159 
2160     z2 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
2161     z3 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
2162     z4 = DEQUANTIZE(inptr[DCTSIZE*6], quantptr[DCTSIZE*6]);
2163 
2164     tmp10 = MULTIPLY(z4, FIX(0.437016024)); /* c12 */
2165     tmp11 = MULTIPLY(z4, FIX(1.144122806)); /* c6 */
2166 
2167     tmp12 = z1 - tmp10;
2168     tmp13 = z1 + tmp11;
2169     z1 -= LEFT_SHIFT(tmp11 - tmp10, 1);     /* c0 = (c6-c12)*2 */
2170 
2171     z4 = z2 - z3;
2172     z3 += z2;
2173     tmp10 = MULTIPLY(z3, FIX(1.337628990)); /* (c2+c4)/2 */
2174     tmp11 = MULTIPLY(z4, FIX(0.045680613)); /* (c2-c4)/2 */
2175     z2 = MULTIPLY(z2, FIX(1.439773946));    /* c4+c14 */
2176 
2177     tmp20 = tmp13 + tmp10 + tmp11;
2178     tmp23 = tmp12 - tmp10 + tmp11 + z2;
2179 
2180     tmp10 = MULTIPLY(z3, FIX(0.547059574)); /* (c8+c14)/2 */
2181     tmp11 = MULTIPLY(z4, FIX(0.399234004)); /* (c8-c14)/2 */
2182 
2183     tmp25 = tmp13 - tmp10 - tmp11;
2184     tmp26 = tmp12 + tmp10 - tmp11 - z2;
2185 
2186     tmp10 = MULTIPLY(z3, FIX(0.790569415)); /* (c6+c12)/2 */
2187     tmp11 = MULTIPLY(z4, FIX(0.353553391)); /* (c6-c12)/2 */
2188 
2189     tmp21 = tmp12 + tmp10 + tmp11;
2190     tmp24 = tmp13 - tmp10 + tmp11;
2191     tmp11 += tmp11;
2192     tmp22 = z1 + tmp11;                     /* c10 = c6-c12 */
2193     tmp27 = z1 - tmp11 - tmp11;             /* c0 = (c6-c12)*2 */
2194 
2195     /* Odd part */
2196 
2197     z1 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
2198     z2 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
2199     z4 = DEQUANTIZE(inptr[DCTSIZE*5], quantptr[DCTSIZE*5]);
2200     z3 = MULTIPLY(z4, FIX(1.224744871));                    /* c5 */
2201     z4 = DEQUANTIZE(inptr[DCTSIZE*7], quantptr[DCTSIZE*7]);
2202 
2203     tmp13 = z2 - z4;
2204     tmp15 = MULTIPLY(z1 + tmp13, FIX(0.831253876));         /* c9 */
2205     tmp11 = tmp15 + MULTIPLY(z1, FIX(0.513743148));         /* c3-c9 */
2206     tmp14 = tmp15 - MULTIPLY(tmp13, FIX(2.176250899));      /* c3+c9 */
2207 
2208     tmp13 = MULTIPLY(z2, - FIX(0.831253876));               /* -c9 */
2209     tmp15 = MULTIPLY(z2, - FIX(1.344997024));               /* -c3 */
2210     z2 = z1 - z4;
2211     tmp12 = z3 + MULTIPLY(z2, FIX(1.406466353));            /* c1 */
2212 
2213     tmp10 = tmp12 + MULTIPLY(z4, FIX(2.457431844)) - tmp15; /* c1+c7 */
2214     tmp16 = tmp12 - MULTIPLY(z1, FIX(1.112434820)) + tmp13; /* c1-c13 */
2215     tmp12 = MULTIPLY(z2, FIX(1.224744871)) - z3;            /* c5 */
2216     z2 = MULTIPLY(z1 + z4, FIX(0.575212477));               /* c11 */
2217     tmp13 += z2 + MULTIPLY(z1, FIX(0.475753014)) - z3;      /* c7-c11 */
2218     tmp15 += z2 - MULTIPLY(z4, FIX(0.869244010)) + z3;      /* c11+c13 */
2219 
2220     /* Final output stage */
2221 
2222     wsptr[8*0]  = (int) RIGHT_SHIFT(tmp20 + tmp10, CONST_BITS-PASS1_BITS);
2223     wsptr[8*14] = (int) RIGHT_SHIFT(tmp20 - tmp10, CONST_BITS-PASS1_BITS);
2224     wsptr[8*1]  = (int) RIGHT_SHIFT(tmp21 + tmp11, CONST_BITS-PASS1_BITS);
2225     wsptr[8*13] = (int) RIGHT_SHIFT(tmp21 - tmp11, CONST_BITS-PASS1_BITS);
2226     wsptr[8*2]  = (int) RIGHT_SHIFT(tmp22 + tmp12, CONST_BITS-PASS1_BITS);
2227     wsptr[8*12] = (int) RIGHT_SHIFT(tmp22 - tmp12, CONST_BITS-PASS1_BITS);
2228     wsptr[8*3]  = (int) RIGHT_SHIFT(tmp23 + tmp13, CONST_BITS-PASS1_BITS);
2229     wsptr[8*11] = (int) RIGHT_SHIFT(tmp23 - tmp13, CONST_BITS-PASS1_BITS);
2230     wsptr[8*4]  = (int) RIGHT_SHIFT(tmp24 + tmp14, CONST_BITS-PASS1_BITS);
2231     wsptr[8*10] = (int) RIGHT_SHIFT(tmp24 - tmp14, CONST_BITS-PASS1_BITS);
2232     wsptr[8*5]  = (int) RIGHT_SHIFT(tmp25 + tmp15, CONST_BITS-PASS1_BITS);
2233     wsptr[8*9]  = (int) RIGHT_SHIFT(tmp25 - tmp15, CONST_BITS-PASS1_BITS);
2234     wsptr[8*6]  = (int) RIGHT_SHIFT(tmp26 + tmp16, CONST_BITS-PASS1_BITS);
2235     wsptr[8*8]  = (int) RIGHT_SHIFT(tmp26 - tmp16, CONST_BITS-PASS1_BITS);
2236     wsptr[8*7]  = (int) RIGHT_SHIFT(tmp27, CONST_BITS-PASS1_BITS);
2237   }
2238 
2239   /* Pass 2: process 15 rows from work array, store into output array. */
2240 
2241   wsptr = workspace;
2242   for (ctr = 0; ctr < 15; ctr++) {
2243     outptr = output_buf[ctr] + output_col;
2244 
2245     /* Even part */
2246 
2247     /* Add fudge factor here for final descale. */
2248     z1 = (INT32) wsptr[0] + (ONE << (PASS1_BITS+2));
2249     z1 = LEFT_SHIFT(z1, CONST_BITS);
2250 
2251     z2 = (INT32) wsptr[2];
2252     z3 = (INT32) wsptr[4];
2253     z4 = (INT32) wsptr[6];
2254 
2255     tmp10 = MULTIPLY(z4, FIX(0.437016024)); /* c12 */
2256     tmp11 = MULTIPLY(z4, FIX(1.144122806)); /* c6 */
2257 
2258     tmp12 = z1 - tmp10;
2259     tmp13 = z1 + tmp11;
2260     z1 -= LEFT_SHIFT(tmp11 - tmp10, 1);     /* c0 = (c6-c12)*2 */
2261 
2262     z4 = z2 - z3;
2263     z3 += z2;
2264     tmp10 = MULTIPLY(z3, FIX(1.337628990)); /* (c2+c4)/2 */
2265     tmp11 = MULTIPLY(z4, FIX(0.045680613)); /* (c2-c4)/2 */
2266     z2 = MULTIPLY(z2, FIX(1.439773946));    /* c4+c14 */
2267 
2268     tmp20 = tmp13 + tmp10 + tmp11;
2269     tmp23 = tmp12 - tmp10 + tmp11 + z2;
2270 
2271     tmp10 = MULTIPLY(z3, FIX(0.547059574)); /* (c8+c14)/2 */
2272     tmp11 = MULTIPLY(z4, FIX(0.399234004)); /* (c8-c14)/2 */
2273 
2274     tmp25 = tmp13 - tmp10 - tmp11;
2275     tmp26 = tmp12 + tmp10 - tmp11 - z2;
2276 
2277     tmp10 = MULTIPLY(z3, FIX(0.790569415)); /* (c6+c12)/2 */
2278     tmp11 = MULTIPLY(z4, FIX(0.353553391)); /* (c6-c12)/2 */
2279 
2280     tmp21 = tmp12 + tmp10 + tmp11;
2281     tmp24 = tmp13 - tmp10 + tmp11;
2282     tmp11 += tmp11;
2283     tmp22 = z1 + tmp11;                     /* c10 = c6-c12 */
2284     tmp27 = z1 - tmp11 - tmp11;             /* c0 = (c6-c12)*2 */
2285 
2286     /* Odd part */
2287 
2288     z1 = (INT32) wsptr[1];
2289     z2 = (INT32) wsptr[3];
2290     z4 = (INT32) wsptr[5];
2291     z3 = MULTIPLY(z4, FIX(1.224744871));                    /* c5 */
2292     z4 = (INT32) wsptr[7];
2293 
2294     tmp13 = z2 - z4;
2295     tmp15 = MULTIPLY(z1 + tmp13, FIX(0.831253876));         /* c9 */
2296     tmp11 = tmp15 + MULTIPLY(z1, FIX(0.513743148));         /* c3-c9 */
2297     tmp14 = tmp15 - MULTIPLY(tmp13, FIX(2.176250899));      /* c3+c9 */
2298 
2299     tmp13 = MULTIPLY(z2, - FIX(0.831253876));               /* -c9 */
2300     tmp15 = MULTIPLY(z2, - FIX(1.344997024));               /* -c3 */
2301     z2 = z1 - z4;
2302     tmp12 = z3 + MULTIPLY(z2, FIX(1.406466353));            /* c1 */
2303 
2304     tmp10 = tmp12 + MULTIPLY(z4, FIX(2.457431844)) - tmp15; /* c1+c7 */
2305     tmp16 = tmp12 - MULTIPLY(z1, FIX(1.112434820)) + tmp13; /* c1-c13 */
2306     tmp12 = MULTIPLY(z2, FIX(1.224744871)) - z3;            /* c5 */
2307     z2 = MULTIPLY(z1 + z4, FIX(0.575212477));               /* c11 */
2308     tmp13 += z2 + MULTIPLY(z1, FIX(0.475753014)) - z3;      /* c7-c11 */
2309     tmp15 += z2 - MULTIPLY(z4, FIX(0.869244010)) + z3;      /* c11+c13 */
2310 
2311     /* Final output stage */
2312 
2313     outptr[0]  = range_limit[(int) RIGHT_SHIFT(tmp20 + tmp10,
2314                                                CONST_BITS+PASS1_BITS+3)
2315                              & RANGE_MASK];
2316     outptr[14] = range_limit[(int) RIGHT_SHIFT(tmp20 - tmp10,
2317                                                CONST_BITS+PASS1_BITS+3)
2318                              & RANGE_MASK];
2319     outptr[1]  = range_limit[(int) RIGHT_SHIFT(tmp21 + tmp11,
2320                                                CONST_BITS+PASS1_BITS+3)
2321                              & RANGE_MASK];
2322     outptr[13] = range_limit[(int) RIGHT_SHIFT(tmp21 - tmp11,
2323                                                CONST_BITS+PASS1_BITS+3)
2324                              & RANGE_MASK];
2325     outptr[2]  = range_limit[(int) RIGHT_SHIFT(tmp22 + tmp12,
2326                                                CONST_BITS+PASS1_BITS+3)
2327                              & RANGE_MASK];
2328     outptr[12] = range_limit[(int) RIGHT_SHIFT(tmp22 - tmp12,
2329                                                CONST_BITS+PASS1_BITS+3)
2330                              & RANGE_MASK];
2331     outptr[3]  = range_limit[(int) RIGHT_SHIFT(tmp23 + tmp13,
2332                                                CONST_BITS+PASS1_BITS+3)
2333                              & RANGE_MASK];
2334     outptr[11] = range_limit[(int) RIGHT_SHIFT(tmp23 - tmp13,
2335                                                CONST_BITS+PASS1_BITS+3)
2336                              & RANGE_MASK];
2337     outptr[4]  = range_limit[(int) RIGHT_SHIFT(tmp24 + tmp14,
2338                                                CONST_BITS+PASS1_BITS+3)
2339                              & RANGE_MASK];
2340     outptr[10] = range_limit[(int) RIGHT_SHIFT(tmp24 - tmp14,
2341                                                CONST_BITS+PASS1_BITS+3)
2342                              & RANGE_MASK];
2343     outptr[5]  = range_limit[(int) RIGHT_SHIFT(tmp25 + tmp15,
2344                                                CONST_BITS+PASS1_BITS+3)
2345                              & RANGE_MASK];
2346     outptr[9]  = range_limit[(int) RIGHT_SHIFT(tmp25 - tmp15,
2347                                                CONST_BITS+PASS1_BITS+3)
2348                              & RANGE_MASK];
2349     outptr[6]  = range_limit[(int) RIGHT_SHIFT(tmp26 + tmp16,
2350                                                CONST_BITS+PASS1_BITS+3)
2351                              & RANGE_MASK];
2352     outptr[8]  = range_limit[(int) RIGHT_SHIFT(tmp26 - tmp16,
2353                                                CONST_BITS+PASS1_BITS+3)
2354                              & RANGE_MASK];
2355     outptr[7]  = range_limit[(int) RIGHT_SHIFT(tmp27,
2356                                                CONST_BITS+PASS1_BITS+3)
2357                              & RANGE_MASK];
2358 
2359     wsptr += 8;         /* advance pointer to next row */
2360   }
2361 }
2362 
2363 
2364 /*
2365  * Perform dequantization and inverse DCT on one block of coefficients,
2366  * producing a 16x16 output block.
2367  *
2368  * Optimized algorithm with 28 multiplications in the 1-D kernel.
2369  * cK represents sqrt(2) * cos(K*pi/32).
2370  */
2371 
2372 GLOBAL(void)
jpeg_idct_16x16(j_decompress_ptr cinfo,jpeg_component_info * compptr,JCOEFPTR coef_block,JSAMPARRAY output_buf,JDIMENSION output_col)2373 jpeg_idct_16x16 (j_decompress_ptr cinfo, jpeg_component_info * compptr,
2374                  JCOEFPTR coef_block,
2375                  JSAMPARRAY output_buf, JDIMENSION output_col)
2376 {
2377   INT32 tmp0, tmp1, tmp2, tmp3, tmp10, tmp11, tmp12, tmp13;
2378   INT32 tmp20, tmp21, tmp22, tmp23, tmp24, tmp25, tmp26, tmp27;
2379   INT32 z1, z2, z3, z4;
2380   JCOEFPTR inptr;
2381   ISLOW_MULT_TYPE * quantptr;
2382   int * wsptr;
2383   JSAMPROW outptr;
2384   JSAMPLE *range_limit = IDCT_range_limit(cinfo);
2385   int ctr;
2386   int workspace[8*16];  /* buffers data between passes */
2387   SHIFT_TEMPS
2388 
2389   /* Pass 1: process columns from input, store into work array. */
2390 
2391   inptr = coef_block;
2392   quantptr = (ISLOW_MULT_TYPE *) compptr->dct_table;
2393   wsptr = workspace;
2394   for (ctr = 0; ctr < 8; ctr++, inptr++, quantptr++, wsptr++) {
2395     /* Even part */
2396 
2397     tmp0 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
2398     tmp0 = LEFT_SHIFT(tmp0, CONST_BITS);
2399     /* Add fudge factor here for final descale. */
2400     tmp0 += 1 << (CONST_BITS-PASS1_BITS-1);
2401 
2402     z1 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
2403     tmp1 = MULTIPLY(z1, FIX(1.306562965));      /* c4[16] = c2[8] */
2404     tmp2 = MULTIPLY(z1, FIX_0_541196100);       /* c12[16] = c6[8] */
2405 
2406     tmp10 = tmp0 + tmp1;
2407     tmp11 = tmp0 - tmp1;
2408     tmp12 = tmp0 + tmp2;
2409     tmp13 = tmp0 - tmp2;
2410 
2411     z1 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
2412     z2 = DEQUANTIZE(inptr[DCTSIZE*6], quantptr[DCTSIZE*6]);
2413     z3 = z1 - z2;
2414     z4 = MULTIPLY(z3, FIX(0.275899379));        /* c14[16] = c7[8] */
2415     z3 = MULTIPLY(z3, FIX(1.387039845));        /* c2[16] = c1[8] */
2416 
2417     tmp0 = z3 + MULTIPLY(z2, FIX_2_562915447);  /* (c6+c2)[16] = (c3+c1)[8] */
2418     tmp1 = z4 + MULTIPLY(z1, FIX_0_899976223);  /* (c6-c14)[16] = (c3-c7)[8] */
2419     tmp2 = z3 - MULTIPLY(z1, FIX(0.601344887)); /* (c2-c10)[16] = (c1-c5)[8] */
2420     tmp3 = z4 - MULTIPLY(z2, FIX(0.509795579)); /* (c10-c14)[16] = (c5-c7)[8] */
2421 
2422     tmp20 = tmp10 + tmp0;
2423     tmp27 = tmp10 - tmp0;
2424     tmp21 = tmp12 + tmp1;
2425     tmp26 = tmp12 - tmp1;
2426     tmp22 = tmp13 + tmp2;
2427     tmp25 = tmp13 - tmp2;
2428     tmp23 = tmp11 + tmp3;
2429     tmp24 = tmp11 - tmp3;
2430 
2431     /* Odd part */
2432 
2433     z1 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
2434     z2 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
2435     z3 = DEQUANTIZE(inptr[DCTSIZE*5], quantptr[DCTSIZE*5]);
2436     z4 = DEQUANTIZE(inptr[DCTSIZE*7], quantptr[DCTSIZE*7]);
2437 
2438     tmp11 = z1 + z3;
2439 
2440     tmp1  = MULTIPLY(z1 + z2, FIX(1.353318001));   /* c3 */
2441     tmp2  = MULTIPLY(tmp11,   FIX(1.247225013));   /* c5 */
2442     tmp3  = MULTIPLY(z1 + z4, FIX(1.093201867));   /* c7 */
2443     tmp10 = MULTIPLY(z1 - z4, FIX(0.897167586));   /* c9 */
2444     tmp11 = MULTIPLY(tmp11,   FIX(0.666655658));   /* c11 */
2445     tmp12 = MULTIPLY(z1 - z2, FIX(0.410524528));   /* c13 */
2446     tmp0  = tmp1 + tmp2 + tmp3 -
2447             MULTIPLY(z1, FIX(2.286341144));        /* c7+c5+c3-c1 */
2448     tmp13 = tmp10 + tmp11 + tmp12 -
2449             MULTIPLY(z1, FIX(1.835730603));        /* c9+c11+c13-c15 */
2450     z1    = MULTIPLY(z2 + z3, FIX(0.138617169));   /* c15 */
2451     tmp1  += z1 + MULTIPLY(z2, FIX(0.071888074));  /* c9+c11-c3-c15 */
2452     tmp2  += z1 - MULTIPLY(z3, FIX(1.125726048));  /* c5+c7+c15-c3 */
2453     z1    = MULTIPLY(z3 - z2, FIX(1.407403738));   /* c1 */
2454     tmp11 += z1 - MULTIPLY(z3, FIX(0.766367282));  /* c1+c11-c9-c13 */
2455     tmp12 += z1 + MULTIPLY(z2, FIX(1.971951411));  /* c1+c5+c13-c7 */
2456     z2    += z4;
2457     z1    = MULTIPLY(z2, - FIX(0.666655658));      /* -c11 */
2458     tmp1  += z1;
2459     tmp3  += z1 + MULTIPLY(z4, FIX(1.065388962));  /* c3+c11+c15-c7 */
2460     z2    = MULTIPLY(z2, - FIX(1.247225013));      /* -c5 */
2461     tmp10 += z2 + MULTIPLY(z4, FIX(3.141271809));  /* c1+c5+c9-c13 */
2462     tmp12 += z2;
2463     z2    = MULTIPLY(z3 + z4, - FIX(1.353318001)); /* -c3 */
2464     tmp2  += z2;
2465     tmp3  += z2;
2466     z2    = MULTIPLY(z4 - z3, FIX(0.410524528));   /* c13 */
2467     tmp10 += z2;
2468     tmp11 += z2;
2469 
2470     /* Final output stage */
2471 
2472     wsptr[8*0]  = (int) RIGHT_SHIFT(tmp20 + tmp0,  CONST_BITS-PASS1_BITS);
2473     wsptr[8*15] = (int) RIGHT_SHIFT(tmp20 - tmp0,  CONST_BITS-PASS1_BITS);
2474     wsptr[8*1]  = (int) RIGHT_SHIFT(tmp21 + tmp1,  CONST_BITS-PASS1_BITS);
2475     wsptr[8*14] = (int) RIGHT_SHIFT(tmp21 - tmp1,  CONST_BITS-PASS1_BITS);
2476     wsptr[8*2]  = (int) RIGHT_SHIFT(tmp22 + tmp2,  CONST_BITS-PASS1_BITS);
2477     wsptr[8*13] = (int) RIGHT_SHIFT(tmp22 - tmp2,  CONST_BITS-PASS1_BITS);
2478     wsptr[8*3]  = (int) RIGHT_SHIFT(tmp23 + tmp3,  CONST_BITS-PASS1_BITS);
2479     wsptr[8*12] = (int) RIGHT_SHIFT(tmp23 - tmp3,  CONST_BITS-PASS1_BITS);
2480     wsptr[8*4]  = (int) RIGHT_SHIFT(tmp24 + tmp10, CONST_BITS-PASS1_BITS);
2481     wsptr[8*11] = (int) RIGHT_SHIFT(tmp24 - tmp10, CONST_BITS-PASS1_BITS);
2482     wsptr[8*5]  = (int) RIGHT_SHIFT(tmp25 + tmp11, CONST_BITS-PASS1_BITS);
2483     wsptr[8*10] = (int) RIGHT_SHIFT(tmp25 - tmp11, CONST_BITS-PASS1_BITS);
2484     wsptr[8*6]  = (int) RIGHT_SHIFT(tmp26 + tmp12, CONST_BITS-PASS1_BITS);
2485     wsptr[8*9]  = (int) RIGHT_SHIFT(tmp26 - tmp12, CONST_BITS-PASS1_BITS);
2486     wsptr[8*7]  = (int) RIGHT_SHIFT(tmp27 + tmp13, CONST_BITS-PASS1_BITS);
2487     wsptr[8*8]  = (int) RIGHT_SHIFT(tmp27 - tmp13, CONST_BITS-PASS1_BITS);
2488   }
2489 
2490   /* Pass 2: process 16 rows from work array, store into output array. */
2491 
2492   wsptr = workspace;
2493   for (ctr = 0; ctr < 16; ctr++) {
2494     outptr = output_buf[ctr] + output_col;
2495 
2496     /* Even part */
2497 
2498     /* Add fudge factor here for final descale. */
2499     tmp0 = (INT32) wsptr[0] + (ONE << (PASS1_BITS+2));
2500     tmp0 = LEFT_SHIFT(tmp0, CONST_BITS);
2501 
2502     z1 = (INT32) wsptr[4];
2503     tmp1 = MULTIPLY(z1, FIX(1.306562965));      /* c4[16] = c2[8] */
2504     tmp2 = MULTIPLY(z1, FIX_0_541196100);       /* c12[16] = c6[8] */
2505 
2506     tmp10 = tmp0 + tmp1;
2507     tmp11 = tmp0 - tmp1;
2508     tmp12 = tmp0 + tmp2;
2509     tmp13 = tmp0 - tmp2;
2510 
2511     z1 = (INT32) wsptr[2];
2512     z2 = (INT32) wsptr[6];
2513     z3 = z1 - z2;
2514     z4 = MULTIPLY(z3, FIX(0.275899379));        /* c14[16] = c7[8] */
2515     z3 = MULTIPLY(z3, FIX(1.387039845));        /* c2[16] = c1[8] */
2516 
2517     tmp0 = z3 + MULTIPLY(z2, FIX_2_562915447);  /* (c6+c2)[16] = (c3+c1)[8] */
2518     tmp1 = z4 + MULTIPLY(z1, FIX_0_899976223);  /* (c6-c14)[16] = (c3-c7)[8] */
2519     tmp2 = z3 - MULTIPLY(z1, FIX(0.601344887)); /* (c2-c10)[16] = (c1-c5)[8] */
2520     tmp3 = z4 - MULTIPLY(z2, FIX(0.509795579)); /* (c10-c14)[16] = (c5-c7)[8] */
2521 
2522     tmp20 = tmp10 + tmp0;
2523     tmp27 = tmp10 - tmp0;
2524     tmp21 = tmp12 + tmp1;
2525     tmp26 = tmp12 - tmp1;
2526     tmp22 = tmp13 + tmp2;
2527     tmp25 = tmp13 - tmp2;
2528     tmp23 = tmp11 + tmp3;
2529     tmp24 = tmp11 - tmp3;
2530 
2531     /* Odd part */
2532 
2533     z1 = (INT32) wsptr[1];
2534     z2 = (INT32) wsptr[3];
2535     z3 = (INT32) wsptr[5];
2536     z4 = (INT32) wsptr[7];
2537 
2538     tmp11 = z1 + z3;
2539 
2540     tmp1  = MULTIPLY(z1 + z2, FIX(1.353318001));   /* c3 */
2541     tmp2  = MULTIPLY(tmp11,   FIX(1.247225013));   /* c5 */
2542     tmp3  = MULTIPLY(z1 + z4, FIX(1.093201867));   /* c7 */
2543     tmp10 = MULTIPLY(z1 - z4, FIX(0.897167586));   /* c9 */
2544     tmp11 = MULTIPLY(tmp11,   FIX(0.666655658));   /* c11 */
2545     tmp12 = MULTIPLY(z1 - z2, FIX(0.410524528));   /* c13 */
2546     tmp0  = tmp1 + tmp2 + tmp3 -
2547             MULTIPLY(z1, FIX(2.286341144));        /* c7+c5+c3-c1 */
2548     tmp13 = tmp10 + tmp11 + tmp12 -
2549             MULTIPLY(z1, FIX(1.835730603));        /* c9+c11+c13-c15 */
2550     z1    = MULTIPLY(z2 + z3, FIX(0.138617169));   /* c15 */
2551     tmp1  += z1 + MULTIPLY(z2, FIX(0.071888074));  /* c9+c11-c3-c15 */
2552     tmp2  += z1 - MULTIPLY(z3, FIX(1.125726048));  /* c5+c7+c15-c3 */
2553     z1    = MULTIPLY(z3 - z2, FIX(1.407403738));   /* c1 */
2554     tmp11 += z1 - MULTIPLY(z3, FIX(0.766367282));  /* c1+c11-c9-c13 */
2555     tmp12 += z1 + MULTIPLY(z2, FIX(1.971951411));  /* c1+c5+c13-c7 */
2556     z2    += z4;
2557     z1    = MULTIPLY(z2, - FIX(0.666655658));      /* -c11 */
2558     tmp1  += z1;
2559     tmp3  += z1 + MULTIPLY(z4, FIX(1.065388962));  /* c3+c11+c15-c7 */
2560     z2    = MULTIPLY(z2, - FIX(1.247225013));      /* -c5 */
2561     tmp10 += z2 + MULTIPLY(z4, FIX(3.141271809));  /* c1+c5+c9-c13 */
2562     tmp12 += z2;
2563     z2    = MULTIPLY(z3 + z4, - FIX(1.353318001)); /* -c3 */
2564     tmp2  += z2;
2565     tmp3  += z2;
2566     z2    = MULTIPLY(z4 - z3, FIX(0.410524528));   /* c13 */
2567     tmp10 += z2;
2568     tmp11 += z2;
2569 
2570     /* Final output stage */
2571 
2572     outptr[0]  = range_limit[(int) RIGHT_SHIFT(tmp20 + tmp0,
2573                                                CONST_BITS+PASS1_BITS+3)
2574                              & RANGE_MASK];
2575     outptr[15] = range_limit[(int) RIGHT_SHIFT(tmp20 - tmp0,
2576                                                CONST_BITS+PASS1_BITS+3)
2577                              & RANGE_MASK];
2578     outptr[1]  = range_limit[(int) RIGHT_SHIFT(tmp21 + tmp1,
2579                                                CONST_BITS+PASS1_BITS+3)
2580                              & RANGE_MASK];
2581     outptr[14] = range_limit[(int) RIGHT_SHIFT(tmp21 - tmp1,
2582                                                CONST_BITS+PASS1_BITS+3)
2583                              & RANGE_MASK];
2584     outptr[2]  = range_limit[(int) RIGHT_SHIFT(tmp22 + tmp2,
2585                                                CONST_BITS+PASS1_BITS+3)
2586                              & RANGE_MASK];
2587     outptr[13] = range_limit[(int) RIGHT_SHIFT(tmp22 - tmp2,
2588                                                CONST_BITS+PASS1_BITS+3)
2589                              & RANGE_MASK];
2590     outptr[3]  = range_limit[(int) RIGHT_SHIFT(tmp23 + tmp3,
2591                                                CONST_BITS+PASS1_BITS+3)
2592                              & RANGE_MASK];
2593     outptr[12] = range_limit[(int) RIGHT_SHIFT(tmp23 - tmp3,
2594                                                CONST_BITS+PASS1_BITS+3)
2595                              & RANGE_MASK];
2596     outptr[4]  = range_limit[(int) RIGHT_SHIFT(tmp24 + tmp10,
2597                                                CONST_BITS+PASS1_BITS+3)
2598                              & RANGE_MASK];
2599     outptr[11] = range_limit[(int) RIGHT_SHIFT(tmp24 - tmp10,
2600                                                CONST_BITS+PASS1_BITS+3)
2601                              & RANGE_MASK];
2602     outptr[5]  = range_limit[(int) RIGHT_SHIFT(tmp25 + tmp11,
2603                                                CONST_BITS+PASS1_BITS+3)
2604                              & RANGE_MASK];
2605     outptr[10] = range_limit[(int) RIGHT_SHIFT(tmp25 - tmp11,
2606                                                CONST_BITS+PASS1_BITS+3)
2607                              & RANGE_MASK];
2608     outptr[6]  = range_limit[(int) RIGHT_SHIFT(tmp26 + tmp12,
2609                                                CONST_BITS+PASS1_BITS+3)
2610                              & RANGE_MASK];
2611     outptr[9]  = range_limit[(int) RIGHT_SHIFT(tmp26 - tmp12,
2612                                                CONST_BITS+PASS1_BITS+3)
2613                              & RANGE_MASK];
2614     outptr[7]  = range_limit[(int) RIGHT_SHIFT(tmp27 + tmp13,
2615                                                CONST_BITS+PASS1_BITS+3)
2616                              & RANGE_MASK];
2617     outptr[8]  = range_limit[(int) RIGHT_SHIFT(tmp27 - tmp13,
2618                                                CONST_BITS+PASS1_BITS+3)
2619                              & RANGE_MASK];
2620 
2621     wsptr += 8;         /* advance pointer to next row */
2622   }
2623 }
2624 
2625 #endif /* IDCT_SCALING_SUPPORTED */
2626 #endif /* DCT_ISLOW_SUPPORTED */
2627