libs/jpeg/jcdctmgr.c

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