jidctint.c

   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)
 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)
 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)
 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)
 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)
 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)
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)
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)
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)
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)
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)
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)
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 */