Mercurial > cgi-bin > hgwebdir.cgi > PR > Applications > VSs > VSs__jpeg_decoder__Proj

annotate VSs_tinyjpeg/jidctflt.c @ 0:a8af8b3fc99d

find changesets by author, revision, files, or words in the commit message
initial commit
author Nina Engelhardt <nengel@mailbox.tu-berlin.de>
date Thu, 05 Jul 2012 11:35:03 +0200
parents
children
rev   line source
nengel@0 1 /*
nengel@0 2 * jidctflt.c
nengel@0 3 *
nengel@0 4 * Copyright (C) 1994-1998, Thomas G. Lane.
nengel@0 5 * This file is part of the Independent JPEG Group's software.
nengel@0 6 *
nengel@0 7 * The authors make NO WARRANTY or representation, either express or implied,
nengel@0 8 * with respect to this software, its quality, accuracy, merchantability, or
nengel@0 9 * fitness for a particular purpose. This software is provided "AS IS", and you,
nengel@0 10 * its user, assume the entire risk as to its quality and accuracy.
nengel@0 11 *
nengel@0 12 * This software is copyright (C) 1991-1998, Thomas G. Lane.
nengel@0 13 * All Rights Reserved except as specified below.
nengel@0 14 *
nengel@0 15 * Permission is hereby granted to use, copy, modify, and distribute this
nengel@0 16 * software (or portions thereof) for any purpose, without fee, subject to these
nengel@0 17 * conditions:
nengel@0 18 * (1) If any part of the source code for this software is distributed, then this
nengel@0 19 * README file must be included, with this copyright and no-warranty notice
nengel@0 20 * unaltered; and any additions, deletions, or changes to the original files
nengel@0 21 * must be clearly indicated in accompanying documentation.
nengel@0 22 * (2) If only executable code is distributed, then the accompanying
nengel@0 23 * documentation must state that "this software is based in part on the work of
nengel@0 24 * the Independent JPEG Group".
nengel@0 25 * (3) Permission for use of this software is granted only if the user accepts
nengel@0 26 * full responsibility for any undesirable consequences; the authors accept
nengel@0 27 * NO LIABILITY for damages of any kind.
nengel@0 28 *
nengel@0 29 * These conditions apply to any software derived from or based on the IJG code,
nengel@0 30 * not just to the unmodified library. If you use our work, you ought to
nengel@0 31 * acknowledge us.
nengel@0 32 *
nengel@0 33 * Permission is NOT granted for the use of any IJG author's name or company name
nengel@0 34 * in advertising or publicity relating to this software or products derived from
nengel@0 35 * it. This software may be referred to only as "the Independent JPEG Group's
nengel@0 36 * software".
nengel@0 37 *
nengel@0 38 * We specifically permit and encourage the use of this software as the basis of
nengel@0 39 * commercial products, provided that all warranty or liability claims are
nengel@0 40 * assumed by the product vendor.
nengel@0 41 *
nengel@0 42 *
nengel@0 43 * This file contains a floating-point implementation of the
nengel@0 44 * inverse DCT (Discrete Cosine Transform). In the IJG code, this routine
nengel@0 45 * must also perform dequantization of the input coefficients.
nengel@0 46 *
nengel@0 47 * This implementation should be more accurate than either of the integer
nengel@0 48 * IDCT implementations. However, it may not give the same results on all
nengel@0 49 * machines because of differences in roundoff behavior. Speed will depend
nengel@0 50 * on the hardware's floating point capacity.
nengel@0 51 *
nengel@0 52 * A 2-D IDCT can be done by 1-D IDCT on each column followed by 1-D IDCT
nengel@0 53 * on each row (or vice versa, but it's more convenient to emit a row at
nengel@0 54 * a time). Direct algorithms are also available, but they are much more
nengel@0 55 * complex and seem not to be any faster when reduced to code.
nengel@0 56 *
nengel@0 57 * This implementation is based on Arai, Agui, and Nakajima's algorithm for
nengel@0 58 * scaled DCT. Their original paper (Trans. IEICE E-71(11):1095) is in
nengel@0 59 * Japanese, but the algorithm is described in the Pennebaker & Mitchell
nengel@0 60 * JPEG textbook (see REFERENCES section in file README). The following code
nengel@0 61 * is based directly on figure 4-8 in P&M.
nengel@0 62 * While an 8-point DCT cannot be done in less than 11 multiplies, it is
nengel@0 63 * possible to arrange the computation so that many of the multiplies are
nengel@0 64 * simple scalings of the final outputs. These multiplies can then be
nengel@0 65 * folded into the multiplications or divisions by the JPEG quantization
nengel@0 66 * table entries. The AA&N method leaves only 5 multiplies and 29 adds
nengel@0 67 * to be done in the DCT itself.
nengel@0 68 * The primary disadvantage of this method is that with a fixed-point
nengel@0 69 * implementation, accuracy is lost due to imprecise representation of the
nengel@0 70 * scaled quantization values. However, that problem does not arise if
nengel@0 71 * we use floating point arithmetic.
nengel@0 72 */
nengel@0 73
nengel@0 74 #include <stdint.h>
nengel@0 75 #include "tinyjpeg-internal.h"
nengel@0 76
nengel@0 77 #define FAST_FLOAT float
nengel@0 78 #define DCTSIZE 8
nengel@0 79 #define DCTSIZE2 (DCTSIZE*DCTSIZE)
nengel@0 80
nengel@0 81 #define DEQUANTIZE(coef,quantval) (((FAST_FLOAT) (coef)) * (quantval))
nengel@0 82
nengel@0 83 static inline unsigned char descale_and_clamp(int x, int shift)
nengel@0 84 {
nengel@0 85 x += (1UL<<(shift-1));
nengel@0 86 if (x<0)
nengel@0 87 x = (x >> shift) | ((~(0UL)) << (32-(shift)));
nengel@0 88 else
nengel@0 89 x >>= shift;
nengel@0 90 x += 128;
nengel@0 91 if (x>255)
nengel@0 92 return 255;
nengel@0 93 else if (x<0)
nengel@0 94 return 0;
nengel@0 95 else
nengel@0 96 return x;
nengel@0 97 }
nengel@0 98
nengel@0 99 /*
nengel@0 100 * Perform dequantization and inverse DCT on one block of coefficients.
nengel@0 101 */
nengel@0 102 void tinyjpeg_idct_float (struct component *compptr, uint8_t *output_buf, int stride)
nengel@0 103 {
nengel@0 104 FAST_FLOAT tmp0, tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7;
nengel@0 105 FAST_FLOAT tmp10, tmp11, tmp12, tmp13;
nengel@0 106 FAST_FLOAT z5, z10, z11, z12, z13;
nengel@0 107 int16_t *inptr;
nengel@0 108 FAST_FLOAT *quantptr;
nengel@0 109 FAST_FLOAT *wsptr;
nengel@0 110 uint8_t *outptr;
nengel@0 111 int ctr;
nengel@0 112 FAST_FLOAT workspace[DCTSIZE2]; /* buffers data between passes */
nengel@0 113
nengel@0 114 /* Pass 1: process columns from input, store into work array. */
nengel@0 115
nengel@0 116 inptr = compptr->DCT;
nengel@0 117 quantptr = compptr->Q_table;
nengel@0 118 wsptr = workspace;
nengel@0 119 for (ctr = DCTSIZE; ctr > 0; ctr--) {
nengel@0 120 /* Due to quantization, we will usually find that many of the input
nengel@0 121 * coefficients are zero, especially the AC terms. We can exploit this
nengel@0 122 * by short-circuiting the IDCT calculation for any column in which all
nengel@0 123 * the AC terms are zero. In that case each output is equal to the
nengel@0 124 * DC coefficient (with scale factor as needed).
nengel@0 125 * With typical images and quantization tables, half or more of the
nengel@0 126 * column DCT calculations can be simplified this way.
nengel@0 127 */
nengel@0 128
nengel@0 129 if (inptr[DCTSIZE*1] == 0 && inptr[DCTSIZE*2] == 0 &&
nengel@0 130 inptr[DCTSIZE*3] == 0 && inptr[DCTSIZE*4] == 0 &&
nengel@0 131 inptr[DCTSIZE*5] == 0 && inptr[DCTSIZE*6] == 0 &&
nengel@0 132 inptr[DCTSIZE*7] == 0) {
nengel@0 133 /* AC terms all zero */
nengel@0 134 FAST_FLOAT dcval = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
nengel@0 135
nengel@0 136 wsptr[DCTSIZE*0] = dcval;
nengel@0 137 wsptr[DCTSIZE*1] = dcval;
nengel@0 138 wsptr[DCTSIZE*2] = dcval;
nengel@0 139 wsptr[DCTSIZE*3] = dcval;
nengel@0 140 wsptr[DCTSIZE*4] = dcval;
nengel@0 141 wsptr[DCTSIZE*5] = dcval;
nengel@0 142 wsptr[DCTSIZE*6] = dcval;
nengel@0 143 wsptr[DCTSIZE*7] = dcval;
nengel@0 144
nengel@0 145 inptr++; /* advance pointers to next column */
nengel@0 146 quantptr++;
nengel@0 147 wsptr++;
nengel@0 148 continue;
nengel@0 149 }
nengel@0 150
nengel@0 151 /* Even part */
nengel@0 152
nengel@0 153 tmp0 = DEQUANTIZE(inptr[DCTSIZE*0], quantptr[DCTSIZE*0]);
nengel@0 154 tmp1 = DEQUANTIZE(inptr[DCTSIZE*2], quantptr[DCTSIZE*2]);
nengel@0 155 tmp2 = DEQUANTIZE(inptr[DCTSIZE*4], quantptr[DCTSIZE*4]);
nengel@0 156 tmp3 = DEQUANTIZE(inptr[DCTSIZE*6], quantptr[DCTSIZE*6]);
nengel@0 157
nengel@0 158 tmp10 = tmp0 + tmp2; /* phase 3 */
nengel@0 159 tmp11 = tmp0 - tmp2;
nengel@0 160
nengel@0 161 tmp13 = tmp1 + tmp3; /* phases 5-3 */
nengel@0 162 tmp12 = (tmp1 - tmp3) * ((FAST_FLOAT) 1.414213562) - tmp13; /* 2*c4 */
nengel@0 163
nengel@0 164 tmp0 = tmp10 + tmp13; /* phase 2 */
nengel@0 165 tmp3 = tmp10 - tmp13;
nengel@0 166 tmp1 = tmp11 + tmp12;
nengel@0 167 tmp2 = tmp11 - tmp12;
nengel@0 168
nengel@0 169 /* Odd part */
nengel@0 170
nengel@0 171 tmp4 = DEQUANTIZE(inptr[DCTSIZE*1], quantptr[DCTSIZE*1]);
nengel@0 172 tmp5 = DEQUANTIZE(inptr[DCTSIZE*3], quantptr[DCTSIZE*3]);
nengel@0 173 tmp6 = DEQUANTIZE(inptr[DCTSIZE*5], quantptr[DCTSIZE*5]);
nengel@0 174 tmp7 = DEQUANTIZE(inptr[DCTSIZE*7], quantptr[DCTSIZE*7]);
nengel@0 175
nengel@0 176 z13 = tmp6 + tmp5; /* phase 6 */
nengel@0 177 z10 = tmp6 - tmp5;
nengel@0 178 z11 = tmp4 + tmp7;
nengel@0 179 z12 = tmp4 - tmp7;
nengel@0 180
nengel@0 181 tmp7 = z11 + z13; /* phase 5 */
nengel@0 182 tmp11 = (z11 - z13) * ((FAST_FLOAT) 1.414213562); /* 2*c4 */
nengel@0 183
nengel@0 184 z5 = (z10 + z12) * ((FAST_FLOAT) 1.847759065); /* 2*c2 */
nengel@0 185 tmp10 = ((FAST_FLOAT) 1.082392200) * z12 - z5; /* 2*(c2-c6) */
nengel@0 186 tmp12 = ((FAST_FLOAT) -2.613125930) * z10 + z5; /* -2*(c2+c6) */
nengel@0 187
nengel@0 188 tmp6 = tmp12 - tmp7; /* phase 2 */
nengel@0 189 tmp5 = tmp11 - tmp6;
nengel@0 190 tmp4 = tmp10 + tmp5;
nengel@0 191
nengel@0 192 wsptr[DCTSIZE*0] = tmp0 + tmp7;
nengel@0 193 wsptr[DCTSIZE*7] = tmp0 - tmp7;
nengel@0 194 wsptr[DCTSIZE*1] = tmp1 + tmp6;
nengel@0 195 wsptr[DCTSIZE*6] = tmp1 - tmp6;
nengel@0 196 wsptr[DCTSIZE*2] = tmp2 + tmp5;
nengel@0 197 wsptr[DCTSIZE*5] = tmp2 - tmp5;
nengel@0 198 wsptr[DCTSIZE*4] = tmp3 + tmp4;
nengel@0 199 wsptr[DCTSIZE*3] = tmp3 - tmp4;
nengel@0 200
nengel@0 201 inptr++; /* advance pointers to next column */
nengel@0 202 quantptr++;
nengel@0 203 wsptr++;
nengel@0 204 }
nengel@0 205
nengel@0 206 /* Pass 2: process rows from work array, store into output array. */
nengel@0 207 /* Note that we must descale the results by a factor of 8 == 2**3. */
nengel@0 208
nengel@0 209 wsptr = workspace;
nengel@0 210 outptr = output_buf;
nengel@0 211 for (ctr = 0; ctr < DCTSIZE; ctr++) {
nengel@0 212 /* Rows of zeroes can be exploited in the same way as we did with columns.
nengel@0 213 * However, the column calculation has created many nonzero AC terms, so
nengel@0 214 * the simplification applies less often (typically 5% to 10% of the time).
nengel@0 215 * And testing floats for zero is relatively expensive, so we don't bother.
nengel@0 216 */
nengel@0 217
nengel@0 218 /* Even part */
nengel@0 219
nengel@0 220 tmp10 = wsptr[0] + wsptr[4];
nengel@0 221 tmp11 = wsptr[0] - wsptr[4];
nengel@0 222
nengel@0 223 tmp13 = wsptr[2] + wsptr[6];
nengel@0 224 tmp12 = (wsptr[2] - wsptr[6]) * ((FAST_FLOAT) 1.414213562) - tmp13;
nengel@0 225
nengel@0 226 tmp0 = tmp10 + tmp13;
nengel@0 227 tmp3 = tmp10 - tmp13;
nengel@0 228 tmp1 = tmp11 + tmp12;
nengel@0 229 tmp2 = tmp11 - tmp12;
nengel@0 230
nengel@0 231 /* Odd part */
nengel@0 232
nengel@0 233 z13 = wsptr[5] + wsptr[3];
nengel@0 234 z10 = wsptr[5] - wsptr[3];
nengel@0 235 z11 = wsptr[1] + wsptr[7];
nengel@0 236 z12 = wsptr[1] - wsptr[7];
nengel@0 237
nengel@0 238 tmp7 = z11 + z13;
nengel@0 239 tmp11 = (z11 - z13) * ((FAST_FLOAT) 1.414213562);
nengel@0 240
nengel@0 241 z5 = (z10 + z12) * ((FAST_FLOAT) 1.847759065); /* 2*c2 */
nengel@0 242 tmp10 = ((FAST_FLOAT) 1.082392200) * z12 - z5; /* 2*(c2-c6) */
nengel@0 243 tmp12 = ((FAST_FLOAT) -2.613125930) * z10 + z5; /* -2*(c2+c6) */
nengel@0 244
nengel@0 245 tmp6 = tmp12 - tmp7;
nengel@0 246 tmp5 = tmp11 - tmp6;
nengel@0 247 tmp4 = tmp10 + tmp5;
nengel@0 248
nengel@0 249 /* Final output stage: scale down by a factor of 8 and range-limit */
nengel@0 250
nengel@0 251 outptr[0] = descale_and_clamp((int)(tmp0 + tmp7), 3);
nengel@0 252 outptr[7] = descale_and_clamp((int)(tmp0 - tmp7), 3);
nengel@0 253 outptr[1] = descale_and_clamp((int)(tmp1 + tmp6), 3);
nengel@0 254 outptr[6] = descale_and_clamp((int)(tmp1 - tmp6), 3);
nengel@0 255 outptr[2] = descale_and_clamp((int)(tmp2 + tmp5), 3);
nengel@0 256 outptr[5] = descale_and_clamp((int)(tmp2 - tmp5), 3);
nengel@0 257 outptr[4] = descale_and_clamp((int)(tmp3 + tmp4), 3);
nengel@0 258 outptr[3] = descale_and_clamp((int)(tmp3 - tmp4), 3);
nengel@0 259
nengel@0 260
nengel@0 261 wsptr += DCTSIZE; /* advance pointer to next row */
nengel@0 262 outptr += stride;
nengel@0 263 }
nengel@0 264 }
nengel@0 265