Libav
|
00001 /* 00002 * jfdctfst.c 00003 * 00004 * This file is part of the Independent JPEG Group's software. 00005 * 00006 * The authors make NO WARRANTY or representation, either express or implied, 00007 * with respect to this software, its quality, accuracy, merchantability, or 00008 * fitness for a particular purpose. This software is provided "AS IS", and 00009 * you, its user, assume the entire risk as to its quality and accuracy. 00010 * 00011 * This software is copyright (C) 1994-1996, Thomas G. Lane. 00012 * All Rights Reserved except as specified below. 00013 * 00014 * Permission is hereby granted to use, copy, modify, and distribute this 00015 * software (or portions thereof) for any purpose, without fee, subject to 00016 * these conditions: 00017 * (1) If any part of the source code for this software is distributed, then 00018 * this README file must be included, with this copyright and no-warranty 00019 * notice unaltered; and any additions, deletions, or changes to the original 00020 * files must be clearly indicated in accompanying documentation. 00021 * (2) If only executable code is distributed, then the accompanying 00022 * documentation must state that "this software is based in part on the work 00023 * of the Independent JPEG Group". 00024 * (3) Permission for use of this software is granted only if the user accepts 00025 * full responsibility for any undesirable consequences; the authors accept 00026 * NO LIABILITY for damages of any kind. 00027 * 00028 * These conditions apply to any software derived from or based on the IJG 00029 * code, not just to the unmodified library. If you use our work, you ought 00030 * to acknowledge us. 00031 * 00032 * Permission is NOT granted for the use of any IJG author's name or company 00033 * name in advertising or publicity relating to this software or products 00034 * derived from it. This software may be referred to only as "the Independent 00035 * JPEG Group's software". 00036 * 00037 * We specifically permit and encourage the use of this software as the basis 00038 * of commercial products, provided that all warranty or liability claims are 00039 * assumed by the product vendor. 00040 * 00041 * This file contains a fast, not so accurate integer implementation of the 00042 * forward DCT (Discrete Cosine Transform). 00043 * 00044 * A 2-D DCT can be done by 1-D DCT on each row followed by 1-D DCT 00045 * on each column. Direct algorithms are also available, but they are 00046 * much more complex and seem not to be any faster when reduced to code. 00047 * 00048 * This implementation is based on Arai, Agui, and Nakajima's algorithm for 00049 * scaled DCT. Their original paper (Trans. IEICE E-71(11):1095) is in 00050 * Japanese, but the algorithm is described in the Pennebaker & Mitchell 00051 * JPEG textbook (see REFERENCES section in file README). The following code 00052 * is based directly on figure 4-8 in P&M. 00053 * While an 8-point DCT cannot be done in less than 11 multiplies, it is 00054 * possible to arrange the computation so that many of the multiplies are 00055 * simple scalings of the final outputs. These multiplies can then be 00056 * folded into the multiplications or divisions by the JPEG quantization 00057 * table entries. The AA&N method leaves only 5 multiplies and 29 adds 00058 * to be done in the DCT itself. 00059 * The primary disadvantage of this method is that with fixed-point math, 00060 * accuracy is lost due to imprecise representation of the scaled 00061 * quantization values. The smaller the quantization table entry, the less 00062 * precise the scaled value, so this implementation does worse with high- 00063 * quality-setting files than with low-quality ones. 00064 */ 00065 00071 #include <stdlib.h> 00072 #include <stdio.h> 00073 #include "libavutil/common.h" 00074 #include "dsputil.h" 00075 00076 #define DCTSIZE 8 00077 #define GLOBAL(x) x 00078 #define RIGHT_SHIFT(x, n) ((x) >> (n)) 00079 00080 /* 00081 * This module is specialized to the case DCTSIZE = 8. 00082 */ 00083 00084 #if DCTSIZE != 8 00085 Sorry, this code only copes with 8x8 DCTs. /* deliberate syntax err */ 00086 #endif 00087 00088 00089 /* Scaling decisions are generally the same as in the LL&M algorithm; 00090 * see jfdctint.c for more details. However, we choose to descale 00091 * (right shift) multiplication products as soon as they are formed, 00092 * rather than carrying additional fractional bits into subsequent additions. 00093 * This compromises accuracy slightly, but it lets us save a few shifts. 00094 * More importantly, 16-bit arithmetic is then adequate (for 8-bit samples) 00095 * everywhere except in the multiplications proper; this saves a good deal 00096 * of work on 16-bit-int machines. 00097 * 00098 * Again to save a few shifts, the intermediate results between pass 1 and 00099 * pass 2 are not upscaled, but are represented only to integral precision. 00100 * 00101 * A final compromise is to represent the multiplicative constants to only 00102 * 8 fractional bits, rather than 13. This saves some shifting work on some 00103 * machines, and may also reduce the cost of multiplication (since there 00104 * are fewer one-bits in the constants). 00105 */ 00106 00107 #define CONST_BITS 8 00108 00109 00110 /* Some C compilers fail to reduce "FIX(constant)" at compile time, thus 00111 * causing a lot of useless floating-point operations at run time. 00112 * To get around this we use the following pre-calculated constants. 00113 * If you change CONST_BITS you may want to add appropriate values. 00114 * (With a reasonable C compiler, you can just rely on the FIX() macro...) 00115 */ 00116 00117 #if CONST_BITS == 8 00118 #define FIX_0_382683433 ((int32_t) 98) /* FIX(0.382683433) */ 00119 #define FIX_0_541196100 ((int32_t) 139) /* FIX(0.541196100) */ 00120 #define FIX_0_707106781 ((int32_t) 181) /* FIX(0.707106781) */ 00121 #define FIX_1_306562965 ((int32_t) 334) /* FIX(1.306562965) */ 00122 #else 00123 #define FIX_0_382683433 FIX(0.382683433) 00124 #define FIX_0_541196100 FIX(0.541196100) 00125 #define FIX_0_707106781 FIX(0.707106781) 00126 #define FIX_1_306562965 FIX(1.306562965) 00127 #endif 00128 00129 00130 /* We can gain a little more speed, with a further compromise in accuracy, 00131 * by omitting the addition in a descaling shift. This yields an incorrectly 00132 * rounded result half the time... 00133 */ 00134 00135 #ifndef USE_ACCURATE_ROUNDING 00136 #undef DESCALE 00137 #define DESCALE(x,n) RIGHT_SHIFT(x, n) 00138 #endif 00139 00140 00141 /* Multiply a DCTELEM variable by an int32_t constant, and immediately 00142 * descale to yield a DCTELEM result. 00143 */ 00144 00145 #define MULTIPLY(var,const) ((DCTELEM) DESCALE((var) * (const), CONST_BITS)) 00146 00147 static av_always_inline void row_fdct(DCTELEM * data){ 00148 int_fast16_t tmp0, tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7; 00149 int_fast16_t tmp10, tmp11, tmp12, tmp13; 00150 int_fast16_t z1, z2, z3, z4, z5, z11, z13; 00151 DCTELEM *dataptr; 00152 int ctr; 00153 00154 /* Pass 1: process rows. */ 00155 00156 dataptr = data; 00157 for (ctr = DCTSIZE-1; ctr >= 0; ctr--) { 00158 tmp0 = dataptr[0] + dataptr[7]; 00159 tmp7 = dataptr[0] - dataptr[7]; 00160 tmp1 = dataptr[1] + dataptr[6]; 00161 tmp6 = dataptr[1] - dataptr[6]; 00162 tmp2 = dataptr[2] + dataptr[5]; 00163 tmp5 = dataptr[2] - dataptr[5]; 00164 tmp3 = dataptr[3] + dataptr[4]; 00165 tmp4 = dataptr[3] - dataptr[4]; 00166 00167 /* Even part */ 00168 00169 tmp10 = tmp0 + tmp3; /* phase 2 */ 00170 tmp13 = tmp0 - tmp3; 00171 tmp11 = tmp1 + tmp2; 00172 tmp12 = tmp1 - tmp2; 00173 00174 dataptr[0] = tmp10 + tmp11; /* phase 3 */ 00175 dataptr[4] = tmp10 - tmp11; 00176 00177 z1 = MULTIPLY(tmp12 + tmp13, FIX_0_707106781); /* c4 */ 00178 dataptr[2] = tmp13 + z1; /* phase 5 */ 00179 dataptr[6] = tmp13 - z1; 00180 00181 /* Odd part */ 00182 00183 tmp10 = tmp4 + tmp5; /* phase 2 */ 00184 tmp11 = tmp5 + tmp6; 00185 tmp12 = tmp6 + tmp7; 00186 00187 /* The rotator is modified from fig 4-8 to avoid extra negations. */ 00188 z5 = MULTIPLY(tmp10 - tmp12, FIX_0_382683433); /* c6 */ 00189 z2 = MULTIPLY(tmp10, FIX_0_541196100) + z5; /* c2-c6 */ 00190 z4 = MULTIPLY(tmp12, FIX_1_306562965) + z5; /* c2+c6 */ 00191 z3 = MULTIPLY(tmp11, FIX_0_707106781); /* c4 */ 00192 00193 z11 = tmp7 + z3; /* phase 5 */ 00194 z13 = tmp7 - z3; 00195 00196 dataptr[5] = z13 + z2; /* phase 6 */ 00197 dataptr[3] = z13 - z2; 00198 dataptr[1] = z11 + z4; 00199 dataptr[7] = z11 - z4; 00200 00201 dataptr += DCTSIZE; /* advance pointer to next row */ 00202 } 00203 } 00204 00205 /* 00206 * Perform the forward DCT on one block of samples. 00207 */ 00208 00209 GLOBAL(void) 00210 fdct_ifast (DCTELEM * data) 00211 { 00212 int_fast16_t tmp0, tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7; 00213 int_fast16_t tmp10, tmp11, tmp12, tmp13; 00214 int_fast16_t z1, z2, z3, z4, z5, z11, z13; 00215 DCTELEM *dataptr; 00216 int ctr; 00217 00218 row_fdct(data); 00219 00220 /* Pass 2: process columns. */ 00221 00222 dataptr = data; 00223 for (ctr = DCTSIZE-1; ctr >= 0; ctr--) { 00224 tmp0 = dataptr[DCTSIZE*0] + dataptr[DCTSIZE*7]; 00225 tmp7 = dataptr[DCTSIZE*0] - dataptr[DCTSIZE*7]; 00226 tmp1 = dataptr[DCTSIZE*1] + dataptr[DCTSIZE*6]; 00227 tmp6 = dataptr[DCTSIZE*1] - dataptr[DCTSIZE*6]; 00228 tmp2 = dataptr[DCTSIZE*2] + dataptr[DCTSIZE*5]; 00229 tmp5 = dataptr[DCTSIZE*2] - dataptr[DCTSIZE*5]; 00230 tmp3 = dataptr[DCTSIZE*3] + dataptr[DCTSIZE*4]; 00231 tmp4 = dataptr[DCTSIZE*3] - dataptr[DCTSIZE*4]; 00232 00233 /* Even part */ 00234 00235 tmp10 = tmp0 + tmp3; /* phase 2 */ 00236 tmp13 = tmp0 - tmp3; 00237 tmp11 = tmp1 + tmp2; 00238 tmp12 = tmp1 - tmp2; 00239 00240 dataptr[DCTSIZE*0] = tmp10 + tmp11; /* phase 3 */ 00241 dataptr[DCTSIZE*4] = tmp10 - tmp11; 00242 00243 z1 = MULTIPLY(tmp12 + tmp13, FIX_0_707106781); /* c4 */ 00244 dataptr[DCTSIZE*2] = tmp13 + z1; /* phase 5 */ 00245 dataptr[DCTSIZE*6] = tmp13 - z1; 00246 00247 /* Odd part */ 00248 00249 tmp10 = tmp4 + tmp5; /* phase 2 */ 00250 tmp11 = tmp5 + tmp6; 00251 tmp12 = tmp6 + tmp7; 00252 00253 /* The rotator is modified from fig 4-8 to avoid extra negations. */ 00254 z5 = MULTIPLY(tmp10 - tmp12, FIX_0_382683433); /* c6 */ 00255 z2 = MULTIPLY(tmp10, FIX_0_541196100) + z5; /* c2-c6 */ 00256 z4 = MULTIPLY(tmp12, FIX_1_306562965) + z5; /* c2+c6 */ 00257 z3 = MULTIPLY(tmp11, FIX_0_707106781); /* c4 */ 00258 00259 z11 = tmp7 + z3; /* phase 5 */ 00260 z13 = tmp7 - z3; 00261 00262 dataptr[DCTSIZE*5] = z13 + z2; /* phase 6 */ 00263 dataptr[DCTSIZE*3] = z13 - z2; 00264 dataptr[DCTSIZE*1] = z11 + z4; 00265 dataptr[DCTSIZE*7] = z11 - z4; 00266 00267 dataptr++; /* advance pointer to next column */ 00268 } 00269 } 00270 00271 /* 00272 * Perform the forward 2-4-8 DCT on one block of samples. 00273 */ 00274 00275 GLOBAL(void) 00276 fdct_ifast248 (DCTELEM * data) 00277 { 00278 int_fast16_t tmp0, tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7; 00279 int_fast16_t tmp10, tmp11, tmp12, tmp13; 00280 int_fast16_t z1; 00281 DCTELEM *dataptr; 00282 int ctr; 00283 00284 row_fdct(data); 00285 00286 /* Pass 2: process columns. */ 00287 00288 dataptr = data; 00289 for (ctr = DCTSIZE-1; ctr >= 0; ctr--) { 00290 tmp0 = dataptr[DCTSIZE*0] + dataptr[DCTSIZE*1]; 00291 tmp1 = dataptr[DCTSIZE*2] + dataptr[DCTSIZE*3]; 00292 tmp2 = dataptr[DCTSIZE*4] + dataptr[DCTSIZE*5]; 00293 tmp3 = dataptr[DCTSIZE*6] + dataptr[DCTSIZE*7]; 00294 tmp4 = dataptr[DCTSIZE*0] - dataptr[DCTSIZE*1]; 00295 tmp5 = dataptr[DCTSIZE*2] - dataptr[DCTSIZE*3]; 00296 tmp6 = dataptr[DCTSIZE*4] - dataptr[DCTSIZE*5]; 00297 tmp7 = dataptr[DCTSIZE*6] - dataptr[DCTSIZE*7]; 00298 00299 /* Even part */ 00300 00301 tmp10 = tmp0 + tmp3; 00302 tmp11 = tmp1 + tmp2; 00303 tmp12 = tmp1 - tmp2; 00304 tmp13 = tmp0 - tmp3; 00305 00306 dataptr[DCTSIZE*0] = tmp10 + tmp11; 00307 dataptr[DCTSIZE*4] = tmp10 - tmp11; 00308 00309 z1 = MULTIPLY(tmp12 + tmp13, FIX_0_707106781); 00310 dataptr[DCTSIZE*2] = tmp13 + z1; 00311 dataptr[DCTSIZE*6] = tmp13 - z1; 00312 00313 tmp10 = tmp4 + tmp7; 00314 tmp11 = tmp5 + tmp6; 00315 tmp12 = tmp5 - tmp6; 00316 tmp13 = tmp4 - tmp7; 00317 00318 dataptr[DCTSIZE*1] = tmp10 + tmp11; 00319 dataptr[DCTSIZE*5] = tmp10 - tmp11; 00320 00321 z1 = MULTIPLY(tmp12 + tmp13, FIX_0_707106781); 00322 dataptr[DCTSIZE*3] = tmp13 + z1; 00323 dataptr[DCTSIZE*7] = tmp13 - z1; 00324 00325 dataptr++; /* advance pointer to next column */ 00326 } 00327 } 00328 00329 00330 #undef GLOBAL 00331 #undef CONST_BITS 00332 #undef DESCALE 00333 #undef FIX_0_541196100 00334 #undef FIX_1_306562965