Dependencies/libjpeg/jcdctmgr.c

   1 /*
   2  * jcdctmgr.c
   3  *
   4  * Copyright (C) 1994-1996, Thomas G. Lane.
   5  * This file is part of the Independent JPEG Group's software.
   6  * For conditions of distribution and use, see the accompanying README file.
   7  *
   8  * This file contains the forward-DCT management logic.
   9  * This code selects a particular DCT implementation to be used,
  10  * and it performs related housekeeping chores including coefficient
  11  * quantization.
  12  */
  13
  14 #define JPEG_INTERNALS
  15 #include "jinclude.h"
  16 #include "jpeglib.h"
  17 #include "jdct.h"               /* Private declarations for DCT subsystem */
  18
  19
  20 /* Private subobject for this module */
  21
  22 typedef struct {
  23   struct jpeg_forward_dct pub;  /* public fields */
  24
  25   /* Pointer to the DCT routine actually in use */
  26   forward_DCT_method_ptr do_dct[MAX_COMPONENTS];
  27
  28   /* The actual post-DCT divisors --- not identical to the quant table
  29    * entries, because of scaling (especially for an unnormalized DCT).
  30    * Each table is given in normal array order.
  31    */
  32   DCTELEM * divisors[NUM_QUANT_TBLS];
  33
  34 #ifdef DCT_FLOAT_SUPPORTED
  35   /* Same as above for the floating-point case. */
  36   float_DCT_method_ptr do_float_dct[MAX_COMPONENTS];
  37   FAST_FLOAT * float_divisors[NUM_QUANT_TBLS];
  38 #endif
  39 } my_fdct_controller;
  40
  41 typedef my_fdct_controller * my_fdct_ptr;
  42
  43
  44 /* The current scaled-DCT routines require ISLOW-style divisor tables,
  45  * so be sure to compile that code if either ISLOW or SCALING is requested.
  46  */
  47 #ifdef DCT_ISLOW_SUPPORTED
  48 #define PROVIDE_ISLOW_TABLES
  49 #else
  50 #ifdef DCT_SCALING_SUPPORTED
  51 #define PROVIDE_ISLOW_TABLES
  52 #endif
  53 #endif
  54
  55
  56 /*
  57  * Perform forward DCT on one or more blocks of a component.
  58  *
  59  * The input samples are taken from the sample_data[] array starting at
  60  * position start_row/start_col, and moving to the right for any additional
  61  * blocks. The quantized coefficients are returned in coef_blocks[].
  62  */
  63
  64 METHODDEF(void)
  65 forward_DCT (j_compress_ptr cinfo, jpeg_component_info * compptr,
  66              JSAMPARRAY sample_data, JBLOCKROW coef_blocks,
  67              JDIMENSION start_row, JDIMENSION start_col,
  68              JDIMENSION num_blocks)
  69 /* This version is used for integer DCT implementations. */
  70 {
  71   /* This routine is heavily used, so it's worth coding it tightly. */
  72   my_fdct_ptr fdct = (my_fdct_ptr) cinfo->fdct;
  73   forward_DCT_method_ptr do_dct = fdct->do_dct[compptr->component_index];
  74   DCTELEM * divisors = fdct->divisors[compptr->quant_tbl_no];
  75   DCTELEM workspace[DCTSIZE2];  /* work area for FDCT subroutine */
  76   JDIMENSION bi;
  77
  78   sample_data += start_row;     /* fold in the vertical offset once */
  79
  80   for (bi = 0; bi < num_blocks; bi++, start_col += compptr->DCT_h_scaled_size) {
  81     /* Perform the DCT */
  82     (*do_dct) (workspace, sample_data, start_col);
  83
  84     /* Quantize/descale the coefficients, and store into coef_blocks[] */
  85     { register DCTELEM temp, qval;
  86       register int i;
  87       register JCOEFPTR output_ptr = coef_blocks[bi];
  88
  89       for (i = 0; i < DCTSIZE2; i++) {
  90         qval = divisors[i];
  91         temp = workspace[i];
  92         /* Divide the coefficient value by qval, ensuring proper rounding.
  93          * Since C does not specify the direction of rounding for negative
  94          * quotients, we have to force the dividend positive for portability.
  95          *
  96          * In most files, at least half of the output values will be zero
  97          * (at default quantization settings, more like three-quarters...)
  98          * so we should ensure that this case is fast.  On many machines,
  99          * a comparison is enough cheaper than a divide to make a special test
 100          * a win.  Since both inputs will be nonnegative, we need only test
 101          * for a < b to discover whether a/b is 0.
 102          * If your machine's division is fast enough, define FAST_DIVIDE.
 103          */
 104 #ifdef FAST_DIVIDE
 105 #define DIVIDE_BY(a,b)  a /= b
 106 #else
 107 #define DIVIDE_BY(a,b)  if (a >= b) a /= b; else a = 0
 108 #endif
 109         if (temp < 0) {
 110           temp = -temp;
 111           temp += qval>>1;      /* for rounding */
 112           DIVIDE_BY(temp, qval);
 113           temp = -temp;
 114         } else {
 115           temp += qval>>1;      /* for rounding */
 116           DIVIDE_BY(temp, qval);
 117         }
 118         output_ptr[i] = (JCOEF) temp;
 119       }
 120     }
 121   }
 122 }
 123
 124
 125 #ifdef DCT_FLOAT_SUPPORTED
 126
 127 METHODDEF(void)
 128 forward_DCT_float (j_compress_ptr cinfo, jpeg_component_info * compptr,
 129                    JSAMPARRAY sample_data, JBLOCKROW coef_blocks,
 130                    JDIMENSION start_row, JDIMENSION start_col,
 131                    JDIMENSION num_blocks)
 132 /* This version is used for floating-point DCT implementations. */
 133 {
 134   /* This routine is heavily used, so it's worth coding it tightly. */
 135   my_fdct_ptr fdct = (my_fdct_ptr) cinfo->fdct;
 136   float_DCT_method_ptr do_dct = fdct->do_float_dct[compptr->component_index];
 137   FAST_FLOAT * divisors = fdct->float_divisors[compptr->quant_tbl_no];
 138   FAST_FLOAT workspace[DCTSIZE2]; /* work area for FDCT subroutine */
 139   JDIMENSION bi;
 140
 141   sample_data += start_row;     /* fold in the vertical offset once */
 142
 143   for (bi = 0; bi < num_blocks; bi++, start_col += compptr->DCT_h_scaled_size) {
 144     /* Perform the DCT */
 145     (*do_dct) (workspace, sample_data, start_col);
 146
 147     /* Quantize/descale the coefficients, and store into coef_blocks[] */
 148     { register FAST_FLOAT temp;
 149       register int i;
 150       register JCOEFPTR output_ptr = coef_blocks[bi];
 151
 152       for (i = 0; i < DCTSIZE2; i++) {
 153         /* Apply the quantization and scaling factor */
 154         temp = workspace[i] * divisors[i];
 155         /* Round to nearest integer.
 156          * Since C does not specify the direction of rounding for negative
 157          * quotients, we have to force the dividend positive for portability.
 158          * The maximum coefficient size is +-16K (for 12-bit data), so this
 159          * code should work for either 16-bit or 32-bit ints.
 160          */
 161         output_ptr[i] = (JCOEF) ((int) (temp + (FAST_FLOAT) 16384.5) - 16384);
 162       }
 163     }
 164   }
 165 }
 166
 167 #endif /* DCT_FLOAT_SUPPORTED */
 168
 169
 170 /*
 171  * Initialize for a processing pass.
 172  * Verify that all referenced Q-tables are present, and set up
 173  * the divisor table for each one.
 174  * In the current implementation, DCT of all components is done during
 175  * the first pass, even if only some components will be output in the
 176  * first scan.  Hence all components should be examined here.
 177  */
 178
 179 METHODDEF(void)
 180 start_pass_fdctmgr (j_compress_ptr cinfo)
 181 {
 182   my_fdct_ptr fdct = (my_fdct_ptr) cinfo->fdct;
 183   int ci, qtblno, i;
 184   jpeg_component_info *compptr;
 185   int method = 0;
 186   JQUANT_TBL * qtbl;
 187   DCTELEM * dtbl;
 188
 189   for (ci = 0, compptr = cinfo->comp_info; ci < cinfo->num_components;
 190        ci++, compptr++) {
 191     /* Select the proper DCT routine for this component's scaling */
 192     switch ((compptr->DCT_h_scaled_size << 8) + compptr->DCT_v_scaled_size) {
 193 #ifdef DCT_SCALING_SUPPORTED
 194     case ((1 << 8) + 1):
 195       fdct->do_dct[ci] = jpeg_fdct_1x1;
 196       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 197       break;
 198     case ((2 << 8) + 2):
 199       fdct->do_dct[ci] = jpeg_fdct_2x2;
 200       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 201       break;
 202     case ((3 << 8) + 3):
 203       fdct->do_dct[ci] = jpeg_fdct_3x3;
 204       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 205       break;
 206     case ((4 << 8) + 4):
 207       fdct->do_dct[ci] = jpeg_fdct_4x4;
 208       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 209       break;
 210     case ((5 << 8) + 5):
 211       fdct->do_dct[ci] = jpeg_fdct_5x5;
 212       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 213       break;
 214     case ((6 << 8) + 6):
 215       fdct->do_dct[ci] = jpeg_fdct_6x6;
 216       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 217       break;
 218     case ((7 << 8) + 7):
 219       fdct->do_dct[ci] = jpeg_fdct_7x7;
 220       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 221       break;
 222     case ((9 << 8) + 9):
 223       fdct->do_dct[ci] = jpeg_fdct_9x9;
 224       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 225       break;
 226     case ((10 << 8) + 10):
 227       fdct->do_dct[ci] = jpeg_fdct_10x10;
 228       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 229       break;
 230     case ((11 << 8) + 11):
 231       fdct->do_dct[ci] = jpeg_fdct_11x11;
 232       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 233       break;
 234     case ((12 << 8) + 12):
 235       fdct->do_dct[ci] = jpeg_fdct_12x12;
 236       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 237       break;
 238     case ((13 << 8) + 13):
 239       fdct->do_dct[ci] = jpeg_fdct_13x13;
 240       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 241       break;
 242     case ((14 << 8) + 14):
 243       fdct->do_dct[ci] = jpeg_fdct_14x14;
 244       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 245       break;
 246     case ((15 << 8) + 15):
 247       fdct->do_dct[ci] = jpeg_fdct_15x15;
 248       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 249       break;
 250     case ((16 << 8) + 16):
 251       fdct->do_dct[ci] = jpeg_fdct_16x16;
 252       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 253       break;
 254     case ((16 << 8) + 8):
 255       fdct->do_dct[ci] = jpeg_fdct_16x8;
 256       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 257       break;
 258     case ((14 << 8) + 7):
 259       fdct->do_dct[ci] = jpeg_fdct_14x7;
 260       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 261       break;
 262     case ((12 << 8) + 6):
 263       fdct->do_dct[ci] = jpeg_fdct_12x6;
 264       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 265       break;
 266     case ((10 << 8) + 5):
 267       fdct->do_dct[ci] = jpeg_fdct_10x5;
 268       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 269       break;
 270     case ((8 << 8) + 4):
 271       fdct->do_dct[ci] = jpeg_fdct_8x4;
 272       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 273       break;
 274     case ((6 << 8) + 3):
 275       fdct->do_dct[ci] = jpeg_fdct_6x3;
 276       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 277       break;
 278     case ((4 << 8) + 2):
 279       fdct->do_dct[ci] = jpeg_fdct_4x2;
 280       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 281       break;
 282     case ((2 << 8) + 1):
 283       fdct->do_dct[ci] = jpeg_fdct_2x1;
 284       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 285       break;
 286     case ((8 << 8) + 16):
 287       fdct->do_dct[ci] = jpeg_fdct_8x16;
 288       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 289       break;
 290     case ((7 << 8) + 14):
 291       fdct->do_dct[ci] = jpeg_fdct_7x14;
 292       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 293       break;
 294     case ((6 << 8) + 12):
 295       fdct->do_dct[ci] = jpeg_fdct_6x12;
 296       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 297       break;
 298     case ((5 << 8) + 10):
 299       fdct->do_dct[ci] = jpeg_fdct_5x10;
 300       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 301       break;
 302     case ((4 << 8) + 8):
 303       fdct->do_dct[ci] = jpeg_fdct_4x8;
 304       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 305       break;
 306     case ((3 << 8) + 6):
 307       fdct->do_dct[ci] = jpeg_fdct_3x6;
 308       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 309       break;
 310     case ((2 << 8) + 4):
 311       fdct->do_dct[ci] = jpeg_fdct_2x4;
 312       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 313       break;
 314     case ((1 << 8) + 2):
 315       fdct->do_dct[ci] = jpeg_fdct_1x2;
 316       method = JDCT_ISLOW;      /* jfdctint uses islow-style table */
 317       break;
 318 #endif
 319     case ((DCTSIZE << 8) + DCTSIZE):
 320       switch (cinfo->dct_method) {
 321 #ifdef DCT_ISLOW_SUPPORTED
 322       case JDCT_ISLOW:
 323         fdct->do_dct[ci] = jpeg_fdct_islow;
 324         method = JDCT_ISLOW;
 325         break;
 326 #endif
 327 #ifdef DCT_IFAST_SUPPORTED
 328       case JDCT_IFAST:
 329         fdct->do_dct[ci] = jpeg_fdct_ifast;
 330         method = JDCT_IFAST;
 331         break;
 332 #endif
 333 #ifdef DCT_FLOAT_SUPPORTED
 334       case JDCT_FLOAT:
 335         fdct->do_float_dct[ci] = jpeg_fdct_float;
 336         method = JDCT_FLOAT;
 337         break;
 338 #endif
 339       default:
 340         ERREXIT(cinfo, JERR_NOT_COMPILED);
 341         break;
 342       }
 343       break;
 344     default:
 345       ERREXIT2(cinfo, JERR_BAD_DCTSIZE,
 346                compptr->DCT_h_scaled_size, compptr->DCT_v_scaled_size);
 347       break;
 348     }
 349     qtblno = compptr->quant_tbl_no;
 350     /* Make sure specified quantization table is present */
 351     if (qtblno < 0 || qtblno >= NUM_QUANT_TBLS ||
 352         cinfo->quant_tbl_ptrs[qtblno] == NULL)
 353       ERREXIT1(cinfo, JERR_NO_QUANT_TABLE, qtblno);
 354     qtbl = cinfo->quant_tbl_ptrs[qtblno];
 355     /* Compute divisors for this quant table */
 356     /* We may do this more than once for same table, but it's not a big deal */
 357     switch (method) {
 358 #ifdef PROVIDE_ISLOW_TABLES
 359     case JDCT_ISLOW:
 360       /* For LL&M IDCT method, divisors are equal to raw quantization
 361        * coefficients multiplied by 8 (to counteract scaling).
 362        */
 363       if (fdct->divisors[qtblno] == NULL) {
 364         fdct->divisors[qtblno] = (DCTELEM *)
 365           (*cinfo->mem->alloc_small) ((j_common_ptr) cinfo, JPOOL_IMAGE,
 366                                       DCTSIZE2 * SIZEOF(DCTELEM));
 367       }
 368       dtbl = fdct->divisors[qtblno];
 369       for (i = 0; i < DCTSIZE2; i++) {
 370         dtbl[i] = ((DCTELEM) qtbl->quantval[i]) << 3;
 371       }
 372       fdct->pub.forward_DCT[ci] = forward_DCT;
 373       break;
 374 #endif
 375 #ifdef DCT_IFAST_SUPPORTED
 376     case JDCT_IFAST:
 377       {
 378         /* For AA&N IDCT method, divisors are equal to quantization
 379          * coefficients scaled by scalefactor[row]*scalefactor[col], where
 380          *   scalefactor[0] = 1
 381          *   scalefactor[k] = cos(k*PI/16) * sqrt(2)    for k=1..7
 382          * We apply a further scale factor of 8.
 383          */
 384 #define CONST_BITS 14
 385         static const INT16 aanscales[DCTSIZE2] = {
 386           /* precomputed values scaled up by 14 bits */
 387           16384, 22725, 21407, 19266, 16384, 12873,  8867,  4520,
 388           22725, 31521, 29692, 26722, 22725, 17855, 12299,  6270,
 389           21407, 29692, 27969, 25172, 21407, 16819, 11585,  5906,
 390           19266, 26722, 25172, 22654, 19266, 15137, 10426,  5315,
 391           16384, 22725, 21407, 19266, 16384, 12873,  8867,  4520,
 392           12873, 17855, 16819, 15137, 12873, 10114,  6967,  3552,
 393            8867, 12299, 11585, 10426,  8867,  6967,  4799,  2446,
 394            4520,  6270,  5906,  5315,  4520,  3552,  2446,  1247
 395         };
 396         SHIFT_TEMPS
 397
 398         if (fdct->divisors[qtblno] == NULL) {
 399           fdct->divisors[qtblno] = (DCTELEM *)
 400             (*cinfo->mem->alloc_small) ((j_common_ptr) cinfo, JPOOL_IMAGE,
 401                                         DCTSIZE2 * SIZEOF(DCTELEM));
 402         }
 403         dtbl = fdct->divisors[qtblno];
 404         for (i = 0; i < DCTSIZE2; i++) {
 405           dtbl[i] = (DCTELEM)
 406             DESCALE(MULTIPLY16V16((INT32) qtbl->quantval[i],
 407                                   (INT32) aanscales[i]),
 408                     CONST_BITS-3);
 409         }
 410       }
 411       fdct->pub.forward_DCT[ci] = forward_DCT;
 412       break;
 413 #endif
 414 #ifdef DCT_FLOAT_SUPPORTED
 415     case JDCT_FLOAT:
 416       {
 417         /* For float AA&N IDCT method, divisors are equal to quantization
 418          * coefficients scaled by scalefactor[row]*scalefactor[col], where
 419          *   scalefactor[0] = 1
 420          *   scalefactor[k] = cos(k*PI/16) * sqrt(2)    for k=1..7
 421          * We apply a further scale factor of 8.
 422          * What's actually stored is 1/divisor so that the inner loop can
 423          * use a multiplication rather than a division.
 424          */
 425         FAST_FLOAT * fdtbl;
 426         int row, col;
 427         static const double aanscalefactor[DCTSIZE] = {
 428           1.0, 1.387039845, 1.306562965, 1.175875602,
 429           1.0, 0.785694958, 0.541196100, 0.275899379
 430         };
 431
 432         if (fdct->float_divisors[qtblno] == NULL) {
 433           fdct->float_divisors[qtblno] = (FAST_FLOAT *)
 434             (*cinfo->mem->alloc_small) ((j_common_ptr) cinfo, JPOOL_IMAGE,
 435                                         DCTSIZE2 * SIZEOF(FAST_FLOAT));
 436         }
 437         fdtbl = fdct->float_divisors[qtblno];
 438         i = 0;
 439         for (row = 0; row < DCTSIZE; row++) {
 440           for (col = 0; col < DCTSIZE; col++) {
 441             fdtbl[i] = (FAST_FLOAT)
 442               (1.0 / (((double) qtbl->quantval[i] *
 443                        aanscalefactor[row] * aanscalefactor[col] * 8.0)));
 444             i++;
 445           }
 446         }
 447       }
 448       fdct->pub.forward_DCT[ci] = forward_DCT_float;
 449       break;
 450 #endif
 451     default:
 452       ERREXIT(cinfo, JERR_NOT_COMPILED);
 453       break;
 454     }
 455   }
 456 }
 457
 458
 459 /*
 460  * Initialize FDCT manager.
 461  */
 462
 463 GLOBAL(void)
 464 jinit_forward_dct (j_compress_ptr cinfo)
 465 {
 466   my_fdct_ptr fdct;
 467   int i;
 468
 469   fdct = (my_fdct_ptr)
 470     (*cinfo->mem->alloc_small) ((j_common_ptr) cinfo, JPOOL_IMAGE,
 471                                 SIZEOF(my_fdct_controller));
 472   cinfo->fdct = (struct jpeg_forward_dct *) fdct;
 473   fdct->pub.start_pass = start_pass_fdctmgr;
 474
 475   /* Mark divisor tables unallocated */
 476   for (i = 0; i < NUM_QUANT_TBLS; i++) {
 477     fdct->divisors[i] = NULL;
 478 #ifdef DCT_FLOAT_SUPPORTED
 479     fdct->float_divisors[i] = NULL;
 480 #endif
 481   }
 482 }