src/SDL2/external/libwebp-0.3.0/src/utils/quant_levels.c

   1 // Copyright 2011 Google Inc. All Rights Reserved.
   2 //
   3 // This code is licensed under the same terms as WebM:
   4 //  Software License Agreement:  http://www.webmproject.org/license/software/
   5 //  Additional IP Rights Grant:  http://www.webmproject.org/license/additional/
   6 // -----------------------------------------------------------------------------
   7 //
   8 // Quantize levels for specified number of quantization-levels ([2, 256]).
   9 // Min and max values are preserved (usual 0 and 255 for alpha plane).
  10 //
  11 // Author: Skal (pascal.massimino@gmail.com)
  12
  13 #include <assert.h>
  14
  15 #include "./quant_levels.h"
  16
  17 #if defined(__cplusplus) || defined(c_plusplus)
  18 extern "C" {
  19 #endif
  20
  21 #define NUM_SYMBOLS     256
  22
  23 #define MAX_ITER  6             // Maximum number of convergence steps.
  24 #define ERROR_THRESHOLD 1e-4    // MSE stopping criterion.
  25
  26 // -----------------------------------------------------------------------------
  27 // Quantize levels.
  28
  29 int QuantizeLevels(uint8_t* const data, int width, int height,
  30                    int num_levels, uint64_t* const sse) {
  31   int freq[NUM_SYMBOLS] = { 0 };
  32   int q_level[NUM_SYMBOLS] = { 0 };
  33   double inv_q_level[NUM_SYMBOLS] = { 0 };
  34   int min_s = 255, max_s = 0;
  35   const size_t data_size = height * width;
  36   int i, num_levels_in, iter;
  37   double last_err = 1.e38, err = 0.;
  38   const double err_threshold = ERROR_THRESHOLD * data_size;
  39
  40   if (data == NULL) {
  41     return 0;
  42   }
  43
  44   if (width <= 0 || height <= 0) {
  45     return 0;
  46   }
  47
  48   if (num_levels < 2 || num_levels > 256) {
  49     return 0;
  50   }
  51
  52   {
  53     size_t n;
  54     num_levels_in = 0;
  55     for (n = 0; n < data_size; ++n) {
  56       num_levels_in += (freq[data[n]] == 0);
  57       if (min_s > data[n]) min_s = data[n];
  58       if (max_s < data[n]) max_s = data[n];
  59       ++freq[data[n]];
  60     }
  61   }
  62
  63   if (num_levels_in <= num_levels) goto End;  // nothing to do!
  64
  65   // Start with uniformly spread centroids.
  66   for (i = 0; i < num_levels; ++i) {
  67     inv_q_level[i] = min_s + (double)(max_s - min_s) * i / (num_levels - 1);
  68   }
  69
  70   // Fixed values. Won't be changed.
  71   q_level[min_s] = 0;
  72   q_level[max_s] = num_levels - 1;
  73   assert(inv_q_level[0] == min_s);
  74   assert(inv_q_level[num_levels - 1] == max_s);
  75
  76   // k-Means iterations.
  77   for (iter = 0; iter < MAX_ITER; ++iter) {
  78     double q_sum[NUM_SYMBOLS] = { 0 };
  79     double q_count[NUM_SYMBOLS] = { 0 };
  80     int s, slot = 0;
  81
  82     // Assign classes to representatives.
  83     for (s = min_s; s <= max_s; ++s) {
  84       // Keep track of the nearest neighbour 'slot'
  85       while (slot < num_levels - 1 &&
  86              2 * s > inv_q_level[slot] + inv_q_level[slot + 1]) {
  87         ++slot;
  88       }
  89       if (freq[s] > 0) {
  90         q_sum[slot] += s * freq[s];
  91         q_count[slot] += freq[s];
  92       }
  93       q_level[s] = slot;
  94     }
  95
  96     // Assign new representatives to classes.
  97     if (num_levels > 2) {
  98       for (slot = 1; slot < num_levels - 1; ++slot) {
  99         const double count = q_count[slot];
 100         if (count > 0.) {
 101           inv_q_level[slot] = q_sum[slot] / count;
 102         }
 103       }
 104     }
 105
 106     // Compute convergence error.
 107     err = 0.;
 108     for (s = min_s; s <= max_s; ++s) {
 109       const double error = s - inv_q_level[q_level[s]];
 110       err += freq[s] * error * error;
 111     }
 112
 113     // Check for convergence: we stop as soon as the error is no
 114     // longer improving.
 115     if (last_err - err < err_threshold) break;
 116     last_err = err;
 117   }
 118
 119   // Remap the alpha plane to quantized values.
 120   {
 121     // double->int rounding operation can be costly, so we do it
 122     // once for all before remapping. We also perform the data[] -> slot
 123     // mapping, while at it (avoid one indirection in the final loop).
 124     uint8_t map[NUM_SYMBOLS];
 125     int s;
 126     size_t n;
 127     for (s = min_s; s <= max_s; ++s) {
 128       const int slot = q_level[s];
 129       map[s] = (uint8_t)(inv_q_level[slot] + .5);
 130     }
 131     // Final pass.
 132     for (n = 0; n < data_size; ++n) {
 133       data[n] = map[data[n]];
 134     }
 135   }
 136  End:
 137   // Store sum of squared error if needed.
 138   if (sse != NULL) *sse = (uint64_t)err;
 139
 140   return 1;
 141 }
 142
 143 #if defined(__cplusplus) || defined(c_plusplus)
 144 }    // extern "C"
 145 #endif