block_encoder.cpp.bak

这是我刚刚完成的关于JPEG2000的C语言实现的部分程序。小波变换是采用97变换

/******************************************************************************
Description:
   Implements the embedded block coding algorithm, including distortion
estimation and R-D covex hull analysis, in addition to the coding
passes themselves.  The low level services offered by the MQ arithmetic coder
appear in "mq_encoder.cpp" and "mq_encoder.h".
******************************************************************************/

#include <math.h>
#include <string.h>
#include "kdu_messaging.h"
#include "kdu_block_coding.h"
#include "block_coding_common.h"
#include "mq_encoder.h"

static kdu_byte *significance_luts[4] =
  {lh_sig_lut, hl_sig_lut, lh_sig_lut, hh_sig_lut};

#define DISTORTION_LSBS 5
#define SIGNIFICANCE_DISTORTIONS (1<<DISTORTION_LSBS)
#define REFINEMENT_DISTORTIONS (1<<(DISTORTION_LSBS+1))

static kdu_int32 significance_distortion_lut[SIGNIFICANCE_DISTORTIONS];
static kdu_int32 significance_distortion_lut_lossless[SIGNIFICANCE_DISTORTIONS];
static kdu_int32 refinement_distortion_lut[REFINEMENT_DISTORTIONS];
static kdu_int32 refinement_distortion_lut_lossless[REFINEMENT_DISTORTIONS];

#define EXTRA_ENCODE_CWORDS 3 // Number of extra context-words between stripes.
#define MAX_POSSIBLE_PASSES (31*3-2)

/* ========================================================================= */
/*                   Local Class and Structure Definitions                   */
/* ========================================================================= */

/*****************************************************************************/
/*                             kd_block_encoder                              */
/*****************************************************************************/

class kd_block_encoder : public kdu_block_encoder_base {
  /* Although we can supply a constructor and a virtual destructor in the
     future, we have no need for these for the moment. */
  protected:
    void encode(kdu_block *block, bool reversible, double msb_wmse,
                kdu_uint16 estimated_slope_threshold);
  private: // Internal implementation
    void reset_states()
      { // See Table 12.1 in the book by Taubman and Marcellin
        for (int n=0; n < 18; n++)
          states[n].init(0,0);
        states[KAPPA_SIG_BASE].init(4,0);
        states[KAPPA_RUN_BASE].init(3,0);
      }
  private: // Data
    mqe_state states[18];
  };

/* ========================================================================= */
/*            Initialization of Distortion Estimation Tables                 */
/* ========================================================================= */

static void initialize_significance_distortion_luts();
static void initialize_refinement_distortion_luts();
static class encoder_local_init {
    public: encoder_local_init()
              { initialize_significance_distortion_luts();
                initialize_refinement_distortion_luts(); }
  } _do_it;

/*****************************************************************************/
/* STATIC          initialize_significance_distortion_luts                   */
/*****************************************************************************/

static void
  initialize_significance_distortion_luts()
{
  double fp_scale = (double)(1<<16);
  for (kdu_int32 n=0; n < SIGNIFICANCE_DISTORTIONS; n++)
    {
      kdu_int32 idx = n | (1<<DISTORTION_LSBS);
      double v_tilde = ((double) idx) / ((double)(1<<DISTORTION_LSBS));
      assert((v_tilde >= 1.0) && (v_tilde < 2.0));
      double sqe_before = v_tilde*v_tilde;
      double sqe_after = (v_tilde-1.5)*(v_tilde-1.5);
      significance_distortion_lut[n] = (int)
        floor(0.5 + fp_scale*(sqe_before-sqe_after));
      significance_distortion_lut_lossless[n] = (int)
        floor(0.5 + fp_scale*sqe_before);
    }
}

/*****************************************************************************/
/* STATIC            initialize_refinement_distortion_luts                   */
/*****************************************************************************/

static void
  initialize_refinement_distortion_luts()
{
  double fp_scale = (double)(1<<16);
  for (kdu_int32 n=0; n < REFINEMENT_DISTORTIONS; n++)
    {
      double v_tilde = ((double) n) / ((double)(1<<DISTORTION_LSBS));
      assert(v_tilde < 2.0);
      double sqe_before = (v_tilde-1.0)*(v_tilde-1.0);
      v_tilde = (n >> DISTORTION_LSBS)?(v_tilde-1.0):v_tilde;
      assert((v_tilde >= 0.0) && (v_tilde < 1.0));
      double sqe_after = (v_tilde-0.5)*(v_tilde-0.5);
      refinement_distortion_lut[n] = (int)
        floor(0.5 + fp_scale*(sqe_before-sqe_after));
      refinement_distortion_lut_lossless[n] = (int)
        floor(0.5 + fp_scale*sqe_before);
    }
}


/* ========================================================================= */
/*             Binding of MQ and Raw Symbol Coding Services                  */
/* ========================================================================= */

#define USE_FAST_MACROS // Comment this out if you want functions instead.

#ifdef USE_FAST_MACROS
#  define _mq_check_out_(coder)                                     \
     register kdu_int32 A; register kdu_int32 C; register kdu_int32 t; \
     kdu_int32 temp; kdu_byte *store;                                 \
     coder.check_out(A,C,t,temp,store)
#  define _mq_check_in_(coder)                                      \
     coder.check_in(A,C,t,temp,store)
#  define _mq_enc_(coder,symbol,state)                              \
     _mq_encode_(symbol,state,A,C,t,temp,store)
#  define _mq_enc_run_(coder,run)                                   \
     _mq_encode_run_(run,A,C,t,temp,store)

#  define _raw_check_out_(coder)                                    \
     register kdu_int32 t; register kdu_int32 temp; kdu_byte *store;   \
     coder.check_out(t,temp,store)
#  define _raw_check_in_(coder)                                     \
     coder.check_in(t,temp,store)
#  define _raw_enc_(coder,symbol)                                   \
     _raw_encode_(symbol,t,temp,store)
#else // Do not use fast macros
#  define _mq_check_out_(coder)
#  define _mq_check_in_(coder)
#  define _mq_enc_(coder,symbol,state) coder.mq_encode(symbol,state)
#  define _mq_enc_run_(coder,run) coder.mq_encode_run(run)

#  define _raw_check_out_(coder)
#  define _raw_check_in_(coder)
#  define _raw_enc_(coder,symbol) coder.raw_encode(symbol)
#endif // USE_FAST_MACROS

 /* The coding pass functions defined below all return a 32-bit integer,
    which represents the normalized reduction in MSE associated with the
    coded symbols.  Specifically, the MSE whose reduction is returned is
    equal to 2^16 * sum_i (x_i/2^p - x_i_hat/2^p)^2 where x_i denotes the
    integer sample values in the `samples' array and x_i_hat denotes the
    quantized representation available from the current coding pass and all
    previous coding passes, assuming a mid-point reconstruction rule.
       The mid-point reconstruction rule satisfies x_i_hat = (q_i+1/2)*Delta
    where q_i denotes the quantization indices and Delta is the quantization
    step size.  This rule is modified only if the `lossless_pass' argument is
    true, which is permitted only when symbols coded in the coding pass
    result in a lossless representation of the corresponding subband samples.
    Of course, this can only happen in the last bit-plane when the reversible
    compression path is being used.  In this case, the function uses the fact
    that all coded symbols have 0 distortion.
       It should be noted that the MSE reduction can be negative, meaning
    that the coding of symbols actually increases distortion. */


/* ========================================================================= */
/*                           Coding pass functions                           */
/* ========================================================================= */

/*****************************************************************************/
/* STATIC                   encode_sig_prop_pass_raw                         */
/*****************************************************************************/

static kdu_int32
  encode_sig_prop_pass_raw(mq_encoder &coder, int p, bool causal,
                           kdu_int32 *samples, kdu_int32 *contexts,
                           int width, int num_stripes, int context_row_gap,
                           bool lossless_pass)
{
  /* Ideally, register storage is available for 9 32-bit integers. Two
     are declared inside the "_raw_check_out_" macro.  The order of priority
     for these registers corresponds roughly to the order in which their
     declarations appear below.  Unfortunately, none of these register
     requests are likely to be honored by the register-starved X86 family
     of processors, but the register declarations may prove useful to
     compilers for other architectures or for hand optimizations of
     assembly code. */
  register kdu_int32 *cp = contexts;
  register int c;
  register kdu_int32 cword;
  _raw_check_out_(coder); // Declares t and temp as registers.
  register kdu_int32 sym;
  register kdu_int32 val;
  register kdu_int32 *sp = samples;
  register kdu_int32 shift = 31-p; assert(shift > 0);
  int r, width_by2=width+width, width_by3=width_by2+width;
  kdu_int32 distortion_change = 0;
  kdu_int32 *distortion_lut = significance_distortion_lut;
  if (lossless_pass)
    distortion_lut = significance_distortion_lut_lossless;

  assert((context_row_gap - width) == EXTRA_ENCODE_CWORDS);
  for (r=num_stripes; r > 0; r--, cp += EXTRA_ENCODE_CWORDS, sp += width_by3)
    for (c=width; c > 0; c--, sp++, cp++)
      {
        if (*cp == 0)
          continue;
        cword = *cp;
        if ((cword & (NBRHD_MASK<<0)) && !(cword & (SIG_PROP_MEMBER_MASK<<0)))
          { // Process first row of stripe column (row 0)
            val = sp[0]<<shift; // Move bit p to sign bit.
            sym = (kdu_int32)(((kdu_uint32) val)>>31); // Move bit into LSB
            _raw_enc_(coder,sym);
            if (val >= 0) // New magnitude bit was 0, so still insignificant
              { cword |= (PI_BIT<<0); goto row_1; }
            // Compute distortion change
            val =  (val>>(31-DISTORTION_LSBS)) & (SIGNIFICANCE_DISTORTIONS-1);
            distortion_change += distortion_lut[val];
            // Encode sign bit
            sym = sp[0];
            sym = (kdu_int32)(((kdu_uint32) sym)>>31); // Move sign into LSB
            _raw_enc_(coder,sym);
            // Broadcast neighbourhood context changes
            if (!causal)
              {
                cp[-context_row_gap-1] |=(SIGMA_BR_BIT<<9);
                cp[-context_row_gap  ] |=(SIGMA_BC_BIT<<9)|(sym<<NEXT_CHI_POS);
                cp[-context_row_gap+1] |=(SIGMA_BL_BIT<<9);
              }
            cp[-1] |= (SIGMA_CR_BIT<<0);
            cp[1]  |= (SIGMA_CL_BIT<<0);
            cword |= (SIGMA_CC_BIT<<0) | (PI_BIT<<0) | (sym<<CHI_POS);
          }
row_1:
        if ((cword & (NBRHD_MASK<<3)) && !(cword & (SIG_PROP_MEMBER_MASK<<3)))
          { // Process second row of stripe column (row 1)
            val = sp[width]<<shift; // Move bit p to sign bit.
            sym = (kdu_int32)(((kdu_uint32) val)>>31); // Move bit into LSB
            _raw_enc_(coder,sym);
            if (val >= 0) // New magnitude bit was 0, so still insignificant
              { cword |= (PI_BIT<<3); goto row_2; }
            // Compute distortion change
            val =  (val>>(31-DISTORTION_LSBS)) & (SIGNIFICANCE_DISTORTIONS-1);
            distortion_change += distortion_lut[val];
            // Encode sign bit
            sym = sp[width];
            sym = (kdu_int32)(((kdu_uint32) sym)>>31); // Move sign into LSB
            _raw_enc_(coder,sym);
            // Broadcast neighbourhood context changes
            cp[-1] |= (SIGMA_CR_BIT<<3);
            cp[1]  |= (SIGMA_CL_BIT<<3);
            cword |= (SIGMA_CC_BIT<<3) | (PI_BIT<<3) | (sym<<(CHI_POS+3));
          }
row_2:
        if ((cword & (NBRHD_MASK<<6)) && !(cword & (SIG_PROP_MEMBER_MASK<<6)))
          { // Process third row of stripe column (row 2)
            val = sp[width_by2]<<shift; // Move bit p to sign bit.
            sym = (kdu_int32)(((kdu_uint32) val)>>31); // Move bit into LSB
            _raw_enc_(coder,sym);
            if (val >= 0) // New magnitude bit was 0, so still insignificant
              { cword |= (PI_BIT<<6); goto row_3; }
            // Compute distortion change
            val =  (val>>(31-DISTORTION_LSBS)) & (SIGNIFICANCE_DISTORTIONS-1);
            distortion_change += distortion_lut[val];
            // Encode sign bit
            sym = sp[width_by2];